source: git/Singular/LIB/presolve.lib @ 2e7859

fieker-DuValspielwiese
Last change on this file since 2e7859 was 38c165, checked in by Kai Krüger <krueger@…>, 24 years ago
Changed excute syntax from command to function git-svn-id: file:///usr/local/Singular/svn/trunk@4318 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 43.4 KB
Line 
1// $Id: presolve.lib,v 1.13 2000-05-12 12:17:17 krueger Exp $
2//(GMG)
3///////////////////////////////////////////////////////////////////////////////
4
5version="$Id: presolve.lib,v 1.13 2000-05-12 12:17:17 krueger Exp $";
6info="
7LIBRARY:  presolve.lib     PROCEDURES FOR PRE-SOLVING POLYNOMIAL EQUATIONS
8AUTHOR:   Gert-Martin Greuel, email: greuel@mathematik.uni-kl.de,
9
10PROCEDURES:
11 degreepart(id,d1,d2);  elements of id of total degree >= d1 and <= d2
12 elimlinearpart(id);    linear part eliminated from id
13 elimpart(id[,n]);      partial elimination of vars [among 1-st n vars]
14 elimpartanyr(i,p);     factors of p partially eliminated from i in any ring
15 fastelim(i,p[..]);     fast elimination of factors of p from i [options]
16 findvars(id[..]);      ideal of variables occuring in id [more information]
17 hilbvec(id[,c,o]);     intvec of Hilberseries of id [in char c and ord o]
18 linearpart(id);        elements of id of total degree <=1
19 tolessvars(id[,]);     maps id to new basering having only vars occuring in id
20 solvelinearpart(id);   reduced std-basis of linear part of id
21 sortandmap(id,s1,s2);  map to new basering with vars sorted w.r.t. complexity
22 sortvars(id[n1,p1..]); sort vars w.r.t. complexity in id [different blocks]
23 valvars(id[..]);       valuation of vars w.r.t. to their complexity in id
24           (parameters in square brackets [] are optional)
25";
26
27LIB "inout.lib";
28LIB "general.lib";
29LIB "matrix.lib";
30LIB "ring.lib";
31LIB "elim.lib";
32///////////////////////////////////////////////////////////////////////////////
33
34proc degreepart (id,int d1,int d2,list #)
35"USAGE:   degreepart(id,d1,d2[,v]);  id=ideal/module, d1,d1=integers, v=intvec
36RETURN:  generators of id of [v-weighted] total degree >= d1 and <= d2
37         (default: v = 1,...,1)
38EXAMPLE: example degreepart; shows an example
39"
40{
41   def dpart = id;
42   int s,ii = size(id),0;
43   if ( size(#)==0 )
44   {
45      for ( ii=1; ii<=s; ii=ii+1 )
46      {
47         dpart[ii] = (jet(id[ii],d1-1)==0)*(id[ii]==jet(id[ii],d2))*id[ii];
48      }
49   }
50   else
51   {
52      for ( ii=1; ii<=s; ii=ii+1 )
53      {
54         dpart[ii]=(jet(id[ii],d1-1,#[1])==0)*(id[ii]==jet(id[ii],d2,#[1]))*id[ii];
55      }
56   }
57   return(simplify(dpart,2));
58}
59example
60{ "EXAMPLE:"; echo = 2;
61   ring r=0,(x,y,z),dp;
62   ideal i=1+x+x2+x3+x4,3,xz+y3+z4z4;
63   degreepart(i,0,4);
64   module m=[x,y,z],x*[x3,y2,z],[1,x2,z3,0,1];
65   intvec v=2,3,6;
66   show(degreepart(m,8,8,v));
67}
68///////////////////////////////////////////////////////////////////////////////
69
70proc elimlinearpart (ideal i,list #)
71"USAGE:   elimlinearpart(i[,n]);  i=ideal, n=integer
72RETURN:  list of of 5 objects:
73         [1]: (interreduced) ideal obtained from i by eliminating (substitute)
74              from the first n variables those which appear in a linear part
75              of i, by putting this part into triangular form
76              (default: n = nvars(basering))
77         [2]: ideal of variables which have been eliminated (= substituted)
78         [3]: ideal, j-th element defines substitution of j-th var in [2]
79         [4]: ideal of variables of basering, eliminated ones are set to 0
80         [5]: ideal, describing the map from the basering to itself such that
81              [1] is the image of i
82NOTE:    the procedure does always interreduces the ideal i internally w.r.t.
83         ordering dp
84         // bei ** spaeter eventuell verbessern
85EXAMPLE: example elimlinearpart; shows an example
86"
87{
88   int ii,n,fi,k;
89   string o, newo;
90   intvec getoption = option(get);
91   option(redSB);
92   def P = basering;
93   n = nvars(P);
94//--------------- replace ordering by dp-ordering if necessary ----------------
95   o = "dp("+string(n)+")";
96   fi = find(ordstr(P),o);
97   if( fi == 0 or find(ordstr(P),"a") != 0 )
98   {
99      execute("ring newP = ("+charstr(P)+"),("+varstr(P)+"),dp;");
100      ideal i = imap(P,i);
101   }
102//---------------------------------- start ------------------------------------
103   if ( size(#)!=0 ) {  n=#[1]; }
104   ideal maxi,rest = maxideal(1),0;
105   if ( n < nvars(P) ) { rest = maxi[n+1..nvars(P)]; }
106   attrib(rest,"isSB",1);
107//-------------------- interreduce and find linear part  ----------------------
108// interred does the only real work. Because of ordering dp the linear part is
109// within the first elements after interreduction
110// **: perhaps Bareiss to constant matrix of linear part instead of interred,
111// and/or for big systems, interred only those generators of id
112// which do not contain elements not to be eliminated
113
114   ideal id = interred(i);
115
116   if(id[1] == 1 )                          //special case of unit ideal
117   {
118     setring P;
119     list L = ideal(1), ideal(0), ideal(0), maxideal(1), maxideal(1);
120     option(set,getoption);
121     return(L);
122   }
123 
124   for ( ii=1; ii<=size(id); ii++ )
125   {
126      if( deg(id[ii]) > 1) { break; }
127      k=ii;
128   }
129   if( k == 0 )       { ideal lin; }
130   else               { ideal lin = id[1..k];}
131   if( k < size(id) ) { id = id[k+1..size(id)]; }
132   else               { id = 0; }
133//----- remove generators from lin containing vars not to be eliminated  ------
134   if ( n < nvars(P) )
135   {
136      for ( ii=1; ii<=size(lin); ii++ )
137      {
138         if ( reduce(lead(lin[ii]),rest) == 0 )
139         {
140            id=lin[ii],id;
141            lin[ii] = 0;
142         }
143      }
144   }
145   lin = simplify(lin,2);
146   attrib(lin,"isSB",1);
147   ideal eva = lead(lin);
148   attrib(eva,"isSB",1);
149   ideal neva = reduce(maxideal(1),eva);
150//------------------ go back to original ring end return  ---------------------
151   if( fi == 0 or find(ordstr(P),"a") != 0 )
152   {
153      setring P;
154      ideal id = imap(newP,id);
155      ideal eva = imap(newP,eva);
156      ideal lin = imap(newP,lin);
157      ideal neva = imap(newP,neva);
158   }
159 
160   eva = eva[ncols(eva)..1];  // sorting according to variables in basering
161   lin = lin[ncols(lin)..1];
162   ideal phi= neva;
163   k = 1;
164   for( ii=1; ii<=n; ii++ )
165   {
166      if( neva[ii] == 0 )
167      {
168         phi[ii] = eva[k]-lin[k];
169         k=k+1;
170      }
171   }
172   list L = id, eva, lin, neva, phi;
173   option(set,getoption);
174   return(L);
175}
176example
177{ "EXAMPLE:"; echo = 2;
178   ring s=0,(x,y,z,t,u,v,w,a,b,c,d,f,e),ds;
179   ideal i = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
180             d2+e2-1, f4+2u, wa+tf, xy+tu+ab, yz+uv+bc,
181             cd+ze, x+y+z+e+1, t+u+v+f-1, w+a+b+c+d;
182  list L= elimlinearpart(i);
183}
184///////////////////////////////////////////////////////////////////////////////
185
186proc elimpart (ideal i,list #)
187"USAGE:   elimpart(i[,n,e]);  i=ideal, n,e=integers
188         consider 1-st n vars for elimination (better: substitution),
189         e =0: substitute from linear part of i (same as elimlinearpart)
190         e!=0: eliminate also by direct substitution
191         (default: n = nvars(basering), e = 1)
192RETURN:  list of of 5 objects:
193         [1]: ideal obtained by substituting from the first n variables those
194              from i which appear in the linear part of i [or, if e!=0, which
195              can be expressed directly in the remaining vars]
196         [2]: ideal, variables which have been substituted
197         [3]: ideal, i-th element defines substitution of i-th var in [2]
198         [4]: ideal of variables of basering, substituted ones are set to 0
199         [5]: ideal, describing the map from the basering, say k[x(1..m)], to
200              itself onto k[..variables fom [4]..] and [1] is the image of i
201         The ideal i is generated by [1] and [3] in k[x(1..m)], the map [5]
202         maps [3] to 0, hence induceds an isomorhism
203                   k[x(1..m)]/i -> k[..variables fom [4]..]/[1]
204NOTE:    If the basering has ordering (c,dp), this is faster for big ideals,
205         since it avoids internal ring change and mapping
206EXAMPLE: example elimpart; shows an example
207"
208{
209   def P = basering;
210   int n,e = nvars(P),1;
211   if ( size(#)==1 ) {  n=#[1]; }
212   if ( size(#)==2 ) {  n=#[1]; e=#[2];}
213//----------- interreduce linear part with proc elimlinearpart ----------------
214// lin = ideal after interreduction of linear part
215// eva = eliminated (substituted) variables
216// sub = polynomials defining substitution
217// neva= not eliminated variables
218// phi = map describing substitution
219   list L = elimlinearpart(i,n);
220   ideal lin, eva, sub, neva, phi = L[1], L[2], L[3], L[4], L[5];
221//-------- direct substitution of variables if possible and if e!=0 -----------
222// first find terms lin1 of pure degree 1 in each poly of lin
223// k1 = pure degree 1 part, k1+k2 = those polys of lin which contained a pure
224// degree 1 part.
225// Then go to ring newP with ordering c,dp(n) and create a matrix with size(k1)
226// colums and 2 rows, such that if [f1,f2] is a column of M then f1+f2 is one
227// of the polys of lin containing a pure degree 1 part and f1 is this part
228// interreduce this matrix (i.e. Gauss elimination on linear part, with rest
229// transformed accordingly).
230   int ii, kk;
231   if ( e!=0 )
232   {
233      ideal k1,k2;
234      int l = size(lin);
235      ideal lin1 = jet(lin,1) - jet(lin,0);  // part of pure degree 1
236      lin = lin - lin1;                      // rest, part of degree 1 substracted
237
238      for( ii=1; ii<=l; ii++ )
239      {
240         if( lin1[ii] != 0 )
241         {
242            kk = kk+1;
243            k1[kk] = lin1[ii];    // part of pure degree 1, renumbered
244            k2[kk] = lin[ii];     // rest of those polys which had a degree 1 part
245            lin[ii] = 0;
246         }
247      }
248
249      if( kk != 0 )
250      {
251         if( ordstr(P) != "c,dp(n)" )
252         {
253            execute("ring newP = ("+charstr(P)+"),("+varstr(P)+"),(c,dp);");
254            ideal k1 = imap(P,k1);
255            ideal k2 = imap(P,k2);
256            ideal lin = imap(P,lin);
257            ideal eva = imap(P,eva);
258            ideal sub = imap(P,sub);
259            ideal neva = imap(P,neva);
260         }
261         ideal k12 = k1,k2;
262         matrix M = matrix(k12,2,kk);
263        // M = interred(M);
264         l = ncols(M);
265         k1 = M[1,1..l];
266         k2 = M[2,1..l];
267         ideal kin = matrix(k1)+matrix(k2);
268         lin = simplify(lin,2);
269
270         l = size(kin);
271         poly p; map phi; ideal max;
272         for ( ii=1; ii<=n; ii++  )
273         {
274            for (kk=1; kk<=l; kk++ )
275            {
276               p = kin[kk]/var(ii);
277               if ( deg(p) == 0 )
278               {
279                  eva = eva+var(ii);
280                  neva[ii] = 0;
281                  sub = sub+kin[kk]/p;
282                  max = maxideal(1);
283                  max[ii] = var(ii) - (kin[kk]/p);
284                  phi = basering,max;
285                  lin = phi(lin);
286                  kin = simplify(phi(kin),2);
287                  l = size(kin);
288                  ii=ii+1;
289                  break;
290               }
291            }
292         }
293         lin = kin+lin;
294      }
295   }
296   for( ii=1; ii<=size(lin); ii++ )
297   {
298      lin[ii] = cleardenom(lin[ii]);
299   }
300   if( defined(newP) )
301   {
302      setring P;
303      lin = imap(newP,lin);
304      eva = imap(newP,eva);
305      sub = imap(newP,sub);
306      neva = imap(newP,neva);
307   }
308   for( ii=1; ii<=n; ii++ )
309   {
310      for( kk=1; kk<=size(eva); kk++ )
311      {
312         if (phi[ii] == eva[kk] )
313         {  phi[ii] = eva[kk]-sub[kk]; break; }
314      }
315   }
316   L = lin, eva, sub, neva, phi;
317  return(L);
318}
319example
320{ "EXAMPLE:"; echo = 2;
321   ring s=0,(x,y,z,t,u,v,w,a,b,c,d,f,e),(c,ds);
322   ideal i = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
323             d2+e2-1, f4+2u, wa+tf, xy+tu+ab, yz+uv+bc,
324             cd+ze, x+y+z+e+1, t+u+v+f-1, w+a+b+c+d;
325   elimpart(i,4);
326}
327
328///////////////////////////////////////////////////////////////////////////////
329
330proc elimpartanyr (ideal i, list #)
331"USAGE:   elimpartanyr(i[,p,e]);  i=ideal, p=product of vars to be eliminated,
332         e=int (default: p=product of all vars, e=1)
333RETURN:  list of of 6 objects:
334         [1]: (interreduced) ideal obtained by substituting from i those vars
335              appearing in p which occur in the linear part of i [or which can
336              be expressed directly in the remaining variables, if e!=0]
337         [2]: ideal, variables which have been substituted
338         [3]: ideal, i-th element defines substitution of i-th var in [2]
339         [4]: ideal of variables of basering, substituted ones are set to 0
340         [5]: ideal, describing the map from the basering, say k[x(1..m)], to
341              itself onto k[..variables fom [4]..] and [1] is the image of i
342         [6]: int, # of vars considered for substitution (= # of factors of p)
343
344         The ideal i is generated by [1] and [3] in k[x(1..m)], the map [5]
345         maps [3] to 0, hence induceds an isomorhism
346                   k[x(1..m)]/i -> k[..variables fom [4]..]/[1]
347NOTE:    the proc uses 'execute' to create a ring with ordering dp and vars
348         placed correctly and then applies 'elimpart';
349EXAMPLE: example elimpartanyr; shows an example
350"
351{
352   def P = basering;
353   int j,n,e = 0,0,1;
354   poly p = product(maxideal(1));
355   if ( size(#)==1 ) { p=#[1]; }
356   if ( size(#)==2 ) { p=#[1]; e=#[2]; }
357   string a,b;
358   for ( j=1; j<=nvars(P); j++ )
359   {
360      if (deg(p/var(j))>=0) { a = a+varstr(j)+","; n = n+1; }
361      else { b = b+varstr(j)+","; }
362   }
363   if ( size(b) != 0 ) { b = b[1,size(b)-1]; }
364   else { a = a[1,size(a)-1]; }
365   execute("ring gnir ="+charstr(P)+",("+a+b+"),dp;");
366   ideal i = imap(P,i);
367   list L = elimpart(i,n,e)+list(n);
368   setring P;
369   list L = imap(gnir,L);
370   return(L);
371}
372example
373{ "EXAMPLE:"; echo = 2;
374   ring s=0,(e,f,x,y,z,t,u,v,w,a,b,c,d),dp;
375   ideal i = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
376             d2+e2-1, f4+2u, wa+tf, xy+tu+ab, yz+uv+bc,
377             cd+ze, x+y+z+e+1, t+u+v+f-1, w+a+b+c+d;
378   show(elimpartanyr(i,xyztuvwabc));"";
379   elimpartanyr(i,xyztuvwabc);
380}
381///////////////////////////////////////////////////////////////////////////////
382
383proc fastelim (ideal i, poly p, list #)
384"USAGE:   fastelim(i,p[h,o,a,b,e,m]); i=ideal, p=product of variables to be
385         eliminated; h,o,a,b,e integers
386         (options for Hilbert-std, 'valvars', elimpart, minimizing)
387         - h !=0: use Hilbert-series driven std-basis computation
388         - o !=0: use proc 'valvars' for a - hopefully - best ordering of vars
389         - a !=0: order vars to be eliminated w.r.t. increasing complexity
390         - b !=0: order vars not to be eliminated w.r.t. increasing complexity
391         - e !=0: use elimpart first to eliminate easy part
392         - m !=0: compute a minimal system of generators
393         replacing '!=' by '=' has the opposite meaning
394         default: h,o,a,b,e,m = 0,1,0,0,0,0
395RETURN:  ideal obtained from i by eliminating those variables which occur in p
396EXAMPLE: example fastelim; shows an example.
397"
398{
399   def P = basering;
400   int h,o,a,b,e,m = 0,1,0,0,0,0;
401   if ( size(#) == 1 ) { h=#[1]; }
402   if ( size(#) == 2 ) { h=#[1]; o=#[2]; }
403   if ( size(#) == 3 ) { h=#[1]; o=#[2]; a=#[3]; }
404   if ( size(#) == 4 ) { h=#[1]; o=#[2]; a=#[3]; b=#[4];}
405   if ( size(#) == 5 ) { h=#[1]; o=#[2]; a=#[3]; b=#[4]; e=#[5]; }
406   if ( size(#) == 6 ) { h=#[1]; o=#[2]; a=#[3]; b=#[4]; e=#[5]; m=#[6]; }
407   list L = elimpartanyr(i,p,e);
408   poly q = product(L[2]);     //product of vars which are already eliminated
409   if ( q==0 ) { q=1; }
410   p = p/q;                    //product of vars which must still be eliminated
411   int nu = size(L[5])-size(L[2]);   //number of vars which must still be eliminated
412   if ( p==1 )                 //ready if no vars are left
413   {                           //compute minbase if 3-rd argument !=0
414      if ( m != 0 ) { L[1]=minbase(L[1]); }
415      return(L);
416   }
417//---------------- create new ring with remaining variables -------------------
418   string newvar = string(L[4]);
419   L = L[1],p;
420   execute("ring r1=("+charstr(P)+"),("+newvar+"),"+"dp;");
421   list L = imap(P,L);
422//------------------- find "best" ordering of variables  ----------------------
423   newvar = string(maxideal(1));
424   if ( o != 0 )
425   {
426      list ordevar = valvars(L[1],a,L[2],b);
427      intvec v = ordevar[1];
428      newvar=string(sort(maxideal(1),v)[1]);
429//------------ create new ring with "best" ordering of variables --------------
430      changevar("r0",newvar);
431      list L = imap(r1,L);
432      kill r1;
433      def r1 = r0;
434      kill r0;
435   }
436//----------------- h==0: eliminate remaining vars directly -------------------
437   if ( h == 0 )
438   {
439      L[1] = eliminate(L[1],L[2]);
440      def r2 = r1;
441   }
442   else
443//------- h!=0: homogenize and compute Hilbert-series using hilbvec ----------
444   {
445      intvec hi = hilbvec(L[1]);         // Hilbert-series of i
446      execute("ring r2=("+charstr(P)+"),("+varstr(basering)+",@homo),dp;");
447      list L = imap(r1,L);
448      L[1] = homog(L[1],@homo);          // @homo = homogenizing var
449//---- use Hilbert-series to eliminate variables with Hilbert-driven std -----
450      L[1] = eliminate(L[1],L[2],hi);
451      L[1]=subst(L[1],@homo,1);          // dehomogenize by setting @homo=1
452   }
453   if ( m != 0 )                         // compute minbase
454   {
455      if ( #[1] != 0 ) { L[1] = minbase(L[1]); }
456   }
457   def id = L[1];
458   setring P;
459   return(imap(r2,id));
460}
461example
462{ "EXAMPLE:"; echo = 2;
463   ring s=31991,(e,f,x,y,z,t,u,v,w,a,b,c,d),dp;
464   ideal i = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
465            d2+e2-1, f4+2u, wa+tf, xy+tu+ab;
466   fastelim(i,xytua);           //with valvars only
467   fastelim(i,xytua,1,1);       //with hilb,valvars,minbase
468   fastelim(i,xytua,1,0);       //with hilb,minbase
469}
470///////////////////////////////////////////////////////////////////////////////
471
472proc faststd (@id,string @s1,string @s2, list #)
473"USAGE:   faststd(id,s1,s2[,\"hilb\",\"sort\",\"dec\",o,\"blocks\"]);
474         id=ideal/module, s1,s2=strings (names for new ring and maped id)
475         o = string (allowed ordstring:\"lp\",\"dp\",\"Dp\",\"ls\",\"ds\",\"Ds\")
476         \"hilb\",\"sort\",\"dec\",\"block\" options for Hilbert-std, sortandmap
477COMPUTE: create a new ring (with \"best\" ordering of vars) and compute a
478         std-basis of id (hopefully faster)
479         - If say, s1=\"R\" and s2=\"j\", the new basering has name R and the
480           std-basis of the image of id in R has name j
481         - \"hilb\"  : use Hilbert-series driven std-basis computation
482         - \"sort\"  : use 'sortandmap' for a best ordering of vars
483         - \"dec\"   : order vars w.r.t. decreasing complexity (with \"sort\")
484         - \"block\" : create blockordering, each block having ordstr=o, s.t.
485                     vars of same complexity are in one block (with \"sort\")
486         - o       : defines the basic ordering of the resulting ring
487         default: o = ordering of 1-st block of basering - if it is allowed,
488                      else o=\"dp\",
489                  \"sort\", if none of the optional parameters is given
490RETURN:  nothing
491NOTE:    This proc is only useful for hard problems where other methods fail.
492         \"hilb\" is useful for hard orderings (as \"lp\") or for characteristic 0,
493         it is correct for \"lp\",\"dp\",\"Dp\" (and for blockorderings combining
494         these) but not for s-orderings or if the vars have different weights.
495         There seem to be only few cases in which \"dec\" is fast
496EXAMPLE: example faststd; shows an example.
497"
498{
499   def @P = basering;
500   int @h,@s,@n,@m,@ii = 0,0,0,0,0;
501   string @o,@va,@c = ordstr(basering),"","";
502//-------------------- prepare ordering and set options -----------------------
503   if ( @o[1]=="c" or @o[1]=="C")
504      {  @o = @o[3,2]; }
505   else
506      { @o = @o[1,2]; }
507   if( @o[1]!="d" and @o[1]!="D" and @o[1]!="l")
508      { @o="dp"; }
509
510   if (size(#) == 0 )
511      { @s = 1; }
512   for ( @ii=1; @ii<=size(#); @ii++ )
513   {
514      if ( typeof(#[@ii]) != "string" )
515      {
516         "// wrong syntax! type: help faststd";
517         return();
518      }
519      else
520      {
521         if ( #[@ii] == "hilb"  ) { @h = 1; }
522         if ( #[@ii] == "dec"   ) { @n = 1; }
523         if ( #[@ii] == "block" ) { @m = 1; }
524         if ( #[@ii] == "sort"  ) { @s = 1; }
525         if ( #[@ii]=="lp" or #[@ii]=="dp" or #[@ii]=="Dp" or #[@ii]=="ls"
526              or #[@ii]=="ds" or #[@ii]=="Ds" ) { @o = #[@ii]; }
527      }
528   }
529   if( voice==2 ) { "// choosen options, hilb sort dec block:",@h,@s,@n,@m; }
530
531//-------------------- nosort: create ring with new name ----------------------
532   if ( @s==0 )
533   {
534      execute("ring "+@s1+"=("+charstr(@P)+"),("+varstr(@P)+"),("+@o+");");
535      def @id = imap(@P,@id);
536      verbose(noredefine);
537      def @P = basering;
538      verbose(redefine);
539      kill `@s1`;
540      if ( @h==0 ) { @id = std(@id); }
541   }
542//--------- sort: create new ring with "best" ordering of variables -----------
543 proc bestorder(@id,string @s1,string @s2,int @n,string @o,int @m,int @l)
544  {
545     intvec getoption = option(get);
546     option(redSB);
547     @id = interred(sort(@id)[1]);
548     poly @p = product(maxideal(1),1..@l);
549     def i,s1,s2,n,p,o,m = @id,@s1,@s2,@n,@p,@o,@m;
550     sortandmap(i,s1,s2,n,p,0,o,m);
551     option(set,getoption);
552     keepring(basering);
553  }
554//---------------------- no hilb: compute SB directly -------------------------
555   if ( @s != 0 and @h == 0 )
556   {
557      bestorder(@id,@s1,@s2,@n,@o,@m,nvars(@P));
558      verbose(noredefine);
559      def @P = basering;
560      verbose(redefine);
561      kill `@s1`;
562      def @id = `@s2`;
563      @id = std(@id);
564   }
565//------- hilb: homogenize and compute Hilbert-series using hilbvec -----------
566// this uses another standardbasis computation
567   if ( @h != 0 )
568   {
569      execute("ring @Q=("+charstr(@P)+"),("+varstr(@P)+",@homo),("+@o+");");
570      def @id = imap(@P,@id);
571      kill @P;
572      @id = homog(@id,@homo);               // @homo = homogenizing var
573      if ( @s != 0 )
574      {
575         bestorder(@id,@s1,@s2,@n,@o,@m,nvars(@Q)-1);
576         verbose(noredefine);
577         def @Q= basering;
578         kill `@s1`;
579         def @id = `@s2`;
580         verbose(redefine);
581      }
582      intvec @hi;                     // encoding of Hilbert-series of i
583      @hi = hilbvec(@id);
584      //if ( @s!=0 ) { @hi = hilbvec(@id,"32003",ordstr(@Q)); }
585      //else { @hi = hilbvec(@id); }
586//-------------------------- use Hilbert-driven std --------------------------
587      @id = std(@id,@hi);
588      @id = subst(@id,@homo,1);             // dehomogenize by setting @homo=1
589      @va = varstr(@Q)[1,size(varstr(@Q))-6];
590      if ( @s!=0 )
591      {
592         @o = ordstr(@Q);
593         if ( @o[1]=="c" or @o[1]=="C") { @o = @o[1,size(@o)-6]; }
594         else { @o = @o[1,size(@o)-8] + @o[size(@o)-1,2]; }
595      }
596      execute("ring @P=("+charstr(@Q)+"),("+@va+"),("+@o+");");
597      def @id = imap(@Q,@id);
598   }
599   def `@s1` = @P;
600   def `@s2` = @id;
601   attrib(`@s2`,"isSB",1);
602   export(`@s2`);
603   kill @P;
604   keepring(basering);
605   if( voice==2 ) { "// basering is now "+@s1+", std-basis has name "+@s2; }
606   return();
607}
608example
609{ "EXAMPLE:"; echo = 2;
610   ring s = 0,(e,f,x,y,z,t,u,v,w,a,b,c,d),(c,lp);
611   ideal i = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
612            d2+e2-1, f4+2u, wa+tf, xy+tu+ab;
613  option(prot); timer=1;
614   int time = timer;
615   ideal j=std(i);
616   timer-time;
617   show(R);dim(j),mult(j);
618   int time = timer;
619   faststd(i,"R","i");                      // use "best" ordering of vars
620   timer-time;
621   show(R);dim(i),mult(i);
622   setring s;time = timer;
623   faststd(i,"R","i","hilb");                // hilb-std only
624   timer-time;
625   show(R);dim(i),mult(i);
626   setring s;time = timer;
627   faststd(i,"R","i","hilb","sort");         // hilb-std,"best" ordering
628   timer-time;
629   show(R);dim(i),mult(i);
630    setring s;time = timer;
631   faststd(i,"R","i","hilb","sort","block","dec"); // hilb-std,"best",blocks
632   timer-time;
633   show(R);dim(i),mult(i);
634  setring s;time = timer;
635   timer-time;time = timer;
636   faststd(i,"R","i","sort","block","Dp"); //"best",decreasing,Dp-blocks
637   timer-time;
638   show(R);dim(i),mult(i);
639}
640///////////////////////////////////////////////////////////////////////////////
641
642proc findvars(id, list #)
643"USAGE:   findvars(id[,any]); id poly/ideal/vector/module/matrix, any=any type
644RETURN:  ideal of variables occuring in id, if no second argument is present
645         list of 4 objects, if a second argument is given (of any type)
646         -[1]: ideal of variables occuring in id
647         -[2]: intvec of variables occuring in id
648         -[3]: ideal of variables not occuring in id
649         -[4]: intvec of variables not occuring in id
650EXAMPLE: example findvars; shows an example
651"
652{
653   int ii,n;
654   ideal found, notfound;
655   intvec f,nf;
656   n = nvars(basering);
657   ideal i = simplify(ideal(matrix(id)),10);
658   matrix M[ncols(i)][1] = i;
659   vector v = module(M)[1];
660   ideal max = maxideal(1);
661
662   for (ii=1; ii<=n; ii++)
663   {
664      if ( v != subst(v,var(ii),0) )
665      {
666         found = found+var(ii);
667         f = f,ii;
668      }
669      else
670      {
671         notfound = notfound+var(ii);
672         nf = nf,ii;
673      }
674   }
675   if ( size(f)>1 ) { f = f[2..size(f)]; }      //intvec of found vars
676   if ( size(nf)>1 ) { nf = nf[2..size(nf)]; }  //intvec of vars not found
677   if( size(#)==0 )  { return(found); }
678   if( size(#)!=0 )  { list L = found,f,notfound,nf; return(L); }
679}
680example
681{ "EXAMPLE:"; echo = 2;
682   ring s  = 0,(e,f,x,y,t,u,v,w,a,d),dp;
683   ideal i = w2+f2-1, x2+t2+a2-1;
684   findvars(i);
685   findvars(i,1);
686}
687///////////////////////////////////////////////////////////////////////////////
688
689proc hilbvec (@id, list #)
690"USAGE:   hilbvec(id[,c,o]); id poly/ideal/vector/module/matrix, c,o=strings
691         c=char, o=ord in which hilb is computed (default: c=\"32003\", o=\"dp\")
692RETURN:  intvec of 1-st Hilbert-series of id, computed in char c and ordering o
693         bei ** aendern falls ringmaps vollstaendig ?
694NOTE:    id must be homogeneous (all vars having weight 1)
695EXAMPLE: example hilbvec; shows an example
696"
697{
698   def @P = basering;
699   string @c,@o = "32003", "dp";
700   if ( size(#) == 1 ) {  @c = #[1]; }
701   if ( size(#) == 2 ) {  @c = #[1]; @o = #[2]; }
702   string @si = typeof(@id)+" @i = "+string(@id)+";";  //** weg
703   execute("ring @r=("+@c+"),("+varstr(basering)+"),("+@o+");");
704   //**def i = imap(P,@id);
705   execute(@si);                   //** weg
706   //show(basering);
707   @i = std(@i);
708   intvec @hi = hilb(@i,1);         // intvec of 1-st Hilbert-series of id
709   return(@hi);
710}
711example
712{ "EXAMPLE:"; echo = 2;
713   ring s   = 0,(e,f,x,y,z,t,u,v,w,a,b,c,d,H),dp;
714   ideal id = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
715              d2+e2-1, f4+2u, wa+tf, xy+tu+ab;
716   id = homog(id,H);
717   hilbvec(id);
718}
719///////////////////////////////////////////////////////////////////////////////
720
721proc linearpart (id)
722"USAGE:   linearpart(id);  id=ideal/module
723RETURN:  generators of id of total degree <= 1
724EXAMPLE: example linearpart; shows an example
725"
726{
727   return(degreepart(id,0,1));
728}
729example
730{ "EXAMPLE:"; echo = 2;
731   ring r=0,(x,y,z),dp;
732   ideal i=1+x+x2+x3,3,x+3y+5z;
733   linearpart(i);
734   module m=[x,y,z],x*[x3,y2,z],[1,x2,z3,0,1];
735   show(linearpart(m));
736}
737///////////////////////////////////////////////////////////////////////////////
738
739proc tolessvars (id ,list #)
740"USAGE:   tolessvars(id,[s1,s2]); id poly/ideal/vector/module/matrix,
741         s1,s2=strings (names of: new ring, new ordering)
742CREATE:  nothing, if id contains all vars of the basering. Else, create
743         a ring with same char as the basering, but with less variables
744         (only those variables which actually occur in id) and map id to the
745         new ring, which will be the basering after the proc has finished.
746         The name of the new ring is by default R(n), where n is the number of
747         variables in the new ring. If, say, s1 = \"newR\" then the new ring has
748         name newR. In s2 a different ordering for the new ring may be given
749         as an allowed ordstring (default is \"dp\" resp. \"ds\", depending whether
750         the first block of the old ordering is a p- resp. an s-ordering).
751DISPLAY: If printlevel >=0, display ideal of vars which have been ommitted from
752         the old ring (default)
753RETURN:  the original ideal id
754NOTE:    You must not type, say, 'ideal id=tolessvars(id);' since the ring
755         to which 'id' would belong will only be defined by the r.h.s.. But you
756         may type 'def id=tolessvars(id);' or 'list id=tolessvars(id);'
757         since then 'id' does not a priory belong to a ring, its type will
758         be defined by the right hand side. Moreover, do not use a name which
759         occurs in the old ring, for the same reason.
760EXAMPLE: example tolessvars; shows an example
761"
762{
763//---------------- initialisation and check occurence of vars -----------------
764   int s,ii,n,fp,fs;
765   string s1,s2,newvar;
766   int pr = printlevel-voice+3;  // p = printlevel+1 (default: p=1)
767   def P = basering;
768   s2 = ordstr(P);
769
770   list L = findvars(id,1);
771   newvar = string(L[1]);    // string of new variables
772   n = size(L[1]);           // number of new variables
773   if( n == 0 )
774   {
775      dbprint( pr,"","// no variable occured in "+typeof(id)+", no change of ring!");
776      return(id);
777   }
778   if( n == nvars(P) )
779   {
780      dbprint( pr,"","// all variables occured in "+typeof(id)+", no change of ring!");
781      return(id);
782   }
783//----------------- prepare new ring, map to it and return --------------------
784   s1 = "R("+string(n)+")";
785   if ( size(#) == 0 )
786   {
787       fp = find(s2,"p");
788       fs = find(s2,"s");
789       if( fs==0 or (fs>=fp && fp!=0) ) { s2="dp"; }
790       else {  s2="ds"; }
791   }
792   if ( size(#) ==1 ) { s1=#[1]; }
793   if ( size(#) ==2 ) { s1=#[1]; s2=#[2]; }
794   //dbprint( pr,"","// variables which did not occur:",simplify(max,2) );
795   dbprint( pr,"","// variables which did not occur:",L[3] );
796
797   execute("ring "+s1+"=("+charstr(P)+"),("+newvar+"),("+s2+");");
798   def id = imap(P,id);
799   export(basering);
800   keepring (basering);
801   dbprint( pr,"// basering is now "+s1 );
802   return(id);
803}
804example
805{ "EXAMPLE:"; echo = 2;
806   ring r  = 0,(x,y,z),dp;
807   ideal i = y2-x3,x-3,y-2x;
808   def j   = tolessvars(i);
809   show(basering);
810   j;
811   setring r;
812   list j = tolessvars(i,"R_r","lp");
813   R_r;
814   kill R_r, R(2);
815}
816///////////////////////////////////////////////////////////////////////////////
817
818proc solvelinearpart (id,list #)
819"USAGE:   solvelinearpart(id[,n]);  id=ideal/module, n=integer
820RETURN:  (interreduced) generators of id of degree <=1 in reduced triangular
821         form  (default) or if n=0 [non-reduced triangular form if n!=0]
822ASSUME:  monomial ordering is a global ordering (p-ordering)
823NOTE:    may be used to solve a system of linear equations
824         see proc 'gauss_row' from 'matrix.lib' for a different method
825WARNING: the result is very likely to be false for 'real' coefficients, use
826         char 0 instead!
827EXAMPLE: example solvelinearpart; shows an example
828"
829{
830   intvec getoption = option(get);
831   option(redSB);
832   if ( size(#)!=0 )
833   {
834      if(#[1]!=0) { option(noredSB); }
835   }
836   def lin = interred(degreepart(id,0,1));
837   if ( size(#)!=0 )
838   {
839      if(#[1]!=0)
840      {
841         return(lin);
842      }
843   }
844   option(set,getoption);
845   return(simplify(lin,1));
846}
847example
848{ "EXAMPLE:"; echo = 2;
849   // Solve the system of linear equations:
850   //         3x +   y +  z -  u = 2
851   //         3x +  8y + 6z - 7u = 1
852   //        14x + 10y + 6z - 7u = 0
853   //         7x +  4y + 3z - 3u = 3
854   ring r = 0,(x,y,z,u),lp;
855   ideal i= 3x +   y +  z -  u,
856           13x +  8y + 6z - 7u,
857           14x + 10y + 6z - 7u,
858            7x +  4y + 3z - 3u;
859   ideal j= 2,1,0,3;
860   j = i-j;                        // difference of 1x4 matrices
861                                   // compute reduced triangular form, setting
862   solvelinearpart(j);             // the RHS equal 0 gives the solutions!
863   solvelinearpart(j,1); "";       // triangular form, not reduced
864   // Solve the same system simultaneously for two RHS's: 2,1,0,3 and 1,2,3,4
865   matrix b[2][size(i)]=2,1,0,3,1,2,3,4;
866   module m = i*[1,1]-b;
867   show(solvelinearpart(m));"";
868   // Solve the same system but with parametric values for the RHS:
869   ring r1 = (0,a,b,c,d),(x,y,z,u),dp;
870   ideal i = 3x +   y +  z -  u - a,
871            13x +  8y + 6z - 7u - b,
872            14x + 10y + 6z - 7u - c,
873             7x +  4y + 3z - 3u - d;
874   solvelinearpart(i);
875}
876///////////////////////////////////////////////////////////////////////////////
877
878proc sortandmap (@id,string @s1,string @s2, list #)
879"USAGE:   sortandmap(id,s1,s2[,n1,p1,n2,p2...,o1,m1,o2,m2...]);
880         id=poly/ideal/vector/module
881         s1,s2=strings (names for new ring and maped id)
882         p1,p2,...= product of vars, n1,n2,...=integers
883         o1,o2,...= allowed ordstrings, m1,m2,...=integers
884         (default: p1=product of all vars, n1=0, o1=\"dp\",m1=0)
885         the last pi (containing the remaining vars) may be omitted
886CREATE:  a new ring and map id into it, the new ring has same char as basering
887         but with new ordering and vars sorted in the following manner:
888         - each block of vars occuring in pi is sorted w.r.t its complexity in
889           id, ni controls the sorting in i-th block (= vars occuring in pi):
890           ni=0 (resp.!=0) means that less (resp. more) complex vars come first
891         - If say, s1=\"R\" and s2=\"j\", the new basering has name R and the image
892           of id in R has name j
893         - oi and mi define the monomial ordering of the i-th block:
894           if mi =0, oi=ordstr(i-th block)
895           if mi!=0, the ordering of the i-th block itself is a blockordering,
896           each subblock having ordstr=oi, such that vars of same complexity
897           are in one block
898           default: oi=\"dp\", mi=0
899         - only simple ordstrings oi are allowed:\"lp\",\"dp\",\"Dp\",\"ls\",\"ds\",\"Ds\"
900RETURN:  nothing
901NOTE:    We define a variable x to be more complex than y (with respect to id)
902         if val(x) > val(y) lexicographically, where val(x) denotes the
903         valuation vector of x: consider id as list of polynomials in x with
904         coefficients in the remaining variables. Then val(x) =
905         (maximal occuring power of x, # of monomials in leading coefficient,
906          # of monomials in coefficient of next smaller power of x,...)
907EXAMPLE: example sortandmap; shows an example
908"
909{
910   def @P = basering;
911   int @ii,@jj;
912   intvec @v;
913   string @o;
914//----------------- find o in # and split # into 2 lists ---------------------
915   # = # +list("dp",0);
916   for ( @ii=1; @ii<=size(#); @ii++)
917   {
918      if ( typeof(#[@ii])=="string" )  break;
919   }
920   if ( @ii==1 ) { list @L1 = list(); }
921   else { list @L1 = #[1..@ii-1]; }
922   list @L2 = #[@ii..size(#)];
923   list @L = sortvars(@id,@L1);
924   string @va = string(@L[1]);
925   list @l = @L[2];   //e.g. @l[4]=intvec describing permutation of 1-st block
926//----------------- construct correct ordering with oi and mi ----------------
927   for ( @ii=4; @ii<=size(@l); @ii=@ii+4 )
928   {
929      @L2=@L2+list("dp",0);
930      if ( @L2[@ii/2] != 0)
931      {
932         @v = @l[@ii];
933         for ( @jj=1; @jj<=size(@v); @jj++ )
934         {
935           @o = @o+@L2[@ii/2 -1]+"("+string(@v[@jj])+"),";
936         }
937      }
938      else
939      {
940         @o = @o+@L2[@ii/2 -1]+"("+string(size(@l[@ii/2]))+"),";
941      }
942   }
943   @o=@o[1..size(@o)-1];
944//------------------ create new ring and make objects global -----------------
945   execute("ring "+@s1+"=("+charstr(@P)+"),("+@va+"),("+@o+");");
946   def @id = imap(@P,@id);
947   execute("def "+ @s2+"=@id;");
948   execute("export("+@s1+");");
949   execute("export("+@s2+");");
950   keepring(basering);
951   return();
952}
953example
954{ "EXAMPLE:"; echo = 2;
955   ring s = 32003,(e,f,x,y,z,t,u,v,w,a,b,c,d),dp;
956   ideal i = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
957             d2+e2-1, f4+2u, wa+tf, xy+tu+ab, yz+uv+bc,
958             cd+ze, x+y+z+e+1, t+u+v+f-1, w+a+b+c+d;
959   sortandmap(i,"R_r","i");
960   // i is now an ideal in the new basering R_r
961   show(R_r);
962   kill R_r; setring s;
963   sortandmap(i,"R_r","i",1,"lp",0);
964   show(R_r);
965   kill R_r; setring s;
966   sortandmap(i,"R_r","i",1,abc,0,xyztuvw,0,"lp",0,"Dp",1);
967   show(R_r);
968   kill R_r;
969}
970///////////////////////////////////////////////////////////////////////////////
971
972proc sortvars (id, list #)
973"USAGE:   sortvars(id[,n1,p1,n2,p2,...]); id=poly/ideal/vector/module,
974         p1,p2,...= product of vars, n1,n2,...=integers
975         (default: p1=product of all vars, n1=0)
976         the last pi (containing the remaining vars) may be omitted
977COMPUTE: sort variables with respect to their complexity in id
978RETURN:  list of two elements, an ideal and a list:
979         [1]: ideal, variables of basering sorted w.r.t their complexity in id
980            ni controls the ordering in i-th block (= vars occuring in pi):
981            ni=0 (resp.!=0) means that less (resp. more) complex vars come
982            first
983         [2]: a list with 4 elements for each pi:
984            ideal ai : vars of pi in correct order,
985            intvec vi: permutation vector describing the ordering in ai,
986            intmat Mi: valuation matrix of ai, the columns of Mi being the
987                       valuation vectors of the vars in ai
988            intvec wi: 1-st,2-nd,...entry = size of 1-st,2-nd,... block of
989                       identically columns of Mi (vars with same valuation)
990NOTE:    We define a variable x to be more complex than y (w.r.t. id)
991         if val(x) > val(y) lexicographically, where val(x) denotes the
992         valuation vector of x: consider id as list of polynomials in x with
993         coefficients in the remaining variables. Then val(x) =
994         (maximal occuring power of x, # of monomials in leading coefficient,
995          # of monomials in coefficient of next smaller power of x,...)
996EXAMPLE: example sortvars; shows an example
997"
998{
999   int ii,jj,n,s;
1000   list L = valvars(id,#);
1001   list L2, L3 = L[2], L[3];
1002   list K; intmat M; intvec v1,v2,w;
1003   ideal i = sort(maxideal(1),L[1])[1];
1004   for ( ii=1; ii<=size(L2); ii++ )
1005   {
1006      M = transpose(L3[2*ii]);
1007      M = M[L2[ii],1..nrows(L3[2*ii])];
1008      w = 0; s = 0;
1009      for ( jj=1; jj<=nrows(M)-1; jj++ )
1010      {
1011         v1 = M[jj,1..ncols(M)];
1012         v2 = M[jj+1,1..ncols(M)];
1013         if ( v1 != v2 ) { n=jj-s; s=s+n; w = w,n; }
1014      }
1015      w=w,nrows(M)-s; w=w[2..size(w)];
1016      K = K+sort(L3[2*ii-1],L2[ii])+list(transpose(M))+list(w);
1017   }
1018   L = i,K;
1019   return(L);
1020}
1021example
1022{ "EXAMPLE:"; echo = 2;
1023   ring r=0,(a,b,c,x,y,z),lp;
1024   poly f=a3+b4+xyz2+xyz+yz+1;
1025   show(sortvars( f,1,abc,1)[1]);"";
1026   ring s=31991,(e,f,x,y,z,t,u,v,w,a,b,c,d),dp;
1027   ideal i = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
1028             d2+e2-1, f4+2u, wa+tf, xy+tu+ab, yz+uv+bc,
1029             cd+ze, x+y+z+e+1, t+u+v+f-1, w+a+b+c+d;
1030   show(sortvars(i,1,uybcazt,0,fewdvx));
1031}
1032///////////////////////////////////////////////////////////////////////////////
1033
1034proc valvars (id, list #)
1035"USAGE:   valvars(id[,n1,p1,n2,p2,...]); id=poly/ideal/vector/module,
1036         p1,p2,...= product of vars, n1,n2,...=integers
1037         ni controls the ordering of vars occuring in pi:
1038         ni=0 (resp.!=0) means that less (resp. more) complex vars come first
1039         (default: p1=product of all vars, n1=0)
1040         the last pi (containing the remaining vars) may be omitted
1041COMPUTE: valuation (complexity) of variables with respect to id
1042RETURN:  list consisting of 3 objects:
1043         [1]: intvec, say v, describing the permutation such that the permuted
1044            ringvariables are ordered with respect to their complexity in id
1045         [2]: list of intvecs, i-th intvec, say v(i) describing prmutation
1046              of vars in a(i) such that v=v(1),v(2),...
1047         [3]: list of ideals and intmat's, say a(i) and M(i), where ideal a(i)
1048            = factors of pi, M(i) = valuation matrix of a(i), such that the
1049            j-th column of M(i) is the valuation vector of j-th generator of a(i)
1050NOTE:    Use proc 'sortvars' for the actual sorting of vars!
1051         We define a variable x to be more complex than y (with respect to id)
1052         if val(x) > val(y) lexicographically, where val(x) denotes the
1053         valuation vector of x: consider id as list of polynomials in x with
1054         coefficients in the remaining variables. Then val(x) =
1055         (maximal occuring power of x, # of monomials in leading coefficient,
1056          # of monomials in coefficient of next smaller power of x,...)
1057EXAMPLE: example valvars; shows an example
1058"
1059{
1060//---------------------------- initialization ---------------------------------
1061   int ii,jj,kk,n;
1062   list L;                    // list of valuation vectors in one block
1063   intvec vec;                // describes permutation of vars (in one block)
1064   list blockvec;             // i-th element = vec of i-th block
1065   intvec varvec;             // result intvector
1066   list Li;                   // result list of ideals
1067   list LM;                   // result list of intmat's
1068   intvec v,w,s;              // w valuation vector for one variable
1069   matrix C;                  // coefficient matrix for different variables
1070   ideal i = simplify(ideal(matrix(id)),10);
1071
1072//---- for each pii in # create ideal a(ii) intvec v(ii) and list L(ii) -------
1073// a(ii) = ideal of vars in product, v(ii)[j]=k <=> a(ii)[j]=var(k)
1074
1075   v = 1..nvars(basering);
1076   int l = size(#);
1077   if ( l >= 2 )
1078   {
1079      ideal m=maxideal(1);
1080      for ( ii=2; ii<=l; ii=ii+2 )
1081      {
1082         int n(ii) = #[ii-1];
1083         ideal a(ii);
1084         intvec v(ii);
1085         for ( jj=1; jj<=nvars(basering); jj++ )
1086         {
1087            if ( #[ii]/var(jj) != 0)
1088            {
1089               a(ii) = a(ii) + var(jj);
1090               v(ii)=v(ii),jj;
1091               m[jj]=0;
1092               v[jj]=0;
1093            }
1094         }
1095         v(ii)=v(ii)[2..size(v(ii))];
1096      }
1097      if ( size(m)!=0 )
1098      {
1099         l = 2*(l/2)+2;
1100         ideal a(l) = simplify(m,2);
1101         intvec v(l) = compress(v);
1102         int n(l);
1103         if ( size(#)==l-1 ) { n(l) = #[l-1]; }
1104      }
1105   }
1106   else
1107   {
1108      l = 2;
1109      ideal a(2) = maxideal(1);
1110      intvec v(2) = v;
1111      int n(2);
1112      if ( size(#)==1 ) { n(2) = #[1]; }
1113   }
1114//------------- start loop to order variables in each a(ii) -------------------
1115
1116   for ( kk=2; kk<=l; kk=kk+2 )
1117   {
1118      L = list();
1119      n = 0;
1120//---------------- get valuation of all variables in a(kk) --------------------
1121      for ( ii=1; ii<=size(a(kk)); ii++ )
1122      {
1123         C = coeffs(i,a(kk)[ii]);
1124         w = nrows(C);
1125         for ( jj=w[1]; jj>1; jj-- )
1126         {
1127            s = size(C[jj,1..ncols(C)]);
1128            w[w[1]-jj+2] = sum(s);
1129         }
1130         L[ii]=w;
1131         n = size(w)*(size(w) > n) + n*(size(w) <= n);
1132      }
1133      intmat M(kk)[size(a(kk))][n];
1134      for ( ii=1; ii<=size(a(kk)); ii++ )
1135      {
1136         if ( n==1 ) { w = L[ii]; M(kk)[ii,1] = w[1]; }
1137         else  { M(kk)[ii,1..n] = L[ii]; }
1138      }
1139      LM[kk-1] = a(kk);
1140      LM[kk] = transpose(compress(M(kk)));
1141//------------------- compare valuation and insert in vec ---------------------
1142      vec = sort(L)[2];
1143      if ( n(kk) != 0 ) { vec = vec[size(vec)..1]; }
1144      blockvec[kk/2] = vec;
1145      vec = sort(v(kk),vec)[1];
1146      varvec = varvec,vec;
1147   }
1148   varvec = varvec[2..size(varvec)];
1149   list result = varvec,blockvec,LM;
1150   return(result);
1151}
1152example
1153{ "EXAMPLE:"; echo = 2;
1154   ring r=0,(a,b,c,x,y,z),lp;
1155   poly f=a3+b4+xyz2+xyz+yz+1;
1156   show(valvars(f,1,abc)[1]);"";
1157   ring s=31991,(e,f,x,y,z,t,u,v,w,a,b,c,d),dp;
1158   ideal i = w2+f2-1, x2+t2+a2-1,  y2+u2+b2-1, z2+v2+c2-1,
1159             d2+e2-1, f4+2u, wa+tf, xy+tu+ab, yz+uv+bc,
1160             cd+ze, x+y+z+e+1, t+u+v+f-1, w+a+b+c+d;
1161   list v6=valvars(i,1,uybcazt,0,efwdvx);
1162   show(v6);
1163}
1164///////////////////////////////////////////////////////////////////////////////
1165/*
1166
1167 ring s=31991,(e,f,x,y,z,t,u,v,w,a,b,c,d),dp;
1168 ring s=31991,(x,y,z,t,u,v,w,a,b,c,d,f,e,h),dp; //standard
1169 ring s1=31991,(y,u,b,c,a,z,t,x,v,d,w,e,f,h),dp; //gut
1170v;
117113,12,11,10,8,7,6,5,4,3,2,1,9,14
1172print(matrix(sort(maxideal(1),v)));
1173f,e,w,d,x,t,z,a,c,b,u,y,v,h
1174print(matrix(maxideal(1)));
1175y,u,b,c,a,z,t,x,v,d,w,e,f,h
1176v0;
117714,9,12,11,10,8,7,6,5,4,3,2,1,13
1178print(matrix(sort(maxideal(1),v0)));
1179h,v,e,w,d,x,t,z,a,c,b,u,y,f
1180v1;v2;
11819,12,11,10,8,7,6,5,4,3,2,1,13,14
118213,12,11,10,8,7,6,5,4,3,2,1,9,14
1183
1184Ev. Gute Ordnung fuer i:
1185========================
1186i=ad*x^d+ad-1*x^(d-1)+...+a1*x+a0, ad!=0
1187mit ar=(ar1,...,ark), k=size(i)
1188    arj in K[..x^..]
1189d=deg_x(i) := max{deg_x(i[k]) | k=1..size(i)}
1190size_x(i,deg_x(i)..0) := size(ad),...,size(a0)
1191x>y  <==
1192  1. deg_x(i)>deg_y(i)
1193  2. "=" in 1. und size_x lexikographisch
1194
1195hier im Beispiel:
1196f: 5,1,0,1,2
1197
1198u: 3,1,4
1199
1200y: 3,1,3
1201b: 3,1,3
1202c: 3,1,3
1203a: 3,1,3
1204z: 3,1,3
1205t: 3,1,3
1206
1207x: 3,1,2
1208v: 3,1,2
1209d: 3,1,2
1210w: 3,1,2
1211e: 3,1,2
1212probier mal:
1213 ring s=31991,(f,u,y,z,t,a,b,c,v,w,d,e,h),dp; //standard
1214
1215*/
Note: See TracBrowser for help on using the repository browser.