Factor the natural number $ n $ into prime factors with $ O (n ^ {\ frac {1} {4}}) $.
When factoring a natural number $ n $ into prime factors, a naive implementation is to divide it sequentially up to the square root of $ n $. This method requires a computational complexity of $ O (\ sqrt {n}) $ and is inefficient for large $ n $.
First, consider performing primality tests on natural numbers at high speed. Here, we use what is called the Miller-Rabin method. First, let's review Fermat's Little Theorem.
According to Fermat's Little Theorem, the following holds for the prime number $ p $ and the coprime natural number $ a
Taking the kinematic pair, about a certain $ a $
If is not true, then $ p $ is not a prime number. Conversely, if we make sure that this is not the case for many $ a $, can we say that $ p $ is a prime number? In fact, unfortunately there is a counterexample of (?) Carmichael number. The Carmichael number $ n $ holds for all natural numbers $ a $ that are composite but $ (2) $ are relatively prime to $ n $. We know that there are an infinite number of Carmichael numbers.
Now, for a prime $ p $ over $ 3 $, in the world of $ \ mathbb {Z} / p \ mathbb {Z} $, there are exactly $ 2 $ square roots of $ 1 $ and $ -1 $. .. This can be seen from the fact that the quotient ring $ \ mathbb {Z} / p \ mathbb {Z} $ becomes the field. If $ n $ is a prime number greater than or equal to $ 3 $, the index $ n-1 $ on the left side of $ (2) $ is even. From $ a ^ {n-1} \ equiv 1 $, $ a ^ {\ frac {n-1} {2}} \ equiv \ pm1 $ should hold. Conversely, if $ a ^ {\ frac {n-1} {2}} \ equiv \ pm1 $ does not hold, then $ n $ is a composite number. If you halve the index of $ a $, it should be $ \ pm1 $ just before $ 1 $. If a number other than $ \ pm1 $ suddenly becomes $ 1 $, you know that $ n $ is a composite number. However, if it becomes $ -1 $ when viewed from the back, it cannot be determined whether $ n $ is a prime number or a composite number by itself, no matter what the front. The probability that such a number, whether prime or not, will be chosen is at most $ \ frac {1} {4} $. In other words, with a lot of trials, the probability of making a correct decision approaches $ 1 $.
Note that if $ n <2 ^ {32} $, it is sufficient to look at $ 2, \ 3, \ 61 $ as $ a $. If $ n <2 ^ {64} $, look for $ 2, \ 3, \ 5, \ 7, \ 11, \ 13, \ 17, \ 19, \ 23, \ 29, \ 31, \ 37 $ is enough. Details of the sufficiency condition can be found in the English version of Wikipedia.
(Japanese) https://ja.wikipedia.org/wiki/%E3%83%9F%E3%83%A9%E3%83%BC%E2%80%93%E3%83%A9%E3%83%93%E3%83%B3%E7%B4%A0%E6%95%B0%E5%88%A4%E5%AE%9A%E6%B3%95
(English) https://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test
Now, let's finally explain the prime factorization. First, I will introduce the Floyd cycle detection method.
Consider the map $ f: A \ rightarrow A $ for a finite set $ A $ consisting of $ n $ elements. Given $ a, \ f (a), \ f (f (a)), \ ... $ for $ a \ in A
For simplicity, let's say $ n $ is represented as $ n = pq $ as the product of $ 2 $ prime factors. Also, by dividing the small prime factors in advance, $ p $ and $ q $ can be reasonably large (for example, $ 100 $ or more). $ F: \ mathbb {Z} / n \ mathbb {Z} \ rightarrow \ mathbb {Z} / n \ mathbb {Z} $ as a pseudo-random number $ f (a) = a ^ {2} + c \ ( It is defined by {\ rm mod} \ n) \ tag {3} $$. From $ n = pq $, this can also be seen as a function on $ \ mathbb {Z} / p \ mathbb {Z} $ or a function on $ \ mathbb {Z} / q \ mathbb {Z} $. .. Now, using the Floyd cycle detection method I mentioned earlier, this circulates somewhere. From the birthday paradox, this should be $ O (\ sqrt p) $ in $ \ mathbb {Z} / p \ mathbb {Z} $ and you should find a cycle. The same is true for $ \ mathbb {Z} / q \ mathbb {Z} $.
How can I find out if I circulated in $ \ mathbb {Z} / p \ mathbb {Z} $? This can be done using $ \ rm GCD
In the above method, taking $ \ rm GCD $ was the bottleneck of computational complexity. Therefore, we will consider speeding up by judging multiple values at the same time. Specifically, $ m $ is decided appropriately, the number of $ \ rm GCD $ to be taken is multiplied by all the $ m $ pieces, and the product is multiplied by $ n $ of $ \ rm GCD. You can judge at once by taking $. If you happen to find $ p $ and $ q $ at the same time, you can go back a bit and move forward by $ 1 $. If you are unlucky and it is detected at the same time, change the pseudo-random number and try again. Specifically, you can try changing $ c $ of $ (3) $ appropriately.
If $ m $ is set to about $ m = n ^ {\ frac {1} {8}} $, the number of times to take $ \ log $ is $ O (n ^ {\ frac {1} {8}} ) $, And the number of multiplications and divisions is $ O (n ^ {\ frac {1} {4}}) $, so $ O (n ^ {\ frac {1} {4}} + n ^ as a whole You can factor in prime factors with {\ frac {1} {8}} \ log n) = O (n ^ {\ frac {1} {4}}) $.
factorize.py
def gcd(a, b):
while b: a, b = b, a % b
return a
def isPrimeMR(n):
d = n - 1
d = d // (d & -d)
L = [2]
for a in L:
t = d
y = pow(a, t, n)
if y == 1: continue
while y != n - 1:
y = (y * y) % n
if y == 1 or t == n - 1: return 0
t <<= 1
return 1
def findFactorRho(n):
m = 1 << n.bit_length() // 8
for c in range(1, 99):
f = lambda x: (x * x + c) % n
y, r, q, g = 2, 1, 1, 1
while g == 1:
x = y
for i in range(r):
y = f(y)
k = 0
while k < r and g == 1:
ys = y
for i in range(min(m, r - k)):
y = f(y)
q = q * abs(x - y) % n
g = gcd(q, n)
k += m
r <<= 1
if g == n:
g = 1
while g == 1:
ys = f(ys)
g = gcd(abs(x - ys), n)
if g < n:
if isPrimeMR(g): return g
elif isPrimeMR(n // g): return n // g
return findFactorRho(g)
def primeFactor(n):
i = 2
ret = {}
rhoFlg = 0
while i*i <= n:
k = 0
while n % i == 0:
n //= i
k += 1
if k: ret[i] = k
i += 1 + i % 2
if i == 101 and n >= 2 ** 20:
while n > 1:
if isPrimeMR(n):
ret[n], n = 1, 1
else:
rhoFlg = 1
j = findFactorRho(n)
k = 0
while n % j == 0:
n //= j
k += 1
ret[j] = k
if n > 1: ret[n] = 1
if rhoFlg: ret = {x: ret[x] for x in sorted(ret)}
return ret
Recommended Posts