[35aab3] | 1 | /**************************************** |
---|
| 2 | * Computer Algebra System SINGULAR * |
---|
| 3 | ****************************************/ |
---|
| 4 | /* |
---|
| 5 | * ABSTRACT: numbers modulo p (<=32003) |
---|
| 6 | */ |
---|
| 7 | |
---|
[16f511] | 8 | #ifdef HAVE_CONFIG_H |
---|
[ba5e9e] | 9 | #include "libpolysconfig.h" |
---|
[16f511] | 10 | #endif /* HAVE_CONFIG_H */ |
---|
[22a09d] | 11 | #include <misc/auxiliary.h> |
---|
| 12 | |
---|
| 13 | #ifdef HAVE_FACTORY |
---|
| 14 | #include <factory/factory.h> |
---|
| 15 | #endif |
---|
| 16 | |
---|
| 17 | #include <string.h> |
---|
[18cb65] | 18 | #include <omalloc/omalloc.h> |
---|
[2d805a] | 19 | #include <coeffs/coeffs.h> |
---|
[31213a4] | 20 | #include <reporter/reporter.h> |
---|
[2d805a] | 21 | #include <coeffs/numbers.h> |
---|
| 22 | #include <coeffs/longrat.h> |
---|
| 23 | #include <coeffs/mpr_complex.h> |
---|
[31213a4] | 24 | #include <misc/mylimits.h> |
---|
[2d805a] | 25 | #include <coeffs/modulop.h> |
---|
[35aab3] | 26 | |
---|
[11206af] | 27 | // int npGen=0; |
---|
| 28 | |
---|
[73a9ffb] | 29 | /// Our Type! |
---|
| 30 | static const n_coeffType ID = n_Zp; |
---|
[35aab3] | 31 | |
---|
[40e88db] | 32 | #ifdef NV_OPS |
---|
| 33 | #pragma GCC diagnostic ignored "-Wlong-long" |
---|
| 34 | static inline number nvMultM(number a, number b, const coeffs r) |
---|
| 35 | { |
---|
| 36 | assume( getCoeffType(r) == ID ); |
---|
| 37 | |
---|
| 38 | #if SIZEOF_LONG == 4 |
---|
| 39 | #define ULONG64 (unsigned long long)(unsigned long) |
---|
| 40 | #else |
---|
| 41 | #define ULONG64 (unsigned long) |
---|
| 42 | #endif |
---|
| 43 | return (number) |
---|
| 44 | (unsigned long)((ULONG64 a)*(ULONG64 b) % (ULONG64 r->ch)); |
---|
| 45 | } |
---|
| 46 | number nvMult (number a, number b, const coeffs r); |
---|
| 47 | number nvDiv (number a, number b, const coeffs r); |
---|
| 48 | number nvInvers (number c, const coeffs r); |
---|
| 49 | void nvPower (number a, int i, number * result, const coeffs r); |
---|
| 50 | #endif |
---|
| 51 | |
---|
| 52 | |
---|
| 53 | |
---|
| 54 | |
---|
[7d90aa] | 55 | BOOLEAN npGreaterZero (number k, const coeffs r) |
---|
[35aab3] | 56 | { |
---|
[0461f0] | 57 | assume( n_Test(k, r) ); |
---|
| 58 | |
---|
[7447d8] | 59 | int h = (int)((long) k); |
---|
[e77676] | 60 | return ((int)h !=0) && (h <= (r->ch>>1)); |
---|
[35aab3] | 61 | } |
---|
| 62 | |
---|
[e77676] | 63 | //unsigned long npMultMod(unsigned long a, unsigned long b, int npPrimeM) |
---|
[35aab3] | 64 | //{ |
---|
| 65 | // unsigned long c = a*b; |
---|
| 66 | // c = c % npPrimeM; |
---|
[e77676] | 67 | // assume(c == (unsigned long) npMultM((number) a, (number) b, npPrimeM)); |
---|
[35aab3] | 68 | // return c; |
---|
| 69 | //} |
---|
| 70 | |
---|
[7d90aa] | 71 | number npMult (number a,number b, const coeffs r) |
---|
[35aab3] | 72 | { |
---|
[0461f0] | 73 | assume( n_Test(a, r) ); |
---|
| 74 | assume( n_Test(b, r) ); |
---|
| 75 | |
---|
[35aab3] | 76 | if (((long)a == 0) || ((long)b == 0)) |
---|
| 77 | return (number)0; |
---|
[0461f0] | 78 | number c = npMultM(a,b, r); |
---|
| 79 | assume( n_Test(c, r) ); |
---|
| 80 | return c; |
---|
[35aab3] | 81 | } |
---|
| 82 | |
---|
| 83 | /*2 |
---|
| 84 | * create a number from int |
---|
| 85 | */ |
---|
[2f3764] | 86 | number npInit (long i, const coeffs r) |
---|
[35aab3] | 87 | { |
---|
[d8b0997] | 88 | long ii=i % (long)r->ch; |
---|
| 89 | if (ii < 0L) ii += (long)r->ch; |
---|
| 90 | |
---|
[0461f0] | 91 | number c = (number)ii; |
---|
| 92 | assume( n_Test(c, r) ); |
---|
| 93 | return c; |
---|
| 94 | |
---|
[35aab3] | 95 | } |
---|
| 96 | |
---|
[0461f0] | 97 | |
---|
[35aab3] | 98 | /*2 |
---|
[0461f0] | 99 | * convert a number to an int in (-p/2 .. p/2] |
---|
| 100 | */ |
---|
[7d90aa] | 101 | int npInt(number &n, const coeffs r) |
---|
[35aab3] | 102 | { |
---|
[0461f0] | 103 | assume( n_Test(n, r) ); |
---|
| 104 | |
---|
[03579d] | 105 | if ((long)n > (((long)r->ch) >>1)) return (int)((long)n -((long)r->ch)); |
---|
| 106 | else return (int)((long)n); |
---|
[35aab3] | 107 | } |
---|
| 108 | |
---|
[7d90aa] | 109 | number npAdd (number a, number b, const coeffs r) |
---|
[35aab3] | 110 | { |
---|
[0461f0] | 111 | assume( n_Test(a, r) ); |
---|
| 112 | assume( n_Test(b, r) ); |
---|
| 113 | |
---|
| 114 | number c = npAddM(a,b, r); |
---|
| 115 | |
---|
| 116 | assume( n_Test(c, r) ); |
---|
| 117 | |
---|
| 118 | return c; |
---|
[35aab3] | 119 | } |
---|
| 120 | |
---|
[7d90aa] | 121 | number npSub (number a, number b, const coeffs r) |
---|
[35aab3] | 122 | { |
---|
[0461f0] | 123 | assume( n_Test(a, r) ); |
---|
| 124 | assume( n_Test(b, r) ); |
---|
| 125 | |
---|
| 126 | number c = npSubM(a,b,r); |
---|
| 127 | |
---|
| 128 | assume( n_Test(c, r) ); |
---|
| 129 | |
---|
| 130 | return c; |
---|
[35aab3] | 131 | } |
---|
| 132 | |
---|
[0461f0] | 133 | BOOLEAN npIsZero (number a, const coeffs r) |
---|
[35aab3] | 134 | { |
---|
[0461f0] | 135 | assume( n_Test(a, r) ); |
---|
| 136 | |
---|
[35aab3] | 137 | return 0 == (long)a; |
---|
| 138 | } |
---|
| 139 | |
---|
[0461f0] | 140 | BOOLEAN npIsOne (number a, const coeffs r) |
---|
[35aab3] | 141 | { |
---|
[0461f0] | 142 | assume( n_Test(a, r) ); |
---|
| 143 | |
---|
[35aab3] | 144 | return 1 == (long)a; |
---|
| 145 | } |
---|
| 146 | |
---|
[7d90aa] | 147 | BOOLEAN npIsMOne (number a, const coeffs r) |
---|
[35aab3] | 148 | { |
---|
[0461f0] | 149 | assume( n_Test(a, r) ); |
---|
| 150 | |
---|
[7d90aa] | 151 | return ((r->npPminus1M == (long)a)&&((long)1!=(long)a)); |
---|
[35aab3] | 152 | } |
---|
| 153 | |
---|
| 154 | #ifdef HAVE_DIV_MOD |
---|
[a74143] | 155 | |
---|
[1fc18b] | 156 | #ifdef USE_NTL_XGCD |
---|
[a74143] | 157 | |
---|
[1fc18b] | 158 | //ifdef HAVE_NTL // in ntl.a |
---|
[35aab3] | 159 | //extern void XGCD(long& d, long& s, long& t, long a, long b); |
---|
| 160 | #include <NTL/ZZ.h> |
---|
| 161 | #ifdef NTL_CLIENT |
---|
| 162 | NTL_CLIENT |
---|
| 163 | #endif |
---|
[a74143] | 164 | |
---|
[1fc18b] | 165 | #endif |
---|
[35aab3] | 166 | |
---|
[7d90aa] | 167 | long InvMod(long a, const coeffs R) |
---|
[1fc18b] | 168 | { |
---|
| 169 | long d, s, t; |
---|
[35aab3] | 170 | |
---|
[1fc18b] | 171 | #ifdef USE_NTL_XGCD |
---|
[e77676] | 172 | XGCD(d, s, t, a, R->ch); |
---|
[1fc18b] | 173 | assume (d == 1); |
---|
| 174 | #else |
---|
| 175 | long u, v, u0, v0, u1, v1, u2, v2, q, r; |
---|
[35aab3] | 176 | |
---|
[1fc18b] | 177 | assume(a>0); |
---|
| 178 | u1=1; u2=0; |
---|
[e77676] | 179 | u = a; v = R->ch; |
---|
[35aab3] | 180 | |
---|
[1fc18b] | 181 | while (v != 0) |
---|
| 182 | { |
---|
[35aab3] | 183 | q = u / v; |
---|
| 184 | r = u % v; |
---|
| 185 | u = v; |
---|
| 186 | v = r; |
---|
| 187 | u0 = u2; |
---|
[1fc18b] | 188 | u2 = u1 - q*u2; |
---|
[35aab3] | 189 | u1 = u0; |
---|
| 190 | } |
---|
| 191 | |
---|
[1fc18b] | 192 | assume(u==1); |
---|
[35aab3] | 193 | s = u1; |
---|
| 194 | #endif |
---|
| 195 | if (s < 0) |
---|
[e77676] | 196 | return s + R->ch; |
---|
[35aab3] | 197 | else |
---|
| 198 | return s; |
---|
| 199 | } |
---|
| 200 | #endif |
---|
| 201 | |
---|
[7d90aa] | 202 | inline number npInversM (number c, const coeffs r) |
---|
[35aab3] | 203 | { |
---|
[0461f0] | 204 | assume( n_Test(c, r) ); |
---|
[35aab3] | 205 | #ifndef HAVE_DIV_MOD |
---|
[0461f0] | 206 | number d = (number)(long)r->npExpTable[r->npPminus1M - r->npLogTable[(long)c]]; |
---|
[35aab3] | 207 | #else |
---|
[7d90aa] | 208 | long inv=(long)r->npInvTable[(long)c]; |
---|
[35aab3] | 209 | if (inv==0) |
---|
| 210 | { |
---|
[7d90aa] | 211 | inv=InvMod((long)c,r); |
---|
| 212 | r->npInvTable[(long)c]=inv; |
---|
[35aab3] | 213 | } |
---|
[0461f0] | 214 | number d = (number)inv; |
---|
[35aab3] | 215 | #endif |
---|
[0461f0] | 216 | assume( n_Test(d, r) ); |
---|
| 217 | return d; |
---|
| 218 | |
---|
[35aab3] | 219 | } |
---|
| 220 | |
---|
[7d90aa] | 221 | number npDiv (number a,number b, const coeffs r) |
---|
[35aab3] | 222 | { |
---|
[0461f0] | 223 | assume( n_Test(a, r) ); |
---|
| 224 | assume( n_Test(b, r) ); |
---|
| 225 | |
---|
[6c15ec] | 226 | //#ifdef NV_OPS |
---|
[e77676] | 227 | // if (r->ch>NV_MAX_PRIME) |
---|
[6c15ec] | 228 | // return nvDiv(a,b); |
---|
| 229 | //#endif |
---|
[35aab3] | 230 | if ((long)a==0) |
---|
| 231 | return (number)0; |
---|
[0461f0] | 232 | number d; |
---|
| 233 | |
---|
[35aab3] | 234 | #ifndef HAVE_DIV_MOD |
---|
[6c15ec] | 235 | if ((long)b==0) |
---|
[35aab3] | 236 | { |
---|
[b7e838] | 237 | WerrorS(nDivBy0); |
---|
[35aab3] | 238 | return (number)0; |
---|
| 239 | } |
---|
[0461f0] | 240 | |
---|
| 241 | int s = r->npLogTable[(long)a] - r->npLogTable[(long)b]; |
---|
| 242 | if (s < 0) |
---|
| 243 | s += r->npPminus1M; |
---|
| 244 | d = (number)(long)r->npExpTable[s]; |
---|
[35aab3] | 245 | #else |
---|
[7d90aa] | 246 | number inv=npInversM(b,r); |
---|
[0461f0] | 247 | d = npMultM(a,inv,r); |
---|
[35aab3] | 248 | #endif |
---|
[0461f0] | 249 | |
---|
| 250 | assume( n_Test(d, r) ); |
---|
| 251 | return d; |
---|
| 252 | |
---|
[35aab3] | 253 | } |
---|
[7d90aa] | 254 | number npInvers (number c, const coeffs r) |
---|
[35aab3] | 255 | { |
---|
[0461f0] | 256 | assume( n_Test(c, r) ); |
---|
| 257 | |
---|
[35aab3] | 258 | if ((long)c==0) |
---|
| 259 | { |
---|
| 260 | WerrorS("1/0"); |
---|
| 261 | return (number)0; |
---|
| 262 | } |
---|
[0461f0] | 263 | number d = npInversM(c,r); |
---|
| 264 | |
---|
| 265 | assume( n_Test(d, r) ); |
---|
| 266 | return d; |
---|
| 267 | |
---|
[35aab3] | 268 | } |
---|
| 269 | |
---|
[7d90aa] | 270 | number npNeg (number c, const coeffs r) |
---|
[35aab3] | 271 | { |
---|
[0461f0] | 272 | assume( n_Test(c, r) ); |
---|
| 273 | |
---|
[35aab3] | 274 | if ((long)c==0) return c; |
---|
[0461f0] | 275 | |
---|
| 276 | #if 0 |
---|
| 277 | number d = npNegM(c,r); |
---|
| 278 | assume( n_Test(d, r) ); |
---|
| 279 | return d; |
---|
| 280 | #else |
---|
| 281 | c = npNegM(c,r); |
---|
| 282 | assume( n_Test(c, r) ); |
---|
| 283 | return c; |
---|
| 284 | #endif |
---|
[35aab3] | 285 | } |
---|
| 286 | |
---|
[0461f0] | 287 | BOOLEAN npGreater (number a,number b, const coeffs r) |
---|
[35aab3] | 288 | { |
---|
[0461f0] | 289 | assume( n_Test(a, r) ); |
---|
| 290 | assume( n_Test(b, r) ); |
---|
| 291 | |
---|
[35aab3] | 292 | //return (long)a != (long)b; |
---|
| 293 | return (long)a > (long)b; |
---|
| 294 | } |
---|
| 295 | |
---|
[7d90aa] | 296 | BOOLEAN npEqual (number a,number b, const coeffs r) |
---|
[35aab3] | 297 | { |
---|
[0461f0] | 298 | assume( n_Test(a, r) ); |
---|
| 299 | assume( n_Test(b, r) ); |
---|
| 300 | |
---|
[35aab3] | 301 | // return (long)a == (long)b; |
---|
[0461f0] | 302 | |
---|
[7d90aa] | 303 | return npEqualM(a,b,r); |
---|
[35aab3] | 304 | } |
---|
| 305 | |
---|
[7d90aa] | 306 | void npWrite (number &a, const coeffs r) |
---|
[35aab3] | 307 | { |
---|
[0461f0] | 308 | assume( n_Test(a, r) ); |
---|
| 309 | |
---|
[493225] | 310 | if ((long)a>(((long)r->ch) >>1)) StringAppend("-%d",(int)(((long)r->ch)-((long)a))); |
---|
| 311 | else StringAppend("%d",(int)((long)a)); |
---|
[35aab3] | 312 | } |
---|
| 313 | |
---|
[7d90aa] | 314 | void npPower (number a, int i, number * result, const coeffs r) |
---|
[35aab3] | 315 | { |
---|
[0461f0] | 316 | assume( n_Test(a, r) ); |
---|
| 317 | |
---|
[35aab3] | 318 | if (i==0) |
---|
| 319 | { |
---|
| 320 | //npInit(1,result); |
---|
| 321 | *(long *)result = 1; |
---|
| 322 | } |
---|
| 323 | else if (i==1) |
---|
| 324 | { |
---|
| 325 | *result = a; |
---|
| 326 | } |
---|
| 327 | else |
---|
| 328 | { |
---|
[7d90aa] | 329 | npPower(a,i-1,result,r); |
---|
| 330 | *result = npMultM(a,*result,r); |
---|
[35aab3] | 331 | } |
---|
| 332 | } |
---|
| 333 | |
---|
[7d90aa] | 334 | static const char* npEati(const char *s, int *i, const coeffs r) |
---|
[35aab3] | 335 | { |
---|
| 336 | |
---|
| 337 | if (((*s) >= '0') && ((*s) <= '9')) |
---|
| 338 | { |
---|
[923b64] | 339 | unsigned long ii=0L; |
---|
[35aab3] | 340 | do |
---|
| 341 | { |
---|
[923b64] | 342 | ii *= 10; |
---|
| 343 | ii += *s++ - '0'; |
---|
[e77676] | 344 | if (ii >= (MAX_INT_VAL / 10)) ii = ii % r->ch; |
---|
[35aab3] | 345 | } |
---|
| 346 | while (((*s) >= '0') && ((*s) <= '9')); |
---|
[e77676] | 347 | if (ii >= r->ch) ii = ii % r->ch; |
---|
[923b64] | 348 | *i=(int)ii; |
---|
[35aab3] | 349 | } |
---|
| 350 | else (*i) = 1; |
---|
| 351 | return s; |
---|
| 352 | } |
---|
| 353 | |
---|
[7d90aa] | 354 | const char * npRead (const char *s, number *a, const coeffs r) |
---|
[35aab3] | 355 | { |
---|
| 356 | int z; |
---|
| 357 | int n=1; |
---|
| 358 | |
---|
[7d90aa] | 359 | s = npEati(s, &z, r); |
---|
[35aab3] | 360 | if ((*s) == '/') |
---|
| 361 | { |
---|
| 362 | s++; |
---|
[7d90aa] | 363 | s = npEati(s, &n, r); |
---|
[35aab3] | 364 | } |
---|
| 365 | if (n == 1) |
---|
| 366 | *a = (number)z; |
---|
[6c15ec] | 367 | else |
---|
[100025] | 368 | { |
---|
| 369 | if ((z==0)&&(n==0)) WerrorS(nDivBy0); |
---|
[35aab3] | 370 | else |
---|
[100025] | 371 | { |
---|
| 372 | #ifdef NV_OPS |
---|
[e77676] | 373 | if (r->ch>NV_MAX_PRIME) |
---|
[7d90aa] | 374 | *a = nvDiv((number)z,(number)n,r); |
---|
[100025] | 375 | else |
---|
[6c15ec] | 376 | #endif |
---|
[7d90aa] | 377 | *a = npDiv((number)z,(number)n,r); |
---|
[100025] | 378 | } |
---|
| 379 | } |
---|
[0461f0] | 380 | assume( n_Test(*a, r) ); |
---|
[35aab3] | 381 | return s; |
---|
| 382 | } |
---|
| 383 | |
---|
| 384 | /*2 |
---|
| 385 | * set the charcteristic (allocate and init tables) |
---|
| 386 | */ |
---|
| 387 | |
---|
[8c484e] | 388 | void npKillChar(coeffs r) |
---|
| 389 | { |
---|
| 390 | #ifdef HAVE_DIV_MOD |
---|
| 391 | if (r->npInvTable!=NULL) |
---|
[e77676] | 392 | omFreeSize( (void *)r->npInvTable, r->ch*sizeof(unsigned short) ); |
---|
[8c484e] | 393 | r->npInvTable=NULL; |
---|
| 394 | #else |
---|
| 395 | if (r->npExpTable!=NULL) |
---|
| 396 | { |
---|
[e77676] | 397 | omFreeSize( (void *)r->npExpTable, r->ch*sizeof(unsigned short) ); |
---|
| 398 | omFreeSize( (void *)r->npLogTable, r->ch*sizeof(unsigned short) ); |
---|
[8c484e] | 399 | r->npExpTable=NULL; r->npLogTable=NULL; |
---|
| 400 | } |
---|
| 401 | #endif |
---|
| 402 | } |
---|
| 403 | |
---|
[aff5ae] | 404 | static BOOLEAN npCoeffsEqual(const coeffs r, n_coeffType n, void * parameter) |
---|
[d2cf18] | 405 | { |
---|
| 406 | /* test, if r is an instance of nInitCoeffs(n,parameter) */ |
---|
[aff5ae] | 407 | return (n==n_Zp) && (r->ch==(int)(long)parameter); |
---|
[d2cf18] | 408 | } |
---|
[eca225] | 409 | #ifdef HAVE_FACTORY |
---|
[abb4787] | 410 | CanonicalForm npConvSingNFactoryN( number n, BOOLEAN setChar, const coeffs r ) |
---|
[eca225] | 411 | { |
---|
[e77676] | 412 | if (setChar) setCharacteristic( r->ch ); |
---|
[eca225] | 413 | CanonicalForm term(npInt( n,r )); |
---|
| 414 | return term; |
---|
| 415 | } |
---|
| 416 | |
---|
| 417 | number npConvFactoryNSingN( const CanonicalForm n, const coeffs r) |
---|
| 418 | { |
---|
| 419 | if (n.isImm()) |
---|
| 420 | { |
---|
[a934fc8] | 421 | return npInit(n.intval(),r); |
---|
[eca225] | 422 | } |
---|
| 423 | else |
---|
| 424 | { |
---|
| 425 | assume(0); |
---|
[c262d5] | 426 | return NULL; |
---|
[eca225] | 427 | } |
---|
| 428 | } |
---|
| 429 | #endif |
---|
| 430 | |
---|
[d2cf18] | 431 | |
---|
[1cce47] | 432 | BOOLEAN npInitChar(coeffs r, void* p) |
---|
[35aab3] | 433 | { |
---|
[73a9ffb] | 434 | assume( getCoeffType(r) == ID ); |
---|
[11206af] | 435 | const int c = (int) (long) p; |
---|
| 436 | |
---|
| 437 | assume( c > 0 ); |
---|
| 438 | |
---|
[35aab3] | 439 | int i, w; |
---|
| 440 | |
---|
[e77676] | 441 | r->ch = c; |
---|
| 442 | r->npPminus1M = c /*r->ch*/ - 1; |
---|
[d2cf18] | 443 | |
---|
| 444 | //r->cfInitChar=npInitChar; |
---|
| 445 | r->cfKillChar=npKillChar; |
---|
| 446 | r->nCoeffIsEqual=npCoeffsEqual; |
---|
| 447 | |
---|
[7bbbef] | 448 | r->cfMult = npMult; |
---|
| 449 | r->cfSub = npSub; |
---|
| 450 | r->cfAdd = npAdd; |
---|
| 451 | r->cfDiv = npDiv; |
---|
| 452 | r->cfIntDiv= npDiv; |
---|
| 453 | //r->cfIntMod= ndIntMod; |
---|
| 454 | r->cfExactDiv= npDiv; |
---|
[d2cf18] | 455 | r->cfInit = npInit; |
---|
[7bbbef] | 456 | //r->cfSize = ndSize; |
---|
| 457 | r->cfInt = npInt; |
---|
[d2cf18] | 458 | #ifdef HAVE_RINGS |
---|
[7bbbef] | 459 | //r->cfDivComp = NULL; // only for ring stuff |
---|
| 460 | //r->cfIsUnit = NULL; // only for ring stuff |
---|
| 461 | //r->cfGetUnit = NULL; // only for ring stuff |
---|
| 462 | //r->cfExtGcd = NULL; // only for ring stuff |
---|
| 463 | // r->cfDivBy = NULL; // only for ring stuff |
---|
[d2cf18] | 464 | #endif |
---|
[7bbbef] | 465 | r->cfNeg = npNeg; |
---|
| 466 | r->cfInvers= npInvers; |
---|
[aff5ae] | 467 | //r->cfCopy = ndCopy; |
---|
[7bbbef] | 468 | //r->cfRePart = ndCopy; |
---|
| 469 | //r->cfImPart = ndReturn0; |
---|
[ce1f78] | 470 | r->cfWriteLong = npWrite; |
---|
[7bbbef] | 471 | r->cfRead = npRead; |
---|
| 472 | //r->cfNormalize=ndNormalize; |
---|
| 473 | r->cfGreater = npGreater; |
---|
| 474 | r->cfEqual = npEqual; |
---|
| 475 | r->cfIsZero = npIsZero; |
---|
| 476 | r->cfIsOne = npIsOne; |
---|
| 477 | r->cfIsMOne = npIsMOne; |
---|
| 478 | r->cfGreaterZero = npGreaterZero; |
---|
| 479 | r->cfPower = npPower; |
---|
[d2cf18] | 480 | r->cfGetDenom = ndGetDenom; |
---|
| 481 | r->cfGetNumerator = ndGetNumerator; |
---|
[7bbbef] | 482 | //r->cfGcd = ndGcd; |
---|
| 483 | //r->cfLcm = ndGcd; |
---|
[aff5ae] | 484 | //r->cfDelete= ndDelete; |
---|
[d2cf18] | 485 | r->cfSetMap = npSetMap; |
---|
[7bbbef] | 486 | //r->cfName = ndName; |
---|
| 487 | r->cfInpMult=ndInpMult; |
---|
[79020f] | 488 | r->cfInit_bigint= nlModP; // npMap0; |
---|
[d2cf18] | 489 | #ifdef NV_OPS |
---|
| 490 | if (c>NV_MAX_PRIME) |
---|
[35aab3] | 491 | { |
---|
[7bbbef] | 492 | r->cfMult = nvMult; |
---|
| 493 | r->cfDiv = nvDiv; |
---|
| 494 | r->cfExactDiv= nvDiv; |
---|
| 495 | r->cfInvers= nvInvers; |
---|
| 496 | r->cfPower= nvPower; |
---|
[d2cf18] | 497 | } |
---|
| 498 | #endif |
---|
[c7e3d7] | 499 | r->cfCoeffWrite=npCoeffWrite; |
---|
[d2cf18] | 500 | #ifdef LDEBUG |
---|
| 501 | // debug stuff |
---|
[7bbbef] | 502 | r->cfDBTest=npDBTest; |
---|
[d2cf18] | 503 | #endif |
---|
[eca225] | 504 | |
---|
| 505 | #ifdef HAVE_FACTORY |
---|
| 506 | r->convSingNFactoryN=npConvSingNFactoryN; |
---|
| 507 | r->convFactoryNSingN=npConvFactoryNSingN; |
---|
| 508 | #endif |
---|
[d2cf18] | 509 | |
---|
| 510 | // the variables: |
---|
| 511 | r->nNULL = (number)0; |
---|
[f0797c] | 512 | r->type = n_Zp; |
---|
[d2cf18] | 513 | r->ch = c; |
---|
| 514 | r->has_simple_Alloc=TRUE; |
---|
| 515 | r->has_simple_Inverse=TRUE; |
---|
| 516 | |
---|
| 517 | // the tables |
---|
[35aab3] | 518 | #ifdef NV_OPS |
---|
[e77676] | 519 | if (r->ch <=NV_MAX_PRIME) |
---|
[35aab3] | 520 | #endif |
---|
[aff5ae] | 521 | { |
---|
[35aab3] | 522 | #if !defined(HAVE_DIV_MOD) || !defined(HAVE_MULT_MOD) |
---|
[e77676] | 523 | r->npExpTable=(unsigned short *)omAlloc( r->ch*sizeof(unsigned short) ); |
---|
| 524 | r->npLogTable=(unsigned short *)omAlloc( r->ch*sizeof(unsigned short) ); |
---|
[aff5ae] | 525 | r->npExpTable[0] = 1; |
---|
| 526 | r->npLogTable[0] = 0; |
---|
[e77676] | 527 | if (r->ch > 2) |
---|
[aff5ae] | 528 | { |
---|
| 529 | w = 1; |
---|
| 530 | loop |
---|
[35aab3] | 531 | { |
---|
[aff5ae] | 532 | r->npLogTable[1] = 0; |
---|
| 533 | w++; |
---|
| 534 | i = 0; |
---|
[35aab3] | 535 | loop |
---|
| 536 | { |
---|
[aff5ae] | 537 | i++; |
---|
| 538 | r->npExpTable[i] =(int)(((long)w * (long)r->npExpTable[i-1]) |
---|
[e77676] | 539 | % r->ch); |
---|
[aff5ae] | 540 | r->npLogTable[r->npExpTable[i]] = i; |
---|
[e77676] | 541 | if (/*(i == r->ch - 1 ) ||*/ (r->npExpTable[i] == 1)) |
---|
[35aab3] | 542 | break; |
---|
| 543 | } |
---|
[e77676] | 544 | if (i == r->ch - 1) |
---|
[aff5ae] | 545 | break; |
---|
[35aab3] | 546 | } |
---|
[aff5ae] | 547 | } |
---|
| 548 | else |
---|
| 549 | { |
---|
| 550 | r->npExpTable[1] = 1; |
---|
| 551 | r->npLogTable[1] = 0; |
---|
| 552 | } |
---|
[35aab3] | 553 | #endif |
---|
| 554 | #ifdef HAVE_DIV_MOD |
---|
[e77676] | 555 | r->npInvTable=(unsigned short*)omAlloc0( r->ch*sizeof(unsigned short) ); |
---|
[35aab3] | 556 | #endif |
---|
[aff5ae] | 557 | } |
---|
[5d594a9] | 558 | return FALSE; |
---|
[35aab3] | 559 | } |
---|
| 560 | |
---|
| 561 | #ifdef LDEBUG |
---|
[2b957a] | 562 | BOOLEAN npDBTest (number a, const char *f, const int l, const coeffs r) |
---|
[35aab3] | 563 | { |
---|
[e77676] | 564 | if (((long)a<0) || ((long)a>r->ch)) |
---|
[35aab3] | 565 | { |
---|
[f1e33bb] | 566 | Print("wrong mod p number %ld at %s,%d\n",(long)a,f,l); |
---|
[35aab3] | 567 | return FALSE; |
---|
| 568 | } |
---|
| 569 | return TRUE; |
---|
| 570 | } |
---|
| 571 | #endif |
---|
| 572 | |
---|
[8c484e] | 573 | number npMapP(number from, const coeffs src, const coeffs dst_r) |
---|
[35aab3] | 574 | { |
---|
| 575 | long i = (long)from; |
---|
[e77676] | 576 | if (i>src->ch/2) |
---|
[35aab3] | 577 | { |
---|
[e77676] | 578 | i-=src->ch; |
---|
| 579 | while (i < 0) i+=dst_r->ch; |
---|
[35aab3] | 580 | } |
---|
[e77676] | 581 | i%=dst_r->ch; |
---|
[35aab3] | 582 | return (number)i; |
---|
| 583 | } |
---|
| 584 | |
---|
[c262d5] | 585 | static number npMapLongR(number from, const coeffs /*src*/, const coeffs dst_r) |
---|
[35aab3] | 586 | { |
---|
| 587 | gmp_float *ff=(gmp_float*)from; |
---|
| 588 | mpf_t *f=ff->_mpfp(); |
---|
| 589 | number res; |
---|
[a604c3] | 590 | mpz_ptr dest,ndest; |
---|
[35aab3] | 591 | int size,i; |
---|
[c262d5] | 592 | int e,al,bl; |
---|
[7447d8] | 593 | long iz; |
---|
[35aab3] | 594 | mp_ptr qp,dd,nn; |
---|
| 595 | |
---|
| 596 | size = (*f)[0]._mp_size; |
---|
| 597 | if (size == 0) |
---|
[7d90aa] | 598 | return npInit(0,dst_r); |
---|
[35aab3] | 599 | if(size<0) |
---|
| 600 | size = -size; |
---|
| 601 | |
---|
| 602 | qp = (*f)[0]._mp_d; |
---|
| 603 | while(qp[0]==0) |
---|
| 604 | { |
---|
| 605 | qp++; |
---|
| 606 | size--; |
---|
| 607 | } |
---|
| 608 | |
---|
[e77676] | 609 | if(dst_r->ch>2) |
---|
[35aab3] | 610 | e=(*f)[0]._mp_exp-size; |
---|
| 611 | else |
---|
| 612 | e=0; |
---|
[896561] | 613 | res = ALLOC_RNUMBER(); |
---|
[35aab3] | 614 | #if defined(LDEBUG) |
---|
| 615 | res->debug=123456; |
---|
| 616 | #endif |
---|
[a604c3] | 617 | dest = res->z; |
---|
[35aab3] | 618 | |
---|
[c262d5] | 619 | int in=0; |
---|
[35aab3] | 620 | if (e<0) |
---|
| 621 | { |
---|
| 622 | al = dest->_mp_size = size; |
---|
| 623 | if (al<2) al = 2; |
---|
| 624 | dd = (mp_ptr)omAlloc(sizeof(mp_limb_t)*al); |
---|
| 625 | for (i=0;i<size;i++) dd[i] = qp[i]; |
---|
| 626 | bl = 1-e; |
---|
| 627 | nn = (mp_ptr)omAlloc(sizeof(mp_limb_t)*bl); |
---|
| 628 | nn[bl-1] = 1; |
---|
| 629 | for (i=bl-2;i>=0;i--) nn[i] = 0; |
---|
[a604c3] | 630 | ndest = res->n; |
---|
[35aab3] | 631 | ndest->_mp_d = nn; |
---|
| 632 | ndest->_mp_alloc = ndest->_mp_size = bl; |
---|
| 633 | res->s = 0; |
---|
[e77676] | 634 | in=mpz_fdiv_ui(ndest,dst_r->ch); |
---|
[35aab3] | 635 | mpz_clear(ndest); |
---|
| 636 | } |
---|
| 637 | else |
---|
| 638 | { |
---|
| 639 | al = dest->_mp_size = size+e; |
---|
| 640 | if (al<2) al = 2; |
---|
| 641 | dd = (mp_ptr)omAlloc(sizeof(mp_limb_t)*al); |
---|
| 642 | for (i=0;i<size;i++) dd[i+e] = qp[i]; |
---|
| 643 | for (i=0;i<e;i++) dd[i] = 0; |
---|
| 644 | res->s = 3; |
---|
| 645 | } |
---|
| 646 | |
---|
| 647 | dest->_mp_d = dd; |
---|
| 648 | dest->_mp_alloc = al; |
---|
[e77676] | 649 | iz=mpz_fdiv_ui(dest,dst_r->ch); |
---|
[35aab3] | 650 | mpz_clear(dest); |
---|
| 651 | if(res->s==0) |
---|
[8c484e] | 652 | iz=(long)npDiv((number)iz,(number)in,dst_r); |
---|
[d663f2] | 653 | FREE_RNUMBER(res); // Q!? |
---|
[35aab3] | 654 | return (number)iz; |
---|
| 655 | } |
---|
| 656 | |
---|
[894f5b1] | 657 | #ifdef HAVE_RINGS |
---|
| 658 | /*2 |
---|
| 659 | * convert from a GMP integer |
---|
| 660 | */ |
---|
[c262d5] | 661 | number npMapGMP(number from, const coeffs /*src*/, const coeffs dst) |
---|
[894f5b1] | 662 | { |
---|
[a604c3] | 663 | int_number erg = (int_number) omAlloc(sizeof(mpz_t)); // evtl. spaeter mit bin |
---|
[894f5b1] | 664 | mpz_init(erg); |
---|
| 665 | |
---|
[e77676] | 666 | mpz_mod_ui(erg, (int_number) from, dst->ch); |
---|
[894f5b1] | 667 | number r = (number) mpz_get_si(erg); |
---|
| 668 | |
---|
| 669 | mpz_clear(erg); |
---|
[7d90aa] | 670 | omFree((void *) erg); |
---|
[894f5b1] | 671 | return (number) r; |
---|
| 672 | } |
---|
| 673 | |
---|
| 674 | /*2 |
---|
| 675 | * convert from an machine long |
---|
| 676 | */ |
---|
[c262d5] | 677 | number npMapMachineInt(number from, const coeffs /*src*/,const coeffs dst) |
---|
[894f5b1] | 678 | { |
---|
[e77676] | 679 | long i = (long) (((unsigned long) from) % dst->ch); |
---|
[894f5b1] | 680 | return (number) i; |
---|
| 681 | } |
---|
| 682 | #endif |
---|
| 683 | |
---|
[e9d796] | 684 | #ifdef HAVE_FACTORY |
---|
[c262d5] | 685 | number npMapCanonicalForm (number a, const coeffs /*src*/, const coeffs dst) |
---|
[e9d796] | 686 | { |
---|
[e77676] | 687 | setCharacteristic (dst ->ch); |
---|
[e9d796] | 688 | CanonicalForm f= CanonicalForm ((InternalCF*)(a)); |
---|
| 689 | return (number) (f.intval()); |
---|
| 690 | } |
---|
| 691 | #endif |
---|
| 692 | |
---|
[7d90aa] | 693 | nMapFunc npSetMap(const coeffs src, const coeffs dst) |
---|
[35aab3] | 694 | { |
---|
[894f5b1] | 695 | #ifdef HAVE_RINGS |
---|
[1cce47] | 696 | if (nCoeff_is_Ring_2toM(src)) |
---|
[894f5b1] | 697 | { |
---|
| 698 | return npMapMachineInt; |
---|
| 699 | } |
---|
[1cce47] | 700 | if (nCoeff_is_Ring_Z(src) || nCoeff_is_Ring_PtoM(src) || nCoeff_is_Ring_ModN(src)) |
---|
[894f5b1] | 701 | { |
---|
| 702 | return npMapGMP; |
---|
| 703 | } |
---|
| 704 | #endif |
---|
[1cce47] | 705 | if (nCoeff_is_Q(src)) |
---|
[35aab3] | 706 | { |
---|
[79020f] | 707 | return nlModP; // npMap0; |
---|
[35aab3] | 708 | } |
---|
[1cce47] | 709 | if ( nCoeff_is_Zp(src) ) |
---|
[35aab3] | 710 | { |
---|
[7d90aa] | 711 | if (n_GetChar(src) == n_GetChar(dst)) |
---|
[35aab3] | 712 | { |
---|
[8c484e] | 713 | return ndCopyMap; |
---|
[35aab3] | 714 | } |
---|
| 715 | else |
---|
| 716 | { |
---|
| 717 | return npMapP; |
---|
| 718 | } |
---|
| 719 | } |
---|
[1cce47] | 720 | if (nCoeff_is_long_R(src)) |
---|
[35aab3] | 721 | { |
---|
| 722 | return npMapLongR; |
---|
| 723 | } |
---|
[e9d796] | 724 | #ifdef HAVE_FACTORY |
---|
[1cce47] | 725 | if (nCoeff_is_CF (src)) |
---|
[e9d796] | 726 | { |
---|
| 727 | return npMapCanonicalForm; |
---|
| 728 | } |
---|
| 729 | #endif |
---|
[35aab3] | 730 | return NULL; /* default */ |
---|
| 731 | } |
---|
| 732 | |
---|
| 733 | // ----------------------------------------------------------- |
---|
| 734 | // operation for very large primes (32003< p < 2^31-1) |
---|
| 735 | // ---------------------------------------------------------- |
---|
| 736 | #ifdef NV_OPS |
---|
| 737 | |
---|
[7d90aa] | 738 | number nvMult (number a,number b, const coeffs r) |
---|
[35aab3] | 739 | { |
---|
[2e1b74] | 740 | //if (((long)a == 0) || ((long)b == 0)) |
---|
| 741 | // return (number)0; |
---|
| 742 | //else |
---|
[4d92d7] | 743 | return nvMultM(a,b,r); |
---|
[35aab3] | 744 | } |
---|
| 745 | |
---|
[4d92d7] | 746 | void nvInpMult(number &a, number b, const coeffs r) |
---|
[2e1b74] | 747 | { |
---|
[4d92d7] | 748 | number n=nvMultM(a,b,r); |
---|
[2e1b74] | 749 | a=n; |
---|
| 750 | } |
---|
| 751 | |
---|
| 752 | |
---|
[a74143] | 753 | inline long nvInvMod(long a, const coeffs R) |
---|
[35aab3] | 754 | { |
---|
[a74143] | 755 | #ifdef HAVE_DIV_MOD |
---|
| 756 | return InvMod(a, R); |
---|
| 757 | #else |
---|
| 758 | /// TODO: use "long InvMod(long a, const coeffs R)"?! |
---|
| 759 | |
---|
[c262d5] | 760 | long s; |
---|
[35aab3] | 761 | |
---|
[a74143] | 762 | long u, u0, u1, u2, q, r; // v0, v1, v2, |
---|
[35aab3] | 763 | |
---|
[a74143] | 764 | u1=1; // v1=0; |
---|
| 765 | u2=0; // v2=1; |
---|
| 766 | u = a; |
---|
| 767 | |
---|
| 768 | long v = R->ch; |
---|
[35aab3] | 769 | |
---|
[9c8927] | 770 | while (v != 0) |
---|
| 771 | { |
---|
[35aab3] | 772 | q = u / v; |
---|
| 773 | r = u % v; |
---|
| 774 | u = v; |
---|
| 775 | v = r; |
---|
| 776 | u0 = u2; |
---|
[a74143] | 777 | // v0 = v2; |
---|
| 778 | u2 = u1 - q*u2; |
---|
| 779 | // v2 = v1 - q*v2; |
---|
[35aab3] | 780 | u1 = u0; |
---|
[a74143] | 781 | // v1 = v0; |
---|
[35aab3] | 782 | } |
---|
| 783 | |
---|
| 784 | s = u1; |
---|
| 785 | //t = v1; |
---|
| 786 | if (s < 0) |
---|
[e77676] | 787 | return s + R->ch; |
---|
[35aab3] | 788 | else |
---|
[a74143] | 789 | return s; |
---|
| 790 | #endif |
---|
[35aab3] | 791 | } |
---|
| 792 | |
---|
[7d90aa] | 793 | inline number nvInversM (number c, const coeffs r) |
---|
[35aab3] | 794 | { |
---|
[7d90aa] | 795 | long inv=nvInvMod((long)c,r); |
---|
[35aab3] | 796 | return (number)inv; |
---|
| 797 | } |
---|
| 798 | |
---|
[7d90aa] | 799 | number nvDiv (number a,number b, const coeffs r) |
---|
[35aab3] | 800 | { |
---|
| 801 | if ((long)a==0) |
---|
| 802 | return (number)0; |
---|
| 803 | else if ((long)b==0) |
---|
| 804 | { |
---|
[b7e838] | 805 | WerrorS(nDivBy0); |
---|
[35aab3] | 806 | return (number)0; |
---|
| 807 | } |
---|
| 808 | else |
---|
| 809 | { |
---|
[7d90aa] | 810 | number inv=nvInversM(b,r); |
---|
| 811 | return nvMultM(a,inv,r); |
---|
[35aab3] | 812 | } |
---|
| 813 | } |
---|
[7d90aa] | 814 | number nvInvers (number c, const coeffs r) |
---|
[35aab3] | 815 | { |
---|
| 816 | if ((long)c==0) |
---|
| 817 | { |
---|
[b7e838] | 818 | WerrorS(nDivBy0); |
---|
[35aab3] | 819 | return (number)0; |
---|
| 820 | } |
---|
[7d90aa] | 821 | return nvInversM(c,r); |
---|
[35aab3] | 822 | } |
---|
[7d90aa] | 823 | void nvPower (number a, int i, number * result, const coeffs r) |
---|
[9c8927] | 824 | { |
---|
| 825 | if (i==0) |
---|
| 826 | { |
---|
| 827 | //npInit(1,result); |
---|
| 828 | *(long *)result = 1; |
---|
| 829 | } |
---|
| 830 | else if (i==1) |
---|
| 831 | { |
---|
| 832 | *result = a; |
---|
| 833 | } |
---|
| 834 | else |
---|
| 835 | { |
---|
[7d90aa] | 836 | nvPower(a,i-1,result,r); |
---|
| 837 | *result = nvMultM(a,*result,r); |
---|
[9c8927] | 838 | } |
---|
| 839 | } |
---|
[35aab3] | 840 | #endif |
---|
[c7e3d7] | 841 | |
---|
[03f7b5] | 842 | void npCoeffWrite (const coeffs r, BOOLEAN /*details*/) |
---|
[c7e3d7] | 843 | { |
---|
[e77676] | 844 | Print("// characteristic : %d\n",r->ch); |
---|
[c7e3d7] | 845 | } |
---|
| 846 | |
---|