qsieve.cc

Go to the documentation of this file.

#include "at_startup.H"

// sanity checks
#ifdef IS_CLIENT
 #error "qsieve is not a client!"
#endif

#ifdef IS_SERVER
 #error "qsieve is not a server!"
#endif

#ifdef USE_NETWORK
 #error "qsieve uses no network features!"
#endif

#ifndef IS_STANDALONE
 #error "qsieve is a standalone application!"
#endif



extern "C"
{
 #include <unistd.h>
}

#include <cstdio>
#include <cstdlib>
#include <fstream>
#include <sstream>
#include <cmath>
#include <gmp.h>
#include <list>
#include <set>
#include <ctime>

#include "qsieve.H"
#include "StaticFactorbase.H"


#include "FactorFound.H"
TFoundFactors FoundFactors;



#include <vector>
#include <stack>
#include <algorithm>


#ifdef NOTIFY_PARENT
 #warning "undefining NOTIFY_PARENT, since this is not a networked server"
 #undef NOTIFY_PARENT
#endif




const string RecoveryFile         = "recovery.dat";
const string FoundFactorsFile     = "FoundFactors.dat";
const string StaticRelationsFile  = "static_relations.dat";
const string SpecialRelationsFile = "special_relations.dat";


// provide streams for storing relations and other data
// *************************************************************************
// **** important:
// **** If "ostream" and "istream" are handled by the same filebuf,
// **** then you *must not* rely on 
// **** "tellp()" and "tellg()" behave independently from each other!
// **** Any write operation may change also the value of "tellg()"
// **** and any read operation may change the value of "tellp()"!
// ************************************************************************* 


// for Recovery
filebuf Recovery_buffer;
ostream Recovery_to_file (&Recovery_buffer);
istream Recovery_from_file (&Recovery_buffer);



#include "mpqsPolynom.H"


mpz_t n, // number to factorize (will be reduced during factorization)
      kN; // input for MPQS (includes a suitable multiplier)


#include "StaticRelations.H"
#include "DynamicRelations.H"
#include "SpecialRelations.H"

#include "Sieving.H"
#include "ConfigFile.cc"


CmpqsPolynom Polynom; // object to manage polynomial computations for multipolynomial sieve (MPQS)
#include "mpqsStatistics.cc" // statistical stuff about found relations, client connections, etc.


#include "modulo.H" // modulo operations for unsigned int
using namespace numtheory;

#include "FactorFound.cc"
#include "polynomial.H" /* for fast polynomial arithmetic & fft-continuation */
#include "fft_param.cc"  /* discrete fast fourier transform */
#include "elliptic_curve.H" /* invariant part of elliptic curve method (ecm) */
#include "elliptic_curve-variant.cc" /* variant part of elliptic curve method (ecm) */

#include "ExitManager.cc" /* controlled termination of program */
#include "StaticRelations.cc"

// this stuff may only be needed by server or stand-alone version
#include "easy_factor.H"


// this is a predicate stating "false"
const bool without_multi_combine_init = false; // true ->multi-init, false ->no multi_init on StaticRelations::insert

#include "CRelation-inc.cc"
#include "SpecialRelations.cc"

#include "parse_term.cc"


void process_ecm()
{
  streampos fpos = Recovery_from_file.tellg();
  int Anzahl_Kurven = 0;
  string s;

  srand(time(NULL)); // setting a random seed (time(NULL) should satisfy for our needs)
  
  unsigned int dezimalstellen = mpz_sizeinbase(n,10);
  tune_parameters(dezimalstellen); // as the name says: tuning of parameters according to our input number

  Recovery_from_file >> Anzahl_Kurven;
  Recovery_from_file.ignore(1,'\n'); // overread "end of line"

  {
   ifstream in("ecm-recover.dat");
   if (in) // does a recovery file exist?
    {
     in.close();
     elliptic_curves elliptic_curve; // create a curve
     elliptic_curve.go(-1,elcu_Phase1,elcu_Phase2); // try to recover curve and retry factorization with elliptic curve
    }
  }

  while (Anzahl_Kurven<elcu_Kurven)
  {

    Anzahl_Kurven++;
    cout << "elliptic curve #" << Anzahl_Kurven << "/" << elcu_Kurven << endl;

    int ecm_sigma = rand(); // choose a random curve
    elliptic_curves elliptic_curve; // create curve
    elliptic_curve.go(ecm_sigma,elcu_Phase1,elcu_Phase2); // try factorization with elliptic curve

    // write number of processed curves to recovery file
    Recovery_to_file.seekp(fpos); 
    Recovery_to_file << Anzahl_Kurven << endl << flush;

    // found a factor?
    if (dezimalstellen==mpz_sizeinbase(n,10)) continue; // no new factor
    
    // digits of n has changed -> new factor

    if (mpz_probab_prime_p(n,probab_prime_checks))
      {
        Factorization_to_file << MAL(n);
        cout << "remaining factor: " << n << endl;
        cout << "remaining factor is most probably prime!" << endl;
        mpz_set_ui(n,1); // n has been factorized
      }

    if ( (mpz_cmp_ui(n,1)==0) || Potenztest(n) ) // factorization complete?
      { 
        cout << "factorization successfully completed." << endl;
        Factorization_to_file << endl;
        ExitManager::StopFactorization(); // terminate program (do finalize work)
      }

    // tune parameters for further ecm-rounds...
    dezimalstellen=mpz_sizeinbase(n,10);
    tune_parameters(dezimalstellen);

    // since the number to factorize has become smaller, the recovery-file need to be updated as well!
    Recovery_buffer.close(); // close file
    Recovery_buffer.open(RecoveryFile.c_str(),ios::in|ios::out|ios::trunc); // re-open (-> create empty file)
    Recovery_to_file << n << endl; // write out number to factorize
    fpos = Recovery_from_file.tellg(); // mark position
    Recovery_to_file << Anzahl_Kurven << endl << flush; // store number of curves done so far
    ecm_curves_processed=Anzahl_Kurven;

  }
  cout << "ECM-Phase (standalone) completed." << endl;
}




