# source:git/libpolys/polys/ext_fields/algext.cc@da0565

spielwiese
Last change on this file since da0565 was da0565, checked in by Oleksandr Motsak <motsak@…>, 12 years ago
fix "equal coeffs test" in preimage.
• Property mode set to `100644`
File size: 24.4 KB
Line
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/**
5  * ABSTRACT: numbers in an algebraic extension field K[a] / < f(a) >
6  *           Assuming that we have a coeffs object cf, then these numbers
7  *           are polynomials in the polynomial ring K[a] represented by
8  *           cf->extRing.
9  *           IMPORTANT ASSUMPTIONS:
10  *           1.) So far we assume that cf->extRing is a valid polynomial
11  *               ring in exactly one variable, i.e., K[a], where K is allowed
12  *               to be any field (representable in SINGULAR and which may
13  *               itself be some extension field, thus allowing for extension
14  *               towers).
15  *           2.) Moreover, this implementation assumes that
16  *               cf->extRing->minideal is not NULL but an ideal with at
17  *               least one non-zero generator which may be accessed by
18  *               cf->extRing->minideal->m and which represents the minimal
19  *               polynomial f(a) of the extension variable 'a' in K[a].
20  *           3.) As soon as an std method for polynomial rings becomes
21  *               availabe, all reduction steps modulo f(a) should be replaced
22  *               by a call to std. Moreover, in this situation one can finally
23  *               move from K[a] / < f(a) > to
24  *                  K[a_1, ..., a_s] / I, with I some zero-dimensional ideal
25  *                                        in K[a_1, ..., a_s] given by a lex
26  *                                        GrÃ¶bner basis.
27  *               The code in algext.h and algext.cc is then capable of
28  *               computing in K[a_1, ..., a_s] / I.
29  **/
30
31#include "config.h"
32#include <misc/auxiliary.h>
33
34#include <omalloc/omalloc.h>
35
36#include <reporter/reporter.h>
37
38#include <coeffs/coeffs.h>
39#include <coeffs/numbers.h>
40#include <coeffs/longrat.h>
41
42#include <polys/monomials/ring.h>
43#include <polys/monomials/p_polys.h>
44#include <polys/simpleideals.h>
45
46#ifdef HAVE_FACTORY
47#include <polys/clapconv.h>
48#include <factory/factory.h>
49#endif
50
51#include "ext_fields/algext.h"
52#define TRANSEXT_PRIVATES 1
53#include "ext_fields/transext.h"
54
55#ifdef LDEBUG
57BOOLEAN  naDBTest(number a, const char *f, const int l, const coeffs r);
58#else
59#define naTest(a) (TRUE)
60#endif
61
62/// Our own type!
63static const n_coeffType ID = n_algExt;
64
65/* polynomial ring in which our numbers live */
66#define naRing cf->extRing
67
68/* coeffs object in which the coefficients of our numbers live;
69 * methods attached to naCoeffs may be used to compute with the
71 * coefficients of our numbers */
72#define naCoeffs cf->extRing->cf
73
74/* minimal polynomial */
75#define naMinpoly naRing->minideal->m
76
77/// forward declarations
78BOOLEAN  naGreaterZero(number a, const coeffs cf);
79BOOLEAN  naGreater(number a, number b, const coeffs cf);
80BOOLEAN  naEqual(number a, number b, const coeffs cf);
81BOOLEAN  naIsOne(number a, const coeffs cf);
82BOOLEAN  naIsMOne(number a, const coeffs cf);
83BOOLEAN  naIsZero(number a, const coeffs cf);
84number   naInit(long i, const coeffs cf);
85int      naInt(number &a, const coeffs cf);
86number   naNeg(number a, const coeffs cf);
87number   naInvers(number a, const coeffs cf);
88number   naAdd(number a, number b, const coeffs cf);
89number   naSub(number a, number b, const coeffs cf);
90number   naMult(number a, number b, const coeffs cf);
91number   naDiv(number a, number b, const coeffs cf);
92void     naPower(number a, int exp, number *b, const coeffs cf);
93number   naCopy(number a, const coeffs cf);
94void     naWrite(number &a, const coeffs cf);
95number   naRePart(number a, const coeffs cf);
96number   naImPart(number a, const coeffs cf);
97number   naGetDenom(number &a, const coeffs cf);
98number   naGetNumerator(number &a, const coeffs cf);
99number   naGcd(number a, number b, const coeffs cf);
100//number   naLcm(number a, number b, const coeffs cf);
101int      naSize(number a, const coeffs cf);
102void     naDelete(number *a, const coeffs cf);
103void     naCoeffWrite(const coeffs cf, BOOLEAN details);
104number   naIntDiv(number a, number b, const coeffs cf);
105const char * naRead(const char *s, number *a, const coeffs cf);
106
107static BOOLEAN naCoeffIsEqual(const coeffs cf, n_coeffType n, void * param);
108
109#ifdef LDEBUG
110BOOLEAN naDBTest(number a, const char *f, const int l, const coeffs cf)
111{
112  assume(getCoeffType(cf) == ID);
113  if (a == NULL) return TRUE;
114  p_Test((poly)a, naRing);
115  if(p_Totaldegree((poly)a, naRing) >= p_Totaldegree(naMinpoly, naRing))
116  {
117    Print("deg >= deg(minpoly) in %s:%d\n",f,l);
118    return FALSE;
119  }
120  return TRUE;
121}
122#endif
123
124void heuristicReduce(poly &p, poly reducer, const coeffs cf);
125void definiteReduce(poly &p, poly reducer, const coeffs cf);
126
127/* returns the bottom field in this field extension tower; if the tower
128   is flat, i.e., if there is no extension, then r itself is returned;
129   as a side-effect, the counter 'height' is filled with the height of
130   the extension tower (in case the tower is flat, 'height' is zero) */
131static coeffs nCoeff_bottom(const coeffs r, int &height)
132{
133  assume(r != NULL);
134  coeffs cf = r;
135  height = 0;
136  while (nCoeff_is_Extension(cf))
137  {
138    assume(cf->extRing != NULL); assume(cf->extRing->cf != NULL);
139    cf = cf->extRing->cf;
140    height++;
141  }
142  return cf;
143}
144
145BOOLEAN naIsZero(number a, const coeffs cf)
146{
147  naTest(a);
148  return (a == NULL);
149}
150
151void naDelete(number * a, const coeffs cf)
152{
153  if (*a == NULL) return;
154  poly aAsPoly = (poly)(*a);
155  p_Delete(&aAsPoly, naRing);
156  *a = NULL;
157}
158
159BOOLEAN naEqual(number a, number b, const coeffs cf)
160{
161  naTest(a); naTest(b);
162
163  /// simple tests
164  if (a == b) return TRUE;
165  if ((a == NULL) && (b != NULL)) return FALSE;
166  if ((b == NULL) && (a != NULL)) return FALSE;
167
168  /// deg test
170  if (a != NULL) aDeg = p_Totaldegree((poly)a, naRing);
171  int bDeg = 0;
172  if (b != NULL) bDeg = p_Totaldegree((poly)b, naRing);
173  if (aDeg != bDeg) return FALSE;
174
175  /// subtraction test
176  number c = naSub(a, b, cf);
177  BOOLEAN result = naIsZero(c, cf);
179  return result;
180}
181
182number naCopy(number a, const coeffs cf)
183{
184  naTest(a);
185  if (a == NULL) return NULL;
186  return (number)p_Copy((poly)a, naRing);
187}
188
189number naGetNumerator(number &a, const coeffs cf)
190{
191  return naCopy(a, cf);
192}
193
194number naGetDenom(number &a, const coeffs cf)
195{
196  naTest(a);
197  return naInit(1, cf);
198}
199
200BOOLEAN naIsOne(number a, const coeffs cf)
201{
202  naTest(a);
203  poly aAsPoly = (poly)a;
204  if (!p_IsConstant(aAsPoly, naRing)) return FALSE;
205  return n_IsOne(p_GetCoeff(aAsPoly, naRing), naCoeffs);
206}
207
208BOOLEAN naIsMOne(number a, const coeffs cf)
209{
210  naTest(a);
211  poly aAsPoly = (poly)a;
212  if (!p_IsConstant(aAsPoly, naRing)) return FALSE;
213  return n_IsMOne(p_GetCoeff(aAsPoly, naRing), naCoeffs);
214}
215
216/// this is in-place, modifies a
217number naNeg(number a, const coeffs cf)
218{
219  naTest(a);
220  if (a != NULL) a = (number)p_Neg((poly)a, naRing);
221  return a;
222}
223
224number naImPart(number a, const coeffs cf)
225{
226  naTest(a);
227  return NULL;
228}
229
230number naInit_bigint(number longratBigIntNumber, const coeffs src, const coeffs cf)
231{
232  assume( cf != NULL );
233
234  const ring A = cf->extRing;
235
236  assume( A != NULL );
237
238  const coeffs C = A->cf;
239
240  assume( C != NULL );
241
242  number n = n_Init_bigint(longratBigIntNumber, src, C);
243
244  if ( n_IsZero(n, C) )
245  {
246    n_Delete(&n, C);
247    return NULL;
248  }
249
250  return (number)p_NSet(n, A);
251}
252
253
254
255number naInit(long i, const coeffs cf)
256{
257  if (i == 0) return NULL;
258  else        return (number)p_ISet(i, naRing);
259}
260
261int naInt(number &a, const coeffs cf)
262{
263  naTest(a);
264  poly aAsPoly = (poly)a;
265  if (!p_IsConstant(aAsPoly, naRing)) return 0;
266  return n_Int(p_GetCoeff(aAsPoly, naRing), naCoeffs);
267}
268
269/* TRUE iff (a != 0 and (b == 0 or deg(a) > deg(b))) */
270BOOLEAN naGreater(number a, number b, const coeffs cf)
271{
272  naTest(a); naTest(b);
273  if (naIsZero(a, cf)) return FALSE;
274  if (naIsZero(b, cf)) return TRUE;
276  if (a != NULL) aDeg = p_Totaldegree((poly)a, naRing);
277  int bDeg = 0;
278  if (b != NULL) bDeg = p_Totaldegree((poly)b, naRing);
280}
281
282/* TRUE iff a != 0 and (LC(a) > 0 or deg(a) > 0) */
283BOOLEAN naGreaterZero(number a, const coeffs cf)
284{
285  naTest(a);
286  if (a == NULL)                                            return FALSE;
287  if (n_GreaterZero(p_GetCoeff((poly)a, naRing), naCoeffs)) return TRUE;
288  if (p_Totaldegree((poly)a, naRing) > 0)                   return TRUE;
289  return FALSE;
290}
291
292void naCoeffWrite(const coeffs cf, BOOLEAN details)
293{
294  assume( cf != NULL );
295
296  const ring A = cf->extRing;
297
298  assume( A != NULL );
299  assume( A->cf != NULL );
300
301  n_CoeffWrite(A->cf, details);
302
303//  rWrite(A);
304
305  const int P = rVar(A);
306  assume( P > 0 );
307
308  Print("//   %d parameter    : ", P);
309
310  for (int nop=0; nop < P; nop ++)
311    Print("%s ", rRingVar(nop, A));
312
313  PrintLn();
314
315  const ideal I = A->minideal;
316
317  assume( I != NULL );
318  assume( IDELEMS(I) == 1 );
319
320
321  if ( details )
322  {
323    PrintS("//   minpoly        : (");
324    p_Write0( I->m, A);
325    PrintS(")");
326  }
327  else
328    PrintS("//   minpoly        : ...");
329
330  PrintLn();
331
332/*
333  char *x = rRingVar(0, A);
334
335  Print("//   Coefficients live in the extension field K[%s]/<f(%s)>\n", x, x);
336  Print("//   with the minimal polynomial f(%s) = %s\n", x,
337        p_String(A->minideal->m, A));
338  PrintS("//   and K: ");
339*/
340}
341
342number naAdd(number a, number b, const coeffs cf)
343{
344  naTest(a); naTest(b);
345  if (a == NULL) return naCopy(b, cf);
346  if (b == NULL) return naCopy(a, cf);
347  poly aPlusB = p_Add_q(p_Copy((poly)a, naRing),
348                        p_Copy((poly)b, naRing), naRing);
349  definiteReduce(aPlusB, naMinpoly, cf);
350  return (number)aPlusB;
351}
352
353number naSub(number a, number b, const coeffs cf)
354{
355  naTest(a); naTest(b);
356  if (b == NULL) return naCopy(a, cf);
357  poly minusB = p_Neg(p_Copy((poly)b, naRing), naRing);
358  if (a == NULL) return (number)minusB;
359  poly aMinusB = p_Add_q(p_Copy((poly)a, naRing), minusB, naRing);
360  definiteReduce(aMinusB, naMinpoly, cf);
361  return (number)aMinusB;
362}
363
364number naMult(number a, number b, const coeffs cf)
365{
366  naTest(a); naTest(b);
367  if (a == NULL) return NULL;
368  if (b == NULL) return NULL;
369  poly aTimesB = p_Mult_q(p_Copy((poly)a, naRing),
370                          p_Copy((poly)b, naRing), naRing);
371  definiteReduce(aTimesB, naMinpoly, cf);
372  return (number)aTimesB;
373}
374
375number naDiv(number a, number b, const coeffs cf)
376{
377  naTest(a); naTest(b);
378  if (b == NULL) WerrorS(nDivBy0);
379  if (a == NULL) return NULL;
380  poly bInverse = (poly)naInvers(b, cf);
381  poly aDivB = p_Mult_q(p_Copy((poly)a, naRing), bInverse, naRing);
384}
385
386/* 0^0 = 0;
387   for |exp| <= 7 compute power by a simple multiplication loop;
388   for |exp| >= 8 compute power along binary presentation of |exp|, e.g.
389      p^13 = p^1 * p^4 * p^8, where we utilise that
390      p^(2^(k+1)) = p^(2^k) * p^(2^k);
391   intermediate reduction modulo the minimal polynomial is controlled by
392   the in-place method heuristicReduce(poly, poly, coeffs); see there.
393*/
394void naPower(number a, int exp, number *b, const coeffs cf)
395{
396  naTest(a);
397
398  /* special cases first */
399  if (a == NULL)
400  {
401    if (exp >= 0) *b = NULL;
402    else          WerrorS(nDivBy0);
403  }
404  else if (exp ==  0) { *b = naInit(1, cf); return; }
405  else if (exp ==  1) { *b = naCopy(a, cf); return; }
406  else if (exp == -1) { *b = naInvers(a, cf); return; }
407
408  int expAbs = exp; if (expAbs < 0) expAbs = -expAbs;
409
410  /* now compute a^expAbs */
411  poly pow; poly aAsPoly = (poly)a;
412  if (expAbs <= 7)
413  {
414    pow = p_Copy(aAsPoly, naRing);
415    for (int i = 2; i <= expAbs; i++)
416    {
417      pow = p_Mult_q(pow, p_Copy(aAsPoly, naRing), naRing);
418      heuristicReduce(pow, naMinpoly, cf);
419    }
420    definiteReduce(pow, naMinpoly, cf);
421  }
422  else
423  {
424    pow = p_ISet(1, naRing);
425    poly factor = p_Copy(aAsPoly, naRing);
426    while (expAbs != 0)
427    {
428      if (expAbs & 1)
429      {
430        pow = p_Mult_q(pow, p_Copy(factor, naRing), naRing);
431        heuristicReduce(pow, naMinpoly, cf);
432      }
433      expAbs = expAbs / 2;
434      if (expAbs != 0)
435      {
436        factor = p_Mult_q(factor, p_Copy(factor, naRing), naRing);
437        heuristicReduce(factor, naMinpoly, cf);
438      }
439    }
440    p_Delete(&factor, naRing);
441    definiteReduce(pow, naMinpoly, cf);
442  }
443
444  /* invert if original exponent was negative */
445  number n = (number)pow;
446  if (exp < 0)
447  {
448    number m = naInvers(n, cf);
450    n = m;
451  }
452  *b = n;
453}
454
455/* may reduce p modulo the reducer by calling definiteReduce;
456   the decision is made based on the following heuristic
457   (which should also only be changed here in this method):
458      if (deg(p) > 10*deg(reducer) then perform reduction;
459   modifies p */
460void heuristicReduce(poly &p, poly reducer, const coeffs cf)
461{
462  #ifdef LDEBUG
463  p_Test((poly)p, naRing);
464  p_Test((poly)reducer, naRing);
465  #endif
466  if (p_Totaldegree(p, naRing) > 10 * p_Totaldegree(reducer, naRing))
467    definiteReduce(p, reducer, cf);
468}
469
470void naWrite(number &a, const coeffs cf)
471{
472  naTest(a);
473  if (a == NULL)
474    StringAppendS("0");
475  else
476  {
477    poly aAsPoly = (poly)a;
478    /* basically, just write aAsPoly using p_Write,
479       but use brackets around the output, if a is not
480       a constant living in naCoeffs = cf->extRing->cf */
481    BOOLEAN useBrackets = !(p_IsConstant(aAsPoly, naRing));
482    if (useBrackets) StringAppendS("(");
483    p_String0(aAsPoly, naRing, naRing);
484    if (useBrackets) StringAppendS(")");
485  }
486}
487
488const char * naRead(const char *s, number *a, const coeffs cf)
489{
490  poly aAsPoly;
491  const char * result = p_Read(s, aAsPoly, naRing);
492  definiteReduce(aAsPoly, naMinpoly, cf);
493  *a = (number)aAsPoly;
494  return result;
495}
496
497#if 0
498/* implemented by the rule lcm(a, b) = a * b / gcd(a, b) */
499number naLcm(number a, number b, const coeffs cf)
500{
501  naTest(a); naTest(b);
502  if (a == NULL) return NULL;
503  if (b == NULL) return NULL;
504  number theProduct = (number)p_Mult_q(p_Copy((poly)a, naRing),
505                                       p_Copy((poly)b, naRing), naRing);
506  /* note that theProduct needs not be reduced w.r.t. naMinpoly;
507     but the final division will take care of the necessary reduction */
508  number theGcd = naGcd(a, b, cf);
510}
511#endif
512
513/* expects *param to be castable to AlgExtInfo */
514static BOOLEAN naCoeffIsEqual(const coeffs cf, n_coeffType n, void * param)
515{
516  if (ID != n) return FALSE;
517  AlgExtInfo *e = (AlgExtInfo *)param;
518  /* for extension coefficient fields we expect the underlying
519     polynomial rings to be IDENTICAL, i.e. the SAME OBJECT;
520     this expectation is based on the assumption that we have properly
521     registered cf and perform reference counting rather than creating
522     multiple copies of the same coefficient field/domain/ring */
523  if (naRing == e->r)
524    return TRUE;
525  /* (Note that then also the minimal ideals will necessarily be
526     the same, as they are attached to the ring.) */
527
528  // NOTE: Q(a)[x] && Q(a)[y] should better share the _same_ Q(a)...
529  if( rEqual(naRing, e->r, TRUE) )
530  {
531    const ideal mi = naRing->minideal;
532    assume( IDELEMS(mi) == 1 );
533    ideal ii = e->i;
534    assume( IDELEMS(ii) == 1 );
535
536    // TODO: the following should be extended for 2 *equal* rings...
537    if( p_EqualPolys(mi->m, ii->m, naRing, e->r) )
538    {
539      id_Delete(&ii, e->r);
540      rDelete(e->r);
541      return TRUE;
542    }
543  }
544
545  return FALSE;
546
547}
548
549int naSize(number a, const coeffs cf)
550{
551  if (a == NULL) return -1;
552  /* this has been taken from the old implementation of field extensions,
553     where we computed the sum of the degree and the number of terms in
554     (poly)a; so we leave it at that, for the time being;
555     maybe, the number of terms alone is a better measure? */
556  poly aAsPoly = (poly)a;
557  int theDegree = 0; int noOfTerms = 0;
558  while (aAsPoly != NULL)
559  {
560    noOfTerms++;
561    int d = p_GetExp(aAsPoly, 1, naRing);
562    if (d > theDegree) theDegree = d;
563    pIter(aAsPoly);
564  }
565  return theDegree + noOfTerms;
566}
567
568/* performs polynomial division and overrides p by the remainder
569   of division of p by the reducer;
570   modifies p */
571void definiteReduce(poly &p, poly reducer, const coeffs cf)
572{
573  #ifdef LDEBUG
574  p_Test((poly)p, naRing);
575  p_Test((poly)reducer, naRing);
576  #endif
577  p_PolyDiv(p, reducer, FALSE, naRing);
578}
579
580void  naNormalize(number &a, const coeffs cf)
581{
582  poly aa=(poly)a;
583  definiteReduce(aa,naMinpoly,cf);
584  a=(number)aa;
585}
586
587#ifdef HAVE_FACTORY
588number naConvFactoryNSingN( const CanonicalForm n, const coeffs cf)
589{
590  if (n.isZero()) return NULL;
591  poly p=convFactoryPSingP(n,naRing);
592  return (number)p;
593}
594CanonicalForm naConvSingNFactoryN( number n, BOOLEAN /*setChar*/, const coeffs cf )
595{
596  naTest(n);
597  if (n==NULL) return CanonicalForm(0);
598
599  return convSingPFactoryP((poly)n,naRing);
600}
601#endif
602
603
604/* IMPORTANT NOTE: Since an algebraic field extension is again a field,
605                   the gcd of two elements is not very interesting. (It
606                   is actually any unit in the field, i.e., any non-
607                   zero element.) Note that the below method does not operate
608                   in this strong sense but rather computes the gcd of
609                   two given elements in the underlying polynomial ring. */
610number naGcd(number a, number b, const coeffs cf)
611{
612  naTest(a); naTest(b);
613  if ((a == NULL) && (b == NULL)) WerrorS(nDivBy0);
614  return (number)p_Gcd((poly)a, (poly)b, naRing);
615}
616
617number naInvers(number a, const coeffs cf)
618{
619  naTest(a);
620  if (a == NULL) WerrorS(nDivBy0);
621  poly aFactor = NULL; poly mFactor = NULL;
622  poly theGcd = p_ExtGcd((poly)a, aFactor, naMinpoly, mFactor, naRing);
623  naTest((number)theGcd); naTest((number)aFactor); naTest((number)mFactor);
624  /* the gcd must be 1 since naMinpoly is irreducible and a != NULL: */
625  assume(naIsOne((number)theGcd, cf));
626  p_Delete(&theGcd, naRing);
627  p_Delete(&mFactor, naRing);
628  return (number)(aFactor);
629}
630
631/* assumes that src = Q, dst = Q(a) */
632number naMap00(number a, const coeffs src, const coeffs dst)
633{
634  if (n_IsZero(a, src)) return NULL;
635  assume(src == dst->extRing->cf);
636  poly result = p_One(dst->extRing);
637  p_SetCoeff(result, n_Copy(a, src), dst->extRing);
638  return (number)result;
639}
640
641/* assumes that src = Z/p, dst = Q(a) */
642number naMapP0(number a, const coeffs src, const coeffs dst)
643{
644  if (n_IsZero(a, src)) return NULL;
645  /* mapping via intermediate int: */
646  int n = n_Int(a, src);
647  number q = n_Init(n, dst->extRing->cf);
648  poly result = p_One(dst->extRing);
649  p_SetCoeff(result, q, dst->extRing);
650  return (number)result;
651}
652
653#if 0
654/* assumes that either src = Q(a), dst = Q(a), or
655                       src = Z/p(a), dst = Z/p(a)     */
656number naCopyMap(number a, const coeffs src, const coeffs dst)
657{
658  return naCopy(a, dst);
659}
660#endif
661
662number naCopyExt(number a, const coeffs src, const coeffs)
663{
664  fraction fa=(fraction)a;
665  return (number)p_Copy(NUM(fa),src->extRing);
666}
667
668/* assumes that src = Q, dst = Z/p(a) */
669number naMap0P(number a, const coeffs src, const coeffs dst)
670{
671  if (n_IsZero(a, src)) return NULL;
672  int p = rChar(dst->extRing);
673
674  number q = nlModP(a, src, dst->extRing->cf);
675
676  poly result = p_NSet(q, dst->extRing);
677
678  return (number)result;
679}
680
681/* assumes that src = Z/p, dst = Z/p(a) */
682number naMapPP(number a, const coeffs src, const coeffs dst)
683{
684  if (n_IsZero(a, src)) return NULL;
685  assume(src == dst->extRing->cf);
686  poly result = p_One(dst->extRing);
687  p_SetCoeff(result, n_Copy(a, src), dst->extRing);
688  return (number)result;
689}
690
691/* assumes that src = Z/u, dst = Z/p(a), where u != p */
692number naMapUP(number a, const coeffs src, const coeffs dst)
693{
694  if (n_IsZero(a, src)) return NULL;
695  /* mapping via intermediate int: */
696  int n = n_Int(a, src);
697  number q = n_Init(n, dst->extRing->cf);
698  poly result = p_One(dst->extRing);
699  p_SetCoeff(result, q, dst->extRing);
700  return (number)result;
701}
702
703nMapFunc naSetMap(const coeffs src, const coeffs dst)
704{
705  /* dst is expected to be an algebraic field extension */
706  assume(getCoeffType(dst) == ID);
707
708  if( src == dst ) return ndCopyMap;
709
710  int h = 0; /* the height of the extension tower given by dst */
711  coeffs bDst = nCoeff_bottom(dst, h); /* the bottom field in the tower dst */
712  coeffs bSrc = nCoeff_bottom(src, h); /* the bottom field in the tower src */
713
714  /* for the time being, we only provide maps if h = 1 and if b is Q or
715     some field Z/pZ: */
716  if (h==0)
717  {
718    if (nCoeff_is_Q(src) && nCoeff_is_Q(bDst))
719      return naMap00;                            /// Q     -->  Q(a)
720    if (nCoeff_is_Zp(src) && nCoeff_is_Q(bDst))
721      return naMapP0;                            /// Z/p   -->  Q(a)
722    if (nCoeff_is_Q(src) && nCoeff_is_Zp(bDst))
723      return naMap0P;                            /// Q      --> Z/p(a)
724    if (nCoeff_is_Zp(src) && nCoeff_is_Zp(bDst))
725    {
726      if (src->ch == dst->ch) return naMapPP;    /// Z/p    --> Z/p(a)
727      else return naMapUP;                       /// Z/u    --> Z/p(a)
728    }
729  }
730  if (h != 1) return NULL;
731  if ((!nCoeff_is_Zp(bDst)) && (!nCoeff_is_Q(bDst))) return NULL;
732  if ((!nCoeff_is_Zp(bSrc)) && (!nCoeff_is_Q(bSrc))) return NULL;
733
734  if (nCoeff_is_Q(bSrc) && nCoeff_is_Q(bDst))
735  {
736    if (strcmp(rRingVar(0, src->extRing),
737               rRingVar(0, dst->extRing)) == 0)
738    {
739      if (src->type==n_algExt)
740         return ndCopyMap; // naCopyMap;         /// Q(a)   --> Q(a)
741      else
742         return naCopyExt;
743    }
744    else
745      return NULL;                               /// Q(b)   --> Q(a)
746  }
747
748  if (nCoeff_is_Zp(bSrc) && nCoeff_is_Zp(bDst))
749  {
750    if (strcmp(rRingVar(0,src->extRing),rRingVar(0,dst->extRing))==0)
751    {
752      if (src->type==n_algExt)
753        return ndCopyMap; // naCopyMap;          /// Z/p(a) --> Z/p(a)
754      else
755         return naCopyExt;
756    }
757    else
758      return NULL;                               /// Z/p(b) --> Z/p(a)
759  }
760
761  return NULL;                                           /// default
762}
763
764int naParDeg(number a, const coeffs cf)
765{
766  if (a == NULL) return -1;
767  poly aa=(poly)a;
768  return cf->extRing->pFDeg(aa,cf->extRing);
769}
770
771BOOLEAN naInitChar(coeffs cf, void * infoStruct)
772{
773  assume( infoStruct != NULL );
774
775  AlgExtInfo *e = (AlgExtInfo *)infoStruct;
776  /// first check whether cf->extRing != NULL and delete old ring???
777
778  assume(e->r                     != NULL);      // extRing;
779  assume(e->r->cf                 != NULL);      // extRing->cf;
780  assume((e->i          != NULL) &&    // minideal has one
781         (IDELEMS(e->i) != 0)    &&    // non-zero generator
782         (e->i->m    != NULL)    ); // at m;
783
784  assume( e->r->minideal == NULL );
785
786  assume( cf != NULL );
787  assume(getCoeffType(cf) == ID);                     // coeff type;
788
789  cf->extRing           = e->r;
790  cf->extRing->ref ++; // increase the ref.counter for the ground poly. ring!
791
792  cf->extRing->minideal = e->i; // make a copy?
793
794  /* propagate characteristic up so that it becomes
795     directly accessible in cf: */
796  cf->ch = cf->extRing->cf->ch;
797
798  #ifdef LDEBUG
799  p_Test((poly)naMinpoly, naRing);
800  #endif
801
802  cf->cfGreaterZero  = naGreaterZero;
803  cf->cfGreater      = naGreater;
804  cf->cfEqual        = naEqual;
805  cf->cfIsZero       = naIsZero;
806  cf->cfIsOne        = naIsOne;
807  cf->cfIsMOne       = naIsMOne;
808  cf->cfInit         = naInit;
809  cf->cfInit_bigint  = naInit_bigint;
810  cf->cfInt          = naInt;
811  cf->cfNeg          = naNeg;
813  cf->cfSub          = naSub;
814  cf->cfMult         = naMult;
817  cf->cfPower        = naPower;
818  cf->cfCopy         = naCopy;
819  cf->cfWrite        = naWrite;
822  cf->cfSetMap       = naSetMap;
823  cf->cfGetDenom     = naGetDenom;
824  cf->cfGetNumerator = naGetNumerator;
825  cf->cfRePart       = naCopy;
826  cf->cfImPart       = naImPart;
827  cf->cfCoeffWrite   = naCoeffWrite;
828  cf->cfNormalize    = naNormalize;
829#ifdef LDEBUG
831#endif
832  cf->cfGcd          = naGcd;
833  //cf->cfLcm          = naLcm;
834  cf->cfSize         = naSize;
835  cf->nCoeffIsEqual  = naCoeffIsEqual;
836  cf->cfInvers       = naInvers;
838#ifdef HAVE_FACTORY
839  cf->convFactoryNSingN=naConvFactoryNSingN;
840  cf->convSingNFactoryN=naConvSingNFactoryN;
841#endif
842  cf->cfParDeg = naParDeg;
843
844  return FALSE;
845}
846
847
848number naParam(const short iParameter, const coeffs cf)
849{
850  assume(getCoeffType(cf) == ID);
851
852  const ring R = cf->extRing;
853  assume( R != NULL );
854  assume( 0 < iParameter && iParameter <= rVar(R) );
855
856  poly p = p_One(R); p_SetExp(p, iParameter, 1, R); p_Setm(p, R);
857
858  return (number) p;
859}
860
861
862/// if m == var(i)/1 => return i,
863int naIsParam(number m, const coeffs cf)
864{
865  assume(getCoeffType(cf) == ID);
866
867  const ring R = cf->extRing;
868  assume( R != NULL );
869
870  return p_Var( (poly)m, R );
871}
Note: See TracBrowser for help on using the repository browser.