Changeset 2203cf in git
 Timestamp:
 Oct 6, 2016, 3:17:04 PM (8 years ago)
 Branches:
 (u'fiekerDuVal', '117eb8c30fc9e991c4decca4832b1d19036c4c65')(u'spielwiese', '4fa496c5b814dbde0a905c54f0796301d03f6dc9')
 Children:
 b61d2e9bb22594df3a38d0289c0394fb3c068301
 Parents:
 17510a828cdc33aa541d40a9bbfb825b72812521
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

Singular/LIB/grobcov.lib
r17510a r2203cf 1 1 // 2 version="version grobcov.lib 4.0. 2.0 Jul_2015"; // $Id$3 // version L; July_2015;2 version="version grobcov.lib 4.0.3.3 Oct_2016 "; // $Id$ 3 // version M; October_2016; 4 4 category="General purpose"; 5 5 info=" 6 LIBRARY: grobcov .libGroebner Cover for parametric ideals.6 LIBRARY: grobcovM.lib Groebner Cover for parametric ideals. 7 7 8 8 Comprehensive Groebner Systems, Groebner Cover, Canonical Forms, … … 16 16 (http://wwwma2.upc.edu/~montes/). 17 17 18 AUTHORS: Antonio Montes , Hans Schoenemann. 18 IMPORTANT: The book, not yet published: 19 20 A. Montes. 21 \"Discussing Parametric Polynomial Systems: The Groebner Cover\" 22 23 can be used as a user manual of all the routines included in this library. 24 It defines and proves all the theoretic results used here, and 25 shows examples of all the routines. 26 27 I expect that it will be published soon. 28 29 AUTHORS: Antonio Montes (Universitat Politecnica de Catalunya), 30 Hans Schoenemann (Technische Universitaet Kaiserslautern). 19 31 20 32 OVERVIEW: 21 33 In 2010, the library was designed to contain MontesWibmer's 22 algorithms for comput ethe canonical Groebner Cover of a34 algorithms for computing the canonical Groebner Cover of a 23 35 parametric ideal as described in the paper: 24 36 … … 27 39 Journal of Symbolic Computation 45 (2010) 13911425. 28 40 29 The central routine is grobcov. Given a parametri c30 ideal,grobcov outputs its Canonical Groebner Cover, consisting41 The central routine is grobcov. Given a parametri ideal, 42 grobcov outputs its Canonical Groebner Cover, consisting 31 43 of a set of pairs of (basis, segment). The basis (after 32 44 normalization) is the reduced Groebner basis for each point … … 45 57 reduced Comprehensive Groebner System with constant lpp. 46 58 For this purpose, in this library, the implemented algorithm is 47 KapurSunWang algorithm, because it is the most efficient48 algorithm known for this purpose.59 KapurSunWang algorithm, because it is actually the most 60 efficient algorithm known for this purpose. 49 61 50 62 D. Kapur, Y. Sun, and D.K. Wang. … … 59 71 union of locally closed sets. It is used in the new version for the 60 72 computation of loci and envelops. It determines the canonical locally closed 61 level sets of a constructuble. They will be described in a forthcoming paper:73 level sets of a constructuble. It is described in: 62 74 63 75 J.M. Brunat, A. Montes, 64 76 \"Computing the canonical representation of constructible sets\". 65 Submited to Mathematics in Computer Science. July 2015.77 Math. Comput. Sci. (2016) 19: 165178. 66 78 67 79 A new set of routines (locus, locusdg, locusto) has been included to … … 73 85 ComputerAided Design 56 (2014) 2233. 74 86 75 Recently also routines for computing the envelop of a family76 of curves (enverlop, envelopdg), to be used in Dynamic Geometry,77 has been included and will be described in a forthcoming paper:78 79 A banades, Botana, Montes, Recio:80 \ ''Envelops in Dynamic Geometry using the Groebner cover\''.81 82 83 This version was finished on 31/07/2015 87 Recently also routines for computing the generalized envelop of a family 88 of hypersurfaces (envelop), to be used in Dynamic Geometry, 89 has been included and is described in the book (not yet published 90 91 A. Montes. 92 \"Discussing Parametric Polynomial Systems: The Groebner Cover\" 93 94 This version was finished on 31/010/2016 95 84 96 85 97 NOTATIONS: All given and determined polynomials and ideals are in the … … 94 106 95 107 PROCEDURES: 108 setglobalrings(); Generates the global rings @R, @P and @RP that are 109 respectively the rings Q[a][x], Q[a], Q[x,a]. 110 It is called inside each of the fundamental routines of the 111 library: grobcov, cgsdr, locus, locusdg and killed before 112 the output. 113 114 If the user want to call the routines Prep or Crep on a 115 parametric ideal Q[a][x] for locally closed sets on Q[a], 116 setglobalrings must be called previously. 117 In the actual version, after the call of setglobalrings on 118 Q[a][x], the public names of the defined ideals generated by 119 setglobalrings are Grobcov::@R, Grobcov::@P, Grobcov::@RP. 120 96 121 grobcov(F); Is the basic routine giving the canonical 97 122 Groebner Cover of the parametric ideal F. … … 103 128 reduced Comprehensive Groebner System that 104 129 is used in grobcov, that can also be used 105 independently if only theCGS is required.130 independently if only a CGS is required. 106 131 It is a more efficient routine than buildtree 107 132 (the own routine of 2010 that is no more used). 108 Now, K SWalgorithm is used.133 Now, KapurSunWang (KSW) algorithm is used. 109 134 110 135 pdivi(f,F); Performs a pseudodivision of a parametric polynomial 111 136 by a parametric ideal. 112 137 138 113 139 pnormalf(f,E,N); Reduces a parametric polynomial f over V(E) \ V(N) 114 ( 140 (E is the null ideal and N the nonnull ideal ) 115 141 over the parameters. 116 142 117 143 Crep(N,M); Computes the canonical Crepresentation of V(N) \ V(M). 144 If the user want to call the routines Crep on a parametric 145 ideal Q[a][x] for locally closed sets on Q[a], setglobalrings() 146 must be called previously. 147 In the actual version, after the call of setglobalrings on 148 Q[a][x], the public names of the defined ideals generated by 149 setglobalrings are Grobcov::@R, Grobcov::@P, Grobcov::@RP. 118 150 119 151 Prep(N,M); Computes the canonical Prepresentation of V(N) \ V(M). 152 If the user want to call the routines Prep on a parametric 153 ideal Q[a][x] for locally closed sets on Q[a], setglobalrings() 154 must be called previously. 155 In the actual version, after the call of setglobalrings on 156 Q[a][x], the public names of the defined ideals generated by 157 setglobalrings are Grobcov::@R, Grobcov::@P, Grobcov::@RP. 120 158 121 159 PtoCrep(L) Starting from the canonical Prep of a locally closed set … … 128 166 full representation of the bases. With the default 129 167 option ('ext',0) only the generic representation of 130 the bases arecomputed, and one can obtain the full168 the bases is computed, and one can obtain the full 131 169 representation using extend. 132 170 133 ConsLevels(L); Given a list L of locally closed sets, it returns the c anonical levels134 of the constructible set of the union of them, as well as the levels135 of the complement.It is described in the paper171 ConsLevels(L); Given a list L of locally closed sets, it returns the closures of the canonical 172 levels of the constructible set and its complements of the union of them. 173 It is described in the paper 136 174 137 175 J.M. Brunat, A. Montes, 138 176 \"Computing the canonical representation of constructible sets\". 139 Submited to Mathematics in Computer Science. July 2015. 177 Math. Comput. Sci. (2016) 19: 165178. 178 179 ConsLevelsToLevels(L): Transforms the output of ConsLevels into the proper 180 Levels of the constructible set. 140 181 141 182 locus(G); Special routine for determining the geometrical locus of points … … 151 192 ComputerAided Design 56 (2014) 2233. 152 193 153 The components can be 'Normal', 'Special', 'Accumulation', 'Degenerate'.154 The output are the components isgiven in Pcanonical form194 The components can be \"Normal\", \"Special\", \"Accumulation\", \"Degenerate\". 195 The output are the components given in Pcanonical form 155 196 It also detects automatically a possible point that is to be 156 avoided by the mover, whose coordinates must be the last two197 avoided by the mover, whose coordinates must be the last 157 198 coordinates in the definition of the ring. If such a point is detected, 158 199 then it eliminates the segments of the grobcov depending on the 159 200 point that is to be avoided. 160 201 161 locusdg(G); Is a special routine that determines the 'Relevant'components202 locusdg(G); Is a special routine that determines the \"Relevant\" components 162 203 of the locus in dynamic geometry. It is to be called to the output of locus 163 204 and selects from it the useful components. 164 205 165 envelop(F,C); Special routine for determining the envelop of a family of curves166 F in Q[x ,y][x_1,..xn] depending on a ideal of constraints C in Q[x_1,..,xn].206 envelop(F,C); Special routine for determining the envelop of a family of hypersurfaces 207 F in Q[x1,..,xn][t1,..,tm] depending on a ideal of constraints C in Q[t1,..,tm]. 167 208 It detemines the different components as well as its type: 168 'Normal', 'Special', 'Accumulation', 'Degenerate'. And169 it also classifies the Specialcomponents, determining the209 \"Normal\", \"Special\", \"Accumulation\", \"Degenerate\". And 210 it also classifies the \"Special\" components, determining the 170 211 zero dimensional antiimage of the component and verifying if 171 212 the component is a special curve of the family or not. … … 173 214 to obtain the complete result. 174 215 The taxonomy that it provides, as well as the algorithms involved 175 will be described in a forthcoming paper: 176 177 Abanades, Botana, Montes, Recio: 178 \''Envelops in Dynamic Geometry using the Gr\"obner cover\''. 179 180 envelopdg(ev); Is a special routine to determine the 'Relevant' components 181 of the envelop of a family of curves to be used in Dynamic Geometry. 182 It must be called to the output of envelop(F,C). 183 184 locusto(L); Transforms the output of locus, locusdg, envelop and envelopdg 216 are described in the book: (not yet published) 217 218 A. Montes. 219 \"Discussing Parametric Polynomial Systems: The Groebner Cover\" 220 221 locusto(L); Transforms the output of locus, locusdg, envelop 185 222 into a string that can be reed from different computational systems. 186 223 224 AssocTanToEnv(F,C,E); Having computed an envelop component E 225 of a family of hypersurfaces F, with constraints C, 226 it returns the parameter values of the unique 227 associated tangent hypersurface of the family 228 passing at one point of the envelop component E. 229 230 FamElemsAtEnvCompPoints(F,C,E) Having computed the an envelop 231 component E of a family of hypersurfaces F, with constraints C, 232 it returns the parameter values of all the hypersurfaes of the family 233 passing at one point of the envelop component E. 234 235 discrim(f,x); Determines the factorized discriminant of a degree 2 polynomial 236 in the variable x. The polynomial can be defined on any ring 237 where x is a variable. The polynomial f can depend on parameters and variables 238 239 WLemma(F,A); Given an ideal F in K[a][x] and a priime ideal A in K[a], 240 it returns the list (lpp,B,S) were B is the reduced Groebner basis 241 of the specialized F over the segment computed in Prepresentation 242 (or optionally in Crepresentation). The basis is given by Iregular 243 functions. 187 244 188 245 SEE ALSO: compregb_lib … … 194 251 // ************ Begin of the grobcov library ********************* 195 252 253 // Development of the library: 196 254 // Library grobcov.lib 197 255 // (Groebner Cover): … … 200 258 // Uses buildtree for cgsdr 201 259 // Final data: 372008 202 // Release 2: (priv ate)260 // Release 2: (prived) 203 261 // Initial data: 692009 204 262 // Last version using buildtree for cgsdr 205 263 // Final data: 25102011 206 // Release B: (priv ate)264 // Release B: (prived) 207 265 // Initial data: 172012 208 266 // Uses both buildtree and KSW for cgsdr … … 213 271 // Final data: 21112013 214 272 // Release L: (public) 215 // New routine ConsLevels: 1072015 216 // Updated locus: 1072015 (uses Conslevels) 217 // New routines for computing the envelop of a family of curves: 2212015 218 // Final data: 2272015 273 // New routine ConsLevels: 2512016 274 // Updated locus: 1072015 (uses ConsLevels) 275 // Release M: (public) 276 // New routines for computing the envelop of a family of 277 // hypersurfaces and associated questions: 2242016: 2092016 278 // New routine WLemma for computing the result of 279 // Wibmer's Lemma: 1992016 280 // Final version October 2016 219 281 220 282 //*************Auxiliary routines************** … … 314 376 } 315 377 316 317 378 // equalideals 318 379 // Input: ideals F and G; … … 325 386 while ((i<=size(F)) and (t)) 326 387 { 327 if (F[i]!=G[i]){t=0;}388 if (F[i]!=G[i]){t=0;} 328 389 i++; 329 390 } … … 372 433 return(t); 373 434 } 374 375 435 376 436 // selectminideals … … 404 464 return(P); 405 465 } 406 407 466 408 467 // Auxiliary routine … … 749 808 //} 750 809 751 static proc setglobalrings() 752 // "USAGE: setglobalrings(); 753 // No arguments 754 // RETURN: After its call the rings Grobcov::@R=Q[a][x], Grobcov::@P=Q[a], 755 // Grobcov::@RP=Q[x,a] are defined as global variables. 756 // (a=parameters, x=variables) 757 // NOTE: It is called internally by many basic routines of the 758 // library grobcov, cgsdr, extend, pdivi, pnormalf, locus, locusdg, 759 // envelop, envelopdg, and killed before the output. 760 // The user does not need to call it, except when it is interested 761 // in using some internal routine of the library that 762 // uses these rings. 763 // The basering R, must be of the form Q[a][x], (a=parameters, 764 // x=variables), and should be defined previously. 765 // KEYWORDS: ring, rings 766 // EXAMPLE: setglobalrings; shows an example" 810 proc setglobalrings() 811 "USAGE: setglobalrings(); 812 No arguments 813 RETURN: After its call the rings Grobcov::@R=Q[a][x], Grobcov::@P=Q[a], 814 Grobcov::@RP=Q[x,a] are defined as global variables. 815 (a=parameters, x=variables) 816 NOTE: It is called internally by many basic routines of the 817 library grobcov, cgsdr, extend, pdivi, pnormalf, locus, locusdg, 818 envelop, WLemma, and killed before the output. 819 The user does not need to call it, except when it is interested 820 in using some internal routine of the library that 821 uses these rings. 822 The basering R, must be of the form Q[a][x], (a=parameters, 823 x=variables), and should be defined previously. 824 When calling Prep or Crep of a locally closed set on a ring 825 R[a][x], if setglobalrings() is called previously, then it is 826 assumed that computing the P or Crepresentation 827 is considered for locally closed sets on Q[a]. 828 KEYWORDS: ring; rings 829 EXAMPLE: setglobalrings; shows an example" 767 830 { 768 831 if (defined(@P)) … … 787 850 setring(RR); 788 851 }; 789 // example 790 // { 791 // "EXAMPLE:"; echo = 2; 792 // ring R=(0,a,b),(x,y,z),dp; 793 // setglobalrings(); 794 // " ";"R=";R; 795 // " ";"Grobcov::@R=";Grobcov::@R; 796 // " ";"Grobcov::@P=";Grobcov::@P; 797 // " ";"Grobcov::@RP=";Grobcov::@RP; 798 // " "; "ringlist(Grobcov::@R)="; ringlist(Grobcov::@R); 799 // " "; "ringlist(Grobcov::@P)="; ringlist(Grobcov::@P); 800 // " "; "ringlist(Grobcov::@RP)="; ringlist(Grobcov::@RP); 801 // } 852 example 853 { 854 "EXAMPLE:"; echo = 2; 855 ring R=(0,a,b),(x,y,z),dp; 856 setglobalrings(); 857 R; 858 859 Grobcov::@R; 860 861 Grobcov::@P; 862 863 Grobcov::@RP; 864 865 ringlist(Grobcov::@R); 866 867 ringlist(Grobcov::@P); 868 869 ringlist(Grobcov::@RP); 870 } 802 871 803 872 // cld : clears denominators of an ideal and normalizes to content 1 … … 859 928 list (r,q,m) such that m*f=r+sum(q.G), and no lpp of a divisor 860 929 divides a pp of r. 861 KEYWORDS: division ,reduce930 KEYWORDS: division; reduce 862 931 EXAMPLE: pdivi; shows an example" 863 932 { … … 919 988 "EXAMPLE:"; echo = 2; 920 989 ring R=(0,a,b,c),(x,y),dp; 990 // Divisor="; 921 991 poly f=(abac)*xy+(ab)*x+(5c); 922 // Divisor="; 923 f; 992 // Dividends="; 924 993 ideal F=ax+b,cy+a; 925 // Dividends=";926 F;927 def r=pdivi(f,F);928 994 // (Remainder, quotients, factor)="; 929 r;995 def r=pdivi(f,F); r; 930 996 // Verifying the division: 931 997 r[3]*f(r[2][1]*F[1]+r[2][2]*F[2]+r[1]); … … 1091 1157 NOTE: Should be called from ring Q[a][x]. 1092 1158 Ideals E and N must be given by polynomials in Q[a]. 1093 KEYWORDS: division , pdivi,reduce1159 KEYWORDS: division; pdivi; reduce 1094 1160 EXAMPLE: pnormalf; shows an example" 1095 1161 { … … 1114 1180 example 1115 1181 { 1116 1182 "EXAMPLE:"; echo = 2; 1117 1183 ring R=(0,a,b,c),(x,y),dp; 1118 1184 short=0; … … 1120 1186 ideal E=(c1); 1121 1187 ideal N=ab; 1122 1123 1188 pnormalf(f,E,N); 1124 1189 } … … 1390 1455 // the list (a,b) of the canonical ideals 1391 1456 // the Crep of V(N) \ V(M) 1392 // Assumed to be called in the ring @R 1393 // Works on the ring @P 1457 // Assumed to be called in the ring @R or the ring @P or a ring ring Q[a] 1394 1458 proc Crep(ideal N, ideal M) 1395 1459 "USAGE: Crep(N,M); 1396 Input:ideal N (null ideal) (not necessarily radical nor maximal)1397 1460 ideal N (null ideal) (not necessarily radical nor maximal) 1461 ideal M (hole ideal) (not necessarily containing N) 1398 1462 RETURN: The canonical Crepresentation of the locally closed set. 1399 1463 [ P,Q ], a pair of radical ideals with P included in Q, 1400 1464 representing the set V(P) \ V(Q) = V(N) \ V(M) 1401 1465 NOTE: Operates in a ring R=Q[a] (a=parameters) 1402 KEYWORDS: locally closed set, canoncial form 1466 If the user want to call the routine Crep on a parametric 1467 ideal Q[a][x] for locally closed sets on Q[a], setglobalrings() 1468 must be called previously. 1469 1470 KEYWORDS: locally closed set; canoncial form 1403 1471 EXAMPLE: Crep; shows an example" 1472 { 1473 int te; 1474 def RR=basering; 1475 if(defined(@P)){te=1; setring(@P); ideal Np=imap(RR,N); ideal Mp=imap(RR,M);} 1476 else {te=0; def Np=N; def Mp=M;} 1477 def La=Crep0(Np,Mp); 1478 if(size(La)==0) 1479 { 1480 if(te==1) {setring(RR); list LL;} 1481 if(te==0){list LL;} 1482 return(LL); 1483 } 1484 else 1485 { 1486 if(te==1) {setring(RR); def LL=imap(@P,La);} 1487 if(te==0){def LL=La;} 1488 return(LL); 1489 } 1490 } 1491 example 1492 { 1493 "EXAMPLE:"; echo = 2; 1494 short=0; 1495 ring R=0,(x,y,z),lp; 1496 ideal E=x*(x^2+y^2+z^225)*(x+y); 1497 ideal N=x*(x3)*(x+y),(y4)*(x+y); 1498 Crep(E,N); 1499 } 1500 1501 // Crep0 1502 // Computes the Crepresentation of V(N) \ V(M). 1503 // input: 1504 // ideal N (null ideal) (not necessarily radical nor maximal) 1505 // ideal M (hole ideal) (not necessarily containing N) 1506 // output: 1507 // the list (a,b) of the canonical ideals 1508 // the Crep0 of V(N) \ V(M) 1509 // Assumed to be called in a ring Q[x] (example @P) 1510 static proc Crep0(ideal N, ideal M) 1404 1511 { 1405 1512 list l; … … 1415 1522 ideal P=ideal(1); 1416 1523 L=minGTZ(Np); 1524 //"T_Np="; Np; 1525 //"T_minGTZ(Np)="; L; 1417 1526 for(i=1;i<=size(L);i++) 1418 1527 { 1528 L[i]=std(L[i]); 1419 1529 if(idcontains(L[i],Q)==0) 1420 1530 { … … 1424 1534 P=std(P); 1425 1535 Q=std(radical(Q+P)); 1536 if(equalideals(P,Q)){return(l);} 1426 1537 list T=P,Q; 1427 1538 return(T); 1428 1539 } 1540 1541 // Prep 1542 // Computes the Prepresentation of V(N) \ V(M). 1543 // input: 1544 // ideal N (null ideal) (not necessarily radical nor maximal) 1545 // ideal M (hole ideal) (not necessarily containing N) 1546 // output: 1547 // the ((p_1,(p_11,p_1k_1)),..,(p_r,(p_r1,p_rk_r))); 1548 // the Prep of V(N) \ V(M) 1549 // Assumed to be called in the ring @R or the ring @P or a ring ring Q[a] 1550 proc Prep(ideal N, ideal M) 1551 "USAGE: Prep(N,M); 1552 ideal N (null ideal) (not necessarily radical nor maximal) 1553 ideal M (hole ideal) (not necessarily containing N) 1554 RETURN: The canonical Prepresentation of the locally closed set V(N) \ V(M) 1555 Output: [ Comp_1, .. , Comp_s ] where 1556 Comp_i=[p_i,[p_i1,..,p_is_i]] 1557 NOTE: To be called from the ring @R or @P or a ring Q[a] 1558 If the user want to call the routine Prep on a parametric 1559 ideal Q[a][x] for locally closed sets on Q[a], setglobalrings() 1560 must be called previously. 1561 KEYWORDS: Locally closed set; Canoncial form 1562 EXAMPLE: Prep; shows an example" 1563 { 1564 int te; 1565 def RR=basering; 1566 if(defined(@P)){te=1; setring(@P); ideal Np=imap(RR,N); ideal Mp=imap(RR,M);} 1567 else {te=0; def Np=N; def Mp=M;} 1568 def La=Prep0(Np,Mp); 1569 if(te==1) {setring(RR); def LL=imap(@P,La);} 1570 if(te==0){def LL=La;} 1571 return(LL); 1572 } 1429 1573 example 1430 1574 { 1431 1575 "EXAMPLE:"; echo = 2; 1432 if(defined(Grobcov::@P)){kill Grobcov::@R; kill Grobcov::@P; kill Grobcov::@RP;} 1576 short=0; 1577 ring R=0,(x,y,z),lp; 1578 ideal E=x*(x^2+y^2+z^225)*(x+y); 1579 ideal N=x*(x3)*(x+y),(y4)*(x+y); 1580 Prep(E,N); 1581 } 1582 1583 // Prep0 1584 // Computes the Prepresentation of V(N) \ V(M). 1585 // input: 1586 // ideal N (null ideal) (not necessarily radical nor maximal) 1587 // ideal M (hole ideal) (not necessarily containing N) 1588 // output: 1589 // the ((p_1,(p_11,p_1k_1)),..,(p_r,(p_r1,p_rk_r))); 1590 // the Prep of V(N) \ V(M) 1591 // Assumed to be called in a ring Q[x] (example @P) 1592 static proc Prep0(ideal N, ideal M) 1593 { 1594 int te; 1595 if (N[1]==1) 1596 { 1597 return(list(list(ideal(1),list(ideal(1))))); 1598 } 1599 int i; int j; list L0; 1600 list Ni=minGTZ(N); 1601 list prep; 1602 for(j=1;j<=size(Ni);j++) 1603 { 1604 option(redSB); 1605 Ni[j]=std(Ni[j]); 1606 } 1607 list Mij; 1608 for (i=1;i<=size(Ni);i++) 1609 { 1610 Mij=minGTZ(Ni[i]+M); 1611 for(j=1;j<=size(Mij);j++) 1612 { 1613 option(redSB); 1614 Mij[j]=std(Mij[j]); 1615 } 1616 if ((size(Mij)==1) and (equalideals(Ni[i],Mij[1])==1)){;} 1617 else 1618 { 1619 prep[size(prep)+1]=list(Ni[i],Mij); 1620 } 1621 } 1622 //"T_before="; prep; 1623 if (size(prep)==0){prep=list(list(ideal(1),list(ideal(1))));} 1624 //"T_Prep="; prep; 1625 //def Lout=CompleteA(prep,prep); 1626 //"T_Lout="; Lout; 1627 return(prep); 1628 } 1629 1630 // PtoCrep 1631 // Computes the Crepresentation from the Prepresentation. 1632 // input: 1633 // list ((p_1,(p_11,p_1k_1)),..,(p_r,(p_r1,p_rk_r))); 1634 // the Prepresentation of V(N) \ V(M) 1635 // output: 1636 // list (ideal ida, ideal idb) 1637 // the Crepresentaion of V(N) \ V(M) = V(ida) \ V(idb) 1638 // Assumed to be called in the ring @R of the ring @P or a ring Q[a] 1639 proc PtoCrep(list L) 1640 "USAGE: PtoCrep(L) 1641 L: list [ Comp_1, .. , Comp_s ] where 1642 Comp_i=[p_i,[p_i1,..,p_is_i] ], is 1643 the Prepresentation of a locally closed set V(N) \ V(M) 1644 RETURN: The canonical Crepresentation of the locally closed set 1645 [ P,Q ], a pair of radical ideals with P included in Q, 1646 representing the set V(P) \ V(Q) 1647 NOTE: To be called from the ring @R or @P or a ring Q[a] 1648 KEYWORDS: locally closed set; canoncial form 1649 EXAMPLE: PtoCrep; shows an example" 1650 { 1651 int te; 1652 def RR=basering; 1653 if(defined(@P)){te=1; setring(@P); list Lp=imap(RR,L);} 1654 else {te=0; def Lp=L;} 1655 def La=PtoCrep0(Lp); 1656 if(te==1) {setring(RR); def LL=imap(@P,La);} 1657 if(te==0){def LL=La;} 1658 return(LL); 1659 } 1660 example 1661 { 1662 "EXAMPLE:"; echo = 2; 1433 1663 ring R=0,(x,y,z),lp; 1434 1664 short=0; … … 1443 1673 } 1444 1674 1445 // Prep1446 // Computes the Prepresentation of V(N) \ V(M).1447 // input:1448 // ideal N (null ideal) (not necessarily radical nor maximal)1449 // ideal M (hole ideal) (not necessarily containing N)1450 // output:1451 // the ((p_1,(p_11,p_1k_1)),..,(p_r,(p_r1,p_rk_r)));1452 // the Prep of V(N) \ V(M)1453 // Assumed to be called in the ring @R1454 // Works on the ring @P1455 proc Prep(ideal N, ideal M)1456 "USAGE: Prep(N,M);1457 Input: ideal N (null ideal) (not necessarily radical nor maximal)1458 ideal M (hole ideal) (not necessarily containing N)1459 RETURN: The canonical Prepresentation of the locally closed set V(N) \ V(M)1460 Output: [ Comp_1, .. , Comp_s ] where1461 Comp_i=[p_i,[p_i1,..,p_is_i]]1462 NOTE: Operates in a ring R=Q[a] (a=parameters)1463 KEYWORDS: Locally closed set, Canoncial form1464 EXAMPLE: Prep; shows an example"1465 {1466 int te;1467 if (N[1]==1)1468 {1469 return(list(list(ideal(1),list(ideal(1)))));1470 }1471 int i; int j; list L0;1472 list Ni=minGTZ(N);1473 list prep;1474 for(j=1;j<=size(Ni);j++)1475 {1476 option(redSB);1477 Ni[j]=std(Ni[j]);1478 }1479 list Mij;1480 for (i=1;i<=size(Ni);i++)1481 {1482 Mij=minGTZ(Ni[i]+M);1483 for(j=1;j<=size(Mij);j++)1484 {1485 option(redSB);1486 Mij[j]=std(Mij[j]);1487 }1488 if ((size(Mij)==1) and (equalideals(Ni[i],Mij[1])==1)){;}1489 else1490 {1491 prep[size(prep)+1]=list(Ni[i],Mij);1492 }1493 }1494 //"T_abans="; prep;1495 if (size(prep)==0){prep=list(list(ideal(1),list(ideal(1))));}1496 //"T_Prep="; prep;1497 //def Lout=CompleteA(prep,prep);1498 //"T_Lout="; Lout;1499 return(prep);1500 }1501 example1502 {1503 "EXAMPLE:"; echo = 2;1504 if(defined(Grobcov::@P)){kill Grobcov::@R; kill Grobcov::@P; kill Grobcov::@RP;}1505 short=0;1506 ring R=0,(x,y,z),lp;1507 ideal E=x*(x^2+y^2+z^225);1508 ideal N=x*(x3),y4;1509 Prep(E,N);1510 }1511 1512 // PtoCrep1513 // Computes the Crepresentation from the Prepresentation.1514 // input:1515 // list ((p_1,(p_11,p_1k_1)),..,(p_r,(p_r1,p_rk_r)));1516 // the Prepresentation of V(N) \ V(M)1517 // output:1518 // list (ideal ida, ideal idb)1519 // the Crepresentaion of V(N) \ V(M) = V(ida) \ V(idb)1520 // Assumed to be called in the ring @R1521 // Works on the ring @P1522 proc PtoCrep(list L)1523 "USAGE: PtoCrep(L)1524 Input L: list [ Comp_1, .. , Comp_s ] where1525 Comp_i=[p_i,[p_i1,..,p_is_i] ], is1526 the Prepresentation of a locally closed set V(N) \ V(M)1527 RETURN: The canonical Crepresentation of the locally closed set1528 [ P,Q ], a pair of radical ideals with P included in Q,1529 representing the set V(P) \ V(Q)1530 NOTE: Operates in a ring R=Q[a] (a=parameters)1531 KEYWORDS: locally closed set, canoncial form1532 EXAMPLE: PtoCrep; shows an example"1533 {1534 int te;1535 def RR=basering;1536 if(defined(@P)){te=1; setring(@P); list Lp=imap(RR,L);}1537 else { te=0; def Lp=L;}1538 def La=PtoCrep0(Lp);1539 if(te==1) {setring(RR); def LL=imap(@P,La);}1540 if(te==0){def LL=La;}1541 return(LL);1542 }1543 example1544 {1545 "EXAMPLE:"; echo = 2;1546 if(defined(Grobcov::@P)){kill Grobcov::@R; kill Grobcov::@P; kill Grobcov::@RP;}1547 short=0;1548 ring R=0,(x,y,z),lp;1549 1550 // (P,Q) represents a locally closed set1551 ideal P=x^3+x*y^2+x*z^225*x;1552 ideal Q=y4,x*z,x^23*x;1553 1554 // Now compute the Prepresentation=1555 def L=Prep(P,Q);1556 L;1557 // Now compute the Crepresentation=1558 def J=PtoCrep(L);1559 J;1560 // Now come back recomputing the Prepresetation of the Crepresentation=1561 Prep(J[1],J[2]);1562 }1563 1564 1675 // PtoCrep0 1565 1676 // Computes the Crepresentation from the Prepresentation. … … 1570 1681 // list (ideal ida, ideal idb) 1571 1682 // the Crepresentation of V(N) \ V(M) = V(ida) \ V(idb) 1572 // Works in a ring Q[u_1,..,u_m] and is called on it.1683 // Assumed to be called in a ring Q[x] (example @P) 1573 1684 static proc PtoCrep0(list L) 1574 1685 { … … 1580 1691 { 1581 1692 option(returnSB); 1693 //"T_Lp[i]="; Lp[i]; 1582 1694 N=Lp[i][1]; 1695 Lb=Lp[i][2]; 1696 //"T_Lb="; Lb; 1583 1697 ida=intersect(ida,N); 1584 Lb=Lp[i][2];1585 1698 for(j=1;j<=size(Lb);j++) 1586 1699 { … … 1594 1707 1595 1708 // input: F a parametric ideal in Q[a][x] 1596 // output: a rComprehensive Groebner System disjoint and reduced.1709 // output: a disjoint and reduced Groebner System. 1597 1710 // It uses KapurSunWang algorithm, and with the options 1598 1711 // can compute the homogenization before (('can',0) or ( 'can',1)) 1599 1712 // and dehomogenize the result. 1600 1713 proc cgsdr(ideal F, list #) 1601 "USAGE: cgsdr(F );1714 "USAGE: cgsdr(F,options]); 1602 1715 F: ideal in Q[a][x] (a=parameters, x=variables) to be discussed. 1603 1716 To compute a disjoint, reduced Comprehensive Groebner System (CGS). … … 1610 1723 option name, value of valid options must be added to the call. 1611 1724 1612 Options:1725 OPTIONS: 1613 1726 \"can\",012: The default value is \"can\",2. In this case no 1614 1727 homogenization is done. With option (\"can\",0) the given … … 1674 1787 x=variables), and should be defined previously, and the ideal 1675 1788 defined on R. 1676 KEYWORDS: CGS , disjoint, reduced,Comprehensive Groebner System1789 KEYWORDS: CGS; disjoint; reduced; Comprehensive Groebner System 1677 1790 EXAMPLE: cgsdr; shows an example" 1678 1791 { … … 1725 1838 if(comment>=1) 1726 1839 { 1727 string("Begin cgsdr with options: ",LL);1840 " "; string("Begin cgsdr with options: ",LL); 1728 1841 } 1729 1842 int ish; … … 1731 1844 if (ish) 1732 1845 { 1733 if(comment>0){ string("The given system is homogneous");}1846 if(comment>0){" "; string("The given system is homogneous");} 1734 1847 def GS=KSW(B,LL); 1735 1848 //can=0; … … 1741 1854 { 1742 1855 // WITHOUT HOMOHGENIZING 1743 if(comment>0){ string("Option of cgsdr: do not homogenize");}1856 if(comment>0){" "; string("Option of cgsdr: do not homogenize");} 1744 1857 def GS=KSW(B,LL); 1745 1858 setglobalrings(); … … 1750 1863 { 1751 1864 // COMPUTING THE HOMOGOENIZED IDEAL 1752 if(comment> =1){string("Homogenizing the whole ideal: option can=1");}1865 if(comment>0){" "; string("Homogenizing the whole ideal: option can=1");} 1753 1866 list RRL=ringlist(RR); 1754 1867 RRL[3][1][1]="dp"; … … 1765 1878 def B1=imap(RR,B); 1766 1879 option(redSB); 1880 if(comment>0){" ";string("Basis before computing its std basis="); B1;} 1767 1881 B1=std(B1); 1882 if(comment>0){" ";string("Basis after computing its std basis="); B1;} 1768 1883 setring(RR); 1769 1884 def B2=imap(RP,B1); … … 1771 1886 else 1772 1887 { // (can=0) 1773 if(comment>0){ string("Homogenizing the basis: option can=0");}1888 if(comment>0){" "; string( "Homogenizing the basis: option can=0");} 1774 1889 def B2=B; 1775 1890 } … … 1786 1901 BH[i]=homog(BH[i],@t); 1787 1902 } 1788 if (comment> =2){string("Homogenized system = "); BH;}1903 if (comment>0){" "; string("Homogenized system = "); BH;} 1789 1904 def GSH=KSW(BH,LH); 1790 1905 setglobalrings(); … … 1812 1927 def GS=imap(RH,GSH); 1813 1928 } 1814 1815 1816 1929 setglobalrings(); 1817 1930 if(out==0) … … 2031 2144 { 2032 2145 //"T_H="; H; 2033 int i; int j; int k; int te; intvec notQ; int l; list sel; int used;2146 int i; int j; int k; int te; intvec notQ; int l; list sel; 2034 2147 intvec jtesC; 2035 2148 if ((size(H)==1) and (equalideals(H[1],ideal(1)))){return(H);} … … 2062 2175 if (te) 2063 2176 { 2064 used++;2065 2177 if ((size(notQ)==1) and (notQ[1]==0)){notQ[1]=i;} 2066 2178 else{notQ[size(notQ)+1]=i;} … … 2100 2212 } 2101 2213 } 2102 //"T_Q1_0="; Q1;2103 2214 sel=selectminideals(Q1); 2104 2215 kill nQ1; list nQ1; … … 2110 2221 } 2111 2222 if(size(Q1)==0){Q1=ideal(1),ideal(1);} 2112 //"T_Q1_1="; Q1;2113 //if(used>0){string("addpartfine was ", used, " times used");}2114 2223 return(Q1); 2115 2224 } 2116 2225 2117 2118 // Auxiliary routine for grobcov: ideal F is assumed to be homogeneous 2226 // Auxiliary rutine for gcover 2227 // Deciding if combine is needed 2228 // input: list LCU=( (basis1, p_1, (p11,..p1s1)), .. (basisr, p_r, (pr1,..prsr)) 2229 // output: (tes); if tes==1 then combine is needed, else not. 2230 static proc needcombine(list LCU,ideal N) 2231 { 2232 //"Deciding if combine is needed";; 2233 ideal BB; 2234 int tes=0; int m=1; int j; int k; poly sp; 2235 while((tes==0) and (m<=size(LCU[1][1]))) 2236 { 2237 j=1; 2238 while((tes==0) and (j<=size(LCU))) 2239 { 2240 k=1; 2241 while((tes==0) and (k<=size(LCU))) 2242 { 2243 if(j!=k) 2244 { 2245 sp=pnormalf(pspol(LCU[j][1][m],LCU[k][1][m]),LCU[k][2],N); 2246 if(sp!=0){tes=1;} 2247 } 2248 k++; 2249 } 2250 j++; 2251 } 2252 if(tes){break;} 2253 m++; 2254 } 2255 return(tes); 2256 } 2257 2258 // Auxiliary routine 2259 // precombine 2260 // input: L: list of ideals (works in @P) 2261 // output: F0: ideal of polys. F0[i] is a poly in the intersection of 2262 // all ideals in L except in the ith one, where it is not. 2263 // L=(p1,..,ps); F0=(f1,..,fs); 2264 // F0[i] \in intersect_{j#i} p_i 2265 static proc precombine(list L) 2266 { 2267 int i; int j; int tes; 2268 def RR=basering; 2269 setring(@P); 2270 list L0; list L1; list L2; list L3; ideal F; 2271 L0=imap(RR,L); 2272 L1[1]=L0[1]; L2[1]=L0[size(L0)]; 2273 for (i=2;i<=size(L0)1;i++) 2274 { 2275 L1[i]=intersect(L1[i1],L0[i]); 2276 L2[i]=intersect(L2[i1],L0[size(L0)i+1]); 2277 } 2278 L3[1]=L2[size(L2)]; 2279 for (i=2;i<=size(L0)1;i++) 2280 { 2281 L3[i]=intersect(L1[i1],L2[size(L0)i]); 2282 } 2283 L3[size(L0)]=L1[size(L1)]; 2284 for (i=1;i<=size(L3);i++) 2285 { 2286 option(redSB); L3[i]=std(L3[i]); 2287 } 2288 for (i=1;i<=size(L3);i++) 2289 { 2290 tes=1; j=0; 2291 while((tes) and (j<size(L3[i]))) 2292 { 2293 j++; 2294 option(redSB); 2295 L0[i]=std(L0[i]); 2296 if(reduce(L3[i][j],L0[i])!=0){tes=0; F[i]=L3[i][j];} 2297 } 2298 if (tes){"ERROR a polynomial in all p_j except p_i was not found";} 2299 } 2300 setring(RR); 2301 def F0=imap(@P,F); 2302 return(F0); 2303 } 2304 2305 // Auxiliary routine 2306 // combine 2307 // input: a list of pairs ((p1,P1),..,(pr,Pr)) where 2308 // ideal pi is a prime component 2309 // poly Pi is the polynomial in Q[a][x] on V(pi)\ V(Mi) 2310 // (p1,..,pr) are the prime decomposition of the lppsegment 2311 // list crep =(ideal ida,ideal idb): the Crep of the segment. 2312 // list Pci of the intersecctions of all pj except the ith one 2313 // output: 2314 // poly P on an open and dense set of V(p_1 int ... p_r) 2315 static proc combine(list L, ideal F) 2316 { 2317 // ATTENTION REVISE AND USE Pci and F 2318 int i; poly f; 2319 f=0; 2320 for(i=1;i<=size(L);i++) 2321 { 2322 f=f+F[i]*L[i][2]; 2323 } 2324 // f=elimconstfac(f); 2325 f=primepartZ(f); 2326 return(f); 2327 } 2328 2329 // Central routine for grobcov: ideal F is assumed to be homogeneous 2119 2330 // gcover 2120 2331 // input: ideal F: a generating set of a homogeneous ideal in Q[a][x] … … 2127 2338 { 2128 2339 int i; int j; int k; ideal lpp; list GPi2; list pairspP; ideal B; int ti; 2129 int i1; int tes; int j1; int selind; int i2; int m;2340 int i1; int tes; int j1; int selind; int i2; //int m; 2130 2341 list prep; list crep; list LCU; poly p; poly lcp; ideal FF; 2131 2342 list lpi; … … 2172 2383 list S; 2173 2384 poly sp; 2174 ideal BB;2175 2385 for (i=1;i<=size(GP);i++) 2176 2386 { … … 2195 2405 } 2196 2406 //"Deciding if combine is needed"; 2197 kill BB;2198 ideal BB;2199 tes=1; m=1;2200 while((tes) and (m<=size(LCU[1][1])))2201 {2202 j=1;2203 while((tes) and (j<=size(LCU)))2204 {2205 k=1;2206 while((tes) and (k<=size(LCU)))2207 {2208 if(j!=k)2209 {2210 sp=pnormalf(pspol(LCU[j][1][m],LCU[k][1][m]),LCU[k][2],N);2211 if(sp!=0){tes=0;}2212 }2213 k++;2214 } //setglobalrings();2215 if(tes)2216 {2217 BB[m]=LCU[j][1][m];2218 }2219 j++;2220 }2221 if(tes==0){break;}2222 m++;2223 }2224 2407 crep=PtoCrep(prep); 2225 if(tes==0) 2408 if(size(LCU)>1){tes=1;} 2409 else 2410 { 2411 tes=0; 2412 for(k=1;k<=size(B);k++){B[k]=pnormalf(B[k],crep[1],crep[2]);} 2413 } 2414 // tes=needcombine(LCU,N); 2415 // if(tes==1){" ";string("combine is needed for segment ",i);" ";} 2416 //crep=PtoCrep(prep); 2417 if(tes) 2226 2418 { 2227 2419 // combine is needed … … 2231 2423 LL[j]=LCU[j][2]; 2232 2424 } 2233 if (size(LCU)>1) 2234 { 2235 FF=precombint(LL); 2236 } 2425 FF=precombine(LL); 2237 2426 for (k=1;k<=size(lpp);k++) 2238 2427 { … … 2242 2431 L[j]=list(LCU[j][2],LCU[j][1][k]); 2243 2432 } 2244 if (size(LCU)>1) 2245 { 2246 B[k]=combine(L,FF); 2247 } 2248 else{B[k]=L[1][2];} 2249 } 2250 } 2251 else{B=BB;} 2433 B[k]=combine(L,FF); 2434 } 2435 } 2252 2436 for(j=1;j<=size(B);j++) 2253 2437 { … … 2277 2461 // N is the null conditions ideal (if desired) 2278 2462 // W is the ideal of nonnull conditions (if desired) 2279 // The value of \"can\" is 1 by default and can be set to 0 if we do not2463 // The value of \"can\" is 1 by default and can be set to 0 if we do not 2280 2464 // need to obtain the canonical GC, but only a GC. 2281 2465 // The value of \"ext\" is 0 by default and so the generic representation … … 2292 2476 // given by the lpp, the basis, and the Prepresentation of the segment 2293 2477 proc grobcov(ideal F,list #) 2294 "USAGE: grobcov(F); This is the fundamental routine of the 2478 "USAGE: grobcov(F,options]); 2479 This is the fundamental routine of the 2295 2480 library. It computes the Groebner cover of a parametric ideal. See 2296 2481 Montes A., Wibmer M., … … 2307 2492 option name, value of valid options must be added to the call. 2308 2493 2309 Options:2494 OPTIONS: 2310 2495 \"null\",ideal E: The default is (\"null\",ideal(0)). 2311 2496 \"nonnull\",ideal N: The default is (\"nonnull\",ideal(1)). … … 2332 2517 \"comment\" higher will provide information about the 2333 2518 development of the computation. 2519 \"showhom\",01: The default is (\"showhom\",0). Setting 2520 \"showhom\",1 will output the set of homogenized lpp of 2521 each segment as last element. 2334 2522 One can give none or whatever of these options. 2335 2523 RETURN: The list 2524 ( 2525 (lpp_1,basis_1,segment_1), 2526 ... 2527 (lpp_s,basis_s,segment_s) 2528 ) 2529 optionally 2336 2530 ( 2337 2531 (lpp_1,basis_1,segment_1,lpph_1), … … 2343 2537 set of lpp of the reduced Groebner basis for each point 2344 2538 of the segment. 2539 With option (\"showhom\",1) the lpph will be shown: 2345 2540 The lpph corresponds to the lpp of the homogenized ideal 2346 2541 and is different for each segment. It is given as a string, … … 2379 2574 x=variables), and should be defined previously. The ideal must 2380 2575 be defined on R. 2381 KEYWORDS: Groebner cover , parametric ideal, canonical,discussion of2382 parametric ideal .2576 KEYWORDS: Groebner cover; parametric ideal; canonical; discussion of 2577 parametric ideal 2383 2578 EXAMPLE: grobcov; shows an example" 2384 2579 { … … 2390 2585 setglobalrings(); 2391 2586 list L0=#; 2587 list Se; 2392 2588 int out=0; 2589 int showhom=0; 2590 int hom; 2393 2591 L0[size(L0)+1]="res"; L0[size(L0)+1]=ideal(1); 2394 2592 // default options … … 2421 2619 { 2422 2620 if (L0[2*i1]=="comment"){comment=L0[2*i];} 2621 else 2622 { 2623 if (L0[2*i1]=="showhom"){showhom=L0[2*i];} 2624 } 2423 2625 } 2424 2626 } … … 2445 2647 LL[11]="ext"; LL[12]=extop; 2446 2648 LL[13]="rep"; LL[14]=repop; 2649 LL[15]="showhom"; LL[16]=showhom; 2447 2650 if (comment>=1) 2448 2651 { … … 2452 2655 def S=gcover(F,LL); 2453 2656 // NOW extend 2657 if (size(S[1])==5){hom=1;} 2454 2658 if(extop) 2455 2659 { 2456 2660 S=extend(S,LL); 2457 2661 } 2662 // NOW representation of the segments by option repop 2663 list Si; list nS; 2664 if (repop==0) 2665 { 2666 for(i=1;i<=size(S);i++) 2667 { 2668 Si=list(S[i][1],S[i][2],S[i][3]); 2669 if(hom){Si[size(Si)+1]=S[i][4];} 2670 nS[size(nS)+1]=Si; 2671 } 2672 S=nS; 2673 } 2458 2674 else 2459 2675 { 2460 // NOW representation of the segments by option repop 2461 list Si; list nS; 2462 if(repop==0) 2676 if (repop==1) 2463 2677 { 2464 2678 for(i=1;i<=size(S);i++) 2465 2679 { 2466 Si=list(S[i][1],S[i][2],S[i][3],S[i][5]); 2680 Si=list(S[i][1],S[i][2],S[i][4]); 2681 if(hom){Si[size(Si)+1]=S[i][5];} 2467 2682 nS[size(nS)+1]=Si; 2468 2683 } 2469 kill S; 2470 def S=nS; 2684 S=nS; 2471 2685 } 2472 2686 else 2473 2687 { 2474 if(repop==1) 2475 { 2476 for(i=1;i<=size(S);i++) 2477 { 2478 Si=list(S[i][1],S[i][2],S[i][4],S[i][5]); 2479 nS[size(nS)+1]=Si; 2480 } 2481 kill S; 2482 def S=nS; 2483 } 2484 else 2485 { 2486 for(i=1;i<=size(S);i++) 2487 { 2488 Si=list(S[i][1],S[i][2],S[i][3],S[i][4],S[i][5]); 2489 nS[size(nS)+1]=Si; 2490 } 2491 kill S; 2492 def S=nS; 2688 for(i=1;i<=size(S);i++) 2689 { 2690 Si=list(S[i][1],S[i][2],S[i][3],S[i][4]); 2691 if(hom){Si[size(Si)+1]=S[i][5];} 2692 nS[size(nS)+1]=Si; 2493 2693 } 2494 2694 } … … 2500 2700 } 2501 2701 if(defined(@P)==1){kill @R; kill @P; kill @RP;} 2702 if (hom and showhom==0) 2703 { 2704 for(i=1;i<=size(S);i++) 2705 { 2706 Se=S[i]; 2707 Se=delete(Se,size(Se)); 2708 S[i]=Se; 2709 } 2710 } 2502 2711 return(S); 2503 2712 } … … 2523 2732 // Can be called from the top 2524 2733 proc extend(list GC, list #); 2525 "USAGE: extend(GC); The default option of grobcov provides 2734 "USAGE: extend(GC[,options]); 2735 The default option of grobcov provides 2526 2736 the bases in generic representation (the Iregular functions 2527 of the bases ar agiven by a single polynomial. It can specialize2737 of the bases are given by a single polynomial. It can specialize 2528 2738 to zero for some points of the segments, but in general, it 2529 2739 is sufficient for many pouposes. Nevertheless the Iregular 2530 functions allow a full representation given b ey a set of2740 functions allow a full representation given by a set of 2531 2741 polynomials specializing to the value of the function (after normalization) 2532 2742 or to zero, but at least one of the polynomials specializes to nonzero. … … 2552 2762 x=variables), and should be defined previously. The ideal must 2553 2763 be defined on R. 2554 KEYWORDS: Groebner cover , parametric ideal,canonical, discussion of2555 parametric ideal , full representation.2764 KEYWORDS: Groebner cover; parametric ideal; canonical, discussion of 2765 parametric ideal; full representation 2556 2766 EXAMPLE: extend; shows an example" 2557 2767 { … … 2567 2777 i++; 2568 2778 } 2779 int hom; 2569 2780 if(m==1){"Warning! grobcov has already extended bases"; return(S);} 2570 if(size(GC[1])!=5){"Warning! extend make sense only when grobcov has been called with options 'rep',2,'ext',0"; " "; return();} 2781 if(not(size(GC[1])==5 or (size(GC[1])==4 and typeof(GC[1][4])=="list"))) 2782 {"Warning! extend make sense only when grobcov has been called with options 'rep',2,'ext',0"; " "; return();} 2783 if(size(GC[1])==5){hom=1;} 2571 2784 int repop=0; 2572 2785 int start3=timer; … … 2634 2847 for(i=1;i<=size(S);i++) 2635 2848 { 2636 Si=list(S[i][1],S[i][2],S[i][3],S[i][5]); 2849 Si=list(S[i][1],S[i][2],S[i][3]); 2850 if(hom){Si[size(Si)+1]=S[i][5];} 2637 2851 nS[size(nS)+1]=Si; 2638 2852 } … … 2645 2859 for(i=1;i<=size(S);i++) 2646 2860 { 2647 Si=list(S[i][1],S[i][2],S[i][4],S[i][5]); 2861 Si=list(S[i][1],S[i][2],S[i][4]); 2862 if(hom){Si[size(Si)+1]=S[i][5];} 2648 2863 nS[size(nS)+1]=Si; 2649 2864 } … … 2654 2869 for(i=1;i<=size(S);i++) 2655 2870 { 2656 Si=list(S[i][1],S[i][2],S[i][3],S[i][4],S[i][5]); 2871 Si=list(S[i][1],S[i][2],S[i][3],S[i][4]); 2872 if(hom){Si[size(Si)+1]=S[i][5];} 2657 2873 nS[size(nS)+1]=Si; 2658 2874 } 2659 2660 2875 } 2661 2876 } … … 2670 2885 short=0; 2671 2886 ideal S=a0*x^2+b0*x+c0, 2672 2887 a1*x^2+b1*x+c1; 2673 2888 def GCS=grobcov(S,"rep",2); 2674 2889 GCS; … … 2718 2933 return(fr); 2719 2934 } 2720 2721 // Auxiliary routine2722 // deltai2723 // input:2724 // int i:2725 // list LPr: (p1,..,pr) of prime components of an ideal in Q[a]2726 // output:2727 // list (fr,fnr) of two polynomials that are equal on V(pi)2728 // and fr=0 on V(P) \ V(pi), and fnr is nonzero on V(pj) for all j.2729 static proc deltai(int i, list LPr)2730 {2731 def RR=basering;2732 setring(@P);2733 def LP=imap(RR,LPr);2734 int j; poly p;2735 def F=ideal(1);2736 poly f;2737 poly fn;2738 ideal LPi;2739 for (j=1;j<=size(LP);j++)2740 {2741 if (j!=i)2742 {2743 F=idint(F,LP[j]);2744 }2745 }2746 p=0; j=1;2747 while ((p==0) and (j<=size(F)))2748 {2749 LPi=LP[i];2750 attrib(LPi,"isSB",1);2751 p=reduce(F[j],LPi);2752 j++;2753 }2754 f=F[j1];2755 fn=nonzerodivisor(f,LP);2756 setring(RR);2757 def fr=imap(@P,f);2758 def fnr=imap(@P,fn);2759 return(list(fr,fnr));2760 }2761 2762 // Auxiliary routine2763 // combine2764 // input: a list of pairs ((p1,P1),..,(pr,Pr)) where2765 // ideal pi is a prime component2766 // poly Pi is the polynomial in Q[a][x] on V(pi)\ V(Mi)2767 // (p1,..,pr) are the prime decomposition of the lppsegment2768 // list crep =(ideal ida,ideal idb): the Crep of the segment.2769 // list Pci of the intersecctions of all pj except the ith one2770 // output:2771 // poly P on an open and dense set of V(p_1 int ... p_r)2772 static proc combine(list L, ideal F)2773 {2774 // ATTENTION REVISE AND USE Pci and F2775 int i; poly f;2776 f=0;2777 for(i=1;i<=size(L);i++)2778 {2779 f=f+F[i]*L[i][2];2780 }2781 // f=elimconstfac(f);2782 f=primepartZ(f);2783 return(f);2784 }2785 2786 2935 2787 2936 //Auxiliary routine … … 3013 3162 M1=intersect(T0[2],T1[2]); 3014 3163 } 3015 T=list(ci,PtoCrep (Prep(N1,M1)));3164 T=list(ci,PtoCrep0(Prep0(N1,M1))); 3016 3165 LL[size(LL)+1]=T; 3017 3166 if(equalideals(T[2][1],ideal(1))){te=1; break;} … … 3046 3195 return(T); 3047 3196 } 3048 3049 3197 3050 3198 // Auxiliary routine … … 3237 3385 } 3238 3386 3239 // // Auxiliary routine3240 // // NonNull: returns 1 if the poly f is nonnull on V(E) \ V(N), 0 otherwise.3241 // // Input:3242 // // f: polynomial3243 // // E: null ideal3244 // // N: nonnull ideal3245 // // Output:3246 // // 1 if f is nonnul in V(P) \ V(Q),3247 // // 0 if it has zeroes inside.3248 // // Works in @P3249 // proc NonNull(poly f, ideal E, ideal N)3250 // {3251 // int te=1; int i;3252 // def RR=basering;3253 // setring(@P);3254 // def fp=imap(RR,f);3255 // def Ep=imap(RR,E);3256 // def Np=imap(RR,N);3257 // ideal H;3258 // ideal Ef=Ep+fp;3259 // for (i=1;i<=size(Np);i++)3260 // {3261 // te=radicalmember(Np[i],Ef);3262 // if (te==0){break;}3263 // }3264 // setring(RR);3265 // return(te);3266 // }3267 3268 3387 // Auxiliary routine 3269 3388 // selectextendcoef … … 3427 3546 } 3428 3547 3429 3430 3548 // Auxiliary routine 3431 3549 // precombint … … 3474 3592 return(F0); 3475 3593 } 3476 3477 3594 3478 3595 // Auxiliary routine … … 3616 3733 int i; 3617 3734 def L=#; 3618 if 3735 if(size(L)>0) 3619 3736 { 3620 3737 for (i=1;i<=size(L) div 2;i++) … … 3658 3775 // Auxiliary routine 3659 3776 // sqf 3660 // This is for releases of Singular before 3513661 // proc sqf(poly f)3662 // {3663 // def RR=basering;3664 // setring(@P);3665 // def ff=imap(RR,f);3666 // def G=sqrfree(ff);3667 // poly fff=1;3668 // int i;3669 // for (i=1;i<=size(G);i++)3670 // {3671 // fff=fff*G[i];3672 // }3673 // setring(RR);3674 // def ffff=imap(@P,fff);3675 // return(ffff);3676 // }3677 3678 // Auxiliary routine3679 // sqf3680 3777 static proc sqf(poly f) 3681 3778 { … … 3688 3785 return(ffff); 3689 3786 } 3690 3691 3787 3692 3788 // Auxiliary routine … … 3803 3899 if (equalideals(Gr,ideal(0))==1){GrHi=H[i];} 3804 3900 else {GrHi=Gr,H[i];} 3805 // else {for (j=1;j<=size(Gr);j++){GrHi[size(GrHi)+1]=Gr[j]*H[i];}}3806 3901 if(comment>1){"Call to KSW with arguments "; GM; GrHi; Nh;} 3807 3902 KS=KSW0(GM,GrHi,Nh,comment); … … 3981 4076 int n=size(L); 3982 4077 if(n<2){return(A);} 4078 intvec w; 3983 4079 for(i=1;i<=size(L);i++) 3984 4080 { … … 3990 4086 { 3991 4087 L[i][2]=L[j][2]; 4088 w[size(w)+1]=j; 3992 4089 } 3993 4090 } 4091 } 4092 } 4093 if(size(w)>0) 4094 { 4095 for(i=1; i<=size(w);i++) 4096 { 4097 j=w[size(w)+1i]; 4098 L=elimfromlist(L, j); 3994 4099 } 3995 4100 } … … 4013 4118 } 4014 4119 if(equalideals(T,ideal(1))==0){L[size(L)+1]=list(std(T),ideal(1))}; 4015 //string("T_n=",n," new n0",size(L));4016 4120 return(L); 4017 4121 } 4018 4122 4019 // Input: list(A) 4020 // A is a list of locally closed sets in Crep. A=[[P1,Q1],[P2,Q2],..,[Pr,Qr]] 4021 // whose union is a constructible set from 4022 // Output list [Lev,C]: 4023 // where Lev is the Crep of the canonical first level of A, and 4024 // C is the complement of the first level Lev wrt to the closure of A. The elements are given in Crep, 4123 // input list A=[[p1,q1],...,[pn,qn]] : 4124 // the list of segments of a constructible set S, where each [pi,qi] is given in Crepresentation 4125 // output list [topA,C] 4126 // where topA is the closure of A 4127 // C is the list of segments of the complement of A given in Crepresentation 4025 4128 static proc FirstLevel(list A) 4026 4129 { … … 4028 4131 list T=zeroone(n); 4029 4132 ideal P; ideal Q; 4030 list Cb; ideal Cc= ideal(1);4133 list Cb; ideal Cc=1; 4031 4134 int i; int j; 4032 4135 intvec t; 4033 ideal Z=ideal(1);4136 ideal topA=1; 4034 4137 list C; 4035 for(i=1;i<=size(A);i++) 4036 { 4037 Z=intersect(Z,A[i][1]); 4038 } 4138 for(i=1;i<=n;i++) 4139 { 4140 topA=intersect(topA,A[i][1]); 4141 } 4142 //topA=std(topA); 4039 4143 for(i=2; i<=size(T);i++) 4040 4144 { 4041 4145 t=T[i]; 4042 P=ideal(1); Q=ideal(0); 4146 //"T_t="; t; 4147 P=0; Q=1; 4043 4148 for(j=1;j<=n;j++) 4044 4149 { 4045 4150 if(t[n+1j]==1) 4046 4151 { 4047 Q=Q+A[j][2];4152 P=P+A[j][2]; 4048 4153 } 4049 4154 else 4050 4155 { 4051 P=intersect(P,A[j][1]); 4052 } 4053 } 4054 //"T_Q="; Q; "T_P="; P; 4055 Cb=Crep(Q,P); 4156 Q=intersect(Q,A[j][1]); 4157 } 4158 } 4159 Cb=Crep0(P,Q); 4056 4160 //"T_Cb="; Cb; 4057 if(Cb[1][1]<>1) 4058 { 4059 C[size(C)+1]=Cb; 4060 Cc=intersect(Cc,Cb[1]); 4061 } 4062 } 4063 list Lev=list(Z,Cc); //Crep(Z,Cc); 4161 if(size(Cb)!=0) 4162 { 4163 if( Cb[1][1]<>1) 4164 { 4165 C[size(C)+1]=Cb; 4166 } 4167 } 4168 } 4064 4169 if(size(C)>1){C=SimplifyUnion(C);} 4065 return(list(Lev,C)); 4066 } 4067 4068 // Input: list(A) 4069 // A is a list of locally closed sets in Crep. A=[[P1,Q1],[P2,Q2],..,[Pr,Qr]] 4070 // whose union is a constructible set from which the algorithm computes its canonical form. 4170 return(list(topA,C)); 4171 } 4172 4173 // Input: 4071 4174 // Output: 4072 // list [L,C] where 4073 // where L is the list of canonical levels of A, 4074 // and C is the list of canonical levels of the complement of A wrt to the closure of A. 4075 // All locally closed sets are given in Crep. 4076 // L=[[1,[p1,p2],[3,[p3,p4],..,[2r1,[p_{2r1},p_2r]]]] is the list of levels of A in Crep. 4077 // C=[[2,p2,p3],[4,[p4,p5],..,[2s,[p_{2s},p_{2s+1}]]] is the list of levels of C, 4078 // the complement of S wrt the closure of A, in Crep. 4079 proc ConsLevels(list A) 4080 "USAGE: ConsLevels(A); 4081 A is a list of locally closed sets in Crep. A=[[P1,Q1],[P2,Q2],..,[Pr,Qr]] 4082 whose union is a constructible set from which the algorithm computes its 4083 canonical form. 4084 RETURN: A list with two elements: [L,C] 4085 where L is the list of canonical levels of A, and C is the list of canonical 4086 levels of the 4087 complement of A wrt to the closure of A. 4088 All locally closed sets are given in Crep. 4089 L=[[1,[p1,p2],[3,[p3,p4],..,[2r1,[p_{2r1},p_2r]]]] 4090 C=[[2,p2,p3],[4,[p4,p5],..,[2s,[p_{2s},p_{2s+1}]]] 4091 NOTE: The basering R, must be of the form Q[a] 4092 KEYWORDS: locally closed set, constructible set 4093 EXAMPLE: ConsLevels; shows an example" 4094 { 4095 list L; list C; int i; 4096 list B; list T; 4175 static proc ConstoPrep(list L) 4176 { 4177 list L1; 4178 int i; int j; 4179 list aux; 4180 for(i=1;i<=size(L);i++) 4181 { 4182 aux=Prep(L[i][2][1],L[i][2][2]); 4183 L1[size(L1)+1]=list(L[i][1],aux); 4184 } 4185 return(L1); 4186 } 4187 4188 // Input: 4189 // list A = [[P1,Q1], .. [Pn,Qn]] 4190 // A constructible set as union of locally closed sets represented by pairs of ideals 4191 // Output: 4192 // list L =[p1,p2,p3,...,pk] 4193 // where pi is the ideal of the closure of level i alternatively of A or its complement 4194 // Note: the levels of A are [p1,p2], [p3,p4], [p5,p6],... 4195 // the levels of C are [p2,p3],[p4,p5], ... 4196 // expressed in Crepresentation 4197 // Assumed to be called in the ring @R or the ring @P or a ring ring Q[a] 4198 proc ConsLevels(list A0) 4199 "USAGE: ConsLevels(L); 4200 where L=[[P1,Q1],...,[Ps,Qs]] 4201 is a list of lists of of pairs of ideals defined over Grobcov::@P, 4202 or over a ring without parameters and represents the constructible set 4203 S=V(P1)\ V(Q1) u ... u V(Ps)\ V(Qs) 4204 RETURN: The list of ideals [a1,a2,...,at] representing the closures of the 4205 canonical levels of S and its complement C wrt to the closure of S. 4206 The canonical levels of S are represented by theirs Crep. So we have: 4207 Levels of S: [a1,a2],[a3,a4],... 4208 Levels of C: [a2,a3],[a4,a5],... 4209 S=V(a1)\ V(a2) u V(a3)\ V(a4) u ... 4210 C=V(a2\ V(a3) u V(a4)\ V(a5) u ... 4211 It is used internally by locus procedure. 4212 NOTE: It is described in the paper 4213 J.M. Brunat, A. Montes 4214 Computing the canonical representation of constructible sets. 4215 Math. Comput. Sci. (2016) 19: 165178. 4216 KEYWORDS: constructible set; locally closed set; canonical form 4217 EXAMPLE: example ConsLevels; shows an example" 4218 { 4219 def RR=basering; 4220 int te; 4221 if(defined(@P)){te=1; setring(@P); list A=imap(RR,A0);} 4222 else {te=0; def A=A0;} 4223 4224 list L; list C; 4225 list B; list T; int i; 4097 4226 for(i=1; i<=size(A);i++) 4098 4227 { 4099 T=Crep (A[i][1],A[i][2]);4228 T=Crep0(A[i][1],A[i][2]); 4100 4229 B[size(B)+1]=T; 4101 4230 } 4102 int level=0;4103 4231 list K; 4104 4232 while(size(B)>0) 4105 4233 { 4106 level++;4107 4234 K=FirstLevel(B); 4108 if(level mod 2==1){L[size(L)+1]=list(level,K[1]);} 4109 else{C[size(C)+1]=list(level,K[1]);} 4110 //"T_L="; L; 4111 //"T_C="; C; 4235 //"T_K="; K; 4236 L[size(L)+1]=K[1]; 4112 4237 B=K[2]; 4113 //"T_size(B)="; size(B); 4114 } 4115 return(list(L,C)); 4238 } 4239 L[size(L)+1]=ideal(1); 4240 if(te==1) {setring(RR); def LL=imap(@P,L);} 4241 if(te==0){def LL=L;} 4242 return(LL); 4116 4243 } 4117 4244 example 4118 { "EXAMPLE:"; echo = 2; 4119 // Example of ConsLevels with nice geometrical interpretetion and 2 levels 4120 4121 if (defined(R)){kill R;} 4122 if (defined(@P)){kill @P; kill @R; kill @RP;} 4123 4124 ring R=0,(x,y,z),lp; 4125 short=0; 4126 ideal P1=x*(x^2+y^2+z^21); 4127 ideal Q1=z,x^2+y^21; 4128 ideal P2=y,x^2+z^21; 4129 ideal Q2=z*(z+1),y,x*(x+1); 4130 4131 list Cr1=Crep(P1,Q1); 4132 list Cr2=Crep(P2,Q2); 4133 4134 list L=list(Cr1,Cr2); 4135 L; 4136 // ConsLevels(L)= 4137 ConsLevels(L); 4138 4139 // 4140 // Begin Problem S93 4141 // Automatic theorem proving 4142 // Generalized SteinerLehmus Theorem 4143 // A.Montes, T.Recio 4144 4145 // Bisector of A(1,0) = Bisector of B(1,0) varying C(a,b) 4146 4147 if(defined(R1)){kill R1;} 4148 ring R1=(0,a,b),(x1,y1,x2,y2,p,r),dp; 4149 ideal S93=(a+1)^2+b^2(p+1)^2, 4150 (a1)^2+b^2(1r)^2, 4151 a*y1b*x1y1+b, 4152 a*y2b*x2+y2b, 4153 2*y1+b*x1(a+p)*y1+b, 4154 2*y2+b*x2(a+r)*y2b, 4155 (x1+1)^2+y1^2(x21)^2y2^2; 4156 4245 { 4246 "EXAMPLE:"; echo = 2; 4247 ring R=0,(x,y,z),lp; 4157 4248 short=0; 4158 def GC93=grobcov(S93,"nonnull",ideal(b),"rep",1); 4159 //"grobcov(S93,'nonnull',ideal(b),"rep",1)="; GC93; 4160 4161 list L0; 4162 for(int i=1;i<=size(GC93);i++) 4163 { 4164 L0[size(L0)+1]=GC93[i][3]; 4165 } 4166 4167 def GC93ab=grobcov(S93,"nonnull",ideal(a*b),"rep",1); 4168 4169 ring RR=0,(a,b),lp; 4170 4171 list L1; 4172 L1=imap(R1,L0); 4173 // L1=Total elements of the grobcov(S93,'nonnull',ideal(b),'rep',1) to be added= 4174 L1; 4175 4176 // Adding all the elements of grobcov(S93,'nonnull',ideal(b),'rep',1)= 4177 ConsLevels(L1); 4249 ideal P1=x*(x^2+y^2+z^21); 4250 ideal Q1=z,x^2+y^21; 4251 ideal P2=y,x^2+z^21; 4252 ideal Q2=z*(z+1),y,x*(x+1); 4253 4254 list Cr1=Crep(P1,Q1); 4255 list Cr2=Crep(P2,Q2); 4256 4257 list L=list(Cr1,Cr2); 4258 L; 4259 ConsLevels(L); 4260 } 4261 4262 // Converts the output of ConsLevels, given by the set of closures of the Levels of the constructible S 4263 // to an expression where the Levels are apparent. 4264 // Input: The ouput of ConsLevels of the form 4265 // [A1,A2,..,Ak], where the Ai's are the closures of the levels. 4266 // Output: An expression of the form 4267 // L1=[[1,[A1,A2]],[3,[A3,A4]],..,[2l1,[A_{2l1},A_{2l}]]] the list of Levels of S 4268 proc ConsLevelsToLevels(list L) 4269 "USAGE: ConsLevelsToLevels(list L); 4270 The input list L must be the output of the call to 4271 the routine 'ConsLevels' of a constructible set: 4272 [A1,A2,..,Ak], where the Ai's are the closures of the levels. 4273 'ConsLevels' outputs the closures of the levels 4274 of the constructible and of its complements. 4275 This routine selects the levels of the 4276 constructible set. 4277 RETURN: The levels of the constructible set: 4278 Lc=[[1,[A1,A2]],[3,[A3,A4]],..,[2l1,[A_{2l1},A_{2l}]]] the list of 4279 Levels of S 4280 NOTE: It must be called to the output of the 'ConsLevels' routine. 4281 KEYWORDS: constructible sets 4282 EXAMPLE: ConsLevelsToLevels shows an example" 4283 { 4284 int n=size(L) div 2; 4285 int i; 4286 list L1; list L2; 4287 for(i=1; i<=n;i++) 4288 { 4289 L1[size(L1)+1]=list(2*i1,list(L[2*i1],L[2*i])); 4290 } 4291 return(L1); 4292 } 4293 example 4294 { 4295 "EXAMPLE:"; echo = 2; 4296 ring R=0,(x,y,z),lp; 4297 short=0; 4298 ideal P1=(x^2+y^2+z^21); 4299 ideal Q1=z,x^2+y^21; 4300 ideal P2=y,x^2+z^21; 4301 ideal Q2=z*(z+1),y,x*(x+1); 4302 ideal P3=x; 4303 ideal Q3=5*z4,5*y3,x; 4304 4305 list Cr1=Crep(P1,Q1); 4306 list Cr2=Crep(P2,Q2); 4307 list Cr3=Crep(P3,Q3); 4308 4309 list L=list(Cr1,Cr2,Cr3); 4310 L; 4311 4312 def LL=ConsLevels(L); 4313 LL; 4314 ConsLevelsToLevels(LL); 4178 4315 } 4179 4316 … … 4247 4384 // ideal N: the holes of the component 4248 4385 // Output: 4249 // int d: the dimension of B on the variables reduced by the holes.4386 // int d: the dimension of B on the variables (antiimage). 4250 4387 // if d>0 then the component is 'Normal' 4251 4388 // if d==0 then the component is 'Singular' … … 4268 4405 } 4269 4406 if(F!=0){Fenv=1;} 4270 //"T_B="; B; "E="; E; "N="; N;4271 4407 list L0; 4272 4408 if(dimP0(E)==0){L0=2,"Normal";} // 2 es fals pero ha de ser >0 encara que sigui 0 … … 4275 4411 if(version==0) 4276 4412 { 4277 // "T_B="; B; //Computing std(B+E,plex(x,y,x1,..xn)) one can detect if there is a first part4413 // Computing std(B+E,plex(x,y,x1,..xn)) one can detect if there is a first part 4278 4414 // independent of parameters giving the variables with dimension 0 4279 4415 dd=indepparameters(B); 4280 if (dd==1){d=0; L0=d,string(B) ,determineF(B,F,E);}4416 if (dd==1){d=0; L0=d,string(B);} // ,determineF(B,F,E) 4281 4417 else{d=1; L0=2,"Normal";} 4282 4418 } … … 4340 4476 { 4341 4477 //" ";string("Antiimage of Special component = ", GGG); 4342 if(Fenv==1)4343 {4344 L0[3]=determineF(GGG,F,E);4345 }4346 4478 } 4347 4479 else … … 4351 4483 } 4352 4484 } 4485 //"T_L0="; L0; 4353 4486 return(L0); 4354 4487 } … … 4381 4514 ideal M=std(AA+FH); 4382 4515 def rh=reduce(EH,M); 4516 //"T_AA="; AA; "T_FH="; FH; "T_EH="; EH; "T_rh="; rh; 4383 4517 if(rh==0){env=1;} else{env=0;} 4384 4518 setring RR; 4385 4519 //L0[3]=env; 4520 //"T_env="; env; 4386 4521 return(env); 4387 4522 } 4388 4523 4389 // DimPar4390 // Auxilliary routine to NorSing determining the dimension of a parametric ideal4391 // Does not use @P and define it directly because is assumes that4392 // variables and parameters have been inverted4393 static proc DimPar(ideal E)4394 {4395 def RRH=basering;4396 def RHx=ringlist(RRH);4397 def Prin=ring(RHx[1]);4398 setring(Prin);4399 def E2=std(imap(RRH,E));4400 def d=dim(E2);4401 setring RRH;4402 return(d);4403 }4524 // DimPar 4525 // Auxilliary routine to NorSing determining the dimension of a parametric ideal 4526 // Does not use @P and define it directly because is assumes that 4527 // variables and parameters have been inverted 4528 static proc DimPar(ideal E) 4529 { 4530 def RRH=basering; 4531 def RHx=ringlist(RRH); 4532 def Prin=ring(RHx[1]); 4533 setring(Prin); 4534 def E2=std(imap(RRH,E)); 4535 def d=dim(E2); 4536 setring RRH; 4537 return(d); 4538 } 4404 4539 4405 4540 // locus0(G): Private routine used by locus (the public routine), that 4406 4541 // builds the diferent components. 4407 4542 // input: The output G of the grobcov (in generic representation, which is the default option) 4408 // Options: The algorithm allows the following options a rpair of arguments:4543 // Options: The algorithm allows the following options as pair of arguments: 4409 4544 // "movdim", d : by default movdim is 2 but it can be set to other values 4545 // when locus is called by envelop then as default is uses d=dim @P 4410 4546 // "version", v : There are two versions of the algorithm. ('version',1) is 4411 4547 // a full algorithm that always distinguishes correctly between 'Normal' … … 4418 4554 // Usually locus problems have mover coordinates, variables and tracer coordinates. 4419 4555 // The mover coordinates are to be placed as the last variables, and by default, 4420 // its number is 2. If one consider locus problems in higer dimensions, the number of 4421 // mover coordinates (placed as the last variables) is to be given as an option. 4556 // its number is dim @P, but as option it can be set to other values. 4422 4557 // output: 4423 4558 // list, the canonical Prepresentation of the Normal and NonNormal locus: … … 4436 4571 int t1=1; int t2=1; int i; 4437 4572 def R=basering; 4438 //if(defined(@P)==1){te=1; kill @P; kill @R; kill @RP; } 4439 //setglobalrings(); 4440 // Options 4573 def RL=ringlist(R); 4574 if(defined(@P)==1){te=1; kill @P; kill @R; kill @RP; } 4575 setglobalrings(); 4576 //Options 4441 4577 list DD=#; 4442 int moverdim=nvars(R); 4578 ideal vmov; 4579 int moverdim=size(ringlist(R)[1][2]); 4580 int dimpar=size(RL[1][2]); 4581 int dimvar=size(RL[2]); 4582 int nv=nvars(R); 4583 if(moverdim>nv){moverdim=nv;} 4584 for(i=1;i<=moverdim;i++) 4585 { 4586 //string("T_moverdim=",moverdim,"T_nv=",nv,"T_var(",i,")=",var(i)); 4587 vmov[size(vmov)+1]=var(i+nvmoverdim); 4588 } 4443 4589 int version=0; 4444 int nv=nvars(R);4445 4590 if(nv<4){version=1;} 4446 4591 int comment=0; … … 4449 4594 for(i=1;i<=(size(DD) div 2);i++) 4450 4595 { 4451 if(DD[2*i1]==" movdim"){moverdim=DD[2*i];}4596 if(DD[2*i1]=="vmov"){vmov=DD[2*i];} 4452 4597 if(DD[2*i1]=="version"){version=DD[2*i];} 4453 4598 if(DD[2*i1]=="system"){Fm=DD[2*i];} … … 4456 4601 } 4457 4602 list HHH; 4458 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);} 4459 list G1; list G2; 4603 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) 4604 {return(HHH);} 4605 list G1; list G2; 4460 4606 def G=GG; 4461 4607 list Q1; list Q2; … … 4465 4611 { 4466 4612 attrib(G[i][1],"IsSB",1); 4467 d= locusdim(G[i][2],moverdim);4613 d=dim(std(G[i][1])); 4468 4614 if(d==0){G1[size(G1)+1]=G[i];} 4469 4615 else … … 4484 4630 for(i=1;i<=size(G1RP);i++) 4485 4631 { 4486 4632 kill B; 4487 4633 ideal B; 4488 4634 for(k=1;k<=size(G1RP[i][3]);k++) … … 4496 4642 P1RP[size(P1RP)+1]=list(G1RP[i][3][k][1],G1RP[i][3][k][2],B); 4497 4643 } 4498 } 4644 } //In P1RP the basis has been reduced wrt the top 4499 4645 setring(R); 4500 4646 ideal h; 4647 list P1; 4501 4648 if(t1) 4502 4649 { 4503 defP1=imap(@RP,P1RP);4650 P1=imap(@RP,P1RP); 4504 4651 for(i=1;i<=size(P1);i++) 4505 4652 { … … 4514 4661 } 4515 4662 } 4516 } 4517 else{list P1;} 4663 } //In P1 factors in the basis independent of the parameters have been obtained 4518 4664 ideal BB; int dd; list NS; 4519 4665 for(i=1;i<=size(P1);i++) 4520 4666 { 4521 4522 4523 if(dd==0){P1[i][3]=NS;}//"Special";4524 4667 NS=NorSing(P1[i][3],P1[i][1],P1[i][2][1],DD); 4668 dd=NS[1]; 4669 if(dd==0){P1[i][3]=NS;} //"Special"; 4670 else{P1[i][3]="Normal";} 4525 4671 } 4526 4672 list P2; … … 4548 4694 def C2=imap(R,P2); 4549 4695 def L2=AddLocus(C2); 4550 }4696 } 4551 4697 else{list L2; list C2; kill P2; list P2;} 4552 4553 4554 4555 d=dim(J); // AQUI4556 if(d==0)4557 4558 L2[i][4]=string("Accumulation",L2[i][4]);4559 4560 else{L2[i][4]=string("Degenerate",L2[i][4]);}4561 4698 for(i=1;i<=size(L2);i++) 4699 { 4700 J=std(L2[i][2]); 4701 d=dim(J); 4702 if(d+1==dimpar) 4703 { 4704 L2[i][4]=string("Degenerate",L2[i][4]); 4705 } 4706 else{L2[i][4]=string("Accumulation",L2[i][4]);} 4707 } 4562 4708 list LN; 4563 4709 if(t1==1) … … 4569 4715 for(i=1;i<=size(L2);i++){LN[size(LN)+1]=L2[i];} 4570 4716 } 4717 int tLN=1; 4718 if(size(LN)==0){tLN=0;} 4571 4719 setring(R); 4572 def L=imap(@P,LN); 4573 for(i=1;i<=size(L);i++){if(size(L[i][2])==0){L[i][2]=ideal(1);}} 4574 //if(te==0){kill @R; kill @RP; kill @P;} 4575 list LL; 4576 for(i=1;i<=size(L);i++) 4577 { 4578 if(typeof(L[i][4])=="list") {L[i][4][1]="Special";} 4579 l[1]=L[i][2]; 4580 l[2]=L[i][3]; 4581 l[3]=L[i][4]; 4582 l[4]=L[i][5]; 4583 L[i]=l; 4584 } 4585 return(L); 4720 if(tLN) 4721 { 4722 def L=imap(@P,LN); 4723 for(i=1;i<=size(L);i++){if(size(L[i][2])==0){L[i][2]=ideal(1);}} 4724 list LL; 4725 for(i=1;i<=size(L);i++) 4726 { 4727 if(typeof(L[i][4])=="list") {L[i][4][1]="Special";} 4728 l[1]=L[i][2]; 4729 l[2]=L[i][3]; 4730 l[3]=L[i][4]; 4731 l[4]=L[i][5]; 4732 L[i]=l; 4733 } 4734 } 4735 else{list L;} 4736 return(L); 4586 4737 } 4587 4738 … … 4593 4744 // The Normal locus has two kind of components: Normal and Special. 4594 4745 // The Nonnormal locus has two kind of components: Accumulation and Degenerate. 4595 // Normal components: for each point in the component, 4596 // the number of solutions in the variables is finite, and 4597 // the solutions depend on the point in the component if the component is not 0dimensional. 4598 // Special components: for each point in the component, 4599 // the number of solutions in the variables is finite, 4600 // the component is not 0dimensional, but the solutions do not depend on the 4601 // values of the parameters in the component. 4602 // Accumlation points: are 0dimensional components for which it exist 4603 // an infinite number of solutions. 4604 // Degenerate components: are components of dimension greater than 0 for which 4605 // for every point in the component there exist infinite solutions. 4746 // Normal component: 4747 //  the component has nonzero antiimage 4748 //  each point in the component has 0dimensional antiimage. 4749 // Special component: 4750 //  Nonzero dimensional component whose antiimage is 0dimensional. 4751 // Accumulation points: 4752 //  Zerodimesnional component whose antiimage is nonzerodimensional. 4753 // Degenerate components: 4754 //  Nonzerodimensional component 4755 //  each point in the component has nonzerodimensional antiimage. 4606 4756 // The output components are given as 4607 4757 // ((p1,(p11,..p1s_1),type_1,level_1),..,(pk,(pk1,..pks_k),type_k,level_k) … … 4610 4760 // gives the depth of the component of the constructible set. 4611 4761 proc locus(list GG, list #) 4612 "USAGE: locus(G) 4762 "USAGE: locus(G[,options]); 4763 Calling sequence: 4764 locus(grobcov(S)); 4613 4765 The input must be the grobcov of a parametrical ideal in Q[a][x], 4614 (a=parameters, x=variables). In fact a must be the tracer coordinates 4615 and x the mover coordinates and other auxiliary ones. 4766 (a=parameters, x=variables). In practice a must be the tracer coordinates 4767 and x the mover coordinates and remaining auxiliary variables. 4768 For a parametric locus computation 'a' can contain some extra parameters. 4616 4769 Special routine for determining the locus of points 4617 4770 of geometrical constructions. Given a parametric ideal J 4618 4771 representing the system determining the locus of points (a) 4619 who verify certain properties, the call to locus on the output of grobcov( J)4772 who verify certain properties, the call to locus on the output of grobcov(J) 4620 4773 determines the different classes of locus components, following 4621 4774 the taxonomy defined in … … 4623 4776 \"An Algebraic Taxonomy for Locus Computation in Dynamic Geometry\". 4624 4777 ComputerAided Design 56 (2014) 2233. 4625 The components can be Normal, Special, Accumulation ,Degenerate.4778 The components can be Normal, Special, Accumulation or Degenerate. 4626 4779 OPTIONS: The algorithm allows the following options as pair of arguments: 4627 \"movdim\", d : by default movdim is 2 but it can be set to other values, 4628 and represents the number of mever variables. they should be given as 4629 the last variables of the ring. 4780 \"vmov\", ideal(mover variables) : by default vmov are the last n variables, 4781 where n is the number of parameters of the ring (tracer variables plus extra 4782 parameters). Thus, if the mover coordinates are not indicated, locus 4783 algorithm will assume that they are the last n ring variables. 4784 When locus is called internallt by envelop, by default, the mover variables 4785 are assumed to be all the ring variables. 4630 4786 \"version\", v : There are two versions of the algorithm. (\"version\",1) is 4631 4787 a full algorithm that always distinguishes correctly between 'Normal' 4632 and 'Special' components, whereas \("version\",0) can dec alre a component4633 as 'Normal' being really 'Special', but is more effective. By default (\"version\",1)4634 is used when the number of variables is less than 4 and 0 if not.4788 and 'Special' components, whereas \("version\",0) can declare a component 4789 to be 'Normal' being in fact 'Special', but it is more effective. By default, 4790 (version,1) is used when the number of variables is less than 4 and 0 if not. 4635 4791 The user can force to use one or other version, but it is not recommended. 4636 \"system\", ideal F: if the initial system is passed as an argument. This is actually not used.4637 4792 \"comments\", c: by default it is 0, but it can be set to 1. 4638 4793 Usually locus problems have mover coordinates, variables and tracer coordinates. 4639 The mover coordinates are to be placed as the last variables, and by default, 4640 its number is 2. If one consider locus problems in higer dimensions, the number of 4641 mover coordinates (placed as the last variables) is to be given as an option. 4642 RETURN: The locus. The output is a list of the components ( C_1,.. C_n ) where 4643 C_i = (p_i,(p_i1,..p_is_i), type_i,l evel_i ) and type_i can be 4644 'Normal', 'Special', Accumulation', 'Degenerate'. The 'Special' components 4645 return more information, namely the antiimage of the component, that 4646 is 0dimensional for these kind of components. 4647 Normal components: for each point in the component, 4648 the number of solutions in the variables is finite, and 4649 the solutions depend on the point in the component. 4650 Special components: for each point in the component, 4651 the number of solutions in the variables is finite. The 4652 antiimage of the component is 0dimensional. 4653 Accumlation points: are 0dimensional components whose 4654 antiimage is not zerodimansional. 4655 Degenerate components: are components of dimension greater than 0 4656 whose antiimage is notzerodiemansional. 4657 The components are given in canonical Prepresentation. 4658 The levels of a class of locus are 1, 4659 because they represent locally closed. sets. 4660 NOTE: It can only be called after computing the grobcov of the 4661 parametrical ideal in generic representation ('ext',0), 4662 which is the default. 4663 The basering R, must be of the form Q[a_1,..,a_m][x_1,..,x_n]. 4664 KEYWORDS: geometrical locus, locus, loci. 4794 Example of option call: 4795 > locus(S,\"version\",1,\"vmov\",ideal(x1,y1)) 4796 RETURN: The output is a list of the components (C_1, .. , C_n ) of the locus. 4797 The locus is divided into two class of subsets: the normal and the nonnormal 4798 locus. 4799 The Normal locus has two kind of components: Normal and Special. 4800 The Nonnormal locus has two kind of components: Accumulation and Degenerate. 4801 Each component is given by 4802 Ci=((pi,(pi1,..pis_i),type_i,level_i) 4803 where 4804 the first element is the canonical Prepresentation of the subset. 4805 The type can be : \"Normal\", \"Special\", \"Accumulation\", \"Degenerate\". 4806 Normal component: 4807  the component has nonzero antiimage 4808  each point in the component has 0dimensional antiimage. 4809 Special component: 4810  Nonzero dimensional component whose antiimage is 0dimensional. 4811 The Special components return more information, namely the antiimage of 4812 the component, that is 0dimensional for these kind of components. 4813 Accumulation points: 4814  Zerodimensional component whose antiimage is nonzerodimensional. 4815 Degenerate components: 4816  Nonzerodimensional component 4817  each point in the component has nonzerodimensional antiimage. 4818 The level is the depth of the segment of the constructible locus subset 4819 (normal and nonnormal subsets). 4820 If all levels of a locus are 1, then both subsets are locally closed. 4821 NOTE: The input must be the grobcov of the locus system in generic 4822 representation ('ext',0), which is the default. 4823 KEYWORDS: geometrical locus; locus; loci. 4665 4824 EXAMPLE: locus; shows an example" 4666 4825 { … … 4669 4828 if(defined(@P)==1){tes=1; kill @P; kill @R; kill @RP;} 4670 4829 setglobalrings(); 4671 // 4830 //Options 4672 4831 list DD=#; 4673 int moverdim=nvars(R); 4832 ideal vmov; 4833 int moverdim=size(ringlist(R)[1][2]); 4834 int nv=nvars(R); 4835 if(moverdim>nv){moverdim=nv;} 4836 for(i=1;i<=moverdim;i++){vmov[size(vmov)+1]=var(i+nvmoverdim);} 4674 4837 int version=0; 4675 int nv=nvars(R);4676 4838 if(nv<4){version=1;} 4677 4839 int comment=0; … … 4679 4841 for(i=1;i<=(size(DD) div 2);i++) 4680 4842 { 4681 if(DD[2*i1]==" movdim"){moverdim=DD[2*i];}4843 if(DD[2*i1]=="vmov"){vmov=DD[2*i];} 4682 4844 if(DD[2*i1]=="version"){version=DD[2*i];} 4683 if(DD[2*i1]=="system"){Fm=DD[2*i];}4684 4845 if(DD[2*i1]=="comment"){comment=DD[2*i];} 4685 if(DD[2*i1]=="family"){poly F=DD[2*i];}4686 }4687 int j; int k; 4846 } 4847 DD=list("vmov",vmov,"version",version,"comment",comment); 4848 int j; int k; int te; 4688 4849 def B0=GG[1][2]; 4689 4850 def H0=GG[1][3][1][1]; 4690 if (equalideals(B0,ideal(1)) or (equalideals(H0,ideal(0))==0))4851 if (equalideals(B0,ideal(1)) ) 4691 4852 { 4692 4853 return(locus0(GG,DD)); … … 4695 4856 { 4696 4857 int n=nvars(R); 4697 ideal v mov=var(n1),var(n);4858 ideal vB; 4698 4859 ideal N; 4699 intvec xw; intvec yw; 4700 for(i=1;i<=n1;i++){xw[i]=0;} 4701 xw[n]=1; 4702 for(i=1;i<=n;i++){yw[i]=0;} 4703 yw[n1]=1; 4704 poly px; poly py; 4705 int te=1; 4706 i=1; 4707 while( te and i<=size(B0)) 4708 { 4709 if((deg(B0[i],xw)==1) and (deg(B0[i])==1)){px=B0[i]; te=0;} 4710 i++; 4711 } 4712 i=1; te=1; 4713 while( te and i<=size(B0)) 4714 { 4715 if((deg(B0[i],yw)==1) and (deg(B0[i])==1)){py=B0[i]; te=0;} 4716 i++; 4717 } 4718 N=px,py; 4860 for(i=1;i<=size(B0);i++) 4861 { 4862 if(subset(variables(B0[i]),vmov)){N[size(N)+1]=B0[i];} 4863 } 4719 4864 te=indepparameters(N); 4720 4865 if(te) 4721 4866 { 4722 string("locus detected that the mover must avoid point (",N,") in order to obtain the correct locus");4723 // 4867 string("locus detected that the mover must avoid points (",N,") in order to obtain the correct locus");" "; 4868 //eliminates segments of GG where N is contained in the basis 4724 4869 list nGP; 4725 4870 def GP=GG; … … 4737 4882 } 4738 4883 } 4739 else{"Warning! Problem with more than one mover. Not able to solve it."; list L; return(L);} 4884 else 4885 { 4886 " ";string("Warning! Problem with more than one mover."); 4887 string("Try option 'vmov',ideal(of mover variables) to avoid some point of the mover"); 4888 " ";"Elements of the basis of the generic segment in mover variables="; N;" "; 4889 list L; return(L); 4890 } 4740 4891 } 4741 4892 def LL=locus0(nGP,DD); … … 4751 4902 // Concoid 4752 4903 ideal S96=x^2+y^24,(b2)*xa*y+2*a,(ax)^2+(by)^21; 4753 // System S96=4754 4904 S96; 4905 4755 4906 locus(grobcov(S96)); 4756 kill R;4757 // ********************************************4758 ring R=(0,a,b),(x4,x3,x2,x1),dp;4759 short=0;4760 ideal S=(x13)^2+(x21)^29,4761 (4x2)*(x33)+(x13)*(x41),4762 (3x1)*(x3x1)+(4x2)*(x4x2),4763 (4x4)*a+(x33)*b+3*x44*x3,4764 (ax1)^2+(bx2)^2(x1x3)^2(x2x4)^2;4765 // System S=4766 S;4767 locus(grobcov(S));4768 kill R;4769 //********************************************4770 4771 ring R=(0,x,y),(x1,x2),dp;4772 short=0;4773 ideal S=(x  5)*(x1  1)  (x2  2)*(y  2),4774 (x1  5)^2 + (x2  2)^2  4,4775 2*(x  5)*(x2  2) + 2*(x1  5)*(y  2);4776 // System S=4777 S;4778 locus(grobcov(S));4779 4907 } 4780 4908 4781 4909 // locusdg(G): Special routine for determining the locus of points 4782 // of geometrical constructions. Given a parametric ideal J with 4783 // parameters (a_1,..a_m) and variables (x_1,..,xn), 4784 // representing the system determining 4785 // the locus of points (a_1,..,a_m) who verify certain 4786 // properties, computing the grobcov G of 4787 // J and applying to it locus, determines the different 4788 // classes of locus components. The components can be 4789 // Normal, Special, Accumulation point, Degenerate. 4790 // The output are the components given in Pcanonical form 4791 // of at most 4 constructible sets: Normal, Special, Accumulation, 4792 // Degenerate. 4793 // The description of the algorithm and definitions is 4794 // given in a forthcoming paper by Abanades, Botana, Montes, Recio. 4795 // Usually locus problems have mover coordinates, variables and tracer coordinates. 4796 // The mover coordinates are to be placed as the last variables, and by default, 4797 // its number is 2. If onw consider locus problems in higer dimensions, the number of 4798 // mover coordinates (placed as the last variables) is to be given as an option. 4799 // 4910 // of geometrical constructions in Dynamic Geometry. 4911 // It is to be applied to the output of locus and selects 4912 // as 'Relevant' the 'Normal' and the 'Accumulation' 4913 // components. 4800 4914 // input: The output of locus(G); 4801 4915 // output: 4802 // list, the canonical Prepresentation of the Normal and NonNormal locus: 4803 // The Normal locus has two kind of components: Normal and Special. 4804 // The Nonnormal locus has two kind of components: Accumulation and Degenerate. 4805 // Normal components: for each point in the component, 4806 // the number of solutions in the variables is finite, and 4807 // the solutions depend on the point in the component if the component is not 0dimensional. 4808 // Special components: for each point in the component, 4809 // the number of solutions in the variables is finite, 4810 // the component is not 0dimensional, but the solutions do not depend on the 4811 // values of the parameters in the component. 4812 // Accumlation points: are 0dimensional components for which it exist 4813 // an infinite number of solutions. 4814 // Degenerate components: are components of dimension greater than 0 for which 4815 // for every point in the component there exist infinite solutions. 4916 // list, the canonical Prepresentation of the 'Relevant' components of the locus. 4816 4917 // The output components are given as 4817 4918 // ((p1,(p11,..p1s_1),type_1,level_1),..,(pk,(pk1,..pks_k),type_k,level_k) … … 4820 4921 // gives the depth of the component of the constructible set. 4821 4922 proc locusdg(list L) 4822 "USAGE: locusdg(L) The call must be locusdg(locus(grobcov(S))). 4923 "USAGE: locusdg(L); 4924 Calling sequence: 4925 locusdg(locus(grobcov(S))). 4823 4926 RETURN: The output is the list of the 'Relevant' components of the locus 4824 4927 in Dynamic Geometry:(C1,..,C:m), where … … 4826 4929 The 'Relevant' components are the 'Normal' and 'Accumulation' components 4827 4930 of the locus. (See help for locus). 4828 4829 NOTE: It can only be called after computing the locus. 4830 Calling sequence: locusdg(locus(grobcov(S))); 4831 KEYWORDS: geometrical locus, locus, loci, dynamic geometry 4832 EXAMPLE: locusdg; shows an example" 4931 KEYWORDS: geometrical locus; locus; loci; dynamic geometry 4932 EXAMPLE: example locusdg; shows an example" 4833 4933 { 4834 4934 list LL; … … 4851 4951 // Concoid 4852 4952 ideal S96=x^2+y^24,(b2)*xa*y+2*a,(ax)^2+(by)^21; 4853 // System S96=4854 4953 S96; 4855 4954 locus(grobcov(S96)); 4856 4955 locusdg(locus(grobcov(S96))); 4857 kill R; 4858 //******************************************** 4859 ring R=(0,a,b),(x4,x3,x2,x1),dp; 4860 ideal S=(x13)^2+(x21)^29, 4861 (4x2)*(x33)+(x13)*(x41), 4862 (3x1)*(x3x1)+(4x2)*(x4x2), 4863 (4x4)*a+(x33)*b+3*x44*x3, 4864 (ax1)^2+(bx2)^2(x1x3)^2(x2x4)^2; 4865 short=0; 4866 locus(grobcov(S)); 4867 locusdg(locus(grobcov(S))); 4868 kill R; 4869 //******************************************** 4870 4871 ring R=(0,x,y),(x1,x2),dp; 4872 short=0; 4873 ideal S=(x  5)*(x1  1)  (x2  2)*(y  2), 4874 (x1  5)^2 + (x2  2)^2  4, 4875 2*(x  5)*(x2  2) + 2*(x1  5)*(y  2); 4876 locus(grobcov(S)); 4877 locusdg(locus(grobcov(S))); 4878 } 4879 4880 // locusto: Transforms the output of locus, locusdg, envelop and envelopdg 4956 } 4957 4958 // locusto: Transforms the output of locus, locusdg, envelop 4881 4959 // into a string that can be reed from different computational systems. 4882 4960 // input: 4883 // list L: The output of locus 4961 // list L: The output of locus or locusdg or envelop. 4884 4962 // output: 4885 // string s: The output of locus convertedto a string readable by other programs4963 // string s: Converts the input into a string readable by other programs 4886 4964 proc locusto(list L) 4887 4965 "USAGE: locusto(L); 4888 4966 The argument must be the output of locus or locusdg or 4889 envelop or envelopdg.4967 envelop. 4890 4968 It transforms the output into a string in standard form 4891 readable in many languages(Geogebra).4969 readable in other languages, not only Singular (Geogebra). 4892 4970 RETURN: The locus in string standard form 4893 4971 NOTE: It can only be called after computing either … … 4895 4973  locusdg(locus(grobcov(F))) > locusto( locusdg(locus(grobcov(F))) ) 4896 4974  envelop(F,C) > locusto( envelop(F,C) ) 4897 envelopdg(envelop(F,C)) > locusto( envelopdg(envelop(F,C)) ) 4898 KEYWORDS: geometrical locus, locus, loci. 4899 EXAMPLE: locusto; shows an example" 4975 KEYWORDS: geometrical locus; locus; loci 4976 EXAMPLE: example locusto; shows an example" 4900 4977 { 4901 4978 int i; int j; int k; … … 4961 5038 locusto(locus(grobcov(S))); 4962 5039 locusto(locusdg(locus(grobcov(S)))); 4963 kill R;4964 //********************************************4965 4966 // 1. Take a fixed line l: x1y1=0 and consider4967 // the family F of a lines parallel to l passing through the mover point M4968 // 2. Consider a circle x1^2+x2^225, and a mover point M(x1,x2) on it.4969 // 3. Compute the envelop of the family of lines.4970 4971 ring R=(0,x,y),(x1,y1),lp;4972 poly F=(yy1)(xx1);4973 ideal C=x1^2+y1^225;4974 short=0;4975 4976 // Curves Family F=4977 F;4978 // Conditions C=4979 C;4980 4981 locusto(envelop(F,C));4982 locusto(envelopdg(envelop(F,C)));4983 kill R;4984 //********************************************4985 4986 // Steiner Deltoid4987 // 1. Consider the circle x1^2+y1^21=0, and a mover point M(x1,y1) on it.4988 // 2. Consider the triangle A(0,1), B(1,0), C(1,0).4989 // 3. Consider lines passing through M perpendicular to two sides of ABC triangle.4990 // 4. Obtain the envelop of the lines above.4991 4992 ring R=(0,x,y),(x1,y1,x2,y2),lp;4993 ideal C=(x1)^2+(y1)^21,4994 x2+y21,4995 x2y2x1+y1;4996 matrix M[3][3]=x,y,1,x2,y2,1,x1,0,1;4997 poly F=det(M);4998 4999 short=0;5000 5001 // Curves Family F=5002 F;5003 // Conditions C=5004 C;5005 5006 locusto(envelop(F,C));5007 locusto(envelopdg(envelop(F,C)));5008 }5009 5010 // Auxiliary routine5011 // locusdim5012 // input:5013 // B: ideal, a basis of a segment of the grobcov5014 // dgdim: int, the dimension of the mover (for locus)5015 // by default dgdim is equal to the number of variables5016 static proc locusdim(ideal B, list #)5017 {5018 def R=basering;5019 int dgdim;5020 int nv=nvars(R);5021 if (size(#)>0){dgdim=#[1];}5022 else {dgdim=nv;}5023 int d;5024 list v;5025 ideal vi;5026 int i;5027 for(i=1;i<=dgdim;i++)5028 {5029 v[size(v)+1]=varstr(nvdgdim+i);5030 vi[size(v)+1]=var(nvdgdim+i);5031 }5032 ideal B0;5033 for(i=1;i<=size(B);i++)5034 {5035 if(subset(variables(B[i]),vi))5036 {5037 B0[size(B0)+1]=B[i];5038 }5039 }5040 def RR=ringlist(R);5041 def RR0=RR;5042 RR0[2]=v;5043 def R0=ring(RR0);5044 setring(R0);5045 def B0r=imap(R,B0);5046 B0r=std(B0r);5047 d=dim(B0r);5048 setring R;5049 return(d);5050 }5051 5052 static proc norspec(ideal F)5053 {5054 def RR=basering;5055 def Rx=ringlist(RR);5056 5057 def Rx=ringlist(RR);5058 def @P=ring(Rx[1]);5059 list Lx;5060 Lx[1]=0;5061 Lx[2]=Rx[2]+Rx[1][2];5062 Lx[3]=Rx[1][3];5063 Lx[4]=Rx[1][4];5064 Rx[1]=0;5065 def D=ring(Rx);5066 def @RP=D+@P;5067 exportto(Top,@R); // global ring Q[a][x]5068 exportto(Top,@P); // global ring Q[a]5069 exportto(Top,@RP); // global ring K[x,a] with product order5070 setring(RR);5071 5040 } 5072 5041 5073 5042 // envelop 5074 // Input: n arguments 5075 // poly F: the polynomial defining the family of curves in ring R=0,(x,y),(x1,..,xn),lp; 5076 // ideal C=g1,..,g_{n1}: the set of constraints 5043 // Input: 5044 // poly F: the polynomial defining the family of hypersurfaces in ring R=0,(x_1,..,x_n),(u_1,..,u_m),lp; 5045 // ideal C=g1,..,g_{n1}: the set of constraints; 5046 // options. 5077 5047 // Output: the components of the envolvent; 5078 5048 proc envelop(poly F, ideal C, list #) 5079 5049 "USAGE: envelop(F,C); 5080 The first argument F must be the family of curves of which 5081 on want to compute the envelop. 5082 The second argument C must be the ideal of conditions 5083 over the variables, and should contain as polynomials 5084 as the number of variables 1. 5050 The first argument poly F must be the family of hypersurfaces for which 5051 on want to compute its envelop. 5052 The second argument C must be the ideal of restrictions on 5053 the variables, and should contain less polynomials 5054 than the number of variables 5055 (x_1,..,x_n) are the variables of the hypersurfaces of F, that are considered 5056 as parameters of the parametric ring. 5057 (u_1,..,u_m) are the parameteres of the hypersurfaces, that are considered 5058 as variables of the parametric ring. 5059 Calling sequence: 5060 ring R=(0,x_1,..,x_n),(u_1,..,u_m),lp; 5061 poly F=F(x_1,..,x_n,u_1,..,u_m); 5062 ideal C=g_1(u_1,..u_m),..,g_s(u_1,..u_m); 5063 envelop(F,C,#options); 5064 where s<m. 5065 OPTIONS: The algorithm allows the following options as pair of arguments: 5066 \"vmov\", ideal(mover variables) : by default vmov are u_1,..,u_m. 5067 But it can be restricted by the user to the more convenient ones. 5068 \"version\", v : There are two versions of the algorithm. (\"version\",1) is 5069 a full algorithm that always distinguishes correctly between 'Normal' 5070 and 'Special' components, whereas \("version\",0) can decalre a component 5071 as 'Normal' being really 'Special', but is more effective. By default (\"version\",1) 5072 is used when the number of variables is less than 4 and 0 if not. 5073 The user can force to use one or other version, but it is not recommended. 5074 \"comments\", c: by default it is 0, but it can be set to 1. 5085 5075 RETURN: The components of the envelop with its taxonomy: 5086 The taxonomy distinguishes 'Normal', 5087 'Special', 'Accumulation', 'Degenerate' components. 5088 In the case of 'Special' components, it also 5089 outputs the antiimage of the component 5090 and an integer (01). If the integer is 0 5091 the component is not a curve of the family and is 5092 not considered as 'Relevant' by the envelopdg routine 5093 applied to it, but is considered as 'Relevant' if the integer is 1. 5076 (see locus help). envelop uses locus. 5077 The taxonomy distinguishes \"Normal\", 5078 \"Special\", \"Accumulation\", \"Degenerate\" components. 5079 In the case of \"Special\" components, it also 5080 outputs the antiimage of the component 5094 5081 NOTE: grobcov is called internally. 5095 5082 The basering R, must be of the form Q[a][x] (a=parameters, x=variables). 5096 KEYWORDS: geometrical locus, locus, loci, envelop 5097 EXAMPLE: envelop; shows an example" 5098 { 5099 int tes=0; int i; int j; 5083 KEYWORDS: geometrical locus; locus; loci; envelop 5084 EXAMPLE: example envelop; shows an example" 5085 { 5100 5086 def R=basering; 5087 int tes=0; int i; int j; int k; int m; 5088 int d; 5089 int dp; 5090 ideal BBB; 5101 5091 if(defined(@P)==1){tes=1; kill @P; kill @R; kill @RP;} 5102 5092 setglobalrings(); 5103 // 5093 //Options 5104 5094 list DD=#; 5105 int moverdim=nvars(R); 5095 ideal vmov; 5096 int nv=nvars(R); 5097 for(i=1;i<=nv;i++){vmov[size(vmov)+1]=var(i);} 5098 int numpars=size(ringlist(R)[1][2]); 5106 5099 int version=0; 5107 int nv=nvars(R);5108 5100 if(nv<4){version=1;} 5109 5101 int comment=0; 5102 int familyinfo; 5110 5103 ideal Fm; 5111 5104 for(i=1;i<=(size(DD) div 2);i++) 5112 5105 { 5113 if(DD[2*i1]==" movdim"){moverdim=DD[2*i];}5106 if(DD[2*i1]=="vmov"){vmov=DD[2*i];} 5114 5107 if(DD[2*i1]=="version"){version=DD[2*i];} 5115 if(DD[2*i1]=="system"){Fm=DD[2*i];}5116 5108 if(DD[2*i1]=="comment"){comment=DD[2*i];} 5117 } 5118 int n=nvars(R); 5119 list v; 5120 for(i=1;i<=n;i++){v[size(v)+1]=var(i);} 5121 def MF=jacob(F); 5122 def TMF=transpose(MF); 5123 def Mg=MF; 5124 def TMg=TMF; 5125 for(i=1;i<=n1;i++) 5126 { 5127 Mg=jacob(C[i]); 5128 TMg=transpose(Mg); 5129 TMF=concat(TMF,TMg); 5130 } 5131 poly J=det(TMF); 5132 ideal S=ideal(F)+C+ideal(J); 5133 DD[size(DD)+1]="family"; 5134 DD[size(DD)+1]=F; 5109 if(DD[2*i1]=="familyinfo"){familyinfo=DD[2*i];} 5110 }; 5111 DD=list("vmov",vmov,"version",version,"comment",comment); 5112 int ng=size(C); 5113 ideal S=F; 5114 for(i=1;i<=size(C);i++){S[size(S)+1]=C[i];} 5115 int s=nvng; 5116 if(s>0) 5117 { 5118 matrix M[ng+1][ng+1]; 5119 def cc=comb(nv,ng+1); 5120 poly J; 5121 for(k=1;k<=size(cc);k++) 5122 { 5123 for(j=1;j<=ng+1;j++) 5124 { 5125 M[1,j]=diff(F,var(cc[k][j])); 5126 } 5127 for(i=1;i<=ng;i++) 5128 { 5129 for(j=1;j<=ng+1;j++) 5130 { 5131 M[i+1,j]=diff(C[i],var(cc[k][j])); 5132 } 5133 } 5134 J=det(M); 5135 S[size(S)+1]=J; 5136 } 5137 } 5138 if(comment>0){"System S before grobcov ="; S;} 5135 5139 def G=grobcov(S,DD); 5136 def L=locus(G, DD); 5137 return(L); 5140 list HHH; 5141 if (G[1][1][1]==1 and G[1][2][1]==1 and G[1][3][1][1][1]==0 and G[1][3][1][2][1]==1) 5142 {return(HHH);} 5143 //DD[size(DD)+1]="vmov"; 5144 //DD[size(DD)+1]=4; 5145 def L=locus(G,DD); 5146 list GL; 5147 ideal fam; ideal env; 5148 5149 def Rx=ringlist(R); 5150 def P=ring(Rx[1]); 5151 list Lx; 5152 Lx[1]=0; 5153 Lx[2]=Rx[2]+Rx[1][2]; 5154 Lx[3]=Rx[1][3]; 5155 Lx[4]=Rx[1][4]; 5156 Rx[1]=0; 5157 def D=ring(Rx); 5158 def RP=P+D; 5159 list LL; 5160 list NormalComp; 5161 ideal Gi; 5162 ideal BBBB; 5163 poly B0; 5164 if(familyinfo==1) 5165 { 5166 for(i=1;i<=size(L);i++) 5167 { 5168 if(typeof(L[i][3])=="string") 5169 { 5170 if(L[i][3]=="Normal") 5171 { 5172 NormalComp[size(NormalComp)+1]=L[i][1]; 5173 Gi=S; 5174 Gi[size(Gi)+1]=L[i][1][1]; 5175 //"T_grobcov(Gi)="; grobcov(Gi); 5176 kill HHH; list HHH; 5177 //"T_L[i]="; L[i]; 5178 //d=DimPar(L[i][1]); 5179 if(defined(SL)){kill SL;} 5180 def SL=C; 5181 SL[size(SL)+1]=F; 5182 for(j=1;j<=size(L[i][1]);j++) 5183 { 5184 SL[size(SL)+1]=L[i][1][j]; 5185 } 5186 setring RP; 5187 if(defined(BBBB)){kill BBBB;} 5188 ideal BBBB; 5189 if(defined(BB)){kill BB;} 5190 def BB=imap(R,SL); 5191 if(defined(B0)){kill B0;} 5192 poly B0; 5193 if(defined(LLL)){kill LLL;} 5194 def LLL=imap(R,L); 5195 //"T_BB="; BB; 5196 BB=std(BB); 5197 for(j=1;j<=size(BB);j++) 5198 { 5199 B0=reduce(BB[j],LLL[i][1]); 5200 if(not(B0==0)){BBBB[size(BBBB)+1]=B0;} 5201 } 5202 setring R; 5203 BBB=imap(RP,BBBB); 5204 L[i][5]=BBB; 5205 } 5206 } 5207 } 5208 LL[1]=L; 5209 LL[2]=NormalComp; 5210 list LLL; list LLLL; 5211 int t; 5212 for(k=1;k<=size(LL[2]);k++) 5213 { 5214 for(i=1;i<=size(G);i++) 5215 { 5216 j=1; t=0; 5217 while(t==0 and j<=size(G[i][3])) 5218 { 5219 //"T_LL[2][k]="; LL[2][k]; 5220 //"T_G[i][3][j][1]="; G[i][3][j][1]; 5221 if(equalideals(LL[2][k],G[i][3][j][1])) 5222 { 5223 LLL[size(LLL)+1]=list(k,i,j); 5224 t=1; 5225 } 5226 j++; 5227 } 5228 } 5229 } 5230 LL[3]=LLL; 5231 5232 for(k=1;k<=size(LLL);k++) 5233 { 5234 for(m=k+1;m<=size(LLL);m++) 5235 { 5236 for(i=1;i<=size(G[LLL[k][2]][3][LLL[k][3]][2]);i++) 5237 { 5238 for(j=1;j<=size(G[LLL[m][2]][3][LLL[m][3]][2]);j++) 5239 { 5240 //string("T_(G[",LLL[k][2],"][3][",LLL[k][3],"][2][",i,"])="); G[LLL[k][2]][3][LLL[k][3]][2][i]; 5241 //string("T_(G[",LLL[m][2],"][3][",LLL[m][3],"][2][",j,"])="); G[LLL[m][2]][3][LLL[m][3]][2][j]; 5242 if(equalideals( G[LLL[k][2]][3][LLL[k][3]][2][i] ,G[LLL[m][2]][3][LLL[m][3]][2][j])) 5243 { 5244 //"T_GGG="; LL[2][LLL[k][1]]; LL[2][LLL[m][1]]; G[LLL[k][2]][3][LLL[k][3]][2][i]; 5245 LLLL[size(LLLL)+1]=list(LL[2][LLL[k][1]],LL[2][LLL[m][1]], G[LLL[k][2]][3][LLL[k][3]][2][i]); 5246 } 5247 } 5248 } 5249 } 5250 } 5251 LL[3]=LLLL; 5252 } 5253 else{LL=L;} 5254 return(LL); 5138 5255 } 5139 5256 example … … 5147 5264 5148 5265 ring R=(0,x,y),(x1,y1,x2,y2),lp; 5266 short=0; 5149 5267 ideal C=(x1)^2+(y1)^21, 5150 5268 x2+y21, … … 5152 5270 matrix M[3][3]=x,y,1,x2,y2,1,x1,0,1; 5153 5271 poly F=det(M); 5154 5155 short=0; 5156 5157 // Curves Family F= 5272 // Curves Family F 5158 5273 F; 5159 5274 // Conditions C= 5160 5275 C; 5276 envelop(F,C); 5277 } 5278 5279 proc AssocTanToEnv(poly F,ideal C, poly E, list #) 5280 "USAGE: AssocTanToEnv(F,C,E); 5281 The first argument poly F must be the family of hypersurfaces for which 5282 on want to compute its envelop. 5283 The second argument C must be the ideal of restrictions on 5284 the variables, and should contain s polynomials 5285 being s<n, 5286 (x_1,..,x_n) are the variables of the hypersurfaces of F, that are considered 5287 as parameters of the parametric ring. 5288 (u_1,..,u_m) are the parameteres of the hypersurfaces, that are considered 5289 as variables of the parametric ring. 5290 The third parameter poly E must be the equation of a component of the 5291 envelop of (F,C) previously determined by envelop. 5292 Calling sequence: 5293 ring R=(0,x_1,..,x_n),(u_1,..,u_m),lp; 5294 poly F=F(x_1,..,x_n,u_1,..,u_m); 5295 ideal C=g_1(u_1,..u_m),..,g_s(u_1,..u_m); 5296 poly E(x_1,..,x_n); 5297 AssocTanToEnv(F,C,E,#options); 5298 5299 OPTIONS: The algorithm allows the following options as pair of arguments: 5300 \"vmov\", ideal(mover variables) : by default vmov are u_1,..,u_m. 5301 But it can be restricted by the user to the more convenient ones. 5302 \"version\", v : There are two versions of the algorithm. (\"version\",1) is 5303 a full algorithm that always distinguishes correctly between 'Normal' 5304 and 'Special' components, whereas \("version\",0) can decalre a component 5305 as 'Normal' being really 'Special', but is more effective. By default (\"version\",1) 5306 is used when the number of variables is less than 4 and 0 if not. 5307 The user can force to use one or other version, but it is not recommended. 5308 \"comments\", c: by default it is 0, but it can be set to 1. 5309 RETURN: The interesting segments of the grobcov each one with (lpp,basis,segment). 5310 Fixing the values of (x_1,..,x_n) inside E, the basis allows to detemine the values of the parameters 5311 u_1,..u_m. 5312 NOTE: grobcov is called internally. 5313 The basering R, must be of the form Q[a][x] (a=parameters, x=variables). 5314 KEYWORDS: geometrical locus; locus; loci; envelop, associated tangent 5315 EXAMPLE: example AssocTanToEnv; shows an example" 5316 { 5317 def R=basering; 5318 int tes=0; int i; int j; int k; int m; 5319 int d; 5320 int dp; 5321 poly EE=E; 5322 int moreinfo=1; 5323 ideal BBB; 5324 if(defined(@P)==1){tes=1; kill @P; kill @R; kill @RP;} 5325 setglobalrings(); 5326 //Options 5327 list DD=#; 5328 ideal vmov; 5329 int nv=nvars(R); 5330 for(i=1;i<=nv;i++){vmov[size(vmov)+1]=var(i);} 5331 int numpars=size(ringlist(R)[1][2]); 5332 int version=0; 5333 if(nv<4){version=1;} 5334 int comment=0; 5335 int familyinfo; 5336 ideal Fm; 5337 for(i=1;i<=(size(DD) div 2);i++) 5338 { 5339 if(DD[2*i1]=="vmov"){vmov=DD[2*i];} 5340 if(DD[2*i1]=="version"){version=DD[2*i];} 5341 if(DD[2*i1]=="comment"){comment=DD[2*i];} 5342 if(DD[2*i1]=="familyinfo"){familyinfo=DD[2*i];} 5343 if(DD[2*i1]=="moreinfo"){moreinfo=DD[2*i];} 5344 }; 5345 DD=list("vmov",vmov,"version",version,"comment",comment); 5346 int ng=size(C); 5347 ideal S=F; 5348 for(i=1;i<=size(C);i++){S[size(S)+1]=C[i];} 5349 int s=nvng; 5350 if(s>0) 5351 { 5352 matrix M[ng+1][ng+1]; 5353 def cc=comb(nv,ng+1); 5354 poly J; 5355 for(k=1;k<=size(cc);k++) 5356 { 5357 for(j=1;j<=ng+1;j++) 5358 { 5359 M[1,j]=diff(F,var(cc[k][j])); 5360 } 5361 for(i=1;i<=ng;i++) 5362 { 5363 for(j=1;j<=ng+1;j++) 5364 { 5365 M[i+1,j]=diff(C[i],var(cc[k][j])); 5366 } 5367 } 5368 J=det(M); 5369 S[size(S)+1]=J; 5370 } 5371 } 5372 S[size(S)+1]=EE; 5373 if(comment>0){"System S before grobcov ="; S;} 5374 def G=grobcov(S,DD); 5375 //"T_G=";G; 5376 list GG; 5377 for(i=2;i<=size(G);i++) 5378 { 5379 GG[size(GG)+1]=G[i]; 5380 } 5381 G=GG; 5382 //"T_G=";G; 5383 if(moreinfo>0){return(G);} 5384 else 5385 { 5386 int t=0; 5387 ideal H; 5388 i=1; 5389 while(t==0 and i<=size(G)) 5390 { 5391 //string("T_G[",i,"][3][1][1][1]="); G[i][3][1][1][1]; 5392 //string("T_EE="); EE; 5393 if(G[i][3][1][1][1]==EE) 5394 { 5395 t=1; 5396 H=G[i][2]; 5397 } 5398 i++; 5399 } 5400 return(H); 5401 } 5402 return(G); 5403 } 5404 example 5405 { 5406 "EXAMPLE:"; echo = 2; 5407 ring R=(0,y0,x,y),(t),dp; 5408 short=0; 5409 poly F=(x5*t)^2+y^23^2*t^2; 5410 F; 5411 ideal C; 5412 C; 5161 5413 5162 5414 def Env=envelop(F,C); 5163 5415 Env; 5164 } 5165 5166 // envelopdg 5167 // Input: list L: the output of envelop(poly F, ideal C, list #) 5168 // Output: the relevant components of the envolvent in dynamic geometry; 5169 proc envelopdg(list L) 5170 "USAGE: envelopdg(L); 5171 The input list L must be the output of the call to 5172 the routine 'envolop' of the family of curves 5173 RETURN: The relevant components of the envelop in Dynamic Geometry. 5174 'Normal' and 'Accumulation' components are always considered 5175 'Relevant'. 'Special' components of the envelop outputs 5176 three objects in its characterization: 'Special', the antiimage ideal, 5177 and the integer 0 or 1, that indicates that the given component is 5178 formed (1) or is not formed (0) by curves of the family. Only if yes, 5179 'envelopdg' considers the component as 'Relevant' . 5180 NOTE: It must be called to the output of the 'envelop' routine. 5181 The basering R, must be of the form Q[a,b,..][x,y,..]. 5182 KEYWORDS: geometrical locus, locus, loci, envelop. 5183 EXAMPLE: envelop; shows an example" 5184 { 5185 list LL; 5186 list Li; 5416 // E is a component of the envelop: 5417 poly E=Env[1][1][1]; 5418 E; 5419 def A=AssocTanToEnv(F,C,E); 5420 A; 5421 // The basis of the parameter values of the associated tangent component is 5422 A[1][2][1]; 5423 // Thus t=(5/12)*y0 and the associated tangent component at (x0,y0) is 5424 subst(F,t,(5/12)*y0); 5425 } 5426 5427 proc FamElemsAtEnvCompPoints(poly F,ideal C, poly E, list #) 5428 "USAGE: FamElemsAtEnvCompPoints(F,C,E[,options]); 5429 The first argument poly F must be the family of hypersurfaces for which 5430 on want to compute its envelop. 5431 The second argument C must be the ideal of restrictions on 5432 the variables, and should contain s polynomials 5433 being s<n, 5434 (x_1,..,x_n) are the variables of the hypersurfaces of F, that are considered 5435 as parameters of the parametric ring. 5436 (u_1,..,u_m) are the parameteres of the hypersurfaces, that are considered 5437 as variables of the parametric ring. 5438 The third parameter poly E must be the equation of a component of the 5439 envelop of (F,C) previously determined by envelop. 5440 Calling sequence: 5441 ring R=(0,x_1,..,x_n),(u_1,..,u_m),lp; 5442 poly F=F(x_1,..,x_n,u_1,..,u_m); 5443 ideal C=g_1(u_1,..u_m),..,g_s(u_1,..u_m); 5444 poly E(x_1,..,x_n); 5445 FamElemsAtEnvCompPoints(F,C,E,#options); 5446 5447 OPTIONS: The algorithm allows the following options as pair of arguments: 5448 \"vmov\", ideal(mover variables) : by default vmov are u_1,..,u_m. 5449 But it can be restricted by the user to the more convenient ones. 5450 \"version\", v : There are two versions of the algorithm. (\"version\",1) is 5451 a full algorithm that always distinguishes correctly between \"Normal\" 5452 and \"Special\" components, whereas (\"version\",0) can declare a component 5453 as \"Normal\" being really \"Special\", but is more effective. By default (\"version\",1) 5454 is used when the number of variables is less than 4 and 0 if not. 5455 The user can force to use one or other version, but it is not recommended. 5456 \"comments\", c: by default it is 0, but it can be set to 1. 5457 RETURN: The interesting segments of the grobcov each one with (lpp,basis,segment). 5458 Fixing the values of (x_1,..,x_n) inside E, the basis allows to detemine the values of the 5459 parameters (u_1,..u_m). 5460 NOTE: grobcov is called internally. 5461 The basering R, must be of the form Q[a][x] (a=parameters, x=variables). 5462 KEYWORDS: geometrical locus; locus; loci; envelop; associated tangent 5463 EXAMPLE: example FamElemsAtEnvCompPoints; shows an example" 5464 { 5465 ideal S=C; 5466 S[size(S)+1]=F; 5467 S[size(S)+1]=E; 5468 def G=grobcov(S); 5469 list GG; 5187 5470 int i; 5188 for(i=1;i<=size(L);i++) 5189 { 5190 if(typeof(L[i][3])=="string") 5191 { 5192 if((L[i][3]=="Normal") or (L[i][3]=="Accumulation")){Li=L[i]; Li[3]="Relevant"; LL[size(LL)+1]=Li;} 5193 } 5471 for(i=2; i<=size(G); i++) 5472 { 5473 GG[size(GG)+1]=G[i]; 5474 } 5475 return(GG); 5476 } 5477 example 5478 { 5479 "EXAMPLE:"; echo = 2; 5480 ring R=(0,y0,x,y),(t),dp; 5481 short=0; 5482 poly F=(x5*t)^2+y^23^2*t^2; 5483 F; 5484 ideal C; 5485 C; 5486 5487 def Env=envelop(F,C); 5488 Env; 5489 5490 // E is a component of the envelop: 5491 poly E=Env[1][1][1]; 5492 E; 5493 def A=AssocTanToEnv(F,C,E); 5494 A; 5495 // The basis of the parameter values of the associated tangent component is 5496 A[1][2][1]; 5497 // Thus t=(5/12)*y0 the assocoated tangent family element at (x0,y0) is 5498 subst(F,t,(5/12)*y0); 5499 5500 FamElemsAtEnvCompPoints(F,C,E); 5501 // Thus (12*t^25*y0)^2=0 and the unique circle of the family passing at (x0,y0) in E 5502 // is the associated tangent circle: 5503 subst(F,t,(5/12)*y0); 5504 } 5505 5506 // discrim 5507 proc discrim(poly F0, poly x0) 5508 "USAGE: discrim(f,x); 5509 poly f: the polynomial in Q[a][x] or Q[x] of degree 2 in x 5510 poly x: a variable in the ring. 5511 RETURN: the factorized discriminant of f wrt x for discussing its sign 5512 KEYWORDS: second degree; solving 5513 EXAMPLE: discrim; shows an example" 5514 { 5515 def RR=basering; 5516 int i; 5517 int te; 5518 int d; int dd; 5519 if(size(ringlist(RR)[1])>0) 5520 { 5521 te=1; 5522 setglobalrings(); 5523 setring @RP; 5524 poly F=imap(RR,F0); 5525 poly X=imap(RR,x0); 5526 } 5527 else 5528 {poly F=F0; poly X=x0;} 5529 matrix M=coef(F,X); 5530 d=deg(M[1,1]); 5531 if(d>2){"Degree is higher than 2. No discriminant"; setring RR; return();} 5532 poly dis=(M[2,2])^24*M[2,1]*M[2,3]; 5533 def disp=factorize(dis,0); 5534 if(te==0){return(disp);} 5194 5535 else 5195 5536 { 5196 if(typeof(L[i][3])=="list") 5197 { 5198 if(L[i][3][3]==1) 5199 { 5200 Li=L[i]; Li[3]="Relevant"; LL[size(LL)+1]=Li; 5201 } 5202 } 5203 } 5204 } 5205 return(LL); 5537 setring RR; 5538 def disp0=imap(@RP,disp); 5539 return(disp0); 5540 } 5206 5541 } 5207 5542 example 5208 5543 { 5209 5544 "EXAMPLE:"; echo = 2; 5210 5211 // 1. Take a fixed line l: x1y1=0 and consider 5212 // the family F of a lines parallel to l passing through the mover point M 5213 // 2. Consider a circle x1^2+x2^225, and a mover point M(x1,x2) on it. 5214 // 3. Compute the envelop of the family of lines. 5215 5216 ring R=(0,x,y),(x1,y1),lp; 5217 short=0; 5218 poly F=(yy1)(xx1); 5219 ideal C=x1^2+y1^225; 5220 short=0; 5221 5222 // Curves Family F= 5223 F; 5224 // Conditions C= 5225 C; 5226 5227 envelop(F,C); 5228 envelopdg(envelop(F,C)); 5545 ring R=(0,a,b,c),(x,y),dp; 5546 poly f=a*x^2*y+b*x*y+c*y; 5547 discrim(f,x); 5229 5548 } 5230 5549 5231 5550 // AddLocus: auxilliary routine for locus0 that computes the components of the constructible: 5232 5551 // Input: the list of locally closed sets to be added, each with its type as third argument 5233 // L=[ [LC[11], ,,LC[1k_1], .., [LC[r1],,,LC[rk_r] ] where5552 // L=[ [LC[11],..,LC[1k_1],.., [LC[r1],..,LC[rk_r] ] where 5234 5553 // LC[1]=[p1,[p11,..,p1k],typ] 5235 // Output: the list of components of the constructible union of theL, with the type of the corresponding top5554 // Output: the list of components of the constructible union of L, with the type of the corresponding top 5236 5555 // and the level of the constructible 5237 5556 // L4= [[v1,p1,[p11,..,p1l],typ_1,level]_1 ,.. [vs,ps,[ps1,..,psl],typ_s,level_s] 5238 5557 static proc AddLocus(list L) 5239 5558 { 5240 // int te0=0;5241 // def RR=basering;5242 // if(defined(@P)){te0=1; def Rx=@R; kill @P; setring RR;}5243 5559 list L1; int i; int j; list L2; list L3; 5244 5560 list l1; list l2; … … 5277 5593 { 5278 5594 v=L2[k0][2]; 5595 l4[1]=v; l4[2]=p1; l4[3]=L3[i][2][j][2]; l4[5]=level; 5596 if(size(L2[k0])>2){l4[4]=L2[k0][3];} 5597 L4[size(L4)+1]=l4; 5279 5598 } 5280 5599 else{"ERROR p1 NOT FOUND";} 5281 l4[1]=v; l4[2]=p1; l4[3]=L3[i][2][j][2]; l4[5]=level;5282 if(size(L2[k0])>2){l4[4]=L2[k0][3];}5283 L4[size(L4)+1]=l4;5284 5600 } 5285 5601 } … … 5305 5621 for(i=1;i<=size(L);i++) 5306 5622 { 5307 Sc=PtoCrep (list(L[i]));5623 Sc=PtoCrep0(list(L[i])); 5308 5624 Lc[size(Lc)+1]=Sc; 5309 5625 } 5310 list S=ConsLevels(Lc)[1]; 5626 list S=ConsLevels(Lc); 5627 S=ConsLevelsToLevels(S); 5311 5628 list Sout; 5312 5629 list Lev; 5313 5630 for(i=1;i<=size(S);i++) 5314 5631 { 5315 Lev=list( i,Prep(S[i][2][1],S[i][2][2]));5632 Lev=list(S[i][1],Prep(S[i][2][1],S[i][2][2])); 5316 5633 Sout[size(Sout)+1]=Lev; 5317 5634 } … … 5321 5638 //********************* End locus **************************** 5322 5639 5640 //********************* Begin WLemma ********************** 5641 5642 // input ideal F in @R 5643 // ideal a in @R but only depending on parameters 5644 // F is a generating ideal in V(a); 5645 // output: ideal b in @R but depending only on parameters 5646 // ideal G=GBasis(F) in V(a) \ V(b) 5647 proc WLemma(ideal F,ideal a, list #) 5648 "USAGE: WLemma(F,A,#); 5649 The first argument ideal F in K[a][x]; 5650 The second argument ideal A in K[a] 5651 ring R=(0,x_1,..,x_n),(u_1,..,u_m),lp; 5652 ideal F=f1(x_1,..,x_n,u_1,..,u_m),..,fs(x_1,..,x_n,u_1,..,u_m); 5653 ideal A=g_1(u_1,..u_m),..,g_s(u_1,..u_m); 5654 list # : Options 5655 Calling sequence: 5656 WLemma(F,A,#); 5657 OPTIONS: either (\"rep\", 0) or (\"rep\",1) the representation of 5658 the resulting segment, by default is 5659 0 =Prepresentation, (default) but can be set to 5660 1=Crepresentation. 5661 RETURN: The list of (lpp,B,S) = (leading power product, basis, segment) 5662 being B the reduced Groebner Basis given by Iregular functions 5663 of the specialized ideal F on the 5664 segment S given in P or Crepresentation 5665 NOTE: The basering R, must be of the form Q[a][x] (a=parameters, x=variables). 5666 KEYWORDS: Wibmer's Lemma 5667 EXAMPLE: WLemma; shows an example" 5668 { 5669 list L=#; 5670 int i; int j; 5671 def RR=basering; 5672 setglobalrings(); 5673 setring(@RP); 5674 ideal FF=imap(RR,F); 5675 FF=std(FF); 5676 ideal AA=imap(RR,a); 5677 AA=std(AA); 5678 ideal FFa; 5679 poly r; 5680 for(i=1; i<=size(FF);i++) 5681 { 5682 r=reduce(FF[i],AA); 5683 if(r!=0){FFa[size(FFa)+1]=r;} 5684 } 5685 setring RR; 5686 ideal Fa=imap(@RP,FFa); 5687 ideal AAA=imap(@RP,AA); 5688 ideal lppFa; 5689 ideal lcFa; 5690 for(i=1;i<=size(Fa);i++) 5691 { 5692 lppFa[size(lppFa)+1]=leadmonom(Fa[i]); 5693 lcFa[size(lcFa)+1]=leadcoef(Fa[i]); 5694 } 5695 setring @RP; 5696 ideal lccr=imap(RR,lppFa); 5697 lccr=std(lccr); 5698 setring RR; 5699 ideal lcc=imap(@RP,lccr); 5700 list J; list Jx; 5701 ideal Jci; 5702 ideal Jxi; 5703 list B; 5704 for(i=1;i<=size(lcc);i++) 5705 { 5706 kill Jci; ideal Jci; kill Jxi; ideal Jxi; 5707 for(j=1;j<=size(Fa);j++) 5708 { 5709 if(lppFa[j]==lcc[i]) 5710 { 5711 Jci[size(Jci)+1]=lcFa[j]; 5712 Jxi[size(Jxi)+1]=Fa[j]; 5713 } 5714 } 5715 J[size(J)+1]=Jci; 5716 B[size(B)+1]=Jxi; 5717 } 5718 setring @P; 5719 list Jp=imap(RR,J); 5720 ideal JL=product(Jp); 5721 setring(RR); 5722 def JLA=imap(@P,JL); 5723 list PR; 5724 if (size(L)>0) 5725 { 5726 if((L[1]=="rep") and (L[2]==1)) 5727 { 5728 PR=Crep(AAA, JLA); 5729 } 5730 else 5731 {PR=Prep(AAA, JLA);} 5732 } 5733 else{PR=Prep(AAA, JLA);} 5734 // setring(RR); 5735 // list PRR=imap(@P,PR); 5736 return(list(lcc,B,PR)); 5737 } 5738 example 5739 { 5740 "EXAMPLE:"; echo = 2; 5741 if(defined(R)){kill R;} 5742 ring R=(0,a,b,c,d,e,f),(x,y),lp; 5743 ideal F=a*x^2+b*x*y+c*y^2,d*x^2+e*x*y+f*y^2; 5744 ideal A=a*eb*d; 5745 WLemma(F,A); 5746 WLemma(F,A,"rep",1); 5747 } 5748 5749 //********************* End WLemma ************************
Note: See TracChangeset
for help on using the changeset viewer.