# source:git/Singular/LIB/freegb.lib@a8052e

spielwiese
Last change on this file since a8052e was a8052e, checked in by Viktor Levandovskyy <levandov@…>, 16 years ago
*levandov: some new freegb functions git-svn-id: file:///usr/local/Singular/svn/trunk@10578 2c84dea3-7e68-4137-9b89-c4e89433aadc
• Property mode set to `100644`
File size: 23.0 KB
Line
1//////////////////////////////////////////////////////////////////////////////
2version="\$Id: freegb.lib,v 1.5 2008-02-15 17:12:41 levandov Exp \$";
3category="Noncommutative";
4info="
5LIBRARY: ratgb.lib  Twosided Noncommutative Groebner bases in Free Algebras
6AUTHOR: Viktor Levandovskyy,     levandov@math.rwth-aachen.de
7
8PROCEDURES:
9freegb(list L, int n);   compute two-sided Groebner basis of ideal, encoded via L, up to degree n
10lst2str(list L);         convert a list (of modules) into polynomials in free algebra
11mod2str(module M);       convert a module into a polynomial in free algebra
12"
13
14// this library computes two-sided GB of an ideal
15// in a free associative algebra
16
17// a monomial is encoded via a vector V
18// where V[1] = coefficient
19// V[1+i] = the corresponding symbol
20
21LIB "qhmoduli.lib"; // for Max
22
23proc vct2mono(vector v)
24{
25  // produces a monomial in new vars
26  // need a ring with new vars!!
27}
28
29proc lshift(module M, int s, string varing, def lpring)
30{
31  // FINALLY IMPLEMENTED AS A PART OT THE CODE
32  // shifts a poly from the ring @R to s positions
33  // M lives in varing, the result in lpring
34  // to be run from varing
35  int i, j, k, sm, sv;
36  vector v;
37  //  execute("setring "+lpring);
38  setring lpring;
39  poly @@p;
40  ideal I;
41  execute("setring "+varing);
42  sm = ncols(M);
43  for (i=1; i<=s; i++)
44  {
45    // modules, e.g. free polynomials
46    for (j=1; j<=sm; j++)
47    {
48      //vectors, e.g. free monomials
49      v  = M[j];
50      sv = size(v);
51      sp = "@@p = @@p + ";
52      for (k=2; k<=sv; k++)
53      {
54        sp = sp + string(v[k])+"("+string(k-1+s)+")*";
55      }
56      sp = sp + string(v[1])+";"; // coef;
57      setring lpring;
58      //      execute("setring "+lpring);
59      execute(sp);
60      execute("setring "+varing);
61    }
62    setring lpring;
63    //    execute("setring "+lpring);
64    I = I,@@p;
65    @@p = 0;
66  }
67  setring lpring;
68  //execute("setring "+lpring);
69  export(I);
70  //  setring varing;
71  execute("setring "+varing);
72}
73
74proc skip0(vector v)
75{
76  // skips zeros in vector
77  int sv = nrows(v);
78  int sw = size(v);
79  if (sv == sw)
80  {
81    return(v);
82  }
83  int i;
84  int j=1;
85  vector w;
86  for (i=1; i<=sv; i++)
87  {
88    if (v[i] != 0)
89    {
90      w = w + v[i]*gen(j);
91      j++;
92    }
93  }
94  return(w);
95}
96
97proc lst2str(list L)
98"USAGE:  lst2str(L);  L a list of modules
99RETURN:  list (of strings)
100PURPOSE: convert a list (of modules) into polynomials in free algebra
101EXAMPLE: example lst2str; shows examples
102"
103{
104  // returns a list of strings
105  // being sentences in words built from L
106  int i;
107  int s    = size(L);
108  list N;
109  for(i=1; i<=s; i++)
110  {
111    if ((typeof(L[i]) == "module") || (typeof(L[i]) == "matrix") )
112    {
113      N[i] = mod2str(L[i]);
114    }
115    else
116    {
117      "module or matrix expected in the list";
118      return(N);
119    }
120  }
121  return(N);
122}
123example
124{
125  "EXAMPLE:"; echo = 2;
126  ring r = 0,(x,y,z),(dp(1),dp(2));
127  module M = [-1,x,y],[-7,y,y],[3,x,x];
128  module N = [1,x,y,x,y],[-2,y,x,y,x],[6,x,y,y,x,y];
129  list L; L[1] = M; L[2] = N;
130  lst2str(L);
131}
132
133
134proc mod2str(module M)
135"USAGE:  mod2str(M);  M a module
136RETURN:  string
137PURPOSE: convert a modules into a polynomial in free algebra
138EXAMPLE: example mod2str; shows examples
139"
140{
141  // returns a string
142  // a sentence in words built from M
143  int i;
144  int s    = ncols(M);
145  string t;
146  string mp;
147  for(i=1; i<=s; i++)
148  {
149    mp = vct2str(M[i]);
150    if (mp[1] == "-")
151    {
152      t = t + mp;
153    }
154    else
155    {
156      t = t + "+" + mp;
157    }
158  }
159  if (t[1]=="+")
160  {
161    t = t[2..size(t)]; // remove first "+"
162  }
163  return(t);
164}
165example
166{
167  "EXAMPLE:"; echo = 2;
168  ring r = 0,(x,y,z),(dp);
169  module M = [1,x,y,x,y],[-2,y,x,y,x],[6,x,y,y,x,y];
170  mod2str(M);
171}
172
173proc vct2str(vector v)
174{
175  int ppl = printlevel-voice+2;
176  // for a word, encoded by v
177  // produces a string for it
178  v = skip0(v);
180  int s = size(v);
181  string vs,vv,vp,err;
182  int i,j,p,q;
183  for (i=1; i<=s-1; i++)
184  {
185    p     = IsVar(v[i+1]);
186    if (p==0)
187    {
188      err = "Error: monomial expected at" + string(i+1);
189      dbprint(ppl,err);
190      return("_");
191    }
192    if (p==1)
193    {
194      vs = vs + string(v[i+1]);
195    }
196    else //power
197    {
198      vv = string(v[i+1]);
199      q = find(vv,"^");
200      if (q==0)
201      {
202        q = find(vv,string(p));
203        if (q==0)
204        {
205          err = "error in find for string "+vv;
206          dbprint(ppl,err);
207          return("_");
208        }
209      }
210      // q>0
211      vp = vv[1..q-1];
212      for(j=1;j<=p;j++)
213      {
214        vs = vs + vp;
215      }
216    }
217  }
218  string scf;
219  if (cf == -1)
220  {
221    scf = "-";
222  }
223  else
224  {
225    scf = string(cf);
226    if (cf == 1)
227    {
228      scf = "";
229    }
230  }
231  vs = scf + vs;
232  return(vs);
233}
234example
235{
236  ring r = (0,a),(x,y3,z(1)),dp;
237  vector v = [-7,x,y3^4,x2,z(1)^3];
238  vct2str(v);
239  vector w = [-7a^5+6a,x,y3,y3,x,z(1),z(1)];
240  vct2str(w);
241}
242
243proc IsVar(poly p)
244{
245  // checks whether p is a variable indeed
246  // if it's a power of a variable, returns the power
247  if (p==0) {  return(0); } //"p=0";
249  if ( (p-lead(p)) !=0 ) { return(0); } // "p-lm(p)>0";
251  int s = size(v);
252  int i=1;
253  int cnt = 0;
254  int pwr = 0;
255  for (i=1; i<=s; i++)
256  {
257    if (v[i] != 0)
258    {
259      cnt++;
260      pwr = v[i];
261    }
262  }
263  //  "cnt:";  cnt;
264  if (cnt==1) { return(pwr); }
265  else { return(0); }
266}
267example
268{
269  ring r = 0,(x,y),dp;
270  poly f = xy+1;
271  IsVar(f);
272  poly g = xy;
273  IsVar(g);
274  poly h = y^3;
275  IsVar(h);
276  poly i = 1;
277  IsVar(i);
278}
279
280// given the element -7xy^2x, it is represented as [-7,x,y^2,x] or as [-7,x,y,y,x]
281// use the orig ord on (x,y,z) and expand it blockwise to (x(i),y(i),z(i))
282
283// the correspondences:
284// monomial in K<x,y,z>    <<--->> vector in R
285// polynomial in K<x,y,z>  <<--->> list of vectors (matrix/module) in R
286// ideal in K<x,y,z>       <<--->> list of matrices/modules in R
287
288
289// 1. form a new ring
290// 2. produce shifted generators
291// 3. compute GB
292// 4. skip shifted elts
293// 5. go back to orig vars, produce strings/modules
294// 6. return the result
295
296proc freegb(list LM, int d)
297"USAGE:  freegb(L, d);  L a list of modules, d an integer
298RETURN:  ring
299PURPOSE: compute the two-sided Groebner basis of an ideal, encoded by L in
300the free associative algebra, up to degree d
301EXAMPLE: example freegb; shows examples
302"
303{
304  // d = up to degree, will be shifted to d+1
305  if (d<1) {"bad d"; return(0);}
306
307  int ppl = printlevel-voice+2;
308  string err = "";
309
310  int i,j,s;
311  def save = basering;
312  // determine max no of places in the input
313  int slm = size(LM); // numbers of polys in the ideal
314  int sm;
315  intvec iv;
316  module M;
317  for (i=1; i<=slm; i++)
318  {
319    // modules, e.g. free polynomials
320    M  = LM[i];
321    sm = ncols(M);
322    for (j=1; j<=sm; j++)
323    {
324      //vectors, e.g. free monomials
325      iv = iv, size(M[j])-1; // 1 place is reserved by the coeff
326    }
327  }
328  int D = Max(iv); // max size of input words
329  if (d<D) {"bad d"; return(LM);}
330  D = D + d-1;
331  //  D = d;
332  list LR  = ringlist(save);
333  list L, tmp;
334  L[1] = LR[1]; // ground field
335  L[4] = LR[4]; // quotient ideal
336  tmp  = LR[2]; // varnames
337  s = size(LR[2]);
338  for (i=1; i<=D; i++)
339  {
340    for (j=1; j<=s; j++)
341    {
342      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
343    }
344  }
345  for (i=1; i<=s; i++)
346  {
347    tmp[i] = string(tmp[i])+"("+string(1)+")";
348  }
349  L[2] = tmp;
350  list OrigNames = LR[2];
351  // ordering: d blocks of the ord on r
352  // try to get whether the ord on r is blockord itself
353  // TODO: make L(2) ordering! exponent is maximally 2
354  s = size(LR[3]);
355  if (s==2)
356  {
357    // not a blockord, 1 block + module ord
358    tmp = LR[3][s]; // module ord
359    for (i=1; i<=D; i++)
360    {
361      LR[3][s-1+i] = LR[3][1];
362    }
363    LR[3][s+D] = tmp;
364  }
365  if (s>2)
366  {
367    // there are s-1 blocks
368    int nb = s-1;
369    tmp = LR[3][s]; // module ord
370    for (i=1; i<=D; i++)
371    {
372      for (j=1; j<=nb; j++)
373      {
374        LR[3][i*nb+j] = LR[3][j];
375      }
376    }
377    //    size(LR[3]);
378    LR[3][nb*(D+1)+1] = tmp;
379  }
380  L[3] = LR[3];
381  def @R = ring(L);
382  setring @R;
383  ideal I;
384  poly @p;
385  s = size(OrigNames);
386  //  "s:";s;
387  // convert LM to canonical vectors (no powers)
388  setring save;
389  kill M; // M was defined earlier
390  module M;
391  slm = size(LM); // numbers of polys in the ideal
392  int sv,k,l;
393  vector v;
394  //  poly p;
395  string sp;
396  setring @R;
397  poly @@p=0;
398  setring save;
399  for (l=1; l<=slm; l++)
400  {
401    // modules, e.g. free polynomials
402    M  = LM[l];
403    sm = ncols(M); // in intvec iv the sizes are stored
404    for (i=0; i<=d-iv[l]; i++)
405    {
406      // modules, e.g. free polynomials
407      for (j=1; j<=sm; j++)
408      {
409        //vectors, e.g. free monomials
410        v  = M[j];
411        sv = size(v);
412        //      "sv:";sv;
413        sp = "@@p = @@p + ";
414        for (k=2; k<=sv; k++)
415        {
416          sp = sp + string(v[k])+"("+string(k-1+i)+")*";
417        }
418        sp = sp + string(v[1])+";"; // coef;
419        setring @R;
420        execute(sp);
421        setring save;
422      }
423      setring @R;
424      //      "@@p:"; @@p;
425      I = I,@@p;
426      @@p = 0;
427      setring save;
428    }
429  }
430  kill sp;
431  // 3. compute GB
432  setring @R;
433  dbprint(ppl,"computing GB");
434  //  ideal J = groebner(I);
435  ideal J = slimgb(I);
436  dbprint(ppl,J);
437  // 4. skip shifted elts
438  ideal K = select1(J,1,s); // s = size(OrigNames)
439  dbprint(ppl,K);
440  dbprint(ppl, "done with GB");
441  // K contains vars x(1),...z(1) = images of originals
442  // 5. go back to orig vars, produce strings/modules
443  if (K[1] == 0)
444  {
445    "no reasonable output, GB gives 0";
446    return(0);
447  }
448  int sk = size(K);
449  int sp, sx, a, b;
450  intvec x;
451  poly p,q;
452  poly pn;
453  // vars in 'save'
454  setring save;
455  module N;
456  list LN;
457  vector V;
458  poly pn;
459  // test and skip exponents >=2
460  setring @R;
461  for(i=1; i<=sk; i++)
462  {
463    p  = K[i];
464    while (p!=0)
465    {
467      //      "processing q:";q;
469      sx = size(x);
470      for(k=1; k<=sx; k++)
471      {
472        if ( x[k] >= 2 )
473        {
474          err = "skip: the value x[k] is " + string(x[k]);
475          dbprint(ppl,err);
476          //        return(0);
477          K[i] = 0;
478          p    = 0;
479          q    = 0;
480          break;
481        }
482      }
483      p  = p - q;
484    }
485  }
486  K  = simplify(K,2);
487  sk = size(K);
488  for(i=1; i<=sk; i++)
489  {
490    //    setring save;
491    //    V  = 0;
492    setring @R;
493    p  = K[i];
494    while (p!=0)
495    {
497      err =  "processing q:" + string(q);
498      dbprint(ppl,err);
500      sx = size(x);
502      setring save;
503      pn = imap(@R,pn);
504      V  = V + leadcoef(pn)*gen(1);
505      for(k=1; k<=sx; k++)
506      {
507        if (x[k] ==1)
508        {
509          a = k / s; // block number=a+1, a!=0
510          b = k % s; // remainder
511          //      printf("a: %s, b: %s",a,b);
512          if (b == 0)
513          {
514            // that is it's the last var in the block
515            b = s;
516            a = a-1;
517          }
518          V = V + var(b)*gen(a+2);
519        }
520//      else
521//      {
522//        printf("error: the value x[k] is %s", x[k]);
523//        return(0);
524//      }
525      }
526      err = "V: " + string(V);
527      dbprint(ppl,err);
528      //      printf("V: %s", string(V));
529      N = N,V;
530      V  = 0;
531      setring @R;
532      p  = p - q;
533      pn = 0;
534    }
535    setring save;
536    LN[i] = simplify(N,2);
537    N     = 0;
538  }
539  setring save;
540  return(LN);
541}
542example
543{
544  "EXAMPLE:"; echo = 2;
545  ring r = 0,(x,y,z),(dp(1),dp(2));
546  module M = [-1,x,y],[-7,y,y],[3,x,x];
547  module N = [1,x,y,x],[-1,y,x,y];
548  list L; L[1] = M; L[2] = N;
549  lst2str(L);
550  def U = freegb(L,5);
551  lst2str(U);
552}
553
554// 1. form a new ring
555// 2. NOP
556// 3. compute GB -> with the kernel stuff
557// 4. skip shifted elts (check that no such exist?)
558// 5. go back to orig vars, produce strings/modules
559// 6. return the result
560
561proc freegbnew(list LM, int d)
562"USAGE:  freegb(L, d);  L a list of modules, d an integer
563RETURN:  ring
564PURPOSE: compute the two-sided Groebner basis of an ideal, encoded by L in
565the free associative algebra, up to degree d
566EXAMPLE: example freegb; shows examples
567"
568{
569  // d = up to degree, will be shifted to d+1
570  if (d<1) {"bad d"; return(0);}
571
572  int ppl = printlevel-voice+2;
573  string err = "";
574
575  int i,j,s;
576  def save = basering;
577  // determine max no of places in the input
578  int slm = size(LM); // numbers of polys in the ideal
579  int sm;
580  intvec iv;
581  module M;
582  for (i=1; i<=slm; i++)
583  {
584    // modules, e.g. free polynomials
585    M  = LM[i];
586    sm = ncols(M);
587    for (j=1; j<=sm; j++)
588    {
589      //vectors, e.g. free monomials
590      iv = iv, size(M[j])-1; // 1 place is reserved by the coeff
591    }
592  }
593  int D = Max(iv); // max size of input words
594  if (d<D) {"bad d"; return(LM);}
595  D = D + d-1;
596  //  D = d;
597  list LR  = ringlist(save);
598  list L, tmp;
599  L[1] = LR[1]; // ground field
600  L[4] = LR[4]; // quotient ideal
601  tmp  = LR[2]; // varnames
602  s = size(LR[2]);
603  for (i=1; i<=D; i++)
604  {
605    for (j=1; j<=s; j++)
606    {
607      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
608    }
609  }
610  for (i=1; i<=s; i++)
611  {
612    tmp[i] = string(tmp[i])+"("+string(1)+")";
613  }
614  L[2] = tmp;
615  list OrigNames = LR[2];
616  // ordering: d blocks of the ord on r
617  // try to get whether the ord on r is blockord itself
618  s = size(LR[3]);
619  if (s==2)
620  {
621    // not a blockord, 1 block + module ord
622    tmp = LR[3][s]; // module ord
623    for (i=1; i<=D; i++)
624    {
625      LR[3][s-1+i] = LR[3][1];
626    }
627    LR[3][s+D] = tmp;
628  }
629  if (s>2)
630  {
631    // there are s-1 blocks
632    int nb = s-1;
633    tmp = LR[3][s]; // module ord
634    for (i=1; i<=D; i++)
635    {
636      for (j=1; j<=nb; j++)
637      {
638        LR[3][i*nb+j] = LR[3][j];
639      }
640    }
641    //    size(LR[3]);
642    LR[3][nb*(D+1)+1] = tmp;
643  }
644  L[3] = LR[3];
645  def @R = ring(L);
646  setring @R;
647  ideal I;
648  poly @p;
649  s = size(OrigNames);
650  //  "s:";s;
651  // convert LM to canonical vectors (no powers)
652  setring save;
653  kill M; // M was defined earlier
654  module M;
655  slm = size(LM); // numbers of polys in the ideal
656  int sv,k,l;
657  vector v;
658  //  poly p;
659  string sp;
660  setring @R;
661  poly @@p=0;
662  setring save;
663  for (l=1; l<=slm; l++)
664  {
665    // modules, e.g. free polynomials
666    M  = LM[l];
667    sm = ncols(M); // in intvec iv the sizes are stored
668    // modules, e.g. free polynomials
669    for (j=1; j<=sm; j++)
670    {
671      //vectors, e.g. free monomials
672      v  = M[j];
673      sv = size(v);
674      //        "sv:";sv;
675      sp = "@@p = @@p + ";
676      for (k=2; k<=sv; k++)
677      {
678        sp = sp + string(v[k])+"("+string(k-1)+")*";
679      }
680      sp = sp + string(v[1])+";"; // coef;
681      setring @R;
682      execute(sp);
683      setring save;
684    }
685    setring @R;
686    //      "@@p:"; @@p;
687    I = I,@@p;
688    @@p = 0;
689    setring save;
690  }
691  kill sp;
692  // 3. compute GB
693  setring @R;
694  dbprint(ppl,"computing GB");
695  //  ideal J = system("",I);
696  ideal J = slimgb(I);
697  dbprint(ppl,J);
698  // 4. skip shifted elts
699  ideal K = select1(J,1,s); // s = size(OrigNames)
700  dbprint(ppl,K);
701  dbprint(ppl, "done with GB");
702  // K contains vars x(1),...z(1) = images of originals
703  // 5. go back to orig vars, produce strings/modules
704  if (K[1] == 0)
705  {
706    "no reasonable output, GB gives 0";
707    return(0);
708  }
709  int sk = size(K);
710  int sp, sx, a, b;
711  intvec x;
712  poly p,q;
713  poly pn;
714  // vars in 'save'
715  setring save;
716  module N;
717  list LN;
718  vector V;
719  poly pn;
720  // test and skip exponents >=2
721  setring @R;
722  for(i=1; i<=sk; i++)
723  {
724    p  = K[i];
725    while (p!=0)
726    {
728      //      "processing q:";q;
730      sx = size(x);
731      for(k=1; k<=sx; k++)
732      {
733        if ( x[k] >= 2 )
734        {
735          err = "skip: the value x[k] is " + string(x[k]);
736          dbprint(ppl,err);
737          //        return(0);
738          K[i] = 0;
739          p    = 0;
740          q    = 0;
741          break;
742        }
743      }
744      p  = p - q;
745    }
746  }
747  K  = simplify(K,2);
748  sk = size(K);
749  for(i=1; i<=sk; i++)
750  {
751    //    setring save;
752    //    V  = 0;
753    setring @R;
754    p  = K[i];
755    while (p!=0)
756    {
758      err =  "processing q:" + string(q);
759      dbprint(ppl,err);
761      sx = size(x);
763      setring save;
764      pn = imap(@R,pn);
765      V  = V + leadcoef(pn)*gen(1);
766      for(k=1; k<=sx; k++)
767      {
768        if (x[k] ==1)
769        {
770          a = k / s; // block number=a+1, a!=0
771          b = k % s; // remainder
772          //      printf("a: %s, b: %s",a,b);
773          if (b == 0)
774          {
775            // that is it's the last var in the block
776            b = s;
777            a = a-1;
778          }
779          V = V + var(b)*gen(a+2);
780        }
781//      else
782//      {
783//        printf("error: the value x[k] is %s", x[k]);
784//        return(0);
785//      }
786      }
787      err = "V: " + string(V);
788      dbprint(ppl,err);
789      //      printf("V: %s", string(V));
790      N = N,V;
791      V  = 0;
792      setring @R;
793      p  = p - q;
794      pn = 0;
795    }
796    setring save;
797    LN[i] = simplify(N,2);
798    N     = 0;
799  }
800  setring save;
801  return(LN);
802}
803example
804{
805  "EXAMPLE:"; echo = 2;
806  ring r = 0,(x,y,z),(dp(1),dp(2));
807  module M = [-1,x,y],[-7,y,y],[3,x,x];
808  module N = [1,x,y,x],[-1,y,x,y];
809  list L; L[1] = M; L[2] = N;
810  lst2str(L);
811  def U = freegbnew(L,5);
812  lst2str(U);
813}
814
815proc crs(list LM, int d)
816"USAGE:  crs(L, d);  L a list of modules, d an integer
817RETURN:  ring
818PURPOSE: create a ring and shift the ideal
819EXAMPLE: example crs; shows examples
820"
821{
822  // d = up to degree, will be shifted to d+1
823  if (d<1) {"bad d"; return(0);}
824
825  int ppl = printlevel-voice+2;
826  string err = "";
827
828  int i,j,s;
829  def save = basering;
830  // determine max no of places in the input
831  int slm = size(LM); // numbers of polys in the ideal
832  int sm;
833  intvec iv;
834  module M;
835  for (i=1; i<=slm; i++)
836  {
837    // modules, e.g. free polynomials
838    M  = LM[i];
839    sm = ncols(M);
840    for (j=1; j<=sm; j++)
841    {
842      //vectors, e.g. free monomials
843      iv = iv, size(M[j])-1; // 1 place is reserved by the coeff
844    }
845  }
846  int D = Max(iv); // max size of input words
847  if (d<D) {"bad d"; return(LM);}
848  D = D + d-1;
849  //  D = d;
850  list LR  = ringlist(save);
851  list L, tmp;
852  L[1] = LR[1]; // ground field
853  L[4] = LR[4]; // quotient ideal
854  tmp  = LR[2]; // varnames
855  s = size(LR[2]);
856  for (i=1; i<=D; i++)
857  {
858    for (j=1; j<=s; j++)
859    {
860      tmp[i*s+j] = string(tmp[j])+"("+string(i)+")";
861    }
862  }
863  for (i=1; i<=s; i++)
864  {
865    tmp[i] = string(tmp[i])+"("+string(0)+")";
866  }
867  L[2] = tmp;
868  list OrigNames = LR[2];
869  // ordering: d blocks of the ord on r
870  // try to get whether the ord on r is blockord itself
871  s = size(LR[3]);
872  if (s==2)
873  {
874    // not a blockord, 1 block + module ord
875    tmp = LR[3][s]; // module ord
876    for (i=1; i<=D; i++)
877    {
878      LR[3][s-1+i] = LR[3][1];
879    }
880    LR[3][s+D] = tmp;
881  }
882  if (s>2)
883  {
884    // there are s-1 blocks
885    int nb = s-1;
886    tmp = LR[3][s]; // module ord
887    for (i=1; i<=D; i++)
888    {
889      for (j=1; j<=nb; j++)
890      {
891        LR[3][i*nb+j] = LR[3][j];
892      }
893    }
894    //    size(LR[3]);
895    LR[3][nb*(D+1)+1] = tmp;
896  }
897  L[3] = LR[3];
898  def @R = ring(L);
899  setring @R;
900  ideal I;
901  poly @p;
902  s = size(OrigNames);
903  //  "s:";s;
904  // convert LM to canonical vectors (no powers)
905  setring save;
906  kill M; // M was defined earlier
907  module M;
908  slm = size(LM); // numbers of polys in the ideal
909  int sv,k,l;
910  vector v;
911  //  poly p;
912  string sp;
913  setring @R;
914  poly @@p=0;
915  setring save;
916  for (l=1; l<=slm; l++)
917  {
918    // modules, e.g. free polynomials
919    M  = LM[l];
920    sm = ncols(M); // in intvec iv the sizes are stored
921    for (i=0; i<=d-iv[l]; i++)
922    {
923      // modules, e.g. free polynomials
924      for (j=1; j<=sm; j++)
925      {
926        //vectors, e.g. free monomials
927        v  = M[j];
928        sv = size(v);
929        //      "sv:";sv;
930        sp = "@@p = @@p + ";
931        for (k=2; k<=sv; k++)
932        {
933          sp = sp + string(v[k])+"("+string(k-2+i)+")*";
934        }
935        sp = sp + string(v[1])+";"; // coef;
936        setring @R;
937        execute(sp);
938        setring save;
939      }
940      setring @R;
941      //      "@@p:"; @@p;
942      I = I,@@p;
943      @@p = 0;
944      setring save;
945    }
946  }
947  setring @R;
948  export I;
949  return(@R);
950}
951example
952{
953  "EXAMPLE:"; echo = 2;
954  ring r = 0,(x,y,z),(dp(1),dp(2));
955  module M = [-1,x,y],[-7,y,y],[3,x,x];
956  module N = [1,x,y,x],[-1,y,x,y];
957  list L; L[1] = M; L[2] = N;
958  lst2str(L);
959  def U = crs(L,5);
960  setring U; U;
961  I;
962}
963
964proc ex_shift()
965{
966  LIB "freegb.lib";
967  ring r = 0,(x,y,z),(dp(1),dp(2));
968  module M = [-1,x,y],[-7,y,y],[3,x,x];
969  module N = [1,x,y,x],[-1,y,x,y];
970  list L; L[1] = M; L[2] = N;
971  lst2str(L);
972  def U = crs(L,5);
973  setring U; U;
974  I;
975  poly p = I[2]; // I[8];
976  p;
977  system("stest",p,7,7,3); // error -> the world is ok
978  poly q1 = system("stest",p,1,7,3); //ok
979  poly q6 = system("stest",p,6,7,3); //ok
980  system("btest",p,3); //ok
981  system("btest",q1,3); //ok
982  system("btest",q6,3); //ok
983}
984
985proc ex2()
986{
987  option(prot);
988  LIB "freegb.lib";
989  ring r = 0,(x,y),dp;
990  module M = [-1,x,y],[3,x,x]; // 3x^2 - xy
991  def U = freegb(M,7);
992  lst2str(U);
993}
994
995proc ex_nonhomog()
996{
997  option(prot);
998  LIB "freegb.lib";
999  ring r = 0,(x,y,h),dp;
1000  list L;
1001  module M;
1002  M = [-1,y,y],[1,x,x,x];  // x3-y2
1003  L[1] = M;
1004  M = [1,x,h],[-1,h,x];  // xh-hx
1005  L[2] = M;
1006  M = [1,y,h],[-1,h,y];  // yh-hy
1007  L[3] = M;
1008  def U = freegb(L,4);
1009  lst2str(U);
1010  // strange elements in the basis
1011}
1012
1013proc ex_nonhomog_comm()
1014{
1015  option(prot);
1016  LIB "freegb.lib";
1017  ring r = 0,(x,y),dp;
1018  module M = [-1,y,y],[1,x,x,x];
1019  def U = freegb(M,5);
1020  lst2str(U);
1021}
1022
1023proc ex_nonhomog_h()
1024{
1025  option(prot);
1026  LIB "freegb.lib";
1027  ring r = 0,(x,y,h),(a(1,1),dp);
1028  module M = [-1,y,y,h],[1,x,x,x]; // x3 - y2h
1029  def U = freegb(M,6);
1030  lst2str(U);
1031}
1032
1033proc ex_nonhomog_h2()
1034{
1035  option(prot);
1036  LIB "freegb.lib";
1037  ring r = 0,(x,y,h),(dp);
1038  list L;
1039  module M;
1040  M = [-1,y,y,h],[1,x,x,x]; // x3 - y2h
1041  L[1] = M;
1042  M = [1,x,h],[-1,h,x]; // xh - hx
1043  L[2] = M;
1044  M = [1,y,h],[-1,h,y]; // yh - hy
1045  L[3] = M;
1046  def U = freegb(L,3);
1047  lst2str(U);
1049}
1050
1051
1052proc ex_nonhomog_3()
1053{
1054  option(prot);
1055  LIB "./freegb.lib";
1056  ring r = 0,(x,y,z),(dp);
1057  list L;
1058  module M;
1059  M = [1,z,y],[-1,x]; // zy - x
1060  L[1] = M;
1061  M = [1,z,x],[-1,y]; // zx - y
1062  L[2] = M;
1063  M = [1,y,x],[-1,z]; // yx - z
1064  L[3] = M;
1065  lst2str(L);
1066  list U = freegb(L,4);
1067  lst2str(U);
1069}
1070
1071
1072
1073proc ex_densep_2()
1074{
1075  option(prot);
1076  LIB "freegb.lib";
1077  ring r = (0,a,b,c),(x,y),(Dp); // deglex
1078  module M = [1,x,x], [a,x,y], [b,y,x], [c,y,y];
1079  lst2str(M);
1080  list U = freegb(M,5);
1081  lst2str(U);
1082  // a=b is important -> finite basis!!!
1083  module M = [1,x,x], [a,x,y], [a,y,x], [c,y,y];
1084  lst2str(M);
1085  list U = freegb(M,5);
1086  lst2str(U);
1087}
1088
1089
1090proc freegbRing(int d)
1091"USAGE:  freegbRing(d); d an integer
1092RETURN:  ring
1093PURPOSE: creates a d-shifted ring
1094EXAMPLE: example freegbRing; shows examples
1095"
1096{
1097  // d = up to degree, will be shifted to d+1
1098  if (d<1) {"bad d"; return(0);}
1099
1100  int ppl = printlevel-voice+2;
1101  string err = "";
1102
1103  int i,j,s;
1104  def save = basering;
1105  int D = d-1;
1106  list LR  = ringlist(save);
1107  list L, tmp;
1108  L[1] = LR[1]; // ground field
1109  L[4] = LR[4]; // quotient ideal
1110  tmp  = LR[2]; // varnames
1111  s = size(LR[2]);
1112  for (i=1; i<=D; i++)
1113  {
1114    for (j=1; j<=s; j++)
1115    {
1116      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
1117    }
1118  }
1119  for (i=1; i<=s; i++)
1120  {
1121    tmp[i] = string(tmp[i])+"("+string(1)+")";
1122  }
1123  L[2] = tmp;
1124  list OrigNames = LR[2];
1125  // ordering: d blocks of the ord on r
1126  // try to get whether the ord on r is blockord itself
1127  // TODO: make L(2) ordering! exponent is maximally 2
1128  s = size(LR[3]);
1129  if (s==2)
1130  {
1131    // not a blockord, 1 block + module ord
1132    tmp = LR[3][s]; // module ord
1133    for (i=1; i<=D; i++)
1134    {
1135      LR[3][s-1+i] = LR[3][1];
1136    }
1137    LR[3][s+D] = tmp;
1138  }
1139  if (s>2)
1140  {
1141    // there are s-1 blocks
1142    int nb = s-1;
1143    tmp = LR[3][s]; // module ord
1144    for (i=1; i<=D; i++)
1145    {
1146      for (j=1; j<=nb; j++)
1147      {
1148        LR[3][i*nb+j] = LR[3][j];
1149      }
1150    }
1151    //    size(LR[3]);
1152    LR[3][nb*(D+1)+1] = tmp;
1153  }
1154  L[3] = LR[3];
1155  def @R = ring(L);
1156  //  setring @R;
1157  return (@R);
1158}
1159example
1160{
1161  "EXAMPLE:"; echo = 2;
1162  ring r = 0,(x,y,z),(dp(1),dp(2));
1163  def A = freegbRing(2);
1164  setring A;
1165  A;
1166}
Note: See TracBrowser for help on using the repository browser.