source: git/libpolys/coeffs/numbers.cc @ 98474f

jengelh-datetimespielwiese
Last change on this file since 98474f was 98474f, checked in by Oleksandr Motsak <motsak@…>, 11 years ago
Abstract iterator/enumerator interfaces + dummies for ClearContent and ClearDenominators
  • Property mode set to 100644
File size: 12.6 KB
Line 
1/*****************************************
2*  Computer Algebra System SINGULAR      *
3*****************************************/
4
5/*
6* ABSTRACT: interface to coefficient aritmetics
7*/
8
9#include <string.h>
10#include <stdlib.h>
11
12#include "config.h"
13#include <misc/auxiliary.h>
14
15#ifdef HAVE_FACTORY
16#include <factory/factory.h>
17#endif
18
19
20#include "coeffs.h"
21#include <coeffs/numbers.h>
22
23#include <reporter/reporter.h>
24#include <omalloc/omalloc.h>
25#include <coeffs/numbers.h>
26#include <coeffs/longrat.h>
27#include <coeffs/modulop.h>
28#include <coeffs/gnumpfl.h>
29#include <coeffs/gnumpc.h>
30#include <coeffs/ffields.h>
31#include <coeffs/shortfl.h>
32
33#ifdef HAVE_RINGS
34#include <coeffs/rmodulo2m.h>
35#include <coeffs/rmodulon.h>
36#include <coeffs/rintegers.h>
37#endif
38
39#ifdef HAVE_POLYEXTENSIONS
40#include <polys/ext_fields/algext.h>
41#include <polys/ext_fields/transext.h>
42#endif
43
44
45
46//static int characteristic = 0;
47extern int IsPrime(int p);
48
49n_Procs_s *cf_root=NULL;
50
51void   nNew(number* d) { *d=NULL; }
52void   ndDelete(number* d, const coeffs) { *d=NULL; }
53void   ndInpMult(number &a, number b, const coeffs r)
54{
55  number n=n_Mult(a,b,r);
56  n_Delete(&a,r);
57  a=n;
58}
59void ndInpAdd(number &a, number b, const coeffs r)
60{
61  number n=n_Add(a,b,r);
62  n_Delete(&a,r);
63  a=n;
64}
65
66#ifdef LDEBUG
67void   nDBDummy1(number* d,char *, int) { *d=NULL; }
68BOOLEAN ndDBTest(number, const char *, const int, const coeffs)
69{
70  return TRUE;
71}
72#endif
73
74number ndFarey(number,number,const coeffs r)
75{
76  Werror("farey not implemented for (c=%d)",getCoeffType(r));
77  return NULL;
78}
79number ndChineseRemainder(number *,number *,int,const coeffs r)
80{
81  Werror("ChineseRemainder not implemented for (c=%d)",getCoeffType(r));
82  return n_Init(0,r); 
83}
84
85static int ndParDeg(number n, const coeffs r)
86{
87  return (-n_IsZero(n,r));
88}
89
90static number ndParameter(const int, const coeffs r)
91{
92  Werror("ndParameter: n_Parameter is not implemented/relevant for (coeff_type = %d)",getCoeffType(r));
93  return NULL;
94}
95
96BOOLEAN n_IsZeroDivisor( number a, const coeffs r)
97{
98  int c = n_GetChar(r);
99  BOOLEAN ret = n_IsZero(a, r);
100  if( (c != 0) && !ret )
101  {
102    number ch = n_Init( c, r ); 
103    number g = n_Gcd( ch, a, r );
104    ret = !n_IsOne (g, r);
105    n_Delete(&ch, r);
106    n_Delete(&g, r);
107  }
108  return ret; 
109}
110
111void   ndNormalize(number&, const coeffs) { }
112
113char * ndName(number, const coeffs) { return NULL; }
114
115number ndReturn0(number, const coeffs r) { return n_Init(0,r); }
116
117number ndGcd(number, number, const coeffs r) { return n_Init(1,r); }
118
119number ndIntMod(number, number, const coeffs r) { return n_Init(0,r); }
120
121number ndGetDenom(number &, const coeffs r) { return n_Init(1,r); }
122number ndGetNumerator(number &a,const coeffs r) { return n_Copy(a,r); }
123
124int ndSize(number a, const coeffs r) { return (int)n_IsZero(a,r)==FALSE; }
125
126/// divide by the first (leading) number and return it, i.e. make monic
127static void ndClearContent(ICoeffsEnumerator& numberCollectionEnumerator, number& c, const coeffs r)
128{
129  assume(r != NULL);
130
131  // no fractions
132  assume(!(  nCoeff_is_Q(r) || nCoeff_is_Q_a(r) ));
133  // all coeffs are given by integers!!!
134  assume( nCoeff_is_Ring(r) || nCoeff_is_Zp(r) || nCoeff_is_numeric(r) || nCoeff_is_GF(r) || nCoeff_is_Zp_a(r) );
135
136  numberCollectionEnumerator.Reset();
137 
138#ifdef HAVE_RINGS
139  /// TODO: move to a separate implementation
140  if (nCoeff_is_Ring(r))
141  {
142    if (nCoeff_has_Units(r))
143    {
144      c = n_GetUnit(numberCollectionEnumerator.Current(), r);
145     
146      if (!n_IsOne(c, r))
147      {
148        number inv = n_Invers(c, r);
149
150        do
151        {
152          n_InpMult(numberCollectionEnumerator.Current(), inv, r);
153        }
154        while( numberCollectionEnumerator.MoveNext() );
155
156        n_Delete(&inv, r);       
157      }     
158    } else c = n_Init(1, r);
159   
160    return;
161  }
162#endif
163
164  assume(!nCoeff_is_Ring(r));
165
166  c = numberCollectionEnumerator.Current();
167 
168  n_Normalize(c, r);
169
170  if (!n_IsOne(c, r))
171  {   
172    numberCollectionEnumerator.Current() = n_Init(1, r); // ???
173   
174    number inv = n_Invers(c, r);
175   
176    while( numberCollectionEnumerator.MoveNext() )
177    {
178      number &n = numberCollectionEnumerator.Current();
179      n_Normalize(n, r); // ?
180      n_InpMult(n, inv, r);
181    }
182   
183    n_Delete(&inv, r);
184  }
185}
186
187/// does nothing (just returns a dummy one number)
188static void ndClearDenominators(ICoeffsEnumerator& /*numberCollectionEnumerator*/, number& d, const coeffs r)
189{
190  assume( r != NULL );
191  assume( !(nCoeff_is_Q(r) || nCoeff_is_transExt(r)) );
192  assume( nCoeff_is_Ring(r) || nCoeff_is_Zp(r) || nCoeff_is_numeric(r) || nCoeff_is_GF(r) || (getCoeffType(r)==n_algExt) );
193
194  d = n_Init(1, r);
195}
196
197number ndCopy(number a, const coeffs) { return a; }
198number ndCopyMap(number a, const coeffs aRing, const coeffs r)
199{
200  assume( getCoeffType(r) == getCoeffType(aRing) );
201  if ( nCoeff_has_simple_Alloc(r) && nCoeff_has_simple_Alloc(aRing) )
202    return a;
203        else
204    return n_Copy(a, r);
205}
206void ndKillChar(coeffs) {}
207void ndSetChar(const coeffs) {}
208
209number nd_Copy(number a, const coeffs r) { return n_Copy(a, r); }
210
211#ifdef HAVE_RINGS
212BOOLEAN ndDivBy(number, number, const coeffs) { return TRUE; } // assume a,b !=0
213int ndDivComp(number, number, const coeffs) { return 2; }
214BOOLEAN ndIsUnit(number a, const coeffs r) { return !n_IsZero(a,r); }
215number  ndExtGcd (number, number, number *, number *, const coeffs r) { return n_Init(1,r); }
216#endif
217
218#ifdef HAVE_FACTORY
219CanonicalForm ndConvSingNFactoryN( number, BOOLEAN /*setChar*/, const coeffs)
220{
221  CanonicalForm term(0);
222  Werror("no conversion to factory");
223  return term;
224}
225
226number ndConvFactoryNSingN( const CanonicalForm, const coeffs)
227{
228  Werror("no conversion from factory");
229  return NULL;
230}
231#endif
232
233number  ndInit_bigint(number, const coeffs, const coeffs)
234{
235  Werror("no conversion from bigint to this field");
236  return NULL;
237}
238
239/**< [in, out] a bigint number >= 0  */
240/**< [out] the GMP equivalent    */
241/// Converts a non-negative bigint number into a GMP number.
242void ndMPZ(mpz_t result, number &n, const coeffs r)
243{
244  mpz_init_set_si( result, n_Int(n, r) );
245}
246
247number ndInitMPZ(mpz_t m, const coeffs r)
248{ 
249  return n_Init( mpz_get_si(m), r);
250}
251
252
253BOOLEAN ndCoeffIsEqual(const coeffs r, n_coeffType n, void *)
254{
255  /* test, if r is an instance of nInitCoeffs(n,parameter) */
256  /* if paramater is not needed */
257  return (n==r->type);
258}
259
260static n_coeffType nLastCoeffs=n_CF;
261cfInitCharProc nInitCharTableDefault[]=
262{ NULL,        /*n_unknown */
263 npInitChar,   /* n_Zp */
264 nlInitChar,   /* n_Q */
265 nrInitChar,   /* n_R */
266 nfInitChar,   /* n_GF */
267 ngfInitChar,  /* n_long_R */
268 #ifdef HAVE_POLYEXTENSIONS
269 naInitChar,  /* n_algExt */
270 ntInitChar,  /* n_transExt */
271 #else
272 NULL,        /* n_algExt */
273 NULL,        /* n_transExt */
274 #endif   
275 ngcInitChar,  /* n_long_C */
276 #ifdef HAVE_RINGS
277 nrzInitChar,  /* n_Z */
278 nrnInitChar,  /* n_Zn */
279 NULL,         /* n_Zpn */
280 nr2mInitChar, /* n_Z2m */
281 #else
282 NULL,         /* n_Z */
283 NULL,         /* n_Zn */
284 NULL,         /* n_Zpn */
285 NULL,         /* n_Z2m */
286 #endif
287 NULL   /* n_CF */
288};
289
290static cfInitCharProc *nInitCharTable=nInitCharTableDefault;
291/*2
292* init operations for coeffs r
293*/
294coeffs nInitChar(n_coeffType t, void * parameter)
295{
296  n_Procs_s *n=cf_root;
297
298  while((n!=NULL) && (n->nCoeffIsEqual!=NULL) && (!n->nCoeffIsEqual(n,t,parameter)))
299      n=n->next;
300
301  if (n==NULL)
302  {
303    n=(n_Procs_s*)omAlloc0(sizeof(n_Procs_s));
304    n->next=cf_root;
305    n->ref=1;
306    n->type=t;
307
308    // default entries (different from NULL) for some routines:
309    n->cfSize = ndSize;
310    n->cfGetDenom= ndGetDenom;
311    n->cfGetNumerator= ndGetNumerator;
312    n->cfName =  ndName;
313    n->cfImPart=ndReturn0;
314    n->cfDelete= ndDelete;
315    n->cfInpMult=ndInpMult;
316    n->cfCopy = ndCopy;
317    n->cfIntMod=ndIntMod; /* dummy !! */
318    n->cfNormalize=ndNormalize;
319    n->cfGcd  = ndGcd;
320    n->cfLcm  = ndGcd; /* tricky, isn't it ?*/
321    n->cfInit_bigint = ndInit_bigint;
322    n->cfInitMPZ = ndInitMPZ;
323    n->cfMPZ = ndMPZ;
324
325    //n->cfKillChar = ndKillChar; /* dummy */
326    n->cfSetChar = ndSetChar; /* dummy */
327    // temp. removed to catch all the coeffs which miss to implement this!
328
329    n->cfChineseRemainder = ndChineseRemainder;
330    n->cfFarey = ndFarey;
331    n->cfParDeg = ndParDeg;
332   
333    n->cfParameter = ndParameter;
334
335    n->cfClearContent = ndClearContent;
336    n->cfClearDenominators = ndClearDenominators;
337
338#ifdef HAVE_RINGS
339    n->cfDivComp = ndDivComp;
340    n->cfDivBy = ndDivBy;
341    n->cfIsUnit = ndIsUnit;
342    n->cfExtGcd = ndExtGcd;
343    //n->cfGetUnit = (nMapFunc)NULL;
344#endif
345
346#ifdef fACTORY
347    n->convSingNFactoryN=ndConvSingNFactoryN;
348    n->convFactoryNSingN=ndConvFactoryNSingN;
349#endif
350   
351    BOOLEAN nOK=TRUE;
352    // init
353    if ((t<=nLastCoeffs) && (nInitCharTable[t]!=NULL))
354      nOK = (nInitCharTable[t])(n,parameter);
355    else
356       Werror("Sorry: the coeff type [%d] was not registered: it is missing in nInitCharTable", (int)t);
357    if (nOK)
358    {
359      omFreeSize(n,sizeof(*n));
360      return NULL;
361    }
362    cf_root=n;
363    // post init settings:
364    if (n->cfRePart==NULL) n->cfRePart=n->cfCopy;
365    if (n->cfIntDiv==NULL) n->cfIntDiv=n->cfDiv;
366   
367#ifdef HAVE_RINGS
368    if (n->cfGetUnit==NULL) n->cfGetUnit=n->cfCopy;
369#endif
370
371    if(n->cfWriteShort==NULL)
372      n->cfWriteShort = n->cfWriteLong;
373
374    assume(n->nCoeffIsEqual!=NULL);
375    assume(n->cfSetChar!=NULL);
376    assume(n->cfMult!=NULL);
377    assume(n->cfSub!=NULL);
378    assume(n->cfAdd!=NULL);
379    assume(n->cfDiv!=NULL);
380    assume(n->cfIntDiv!=NULL);
381    assume(n->cfIntMod!=NULL);
382    assume(n->cfExactDiv!=NULL);
383    assume(n->cfInit!=NULL);
384    assume(n->cfInitMPZ!=NULL);
385    assume(n->cfSize!=NULL);
386    assume(n->cfInt!=NULL);
387    assume(n->cfMPZ!=NULL);
388    //assume(n->n->cfDivComp!=NULL);
389    //assume(n->cfIsUnit!=NULL);
390    //assume(n->cfGetUnit!=NULL);
391    //assume(n->cfExtGcd!=NULL);
392    assume(n->cfNeg!=NULL);
393    assume(n->cfCopy!=NULL);
394    assume(n->cfRePart!=NULL);
395    assume(n->cfImPart!=NULL);
396
397    assume(n->cfWriteLong!=NULL);
398    assume(n->cfWriteShort!=NULL);
399
400    assume(n->iNumberOfParameters>= 0);
401
402    assume( (n->iNumberOfParameters == 0 && n->pParameterNames == NULL) ||
403            (n->iNumberOfParameters >  0 && n->pParameterNames != NULL) );           
404
405    assume(n->cfParameter!=NULL);
406    assume(n->cfParDeg!=NULL);
407     
408    assume(n->cfRead!=NULL);
409    assume(n->cfNormalize!=NULL);
410    assume(n->cfGreater!=NULL);
411    //assume(n->cfDivBy!=NULL);
412    assume(n->cfEqual!=NULL);
413    assume(n->cfIsZero!=NULL);
414    assume(n->cfIsOne!=NULL);
415    assume(n->cfIsMOne!=NULL);
416    assume(n->cfGreaterZero!=NULL);
417    assume(n->cfPower!=NULL);
418    assume(n->cfGetDenom!=NULL);
419    assume(n->cfGetNumerator!=NULL);
420    assume(n->cfGcd!=NULL);
421    assume(n->cfLcm!=NULL);
422    assume(n->cfDelete!=NULL);
423    assume(n->cfSetMap!=NULL);
424    assume(n->cfName!=NULL);
425    assume(n->cfInpMult!=NULL);
426//    assume(n->cfInit_bigint!=NULL);
427    assume(n->cfCoeffWrite != NULL);
428
429    assume(n->cfClearContent != NULL);
430    assume(n->cfClearDenominators != NULL);
431   
432#ifdef LDEBUG
433    assume(n->cfDBTest!=NULL);
434#endif
435    assume(n->type==t);
436     
437#ifndef NDEBUG
438    if(n->cfKillChar==NULL) Warn("cfKillChar is NULL for coeff %d",t);
439    if(n->cfWriteLong==NULL) Warn("cfWrite is NULL for coeff %d",t);
440    if(n->cfWriteShort==NULL) Warn("cfWriteShort is NULL for coeff %d",t);
441#endif
442     
443   if( n->nNULL == NULL )
444     n->nNULL = n_Init(0, n); // may still remain NULL
445  }
446  else
447  {
448    n->ref++;
449  }
450  return n;
451}
452
453void nKillChar(coeffs r)
454{
455  if (r!=NULL)
456  {
457    r->ref--;
458    if (r->ref<=0)
459    {
460      n_Procs_s tmp;
461      n_Procs_s* n=&tmp;
462      tmp.next=cf_root;
463      while((n->next!=NULL) && (n->next!=r)) n=n->next;
464      if (n->next==r)
465      {
466        n->next=n->next->next;
467        if (cf_root==r) cf_root=n->next;
468        r->cfDelete(&(r->nNULL),r);
469        if (r->cfKillChar!=NULL) r->cfKillChar(r);
470        omFreeSize((void *)r, sizeof(n_Procs_s));
471        r=NULL;
472      }
473      else
474      {
475        WarnS("cf_root list destroyed");
476      }
477    }
478  }
479}
480
481
482n_coeffType nRegister(n_coeffType n, cfInitCharProc p)
483{
484  if (n==n_unknown)
485  {
486    nLastCoeffs=(n_coeffType)(int(nLastCoeffs)+1);
487    if (nInitCharTable==nInitCharTableDefault)
488    {
489      nInitCharTable=(cfInitCharProc*)omAlloc0(
490                                          nLastCoeffs*sizeof(cfInitCharProc));
491      memcpy(nInitCharTable,nInitCharTableDefault,
492              (nLastCoeffs-1)*sizeof(cfInitCharProc));
493    }
494    else
495    {
496      nInitCharTable=(cfInitCharProc*)omReallocSize(nInitCharTable,
497                                          (((int)nLastCoeffs)-1)*sizeof(cfInitCharProc),
498                                          ((int)nLastCoeffs)*sizeof(cfInitCharProc));
499    }
500
501    nInitCharTable[nLastCoeffs]=p;
502    return nLastCoeffs;
503  }
504  else
505  {
506    if (nInitCharTable[n]!=NULL) Print("coeff %d already initialized\n",n);
507    nInitCharTable[n]=p;
508    return n;
509  }
510}
511
Note: See TracBrowser for help on using the repository browser.