# Changeset e90fd1 in git

Ignore:
Timestamp:
Jan 30, 2015, 6:00:54 PM (9 years ago)
Branches:
Children:
10e299defc5ae0c54c970239c42f261ddbfaa5d4
Parents:
c9fcd5e3b1b695444030076b43c81e910fb047db
Message:
new version: grobcov.lib (grobcovK)
Files:
2 edited

### Legend:

Unmodified
 rc9fcd5 ////////////////////////////////////////////////////////////////////////////// version="version grobcov.lib 4-0-0-0 Dec_2013 "; // version="version grobcov.lib 4.0.1.2 Jan_2015 "; // $Id$ category="General purpose"; info=" LIBRARY:  grobcov.lib   Groebner Cover for parametric ideals. PURPOSE:  Comprehensive Groebner Systems, Groebner Cover, Canonical Forms, Parametric Polynomial Systems. The library contains Montes-Wibmer's algorithms to compute the canonical Groebner cover of a parametric ideal as described in the paper: Montes A., Wibmer M., \"Groebner Bases for Polynomial Systems with parameters\". Journal of Symbolic Computation 45 (2010) 1391-1425. The locus algorithm and definitions will be published in a forthcoming paper by Abanades, Botana, Montes, Recio: \''An Algebraic Taxonomy for Locus Computation in Dynamic Geometry\''. The central routine is grobcov. Given a parametric ideal, grobcov outputs its Canonical Groebner Cover, consisting of a set of pairs of (basis, segment). The basis (after normalization) is the reduced Groebner basis for each point of the segment. The segments are disjoint, locally closed and correspond to constant lpp (leading power product) of the basis, and are represented in canonical prime representation. The segments are disjoint and cover the whole parameter space. The output is canonical, it only depends on the given parametric ideal and the monomial order. This is much more than a simple Comprehensive Groebner System. The algorithm grobcov allows options to solve partially the problem when the whole automatic algorithm does not finish in reasonable time. grobcov uses a first algorithm cgsdr that outputs a disjoint reduced Comprehensive Groebner System with constant lpp. For this purpose, in this library, the implemented algorithm is Kapur-Sun-Wang algorithm, because it is the most efficient algorithm known for this purpose. D. Kapur, Y. Sun, and D.K. Wang. \"A New Algorithm for Computing Comprehensive Groebner Systems\". Proceedings of ISSAC'2010, ACM Press, (2010), 29-36. cgsdr can be called directly if only a disjoint reduced Comprehensive Groebner System (CGS) is required. AUTHORS:  Antonio Montes , Hans Schoenemann. OVERVIEW: see \"Groebner Bases for Polynomial Systems with parameters\" Montes A., Wibmer M., Comprehensive Groebner Systems, Groebner Cover, Canonical Forms, Parametric Polynomial Systems, Dynamic Geometry, Loci, Envelop, Constructible sets. See A. Montes A, M. Wibmer, \"Groebner Bases for Polynomial Systems with parameters\", Journal of Symbolic Computation 45 (2010) 1391-1425. (http://www-ma2.upc.edu/~montes/). AUTHORS:  Antonio Montes , Hans Schoenemann. OVERVIEW: In 2010, the library was designed to contain Montes-Wibmer's algorithms for compute the canonical Groebner Cover of a parametric ideal as described in the paper: Montes A., Wibmer M., \"Groebner Bases for Polynomial Systems with parameters\". Journal of Symbolic Computation 45 (2010) 1391-1425. The central routine is grobcov. Given a parametric ideal, grobcov outputs its Canonical Groebner Cover, consisting of a set of pairs of (basis, segment). The basis (after normalization) is the reduced Groebner basis for each point of the segment. The segments are disjoint, locally closed and correspond to constant lpp (leading power product) of the basis, and are represented in canonical prime representation. The segments are disjoint and cover the whole parameter space. The output is canonical, it only depends on the given parametric ideal and the monomial order. This is much more than a simple Comprehensive Groebner System. The algorithm grobcov allows options to solve partially the problem when the whole automatic algorithm does not finish in reasonable time. grobcov uses a first algorithm cgsdr that outputs a disjoint reduced Comprehensive Groebner System with constant lpp. For this purpose, in this library, the implemented algorithm is Kapur-Sun-Wang algorithm, because it is the most efficient algorithm known for this purpose. D. Kapur, Y. Sun, and D.K. Wang. \"A New Algorithm for Computing Comprehensive Groebner Systems\". Proceedings of ISSAC'2010, ACM Press, (2010), 29-36. The library has evolved to include new applications of the Groebner Cover, and new theoretical developments have been done. A new set of routines (locus, locusdg, locusto) has been included to compute loci of points. The routines are used in the Dynamic Geometry software Geogebra. They are described in: Abanades, Botana, Montes, Recio: \''An Algebraic Taxonomy for Locus Computation in Dynamic Geometry\''. Computer-Aided Design 56 (2014) 22-33. Recently also routines for computing the envelop of a family of curves (enverlop, envelopdg), to be used in Dynamic Geometry, has been included and will be described in a forthcoming paper: Abanades, Botana, Montes, Recio: \''Envelops in Dynamic Geometry using the Groebner cover\''. The actual version also includes two routines (AddCons and AddconsP) for computing the canonical form of a constructible set, given as a union of locally closed sets. They are used in the new version for the computation of loci and envelops. It determines the canonical locally closed level sets of a constructuble. They will be described in a forthcoming paper: A. Montes, J.M. Brunat, \"Canonical representations of constructible sets\". This version was finished on 31/01/2015 NOTATIONS: All given and determined polynomials and ideals are in the @*         If you want to use some internal routine you must @*         create before the above rings by calling setglobalrings(); @*         because most of the internal routines use these rings. @*         because some of the internal routines use these rings. @*         The call to the basic routines grobcov, cgsdr will @*         kill these rings. PROCEDURES: grobcov(F);  Is the basic routine giving the canonical Groebner cover of the parametric ideal F. This routine accepts many options, that allow to obtain results even when the canonical computation does not finish in reasonable time. cgsdr(F);      Is the procedure for obtaining a first disjoint, reduced Comprehensive Groebner System that is used in grobcov, that can also be used independently if only the CGS is required. It is a more efficient routine than buildtree (the own routine that is no more used). It uses now KSW algorithm. setglobalrings();  Generates the global rings @R, @P and @PR that are respectively the rings Q[a][x], Q[a], Q[x,a]. It is called inside each of the fundamental routines of the library: grobcov, cgsdr, locus, locusdg and killed before the output. So, if the user want to use some other internal routine, then setglobalrings() is to be called before, as most of them use these rings. pdivi(f,F);     Performs a pseudodivision of a parametric polynomial by a parametric ideal. Can be used without calling presiouly setglobalrings(), Groebner Cover of the parametric ideal F. This routine accepts many options, that allow to obtain results even when the canonical computation does not finish in reasonable time. cgsdr(F); Is the procedure for obtaining a first disjoint, reduced Comprehensive Groebner System that is used in grobcov, that can also be used independently if only the CGS is required. It is a more efficient routine than buildtree (the own routine of 2010 that is no more used). Now, KSW algorithm is used. setglobalrings();  Generates the global rings @R, @P and @PR that are respectively the rings Q[a][x], Q[a], Q[x,a].  (a=parameters, x=variables) It is called inside each of the fundamental routines of the library: grobcov, cgsdr, locus, locusdg and killed before the output. If the user want to use some other internal routine, then setglobalrings() is to be called before, as most of them use these rings. pdivi(f,F); Performs a pseudodivision of a parametric polynomial by a parametric ideal. pnormalf(f,E,N);   Reduces a parametric polynomial f over V(E) \ V(N) E is the null ideal and N the non-null ideal over the parameters. Can be used without calling presiouly setglobalrings(), Prep(N,M);   Computes the P-representation of V(N) \ V(M). Can be used without calling previously setglobalrings(). extend(GC); When the grobcov of an ideal has been computed with the default option ('ext',0) and the explicit option ('rep',2) (which is not the default), then one can call extend (GC) (and options) to obtain the full representation of the bases. With the default option ('ext',0) only the generic representation of the bases are computed, and one can obtain the full representation using extend. Can be used without calling presiouly setglobalrings(), locus(G):     Special routine for determining the locus of points of  objects. Given a parametric ideal J with parameters (a_1,..a_m) and variables (x_1,..,xn), representing the system determining the locus of points (a_1,..,a_m)) who verify certain properties, computing the grobcov G of J and applying to it locus, determines the different classes of locus components. They can be Normal, Special, Accumulation point, Degenerate. The output are the components given in P-canonical form of two constructible sets: Normal, and Non-Normal The Normal Set has Normal and Special components The Non-Normal set has Accumulation and Degenerate components. The description of the algorithm and definitions will be given in a forthcoming paper by Abanades, Botana, Montes, Recio: ''An Algebraic Taxonomy for Locus Computation in Dynamic Geometry''. Can be used without calling presiouly setglobalrings(), locusdg(G):  Is a special routine for computing the locus in dinamical geometry. It detects automatically a possible point that is to be avoided by the mover, whose coordinates must be the last two coordinates in the definition of the ring. If such a point is detected, then it eliminates the segments of the grobcov depending on the point that is to be avoided. Then it calls locus. Can be used without calling presiouly setglobalrings(), locusto(L):  Transforms the output of locus to a string that can be reed from different computational systems. Can be used without calling presiouly setglobalrings(), addcons(L): Given a disjoint set of locally closed subsets in P-representation, it returns the canonical P-representation of the constructible set formed by the union of them. Can be used without calling presiouly setglobalrings(), ( E is the null ideal and N the non-null ideal ) over the parameters. Crep(N,M);  Computes the canonical C-representation of V(N) \ V(M). Prep(N,M);  Computes the canonical P-representation of V(N) \ V(M). PtoCrep(L)  Starting from the canonical Prep of a locally closed set computes its Crep. extend(GC);  When the grobcov of an ideal has been computed with the default option ('ext',0) and the explicit option ('rep',2) (which is not the default), then one can call extend (GC) (and options) to obtain the full representation of the bases. With the default option ('ext',0) only the generic representation of the bases are computed, and one can obtain the full representation using extend. locus(G);    Special routine for determining the geometrical locus of points verifying given conditions. Given a parametric ideal J with parameters (x,y) and variables (x_1,..,xn), representing the system determining the locus of points (x,y) who verify certain properties, one can apply locus to the output of  grobcov(J), locus determines the different classes of locus components. described in the paper: \"An Algebraic Taxonomy for Locus Computation in Dynamic Geometry\", M. Abanades, F. Botana, A. Montes, T. Recio, Computer-Aided Design 56 (2014) 22-33. The components can be 'Normal', 'Special', 'Accumulation', 'Degenerate'. The output are the components is given in P-canonical form It also detects automatically a possible point that is to be avoided by the mover, whose coordinates must be the last two coordinates in the definition of the ring. If such a point is detected, then it eliminates the segments of the grobcov depending on the point that is to be avoided. locusdg(G);  Is a special routine that determines the  'Relevant' components of the locus in dynamic geometry. It is to be called to the output of locus and selects from it the useful components. envelop(F,C); Special routine for determining the envelop of a family of curves F in Q[x,y][x_1,..xn] depending on a ideal of constraints C in Q[x_1,..,xn]. It detemines the different components as well as its type: 'Normal', 'Special', 'Accumulation', 'Degenerate'. And it also classifies the Special components, determining the zero dimensional antiimage of the component and verifying if the component is a special curve of the family or not. It calls internally first grobcov and then locus with special options to obtain the complete result. The taxonomy that it provides, as well as the algorithms involved will be described in a forthcoming paper: Abanades, Botana, Montes, Recio: \''Envelops in Dynamic Geometry using the Gr\"obner cover\''. envelopdg(ev); Is a special routine to determine the 'Relevant' components of the envelop of a family of curves to be used in Dynamic Geometry. It must be called to the output of envelop(F,C). locusto(L); Transforms the output of locus, locusdg, envelop and  envelopdg into a string that can be reed from different computational systems. AddCons(L); Uses the routine AddConsP. Given a set of locally closed sets as difference of of varieties (does not need to be in C-representation) it returns the canonical P-representation of the constructible set formed by the union of them. The result is formed by a set of embedded, disjoint, locally closed sets (levels). AddConsP(L);  Given a set of locally closed P-components, it returns the canonical P-representation of the constructible set formed by the union of them, consisting of a set of embedded, disjoint locally closed sets (levels). The routines AddCons and AddConsP and the canonical structure of the constructible sets will be described in a forthcoming paper. A. Montes, J.M. Brunat, \"Canonical representations of constructible sets\". SEE ALSO: compregb_lib // Library grobcov.lib // (Groebner cover): // Release 1: (public) // (Groebner Cover): // Release 0: (public) // Initial data: 21-1-2008 // Uses buildtree for cgsdr // Final data: 3-7-2008 // Release 2: (private) // Initial data: 6-9-2009 // Last version using buildtree for cgsdr // Final data: 25-10-2011 // Release 3: // Release B: (private) // Initial data: 1-7-2012 // Uses both buildtree and KSW for cgsdr // Final data: 4-9-2012 // Release 4: (this release, public) // Release G: (public) // Initial data: 4-9-2012 // Uses KSW algorithm for cgsdr // Final data: 21-11-2013 // basering Q[a][x]; // ************ Begin of buildtree ****************************** proc setglobalrings() "USAGE:   setglobalrings(); No arguments RETURN:   After its call the rings @R=Q[a][x], @P=Q[a], @RP=Q[x,a] are defined as global variables. NOTE: It is called internally by the fundamental routines of the library grobcov, cgsdr, extend, pdivi, pnormalf, locus, locusdg, locusto, and killed before the output. The user does not need to call it, except when it is interested in using some internal routine of the library that uses these rings. The basering R, must be of the form Q[a][x], a=parameters, x=variables, and should be defined previously. KEYWORDS: ring, rings EXAMPLE:  setglobalrings; shows an example" { if (defined(@P)) { kill @P; kill @R; kill @RP; } // Release K: (public) // Updated locus: 8-1-2015 // Updated AddConsP and AddCons: 8-1-2015 // Reformed many routines, examples and helps: 8-1-2015 // New routines for computing the envelop of a family of curves: 22-1-2015 // Final data: 22-1-2015 //*************Auxiliary routines************** // elimintfromideal: elimine the constant numbers from an ideal //        (designed for W, nonnull conditions) // Input: ideal J // Output:ideal K with the elements of J that are non constants, in the //        ring K[x1,..,xm] static proc elimintfromideal(ideal J) { int i; int j=0; ideal K; if (size(J)==0){return(ideal(0));} for (i=1;i<=ncols(J);i++){if (size(variables(J[i])) !=0){j++; K[j]=J[i];}} return(K); } // delfromideal: deletes the i-th polynomial from the ideal F //    Works in any kind of ideal static proc delfromideal(ideal F, int i) { int j; ideal G; if (size(F)0) { i=1; while ((t) and (i<=size(L))) { if (equalideals(L[i],M[i])==0){t=0;} i++; } } return(t); } } // idcontains // Input: ideal p, ideal q // Output: 1 if p contains q,  0 otherwise // If the routine is to be called from the top, a previous call to // setglobalrings() is needed. static proc idcontains(ideal p, ideal q) { int t; int i; t=1; i=1; def P=p; def Q=q; attrib(P,"isSB",1); poly r; while ((t) and (i<=size(Q))) { r=reduce(Q[i],P); if (r!=0){t=0;} i++; } return(t); } // selectminideals //   given a list of ideals returns the list of integers corresponding //   to the minimal ideals in the list // Input: L (list of ideals) // Output: the list of integers corresponding to the minimal ideals in L //   Works in Q[u_1,..,u_m] static proc selectminideals(list L) { list P; int i; int j; int t; if(size(L)==0){return(L)}; if(size(L)==1){P[1]=1; return(P);} for (i=1;i<=size(L);i++) { t=1; j=1; while ((t) and (j<=size(L))) { if (i!=j) { if(idcontains(L[i],L[j])==1) { t=0; } } j++; } if (t){P[size(P)+1]=i;} } return(P); } // Auxiliary routine // elimconstfac: eliminate the factors in the polynom f that are in Q[a] // Input: //   poly f: //   list L: of components of the segment // Output: //   poly f2  where the factors of f in Q[a] that are non-null on any component //   have been dropped from f static proc elimconstfac(poly f) { int cond; int i; int j; int t; if (f==0){return(f);} def RR=basering; def @R=basering;  // must be of the form K[a][x], a=parameters, x=variables list Rx=ringlist(RR); def @P=ring(Rx[1]); list Lx; Lx[1]=0; Lx[2]=Rx[2]+Rx[1][2]; Lx[3]=Rx[1][3]; Lx[4]=Rx[1][4]; Rx[1]=0; def D=ring(Rx); def @RP=D+@P; export(@R);      // global ring K[a][x] export(@P);      // global ring K[a] export(@RP);     // global ring K[x,a] with product order setring(@R); def ff=imap(RR,f); def l=factorize(ff,0); poly f1=1; for(i=2;i<=size(l[1]);i++) { f1=f1*(l[1][i])^(l[2][i]); } setring(RR); } example { "EXAMPLE:"; echo = 2; ring R=(0,a,b),(x,y,z),dp; setglobalrings(); Grobcov::@R; Grobcov::@P; Grobcov::@RP; ringlist(Grobcov::@R); ringlist(Grobcov::@P); ringlist(Grobcov::@RP); } //*************Auxiliary routines************** // cld : clears denominators of an ideal and normalizes to content 1 //        can be used in @R or @P or @RP // input: //        ideal J (J can be also poly), but the output is an ideal; // output: //        ideal Jc (the new form of ideal J without denominators and //        normalized to content 1) proc cld(ideal J) { if (size(J)==0){return(ideal(0));} int te=0; def RR=basering; if(not(defined(@RP))) { te=1; setglobalrings(); } setring(@RP); def Ja=imap(RR,J); ideal Jb; if (size(Ja)==0){setring(RR); return(ideal(0));} int i; def j=0; for (i=1;i<=ncols(Ja);i++){if (size(Ja[i])!=0){j++; Jb[j]=cleardenom(Ja[i]);}} setring(RR); def Jc=imap(@RP,Jb); if(te){kill @R; kill @RP; kill @P;} return(Jc); } proc memberpos(f,J) def f2=imap(@R,f1); return(f2); }; static proc memberpos(f,J) //"USAGE:  memberpos(f,J); //         (f,J) expected (polynomial,ideal) //               or       (ideal,list(ideal)) //               or       (list(intvec),  list(list(intvec))). //         The ring can be @R or @P or @RP or any other. //         Works in any kind of ideals //RETURN:  The list (t,pos) t int; pos int; //         t is 1 if f belongs to J and 0 if not. //} proc subset(J,K) // Auxiliary routine // pos // Input:  intvec p of zeros and ones // Output: intvec W of the positions where p has ones. static proc pos(intvec p) { int i; intvec W; int j=1; for (i=1; i<=size(p); i++) { if (p[i]==1){W[j]=i; j++;} } return(W); } // Input: //  A,B: lists of ideals // Output: //   1 if both lists of ideals are equal, or 0 if not static proc equallistsofideals(list A, list B) { int i; int tes=0; if (size(A)!=size(B)){return(tes);} tes=1; i=1; while(tes==1 and i<=size(A)) { if (equalideals(A[i],B[i])==0){tes=0; return(tes);} i++; } return(tes); } // Input: //  A,B:  lists of P-rep, i.e. of the form [p_i,[p_{i1},..,p_{ij_i}]] // Output: //   1 if both lists of P-reps are equal, or 0 if not static proc equallistsA(list A, list B) { int tes=0; if(equalideals(A[1],B[1])==0){return(tes);} tes=equallistsofideals(A[2],B[2]); return(tes); } // Input: //  A,B:  lists lists of of P-rep, i.e. of the form [[p_1,[p_{11},..,p_{1j_1}]] .. [p_i,[p_{i1},..,p_{ij_i}]] // Output: //   1 if both lists of lists of P-rep are equal, or 0 if not static proc equallistsAall(list A,list B) { int i; int tes; if(size(A)!=size(B)){return(tes);} tes=1; i=1; while(tes and i<=size(A)) { tes=equallistsA(A[i],B[i]); i++; } return(tes); } // idint: ideal intersection //        in the ring @P. //        it works in an extended ring // input: two ideals in the ring @P // output the intersection of both (is not a GB) static proc idint(ideal I, ideal J) { def RR=basering; ring T=0,t,lp; def K=T+RR; setring(K); def Ia=imap(RR,I); def Ja=imap(RR,J); ideal IJ; int i; for(i=1;i<=size(Ia);i++){IJ[i]=t*Ia[i];} for(i=1;i<=size(Ja);i++){IJ[size(Ia)+i]=(1-t)*Ja[i];} ideal eIJ=eliminate(IJ,t); setring(RR); return(imap(K,eIJ)); } //purpose ideal intersection called in @R and computed in @P static proc idintR(ideal N, ideal M) { def RR=basering; setring(@P); def Np=imap(RR,N); def Mp=imap(RR,M); def Jp=idint(Np,Mp); setring(RR); return(imap(@P,Jp)); } // Auxiliary routine // comb: the list of combinations of elements (1,..n) of order p static proc comb(int n, int p) { list L; list L0; intvec c; intvec d; int i; int j; int last; if ((n<0) or (n=g proc lesspol(poly f, poly g) static proc lesspol(poly f, poly g) { if (leadmonom(f) proc incquotient(ideal N, poly f) static proc incquotient(ideal N, poly f) { poly g; int i; } // eliminates the ith element from a list or an intvec proc elimfromlist(def l, int i) { if(typeof(l)=="list"){list L;} if (typeof(l)=="intvec"){intvec L;} if (typeof(l)=="ideal"){ideal L;} int j; if((size(l)==0) or (size(l)==1 and i!=1)){return(l);} if (size(l)==1 and i==1){return(L);} // L=l[1]; if(i==1) { for(j=2;j<=size(l);j++) { L[j-1]=l[j]; } } else { for(j=1;j<=i-1;j++) { L[j]=l[j]; } for(j=i+1;j<=size(l);j++) { L[j-1]=l[j]; } } return(L); } // Auxiliary routine to define an order for ideals // Returns 1 if the ideal a is shoud precede ideal b by sorting them in idbefid order //             2 if the the contrary happen //             0 if the order is not relevant proc idbefid(ideal a, ideal b) static proc idbefid(ideal a, ideal b) { poly fa; poly fb; poly la; poly lb; // sort a list of ideals using idbefid proc sortlistideals(list L) static proc sortlistideals(list L) { int i; int j; int n; } // returns 1 if the two lists of ideals are equal and 0 if not proc equallistideals(list L, list M) { int t; int i; if (size(L)!=size(M)){return(0);} else { t=1; if (size(L)>0) { i=1; while ((t) and (i<=size(L))) { if (equalideals(L[i],M[i])==0){t=0;} i++; } } return(t); } // Crep // Computes the C-representation of V(N) \ V(M). // input: //    ideal N (null ideal) (not necessarily radical nor maximal) //    ideal M (hole ideal) (not necessarily containing N) // output: //    the list (a,b) of the canonical ideals //    the Crep of V(N) \ V(M) // Assumed to be called in the ring @R // Works on the ring @P proc Crep(ideal N, ideal M) "USAGE:  Crep(N,M); Input: ideal N (null ideal) (not necessarily radical nor maximal) ideal M (hole ideal) (not necessarily containing N) RETURN:  The canonical C-representation of the locally closed set. [ P,Q ], a pair of radical ideals with P included in Q, representing the set V(P) \ V(Q) = V(N) \ V(M) NOTE:   Operates in a ring R=Q[a] (a=parameters) KEYWORDS: locally closed set, canoncial form EXAMPLE:  Crep; shows an example" { list l; ideal Np=std(N); if (equalideals(Np,ideal(1))) { l=ideal(1),ideal(1); return(l); } int i; list L; ideal Q=Np+M; ideal P=ideal(1); L=minGTZ(Np); for(i=1;i<=size(L);i++) { if(idcontains(L[i],Q)==0) { P=intersect(P,L[i]); } } P=std(P); Q=std(radical(Q+P)); list T=P,Q; return(T); } example { "EXAMPLE:"; echo = 2; if(defined(Grobcov::@P)){kill Grobcov::@R; kill Grobcov::@P; kill Grobcov::@RP;} ring R=0,(x,y,z),lp; short=0; ideal E=x*(x^2+y^2+z^2-25); ideal N=x*(x-3),y-4; def Cr=Crep(E,N); Cr; def L=Prep(E,N); L; def Cr1=PtoCrep(L); Cr1; } // output: //    the ((p_1,(p_11,p_1k_1)),..,(p_r,(p_r1,p_rk_r))); //    the Prep of V(N)\V(M) // Assumed to work in the ring @P of the parameters // Can be used without calling previously setglobalrings(). //    the Prep of V(N) \ V(M) // Assumed to be called in the ring @R // Works on the ring @P proc Prep(ideal N, ideal M) "USAGE:  Prep(N,M); Input: ideal N (null ideal) (not necessarily radical nor maximal) ideal M (hole ideal) (not necessarily containing N) RETURN: The canonical P-representation of the locally closed set V(N) \ V(M) Output: [ Comp_1, .. , Comp_s ] where Comp_i=[p_i,[p_i1,..,p_is_i]] NOTE:   Operates in a ring R=Q[a]  (a=parameters) KEYWORDS: Locally closed set, Canoncial form EXAMPLE:  Prep; shows an example" { int te; return(list(list(ideal(1),list(ideal(1))))); } def RR=basering; if(defined(@P)){te=1;} else{te=0; setglobalrings();} setring(@P); ideal Np=imap(RR,N); ideal Mp=imap(RR,M); int i; int j; list L0; list Ni=minGTZ(Np); list Ni=minGTZ(N); list prep; for(j=1;j<=size(Ni);j++) for (i=1;i<=size(Ni);i++) { Mij=minGTZ(Ni[i]+Mp); Mij=minGTZ(Ni[i]+M); for(j=1;j<=size(Mij);j++) { } if (size(prep)==0){prep=list(list(ideal(1),list(ideal(1))));} setring(RR); def pprep=imap(@P,prep); if(te==0){kill @P; kill @R; kill @RP;} return(pprep); def Lout=CompleteA(prep,prep); return(Lout); } example { "EXAMPLE:"; echo = 2; if(defined(Grobcov::@P)){kill Grobcov::@R; kill Grobcov::@P; kill Grobcov::@RP;} short=0; ring R=0,(x,y,z),lp; ideal E=x*(x^2+y^2+z^2-25); ideal N=x*(x-3),y-4; Prep(E,N); } // input: //    list ((p_1,(p_11,p_1k_1)),..,(p_r,(p_r1,p_rk_r))); //         the P-representation of V(N)\V(M) //         the P-representation of V(N) \ V(M) // output: //    list (ideal ida, ideal idb) //    the C-representaion of V(N)\V(M) = V(ida)\V(idb) // Assumed to work in the ring @P of the parameters // If the routine is to be called from the top, a previous call to // setglobalrings() is needed. //    the C-representaion of V(N) \ V(M) = V(ida) \ V(idb) // Assumed to be called in the ring @R // Works on the ring @P proc PtoCrep(list L) "USAGE: PtoCrep(L) Input L: list  [ Comp_1, .. , Comp_s ] where Comp_i=[p_i,[p_i1,..,p_is_i] ], is the P-representation of a locally closed set V(N) \ V(M) RETURN:  The canonical C-representation of the locally closed set [ P,Q ], a pair of radical ideals with P included in Q, representing the set V(P) \ V(Q) NOTE: Operates in a ring R=Q[a]  (a=parameters) KEYWORDS: locally closed set, canoncial form EXAMPLE:  PtoCrep; shows an example" { int te; def RR=basering; if(defined(@P)){te=1; setring(@P); list Lp=imap(RR,L);} else {  te=0; def Lp=L;} def La=PtoCrep0(Lp); if(te==1) {setring(RR); def LL=imap(@P,La);} if(te==0){def LL=La;} return(LL); } example { "EXAMPLE:"; echo = 2; if(defined(Grobcov::@P)){kill Grobcov::@R; kill Grobcov::@P; kill Grobcov::@RP;} short=0; ring R=0,(x,y,z),lp; // (P,Q) represents a locally closed set ideal P=x^3+x*y^2+x*z^2-25*x; ideal Q=y-4,x*z,x^2-3*x; // Now compute the P-representation= def L=Prep(P,Q); L; // Now compute the C-representation= def J=PtoCrep(L); J; // Now come back recomputing the P-represetation of the C-representation= Prep(J[1],J[2]); } // PtoCrep0 // Computes the C-representation from the P-representation. // input: //    list ((p_1,(p_11,p_1k_1)),..,(p_r,(p_r1,p_rk_r))); //         the P-representation of V(N) \ V(M) // output: //    list (ideal ida, ideal idb) //    the C-representation of V(N) \ V(M) = V(ida) \ V(idb) // Works in a ring Q[u_1,..,u_m] and is called on it. static proc PtoCrep0(list L) { int te=0; def RR=basering; if(defined(@P)){te=1;} else{te=0; setglobalrings();} setring(@P); def Lp=imap(RR,L); def Lp=L; int i; int j; ideal ida=ideal(1); ideal idb=ideal(1); list Lb; ideal N; } } //idb=radical(idb); def La=list(ida,idb); setring(RR); def LL=imap(@P,La); if(te==0){kill @P; kill @R; kill @RP;} return(LL); return(La); } //      and dehomogenize the result. proc cgsdr(ideal F, list #) "USAGE: cgsdr(F); To compute a disjoint, reduced CGS. "USAGE: cgsdr(F); F: ideal in Q[a][x] (a=parameters, x=variables) to be discussed. To compute a disjoint, reduced Comprehensive Groebner System (CGS). cgsdr is the starting point of the fundamental routine grobcov. Inside grobcov it is used only with options 'can' set to 0,1 and Inside grobcov it is used with options 'can' set to 0,1 and not with options ('can',2). It is to be used if only a disjoint reduced CGS is required. F: ideal in Q[a][x] (parameters and variables) to be discussed. Options: To modify the default options, pairs of arguments \"null\",ideal E: The default is (\"null\",ideal(0)). \"nonnull\",ideal N: The default (\"nonnull\",ideal(1)). When options 'null' and/or 'nonnull' are given, then the parameter space is restricted to V(E)\V(N). When options \"null\" and/or \"nonnull\" are given, then the parameter space is restricted to V(E) \ V(N). \"comment\",0-1: The default is (\"comment\",0). Setting (\"comment\",1) will provide information about the development of the and the segments grouped by lpp With options (\"can\",0) and (\"can\",1) the option (\"out\",1) is set to (out,0) because it is not compatible. is set to (\"out\",0) because it is not compatible. One can give none or whatever of these options. With the default options (\"can\",2,\"out\",1), only the Kapur-Sun-Wang algorithm is computed. This is very effectif but is only the starting point for the computation. Kapur-Sun-Wang algorithm is computed. This is very efficient but is only the starting point for the computation of grobcov. When grobcov is computed, the call to cgsdr inside uses specific options that are more expensive ("can",0-1,"out",0). RETURN:   Returns a list T describing a reduced and disjoint RETURN: Returns a list T describing a reduced and disjoint Comprehensive Groebner System (CGS), With option (\"out\",0) the segments are grouped by leading power products (lpp) of the reduced Groebner basis and given in P-representation. The returned list is of the form: the segments are grouped by leading power products (lpp) of the reduced Groebner basis and given in P-representation. The returned list is of the form: ( (lpp, (num,basis,segment),...,(num,basis,segment),lpp), (lpp, (num,basis,segment),...,(num,basis,segment),lpp) ) The bases are the reduced Groebner bases (after normalization) for each point of the corresponding segment. The third element of each lpp segment is the lpp of the used ideal in the CGS as a string: The bases are the reduced Groebner bases (after normalization) for each point of the corresponding segment. The third element of each lpp segment is the lpp of the used ideal in the CGS as a string: with option (\"can\",0) the homogenized basis is used with option (\"can\",1) the homogenized ideal is used With option (\"out\",1) (default) only KSW is applied and segments are given as difference of varieties and are not grouped The returned list is of the form: only KSW is applied and segments are given as difference of varieties and are not grouped The returned list is of the form: ( (E,N,B),..(E,N,B) ) E is the null variety N is the nonnull variety segment = V(E)\V(N) B is the reduced Groebner basis NOTE:     The basering R, must be of the form Q[a][x], a=parameters, x=variables, and should be defined previously, and the ideal E is the null variety N is the nonnull variety segment = V(E) \ V(N) B is the reduced Groebner basis NOTE:  The basering R, must be of the form Q[a][x], (a=parameters, x=variables), and should be defined previously, and the ideal defined on R. KEYWORDS: CGS, disjoint, reduced, Comprehensive Groebner System // INITIALIZING OPTIONS int i; int j; def E=ideal(0); def N=ideal(1); int comment=0; int can=2; int out=1; poly f; ideal B; def E=ideal(0); def N=ideal(1); int comment=0; int start=timer; list L=#; { // COMPUTING THE HOMOGOENIZED IDEAL if(comment>0){string("Homogenizing the whole ideal: option can=1");} if(comment>=1){string("Homogenizing the whole ideal: option can=1");} list RRL=ringlist(RR); RRL[3][1][1]="dp"; BH[i]=homog(BH[i],@t); } if (comment>=1){string("Homogenized system = "); BH;} if (comment>=2){string("Homogenized system = "); BH;} def GSH=KSW(BH,LH); setglobalrings(); GS[i][3]=postredgb(mingb(GS[i][3])); if (typeof(GS[i][7])=="ideal") { GS[i][7]=postredgb(mingb(GS[i][7])); } { GS[i][7]=postredgb(mingb(GS[i][7]));} } } } example { "EXAMPLE:"; echo = 2; "Casas conjecture for degree 4"; { "EXAMPLE:"; echo = 2; // Casas conjecture for degree 4: ring R=(0,a0,a1,a2,a3,a4),(x1,x2,x3),dp; short=0; ideal F=x1^4+(4*a3)*x1^3+(6*a2)*x1^2+(4*a1)*x1+(a0), x1^3+(3*a3)*x1^2+(3*a2)*x1+(a1), //            lpp segments and improve the output // output: grouped segments by lpp obtained in cgsdr proc grsegments(list T) static proc grsegments(list T) { int i; } // idcontains // input: ideal p, ideal q // output: 1 if p contains q,  0 otherwise // If the routine is to be called from the top, a previous call to // setglobalrings() is needed. proc idcontains(ideal p, ideal q) { int t; int i; t=1; i=1; def RR=basering; setring(@P); def P=imap(RR,p); def Q=imap(RR,q); attrib(P,"isSB",1); poly r; while ((t) and (i<=size(Q))) { r=reduce(Q[i],P); if (r!=0){t=0;} i++; } setring(RR); return(t); } // selectminindeals //   given a list of ideals returns the list of integers corresponding //   to the minimal ideals in the list // input: L (list of ideals) // output: the list of integers corresponding to the minimal ideals in L // If the routine is to be called from the top, a previous call to // setglobalrings() is needed. proc selectminideals(list L) { if (size(L)==0){return(L)}; def RR=basering; setring(@P); def Lp=imap(RR,L); int i; int j; int t; intvec notsel; list P; for (i=1;i<=size(Lp);i++) { if(memberpos(i,notsel)[1]) { i++; if(i>size(Lp)){break;} } t=1; j=1; while ((t) and (j<=size(Lp))) { if (i==j){j++;} if ((j<=size(Lp)) and (memberpos(j,notsel)[1]==0)) { if (idcontains(Lp[i],Lp[j])) { notsel[size(notsel)+1]=i; t=0; } } j++; } if (t){P[size(P)+1]=i;} } setring(RR); return(P); } // LCUnion // Given a list of the P-representations of disjoint locally closed segments // Given a list of the P-representations of locally closed segments // for which we know that the union is also locally closed // it returns the P-representation of its union // output: P-representation of the union //       ((P_j,(P_j1,...,P_jk_j | j=1..t))) // If the routine is to be called from the top, a previous call to // setglobalrings() is needed. // Auxiliary routine called by grobcov and addcons // A previous call to setglobarings() is needed if it is to be used at the top. proc LCUnion(list LL) static proc LCUnion(list LL) { def RR=basering; setring(@P); int care; def L=imap(RR,LL); int i; int j; int k; list H; list C; list T; list L0; list P0; list P; list Q0; list Q;  intvec Qi; if(not(defined(@Q2))){list @Q2;} else{kill @Q2; list @Q2;} exportto(Top,@Q2); list L0; list P0; list P; list Q0; list Q; for (i=1;i<=size(L);i++) { } Q0=selectminideals(P0); //"T_3Q0="; Q0; kill P; list P; for (i=1;i<=size(Q0);i++) { Qi=L0[Q0[i]]; Q[size(Q)+1]=Qi; P[size(P)+1]=L[Qi[1]][Qi[2]]; } @Q2=Q; if(size(P)==0) { setring(RR); list TT; return(TT); Q[i]=L0[Q0[i]]; P[i]=L[Q[i][1]][Q[i][2]]; } // P is the list of the maximal components of the union { C[size(C)+1]=L[i][j]; C[size(C)][3]=intvec(i,j); } } } if(size(P[k])>=3){T[size(T)+1]=list(Q[k],P[k][1],addpart(H,C),P[k][3]);} else{T[size(T)+1]=list(Q[k],P[k][1],addpart(H,C));} } @Q2=sortpairs(@Q2); T[size(T)+1]=list(Q[k],P[k][1],addpart(H,C)); } setring(RR); def TT=imap(@P,T); // posi=(i,j) position of the component //        the list of segments to be added to the holes proc addpart(list H, list C) static proc addpart(list H, list C) { list Q; int i; int j; int k; int l; int t; int t1; if (equalideals(q,C[j][1])) { @Q2[size(@Q2)+1]=C[j][3]; // \\ @Q2[size(@Q2)+1]=C[j][3]; t=0; for (k=1;k<=size(C[j][2]);k++) { t1=1; //kill addq; //list addq; l=1; while((t1) and (l<=size(Q))) { addq[size(addq)+1]=C[j][2][k]; @Q2[size(@Q2)+1]=C[j][3]; // \\ @Q2[size(@Q2)+1]=C[j][3]; } } list addq; } //print("Q="); Q; print("notQ="); notQ; } i++; // that part. // Works on @P ring. proc addpartfine(list H, list C0) { static proc addpartfine(list H, list C0) { //"T_H="; H; int i; int j; int k; int te; intvec notQ; int l; list sel; int used; intvec jtesC; if ((size(H)==1) and (equalideals(H[1],ideal(1)))){return(H);} if (size(C0)==0){return(H);} def RR=basering; setring(@P); list newQ; list nQ; list Q; list nQ1; list Q0; def Q1=imap(RR,H); //Q1=sortlistideals(Q1); def C=imap(RR,C0); def Q1=H; //Q1=sortlistideals(Q1,idbefid); def C=C0; while(equallistideals(Q0,Q1)==0) { } } //"T_Q1_0="; Q1; sel=selectminideals(Q1); kill nQ1; list nQ1; Q1=nQ1; } setring(RR); if(size(Q1)==0){Q1=ideal(1),ideal(1);} //"T_Q1_1="; Q1; //if(used>0){string("addpartfine was ", used, " times used");} return(imap(@P,Q1)); } return(Q1); } // Auxiliary routine for grobcov: ideal F is assumed to be homogeneous //      where a Prep is ( (p1,(p11,..,p1k_1)),..,(pj,(pj1,..,p1k_j)) ) //            a Crep is ( ida, idb ) proc gcover(ideal F,list #) static proc gcover(ideal F,list #) { int i; int j; int k; ideal lpp; list GPi2; list pairspP; ideal B; int ti; // grobcov // input: //    ideal F: a parametric ideal in Q[a][x], where a are the parameters //             and x the variables //    ideal F: a parametric ideal in Q[a][x], (a=parameters, x=variables). //    list #: (options) list("null",N,"nonnull",W,"can",0-1,ext",0-1, "rep",0-1-2) //            where proc grobcov(ideal F,list #) "USAGE:   grobcov(F); This is the fundamental routine of the library. It computes the Groebner cover of a parametric ideal (see (*) Montes A., Wibmer M., Groebner Bases for Polynomial Systems with parameters. JSC 45 (2010) 1391-1425.) The Groebner cover of a parametric ideal consist of a set of library. It computes the Groebner cover of a parametric ideal. See Montes A., Wibmer M., \"Groebner Bases for Polynomial Systems with parameters\". JSC 45 (2010) 1391-1425.) The Groebner Cover of a parametric ideal consist of a set of pairs(S_i,B_i), where the S_i are disjoint locally closed segments of the parameter space, and the B_i are the reduced The ideal F must be defined on a parametric ring Q[a][x]. (a=parameters, x=variables) Options: To modify the default options, pair of arguments -option name, value- of valid options must be added to the call. Options: \"null\",ideal E: The default is (\"null\",ideal(0)). \"nonnull\",ideal N: The default (\"nonnull\",ideal(1)). \"nonnull\",ideal N: The default is (\"nonnull\",ideal(1)). When options \"null\" and/or \"nonnull\" are given, then the parameter space is restricted to V(E)\V(N). the parameter space is restricted to V(E) \ V(N). \"can\",0-1: The default is (\"can\",1). With the default option the homogenized ideal is computed before obtaining the Groebner cover, so that the result is the canonical Groebner cover. Setting (\"can\",0) only homogenizes the Groebner Cover, so that the result is the canonical Groebner Cover. Setting (\"can\",0) only homogenizes the basis so the result is not exactly canonical, but the computation is shorter. \"ext\",0-1: The default is (\"ext\",0). With the default (\"ext\",0), only the generic representation is computed (single polynomials, but not specializing to non-zero at each point of the segment. With option (\"ext\",1) the (\"ext\",0), only the generic representation of the bases is computed (single polynomials, but not specializing to non-zero for every point of the segment. With option (\"ext\",1) the full representation of the bases is computed (possible sheaves) and sometimes a simpler result is obtained. sheaves) and sometimes a simpler result is obtained, but the computation is more time consuming. \"rep\",0-1-2: The default is (\"rep\",0) and then the segments are given in canonical P-representation. Option (\"rep\",1) of the segment. The lpph corresponds to the lpp of the homogenized ideal and is different for each segment. It is given as a string. Basis: to each element of lpp corresponds an I-regular function given in full representation (by option (\"ext\",1)) or in and is different for each segment. It is given as a string, and shown only for information. With the default option \"can\",1, the segments have different lpph. Basis: to each element of lpp corresponds an I-regular function given in full representation (by option (\"ext\",1)) or in generic representation (default option (\"ext\",0)). The I-regular function is the corresponding element of the reduced The P-representation of a segment is of the form ((p_1,(p_11,..,p_1k1)),..,(p_r,(p_r1,..,p_rkr)) representing the segment U_i (V(p_i) \ U_j (V(p_ij))), representing the segment Union_i ( V(p_i) \ ( Union_j V(p_ij) ) ), where the p's are prime ideals. The C-representation of a segment is of the form (E,N) representing V(E)\V(N), and the ideals E and N are (E,N) representing V(E) \ V(N), and the ideals E and N are radical and N contains E. NOTE: The basering R, must be of the form Q[a][x], a=parameters, x=variables, and should be defined previously. The ideal must NOTE: The basering R, must be of the form Q[a][x], (a=parameters, x=variables), and should be defined previously. The ideal must be defined on R. KEYWORDS: Groebner cover, parametric ideal, canonical, discussion of } example { "EXAMPLE:"; echo = 2; "Casas conjecture for degree 4"; { "EXAMPLE:"; echo = 2; // Casas conjecture for degree 4: ring R=(0,a0,a1,a2,a3,a4),(x1,x2,x3),dp; short=0; ideal F=x1^4+(4*a3)*x1^3+(6*a2)*x1^2+(4*a1)*x1+(a0), x1^3+(3*a3)*x1^2+(3*a2)*x1+(a1), x2^4+(4*a3)*x2^3+(6*a2)*x2^2+(4*a1)*x2+(a0), x2^2+(2*a3)*x2+(a2), x3^4+(4*a3)*x3^3+(6*a2)*x3^2+(4*a1)*x3+(a0), x3+(a3); x1^3+(3*a3)*x1^2+(3*a2)*x1+(a1), x2^4+(4*a3)*x2^3+(6*a2)*x2^2+(4*a1)*x2+(a0), x2^2+(2*a3)*x2+(a2), x3^4+(4*a3)*x3^3+(6*a2)*x3^2+(4*a1)*x3+(a0), x3+(a3); grobcov(F); } // Can be called from the top proc extend(list GC, list #); "USAGE:   extend(GC); When the grobcov of an ideal has been computed with the default option (\"ext\",0) and the explicit option (\"rep\",2) (which is not the default), then one can call extend (GC) (and options) to obtain the full representation of the bases. With the default option (\"ext\",0) only the generic representation of the bases are computed, and one can obtain the full representation using extend. \"rep\",0-1-2: The default is (\"rep\",0) and then the segments are given in canonical P-representation. Option (\"rep\",1) represents the segments in canonical C-representation, and option (\"rep\",2) gives both representations. \"comment\",0-1: The default is (\"comment\",0). Setting \"comment\" higher will provide information about the time used in the computation. One can give none or whatever of these options. RETURN:   The list ( (lpp_1,basis_1,segment_1,lpph_1), ... (lpp_s,basis_s,segment_s,lpph_s) ) The lpp are constant over a segment and correspond to the set of lpp of the reduced Groebner basis for each point of the segment. The lpph corresponds to the lpp of the homogenized ideal and is different for each segment. It is given as a string. Basis: to each element of lpp corresponds an I-regular function given in full representation. The I-regular function is the corresponding element of the reduced Groebner basis for each point of the segment with the given lpp. For each point in the segment, the polynomial or the set of polynomials representing it, if they do not specialize to 0, then after normalization, specializes to the corresponding element of the reduced Groebner basis. In the full representation at least one of the polynomials representing the I-regular function specializes to non-zero. With the default option (\"rep\",0) the segments are given in P-representation. With option (\"rep\",1) the segments are given in C-representation. With option (\"rep\",2) both representations of the segments are given. The P-representation of a segment is of the form ((p_1,(p_11,..,p_1k1)),..,(p_r,(p_r1,..,p_rkr)) representing the segment U_i (V(p_i) \ U_j (V(p_ij))), where the p's are prime ideals. The C-representation of a segment is of the form (E,N) representing V(E)\V(N), and the ideals E and N are radical and N contains E. NOTE: The basering R, must be of the form Q[a][x], a=parameters, x=variables, and should be defined previously. The ideal must "USAGE: extend(GC); The default option of grobcov provides the bases in generic representation (the I-regular functions of the bases ara given by a single polynomial. It can specialize to zero for some points of the segments, but in general, it is sufficient for many pouposes. Nevertheless the I-regular functions allow a full representation given bey a set of polynomials specializing to the value of the function (after normalization) or to zero, but at least one of the polynomials specializes to non-zero. The full representation can be obtained by computing the grobcov with option \"ext\",1. The default option is \"ext\",0. With option \"ext\",1 the computation can be much more time consuming, even if the result can be simpler. Alternatively, one can compute the full representation of the bases after computing grobcov with the defaoult option \"ext\",0 and the option \"rep\",2, that outputs both the Prep and the Crep of the segments and then call \"extend\" to the output. RETURN: When calling extend(grobcov(S,\"rep\",2)) the result is of the form ( (lpp_1,basis_1,segment_1,lpph_1), ... (lpp_s,basis_s,segment_s,lpph_s) ) where each function of the basis can be given by an ideal of representants. NOTE: The basering R, must be of the form Q[a][x], (a=parameters, x=variables), and should be defined previously. The ideal must be defined on R. KEYWORDS: Groebner cover, parametric ideal, canonical, discussion of } if(m==1){"Warning! grobcov has already extended bases"; return(S);} if(size(GC[1])!=5){"Warning! extend make sense only when grobcov has been called with options 'rep',2,'ext',0"; " "; return();} if(size(GC[1])!=5){"Warning! extend make sense only when grobcov has been called with options  'rep',2,'ext',0"; " "; return();} int repop=0; int start3=timer; example { ring R=(0,a0,b0,c0,a1,b1,c1,a2,b2,c2),(x), dp; "EXAMPLE:"; echo = 2; ring R=(0,a0,b0,c0,a1,b1,c1),(x), dp; short=0; ideal S=a0*x^2+b0*x+c0, a1*x^2+b1*x+c1, a2*x^2+b2*x+c2; "System S="; S; def GCS=grobcov(S,"rep",2,"comment",1); "grobcov(S,'rep',2,'comment',1)="; GCS; def FGC=extend(GCS,"rep",0,"comment",1); "Full representation="; FGC; a1*x^2+b1*x+c1; def GCS=grobcov(S,"rep",2); GCS; def FGC=extend(GCS,"rep",0); // Full representation= FGC; } // nonzerodivisor // input: //    poly g in K[a], //    poly g in Q[a], //    list P=(p_1,..p_r) representing a minimal prime decomposition // output //    poly f such that f notin p_i for all i and //           g-f in p_i for all i such that g notin p_i proc nonzerodivisor(poly gr, list Pr) static proc nonzerodivisor(poly gr, list Pr) { def RR=basering; // input: //   int i: //   list LPr: (p1,..,pr) of prime components of an ideal in K[a] //   list LPr: (p1,..,pr) of prime components of an ideal in Q[a] // output: //   list (fr,fnr) of two polynomials that are equal on V(pi) //       and fr=0 on V(P) \ V(pi), and fnr is nonzero on V(pj) for all j. proc deltai(int i, list LPr) static proc deltai(int i, list LPr) { def RR=basering; def LP=imap(RR,LPr); int j; poly p; ideal F=ideal(1); def F=ideal(1); poly f; poly fn; // output: //    poly P on an open and dense set of V(p_1 int ... p_r) proc combine(list L, ideal F) static proc combine(list L, ideal F) { // ATTENTION REVISE AND USE Pci and F } // Auxiliary routine // elimconstfac: eliminate the factors in the polynom f that are in K[a] // input: //   poly f: //   list L: of components of the segment // output: //   poly f2  where the factors of f in K[a] that are non-null on any component //   have been dropped from f proc elimconstfac(poly f) { int cond; int i; int j; int t; if (f==0){return(f);} def RR=basering; setring(@R); def ff=imap(RR,f); list l=factorize(ff,0); poly f1=1; for(i=2;i<=size(l[1]);i++) { f1=f1*(l[1][i])^(l[2][i]); } setring(RR); def f2=imap(@R,f1); return(f2); } //Auxiliary routine // output: //   t:  with value 1 if f reduces modulo P, 0 if not. proc nullin(poly f,ideal P) static proc nullin(poly f,ideal P) { int t; // Input: A polynomial f // Output: The list of leading terms proc monoms(poly f) static proc monoms(poly f) { list L; //   poly f: a generic polynomial in the basis //   ideal idp: such that ideal(S)=idp //   ideal idq: such that S=V(idp)\V(idq) //   ideal idq: such that S=V(idp) \ V(idq) ////   NW the list of ((N1,W1),..,(Ns,Ws)) of red-rep of the grouped ////      segments in the lpp-segment  NO MORE USED // output: proc extend0(poly f, ideal idp, ideal idq) static proc extend0(poly f, ideal idp, ideal idq) { matrix CC; poly Q; list NewMonoms; //        each intvec v=(i_1,..,is) corresponds to a polynomial in the sheaf //        that will be built from it in extend procedure. proc findindexpolys(list coefs) static proc findindexpolys(list coefs) { int i; int j; intvec numdens; // Auxiliary routine // extendcoef: given Q,P in K[a] where P/Q specializes on an open and dense subset // extendcoef: given Q,P in Q[a] where P/Q specializes on an open and dense subset //      of the whole V(p1 int...int pr), it returns a basis of the module //      of all syzygies equivalent to P/Q, proc extendcoef(poly P, poly Q, ideal idp, ideal idq) static proc extendcoef(poly P, poly Q, ideal idp, ideal idq) { def RR=basering; // input: //   list L of the polynomials matrix CC //      (we assume that one of them is non-null on V(N)\V(M)) //   ideal N, ideal M: ideals representing the locally closed set V(N)\V(M) //      (we assume that one of them is non-null on V(N) \ V(M)) //   ideal N, ideal M: ideals representing the locally closed set V(N) \ V(M) // assume to work in @P proc selectregularfun(matrix CC, ideal NN, ideal MM) static proc selectregularfun(matrix CC, ideal NN, ideal MM) { int numcombused; // output: //   object T with index c proc searchinlist(intvec c,list L) static proc searchinlist(intvec c,list L) { int i; list T; } // Auxiliary routine // comb: the list of combinations of elements (1,..n) of order p proc comb(int n, int p) { list L; list L0; intvec c; intvec d; int i; int j; int last; if ((n<0) or (n0 encara que sigui 0 else { if(version==0) { //"T_B="; B;  // Computing std(B+E,plex(x,y,x1,..xn)) one can detect if there is a first part // independent of parameters giving the variables with dimension 0 dd=indepparameters(B); if (dd==1){d=0; L0=d,string(B),determineF(B,F,E);} else{d=1; L0=2,"Normal";} } else { def RH=ringlist(RR); //"T_RH="; RH; def H=RH; H[1]=0; H[2]=RH[1][2]+RH[2]; int n=size(H[2]); intvec ll; for(i=1;i<=n;i++) { ll[i]=1; } H[3][1][1]="lp"; H[3][1][2]=ll; def RRH=ring(H); setring(RRH); ideal BH=imap(RR,B); ideal EH=imap(RR,E); ideal NH=imap(RR,N); if(Fenv==1){poly FH=imap(RR,F);} for(i=1;i<=size(EH);i++){BH[size(BH)+1]=EH[i];} BH=std(BH);  // MOLT COSTOS!!! ideal G; ideal r; poly q; for(i=1;i<=size(BH);i++) { r=factorize(BH[i],1); q=1; for(j=1;j<=size(r);j++) { if((pdivi(r[j],NH)[1] != 0) or (equalideals(ideal(NH),ideal(1)))) { q=q*r[j]; } } if(q!=1){G[size(G)+1]=q;} } setring RR; def GG=imap(RRH,G); ideal GGG; if(defined(L0)){kill L0; list L0;} for(i=1;i<=size(GG);i++) { if(indepparameters(GG[i])){GGG[size(GGG)+1]=GG[i];} } GGG=std(GGG); ideal GLM; for(i=1;i<=size(GGG);i++) { GLM[i]=leadmonom(GGG[i]); } attrib(GLM,"IsSB",1); d=dim(std(GLM)); string antiim=string(GGG); L0=d,antiim; if(d==0) { //" ";string("Antiimage of Special component = ", GGG); if(Fenv==1) { L0[3]=determineF(GGG,F,E); } } else { L0[2]="Normal"; } } } return(L0); } static proc determineF(ideal A,poly F,ideal E) { int env; int i; def RR=basering; def RH=ringlist(RR); def H=RH; H[1]=0; H[2]=RH[1][2]+RH[2]; int n=size(H[2]); intvec ll; for(i=1;i<=n;i++) { ll[i]=1; } H[3][1][1]="lp"; H[3][1][2]=ll; def RRH=ring(H); //" ";string("Antiimage of Special component = ", GGG); setring(RRH); list LL; def AA=imap(RR,A); def FH=imap(RR,F); def EH=imap(RR,E); ideal M=std(AA+FH); def rh=reduce(EH,M); if(rh==0){env=1;} else{env=0;} setring RR; //L0[3]=env; return(env); } // DimPar // Auxilliary routine to NorSing determining the dimension of a parametric ideal // Does not use @P and define it directly because is assumes that //                             variables and parameters have been inverted static proc DimPar(ideal E) { def RRH=basering; def RHx=ringlist(RRH); def Prin=ring(RHx[1]); setring(Prin); def E2=std(imap(RRH,E)); def d=dim(E2); setring RRH; return(d); } // locus0(G): Private routine used by locus (the public routine), that //                builds the diferent components. // input:      The output G of the grobcov (in generic representation, which is the default option) //       Options: The algorithm allows the following options ar pair of arguments: //                "movdim", d  : by default movdim is 2 but it can be set to other values //                "version", v   :  There are two versions of the algorithm. ('version',1) is //                 a full algorithm that always distinguishes correctly between 'Normal' //                 and 'Special' components, whereas ('version',0) can decalre a component //                 as 'Normal' being really 'Special', but is more effective. By default ('version',1) //                 is used when the number of variables is less than 4 and 0 if not. //                 The user can force to use one or other version, but it is not recommended. //                 "system", ideal F: if the initial systrem is passed as an argument. This is actually not used. //                 "comments", c: by default it is 0, but it can be set to 1. //                 Usually locus problems have mover coordinates, variables and tracer coordinates. //                 The mover coordinates are to be placed as the last variables, and by default, //                 its number is 2. If one consider locus problems in higer dimensions, the number of //                 mover coordinates (placed as the last variables) is to be given as an option. // output: //         list, the canonical P-representation of the Normal and Non-Normal locus: //              If all levels of a class of locus are 1, then the set is locally closed. Otherwise the level //              gives the depth of the component. proc locus0(list GG, list #) { int t1=1; int t2=1; static proc locus0(list GG, list #) { int te=0; int t1=1; int t2=1; int i; def R=basering; //if(defined(@P)==1){te=1; kill @P; kill @R; kill @RP; } //setglobalrings(); // Options list DD=#; int moverdim=nvars(R); int version=0; int nv=nvars(R); if(nv<4){version=1;} int comment=0; ideal Fm; poly F; for(i=1;i<=(size(DD) div 2);i++) { if(DD[2*i-1]=="movdim"){moverdim=DD[2*i];} if(DD[2*i-1]=="version"){version=DD[2*i];} if(DD[2*i-1]=="system"){Fm=DD[2*i];} if(DD[2*i-1]=="comment"){comment=DD[2*i];} if(DD[2*i-1]=="family"){F=DD[2*i];} } list HHH; if (GG[1][1][1]==1 and GG[1][2][1]==1 and GG[1][3][1][1][1]==0 and GG[1][3][1][2][1]==1){return(HHH);} list Lop=#; if(size(Lop)>0){moverdim=Lop[1];} setglobalrings(); list G1; list G2; list G1; list G2; def G=GG; list Q1; list Q2; int i; int d; int j; int k; int d; int j; int k; t1=1; for(i=1;i<=size(G);i++) } else{list P1;} ideal BB; ideal BB; int dd; list NS; for(i=1;i<=size(P1);i++) { if (indepparameters(P1[i][3])==1){P1[i][3]="Special";} else{P1[i][3]="Normal";} NS=NorSing(P1[i][3],P1[i][1],P1[i][2][1],DD); dd=NS[1]; if(dd==0){P1[i][3]=NS;} //"Special"; else{P1[i][3]="Normal";} } list P2; for(i=1;i<=size(P2);i++){Q2[i]=l; Q2[i][1]=P2[i];} P2=Q2; setring @P; setring(@P); ideal J; if(t1==1) { def C1=imap(R,P1); def L1=addcons(C1); def L1=AddLocus(C1); } else{list C1; list L1; kill P1; list P1;} { def C2=imap(R,P2); def L2=addcons(C2); } def L2=AddLocus(C2); } else{list L2; list C2; kill P2; list P2;} for(i=1;i<=size(L2);i++) def L=imap(@P,LN); for(i=1;i<=size(L);i++){if(size(L[i][2])==0){L[i][2]=ideal(1);}} kill @R; kill @RP; kill @P; //if(te==0){kill @R; kill @RP; kill @P;} list LL; for(i=1;i<=size(L);i++) { LL[i]=l; LL[i]=elimfromlist(L[i],1); } return(LL); } // locus0dg(G): Private routine used by locusdg (the public routine), that //                builds the diferent components used in Dynamical Geometry // input:      The output G of the grobcov (in generic representation, which is the default option) // output: //         list, the canonical P-representation of the Relevant components of the locus in Dynamical Geometry, //              i.e. the Normal and Accumulation components. //              This routine is compemented by locusdg that calls it in order to eliminate problems //              with degenerate points of the mover. //         The output components are given as //              ((p1,(p11,..p1s_1),type_1,level_1),..,(pk,(pk1,..pks_k),type_k,level_k) //              whwre type_i is always "Relevant". //         The components are given in canonical P-representation of the subset. //              If all levels of a class of locus are 1, then the set is locally closed. Otherwise the level //              gives the depth of the component. proc locus0dg(list GG, list #) { int t1=1; int t2=1; int ojo; def R=basering; int moverdim=nvars(R); list HHH; if (GG[1][1][1]==1 and GG[1][2][1]==1 and GG[1][3][1][1][1]==0 and GG[1][3][1][2][1]==1){return(HHH);} list Lop=#; if(size(Lop)>0){moverdim=Lop[1];} setglobalrings(); list G1; list G2; def G=GG; list Q1; list Q2; int i; int d; int j; int k; t1=1; for(i=1;i<=size(G);i++) { attrib(G[i][1],"IsSB",1); d=locusdim(G[i][2],moverdim); if(d==0){G1[size(G1)+1]=G[i];} else { if(d>0){G2[size(G2)+1]=G[i];} } } if(size(G1)==0){t1=0;} if(size(G2)==0){t2=0;} setring(@RP); if(t1) { list G1RP=imap(R,G1); } else {list G1RP;} list P1RP; ideal B; for(i=1;i<=size(G1RP);i++) { kill B; ideal B; for(k=1;k<=size(G1RP[i][3]);k++) { attrib(G1RP[i][3][k][1],"IsSB",1); G1RP[i][3][k][1]=std(G1RP[i][3][k][1]); for(j=1;j<=size(G1RP[i][2]);j++) { B[j]=reduce(G1RP[i][2][j],G1RP[i][3][k][1]); } P1RP[size(P1RP)+1]=list(G1RP[i][3][k][1],G1RP[i][3][k][2],B); } } setring(R); ideal h; if(t1) { def P1=imap(@RP,P1RP); for(i=1;i<=size(P1);i++) { for(j=1;j<=size(P1[i][3]);j++) { h=factorize(P1[i][3][j],1); P1[i][3][j]=h[1]; for(k=2;k<=size(h);k++) { P1[i][3][j]=P1[i][3][j]*h[k]; } } } } else{list P1;} ideal BB; for(i=1;i<=size(P1);i++) { if (indepparameters(P1[i][3])==1){P1[i][3]="Special";} else{P1[i][3]="Relevant";} } list P2; for(i=1;i<=size(G2);i++) { for(k=1;k<=size(G2[i][3]);k++) { P2[size(P2)+1]=list(G2[i][3][k][1],G2[i][3][k][2]); } } list l; for(i=1;i<=size(P1);i++){Q1[i]=l; Q1[i][1]=P1[i];} P1=Q1; for(i=1;i<=size(P2);i++){Q2[i]=l; Q2[i][1]=P2[i];} P2=Q2; setring @P; ideal J; if(t1==1) { def C1=imap(R,P1); def L1=addcons(C1); } else{list C1; list L1; kill P1; list P1;} if(t2==1) { def C2=imap(R,P2); def L2=addcons(C2); } else{list L2; list C2; kill P2; list P2;} for(i=1;i<=size(L2);i++) { J=std(L2[i][2]); d=dim(J); // AQUI if(d==0) { L2[i][4]=string("Relevant",L2[i][4]); } else{L2[i][4]=string("Degenerate",L2[i][4]);} } list LN; list ll; list l0; if(t1==1) { for(i=1;i<=size(L1);i++) { if(L1[i][4]=="Relevant") { LN[size(LN)+1]=l0; LN[size(LN)][1]=L1[i]; } } } if(t2==1) { for(i=1;i<=size(L2);i++) { if(L2[i][4]=="Relevant") { LN[size(LN)+1]=l0; LN[size(LN)][1]=L2[i]; } } } list LNN; kill ll; list ll; list lll; list llll; for(i=1;i<=size(LN);i++) { LNN[size(LNN)+1]=l0; ll=LN[i][1]; lll[1]=ll[2]; lll[2]=ll[3]; lll[3]=ll[4]; LNN[size(LNN)][1]=lll; } kill LN; list LN; LN=addcons(LNN); if(size(LN)==0){ojo=1;} setring(R); if(ojo==1){list L;} else { def L=imap(@P,LN); for(i=1;i<=size(L);i++){if(size(L[i][2])==0){L[i][2]=ideal(1);}} } kill @R; kill @RP; kill @P; list LL; for(i=1;i<=size(L);i++) { LL[i]=l; LL[i]=elimfromlist(L[i],1); } return(LL); } if(typeof(L[i][4])=="list") {L[i][4][1]="Special";} l[1]=L[i][2]; l[2]=L[i][3]; l[3]=L[i][4]; l[4]=L[i][5]; L[i]=l; } return(L); } //  locus(G):  Special routine for determining the locus of points //                 of  geometrical constructions. Given a parametric ideal J with //                 parameters (a_1,..a_m) and variables (x_1,..,xn), //                 representing the system determining //                 the locus of points (a_1,..,a_m) who verify certain //                 properties, computing the grobcov G of //                 J and applying to it locus, determines the different //                 classes of locus components. The components can be //                 Normal, Special, Accumulation, Degenerate. //                 The output are the components given in P-canonical form //                 of at most 4 constructible sets: Normal, Special, Accumulation, //                 Degenerate. //                 The description of the algorithm and definitions is //                 given in a forthcoming paper by Abanades, Botana, Montes, Recio. //                 Usually locus problems have mover coordinates, variables and tracer coordinates. //                 The mover coordinates are to be placed as the last variables, and by default, //                 its number is 2. If one consider locus problems in higer dimensions, the number of //                 mover coordinates (placed as the last variables) is to be given as an option. // //                 of  geometrical constructions. //  input:      The output G of the grobcov (in generic representation, which is the default option for grobcov) //  output: //               gives the depth of the component of the constructible set. proc locus(list GG, list #) "USAGE:   locus(G); The input must be the grobcov  of a parametrical ideal RETURN:  The  locus. The output are the components of four constructible subsets of the locus of the parametrical system: Normal , Special, Accumulation and Degenerate. These components are given as a list of  (pi,(pi1,..pis_i),type_i,level_i) varying i, where the p's are prime ideals, the type can be: Normal, Special, Accumulation, Degenerate. NOTE:      It can only be called after computing the grobcov of the parametrical ideal in generic representation ('ext',0), which is the default. The basering R, must be of the form Q[a_1,..,a_m][x_1,..,x_n]. KEYWORDS: geometrical locus, locus, loci. EXAMPLE:  locus; shows an example" { "USAGE: locus(G) The input must be the grobcov  of a parametrical ideal in Q[a][x], (a=parameters, x=variables). In fact a must be the tracer coordinates and x the mover coordinates and other auxiliary ones. Special routine for determining the locus of points of  geometrical constructions. Given a parametric ideal J representing the system determining the locus of points (a) who verify certain properties, the call to locus on the output of grobcov( J ) determines the different classes of locus components, following the taxonomy defined in Abanades, Botana, Montes, Recio: \"An Algebraic Taxonomy for Locus Computation in Dynamic Geometry\". Computer-Aided Design 56 (2014) 22-33. The components can be Normal, Special, Accumulation, Degenerate. OPTIONS: The algorithm allows the following options as pair of arguments: \"movdim\", d  : by default movdim is 2 but it can be set to other values, and represents the number of mever variables. they should be given as the last variables of the ring. \"version\", v   :  There are two versions of the algorithm. (\"version\",1) is a full algorithm that always distinguishes correctly between 'Normal' and 'Special' components, whereas \("version\",0) can decalre a component as 'Normal' being really 'Special', but is more effective. By default (\"version\",1) is used when the number of variables is less than 4 and 0 if not. The user can force to use one or other version, but it is not recommended. \"system\", ideal F: if the initial system is passed as an argument. This is actually not used. \"comments\", c: by default it is 0, but it can be set to 1. Usually locus problems have mover coordinates, variables and tracer coordinates. The mover coordinates are to be placed as the last variables, and by default, its number is 2. If one consider locus problems in higer dimensions, the number of mover coordinates (placed as the last variables) is to be given as an option. RETURN: The locus. The output is a list of the components ( C_1,.. C_n ) where C_i =  (p_i,(p_i1,..p_is_i), type_i,l evel_i )   and type_i can be 'Normal', 'Special', Accumulation', 'Degenerate'. The 'Special' components return more information, namely the antiimage of the component, that is 0-dimensional for these kind of components. Normal components: for each point in the component, the number of solutions in the variables is finite, and the solutions depend on the point in the component. Special components:  for each point in the component, the number of solutions in the variables is finite. The antiimage of the component is 0-dimensional. Accumlation points: are 0-dimensional components whose antiimage is not zero-dimansional. Degenerate components: are components of dimension greater than 0 whose antiimage is not-zero-diemansional. The components are given in canonical P-representation. The levels of a class of locus are 1, because they represent locally closed. sets. NOTE: It can only be called after computing the grobcov of the parametrical ideal in generic representation ('ext',0), which is the default. The basering R, must be of the form Q[a_1,..,a_m][x_1,..,x_n]. KEYWORDS: geometrical locus, locus, loci. EXAMPLE: locus; shows an example" { int tes=0; int i; def R=basering; int dd=2; int comment=0; if(defined(@P)==1){tes=1; kill @P; kill @R; kill @RP;} setglobalrings(); // Options list DD=#; if (size(DD)>1){comment=1;} if(size(DD)>0){dd=DD[1];} int i; int j; int k; int moverdim=nvars(R); int version=0; int nv=nvars(R); if(nv<4){version=1;} int comment=0; ideal Fm; for(i=1;i<=(size(DD) div 2);i++) { if(DD[2*i-1]=="movdim"){moverdim=DD[2*i];} if(DD[2*i-1]=="version"){version=DD[2*i];} if(DD[2*i-1]=="system"){Fm=DD[2*i];} if(DD[2*i-1]=="comment"){comment=DD[2*i];} if(DD[2*i-1]=="family"){poly F=DD[2*i];} } int j; int k; def B0=GG[1][2]; if (equalideals(B0,ideal(1))) { return(locus0(GG)); def H0=GG[1][3][1][1]; if (equalideals(B0,ideal(1)) or (equalideals(H0,ideal(0))==0)) { return(locus0(GG,DD)); } else } N=px,py; setglobalrings(); te=indepparameters(N); if(te) if(te==0){nGP[size(nGP)+1]=GP[j];} } } } } else{"Warning! Problem with more than one mover. Not able to solve it."; list L; return(L);} } def LL=locus0(nGP,DD); kill @RP; kill @P; kill @R; return(locus0(nGP,dd)); return(LL); } example {"EXAMPLE:"; echo = 2; ring R=(0,a,b),(x4,x3,x2,x1),dp; ideal S=(x1-3)^2+(x2-1)^2-9, (4-x2)*(x3-3)+(x1-3)*(x4-1), (3-x1)*(x3-x1)+(4-x2)*(x4-x2), (4-x4)*a+(x3-3)*b+3*x4-4*x3, (a-x1)^2+(b-x2)^2-(x1-x3)^2-(x2-x4)^2; short=0; locus(grobcov(S)); " "; kill R; { "EXAMPLE:"; echo = 2; ring R=(0,a,b),(x,y),dp; short=0; "Concoid"; // Concoid ideal S96=x^2+y^2-4,(b-2)*x-a*y+2*a,(a-x)^2+(b-y)^2-1; "System="; S96; " "; // System S96= S96; locus(grobcov(S96)); kill R; ring R=(0,x,y),(x1,x2),dp; kill R; // ******************************************** ring R=(0,a,b),(x4,x3,x2,x1),dp; short=0; ideal S=(x1-3)^2+(x2-1)^2-9, (4-x2)*(x3-3)+(x1-3)*(x4-1), (3-x1)*(x3-x1)+(4-x2)*(x4-x2), (4-x4)*a+(x3-3)*b+3*x4-4*x3, (a-x1)^2+(b-x2)^2-(x1-x3)^2-(x2-x4)^2; // System S= S; locus(grobcov(S)); kill R; "********************************************"; ring R=(0,x,y),(x1,x2),dp; short=0; ideal S=-(x - 5)*(x1 - 1) - (x2 - 2)*(y - 2), (x1 - 5)^2 + (x2 - 2)^2 - 4, -2*(x - 5)*(x2 - 2) + 2*(x1 - 5)*(y - 2); "System="; S; locus(grobcov(S)); " "; ideal S1=-(x - x1)*(x1 - 4) + (x2 - y)*(x2 - 4), (x1 - 4)^2 + (x2 - 2)^2 - 4, -2*(x1 - 4)*(2*x2 - y - 4) - 2*(x - 2*x1 + 4)*(x2 - 2); "System="; S1; locus(grobcov(S1)); (x1 - 5)^2 + (x2 - 2)^2 - 4, -2*(x - 5)*(x2 - 2) + 2*(x1 - 5)*(y - 2); // System S= S; locus(grobcov(S)); } //                 mover coordinates (placed as the last variables) is to be given as an option. // //  input:      The output G of the grobcov (in generic representation, which is the default option for grobcov) //  input:      The output of locus(G); //  output: //          list, the canonical P-representation of the Normal and Non-Normal locus: //               If all levels of a class of locus are 1, then the set is locally closed. Otherwise the level //               gives the depth of the component of the constructible set. proc locusdg(list GG, list #) "USAGE:   locus(G); The input must be the grobcov  of a parametrical ideal RETURN:  The  locus. The output are the components of two constructible subsets of the locus of the parametrical system.: Normal and Non-normal. These components are given as a list of  (pi,(pi1,..pis_i),type_i,level_i) varying i, where the p's are prime ideals, the type can be: Normal, Special, Accumulation, Degenerate. NOTE:      It can only be called after computing the grobcov of the parametrical ideal in generic representation ('ext',0), which is the default. The basering R, must be of the form Q[a_1,..,a_m][x_1,..,x_n]. KEYWORDS: geometrical locus, locus, loci. EXAMPLE:  locusdg; shows an example" { def R=basering; int dd=2; int comment=0; list DD=#; if (size(DD)>1){comment=1;} if(size(DD)>0){dd=DD[1];} int i; int j; int k; def B0=GG[1][2]; if (equalideals(B0,ideal(1))) { return(locus0dg(GG)); } else { int n=nvars(R); ideal vmov=var(n-1),var(n); ideal N; intvec xw; intvec yw; for(i=1;i<=n-1;i++){xw[i]=0;} xw[n]=1; for(i=1;i<=n;i++){yw[i]=0;} yw[n-1]=1; poly px; poly py; int te=1; i=1; while( te and i<=size(B0)) { if((deg(B0[i],xw)==1) and (deg(B0[i])==1)){px=B0[i]; te=0;} i++; } i=1; te=1; while( te and i<=size(B0)) { if((deg(B0[i],yw)==1) and (deg(B0[i])==1)){py=B0[i]; te=0;} i++; } N=px,py; setglobalrings(); te=indepparameters(N); if(te) { string("locus detected that the mover must avoid point (",N,") in order to obtain the correct locus"); // eliminates segments of GG where N is contained in the basis list nGP; def GP=GG; ideal BP; for(j=1;j<=size(GP);j++) { te=1; k=1; BP=GP[j][2]; while((te==1) and (k<=size(N))) { if(pdivi(N[k],BP)[1]!=0){te=0;} k++; } if(te==0){nGP[size(nGP)+1]=GP[j];} } } } //"T_nGP="; nGP; kill @RP; kill @P; kill @R; return(locus0dg(nGP,dd)); proc locusdg(list L) "USAGE: locusdg(L)   The call must be locusdg(locus(grobcov(S))). RETURN: The output is the list of the 'Relevant' components of the locus in Dynamic Geometry:(C1,..,C:m), where C_i= ( p_i,(p_i1,..p_is_i), 'Relevant', level_i ) The 'Relevant' components are the 'Normal' and 'Accumulation' components of the locus. (See help for locus). NOTE: It can only be called after computing the locus. Calling sequence:    locusdg(locus(grobcov(S))); KEYWORDS: geometrical locus, locus, loci, dynamic geometry EXAMPLE: locusdg; shows an example" { list LL; int i; for(i=1;i<=size(L);i++) { if(typeof(L[i][3])=="string") { if((L[i][3]=="Normal") or (L[i][3]=="Accumulation")){L[i][3]="Relevant"; LL[size(LL)+1]=L[i];} } } return(LL); } example {"EXAMPLE:"; echo = 2; { "EXAMPLE:"; echo = 2; ring R=(0,a,b),(x,y),dp; short=0; // Concoid ideal S96=x^2+y^2-4,(b-2)*x-a*y+2*a,(a-x)^2+(b-y)^2-1; // System S96= S96; locus(grobcov(S96)); locusdg(locus(grobcov(S96))); kill R; "********************************************"; ring R=(0,a,b),(x4,x3,x2,x1),dp; ideal S=(x1-3)^2+(x2-1)^2-9, (4-x2)*(x3-3)+(x1-3)*(x4-1), (3-x1)*(x3-x1)+(4-x2)*(x4-x2), (4-x4)*a+(x3-3)*b+3*x4-4*x3, (a-x1)^2+(b-x2)^2-(x1-x3)^2-(x2-x4)^2; (4-x2)*(x3-3)+(x1-3)*(x4-1), (3-x1)*(x3-x1)+(4-x2)*(x4-x2), (4-x4)*a+(x3-3)*b+3*x4-4*x3, (a-x1)^2+(b-x2)^2-(x1-x3)^2-(x2-x4)^2; short=0; locus(grobcov(S)); " "; locus(grobcov(S)); locusdg(locus(grobcov(S))); kill R; ring R=(0,a,b),(x,y),dp; "********************************************"; ring R=(0,x,y),(x1,x2),dp; short=0; "Concoid"; ideal S96=x^2+y^2-4,(b-2)*x-a*y+2*a,(a-x)^2+(b-y)^2-1; "System="; S96; " "; locusdg(grobcov(S96)); kill R; ring R=(0,x,y),(x1,x2),dp; ideal S=-(x - 5)*(x1 - 1) - (x2 - 2)*(y - 2), (x1 - 5)^2 + (x2 - 2)^2 - 4, -2*(x - 5)*(x2 - 2) + 2*(x1 - 5)*(y - 2); "System="; S; locusdg(grobcov(S)); " "; ideal S1=-(x - x1)*(x1 - 4) + (x2 - y)*(x2 - 4), (x1 - 4)^2 + (x2 - 2)^2 - 4, -2*(x1 - 4)*(2*x2 - y - 4) - 2*(x - 2*x1 + 4)*(x2 - 2); "System="; S1; locusdg(grobcov(S1)); (x1 - 5)^2 + (x2 - 2)^2 - 4, -2*(x - 5)*(x2 - 2) + 2*(x1 - 5)*(y - 2); locus(grobcov(S)); locusdg(locus(grobcov(S))); } //     string s: The output of locus converted to a string readable by other programs proc locusto(list L) "USAGE:   locusto(L); The argument must be the output of locus of a parametrical ideal "USAGE: locusto(L); The argument must be the output of locus or locusdg or envelop or envelopdg. It transforms the output into a string in standard form readable in many languages (Geogebra). RETURN: The locus in string standard form NOTE:     It can only be called after computing the locus(grobcov(F)) of the parametrical ideal. The basering R, must be of the form Q[a,b,..][x,y,..]. NOTE: It can only be called after computing either - locus(grobcov(F))                -> locusto( locus(grobcov(F)) ) - locusdg(locus(grobcov(F)))  -> locusto( locusdg(locus(grobcov(F))) ) - envelop(F,C)                       -> locusto( envelop(F,C) ) -envelopdg(envelop(F,C))      -> locusto( envelopdg(envelop(F,C)) ) KEYWORDS: geometrical locus, locus, loci. EXAMPLE:  locusto; shows an example" if(size(L[i])>=3) { s[size(s)]=","; s=string(s,string(L[i][3]),"]"); s=string(s,",["); if(typeof(L[i][3])=="string") { s=string(s,string(L[i][3]),"]]"); } else { for(k=1;k<=size(L[i][3]);k++) { s=string(s,"[",L[i][3][k],"],"); } s[size(s)]="]"; s=string(s,"]"); } } if(size(L[i])>=4) } example {"EXAMPLE:"; echo = 2; ring R=(0,a,b),(x,y),dp; { "EXAMPLE:"; echo = 2; ring R=(0,x,y),(x1,y1),dp; short=0; ideal S96=x^2+y^2-4,(b-2)*x-a*y+2*a,(a-x)^2+(b-y)^2-1; "System="; S96; " "; locusto(locus(grobcov(S96))); } // Auxiliary routione ideal S=x1^2+y1^2-4,(y-2)*x1-x*y1+2*x,(x-x1)^2+(y-y1)^2-1; locusto(locus(grobcov(S))); locusto(locusdg(locus(grobcov(S)))); kill R; "********************************************"; // 1. Take a fixed line l: x1-y1=0  and consider //    the family F of a lines parallel to l passing through the mover point M // 2. Consider a circle x1^2+x2^2-25, and a mover point M(x1,x2) on it. // 3. Compute the envelop of the family of lines. ring R=(0,x,y),(x1,y1),lp; poly F=(y-y1)-(x-x1); ideal C=x1^2+y1^2-25; short=0; // Curves Family F= F; // Conditions C= C; locusto(envelop(F,C)); locusto(envelopdg(envelop(F,C))); kill R; "********************************************"; // Steiner Deltoid // 1. Consider the circle x1^2+y1^2-1=0, and a mover point M(x1,y1) on it. // 2. Consider the triangle A(0,1), B(-1,0), C(1,0). // 3. Consider lines passing through M perpendicular to two sides of ABC triangle. // 4. Obtain the envelop of the lines above. ring R=(0,x,y),(x1,y1,x2,y2),lp; ideal C=(x1)^2+(y1)^2-1, x2+y2-1, x2-y2-x1+y1; matrix M[3][3]=x,y,1,x2,y2,1,x1,0,1; poly F=det(M); short=0; // Curves Family F= F; // Conditions C= C; locusto(envelop(F,C)); locusto(envelopdg(envelop(F,C))); } // Auxiliary routine // locusdim // input: //    B:         ideal, a basis of a segment of the grobcov //    dgdim: int, the dimension of the mover (for locusdg) //    dgdim: int, the dimension of the mover (for locus) //                by default dgdim is equal to the number of variables proc locusdim(ideal B, list #) static proc locusdim(ideal B, list #) { def R=basering; } } //"T_B0="; B0; //"T_v="; v; def RR=ringlist(R); def RR0=RR; RR0[2]=v; def R0=ring(RR0); //ringlist(R0); setring(R0); def B0r=imap(R,B0); } // locusdgto: Transforms the output of locus to a string that //      can be read by different computational systems. // input: //     list L: The output of locus // output: //     string s: The output of locus converted to a string readable by other programs //                   Outputs only the relevant dynamical geometry components. //                   Without unnecessary parenthesis proc locusdgto(list LL) // // locusdgto: Transforms the output of locusdg to a string that // //      can be read by different computational systems. // // input: // //     list L: The output of locus // // output: // //     string s: The output of locus converted to a string readable by other programs // //                   Outputs only the relevant dynamical geometry components. // //                   Without unnecessary parenthesis // proc locusdgto(list LL) // "USAGE: locusdgto(L); //           The argument must be the output of locusdg of a parametrical ideal //           It transforms the output into a string in standard form //           readable in many languages (Geogebra). // RETURN: The locusdg in string standard form // NOTE: It can only be called after computing the locusdg(grobcov(F)) of the //           parametrical ideal. //           The basering R, must be of the form Q[a,b,..][x,y,..]. // KEYWORDS: geometrical locus, locus, loci. // EXAMPLE:  locusdgto; shows an example" // { //   int i; int j; int k; int short0=short; int ojo; //   int te=0; //   short=0; //   if(size(LL)==0){ojo=1; list L;} //   else //   { //     def L=LL; //   } //   string s="["; string sf="]"; string st=s+sf; //   if(size(L)==0){return(st);} //   ideal p; //   ideal q; //   for(i=1;i<=size(L);i++) //   { //     if(L[i][3]=="Relevant") //     { //       s=string(s,"[["); //       for (j=1;j<=size(L[i][1]);j++) //       { //         s=string(s,L[i][1][j],","); //       } //       s[size(s)]="]"; //       s=string(s,",["); //       for(j=1;j<=size(L[i][2]);j++) //       { //         s=string(s,"["); //         for(k=1;k<=size(L[i][2][j]);k++) //         { //           s=string(s,L[i][2][j][k],","); //         } //         s[size(s)]="]"; //         s=string(s,","); //       } //       s[size(s)]="]"; //       s=string(s,"]"); //       s[size(s)]="]"; //       s=string(s,","); //     } //   } //   if(s=="["){s="[]";} //   else{s[size(s)]="]";} //   short=short0; //   return(s); // } // example // {"EXAMPLE:"; echo = 2; //   ring R=(0,a,b),(x,y),dp; //   short=0; //   ideal S96=x^2+y^2-4,(b-2)*x-a*y+2*a,(a-x)^2+(b-y)^2-1; //   "System="; S96; " "; //   "locusdgto(locusdg(grobcov(S96)))="; //   locusdgto(locusdg(grobcov(S96))); // } static proc norspec(ideal F) { def RR=basering; int i; int j; int k; int short0=short; int ojo; int te=0; if(defined(@P)){te=1;} else{setglobalrings();} setring @P; def Rx=ringlist(RR); def Rx=ringlist(RR); def @P=ring(Rx[1]); list Lx; Lx[1]=0; Lx[2]=Rx[2]+Rx[1][2]; Lx[3]=Rx[1][3]; Lx[4]=Rx[1][4]; Rx[1]=0; def D=ring(Rx); def @RP=D+@P; exportto(Top,@R);      // global ring Q[a][x] exportto(Top,@P);      // global ring Q[a] exportto(Top,@RP);     // global ring K[x,a] with product order setring(RR); } // envelop // Input: n arguments //   poly F: the polynomial defining the family of curves in ring R=0,(x,y),(x1,..,xn),lp; //   ideal C=g1,..,g_{n-1}:  the set of constraints // Output: the components of the envolvent; proc envelop(poly F, ideal C, list #) "USAGE: envelop(F,C); The first argument F must be the family of curves of which on want to compute the envelop. The second argument C must be the ideal of conditions over the variables, and should contain as polynomials as the number of variables -1. RETURN: The components of the envelop with its taxonomy: The taxonomy distinguishes 'Normal', 'Special', 'Accumulation', 'Degenerate' components. In the case of 'Special' components, it also outputs the antiimage of the component and an integer (0-1). If the integer is 0 the component is not a curve of the family and is not considered as 'Relevant' by the envelopdg routine applied to it, but is considered as 'Relevant' if the integer is 1. NOTE: grobcov is called internally. The basering R, must be of the form Q[a][x] (a=parameters, x=variables). KEYWORDS: geometrical locus, locus, loci, envelop EXAMPLE:  envelop; shows an example" { int tes=0; int i;   int j; def R=basering; if(defined(@P)==1){tes=1; kill @P; kill @R; kill @RP;} setglobalrings(); // Options list DD=#; int moverdim=nvars(R); int version=0; int nv=nvars(R); if(nv<4){version=1;} int comment=0; ideal Fm; for(i=1;i<=(size(DD) div 2);i++) { if(DD[2*i-1]=="movdim"){moverdim=DD[2*i];} if(DD[2*i-1]=="version"){version=DD[2*i];} if(DD[2*i-1]=="system"){Fm=DD[2*i];} if(DD[2*i-1]=="comment"){comment=DD[2*i];} } int n=nvars(R); list v; for(i=1;i<=n;i++){v[size(v)+1]=var(i);} def MF=jacob(F); def TMF=transpose(MF); def Mg=MF; def TMg=TMF; for(i=1;i<=n-1;i++) { Mg=jacob(C[i]); TMg=transpose(Mg); TMF=concat(TMF,TMg); } poly J=det(TMF); ideal S=ideal(F)+C+ideal(J); DD[size(DD)+1]="family"; DD[size(DD)+1]=F; def G=grobcov(S,DD); def L=locus(G, DD); return(L); } example { "EXAMPLE:"; echo = 2; // Steiner Deltoid // 1. Consider the circle x1^2+y1^2-1=0, and a mover point M(x1,y1) on it. // 2. Consider the triangle A(0,1), B(-1,0), C(1,0). // 3. Consider lines passing through M perpendicular to two sides of ABC triangle. // 4. Obtain the envelop of the lines above. ring R=(0,x,y),(x1,y1,x2,y2),lp; ideal C=(x1)^2+(y1)^2-1, x2+y2-1, x2-y2-x1+y1; matrix M[3][3]=x,y,1,x2,y2,1,x1,0,1; poly F=det(M); short=0; if(size(LL)==0){ojo=1; list L;} else { def L=imap(RR,LL); } string s="["; string sf="]"; string st=s+sf; if(size(L)==0){return(st);} ideal p; ideal q; // Curves Family F= F; // Conditions C= C; def Env=envelop(F,C); Env; } // envelopdg // Input: list L: the output of envelop(poly F, ideal C, list #) // Output: the relevant components of the envolvent in dynamic geometry; proc envelopdg(list L) "USAGE: envelopdg(L); The input list L must be the output of the call to the routine 'envolop' of the family of curves RETURN: The relevant components of the envelop in Dynamic Geometry. 'Normal' and 'Accumulation' components are always considered 'Relevant'. 'Special' components of the envelop outputs three objects in its characterization: 'Special', the antiimage ideal, and the integer 0 or 1, that indicates that the given component is formed (1) or is not formed (0) by curves of the family. Only if yes, 'envelopdg' considers the component as 'Relevant' . NOTE: It must be called to the output of the 'envelop' routine. The basering R, must be of the form Q[a,b,..][x,y,..]. KEYWORDS: geometrical locus, locus, loci, envelop. EXAMPLE:  envelop; shows an example" { list LL; list Li; int i; for(i=1;i<=size(L);i++) { if(L[i][3]=="Relevant") { s=string(s,"[["); for (j=1;j<=size(L[i][1]);j++) { s=string(s,L[i][1][j],","); } s[size(s)]="]"; s=string(s,",["); for(j=1;j<=size(L[i][2]);j++) { s=string(s,"["); for(k=1;k<=size(L[i][2][j]);k++) if(typeof(L[i][3])=="string") { if((L[i][3]=="Normal") or (L[i][3]=="Accumulation")){Li=L[i]; Li[3]="Relevant"; LL[size(LL)+1]=Li;} } else { if(typeof(L[i][3])=="list") { if(L[i][3][3]==1) { s=string(s,L[i][2][j][k],","); Li=L[i]; Li[3]="Relevant"; LL[size(LL)+1]=Li; } s[size(s)]="]"; s=string(s,","); } s[size(s)]="]"; s=string(s,"]"); s[size(s)]="]"; s=string(s,","); } } if(s=="["){s="[]";} else{s[size(s)]="]";} setring(RR); short=short0; if(te==0){kill @P; kill @R; kill @RP;} return(s); } } } return(LL); } example { "EXAMPLE:"; echo = 2; // 1. Take a fixed line l: x1-y1=0  and consider //    the family F of a lines parallel to l passing through the mover point M // 2. Consider a circle x1^2+x2^2-25, and a mover point M(x1,x2) on it. // 3. Compute the envelop of the family of lines. ring R=(0,x,y),(x1,y1),lp; short=0; poly F=(y-y1)-(x-x1); ideal C=x1^2+y1^2-25; short=0; // Curves Family F= F; // Conditions C= C; envelop(F,C); envelopdg(envelop(F,C)); } //********************* End locus **************************** //********************* Begin AddCons ********************** // Input: L1,L2: lists of components with common top // Output L: list of the union of L1 and L2. static proc Add2ComWithCommonTop(list L1, list L2) { int i; int j; ideal pij; list L; list Lp; list PR; int k; for(i=1;i<=size(L1[2]);i++) { for(j=1;j<=size(L2[2]);j++) { pij=std(L1[2][i]+L2[2][j]); PR=minGTZ(pij); for(k=1;k<=size(PR);k++) { Lp[size(Lp)+1]=PR[k]; } } } for(i=1; i<=size(Lp);i++) { for(j=i+1;j<=size(Lp);j++) { if(idcontains(Lp[i],Lp[j])) {Lp=delete(Lp,j);} } for(j=1;j... //"T_te="; te; pause(); if(te==0) { modif=1; kill Pij; ideal Pij=1; for(l=1; l<=size(L1[k][2]);l++) { pkl=L1[k][2][l]; pijkl=std(pij+pkl); Pij=intersect(Pij,pijkl); } kill Lp; list Lp; Lp=minGTZ(Pij); for(m=1;m<=size(Lp);m++) { nl[size(nl)+1]=Lp[m]; } A[i][2]=delete(A[i][2], j); for(n=1;n<=size(nl);n++) { A[i][2][size(A[i][2])+1]=nl[n]; } } // end if(te==0) } // end if(idcontains(pij,pk)==1) }  // end if (k!=i) } // end for k }  // end for j kill vvv; list vvv; if(modif==1) // Select the maximal ideals of the set A[I][2][j] { kill nl; list nl; nl=selectminideals(A[i][2]); kill AAA; list AAA; for(j=1;j<=size(nl);j++) { AAA[size(AAA)+1]=A[i][2][nl[j]]; } A[i][2]=AAA; } // end if(modif=1) }  // end for i modif=1-equallistsAall(A,A0); } // end while(modif==1) return(A); } // Input: //   A: list of the top P-reps of one level //   B: list of remaining lower P-reps that have not yeen be possible to add // Output: //   Bn: list B where the elements that are completely included in A are removed and the parts that are //         included in A also. static proc ReduceB(list A,list B) { int i; int j; list nl; list Bn; int te; int k; int elim; ideal pC; ideal pD; ideal pCD; ideal pBC; list nB; int j0; for(i=1;i<=size(B);i++) { j=1; te=0; while((te==0) and (j<=size(A))) { if(idcontains(B[i][1],A[j][1])){te=1; j0=j;} else{j++;} } pD=B[i][2][1]; for(k=2;k<=size(B[i][2]);k++){pD=intersect(pD,B[i][2][k]);} pC=A[j0][2][1]; for(k=2;k<=size(A[j0][2]);k++) {pC=intersect(pC,A[j0][2][k]);} pCD=std(pD+pC); pBC=std(B[i][1]+pC); elim=0; if(idcontains(pBC,pCD)){elim=1;}   // B=delfromlist(B,i);i=i-1; else { nB=Prep(pBC,pCD); if(equalideals(nB[1][1],ideal(1))==0) { for(k=1;k<=size(nB);k++){Bn[size(Bn)+1]=nB[k];} } } }   // end for i return(Bn); } // AddConsP: given a set of components of locally closed sets in P-representation, it builds the //       canonical P-representation of the corresponding constructible set, of its union, //       including levels it they are. proc AddConsP(list L) "USAGE:   AddConsP(L) Input L: list of components in P-rep to be added [  [[p_1,[p_11,..,p_1,r1]],..[p_k,[p_k1,..,p_kr_k]]  ] RETURN: list of lists of levels of the different locally closed sets of the canonical P-rep of the constructible. [  [level_1,[ [Comp_11,..Comp_1r_1] ] ], .. , [level_s,[ [Comp_s1,..Comp_sr_1] ] ] where level_i=i,   Comp_ij=[ p_i,[p_i1,..,p_it_i] ] is a prime component. NOTE:     Operates in a ring R=Q[u_1,..,u_m] KEYWORDS: Constructible set, Canoncial form EXAMPLE:  AddConsP; shows an example" { list LL; list A; list B; list L1; list T; int level=0; list h; LL=L; int i; while(size(LL)!=0) { level++; L1=SepareAB(LL); A=L1[1]; B=L1[2]; LL=L1[3]; A=CompleteA(A,LL); for(i=1;i<=size(A);i++) { LL[i]=A[i]; } h[1]=level; h[2]=A; T[size(T)+1]=h; LL=ReduceB(A,B); } return(T); } example { "EXAMPLE:"; echo = 2; if (defined(Grobcov::@P)){kill Grobcov::@P; kill Grobcov::@R; kill Grobcov::@RP;} ring R=0,(x,y,z),lp; short=0; ideal P1=x; ideal Q1=x,y; ideal P2=y; ideal Q2=y,z; list L=list(Prep(P1,Q1)[1],Prep(P2,Q2)[1]); L; AddConsP(L); } // AddCons:  given a set of  locally closed sets by pairs of ideal, it builds the //       canonical P-representation of the corresponding constructible set, of its union, //       including levels it they are. //       It makes a call to AddConsP after transforming the input. // Input list of lists of pairs of ideals representing locally colsed sets: //     L=  [ [E1,N1], .. , [E_s,N_s] ] // Output: The canonical frepresentation of the constructible set union of the V(E_i) \ V(N_i) //     T=[  [level_1,[ [p_1,[p_11,..,p_1s_1]],.., [p_k,[p_k1,..,p_ks_k]] ]],, .. , [level_r,[..       ]]  ] proc AddCons(list L) "USAGE:   AddCons(L) Calls internally AddConsP and allows a different input. Input L: list of pairs of of ideals [ [P_1,Q_1], .., [Pr,Qr] ] representing a set of locally closed setsV(P_i) \ V(Q_i) to be added. RETURN: list of lists of levels of the different locally closed sets of the canonical P-rep of the constructible. [  [level_1,[ [Comp_11,..Comp_1r_1] ] ], .. , [level_s,[ [Comp_s1,..Comp_sr_1] ] ] where level_i=i,   Comp_ij=[ p_i,[p_i1,..,p_it_i] ] is a prime component. NOTE: Operates in a ring R=Q[u_1,..,u_m] KEYWORDS: Constructible set, Canoncial form EXAMPLE: AddCons; shows an example" { int i; list H; list P; int j; for(i=1;i<=size(L);i++) { P=Prep(L[i][1],L[i][2]); for(j=1;j<=size(P);j++) { H[size(H)+1]=P[j]; } } return(AddConsP(H)); } example { "EXAMPLE:"; echo = 2; if (defined(Grobcov::@P)){kill Grobcov::@P; kill Grobcov::@R; kill Grobcov::@RP;} ring R=0,(x,y,z),lp; short=0; ideal P1=x; ideal Q1=x,y; ideal P2=y; ideal Q2=y,z; list L=list(list(P1,Q1), list(P2,Q2) ); L; AddCons(L); } // AddLocus: auxilliary routine for locus0 that computes the components of the constructible: // Input:  the list of locally closed sets to be added, each with its type as third argument //     L=[ [LC[11],,,LC[1k_1],  .., [LC[r1],,,LC[rk_r] ] where //            LC[1]=[p1,[p11,..,p1k],typ] // Output:  the list of components of the constructible union of the L, with the type of the corresponding top //               and the level of the constructible //     L4= [[v1,p1,[p11,..,p1l],typ_1,level]_1 ,.. [vs,ps,[ps1,..,psl],typ_s,level_s] static proc AddLocus(list L) { //  int te0=0; //  def RR=basering; //  if(defined(@P)){te0=1;  def Rx=@R;  kill @P; setring RR;} list L1; int i; int j;  list L2; list L3; list l1; list l2; intvec v; for(i=1; i<=size(L); i++) { for(j=1;j<=size(L[i]);j++) { l1[1]=L[i][j][1]; l1[2]=L[i][j][2]; l2[1]=l1[1]; if(size(L[i][j])>2){l2[3]=L[i][j][3];} v[1]=i; v[2]=j; l2[2]=v; L1[size(L1)+1]=l1; L2[size(L2)+1]=l2; } } L3=AddConsP(L1); list L4; int level; ideal p1; ideal pp1; int t; int k; int k0; string typ; list l4; for(i=1;i<=size(L3);i++) { level=L3[i][1]; for(j=1;j<=size(L3[i][2]);j++) { p1=L3[i][2][j][1]; t=1; k=1; while((t==1) and (k<=size(L2))) { pp1=L2[k][1]; if(equalideals(p1,pp1)){t=0; k0=k;} k++; } if(t==0) { v=L2[k0][2]; } else{"ERROR p1 NOT FOUND";} l4[1]=v; l4[2]=p1; l4[3]=L3[i][2][j][2];  l4[5]=level; if(size(L2[k0])>2){l4[4]=L2[k0][3];} L4[size(L4)+1]=l4; } } return(L4); } //********************* End AddCons ********************** ;