source: git/libpolys/polys/ext_fields/algext.cc @ 4d94c97

spielwiese
Last change on this file since 4d94c97 was 4d94c97, checked in by Oleksandr Motsak <motsak@…>, 12 years ago
improved n[at]CoeffIsEqual chg: they should test the *equality* of underlying rings, e.g. in order to assure the coeffs of the following rings are the same:Q(P)[X],Q(P)[Y]
  • Property mode set to 100644
File size: 24.3 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  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[0], ii->m[0], 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  int h = 0; /* the height of the extension tower given by dst */
709  coeffs bDst = nCoeff_bottom(dst, h); /* the bottom field in the tower dst */
710  coeffs bSrc = nCoeff_bottom(src, h); /* the bottom field in the tower src */
711 
712  /* for the time being, we only provide maps if h = 1 and if b is Q or
713     some field Z/pZ: */
714  if (h==0)
715  {
716    if (nCoeff_is_Q(src) && nCoeff_is_Q(bDst))
717      return naMap00;                            /// Q     -->  Q(a)
718    if (nCoeff_is_Zp(src) && nCoeff_is_Q(bDst))
719      return naMapP0;                            /// Z/p   -->  Q(a)
720    if (nCoeff_is_Q(src) && nCoeff_is_Zp(bDst))
721      return naMap0P;                            /// Q      --> Z/p(a)
722    if (nCoeff_is_Zp(src) && nCoeff_is_Zp(bDst))
723    {
724      if (src->ch == dst->ch) return naMapPP;    /// Z/p    --> Z/p(a)
725      else return naMapUP;                       /// Z/u    --> Z/p(a)
726    }
727  }
728  if (h != 1) return NULL;
729  if ((!nCoeff_is_Zp(bDst)) && (!nCoeff_is_Q(bDst))) return NULL;
730  if ((!nCoeff_is_Zp(bSrc)) && (!nCoeff_is_Q(bSrc))) return NULL;
731 
732  if (nCoeff_is_Q(bSrc) && nCoeff_is_Q(bDst))
733  {
734    if (strcmp(rRingVar(0, src->extRing),
735               rRingVar(0, dst->extRing)) == 0)
736    {
737      if (src->type==n_algExt)
738         return ndCopyMap; // naCopyMap;         /// Q(a)   --> Q(a)
739      else
740         return naCopyExt;
741    }
742    else
743      return NULL;                               /// Q(b)   --> Q(a)
744  }
745 
746  if (nCoeff_is_Zp(bSrc) && nCoeff_is_Zp(bDst))
747  {
748    if (strcmp(rRingVar(0,src->extRing),rRingVar(0,dst->extRing))==0)
749    {
750      if (src->type==n_algExt)
751        return ndCopyMap; // naCopyMap;          /// Z/p(a) --> Z/p(a)
752      else
753         return naCopyExt;
754    }
755    else
756      return NULL;                               /// Z/p(b) --> Z/p(a)
757  }
758 
759  return NULL;                                           /// default
760}
761
762int naParDeg(number a, const coeffs cf)
763{
764  if (a == NULL) return -1;
765  poly aa=(poly)a;
766  return cf->extRing->pFDeg(aa,cf->extRing);
767}
768
769BOOLEAN naInitChar(coeffs cf, void * infoStruct)
770{ 
771  assume( infoStruct != NULL );
772
773  AlgExtInfo *e = (AlgExtInfo *)infoStruct;
774  /// first check whether cf->extRing != NULL and delete old ring???
775
776  assume(e->r                     != NULL);      // extRing;
777  assume(e->r->cf                 != NULL);      // extRing->cf;
778  assume((e->i          != NULL) &&    // minideal has one
779         (IDELEMS(e->i) != 0)    &&    // non-zero generator
780         (e->i->m[0]    != NULL)    ); // at m[0];
781
782  assume( e->r->minideal == NULL );
783
784  assume( cf != NULL );
785  assume(getCoeffType(cf) == ID);                     // coeff type;
786 
787  cf->extRing           = e->r;
788  cf->extRing->ref ++; // increase the ref.counter for the ground poly. ring!
789
790  cf->extRing->minideal = e->i; // make a copy?
791
792  /* propagate characteristic up so that it becomes
793     directly accessible in cf: */
794  cf->ch = cf->extRing->cf->ch;
795 
796  #ifdef LDEBUG
797  p_Test((poly)naMinpoly, naRing);
798  #endif
799 
800  cf->cfGreaterZero  = naGreaterZero;
801  cf->cfGreater      = naGreater;
802  cf->cfEqual        = naEqual;
803  cf->cfIsZero       = naIsZero;
804  cf->cfIsOne        = naIsOne;
805  cf->cfIsMOne       = naIsMOne;
806  cf->cfInit         = naInit;
807  cf->cfInit_bigint  = naInit_bigint; 
808  cf->cfInt          = naInt;
809  cf->cfNeg          = naNeg;
810  cf->cfAdd          = naAdd;
811  cf->cfSub          = naSub;
812  cf->cfMult         = naMult;
813  cf->cfDiv          = naDiv;
814  cf->cfExactDiv     = naDiv;
815  cf->cfPower        = naPower;
816  cf->cfCopy         = naCopy;
817  cf->cfWrite        = naWrite;
818  cf->cfRead         = naRead;
819  cf->cfDelete       = naDelete;
820  cf->cfSetMap       = naSetMap;
821  cf->cfGetDenom     = naGetDenom;
822  cf->cfGetNumerator = naGetNumerator;
823  cf->cfRePart       = naCopy;
824  cf->cfImPart       = naImPart;
825  cf->cfCoeffWrite   = naCoeffWrite;
826  cf->cfNormalize    = naNormalize;
827#ifdef LDEBUG
828  cf->cfDBTest       = naDBTest;
829#endif
830  cf->cfGcd          = naGcd;
831  //cf->cfLcm          = naLcm;
832  cf->cfSize         = naSize;
833  cf->nCoeffIsEqual  = naCoeffIsEqual;
834  cf->cfInvers       = naInvers;
835  cf->cfIntDiv       = naDiv;
836#ifdef HAVE_FACTORY
837  cf->convFactoryNSingN=naConvFactoryNSingN;
838  cf->convSingNFactoryN=naConvSingNFactoryN;
839#endif
840  cf->cfParDeg = naParDeg;
841 
842  return FALSE;
843}
844
845
846number naParam(const short iParameter, const coeffs cf)
847{
848  assume(getCoeffType(cf) == ID);
849 
850  const ring R = cf->extRing;
851  assume( R != NULL ); 
852  assume( 0 < iParameter && iParameter <= rVar(R) );
853 
854  poly p = p_One(R); p_SetExp(p, iParameter, 1, R); p_Setm(p, R);
855 
856  return (number) p; 
857}
858
859
860/// if m == var(i)/1 => return i,
861int naIsParam(number m, const coeffs cf)
862{
863  assume(getCoeffType(cf) == ID);
864
865  const ring R = cf->extRing;
866  assume( R != NULL ); 
867
868  return p_Var( (poly)m, R ); 
869}
Note: See TracBrowser for help on using the repository browser.