source: git/libpolys/polys/ext_fields/transext.cc @ 96847c

spielwiese
Last change on this file since 96847c was 96847c, checked in by Janko Boehm <boehm@…>, 11 years ago
added tests for ntDiff and ellipticcovers.lib
  • Property mode set to 100644
File size: 61.8 KB
RevLine 
[6ccdd3a]1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
[2c7f28]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
[06df101]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. (This definite cancel-
24*           lation will also be performed at the beginning of ntWrite,
25*           ensuring that any output is free of common factors.
26*           For the special case of K = Q (i.e., when computing over the
27*           rationals), this definite cancellation procedure will also take
28*           care of nested fractions: If there are fractional coefficients
29*           in the numerator or denominator of a number, then this number
30*           is being replaced by a quotient of two polynomials over Z, or
31*           - if the denominator is a constant - by a polynomial over Q.
[dc79bd]32*
33*           TODO: the description above needs a major update!!!
[6ccdd3a]34*/
[12f6ce0]35
36
37
38
39
[1f414c8]40#define TRANSEXT_PRIVATES
[6ccdd3a]41
[16f511]42#ifdef HAVE_CONFIG_H
[ba5e9e]43#include "libpolysconfig.h"
[16f511]44#endif /* HAVE_CONFIG_H */
[6ccdd3a]45#include <misc/auxiliary.h>
46
47#include <omalloc/omalloc.h>
48
49#include <reporter/reporter.h>
50
51#include <coeffs/coeffs.h>
52#include <coeffs/numbers.h>
53#include <coeffs/longrat.h>
54
55#include <polys/monomials/ring.h>
56#include <polys/monomials/p_polys.h>
57#include <polys/simpleideals.h>
58
[e5d267]59#ifdef HAVE_FACTORY
60#include <polys/clapsing.h>
[fc4977]61#include <polys/clapconv.h>
62#include <factory/factory.h>
[e5d267]63#endif
64
[dc79bd]65#include <polys/ext_fields/transext.h>
66#include <polys/prCopy.h>
67
68#include <polys/PolyEnumerator.h>
69
[1f414c8]70
71/* constants for controlling the complexity of numbers */
72#define ADD_COMPLEXITY 1   /**< complexity increase due to + and - */
73#define MULT_COMPLEXITY 2   /**< complexity increase due to * and / */
[12f6ce0]74#define DIFF_COMPLEXITY 2   /**< complexity increase due to * and / */
[1f414c8]75#define BOUND_COMPLEXITY 10   /**< maximum complexity of a number */
76
[dc79bd]77/// TRUE iff num. represents 1
78#define NUMIS1(f) (p_IsOne(NUM(f), cf->extRing))
79
[1f414c8]80#define COM(f) f->complexity
81
82
83#ifdef LDEBUG
[dc79bd]84#define ntTest(a) assume(ntDBTest(a,__FILE__,__LINE__,cf))
[1f414c8]85BOOLEAN  ntDBTest(number a, const char *f, const int l, const coeffs r);
86#else
[d101b1]87#define ntTest(a) do {} while (0)
[1f414c8]88#endif
89
90/// Our own type!
91static const n_coeffType ID = n_transExt;
92
93/* polynomial ring in which the numerators and denominators of our
94   numbers live */
95#define ntRing cf->extRing
96
97/* coeffs object in which the coefficients of our numbers live;
98 * methods attached to ntCoeffs may be used to compute with the
99 * coefficients of our numbers, e.g., use ntCoeffs->nAdd to add
100 * coefficients of our numbers */
101#define ntCoeffs cf->extRing->cf
102
[6ccdd3a]103
[fc4977]104omBin fractionObjectBin = omGetSpecBin(sizeof(fractionObject));
[de90c01]105
[6ccdd3a]106/// forward declarations
[c14846c]107BOOLEAN  ntGreaterZero(number a, const coeffs cf);
[6ccdd3a]108BOOLEAN  ntGreater(number a, number b, const coeffs cf);
109BOOLEAN  ntEqual(number a, number b, const coeffs cf);
110BOOLEAN  ntIsOne(number a, const coeffs cf);
111BOOLEAN  ntIsMOne(number a, const coeffs cf);
112BOOLEAN  ntIsZero(number a, const coeffs cf);
[2f3764]113number   ntInit(long i, const coeffs cf);
[6ccdd3a]114int      ntInt(number &a, const coeffs cf);
115number   ntNeg(number a, const coeffs cf);
116number   ntInvers(number a, const coeffs cf);
117number   ntAdd(number a, number b, const coeffs cf);
118number   ntSub(number a, number b, const coeffs cf);
119number   ntMult(number a, number b, const coeffs cf);
120number   ntDiv(number a, number b, const coeffs cf);
121void     ntPower(number a, int exp, number *b, const coeffs cf);
122number   ntCopy(number a, const coeffs cf);
[ce1f78]123void     ntWriteLong(number &a, const coeffs cf);
124void     ntWriteShort(number &a, const coeffs cf);
[6ccdd3a]125number   ntRePart(number a, const coeffs cf);
126number   ntImPart(number a, const coeffs cf);
127number   ntGetDenom(number &a, const coeffs cf);
128number   ntGetNumerator(number &a, const coeffs cf);
129number   ntGcd(number a, number b, const coeffs cf);
130number   ntLcm(number a, number b, const coeffs cf);
[2c7f28]131int      ntSize(number a, const coeffs cf);
[6ccdd3a]132void     ntDelete(number * a, const coeffs cf);
[03f7b5]133void     ntCoeffWrite(const coeffs cf, BOOLEAN details);
[6ccdd3a]134number   ntIntDiv(number a, number b, const coeffs cf);
135const char * ntRead(const char *s, number *a, const coeffs cf);
136static BOOLEAN ntCoeffIsEqual(const coeffs cf, n_coeffType n, void * param);
137
[2c7f28]138void heuristicGcdCancellation(number a, const coeffs cf);
[010f3b]139void definiteGcdCancellation(number a, const coeffs cf,
[06df101]140                             BOOLEAN simpleTestsHaveAlreadyBeenPerformed);
141void handleNestedFractionsOverQ(fraction f, const coeffs cf);
[2c7f28]142
[6ccdd3a]143#ifdef LDEBUG
144BOOLEAN ntDBTest(number a, const char *f, const int l, const coeffs cf)
145{
[1f414c8]146  assume(getCoeffType(cf) == ID);
[dc79bd]147
148  if (IS0(a)) return TRUE;
149
150  const fraction t = (fraction)a;
151
152  const poly num = NUM(t);
153  assume(num != NULL);   /**< t != 0 ==> numerator(t) != 0 */
[f7286a]154  assume( _p_Test(num, ntRing,1) );
[dc79bd]155
156  const poly den = DEN(t);
[a355723]157
[dc79bd]158  if (den != NULL) // !DENIS1(f)
[1374bc]159  {
[f7286a]160    assume( _p_Test(den, ntRing,1) );
[a355723]161
[e5c9e5]162    if(p_IsConstant(den, ntRing) && (n_IsOne(pGetCoeff(den), ntCoeffs)))
[1374bc]163    {
164      Print("?/1 in %s:%d\n",f,l);
[dc79bd]165      return FALSE;
166    }
[a355723]167
[e5c9e5]168    if( !n_GreaterZero(pGetCoeff(den), ntCoeffs) )
[dc79bd]169    {
170      Print("negative sign of DEN. of a fraction in %s:%d\n",f,l);
171      return FALSE;
[1374bc]172    }
[a355723]173
[dc79bd]174    // test that den is over integers!?
[a355723]175
[dc79bd]176  } else
177  {  // num != NULL // den == NULL
[a355723]178
[dc79bd]179//    if( COM(t) != 0 )
180//    {
[a355723]181//      Print("?//NULL with non-zero complexity: %d in %s:%d\n", COM(t), f, l);
[dc79bd]182//      return FALSE;
183//    }
184    // test that nume is over integers!?
[1374bc]185  }
[a44bcf]186  if (getCoeffType(ntCoeffs)==n_Q)
187  {
188    poly p=num; // !=NULL
189    do
190    {
191      number n=pGetCoeff(p);
192      n_Test(n,ntCoeffs);
[2e7132]193      if ((!(SR_HDL(n) & SR_INT))&&(n->s==0))
194      /* not normalized, just do for the following test*/
195      {
196        n_Normalize(pGetCoeff(p),ntCoeffs);
197        n=pGetCoeff(p);
198      }
[a44bcf]199      if (!(SR_HDL(n) & SR_INT))
200      {
201        if (n->s<2)
202          Print("rational coeff in num: %s:%d\n",f,l);
203      }
204      pIter(p);
205    } while(p!=NULL);
206    p=den;
207    while(p!=NULL)
208    {
209      number n=pGetCoeff(p);
210      if (!(SR_HDL(n) & SR_INT))
211      {
212        if (n->s!=3)
213          Print("rational coeff in den.:%s:%d\n",f,l);
214      }
215      pIter(p);
216    }
217  }
[6ccdd3a]218  return TRUE;
219}
220#endif
221
222/* returns the bottom field in this field extension tower; if the tower
223   is flat, i.e., if there is no extension, then r itself is returned;
224   as a side-effect, the counter 'height' is filled with the height of
225   the extension tower (in case the tower is flat, 'height' is zero) */
226static coeffs nCoeff_bottom(const coeffs r, int &height)
227{
228  assume(r != NULL);
229  coeffs cf = r;
230  height = 0;
231  while (nCoeff_is_Extension(cf))
232  {
233    assume(cf->extRing != NULL); assume(cf->extRing->cf != NULL);
234    cf = cf->extRing->cf;
235    height++;
236  }
237  return cf;
238}
239
[2c7f28]240BOOLEAN ntIsZero(number a, const coeffs cf)
[6ccdd3a]241{
[dc79bd]242  ntTest(a); // !!!
[e5d267]243  return (IS0(a));
[6ccdd3a]244}
245
[2c7f28]246void ntDelete(number * a, const coeffs cf)
[6ccdd3a]247{
[dc79bd]248  ntTest(*a); // !!!
[2c7f28]249  fraction f = (fraction)(*a);
[e5d267]250  if (IS0(f)) return;
251  p_Delete(&NUM(f), ntRing);
252  if (!DENIS1(f)) p_Delete(&DEN(f), ntRing);
[2c7f28]253  omFreeBin((ADDRESS)f, fractionObjectBin);
[6ccdd3a]254  *a = NULL;
255}
256
[2c7f28]257BOOLEAN ntEqual(number a, number b, const coeffs cf)
[6ccdd3a]258{
[dc79bd]259  ntTest(a);
260  ntTest(b);
[c14846c]261
[6ccdd3a]262  /// simple tests
263  if (a == b) return TRUE;
[e5d267]264  if ((IS0(a)) && (!IS0(b))) return FALSE;
265  if ((IS0(b)) && (!IS0(a))) return FALSE;
[c14846c]266
267  /// cheap test if gcd's have been cancelled in both numbers
[2c7f28]268  fraction fa = (fraction)a;
269  fraction fb = (fraction)b;
[e5d267]270  if ((COM(fa) == 1) && (COM(fb) == 1))
[2c7f28]271  {
[e5d267]272    poly f = p_Add_q(p_Copy(NUM(fa), ntRing),
273                     p_Neg(p_Copy(NUM(fb), ntRing), ntRing),
[2c7f28]274                     ntRing);
275    if (f != NULL) { p_Delete(&f, ntRing); return FALSE; }
[e5d267]276    if (DENIS1(fa) && DENIS1(fb))  return TRUE;
277    if (DENIS1(fa) && !DENIS1(fb)) return FALSE;
278    if (!DENIS1(fa) && DENIS1(fb)) return FALSE;
279    f = p_Add_q(p_Copy(DEN(fa), ntRing),
280                p_Neg(p_Copy(DEN(fb), ntRing), ntRing),
[2c7f28]281                ntRing);
282    if (f != NULL) { p_Delete(&f, ntRing); return FALSE; }
283    return TRUE;
284  }
[c14846c]285
[2c7f28]286  /* default: the more expensive multiplication test
287              a/b = c/d  <==>  a*d = b*c */
[e5d267]288  poly f = p_Copy(NUM(fa), ntRing);
289  if (!DENIS1(fb)) f = p_Mult_q(f, p_Copy(DEN(fb), ntRing), ntRing);
290  poly g = p_Copy(NUM(fb), ntRing);
291  if (!DENIS1(fa)) g = p_Mult_q(g, p_Copy(DEN(fa), ntRing), ntRing);
[2c7f28]292  poly h = p_Add_q(f, p_Neg(g, ntRing), ntRing);
293  if (h == NULL) return TRUE;
294  else
295  {
296    p_Delete(&h, ntRing);
297    return FALSE;
298  }
[6ccdd3a]299}
300
[2c7f28]301number ntCopy(number a, const coeffs cf)
[6ccdd3a]302{
[dc79bd]303  ntTest(a); // !!!
[e5d267]304  if (IS0(a)) return NULL;
[2c7f28]305  fraction f = (fraction)a;
[e5d267]306  poly g = p_Copy(NUM(f), ntRing);
307  poly h = NULL; if (!DENIS1(f)) h = p_Copy(DEN(f), ntRing);
[2c7f28]308  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[e5d267]309  NUM(result) = g;
310  DEN(result) = h;
311  COM(result) = COM(f);
[dc79bd]312  ntTest((number)result);
[2c7f28]313  return (number)result;
[6ccdd3a]314}
315
[dc79bd]316/// TODO: normalization of a!?
[2c7f28]317number ntGetNumerator(number &a, const coeffs cf)
[6ccdd3a]318{
[2c7f28]319  ntTest(a);
[e5d267]320  if (IS0(a)) return NULL;
[dc79bd]321
[9d68fd]322  definiteGcdCancellation(a, cf, FALSE);
323
[2c7f28]324  fraction f = (fraction)a;
325  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[dc79bd]326
327  const BOOLEAN denis1= DENIS1 (f);
328
[28216f]329  if (getCoeffType (ntCoeffs) == n_Q && !denis1)
330    handleNestedFractionsOverQ (f, cf);
[dc79bd]331
[28216f]332  if (getCoeffType (ntCoeffs) == n_Q && denis1)
333  {
[dc79bd]334    assume( DEN (f) == NULL );
[a355723]335
[dc79bd]336    number g;
337    // TODO/NOTE: the following should not be necessary (due to
[a355723]338    // Hannes!) as NUM (f) should be over Z!!!
[dc79bd]339    CPolyCoeffsEnumerator itr(NUM(f));
340
341
[e5c9e5]342    n_ClearDenominators(itr, g, ntCoeffs);
[dc79bd]343
[e5c9e5]344    if( !n_GreaterZero(g, ntCoeffs) )
[28216f]345    {
[dc79bd]346      NUM (f) = p_Neg(NUM (f), ntRing); // Ugly :(((
[e5c9e5]347      g = n_Neg(g, ntCoeffs);
[28216f]348    }
[dc79bd]349
350    // g should be a positive integer now!
[e5c9e5]351    assume( n_GreaterZero(g, ntCoeffs) );
[a355723]352
[e5c9e5]353    if( !n_IsOne(g, ntCoeffs) )
[dc79bd]354    {
355      DEN (f) = p_NSet(g, ntRing); // update COM(f)???
356      COM (f) ++;
[a355723]357      assume( DEN (f) != NULL );
[dc79bd]358    }
359    else
[e5c9e5]360      n_Delete(&g, ntCoeffs);
[dc79bd]361
362    ntTest(a);
[28216f]363  }
[dc79bd]364
365  // Call ntNormalize instead of above?!?
[a355723]366
[dc79bd]367  NUM (result) = p_Copy (NUM (f), ntRing); // ???
368  DEN (result) = NULL;
369  COM (result) = 0;
[a355723]370
[dc79bd]371  ntTest((number)result);
[2c7f28]372  return (number)result;
[6ccdd3a]373}
374
[dc79bd]375/// TODO: normalization of a!?
[2c7f28]376number ntGetDenom(number &a, const coeffs cf)
[6ccdd3a]377{
[2c7f28]378  ntTest(a);
[dc79bd]379
[2c7f28]380  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[dc79bd]381  DEN (result)= NULL;
382  COM (result)= 0;
383
[9501be]384  if (IS0(a)) 
[6b32a6b]385  {
[9501be]386    NUM (result) = p_One(ntRing);
[6b32a6b]387    return (number)result;
[9501be]388  }
389     
390  definiteGcdCancellation(a, cf, FALSE);
391   
392  fraction f = (fraction)a;
393   
394  assume( !IS0(f) );
395
[dc79bd]396  const BOOLEAN denis1 = DENIS1 (f);
[6b32a6b]397 
[9501be]398  if( denis1 && (getCoeffType (ntCoeffs) != n_Q) ) // */1 or 0
[28216f]399  {
400    NUM (result)= p_One(ntRing);
[dc79bd]401    ntTest((number)result);
402    return (number)result;
[28216f]403  }
[dc79bd]404
405  if (!denis1) // */* / Q
[28216f]406  {
[dc79bd]407    assume( DEN (f) != NULL );
408
[28216f]409    if (getCoeffType (ntCoeffs) == n_Q)
410      handleNestedFractionsOverQ (f, cf);
[dc79bd]411
412    ntTest(a);
413
414    if( DEN (f) != NULL ) // is it ?? // 1 now???
415    {
416      assume( !p_IsOne(DEN (f), ntRing) );
417
418      NUM (result) = p_Copy (DEN (f), ntRing);
419      ntTest((number)result);
420      return (number)result;
421    }
422//    NUM (result) = p_One(ntRing); // NOTE: just in order to be sure...
423  }
[a355723]424
[dc79bd]425  // */1 / Q
426  assume( getCoeffType (ntCoeffs) == n_Q );
[a355723]427  assume( DEN (f) == NULL );
428
429  number g;
[dc79bd]430//    poly num= p_Copy (NUM (f), ntRing); // ???
431
432
433  // TODO/NOTE: the following should not be necessary (due to
434  // Hannes!) as NUM (f) should be over Z!!!
435  CPolyCoeffsEnumerator itr(NUM(f));
[a355723]436
[e5c9e5]437  n_ClearDenominators(itr, g, ntCoeffs); // may return -1 :(((
[a355723]438
[e5c9e5]439  if( !n_GreaterZero(g, ntCoeffs) )
[dc79bd]440  {
441//     NUM (f) = p_Neg(NUM (f), ntRing); // Ugly :(((
[e5c9e5]442//     g = n_Neg(g, ntCoeffs);
[dc79bd]443    NUM (f) = p_Neg(NUM (f), ntRing); // Ugly :(((
[e5c9e5]444    g = n_Neg(g, ntCoeffs);
[dc79bd]445  }
446
447  // g should be a positive integer now!
[e5c9e5]448  assume( n_GreaterZero(g, ntCoeffs) );
[dc79bd]449
[e5c9e5]450  if( !n_IsOne(g, ntCoeffs) )
[dc79bd]451  {
[e5c9e5]452    assume( n_GreaterZero(g, ntCoeffs) );
453    assume( !n_IsOne(g, ntCoeffs) );
[a355723]454
[dc79bd]455    DEN (f) = p_NSet(g, ntRing); // update COM(f)???
456    assume( DEN (f) != NULL );
457    COM (f) ++;
[a355723]458
[28216f]459    NUM (result)= p_Copy (DEN (f), ntRing);
460  }
[dc79bd]461  else
462  { // common denom == 1?
463    NUM (result)= p_NSet(g, ntRing); // p_Copy (DEN (f), ntRing);
[e5c9e5]464//  n_Delete(&g, ntCoeffs);
[a355723]465  }
466
[dc79bd]467//    if (!p_IsConstant (num, ntRing) && pNext(num) != NULL)
468//    else
469//      g= p_GetAllDenom (num, ntRing);
[e5c9e5]470//    result= (fraction) ntSetMap (ntCoeffs, cf) (g, ntCoeffs, cf);
[dc79bd]471
472  ntTest((number)result);
[2c7f28]473  return (number)result;
[6ccdd3a]474}
475
[2c7f28]476BOOLEAN ntIsOne(number a, const coeffs cf)
[6ccdd3a]477{
[a355723]478  ntTest(a); // !!!
[010f3b]479  definiteGcdCancellation(a, cf, FALSE);
[2c7f28]480  fraction f = (fraction)a;
[31c731]481  return (f!=NULL) && DENIS1(f) && NUMIS1(f);
[6ccdd3a]482}
483
[2c7f28]484BOOLEAN ntIsMOne(number a, const coeffs cf)
[6ccdd3a]485{
[2c7f28]486  ntTest(a);
[010f3b]487  definiteGcdCancellation(a, cf, FALSE);
[2c7f28]488  fraction f = (fraction)a;
[31c731]489  if ((f==NULL) || (!DENIS1(f))) return FALSE;
[e5d267]490  poly g = NUM(f);
[2c7f28]491  if (!p_IsConstant(g, ntRing)) return FALSE;
492  return n_IsMOne(p_GetCoeff(g, ntRing), ntCoeffs);
[6ccdd3a]493}
494
495/// this is in-place, modifies a
[2c7f28]496number ntNeg(number a, const coeffs cf)
[6ccdd3a]497{
[2c7f28]498  ntTest(a);
[e5d267]499  if (!IS0(a))
[2c7f28]500  {
501    fraction f = (fraction)a;
[e5d267]502    NUM(f) = p_Neg(NUM(f), ntRing);
[2c7f28]503  }
[dc79bd]504  ntTest(a);
[6ccdd3a]505  return a;
506}
507
[2c7f28]508number ntImPart(number a, const coeffs cf)
[6ccdd3a]509{
[2c7f28]510  ntTest(a);
[6ccdd3a]511  return NULL;
512}
513
[61b2e16]514number ntInit_bigint(number longratBigIntNumber, const coeffs src, const coeffs cf)
515{
516  assume( cf != NULL );
517
518  const ring A = cf->extRing;
519
520  assume( A != NULL );
521
522  const coeffs C = A->cf;
523
524  assume( C != NULL );
525
526  number n = n_Init_bigint(longratBigIntNumber, src, C);
527
528  if ( n_IsZero(n, C) )
529  {
530    n_Delete(&n, C);
531    return NULL;
[a355723]532  }
[61b2e16]533
534  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[a355723]535
[dc79bd]536  number den = n_GetDenom(n, C);
[a355723]537
[dc79bd]538  assume( n_GreaterZero(den, C) );
[a355723]539
[dc79bd]540  if( n_IsOne(den, C) )
541  {
542     NUM(result) = p_NSet(n, A);
[a355723]543     DEN(result) = NULL;
544     n_Delete(&den, C);
[dc79bd]545  } else
546  {
[a355723]547     DEN(result) = p_NSet(den, A);
548     NUM(result) = p_NSet(n_GetNumerator(n, C), A);
549     n_Delete(&n, C);
[dc79bd]550  }
[61b2e16]551
552  COM(result) = 0;
[dc79bd]553
554  ntTest((number)result);
555
[61b2e16]556  return (number)result;
557}
558
[a355723]559
[2f3764]560number ntInit(long i, const coeffs cf)
[6ccdd3a]561{
[6ac003]562  if (i != 0)
[2c7f28]563  {
[6ac003]564    poly p=p_ISet(i, ntRing);
565    if (p!=NULL)
566    {
567      fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
568      NUM(result) = p;
569      //DEN(result) = NULL; // done by omAlloc0Bin
570      //COM(result) = 0; // done by omAlloc0Bin
571      ntTest((number)result);
572      return (number)result;
573    }
[2c7f28]574  }
[6ac003]575  return NULL;
[6ccdd3a]576}
577
[dc79bd]578
579/// takes over p!
[b38d70]580number ntInit(poly p, const coeffs cf)
581{
[6ac003]582  if (p == NULL) return NULL;
[a355723]583
[e5c9e5]584  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
585
586  if (nCoeff_is_Q(ntCoeffs))
587  {
[dc79bd]588    number g;
589    // TODO/NOTE: the following should not be necessary (due to
[a355723]590    // Hannes!) as NUM (f) should be over Z!!!
[e5c9e5]591    // but it is not: normalize it
[dc79bd]592    CPolyCoeffsEnumerator itr(p);
593
[e5c9e5]594    n_ClearDenominators(itr, g, ntCoeffs);
[dc79bd]595
[e5c9e5]596    if( !n_GreaterZero(g, ntCoeffs) )
[dc79bd]597    {
[e5c9e5]598      p = p_Neg(p, ntRing);
599      g = n_Neg(g, ntCoeffs);
[dc79bd]600    }
601
602    // g should be a positive integer now!
[e5c9e5]603    assume( n_GreaterZero(g, ntCoeffs) );
[a355723]604
[e5c9e5]605    if( !n_IsOne(g, ntCoeffs) )
[dc79bd]606    {
[a355723]607      DEN (f) = p_NSet(g, ntRing);
608      assume( DEN (f) != NULL );
[dc79bd]609    }
610    else
611    {
[e5c9e5]612      //DEN(f) = NULL; // done by omAlloc0
613      n_Delete(&g, ntCoeffs);
[a355723]614    }
[e5c9e5]615  }
[a355723]616
[e5c9e5]617  NUM(f) = p;
618  COM(f) = 0;
[dc79bd]619
[e5c9e5]620  ntTest((number)f);
621  return (number)f;
[b38d70]622}
623
[2c7f28]624int ntInt(number &a, const coeffs cf)
[6ccdd3a]625{
[2c7f28]626  ntTest(a);
[e5d267]627  if (IS0(a)) return 0;
[010f3b]628  definiteGcdCancellation(a, cf, FALSE);
[2c7f28]629  fraction f = (fraction)a;
[e5d267]630  if (!DENIS1(f)) return 0;
[a355723]631
[1090a98]632  const poly aAsPoly = NUM(f);
633
634  if(aAsPoly == NULL)
635    return 0;
[a355723]636
[1090a98]637  if (!p_IsConstant(aAsPoly, ntRing))
638    return 0;
639
640  assume( aAsPoly != NULL );
[a355723]641
[1090a98]642  return n_Int(p_GetCoeff(aAsPoly, ntRing), ntCoeffs);
[2c7f28]643}
644
645/* This method will only consider the numerators of a and b, without
646   cancelling gcd's before.
647   Moreover it may return TRUE only if one or both numerators
648   are zero or if their degrees are equal. Then TRUE is returned iff
649   coeff(numerator(a)) > coeff(numerator(b));
650   In all other cases, FALSE will be returned. */
651BOOLEAN ntGreater(number a, number b, const coeffs cf)
652{
[dc79bd]653  ntTest(a);
654  ntTest(b);
[c5fbc4]655  number aNumCoeff = NULL; int aNumDeg = 0;
656  number aDenCoeff = NULL; int aDenDeg = 0;
657  number bNumCoeff = NULL; int bNumDeg = 0;
658  number bDenCoeff = NULL; int bDenDeg = 0;
[e5d267]659  if (!IS0(a))
[2c7f28]660  {
661    fraction fa = (fraction)a;
[e5d267]662    aNumDeg = p_Totaldegree(NUM(fa), ntRing);
663    aNumCoeff = p_GetCoeff(NUM(fa), ntRing);
[fd60be]664    if (DEN(fa)!=NULL)
665    {
666      aDenDeg = p_Totaldegree(DEN(fa), ntRing);
667      aDenCoeff=p_GetCoeff(DEN(fa),ntRing);
668    }
[2c7f28]669  }
[95c826]670  else return !(ntGreaterZero (b,cf));
[e5d267]671  if (!IS0(b))
[2c7f28]672  {
673    fraction fb = (fraction)b;
[e5d267]674    bNumDeg = p_Totaldegree(NUM(fb), ntRing);
675    bNumCoeff = p_GetCoeff(NUM(fb), ntRing);
[fd60be]676    if (DEN(fb)!=NULL)
677    {
678      bDenDeg = p_Totaldegree(DEN(fb), ntRing);
679      bDenCoeff=p_GetCoeff(DEN(fb),ntRing);
680    }
[2c7f28]681  }
[95c826]682  else return ntGreaterZero(a,cf);
[fd60be]683  if (aNumDeg-aDenDeg > bNumDeg-bDenDeg) return TRUE;
684  if (aNumDeg-aDenDeg < bNumDeg-bDenDeg) return FALSE;
685  number aa;
686  number bb;
[e5c9e5]687  if (bDenCoeff==NULL) aa=n_Copy(aNumCoeff,ntCoeffs);
688  else                 aa=n_Mult(aNumCoeff,bDenCoeff,ntCoeffs);
689  if (aDenCoeff==NULL) bb=n_Copy(bNumCoeff,ntCoeffs);
690  else                 bb=n_Mult(bNumCoeff,aDenCoeff,ntCoeffs);
[fd60be]691  BOOLEAN rr= n_Greater(aa, bb, ntCoeffs);
[e5c9e5]692  n_Delete(&aa,ntCoeffs);
693  n_Delete(&bb,ntCoeffs);
[fd60be]694  return rr;
[6ccdd3a]695}
696
[2c7f28]697/* this method will only consider the numerator of a, without cancelling
698   the gcd before;
699   returns TRUE iff the leading coefficient of the numerator of a is > 0
700                    or the leading term of the numerator of a is not a
701                    constant */
702BOOLEAN ntGreaterZero(number a, const coeffs cf)
[6ccdd3a]703{
[2c7f28]704  ntTest(a);
[e5d267]705  if (IS0(a)) return FALSE;
[2c7f28]706  fraction f = (fraction)a;
[e5d267]707  poly g = NUM(f);
[7964f9c]708  return (!p_LmIsConstant(g,ntRing)|| n_GreaterZero(pGetCoeff(g), ntCoeffs));
[6ccdd3a]709}
710
[03f7b5]711void ntCoeffWrite(const coeffs cf, BOOLEAN details)
[6ccdd3a]712{
[a55ef0]713  assume( cf != NULL );
714
715  const ring A = cf->extRing;
716
717  assume( A != NULL );
718  assume( A->cf != NULL );
719
[03f7b5]720  n_CoeffWrite(A->cf, details);
[a55ef0]721
722//  rWrite(A);
723
724  const int P = rVar(A);
725  assume( P > 0 );
[c14846c]726
[a55ef0]727  Print("//   %d parameter    : ", P);
[c14846c]728
[a55ef0]729  for (int nop=0; nop < P; nop ++)
730    Print("%s ", rRingVar(nop, A));
731
[dd668f]732  assume( A->qideal == NULL );
[c14846c]733
[a55ef0]734  PrintS("\n//   minpoly        : 0\n");
735
736/*
[2c7f28]737  PrintS("//   Coefficients live in the rational function field\n");
738  Print("//   K(");
739  for (int i = 0; i < rVar(ntRing); i++)
740  {
[a55ef0]741    if (i > 0) PrintS(" ");
[2c7f28]742    Print("%s", rRingVar(i, ntRing));
743  }
744  PrintS(") with\n");
745  PrintS("//   K: "); n_CoeffWrite(cf->extRing->cf);
[a55ef0]746*/
[6ccdd3a]747}
748
[12f6ce0]749number ntDiff(number a, number d, const coeffs cf)
750{
751  ntTest(a);
752  ntTest(d);
753
[0b8402]754  if (IS0(d)) 
755  {
756    WerrorS("ringvar expected");
757    return NULL;
758  }
[96847c]759  fraction t = (fraction) d;
[12f6ce0]760  if (!DENIS1(t))
761  {
762    WerrorS("expected differentiation by a variable");
[96847c]763    return NULL;
[12f6ce0]764  }
765  int k=p_Var(NUM(t),ntRing);
766  if (k==0)
767  {
768    WerrorS("expected differentiation by a variable");
[96847c]769    return NULL;
[12f6ce0]770  }
771
772  if (IS0(a)) return ntCopy(a, cf);
773
774  fraction fa = (fraction)a;
775  if (DENIS1(fa)) {
776
777     fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[1e8006]778     NUM(result) = p_Diff(NUM(fa),k,ntRing);
[12f6ce0]779     DEN(result) = NULL;
780     COM(result) = COM(fa);
781     return (number)result;
782  }
783
784  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[1e8006]785  poly fg = p_Mult_q(p_Copy(DEN(fa),ntRing),p_Diff(NUM(fa),k,ntRing),ntRing);
786  poly gf = p_Mult_q(p_Copy(NUM(fa),ntRing),p_Diff(DEN(fa),k,ntRing),ntRing);
787  NUM(result) = p_Sub(fg,gf,ntRing);
[5d258e9]788  if (NUM(result)==NULL) return(NULL);
[1e8006]789  DEN(result) = pp_Mult_qq(DEN(fa), DEN(fa), ntRing);
790  COM(result) = COM(fa) + COM(fa) + DIFF_COMPLEXITY;
[12f6ce0]791  heuristicGcdCancellation((number)result, cf);
792
793  return (number)result;
794}
795
796
[2c7f28]797number ntAdd(number a, number b, const coeffs cf)
[6ccdd3a]798{
[dc79bd]799  ntTest(a);
800  ntTest(b);
[e5d267]801  if (IS0(a)) return ntCopy(b, cf);
802  if (IS0(b)) return ntCopy(a, cf);
[c14846c]803
[2c7f28]804  fraction fa = (fraction)a;
805  fraction fb = (fraction)b;
[c14846c]806
[e5d267]807  poly g = p_Copy(NUM(fa), ntRing);
808  if (!DENIS1(fb)) g = p_Mult_q(g, p_Copy(DEN(fb), ntRing), ntRing);
809  poly h = p_Copy(NUM(fb), ntRing);
810  if (!DENIS1(fa)) h = p_Mult_q(h, p_Copy(DEN(fa), ntRing), ntRing);
[de90c01]811  g = p_Add_q(g, h, ntRing);
[c14846c]812
[de90c01]813  if (g == NULL) return NULL;
[c14846c]814
[2c7f28]815  poly f;
[e5d267]816  if      (DENIS1(fa) && DENIS1(fb))  f = NULL;
817  else if (!DENIS1(fa) && DENIS1(fb)) f = p_Copy(DEN(fa), ntRing);
818  else if (DENIS1(fa) && !DENIS1(fb)) f = p_Copy(DEN(fb), ntRing);
819  else /* both denom's are != 1 */    f = p_Mult_q(p_Copy(DEN(fa), ntRing),
820                                                   p_Copy(DEN(fb), ntRing),
[2c7f28]821                                                   ntRing);
[c14846c]822
[2c7f28]823  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[e5d267]824  NUM(result) = g;
825  DEN(result) = f;
826  COM(result) = COM(fa) + COM(fb) + ADD_COMPLEXITY;
[2c7f28]827  heuristicGcdCancellation((number)result, cf);
[dc79bd]828
829//  ntTest((number)result);
[a355723]830
[2c7f28]831  return (number)result;
[6ccdd3a]832}
833
[2c7f28]834number ntSub(number a, number b, const coeffs cf)
[6ccdd3a]835{
[dc79bd]836  ntTest(a);
837  ntTest(b);
[e5d267]838  if (IS0(a)) return ntNeg(ntCopy(b, cf), cf);
839  if (IS0(b)) return ntCopy(a, cf);
[c14846c]840
[2c7f28]841  fraction fa = (fraction)a;
842  fraction fb = (fraction)b;
[c14846c]843
[e5d267]844  poly g = p_Copy(NUM(fa), ntRing);
845  if (!DENIS1(fb)) g = p_Mult_q(g, p_Copy(DEN(fb), ntRing), ntRing);
846  poly h = p_Copy(NUM(fb), ntRing);
847  if (!DENIS1(fa)) h = p_Mult_q(h, p_Copy(DEN(fa), ntRing), ntRing);
[de90c01]848  g = p_Add_q(g, p_Neg(h, ntRing), ntRing);
[c14846c]849
[de90c01]850  if (g == NULL) return NULL;
[c14846c]851
[2c7f28]852  poly f;
[e5d267]853  if      (DENIS1(fa) && DENIS1(fb))  f = NULL;
854  else if (!DENIS1(fa) && DENIS1(fb)) f = p_Copy(DEN(fa), ntRing);
855  else if (DENIS1(fa) && !DENIS1(fb)) f = p_Copy(DEN(fb), ntRing);
856  else /* both den's are != 1 */      f = p_Mult_q(p_Copy(DEN(fa), ntRing),
857                                                   p_Copy(DEN(fb), ntRing),
[2c7f28]858                                                   ntRing);
[c14846c]859
[2c7f28]860  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[e5d267]861  NUM(result) = g;
862  DEN(result) = f;
863  COM(result) = COM(fa) + COM(fb) + ADD_COMPLEXITY;
[2c7f28]864  heuristicGcdCancellation((number)result, cf);
[dc79bd]865//  ntTest((number)result);
[2c7f28]866  return (number)result;
[6ccdd3a]867}
868
[2c7f28]869number ntMult(number a, number b, const coeffs cf)
[6ccdd3a]870{
[dc79bd]871  ntTest(a); // !!!?
872  ntTest(b); // !!!?
[a355723]873
[e5d267]874  if (IS0(a) || IS0(b)) return NULL;
[c14846c]875
[2c7f28]876  fraction fa = (fraction)a;
877  fraction fb = (fraction)b;
[c14846c]878
[dc79bd]879  const poly g = pp_Mult_qq(NUM(fa), NUM(fb), ntRing);
[c14846c]880
[dc79bd]881  if (g == NULL) return NULL; // may happen due to zero divisors???
[c14846c]882
[2c7f28]883  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[dc79bd]884
[e5d267]885  NUM(result) = g;
[dc79bd]886
887  const poly da = DEN(fa);
888  const poly db = DEN(fb);
889
[a355723]890
[dc79bd]891  if (db == NULL)
892  {
893    // b = ? // NULL
[a355723]894
[dc79bd]895    if(da == NULL)
896    { // both fa && fb are ?? // NULL!
897      assume (da == NULL && db == NULL);
898      DEN(result) = NULL;
899      COM(result) = 0;
900    }
[a355723]901    else
[dc79bd]902    {
903      assume (da != NULL && db == NULL);
904      DEN(result) = p_Copy(da, ntRing);
905      COM(result) = COM(fa) + MULT_COMPLEXITY;
906      heuristicGcdCancellation((number)result, cf);
907    }
908  } else
909  { // b = ?? / ??
910    if (da == NULL)
911    { // a == ? // NULL
912      assume( db != NULL && da == NULL);
913      DEN(result) = p_Copy(db, ntRing);
914      COM(result) = COM(fb) + MULT_COMPLEXITY;
915      heuristicGcdCancellation((number)result, cf);
916    }
917    else /* both den's are != 1 */
918    {
919      assume (da != NULL && db != NULL);
920      DEN(result) = pp_Mult_qq(da, db, ntRing);
921      COM(result) = COM(fa) + COM(fb) + MULT_COMPLEXITY;
922      heuristicGcdCancellation((number)result, cf);
923    }
924  }
925
926//  ntTest((number)result);
[a355723]927
[2c7f28]928  return (number)result;
[6ccdd3a]929}
930
[2c7f28]931number ntDiv(number a, number b, const coeffs cf)
[6ccdd3a]932{
[dc79bd]933  ntTest(a);
934  ntTest(b);
[e5d267]935  if (IS0(a)) return NULL;
936  if (IS0(b)) WerrorS(nDivBy0);
[c14846c]937
[2c7f28]938  fraction fa = (fraction)a;
939  fraction fb = (fraction)b;
[c14846c]940
[e5d267]941  poly g = p_Copy(NUM(fa), ntRing);
942  if (!DENIS1(fb)) g = p_Mult_q(g, p_Copy(DEN(fb), ntRing), ntRing);
[c14846c]943
[de90c01]944  if (g == NULL) return NULL;   /* may happen due to zero divisors */
[c14846c]945
[e5d267]946  poly f = p_Copy(NUM(fb), ntRing);
947  if (!DENIS1(fa)) f = p_Mult_q(f, p_Copy(DEN(fa), ntRing), ntRing);
[c14846c]948
[2c7f28]949  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[e5d267]950  NUM(result) = g;
[e5c9e5]951  if (!n_GreaterZero(pGetCoeff(f),ntCoeffs))
[644b31]952  {
953    g=p_Neg(g,ntRing);
954    f=p_Neg(f,ntRing);
955    NUM(result) = g;
956  }
[e5c9e5]957  if (!p_IsConstant(f,ntRing) || !n_IsOne(pGetCoeff(f),ntCoeffs))
[644b31]958  {
[1374bc]959    DEN(result) = f;
[644b31]960  }
[e5d267]961  COM(result) = COM(fa) + COM(fb) + MULT_COMPLEXITY;
[2c7f28]962  heuristicGcdCancellation((number)result, cf);
[dc79bd]963//  ntTest((number)result);
[2c7f28]964  return (number)result;
[6ccdd3a]965}
966
967/* 0^0 = 0;
968   for |exp| <= 7 compute power by a simple multiplication loop;
969   for |exp| >= 8 compute power along binary presentation of |exp|, e.g.
970      p^13 = p^1 * p^4 * p^8, where we utilise that
[2c7f28]971      p^(2^(k+1)) = p^(2^k) * p^(2^k);
972   intermediate cancellation is controlled by the in-place method
973   heuristicGcdCancellation; see there.
[6ccdd3a]974*/
[2c7f28]975void ntPower(number a, int exp, number *b, const coeffs cf)
[6ccdd3a]976{
[2c7f28]977  ntTest(a);
[c14846c]978
[6ccdd3a]979  /* special cases first */
[e5d267]980  if (IS0(a))
[6ccdd3a]981  {
982    if (exp >= 0) *b = NULL;
983    else          WerrorS(nDivBy0);
984  }
[35e86e]985  else if (exp ==  0) { *b = ntInit(1, cf); return;}
986  else if (exp ==  1) { *b = ntCopy(a, cf); return;}
987  else if (exp == -1) { *b = ntInvers(a, cf); return;}
[c14846c]988
[6ccdd3a]989  int expAbs = exp; if (expAbs < 0) expAbs = -expAbs;
[c14846c]990
[2c7f28]991  /* now compute a^expAbs */
992  number pow; number t;
[6ccdd3a]993  if (expAbs <= 7)
994  {
[2c7f28]995    pow = ntCopy(a, cf);
[6ccdd3a]996    for (int i = 2; i <= expAbs; i++)
997    {
[2c7f28]998      t = ntMult(pow, a, cf);
999      ntDelete(&pow, cf);
1000      pow = t;
1001      heuristicGcdCancellation(pow, cf);
[6ccdd3a]1002    }
1003  }
1004  else
1005  {
[2c7f28]1006    pow = ntInit(1, cf);
1007    number factor = ntCopy(a, cf);
[6ccdd3a]1008    while (expAbs != 0)
1009    {
1010      if (expAbs & 1)
1011      {
[2c7f28]1012        t = ntMult(pow, factor, cf);
1013        ntDelete(&pow, cf);
1014        pow = t;
1015        heuristicGcdCancellation(pow, cf);
[6ccdd3a]1016      }
1017      expAbs = expAbs / 2;
1018      if (expAbs != 0)
1019      {
[2c7f28]1020        t = ntMult(factor, factor, cf);
1021        ntDelete(&factor, cf);
1022        factor = t;
1023        heuristicGcdCancellation(factor, cf);
[6ccdd3a]1024      }
1025    }
[2c7f28]1026    ntDelete(&factor, cf);
[6ccdd3a]1027  }
[c14846c]1028
[6ccdd3a]1029  /* invert if original exponent was negative */
1030  if (exp < 0)
1031  {
[2c7f28]1032    t = ntInvers(pow, cf);
1033    ntDelete(&pow, cf);
1034    pow = t;
[6ccdd3a]1035  }
[2c7f28]1036  *b = pow;
[dc79bd]1037  ntTest(*b);
[6ccdd3a]1038}
1039
[06df101]1040/* assumes that cf represents the rationals, i.e. Q, and will only
1041   be called in that case;
1042   assumes furthermore that f != NULL and that the denominator of f != 1;
1043   generally speaking, this method removes denominators in the rational
1044   coefficients of the numerator and denominator of 'a';
1045   more concretely, the following normalizations will be performed,
1046   where t^alpha denotes a monomial in the transcendental variables t_k
1047   (1) if 'a' is of the form
1048          (sum_alpha a_alpha/b_alpha * t^alpha)
1049          -------------------------------------
1050            (sum_beta c_beta/d_beta * t^beta)
1051       with integers a_alpha, b_alpha, c_beta, d_beta, then both the
1052       numerator and the denominator will be multiplied by the LCM of
1053       the b_alpha's and the d_beta's (if this LCM is != 1),
1054   (2) if 'a' is - e.g. after having performed step (1) - of the form
1055          (sum_alpha a_alpha * t^alpha)
1056          -----------------------------
1057            (sum_beta c_beta * t^beta)
1058       with integers a_alpha, c_beta, and with a non-constant denominator,
1059       then both the numerator and the denominator will be divided by the
1060       GCD of the a_alpha's and the c_beta's (if this GCD is != 1),
1061   this procedure does not alter COM(f) (this has to be done by the
1062   calling procedure);
1063   modifies f */
1064void handleNestedFractionsOverQ(fraction f, const coeffs cf)
1065{
1066  assume(nCoeff_is_Q(ntCoeffs));
1067  assume(!IS0(f));
1068  assume(!DENIS1(f));
[c14846c]1069
[06df101]1070  { /* step (1); see documentation of this procedure above */
[d12f186]1071    p_Normalize(NUM(f), ntRing);
1072    p_Normalize(DEN(f), ntRing);
[06df101]1073    number lcmOfDenominators = n_Init(1, ntCoeffs);
1074    number c; number tmp;
1075    poly p = NUM(f);
1076    /* careful when using n_Lcm!!! It computes the lcm of the numerator
1077       of the 1st argument and the denominator of the 2nd!!! */
1078    while (p != NULL)
1079    {
1080      c = p_GetCoeff(p, ntRing);
1081      tmp = n_Lcm(lcmOfDenominators, c, ntCoeffs);
1082      n_Delete(&lcmOfDenominators, ntCoeffs);
1083      lcmOfDenominators = tmp;
1084      pIter(p);
1085    }
1086    p = DEN(f);
1087    while (p != NULL)
1088    {
1089      c = p_GetCoeff(p, ntRing);
1090      tmp = n_Lcm(lcmOfDenominators, c, ntCoeffs);
1091      n_Delete(&lcmOfDenominators, ntCoeffs);
1092      lcmOfDenominators = tmp;
1093      pIter(p);
1094    }
1095    if (!n_IsOne(lcmOfDenominators, ntCoeffs))
1096    { /* multiply NUM(f) and DEN(f) with lcmOfDenominators */
1097      NUM(f) = p_Mult_nn(NUM(f), lcmOfDenominators, ntRing);
[d12f186]1098      p_Normalize(NUM(f), ntRing);
[06df101]1099      DEN(f) = p_Mult_nn(DEN(f), lcmOfDenominators, ntRing);
[d12f186]1100      p_Normalize(DEN(f), ntRing);
[06df101]1101    }
1102    n_Delete(&lcmOfDenominators, ntCoeffs);
[2e7132]1103    if (DEN(f)!=NULL)
[06df101]1104    { /* step (2); see documentation of this procedure above */
1105      p = NUM(f);
1106      number gcdOfCoefficients = n_Copy(p_GetCoeff(p, ntRing), ntCoeffs);
1107      pIter(p);
1108      while ((p != NULL) && (!n_IsOne(gcdOfCoefficients, ntCoeffs)))
1109      {
1110        c = p_GetCoeff(p, ntRing);
1111        tmp = n_Gcd(c, gcdOfCoefficients, ntCoeffs);
1112        n_Delete(&gcdOfCoefficients, ntCoeffs);
1113        gcdOfCoefficients = tmp;
1114        pIter(p);
1115      }
1116      p = DEN(f);
1117      while ((p != NULL) && (!n_IsOne(gcdOfCoefficients, ntCoeffs)))
1118      {
1119        c = p_GetCoeff(p, ntRing);
1120        tmp = n_Gcd(c, gcdOfCoefficients, ntCoeffs);
1121        n_Delete(&gcdOfCoefficients, ntCoeffs);
1122        gcdOfCoefficients = tmp;
1123        pIter(p);
1124      }
1125      if (!n_IsOne(gcdOfCoefficients, ntCoeffs))
1126      { /* divide NUM(f) and DEN(f) by gcdOfCoefficients */
1127        number inverseOfGcdOfCoefficients = n_Invers(gcdOfCoefficients,
1128                                                     ntCoeffs);
1129        NUM(f) = p_Mult_nn(NUM(f), inverseOfGcdOfCoefficients, ntRing);
[d12f186]1130        p_Normalize(NUM(f), ntRing);
[06df101]1131        DEN(f) = p_Mult_nn(DEN(f), inverseOfGcdOfCoefficients, ntRing);
[d12f186]1132        p_Normalize(DEN(f), ntRing);
[06df101]1133        n_Delete(&inverseOfGcdOfCoefficients, ntCoeffs);
1134      }
1135      n_Delete(&gcdOfCoefficients, ntCoeffs);
1136    }
1137  }
[c14846c]1138
[06df101]1139  /* Now, due to the above computations, DEN(f) may have become the
1140     1-polynomial which needs to be represented by NULL: */
1141  if ((DEN(f) != NULL) &&
1142      p_IsConstant(DEN(f), ntRing) &&
1143      n_IsOne(p_GetCoeff(DEN(f), ntRing), ntCoeffs))
1144  {
1145    p_Delete(&DEN(f), ntRing); DEN(f) = NULL;
1146  }
[a355723]1147
[dc79bd]1148  if( DEN(f) != NULL )
1149    if( !n_GreaterZero(pGetCoeff(DEN(f)), ntCoeffs) )
1150    {
1151      NUM(f) = p_Neg(NUM(f), ntRing);
1152      DEN(f) = p_Neg(DEN(f), ntRing);
[a355723]1153    }
1154
[dc79bd]1155  ntTest((number)f); // TODO!
[06df101]1156}
1157
[2c7f28]1158/* modifies a */
1159void heuristicGcdCancellation(number a, const coeffs cf)
[6ccdd3a]1160{
[dc79bd]1161//  ntTest(a); // !!!!????
[e5d267]1162  if (IS0(a)) return;
[c14846c]1163
[2c7f28]1164  fraction f = (fraction)a;
[5a4e17]1165  if (COM(f)!=0) p_Normalize(NUM(f), ntRing);
1166  if (DENIS1(f) || NUMIS1(f)) { COM(f) = 0; return; }
1167
1168  p_Normalize(DEN(f), ntRing);
[dc79bd]1169
1170  assume( DEN(f) != NULL );
[c14846c]1171
[e5d267]1172  /* check whether NUM(f) = DEN(f), and - if so - replace 'a' by 1 */
1173  if (p_EqualPolys(NUM(f), DEN(f), ntRing))
1174  { /* numerator and denominator are both != 1 */
1175    p_Delete(&NUM(f), ntRing); NUM(f) = p_ISet(1, ntRing);
1176    p_Delete(&DEN(f), ntRing); DEN(f) = NULL;
1177    COM(f) = 0;
[603ac7a]1178  }
1179  else
[dc79bd]1180  {
[a355723]1181    if (COM(f) > BOUND_COMPLEXITY)
[dc79bd]1182      definiteGcdCancellation(a, cf, TRUE);
[c14846c]1183
[7964f9c]1184    // TODO: check if it is enough to put the following into definiteGcdCancellation?!
1185    if( DEN(f) != NULL )
[dc79bd]1186    {
[7964f9c]1187      if( !n_GreaterZero(pGetCoeff(DEN(f)), ntCoeffs) )
1188      {
1189        NUM(f) = p_Neg(NUM(f), ntRing);
1190        DEN(f) = p_Neg(DEN(f), ntRing);
1191      }
1192      if (ntCoeffs->has_simple_Inverse)
1193      {
1194        if (!n_IsOne(pGetCoeff(DEN(f)),ntCoeffs))
1195        {
1196          number inv=n_Invers(pGetCoeff(DEN(f)),ntCoeffs);
1197          DEN(f)=p_Mult_nn(DEN(f),inv,ntRing);
1198          NUM(f)=p_Mult_nn(NUM(f),inv,ntRing);
1199        }
1200        if(p_LmIsConstant(DEN(f),ntRing))
1201        {
1202          p_Delete(&DEN(f),ntRing);
1203          COM(f)=0;
1204        }
1205      }
[a355723]1206    }
[dc79bd]1207  }
[a355723]1208
[603ac7a]1209  ntTest(a);
[6ccdd3a]1210}
1211
[dc79bd]1212/// modifies a
[010f3b]1213void definiteGcdCancellation(number a, const coeffs cf,
[06df101]1214                             BOOLEAN simpleTestsHaveAlreadyBeenPerformed)
[6ccdd3a]1215{
[dc79bd]1216  ntTest(a); // !!!!
[c14846c]1217
[2c7f28]1218  fraction f = (fraction)a;
[c14846c]1219
[603ac7a]1220  if (IS0(a)) return;
[14a0ca]1221  if (NUM(f)!=NULL) p_Normalize(NUM(f), ntRing);
1222  if (DEN(f)!=NULL) p_Normalize(DEN(f), ntRing);
[06df101]1223  if (!simpleTestsHaveAlreadyBeenPerformed)
[2c7f28]1224  {
[14a0ca]1225    //p_Normalize(NUM(f), ntRing);
1226    //if (DEN(f)!=NULL) p_Normalize(DEN(f), ntRing);
[e5d267]1227    if (DENIS1(f) || NUMIS1(f)) { COM(f) = 0; return; }
[c14846c]1228
[e5d267]1229    /* check whether NUM(f) = DEN(f), and - if so - replace 'a' by 1 */
1230    if (p_EqualPolys(NUM(f), DEN(f), ntRing))
1231    { /* numerator and denominator are both != 1 */
1232      p_Delete(&NUM(f), ntRing); NUM(f) = p_ISet(1, ntRing);
1233      p_Delete(&DEN(f), ntRing); DEN(f) = NULL;
1234      COM(f) = 0;
[dc79bd]1235      ntTest(a); // !!!!
[010f3b]1236      return;
1237    }
[2c7f28]1238  }
[14a0ca]1239  /*if (rField_is_Q(ntRing))
[5a4e17]1240  {
1241    number c=n_Copy(pGetCoeff(NUM(f)),ntCoeffs);
1242    poly p=pNext(NUM(f));
1243    while((p!=NULL)&&(!n_IsOne(c,ntCoeffs)))
1244    {
1245      number cc=n_Gcd(c,pGetCoeff(p),ntCoeffs);
1246      n_Delete(&c,ntCoeffs);
1247      c=cc;
1248      pIter(p);
1249    };
1250    p=DEN(f);
1251    while((p!=NULL)&&(!n_IsOne(c,ntCoeffs)))
1252    {
1253      number cc=n_Gcd(c,pGetCoeff(p),ntCoeffs);
1254      n_Delete(&c,ntCoeffs);
1255      c=cc;
1256      pIter(p);
1257    };
1258    if(!n_IsOne(c,ntCoeffs))
1259    {
1260      p=NUM(f);
1261      do
1262      {
1263        number cc=n_Div(pGetCoeff(p),c,ntCoeffs);
1264        n_Normalize(cc,ntCoeffs);
1265        p_SetCoeff(p,cc,ntRing);
1266        pIter(p);
1267      } while(p!=NULL);
1268      p=DEN(f);
1269      do
1270      {
1271        number cc=n_Div(pGetCoeff(p),c,ntCoeffs);
1272        n_Normalize(cc,ntCoeffs);
1273        p_SetCoeff(p,cc,ntRing);
1274        pIter(p);
1275      } while(p!=NULL);
1276      n_Delete(&c,ntCoeffs);
1277      if(pNext(DEN(f))==NULL)
1278      {
1279        if (p_IsOne(DEN(f),ntRing))
1280        {
1281          p_LmDelete(&DEN(f),ntRing);
1282          COM(f)=0;
1283          return;
1284        }
1285        else
1286        {
1287          return;
1288        }
1289      }
1290    }
[14a0ca]1291  }*/
[c14846c]1292
1293#ifdef HAVE_FACTORY
[5a4e17]1294  poly pGcd;
1295  /* here we assume: NUM(f), DEN(f) !=NULL, in Z_a reqp. Z/p_a */
1296    pGcd = singclap_gcd_r(NUM(f), DEN(f), ntRing);
1297  if (p_IsConstant(pGcd, ntRing)
[14a0ca]1298  && n_IsOne(p_GetCoeff(pGcd, ntRing), ntCoeffs)
[5a4e17]1299  )
[06df101]1300  { /* gcd = 1; nothing to cancel;
1301       Suppose the given rational function field is over Q. Although the
1302       gcd is 1, we may have produced fractional coefficients in NUM(f),
1303       DEN(f), or both, due to previous arithmetics. The next call will
1304       remove those nested fractions, in case there are any. */
[9053f02]1305    if (nCoeff_is_Zp(ntCoeffs) && p_IsConstant (DEN (f), ntRing))
1306    {
1307      NUM (f) = p_Div_nn (NUM (f), p_GetCoeff (DEN(f),ntRing), ntRing);
1308      //poly newNum= singclap_pdivide (NUM(f), DEN (f), ntRing);
1309      //p_Delete(&NUM (f), ntRing);
1310      //NUM (f)= newNum;
1311      p_Delete(&DEN (f), ntRing);
1312      DEN (f) = NULL;
1313      COM (f) = 0;
1314    } else if (nCoeff_is_Q(ntCoeffs)) handleNestedFractionsOverQ(f, cf);
[e5d267]1315  }
1316  else
[06df101]1317  { /* We divide both NUM(f) and DEN(f) by the gcd which is known
1318       to be != 1. */
1319    poly newNum = singclap_pdivide(NUM(f), pGcd, ntRing);
1320    p_Delete(&NUM(f), ntRing);
1321    NUM(f) = newNum;
1322    poly newDen = singclap_pdivide(DEN(f), pGcd, ntRing);
1323    p_Delete(&DEN(f), ntRing);
1324    DEN(f) = newDen;
1325    if (p_IsConstant(DEN(f), ntRing) &&
1326        n_IsOne(p_GetCoeff(DEN(f), ntRing), ntCoeffs))
1327    {
1328      /* DEN(f) = 1 needs to be represented by NULL! */
[e5d267]1329      p_Delete(&DEN(f), ntRing);
[06df101]1330      newDen = NULL;
1331    }
1332    else
1333    { /* Note that over Q, by cancelling the gcd, we may have produced
1334         fractional coefficients in NUM(f), DEN(f), or both. The next
1335         call will remove those nested fractions, in case there are
1336         any. */
1337      if (nCoeff_is_Q(ntCoeffs)) handleNestedFractionsOverQ(f, cf);
1338    }
[e5d267]1339  }
1340  COM(f) = 0;
1341  p_Delete(&pGcd, ntRing);
[dc79bd]1342
1343  if( DEN(f) != NULL )
1344    if( !n_GreaterZero(pGetCoeff(DEN(f)), ntCoeffs) )
1345    {
1346      NUM(f) = p_Neg(NUM(f), ntRing);
1347      DEN(f) = p_Neg(DEN(f), ntRing);
[a355723]1348    }
[e5d267]1349#endif /* HAVE_FACTORY */
[a355723]1350
[dc79bd]1351  ntTest(a); // !!!!
[2c7f28]1352}
1353
[ce1f78]1354// NOTE: modifies a
1355void ntWriteLong(number &a, const coeffs cf)
1356{
1357  ntTest(a);
1358  definiteGcdCancellation(a, cf, FALSE);
1359  if (IS0(a))
1360    StringAppendS("0");
1361  else
1362  {
1363    fraction f = (fraction)a;
1364    // stole logic from napWrite from kernel/longtrans.cc of legacy singular
1365    BOOLEAN omitBrackets = p_IsConstant(NUM(f), ntRing);
1366    if (!omitBrackets) StringAppendS("(");
1367    p_String0Long(NUM(f), ntRing, ntRing);
1368    if (!omitBrackets) StringAppendS(")");
1369    if (!DENIS1(f))
1370    {
1371      StringAppendS("/");
1372      omitBrackets = p_IsConstant(DEN(f), ntRing);
1373      if (!omitBrackets) StringAppendS("(");
1374      p_String0Long(DEN(f), ntRing, ntRing);
1375      if (!omitBrackets) StringAppendS(")");
[a355723]1376    }
[ce1f78]1377  }
[dc79bd]1378  ntTest(a); // !!!!
[ce1f78]1379}
1380
1381// NOTE: modifies a
1382void ntWriteShort(number &a, const coeffs cf)
[2c7f28]1383{
1384  ntTest(a);
[010f3b]1385  definiteGcdCancellation(a, cf, FALSE);
[e5d267]1386  if (IS0(a))
[6ccdd3a]1387    StringAppendS("0");
1388  else
1389  {
[2c7f28]1390    fraction f = (fraction)a;
[a5071b9]1391    // stole logic from napWrite from kernel/longtrans.cc of legacy singular
[a0a9f0]1392    BOOLEAN omitBrackets = p_IsConstant(NUM(f), ntRing);
[a5071b9]1393    if (!omitBrackets) StringAppendS("(");
[ce1f78]1394    p_String0Short(NUM(f), ntRing, ntRing);
[a5071b9]1395    if (!omitBrackets) StringAppendS(")");
[e5d267]1396    if (!DENIS1(f))
[2c7f28]1397    {
1398      StringAppendS("/");
[a0a9f0]1399      omitBrackets = p_IsConstant(DEN(f), ntRing);
1400      if (!omitBrackets) StringAppendS("(");
[ce1f78]1401      p_String0Short(DEN(f), ntRing, ntRing);
[a0a9f0]1402      if (!omitBrackets) StringAppendS(")");
[2c7f28]1403    }
[6ccdd3a]1404  }
[dc79bd]1405  ntTest(a);
[6ccdd3a]1406}
1407
[2c7f28]1408const char * ntRead(const char *s, number *a, const coeffs cf)
[6ccdd3a]1409{
[2c7f28]1410  poly p;
1411  const char * result = p_Read(s, p, ntRing);
[dc79bd]1412  if (p == NULL) *a = NULL;
1413  else *a = ntInit(p, cf);
1414  return result;
[6ccdd3a]1415}
1416
[237b4dd]1417void ntNormalize (number &a, const coeffs cf)
1418{
[e56eb1a]1419  if ((a!=NULL))
[6ad7e6]1420  {
1421    definiteGcdCancellation(a, cf, FALSE);
1422    if ((DEN(a)!=NULL)
[e5c9e5]1423    &&(!n_GreaterZero(pGetCoeff(DEN(a)),ntCoeffs)))
[6ad7e6]1424    {
1425      NUM(a)=p_Neg(NUM(a),ntRing);
1426      DEN(a)=p_Neg(DEN(a),ntRing);
1427    }
1428  }
[dc79bd]1429  ntTest(a); // !!!!
[237b4dd]1430}
1431
[2c7f28]1432/* expects *param to be castable to TransExtInfo */
1433static BOOLEAN ntCoeffIsEqual(const coeffs cf, n_coeffType n, void * param)
[6ccdd3a]1434{
[1f414c8]1435  if (ID != n) return FALSE;
[2c7f28]1436  TransExtInfo *e = (TransExtInfo *)param;
1437  /* for rational function fields we expect the underlying
1438     polynomial rings to be IDENTICAL, i.e. the SAME OBJECT;
[6ccdd3a]1439     this expectation is based on the assumption that we have properly
1440     registered cf and perform reference counting rather than creating
1441     multiple copies of the same coefficient field/domain/ring */
[4d94c97]1442  if (ntRing == e->r)
1443    return TRUE;
1444
1445  // NOTE: Q(a)[x] && Q(a)[y] should better share the _same_ Q(a)...
1446  if( rEqual(ntRing, e->r, TRUE) )
1447  {
1448    rDelete(e->r);
1449    return TRUE;
1450  }
1451
[a355723]1452  return FALSE;
[6ccdd3a]1453}
1454
[2c7f28]1455number ntLcm(number a, number b, const coeffs cf)
[6ccdd3a]1456{
[dc79bd]1457  ntTest(a);
1458  ntTest(b);
[1577ebd]1459  fraction fb = (fraction)b;
1460  if ((b==NULL)||(DEN(fb)==NULL)) return ntCopy(a,cf);
[c14846c]1461#ifdef HAVE_FACTORY
[1577ebd]1462  fraction fa = (fraction)a;
1463  /* singclap_gcd destroys its arguments; we hence need copies: */
1464  poly pa = p_Copy(NUM(fa), ntRing);
1465  poly pb = p_Copy(DEN(fb), ntRing);
[c14846c]1466
[0635d51]1467  poly pGcd;
1468  if (nCoeff_is_Q(ntCoeffs))
1469  {
1470    if (p_IsConstant(pa,ntRing) && p_IsConstant(pb,ntRing))
1471    {
1472      pGcd = pa;
1473      p_SetCoeff (pGcd, n_Gcd (pGetCoeff(pGcd), pGetCoeff(pb), ntCoeffs), ntRing);
1474    }
1475    else
1476    {
1477      number contentpa, contentpb, tmp;
1478
1479      contentpb= p_GetCoeff(pb, ntRing);
1480      pIter(pb);
1481      while (pb != NULL)
1482      {
1483        tmp = n_Gcd(contentpb, p_GetCoeff(pb, ntRing) , ntCoeffs);
1484        n_Delete(&contentpb, ntCoeffs);
1485        contentpb = tmp;
1486        pIter(pb);
1487      }
1488
1489      contentpa= p_GetCoeff(pa, ntRing);
1490      pIter(pa);
1491      while (pa != NULL)
1492      {
1493        tmp = n_Gcd(contentpa, p_GetCoeff(pa, ntRing), ntCoeffs);
1494        n_Delete(&contentpa, ntCoeffs);
1495        contentpa = tmp;
1496        pIter(pa);
1497      }
1498
1499      tmp= n_Gcd (contentpb, contentpa, ntCoeffs);
1500      n_Delete(&contentpa, ntCoeffs);
1501      n_Delete(&contentpb, ntCoeffs);
1502      contentpa= tmp;
1503      p_Delete(&pb, ntRing);
1504      p_Delete(&pa, ntRing);
1505
1506      /* singclap_gcd destroys its arguments; we hence need copies: */
1507      pGcd = singclap_gcd(p_Copy(NUM(fa),ntRing), p_Copy(DEN(fb),ntRing), ntRing);
1508      pGcd= p_Mult_nn (pGcd, contentpa, ntRing);
1509      n_Delete(&contentpa, ntCoeffs);
1510    }
1511  }
1512  else
1513    pGcd = singclap_gcd(pa, pb, cf->extRing);
1514
[1577ebd]1515  /* Note that, over Q, singclap_gcd will remove the denominators in all
1516     rational coefficients of pa and pb, before starting to compute
1517     the gcd. Thus, we do not need to ensure that the coefficients of
1518     pa and pb live in Z; they may well be elements of Q\Z. */
[0635d51]1519
[1577ebd]1520  if (p_IsConstant(pGcd, ntRing) &&
1521      n_IsOne(p_GetCoeff(pGcd, ntRing), ntCoeffs))
1522  { /* gcd = 1; return pa*pb*/
1523    p_Delete(&pGcd,ntRing);
1524    fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1525    NUM(result) = pp_Mult_qq(NUM(fa),DEN(fb),ntRing);
[dc79bd]1526
1527    ntTest((number)result); // !!!!
[a355723]1528
[1577ebd]1529    return (number)result;
1530  }
[a355723]1531
[dc79bd]1532
1533  /* return pa*pb/gcd */
[1577ebd]1534    poly newNum = singclap_pdivide(NUM(fa), pGcd, ntRing);
[331fd0]1535    p_Delete(&pGcd,ntRing);
[1577ebd]1536    fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1537    NUM(result) = p_Mult_q(p_Copy(DEN(fb),ntRing),newNum,ntRing);
[dc79bd]1538    ntTest((number)result); // !!!!
[1577ebd]1539    return (number)result;
[a355723]1540
[1577ebd]1541#else
1542  Print("// factory needed: transext.cc:ntLcm\n");
1543  return NULL;
1544#endif /* HAVE_FACTORY */
[2c7f28]1545  return NULL;
[6ccdd3a]1546}
1547
[2c7f28]1548number ntGcd(number a, number b, const coeffs cf)
[6ccdd3a]1549{
[dc79bd]1550  ntTest(a);
1551  ntTest(b);
[1577ebd]1552  if (a==NULL) return ntCopy(b,cf);
1553  if (b==NULL) return ntCopy(a,cf);
[c14846c]1554#ifdef HAVE_FACTORY
[1577ebd]1555  fraction fa = (fraction)a;
1556  fraction fb = (fraction)b;
[0635d51]1557
[1577ebd]1558  poly pa = p_Copy(NUM(fa), ntRing);
1559  poly pb = p_Copy(NUM(fb), ntRing);
[c14846c]1560
[0635d51]1561  poly pGcd;
1562  if (nCoeff_is_Q(ntCoeffs))
1563  {
1564    if (p_IsConstant(pa,ntRing) && p_IsConstant(pb,ntRing))
1565    {
1566      pGcd = pa;
1567      p_SetCoeff (pGcd, n_Gcd (pGetCoeff(pGcd), pGetCoeff(pb), ntCoeffs), ntRing);
1568    }
1569    else
1570    {
1571      number contentpa, contentpb, tmp;
1572
1573      contentpb= p_GetCoeff(pb, ntRing);
1574      pIter(pb);
1575      while (pb != NULL)
1576      {
1577        tmp = n_Gcd(contentpb, p_GetCoeff(pb, ntRing) , ntCoeffs);
1578        n_Delete(&contentpb, ntCoeffs);
1579        contentpb = tmp;
1580        pIter(pb);
1581      }
1582
1583      contentpa= p_GetCoeff(pa, ntRing);
1584      pIter(pa);
1585      while (pa != NULL)
1586      {
1587        tmp = n_Gcd(contentpa, p_GetCoeff(pa, ntRing), ntCoeffs);
1588        n_Delete(&contentpa, ntCoeffs);
1589        contentpa = tmp;
1590        pIter(pa);
1591      }
1592
1593      tmp= n_Gcd (contentpb, contentpa, ntCoeffs);
1594      n_Delete(&contentpa, ntCoeffs);
1595      n_Delete(&contentpb, ntCoeffs);
1596      contentpa= tmp;
1597      p_Delete(&pb, ntRing);
1598      p_Delete(&pa, ntRing);
1599
1600      /* singclap_gcd destroys its arguments; we hence need copies: */
1601      pGcd = singclap_gcd(p_Copy(NUM(fa),ntRing), p_Copy(NUM(fb),ntRing), ntRing);
1602      pGcd= p_Mult_nn (pGcd, contentpa, ntRing);
1603      n_Delete(&contentpa, ntCoeffs);
1604    }
1605  }
1606  else
1607    pGcd = singclap_gcd(pa, pb, cf->extRing);
[1577ebd]1608  /* Note that, over Q, singclap_gcd will remove the denominators in all
1609     rational coefficients of pa and pb, before starting to compute
1610     the gcd. Thus, we do not need to ensure that the coefficients of
1611     pa and pb live in Z; they may well be elements of Q\Z. */
[0635d51]1612
[1577ebd]1613  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1614  NUM(result) = pGcd;
[dc79bd]1615  ntTest((number)result); // !!!!
[1577ebd]1616  return (number)result;
1617#else
1618  Print("// factory needed: transext.cc:ntGcd\n");
[2c7f28]1619  return NULL;
[1577ebd]1620#endif /* HAVE_FACTORY */
[6ccdd3a]1621}
1622
[2c7f28]1623int ntSize(number a, const coeffs cf)
[6ccdd3a]1624{
[2c7f28]1625  ntTest(a);
[e5d267]1626  if (IS0(a)) return -1;
[2c7f28]1627  /* this has been taken from the old implementation of field extensions,
1628     where we computed the sum of the degrees and the numbers of terms in
1629     the numerator and denominator of a; so we leave it at that, for the
1630     time being */
1631  fraction f = (fraction)a;
[e5d267]1632  poly p = NUM(f);
[2c7f28]1633  int noOfTerms = 0;
1634  int numDegree = 0;
1635  while (p != NULL)
1636  {
1637    noOfTerms++;
1638    int d = 0;
1639    for (int i = 1; i <= rVar(ntRing); i++)
1640      d += p_GetExp(p, i, ntRing);
1641    if (d > numDegree) numDegree = d;
1642    pIter(p);
1643  }
1644  int denDegree = 0;
[e5d267]1645  if (!DENIS1(f))
[2c7f28]1646  {
[e5d267]1647    p = DEN(f);
[2c7f28]1648    while (p != NULL)
1649    {
1650      noOfTerms++;
1651      int d = 0;
1652      for (int i = 1; i <= rVar(ntRing); i++)
1653        d += p_GetExp(p, i, ntRing);
1654      if (d > denDegree) denDegree = d;
1655      pIter(p);
1656    }
1657  }
[dc79bd]1658  ntTest(a); // !!!!
[2c7f28]1659  return numDegree + denDegree + noOfTerms;
[6ccdd3a]1660}
1661
[2c7f28]1662number ntInvers(number a, const coeffs cf)
[6ccdd3a]1663{
[2c7f28]1664  ntTest(a);
[dc79bd]1665  if (IS0(a))
[a355723]1666  {
[dc79bd]1667    WerrorS(nDivBy0);
1668    return NULL;
[a355723]1669  }
[2c7f28]1670  fraction f = (fraction)a;
[dc79bd]1671  assume( f != NULL );
1672
[2c7f28]1673  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[dc79bd]1674
1675  assume( NUM(f) != NULL );
1676  const poly den = DEN(f);
[a355723]1677
[dc79bd]1678  if (den == NULL)
1679    NUM(result) = p_One(ntRing);
1680  else
1681    NUM(result) = p_Copy(den, ntRing);
1682
1683  if( !NUMIS1(f) )
1684  {
[a355723]1685    poly num_f=NUM(f);
[e5c9e5]1686    BOOLEAN neg= !n_GreaterZero(pGetCoeff(num_f),ntCoeffs);
[a355723]1687    if (neg)
1688    {
1689      num_f=p_Neg(p_Copy(num_f, ntRing), ntRing);
1690      NUM(result)=p_Neg(NUM(result), ntRing);
1691    }
1692    else
1693    {
1694      num_f=p_Copy(num_f, ntRing);
1695    }
1696    DEN(result) = num_f;
[dc79bd]1697    COM(result) = COM(f);
[a355723]1698    if (neg)
1699    {
1700      if (p_IsOne(num_f, ntRing))
1701      {
1702        DEN(result)=NULL;
1703        COM(result) = 0;
1704        p_Delete(&num_f,ntRing);
1705      }
1706    }
[dc79bd]1707  }
1708  else
1709  {
1710    DEN(result) = NULL;
1711    COM(result) = 0;
1712  }
1713  ntTest((number)result); // !!!!
[2c7f28]1714  return (number)result;
[6ccdd3a]1715}
1716
[2c7f28]1717/* assumes that src = Q, dst = Q(t_1, ..., t_s) */
1718number ntMap00(number a, const coeffs src, const coeffs dst)
[6ccdd3a]1719{
[2c7f28]1720  if (n_IsZero(a, src)) return NULL;
[dc79bd]1721  assume(n_Test(a, src));
[6ccdd3a]1722  assume(src == dst->extRing->cf);
[0c516e]1723  if ((SR_HDL(a) & SR_INT) || (a->s==3))
[9d68fd]1724  {
1725    number res=ntInit(p_NSet(n_Copy(a, src), dst->extRing), dst);
1726    n_Test(res,dst);
1727    return res;
1728  }
[0c516e]1729  number nn=nlGetDenom(a,src);
1730  number zz=nlGetNumerator(a,src);
1731  number res=ntInit(p_NSet(zz,dst->extRing), dst);
1732  fraction ff=(fraction)res;
[9d68fd]1733  if (n_IsOne(nn,src)) DEN(ff)=NULL;
1734  else                 DEN(ff)=p_NSet(nn,dst->extRing);
1735  n_Test((number)ff,dst);
[0c516e]1736  return (number)ff;
[6ccdd3a]1737}
1738
[2c7f28]1739/* assumes that src = Z/p, dst = Q(t_1, ..., t_s) */
1740number ntMapP0(number a, const coeffs src, const coeffs dst)
[6ccdd3a]1741{
[2c7f28]1742  if (n_IsZero(a, src)) return NULL;
[dc79bd]1743  assume(n_Test(a, src));
[6ccdd3a]1744  /* mapping via intermediate int: */
1745  int n = n_Int(a, src);
1746  number q = n_Init(n, dst->extRing->cf);
[2c7f28]1747  if (n_IsZero(q, dst->extRing->cf))
1748  {
1749    n_Delete(&q, dst->extRing->cf);
1750    return NULL;
1751  }
[dc79bd]1752  return ntInit(p_NSet(q, dst->extRing), dst);
[6ccdd3a]1753}
1754
[2c7f28]1755/* assumes that either src = Q(t_1, ..., t_s), dst = Q(t_1, ..., t_s), or
1756                       src = Z/p(t_1, ..., t_s), dst = Z/p(t_1, ..., t_s) */
[b38d70]1757number ntCopyMap(number a, const coeffs cf, const coeffs dst)
[6ccdd3a]1758{
[b38d70]1759//  if (n_IsZero(a, cf)) return NULL;
[a355723]1760
[b38d70]1761  ntTest(a);
1762
1763  if (IS0(a)) return NULL;
[a355723]1764
[b38d70]1765  const ring rSrc = cf->extRing;
1766  const ring rDst = dst->extRing;
[a355723]1767
[b38d70]1768  if( rSrc == rDst )
1769    return ntCopy(a, dst); // USUALLY WRONG!
[a355723]1770
[b38d70]1771  fraction f = (fraction)a;
1772  poly g = prCopyR(NUM(f), rSrc, rDst);
[a355723]1773
[b38d70]1774  poly h = NULL;
[a355723]1775
[b38d70]1776  if (!DENIS1(f))
1777     h = prCopyR(DEN(f), rSrc, rDst);
[a355723]1778
[b38d70]1779  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
[a355723]1780
[b38d70]1781  NUM(result) = g;
1782  DEN(result) = h;
1783  COM(result) = COM(f);
[dc79bd]1784  assume(n_Test((number)result, dst));
[a355723]1785  return (number)result;
[6ccdd3a]1786}
1787
[b38d70]1788number ntCopyAlg(number a, const coeffs cf, const coeffs dst)
[31c731]1789{
[dc79bd]1790  assume( n_Test(a, cf) );
[b38d70]1791  if (n_IsZero(a, cf)) return NULL;
[a355723]1792
1793  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
[b38d70]1794  // DEN(f) = NULL; COM(f) = 0;
1795  NUM(f) = prCopyR((poly)a, cf->extRing, dst->extRing);
[dc79bd]1796  assume(n_Test((number)f, dst));
[31c731]1797  return (number)f;
1798}
1799
[2c7f28]1800/* assumes that src = Q, dst = Z/p(t_1, ..., t_s) */
1801number ntMap0P(number a, const coeffs src, const coeffs dst)
[6ccdd3a]1802{
[dc79bd]1803  assume( n_Test(a, src) );
[2c7f28]1804  if (n_IsZero(a, src)) return NULL;
[6909cfb]1805  // int p = rChar(dst->extRing);
[79020f]1806  number q = nlModP(a, src, dst->extRing->cf);
1807
[2c7f28]1808  if (n_IsZero(q, dst->extRing->cf))
1809  {
1810    n_Delete(&q, dst->extRing->cf);
1811    return NULL;
1812  }
[a355723]1813
[79020f]1814  poly g = p_NSet(q, dst->extRing);
1815
[2c7f28]1816  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
[31c731]1817  NUM(f) = g; // DEN(f) = NULL; COM(f) = 0;
[dc79bd]1818  assume(n_Test((number)f, dst));
[2c7f28]1819  return (number)f;
[6ccdd3a]1820}
1821
[2c7f28]1822/* assumes that src = Z/p, dst = Z/p(t_1, ..., t_s) */
1823number ntMapPP(number a, const coeffs src, const coeffs dst)
[6ccdd3a]1824{
[dc79bd]1825  assume( n_Test(a, src) );
[2c7f28]1826  if (n_IsZero(a, src)) return NULL;
[6ccdd3a]1827  assume(src == dst->extRing->cf);
[2c7f28]1828  poly p = p_One(dst->extRing);
[a0acbc]1829  p_SetCoeff(p, n_Copy(a, src), dst->extRing);
[2c7f28]1830  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
[e5d267]1831  NUM(f) = p; DEN(f) = NULL; COM(f) = 0;
[dc79bd]1832  assume(n_Test((number)f, dst));
[2c7f28]1833  return (number)f;
[6ccdd3a]1834}
1835
[2c7f28]1836/* assumes that src = Z/u, dst = Z/p(t_1, ..., t_s), where u != p */
1837number ntMapUP(number a, const coeffs src, const coeffs dst)
[6ccdd3a]1838{
[dc79bd]1839  assume( n_Test(a, src) );
[2c7f28]1840  if (n_IsZero(a, src)) return NULL;
[6ccdd3a]1841  /* mapping via intermediate int: */
1842  int n = n_Int(a, src);
1843  number q = n_Init(n, dst->extRing->cf);
[2c7f28]1844  poly p;
1845  if (n_IsZero(q, dst->extRing->cf))
1846  {
1847    n_Delete(&q, dst->extRing->cf);
1848    return NULL;
1849  }
1850  p = p_One(dst->extRing);
1851  p_SetCoeff(p, q, dst->extRing);
1852  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
[e5d267]1853  NUM(f) = p; DEN(f) = NULL; COM(f) = 0;
[dc79bd]1854  assume(n_Test((number)f, dst));
[2c7f28]1855  return (number)f;
[6ccdd3a]1856}
1857
[2c7f28]1858nMapFunc ntSetMap(const coeffs src, const coeffs dst)
[6ccdd3a]1859{
[2c7f28]1860  /* dst is expected to be a rational function field */
[1f414c8]1861  assume(getCoeffType(dst) == ID);
[c14846c]1862
[da0565b]1863  if( src == dst ) return ndCopyMap;
1864
[6ccdd3a]1865  int h = 0; /* the height of the extension tower given by dst */
1866  coeffs bDst = nCoeff_bottom(dst, h); /* the bottom field in the tower dst */
[331fd0]1867  coeffs bSrc = nCoeff_bottom(src, h); /* the bottom field in the tower src */
[c14846c]1868
[6ccdd3a]1869  /* for the time being, we only provide maps if h = 1 and if b is Q or
1870     some field Z/pZ: */
[31c731]1871  if (h==0)
1872  {
1873    if (nCoeff_is_Q(src) && nCoeff_is_Q(bDst))
1874      return ntMap00;                                 /// Q       -->  Q(T)
1875    if (nCoeff_is_Zp(src) && nCoeff_is_Q(bDst))
1876      return ntMapP0;                                 /// Z/p     -->  Q(T)
1877    if (nCoeff_is_Q(src) && nCoeff_is_Zp(bDst))
1878      return ntMap0P;                                 /// Q       --> Z/p(T)
1879    if (nCoeff_is_Zp(src) && nCoeff_is_Zp(bDst))
1880    {
1881      if (src->ch == dst->ch) return ntMapPP;         /// Z/p     --> Z/p(T)
1882      else return ntMapUP;                            /// Z/u     --> Z/p(T)
1883    }
1884  }
[6ccdd3a]1885  if (h != 1) return NULL;
1886  if ((!nCoeff_is_Zp(bDst)) && (!nCoeff_is_Q(bDst))) return NULL;
[c14846c]1887
[2c7f28]1888  /* Let T denote the sequence of transcendental extension variables, i.e.,
1889     K[t_1, ..., t_s] =: K[T];
1890     Let moreover, for any such sequence T, T' denote any subsequence of T
1891     of the form t_1, ..., t_w with w <= s. */
[c14846c]1892
[6ccdd3a]1893  if ((!nCoeff_is_Zp(bSrc)) && (!nCoeff_is_Q(bSrc))) return NULL;
[c14846c]1894
[6ccdd3a]1895  if (nCoeff_is_Q(bSrc) && nCoeff_is_Q(bDst))
1896  {
[a355723]1897    if (rVar(src->extRing) > rVar(dst->extRing))
[b38d70]1898       return NULL;
[a355723]1899
[2c7f28]1900    for (int i = 0; i < rVar(src->extRing); i++)
[a355723]1901      if (strcmp(rRingVar(i, src->extRing), rRingVar(i, dst->extRing)) != 0)
1902         return NULL;
1903
[b38d70]1904    if (src->type==n_transExt)
1905       return ntCopyMap;          /// Q(T')   --> Q(T)
1906    else
1907       return ntCopyAlg;
[6ccdd3a]1908  }
[c14846c]1909
[6ccdd3a]1910  if (nCoeff_is_Zp(bSrc) && nCoeff_is_Zp(bDst))
1911  {
[a355723]1912    if (rVar(src->extRing) > rVar(dst->extRing))
[b38d70]1913       return NULL;
[a355723]1914
[2c7f28]1915    for (int i = 0; i < rVar(src->extRing); i++)
[a355723]1916      if (strcmp(rRingVar(i, src->extRing), rRingVar(i, dst->extRing)) != 0)
1917         return NULL;
1918
[b38d70]1919    if (src->type==n_transExt)
1920       return ntCopyMap;         /// Z/p(T') --> Z/p(T)
1921    else
1922       return ntCopyAlg;
[6ccdd3a]1923  }
[c14846c]1924
1925  return NULL;                                 /// default
[6ccdd3a]1926}
[c8e030]1927#if 0
1928nMapFunc ntSetMap_T(const coeffs src, const coeffs dst)
1929{
1930  nMapFunc n=ntSetMap(src,dst);
1931  if (n==ntCopyAlg) printf("n=ntCopyAlg\n");
1932  else if (n==ntCopyMap) printf("n=ntCopyMap\n");
1933  else if (n==ntMapUP) printf("n=ntMapUP\n");
1934  else if (n==ntMap0P) printf("n=ntMap0P\n");
1935  else if (n==ntMapP0) printf("n=ntMapP0\n");
1936  else if (n==ntMap00) printf("n=ntMap00\n");
1937  else if (n==NULL) printf("n=NULL\n");
1938  else printf("n=?\n");
1939  return n;
1940}
1941#endif
[6ccdd3a]1942
[31c731]1943void ntKillChar(coeffs cf)
1944{
[a195ed]1945  if ((--cf->extRing->ref) == 0)
1946    rDelete(cf->extRing);
1947}
[fc4977]1948#ifdef HAVE_FACTORY
1949number ntConvFactoryNSingN( const CanonicalForm n, const coeffs cf)
1950{
1951  if (n.isZero()) return NULL;
1952  poly p=convFactoryPSingP(n,ntRing);
1953  fraction result = (fraction)omAlloc0Bin(fractionObjectBin);
1954  NUM(result) = p;
1955  //DEN(result) = NULL; // done by omAlloc0Bin
1956  //COM(result) = 0; // done by omAlloc0Bin
[dc79bd]1957  ntTest((number)result);
[fc4977]1958  return (number)result;
1959}
[2e4ec14]1960CanonicalForm ntConvSingNFactoryN( number n, BOOLEAN /*setChar*/, const coeffs cf )
[fc4977]1961{
1962  ntTest(n);
1963  if (IS0(n)) return CanonicalForm(0);
1964
1965  fraction f = (fraction)n;
1966  return convSingPFactoryP(NUM(f),ntRing);
1967}
1968#endif
[a195ed]1969
[da5d77]1970static int ntParDeg(number a, const coeffs cf)
[48a41a]1971{
[dc79bd]1972  ntTest(a);
[48a41a]1973  if (IS0(a)) return -1;
1974  fraction fa = (fraction)a;
1975  return cf->extRing->pFDeg(NUM(fa),cf->extRing);
1976}
1977
[7fee876]1978/// return the specified parameter as a number in the given trans.ext.
1979static number ntParameter(const int iParameter, const coeffs cf)
1980{
1981  assume(getCoeffType(cf) == ID);
1982
1983  const ring R = cf->extRing;
1984  assume( R != NULL );
1985  assume( 0 < iParameter && iParameter <= rVar(R) );
1986
1987  poly p = p_One(R); p_SetExp(p, iParameter, 1, R); p_Setm(p, R);
[9d68fd]1988  p_Test(p,R);
[7fee876]1989
1990  fraction f = (fraction)omAlloc0Bin(fractionObjectBin);
1991  NUM(f) = p;
1992  DEN(f) = NULL;
1993  COM(f) = 0;
1994
[dc79bd]1995  ntTest((number)f);
1996
[7fee876]1997  return (number)f;
1998}
1999
2000/// if m == var(i)/1 => return i,
2001int ntIsParam(number m, const coeffs cf)
2002{
[dc79bd]2003  ntTest(m);
[7fee876]2004  assume(getCoeffType(cf) == ID);
2005
2006  const ring R = cf->extRing;
2007  assume( R != NULL );
2008
2009  fraction f = (fraction)m;
2010
2011  if( DEN(f) != NULL )
2012    return 0;
2013
2014  return p_Var( NUM(f), R );
2015}
2016
[dc79bd]2017struct NTNumConverter
2018{
2019  static inline poly convert(const number& n)
2020  {
2021    // suitable for trans. ext. numbers that are fractions of polys
2022    return NUM((fraction)n); // return the numerator
2023  }
2024};
2025
2026
2027static void ntClearContent(ICoeffsEnumerator& numberCollectionEnumerator, number& c, const coeffs cf)
[de88371]2028{
2029  assume(cf != NULL);
2030  assume(getCoeffType(cf) == ID);
[dc79bd]2031  // all coeffs are given by fractions of polynomails over integers!!!
2032  // without denominators!!!
2033
2034  const ring   R = cf->extRing;
2035  assume(R != NULL);
[a355723]2036  const coeffs Q = R->cf;
2037  assume(Q != NULL);
2038  assume(nCoeff_is_Q(Q));
2039
[dc79bd]2040
2041  numberCollectionEnumerator.Reset();
2042
2043  if( !numberCollectionEnumerator.MoveNext() ) // empty zero polynomial?
2044  {
2045    c = ntInit(1, cf);
2046    return;
2047  }
2048
2049  // all coeffs are given by integers after returning from this routine
2050
2051  // part 1, collect product of all denominators /gcds
2052  poly cand = NULL;
2053
2054  do
2055  {
2056    number &n = numberCollectionEnumerator.Current();
2057
2058    ntNormalize(n, cf);
2059
2060    fraction f = (fraction)n;
2061
2062    assume( f != NULL );
2063
2064    const poly den = DEN(f);
2065
2066    assume( den == NULL ); // ?? / 1 ?
2067
2068    const poly num = NUM(f);
2069
2070    if( cand == NULL )
2071      cand = p_Copy(num, R);
2072    else
2073      cand = singclap_gcd(cand, p_Copy(num, R), R); // gcd(cand, num)
2074
2075    if( p_IsConstant(cand, R) )
2076      break;
2077  }
2078  while( numberCollectionEnumerator.MoveNext() ) ;
2079
[a355723]2080
[dc79bd]2081  // part2: all coeffs = all coeffs * cand
2082  if( cand != NULL )
2083  {
2084  if( !p_IsConstant(cand, R) )
2085  {
[a355723]2086    c = ntInit(cand, cf);
[dc79bd]2087    numberCollectionEnumerator.Reset();
2088    while (numberCollectionEnumerator.MoveNext() )
2089    {
2090      number &n = numberCollectionEnumerator.Current();
2091      const number t = ntDiv(n, c, cf); // TODO: rewrite!?
2092      ntDelete(&n, cf);
2093      n = t;
2094    }
2095  } // else NUM (result) = p_One(R);
2096  else { p_Delete(&cand, R); cand = NULL; }
2097  }
[a355723]2098
[dc79bd]2099  // Quick and dirty fix for constant content clearing: consider numerators???
2100  CRecursivePolyCoeffsEnumerator<NTNumConverter> itr(numberCollectionEnumerator); // recursively treat the NUM(numbers) as polys!
2101  number cc;
[a355723]2102
[558f3cc]2103  n_ClearContent(itr, cc, Q);
[dc79bd]2104  number g = ntInit(p_NSet(cc, R), cf);
[a355723]2105
[dc79bd]2106  if( cand != NULL )
2107  {
2108    number gg = ntMult(g, c, cf);
[a355723]2109    ntDelete(&g, cf);
[dc79bd]2110    ntDelete(&c, cf); c = gg;
[a355723]2111  } else
[dc79bd]2112    c = g;
2113  ntTest(c);
[de88371]2114}
2115
[dc79bd]2116static void ntClearDenominators(ICoeffsEnumerator& numberCollectionEnumerator, number& c, const coeffs cf)
[de88371]2117{
2118  assume(cf != NULL);
2119  assume(getCoeffType(cf) == ID); // both over Q(a) and Zp(a)!
[dc79bd]2120  // all coeffs are given by fractions of polynomails over integers!!!
2121
2122  numberCollectionEnumerator.Reset();
2123
2124  if( !numberCollectionEnumerator.MoveNext() ) // empty zero polynomial?
2125  {
2126    c = ntInit(1, cf);
2127    return;
2128  }
2129
2130  // all coeffs are given by integers after returning from this routine
2131
2132  // part 1, collect product of all denominators /gcds
2133  poly cand = NULL;
2134
2135  const ring R = cf->extRing;
[a355723]2136  assume(R != NULL);
[dc79bd]2137
[a355723]2138  const coeffs Q = R->cf;
2139  assume(Q != NULL);
2140//  assume(nCoeff_is_Q(Q));
[dc79bd]2141
2142  do
2143  {
2144    number &n = numberCollectionEnumerator.Current();
[a355723]2145
[dc79bd]2146    ntNormalize(n, cf);
2147
[e56eb1a]2148    fraction f = (fraction)ntGetDenom (n, cf);
[dc79bd]2149
2150    assume( f != NULL );
2151
[e56eb1a]2152    const poly den = NUM(f);
[dc79bd]2153
2154    if( den == NULL ) // ?? / 1 ?
2155      continue;
2156
2157    if( cand == NULL )
2158      cand = p_Copy(den, R);
2159    else
2160    {
2161      // cand === LCM( cand, den )!!!!
2162      // NOTE: maybe it's better to make the product and clearcontent afterwards!?
2163      // TODO: move the following to factory?
2164      poly gcd = singclap_gcd(p_Copy(cand, R), p_Copy(den, R), R); // gcd(cand, den) is monic no mater leading coeffs! :((((
[e56eb1a]2165      if (nCoeff_is_Q (Q))
2166      {
2167        number LcGcd= n_Gcd (p_GetCoeff (cand, R), p_GetCoeff(den, R), Q);
2168        gcd = p_Mult_nn(gcd, LcGcd, R);
[20c99b]2169        n_Delete(&LcGcd,Q);
[e56eb1a]2170      }
[dc79bd]2171//      assume( n_IsOne(pGetCoeff(gcd), Q) ); // TODO: this may be wrong...
[a355723]2172      cand = p_Mult_q(cand, p_Copy(den, R), R); // cand *= den
[dc79bd]2173      const poly t = singclap_pdivide( cand, gcd, R ); // cand' * den / gcd(cand', den)
2174      p_Delete(&cand, R);
2175      p_Delete(&gcd, R);
2176      cand = t;
2177    }
2178  }
2179  while( numberCollectionEnumerator.MoveNext() );
2180
2181  if( cand == NULL )
2182  {
2183    c = ntInit(1, cf);
2184    return;
[a355723]2185  }
[dc79bd]2186
[a355723]2187  c = ntInit(cand, cf);
[dc79bd]2188
2189  numberCollectionEnumerator.Reset();
[a355723]2190
[dc79bd]2191  number d = NULL;
2192
2193  while (numberCollectionEnumerator.MoveNext() )
2194  {
2195    number &n = numberCollectionEnumerator.Current();
2196    number t = ntMult(n, c, cf); // TODO: rewrite!?
2197    ntDelete(&n, cf);
2198
2199    ntNormalize(t, cf); // TODO: needed?
2200    n = t;
[a355723]2201
[dc79bd]2202    fraction f = (fraction)t;
2203    assume( f != NULL );
[de88371]2204
[dc79bd]2205    const poly den = DEN(f);
2206
2207    if( den != NULL ) // ?? / ?? ?
2208    {
2209      assume( p_IsConstant(den, R) );
2210      assume( pNext(den) == NULL );
[a355723]2211
[dc79bd]2212      if( d == NULL )
[a355723]2213        d = n_Copy(pGetCoeff(den), Q);
[dc79bd]2214      else
2215      {
[a355723]2216        number g = n_Lcm(d, pGetCoeff(den), Q);
[dc79bd]2217        n_Delete(&d, Q); d = g;
2218      }
2219    }
2220  }
[a355723]2221
[dc79bd]2222  if( d != NULL )
2223  {
2224    numberCollectionEnumerator.Reset();
2225    while (numberCollectionEnumerator.MoveNext() )
2226    {
2227      number &n = numberCollectionEnumerator.Current();
2228      fraction f = (fraction)n;
2229
2230      assume( f != NULL );
2231
2232      const poly den = DEN(f);
2233
2234      if( den == NULL ) // ?? / 1 ?
2235        NUM(f) = p_Mult_nn(NUM(f), d, R);
2236      else
2237      {
2238        assume( p_IsConstant(den, R) );
2239        assume( pNext(den) == NULL );
[a355723]2240
2241        number ddd = n_Div(d, pGetCoeff(den), Q); // but be an integer now!!!
[dc79bd]2242        NUM(f) = p_Mult_nn(NUM(f), ddd, R);
[a355723]2243        n_Delete(&ddd, Q);
2244
2245        p_Delete(&DEN(f), R);
2246        DEN(f) = NULL; // TODO: check if this is needed!?
2247      }
2248
[dc79bd]2249      assume( DEN(f) == NULL );
2250    }
[a355723]2251
[dc79bd]2252    NUM(c) = p_Mult_nn(NUM(c), d, R);
2253    n_Delete(&d, Q);
2254  }
[a355723]2255
2256
[dc79bd]2257  ntTest(c);
[de88371]2258}
2259
[2c7f28]2260BOOLEAN ntInitChar(coeffs cf, void * infoStruct)
[a55ef0]2261{
2262
2263  assume( infoStruct != NULL );
[c14846c]2264
[2c7f28]2265  TransExtInfo *e = (TransExtInfo *)infoStruct;
[c14846c]2266
[ec5ec8]2267  assume( e->r                != NULL);      // extRing;
2268  assume( e->r->cf            != NULL);      // extRing->cf;
[dd668f]2269  assume( e->r->qideal == NULL );
[2c7f28]2270
[ec5ec8]2271  assume( cf != NULL );
[1f414c8]2272  assume(getCoeffType(cf) == ID);                // coeff type;
[ec5ec8]2273
[7fee876]2274  ring R = e->r;
2275  assume(R != NULL);
[a355723]2276
[7fee876]2277  R->ref ++; // increase the ref.counter for the ground poly. ring!
[ec5ec8]2278
[7fee876]2279  cf->extRing           = R;
[6ccdd3a]2280  /* propagate characteristic up so that it becomes
2281     directly accessible in cf: */
[7fee876]2282  cf->ch = R->cf->ch;
2283  cf->factoryVarOffset = R->cf->factoryVarOffset + rVar(R);
[a355723]2284
[2c7f28]2285  cf->cfGreaterZero  = ntGreaterZero;
2286  cf->cfGreater      = ntGreater;
2287  cf->cfEqual        = ntEqual;
2288  cf->cfIsZero       = ntIsZero;
2289  cf->cfIsOne        = ntIsOne;
2290  cf->cfIsMOne       = ntIsMOne;
2291  cf->cfInit         = ntInit;
[61b2e16]2292  cf->cfInit_bigint  = ntInit_bigint;
[2c7f28]2293  cf->cfInt          = ntInt;
2294  cf->cfNeg          = ntNeg;
2295  cf->cfAdd          = ntAdd;
2296  cf->cfSub          = ntSub;
2297  cf->cfMult         = ntMult;
2298  cf->cfDiv          = ntDiv;
2299  cf->cfExactDiv     = ntDiv;
2300  cf->cfPower        = ntPower;
2301  cf->cfCopy         = ntCopy;
[ce1f78]2302  cf->cfWriteLong    = ntWriteLong;
[2c7f28]2303  cf->cfRead         = ntRead;
[237b4dd]2304  cf->cfNormalize    = ntNormalize;
[2c7f28]2305  cf->cfDelete       = ntDelete;
2306  cf->cfSetMap       = ntSetMap;
2307  cf->cfGetDenom     = ntGetDenom;
2308  cf->cfGetNumerator = ntGetNumerator;
2309  cf->cfRePart       = ntCopy;
2310  cf->cfImPart       = ntImPart;
2311  cf->cfCoeffWrite   = ntCoeffWrite;
[dbcf787]2312#ifdef LDEBUG
[2c7f28]2313  cf->cfDBTest       = ntDBTest;
[dbcf787]2314#endif
[2c7f28]2315  cf->cfGcd          = ntGcd;
2316  cf->cfLcm          = ntLcm;
2317  cf->cfSize         = ntSize;
2318  cf->nCoeffIsEqual  = ntCoeffIsEqual;
2319  cf->cfInvers       = ntInvers;
2320  cf->cfIntDiv       = ntDiv;
[a195ed]2321  cf->cfKillChar     = ntKillChar;
[c14846c]2322
[ce1f78]2323  if( rCanShortOut(ntRing) )
2324    cf->cfWriteShort = ntWriteShort;
2325  else
2326    cf->cfWriteShort = ntWriteLong;
2327
[e5d267]2328#ifndef HAVE_FACTORY
2329  PrintS("// Warning: The 'factory' module is not available.\n");
2330  PrintS("//          Hence gcd's cannot be cancelled in any\n");
2331  PrintS("//          computed fraction!\n");
[fc4977]2332#else
2333  cf->convFactoryNSingN =ntConvFactoryNSingN;
2334  cf->convSingNFactoryN =ntConvSingNFactoryN;
[e5d267]2335#endif
[48a41a]2336  cf->cfParDeg = ntParDeg;
[c14846c]2337
[7fee876]2338  cf->iNumberOfParameters = rVar(R);
2339  cf->pParameterNames = R->names;
2340  cf->cfParameter = ntParameter;
[6637ee]2341
[de88371]2342  if( nCoeff_is_Q(R->cf) )
2343    cf->cfClearContent = ntClearContent;
2344
2345  cf->cfClearDenominators = ntClearDenominators;
[6637ee]2346
[7fee876]2347  return FALSE;
[e82417]2348}
[35564a5]2349
2350template class CRecursivePolyCoeffsEnumerator<NTNumConverter>;
2351template class IEnumerator<snumber*>;
Note: See TracBrowser for help on using the repository browser.