source: git/Singular/LIB/freegb.lib @ 0be039

spielwiese
Last change on this file since 0be039 was 0be039, checked in by Karim Abou Zeid <karim23697@…>, 6 years ago
Merge branch 'lpDivision' into stable
  • Property mode set to 100644
File size: 99.2 KB
Line 
1/////////////////////////////////////////////////////////////////////////////
2version="version freegb.lib 4.1.1.0 Dec_2017 "; // $Id$
3category="Noncommutative";
4info="
5LIBRARY: freegb.lib   Compute two-sided Groebner bases in free algebras via letterplace approach
6AUTHORS: Viktor Levandovskyy,     viktor.levandovskyy@math.rwth-aachen.de
7       Grischa Studzinski,      grischa.studzinski@math.rwth-aachen.de
8
9OVERVIEW: For the theory, see chapter 'Letterplace' in the @sc{Singular} Manual
10
11Support: Joint projects LE 2697/2-1 and KR 1907/3-1 of the Priority Programme SPP 1489:
12'Algorithmische und Experimentelle Methoden in Algebra, Geometrie und Zahlentheorie'
13of the German DFG
14and Project II.6 of the transregional collaborative research centre
15SFB-TRR 195 'Symbolic Tools in Mathematics and their Application' of the German DFG
16
17PROCEDURES:
18makeLetterplaceRing(d); creates a ring with d blocks of shifted original variables
19letplaceGBasis(I); computes two-sided Groebner basis of a letterplace ideal I up to a degree bound
20lpNF(f,I); two-sided normal form of f with respect to ideal I
21setLetterplaceAttributes(R,d,b); supplies ring R with the letterplace structure
22freeGBasis(L, n);  computes two-sided Groebner basis of an ideal, encoded via list L, up to degree n
23
24lpMult(f,g);    letterplace multiplication of letterplace polynomials
25shiftPoly(p,i); compute the i-th shift of letterplace polynomial p
26lpPower(f,n);   natural power of a letterplace polynomial
27lieBracket(a,b[, N]);  compute Lie bracket ab-ba of two letterplace polynomials
28
29lp2lstr(K, s);  convert a letterplace ideal into a list of modules
30lst2str(L[, n]); convert a list (of modules) into polynomials in free algebra via strings
31mod2str(M[, n]); convert a module into a polynomial in free algebra via strings
32vct2str(M[, n]);   convert a vector into a word in free algebra
33
34serreRelations(A,z); compute the homogeneous part of Serre's relations associated to a generalized Cartan matrix A
35fullSerreRelations(A,N,C,P,d); compute the ideal of all Serre's relations associated to a generalized Cartan matrix A
36isVar(p);              check whether p is a power of a single variable
37ademRelations(i,j);    compute the ideal of Adem relations for i<2j in char 0
38
39lpPrint(ideal I, def @r); print a letterplace ideal as an easily readable string
40
41SEE ALSO: fpadim_lib, LETTERPLACE
42";
43
44// this library computes two-sided GB of an ideal
45// in a free associative algebra
46
47// a monomial is encoded via a vector V
48// where V[1] = coefficient
49// V[1+i] = the corresponding symbol
50
51LIB "qhmoduli.lib"; // for Max
52LIB "bfun.lib"; // for inForm
53LIB "fpadim.lib"; // for intvec conversion
54
55proc tstfreegb()
56{
57    /* tests all procs for consistency */
58  /* adding the new proc, add it here */
59  example makeLetterplaceRing;
60  example letplaceGBasis;
61  example lpNF;
62  example freeGBasis;
63  example setLetterplaceAttributes;
64  /* secondary */
65  example   lpMult;
66  example   shiftPoly;
67  example   lpPower;
68  example   lieBracket;
69  example   lp2lstr;
70  example   lst2str;
71  example   mod2str;
72  example   vct2str;
73  example   serreRelations;
74  example   fullSerreRelations;
75  example   isVar;
76  example   ademRelations;
77}
78
79proc setLetterplaceAttributes(def R, int uptodeg, int lV)
80"USAGE: setLetterplaceAttributes(R, d, b); R a ring, b,d integers
81RETURN: ring with special attributes set
82PURPOSE: sets attributes for a letterplace ring:
83@*      'isLetterplaceRing' = true, 'uptodeg' = d, 'lV' = b, where
84@*      'uptodeg' stands for the degree bound,
85@*      'lV' for the number of variables in the block 0.
86NOTE: Activate the resulting ring by using @code{setring}
87"
88{
89  if (uptodeg*lV != nvars(R))
90  {
91    ERROR("uptodeg and lV do not agree on the basering!");
92  }
93
94    // Set letterplace-specific attributes for the output ring!
95  attrib(R, "uptodeg", uptodeg);
96  attrib(R, "lV", lV);
97  attrib(R, "isLetterplaceRing", 1);
98  return (R);
99}
100example
101{
102  "EXAMPLE:"; echo = 2;
103  ring r = 0,(x(1),y(1),x(2),y(2),x(3),y(3),x(4),y(4)),dp;
104  def R = setLetterplaceAttributes(r, 4, 2); setring R;
105  attrib(R,"isLetterplaceRing");
106  lieBracket(x(1),y(1),2);
107}
108
109
110// obsolete?
111
112static proc lshift(module M, int s, string varing, def lpring)
113{
114  // FINALLY IMPLEMENTED AS A PART OT THE CODE
115  // shifts a polynomial from the ring R to s positions
116  // M lives in varing, the result in lpring
117  // to be run from varing
118  int i, j, k, sm, sv;
119  vector v;
120  //  execute("setring "+lpring);
121  setring lpring;
122  poly @@p;
123  ideal I;
124  execute("setring "+varing);
125  sm = ncols(M);
126  for (i=1; i<=s; i++)
127  {
128    // modules, e.g. free polynomials
129    for (j=1; j<=sm; j++)
130    {
131      //vectors, e.g. free monomials
132      v  = M[j];
133      sv = size(v);
134      sp = "@@p = @@p + ";
135      for (k=2; k<=sv; k++)
136      {
137        sp = sp + string(v[k])+"("+string(k-1+s)+")*";
138      }
139      sp = sp + string(v[1])+";"; // coef;
140      setring lpring;
141      //      execute("setring "+lpring);
142      execute(sp);
143      execute("setring "+varing);
144    }
145    setring lpring;
146    //    execute("setring "+lpring);
147    I = I,@@p;
148    @@p = 0;
149  }
150  setring lpring;
151  //execute("setring "+lpring);
152  export(I);
153  //  setring varing;
154  execute("setring "+varing);
155}
156
157static proc skip0(vector v)
158{
159  // skips zeros in a vector, producing another vector
160  if ( (v[1]==0) || (v==0) ) { return(vector(0)); }
161  int sv = nrows(v);
162  int sw = size(v);
163  if (sv == sw)
164  {
165    return(v);
166  }
167  int i;
168  int j=1;
169  vector w;
170  for (i=1; i<=sv; i++)
171  {
172    if (v[i] != 0)
173    {
174      w = w + v[i]*gen(j);
175      j++;
176    }
177  }
178  return(w);
179}
180
181proc lst2str(list L, list #)
182"USAGE:  lst2str(L[,n]);  L a list of modules, n an optional integer
183RETURN:  list (of strings)
184PURPOSE: convert a list (of modules) into polynomials in free algebra
185EXAMPLE: example lst2str; shows examples
186NOTE: if an optional integer is not 0, stars signs are used in multiplication
187"
188{
189  // returns a list of strings
190  // being sentences in words built from L
191  // if #[1] = 1, use * between generators
192  int useStar = 0;
193  if ( size(#)>0 )
194  {
195    if ( typeof(#[1]) != "int")
196    {
197      ERROR("Second argument of type int expected");
198    }
199    if (#[1])
200    {
201      useStar = 1;
202    }
203  }
204  int i;
205  int s    = size(L);
206  if (s<1) { return(list(""));}
207  list N;
208  for(i=1; i<=s; i++)
209  {
210    if ((typeof(L[i]) == "module") || (typeof(L[i]) == "matrix") )
211    {
212      N[i] = mod2str(L[i],useStar);
213    }
214    else
215    {
216      "module or matrix expected in the list";
217      return(N);
218    }
219  }
220  return(N);
221}
222example
223{
224  "EXAMPLE:"; echo = 2;
225  ring r = 0,(x,y,z),(dp(1),dp(2));
226  module M = [-1,x,y],[-7,y,y],[3,x,x];
227  module N = [1,x,y,x,y],[-2,y,x,y,x],[6,x,y,y,x,y];
228  list L; L[1] = M; L[2] = N;
229  lst2str(L);
230  lst2str(L[1],1);
231}
232
233
234proc mod2str(module M, list #)
235"USAGE:  mod2str(M[,n]);  M a module, n an optional integer
236RETURN:  string
237PURPOSE: convert a module into a polynomial in free algebra
238EXAMPLE: example mod2str; shows examples
239NOTE: if an optional integer is not 0, stars signs are used in multiplication
240"
241{
242  if (size(M)==0) { return(""); }
243  // returns a string
244  // a sentence in words built from M
245  // if #[1] = 1, use * between generators
246  int useStar = 0;
247  if ( size(#)>0 )
248  {
249    if ( typeof(#[1]) != "int")
250    {
251      ERROR("Second argument of type int expected");
252    }
253    if (#[1])
254    {
255      useStar = 1;
256    }
257  }
258  int i;
259  int s    = ncols(M);
260  string t;
261  string mp;
262  for(i=1; i<=s; i++)
263  {
264    mp = vct2str(M[i],useStar);
265    if (mp[1] == "-")
266    {
267      t = t + mp;
268    }
269    else
270    {
271      if (mp != "")
272      {
273         t = t + "+" + mp;
274      }
275    }
276  }
277  if (t[1]=="+")
278  {
279    t = t[2..size(t)]; // remove first "+"
280  }
281  return(t);
282}
283example
284{
285  "EXAMPLE:"; echo = 2;
286  ring r = 0,(x,y,z),(dp);
287  module M = [1,x,y,x,y],[-2,y,x,y,x],[6,x,y,y,x,y];
288  mod2str(M);
289  mod2str(M,1);
290}
291
292proc vct2str(vector v, list #)
293"USAGE:  vct2str(v[,n]);  v a vector, n an optional integer
294RETURN:  string
295PURPOSE: convert a vector into a word in free algebra
296EXAMPLE: example vct2str; shows examples
297NOTE: if an optional integer is not 0, stars signs are used in multiplication
298"
299{
300  if (v==0) { return(""); }
301  // if #[1] = 1, use * between generators
302  int useStar = 0;
303  if ( size(#)>0 )
304  {
305    if (#[1])
306    {
307      useStar = 1;
308    }
309  }
310  int ppl = printlevel-voice+2;
311  // for a word, encoded by v
312  // produces a string for it
313  v = skip0(v);
314  if (v==0) { return(string(""));}
315  number cf = leadcoef(v[1]);
316  int s = size(v);
317  string vs,vv,vp,err;
318  int i,j,p,q;
319  for (i=1; i<=s-1; i++)
320  {
321    p     = isVar(v[i+1]);
322    if (p==0)
323    {
324      err = "Error: monomial expected at nonzero position " + string(i+1);
325      ERROR(err+" in vct2str");
326      //      dbprint(ppl,err);
327      //      return("_");
328    }
329    if (p==1)
330    {
331      if (useStar && (size(vs) >0))       {   vs = vs + "*"; }
332        vs = vs + string(v[i+1]);
333    }
334    else //power
335    {
336      vv = string(v[i+1]);
337      q = find(vv,"^");
338      if (q==0)
339      {
340        q = find(vv,string(p));
341        if (q==0)
342        {
343          err = "error in find for string "+vv;
344          dbprint(ppl,err);
345          return("_");
346        }
347      }
348      // q>0
349      vp = vv[1..q-1];
350      for(j=1;j<=p;j++)
351      {
352         if (useStar && (size(vs) >0))       {   vs = vs + "*"; }
353         vs = vs + vp;
354      }
355    }
356  }
357  string scf;
358  if (cf == -1)
359  {
360    scf = "-";
361  }
362  else
363  {
364    scf = string(cf);
365    if ( (cf == 1) && (size(vs)>0) )
366    {
367      scf = "";
368    }
369  }
370  if (useStar && (size(scf) >0) && (scf!="-") )       {   scf = scf + "*"; }
371  vs = scf + vs;
372  return(vs);
373}
374example
375{
376  "EXAMPLE:"; echo = 2;
377  ring r = (0,a),(x,y3,z(1)),dp;
378  vector v = [-7,x,y3^4,x2,z(1)^3];
379  vct2str(v);
380  vct2str(v,1);
381  vector w = [-7a^5+6a,x,y3,y3,x,z(1),z(1)];
382  vct2str(w);
383  vct2str(w,1);
384}
385
386proc isVar(poly p)
387"USAGE:  isVar(p);  poly p
388RETURN:  int
389PURPOSE: check, whether leading monomial of p is a power of a single variable
390@* from the basering. Returns the exponent or 0 if p is multivariate.
391EXAMPLE: example isVar; shows examples
392"
393{
394  // checks whether p is a variable indeed
395  // if it's a power of a variable, returns the power
396  if (p==0) {  return(0); } //"p=0";
397  poly q   = leadmonom(p);
398  if ( (p-lead(p)) !=0 ) { return(0); } // "p-lm(p)>0";
399  intvec v = leadexp(p);
400  int s = size(v);
401  int i=1;
402  int cnt = 0;
403  int pwr = 0;
404  for (i=1; i<=s; i++)
405  {
406    if (v[i] != 0)
407    {
408      cnt++;
409      pwr = v[i];
410    }
411  }
412  //  "cnt:";  cnt;
413  if (cnt==1) { return(pwr); }
414  else { return(0); }
415}
416example
417{
418  "EXAMPLE:"; echo = 2;
419  ring r = 0,(x,y),dp;
420  poly f = xy+1;
421  isVar(f);
422  poly g = y^3;
423  isVar(g);
424  poly h = 7*x^3;
425  isVar(h);
426  poly i = 1;
427  isVar(i);
428}
429
430// new conversion routines
431
432static proc id2words(ideal I, int d)
433{
434  // NOT FINISHED
435  // input: ideal I of polys in letter-place notation
436  // in the ring with d real vars
437  // output: the list of strings: associative words
438  // extract names of vars
439  int i,m,n; string s; string place = "(1)";
440  list lv;
441  for(i=1; i<=d; i++)
442  {
443    s = string(var(i));
444    // get rid of place
445    n = find(s, place);
446    if (n>0)
447    {
448      s = s[1..n-1];
449    }
450    lv[i] = s;
451  }
452  poly p,q;
453  for (i=1; i<=ncols(I); i++)
454  {
455    if (I[i] != 0)
456    {
457      p = I[i];
458      while (p!=0)
459      {
460         q = leadmonom(p);
461      }
462    }
463  }
464
465  return(lv);
466}
467example
468{
469  "EXAMPLE:"; echo = 2;
470  ring r = 0,(x(1),y(1),z(1),x(2),y(2),z(2)),dp;
471  ideal I = x(1)*y(2) -z(1)*x(2);
472  id2words(I,3);
473}
474
475static proc mono2word(poly p, int d)
476{
477}
478
479proc letplaceGBasis(ideal I)
480"USAGE: letplaceGBasis(I);  I an ideal
481RETURN: ideal
482ASSUME: basering is a Letterplace ring, an ideal consists of Letterplace
483@*      polynomials
484PURPOSE: compute the two-sided Groebner basis of an ideal I via Letterplace
485@*       algorithm
486NOTE: the degree bound for this computation is read off the letterplace
487@*    structure of basering
488EXAMPLE: example letplaceGBasis; shows examples
489"
490{
491  if (lpAssumeViolation())
492  {
493    ERROR("Incomplete Letterplace structure on the basering!");
494  }
495  int ppl = printlevel-voice+2;
496  def save = basering;
497  // assumes of the ring have been checked
498  // run the computation - it will test assumes on the ideal
499  int uptodeg = attrib(save,"uptodeg");
500  int lV = attrib(save,"lV");
501  dbprint(ppl,"start computing GB");
502  ideal J = system("freegb",I,uptodeg,lV);
503  dbprint(ppl,"finished computing GB");
504  dbprint(ppl-1,"the result is:");
505  dbprint(ppl-1,J);
506  return(J);
507}
508example
509{
510  "EXAMPLE:"; echo = 2;
511  ring r = 0,(x,y,z),(dp(1),dp(2));
512  int degree_bound = 5;
513  def R = makeLetterplaceRing(5);
514  setring R;
515  ideal I = -x(1)*y(2)-7*y(1)*y(2)+3*x(1)*x(2), x(1)*y(2)*x(3)-y(1)*x(2)*y(3);
516  ideal J = letplaceGBasis(I);
517  J;
518  // now transfom letterplace polynomials into strings of words
519  lp2lstr(J,r); // export an object called @code{@LN} to the ring r
520  setring r;  // change to the ring r
521  lst2str(@LN,1);
522}
523
524// given the element -7xy^2x, it is represented as [-7,x,y^2,x] or as [-7,x,y,y,x]
525// use the orig ord on (x,y,z) and expand it blockwise to (x(i),y(i),z(i))
526
527// the correspondences:
528// monomial in K<x,y,z>    <<--->> vector in R
529// polynomial in K<x,y,z>  <<--->> list of vectors (matrix/module) in R
530// ideal in K<x,y,z>       <<--->> list of matrices/modules in R
531
532
533// 1. form a new ring
534// 2. NOP
535// 3. compute GB -> with the kernel stuff
536// 4. skip shifted elts (check that no such exist?)
537// 5. go back to orig vars, produce strings/modules
538// 6. return the result
539
540proc freeGBasis(list LM, int d)
541"USAGE:  freeGBasis(L, d);  L a list of modules, d an integer
542RETURN:  ring
543ASSUME: L has a special form. Namely, it is a list of modules, where
544
545 - each generator of every module stands for a monomial times coefficient in
546@* free algebra,
547
548 - in such a vector generator, the 1st entry is a nonzero coefficient from the
549@* ground field
550
551 - and each next entry hosts a variable from the basering.
552PURPOSE: compute the two-sided Groebner basis of an ideal, encoded by L
553@* in the free associative algebra, up to degree d
554NOTE: Apply @code{lst2str} to the output in order to obtain a better readable
555@*    presentation
556EXAMPLE: example freeGBasis; shows examples
557"
558{
559  // d = up to degree, will be shifted to d+1
560  if (d<1) {"bad d"; return(0);}
561
562  int ppl = printlevel-voice+2;
563  string err = "";
564
565  int i,j,s;
566  def save = basering;
567  // determine max no of places in the input
568  int slm = size(LM); // numbers of polys in the ideal
569  int sm;
570  intvec iv;
571  module M;
572  for (i=1; i<=slm; i++)
573  {
574    // modules, e.g. free polynomials
575    M  = LM[i];
576    sm = ncols(M);
577    for (j=1; j<=sm; j++)
578    {
579      //vectors, e.g. free monomials
580      iv = iv, size(M[j])-1; // 1 place is reserved by the coeff
581    }
582  }
583  int D = Max(iv); // max size of input words
584  if (d<D) {"bad d"; return(LM);}
585  D = D + d-1;
586  //  D = d;
587  list LR  = ringlist(save);
588  list L, tmp;
589  L[1] = LR[1]; // ground field
590  L[4] = LR[4]; // quotient ideal
591  tmp  = LR[2]; // varnames
592  s = size(LR[2]);
593  for (i=1; i<=D; i++)
594  {
595    for (j=1; j<=s; j++)
596    {
597      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
598    }
599  }
600  for (i=1; i<=s; i++)
601  {
602    tmp[i] = string(tmp[i])+"("+string(1)+")";
603  }
604  L[2] = tmp;
605  list OrigNames = LR[2];
606  // ordering: d blocks of the ord on r
607  // try to get whether the ord on r is blockord itself
608  s = size(LR[3]);
609  if (s==2)
610  {
611    // not a blockord, 1 block + module ord
612    tmp = LR[3][s]; // module ord
613    for (i=1; i<=D; i++)
614    {
615      LR[3][s-1+i] = LR[3][1];
616    }
617    LR[3][s+D] = tmp;
618  }
619  if (s>2)
620  {
621    // there are s-1 blocks
622    int nb = s-1;
623    tmp = LR[3][s]; // module ord
624    for (i=1; i<=D; i++)
625    {
626      for (j=1; j<=nb; j++)
627      {
628        LR[3][i*nb+j] = LR[3][j];
629      }
630    }
631    //    size(LR[3]);
632    LR[3][nb*(D+1)+1] = tmp;
633  }
634  L[3] = LR[3];
635  def @R = ring(L);
636  setring @R;
637  ideal I;
638  poly @p;
639  s = size(OrigNames);
640  //  "s:";s;
641  // convert LM to canonical vectors (no powers)
642  setring save;
643  kill M; // M was defined earlier
644  module M;
645  slm = size(LM); // numbers of polys in the ideal
646  int sv,k,l;
647  vector v;
648  //  poly p;
649  string sp;
650  setring @R;
651  poly @@p=0;
652  setring save;
653  for (l=1; l<=slm; l++)
654  {
655    // modules, e.g. free polynomials
656    M  = LM[l];
657    sm = ncols(M); // in intvec iv the sizes are stored
658    // modules, e.g. free polynomials
659    for (j=1; j<=sm; j++)
660    {
661      //vectors, e.g. free monomials
662      v  = M[j];
663      sv = size(v);
664      //        "sv:";sv;
665      sp = "@@p = @@p + ";
666      for (k=2; k<=sv; k++)
667      {
668        sp = sp + string(v[k])+"("+string(k-1)+")*";
669      }
670      sp = sp + string(v[1])+";"; // coef;
671      setring @R;
672      execute(sp);
673      setring save;
674    }
675    setring @R;
676    //      "@@p:"; @@p;
677    I = I,@@p;
678    @@p = 0;
679    setring save;
680  }
681  kill sp;
682  // 3. compute GB
683  setring @R;
684  dbprint(ppl,"computing GB");
685  ideal J = system("freegb",I,d,nvars(save));
686  //  ideal J = slimgb(I);
687  dbprint(ppl,J);
688  // 4. skip shifted elts
689  ideal K = select1(J,1..s); // s = size(OrigNames)
690  dbprint(ppl,K);
691  dbprint(ppl, "done with GB");
692  // K contains vars x(1),...z(1) = images of originals
693  // 5. go back to orig vars, produce strings/modules
694  if (K[1] == 0)
695  {
696    "no reasonable output, GB gives 0";
697    return(0);
698  }
699  int sk = size(K);
700  int sp, sx, a, b;
701  intvec x;
702  poly p,q;
703  poly pn;
704  // vars in 'save'
705  setring save;
706  module N;
707  list LN;
708  vector V;
709  poly pn;
710  // test and skip exponents >=2
711  setring @R;
712  for(i=1; i<=sk; i++)
713  {
714    p  = K[i];
715    while (p!=0)
716    {
717      q  = lead(p);
718      //      "processing q:";q;
719      x  = leadexp(q);
720      sx = size(x);
721      for(k=1; k<=sx; k++)
722      {
723        if ( x[k] >= 2 )
724        {
725          err = "skip: the value x[k] is " + string(x[k]);
726          dbprint(ppl,err);
727          //            return(0);
728          K[i] = 0;
729          p    = 0;
730          q    = 0;
731          break;
732        }
733      }
734      p  = p - q;
735    }
736  }
737  K  = simplify(K,2);
738  sk = size(K);
739  for(i=1; i<=sk; i++)
740  {
741    //    setring save;
742    //    V  = 0;
743    setring @R;
744    p  = K[i];
745    while (p!=0)
746    {
747      q  = lead(p);
748      err =  "processing q:" + string(q);
749      dbprint(ppl,err);
750      x  = leadexp(q);
751      sx = size(x);
752      pn = leadcoef(q);
753      setring save;
754      pn = imap(@R,pn);
755      V  = V + leadcoef(pn)*gen(1);
756      for(k=1; k<=sx; k++)
757      {
758        if (x[k] ==1)
759        {
760          a = k div s; // block number=a+1, a!=0
761          b = k % s; // remainder
762          //          printf("a: %s, b: %s",a,b);
763          if (b == 0)
764          {
765            // that is it's the last var in the block
766            b = s;
767            a = a-1;
768          }
769          V = V + var(b)*gen(a+2);
770        }
771//         else
772//         {
773//           printf("error: the value x[k] is %s", x[k]);
774//           return(0);
775//         }
776      }
777      err = "V: " + string(V);
778      dbprint(ppl,err);
779      //      printf("V: %s", string(V));
780      N = N,V;
781      V  = 0;
782      setring @R;
783      p  = p - q;
784      pn = 0;
785    }
786    setring save;
787    LN[i] = simplify(N,2);
788    N     = 0;
789  }
790  setring save;
791  return(LN);
792}
793example
794{
795  "EXAMPLE:"; echo = 2;
796  ring r = 0,(x,y,z),(dp(1),dp(2)); //  ring r = 0,(x,y,z),(a(3,0,2), dp(2));
797  module M = [-1,x,y],[-7,y,y],[3,x,x]; // stands for free poly -xy - 7yy - 3xx
798  module N = [1,x,y,x],[-1,y,x,y]; // stands for free poly xyx - yxy
799  list L; L[1] = M; L[2] = N; // list of modules stands for an ideal in free algebra
800  lst2str(L); // list to string conversion of input polynomials
801  def U = freeGBasis(L,5); // 5 is the degree bound
802  lst2str(U);
803}
804
805static proc crs(list LM, int d)
806"USAGE:  crs(L, d);  L a list of modules, d an integer
807RETURN:  ring
808PURPOSE: create a ring and shift the ideal
809EXAMPLE: example crs; shows examples
810"
811{
812  // d = up to degree, will be shifted to d+1
813  if (d<1) {"bad d"; return(0);}
814
815  int ppl = printlevel-voice+2;
816  string err = "";
817
818  int i,j,s;
819  def save = basering;
820  // determine max no of places in the input
821  int slm = size(LM); // numbers of polys in the ideal
822  int sm;
823  intvec iv;
824  module M;
825  for (i=1; i<=slm; i++)
826  {
827    // modules, e.g. free polynomials
828    M  = LM[i];
829    sm = ncols(M);
830    for (j=1; j<=sm; j++)
831    {
832      //vectors, e.g. free monomials
833      iv = iv, size(M[j])-1; // 1 place is reserved by the coeff
834    }
835  }
836  int D = Max(iv); // max size of input words
837  if (d<D) {"bad d"; return(LM);}
838  D = D + d-1;
839  //  D = d;
840  list LR  = ringlist(save);
841  list L, tmp;
842  L[1] = LR[1]; // ground field
843  L[4] = LR[4]; // quotient ideal
844  tmp  = LR[2]; // varnames
845  s = size(LR[2]);
846  for (i=1; i<=D; i++)
847  {
848    for (j=1; j<=s; j++)
849    {
850      tmp[i*s+j] = string(tmp[j])+"("+string(i)+")";
851    }
852  }
853  for (i=1; i<=s; i++)
854  {
855    tmp[i] = string(tmp[i])+"("+string(0)+")";
856  }
857  L[2] = tmp;
858  list OrigNames = LR[2];
859  // ordering: d blocks of the ord on r
860  // try to get whether the ord on r is blockord itself
861  s = size(LR[3]);
862  if (s==2)
863  {
864    // not a blockord, 1 block + module ord
865    tmp = LR[3][s]; // module ord
866    for (i=1; i<=D; i++)
867    {
868      LR[3][s-1+i] = LR[3][1];
869    }
870    LR[3][s+D] = tmp;
871  }
872  if (s>2)
873  {
874    // there are s-1 blocks
875    int nb = s-1;
876    tmp = LR[3][s]; // module ord
877    for (i=1; i<=D; i++)
878    {
879      for (j=1; j<=nb; j++)
880      {
881        LR[3][i*nb+j] = LR[3][j];
882      }
883    }
884    //    size(LR[3]);
885    LR[3][nb*(D+1)+1] = tmp;
886  }
887  L[3] = LR[3];
888  def @R = ring(L);
889  setring @R;
890  ideal I;
891  poly @p;
892  s = size(OrigNames);
893  //  "s:";s;
894  // convert LM to canonical vectors (no powers)
895  setring save;
896  kill M; // M was defined earlier
897  module M;
898  slm = size(LM); // numbers of polys in the ideal
899  int sv,k,l;
900  vector v;
901  //  poly p;
902  string sp;
903  setring @R;
904  poly @@p=0;
905  setring save;
906  for (l=1; l<=slm; l++)
907  {
908    // modules, e.g. free polynomials
909    M  = LM[l];
910    sm = ncols(M); // in intvec iv the sizes are stored
911    for (i=0; i<=d-iv[l]; i++)
912    {
913      // modules, e.g. free polynomials
914      for (j=1; j<=sm; j++)
915      {
916        //vectors, e.g. free monomials
917        v  = M[j];
918        sv = size(v);
919        //        "sv:";sv;
920        sp = "@@p = @@p + ";
921        for (k=2; k<=sv; k++)
922        {
923          sp = sp + string(v[k])+"("+string(k-2+i)+")*";
924        }
925        sp = sp + string(v[1])+";"; // coef;
926        setring @R;
927        execute(sp);
928        setring save;
929      }
930      setring @R;
931      //      "@@p:"; @@p;
932      I = I,@@p;
933      @@p = 0;
934      setring save;
935    }
936  }
937  setring @R;
938  export I;
939  return(@R);
940}
941example
942{
943  "EXAMPLE:"; echo = 2;
944  ring r = 0,(x,y,z),(dp(1),dp(2));
945  module M = [-1,x,y],[-7,y,y],[3,x,x];
946  module N = [1,x,y,x],[-1,y,x,y];
947  list L; L[1] = M; L[2] = N;
948  lst2str(L);
949  def U = crs(L,5);
950  setring U; U;
951  I;
952}
953
954static proc polylen(ideal I)
955{
956  // returns the ideal of length of polys
957  int i;
958  intvec J;
959  number s = 0;
960  for(i=1;i<=ncols(I);i++)
961  {
962    J[i] = size(I[i]);
963    s = s + J[i];
964  }
965  printf("the sum of length %s",s);
966  //  print(s);
967  return(J);
968}
969
970// new: uniting both mLR1 (homog) and mLR2 (nonhomog)
971proc makeLetterplaceRing(int d, list #)
972"USAGE:  makeLetterplaceRing(d [,h]); d an integer, h an optional integer
973RETURN:  ring
974PURPOSE: creates a ring with the ordering, used in letterplace computations
975NOTE: h = 0 (default) : Dp ordering will be used
976h = 2 : weights 1 used for all the variables, a tie breaker is a list of block of original ring
977h = 1 : the pure homogeneous letterplace block ordering (applicable in the situation of homogeneous input ideals) will be used.
978EXAMPLE: example makeLetterplaceRing; shows examples
979"
980{
981  int alternativeVersion = 0;
982  if ( size(#)>0 )
983  {
984    if (typeof(#[1]) == "int")
985    {
986      alternativeVersion = #[1];
987    }
988  }
989  if (alternativeVersion == 1)
990  {
991    def @A = makeLetterplaceRing1(d);
992  }
993  else {
994    if (alternativeVersion == 2)
995    {
996      def @A = makeLetterplaceRing2(d);
997    }
998    else {
999      def @A = makeLetterplaceRing4(d);
1000    }
1001  }
1002  return(@A);
1003}
1004example
1005{
1006  "EXAMPLE:"; echo = 2;
1007  ring r = 0,(x,y,z),(dp(1),dp(2));
1008  def A = makeLetterplaceRing(2);
1009  setring A;  A;
1010  attrib(A,"isLetterplaceRing");
1011  attrib(A,"uptodeg");  // degree bound
1012  attrib(A,"lV"); // number of variables in the main block
1013  setring r; def B = makeLetterplaceRing(2,1); // to compare:
1014  setring B;  B;
1015}
1016
1017//proc freegbRing(int d)
1018static proc makeLetterplaceRing1(int d)
1019"USAGE:  makeLetterplaceRing1(d); d an integer
1020RETURN:  ring
1021PURPOSE: creates a ring with a special ordering, suitable for
1022@* the use of homogeneous letterplace (d blocks of shifted original variables)
1023EXAMPLE: example makeLetterplaceRing1; shows examples
1024"
1025{
1026  // d = up to degree, will be shifted to d+1
1027  if (d<1) {"bad d"; return(0);}
1028
1029  int uptodeg = d; int lV = nvars(basering);
1030
1031  int ppl = printlevel-voice+2;
1032  string err = "";
1033
1034  int i,j,s;
1035  def save = basering;
1036  int D = d-1;
1037  list LR  = ringlist(save);
1038  list L, tmp;
1039  L[1] = LR[1]; // ground field
1040  L[4] = LR[4]; // quotient ideal
1041  tmp  = LR[2]; // varnames
1042  s = size(LR[2]);
1043  for (i=1; i<=D; i++)
1044  {
1045    for (j=1; j<=s; j++)
1046    {
1047      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
1048    }
1049  }
1050  for (i=1; i<=s; i++)
1051  {
1052    tmp[i] = string(tmp[i])+"("+string(1)+")";
1053  }
1054  L[2] = tmp;
1055  list OrigNames = LR[2];
1056  // ordering: d blocks of the ord on r
1057  // try to get whether the ord on r is blockord itself
1058  // TODO: make L(2) ordering! exponent is maximally 2
1059  s = size(LR[3]);
1060  if (s==2)
1061  {
1062    // not a blockord, 1 block + module ord
1063    tmp = LR[3][s]; // module ord
1064    for (i=1; i<=D; i++)
1065    {
1066      LR[3][s-1+i] = LR[3][1];
1067    }
1068    LR[3][s+D] = tmp;
1069  }
1070  if (s>2)
1071  {
1072    // there are s-1 blocks
1073    int nb = s-1;
1074    tmp = LR[3][s]; // module ord
1075    for (i=1; i<=D; i++)
1076    {
1077      for (j=1; j<=nb; j++)
1078      {
1079        LR[3][i*nb+j] = LR[3][j];
1080      }
1081    }
1082    //    size(LR[3]);
1083    LR[3][nb*(D+1)+1] = tmp;
1084  }
1085  L[3] = LR[3];
1086  attrib(L,"maxExp",1);
1087  def @R = ring(L);
1088  //  setring @R;
1089  //  int uptodeg = d; int lV = nvars(basering); // were defined before
1090  def @@R = setLetterplaceAttributes(@R,uptodeg,lV);
1091  return (@@R);
1092}
1093example
1094{
1095  "EXAMPLE:"; echo = 2;
1096  ring r = 0,(x,y,z),(dp(1),dp(2));
1097  def A = makeLetterplaceRing1(2);
1098  setring A;
1099  A;
1100  attrib(A,"isLetterplaceRing");
1101  attrib(A,"uptodeg");  // degree bound
1102  attrib(A,"lV"); // number of variables in the main block
1103}
1104
1105static proc makeLetterplaceRing2(int d)
1106"USAGE:  makeLetterplaceRing2(d); d an integer
1107RETURN:  ring
1108PURPOSE: creates a ring with a special ordering, suitable for
1109@* the use of non-homogeneous letterplace
1110NOTE: the matrix for the ordering looks as follows: first row is 1,1,...,1
1111@* then there come 'd' blocks of shifted original variables
1112EXAMPLE: example makeLetterplaceRing2; shows examples
1113"
1114{
1115
1116  // ToDo future: inherit positive weights in the orig ring
1117  // complain on nonpositive ones
1118
1119  // d = up to degree, will be shifted to d+1
1120  if (d<1) {"bad d"; return(0);}
1121
1122  int uptodeg = d; int lV = nvars(basering);
1123
1124  int ppl = printlevel-voice+2;
1125  string err = "";
1126
1127  int i,j,s;
1128  def save = basering;
1129  int D = d-1;
1130  list LR  = ringlist(save);
1131  list L, tmp, tmp2, tmp3;
1132  L[1] = LR[1]; // ground field
1133  L[4] = LR[4]; // quotient ideal
1134  tmp  = LR[2]; // varnames
1135  s = size(LR[2]);
1136  for (i=1; i<=D; i++)
1137  {
1138    for (j=1; j<=s; j++)
1139    {
1140      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
1141    }
1142  }
1143  for (i=1; i<=s; i++)
1144  {
1145    tmp[i] = string(tmp[i])+"("+string(1)+")";
1146  }
1147  L[2] = tmp;
1148  list OrigNames = LR[2];
1149  // ordering: one 1..1 a above
1150  // ordering: d blocks of the ord on r
1151  // try to get whether the ord on r is blockord itself
1152  // TODO: make L(2) ordering! exponent is maximally 2
1153  s = size(LR[3]);
1154  if (s==2)
1155  {
1156    // not a blockord, 1 block + module ord
1157    tmp = LR[3][s]; // module ord
1158    for (i=1; i<=d; i++)
1159    {
1160      LR[3][s-1+i] = LR[3][1];
1161    }
1162    //    LR[3][s+D] = tmp;
1163    LR[3][s+1+D] = tmp;
1164    LR[3][1] = list("a",intvec(1: int(d*lV))); // deg-ord
1165  }
1166  if (s>2)
1167  {
1168    // there are s-1 blocks
1169    int nb = s-1;
1170    tmp = LR[3][s]; // module ord to place at the very end
1171    tmp2 = LR[3]; tmp2 = tmp2[1..nb];
1172    tmp3[1] = list("a",intvec(1: int(d*lV))); // deg-ord, insert as the 1st
1173    for (i=1; i<=d; i++)
1174    {
1175      tmp3 = tmp3 + tmp2;
1176    }
1177    tmp3 = tmp3 + list(tmp);
1178    LR[3] = tmp3;
1179//     for (i=1; i<=d; i++)
1180//     {
1181//       for (j=1; j<=nb; j++)
1182//       {
1183//         //        LR[3][i*nb+j+1]= LR[3][j];
1184//         LR[3][i*nb+j+1]= tmp2[j];
1185//       }
1186//     }
1187//     //    size(LR[3]);
1188//     LR[3][(s-1)*d+2] = tmp;
1189//     LR[3] = list("a",intvec(1: int(d*lV))) + LR[3]; // deg-ord, insert as the 1st
1190    // remove everything behind nb*(D+1)+1 ?
1191    //    tmp = LR[3];
1192    //    LR[3] = tmp[1..size(tmp)-1];
1193  }
1194  L[3] = LR[3];
1195  attrib(L,"maxExp",1);
1196  def @R = ring(L);
1197  //  setring @R;
1198  //  int uptodeg = d; int lV = nvars(basering); // were defined before
1199  def @@R = setLetterplaceAttributes(@R,uptodeg,lV);
1200  return (@@R);
1201}
1202example
1203{
1204  "EXAMPLE:"; echo = 2;
1205  ring r = 0,(x,y,z),(dp(1),dp(2));
1206  def A = makeLetterplaceRing2(2);
1207  setring A;
1208  A;
1209  attrib(A,"isLetterplaceRing");
1210  attrib(A,"uptodeg");  // degree bound
1211  attrib(A,"lV"); // number of variables in the main block
1212}
1213
1214static proc makeLetterplaceRing4(int d)
1215"USAGE:  makeLetterplaceRing2(d); d an integer
1216RETURN:  ring
1217PURPOSE: creates a Letterplace ring with a Dp ordering, suitable for
1218@* the use of non-homogeneous letterplace
1219NOTE: the matrix for the ordering looks as follows: first row is 1,1,...,1
1220EXAMPLE: example makeLetterplaceRing2; shows examples
1221"
1222{
1223
1224  // ToDo future: inherit positive weights in the orig ring
1225  // complain on nonpositive ones
1226
1227  // d = up to degree, will be shifted to d+1
1228  if (d<1) {"bad d"; return(0);}
1229
1230  int uptodeg = d; int lV = nvars(basering);
1231
1232  int ppl = printlevel-voice+2;
1233  string err = "";
1234
1235  int i,j,s;
1236  def save = basering;
1237  int D = d-1;
1238  list LR  = ringlist(save);
1239  list L, tmp, tmp2, tmp3;
1240  L[1] = LR[1]; // ground field
1241  L[4] = LR[4]; // quotient ideal
1242  tmp  = LR[2]; // varnames
1243  s = size(LR[2]);
1244  for (i=1; i<=D; i++)
1245  {
1246    for (j=1; j<=s; j++)
1247    {
1248      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
1249    }
1250  }
1251  for (i=1; i<=s; i++)
1252  {
1253    tmp[i] = string(tmp[i])+"("+string(1)+")";
1254  }
1255  L[2] = tmp;
1256  list OrigNames = LR[2];
1257
1258  s = size(LR[3]);
1259  list ordering;
1260  ordering[1] = list("Dp",intvec(1: int(d*lV)));
1261  ordering[2] = LR[3][s]; // module ord to place at the very end
1262  LR[3] = ordering;
1263
1264  L[3] = LR[3];
1265  attrib(L,"maxExp",1);
1266  def @R = ring(L);
1267  def @@R = setLetterplaceAttributes(@R,uptodeg,lV);
1268  return (@@R);
1269}
1270example
1271{
1272  "EXAMPLE:"; echo = 2;
1273  ring r = 0,(x,y,z),(dp(1),dp(2));
1274  def A = makeLetterplaceRing2(2);
1275  setring A;
1276  A;
1277  attrib(A,"isLetterplaceRing");
1278  attrib(A,"uptodeg");  // degree bound
1279  attrib(A,"lV"); // number of variables in the main block
1280}
1281
1282// P[s;sigma] approach
1283static proc makeLetterplaceRing3(int d)
1284"USAGE:  makeLetterplaceRing1(d); d an integer
1285RETURN:  ring
1286PURPOSE: creates a ring with a special ordering, representing
1287@* the original P[s,sigma] (adds d blocks of shifted s)
1288ASSUME: basering is a letterplace ring
1289NOTE: experimental status
1290EXAMPLE: example makeLetterplaceRing1; shows examples
1291"
1292{
1293  // d = up to degree, will be shifted to d+1
1294  if (d<1) {"bad d"; return(0);}
1295
1296  int uptodeg = d; int lV = nvars(basering);
1297
1298  int ppl = printlevel-voice+2;
1299  string err = "";
1300
1301  int i,j,s;
1302  def save = basering;
1303  int D = d-1;
1304  list LR  = ringlist(save);
1305  list L, tmp;
1306  L[1] = LR[1]; // ground field
1307  L[4] = LR[4]; // quotient ideal
1308  tmp  = LR[2]; // varnames
1309  tmp[size(tmp)+1] = "s";
1310  // add s's
1311  //  string newSname = "@s";
1312  s = size(LR[2]);
1313  for (i=1; i<=D; i++)
1314  {
1315    for (j=1; j<=s; j++)
1316    {
1317      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
1318    }
1319  }
1320  // the final index is D*s+s = (D+1)*s = degBound*s
1321  for (i=1; i<=d; i++)
1322  {
1323    tmp[FIndex + i] =  string(newSname)+"("+string(i)+")";
1324  }
1325  L[2] = tmp;
1326  list OrigNames = LR[2];
1327  // ordering: d blocks of the MODIFIED ord on r
1328  // try to get whether the ord on r is blockord itself
1329  // TODO: make L(2) ordering! exponent is maximally 2
1330  s = size(LR[3]);
1331
1332  // assume: basering was a letterplace, so get its block
1333  tmp = LR[3][1]; // ASSUME: it's a nice block
1334  // modify it
1335  // add (0,..,0,1) ... as antiblock part
1336  intvec iv; list ttmp, tmp1;
1337  for (i=1; i<=d; i++)
1338  {
1339    // the position to hold 1:
1340    iv = intvec( gen( i*(lV+1)-1 ) );
1341    ttmp[1] = "a";
1342    ttmp[2] = iv;
1343    tmp1[i] = ttmp;
1344  }
1345  // finished: antiblock part //TOCONTINUE
1346
1347  if (s==2)
1348  {
1349    // not a blockord, 1 block + module ord
1350    tmp = LR[3][s]; // module ord
1351    for (i=1; i<=D; i++)
1352    {
1353      LR[3][s-1+i] = LR[3][1];
1354    }
1355    LR[3][s+D] = tmp;
1356  }
1357  if (s>2)
1358  {
1359    // there are s-1 blocks
1360    int nb = s-1;
1361    tmp = LR[3][s]; // module ord
1362    for (i=1; i<=D; i++)
1363    {
1364      for (j=1; j<=nb; j++)
1365      {
1366        LR[3][i*nb+j] = LR[3][j];
1367      }
1368    }
1369    //    size(LR[3]);
1370    LR[3][nb*(D+1)+1] = tmp;
1371  }
1372  L[3] = LR[3];
1373  attrib(L,"maxExp",1);
1374  def @R = ring(L);
1375  //  setring @R;
1376  //  int uptodeg = d; int lV = nvars(basering); // were defined before
1377  def @@R = setLetterplaceAttributes(@R,uptodeg,lV);
1378  return (@@R);
1379}
1380example
1381{
1382  "EXAMPLE:"; echo = 2;
1383  ring r = 0,(x,y,z),(dp(1),dp(2));
1384  def A = makeLetterplaceRing3(2);
1385  setring A;
1386  A;
1387  attrib(A,"isLetterplaceRing");
1388  attrib(A,"uptodeg");  // degree bound
1389  attrib(A,"lV"); // number of variables in the main block
1390}
1391
1392/* EXAMPLES:
1393
1394//static proc ex_shift()
1395{
1396  LIB "freegb.lib";
1397  ring r = 0,(x,y,z),(dp(1),dp(2));
1398  module M = [-1,x,y],[-7,y,y],[3,x,x];
1399  module N = [1,x,y,x],[-1,y,x,y];
1400  list L; L[1] = M; L[2] = N;
1401  lst2str(L);
1402  def U = crs(L,5);
1403  setring U; U;
1404  I;
1405  poly p = I[2]; // I[8];
1406  p;
1407  system("stest",p,7,7,3); // error -> the world is ok
1408  poly q1 = system("stest",p,1,7,3); //ok
1409  poly q6 = system("stest",p,6,7,3); //ok
1410  system("btest",p,3); //ok
1411  system("btest",q1,3); //ok
1412  system("btest",q6,3); //ok
1413}
1414
1415//static proc test_shrink()
1416{
1417  LIB "freegb.lib";
1418  ring r =0,(x,y,z),dp;
1419  int d = 5;
1420  def R = makeLetterplaceRing(d);
1421  setring R;
1422  poly p1 = x(1)*y(2)*z(3);
1423  poly p2 = x(1)*y(4)*z(5);
1424  poly p3 = x(1)*y(1)*z(3);
1425  poly p4 = x(1)*y(2)*z(2);
1426  poly p5 = x(3)*z(5);
1427  poly p6 = x(1)*y(1)*x(3)*z(5);
1428  poly p7 = x(1)*y(2)*x(3)*y(4)*z(5);
1429  poly p8 = p1+p2+p3+p4+p5 + p6 + p7;
1430  p1; system("shrinktest",p1,3);
1431  p2; system("shrinktest",p2,3);
1432  p3; system("shrinktest",p3,3);
1433  p4; system("shrinktest",p4,3);
1434  p5; system("shrinktest",p5,3);
1435  p6; system("shrinktest",p6,3);
1436  p7; system("shrinktest",p7,3);
1437  p8; system("shrinktest",p8,3);
1438  poly p9 = p1 + 2*p2 + 5*p5 + 7*p7;
1439  p9; system("shrinktest",p9,3);
1440}
1441
1442//static proc ex2()
1443{
1444  option(prot);
1445  LIB "freegb.lib";
1446  ring r = 0,(x,y),dp;
1447  module M = [-1,x,y],[3,x,x]; // 3x^2 - xy
1448  def U = freegb(M,7);
1449  lst2str(U);
1450}
1451
1452//static proc ex_nonhomog()
1453{
1454  option(prot);
1455  LIB "freegb.lib";
1456  ring r = 0,(x,y,h),dp;
1457  list L;
1458  module M;
1459  M = [-1,y,y],[1,x,x,x];  // x3-y2
1460  L[1] = M;
1461  M = [1,x,h],[-1,h,x];  // xh-hx
1462  L[2] = M;
1463  M = [1,y,h],[-1,h,y];  // yh-hy
1464  L[3] = M;
1465  def U = freegb(L,4);
1466  lst2str(U);
1467  // strange elements in the basis
1468}
1469
1470//static proc ex_nonhomog_comm()
1471{
1472  option(prot);
1473  LIB "freegb.lib";
1474  ring r = 0,(x,y),dp;
1475  module M = [-1,y,y],[1,x,x,x];
1476  def U = freegb(M,5);
1477  lst2str(U);
1478}
1479
1480//static proc ex_nonhomog_h()
1481{
1482  option(prot);
1483  LIB "freegb.lib";
1484  ring r = 0,(x,y,h),(a(1,1),dp);
1485  module M = [-1,y,y,h],[1,x,x,x]; // x3 - y2h
1486  def U = freegb(M,6);
1487  lst2str(U);
1488}
1489
1490//static proc ex_nonhomog_h2()
1491{
1492  option(prot);
1493  LIB "freegb.lib";
1494  ring r = 0,(x,y,h),(dp);
1495  list L;
1496  module M;
1497  M = [-1,y,y,h],[1,x,x,x]; // x3 - y2h
1498  L[1] = M;
1499  M = [1,x,h],[-1,h,x]; // xh - hx
1500  L[2] = M;
1501  M = [1,y,h],[-1,h,y]; // yh - hy
1502  L[3] = M;
1503  def U = freeGBasis(L,3);
1504  lst2str(U);
1505  // strange answer CHECK
1506}
1507
1508
1509//static proc ex_nonhomog_3()
1510{
1511  option(prot);
1512  LIB "./freegb.lib";
1513  ring r = 0,(x,y,z),(dp);
1514  list L;
1515  module M;
1516  M = [1,z,y],[-1,x]; // zy - x
1517  L[1] = M;
1518  M = [1,z,x],[-1,y]; // zx - y
1519  L[2] = M;
1520  M = [1,y,x],[-1,z]; // yx - z
1521  L[3] = M;
1522  lst2str(L);
1523  list U = freegb(L,4);
1524  lst2str(U);
1525  // strange answer CHECK
1526}
1527
1528//static proc ex_densep_2()
1529{
1530  option(prot);
1531  LIB "freegb.lib";
1532  ring r = (0,a,b,c),(x,y),(Dp); // deglex
1533  module M = [1,x,x], [a,x,y], [b,y,x], [c,y,y];
1534  lst2str(M);
1535  list U = freegb(M,5);
1536  lst2str(U);
1537  // a=b is important -> finite basis!!!
1538  module M = [1,x,x], [a,x,y], [a,y,x], [c,y,y];
1539  lst2str(M);
1540  list U = freegb(M,5);
1541  lst2str(U);
1542}
1543
1544// END COMMENTED EXAMPLES
1545
1546*/
1547
1548// 1. form a new ring
1549// 2. produce shifted generators
1550// 3. compute GB
1551// 4. skip shifted elts
1552// 5. go back to orig vars, produce strings/modules
1553// 6. return the result
1554
1555static proc freegbold(list LM, int d)
1556"USAGE:  freegbold(L, d);  L a list of modules, d an integer
1557RETURN:  ring
1558PURPOSE: compute the two-sided Groebner basis of an ideal, encoded by L in
1559the free associative algebra, up to degree d
1560EXAMPLE: example freegbold; shows examples
1561"
1562{
1563  // d = up to degree, will be shifted to d+1
1564  if (d<1) {"bad d"; return(0);}
1565
1566  int ppl = printlevel-voice+2;
1567  string err = "";
1568
1569  int i,j,s;
1570  def save = basering;
1571  // determine max no of places in the input
1572  int slm = size(LM); // numbers of polys in the ideal
1573  int sm;
1574  intvec iv;
1575  module M;
1576  for (i=1; i<=slm; i++)
1577  {
1578    // modules, e.g. free polynomials
1579    M  = LM[i];
1580    sm = ncols(M);
1581    for (j=1; j<=sm; j++)
1582    {
1583      //vectors, e.g. free monomials
1584      iv = iv, size(M[j])-1; // 1 place is reserved by the coeff
1585    }
1586  }
1587  int D = Max(iv); // max size of input words
1588  if (d<D) {"bad d"; return(LM);}
1589  D = D + d-1;
1590  //  D = d;
1591  list LR  = ringlist(save);
1592  list L, tmp;
1593  L[1] = LR[1]; // ground field
1594  L[4] = LR[4]; // quotient ideal
1595  tmp  = LR[2]; // varnames
1596  s = size(LR[2]);
1597  for (i=1; i<=D; i++)
1598  {
1599    for (j=1; j<=s; j++)
1600    {
1601      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
1602    }
1603  }
1604  for (i=1; i<=s; i++)
1605  {
1606    tmp[i] = string(tmp[i])+"("+string(1)+")";
1607  }
1608  L[2] = tmp;
1609  list OrigNames = LR[2];
1610  // ordering: d blocks of the ord on r
1611  // try to get whether the ord on r is blockord itself
1612  // TODO: make L(2) ordering! exponent is maximally 2
1613  s = size(LR[3]);
1614  if (s==2)
1615  {
1616    // not a blockord, 1 block + module ord
1617    tmp = LR[3][s]; // module ord
1618    for (i=1; i<=D; i++)
1619    {
1620      LR[3][s-1+i] = LR[3][1];
1621    }
1622    LR[3][s+D] = tmp;
1623  }
1624  if (s>2)
1625  {
1626    // there are s-1 blocks
1627    int nb = s-1;
1628    tmp = LR[3][s]; // module ord
1629    for (i=1; i<=D; i++)
1630    {
1631      for (j=1; j<=nb; j++)
1632      {
1633        LR[3][i*nb+j] = LR[3][j];
1634      }
1635    }
1636    //    size(LR[3]);
1637    LR[3][nb*(D+1)+1] = tmp;
1638  }
1639  L[3] = LR[3];
1640  def @R = ring(L);
1641  setring @R;
1642  ideal I;
1643  poly @p;
1644  s = size(OrigNames);
1645  //  "s:";s;
1646  // convert LM to canonical vectors (no powers)
1647  setring save;
1648  kill M; // M was defined earlier
1649  module M;
1650  slm = size(LM); // numbers of polys in the ideal
1651  int sv,k,l;
1652  vector v;
1653  //  poly p;
1654  string sp;
1655  setring @R;
1656  poly @@p=0;
1657  setring save;
1658  for (l=1; l<=slm; l++)
1659  {
1660    // modules, e.g. free polynomials
1661    M  = LM[l];
1662    sm = ncols(M); // in intvec iv the sizes are stored
1663    for (i=0; i<=d-iv[l]; i++)
1664    {
1665      // modules, e.g. free polynomials
1666      for (j=1; j<=sm; j++)
1667      {
1668        //vectors, e.g. free monomials
1669        v  = M[j];
1670        sv = size(v);
1671        //        "sv:";sv;
1672        sp = "@@p = @@p + ";
1673        for (k=2; k<=sv; k++)
1674        {
1675          sp = sp + string(v[k])+"("+string(k-1+i)+")*";
1676        }
1677        sp = sp + string(v[1])+";"; // coef;
1678        setring @R;
1679        execute(sp);
1680        setring save;
1681      }
1682      setring @R;
1683      //      "@@p:"; @@p;
1684      I = I,@@p;
1685      @@p = 0;
1686      setring save;
1687    }
1688  }
1689  kill sp;
1690  // 3. compute GB
1691  setring @R;
1692  dbprint(ppl,"computing GB");
1693  //  ideal J = groebner(I);
1694  ideal J = slimgb(I);
1695  dbprint(ppl,J);
1696  // 4. skip shifted elts
1697  ideal K = select1(J,1..s); // s = size(OrigNames)
1698  dbprint(ppl,K);
1699  dbprint(ppl, "done with GB");
1700  // K contains vars x(1),...z(1) = images of originals
1701  // 5. go back to orig vars, produce strings/modules
1702  if (K[1] == 0)
1703  {
1704    "no reasonable output, GB gives 0";
1705    return(0);
1706  }
1707  int sk = size(K);
1708  int sp, sx, a, b;
1709  intvec x;
1710  poly p,q;
1711  poly pn;
1712  // vars in 'save'
1713  setring save;
1714  module N;
1715  list LN;
1716  vector V;
1717  poly pn;
1718  // test and skip exponents >=2
1719  setring @R;
1720  for(i=1; i<=sk; i++)
1721  {
1722    p  = K[i];
1723    while (p!=0)
1724    {
1725      q  = lead(p);
1726      //      "processing q:";q;
1727      x  = leadexp(q);
1728      sx = size(x);
1729      for(k=1; k<=sx; k++)
1730      {
1731        if ( x[k] >= 2 )
1732        {
1733          err = "skip: the value x[k] is " + string(x[k]);
1734          dbprint(ppl,err);
1735          //            return(0);
1736          K[i] = 0;
1737          p    = 0;
1738          q    = 0;
1739          break;
1740        }
1741      }
1742      p  = p - q;
1743    }
1744  }
1745  K  = simplify(K,2);
1746  sk = size(K);
1747  for(i=1; i<=sk; i++)
1748  {
1749    //    setring save;
1750    //    V  = 0;
1751    setring @R;
1752    p  = K[i];
1753    while (p!=0)
1754    {
1755      q  = lead(p);
1756      err =  "processing q:" + string(q);
1757      dbprint(ppl,err);
1758      x  = leadexp(q);
1759      sx = size(x);
1760      pn = leadcoef(q);
1761      setring save;
1762      pn = imap(@R,pn);
1763      V  = V + leadcoef(pn)*gen(1);
1764      for(k=1; k<=sx; k++)
1765      {
1766        if (x[k] ==1)
1767        {
1768          a = k div s; // block number=a+1, a!=0
1769          b = k % s; // remainder
1770          //          printf("a: %s, b: %s",a,b);
1771          if (b == 0)
1772          {
1773            // that is it's the last var in the block
1774            b = s;
1775            a = a-1;
1776          }
1777          V = V + var(b)*gen(a+2);
1778        }
1779//         else
1780//         {
1781//           printf("error: the value x[k] is %s", x[k]);
1782//           return(0);
1783//         }
1784      }
1785      err = "V: " + string(V);
1786      dbprint(ppl,err);
1787      //      printf("V: %s", string(V));
1788      N = N,V;
1789      V  = 0;
1790      setring @R;
1791      p  = p - q;
1792      pn = 0;
1793    }
1794    setring save;
1795    LN[i] = simplify(N,2);
1796    N     = 0;
1797  }
1798  setring save;
1799  return(LN);
1800}
1801example
1802{
1803  "EXAMPLE:"; echo = 2;
1804  ring r = 0,(x,y,z),(dp(1),dp(2));
1805  module M = [-1,x,y],[-7,y,y],[3,x,x];
1806  module N = [1,x,y,x],[-1,y,x,y];
1807  list L; L[1] = M; L[2] = N;
1808  lst2str(L);
1809  def U = freegbold(L,5);
1810  lst2str(U);
1811}
1812
1813/* begin older procs and tests
1814
1815static proc sgb(ideal I, int d)
1816{
1817  // new code
1818  // map x_i to x_i(1) via map()
1819  //LIB "freegb.lib";
1820  def save = basering;
1821  //int d =7;// degree
1822  int nv = nvars(save);
1823  def R = makeLetterplaceRing(d);
1824  setring R;
1825  int i;
1826  ideal Imap;
1827  for (i=1; i<=nv; i++)
1828  {
1829    Imap[i] = var(i);
1830  }
1831  //ideal I = x(1)*y(2), y(1)*x(2)+z(1)*z(2);
1832  ideal I = x(1)*x(2),x(1)*y(2) + z(1)*x(2);
1833  option(prot);
1834  //option(teach);
1835  ideal J = system("freegb",I,d,nv);
1836}
1837
1838static proc checkCeq()
1839{
1840  ring r = 0,(x,y),Dp;
1841  def A = makeLetterplaceRing(4);
1842  setring A;
1843  A;
1844  // I = x2-xy
1845  ideal I = x(1)*x(2) - x(1)*y(2), x(2)*x(3) - x(2)*y(3), x(3)*x(4) - x(3)*y(4);
1846  ideal C = x(2)-x(1),x(3)-x(2),x(4)-x(3),y(2)-y(1),y(3)-y(2),y(4)-y(3);
1847  ideal K = I,C;
1848  groebner(K);
1849}
1850
1851static proc exHom1()
1852{
1853  // we start with
1854  // z*y - x, z*x - y, y*x - z
1855  LIB "freegb.lib";
1856  LIB "elim.lib";
1857  ring r = 0,(x,y,z,h),dp;
1858  list L;
1859  module M;
1860  M = [1,z,y],[-1,x,h]; // zy - xh
1861  L[1] = M;
1862  M = [1,z,x],[-1,y,h]; // zx - yh
1863  L[2] = M;
1864  M = [1,y,x],[-1,z,h]; // yx - zh
1865  L[3] = M;
1866  lst2str(L);
1867  def U = crs(L,4);
1868  setring U;
1869  I = I,
1870    y(2)*h(3)+z(2)*x(3),     y(3)*h(4)+z(3)*x(4),
1871    y(2)*x(3)-z(2)*h(3),     y(3)*x(4)-z(3)*h(4);
1872  I = simplify(I,2);
1873  ring r2 = 0,(x(0..4),y(0..4),z(0..4),h(0..4)),dp;
1874  ideal J = imap(U,I);
1875  //  ideal K = homog(J,h);
1876  option(redSB);
1877  option(redTail);
1878  ideal L = groebner(J); //(K);
1879  ideal LL = sat(L,ideal(h))[1];
1880  ideal M = subst(LL,h,1);
1881  M = simplify(M,2);
1882  setring U;
1883  ideal M = imap(r2,M);
1884  lst2str(U);
1885}
1886
1887static proc test1()
1888{
1889  LIB "freegb.lib";
1890  ring r = 0,(x,y),Dp;
1891  int d = 10; // degree
1892  def R = makeLetterplaceRing(d);
1893  setring R;
1894  ideal I = x(1)*x(2) - y(1)*y(2);
1895  option(prot);
1896  option(teach);
1897  ideal J = system("freegb",I,d,2);
1898  J;
1899}
1900
1901static proc test2()
1902{
1903  LIB "freegb.lib";
1904  ring r = 0,(x,y),Dp;
1905  int d = 10; // degree
1906  def R = makeLetterplaceRing(d);
1907  setring R;
1908  ideal I = x(1)*x(2) - x(1)*y(2);
1909  option(prot);
1910  option(teach);
1911  ideal J = system("freegb",I,d,2);
1912  J;
1913}
1914
1915static proc test3()
1916{
1917  LIB "freegb.lib";
1918  ring r = 0,(x,y,z),dp;
1919  int d =5; // degree
1920  def R = makeLetterplaceRing(d);
1921  setring R;
1922  ideal I = x(1)*y(2), y(1)*x(2)+z(1)*z(2);
1923  option(prot);
1924  option(teach);
1925  ideal J = system("freegb",I,d,3);
1926}
1927
1928static proc schur2-3()
1929{
1930  // nonhomog:
1931  //  h^4-10*h^2+9,f*e-e*f+h, h*2-e*h-2*e,h*f-f*h+2*f
1932  // homogenized with t
1933  //  h^4-10*h^2*t^2+9*t^4,f*e-e*f+h*t, h*2-e*h-2*e*t,h*f-f*h+2*f*t,
1934  // t*h - h*t, t*f - f*t, t*e - e*t
1935}
1936
1937end older procs and tests */
1938
1939proc ademRelations(int i, int j)
1940"USAGE:  ademRelations(i,j); i,j int
1941RETURN:  ring (and exports ideal)
1942ASSUME: there are at least i+j variables in the basering
1943PURPOSE: compute the ideal of Adem relations for i<2j in characteristic 0
1944@*  the ideal is exported under the name AdemRel in the output ring
1945EXAMPLE: example ademRelations; shows examples
1946"
1947{
1948  // produces Adem relations for i<2j in char 0
1949  // assume: 0<i<2j
1950  // requires presence of vars up to i+j
1951  if ( (i<0) || (i >= 2*j) )
1952  {
1953    ERROR("arguments out of range"); return(0);
1954  }
1955  ring @r = 0,(s(i+j..0)),lp;
1956  poly p,q;
1957  number n;
1958  int ii = i div 2; int k;
1959  // k=0 => s(0)=1
1960  n = binomial(j-1,i);
1961  q = n*s(i+j)*s(0);
1962  //  printf("k=0, term=%s",q);
1963  p = p + q;
1964  for (k=1; k<= ii; k++)
1965  {
1966    n = binomial(j-k-1,i-2*k);
1967    q = n*s(i+j-k)*s(k);;
1968    //    printf("k=%s, term=%s",k,q);
1969    p = p + q;
1970  }
1971  poly AdemRel = p;
1972  export AdemRel;
1973  return(@r);
1974}
1975example
1976{
1977  "EXAMPLE:"; echo = 2;
1978  def A = ademRelations(2,5);
1979  setring A;
1980  AdemRel;
1981}
1982
1983/*
19841,1: 0
19851,2: s(3)*s(0) == s(3) -> def for s(3):=s(1)s(2)
19862,1: adm
19872,2: s(3)*s(1) == s(1)s(2)s(1)
19881,3: 0 ( since 2*s(4)*s(0) = 0 mod 2)
19893,1: adm
19902,3: s(5)*s(0)+s(4)*s(1) == s(5)+s(4)*s(1)
19913,2: 0
19923,3: s(5)*s(1)
19931,4: 3*s(5)*s(0) == s(5)  -> def for s(5):=s(1)*s(4)
19944,1: adm
19952,4: 3*s(6)*s(0)+s(5)*s(1) == s(6) + s(5)*s(1) == s(6) + s(1)*s(4)*s(1)
19964,2: adm
19974,3: s(5)*s(2)
19983,4: s(7)*s(0)+2*s(6)*s(1) == s(7) -> def for s(7):=s(3)*s(4)
19994,4: s(7)*s(1)+s(6)*s(2)
2000*/
2001
2002/* s1,s2:
2003s1*s1 =0, s2*s2 = s1*s2*s1
2004*/
2005
2006/*
2007try char 0:
2008s1,s2:
2009s1*s1 =0, s2*s2 = s1*s2*s1, s(1)*s(3)== s(1)*s(1)*s(3) == 0 = 2*s(4) ->def for s(4)
2010hence 2==0! only in char 2
2011 */
2012
2013  // Adem rels modulo 2 are interesting
2014
2015static proc stringpoly2lplace(string s)
2016{
2017  // decomposes sentence into terms
2018  s = replace(s,newline,""); // get rid of newlines
2019  s = replace(s," ",""); // get rid of empties
2020  //arith symbols: +,-
2021  // decompose into words with coeffs
2022  list LS;
2023  int i,j,ie,je,k,cnt;
2024  // s[1]="-" situation
2025  if (s[1]=="-")
2026  {
2027    LS = stringpoly2lplace(string(s[2..size(s)]));
2028    LS[1] = string("-"+string(LS[1]));
2029    return(LS);
2030  }
2031  i = find(s,"-",2);
2032  // i==1 might happen if the 1st symbol coeff is negative
2033  j = find(s,"+");
2034  list LL;
2035  if (i==j)
2036  {
2037    "return a monomial";
2038    // that is both are 0 -> s is a monomial
2039    LS[1] = s;
2040    return(LS);
2041  }
2042  if (i==0)
2043  {
2044    "i==0 situation";
2045    // no minuses at all => pluses only
2046    cnt++;
2047    LS[cnt] = string(s[1..j-1]);
2048    s = s[j+1..size(s)];
2049    while (s!= "")
2050    {
2051      j = find(s,"+");
2052      cnt++;
2053      if (j==0)
2054      {
2055        LS[cnt] = string(s);
2056        s = "";
2057      }
2058      else
2059      {
2060        LS[cnt] = string(s[1..j-1]);
2061        s = s[j+1..size(s)];
2062      }
2063    }
2064    return(LS);
2065  }
2066  if (j==0)
2067  {
2068    "j==0 situation";
2069    // no pluses at all except the lead coef => the rest are minuses only
2070    cnt++;
2071    LS[cnt] = string(s[1..i-1]);
2072    s = s[i..size(s)];
2073    while (s!= "")
2074    {
2075      i = find(s,"-",2);
2076      cnt++;
2077      if (i==0)
2078      {
2079        LS[cnt] = string(s);
2080        s = "";
2081      }
2082      else
2083      {
2084        LS[cnt] = string(s[1..i-1]);
2085        s = s[i..size(s)];
2086      }
2087    }
2088    return(LS);
2089  }
2090  // now i, j are nonzero
2091  if (i>j)
2092  {
2093    "i>j situation";
2094    // + comes first, at place j
2095    cnt++;
2096    //    "cnt:"; cnt; "j:"; j;
2097    LS[cnt] = string(s[1..j-1]);
2098    s = s[j+1..size(s)];
2099    LL = stringpoly2lplace(s);
2100    LS = LS + LL;
2101    kill LL;
2102    return(LS);
2103  }
2104  else
2105  {
2106    "j>i situation";
2107    // - might come first, at place i
2108    if (i>1)
2109    {
2110      cnt++;
2111      LS[cnt] = string(s[1..i-1]);
2112      s = s[i..size(s)];
2113    }
2114    else
2115    {
2116      // i==1->  minus at leadcoef
2117      ie = find(s,"-",i+1);
2118      je = find(s,"+",i+1);
2119      if (je == ie)
2120      {
2121         "ie=je situation";
2122        //monomial
2123        cnt++;
2124        LS[cnt] = s;
2125        return(LS);
2126      }
2127      if (je < ie)
2128      {
2129         "je<ie situation";
2130        // + comes first
2131        cnt++;
2132        LS[cnt] = s[1..je-1];
2133        s = s[je+1..size(s)];
2134      }
2135      else
2136      {
2137        // ie < je
2138         "ie<je situation";
2139        cnt++;
2140        LS[cnt] = s[1..ie-1];
2141        s = s[ie..size(s)];
2142      }
2143    }
2144    "going into recursion with "+s;
2145    LL = stringpoly2lplace(s);
2146    LS = LS + LL;
2147    return(LS);
2148  }
2149}
2150example
2151{
2152  "EXAMPLE:"; echo = 2;
2153  string s = "x*y+y*z+z*t"; // + only
2154  stringpoly2lplace(s);
2155  string s2 = "x*y - y*z-z*t*w*w"; // +1, - only
2156  stringpoly2lplace(s2);
2157  string s3 = "-x*y + y - 2*x +7*w*w*w";
2158  stringpoly2lplace(s3);
2159}
2160
2161static proc addplaces(list L)
2162{
2163  // adds places to the list of strings
2164  // according to their order in the list
2165  int s = size(L);
2166  int i;
2167  for (i=1; i<=s; i++)
2168  {
2169    if (typeof(L[i]) == "string")
2170    {
2171      L[i] = L[i] + "(" + string(i) + ")";
2172    }
2173    else
2174    {
2175      ERROR("entry of type string expected");
2176      return(0);
2177    }
2178  }
2179  return(L);
2180}
2181example
2182{
2183  "EXAMPLE:"; echo = 2;
2184  string a = "f1";   string b = "f2";
2185  list L = a,b,a;
2186  addplaces(L);
2187}
2188
2189static proc sent2lplace(string s)
2190{
2191  // SENTence of words TO LetterPLACE
2192  list L =   stringpoly2lplace(s);
2193  int i; int ss = size(L);
2194  for(i=1; i<=ss; i++)
2195  {
2196    L[i] = str2lplace(L[i]);
2197  }
2198  return(L);
2199}
2200example
2201{
2202  "EXAMPLE:"; echo = 2;
2203  ring r = 0,(f2,f1),dp;
2204  string s = "f2*f1*f1 - 2*f1*f2*f1+ f1*f1*f2";
2205  sent2lplace(s);
2206}
2207
2208static proc testnumber(string s)
2209{
2210  string t;
2211  if (s[1]=="-")
2212  {
2213    // two situations: either there's a negative number
2214    t = s[2..size(s)];
2215    if (testnumber(t))
2216    {
2217      //a negative number
2218    }
2219    else
2220    {
2221      // a variable times -1
2222    }
2223    // or just a "-" for -1
2224  }
2225  t = "ring @r=(";
2226  t = t + charstr(basering)+"),";
2227  t = t + string(var(1))+",dp;";
2228  //  write(":w tstnum.tst",t);
2229  t = t+ "number @@Nn = " + s + ";"+"$";
2230  write(":w tstnum.tst",t);
2231  string runsing = system("Singular");
2232  int k;
2233  t = runsing+ " -teq <tstnum.tst >tstnum.out";
2234  k = system("sh",t);
2235  if (k!=0)
2236  {
2237    ERROR("Problems running Singular");
2238  }
2239  int i = system("sh", "grep error tstnum.out > /dev/NULL");
2240  if (i!=0)
2241  {
2242    // no error: s is a number
2243    i = 1;
2244  }
2245  k = system("sh","rm tstnum.tst tstnum.out > /dev/NULL");
2246  return(i);
2247}
2248example
2249{
2250  "EXAMPLE:"; echo = 2;
2251  ring r = (0,a),x,dp;
2252  string s = "a^2+7*a-2";
2253  testnumber(s);
2254  s = "b+a";
2255  testnumber(s);
2256}
2257
2258static proc str2lplace(string s)
2259{
2260  // converts a word (monomial) with coeff into letter-place
2261  // string: coef*var1^exp1*var2^exp2*...varN^expN
2262  s = strpower2rep(s); // expand powers
2263  if (size(s)==0) { return(0); }
2264  int i,j,k,insC;
2265  string a,b,c,d,t;
2266  // 1. get coeff
2267  i = find(s,"*");
2268  if (i==0) { return(s); }
2269  list VN;
2270  c = s[1..i-1]; // incl. the case like (-a^2+1)
2271  int tn = testnumber(c);
2272  if (tn == 0)
2273  {
2274    // failed test
2275    if (c[1]=="-")
2276    {
2277      // two situations: either there's a negative number
2278      t = c[2..size(c)];
2279      if (testnumber(t))
2280      {
2281         //a negative number
2282        // nop here
2283      }
2284      else
2285      {
2286         // a variable times -1
2287          c = "-1";
2288          j++; VN[j] = t; //string(c[2..size(c)]);
2289          insC = 1;
2290      }
2291    }
2292    else
2293    {
2294      // just a variable with coeff 1
2295          j++; VN[j] = string(c);
2296          c = "1";
2297          insC = 1;
2298    }
2299  }
2300 // get vars
2301  t = s;
2302  //  t = s[i+1..size(s)];
2303  k = size(t); //j = 0;
2304  while (k>0)
2305  {
2306    t = t[i+1..size(t)]; //next part
2307    i = find(t,"*"); // next *
2308    if (i==0)
2309    {
2310      // last monomial
2311      j++;
2312      VN[j] = t;
2313      k = size(t);
2314      break;
2315    }
2316    b = t[1..i-1];
2317    //    print(b);
2318    j++;
2319    VN[j] = b;
2320    k = size(t);
2321  }
2322  VN = addplaces(VN);
2323  VN[size(VN)+1] = string(c);
2324  return(VN);
2325}
2326example
2327{
2328  "EXAMPLE:"; echo = 2;
2329  ring r = (0,a),(f2,f1),dp;
2330  str2lplace("-2*f2^2*f1^2*f2");
2331  str2lplace("-f1*f2");
2332  str2lplace("(-a^2+7a)*f1*f2");
2333}
2334
2335static proc strpower2rep(string s)
2336{
2337  // makes x*x*x*x out of x^4 ., rep statys for repetitions
2338  // looks for "-" problem
2339  // exception: "-" as coeff
2340  string ex,t;
2341  int i,j,k;
2342
2343  i = find(s,"^"); // first ^
2344  if (i==0) { return(s); } // no ^ signs
2345
2346  if (s[1] == "-")
2347  {
2348    // either -coef or -1
2349    // got the coeff:
2350    i = find(s,"*");
2351    if (i==0)
2352    {
2353      // no *'s   => coef == -1 or s == -23
2354      i = size(s)+1;
2355    }
2356    t = string(s[2..i-1]); // without "-"
2357    if ( testnumber(t) )
2358    {
2359      // a good number
2360      t = strpower2rep(string(s[2..size(s)]));
2361      t = "-"+t;
2362      return(t);
2363    }
2364    else
2365    {
2366      // a variable
2367      t = strpower2rep(string(s[2..size(s)]));
2368      t = "-1*"+ t;
2369      return(t);
2370    }
2371  }
2372  // the case when leadcoef is a number in ()
2373  if (s[1] == "(")
2374  {
2375    i = find(s,")",2);    // must be nonzero
2376    t = s[2..i-1];
2377    if ( testnumber(t) )
2378    {
2379      // a good number
2380    }
2381    else {"strpower2rep: bad number as coef";}
2382    ex = string(s[i+2..size(s)]); // 2 because of *
2383    ex =  strpower2rep(ex);
2384    t = "("+t+")*"+ex;
2385    return(t);
2386  }
2387
2388  i = find(s,"^"); // first ^
2389  j = find(s,"*",i+1); // next * == end of ^
2390  if (j==0)
2391  {
2392    ex = s[i+1..size(s)];
2393  }
2394  else
2395  {
2396    ex = s[i+1..j-1];
2397  }
2398  execute("int @exp = " + ex + ";"); //@exp = exponent
2399  // got varname
2400  for (k=i-1; k>0; k--)
2401  {
2402    if (s[k] == "*") break;
2403  }
2404  string varn = s[k+1..i-1];
2405  //  "varn:";  varn;
2406  string pref;
2407  if (k>0)
2408  {
2409    pref = s[1..k]; // with * on the k-th place
2410  }
2411  //  "pref:";  pref;
2412  string suf;
2413  if ( (j>0) && (j+1 <= size(s)) )
2414  {
2415    suf = s[j+1..size(s)]; // without * on the 1st place
2416  }
2417  //  "suf:"; suf;
2418  string toins;
2419  for (k=1; k<=@exp; k++)
2420  {
2421    toins = toins + varn+"*";
2422  }
2423  //  "toins: ";  toins;
2424  if (size(suf) == 0)
2425  {
2426    toins = toins[1..size(toins)-1]; // get rid of trailing *
2427  }
2428  else
2429  {
2430    suf = strpower2rep(suf);
2431  }
2432  ex = pref + toins + suf;
2433  return(ex);
2434  //  return(strpower2rep(ex));
2435}
2436example
2437{
2438  "EXAMPLE:"; echo = 2;
2439  ring r = (0,a),(x,y,z,t),dp;
2440  strpower2rep("-x^4");
2441  strpower2rep("-2*x^4*y^3*z*t^2");
2442  strpower2rep("-a^2*x^4");
2443}
2444
2445proc lieBracket(poly a, poly b, list #)
2446"USAGE:  lieBracket(a,b[,N]); a,b letterplace polynomials, N an optional integer
2447RETURN:  poly
2448ASSUME: basering has a letterplace ring structure
2449PURPOSE:compute the Lie bracket [a,b] = ab - ba between letterplace polynomials
2450NOTE: if N>1 is specified, then the left normed bracket [a,[...[a,b]]]] is
2451@*    computed.
2452EXAMPLE: example lieBracket; shows examples
2453"
2454{
2455  if (lpAssumeViolation())
2456  {
2457    //    ERROR("Either 'uptodeg' or 'lV' global variables are not set!");
2458    ERROR("Incomplete Letterplace structure on the basering!");
2459  }
2460  // alias ppLiebr;
2461  //if int N is given compute [a,[...[a,b]]]] left normed bracket
2462  poly q;
2463  int N=1;
2464  if (size(#)>0)
2465  {
2466    if (typeof(#[1])=="int")
2467    {
2468      N = int(#[1]);
2469    }
2470  }
2471  if (N<=0) { return(q); }
2472  while (b!=0)
2473  {
2474    q = q + pmLiebr(a,lead(b));
2475    b = b - lead(b);
2476  }
2477  int i;
2478  if (N >1)
2479  {
2480    for(i=1; i<=N; i++)
2481    {
2482      q = lieBracket(a,q);
2483    }
2484  }
2485  return(q);
2486}
2487example
2488{
2489  "EXAMPLE:"; echo = 2;
2490  ring r = 0,(x(1),y(1),x(2),y(2),x(3),y(3),x(4),y(4)),dp;
2491  def R = setLetterplaceAttributes(r,4,2); // supply R with letterplace structure
2492  setring R;
2493  poly a = x(1)*y(2); poly b = y(1);
2494  lieBracket(a,b);
2495  lieBracket(x(1),y(1),2);
2496}
2497
2498static proc pmLiebr(poly a, poly b)
2499{
2500  //  a poly, b mono
2501  poly s;
2502  while (a!=0)
2503  {
2504    s = s + mmLiebr(lead(a),lead(b));
2505    a = a - lead(a);
2506  }
2507  return(s);
2508}
2509
2510proc shiftPoly(poly a, int i)
2511"USAGE:  shiftPoly(p,i); p letterplace poly, i int
2512RETURN: poly
2513ASSUME: basering has letterplace ring structure
2514PURPOSE: compute the i-th shift of letterplace polynomial p
2515EXAMPLE: example shiftPoly; shows examples
2516"
2517{
2518  // shifts a monomial a by i
2519  // calls pLPshift(p,sh,uptodeg,lVblock);
2520  if (lpAssumeViolation())
2521  {
2522    ERROR("Incomplete Letterplace structure on the basering!");
2523  }
2524  int uptodeg = attrib(basering,"uptodeg");
2525  int lV = attrib(basering,"lV");
2526  if (deg(a) + i > uptodeg)
2527  {
2528    ERROR("degree bound violated by the shift!");
2529  }
2530  return(system("stest",a,i,uptodeg,lV));
2531}
2532example
2533{
2534  "EXAMPLE:"; echo = 2;
2535  ring r = 0,(x,y,z),dp;
2536  int uptodeg = 5; int lV = 3;
2537  def R = makeLetterplaceRing(uptodeg);
2538  setring R;
2539  poly f = x(1)*z(2)*y(3) - 2*z(1)*y(2) + 3*x(1);
2540  shiftPoly(f,1);
2541  shiftPoly(f,2);
2542}
2543
2544
2545static proc mmLiebr(poly a, poly b)
2546{
2547  // a,b, monomials
2548  a = lead(a);
2549  b = lead(b);
2550  int sa = deg(a);
2551  int sb = deg(b);
2552  poly v = a*shiftPoly(b,sa) - b*shiftPoly(a,sb);
2553  return(v);
2554}
2555
2556static proc test_shift()
2557{
2558  LIB "freegb.lib";
2559  ring r = 0,(a,b),dp;
2560  int d =5;
2561  def R = makeLetterplaceRing(d);
2562  setring R;
2563  int uptodeg = d;
2564  int lV = 2;
2565  def R = setLetterplaceAttributes(r,uptodeg,2); // supply R with letterplace structure
2566  setring R;
2567  poly p = mmLiebr(a(1),b(1));
2568  poly p = lieBracket(a(1),b(1));
2569}
2570
2571proc serreRelations(intmat A, int zu)
2572"USAGE:  serreRelations(A,z); A an intmat, z an int
2573RETURN:  ideal
2574ASSUME: basering has a letterplace ring structure and
2575@*          A is a generalized Cartan matrix with integer entries
2576PURPOSE: compute the ideal of Serre's relations associated to A
2577EXAMPLE: example serreRelations; shows examples
2578"
2579{
2580  // zu = 1 -> with commutators [f_i,f_j]; zu == 0 without them
2581  // suppose that A is cartan matrix
2582  // then Serre's relations are
2583  // (ad f_j)^{1-A_{ij}} ( f_i)
2584  int ppl = printlevel-voice+2;
2585  int n = ncols(A); // hence n variables
2586  int i,j,k,el;
2587  poly p,q;
2588  ideal I;
2589  for (i=1; i<=n; i++)
2590  {
2591    for (j=1; j<=n; j++)
2592    {
2593      el = 1 - A[i,j];
2594      //     printf("i:%s, j: %s, l: %s",i,j,l);
2595      dbprint(ppl,"i, j, l: ",i,j,el);
2596      //      if ((i!=j) && (l >0))
2597      //      if ( (i!=j) &&  ( ((zu ==0) &&  (l >=2)) || ((zu ==1) &&  (l >=1)) ) )
2598      if ((i!=j) && (el >0))
2599      {
2600        q = lieBracket(var(j),var(i));
2601        dbprint(ppl,"first bracket: ",q);
2602        //        if (l >=2)
2603        //        {
2604          for (k=1; k<=el-1; k++)
2605          {
2606            q = lieBracket(var(j),q);
2607            dbprint(ppl,"further bracket:",q);
2608          }
2609          //        }
2610      }
2611      if (q!=0) { I = I,q; q=0;}
2612    }
2613  }
2614  I = simplify(I,2);
2615  return(I);
2616}
2617example
2618{
2619  "EXAMPLE:"; echo = 2;
2620  intmat A[3][3] =
2621    2, -1, 0,
2622    -1, 2, -3,
2623    0, -1, 2; // G^1_2 Cartan matrix
2624  ring r = 0,(f1,f2,f3),dp;
2625  int uptodeg = 5;
2626  def R = makeLetterplaceRing(uptodeg);
2627  setring R;
2628  ideal I = serreRelations(A,1); I = simplify(I,1+2+8);
2629  I;
2630}
2631
2632/* setup for older example:
2633  intmat A[2][2] = 2, -1, -1, 2; // sl_3 == A_2
2634  ring r = 0,(f1,f2),dp;
2635  int uptodeg = 5; int lV = 2;
2636*/
2637
2638proc fullSerreRelations(intmat A, ideal rNegative, ideal rCartan, ideal rPositive, int uptodeg)
2639"USAGE:  fullSerreRelations(A,N,C,P,d); A an intmat, N,C,P ideals, d an int
2640RETURN:  ring (and ideal)
2641PURPOSE: compute the inhomogeneous Serre's relations associated to A in given
2642@*       variable names
2643ASSUME: three ideals in the input are of the same sizes and contain merely
2644@* variables which are interpreted as follows: N resp. P stand for negative
2645@* resp. positive roots, C stand for Cartan elements. d is the degree bound for
2646@* letterplace ring, which will be returned.
2647@* The matrix A is a generalized Cartan matrix with integer entries
2648@* The result is the ideal called 'fsRel' in the returned ring.
2649EXAMPLE: example fullSerreRelations; shows examples
2650"
2651{
2652  /* SerreRels on rNeg and rPos plus Cartans etc. */
2653  int ppl = printlevel -voice+2;
2654  /* ideals must be written in variables: assume each term is of degree 1 */
2655  int i,j,k;
2656  int N = nvars(basering);
2657  def save = basering;
2658  int comFlag = 0;
2659  /* assume:  (size(rNegative) == size(rPositive)) */
2660  /* assume:  (size(rNegative) == size(rCartan)) i.e. nonsimple Cartans */
2661  if ( (size(rNegative) != size(rPositive)) || (size(rNegative) != size(rCartan)) )
2662  {
2663    ERROR("All input ideals must be of the same size");
2664  }
2665
2666//   if (size(rNegative) != size(rPositive))
2667//   {
2668//     ERROR("The 1st and the 3rd input ideals must be of the same size");
2669//   }
2670
2671  /* assume:  2*size(rNegative) + size(rCartan) >= nvars(basering) */
2672  i = 2*size(rNegative) + size(rCartan);
2673  if (i>N)
2674  {
2675    string s1="The total number of elements in input ideals";
2676    string s2="must not exceed the dimension of the ground ring";
2677    ERROR(s1+s2);
2678  }
2679  if (i < N)
2680  {
2681    comFlag = N-i; // so many elements will commute
2682    "Warning: some elements will be treated as mutually commuting";
2683  }
2684  /* extract varnames from input ideals */
2685  intvec iNeg = varIdeal2intvec(rNegative);
2686  intvec iCartan = varIdeal2intvec(rCartan);
2687  intvec iPos = varIdeal2intvec(rPositive);
2688  /* for each vector in rNeg and rPositive, go into the corr. ring and create SerreRels */
2689  /* rNegative: */
2690  list L = ringlist(save);
2691  def LPsave = makeLetterplaceRing2(uptodeg); setring save;
2692  list LNEG = L; list tmp;
2693  /* L[1] field as is; L[2] vars: a subset; L[3] ordering: dp, L[4] as is */
2694  for (i=1; i<=size(iNeg); i++)
2695  {
2696    tmp[i] = string(var(iNeg[i]));
2697  }
2698  LNEG[2] = tmp; LNEG[3] = list(list("dp",intvec(1:size(iNeg))), list("C",0));
2699  def RNEG = ring(LNEG); setring RNEG;
2700  def RRNEG = makeLetterplaceRing2(uptodeg);
2701  setring RRNEG;
2702  ideal I = serreRelations(A,1); I = simplify(I,1+2+8);
2703  setring LPsave;
2704  ideal srNeg = imap(RRNEG,I);
2705  dbprint(ppl,"0-1 ideal of negative relations is ready");
2706  dbprint(ppl-1,srNeg);
2707  setring save; kill L,tmp,RRNEG,RNEG, LNEG;
2708  /* rPositive: */
2709  list L = ringlist(save);
2710  list LPOS = L; list tmp;
2711  /* L[1] field as is; L[2] vars: a subset; L[3] ordering: dp, L[4] as is */
2712  for (i=1; i<=size(iPos); i++)
2713  {
2714    tmp[i] = string(var(iPos[i]));
2715  }
2716  LPOS[2] = tmp; LPOS[3] = list(list("dp",intvec(1:size(iPos))), list("C",0));
2717  def RPOS = ring(LPOS); setring RPOS;
2718  def RRPOS = makeLetterplaceRing2(uptodeg);
2719  setring RRPOS;
2720  ideal I = serreRelations(A,1); I = simplify(I,1+2+8);
2721  setring LPsave;
2722  ideal srPos = imap(RRPOS,I);
2723  dbprint(ppl,"0-2 ideal of positive relations is ready");
2724  dbprint(ppl-1,srPos);
2725  setring save; kill L,tmp,RRPOS,RPOS, LPOS;
2726  string sMap = "ideal Mmap =";
2727  for (i=1; i<=nvars(save); i++)
2728  {
2729    sMap = sMap + string(var(i)) +"(1),";
2730  }
2731  sMap[size(sMap)] = ";";
2732  /* cartans: h_j h_i = h_i h_j */
2733  setring LPsave;
2734  ideal ComCartan;
2735  for (i=1; i<size(iCartan); i++)
2736  {
2737    for (j=i+1; j<=size(iCartan); j++)
2738    {
2739      ComCartan =  ComCartan + lieBracket(var(iCartan[j]),var(iCartan[i]));
2740    }
2741  }
2742  ComCartan = simplify(ComCartan,1+2+8);
2743  execute(sMap); // defines an ideal Mmap
2744  map F = save, Mmap;
2745  dbprint(ppl,"1. commuting Cartans: ");
2746  dbprint(ppl-1,ComCartan);
2747  /* [e_i, f_j] =0 if i<>j */
2748  ideal ComPosNeg; // assume: #Neg=#Pos
2749  for (i=1; i<size(iPos); i++)
2750  {
2751    for (j=1; j<=size(iPos); j++)
2752    {
2753      if (j !=i)
2754      {
2755        ComPosNeg =  ComPosNeg + lieBracket(var(iPos[i]),var(iNeg[j]));
2756        ComPosNeg =  ComPosNeg + lieBracket(var(iPos[j]),var(iNeg[i]));
2757      }
2758    }
2759  }
2760  ComPosNeg = simplify(ComPosNeg,1+2+8);
2761  dbprint(ppl,"2. commuting Positive and Negative:");
2762  dbprint(ppl-1,ComPosNeg);
2763  /* [e_i, f_i] = h_i */
2764  poly tempo;
2765  for (i=1; i<=size(iCartan); i++)
2766  {
2767    tempo = lieBracket(var(iPos[i]),var(iNeg[i])) - var(iCartan[i]);
2768    ComPosNeg =  ComPosNeg + tempo;
2769  }
2770  //  ComPosNeg = simplify(ComPosNeg,1+2+8);
2771  dbprint(ppl,"3. added sl2 triples [e_i,f_i]=h_i");
2772  dbprint(ppl-1,ComPosNeg);
2773
2774  /* [h_i, e_j] = A_ij e_j */
2775  /* [h_i, f_j] = -A_ij f_j */
2776  ideal ActCartan; // assume: #Neg=#Pos
2777  for (i=1; i<=size(iCartan); i++)
2778  {
2779    for (j=1; j<=size(iCartan); j++)
2780    {
2781      tempo = lieBracket(var(iCartan[i]),var(iPos[j])) - A[i,j]*var(iPos[j]);
2782      ActCartan = ActCartan + tempo;
2783      tempo = lieBracket(var(iCartan[i]),var(iNeg[j])) + A[i,j]*var(iNeg[j]);
2784      ActCartan = ActCartan + tempo;
2785    }
2786  }
2787  ActCartan = simplify(ActCartan,1+2+8);
2788  dbprint(ppl,"4. actions of Cartan:");
2789  dbprint(ppl-1, ActCartan);
2790
2791  /* final part: prepare the output */
2792  setring LPsave;
2793  ideal fsRel = srNeg, srPos, ComPosNeg, ComCartan, ActCartan;
2794  export fsRel;
2795  setring save;
2796  return(LPsave);
2797}
2798example
2799{
2800  "EXAMPLE:"; echo = 2;
2801  intmat A[2][2] =
2802    2, -1,
2803    -1, 2; // A_2 = sl_3 Cartan matrix
2804  ring r = 0,(f1,f2,h1,h2,e1,e2),dp;
2805  ideal negroots = f1,f2; ideal cartans = h1,h2; ideal posroots = e1,e2;
2806  int uptodeg = 5;
2807  def RS = fullSerreRelations(A,negroots,cartans,posroots,uptodeg);
2808  setring RS; fsRel;
2809}
2810
2811static proc varIdeal2intvec(ideal I)
2812{
2813  // used in SerreRelations
2814  /* assume1:  input ideal is a list of variables of the ground ring */
2815  int i,j; intvec V;
2816  for (i=1; i<= size(I); i++)
2817  {
2818    j = univariate(I[i]);
2819    if (j<=0)
2820    {
2821      ERROR("input ideal must contain only variables");
2822    }
2823    V[i] = j;
2824  }
2825  dbprint(printlevel-voice+2,V);
2826  /* now we make a smaller list of non-repeating entries */
2827  ideal iW = simplify(ideal(V),2+4); // no zeros, no repetitions
2828  if (size(iW) < size(V))
2829  {
2830    /* extract intvec from iW */
2831    intvec inW;
2832    for(j=1; j<=size(iW); j++)
2833    {
2834      inW[j] = int(leadcoef(iW[j]));
2835    }
2836    return(inW);
2837  }
2838  return(V);
2839}
2840example
2841{
2842  "EXAMPLE:"; echo = 2;
2843  ring r = 0,(x,y,z),dp;
2844  ideal I = x,z;
2845  varIdeal2intvec(I);
2846  varIdeal2intvec(ideal(x2,y^3,x+1));
2847  varIdeal2intvec(ideal(x*y,y,x+1));
2848}
2849
2850proc lp2lstr(ideal K, def save)
2851"USAGE:  lp2lstr(K,s); K an ideal, s a ring name
2852RETURN:  nothing (exports object @LN into the ring named s)
2853ASSUME: basering has a letterplace ring structure
2854PURPOSE: converts letterplace ideal to list of modules
2855NOTE: useful as preprocessing to 'lst2str'
2856EXAMPLE: example lp2lstr; shows examples
2857"
2858{
2859  def @R = basering;
2860  string err;
2861  int s = nvars(save);
2862  int i,j,k;
2863    // K contains vars x(1),...z(1) = images of originals
2864  // 5. go back to orig vars, produce strings/modules
2865  int sk = size(K);
2866  int sp, sx, a, b;
2867  intvec x;
2868  poly p,q;
2869  poly pn;
2870  // vars in 'save'
2871  setring save;
2872  module N;
2873  list LN;
2874  vector V;
2875  poly pn;
2876  // test and skip exponents >=2
2877  setring @R;
2878  for(i=1; i<=sk; i++)
2879  {
2880    p  = K[i];
2881    while (p!=0)
2882    {
2883      q  = lead(p);
2884      //      "processing q:";q;
2885      x  = leadexp(q);
2886      sx = size(x);
2887      for(k=1; k<=sx; k++)
2888      {
2889        if ( x[k] >= 2 )
2890        {
2891          err = "skip: the value x[k] is " + string(x[k]);
2892          dbprint(ppl,err);
2893          //            return(0);
2894          K[i] = 0;
2895          p    = 0;
2896          q    = 0;
2897          break;
2898        }
2899      }
2900      p  = p - q;
2901    }
2902  }
2903  K  = simplify(K,2);
2904  sk = size(K);
2905  for(i=1; i<=sk; i++)
2906  {
2907    //    setring save;
2908    //    V  = 0;
2909    setring @R;
2910    p  = K[i];
2911    while (p!=0)
2912    {
2913      q  = lead(p);
2914      err =  "processing q:" + string(q);
2915      dbprint(ppl,err);
2916      x  = leadexp(q);
2917      sx = size(x);
2918      pn = leadcoef(q);
2919      setring save;
2920      pn = imap(@R,pn);
2921      V  = V + leadcoef(pn)*gen(1);
2922      for(k=1; k<=sx; k++)
2923      {
2924        if (x[k] ==1)
2925        {
2926          a = k div s; // block number=a+1, a!=0
2927          b = k % s; // remainder
2928          //          printf("a: %s, b: %s",a,b);
2929          if (b == 0)
2930          {
2931            // that is it's the last var in the block
2932            b = s;
2933            a = a-1;
2934          }
2935          V = V + var(b)*gen(a+2);
2936        }
2937      }
2938      err = "V: " + string(V);
2939      dbprint(ppl,err);
2940      //      printf("V: %s", string(V));
2941      N = N,V;
2942      V  = 0;
2943      setring @R;
2944      p  = p - q;
2945      pn = 0;
2946    }
2947    setring save;
2948    LN[i] = simplify(N,2);
2949    N     = 0;
2950  }
2951  setring save;
2952  list @LN = LN;
2953  export @LN;
2954  //  return(LN);
2955}
2956example
2957{
2958  "EXAMPLE:"; echo = 2;
2959  intmat A[2][2] = 2, -1, -1, 2; // sl_3 == A_2
2960  ring r = 0,(f1,f2),dp;
2961  def R = makeLetterplaceRing(3);
2962  setring R;
2963  ideal I = serreRelations(A,1);
2964  lp2lstr(I,r);
2965  setring r;
2966  lst2str(@LN,1);
2967}
2968
2969static proc strList2poly(list L)
2970{
2971  //  list L comes from sent2lplace (which takes a polynomial as input)
2972  // each entry of L is a sublist with the coef on the last place
2973  int s = size(L); int t;
2974  int i,j;
2975  list M;
2976  poly p,q;
2977  string Q;
2978  for(i=1; i<=s; i++)
2979  {
2980    M = L[i];
2981    t = size(M);
2982    //    q = M[t]; // a constant
2983    Q = string(M[t]);
2984    for(j=1; j<t; j++)
2985    {
2986      //      q = q*M[j];
2987      Q = Q+"*"+string(M[j]);
2988    }
2989    execute("q="+Q+";");
2990    //    q;
2991    p = p + q;
2992  }
2993  kill Q;
2994  return(p);
2995}
2996example
2997{
2998  "EXAMPLE:"; echo = 2;
2999  ring r =0,(x,y,z,t),Dp;
3000  def A = makeLetterplaceRing(4);
3001  setring A;
3002  string t = "-2*y*z*y*z + y*t*z*z - z*x*x*y  + 2*z*y*z*y";
3003  list L = sent2lplace(t);
3004  L;
3005  poly p = strList2poly(L);
3006  p;
3007}
3008
3009static proc file2lplace(string fname)
3010"USAGE:  file2lplace(fnm);  fnm a string
3011RETURN:  ideal
3012PURPOSE: convert the contents of the file fnm into ideal of polynomials in free algebra
3013EXAMPLE: example file2lplace; shows examples
3014"
3015{
3016  // format: from the usual string to letterplace
3017  string s = read(fname);
3018  // assume: file is a comma-sep list of polys
3019  // the vars are declared before
3020  // the file ends with ";"
3021  string t; int i;
3022  ideal I;
3023  list tst;
3024  while (s!="")
3025  {
3026    i = find(s,",");
3027    "i"; i;
3028    if (i==0)
3029    {
3030      i = find(s,";");
3031      if (i==0)
3032      {
3033        // no ; ??
3034         "no colon or semicolon found anymore";
3035         return(I);
3036      }
3037      // no "," but ";" on the i-th place
3038      t = s[1..i-1];
3039      s = "";
3040      "processing: "; t;
3041      tst = sent2lplace(t);
3042      tst;
3043      I = I, strList2poly(tst);
3044      return(I);
3045    }
3046    // here i !=0
3047    t = s[1..i-1];
3048    s = s[i+1..size(s)];
3049    "processing: "; t;
3050    tst = sent2lplace(t);
3051    tst;
3052    I = I, strList2poly(tst);
3053  }
3054  return(I);
3055}
3056example
3057{
3058  "EXAMPLE:"; echo = 2;
3059  ring r =0,(x,y,z,t),dp;
3060  def A = makeLetterplaceRing(4);
3061  setring A;
3062  string fn = "myfile";
3063  string s1 = "z*y*y*y - 3*y*z*x*y  + 3*y*y*z*y - y*x*y*z,";
3064  string s2 = "-2*y*x*y*z + y*y*z*z - z*z*y*y + 2*z*y*z*y,";
3065  string s3 = "z*y*x*t - 2*y*z*y*t + y*y*z*t - t*z*y*y + 2*t*y*z*y - t*x*y*z;";
3066  write(":w "+fn,s1);  write(":a "+fn,s2);   write(":a "+fn,s3);
3067  read(fn);
3068  ideal I = file2lplace(fn);
3069  I;
3070}
3071
3072/* EXAMPLES AGAIN:
3073//static proc get_ls3nilp()
3074{
3075//first app of file2lplace
3076  ring r =0,(x,y,z,t),dp;
3077  int d = 10;
3078  def A = makeLetterplaceRing(d);
3079  setring A;
3080  ideal I = file2lplace("./ls3nilp.bg");
3081  // and now test the correctness: go back from lplace to strings
3082  lp2lstr(I,r);
3083  setring r;
3084  lst2str(@LN,1); // agree!
3085}
3086
3087//static proc doc_example()
3088{
3089  LIB "freegb.lib";
3090  ring r = 0,(x,y,z),dp;
3091  int d =4; // degree bound
3092  def R = makeLetterplaceRing(d);
3093  setring R;
3094  ideal I = x(1)*y(2) + y(1)*z(2), x(1)*x(2) + x(1)*y(2) - y(1)*x(2) - y(1)*y(2);
3095  option(redSB);option(redTail);
3096  ideal J = system("freegb",I,d,nvars(r));
3097  J;
3098  // visualization:
3099  lp2lstr(J,r); // export an object called @LN to the ring r
3100  setring r;  // change to the ring r
3101  lst2str(@LN,1); // output the strings
3102}
3103
3104*/
3105
3106static proc lpMultX(poly f, poly g)
3107{
3108  /* multiplies two polys in a very general setting correctly */
3109  /* alternative to lpMult, possibly better at non-positive orderings */
3110
3111  if (lpAssumeViolation())
3112  {
3113    ERROR("Incomplete Letterplace structure on the basering!");
3114  }
3115  // decompose f,g into graded pieces with inForm: need dmodapp.lib
3116  int b = attrib(basering,"lV");  // the length of the block
3117  intvec w; // inherit the graded on the oridinal ring
3118  int i;
3119  for(i=1; i<=b; i++)
3120  {
3121    w[i] = deg(var(i));
3122  }
3123  intvec v = w;
3124  for(i=1; i< attrib(basering,"uptodeg"); i++)
3125  {
3126    v = v,w;
3127  }
3128  w = v;
3129  poly p,q,s, result;
3130  s = g;
3131  while (f!=0)
3132  {
3133    p = inForm(f,w)[1];
3134    f = f - p;
3135    s = g;
3136    while (s!=0)
3137    {
3138      q = inForm(s,w)[1];
3139      s = s - q;
3140      result = result + lpMult(p,q);
3141    }
3142  }
3143  // shrinking
3144  //  result;
3145  return( system("shrinktest",result,attrib(basering, "lV")) );
3146}
3147example
3148{
3149  "EXAMPLE:"; echo = 2;
3150  // define a ring in letterplace form as follows:
3151  ring r = 0,(x(1),y(1),x(2),y(2),x(3),y(3),x(4),y(4)),dp;
3152  def R = setLetterplaceAttributes(r,4,2); // supply R with letterplace structure
3153  setring R;
3154  poly a = x(1)*y(2)+x(1)+y(1); poly b = y(1)+3;
3155  lpMultX(b,a);
3156  lpMultX(a,b);
3157}
3158
3159// multiply two letterplace polynomials, lpMult: done
3160// reduction/ Normalform? needs kernel stuff
3161
3162
3163proc lpMult(poly f, poly g)
3164"USAGE:  lpMult(f,g); f,g letterplace polynomials
3165RETURN:  poly
3166ASSUME: basering has a letterplace ring structure
3167PURPOSE: compute the letterplace form of f*g
3168EXAMPLE: example lpMult; shows examples
3169"
3170{
3171
3172  // changelog:
3173  // VL oct 2010: deg -> deg(_,w) for the length
3174  // shrink the result => don't need to decompose polys
3175  // since the shift is big enough
3176
3177  // indeed it's better to have that
3178  // ASSUME: both f and g are quasi-homogeneous
3179
3180  if (lpAssumeViolation())
3181  {
3182    ERROR("Incomplete Letterplace structure on the basering!");
3183  }
3184  intvec w = 1:nvars(basering);
3185  int sf = deg(f,w); // VL Oct 2010: we need rather length than degree
3186  int sg = deg(g,w); // esp. in the case of weighted ordering
3187  int uptodeg = attrib(basering, "uptodeg");
3188  if (sf+sg > uptodeg)
3189  {
3190    ERROR("degree bound violated by the product!");
3191  }
3192  //  if (sf>1) { sf = sf -1; }
3193  poly v = f*shiftPoly(g,sf);
3194  // bug, reported by Simon King: in nonhomog case [solved]
3195  // we need to shrink
3196  return( system("shrinktest",v,attrib(basering, "lV")) );
3197}
3198example
3199{
3200  "EXAMPLE:"; echo = 2;
3201  // define a ring in letterplace form as follows:
3202  ring r = 0,(x(1),y(1),x(2),y(2),x(3),y(3),x(4),y(4)),dp;
3203  def R = setLetterplaceAttributes(r,4,2); // supply R with letterplace structure
3204  setring R;
3205  poly a = x(1)*y(2)+x(1)+y(1); poly b = y(1)+3;
3206  lpMult(b,a);
3207  lpMult(a,b);
3208}
3209
3210proc lpPower(poly f, int n)
3211"USAGE:  lpPower(f,n); f letterplace polynomial, int n
3212RETURN:  poly
3213ASSUME: basering has a letterplace ring structure
3214PURPOSE: compute the letterplace form of f^n
3215EXAMPLE: example lpPower; shows examples
3216"
3217{
3218  if (n<0) { ERROR("the power must be a natural number!"); }
3219  if (n==0) { return(poly(1)); }
3220  if (n==1) { return(f); }
3221  int i;
3222  poly p = 1;
3223  for(i=1; i<= n; i++)
3224  {
3225    p = lpMult(p,f);
3226  }
3227  return(p);
3228}
3229example
3230{
3231  "EXAMPLE:"; echo = 2;
3232  // define a ring in letterplace form as follows:
3233  ring r = 0,(x(1),y(1),x(2),y(2),x(3),y(3),x(4),y(4)),dp;
3234  def R = setLetterplaceAttributes(r,4,2); // supply R with letterplace structure
3235  setring R;
3236  poly a = x(1)*y(2) + y(1); poly b = y(1) - 1;
3237  lpPower(a,2);
3238  lpPower(b,4);
3239}
3240
3241// new: lp normal from by using shift-invariant data by Grischa Studzinski
3242
3243/////////////////////////////////////////////////////////
3244//   ASSUMPTIONS: every polynomial is an element of V',
3245//@* else there wouldn't be an dvec representation
3246
3247//Main procedure for the user
3248
3249proc lpNF(poly p, ideal G)
3250"USAGE: lpNF(p,G); f letterplace polynomial, ideal I
3251RETURN: poly
3252PURPOSE: computation of the normal form of p with respect to G
3253ASSUME: p is a Letterplace polynomial, G is a set Letterplace polynomials,
3254being a Letterplace Groebner basis (no check for this will be done)
3255NOTE: Strategy: take the smallest monomial wrt ordering for reduction
3256-     For homogenous ideals the shift does not matter
3257-     For non-homogenous ideals the first shift will be the smallest monomial
3258EXAMPLE: example lpNF; shows examples
3259"
3260{if ((p==0) || (size(G) == 0)){return(p);}
3261 checkAssumptions(p,G);
3262 G = sort(G)[1];
3263 list L = makeDVecI(G);
3264 return(lpNormalForm2(p,G,L));
3265}
3266example
3267{
3268  "EXAMPLE:"; echo = 2;
3269ring r = 0,(x,y,z),dp;
3270int d =5; // degree
3271def R = makeLetterplaceRing(d);
3272setring R;
3273ideal I = y(1)*x(2)*y(3) - z(1)*y(2)*z(3), x(1)*y(2)*x(3) - z(1)*x(2)*y(3), z(1)*x(2)*z(3) - y(1)*z(2)*x(3), x(1)*x(2)*x(3) + y(1)*y(2)*y(3) + z(1)*z(2)*z(3) + x(1)*y(2)*z(3);
3274ideal J = letplaceGBasis(I); // compute a Letterplace Groebner basis
3275poly p = y(1)*x(2)*y(3)*z(4)*y(5) - y(1)*z(2)*z(3)*y(4) + z(1)*y(2)*z(3);
3276poly q = z(1)*x(2)*z(3)*y(4)*z(5) - y(1)*z(2)*x(3)*y(4)*z(5);
3277lpNF(p,J);
3278lpNF(q,J);
3279}
3280
3281proc lpDivision(poly p, ideal I)
3282"ASSUME: I is a Groebner basis G = {g1,...,gN}, the original ring of the Letterplace ring has the name 'r' and no variable is called 'tag_i' for i in 1...N
3283RETURN: list L
3284NOTE: - L[1] NF(p,I)
3285      - L[2] list of expressions [i,l_{ij},r_{ij}] with \sum_{i,j} l_{ij} g_i r_{ij} = p - NF(p,I)
3286"
3287{
3288  if (p == 0 || size(I) == 0) {
3289    list L = 0;
3290    list empty;
3291    L[2] = empty;
3292    return (L);
3293  }
3294  poly pNF = lpNF(p,I);
3295  p = p - pNF;
3296
3297  // make new ring
3298  def save = basering;
3299  int norigvars = nvars(r);
3300  string tagvarstr = "(tag_1";
3301  for (int i = 2; i <= size(I); i++) {
3302    tagvarstr = tagvarstr + ",tag_" + string(i);
3303  } kill i;
3304  tagvarstr = tagvarstr + ")";
3305  string field = string(ringlist(r)[1]);
3306  execute("ring @tags = " + field + "," + tagvarstr + ",dp;");
3307  ring @tagged = (r + @tags);
3308  def @R = makeLetterplaceRing(attrib(save,"uptodeg")); setring @R;
3309
3310  // restore vars
3311  poly p = imap(save, p);
3312  poly pNF = imap(save, pNF);
3313  ideal I = imap(save, I);
3314  for (int i = 1; i <= size(I); i++) {
3315    I[i] = I[i] - var(norigvars + i);
3316  } kill i;
3317
3318  list summands;
3319  list L = pNF;
3320  poly pTaggedNF = lpNF(p,I);
3321  pTaggedNF = pTaggedNF * leadcoef(p); // somehow pTaggedNF gets normalized
3322  for (int i = 1; i <= size(pTaggedNF); i++) {
3323    intvec iv = lp2iv(pTaggedNF[i]);
3324    for (int j = 1; j <= size(iv); j++) {
3325      if (iv[j] > norigvars) {
3326        intvec left;
3327        intvec right;
3328        if (j > 1) {
3329          left = iv[1..(j-1)];
3330        }
3331        if (j < size(iv)) {
3332          right = iv[(j+1)..size(iv)];
3333        }
3334        list summand = (iv[j] - norigvars), leadcoef(pTaggedNF[i])*iv2lp(left), iv2lp(right);
3335        summands = insert(summands, summand, size(summands));
3336
3337        kill left;
3338        kill right;
3339        kill summand;
3340        break;
3341      }
3342    } kill j;
3343    kill iv;
3344  } kill i;
3345
3346  L[2] = summands;
3347
3348  setring save;
3349  list L = imap(@R,L);
3350  return (L);
3351}
3352
3353proc lpGBPres2Poly(list L, ideal I)
3354"ASSUME: L is a valid Groebner presentation as for example the result of lpDivision
3355RETURN: poly
3356NOTE: computes p = \sum_{i,j} l_{ij} g_i r_{ij} + NF(p,I) = \sum_{i} L[2][i][2] I[L[2][i][1]] L[2][i][3] + L[1]
3357"
3358{
3359  poly p;
3360  for (int i = 1; i <= size(L[2]); i++) {
3361    p = p + lpMult(lpMult(L[2][i][2], I[L[2][i][1]]), L[2][i][3]);
3362  }
3363  p = p + L[1];
3364  return (p);
3365}
3366
3367//procedures to convert monomials into the DVec representation, all static
3368////////////////////////////////////////////////////////
3369
3370
3371static proc getExpVecs(ideal G)
3372"USUAGE: getExpVecs(G);
3373RETURN: list of intvecs
3374PURPOSE: convert G into a list of intvecs, corresponding to the exponent vector
3375 of the leading monomials of G
3376"
3377{int i; list L;
3378 for (i = 1; i <= size(G); i++) {L[i] = leadexp(G[i]); }
3379 return(L);
3380}
3381
3382static proc delSupZero(intvec I)
3383"USUAGE:delSupZero(I);
3384RETURN: intvec
3385PURPOSE: Deletes superfluous zero blocks of an exponent vector
3386ASSUME: Intvec is an exponent vector of a letterplace monomial contained in V'
3387"
3388{if (I==intvec(0)) {return(intvec(0));}
3389 int j,k,l;
3390 int n = attrib(basering,"lV"); int d = attrib(basering,"uptodeg");
3391 intvec w; j = 1;
3392 while (j <= d)
3393 {w = I[1..n];
3394  if (w<>intvec(0)){break;}
3395   else {I = I[(n+1)..(n*d)]; d = d-1; j++;}
3396 }
3397 for (j = 1; j <= d; j++)
3398  {l=(j-1)*n+1; k= j*n;
3399   w = I[l..k];
3400   if (w==intvec(0)){w = I[1..(l-1)]; return(w);}//if a zero block is found there are only zero blocks left,
3401                                                 //otherwise there would be a hole in the monomial
3402                                                 // shrink should take care that this will not happen
3403  }
3404 return(I);
3405}
3406
3407static proc delSupZeroList(list L)
3408"USUAGE:delSupZeroList(L); L a list, containing intvecs
3409RETURN: list, containing intvecs
3410PURPOSE: Deletes all superfluous zero blocks for a list of exponent vectors
3411ASSUME: All intvecs are exponent vectors of letterplace monomials contained in V'
3412"
3413{int i;
3414 for (i = size(L); 0 < i; i--){L[i] = delSupZero(L[i]);}
3415 return(L);
3416}
3417
3418
3419static proc makeDVec(intvec V)
3420"USUAGE:makeDVec(V);
3421RETURN: intvec
3422PURPOSE: Converts an modified exponent vector into an Dvec
3423NOTE: Superfluos zero blocks must have been deleted befor using this procedure
3424"
3425{int i,j,k,r1,r2; intvec D;
3426 int n = attrib(basering,"lV");
3427 k = size(V) div n; r1 = 0; r2 = 0;
3428 for (i=1; i<= k; i++)
3429  {for (j=(1+((i-1)*n)); j <= (i*n); j++)
3430   {if (V[j]>0){r2 = j - ((i-1)*n); j = (j mod n); break;}
3431   }
3432   D[size(D)+1] = r1+r2;
3433   if (j == 0) {r1 = 0;} else{r1= n-j;}
3434  }
3435 D = D[2..size(D)];
3436 return(D);
3437}
3438
3439static proc makeDVecL(list L)
3440"USUAGE:makeDVecL(L); L, a list containing intvecs
3441RETURN: list, containing intvecs
3442ASSUME:
3443"
3444{int i; list R;
3445 for (i=1; i <= size(L); i++) {R[i] = makeDVec(L[i]);}
3446 return(R);
3447}
3448
3449static proc makeDVecI(ideal G)
3450"USUAGE:makeDVecI(G);
3451RETURN:list, containing intvecs
3452PURPOSE:computing the DVec representation for lead(G)
3453ASSUME:
3454"
3455{list L = delSupZeroList(getExpVecs(G));
3456 return(makeDVecL(L));
3457}
3458
3459
3460//procedures, which are dealing with the DVec representation, all static
3461
3462static proc dShiftDiv(intvec V, intvec W)
3463"USUAGE: dShiftDiv(V,W);
3464RETURN: a list,containing integers, or -1, if no shift of W divides V
3465PURPOSE: find all possible shifts s, such that s.W|V
3466ASSUME: V,W are DVecs of monomials contained in V'
3467"
3468{if(size(V)<size(W)){return(list(-1));}
3469
3470 int i,j,r; intvec T; list R;
3471 int n = attrib(basering,"lV");
3472 int k = size(V) - size(W) + 1;
3473 if (intvec(V[1..size(W)])-W == 0){R[1]=0;}
3474 for (i =2; i <=k; i++)
3475 {r = 0; kill T; intvec T;
3476  for (j =1; j <= i; j++) {r = r + V[j];}
3477  //if (i==1) {T[1] = r-(i-1)*n;} else
3478  T[1] = r-(i-1)*n; if (size(W)>1) {T[2..size(W)] = V[(i+1)..(size(W)+i-1)];}
3479  if (T-W == 0) {R[size(R)+1] = i-1;}
3480 }
3481 if (size(R)>0) {return(R);}
3482 else {return(list(-1));}
3483}
3484
3485//the first normal form procedure, if a user want not to presort the ideal, just make it not static
3486
3487static proc lpNormalForm1(poly p, ideal G, list L)
3488"USUAGE:lpNormalForm1(p,G);
3489RETURN:poly
3490PURPOSE:computation of the normalform of p w.r.t. G
3491ASSUME: p is a Letterplace polynomial, G is a set of Letterplace polynomials
3492NOTE: Taking the first possible reduction
3493"
3494{
3495 if (deg(p) <1) {return(p);}
3496 else
3497  {
3498  int i; int s;
3499  intvec V = makeDVec(delSupZero(leadexp(p)));
3500  for (i = 1; i <= size(L); i++)
3501  {s = dShiftDiv(V, L[i])[1];
3502   if (s <> -1)
3503   {p = lpReduce(p,G[i],s);
3504    p = lpNormalForm1(p,G,L);
3505    break;
3506   }
3507  }
3508  p = p[1] + lpNormalForm1(p-p[1],G,L);
3509  return(p);
3510 }
3511}
3512
3513
3514// new VL; called from lpNF
3515static proc lpNormalForm2(poly pp, ideal G, list L)
3516"USUAGE:lpNormalForm2(p,G);
3517RETURN:poly
3518PURPOSE:computation of the normal form of p w.r.t. G
3519ASSUME: p is a Letterplace polynomial, G is a set of Letterplace polynomials
3520NOTE: Taking the first possible reduction
3521"
3522{
3523 poly one = 1;
3524 if ( (pp == 0) || (leadmonom(pp) == one) ) { return(pp); }
3525 poly p = pp; poly q;
3526 int i; int s; intvec V;
3527 while ( (p != 0) && (leadmonom(p) != one) )
3528 {
3529   //"entered while with p="; p;
3530   V = makeDVec(delSupZero(leadexp(p)));
3531   i = 0;
3532   s = -1;
3533   //"look for divisor";
3534   while ( (s == -1) && (i<size(L)) )
3535   {
3536     i = i+1;
3537     s = dShiftDiv(V, L[i])[1];
3538   }
3539 // now, out of here: either i=size(L) and s==-1 => no reduction
3540 // otherwise: i<=size(L) and s!= -1 => reduction
3541    //"out of divisor search: s="; s; "i="; i;
3542    if (s != -1)
3543    {
3544    //"start reducing with G[i]:";
3545      p = lpReduce(p,G[i],s); // lm-reduction
3546      //"reduced to p="; p;
3547    }
3548    else
3549    {
3550      // ie no lm-reduction possible; proceed with the tail reduction
3551      q = p-lead(p);
3552      p = lead(p);
3553      if (q!=0)
3554      {
3555        p = p + lpNormalForm2(q,G,L);
3556      }
3557      return(p);
3558    }
3559 }
3560 // out of while when p==0 or p == const
3561 return(p);
3562}
3563
3564
3565
3566
3567//secundary procedures, all static
3568
3569static proc getlpCoeffs(poly q, poly p)
3570{list R; poly m; intvec cq,t,lv,rv,bla;
3571 int n = attrib(basering,"lV"); int d = attrib(basering,"uptodeg");
3572 int i;
3573 m = p/q;
3574 cq = leadexp(m);
3575 for (i = 1; i<= d; i++)
3576 {bla = cq[((i-1)*n+1)..(i*n)];
3577  if (bla == 0) {lv = cq[1..i*n]; cq = cq[(i*n+1)..(d*n)]; break;}
3578 }
3579
3580 d = size(cq) div n;
3581 for (i = 1; i<= d; i++)
3582 {bla = cq[((i-1)*n+1)..(i*n)];
3583  if (bla <> 0){rv = cq[((i-1)*n+1)..(d*n)]; break;}
3584 }
3585 return(list(monomial(lv),monomial(rv)));
3586}
3587
3588static proc lpReduce(poly p, poly g, int s)
3589"NOTE: shift can not exceed the degree bound, because s*g | p
3590"
3591{poly l,r,qt; int i;
3592 list K = getlpCoeffs(lead(shiftPoly(g,s)), lead(p));
3593 l = K[1]; r = K[2];
3594 kill K;
3595 for (i = 1; i <= size(g); i++)
3596 {
3597   qt = qt + lpMult(lpMult(l,g[i]),r);
3598 }
3599 return(p - leadcoef(p)*normalize(qt));
3600}
3601
3602
3603static proc lpShrink(poly p)
3604"
3605"
3606{int n;
3607 if (typeof(attrib(basering,"isLetterplaceRing"))=="int")
3608 {n = attrib(basering,"lV");
3609  return(system("shrinktest",p,n));
3610 }
3611 else {ERROR("Basering is not a Letterplace ring!");}
3612}
3613
3614static proc checkAssumptions(poly p, ideal G)
3615"
3616"
3617{checkLPRing();
3618 checkAssumptionPoly(p);
3619 checkAssumptionIdeal(G);
3620 return();
3621}
3622
3623static proc checkLPRing();
3624"
3625"
3626{if (typeof(attrib(basering,"isLetterplaceRing"))=="string") {ERROR("Basering is not a Letterplace ring!");}
3627 return();
3628}
3629
3630static proc checkAssumptionIdeal(ideal G)
3631"PURPOSE:Check if all elements of ideal are elements of V'
3632"
3633{ideal L = lead(normalize(G));
3634 int i;
3635 for (i = 1; i <= ncols(G); i++) {if (!isContainedInVp(G[i])) {ERROR("Ideal containes elements not contained in V'");}}
3636 return();
3637}
3638
3639static proc checkAssumptionPoly(poly p)
3640"PURPOSE:Check if p is an element of V'
3641"
3642{poly l = lead(normalize(p));
3643 if (!isContainedInVp(l)) {ERROR("Polynomial is not contained in V'");}
3644 return();
3645}
3646
3647static proc isContainedInVp(poly p)
3648"PURPOSE: Check monomial for holes in the places
3649"
3650{int r = 0; intvec w;
3651 intvec l = leadexp(p);
3652 int n = attrib(basering,"lV"); int d = attrib(basering,"uptodeg");
3653 int i,j,c,c1;
3654 while (1 <= d)
3655 {w = l[1..n];
3656  if (w<>intvec(0)){break;}
3657   else {l = l[(n+1)..(n*d)]; d = d-1;}
3658 }
3659
3660 while (1 <= d)
3661  {for (j = 1; j <= n; j++)
3662   {if (l[j]<>0)
3663    {if (c1<>0){return(0);}
3664     if (c<>0){return(0);}
3665     if (l[j]<>1){return(0);}
3666     c=1;
3667    }
3668   }
3669   if (c == 0){c1=1;if (1 < d){l = l[(n+1)..(n*d)]; d = d-1;} else {d = d -1;}}
3670    else {c = 0; if (1 < d){l = l[(n+1)..(n*d)]; d = d-1;} else {d = d -1;}}
3671  }
3672 return(1);
3673}
3674
3675// under development for Roberto
3676static proc extractLinearPart(module M)
3677{
3678  /* returns vectors from a module whose max leadexp is 1 */
3679  /* does not take place nonlinearity into account yet */
3680  /* use rather kernel function isinV to get really nonlinear things */
3681  int i; int s = ncols(M);
3682  int answer = 1;
3683  vector v; module Ret;
3684  for(i=1; i<=s; i++)
3685  {
3686    if ( isLinearVector(M[i]) )
3687    {
3688      Ret = Ret, M[i];
3689    }
3690  }
3691  Ret = simplify(Ret,2);
3692  return(Ret);
3693}
3694
3695// under development for Roberto
3696static proc isLinearVector(vector v)
3697{
3698  /* returns true iff max leadexp is 1 */
3699  int i,j,k;
3700  intvec w;
3701  int s = size(v);
3702  poly p;
3703  int answer = 1;
3704  for(i=1; i<=s; i++)
3705  {
3706    p = v[i];
3707    while (p != 0)
3708    {
3709      w = leadexp(p);
3710      j = Max(w);
3711      if (j >=2)
3712      {
3713        answer = 0;
3714        return(answer);
3715      }
3716      p = p-lead(p);
3717    }
3718  }
3719  return(answer);
3720}
3721
3722
3723// // the following is to determine a shift of a mono/poly from the
3724// // interface
3725
3726// static proc whichshift(poly p, int numvars)
3727// {
3728// // numvars = number of vars of the orig free algebra
3729// // assume: we are in the letterplace ring
3730// // takes  monomial on the input
3731// poly q = lead(p);
3732// intvec v = leadexp(v);
3733// if (v==0) { return(int(0)); }
3734// int sv = size(v);
3735// int i=1;
3736// while ( (v[i]==0) && (i<sv) ) { i++; }
3737// i = sv div i;
3738// return(i);
3739// }
3740
3741
3742// LIB "qhmoduli.lib";
3743// static proc polyshift(poly p,  int numvars)
3744// {
3745//   poly q = p; int i = 0;
3746//   while (q!=0)
3747//   {
3748//     i = Max(i, whichshift(q,numvars));
3749//     q = q - lead(q);
3750//   }
3751//   return(q);
3752// }
3753
3754static proc lpAssumeViolation()
3755{
3756  // checks whether the global vars
3757  // uptodeg and lV are defined
3758  // returns Boolean : yes/no [for assume violation]
3759  def lpring = attrib(basering,"isLetterplaceRing");
3760  if ( typeof(lpring)!="int" )
3761  {
3762    //  if ( typeof(lpring)=="string" ) ??
3763    // basering is NOT lp Ring
3764    return(1);
3765  }
3766  def uptodeg = attrib(basering,"uptodeg");
3767  if ( typeof(uptodeg)!="int" )
3768  {
3769    return(1);
3770  }
3771  def lV = attrib(basering,"lV");
3772  if ( typeof(lV)!="int" )
3773  {
3774    return(1);
3775  }
3776  //  int i = ( defined(uptodeg) && (defined(lV)) );
3777  //  return ( !i );
3778  return(0);
3779}
3780
3781static proc bugSKing()
3782{
3783  LIB "freegb.lib";
3784  ring r=0,(a,b),dp;
3785  def R = makeLetterplaceRing(5);
3786  setring R;
3787  poly p = a(1);
3788  poly q = b(1);
3789  poly p2 = lpPower(p,2);
3790  lpMult(p2+q,q)-lpMult(p2,q)-lpMult(q,q); // now its 0
3791}
3792
3793static proc bugRucker()
3794{
3795  // needs unstatic lpMultX
3796  LIB "freegb.lib";
3797  ring r=0,(a,b,c,d,p,q,r,s,t,u,v,w),(a(7,1,1,7),dp);
3798  def R=makeLetterplaceRing(20,1);
3799  setring R;
3800  option(redSB); option(redTail);
3801  ideal I=a(1)*b(2)*c(3)-p(1)*q(2)*r(3)*s(4)*t(5)*u(6),b(1)*c(2)*d(3)-v(1)*w(2);
3802  poly ttt = a(1)*v(2)*w(3)-p(1)*q(2)*r(3)*s(4)*t(5)*u(6)*d(7);
3803  // with lpMult
3804  lpMult(I[1],d(1)) - lpMult(a(1),I[2]); // spoly; has been incorrect before
3805  _ - ttt;
3806  // with lpMultX
3807  lpMultX(I[1],d(1)) - lpMultX(a(1),I[2]); // spoly; has been incorrect before
3808  _ - ttt;
3809}
3810
3811static proc checkWeightedExampleLP()
3812{
3813  ring r = 0,(x(1),y(1),x(2),y(2),x(3),y(3),x(4),y(4)),wp(2,1,2,1,2,1,2,1);
3814  def R = setLetterplaceAttributes(r,4,2); // supply R with letterplace structure
3815  setring R;
3816  poly a = x(1)*y(2)+x(1)+y(1); poly b = y(1)+3;
3817  lpMultX(b,a);
3818  lpMultX(a,b); // seems to work properly
3819}
3820
3821proc lpPrint(ideal I, def @r)
3822"USAGE: lpPrint(I, r); I an ideal, r a ring
3823RETURN: list of strings
3824PURPOSE: represent Letterplace ideal in the form of words
3825ASSUME: - basering is a Letterplace ring, r is the commutative ring
3826from which basering has been built
3827EXAMPLE: example lpPrint; shows example
3828"
3829{
3830        def save = basering;
3831        lp2lstr(I,@r); // export an object called @code{@LN} to the ring r
3832        setring @r;  // change to the ring r
3833        list @L = lst2str(@LN,1);
3834        export @L;
3835        setring save;
3836        list @@L = @L;
3837        setring @r;
3838        kill @L;
3839        kill @LN;
3840        setring save;
3841        return(@@L);
3842}
3843example
3844{
3845 "EXAMPLE:"; echo = 2;
3846 ring r = (0,a,b,g),(x,y),Dp;
3847 def R = makeLetterplaceRing(4); // constructs a Letterplace ring
3848 setring R; // downup algebra A
3849 ideal J = x(1)*x(2)*y(3)-a*x(1)*y(2)*x(3) - b*y(1)*x(2)*x(3) - g*x(1),
3850 x(1)*y(2)*y(3)-a*y(1)*x(2)*y(3) - b*y(1)*y(2)*x(3) - g*y(1);
3851 list L = lpPrint(J,r);
3852 L;
3853}
3854
3855/* THE FOLLOWING ARE UNDER DEVELOPMENT
3856// copied following from freegb_wrkcp.lib by Karim Abou Zeid on 07.04.2017:
3857// makeLetterplaceRingElim(int d)
3858// makeLetterplaceRingNDO(int d)
3859// setLetterplaceAttributesElim(def R, int uptodeg, int lV)
3860// lpElimIdeal(ideal I)
3861// makeLetterplaceRingWt(int d, intvec W)
3862
3863static proc makeLetterplaceRingElim(int d)
3864"USAGE:  makeLetterplaceRingElim(d); d integers
3865RETURN:  ring
3866PURPOSE: creates a ring with an elimination ordering
3867NOTE: the matrix for the ordering looks as follows: first row is 1,..,0,1,0,..
3868@* then 0,1,0,...,0,0,1,0... and so on, lastly its lp
3869@* this ordering is only correct if only polys with same shift are compared
3870EXAMPLE: example makeLetterplaceRingElim; shows examples
3871"
3872{
3873
3874  // ToDo future: inherit positive weights in the orig ring
3875  // complain on nonpositive ones
3876
3877  // d = up to degree, will be shifted to d+1
3878  if (d<1) {"bad d"; return(0);}
3879
3880  int uptodeg = d; int lV = nvars(basering);
3881
3882  int ppl = printlevel-voice+2;
3883  string err = "";
3884
3885  int i,j,s; intvec iV,iVl;
3886  def save = basering;
3887  int D = d-1;
3888  list LR  = ringlist(save);
3889  list L, tmp, tmp2, tmp3;
3890  L[1] = LR[1]; // ground field
3891  L[4] = LR[4]; // quotient ideal
3892  tmp  = LR[2]; // varnames
3893  s = size(LR[2]);
3894  for (i=1; i<=D; i++)
3895  {
3896    for (j=1; j<=s; j++)
3897    {
3898      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
3899    }
3900  }
3901  for (i=1; i<=s; i++)
3902  {
3903    tmp[i] = string(tmp[i])+"("+string(1)+")";
3904  }
3905  L[2] = tmp;
3906  L[3] = list();
3907  list OrigNames = LR[2];
3908  s = size(LR[3]);
3909  //creation of first block
3910
3911  if (s==2)
3912  {
3913    // not a blockord, 1 block + module ord
3914    tmp = LR[3][s]; // module ord
3915    for (i = 1; i <= lV;  i++)
3916    {
3917      iV = (0: lV);
3918      iV[i] = 1;
3919      iVl = iV;
3920      for (j = 1; j <= D; j++)
3921       { iVl = iVl,iV; }
3922      L[3][i] = list("a",iVl);
3923    }
3924//    for (i=1; i<=d; i++)
3925//    {
3926//      LR[3][s-1+i] = LR[3][1];
3927//    }
3928    //    LR[3][s+D] = tmp;
3929    //iV = (1:(d*lV));
3930    L[3][lV+1] = list("lp",(1:(d*lV)));
3931    L[3][lV+2] = tmp;
3932  }
3933  else {ERROR("Please set the ordering of basering to dp");}
3934//  if (s>2)
3935//  {
3936//    // there are s-1 blocks
3937//    int nb = s-1;
3938//    tmp = LR[3][s]; // module ord to place at the very end
3939//   tmp2 = LR[3]; tmp2 = tmp2[1..nb];
3940//    LR[3][1] = list("a",LTO);
3941//    //tmp3[1] = list("a",intvec(1: int(d*lV))); // deg-ord, insert as the 1st
3942//    for (i=1; i<=d; i++)
3943//    {
3944//      tmp3 = tmp3 + tmp2;
3945//    }
3946//    tmp3 = tmp3 + list(tmp);
3947//    LR[3] = tmp3;
3948//     for (i=1; i<=d; i++)
3949//     {
3950//       for (j=1; j<=nb; j++)
3951//       {
3952//         //        LR[3][i*nb+j+1]= LR[3][j];
3953//         LR[3][i*nb+j+1]= tmp2[j];
3954//       }
3955//     }
3956//     //    size(LR[3]);
3957//     LR[3][(s-1)*d+2] = tmp;
3958//     LR[3] = list("a",intvec(1: int(d*lV))) + LR[3]; // deg-ord, insert as the 1st
3959    // remove everything behind nb*(D+1)+1 ?
3960    //    tmp = LR[3];
3961    //    LR[3] = tmp[1..size(tmp)-1];
3962 // }
3963 // L[3] = LR[3];
3964  def @R = ring(L);
3965  //  setring @R;
3966  //  int uptodeg = d; int lV = nvars(basering); // were defined before
3967  def @@R = setLetterplaceAttributesElim(@R,uptodeg,lV);
3968  return (@@R);
3969}
3970example
3971{
3972  "EXAMPLE:"; echo = 2;
3973  ring r = 0,(x,y,z),lp;
3974  def A = makeLetterplaceRingElim(2);
3975  setring A;
3976  A;
3977  attrib(A,"isLetterplaceRing");
3978  attrib(A,"uptodeg");  // degree bound
3979  attrib(A,"lV"); // number of variables in the main block
3980}
3981
3982
3983
3984static proc makeLetterplaceRingNDO(int d)
3985"USAGE:  makeLetterplaceRingNDO(d); d an integer
3986RETURN:  ring
3987PURPOSE: creates a ring with a non-degree first ordering, suitable for
3988@* the use of non-homogeneous letterplace
3989NOTE: the matrix for the ordering looks as follows:
3990@*    'd' blocks of shifted original variables
3991EXAMPLE: example makeLetterplaceRingNDO; shows examples
3992"
3993{
3994
3995  // ToDo future: inherit positive weights in the orig ring
3996  // complain on nonpositive ones
3997
3998  // d = up to degree, will be shifted to d+1
3999  if (d<1) {"bad d"; return(0);}
4000
4001  int uptodeg = d; int lV = nvars(basering);
4002
4003  int ppl = printlevel-voice+2;
4004  string err = "";
4005
4006  int i,j,s;
4007  def save = basering;
4008  int D = d-1;
4009  list LR  = ringlist(save);
4010  list L, tmp, tmp2, tmp3;
4011  L[1] = LR[1]; // ground field
4012  L[4] = LR[4]; // quotient ideal
4013  tmp  = LR[2]; // varnames
4014  s = size(LR[2]);
4015  for (i=1; i<=D; i++)
4016  {
4017    for (j=1; j<=s; j++)
4018    {
4019      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
4020    }
4021  }
4022  for (i=1; i<=s; i++)
4023  {
4024    tmp[i] = string(tmp[i])+"("+string(1)+")";
4025  }
4026  L[2] = tmp;
4027  list OrigNames = LR[2];
4028  // ordering: one 1..1 a above
4029  // ordering: d blocks of the ord on r
4030  // try to get whether the ord on r is blockord itself
4031  // TODO: make L(2) ordering! exponent is maximally 2
4032  s = size(LR[3]);
4033  if (s==2)
4034  {
4035    // not a blockord, 1 block + module ord
4036    tmp = LR[3][s]; // module ord
4037    for (i=1; i<=d; i++)
4038    {
4039      LR[3][i] = LR[3][1];
4040    }
4041    //    LR[3][s+D] = tmp;
4042    LR[3][d+1] = tmp;
4043    //LR[3][1] = list("a",intvec(1: int(d*lV))); // deg-ord not wanted here
4044  }
4045  if (s>2)
4046  {
4047    // there are s-1 blocks
4048    int nb = s-1;
4049    tmp = LR[3][s]; // module ord to place at the very end
4050    tmp2 = LR[3]; tmp2 = tmp2[1..nb];
4051    //tmp3[1] = list("a",intvec(1: int(d*lV))); // deg-ord not wanted here
4052    for (i=1; i<=d; i++)
4053    {
4054      tmp3 = tmp3 + tmp2;
4055    }
4056    tmp3 = tmp3 + list(tmp);
4057    LR[3] = tmp3;
4058//     for (i=1; i<=d; i++)
4059//     {
4060//       for (j=1; j<=nb; j++)
4061//       {
4062//         //        LR[3][i*nb+j+1]= LR[3][j];
4063//         LR[3][i*nb+j+1]= tmp2[j];
4064//       }
4065//     }
4066//     //    size(LR[3]);
4067//     LR[3][(s-1)*d+2] = tmp;
4068//     LR[3] = list("a",intvec(1: int(d*lV))) + LR[3]; // deg-ord, insert as the 1st
4069    // remove everything behind nb*(D+1)+1 ?
4070    //    tmp = LR[3];
4071    //    LR[3] = tmp[1..size(tmp)-1];
4072  }
4073  L[3] = LR[3];
4074  def @R = ring(L);
4075  //  setring @R;
4076  //  int uptodeg = d; int lV = nvars(basering); // were defined before
4077  def @@R = setLetterplaceAttributes(@R,uptodeg,lV);
4078  return (@@R);
4079}
4080example
4081{
4082  "EXAMPLE:"; echo = 2;
4083  ring r = 0,(x,y,z),lp;
4084  def A = makeLetterplaceRingNDO(2);
4085  setring A;
4086  A;
4087  attrib(A,"isLetterplaceRing");
4088  attrib(A,"uptodeg");  // degree bound
4089  attrib(A,"lV"); // number of variables in the main block
4090}
4091
4092static proc setLetterplaceAttributesElim(def R, int uptodeg, int lV)
4093"USAGE: setLetterplaceAttributesElim(R, d, b, eV); R a ring, b,d, eV integers
4094RETURN: ring with special attributes set
4095PURPOSE: sets attributes for a letterplace ring:
4096@*      'isLetterplaceRing' = true, 'uptodeg' = d, 'lV' = b, 'eV' = eV, where
4097@*      'uptodeg' stands for the degree bound,
4098@*      'lV' for the number of variables in the block 0
4099@*      'eV' for the number of elimination variables
4100NOTE: Activate the resulting ring by using @code{setring}
4101"
4102{
4103  if (uptodeg*lV != nvars(R))
4104  {
4105    ERROR("uptodeg and lV do not agree on the basering!");
4106  }
4107
4108
4109    // Set letterplace-specific attributes for the output ring!
4110  attrib(R, "uptodeg", uptodeg);
4111  attrib(R, "lV", lV);
4112  attrib(R, "isLetterplaceRing", 1);
4113  attrib(R, "HasElimOrd", 1);
4114  return (R);
4115}
4116example
4117{
4118  "EXAMPLE:"; echo = 2;
4119  ring r = 0,(x(1),y(1),x(2),y(2),x(3),y(3),x(4),y(4)),dp;
4120  def R = setLetterplaceAttributesElim(r, 4, 2, 1); setring R;
4121  attrib(R,"isLetterplaceRing");
4122  lieBracket(x(1),y(1),2);
4123}
4124
4125
4126static proc lpElimIdeal(ideal I)
4127"
4128does not work for degree reasons (deg function does not work for lp rings -> newone!)
4129"
4130{
4131  def lpring = attrib(basering,"isLetterplaceRing");
4132  def lpEO =  attrib(basering,"HasElimOrd");
4133  if ( typeof(lpring)!="int" && typeof(lpEO)!="int")
4134  {
4135    ERROR("Ring is not a lp-ring with an elimination ordering");
4136  }
4137
4138  //int nE = attrib(basering, "eV");
4139
4140  return(letplaceGBasis(I));
4141}
4142
4143
4144static proc makeLetterplaceRingWt(int d, intvec W)
4145"USAGE:  makeLetterplaceRingWt(d,W); d an integer, W a vector of positive integers
4146RETURN:  ring
4147PURPOSE: creates a ring with a special ordering, suitable for
4148@* the use of non-homogeneous letterplace
4149NOTE: the matrix for the ordering looks as follows: first row is W,W,W,...
4150@* then there come 'd' blocks of shifted original variables
4151EXAMPLE: example makeLetterplaceRing2; shows examples
4152"
4153{
4154
4155  // ToDo future: inherit positive weights in the orig ring
4156  // complain on nonpositive ones
4157
4158  // d = up to degree, will be shifted to d+1
4159  if (d<1) {"bad d"; return(0);}
4160
4161  int uptodeg = d; int lV = nvars(basering);
4162
4163  //check weightvector
4164  if (size(W) <> lV) {"bad weights"; return(0);}
4165
4166  int i;
4167  for (i = 1; i <= size(W); i++) {if (W[i] < 0) {"bad weights"; return(0);}}
4168  intvec Wt = W;
4169  for (i = 2; i <= d; i++) {Wt = Wt, W;}
4170  kill i;
4171
4172  int ppl = printlevel-voice+2;
4173  string err = "";
4174
4175  int i,j,s;
4176  def save = basering;
4177  int D = d-1;
4178  list LR  = ringlist(save);
4179  list L, tmp, tmp2, tmp3;
4180  L[1] = LR[1]; // ground field
4181  L[4] = LR[4]; // quotient ideal
4182  tmp  = LR[2]; // varnames
4183  s = size(LR[2]);
4184  for (i=1; i<=D; i++)
4185  {
4186    for (j=1; j<=s; j++)
4187    {
4188      tmp[i*s+j] = string(tmp[j])+"("+string(i+1)+")";
4189    }
4190  }
4191  for (i=1; i<=s; i++)
4192  {
4193    tmp[i] = string(tmp[i])+"("+string(1)+")";
4194  }
4195  L[2] = tmp;
4196  list OrigNames = LR[2];
4197  // ordering: one 1..1 a above
4198  // ordering: d blocks of the ord on r
4199  // try to get whether the ord on r is blockord itself
4200  // TODO: make L(2) ordering! exponent is maximally 2
4201  s = size(LR[3]);
4202  if (s==2)
4203  {
4204    // not a blockord, 1 block + module ord
4205    tmp = LR[3][s]; // module ord
4206    for (i=1; i<=d; i++)
4207    {
4208      LR[3][s-1+i] = LR[3][1];
4209    }
4210    //    LR[3][s+D] = tmp;
4211    LR[3][s+1+D] = tmp;
4212    LR[3][1] = list("a",Wt); // deg-ord
4213  }
4214  if (s>2)
4215  {
4216    // there are s-1 blocks
4217    int nb = s-1;
4218    tmp = LR[3][s]; // module ord to place at the very end
4219    tmp2 = LR[3]; tmp2 = tmp2[1..nb];
4220    tmp3[1] = list("a",Wt); // deg-ord, insert as the 1st
4221    for (i=1; i<=d; i++)
4222    {
4223      tmp3 = tmp3 + tmp2;
4224    }
4225    tmp3 = tmp3 + list(tmp);
4226    LR[3] = tmp3;
4227
4228  }
4229  L[3] = LR[3];
4230  def @R = ring(L);
4231  //  setring @R;
4232  //  int uptodeg = d; int lV = nvars(basering); // were defined before
4233  def @@R = setLetterplaceAttributes(@R,uptodeg,lV);
4234  return (@@R);
4235}
4236example
4237{
4238  "EXAMPLE:"; echo = 2;
4239  ring r = 0,(x,y,z),(dp(1),dp(2));
4240  def A = makeLetterplaceRingWt(2,intvec(1,2,3));
4241  setring A;
4242  A;
4243  attrib(A,"isLetterplaceRing");
4244  attrib(A,"uptodeg");  // degree bound
4245  attrib(A,"lV"); // number of variables in the main block
4246}
4247*/
Note: See TracBrowser for help on using the repository browser.