Changeset 3517649 in git
- Timestamp:
- Jun 12, 2009, 10:01:11 AM (14 years ago)
- Branches:
- (u'jengelh-datetime', 'ceac47cbc86fe4a15902392bdbb9bd2ae0ea02c6')(u'spielwiese', '00e2e9c41af3fde1273eb3633f4c0c7c3db2579d')
- Children:
- 008e8143a13845539fcd39e47a9b1267e6d78265
- Parents:
- b5d4d19ebfd46b641a2428d8d8f86d9d6dafa6da
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
Singular/LIB/modstd.lib
rb5d4d1 r3517649 1 1 //GP, last modified 23.10.06 2 2 /////////////////////////////////////////////////////////////////////////////// 3 version="$Id: modstd.lib,v 1.1 5 2009-04-15 11:14:36 seelischExp $";3 version="$Id: modstd.lib,v 1.16 2009-06-12 08:01:11 Singular Exp $"; 4 4 category="Commutative Algebra"; 5 5 info=" … … 28 28 LIB "crypto.lib"; 29 29 /////////////////////////////////////////////////////////////////////////////// 30 proc pTestSB(ideal I, ideal J, list L) 31 "USAGE: pTestSB(I,J,L); I,J ideals, L intvec of primes 32 RETURN: 1 (resp. 0) if for a randomly chosen prime p not in L 33 J mod p is (resp. is not) a standard basis of I mod p 34 EXAMPLE:example primList; shows an example 35 " 36 { 37 int i,j,k,p; 38 def R=basering; 39 list r= ringlist(R); 40 41 while(!j) 42 { 43 j=1; 44 p=prime(random(1000000000,2134567879)); 45 for(i=1;i<=size(L);i++) 46 { 47 if(p==L[i]){j=0;break} 48 } 49 if(j) 50 { 51 for(i=1;i<=ncols(J);i++) 52 { 53 for(k=2;k<=size(J[i]);k++) 54 { 55 if((denominator(leadcoef(J[i][k])) mod 56 p)==0){j=0;break} 57 } 58 if(!j){break;} 59 } 60 } 61 } 62 r[1]=p; 63 def @R=ring(r); 64 setring @R; 65 ideal I=imap(R,I); 66 ideal J=imap(R,J); 67 attrib(J,"isSB",1); 68 j=1; 69 if(size(reduce(I,J))!=0){j=0;} 70 if(j) 71 { 72 ideal K=std(J); 73 if(size(reduce(K,J))!=0){j=0;} 74 } 75 setring R; 76 return(j); 77 } 78 example 79 { "EXAMPLE:"; echo = 2; 80 intvec L=2,3,5; 81 ring r=0,(x,y,z),dp; 82 ideal i=x+1,x+y+1; 83 ideal j=x+1,y; 84 pTestSB(i,i,L); 85 pTestSB(i,j,L); 86 } 87 88 proc primeList(int n, list #) 89 "USAGE: primeList(n); (resp. primeList(n,L); ) 90 RETURN: the intvec of n greatest primes <= 2147483647 (resp. n greatest 91 primes < L[size(L)] union with L) 92 EXAMPLE:example primList; shows an example 93 " 94 { 95 intvec L; 96 int i,p; 97 if(size(#)==0) 98 { 99 p=2147483647; 100 L[1]=p; 101 } 102 else 103 { 104 L=#[1]; 105 p=prime(L[size(L)]-1); 106 L[size(L)+1]=p; 107 } 108 if(p==2){ERROR("no more primes");} 109 for(i=2;i<=n;i++) 110 { 111 p=prime(p-1); 112 L[size(L)+1]=p; 113 } 114 return(L); 115 } 116 example 117 { "EXAMPLE:"; echo = 2; 118 intvec L=primeList(10); 119 size(L); 120 L[size(L)]; 121 L=primeList(5,L); 122 size(L); 123 L[size(L)]; 124 } 125 126 30 127 proc modStd(ideal I) 31 128 "USAGE: modStd(I); … … 33 130 NOTE: the procedure computes a standard basis of I (over the 34 131 rational numbers) by using modular methods. If a 35 warning appears then the result is a standard basis with no defined36 relation to I; this is a sign that not enough prime numbers have37 been used.For further experiments see procedure modS.132 warning appears then the result is a standard basis 133 containing I and with high probability a standard basis of I. 134 For further experiments see procedure modS. 38 135 EXAMPLE: example modStd; shows an example 39 136 " … … 43 140 if((npars(R0)>0)||(rl[1]>0)) 44 141 { 45 ERROR("characteristic of basering should be zero"); 46 } 47 int l,j,k,q; 142 ERROR("characteristic of basering should be zero, basering should have 143 no parameters"); 144 } 145 146 int k,c; 147 int pd=printlevel-voice+2; 148 int j=1; 149 int h=homog(I); 48 150 int en=2134567879; 49 151 int an=1000000000; 50 intvec hi,hl,hc,hpl,hpc; 51 list T,TT; 52 intvec L=primeList(5); 53 L[6]=prime(random(an,en)); 54 ideal J,cT,lT,K; 55 ideal I0=I; 56 int h=homog(I); 57 if((!h)&&(ord_test(R0)==0)) 58 { 59 ERROR("input is not homogeneous and ordering is not local"); 60 } 61 if(h) 62 { 63 execute("ring gn="+string(L[6])+",x(1.."+string(nvars(R0))+"),dp;"); 64 ideal I=fetch(R0,I); 65 ideal J=std(I); 66 hi=hilb(J,1); 67 setring R0; 68 } 69 for (j=1;j<=size(L);j++) 70 { 71 rl[1]=L[j]; 72 def oro=ring(rl); 73 setring oro; 74 ideal I=fetch(R0,I); 75 option(redSB); 76 if(h) 77 { 78 ideal I1=std(I,hi); 79 } 80 else 81 { 82 if(ord_test(R0)==-1) 83 { 84 ideal I1=std(I); 85 } 86 else 87 { 88 matrix M; 89 ideal I1=liftstd(I,M); 90 } 91 } 92 setring R0; 93 T[j]=fetch(oro,I1); 94 kill oro; 95 } 152 153 intvec opt = option(get); // Save current options 154 intvec L=primeList(10); 155 L[5]=prime(random(an,en)); 156 list T,TT,TH,LL; 157 158 159 ideal J,K; 160 161 option(redSB); 162 163 while(1) 164 { 165 while(j<=size(L)) 166 { 167 if(pd>2){c++;c;} 168 rl[1]=L[j]; 169 def @r=ring(rl); 170 setring @r; 171 ideal i=fetch(R0,I); 172 i=groebner(i); 173 setring R0; 174 T[size(T)+1]=fetch(@r,i); 175 kill @r; 176 j++; 177 } 178 if(pd>2){"lifting";} 96 179 //================= delete unlucky primes ==================== 97 180 // unlucky if and only if the leading ideal is wrong 98 list LL=deleteUnluckyPrimes(T,L); 99 T=LL[1]; 100 L=LL[2]; 101 lT=LL[3]; 181 LL=deleteUnluckyPrimes(T,L); 182 TH=LL[1]; 183 L=LL[2]; 102 184 //============ now all leading ideals are the same ============ 103 for(j=1;j<=ncols(T[1]);j++) 104 { 105 for(k=1;k<=size(L);k++) 106 { 107 TT[k]=T[k][j]; 108 } 109 J[j]=liftPoly(TT,L); 110 } 111 //=========== chooses more primes up to the moment the result becomes stable 112 while(1) 113 { 114 k=0; 115 q=prime(random(an,en)); 116 while(k<size(L)) 185 for(j=1;j<=ncols(TH[1]);j++) 117 186 { 118 k++; 119 if(L[k]==q) 187 for(k=1;k<=size(L);k++) 188 { 189 TT[k]=TH[k][j]; 190 } 191 J[j]=liftPoly(TT,L); 192 } 193 if(pd>2){"list of primes:";L;"pTest";} 194 if(pTestSB(I,J,L)) 195 { 196 attrib(J,"isSB",1); 197 K=std(J); 198 if(size(reduce(I,J))==0) 120 199 { 121 k=0; 122 q=prime(random(an,en)); 123 } 124 } 125 L[size(L)+1]=q; 126 rl[1]=L[size(L)]; 127 def @r=ring(rl); 128 setring @r; 129 ideal i=fetch(R0,I); 130 option(redSB); 131 if(h) 132 { 133 i=std(i,hi); 134 } 135 else 136 { 137 if(ord_test(R0)==-1) 138 { 139 i=std(i); 140 } 141 else 142 { 143 matrix M; 144 i=liftstd(i,M); 145 } 146 } 147 setring R0; 148 T[size(T)+1]=fetch(@r,i); 149 kill @r; 150 cT=lead(T[size(T)]); 151 attrib(cT,"isSB",1); 152 if((size(reduce(cT,lT))!=0)||(size(reduce(lT,cT))!=0)) 153 { 154 T=delete(T,size(T)); 155 L=L[1..size(L)-1]; 156 k=0; 157 } 158 else 159 { 160 for(j=1;j<=ncols(T[1]);j++) 161 { 162 for(k=1;k<=size(L);k++) 200 if(size(reduce(K,J))==0) 163 201 { 164 TT[k]=T[k][j]; 202 if(!((h)||(ord_test(R0)==-1))) 203 { 204 205 "=================================================================="; 206 " The input is not homogeneous and the ordering is not 207 local. "; 208 "WARNING: ideal generated by output may be greater then 209 input ideal"; 210 211 "=================================================================="; 212 } 213 option(set, opt); 214 return(J); 165 215 } 166 K[j]=liftPoly(TT,L); 167 } 168 k=1; 169 for(j=1;j<=size(K);j++) 170 { 171 if(K[j]-J[j]!=0) 172 { 173 k=0; 174 J=K; 175 break; 176 } 177 } 178 } 179 if(k){break;} 180 } 181 //============ test for standard basis and I=J ======= 182 J=std(J); 183 I0=reduce(I0,J); 184 if(size(I0)>0) 185 { 186 "WARNING: The input ideal is not contained 187 in the ideal generated by the standardbasis"; 188 "list of primes used:"; 189 L; 190 } 191 attrib(J,"isSB",1); 192 return(J); 193 } 194 example 195 { "EXAMPLE:"; echo = 2; 196 ring r=0,(x,y,z),dp; 216 if(pd>2){"pTest o.k. but result wrong";} 217 } 218 if(pd>2){"pTest o.k. but result wrong";} 219 } 220 j=size(L)+1; 221 L=primeList(10,L); 222 } 223 } 224 example 225 { "EXAMPLE:"; echo = 2; 226 ring r=0,(x,y,z,t),dp; 197 227 ideal I=3x3+x2+1,11y5+y3+2,5z4+z2+4; 228 ideal J=modStd(I); 229 J; 230 I=homog(I,t); 231 J=modStd(I); 232 J; 233 ring s=0,(x,y,z),ds; 234 ideal I=jacob(x5+y6+z7+xyz); 198 235 ideal J=modStd(I); 199 236 J; … … 517 554 } 518 555 519 ///////////////////////////////////////////////////////////////////////////////520 proc primeList(int n)521 "USAGE: primeList(n);522 RETURN: the intvec of n greatest primes <= 2134567879523 EXAMPLE:example primList; shows an example524 "525 {526 intvec L=0:n;527 int i;528 int p=2134567879;529 for(i=1;i<=n;i++)530 {531 L[i]=p;532 p=prime(p-1);533 }534 return(L);535 }536 example537 { "EXAMPLE:"; echo = 2;538 intvec L=primeList(10);539 size(L);540 L[size(L)];541 }542 556 /////////////////////////////////////////////////////////////////////////////// 543 557 proc pStd(int p,ideal i)
Note: See TracChangeset
for help on using the changeset viewer.