Changeset 3bb38d3 in git
- Timestamp:
- Mar 5, 2007, 8:39:23 PM (17 years ago)
- Branches:
- (u'spielwiese', '4a9821a93ffdc22a6696668bd4f6b8c9de3e6c5f')
- Children:
- f3587d92249476f601fedf863d711cd00d783635
- Parents:
- f097f1d706a913c97d6b90901ad7dd0e830cb1dd
- File:
-
- 1 edited
-
Singular/LIB/finvar.lib (modified) (15 diffs)
Legend:
- Unmodified
- Added
- Removed
-
Singular/LIB/finvar.lib
rf097f1 r3bb38d3 1 1 /////////////////////////////////////////////////////////////////////////////// 2 version="$Id: finvar.lib,v 1.6 5 2007-02-24 15:08:04 SingularExp $"2 version="$Id: finvar.lib,v 1.66 2007-03-05 19:39:23 king Exp $" 3 3 category="Invariant theory"; 4 4 info=" … … 5040 5040 int LastNewIS=-10; // if the Last New Irr. Sec. was found in degree i-2 then 5041 5041 // we compute TotalReductor again 5042 int TotalSet = 1; // It TotalSet is 0, we eventually need to re-compute TotalReductor.5042 int NoPP; // It NoPP==1, we must not use Power Products 5043 5043 ideal Reductor,SaveRed; // sP union Reductor is a Groebner basis up to degree i-1 5044 5044 int SizeSave; … … 5067 5067 //--------------------- generating secondary invariants ---------------------- 5068 5068 for (i=2;i<=m;i++) // going through dimmat - 5069 { if ((LastNewIS == i-3) and (UsePP==0)) // We guess that all irreducibles have been found! 5069 { // Case 1: PP are not required by the user, and we guess that all irreducibles have 5070 // been found. 5071 // Afterwards, the use of PP is impossible. 5072 if ((LastNewIS == i-3) and (UsePP==0)) 5070 5073 { if (v) {"Computing Groebner basis for primaries and previously found irreducibles...";} 5074 degBound = 0; 5071 5075 TotalReductor = groebner(sP+IS); 5072 TotalSet = 1; 5076 NoPP = 1; 5077 } 5078 // Case 2: The use of PP is impossible, but we did find a new irr. sec. inv. 5079 // in the preceding degree 5080 if ((LastNewIS == i-2) and (NoPP==1)) 5081 { degBound = i-1; 5082 TotalReductor = groebner(TotalReductor+SaveRed); // lift from degree i-2 to degree i-1 5083 attrib(TotalReductor, "isSB",1); 5073 5084 } 5074 5085 if (int(dimmat[i,1])<>0) // when it is == 0 we need to find no … … 5085 5096 attrib(Reductor,"isSB",1); 5086 5097 attrib(SaveRed,"isSB",1); 5087 5088 if ( TotalSet==0 ) 5098 5099 // Case A: We use PP 5100 if ( NoPP==0 ) 5089 5101 { if (v) 5090 5102 { 5091 5103 " Looking for Power Products..."; 5092 5104 } 5093 sP_Reductor = sP;5105 sP_Reductor = sP; 5094 5106 // We start searching for reducible secondary invariants in degree i-1, i.e., those 5095 5107 // that are power products of irreducible secondary invariants. … … 5102 5114 // work with their reduction modulo sP --- this allows to detect a secondary invariant 5103 5115 // without to actually compute it! 5104 for (k=1;k<i-1;k++) 5105 { if (int(dimmat[i,1])==counter) 5106 { break; 5107 } 5108 if (typeof(RedISSort[k])<>"none") 5109 { Mul1 = RedISSort[k]; 5110 } 5111 else 5112 { Mul1 = ideal(0); 5113 } 5114 if ((int(dimmat[i-k,1])>0) && (Mul1[1]<>0)) 5115 { for (minD=i-k-1;minD>0;minD--) 5116 { if (int(dimmat[i,1])==counter) 5117 { break; 5118 } 5119 for (k2=1;k2 <= ((attrib(RedSSort[i-k-1][minD],"size")-1) div pieces)+1; k2++) 5116 for (k=1;k<i-1;k++) 5117 { if (int(dimmat[i,1])==counter) 5118 { break; 5119 } 5120 if (typeof(RedISSort[k])<>"none") 5121 { Mul1 = RedISSort[k]; 5122 } 5123 else 5124 { Mul1 = ideal(0); 5125 } 5126 if ((int(dimmat[i-k,1])>0) && (Mul1[1]<>0)) 5127 { for (minD=i-k-1;minD>0;minD--) 5120 5128 { if (int(dimmat[i,1])==counter) 5121 5129 { break; 5122 5130 } 5123 Mul2=ideal(0); 5124 if (attrib(RedSSort[i-k-1][minD],"size")>=k2*pieces) 5125 { for (k3=1;k3<=pieces;k3++) 5126 { Mul2[k3] = RedSSort[i-k-1][minD][((k2-1)*pieces)+k3]; 5131 for (k2=1;k2 <= ((attrib(RedSSort[i-k-1][minD],"size")-1) div pieces)+1; k2++) 5132 { if (int(dimmat[i,1])==counter) 5133 { break; 5127 5134 } 5135 Mul2=ideal(0); 5136 if (attrib(RedSSort[i-k-1][minD],"size")>=k2*pieces) 5137 { for (k3=1;k3<=pieces;k3++) 5138 { Mul2[k3] = RedSSort[i-k-1][minD][((k2-1)*pieces)+k3]; 5139 } 5140 } 5141 else 5142 { for (k3=1;k3<=(attrib(RedSSort[i-k-1][minD],"size") mod pieces);k3++) 5143 { Mul2[k3] = RedSSort[i-k-1][minD][((k2-1)*pieces)+k3]; 5144 } 5145 } 5146 ProdCand = simplify(Mul1*Mul2,4); 5147 ReducedCandidates = reduce(ProdCand,sP); 5148 // sP union SaveRed union Reductor is a homogeneous Groebner basis 5149 // up to degree i-1. 5150 // We first reduce by sP (which is fixed, so we can do it once for all), 5151 // then by SaveRed resp. by Reductor (which is modified during 5152 // the computations). 5153 Indicator = reduce(ReducedCandidates,SaveRed); 5154 // If Indicator[ii]==0 then ReducedCandidates it the reduction 5155 // of an invariant that is in the algebra generated by primary 5156 // invariants and previously computed secondary invariants. 5157 // Otherwise ReducedCandidates[ii] is the reduction of an invariant 5158 // that we can take as secondary invariant. 5159 if (size(Indicator)<>0) 5160 { for (ii=1;ii<=ncols(ProdCand);ii++) // going through all the power products 5161 { helpP = Indicator[ii]; 5162 if (helpP <> 0) 5163 { counter++; 5164 saveAttr = attrib(RedSSort[i-1][k],"size")+1; 5165 RedSSort[i-1][k][saveAttr] = ReducedCandidates[ii]; 5166 // By construction, this is the reduction of a reducible s.i. 5167 // of degree i-1. 5168 attrib(RedSSort[i-1][k],"size",saveAttr); 5169 if (v) 5170 { " We found reducible sec. inv. number "+string(counter); 5171 } 5172 if (int(dimmat[i,1])<>counter) 5173 { Reductor = ideal(helpP); 5174 attrib(Reductor, "isSB",1); 5175 Indicator=reduce(Indicator,Reductor); 5176 SizeSave++; 5177 SaveRed[SizeSave] = helpP; 5178 sP_Reductor[sPOffset+SizeSave] = helpP; 5179 attrib(sP_Reductor, "isSB",1); 5180 // Lemma: If G is a homogeneous Groebner basis up to degree i-1 and p is a 5181 // homogeneous polynomial of degree i-1 then G union NF(p,G) is 5182 // a homogeneous Groebner basis up to degree i-1. Hence, sP_Reductor 5183 // is a homog. GB up to degree i-1. 5184 attrib(SaveRed, "isSB",1); 5185 attrib(Reductor, "isSB",1); 5186 } 5187 else 5188 { break; 5189 } 5190 } // new reducible sec. inv. found 5191 } // loop through ProdCand 5192 } // if there is some reducible sec. inv. 5193 attrib(SaveRed, "isSB",1); 5194 attrib(Reductor, "isSB",1); 5128 5195 } 5129 else5130 { for (k3=1;k3<=(attrib(RedSSort[i-k-1][minD],"size") mod pieces);k3++)5131 { Mul2[k3] = RedSSort[i-k-1][minD][((k2-1)*pieces)+k3];5132 }5133 }5134 ProdCand = simplify(Mul1*Mul2,4);5135 ReducedCandidates = reduce(ProdCand,sP);5136 // sP union SaveRed union Reductor is a homogeneous Groebner basis5137 // up to degree i-1.5138 // We first reduce by sP (which is fixed, so we can do it once for all),5139 // then by SaveRed resp. by Reductor (which is modified during5140 // the computations).5141 Indicator = reduce(ReducedCandidates,SaveRed);5142 // If Indicator[ii]==0 then ReducedCandidates it the reduction5143 // of an invariant that is in the algebra generated by primary5144 // invariants and previously computed secondary invariants.5145 // Otherwise ReducedCandidates[ii] is the reduction of an invariant5146 // that we can take as secondary invariant.5147 if (size(Indicator)<>0)5148 { for (ii=1;ii<=ncols(ProdCand);ii++) // going through all the power products5149 { helpP = Indicator[ii];5150 if (helpP <> 0)5151 { counter++;5152 saveAttr = attrib(RedSSort[i-1][k],"size")+1;5153 RedSSort[i-1][k][saveAttr] = ReducedCandidates[ii];5154 // By construction, this is the reduction of a reducible s.i.5155 // of degree i-1.5156 attrib(RedSSort[i-1][k],"size",saveAttr);5157 if (v)5158 {5159 " We found reducible sec. inv. number "+string(counter);5160 }5161 if (int(dimmat[i,1])<>counter)5162 { Reductor = ideal(helpP);5163 attrib(Reductor, "isSB",1);5164 Indicator=reduce(Indicator,Reductor);5165 SizeSave++;5166 SaveRed[SizeSave] = helpP;5167 sP_Reductor[sPOffset+SizeSave] = helpP;5168 attrib(sP_Reductor, "isSB",1);5169 // Lemma: If G is a homogeneous Groebner basis up to degree i-1 and p is a5170 // homogeneous polynomial of degree i-1 then G union NF(p,G) is5171 // a homogeneous Groebner basis up to degree i-1. Hence, sP_Reductor5172 // is a homog. GB up to degree i-1.5173 attrib(SaveRed, "isSB",1);5174 attrib(Reductor, "isSB",1);5175 }5176 else5177 { break;5178 }5179 } // new reducible sec. inv. found5180 } // loop through ProdCand5181 } // if there is some reducible sec. inv.5182 attrib(SaveRed, "isSB",1);5183 attrib(Reductor, "isSB",1);5184 5196 } 5185 5197 } 5186 5198 } 5187 }5188 } // if TotalSet==05189 else5190 { if (v) {" We don't compute Power Products!"; }5199 TotalReductor = sP_Reductor; 5200 } // if NoPP==0 5201 else 5202 { if (v) {" We don't compute Power Products!"; } 5191 5203 // Instead of computing Power Products, we can compute a Groebner basis of the ideal 5192 5204 // generated by the primary and previously found irr. sec. invariants up to degree i-2. An invariant … … 5195 5207 // Hence, if we find an invariant outside this ideal, it is an irreducible secondary 5196 5208 // invariant of degree i-1. 5197 } 5198 // The remaining sec. inv. are irreducible! 5199 if (int(dimmat[i,1])<>counter) // need more than all the power products 5209 } // dealing with reducible sec. inv. 5210 // The remaining sec. inv. are irreducible! 5211 // Reason: If we use PP, TotalReductor detects the <<P>>-module generated by 5212 // all secondaries 5213 // If we don't, TotalReductor detects the algebra <<P union IS>> 5214 if (int(dimmat[i,1])<>counter) // need more than all the power products 5200 5215 // or work without power products 5201 5216 { if (v) 5202 5217 { " Looking for irreducible secondary invariants in degree ", i-1; 5203 5218 } 5204 if (TotalSet==0) 5205 { mon = kbase(sP_Reductor,i-1); 5206 } 5207 else 5208 { mon = kbase(TotalReductor,i-1); 5209 } 5219 mon = kbase(TotalReductor,i-1); 5210 5220 ii = ncols(mon); 5211 5221 if (mon[1]<>0) 5212 5222 { if ((i>IrrSwitch) or (ii>200)) // use sparse algorithm 5213 { if ( TotalSet==0 && ii+counter==int(dimmat[i,1])) // then we can use all of mon5223 { if (NoPP==0 && ii+counter==int(dimmat[i,1])) // then we can use all of mon 5214 5224 { B=normalize(evaluate_reynolds(REY,mon)); 5215 5225 IS=IS+B; … … 5223 5233 { RedISSort[i-1] = NF(B,sP); 5224 5234 } 5225 TotalSet=0; // i.e., TotalReductor needs to be re-computed5226 5235 if (v) {" We found "+string(size(B))+" irred. sec. inv.";} 5227 5236 } … … 5235 5244 { ideal tmp = normalize(evaluate_reynolds(REY,ideal(mon[j+1..j+MonStep]))); 5236 5245 B[j+1..j+MonStep] = tmp[1..MonStep]; 5237 tmp = reduce(tmp,sP); 5238 ReducedCandidates[j+1..j+MonStep] = tmp[1..MonStep]; 5239 if (TotalSet == 0) 5240 { tmp = reduce(tmp,SaveRed); 5246 if (NoPP == 0) // we are still working with PowerProducts, 5247 // hence, we need to store NF(sec.inv.,sP) 5248 { tmp = reduce(tmp,sP); 5249 ReducedCandidates[j+1..j+MonStep] = tmp[1..MonStep]; 5250 tmp = reduce(tmp,SaveRed); 5241 5251 } 5242 5252 else … … 5249 5259 { ideal tmp = normalize(evaluate_reynolds(REY,ideal(mon[j+1..ii]))); 5250 5260 B[j+1..ii] = tmp[1..ii-j]; 5251 tmp = reduce(tmp,sP);5252 ReducedCandidates[j+1..ii] = tmp[1..ii-j];5253 if (TotalSet == 0)5254 {tmp = reduce(tmp,SaveRed);5261 if (NoPP == 0) 5262 { tmp = reduce(tmp,sP); 5263 ReducedCandidates[j+1..ii] = tmp[1..ii-j]; 5264 tmp = reduce(tmp,SaveRed); 5255 5265 } 5256 5266 else … … 5264 5274 helpP = Indicator[j]; 5265 5275 if (helpP <>0) // B[j] should be added 5266 { TotalSet=0; // i.e., TotalReductor needs to be re-computed 5267 counter++; irrcounter++; 5276 { counter++; irrcounter++; 5268 5277 IS=IS,B[j]; 5269 5278 saveAttr = attrib(RedSSort[i-1][i-1],"size")+1; 5270 RedSSort[i-1][i-1][saveAttr] = ReducedCandidates[j]; 5271 attrib(RedSSort[i-1][i-1],"size",saveAttr); 5272 if (typeof(RedISSort[i-1]) <> "none") 5273 { RedISSort[i-1][irrcounter] = ReducedCandidates[j]; 5274 } 5275 else 5276 { RedISSort[i-1] = ideal(ReducedCandidates[j]); 5279 if (NoPP==0) // we will still be working with Power Products 5280 { RedSSort[i-1][i-1][saveAttr] = ReducedCandidates[j]; 5281 attrib(RedSSort[i-1][i-1],"size",saveAttr); 5282 if (typeof(RedISSort[i-1]) <> "none") 5283 { RedISSort[i-1][irrcounter] = ReducedCandidates[j]; 5284 } 5285 else 5286 { RedISSort[i-1] = ideal(ReducedCandidates[j]); 5287 } 5277 5288 } 5278 5289 if (v) … … 5297 5308 } // if i>IrrSwitch 5298 5309 else // use fast algorithm 5299 { if (TotalSet == 0) 5300 { B=sort_of_invariant_basis(sP_Reductor,REY,i-1,int(dimmat[i,1])*6); 5301 } 5302 else 5303 { B=sort_of_invariant_basis(TotalReductor,REY,i-1,int(dimmat[i,1])*6); 5304 } 5310 { B=sort_of_invariant_basis(TotalReductor,REY,i-1,int(dimmat[i,1])*6); 5305 5311 // B is a linearly independent set of reynolds 5306 5312 // images of monomials that do not occur … … 5308 5314 // by primary and previously found irreducible 5309 5315 // secondary invariants. 5310 if ( TotalSet==0 && ncols(B)+counter==int(dimmat[i,1])) // then we can take all of B5316 if (NoPP==0 && ncols(B)+counter==int(dimmat[i,1])) // then we can take all of B 5311 5317 { IS=IS+B; 5312 5318 saveAttr = attrib(RedSSort[i-1][i-1],"size")+int(dimmat[i,1]); … … 5332 5338 helpP = Indicator[j]; 5333 5339 if (helpP <>0) // B[j] should be added 5334 { TotalSet=0; // i.e., TotalReductor needs to be re-computed5340 { NoPP=0; // i.e., TotalReductor needs to be re-computed 5335 5341 counter++; irrcounter++; 5336 5342 IS=IS,B[j]; … … 6845 6851 int Max = vdim(sP); 6846 6852 if (Max<0) 6847 { "ERROR: The second parameter ought to be the Reynolds operator.";6853 { "ERROR: The first parameter ought to be a matrix of primary invariants."; 6848 6854 return(); 6849 6855 } … … 8220 8226 counter = 0; // counts generators in degree d 8221 8227 OrbSums = orbit_sums(kbase(sG,d),M); 8222 RedSums = reduce(OrbSums,sG);8223 sGoffset=ncols(sG);8224 8228 OrbSize = ncols(OrbSums); 8225 8229 if (v) 8226 8230 { " We have "+string(OrbSize)+" orbit sums of degree "+string(d); 8227 8231 } 8232 RedSums = reduce(OrbSums,sG); 8233 sGoffset=ncols(sG); 8228 8234 // loop through orbit_sum, and select generators among them 8229 8235 for (i=1;i<=OrbSize;i++)
Note: See TracChangeset
for help on using the changeset viewer.
