#include <iostream.h>
#include <gmpxx.h>

class Tester {
public:
  mpz_class pow(mpz_class x, mpz_class e, mpz_class m) ;
  mpz_class exp (mpz_class x, mpz_class e) ;
  int fermat_test (mpz_class n, mpz_class b) ;
  int miller_rabin (mpz_class n, mpz_class b) ;
  mpz_class gcd(mpz_class m, mpz_class n) ;
  mpz_class trial_division(mpz_class n) ;
  mpz_class next_prime(mpz_class n) ;
  mpz_class pollard_rho (mpz_class n) ;
};

/* to compile, do   g++ -c -o <file>.o <file>.C -lgmpxx   */

mpz_class Tester::pow(mpz_class x, mpz_class e, mpz_class m) {
  mpz_class X = x;
  mpz_class E = e;
  mpz_class Y = 1;
  while (E>0) {
	 if (E % 2 == 0) {
		X = (X * X) % m;
		E = E/2;
	 }
	 else {
		Y = (X * Y) % m;
		E = E - 1;
	 }
  }
  return Y;
}

mpz_class Tester::exp (mpz_class x, mpz_class e) {
  mpz_class X = x;
  mpz_class E = e;
  mpz_class Y = 1;
  while (E>0) {
	 if (E % 2 == 0) {
		X = X * X;
		E = E/2;
	 }
	 else {
		Y = X * Y;
		E = E - 1;
	 }
  }
  return Y;
}

int Tester::fermat_test (mpz_class n, mpz_class b) {
  if ( pow(b,n,n) == b ) {
	 return 1;
  }
  else {
	 return 0;
  }
}

int Tester::miller_rabin (mpz_class n, mpz_class b) {
  mpz_class m = n-1;
  long int s = 0;
  while ( m % 2 == 0) {
	 s++;
	 m /= 2;
  }
  mpz_class B = pow(b,m,n);
  if (B == 1 || B == n-1) {
	 return 1;
  }
  long int t = 0;
  while (t < s) {
	 t++;
	 B = (B * B) % n;
	 if (B == 1) {
		return 0;
	 }
	 else if (B == n-1) {
		return 1;
	 }
  }
  return 0;
}

mpz_class Tester::gcd(mpz_class m, mpz_class n) {
  
  mpz_class r;

  if (m < 0) {
	 m = -m;
  }
  if (n < 0) {
	 n = -n;
  }

  r = m % n;

  if (r == 0) {
	 return n;
  }

  while(r > 0) {
    m = n;
	 n = r;
	 r = m % n;
  }
  return n;
}

/******************************/

mpz_class Tester::trial_division(mpz_class n) {
  if (n % 2 == 0) {
	 return 2;
  }
  mpz_class d = 3;
  while (d * d <= n) {
	 if (n % d == 0) {
		return d;
	 }
	 d += 2;
  }
  return 1;
}

/************************************************/

mpz_class Tester::next_prime(mpz_class n) {
  if (n % 2 == 0) {
	 n++;
  }
  while ( fermat_test(n,2) == 0 ||
			 miller_rabin(n,2) == 0 || 
			 miller_rabin(n,3) == 0 || 
			 miller_rabin(n,5) == 0) {
	 n += 2;
  }
  return n;
}

/**************************************************/

mpz_class Tester::pollard_rho (mpz_class n) {
  mpz_class x,y;
  x = 2;
  y = 2;

  mpz_class g = 1;
  
  int ct = 0;

  while (g == 1) {
	 ct++;
	 x = (x * x + 2) % n;
	 y = (y * y + 2) % n;
	 y = (y * y + 2) % n;
	 g = gcd(x-y,n);
  }

  cout << "... now on step " << ct << " in Pollard rho\n";
  return g;
}

/**********************************************************************/