Changeset ea3c28a in git
- Timestamp:
- Jan 8, 2007, 1:59:38 PM (17 years ago)
- Branches:
- (u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
- Children:
- d1b006571ed1b6586799448402366ab672020d06
- Parents:
- f183dd53eec7359ae570b06c0f66b28fb8a64b3b
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
Singular/LIB/noether.lib
rf183dd rea3c28a 1 // AH /GP last modified: 12.05.20061 // AH last modified: 01.07.2007 2 2 ////////////////////////////////////////////////////////////////////////////// 3 version = "$Id: noether.lib,v 1. 3 2006-12-11 18:57:42 Singular Exp $";3 version = "$Id: noether.lib,v 1.4 2007-01-08 12:59:38 pfister Exp $"; 4 4 category="Commutative Algebra"; 5 5 info=" 6 LIBRARY: noether.lib Noether normalization of an ideal (not nessecary homogeneous) 6 LIBRARY: noether.lib Noether normalization of an ideal (not nessecary 7 homogeneous) 7 8 AUTHORS: A. Hashemi, Amir.Hashemi@lip6.fr 8 9 9 10 10 11 OVERVIEW: 11 A library for computing the Noether normalization of an ideal that DOES NOT require the computation of the dimension of the ideal. 12 It checks whether an ideal is in Noether position. A modular version of these algorithms is also provided. 13 The procedures are based on a paper of Amir Hashemi 'Efficient Algorithms for Computing Noether Normalization' 14 Submitted to: Special Issue of Mathematics in Computer Science on Symbolic and Numeric Computation. 15 16 This library computes also Castelnuovo-Mumford regularity and satiety of an ideal. A modular version of these algorithms is also provided. 17 The procedures are based on a paper of Amir Hashemi 'Computation of Castelnuovo-Mumford regularity and satiety' 12 A library for computing the Noether normalization of an ideal that DOES NOT 13 require the computation of the dimension of the ideal. 14 It checks whether an ideal is in Noether position. A modular version of 15 these algorithms is also provided. 16 The procedures are based on a paper of Amir Hashemi 'Efficient Algorithms for 17 Computing Noether Normalization' 18 Submitted to: Special Issue of Mathematics in Computer Science on Symbolic 19 and Numeric Computation. 20 21 This library computes also Castelnuovo-Mumford regularity and satiety of an 22 ideal. A modular version of these algorithms is also provided. 23 The procedures are based on a paper of Amir Hashemi 'Computation of 24 Castelnuovo-Mumford regularity and satiety' 18 25 Submitted to: ISSAC 2007. 19 26 … … 21 28 PROCEDURES: 22 29 NP_test(id); checks whether monomial ideal id is in Noether position 23 MNP_test(id); checks whether ideal id is in Noether position bymodular methods24 NP(id); Noether normalization of ideal id25 MNP(id); Noether normalization of ideal id by modular methods26 nsatiety(id); 30 MNP_test(id); the same as above using modular methods 31 NP(id); Noether normalization of ideal id 32 MNP(id); Noether normalization of ideal id by modular methods 33 nsatiety(id); Satiety of ideal id 27 34 msatiety(id) Satiety of ideal id by modular methods 28 regCM(id); Castelnuovo-Mumford regularity of ideal id29 mregCM(id); Castelnuovo-Mumford regularity of ideal id by modular methods35 regCM(id); Castelnuovo-Mumford regularity of ideal id 36 mregCM(id); Castelnuovo-Mumford regularity of ideal id by modular methods 30 37 "; 31 38 LIB "elim.lib"; … … 36 43 37 44 /////////////////////////////////////////////////////////////////////////////// 45 38 46 proc NP_test (ideal I) 39 47 " 40 USAGE: 41 RETURN: 42 0 otherwise. The second element of this list is the list of variable which 43 44 45 ASSUME: 48 USAGE: NP-test (I); I monomial ideal 49 RETURN: A list which first element is 1 if i is in Noether position 50 0 otherwise. The second element of this list is the list of variable which 51 its first part is the variable such that a power of this varaibles belong to the initial of i. 52 It return also the dimension of i if i is in Noether position 53 ASSUME: i is a nonzero monomial ideal. 46 54 " 47 55 { … … 54 62 if (I[1]==1) 55 63 { 56 64 print("The ideal is 1");return(1); 57 65 } 58 66 for ( ii = 1; ii <= n+1; ii++ ) 59 67 { 60 68 L[ii]=0; 61 69 } 62 70 for ( ii = 1; ii <= size(I); ii++ ) 63 71 { 64 65 66 if (size(Y[1])==1) 67 68 69 70 71 72 73 74 72 Y=findvars(I[ii],1)[1]; 73 l=rvar(Y[1][1]); 74 if (size(Y[1])==1) 75 { 76 L[l]=1; 77 P1=insert(P1,Y[1][1]); 78 } 79 if (L[l]==0) 80 { 81 L[l]=-1; 82 } 75 83 } 76 84 t=size(P1); … … 79 87 for ( jj = 1; jj <= n+1; jj++ ) 80 88 { 81 82 } 89 P3=insert(P3,varstr(jj)); 90 } 83 91 } 84 92 else … … 87 95 for ( jj = 1; jj <= size(P2[1]); jj++ ) 88 96 { 89 97 P3=insert(P3,P2[1][jj]); 90 98 } 91 99 } 92 100 if (L[n+1]==-1) 93 101 { 94 return(list(0,P1+P3)); 95 102 return(list(0,P1+P3)); 103 } 96 104 for ( ii = 1; ii <= n; ii++ ) 97 105 { 98 106 if (L[ii]==-1) 99 107 { 100 return(list(0,P1+P3)); 101 102 108 return(list(0,P1+P3)); 109 } 110 if (L[ii]==0 and L[ii+1]==1) 103 111 { 104 return(list(0,P1+P3)); 105 112 return(list(0,P1+P3)); 113 } 106 114 } 107 115 d=n+1-sum(L); … … 112 120 ////////////////////////////////////////// 113 121 proc MNP_test (ideal i) 114 "USAGE: 115 RETURN: 116 NOTE: This test is a probabilistic test, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31).122 "USAGE: MNP-test (i); i an ideal 123 RETURN: 1 if i is in Noether position 0 otherwise. 124 NOTE: This test is a probabilistic test, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31). 117 125 " 118 126 { … … 145 153 " 146 154 USAGE: NP (i); i ideal 147 RETURN: A linear map phi such that phi(i) is in Noether position 155 RETURN: A linear map phi such that phi(i) is in Noether position 148 156 " 149 157 { … … 168 176 if ( #[1]== 1 ) 169 177 { 170 178 return ("The ideal is in Noether position and the time of this computation is:",rtimer-time,"/10 sec."); 171 179 } 172 180 else 173 181 { 174 175 176 for ( ii = 1; ii<=n+1; ii++ )177 178 179 180 182 L=maxideal(1); 183 chcoord=maxideal(1); 184 for ( ii = 1; ii<=n+1; ii++ ) 185 { 186 chcoord[rvar(#[2][ii])]=L[ii]; 187 } 188 phi=r1,chcoord; 181 189 sbi=phi(sbi); 182 190 if ( NP_test(sbi)[1] == 1 ) 183 184 185 186 191 { 192 setring r0; 193 chcoord=fetch(r1,chcoord); 194 return (chcoord,"and the time of this computation is:",rtimer-time,"/10 sec."); 187 195 } 188 196 } … … 190 198 { 191 199 nl=nl+1; 192 I=i; 200 I=i; 193 201 L=maxideal(1); 194 202 for ( ii = n; ii>=0; ii-- ) … … 201 209 for ( jj = 1; jj<=ii+1; jj++ ) 202 210 { 203 211 P=P+ran[1,jj]*m[jj]; 204 212 } 205 213 chcoord[ii+1]=P; 206 214 L[ii+1]=P; 207 215 P=0; 208 216 phi=r1,chcoord; 209 217 I=phi(I); 210 218 if ( NP_test(sort(lead(std(I)))[1])[1] == 1 ) 211 212 213 214 215 ideal L=fetch(r1,L); 216 219 { 220 K=x(ii..n); 221 setring r0; 222 K=fetch(r1,K); 223 ideal L=fetch(r1,L); 224 return (L,"and the time of this computation is:",rtimer-time,"/10 sec."); 217 225 } 218 226 … … 226 234 /////////////////////////////////////////////////////////////////////////////// 227 235 proc MNP (ideal i) 228 "USAGE: 229 RETURN: 230 NOTE: It uses the procedure MNP_test to test Noether position.236 "USAGE: MNP (i); i ideal 237 RETURN: A linear map phi such that phi(i) is in Noether position 238 NOTE: It uses the procedure MNP_test to test Noether position. 231 239 " 232 240 { … … 246 254 time=rtimer; 247 255 system("--ticks-per-sec",10); 248 256 #=MNP_test(i); 249 257 if ( #[1]== 1 ) 250 258 { 251 259 return ("The ideal is in Noether position and the time of this computation is:",rtimer-time,"/10 sec."); 252 260 } 253 261 else 254 262 { 255 256 257 for ( ii = 1; ii<=n+1; ii++ )258 259 260 261 263 L=maxideal(1); 264 chcoord=maxideal(1); 265 for ( ii = 1; ii<=n+1; ii++ ) 266 { 267 chcoord[rvar(#[2][ii])]=L[ii]; 268 } 269 phi=r1,chcoord; 262 270 I=phi(i); 263 271 if ( MNP_test(I)[1] == 1 ) 264 265 266 267 272 { 273 setring r0; 274 chcoord=fetch(r1,chcoord); 275 return (chcoord,"and the time of this computation is:",rtimer-time,"/10 sec."); 268 276 } 269 277 } … … 271 279 { 272 280 nl=nl+1; 273 I=i; 281 I=i; 274 282 L=maxideal(1); 275 283 for ( ii = n; ii>=0; ii-- ) … … 282 290 for ( jj = 1; jj<=ii+1; jj++ ) 283 291 { 284 292 P=P+ran[1,jj]*m[jj]; 285 293 } 286 294 chcoord[ii+1]=P; … … 290 298 I=phi(I); 291 299 if ( MNP_test(I)[1] == 1 ) 292 293 294 295 296 ideal L=fetch(r1,L); 297 300 { 301 K=x(ii..n); 302 setring r0; 303 K=fetch(r1,K); 304 ideal L=fetch(r1,L); 305 return (L,"and the time of this computation is:",rtimer-time,"/10 sec."); 298 306 } 299 307 … … 309 317 proc Test (ideal i) 310 318 " 311 USAGE: Test (i); i a monomial ideal, 319 USAGE: Test (i); i a monomial ideal, 312 320 RETURN: 1 if the last variable is in generic position for i and 0 otherwise. 313 THEORY: The last variable is in generic position if the quotient of the ideal 314 321 THEORY: The last variable is in generic position if the quotient of the ideal 322 with respect to this variable is equal to the quotient of the ideal with respect to the maximal ideal. 315 323 " 316 324 { … … 329 337 if (size(reduce(I,i)) <> 0) 330 338 { 331 339 ret=0; 332 340 } 333 341 return(ret); … … 338 346 proc nsatiety (ideal i) 339 347 " 340 USAGE: nsatiety (i); i ideal, 348 USAGE: nsatiety (i); i ideal, 341 349 RETURN: an integer, the satiety of i. 342 350 (returns -1 if i is not homogeneous) … … 344 352 THEORY: The satiety, or saturation index, of a homogeneous ideal i is the 345 353 least integer s such that, for all d>=s, the degree d part of the 346 ideals i and isat iety=satiety(i,maxideal(1))[1] coincide.354 ideals i and isat=sat(i,maxideal(1))[1] coincide. 347 355 " 348 356 { … … 365 373 { 366 374 dbprint(2,"The ideal is not homogeneous, and time for this test is: " + string(rtimer-time) + "/100sec."); 367 375 return (); 368 376 } 369 377 I=simplify(lead(sbi),1); … … 372 380 if (size(K) == 0) 373 381 { 374 dbprint(2,"satiety(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");375 382 dbprint(2,"sat(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec."); 383 return(); 376 384 } 377 385 if (Test(I) == 1 ) 378 386 { 379 dbprint(2,"satiety(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");380 387 dbprint(2,"sat(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec."); 388 return(); 381 389 } 382 390 while ( nl < 5 ) 383 391 { 384 nl=nl+1; 392 nl=nl+1; 385 393 chcoord=select1(maxideal(1),1,(n)); 386 394 ran=random(100,1,n); … … 390 398 for ( jj = 1; jj<=n+1; jj++ ) 391 399 { 392 400 P=P+ran[1,jj]*m[jj]; 393 401 } 394 402 chcoord[n+1]=P; … … 398 406 I=simplify(lead(L),1); 399 407 attrib(I,"isSB",1); 400 K=select(I,n+1); 408 K=select(I,n+1); 401 409 if (size(K) == 0) 402 410 { 403 dbprint(2,"satiety(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");404 411 dbprint(2,"sat(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec."); 412 return(); 405 413 } 406 414 if (Test(I) == 1 ) 407 415 { 408 dbprint(2,"satiety(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");409 416 dbprint(2,"sat(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec."); 417 return(); 410 418 } 411 419 } … … 416 424 proc msatiety (ideal i) 417 425 " 418 USAGE: msatiety (i); i ideal, 426 USAGE: msatiety (i); i ideal, 419 427 RETURN: an integer, the satiety of i. 420 428 (returns -1 if i is not homogeneous) … … 422 430 THEORY: The satiety, or saturation index, of a homogeneous ideal i is the 423 431 least integer s such that, for all d>=s, the degree d part of the 424 ideals i and isat iety=satiety(i,maxideal(1))[1] coincide.425 NOTE: 432 ideals i and isat=sat(i,maxideal(1))[1] coincide. 433 NOTE: This is a probabilistic procedure, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31). 426 434 " 427 435 { 428 436 //--------------------------- initialisation --------------------------------- 429 "// WARNING: 430 // The procedure is probabilistic and it computes the initial of the ideal modulo the prime number 2147483647"; 437 "// WARNING: The characteristic of base field must be zero. 438 // The procedure is probabilistic and it computes the 439 //initial ideals modulo the prime number 2147483647."; 431 440 int e,ii,jj,h,d,time,lastv,nl,ret,s1,d1,siz,j,si,u,k,p; 432 441 intvec v1; … … 449 458 "// WARNING: The ideal is not homogeneous."; 450 459 dbprint(2,"Time for this test is: " + string(rtimer-time) + "/100sec."); 451 460 return (); 452 461 } 453 462 option(redSB); … … 467 476 if (size(K) == 0) 468 477 { 469 dbprint(2,"msatiety(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");470 478 dbprint(2,"msat(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec."); 479 return(); 471 480 } 472 481 if (Test(I) == 1 ) 473 482 { 474 dbprint(2,"msatiety(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");475 483 dbprint(2,"msat(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec."); 484 return(); 476 485 } 477 486 while ( nl < 30 ) 478 487 { 479 nl=nl+1; 488 nl=nl+1; 480 489 chcoord=select1(maxideal(1),1,(n)); 481 490 ran=random(100,1,n); … … 485 494 for ( jj = 1; jj<=n+1; jj++ ) 486 495 { 487 496 P=P+ran[1,jj]*m[jj]; 488 497 } 489 498 chcoord[n+1]=P; … … 502 511 lsbi1=lead(sbi); 503 512 attrib(lsbi1,"isSB",1); 504 K=select(lsbi1,n+1); 513 K=select(lsbi1,n+1); 505 514 if (size(K) == 0) 506 515 { 507 dbprint(2,"msatiety(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");508 516 dbprint(2,"msat(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec."); 517 return(); 509 518 } 510 519 if (Test(lsbi1) == 1 ) 511 520 { 512 dbprint(2,"msatiety(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");513 521 dbprint(2,"msat(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec."); 522 return(); 514 523 } 515 524 } … … 518 527 ////////////////////////////////////////////////////////////////////////////// 519 528 // 520 proc reg CM(ideal i)521 " 522 USAGE: regCM(i); i ideal523 RETURN: 529 proc reg (ideal i) 530 " 531 USAGE: reg (i); i ideal 532 RETURN: the Castelnuovo-Mumford regularity of i. 524 533 (returns -1 if i is not homogeneous) 525 534 ASSUME: i is a homogeneous ideal. … … 527 536 { 528 537 //--------------------------- initialisation --------------------------------- 529 int e,ii,jj,H,h,d,time,lastv,nv,nl; 530 int lastind,ch,nesttest,NPtest,nl,N,acc; 531 intmat ran; 538 int e,ii,jj,H,h,d,time,nl; 532 539 def r0 = basering; 533 540 int n = nvars(r0)-1; 534 541 string s = "ring r1 = ",charstr(r0),",x(0..n),dp;"; 535 542 execute(s); 536 ideal i,sbi,I,J,K, chcoord,m,L;543 ideal i,sbi,I,J,K,L; 537 544 list #; 538 545 poly P; … … 540 547 i = fetch(r0,i); 541 548 time=rtimer; 549 system("--ticks-per-sec",100); 542 550 sbi=std(i); 543 551 //----- Check ideal homogeneous … … 550 558 attrib(I,"isSB",1); 551 559 d=dim(I); 560 if (char(r1) > 0 and d == 0) 561 { 562 def r2=changechar("0",r1); 563 setring r2; 564 ideal sbi,I,i,K,T; 565 map phi; 566 I = fetch(r1,I); 567 i=I; 568 attrib(I,"isSB",1); 569 } 570 else 571 { 572 def r2=changechar(charstr(r1),r1); 573 setring r2; 574 ideal sbi,I,i,K,T,ic,Ic; 575 map phi; 576 I = imap(r1,I); 577 Ic=I; 578 attrib(I,"isSB",1); 579 i = imap(r1,i); 580 ic=i; 581 } 552 582 K=select(I,n+1); 553 583 if (size(K) == 0) 554 584 { 555 h=0; 556 } 557 if (Test(I) == 1) 558 { 559 h=maxdeg1(K); 585 h=0; 560 586 } 561 587 else 562 588 { 563 while ( nl < 5 ) 564 { 565 nl=nl+1; 566 chcoord=select1(maxideal(1),1,(n)); 567 ran=random(100,1,n); 568 ran=intmat(ran,1,n+1); 569 ran[1,n+1]=1; 570 m=select1(maxideal(1),1,(n+1)); 571 for ( jj = 1; jj<=n+1; jj++ ) 572 { 573 P=P+ran[1,jj]*m[jj]; 574 } 575 chcoord[n+1]=P; 576 P=0; 577 phi=r1,chcoord; 578 i=phi(i); 579 I=simplify(lead(std(i)),1); 580 attrib(I,"isSB",1); 581 K=select(I,n+1); 582 if (size(K) == 0) 583 { 584 h=0;break 585 } 586 if (Test(I) == 1 ) 587 { 588 h=maxdeg1(K);break; 589 } 590 } 591 } 589 if (Test(I) == 1) 590 { 591 h=maxdeg1(K); 592 } 593 else 594 { 595 while ( nl < 30 ) 596 { 597 nl=nl+1; 598 phi=r2,randomLast(100); 599 T=phi(i); 600 I=simplify(lead(std(T)),1); 601 attrib(I,"isSB",1); 602 K=select(I,n+1); 603 if (size(K) == 0) 604 { 605 h=0;break 606 } 607 if (Test(I) == 1 ) 608 { 609 h=maxdeg1(K);break; 610 } 611 } 612 i=T; 613 } 614 } 592 615 for ( ii = n; ii>=n-d+1; ii-- ) 593 616 { 594 i=subst(i,x(ii),0); 595 s = "ring mr = ",charstr(r0),",x(0..ii-1),dp;"; 596 execute(s); 597 ideal i,sbi,I,J,K,chcoord,m,L; 598 poly P; 599 map phi; 600 i=imap(r1,i); 601 sbi=std(i); 602 I=simplify(lead(sbi),1); 603 attrib(I,"isSB",1); 604 K=select(I,ii); 605 if (size(K) == 0) 606 { 607 H=0; 608 } 609 if (Test(I) == 1) 610 { 611 H=maxdeg1(K); 612 } 613 else 614 { 615 while ( nl < 5 ) 616 { 617 nl=nl+1; 618 chcoord=select1(maxideal(1),1,(ii-1)); 619 ran=random(100,1,ii-1); 620 ran=intmat(ran,1,ii); 621 ran[1,ii]=1; 622 m=select1(maxideal(1),1,(ii)); 623 for ( jj = 1; jj<=ii; jj++ ) 624 { 625 P=P+ran[1,jj]*m[jj]; 626 } 627 chcoord[ii]=P; 628 P=0; 629 phi=mr,chcoord; 630 i=phi(i); 631 I=simplify(lead(std(i)),1); 632 attrib(I,"isSB",1); 633 K=select(I,ii); 634 if (size(K) == 0) 635 { 636 H=0;break; 637 } 638 if (Test(I) == 1 ) 639 { 640 H=maxdeg1(K);break; 641 } 642 } 643 setring r1; 644 i=imap(mr,i); 645 kill mr; 646 } 647 if (H > h) 648 { 649 h=H; 650 } 651 } 652 dbprint(2,"regCM(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + "/100sec."); 653 return(); 654 } 617 i=subst(i,x(ii),0); 618 s = "ring mr = ",charstr(r1),",x(0..ii-1),dp;"; 619 execute(s); 620 ideal i,sbi,I,J,K,L,T; 621 poly P; 622 map phi; 623 i=imap(r2,i); 624 I=simplify(lead(std(i)),1); 625 attrib(I,"isSB",1); 626 K=select(I,ii); 627 if (size(K) == 0) 628 { 629 H=0; 630 } 631 else 632 { 633 if (Test(I) == 1) 634 { 635 H=maxdeg1(K); 636 } 637 else 638 { 639 while ( nl < 30 ) 640 { 641 nl=nl+1; 642 phi=mr,randomLast(100); 643 T=phi(i); 644 I=simplify(lead(std(T)),1); 645 attrib(I,"isSB",1); 646 K=select(I,ii); 647 if (size(K) == 0) 648 { 649 H=0;break; 650 } 651 if (Test(I) == 1 ) 652 { 653 H=maxdeg1(K);break; 654 } 655 } 656 setring r2; 657 i=imap(mr,T); 658 kill mr; 659 } 660 } 661 if (H > h) 662 { 663 h=H; 664 } 665 } 666 if (nl < 30) 667 { 668 dbprint(2,"reg(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + " sec./100"); 669 return(); 670 } 671 else 672 { 673 I=Ic; 674 attrib(I,"isSB",1); 675 i=ic; 676 K=subst(select(I,n+1),x(n),1); 677 K=K*maxideal(maxdeg1(I)); 678 if (size(reduce(K,I)) <> 0) 679 { 680 nl=0; 681 while ( nl < 30 ) 682 { 683 nl=nl+1; 684 phi=r1,randomLast(100); 685 sbi=phi(i); 686 I=simplify(lead(std(sbi)),1); 687 attrib(I,"isSB",1); 688 K=subst(select(I,n+1),x(n),1); 689 K=K*maxideal(maxdeg1(I)); 690 if (size(reduce(K,I)) == 0) 691 { 692 break; 693 } 694 } 695 } 696 h=maxdeg1(simplify(reduce(quotient(I,maxideal(1)),I),2))+1; 697 for ( ii = n; ii> n-d+1; ii-- ) 698 { 699 sbi=subst(sbi,x(ii),0); 700 s = "ring mr = ",charstr(r0),",x(0..ii-1),dp;"; 701 execute(s); 702 ideal sbi,I,L,K,T; 703 map phi; 704 sbi=imap(r1,sbi); 705 I=simplify(lead(std(sbi)),1); 706 attrib(I,"isSB",1); 707 K=subst(select(I,ii),x(ii-1),1); 708 K=K*maxideal(maxdeg1(I)); 709 if (size(reduce(K,I)) <> 0) 710 { 711 nl=0; 712 while ( nl < 30 ) 713 { 714 nl=nl+1; 715 L=randomLast(100); 716 phi=mr,L; 717 T=phi(sbi); 718 I=simplify(lead(std(T)),1); 719 attrib(I,"isSB",1); 720 K=subst(select(I,ii),x(ii-1),1); 721 K=K*maxideal(maxdeg1(I)); 722 if (size(reduce(K,I)) == 0) 723 { 724 sbi=T; 725 break; 726 } 727 } 728 } 729 H=maxdeg1(simplify(reduce(quotient(I,maxideal(1)),I),2))+1; 730 if (H > h) 731 { 732 h=H; 733 } 734 setring r1; 735 sbi=fetch(mr,sbi); 736 kill mr; 737 } 738 sbi=subst(sbi,x(n-d+1),0); 739 s = "ring mr = ",charstr(r0),",x(0..n-d),dp;"; 740 execute(s); 741 ideal sbi,I,L,K,T; 742 map phi; 743 sbi=imap(r1,sbi); 744 I=simplify(lead(std(sbi)),1); 745 attrib(I,"isSB",1); 746 H=maxdeg1(simplify(reduce(quotient(I,maxideal(1)),I),2))+1; 747 if (H > h) 748 { 749 h=H; 750 } 751 dbprint(2,"reg(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + " sec./100"); 752 return(); 753 } 754 } 755 655 756 ////////////////////////////////////////////////////////////////////////////// 656 757 // 657 proc mreg CM(ideal i)658 " 659 USAGE: mregCM(i); i ideal660 RETURN: 758 proc mreg (ideal i) 759 " 760 USAGE: mreg (i); i ideal 761 RETURN: an integer, the Castelnuovo-Mumford regularity of i. 661 762 (returns -1 if i is not homogeneous) 662 ASSUME: i is a homogeneous ideal .663 NOTE: 763 ASSUME: i is a homogeneous ideal and the characteristic of base field is zero.. 764 NOTE: This is a probabilistic procedure, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31). 664 765 " 665 766 { 666 767 //--------------------------- initialisation --------------------------------- 667 "// WARNING: 668 // The procedure is probabilistic and it computes the initial of the ideal modulo the prime number 2147483647"; 669 int e,ii,jj,H,h,d,time,lastv,sat,firstind,nv,nl,p; 670 int lastind,ch,nesttest,NPtest,nl,N,acc; 671 intmat ran; 768 "// WARNING: The characteristic of base field musr be zero. 769 // This procedure is probabilistic and it computes the initial 770 //ideals modulo the prime number 2147483647"; 771 int e,ii,jj,H,h,d,time,p,nl; 672 772 def r0 = basering; 673 773 int n = nvars(r0)-1; 674 774 string s = "ring r1 = ",charstr(r0),",x(0..n),dp;"; 675 775 execute(s); 676 ideal i,sbi,I,J,K, chcoord,m,L,lsbi1,lsbi2;776 ideal i,sbi,I,J,K,L,lsbi1,lsbi2; 677 777 list #; 678 778 poly P; … … 680 780 i = fetch(r0,i); 681 781 time=rtimer; 782 system("--ticks-per-sec",100); 682 783 //----- Check ideal homogeneous 683 784 if ( homog(i) == 0 ) … … 702 803 I=lsbi1; 703 804 d=dim(I); 704 K=select(I,n+1); 805 K=select(I,n+1); 705 806 if (size(K) == 0) 706 807 { 707 h=0; 708 } 709 if (Test(I) == 1) 710 { 711 h=maxdeg1(K); 808 h=0; 712 809 } 713 810 else 714 811 { 715 while ( nl < 5 ) 716 { 717 nl=nl+1; 718 chcoord=select1(maxideal(1),1,(n)); 719 ran=random(100,1,n); 720 ran=intmat(ran,1,n+1); 721 ran[1,n+1]=1; 722 m=select1(maxideal(1),1,(n+1)); 723 for ( jj = 1; jj<=n+1; jj++ ) 724 { 725 P=P+ran[1,jj]*m[jj]; 726 } 727 chcoord[n+1]=P; 728 P=0; 729 phi=r1,chcoord; 730 i=phi(i); 731 #=ringlist(r1); 732 #[1]=p; 733 def oro=ring(#); 734 setring oro; 735 ideal sbi,lsbi; 736 sbi=fetch(r1,i); 737 lsbi=lead(std(sbi)); 738 setring r1; 739 lsbi1=fetch(oro,lsbi); 740 lsbi1=simplify(lsbi1,1); 741 attrib(lsbi1,"isSB",1); 742 kill oro; 743 I=lsbi1; 744 K=select(I,n+1); 745 if (size(K) == 0) 746 { 747 h=0;break 748 } 749 if (Test(I) == 1 ) 750 { 751 h=maxdeg1(K);break; 752 } 753 } 754 } 812 if (Test(I) == 1) 813 { 814 h=maxdeg1(K); 815 } 816 else 817 { 818 while ( nl < 30 ) 819 { 820 nl=nl+1; 821 phi=r1,randomLast(100); 822 sbi=phi(i); 823 #=ringlist(r1); 824 #[1]=p; 825 def oro=ring(#); 826 setring oro; 827 ideal sbi,lsbi; 828 sbi=fetch(r1,sbi); 829 lsbi=lead(std(sbi)); 830 setring r1; 831 lsbi1=fetch(oro,lsbi); 832 lsbi1=simplify(lsbi1,1); 833 attrib(lsbi1,"isSB",1); 834 kill oro; 835 I=lsbi1; 836 K=select(I,n+1); 837 if (size(K) == 0) 838 { 839 h=0;break 840 } 841 if (Test(I) == 1 ) 842 { 843 h=maxdeg1(K);break; 844 } 845 } 846 i=sbi; 847 } 848 } 755 849 for ( ii = n; ii>=n-d+1; ii-- ) 756 850 { 757 758 s = "ring mr = ",charstr(r0),",x(0..ii-1),dp;"; 759 760 ideal i,sbi,I,J,K,chcoord,m,L,lsbi1;761 851 i=subst(i,x(ii),0); 852 s = "ring mr = ","0",",x(0..ii-1),dp;"; 853 execute(s); 854 ideal i,sbi,I,J,K,L,lsbi1; 855 poly P; 762 856 list #; 763 map phi; 764 i=imap(r1,i); 765 #=ringlist(mr); 766 #[1]=p; 767 def oro=ring(#); 768 setring oro; 769 ideal sbi,lsbi; 770 sbi=fetch(mr,i); 771 lsbi=lead(std(sbi)); 772 setring mr; 773 lsbi1=fetch(oro,lsbi); 774 lsbi1=simplify(lsbi1,1); 775 attrib(lsbi1,"isSB",1); 776 kill oro; 777 I=lsbi1; 778 K=select(I,ii); 779 if (size(K) == 0) 780 { 781 H=0; 782 } 783 if (Test(I) == 1) 784 { 785 H=maxdeg1(K); 786 } 787 else 788 { 789 while ( nl < 5 ) 790 { 791 nl=nl+1; 792 chcoord=select1(maxideal(1),1,(ii-1)); 793 ran=random(100,1,ii-1); 794 ran=intmat(ran,1,ii); 795 ran[1,ii]=1; 796 m=select1(maxideal(1),1,(ii)); 797 for ( jj = 1; jj<=ii; jj++ ) 798 { 799 P=P+ran[1,jj]*m[jj]; 800 } 801 chcoord[ii]=P; 802 P=0; 803 phi=mr,chcoord; 804 i=phi(i); 805 #=ringlist(mr); 806 #[1]=p; 807 def oro=ring(#); 808 setring oro; 809 ideal sbi,lsbi; 810 sbi=fetch(mr,i); 811 lsbi=lead(std(sbi)); 812 setring mr; 813 lsbi1=fetch(oro,lsbi); 814 lsbi1=simplify(lsbi1,1); 815 attrib(lsbi1,"isSB",1); 816 kill oro; 817 I=lsbi1; 818 K=select(I,ii); 819 if (size(K) == 0) 820 { 821 H=0;break; 822 } 823 if (Test(I) == 1 ) 824 { 825 H=maxdeg1(K);break; 826 } 827 } 828 setring r1; 829 i=imap(mr,i); 830 kill mr; 831 } 832 if (H > h) 833 { 834 h=H; 835 } 836 } 837 dbprint(2,"mregCM(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + "/100sec."); 857 map phi; 858 i=imap(r1,i); 859 #=ringlist(mr); 860 #[1]=p; 861 def oro=ring(#); 862 setring oro; 863 ideal sbi,lsbi; 864 sbi=fetch(mr,i); 865 lsbi=lead(std(sbi)); 866 setring mr; 867 lsbi1=fetch(oro,lsbi); 868 lsbi1=simplify(lsbi1,1); 869 attrib(lsbi1,"isSB",1); 870 kill oro; 871 I=lsbi1; 872 K=select(I,ii); 873 if (size(K) == 0) 874 { 875 H=0; 876 } 877 else 878 { 879 if (Test(I) == 1) 880 { 881 H=maxdeg1(K); 882 } 883 else 884 { 885 nl=0; 886 while ( nl < 30 ) 887 { 888 nl=nl+1; 889 phi=mr,randomLast(100); 890 sbi=phi(i); 891 #=ringlist(mr); 892 #[1]=p; 893 def oro=ring(#); 894 setring oro; 895 ideal sbi,lsbi; 896 sbi=fetch(mr,sbi); 897 lsbi=lead(std(sbi)); 898 setring mr; 899 lsbi1=fetch(oro,lsbi); 900 lsbi1=simplify(lsbi1,1); 901 kill oro; 902 I=lsbi1; 903 attrib(I,"isSB",1); 904 K=select(I,ii); 905 if (size(K) == 0) 906 { 907 H=0;break; 908 } 909 if (Test(I) == 1 ) 910 { 911 H=maxdeg1(K);break; 912 } 913 } 914 setring r1; 915 i=imap(mr,sbi); 916 kill mr; 917 } 918 } 919 if (H > h) 920 { 921 h=H; 922 } 923 } 924 dbprint(2,"mreg(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + "sec./100"); 838 925 return(); 839 926 } 927 928 929 930 931 932 933 ////////////////////////////////////////////////////////////// 840 934 example 841 935 { "EXAMPLE:"; echo = 2; … … 917 1011 } 918 1012 919
Note: See TracChangeset
for help on using the changeset viewer.