]> www.ginac.de Git - cln.git/blobdiff - src/modinteger/cl_MI_std.h
Convert complex numbers to real numbers if imaginary part is floating-point 0.0.
[cln.git] / src / modinteger / cl_MI_std.h
index 4cd56e07ed86be282d67fccda21264171c974853..ede3818d5918cb57e54cba87beaf726ff15daf62 100644 (file)
@@ -1,6 +1,8 @@
 // m > 1, standard representation, no tricks
 
-static void std_fprint (cl_heap_modint_ring* R, cl_ostream stream, const _cl_MI &x)
+namespace cln {
+
+static void std_fprint (cl_heap_modint_ring* R, std::ostream& stream, const _cl_MI &x)
 {
        fprint(stream,R->_retract(x));
        fprint(stream," mod ");
@@ -23,7 +25,7 @@ static const cl_I std_retract (cl_heap_modint_ring* R, const _cl_MI& x)
        return x.rep;
 }
 
-static const _cl_MI std_random (cl_heap_modint_ring* R, cl_random_state& randomstate)
+static const _cl_MI std_random (cl_heap_modint_ring* R, random_state& randomstate)
 {
        return _cl_MI(R, random_I(randomstate,R->modulus));
 }
@@ -33,7 +35,7 @@ static const _cl_MI std_zero (cl_heap_modint_ring* R)
        return _cl_MI(R, 0);
 }
 
-static cl_boolean std_zerop (cl_heap_modint_ring* R, const _cl_MI& x)
+static bool std_zerop (cl_heap_modint_ring* R, const _cl_MI& x)
 {
        unused R;
        return zerop(x.rep);
@@ -80,7 +82,7 @@ static const cl_MI_x std_recip (cl_heap_modint_ring* R, const _cl_MI& x)
        if (eq(g,1))
                return cl_MI(R, (minusp(u) ? u + R->modulus : u));
        if (zerop(xr))
-               cl_error_division_by_0();
+               throw division_by_0_exception();
        return cl_notify_composite(R,xr);
 }
 
@@ -93,7 +95,7 @@ static const cl_MI_x std_div (cl_heap_modint_ring* R, const _cl_MI& x, const _cl
        if (eq(g,1))
                return cl_MI(R, mod(x.rep * (minusp(u) ? u + R->modulus : u), R->modulus));
        if (zerop(yr))
-               cl_error_division_by_0();
+               throw division_by_0_exception();
        return cl_notify_composite(R,yr);
 }
 
