source: git/libpolys/coeffs/rintegers.cc @ fe99ed

spielwiese
Last change on this file since fe99ed was fe99ed, checked in by Claus Fieker <fieker@…>, 9 years ago
"fix" the SI==3 case in rintegers and add the missing bits in rmodulon again Try a short patch to keep it simple
  • Property mode set to 100644
File size: 36.6 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5* ABSTRACT: numbers (integers)
6*/
7
8
9
10#include <misc/auxiliary.h>
11
12#ifdef HAVE_RINGS
13
14#include <string.h>
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/rintegers.h>
23#include <coeffs/rmodulon.h>
24#include "si_gmp.h"
25
26/// Our Type!
27static const n_coeffType ID = n_Z;
28
29omBin gmp_nrz_bin = omGetSpecBin(sizeof(mpz_t));
30
31#if SI_INTEGER_VARIANT == 2
32/*
33 * Multiply two numbers
34 */
35number nrzMult (number a, number b, const coeffs)
36{
37  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
38  mpz_init(erg);
39  mpz_mul(erg, (int_number) a, (int_number) b);
40  return (number) erg;
41}
42
43/*
44 * Give the smallest non unit k, such that a * x = k = b * y has a solution
45 */
46number nrzLcm (number a,number b,const coeffs)
47{
48  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
49  mpz_init(erg);
50  mpz_lcm(erg, (int_number) a, (int_number) b);
51  return (number) erg;
52}
53
54/*
55 * Give the largest non unit k, such that a = x * k, b = y * k has
56 * a solution.
57 */
58number nrzGcd (number a,number b,const coeffs)
59{
60  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
61  mpz_init(erg);
62  mpz_gcd(erg, (int_number) a, (int_number) b);
63  return (number) erg;
64}
65
66/*
67 * Give the largest non unit k, such that a = x * k, b = y * k has
68 * a solution and r, s, s.t. k = s*a + t*b
69 */
70number  nrzExtGcd (number a, number b, number *s, number *t, const coeffs)
71{
72  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
73  int_number bs = (int_number) omAllocBin(gmp_nrz_bin);
74  int_number bt = (int_number) omAllocBin(gmp_nrz_bin);
75  mpz_init(erg);
76  mpz_init(bs);
77  mpz_init(bt);
78  mpz_gcdext(erg, bs, bt, (int_number) a, (int_number) b);
79  *s = (number) bs;
80  *t = (number) bt;
81  return (number) erg;
82}
83
84void nrzPower (number a, int i, number * result, const coeffs)
85{
86  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
87  mpz_init(erg);
88  mpz_pow_ui(erg, (int_number) a, i);
89  *result = (number) erg;
90}
91
92/*
93 * create a number from int
94 */
95number nrzInit (long i, const coeffs)
96{
97  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
98  mpz_init_set_si(erg, i);
99  return (number) erg;
100}
101
102void nrzDelete(number *a, const coeffs)
103{
104  if (*a == NULL) return;
105  mpz_clear((int_number) *a);
106  omFreeBin((ADDRESS) *a, gmp_nrz_bin);
107  *a = NULL;
108}
109
110number nrzCopy(number a, const coeffs)
111{
112  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
113  mpz_init_set(erg, (int_number) a);
114  return (number) erg;
115}
116
117#if 0
118number nrzCopyMap(number a, const coeffs /*src*/, const coeffs dst)
119{
120  return nrzCopy(a,dst);
121}
122#endif
123
124int nrzSize(number a, const coeffs)
125{
126  if (a == NULL) return 0;
127  return sizeof(mpz_t);
128}
129
130/*
131 * convert a number to int
132 */
133int nrzInt(number &n, const coeffs)
134{
135  return (int) mpz_get_si( (int_number)n);
136}
137
138number nrzAdd (number a, number b, const coeffs)
139{
140  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
141  mpz_init(erg);
142  mpz_add(erg, (int_number) a, (int_number) b);
143  return (number) erg;
144}
145
146number nrzSub (number a, number b, const coeffs)
147{
148  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
149  mpz_init(erg);
150  mpz_sub(erg, (int_number) a, (int_number) b);
151  return (number) erg;
152}
153
154number  nrzGetUnit (number, const coeffs r)
155{
156  return nrzInit(1, r);
157}
158
159BOOLEAN nrzIsUnit (number a, const coeffs)
160{
161  return 0 == mpz_cmpabs_ui((int_number) a, 1);
162}
163
164BOOLEAN nrzIsZero (number  a, const coeffs)
165{
166  return 0 == mpz_cmpabs_ui((int_number) a, 0);
167}
168
169BOOLEAN nrzIsOne (number a, const coeffs)
170{
171  return (a!=NULL) && (0 == mpz_cmp_si((int_number) a, 1));
172}
173
174BOOLEAN nrzIsMOne (number a, const coeffs)
175{
176  return (a!=NULL) && (0 == mpz_cmp_si((int_number) a, -1));
177}
178
179BOOLEAN nrzEqual (number a,number b, const coeffs)
180{
181  return 0 == mpz_cmp((int_number) a, (int_number) b);
182}
183
184BOOLEAN nrzGreater (number a,number b, const coeffs)
185{
186  return 0 < mpz_cmp((int_number) a, (int_number) b);
187}
188
189BOOLEAN nrzGreaterZero (number k, const coeffs)
190{
191  return 0 < mpz_cmp_si((int_number) k, 0);
192}
193
194int nrzDivComp(number a, number b, const coeffs r)
195{
196  if (nrzDivBy(a, b, r))
197  {
198    if (nrzDivBy(b, a, r)) return 2;
199    return -1;
200  }
201  if (nrzDivBy(b, a, r)) return 1;
202  return 0;
203}
204
205BOOLEAN nrzDivBy (number a,number b, const coeffs)
206{
207  return mpz_divisible_p((int_number) a, (int_number) b) != 0;
208}
209
210number nrzDiv (number a,number b, const coeffs R)
211{
212  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
213  mpz_init(erg);
214  int_number r = (int_number) omAllocBin(gmp_nrz_bin);
215  mpz_init(r);
216  mpz_tdiv_qr(erg, r, (int_number) a, (int_number) b);
217  //if (!nrzIsZero((number) r, R))
218  //{
219  //  WerrorS("Division by non divisible element.");
220  //  WerrorS("Result is without remainder.");
221  //}
222  mpz_clear(r);
223  omFreeBin(r, gmp_nrz_bin);
224  return (number) erg;
225}
226
227number nrzExactDiv (number a,number b, const coeffs)
228{
229  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
230  mpz_init(erg);
231  mpz_tdiv_q(erg, (int_number) a, (int_number) b);
232  return (number) erg;
233}
234
235number nrzIntMod (number a,number b, const coeffs)
236{
237  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
238  mpz_init(erg);
239  int_number r = (int_number) omAllocBin(gmp_nrz_bin);
240  mpz_init(r);
241  mpz_tdiv_qr(erg, r, (int_number) a, (int_number) b);
242  mpz_clear(erg);
243  return (number) r;
244}
245
246number  nrzInvers (number c, const coeffs r)
247{
248  if (!nrzIsUnit((number) c, r))
249  {
250    WerrorS("Non invertible element.");
251    return (number)0; //TODO
252  }
253  return nrzCopy(c,r);
254}
255
256number nrzNeg (number c, const coeffs)
257{
258// nNeg inplace !!!
259  mpz_mul_si((int_number) c, (int_number) c, -1);
260  return c;
261}
262
263number nrzMapMachineInt(number from, const coeffs /*src*/, const coeffs /*dst*/)
264{
265  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
266  mpz_init_set_ui(erg, (NATNUMBER) from);
267  return (number) erg;
268}
269
270number nrzMapZp(number from, const coeffs /*src*/, const coeffs /*dst*/)
271{
272  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
273  mpz_init_set_si(erg, (long) from);
274  return (number) erg;
275}
276
277number nrzMapQ(number from, const coeffs src, const coeffs /*dst*/)
278{
279  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
280  mpz_init(erg);
281  nlGMP(from, (number) erg, src);
282  return (number) erg;
283}
284
285nMapFunc nrzSetMap(const coeffs src, const coeffs /*dst*/)
286{
287  /* dst = currRing */
288  if (nCoeff_is_Ring_Z(src) || nCoeff_is_Ring_ModN(src) || nCoeff_is_Ring_PtoM(src))
289  {
290    return ndCopyMap; //nrzCopyMap;
291  }
292  if (nCoeff_is_Ring_2toM(src))
293  {
294    return nrzMapMachineInt;
295  }
296  if (nCoeff_is_Zp(src))
297  {
298    return nrzMapZp;
299  }
300  if (getCoeffType(src)==n_Q /*nCoeff_is_Q(src) or coeffs_BIGINT*/)
301  {
302    return nrzMapQ;
303  }
304  return NULL;      // default
305}
306
307
308/*
309 * set the exponent (allocate and init tables) (TODO)
310 */
311
312void nrzSetExp(int, coeffs)
313{
314}
315
316void nrzInitExp(int, coeffs)
317{
318}
319
320#ifdef LDEBUG
321BOOLEAN nrzDBTest (number, const char *, const int, const coeffs)
322{
323  return TRUE;//TODO
324}
325#endif
326
327void nrzWrite (number &a, const coeffs)
328{
329  char *s,*z;
330  if (a==NULL)
331  {
332    StringAppendS("o");
333  }
334  else
335  {
336    int l=mpz_sizeinbase((int_number) a, 10) + 2;
337    s=(char*)omAlloc(l);
338    z=mpz_get_str(s,10,(int_number) a);
339    StringAppendS(z);
340    omFreeSize((ADDRESS)s,l);
341  }
342}
343
344/*2
345* extracts a long integer from s, returns the rest    (COPY FROM longrat0.cc)
346*/
347static const char * nlEatLongC(char *s, mpz_ptr i)
348{
349  const char * start=s;
350
351  if (*s<'0' || *s>'9')
352  {
353    mpz_set_si(i,1);
354    return s;
355  }
356  while (*s >= '0' && *s <= '9') s++;
357  if (*s=='\0')
358  {
359    mpz_set_str(i,start,10);
360  }
361  else
362  {
363    char c=*s;
364    *s='\0';
365    mpz_set_str(i,start,10);
366    *s=c;
367  }
368  return s;
369}
370
371const char * nrzRead (const char *s, number *a, const coeffs)
372{
373  int_number z = (int_number) omAllocBin(gmp_nrz_bin);
374  {
375    mpz_init(z);
376    s = nlEatLongC((char *) s, z);
377  }
378  *a = (number) z;
379  return s;
380}
381
382void    nrzCoeffWrite  (const coeffs, BOOLEAN /*details*/)
383{
384  PrintS("//   coeff. ring is : Integers\n");
385}
386
387static char* nrzCoeffString(const coeffs)
388{
389  return omStrDup("integer");
390}
391
392coeffs nrzQuot1(number c, const coeffs r)
393{
394    int ch = r->cfInt(c, r);
395    int_number dummy;
396    dummy = (int_number) omAlloc(sizeof(mpz_t));
397    mpz_init_set_ui(dummy, ch);
398    ZnmInfo info;
399    info.base = dummy;
400    info.exp = (unsigned long) 1;
401    coeffs rr = nInitChar(n_Zn, (void*)&info);
402    return(rr);
403}
404
405BOOLEAN nrzInitChar(coeffs r,  void *)
406{
407  assume( getCoeffType(r) == ID );
408
409  r->is_field=FALSE;
410  r->is_domain=TRUE;
411  r->rep=n_rep_gmp;
412
413  r->nCoeffIsEqual = ndCoeffIsEqual;
414  r->cfCoeffString = nrzCoeffString;
415  r->cfKillChar = ndKillChar;
416  r->cfMult  = nrzMult;
417  r->cfSub   = nrzSub;
418  r->cfAdd   = nrzAdd;
419  r->cfDiv   = nrzDiv;
420  r->cfIntMod= nrzIntMod;
421  r->cfExactDiv= nrzExactDiv;
422  r->cfInit = nrzInit;
423  r->cfSize  = nrzSize;
424  r->cfInt  = nrzInt;
425  //#ifdef HAVE_RINGS
426  r->cfDivComp = nrzDivComp; // only for ring stuff
427  r->cfIsUnit = nrzIsUnit; // only for ring stuff
428  r->cfGetUnit = nrzGetUnit; // only for ring stuff
429  r->cfExtGcd = nrzExtGcd; // only for ring stuff
430  r->cfDivBy = nrzDivBy; // only for ring stuff
431  //#endif
432  r->cfInpNeg   = nrzNeg;
433  r->cfInvers= nrzInvers;
434  r->cfCopy  = nrzCopy;
435  r->cfWriteLong = nrzWrite;
436  r->cfRead = nrzRead;
437  r->cfGreater = nrzGreater;
438  r->cfEqual = nrzEqual;
439  r->cfIsZero = nrzIsZero;
440  r->cfIsOne = nrzIsOne;
441  r->cfIsMOne = nrzIsMOne;
442  r->cfGreaterZero = nrzGreaterZero;
443  r->cfPower = nrzPower;
444  r->cfGcd  = nrzGcd;
445  r->cfLcm  = nrzLcm;
446  r->cfDelete= nrzDelete;
447  r->cfSetMap = nrzSetMap;
448  r->cfCoeffWrite = nrzCoeffWrite;
449  r->cfQuot1 = nrzQuot1;
450  // debug stuff
451
452#ifdef LDEBUG
453  r->cfDBTest=nrzDBTest;
454#endif
455
456  r->nNULL = 0;
457  r->ch = 0;
458  r->has_simple_Alloc=FALSE;
459  r->has_simple_Inverse=FALSE;
460  return FALSE;
461}
462
463#elif SI_INTEGER_VARIANT == 3
464
465//make sure that a small number is an immediate integer
466//bascially coped from longrat.cc nlShort3
467//TODO: is there any point in checking 0 first???
468//TODO: it is not clear that this works in 32/64 bit everywhere.
469//      too many hacks.
470#ifdef LDEBUG
471#define nrzTest(A) nrzDBTest(A,__FILE__,__LINE__,NULL)
472BOOLEAN nrzDBTest (number x, const char *f, const int l, const coeffs);
473#else
474#define nrzTest(A)
475#endif
476
477#define CF_DEBUG 0
478static inline number nrz_short(number x)
479{
480#if CF_DEBUG
481  StringAppendS("short(");
482  nrzWrite(x, NULL);
483#endif
484  if (mpz_cmp_ui((int_number) x,(long)0)==0)
485  {
486    mpz_clear((int_number)x);
487    omFreeBin(x, gmp_nrz_bin);
488#if CF_DEBUG
489    StringAppendS(")=0");
490#endif
491    return INT_TO_SR(0);
492  }
493  if (mpz_size1((int_number)x)<=MP_SMALL)
494  {
495    int ui=mpz_get_si((int_number)x);
496    if ((((ui<<3)>>3)==ui)
497    && (mpz_cmp_si((int_number)x,(long)ui)==0))
498    {
499      mpz_clear((int_number)x);
500      omFreeBin(x, gmp_nrz_bin);
501#if CF_DEBUG
502    StringAppendS(")=imm");
503#endif
504      return INT_TO_SR(ui);
505    }
506  }
507#if CF_DEBUG
508  StringAppendS(")");
509#endif
510  return x;
511}
512
513
514
515
516/*
517 * Multiply two numbers
518 * check for 0, 1, -1 maybe
519 */
520#if CF_DEBUG
521number _nrzMult(number, number, const coeffs);
522number nrzMult(number a, number b, const coeffs R)
523{
524  StringSetS("Mult: ");
525  nrzWrite(a, R);
526  StringAppendS(" by ");
527  nrzWrite(b, R);
528  number c = _nrzMult(a, b, R);
529  StringAppendS(" is ");
530  nrzWrite(c, R);
531  char * s = StringEndS();
532  Print("%s\n", s);
533  omFree(s);
534  return c;
535}
536number _nrzMult (number a, number b, const coeffs R)
537#else
538number nrzMult (number a, number b, const coeffs R)
539#endif
540{
541  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b)) {
542  //from longrat.cc
543    if (SR_TO_INT(a)==0)
544      return a;
545    if (SR_TO_INT(b)==0)
546      return b;
547    long r=(long)((unsigned long)(SR_HDL(a)-1L))*((unsigned long)(SR_HDL(b)>>1));
548    if ((r/(SR_HDL(b)>>1))==(SR_HDL(a)-1L))
549    {
550      number u=((number) ((r>>1)+SR_INT));
551    //  if (((((long)SR_HDL(u))<<1)>>1)==SR_HDL(u)) return (u);
552      return nrzInit(SR_HDL(u)>>2, R);
553    }
554    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
555    mpz_init(erg);
556    mpz_set_si(erg, SR_TO_INT(a));
557    mpz_mul_si(erg, erg, SR_TO_INT(b));
558    nrzTest((number)erg);
559    return (number) erg;
560  }
561  else if (n_Z_IS_SMALL(a))
562  {
563    if (SR_TO_INT(a)==0)
564      return a;
565    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
566    mpz_init_set(erg, (int_number) b);
567    mpz_mul_si(erg, erg, SR_TO_INT(a));
568    nrzTest((number)erg);
569    return (number) erg;
570  }
571  else if (n_Z_IS_SMALL(b))
572  {
573    if (SR_TO_INT(b)==0)
574      return b;
575    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
576    mpz_init_set(erg, (int_number) a);
577    mpz_mul_si(erg, erg, SR_TO_INT(b));
578    nrzTest((number)erg);
579    return (number) erg;
580  }
581  else
582  {
583    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
584    mpz_init(erg);
585    mpz_mul(erg, (int_number) a, (int_number) b);
586    nrzTest((number)erg);
587    return (number) erg;
588  }
589}
590
591
592static int int_gcd(int a, int b)
593{
594  int r;
595  a = ABS(a);
596  b = ABS(b);
597  if (!a) return b;
598  if (!b) return a;
599  do
600  {
601    r = a % b;
602    a = b;
603    b = r;
604  } while (b);
605  return ABS(a); // % in c doeas not imply a signn
606                 // it would be unlikely to see a negative here
607                 // but who knows
608}
609
610/*
611 * Give the smallest non unit k, such that a * x = k = b * y has a solution
612 */
613number nrzLcm (number a, number b, const coeffs R)
614{
615  PrintS("nrzLcm\n");
616  int_number erg;
617  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
618  {
619    int g = int_gcd(SR_TO_INT(a), SR_TO_INT(b));
620    return nrzMult(a, INT_TO_SR(SR_TO_INT(b)/g), R);
621  }
622  else if (n_Z_IS_SMALL(a))
623  {
624    erg = (int_number) omAllocBin(gmp_nrz_bin);
625    mpz_init_set(erg, (int_number) b);
626    unsigned long g = mpz_gcd_ui(NULL, erg, (unsigned long) ABS(SR_TO_INT(a)));
627    mpz_mul_si(erg, erg, SR_TO_INT(a)/g);
628  }
629  else if (n_Z_IS_SMALL(b))
630  {
631    erg = (int_number) omAllocBin(gmp_nrz_bin);
632    mpz_init_set(erg, (int_number) a);
633    unsigned long g = mpz_gcd_ui(NULL, erg, (unsigned long) ABS(SR_TO_INT(b)));
634    mpz_mul_si(erg, erg, SR_TO_INT(b)/g);
635  }
636  else
637  {
638    erg = (int_number) omAllocBin(gmp_nrz_bin);
639    mpz_init(erg);
640    mpz_lcm(erg, (int_number) a, (int_number) b);
641  }
642  return (number) erg;
643}
644
645/*
646 * Give the largest non unit k, such that a = x * k, b = y * k has
647 * a solution.
648 */
649number nrzGcd (number a,number b,const coeffs R)
650{
651  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
652  {
653    int g = int_gcd(SR_TO_INT(a), SR_TO_INT(b));
654    return INT_TO_SR(g);
655  }
656  else if (n_Z_IS_SMALL(a))
657  {
658    if (a==INT_TO_SR(0))
659      return nrzCopy(b, R);
660    unsigned long g = mpz_gcd_ui(NULL, (int_number)b, (unsigned long) ABS(SR_TO_INT(a)));
661    return INT_TO_SR( g);
662  }
663  else if (n_Z_IS_SMALL(b))
664  {
665    if (b==INT_TO_SR(0))
666      return nrzCopy(a, R);
667    unsigned long g = mpz_gcd_ui(NULL, (int_number)a, (unsigned long) ABS(SR_TO_INT(b)));
668    return INT_TO_SR(g);
669  }
670  else
671  {
672    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
673    mpz_init(erg);
674    mpz_gcd(erg, (int_number) a, (int_number) b);
675    return (number) erg;
676  }
677}
678
679/*
680 * Give the largest non unit k, such that a = x * k, b = y * k has
681 * a solution and r, s, s.t. k = s*a + t*b
682 */
683static int int_extgcd(int a, int b, int * u, int* x, int * v, int* y)
684{
685  int q, r;
686  if (!a)
687  {
688    *u = 0;
689    *v = 1;
690    *x = -1;
691    *y = 0;
692    return b;
693  }
694  if (!b)
695  {
696    *u = 1;
697    *v = 0;
698    *x = 0;
699    *y = 1;
700    return a;
701  }
702  *u=1;
703  *v=0;
704  *x=0;
705  *y=1;
706  do
707  {
708    q = a/b;
709    r = a%b;
710    assume (q*b+r == a);
711    a = b;
712    b = r;
713
714    r = -(*v)*q+(*u);
715    (*u) =(*v);
716    (*v) = r;
717
718    r = -(*y)*q+(*x);
719    (*x) = (*y);
720    (*y) = r;
721  } while (b);
722
723  return a;
724}
725
726number  nrzExtGcd (number a, number b, number *s, number *t, const coeffs)
727{
728  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
729  {
730    int u, v, x, y;
731    int g = int_extgcd(SR_TO_INT(a), SR_TO_INT(b), &u, &v, &x, &y);
732    *s = INT_TO_SR(u);
733    *t = INT_TO_SR(v);
734    return INT_TO_SR(g);
735  }
736  else
737  {
738    mpz_t aa, bb;
739    if (n_Z_IS_SMALL(a))
740    {
741      mpz_init_set_si(aa, SR_TO_INT(a));
742    }
743    else
744    {
745      mpz_init_set(aa, (int_number) a);
746    }
747    if (n_Z_IS_SMALL(b))
748    {
749      mpz_init_set_si(bb, SR_TO_INT(b));
750    }
751    else
752    {
753      mpz_init_set(bb, (int_number) b);
754    }
755    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
756    int_number bs = (int_number) omAllocBin(gmp_nrz_bin);
757    int_number bt = (int_number) omAllocBin(gmp_nrz_bin);
758    mpz_init(erg);
759    mpz_init(bs);
760    mpz_init(bt);
761    mpz_gcdext(erg, bs, bt, aa, bb);
762    *s = nrz_short((number) bs);
763    *t = nrz_short((number) bt);
764    mpz_clear(aa);
765    mpz_clear(bb);
766    return nrz_short((number) erg);
767  }
768}
769#if CF_DEBUG
770number _nrzXExtGcd(number, number, number *, number *, number *, number *, const coeffs);
771number nrzXExtGcd(number a, number b, number *x, number * y, number * u, number * v, const coeffs R)
772{
773  char * s;
774  StringSetS("XExtGcd: ");
775  nrzWrite(a, R);
776  StringAppendS(" by ");
777  nrzWrite(b, R);
778  number c = _nrzXExtGcd(a, b, x, y, u, v, R);
779  StringAppendS(" is ");
780  nrzWrite(c, R);
781  StringAppendS("[[");
782  nrzWrite(*x, R);
783  StringAppendS(", ");
784  nrzWrite(*y, R);
785  StringAppendS("], ");
786  nrzWrite(*u, R);
787  StringAppendS(", ");
788  nrzWrite(*v, R);
789  s=StringEndS();
790  Print("%s]]\n", s);
791  omFree(s);
792  return c;
793}
794number  _nrzXExtGcd (number a, number b, number *s, number *t, number *u, number *v, const coeffs )
795#else
796number  nrzXExtGcd (number a, number b, number *s, number *t, number *u, number *v, const coeffs )
797#endif
798{
799  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
800  {
801    int uu, vv, x, y;
802    int g = int_extgcd(SR_TO_INT(a), SR_TO_INT(b), &uu, &vv, &x, &y);
803    *s = INT_TO_SR(uu);
804    *t = INT_TO_SR(vv);
805    *u = INT_TO_SR(x);
806    *v = INT_TO_SR(y);
807    return INT_TO_SR(g);
808  }
809  else
810  {
811    mpz_t aa, bb;
812    if (n_Z_IS_SMALL(a))
813    {
814      mpz_init_set_si(aa, SR_TO_INT(a));
815    }
816    else
817    {
818      mpz_init_set(aa, (int_number) a);
819    }
820    if (n_Z_IS_SMALL(b))
821    {
822      mpz_init_set_si(bb, SR_TO_INT(b));
823    }
824    else
825    {
826      mpz_init_set(bb, (int_number) b);
827    }
828    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
829    int_number bs = (int_number) omAllocBin(gmp_nrz_bin);
830    int_number bt = (int_number) omAllocBin(gmp_nrz_bin);
831    mpz_init(erg);
832    mpz_init(bs);
833    mpz_init(bt);
834
835    mpz_gcdext(erg, bs, bt, aa, bb);
836
837    int_number bu = (int_number) omAllocBin(gmp_nrz_bin);
838    int_number bv = (int_number) omAllocBin(gmp_nrz_bin);
839
840    mpz_init_set(bu, (int_number) bb);
841    mpz_init_set(bv, (int_number) aa);
842
843    mpz_clear(aa);
844    mpz_clear(bb);
845    assume(mpz_cmp_si(erg, 0));
846
847    mpz_div(bu, bu, erg);
848    mpz_div(bv, bv, erg);
849
850    mpz_mul_si(bu, bu, -1);
851    *u = nrz_short((number) bu);
852    *v = nrz_short((number) bv);
853
854    *s = nrz_short((number) bs);
855    *t = nrz_short((number) bt);
856    return nrz_short((number) erg);
857  }
858}
859#if CF_DEBUG
860number _nrzQuotRem(number, number, number *, const coeffs);
861number nrzQuotRem(number a, number b, number * r, const coeffs R)
862{
863  StringSetS("QuotRem: ");
864  nrzWrite(a, R);
865  StringAppendS(" by ");
866  nrzWrite(b, R);
867  number c = _nrzQuotRem(a, b, r, R);
868  StringAppendS(" is ");
869  nrzWrite(c, R);
870  if (r) {
871    StringAppendS("+R(");
872    nrzWrite(*r, R);
873    StringAppendS(")");
874  }
875  char * s = StringEndS();
876  Print("%s\n", s);
877  omFree(s);
878  return c;
879}
880number _nrzQuotRem (number a, number b, number * r, const coeffs )
881#else
882number nrzQuotRem (number a, number b, number * r, const coeffs )
883#endif
884{
885  assume(SR_TO_INT(b));
886  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
887  {
888    if (r)
889      *r = INT_TO_SR(SR_TO_INT(a) % SR_TO_INT(b));
890    return INT_TO_SR(SR_TO_INT(a)/SR_TO_INT(b));
891  }
892  else if (n_Z_IS_SMALL(a))
893  {
894    //a is small, b is not, so q=0, r=a
895    if (r)
896      *r = a;
897    return INT_TO_SR(0);
898  }
899  else if (n_Z_IS_SMALL(b))
900  {
901    unsigned long rr;
902    int_number qq = (int_number) omAllocBin(gmp_nrz_bin);
903    mpz_init(qq);
904    mpz_t rrr;
905    mpz_init(rrr);
906    rr = mpz_divmod_ui(qq, rrr, (int_number) a, (unsigned long)ABS(SR_TO_INT(b)));
907    mpz_clear(rrr);
908
909    if (r)
910      *r = INT_TO_SR(rr);
911    if (SR_TO_INT(b)<0)
912    {
913      mpz_mul_si(qq, qq, -1);
914    }
915    return nrz_short((number)qq);
916  }
917  int_number qq = (int_number) omAllocBin(gmp_nrz_bin),
918             rr = (int_number) omAllocBin(gmp_nrz_bin);
919  mpz_init(qq);
920  mpz_init(rr);
921  mpz_divmod(qq, rr, (int_number)a, (int_number)b);
922  if (r)
923    *r = (number) rr;
924  else
925  {
926    mpz_clear(rr);
927  }
928  nrzTest((number)qq);
929  return (number) qq;
930}
931
932
933void nrzPower (number a, int i, number * result, const coeffs)
934{
935  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
936  mpz_init(erg);
937  mpz_t aa;
938  if (n_Z_IS_SMALL(a))
939    mpz_init_set_si(aa, SR_TO_INT(a));
940  else
941    mpz_init_set(aa, (int_number) a);
942  mpz_pow_ui(erg, aa, i);
943  *result = nrz_short((number) erg);
944}
945
946/*
947 * create a number from int
948 * TODO: do not create an mpz if not necessary
949 */
950number nrzInit (long i, const coeffs)
951{
952  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
953  mpz_init_set_si(erg, i);
954  return nrz_short((number) erg);
955}
956
957void nrzDelete(number *a, const coeffs)
958{
959  if (*a == NULL) return;
960  if (n_Z_IS_SMALL(*a)==0)
961  {
962    mpz_clear((int_number) *a);
963    omFreeBin((ADDRESS) *a, gmp_nrz_bin);
964  }
965  *a = NULL;
966}
967
968number nrzCopy(number a, const coeffs)
969{
970  if (n_Z_IS_SMALL(a)) return a;
971  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
972  mpz_init_set(erg, (int_number) a);
973  return (number) erg;
974}
975
976int nrzSize(number a, const coeffs)
977{
978  if (a == NULL) return 0;
979  if (n_Z_IS_SMALL(a)) return 1;
980  return mpz_size1((int_number)a)+1;
981}
982
983/*
984 * convert a number to int
985 */
986int nrzInt(number &n, const coeffs)
987{
988  if (n_Z_IS_SMALL(n)) return SR_TO_INT(n);
989  return (int) mpz_get_si( (int_number)n);
990}
991#if CF_DEBUG
992number _nrzAdd(number, number, const coeffs);
993number nrzAdd(number a, number b, const coeffs R)
994{
995  StringSetS("Add: ");
996  nrzWrite(a, R);
997  StringAppendS(" to ");
998  nrzWrite(b, R);
999  number c = _nrzAdd(a, b, R);
1000  StringAppendS(" is ");
1001  nrzWrite(c, R);
1002  char * s = StringEndS();
1003  Print("%s\n", s);
1004  omFree(s);
1005  return c;
1006}
1007number _nrzAdd (number a, number b, const coeffs )
1008#else
1009number nrzAdd (number a, number b, const coeffs )
1010#endif
1011{
1012  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
1013  {
1014    int c = SR_TO_INT(a) + SR_TO_INT(b);
1015    if (INT_IS_SMALL(c))
1016      return INT_TO_SR(c);
1017    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1018    mpz_init_set_si(erg, c);
1019
1020    nrzTest((number)erg);
1021    return (number) erg;
1022  }
1023  else if (n_Z_IS_SMALL(a))
1024  {
1025    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1026    mpz_init(erg);
1027    if (SR_TO_INT(a)>0)
1028      mpz_add_ui(erg, (int_number) b, (unsigned long)SR_TO_INT(a));
1029    else
1030      mpz_sub_ui(erg, (int_number) b, (unsigned long)-(SR_TO_INT(a)));
1031    return nrz_short((number) erg);
1032  }
1033  else if (n_Z_IS_SMALL(b))
1034  {
1035    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1036    mpz_init(erg);
1037    if (SR_TO_INT(b)>0)
1038      mpz_add_ui(erg, (int_number) a, (unsigned long)SR_TO_INT(b));
1039    else
1040      mpz_sub_ui(erg, (int_number) a, (unsigned long)-(SR_TO_INT(b)));
1041    return nrz_short((number) erg);
1042  }
1043  else
1044  {
1045    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1046    mpz_init(erg);
1047    mpz_add(erg, (int_number) a, (int_number) b);
1048    return nrz_short((number) erg);
1049  }
1050}
1051
1052number nrzSub (number a, number b,  const coeffs )
1053{
1054  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
1055  {
1056    int c = SR_TO_INT(a) - SR_TO_INT(b);
1057    if (INT_IS_SMALL(c))
1058      return INT_TO_SR(c);
1059    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1060    mpz_init_set_si(erg, c);
1061    nrzTest((number)erg);
1062    return (number) erg;
1063  }
1064  else if (n_Z_IS_SMALL(a))
1065  {
1066    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1067    mpz_init(erg);
1068
1069    if (SR_TO_INT(a)>0)
1070      mpz_ui_sub(erg, (unsigned long)SR_TO_INT(a), (int_number) b);
1071    else
1072    {
1073      mpz_add_ui(erg, (int_number) b, (unsigned long)-SR_TO_INT(a));
1074      mpz_neg(erg, erg);
1075    }
1076    return nrz_short((number) erg);
1077  }
1078  else if (n_Z_IS_SMALL(b))
1079  {
1080    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1081    mpz_init(erg);
1082    if (SR_TO_INT(b)>0)
1083      mpz_sub_ui(erg, (int_number) a, (unsigned long)SR_TO_INT(b));
1084    else
1085      mpz_add_ui(erg, (int_number) a, (unsigned long)-SR_TO_INT(b));
1086    return nrz_short((number) erg);
1087  }
1088  else
1089  {
1090    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1091    mpz_init(erg);
1092    mpz_sub(erg, (int_number) a, (int_number) b);
1093    return nrz_short((number) erg);
1094  }
1095}
1096
1097number  nrzGetUnit (number n, const coeffs r)
1098{
1099  if (nrzGreaterZero(n, r))
1100    return INT_TO_SR(1);
1101  else
1102    return INT_TO_SR(-1);
1103}
1104
1105number nrzAnn(number n, const coeffs)
1106{
1107  if (SR_TO_INT(n))  // in Z: the annihilator of !=0 is 0
1108    return INT_TO_SR(0);
1109  else
1110    return INT_TO_SR(1);
1111}
1112
1113BOOLEAN nrzIsUnit (number a, const coeffs)
1114{
1115  return ABS(SR_TO_INT(a))==1;
1116}
1117
1118BOOLEAN nrzIsZero (number  a, const coeffs)
1119{
1120  return a==INT_TO_SR(0);
1121}
1122
1123BOOLEAN nrzIsOne (number a, const coeffs)
1124{
1125  return a==INT_TO_SR(1);
1126}
1127
1128BOOLEAN nrzIsMOne (number a, const coeffs)
1129{
1130  return a==INT_TO_SR(-1);
1131}
1132
1133BOOLEAN nrzEqual (number a,number b, const coeffs)
1134{
1135  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
1136    return a==b;
1137  else if (n_Z_IS_SMALL(a) || n_Z_IS_SMALL(b))
1138    return FALSE;
1139  else
1140    return 0 == mpz_cmp((int_number) a, (int_number) b);
1141}
1142
1143BOOLEAN nrzGreater (number a,number b, const coeffs)
1144{
1145  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
1146    return ((long)a)>((long)b);
1147  else if (n_Z_IS_SMALL(a))
1148    return 0 > mpz_cmp_si((int_number)b,SR_TO_INT(a));
1149  else if (n_Z_IS_SMALL(b))
1150    return 0 < mpz_cmp_si((int_number)a,SR_TO_INT(b));
1151  return 0 < mpz_cmp((int_number) a, (int_number) b);
1152}
1153
1154BOOLEAN nrzGreaterZero (number k, const coeffs C)
1155{
1156  return nrzGreater(k, INT_TO_SR(0), C);
1157}
1158
1159int nrzDivComp(number a, number b, const coeffs r)
1160{
1161  if (nrzDivBy(a, b, r))
1162  {
1163    if (nrzDivBy(b, a, r)) return 2;
1164    return -1;
1165  }
1166  if (nrzDivBy(b, a, r)) return 1;
1167  return 0;
1168}
1169
1170BOOLEAN nrzDivBy (number a,number b, const coeffs)
1171{
1172  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
1173  {
1174    return SR_TO_INT(a) %SR_TO_INT(b) ==0;
1175  }
1176  else if (n_Z_IS_SMALL(a))
1177  {
1178    return a==INT_TO_SR(0);
1179  }
1180  else if (n_Z_IS_SMALL(b))
1181  {
1182    return mpz_divisible_ui_p((int_number)a, (unsigned long)ABS(SR_TO_INT(b))) != 0;
1183  }
1184  else
1185    return mpz_divisible_p((int_number) a, (int_number) b) != 0;
1186}
1187
1188number nrzDiv (number a,number b, const coeffs)
1189{
1190  assume(SR_TO_INT(b));
1191  if (n_Z_IS_SMALL(a) && n_Z_IS_SMALL(b))
1192  {
1193    //if (SR_TO_INT(a) % SR_TO_INT(b))
1194    //{
1195    //  WerrorS("1:Division by non divisible element.");
1196    //  WerrorS("Result is without remainder.");
1197    //}
1198    return INT_TO_SR(SR_TO_INT(a)/SR_TO_INT(b));
1199  }
1200  else if (n_Z_IS_SMALL(a))
1201  {
1202    //if (SR_TO_INT(a))
1203    //{
1204    //  WerrorS("2:Division by non divisible element.");
1205    //  WerrorS("Result is without remainder.");
1206    //}
1207    return INT_TO_SR(0);
1208  }
1209  else if (n_Z_IS_SMALL(b))
1210  {
1211    int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1212    mpz_t r;
1213    mpz_init(r);
1214    mpz_init(erg);
1215    if (mpz_divmod_ui(erg, r, (int_number) a, (unsigned long)ABS(SR_TO_INT(b)))) {
1216    //  WerrorS("3:Division by non divisible element.");
1217    //  WerrorS("Result is without remainder.");
1218    }
1219    mpz_clear(r);
1220    if (SR_TO_INT(b)<0)
1221      mpz_neg(erg, erg);
1222    return nrz_short((number) erg);
1223  }
1224  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1225  mpz_init(erg);
1226  mpz_t r;
1227  mpz_init(r);
1228  mpz_tdiv_qr(erg, r, (int_number) a, (int_number) b);
1229#if CF_DEBUG
1230  StringSetS("division of");
1231  nrzWrite(a, R);
1232  StringAppendS(" by ");
1233  nrzWrite(b, R);
1234  StringAppendS(" is ");
1235  number du;
1236  nrzWrite(du = (number)erg, R);
1237  StringAppendS(" rest ");
1238  nrzWrite(du = (number)r, R);
1239  char * s = StringEndS();
1240  Print("%s\n", s);
1241  omFree(s);
1242#endif
1243
1244  if (mpz_cmp_si(r, 0)!=0)
1245  {
1246    //WerrorS("4:Division by non divisible element.");
1247    //WerrorS("Result is without remainder.");
1248  }
1249  mpz_clear(r);
1250  return nrz_short((number) erg);
1251}
1252
1253number nrzExactDiv (number a,number b, const coeffs)
1254{
1255  assume(SR_TO_INT(b));
1256  mpz_t aa, bb;
1257  if (n_Z_IS_SMALL(a))
1258    mpz_init_set_si(aa, SR_TO_INT(a));
1259  else
1260    mpz_init_set(aa, (int_number) a);
1261  if (n_Z_IS_SMALL(b))
1262    mpz_init_set_si(bb, SR_TO_INT(b));
1263  else
1264    mpz_init_set(bb, (int_number) b);
1265  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1266  mpz_init(erg);
1267  mpz_tdiv_q(erg, (int_number) aa, (int_number) bb);
1268  mpz_clear(aa);
1269  mpz_clear(bb);
1270  nrzTest((number)erg);
1271  return (number) erg;
1272}
1273
1274number nrzIntMod (number a,number b, const coeffs)
1275{
1276  mpz_t aa, bb;
1277  assume(SR_TO_INT(b));
1278  if (n_Z_IS_SMALL(a))
1279    mpz_init_set_si(aa, SR_TO_INT(a));
1280  else
1281    mpz_init_set(aa, (int_number) a);
1282  if (n_Z_IS_SMALL(b))
1283    mpz_init_set_si(bb, SR_TO_INT(b));
1284  else
1285    mpz_init_set(bb, (int_number) b);
1286
1287  mpz_t erg;
1288  mpz_init(erg);
1289  int_number r = (int_number) omAllocBin(gmp_nrz_bin);
1290  mpz_init(r);
1291  mpz_tdiv_qr(erg, r, (int_number) aa, (int_number) bb);
1292  mpz_clear(erg);
1293  mpz_clear(aa);
1294  mpz_clear(bb);
1295 
1296  return nrz_short((number) r);
1297}
1298
1299number  nrzInvers (number c, const coeffs r)
1300{
1301  if (!nrzIsUnit((number) c, r))
1302  {
1303    WerrorS("Non invertible element.");
1304    return (number)0; //TODO
1305  }
1306  return c; // has to be 1 or -1....
1307}
1308
1309number nrzNeg (number c, const coeffs)
1310{
1311// nNeg inplace !!!
1312  if (n_Z_IS_SMALL(c))
1313    return INT_TO_SR(-SR_TO_INT(c));
1314  mpz_mul_si((int_number) c, (int_number) c, -1);
1315  return c;
1316}
1317
1318static number nrzFarey(number r, number N, const coeffs R)
1319{
1320  number a0 = nrzCopy(N, R);
1321  number b0 = nrzInit(0, R);
1322  number a1 = nrzCopy(r, R);
1323  number b1 = nrzInit(1, R);
1324  number two = nrzInit(2, R);
1325#if 0
1326  Print("Farey start with ");
1327  n_Print(r, R);
1328  Print(" mod ");
1329  n_Print(N, R);
1330  Print("\n");
1331#endif
1332  while (1)
1333  {
1334    number as = nrzMult(a1, a1, R);
1335    n_InpMult(as, two, R);
1336    if (nrzGreater(N, as, R))
1337    {
1338      nrzDelete(&as, R);
1339      break;
1340    }
1341    nrzDelete(&as, R);
1342    number q = nrzDiv(a0, a1, R);
1343    number t = nrzMult(a1, q, R),
1344           s = nrzSub(a0, t, R);
1345    nrzDelete(&a0, R);
1346    a0 = a1;
1347    a1 = s;
1348    nrzDelete(&t, R);
1349
1350    t = nrzMult(b1, q, R);
1351    s = nrzSub(b0, t, R);
1352    nrzDelete(&b0, R);
1353    b0 = b1;
1354    b1 = s;
1355    nrzDelete(&t, R);
1356    nrzDelete(&q, R);
1357  }
1358  number as = nrzMult(b1, b1, R);
1359  n_InpMult(as, two, R);
1360  nrzDelete(&two, R);
1361  if (nrzGreater(as, N, R))
1362  {
1363    nrzDelete(&a0, R);
1364    nrzDelete(&a1, R);
1365    nrzDelete(&b0, R);
1366    nrzDelete(&b1, R);
1367    nrzDelete(&as, R);
1368    return NULL;
1369  }
1370  nrzDelete(&as, R);
1371  nrzDelete(&a0, R);
1372  nrzDelete(&b0, R);
1373
1374  number a, b, ab;
1375  coeffs Q = nInitChar(n_Q, 0);
1376  nMapFunc f = n_SetMap(R, Q);
1377  a = f(a1, R, Q);
1378  b = f(b1, R, Q);
1379  ab = n_Div(a, b, Q);
1380  n_Delete(&a, Q);
1381  n_Delete(&b, Q);
1382  nKillChar(Q);
1383
1384  nrzDelete(&a1, R);
1385  nrzDelete(&b1, R);
1386  return ab;
1387}
1388
1389number nrzMapMachineInt(number from, const coeffs /*src*/, const coeffs /*dst*/)
1390{
1391  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1392  mpz_init_set_ui(erg, (NATNUMBER) from);
1393  return nrz_short((number) erg);
1394}
1395
1396number nrzMapZp(number from, const coeffs /*src*/, const coeffs /*dst*/)
1397{
1398  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1399  mpz_init_set_si(erg, (long) from);
1400  return nrz_short((number) erg);
1401}
1402
1403number nrzModNMap(number from, const coeffs /* src */, const coeffs /*dst*/)
1404{
1405  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1406  mpz_init_set(erg, (int_number) from);
1407  return nrz_short((number) erg);
1408}
1409
1410number nrzMapQ(number from, const coeffs /* src */, const coeffs dst)
1411{
1412  if (SR_HDL(from) & SR_INT)
1413    return nrzInit(SR_TO_INT(from),dst);
1414  if (from->s!=3)
1415  {
1416    WerrorS("rational in map to integer");
1417    return NULL;
1418  }
1419  int_number erg = (int_number) omAllocBin(gmp_nrz_bin);
1420  mpz_init_set(erg, from->z);
1421  return nrz_short((number) erg);
1422}
1423
1424nMapFunc nrzSetMap(const coeffs src, const coeffs /*dst*/)
1425{
1426  /* dst = rintegers */
1427  if (src->rep==n_rep_gmp) //nCoeff_is_Ring_ModN(src) || nCoeff_is_Ring_PtoM(src))
1428    return nrzModNMap;
1429
1430  if ((src->rep==n_rep_gap_gmp) && nCoeff_is_Ring_Z(src))
1431  {
1432    return ndCopyMap; //nrzCopyMap;
1433  }
1434  if (src->rep==n_rep_gap_rat) /*&& nCoeff_is_Ring_Z(src)) Q, bigint*/
1435  {
1436    return nrzMapQ;
1437  }
1438  if ((src->rep=n_rep_int) && nCoeff_is_Ring_2toM(src))
1439  {
1440    return nrzMapMachineInt;
1441  }
1442  if ((src->rep=n_rep_int) && nCoeff_is_Zp(src))
1443  {
1444    return nrzMapZp;
1445  }
1446  return NULL;      // default
1447}
1448
1449
1450/*
1451 * set the exponent (allocate and init tables) (TODO)
1452 */
1453
1454void nrzSetExp(int, coeffs)
1455{
1456}
1457
1458void nrzInitExp(int, coeffs)
1459{
1460}
1461
1462#ifdef LDEBUG
1463BOOLEAN nrzDBTest (number x, const char *f, const int l, const coeffs)
1464{
1465  if (SR_HDL(x) & SR_INT) return TRUE;
1466  if (mpz_cmp_ui((int_number) x,(long)0)==0)
1467  {
1468    Print("gmp-0 %s:%d\n",f,l);
1469    return FALSE;
1470  }
1471  if (mpz_size1((int_number)x)<=MP_SMALL)
1472  {
1473    int ui=mpz_get_si((int_number)x);
1474    if ((((ui<<3)>>3)==ui)
1475    && (mpz_cmp_si((int_number)x,(long)ui)==0))
1476    {
1477      Print("gmp-small %s:%d\n",f,l);
1478      return FALSE;
1479    }
1480  }
1481  return TRUE;
1482}
1483#endif
1484
1485void nrzWrite (number &a, const coeffs)
1486{
1487  char *s,*z;
1488  if (a==NULL)
1489  {
1490    StringAppendS("o");
1491  }
1492  else
1493  {
1494    if (n_Z_IS_SMALL(a))
1495    {
1496      StringAppend("%d", SR_TO_INT(a));
1497    }
1498    else
1499    {
1500      int l=mpz_sizeinbase((int_number) a, 10) + 2;
1501      s=(char*)omAlloc(l);
1502      z=mpz_get_str(s,10,(int_number) a);
1503      StringAppendS(z);
1504      omFreeSize((ADDRESS)s,l);
1505    }
1506  }
1507}
1508
1509/*2
1510* extracts a long integer from s, returns the rest    (COPY FROM longrat0.cc)
1511*/
1512static const char * nlEatLongC(char *s, mpz_ptr i)
1513{
1514  const char * start=s;
1515
1516  if (*s<'0' || *s>'9')
1517  {
1518    mpz_set_si(i,1);
1519    return s;
1520  }
1521  while (*s >= '0' && *s <= '9') s++;
1522  if (*s=='\0')
1523  {
1524    mpz_set_str(i,start,10);
1525  }
1526  else
1527  {
1528    char c=*s;
1529    *s='\0';
1530    mpz_set_str(i,start,10);
1531    *s=c;
1532  }
1533  return s;
1534}
1535
1536const char * nrzRead (const char *s, number *a, const coeffs)
1537{
1538  int_number z = (int_number) omAllocBin(gmp_nrz_bin);
1539  {
1540    mpz_init(z);
1541    s = nlEatLongC((char *) s, z);
1542  }
1543  *a = nrz_short((number) z);
1544  return s;
1545}
1546
1547void    nrzCoeffWrite  (const coeffs, BOOLEAN /*details*/)
1548{
1549  //PrintS("// ZZ\n");
1550  PrintS("//   coeff. ring is : Integers\n");
1551}
1552
1553static char* nrzCoeffString(const coeffs)
1554{
1555  return omStrDup("integer");
1556}
1557
1558#include "factory/canonicalform.h"
1559#include <factory/cf_gmp.h>
1560#include "factory/singext.h"
1561
1562static CanonicalForm nrzConvSingNFactoryN( number n, BOOLEAN setChar, const coeffs /*r*/ )
1563{
1564  if (setChar) setCharacteristic( 0 );
1565
1566  CanonicalForm term;
1567  if ( n_Z_IS_SMALL(n))
1568  {
1569    term = SR_TO_INT(n);
1570  }
1571  else
1572  {
1573    mpz_t dummy;
1574    mpz_init_set( dummy,n->z );
1575    term = make_cf( dummy );
1576  }
1577  return term;
1578}
1579
1580static number nrzConvFactoryNSingN( const CanonicalForm n, const coeffs r)
1581{
1582  if (n.isImm())
1583  {
1584    return nrzInit(n.intval(),r);
1585  }
1586  else
1587  {
1588    if ( !n.den().isOne() )
1589    {
1590     WerrorS("rational in conversion to integer");
1591     return NULL;
1592    }
1593    int_number z = (int_number) omAlloc0Bin(gmp_nrz_bin);
1594    gmp_numerator( n,z);
1595    return nrz_short((number)z);
1596  }
1597}
1598
1599static void nrzMPZ(mpz_t res, number &a, const coeffs)
1600{
1601  if (n_Z_IS_SMALL(a))
1602    mpz_init_set_si(res, SR_TO_INT(a));
1603  else
1604    mpz_init_set(res, (int_number) a);
1605}
1606
1607coeffs nrzQuot1(number c, const coeffs r)
1608{
1609    int_number dummy;
1610    dummy = (int_number) omAlloc(sizeof(mpz_t));
1611    if(n_Z_IS_SMALL(c))
1612    {
1613      int ch = r->cfInt(c, r);
1614      mpz_init_set_ui(dummy, ch);
1615    }
1616    else
1617    {
1618      mpz_init_set(dummy, (int_number)c);
1619    }
1620    ZnmInfo info;
1621    info.base = dummy;
1622    info.exp = (unsigned long) 1;
1623    coeffs rr = nInitChar(n_Zn, (void*)&info);
1624    return(rr);
1625}
1626
1627BOOLEAN nrzInitChar(coeffs r,  void *)
1628{
1629  assume( getCoeffType(r) == ID );
1630
1631  r->is_field=FALSE;
1632  r->is_domain=TRUE;
1633  r->rep=n_rep_gap_gmp;
1634
1635  r->nCoeffIsEqual = ndCoeffIsEqual;
1636  r->cfCoeffString = nrzCoeffString;
1637  r->cfKillChar = ndKillChar;
1638  r->cfMult  = nrzMult;
1639  r->cfSub   = nrzSub;
1640  r->cfAdd   = nrzAdd;
1641  r->cfDiv   = nrzDiv;
1642  r->cfIntMod= nrzIntMod;
1643  r->cfExactDiv= nrzExactDiv;
1644  r->cfInit = nrzInit;
1645  r->cfSize  = nrzSize;
1646  r->cfInt  = nrzInt;
1647  //#ifdef HAVE_RINGS
1648  r->cfDivComp = nrzDivComp; // only for ring stuff
1649  r->cfIsUnit = nrzIsUnit; // only for ring stuff
1650  r->cfGetUnit = nrzGetUnit; // only for ring stuff
1651  r->cfAnn = nrzAnn;
1652  r->cfExtGcd = nrzExtGcd; // only for ring stuff
1653  r->cfXExtGcd = nrzXExtGcd; // only for ring stuff
1654  r->cfQuotRem = nrzQuotRem;
1655  r->cfDivBy = nrzDivBy; // only for ring stuff
1656  //#endif
1657  r->cfInpNeg   = nrzNeg;
1658  r->cfInvers= nrzInvers;
1659  r->cfCopy  = nrzCopy;
1660  r->cfWriteLong = nrzWrite;
1661  r->cfRead = nrzRead;
1662  r->cfGreater = nrzGreater;
1663  r->cfEqual = nrzEqual;
1664  r->cfIsZero = nrzIsZero;
1665  r->cfIsOne = nrzIsOne;
1666  r->cfIsMOne = nrzIsMOne;
1667  r->cfGreaterZero = nrzGreaterZero;
1668  r->cfPower = nrzPower;
1669  r->cfGcd  = nrzGcd;
1670  //r->cfLcm  = nrzLcm;
1671  r->cfDelete= nrzDelete;
1672  r->cfSetMap = nrzSetMap;
1673  r->cfCoeffWrite = nrzCoeffWrite;
1674  r->convSingNFactoryN = nrzConvSingNFactoryN;
1675  r->convFactoryNSingN = nrzConvFactoryNSingN;
1676  r->cfMPZ = nrzMPZ;
1677  r->cfFarey = nrzFarey;
1678
1679  r->cfQuot1 = nrzQuot1;
1680  // debug stuff
1681
1682#ifdef LDEBUG
1683  r->cfDBTest=nrzDBTest;
1684#endif
1685
1686  r->nNULL = 0;
1687  r->ch = 0;
1688  r->has_simple_Alloc=FALSE;
1689  r->has_simple_Inverse=FALSE;
1690  return FALSE;
1691}
1692
1693#else
1694#error set SI_INTEGER_VARIANT
1695#endif
1696#endif
Note: See TracBrowser for help on using the repository browser.