check_polynomial.cc

Go to the documentation of this file.

// only for testing the functions

extern "C"
{
  #include <sys/times.h>
}

#define USER_GMP_WRAP
#include "modulo.cc"
#include "polynomial.cc"

// to check the stuff do:  g++ -O2 -g -Wall check_polynomial.cc -lgmp
// and for  dft-check do:  g++ -O2 -g -Wall -DUSE_DFT check_polynomial.cc -lgmp

using namespace polynomial;


double difftimes(clock_t t1, clock_t t0)
{
  // refer "man 2 times"
  return static_cast<double>(t1-t0)/static_cast<double>(sysconf(_SC_CLK_TCK));
}


bool do_check(const int maxk, const int anzPoints, const mpz_t m)
{
  cout << endl << endl << "====================================" << endl;
  cout << " MULTIPOINT POLYNOM-EVAL SANITY " << maxk << ", " << anzPoints << endl; 

  int i;

  mpz_t Polynom [maxk];
  mpz_t Points [anzPoints];

  {
   const size_t bits = mpz_sizeinbase(m,2);

   gmp_randstate_t randstate;
   gmp_randinit_default(randstate);

   for (i=0; i<maxk; ++i) { mpz_init(Polynom[i]); mpz_rrandomb(Polynom[i],randstate,bits); mpz_mod(Polynom[i],Polynom[i],m); }
   for (i=0; i<anzPoints; ++i) { mpz_init(Points[i]); mpz_rrandomb(Points[i],randstate,bits); mpz_mod(Points[i],Points[i],m); }

   gmp_randclear(randstate);
  }

  mpz_t Res_Horner [anzPoints];
  for (i=0; i<anzPoints; ++i) mpz_init(Res_Horner[i]);
 
  mpz_t Res_Multipoint [anzPoints];
  for (i=0; i<anzPoints; ++i) mpz_init(Res_Multipoint[i]);

  // evaluate using Horner scheme
  cout << "Starting evaluation using Horner scheme" << endl;
  clock_t Start = times(NULL);
  for (i=0; i<anzPoints; ++i) eval(Res_Horner[i],Polynom,maxk,Points[i],m);
  clock_t Stop = times(NULL);
  cout << difftimes(Stop,Start) << " seconds (user time)." << endl;

  // now fast multipoint evaluation
  cout << "Starting fast multipoint evaluation" << endl;
  Start = times(NULL);
  multipoint_eval(Res_Multipoint,Polynom,maxk,Points,anzPoints,m);
  Stop = times(NULL);
  cout << difftimes(Stop,Start) << " seconds (user time)." << endl;

  cout << "verifying points for Horner <-> Multipoint..." << endl;
  bool error = false;
  for (int i=0; i<anzPoints; ++i)
   {
#ifdef DEBUG
     mpz_out_str(stdout,10,Points[i]); cout << " -> ";
     mpz_out_str(stdout,10,Res_Horner[i]); cout << " (Horner)  ";
     mpz_out_str(stdout,10,Res_Multipoint[i]); cout << " (Multipoint)" << endl;
#endif
     if (mpz_cmp(Res_Horner[i],Res_Multipoint[i]))
       {
         error=true;
         cerr << "values differ!!" << i << endl;
         mpz_out_str(stdout,10,Points[i]); cout << " -> ";
         mpz_out_str(stdout,10,Res_Horner[i]); cout << " (Horner)  ";
         mpz_out_str(stdout,10,Res_Multipoint[i]); cout << " (Multipoint)" << endl;
         break;
       }
   }
  if (error)
   {
     cout << "FAILED!" << endl;
     exit(1);
   }

  cout << "conforming; test passed." << endl;
  cout << "PASSED." << endl;

  for (i=0; i<maxk; ++i) { mpz_clear(Polynom[i]); }
  for (i=0; i<anzPoints; ++i) mpz_clear(Points[i]);
  for (i=0; i<anzPoints; ++i) mpz_clear(Res_Horner[i]);
  for (i=0; i<anzPoints; ++i) mpz_clear(Res_Multipoint[i]);
  return true;
}

bool performance_check(const int maxk, const int anzPoints, const mpz_t m)
{
  cout << endl << endl << "====================================" << endl;
  cout << " MULTIPOINT POLYNOM-EVAL PERFORMANCE " << maxk << ", " << anzPoints << endl; 

  int i;

  mpz_t Polynom [maxk];
  mpz_t Points [anzPoints];

  {
   const size_t bits = mpz_sizeinbase(m,2);

   gmp_randstate_t randstate;
   gmp_randinit_default(randstate);

   for (i=0; i<maxk; ++i) { mpz_init(Polynom[i]); mpz_rrandomb(Polynom[i],randstate,bits); mpz_mod(Polynom[i],Polynom[i],m); }
   for (i=0; i<anzPoints; ++i) { mpz_init(Points[i]); mpz_rrandomb(Points[i],randstate,bits); mpz_mod(Points[i],Points[i],m); }

   gmp_randclear(randstate);
  }


  mpz_t Res_Multipoint [anzPoints];
  for (i=0; i<anzPoints; ++i) mpz_init(Res_Multipoint[i]);

  // now fast multipoint evaluation
  cout << "Starting fast multipoint evaluation" << endl;
  clock_t Start = times(NULL);
  multipoint_eval(Res_Multipoint,Polynom,maxk,Points,anzPoints,m);
  clock_t Stop = times(NULL);
  cout << difftimes(Stop,Start) << " seconds (user time)." << endl;

  for (i=0; i<maxk; ++i) { mpz_clear(Polynom[i]); }
  for (i=0; i<anzPoints; ++i) mpz_clear(Points[i]);
  for (i=0; i<anzPoints; ++i) mpz_clear(Res_Multipoint[i]);
  return true;
}


int main()
{
  cout << endl << endl
       << "===============================================" << endl;
  cout << "POLYNOMIAL-MULTIPOINT-EVAL doing some checks..." << endl;
  cout << "===============================================" << endl;

  mpz_t m;
  mpz_init(m);

  // generate a 120-digit number
  mpz_set_ui(m,10); mpz_pow_ui(m,m,120); mpz_add_ui(m,m,3);
//    mpz_set_ui(m,10); mpz_pow_ui(m,m,250); mpz_add_ui(m,m,3);

#if 0
  // torture check...
  for (int maxk=4; maxk<(1<<14); maxk*=2)
   {
     for (int i=3; i<maxk; ++i) do_check(maxk,i,m);
     for (int i=3; i<maxk+1; ++i) do_check(maxk+1,i,m);
   }
#endif


#if 1
  for (int maxk=1024; maxk<(1<<19); maxk*=2)
   {
     performance_check(maxk+1,maxk,m);
   }
#endif

  //performance_check(32769,32768,m);

#if 0
  // now a normal check
  for (int maxk=4; maxk<(1<<16); maxk*=2)
   {
     do_check(maxk,maxk-1,m);
     do_check(maxk+1,maxk,m);
   }

  mpz_set_str(m,"3923385745693995079670229419275984584311007321932374190635656246740175165573932140787529348954892963218868359081838772941945556717",10);
  for (int maxk=4; maxk<2048; maxk*=2)
   {
     do_check(maxk,maxk-1,m);
     do_check(maxk+1,maxk,m);
   }
#endif

  mpz_clear(m);
}

---