source: git/dyn_modules/pcv.mod @ c1ec9a

spielwiese
Last change on this file since c1ec9a was 6ce030f, checked in by Oleksandr Motsak <motsak@…>, 12 years ago
removal of the $Id$ svn tag from everywhere NOTE: the git SHA1 may be used instead (only on special places) NOTE: the libraries Singular/LIB/*.lib still contain the marker due to our current use of svn
  • Property mode set to 100644
File size: 10.3 KB
Line 
1%{
2/*
3 *
4 *  Test mod fuer modgen
5 */
6
7#include <stdio.h>
8#include <kernel/mod2.h>
9#include <tok.h>
10#include <ipid.h>
11#include <numbers.h>
12#include <polys.h>
13#include <ideals.h>
14#include <lists.h>
15#include <matpol.h>
16#include <febase.h>
17
18static int pcvMaxDegree;
19static int pcvTableSize;
20static int pcvIndexSize;
21static unsigned* pcvTable=NULL;
22static unsigned** pcvIndex=NULL;
23
24#ifndef PCV_H
25#define PCV_H
26
27lists pcvLAddL(lists l1,lists l2);
28lists pcvPMulL(poly p,lists l1);
29BOOLEAN pcvLAddL(leftv res,leftv h);
30BOOLEAN pcvPMulL(leftv res,leftv h);
31int pcvDeg(poly p);
32int pcvMinDeg(poly p);
33int pcvMinDeg(matrix m);
34BOOLEAN pcvMinDeg(leftv res,leftv h);
35void pcvInit(int d);
36void pcvClean();
37int pcvM2N(poly m);
38poly pcvN2M(int n);
39poly pcvP2CV(poly p,int d0,int d1);
40poly pcvCV2P(poly cv,int d0,int d1);
41lists pcvP2CV(lists pl,int d0,int d1);
42ideal pcvP2CV(ideal p,int d0,int d1);
43lists pcvCV2P(lists cvl,int d0,int d1);
44ideal pcvCV2P(ideal cv,int d0,int d1);
45BOOLEAN pcvP2CV(leftv res,leftv h);
46BOOLEAN pcvCV2P(leftv res,leftv h);
47int pcvDim(int d0,int d1);
48BOOLEAN pcvDim(leftv res,leftv h);
49int pcvBasis(lists b,int i,poly m,int d,int n);
50lists pcvBasis(int d0,int d1);
51BOOLEAN pcvBasis(leftv res,leftv h);
52
53#endif
54%}
55
56
57// module="pcv";
58package="pcv";
59
60version="$Id$";
61info="
62LIBRARY: pcv.so  CONVERSION BETWEEN POLYS AND COEF VECTORS
63AUTHOR:  Mathias Schulze, email: mschulze@mathematik.uni-kl.de
64
65 MinDeg(p);      min deg of monomials of poly p
66 P2CV(l,d0,d1);  list of coef vectors from deg d0 to d1 of polys in list l
67 CV2P(l,d0,d1);  list of polys with coef vectors from deg d0 to d1 in list l
68 Dim(d0,d1);     number of monomials from deg d0 to d1
69 Basis(d0,d1);   list of monomials from deg d0 to d1
70";
71
72%modinitial
73// no commands to be run upon loading the module
74%endinitial
75
76%procedures
77
78
79int MinDeg(poly p) {
80   %declaration;
81   %typecheck;
82   %return(pcvMinDeg(p));
83}
84
85list P2CV(list pl,int d0,int d1) {
86  %declaration;
87
88  /* check if current RingHandle is set */
89  if(currRingHdl == NULL)
90  {
91    WerrorS("no ring active");
92    return TRUE;
93  }
94 
95  %typecheck;
96  %return(pcvP2CV(pl, d0, d1));
97}
98
99list CV2P(list pl,int d0,int d1)
100{
101  %declaration;
102
103  /* check if current RingHandle is set */
104  if(currRingHdl == NULL)
105  {
106    WerrorS("no ring active");
107    return TRUE;
108  }
109 
110  %typecheck;
111  %return(pcvCV2P(pl, d0, d1));
112}
113
114int Dim(int d0,int d1)
115{
116  %declaration;
117
118  /* check if current RingHandle is set */
119  if(currRingHdl == NULL)
120  {
121    WerrorS("no ring active");
122    return TRUE;
123  }
124 
125  %typecheck;
126  %return(pcvDim);
127}
128
129list Basis(int d0,int d1)
130{
131  %declaration;
132
133  /* check if current RingHandle is set */
134  if(currRingHdl == NULL)
135  {
136    WerrorS("no ring active");
137    return TRUE;
138  }
139 
140  %typecheck;
141  %return(pcvBasis);
142}
143
144%C
145lists pcvLAddL(lists l1,lists l2)
146{
147  lists l0=(lists)omAllocBin(slists_bin);
148  int i=l1->nr;
149  if(l1->nr<l2->nr) i=l2->nr;
150  l0->Init(i+1);
151  for(;i>=0;i--)
152  {
153    if(i<=l1->nr&&(l1->m[i].rtyp==POLY_CMD||l1->m[i].rtyp==VECTOR_CMD))
154    {
155      l0->m[i].rtyp=l1->m[i].rtyp;
156      l0->m[i].data=pCopy((poly)l1->m[i].data);
157      if(i<=l2->nr&&l2->m[i].rtyp==l1->m[i].rtyp)
158        l0->m[i].data=pAdd((poly)l0->m[i].data,pCopy((poly)l2->m[i].data));
159    }
160    else
161    if(i<=l2->nr&&(l2->m[i].rtyp==POLY_CMD||l2->m[i].rtyp==VECTOR_CMD))
162    {
163      l0->m[i].rtyp=l2->m[i].rtyp;
164      l0->m[i].data=pCopy((poly)l2->m[i].data);
165    }
166  }
167  return(l0);
168}
169
170lists pcvPMulL(poly p,lists l1)
171{
172  lists l0=(lists)omAllocBin(slists_bin);
173  l0->Init(l1->nr+1);
174  for(int i=l1->nr;i>=0;i--)
175  {
176    if(l1->m[i].rtyp==POLY_CMD)
177    {
178      l0->m[i].rtyp=POLY_CMD;
179      l0->m[i].data=ppMult_qq(p,(poly)l1->m[i].data);
180    }
181  }
182  return(l0);
183}
184
185BOOLEAN pcvLAddL(leftv res,leftv h)
186{
187  if(h&&h->Typ()==LIST_CMD)
188  {
189    lists l1=(lists)h->Data();
190    h=h->next;
191    if(h&&h->Typ()==LIST_CMD)
192    {
193      lists l2=(lists)h->Data();
194      res->rtyp=LIST_CMD;
195      res->data=(void*)pcvLAddL(l1,l2);
196      return FALSE;
197    }
198  }
199  WerrorS("<list>,<list> expected");
200  return TRUE;
201}
202
203BOOLEAN pcvPMulL(leftv res,leftv h)
204{
205  if(h&&h->Typ()==POLY_CMD)
206  {
207    poly p=(poly)h->Data();
208    h=h->next;
209    if(h&&h->Typ()==LIST_CMD)
210    {
211      lists l=(lists)h->Data();
212      res->rtyp=LIST_CMD;
213      res->data=(void*)pcvPMulL(p,l);
214      return FALSE;
215    }
216  }
217  WerrorS("<poly>,<list> expected");
218  return TRUE;
219}
220
221int pcvDeg(poly p)
222{
223  int d=0;
224  for(int i=pVariables;i>=1;i--) d+=pGetExp(p,i);
225  return d;
226}
227
228int pcvMinDeg(poly p)
229{
230  if(!p) return -1;
231  int md=pcvDeg(p);
232  pIter(p);
233  while(p)
234  {
235    int d=pcvDeg(p);
236    if(d<md) md=d;
237    pIter(p);
238  }
239  return md;
240}
241
242int pcvMinDeg(matrix m)
243{
244  int i,j,d;
245  int md=-1;
246  for(i=1;i<=MATROWS(m);i++)
247  {
248    for(j=1;j<=MATCOLS(m);j++)
249    {
250      d=pcvMinDeg(MATELEM(m,i,j));
251      if((d>=0&&md>d)||md==-1) md=d;
252    }
253  }
254  return(md);
255}
256
257BOOLEAN pcvMinDeg(leftv res,leftv h)
258{
259  if(h)
260  {
261    if(h->Typ()==POLY_CMD)
262    {
263      res->rtyp=INT_CMD;
264      res->data=(void*)pcvMinDeg((poly)h->Data());
265      return FALSE;
266    }
267    else
268    if(h->Typ()==MATRIX_CMD)
269    {
270      res->rtyp=INT_CMD;
271      res->data=(void*)pcvMinDeg((matrix)h->Data());
272      return FALSE;
273    }
274  }
275  WerrorS("<poly> expected");
276  return TRUE;
277}
278
279void pcvInit(int d)
280{
281  if(d<0) d=1;
282  pcvMaxDegree=d+1;
283  pcvTableSize=pVariables*pcvMaxDegree*sizeof(unsigned);
284  pcvTable=(unsigned*)omAlloc0(pcvTableSize);
285  pcvIndexSize=pVariables*sizeof(unsigned*);
286  pcvIndex=(unsigned**)omAlloc(pcvIndexSize);
287  for(int i=0;i<pVariables;i++)
288    pcvIndex[i]=pcvTable+i*pcvMaxDegree;
289  for(int i=0;i<pcvMaxDegree;i++)
290    pcvIndex[0][i]=i;
291  unsigned k,l;
292  for(int i=1;i<pVariables;i++)
293  {
294    k=0;
295    for(int j=0;j<pcvMaxDegree;j++)
296    {
297      l=pcvIndex[i-1][j];
298      if(l>unsigned(~0)-k)
299      {
300        j=pcvMaxDegree;
301        i=pVariables;
302        WerrorS("unsigned overflow");
303      }
304      else pcvIndex[i][j]=k+=l;
305    }
306  }
307}
308
309void pcvClean()
310{
311  if(pcvTable)
312  {
313    omFreeSize(pcvTable,pcvTableSize);
314    pcvTable=NULL;
315  }
316  if(pcvIndex)
317  {
318    omFreeSize(pcvIndex,pcvIndexSize);
319    pcvIndex=NULL;
320  }
321}
322
323int pcvM2N(poly m)
324{
325  unsigned n=0,dn,d=0;
326  for(int i=0;i<pVariables;i++)
327  {
328    d+=pGetExp(m,i+1);
329    dn=pcvIndex[i][d];
330    if(dn>MAX_INT_VAL-n)
331    {
332      i=pVariables;
333      WerrorS("component overflow");
334    }
335    else n+=dn;
336  }
337  return n+1;
338}
339
340poly pcvN2M(int n)
341{
342  n--;
343  poly m=pOne();
344  int i,j,k;
345  for(i=pVariables-1;i>=0;i--)
346  {
347    k=j;
348    for(j=0; (j<pcvMaxDegree) && (pcvIndex[i][j]<=(unsigned)n); j++);
349    j--;
350    n-=pcvIndex[i][j];
351    if(i<pVariables-1) pSetExp(m,i+2,k-j);
352  }
353  if(n==0)
354  {
355    pSetExp(m,1,j);
356    pSetm(m);
357    return m;
358  }
359  else
360  {
361    pDeleteLm(&m);
362    return NULL;
363  }
364}
365
366poly pcvP2CV(poly p,int d0,int d1)
367{
368  poly cv=NULL;
369  while(p)
370  {
371    int d=pcvDeg(p);
372    if(d0<=d&&d<d1)
373    {
374      poly c=pNSet(nCopy(pGetCoeff(p)));
375      pSetComp(c,pcvM2N(p));
376      cv=pAdd(cv,c);
377    }
378    pIter(p);
379  }
380  return cv;
381}
382
383poly pcvCV2P(poly cv,int d0,int d1)
384{
385  poly p=NULL;
386  while(cv)
387  {
388    poly m=pcvN2M(pGetComp(cv));
389    if(m)
390    {
391      int d=pcvDeg(m);
392      if(d0<=d&&d<d1)
393      {
394        pSetCoeff(m,nCopy(pGetCoeff(cv)));
395        p=pAdd(p,m);
396      }
397    }
398    pIter(cv);
399  }
400  return p;
401}
402
403lists pcvP2CV(lists pl,int d0,int d1)
404{
405  lists cvl=(lists)omAllocBin(slists_bin);
406  cvl->Init(pl->nr+1);
407  pcvInit(d1);
408  for(int i=pl->nr;i>=0;i--)
409  {
410    if(pl->m[i].rtyp==POLY_CMD)
411    {
412      cvl->m[i].rtyp=VECTOR_CMD;
413      cvl->m[i].data=pcvP2CV((poly)pl->m[i].data,d0,d1);
414    }
415  }
416  pcvClean();
417  return cvl;
418}
419
420lists pcvCV2P(lists cvl,int d0,int d1)
421{
422  lists pl=(lists)omAllocBin(slists_bin);
423  pl->Init(cvl->nr+1);
424  pcvInit(d1);
425  for(int i=cvl->nr;i>=0;i--)
426  {
427    if(cvl->m[i].rtyp==VECTOR_CMD)
428    {
429      pl->m[i].rtyp=POLY_CMD;
430      pl->m[i].data=pcvCV2P((poly)cvl->m[i].data,d0,d1);
431    }
432  }
433  pcvClean();
434  return pl;
435}
436
437BOOLEAN pcvP2CV(leftv res,leftv h)
438{
439  if(currRingHdl)
440  {
441    if(h&&h->Typ()==LIST_CMD)
442    {
443      lists p=(lists)h->Data();
444      h=h->next;
445      if(h&&h->Typ()==INT_CMD)
446      {
447        int d0=(int)(long)h->Data();
448        h=h->next;
449        if(h&&h->Typ()==INT_CMD)
450        {
451          int d1=(int)(long)h->Data();
452          res->rtyp=LIST_CMD;
453          res->data=(void*)pcvP2CV(p,d0,d1);
454          return FALSE;
455        }
456      }
457    }
458    WerrorS("<list>,<int>,<int> expected");
459    return TRUE;
460  }
461  WerrorS("no ring active");
462  return TRUE;
463}
464
465BOOLEAN pcvCV2P(leftv res,leftv h)
466{
467  if(currRingHdl)
468  {
469    if(h&&h->Typ()==LIST_CMD)
470    {
471      lists pl=(lists)h->Data();
472      h=h->next;
473      if(h&&h->Typ()==INT_CMD)
474      {
475        int d0=(int)(long)h->Data();
476        h=h->next;
477        if(h&&h->Typ()==INT_CMD)
478        {
479          int d1=(int)(long)h->Data();
480          res->rtyp=LIST_CMD;
481          res->data=(void*)pcvCV2P(pl,d0,d1);
482          return FALSE;
483        }
484      }
485    }
486    WerrorS("<list>,<int>,<int> expected");
487    return TRUE;
488  }
489  WerrorS("no ring active");
490  return TRUE;
491}
492
493int pcvDim(int d0,int d1)
494{
495  if(d0<0) d0=0;
496  if(d1<0) d1=0;
497  pcvInit(d1);
498  int d=pcvIndex[pVariables-1][d1]-pcvIndex[pVariables-1][d0];
499  pcvClean();
500  return d;
501}
502
503BOOLEAN pcvDim(leftv res,leftv h)
504{
505  if(currRingHdl)
506  {
507    if(h&&h->Typ()==INT_CMD)
508    {
509      int d0=(int)(long)h->Data();
510      h=h->next;
511      if(h&&h->Typ()==INT_CMD)
512      {
513        int d1=(int)(long)h->Data();
514        res->rtyp=INT_CMD;
515        res->data=(void*)pcvDim(d0,d1);
516        return FALSE;
517      }
518    }
519    WerrorS("<int>,<int> expected");
520    return TRUE;
521  }
522  WerrorS("no ring active");
523  return TRUE;
524}
525
526int pcvBasis(lists b,int i,poly m,int d,int n)
527{
528  if(n<pVariables)
529  {
530    for(int k=0,l=d;k<=l;k++,d--)
531    {
532      pSetExp(m,n,k);
533      i=pcvBasis(b,i,m,d,n+1);
534    }
535  }
536  else
537  {
538    pSetExp(m,n,d);
539    pSetm(m);
540    b->m[i].rtyp=POLY_CMD;
541    b->m[i++].data=pCopy(m);
542  }
543  return i;
544}
545
546lists pcvBasis(int d0,int d1)
547{
548  if(d0<0) d0=0;
549  if(d1<0) d1=0;
550  lists b=(lists)omAllocBin(slists_bin);
551  b->Init(pcvDim(d0,d1));
552  poly m=pOne();
553  for(int d=d0,i=0;d<d1;d++)
554    i=pcvBasis(b,i,m,d,1);
555  pDeleteLm(&m);
556  return b;
557}
558
559BOOLEAN pcvBasis(leftv res,leftv h)
560{
561  if(currRingHdl)
562  {
563    if(h&&h->Typ()==INT_CMD)
564    {
565      int d0=(int)(long)h->Data();
566      h=h->next;
567      if(h&&h->Typ()==INT_CMD)
568      {
569        int d1=(int)(long)h->Data();
570        res->rtyp=LIST_CMD;
571        res->data=(void*)pcvBasis(d0,d1);
572        return FALSE;
573      }
574    }
575    WerrorS("<int>,<int> expected");
576    return TRUE;
577  }
578  WerrorS("no ring active");
579  return TRUE;
580}
Note: See TracBrowser for help on using the repository browser.