Untitled

"""A partial list of built-in functions in SAGE we shall use: 1. gcd(a,b) returns the greatest common divisor of the integers a and b. 2. xgcd(a,b) returns three integers (c, n, m) where c = n*a + m*b 3. a^b returns a to the power b. 4. a % b or mod(a, b) returns a mod b. 5. a ^ (-1) % b returns n such that a*n mod b = 1. It will return an error if gcd(a,b) > 1. 6. crt(a1,a2,m1,m2) returns an integer x such that x mod ai = bi Note: this is the Chinese Remainder Theorem. Even though it is valid for any number of arguments SAGE only accepts 2 terms. To calculate any number of terms write your own function. 7. len(str(n)) returns the number of digits of the integer n. 8. n // m = floor(n/m). n/m is a fraction, not Integer. 9. factor(n) returns the factors of n (if it can execute it). 10. factor(n, algorithm = 'kash') Uses a sophisticated algorithm to factor n. You need to load the package kash if you wish to use it. 11. is_prime(n) returns True or False 12. primes(n,m) returns a list of all primes between n and m. 13. pow(a,b,c) returns a^b mod c. 14. fctorial(n) returns n! 15. next_prime(n) returns the smallest prime p >= n. This list is not complete. There are other functions and we can also implement other functions. Examples follow: """ 
       
#Getting help #topic? primes? 
       

File: /home/sage/sage/local/lib/python2.6/site-packages/sage/rings/arith.py

Type: <type ‘function’>

Definition: primes(start, stop=None, proof=None)

Docstring:

Returns an iterator over all primes between start and stop-1, inclusive. This is much slower than prime_range, but potentially uses less memory. As with next_prime(), the optional argument proof controls whether the numbers returned are guaranteed to be prime or not.

This command is like the xrange command, except it only iterates over primes. In some cases it is better to use primes thanprime_range, because primes does not build a list of all primes in the range in memory all at once. However, it is potentially much slower since it simply calls the next_prime() function repeatedly, and next_prime() is slow.

INPUT:

  • start - an integer - lower bound for the primes
  • stop - an integer (or infinity) optional argument - giving upper (open) bound for the primes
  • proof - bool or None (default: None) If True, the function yields only proven primes. If False, the function uses a pseudo-primality test, which is much faster for really big numbers but does not provide a proof of primality. If None, uses the global default (see sage.structure.proof.proof)

OUTPUT:

  • an iterator over primes from start to stop-1, inclusive

EXAMPLES:

sage: for p in primes(5,10):
...     print p
...
5
7
sage: list(primes(13))
[2, 3, 5, 7, 11]
sage: list(primes(10000000000, 10000000100))
[10000000019, 10000000033, 10000000061, 10000000069, 10000000097]
sage: max(primes(10^100, 10^100+10^4, proof=False))
10000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000009631
sage: next(p for p in primes(10^20, infinity) if is_prime(2*p+1))
100000000000000001243

TESTS:

sage: for a in range(-10, 50):
...    for b in range(-10, 50):
...        assert list(primes(a,b)) == list(filter(is_prime, xrange(a,b)))
...
sage: sum(primes(-10, 9973, proof=False)) == sum(filter(is_prime, range(-10, 9973)))
True
sage: for p in primes(10, infinity):
...    if p > 20: break
...    print p
...
11
13
17
19
sage: next(p for p in primes(10,oo)) # checks alternate infinity notation
11
gcd(3456, 3906) 
       
18
xgcd(3456, 3906) 
       
(18, 26, -23)
#Check: 26*3456 -23 *3906 
       
18
105^2, 11^4 
       
(11025, 14641)
pow(105, 7, 3906) 
       
3465
 
       
107^(-1) % 3906 
       
3833
3833*107 % 3906 
       
1
crt(34, 71, 119, 235) 
       
12291
#testing: 12291 % 119, 12291 % 235 
       
(34, 71)
len(str(2^256 + 1)) 
       
78
96570387/ 77, 96570387 // 77 
       
(96570387/77, 1254160)
factor(factorial(23) + 1) 
       
47^2 * 79 * 148139754736864591
n = 2^256 + 1 n 
       
11579208923731619542357098500868790785326998466564056403945758400791\
3129639937
is_prime(n) 
       
False
#factor(n) will take time. 
       
for i in primes(10001, 12000): print i, 
       
10007 10009 10037 10039 10061 10067 10069 10079 10091 10093 10099
10103 10111 10133 10139 10141 10151 10159 10163 10169 10177 10181
10193 10211 10223 10243 10247 10253 10259 10267 10271 10273 10289
10301 10303 10313 10321 10331 10333 10337 10343 10357 10369 10391
10399 10427 10429 10433 10453 10457 10459 10463 10477 10487 10499
10501 10513 10529 10531 10559 10567 10589 10597 10601 10607 10613
10627 10631 10639 10651 10657 10663 10667 10687 10691 10709 10711
10723 10729 10733 10739 10753 10771 10781 10789 10799 10831 10837
10847 10853 10859 10861 10867 10883 10889 10891 10903 10909 10937
10939 10949 10957 10973 10979 10987 10993 11003 11027 11047 11057
11059 11069 11071 11083 11087 11093 11113 11117 11119 11131 11149
11159 11161 11171 11173 11177 11197 11213 11239 11243 11251 11257
11261 11273 11279 11287 11299 11311 11317 11321 11329 11351 11353
11369 11383 11393 11399 11411 11423 11437 11443 11447 11467 11471
11483 11489 11491 11497 11503 11519 11527 11549 11551 11579 11587
11593 11597 11617 11621 11633 11657 11677 11681 11689 11699 11701
11717 11719 11731 11743 11777 11779 11783 11789 11801 11807 11813
11821 11827 11831 11833 11839 11863 11867 11887 11897 11903 11909
11923 11927 11933 11939 11941 11953 11959 11969 11971 11981 11987
pow(2, n-1, n) 
       
