easy_factor.cc

Go to the documentation of this file.

#include "easy_factor.H"
#include <cmath>
#include "modulo.H"
#include "elliptic_curve.H"
#include <sstream>

extern std::string MAL(const mpz_t factor, const std::ostringstream &comment);
extern std::string MAL(const mpz_t factor, const unsigned int exponent, const std::ostringstream &comment);
extern std::string MAL(const mpz_t factor, const unsigned int exponent=1);
extern std::string MAL(const unsigned long int factor, const unsigned int exponent=1);
extern void tune_parameters(const unsigned int Dezimalstellen);

#include <fstream>
extern std::ofstream Factorization_to_file;

#ifdef REACT_ON_SIGUSR
 #include "usr_signals.H"
 extern Cusr_signal_proxy USRSignalHandler;
#endif

class ExitManager
{
 public:
  static void StopFactorization();
};

// propagate DEFINE condition for IS_STANDALONE from compiletime
// to runtime...
extern bool compiled_IS_STANDALONE();



void pollard(const int runden) /* Pollard-rho-Methode */
{
#ifdef VERBOSE_INFO
  cout << "Trying Pollard-rho up to " << runden << endl;
#endif
  mpz_t a,a2,x;
  mpz_init(a); mpz_init(a2); mpz_init(x);
  mpz_set_ui(a,1); mpz_set(a2,a);
  
  for (int i=1; i<=runden; i++)
   {
#ifdef REACT_ON_SIGUSR
     if (USRSignalHandler.got_SIGUSR1()) break;
     if (USRSignalHandler.got_SIGUSR2()) break;
#endif
     mpz_mul(a,a,a); mpz_add_ui(a,a,1); mpz_mod(a,a,n);
     mpz_mul(a2,a2,a2); mpz_add_ui(a2,a2,1); mpz_mod(a2,a2,n);
     mpz_mul(a2,a2,a2); mpz_add_ui(a2,a2,1); mpz_mod(a2,a2,n);
     mpz_sub(x,a2,a);
     mpz_gcd(x,x,n);
     while (mpz_cmp_ui(x,1)!=0)
      {
        if (mpz_cmp(x,n)==0)
          {
#ifdef VERBOSE_INFO
            cout << "trivial factor in round " << i << endl;
#endif
            continue;
          }
        unsigned int exp=mpz_remove(n,n,x); // divide out all occurrences of x in n
        if (mpz_probab_prime_p(x,probab_prime_checks))
         {
           cout << x << " is factor. [round " << i <<"]" << endl;
           Factorization_to_file << MAL(x,exp) << flush;
         }
        else
         {
           cout << x << " is a composite factor. [round " << i <<"]" << endl;
           std::ostringstream comment;
           comment << " [rho/composite]";
           Factorization_to_file << MAL(x,exp,comment) <<  flush;
         }
        if (mpz_probab_prime_p(n,probab_prime_checks)) goto done; // remaining number is most probably prime
        mpz_gcd(x,x,n); // possible case: x=p*q and n=p*q*q*r -> second q would not be found!
      }
   };
done:
  mpz_clear(a); mpz_clear(a2); mpz_clear(x);
}

bool try_memorized_factors(mpz_t n, const std::string& memfilename)
{
  // trial division by some factors given in a file...
  bool reduced = false;

  std::ifstream in(memfilename.c_str());
  if (in)
    {
#ifdef VERBOSE_INFO
      cout << "Trial division by factors in " << memfilename << endl; 
#endif
      mpz_t x;
      mpz_init(x);
      while (!in.eof())
        {
          while (isspace(in.peek())) in.get(); // ignore all whitespaces
          if (in.peek()=='#')
           {
             in.ignore(1000,'\n'); // ignore comment lines
             continue;
           }

          // it seems that gcc-3.2.3 has some problems detecting eof,
          // I hope this helps...
          in.get();
          if (in.eof()) break;
          in.unget();

          in >> x;
          //cout << "try " << x << endl;
          if (mpz_cmp(n,x)==0) return reduced; // remaining factor should reside in n!
          while (mpz_divisible_p(n,x))
           {
             if (mpz_probab_prime_p(x,probab_prime_checks))
               {
                 cout << x << " is factor. [memorized]" << endl;
                 std::ostringstream comment;
                 comment << " [memorized]";
                 Factorization_to_file << MAL(x,comment) << flush;
               }
             else
               {
                 cout << x << " is a composite factor. [memorized]" << endl;
                 std::ostringstream comment;
                 comment << " [composite/memorized]";
                 Factorization_to_file << MAL(x,comment) <<  flush;
               }
             mpz_divexact(n,n,x); reduced=true;
           }
          while (isspace(in.peek())) in.get(); // ignore all whitespaces
        }
      mpz_clear(x);
    }
  else
    {
#ifdef VERBOSE_WARN
      cerr << memfilename << " inaccessible!" << endl;
#endif
    }
 return reduced;
}

