Changeset 27e935 in git for IntegerProgramming/toric.lib
- Timestamp:
- May 11, 2000, 11:04:05 AM (24 years ago)
- Branches:
- (u'spielwiese', '5b153614cbc72bfa198d75b1e9e33dab2645d9fe')
- Children:
- 54b3568dedfbe0f7886ed616396038a36933c3e3
- Parents:
- 0658d94b02b071f2a019543baea28fe2c95d9226
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
IntegerProgramming/toric.lib
r0658d9 r27e935 1 // $Id: toric.lib,v 1. 1.1.1 2000-05-02 12:58:31 Singular Exp $1 // $Id: toric.lib,v 1.2 2000-05-11 09:03:41 Singular Exp $ 2 2 // 3 3 // author : Christine Theis 4 4 // 5 6 //version="$Id: toric.lib,v 1.1.1.1 2000-05-02 12:58:31 Singular Exp $"; 5 //version="$Id: toric.lib,v 1.2 2000-05-11 09:03:41 Singular Exp $"; 7 6 8 7 /////////////////////////////////////////////////////////////////////////////// 9 8 10 9 info=" 11 LIBRARY: toric.lib PROCEDURES FOR COMPUTING TORIC IDEALS12 13 14 Let A an integral (mxn)-matrix. The toric ideal I_A of A is defined 10 LIBRARY: toric.lib PROCEDURES FOR COMPUTING TORIC IDEALS 11 12 13 Let A an integral (mxn)-matrix. The toric ideal I_A of A is defined 15 14 as the ideal 16 15 … … 18 17 19 18 in the ring of n variables x:=x1,...,xn. 20 Toric ideals play an important role in polyhedral geometry and may also be 21 used for integer programming. They are generated by binomials with 22 relatively prime monomials. Buchberger's algorithm can be specialized to 23 these structures in a way that considerably speeds up computation. 19 Toric ideals play an important role in polyhedral geometry and may also be 20 used for integer programming. They are generated by binomials with 21 relatively prime monomials. Buchberger's algorithm can be specialized to 22 these structures in a way that considerably speeds up computation. 24 23 25 24 26 25 toric_ideal(intmat A, string alg); 27 26 toric_ideal(intmat A, string alg, intvec prsv); 28 27 29 28 procedures for computing the toric ideal of A 30 They return the standard basis of the toric ideal of A with respect 29 They return the standard basis of the toric ideal of A with respect 31 30 to the term ordering in the actual basering. 32 When calling this procedure, a ring with n variables should be active. 31 When calling this procedure, a ring with n variables should be active. 33 32 Not all term orderings are supported: The usual global term orderings 34 33 may be used, but no block orderings combining them. … … 42 41 43 42 The last two seem to be the fastest in the actual implementation. 44 `alg' should be the abbreviation (in brackets) for an algorithm 45 as above. 43 `alg' should be the abbreviation (in brackets) for an algorithm 44 as above. 46 45 If `alg' is chosen to be `blr' or `hs', the algorithm needs a 47 46 vector with positive coefficcients in the row space of A. If 48 47 no row of A contains only positive entries, one must use the 49 second version of toric_ideal which takes such a vector in the 50 third argument. 51 48 second version of toric_ideal which takes such a vector in the 49 third argument. 50 52 51 53 52 toric_std(ideal I); 54 53 55 computes the standard basis of I using the specialized Buchberger 56 algorithm. The generating system by which I is given has to consist 54 computes the standard basis of I using the specialized Buchberger 55 algorithm. The generating system by which I is given has to consist 57 56 of binomials of the form x^u-x^v (although there are other generating 58 57 systems of toric ideals). There is no real check if I is toric. 59 If the generator list of I contains a binomial whose monomials are not 58 If the generator list of I contains a binomial whose monomials are not 60 59 relatively prime, the procedure outputs a warning. If I is generated by 61 binomials of the above form, but not toric, toric_std computes an ideal 60 binomials of the above form, but not toric, toric_std computes an ideal 62 61 `between' I and its saturation with respect to all variables. 63 62 64 63 "; 65 64 66 65 /////////////////////////////////////////////////////////////////////////////// 67 66 68 69 70 67 static proc toric_ideal_1(intmat A, string alg) 71 { 68 { 72 69 ideal I; 73 70 // to be returned … … 78 75 "ERROR: number of matrix columns must be greater or equal number of ring variables"; 79 76 return(I); 80 } 77 } 81 78 82 79 // check suitability of actual term ordering … … 167 164 "Warning: block orderings are not supported; Dp used for computation"; 168 165 } 169 } 166 } 170 167 171 168 int pos; … … 263 260 return(I); 264 261 } 265 262 266 263 // check algorithm 267 264 if(alg=="ct" || alg=="pct") … … 274 271 } 275 272 276 // create temporary file with that the external program is called 273 // create temporary file with that the external program is called 277 274 278 275 int dummy; 279 276 int process=system("pid"); 280 string matrixfile="temp_MATRIX"+string(process); 277 string matrixfile="temp_MATRIX"+string(process); 281 278 link MATRIX=":w "+matrixfile; 282 279 open(MATRIX); … … 295 292 } 296 293 } 297 294 298 295 // search for positive row space vector, if required by the 299 // algorithm 296 // algorithm 300 297 int found=0; 301 298 if((alg=="blr") || (alg=="hs")) 302 { 299 { 303 300 for(i=1;i<=nrows(A);i++) 304 301 { … … 318 315 if(found==0) 319 316 { 320 "ERROR: algorithm needs positive vector in the row space of the matrix"; 317 "ERROR: algorithm needs positive vector in the row space of the matrix"; 321 318 close(MATRIX); 322 319 dummy=system("sh","rm -f "+matrixfile); 323 320 return(I); 324 } 321 } 325 322 write(MATRIX,"positive row space vector:"); 326 323 for(j=1;j<=ncols(A);j++) 327 { 324 { 328 325 write(MATRIX,A[found,j]); 329 326 } … … 343 340 pos=find(toric_ideal,":",pos); 344 341 pos++; 345 342 346 343 while(toric_ideal[pos]==" " || toric_ideal[pos]==newline) 347 344 { … … 355 352 } 356 353 execute("generators="+number_string+";"); 357 354 358 355 intvec v; 359 356 poly head; … … 388 385 if(v[j]>0) 389 386 { 390 head=head*var(j)^v[j]; 387 head=head*var(j)^v[j]; 391 388 } 392 389 } 393 390 I[i]=head-tail; 394 391 } 395 392 396 393 // delete all created files 397 394 dummy=system("sh","rm -f "+matrixfile); … … 401 398 } 402 399 403 404 405 400 static proc toric_ideal_2(intmat A, string alg, intvec prsv) 406 { 401 { 407 402 ideal I; 408 403 // to be returned … … 413 408 "ERROR: number of matrix columns must equal size of positive row space vector"; 414 409 return(I); 415 } 410 } 416 411 417 412 // check suitability of actual basering … … 420 415 "ERROR: number of matrix columns must be greater or equal number of ring variables"; 421 416 return(I); 422 } 417 } 423 418 424 419 // check suitability of actual term ordering … … 509 504 "Warning: block orderings are not supported; Dp used for computation"; 510 505 } 511 } 506 } 512 507 513 508 int pos; … … 605 600 return(I); 606 601 } 607 602 608 603 // check algorithm 609 604 if(alg=="ct" || alg=="pct") … … 616 611 } 617 612 618 // create temporary file with that the external program is called 613 // create temporary file with that the external program is called 619 614 620 615 int dummy; 621 616 int process=system("pid"); 622 string matrixfile="temp_MATRIX"+string(process); 617 string matrixfile="temp_MATRIX"+string(process); 623 618 link MATRIX=":w "+matrixfile; 624 619 open(MATRIX); … … 640 635 // enter positive row space vector, if required by the algorithm 641 636 if((alg=="blr") || (alg=="hs")) 642 { 637 { 643 638 write(MATRIX,"positive row space vector:"); 644 639 for(j=1;j<=ncols(A);j++) 645 { 640 { 646 641 write(MATRIX,prsv[j]); 647 642 } … … 660 655 pos=find(toric_ideal,":",pos); 661 656 pos++; 662 657 663 658 while(toric_ideal[pos]==" " || toric_ideal[pos]==newline) 664 659 { … … 672 667 } 673 668 execute("generators="+number_string+";"); 674 669 675 670 intvec v; 676 671 poly head; … … 705 700 if(v[j]>0) 706 701 { 707 head=head*var(j)^v[j]; 702 head=head*var(j)^v[j]; 708 703 } 709 704 } 710 705 I[i]=head-tail; 711 706 } 712 707 713 708 // delete all created files 714 709 dummy=system("sh","rm -f "+matrixfile); … … 718 713 } 719 714 720 721 722 715 proc toric_ideal 723 "USAGE: 716 "USAGE: 724 717 toric_ideal(A,alg); A intmat, alg string 725 toric_ideal(A,alg,prsv); A intmat, alg string, prsv intvec 718 toric_ideal(A,alg,prsv); A intmat, alg string, prsv intvec 726 719 RETURN: toric ideal of A as explained in toric.lib 727 720 return type = ideal … … 737 730 } 738 731 } 739 740 741 742 732 example 743 733 { 744 "EXAMPLE"; echo=2;745 746 ring r=0,(x,y,z),dp; 747 748 // call with two arguments749 intmat A[2][3]=1,1,0,0,1,1;750 A;751 752 ideal I=toric_ideal(A,"du");753 I;754 755 I=toric_ideal(A,"blr");756 I;757 758 // call with three arguments759 intvec prsv=1,2,1;760 I=toric_ideal(A,"blr",prsv);761 I;762 734 "EXAMPLE"; echo=2; 735 736 ring r=0,(x,y,z),dp; 737 738 // call with two arguments 739 intmat A[2][3]=1,1,0,0,1,1; 740 A; 741 742 ideal I=toric_ideal(A,"du"); 743 I; 744 745 I=toric_ideal(A,"blr"); 746 I; 747 748 // call with three arguments 749 intvec prsv=1,2,1; 750 I=toric_ideal(A,"blr",prsv); 751 I; 752 763 753 } 764 754 765 766 767 755 proc toric_std(ideal I) 768 "USAGE: toric_std(I); I ideal 756 "USAGE: toric_std(I); I ideal 769 757 RETURN: standard basis of I as explained in toric.lib 770 758 return type = ideal … … 861 849 "Warning: block orderings are not supported; Dp used for computation"; 862 850 } 863 } 851 } 864 852 865 853 int pos; … … 887 875 } 888 876 } 889 877 890 878 if(external_ord=="" && find(singular_ord,"wp")==3) 891 879 { … … 958 946 return(I); 959 947 } 960 948 961 949 // create first temporary file with which the external program is called 962 950 963 951 int dummy; 964 952 int process=system("pid"); 965 string groebnerfile="temp_GROEBNER"+string(process); 953 string groebnerfile="temp_GROEBNER"+string(process); 966 954 link GROEBNER=":w "+groebnerfile; 967 955 open(GROEBNER); … … 969 957 write(GROEBNER,"GROEBNER","computed with algorithm:","pt","term ordering:","elimination block",0,"weighted block",nvars(basering),external_ord); 970 958 // algorithm is totally unimportant, only required by the external program 971 959 972 960 for(i=1;i<=nvars(basering);i++) 973 961 { 974 962 write(GROEBNER,weightvec[i]); 975 963 } 976 964 977 965 write(GROEBNER,"size:",size(I),"Groebner basis:"); 978 966 poly head; … … 983 971 for(j=1;j<=size(I);j++) 984 972 { 985 // test suitability of generator j 973 // test suitability of generator j 986 974 rest=I[j]; 987 975 head=lead(rest); … … 989 977 tail=lead(rest); 990 978 rest=rest-tail; 991 979 992 980 if(head==0 && tail==0 && rest!=0) 993 981 { … … 997 985 return(J); 998 986 } 999 987 1000 988 if(leadcoef(tail)!=-leadcoef(head)) 1001 989 // generator no difference of monomials (or a constant multiple) … … 1006 994 return(J); 1007 995 } 1008 996 1009 997 if(gcd(head,tail)!=1) 1010 998 { 1011 999 "Warning: monomials of generator "+string(j)+" of input ideal are not relatively prime"; 1012 1000 } 1013 1001 1014 1002 // write vector representation of generator j into the file 1015 1003 v=leadexp(head)-leadexp(tail); … … 1020 1008 } 1021 1009 close(GROEBNER); 1022 1010 1023 1011 // create second temporary file 1024 1012 1025 string newcostfile="temp_NEW_COST"+string(process); 1013 string newcostfile="temp_NEW_COST"+string(process); 1026 1014 link NEW_COST=":w "+newcostfile; 1027 1015 open(NEW_COST); 1028 1016 1029 write(NEW_COST,"NEW_COST","variables:",nvars(basering),"cost vector:"); 1017 write(NEW_COST,"NEW_COST","variables:",nvars(basering),"cost vector:"); 1030 1018 for(i=1;i<=nvars(basering);i++) 1031 1019 { 1032 1020 write(NEW_COST,weightvec[i]); 1033 1021 } 1034 1022 1035 1023 // call external program 1036 1024 dummy=system("sh","change_cost "+groebnerfile+" "+newcostfile); 1037 1025 1038 1026 // read toric standard basis from created file 1039 1027 link TORIC_IDEAL=":r "+newcostfile+".GB.pt"; … … 1044 1032 pos=find(toric_ideal,":",pos); 1045 1033 pos++; 1046 1034 1047 1035 while(toric_ideal[pos]==" " || toric_ideal[pos]==newline) 1048 1036 { … … 1085 1073 if(v[i]>0) 1086 1074 { 1087 head=head*var(i)^v[i]; 1075 head=head*var(i)^v[i]; 1088 1076 } 1089 1077 } 1090 1078 J[j]=head-tail; 1091 1079 } 1092 1080 1093 1081 // delete all created files 1094 1082 dummy=system("sh","rm -f "+groebnerfile); 1095 1083 dummy=system("sh","rm -f "+groebnerfile+".GB.pt"); 1096 dummy=system("sh","rm -f "+newcostfile); 1084 dummy=system("sh","rm -f "+newcostfile); 1097 1085 1098 1086 return(J); 1099 1087 } 1100 1101 1102 1103 1088 example 1104 1089 { 1105 "EXAMPLE"; echo=2;1106 1107 ring r=0,(x,y,z),wp(3,2,1);1108 1109 // call with toric ideal (of the matrix A=(1,1,1))1110 ideal I=x-y,x-z;1111 ideal J=toric_std(I);1112 J;1113 1114 // call with the same ideal, but badly chosen generators:1115 // not only binomials 1116 I=x-y,2x-y-z;1117 J=toric_std(I);1118 // binomials whose monomials are not relatively prime1119 I=x-y,xy-yz,y-z;1120 J=toric_std(I);1121 J;1122 1123 // call with a non-toric ideal that seems to be toric1124 I=x-yz,xy-z;1125 J=toric_std(I);1126 J;1127 // comparison with real standard basis and saturation 1128 ideal H=std(I);1129 H;1130 LIB "elim.lib";1131 sat(H,xyz);1090 "EXAMPLE"; echo=2; 1091 1092 ring r=0,(x,y,z),wp(3,2,1); 1093 1094 // call with toric ideal (of the matrix A=(1,1,1)) 1095 ideal I=x-y,x-z; 1096 ideal J=toric_std(I); 1097 J; 1098 1099 // call with the same ideal, but badly chosen generators: 1100 // not only binomials 1101 I=x-y,2x-y-z; 1102 J=toric_std(I); 1103 // binomials whose monomials are not relatively prime 1104 I=x-y,xy-yz,y-z; 1105 J=toric_std(I); 1106 J; 1107 1108 // call with a non-toric ideal that seems to be toric 1109 I=x-yz,xy-z; 1110 J=toric_std(I); 1111 J; 1112 // comparison with real standard basis and saturation 1113 ideal H=std(I); 1114 H; 1115 LIB "elim.lib"; 1116 sat(H,xyz); 1132 1117 } 1133 1118 1134 1135 1136 1119 ////////////////////////////////////////////////////////////////////////////// 1137 1138 1139 // toric_std(ideal I); ring term ordering
Note: See TracChangeset
for help on using the changeset viewer.