source: git/libpolys/coeffs/flintcf_Zn.cc @ 96694a

spielwiese
Last change on this file since 96694a was 96694a, checked in by Hans Schoenemann <hannes@…>, 5 years ago
fix: sscanf
  • Property mode set to 100644
File size: 14.3 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5* ABSTRACT: flint: nmod_poly_t
6*/
7#include <ctype.h> /* isdigit*/
8
9#include "misc/auxiliary.h"
10
11#ifdef HAVE_FLINT
12
13#include <flint/flint.h>
14#include <flint/nmod_poly.h>
15#include "factory/factory.h"
16
17#include "omalloc/omalloc.h"
18#include "coeffs/coeffs.h"
19
20#include "coeffs/numbers.h"
21#include "coeffs/longrat.h"
22#include "coeffs/modulop.h"
23#include "coeffs/flintcf_Zn.h"
24
25typedef nmod_poly_struct *nmod_poly_ptr;
26
27/*2
28* extracts a long integer from s, returns the rest
29*/
30static const char* Eati(const char *s, int *i)
31{
32
33  if (((*s) >= '0') && ((*s) <= '9'))
34  {
35    unsigned long ii=0L;
36    do
37    {
38      ii *= 10;
39      ii += *s++ - '0';
40    }
41    while (((*s) >= '0') && ((*s) <= '9'));
42    *i=(int)ii;
43  }
44  else (*i) = 1;
45  return s;
46}
47
48static void CoeffWrite(const coeffs r, BOOLEAN details)
49{
50  Print("flint:Z/%d[%s]",r->ch,r->pParameterNames[0]);
51}
52
53static BOOLEAN CoeffIsEqual(const coeffs r, n_coeffType n, void * parameter)
54{
55  flintZn_struct *pp=(flintZn_struct*)parameter;
56  return (r->type==n) &&(r->ch==pp->ch)
57          &&(r->pParameterNames!=NULL)
58          &&(strcmp(r->pParameterNames[0],pp->name)==0);
59}
60static void KillChar(coeffs r)
61{
62  // not yet
63}
64static void SetChar(const coeffs r)
65{
66  // dummy
67}
68static number Mult(number a, number b, const coeffs c)
69{
70  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(*res));
71  nmod_poly_init(res,c->ch);
72  nmod_poly_mul(res,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
73  return (number)res;
74}
75static number Sub(number a, number b, const coeffs c)
76{
77  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
78  nmod_poly_init(res,c->ch);
79  nmod_poly_sub(res,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
80  return (number)res;
81}
82static number Add(number a, number b, const coeffs c)
83{
84  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
85  nmod_poly_init(res,c->ch);
86  nmod_poly_add(res,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
87  return (number)res;
88}
89static number Div(number a, number b, const coeffs c)
90{
91  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
92  nmod_poly_init(res,c->ch);
93  if(nmod_poly_is_zero((nmod_poly_ptr)b))
94  {
95     WerrorS(nDivBy0);
96  }
97  else
98  {
99    nmod_poly_div(res,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
100    nmod_poly_t mod;
101    nmod_poly_init(mod,c->ch);
102    nmod_poly_rem(mod,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
103    if (!nmod_poly_is_zero((nmod_poly_ptr)mod))
104    {
105      WerrorS("cannot divide");
106    }
107    nmod_poly_clear(mod);
108  }
109  return (number)res;
110}
111static number ExactDiv(number a, number b, const coeffs c)
112{
113  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
114  nmod_poly_init(res,c->ch);
115  if(nmod_poly_is_zero((nmod_poly_ptr)b))
116  {
117     WerrorS(nDivBy0);
118  }
119  else
120    nmod_poly_div(res,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
121  return (number)res;
122}
123static number IntMod(number a, number b, const coeffs c)
124{
125  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
126  nmod_poly_init(res,c->ch);
127  nmod_poly_rem(res,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
128  return (number)res;
129}
130static number Init (long i, const coeffs r)
131{
132  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
133  nmod_poly_init(res,r->ch);
134  i= i%r->ch;
135  if (i<0) i+=r->ch;
136  nmod_poly_set_coeff_ui(res,0,i);
137  return (number)res;
138}
139static number InitMPZ (mpz_t i, const coeffs r)
140{
141  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
142  nmod_poly_init(res,r->ch);
143  mpz_t tmp;
144  mpz_init(tmp);
145  slong ii=mpz_mod_ui(tmp,i,r->ch);
146  mpz_clear(tmp);
147  nmod_poly_set_coeff_ui(res,0,ii);
148  return (number)res;
149}
150static int Size (number n, const coeffs r)
151{
152  return nmod_poly_degree((nmod_poly_ptr)n);
153}
154static long Int (number &n, const coeffs r)
155{
156  if (nmod_poly_degree((nmod_poly_ptr)n)==0)
157  {
158    slong m;
159    m=nmod_poly_get_coeff_ui((nmod_poly_ptr)n,0);
160    return (long)m;
161  }
162  return 0;
163}
164static void MPZ(mpz_t result, number &n, const coeffs r)
165{
166  mpz_init(result);
167  if (nmod_poly_degree((nmod_poly_ptr)n)==0)
168  {
169    slong m;
170    m=nmod_poly_get_coeff_ui((nmod_poly_ptr)n,0);
171    mpz_set_ui(result,m);
172  }
173}
174static number Neg(number a, const coeffs r)
175{
176  nmod_poly_neg((nmod_poly_ptr)a,(nmod_poly_ptr)a);
177  return a;
178}
179static number Invers(number a, const coeffs r)
180{
181  if(nmod_poly_is_zero((nmod_poly_ptr)a))
182  {
183    WerrorS(nDivBy0);
184    return NULL;
185  }
186  if (nmod_poly_degree((nmod_poly_ptr)a)==0)
187  {
188    nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
189    nmod_poly_init(res,r->ch);
190    slong c=nmod_poly_get_coeff_ui((nmod_poly_ptr)a,0);
191    extern number  nvInvers (number c, const coeffs r);
192    c=(slong)nvInvers((number)c,r);
193    nmod_poly_set_coeff_ui((nmod_poly_ptr)a,0,c);
194    return (number)res;
195  }
196  else
197  {
198    WerrorS("not invertable");
199    return NULL;
200  }
201}
202static number Copy(number a, const coeffs r)
203{
204 nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
205 nmod_poly_init(res,r->ch);
206 nmod_poly_set(res,(nmod_poly_ptr)a);
207 return (number)res;
208}
209//static number RePart(number a, const coeffs r)
210//{
211//}
212//static number ImPart(number a, const coeffs r)
213//{
214//}
215static BOOLEAN IsOne (number a, const coeffs r);
216static BOOLEAN IsZero (number a, const coeffs r);
217//static void WriteLong(number &a, const coeffs r)
218//{
219//}
220static void WriteShort(number a, const coeffs r)
221{
222  //nmod_poly_print_pretty((nmod_poly_ptr)a,r->pParameterNames[0]);
223  if (IsOne(a,r)) StringAppendS("1");
224  else if (IsZero(a,r)) StringAppendS("0");
225  else
226  {
227    StringAppendS("(");
228    BOOLEAN need_plus=FALSE;
229    for(int i=nmod_poly_length((nmod_poly_ptr)a);i>=0;i--)
230    {
231      slong m=nmod_poly_get_coeff_ui((nmod_poly_ptr)a,i);
232      if (m!=0)
233      {
234        if (need_plus) StringAppendS("+");
235        need_plus=TRUE;
236        if (i>0)
237        {
238          if (m!=1) StringAppend("%d*",(int)m);
239          if (i>1)
240            StringAppend("%s^%d",r->pParameterNames[0],i);
241          else if (i==1)
242            StringAppend("%s",r->pParameterNames[0]);
243        }
244        else StringAppend("%d",(int)m);
245      }
246    }
247    StringAppendS(")");
248  }
249}
250static const char* Read(const char * st, number * a, const coeffs r)
251{
252// we only read "monomials" (i.e. [-][digits][parameter]),
253// everythings else (+,*,^,()) is left to the singular interpreter
254  const char *s=st;
255  *a=(number)omAlloc(sizeof(nmod_poly_t));
256  nmod_poly_init((nmod_poly_ptr)(*a),r->ch);
257  BOOLEAN neg=FALSE;
258  if (*s=='-') { neg=TRUE; s++;}
259  if (isdigit(*s))
260  {
261    int z;
262    s=Eati((char *)s, &z);
263    nmod_poly_set_coeff_ui((nmod_poly_ptr)(*a),0,z);
264  }
265  else if(strncmp(s,r->pParameterNames[0],strlen(r->pParameterNames[0]))==0)
266  {
267    nmod_poly_set_coeff_ui((nmod_poly_ptr)(*a),1,1);
268    s+=strlen(r->pParameterNames[0]);
269    if(isdigit(*s))
270    {
271      int i=1;
272      s=Eati(s,&i);
273      if (i!=1)
274      {
275        nmod_poly_set_coeff_ui((nmod_poly_ptr)(*a),1,0);
276        nmod_poly_set_coeff_ui((nmod_poly_ptr)(*a),i,1);
277      }
278    }
279  }
280  if (neg)
281    nmod_poly_neg((nmod_poly_ptr)(*a),(nmod_poly_ptr)(*a));
282  return s;
283}
284static void Normalize(number &a, const coeffs r)
285{
286}
287static BOOLEAN Greater (number a, number b, const coeffs r)
288{
289  if (nmod_poly_length((nmod_poly_ptr)a)>nmod_poly_length((nmod_poly_ptr)b))
290    return TRUE;
291   else if (nmod_poly_length((nmod_poly_ptr)a)<nmod_poly_length((nmod_poly_ptr)b))
292     return FALSE;
293   for(int i=nmod_poly_length((nmod_poly_ptr)a);i>=0;i--)
294   {
295     slong ac=nmod_poly_get_coeff_ui((nmod_poly_ptr)a,i);
296     slong bc=nmod_poly_get_coeff_ui((nmod_poly_ptr)b,i);
297     if (ac>bc) return TRUE;
298     else if (ac<bc) return FALSE;
299   }
300   return FALSE;
301}
302static BOOLEAN Equal (number a, number b, const coeffs r)
303{
304  return (nmod_poly_equal((nmod_poly_ptr)a,(nmod_poly_ptr)b));
305}
306static BOOLEAN IsZero (number a, const coeffs r)
307{
308  return nmod_poly_is_zero((nmod_poly_ptr)a);
309}
310static BOOLEAN IsOne (number a, const coeffs r)
311{
312  return nmod_poly_is_one((nmod_poly_ptr)a);
313}
314static BOOLEAN IsMOne (number k, const coeffs r)
315{
316  if (nmod_poly_length((nmod_poly_ptr)k)>0) return FALSE;
317  slong m=nmod_poly_get_coeff_ui((nmod_poly_ptr)k,0);
318  return (m+1==r->ch);
319}
320static BOOLEAN GreaterZero (number k, const coeffs r)
321{
322  // does it have a leading sign?
323  // no: 0 and 1 do not have, everything else is in (...)
324  return TRUE;
325}
326static void Power(number a, int i, number * result, const coeffs r)
327{
328  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
329  nmod_poly_init(res,r->ch);
330  *result=(number)res;
331  nmod_poly_pow((nmod_poly_ptr)(*result),(nmod_poly_ptr)a,i);
332}
333static number Gcd(number a, number b, const coeffs r)
334{
335  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
336  nmod_poly_init(res,r->ch);
337  nmod_poly_gcd(res,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
338  return (number)res;
339}
340static number ExtGcd(number a, number b, number *s, number *t,const coeffs r)
341{
342  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
343  nmod_poly_init(res,r->ch);
344  nmod_poly_init((nmod_poly_ptr)*s,r->ch);
345  nmod_poly_init((nmod_poly_ptr)*t,r->ch);
346  nmod_poly_xgcd(res,(nmod_poly_ptr)*s,(nmod_poly_ptr)*t,(nmod_poly_ptr)a,(nmod_poly_ptr)b);
347  return (number)res;
348}
349static number Lcm(number a, number b, const coeffs r)
350{
351  WerrorS("not yet: Lcm");
352}
353static void Delete(number * a, const coeffs r)
354{
355  if ((*a)!=NULL)
356  {
357    nmod_poly_clear((nmod_poly_ptr)*a);
358    omFree(*a);
359    *a=NULL;
360  }
361}
362static nMapFunc SetMap(const coeffs src, const coeffs dst)
363{
364  WerrorS("not yet: SetMap");
365  return NULL;
366}
367//static void InpMult(number &a, number b, const coeffs r)
368//{
369//}
370//static void InpAdd(number &a, number b, const coeffs r)
371//{
372//}
373static number Init_bigint(number i, const coeffs dummy, const coeffs dst)
374{
375  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
376  nmod_poly_init(res,dst->ch);
377  long ii;
378  if (SR_HDL(i) & SR_INT)
379  {
380    ii=SR_TO_INT(i) % dst->ch;
381  }
382  else
383  {
384    mpz_t tmp;
385    mpz_init(tmp);
386    ii=mpz_mod_ui(tmp,i->z,dst->ch);
387    mpz_clear(tmp);
388  }
389  if (ii<0) ii+=dst->ch;
390  nmod_poly_set_coeff_ui(res,0,ii);
391  return (number)res;
392}
393static number Farey(number p, number n, const coeffs)
394{
395  WerrorS("not yet: Farey");
396}
397static number ChineseRemainder(number *x, number *q,int rl, BOOLEAN sym,CFArray &inv_cache,const coeffs)
398{
399  WerrorS("not yet: ChineseRemainder");
400}
401static int ParDeg(number x,const coeffs r)
402{
403  return nmod_poly_degree((nmod_poly_ptr)x);
404}
405static number Parameter(const int i, const coeffs r)
406{
407  nmod_poly_ptr res=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
408  nmod_poly_init(res,r->ch);
409  nmod_poly_set_coeff_ui(res,1,1);
410  return (number)res;
411}
412// cfClearContent
413// cfClearDenominators
414static number ConvFactoryNSingN( const CanonicalForm n, const coeffs r)
415{
416}
417static CanonicalForm ConvSingNFactoryN( number n, BOOLEAN setChar, const coeffs r )
418{
419  WerrorS("not yet: ConvSingNFactoryN");
420}
421static char * CoeffName(const coeffs r)
422{
423  static char CoeffName_flint_Zn[20];
424  sprintf(CoeffName_flint_Zn,"flint:Z/%d[%s]",r->ch,r->pParameterNames[0]);
425  return (char*)CoeffName_flint_Zn;
426}
427static char* CoeffString(const coeffs r)
428{
429  char *buf=(char*)omAlloc(12+10 /*ch*/+strlen(r->pParameterNames[0]));
430  sprintf(buf,"flintZn(%d,\"%s\")",r->ch,r->pParameterNames[0]);
431  return buf;
432}
433coeffs flintZnInitCfByName(char *s,n_coeffType n)
434{
435  const char start[]="flint:Z/";
436  const int start_len=strlen(start);
437  if (strncmp(s,start,start_len)==0)
438  {
439    s+=start_len;
440    int p;
441    char st[10];
442    int l=sscanf(s,"%d[%s",&p,st);
443    if (l==2)
444    {
445      flintZn_struct info;
446      info.ch=p;
447      while (st[strlen(st)-1]==']') st[strlen(st)-1]='\0';
448      info.name=st;
449      return nInitChar(n,(void*)&info);
450    }
451  }
452  return NULL;
453}
454static void WriteFd(number a, const ssiInfo *d, const coeffs)
455{
456  // format: len a_len .. a_0
457  nmod_poly_ptr aa=(nmod_poly_ptr)a;
458  int l=nmod_poly_length(aa);
459  fprintf(d->f_write,"%d ",l);
460  for(int i=l; i>=0; i--)
461  {
462    ulong ul=nmod_poly_get_coeff_ui(aa,i);
463    fprintf(d->f_write,"%lu ", ul);
464  }
465}
466
467static number ReadFd(const ssiInfo *d, const coeffs r)
468{
469  // format: len a_len .. a_0
470  nmod_poly_ptr aa=(nmod_poly_ptr)omAlloc(sizeof(nmod_poly_t));
471  nmod_poly_init(aa,r->ch);
472  int l=s_readint(d->f_read);
473  unsigned long ul;
474  for (int i=l;i>=0;i--)
475  {
476    unsigned long ul=s_readlong(d->f_read);
477    nmod_poly_set_coeff_ui(aa,i,ul);
478  }
479  return (number)aa;
480}
481#ifdef LDEBUG
482static BOOLEAN DBTest(number a, const char *f, const int l, const coeffs r)
483{
484  return TRUE;
485}
486#endif
487BOOLEAN flintZn_InitChar(coeffs cf, void * infoStruct)
488{
489  flintZn_struct *pp=(flintZn_struct*)infoStruct;
490  cf->ch=pp->ch;
491
492  cf->cfCoeffString  = CoeffString;
493  cf->cfCoeffName    = CoeffName;
494  cf->cfCoeffWrite   = CoeffWrite;
495  cf->nCoeffIsEqual  = CoeffIsEqual;
496  cf->cfKillChar     = KillChar;
497  cf->cfSetChar      = SetChar;
498  cf->cfMult         = Mult;
499  cf->cfSub          = Sub;
500  cf->cfAdd          = Add;
501  cf->cfDiv          = Div;
502  cf->cfExactDiv     = ExactDiv; // ???
503  cf->cfInit         = Init;
504  cf->cfInitMPZ      = InitMPZ;
505  cf->cfSize         = Size;
506  cf->cfInt          = Int;
507  cf->cfMPZ          = MPZ;
508  cf->cfInpNeg       = Neg;
509  cf->cfInvers       = Invers;
510  cf->cfCopy         = Copy;
511  cf->cfRePart       = Copy;
512  // default: cf->cfImPart       = ndReturn0;
513  cf->cfWriteLong    = WriteShort; //WriteLong;
514  cf->cfWriteShort = WriteShort;
515  cf->cfRead         = Read;
516  cf->cfNormalize    = Normalize;
517
518  //cf->cfDivComp=
519  //cf->cfIsUnit=
520  //cf->cfGetUnit=
521  //cf->cfDivBy=
522
523  cf->cfGreater=Greater;
524  cf->cfEqual  =Equal;
525  cf->cfIsZero =IsZero;
526  cf->cfIsOne  =IsOne;
527  cf->cfIsMOne =IsMOne;
528  cf->cfGreaterZero=GreaterZero;
529
530  cf->cfPower        = Power;
531  //default: cf->cfGetDenom     = GetDenom;
532  //default: cf->cfGetNumerator = GetNumerator;
533  cf->cfGcd          = Gcd;
534  cf->cfExtGcd         = ExtGcd;
535  cf->cfLcm          = Lcm;
536  cf->cfDelete       = Delete;
537  cf->cfSetMap       = SetMap;
538  // default: cf->cfInpMult
539  // default: cf->cfInpAdd
540  cf->cfFarey        =Farey;
541  cf->cfChineseRemainder=ChineseRemainder;
542  cf->cfParDeg = ParDeg;
543  cf->cfParameter = Parameter;
544  //  cf->cfClearContent = ClearContent;
545  //  cf->cfClearDenominators = ClearDenominators;
546  cf->convFactoryNSingN=ConvFactoryNSingN;
547  cf->convSingNFactoryN=ConvSingNFactoryN;
548  cf->cfWriteFd = WriteFd;
549  cf->cfReadFd  = ReadFd;
550#ifdef LDEBUG
551  cf->cfDBTest       = DBTest;
552#endif
553
554  cf->iNumberOfParameters = 1;
555  char **pn=(char**)omAlloc0(sizeof(char*));
556  pn[0]=(char*)omStrDup(pp->name);
557  cf->pParameterNames = (const char **)pn;
558  cf->has_simple_Inverse= FALSE;
559  cf->has_simple_Alloc= FALSE;
560  cf->is_field=FALSE;
561
562  return FALSE;
563}
564#endif
Note: See TracBrowser for help on using the repository browser.