00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #ifndef __MMX__MODULAR__HPP
00014 #define __MMX__MODULAR__HPP
00015 #include <basix/generic.hpp>
00016 #include <basix/wrap.hpp>
00017 #include <basix/operators.hpp>
00018 #include <basix/int.hpp>
00019 #include <numerix/integer.hpp>
00020 #include <numerix/modulus.hpp>
00021
00022 namespace mmx {
00023 #define Modular_variant(M) modular_variant_helper<M>::MV
00024 #define TMPL_DEF \
00025 template<typename M, typename V= typename Modular_variant(M)>
00026 #define TMPL template<typename M, typename V>
00027 #define Modulus modulus<C,V>
00028 #define Modular modular<M,V>
00029
00030 template<typename M>
00031 struct modular_variant_helper;
00032 TMPL class modular;
00033
00034
00035
00036
00037
00038
00039 struct modular_local {
00040 template<typename M>
00041 class modulus_storage {
00042 static inline M& dyn_modulus () {
00043 static M modulus = M ();
00044 return modulus; }
00045 public:
00046 static inline void set_modulus (const M& p) {
00047 if (! is_exact_zero (*p)) dyn_modulus () = p; }
00048 static inline M get_modulus () { return dyn_modulus (); }
00049 };
00050 };
00051
00052
00053 struct modular_global {
00054 template<typename M>
00055 class modulus_storage {
00056 static inline M& dyn_modulus () {
00057 static M modulus = M ();
00058 return modulus; }
00059 public:
00060 static inline void set_modulus (const M& p) { dyn_modulus () = p; }
00061 static inline M get_modulus () { return dyn_modulus (); }
00062 };
00063 };
00064
00065
00066 template<typename C, C p>
00067 struct modular_fixed {
00068 template<typename M>
00069 class modulus_storage {
00070 static inline M& dyn_modulus () {
00071 static fixed_value <C, p> _p;
00072 static M modulus = M (_p);
00073 return modulus; }
00074 public:
00075 static inline void set_modulus (const M& q) {
00076 ASSERT (q == dyn_modulus (),
00077 "modular_fixed: can not change the modulus");
00078 }
00079 static inline M get_modulus () { return dyn_modulus (); }
00080 };
00081 };
00082
00083
00084
00085
00086
00087 template<typename M>
00088 struct modular_variant_helper {
00089 typedef modular_global MV; };
00090
00091
00092
00093
00094
00095 #define DEFINE_MODULAR_GLOBAL(name) \
00096 struct modular_global_##name { \
00097 template<typename P> \
00098 class modulus_storage { \
00099 static inline P& dyn_modulus () { \
00100 static P modulus = P (); \
00101 return modulus; } \
00102 public: \
00103 static inline void set_modulus (const P& p) { dyn_modulus () = p; } \
00104 static inline P get_modulus () { return dyn_modulus (); } \
00105 }; \
00106 };
00107
00108 #define Modular_global(name) modular_global_##name
00109
00110
00111
00112
00113
00114 TMPL_DEF
00115 class modular {
00116 MMX_ALLOCATORS
00117 public:
00118 typedef M modulus;
00119 typedef typename M::base C;
00120 typedef typename V::template modulus_storage<M> S;
00121
00122 C rep;
00123 static inline M get_modulus () { return S::get_modulus (); }
00124 static inline void set_modulus (const M& p) { S::set_modulus (p); }
00125
00126 inline C operator * () const {
00127 C dest;
00128 decode_mod (dest, rep, get_modulus ());
00129 return dest; }
00130
00131 inline modular () { encode_mod (rep, C (), get_modulus ()); }
00132
00133 inline modular (const Modular& s) : rep (s.rep) {}
00134
00135 inline modular (const C& a, bool reduced) {
00136 if (reduced) rep = a;
00137 else encode_mod (rep, a, get_modulus ()); }
00138
00139 inline modular (const C& s, const M& p, bool reduced=false) {
00140 VERIFY (p == get_modulus (), "wrong modulus");
00141 if (reduced) rep= s;
00142 else encode_mod (rep, s, p); }
00143
00144 template<typename O>
00145 inline modular (const O& a) {
00146 encode_mod (rep, as<C> (a), get_modulus ()); }
00147
00148 inline Modular& operator = (const Modular& a) {
00149 rep = a.rep;
00150 return *this; }
00151
00152 inline bool operator == (const Modular& a) const {
00153 return rep == a.rep; }
00154
00155 inline bool operator != (const Modular& a) const {
00156 return rep != a.rep; }
00157 };
00158
00159 TMPL inline M get_modulus (const Modular& x) {
00160 (void) x; return Modular::get_modulus (); }
00161 TMPL inline void set_modulus (Modular& x, const M& p) {
00162 (void) x; Modular::set_modulus (p); }
00163
00164 template<typename TM, typename TV, typename SM, typename SV>
00165 struct as_helper<modular<TM,TV>,modular<SM,SV> > {
00166 static inline modular<TM,TV>
00167 cv (const modular<SM,SV>& p) { return modular<TM,TV> (*p); }
00168 };
00169
00170 template<typename C, typename V, typename W>
00171 struct as_helper<C,modular<modulus<C,V>,W> > {
00172 static inline C
00173 cv (const modular<modulus<C,V>,W>& p) { return *p; }
00174 };
00175
00176 template<typename C, typename V, typename W>
00177 struct as_helper<generic, modular<modulus<C,V>,W> > {
00178 static inline generic cv (const modular<modulus<C,V>,W> & x) {
00179 return new generic_concrete_rep<modular<modulus<C,V>,W> > (x); }
00180 };
00181
00182 WRAP_WRAPPED_IMPL(TMPL inline,Modular)
00183
00184
00185
00186
00187
00188 UNARY_RETURN_TYPE(TMPL,lift_op,Modular,typename Modular::C);
00189
00190 TMPL
00191 struct lift_helper<Modular> {
00192 template<typename R> static inline void
00193 set_op (R& y, const Modular& x) { set_as (y, *x); }
00194 static inline Lift_type(Modular)
00195 op (const Modular& x) { return *x; }
00196 };
00197
00198
00199
00200
00201
00202 UNARY_RETURN_TYPE(template<typename C>,project_op,C,modular<modulus<C> >);
00203
00204 template<typename C> template<typename R> void
00205 project_helper<C>::set_op (R& y, const C& x) { y= R (x); }
00206
00207 template<typename C> Project_type(C)
00208 project_helper<C>::op (const C& x) { return Project_type(C) (x); }
00209
00210
00211
00212
00213
00214 TMPL inline Modular&
00215 operator += (Modular& dest, const Modular& s) {
00216 add_mod (dest.rep, s.rep, Modular::get_modulus ());
00217 return dest; }
00218
00219 TMPL inline Modular&
00220 operator -= (Modular& dest, const Modular& s) {
00221 sub_mod (dest.rep, s.rep, Modular::get_modulus ());
00222 return dest; }
00223
00224 TMPL inline Modular&
00225 operator *= (Modular& dest, const Modular& s) {
00226 mul_mod (dest.rep, s.rep, Modular::get_modulus ());
00227 return dest; }
00228
00229 TMPL inline Modular&
00230 operator /= (Modular& dest, const Modular& s) {
00231 div_mod (dest.rep, s.rep, Modular::get_modulus ());
00232 return dest; }
00233
00234 TMPL inline Modular
00235 operator + (const Modular& a, const Modular& b) {
00236 Modular c;
00237 add_mod (c.rep, a.rep, b.rep, Modular::get_modulus ());
00238 return c; }
00239
00240 TMPL inline Modular
00241 operator + (const Modular& a, const typename M::base& b) {
00242 return a + Modular (b); }
00243
00244 TMPL inline Modular
00245 operator + (const typename M::base& a, const Modular& b) {
00246 return Modular (a) + b; }
00247
00248 TMPL inline Modular
00249 operator - (const Modular& a) {
00250 Modular c;
00251 neg_mod (c.rep, a.rep, Modular::get_modulus ());
00252 return c; }
00253
00254 TMPL inline Modular
00255 operator - (const Modular& a, const Modular& b) {
00256 Modular c;
00257 sub_mod (c.rep, a.rep, b.rep, Modular::get_modulus ());
00258 return c; }
00259
00260 TMPL inline Modular
00261 operator - (const Modular& a, const typename M::base& b) {
00262 return a - Modular (b); }
00263
00264 TMPL inline Modular
00265 operator - (const typename M::base& a, const Modular& b) {
00266 return Modular (a) - b; }
00267
00268 TMPL inline Modular
00269 operator * (const Modular& a, const Modular& b) {
00270 Modular c;
00271 mul_mod (c.rep, a.rep, b.rep, Modular::get_modulus ());
00272 return c; }
00273
00274 TMPL inline Modular
00275 operator * (const Modular& a, const typename M::base& b) {
00276 return a * Modular (b); }
00277
00278 TMPL inline Modular
00279 operator * (const typename M::base& a, const Modular& b) {
00280 return Modular (a) * b; }
00281
00282 TMPL inline Modular
00283 operator / (const Modular& a, const Modular& b) {
00284 Modular c;
00285 div_mod (c.rep, a.rep, b.rep, Modular::get_modulus ());
00286 return c; }
00287
00288 TMPL inline Modular
00289 operator / (const Modular& a, const typename M::base& b) {
00290 return a / Modular (b); }
00291
00292 TMPL inline Modular
00293 operator / (const typename M::base& a, const Modular& b) {
00294 return Modular (a) / b; }
00295
00296 TMPL inline Modular
00297 invert (const Modular& a) {
00298 Modular c;
00299 inv_mod (c.rep, a.rep, Modular::get_modulus ());
00300 return c; }
00301
00302 TMPL inline bool
00303 is_invertible (const Modular& a) {
00304 return is_invertible_mod (a.rep, Modular::get_modulus ()); }
00305
00306
00307
00308
00309
00310 TMPL inline Modular
00311 gcd (const Modular& a, const Modular& b) {
00312 if (a == 0) return b;
00313 return a; }
00314
00315 TMPL inline Modular
00316 lcm (const Modular& a, const Modular& b) {
00317 return a * b; }
00318
00319
00320
00321
00322
00323 template<typename C, typename V, typename W> inline bool
00324 operator == (const modular<modulus<C,W>, V>& m, const C& c) {
00325 return m == modular<modulus<C,W>, V> (c); }
00326 template<typename C, typename V, typename W> inline bool
00327 operator != (const modular<modulus<C,W>, V>& m, const C& c) {
00328 return m != modular<modulus<C,W>, V> (c); }
00329 template<typename C, typename V, typename W> inline bool
00330 operator == (const C& c, const modular<modulus<C,W>, V>& m) {
00331 return modular<modulus<C,W>, V> (c) == m; }
00332 template<typename C, typename V, typename W> inline bool
00333 operator != (const C& c, const modular<modulus<C,W>, V>& m) {
00334 return modular<modulus<C,W>, V> (c) != m; }
00335
00336
00337
00338
00339
00340 WRAP_PRINT_IMPL(TMPL inline,Modular)
00341
00342
00343
00344
00345
00346 ARITH_INT_SUGAR(TMPL,Modular)
00347 POOR_MAN_SQRT_SUGAR(TMPL,Modular)
00348 POOR_MAN_ELEMENTARY_SUGAR(TMPL,Modular)
00349
00350 #undef TMPL_DEF
00351 #undef TMPL
00352 #undef Modulus
00353 #undef Modular
00354
00355
00356
00357
00358
00359 #define TMPL template<typename M>
00360 #define Modulus modulus<C,V>
00361 #define Modular modular<M,modular_local>
00362
00363 TMPL
00364 class modular<M, modular_local> {
00365 MMX_ALLOCATORS
00366 public:
00367 typedef M modulus;
00368 typedef modular_local V;
00369 typedef typename M::base C;
00370 typedef typename V::template modulus_storage<M> S;
00371
00372 C rep;
00373 M mod;
00374 static inline M get_modulus () { return S::get_modulus (); }
00375 static inline void set_modulus (const M& p) { S::set_modulus (p); }
00376
00377 inline M get_local_modulus () const { return mod; }
00378 inline void set_local_modulus (const M& p) { mod= p; S::set_modulus (p); }
00379
00380 inline C operator * () const {
00381 C dest;
00382 decode_mod (dest, rep, mod);
00383 return dest; }
00384
00385 inline modular () {
00386 mod= get_modulus ();
00387 encode_mod (rep, C (), mod); }
00388
00389 inline modular (const Modular& s) : rep (s.rep), mod (s.mod) {}
00390
00391 inline modular (const C& s, bool reduced) {
00392 mod= get_modulus ();
00393 if (reduced) rep= s;
00394 else encode_mod (rep, s, mod); }
00395
00396 inline modular (const C& s, const M& p, bool reduced=false) {
00397 mod= p;
00398 set_modulus (p);
00399 if (reduced) rep= s;
00400 else encode_mod (rep, s, mod); }
00401
00402 template<typename O>
00403 inline modular (const O& a) {
00404 mod= get_modulus ();
00405 encode_mod (rep, as<C> (a), mod); }
00406
00407 inline Modular& operator = (const Modular& a) {
00408 rep = a.rep;
00409 mod = a.mod;
00410 return *this; }
00411
00412 inline bool operator == (const Modular& a) const {
00413 if (mod == a.mod) return rep == a.rep;
00414 if (mod == C(0)) return Modular (* (*this), a.mod) == a;
00415 if (a.mod == C(0)) return Modular (* a, mod) == *this;
00416 return false; }
00417
00418 inline bool operator != (const Modular& a) const {
00419 if (mod == a.mod) return rep != a.rep;
00420 if (mod == C(0)) return Modular (* (*this), a.mod) != a;
00421 if (a.mod == C(0)) return Modular (* a, mod) != *this;
00422 return true; }
00423 };
00424
00425 template<typename M> inline M
00426 get_modulus (const modular<M, modular_local>& x) {
00427 return x.get_local_modulus (); }
00428
00429 template<typename M> inline void
00430 set_modulus (Modular& x, const M& p) {
00431 x.set_local_modulus (p); }
00432
00433 template<typename TM, typename SM, typename SV>
00434 struct as_helper<modular<TM,modular_local>,modular<SM,SV> > {
00435 static inline modular<TM,modular_local>
00436 cv (const modular<SM,SV>& p) {
00437 return modular<TM,modular_local> (* p, get_modulus (p)); }
00438 };
00439
00440 TMPL nat
00441 hash (const Modular& c) {
00442 nat h= hash (*c);
00443 return (h<<1) ^ (h<<5) ^ (h>>29) ^ hash (get_modulus (c)); }
00444
00445 TMPL nat
00446 exact_hash (const Modular& c) {
00447 nat h= exact_hash (*c);
00448 return (h<<1) ^ (h<<5) ^ (h>>29) ^ exact_hash (get_modulus (c)); }
00449
00450 TMPL nat hard_hash (const Modular& c) {
00451 nat h= hard_hash (*c);
00452 return (h<<1) ^ (h<<5) ^ (h>>29) ^ hard_hash (get_modulus (c)); }
00453
00454 TMPL bool exact_eq (const Modular& c1, const Modular& c2) {
00455 return exact_eq (*c1, *c2) &&
00456 exact_eq (get_modulus (c1), get_modulus (c2)); }
00457
00458 TMPL bool exact_neq (const Modular& c1, const Modular& c2) {
00459 return ! exact_eq (c1, c2); }
00460
00461 TMPL bool hard_eq (const Modular& c1, const Modular& c2) {
00462 return hard_eq (*c1, *c2) &&
00463 hard_eq (get_modulus (c1), get_modulus (c2)); }
00464
00465 TMPL bool hard_neq (const Modular& c1, const Modular& c2) {
00466 return ! hard_eq (c1, c2); }
00467
00468 TMPL inline Modular&
00469 operator += (Modular& dest, const Modular& s) {
00470 M md= get_modulus (dest), ms= get_modulus (s);
00471 M mod= md == 0 ? ms : md;
00472 VERIFY (md == 0 || ms == 0 || md == ms, "incompatible moduli");
00473 add_mod (dest.rep, s.rep, mod);
00474 return dest; }
00475
00476 TMPL inline Modular&
00477 operator -= (Modular& dest, const Modular& s) {
00478 M md= get_modulus (dest), ms= get_modulus (s);
00479 M mod= md == 0 ? ms : md;
00480 VERIFY (md == 0 || ms == 0 || md == ms, "incompatible moduli");
00481 sub_mod (dest.rep, s.rep, mod);
00482 return dest; }
00483
00484 TMPL inline Modular&
00485 operator *= (Modular& dest, const Modular& s) {
00486 M md= get_modulus (dest), ms= get_modulus (s);
00487 M mod= md == 0 ? ms : md;
00488 VERIFY (md == 0 || ms == 0 || md == ms, "incompatible moduli");
00489 mul_mod (dest.rep, s.rep, mod);
00490 return dest; }
00491
00492 TMPL inline Modular&
00493 operator /= (Modular& dest, const Modular& s) {
00494 M md= get_modulus (dest), ms= get_modulus (s);
00495 M mod= md == 0 ? ms : md;
00496 VERIFY (md == 0 || ms == 0 || md == ms, "incompatible moduli");
00497 div_mod (dest.rep, s.rep, mod);
00498 return dest; }
00499
00500 TMPL inline Modular
00501 operator + (const Modular& a, const Modular& b) {
00502 M ma= get_modulus (a), mb= get_modulus (b);
00503 M mod= ma == 0 ? mb : ma;
00504 VERIFY (ma == 0 || mb == 0 || ma == mb, "incompatible moduli");
00505 Modular c;
00506 add_mod (c.rep, a.rep, b.rep, mod);
00507 c.mod= mod;
00508 return c; }
00509
00510 TMPL inline Modular
00511 operator - (const Modular& a) {
00512 M mod= get_modulus (a);
00513 Modular c;
00514 neg_mod (c.rep, a.rep, mod);
00515 c.mod= mod;
00516 return c; }
00517
00518 TMPL inline Modular
00519 operator - (const Modular& a, const Modular& b) {
00520 M ma= get_modulus (a), mb= get_modulus (b);
00521 M mod= ma == 0 ? mb : ma;
00522 VERIFY (ma == 0 || mb == 0 || ma == mb, "incompatible moduli");
00523 Modular c;
00524 sub_mod (c.rep, a.rep, b.rep, mod);
00525 c.mod= mod;
00526 return c; }
00527
00528 TMPL inline Modular
00529 operator * (const Modular& a, const Modular& b) {
00530 M ma= get_modulus (a), mb= get_modulus (b);
00531 M mod= ma == 0 ? mb : ma;
00532 VERIFY (ma == 0 || mb == 0 || ma == mb, "incompatible moduli");
00533 Modular c;
00534 mul_mod (c.rep, a.rep, b.rep, mod);
00535 c.mod= mod;
00536 return c; }
00537
00538 TMPL inline Modular
00539 operator / (const Modular& a, const Modular& b) {
00540 M ma= get_modulus (a), mb= get_modulus (b);
00541 M mod= ma == 0 ? mb : ma;
00542 VERIFY (ma == 0 || mb == 0 || ma == mb, "incompatible moduli");
00543 Modular c;
00544 div_mod (c.rep, a.rep, b.rep, mod);
00545 c.mod= mod;
00546 return c; }
00547
00548 TMPL inline Modular
00549 invert (const Modular& a) {
00550 M mod= get_modulus (a);
00551 Modular c;
00552 inv_mod (c.rep, a.rep, mod);
00553 c.mod= mod;
00554 return c; }
00555
00556 TMPL inline bool
00557 is_invertible (const Modular& a) {
00558 return is_invertible_mod (a.rep, get_modulus (a)); }
00559
00560
00561
00562
00563
00564 #define mmx_modular(R) modular<modulus<R >,modular_local>
00565
00566 #undef TMPL
00567 #undef Modulus
00568 #undef Modular
00569 }
00570 #endif //__MMX__MODULAR__HPP