ginac/polynomial/primes_factory.h

   1 /** @file primes_factory.h
   2  *
   3  *  Factory for prime numbers. */
   4
   5 /*
   6  *  GiNaC Copyright (C) 1999-2022 Johannes Gutenberg University Mainz, Germany
   7  *
   8  *  This program is free software; you can redistribute it and/or modify
   9  *  it under the terms of the GNU General Public License as published by
  10  *  the Free Software Foundation; either version 2 of the License, or
  11  *  (at your option) any later version.
  12  *
  13  *  This program is distributed in the hope that it will be useful,
  14  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
  15  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  16  *  GNU General Public License for more details.
  17  *
  18  *  You should have received a copy of the GNU General Public License
  19  *  along with this program; if not, write to the Free Software
  20  *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
  21  */
  22
  23 #ifndef GINAC_CHINREM_GCD_PRIMES_FACTORY_H
  24 #define GINAC_CHINREM_GCD_PRIMES_FACTORY_H
  25
  26 #include "smod_helpers.h"
  27 #include "debug.h"
  28
  29 #include <cln/integer.h>
  30 #include <cln/numtheory.h>
  31 #include <limits>
  32
  33 namespace GiNaC {
  34
  35 /**
  36  * Find a `big' prime p such that lc mod p != 0. Helper class used by modular
  37  * GCD algorithm.
  38  */
  39 class primes_factory
  40 {
  41 private:
  42         // These primes need to be large enough, so that the number of images
  43         // we need to reconstruct the GCD (in Z) is reasonable. On the other
  44         // hand, they should be as small as possible, so that operations on
  45         // coefficients are efficient. Practically this means we coefficients
  46         // should be native integers. (N.B.: as of now chinrem_gcd uses cl_I
  47         // or even numeric. Eventually this will be fixed).
  48         cln::cl_I last;
  49         // This ensures coefficients are immediate.
  50         static const int immediate_bits = 8*sizeof(void *) - alignof(void *);
  51         static const long opt_hint = (1L << (immediate_bits >> 1)) - 1;
  52 public:
  53         primes_factory()
  54         {
  55                 last = cln::nextprobprime(cln::cl_I(opt_hint));
  56         }
  57
  58         bool operator()(long& p, const cln::cl_I& lc)
  59         {
  60                 static const cln::cl_I maxval(std::numeric_limits<long>::max());
  61                 while (last < maxval) {
  62                         long p_ = cln::cl_I_to_long(last);
  63                         last = cln::nextprobprime(last + 1);
  64
  65                         if (!zerop(smod(lc, p_))) {
  66                                 p = p_;
  67                                 return true;
  68                         }
  69                 }
  70                 return false;
  71         }
  72
  73         bool has_primes() const
  74         {
  75                 static const cln::cl_I maxval(std::numeric_limits<long>::max());
  76                 return last < maxval;
  77         }
  78 };
  79
  80 } // namespace GiNaC
  81
  82 #endif // ndef GINAC_CHINREM_GCD_PRIMES_FACTORY_H