source: git/libpolys/coeffs/flintcf_Zn.cc @ ce5acd0

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