source: git/libpolys/polys/ext_fields/algext.cc @ 28216f

fieker-DuValspielwiese
Last change on this file since 28216f was 2f3764, checked in by Oleksandr Motsak <motsak@…>, 12 years ago
n_Init and p_ISet take a "long" (instead of "int") argument
  • Property mode set to 100644
File size: 23.9 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[0] 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
56#define naTest(a) naDBTest(a,__FILE__,__LINE__,cf)
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
70 * coefficients of our numbers, e.g., use naCoeffs->nAdd to add
71 * coefficients of our numbers */
72#define naCoeffs cf->extRing->cf
73
74/* minimal polynomial */
75#define naMinpoly naRing->minideal->m[0]
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
169  int aDeg = 0;
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);
178  naDelete(&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;
275  int aDeg = 0;
276  if (a != NULL) aDeg = p_Totaldegree((poly)a, naRing);
277  int bDeg = 0;
278  if (b != NULL) bDeg = p_Totaldegree((poly)b, naRing);
279  return (aDeg > bDeg);
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[0], 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[0], 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);
382  definiteReduce(aDivB, naMinpoly, cf);
383  return (number)aDivB;
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);
449    naDelete(&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);
509  return naDiv(theProduct, theGcd, 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  return (naRing == e->r);
524  /* (Note that then also the minimal ideals will necessarily be
525     the same, as they are attached to the ring.) */
526}
527
528int naSize(number a, const coeffs cf)
529{
530  if (a == NULL) return -1;
531  /* this has been taken from the old implementation of field extensions,
532     where we computed the sum of the degree and the number of terms in
533     (poly)a; so we leave it at that, for the time being;
534     maybe, the number of terms alone is a better measure? */
535  poly aAsPoly = (poly)a;
536  int theDegree = 0; int noOfTerms = 0;
537  while (aAsPoly != NULL)
538  {
539    noOfTerms++;
540    int d = p_GetExp(aAsPoly, 1, naRing);
541    if (d > theDegree) theDegree = d;
542    pIter(aAsPoly);
543  }
544  return theDegree + noOfTerms;
545}
546
547/* performs polynomial division and overrides p by the remainder
548   of division of p by the reducer;
549   modifies p */
550void definiteReduce(poly &p, poly reducer, const coeffs cf)
551{
552  #ifdef LDEBUG
553  p_Test((poly)p, naRing);
554  p_Test((poly)reducer, naRing);
555  #endif
556  p_PolyDiv(p, reducer, FALSE, naRing);
557}
558
559void  naNormalize(number &a, const coeffs cf)
560{
561  poly aa=(poly)a;
562  definiteReduce(aa,naMinpoly,cf); 
563  a=(number)aa;
564}
565
566#ifdef HAVE_FACTORY
567number naConvFactoryNSingN( const CanonicalForm n, const coeffs cf)
568{
569  if (n.isZero()) return NULL;
570  poly p=convFactoryPSingP(n,naRing);
571  return (number)p;
572}
573CanonicalForm naConvSingNFactoryN( number n, BOOLEAN /*setChar*/, const coeffs cf )
574{
575  naTest(n);
576  if (n==NULL) return CanonicalForm(0);
577
578  return convSingPFactoryP((poly)n,naRing);
579}
580#endif
581
582
583/* IMPORTANT NOTE: Since an algebraic field extension is again a field,
584                   the gcd of two elements is not very interesting. (It
585                   is actually any unit in the field, i.e., any non-
586                   zero element.) Note that the below method does not operate
587                   in this strong sense but rather computes the gcd of
588                   two given elements in the underlying polynomial ring. */
589number naGcd(number a, number b, const coeffs cf)
590{
591  naTest(a); naTest(b);
592  if ((a == NULL) && (b == NULL)) WerrorS(nDivBy0);
593  return (number)p_Gcd((poly)a, (poly)b, naRing);
594}
595
596number naInvers(number a, const coeffs cf)
597{
598  naTest(a);
599  if (a == NULL) WerrorS(nDivBy0);
600  poly aFactor = NULL; poly mFactor = NULL;
601  poly theGcd = p_ExtGcd((poly)a, aFactor, naMinpoly, mFactor, naRing);
602  naTest((number)theGcd); naTest((number)aFactor); naTest((number)mFactor);
603  /* the gcd must be 1 since naMinpoly is irreducible and a != NULL: */
604  assume(naIsOne((number)theGcd, cf));     
605  p_Delete(&theGcd, naRing);
606  p_Delete(&mFactor, naRing);
607  return (number)(aFactor);
608}
609
610/* assumes that src = Q, dst = Q(a) */
611number naMap00(number a, const coeffs src, const coeffs dst)
612{
613  if (n_IsZero(a, src)) return NULL;
614  assume(src == dst->extRing->cf);
615  poly result = p_One(dst->extRing);
616  p_SetCoeff(result, n_Copy(a, src), dst->extRing);
617  return (number)result;
618}
619
620/* assumes that src = Z/p, dst = Q(a) */
621number naMapP0(number a, const coeffs src, const coeffs dst)
622{
623  if (n_IsZero(a, src)) return NULL;
624  /* mapping via intermediate int: */
625  int n = n_Int(a, src);
626  number q = n_Init(n, dst->extRing->cf);
627  poly result = p_One(dst->extRing);
628  p_SetCoeff(result, q, dst->extRing);
629  return (number)result;
630}
631
632#if 0
633/* assumes that either src = Q(a), dst = Q(a), or
634                       src = Z/p(a), dst = Z/p(a)     */
635number naCopyMap(number a, const coeffs src, const coeffs dst)
636{
637  return naCopy(a, dst);
638}
639#endif
640
641number naCopyExt(number a, const coeffs src, const coeffs)
642{
643  fraction fa=(fraction)a;
644  return (number)p_Copy(NUM(fa),src->extRing);
645}
646
647/* assumes that src = Q, dst = Z/p(a) */
648number naMap0P(number a, const coeffs src, const coeffs dst)
649{
650  if (n_IsZero(a, src)) return NULL;
651  int p = rChar(dst->extRing);
652
653  number q = nlModP(a, src, dst->extRing->cf);
654
655  poly result = p_NSet(q, dst->extRing);
656 
657  return (number)result;
658}
659
660/* assumes that src = Z/p, dst = Z/p(a) */
661number naMapPP(number a, const coeffs src, const coeffs dst)
662{
663  if (n_IsZero(a, src)) return NULL;
664  assume(src == dst->extRing->cf);
665  poly result = p_One(dst->extRing);
666  p_SetCoeff(result, n_Copy(a, src), dst->extRing);
667  return (number)result;
668}
669
670/* assumes that src = Z/u, dst = Z/p(a), where u != p */
671number naMapUP(number a, const coeffs src, const coeffs dst)
672{
673  if (n_IsZero(a, src)) return NULL;
674  /* mapping via intermediate int: */
675  int n = n_Int(a, src);
676  number q = n_Init(n, dst->extRing->cf);
677  poly result = p_One(dst->extRing);
678  p_SetCoeff(result, q, dst->extRing);
679  return (number)result;
680}
681
682nMapFunc naSetMap(const coeffs src, const coeffs dst)
683{
684  /* dst is expected to be an algebraic field extension */
685  assume(getCoeffType(dst) == ID);
686 
687  int h = 0; /* the height of the extension tower given by dst */
688  coeffs bDst = nCoeff_bottom(dst, h); /* the bottom field in the tower dst */
689  coeffs bSrc = nCoeff_bottom(src, h); /* the bottom field in the tower src */
690 
691  /* for the time being, we only provide maps if h = 1 and if b is Q or
692     some field Z/pZ: */
693  if (h==0)
694  {
695    if (nCoeff_is_Q(src) && nCoeff_is_Q(bDst))
696      return naMap00;                            /// Q     -->  Q(a)
697    if (nCoeff_is_Zp(src) && nCoeff_is_Q(bDst))
698      return naMapP0;                            /// Z/p   -->  Q(a)
699    if (nCoeff_is_Q(src) && nCoeff_is_Zp(bDst))
700      return naMap0P;                            /// Q      --> Z/p(a)
701    if (nCoeff_is_Zp(src) && nCoeff_is_Zp(bDst))
702    {
703      if (src->ch == dst->ch) return naMapPP;    /// Z/p    --> Z/p(a)
704      else return naMapUP;                       /// Z/u    --> Z/p(a)
705    }
706  }
707  if (h != 1) return NULL;
708  if ((!nCoeff_is_Zp(bDst)) && (!nCoeff_is_Q(bDst))) return NULL;
709  if ((!nCoeff_is_Zp(bSrc)) && (!nCoeff_is_Q(bSrc))) return NULL;
710 
711  if (nCoeff_is_Q(bSrc) && nCoeff_is_Q(bDst))
712  {
713    if (strcmp(rRingVar(0, src->extRing),
714               rRingVar(0, dst->extRing)) == 0)
715    {
716      if (src->type==n_algExt)
717         return ndCopyMap; // naCopyMap;         /// Q(a)   --> Q(a)
718      else
719         return naCopyExt;
720    }
721    else
722      return NULL;                               /// Q(b)   --> Q(a)
723  }
724 
725  if (nCoeff_is_Zp(bSrc) && nCoeff_is_Zp(bDst))
726  {
727    if (strcmp(rRingVar(0,src->extRing),rRingVar(0,dst->extRing))==0)
728    {
729      if (src->type==n_algExt)
730        return ndCopyMap; // naCopyMap;          /// Z/p(a) --> Z/p(a)
731      else
732         return naCopyExt;
733    }
734    else
735      return NULL;                               /// Z/p(b) --> Z/p(a)
736  }
737 
738  return NULL;                                           /// default
739}
740
741int naParDeg(number a, const coeffs cf)
742{
743  if (a == NULL) return -1;
744  poly aa=(poly)a;
745  return cf->extRing->pFDeg(aa,cf->extRing);
746}
747
748BOOLEAN naInitChar(coeffs cf, void * infoStruct)
749{ 
750  assume( infoStruct != NULL );
751
752  AlgExtInfo *e = (AlgExtInfo *)infoStruct;
753  /// first check whether cf->extRing != NULL and delete old ring???
754
755  assume(e->r                     != NULL);      // extRing;
756  assume(e->r->cf                 != NULL);      // extRing->cf;
757  assume((e->i          != NULL) &&    // minideal has one
758         (IDELEMS(e->i) != 0)    &&    // non-zero generator
759         (e->i->m[0]    != NULL)    ); // at m[0];
760
761  assume( e->r->minideal == NULL );
762
763  assume( cf != NULL );
764  assume(getCoeffType(cf) == ID);                     // coeff type;
765 
766  cf->extRing           = e->r;
767  cf->extRing->ref ++; // increase the ref.counter for the ground poly. ring!
768
769  cf->extRing->minideal = e->i; // make a copy?
770
771  /* propagate characteristic up so that it becomes
772     directly accessible in cf: */
773  cf->ch = cf->extRing->cf->ch;
774 
775  #ifdef LDEBUG
776  p_Test((poly)naMinpoly, naRing);
777  #endif
778 
779  cf->cfGreaterZero  = naGreaterZero;
780  cf->cfGreater      = naGreater;
781  cf->cfEqual        = naEqual;
782  cf->cfIsZero       = naIsZero;
783  cf->cfIsOne        = naIsOne;
784  cf->cfIsMOne       = naIsMOne;
785  cf->cfInit         = naInit;
786  cf->cfInit_bigint  = naInit_bigint; 
787  cf->cfInt          = naInt;
788  cf->cfNeg          = naNeg;
789  cf->cfAdd          = naAdd;
790  cf->cfSub          = naSub;
791  cf->cfMult         = naMult;
792  cf->cfDiv          = naDiv;
793  cf->cfExactDiv     = naDiv;
794  cf->cfPower        = naPower;
795  cf->cfCopy         = naCopy;
796  cf->cfWrite        = naWrite;
797  cf->cfRead         = naRead;
798  cf->cfDelete       = naDelete;
799  cf->cfSetMap       = naSetMap;
800  cf->cfGetDenom     = naGetDenom;
801  cf->cfGetNumerator = naGetNumerator;
802  cf->cfRePart       = naCopy;
803  cf->cfImPart       = naImPart;
804  cf->cfCoeffWrite   = naCoeffWrite;
805  cf->cfNormalize    = naNormalize;
806#ifdef LDEBUG
807  cf->cfDBTest       = naDBTest;
808#endif
809  cf->cfGcd          = naGcd;
810  //cf->cfLcm          = naLcm;
811  cf->cfSize         = naSize;
812  cf->nCoeffIsEqual  = naCoeffIsEqual;
813  cf->cfInvers       = naInvers;
814  cf->cfIntDiv       = naDiv;
815#ifdef HAVE_FACTORY
816  cf->convFactoryNSingN=naConvFactoryNSingN;
817  cf->convSingNFactoryN=naConvSingNFactoryN;
818#endif
819  cf->cfParDeg = naParDeg;
820 
821  return FALSE;
822}
823
824
825number naParam(const short iParameter, const coeffs cf)
826{
827  assume(getCoeffType(cf) == ID);
828 
829  const ring R = cf->extRing;
830  assume( R != NULL ); 
831  assume( 0 < iParameter && iParameter <= rVar(R) );
832 
833  poly p = p_One(R); p_SetExp(p, iParameter, 1, R); p_Setm(p, R);
834 
835  return (number) p; 
836}
837
838
839/// if m == var(i)/1 => return i,
840int naIsParam(number m, const coeffs cf)
841{
842  assume(getCoeffType(cf) == ID);
843
844  const ring R = cf->extRing;
845  assume( R != NULL ); 
846
847  return p_Var( (poly)m, R ); 
848}
Note: See TracBrowser for help on using the repository browser.