source: git/Singular/LIB/nctools.lib @ a57b65

spielwiese
Last change on this file since a57b65 was 3686937, checked in by Oleksandr Motsak <motsak@…>, 11 years ago
Added '$Id$' as a comment to all libs (LIB/*.lib)
  • Property mode set to 100644
File size: 46.8 KB
Line 
1//////////////////////////////////////////////////////////////////////////////
2version="version nctools.lib 4.0.0.0 Jun_2013 "; // $Id$
3category="Noncommutative";
4info="
5LIBRARY: nctools.lib     General tools for noncommutative algebras
6AUTHORS:   Levandovskyy V.,     levandov@mathematik.uni-kl.de,
7@*         Lobillo, F.J.,       jlobillo@ugr.es,
8@*         Rabelo, C.,          crabelo@ugr.es,
9@*         Motsak, O.,          U@D, where U={motsak}, D={mathematik.uni-kl.de}
10
11
12OVERVIEW:
13Support: DFG (Deutsche Forschungsgesellschaft) and Metodos algebraicos y efectivos
14en grupos cuanticos, BFM2001-3141, MCYT, Jose Gomez-Torrecillas (Main researcher).
15
16PROCEDURES:
17Gweights(r);              compute weights for a compatible ordering in a G-algebra,
18weightedRing(r);          change the ordering of a ring to a weighted one,
19ndcond();                 the ideal of non-degeneracy conditions in G-algebra,
20Weyl([p]);                create Weyl algebra structure in a basering (two different realizations),
21makeWeyl(n, [p]);         return n-th Weyl algebra in (x(i),D(i)) presentation,
22makeHeisenberg(N, [p,d]); return n-th Heisenberg algebra in (x(i),y(i),h) realization,
23Exterior();               return qring, the exterior algebra of a basering,
24findimAlgebra(M,[r]);     create finite dimensional algebra structure from the basering and the multiplication matrix M,
25superCommutative([b,e,Q]);  return qring, a super-commutative algebra over a basering,
26rightStd(I);              compute right Groebner basis of an ideal,
27rightNF(f,I);             compute right normal form wrt a submodule,
28rightModulo(M,N);         compute kernel of a homomorphism of right modules,
29moduloSlim(A,B);     compute modulo command via slimgb
30ncRelations(r);      recover the non-commutative relations of a G-algebra,
31isCentral(p);        check for the commutativity of a polynomial in the G-algebra,
32isNC();              check whether basering is noncommutative,
33isCommutative();     check whether basering is commutative
34isWeyl();            check whether basering is a Weyl algebra
35UpOneMatrix();       return NxN matrix with 1's in the whole upper triagle,
36AltVarStart();       return first alternating variable of a super-commutative algebra,
37AltVarEnd();         return last alternating variable of a super-commutative algebra,
38IsSCA();             check whether current ring is a super-commutative algebra,
39makeModElimRing(R);  equip a ring with module elimination ordering,
40embedMat(M,m,n);     embeds matrix M in a left upper corner of m times n matrix
41";
42
43
44LIB "ring.lib"; // for rootofUnity
45LIB "poly.lib"; // for newtonDiag
46LIB "matrix.lib"; // for submat
47
48///////////////////////////////////////////////////////////////////////////////
49
50// This procedure computes a weights vector for a G-algebra r
51
52proc Gweights(def r)
53"USAGE:   Gweights(r); r a ring or a square matrix
54RETURN:   intvec
55PURPOSE: compute an appropriate weight int vector for a G-algebra, i.e., such that
56\foral\;i<j\;\;lm_w(d_{ij}) <_w x_i x_j.
57@*       the polynomials d_{ij} are taken from r itself, if it is of the type ring
58@*       or defined by the given square polynomial matrix
59THEORY:   @code{Gweights} returns an integer vector, whose weighting should be used to redefine the G-algebra in order
60to get the same non-commutative structure w.r.t. a weighted ordering. If the input is a matrix and the output is the zero
61vector then there is not a G-algebra structure associated to these relations with respect to the given variables.
62@*Another possibility is to use @code{weightedRing} to obtain directly a G-algebra with the new appropriate (weighted) ordering.
63EXAMPLE: example Gweights; shows examples
64SEE ALSO: weightedRing
65"{
66  int novalid=0;
67  if (typeof(r)=="ring") //a ring is admissible as input
68  {
69    setring r;
70    matrix tails;
71    def l = ncRelations(r);
72    tails = l[2]; // l=C,D we need D, the tails of the relations
73  }
74  else
75  {
76    matrix tails;
77    if ( (typeof(r)=="matrix") || (typeof(r)=="intmat") )
78    {
79      if ( nrows(r)==ncols(r) ) //the input is a square matrix
80      {
81        tails = matrix(r);
82      }
83      else
84      {
85        novalid = 1;
86      }
87    }
88    else
89    {
90      novalid=1;
91    }
92  }
93  if (novalid==0)
94  {
95    intmat IM = SimplMat(tails);
96    if ( size(IM)>1 )
97    {
98      int n  = ncols(tails);
99      int m  = nrows(IM)-1;
100      int m1 = 0;
101      int m2 = m;
102      int m3 = 0;
103      ring simplexring=(real,10),(x),lp;// The simplex procedure requires a basering of this type
104      matrix M = IM;
105      list sol = simplex (M,m,n,m1,m2,m3);
106      return(weightvector(sol));
107    }
108    else
109    {
110      "Invalid input"; //usually because the input is a one variable ring
111      return();
112    }
113  }
114  else
115  {
116    "The input must be a ring or a square matrix";
117    return();
118  }
119}
120example
121{
122  "EXAMPLE:";echo=2;
123  ring r = (0,q),(a,b,c,d),lp;
124  matrix C[4][4];
125  C[1,2]=q; C[1,3]=q; C[1,4]=1; C[2,3]=1; C[2,4]=q; C[3,4]=q;
126  matrix D[4][4];
127  D[1,4]=(q-1/q)*b*c;
128  def S = nc_algebra(C,D); setring S; S;
129  Gweights(S);
130  def D=fetch(r,D);
131  Gweights(D);
132}
133
134///////////////////////////////////////////////////////////////////////////////
135
136// This procedure take a ring r, call to Gweights(r) and use the output
137// of Gweights(r) to make a change of order in r
138// The output is a new ring, equal to r but the order
139// r must be a G-algebra
140
141proc weightedRing(def r)
142"USAGE:   weightedRing(r); r a ring
143RETURN:  ring
144PURPOSE:  equip the variables of the given ring with weights such that the relations of new ring (with weighted variables) satisfies the ordering condition for G-algebras:
145e.g. \forall\;i<j\;\;lm_w(d_{ij})<_w x_i x_j.
146NOTE:    activate this ring with the \"setring\" command
147EXAMPLE: example weightedRing; shows examples
148SEE ALSO: Gweights
149"{
150  def wv=Gweights(r);
151  if (typeof(wv)=="intvec")
152  {
153    setring r;
154    int n=nvars(r);
155    // Generating an nxn-intmat order
156    intmat m[n][n];
157    m[1,1]=wv[1];
158    int i;
159    for (i=2; i<=n; i++)
160    {
161      m[1,i]=wv[i];
162      m[i,n+2-i]=1;
163    }
164    // End of generation.
165    def lr=ncRelations(r);
166    string newringstring="ring newring=("+charstr(r)+"),("+varstr(r)+"),M("+string(m)+")";
167    execute (newringstring);
168    def lnewring=imap(r,lr);
169    return( nc_algebra(lnewring[1],lnewring[2]) );
170  }
171  else
172  {
173    "Invalid input.";//usually because the input is a one variable ring
174    return();
175  }
176}
177example
178{
179  "EXAMPLE:";echo=2;
180  ring r = (0,q),(a,b,c,d),lp;
181  matrix C[4][4];
182  C[1,2]=q; C[1,3]=q; C[1,4]=1; C[2,3]=1; C[2,4]=q; C[3,4]=q;
183  matrix D[4][4];
184  D[1,4]=(q-1/q)*b*c;
185  def S = nc_algebra(C,D); setring S; S;
186  def t=weightedRing(S);
187  setring t; t;
188}
189
190///////////////////////////////////////////////////////////////////////////////
191
192// This procedure computes ei+ej-f with f running in Newton(pij) and deletes the zero rows
193
194static proc Cij(intmat M, int i,j)
195{
196  M=(-1)*M;
197  int nc=ncols(M);
198  intvec N;
199  int k;
200  for (k=1; k<=nrows(M); k++)
201  {
202    M[k,i]=M[k,i]+1;
203    M[k,j]=M[k,j]+1;
204    if (intvec(M[k,1..nc])!=0)
205    {
206      N=N,intvec(M[k,1..nc]);
207    } // we only want non-zero rows
208  }
209  if (size(N)>1)
210  {
211    N=N[2..size(N)]; // Deleting the zero added in the definition of N
212    M=intmat(N,size(N) div nc,nc); // Conversion from vector to matrix
213  }
214  else
215  {
216    intmat M[1][1]=0;
217  }
218  return (M);
219}
220
221///////////////////////////////////////////////////////////////////////////////
222
223// This procedure run over the matrix of pij calculating Cij
224
225static proc Ct(matrix P)
226{
227  int    k = ncols(P);
228  intvec T = 0;
229  int    i,j;
230//  int notails=1;
231  def S;
232  for (j=2; j<=k; j++)
233  {
234    for (i=1; i<j; i++)
235    {
236      if ( P[i,j] != 0 )
237      {
238//        notails=0;
239        S = newtonDiag(P[i,j]);
240        S = Cij(S,i,j);
241        if ( size(S)>1 )
242        {
243          T = T,S;
244        }
245      }
246    }
247  }
248  if ( size(T)==1 )
249  {
250    intmat C[1][1] = 0;
251  }
252  else
253  {
254    T=T[2..size(T)]; // Deleting the zero added in the definition of T
255    intmat C = intmat(T,size(T) div k,k); // Conversion from vector to matrix
256  }
257  return (C);
258}
259
260///////////////////////////////////////////////////////////////////////////////
261
262// The purpose of this procedure is to produce the input matrix required by simplex procedure
263
264static proc SimplMat(matrix P)
265{
266  intmat C=Ct(P);
267  if (size(C)>1)
268  {
269    int r = nrows(C);
270    int n = ncols(C);
271    int f = 1+n+r;
272    intmat M[f][n+1]=0;
273    int i;
274    for (i=2; i<=(n+1); i++)
275    {
276      M[1,i]=-1; // (0,-1,-1,-1,...) objective function in the first row
277    }
278    for (i=2; i<=f; i++) {M[i,1]=1;} // All the independent terms are 1
279    for (i=2; i<=(n+1); i++) {M[i,i]=-1;} // wi>=1 is an identity matrix
280    M[(n+2)..f,2..(n+1)]=(-1)*intvec(C); // <wi,a> >= 1, a in C ...
281  }
282  else
283  {
284    int n = ncols(P);
285    int f = 1+n;
286    intmat M[f][n+1]=0;
287    int i;
288    for (i=2; i<=(n+1); i++) {M[1,i]=-1;} // (0,-1,-1,-1,...) objective function in the first row
289    for (i=2; i<=f; i++) {M[i,1]=1;} // All the independent terms are 1
290    for (i=2; i<=(n+1); i++) {M[i,i]=-1;} // wi>=1 is an identity matrix
291  }
292  return (M);
293}
294
295///////////////////////////////////////////////////////////////////////////////
296
297// This procedure generates a nice output of the simplex method consisting of a vector
298// with the solutions. The vector is ordered.
299
300static proc weightvector(list l)
301"ASSUME:  l is the output of simplex.
302RETURN: if there is a solution, an intvec with it will be returned"
303{
304  matrix m=l[1];
305  intvec nv=l[3];
306  int sol=l[2];
307  int rows=nrows(m);
308  int N=l[6];
309  intmat wv[1][N]=0;
310  int i;
311  if (sol)
312  {
313    "no solution satisfies the given constraints";
314  }
315  else
316  {
317    for ( i = 2; i <= rows; i++ )
318    {
319      if ( nv[i-1] <= N )
320      {
321        wv[1,nv[i-1]]=int(m[i,1]);
322      }
323    }
324  }
325  return (intvec(wv));
326}
327
328
329
330///////////////////////////////////////////////////////////////////////////////
331
332// This procedure recover the non-conmutative relations (matrices C and D)
333
334proc ncRelations(def r)
335"USAGE:   ncRelations(r); r a ring
336RETURN:  list L with two elements, both elements are of type matrix:
337@*         L[1] = matrix of coefficients C,
338@*         L[2] = matrix of polynomials D
339PURPOSE: recover the noncommutative relations via matrices C and D from
340a noncommutative ring
341SEE ALSO: ringlist, G-algebras
342EXAMPLE: example ncRelations; shows examples
343"{
344  list l;
345  if (typeof(r)=="ring")
346  {
347    int n=nvars(r);
348    matrix C[n][n]=0;
349    matrix D[n][n]=0;
350    poly f; poly g;
351    if (n>1)
352    {
353      int i,j;
354      for (i=2; i<=n; i++)
355      {
356        for (j=1; j<i; j++)
357        {
358          f=var(i)*var(j); // yx=c*xy+...
359          g=var(j)*var(i); // xy
360          while (C[j,i]==0)
361          {
362            if (leadmonom(f)==leadmonom(g))
363            {
364              C[j,i]=leadcoef(f);
365              D[j,i]=D[j,i]+f-lead(f);
366            }
367            else
368            {
369              D[j,i]=D[j,i]+lead(f);
370              f=f-lead(f);
371            }
372          }
373        }
374      }
375      l=C,D;
376    }
377    else { "The ring must have two or more variables"; }
378  }
379  else { "The input must be of a type ring";}
380  return (l);
381}
382example
383{
384  "EXAMPLE:";echo=2;
385  ring r = 0,(x,y,z),dp;
386  matrix C[3][3]=0,1,2,0,0,-1,0,0,0;
387  print(C);
388  matrix D[3][3]=0,1,2y,0,0,-2x+y+1;
389  print(D);
390  def S=nc_algebra(C,D);setring S; S;
391  def l=ncRelations(S);
392  print (l[1]);
393  print (l[2]);
394}
395
396///////////////////////////////////////////////////////////////////////////////
397
398proc findimAlgebra(matrix M, list #)
399"USAGE:   findimAlgebra(M,[r]); M a matrix, r an optional ring
400RETURN:  ring
401PURPOSE: define a finite dimensional algebra structure on a ring
402NOTE:  the matrix M is used to define the relations x(i)*x(j) = M[i,j] in the
403basering (by default) or in the optional ring r.
404@* The procedure equips the ring with the noncommutative structure.
405@* The procedure exports the ideal (not a two-sided Groebner basis!), called @code{fdQuot}, for further qring definition.
406THEORY: finite dimensional algebra can be represented as a factor algebra
407of a G-algebra modulo certain two-sided ideal. The relations of a f.d. algebra are thus naturally divided into two groups: firstly, the relations
408on the variables of the ring, making it into G-algebra and the rest of them, which constitute the ideal which will be factored out.
409EXAMPLE: example findimAlgebra; shows examples
410"
411{
412  if (size(#) >0)
413  {
414    if ( typeof(#[1])!="ring" ) { return();}
415    else
416    {
417      def @R1 = #[1];
418      setring @R1;
419    }
420  }
421  int i,j;
422  int n=nvars(basering);
423  poly p;
424  ideal I;
425  number c;
426  matrix C[n][n];
427  matrix D[n][n];
428  for (i=1; i<=n; i++)
429  {
430    for (j=i; j<=n; j++)
431    {
432      p=var(i)*var(j)-M[i,j];
433      if ( (ncols(I)==1) && (I[1]==0) )   { I=p; }
434      else { I=I,p; }
435      if (j>i)
436      {
437        if ((M[i,j]!=0) && (M[j,i]!=0))
438        {
439          c = leadcoef(M[j,i])/leadcoef(M[i,j]);
440        }
441        else
442        {
443          c = 1;
444        }
445        C[i,j]=c;
446        D[i,j]= M[j,i] -c*M[i,j];
447      }
448    }
449  }
450  def save = basering;
451  def S = nc_algebra(C,D); setring S;
452  ideal fdQuot = fetch(save,I);
453  export fdQuot;
454  return(S);
455}
456example
457{
458  "EXAMPLE:";echo=2;
459  ring r=(0,a,b),(x(1..3)),dp;
460  matrix S[3][3];
461  S[2,3]=a*x(1); S[3,2]=-b*x(1);
462  def A=findimAlgebra(S); setring A;
463  fdQuot = twostd(fdQuot);
464  qring Qr = fdQuot;
465  Qr;
466}
467
468///////////////////////////////////////////////////////////////////////////////
469
470proc isCentral(poly p, list #)
471"USAGE:   isCentral(p); p poly
472RETURN:  int, 1 if p commutes with all variables and 0 otherwise
473PURPOSE: check whether p is central in a basering (that is, commutes with every generator of the ring)
474NOTE: if @code{printlevel} > 0, the procedure displays intermediate information (by default, @code{printlevel}=0 )
475EXAMPLE: example isCentral; shows examples
476"{
477  //v an integer (with v!=0, procedure will be verbose)
478  int N = nvars(basering);
479  int in;
480  int flag = 1;
481  poly   q = 0;
482  for (in=1; in<=N; in++)
483  {
484    q = p*var(in)-var(in)*p;
485    if (q!=0)
486    {
487      if ( (size(#) >0 ) || (printlevel>0) )
488      {
489        "Non-central at:", var(in);
490      }
491      flag = 0;
492    }
493  }
494  return(flag);
495}
496example
497{
498  "EXAMPLE:";echo=2;
499  ring r=0,(x,y,z),dp;
500  matrix D[3][3]=0;
501  D[1,2]=-z;
502  D[1,3]=2*x;
503  D[2,3]=-2*y;
504  def S = nc_algebra(1,D); setring S;
505  S; // this is U(sl_2)
506  poly c = 4*x*y+z^2-2*z;
507  printlevel = 0;
508  isCentral(c);
509  poly h = x*c;
510  printlevel = 1;
511  isCentral(h);
512}
513
514///////////////////////////////////////////////////////////////////////////////
515
516proc UpOneMatrix(int N)
517"USAGE:   UpOneMatrix(n); n an integer
518RETURN:  intmat
519PURPOSE: compute an  n x n matrix with 1's in the whole upper triangle
520NOTE: helpful for setting noncommutative algebras with complicated
521coefficient matrices
522EXAMPLE: example UpOneMatrix; shows examples
523"{
524  int ii,jj;
525  intmat U[N][N]=0;
526  for (ii=1;ii<N;ii++)
527  {
528    for (jj=ii+1;jj<=N;jj++)
529    {
530      U[ii,jj]=1;
531    }
532  }
533  return(U);
534}
535example
536{
537  "EXAMPLE:";echo=2;
538  ring   r = (0,q),(x,y,z),dp;
539  matrix C = UpOneMatrix(3);
540  C[1,3]   = q;
541  print(C);
542  def S = nc_algebra(C,0); setring S;
543  S;
544}
545
546///////////////////////////////////////////////////////////////////////////////
547proc ndcond(list #)
548"USAGE:   ndcond();
549RETURN:  ideal
550PURPOSE: compute the non-degeneracy conditions of the basering
551NOTE: if @code{printlevel} > 0, the procedure displays intermediate information (by default, @code{printlevel}=0 )
552EXAMPLE: example ndcond; shows examples
553"
554{
555  // internal documentation, for tests etc
556  // 1st arg: v an optional integer (if v!=0, will be verbose)
557  // if the second argument is given, produces ndc w.r.t. powers x^N
558  int N = 1;
559  int Verbose = 0;
560  if ( size(#)>=1 ) { Verbose = int(#[1]); }
561  if ( size(#)>=2 ) { N = int(#[2]); }
562  Verbose = ((Verbose) || (printlevel>0));
563  int cnt = 1;
564  int numvars = nvars(basering);
565  int a,b,c;
566  poly p = 1;
567  ideal res = 0;
568  for (cnt=1; cnt<=N; cnt++)
569  {
570    if (Verbose) { "Processing degree :",cnt;}
571    for (a=1; a<=numvars-2; a++)
572    {
573      for (b=a+1; b<=numvars-1; b++)
574      {
575        for(c=b+1; c<=numvars; c++)
576        {
577          p = (var(c)^cnt)*(var(b)^cnt);
578          p = p*(var(a)^cnt);
579          p = p-(var(c)^cnt)*((var(b)^cnt)*(var(a)^cnt));
580          if (Verbose) {a,".",b,".",c,".";}
581          if (p!=0)
582          {
583            if ( res==0 )
584            {
585              res[1] = p;
586            }
587            else
588            {
589              res = res,p;
590            }
591            if (Verbose) { "failed:",p; }
592          }
593        }
594      }
595    }
596    if (Verbose) { "done"; }
597  }
598  return(res);
599}
600example
601{
602  "EXAMPLE:";echo=2;
603  ring r = (0,q1,q2),(x,y,z),dp;
604  matrix C[3][3];
605  C[1,2]=q2; C[1,3]=q1; C[2,3]=1;
606  matrix D[3][3];
607  D[1,2]=x; D[1,3]=z;
608  def S = nc_algebra(C,D); setring S;
609  S;
610  ideal j=ndcond(); // the silent version
611  j;
612  printlevel=1;
613  ideal i=ndcond(); // the verbose version
614  i;
615}
616
617
618///////////////////////////////////////////////////////////////////////////////
619proc Weyl(list #)
620"USAGE:   Weyl()
621RETURN:  ring
622PURPOSE: create a Weyl algebra structure on the basering
623NOTE: Activate this ring using the command @code{setring}.
624@*Assume the number of variables of a basering is 2k.
625(if the number of variables is odd, an error message will be returned)
626@*    by default, the procedure treats first k variables as coordinates x_i and the last k as differentials d_i
627@*    if a non-zero optional argument is given, the procedure treats 2k variables of a basering as k pairs (x_i,d_i), i.e. variables with odd numbers are treated as coordinates and with even numbers as differentials
628SEE ALSO: makeWeyl
629EXAMPLE: example Weyl; shows examples
630"
631{
632  //there are two possibilities for choosing the PBW basis.
633  //The variables have names x(i) for coordinates and d(i) for partial
634  // differentiations. By default, the procedure
635  //creates a ring, where the variables are ordered as x(1..n),d(1..n).  the
636  // tensor product-like realization x(1),d(1),x(2),d(2),... is used.
637  string rname=nameof(basering);
638  if ( rname == "basering") // i.e. no ring has been set yet
639  {
640    "You have to call the procedure from the ring";
641    return();
642  }
643  int @chr = 0;
644  if ( size(#) > 0 )
645  {
646    if ( typeof( #[1] ) == "int" )
647    {
648      @chr = #[1];
649    }
650  }
651  int nv = nvars(basering);
652  int N = nv div 2;
653  if ((nv % 2) != 0)
654  {
655    "Cannot create Weyl structure for an odd number of generators";
656    return();
657  }
658  matrix @D[nv][nv];
659  int i;
660  for ( i=1; i<=N; i++ )
661  {
662    if ( @chr==0 ) // default
663    {
664      @D[i,N+i]=1;
665    }
666    else
667    {
668      @D[2*i-1,2*i]=1;
669    }
670  }
671  def @R = nc_algebra(1,@D);
672  return(@R);
673}
674example
675{
676  "EXAMPLE:";echo=2;
677  ring A1=0,(x(1..2),d(1..2)),dp;
678  def S=Weyl();
679  setring S;  S;
680  kill A1,S;
681  ring B1=0,(x1,d1,x2,d2),dp;
682  def S=Weyl(1);
683  setring S;  S;
684}
685
686///////////////////////////////////////////////////////////////////////////////
687proc makeHeisenberg(int N, list #)
688"USAGE:  makeHeisenberg(n, [p,d]); int n (setting 2n+1 variables), optional int p (field characteristic), optional int d (power of h in the commutator)
689RETURN: ring
690PURPOSE: create the n-th Heisenberg algebra in the variables x(1),y(1),...,x(n),y(n),h over the rationals Q or F_p with the relations
691\forall\;i\in\{1,2,\ldots,n\}\;\;y(j)x(i) = x(i)y(j)+h^d.
692SEE ALSO: makeWeyl
693NOTE: activate this ring with the @code{setring} command
694@*       If p is not prime, the next larger prime number will be used.
695EXAMPLE: example makeHeisenberg; shows examples
696"
697{
698  int @chr = 0;
699  int @deg = 1;
700  if ( size(#) > 0 )
701  {
702    if ( typeof( #[1] ) == "int" )
703    {
704      @chr = #[1];
705    }
706  }
707  if ( size(#) > 1 )
708  {
709    if ( typeof( #[2] ) == "int" )
710    {
711      @deg = #[2];
712      if (@deg <1) { @deg = 1; }
713    }
714  }
715  ring @@r=@chr,(x(1..N),y(1..N),h),lp;
716  matrix D[2*N+1][2*N+1];
717  int i;
718  for (i=1;i<=N;i++)
719  {
720    D[i,N+i]=h^@deg;
721  }
722  return(nc_algebra(1,D));
723}
724example
725{
726  "EXAMPLE:";echo=2;
727  def a = makeHeisenberg(2);
728  setring a;   a;
729  def H3 = makeHeisenberg(3, 7, 2);
730  setring H3;  H3;
731}
732
733
734///////////////////////////////////////////////////////////////////////////////
735proc superCommutative(list #)
736"USAGE:   superCommutative([b,[e, [Q]]]);
737RETURN:  qring
738PURPOSE:  create a super-commutative algebra (as a GR-algebra) over a basering,
739NOTE: activate this qring with the \"setring\" command.
740NOTE: if b==e then the resulting ring is commutative.
741@* By default, @code{b=1, e=nvars(basering), Q=0}.
742THEORY: given a basering, this procedure introduces the anti-commutative relations
743@* var(j)var(i)=-var(i)var(j) for all e>=j>i>=b and creates the quotient
744@* of the anti-commutative algebra modulo the two-sided ideal, generated by
745@* x(b)^2, ..., x(e)^2[ + Q]
746DISPLAY: If @code{printlevel} > 1, warning debug messages will be printed
747EXAMPLE: example superCommutative; shows examples
748"
749{
750  int fprot = (printlevel > 1); // (find(option(),"prot") != 0);
751
752  string rname=nameof(basering);
753
754  if ( rname == "basering") // i.e. no ring has been set yet
755  {
756    ERROR("You have to call the procedure from the ring");
757    return();
758  }
759
760  def saveRing = basering;
761
762  int N = nvars(saveRing);
763  int b = 1;
764  int e = N;
765  int flag = 0;
766
767  ideal Q = 0;
768
769  if(size(#)>0)
770  {
771    if(typeof(#[1]) != "int")
772    {
773      ERROR("The argument 'b' must be an integer!");
774      return();
775    }
776    b = #[1];
777
778    if((b < 1)||(b > N))
779    {
780      ERROR("The argument 'b' must within [1..nvars(basering)]!");
781      return();
782    }
783
784  }
785
786  if(size(#)>1)
787  {
788    if(typeof(#[2]) != "int")
789    {
790      ERROR("The argument 'e' must be an integer!");
791      return();
792    }
793    e = #[2];
794
795    if((e < 1)||(e > N))
796    {
797      ERROR("The argument 'e' must within [1..nvars(basering)]!");
798      return();
799    }
800
801    if(e < b)
802    {
803      ERROR("The argument 'e' must be bigger or equal to 'b'!");
804      return();
805    }
806  }
807
808  if(size(#)>2)
809  {
810    if(typeof(#[3]) != "ideal")
811    {
812      ERROR("The argument 'Q' must be an ideal!");
813      return();
814    }
815    Q = #[3];
816  }
817
818/*  if(size(#)>3)
819  {
820    if(typeof(#[4]) != "int")
821    {
822      ERROR("The argument 'flag' must be an integer!");
823      return();
824    }
825    flag = #[4];
826  }
827*/
828
829  int iSavedDegBoung = degBound;
830
831  if( (b == e) && (flag == 0) ) // commutative ring!!!
832  {
833    if( fprot == 1)
834    {
835      print("Warning: (b==e) means that the resulting ring will be commutative!");
836    }
837
838    degBound=0;
839    Q = std(Q + (var(b)^2));
840    degBound = iSavedDegBoung;
841
842    qring @EA = Q; // and it will be internally commutative as well!!!
843
844    return(@EA);
845  }
846
847/*
848  // Singular'(H.S.) politics: no ring copies!
849  // in future nc_algebra() should return a new ring!!!
850  list CurrRing = ringlist(basering);
851  def @R = ring(CurrRing);
852  setring @R; // @R;
853*/
854  int i, j;
855
856  if( (char(basering)==2) && (flag == 0) )// commutative ring!!!
857  {
858    if( fprot == 1)
859    {
860      print("Warning: (char == 2) means that the resulting ring will be commutative!");
861    }
862
863    ideal I;
864
865    for (i = e - b + 1; i > 0; i--)
866    {
867      I[i] = var(i + b - 1)^2;
868    }
869
870    degBound=0;
871    Q = std(I + Q);
872    degBound = iSavedDegBoung;
873
874    qring @EA = Q; // and it will be internally commutative as well!!!
875    return(@EA);
876  }
877
878
879
880  if( (b == 1) && (e == N) ) // just an exterior algebra?
881  {
882    def S = nc_algebra(-1, 0); // define ground G-algebra!
883    setring S;
884  } else
885  {
886    matrix @E = UpOneMatrix(N);
887
888    for ( i = b; i < e; i++ )
889    {
890      for ( j = i+1; j <= e; j++ )
891      {
892        @E[i, j] = -1;
893      }
894    }
895    def S = nc_algebra(@E, 0); // define ground G-algebra!
896    setring S;
897  }
898
899  ideal @I;
900
901  for (i = e - b + 1; i > 0; i--)
902  {
903    @I[i] = var(i + b - 1)^2;
904  }
905
906
907  degBound=0;
908  @I = twostd(@I); // must be computed within the ground G-algebra => problems with local orderings!
909  degBound = iSavedDegBoung;
910
911  qring @EA = @I;
912
913  ideal @Q = twostd(fetch(saveRing, Q));
914
915  if( size(@Q) > 0 )
916  {
917    qring @EA2 = @Q;
918  }
919
920  attrib(basering, "isSCA", 1==1);
921  attrib(basering, "iAltVarStart", b);
922  attrib(basering, "iAltVarEnd", e);
923
924//   "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
925  return(basering);
926}
927example
928{
929  "EXAMPLE:";echo=2;
930  ring R = 0,(x(1..4)),dp; // global!
931  def ER = superCommutative(); // the same as Exterior (b = 1, e = N)
932  setring ER; ER;
933  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
934  kill R; kill ER;
935  ring R = 0,(x(1..4)),(lp(1), dp(3)); // global!
936  def ER = superCommutative(2); // b = 2, e = N
937  setring ER; ER;
938  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
939  kill R; kill ER;
940  ring R = 0,(x, y, z),(ds(1), dp(2)); // mixed!
941  def ER = superCommutative(2,3); // b = 2, e = 3
942  setring ER; ER;
943  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
944  x + 1 + z + y; // ordering on variables: y > z > 1 > x
945  std(x - x*x*x);
946  std(ideal(x - x*x*x, x*x*z + y, z + y*x*x));
947  kill R; kill ER;
948  ring R = 0,(x, y, z),(ds(1), dp(2)); // mixed!
949  def ER = superCommutative(2, 3, ideal(x - x*x, x*x*z + y, z + y*x*x )); // b = 2, e = 3
950  setring ER; ER;
951  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
952}
953
954// Please, don't throw this away!!! Needed for backward compatibility.
955proc SuperCommutative(list #)
956"USAGE:   please use @code{superCommutative} instead
957"
958{
959  "// This procedure is deprecated. Please use superCommutative instead";
960  return( superCommutative(#) );
961}
962example
963{
964  "EXAMPLE:";
965  "Procedure is deprecated. Please use superCommutative instead";
966}
967
968static proc ParseSCA()
969"
970RETURN: list {AltVarStart, AltVarEnd} is currRing is SCA, returns undef otherwise.
971NOTE: rings with only one non-commutative variable are commutative rings which are super-sommutative itself!
972"
973{
974  if(typeof(attrib(basering, "isSCA"))=="int") // workaround, if(defined()) doesn't work!!!!
975  {
976    if(typeof(attrib(basering, "iAltVarStart"))=="int")
977    {
978      if(typeof(attrib(basering, "iAltVarEnd"))=="int")
979      {
980        if(attrib(basering, "isSCA"))
981        {
982          return(list(
983            attrib(basering, "iAltVarStart"),
984            attrib(basering, "iAltVarEnd")
985                 ));
986        }
987      }
988    }
989  }
990
991  def saveRing = basering;
992
993  int i, j;
994  int N = nvars(saveRing);
995
996  int b = N+1;
997  int e =  -1;
998
999  int fprot = 0; // (find(option(),"prot") != 0);
1000
1001
1002  if( size(ideal(saveRing)) == 0 )
1003  {
1004    return("SCA rings are factors by (at least) squares!"); // no squares in the factor ideal!
1005  }
1006
1007  list L = ringlist(saveRing);
1008
1009  if( size(L)!=6 )
1010  {
1011    if(fprot)
1012    {
1013      print("// Warning: The current ring is internally commutative!");
1014    }
1015
1016    for( i = N; i > 0; i-- )
1017    {
1018      if( NF(var(i)^2, std(0)) == 0 )
1019      {
1020        if( (fprot == 1) and (i > 1) )
1021        {
1022          print("// Warning: the SCA representation of the current commutative factor ring may be ambiguous!");
1023        }
1024
1025        return( list(i, i) ); // this is not unique in this case! there may be other squares in the factor ideal!
1026      }
1027    }
1028
1029    return("The current commutative ring is not SCA! (Wrong quotient ideal)"); // no squares in the factor ideal!
1030  }
1031
1032  module D = simplify(L[6], 2 + 4);
1033
1034  if( size(D)>0 )
1035  {
1036    return("The current ring is not SCA! (D!=0)");
1037  }
1038
1039  matrix C = L[5];
1040  poly c;
1041
1042  for( i = 1; i < N; i++ )
1043  {
1044    for( j = i+1; j <= N; j++ )
1045    {
1046      c = C[i, j];
1047
1048      if( c == -1 )
1049      {
1050        if(i < b)
1051        {
1052          b = i;
1053        }
1054
1055        if(j > e)
1056        {
1057          e = j;
1058        }
1059      } else
1060      { // should commute
1061        if( c!=1 )
1062        {
1063          return("The current ring is not SCA! (C["+ string(i)+"," + string(j)+"]!=1)");
1064        }
1065      }
1066    }
1067  }
1068
1069  if( (b > N) || (e < 1))
1070  {
1071    if(fprot)
1072    {
1073      print("Warning: The current ring is a commutative GR-algebra!");
1074    }
1075
1076    for( i = N; i > 0; i-- )
1077    {
1078      if( NF(var(i)^2, std(0)) == 0 )
1079      {
1080        if( (fprot == 1) and (i > 1) )
1081        {
1082          print("Warning: the SCA representation of the current factor ring may be ambiguous!");
1083        }
1084
1085        return( list(i, i) ); // this is not unique in this case! there may be other squares in the factor ideal!
1086      }
1087    }
1088
1089    return("The current commutative GR-algebra is not SCA! (Wrong quotient ideal)"); // no squares in the factor ideal!
1090  }
1091
1092  for( i = 1; i < N; i++ )
1093  {
1094    for( j = i+1; j <= N; j++ )
1095    {
1096      c = C[i, j];
1097
1098      if( (b <= i) && (j <= e) ) // S <= i < j <= E
1099      { // anticommutative part
1100        if( c!= -1 )
1101        {
1102          return("The current ring is not SCA! (C["+ string(i)+"," + string(j)+"]!=-1)");
1103        }
1104      } else
1105      { // should commute
1106        if( c!=1 )
1107        {
1108          return("The current ring is not SCA! (C["+ string(i)+"," + string(j)+"]!=1)");
1109        }
1110      }
1111    }
1112  }
1113
1114  for( i = b; i <= e; i++ )
1115  {
1116    if( NF(var(i)^2, std(0)) != 0 )
1117    {
1118      return("The current ring is not SCA! (Wrong quotient ideal)");
1119    }
1120  }
1121
1122  ////////////////////////////////////////////////////////////////////////
1123  // ok. this is a SCA!!!
1124
1125  return(list(b, e));
1126}
1127
1128///////////////////////////////////////////////////////////////////////////////
1129proc AltVarStart()
1130"USAGE:   AltVarStart();
1131RETURN:  int
1132PURPOSE:  returns the number of the first alternating variable of basering
1133NOTE:  basering should be a super-commutative algebra constructed by
1134@*     the procedure @code{superCommutative}, emits an error otherwise
1135EXAMPLE: example AltVarStart; shows examples
1136"
1137{
1138  def l = ParseSCA();
1139
1140  if( typeof(l) != "string" )
1141  {
1142    return(l[1]);
1143  }
1144
1145  ERROR(l);
1146  return();
1147}
1148example
1149{
1150  "EXAMPLE:";echo=2;
1151  ring R = 0,(x(1..4)),dp; // global!
1152  def ER = superCommutative(2); // (b = 2, e = N)
1153  setring ER; ER;
1154  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
1155  setring R;
1156  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
1157  kill R, ER;
1158  //////////////////////////////////////////////////////////////////
1159  ring R = 2,(x(1..4)),dp; // the same in char. = 2!
1160  def ER = superCommutative(2); // (b = 2, e = N)
1161  setring ER; ER;
1162  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
1163  setring R;
1164  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
1165}
1166
1167///////////////////////////////////////////////////////////////////////////////
1168proc AltVarEnd()
1169"USAGE:   AltVarStart();
1170RETURN:  int
1171PURPOSE:  returns the number of the last alternating variable of basering
1172NOTE:  basering should be a super-commutative algebra constructed by
1173@*     the procedure @code{superCommutative}, emits an error otherwise
1174EXAMPLE: example AltVarEnd; shows examples
1175"
1176{
1177  def l = ParseSCA();
1178
1179  if( typeof(l) != "string" )
1180  {
1181    return(l[2]);
1182  }
1183
1184  ERROR(l);
1185  return();
1186}
1187example
1188{
1189  "EXAMPLE:";echo=2;
1190  ring R = 0,(x(1..4)),dp; // global!
1191  def ER = superCommutative(2); // (b = 2, e = N)
1192  setring ER; ER;
1193  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
1194  setring R;
1195  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
1196  kill R, ER;
1197  //////////////////////////////////////////////////////////////////
1198  ring R = 2,(x(1..4)),dp; // the same in char. = 2!
1199  def ER = superCommutative(2); // (b = 2, e = N)
1200  setring ER; ER;
1201  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
1202  setring R;
1203  "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "].";
1204}
1205
1206///////////////////////////////////////////////////////////////////////////////
1207proc IsSCA()
1208"USAGE:   IsSCA();
1209RETURN:  int
1210PURPOSE:  returns 1 if basering is a super-commutative algebra and 0 otherwise
1211EXAMPLE: example IsSCA; shows examples
1212"
1213{
1214  def l = ParseSCA();
1215
1216  if( typeof(l) != "string" )
1217  {
1218    return(1);
1219  }
1220
1221  if( find(option(),"prot") != 0 )
1222  {
1223    print(l);
1224  }
1225
1226  return(0);
1227}
1228example
1229{
1230  "EXAMPLE:";echo=2;
1231/////////////////////////////////////////////////////////////////////
1232  ring R = 0,(x(1..4)),dp; // commutative
1233  if(IsSCA())
1234    { "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "]."; }
1235  else
1236    { "Not a super-commutative algebra!!!"; }
1237  kill R;
1238/////////////////////////////////////////////////////////////////////
1239  ring R = 0,(x(1..4)),dp;
1240  def S = nc_algebra(1, 0); setring S; S; // still commutative!
1241  if(IsSCA())
1242    { "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "]."; }
1243  else
1244    { "Not a super-commutative algebra!!!"; }
1245  kill R, S;
1246/////////////////////////////////////////////////////////////////////
1247  ring R = 0,(x(1..4)),dp;
1248  list CurrRing = ringlist(R);
1249  def ER = ring(CurrRing);
1250  setring ER; // R;
1251
1252  matrix E = UpOneMatrix(nvars(R));
1253
1254  int i, j; int b = 2; int e = 3;
1255
1256  for ( i = b; i < e; i++ )
1257  {
1258    for ( j = i+1; j <= e; j++ )
1259    {
1260      E[i, j] = -1;
1261    }
1262  }
1263
1264  def S = nc_algebra(E,0); setring S; S;
1265
1266  if(IsSCA())
1267    { "Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "]."; }
1268  else
1269    { "Not a super-commutative algebra!!!"; }
1270  kill R, ER, S;
1271/////////////////////////////////////////////////////////////////////
1272  ring R = 0,(x(1..4)),dp;
1273  def ER = superCommutative(2); // (b = 2, e = N)
1274  setring ER; ER;
1275  if(IsSCA())
1276    { "This is a SCA! Alternating variables: [", AltVarStart(), ",", AltVarEnd(), "]."; }
1277  else
1278    { "Not a super-commutative algebra!!!"; }
1279  kill R, ER;
1280}
1281
1282
1283
1284///////////////////////////////////////////////////////////////////////////////
1285proc Exterior(list #)
1286"USAGE:   Exterior();
1287RETURN:  qring
1288PURPOSE:  create the exterior algebra of a basering
1289NOTE:  activate this qring with the \"setring\" command
1290THEORY: given a basering, this procedure introduces the anticommutative relations x(j)x(i)=-x(i)x(j) for all j>i,
1291@* moreover, creates a factor algebra modulo the two-sided ideal, generated by x(i)^2 for all i
1292EXAMPLE: example Exterior; shows examples
1293"
1294{
1295  string rname=nameof(basering);
1296  if ( rname == "basering") // i.e. no ring has been set yet
1297  {
1298    "You have to call the procedure from the ring";
1299    return();
1300  }
1301  int N = nvars(basering);
1302  string NewRing = "ring @R=("+charstr(basering)+"),("+varstr(basering)+"),("+ordstr(basering)+");";
1303  execute(NewRing);
1304  matrix @E = UpOneMatrix(N);
1305  @E = -1*(@E);
1306  def @@RR = nc_algebra(@E,0); setring @@RR;
1307  int i;
1308  ideal Q;
1309  for ( i=1; i<=N; i++ )
1310  {
1311    Q[i] = var(i)^2;
1312  }
1313  Q = twostd(Q);
1314  qring @EA = Q;
1315  return(@EA);
1316}
1317example
1318{
1319  "EXAMPLE:";echo=2;
1320  ring R = 0,(x(1..3)),dp;
1321  def ER = Exterior();
1322  setring ER;
1323  ER;
1324}
1325
1326///////////////////////////////////////////////////////////////////////////////
1327proc makeWeyl(int n, list #)
1328"USAGE:  makeWeyl(n,[p]); n an integer, n>0; p an optional integer (field characteristic)
1329RETURN:  ring
1330PURPOSE: create the n-th Weyl algebra over the rationals Q or F_p
1331NOTE:    activate this ring with the \"setring\" command.
1332@*       The presentation of an n-th Weyl algebra is classical: D(i)x(i)=x(i)D(i)+1,
1333@*       where x(i) correspond to coordinates and D(i) to partial differentiations, i=1,...,n.
1334@*       If p is not prime, the next larger prime number will be used.
1335SEE ALSO: Weyl
1336EXAMPLE: example makeWeyl; shows examples
1337"{
1338  if (n<1)
1339  {
1340    print("Incorrect input");
1341    return();
1342  }
1343  int @p = 0;
1344  if ( size(#) > 0 )
1345  {
1346    if ( typeof( #[1] ) == "int" )
1347    {
1348      @p = #[1];
1349    }
1350  }
1351  if (n ==1)
1352  {
1353    ring @rr = @p,(x,D),dp;
1354  }
1355  else
1356  {
1357    ring @rr = @p,(x(1..n),D(1..n)),dp;
1358  }
1359  setring @rr;
1360  def @rrr = Weyl();
1361  return(@rrr);
1362}
1363example
1364{ "EXAMPLE:"; echo = 2;
1365   def a = makeWeyl(3);
1366   setring a;
1367   a;
1368}
1369
1370//////////////////////////////////////////////////////////////////////
1371proc isNC()
1372"USAGE:   isNC();
1373PURPOSE: check whether a basering is commutative or not
1374RETURN:   int, 1 if basering is noncommutative and 0 otherwise
1375EXAMPLE: example isNC; shows examples
1376"{
1377  string rname=nameof(basering);
1378  if ( rname == "basering") // i.e. no ring has been set yet
1379  {
1380    "You have to call the procedure from the ring";
1381    return();
1382  }
1383  int n = nvars(basering);
1384  int i,j;
1385  poly p;
1386  for (i=1; i<n; i++)
1387  {
1388    for (j=i+1; j<=n; j++)
1389    {
1390      p = var(j)*var(i) - var(i)*var(j);
1391      if (p!=0) { return(1);}
1392    }
1393  }
1394  return(0);
1395}
1396example
1397{ "EXAMPLE:"; echo = 2;
1398   def a = makeWeyl(2);
1399   setring a;
1400   isNC();
1401   kill a;
1402   ring r = 17,(x(1..7)),dp;
1403   isNC();
1404   kill r;
1405}
1406
1407///////////////////////////////////////////////////////////////////////////////
1408proc rightStd(def I)
1409"USAGE:  rightStd(I); I an ideal/ module
1410PURPOSE: compute a right Groebner basis of I
1411RETURN:  the same type as input
1412EXAMPLE: example rightStd; shows examples
1413"
1414{
1415  def A = basering;
1416  def Aopp = opposite(A);
1417  setring Aopp;
1418  def Iopp = oppose(A,I);
1419  def Jopp = groebner(Iopp);
1420  setring A;
1421  def J = oppose(Aopp,Jopp);
1422  return(J);
1423}
1424example
1425{ "EXAMPLE:"; echo = 2;
1426  LIB "ncalg.lib";
1427  def A = makeUsl(2);
1428  setring A;
1429  ideal I = e2,f;
1430  option(redSB);
1431  option(redTail);
1432  ideal LI = std(I);
1433  LI;
1434  ideal RI = rightStd(I);
1435  RI;
1436}
1437
1438///////////////////////////////////////////////////////////////////////////////
1439proc rightSyz(def I)
1440"USAGE:  rightSyz(I); I an ideal/ module
1441PURPOSE: compute a right syzygy module of I
1442RETURN:  the same type as input
1443EXAMPLE: example rightSyz; shows examples
1444"
1445{
1446  def A = basering;
1447  def Aopp = opposite(A);
1448  setring Aopp;
1449  def Iopp = oppose(A,I);
1450  def Jopp = syz(Iopp);
1451  setring A;
1452  def J = oppose(Aopp,Jopp);
1453  return(J);
1454}
1455example
1456{ "EXAMPLE:"; echo = 2;
1457  ring r = 0,(x,d),dp;
1458  def S = nc_algebra(1,1); setring S; // the first Weyl algebra
1459  ideal I = x,d;
1460  module LS = syz(I);
1461  print(LS);
1462  module RS = rightSyz(I);
1463  print(RS);
1464}
1465
1466///////////////////////////////////////////////////////////////////////////////
1467proc rightNF(def v, def M)
1468"USAGE:  rightNF(I); v a poly/vector, M an ideal/module
1469PURPOSE: compute a right normal form of v w.r.t. M
1470RETURN:  poly/vector (as of the 1st argument)
1471EXAMPLE: example rightNF; shows examples
1472"
1473{
1474  def A = basering;
1475  def Aopp = opposite(A);
1476  setring Aopp;
1477  def vopp = oppose(A,v);
1478  def Mopp = oppose(A,M);
1479  Mopp = std(Mopp);
1480  def wopp = NF(vopp,Mopp);
1481  setring A;
1482  def w    = oppose(Aopp,wopp);
1483  w = simplify(w,2); // skip zeros in ideal/module
1484  return(w);
1485}
1486example
1487{ "EXAMPLE:"; echo = 2;
1488  LIB "ncalg.lib";
1489  ring r = 0,(x,d),dp;
1490  def S = nc_algebra(1,1); setring S; // Weyl algebra
1491  ideal I = x; I = std(I);
1492  poly  p = x*d+1;
1493  NF(p,I); // left normal form
1494  rightNF(p,I); // right normal form
1495}
1496
1497// **********************************
1498// * NF: Example for vector/module: *
1499// **********************************
1500// module M = [x,0],[0,d]; M = std(M);
1501// vector v = (x*d+1)*[1,1];
1502// print(NF(v,M));
1503// print(rightNF(v,M));
1504
1505///////////////////////////////////////////////////////////////////////////////
1506proc rightModulo(def M, def N)
1507"USAGE:  rightModulo(M,N); M,N are ideals/modules
1508PURPOSE: compute a right representation of the module (M+N)/N
1509RETURN:  module
1510ASSUME:  M,N are presentation matrices for right modules
1511EXAMPLE: example rightModulo; shows examples
1512"
1513{
1514  def A = basering;
1515  def Aopp = opposite(A);
1516  setring Aopp;
1517  def Mopp = oppose(A,M);
1518  def Nopp = oppose(A,N);
1519  def Kopp = modulo(Mopp,Nopp);
1520  setring A;
1521  def K = oppose(Aopp,Kopp);
1522  return(K);
1523}
1524example
1525{ "EXAMPLE:"; echo = 2;
1526  LIB "ncalg.lib";
1527  def A = makeUsl(2);
1528  setring A;
1529  option(redSB);
1530  option(redTail);
1531  ideal I = e2,f2,h2-1;
1532  I = twostd(I);
1533  print(matrix(I));
1534  ideal E  = std(e);
1535  ideal TL = e,h-1; // the result of left modulo
1536  TL;
1537  ideal T = rightModulo(E,I);
1538  T = rightStd(T+I);
1539  T = rightStd(rightNF(T,I)); // make the output canonic
1540  T;
1541}
1542
1543//////////////////////////////////////////////////////////////////////
1544
1545proc isCommutative ()
1546"USAGE:  isCommutative();
1547RETURN:  int, 1 if basering is commutative, or 0 otherwise
1548PURPOSE: check whether basering is commutative
1549EXAMPLE: example isCommutative; shows an example
1550"
1551{
1552  int iscom = 1;
1553  list L = ringlist(basering);
1554  if (size(L) > 4) // basering is nc_algebra
1555  {
1556    matrix C = L[5];
1557    matrix D = L[6];
1558    if (size(module(D)) <> 0) { iscom = 0; }
1559    else
1560    {
1561      matrix U = UpOneMatrix(nvars(basering));
1562      if (size(module(C-U)) <> 0) { iscom = 0; }
1563    }
1564  }
1565  return(iscom);
1566}
1567example
1568{
1569  "EXAMPLE:"; echo = 2;
1570  ring r = 0,(x,y),dp;
1571  isCommutative();
1572  def D = Weyl(); setring D;
1573  isCommutative();
1574  setring r;
1575  def R = nc_algebra(1,0); setring R;
1576  isCommutative();
1577}
1578
1579//////////////////////////////////////////////////////////////////////
1580
1581proc isWeyl ()
1582"USAGE:  isWeyl();
1583RETURN:  int, 1 if basering is a Weyl algebra, or 0 otherwise
1584PURPOSE: check whether basering is a Weyl algebra
1585EXAMPLE: example isWeyl; shows an example
1586"
1587{
1588  int i,j;
1589  int notW = 0;
1590  int N = nvars(basering);
1591  if (N mod 2 <> 0) { return(notW); } // odd number of generators
1592  int n = N div 2;
1593  list L = ringlist(basering);
1594  if (size(L) < 6) { return(notW); } // basering is commutative
1595  matrix C = L[5];
1596  matrix D = L[6];
1597  matrix U = UpOneMatrix(N);
1598  if (size(ideal(C-U)) <> 0) { return(notW); } // lt(xy)<>lt(yx)
1599  ideal I = D;
1600  if (size(I) <> n) { return(notW); } // not n entries<>0
1601  I = simplify(I,4+2);
1602  int sI = size(I);
1603  if (sI > 2) { return(notW); }  // more than 2 distinct entries
1604  for (i=1; i<=sI; i++)
1605  {
1606    if (I[i]<>1 && I[i]<>-1) { return (notW); } // other values apart from 1,-1
1607  }
1608  ideal Ro,Co;
1609  for (i=1; i<=N; i++)
1610  {
1611    Ro = D[1..N,i];
1612    Co = D[i,1..N];
1613    if (size(Ro)>1 || size(Co)>1)
1614    {
1615      return(int(0)); // var(i) doesn't commute with more than 1 other vars
1616    }
1617  }
1618  return(int(1)); // all tests passed: basering is Weyl algebra
1619}
1620example
1621{
1622  "EXAMPLE:"; echo = 2;
1623  ring r = 0,(a,b,c,d),dp;
1624  isWeyl();
1625  def D = Weyl(1); setring D; //make from r a Weyl algebra
1626  b*a;
1627  isWeyl();
1628  ring t = 0,(Dx,x,y,Dy),dp;
1629  matrix M[4][4]; M[1,2]=-1; M[3,4]=1;
1630  def T = nc_algebra(1,M); setring T;
1631  isWeyl();
1632}
1633
1634//////////////////////////////////////////////////////////////////////
1635
1636proc embedMat(matrix A, int m, int n)
1637"USAGE:  embedMat(A,m,n); A,B matrix/module
1638RETURN:  matrix
1639PURPOSE: embed A in the left upper corner of mxn matrix
1640EXAMPLE: example embedMat; shows an example
1641"
1642{
1643  // returns A embedded in the left upper corner of mxn matrix
1644  int rA = nrows(A);
1645  int cA = ncols(A);
1646  if ((rA >m) || (cA>n))
1647  {
1648    ERROR("wrong dimensions of the new matrix");
1649  }
1650  matrix @M[m][n];
1651  int i,j;
1652  for(i=1;i<=rA; i++)
1653  {
1654    for(j=1;j<=cA; j++)
1655    {
1656      @M[i,j]=A[i,j];
1657    }
1658  }
1659  return(@M);
1660}
1661example
1662{
1663  "EXAMPLE:"; echo = 2;
1664  ring r = 0,(a,b,c,d),dp;
1665  matrix M[2][3]; M[1,1]=a; M[1,2]=b;M[2,2]=d;M[1,3]=c;
1666  print(M);
1667  print(embedMat(M,3,4));
1668  matrix N = M; N[2,2]=0;
1669  print(embedMat(N,3,4));
1670}
1671
1672//proc moduloSlim (matrix A, matrix B)
1673proc moduloSlim (module A, module B)
1674"USAGE:  moduloSlim(A,B); A,B module/matrix/ideal
1675RETURN:  module
1676PURPOSE: compute @code{modulo} with slimgb as engine
1677EXAMPLE: example moduloSlim; shows an example
1678"
1679{
1680  def save  = basering;
1681  int rA = nrows(A);  int rB = nrows(B);
1682  int cA = ncols(A);  int cB = ncols(B);
1683  int j;
1684  int dab; // difference a,b
1685  dab = rA - rB;
1686  if (dab <0)
1687  {
1688    // rA<rB: add zero rows to A
1689    dab = -dab;
1690    A = embedMat(A,rB,cA);
1691  }
1692  else
1693  {
1694    // rA>rB: add zero rows to B
1695    B = embedMat(B,rA,cB);
1696  }
1697  def mering = makeModElimRing(save);
1698  setring mering;
1699  module A = imap(save, A);
1700  module B = imap(save, B);
1701  // create matrix C
1702  //  matrix C[2*rA][cA+cB];
1703  module C;
1704  int i;
1705  for(i=1; i<= cA; i++)
1706  {
1707    C = C, A[i] + gen(rA + i);
1708  }
1709  C = C,B;
1710//   for(i=1; i<=cB; i++)
1711//   {
1712//     C = C, B[i];
1713//   }
1714  C = C[2..ncols(C)];
1715//  print(C);
1716  matrix D = slimgb(C);
1717  module E; int k;
1718  // TODO: why only first row? need smth like rA rows...
1719  for(i=1; i<= ncols(D); i++)
1720  {
1721    k=1;
1722    // determine first zero in the column
1723    while ( (D[k,i]==0) && (k<= cA+rA) )
1724    {
1725      k++;
1726    }
1727    // what can that be: k = cA+rA+1=> zero column
1728    // k<=rA => column not in ker
1729    // rA+1 <= k <= rA+cA => column in ker
1730    if ( ( k>=rA+1) && (k<=rA+cA) )
1731    {
1732      E = E,D[i];
1733    }
1734  }
1735//   for(i=1; i<= ncols(D); i++)
1736//   {
1737//     if (D[1,i]==0)
1738//     {
1739//       E = E,D[i];
1740//     }
1741//   }
1742//  // this E has 1st column and 1st row zero
1743  // use submat@matrix.lib
1744  //  E = submat(E,intvec(2..nrows(E)),intvec(2..ncols(E)));
1745  E = submat(E,intvec(rA+1..nrows(E)),intvec(2..ncols(E)));
1746  setring save;
1747  module E = imap(mering,E);
1748  kill mering;
1749  // TODO: clean components!
1750  return(E);
1751}
1752example
1753{
1754  "EXAMPLE:"; echo = 2;
1755  LIB "ncalg.lib";
1756  ring r; // first classical example for modulo
1757  ideal h1=x,y,z;    ideal h2=x;
1758  module m=moduloSlim(h1,h2);
1759  print(m);
1760  // now, a noncommutative example
1761  def A = makeUsl2(); setring A; // this algebra is U(sl_2)
1762  ideal H2 = e2,f2,h2-1; H2 = twostd(H2);
1763  print(matrix(H2)); // print H2 in a compact form
1764  ideal H1 = std(e);
1765  ideal T = moduloSlim(H1,H2);
1766  T = std( NF(std(H2+T),H2) );
1767  T;
1768  // now, a matrix example:
1769  ring r2 = 0,(x,d), (dp);
1770  def R = nc_algebra(1,1); setring R;
1771  matrix M[2][2] = d, 0, 0, d*(x*d);
1772  matrix P[2][1] = (8x+7)*d+9x, (x2+1)*d + 5*x;
1773  module X = moduloSlim(P,M);
1774  print(X);
1775}
1776
1777//////////////////////////////////////////////////////////////////////
1778
1779proc makeModElimRing(list #)
1780"USAGE:  makeModElimRing(L); L a list
1781RETURN:  ring
1782PURPOSE: create a copy of a given ring equipped with the
1783@* elimination ordering for module components @code{(c,<)}
1784NOTE: usually the list argument contains a ring to work with
1785EXAMPLE: example makeModElimRing; shows an example
1786"
1787{
1788  // supports qring;
1789  // can be extended to handle C istead of c
1790  /* input/basering business */
1791  def save; int Noinput = 0;
1792  if ( size(#)>0 )
1793  {
1794    if ( (typeof(#[1]) == "ring" ) || (typeof(#[1]) == "qring" ) )
1795    {
1796      save = #[1];
1797    }
1798    else
1799    {
1800      print("unsupported input type, proceeding with basering");
1801      Noinput = 1;
1802    }
1803  }
1804  if (Noinput)
1805  {
1806    if (nameof(basering)=="basering")
1807    {
1808      ERROR("no rings are given");
1809    }
1810    else
1811    {
1812      save = basering;
1813    }
1814  }
1815  /* END input/basering business */
1816  list L = ringlist(save);
1817  list Ord = L[3];
1818  int s = size(Ord); int done;
1819  // detect where module ordering is located: either 1st or last entry
1820  int i,j;
1821  for(i=1; i<=s; i++)
1822  {
1823    if ( (Ord[i][1] == "C") || (Ord[i][1] == "c") )
1824    {
1825      Ord[i][1] = "c";
1826      j = i; i=s;
1827    }
1828  }
1829  if (j==0) { ERROR("no component entry found in the ringlist"); }
1830  list N;
1831  N[1] = Ord[j];
1832  for(i=2; i<=j; i++)
1833  {
1834    N[i] = Ord[i-1];
1835  }
1836  for(i=j+1; i<=s; i++)
1837  {
1838    N[i] = Ord[i];
1839  }
1840  L[3] = N; def NR = ring(L);
1841  return(NR);
1842}
1843example
1844{
1845  "EXAMPLE:"; echo = 2;
1846  ring r1 = 0,(x,y,z),(C,Dp);
1847  def r2 = makeModElimRing(r1); setring r2; r2;   kill r2;
1848  ring r3 = 0,(z,t),(wp(2,3),c);
1849  def r2 = makeModElimRing(r3); setring r2; r2; kill r2;
1850  ring r4 = 0,(z,t,u,w),(a(1,2),C,wp(2,3,4,5));
1851  def r2 = makeModElimRing(r4); setring r2; r2;
1852}
1853
1854proc isLieType()
1855"USAGE:  isLieType();
1856RETURN:  int, 1 if basering is a G-algebra of Lie type, 0 otherwise
1857PURPOSE: G-algebra of Lie type has relations of the kind Y*X=X*Y+D
1858EXAMPLE: example isLieType; shows an example
1859"
1860{
1861  def @B    = basering; //save the name of basering
1862  int NVars = nvars(@B); //number of variables in basering
1863  int i, j;
1864
1865  int answer = 1;
1866
1867  // check basering is of Lie type:
1868  matrix @@CC[NVars][NVars];
1869  for(i=1; i<NVars; i++)
1870  {
1871    for(j=i+1; j<=NVars; j++)
1872    {
1873      @@CC[i,j]=leadcoef(var(j)*var(i));
1874    }
1875  }
1876  ideal @C@ = simplify(ideal(@@CC),2+4);// skip zeroes and repeated entries
1877  if (  (size(@C@) >1 ) || ( (size(@C@)==1) && (@C@[1]!=1) )  )
1878  {
1879    answer = 0;
1880  }
1881  return(answer);
1882}
1883example
1884{
1885  "EXAMPLE:"; echo = 2;
1886  ring r = 0,(x,y),dp;
1887  y*x;
1888  isLieType(); //yes
1889  def D = Weyl(); setring D;
1890  y*x;
1891  isLieType(); //yes
1892  setring r;
1893  def R = nc_algebra(-3,0); setring R;
1894  y*x;
1895  isLieType(); // no
1896  kill R; kill r;
1897  ring s = (0,q),(x,y),dp;
1898  def S = nc_algebra(q,0); setring S;
1899  y*x;
1900  isLieType(); //no
1901  kill S; setring s;
1902  def S = nc_algebra(q,y^2); setring S;
1903  y*x;
1904  isLieType(); //no
1905}
Note: See TracBrowser for help on using the repository browser.