Pollard Rho factorization method implementation in C

Question

Can anyone help me out with the pollard rho implementation? I have implemented this in C. It's working fine for numbers upto 10 digits but it's not able to handle greater numbers.

Please help me out to improve it to carry out factorization of numbers upto 18 digits . My code is this:

#include<stdio.h>
#include<math.h>
int gcd(int a, int b)
{
    if(b==0) return a ;
    else
    return(gcd(b,a%b)) ;
}

long long int mod(long long int a , long long int b , long long int n )
{    
     long long int x=1 , y=a ;
     while(b>0)
     {
        if(b%2==1)  x = ((x%n)*(y%n))%n ;
        y = ((y%n)*(y%n))%n ;
         b/=2 ;
     }
     return x%n ;
}

int isprimes(long long int u)
{  
    if(u==3)
    return 1 ;
     int a = 2 , i ;
     long long int k , t = 0 , r , p ;
     k = u-1 ;
     while(k%2==0)
     { k/=2 ; t++ ; }

         while(a<=3)                                                              /*der are no strong pseudoprimes common in base 2 and base 3*/
         {   
         r = mod(a,k,u) ;
             for(i = 1 ; i<=t ; i++)
             {
                  p = ((r%u)*(r%u))%u ;
                  if((p==1)&&(r!=1)&&(r!=(u-1)))
                  {  return 0 ; }
                  r = p ;
             }
          if(p!=1)
          return 0 ;
         else
          a++ ;
         } 

          if(a==4)
          return 1 ;

}

long long int pol(long long int u)
{ 
  long long int x = 2 , k , i , a , y , c , s;
  int d = 1 ;
   k = 2 ;
   i = 1 ;
   y = x ;
   a = u ;
   if(isprimes(u)==1)
   { 
     return 1;
   }
   c=-1 ;
   s = 2 ;
   while(1)
   {
     i++;
     x=((x%u)*(x%u)-1)% u ;

     d = gcd(abs(y-x),u) ;

     if(d!=1&&d!=u)
     { printf("%d ",d);
       while(a%d==0) { a=a/d;  }

        x = 2 ;
        k = 2 ;
        i = 1 ;
        y = x ;
        if(a==1)
        { return 0 ; }
        if(isprimes(a)!=0)
        { return a ; }
        u=a ;

     }
     if(i==k)
     {y = x ; k*=2 ; c = x ;}                                                       /*floyd cycle detection*/
        if(c==x)                                                                 
     { x = ++s ; }
    }
    return ;

}

int main()
{
   long long int t ;
    long long int i , n , j , k , a , b  , u ;
    while(scanf("%lld",&n)&&n!=0)
    { u = n ; k = 0 ;
    while(u%2==0)
       {  u/=2 ; k = 1 ; }
      if(k==1) printf("2 ") ;
      if(u!=1)
      t = pol(u) ;
        if(u!=1) 
      {
           if(t==1)
           { printf("%lld",u) ; }
           else
           if(t!=0)
           { printf("%lld",t) ; }
      }
          printf("\n");
    }
    return 0;
}

sorry for the long code ..... I am a new coder.

Answer 1

You can try this ~100 lines C implementation of Pollard Rho :

There is some helpers :

#include <stdlib.h>
#include <stdint.h>

typedef uint_fast64_t num ;

static inline num mul_mod(num a, num b, const num mod) {
    // Return (a * b) % mod, avoiding overflow while doing modular multiplication.
    num res = 0, tmp;
    for (b %= mod; a; a & 1 ? b >= mod - res ? res -= mod : 0, res += b : 0, a >>= 1, (tmp = b) >= mod - b ? tmp -= mod : 0, b += tmp);
    return res % mod;
}

static inline num square_root(num n) {
    // Return the number that was multiplied by itself to reach N.
    num a = 0, b, c;
    for (b = 1ULL << 62; b; c = a + b, n -= c &= -(n >= c), a = (a >> 1) | (c & b), b >>= 2);
    // Variable n contains the remainder.
    return a;
}

There is a required prime checker :

static  int is_prime(const num n, size_t iterations) {
    // Perform a Miller-Rabin test.
    num a = 0, b, c, d, e, f; int h, i;
    if ((n == 1) == (n & 1)) return n == 2;
    for (b = c = n - 1, h = 0; !(b & 1); b >>= 1, ++h);
    for (; iterations--;) {
        for (size_t g = 0; g < sizeof(a); ((char*)&a)[g++] = rand()); // rand input.
        do for (d = e = 1 + a % c, f = n; (d %= f) && (f %= d););
        while (d > 1 && f > 1);
        for (d = f = 1; f <= b; f <<= 1);
        for (; f >>= 1; d = mul_mod(d, d, n), f & b && (d = mul_mod(e, d, n)));
        if (d == 1) continue;
        for (i = h; i-- && d != c; d = mul_mod(d, d, n));
        if (d != c) return 0;
    }
    return 1;
}

There is two factorization functions :

static inline num factor_worker_2(const num n, size_t limit) {
    // Perform a Pollard's Rho test.
    size_t a = -1, b = 2;
    num c, d = 1 + rand(), e = 1, f = 1;
    for (c = d %= n; f == 1 && --limit; d = c, b <<= 1, a = -1) {
        for (; f |= e, f == 1 && ++a != b;) {
            c = mul_mod(c, c, n);
            for (++c, c *= c != n, e = n, f = c > d ? c - d : d - c; (f %= e) && (e %= f););
        }
    }
    return f;
}

static inline num factor_worker_1(const num n) {
    // Perform a trial divisions test on N.
    static const num list[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 1};
    size_t i;
    for (i = -1; n % list[++i];);
    return list[i];
}

There is a factorizarion manager :

num factor(const num n) {
// Factorization manager, detect primes, perfect squares, execute workers.
    num res;
    switch (n) {
        case 0: case 1: case 2: case 3:
            res = 1; break;
        default:
            res = factor_worker_1(n);
            if (res == 1 && !is_prime(n, 20)) {
                res = square_root(n);
                if (res * res != n)
                    for(;res = factor_worker_2(n, -1), res == 1 || res == n;);
            }
    }
    return res;
}

There is a main :

#include <assert.h>
#include <stdio.h>

int main(void) {
    num N;
    N = 14898780622628181263U;
    printf("factor is %zu\n", factor(N));
}

To try it :

// You can put it into a main.c file then compile + execute :
// gcc -O3 -std=c99 -Wall -pedantic main.c ; ./a.out ;

Here is the source for more informations, Thank You.