Changeset 09830b6 in git
 Timestamp:
 Jul 25, 2013, 12:17:25 PM (10 years ago)
 Branches:
 (u'jengelhdatetime', 'ceac47cbc86fe4a15902392bdbb9bd2ae0ea02c6')(u'spielwiese', 'a800fe4b3e9d37a38c5a10cc0ae9dfa0c15a4ee6')
 Children:
 2ea18609e14482a4ddbe0fd0d68b8d5885268a0c
 Parents:
 0ab2b929de2eb4a16c8561981bbcef044441f401
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

Singular/LIB/primdecint.lib
r0ab2b92 r09830b6 17 17 18 18 PROCEDURES: 19 primdecZ(I); compute the primary decomposition of I 19 primdecZ(I); compute the primary decomposition of ideal I 20 primdecZM(I); compute the primary decomposition of module I 20 21 minAssZ(I); compute the minimal associated primes of I 21 22 radicalZ(I); compute the radical of I … … 224 225 J=simplify(J,2); 225 226 attrib(J,"isSB",1); 226 IS= maxIndependSet(J);227 IS=Primdec::maxIndependSet(J); 227 228 setring R; 228 229 //=== computing the pseudo primary and extract it … … 809 810 J=simplify(J,2); 810 811 attrib(J,"isSB",1); 811 IS= maxIndependSet(J);812 IS=Primdec::maxIndependSet(J); 812 813 setring R; 813 814 //=== computing the pseudo primary and extract it … … 1099 1100 //=== this is our way to obtain the coefficients in Z[u] of the 1100 1101 //=== leading terms of the Groebner basis above 1101 string quotring= prepareQuotientring(nvars(basering)L[1][3]);1102 string quotring=Primdec::prepareQuotientring(nvars(basering)L[1][3]); 1102 1103 execute(quotring); 1103 1104 ideal I=imap(Shelp,I); … … 1184 1185 } 1185 1186 return(J); 1186 }1187 1188 ////////////////////////////////////////////////////////////////////////////////1189 1190 static proc prepareQuotientring (int nnp)1191 {1192 //=== this is from primdec.lib, it is static there, should be imported later1193 //=== if it is no more static1194 ideal @ih,@jh;1195 int npar=npars(basering);1196 int @n;1197 1198 string quotring= "ring quring = ("+charstr(basering);1199 for(@n=nnp+1;@n<=nvars(basering);@n++)1200 {1201 quotring=quotring+",var("+string(@n)+")";1202 @ih=@ih+var(@n);1203 }1204 1205 quotring=quotring+"),(var(1)";1206 @jh=@jh+var(1);1207 for(@n=2;@n<=nnp;@n++)1208 {1209 quotring=quotring+",var("+string(@n)+")";1210 @jh=@jh+var(@n);1211 }1212 quotring=quotring+"),(C,lp);";1213 1214 return(quotring);1215 }1216 1217 ////////////////////////////////////////////////////////////////////////////////1218 1219 static proc maxIndependSet (ideal j)1220 {1221 //=== this is from primdec.lib, it is static there, should be imported later1222 //=== if it is no more static1223 int n,k,di;1224 1225 list resu,hilf;1226 if(size(j)==0)1227 {1228 resu[1]=varstr(basering);1229 resu[2]=ordstr(basering);1230 resu[3]=0;1231 return(list(resu));1232 }1233 string var1,var2;1234 list v=indepSet(j,0);1235 1236 for(n=1;n<=size(v);n++)1237 {1238 di=0;1239 var1="";1240 var2="";1241 for(k=1;k<=size(v[n]);k++)1242 {1243 if(v[n][k]!=0)1244 {1245 di++;1246 var2=var2+"var("+string(k)+"),";1247 }1248 else1249 {1250 var1=var1+"var("+string(k)+"),";1251 }1252 }1253 if(di>0)1254 {1255 var1=var1+var2;1256 var1=var1[1..size(var1)1];1257 hilf[1]=var1;1258 hilf[2]="lp";1259 hilf[3]=di;1260 resu[n]=hilf;1261 }1262 else1263 {1264 resu[n]=varstr(basering),ordstr(basering),0;1265 }1266 }1267 return(resu);1268 1187 } 1269 1188 … … 1402 1321 return(0); 1403 1322 } 1323 1324 //////////////////////////////////////////////////////////////////////////////// 1325 1326 static proc pseudo_primdecZM(module N) 1327 { 1328 ideal I=quotient(N,freemodule(nrows(N))); 1329 if(size(I)==0){return(list(list(N,I)));} 1330 1331 list B=minAssZ(I); 1332 list S,R,L; 1333 ideal K; 1334 if(size(B)==0){return(S);} 1335 for(int i=1;i<=size(B);i++) 1336 { 1337 S[i]=separatorsZ(i,B); 1338 } 1339 for(i=1;i<=size(B);i++) 1340 { 1341 L=sat(N,S[i]); 1342 K[i]=S[i]^L[2]; 1343 R[i]=list(L[1],B[i]); 1344 } 1345 L=pseudo_primdecZM(N+K*freemodule(nrows(N))); 1346 for(i=1;i<=size(L);i++) 1347 { 1348 R[size(R)+1]=L[i]; 1349 } 1350 return(R); 1351 } 1352 1353 1354 1355 static proc prepare_extractZM(list L) 1356 { 1357 def R=basering; 1358 module N=L[1]; 1359 ideal I=quotient(N,freemodule(nrows(N))); 1360 list B=primdecZ(I); 1361 list M; 1362 if(size(B)==1){return(M);} 1363 I=std(I); 1364 list rl=ringlist(R); 1365 if(deg(I[1])==0) 1366 { 1367 execute("int p="+string(I[1])+";"); 1368 rl[1]=p; 1369 } 1370 else 1371 { 1372 rl[1]=0; 1373 } 1374 def Shelp =ring(rl); 1375 setring Shelp; 1376 ideal I=imap(R,I); 1377 I=std(I); 1378 M=Primdec::maxIndependSet(I); 1379 setring R; 1380 return(M); 1381 } 1382 1383 1384 static proc extractZM(list M, list L) 1385 { 1386 //=== M is a list of a pseudo primary module and the corresponding prime 1387 //=== L is a list of maximal independent sets for P 1388 def R=basering; 1389 ideal P=M[2]; 1390 module I=M[1]; 1391 poly h=1; 1392 1393 //=== size(L)=0 means P is maximal ideal and I is primary 1394 if(size(L)>0) 1395 { 1396 if(L[1][3]!=0) 1397 { 1398 //=== if u in x is an independent set of L then we compute a Groebner 1399 //=== Basis in Z[u][xu] 1400 execute("ring S=integer,("+L[1][1]+"),lp;"); 1401 module I=imap(R,I); 1402 I=std(I); 1403 list rl=ringlist(S); 1404 rl[1]=0; 1405 def Shelp =ring(rl); 1406 setring Shelp; 1407 module I=imap(S,I); 1408 //=== this is our way to obtain the coefficients in Z[u] of the 1409 //=== leading terms of the Groebner basis above 1410 string quotring=Primdec::prepareQuotientring(nvars(basering)L[1][3]); 1411 execute(quotring); 1412 module I=imap(Shelp,I); 1413 list C; 1414 int i; 1415 for(i=1;i<=size(I);i++) 1416 { 1417 C[i]=leadcoef(I[i]); 1418 } 1419 setring Shelp; 1420 list C=imap(quring,C); 1421 setring R; 1422 list C=imap(Shelp,C); 1423 } 1424 else 1425 { 1426 // this is the case that P=<p>, p prime 1427 I=std(I); 1428 ideal IC=simplify(flatten(lead(I)),2); 1429 list C; 1430 int i; 1431 for(i=1;i<=size(IC);i++) 1432 { 1433 C[i]=I[i]; 1434 } 1435 list rl=ringlist(R); 1436 rl[1]=0; 1437 def Shelp =ring(rl); 1438 } 1439 for(i=1;i<=size(C);i++) 1440 { 1441 if(deg(C[i])>0){h=h*C[i];} // das muss noch besser gemacht werden, 1442 // nicht ausmultiplizieren! 1443 } 1444 setring Shelp; 1445 poly h=imap(R,h); 1446 ideal fac=factorize(h,1); 1447 setring R; 1448 list II; 1449 h=1; 1450 ideal fac=imap(Shelp,fac); 1451 for(i=1;i<=size(fac);i++) 1452 { 1453 II=sat(I,fac[i]); 1454 I=II[1]; 1455 h=h*fac[i]^II[2]; 1456 } 1457 } 1458 I=std(I); 1459 return(list(I,h)); 1460 } 1461 1462 1463 proc primdecZM(module N) 1464 "USAGE: primdecZM(N); N module 1465 RETURN: a list pr of primary modules and their associated primes: 1466 @format 1467 pr[i][1] the ith primary component, 1468 pr[i][2] the ith prime component. 1469 @end format 1470 EXAMPLE: example primdecZM; shows an example 1471 " 1472 { 1473 list P,K,S; 1474 int i,j; 1475 list L=pseudo_primdecZM(N); 1476 list M,O; 1477 for(i=1;i<=size(L);i++) 1478 { 1479 if(size(L[i][2])!=0) 1480 { 1481 M=prepare_extractZM(L[i]); 1482 O=extractZM(L[i],M); 1483 P[size(P)+1]=list(O[1],L[i][2]); 1484 K[size(K)+1]=L[i][1]+O[2]*freemodule(nrows(L[i][1])); 1485 } 1486 else 1487 { 1488 P[size(P)+1]=L[i]; 1489 } 1490 } 1491 for(j=1;j<=size(K);j++) 1492 { 1493 S=primdecZM(K[j]); 1494 for(i=1;i<=size(S);i++) 1495 { 1496 P[size(P)+1]=S[i]; 1497 } 1498 } 1499 return(P); 1500 } 1501 example 1502 { "EXAMPLE:"; echo = 2; 1503 ring R=integer,(x,y),(c,lp); 1504 module N=[0,0,xy2x2xy],[0,y,x],[0,x,2xyx],[x,0,xy],[0,0,18x]; 1505 primdecZM(N); 1506 } 1507 1404 1508 1405 1509 //////////////////////////////////////////////////////////////////////////////// … … 1602 1706 p^3*P^21, 1603 1707 F^2*f^31; 1708 1709 ring R=integer,(x,y),(c,lp); 1710 module N=[0,0,xy2x2xy],[0,y,x],[0,x,2xyx],[x,0,xy],[0,0,18x]; 1711 1712 ring R=integer,(x,y),(c,lp); 1713 module N=[0,0,xy2x2xy],[0,y,x],[0,x,2xyx],[x,0,xy],[0,0,18]; 1714 1715 ring R=integer,(x,y),(c,lp); 1716 module N=[y,7,0],[2y3y2],[3x,y2],[2yy2,x],[4,5x3]; 1717 1718 ring r=integer,(x,y),(c,lp); 1719 module N=[0,0,xy2x2xy],[0,y,x],[0,x,xyx],[x,0,xy],[5x,0,0]; 1720 1721 ring R2=integer,(a(1),a(2),a(3),b(1),b(2),b(3)),(c,lp); 1722 module N=[a(1)*b(1),a(2)*b(1),a(3)*b(1)],[a(1)*b(2),a(2)*b(2),a(3)*b(2)],[a(1)*b(3),a(2)*b(3),a(3)*b(3)]; 1723 1724 ring R3=integer,(x,y,z),(c,lp); 1725 module N=[y2+z2,xy,xz],[xy,x2+z2,yz],[xz,yz,x2+y2]; 1726 1727 ring R4=integer,(x,y,z,a,b,c),(c,lp); 1728 module N=[x3y2z2c,x2y3z2c,x2y2z3c],[x3y2z2b,x2y3z2b,x2y2z3b],[x3y2z2a,x2y3z2a,x2y2z3a]; 1729 1604 1730 */
Note: See TracChangeset
for help on using the changeset viewer.