source: git/libpolys/coeffs/rmodulo2m.cc @ 0acf3e

spielwiese
Last change on this file since 0acf3e was 45cc512, checked in by Hans Schoenemann <hannes@…>, 10 years ago
chg: rCharstr is now a wrapper for r->cf->cfCoeffString fixes also: charstr for integer,2,3
  • Property mode set to 100644
File size: 17.5 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5* ABSTRACT: numbers modulo 2^m
6*/
7
8#ifdef HAVE_CONFIG_H
9#include "libpolysconfig.h"
10#endif /* HAVE_CONFIG_H */
11#include <misc/auxiliary.h>
12
13#ifdef HAVE_RINGS
14
15#include <misc/mylimits.h>
16#include <coeffs/coeffs.h>
17#include <reporter/reporter.h>
18#include <omalloc/omalloc.h>
19#include <coeffs/numbers.h>
20#include <coeffs/longrat.h>
21#include <coeffs/mpr_complex.h>
22#include <coeffs/rmodulo2m.h>
23#include "si_gmp.h"
24
25#include <string.h>
26
27/// Our Type!
28static const n_coeffType ID = n_Z2m;
29
30extern omBin gmp_nrz_bin; /* init in rintegers*/
31
32void    nr2mCoeffWrite  (const coeffs r, BOOLEAN /*details*/)
33{
34  PrintS("//   coeff. ring is : ");
35  Print("Z/2^%lu\n", r->modExponent);
36}
37
38BOOLEAN nr2mCoeffIsEqual(const coeffs r, n_coeffType n, void * p)
39{
40  if (n==n_Z2m)
41  {
42    int m=(int)(long)(p);
43    unsigned long mm=r->mod2mMask;
44    if ((mm>>m)==1L) return TRUE;
45  }
46  return FALSE;
47}
48
49static char* nr2mCoeffString(const coeffs r)
50{
51  char* s = (char*) omAlloc(11+11);
52  sprintf(s,"integer,2,%lu",r->modExponent);
53  return s;
54}
55
56/* for initializing function pointers */
57BOOLEAN nr2mInitChar (coeffs r, void* p)
58{
59  assume( getCoeffType(r) == ID );
60  nr2mInitExp((int)(long)(p), r);
61  r->cfKillChar    = ndKillChar; /* dummy*/
62  r->nCoeffIsEqual = nr2mCoeffIsEqual;
63  r->cfCoeffString = nr2mCoeffString;
64
65  r->modBase = (int_number) omAllocBin (gmp_nrz_bin);
66  mpz_init_set_si (r->modBase, 2L);
67  r->modNumber= (int_number) omAllocBin (gmp_nrz_bin);
68  mpz_init (r->modNumber);
69  mpz_pow_ui (r->modNumber, r->modBase, r->modExponent);
70
71  /* next cast may yield an overflow as mod2mMask is an unsigned long */
72  r->ch = (int)r->mod2mMask + 1;
73
74  r->cfInit        = nr2mInit;
75  r->cfCopy        = ndCopy;
76  r->cfInt         = nr2mInt;
77  r->cfAdd         = nr2mAdd;
78  r->cfSub         = nr2mSub;
79  r->cfMult        = nr2mMult;
80  r->cfDiv         = nr2mDiv;
81  r->cfIntDiv      = nr2mIntDiv;
82  r->cfIntMod      = nr2mMod;
83  r->cfExactDiv    = nr2mDiv;
84  r->cfNeg         = nr2mNeg;
85  r->cfInvers      = nr2mInvers;
86  r->cfDivBy       = nr2mDivBy;
87  r->cfDivComp     = nr2mDivComp;
88  r->cfGreater     = nr2mGreater;
89  r->cfEqual       = nr2mEqual;
90  r->cfIsZero      = nr2mIsZero;
91  r->cfIsOne       = nr2mIsOne;
92  r->cfIsMOne      = nr2mIsMOne;
93  r->cfGreaterZero = nr2mGreaterZero;
94  r->cfWriteLong       = nr2mWrite;
95  r->cfRead        = nr2mRead;
96  r->cfPower       = nr2mPower;
97  r->cfSetMap      = nr2mSetMap;
98  r->cfNormalize   = ndNormalize;
99  r->cfLcm         = nr2mLcm;
100  r->cfGcd         = nr2mGcd;
101  r->cfIsUnit      = nr2mIsUnit;
102  r->cfGetUnit     = nr2mGetUnit;
103  r->cfExtGcd      = nr2mExtGcd;
104  r->cfName        = ndName;
105  r->cfCoeffWrite  = nr2mCoeffWrite;
106  r->cfInit_bigint = nr2mMapQ;
107#ifdef LDEBUG
108  r->cfDBTest      = nr2mDBTest;
109#endif
110  r->has_simple_Alloc=TRUE;
111  return FALSE;
112}
113
114/*
115 * Multiply two numbers
116 */
117number nr2mMult(number a, number b, const coeffs r)
118{
119  if (((NATNUMBER)a == 0) || ((NATNUMBER)b == 0))
120    return (number)0;
121  else
122    return nr2mMultM(a, b, r);
123}
124
125/*
126 * Give the smallest k, such that a * x = k = b * y has a solution
127 */
128number nr2mLcm(number a, number b, const coeffs)
129{
130  NATNUMBER res = 0;
131  if ((NATNUMBER)a == 0) a = (number) 1;
132  if ((NATNUMBER)b == 0) b = (number) 1;
133  while ((NATNUMBER)a % 2 == 0)
134  {
135    a = (number)((NATNUMBER)a / 2);
136    if ((NATNUMBER)b % 2 == 0) b = (number)((NATNUMBER)b / 2);
137    res++;
138  }
139  while ((NATNUMBER)b % 2 == 0)
140  {
141    b = (number)((NATNUMBER)b / 2);
142    res++;
143  }
144  return (number)(1L << res);  // (2**res)
145}
146
147/*
148 * Give the largest k, such that a = x * k, b = y * k has
149 * a solution.
150 */
151number nr2mGcd(number a, number b, const coeffs)
152{
153  NATNUMBER res = 0;
154  if ((NATNUMBER)a == 0 && (NATNUMBER)b == 0) return (number)1;
155  while ((NATNUMBER)a % 2 == 0 && (NATNUMBER)b % 2 == 0)
156  {
157    a = (number)((NATNUMBER)a / 2);
158    b = (number)((NATNUMBER)b / 2);
159    res++;
160  }
161//  if ((NATNUMBER)b % 2 == 0)
162//  {
163//    return (number)((1L << res)); // * (NATNUMBER) a);  // (2**res)*a    a is a unit
164//  }
165//  else
166//  {
167    return (number)((1L << res)); // * (NATNUMBER) b);  // (2**res)*b    b is a unit
168//  }
169}
170
171/*
172 * Give the largest k, such that a = x * k, b = y * k has
173 * a solution.
174 */
175number nr2mExtGcd(number a, number b, number *s, number *t, const coeffs r)
176{
177  NATNUMBER res = 0;
178  if ((NATNUMBER)a == 0 && (NATNUMBER)b == 0) return (number)1;
179  while ((NATNUMBER)a % 2 == 0 && (NATNUMBER)b % 2 == 0)
180  {
181    a = (number)((NATNUMBER)a / 2);
182    b = (number)((NATNUMBER)b / 2);
183    res++;
184  }
185  if ((NATNUMBER)b % 2 == 0)
186  {
187    *t = NULL;
188    *s = nr2mInvers(a,r);
189    return (number)((1L << res)); // * (NATNUMBER) a);  // (2**res)*a    a is a unit
190  }
191  else
192  {
193    *s = NULL;
194    *t = nr2mInvers(b,r);
195    return (number)((1L << res)); // * (NATNUMBER) b);  // (2**res)*b    b is a unit
196  }
197}
198
199void nr2mPower(number a, int i, number * result, const coeffs r)
200{
201  if (i == 0)
202  {
203    *(NATNUMBER *)result = 1;
204  }
205  else if (i == 1)
206  {
207    *result = a;
208  }
209  else
210  {
211    nr2mPower(a, i-1, result, r);
212    *result = nr2mMultM(a, *result, r);
213  }
214}
215
216/*
217 * create a number from int
218 */
219number nr2mInit(long i, const coeffs r)
220{
221  if (i == 0) return (number)(NATNUMBER)i;
222
223  long ii = i;
224  NATNUMBER j = (NATNUMBER)1;
225  if (ii < 0) { j = r->mod2mMask; ii = -ii; }
226  NATNUMBER k = (NATNUMBER)ii;
227  k = k & r->mod2mMask;
228  /* now we have: i = j * k mod 2^m */
229  return (number)nr2mMult((number)j, (number)k, r);
230}
231
232/*
233 * convert a number to an int in ]-k/2 .. k/2],
234 * where k = 2^m; i.e., an int in ]-2^(m-1) .. 2^(m-1)];
235 * note that the code computes a long which will then
236 * automatically casted to int
237 */
238static long nr2mLong(number &n, const coeffs r)
239{
240  NATNUMBER nn = (unsigned long)(NATNUMBER)n & r->mod2mMask;
241  unsigned long l = r->mod2mMask >> 1; l++; /* now: l = 2^(m-1) */
242  if ((NATNUMBER)nn > l)
243    return (long)((NATNUMBER)nn - r->mod2mMask - 1);
244  else
245    return (long)((NATNUMBER)nn);
246}
247int nr2mInt(number &n, const coeffs r)
248{
249  return (int)nr2mLong(n,r);
250}
251
252number nr2mAdd(number a, number b, const coeffs r)
253{
254  return nr2mAddM(a, b, r);
255}
256
257number nr2mSub(number a, number b, const coeffs r)
258{
259  return nr2mSubM(a, b, r);
260}
261
262BOOLEAN nr2mIsUnit(number a, const coeffs)
263{
264  return ((NATNUMBER)a % 2 == 1);
265}
266
267number nr2mGetUnit(number k, const coeffs)
268{
269  if (k == NULL) return (number)1;
270  NATNUMBER erg = (NATNUMBER)k;
271  while (erg % 2 == 0) erg = erg / 2;
272  return (number)erg;
273}
274
275BOOLEAN nr2mIsZero(number a, const coeffs)
276{
277  return 0 == (NATNUMBER)a;
278}
279
280BOOLEAN nr2mIsOne(number a, const coeffs)
281{
282  return 1 == (NATNUMBER)a;
283}
284
285BOOLEAN nr2mIsMOne(number a, const coeffs r)
286{
287  return (r->mod2mMask  == (NATNUMBER)a);
288}
289
290BOOLEAN nr2mEqual(number a, number b, const coeffs)
291{
292  return (a == b);
293}
294
295BOOLEAN nr2mGreater(number a, number b, const coeffs r)
296{
297  return nr2mDivBy(a, b,r);
298}
299
300/* Is 'a' divisible by 'b'? There are two cases:
301   1) a = 0 mod 2^m; then TRUE iff b = 0 or b is a power of 2
302   2) a, b <> 0; then TRUE iff b/gcd(a, b) is a unit mod 2^m */
303BOOLEAN nr2mDivBy (number a, number b, const coeffs r)
304{
305  if (a == NULL)
306  {
307    NATNUMBER c = r->mod2mMask + 1;
308    if (c != 0) /* i.e., if no overflow */
309      return (c % (NATNUMBER)b) == 0;
310    else
311    {
312      /* overflow: we need to check whether b
313         is zero or a power of 2: */
314      c = (NATNUMBER)b;
315      while (c != 0)
316      {
317        if ((c % 2) != 0) return FALSE;
318        c = c >> 1;
319      }
320      return TRUE;
321    }
322  }
323  else
324  {
325    number n = nr2mGcd(a, b, r);
326    n = nr2mDiv(b, n, r);
327    return nr2mIsUnit(n, r);
328  }
329}
330
331int nr2mDivComp(number as, number bs, const coeffs)
332{
333  NATNUMBER a = (NATNUMBER)as;
334  NATNUMBER b = (NATNUMBER)bs;
335  assume(a != 0 && b != 0);
336  while (a % 2 == 0 && b % 2 == 0)
337  {
338    a = a / 2;
339    b = b / 2;
340  }
341  if (a % 2 == 0)
342  {
343    return -1;
344  }
345  else
346  {
347    if (b % 2 == 1)
348    {
349      return 2;
350    }
351    else
352    {
353      return 1;
354    }
355  }
356}
357
358/* TRUE iff 0 < k <= 2^m / 2 */
359BOOLEAN nr2mGreaterZero(number k, const coeffs r)
360{
361  if ((NATNUMBER)k == 0) return FALSE;
362  if ((NATNUMBER)k > ((r->mod2mMask >> 1) + 1)) return FALSE;
363  return TRUE;
364}
365
366/* assumes that 'a' is odd, i.e., a unit in Z/2^m, and computes
367   the extended gcd of 'a' and 2^m, in order to find some 's'
368   and 't' such that a * s + 2^m * t = gcd(a, 2^m) = 1;
369   this code will always find a positive 's' */
370void specialXGCD(unsigned long& s, unsigned long a, const coeffs r)
371{
372  int_number u = (int_number)omAlloc(sizeof(mpz_t));
373  mpz_init_set_ui(u, a);
374  int_number u0 = (int_number)omAlloc(sizeof(mpz_t));
375  mpz_init(u0);
376  int_number u1 = (int_number)omAlloc(sizeof(mpz_t));
377  mpz_init_set_ui(u1, 1);
378  int_number u2 = (int_number)omAlloc(sizeof(mpz_t));
379  mpz_init(u2);
380  int_number v = (int_number)omAlloc(sizeof(mpz_t));
381  mpz_init_set_ui(v, r->mod2mMask);
382  mpz_add_ui(v, v, 1); /* now: v = 2^m */
383  int_number v0 = (int_number)omAlloc(sizeof(mpz_t));
384  mpz_init(v0);
385  int_number v1 = (int_number)omAlloc(sizeof(mpz_t));
386  mpz_init(v1);
387  int_number v2 = (int_number)omAlloc(sizeof(mpz_t));
388  mpz_init_set_ui(v2, 1);
389  int_number q = (int_number)omAlloc(sizeof(mpz_t));
390  mpz_init(q);
391  int_number rr = (int_number)omAlloc(sizeof(mpz_t));
392  mpz_init(rr);
393
394  while (mpz_cmp_ui(v, 0) != 0) /* i.e., while v != 0 */
395  {
396    mpz_div(q, u, v);
397    mpz_mod(rr, u, v);
398    mpz_set(u, v);
399    mpz_set(v, rr);
400    mpz_set(u0, u2);
401    mpz_set(v0, v2);
402    mpz_mul(u2, u2, q); mpz_sub(u2, u1, u2); /* u2 = u1 - q * u2 */
403    mpz_mul(v2, v2, q); mpz_sub(v2, v1, v2); /* v2 = v1 - q * v2 */
404    mpz_set(u1, u0);
405    mpz_set(v1, v0);
406  }
407
408  while (mpz_cmp_ui(u1, 0) < 0) /* i.e., while u1 < 0 */
409  {
410    /* we add 2^m = (2^m - 1) + 1 to u1: */
411    mpz_add_ui(u1, u1, r->mod2mMask);
412    mpz_add_ui(u1, u1, 1);
413  }
414  s = mpz_get_ui(u1); /* now: 0 <= s <= 2^m - 1 */
415
416  mpz_clear(u);  omFree((ADDRESS)u);
417  mpz_clear(u0); omFree((ADDRESS)u0);
418  mpz_clear(u1); omFree((ADDRESS)u1);
419  mpz_clear(u2); omFree((ADDRESS)u2);
420  mpz_clear(v);  omFree((ADDRESS)v);
421  mpz_clear(v0); omFree((ADDRESS)v0);
422  mpz_clear(v1); omFree((ADDRESS)v1);
423  mpz_clear(v2); omFree((ADDRESS)v2);
424  mpz_clear(q); omFree((ADDRESS)q);
425  mpz_clear(rr); omFree((ADDRESS)rr);
426}
427
428NATNUMBER InvMod(NATNUMBER a, const coeffs r)
429{
430  assume((NATNUMBER)a % 2 != 0);
431  unsigned long s;
432  specialXGCD(s, a, r);
433  return s;
434}
435//#endif
436
437inline number nr2mInversM(number c, const coeffs r)
438{
439  assume((NATNUMBER)c % 2 != 0);
440  // Table !!!
441  NATNUMBER inv;
442  inv = InvMod((NATNUMBER)c,r);
443  return (number)inv;
444}
445
446number nr2mDiv(number a, number b, const coeffs r)
447{
448  if ((NATNUMBER)a == 0) return (number)0;
449  else if ((NATNUMBER)b % 2 == 0)
450  {
451    if ((NATNUMBER)b != 0)
452    {
453      while (((NATNUMBER)b % 2 == 0) && ((NATNUMBER)a % 2 == 0))
454      {
455        a = (number)((NATNUMBER)a / 2);
456        b = (number)((NATNUMBER)b / 2);
457      }
458    }
459    if ((NATNUMBER)b % 2 == 0)
460    {
461      WerrorS("Division not possible, even by cancelling zero divisors.");
462      WerrorS("Result is integer division without remainder.");
463      return (number) ((NATNUMBER) a / (NATNUMBER) b);
464    }
465  }
466  return (number)nr2mMult(a, nr2mInversM(b,r),r);
467}
468
469number nr2mMod(number a, number b, const coeffs r)
470{
471  /*
472    We need to return the number rr which is uniquely determined by the
473    following two properties:
474      (1) 0 <= rr < |b| (with respect to '<' and '<=' performed in Z x Z)
475      (2) There exists some k in the integers Z such that a = k * b + rr.
476    Consider g := gcd(2^m, |b|). Note that then |b|/g is a unit in Z/2^m.
477    Now, there are three cases:
478      (a) g = 1
479          Then |b| is a unit in Z/2^m, i.e. |b| (and also b) divides a.
480          Thus rr = 0.
481      (b) g <> 1 and g divides a
482          Then a = (a/g) * (|b|/g)^(-1) * b (up to sign), i.e. again rr = 0.
483      (c) g <> 1 and g does not divide a
484          Let's denote the division with remainder of a by g as follows:
485          a = s * g + t. Then t = a - s * g = a - s * (|b|/g)^(-1) * |b|
486          fulfills (1) and (2), i.e. rr := t is the correct result. Hence
487          in this third case, rr is the remainder of division of a by g in Z.
488    This algorithm is the same as for the case Z/n, except that we may
489    compute the gcd of |b| and 2^m "by hand": We just extract the highest
490    power of 2 (<= 2^m) that is contained in b.
491  */
492  assume((NATNUMBER) b != 0);
493  NATNUMBER g = 1;
494  NATNUMBER b_div = (NATNUMBER) b;
495
496  /*
497   * b_div is unsigned, so that (b_div < 0) evaluates false at compile-time
498   *
499  if (b_div < 0) b_div = -b_div; // b_div now represents |b|, BUT b_div is unsigned!
500  */
501
502  NATNUMBER rr = 0;
503  while ((g < r->mod2mMask ) && (b_div > 0) && (b_div % 2 == 0))
504  {
505    b_div = b_div >> 1;
506    g = g << 1;
507  } // g is now the gcd of 2^m and |b|
508
509  if (g != 1) rr = (NATNUMBER)a % g;
510  return (number)rr;
511}
512
513number nr2mIntDiv(number a, number b, const coeffs r)
514{
515  if ((NATNUMBER)a == 0)
516  {
517    if ((NATNUMBER)b == 0)
518      return (number)1;
519    if ((NATNUMBER)b == 1)
520      return (number)0;
521    NATNUMBER c = r->mod2mMask + 1;
522    if (c != 0) /* i.e., if no overflow */
523      return (number)(c / (NATNUMBER)b);
524    else
525    {
526      /* overflow: c = 2^32 resp. 2^64, depending on platform */
527      int_number cc = (int_number)omAlloc(sizeof(mpz_t));
528      mpz_init_set_ui(cc, r->mod2mMask); mpz_add_ui(cc, cc, 1);
529      mpz_div_ui(cc, cc, (unsigned long)(NATNUMBER)b);
530      unsigned long s = mpz_get_ui(cc);
531      mpz_clear(cc); omFree((ADDRESS)cc);
532      return (number)(NATNUMBER)s;
533    }
534  }
535  else
536  {
537    if ((NATNUMBER)b == 0)
538      return (number)0;
539    return (number)((NATNUMBER) a / (NATNUMBER) b);
540  }
541}
542
543number nr2mInvers(number c, const coeffs r)
544{
545  if ((NATNUMBER)c % 2 == 0)
546  {
547    WerrorS("division by zero divisor");
548    return (number)0;
549  }
550  return nr2mInversM(c, r);
551}
552
553number nr2mNeg(number c, const coeffs r)
554{
555  if ((NATNUMBER)c == 0) return c;
556  return nr2mNegM(c, r);
557}
558
559number nr2mMapMachineInt(number from, const coeffs /*src*/, const coeffs dst)
560{
561  NATNUMBER i = ((NATNUMBER)from) % dst->mod2mMask ;
562  return (number)i;
563}
564
565number nr2mMapZp(number from, const coeffs /*src*/, const coeffs dst)
566{
567  NATNUMBER j = (NATNUMBER)1;
568  long ii = (long)from;
569  if (ii < 0) { j = dst->mod2mMask; ii = -ii; }
570  NATNUMBER i = (NATNUMBER)ii;
571  i = i & dst->mod2mMask;
572  /* now we have: from = j * i mod 2^m */
573  return (number)nr2mMult((number)i, (number)j, dst);
574}
575
576number nr2mMapQ(number from, const coeffs src, const coeffs dst)
577{
578  int_number erg = (int_number)omAllocBin(gmp_nrz_bin);
579  mpz_init(erg);
580  int_number k = (int_number)omAlloc(sizeof(mpz_t));
581  mpz_init_set_ui(k, dst->mod2mMask);
582
583  nlGMP(from, (number)erg, src);
584  mpz_and(erg, erg, k);
585  number res = (number)mpz_get_ui(erg);
586
587  mpz_clear(erg); omFree((ADDRESS)erg);
588  mpz_clear(k);   omFree((ADDRESS)k);
589
590  return (number)res;
591}
592
593number nr2mMapGMP(number from, const coeffs /*src*/, const coeffs dst)
594{
595  int_number erg = (int_number)omAllocBin(gmp_nrz_bin);
596  mpz_init(erg);
597  int_number k = (int_number)omAlloc(sizeof(mpz_t));
598  mpz_init_set_ui(k, dst->mod2mMask);
599
600  mpz_and(erg, (int_number)from, k);
601  number res = (number) mpz_get_ui(erg);
602
603  mpz_clear(erg); omFree((ADDRESS)erg);
604  mpz_clear(k);   omFree((ADDRESS)k);
605
606  return (number)res;
607}
608
609nMapFunc nr2mSetMap(const coeffs src, const coeffs dst)
610{
611  if (nCoeff_is_Ring_2toM(src)
612     && (src->mod2mMask == dst->mod2mMask))
613  {
614    return ndCopyMap;
615  }
616  if (nCoeff_is_Ring_2toM(src)
617     && (src->mod2mMask < dst->mod2mMask))
618  { /* i.e. map an integer mod 2^s into Z mod 2^t, where t < s */
619    return nr2mMapMachineInt;
620  }
621  if (nCoeff_is_Ring_2toM(src)
622     && (src->mod2mMask > dst->mod2mMask))
623  { /* i.e. map an integer mod 2^s into Z mod 2^t, where t > s */
624    // to be done
625  }
626  if (nCoeff_is_Ring_Z(src))
627  {
628    return nr2mMapGMP;
629  }
630  if (nCoeff_is_Q(src))
631  {
632    return nr2mMapQ;
633  }
634  if (nCoeff_is_Zp(src) && (src->ch == 2))
635  {
636    return nr2mMapZp;
637  }
638  if (nCoeff_is_Ring_PtoM(src) || nCoeff_is_Ring_ModN(src))
639  {
640    if (mpz_divisible_2exp_p(src->modNumber,dst->modExponent))
641      return nr2mMapGMP;
642  }
643  return NULL;      // default
644}
645
646/*
647 * set the exponent
648 */
649
650void nr2mSetExp(int m, coeffs r)
651{
652  if (m > 1)
653  {
654    /* we want mod2mMask to be the bit pattern
655       '111..1' consisting of m one's: */
656    r->modExponent= m;
657    r->mod2mMask = 1;
658    for (int i = 1; i < m; i++) r->mod2mMask = (r->mod2mMask << 1) + 1;
659  }
660  else
661  {
662    r->modExponent= 2;
663    /* code unexpectedly called with m = 1; we continue with m = 2: */
664    r->mod2mMask = 3; /* i.e., '11' in binary representation */
665  }
666}
667
668void nr2mInitExp(int m, coeffs r)
669{
670  nr2mSetExp(m, r);
671  if (m < 2)
672    WarnS("nr2mInitExp unexpectedly called with m = 1 (we continue with Z/2^2");
673}
674
675#ifdef LDEBUG
676BOOLEAN nr2mDBTest (number a, const char *, const int, const coeffs r)
677{
678  //if ((NATNUMBER)a < 0) return FALSE; // is unsigned!
679  if (((NATNUMBER)a & r->mod2mMask) != (NATNUMBER)a) return FALSE;
680  return TRUE;
681}
682#endif
683
684void nr2mWrite (number &a, const coeffs r)
685{
686  long i = nr2mLong(a, r);
687  StringAppend("%ld", i);
688}
689
690static const char* nr2mEati(const char *s, int *i, const coeffs r)
691{
692
693  if (((*s) >= '0') && ((*s) <= '9'))
694  {
695    (*i) = 0;
696    do
697    {
698      (*i) *= 10;
699      (*i) += *s++ - '0';
700      if ((*i) >= (MAX_INT_VAL / 10)) (*i) = (*i) & r->mod2mMask;
701    }
702    while (((*s) >= '0') && ((*s) <= '9'));
703    (*i) = (*i) & r->mod2mMask;
704  }
705  else (*i) = 1;
706  return s;
707}
708
709const char * nr2mRead (const char *s, number *a, const coeffs r)
710{
711  int z;
712  int n=1;
713
714  s = nr2mEati(s, &z,r);
715  if ((*s) == '/')
716  {
717    s++;
718    s = nr2mEati(s, &n,r);
719  }
720  if (n == 1)
721    *a = (number)(long)z;
722  else
723      *a = nr2mDiv((number)(long)z,(number)(long)n,r);
724  return s;
725}
726#endif
727/* #ifdef HAVE_RINGS */
Note: See TracBrowser for help on using the repository browser.