]> www.ginac.de Git - cln.git/blob - src/modinteger/cl_MI.cc
a3d55116e7aecc6f0dd7a0c1e32a7c8845100c29
[cln.git] / src / modinteger / cl_MI.cc
1 // Modular integer operations.
2
3 // General includes.
4 #include "cl_sysdep.h"
5
6 CL_PROVIDE(cl_MI)
7
8 // Specification.
9 #include "cln/modinteger.h"
10
11
12 // Implementation.
13
14 #include "cl_I.h"
15 #include "cl_DS.h"
16 #include "cl_2DS.h"
17 #include "cln/io.h"
18 #include "cln/integer_io.h"
19 #include "cl_N.h"
20 #include "cl_MI.h"
21 #include "cln/exception.h"
22 #include "cl_alloca.h"
23
24 // MacOS X does "#define _R 0x00040000L"
25 // Grr...
26 #undef _R
27
28 namespace cln {
29
30 static void cl_modint_ring_destructor (cl_heap* pointer)
31 {
32         (*(cl_heap_modint_ring*)pointer).~cl_heap_modint_ring();
33 }
34
35 cl_class cl_class_modint_ring = {
36         cl_modint_ring_destructor,
37         cl_class_flags_modint_ring
38 };
39
40 cl_heap_modint_ring::cl_heap_modint_ring (cl_I m, cl_modint_setops* setopv, cl_modint_addops* addopv, cl_modint_mulops* mulopv)
41         : setops (setopv), addops (addopv), mulops (mulopv), modulus (m)
42 {
43         refcount = 0; // will be incremented by the `cl_modint_ring' constructor
44         type = &cl_class_modint_ring;
45         if (minusp(m)) throw runtime_exception();
46         if (!cln::zerop(m)) {
47                 var uintC b = integer_length(m-1);
48                 // m <= 2^b, hence one needs b bits for a representative mod m.
49                 if (b <= 1) {
50                         log2_bits = 0; bits = 1;
51                 } else if (b <= cl_word_size) {
52                         var uintL bb;
53                         integerlengthC(b-1,bb=); // b <= 2^bb with bb minimal
54                         log2_bits = bb; bits = 1<<bb;
55                 } else {
56                         log2_bits = -1; bits = -1;
57                 }
58         } else {
59                 log2_bits = -1; bits = -1;
60         }
61 }
62
63
64 static bool modint_equal (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
65 {
66         unused R;
67         return equal(x.rep,y.rep);
68 }
69
70 }  // namespace cln
71 #include "cl_MI_int.h"
72 #include "cl_MI_std.h"
73 #include "cl_MI_fix16.h"
74 #if (cl_value_len <= 32)
75 #include "cl_MI_fix29.h"
76 #include "cl_MI_int32.h"
77 #else
78 #include "cl_MI_fix32.h"
79 #endif
80 #include "cl_MI_pow2.h"
81 #include "cl_MI_pow2m1.h"
82 #include "cl_MI_pow2p1.h"
83 #include "cl_MI_montgom.h"
84 namespace cln {
85
86
87 static inline cl_heap_modint_ring* make_modint_ring (const cl_I& m) // m >= 0
88 {
89         if (m == 0)
90                 return new cl_heap_modint_ring_int();
91         // Now m > 0.
92         {
93                 var uintC log2_m = power2p(m);
94                 if (log2_m)
95                         return new cl_heap_modint_ring_pow2(m,log2_m-1);
96         }
97         // Now m > 1.
98         {
99                 var uintC log2_m = integer_length(m); // = integerlength(m-1)
100                 if (log2_m < 16) // m < 0x10000
101                         return new cl_heap_modint_ring_fix16(m);
102                 #if (cl_value_len <= 32)
103                 if (log2_m < cl_value_len)
104                         return new cl_heap_modint_ring_fix29(m);
105                 if (log2_m < 32) // m < 0x100000000
106                         return new cl_heap_modint_ring_int32(m);
107                 #else
108                 if (log2_m < 32) // m < 0x100000000
109                         return new cl_heap_modint_ring_fix32(m);
110                 #endif
111         }
112         {
113                 var uintC m1 = power2p(m+1);
114                 if (m1)
115                         return new cl_heap_modint_ring_pow2m1(m,m1-1);
116         }
117         {
118                 var uintC m1 = power2p(m-1);
119                 if (m1)
120                         return new cl_heap_modint_ring_pow2p1(m,m1-1);
121         }
122         {
123                 var cl_heap_modint_ring* R = try_make_modint_ring_montgom(m);
124                 if (R)
125                         return R;
126         }
127         return new cl_heap_modint_ring_std(m);
128 }
129
130
131 // The table of modular integer rings.
132 // A weak hash table cl_I -> cl_modint_ring.
133 // (It could also be a weak hashuniq table cl_I -> cl_modint_ring.)
134
135 }  // namespace cln
136 #include "cl_I_hashweak_rcpointer.h"
137 namespace cln {
138
139 // An entry can be collected when the value (the ring) isn't referenced any more
140 // except from the hash table, and when the key (the modulus) isn't referenced
141 // any more except from the hash table and the ring. Note that the ring contains
142 // exactly one reference to the modulus.
143
144 static bool maygc_htentry (const cl_htentry_from_integer_to_rcpointer& entry)
145 {
146         if (!entry.key.pointer_p() || (entry.key.heappointer->refcount == 2))
147                 if (!entry.val.pointer_p() || (entry.val.heappointer->refcount == 1))
148                         return true;
149         return false;
150 }
151
152 static const cl_wht_from_integer_to_rcpointer modint_ring_table = cl_wht_from_integer_to_rcpointer(maygc_htentry);
153
154 static inline cl_modint_ring* get_modint_ring (const cl_I& m)
155 {
156         return (cl_modint_ring*) modint_ring_table.get(m);
157 }
158
159 static inline void store_modint_ring (const cl_modint_ring& R)
160 {
161         modint_ring_table.put(R->modulus,R);
162 }
163
164
165 const cl_modint_ring find_modint_ring (const cl_I& m)
166 {
167  {      Mutable(cl_I,m);
168         m = abs(m);
169         var cl_modint_ring* ring_in_table = get_modint_ring(m);
170         if (!ring_in_table) {
171                 var cl_modint_ring R = make_modint_ring(m);
172                 store_modint_ring(R);
173                 ring_in_table = get_modint_ring(m);
174                 if (!ring_in_table)
175                         throw runtime_exception();
176         }
177         return *ring_in_table;
178 }}
179
180
181 const cl_modint_ring cl_modint0_ring = find_modint_ring(0);
182
183 }  // namespace cln
184
185 CL_PROVIDE_END(cl_MI)