COME ON CODE ON

A blog about programming and more programming.

Archive for September 2010

Pollard Rho Brent Integer Factorization

with 4 comments

Pollard Rho is an integer factorization algorithm, which is quite fast for large numbers. It is based on Floyd’s cycle-finding algorithm and on the observation that two numbers x and y are congruent modulo p with probability 0.5 after 1.177\sqrt{p} numbers have been randomly chosen.

Algorithm
Input : A number N to be factorized
Output : A divisor of N
If x mod 2 is 0
	return 2

Choose random x and c
y = x
g = 1
while g=1
	x = f(x)
	y = f(f(y))
	g = gcd(x-y,N)
return g

Note that this algorithm may not find the factors and will return failure for composite n. In that case, use a different f(x) and try again. Note, as well, that this algorithm does not work when n is a prime number, since, in this case, d will be always 1. We choose f(x) = x*x + c. Here’s a python implementation :

def pollardRho(N):
        if N%2==0:
                return 2
        x = random.randint(1, N-1)
        y = x
        c = random.randint(1, N-1)
        g = 1
        while g==1:             
                x = ((x*x)%N+c)%N
                y = ((y*y)%N+c)%N
                y = ((y*y)%N+c)%N
                g = gcd(abs(x-y),N)
        return g

In 1980, Richard Brent published a faster variant of the rho algorithm. He used the same core ideas as Pollard but a different method of cycle detection, replacing Floyd’s cycle-finding algorithm with the related Brent’s cycle finding method. It is quite faster than pollard rho. Here’s a python implementation :
def brent(N):
        if N%2==0:
                return 2
        y,c,m = random.randint(1, N-1),random.randint(1, N-1),random.randint(1, N-1)
        g,r,q = 1,1,1
        while g==1:             
                x = y
                for i in range(r):
                        y = ((y*y)%N+c)%N
                k = 0
                while (k<r and g==1):
                        ys = y
                        for i in range(min(m,r-k)):
                                y = ((y*y)%N+c)%N
                                q = q*(abs(x-y))%N
                        g = gcd(q,N)
                        k = k + m
                r = r*2
        if g==N:
                while True:
                        ys = ((ys*ys)%N+c)%N
                        g = gcd(abs(x-ys),N)
                        if g>1:
                                break
        
        return g    

-fR0DDY

Written by fR0DDY

September 18, 2010 at 11:51 PM

Miller Rabin Primality Test

with 13 comments

Miller Rabin Primality Test is a probabilistic test to check whether a number is a prime or not. It relies on an equality or set of equalities that hold true for prime values, then checks whether or not they hold for a number that we want to test for primality.

Theory

1> Fermat’s little theorem states that if p is a prime and 1 ≤ a < p then a^{p-1}\equiv1\pmod{p}.
2> If p is a prime and x^2\equiv1\pmod{p} or \left(x-1\right)\left(x+1\right)\equiv0\pmod{p} then x\equiv1\pmod{p} or x\equiv-1\pmod{p}.
3> If n is an odd prime then n-1 is an even number and can be written as 2^{s}.d. By Fermat’s Little Theorem either a^{d}\equiv1\pmod{n} or a^{2^r\cdot d}\equiv -1\pmod{n} for some 0 ≤ r ≤  s-1.
4> The Miller–Rabin primality test is based on the contrapositive of the above claim. That is, if we can find an a such that a^{d}\not\equiv1\pmod{n} and a^{2^r\cdot d}\not\equiv -1\pmod{n} for all 0 ≤ r ≤  s-1 then a is witness of compositeness of n and we can say n is not a prime. Otherwise, n may be a prime.
5> We test our number N for some random a and either declare that N is definitely a composite or probably a prime. The probably that a composite number is returned as prime after k itereations is 4^{-k}.

Algorithm

Input :A number N to be tested and a variable iteration-the number 
of 'a' for which algorithm will test N.
Output :0 if N is definitely a composite and 1 if N is probably a prime.

Write N as 2^{s}.d
For each iteration
	Pick a random a in [1,N-1]
	x = a^{d} mod n
	if x =1 or x = n-1 
		Next iteration
	for r = 1 to s-1
		x  = x^{2} mod n
		if x = 1
			return false
		if x = N-1
			Next iteration
	return false
return true

Here’s a python implementation :

import random
def modulo(a,b,c):
        x = 1
        y = a
        while b>0:
                if b%2==1:
                        x = (x*y)%c
                y = (y*y)%c
                b = b/2
        return x%c
        
def millerRabin(N,iteration):
        if N<2:
                return False
        if N!=2 and N%2==0:
                return False
        
        d=N-1
        while d%2==0:
                d = d/2
        
        for i in range(iteration):
                a = random.randint(1, N-1)
                temp = d
                x = modulo(a,temp,N)
                while (temp!=N-1 and x!=1 and x!=N-1):
                        x = (x*x)%N
                        temp = temp*2
                
                if (x!=N-1 and temp%2==0):
                        return False
        
        return True

-fR0DDY

Written by fR0DDY

September 18, 2010 at 12:23 PM

Follow

Get every new post delivered to your Inbox.

Join 99 other followers