# source:git/Singular/LIB/invar.lib@5480da

spielwiese
Last change on this file since 5480da was 5480da, checked in by Kai Krüger <krueger@…>, 26 years ago
• Property mode set to `100644`
File size: 14.9 KB
Line
1// \$Id: invar.lib,v 1.5 1998-04-03 22:47:07 krueger Exp \$
2///////////////////////////////////////////////////////
3// invar.lib
4// algorithm for computing the ring of invariants under
5// the action of the additive group (C,+)
6// written by Gerhard Pfister
7//////////////////////////////////////////////////////
8
9version="\$Id: invar.lib,v 1.5 1998-04-03 22:47:07 krueger Exp \$";
10info="
11LIBRARY: invar.lib PROCEDURES FOR COMPUTING INVARIANTS OF (C,+)-ACTIONS
12
13  invariantRing(matrix m,poly p,poly q,int choose)
14  // ring of invariants of the action of the additive group
15  // defined by the vectorfield corresponding to the matrix m
16  // (m[i,1] are the coefficients of d/dx(i))
17  // the polys p and q are assumed to be variables x(i) and x(j)
18  // such that m[j,1]=0 and m[i,1]=x(j)
19  // if choose=0 the computation stops if generators of the ring
20  // of invariants are computed (to be used only if you know that
21  // the ring of invariants is finitey generated)
22  // if choose<>0 it computes invariants up to degree choose
23
24  actionIsProper(matrix m)
25  // returns 1 if the action of the additive group defined by the
26  // matrix m as above i proper and 0 if not.
27";
28
29///////////////////////////////////////////////////////////////////////////////
30
31LIB "inout.lib";
32LIB "general.lib";
33
34///////////////////////////////////////////////////////////////////////////////
35
36
37proc sortier(ideal id)
38{
39  if(size(id)==0)
40  {
41     return(id);
42  }
43  intvec i=sortvec(id);
44  int j;
45  ideal m;
46  for (j=1;j<=size(i);j++)
47  {
48    m[j] = id[i[j]];
49  }
50  return(m);
51}
52example
53{ "EXAMPLE:"; echo = 2;
54   ring q=0,(x,y,z,u,v,w),dp;
55   ideal i=w,x,z,y,v;
56   ideal j=sortier(i);
57   j;
58}
59
60
61///////////////////////////////////////////////////////////////////////////////
62
63
64proc der (matrix m, poly f)
65USAGE:   der(m,f);  m matrix, f poly
66RETURN:  poly= application of the vectorfield m befined by the matrix m
67         (m[i,1] are the coefficients of d/dx(i)) to f
68NOTE:
69EXAMPLE: example der; shows an example
70{
71  matrix mh=matrix(jacob(f))*m;
72  return(mh[1,1]);
73}
74example
75{ "EXAMPLE:"; echo = 2;
76   ring q=0,(x,y,z,u,v,w),dp;
77   poly f=2xz-y2;
78   matrix m[6][1];
79   m[2,1]=x;
80   m[3,1]=y;
81   m[5,1]=u;
82   m[6,1]=v;
83   der(m,f);
84}
85
86///////////////////////////////////////////////////////////////////////////////
87
88
89proc actionIsProper(matrix m)
90USAGE:   actionIsProper(m); m matrix
91RETURN:  int= 1 if is proper, 0 else
92NOTE:
93EXAMPLE: example actionIsProper; shows an example
94{
95  int i;
96  ideal id=maxideal(1);
97  def bsr=basering;
98
99  //changes the basering bsr to bsr[@t]
100  execute "ring s="+charstr(basering)+",("+varstr(basering)+",@t),dp;";
101  poly inv,delta,tee,j;
102  ideal id=imap(bsr,id);
103  matrix @m[size(id)+1][1];
104  @m=imap(bsr,m),0;
105
106  //computes the exp(@t*m)(var(i)) for all i
107  for(i=1;i<=nvars(basering)-1;i++)
108  {
109     inv=var(i);
110     delta=der(@m,inv);
111     j=1;
112     tee=@t;
113     while(delta!=0)
114     {
115        inv=inv+1/j*delta*tee;
116        j=j*(j+1);
117        tee=tee*@t;
118        delta=der(@m,delta);
119     }
120     id=id+ideal(inv);
121  }
122  i=inSubring(@t,id)[1];
123  setring(bsr);
124  return(i);
125}
126example
127{ "EXAMPLE:"; echo = 2;
128
129  ring rf=0,(x(1..7)),dp;
130  matrix m[7][1];
131  m[4,1]=x(1)^3;
132  m[5,1]=x(2)^3;
133  m[6,1]=x(3)^3;
134  m[7,1]=(x(1)*x(2)*x(3))^2;
135  actionIsProper(m);
136
137  ring rd=0,(x(1..5)),dp;
138  matrix m[5][1];
139  m[3,1]=x(1);
140  m[4,1]=x(2);
141  m[5,1]=1+x(1)*x(4)^2
142  actionIsProper(m);
143}
144///////////////////////////////////////////////////////////////////////////////
145
146
147proc reduction(poly p, ideal dom, list #)
148USAGE:   reduction(p,dom,q); p poly, dom ideal, q (optional) monomial
149RETURN:  poly= (p-H(f1,...,fr))/q^a, if Lt(p)=H(Lt(f1),...,Lt(fr)) for
150               some polynomial H in r variables over the base field,
151               a maximal such that q^a devides p-H(f1,...,fr),
152               dom =(f1,...,fr)
153NOTE:
154EXAMPLE: example reduction; shows an example
155{
156  int i;
157  int z=size(dom);
158  def bsr=basering;
159
160  //arranges the monomial v for elimination
161  poly v=var(1);
162  for (i=2;i<=nvars(basering);i=i+1)
163  {
164    v=v*var(i);
165  }
166
167  //changes the basering bsr to bsr[@(0),...,@(z)]
168  execute "ring s="+charstr(basering)+",("+varstr(basering)+",@(0..z)),dp;";
169
170  //costructes the ideal dom=(p-@(0),dom[1]-@(1),...,dom[z]-@(z))
171  ideal dom=imap(bsr,dom);
172  for (i=1;i<=z;i++)
173  {
175  }
177
178  //eliminates the variables of the basering bsr
179  //i.e. computes dom intersected with K[@(0),...,@(z)]
180  ideal kern=eliminate(dom,imap(bsr,v));
181
182  // test wether @(0)-h(@(1),...,@(z)) is in ker for some poly h
183  poly h;
184  z=size(kern);
185  for (i=1;i<=z;i++)
186  {
187     h=kern[i]/@(0);
188     if (deg(h)==0)
189     {
190        h=(1/h)*kern[i];
191        // defines the map psi : s ---> bsr defined by @(i) ---> p,dom[i]
192        setring bsr;
193        map psi=s,maxideal(1),p,dom;
194        poly re=psi(h);
195
196        // devides by the maximal power of #[1]
197        if (size(#)>0)
198        {
199           while ((re!=0) && (re!=#[1]) &&(subst(re,#[1],0)==0))
200           {
201             re=re/#[1];
202           }
203        }
204
205        return(re);
206     }
207  }
208  setring bsr;
209  return(p);
210}
211example
212{ "EXAMPLE:"; echo = 2;
213   ring q=0,(x,y,z,u,v,w),dp;
214   poly p=x2yz-x2v;
215   ideal dom =x-w,u2w+1,yz-v;
216   reduction(p,dom);
217   reduction(p,dom,w);
218}
219
220///////////////////////////////////////////////////////////////////////////////
221
222proc completeReduction(poly p, ideal dom, list #)
223USAGE:   completeReduction(p,dom,q); p poly, dom ideal,
224                                     q (optional) monomial
225RETURN:  poly= the polynomial p reduced with dom via the procedure
226               reduction as long as possible
227NOTE:
228EXAMPLE: example completeReduction; shows an example
229{
230  poly p1=p;
231  poly p2=reduction(p,dom,#);
232  while (p1!=p2)
233  {
234    p1=p2;
235    p2=reduction(p1,dom,#);
236  }
237  return(p2);
238}
239example
240{ "EXAMPLE:"; echo = 2;
241   ring q=0,(x,y,z,u,v,w),dp;
242   poly p=x2yz-x2v;
243   ideal dom =x-w,u2w+1,yz-v;
244   completeReduction(p,dom);
245   completeReduction(p,dom,w);
246}
247///////////////////////////////////////////////////////////////////////////////
248
249proc inSubring(poly p, ideal dom)
250USAGE:   inSubring(p,dom); p poly, dom ideal
251RETURN:  list= 1,string(@(0)-h(@(1),...,@(size(dom)))) :if p = h(dom[1],...,dom[size(dom)])
252              0,string(h(@(0),...,@(size(dom)))) :if there is only a nonlinear relation
253              h(p,dom[1],...,dom[size(dom)])=0.
254NOTE:
255EXAMPLE: example inSubring; shows an example
256{
257  int z=size(dom);
258  int i;
259  def gnir=basering;
260  list l;
261  poly mile=var(1);
262  for (i=2;i<=nvars(basering);i++)
263  {
264    mile=mile*var(i);
265  }
266  string eli=string(mile);
267  // the intersection of ideal nett=(p-@(0),dom[1]-@(1),...)
268  // with the ring k[@(0),...,@(n)] is computed, the result is ker
269  execute "ring r1="+charstr(basering)+",("+varstr(basering)+",@(0..z)),dp;";
270  ideal nett=imap(gnir,dom);
271  poly p;
272  for (i=1;i<=z;i++)
273  {
274    execute "p=@("+string(i)+");";
275    nett[i]=nett[i]-p;
276  }
277  nett=imap(gnir,p)-@(0),nett;
278  execute "ideal ker=eliminate(nett,"+eli+");";
279  // test wether @(0)-h(@(1),...,@(z)) is in ker
280  l[1]=0;
281  l[2]="";
282  for (i=1;i<=size(ker);i++)
283  {
284     if (deg(ker[i]/@(0))==0)
285     {
286        string str=string(ker[i]);
287        setring gnir;
288        l[1]=1;
289        l[2]=str;
290        return(l);
291     }
292     if (deg(ker[i]/@(0))>0)
293     {
294        l[2]=l[2]+string(ker[i]);
295     }
296  }
297  setring gnir;
298  return(l);
299}
300example
301{ "EXAMPLE:"; echo = 2;
302   ring q=0,(x,y,z,u,v,w),dp;
303   poly p=xyzu2w-1yzu2w2+u4w2-1xu2vw+u2vw2+xyz-1yzw+2u2w-1xv+vw+2;
304   ideal dom =x-w,u2w+1,yz-v;
305   inSubring(p,dom);
306}
307
308///////////////////////////////////////////////////////////////////////////////
309
310proc localInvar(matrix m, poly p, poly q, poly h)
311USAGE:   localInvar(m,p,q,h); m matrix, p,q,h poly
312RETURN:  poly= the invariant of the vectorfield m=Sum m[i,1]*d/dx(i) with respect
313               to p,q,h, i.e.
314               Sum (-1)^v*(1/v!)*m^v(p)*(q/m(q))^v)*m(q)^N, m^N(q)=0, m^(N-1)(q)<>0
315               it is assumed that m(q) and h are invariant
316               the sum above is divided by h as much as possible
317NOTE:
318EXAMPLE: example localInvar; shows an example
319{
320  if ((der(m,h) !=0) || (der(m,der(m,q)) !=0))
321  {
322    "the last variable defines not an invariant function ";
323    return(q);
324  }
325  poly inv=p;
326  poly dif= der(m,inv);
327  poly a=der(m,q);
328  poly sgn=-1;
329  poly coeff=sgn*q;
330  int k=1;
331  if (dif==0)
332  {
333    return(inv);
334  }
335  while (dif!=0)
336  {
337    inv=(a*inv)+(coeff*dif);
338    dif=der(m,dif);
339    k=k+1;
340    coeff=q*coeff*sgn/k;
341  }
342  while ((inv!=0) && (inv!=h) &&(subst(inv,h,0)==0))
343 {
344   inv=inv/h;
345  }
346  return(inv);
347}
348example
349{ "EXAMPLE:"; echo = 2;
350   ring q=0,(x,y,z),dp;
351   matrix m[3][1];
352   m[2,1]=x;
353   m[3,1]=y;
354   poly in=localInvar(m,z,y,x);
355   in;
356}
357///////////////////////////////////////////////////////////////////////////////
358
359
360
361proc furtherInvar(matrix m, ideal id, ideal karl, poly q)
362USAGE:   furtherInvar(m,id,karl,q); m matrix, id,karl ideal,q poly
363RETURN:  ideal= further invariants of the vectorfield m=Sum m[i,1]*d/dx(i) with respect
364               to id,p,q, i.e.
365               the ideal id contains invariants of m and we are looking for elements
366               in the subring generated by id which are divisible by q
367               it is assumed that m(p) and q are invariant
368               the elements mentioned  above are computed and divided by q
369               as much as possible
370               the ideal karl contains all invariants computed yet
371NOTE:
372EXAMPLE: example furtherInvar; shows an example
373{
374  int i;
375  ideal null;
376  int z=size(id);
377  intvec v;
378  def @r=basering;
379  ideal su;
380  for (i=1;i<=z;i++)
381  {
382    su[i]=subst(id[i],q,0);
383  }
384  // defines the map phi : r1 ---> @r defined by
385  // y(i) ---> id[i](q=0)
386  execute "ring r1="+charstr(basering)+",(y(1..z)),dp";
387  setring @r;
388  map phi=r1,su;
389  setring r1;
390  // computes the kernel of phi
391  execute "ideal ker=preimage(@r,phi,null)";
392  // defines the map psi : r1 ---> @r defined by y(i) ---> id[i]
393  setring @r;
394  map psi=r1,id;
395  // computes psi(ker(phi))
396  ideal rel=psi(ker);
397  // devides by the maximal power of q
398  // and tests wether we really obtain invariants
399  for (i=1;i<=size(rel);i++)
400  {
401    while ((rel[i]!=0) && (rel[i]!=q) &&(subst(rel[i],q,0)==0))
402    {
403      rel[i]=rel[i]/q;
404      if (der(m,rel[i])!=0)
405      {
406         "error in furtherInvar, function not invariant";
407         rel[i];
408      }
409    }
410    rel[i]=simplify(rel[i],1);
411  }
412  // test whether some variables occur linearly
413  // and deletes the corresponding invariant function
414  setring r1;
415  int j;
416  for (i=1;i<=size(ker);i=i+1)
417  {
418     for (j=1;j<=z;j++)
419     {
420        if (deg(ker[i]/y(j))==0)
421        {
422           setring @r;
423           rel[i]= completeReduction(rel[i],karl,q);
424           if(rel[i]!=0)
425           {
426              karl[j+1]=rel[i];
427              rel[i]=0;
428           }
429           setring r1;
430        }
431     }
432
433  }
434  setring @r;
435  list l=rel+null,karl;
436  return(l);
437}
438example
439{ "EXAMPLE:"; echo = 2;
440   ring r=0,(x,y,z,u),dp;
441   matrix m[4][1];
442   m[2,1]=x;
443   m[3,1]=y;
444   m[4,1]=z;
445   ideal id=localInvar(m,z,y,x),localInvar(m,u,y,x);
446   ideal karl=id,x;
447   list in=furtherInvar(m,id,karl,x);
448   in;
449}
450///////////////////////////////////////////////////////////////////////////////
451
452
453
454proc invariantRing(matrix m, poly p, poly q,list #)
455USAGE:   invariantRing(m,p,q); m matrix, p,q poly
456RETURN:  ideal= the invariants of the vectorfield m=Sum m[i,1]*d/dx(i)
457                p,q variables with m(p)=q invariant
458NOTE:
459EXAMPLE: example furtherInvar; shows an example
460{
461  ideal j;
462  int i,it;
463  int bou=-1;
464  if(size(#)>0)
465  {
466     bou=#[1];
467  }
468  int z;
469  ideal karl;
470  ideal k1=1;
471  list k2;
472  // computation of local invariants
473  for (i=1;i<=nvars(basering);i++)
474  {
475    karl=karl+localInvar(m,var(i),p,q);
476  }
477  if(bou==0)
478  {
479     "                     ";
480     "the local invariants:";
481     "                     ";
482     karl;
483    // pause;
484     "                     ";
485  }
486  // computation of further invariants
487  it=0;
488  while (size(k1)!=0)
489 {
490    // test if the new invariants are already in the ring generated
491    // by the invariants we constructed already
492    it++;
493    karl=sortier(karl);
494    j=q;
495    for (i=1;i<=size(karl);i++)
496    {
497       j=j+ simplify(completeReduction(karl[i],j,q),1);
498    }
499    karl=j;
500    j[1]=0;
501    j=simplify(j,2);
502    k2=furtherInvar(m,j,karl,q);
503    k1=k2[1];
504    karl=k2[2];
505    k1=sortier(k1);
506    z=size(k1);
507    for (i=1;i<=z;i++)
508    {
509      k1[i]= completeReduction(k1[i],karl,q);
510      if (k1[i]!=0)
511      {
512        karl=karl+simplify(k1[i],1);
513      }
514    }
515    if(bou==0)
516    {
517       "                                  ";
518       "the invariants after the iteration";
519       it;
520       "                                  ";
521       karl;
522      // pause;
523       "                                  ";
524    }
525    if((bou>0) && (size(k1)>0))
526    {
527      if(deg(k1[size(k1)])>bou)
528      {
529         return(karl);
530      }
531    }
532  }
533  return(karl);
534}
535example
536{ "EXAMPLE:"; echo = 2;
537
538  //Winkelmann: free action but Spec k[x(1),...,x(5)]---> Spec In-
539  //variantring is not surjective
540
541  ring rw=0,(x(1..5)),dp;
542  matrix m[5][1];
543  m[3,1]=x(1);
544  m[4,1]=x(2);
545  m[5,1]=1+x(1)*x(4)+x(2)*x(3);
546  ideal in=invariantRing(m,x(3),x(1),0);
547  in;
548
549  //Deveney/Finston: The ring of invariants is not finitely generated
550
551  ring rf=0,(x(1..7)),dp;
552  matrix m[7][1];
553  m[4,1]=x(1)^3;
554  m[5,1]=x(2)^3;
555  m[6,1]=x(3)^3;
556  m[7,1]=(x(1)*x(2)*x(3))^2;
557  ideal in=invariantRing(m,x(4),x(1),6);
558  in;
559
560
561  //Deveney/Finston:Proper Ga-action which is not locally trivial
562  //r[x(1),...,x(5)] is not flat over the ring of invariants
563
564  ring rd=0,(x(1..5)),dp;
565  matrix m[5][1];
566  m[3,1]=x(1);
567  m[4,1]=x(2);
568  m[5,1]=1+x(1)*x(4)^2;
569  ideal in=invariantRing(m,x(3),x(1));
570  in;
571
572  actionIsProper(m);
573
574  //computes the relations between the invariants
575  int z=size(in);
576  ideal null;
577  ring r1=0,(y(1..z)),dp;
578  setring rd;
579  map phi=r1,in;
580  setring r1;
581  ideal ker=preimage(rd,phi,null);
582  ker;
583
584  //the discriminant
585
586  ring r=0,(x(1..2),y(1..2),z,t),dp;
587  poly p=z+(1+x(1)*y(2)^2)*t+x(1)*y(1)*y(2)*t^2+(1/3)*x(1)*y(1)^2*t^3;
588
589  matrix m[5][5];
590  m[1,1]=z;
591  m[1,2]=x(1)*y(2)^2+1;
592  m[1,3]=x(1)*y(1)*y(2);
593  m[1,4]=1/3*x(1)*y(1)^2;
594  m[1,5]=0;
595  m[2,1]=0;
596  m[2,2]=z;
597  m[2,3]=x(1)*y(2)^2+1;
598  m[2,4]=x(1)*y(1)*y(2);
599  m[2,5]=1/3*x(1)*y(1)^2;
600  m[3,1]=x(1)*y(2)^2+1;
601  m[3,2]=2*x(1)*y(1)*y(2);
602  m[3,3]=x(1)*y(1)^2;
603  m[3,4]=0;
604  m[3,5]=0;
605  m[4,1]=0;
606  m[4,2]=x(1)*y(2)^2+1;
607  m[4,3]=2*x(1)*y(1)*y(2);
608  m[4,4]=x(1)*y(1)^2;
609  m[4,5]=0;
610  m[5,1]=0;
611  m[5,2]=0;
612  m[5,3]=x(1)*y(2)^2+1;
613  m[5,4]=2*x(1)*y(1)*y(2);
614  m[5,5]=x(1)*y(1)^2;
615
616  poly disc=9*det(m)/(x(1)^2*y(1)^4);
617
618  LIB "invar.lib";
619  matrix n[6][1];
620  n[2,1]=x(1);
621  n[4,1]=y(1);
622  n[5,1]=1+x(1)*y(2)^2;
623
624  der(n,disc);
625
626  //x(1)^3*y(2)^6-6*x(1)^2*y(1)*y(2)^3*z+6*x(1)^2*y(2)^4+9*x(1)*y(1)^2*z^2-18*x(1)*y(1)*y(2)*z+9*x(1)*y(2)^2+4
627
628
629  //constructive approach to Weizenbcks theorem
630
631  int n=5;
632
633  ring w=0,(x(1..n)),wp(1..n);
634
635  // definition of the vectorfield m=sum m[i]*d/dx(i)
636  matrix m[n][1];
637  int i;
638  for (i=1;i<=n-1;i=i+1)
639  {
640    m[i+1,1]=x(i);
641  }
642  ideal in=invariantRing(m,x(2),x(1),0);
643  in;
644
645
646
647}
Note: See TracBrowser for help on using the repository browser.