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

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