My Project
transext.cc
Go to the documentation of this file.
1 /****************************************
2 * Computer Algebra System SINGULAR *
3 ****************************************/
4 /*
5 * ABSTRACT: numbers in a rational function field K(t_1, .., t_s) with
6 * transcendental variables t_1, ..., t_s, where s >= 1.
7 * Denoting the implemented coeffs object by cf, then these numbers
8 * are represented as quotients of polynomials living in the
9 * polynomial ring K[t_1, .., t_s] represented by cf->extRing.
10 *
11 * An element of K(t_1, .., t_s) may have numerous representations,
12 * due to the possibility of common polynomial factors in the
13 * numerator and denominator. This problem is handled by a
14 * cancellation heuristic: Each number "knows" its complexity
15 * which is 0 if and only if common factors have definitely been
16 * cancelled, and some positive integer otherwise.
17 * Each arithmetic operation of two numbers with complexities c1
18 * and c2 will result in a number of complexity c1 + c2 + some
19 * penalty (specific for each arithmetic operation; see constants
20 * in the *.h file). Whenever the resulting complexity exceeds a
21 * certain threshold (see constant in the *.h file), then the
22 * cancellation heuristic will call 'factory' to compute the gcd
23 * and cancel it out in the given number.
24 * For the special case of K = Q (i.e., when computing over the
25 * rationals), this definite cancellation procedure will also take
26 * care of nested fractions: If there are fractional coefficients
27 * in the numerator or denominator of a number, then this number
28 * is being replaced by a quotient of two polynomials over Z, or
29 * - if the denominator is a constant - by a polynomial over Q.
30 *
31 * TODO: the description above needs a major update!!!
32 */
33 #define TRANSEXT_PRIVATES
34 
35 #include "misc/auxiliary.h"
36 
37 #include "factory/factory.h"
38 
39 #include "reporter/reporter.h"
40 
41 #include "coeffs/coeffs.h"
42 #include "coeffs/numbers.h"
43 
44 #include "coeffs/longrat.h"
45 
46 #include "polys/monomials/ring.h"
48 #include "polys/simpleideals.h"
49 
50 #include "polys/clapsing.h"
51 #include "polys/clapconv.h"
52 
53 #include "polys/prCopy.h"
54 #include "transext.h"
55 #include "algext.h"
56 
57 #include "polys/PolyEnumerator.h"
58 
59 
60 /* constants for controlling the complexity of numbers */
61 #define ADD_COMPLEXITY 1 /**< complexity increase due to + and - */
62 #define MULT_COMPLEXITY 2 /**< complexity increase due to * and / */
63 #define DIFF_COMPLEXITY 2 /**< complexity increase due to diff */
64 #define BOUND_COMPLEXITY 10 /**< maximum complexity of a number */
65 
66 /// TRUE iff num. represents 1
67 #define NUMIS1(f) (p_IsOne(NUM(f), cf->extRing))
68 
69 #define COM(f) (f)->complexity
70 
71 
72 #ifdef LDEBUG
73 static BOOLEAN ntDBTest(number a, const char *f, const int l, const coeffs r);
74 #endif
75 
76 #define ntTest(a) n_Test(a, cf)
77 
78 /* polynomial ring in which the numerators and denominators of our
79  numbers live */
80 #define ntRing cf->extRing
81 
82 /* coeffs object in which the coefficients of our numbers live;
83  * methods attached to ntCoeffs may be used to compute with the
84  * coefficients of our numbers, e.g., use ntCoeffs->nAdd to add
85  * coefficients of our numbers */
86 #define ntCoeffs cf->extRing->cf
87 
88 
89 VAR omBin fractionObjectBin = omGetSpecBin(sizeof(fractionObject));
90 
91 /// forward declarations
92 static void heuristicGcdCancellation(number a, const coeffs cf);
93 static void definiteGcdCancellation(number a, const coeffs cf,
94  BOOLEAN simpleTestsHaveAlreadyBeenPerformed);
95 
96 /* test routine, usualy disabled *
97  * if want to activate it, activate also the calls to check_N *
98  *
99 void check_normalized(number t,const coeffs cf, const char *f, int l)
100 {
101  if (IS0(t)) return;
102  if(rField_is_Q(ntRing))
103  {
104  poly pp=NUM(t);
105  while(pp!=NULL)
106  {
107  if (((SR_HDL(pGetCoeff(pp)) & SR_INT)==0)&&(SR_HDL(pGetCoeff(pp))!=NULL))
108  {
109  if (pGetCoeff(pp)->s==0)
110  {
111  Print("NUM not normalized in %s:%d\n",f,l);
112  p_Normalize(pp,ntRing);
113  }
114  else if (pGetCoeff(pp)->s==1)
115  Print("NUM is rational in %s:%d\n",f,l);
116  }
117  pIter(pp);
118  }
119  pp=DEN(t);
120  while(pp!=NULL)
121  {
122  if (((SR_HDL(pGetCoeff(pp)) & SR_INT)==0)&&(SR_HDL(pGetCoeff(pp))!=NULL))
123  {
124  if (pGetCoeff(pp)->s==0)
125  {
126  Print("NUM not normalized in %s:%d\n",f,l);
127  p_Normalize(pp,ntRing);
128  }
129  else if (pGetCoeff(pp)->s==1)
130  Print("DEN is rational in %s:%d\n",f,l);
131  }
132  pIter(pp);
133  }
134  }
135 }
136 #define check_N(A,B) check_normalized(A,B,__FILE__,__LINE__)
137 */
138 
139 #ifdef LDEBUG
140 static BOOLEAN ntDBTest(number a, const char *f, const int l, const coeffs cf)
141 {
143 
144  if (IS0(a)) return TRUE;
145 
146  const fraction t = (fraction)a;
147 
148  //check_N(a,cf);
149  const poly num = NUM(t);
150  assume(num != NULL); ///< t != 0 ==> numerator(t) != 0
151 
152  p_Test(num, ntRing);
153 
154  if (getCoeffType(ntCoeffs)==n_Q)
155  for( poly p = num; p != NULL; pIter(p) )
156  if (! nlIsInteger( p_GetCoeff(p, ntRing), ntCoeffs) )
157  {
158  Print("ERROR in %s:%d: non-integer Q coeff in num. poly\n",f,l);
159  Print("TERM: "); p_wrp(p, ntRing); PrintLn();
160  return FALSE;
161  }
162 
163  const poly den = DEN(t);
164 
165  if (den != NULL) // !DENIS1(f)
166  {
167  p_Test(den, ntRing);
168 
169  if (getCoeffType(ntCoeffs)==n_Q)
170  for( poly p = den; p != NULL; pIter(p) )
171  if (! nlIsInteger( p_GetCoeff(p, ntRing), ntCoeffs) )
172  {
173  Print("ERROR in %s:%d: non-integer Q coeff in den. poly\n",f,l);
174  Print("TERM: "); p_wrp(p, ntRing); PrintLn();
175  return FALSE;
176  }
177 
178  if (getCoeffType(ntCoeffs)==n_Zp)
179  {
180  if( p_IsConstant(den, ntRing) )
181  {
182  Print("ERROR in %s:%d: constant den. poly / Zp\n",f,l);
183  PrintS("NUM: "); p_Write(num, ntRing);
184  PrintS("DEN: "); p_Write(den, ntRing);
185  return FALSE;
186  }
187 
188  if( !n_IsOne(pGetCoeff(den), ntCoeffs) )
189  {
190  Print("ERROR in %s:%d: non-monic den. poly / Zp\n",f,l);
191  PrintS("NUM: "); p_Write(num, ntRing);
192  PrintS("DEN: "); p_Write(den, ntRing);
193  return FALSE;
194  }
195  }
196 
197  if (COM(t)==0)
198  {
199  poly gcd = singclap_gcd_r( num, den, ntRing );
200  if(gcd!=NULL)
201  {
202  if( !p_IsOne(gcd, ntRing) )
203  {
204  Print("ERROR in %s:%d: 1 != GCD between num. & den. poly\n",f,l);
205  PrintS("GCD: "); p_Write(gcd, ntRing);
206  PrintS("NUM: "); p_Write(num, ntRing);
207  PrintS("DEN: "); p_Write(den, ntRing);
208  return FALSE;
209  }
210  p_Delete( &gcd, ntRing );
211  }
212  }
213  return TRUE;
214 
216  {
217  Print("?/1 in %s:%d\n",f,l);
218  return FALSE;
219  }
221  {
222  Print("negative sign of DEN. of a fraction in %s:%d\n",f,l);
223  return FALSE;
224  }
225  // test that den is over integers!?
226  }
227  else
228  {
229  return TRUE;
230 
231  // num != NULL // den == NULL
232 // if( COM(t) != 0 )
233 // {
234 // Print("?//NULL with non-zero complexity: %d in %s:%d\n", COM(t), f, l);
235 // return FALSE;
236 // }
237  // test that nume is over integers!?
238  }
239  if (getCoeffType(ntCoeffs)==n_Q)
240  {
241  poly p=num; // !=NULL
242  do
243  {
244  number n=pGetCoeff(p);
245  n_Test(n,ntCoeffs);
246  if ((!(SR_HDL(n) & SR_INT))&&(n->s==0))
247  /* not normalized, just do for the following test*/
248  {
250  n=pGetCoeff(p);
251  }
252  if (!(SR_HDL(n) & SR_INT))
253  {
254  if (n->s<2)
255  Print("rational coeff in num: %s:%d\n",f,l);
256  }
257  pIter(p);
258  } while(p!=NULL);
259  p=den;
260  while(p!=NULL)
261  {
262  number n=pGetCoeff(p);
263  if (!(SR_HDL(n) & SR_INT))
264  {
265  if (n->s!=3)
266  Print("rational coeff in den.:%s:%d\n",f,l);
267  }
268  pIter(p);
269  }
270  }
271  return TRUE;
272 }
273 #endif
274 
275 poly gcd_over_Q ( poly f, poly g, const ring r)
276 {
277  poly res;
278  f=p_Copy(f,r);
279  p_Cleardenom(f, r);
280  g=p_Copy(g,r);
281  p_Cleardenom(g, r);
282  res=singclap_gcd_r(f,g,r);
283  p_Delete(&f, r);
284  p_Delete(&g, r);
285  return res;
286 }
287 
288 /* returns the bottom field in this field extension tower; if the tower
289  is flat, i.e., if there is no extension, then r itself is returned;
290  as a side-effect, the counter 'height' is filled with the height of
291  the extension tower (in case the tower is flat, 'height' is zero) */
292 static coeffs nCoeff_bottom(const coeffs r, int &height)
293 {
294  assume(r != NULL);
295  coeffs cf = r;
296  height = 0;
297  while (nCoeff_is_Extension(cf))
298  {
299  assume(cf->extRing != NULL); assume(cf->extRing->cf != NULL);
300  cf = cf->extRing->cf;
301  height++;
302  }
303  return cf;
304 }
305 
306 static BOOLEAN ntIsZero(number a, const coeffs cf)
307 {
308  //check_N(a,cf);
309  ntTest(a); // !!!
310  return (IS0(a));
311 }
312 
313 static void ntDelete(number * a, const coeffs cf)
314 {
315  //check_N(*a,cf);
316  ntTest(*a); // !!!
317 
318  fraction f = (fraction)(*a);
319  if (IS0(f)) return;
320  p_Delete(&NUM(f), ntRing);
321  if (!DENIS1(f)) p_Delete(&DEN(f), ntRing);
323  *a = NULL;
324 }
325 
326 static BOOLEAN ntEqual(number a, number b, const coeffs cf)
327 {
328  //check_N(a,cf);
329  //check_N(b,cf);
330  ntTest(a);
331  ntTest(b);
332 
333  /// simple tests
334  if (a == b) return TRUE;
335  if ((IS0(a)) && (!IS0(b))) return FALSE;
336  if ((IS0(b)) && (!IS0(a))) return FALSE;
337 
338  /// cheap test if gcd's have been cancelled in both numbers
339  fraction fa = (fraction)a;
340  fraction fb = (fraction)b;
341  if ((COM(fa) == 1) && (COM(fb) == 1))
342  {
343  poly f = p_Add_q(p_Copy(NUM(fa), ntRing),
345  ntRing);
346  if (f != NULL) { p_Delete(&f, ntRing); return FALSE; }
347  if (DENIS1(fa) && DENIS1(fb)) return TRUE;
348  if (DENIS1(fa) && !DENIS1(fb)) return FALSE;
349  if (!DENIS1(fa) && DENIS1(fb)) return FALSE;
350  f = p_Add_q(p_Copy(DEN(fa), ntRing),
351  p_Neg(p_Copy(DEN(fb), ntRing), ntRing),
352  ntRing);
353  if (f != NULL) { p_Delete(&f, ntRing); return FALSE; }
354  return TRUE;
355  }
356 
357  /* default: the more expensive multiplication test
358  a/b = c/d <==> a*d = b*c */
359  poly f = p_Copy(NUM(fa), ntRing);
360  if (!DENIS1(fb)) f = p_Mult_q(f, p_Copy(DEN(fb), ntRing), ntRing);
361  poly g = p_Copy(NUM(fb), ntRing);
362  if (!DENIS1(fa)) g = p_Mult_q(g, p_Copy(DEN(fa), ntRing), ntRing);
363  poly h = p_Add_q(f, p_Neg(g, ntRing), ntRing);
364  if (h == NULL) return TRUE;
365  else
366  {
367  p_Delete(&h, ntRing);
368  return FALSE;
369  }
370 }
371 
372 static number ntCopy(number a, const coeffs cf)
373 {
374  //check_N(a,cf);
375  ntTest(a); // !!!
376  if (IS0(a)) return NULL;
377  fraction f = (fraction)a;
378  poly g = NUM(f);
379  poly h = NULL;
380  h =DEN(f);
381  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
382  NUM(result) = p_Copy(g,cf->extRing);
383  DEN(result) = p_Copy(h,cf->extRing);
384  COM(result) = COM(f);
385  ntTest((number)result);
386  return (number)result;
387 }
388 
389 /* assumes that cf represents the rationals, i.e. Q, and will only
390  be called in that case;
391  assumes furthermore that f != NULL and that the denominator of f != 1;
392  generally speaking, this method removes denominators in the rational
393  coefficients of the numerator and denominator of 'a';
394  more concretely, the following normalizations will be performed,
395  where t^alpha denotes a monomial in the transcendental variables t_k
396  (1) if 'a' is of the form
397  (sum_alpha a_alpha/b_alpha * t^alpha)
398  -------------------------------------
399  (sum_beta c_beta/d_beta * t^beta)
400  with integers a_alpha, b_alpha, c_beta, d_beta, then both the
401  numerator and the denominator will be multiplied by the LCM of
402  the b_alpha's and the d_beta's (if this LCM is != 1),
403  (2) if 'a' is - e.g. after having performed step (1) - of the form
404  (sum_alpha a_alpha * t^alpha)
405  -----------------------------
406  (sum_beta c_beta * t^beta)
407  with integers a_alpha, c_beta, and with a non-constant denominator,
408  then both the numerator and the denominator will be divided by the
409  GCD of the a_alpha's and the c_beta's (if this GCD is != 1),
410  this procedure does not alter COM(f) (this has to be done by the
411  calling procedure);
412  modifies f */
413 static void handleNestedFractionsOverQ(fraction f, const coeffs cf)
414 {
416  assume(!IS0(f));
417  assume(!DENIS1(f));
418 
419  { /* step (1); see documentation of this procedure above */
420  number lcmOfDenominators = n_Init(1, ntCoeffs);
421  number c; number tmp;
422  poly p = NUM(f);
423  /* careful when using n_NormalizeHelper!!! It computes the lcm of the numerator
424  of the 1st argument and the denominator of the 2nd!!! */
425  while (p != NULL)
426  {
427  c = p_GetCoeff(p, ntRing);
428  tmp = n_NormalizeHelper(lcmOfDenominators, c, ntCoeffs);
429  n_Delete(&lcmOfDenominators, ntCoeffs);
430  lcmOfDenominators = tmp;
431  pIter(p);
432  }
433  p = DEN(f);
434  while (p != NULL)
435  {
436  c = p_GetCoeff(p, ntRing);
437  tmp = n_NormalizeHelper(lcmOfDenominators, c, ntCoeffs);
438  n_Delete(&lcmOfDenominators, ntCoeffs);
439  lcmOfDenominators = tmp;
440  pIter(p);
441  }
442  if (!n_IsOne(lcmOfDenominators, ntCoeffs))
443  { /* multiply NUM(f) and DEN(f) with lcmOfDenominators */
444  NUM(f) = __p_Mult_nn(NUM(f), lcmOfDenominators, ntRing);
445  p_Normalize(NUM(f), ntRing);
446  DEN(f) = __p_Mult_nn(DEN(f), lcmOfDenominators, ntRing);
447  p_Normalize(DEN(f), ntRing);
448  }
449  n_Delete(&lcmOfDenominators, ntCoeffs);
450  if (DEN(f)!=NULL)
451  { /* step (2); see documentation of this procedure above */
452  p = NUM(f);
453  number gcdOfCoefficients = n_Copy(p_GetCoeff(p, ntRing), ntCoeffs);
454  pIter(p);
455  while ((p != NULL) && (!n_IsOne(gcdOfCoefficients, ntCoeffs)))
456  {
457  c = p_GetCoeff(p, ntRing);
458  tmp = n_Gcd(c, gcdOfCoefficients, ntCoeffs);
459  n_Delete(&gcdOfCoefficients, ntCoeffs);
460  gcdOfCoefficients = tmp;
461  pIter(p);
462  }
463  p = DEN(f);
464  while ((p != NULL) && (!n_IsOne(gcdOfCoefficients, ntCoeffs)))
465  {
466  c = p_GetCoeff(p, ntRing);
467  tmp = n_Gcd(c, gcdOfCoefficients, ntCoeffs);
468  n_Delete(&gcdOfCoefficients, ntCoeffs);
469  gcdOfCoefficients = tmp;
470  pIter(p);
471  }
472  if (!n_IsOne(gcdOfCoefficients, ntCoeffs))
473  { /* divide NUM(f) and DEN(f) by gcdOfCoefficients */
474  number inverseOfGcdOfCoefficients = n_Invers(gcdOfCoefficients,
475  ntCoeffs);
476  NUM(f) = __p_Mult_nn(NUM(f), inverseOfGcdOfCoefficients, ntRing);
477  p_Normalize(NUM(f), ntRing);
478  DEN(f) = __p_Mult_nn(DEN(f), inverseOfGcdOfCoefficients, ntRing);
479  p_Normalize(DEN(f), ntRing);
480  n_Delete(&inverseOfGcdOfCoefficients, ntCoeffs);
481  }
482  n_Delete(&gcdOfCoefficients, ntCoeffs);
483  }
484  }
485 
486  /* Now, due to the above computations, DEN(f) may have become the
487  1-polynomial which needs to be represented by NULL: */
488  if ((DEN(f) != NULL) &&
489  p_IsConstant(DEN(f), ntRing) &&
491  {
492  p_Delete(&DEN(f), ntRing); DEN(f) = NULL;
493  }
494 
495  if( DEN(f) != NULL )
496  if( !n_GreaterZero(pGetCoeff(DEN(f)), ntCoeffs) )
497  {
498  NUM(f) = p_Neg(NUM(f), ntRing);
499  DEN(f) = p_Neg(DEN(f), ntRing);
500  }
502  ntTest((number)f); // TODO!
503 }
504 
505 /// TODO: normalization of a!?
506 static number ntGetNumerator(number &a, const coeffs cf)
507 {
508  //check_N(a,cf);
509  ntTest(a);
510  if (IS0(a)) return NULL;
511 
513 
514  fraction f = (fraction)a;
515  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
516 
517  const BOOLEAN denis1= DENIS1 (f);
518 
519  if (getCoeffType (ntCoeffs) == n_Q && !denis1)
521 
522  if (getCoeffType (ntCoeffs) == n_Q && denis1)
523  {
524  assume( DEN (f) == NULL );
525 
526  number g;
527  // TODO/NOTE: the following should not be necessary (due to
528  // Hannes!) as NUM (f) should be over Z!!!
530 
531 
533 
534  if( !n_GreaterZero(g, ntCoeffs) )
535  {
536  NUM (f) = p_Neg(NUM (f), ntRing);
537  g = n_InpNeg(g, ntCoeffs);
538  }
539 
540  // g should be a positive integer now!
542 
543  if( !n_IsOne(g, ntCoeffs) )
544  {
545  DEN (f) = p_NSet(g, ntRing);
546  COM (f) ++;
547  assume( DEN (f) != NULL );
548  }
549  else
550  n_Delete(&g, ntCoeffs);
551 
552  ntTest(a);
553  }
554 
555  // Call ntNormalize instead of above?!?
556 
557  NUM (result) = p_Copy (NUM (f), ntRing); // ???
558  //DEN (result) = NULL; // done by ..Alloc0..
559  //COM (result) = 0; // done by ..Alloc0..
560 
561  ntTest((number)result);
562  //check_N((number)result,cf);
563  return (number)result;
564 }
565 
566 /// TODO: normalization of a!?
567 static number ntGetDenom(number &a, const coeffs cf)
568 {
569  //check_N(a,cf);
570  ntTest(a);
571 
572  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
573  //DEN (result)= NULL; // done by ..Alloc0..
574  //COM (result)= 0; // done by ..Alloc0..
575 
576  if (IS0(a))
577  {
578  NUM (result) = p_One(ntRing);
579  return (number)result;
580  }
581 
583 
584  fraction f = (fraction)a;
585 
586  assume( !IS0(f) );
587 
588  const BOOLEAN denis1 = DENIS1 (f);
589 
590  if( denis1 && (getCoeffType (ntCoeffs) != n_Q) ) // */1 or 0
591  {
592  NUM (result)= p_One(ntRing);
593  ntTest((number)result);
594  return (number)result;
595  }
596 
597  if (!denis1) // */* / Q
598  {
599  assume( DEN (f) != NULL );
600 
601  if (getCoeffType (ntCoeffs) == n_Q)
603 
604  ntTest(a);
605 
606  if( DEN (f) != NULL ) // is it ?? // 1 now???
607  {
608  assume( !p_IsOne(DEN (f), ntRing) );
609 
610  NUM (result) = p_Copy (DEN (f), ntRing);
611  ntTest((number)result);
612  return (number)result;
613  }
614 // NUM (result) = p_One(ntRing); // NOTE: just in order to be sure...
615  }
616 
617  // */1 / Q
619  assume( DEN (f) == NULL );
620 
621  number g;
622 // poly num= p_Copy (NUM (f), ntRing); // ???
623 
624 
625  // TODO/NOTE: the following should not be necessary (due to
626  // Hannes!) as NUM (f) should be over Z!!!
628 
629  n_ClearDenominators(itr, g, ntCoeffs); // may return -1 :(((
630 
631  if( !n_GreaterZero(g, ntCoeffs) )
632  {
633 // NUM (f) = p_Neg(NUM (f), ntRing); // Ugly :(((
634 // g = n_InpNeg(g, ntCoeffs);
635  NUM (f) = p_Neg(NUM (f), ntRing); // Ugly :(((
636  g = n_InpNeg(g, ntCoeffs);
637  }
638 
639  // g should be a positive integer now!
641 
642  if( !n_IsOne(g, ntCoeffs) )
643  {
645  assume( !n_IsOne(g, ntCoeffs) );
646 
647  DEN (f) = p_NSet(g, ntRing); // update COM(f)???
648  assume( DEN (f) != NULL );
649  COM (f) ++;
650 
651  NUM (result)= p_Copy (DEN (f), ntRing);
652  }
653  else
654  { // common denom == 1?
655  NUM (result)= p_NSet(g, ntRing); // p_Copy (DEN (f), ntRing);
656 // n_Delete(&g, ntCoeffs);
657  }
658 
659 // if (!p_IsConstant (num, ntRing) && pNext(num) != NULL)
660 // else
661 // g= p_GetAllDenom (num, ntRing);
662 // result= (fraction) ntSetMap (ntCoeffs, cf) (g, ntCoeffs, cf);
663 
664  ntTest((number)result);
665  //check_N((number)result,cf);
666  return (number)result;
667 }
668 
669 static BOOLEAN ntIsOne(number a, const coeffs cf)
670 {
671  //check_N(a,cf);
672  ntTest(a); // !!!
674  fraction f = (fraction)a;
675  return (f!=NULL) && DENIS1(f) && NUMIS1(f);
676 }
677 
678 static BOOLEAN ntIsMOne(number a, const coeffs cf)
679 {
680  //check_N(a,cf);
681  ntTest(a);
683  fraction f = (fraction)a;
684  if ((f==NULL) || (!DENIS1(f))) return FALSE;
685  poly g = NUM(f);
686  if (!p_IsConstant(g, ntRing)) return FALSE;
687  return n_IsMOne(p_GetCoeff(g, ntRing), ntCoeffs);
688 }
689 
690 /// this is in-place, modifies a
691 static number ntNeg(number a, const coeffs cf)
692 {
693  //check_N(a,cf);
694  ntTest(a);
695  if (!IS0(a))
696  {
697  fraction f = (fraction)a;
698  NUM(f) = p_Neg(NUM(f), ntRing);
699  }
700  ntTest(a);
701  return a;
702 }
703 
704 number ntInit(long i, const coeffs cf)
705 {
706  if (i != 0)
707  {
708  poly p=p_ISet(i, ntRing);
709  if (p!=NULL)
710  {
711  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
712  NUM(result) = p;
713  //DEN(result) = NULL; // done by omAlloc0Bin
714  //COM(result) = 0; // done by omAlloc0Bin
715  ntTest((number)result);
716  //check_N((number)result,cf);
717  return (number)result;
718  }
719  }
720  return NULL;
721 }
722 
723 
724 // takes over p!
725 number ntInit(poly p, const coeffs cf)
726 {
727  if (p == NULL) return NULL;
728 
729  p_Test( p, ntRing);
730  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
731 
732  if (nCoeff_is_Q(ntCoeffs))
733  {
734  number g;
735  // the following is necessary because
736  // NUM (f) should be over Z,
737  // while p may be over Q
739 
741 
742  if( !n_GreaterZero(g, ntCoeffs) )
743  {
744  p = p_Neg(p, ntRing);
745  g = n_InpNeg(g, ntCoeffs);
746  }
747 
748  // g should be a positive integer now!
750 
751  if( !n_IsOne(g, ntCoeffs) )
752  {
753  DEN (f) = p_NSet(g, ntRing);
754  p_Normalize(DEN(f), ntRing);
755  assume( DEN (f) != NULL );
756  }
757  else
758  {
759  //DEN(f) = NULL; // done by omAlloc0
760  n_Delete(&g, ntCoeffs);
761  }
762  }
763 
764  p_Normalize(p, ntRing);
765  NUM(f) = p;
766  //COM(f) = 0; // done by omAlloc0
767 
768  //check_N((number)f,cf);
769  ntTest((number)f);
770  return (number)f;
771 }
772 
773 static long ntInt(number &a, const coeffs cf)
774 {
775  //check_N(a,cf);
776  ntTest(a);
777  if (IS0(a)) return 0;
779  fraction f = (fraction)a;
780  if (!DENIS1(f)) return 0;
781 
782  const poly aAsPoly = NUM(f);
783 
784  if(aAsPoly == NULL)
785  return 0;
786 
787  if (!p_IsConstant(aAsPoly, ntRing))
788  return 0;
789 
790  assume( aAsPoly != NULL );
791 
792  return n_Int(p_GetCoeff(aAsPoly, ntRing), ntCoeffs);
793 }
794 
795 /* return FALSE, if a is a constant <=0 */
796 static BOOLEAN ntGreaterZero(number a, const coeffs cf)
797 {
798  //check_N(a,cf);
799  ntTest(a);
800  if (IS0(a)) return FALSE;
801  fraction f = (fraction)a;
802  poly g = NUM(f);
804 }
805 
806 static BOOLEAN ntGreater(number a, number b, const coeffs cf)
807 {
808  //check_N(a,cf);
809  //check_N(b,cf);
810  ntTest(a);
811  ntTest(b);
812  if (IS0(a))
813  {
814  if (IS0(b)) return FALSE;
815  fraction fb = (fraction)b;
816  return (!n_GreaterZero(pGetCoeff(NUM(fb)), ntCoeffs));
817  }
818  if (IS0(b))
819  {
820  fraction fa = (fraction)a;
821  return (n_GreaterZero(pGetCoeff(NUM(fa)), ntCoeffs));
822  }
823  // now: a!=0, b!=0
824  fraction fa = (fraction)a;
825  number aNumCoeff = p_GetCoeff(NUM(fa), ntRing);
826  int aNumDeg = p_Totaldegree(NUM(fa), ntRing);
827  number aDenCoeff = NULL; int aDenDeg = 0;
828  if (DEN(fa)!=NULL)
829  {
830  aDenCoeff=p_GetCoeff(DEN(fa),ntRing);
831  aDenDeg = p_Totaldegree(DEN(fa), ntRing);
832  }
833  fraction fb = (fraction)b;
834  number bNumCoeff = p_GetCoeff(NUM(fb), ntRing);
835  int bNumDeg = p_Totaldegree(NUM(fb), ntRing);
836  number bDenCoeff = NULL; int bDenDeg = 0;
837  if (DEN(fb)!=NULL)
838  {
839  bDenCoeff=p_GetCoeff(DEN(fb),ntRing);
840  bDenDeg = p_Totaldegree(DEN(fb), ntRing);
841  }
842  if (aNumDeg-aDenDeg > bNumDeg-bDenDeg) return TRUE;
843  if (aNumDeg-aDenDeg < bNumDeg-bDenDeg) return FALSE;
844  number aa;
845  number bb;
846  if (bDenCoeff==NULL) aa=n_Copy(aNumCoeff,ntCoeffs);
847  else aa=n_Mult(aNumCoeff,bDenCoeff,ntCoeffs);
848  if (aDenCoeff==NULL) bb=n_Copy(bNumCoeff,ntCoeffs);
849  else bb=n_Mult(bNumCoeff,aDenCoeff,ntCoeffs);
850  BOOLEAN rr= n_Greater(aa, bb, ntCoeffs);
851  n_Delete(&aa,ntCoeffs);
852  n_Delete(&bb,ntCoeffs);
853  return rr;
854 }
855 
856 static void ntCoeffWrite(const coeffs cf, BOOLEAN details)
857 {
858  assume( cf != NULL );
859 
860  const ring A = cf->extRing;
861 
862  assume( A != NULL );
863  assume( A->cf != NULL );
864 
865  n_CoeffWrite(A->cf, details);
866 
867 // rWrite(A);
868 
869  const int P = rVar(A);
870  assume( P > 0 );
871 
872  PrintS("(");
873 
874  for (int nop=0; nop < P; nop ++)
875  {
876  Print("%s", rRingVar(nop, A));
877  if (nop!=P-1) PrintS(", ");
878  }
879 
880  PrintS(")");
881 
882  assume( A->qideal == NULL );
883 
884 /*
885  PrintS("// Coefficients live in the rational function field\n");
886  Print("// K(");
887  for (int i = 0; i < rVar(ntRing); i++)
888  {
889  if (i > 0) PrintS(" ");
890  Print("%s", rRingVar(i, ntRing));
891  }
892  PrintS(") with\n");
893  PrintS("// K: "); n_CoeffWrite(cf->extRing->cf);
894 */
895 }
896 
897 number ntDiff(number a, number d, const coeffs cf)
898 {
899  //check_N(a,cf);
900  //check_N(d,cf);
901  ntTest(a);
902  ntTest(d);
903 
904  if (IS0(d))
905  {
906  WerrorS("ringvar expected");
907  return NULL;
908  }
909  fraction t = (fraction) d;
910  if (!DENIS1(t))
911  {
912  WerrorS("expected differentiation by a variable");
913  return NULL;
914  }
915  int k=p_Var(NUM(t),ntRing);
916  if (k==0)
917  {
918  WerrorS("expected differentiation by a variable");
919  return NULL;
920  }
921 
922  if (IS0(a)) return ntCopy(a, cf);
923 
924  fraction fa = (fraction)a;
925  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
926  if (DENIS1(fa))
927  {
928  NUM(result) = p_Diff(NUM(fa),k,ntRing);
929  //DEN(result) = NULL; // done by ..Alloc0..
930  if (NUM(result)==NULL)
931  {
933  return(NULL);
934  }
936  //check_N((number)result,cf);
937  ntTest((number)result);
938  return (number)result;
939  }
940 
941  poly fg = p_Mult_q(p_Copy(DEN(fa),ntRing),p_Diff(NUM(fa),k,ntRing),ntRing);
942  poly gf = p_Mult_q(p_Copy(NUM(fa),ntRing),p_Diff(DEN(fa),k,ntRing),ntRing);
943  NUM(result) = p_Sub(fg,gf,ntRing);
944  if (NUM(result)==NULL) return(NULL);
945  DEN(result) = pp_Mult_qq(DEN(fa), DEN(fa), ntRing);
946  COM(result) = COM(fa) + COM(fa) + DIFF_COMPLEXITY;
948 
949  //check_N((number)result,cf);
950  ntTest((number)result);
951  return (number)result;
952 }
953 
954 static number ntAdd(number a, number b, const coeffs cf)
955 {
956  //check_N(a,cf);
957  //check_N(b,cf);
958  ntTest(a);
959  ntTest(b);
960  if (IS0(a)) return ntCopy(b, cf);
961  if (IS0(b)) return ntCopy(a, cf);
962 
963  fraction fa = (fraction)a;
964  fraction fb = (fraction)b;
965 
966  poly g = p_Copy(NUM(fa), ntRing);
967  if (!DENIS1(fb)) g = p_Mult_q(g, p_Copy(DEN(fb), ntRing), ntRing);
968  poly h = p_Copy(NUM(fb), ntRing);
969  if (!DENIS1(fa)) h = p_Mult_q(h, p_Copy(DEN(fa), ntRing), ntRing);
970  g = p_Add_q(g, h, ntRing);
971 
972  if (g == NULL) return NULL;
973 
974  poly f;
975  if (DENIS1(fa) && DENIS1(fb)) f = NULL;
976  else if (!DENIS1(fa) && DENIS1(fb)) f = p_Copy(DEN(fa), ntRing);
977  else if (DENIS1(fa) && !DENIS1(fb)) f = p_Copy(DEN(fb), ntRing);
978  else /* both denom's are != 1 */ f = p_Mult_q(p_Copy(DEN(fa), ntRing),
979  p_Copy(DEN(fb), ntRing),
980  ntRing);
981 
982  fraction result = (fraction)omAllocBin(fractionObjectBin);
983  NUM(result) = g;
984  DEN(result) = f;
985  COM(result) = COM(fa) + COM(fb) + ADD_COMPLEXITY;
987 
988 // ntTest((number)result);
989 
990  //check_N((number)result,cf);
991  ntTest((number)result);
992  return (number)result;
993 }
994 
995 static number ntSub(number a, number b, const coeffs cf)
996 {
997  //check_N(a,cf);
998  //check_N(b,cf);
999  ntTest(a);
1000  ntTest(b);
1001  if (IS0(a)) return ntNeg(ntCopy(b, cf), cf);
1002  if (IS0(b)) return ntCopy(a, cf);
1003 
1004  fraction fa = (fraction)a;
1005  fraction fb = (fraction)b;
1006 
1007  poly g = p_Copy(NUM(fa), ntRing);
1008  if (!DENIS1(fb)) g = p_Mult_q(g, p_Copy(DEN(fb), ntRing), ntRing);
1009  poly h = p_Copy(NUM(fb), ntRing);
1010  if (!DENIS1(fa)) h = p_Mult_q(h, p_Copy(DEN(fa), ntRing), ntRing);
1011  g = p_Add_q(g, p_Neg(h, ntRing), ntRing);
1012 
1013  if (g == NULL) return NULL;
1014 
1015  poly f;
1016  if (DENIS1(fa) && DENIS1(fb)) f = NULL;
1017  else if (!DENIS1(fa) && DENIS1(fb)) f = p_Copy(DEN(fa), ntRing);
1018  else if (DENIS1(fa) && !DENIS1(fb)) f = p_Copy(DEN(fb), ntRing);
1019  else /* both den's are != 1 */ f = p_Mult_q(p_Copy(DEN(fa), ntRing),
1020  p_Copy(DEN(fb), ntRing),
1021  ntRing);
1022 
1023  fraction result = (fraction)omAllocBin(fractionObjectBin);
1024  NUM(result) = g;
1025  DEN(result) = f;
1026  COM(result) = COM(fa) + COM(fb) + ADD_COMPLEXITY;
1028 // ntTest((number)result);
1029  //check_N((number)result,cf);
1030  ntTest((number)result);
1031  return (number)result;
1032 }
1033 
1034 static number ntMult(number a, number b, const coeffs cf)
1035 {
1036  //check_N(a,cf);
1037  //check_N(b,cf);
1038  ntTest(a); // !!!?
1039  ntTest(b); // !!!?
1040 
1041  if (IS0(a) || IS0(b)) return NULL;
1042 
1043  fraction fa = (fraction)a;
1044  fraction fb = (fraction)b;
1045 
1046  const poly g = pp_Mult_qq(NUM(fa), NUM(fb), ntRing);
1047 
1048  if (g == NULL) return NULL; // may happen due to zero divisors???
1049 
1050  fraction result = (fraction)omAllocBin(fractionObjectBin);
1051 
1052  NUM(result) = g;
1053 
1054  const poly da = DEN(fa);
1055  const poly db = DEN(fb);
1056 
1057 
1058  //check_N((number)result,cf);
1059  if (db == NULL)
1060  {
1061  // b = ? // NULL
1062 
1063  if(da == NULL)
1064  { // both fa && fb are ?? // NULL!
1065  DEN(result) = NULL;
1066  COM(result) = 0;
1067  p_Normalize(g,ntRing);
1068  }
1069  else
1070  {
1071  DEN(result) = p_Copy(da, ntRing);
1074  //check_N((number)result,cf);
1075  }
1076  }
1077  else
1078  { // b = ?? / ??
1079  if (da == NULL)
1080  { // a == ? // NULL
1081  DEN(result) = p_Copy(db, ntRing);
1084  //check_N((number)result,cf);
1085  }
1086  else /* both den's are != 1 */
1087  {
1088  DEN(result) = pp_Mult_qq(da, db, ntRing);
1089  COM(result) = COM(fa) + COM(fb) + MULT_COMPLEXITY;
1091  //check_N((number)result,cf);
1092  }
1093  }
1094 
1095 // ntTest((number)result);
1096 
1097  //check_N((number)result,cf);
1098  ntTest((number)result);
1099  return (number)result;
1100 }
1101 
1102 static void ntNormalizeDen(fraction result, const ring R)
1103 {
1104  if ((nCoeff_has_simple_inverse(R->cf))
1105  && (result!=NULL)
1106  && (DEN(result)!=NULL))
1107  {
1108  poly n=DEN(result);
1109  if (!n_IsOne(pGetCoeff(n),R->cf))
1110  {
1111  number inv=n_Invers(pGetCoeff(n),R->cf);
1112  DEN(result)=__p_Mult_nn(n,inv,R);
1113  NUM(result)=__p_Mult_nn(NUM(result),inv,R);
1114  n_Delete(&inv,R->cf);
1115  if (p_IsOne(DEN(result), R))
1116  {
1117  n=DEN(result);
1118  DEN(result)=NULL;
1119  COM(result) = 0;
1120  p_Delete(&n,R);
1121  }
1122  }
1123  }
1124 }
1125 
1126 static number ntDiv(number a, number b, const coeffs cf)
1127 {
1128  //check_N(a,cf);
1129  //check_N(b,cf);
1130  ntTest(a);
1131  ntTest(b);
1132  if (IS0(a)) return NULL;
1133  if (IS0(b)) WerrorS(nDivBy0);
1134 
1135  fraction fa = (fraction)a;
1136  fraction fb = (fraction)b;
1137 
1138  poly g = p_Copy(NUM(fa), ntRing);
1139  if (!DENIS1(fb)) g = p_Mult_q(g, p_Copy(DEN(fb), ntRing), ntRing);
1140 
1141  if (g == NULL) return NULL; /* may happen due to zero divisors */
1142 
1143  poly f = p_Copy(NUM(fb), ntRing);
1144  if (!DENIS1(fa)) f = p_Mult_q(f, p_Copy(DEN(fa), ntRing), ntRing);
1145 
1146  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1147  NUM(result) = g;
1149  {
1150  g=p_Neg(g,ntRing);
1151  f=p_Neg(f,ntRing);
1152  NUM(result) = g;
1153  }
1155  {
1156  DEN(result) = f;
1157  }
1158  else
1159  {
1160  p_Delete(&f, ntRing);
1161  }
1162  COM(result) = COM(fa) + COM(fb) + MULT_COMPLEXITY;
1163 // definiteGcdCancellation((number)result, cf,FALSE);
1165 // ntTest((number)result);
1166  //check_N((number)result,cf);
1168  ntTest((number)result);
1169  return (number)result;
1170 }
1171 
1172 static number ntInvers(number a, const coeffs cf)
1173 {
1174  //check_N(a,cf);
1175  ntTest(a);
1176  if (IS0(a))
1177  {
1178  WerrorS(nDivBy0);
1179  return NULL;
1180  }
1181  fraction f = (fraction)a;
1182  assume( f != NULL );
1183 
1184  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1185 
1186  assume( NUM(f) != NULL );
1187  const poly den = DEN(f);
1188 
1189  if (den == NULL)
1190  NUM(result) = p_One(ntRing);
1191  else
1192  NUM(result) = p_Copy(den, ntRing);
1193 
1194  if( !NUMIS1(f) )
1195  {
1196  poly num_f=NUM(f);
1197  BOOLEAN neg= !n_GreaterZero(pGetCoeff(num_f),ntCoeffs);
1198  if (neg)
1199  {
1200  num_f=p_Neg(p_Copy(num_f, ntRing), ntRing);
1202  }
1203  else
1204  {
1205  num_f=p_Copy(num_f, ntRing);
1206  }
1207  DEN(result) = num_f;
1208  COM(result) = COM(f);
1209  if (neg)
1210  {
1211  if (p_IsOne(num_f, ntRing))
1212  {
1213  DEN(result)=NULL;
1214  //COM(result) = 0;
1215  p_Delete(&num_f,ntRing);
1216  }
1217  }
1218  }
1219  //else// Alloc0
1220  //{
1221  // DEN(result) = NULL;
1222  // COM(result) = 0;
1223  //}
1225  ntTest((number)result); // !!!!
1226  //check_N((number)result,cf);
1227  return (number)result;
1228 }
1229 
1230 /* 0^0 = 0;
1231  for |exp| <= 7 compute power by a simple multiplication loop;
1232  for |exp| >= 8 compute power along binary presentation of |exp|, e.g.
1233  p^13 = p^1 * p^4 * p^8, where we utilise that
1234  p^(2^(k+1)) = p^(2^k) * p^(2^k);
1235  intermediate cancellation is controlled by the in-place method
1236  heuristicGcdCancellation; see there.
1237 */
1238 static void ntPower(number a, int exp, number *b, const coeffs cf)
1239 {
1240  ntTest(a);
1241 
1242  /* special cases first */
1243  if (IS0(a))
1244  {
1245  if (exp >= 0) *b = NULL;
1246  else WerrorS(nDivBy0);
1247  }
1248  else if (exp == 0) { *b = ntInit(1, cf); return;}
1249  else if (exp == 1) { *b = ntCopy(a, cf); return;}
1250  else if (exp == -1) { *b = ntInvers(a, cf); return;}
1251 
1252  int expAbs = exp; if (expAbs < 0) expAbs = -expAbs;
1253 
1254  /* now compute a^expAbs */
1255  number pow; number t;
1256  if (expAbs <= 7)
1257  {
1258  pow = ntCopy(a, cf);
1259  for (int i = 2; i <= expAbs; i++)
1260  {
1261  t = ntMult(pow, a, cf);
1262  ntDelete(&pow, cf);
1263  pow = t;
1265  }
1266  }
1267  else
1268  {
1269  pow = ntInit(1, cf);
1270  number factor = ntCopy(a, cf);
1271  while (expAbs != 0)
1272  {
1273  if (expAbs & 1)
1274  {
1275  t = ntMult(pow, factor, cf);
1276  ntDelete(&pow, cf);
1277  pow = t;
1279  }
1280  expAbs = expAbs / 2;
1281  if (expAbs != 0)
1282  {
1283  t = ntMult(factor, factor, cf);
1284  ntDelete(&factor, cf);
1285  factor = t;
1287  }
1288  }
1289  ntDelete(&factor, cf);
1290  }
1291 
1292  /* invert if original exponent was negative */
1293  if (exp < 0)
1294  {
1295  t = ntInvers(pow, cf);
1296  ntDelete(&pow, cf);
1297  pow = t;
1298  }
1299  *b = pow;
1300  ntTest(*b);
1301  //check_N(*b,cf);
1302 }
1303 
1304 /* modifies a */
1305 /* this is an intermediate simplification routine - not a comple "normalize" */
1306 static void heuristicGcdCancellation(number a, const coeffs cf)
1307 {
1308  if (IS0(a)) return;
1309 
1310  fraction f = (fraction)a;
1311  p_Normalize(NUM(f),ntRing);
1312  if (DENIS1(f) || NUMIS1(f)) { COM(f) = 0; return; }
1313 
1314  assume( DEN(f) != NULL );
1315  p_Normalize(DEN(f),ntRing);
1316 
1317  /* check whether NUM(f) = DEN(f), and - if so - replace 'a' by 1 */
1318  if (p_EqualPolys(NUM(f), DEN(f), ntRing))
1319  { /* numerator and denominator are both != 1 */
1320  p_Delete(&NUM(f), ntRing); NUM(f) = p_ISet(1, ntRing);
1321  p_Delete(&DEN(f), ntRing); DEN(f) = NULL;
1322  COM(f) = 0;
1323  }
1324  else
1325  {
1326  if (COM(f) > BOUND_COMPLEXITY)
1328 
1329  // TODO: check if it is enough to put the following into definiteGcdCancellation?!
1330  if( DEN(f) != NULL )
1331  {
1332  if( !n_GreaterZero(pGetCoeff(DEN(f)), ntCoeffs) )
1333  {
1334  NUM(f) = p_Neg(NUM(f), ntRing);
1335  DEN(f) = p_Neg(DEN(f), ntRing);
1336  }
1337  if (ntCoeffs->has_simple_Inverse)
1338  {
1339  if (!n_IsOne(pGetCoeff(DEN(f)),ntCoeffs))
1340  {
1341  number inv=n_Invers(pGetCoeff(DEN(f)),ntCoeffs);
1342  DEN(f)=__p_Mult_nn(DEN(f),inv,ntRing);
1343  NUM(f)=__p_Mult_nn(NUM(f),inv,ntRing);
1344  }
1345  if(p_LmIsConstant(DEN(f),ntRing))
1346  {
1347  p_Delete(&DEN(f),ntRing);
1348  COM(f)=0;
1349  }
1350  }
1351  if ((DEN(f)!=NULL)
1352  && (pNext(DEN(f))==NULL))
1353  {
1354  poly den_f=DEN(f);
1355  poly h=NUM(f);
1356  loop
1357  {
1358  if (h==NULL)
1359  {
1360  h=NUM(f);
1361  do
1362  {
1363  p_ExpVectorDiff(h,h,den_f,ntRing);
1364  pIter(h);
1365  } while(h!=NULL);
1366  p_ExpVectorDiff(den_f,den_f,den_f,ntRing);
1367  break;
1368  }
1369  int i=0;
1370  do
1371  {
1372  i++;
1373  if (p_GetExp(den_f,i,ntRing) > p_GetExp(h,i,ntRing)) return;
1374  } while(i<ntRing->N);
1375  pIter(h);
1376  }
1377  }
1378  }
1379  }
1380  if ((DEN(f)!=NULL)
1381  && (pNext(DEN(f))==NULL)
1382  && (p_LmIsConstantComp(DEN(f),ntRing))
1383  && (n_IsOne(pGetCoeff(DEN(f)),ntCoeffs)))
1384  {
1385  p_Delete(&DEN(f),ntRing);
1386  COM(f)=0;
1387  }
1388 }
1389 
1390 /// modifies a
1391 static void definiteGcdCancellation(number a, const coeffs cf,
1392  BOOLEAN simpleTestsHaveAlreadyBeenPerformed)
1393 {
1394 // ntTest(a); // !!!!
1395 
1396  fraction f = (fraction)a;
1397 
1398  if (IS0(a)) return;
1399  if (COM(f)==0) return;
1400  if (DENIS1(f) || NUMIS1(f)) { COM(f) = 0; ntTest(a); return; }
1401  if (!simpleTestsHaveAlreadyBeenPerformed)
1402  {
1403 
1404  /* check whether NUM(f) = DEN(f), and - if so - replace 'a' by 1 */
1405  if (p_EqualPolys(NUM(f), DEN(f), ntRing))
1406  { /* numerator and denominator are both != 1 */
1407  p_Delete(&NUM(f), ntRing); NUM(f) = p_ISet(1, ntRing);
1408  p_Delete(&DEN(f), ntRing); DEN(f) = NULL;
1409  COM(f) = 0;
1410  ntTest(a);
1411  return;
1412  }
1413  }
1414  /*if (rField_is_Q(ntRing))
1415  {
1416  number c=n_Copy(pGetCoeff(NUM(f)),ntCoeffs);
1417  poly p=pNext(NUM(f));
1418  while((p!=NULL)&&(!n_IsOne(c,ntCoeffs)))
1419  {
1420  number cc=n_Gcd(c,pGetCoeff(p),ntCoeffs);
1421  n_Delete(&c,ntCoeffs);
1422  c=cc;
1423  pIter(p);
1424  };
1425  p=DEN(f);
1426  while((p!=NULL)&&(!n_IsOne(c,ntCoeffs)))
1427  {
1428  number cc=n_Gcd(c,pGetCoeff(p),ntCoeffs);
1429  n_Delete(&c,ntCoeffs);
1430  c=cc;
1431  pIter(p);
1432  };
1433  if(!n_IsOne(c,ntCoeffs))
1434  {
1435  p=NUM(f);
1436  do
1437  {
1438  number cc=n_Div(pGetCoeff(p),c,ntCoeffs);
1439  n_Normalize(cc,ntCoeffs);
1440  p_SetCoeff(p,cc,ntRing);
1441  pIter(p);
1442  } while(p!=NULL);
1443  p=DEN(f);
1444  do
1445  {
1446  number cc=n_Div(pGetCoeff(p),c,ntCoeffs);
1447  n_Normalize(cc,ntCoeffs);
1448  p_SetCoeff(p,cc,ntRing);
1449  pIter(p);
1450  } while(p!=NULL);
1451  n_Delete(&c,ntCoeffs);
1452  if(pNext(DEN(f))==NULL)
1453  {
1454  if (p_IsOne(DEN(f),ntRing))
1455  {
1456  p_LmDelete(&DEN(f),ntRing);
1457  COM(f)=0;
1458  return;
1459  }
1460  else
1461  {
1462  return;
1463  }
1464  }
1465  }
1466  }*/
1467 
1468  /* here we assume: NUM(f), DEN(f) !=NULL, in Z_a reqp. Z/p_a */
1469  //StringSetS("");ntWriteLong(a,cf);
1470  poly pGcd = singclap_gcd_and_divide(NUM(f), DEN(f), ntRing);
1471  //PrintS("gcd= ");p_wrp(pGcd,ntRing);PrintLn();
1472  if (p_IsConstant(pGcd, ntRing)
1473  && n_IsOne(p_GetCoeff(pGcd, ntRing), ntCoeffs)
1474  )
1475  { /* gcd = 1; nothing to cancel;
1476  Suppose the given rational function field is over Q. Although the
1477  gcd is 1, we may have produced fractional coefficients in NUM(f),
1478  DEN(f), or both, due to previous arithmetics. The next call will
1479  remove those nested fractions, in case there are any. */
1480  if (nCoeff_is_Zp(ntCoeffs))
1481  {
1482  number d=p_GetCoeff (DEN(f),ntRing);
1483  BOOLEAN d_not_1=FALSE;
1484  if (!n_IsOne(d,ntCoeffs))
1485  {
1486  NUM (f) = p_Div_nn (NUM (f), d, ntRing);
1487  d_not_1=TRUE;
1488  }
1489  if (p_IsConstant (DEN (f), ntRing))
1490  {
1491  p_Delete(&DEN (f), ntRing);
1492  DEN (f) = NULL;
1493  }
1494  else if (d_not_1)
1495  {
1496  DEN (f) = p_Div_nn (DEN (f), d, ntRing);
1497  }
1499  }
1500  else
1501  { /* We divide both NUM(f) and DEN(f) by the gcd which is known
1502  to be != 1. */
1503  if (p_IsConstant(DEN(f), ntRing) &&
1504  n_IsOne(p_GetCoeff(DEN(f), ntRing), ntCoeffs))
1505  {
1506  /* DEN(f) = 1 needs to be represented by NULL! */
1507  p_Delete(&DEN(f), ntRing);
1508  DEN(f) = NULL;
1509  }
1510  else
1511  {
1512  if (nCoeff_is_Zp(ntCoeffs))
1513  {
1514  NUM (f) = p_Div_nn (NUM (f), p_GetCoeff (DEN(f),ntRing), ntRing);
1515  if (p_IsConstant (DEN (f), ntRing))
1516  {
1517  p_Delete(&DEN (f), ntRing);
1518  DEN (f) = NULL;
1519  }
1520  else
1521  {
1522  p_Norm (DEN (f),ntRing);
1523  }
1524  }
1525  }
1526  }
1527  p_Delete(&pGcd, ntRing);
1528 // StringAppendS(" -> ");ntWriteLong(a,cf);StringAppendS("\n");{ char* s = StringEndS(); Print("%s", s); omFree(s); }
1529  COM(f) = 0;
1530 
1531  if( DEN(f) != NULL )
1532  {
1533  if( !n_GreaterZero(pGetCoeff(DEN(f)), ntCoeffs) )
1534  {
1535  NUM(f) = p_Neg(NUM(f), ntRing);
1536  DEN(f) = p_Neg(DEN(f), ntRing);
1537  if (p_IsConstant(DEN(f), ntRing) &&
1538  n_IsOne(p_GetCoeff(DEN(f), ntRing), ntCoeffs))
1539  {
1540  /* DEN(f) = 1 needs to be represented by NULL! */
1541  p_Delete(&DEN(f), ntRing);
1542  DEN (f) = NULL;
1543  }
1544  }
1545  }
1546  ntTest(a); // !!!!
1547 }
1548 
1549 static void ntWriteLong(number a, const coeffs cf)
1550 {
1551  ntTest(a);
1552  if (IS0(a))
1553  StringAppendS("0");
1554  else
1555  {
1556  fraction f = (fraction)a;
1557  // stole logic from napWrite from kernel/longtrans.cc of legacy singular
1558  BOOLEAN omitBrackets = p_IsConstant(NUM(f), ntRing);
1559  if (!omitBrackets) StringAppendS("(");
1561  if (!omitBrackets) StringAppendS(")");
1562  if (!DENIS1(f))
1563  {
1564  StringAppendS("/");
1565  omitBrackets = p_IsConstant(DEN(f), ntRing);
1566  if (!omitBrackets) StringAppendS("(");
1567  p_String0Long(DEN(f), ntRing, ntRing);
1568  if (!omitBrackets) StringAppendS(")");
1569  }
1570  }
1571  ntTest(a); // !!!!
1572 }
1573 
1574 static void ntWriteShort(number a, const coeffs cf)
1575 {
1576  ntTest(a);
1577  if (IS0(a))
1578  StringAppendS("0");
1579  else
1580  {
1581  fraction f = (fraction)a;
1582  // stole logic from napWrite from kernel/longtrans.cc of legacy singular
1583  BOOLEAN omitBrackets = p_IsConstant(NUM(f), ntRing);
1584  if (!omitBrackets) StringAppendS("(");
1586  if (!omitBrackets) StringAppendS(")");
1587  if (!DENIS1(f))
1588  {
1589  StringAppendS("/");
1590  omitBrackets = p_IsConstant(DEN(f), ntRing);
1591  if (!omitBrackets) StringAppendS("(");
1592  p_String0Short(DEN(f), ntRing, ntRing);
1593  if (!omitBrackets) StringAppendS(")");
1594  }
1595  }
1596  ntTest(a);
1597 }
1598 
1599 static const char * ntRead(const char *s, number *a, const coeffs cf)
1600 {
1601  poly p;
1602  const char * result = p_Read(s, p, ntRing);
1603  if (p == NULL) *a = NULL;
1604  else *a = ntInit(p, cf);
1605  ntTest(*a);
1606  return result;
1607 }
1608 
1609 static void ntNormalize (number &a, const coeffs cf)
1610 {
1611  if ( /*(*/ a!=NULL /*)*/ )
1612  {
1613  //PrintS("num=");p_wrp(NUM(a),ntRing);
1614  //PrintS(" den=");p_wrp(DEN(a),ntRing);PrintLn();
1615  if (COM((fraction)a)>0) definiteGcdCancellation(a, cf, FALSE);
1616  if ((DEN((fraction)a)!=NULL)
1617  &&(!n_GreaterZero(pGetCoeff(DEN((fraction)a)),ntCoeffs)))
1618  {
1619  NUM((fraction)a)=p_Neg(NUM((fraction)a),ntRing);
1620  DEN((fraction)a)=p_Neg(DEN((fraction)a),ntRing);
1621  }
1622  }
1623  ntNormalizeDen((fraction)a,ntRing);
1624  ntTest(a); // !!!!
1625 }
1626 
1627 static number ntExactDiv(number a, number b, const coeffs cf)
1628 {
1629  number r=ntDiv(a,b,cf);
1630  ntNormalize(r,cf);
1631  return r;
1632 }
1633 
1634 /* expects *param to be castable to TransExtInfo */
1635 static BOOLEAN ntCoeffIsEqual(const coeffs cf, n_coeffType n, void * param)
1636 {
1637  if (n_transExt != n) return FALSE;
1638  TransExtInfo *e = (TransExtInfo *)param;
1639  /* for rational function fields we expect the underlying
1640  polynomial rings to be IDENTICAL, i.e. the SAME OBJECT;
1641  this expectation is based on the assumption that we have properly
1642  registered cf and perform reference counting rather than creating
1643  multiple copies of the same coefficient field/domain/ring */
1644  if (ntRing == e->r)
1645  return TRUE;
1646 
1647  // NOTE: Q(a)[x] && Q(a)[y] should better share the _same_ Q(a)...
1648  if( rEqual(ntRing, e->r, TRUE) )
1649  {
1650  rDelete(e->r);
1651  return TRUE;
1652  }
1653 
1654  return FALSE;
1655 }
1656 
1657 static number ntNormalizeHelper(number a, number b, const coeffs cf)
1658 {
1659  ntTest(a);
1660  ntTest(b);
1661  fraction fb = (fraction)b;
1662  if ((b==NULL)||(DEN(fb)==NULL)) return ntCopy(a,cf);
1663  fraction fa = (fraction)a;
1664 
1665  poly pGcd;
1666  if (nCoeff_is_Q(ntCoeffs))
1667  {
1668  poly pa = NUM(fa);
1669  poly pb = DEN(fb);
1671  {
1672  pGcd = p_Copy(pa,ntRing);
1673  p_SetCoeff (pGcd, n_Gcd (pGetCoeff(pGcd), pGetCoeff(pb), ntCoeffs), ntRing);
1674  }
1675  else
1676  {
1677  number contentpa, contentpb, tmp;
1678 
1679  contentpb= n_Copy(p_GetCoeff(pb, ntRing),ntCoeffs);
1680  pIter(pb);
1681  while (pb != NULL)
1682  {
1683  tmp = n_SubringGcd(contentpb, p_GetCoeff(pb, ntRing) , ntCoeffs);
1684  n_Delete(&contentpb, ntCoeffs);
1685  contentpb = tmp;
1686  pIter(pb);
1687  }
1688 
1689  contentpa= n_Copy(p_GetCoeff(pa, ntRing),ntCoeffs);
1690  pIter(pa);
1691  while (pa != NULL)
1692  {
1693  tmp = n_SubringGcd(contentpa, p_GetCoeff(pa, ntRing), ntCoeffs);
1694  n_Delete(&contentpa, ntCoeffs);
1695  contentpa = tmp;
1696  pIter(pa);
1697  }
1698 
1699  tmp= n_SubringGcd (contentpb, contentpa, ntCoeffs);
1700  n_Delete(&contentpa, ntCoeffs);
1701  n_Delete(&contentpb, ntCoeffs);
1702  contentpa= tmp;
1703 
1704  pGcd = gcd_over_Q(NUM(fa), DEN(fb), ntRing);
1705  pGcd= __p_Mult_nn (pGcd, contentpa, ntRing);
1706  n_Delete(&contentpa, ntCoeffs);
1707  }
1708  }
1709  else
1710  pGcd = singclap_gcd_r(NUM(fa), DEN(fb), ntRing);
1711 
1712  /* Note that, over Q, singclap_gcd will remove the denominators in all
1713  rational coefficients of pa and pb, before starting to compute
1714  the gcd. Thus, we do not need to ensure that the coefficients of
1715  pa and pb live in Z; they may well be elements of Q\Z. */
1716 
1717  if (p_IsConstant(pGcd, ntRing) &&
1718  n_IsOne(p_GetCoeff(pGcd, ntRing), ntCoeffs))
1719  { /* gcd = 1; return pa*pb*/
1720  p_Delete(&pGcd,ntRing);
1721  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1722  NUM(result) = pp_Mult_qq(NUM(fa),DEN(fb),ntRing);
1723 
1724  ntTest((number)result); // !!!!
1725 
1726  return (number)result;
1727  }
1728 
1729  /* return pa*pb/gcd */
1730  poly newNum = singclap_pdivide(NUM(fa), pGcd, ntRing);
1731  p_Delete(&pGcd,ntRing);
1732  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1733  NUM(result) = p_Mult_q(p_Copy(DEN(fb),ntRing),newNum,ntRing);
1734  ntTest((number)result); // !!!!
1735  return (number)result;
1736 }
1737 
1738 static number ntGcd(number a, number b, const coeffs cf)
1739 {
1740  ntTest(a);
1741  ntTest(b);
1742  if (a==NULL) return ntCopy(b,cf);
1743  if (b==NULL) return ntCopy(a,cf);
1744  fraction fa = (fraction)a;
1745  fraction fb = (fraction)b;
1746 
1747 
1748  poly pGcd;
1749  if (nCoeff_is_Q(ntCoeffs))
1750  {
1751  poly pa = NUM(fa);
1752  poly pb = NUM(fb);
1754  {
1755  pGcd = p_Copy(pa,ntRing);
1757  }
1758  else
1759  {
1760  number contentpa, contentpb, tmp;
1761 
1762  contentpb= n_Copy(p_GetCoeff(pb, ntRing),ntCoeffs);
1763  pIter(pb);
1764  while (pb != NULL)
1765  {
1766  tmp = n_SubringGcd(contentpb, p_GetCoeff(pb, ntRing) , ntCoeffs);
1767  n_Delete(&contentpb, ntCoeffs);
1768  contentpb = tmp;
1769  pIter(pb);
1770  }
1771 
1772  contentpa= n_Copy(p_GetCoeff(pa, ntRing),ntCoeffs);
1773  pIter(pa);
1774  while (pa != NULL)
1775  {
1776  tmp = n_SubringGcd(contentpa, p_GetCoeff(pa, ntRing), ntCoeffs);
1777  n_Delete(&contentpa, ntCoeffs);
1778  contentpa = tmp;
1779  pIter(pa);
1780  }
1781 
1782  tmp= n_SubringGcd (contentpb, contentpa, ntCoeffs);
1783  n_Delete(&contentpa, ntCoeffs);
1784  n_Delete(&contentpb, ntCoeffs);
1785  contentpa= tmp;
1786 
1787  pGcd = gcd_over_Q(NUM(fa), NUM(fb), ntRing);
1788  pGcd= __p_Mult_nn (pGcd, contentpa, ntRing);
1789  n_Delete(&contentpa, ntCoeffs);
1790  }
1791  }
1792  else
1793  pGcd = singclap_gcd_r(NUM(fa), NUM(fb), ntRing);
1794  /* Note that, over Q, singclap_gcd will remove the denominators in all
1795  rational coefficients of pa and pb, before starting to compute
1796  the gcd. Thus, we do not need to ensure that the coefficients of
1797  pa and pb live in Z; they may well be elements of Q\Z. */
1798 
1799  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1800  NUM(result) = pGcd;
1801  ntTest((number)result); // !!!!
1802  return (number)result;
1803 }
1804 //number ntGcd_dummy(number a, number b, const coeffs cf)
1805 //{
1806 // extern char my_yylinebuf[80];
1807 // Print("ntGcd in >>%s<<\n",my_yylinebuf);
1808 // return ntGcd(a,b,cf);
1809 //}
1810 
1811 static int ntSize(number a, const coeffs cf)
1812 {
1813  ntTest(a);
1814  if (IS0(a)) return 0;
1815  fraction f = (fraction)a;
1816  poly p = NUM(f);
1817  unsigned long noOfTerms = 0;
1818  unsigned long numDegree = 0;
1819  if (p!=NULL)
1820  {
1821  numDegree = p_Totaldegree(p,ntRing);
1822  noOfTerms = pLength(p);
1823  }
1824  unsigned long denDegree = 0;
1825  if (!DENIS1(f))
1826  {
1827  denDegree = p_Totaldegree(DEN(f),ntRing);
1828  noOfTerms += pLength(DEN(f));
1829  }
1830  ntTest(a); // !!!!
1831  // avoid int overflow:
1832  unsigned long t= ((numDegree + denDegree)*(numDegree + denDegree) + 1) * noOfTerms; // must be >0
1833  if (t>INT_MAX) return INT_MAX;
1834  else return (int)t;
1835 }
1836 
1837 /* assumes that src = Q or Z, dst = Q(t_1, ..., t_s) */
1838 static number ntMap00(number a, const coeffs src, const coeffs dst)
1839 {
1840  n_Test(a, src);
1841 
1842  if (n_IsZero(a, src)) return NULL;
1843  assume(src->rep == dst->extRing->cf->rep);
1844  if ((SR_HDL(a) & SR_INT) || (a->s==3))
1845  {
1846  number res=ntInit(p_NSet(n_Copy(a, src), dst->extRing), dst);
1847  n_Test(res, dst);
1848  return res;
1849  }
1850  number nn=n_GetDenom(a,src);
1851  number zz=n_GetNumerator(a,src);
1852  number res=ntInit(p_NSet(zz,dst->extRing), dst);
1853  fraction ff=(fraction)res;
1854  if (n_IsOne(nn,src)) DEN(ff)=NULL;
1855  else DEN(ff)=p_NSet(nn,dst->extRing);
1856 
1857  n_Test((number)ff,dst);
1858  //check_N((number)ff,dst);
1859  return (number)ff;
1860 }
1861 
1862 static number ntMapZ0(number a, const coeffs src, const coeffs dst)
1863 {
1864  n_Test(a, src);
1865  if (n_IsZero(a, src)) return NULL;
1866  nMapFunc nMap=n_SetMap(src,dst->extRing->cf);
1867  poly p=p_NSet(nMap(a, src,dst->extRing->cf), dst->extRing);
1868  if (n_IsZero(pGetCoeff(p),dst->extRing->cf))
1869  p_Delete(&p,dst->extRing);
1870  number res=ntInit(p, dst);
1871  n_Test(res,dst);
1872  return res;
1873 }
1874 
1875 /* assumes that src = Z/p, dst = Q(t_1, ..., t_s) */
1876 static number ntMapP0(number a, const coeffs src, const coeffs dst)
1877 {
1878  n_Test(a, src);
1879  if (n_IsZero(a, src)) return NULL;
1880  /* mapping via intermediate int: */
1881  int n = n_Int(a, src);
1882  number q = n_Init(n, dst->extRing->cf);
1883  if (n_IsZero(q, dst->extRing->cf))
1884  {
1885  n_Delete(&q, dst->extRing->cf);
1886  return NULL;
1887  }
1888  return ntInit(p_NSet(q, dst->extRing), dst);
1889 }
1890 
1891  /* assumes that either src = K(t_1, ..., t_s), dst = K(t_1, ..., t_s) */
1892 static number ntCopyMap(number a, const coeffs cf, const coeffs dst)
1893 {
1894  ntTest(a);
1895  if (IS0(a)) return NULL;
1896 
1897  const ring rSrc = cf->extRing;
1898  const ring rDst = dst->extRing;
1899 
1900  if( rSrc == rDst )
1901  return ntCopy(a, dst); // USUALLY WRONG!
1902 
1903  fraction f = (fraction)a;
1904  poly g = prCopyR(NUM(f), rSrc, rDst);
1905 
1906  poly h = NULL;
1907 
1908  if (!DENIS1(f))
1909  h = prCopyR(DEN(f), rSrc, rDst);
1910 
1911  fraction result = (fraction)omAllocBin(fractionObjectBin);
1912 
1913  NUM(result) = g;
1914  DEN(result) = h;
1915  COM(result) = COM(f);
1916  //check_N((number)result,dst);
1917  n_Test((number)result, dst);
1918  return (number)result;
1919 }
1920 
1921 static number ntGenMap(number a, const coeffs cf, const coeffs dst)
1922 {
1923  ntTest(a);
1924  if (IS0(a)) return NULL;
1925 
1926  const ring rSrc = cf->extRing;
1927  const ring rDst = dst->extRing;
1928 
1929  const nMapFunc nMap=n_SetMap(rSrc->cf,rDst->cf);
1930  fraction f = (fraction)a;
1931  poly g = prMapR(NUM(f), nMap, rSrc, rDst);
1932  /* g may contain summands with coeff 0 */
1933  poly hh=g;
1934  poly prev=NULL;
1935  while(hh!=NULL)
1936  {
1937  if (n_IsZero(pGetCoeff(hh),rDst->cf))
1938  {
1939  if (prev==NULL)
1940  {
1941  g=p_LmFreeAndNext(g,rDst);
1942  hh=g;
1943  }
1944  else
1945  {
1946  prev->next=p_LmFreeAndNext(prev->next,rDst);
1947  hh=prev->next;
1948  }
1949  }
1950  else
1951  {
1952  prev=hh;
1953  pIter(hh);
1954  }
1955  }
1956  if (g==NULL) return NULL;
1957 
1958  poly h = NULL;
1959 
1960  if (!DENIS1(f))
1961  {
1962  h = prMapR(DEN(f), nMap, rSrc, rDst);
1963  /* h may contain summands with coeff 0 */
1964  hh=h;
1965  prev=NULL;
1966  while(hh!=NULL)
1967  {
1968  if (n_IsZero(pGetCoeff(hh),rDst->cf))
1969  {
1970  if (prev==NULL)
1971  {
1972  h=p_LmFreeAndNext(h,rDst);
1973  hh=h;
1974  }
1975  else
1976  {
1977  prev->next=p_LmFreeAndNext(prev->next,rDst);
1978  hh=prev->next;
1979  }
1980  }
1981  else
1982  {
1983  prev=hh;
1984  pIter(hh);
1985  }
1986  }
1987  if (h==NULL) WerrorS("mapping to */0");
1988  }
1989 
1990  fraction result = (fraction)omAllocBin(fractionObjectBin);
1991 
1992  NUM(result) = g;
1993  DEN(result) = h;
1994  COM(result) = COM(f);
1995  //check_N((number)result,dst);
1996  n_Test((number)result, dst);
1997  return (number)result;
1998 }
1999 
2000 static number ntCopyAlg(number a, const coeffs cf, const coeffs dst)
2001 {
2002  n_Test(a, cf) ;
2003  if (n_IsZero(a, cf)) return NULL;
2004  return ntInit(prCopyR((poly)a, cf->extRing, dst->extRing),dst);
2005 }
2006 
2007 static number ntGenAlg(number a, const coeffs cf, const coeffs dst)
2008 {
2009  n_Test(a, cf) ;
2010  if (n_IsZero(a, cf)) return NULL;
2011 
2012  const nMapFunc nMap=n_SetMap(cf->extRing->cf,dst->extRing->cf);
2013  return ntInit(prMapR((poly)a, nMap, cf->extRing, dst->extRing),dst);
2014 }
2015 
2016 /* assumes that src = Q, dst = Z/p(t_1, ..., t_s) */
2017 static number ntMap0P(number a, const coeffs src, const coeffs dst)
2018 {
2019  n_Test(a, src) ;
2020  if (n_IsZero(a, src)) return NULL;
2021  // int p = rChar(dst->extRing);
2022 
2023  number q = nlModP(a, src, dst->extRing->cf); // FIXME? TODO? // extern number nlModP(number q, const coeffs Q, const coeffs Zp); // Map q \in QQ \to Zp
2024 
2025  if (n_IsZero(q, dst->extRing->cf))
2026  {
2027  n_Delete(&q, dst->extRing->cf);
2028  return NULL;
2029  }
2030 
2031  poly g = p_NSet(q, dst->extRing);
2032 
2033  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
2034  NUM(f) = g; // DEN(f) = NULL; COM(f) = 0;
2035  n_Test((number)f, dst);
2036  //check_N((number)f,dst);
2037  return (number)f;
2038 }
2039 
2040 /* assumes that src = Z/p, dst = Z/p(t_1, ..., t_s) */
2041 static number ntMapPP(number a, const coeffs src, const coeffs dst)
2042 {
2043  n_Test(a, src) ;
2044  if (n_IsZero(a, src)) return NULL;
2045  assume(src == dst->extRing->cf);
2046  poly p = p_One(dst->extRing);
2047  p_SetCoeff(p, n_Copy(a, src), dst->extRing);
2048  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
2049  NUM(f) = p; // DEN(f) = NULL; COM(f) = 0;
2050  n_Test((number)f, dst);
2051  //check_N((number)f,dst);
2052  return (number)f;
2053 }
2054 
2055 /* assumes that src = Z/u, dst = Z/p(t_1, ..., t_s), where u != p */
2056 static number ntMapUP(number a, const coeffs src, const coeffs dst)
2057 {
2058  n_Test(a, src) ;
2059  if (n_IsZero(a, src)) return NULL;
2060  /* mapping via intermediate int: */
2061  int n = n_Int(a, src);
2062  number q = n_Init(n, dst->extRing->cf);
2063  poly p;
2064  if (n_IsZero(q, dst->extRing->cf))
2065  {
2066  n_Delete(&q, dst->extRing->cf);
2067  return NULL;
2068  }
2069  p = p_One(dst->extRing);
2070  p_SetCoeff(p, q, dst->extRing);
2071  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
2072  NUM(f) = p; // DEN(f) = NULL; COM(f) = 0;
2073  n_Test((number)f, dst);
2074  //check_N((number)f,dst);
2075  return (number)f;
2076 }
2077 
2078 nMapFunc ntSetMap(const coeffs src, const coeffs dst)
2079 {
2080  /* dst is expected to be a rational function field */
2081  assume(getCoeffType(dst) == n_transExt);
2082 
2083  int h = 0; /* the height of the extension tower given by dst */
2084  coeffs bDst = nCoeff_bottom(dst, h); /* the bottom field in the tower dst */
2085  coeffs bSrc = nCoeff_bottom(src, h); /* the bottom field in the tower src */
2086 
2087  /* for the time being, we only provide maps if h = 1 and if b is Q or
2088  some field Z/pZ: */
2089  if (h==0)
2090  {
2091  if (((src->rep==n_rep_gap_rat) || (src->rep==n_rep_gap_gmp))
2092  && (nCoeff_is_Q(dst->extRing->cf) || nCoeff_is_Z(dst->extRing->cf)))
2093  return ntMap00; /// Q or Z --> Q(T)
2094  if (src->rep==n_rep_gmp)
2095  return ntMapZ0; /// Z --> K(T)
2096  if (nCoeff_is_Zp(src) && nCoeff_is_Q(bDst))
2097  return ntMapP0; /// Z/p --> Q(T)
2098  if (nCoeff_is_Q_or_BI(src) && nCoeff_is_Zp(bDst))
2099  return ntMap0P; /// Q --> Z/p(T)
2100  if (nCoeff_is_Zp(src) && nCoeff_is_Zp(bDst))
2101  {
2102  if (src->ch == dst->ch) return ntMapPP; /// Z/p --> Z/p(T)
2103  else return ntMapUP; /// Z/u --> Z/p(T)
2104  }
2105  if (nCoeff_is_Zn(src) && nCoeff_is_Zn(bDst))
2106  {
2107  if (mpz_cmp(src->modNumber,bDst->modNumber)==0) return ntMapPP; /// Z/p --> Z/p(T)
2108  }
2109  }
2110  if (h != 1) return NULL;
2111  //if ((!nCoeff_is_Zp(bDst)) && (!nCoeff_is_Q(bDst))) return NULL;
2112 
2113  /* Let T denote the sequence of transcendental extension variables, i.e.,
2114  K[t_1, ..., t_s] =: K[T];
2115  Let moreover, for any such sequence T, T' denote any subsequence of T
2116  of the form t_1, ..., t_w with w <= s. */
2117 
2118  if (rVar(src->extRing) > rVar(dst->extRing))
2119  return NULL;
2120 
2121  for (int i = 0; i < rVar(src->extRing); i++)
2122  if (strcmp(rRingVar(i, src->extRing), rRingVar(i, dst->extRing)) != 0)
2123  return NULL;
2124 
2125  if (src->type==n_transExt)
2126  {
2127  if (src->extRing->cf==dst->extRing->cf)
2128  return ntCopyMap; /// K(T') --> K(T)
2129  else
2130  return ntGenMap; /// K(T') --> K'(T)
2131  }
2132  else
2133  {
2134  if (src->extRing->cf==dst->extRing->cf)
2135  return ntCopyAlg; /// K(T') --> K(T)
2136  else
2137  return ntGenAlg; /// K(T') --> K'(T)
2138  }
2139 
2140  return NULL; /// default
2141 }
2142 #if 0
2143 nMapFunc ntSetMap_T(const coeffs src, const coeffs dst)
2144 {
2145  nMapFunc n=ntSetMap(src,dst);
2146  if (n==ntCopyAlg) printf("n=ntCopyAlg\n");
2147  else if (n==ntCopyMap) printf("n=ntCopyMap\n");
2148  else if (n==ntMapUP) printf("n=ntMapUP\n");
2149  else if (n==ntMap0P) printf("n=ntMap0P\n");
2150  else if (n==ntMapP0) printf("n=ntMapP0\n");
2151  else if (n==ntMap00) printf("n=ntMap00\n");
2152  else if (n==NULL) printf("n=NULL\n");
2153  else printf("n=?\n");
2154  return n;
2155 }
2156 #endif
2157 
2158 static void ntKillChar(coeffs cf)
2159 {
2160  rDecRefCnt(cf->extRing);
2161  if (cf->extRing->ref < 0)
2162  rDelete(cf->extRing);
2163 }
2164 static number ntConvFactoryNSingN( const CanonicalForm n, const coeffs cf)
2165 {
2166  if (n.isZero()) return NULL;
2167  poly p=convFactoryPSingP(n,ntRing);
2168  p_Normalize(p,ntRing);
2169  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
2170  NUM(result) = p;
2171  //DEN(result) = NULL; // done by omAlloc0Bin
2172  //COM(result) = 0; // done by omAlloc0Bin
2173  ntTest((number)result);
2174  return (number)result;
2175 }
2176 static CanonicalForm ntConvSingNFactoryN( number n, BOOLEAN /*setChar*/, const coeffs cf )
2177 {
2178  ntTest(n);
2179  if (IS0(n)) return CanonicalForm(0);
2180 
2181  fraction f = (fraction)n;
2182  return convSingPFactoryP(NUM(f),ntRing);
2183 }
2184 
2185 static int ntParDeg(number a, const coeffs cf)
2186 {
2187  ntTest(a);
2188  if (IS0(a)) return -1;
2189  fraction fa = (fraction)a;
2190  return cf->extRing->pFDeg(NUM(fa),cf->extRing);
2191 }
2192 
2193 /// return the specified parameter as a number in the given trans.ext.
2194 static number ntParameter(const int iParameter, const coeffs cf)
2195 {
2197 
2198  const ring R = cf->extRing;
2199  assume( R != NULL );
2200  assume( 0 < iParameter && iParameter <= rVar(R) );
2201 
2202  poly p = p_One(R); p_SetExp(p, iParameter, 1, R); p_Setm(p, R);
2203  p_Test(p,R);
2204 
2205  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
2206  NUM(f) = p;
2207  //DEN(f) = NULL;
2208  //COM(f) = 0;
2209 
2210  ntTest((number)f);
2211 
2212  return (number)f;
2213 }
2214 
2215 /// if m == var(i)/1 => return i,
2216 int ntIsParam(number m, const coeffs cf)
2217 {
2218  ntTest(m);
2220 
2221  const ring R = cf->extRing;
2222  assume( R != NULL );
2223 
2224  fraction f = (fraction)m;
2225 
2226  if( DEN(f) != NULL )
2227  return 0;
2228 
2229  return p_Var( NUM(f), R );
2230 }
2231 
2233 {
2234  static inline poly convert(const number& n)
2235  {
2236  // suitable for trans. ext. numbers that are fractions of polys
2237  return NUM((fraction)n); // return the numerator
2238  }
2239 };
2240 
2241 
2242 static void ntClearContent(ICoeffsEnumerator& numberCollectionEnumerator, number& c, const coeffs cf)
2243 {
2244  assume(cf != NULL);
2246  // all coeffs are given by fractions of polynomails over integers!!!
2247  // without denominators!!!
2248 
2249  const ring R = cf->extRing;
2250  assume(R != NULL);
2251  const coeffs Q = R->cf;
2252  assume(Q != NULL);
2253  assume(nCoeff_is_Q(Q));
2254 
2255 
2256  numberCollectionEnumerator.Reset();
2257 
2258  if( !numberCollectionEnumerator.MoveNext() ) // empty zero polynomial?
2259  {
2260  c = ntInit(1, cf);
2261  return;
2262  }
2263 
2264  // all coeffs are given by integers after returning from this routine
2265 
2266  // part 1, collect product of all denominators /gcds
2267  poly cand = NULL;
2268 
2269  do
2270  {
2271  number &n = numberCollectionEnumerator.Current();
2272 
2273  ntNormalize(n, cf);
2274 
2275  fraction f = (fraction)n;
2276 
2277  assume( f != NULL );
2278 
2279  const poly den = DEN(f);
2280 
2281  assume( den == NULL ); // ?? / 1 ?
2282 
2283  const poly num = NUM(f);
2284 
2285  if( cand == NULL )
2286  cand = p_Copy(num, R);
2287  else
2288  {
2289  poly tmp = singclap_gcd_r(cand, num, R); // gcd(cand, num)
2290  p_Delete(&cand,R);
2291  cand=tmp;
2292  }
2293 
2294  if( p_IsConstant(cand, R) )
2295  break;
2296  }
2297  while( numberCollectionEnumerator.MoveNext() ) ;
2298 
2299 
2300  // part2: all coeffs = all coeffs * cand
2301  if( cand != NULL )
2302  {
2303  if( !p_IsConstant(cand, R) )
2304  {
2305  c = ntInit(cand, cf);
2306  numberCollectionEnumerator.Reset();
2307  while (numberCollectionEnumerator.MoveNext() )
2308  {
2309  number &n = numberCollectionEnumerator.Current();
2310  const number t = ntDiv(n, c, cf); // TODO: rewrite!?
2311  ntDelete(&n, cf);
2312  n = t;
2313  }
2314  } // else NUM (result) = p_One(R);
2315  else { p_Delete(&cand, R); cand = NULL; }
2316  }
2317 
2318  // Quick and dirty fix for constant content clearing: consider numerators???
2319  CRecursivePolyCoeffsEnumerator<NTNumConverter> itr(numberCollectionEnumerator); // recursively treat the NUM(numbers) as polys!
2320  number cc;
2321 
2322  n_ClearContent(itr, cc, Q);
2323  number g = ntInit(p_NSet(cc, R), cf);
2324 
2325  if( cand != NULL )
2326  {
2327  number gg = ntMult(g, c, cf);
2328  ntDelete(&g, cf);
2329  ntDelete(&c, cf); c = gg;
2330  } else
2331  c = g;
2332  ntTest(c);
2333 }
2334 
2335 static void ntClearDenominators(ICoeffsEnumerator& numberCollectionEnumerator, number& c, const coeffs cf)
2336 {
2337  assume(cf != NULL);
2338  assume(getCoeffType(cf) == n_transExt); // both over Q(a) and Zp(a)!
2339  // all coeffs are given by fractions of polynomails over integers!!!
2340 
2341  numberCollectionEnumerator.Reset();
2342 
2343  if( !numberCollectionEnumerator.MoveNext() ) // empty zero polynomial?
2344  {
2345  c = ntInit(1, cf);
2346  return;
2347  }
2348 
2349  // all coeffs are given by integers after returning from this routine
2350 
2351  // part 1, collect product of all denominators /gcds
2352  poly cand = NULL;
2353 
2354  const ring R = cf->extRing;
2355  assume(R != NULL);
2356 
2357  const coeffs Q = R->cf;
2358  assume(Q != NULL);
2359 // assume(nCoeff_is_Q(Q));
2360 
2361  do
2362  {
2363  number &n = numberCollectionEnumerator.Current();
2364 
2365  ntNormalize(n, cf);
2366 
2367  fraction f = (fraction)ntGetDenom (n, cf);
2368 
2369  assume( f != NULL );
2370 
2371  const poly den = NUM(f);
2372 
2373  if( den == NULL ) // ?? / 1 ?
2374  continue;
2375 
2376  if( cand == NULL )
2377  cand = p_Copy(den, R);
2378  else
2379  {
2380  // cand === LCM( cand, den )!!!!
2381  // NOTE: maybe it's better to make the product and clearcontent afterwards!?
2382  // TODO: move the following to factory?
2383  poly gcd = singclap_gcd_r(cand, den, R); // gcd(cand, den) is monic no mater leading coeffs! :((((
2384  if (nCoeff_is_Q (Q))
2385  {
2386  number LcGcd= n_SubringGcd (p_GetCoeff (cand, R), p_GetCoeff(den, R), Q);
2387  gcd = __p_Mult_nn(gcd, LcGcd, R);
2388  n_Delete(&LcGcd,Q);
2389  }
2390 // assume( n_IsOne(pGetCoeff(gcd), Q) ); // TODO: this may be wrong...
2391  cand = p_Mult_q(cand, p_Copy(den, R), R); // cand *= den
2392  const poly t = singclap_pdivide( cand, gcd, R ); // cand' * den / gcd(cand', den)
2393  p_Delete(&cand, R);
2394  p_Delete(&gcd, R);
2395  cand = t;
2396  }
2397  }
2398  while( numberCollectionEnumerator.MoveNext() );
2399 
2400  if( cand == NULL )
2401  {
2402  c = ntInit(1, cf);
2403  return;
2404  }
2405 
2406  c = ntInit(cand, cf);
2407 
2408  numberCollectionEnumerator.Reset();
2409 
2410  number d = NULL;
2411 
2412  while (numberCollectionEnumerator.MoveNext() )
2413  {
2414  number &n = numberCollectionEnumerator.Current();
2415  number t = ntMult(n, c, cf); // TODO: rewrite!?
2416  ntDelete(&n, cf);
2417 
2418  ntNormalize(t, cf); // TODO: needed?
2419  n = t;
2420 
2421  fraction f = (fraction)t;
2422  assume( f != NULL );
2423 
2424  const poly den = DEN(f);
2425 
2426  if( den != NULL ) // ?? / ?? ?
2427  {
2428  assume( p_IsConstant(den, R) );
2429  assume( pNext(den) == NULL );
2430 
2431  if( d == NULL )
2432  d = n_Copy(pGetCoeff(den), Q);
2433  else
2434  {
2435  number g = n_NormalizeHelper(d, pGetCoeff(den), Q);
2436  n_Delete(&d, Q); d = g;
2437  }
2438  }
2439  }
2440 
2441  if( d != NULL )
2442  {
2443  numberCollectionEnumerator.Reset();
2444  while (numberCollectionEnumerator.MoveNext() )
2445  {
2446  number &n = numberCollectionEnumerator.Current();
2447  fraction f = (fraction)n;
2448 
2449  assume( f != NULL );
2450 
2451  const poly den = DEN(f);
2452 
2453  if( den == NULL ) // ?? / 1 ?
2454  NUM(f) = __p_Mult_nn(NUM(f), d, R);
2455  else
2456  {
2457  assume( p_IsConstant(den, R) );
2458  assume( pNext(den) == NULL );
2459 
2460  number ddd = n_Div(d, pGetCoeff(den), Q); // but be an integer now!!!
2461  NUM(f) = __p_Mult_nn(NUM(f), ddd, R);
2462  n_Delete(&ddd, Q);
2463 
2464  p_Delete(&DEN(f), R);
2465  DEN(f) = NULL; // TODO: check if this is needed!?
2466  }
2467 
2468  assume( DEN(f) == NULL );
2469  }
2470 
2471  NUM((fraction)c) = __p_Mult_nn(NUM((fraction)c), d, R);
2472  n_Delete(&d, Q);
2473  }
2474 
2475 
2476  ntTest(c);
2477 }
2478 
2479 static number ntChineseRemainder(number *x, number *q,int rl, BOOLEAN /*sym*/,CFArray &inv_cache,const coeffs cf)
2480 {
2481  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
2482 
2483  poly *P=(poly*)omAlloc(rl*sizeof(poly*));
2484  number *X=(number *)omAlloc(rl*sizeof(number));
2485 
2486  int i;
2487 
2488  for(i=0;i<rl;i++) P[i]=p_Copy(NUM((fraction)(x[i])),cf->extRing);
2489  NUM(result)=p_ChineseRemainder(P,X,q,rl,inv_cache,cf->extRing);
2490 
2491  for(i=0;i<rl;i++)
2492  {
2493  P[i]=p_Copy(DEN((fraction)(x[i])),cf->extRing);
2494  if (P[i]==NULL) P[i]=p_One(cf->extRing);
2495  }
2496  DEN(result)=p_ChineseRemainder(P,X,q,rl,inv_cache,cf->extRing);
2497 
2498  omFreeSize(X,rl*sizeof(number));
2499  omFreeSize(P,rl*sizeof(poly*));
2500  if (p_IsConstant(DEN(result), ntRing)
2501  && n_IsOne(pGetCoeff(DEN(result)), ntCoeffs))
2502  {
2503  p_Delete(&DEN(result),ntRing);
2504  }
2505  ntTest((number)result);
2506  return ((number)result);
2507 }
2508 
2509 static number ntFarey(number p, number n, const coeffs cf)
2510 {
2511  // n is really a bigint
2512  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
2513  NUM(result)=p_Farey(p_Copy(NUM((fraction)p),cf->extRing),n,cf->extRing);
2514  DEN(result)=p_Farey(p_Copy(DEN((fraction)p),cf->extRing),n,cf->extRing);
2515  ntTest((number)result);
2516  return ((number)result);
2517 }
2518 
2519 static number ntInitMPZ(mpz_t m, const coeffs r)
2520 {
2521  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
2522  number n=n_InitMPZ(m,r->extRing->cf);
2523  NUM(result)=p_NSet(n,r->extRing);
2524  return ((number)result);
2525 }
2526 
2527 static void ntMPZ(mpz_t m, number &n, const coeffs r)
2528 {
2529  mpz_init(m);
2530  if (n!=NULL)
2531  {
2532  fraction nn=(fraction)n;
2533  if (DENIS1(nn))
2534  {
2535  if (p_IsConstant(NUM(nn),r->extRing))
2536  {
2537  n_MPZ(m,pGetCoeff(NUM(nn)),r->extRing->cf);
2538  return;
2539  }
2540  }
2541  }
2542 }
2543 
2544 BOOLEAN ntInitChar(coeffs cf, void * infoStruct)
2545 {
2546 
2547  assume( infoStruct != NULL );
2548 
2549  TransExtInfo *e = (TransExtInfo *)infoStruct;
2550 
2551  assume( e->r != NULL); // extRing;
2552  assume( e->r->cf != NULL); // extRing->cf;
2553  assume( e->r->qideal == NULL );
2554 
2555  assume( cf != NULL );
2556  assume(getCoeffType(cf) == n_transExt); // coeff type;
2557 
2558  ring R = e->r;
2559  assume(R != NULL);
2560 
2561  rIncRefCnt(R); // increase the ref.counter for the ground poly. ring!
2562 
2563  cf->extRing = R;
2564  /* propagate characteristic up so that it becomes
2565  directly accessible in cf: */
2566  cf->ch = R->cf->ch;
2567 
2568  cf->is_field=TRUE;
2569  cf->is_domain=TRUE;
2570  cf->rep=n_rep_rat_fct;
2571 
2572  cf->factoryVarOffset = R->cf->factoryVarOffset + rVar(R);
2573 
2574  cf->cfCoeffName = naCoeffName; // FIXME? TODO? // extern char* naCoeffString(const coeffs r);
2575 
2576  cf->cfGreaterZero = ntGreaterZero;
2577  cf->cfGreater = ntGreater;
2578  cf->cfEqual = ntEqual;
2579  cf->cfIsZero = ntIsZero;
2580  cf->cfIsOne = ntIsOne;
2581  cf->cfIsMOne = ntIsMOne;
2582  cf->cfInit = ntInit;
2583  cf->cfFarey = ntFarey;
2584  cf->cfChineseRemainder = ntChineseRemainder;
2585  cf->cfInt = ntInt;
2586  cf->cfAdd = ntAdd;
2587  cf->cfInpNeg = ntNeg;
2588  cf->cfSub = ntSub;
2589  cf->cfMult = ntMult;
2590  cf->cfDiv = ntDiv;
2591  cf->cfExactDiv = ntExactDiv;
2592  cf->cfPower = ntPower;
2593  cf->cfCopy = ntCopy;
2594  cf->cfWriteLong = ntWriteLong;
2595  cf->cfRead = ntRead;
2596  cf->cfNormalize = ntNormalize;
2597  cf->cfDelete = ntDelete;
2598  cf->cfSetMap = ntSetMap;
2599  cf->cfGetDenom = ntGetDenom;
2600  cf->cfGetNumerator = ntGetNumerator;
2601  //cf->cfRePart = ntCopy;
2602  //cf->cfImPart = ntImPart;
2603  cf->cfCoeffWrite = ntCoeffWrite;
2604 #ifdef LDEBUG
2605  cf->cfDBTest = ntDBTest;
2606 #endif
2607  //cf->cfGcd = ntGcd_dummy;
2608  cf->cfSubringGcd = ntGcd;
2609  cf->cfNormalizeHelper = ntNormalizeHelper;
2610  cf->cfSize = ntSize;
2611  cf->nCoeffIsEqual = ntCoeffIsEqual;
2612  cf->cfInvers = ntInvers;
2613  cf->cfKillChar = ntKillChar;
2614  cf->cfInitMPZ = ntInitMPZ;
2615  cf->cfMPZ = ntMPZ;
2616 
2617  if( rCanShortOut(ntRing) )
2618  cf->cfWriteShort = ntWriteShort;
2619  else
2620  cf->cfWriteShort = ntWriteLong;
2621 
2622  cf->convFactoryNSingN =ntConvFactoryNSingN;
2623  cf->convSingNFactoryN =ntConvSingNFactoryN;
2624  cf->cfParDeg = ntParDeg;
2625 
2626  cf->iNumberOfParameters = rVar(R);
2627  cf->pParameterNames = (const char**)R->names;
2628  cf->cfParameter = ntParameter;
2629  cf->has_simple_Inverse= FALSE;
2630  /* cf->has_simple_Alloc= FALSE; */
2631 
2632 
2633  if( nCoeff_is_Q(R->cf) )
2634  cf->cfClearContent = ntClearContent;
2635 
2636  cf->cfClearDenominators = ntClearDenominators;
2637 
2638  return FALSE;
2639 }
2640 
2642 template class IEnumerator<snumber*>;
Rational pow(const Rational &a, int e)
Definition: GMPrat.cc:411
Concrete implementation of enumerators over polynomials.
char * naCoeffName(const coeffs r)
Definition: algext.cc:1335
All the auxiliary stuff.
int BOOLEAN
Definition: auxiliary.h:87
#define TRUE
Definition: auxiliary.h:100
#define FALSE
Definition: auxiliary.h:96
void * ADDRESS
Definition: auxiliary.h:119
CanonicalForm num(const CanonicalForm &f)
CanonicalForm den(const CanonicalForm &f)
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:56
int l
Definition: cfEzgcd.cc:100
int m
Definition: cfEzgcd.cc:128
int i
Definition: cfEzgcd.cc:132
int k
Definition: cfEzgcd.cc:99
const CanonicalForm const CanonicalForm const CanonicalForm const CanonicalForm & cand
Definition: cfModGcd.cc:69
Variable x
Definition: cfModGcd.cc:4082
int p
Definition: cfModGcd.cc:4078
g
Definition: cfModGcd.cc:4090
CanonicalForm cf
Definition: cfModGcd.cc:4083
CanonicalForm b
Definition: cfModGcd.cc:4103
FILE * f
Definition: checklibs.c:9
CanonicalForm convSingPFactoryP(poly p, const ring r)
Definition: clapconv.cc:136
poly convFactoryPSingP(const CanonicalForm &f, const ring r)
Definition: clapconv.cc:40
poly singclap_pdivide(poly f, poly g, const ring r)
Definition: clapsing.cc:601
poly singclap_gcd_and_divide(poly &f, poly &g, const ring r)
clears denominators of f and g, divides by gcd(f,g)
Definition: clapsing.cc:147
poly singclap_gcd_r(poly f, poly g, const ring r)
Definition: clapsing.cc:45
This is a polynomial enumerator for simple iteration over coefficients of polynomials.
go into polynomials over an alg. extension recursively
factory's main class
Definition: canonicalform.h:86
CF_NO_INLINE bool isZero() const
virtual reference Current()=0
Gets the current element in the collection (read and write).
virtual void Reset()=0
Sets the enumerator to its initial position: -1, which is before the first element in the collection.
virtual bool MoveNext()=0
Advances the enumerator to the next element of the collection. returns true if the enumerator was suc...
Templated enumerator interface for simple iteration over a generic collection of T's.
Definition: Enumerator.h:125
Coefficient rings, fields and other domains suitable for Singular polynomials.
static FORCE_INLINE number n_Mult(number a, number b, const coeffs r)
return the product of 'a' and 'b', i.e., a*b
Definition: coeffs.h:636
static FORCE_INLINE long n_Int(number &n, const coeffs r)
conversion of n to an int; 0 if not possible in Z/pZ: the representing int lying in (-p/2 ....
Definition: coeffs.h:547
static FORCE_INLINE number n_Copy(number n, const coeffs r)
return a copy of 'n'
Definition: coeffs.h:451
static FORCE_INLINE number n_NormalizeHelper(number a, number b, const coeffs r)
assume that r is a quotient field (otherwise, return 1) for arguments (a1/a2,b1/b2) return (lcm(a1,...
Definition: coeffs.h:695
static FORCE_INLINE number n_GetDenom(number &n, const coeffs r)
return the denominator of n (if elements of r are by nature not fractional, result is 1)
Definition: coeffs.h:603
static FORCE_INLINE void n_CoeffWrite(const coeffs r, BOOLEAN details=TRUE)
output the coeff description
Definition: coeffs.h:719
static FORCE_INLINE BOOLEAN nCoeff_is_Z(const coeffs r)
Definition: coeffs.h:816
static FORCE_INLINE BOOLEAN nCoeff_is_Extension(const coeffs r)
Definition: coeffs.h:846
#define n_Test(a, r)
BOOLEAN n_Test(number a, const coeffs r)
Definition: coeffs.h:712
n_coeffType
Definition: coeffs.h:27
@ n_Q
rational (GMP) numbers
Definition: coeffs.h:30
@ n_Zp
\F{p < 2^31}
Definition: coeffs.h:29
@ n_transExt
used for all transcendental extensions, i.e., the top-most extension in an extension tower is transce...
Definition: coeffs.h:38
static FORCE_INLINE number n_Gcd(number a, number b, const coeffs r)
in Z: return the gcd of 'a' and 'b' in Z/nZ, Z/2^kZ: computed as in the case Z in Z/pZ,...
Definition: coeffs.h:664
static FORCE_INLINE number n_Invers(number a, const coeffs r)
return the multiplicative inverse of 'a'; raise an error if 'a' is not invertible
Definition: coeffs.h:564
static FORCE_INLINE BOOLEAN n_GreaterZero(number n, const coeffs r)
ordered fields: TRUE iff 'n' is positive; in Z/pZ: TRUE iff 0 < m <= roundedBelow(p/2),...
Definition: coeffs.h:494
static FORCE_INLINE void n_MPZ(mpz_t result, number &n, const coeffs r)
conversion of n to a GMP integer; 0 if not possible
Definition: coeffs.h:551
static FORCE_INLINE BOOLEAN n_IsMOne(number n, const coeffs r)
TRUE iff 'n' represents the additive inverse of the one element, i.e. -1.
Definition: coeffs.h:472
static FORCE_INLINE nMapFunc n_SetMap(const coeffs src, const coeffs dst)
set the mapping function pointers for translating numbers from src to dst
Definition: coeffs.h:700
static FORCE_INLINE number n_InpNeg(number n, const coeffs r)
in-place negation of n MUST BE USED: n = n_InpNeg(n) (no copy is returned)
Definition: coeffs.h:557
static FORCE_INLINE BOOLEAN nCoeff_has_simple_inverse(const coeffs r)
TRUE, if the computation of the inverse is fast, i.e. prefer leading coeff. 1 over content.
Definition: coeffs.h:902
static FORCE_INLINE number n_Div(number a, number b, const coeffs r)
return the quotient of 'a' and 'b', i.e., a/b; raises an error if 'b' is not invertible in r exceptio...
Definition: coeffs.h:615
static FORCE_INLINE BOOLEAN nCoeff_is_Q(const coeffs r)
Definition: coeffs.h:806
static FORCE_INLINE BOOLEAN n_Greater(number a, number b, const coeffs r)
ordered fields: TRUE iff 'a' is larger than 'b'; in Z/pZ: TRUE iff la > lb, where la and lb are the l...
Definition: coeffs.h:511
static FORCE_INLINE BOOLEAN n_IsZero(number n, const coeffs r)
TRUE iff 'n' represents the zero element.
Definition: coeffs.h:464
static FORCE_INLINE void n_ClearDenominators(ICoeffsEnumerator &numberCollectionEnumerator, number &d, const coeffs r)
(inplace) Clears denominators on a collection of numbers number d is the LCM of all the coefficient d...
Definition: coeffs.h:935
static FORCE_INLINE BOOLEAN nCoeff_is_Q_or_BI(const coeffs r)
Definition: coeffs.h:829
static FORCE_INLINE n_coeffType getCoeffType(const coeffs r)
Returns the type of coeffs domain.
Definition: coeffs.h:421
static FORCE_INLINE void n_Delete(number *p, const coeffs r)
delete 'p'
Definition: coeffs.h:455
static FORCE_INLINE BOOLEAN nCoeff_is_Zn(const coeffs r)
Definition: coeffs.h:826
static FORCE_INLINE number n_InitMPZ(mpz_t n, const coeffs r)
conversion of a GMP integer to number
Definition: coeffs.h:542
static FORCE_INLINE BOOLEAN nCoeff_is_Zp(const coeffs r)
Definition: coeffs.h:800
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:538
static FORCE_INLINE void n_ClearContent(ICoeffsEnumerator &numberCollectionEnumerator, number &c, const coeffs r)
Computes the content and (inplace) divides it out on a collection of numbers number c is the content ...
Definition: coeffs.h:928
static FORCE_INLINE number n_GetNumerator(number &n, const coeffs r)
return the numerator of n (if elements of r are by nature not fractional, result is n)
Definition: coeffs.h:608
@ n_rep_gap_rat
(number), see longrat.h
Definition: coeffs.h:111
@ n_rep_gap_gmp
(), see rinteger.h, new impl.
Definition: coeffs.h:112
@ n_rep_gmp
(mpz_ptr), see rmodulon,h
Definition: coeffs.h:115
@ n_rep_rat_fct
(fraction), see transext.h
Definition: coeffs.h:114
static FORCE_INLINE number n_SubringGcd(number a, number b, const coeffs r)
Definition: coeffs.h:666
number(* nMapFunc)(number a, const coeffs src, const coeffs dst)
maps "a", which lives in src, into dst
Definition: coeffs.h:73
static FORCE_INLINE void n_Normalize(number &n, const coeffs r)
inplace-normalization of n; produces some canonical representation of n;
Definition: coeffs.h:578
static FORCE_INLINE BOOLEAN n_IsOne(number n, const coeffs r)
TRUE iff 'n' represents the one element.
Definition: coeffs.h:468
BOOLEAN pa(leftv res, leftv args)
Definition: cohomo.cc:4344
BOOLEAN pb(leftv res, leftv args)
Definition: cohomo.cc:4371
BOOLEAN fa(leftv res, leftv args)
Definition: cohomo.cc:4390
BOOLEAN fb(leftv res, leftv args)
Definition: cohomo.cc:4331
#define Print
Definition: emacs.cc:80
return result
Definition: facAbsBiFact.cc:75
const CanonicalForm int s
Definition: facAbsFact.cc:51
CanonicalForm res
Definition: facAbsFact.cc:60
CanonicalForm factor
Definition: facAbsFact.cc:97
void WerrorS(const char *s)
Definition: feFopen.cc:24
#define VAR
Definition: globaldefs.h:5
STATIC_VAR Poly * h
Definition: janet.cc:971
STATIC_VAR jList * Q
Definition: janet.cc:30
poly p_ChineseRemainder(poly *xx, mpz_ptr *x, mpz_ptr *q, int rl, mpz_ptr *C, const ring R)
number nlModP(number q, const coeffs, const coeffs Zp)
Definition: longrat.cc:1577
#define SR_INT
Definition: longrat.h:67
static FORCE_INLINE BOOLEAN nlIsInteger(number q, const coeffs r)
Definition: longrat.h:94
#define assume(x)
Definition: mod2.h:387
#define pIter(p)
Definition: monomials.h:37
#define pNext(p)
Definition: monomials.h:36
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy
Definition: monomials.h:44
#define p_GetCoeff(p, r)
Definition: monomials.h:50
gmp_float exp(const gmp_float &a)
Definition: mpr_complex.cc:357
The main handler for Singular numbers which are suitable for Singular polynomials.
const char *const nDivBy0
Definition: numbers.h:87
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260
#define omAlloc(size)
Definition: omAllocDecl.h:210
#define omAllocBin(bin)
Definition: omAllocDecl.h:205
#define omAlloc0Bin(bin)
Definition: omAllocDecl.h:206
#define omFreeBin(addr, bin)
Definition: omAllocDecl.h:259
#define omGetSpecBin(size)
Definition: omBin.h:11
#define NULL
Definition: omList.c:12
omBin_t * omBin
Definition: omStructs.h:12
poly p_Diff(poly a, int k, const ring r)
Definition: p_polys.cc:1890
poly p_Farey(poly p, number N, const ring r)
Definition: p_polys.cc:54
poly p_ISet(long i, const ring r)
returns the poly representing the integer i
Definition: p_polys.cc:1293
poly p_Div_nn(poly p, const number n, const ring r)
Definition: p_polys.cc:1497
void p_Normalize(poly p, const ring r)
Definition: p_polys.cc:3847
void p_Norm(poly p1, const ring r)
Definition: p_polys.cc:3793
const char * p_Read(const char *st, poly &rc, const ring r)
Definition: p_polys.cc:1366
poly p_Cleardenom(poly p, const ring r)
Definition: p_polys.cc:2906
int p_Var(poly m, const ring r)
Definition: p_polys.cc:4689
poly p_One(const ring r)
Definition: p_polys.cc:1309
poly p_Sub(poly p1, poly p2, const ring r)
Definition: p_polys.cc:1982
poly p_NSet(number n, const ring r)
returns the poly representing the number n, destroys n
Definition: p_polys.cc:1465
BOOLEAN p_EqualPolys(poly p1, poly p2, const ring r)
Definition: p_polys.cc:4545
static poly p_Neg(poly p, const ring r)
Definition: p_polys.h:1079
static poly p_Add_q(poly p, poly q, const ring r)
Definition: p_polys.h:908
static poly p_Mult_q(poly p, poly q, const ring r)
Definition: p_polys.h:1086
void p_Write(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:342
static unsigned long p_SetExp(poly p, const unsigned long e, const unsigned long iBitmask, const int VarOffset)
set a single variable exponent @Note: VarOffset encodes the position in p->exp
Definition: p_polys.h:488
static void p_ExpVectorDiff(poly pr, poly p1, poly p2, const ring r)
Definition: p_polys.h:1446
void p_String0Long(const poly p, ring lmRing, ring tailRing)
print p in a long way
Definition: polys0.cc:203
void p_String0Short(const poly p, ring lmRing, ring tailRing)
print p in a short way, if possible
Definition: polys0.cc:184
static void p_Setm(poly p, const ring r)
Definition: p_polys.h:233
static number p_SetCoeff(poly p, number n, ring r)
Definition: p_polys.h:412
static BOOLEAN p_LmIsConstantComp(const poly p, const ring r)
Definition: p_polys.h:978
static long p_GetExp(const poly p, const unsigned long iBitmask, const int VarOffset)
get a single variable exponent @Note: the integer VarOffset encodes:
Definition: p_polys.h:469
static BOOLEAN p_LmIsConstant(const poly p, const ring r)
Definition: p_polys.h:995
static BOOLEAN p_IsOne(const poly p, const ring R)
either poly(1) or gen(k)?!
Definition: p_polys.h:1990
static BOOLEAN p_IsConstant(const poly p, const ring r)
Definition: p_polys.h:1983
static void p_Delete(poly *p, const ring r)
Definition: p_polys.h:873
static unsigned pLength(poly a)
Definition: p_polys.h:191
static poly pp_Mult_qq(poly p, poly q, const ring r)
Definition: p_polys.h:1123
static poly p_LmFreeAndNext(poly p, ring)
Definition: p_polys.h:703
static poly p_Copy(poly p, const ring r)
returns a copy of p
Definition: p_polys.h:818
static long p_Totaldegree(poly p, const ring r)
Definition: p_polys.h:1479
#define p_Test(p, r)
Definition: p_polys.h:162
#define __p_Mult_nn(p, n, r)
Definition: p_polys.h:943
void p_wrp(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:373
poly prMapR(poly src, nMapFunc nMap, ring src_r, ring dest_r)
Definition: prCopy.cc:45
poly prCopyR(poly p, ring src_r, ring dest_r)
Definition: prCopy.cc:34
@ NUM
Definition: readcf.cc:170
void StringAppendS(const char *st)
Definition: reporter.cc:107
void PrintS(const char *s)
Definition: reporter.cc:284
void PrintLn()
Definition: reporter.cc:310
void rDelete(ring r)
unconditionally deletes fields in r
Definition: ring.cc:449
BOOLEAN rEqual(ring r1, ring r2, BOOLEAN qr)
returns TRUE, if r1 equals r2 FALSE, otherwise Equality is determined componentwise,...
Definition: ring.cc:1660
static char * rRingVar(short i, const ring r)
Definition: ring.h:578
static ring rIncRefCnt(ring r)
Definition: ring.h:843
static void rDecRefCnt(ring r)
Definition: ring.h:844
static BOOLEAN rCanShortOut(const ring r)
Definition: ring.h:587
static short rVar(const ring r)
#define rVar(r) (r->N)
Definition: ring.h:593
#define R
Definition: sirandom.c:27
#define A
Definition: sirandom.c:24
static poly convert(const number &n)
Definition: transext.cc:2234
#define loop
Definition: structs.h:79
#define SR_HDL(A)
Definition: tgb.cc:35
static int ntSize(number a, const coeffs cf)
Definition: transext.cc:1811
static BOOLEAN ntDBTest(number a, const char *f, const int l, const coeffs r)
Definition: transext.cc:140
static number ntParameter(const int iParameter, const coeffs cf)
return the specified parameter as a number in the given trans.ext.
Definition: transext.cc:2194
static void ntWriteLong(number a, const coeffs cf)
Definition: transext.cc:1549
static BOOLEAN ntCoeffIsEqual(const coeffs cf, n_coeffType n, void *param)
Definition: transext.cc:1635
static void heuristicGcdCancellation(number a, const coeffs cf)
forward declarations
Definition: transext.cc:1306
number ntDiff(number a, number d, const coeffs cf)
Definition: transext.cc:897
nMapFunc ntSetMap(const coeffs src, const coeffs dst)
Get a mapping function from src into the domain of this type (n_transExt)
Definition: transext.cc:2078
static number ntSub(number a, number b, const coeffs cf)
Definition: transext.cc:995
static number ntCopyAlg(number a, const coeffs cf, const coeffs dst)
Definition: transext.cc:2000
static long ntInt(number &a, const coeffs cf)
Definition: transext.cc:773
static number ntMapPP(number a, const coeffs src, const coeffs dst)
Definition: transext.cc:2041
static void ntDelete(number *a, const coeffs cf)
Definition: transext.cc:313
static BOOLEAN ntIsOne(number a, const coeffs cf)
Definition: transext.cc:669
static void ntNormalizeDen(fraction result, const ring R)
Definition: transext.cc:1102
poly gcd_over_Q(poly f, poly g, const ring r)
helper routine for calling singclap_gcd_r
Definition: transext.cc:275
BOOLEAN ntInitChar(coeffs cf, void *infoStruct)
Initialize the coeffs object.
Definition: transext.cc:2544
static BOOLEAN ntIsZero(number a, const coeffs cf)
Definition: transext.cc:306
static number ntGetNumerator(number &a, const coeffs cf)
TODO: normalization of a!?
Definition: transext.cc:506
static CanonicalForm ntConvSingNFactoryN(number n, BOOLEAN, const coeffs cf)
Definition: transext.cc:2176
static void ntClearDenominators(ICoeffsEnumerator &numberCollectionEnumerator, number &c, const coeffs cf)
Definition: transext.cc:2335
#define ntRing
Definition: transext.cc:80
static BOOLEAN ntIsMOne(number a, const coeffs cf)
Definition: transext.cc:678
static number ntMapZ0(number a, const coeffs src, const coeffs dst)
Definition: transext.cc:1862
static void ntMPZ(mpz_t m, number &n, const coeffs r)
Definition: transext.cc:2527
static void ntCoeffWrite(const coeffs cf, BOOLEAN details)
Definition: transext.cc:856
static number ntMap0P(number a, const coeffs src, const coeffs dst)
Definition: transext.cc:2017
#define NUMIS1(f)
TRUE iff num. represents 1.
Definition: transext.cc:67
static number ntGcd(number a, number b, const coeffs cf)
Definition: transext.cc:1738
static number ntConvFactoryNSingN(const CanonicalForm n, const coeffs cf)
Definition: transext.cc:2164
static void ntClearContent(ICoeffsEnumerator &numberCollectionEnumerator, number &c, const coeffs cf)
Definition: transext.cc:2242
#define ntCoeffs
Definition: transext.cc:86
static number ntInvers(number a, const coeffs cf)
Definition: transext.cc:1172
static number ntNormalizeHelper(number a, number b, const coeffs cf)
Definition: transext.cc:1657
static void ntPower(number a, int exp, number *b, const coeffs cf)
Definition: transext.cc:1238
number ntInit(long i, const coeffs cf)
Definition: transext.cc:704
#define COM(f)
Definition: transext.cc:69
static number ntMapP0(number a, const coeffs src, const coeffs dst)
Definition: transext.cc:1876
static number ntGetDenom(number &a, const coeffs cf)
TODO: normalization of a!?
Definition: transext.cc:567
static const char * ntRead(const char *s, number *a, const coeffs cf)
Definition: transext.cc:1599
static number ntMapUP(number a, const coeffs src, const coeffs dst)
Definition: transext.cc:2056
#define MULT_COMPLEXITY
complexity increase due to * and /
Definition: transext.cc:62
static void ntKillChar(coeffs cf)
Definition: transext.cc:2158
static void ntWriteShort(number a, const coeffs cf)
Definition: transext.cc:1574
static number ntDiv(number a, number b, const coeffs cf)
Definition: transext.cc:1126
static number ntCopy(number a, const coeffs cf)
Definition: transext.cc:372
static int ntParDeg(number a, const coeffs cf)
Definition: transext.cc:2185
static BOOLEAN ntGreaterZero(number a, const coeffs cf)
Definition: transext.cc:796
static number ntChineseRemainder(number *x, number *q, int rl, BOOLEAN, CFArray &inv_cache, const coeffs cf)
Definition: transext.cc:2479
static number ntGenMap(number a, const coeffs cf, const coeffs dst)
Definition: transext.cc:1921
int ntIsParam(number m, const coeffs cf)
if m == var(i)/1 => return i,
Definition: transext.cc:2216
static number ntMap00(number a, const coeffs src, const coeffs dst)
Definition: transext.cc:1838
static number ntMult(number a, number b, const coeffs cf)
Definition: transext.cc:1034
#define ADD_COMPLEXITY
complexity increase due to + and -
Definition: transext.cc:61
#define ntTest(a)
Definition: transext.cc:76
static number ntGenAlg(number a, const coeffs cf, const coeffs dst)
Definition: transext.cc:2007
static number ntNeg(number a, const coeffs cf)
this is in-place, modifies a
Definition: transext.cc:691
static void definiteGcdCancellation(number a, const coeffs cf, BOOLEAN simpleTestsHaveAlreadyBeenPerformed)
modifies a
Definition: transext.cc:1391
static number ntAdd(number a, number b, const coeffs cf)
Definition: transext.cc:954
static number ntFarey(number p, number n, const coeffs cf)
Definition: transext.cc:2509
static number ntCopyMap(number a, const coeffs cf, const coeffs dst)
Definition: transext.cc:1892
static void ntNormalize(number &a, const coeffs cf)
Definition: transext.cc:1609
static BOOLEAN ntEqual(number a, number b, const coeffs cf)
Definition: transext.cc:326
VAR omBin fractionObjectBin
Definition: transext.cc:89
static coeffs nCoeff_bottom(const coeffs r, int &height)
Definition: transext.cc:292
#define DIFF_COMPLEXITY
complexity increase due to diff
Definition: transext.cc:63
static BOOLEAN ntGreater(number a, number b, const coeffs cf)
Definition: transext.cc:806
static void handleNestedFractionsOverQ(fraction f, const coeffs cf)
Definition: transext.cc:413
static number ntExactDiv(number a, number b, const coeffs cf)
Definition: transext.cc:1627
#define BOUND_COMPLEXITY
maximum complexity of a number
Definition: transext.cc:64
static number ntInitMPZ(mpz_t m, const coeffs r)
Definition: transext.cc:2519
struct for passing initialization parameters to naInitChar
Definition: transext.h:88
int gcd(int a, int b)
Definition: walkSupport.cc:836