void cleanup_files()
{
#ifdef VERBOSE_INFO
  cout << "cleaning files..." << endl;
#endif
  StaticRelations::cleanup_files();
  SpecialRelations::cleanup_files();
  DynamicRelations::cleanup_files();
  Recovery_buffer.close();
  remove(RecoveryFile.c_str());
}

void cleanup_memory()
{
#ifdef VERBOSE_INFO
  cout << "cleanup allocated memory" << endl;
#endif
  StaticRelations::cleanup_memory();
  CSieveStaticFactors::destruct();
  mpz_clear(CmpqsFactor::DLP_Threshold);
  mpz_clear(kN); mpz_clear(n);
}




int main(const int argc, const char* const argv[])
{

#ifdef USE_NCURSES
  new Cncursed(); // trigger activation of ncursed streams
#endif

  PrintHeader("Qsieve standalone");
    
  if (argc>2)
    {
      cerr << "number for factorization expected!" << endl;
      exit(1);
    }
  
  const bool recover = (argc==1); // without argument: Recovery-Mode

  cout.setf(ios::fixed); // fixed decimal notation, not scientific notation!

  mpz_init(n); // our number to factorize
  mpz_init(kN);
  ExitManager::register_exithandler(cleanup_memory); // on successful exit free allocated data

  Read_ConfigFile();

  if (!recover) // not in recovery-mode -> evaluate arguments
    { 
      // get and evaluate given number
      char* const neuer_str = strdup(argv[1]);
      char* str = neuer_str; 
      if (!parse_term::get_number(n, str))
        {
          cout << "Wrong input at: '" << str << "'" << endl;
          exit(1);
        }
      else
        if (str[0]!='\0')
          {    
            cout << "Syntax error in input term. (parenthesis?)" << endl;
            exit(1);
          }
        else
          if (mpz_cmp_ui(n,0)<=0)
            {
              cout << "Input must be positive natural number!" << endl;
              exit(1);
            }
      free(neuer_str); // don't call "delete []" because of "stdup/malloc"
      FoundFactors.regarding=argv[1];
      mpz_set(FoundFactors.regarding_numeric_value,n);
    }

  if (recover)
    {
      // Recovery: continue a factorization, which was previously aborted
      // (try to) open the necessary file streams
      StaticRelations::FileBuffer.open(StaticRelationsFile.c_str(),ios::in|ios::out);
      SpecialRelations::FileBuffer.open(SpecialRelationsFile.c_str(),ios::in|ios::out);
      if (!DynamicRelations::FileBuffer.open(DynamicRelationsFile.c_str(),ios::in|ios::out))
       {
         cerr << "Cannot access " << DynamicRelationsFile << endl;
         exit(1);
       }
      Recovery_buffer.open(RecoveryFile.c_str(),ios::in|ios::out|ios::ate);

      // offenbar wird trotz ios::ate nicht (immer) ans Ende positioniert
      // deshalb wird die Testabfrage modifiziert:
      if (Recovery_from_file) Recovery_from_file.seekg(0,ios::end);
#ifdef STL_STREAM_workaround
      if ( (!Recovery_from_file) || (Recovery_from_file.tellg()==std::streampos(0)) || (Recovery_from_file.tellg()==std::streampos(-1)) )
// tellg()==0 indicates empty file -> we cannot recover
// tellg()==-1 indicates failure -> we cannot recover
// remark:
//  in gcc 3.4 (cvs-experimental-2003-10-17 we cannot compare with "<" !!)
//  do we really need a workaround to check this condition? 
#else
      if ( (!Recovery_from_file) || (Recovery_from_file.tellg()<1) )
#endif /* STL_STREAM_workaround */
        {
          cerr << "Recovery not possible!" << endl;
          exit(1);
        }
      
      Recovery_from_file.seekg(0,ios::beg);
      Recovery_to_file.seekp(0,ios::beg);
      Recovery_from_file >> n; // retrieve number
      Recovery_from_file.ignore(1,'\n');

      cout << "Continue factorization for " << endl;
      cout << n << endl;

      // mark recovery-mode in Factorization-File! One never knows, what has
      // happened to that file in the mean time!
      Factorization_to_file << " [RECOVERY] "; // flush not necessary... 
      FoundFactors.AutoLoad();
    }
  else
    {
      // First start of this factorization (i.e. not a recovery)
      // open streams
      StaticRelations::FileBuffer.open(StaticRelationsFile.c_str(),ios::out|ios::trunc);
      SpecialRelations::FileBuffer.open(SpecialRelationsFile.c_str(),ios::in|ios::out|ios::trunc);
      if (!DynamicRelations::FileBuffer.open(DynamicRelationsFile.c_str(),ios::in|ios::out|ios::trunc))
       {
         cerr << "Cannot access " << DynamicRelationsFile << endl;
         exit(1);
       }

      Recovery_buffer.open(RecoveryFile.c_str(),ios::in|ios::out|ios::trunc);
      
      cout_status << "Starting factorization for" << endl;
      cout_status << n << endl;
      
      Factorization_to_file << endl << "Factorization of " << argv[1] << ":" << endl;
      Factorization_to_file << n << " = " << endl;
      
      easy_factor(); // remove "easy" detectable factors
      
#ifdef VERBOSE_NOTICE
      cout_status << "Starting factorization with ECM and MPQS for" << endl;
      cout_status << n << endl;
#endif
      Recovery_to_file << n << endl;
      unlink("ecm-recover.dat"); // remove any existent ecm recovery file
    }
  ExitManager::register_exithandler(cleanup_files); // on successful exit delete the files



  // The standalone-version should use elliptic curves, too.
  if (!SkipECM) process_ecm(); // start ECM...

#if defined(USE_DFT) || (CONTINUATION_METHOD==2)
  // at this place, all possible functions that made possible use of
  // discrete fast fourier transform have been finished.
  //  --> we can release the DFT-object and associated auxiliary memory...
  polynomial::clear_dft_tempmemory();
#endif /* USE_DFT */

  // now prepare everything for the Quadratic Sieve...

  determine_best_MPQS_Multiplier(n,kN,MPQS_Multiplier); // calculate a good/optimal value for fastest possible sieving
  tune_parameters(mpz_sizeinbase(n,10)); // tune parameters for n

  
  if ( sqrt(mpz_get_d(kN)) < PhysicalSieveSize )
    {
      cerr << "Sieve size too big (you may want to reduce its size)!" << endl;
      exit(1);
    }

  // Set a threshold for Double-Large-Primes,
  // this is the square of the maximal Single Large Prime...
  mpz_init(CmpqsFactor::DLP_Threshold);
  mpz_set_ui(CmpqsFactor::DLP_Threshold,SingleLargePrime_Threshold);
  mpz_mul(CmpqsFactor::DLP_Threshold,CmpqsFactor::DLP_Threshold,CmpqsFactor::DLP_Threshold);
  
  StaticFactorbase::compute_StaticFactorbase();
  CSieveStaticFactors::initialize();
  SieveControl::compute_SieveThreshold(); // Threshold for detecting relations during sieving phase
  for (int i=0; i<64; ++i) SieveArray_[i]=-1; // initialize array underflow trigger values
  Polynom.compute_first_polynomial(kN,LogicalSieveSize); // compute the first MPQS polynomial to sieve with
  
  if (recover)
    {
      StaticRelations::Load(); // read in Static Relations
      DynamicRelations::Load(); // read in Dynamic Relations
      SpecialRelations::Load();
      statistical_data::ProgressStats.take_sample(); // first statistical sample right here
        // CycleSearch triggers taking a sample when finding a cycle,
        // but this sample wouldn't contain the full DLP size (see CycleSearch)!
        // So we need the first sample taken before the first CycleSearch...
      SpecialRelations::CycleSearch();

#ifdef SAFEMODE
      // are all the relations valid?
      CRelation::SanityCheckRelationsFile(StaticRelationsFile);
      CRelation::SanityCheckRelationsFile(DynamicRelationsFile);
      CRelation::SanityCheckRelationsFile(SpecialRelationsFile);
#endif      
      streampos fpos = Recovery_from_file.tellg();
      Polynom.load_if_available(Recovery_from_file); // get current polynomial (to continue factorization)
      Recovery_to_file.seekp(fpos);
    }


  display_StatusLegend();
  
  // The ultimate loop:
  while(true) // sieve relations as long no factorization has been found.
   {
     // for Recovery:
     streampos fpos = Recovery_to_file.tellp();
     Polynom.save(Recovery_to_file); // store current MPQS polynomial (for recovery)
     Recovery_to_file.seekp(fpos);
     do_sieving();
     Polynom.compute_next_polynomial();
   }

#ifdef VERBOSE_INFO  
  cout << endl << "Session ended successfully." << endl;
#endif
  return 0;
}

---