Changeset f0daaa2 in git
 Timestamp:
 May 15, 2006, 12:59:17 PM (17 years ago)
 Branches:
 (u'spielwiese', '0d6b7fcd9813a1ca1ed4220cfa2b104b97a0a003')
 Children:
 45c9af467659b9010e2fc67a7cb8da665abaf595
 Parents:
 181148eb58afcf5d585ec06f89fbab21c1ab9c99
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

Singular/LIB/primdec.lib
r181148 rf0daaa2 1 1 /////////////////////////////////////////////////////////////////////////////// 2 version="$Id: primdec.lib,v 1.11 4 20060512 12:12:05Singular Exp $";2 version="$Id: primdec.lib,v 1.115 20060515 10:59:17 Singular Exp $"; 3 3 category="Commutative Algebra"; 4 4 info=" … … 5509 5509 5510 5510 /////////////////////////////////////////////////////////////////////////////// 5511 proc radical(ideal i,list #) 5512 "USAGE: radical(i); i ideal. 5511 proc radical(ideal i, list #) 5512 "USAGE: radical(i); i ideal. 5513 Optional parameters in list #: (can be entered in any order) 5514 1, "equiRad" > equiRadical is computed 5515 0, "fullRad" > full radical is computed (default) 5516 "SL" > the new algorithm is used (default) 5517 "KL" > the old algorithm is used 5518 "KLdp" > the old algorithm, with modifications 5519 "facstd" > uses facstd to first decompose the ideal (default for non homogeneous ideals) 5520 "noFacstd" > does not use facstd (default for homogeneous ideals) 5513 5521 RETURN: ideal, the radical of i. 5514 5522 NOTE: A combination of the algorithms of Krick/Logar and Kemper is used. … … 5517 5525 " 5518 5526 { 5527 dbprint(printlevel  voice, "Radical, version 2006.05.08"); 5519 5528 if(ord_test(basering)!=1) 5520 5529 { … … 5524 5533 } 5525 5534 def @P=basering; 5526 int j,il; 5527 if(size(#)>0){il=#[1];} 5528 if(size(i)==0){return(ideal(0));} 5529 ideal re=1; 5535 int j; 5536 int il; 5537 string algorithm; 5538 int useFac; 5539 5540 // Set input parameters 5541 algorithm = "SL"; // Default: SL algorithm 5542 il = 0; // Default: Full radical (not only equidim part) 5543 if (homog(i) == 1) { // Default: facStd is used, except if the ideal is homogeneous. 5544 useFac = 0; 5545 } else { 5546 useFac = 1; 5547 }; 5548 if(size(#) > 0){ 5549 int valid; 5550 for(j = 1; j <= size(#); j++){ 5551 valid = 0; 5552 if((typeof(#[j]) == "int") or (typeof(#[j]) == "number")) { 5553 il = #[j]; // If il == 1, equiRadical is computed 5554 valid = 1; 5555 }; 5556 if(typeof(#[j]) == "string"){ 5557 if(#[j] == "KL") { 5558 algorithm = "KL"; 5559 valid = 1}; 5560 if(#[j] == "KLdp") { 5561 algorithm = "KLdp"; 5562 valid = 1}; 5563 if(#[j] == "SL") { 5564 algorithm = "SL"; 5565 valid = 1}; 5566 if(#[j] == "noFacstd") { 5567 useFac = 0; 5568 valid = 1}; 5569 if(#[j] == "facstd") { 5570 useFac = 1; 5571 valid = 1}; 5572 if(#[j] == "equiRad") { 5573 il = 1; 5574 valid = 1}; 5575 if(#[j] == "fullRad") { 5576 il = 0; 5577 valid = 1}; 5578 }; 5579 if(valid == 0) { 5580 dbprint(1, "Warning! The following input parameter was not recognized:", #[j]); 5581 }; 5582 }; 5583 }; 5584 5585 if(size(i) == 0){return(ideal(0));} 5586 ideal rad = 1; 5530 5587 intvec op = option(get); 5531 list qr=simplifyIdeal(i); 5532 ideal isave=i; 5533 map phi=@P,qr[2]; 5588 list qr = simplifyIdeal(i); 5589 // SL  Line removed. The ideal isave was never used. 5590 //ideal isave=i; 5591 map phi = @P, qr[2]; 5534 5592 5535 5593 option(redSB); 5536 i =groebner(qr[1]);5537 option(set, op);5538 int di =dim(i);5539 5540 if(di ==0)5541 { 5542 i =zeroRad(i,qr[1]);5594 i = groebner(qr[1]); 5595 option(set, op); 5596 int di = dim(i); 5597 5598 if(di == 0) 5599 { 5600 i = zeroRad(i, qr[1]); 5543 5601 return(interred(phi(i))); 5544 5602 } 5545 5603 5546 5604 option(redSB); 5547 list pr=i; 5548 if (!homog(i)) 5549 { 5550 pr=facstd(i); 5551 } 5552 option(set,op); 5553 int s=size(pr); 5554 5555 for(j=1;j<=s;j++) 5556 { 5557 attrib(pr[s+1j],"isSB",1); 5558 if((size(reduce(re,pr[s+1j],1))!=0)&&((dim(pr[s+1j])==di)!il)) 5559 { 5560 re=intersect(re,radicalKL(pr[s+1j],re,il)); 5561 } 5562 } 5563 return(interred(phi(re))); 5605 list pr; 5606 if(useFac == 1) 5607 { 5608 pr = facstd(i); 5609 } else { 5610 pr = i 5611 }; 5612 option(set, op); 5613 int s = size(pr); 5614 if(useFac == 1) { 5615 dbprint(printlevel  voice, "Number of components returned by facstd: ", s); 5616 }; 5617 for(j = 1; j <= s; j++) 5618 { 5619 attrib(pr[s + 1  j], "isSB", 1); 5620 if((size(reduce(rad, pr[s + 1  j], 1)) != 0) && ((dim(pr[s + 1  j]) == di)  !il)) 5621 { 5622 // SL Debug messages 5623 dbprint(printlevelvoice, "We shall compute the radical of ", pr[s + 1  j]); 5624 dbprint(printlevelvoice, "The dimension is: ", dim(pr[s+1j])); 5625 5626 if(algorithm == "KL") { 5627 rad = intersect(rad, radicalKL(pr[s + 1  j], rad, il)); 5628 } 5629 if(algorithm == "KLdp") { 5630 rad = intersect(rad, newRadicalKL(pr[s + 1  j], rad, il)); 5631 } 5632 if(algorithm == "SL") { 5633 rad = intersect(rad, newRadicalSL(pr[s + 1  j], il)); 5634 }; 5635 } else 5636 { 5637 // SL Debug 5638 dbprint(printlevelvoice, "The radical of this component is not needed."); 5639 dbprint(printlevelvoice, "size(reduce(rad, pr[s + 1  j], 1))", size(reduce(rad, pr[s + 1  j], 1))); 5640 dbprint(printlevelvoice, "dim(pr[s + 1  j])", dim(pr[s + 1  j])); 5641 dbprint(printlevelvoice, "il", il); 5642 }; 5643 } 5644 return(interred(phi(rad))); 5564 5645 } 5565 5646 example … … 5569 5650 poly q = z3+2; 5570 5651 ideal i = p*q^2,yz2; 5571 ideal pr = radical(i);5652 ideal pr = radical(i); 5572 5653 pr; 5573 5654 } 5655 5656 /////////////////////////////////////////////////////////////////////////////// 5657 // 5658 // Computes the radical of I using KL algorithm. 5659 // The only difference with the the previous implementation of KL algorithm is 5660 // that now it uses block dp instead of lp ordering for the reduction to the 5661 // zerodimensional case. 5662 // The reduction step has been moved to the new routine newRadicalReduction, so that it can be 5663 // used also by the new algorithm. 5664 // This algorithm should be left in the final version of the library, 5665 // and radicalKL removed. 5666 // 5667 static proc newRadicalKL(ideal I, ideal ser, list #) { 5668 // I The ideal for which the radical is computed 5669 // ideal ser Same as in radicalKL (used to reduce components already obtained) 5670 // list # If #[1] = 1, equiradical is computed. 5671 5672 // I needs to be a Groebner basis. 5673 if (attrib(I, "isSB") != 1) { 5674 I = groebner(I); 5675 }; 5676 5677 ideal rad; // The radical 5678 int allIndep = 1; // All max independent sets are used 5679 5680 list result = newRadicalReduction(I, ser, allIndep, #); 5681 int done = result[3]; 5682 rad = result[1]; 5683 if (done == 0) { 5684 rad = intersect(rad, newRadicalKL(result[2], ideal(1), #)); 5685 }; 5686 return(rad); 5687 }; 5688 5689 5690 /////////////////////////////////////////////////////////////////////////////// 5691 // 5692 // Computes the radical of I via the new algorithm, using zerodimensional radical in 5693 // the zero dimensional case. 5694 // For the reduction to the zerodimensional case, it uses the procedure 5695 // radical of primdec.lib, with some modifications to avoid the recursion. 5696 // 5697 static proc newRadicalSL(ideal I, list #) 5698 // Input = I, ideal 5699 // #, list. If #[1] = 1, then computes only the equiradical. 5700 // Output = (P, primaryDec) where P = rad(I) and primaryDec is the list of the radicals 5701 // obtained in intermediate steps. 5702 { 5703 ideal rad = 1; 5704 ideal equiRad = 1; 5705 list primes; 5706 int k; // Counter 5707 int il; // If il = 1, only the equiradical is required. 5708 int iDim; // The dimension of I 5709 int stop = 0; // Checks if the radical has been obtained 5710 5711 if (attrib(I, "isSB") != 1) { 5712 I = groebner(I); 5713 }; 5714 iDim = dim(I); 5715 5716 // Checks if only equiradical is required 5717 if (size(#) > 0) { 5718 il = #[1]; 5719 }; 5720 5721 while(stop == 0) { 5722 dbprint (printlevelvoice, "// We call radLoopR to find new prime ideals."); 5723 primes = newRadicalSLIteration(I, rad); // A list of primes or intersections of primes, not included in P 5724 dbprint (printlevel  voice, "// Output of Iteration Step:"); 5725 dbprint (printlevel  voice, primes); 5726 if (size(primes) > 0) { 5727 dbprint (printlevel  voice, "// We intersect P with the ideal just obtained."); 5728 for(k = 1; k <= size(primes); k++) { 5729 rad = intersect(rad, primes[k]); 5730 if (il == 1){ 5731 if (attrib(primes[k], "isSB") != 1) { 5732 primes[k] = groebner(primes[k]); 5733 }; 5734 if (iDim == dim(primes[k])) { 5735 equiRad = intersect(equiRad, primes[k]); 5736 }; 5737 }; 5738 }; 5739 } else { 5740 stop = 1 5741 }; 5742 }; 5743 if (il == 0) { 5744 return(rad); 5745 } else { 5746 return(equiRad); 5747 }; 5748 }; 5749 5750 ////////////////////////////////////////////////////////////////////////// 5751 // Based on radicalKL. 5752 // It contains all of proc radicalKL except the recursion call. 5753 5754 // Output: 5755 // #1 > output ideal, the part of the radical that has been computed 5756 // #2 > complementary ideal, the part of the ideal I whose radical remains to be computed 5757 // = (I, h) in KL algorithm 5758 // This is not used in the new algorithm. It is part of KL algorithm 5759 // #3 > done, 1: output = radical, there is no need to continue 5760 // 0: radical = output \cap \sqrt{complementary ideal} 5761 // This is not used in the new algorithm. It is part of KL algorithm 5762 5763 static proc newRadicalReduction(ideal I, ideal ser, int allIndep, list #) { 5764 // allMaximal 1 > Indicates that the reduction to the zerodim case 5765 // must be done for all indep set of the leading terms ideal 5766 // 0 > Otherwise 5767 // ideal ser Only for radicalKL. (Same as in radicalKL) 5768 // list # Only for radicalKL (If #[1] = 1, only equiradical is required. 5769 // It is used to set the value of done.) 5770 5771 attrib(I, "isSB", 1); // I needs to be a reduced standard basis 5772 list indep, fett; 5773 intvec @w, @hilb, op; 5774 int @wr, @n, @m, lauf, di; 5775 ideal fac, @h, collectrad, lsau; 5776 poly @q; 5777 string @va, quotring; 5778 5779 def @P = basering; 5780 int jdim = dim(I); // Computes the dimension of I 5781 int homo = homog(I); // Finds out if I is homogeneous 5782 ideal rad = ideal(1); // The unit ideal 5783 ideal te = ser; 5784 if(size(#) > 0) 5785 { 5786 @wr = #[1]; 5787 } 5788 if(homo == 1) 5789 { 5790 for(@n = 1; @n <= nvars(basering); @n++) 5791 { 5792 @w[@n] = ord(var(@n)); 5793 } 5794 @hilb = hilb(I, 1, @w); 5795 } 5796 5797 // SL 2006.04.11 1 Debug messages 5798 dbprint(printlevelvoice, "//Computes the radical of the ideal:", I); 5799 // SL 2006.04.11 2 Debug messages 5800 5801 // 5802 //j is the ring 5803 // 5804 5805 if (jdim==1) 5806 { 5807 return(ideal(1), ideal(1), 1); 5808 } 5809 5810 // 5811 //the zerodimensional case 5812 // 5813 5814 if (jdim==0) 5815 { 5816 return(zeroRad(I), ideal(1), 1); 5817 } 5818 5819 // 5820 //search for a maximal independent set indep,i.e. 5821 //look for subring such that the intersection with the ideal is zero 5822 //j intersected with K[var(indep[3]+1),...,var(nvar)] is zero, 5823 //indep[1] is the new varstring, indep[2] the string for the blockordering 5824 // 5825 5826 // SL 20060424 1 If allIndep = 0, then it only computes one maximal 5827 // independent set. 5828 // This looks better for the new algorithm but not for KL 5829 // algorithm 5830 list parameters = allIndep; 5831 indep = newMaxIndependSet(I, parameters); 5832 // SL 20060424 2 5833 5834 for(@m = 1; @m <= size(indep); @m++) 5835 { 5836 if((indep[@m][1] == varstr(basering)) && (@m == 1)) 5837 //this is the good case, nothing to do, just to have the same notations 5838 //change the ring 5839 { 5840 execute("ring gnir1 = ("+charstr(basering)+"),("+varstr(basering)+"),(" 5841 +ordstr(basering)+");"); 5842 ideal @j = fetch(@P, I); 5843 attrib(@j, "isSB", 1); 5844 } 5845 else 5846 { 5847 @va = string(maxideal(1)); 5848 5849 execute("ring gnir1 = (" + charstr(basering) + "), (" + indep[@m][1] + "),(" 5850 + indep[@m][2] + ");"); 5851 execute("map phi = @P," + @va + ";"); 5852 if(homo == 1) 5853 { 5854 ideal @j = std(phi(I), @hilb, @w); 5855 } 5856 else 5857 { 5858 ideal @j = groebner(phi(I)); 5859 } 5860 } 5861 if((deg(@j[1]) == 0)  (dim(@j) < jdim)) 5862 { 5863 setring @P; 5864 break; 5865 } 5866 for (lauf = 1; lauf <= size(@j); lauf++) 5867 { 5868 fett[lauf] = size(@j[lauf]); 5869 } 5870 // 5871 // We have now the following situation: 5872 // j intersected with K[var(nnp+1),..,var(nva)] is zero so we may pass 5873 // to this quotientring, j is there still a standardbasis, the 5874 // leading coefficients of the polynomials there (polynomials in 5875 // K[var(nnp+1),..,var(nva)]) are collected in the list h, 5876 // we need their LCM, gh, because of the following: 5877 // let (j:gh^n)=(j:gh^infinity) then j*K(var(nnp+1),..,var(nva))[..rest..] 5878 // intersected with K[var(1),...,var(nva)] is (j:gh^n) 5879 // on the other hand j = ((j, gh^n) intersected with (j : gh^n)) 5880 5881 // 5882 // The arrangement for the quotientring K(var(nnp+1),..,var(nva))[..rest..] 5883 // and the map phi:K[var(1),...,var(nva)] > 5884 // K(var(nnpr+1),..,var(nva))[..the rest..] 5885 // 5886 quotring = newPrepareQuotientRing(nvars(basering)  indep[@m][3]); 5887 // 5888 // We pass to the quotientring K(var(nnp+1),..,var(nva))[..the rest..] 5889 // 5890 5891 execute(quotring); 5892 5893 // @j considered in the quotientring 5894 ideal @j = imap(gnir1, @j); 5895 5896 kill gnir1; 5897 5898 // j is a standardbasis in the quotientring but usually not minimal 5899 // here it becomes minimal 5900 5901 @j = clearSB(@j, fett); 5902 5903 // We need later LCM(h[1],...) = gh for saturation 5904 ideal @h; 5905 if(deg(@j[1]) > 0) 5906 { 5907 for(@n = 1; @n <= size(@j); @n++) 5908 { 5909 @h[@n] = leadcoef(@j[@n]); 5910 } 5911 op = option(get); 5912 option(redSB); 5913 @j = interred(@j); //to obtain a reduced standardbasis 5914 attrib(@j, "isSB", 1); 5915 option(set, op); 5916 5917 // SL 1 Debug messages 5918 dbprint(printlevel  voice, "zero_rad", basering, @j, dim(groebner(@j))); 5919 ideal zero_rad = zeroRad(@j); 5920 dbprint(printlevel  voice, "zero_rad passed"); 5921 // SL 2 5922 } 5923 else 5924 { 5925 ideal zero_rad = ideal(1); 5926 } 5927 5928 // We need the intersection of the ideals in the list quprimary with the 5929 // polynomialring, i.e. let q=(f1,...,fr) in the quotientring such an ideal 5930 // but fi polynomials, then the intersection of q with the polynomialring 5931 // is the saturation of the ideal generated by f1,...,fr with respect to 5932 // h which is the lcm of the leading coefficients of the fi considered in 5933 // the quotientring: this is coded in saturn 5934 5935 zero_rad = std(zero_rad); 5936 5937 ideal hpl; 5938 5939 for(@n = 1; @n <= size(zero_rad); @n++) 5940 { 5941 hpl = hpl, leadcoef(zero_rad[@n]); 5942 } 5943 5944 // 5945 // We leave the quotientring K(var(nnp+1),..,var(nva))[..the rest..] 5946 // back to the polynomialring 5947 // 5948 setring @P; 5949 5950 collectrad = imap(quring, zero_rad); 5951 lsau = simplify(imap(quring, hpl), 2); 5952 @h = imap(quring, @h); 5953 5954 kill quring; 5955 5956 // Here the intersection with the polynomialring 5957 // mentioned above is really computed 5958 5959 collectrad = sat2(collectrad, lsau)[1]; 5960 if(deg(@h[1])>=0) 5961 { 5962 fac = ideal(0); 5963 for(lauf = 1; lauf <= ncols(@h); lauf++) 5964 { 5965 if(deg(@h[lauf]) > 0) 5966 { 5967 fac = fac + factorize(@h[lauf], 1); 5968 } 5969 } 5970 fac = simplify(fac, 4); 5971 @q = 1; 5972 for(lauf = 1; lauf <= size(fac); lauf++) 5973 { 5974 @q = @q * fac[lauf]; 5975 } 5976 op = option(get); 5977 option(returnSB); 5978 option(redSB); 5979 I = quotient(I + ideal(@q), rad); 5980 attrib(I, "isSB", 1); 5981 option(set, op); 5982 } 5983 if((deg(rad[1]) > 0) && (deg(collectrad[1]) > 0)) 5984 { 5985 rad = intersect(rad, collectrad); 5986 te = intersect(te, collectrad); 5987 te = simplify(reduce(te, I, 1), 2); 5988 } 5989 else 5990 { 5991 if(deg(collectrad[1]) > 0) 5992 { 5993 rad = collectrad; 5994 te = intersect(te, collectrad); 5995 te = simplify(reduce(te, I, 1), 2); 5996 } 5997 } 5998 5999 if((dim(I) < jdim)(size(te) == 0)) 6000 { 6001 break; 6002 } 6003 if(homo==1) 6004 { 6005 @hilb = hilb(I, 1, @w); 6006 } 6007 } 6008 6009 // SL 2006.04.11 1 Debug messages 6010 dbprint (printlevelvoice, "// Part of the Radical already computed:", rad); 6011 dbprint (printlevelvoice, "// Dimension:", dim(groebner(rad))); 6012 // SL 2006.04.11 2 Debug messages 6013 6014 // SL 2006.04.21 1 New variable "done". 6015 // It tells if the radical is already computed or 6016 // if it still has to be computed the radical of the new ideal I 6017 int done; 6018 if(((@wr == 1) && (dim(I)<jdim))  (deg(I[1])==0)  (size(te) == 0)) 6019 { 6020 done = 1; 6021 } else { 6022 done = 0; 6023 }; 6024 // SL 2006.04.21 2 6025 6026 // SL 2006.04.21 1 See details of the output at the beggining of this proc. 6027 list result = rad, I, done; 6028 return(result); 6029 // SL 2006.04.21 2 6030 }; 6031 6032 6033 6034 /////////////////////////////////////////////////////////////////////////////// 6035 // Given an ideal I and an ideal P (intersection of some minimal prime ideals 6036 // associated to I), it calculates the intersection of new minimal prime ideals 6037 // associated to I which where not used to calculate P. 6038 // This version uses ZD Radical in the zerodimensional case. 6039 static proc newRadicalSLIteration (ideal I, ideal P); 6040 // Input: I, ideal. The ideal from which new prime components will be obtained. 6041 // P, ideal. Intersection of some prime ideals of I. 6042 // Output: ideal. Intersection of some primes of I different from the ones in P. 6043 { 6044 int k = 1; // Counter 6045 int good = 0; // Checks if an element of P is in rad(I) 6046 6047 dbprint (printlevelvoice, "// We search for an element in P  sqrt(I)."); 6048 while ((k <= size(P)) and (good == 0)) { 6049 dbprint (printlevelvoice, "// We try with:", P[k]); 6050 good = 1  rad_con(P[k], I); 6051 k++; 6052 }; 6053 k; 6054 if (good == 0) { 6055 dbprint (printlevelvoice, "// No element was found, P = sqrt(I)."); 6056 list emptyList = list(); 6057 return (emptyList); 6058 }; 6059 dbprint(printlevel  voice, "// That one was good!"); 6060 dbprint(printlevel  voice, "// We saturate I with respect to this element."); 6061 if (P[k] != 1) { 6062 ideal J = sat(I, P[k])[1]; 6063 } else { 6064 dbprint(printlevel  voice, "// The polynomial is 1, the saturation in not actually computed."); 6065 ideal J = I; 6066 }; 6067 6068 // We now call proc radicalNew; 6069 dbprint(printlevel  voice, "// We do the reduction to the zerodimensional case, via radical."); 6070 dbprint(printlevel  voice, "// The ideal is ", J); 6071 dbprint(printlevel  voice, "// The dimension is ", dim(groebner(J))); 6072 6073 int allMaximal = 0; // Compute the zerodim reduction for only one indep set. 6074 ideal re = 1; // No reduction is need, there are not redundant components. 6075 list emptyList = list(); // Look for primes of any dimension, not only of max dimension. 6076 list result = newRadicalReduction(J, re, allMaximal, emptyList); 6077 6078 return(result[1]); 6079 }; 6080 6081 6082 6083 /////////////////////////////////////////////////////////////////////////////////// 6084 // Based on maxIndependSet 6085 // Added list # as parameter 6086 // If the first element of # is 0, the output is only 1 max indep set. 6087 // If no list is specified or #[1] = 1, the output is all the max indep set of the 6088 // leading terms ideal. This is the original output of maxIndependSet 6089 6090 // The ordering given in the output has been changed to block dp instead of lp. 6091 6092 static proc newMaxIndependSet(ideal j, list #) 6093 // #[1]=0 > computes only one max indep set 6094 "USAGE: maxIndependentSet(i); i ideal 6095 RETURN: list = #1. new varstring with the maximal independent set at the end, 6096 #2. ordstring with the corresponding block ordering, 6097 #3. the number of independent variables 6098 // SL #3 explanation modified, 6099 // it was: integer where the independent set starts in the varstring 6100 NOTE: 6101 EXAMPLE: example maxIndependentSet; shows an example 6102 " 6103 { 6104 int n, k, di; 6105 list resu, hilf; 6106 string var1, var2; 6107 list v = indepSet(j, 0); 6108 6109 // SL 2006.04.21 1 Lines modified to use only one independent Set 6110 int allMaximal; 6111 if (size(#) > 0) { 6112 allMaximal = #[1]; 6113 } else { 6114 allMaximal = 1; 6115 }; 6116 6117 int nMax; 6118 if (allMaximal == 1) { 6119 nMax = size(v); 6120 } else { 6121 nMax = 1; 6122 }; 6123 6124 for(n = 1; n <= nMax; n++) 6125 // SL 2006.04.21 2 6126 { 6127 di = 0; 6128 var1 = ""; 6129 var2 = ""; 6130 for(k = 1; k <= size(v[n]); k++) 6131 { 6132 if(v[n][k] != 0) 6133 { 6134 di++; 6135 var2 = var2 + "var(" + string(k) + "), "; 6136 } 6137 else 6138 { 6139 var1 = var1 + "var(" + string(k) + "), "; 6140 } 6141 } 6142 if(di > 0) 6143 { 6144 var1 = var1 + var2; 6145 var1 = var1[1..size(var1)  2]; // The " 2" removes the trailer comma 6146 hilf[1] = var1; 6147 // SL 2006.21.04 1 The order is now block dp instead of lp 6148 hilf[2] = "dp(" + string(nvars(basering)  di) + "), dp(" + string(di) + ")"; 6149 // SL 2006.21.04 2 6150 hilf[3] = di; 6151 resu[n] = hilf; 6152 } 6153 else 6154 { 6155 resu[n] = varstr(basering), ordstr(basering), 0; 6156 } 6157 } 6158 return(resu); 6159 } 6160 example 6161 { "EXAMPLE:"; echo = 2; 6162 ring s1 = (0, x, y), (a, b, c, d, e, f, g), lp; 6163 ideal i = ea  fbg, fa + be, ec  fdg, fc + de; 6164 i = std(i); 6165 list l = newMaxIndependSet(i); 6166 l; 6167 i = i, g; 6168 l = newMaxIndependSet(i); 6169 l; 6170 6171 ring s = 0, (x, y, z), lp; 6172 ideal i = z, yx; 6173 list l = newMaxIndependSet(i); 6174 l; 6175 } 6176 6177 6178 /////////////////////////////////////////////////////////////////////////////// 6179 // based on prepareQuotientring 6180 // The order returned is now (C, dp) instead of (C, lp) 6181 6182 static proc newPrepareQuotientRing (int nnp) 6183 "USAGE: newPrepareQuotientRing(nnp); nnp int 6184 RETURN: string = to define Kvar(nnp+1),...,var(nvars)[..rest ] 6185 NOTE: 6186 EXAMPLE: example newPrepareQuotientRing; shows an example 6187 " 6188 { 6189 ideal @ih,@jh; 6190 int npar=npars(basering); 6191 int @n; 6192 6193 string quotring= "ring quring = ("+charstr(basering); 6194 for(@n = nnp + 1; @n <= nvars(basering); @n++) 6195 { 6196 quotring = quotring + ", var(" + string(@n) + ")"; 6197 @ih = @ih + var(@n); 6198 } 6199 6200 quotring = quotring+"),(var(1)"; 6201 @jh = @jh + var(1); 6202 for(@n = 2; @n <= nnp; @n++) 6203 { 6204 quotring = quotring + ", var(" + string(@n) + ")"; 6205 @jh = @jh + var(@n); 6206 } 6207 // SL 20060421 1 The order returned is now (C, dp) instead of (C, lp) 6208 quotring = quotring + "), (C, dp);"; 6209 // SL 20060421 2 6210 6211 return(quotring); 6212 } 6213 example 6214 { "EXAMPLE:"; echo = 2; 6215 ring s1=(0,x),(a,b,c,d,e,f,g),lp; 6216 def @Q=basering; 6217 list l= newPrepareQuotientRing(3); 6218 l; 6219 execute(l[1]); 6220 execute(l[2]); 6221 basering; 6222 phi; 6223 setring @Q; 6224 6225 } 6226 5574 6227 /////////////////////////////////////////////////////////////////////////////// 5575 6228 proc prepareAss(ideal i)
Note: See TracChangeset
for help on using the changeset viewer.