longrat.cc
Go to the documentation of this file.
1 /****************************************
2 * Computer Algebra System SINGULAR *
3 ****************************************/
4 /*
5 * ABSTRACT: computation with long rational numbers (Hubert Grassmann)
6 */
7 
8 #include <misc/auxiliary.h>
9 #include <omalloc/omalloc.h>
10 
11 #include <factory/factory.h>
12 
13 #include <misc/sirandom.h>
14 #include <misc/prime.h>
15 #include <reporter/reporter.h>
16 
17 #include "rmodulon.h" // ZnmInfo
18 #include "longrat.h"
19 #include "shortfl.h"
20 #include "modulop.h"
21 
22 // allow inlining only from p_Numbers.h and if ! LDEBUG
23 #if defined(DO_LINLINE) && defined(P_NUMBERS_H) && !defined(LDEBUG)
24 #define LINLINE static FORCE_INLINE
25 #else
26 #define LINLINE
27 #undef DO_LINLINE
28 #endif // DO_LINLINE
29 
30 LINLINE BOOLEAN nlEqual(number a, number b, const coeffs r);
31 LINLINE number nlInit(long i, const coeffs r);
32 LINLINE BOOLEAN nlIsOne(number a, const coeffs r);
33 LINLINE BOOLEAN nlIsZero(number za, const coeffs r);
34 LINLINE number nlCopy(number a, const coeffs r);
35 LINLINE number nl_Copy(number a, const coeffs r);
36 LINLINE void nlDelete(number *a, const coeffs r);
37 LINLINE number nlNeg(number za, const coeffs r);
38 LINLINE number nlAdd(number la, number li, const coeffs r);
39 LINLINE number nlSub(number la, number li, const coeffs r);
40 LINLINE number nlMult(number a, number b, const coeffs r);
41 LINLINE void nlInpAdd(number &a, number b, const coeffs r);
42 LINLINE void nlInpMult(number &a, number b, const coeffs r);
43 
44 number nlRInit (long i);
45 
46 
47 // number nlInitMPZ(mpz_t m, const coeffs r);
48 // void nlMPZ(mpz_t m, number &n, const coeffs r);
49 
50 void nlNormalize(number &x, const coeffs r);
51 
52 number nlGcd(number a, number b, const coeffs r);
53 number nlExtGcd(number a, number b, number *s, number *t, const coeffs);
54 number nlNormalizeHelper(number a, number b, const coeffs r); /*special routine !*/
55 BOOLEAN nlGreater(number a, number b, const coeffs r);
56 BOOLEAN nlIsMOne(number a, const coeffs r);
57 long nlInt(number &n, const coeffs r);
58 number nlBigInt(number &n);
59 
60 #ifdef HAVE_RINGS
61 number nlMapGMP(number from, const coeffs src, const coeffs dst);
62 #endif
63 
64 BOOLEAN nlGreaterZero(number za, const coeffs r);
65 number nlInvers(number a, const coeffs r);
66 number nlDiv(number a, number b, const coeffs r);
67 number nlExactDiv(number a, number b, const coeffs r);
68 number nlIntDiv(number a, number b, const coeffs r);
69 number nlIntMod(number a, number b, const coeffs r);
70 void nlPower(number x, int exp, number *lu, const coeffs r);
71 const char * nlRead (const char *s, number *a, const coeffs r);
72 void nlWrite(number a, const coeffs r);
73 
74 void nlCoeffWrite(const coeffs r, BOOLEAN details);
75 number nlChineseRemainder(number *x, number *q,int rl, const coeffs C);
76 number nlFarey(number nN, number nP, const coeffs CF);
77 
78 #ifdef LDEBUG
79 BOOLEAN nlDBTest(number a, const char *f, const int l);
80 #endif
81 
82 nMapFunc nlSetMap(const coeffs src, const coeffs dst);
83 
84 // in-place operations
85 void nlInpIntDiv(number &a, number b, const coeffs r);
86 
87 #ifdef LDEBUG
88 #define nlTest(a, r) nlDBTest(a,__FILE__,__LINE__, r)
89 BOOLEAN nlDBTest(number a, const char *f,int l, const coeffs r);
90 #else
91 #define nlTest(a, r) do {} while (0)
92 #endif
93 
94 
95 // 64 bit version:
96 //#if SIZEOF_LONG == 8
97 #if 0
98 #define MAX_NUM_SIZE 60
99 #define POW_2_28 (1L<<60)
100 #define POW_2_28_32 (1L<<28)
101 #define LONG long
102 #else
103 #define MAX_NUM_SIZE 28
104 #define POW_2_28 (1L<<28)
105 #define POW_2_28_32 (1L<<28)
106 #define LONG int
107 #endif
108 
109 
110 static inline number nlShort3(number x) // assume x->s==3
111 {
112  assume(x->s==3);
113  if (mpz_cmp_ui(x->z,0L)==0)
114  {
115  mpz_clear(x->z);
116  FREE_RNUMBER(x);
117  return INT_TO_SR(0);
118  }
119  if (mpz_size1(x->z)<=MP_SMALL)
120  {
121  LONG ui=mpz_get_si(x->z);
122  if ((((ui<<3)>>3)==ui)
123  && (mpz_cmp_si(x->z,(long)ui)==0))
124  {
125  mpz_clear(x->z);
126  FREE_RNUMBER(x);
127  return INT_TO_SR(ui);
128  }
129  }
130  return x;
131 }
132 
133 #ifndef LONGRAT_CC
134 #define LONGRAT_CC
135 
136 #include <string.h>
137 #include <float.h>
138 
139 #include <coeffs/coeffs.h>
140 #include <reporter/reporter.h>
141 #include <omalloc/omalloc.h>
142 
143 #include <coeffs/numbers.h>
144 #include <coeffs/mpr_complex.h>
145 
146 #ifndef BYTES_PER_MP_LIMB
147 #define BYTES_PER_MP_LIMB sizeof(mp_limb_t)
148 #endif
149 
150 //#define SR_HDL(A) ((long)(A))
151 /*#define SR_INT 1L*/
152 /*#define INT_TO_SR(INT) ((number) (((long)INT << 2) + SR_INT))*/
153 // #define SR_TO_INT(SR) (((long)SR) >> 2)
154 
155 #define MP_SMALL 1
156 //#define mpz_isNeg(A) (mpz_cmp_si(A,0L)<0)
157 #define mpz_isNeg(A) ((A)->_mp_size<0)
158 #define mpz_limb_size(A) ((A)->_mp_size)
159 #define mpz_limb_d(A) ((A)->_mp_d)
160 
161 void _nlDelete_NoImm(number *a);
162 
163 /***************************************************************
164  *
165  * Routines which are never inlined by p_Numbers.h
166  *
167  *******************************************************************/
168 #ifndef P_NUMBERS_H
169 
170 number nlShort3_noinline(number x) // assume x->s==3
171 {
172  return nlShort3(x);
173 }
174 
175 
176 #if (__GNU_MP_VERSION*10+__GNU_MP_VERSION_MINOR < 31)
177 void mpz_mul_si (mpz_ptr r, mpz_srcptr s, long int si)
178 {
179  if (si>=0)
180  mpz_mul_ui(r,s,si);
181  else
182  {
183  mpz_mul_ui(r,s,-si);
184  mpz_neg(r,r);
185  }
186 }
187 #endif
188 
189 static number nlMapP(number from, const coeffs src, const coeffs dst)
190 {
191  assume( getCoeffType(src) == n_Zp );
192 
193  number to = nlInit(npInt(from,src), dst); // FIXME? TODO? // extern long npInt (number &n, const coeffs r);
194 
195  return to;
196 }
197 
198 static number nlMapLongR(number from, const coeffs src, const coeffs dst);
199 static number nlMapR(number from, const coeffs src, const coeffs dst);
200 
201 
202 #ifdef HAVE_RINGS
203 /*2
204 * convert from a GMP integer
205 */
206 number nlMapGMP(number from, const coeffs /*src*/, const coeffs /*dst*/)
207 {
208  number z=ALLOC_RNUMBER();
209 #if defined(LDEBUG)
210  z->debug=123456;
211 #endif
212  mpz_init_set(z->z,(mpz_ptr) from);
213  //mpz_init_set_ui(&z->n,1);
214  z->s = 3;
215  z=nlShort3(z);
216  return z;
217 }
218 
219 number nlMapZ(number from, const coeffs src, const coeffs dst)
220 {
221  if (SR_HDL(from) & SR_INT)
222  {
223  return from;
224  }
225  return nlMapGMP(from,src,dst);
226 }
227 
228 /*2
229 * convert from an machine long
230 */
231 number nlMapMachineInt(number from, const coeffs /*src*/, const coeffs /*dst*/)
232 {
233  number z=ALLOC_RNUMBER();
234 #if defined(LDEBUG)
235  z->debug=123456;
236 #endif
237  mpz_init_set_ui(z->z,(unsigned long) from);
238  z->s = 3;
239  z=nlShort3(z);
240  return z;
241 }
242 #endif
243 
244 
245 #ifdef LDEBUG
246 BOOLEAN nlDBTest(number a, const char *f,const int l, const coeffs /*r*/)
247 {
248  if (a==NULL)
249  {
250  Print("!!longrat: NULL in %s:%d\n",f,l);
251  return FALSE;
252  }
253  //if ((int)a==1) Print("!! 0x1 as number ? %s %d\n",f,l);
254  if ((((long)a)&3L)==3L)
255  {
256  Print(" !!longrat:ptr(3) in %s:%d\n",f,l);
257  return FALSE;
258  }
259  if ((((long)a)&3L)==1L)
260  {
261  if (((((LONG)(long)a)<<1)>>1)!=((LONG)(long)a))
262  {
263  Print(" !!longrat:arith:%lx in %s:%d\n",(long)a, f,l);
264  return FALSE;
265  }
266  return TRUE;
267  }
268  /* TODO: If next line is active, then computations in algebraic field
269  extensions over Q will throw a lot of assume violations although
270  everything is computed correctly and no seg fault appears.
271  Maybe the test is not appropriate in this case. */
272  omCheckIf(omCheckAddrSize(a,sizeof(*a)), return FALSE);
273  if (a->debug!=123456)
274  {
275  Print("!!longrat:debug:%d in %s:%d\n",a->debug,f,l);
276  a->debug=123456;
277  return FALSE;
278  }
279  if ((a->s<0)||(a->s>4))
280  {
281  Print("!!longrat:s=%d in %s:%d\n",a->s,f,l);
282  return FALSE;
283  }
284  /* TODO: If next line is active, then computations in algebraic field
285  extensions over Q will throw a lot of assume violations although
286  everything is computed correctly and no seg fault appears.
287  Maybe the test is not appropriate in this case. */
288  //omCheckAddrSize(a->z[0]._mp_d,a->z[0]._mp_alloc*BYTES_PER_MP_LIMB);
289  if (a->z[0]._mp_alloc==0)
290  Print("!!longrat:z->alloc=0 in %s:%d\n",f,l);
291 
292  if (a->s<2)
293  {
294  if ((a->n[0]._mp_d[0]==0)&&(a->n[0]._mp_alloc<=1))
295  {
296  Print("!!longrat: n==0 in %s:%d\n",f,l);
297  return FALSE;
298  }
299  /* TODO: If next line is active, then computations in algebraic field
300  extensions over Q will throw a lot of assume violations although
301  everything is computed correctly and no seg fault appears.
302  Maybe the test is not appropriate in this case. */
303  //omCheckIf(omCheckAddrSize(a->n[0]._mp_d,a->n[0]._mp_alloc*BYTES_PER_MP_LIMB), return FALSE);
304  if (a->z[0]._mp_alloc==0)
305  Print("!!longrat:n->alloc=0 in %s:%d\n",f,l);
306  if ((mpz_size1(a->n) ==1) && (mpz_cmp_si(a->n,1L)==0))
307  {
308  Print("!!longrat:integer as rational in %s:%d\n",f,l);
309  mpz_clear(a->n); a->s=3;
310  return FALSE;
311  }
312  else if (mpz_isNeg(a->n))
313  {
314  Print("!!longrat:div. by negative in %s:%d\n",f,l);
315  mpz_neg(a->z,a->z);
316  mpz_neg(a->n,a->n);
317  return FALSE;
318  }
319  return TRUE;
320  }
321  //if (a->s==2)
322  //{
323  // Print("!!longrat:s=2 in %s:%d\n",f,l);
324  // return FALSE;
325  //}
326  if (mpz_size1(a->z)>MP_SMALL) return TRUE;
327  LONG ui=(LONG)mpz_get_si(a->z);
328  if ((((ui<<3)>>3)==ui)
329  && (mpz_cmp_si(a->z,(long)ui)==0))
330  {
331  Print("!!longrat:im int %d in %s:%d\n",ui,f,l);
332  return FALSE;
333  }
334  return TRUE;
335 }
336 #endif
337 
338 static CanonicalForm nlConvSingNFactoryN( number n, const BOOLEAN setChar, const coeffs /*r*/ )
339 {
340  if (setChar) setCharacteristic( 0 );
341 
343  if ( SR_HDL(n) & SR_INT )
344  {
345  long nn=SR_TO_INT(n);
346  term = nn;
347  }
348  else
349  {
350  if ( n->s == 3 )
351  {
352  mpz_t dummy;
353  long lz=mpz_get_si(n->z);
354  if (mpz_cmp_si(n->z,lz)==0) term=lz;
355  else
356  {
357  mpz_init_set( dummy,n->z );
358  term = make_cf( dummy );
359  }
360  }
361  else
362  {
363  // assume s==0 or s==1
364  mpz_t num, den;
365  On(SW_RATIONAL);
366  mpz_init_set( num, n->z );
367  mpz_init_set( den, n->n );
368  term = make_cf( num, den, ( n->s != 1 ));
369  }
370  }
371  return term;
372 }
373 
374 number nlRInit (long i);
375 
376 static number nlConvFactoryNSingN( const CanonicalForm f, const coeffs r)
377 {
378  if (f.isImm())
379  {
380  return nlInit(f.intval(),r);
381  }
382  else
383  {
384  number z = ALLOC_RNUMBER();
385 #if defined(LDEBUG)
386  z->debug=123456;
387 #endif
388  gmp_numerator( f, z->z );
389  if ( f.den().isOne() )
390  {
391  z->s = 3;
392  z=nlShort3(z);
393  }
394  else
395  {
396  gmp_denominator( f, z->n );
397  z->s = 1;
398  }
399  return z;
400  }
401 }
402 
403 static number nlMapR(number from, const coeffs src, const coeffs dst)
404 {
405  assume( getCoeffType(src) == n_R );
406 
407  double f=nrFloat(from); // FIXME? TODO? // extern float nrFloat(number n);
408  if (f==0.0) return INT_TO_SR(0);
409  int f_sign=1;
410  if (f<0.0)
411  {
412  f_sign=-1;
413  f=-f;
414  }
415  int i=0;
416  mpz_t h1;
417  mpz_init_set_ui(h1,1);
418  while((FLT_RADIX*f) < DBL_MAX && i<DBL_MANT_DIG)
419  {
420  f*=FLT_RADIX;
421  mpz_mul_ui(h1,h1,FLT_RADIX);
422  i++;
423  }
424  number re=nlRInit(1);
425  mpz_set_d(re->z,f);
426  memcpy(&(re->n),&h1,sizeof(h1));
427  re->s=0; /* not normalized */
428  if(f_sign==-1) re=nlNeg(re,dst);
429  nlNormalize(re,dst);
430  return re;
431 }
432 
433 static number nlMapLongR(number from, const coeffs src, const coeffs dst)
434 {
435  assume( getCoeffType(src) == n_long_R );
436 
437  gmp_float *ff=(gmp_float*)from;
438  mpf_t *f=ff->_mpfp();
439  number res;
440  mpz_ptr dest,ndest;
441  int size, i,negative;
442  int e,al,bl;
443  mp_ptr qp,dd,nn;
444 
445  size = (*f)[0]._mp_size;
446  if (size == 0)
447  return INT_TO_SR(0);
448  if(size<0)
449  {
450  negative = 1;
451  size = -size;
452  }
453  else
454  negative = 0;
455 
456  qp = (*f)[0]._mp_d;
457  while(qp[0]==0)
458  {
459  qp++;
460  size--;
461  }
462 
463  e=(*f)[0]._mp_exp-size;
464  res = ALLOC_RNUMBER();
465 #if defined(LDEBUG)
466  res->debug=123456;
467 #endif
468  dest = res->z;
469 
470  void* (*allocfunc) (size_t);
471  mp_get_memory_functions (&allocfunc,NULL, NULL);
472  if (e<0)
473  {
474  al = dest->_mp_size = size;
475  if (al<2) al = 2;
476  dd = (mp_ptr)allocfunc(sizeof(mp_limb_t)*al);
477  for (i=0;i<size;i++) dd[i] = qp[i];
478  bl = 1-e;
479  nn = (mp_ptr)allocfunc(sizeof(mp_limb_t)*bl);
480  nn[bl-1] = 1;
481  ndest = res->n;
482  ndest->_mp_d = nn;
483  ndest->_mp_alloc = ndest->_mp_size = bl;
484  res->s = 0;
485  }
486  else
487  {
488  al = dest->_mp_size = size+e;
489  if (al<2) al = 2;
490  dd = (mp_ptr)allocfunc(sizeof(mp_limb_t)*al);
491  for (i=0;i<size;i++) dd[i+e] = qp[i];
492  for (i=0;i<e;i++) dd[i] = 0;
493  res->s = 3;
494  }
495 
496  dest->_mp_d = dd;
497  dest->_mp_alloc = al;
498  if (negative) mpz_neg(dest,dest);
499 
500  if (res->s==0)
501  nlNormalize(res,dst);
502  else if (mpz_size1(res->z)<=MP_SMALL)
503  {
504  // res is new, res->ref is 1
505  res=nlShort3(res);
506  }
507  nlTest(res, dst);
508  return res;
509 }
510 
511 //static number nlMapLongR(number from)
512 //{
513 // gmp_float *ff=(gmp_float*)from;
514 // const mpf_t *f=ff->mpfp();
515 // int f_size=ABS((*f)[0]._mp_size);
516 // if (f_size==0)
517 // return nlInit(0);
518 // int f_sign=1;
519 // number work=ngcCopy(from);
520 // if (!ngcGreaterZero(work))
521 // {
522 // f_sign=-1;
523 // work=ngcNeg(work);
524 // }
525 // int i=0;
526 // mpz_t h1;
527 // mpz_init_set_ui(h1,1);
528 // while((FLT_RADIX*f) < DBL_MAX && i<DBL_MANT_DIG)
529 // {
530 // f*=FLT_RADIX;
531 // mpz_mul_ui(h1,h1,FLT_RADIX);
532 // i++;
533 // }
534 // number r=nlRInit(1);
535 // mpz_set_d(&(r->z),f);
536 // memcpy(&(r->n),&h1,sizeof(h1));
537 // r->s=0; /* not normalized */
538 // nlNormalize(r);
539 // return r;
540 //
541 //
542 // number r=nlRInit(1);
543 // int f_shift=f_size+(*f)[0]._mp_exp;
544 // if ( f_shift > 0)
545 // {
546 // r->s=0;
547 // mpz_init(&r->n);
548 // mpz_setbit(&r->n,f_shift*BYTES_PER_MP_LIMB*8);
549 // mpz_setbit(&r->z,f_size*BYTES_PER_MP_LIMB*8-1);
550 // // now r->z has enough space
551 // memcpy(mpz_limb_d(&r->z),((*f)[0]._mp_d),f_size*BYTES_PER_MP_LIMB);
552 // nlNormalize(r);
553 // }
554 // else
555 // {
556 // r->s=3;
557 // if (f_shift==0)
558 // {
559 // mpz_setbit(&r->z,f_size*BYTES_PER_MP_LIMB*8-1);
560 // // now r->z has enough space
561 // memcpy(mpz_limb_d(&r->z),((*f)[0]._mp_d),f_size*BYTES_PER_MP_LIMB);
562 // }
563 // else /* f_shift < 0 */
564 // {
565 // mpz_setbit(&r->z,(f_size-f_shift)*BYTES_PER_MP_LIMB*8-1);
566 // // now r->z has enough space
567 // memcpy(mpz_limb_d(&r->z)-f_shift,((*f)[0]._mp_d),
568 // f_size*BYTES_PER_MP_LIMB);
569 // }
570 // }
571 // if ((*f)[0]._mp_size<0);
572 // r=nlNeg(r);
573 // return r;
574 //}
575 
576 int nlSize(number a, const coeffs)
577 {
578  if (a==INT_TO_SR(0))
579  return 0; /* rational 0*/
580  if (SR_HDL(a) & SR_INT)
581  return 1; /* immediate int */
582  int s=a->z[0]._mp_alloc;
583 // while ((s>0) &&(a->z._mp_d[s]==0L)) s--;
584 //#if SIZEOF_LONG == 8
585 // if (a->z._mp_d[s] < (unsigned long)0x100000000L) s=s*2-1;
586 // else s *=2;
587 //#endif
588 // s++;
589  if (a->s<2)
590  {
591  int d=a->n[0]._mp_alloc;
592 // while ((d>0) && (a->n._mp_d[d]==0L)) d--;
593 //#if SIZEOF_LONG == 8
594 // if (a->n._mp_d[d] < (unsigned long)0x100000000L) d=d*2-1;
595 // else d *=2;
596 //#endif
597  s+=d;
598  }
599  return s;
600 }
601 
602 /*2
603 * convert number to int
604 */
605 long nlInt(number &i, const coeffs r)
606 {
607  nlTest(i, r);
608  nlNormalize(i,r);
609  if (SR_HDL(i) & SR_INT)
610  {
611  return SR_TO_INT(i);
612  }
613  if (i->s==3)
614  {
615  if(mpz_size1(i->z)>MP_SMALL) return 0;
616  long ul=mpz_get_si(i->z);
617  if (mpz_cmp_si(i->z,ul)!=0) return 0;
618  return ul;
619  }
620  mpz_t tmp;
621  long ul;
622  mpz_init(tmp);
623  mpz_tdiv_q(tmp,i->z,i->n);
624  if(mpz_size1(tmp)>MP_SMALL) ul=0;
625  else
626  {
627  ul=mpz_get_si(tmp);
628  if (mpz_cmp_si(tmp,ul)!=0) ul=0;
629  }
630  mpz_clear(tmp);
631  return ul;
632 }
633 
634 /*2
635 * convert number to bigint
636 */
637 number nlBigInt(number &i, const coeffs r)
638 {
639  nlTest(i, r);
640  nlNormalize(i,r);
641  if (SR_HDL(i) & SR_INT) return (i);
642  if (i->s==3)
643  {
644  return nlCopy(i,r);
645  }
646  number tmp=nlRInit(1);
647  mpz_tdiv_q(tmp->z,i->z,i->n);
648  tmp=nlShort3(tmp);
649  return tmp;
650 }
651 
652 /*
653 * 1/a
654 */
655 number nlInvers(number a, const coeffs r)
656 {
657  nlTest(a, r);
658  number n;
659  if (SR_HDL(a) & SR_INT)
660  {
661  if ((a==INT_TO_SR(1L)) || (a==INT_TO_SR(-1L)))
662  {
663  return a;
664  }
665  if (nlIsZero(a,r))
666  {
667  WerrorS(nDivBy0);
668  return INT_TO_SR(0);
669  }
670  n=ALLOC_RNUMBER();
671 #if defined(LDEBUG)
672  n->debug=123456;
673 #endif
674  n->s=1;
675  if (((long)a)>0L)
676  {
677  mpz_init_set_si(n->z,1L);
678  mpz_init_set_si(n->n,(long)SR_TO_INT(a));
679  }
680  else
681  {
682  mpz_init_set_si(n->z,-1L);
683  mpz_init_set_si(n->n,(long)-SR_TO_INT(a));
684  }
685  nlTest(n, r);
686  return n;
687  }
688  n=ALLOC_RNUMBER();
689 #if defined(LDEBUG)
690  n->debug=123456;
691 #endif
692  {
693  mpz_init_set(n->n,a->z);
694  switch (a->s)
695  {
696  case 0:
697  case 1:
698  n->s=a->s;
699  mpz_init_set(n->z,a->n);
700  if (mpz_isNeg(n->n)) /* && n->s<2*/
701  {
702  mpz_neg(n->z,n->z);
703  mpz_neg(n->n,n->n);
704  }
705  if (mpz_cmp_ui(n->n,1L)==0)
706  {
707  mpz_clear(n->n);
708  n->s=3;
709  n=nlShort3(n);
710  }
711  break;
712  case 3:
713  // i.e. |a| > 2^...
714  n->s=1;
715  if (mpz_isNeg(n->n)) /* && n->s<2*/
716  {
717  mpz_neg(n->n,n->n);
718  mpz_init_set_si(n->z,-1L);
719  }
720  else
721  {
722  mpz_init_set_si(n->z,1L);
723  }
724  break;
725  }
726  }
727  nlTest(n, r);
728  return n;
729 }
730 
731 
732 /*2
733 * u := a / b in Z, if b | a (else undefined)
734 */
735 number nlExactDiv(number a, number b, const coeffs r)
736 {
737  if (b==INT_TO_SR(0))
738  {
739  WerrorS(nDivBy0);
740  return INT_TO_SR(0);
741  }
742  if (a==INT_TO_SR(0))
743  return INT_TO_SR(0);
744  number u;
745  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
746  {
747  /* the small int -(1<<28) divided by -1 is the large int (1<<28) */
748  if ((a==INT_TO_SR(-(POW_2_28)))&&(b==INT_TO_SR(-1L)))
749  {
750  return nlRInit(POW_2_28);
751  }
752  long aa=SR_TO_INT(a);
753  long bb=SR_TO_INT(b);
754  return INT_TO_SR(aa/bb);
755  }
756  number bb=NULL;
757  if (SR_HDL(b) & SR_INT)
758  {
759  bb=nlRInit(SR_TO_INT(b));
760  b=bb;
761  }
762  u=ALLOC_RNUMBER();
763 #if defined(LDEBUG)
764  u->debug=123456;
765 #endif
766  mpz_init(u->z);
767  /* u=a/b */
768  u->s = 3;
769  mpz_divexact(u->z,a->z,b->z);
770  if (bb!=NULL)
771  {
772  mpz_clear(bb->z);
773 #if defined(LDEBUG)
774  bb->debug=654324;
775 #endif
776  FREE_RNUMBER(bb); // omFreeBin((void *)bb, rnumber_bin);
777  }
778  u=nlShort3(u);
779  nlTest(u, r);
780  return u;
781 }
782 
783 /*2
784 * u := a / b in Z
785 */
786 number nlIntDiv (number a, number b, const coeffs r)
787 {
788  if (b==INT_TO_SR(0))
789  {
790  WerrorS(nDivBy0);
791  return INT_TO_SR(0);
792  }
793  if (a==INT_TO_SR(0))
794  return INT_TO_SR(0);
795  number u;
796  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
797  {
798  /* the small int -(1<<28) divided by -1 is the large int (1<<28) */
799  if ((a==INT_TO_SR(-(POW_2_28)))&&(b==INT_TO_SR(-1L)))
800  {
801  return nlRInit(POW_2_28);
802  }
803  LONG aa=SR_TO_INT(a);
804  LONG bb=SR_TO_INT(b);
805  LONG rr=aa%bb;
806  if (rr<0) rr+=ABS(bb);
807  LONG cc=(aa-rr)/bb;
808  return INT_TO_SR(cc);
809  }
810  number aa=NULL;
811  if (SR_HDL(a) & SR_INT)
812  {
813  /* the small int -(1<<28) divided by 2^28 is 1 */
814  if (a==INT_TO_SR(-(POW_2_28)))
815  {
816  if(mpz_cmp_si(b->z,(POW_2_28))==0)
817  {
818  return INT_TO_SR(-1);
819  }
820  }
821  aa=nlRInit(SR_TO_INT(a));
822  a=aa;
823  }
824  number bb=NULL;
825  if (SR_HDL(b) & SR_INT)
826  {
827  bb=nlRInit(SR_TO_INT(b));
828  b=bb;
829  }
830  u=ALLOC_RNUMBER();
831 #if defined(LDEBUG)
832  u->debug=123456;
833 #endif
834  assume(a->s==3);
835  assume(b->s==3);
836  mpz_init_set(u->z,a->z);
837  /* u=u/b */
838  u->s = 3;
839  number rr=nlIntMod(a,b,r);
840  if (SR_HDL(rr) & SR_INT) mpz_sub_ui(u->z,u->z,SR_TO_INT(rr));
841  else mpz_sub(u->z,u->z,rr->z);
842  mpz_divexact(u->z,u->z,b->z);
843  if (aa!=NULL)
844  {
845  mpz_clear(aa->z);
846 #if defined(LDEBUG)
847  aa->debug=654324;
848 #endif
849  FREE_RNUMBER(aa);
850  }
851  if (bb!=NULL)
852  {
853  mpz_clear(bb->z);
854 #if defined(LDEBUG)
855  bb->debug=654324;
856 #endif
857  FREE_RNUMBER(bb);
858  }
859  u=nlShort3(u);
860  nlTest(u,r);
861  return u;
862 }
863 
864 /*2
865 * u := a mod b in Z, u>=0
866 */
867 number nlIntMod (number a, number b, const coeffs r)
868 {
869  if (b==INT_TO_SR(0))
870  {
871  WerrorS(nDivBy0);
872  return INT_TO_SR(0);
873  }
874  if (a==INT_TO_SR(0))
875  return INT_TO_SR(0);
876  number u;
877  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
878  {
879  LONG aa=SR_TO_INT(a);
880  LONG bb=SR_TO_INT(b);
881  LONG c=aa % bb;
882  if (c<0) c+=ABS(bb);
883  return INT_TO_SR(c);
884  }
885  if (SR_HDL(a) & SR_INT)
886  {
887  LONG ai=SR_TO_INT(a);
888  mpz_t aa;
889  mpz_init_set_si(aa, ai);
890  u=ALLOC_RNUMBER();
891 #if defined(LDEBUG)
892  u->debug=123456;
893 #endif
894  u->s = 3;
895  mpz_init(u->z);
896  mpz_mod(u->z, aa, b->z);
897  mpz_clear(aa);
898  u=nlShort3(u);
899  nlTest(u,r);
900  return u;
901  }
902  number bb=NULL;
903  if (SR_HDL(b) & SR_INT)
904  {
905  bb=nlRInit(SR_TO_INT(b));
906  b=bb;
907  }
908  u=ALLOC_RNUMBER();
909 #if defined(LDEBUG)
910  u->debug=123456;
911 #endif
912  mpz_init(u->z);
913  u->s = 3;
914  mpz_mod(u->z, a->z, b->z);
915  if (bb!=NULL)
916  {
917  mpz_clear(bb->z);
918 #if defined(LDEBUG)
919  bb->debug=654324;
920 #endif
921  FREE_RNUMBER(bb);
922  }
923  u=nlShort3(u);
924  nlTest(u,r);
925  return u;
926 }
927 
928 BOOLEAN nlDivBy (number a,number b, const coeffs)
929 {
930  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
931  {
932  return ((SR_TO_INT(a) % SR_TO_INT(b))==0);
933  }
934  if (SR_HDL(b) & SR_INT)
935  {
936  return (mpz_divisible_ui_p(a->z,SR_TO_INT(b))!=0);
937  }
938  if (SR_HDL(a) & SR_INT) return FALSE;
939  return mpz_divisible_p(a->z, b->z) != 0;
940 }
941 
942 int nlDivComp(number a, number b, const coeffs r)
943 {
944  if (nlDivBy(a, b, r))
945  {
946  if (nlDivBy(b, a, r)) return 2;
947  return -1;
948  }
949  if (nlDivBy(b, a, r)) return 1;
950  return 0;
951 }
952 
953 number nlGetUnit (number, const coeffs)
954 {
955  return INT_TO_SR(1);
956 }
957 
958 coeffs nlQuot1(number c, const coeffs r)
959 {
960  long ch = r->cfInt(c, r);
961  int p=IsPrime(ch);
962  coeffs rr=NULL;
963  if (((long)p)==ch)
964  {
965  rr = nInitChar(n_Zp,(void*)ch);
966  }
967  #ifdef HAVE_RINGS
968  else
969  {
970  mpz_ptr dummy;
971  dummy = (mpz_ptr) omAlloc(sizeof(mpz_t));
972  mpz_init_set_ui(dummy, ch);
973  ZnmInfo info;
974  info.base = dummy;
975  info.exp = (unsigned long) 1;
976  rr = nInitChar(n_Zn, (void*)&info);
977  }
978  #endif
979  return(rr);
980 }
981 
982 
983 BOOLEAN nlIsUnit (number a, const coeffs)
984 {
985  return ((SR_HDL(a) & SR_INT) && (ABS(SR_TO_INT(a))==1));
986 }
987 
988 
989 /*2
990 * u := a / b
991 */
992 number nlDiv (number a, number b, const coeffs r)
993 {
994  if (nlIsZero(b,r))
995  {
996  WerrorS(nDivBy0);
997  return INT_TO_SR(0);
998  }
999  number u;
1000 // ---------- short / short ------------------------------------
1001  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
1002  {
1003  LONG i=SR_TO_INT(a);
1004  LONG j=SR_TO_INT(b);
1005  if (j==1L) return a;
1006  if ((i==-POW_2_28) && (j== -1L))
1007  {
1008  return nlRInit(POW_2_28);
1009  }
1010  LONG r=i%j;
1011  if (r==0)
1012  {
1013  return INT_TO_SR(i/j);
1014  }
1015  u=ALLOC_RNUMBER();
1016  u->s=0;
1017  #if defined(LDEBUG)
1018  u->debug=123456;
1019  #endif
1020  mpz_init_set_si(u->z,(long)i);
1021  mpz_init_set_si(u->n,(long)j);
1022  }
1023  else
1024  {
1025  u=ALLOC_RNUMBER();
1026  u->s=0;
1027  #if defined(LDEBUG)
1028  u->debug=123456;
1029  #endif
1030  mpz_init(u->z);
1031 // ---------- short / long ------------------------------------
1032  if (SR_HDL(a) & SR_INT)
1033  {
1034  // short a / (z/n) -> (a*n)/z
1035  if (b->s<2)
1036  {
1037  mpz_mul_si(u->z,b->n,SR_TO_INT(a));
1038  }
1039  else
1040  // short a / long z -> a/z
1041  {
1042  mpz_set_si(u->z,SR_TO_INT(a));
1043  }
1044  if (mpz_cmp(u->z,b->z)==0)
1045  {
1046  mpz_clear(u->z);
1047  FREE_RNUMBER(u);
1048  return INT_TO_SR(1);
1049  }
1050  mpz_init_set(u->n,b->z);
1051  }
1052 // ---------- long / short ------------------------------------
1053  else if (SR_HDL(b) & SR_INT)
1054  {
1055  mpz_set(u->z,a->z);
1056  // (z/n) / b -> z/(n*b)
1057  if (a->s<2)
1058  {
1059  mpz_init_set(u->n,a->n);
1060  if (((long)b)>0L)
1061  mpz_mul_ui(u->n,u->n,SR_TO_INT(b));
1062  else
1063  {
1064  mpz_mul_ui(u->n,u->n,-SR_TO_INT(b));
1065  mpz_neg(u->z,u->z);
1066  }
1067  }
1068  else
1069  // long z / short b -> z/b
1070  {
1071  //mpz_set(u->z,a->z);
1072  mpz_init_set_si(u->n,SR_TO_INT(b));
1073  }
1074  }
1075 // ---------- long / long ------------------------------------
1076  else
1077  {
1078  mpz_set(u->z,a->z);
1079  mpz_init_set(u->n,b->z);
1080  if (a->s<2) mpz_mul(u->n,u->n,a->n);
1081  if (b->s<2) mpz_mul(u->z,u->z,b->n);
1082  }
1083  }
1084  if (mpz_isNeg(u->n))
1085  {
1086  mpz_neg(u->z,u->z);
1087  mpz_neg(u->n,u->n);
1088  }
1089  if (mpz_cmp_si(u->n,1L)==0)
1090  {
1091  mpz_clear(u->n);
1092  u->s=3;
1093  u=nlShort3(u);
1094  }
1095  nlTest(u, r);
1096  return u;
1097 }
1098 
1099 /*2
1100 * u:= x ^ exp
1101 */
1102 void nlPower (number x,int exp,number * u, const coeffs r)
1103 {
1104  *u = INT_TO_SR(0); // 0^e, e!=0
1105  if (exp==0)
1106  *u= INT_TO_SR(1);
1107  else if (!nlIsZero(x,r))
1108  {
1109  nlTest(x, r);
1110  number aa=NULL;
1111  if (SR_HDL(x) & SR_INT)
1112  {
1113  aa=nlRInit(SR_TO_INT(x));
1114  x=aa;
1115  }
1116  else if (x->s==0)
1117  nlNormalize(x,r);
1118  *u=ALLOC_RNUMBER();
1119 #if defined(LDEBUG)
1120  (*u)->debug=123456;
1121 #endif
1122  mpz_init((*u)->z);
1123  mpz_pow_ui((*u)->z,x->z,(unsigned long)exp);
1124  if (x->s<2)
1125  {
1126  if (mpz_cmp_si(x->n,1L)==0)
1127  {
1128  x->s=3;
1129  mpz_clear(x->n);
1130  }
1131  else
1132  {
1133  mpz_init((*u)->n);
1134  mpz_pow_ui((*u)->n,x->n,(unsigned long)exp);
1135  }
1136  }
1137  (*u)->s = x->s;
1138  if ((*u)->s==3) *u=nlShort3(*u);
1139  if (aa!=NULL)
1140  {
1141  mpz_clear(aa->z);
1142  FREE_RNUMBER(aa);
1143  }
1144  }
1145 #ifdef LDEBUG
1146  if (exp<0) Print("nlPower: neg. exp. %d\n",exp);
1147  nlTest(*u, r);
1148 #endif
1149 }
1150 
1151 
1152 /*2
1153 * za >= 0 ?
1154 */
1155 BOOLEAN nlGreaterZero (number a, const coeffs r)
1156 {
1157  nlTest(a, r);
1158  if (SR_HDL(a) & SR_INT) return SR_HDL(a)>1L /* represents number(0) */;
1159  return (!mpz_isNeg(a->z));
1160 }
1161 
1162 /*2
1163 * a > b ?
1164 */
1165 BOOLEAN nlGreater (number a, number b, const coeffs r)
1166 {
1167  nlTest(a, r);
1168  nlTest(b, r);
1169  number re;
1170  BOOLEAN rr;
1171  re=nlSub(a,b,r);
1172  rr=(!nlIsZero(re,r)) && (nlGreaterZero(re,r));
1173  nlDelete(&re,r);
1174  return rr;
1175 }
1176 
1177 /*2
1178 * a == -1 ?
1179 */
1180 BOOLEAN nlIsMOne (number a, const coeffs r)
1181 {
1182 #ifdef LDEBUG
1183  if (a==NULL) return FALSE;
1184  nlTest(a, r);
1185 #endif
1186  return (a==INT_TO_SR(-1L));
1187 }
1188 
1189 /*2
1190 * result =gcd(a,b)
1191 */
1192 number nlGcd(number a, number b, const coeffs r)
1193 {
1194  number result;
1195  nlTest(a, r);
1196  nlTest(b, r);
1197  //nlNormalize(a);
1198  //nlNormalize(b);
1199  if ((a==INT_TO_SR(1L))||(a==INT_TO_SR(-1L))
1200  || (b==INT_TO_SR(1L))||(b==INT_TO_SR(-1L)))
1201  return INT_TO_SR(1L);
1202  if (a==INT_TO_SR(0)) /* gcd(0,b) ->b */
1203  return nlCopy(b,r);
1204  if (b==INT_TO_SR(0)) /* gcd(a,0) -> a */
1205  return nlCopy(a,r);
1206  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
1207  {
1208  long i=SR_TO_INT(a);
1209  long j=SR_TO_INT(b);
1210  if((i==0L)||(j==0L))
1211  return INT_TO_SR(1);
1212  long l;
1213  i=ABS(i);
1214  j=ABS(j);
1215  do
1216  {
1217  l=i%j;
1218  i=j;
1219  j=l;
1220  } while (l!=0L);
1221  if (i==POW_2_28)
1222  result=nlRInit(POW_2_28);
1223  else
1224  result=INT_TO_SR(i);
1225  nlTest(result,r);
1226  return result;
1227  }
1228  if (((!(SR_HDL(a) & SR_INT))&&(a->s<2))
1229  || ((!(SR_HDL(b) & SR_INT))&&(b->s<2))) return INT_TO_SR(1);
1230  if (SR_HDL(a) & SR_INT)
1231  {
1232  LONG aa=ABS(SR_TO_INT(a));
1233  unsigned long t=mpz_gcd_ui(NULL,b->z,(long)aa);
1234  if (t==POW_2_28)
1235  result=nlRInit(POW_2_28);
1236  else
1237  result=INT_TO_SR(t);
1238  }
1239  else
1240  if (SR_HDL(b) & SR_INT)
1241  {
1242  LONG bb=ABS(SR_TO_INT(b));
1243  unsigned long t=mpz_gcd_ui(NULL,a->z,(long)bb);
1244  if (t==POW_2_28)
1245  result=nlRInit(POW_2_28);
1246  else
1247  result=INT_TO_SR(t);
1248  }
1249  else
1250  {
1251  result=ALLOC0_RNUMBER();
1252  result->s = 3;
1253  #ifdef LDEBUG
1254  result->debug=123456;
1255  #endif
1256  mpz_init(result->z);
1257  mpz_gcd(result->z,a->z,b->z);
1258  result=nlShort3(result);
1259  }
1260  nlTest(result, r);
1261  return result;
1262 }
1263 static int int_extgcd(int a, int b, int * u, int* x, int * v, int* y)
1264 {
1265  int q, r;
1266  if (a==0)
1267  {
1268  *u = 0;
1269  *v = 1;
1270  *x = -1;
1271  *y = 0;
1272  return b;
1273  }
1274  if (b==0)
1275  {
1276  *u = 1;
1277  *v = 0;
1278  *x = 0;
1279  *y = 1;
1280  return a;
1281  }
1282  *u=1;
1283  *v=0;
1284  *x=0;
1285  *y=1;
1286  do
1287  {
1288  q = a/b;
1289  r = a%b;
1290  assume (q*b+r == a);
1291  a = b;
1292  b = r;
1293 
1294  r = -(*v)*q+(*u);
1295  (*u) =(*v);
1296  (*v) = r;
1297 
1298  r = -(*y)*q+(*x);
1299  (*x) = (*y);
1300  (*y) = r;
1301  } while (b);
1302 
1303  return a;
1304 }
1305 
1306 //number nlGcd_dummy(number a, number b, const coeffs r)
1307 //{
1308 // extern char my_yylinebuf[80];
1309 // Print("nlGcd in >>%s<<\n",my_yylinebuf);
1310 // return nlGcd(a,b,r);;
1311 //}
1312 
1313 number nlShort1(number x) // assume x->s==0/1
1314 {
1315  assume(x->s<2);
1316  if (mpz_cmp_ui(x->z,0L)==0)
1317  {
1318  _nlDelete_NoImm(&x);
1319  return INT_TO_SR(0);
1320  }
1321  if (x->s<2)
1322  {
1323  if (mpz_cmp(x->z,x->n)==0)
1324  {
1325  _nlDelete_NoImm(&x);
1326  return INT_TO_SR(1);
1327  }
1328  }
1329  return x;
1330 }
1331 /*2
1332 * simplify x
1333 */
1334 void nlNormalize (number &x, const coeffs r)
1335 {
1336  if ((SR_HDL(x) & SR_INT) ||(x==NULL))
1337  return;
1338  if (x->s==3)
1339  {
1340  x=nlShort3_noinline(x);
1341  nlTest(x,r);
1342  return;
1343  }
1344  else if (x->s==0)
1345  {
1346  if (mpz_cmp_si(x->n,1L)==0)
1347  {
1348  mpz_clear(x->n);
1349  x->s=3;
1350  x=nlShort3(x);
1351  }
1352  else
1353  {
1354  mpz_t gcd;
1355  mpz_init(gcd);
1356  mpz_gcd(gcd,x->z,x->n);
1357  x->s=1;
1358  if (mpz_cmp_si(gcd,1L)!=0)
1359  {
1360  mpz_divexact(x->z,x->z,gcd);
1361  mpz_divexact(x->n,x->n,gcd);
1362  if (mpz_cmp_si(x->n,1L)==0)
1363  {
1364  mpz_clear(x->n);
1365  x->s=3;
1366  x=nlShort3_noinline(x);
1367  }
1368  }
1369  mpz_clear(gcd);
1370  }
1371  }
1372  nlTest(x, r);
1373 }
1374 
1375 /*2
1376 * returns in result->z the lcm(a->z,b->n)
1377 */
1378 number nlNormalizeHelper(number a, number b, const coeffs r)
1379 {
1380  number result;
1381  nlTest(a, r);
1382  nlTest(b, r);
1383  if ((SR_HDL(b) & SR_INT)
1384  || (b->s==3))
1385  {
1386  // b is 1/(b->n) => b->n is 1 => result is a
1387  return nlCopy(a,r);
1388  }
1389  result=ALLOC_RNUMBER();
1390 #if defined(LDEBUG)
1391  result->debug=123456;
1392 #endif
1393  result->s=3;
1394  mpz_t gcd;
1395  mpz_init(gcd);
1396  mpz_init(result->z);
1397  if (SR_HDL(a) & SR_INT)
1398  mpz_gcd_ui(gcd,b->n,ABS(SR_TO_INT(a)));
1399  else
1400  mpz_gcd(gcd,a->z,b->n);
1401  if (mpz_cmp_si(gcd,1L)!=0)
1402  {
1403  mpz_t bt;
1404  mpz_init_set(bt,b->n);
1405  mpz_divexact(bt,bt,gcd);
1406  if (SR_HDL(a) & SR_INT)
1407  mpz_mul_si(result->z,bt,SR_TO_INT(a));
1408  else
1409  mpz_mul(result->z,bt,a->z);
1410  mpz_clear(bt);
1411  }
1412  else
1413  if (SR_HDL(a) & SR_INT)
1414  mpz_mul_si(result->z,b->n,SR_TO_INT(a));
1415  else
1416  mpz_mul(result->z,b->n,a->z);
1417  mpz_clear(gcd);
1418  result=nlShort3(result);
1419  nlTest(result, r);
1420  return result;
1421 }
1422 
1423 // Map q \in QQ or ZZ \to Zp or an extension of it
1424 // src = Q or Z, dst = Zp (or an extension of Zp)
1425 number nlModP(number q, const coeffs /*Q*/, const coeffs Zp)
1426 {
1427  const int p = n_GetChar(Zp);
1428  assume( p > 0 );
1429 
1430  const long P = p;
1431  assume( P > 0 );
1432 
1433  // embedded long within q => only long numerator has to be converted
1434  // to int (modulo char.)
1435  if (SR_HDL(q) & SR_INT)
1436  {
1437  long i = SR_TO_INT(q);
1438  return n_Init( i, Zp );
1439  }
1440 
1441  const unsigned long PP = p;
1442 
1443  // numerator modulo char. should fit into int
1444  number z = n_Init( static_cast<long>(mpz_fdiv_ui(q->z, PP)), Zp );
1445 
1446  // denominator != 1?
1447  if (q->s!=3)
1448  {
1449  // denominator modulo char. should fit into int
1450  number n = n_Init( static_cast<long>(mpz_fdiv_ui(q->n, PP)), Zp );
1451 
1452  number res = n_Div( z, n, Zp );
1453 
1454  n_Delete(&z, Zp);
1455  n_Delete(&n, Zp);
1456 
1457  return res;
1458  }
1459 
1460  return z;
1461 }
1462 
1463 #ifdef HAVE_RINGS
1464 /*2
1465 * convert number i (from Q) to GMP and warn if denom != 1
1466 */
1467 void nlGMP(number &i, number n, const coeffs r)
1468 {
1469  // Hier brauche ich einfach die GMP Zahl
1470  nlTest(i, r);
1471  nlNormalize(i, r);
1472  if (SR_HDL(i) & SR_INT)
1473  {
1474  mpz_set_si((mpz_ptr) n, SR_TO_INT(i));
1475  return;
1476  }
1477  if (i->s!=3)
1478  {
1479  WarnS("Omitted denominator during coefficient mapping !");
1480  }
1481  mpz_set((mpz_ptr) n, i->z);
1482 }
1483 #endif
1484 
1485 /*2
1486 * acces to denominator, other 1 for integers
1487 */
1488 number nlGetDenom(number &n, const coeffs r)
1489 {
1490  if (!(SR_HDL(n) & SR_INT))
1491  {
1492  if (n->s==0)
1493  {
1494  nlNormalize(n,r);
1495  }
1496  if (!(SR_HDL(n) & SR_INT))
1497  {
1498  if (n->s!=3)
1499  {
1500  number u=ALLOC_RNUMBER();
1501  u->s=3;
1502 #if defined(LDEBUG)
1503  u->debug=123456;
1504 #endif
1505  mpz_init_set(u->z,n->n);
1506  u=nlShort3_noinline(u);
1507  return u;
1508  }
1509  }
1510  }
1511  return INT_TO_SR(1);
1512 }
1513 
1514 /*2
1515 * acces to Nominator, nlCopy(n) for integers
1516 */
1517 number nlGetNumerator(number &n, const coeffs r)
1518 {
1519  if (!(SR_HDL(n) & SR_INT))
1520  {
1521  if (n->s==0)
1522  {
1523  nlNormalize(n,r);
1524  }
1525  if (!(SR_HDL(n) & SR_INT))
1526  {
1527  number u=ALLOC_RNUMBER();
1528 #if defined(LDEBUG)
1529  u->debug=123456;
1530 #endif
1531  u->s=3;
1532  mpz_init_set(u->z,n->z);
1533  if (n->s!=3)
1534  {
1535  u=nlShort3_noinline(u);
1536  }
1537  return u;
1538  }
1539  }
1540  return n; // imm. int
1541 }
1542 
1543 /***************************************************************
1544  *
1545  * routines which are needed by Inline(d) routines
1546  *
1547  *******************************************************************/
1549 {
1550  assume(! (SR_HDL(a) & SR_HDL(b) & SR_INT));
1551 // long - short
1552  BOOLEAN bo;
1553  if (SR_HDL(b) & SR_INT)
1554  {
1555  if (a->s!=0) return FALSE;
1556  number n=b; b=a; a=n;
1557  }
1558 // short - long
1559  if (SR_HDL(a) & SR_INT)
1560  {
1561  if (b->s!=0)
1562  return FALSE;
1563  if ((((long)a) > 0L) && (mpz_isNeg(b->z)))
1564  return FALSE;
1565  if ((((long)a) < 0L) && (!mpz_isNeg(b->z)))
1566  return FALSE;
1567  mpz_t bb;
1568  mpz_init_set(bb,b->n);
1569  mpz_mul_si(bb,bb,(long)SR_TO_INT(a));
1570  bo=(mpz_cmp(bb,b->z)==0);
1571  mpz_clear(bb);
1572  return bo;
1573  }
1574 // long - long
1575  if (((a->s==1) && (b->s==3))
1576  || ((b->s==1) && (a->s==3)))
1577  return FALSE;
1578  if (mpz_isNeg(a->z)&&(!mpz_isNeg(b->z)))
1579  return FALSE;
1580  if (mpz_isNeg(b->z)&&(!mpz_isNeg(a->z)))
1581  return FALSE;
1582  mpz_t aa;
1583  mpz_t bb;
1584  mpz_init_set(aa,a->z);
1585  mpz_init_set(bb,b->z);
1586  if (a->s<2) mpz_mul(bb,bb,a->n);
1587  if (b->s<2) mpz_mul(aa,aa,b->n);
1588  bo=(mpz_cmp(aa,bb)==0);
1589  mpz_clear(aa);
1590  mpz_clear(bb);
1591  return bo;
1592 }
1593 
1594 // copy not immediate number a
1595 number _nlCopy_NoImm(number a)
1596 {
1597  assume(!((SR_HDL(a) & SR_INT)||(a==NULL)));
1598  //nlTest(a, r);
1599  number b=ALLOC_RNUMBER();
1600 #if defined(LDEBUG)
1601  b->debug=123456;
1602 #endif
1603  switch (a->s)
1604  {
1605  case 0:
1606  case 1:
1607  mpz_init_set(b->n,a->n);
1608  case 3:
1609  mpz_init_set(b->z,a->z);
1610  break;
1611  }
1612  b->s = a->s;
1613  return b;
1614 }
1615 
1616 void _nlDelete_NoImm(number *a)
1617 {
1618  {
1619  switch ((*a)->s)
1620  {
1621  case 0:
1622  case 1:
1623  mpz_clear((*a)->n);
1624  case 3:
1625  mpz_clear((*a)->z);
1626 #ifdef LDEBUG
1627  (*a)->s=2;
1628 #endif
1629  }
1630  FREE_RNUMBER(*a); // omFreeBin((void *) *a, rnumber_bin);
1631  }
1632 }
1633 
1634 number _nlNeg_NoImm(number a)
1635 {
1636  {
1637  mpz_neg(a->z,a->z);
1638  if (a->s==3)
1639  {
1640  a=nlShort3(a);
1641  }
1642  }
1643  return a;
1644 }
1645 
1646 // conditio to use nlNormalize_Gcd in intermediate computations:
1647 #define GCD_NORM_COND(OLD,NEW) (mpz_size1(NEW->z)>mpz_size1(OLD->z))
1648 
1649 static void nlNormalize_Gcd(number &x)
1650 {
1651  mpz_t gcd;
1652  mpz_init(gcd);
1653  mpz_gcd(gcd,x->z,x->n);
1654  x->s=1;
1655  if (mpz_cmp_si(gcd,1L)!=0)
1656  {
1657  mpz_divexact(x->z,x->z,gcd);
1658  mpz_divexact(x->n,x->n,gcd);
1659  if (mpz_cmp_si(x->n,1L)==0)
1660  {
1661  mpz_clear(x->n);
1662  x->s=3;
1663  x=nlShort3_noinline(x);
1664  }
1665  }
1666  mpz_clear(gcd);
1667 }
1668 
1669 number _nlAdd_aNoImm_OR_bNoImm(number a, number b)
1670 {
1671  number u=ALLOC_RNUMBER();
1672 #if defined(LDEBUG)
1673  u->debug=123456;
1674 #endif
1675  mpz_init(u->z);
1676  if (SR_HDL(b) & SR_INT)
1677  {
1678  number x=a;
1679  a=b;
1680  b=x;
1681  }
1682  if (SR_HDL(a) & SR_INT)
1683  {
1684  switch (b->s)
1685  {
1686  case 0:
1687  case 1:/* a:short, b:1 */
1688  {
1689  mpz_t x;
1690  mpz_init(x);
1691  mpz_mul_si(x,b->n,SR_TO_INT(a));
1692  mpz_add(u->z,b->z,x);
1693  mpz_clear(x);
1694  if (mpz_cmp_ui(u->z,0L)==0)
1695  {
1696  mpz_clear(u->z);
1697  FREE_RNUMBER(u);
1698  return INT_TO_SR(0);
1699  }
1700  if (mpz_cmp(u->z,b->n)==0)
1701  {
1702  mpz_clear(u->z);
1703  FREE_RNUMBER(u);
1704  return INT_TO_SR(1);
1705  }
1706  mpz_init_set(u->n,b->n);
1707  u->s = 0;
1708  if (GCD_NORM_COND(b,u)) { nlNormalize_Gcd(u); }
1709  break;
1710  }
1711  case 3:
1712  {
1713  if (((long)a)>0L)
1714  mpz_add_ui(u->z,b->z,SR_TO_INT(a));
1715  else
1716  mpz_sub_ui(u->z,b->z,-SR_TO_INT(a));
1717  u->s = 3;
1718  u=nlShort3(u);
1719  break;
1720  }
1721  }
1722  }
1723  else
1724  {
1725  switch (a->s)
1726  {
1727  case 0:
1728  case 1:
1729  {
1730  switch(b->s)
1731  {
1732  case 0:
1733  case 1:
1734  {
1735  mpz_t x;
1736  mpz_init(x);
1737 
1738  mpz_mul(x,b->z,a->n);
1739  mpz_mul(u->z,a->z,b->n);
1740  mpz_add(u->z,u->z,x);
1741  mpz_clear(x);
1742 
1743  if (mpz_cmp_ui(u->z,0L)==0)
1744  {
1745  mpz_clear(u->z);
1746  FREE_RNUMBER(u);
1747  return INT_TO_SR(0);
1748  }
1749  mpz_init(u->n);
1750  mpz_mul(u->n,a->n,b->n);
1751  if (mpz_cmp(u->z,u->n)==0)
1752  {
1753  mpz_clear(u->z);
1754  mpz_clear(u->n);
1755  FREE_RNUMBER(u);
1756  return INT_TO_SR(1);
1757  }
1758  u->s = 0;
1759  if (GCD_NORM_COND(b,u)) { nlNormalize_Gcd(u); }
1760  break;
1761  }
1762  case 3: /* a:1 b:3 */
1763  {
1764  mpz_mul(u->z,b->z,a->n);
1765  mpz_add(u->z,u->z,a->z);
1766  if (mpz_cmp_ui(u->z,0L)==0)
1767  {
1768  mpz_clear(u->z);
1769  FREE_RNUMBER(u);
1770  return INT_TO_SR(0);
1771  }
1772  if (mpz_cmp(u->z,a->n)==0)
1773  {
1774  mpz_clear(u->z);
1775  FREE_RNUMBER(u);
1776  return INT_TO_SR(1);
1777  }
1778  mpz_init_set(u->n,a->n);
1779  u->s = 0;
1780  if (GCD_NORM_COND(a,u)) { nlNormalize_Gcd(u); }
1781  break;
1782  }
1783  } /*switch (b->s) */
1784  break;
1785  }
1786  case 3:
1787  {
1788  switch(b->s)
1789  {
1790  case 0:
1791  case 1:/* a:3, b:1 */
1792  {
1793  mpz_mul(u->z,a->z,b->n);
1794  mpz_add(u->z,u->z,b->z);
1795  if (mpz_cmp_ui(u->z,0L)==0)
1796  {
1797  mpz_clear(u->z);
1798  FREE_RNUMBER(u);
1799  return INT_TO_SR(0);
1800  }
1801  if (mpz_cmp(u->z,b->n)==0)
1802  {
1803  mpz_clear(u->z);
1804  FREE_RNUMBER(u);
1805  return INT_TO_SR(1);
1806  }
1807  mpz_init_set(u->n,b->n);
1808  u->s = 0;
1809  if (GCD_NORM_COND(b,u)) { nlNormalize_Gcd(u); }
1810  break;
1811  }
1812  case 3:
1813  {
1814  mpz_add(u->z,a->z,b->z);
1815  u->s = 3;
1816  u=nlShort3(u);
1817  break;
1818  }
1819  }
1820  break;
1821  }
1822  }
1823  }
1824  return u;
1825 }
1826 
1827 void _nlInpAdd_aNoImm_OR_bNoImm(number &a, number b)
1828 {
1829  if (SR_HDL(b) & SR_INT)
1830  {
1831  switch (a->s)
1832  {
1833  case 0:
1834  case 1:/* b:short, a:1 */
1835  {
1836  mpz_t x;
1837  mpz_init(x);
1838  mpz_mul_si(x,a->n,SR_TO_INT(b));
1839  mpz_add(a->z,a->z,x);
1840  mpz_clear(x);
1841  nlNormalize_Gcd(a);
1842  break;
1843  }
1844  case 3:
1845  {
1846  if (((long)b)>0L)
1847  mpz_add_ui(a->z,a->z,SR_TO_INT(b));
1848  else
1849  mpz_sub_ui(a->z,a->z,-SR_TO_INT(b));
1850  a->s = 3;
1851  a=nlShort3_noinline(a);
1852  break;
1853  }
1854  }
1855  return;
1856  }
1857  else if (SR_HDL(a) & SR_INT)
1858  {
1859  number u=ALLOC_RNUMBER();
1860  #if defined(LDEBUG)
1861  u->debug=123456;
1862  #endif
1863  mpz_init(u->z);
1864  switch (b->s)
1865  {
1866  case 0:
1867  case 1:/* a:short, b:1 */
1868  {
1869  mpz_t x;
1870  mpz_init(x);
1871 
1872  mpz_mul_si(x,b->n,SR_TO_INT(a));
1873  mpz_add(u->z,b->z,x);
1874  mpz_clear(x);
1875  // result cannot be 0, if coeffs are normalized
1876  mpz_init_set(u->n,b->n);
1877  u->s=0;
1878  if (GCD_NORM_COND(b,u)) { nlNormalize_Gcd(u); }
1879  else { u=nlShort1(u); }
1880  break;
1881  }
1882  case 3:
1883  {
1884  if (((long)a)>0L)
1885  mpz_add_ui(u->z,b->z,SR_TO_INT(a));
1886  else
1887  mpz_sub_ui(u->z,b->z,-SR_TO_INT(a));
1888  // result cannot be 0, if coeffs are normalized
1889  u->s = 3;
1890  u=nlShort3_noinline(u);
1891  break;
1892  }
1893  }
1894  a=u;
1895  }
1896  else
1897  {
1898  switch (a->s)
1899  {
1900  case 0:
1901  case 1:
1902  {
1903  switch(b->s)
1904  {
1905  case 0:
1906  case 1: /* a:1 b:1 */
1907  {
1908  mpz_t x;
1909  mpz_t y;
1910  mpz_init(x);
1911  mpz_init(y);
1912  mpz_mul(x,b->z,a->n);
1913  mpz_mul(y,a->z,b->n);
1914  mpz_add(a->z,x,y);
1915  mpz_clear(x);
1916  mpz_clear(y);
1917  mpz_mul(a->n,a->n,b->n);
1918  a->s=0;
1919  if (GCD_NORM_COND(b,a)) { nlNormalize_Gcd(a); }
1920  else { a=nlShort1(a);}
1921  break;
1922  }
1923  case 3: /* a:1 b:3 */
1924  {
1925  mpz_t x;
1926  mpz_init(x);
1927  mpz_mul(x,b->z,a->n);
1928  mpz_add(a->z,a->z,x);
1929  mpz_clear(x);
1930  a->s=0;
1931  if (GCD_NORM_COND(b,a)) { nlNormalize_Gcd(a); }
1932  else { a=nlShort1(a);}
1933  break;
1934  }
1935  } /*switch (b->s) */
1936  break;
1937  }
1938  case 3:
1939  {
1940  switch(b->s)
1941  {
1942  case 0:
1943  case 1:/* a:3, b:1 */
1944  {
1945  mpz_t x;
1946  mpz_init(x);
1947  mpz_mul(x,a->z,b->n);
1948  mpz_add(a->z,b->z,x);
1949  mpz_clear(x);
1950  mpz_init_set(a->n,b->n);
1951  a->s=0;
1952  if (GCD_NORM_COND(b,a)) { nlNormalize_Gcd(a); }
1953  else { a=nlShort1(a);}
1954  break;
1955  }
1956  case 3:
1957  {
1958  mpz_add(a->z,a->z,b->z);
1959  a->s = 3;
1960  a=nlShort3_noinline(a);
1961  break;
1962  }
1963  }
1964  break;
1965  }
1966  }
1967  }
1968 }
1969 
1970 number _nlSub_aNoImm_OR_bNoImm(number a, number b)
1971 {
1972  number u=ALLOC_RNUMBER();
1973 #if defined(LDEBUG)
1974  u->debug=123456;
1975 #endif
1976  mpz_init(u->z);
1977  if (SR_HDL(a) & SR_INT)
1978  {
1979  switch (b->s)
1980  {
1981  case 0:
1982  case 1:/* a:short, b:1 */
1983  {
1984  mpz_t x;
1985  mpz_init(x);
1986  mpz_mul_si(x,b->n,SR_TO_INT(a));
1987  mpz_sub(u->z,x,b->z);
1988  mpz_clear(x);
1989  if (mpz_cmp_ui(u->z,0L)==0)
1990  {
1991  mpz_clear(u->z);
1992  FREE_RNUMBER(u);
1993  return INT_TO_SR(0);
1994  }
1995  if (mpz_cmp(u->z,b->n)==0)
1996  {
1997  mpz_clear(u->z);
1998  FREE_RNUMBER(u);
1999  return INT_TO_SR(1);
2000  }
2001  mpz_init_set(u->n,b->n);
2002  u->s=0;
2003  if (GCD_NORM_COND(b,u)) { nlNormalize_Gcd(u); }
2004  break;
2005  }
2006  case 3:
2007  {
2008  if (((long)a)>0L)
2009  {
2010  mpz_sub_ui(u->z,b->z,SR_TO_INT(a));
2011  mpz_neg(u->z,u->z);
2012  }
2013  else
2014  {
2015  mpz_add_ui(u->z,b->z,-SR_TO_INT(a));
2016  mpz_neg(u->z,u->z);
2017  }
2018  u->s = 3;
2019  u=nlShort3(u);
2020  break;
2021  }
2022  }
2023  }
2024  else if (SR_HDL(b) & SR_INT)
2025  {
2026  switch (a->s)
2027  {
2028  case 0:
2029  case 1:/* b:short, a:1 */
2030  {
2031  mpz_t x;
2032  mpz_init(x);
2033  mpz_mul_si(x,a->n,SR_TO_INT(b));
2034  mpz_sub(u->z,a->z,x);
2035  mpz_clear(x);
2036  if (mpz_cmp_ui(u->z,0L)==0)
2037  {
2038  mpz_clear(u->z);
2039  FREE_RNUMBER(u);
2040  return INT_TO_SR(0);
2041  }
2042  if (mpz_cmp(u->z,a->n)==0)
2043  {
2044  mpz_clear(u->z);
2045  FREE_RNUMBER(u);
2046  return INT_TO_SR(1);
2047  }
2048  mpz_init_set(u->n,a->n);
2049  u->s=0;
2050  if (GCD_NORM_COND(a,u)) { nlNormalize_Gcd(u); }
2051  break;
2052  }
2053  case 3:
2054  {
2055  if (((long)b)>0L)
2056  {
2057  mpz_sub_ui(u->z,a->z,SR_TO_INT(b));
2058  }
2059  else
2060  {
2061  mpz_add_ui(u->z,a->z,-SR_TO_INT(b));
2062  }
2063  u->s = 3;
2064  u=nlShort3(u);
2065  break;
2066  }
2067  }
2068  }
2069  else
2070  {
2071  switch (a->s)
2072  {
2073  case 0:
2074  case 1:
2075  {
2076  switch(b->s)
2077  {
2078  case 0:
2079  case 1:
2080  {
2081  mpz_t x;
2082  mpz_t y;
2083  mpz_init(x);
2084  mpz_init(y);
2085  mpz_mul(x,b->z,a->n);
2086  mpz_mul(y,a->z,b->n);
2087  mpz_sub(u->z,y,x);
2088  mpz_clear(x);
2089  mpz_clear(y);
2090  if (mpz_cmp_ui(u->z,0L)==0)
2091  {
2092  mpz_clear(u->z);
2093  FREE_RNUMBER(u);
2094  return INT_TO_SR(0);
2095  }
2096  mpz_init(u->n);
2097  mpz_mul(u->n,a->n,b->n);
2098  if (mpz_cmp(u->z,u->n)==0)
2099  {
2100  mpz_clear(u->z);
2101  mpz_clear(u->n);
2102  FREE_RNUMBER(u);
2103  return INT_TO_SR(1);
2104  }
2105  u->s=0;
2106  if (GCD_NORM_COND(a,u)) { nlNormalize_Gcd(u); }
2107  break;
2108  }
2109  case 3: /* a:1, b:3 */
2110  {
2111  mpz_t x;
2112  mpz_init(x);
2113  mpz_mul(x,b->z,a->n);
2114  mpz_sub(u->z,a->z,x);
2115  mpz_clear(x);
2116  if (mpz_cmp_ui(u->z,0L)==0)
2117  {
2118  mpz_clear(u->z);
2119  FREE_RNUMBER(u);
2120  return INT_TO_SR(0);
2121  }
2122  if (mpz_cmp(u->z,a->n)==0)
2123  {
2124  mpz_clear(u->z);
2125  FREE_RNUMBER(u);
2126  return INT_TO_SR(1);
2127  }
2128  mpz_init_set(u->n,a->n);
2129  u->s=0;
2130  if (GCD_NORM_COND(a,u)) { nlNormalize_Gcd(u); }
2131  break;
2132  }
2133  }
2134  break;
2135  }
2136  case 3:
2137  {
2138  switch(b->s)
2139  {
2140  case 0:
2141  case 1: /* a:3, b:1 */
2142  {
2143  mpz_t x;
2144  mpz_init(x);
2145  mpz_mul(x,a->z,b->n);
2146  mpz_sub(u->z,x,b->z);
2147  mpz_clear(x);
2148  if (mpz_cmp_ui(u->z,0L)==0)
2149  {
2150  mpz_clear(u->z);
2151  FREE_RNUMBER(u);
2152  return INT_TO_SR(0);
2153  }
2154  if (mpz_cmp(u->z,b->n)==0)
2155  {
2156  mpz_clear(u->z);
2157  FREE_RNUMBER(u);
2158  return INT_TO_SR(1);
2159  }
2160  mpz_init_set(u->n,b->n);
2161  u->s=0;
2162  if (GCD_NORM_COND(b,u)) { nlNormalize_Gcd(u); }
2163  break;
2164  }
2165  case 3: /* a:3 , b:3 */
2166  {
2167  mpz_sub(u->z,a->z,b->z);
2168  u->s = 3;
2169  u=nlShort3(u);
2170  break;
2171  }
2172  }
2173  break;
2174  }
2175  }
2176  }
2177  return u;
2178 }
2179 
2180 // a and b are intermediate, but a*b not
2181 number _nlMult_aImm_bImm_rNoImm(number a, number b)
2182 {
2183  number u=ALLOC_RNUMBER();
2184 #if defined(LDEBUG)
2185  u->debug=123456;
2186 #endif
2187  u->s=3;
2188  mpz_init_set_si(u->z,SR_TO_INT(a));
2189  mpz_mul_si(u->z,u->z,SR_TO_INT(b));
2190  return u;
2191 }
2192 
2193 // a or b are not immediate
2194 number _nlMult_aNoImm_OR_bNoImm(number a, number b)
2195 {
2196  assume(! (SR_HDL(a) & SR_HDL(b) & SR_INT));
2197  number u=ALLOC_RNUMBER();
2198 #if defined(LDEBUG)
2199  u->debug=123456;
2200 #endif
2201  mpz_init(u->z);
2202  if (SR_HDL(b) & SR_INT)
2203  {
2204  number x=a;
2205  a=b;
2206  b=x;
2207  }
2208  if (SR_HDL(a) & SR_INT)
2209  {
2210  u->s=b->s;
2211  if (u->s==1) u->s=0;
2212  if (((long)a)>0L)
2213  {
2214  mpz_mul_ui(u->z,b->z,(unsigned long)SR_TO_INT(a));
2215  }
2216  else
2217  {
2218  if (a==INT_TO_SR(-1))
2219  {
2220  mpz_set(u->z,b->z);
2221  mpz_neg(u->z,u->z);
2222  u->s=b->s;
2223  }
2224  else
2225  {
2226  mpz_mul_ui(u->z,b->z,(unsigned long)-SR_TO_INT(a));
2227  mpz_neg(u->z,u->z);
2228  }
2229  }
2230  if (u->s<2)
2231  {
2232  if (mpz_cmp(u->z,b->n)==0)
2233  {
2234  mpz_clear(u->z);
2235  FREE_RNUMBER(u);
2236  return INT_TO_SR(1);
2237  }
2238  mpz_init_set(u->n,b->n);
2239  if (GCD_NORM_COND(b,u)) { nlNormalize_Gcd(u); }
2240  }
2241  else //u->s==3
2242  {
2243  u=nlShort3(u);
2244  }
2245  }
2246  else
2247  {
2248  mpz_mul(u->z,a->z,b->z);
2249  u->s = 0;
2250  if(a->s==3)
2251  {
2252  if(b->s==3)
2253  {
2254  u->s = 3;
2255  }
2256  else
2257  {
2258  if (mpz_cmp(u->z,b->n)==0)
2259  {
2260  mpz_clear(u->z);
2261  FREE_RNUMBER(u);
2262  return INT_TO_SR(1);
2263  }
2264  mpz_init_set(u->n,b->n);
2265  if (GCD_NORM_COND(b,u)) { nlNormalize_Gcd(u); }
2266  }
2267  }
2268  else
2269  {
2270  if(b->s==3)
2271  {
2272  if (mpz_cmp(u->z,a->n)==0)
2273  {
2274  mpz_clear(u->z);
2275  FREE_RNUMBER(u);
2276  return INT_TO_SR(1);
2277  }
2278  mpz_init_set(u->n,a->n);
2279  if (GCD_NORM_COND(a,u)) { nlNormalize_Gcd(u); }
2280  }
2281  else
2282  {
2283  mpz_init(u->n);
2284  mpz_mul(u->n,a->n,b->n);
2285  if (mpz_cmp(u->z,u->n)==0)
2286  {
2287  mpz_clear(u->z);
2288  mpz_clear(u->n);
2289  FREE_RNUMBER(u);
2290  return INT_TO_SR(1);
2291  }
2292  if (GCD_NORM_COND(a,u)) { nlNormalize_Gcd(u); }
2293  }
2294  }
2295  }
2296  return u;
2297 }
2298 
2299 /*2
2300 * copy a to b for mapping
2301 */
2302 number nlCopyMap(number a, const coeffs /*src*/, const coeffs /*dst*/)
2303 {
2304  if ((SR_HDL(a) & SR_INT)||(a==NULL))
2305  {
2306  return a;
2307  }
2308  return _nlCopy_NoImm(a);
2309 }
2310 
2311 nMapFunc nlSetMap(const coeffs src, const coeffs /*dst*/)
2312 {
2313  if (src->rep==n_rep_gap_rat) /*Q, coeffs_BIGINT */
2314  {
2315  return ndCopyMap;
2316  }
2317  if ((src->rep==n_rep_int) && nCoeff_is_Zp(src))
2318  {
2319  return nlMapP;
2320  }
2321  if ((src->rep==n_rep_float) && nCoeff_is_R(src))
2322  {
2323  return nlMapR;
2324  }
2325  if ((src->rep==n_rep_gmp_float) && nCoeff_is_long_R(src))
2326  {
2327  return nlMapLongR; /* long R -> Q */
2328  }
2329 #ifdef HAVE_RINGS
2330  if (src->rep==n_rep_gmp) // nCoeff_is_Ring_Z(src) || nCoeff_is_Ring_PtoM(src) || nCoeff_is_Ring_ModN(src))
2331  {
2332  return nlMapGMP;
2333  }
2334  if (src->rep==n_rep_gap_gmp)
2335  {
2336  return nlMapZ;
2337  }
2338  if ((src->rep==n_rep_int) && nCoeff_is_Ring_2toM(src))
2339  {
2340  return nlMapMachineInt;
2341  }
2342 #endif
2343  return NULL;
2344 }
2345 /*2
2346 * z := i
2347 */
2348 number nlRInit (long i)
2349 {
2350  number z=ALLOC_RNUMBER();
2351 #if defined(LDEBUG)
2352  z->debug=123456;
2353 #endif
2354  mpz_init_set_si(z->z,i);
2355  z->s = 3;
2356  return z;
2357 }
2358 
2359 /*2
2360 * z := i/j
2361 */
2362 number nlInit2 (int i, int j, const coeffs r)
2363 {
2364  number z=ALLOC_RNUMBER();
2365 #if defined(LDEBUG)
2366  z->debug=123456;
2367 #endif
2368  mpz_init_set_si(z->z,(long)i);
2369  mpz_init_set_si(z->n,(long)j);
2370  z->s = 0;
2371  nlNormalize(z,r);
2372  return z;
2373 }
2374 
2375 number nlInit2gmp (mpz_t i, mpz_t j, const coeffs r)
2376 {
2377  number z=ALLOC_RNUMBER();
2378 #if defined(LDEBUG)
2379  z->debug=123456;
2380 #endif
2381  mpz_init_set(z->z,i);
2382  mpz_init_set(z->n,j);
2383  z->s = 0;
2384  nlNormalize(z,r);
2385  return z;
2386 }
2387 
2388 #else // DO_LINLINE
2389 
2390 // declare immedate routines
2391 number nlRInit (long i);
2392 BOOLEAN _nlEqual_aNoImm_OR_bNoImm(number a, number b);
2393 number _nlCopy_NoImm(number a);
2394 number _nlNeg_NoImm(number a);
2395 number _nlAdd_aNoImm_OR_bNoImm(number a, number b);
2396 void _nlInpAdd_aNoImm_OR_bNoImm(number &a, number b);
2397 number _nlSub_aNoImm_OR_bNoImm(number a, number b);
2398 number _nlMult_aNoImm_OR_bNoImm(number a, number b);
2399 number _nlMult_aImm_bImm_rNoImm(number a, number b);
2400 
2401 #endif
2402 
2403 /***************************************************************
2404  *
2405  * Routines which might be inlined by p_Numbers.h
2406  *
2407  *******************************************************************/
2408 #if defined(DO_LINLINE) || !defined(P_NUMBERS_H)
2409 
2410 // routines which are always inlined/static
2411 
2412 /*2
2413 * a = b ?
2414 */
2415 LINLINE BOOLEAN nlEqual (number a, number b, const coeffs r)
2416 {
2417  nlTest(a, r);
2418  nlTest(b, r);
2419 // short - short
2420  if (SR_HDL(a) & SR_HDL(b) & SR_INT) return a==b;
2421  return _nlEqual_aNoImm_OR_bNoImm(a, b);
2422 }
2423 
2424 LINLINE number nlInit (long i, const coeffs r)
2425 {
2426  number n;
2427  #if MAX_NUM_SIZE == 60
2428  if (((i << 3) >> 3) == i) n=INT_TO_SR(i);
2429  else n=nlRInit(i);
2430  #else
2431  LONG ii=(LONG)i;
2432  if ( ((((long)ii)==i) && ((ii << 3) >> 3) == ii )) n=INT_TO_SR(ii);
2433  else n=nlRInit(i);
2434  #endif
2435  nlTest(n, r);
2436  return n;
2437 }
2438 
2439 /*2
2440 * a == 1 ?
2441 */
2442 LINLINE BOOLEAN nlIsOne (number a, const coeffs r)
2443 {
2444 #ifdef LDEBUG
2445  if (a==NULL) return FALSE;
2446  nlTest(a, r);
2447 #endif
2448  return (a==INT_TO_SR(1));
2449 }
2450 
2452 {
2453  #if 0
2454  if (a==INT_TO_SR(0)) return TRUE;
2455  if ((SR_HDL(a) & SR_INT)||(a==NULL)) return FALSE;
2456  if (mpz_cmp_si(a->z,0L)==0)
2457  {
2458  printf("gmp-0 in nlIsZero\n");
2459  dErrorBreak();
2460  return TRUE;
2461  }
2462  return FALSE;
2463  #else
2464  return (a==INT_TO_SR(0));
2465  #endif
2466 }
2467 
2468 /*2
2469 * copy a to b
2470 */
2471 LINLINE number nlCopy(number a, const coeffs)
2472 {
2473  if ((SR_HDL(a) & SR_INT)||(a==NULL))
2474  {
2475  return a;
2476  }
2477  return _nlCopy_NoImm(a);
2478 }
2479 
2480 
2481 /*2
2482 * delete a
2483 */
2484 LINLINE void nlDelete (number * a, const coeffs r)
2485 {
2486  if (*a!=NULL)
2487  {
2488  nlTest(*a, r);
2489  if ((SR_HDL(*a) & SR_INT)==0)
2490  {
2491  _nlDelete_NoImm(a);
2492  }
2493  *a=NULL;
2494  }
2495 }
2496 
2497 /*2
2498 * za:= - za
2499 */
2500 LINLINE number nlNeg (number a, const coeffs R)
2501 {
2502  nlTest(a, R);
2503  if(SR_HDL(a) &SR_INT)
2504  {
2505  LONG r=SR_TO_INT(a);
2506  if (r==(-(POW_2_28))) a=nlRInit(POW_2_28);
2507  else a=INT_TO_SR(-r);
2508  return a;
2509  }
2510  a = _nlNeg_NoImm(a);
2511  nlTest(a, R);
2512  return a;
2513 
2514 }
2515 
2516 /*2
2517 * u:= a + b
2518 */
2519 LINLINE number nlAdd (number a, number b, const coeffs R)
2520 {
2521  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
2522  {
2523  LONG r=SR_HDL(a)+SR_HDL(b)-1L;
2524  if ( ((r << 1) >> 1) == r )
2525  return (number)(long)r;
2526  else
2527  return nlRInit(SR_TO_INT(r));
2528  }
2529  number u = _nlAdd_aNoImm_OR_bNoImm(a, b);
2530  nlTest(u, R);
2531  return u;
2532 }
2533 
2534 number nlShort1(number a);
2535 number nlShort3_noinline(number x);
2536 
2537 LINLINE void nlInpAdd(number &a, number b, const coeffs r)
2538 {
2539  // a=a+b
2540  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
2541  {
2542  LONG r=SR_HDL(a)+SR_HDL(b)-1L;
2543  if ( ((r << 1) >> 1) == r )
2544  a=(number)(long)r;
2545  else
2546  a=nlRInit(SR_TO_INT(r));
2547  }
2548  else
2549  {
2551  nlTest(a,r);
2552  }
2553 }
2554 
2555 LINLINE number nlMult (number a, number b, const coeffs R)
2556 {
2557  nlTest(a, R);
2558  nlTest(b, R);
2559  if (a==INT_TO_SR(0)) return INT_TO_SR(0);
2560  if (b==INT_TO_SR(0)) return INT_TO_SR(0);
2561  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
2562  {
2563  LONG r=(LONG)((unsigned LONG)(SR_HDL(a)-1L))*((unsigned LONG)(SR_HDL(b)>>1));
2564  if ((r/(SR_HDL(b)>>1))==(SR_HDL(a)-1L))
2565  {
2566  number u=((number) ((r>>1)+SR_INT));
2567  if (((((LONG)SR_HDL(u))<<1)>>1)==SR_HDL(u)) return (u);
2568  return nlRInit(SR_HDL(u)>>2);
2569  }
2570  number u = _nlMult_aImm_bImm_rNoImm(a, b);
2571  nlTest(u, R);
2572  return u;
2573 
2574  }
2575  number u = _nlMult_aNoImm_OR_bNoImm(a, b);
2576  nlTest(u, R);
2577  return u;
2578 
2579 }
2580 
2581 
2582 /*2
2583 * u:= a - b
2584 */
2585 LINLINE number nlSub (number a, number b, const coeffs r)
2586 {
2587  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
2588  {
2589  LONG r=SR_HDL(a)-SR_HDL(b)+1;
2590  if ( ((r << 1) >> 1) == r )
2591  {
2592  return (number)(long)r;
2593  }
2594  else
2595  return nlRInit(SR_TO_INT(r));
2596  }
2597  number u = _nlSub_aNoImm_OR_bNoImm(a, b);
2598  nlTest(u, r);
2599  return u;
2600 
2601 }
2602 
2603 LINLINE void nlInpMult(number &a, number b, const coeffs r)
2604 {
2605  number aa=a;
2606  if (((SR_HDL(b)|SR_HDL(aa))&SR_INT))
2607  {
2608  number n=nlMult(aa,b,r);
2609  nlDelete(&a,r);
2610  a=n;
2611  }
2612  else
2613  {
2614  mpz_mul(aa->z,a->z,b->z);
2615  if (aa->s==3)
2616  {
2617  if(b->s!=3)
2618  {
2619  mpz_init_set(a->n,b->n);
2620  a->s=0;
2621  }
2622  }
2623  else
2624  {
2625  if(b->s!=3)
2626  {
2627  mpz_mul(a->n,a->n,b->n);
2628  }
2629  a->s=0;
2630  }
2631  }
2632 }
2633 #endif // DO_LINLINE
2634 
2635 #ifndef P_NUMBERS_H
2636 
2637 static void nlMPZ(mpz_t m, number &n, const coeffs r)
2638 {
2639  nlTest(n, r);
2640  nlNormalize(n, r);
2641  if (SR_HDL(n) & SR_INT) mpz_init_set_si(m, SR_TO_INT(n)); /* n fits in an int */
2642  else mpz_init_set(m, (mpz_ptr)n->z);
2643 }
2644 
2645 
2646 static number nlInitMPZ(mpz_t m, const coeffs)
2647 {
2648  number z = ALLOC_RNUMBER();
2649  z->s = 3;
2650  #ifdef LDEBUG
2651  z->debug=123456;
2652  #endif
2653  mpz_init_set(z->z, m);
2654  z=nlShort3(z);
2655  return z;
2656 }
2657 
2658 number nlXExtGcd (number a, number b, number *s, number *t, number *u, number *v, const coeffs r)
2659 {
2660  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
2661  {
2662  int uu, vv, x, y;
2663  int g = int_extgcd(SR_TO_INT(a), SR_TO_INT(b), &uu, &vv, &x, &y);
2664  *s = INT_TO_SR(uu);
2665  *t = INT_TO_SR(vv);
2666  *u = INT_TO_SR(x);
2667  *v = INT_TO_SR(y);
2668  return INT_TO_SR(g);
2669  }
2670  else
2671  {
2672  mpz_t aa, bb;
2673  if (SR_HDL(a) & SR_INT)
2674  {
2675  mpz_init_set_si(aa, SR_TO_INT(a));
2676  }
2677  else
2678  {
2679  mpz_init_set(aa, a->z);
2680  }
2681  if (SR_HDL(b) & SR_INT)
2682  {
2683  mpz_init_set_si(bb, SR_TO_INT(b));
2684  }
2685  else
2686  {
2687  mpz_init_set(bb, b->z);
2688  }
2689  mpz_t erg; mpz_t bs; mpz_t bt;
2690  mpz_init(erg);
2691  mpz_init(bs);
2692  mpz_init(bt);
2693 
2694  mpz_gcdext(erg, bs, bt, aa, bb);
2695 
2696  mpz_div(aa, aa, erg);
2697  *u=nlInitMPZ(bb,r);
2698  *u=nlNeg(*u,r);
2699  *v=nlInitMPZ(aa,r);
2700 
2701  mpz_clear(aa);
2702  mpz_clear(bb);
2703 
2704  *s = nlInitMPZ(bs,r);
2705  *t = nlInitMPZ(bt,r);
2706  return nlInitMPZ(erg,r);
2707  }
2708 }
2709 
2710 number nlQuotRem (number a, number b, number * r, const coeffs R)
2711 {
2712  assume(SR_TO_INT(b)!=0);
2713  if (SR_HDL(a) & SR_HDL(b) & SR_INT)
2714  {
2715  if (r!=NULL)
2716  *r = INT_TO_SR(SR_TO_INT(a) % SR_TO_INT(b));
2717  return INT_TO_SR(SR_TO_INT(a)/SR_TO_INT(b));
2718  }
2719  else if (SR_HDL(a) & SR_INT)
2720  {
2721  // -2^xx / 2^xx
2722  if ((a==INT_TO_SR(-(POW_2_28)))&&(b==INT_TO_SR(-1L)))
2723  {
2724  if (r!=NULL) *r=INT_TO_SR(0);
2725  return nlRInit(POW_2_28);
2726  }
2727  //a is small, b is not, so q=0, r=a
2728  if (r!=NULL)
2729  *r = a;
2730  return INT_TO_SR(0);
2731  }
2732  else if (SR_HDL(b) & SR_INT)
2733  {
2734  unsigned long rr;
2735  mpz_t qq;
2736  mpz_init(qq);
2737  mpz_t rrr;
2738  mpz_init(rrr);
2739  rr = mpz_divmod_ui(qq, rrr, a->z, (unsigned long)ABS(SR_TO_INT(b)));
2740  mpz_clear(rrr);
2741 
2742  if (r!=NULL)
2743  *r = INT_TO_SR(rr);
2744  if (SR_TO_INT(b)<0)
2745  {
2746  mpz_mul_si(qq, qq, -1);
2747  }
2748  return nlInitMPZ(qq,R);
2749  }
2750  mpz_t qq,rr;
2751  mpz_init(qq);
2752  mpz_init(rr);
2753  mpz_divmod(qq, rr, a->z, b->z);
2754  if (r!=NULL)
2755  *r = nlInitMPZ(rr,R);
2756  else
2757  {
2758  mpz_clear(rr);
2759  }
2760  return nlInitMPZ(qq,R);
2761 }
2762 
2763 void nlInpGcd(number &a, number b, const coeffs r)
2764 {
2765  if ((SR_HDL(b)|SR_HDL(a))&SR_INT)
2766  {
2767  number n=nlGcd(a,b,r);
2768  nlDelete(&a,r);
2769  a=n;
2770  }
2771  else
2772  {
2773  mpz_gcd(a->z,a->z,b->z);
2774  a=nlShort3_noinline(a);
2775  }
2776 }
2777 
2778 void nlInpIntDiv(number &a, number b, const coeffs r)
2779 {
2780  if ((SR_HDL(b)|SR_HDL(a))&SR_INT)
2781  {
2782  number n=nlIntDiv(a,b, r);
2783  nlDelete(&a,r);
2784  a=n;
2785  }
2786  else
2787  {
2788  number rr=nlIntMod(a,b,r);
2789  if (SR_HDL(rr) & SR_INT) mpz_sub_ui(a->z,a->z,SR_TO_INT(rr));
2790  else mpz_sub(a->z,a->z,rr->z);
2791  mpz_divexact(a->z,a->z,b->z);
2792  a=nlShort3_noinline(a);
2793  }
2794 }
2795 
2796 number nlFarey(number nN, number nP, const coeffs r)
2797 {
2798  mpz_t tmp; mpz_init(tmp);
2799  mpz_t A,B,C,D,E,N,P;
2800  if (SR_HDL(nN) & SR_INT) mpz_init_set_si(N,SR_TO_INT(nN));
2801  else mpz_init_set(N,nN->z);
2802  if (SR_HDL(nP) & SR_INT) mpz_init_set_si(P,SR_TO_INT(nP));
2803  else mpz_init_set(P,nP->z);
2804  assume(!mpz_isNeg(P));
2805  if (mpz_isNeg(N)) mpz_add(N,N,P);
2806  mpz_init_set_si(A,0L);
2807  mpz_init_set_ui(B,(unsigned long)1);
2808  mpz_init_set_si(C,0L);
2809  mpz_init(D);
2810  mpz_init_set(E,P);
2811  number z=INT_TO_SR(0);
2812  while(mpz_cmp_si(N,0L)!=0)
2813  {
2814  mpz_mul(tmp,N,N);
2815  mpz_add(tmp,tmp,tmp);
2816  if (mpz_cmp(tmp,P)<0)
2817  {
2818  if (mpz_isNeg(B))
2819  {
2820  mpz_neg(B,B);
2821  mpz_neg(N,N);
2822  }
2823  // check for gcd(N,B)==1
2824  mpz_gcd(tmp,N,B);
2825  if (mpz_cmp_ui(tmp,1)==0)
2826  {
2827  // return N/B
2828  z=ALLOC_RNUMBER();
2829  #ifdef LDEBUG
2830  z->debug=123456;
2831  #endif
2832  mpz_init_set(z->z,N);
2833  mpz_init_set(z->n,B);
2834  z->s = 0;
2835  nlNormalize(z,r);
2836  }
2837  else
2838  {
2839  // return nN (the input) instead of "fail"
2840  z=nlCopy(nN,r);
2841  }
2842  break;
2843  }
2844  //mpz_mod(D,E,N);
2845  //mpz_div(tmp,E,N);
2846  mpz_divmod(tmp,D,E,N);
2847  mpz_mul(tmp,tmp,B);
2848  mpz_sub(C,A,tmp);
2849  mpz_set(E,N);
2850  mpz_set(N,D);
2851  mpz_set(A,B);
2852  mpz_set(B,C);
2853  }
2854  mpz_clear(tmp);
2855  mpz_clear(A);
2856  mpz_clear(B);
2857  mpz_clear(C);
2858  mpz_clear(D);
2859  mpz_clear(E);
2860  mpz_clear(N);
2861  mpz_clear(P);
2862  return z;
2863 }
2864 
2865 number nlExtGcd(number a, number b, number *s, number *t, const coeffs)
2866 {
2867  mpz_ptr aa,bb;
2868  *s=ALLOC_RNUMBER();
2869  mpz_init((*s)->z); (*s)->s=3;
2870  (*t)=ALLOC_RNUMBER();
2871  mpz_init((*t)->z); (*t)->s=3;
2872  number g=ALLOC_RNUMBER();
2873  mpz_init(g->z); g->s=3;
2874  #ifdef LDEBUG
2875  g->debug=123456;
2876  (*s)->debug=123456;
2877  (*t)->debug=123456;
2878  #endif
2879  if (SR_HDL(a) & SR_INT)
2880  {
2881  aa=(mpz_ptr)omAlloc(sizeof(mpz_t));
2882  mpz_init_set_si(aa,SR_TO_INT(a));
2883  }
2884  else
2885  {
2886  aa=a->z;
2887  }
2888  if (SR_HDL(b) & SR_INT)
2889  {
2890  bb=(mpz_ptr)omAlloc(sizeof(mpz_t));
2891  mpz_init_set_si(bb,SR_TO_INT(b));
2892  }
2893  else
2894  {
2895  bb=b->z;
2896  }
2897  mpz_gcdext(g->z,(*s)->z,(*t)->z,aa,bb);
2898  g=nlShort3(g);
2899  (*s)=nlShort3((*s));
2900  (*t)=nlShort3((*t));
2901  if (SR_HDL(a) & SR_INT)
2902  {
2903  mpz_clear(aa);
2904  omFreeSize(aa, sizeof(mpz_t));
2905  }
2906  if (SR_HDL(b) & SR_INT)
2907  {
2908  mpz_clear(bb);
2909  omFreeSize(bb, sizeof(mpz_t));
2910  }
2911  return g;
2912 }
2913 
2914 void nlCoeffWrite (const coeffs r, BOOLEAN /*details*/)
2915 {
2916  if (r->is_field)
2917  PrintS("// characteristic : 0\n");
2918  else
2919  PrintS("// coeff. ring is : Integers\n");
2920 }
2921 
2923 number nlChineseRemainderSym(number *x, number *q,int rl, BOOLEAN sym, CFArray &inv_cache,const coeffs CF)
2924 // elemenst in the array are x[0..(rl-1)], q[0..(rl-1)]
2925 {
2926  setCharacteristic( 0 ); // only in char 0
2927  Off(SW_RATIONAL);
2928  CFArray X(rl), Q(rl);
2929  int i;
2930  for(i=rl-1;i>=0;i--)
2931  {
2932  X[i]=CF->convSingNFactoryN(x[i],FALSE,CF); // may be larger MAX_INT
2933  Q[i]=CF->convSingNFactoryN(q[i],FALSE,CF); // may be larger MAX_INT
2934  }
2935  CanonicalForm xnew,qnew;
2936  if (n_SwitchChinRem)
2937  chineseRemainder(X,Q,xnew,qnew);
2938  else
2939  chineseRemainderCached(X,Q,xnew,qnew,inv_cache);
2940  number n=CF->convFactoryNSingN(xnew,CF);
2941  if (sym)
2942  {
2943  number p=CF->convFactoryNSingN(qnew,CF);
2944  number p2=nlIntDiv(p,nlInit(2, CF),CF);
2945  if (nlGreater(n,p2,CF))
2946  {
2947  number n2=nlSub(n,p,CF);
2948  nlDelete(&n,CF);
2949  n=n2;
2950  }
2951  nlDelete(&p2,CF);
2952  nlDelete(&p,CF);
2953  }
2954  nlNormalize(n,CF);
2955  return n;
2956 }
2957 number nlChineseRemainder(number *x, number *q,int rl, const coeffs C)
2958 {
2959  CFArray inv(rl);
2960  return nlChineseRemainderSym(x,q,rl,TRUE,inv,C);
2961 }
2962 
2963 static void nlClearContent(ICoeffsEnumerator& numberCollectionEnumerator, number& c, const coeffs cf)
2964 {
2965  assume(cf != NULL);
2966 
2967  numberCollectionEnumerator.Reset();
2968 
2969  if( !numberCollectionEnumerator.MoveNext() ) // empty zero polynomial?
2970  {
2971  c = nlInit(1, cf);
2972  return;
2973  }
2974 
2975  // all coeffs are given by integers!!!
2976 
2977  // part 1, find a small candidate for gcd
2978  number cand1,cand;
2979  int s1,s;
2980  s=2147483647; // max. int
2981 
2982  const BOOLEAN lc_is_pos=nlGreaterZero(numberCollectionEnumerator.Current(),cf);
2983 
2984  int normalcount = 0;
2985  do
2986  {
2987  number& n = numberCollectionEnumerator.Current();
2988  nlNormalize(n, cf); ++normalcount;
2989  cand1 = n;
2990 
2991  if (SR_HDL(cand1)&SR_INT) { cand=cand1; break; }
2992  assume(cand1->s==3); // all coeffs should be integers // ==0?!! after printing
2993  s1=mpz_size1(cand1->z);
2994  if (s>s1)
2995  {
2996  cand=cand1;
2997  s=s1;
2998  }
2999  } while (numberCollectionEnumerator.MoveNext() );
3000 
3001 // assume( nlGreaterZero(cand,cf) ); // cand may be a negative integer!
3002 
3003  cand=nlCopy(cand,cf);
3004  // part 2: compute gcd(cand,all coeffs)
3005 
3006  numberCollectionEnumerator.Reset();
3007 
3008  while (numberCollectionEnumerator.MoveNext() )
3009  {
3010  number& n = numberCollectionEnumerator.Current();
3011 
3012  if( (--normalcount) <= 0)
3013  nlNormalize(n, cf);
3014 
3015  nlInpGcd(cand, n, cf);
3016  assume( nlGreaterZero(cand,cf) );
3017 
3018  if(nlIsOne(cand,cf))
3019  {
3020  c = cand;
3021 
3022  if(!lc_is_pos)
3023  {
3024  // make the leading coeff positive
3025  c = nlNeg(c, cf);
3026  numberCollectionEnumerator.Reset();
3027 
3028  while (numberCollectionEnumerator.MoveNext() )
3029  {
3030  number& nn = numberCollectionEnumerator.Current();
3031  nn = nlNeg(nn, cf);
3032  }
3033  }
3034  return;
3035  }
3036  }
3037 
3038  // part3: all coeffs = all coeffs / cand
3039  if (!lc_is_pos)
3040  cand = nlNeg(cand,cf);
3041 
3042  c = cand;
3043  numberCollectionEnumerator.Reset();
3044 
3045  while (numberCollectionEnumerator.MoveNext() )
3046  {
3047  number& n = numberCollectionEnumerator.Current();
3048  number t=nlIntDiv(n, cand, cf); // simple integer exact division, no ratios to remain
3049  nlDelete(&n, cf);
3050  n = t;
3051  }
3052 }
3053 
3054 static void nlClearDenominators(ICoeffsEnumerator& numberCollectionEnumerator, number& c, const coeffs cf)
3055 {
3056  assume(cf != NULL);
3057 
3058  numberCollectionEnumerator.Reset();
3059 
3060  if( !numberCollectionEnumerator.MoveNext() ) // empty zero polynomial?
3061  {
3062  c = nlInit(1, cf);
3063 // assume( n_GreaterZero(c, cf) );
3064  return;
3065  }
3066 
3067  // all coeffs are given by integers after returning from this routine
3068 
3069  // part 1, collect product of all denominators /gcds
3070  number cand;
3071  cand=ALLOC_RNUMBER();
3072 #if defined(LDEBUG)
3073  cand->debug=123456;
3074 #endif
3075  cand->s=3;
3076 
3077  int s=0;
3078 
3079  const BOOLEAN lc_is_pos=nlGreaterZero(numberCollectionEnumerator.Current(),cf);
3080 
3081  do
3082  {
3083  number& cand1 = numberCollectionEnumerator.Current();
3084 
3085  if (!(SR_HDL(cand1)&SR_INT))
3086  {
3087  nlNormalize(cand1, cf);
3088  if ((!(SR_HDL(cand1)&SR_INT)) // not a short int
3089  && (cand1->s==1)) // and is a normalised rational
3090  {
3091  if (s==0) // first denom, we meet
3092  {
3093  mpz_init_set(cand->z, cand1->n); // cand->z = cand1->n
3094  s=1;
3095  }
3096  else // we have already something
3097  {
3098  mpz_lcm(cand->z, cand->z, cand1->n);
3099  }
3100  }
3101  }
3102  }
3103  while (numberCollectionEnumerator.MoveNext() );
3104 
3105 
3106  if (s==0) // nothing to do, all coeffs are already integers
3107  {
3108 // mpz_clear(tmp);
3109  FREE_RNUMBER(cand);
3110  if (lc_is_pos)
3111  c=nlInit(1,cf);
3112  else
3113  {
3114  // make the leading coeff positive
3115  c=nlInit(-1,cf);
3116 
3117  // TODO: incorporate the following into the loop below?
3118  numberCollectionEnumerator.Reset();
3119  while (numberCollectionEnumerator.MoveNext() )
3120  {
3121  number& n = numberCollectionEnumerator.Current();
3122  n = nlNeg(n, cf);
3123  }
3124  }
3125 // assume( n_GreaterZero(c, cf) );
3126  return;
3127  }
3128 
3129  cand = nlShort3(cand);
3130 
3131  // part2: all coeffs = all coeffs * cand
3132  // make the lead coeff positive
3133  numberCollectionEnumerator.Reset();
3134 
3135  if (!lc_is_pos)
3136  cand = nlNeg(cand, cf);
3137 
3138  c = cand;
3139 
3140  while (numberCollectionEnumerator.MoveNext() )
3141  {
3142  number &n = numberCollectionEnumerator.Current();
3143  nlInpMult(n, cand, cf);
3144  }
3145 
3146 }
3147 
3148 char * nlCoeffName(const coeffs r)
3149 {
3150  if (r->cfDiv==nlDiv) return (char*)"QQ";
3151  else return (char*)"ZZ";
3152 }
3153 
3154 static char* nlCoeffString(const coeffs r)
3155 {
3156  //return omStrDup(nlCoeffName(r));
3157  if (r->cfDiv==nlDiv) return omStrDup("QQ");
3158  else return omStrDup("ZZ");
3159 }
3160 
3161 static void nlWriteFd(number n,FILE* f, const coeffs)
3162 {
3163  if(SR_HDL(n) & SR_INT)
3164  {
3165  #if SIZEOF_LONG == 4
3166  fprintf(f,"4 %ld ",SR_TO_INT(n));
3167  #else
3168  long nn=SR_TO_INT(n);
3169  if ((nn<POW_2_28_32)&&(nn>= -POW_2_28_32))
3170  {
3171  int nnn=(int)nn;
3172  fprintf(f,"4 %d ",nnn);
3173  }
3174  else
3175  {
3176  mpz_t tmp;
3177  mpz_init_set_si(tmp,nn);
3178  fputs("8 ",f);
3179  mpz_out_str (f,SSI_BASE, tmp);
3180  fputc(' ',f);
3181  mpz_clear(tmp);
3182  }
3183  #endif
3184  }
3185  else if (n->s<2)
3186  {
3187  //gmp_fprintf(f,"%d %Zd %Zd ",n->s,n->z,n->n);
3188  fprintf(f,"%d ",n->s+5);
3189  mpz_out_str (f,SSI_BASE, n->z);
3190  fputc(' ',f);
3191  mpz_out_str (f,SSI_BASE, n->n);
3192  fputc(' ',f);
3193 
3194  //if (d->f_debug!=NULL) gmp_fprintf(d->f_debug,"number: s=%d gmp/gmp \"%Zd %Zd\" ",n->s,n->z,n->n);
3195  }
3196  else /*n->s==3*/
3197  {
3198  //gmp_fprintf(d->f_write,"3 %Zd ",n->z);
3199  fputs("8 ",f);
3200  mpz_out_str (f,SSI_BASE, n->z);
3201  fputc(' ',f);
3202 
3203  //if (d->f_debug!=NULL) gmp_fprintf(d->f_debug,"number: gmp \"%Zd\" ",n->z);
3204  }
3205 }
3206 
3207 static number nlReadFd(s_buff f, const coeffs)
3208 {
3209  int sub_type=-1;
3210  sub_type=s_readint(f);
3211  switch(sub_type)
3212  {
3213  case 0:
3214  case 1:
3215  {// read mpz_t, mpz_t
3216  number n=nlRInit(0);
3217  mpz_init(n->n);
3218  s_readmpz(f,n->z);
3219  s_readmpz(f,n->n);
3220  n->s=sub_type;
3221  return n;
3222  }
3223 
3224  case 3:
3225  {// read mpz_t
3226  number n=nlRInit(0);
3227  s_readmpz(f,n->z);
3228  n->s=3; /*sub_type*/
3229  #if SIZEOF_LONG == 8
3230  n=nlShort3(n);
3231  #endif
3232  return n;
3233  }
3234  case 4:
3235  {
3236  LONG dd=s_readlong(f);
3237  //#if SIZEOF_LONG == 8
3238  return INT_TO_SR(dd);
3239  //#else
3240  //return nlInit(dd,NULL);
3241  //#endif
3242  }
3243  case 5:
3244  case 6:
3245  {// read raw mpz_t, mpz_t
3246  number n=nlRInit(0);
3247  mpz_init(n->n);
3248  s_readmpz_base (f,n->z, SSI_BASE);
3249  s_readmpz_base (f,n->n, SSI_BASE);
3250  n->s=sub_type-5;
3251  return n;
3252  }
3253  case 8:
3254  {// read raw mpz_t
3255  number n=nlRInit(0);
3256  s_readmpz_base (f,n->z, SSI_BASE);
3257  n->s=sub_type=3; /*subtype-5*/
3258  #if SIZEOF_LONG == 8
3259  n=nlShort3(n);
3260  #endif
3261  return n;
3262  }
3263 
3264  default: Werror("error in reading number: invalid subtype %d",sub_type);
3265  return NULL;
3266  }
3267  return NULL;
3268 }
3270 {
3271  /* test, if r is an instance of nInitCoeffs(n,parameter) */
3272  /* if parameter is not needed */
3273  if (n==r->type)
3274  {
3275  if ((p==NULL)&&(r->cfDiv==nlDiv)) return TRUE;
3276  if ((p!=NULL)&&(r->cfDiv!=nlDiv)) return TRUE;
3277  }
3278  return FALSE;
3279 }
3280 
3281 static number nlLcm(number a,number b,const coeffs r)
3282 {
3283  number g=nlGcd(a,b,r);
3284  number n1=nlMult(a,b,r);
3285  number n2=nlIntDiv(n1,g,r);
3286  nlDelete(&g,r);
3287  nlDelete(&n1,r);
3288  return n2;
3289 }
3290 
3291 static number nlRandom(siRandProc p, number v2, number, const coeffs cf)
3292 {
3293  number a=nlInit(p(),cf);
3294  if (v2!=NULL)
3295  {
3296  number b=nlInit(p(),cf);
3297  number c=nlDiv(a,b,cf);
3298  nlDelete(&b,cf);
3299  nlDelete(&a,cf);
3300  a=c;
3301  }
3302  return a;
3303 }
3304 
3306 {
3307  r->is_domain=TRUE;
3308  r->rep=n_rep_gap_rat;
3309 
3310  //const int ch = (int)(long)(p);
3311 
3312  r->nCoeffIsEqual=nlCoeffIsEqual;
3313  //r->cfKillChar = ndKillChar; /* dummy */
3314  r->cfCoeffString=nlCoeffString;
3315  r->cfCoeffName=nlCoeffName;
3316 
3317  r->cfInitMPZ = nlInitMPZ;
3318  r->cfMPZ = nlMPZ;
3319 
3320  r->cfMult = nlMult;
3321  r->cfSub = nlSub;
3322  r->cfAdd = nlAdd;
3323  if (p==NULL) /* Q */
3324  {
3325  r->is_field=TRUE;
3326  r->cfDiv = nlDiv;
3327  //r->cfGcd = ndGcd_dummy;
3328  r->cfSubringGcd = nlGcd;
3329  }
3330  else /* Z: coeffs_BIGINT */
3331  {
3332  r->is_field=FALSE;
3333  r->cfDiv = nlIntDiv;
3334  r->cfIntMod= nlIntMod;
3335  r->cfGcd = nlGcd;
3336  r->cfDivBy=nlDivBy;
3337  r->cfDivComp = nlDivComp;
3338  r->cfIsUnit = nlIsUnit;
3339  r->cfGetUnit = nlGetUnit;
3340  r->cfQuot1 = nlQuot1;
3341  r->cfLcm = nlLcm;
3342  r->cfXExtGcd=nlXExtGcd;
3343  r->cfQuotRem=nlQuotRem;
3344  }
3345  r->cfExactDiv= nlExactDiv;
3346  r->cfInit = nlInit;
3347  r->cfSize = nlSize;
3348  r->cfInt = nlInt;
3349 
3350  r->cfChineseRemainder=nlChineseRemainderSym;
3351  r->cfFarey=nlFarey;
3352  r->cfInpNeg = nlNeg;
3353  r->cfInvers= nlInvers;
3354  r->cfCopy = nlCopy;
3355  r->cfRePart = nlCopy;
3356  //r->cfImPart = ndReturn0;
3357  r->cfWriteLong = nlWrite;
3358  r->cfRead = nlRead;
3359  r->cfNormalize=nlNormalize;
3360  r->cfGreater = nlGreater;
3361  r->cfEqual = nlEqual;
3362  r->cfIsZero = nlIsZero;
3363  r->cfIsOne = nlIsOne;
3364  r->cfIsMOne = nlIsMOne;
3365  r->cfGreaterZero = nlGreaterZero;
3366  r->cfPower = nlPower;
3367  r->cfGetDenom = nlGetDenom;
3368  r->cfGetNumerator = nlGetNumerator;
3369  r->cfExtGcd = nlExtGcd; // only for ring stuff and Z
3370  r->cfNormalizeHelper = nlNormalizeHelper;
3371  r->cfDelete= nlDelete;
3372  r->cfSetMap = nlSetMap;
3373  //r->cfName = ndName;
3374  r->cfInpMult=nlInpMult;
3375  r->cfInpAdd=nlInpAdd;
3376  r->cfCoeffWrite=nlCoeffWrite;
3377 
3378  r->cfClearContent = nlClearContent;
3379  r->cfClearDenominators = nlClearDenominators;
3380 
3381 #ifdef LDEBUG
3382  // debug stuff
3383  r->cfDBTest=nlDBTest;
3384 #endif
3385  r->convSingNFactoryN=nlConvSingNFactoryN;
3386  r->convFactoryNSingN=nlConvFactoryNSingN;
3387 
3388  r->cfRandom=nlRandom;
3389 
3390  // io via ssi
3391  r->cfWriteFd=nlWriteFd;
3392  r->cfReadFd=nlReadFd;
3393 
3394  // the variables: general stuff (required)
3395  r->nNULL = INT_TO_SR(0);
3396  //r->type = n_Q;
3397  r->ch = 0;
3398  r->has_simple_Alloc=FALSE;
3399  r->has_simple_Inverse=FALSE;
3400 
3401  // variables for this type of coeffs:
3402  // (none)
3403  return FALSE;
3404 }
3405 #if 0
3406 number nlMod(number a, number b)
3407 {
3408  if (((SR_HDL(b)&SR_HDL(a))&SR_INT)
3409  {
3410  int bi=SR_TO_INT(b);
3411  int ai=SR_TO_INT(a);
3412  int bb=ABS(bi);
3413  int c=ai%bb;
3414  if (c<0) c+=bb;
3415  return (INT_TO_SR(c));
3416  }
3417  number al;
3418  number bl;
3419  if (SR_HDL(a))&SR_INT)
3420  al=nlRInit(SR_TO_INT(a));
3421  else
3422  al=nlCopy(a);
3423  if (SR_HDL(b))&SR_INT)
3424  bl=nlRInit(SR_TO_INT(b));
3425  else
3426  bl=nlCopy(b);
3427  number r=nlRInit(0);
3428  mpz_mod(r->z,al->z,bl->z);
3429  nlDelete(&al);
3430  nlDelete(&bl);
3431  if (mpz_size1(&r->z)<=MP_SMALL)
3432  {
3433  LONG ui=(int)mpz_get_si(&r->z);
3434  if ((((ui<<3)>>3)==ui)
3435  && (mpz_cmp_si(x->z,(long)ui)==0))
3436  {
3437  mpz_clear(&r->z);
3438  FREE_RNUMBER(r); // omFreeBin((void *)r, rnumber_bin);
3439  r=INT_TO_SR(ui);
3440  }
3441  }
3442  return r;
3443 }
3444 #endif
3445 #endif // not P_NUMBERS_H
3446 #endif // LONGRAT_CC
mpz_ptr base
Definition: rmodulon.h:19
mpz_t z
Definition: longrat.h:51
long intval() const
conversion functions
LINLINE number nlSub(number la, number li, const coeffs r)
Definition: longrat.cc:2585
static void nlClearDenominators(ICoeffsEnumerator &numberCollectionEnumerator, number &c, const coeffs cf)
Definition: longrat.cc:3054
const CanonicalForm int s
Definition: facAbsFact.cc:55
const CanonicalForm int const CFList const Variable & y
Definition: facAbsFact.cc:57
CF_NO_INLINE bool isOne() const
CF_INLINE bool CanonicalForm::isOne, isZero () const.
Definition: cf_inline.cc:354
void mpz_mul_si(mpz_ptr r, mpz_srcptr s, long int si)
Definition: longrat.cc:177
#define omCheckAddrSize(addr, size)
Definition: omAllocDecl.h:327
#define D(A)
Definition: gentable.cc:123
#define INT_TO_SR(INT)
Definition: longrat.h:69
const poly a
Definition: syzextra.cc:212
BOOLEAN nlCoeffIsEqual(const coeffs r, n_coeffType n, void *p)
Definition: longrat.cc:3269
#define Print
Definition: emacs.cc:83
only used if HAVE_RINGS is defined
Definition: coeffs.h:44
static FORCE_INLINE BOOLEAN nCoeff_is_Zp(const coeffs r)
Definition: coeffs.h:834
static int int_extgcd(int a, int b, int *u, int *x, int *v, int *y)
Definition: longrat.cc:1263
char * nlCoeffName(const coeffs r)
Definition: longrat.cc:3148
void Off(int sw)
switches
CanonicalForm num(const CanonicalForm &f)
long npInt(number &n, const coeffs r)
Definition: modulop.cc:140
Definition: int_poly.h:33
static number nlConvFactoryNSingN(const CanonicalForm f, const coeffs r)
Definition: longrat.cc:376
BOOLEAN nlGreaterZero(number za, const coeffs r)
Definition: longrat.cc:1155
number nlModP(number q, const coeffs, const coeffs Zp)
Definition: longrat.cc:1425
#define FALSE
Definition: auxiliary.h:94
static void nlMPZ(mpz_t m, number &n, const coeffs r)
Definition: longrat.cc:2637
mpf_t * _mpfp()
Definition: mpr_complex.h:134
return P p
Definition: myNF.cc:203
number _nlMult_aNoImm_OR_bNoImm(number a, number b)
Definition: longrat.cc:2194
f
Definition: cfModGcd.cc:4022
number nlShort1(number x)
Definition: longrat.cc:1313
number nlNormalizeHelper(number a, number b, const coeffs r)
Definition: longrat.cc:1378
number ndCopyMap(number a, const coeffs aRing, const coeffs r)
Definition: numbers.cc:244
bool isImm() const
static FORCE_INLINE BOOLEAN nCoeff_is_long_R(const coeffs r)
Definition: coeffs.h:908
LINLINE void nlInpAdd(number &a, number b, const coeffs r)
Definition: longrat.cc:2537
number nlGetDenom(number &n, const coeffs r)
Definition: longrat.cc:1488
void nlWrite(number a, const coeffs r)
Definition: longrat0.cc:116
int nlSize(number a, const coeffs)
Definition: longrat.cc:576
int n_SwitchChinRem
Definition: longrat.cc:2922
static FORCE_INLINE number n_Init(long i, const coeffs r)
a number representing i in the given coeff field/ring r
Definition: coeffs.h:542
static FORCE_INLINE BOOLEAN nCoeff_is_R(const coeffs r)
Definition: coeffs.h:853
#define omCheckIf(cond, test)
Definition: omAllocDecl.h:323
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260
LINLINE number nlAdd(number la, number li, const coeffs r)
Definition: longrat.cc:2519
{p < 2^31}
Definition: coeffs.h:30
void nlInpGcd(number &a, number b, const coeffs r)
Definition: longrat.cc:2763
BOOLEAN nlIsMOne(number a, const coeffs r)
Definition: longrat.cc:1180
static FORCE_INLINE BOOLEAN nCoeff_is_Ring_2toM(const coeffs r)
Definition: coeffs.h:750
(), see rinteger.h, new impl.
Definition: coeffs.h:112
void nlCoeffWrite(const coeffs r, BOOLEAN details)
Definition: longrat.cc:2914
number nlIntDiv(number a, number b, const coeffs r)
Definition: longrat.cc:786
static FORCE_INLINE int n_GetChar(const coeffs r)
Return the characteristic of the coeff. domain.
Definition: coeffs.h:448
number nlGcd(number a, number b, const coeffs r)
Definition: longrat.cc:1192
factory&#39;s main class
Definition: canonicalform.h:75
#define TRUE
Definition: auxiliary.h:98
coeffs nlQuot1(number c, const coeffs r)
Definition: longrat.cc:958
#define POW_2_28
Definition: longrat.cc:104
number nlRInit(long i)
Definition: longrat.cc:2348
#define FREE_RNUMBER(x)
Definition: coeffs.h:86
number nlInit2(int i, int j, const coeffs r)
create a rational i/j (implicitly) over Q NOTE: make sure to use correct Q in debug mode ...
Definition: longrat.cc:2362
g
Definition: cfModGcd.cc:4031
const char * nlRead(const char *s, number *a, const coeffs r)
Definition: longrat0.cc:57
void WerrorS(const char *s)
Definition: feFopen.cc:24
bool negative(N n)
Definition: ValueTraits.h:119
void s_readmpz_base(s_buff F, mpz_ptr a, int base)
Definition: s_buff.cc:207
CanonicalForm make_cf(const mpz_ptr n)
Definition: singext.cc:67
void nlGMP(number &i, number n, const coeffs r)
Definition: longrat.cc:1467
#define Q
Definition: sirandom.c:25
number nlIntMod(number a, number b, const coeffs r)
Definition: longrat.cc:867
#define POW_2_28_32
Definition: longrat.cc:105
number nlInit2gmp(mpz_t i, mpz_t j, const coeffs r)
create a rational i/j (implicitly) over Q NOTE: make sure to use correct Q in debug mode ...
Definition: longrat.cc:2375
#define WarnS
Definition: emacs.cc:81
BOOLEAN nlInitChar(coeffs r, void *p)
Definition: longrat.cc:3305
#define omAlloc(size)
Definition: omAllocDecl.h:210
void setCharacteristic(int c)
Definition: cf_char.cc:23
BOOLEAN nlGreater(number a, number b, const coeffs r)
Definition: longrat.cc:1165
LINLINE number nl_Copy(number a, const coeffs r)
static number nlInitMPZ(mpz_t m, const coeffs)
Definition: longrat.cc:2646
real floating point (GMP) numbers
Definition: coeffs.h:34
number nlMapZ(number from, const coeffs src, const coeffs dst)
Definition: longrat.cc:219
LINLINE BOOLEAN nlIsOne(number a, const coeffs r)
Definition: longrat.cc:2442
poly res
Definition: myNF.cc:322
virtual void Reset()=0
Sets the enumerator to its initial position: -1, which is before the first element in the collection...
single prescision (6,6) real numbers
Definition: coeffs.h:32
#define MP_SMALL
Definition: longrat.cc:155
static number nlLcm(number a, number b, const coeffs r)
Definition: longrat.cc:3281
mpz_t n
Definition: longrat.h:52
const ring r
Definition: syzextra.cc:208
static number nlMapP(number from, const coeffs src, const coeffs dst)
Definition: longrat.cc:189
LINLINE number nlNeg(number za, const coeffs r)
Definition: longrat.cc:2500
number nlXExtGcd(number a, number b, number *s, number *t, number *u, number *v, const coeffs r)
Definition: longrat.cc:2658
float nrFloat(number n)
Converts a n_R number into a float. Needed by Maps.
Definition: shortfl.cc:75
Coefficient rings, fields and other domains suitable for Singular polynomials.
void s_readmpz(s_buff F, mpz_t a)
Definition: s_buff.cc:182
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:49
int s_readint(s_buff F)
Definition: s_buff.cc:110
int j
Definition: myNF.cc:70
number nlDiv(number a, number b, const coeffs r)
Definition: longrat.cc:992
LINLINE number nlMult(number a, number b, const coeffs r)
Definition: longrat.cc:2555
number nlInvers(number a, const coeffs r)
Definition: longrat.cc:655
number _nlCopy_NoImm(number a)
Definition: longrat.cc:1595
#define assume(x)
Definition: mod2.h:394
void _nlInpAdd_aNoImm_OR_bNoImm(number &a, number b)
Definition: longrat.cc:1827
The main handler for Singular numbers which are suitable for Singular polynomials.
Templated enumerator interface for simple iteration over a generic collection of T&#39;s.
Definition: Enumerator.h:124
int nlDivComp(number a, number b, const coeffs r)
Definition: longrat.cc:942
LINLINE void nlInpMult(number &a, number b, const coeffs r)
Definition: longrat.cc:2603
static void nlWriteFd(number n, FILE *f, const coeffs)
Definition: longrat.cc:3161
const ExtensionInfo & info
< [in] sqrfree poly
#define A
Definition: sirandom.c:23
number(* nMapFunc)(number a, const coeffs src, const coeffs dst)
maps "a", which lives in src, into dst
Definition: coeffs.h:73
const ring R
Definition: DebugPrint.cc:36
void _nlDelete_NoImm(number *a)
Definition: longrat.cc:1616
virtual reference Current()=0
Gets the current element in the collection (read and write).
#define LINLINE
Definition: longrat.cc:26
static number nlMapLongR(number from, const coeffs src, const coeffs dst)
Definition: longrat.cc:433
static const int SW_RATIONAL
set to 1 for computations over Q
Definition: cf_defs.h:28
All the auxiliary stuff.
int m
Definition: cfEzgcd.cc:119
#define mpz_isNeg(A)
Definition: longrat.cc:157
const char *const nDivBy0
Definition: numbers.h:83
void On(int sw)
switches
LINLINE BOOLEAN nlEqual(number a, number b, const coeffs r)
Definition: longrat.cc:2415
unsigned long exp
Definition: rmodulon.h:19
void nlPower(number x, int exp, number *lu, const coeffs r)
Definition: longrat.cc:1102
#define SSI_BASE
Definition: auxiliary.h:131
int i
Definition: cfEzgcd.cc:123
void PrintS(const char *s)
Definition: reporter.cc:284
number nlQuotRem(number a, number b, number *r, const coeffs R)
Definition: longrat.cc:2710
static number nlReadFd(s_buff f, const coeffs)
Definition: longrat.cc:3207
void nlInpIntDiv(number &a, number b, const coeffs r)
Definition: longrat.cc:2778
number nlExtGcd(number a, number b, number *s, number *t, const coeffs)
Definition: longrat.cc:2865
LINLINE BOOLEAN nlIsZero(number za, const coeffs r)
Definition: longrat.cc:2451
int IsPrime(int p)
Definition: prime.cc:61
(mpz_ptr), see rmodulon,h
Definition: coeffs.h:115
static void nlNormalize_Gcd(number &x)
Definition: longrat.cc:1649
BOOLEAN _nlEqual_aNoImm_OR_bNoImm(number a, number b)
Definition: longrat.cc:1548
number nlShort3_noinline(number x)
Definition: longrat.cc:170
static void nlClearContent(ICoeffsEnumerator &numberCollectionEnumerator, number &c, const coeffs cf)
Definition: longrat.cc:2963
static CanonicalForm nlConvSingNFactoryN(number n, const BOOLEAN setChar, const coeffs)
Definition: longrat.cc:338
static FORCE_INLINE n_coeffType getCoeffType(const coeffs r)
Returns the type of coeffs domain.
Definition: coeffs.h:425
int(* siRandProc)()
Definition: sirandom.h:9
number nlGetUnit(number, const coeffs)
Definition: longrat.cc:953
number nlChineseRemainder(number *x, number *q, int rl, const coeffs C)
Definition: longrat.cc:2957
int size(const CanonicalForm &f, const Variable &v)
int size ( const CanonicalForm & f, const Variable & v )
Definition: cf_ops.cc:600
number nlBigInt(number &n)
void chineseRemainderCached(CFArray &a, CFArray &n, CanonicalForm &xnew, CanonicalForm &prod, CFArray &inv)
Definition: cf_chinese.cc:264
#define SR_TO_INT(SR)
Definition: longrat.h:70
void gmp_numerator(const CanonicalForm &f, mpz_ptr result)
Definition: singext.cc:20
(number), see longrat.h
Definition: coeffs.h:111
void nlNormalize(number &x, const coeffs r)
Definition: longrat.cc:1334
void chineseRemainder(const CanonicalForm &x1, const CanonicalForm &q1, const CanonicalForm &x2, const CanonicalForm &q2, CanonicalForm &xnew, CanonicalForm &qnew)
void chineseRemainder ( const CanonicalForm & x1, const CanonicalForm & q1, const CanonicalForm & x2...
Definition: cf_chinese.cc:52
number nlChineseRemainderSym(number *x, number *q, int rl, BOOLEAN sym, CFArray &inv_cache, const coeffs CF)
Definition: longrat.cc:2923
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:37
#define mpz_size1(A)
Definition: si_gmp.h:12
n_coeffType
Definition: coeffs.h:27
CanonicalForm cf
Definition: cfModGcd.cc:4024
number _nlNeg_NoImm(number a)
Definition: longrat.cc:1634
#define NULL
Definition: omList.c:10
static number nlShort3(number x)
Definition: longrat.cc:110
static number nlRandom(siRandProc p, number v2, number, const coeffs cf)
Definition: longrat.cc:3291
CanonicalForm den() const
den() returns the denominator of CO if CO is a rational number, 1 (from the current domain!) otherwis...
REvaluation E(1, terms.length(), IntRandom(25))
CanonicalForm den(const CanonicalForm &f)
LINLINE void nlDelete(number *a, const coeffs r)
Definition: longrat.cc:2484
static FORCE_INLINE number n_Div(number a, number b, const coeffs r)
return the quotient of &#39;a&#39; and &#39;b&#39;, i.e., a/b; raises an error if &#39;b&#39; is not invertible in r exceptio...
Definition: coeffs.h:619
BOOLEAN nlDivBy(number a, number b, const coeffs)
Definition: longrat.cc:928
(gmp_float), see
Definition: coeffs.h:117
b *CanonicalForm B
Definition: facBivar.cc:51
#define ABS(x)
Definition: auxiliary.h:111
virtual bool MoveNext()=0
Advances the enumerator to the next element of the collection. returns true if the enumerator was suc...
int gcd(int a, int b)
Definition: walkSupport.cc:839
BOOLEAN nlIsUnit(number a, const coeffs)
Definition: longrat.cc:983
long nlInt(number &n, const coeffs r)
Definition: longrat.cc:605
#define SR_INT
Definition: longrat.h:68
number _nlMult_aImm_bImm_rNoImm(number a, number b)
Definition: longrat.cc:2181
#define ALLOC_RNUMBER()
Definition: coeffs.h:87
nMapFunc nlSetMap(const coeffs src, const coeffs dst)
Definition: longrat.cc:2311
Variable x
Definition: cfModGcd.cc:4023
number _nlAdd_aNoImm_OR_bNoImm(number a, number b)
Definition: longrat.cc:1669
#define GCD_NORM_COND(OLD, NEW)
Definition: longrat.cc:1647
long s_readlong(s_buff F)
Definition: s_buff.cc:138
number nlExactDiv(number a, number b, const coeffs r)
Definition: longrat.cc:735
number _nlSub_aNoImm_OR_bNoImm(number a, number b)
Definition: longrat.cc:1970
number nlMapGMP(number from, const coeffs src, const coeffs dst)
Definition: longrat.cc:206
BOOLEAN nlDBTest(number a, const char *f, const int l)
p exp[i]
Definition: DebugPrint.cc:39
LINLINE number nlInit(long i, const coeffs r)
Definition: longrat.cc:2424
(int), see modulop.h
Definition: coeffs.h:110
#define SR_HDL(A)
Definition: tgb.cc:35
END_NAMESPACE const void * p2
Definition: syzextra.cc:202
static char * nlCoeffString(const coeffs r)
Definition: longrat.cc:3154
static FORCE_INLINE void n_Delete(number *p, const coeffs r)
delete &#39;p&#39;
Definition: coeffs.h:459
number nlMapMachineInt(number from, const coeffs, const coeffs)
Definition: longrat.cc:231
kBucketDestroy & P
Definition: myNF.cc:191
number nlCopyMap(number a, const coeffs, const coeffs)
Definition: longrat.cc:2302
int BOOLEAN
Definition: auxiliary.h:85
static number nlMapR(number from, const coeffs src, const coeffs dst)
Definition: longrat.cc:403
const poly b
Definition: syzextra.cc:213
number nlGetNumerator(number &n, const coeffs r)
Definition: longrat.cc:1517
LINLINE number nlCopy(number a, const coeffs r)
Definition: longrat.cc:2471
void dErrorBreak()
Definition: dError.cc:141
void Werror(const char *fmt,...)
Definition: reporter.cc:189
#define LONG
Definition: longrat.cc:106
return result
Definition: facAbsBiFact.cc:76
int l
Definition: cfEzgcd.cc:94
#define ALLOC0_RNUMBER()
Definition: coeffs.h:88
#define nlTest(a, r)
Definition: longrat.cc:88
const CanonicalForm const CanonicalForm const CanonicalForm const CanonicalForm & cand
Definition: cfModGcd.cc:69
number nlFarey(number nN, number nP, const coeffs CF)
Definition: longrat.cc:2796
void gmp_denominator(const CanonicalForm &f, mpz_ptr result)
Definition: singext.cc:40
coeffs nInitChar(n_coeffType t, void *parameter)
one-time initialisations for new coeffs in case of an error return NULL
Definition: numbers.cc:334
(float), see shortfl.h
Definition: coeffs.h:116
#define omStrDup(s)
Definition: omAllocDecl.h:263