Let $p,q$ be two large primes of comparable size. Put $n = pq$. This computation takes about $O(\log_2(n)^2)$ many flops if you use elementary school method of multiplying these numbers (it is actually $O(\log_2(n)\log\log_2(n))$ many flops if you employ Fast Fourier Transform multiplication algorithm).
However trying to factor $n$ into $p$ and $q$ using the brute force method will take $O(\max\{p,q\})$ modular divisions. which is on the order of $O(\sqrt{n})$ if $p$ and $q$ are close enough to one another.
163298476128374618237491823746912837468127467125648723654817356418765*298734691328746192837461923874619287461983278
48782919860664620234782399461351856881494516091785396502382566480860785545607256493096155801868780530190395411670
def BruteForceFactor(n):
primelist = []
d = 2
while n!= 1:
if d*d>n:
primelist.append(n)
return primelist
elif n%d == 0:
primelist.append(d)
n = n//d
else:
d = d + 1
return primelist
BruteForceFactor(156)
[2, 2, 3, 13]
BruteForceFactor(123456789)
[3, 3, 3607, 3803]
BruteForceFactor(1234567891)
[1234567891]
BruteForceFactor(1234567891011121314151617181920)
[2, 2, 2, 2, 2, 3, 5, 323339, 3347983, 2375923237887317]
We are interested in the order of the group $(\mathbb{Z}/n\mathbb{Z})^*$. That is, what is the number $$\phi(n) := \#\{1\leq a \leq n: \operatorname{gcd}(a,n) = 1\}.$$
If $n = p$ is a prime, we know that every integer less than $p$ is relatively prime to $p$, and hence $\phi(p) = p-1$.
If $n = p^k$ a power of a prime, then the only integers that are not relatively prime to $n$ are multiples of $p$, and there are $n/p = p^{k-1}$ many of them. Hence $\phi(p^k) = p^k - p^{k-1} = p^k(1 - 1/p)$.
Lastly by Chinese Remainder theorem we note that if $n = ab$ with $\operatorname{gcd}(a,b) = 1$, then there is an isomorphism of groups $$ (\mathbb{Z}/n\mathbb{Z})^* \cong (\mathbb{Z}/a\mathbb{Z})^* \times (\mathbb{Z}/b\mathbb{Z})^*. $$ There are $\phi(n)$ elements on the left, and $\phi(a)\phi(b)$ elements on the right.
Because of this we say that $\phi$ is a multiplicative function and for a general $n$ of the form, $$ n = p_1^{k_1} p_2^{k_2} \cdots p_r^{k_r} $$ we deduce that $$ \phi(n) = (p_1^{k_1} - p_1^{k_1-1}) \cdots (p_r^{k_r} - p_r^{k_r-1}) = n \prod_{p | n} \left(1 - \frac 1p\right) $$
def EulerPhi(n, factors = None):
if factors == None:
primelist = BruteForceFactor(n)
else:
primelist= factors
result = 1
while len(primelist)>1:
ult = primelist[-1]
penult = primelist[-2]
if ult == penult:
result = result*ult
primelist.pop()
else:
result = result*(ult -1)
primelist.pop()
return result*(primelist[-1]-1)
EulerPhi(101)
100
EulerPhi(25)
20
EulerPhi(1234567891)
1234567890
EulerPhi(6)
2
def fastexp(x,n,m):
z = 1
while n!=0:
if n%2==1:
z = (z*x)%m
x = (x*x)%m
n = n//2
return z
Let $p = 1234567891$ be a prime, and choose a number $m$ as the message to be sent.
Computing $m^{123876} \pmod{p}$ gives:
p = 1234567891
m = 58008 #a message for the other side
c = fastexp(m,1234567, p)
c
1010044275
import matplotlib.pyplot as plt
ks = range(2,640)
y = [fastexp(m,k,p) for k in ks]
s = [5 for i in ks]
plt.scatter(ks,y,s)
plt.show()
<matplotlib.figure.Figure at 0x1f9a79f7710>
Actually it is not very hard to solve $m$ given $c$ using the fact that $p$ is prime, and we have $$|(\mathbb{Z}/p\mathbb{Z})^*| = p-1$$ for a prime number $p$. Hence we have that $$ a^{p-1} \equiv 1 \pmod{p}.$$ This is because the order of any element $a$ divides the order of the group (that is $p-1$). So $a^p \equiv a \pmod{p}$ for all integers $a$, or in fact $a^{e} \equiv a \pmod p$ for any $e \equiv 1 \pmod {p-1}$.
#Let us bring the things we have written again.
def bezout(n,m):
(g,u,v) = (n,1,0)
(e,s,t) = (m,0,1)
r = g%e
while r != 0:
q = g//e
(eprime, sprime,tprime) = (e,s,t)
(e,s,t) = (r,u - q*s,v - t*q)
(g,u,v) = (eprime,sprime,tprime)
r = g%e
return (e,s,t)
def gcd(n,m):
(g,u,v) = bezout(n,m)
return g
def inverse(x,n):
(g,u,v) = bezout(x,n)
if g!=1:
return None
else:
return (u%n)
e = 1234567
d = inverse(e, p-1)
d
18033013
(e*d)%1234567890
1
n = 10403
fastexp(fastexp(fastexp(fastexp(fastexp(fastexp(7580,5,10403),6,n),7,n),8,n),9,n),10,n)
9798
fastexp(c,d,p)
58008
gcd(9797,10403)
101
for k in range(15):
print("2^", k," =", fastexp(2,k,7), 'mod 7')
2^ 0 = 1 mod 7 2^ 1 = 2 mod 7 2^ 2 = 4 mod 7 2^ 3 = 1 mod 7 2^ 4 = 2 mod 7 2^ 5 = 4 mod 7 2^ 6 = 1 mod 7 2^ 7 = 2 mod 7 2^ 8 = 4 mod 7 2^ 9 = 1 mod 7 2^ 10 = 2 mod 7 2^ 11 = 4 mod 7 2^ 12 = 1 mod 7 2^ 13 = 2 mod 7 2^ 14 = 4 mod 7
We got $m$ back!!! Why was that, well here is the reason $$ c^{d} \equiv (m^e)^d \equiv m^{ed} \pmod p$$ and $ed \equiv 1 \pmod {p-1}$ and therefore $$c^{d} \equiv m^{ed} \equiv m \pmod {p}.$$
Then we can do calculation in $G = (\mathbb{Z}/n\mathbb{Z})^*$ however we don't necessarily know $|G| = \phi(n)$. If we do, then the above method also solves the equation for $m$: $$ m^e = \equiv c \pmod{n}.$$
Knowing $\phi(n)$ is easy if you know the factorization of $n$, which means that we are looking for numbers that are hard to factor. The smallest such numbers are of the form $$ n = pq$$ for primes $p$ and $q$.
q = 12345678910111
BruteForceFactor(q)
[12345678910111]
n = p*q
n
15241578775018915845901
Dont even try BruteForceFactor(n), it takes too long!!!
from random import randint as rdint
def RSAkeys(n,factors):
e = rdint(n//4,n//2)
while not all([gcd(e, p-1)==1 for p in factors]):
e = rdint(n//4,n//2)
order = EulerPhi(n,factors)
d = inverse(e,order)
return (e,d)
def RSAencrypt(m, n, e):
c = fastexp(m,e,n)
return c
def RSAdecrypt(c,n, d):
m = fastexp(c,d,n)
return m
(e,d) = RSAkeys(n,[p,q])
(e,d)
(3937139939893578395651, 15184922774872748407451)
(p-1)*(q-1)
15241578762672002367900
(e*d)%((p-1)*(q-1))
1
c = RSAencrypt(58008, n, e)
c
7399584266783735548498
RSAdecrypt(c,n,d)
58008
Alice chooses two primes $p,q$ and computes $n = pq$. Choses an $e$ such that $\operatorname{gcd}(e,(p-1)(q-1))$. Computes $d$ such that $ed \equiv 1 \pmod{\phi(n)}$.
Alice publishes $e, n$ and keeps $d$ to herself.
Bob has a message $m$, and computes $c \equiv m^e \pmod{n}$ at the privacy of his own office and sends Alice the cipher $c$.
Alice (since only she has knowledge of the secret key $e$) computes $c^e \pmod{n}$ which is equivalent to $m$.
Note that if $n = pq$ then $\phi(n) = (p-1)(q-1) = pq - (p + q) + 1= n - (p + q) + 1$. Thus if we know $n$ and $\phi(n)$ then we can easily calculate $p + q = n- \phi(n) +1$, and knowing both $pq$ and $p+q$ we solve the quadratic equation $$X^2 + (p + q) X + pq = 0$$ in order to get $p$ and $q$.
We are given $c \equiv m^{e} \pmod{n}$ and for $\operatorname{gcd}(e, (p-1)(q-1))=1$. So we are trying to find the $e$-th root of $c$ modulo $n$. We said that with $d$ given as $ed \equiv 1 \pmod \phi(n)$ then $x \equiv c^d\pmod{n}$ is a solution to $x^e \equiv c \pmod{n}$. It is also the only solution:
Assume $x$ were a solution, and that $de = 1 + k\phi(n) $ for some integer $k$, $$ \begin{align*} x &\equiv x^{de - k\phi(n) }\pmod{n} \\ &\equiv (x^e)^d (x^{\phi(n)})^{-k} \pmod{n}\\ &\equiv c^d 1^{-k} \equiv c^d \pmod{n}\end{align*} $$
(Note that if the exponent $e$ is not relatively prime to the order of the group, it is not necessarily true that $x^e \equiv c \pmod{n}$ is solvable for generic $c$. For example not every integer modulo $p$ is a square.)
It is generally assumed (especially by public) that the difficulty of breaking RSA is equivalent to factorization. Technically what one needs is to solve the equation $x^e \equiv c \pmod{n}$, which may or may not be strictly easier than the problem of factorizing $n$. The expert opinion is a bit more complicated:
Divesh Aggarwal, Ueli Maurer, <i> <b>Breaking RSA Generically Is Equivalent to Factoring </b></i> Joux A. (eds) Advances in Cryptology - EUROCRYPT 2009. EUROCRYPT 2009. Lecture Notes in Computer Science, vol 5479. Springer, Berlin, Heidelberg
If $p$ is a prime number then $(\mathbb{Z}/p\mathbb{Z})^*$ has $p-1$ integers and therefore the order of any nonzero element $a \pmod{p}$ divides $p-1$. In particular, $$ a^{p-1} \equiv 1 \pmod{p}.$$ Equivalently $$ a^{p}\equiv a \pmod{p}$$ for all integers $a$.
n = 329487987
fastexp(2,n, n)
211367123
BruteForceFactor(n)
[3, 19, 5780491]
Since $2^n \not\equiv 2 \pmod{n}$ we know that $n$ is not a prime.
n = 12801
fastexp(2,n,n)
2
BruteForceFactor(n)
[3, 17, 251]
fastexp(3,n,n)
9693
As can be seen even though $n = 12801$ is not a prime, $2^n \equiv 2 \pmod{n}$. Yet the base $3$ comes to save the day, and is a witness that $n$ is composite.
def CompositeWitness(a,n):
return fastexp(a,n,n)!=(a%n)
CompositeWitness(2,n)
False
CompositeWitness(3,n)
True
Let us see of all the composite numbers, how many of them are caught by various Fermat witnesses.
def CompositeList(N):
return list(filter(lambda n: len(BruteForceFactor(n))>1, range(2,N)))
comp = CompositeList(10000)
len(comp)
comp[:20]
[4, 6, 8, 9, 10, 12, 14, 15, 16, 18, 20, 21, 22, 24, 25, 26, 27, 28, 30, 32]
list(filter(lambda n: not CompositeWitness(2,n),comp))
[341, 561, 645, 1105, 1387, 1729, 1905, 2047, 2465, 2701, 2821, 3277, 4033, 4369, 4371, 4681, 5461, 6601, 7957, 8321, 8481, 8911]
list(filter(lambda n: not (CompositeWitness(2,n) or CompositeWitness(3,n)),comp))
[561, 1105, 1729, 2465, 2701, 2821, 6601, 8911]
def CarmichaelNumbers(N):
#Returns all the composite numbers n less than N such that
#All a < n are non-witnesses for the compositeness of n
comp = CompositeList(N)
a = 2
while (a < comp[-1]) and len(comp)>0:
comp = list(filter(lambda n: not CompositeWitness(a,n), comp))
a = a + 1
return list(comp)
CarmichaelNumbers(10000)
[561, 1105, 1729, 2465, 2821, 6601, 8911]
Even though $561 = 3\times 11\times 17$, it is true that $a^n \equiv a \pmod{n}$ for all integers $a$. So we have a sufficient method for testing primality.
Theorem. All Charmichael numbers are of the form $n = p_1 p_2 \cdots p_k$ with $k>3$ such that the prime divisors $p_i$ are all distinct and $(p_i -1) | (n-1)$ for all $i$.
BruteForceFactor(561)
[3, 11, 17]
X = range(10**8,10**8 + 100000)
X = filter(lambda n: not(CompositeWitness(2,n) or CompositeWitness(3,n)),X)
X = list(X)
print(len(X))
5412
for p in X:
if len(BruteForceFactor(p))>1:
print(str(p) + ' is not a prime, but pretends to be')
100017223 is not a prime, but pretends to be
Let us test whether $n$ is a prime using the following algorithm. First factorize as $$n -1 = 2^k q$$ with $q$ an odd integer. Let $L_n$ be the set of nonzero integers $a$ less than $n$ such that either $$a^{2^kq} \equiv 1 \pmod{n}$$ or in the set $$a^q, a^{2q}, a^{2^2q}, \ldots, a^{2^{k-1}q}$$ of successive squares they are all congruent to $1$ or one element is congruent to $-1$ modulo $n$.
If $n$ is prime then $L_n = (\mathbb{Z}/N\mathbb{Z})^*$.
If $n$ is not prime then $|L_n| \leq (n-1)/4$. This means that for $n$ composite $\%75$ of the bases $a$ are Witnesses for the compositeness of $n$ using the Miller Rabin primality test. This means that a probabilistic argument can be made.
Partial Proof: If $n=p$ is a prime then since $a^{n-1} \equiv 1 \pmod {n}$ then eitber all integers in the sequence are equivalent to $1$ or in this list of succesive squares one integer $b \not \equiv 1 \pmod{p}$ and $b^2 \equiv 1 \pmod{p}$. This implies that $p | b^2 -1$ hence $p | (b-1)$ or $p | b + 1$. The first case would mean $b \equiv 1 \pmod{p}$, and hence we have $b \equiv -1 \pmod{p}$.
def powersOfTwo(n):
"""
Returns (k,q) such that n = 2^k q
"""
q = n
k = 0
while q%2==0:
k = k + 1
q = q//2
return (k,q)
def MillerRabinWitness(a,n):
"""
Tests whether a is a witness for the compositeness of n
"""
g = gcd(n,a)
if g > 1 and g < n:
return True
(k,q) = powersOfTwo(n-1)
b = fastexp(a,q,n)
if b==1:
return False
for i in range(k):
if b==n-1:
return False
b = (b*b)%n
return True
from math import log as log
from random import randint as randint
def MillerRabinPrimalityTest(n, precision = None):
if precision == None:
precision = min(1/n,1/1000)
if n==2:
return True
if n%2==0:
return False
A = int(log(precision,3/4))+1
for trial in range(A):
a = randint(2,n-1)
if MillerRabinWitness(a,n):
return False
return True
MillerRabinPrimalityTest(105,0.001)
False
def FindPrimes(N,h):
#Returns the list of (Miller-Rabin)--primes in the interval [N, N + h)
X = range(N, N + h)
X = list(filter(MillerRabinPrimalityTest, X))
return X
FindPrimes(1000000000000000000000000000000,100)
[1000000000000000000000000000057, 1000000000000000000000000000099]
It can be shown that if one assumes the Riemann hypothesis, then if $n$ is not a prime, then an integer $a$ will be a witness for the compositeness of $a$ such that $$a \leq 2\ln(n).$$ This gives us a deterministic test for primality (dependent on Generalized Riemann Hypothesis).
def Pollardp(n, B):
a = 2
k = max(1,B//1000) #this is a optimization technique, we don't check the gcd everytime
for j in range(2,B):
a = fastexp(a,j,n)
if j%k ==0:
g = gcd(a-1,n)
if g >1 and g<n:
return (g,n//g)
return None
Pollardp(101*103,100)
(101, 103)
B = 10**7
def BSmooth(n,B):
#Checks whether n is B smooth or not
for prime in range(2,B):
while n%prime==0:
n=n//prime
if n ==1:
return True
else:
return False
def BSmoothFinder(alist, B):
for n in alist:
if BSmooth(n-1,B):
return n
return None
primes = FindPrimes(223456789101,1000)
q = primes.pop()
p = BSmoothFinder(primes,B)
print('p = ',p, 'q = ', q)
n = p*q
n
p = 223456789121 q = 223456790033
49932936808059655630993
Pollardp(n,B)
(223456789121, 223456790033)
BruteForceFactor(p-1)
[2, 2, 2, 2, 2, 2, 2, 5, 31, 71, 158633]
BruteForceFactor(q-1)
[2, 2, 2, 2, 7, 643, 3102877]
If we find integers such that $a^2 \equiv b^2 \pmod n$ with $n = pq$ then we can use that since $$pq | (a-b) (a + b) $$ we may have that $p | (a + b)$ and $b | (a - b) $.
How to find such numbers.
Let $p$ and $q$ be as above.
n = 914387
from math import sqrt as sqrt
a = int(sqrt(n)) + 1
sample =0
while sample <5:
x = (a*a)%n
if BSmooth(x,13):
print(a,"^2 =", BruteForceFactor(x), "mod n")
sample = sample +1
a = a + 1
1869 ^2 = [2, 2, 2, 2, 3, 5, 5, 5, 5, 5, 5] mod n 1909 ^2 = [2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 5, 11] mod n 3387 ^2 = [3, 5, 5, 5, 11, 11, 11] mod n 4586 ^2 = [3, 3, 5, 11] mod n 5023 ^2 = [2, 2, 2, 2, 2, 2, 2, 5, 7, 11, 11] mod n
a = (1869*1909*3387)%n
b = (2**9*3*5**5*11**2)%n
print(a,b)
9835 164255
gcd(a - b,n)
1103
gcd(a + b,n)
829
1103 * 829 ==n
True
n = 9873417
from math import sqrt as sqrt
a = int(sqrt(n)) + 1
sample =0
while sample <10:
x = (a*a)%n
if BSmooth(x,20):
print(a,"^2 =", BruteForceFactor(x), "mod n")
sample = sample +1
a = a + 1
3144 ^2 = [3, 7, 7, 7, 11] mod n 3177 ^2 = [2, 2, 2, 3, 7, 7, 11, 17] mod n 3487 ^2 = [2, 2, 2, 7, 7, 7, 7, 7, 17] mod n 3555 ^2 = [2, 2, 2, 2, 2, 2, 3, 7, 11, 11, 17] mod n 4578 ^2 = [2, 3, 5, 5, 5, 5, 17, 19] mod n 4806 ^2 = [2, 3, 7, 13, 17, 19, 19] mod n 4833 ^2 = [3, 5, 7, 7, 17, 17, 17] mod n 5499 ^2 = [2, 3, 3, 5, 5, 5, 5, 5, 11] mod n 5655 ^2 = [2, 3, 3, 3, 11, 11, 19, 19] mod n 6288 ^2 = [2, 2, 3, 7, 7, 7, 11] mod n
a = (3144*3177*3487)%n
b = (2**3*3**7**5*11*17)%n
print(a, b)
6315897 3814236
gcd(a + b, n)
3
BruteForceFactor(n)
[3, 23, 143093]
n = 81407
a = int(sqrt(n)) + 1
sample =0
while sample <10:
x = (a*a)%n
if BSmooth(x,12):
print(a,"^2 =", BruteForceFactor(x), "mod n")
sample = sample +1
a = a + 1
407 ^2 = [3, 3, 3, 3, 5, 7] mod n 412 ^2 = [2, 3, 3, 5, 7, 11] mod n 435 ^2 = [7, 7, 7, 7, 11] mod n 484 ^2 = [2, 3, 3, 3, 3, 3, 3, 7, 7] mod n 638 ^2 = [3, 3] mod n 755 ^2 = [2, 2, 2, 2, 11] mod n 757 ^2 = [2, 2, 2, 2, 2, 2, 2, 5, 5] mod n 763 ^2 = [2, 2, 2, 2, 2, 5, 7, 11] mod n 777 ^2 = [2, 2, 2, 5, 7, 11, 11] mod n 814 ^2 = [2, 2, 3, 3, 3, 3, 5, 7] mod n
(503**2)%n
8788
(313*484)%n
70085
(2*3**3*7**2*13)%n
34398
a = 78005
b = 34398
n
81407
a + b
112403
gcd(a - b,n)
1
a = (2*3**3*5*7**4*11*13)%n
b = (313*407*412*435)%n
print(a,b)
61444 14835
gcd(a + b,n)
641