1
pow(3, n-1, n) 
       
11308059312705222464474529196106459540324134768955225107825802801824\
6279223993
""" Finding sqrt(a) mod p when p is a prime number. 1. Some facts about the finite field GF(p). a. GF(p) has primitive elements. a is primitive if {1,a,a^2, a^3, ..., a^{p-2}} = {1,2,...,p-1} Example: 1, 2, 2^2 = 4, 2^3 = 3 mod 5, 2 is a primitive element mod 5. b. q is a quadratic residue mod p if q = x^2 mod p. c. (p-1)/2 numbers in GF(p) are quadratic residues. d. q is a quadratic residue mod p if and only if q^{(p-1)/2} = 1 mod p. 2. Given q, a quadratic residue mod p, we wish to find its square root, that is find x such that x^2 mod p = q. 3. Here are a couple of calculations with explanations. """ 
       
""" Does 111111 have a square root mod z? Be careful. Even though (z-1)/2 is an integer, SAGE casts it as rational, so you cannot execute powmod(111111, (z-1)/2, z). Remedy cast it as Integer. There are two ways to do it: Integer((z-1)/2) or (z-1)//2 """ z = next_prime(3792643488006829893221399440992214604544353) pow(111111, Integer((z-1)/2),z) 
       
3792643488006829893221399440992214604544370
#We are not lucky. 111111 is not a quadratic residue mod z z 
       
3792643488006829893221399440992214604544371
#Let us try 1111111 pow(1111111, (z-1)//2, z) 
       
1
#We are lucky. m = (z-1)//2 m 
       
1896321744003414946610699720496107302272185
#We have: 1111111^m = 1 or 1111111^(m+1) = 1111111. b = pow(1111111, (m+1)//2, z) b 
       
889727700437799721649009803358653820921488
#Testing: b^2 % z 
       
1111111
#We were lucky as z = 3 mod 4 and (z-1)//2 was odd. #Let us see how thos works when p = 1 mod 4. p = next_prime(z) p 
       
3792643488006829893221399440992214604544473
#Good. p = 1 mod 4. pow(1111, (p-1)//2, p) 
       
1
#1111 is a quadratic residue mod p. n1 = (p-1)//2 n1 
       
1896321744003414946610699720496107302272236
#We know that 1111^n1 mod p = 1. So 1111^(n1/2) is either 1 or -1. a1 = pow(1111, n1//2, p) a1 
       
1
n2 = n1//2 n2 
       
948160872001707473305349860248053651136118
pow(1111, n2//2, p) 
       
1
#Note: n2//2 is an odd integer. Also 1111^(n2//2 + 1) = 1111. So: rt = pow(1111, (n2//2 +1)//2, p) rt 
       
3728998878451604586953703188757611624985334
rt^2 % p 
       
1111
#Quite amazing! a small number like 1111 has such a huge square root. You probably wonder what is the # other square root: st = (p - rt) st 
       
63644609555225306267696252234602979559139
st^2 % p 
       
1111
"""Finding square roots of n mod pq. An example. Note: we assume that we know the primes p and q. """ 
       
#WE start by choosing two very large integers: p = next_prime(73 ^ 54) p 
       
41632687486237173187939287634191588382353880270019124010730946361508\
962534646068310417373812845345157
q = next_prime(5*p) q 
       
20816343743118586593969643817095794191176940135009562005365473180754\
4812673230341552086869064226725861
len(str(p*q)) 
       
202
"""We have a 202 decimal digit long integer. Do not bother to ask SAGE to factor it unless you plan to live to about 250 years of age. """ p*q 
       
86664033366334465731156283197630132969913941415627514510407867309397\
67297782967082791305269130056723358212019627931603556954213403413843\
385561078736730894541271046896698119880843244321621736788563005177
#We do not know whether n has a square root mod pq. But it's square definitely does. key = p*q n 
       
1111111111
m = n^2 key = p*q m 
       
1234567900987654321
 
       
x1 = (n*ap*p - n*aq*q)%key x1 
       
84383400909325664001389012587166182102284627167847843075923449748624\
05000472889001665218288357389633053160547702150425347063947609327701\
098838164971184023930472120531215541863956614064356309828730138742
x1^2 % key 
       
1234567900987654321
#so 1111111111 and x1 are two different square roots of (1111111111)^2 
       
len(str(x1)) 
       
202
#And here are the other roots: x3 = (key - x1); x4 =(key - 1111111111) x3, x4 
       
(2280632457008801729767270610463950867629314247779671434484417560773\
62297310078081126086980772667090305051471925781178209890265794086142\
286722913765546870610798926365482578016886630257265426959832866435,
86664033366334465731156283197630132969913941415627514510407867309397\
67297782967082791305269130056723358212019627931603556954213403413843\
385561078736730894541271046896698119880843244321621736787451894066)
#Check x3^2 % key 
       
1234567900987654321
#Suppose you "discovered" the two roots 1111111111 and x3. We can try to use them to factor the key. y = gcd(x3 - 1111111111, key) y 
       
41632687486237173187939287634191588382353880270019124010730946361508\
962534646068310417373812845345157
#Looks like y = p. Let us check: y - p 
       
0
#Indeed we managed to factor the key. Finding the second factor is now easy: z = key/y z 
       
20816343743118586593969643817095794191176940135009562005365473180754\
4812673230341552086869064226725861
z - q 
       
0