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