@@ -169,11 +171,11 @@ static const _cl_MI std_expt_pos (cl_heap_modint_ring* R, const _cl_MI& x, const
        //   k = 3 for nn <= 69.8
        //   ... k for nn <= k*(k+1)*2^(2k)/(2^(k+1)-k-2)
        var cl_I n = y;
-       var uintL nn = integer_length(n);
+       var uintC nn = integer_length(n);
        // n has nn bits.
        if (nn <= 8) {
                // k = 1, normal Left-Right Binary algorithm.
-               var uintL _n = FN_to_UL(n);
+               var uintL _n = FN_to_UV(n);
                var _cl_MI a = x;
                for (var int i = nn-2; i >= 0; i--) {
                        a = R->_square(a);
@@ -200,17 +202,17 @@ static const _cl_MI std_expt_pos (cl_heap_modint_ring* R, const _cl_MI& x, const
                else if (nn <= 1721160) k = 14;
                else if (nn <= 3933180) k = 15;
                else /* if (nn <= 8914120) */ k = 16;
-               var uintL nnk = ceiling(nn,k); // number of base-2^k digits in n
+               var uintC nnk = ceiling(nn,k); // number of base-2^k digits in n
                var uint16* n_digits = cl_alloc_array(uint16,nnk);
                // Split n into base-2^k digits.
                {
                        var const uintD* n_LSDptr;
                        var const uintD* n_MSDptr;
-                       I_to_NDS_nocopy(n, n_MSDptr=,,n_LSDptr=,cl_false,);
+                       I_to_NDS_nocopy(n, n_MSDptr=,,n_LSDptr=,false,);
                        var const uintL k_mask = bit(k)-1;
                        var uintD carry = 0;
                        var unsigned int carrybits = 0;
-                       for (var uintL i = 0; i < nnk; i++) {
+                       for (var uintC i = 0; i < nnk; i++) {
                                if (carrybits >= k) {
                                        n_digits[i] = carry & k_mask;
                                        carry = carry >> k;
@@ -227,7 +229,7 @@ static const _cl_MI std_expt_pos (cl_heap_modint_ring* R, const _cl_MI& x, const
                // Compute maximum odd base-2^k digit.
                var uintL maxodd = 1;
                if (k <= 8) {
-                       for (var uintL i = 0; i < nnk; i++) {
+                       for (var uintC i = 0; i < nnk; i++) {
                                var uintL d = n_digits[i];
                                if (d > 0) {
                                        d = odd_table[d];
@@ -235,7 +237,7 @@ static const _cl_MI std_expt_pos (cl_heap_modint_ring* R, const _cl_MI& x, const
                                }
                        }
                } else {
-                       for (var uintL i = 0; i < nnk; i++) {
+                       for (var uintC i = 0; i < nnk; i++) {
                                var uintL d = n_digits[i];
                                if (d > 0) {
                                        var uintL d2; ord2_32(d,d2=);
@@ -256,7 +258,7 @@ static const _cl_MI std_expt_pos (cl_heap_modint_ring* R, const _cl_MI& x, const
                // Compute a = x^n_digits[nnk-1].
                {
                        var uintL d = n_digits[nnk-1];
-                       if (d == 0) cl_abort();
+                       if (d == 0) throw runtime_exception();
                        var uintL d2;
                        if (k <= 8)
                                d2 = ord2_table[d];
@@ -269,11 +271,11 @@ static const _cl_MI std_expt_pos (cl_heap_modint_ring* R, const _cl_MI& x, const
                        if (d==0 && maxodd > 1 && d2>0) {
                                a = x2; d2--;
                        }
-                       if (!(d2 < k)) cl_abort();
+                       if (!(d2 < k)) throw runtime_exception();
                        for ( ; d2>0; d2--)
                                a = R->_square(a);
                }
-               for (var sintL i = nnk-2; i >= 0; i--) {
+               for (var sintC i = nnk-2; i >= 0; i--) {
                        // Compute a := a^(2^k) * x^n_digits[i].
                        var uintL d = n_digits[i];
                        var uintL d2;
@@ -290,7 +292,7 @@ static const _cl_MI std_expt_pos (cl_heap_modint_ring* R, const _cl_MI& x, const
                        } else
                                d2 = k;
                        // Square d2 times.
-                       if (!(d2 <= k)) cl_abort();
+                       if (!(d2 <= k)) throw runtime_exception();
                        for ( ; d2>0; d2--)
                                a = R->_square(a);
                }
@@ -343,8 +345,26 @@ class cl_heap_modint_ring_std : public cl_heap_modint_ring {
        SUBCLASS_cl_heap_modint_ring()
 public:
        // Constructor.
-       cl_heap_modint_ring_std (const cl_I& m)
-               : cl_heap_modint_ring (m, &std_setops, &std_addops, &std_mulops) {}
+       cl_heap_modint_ring_std (const cl_I& m);
        // Virtual destructor.
        ~cl_heap_modint_ring_std () {}
 };
+
+static void cl_heap_modint_ring_std_destructor (cl_heap* pointer)
+{
+       (*(cl_heap_modint_ring_std*)pointer).~cl_heap_modint_ring_std();
+}
+
+cl_class cl_class_modint_ring_std = {
+       cl_heap_modint_ring_std_destructor,
+       cl_class_flags_modint_ring
+};
+
+// Constructor.
+inline cl_heap_modint_ring_std::cl_heap_modint_ring_std (const cl_I& m)
+       : cl_heap_modint_ring (m, &std_setops, &std_addops, &std_mulops)
+{
+       type = &cl_class_modint_ring_std;
+}
+
+}  // namespace cln