modular_mult.cc

Go to the documentation of this file.

// Modular Multiplication Without Trial Division
/* reference:
   Peter L. Montgomery: "Modular Multiplication Without Trial Division",
   Mathematics Of Computation, Vol. 44 (170), April 1985, p. 519-521
*/

// this module written by Thorsten Reinecke, 1999-11-16
// last change: 2004-08-01

// not efficient enough (-> should use mpn instead of mpz)

class CN_Residue
{
 private:
  mpz_t N;
  mpz_t R,Ri,Ns;
  mpz_t m,t;
  unsigned int k;
 public:
  inline CN_Residue() /* Konstruktor */
    {
      cout << "Modular Multiplication Without Trial Division activated." << endl;
      mpz_init(N);
      mpz_init(R); mpz_init(Ri); mpz_init(Ns);
      mpz_init(m); mpz_init(t);
      k=0; // Dummywert-Initialisierung
    };
  inline void init(const mpz_t M) // Initialization
    {
      mpz_set(N,M); // NResidue --> mod N
      // select R=2^k, such that x div/mod R can be computed fast
      // compute Ri:=R^(-1) mod N; Ns:=(R*Ri-1)/N;
      k=mpz_sizeinbase(N,2)+1; mpz_set_ui(R,1); mpz_mul_2exp(R,R,k);
      mpz_invert(Ri,R,N); mpz_mod(Ri,Ri,N);
      mpz_mul(Ns,R,Ri); mpz_sub_ui(Ns,Ns,1); mpz_div(Ns,Ns,N);
    };
  inline explicit CN_Residue(const mpz_t M) // constructor
    { 
      CN_Residue();
      init(M);
    };
  inline ~CN_Residue() // destructor
    {
      mpz_clear(N);
      mpz_clear(R); mpz_clear(Ri); mpz_clear(Ns);
      mpz_clear(m); mpz_clear(t);
    };
  inline void redc(mpz_t res, const mpz_t T) 
    {
      //mpz_mod(m,T,R); // isn't efficient!
      mpz_tdiv_r_2exp(m,T,k); // this should be better...
      mpz_mul(m,m,Ns);
      //mpz_mod(m,m,R); // isn't efficient!
      mpz_tdiv_r_2exp(m,m,k); // this should be better...
      mpz_mul(t,m,N); mpz_add(t,t,T);
      //mpz_div(t,t,R); // isn't efficient!
      mpz_tdiv_q_2exp(res,t,k); // this should be better...
      if (mpz_cmp(res,N)>=0)
       {
         mpz_sub(res,res,N);
         if (mpz_cmp(res,N)>=0)
          {
            // overflow! arguments are not normalized!
            //MARK;
            mpz_mod(res,res,N); // fallback: standard mod
          }
       }
    };
  inline void convert(mpz_t x) const 
   {
     //mpz_mul(x,x,R); // isn't efficient!
     mpz_mul_2exp(x,x,k); // this should be better...
     mpz_mod(x,x,N);
   };
  inline void convert_back(mpz_t x) const 
   {
     mpz_mul(x,x,Ri); mpz_mod(x,x,N);
   };
  inline int invert(mpz_t res, const mpz_t T) const 
   {
     mpz_set(res,T);
     convert_back(res);
     int i=mpz_invert(res,res,N);
     convert(res);
     return i;
   };
};

---