#include "pollard_phi.cc"
#include "fibonacci_ppm1.cc"

bool Potenztest(const mpz_t n)
{
  // check whether n can be written as a square, a cube or
  // any other x-th power of a natural number.
  mpz_t y;
  mpz_init(y);
  for (unsigned int k=mpz_sizeinbase(n,256); k>1; k--)
    {
      //cout << "k=" << k << "?" << endl;
      if (mpz_root(y,n,k))
        {
          cout << y << "^" << k <<  " is ";
          Factorization_to_file << MAL(y,k);
          if (!mpz_probab_prime_p(y,probab_prime_checks))
            {
              cout << "composite ";
              Factorization_to_file << " [composite power]";
            }
          Factorization_to_file << flush;
          cout << "factor." << endl;
          mpz_clear(y);
          return true; 
        }
    }
  mpz_clear(y);
  return false;
}


#include "fermat.cc"

void easy_factor()
{
  /* einfache Faktorisierung versuchen:
     vielleicht kann die Zahl dadurch komplett faktorisiert werden
     oder zumindest um einige Faktoren reduziert werden... */

  if (!SkipFermat) fermat_like_method();
  /* try a Fermat-like method to check for cases, where the
     number can be written as (a-d)*(a+d) using a small d... */
  
  const unsigned long int UpperLimit_Easy = 1000000;
#ifdef VERBOSE_INFO
  cout << "Trialdivision up to " << UpperLimit_Easy << endl;
#endif
  {
    const unsigned int some_primes [] = { 2,3,5,7,11,13,0 };
    unsigned long int p = 0;

    // trial division using small prime numbers
    for (unsigned int index=0; some_primes[index]; ++index)
     {
       p=some_primes[index];
       //cout << "Trial Division with " << p << endl;
       int exp=0;
       while (mpz_divisible_ui_p(n,p)) { mpz_divexact_ui(n,n,p); exp++; }
       if (exp>0)
        {
          cout << p << "^" << exp << " is factor." << endl;
          Factorization_to_file << MAL(p,exp);
        }
     }
    
    /* Trial division up to UpperLimit_Easy */
    unsigned int delta;
    if (p%6==1) delta=4;
     else 
      if (p%6==5) delta=2;
      else
       {
         MARK;
         cerr << "Something is wrong with some_primes!" << endl;
         exit(1);
       }
      
    while (p<UpperLimit_Easy)
     {
       p+=delta; delta=6-delta;
       if ( p%5==0 || p%7==0 || p%11==0 || p%13==0 ) continue;
       //cout << "Trial Division with " << p << endl;
       int exp=0;
       while (mpz_divisible_ui_p(n,p)) { mpz_divexact_ui(n,n,p); exp++; }
       if (exp>0)
         {
           cout << p << "^" << exp << " is factor." << endl;
           Factorization_to_file << MAL(p,exp);
           if (mpz_cmp_ui(n,1)==0) break; // work is done!
         }
     } 
    
    Factorization_to_file.flush();
  }

  if (Potenztest(n)) mpz_set_ui(n,1);
  if (mpz_cmp_ui(n,1)==0) goto done; 
  if (mpz_probab_prime_p(n,probab_prime_checks)) goto done;

  { 
    std::string cached_factorfile;
 #if defined(PKGDATADIR)
    cached_factorfile = PKGDATADIR; // use this prefix (if defined)
    cached_factorfile+="/";
 #endif
   cached_factorfile+= "fibonacci.factors";

   if (try_memorized_factors(n,cached_factorfile)) // "cached" factors
    {
      if (mpz_probab_prime_p(n,probab_prime_checks)) goto done;
    }
  }

  {
    unsigned int dezimalstellen = mpz_sizeinbase(n,10);
    tune_parameters(dezimalstellen);
    pollard(rho_Phase); // try the pollard-rho factorization method

    if (Potenztest(n)) mpz_set_ui(n,1);
    if (mpz_cmp_ui(n,1)==0) goto done;
    if (mpz_probab_prime_p(n,probab_prime_checks)) goto done;

    if (dezimalstellen != mpz_sizeinbase(n,10))
     {
       dezimalstellen=mpz_sizeinbase(n,10);
       tune_parameters(dezimalstellen);
     }

  if (UsePhimat)
   {
     // this is an experimental method!
     // Don't use it by default since it takes very long
     // (except for certain numbers)
     phimat2(); if (mpz_cmp_ui(n,1)==0) goto done;
   }

#if 1
   if (dezimalstellen>200)
    {
      // number seems to be very large
      // if it partly contains fibonacci factors, we should try to split it
      // here to save time.  (Since I'm trying to factor some fibonacci &
      // lucas numbers, this test is implemented.)
      fibonacci_ppm1(100000,0); /* Fibonacci (p+-1)-Methode */
      if (dezimalstellen != mpz_sizeinbase(n,10))
       {
         if (Potenztest(n)) mpz_set_ui(n,1);
         if (mpz_cmp_ui(n,1)==0) goto done;
         if (mpz_probab_prime_p(n,probab_prime_checks)) goto done;
         dezimalstellen=mpz_sizeinbase(n,10);
         tune_parameters(dezimalstellen);
       }
    }
#endif


#if 1
    // in the standalone version the mpqs will be quite fast in finding
    // small factors, so we do not need to run ecm on small composite numbers.
    // However, in a networked client-server system, factoring small composite
    // numbers causes too much overhead and is therefore annoying...

    // to make it separately compileable, we must not evaluate the DEFINE "IS_STANDALONE";
    // we have it to evaluate on runtime, not at compiletime!
    if ( !compiled_IS_STANDALONE() || (compiled_IS_STANDALONE() && mpz_sizeinbase(n,10)>70))
     if (!SkipEasyECM)
      {
        // easy-ecm: try a few elliptic curves just to eliminate some small numbers
        // in the range up to let's say 10-20 digits.
        // This is especially advantageous for composite Fibonacci-
        // and Lucasnumbers, because fibfactor has the chance to find the
        // remaining cofactors!
        for (int i=0; i<10; ++i)
         {
           elliptic_curves elliptic_curve; // instantiate an elliptic curve
           elliptic_curve.go(7+i,60000,15000000); // and try a factorization using this curve

           if (dezimalstellen != mpz_sizeinbase(n,10))
            {
              if (Potenztest(n)) mpz_set_ui(n,1);
              if (mpz_cmp_ui(n,1)==0) goto done;
              if (mpz_probab_prime_p(n,probab_prime_checks)) goto done;
              dezimalstellen=mpz_sizeinbase(n,10);
              tune_parameters(dezimalstellen);
            }
         }
      }
#endif


#if 1
    if (SkipPhi || (phi_Phase1<=0 && phi_Phase2<=0))
     {
#ifdef VERBOSE_NOTICE
       cout << "Phi-Phase is deactivated!" << endl;
#endif
     }
    else
     {
       pollard_phi(phi_Phase1,phi_Phase2); /* Faktorisierung nach Pollard-Phi-Methode versuchen */
       if (dezimalstellen != mpz_sizeinbase(n,10))
        {
          if (Potenztest(n)) mpz_set_ui(n,1);
          if (mpz_cmp_ui(n,1)==0) goto done;
          if (mpz_probab_prime_p(n,probab_prime_checks)) goto done;
          dezimalstellen=mpz_sizeinbase(n,10);
          tune_parameters(dezimalstellen);
        }
     }
#endif

    if (SkipFibonacci || (phi_Phase1<=0 && phi_Phase2<=0))
     {
#ifdef VERBOSE_NOTICE
       cout << "Fibonacci-Phase is deactivated!" << endl;
#endif
     }
    else
     {
       fibonacci_ppm1(phi_Phase1,phi_Phase2); /* Fibonacci (p+-1)-Methode */
       if (dezimalstellen != mpz_sizeinbase(n,10))
        {
          if (Potenztest(n)) mpz_set_ui(n,1);
          if (mpz_cmp_ui(n,1)==0) goto done;
          if (mpz_probab_prime_p(n,probab_prime_checks)) goto done;
          dezimalstellen=mpz_sizeinbase(n,10);
          tune_parameters(dezimalstellen);
        }
     }
  }

done:
#ifdef VERBOSE_NOTICE  
  cout << "remaining number:" << endl << n << endl;
#endif
  
  if ( (mpz_cmp_ui(n,1)==0) || Potenztest(n) )
    {
      Factorization_to_file << endl;
      ExitManager::StopFactorization(); // factorization complete!
    }  

  if (mpz_probab_prime_p(n,probab_prime_checks))
    {
      Factorization_to_file << MAL(n) << endl;
      cout << "remaining number is most probably a prime!" << endl;
      mpz_set_ui(n,1);
      ExitManager::StopFactorization();
    }
//Factorization_to_file << " * " << n << " [REST]" << endl; exit(0);
}

---