Changeset 2abb041 in git
- Timestamp:
- Feb 24, 2007, 4:08:04 PM (16 years ago)
- Branches:
- (u'spielwiese', '0d6b7fcd9813a1ca1ed4220cfa2b104b97a0a003')
- Children:
- bd7468e27c17b2493beec0180f366c438528745d
- Parents:
- b863e9dd32456f975a36c9963618660981b3bdb9
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
Singular/LIB/finvar.lib
rb863e9 r2abb041 1 1 /////////////////////////////////////////////////////////////////////////////// 2 version="$Id: finvar.lib,v 1.6 4 2007-02-21 13:39:21 kingExp $"2 version="$Id: finvar.lib,v 1.65 2007-02-24 15:08:04 Singular Exp $" 3 3 category="Invariant theory"; 4 4 info=" … … 4769 4769 { ISSort[i-1] = ideal(B[j]); 4770 4770 } 4771 if (v) 4772 { " We found", irrcounter,"of", NrIS ,"irr. sec. inv. in degree ",i-1; 4771 if (v) 4772 { " We found", irrcounter,"of", NrIS ,"irr. sec. inv. in degree ",i-1; 4773 4773 } 4774 4774 Reductor = ideal(helpP); … … 4788 4788 else // use fast algorithm 4789 4789 { B=sort_of_invariant_basis(sP_Reductor,REY,i-1,int(dimmat[i,1])*6); 4790 // B contains 4790 // B contains 4791 4791 // images of kbase(sP,i-1) under the 4792 4792 // Reynolds operator that are linearly … … 4825 4825 { ISSort[i-1] = ideal(B[j]); 4826 4826 } 4827 if (v) 4828 { " We found", irrcounter,"of", NrIS,"irr. sec. inv. in degree ",i-1; 4827 if (v) 4828 { " We found", irrcounter,"of", NrIS,"irr. sec. inv. in degree ",i-1; 4829 4829 } 4830 4830 Reductor = ideal(helpP); … … 4900 4900 @* \"PP\": if this string occurs as (optional) parameter, then in 4901 4901 all degrees power products of irr. sec. inv. will be computed. 4902 RETURN: Irreducible homogeneous secondary invariants of the invariant ring 4902 RETURN: Irreducible homogeneous secondary invariants of the invariant ring 4903 4903 (type <matrix>) 4904 4904 ASSUME: We are in the non-modular case, i.e., the characteristic of the basering … … 4915 4915 independent modulo the primary invariants using Groebner basis techniques 4916 4916 (see paper \"Fast Computation of Secondary Invariants\" by S. King). 4917 The size of this set can be read off from the Molien series. Here, only 4918 irreducible secondary invariants are explicitly computed, which saves time and 4917 The size of this set can be read off from the Molien series. Here, only 4918 irreducible secondary invariants are explicitly computed, which saves time and 4919 4919 memory. 4920 4920 @* Moreover, if no irr. sec. inv. in degree d-1 have been found and unless the last 4921 optional paramter \"PP\" is used, a Groebner basis of primary invariants and 4922 irreducible secondary invariants up to degree d-2 is computed, which allows to 4923 detect irr. sec. inv. in degree d without computing power products. 4921 optional paramter \"PP\" is used, a Groebner basis of primary invariants and 4922 irreducible secondary invariants up to degree d-2 is computed, which allows to 4923 detect irr. sec. inv. in degree d without computing power products. 4924 4924 @* There are three internal parameters \"pieces\", \"MonStep\" and \"IrrSwitch\". 4925 4925 The default values of the parameters should be fine in most cases. However, … … 4937 4937 // If this parameter is small, the memory consumption will decrease. 4938 4938 int IrrSwitch = 7; // Up to degree #IrrSwitch#-1, or if there are few monomials in kbase(sP,i-1), 4939 // we use a method for generating irreducible secondaries that tends to 4940 // produce a smaller output (which is good for subsequent computations). 4941 // However, this method needs to much memory in high degrees, and so we 4939 // we use a method for generating irreducible secondaries that tends to 4940 // produce a smaller output (which is good for subsequent computations). 4941 // However, this method needs to much memory in high degrees, and so we 4942 4942 // use a sparser method from degree #IrrSwitch# upwards. 4943 4943 def br=basering; … … 5010 5010 // of a certain degree - 5011 5011 m=nrows(dimmat); // m-1 is the highest degree 5012 if (v) 5012 if (v) 5013 5013 { "There are"; 5014 5014 for (j=1;j<=m;j++) … … 5069 5069 { if ((LastNewIS == i-3) and (UsePP==0)) // We guess that all irreducibles have been found! 5070 5070 { if (v) {"Computing Groebner basis for primaries and previously found irreducibles...";} 5071 TotalReductor = groebner(sP+IS); 5072 TotalSet = 1; 5071 TotalReductor = groebner(sP+IS); 5072 TotalSet = 1; 5073 5073 } 5074 5074 if (int(dimmat[i,1])<>0) // when it is == 0 we need to find no … … 5156 5156 attrib(RedSSort[i-1][k],"size",saveAttr); 5157 5157 if (v) 5158 { 5158 { 5159 5159 " We found reducible sec. inv. number "+string(counter); 5160 5160 } … … 5187 5187 } 5188 5188 } // if TotalSet==0 5189 else 5189 else 5190 5190 { if (v) {" We don't compute Power Products!"; } 5191 5191 // Instead of computing Power Products, we can compute a Groebner basis of the ideal … … 5193 5193 // polynomial of degree i-1 belongs to this ideal if and only if it belongs to the sub-algebra 5194 5194 // generated by primary and irreducible secondary invariants up to degree i-2. 5195 // Hence, if we find an invariant outside this ideal, it is an irreducible secondary 5195 // Hence, if we find an invariant outside this ideal, it is an irreducible secondary 5196 5196 // invariant of degree i-1. 5197 5197 } … … 5231 5231 // Compare the comments on the computation of reducible sec. inv.! 5232 5232 while (j < ii) // we will break the loop if counter==int(dimmat[i,1]) 5233 { if ((j mod MonStep) == 0) 5233 { if ((j mod MonStep) == 0) 5234 5234 { if ((j+MonStep) <= ii) 5235 5235 { ideal tmp = normalize(evaluate_reynolds(REY,ideal(mon[j+1..j+MonStep]))); … … 5260 5260 kill tmp; 5261 5261 } 5262 } 5262 } 5263 5263 j++; 5264 5264 helpP = Indicator[j]; … … 5303 5303 { B=sort_of_invariant_basis(TotalReductor,REY,i-1,int(dimmat[i,1])*6); 5304 5304 } 5305 // B is a linearly independent set of reynolds 5305 // B is a linearly independent set of reynolds 5306 5306 // images of monomials that do not occur 5307 5307 // as leading monomials of the ideal spanned 5308 // by primary and previously found irreducible 5308 // by primary and previously found irreducible 5309 5309 // secondary invariants. 5310 5310 if (TotalSet==0 && ncols(B)+counter==int(dimmat[i,1])) // then we can take all of B … … 5323 5323 else 5324 5324 { irrcounter=0; 5325 j=0; // j goes through all of B 5325 j=0; // j goes through all of B 5326 5326 // Compare the comments on the computation of reducible sec. inv.! 5327 5327 ReducedCandidates = reduce(B,sP); … … 5960 5960 { ISSort[i-1] = ideal(B[j]); 5961 5961 } 5962 if (v) 5963 { " We found", irrcounter,"of", NrIS ,"irr. sec. inv. in degree ",i-1; 5962 if (v) 5963 { " We found", irrcounter,"of", NrIS ,"irr. sec. inv. in degree ",i-1; 5964 5964 } 5965 5965 Reductor = ideal(helpP); … … 6016 6016 { ISSort[i-1] = ideal(B[j]); 6017 6017 } 6018 if (v) 6019 { " We found", irrcounter,"of", NrIS ,"irr. sec. inv. in degree ",i-1; 6018 if (v) 6019 { " We found", irrcounter,"of", NrIS ,"irr. sec. inv. in degree ",i-1; 6020 6020 } 6021 6021 Reductor = ideal(helpP); … … 6099 6099 hence the number of secondary invariants is the product of the degrees of 6100 6100 primary invariants divided by the group order. 6101 <secondary_and_irreducible _no_molien> should usually be faster and of more useful6101 <secondary_and_irreducibles_no_molien> should usually be faster and of more useful 6102 6102 functionality. 6103 SEE ALSO: secondary_and_irreducible _no_molien6103 SEE ALSO: secondary_and_irreducibles_no_molien 6104 6104 EXAMPLE: example secondary_no_molien; shows an example 6105 6105 " … … 6759 6759 THEORY: Irred. secondary invariants are calculated by finding a basis (in terms of 6760 6760 monomials) of the basering modulo primary and previously found secondary 6761 invariants, mapping those to invariants with the Reynolds operator. Among 6762 these images we pick a maximal subset that is linearly independent modulo 6763 the primary and previously found secondary invariants, using Groebner basis 6761 invariants, mapping those to invariants with the Reynolds operator. Among 6762 these images we pick a maximal subset that is linearly independent modulo 6763 the primary and previously found secondary invariants, using Groebner basis 6764 6764 techniques. 6765 6765 EXAMPLE: example irred_secondary_no_molien; shows an example … … 6815 6815 return(); 6816 6816 } 6817 if (ncols(REY)<>n) 6817 if (ncols(REY)<>n) 6818 6818 { "ERROR: The second parameter ought to be the Reynolds operator."; 6819 6819 return(); … … 6838 6838 degBound = 0; 6839 6839 int d,block, counter, totalcount, monsize,i,j, fin,hlp; 6840 poly helpP,lmp; 6840 poly helpP,lmp; 6841 6841 block = 1; 6842 6842 ideal mon, MON, Inv, IS,RedInv, NewIS, NewReductor; … … 6844 6844 ideal sP = slimgb(ideal(P)); 6845 6845 int Max = vdim(sP); 6846 if (Max<0) 6846 if (Max<0) 6847 6847 { "ERROR: The second parameter ought to be the Reynolds operator."; 6848 6848 return(); … … 6864 6864 block=0; // Loops through the monomials 6865 6865 while (block<monsize) 6866 { if ((block mod MonStep) == 0) 6866 { if ((block mod MonStep) == 0) 6867 6867 { if ((block+MonStep) <= monsize) 6868 6868 { ideal tmp = normalize(evaluate_reynolds(REY,ideal(mon[block+1..block+MonStep]))); … … 6879 6879 kill tmp; 6880 6880 } 6881 } 6881 } 6882 6882 block++; 6883 6883 helpP = RedInv[block]; … … 6891 6891 NewIS[NewIScounter] = helpP; 6892 6892 mon[block]=0; 6893 Inv[block]=0; 6893 Inv[block]=0; 6894 6894 RedInv[block]=0; 6895 6895 NewReductor[1]=helpP; … … 6904 6904 { degBound = 0; 6905 6905 sIS = groebner(sIS); 6906 fin=1; 6906 fin=1; 6907 6907 } 6908 6908 } 6909 6909 else 6910 6910 { degBound = d+1; 6911 sIS = groebner(sIS+NewIS); 6911 sIS = groebner(sIS+NewIS); 6912 6912 degBound = 0; 6913 6913 fin = 0; … … 6915 6915 } 6916 6916 attrib(sIS,"isSB",1); 6917 } // if degree d is not excluded by deg_vec 6917 } // if degree d is not excluded by deg_vec 6918 6918 else 6919 6919 { if (size(deg_vec)==l) … … 7963 7963 ASSUME: mon contains at least the minimal monomial in each orbit sum of degree d 7964 7964 under the permutation group given by Per. 7965 RETURN: an <ideal> whose generators are the orbit sums of degree d under the permutation 7965 RETURN: an <ideal> whose generators are the orbit sums of degree d under the permutation 7966 7966 group given by Per. 7967 7967 SEE ALSO: invariant_algebra_perm … … 7969 7969 EXAMPLE: example orbit_sums; shows an example 7970 7970 " 7971 { 7971 { 7972 7972 //----------------- checking input ------------------ 7973 7973 mon=simplify(mon,2); // remove 0 … … 7975 7975 { "ERROR: Monomials in the first parameter must be of the same degree."; 7976 7976 return(); 7977 } 7977 } 7978 7978 int k,k2; 7979 7979 for (k=1;k<=ncols(mon);k++) … … 7989 7989 return(); 7990 7990 } 7991 if ((ncols(M[k2])!=nvars(basering)) 7991 if ((ncols(M[k2])!=nvars(basering)) 7992 7992 or (ncols(M[k2])!=ncols(simplify(M[k2],6))) 7993 7993 or (M[k2][1]==0)) … … 8029 8029 { if ( mon[k]!=0) 8030 8030 { if (tmpI[k] < mon[k]) 8031 { 8031 { 8032 8032 mon[k]=0; 8033 8033 } 8034 8034 else 8035 8035 { if (tmpI[k]>mon[k]) 8036 { 8036 { 8037 8037 mon=NF(mon,std(tmpI[k])); 8038 8038 } … … 8052 8052 poly tmp; 8053 8053 int GenOK, j,i; 8054 poly OrbitS=mon[1]; 8054 poly OrbitS=mon[1]; 8055 8055 mon[1]=0; 8056 8056 int countOrbits; … … 8073 8073 newL[k] = simplify(NF(newL[k],L[k]),2); 8074 8074 attrib(newL[k],"isSB",1); 8075 if (newL[k][1] != 0) // i.e., the Orbit is not complete under all generators yet 8076 { GenOK=-1; 8075 if (newL[k][1] != 0) // i.e., the Orbit is not complete under all generators yet 8076 { GenOK=-1; 8077 8077 L[k] = L[k],newL[k]; 8078 8078 attrib(L[k],"isSB",1); … … 8099 8099 if (GenOK == 0) // i.e., all orbits are complete, so we can break the iteration 8100 8100 { break; 8101 } 8101 } 8102 8102 } // iterate 8103 8103 degBound = dgb; … … 8112 8112 orb; 8113 8113 } 8114 8114 8115 8115 /////////////////////////////////////////////////////////////////////////////// 8116 8116 8117 8117 proc invariant_algebra_perm (list #) 8118 8118 "USAGE: invariant_algebra_perm(GEN[,v]); 8119 @* GEN: a list of generators of a permutation group. It is given in disjoint cycle 8120 form, where trivial cycles can be omitted; e.g., the generator (1,2)(3,4)(5) is 8119 @* GEN: a list of generators of a permutation group. It is given in disjoint cycle 8120 form, where trivial cycles can be omitted; e.g., the generator (1,2)(3,4)(5) is 8121 8121 given by <list(list(1,2),list(3,4))>. 8122 8122 @* v: an optional <int> … … 8124 8124 by GEN, type <matrix> 8125 8125 ASSUME: We are in the non-modular case, i.e., the characteristic of the basering 8126 does not divide the group order. Note that the function does not verify whether 8126 does not divide the group order. Note that the function does not verify whether 8127 8127 the assumption holds or not 8128 8128 DISPLAY: Information on the progress of computations if v does not equal 0 … … 8134 8134 EXAMPLE: example invariant_algebra_perm; shows an example 8135 8135 " 8136 { 8136 { 8137 8137 //----------------- checking input and setting verbose mode ------------------ 8138 8138 if (size(#)==0) … … 8192 8192 NextCand, sGoffset, fin,MaxD,LastGDeg; 8193 8193 ideal OrbSums, RedSums, NewReductor; 8194 poly helpP,lmp; 8194 poly helpP,lmp; 8195 8195 8196 8196 … … 8216 8216 while ((fin==0) or (d<=MaxD)) 8217 8217 { if (v) 8218 { "Searching generators in degree",d; 8218 { "Searching generators in degree",d; 8219 8219 } 8220 8220 counter = 0; // counts generators in degree d … … 8223 8223 sGoffset=ncols(sG); 8224 8224 OrbSize = ncols(OrbSums); 8225 if (v) 8225 if (v) 8226 8226 { " We have "+string(OrbSize)+" orbit sums of degree "+string(d); 8227 8227 } 8228 8228 // loop through orbit_sum, and select generators among them 8229 for (i=1;i<=OrbSize;i++) 8229 for (i=1;i<=OrbSize;i++) 8230 8230 { helpP = RedSums[i]; 8231 8231 if (helpP<>0) 8232 8232 { counter++; // found new inv. algebra generator! 8233 8233 totalcount++; 8234 if (v) 8234 if (v) 8235 8235 { " We found generator number "+string(totalcount)+ 8236 " in degree "+string(d); 8236 " in degree "+string(d); 8237 8237 } 8238 8238 G[totalcount] = OrbSums[i]; … … 8256 8256 { degBound=0; 8257 8257 if (v) {"Presumably all generators have been found";} 8258 sG = groebner(sG); 8258 sG = groebner(sG); 8259 8259 fin = 0; 8260 8260 } … … 8268 8268 if ((fin==0) and (dim(sG)==0)) 8269 8269 { MaxD = maxdeg1(kbase(sG)); 8270 if (v) 8270 if (v) 8271 8271 { if (MaxD>=d) {"We found the degree bound",MaxD; } 8272 8272 else {"We found the degree bound",d-1;} … … 8277 8277 } 8278 8278 fin = 1; 8279 } // computing degree bound 8279 } // computing degree bound 8280 8280 } // for increasing degrees 8281 8281 kill sG; … … 8295 8295 proc invariant_algebra_reynolds (matrix REY, list #) 8296 8296 "USAGE: invariant_algebra_reynolds(REY[,v]); 8297 @* REY: a gxn <matrix> representing the Reynolds operator of a finite matrix group, 8297 @* REY: a gxn <matrix> representing the Reynolds operator of a finite matrix group, 8298 8298 where g ist the group order and n is the number of variables of the basering; 8299 8299 @* v: an optional <int> … … 8305 8305 DISPLAY: Information on the progress of computations if v does not equal 0 8306 8306 THEORY: We do an incremental search in increasing degree d. Generators of the invariant 8307 ring are found among the Reynolds images of monomials of degree d. The generators are 8307 ring are found among the Reynolds images of monomials of degree d. The generators are 8308 8308 chosen by Groebner basis techniques. 8309 8309 SEE ALSO: invariant_algebra_perm … … 8318 8318 int block, counter, totalcount, monsize,i,j, 8319 8319 sGoffset, fin,MaxD,LastGDeg; 8320 poly helpP,lmp; 8320 poly helpP,lmp; 8321 8321 block = 1; 8322 8322 ideal mon, MON, Inv, G; // G will contain homog. alg.generators … … 8339 8339 block=0; // Loops through the monomials 8340 8340 while (block<monsize) 8341 { if ((block mod MonStep) == 0) 8341 { if ((block mod MonStep) == 0) 8342 8342 { if ((block+MonStep) <= monsize) 8343 8343 { ideal tmp = normalize(evaluate_reynolds(REY,ideal(mon[block+1..block+MonStep]))); … … 8354 8354 kill tmp; 8355 8355 } 8356 } 8356 } 8357 8357 block++; 8358 8358 helpP = RedInv[block]; … … 8367 8367 NewG[NewGcounter] = helpP; 8368 8368 mon[block]=0; 8369 Inv[block]=0; 8369 Inv[block]=0; 8370 8370 RedInv[block]=0; 8371 8371 NewReductor[1]=helpP; … … 8387 8387 { degBound=0; 8388 8388 if (v) {"Presumably all generators have been found";} 8389 sG = groebner(sG); 8389 sG = groebner(sG); 8390 8390 fin = 0; 8391 8391 } … … 8399 8399 if ((fin==0) and (dim(sG)==0)) 8400 8400 { MaxD = maxdeg1(kbase(sG)); 8401 if (v) 8401 if (v) 8402 8402 { if (MaxD>=d) {"We found the degree bound",MaxD; } 8403 8403 else {"We found the degree bound",d-1;} … … 8408 8408 } 8409 8409 fin = 1; 8410 } // computing degree bound 8410 } // computing degree bound 8411 8411 } // for increasing degrees 8412 8412 kill sG;
Note: See TracChangeset
for help on using the changeset viewer.