X-Git-Url: https://ginac.de/CLN/cln.git//cln.git?a=blobdiff_plain;f=src%2Fnumtheory%2Fcl_nt_sqrtmodp.cc;h=8c0741e9c7524cb7c6b7316a79bcf887b431e987;hb=8b3d91dec77438c0fe679b10869ab29e6cdeba58;hp=6a74c5f26805af6bcb4dd89f11004cf453120816;hpb=850abfde7f0d985ba01526c346bcd0d733562943;p=cln.git diff --git a/src/numtheory/cl_nt_sqrtmodp.cc b/src/numtheory/cl_nt_sqrtmodp.cc index 6a74c5f..8c0741e 100644 --- a/src/numtheory/cl_nt_sqrtmodp.cc +++ b/src/numtheory/cl_nt_sqrtmodp.cc @@ -10,12 +10,15 @@ // Implementation. #include "cl_I.h" -#include "cln/abort.h" +#include "cln/exception.h" #undef floor #include #define floor cln_floor +// MacOS X does "#define _R 0x00040000L". Grr... +#undef _R + namespace cln { // Algorithm 1 (for very small p only): @@ -173,6 +176,16 @@ static const sqrt_mod_p_t cantor_zassenhaus_sqrt (const cl_modint_ring& R, const } } +#if defined(__GNUC__) && defined(__s390__) && (__GNUC__ == 2) // Workaround GCC-bug (see below) + struct cl_sylow2gen_property : public cl_property { + SUBCLASS_cl_property(); + public: + cl_I h_rep; + // Constructor. + cl_sylow2gen_property (const cl_symbol& k, const cl_MI& h) : cl_property (k), h_rep (h.rep) {} + }; +#endif + // Algorithm 3 (for p > 2 only): // Tonelli-Shanks. // [Cohen, A Course in Computational Algebraic Number Theory, @@ -196,11 +209,11 @@ static const sqrt_mod_p_t tonelli_shanks_sqrt (const cl_modint_ring& R, const cl // take h in G_(j+1) \ G_j, so that h^2 in G_j \ G_(j-1), and // a^-1*b^2*h^2 is in G_(j-1). So multiply b with h. var cl_I& p = R->modulus; - var uintL e = ord2(p-1); + var uintC e = ord2(p-1); var cl_I m = (p-1) >> e; // p-1 = 2^e*m, m odd. // We will have the invariant c = a^-1*b^2 in G/G_j. - var uintL j = e; + var uintC j = e; // Initialize b = a^((m+1)/2), c = a^m, but avoid to divide by a. var cl_MI c = R->expt_pos(a,(m-1)>>1); var cl_MI b = R->mul(a,c); @@ -213,6 +226,7 @@ static const sqrt_mod_p_t tonelli_shanks_sqrt (const cl_modint_ring& R, const cl // Since this computation is a bit costly, we cache its result // on the ring's property list. static const cl_symbol key = (cl_symbol)(cl_string)"generator of 2-Sylow subgroup of (Z/pZ)^*"; +#if !(defined(__GNUC__) && defined(__s390__) && (__GNUC__ == 2)) // Workaround GCC-bug (see above) struct cl_sylow2gen_property : public cl_property { SUBCLASS_cl_property(); public: @@ -220,6 +234,7 @@ static const sqrt_mod_p_t tonelli_shanks_sqrt (const cl_modint_ring& R, const cl // Constructor. cl_sylow2gen_property (const cl_symbol& k, const cl_MI& h) : cl_property (k), h_rep (h.rep) {} }; +#endif var cl_sylow2gen_property* prop = (cl_sylow2gen_property*) R->get_property(key); if (prop) h = cl_MI(R,prop->h_rep); @@ -233,7 +248,7 @@ static const sqrt_mod_p_t tonelli_shanks_sqrt (const cl_modint_ring& R, const cl do { // Now c = a^-1*b^2 in G_j, h in G_j \ G_(j-1). // Determine the smallest i such that c in G_i. - var uintL i = 0; + var uintC i = 0; var cl_MI ci = c; // c_i = c^(2^i) for ( ; i < j; i++, ci = R->square(ci)) if (ci == R->one()) @@ -244,7 +259,7 @@ static const sqrt_mod_p_t tonelli_shanks_sqrt (const cl_modint_ring& R, const cl // Indicates that p is not prime. return new cl_composite_condition(p); // OK, i < j. - for (var uintL count = j-i-1; count > 0; count--) + for (var uintC count = j-i-1; count > 0; count--) h = R->square(h); // Now h in G_(i+1) \ G_i. b = R->mul(b,h); @@ -272,7 +287,7 @@ static const sqrt_mod_p_t tonelli_shanks_sqrt (const cl_modint_ring& R, const cl const sqrt_mod_p_t sqrt_mod_p (const cl_modint_ring& R, const cl_MI& a) { - if (!(a.ring() == R)) cl_abort(); + if (!(a.ring() == R)) throw runtime_exception(); var cl_I& p = R->modulus; var cl_I aa = R->retract(a); switch (jacobi(aa,p)) { @@ -297,9 +312,9 @@ const sqrt_mod_p_t sqrt_mod_p (const cl_modint_ring& R, const cl_MI& a) else return sqrt_mod_p_t(2,R->canonhom(x1),R->canonhom(x2)); } - var uintL l = integer_length(p); - var uintL e = ord2(p-1); - //if (e > 30 && e > l/(log((double)l)*0.72-1)) + var uintC l = integer_length(p); + var uintC e = ord2(p-1); + //if (e > 30 && e > l/(::log((double)l)*0.72-1)) if (e > 30 && e > l/(::log((double)l)*0.92-2.41)) // Algorithm 2. return cantor_zassenhaus_sqrt(R,a);