source: git/kernel/numbers.cc @ 994445

spielwiese
Last change on this file since 994445 was 994445, checked in by Oliver Wienand <wienand@…>, 16 years ago
Neuer Koeffizientenbereich Z/n Neuer Makrodatentyp NATNUMBER = unsigned long Bisschen aufgeraeumt git-svn-id: file:///usr/local/Singular/svn/trunk@10015 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 15.1 KB
Line 
1/*****************************************
2*  Computer Algebra System SINGULAR      *
3*****************************************/
4/* $Id: numbers.cc,v 1.6 2007-05-03 13:50:09 wienand Exp $ */
5
6/*
7* ABSTRACT: interface to coefficient aritmetics
8*/
9
10#include <string.h>
11#include <stdlib.h>
12#include "mod2.h"
13#include "structs.h"
14#include "febase.h"
15#include "kstd1.h"
16#include "numbers.h"
17#include "longrat.h"
18#include "longalg.h"
19#include "modulop.h"
20#include "gnumpfl.h"
21#include "gnumpc.h"
22#include "ring.h"
23#include "ffields.h"
24#include "shortfl.h"
25#ifdef HAVE_RING2TOM
26#include "rmodulo2m.h"
27#endif
28#ifdef HAVE_RINGMODN
29#include "rmodulon.h"
30#endif
31
32//static int characteristic = 0;
33extern int IsPrime(int p);
34
35void   (*nNew)(number *a);
36number (*nInit)(int i);
37number (*nPar)(int i);
38int    (*nParDeg)(number n);
39int    (*nSize)(number n);
40int    (*nInt)(number &n);
41numberfunc nMult, nSub, nAdd, nDiv, nIntDiv, nIntMod, nExactDiv;
42number (*nNeg)(number a);
43number (*nInvers)(number a);
44void   (*nNormalize)(number &a);
45number (*nCopy)(number a);
46number (*nRePart)(number a);
47number (*nImPart)(number a);
48BOOLEAN (*nGreater)(number a,number b);
49BOOLEAN (*nEqual)(number a,number b);
50BOOLEAN (*nIsZero)(number a);
51BOOLEAN (*nIsOne)(number a);
52BOOLEAN (*nIsMOne)(number a);
53BOOLEAN (*nGreaterZero)(number a);
54void    (*nWrite)(number &a);
55char *  (*nRead)(char *s,number *a);
56void    (*nPower)(number a, int i, number * result);
57number  (*nGcd)(number a, number b, const ring r);
58number  (*nLcm)(number a, number b, const ring r);
59char * (*nName)(number n);
60void   (*n__Delete)(number *a, const ring r);
61
62/*0 implementation*/
63number nNULL; /* the 0 as constant */
64
65
66n_Procs_s *cf_root=NULL;
67
68void   nDummy1(number* d) { *d=NULL; }
69void   ndDelete(number* d, const ring r) { *d=NULL; }
70void   ndInpMult(number &a, number b, const ring r)
71{
72  number n=n_Mult(a,b,r);
73  n_Delete(&a,r);
74  a=n;
75}
76
77#ifdef LDEBUG
78void   nDBDummy1(number* d,char *f, int l) { *d=NULL; }
79#endif
80
81void   nDummy2(number& d) { }
82
83char * ndName(number n) { return NULL; }
84
85number ndPar(int i) { return nInit(0); }
86
87number ndReturn0(number n) { return nInit(0); }
88
89int    ndParDeg(number n) { return 0; }
90
91number ndGcd(number a, number b, const ring r) { return r->cf->nInit(1); }
92
93number ndIntMod(number a, number b) { return nInit(0); }
94
95number ndGetDenom(number &n, const ring r) { return n_Init(1,r); }
96
97int ndSize(number a) { return (int)nIsZero(a)==FALSE; }
98
99number ndCopy(number a) { return a; }
100number nd_Copy(number a,const ring r) { return r->cf->nCopy(a); }
101
102/*2
103* init operations for characteristic c (complete==TRUE)
104* init nDelete    for characteristic c (complete==FALSE)
105*/
106void nSetChar(ring r)
107{
108  int c=rInternalChar(r);
109
110  n__Delete= r->cf->cfDelete;
111  if (rField_is_Extension(r))
112  {
113    naSetChar(c,r);
114  }
115#ifdef HAVE_RING2TOM
116  /*----------------------ring Z / 2^m----------------*/
117  else if (rField_is_Ring_2toM(r))
118  {
119    nr2mSetExp(c, r);
120  }
121#endif 
122#ifdef HAVE_RINGMODN
123  /*----------------------ring Z / n----------------*/
124  else if (rField_is_Ring_ModN(r))
125  {
126    nrnSetExp(c, r);
127  }
128#endif
129  else if (rField_is_Zp(r))
130  /*----------------------char. p----------------*/
131  {
132    npSetChar(c, r);
133  }
134  /* -------------- GF(p^m) -----------------------*/
135  else if (rField_is_GF(r))
136  {
137    nfSetChar(c,r->parameter);
138  }
139  /* -------------- R -----------------------*/
140  //if (c==(-1))
141  else if (rField_is_R(r))
142  {
143  }
144  /* -------------- long R -----------------------*/
145  /* -------------- long C -----------------------*/
146  else if ((rField_is_long_R(r))
147  || (rField_is_long_C(r)))
148  {
149    setGMPFloatDigits(r->float_len,r->float_len2);
150  }
151#ifdef TEST
152  /* -------------- R -----------------------*/
153  //if (c==(-1))
154  else if (!rField_is_R(r) && !rField_is_Q(r))
155  {
156    WerrorS("unknown field");
157  }
158#endif
159  nNew   = r->cf->nNew;
160  nNormalize=r->cf->nNormalize;
161  nInit  = r->cf->nInit;
162  nPar   = r->cf->nPar;
163  nParDeg= r->cf->nParDeg;
164  nInt   = r->cf->nInt;
165  nAdd   = r->cf->nAdd;
166  nSub   = r->cf->nSub;
167  nMult  = r->cf->nMult;
168  nDiv   = r->cf->nDiv;
169  nExactDiv= r->cf->nExactDiv;
170  nIntDiv= r->cf->nIntDiv;
171  nIntMod= r->cf->nIntMod;
172  nNeg   = r->cf->nNeg;
173  nInvers= r->cf->nInvers;
174  nCopy  = r->cf->nCopy;
175  nGreater = r->cf->nGreater;
176  nEqual = r->cf->nEqual;
177  nIsZero = r->cf->nIsZero;
178  nIsOne = r->cf->nIsOne;
179  nIsMOne = r->cf->nIsMOne;
180  nGreaterZero = r->cf->nGreaterZero;
181  nWrite = r->cf->nWrite;
182  nRead = r->cf->nRead;
183  nPower = r->cf->nPower;
184  nGcd  = r->cf->nGcd;
185  nLcm  = r->cf->nLcm;
186  nName= r->cf->nName;
187  nSize  = r->cf->nSize;
188  nRePart = r->cf->nRePart;
189  nImPart = r->cf->nImPart;
190  nNULL=r->cf->nNULL;
191}
192
193/*2
194* init operations for ring r
195*/
196void nInitChar(ring r)
197{
198  int c=rInternalChar(r);
199  n_coeffType t=rFieldType(r);
200
201  if (rField_is_Extension(r))
202  {
203    if (r->algring==NULL)
204    {
205      int ch=-c;
206      if (c==1) ch=0;
207      r->algring=(ring) rDefault(ch,r->P,r->parameter);
208      //r->algring->ShortOut=r->ShortOut;
209      // includes: nInitChar(r->algring);
210    }
211  }
212
213  n_Procs_s *n=cf_root;
214  while((n!=NULL)
215    && ((n->nChar!=c) || (n->type!=t)))
216      n=n->next;
217  if (n==NULL)
218  {
219    n=(n_Procs_s*)omAlloc0(sizeof(n_Procs_s));
220    n->next=cf_root;
221    n->ref=1;
222    n->nChar=c;
223    n->type=t;
224    cf_root=n;
225  }
226  else if ((n->nChar==c) && (n->type==t))
227  {
228    n->ref++;
229    r->cf=n;
230    return;
231  }
232  else
233  {
234    WerrorS("nInitChar failed");
235    return;
236  }
237  r->cf=n;
238  n->nChar = c;
239  n->nPar  = ndPar;
240  n->nParDeg=ndParDeg;
241  n->nSize = ndSize;
242  n->n_GetDenom= ndGetDenom;
243  n->nName =  ndName;
244  n->nImPart=ndReturn0;
245  n->cfDelete= ndDelete;
246  n->nNew=nDummy1;
247  n->nInpMult=ndInpMult;
248  n->cfCopy=nd_Copy;
249  n->nIntMod=ndIntMod; /* dummy !! */
250  n->nNormalize=nDummy2;
251  n->nGcd  = ndGcd;
252  n->nLcm  = ndGcd; /* tricky, isn't it ?*/
253  if (rField_is_Extension(r))
254  {
255    //naInitChar(c,TRUE,r);
256    n->cfDelete = naDelete;
257    n-> nNew       = naNew;
258    n-> nNormalize = naNormalize;
259    n->nInit       = naInit;
260    n->nPar        = naPar;
261    n->nParDeg     = naParDeg;
262    n->nInt        = naInt;
263    n->nAdd        = naAdd;
264    n->nSub        = naSub;
265    n->nMult       = naMult;
266    n->nDiv        = naDiv;
267    n->nExactDiv   = naDiv;
268    n->nIntDiv     = naIntDiv;
269    n->nNeg        = naNeg;
270    n->nInvers     = naInvers;
271    n->nCopy       = naCopy;
272    n->cfCopy      = na_Copy;
273    n->nGreater    = naGreater;
274    n->nEqual      = naEqual;
275    n->nIsZero     = naIsZero;
276    n->nIsOne      = naIsOne;
277    n->nIsMOne     = naIsMOne;
278    n->nGreaterZero= naGreaterZero;
279    n->nWrite      = naWrite;
280    n->nRead       = naRead;
281    n->nPower      = naPower;
282    n->nGcd        = naGcd;
283    n->nLcm        = naLcm;
284    n->cfSetMap    = naSetMap;
285    n->nName       = naName;
286    n->nSize       = naSize;
287    n->n_GetDenom   = naGetDenom;
288#ifdef LDEBUG
289    //n->nDBTest     = naDBTest;
290#endif
291  }
292#ifdef HAVE_RING2TOM
293  /* -------------- Z/2^m ----------------------- */
294  else if (rField_is_Ring_2toM(r))
295  {
296     nr2mInitExp(c,r);
297     n->nInit  = nr2mInit;
298     n->nCopy  = ndCopy;
299     n->nInt   = nr2mInt;
300     n->nAdd   = nr2mAdd;
301     n->nSub   = nr2mSub;
302     n->nMult  = nr2mMult;
303     n->nDiv   = nr2mDiv;
304     n->nIntDiv       = nr2mIntDiv;
305     n->nExactDiv= nr2mDiv;
306     n->nNeg   = nr2mNeg;
307     n->nInvers= nr2mInvers;
308     n->nGreater = nr2mGreater;
309     n->nEqual = nr2mEqual;
310     n->nIsZero = nr2mIsZero;
311     n->nIsOne = nr2mIsOne;
312     n->nIsMOne = nr2mIsMOne;
313     n->nGreaterZero = nr2mGreaterZero;
314     n->nWrite = nr2mWrite;
315     n->nRead = nr2mRead;
316     n->nPower = nr2mPower;
317     n->cfSetMap = nr2mSetMap;
318     n->nNormalize = nDummy2;
319     n->nLcm          = nr2mLcm;
320     n->nGcd          = nr2mGcd;
321//     n->nGetUnit = nr2mGetUnit; //TODO OLIVER
322     n->nName= ndName;
323#ifdef LDEBUG
324//     n->nDBTest=nr2mDBTest;
325#endif
326  }
327#endif
328#ifdef HAVE_RINGMODN
329  /* -------------- Z/n ----------------------- */
330  else if (rField_is_Ring_ModN(r))
331  {
332     nrnInitExp(c,r);
333     n->nInit  = nrnInit;
334     n->nCopy  = ndCopy;
335     n->nInt   = nrnInt;
336     n->nAdd   = nrnAdd;
337     n->nSub   = nrnSub;
338     n->nMult  = nrnMult;
339     n->nDiv   = nrnDiv;
340     n->nIntDiv       = nrnIntDiv;
341     n->nExactDiv= nrnDiv;
342     n->nNeg   = nrnNeg;
343     n->nInvers= nrnInvers;
344     n->nGreater = nrnGreater;
345     n->nEqual = nrnEqual;
346     n->nIsZero = nrnIsZero;
347     n->nIsOne = nrnIsOne;
348     n->nIsMOne = nrnIsMOne;
349     n->nGreaterZero = nrnGreaterZero;
350     n->nWrite = nrnWrite;
351     n->nRead = nrnRead;
352     n->nPower = nrnPower;
353     n->cfSetMap = nrnSetMap;
354     n->nNormalize = nDummy2;
355     n->nLcm          = nrnLcm;
356     n->nGcd          = nrnGcd;
357//     n->nGetUnit = nrnGetUnit; //TODO OLIVER
358     n->nName= ndName;
359#ifdef LDEBUG
360//     n->nDBTest=nrnDBTest;
361#endif
362  }
363#endif
364  else if (rField_is_Q(r))
365  {
366    n->cfDelete= nlDelete;
367    n->nNew   = nlNew;
368    n->nNormalize=nlNormalize;
369    n->nInit  = nlInit;
370    n->nInt   = nlInt;
371    n->nAdd   = nlAdd;
372    n->nSub   = nlSub;
373    n->nMult  = nlMult;
374    n->nInpMult=nlInpMult;
375    n->nDiv   = nlDiv;
376    n->nExactDiv= nlExactDiv;
377    n->nIntDiv= nlIntDiv;
378    n->nIntMod= nlIntMod;
379    n->nNeg   = nlNeg;
380    n->nInvers= nlInvers;
381    n->nCopy  = nlCopy;
382    n->nGreater = nlGreater;
383    n->nEqual = nlEqual;
384    n->nIsZero = nlIsZero;
385    n->nIsOne = nlIsOne;
386    n->nIsMOne = nlIsMOne;
387    n->nGreaterZero = nlGreaterZero;
388    n->nWrite = nlWrite;
389    n->nRead = nlRead;
390    n->nPower = nlPower;
391    n->nGcd  = nlGcd;
392    n->nLcm  = nlLcm;
393    n->cfSetMap = nlSetMap;
394    n->nSize  = nlSize;
395    n->n_GetDenom = nlGetDenom;
396#ifdef LDEBUG
397    //n->nDBTest=nlDBTest;
398#endif
399  }
400  else if (rField_is_Zp(r))
401  /*----------------------char. p----------------*/
402  {
403    npInitChar(c,r);
404    n->nInit  = npInit;
405    n->nInt   = npInt;
406    n->nAdd   = npAdd;
407    n->nSub   = npSub;
408    n->nMult  = npMult;
409    n->nDiv   = npDiv;
410    n->nExactDiv= npDiv;
411    n->nNeg   = npNeg;
412    n->nInvers= npInvers;
413    n->nCopy  = ndCopy;
414    n->nGreater = npGreater;
415    n->nEqual = npEqual;
416    n->nIsZero = npIsZero;
417    n->nIsOne = npIsOne;
418    n->nIsMOne = npIsMOne;
419    n->nGreaterZero = npGreaterZero;
420    n->nWrite = npWrite;
421    n->nRead = npRead;
422    n->nPower = npPower;
423    n->cfSetMap = npSetMap;
424    /* nName= ndName; */
425    /*nSize  = ndSize;*/
426#ifdef LDEBUG
427    //n->nDBTest=npDBTest;
428#endif
429#ifdef NV_OPS
430    if (c>NV_MAX_PRIME)
431    {
432      n->nMult  = nvMult;
433      n->nDiv   = nvDiv;
434      n->nExactDiv= nvDiv;
435      n->nInvers= nvInvers;
436    }
437#endif
438  }
439  /* -------------- GF(p^m) -----------------------*/
440  else if (rField_is_GF(r))
441  {
442    //nfSetChar(c,r->parameter);
443    n->nInit  = nfInit;
444    n->nPar   = nfPar;
445    n->nParDeg= nfParDeg;
446    n->nInt   = nfInt;
447    n->nAdd   = nfAdd;
448    n->nSub   = nfSub;
449    n->nMult  = nfMult;
450    n->nDiv   = nfDiv;
451    n->nExactDiv= nfDiv;
452    n->nNeg   = nfNeg;
453    n->nInvers= nfInvers;
454    n->nCopy  = ndCopy;
455    n->nGreater = nfGreater;
456    n->nEqual = nfEqual;
457    n->nIsZero = nfIsZero;
458    n->nIsOne = nfIsOne;
459    n->nIsMOne = nfIsMOne;
460    n->nGreaterZero = nfGreaterZero;
461    n->nWrite = nfWrite;
462    n->nRead = nfRead;
463    n->nPower = nfPower;
464    n->cfSetMap = nfSetMap;
465    n->nName= nfName;
466    /*nSize  = ndSize;*/
467#ifdef LDEBUG
468    //n->nDBTest=nfDBTest;
469#endif
470  }
471  /* -------------- R -----------------------*/
472  //if (c==(-1))
473  else if (rField_is_R(r))
474  {
475    n->nInit  = nrInit;
476    n->nInt   = nrInt;
477    n->nAdd   = nrAdd;
478    n->nSub   = nrSub;
479    n->nMult  = nrMult;
480    n->nDiv   = nrDiv;
481    n->nExactDiv= nrDiv;
482    n->nNeg   = nrNeg;
483    n->nInvers= nrInvers;
484    n->nCopy  = ndCopy;
485    n->nGreater = nrGreater;
486    n->nEqual = nrEqual;
487    n->nIsZero = nrIsZero;
488    n->nIsOne = nrIsOne;
489    n->nIsMOne = nrIsMOne;
490    n->nGreaterZero = nrGreaterZero;
491    n->nWrite = nrWrite;
492    n->nRead = nrRead;
493    n->nPower = nrPower;
494    n->cfSetMap=nrSetMap;
495    /* nName= ndName; */
496    /*nSize  = ndSize;*/
497#ifdef LDEBUG
498    //n->nDBTest=nrDBTest;
499#endif
500  }
501  /* -------------- long R -----------------------*/
502  else if (rField_is_long_R(r))
503  {
504    n->cfDelete= ngfDelete;
505    n->nNew=ngfNew;
506    n->nInit  = ngfInit;
507    n->nInt   = ngfInt;
508    n->nAdd   = ngfAdd;
509    n->nSub   = ngfSub;
510    n->nMult  = ngfMult;
511    n->nDiv   = ngfDiv;
512    n->nExactDiv= ngfDiv;
513    n->nNeg   = ngfNeg;
514    n->nInvers= ngfInvers;
515    n->nCopy  = ngfCopy;
516    n->nGreater = ngfGreater;
517    n->nEqual = ngfEqual;
518    n->nIsZero = ngfIsZero;
519    n->nIsOne = ngfIsOne;
520    n->nIsMOne = ngfIsMOne;
521    n->nGreaterZero = ngfGreaterZero;
522    n->nWrite = ngfWrite;
523    n->nRead = ngfRead;
524    n->nPower = ngfPower;
525    n->cfSetMap=ngfSetMap;
526    n->nName= ndName;
527    n->nSize  = ndSize;
528#ifdef LDEBUG
529    //n->nDBTest=ngfDBTest;
530#endif
531  }
532  /* -------------- long C -----------------------*/
533  else if (rField_is_long_C(r))
534  {
535    n->cfDelete= ngcDelete;
536    n->nNew=ngcNew;
537    n->nNormalize=nDummy2;
538    n->nInit  = ngcInit;
539    n->nInt   = ngcInt;
540    n->nAdd   = ngcAdd;
541    n->nSub   = ngcSub;
542    n->nMult  = ngcMult;
543    n->nDiv   = ngcDiv;
544    n->nExactDiv= ngcDiv;
545    n->nNeg   = ngcNeg;
546    n->nInvers= ngcInvers;
547    n->nCopy  = ngcCopy;
548    n->nGreater = ngcGreater;
549    n->nEqual = ngcEqual;
550    n->nIsZero = ngcIsZero;
551    n->nIsOne = ngcIsOne;
552    n->nIsMOne = ngcIsMOne;
553    n->nGreaterZero = ngcGreaterZero;
554    n->nWrite = ngcWrite;
555    n->nRead = ngcRead;
556    n->nPower = ngcPower;
557    n->cfSetMap=ngcSetMap;
558    n->nPar=ngcPar;
559    n->nRePart=ngcRePart;
560    n->nImPart=ngcImPart;
561    /*nSize  = ndSize;*/
562#ifdef LDEBUG
563    //n->nDBTest=ngcDBTest;
564#endif
565  }
566#ifdef TEST
567  else
568  {
569    WerrorS("unknown field");
570  }
571#endif
572  if (!errorreported)
573  {
574    n->nNULL=n->nInit(0);
575    if (n->nRePart==NULL)
576      n->nRePart=n->nCopy;
577    if (n->nIntDiv==NULL)
578      n->nIntDiv=n->nDiv;
579  }
580}
581
582void nKillChar(ring r)
583{
584  if (r!=NULL)
585  {
586    if (r->cf!=NULL)
587    {
588      r->cf->ref--;
589      if (r->cf->ref<=0)
590      {
591        n_Procs_s tmp;
592        n_Procs_s* n=&tmp;
593        tmp.next=cf_root;
594        while((n->next!=NULL) && (n->next!=r->cf)) n=n->next;
595        if (n->next==r->cf)
596        {
597          n->next=n->next->next;
598          if (cf_root==r->cf) cf_root=n->next;
599          r->cf->cfDelete(&(r->cf->nNULL),r);
600          switch(r->cf->type)
601          {
602            case n_Zp:
603                 #ifdef HAVE_DIV_MOD
604                 if (r->cf->npInvTable!=NULL)
605                 omFreeSize( (ADDRESS)r->cf->npInvTable,
606                             r->cf->npPrimeM*sizeof(CARDINAL) );
607                 #else
608                 if (r->cf->npExpTable!=NULL)
609                 {
610                   omFreeSize( (ADDRESS)r->cf->npExpTable,
611                               r->cf->npPrimeM*sizeof(CARDINAL) );
612                   omFreeSize( (ADDRESS)r->cf->npLogTable,
613                               r->cf->npPrimeM*sizeof(CARDINAL) );
614                 }
615                 #endif
616                 break;
617            case n_Zp_a:
618            case n_Q_a:
619                 {
620                   number n=r->minpoly;
621                   if (n!=NULL)
622                   {
623                     r->minpoly=NULL;
624                     if (r==currRing) naMinimalPoly=NULL;
625                     naDelete(&n,r);
626                   }
627                 }
628                 break;
629
630            default:
631                 break;
632          }
633          omFreeSize((ADDRESS)r->cf, sizeof(n_Procs_s));
634          r->cf=NULL;
635        }
636        else
637        {
638          WarnS("cf_root list destroyed");
639        }
640      }
641    }
642    if (r->algring!=NULL)
643    {
644      rKill(r->algring);
645      r->algring=NULL;
646    }
647  }
648}
Note: See TracBrowser for help on using the repository browser.