Changeset 2653d3 in git for Singular/LIB/bfct.lib
- Timestamp:
- Jan 15, 2009, 8:19:09 PM (15 years ago)
- Branches:
- (u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
- Children:
- 3ab3910724146bbcc35aa79fd9d30888ef6226e6
- Parents:
- a0350e9a597a2d719cdc7e27752cc38d0a6147b8
- File:
-
- 1 edited
-
Singular/LIB/bfct.lib (modified) (40 diffs)
Legend:
- Unmodified
- Added
- Removed
-
Singular/LIB/bfct.lib
ra0350e9 r2653d3 1 1 ////////////////////////////////////////////////////////////////////////////// 2 version="$Id: bfct.lib,v 1.1 1 2009-01-14 16:07:03 SingularExp $";2 version="$Id: bfct.lib,v 1.12 2009-01-15 19:19:09 levandov Exp $"; 3 3 category="Noncommutative"; 4 4 info=" … … 16 16 @* - the list of its roots, which are known to be rational 17 17 @* - the multiplicities of the roots. 18 @* References: Saito, St rurmfels, Takayama: Groebner Deformations of Hypergeometric18 @* References: Saito, Sturmfels, Takayama: Groebner Deformations of Hypergeometric 19 19 @* Differential Equations (2000), Noro: An Efficient Modular Algorithm for Computing 20 20 @* the Global b-function, (2002). … … 29 29 bfctIdeal(I,w[,s,t]); computes the global b-function of ideal I w.r.t. the weight w 30 30 pIntersect(f,I); intersection of the ideal I with the subalgebra K[f] for a poly f 31 pIntersectSyz(f,I[,p,s,t]); 31 pIntersectSyz(f,I[,p,s,t]); intersection of the ideal I with the subalgebra K[f] for a poly f 32 32 linReduce(f,I[,s]); reduces a poly by linear reductions only 33 linReduceIdeal(I,[s,t]); reduces generators of ideal by linear reductions only 33 34 linSyzSolve(I[,s]); computes a linear dependency of the elements of ideal I 34 35 35 36 AUXILIARY PROCEDURES: 36 37 37 allPositive(v); checks whether all entries of an intvec are positive 38 allPositive(v); checks whether all entries of an intvec are positive 38 39 scalarProd(v,w); computes the standard scalar product of two intvecs 39 40 vec2poly(v[,i]); constructs an univariate poly with given coefficients … … 79 80 EXAMPLE: example gradedWeyl; shows examples 80 81 NOTE: u[i] is the weight of x(i), v[i] the weight of D(i). 81 @* u+v has to be a non-negative intvec. 82 @* u+v has to be a non-negative intvec. 82 83 " 83 84 { … … 88 89 { 89 90 ERROR("weight vectors have wrong dimension"); 90 } 91 } 91 92 intvec uv,gr; 92 93 uv = u+v; … … 97 98 if (uv[i]==0) 98 99 { 99 100 gr[i] = 0; 100 101 } 101 102 else 102 103 { 103 104 gr[i] = 1; 104 105 } 105 106 } … … 114 115 for (i=1; i<=n; i++) 115 116 { 116 if (gr[i] == 1) 117 if (gr[i] == 1) 117 118 { 118 119 l2[n+i] = "xi("+string(i)+")"; … … 229 230 proc linReduceIdeal(ideal I, list #) 230 231 "USAGE: linReduceIdeal(I [,s,t]); I an ideal, s,t optional ints 231 RETURN: ideal/list, linear reductum (over field) of f by elements from I232 PURPOSE: reduce a polyonly by linear reductions (no monomial multiplications)233 NOTE: If s<>0, a list consisting of the reduced ideal and the coefficient 232 RETURN: ideal/list, linear reductum (over field) of I by its elements 233 PURPOSE: reduce a list of polys only by linear reductions (no monomial multiplications) 234 NOTE: If s<>0, a list consisting of the reduced ideal and the coefficient 234 235 @* vectors of the used reductions given as module is returned. 235 236 @* Otherwise (and by default), only the reduced ideal is returned. … … 255 256 if (typeof(#[2])=="int" || typeof(#[2])=="number") 256 257 { 257 258 redtail = #[2]; 258 259 } 259 260 } … … 273 274 for (i=1; i<=sI; i++) 274 275 { 275 276 277 278 279 280 281 282 283 284 285 276 if (I[i] == 0) 277 { 278 j++; 279 J[j] = 0; 280 ordJ[j] = -1; 281 M[j] = gen(i); 282 } 283 else 284 { 285 M[i+sZeros-j] = gen(lJ[2][i-j]+j); 286 } 286 287 } 287 288 } … … 290 291 for (i=1; i<=sZeros; i++) 291 292 { 292 293 293 J[i] = 0; 294 ordJ[i] = -1; 294 295 } 295 296 } … … 327 328 for (j=sZeros+1; j<i; j++) 328 329 { 329 330 331 332 { 333 334 335 336 337 338 339 340 341 342 343 344 345 346 330 if (lm == 0) { break; } 331 if (ordlm > maxordJ) { break; } 332 if (ordlm == ordJ[j]) 333 { 334 if (lm > maxlmJ) { break; } 335 if (lm == lmJ[j]) 336 { 337 dbprint(ppl,"reducing " + string(redpoly)); 338 dbprint(ppl," with " + string(J[j])); 339 c = leadcoef(redpoly)/leadcoef(J[j]); 340 redpoly = redpoly - c*J[j]; 341 dbprint(ppl," to " + string(redpoly)); 342 lm = leadmonom(redpoly); 343 ordlm = ord(lm); 344 if (remembercoeffs <> 0) { M[i] = M[i] - c * M[j]; } 345 reduction = 1; 346 } 347 } 347 348 } 348 349 } … … 354 355 lmJ = insertGenerator(lmJ,lm,j); 355 356 ordJ = insertGenerator(ordJ,poly(ordlm),j); 356 if (remembercoeffs <> 0) 357 if (remembercoeffs <> 0) 357 358 { 358 359 v = M[i]; … … 373 374 for (j=i+1; j<=sI; j++) 374 375 { 375 376 377 378 { 379 380 381 382 383 384 385 386 387 388 389 376 for (k=2; k<=size(J[j]); k++) // run over all terms in J[j] 377 { 378 if (ordJ[i] == ord(J[j][k])) 379 { 380 if (lm == normalize(J[j][k])) 381 { 382 c = leadcoef(J[j][k])/leadcoef(J[i]); 383 dbprint(ppl,"reducing " + string(J[j])); 384 dbprint(ppl," with " + string(J[i])); 385 J[j] = J[j] - c*J[i]; 386 dbprint(ppl," to " + string(J[j])); 387 if (remembercoeffs <> 0) { M[j] = M[j] - c * M[i]; } 388 } 389 } 390 } 390 391 } 391 392 } … … 410 411 RETURN: poly/list, linear reductum (over field) of f by elements from I 411 412 PURPOSE: reduce a poly only by linear reductions (no monomial multiplications) 412 NOTE: If s<>0, a list consisting of the reduced poly and the coefficient 413 NOTE: If s<>0, a list consisting of the reduced poly and the coefficient 413 414 @* vector of the used reductions is returned, otherwise (and by default) 414 415 @* only reduced poly is returned. … … 433 434 if (typeof(#[2])=="int" || typeof(#[2])=="number") 434 435 { 435 436 redtail = #[2]; 436 437 } 437 438 if (size(#)>2) 438 439 { 439 440 if (typeof(#[3])=="int" || typeof(#[3])=="number") 440 441 { 441 442 442 prepareideal = #[3]; 443 } 443 444 } 444 445 } … … 469 470 for (i=1; i<=sI; i++) 470 471 { 471 472 M[i] = gen(i); 472 473 } 473 474 } … … 495 496 if (ordf == ordI[i]) 496 497 { 497 498 499 500 501 502 503 504 505 506 507 498 if (lm == lmI[i]) 499 { 500 c = leadcoef(f)/lcI[i]; 501 f = f - c*I[i]; 502 lm = leadmonom(f); 503 ordf = ord(lm); 504 if (remembercoeffs <> 0) 505 { 506 v = v - c * M[i]; 507 } 508 } 508 509 } 509 510 } … … 518 519 for (j=1; j<=size(f); j++) 519 520 { 520 521 522 523 524 525 526 527 528 529 530 531 532 521 if (ord(f[j]) == ordI[i]) 522 { 523 if (normalize(f[j]) == lmI[i]) 524 { 525 c = leadcoef(f[j])/lcI[i]; 526 f = f - c*I[i]; 527 dbprint(ppl,"reducing poly to ",f); 528 if (remembercoeffs <> 0) 529 { 530 v = v - c * M[i]; 531 } 532 } 533 } 533 534 } 534 535 } … … 556 557 f = x3 + y2 + x2 + y + x; 557 558 I = x3 - y3, y3 - x2, x3 - y2, x2 - y, y2-x; 558 list l = linReduce(f, I, 1); 559 list l = linReduce(f, I, 1); 559 560 l; 560 561 module M = I; … … 601 602 if (p <> 0) 602 603 { 603 604 whichengine = 1; 604 605 } 605 606 } … … 668 669 proc pIntersect (poly s, ideal I) 669 670 "USAGE: pIntersect(f, I); f a poly, I an ideal 670 RETURN: vector, coefficient vector of the monic polynomial 671 RETURN: vector, coefficient vector of the monic polynomial 671 672 PURPOSE: compute the intersection of ideal I with the subalgebra K[f] 672 673 ASSUME: I is given as Groebner basis. … … 684 685 attrib(I,"isSB",1); // set attribute for suppressing NF messages 685 686 } 686 int ppl = printlevel-voice+ 2;687 int ppl = printlevel-voice+1; 687 688 // ---case 1: I = basering--- 688 689 if (size(I) == 1) … … 712 713 if (degs[j] == 0) 713 714 { 714 715 716 717 715 if (degI[i][j] <> 0) 716 { 717 break; 718 } 718 719 } 719 720 if (j == n) 720 721 { 721 722 722 k++; 723 possdegbounds[k] = Max(degI[i]); 723 724 } 724 725 } … … 749 750 if (tobracket==0) // todo bug in bracket? 750 751 { 751 752 toNF = 0; 752 753 } 753 754 else 754 755 { 755 756 toNF = bracket(tobracket,secNF); 756 757 } 757 758 newNF = NF(toNF+oldNF*secNF,I); // = NF(s^i,I) … … 802 803 v = v + m[j,1]*gen(j); 803 804 } 804 setring save; 805 setring save; 805 806 v = imap(@R,v); 806 807 kill @R; … … 864 865 if (size(#)>2) 865 866 { 866 867 868 869 867 if (typeof(#[3])=="int" || typeof(#[3])=="number") 868 { 869 modengine = int(#[3]); 870 } 870 871 } 871 872 } … … 909 910 if (tobracket!=0) 910 911 { 911 912 toNF = bracket(tobracket,NI[2]) + NI[i]*NI[2]; 912 913 } 913 914 else 914 915 { 915 916 toNF = NI[i]*NI[2]; 916 917 } 917 918 newNF = NF(toNF,I); … … 927 928 if (v!=0) // there is a modular solution 928 929 { 929 930 931 932 933 934 935 930 dbprint(ppl,"got solution in char ",solveincharp," of degree " ,i); 931 setring save; 932 v = linSyzSolve(NI,whichengine); 933 if (v==0) 934 { 935 break; 936 } 936 937 } 937 938 else // no modular solution 938 939 { 939 940 940 setring save; 941 v = 0; 941 942 } 942 943 } … … 965 966 { 966 967 dbprint(ppl,"linSyzSolve: got solution!"); 967 968 // "got solution!"; 968 969 break; 969 970 } … … 1015 1016 if (typeof(#[2])=="int" || typeof(#[2])=="number") 1016 1017 { 1017 1018 ringvar = int(#[2]); 1018 1019 } 1019 1020 } … … 1094 1095 def save = basering; 1095 1096 int n = nvars(save); 1097 if (size(variables(f)) == 0) // f is constant 1098 { 1099 return(list(ideal(0),intvec(0))); 1100 } 1096 1101 if (inorann == 0) // bfct using initial ideal 1097 1102 { … … 1121 1126 while (b == 0) 1122 1127 { 1123 1124 1125 1128 dbprint(ppl,"number of run in the loop: "+string(i)); 1129 int q = prime(random(lb,ub)); 1130 if (findFirst(usedprimes,q)==0) // if q was not already used 1126 1131 { 1127 1128 1129 1130 1131 1132 } 1133 } 1134 else // pIntersectSyz::non-modular 1132 usedprimes = usedprimes,q; 1133 dbprint(ppl,"used prime is: "+string(q)); 1134 b = pIntersectSyz(s,J,q,whichengine,modengine); 1135 } 1136 i++; 1137 } 1138 } 1139 else // pIntersectSyz::non-modular 1135 1140 { 1136 1141 b = pIntersectSyz(s,J,0,whichengine); … … 1158 1163 RETURN: list of ideal and intvec 1159 1164 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) 1160 @* for the hypersurface defined by f. 1165 @* for the hypersurface defined by f. 1161 1166 ASSUME: The basering is a commutative polynomial ring in char 0. 1162 BACKGROUND: In this proc, the initial Malgrange ideal is computed according to the algorithm 1167 BACKGROUND: In this proc, the initial Malgrange ideal is computed according to the algorithm 1163 1168 @* by Masayuki Noro and then a system of linear equations is solved by linear reductions. 1164 NOTE: In the output list, the ideal contains all the roots 1169 NOTE: In the output list, the ideal contains all the roots 1165 1170 @* and the intvec their multiplicities. 1166 1171 @* If s<>0, @code{std} is used for GB computations, 1167 @* otherwise, and by default, @code{slimgb} is used. 1172 @* otherwise, and by default, @code{slimgb} is used. 1168 1173 @* If t<>0, a matrix ordering is used for Groebner basis computations, 1169 1174 @* otherwise, and by default, a block ordering is used. 1170 @* If v is a positive weight vector, v is used for homogenization computations, 1175 @* If v is a positive weight vector, v is used for homogenization computations, 1171 1176 @* otherwise and by default, no weights are used. 1172 1177 DISPLAY: If printlevel=1, progress debug messages will be printed, … … 1199 1204 if (size(#)>2) 1200 1205 { 1201 1206 if (typeof(#[3])=="intvec" && size(#[3])==n && allPositive(#[3])==1) 1202 1207 { 1203 1204 1208 u0 = #[3]; 1209 } 1205 1210 } 1206 1211 } … … 1222 1227 "USAGE: bfctSyz(f [,r,s,t,u,v]); f a poly, r,s,t,u optional ints, v an optional intvec 1223 1228 RETURN: list of ideal and intvec 1224 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) 1229 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) 1225 1230 @* for the hypersurface defined by f 1226 1231 ASSUME: The basering is a commutative polynomial ring in char 0. 1227 BACKGROUND: In this proc, the initial Malgrange ideal is computed according to the algorithm 1232 BACKGROUND: In this proc, the initial Malgrange ideal is computed according to the algorithm 1228 1233 @* by Masayuki Noro and then a system of linear equations is solved by computing syzygies. 1229 NOTE: In the output list, the ideal contains all the roots 1234 NOTE: In the output list, the ideal contains all the roots 1230 1235 @* and the intvec their multiplicities. 1231 @* If r<>0, @code{std} is used for GB computations in characteristic 0, 1232 @* otherwise, and by default, @code{slimgb} is used. 1233 @* If s<>0, a matrix ordering is used for GB computations, otherwise, 1236 @* If r<>0, @code{std} is used for GB computations in characteristic 0, 1237 @* otherwise, and by default, @code{slimgb} is used. 1238 @* If s<>0, a matrix ordering is used for GB computations, otherwise, 1234 1239 @* and by default, a block ordering is used. 1235 @* If t<>0, the computation of the intersection is solely performed over 1240 @* If t<>0, the computation of the intersection is solely performed over 1236 1241 @* charasteristic 0, otherwise and by default, a modular method is used. 1237 @* If u<>0 and by default, @code{std} is used for GB computations in 1238 @* characteristic >0, otherwise, @code{slimgb} is used. 1239 @* If v is a positive weight vector, v is used for homogenization 1242 @* If u<>0 and by default, @code{std} is used for GB computations in 1243 @* characteristic >0, otherwise, @code{slimgb} is used. 1244 @* If v is a positive weight vector, v is used for homogenization 1240 1245 @* computations, otherwise and by default, no weights are used. 1241 1246 DISPLAY: If printlevel=1, progress debug messages will be printed, 1242 1247 @* if printlevel>=2, all the debug messages will be printed. 1243 EXAMPLE: example bfct ; shows examples1248 EXAMPLE: example bfctSyz; shows examples 1244 1249 " 1245 1250 { … … 1275 1280 if (typeof(#[3])=="int" || typeof(#[3])=="number") 1276 1281 { 1277 1278 1279 1282 pIntersectchar = int(#[3]); 1283 } 1284 if (size(#)>3) 1280 1285 { 1281 1286 if (typeof(#[4])=="int" || typeof(#[4])=="number") 1282 1287 { 1283 1284 1285 1286 1287 1288 1289 1290 1291 1292 1288 modengine = int(#[4]); 1289 } 1290 if (size(#)>4) 1291 { 1292 if (typeof(#[5])=="intvec" && size(#[5])==n && allPositive(#[5])==1) 1293 { 1294 u0 = #[5]; 1295 } 1296 } 1297 } 1293 1298 } 1294 1299 } … … 1315 1320 BACKGROUND: In this proc, the initial ideal of I is computed according to the algorithm by 1316 1321 @* Masayuki Noro and then a system of linear equations is solved by linear reductions. 1317 NOTE: In the output list, the ideal contains all the roots 1322 NOTE: In the output list, the ideal contains all the roots 1318 1323 @* and the intvec their multiplicities. 1319 1324 @* If s<>0, @code{std} is used for GB computations in characteristic 0, 1320 @* otherwise, and by default, @code{slimgb} is used. 1325 @* otherwise, and by default, @code{slimgb} is used. 1321 1326 @* If t<>0, a matrix ordering is used for GB computations, otherwise, 1322 1327 @* and by default, a block ordering is used. 1323 1328 DISPLAY: If printlevel=1, progress debug messages will be printed, 1324 1329 @* if printlevel>=2, all the debug messages will be printed. 1325 EXAMPLE: example bfct ideal; shows examples1330 EXAMPLE: example bfctIdeal; shows examples 1326 1331 " 1327 1332 { … … 1374 1379 "USAGE: bfctOneGB(f [,s,t]); f a poly, s,t optional ints 1375 1380 RETURN: list of ideal and intvec 1376 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) for the 1381 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) for the 1377 1382 @* hypersurface defined by f, using only one GB computation 1378 1383 ASSUME: The basering is a commutative polynomial ring in char 0. 1379 BACKGROUND: In this proc, the initial Malgrange ideal is computed based on the 1384 BACKGROUND: In this proc, the initial Malgrange ideal is computed based on the 1380 1385 @* algorithm by Masayuki Noro and combined with an elimination ordering. 1381 NOTE: In the output list, the ideal contains all the roots 1386 NOTE: In the output list, the ideal contains all the roots 1382 1387 @* and the intvec their multiplicities. 1383 1388 @* If s<>0, @code{std} is used for the GB computation, otherwise, 1384 @* and by default, @code{slimgb} is used. 1389 @* and by default, @code{slimgb} is used. 1385 1390 @* If t<>0, a matrix ordering is used for GB computations, 1386 1391 @* otherwise, and by default, a block ordering is used. … … 1492 1497 "USAGE: bfctAnn(f [,r]); f a poly, r an optional int 1493 1498 RETURN: list of ideal and intvec 1494 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) 1499 PURPOSE: computes the roots of the Bernstein-Sato polynomial b(s) 1495 1500 @* for the hypersurface defined by f 1496 1501 ASSUME: The basering is a commutative polynomial ring in char 0. 1497 1502 BACKGROUND: In this proc, ann(f^s) is computed and then a system of linear 1498 1503 @* equations is solved by linear reductions. 1499 NOTE: In the output list, the ideal contains all the roots 1504 NOTE: In the output list, the ideal contains all the roots 1500 1505 @* and the intvec their multiplicities. 1501 1506 @* If r<>0, @code{std} is used for GB computations, 1502 @* otherwise, and by default, @code{slimgb} is used. 1507 @* otherwise, and by default, @code{slimgb} is used. 1503 1508 DISPLAY: If printlevel=1, progress debug messages will be printed, 1504 1509 @* if printlevel>=2, all the debug messages will be printed. 1505 EXAMPLE: example bfct ann; shows examples1510 EXAMPLE: example bfctAnn; shows examples 1506 1511 " 1507 1512 {
Note: See TracChangeset
for help on using the changeset viewer.
