1 | //////////////////////////////////////////////////////////////////////////////// |
---|
2 | version="version autgradalg.lib 4.1.2.0 Feb_2019 "; // $Id$ |
---|
3 | category="Commutative Algebra, Algebraic Geometry"; |
---|
4 | info=" |
---|
5 | LIBRARY: autgradalg.lib Compute automorphism groups of pointedly graded algebras and of Mori dream spaces. |
---|
6 | |
---|
7 | AUTHORS: Simon Keicher |
---|
8 | |
---|
9 | OVERVIEW: This library provides a framework for computing automorphisms of integral, finitely generated algebras that are graded by a finitely generated abelian group. This library also contains functions to compute automorphism groups of Mori dream spaces. The results are ideals I such that the respective automorphism group is isomorphic to the subgroup V(I) in some GL(n). |
---|
10 | |
---|
11 | ASSUMPTIONS: |
---|
12 | * the algebra R is given as factor algebra S/I with a graded polynomial ring S = KK[T_1,...,T_r]. We will always assume that the basering is S and it is given over the rationals QQ or a number field QQ(a). |
---|
13 | * R must be minimally presented, i.e., I is contained in <T_1,...,T_r>^2. |
---|
14 | * S (and hence R) are graded via 'setBaseMultigrading(Q)' from 'multigrading.lib'. The last rows of the matrix Q are interpreted over ZZ/a_iZZ if the respective entry of the list TOR is a_i and has been provided as parameter to the respective function. (See the functions for more details.) |
---|
15 | * For all 1 <= i <= r: I_{w_i} = 0 where w_i := deg(T_i). |
---|
16 | * the grading is pointed, i.e., no generator has degree 0 and the cone over all generator degrees is pointed. |
---|
17 | * For Mori dream spaces X, we assume them to be given as X = X(R,w) with the Cox ring R of X (given as the algebra R before) and an ample class w in the grading group K with the torsion entries removed. |
---|
18 | |
---|
19 | NOTE: we require the user to execute 'LIB'gfanlib.so'' before using this library. |
---|
20 | |
---|
21 | PROCEDURES: |
---|
22 | autKS(): compute the automorphism group of the basering (must be a polynomial ring) as an algebraic subgroup V(I) of some GL(n) |
---|
23 | autGradAlg(I0, TOR): compute the automorphism group of R as an algebraic subgroup V(I) of some GL(n). |
---|
24 | autGenWeights(): compute the automorphisms of the grading group that fix the generator degrees. |
---|
25 | stabilizer(I0, TOR): compute the stabilizer of the given ideal |
---|
26 | autXhat(I0, w, TOR): compute the automorphism group of \widehat X as an algebraic subgroup V(I) of some GL(n). |
---|
27 | autX(I0, w, TOR): compute the automorphism group of X=X(R,w) as an algebraic subgroup V(I) of some GL(n). |
---|
28 | |
---|
29 | NOTE: the following functions were taken from 'quotsingcox.lib' by M.Donten-Bury and S.Keicher: 'hilbertBas'. |
---|
30 | NOTE: This library comes without any warranty whatsoever. Use it at your own risk. |
---|
31 | |
---|
32 | KEYWORDS: automorphisms; graded algebra; Mori dream spaces; automorphism group; symmetries |
---|
33 | "; |
---|
34 | |
---|
35 | /* Printlevel settings: |
---|
36 | * 0: Silent mode |
---|
37 | * >0: Show status information |
---|
38 | */ |
---|
39 | |
---|
40 | //Included libraries |
---|
41 | LIB "multigrading.lib"; |
---|
42 | LIB "gitfan.lib"; |
---|
43 | LIB "normaliz.lib"; |
---|
44 | |
---|
45 | //////////////////////// |
---|
46 | |
---|
47 | static proc col(intmat Q, int i) |
---|
48 | "USAGE: col(Q, i): Q is an intmat, i an integer. |
---|
49 | PURPOSE: return the i-th column as an intvec. |
---|
50 | RETURN: an intvec. |
---|
51 | EXAMPLE: shows an example." |
---|
52 | { |
---|
53 | return(intvec(Q[1..nrows(Q),i])); |
---|
54 | } |
---|
55 | example |
---|
56 | { |
---|
57 | echo=2; |
---|
58 | |
---|
59 | intmat Q[2][4] = |
---|
60 | -2, 2, -1, 1, |
---|
61 | 1, 1, 1, 1; |
---|
62 | |
---|
63 | col(Q, 2); |
---|
64 | } |
---|
65 | |
---|
66 | /////////////////////// |
---|
67 | |
---|
68 | static proc gradedCompIdeal(ideal RL, intvec w, int dd, list #) |
---|
69 | "USAGE: gradedCompIdeal(I, w, dd): I is an ideal, w an intvec, dd an integer (0 or 1). Optional input: a list of integers representing the torsion. |
---|
70 | ASSUME: the basering must be graded (see setBaseMultigrading()) and the cone over the degrees of the variables must be pointed; there must not be 0-degrees. The vector w must be an element of the cone over the degrees of the variables. |
---|
71 | PURPOSE: returns a basis for the component I_w of the ideal; if dd = 0, a list of basis vectors is returned, if dd = 1, a list of monomials is returned. |
---|
72 | RETURN: a list L where L[1] is a list of all monomials of degree w and L[2] is list of polynomials. |
---|
73 | EXAMPLE: shows an example." |
---|
74 | { |
---|
75 | list TOR; |
---|
76 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
77 | TOR = #[1]; |
---|
78 | } |
---|
79 | |
---|
80 | // columns of Q are the degrees of the variables of the basering |
---|
81 | intmat Q = getVariableWeights(); |
---|
82 | |
---|
83 | // compute all monomials of degree w: |
---|
84 | dbprint(printlevel, "wmonomials..."); |
---|
85 | list MON = wmonomials(w, 0, TOR); |
---|
86 | dbprint(printlevel, "...done"); |
---|
87 | |
---|
88 | list Iw0 = gradedCompIdeal2basevecs(w, Q, RL, MON, TOR); |
---|
89 | |
---|
90 | if(size(Iw0) == 0){ |
---|
91 | list L; |
---|
92 | L[1] = MON; |
---|
93 | L[2] = list(); |
---|
94 | return(L); |
---|
95 | } |
---|
96 | |
---|
97 | // Choose a maximal linearly independent |
---|
98 | // subset out of the elements of Iw0: |
---|
99 | // Iteratively add columns from Tmp: |
---|
100 | int rg1 = 0; |
---|
101 | int rg2 = 0; |
---|
102 | |
---|
103 | list Tmp; |
---|
104 | Tmp[1] = Iw0[1]; |
---|
105 | |
---|
106 | for(int m = 2; m <= size(Iw0); m++){ |
---|
107 | matrix IwTmp[size(Tmp)+1][size(MON)]; |
---|
108 | |
---|
109 | for(int j = 1; j <= size(Tmp); j++){ |
---|
110 | IwTmp[j,1..size(MON)] = Tmp[j]; |
---|
111 | } |
---|
112 | |
---|
113 | // add the new column: |
---|
114 | IwTmp[size(Tmp)+1,1..size(MON)] = Iw0[m]; |
---|
115 | |
---|
116 | rg1 = rg2; |
---|
117 | rg2 = rank(IwTmp); |
---|
118 | |
---|
119 | if(rg1 < rg2){ |
---|
120 | Tmp[size(Tmp)+1] = Iw0[m]; |
---|
121 | } |
---|
122 | |
---|
123 | kill IwTmp, j; |
---|
124 | } |
---|
125 | |
---|
126 | list Iw = Tmp; |
---|
127 | list L; // will be returned |
---|
128 | L[1] = MON; |
---|
129 | |
---|
130 | // translate back to polynomials unless 'basevecs' specified: |
---|
131 | if(dd == 0){ |
---|
132 | L[2] = Iw; |
---|
133 | return(L); |
---|
134 | } |
---|
135 | |
---|
136 | matrix MONmat[1][size(MON)] = MON[1..size(MON)]; |
---|
137 | matrix Iwmat[size(Iw)][size(MON)]; |
---|
138 | |
---|
139 | for(m = 1; m <= size(Iw); m++){ |
---|
140 | matrix BB = Iw[m]; // has only one column |
---|
141 | Iwmat[m, 1..size(MON)] = BB[1..nrows(BB), 1]; |
---|
142 | |
---|
143 | kill BB; |
---|
144 | } |
---|
145 | |
---|
146 | matrix Iwmons = Iwmat * transpose(MONmat); |
---|
147 | list Iwmonslist = Iwmons[1..nrows(Iwmons),1]; |
---|
148 | L[2] = Iwmonslist; |
---|
149 | |
---|
150 | return(L); |
---|
151 | } |
---|
152 | example |
---|
153 | { |
---|
154 | echo = 2; |
---|
155 | |
---|
156 | // 1st example ///////// |
---|
157 | // example with torsion: ZZ + ZZ/5ZZ: |
---|
158 | intmat Q[2][5] = |
---|
159 | 1,1,1,1,1, |
---|
160 | 2,3,1,4,0; |
---|
161 | |
---|
162 | list TOR = 5; |
---|
163 | ring R = 0,T(1..5),dp; |
---|
164 | |
---|
165 | // attach degree matrix Q to R: |
---|
166 | setBaseMultigrading(Q); |
---|
167 | |
---|
168 | ideal I = T(1)*T(2) + T(3)*T(4) + T(5)^2; |
---|
169 | intvec w = 3,1; |
---|
170 | |
---|
171 | // as list of polynomials: |
---|
172 | list Lpol = gradedCompIdeal(I, w, 1, TOR); |
---|
173 | Lpol; |
---|
174 | |
---|
175 | kill R, w, Q, TOR; |
---|
176 | |
---|
177 | // 2nd example ///////// |
---|
178 | ring R = 0, T(1..4), dp; |
---|
179 | |
---|
180 | // Weights of variables |
---|
181 | intmat Q[2][4] = |
---|
182 | -2, 2, -1, 1, |
---|
183 | 1, 1, 1, 1; |
---|
184 | |
---|
185 | // attach degree matrix Q to R: |
---|
186 | setBaseMultigrading(Q); |
---|
187 | |
---|
188 | ideal I = T(1)*T(2) + T(3)*T(4); |
---|
189 | intvec w = 0,2; |
---|
190 | |
---|
191 | // as basis vectors in \KK^n: |
---|
192 | list Lbas = gradedCompIdeal(I, w, 0); |
---|
193 | Lbas; |
---|
194 | |
---|
195 | // as list of polynomials: |
---|
196 | list Lpol = gradedCompIdeal(I, w, 1); |
---|
197 | Lpol; |
---|
198 | |
---|
199 | kill w, I, R, Q; |
---|
200 | |
---|
201 | // 3rd example ///////// |
---|
202 | // torsion, parameter and 2 generators: |
---|
203 | intmat Q[3][6] = |
---|
204 | 1,1,1,1,1,1, |
---|
205 | -1,1,1,-1,0,0, |
---|
206 | 0,0,1,1,1,0; |
---|
207 | |
---|
208 | list TOR = 2; // means CL(X) = ZZ + ZZ/TOR[1]*ZZ + ... |
---|
209 | |
---|
210 | //int lam = 13; |
---|
211 | ring R = (0,lam),T(1..6),dp; |
---|
212 | |
---|
213 | // attach degree matrix Q to R: |
---|
214 | setBaseMultigrading(Q); |
---|
215 | |
---|
216 | ideal I = |
---|
217 | T(1)*T(2) + T(3)*T(4) + T(5)^2, |
---|
218 | lam*T(3)*T(4) + T(5)^2 + T(6)^2; |
---|
219 | |
---|
220 | intvec w = 2,0,0; |
---|
221 | gradedCompIdeal(I, w, 1, TOR); |
---|
222 | } |
---|
223 | |
---|
224 | /////////////////////////////////// |
---|
225 | |
---|
226 | // helper |
---|
227 | // returns a list of matrices over the coeff. ring of the basering |
---|
228 | // see also the subsequent example. |
---|
229 | static proc gradedCompIdeal2basevecs(intvec w, intmat Q, ideal RL, list MON, list #){ |
---|
230 | list TOR; |
---|
231 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
232 | TOR = #[1]; |
---|
233 | } |
---|
234 | |
---|
235 | list RLw = gradedCompIdealhelper(w, Q, RL, TOR); |
---|
236 | list EL; |
---|
237 | |
---|
238 | for(int j = 1; j <= size(RLw); j++){ |
---|
239 | // translate e.g. T[1]^3-T[1]*T[2]^2 to e_1-e_3 |
---|
240 | // if T[1]^3 is the 1st entry and T[1]*T[2]^2 the 3rd entry |
---|
241 | // store then 1 in v[...] and -1 in v[...]: |
---|
242 | matrix v[1][size(MON)] = 0:size(MON); // Elements are in the coef-ring |
---|
243 | |
---|
244 | for(int k = 1; k <= size(MON); k++){ |
---|
245 | v[k,1] = moncoef(RLw[j], MON[k]); // Note: it is v[k,1], not v[1,k]. |
---|
246 | } |
---|
247 | |
---|
248 | EL[size(EL)+1] = v; |
---|
249 | |
---|
250 | kill v, k; |
---|
251 | } |
---|
252 | |
---|
253 | return(EL); |
---|
254 | } |
---|
255 | example |
---|
256 | { |
---|
257 | echo=2; |
---|
258 | |
---|
259 | ring R = 0, T(1..4), dp; |
---|
260 | |
---|
261 | // Weights of variables |
---|
262 | intmat Q[2][4] = |
---|
263 | -2, 2, -1, 1, |
---|
264 | 1, 1, 1, 1; |
---|
265 | |
---|
266 | // attach degree matrix Q to R: |
---|
267 | setBaseMultigrading(Q); |
---|
268 | |
---|
269 | ideal I = T(1)*T(2) + T(3)*T(4); |
---|
270 | |
---|
271 | intvec w = 0,4; |
---|
272 | list MON = wmonomials(w, 0); |
---|
273 | MON; |
---|
274 | gradedCompIdeal2basevecs(w, Q, I, MON); |
---|
275 | } |
---|
276 | |
---|
277 | ///////////////////////////////////// |
---|
278 | |
---|
279 | // helper |
---|
280 | // |
---|
281 | // Returns vector space generators for all polynomials of degree w |
---|
282 | // by shifting the generators in RL to degree w by a monomial. |
---|
283 | // |
---|
284 | // Returns a list of polynomials. |
---|
285 | // |
---|
286 | static proc gradedCompIdealhelper(intvec w, intmat Q, ideal RL, list #){ |
---|
287 | |
---|
288 | list TOR; |
---|
289 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
290 | TOR = #[1]; |
---|
291 | } |
---|
292 | |
---|
293 | list RLw; |
---|
294 | cone cQ = coneViaPoints(transpose(Q)); |
---|
295 | |
---|
296 | for(int i = 1; i <= size(RL); i++){ |
---|
297 | // shift RL[i] to deg w by multiplication with T^delta |
---|
298 | intvec delta = reduceIntvec(w - multiDeg(RL[i]), TOR); |
---|
299 | |
---|
300 | if(delta == intvec(0:size(w))){ |
---|
301 | RLw[size(RLw)+1] = RL[i]; |
---|
302 | } else { |
---|
303 | // if delta is not in cQ, then there is no monomial of this degree: |
---|
304 | if(containsInSupport(cQ, delta)){ |
---|
305 | // calculate all monomials of degree delta |
---|
306 | list deltaMON = wmonomials(delta, 0, TOR); |
---|
307 | |
---|
308 | for(int k = 1; k <= size(deltaMON); k++){ |
---|
309 | poly m = deltaMON[k]; |
---|
310 | RLw[size(RLw)+1] = RL[i] * m; |
---|
311 | |
---|
312 | kill m; |
---|
313 | } |
---|
314 | |
---|
315 | kill k; |
---|
316 | } |
---|
317 | } |
---|
318 | |
---|
319 | kill delta; |
---|
320 | } |
---|
321 | |
---|
322 | return(RLw); |
---|
323 | } |
---|
324 | example |
---|
325 | { |
---|
326 | echo=2; |
---|
327 | |
---|
328 | ring R = 0, T(1..4), dp; |
---|
329 | |
---|
330 | // Weights of variables |
---|
331 | intmat Q[2][4] = |
---|
332 | -2, 2, -1, 1, |
---|
333 | 1, 1, 1, 1; |
---|
334 | |
---|
335 | // attach degree matrix Q to R: |
---|
336 | setBaseMultigrading(Q); |
---|
337 | |
---|
338 | ideal I = T(1)*T(2) + T(3)*T(4); |
---|
339 | |
---|
340 | intvec w = 0,4; |
---|
341 | gradedCompIdealhelper(w, Q, I); |
---|
342 | } |
---|
343 | |
---|
344 | /////////////////////////////////// |
---|
345 | |
---|
346 | static proc wmonomials(intvec w1, int dd, list #) |
---|
347 | "USAGE: wmonomials(w, dd): w is an intvec and dd either 1 or 0; 1 means that the function should compute a list of exponent vectors and 0 means the function should compute a list of monomials. As an optional third argument a list of integers may be provided to indicate that the last rows of the degree matrix are residue classes in ZZ/aZZ. For instance [2,3] means that the second to last row consists of elements in ZZ/2ZZ and the last one of elements in ZZ/3ZZ. |
---|
348 | ASSUME: the basering must be graded (see setBAseMultigrading()) and the cone over |
---|
349 | the degrees of the variables must be pointed; there mustn't be 0-degrees. |
---|
350 | PURPOSE: computes all monomials T^\nu with deg(T^\nu) = w where |
---|
351 | the degree of the i-th variable of the basering is assumed to have |
---|
352 | RETURN: a list of monomials or a list of their exponent vectors if a second argument is provided (an integer). |
---|
353 | EXAMPLE: shows an example." |
---|
354 | { |
---|
355 | intmat Qfull = getVariableWeights(); |
---|
356 | |
---|
357 | // check whether a torsion list is given: |
---|
358 | list TOR; |
---|
359 | |
---|
360 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
361 | TOR = #[1]; |
---|
362 | |
---|
363 | // cut off the degree matrix (i.e., the free part): |
---|
364 | intmat Q[nrows(Qfull) - size(TOR)][ncols(Qfull)]; |
---|
365 | |
---|
366 | for(int i = 1; i <= nrows(Q); i++){ |
---|
367 | Q[i, 1..ncols(Q)] = Qfull[i, 1..ncols(Q)]; |
---|
368 | } |
---|
369 | |
---|
370 | // cut off the free part w0 of w and store as bigintmat |
---|
371 | // for easier multiplication: |
---|
372 | intvec w0 = w1[1..size(w1)-size(TOR)]; |
---|
373 | bigintmat w[1][size(w0)] = w0[1..size(w1)]; |
---|
374 | kill i; |
---|
375 | |
---|
376 | } else { |
---|
377 | // deg(T_i) =... |
---|
378 | intmat Q = Qfull; |
---|
379 | |
---|
380 | // for easier multiplication: |
---|
381 | intvec w0 = w1; |
---|
382 | bigintmat w[1][size(w0)] = w0[1..size(w0)]; |
---|
383 | } |
---|
384 | |
---|
385 | // look for a bound for monomials: |
---|
386 | cone c = coneViaPoints(transpose(Q)); |
---|
387 | cone cd = dualCone(c); |
---|
388 | bigintmat v = transpose(relativeInteriorPoint(cd)); // 1 column |
---|
389 | |
---|
390 | // if dotprod(a,w) > dotprod(v,w), then |
---|
391 | // we need not consider T^a any more |
---|
392 | // --> gives a degree bound bnd: |
---|
393 | bigint bnd = (w * v)[1,1]; |
---|
394 | |
---|
395 | // traverse all possible exponents v_i of T_i: |
---|
396 | // i-th entry of BNDlist will be maximum exponent for T_i |
---|
397 | // such that the bound can still be satisfied. |
---|
398 | int r = nvars(basering); |
---|
399 | intvec BNDlist = 0:r; |
---|
400 | |
---|
401 | for(int i = 1; i <= r; i++){ |
---|
402 | int bndi = 1; |
---|
403 | intvec coli = col(Q, i); |
---|
404 | bigintmat vi[1][size(coli)] = coli[1..nrows(Q)]; |
---|
405 | |
---|
406 | while(((bndi * vi) * v)[1,1] < bnd){ |
---|
407 | bndi++; |
---|
408 | } |
---|
409 | |
---|
410 | BNDlist[i] = bndi; |
---|
411 | |
---|
412 | kill bndi, vi, coli; |
---|
413 | } |
---|
414 | |
---|
415 | // build up solutions recursively: |
---|
416 | // traverse all possible exponents |
---|
417 | // that still satisfy the bound condition |
---|
418 | // |
---|
419 | // NOTE: this is rather brute-force, |
---|
420 | // but its running time is dominated by the |
---|
421 | // other steps. |
---|
422 | list sol; |
---|
423 | list WL0 = wmonomialsRec(v, bnd, Q, BNDlist, sol); // candidates |
---|
424 | |
---|
425 | // - in case of no torsion: |
---|
426 | // select those ww out of WL0 with Q*ww = w0 |
---|
427 | // - if case of torsion: |
---|
428 | // select those ww out of WL0 with reduce(Q*ww) = w1 |
---|
429 | // where reduce() means to take the last elements mod TOR[i]: |
---|
430 | list WL; |
---|
431 | |
---|
432 | for(int m = 1; m <= size(WL0); m++){ |
---|
433 | list v0 = WL0[m]; |
---|
434 | intvec vv = v0[1..size(v0)]; |
---|
435 | //bigintmat vv[1][size(v0)] = v0[1..size(v0)]; |
---|
436 | |
---|
437 | if(size(TOR) == 0){ // free case |
---|
438 | if(Q * vv == w0){ |
---|
439 | WL[size(WL) + 1] = vv; |
---|
440 | } |
---|
441 | } else { // torsion case |
---|
442 | // reduce the last entries mod p: |
---|
443 | intvec Qvvred = Qfull * vv; |
---|
444 | intvec w1red = w1; |
---|
445 | int pos; |
---|
446 | |
---|
447 | for(int l = 1; l <= size(TOR); l++){ |
---|
448 | pos = nrows(Qfull) - size(TOR) + l; |
---|
449 | Qvvred[pos] = Qvvred[pos] mod TOR[l]; |
---|
450 | w1red[pos] = w1red[pos] mod TOR[l]; |
---|
451 | |
---|
452 | // make positive: e.g. -1 becomes 4 in ZZ/5ZZ: |
---|
453 | if(Qvvred[pos] < 0){ |
---|
454 | Qvvred[pos] = Qvvred[pos] + TOR[l]; |
---|
455 | } |
---|
456 | if(w1red[pos] < 0){ |
---|
457 | w1red[pos] = w1red[pos] + TOR[l]; |
---|
458 | } |
---|
459 | } |
---|
460 | |
---|
461 | if(Qvvred == w1red){ |
---|
462 | WL[size(WL) + 1] = vv; |
---|
463 | } |
---|
464 | |
---|
465 | kill pos, l, Qvvred, w1red; |
---|
466 | } |
---|
467 | |
---|
468 | kill vv, v0; |
---|
469 | } |
---|
470 | |
---|
471 | // output the list of exponentvectors if 2nd argument provided: |
---|
472 | if(dd == 1){ |
---|
473 | return(WL); |
---|
474 | } |
---|
475 | |
---|
476 | // otherwise: as monomials: |
---|
477 | list MONs; |
---|
478 | |
---|
479 | for(i = 1; i <= size(WL); i++){ |
---|
480 | MONs[size(MONs)+1] = vec2monomial(WL[i]); |
---|
481 | } |
---|
482 | |
---|
483 | return(MONs); |
---|
484 | } |
---|
485 | example |
---|
486 | { |
---|
487 | echo=2; |
---|
488 | |
---|
489 | // example with torsion: ZZ + ZZ/5ZZ: |
---|
490 | intmat Q[2][5] = |
---|
491 | 1,1,1,1,1, |
---|
492 | 2,3,1,4,0; |
---|
493 | |
---|
494 | list TOR = 5; |
---|
495 | ring R = 0,T(1..5),dp; |
---|
496 | |
---|
497 | // attach degree matrix Q to R: |
---|
498 | setBaseMultigrading(Q); |
---|
499 | |
---|
500 | ideal I = T(1)*T(2) + T(3)*T(4) + T(5)^2; |
---|
501 | intvec w = 2,0; |
---|
502 | wmonomials(w, 0, TOR); |
---|
503 | |
---|
504 | kill R, w, Q; |
---|
505 | |
---|
506 | ////////////////////// |
---|
507 | // free example: |
---|
508 | ring R = 0, T(1..4), dp; |
---|
509 | |
---|
510 | // Weights of variables |
---|
511 | intmat Q[2][4] = |
---|
512 | -2, 2, -1, 1, |
---|
513 | 1, 1, 1, 1; |
---|
514 | |
---|
515 | // attach degree matrix Q to R: |
---|
516 | setBaseMultigrading(Q); |
---|
517 | |
---|
518 | intvec w = 0,2; |
---|
519 | |
---|
520 | // as list of exponent vectors |
---|
521 | wmonomials(w, 1); |
---|
522 | |
---|
523 | // as list of monomials |
---|
524 | wmonomials(w, 0); |
---|
525 | } |
---|
526 | |
---|
527 | /////////////////////// |
---|
528 | |
---|
529 | // helper: |
---|
530 | // returns all vectors (without looking at w0, Q) with i-th entry at most BND[i] |
---|
531 | // but only those v0 with <v,v0> <= bnd for all i |
---|
532 | static proc wmonomialsRec(bigintmat v, bigint bnd, intmat Q, intvec BNDs, list sol){ |
---|
533 | int pos = size(sol) + 1; |
---|
534 | list Res; |
---|
535 | |
---|
536 | if(size(sol) == nvars(basering)){ |
---|
537 | return(list(sol)); |
---|
538 | } |
---|
539 | |
---|
540 | for(int i = 0; i <= BNDs[pos]; i++){ |
---|
541 | list sol0 = sol; |
---|
542 | sol0[pos] = i; |
---|
543 | |
---|
544 | list sol0tmp = sol; |
---|
545 | |
---|
546 | for(int k = 1; k <= nvars(basering)-size(sol); k++){ |
---|
547 | sol0tmp[size(sol)+k] = 0; |
---|
548 | } |
---|
549 | |
---|
550 | intvec v0tmp = sol0tmp[1 .. size(sol0tmp)]; |
---|
551 | intvec ivQv0tmp = Q*v0tmp; |
---|
552 | bigintmat Qv0tmp[1][size(ivQv0tmp)] = ivQv0tmp[1..size(ivQv0tmp)]; |
---|
553 | |
---|
554 | // this monmial will be too big to have degree w0: |
---|
555 | if((Qv0tmp * v)[1,1] > bnd){ |
---|
556 | return(list(sol0tmp)); |
---|
557 | } |
---|
558 | |
---|
559 | // otherwise: recurse |
---|
560 | list Res0 = wmonomialsRec(v, bnd, Q, BNDs, sol0); |
---|
561 | |
---|
562 | for(k = 1; k <= size(Res0); k++){ |
---|
563 | Res[size(Res)+1] = Res0[k]; |
---|
564 | } |
---|
565 | |
---|
566 | kill k, sol0, Res0, sol0tmp, v0tmp, ivQv0tmp, Qv0tmp; |
---|
567 | } |
---|
568 | |
---|
569 | return(Res); |
---|
570 | } |
---|
571 | |
---|
572 | ////////////////////////// |
---|
573 | |
---|
574 | // helper |
---|
575 | static proc vec2monomial(intvec v){ |
---|
576 | poly f = 1; |
---|
577 | |
---|
578 | for(int i = 1; i <= size(v); i++){ |
---|
579 | f = f * var(i)^v[i]; |
---|
580 | } |
---|
581 | |
---|
582 | return(f); |
---|
583 | } |
---|
584 | |
---|
585 | ////////////////////////// |
---|
586 | |
---|
587 | static proc moncoef(poly f, poly mon) |
---|
588 | "USAGE: moncoef(f, mon): f is a polynomial, mon a monomial. |
---|
589 | PURPOSE: returns the (constant) coefficient of mon in f. |
---|
590 | NOTE: the coefficient of T(1)^2 in T(1)^2*T(2) is assumed to be 1 not T(2). See moncoef2 of the latter is wanted. |
---|
591 | RETURN: a ring element. |
---|
592 | EXAMPLE: shows an example." |
---|
593 | { |
---|
594 | |
---|
595 | poly p = f; |
---|
596 | |
---|
597 | while (p != 0) { |
---|
598 | if( leadmonom(p) / mon != 0){ |
---|
599 | return(leadcoef(p)); |
---|
600 | } |
---|
601 | |
---|
602 | p = p - lead(p); |
---|
603 | } |
---|
604 | |
---|
605 | return(0); |
---|
606 | } |
---|
607 | example |
---|
608 | { |
---|
609 | echo=2; |
---|
610 | |
---|
611 | ring R = (0,a),T(1..4),dp; |
---|
612 | |
---|
613 | poly f = T(1)^2*T(2)^3 + 7*a*T(3)^3 -8*T(3)^2 +7; |
---|
614 | poly mon = T(3)^3; |
---|
615 | |
---|
616 | moncoef(f, mon); |
---|
617 | kill f, mon, R; |
---|
618 | |
---|
619 | // 2nd example |
---|
620 | ring R = 0,(T(1..3), Y(1..12)),dp; |
---|
621 | poly f = T(1)^2*Y(1)^2; |
---|
622 | poly m = T(1)^2; |
---|
623 | moncoef(f,m); |
---|
624 | } |
---|
625 | |
---|
626 | ///////////////////////////////// |
---|
627 | |
---|
628 | // helper |
---|
629 | // for autGenWeights, the case of torsion. |
---|
630 | static proc autGenWeightsTorsion(intmat Q, list TOR, list Origs, list OrigDim){ |
---|
631 | // cut off the torsion rows: |
---|
632 | intmat Q0[nrows(Q)-size(TOR)][ncols(Q)]; |
---|
633 | |
---|
634 | for(int i = 1; i <= nrows(Q) - size(TOR); i++){ |
---|
635 | Q0[i, 1..ncols(Q)] = Q[i, 1..ncols(Q)]; |
---|
636 | } |
---|
637 | |
---|
638 | // compute automorphisms B in AUT0 fixing the free part of |
---|
639 | // the generator weights: |
---|
640 | list AUT0 = autGenWeights(Q0); |
---|
641 | list AUT; |
---|
642 | |
---|
643 | // each element of AUT should be of shape |
---|
644 | // |
---|
645 | // B | 0 |
---|
646 | // C | D |
---|
647 | // |
---|
648 | // where D = diag(d1,..,dn) and C has entries |
---|
649 | // 0 <= c_ij < TOR[i] |
---|
650 | |
---|
651 | // form all possible C-matrices and D-matrices: |
---|
652 | list DL = getDmats(TOR); |
---|
653 | list CL = getCmats(TOR, nrows(Q0)); |
---|
654 | |
---|
655 | int j; |
---|
656 | int k; |
---|
657 | |
---|
658 | for(i = 1; i <= size(AUT0); i++){ |
---|
659 | for(j = 1; j <= size(CL); j++){ |
---|
660 | for(k = 1; k <= size(DL); k++){ |
---|
661 | intmat A = concatBCD(AUT0[i], CL[j], DL[k]); |
---|
662 | |
---|
663 | // M * Origs = Origs?: |
---|
664 | if(fixesOrig(A, Origs, OrigDim, TOR)){ |
---|
665 | AUT[size(AUT)+1] = A; |
---|
666 | } |
---|
667 | |
---|
668 | kill A; |
---|
669 | } |
---|
670 | } |
---|
671 | } |
---|
672 | |
---|
673 | return(AUT); |
---|
674 | } |
---|
675 | |
---|
676 | ///////////////////////////////// |
---|
677 | |
---|
678 | // helper |
---|
679 | // returns all lattice bases among QQ of size size(QQ[1]) |
---|
680 | static proc allLatticeBases(list QQ){ |
---|
681 | int d = size(QQ[1]); |
---|
682 | int r = size(QQ); |
---|
683 | list Bases; |
---|
684 | int j; |
---|
685 | |
---|
686 | // compute all lattice bases of size d among QQ: |
---|
687 | // traverse subsets (given as binary vectors J): |
---|
688 | for(int i = 1; i <= 2^r; i++){ |
---|
689 | intvec J = int2face(i, r); |
---|
690 | list QJ; |
---|
691 | |
---|
692 | // status info since this may run for a longer time: |
---|
693 | if(i mod 300 == 0){ |
---|
694 | dbprint(printlevel, "i=" + string(i) + " of " + string(2^r)); |
---|
695 | } |
---|
696 | |
---|
697 | // translate J to a 'subset' QJ (given as a list) of QQ: |
---|
698 | for(j = 1; j <= size(J); j++){ |
---|
699 | if(J[j] == 1){ |
---|
700 | QJ[size(QJ) + 1] = QQ[j]; |
---|
701 | } |
---|
702 | } |
---|
703 | |
---|
704 | // sizes are invariant --> look for size d: |
---|
705 | if(size(QJ) == d){ |
---|
706 | intmat QJmat[size(QQ[1])][d]; |
---|
707 | |
---|
708 | for(int k = 1; k <= d; k++){ |
---|
709 | QJmat[1..size(QQ[1]), k] = QJ[k]; |
---|
710 | } |
---|
711 | if(absValue(det(QJmat)) == 1){ |
---|
712 | Bases[size(Bases)+1] = QJmat; |
---|
713 | } |
---|
714 | |
---|
715 | kill QJmat, k; |
---|
716 | } |
---|
717 | |
---|
718 | kill QJ, J; |
---|
719 | } |
---|
720 | |
---|
721 | if(size(Bases) < 1){ |
---|
722 | // for only one row, it is ok to return simply [1]: |
---|
723 | if(d == 1){ |
---|
724 | intmat QJmat[1][1] = 1; |
---|
725 | Bases[1] = QJmat; |
---|
726 | return(Bases); |
---|
727 | } |
---|
728 | |
---|
729 | ERROR("Sorry. Could not find a basis in the columns."); |
---|
730 | } |
---|
731 | |
---|
732 | return(Bases); |
---|
733 | } |
---|
734 | |
---|
735 | //////////////////////////////// |
---|
736 | |
---|
737 | // helper |
---|
738 | // for autGenWeights, the case of no torsion. |
---|
739 | static proc autGenWeightsFree(list QQ, list Origs, list OrigDim) |
---|
740 | { |
---|
741 | dbprint(printlevel, "autGenWeights: entering free case..."); |
---|
742 | list Bases = allLatticeBases(QQ); |
---|
743 | |
---|
744 | dbprint(printlevel, "Number of bases: " + string(size(Bases))); |
---|
745 | |
---|
746 | // take only those A, with A*Origs = Origs: |
---|
747 | list G; |
---|
748 | |
---|
749 | // Any basis is mapped to another basis, |
---|
750 | // so it suffices to fix one basis A and |
---|
751 | // see where it is mapped to. |
---|
752 | intmat A = Bases[1]; |
---|
753 | list AA; |
---|
754 | |
---|
755 | for(int j = 1; j <= ncols(A); j++){ |
---|
756 | AA[size(AA)+1] = col(A, j); |
---|
757 | } |
---|
758 | |
---|
759 | // does not work for 1x1: |
---|
760 | // workaround: |
---|
761 | if(nrows(A) == 1){ |
---|
762 | if(A[1,1] == 1){ |
---|
763 | intmat Ainv = A; |
---|
764 | } else{ // = -1 |
---|
765 | intmat Ainv = -A; |
---|
766 | } |
---|
767 | } else { |
---|
768 | intmat Ainv = intInverse(A); |
---|
769 | } |
---|
770 | |
---|
771 | // we also have to consider permutations of each |
---|
772 | // "image basis": |
---|
773 | int d = size(QQ[1]); |
---|
774 | intvec vd = 1..d; |
---|
775 | list Ld = vd[1..size(vd)]; |
---|
776 | list PER = permutations(Ld, list()); |
---|
777 | |
---|
778 | dbprint(printlevel, "Number of permutations: " + string(size(PER))); |
---|
779 | |
---|
780 | int b; |
---|
781 | list RES; |
---|
782 | |
---|
783 | for(j = 1; j <= size(Bases); j++){ |
---|
784 | dbprint(printlevel, "Checking Basis no. " + string(j) + " of " + string(size(Bases))); |
---|
785 | |
---|
786 | intmat B0 = Bases[j]; // note: d == ncols(B0) >= nrows(B0) |
---|
787 | |
---|
788 | for(b = 1; b <= size(PER); b++){ |
---|
789 | intmat B[nrows(B0)][ncols(B0)]; |
---|
790 | list sigma = PER[b]; |
---|
791 | |
---|
792 | for(int k = 1; k <= d; k++){ |
---|
793 | // the k-th col of B is the |
---|
794 | // sigma(k)-th one of B0: |
---|
795 | intvec v = col(B0, sigma[k]); |
---|
796 | B[1..nrows(B0), k] = v[1..size(v)]; |
---|
797 | |
---|
798 | kill v; |
---|
799 | } |
---|
800 | |
---|
801 | // init M: note that det(M) = +-1 |
---|
802 | // since det(B)=+-1 and det(Ainv)=+-1. |
---|
803 | intmat M = B * Ainv; |
---|
804 | |
---|
805 | // M * Origs = Origs?: |
---|
806 | int takeit = fixesOrig(M, Origs, OrigDim, list()); |
---|
807 | |
---|
808 | if(takeit){ |
---|
809 | // M already present?: |
---|
810 | for(int u = 1; u <= size(RES); u++){ |
---|
811 | if(RES[u] == M){ |
---|
812 | takeit = 0; |
---|
813 | break; |
---|
814 | } |
---|
815 | } |
---|
816 | |
---|
817 | if(takeit){ |
---|
818 | RES[size(RES)+1] = M; |
---|
819 | } |
---|
820 | |
---|
821 | kill u; |
---|
822 | } |
---|
823 | |
---|
824 | kill k, sigma, B, M, takeit; |
---|
825 | } |
---|
826 | |
---|
827 | kill B0; |
---|
828 | } |
---|
829 | |
---|
830 | return(RES); |
---|
831 | } |
---|
832 | |
---|
833 | ///////////////////////////////// |
---|
834 | |
---|
835 | proc autGenWeights(intmat Q, list #) |
---|
836 | "USAGE: autGenWeights(Q): Q is an intmat (columns must contain a lattice basis). |
---|
837 | ASSUME: the cone over Q must be pointed and the columns of Q contain a lattice basis; there must be no 0-columns in Q. We assume that, in the torsion case, the torsion rows of Q are reduced (for example, a row of Q standing for entries in ZZ/5ZZ must not contain elements > 5 or < 0). |
---|
838 | PURPOSE: computes generators for the subgroup aut(Omega_S) of GL(n, \ZZ) that consists of all invertible integer kxk matrices which fix the set Omega_S of degrees of the variables of the basering S. The set of columns of Q equals Omega_S. |
---|
839 | REFERENCE: Remark 3.1. |
---|
840 | RETURN: a list of integral matrices A with |det A| = 1 such that A*{columns of Q} = {columns of Q}. |
---|
841 | EXAMPLE:" |
---|
842 | { |
---|
843 | dbprint(printlevel, "Entering autGenWeights."); |
---|
844 | |
---|
845 | // if S is the basering, and w1, w2 are weights of some generators |
---|
846 | // the component S_w1 can only be mapped to S_w2 if the dimensions |
---|
847 | // coincide. --> store these in OrigDim. |
---|
848 | // Store in Origs the generator weights without repetition (also known as Omega_S). |
---|
849 | list Origs = Q2Origs(Q); |
---|
850 | list OrigDim = Origs[2]; |
---|
851 | list QQ = Origs[3]; // store columns of Q here |
---|
852 | Origs = Origs[1]; |
---|
853 | |
---|
854 | dbprint(printlevel, "--- Originary degrees ---"); |
---|
855 | dbprint(printlevel, Origs); |
---|
856 | dbprint(printlevel, "--- Dimension of the respective component ---"); |
---|
857 | dbprint(printlevel, OrigDim); |
---|
858 | |
---|
859 | // distinguish between torsion and free cases. |
---|
860 | list TOR; |
---|
861 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
862 | // this treats the torsion case |
---|
863 | dbprint(printlevel, "autGenWeights: entering torsion case..."); |
---|
864 | |
---|
865 | TOR = #[1]; |
---|
866 | |
---|
867 | return(autGenWeightsTorsion(Q, TOR, Origs, OrigDim)); |
---|
868 | } |
---|
869 | |
---|
870 | // now, we treat the free case: |
---|
871 | return(autGenWeightsFree(QQ, Origs, OrigDim)); |
---|
872 | } |
---|
873 | example |
---|
874 | { |
---|
875 | echo=2; |
---|
876 | |
---|
877 | // torsion example |
---|
878 | // ZZ + ZZ/5ZZ: |
---|
879 | intmat Q[2][5] = |
---|
880 | 1,1,1,1,1, |
---|
881 | 2,3,1,4,0; |
---|
882 | |
---|
883 | list TOR = 5; |
---|
884 | |
---|
885 | autGenWeights(Q, TOR); |
---|
886 | kill Q, TOR; |
---|
887 | |
---|
888 | // another free example |
---|
889 | intmat Q[2][6] = |
---|
890 | -2,2,-1,1,-1,1, |
---|
891 | 1,1,1,1,1,1; |
---|
892 | |
---|
893 | autGenWeights(Q); |
---|
894 | kill Q; |
---|
895 | |
---|
896 | //---------------- |
---|
897 | // 2nd free example |
---|
898 | intmat Q[2][4] = |
---|
899 | 1,0,1,1, |
---|
900 | 0,1,1,1; |
---|
901 | |
---|
902 | autGenWeights(Q); |
---|
903 | kill Q; |
---|
904 | } |
---|
905 | |
---|
906 | /////////////////////////////////// |
---|
907 | |
---|
908 | // converts e.g. n=5 to its binary representation, i.e. 0,1,0,1 |
---|
909 | // if r = 3. |
---|
910 | // and stores it in an intvec. |
---|
911 | // r gives the bound for n <= 2^r: |
---|
912 | static proc int2face(int n, int r) |
---|
913 | { |
---|
914 | int k = r; |
---|
915 | list v; |
---|
916 | int n0 = n; |
---|
917 | |
---|
918 | while(k >= 0){ |
---|
919 | if(2^k > n0){ |
---|
920 | v[size(v)+1] = 0; |
---|
921 | } else { |
---|
922 | v[size(v)+1] = 1; |
---|
923 | n0 = n0 - 2^k; |
---|
924 | } |
---|
925 | |
---|
926 | k--; |
---|
927 | } |
---|
928 | |
---|
929 | if(size(v)>=2){ |
---|
930 | v = v[2..size(v)]; |
---|
931 | } |
---|
932 | |
---|
933 | return(intvec(v[1..size(v)])); |
---|
934 | } |
---|
935 | example |
---|
936 | { |
---|
937 | echo = 2; |
---|
938 | int n = 5; |
---|
939 | int r = 4; |
---|
940 | int2face(n, r); |
---|
941 | |
---|
942 | n = 1; |
---|
943 | r = 1; |
---|
944 | int2face(n,r); |
---|
945 | } |
---|
946 | |
---|
947 | ////////////////////////////// |
---|
948 | |
---|
949 | // helper |
---|
950 | // returns all permutations of the input list L. |
---|
951 | static proc permutations(list Left, list sol){ |
---|
952 | list Res; |
---|
953 | |
---|
954 | if(size(Left) == 0){ |
---|
955 | return(list(sol)); |
---|
956 | } |
---|
957 | |
---|
958 | for(int i = 1; i <= size(Left); i++){ |
---|
959 | list sol0 = sol; |
---|
960 | sol0[size(sol0)+1] = Left[i]; |
---|
961 | |
---|
962 | if(i > 1){ |
---|
963 | if(i < size(Left)){ |
---|
964 | list Left0 = Left[1..i-1], Left[i+1..size(Left)]; |
---|
965 | } else { // i = size(Left) |
---|
966 | list Left0 = Left[1..i-1]; |
---|
967 | } |
---|
968 | } else{ // i = 1 |
---|
969 | if(1 == size(Left)){ |
---|
970 | list Left0; |
---|
971 | } else { |
---|
972 | list Left0 = Left[i+1..size(Left)]; |
---|
973 | } |
---|
974 | } |
---|
975 | |
---|
976 | list Res0 = permutations(Left0, sol0); |
---|
977 | |
---|
978 | for(int k = 1; k <= size(Res0); k++){ |
---|
979 | Res[size(Res)+1] = Res0[k]; |
---|
980 | } |
---|
981 | |
---|
982 | kill k, sol0, Res0, Left0; |
---|
983 | } |
---|
984 | |
---|
985 | return(Res); |
---|
986 | } |
---|
987 | example |
---|
988 | { |
---|
989 | echo=2; |
---|
990 | |
---|
991 | list L = 1,2,3; |
---|
992 | permutations(L, list()); |
---|
993 | |
---|
994 | kill L; |
---|
995 | } |
---|
996 | |
---|
997 | //////////////////////////// |
---|
998 | |
---|
999 | static proc fixesOrig(intmat M, list Origs, list OrigDim, list TOR){ |
---|
1000 | int takeit = 1; |
---|
1001 | intvec found = 0:size(Origs); |
---|
1002 | |
---|
1003 | for(int l = 1; l <= size(Origs); l++){ |
---|
1004 | intmat W = Origs[l]; // 1 col |
---|
1005 | intmat MW = M * W; // 1 col |
---|
1006 | |
---|
1007 | // test whether M*W = some orig |
---|
1008 | // (and also not yet found)?: |
---|
1009 | for(int m = 1; m <= size(Origs); m++){ |
---|
1010 | // reduce torsion rows: |
---|
1011 | intmat MWred = MW; // 1 col |
---|
1012 | |
---|
1013 | for(int i = 1; i <= size(TOR); i++){ |
---|
1014 | MWred[nrows(MWred) - size(TOR) + i,1] = MWred[nrows(MWred) - size(TOR) + i,1] mod TOR[i]; |
---|
1015 | |
---|
1016 | if(MWred[nrows(MWred) - size(TOR) + i,1] < 0){ |
---|
1017 | MWred[nrows(MWred) - size(TOR) + i,1] = MWred[nrows(MWred) - size(TOR) + i,1] + TOR[i]; |
---|
1018 | } |
---|
1019 | } |
---|
1020 | |
---|
1021 | if(Origs[m] == MWred){ |
---|
1022 | if(found[m] == OrigDim[m]){ |
---|
1023 | takeit = 0; |
---|
1024 | break; |
---|
1025 | } |
---|
1026 | |
---|
1027 | found[m] = 1; |
---|
1028 | } |
---|
1029 | |
---|
1030 | kill MWred, i; |
---|
1031 | } |
---|
1032 | |
---|
1033 | kill m, W, MW; |
---|
1034 | |
---|
1035 | if(!takeit){ |
---|
1036 | return(0); // false |
---|
1037 | } |
---|
1038 | } |
---|
1039 | |
---|
1040 | if(takeit and found == 1:size(Origs)){ |
---|
1041 | return(1); // true |
---|
1042 | } |
---|
1043 | |
---|
1044 | return(0); // false |
---|
1045 | } |
---|
1046 | |
---|
1047 | ///////////////////////////////////// |
---|
1048 | |
---|
1049 | // helper: |
---|
1050 | // returns a list of all matrices of shape |
---|
1051 | // diag(d_1,...,d_size(TOR)) |
---|
1052 | // where gcd(d_i, TOR[i]) = 1. |
---|
1053 | static proc getDmats(list TOR){ |
---|
1054 | list sol; |
---|
1055 | return(getDmatsrec(TOR, sol)); |
---|
1056 | } |
---|
1057 | example |
---|
1058 | { |
---|
1059 | list TOR = 3,5; |
---|
1060 | TOR; |
---|
1061 | getDmats(TOR); |
---|
1062 | } |
---|
1063 | |
---|
1064 | /////////////////////////// |
---|
1065 | |
---|
1066 | // helper |
---|
1067 | // for getDmats |
---|
1068 | static proc getDmatsrec(list TOR, list sol){ |
---|
1069 | int pos = size(sol) + 1; |
---|
1070 | list Res; |
---|
1071 | |
---|
1072 | if(size(sol) == size(TOR)){ |
---|
1073 | intmat D[size(TOR)][size(TOR)]; |
---|
1074 | |
---|
1075 | for(int k = 1; k <= size(TOR); k++){ |
---|
1076 | list L0 = sol[k]; |
---|
1077 | D[k, 1..ncols(D)] = L0[1..size(L0)]; |
---|
1078 | |
---|
1079 | kill L0; |
---|
1080 | } |
---|
1081 | |
---|
1082 | return(list(D)); |
---|
1083 | } |
---|
1084 | |
---|
1085 | for(int i = 1; i < TOR[pos]; i++){ |
---|
1086 | list sol0 = sol; |
---|
1087 | |
---|
1088 | if(gcd(i, TOR[pos]) == 1){ |
---|
1089 | sol0[pos] = i; |
---|
1090 | |
---|
1091 | list Res0 = getDmatsrec(TOR, sol0); |
---|
1092 | for(int k = 1; k <= size(Res0); k++){ |
---|
1093 | Res[size(Res)+1] = Res0[k]; |
---|
1094 | } |
---|
1095 | |
---|
1096 | kill k; |
---|
1097 | kill Res0; |
---|
1098 | } |
---|
1099 | |
---|
1100 | kill sol0; |
---|
1101 | } |
---|
1102 | |
---|
1103 | return(Res); |
---|
1104 | } |
---|
1105 | |
---|
1106 | /////////////////////////// |
---|
1107 | |
---|
1108 | // helper: returns all intvecs |
---|
1109 | // of length n with entries 0 <= c < m: |
---|
1110 | static proc getCvecsrec(int m, int n, list sol){ |
---|
1111 | int pos = size(sol) + 1; |
---|
1112 | list Res; |
---|
1113 | |
---|
1114 | if(size(sol) == n){ |
---|
1115 | return(list(sol)); |
---|
1116 | } |
---|
1117 | |
---|
1118 | for(int i = 0; i < m; i++){ |
---|
1119 | list sol0 = sol; |
---|
1120 | sol0[pos] = i; |
---|
1121 | |
---|
1122 | list Res0 = getCvecsrec(m, n, sol0); |
---|
1123 | for(int k = 1; k <= size(Res0); k++){ |
---|
1124 | Res[size(Res)+1] = Res0[k]; |
---|
1125 | } |
---|
1126 | |
---|
1127 | kill k, Res0, sol0; |
---|
1128 | } |
---|
1129 | |
---|
1130 | return(Res); |
---|
1131 | } |
---|
1132 | |
---|
1133 | ////////////////////////////////////////// |
---|
1134 | |
---|
1135 | // helper: |
---|
1136 | // returns a list of all size(TOR) x n |
---|
1137 | // matrices C with entries |
---|
1138 | // 0 <= c_ij < TOR[i] |
---|
1139 | static proc getCmats(list TOR, int n){ |
---|
1140 | list ROWS; |
---|
1141 | |
---|
1142 | for(int i = 1; i <= size(TOR); i++){ |
---|
1143 | list sol; |
---|
1144 | ROWS[i] = getCvecsrec(TOR[i], n, sol); |
---|
1145 | |
---|
1146 | kill sol; |
---|
1147 | } |
---|
1148 | |
---|
1149 | // choose one from each element of ROWs: |
---|
1150 | // i.e., form the cartesian product: |
---|
1151 | list CL = getCmatsrec(ROWS, list()); |
---|
1152 | |
---|
1153 | return(CL); |
---|
1154 | } |
---|
1155 | example |
---|
1156 | { |
---|
1157 | list TOR = 3,5; |
---|
1158 | |
---|
1159 | TOR; |
---|
1160 | |
---|
1161 | "ZZ^2 + ZZ/3ZZ + ZZ/5ZZ:"; |
---|
1162 | getCmats(TOR, 2); |
---|
1163 | |
---|
1164 | } |
---|
1165 | |
---|
1166 | /////////////////////////// |
---|
1167 | |
---|
1168 | // helper: chooses one element of each |
---|
1169 | // element of ROWs and forms the matrix |
---|
1170 | static proc getCmatsrec(list ROWS, list sol){ |
---|
1171 | int pos = size(sol) + 1; |
---|
1172 | list Res; |
---|
1173 | |
---|
1174 | if(size(sol) == size(ROWS)){ |
---|
1175 | intmat C[size(ROWS)][size(ROWS[1][1])]; |
---|
1176 | |
---|
1177 | for(int k = 1; k <= size(ROWS); k++){ |
---|
1178 | list ROWSk = sol[k]; |
---|
1179 | C[k, 1..ncols(C)] = ROWSk[1..size(ROWSk)]; |
---|
1180 | |
---|
1181 | kill ROWSk; |
---|
1182 | } |
---|
1183 | |
---|
1184 | return(list(C)); //return(list(sol)); |
---|
1185 | } |
---|
1186 | |
---|
1187 | for(int i = 1; i <= size(ROWS[pos]); i++){ |
---|
1188 | list sol0 = sol; |
---|
1189 | sol0[pos] = ROWS[pos][i]; |
---|
1190 | |
---|
1191 | list Res0 = getCmatsrec(ROWS, sol0); |
---|
1192 | for(int k = 1; k <= size(Res0); k++){ |
---|
1193 | Res[size(Res)+1] = Res0[k]; |
---|
1194 | } |
---|
1195 | |
---|
1196 | kill k, Res0, sol0; |
---|
1197 | } |
---|
1198 | |
---|
1199 | return(Res); |
---|
1200 | } |
---|
1201 | |
---|
1202 | ////////////////////////////////////////// |
---|
1203 | |
---|
1204 | // helper: concats the given matrices to |
---|
1205 | // B | 0 |
---|
1206 | // C | D |
---|
1207 | static proc concatBCD(intmat B, intmat C, intmat D){ |
---|
1208 | intmat A[nrows(B) + nrows(C)][ncols(B) + ncols(D)]; |
---|
1209 | |
---|
1210 | for(int i = 1; i <= nrows(B); i++){ |
---|
1211 | A[i, 1..ncols(B)] = B[i,1..ncols(B)]; //, 0:ncols(D); |
---|
1212 | } |
---|
1213 | for(i = nrows(B)+1; i <= nrows(B) + nrows(C); i++){ |
---|
1214 | //A[i, 1..ncols(A)] = C[i - nrows(B), 1..ncols(C)], D[i - nrows(B), 1..ncols(D)]; |
---|
1215 | A[i, 1..ncols(C)] = C[i - nrows(B), 1..ncols(C)]; |
---|
1216 | A[i, ncols(C)+1..ncols(A)] = D[i - nrows(B), 1..ncols(D)]; |
---|
1217 | } |
---|
1218 | |
---|
1219 | return(A); |
---|
1220 | } |
---|
1221 | |
---|
1222 | ///////////////////////////////// |
---|
1223 | |
---|
1224 | // helper |
---|
1225 | // assume: BL of shape [[d, w, Rw],...] where d is the dimension of R_w, |
---|
1226 | // Rw a basis of R_w and w the weight vector. |
---|
1227 | // Returns a ring S and exports a matrix of name Aexported over S |
---|
1228 | // and a list of monomials MONSexported. |
---|
1229 | static proc createAutMatrix(list BL) |
---|
1230 | { |
---|
1231 | // A will be nxn; compute n: |
---|
1232 | list MON = sortMons(BL); |
---|
1233 | int n = size(MON); |
---|
1234 | |
---|
1235 | // the first rows of A correspond to the variables T(i): |
---|
1236 | def R = basering; |
---|
1237 | int n0 = nvars(basering); |
---|
1238 | intmat Q = getVariableWeights(); |
---|
1239 | intmat QS[nrows(Q)][ncols(Q) + n0*n]; |
---|
1240 | |
---|
1241 | for(int l = 1; l <= ncols(Q); l++){ |
---|
1242 | QS[1..nrows(Q),l] = Q[1..nrows(Q),l]; |
---|
1243 | } |
---|
1244 | |
---|
1245 | int ii; |
---|
1246 | list l3; |
---|
1247 | for (ii = 1; ii <= n0; ii++) |
---|
1248 | { |
---|
1249 | l3[ii] = "T("+string(ii)+")"; |
---|
1250 | } |
---|
1251 | for (ii = 1; ii <= n0*n; ii++) |
---|
1252 | { |
---|
1253 | l3[n0+ii] = "Y("+string(ii)+")"; |
---|
1254 | } |
---|
1255 | ring S = create_ring(ring_list(basering)[1], l3, "dp("+string(n0+n0*n)+")", "no_minpoly"); |
---|
1256 | setring S; |
---|
1257 | setBaseMultigrading(QS); |
---|
1258 | |
---|
1259 | matrix A[n][n]; |
---|
1260 | map f = R, T(1..n0); |
---|
1261 | list MONS = f(MON); |
---|
1262 | list BLS = f(BL); |
---|
1263 | int k; |
---|
1264 | |
---|
1265 | for(int i = 1; i <= n0; i++){ |
---|
1266 | // find out j such that BLS[j] contains var(i): |
---|
1267 | int ind = 0; |
---|
1268 | |
---|
1269 | for(int j = 1; j <= size(BLS); j++){ |
---|
1270 | for(k = 1; k <= size(BLS[j][3]); k++){ |
---|
1271 | if(var(i) == BLS[j][3][k]){ |
---|
1272 | ind = j; |
---|
1273 | break; |
---|
1274 | } |
---|
1275 | } |
---|
1276 | |
---|
1277 | if(ind > 0){ |
---|
1278 | break; |
---|
1279 | } |
---|
1280 | } |
---|
1281 | |
---|
1282 | list MONi = BLS[ind][3]; |
---|
1283 | |
---|
1284 | for(j = 1; j <= n; j++){ |
---|
1285 | // test if MON[j] is in the set {MONi}: |
---|
1286 | for(int m = 1; m <= size(MONi); m++){ |
---|
1287 | if(MONi[m] == MONS[j]){ |
---|
1288 | A[i,j] = Y((i-1)*n+j); // the Y-variables are ordered row-wise |
---|
1289 | break; |
---|
1290 | } |
---|
1291 | } |
---|
1292 | |
---|
1293 | kill m; |
---|
1294 | } |
---|
1295 | |
---|
1296 | kill MONi, ind, j; |
---|
1297 | } |
---|
1298 | |
---|
1299 | // Form the A*T[i] first for the lower rows: |
---|
1300 | // store A*T[i] in H[1,i]: |
---|
1301 | matrix H[1][n0]; |
---|
1302 | |
---|
1303 | for(k = 1; k <= n0; k++){ |
---|
1304 | for(i = 1; i <= size(MONS); i++){ |
---|
1305 | H[1,k] = H[1,k] + A[k,i] * MONS[i]; |
---|
1306 | } |
---|
1307 | } |
---|
1308 | |
---|
1309 | // the lower rows of of A |
---|
1310 | // are determined by the first ones: |
---|
1311 | // e.g. A*(T_1^2) = (A*T_1)^2: |
---|
1312 | for(i = n0 + 1; i <= size(MONS); i++){ |
---|
1313 | // e.g. for MON[i] = T[1]*T[2]^2, |
---|
1314 | // we have v = [1,2,0,0] and |
---|
1315 | // A * (T[1]*T[2]^2) = (A*T[1])(A*T[2])^2: |
---|
1316 | poly Av = 1; |
---|
1317 | intvec v = leadexp(MONS[i]); // the last entries are all 0s |
---|
1318 | |
---|
1319 | for(int m = 1; m <= n0; m++){ |
---|
1320 | Av = Av * H[1,m]^v[m]; |
---|
1321 | } |
---|
1322 | |
---|
1323 | // if in Av there is e.g. the monomial Y(22)^2*T(1)^2 |
---|
1324 | // then A[i,j] = Y(22)^2 where MONS[j] = T(1)^2: |
---|
1325 | for(int j = 1; j <= size(MONS); j++){ |
---|
1326 | A[i,j] = moncoef2(Av, MONS[j], n0); |
---|
1327 | } |
---|
1328 | |
---|
1329 | kill Av, v, m, j; |
---|
1330 | } |
---|
1331 | |
---|
1332 | matrix Aexported = A; |
---|
1333 | list MONSexported = MONS; |
---|
1334 | |
---|
1335 | // ring-dependent: |
---|
1336 | export(Aexported); |
---|
1337 | export(MONSexported); |
---|
1338 | |
---|
1339 | return(S); |
---|
1340 | } |
---|
1341 | |
---|
1342 | //////////////////////////// |
---|
1343 | |
---|
1344 | // helper |
---|
1345 | // returns bases for <RL>_w for all w in Orig. |
---|
1346 | static proc allMonomials(ideal RL, list #){ |
---|
1347 | list TOR; |
---|
1348 | |
---|
1349 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
1350 | TOR = #[1]; |
---|
1351 | } |
---|
1352 | |
---|
1353 | list BL; |
---|
1354 | intmat Q = getVariableWeights(); |
---|
1355 | |
---|
1356 | // by assumption: Origs = degrees (no repetition) |
---|
1357 | list Origs = Q2Origs(Q)[1]; |
---|
1358 | |
---|
1359 | for(int i = 1; i <= size(Origs); i++){ |
---|
1360 | intvec w = Origs[i]; |
---|
1361 | list L = gradedCompRing(RL, w, TOR); |
---|
1362 | list Rw = L[1]; |
---|
1363 | |
---|
1364 | BL[size(BL)+1] = list(size(Rw), w, Rw); |
---|
1365 | |
---|
1366 | kill L, w, Rw; |
---|
1367 | } |
---|
1368 | |
---|
1369 | return(BL); |
---|
1370 | } |
---|
1371 | |
---|
1372 | //////////////////////////// |
---|
1373 | |
---|
1374 | // helper |
---|
1375 | // returns the columns of Q without repetition |
---|
1376 | // and the list of dimensions. |
---|
1377 | static proc Q2Origs(intmat Q){ |
---|
1378 | list Origs; |
---|
1379 | list OrigDim; |
---|
1380 | list QQ; |
---|
1381 | |
---|
1382 | for(int i = 1; i <= ncols(Q); i++){ |
---|
1383 | int takeit = 1; |
---|
1384 | intvec v = col(Q, i); |
---|
1385 | QQ[i] = v; |
---|
1386 | |
---|
1387 | for(int j = 1; j <= size(Origs); j++){ |
---|
1388 | if(Origs[j] == v){ |
---|
1389 | OrigDim[j] = OrigDim[j] + 1; // OrigDim[j] is already defined |
---|
1390 | takeit = 0; |
---|
1391 | break; |
---|
1392 | } |
---|
1393 | } |
---|
1394 | |
---|
1395 | if(takeit == 1){ |
---|
1396 | Origs[size(Origs)+1] = v; |
---|
1397 | OrigDim[size(OrigDim)+1] = 1; |
---|
1398 | } |
---|
1399 | |
---|
1400 | kill j, v, takeit; |
---|
1401 | } |
---|
1402 | |
---|
1403 | return(list(Origs, OrigDim, QQ)); |
---|
1404 | } |
---|
1405 | |
---|
1406 | ////////////////////////////////////////// |
---|
1407 | |
---|
1408 | static proc gradedCompRing(ideal RL, intvec w, list #) |
---|
1409 | "USAGE: gradedCompRing(I, w): I is an ideal, w an intvec. Optional 3rd argument: a list of integers representing the torsion part. |
---|
1410 | ASSUME: the basering must be graded (see setBAseMultigrading()) and the cone over |
---|
1411 | the degrees of the variables must be pointed; there mustn't be 0-degrees. The vector |
---|
1412 | w must be an element of the cone over the degrees of the variables. |
---|
1413 | PURPOSE: returns a basis for the vector space R_w where R is the factor ring of the basering by I. |
---|
1414 | RETURN: a list L where L[1] = is the (monomial) Basis, L[2] = is a matrix with columns representing the basis vectors, L[3] = Monomials of degree w; |
---|
1415 | EXAMPLE: shows an example." |
---|
1416 | { |
---|
1417 | list TOR; |
---|
1418 | |
---|
1419 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
1420 | TOR = #[1]; |
---|
1421 | } |
---|
1422 | |
---|
1423 | // return a basis for I_w as a list of basis vectors: |
---|
1424 | list L = gradedCompIdeal(RL, w, 0, TOR); |
---|
1425 | list MON = L[1]; |
---|
1426 | list RLw = L[2]; |
---|
1427 | |
---|
1428 | list Bas; |
---|
1429 | int r; |
---|
1430 | |
---|
1431 | if(size(RLw) == 0){ |
---|
1432 | r = 0; |
---|
1433 | } else { |
---|
1434 | matrix B[size(RLw)][size(MON)]; |
---|
1435 | |
---|
1436 | for(int i = 1; i <= size(RLw); i++){ |
---|
1437 | matrix RLwi = RLw[i]; |
---|
1438 | B[i, 1..size(MON)] = RLwi[1..nrows(RLwi), 1]; |
---|
1439 | |
---|
1440 | kill RLwi; |
---|
1441 | } |
---|
1442 | kill i; |
---|
1443 | |
---|
1444 | r = rank(B); |
---|
1445 | } |
---|
1446 | |
---|
1447 | // Add Elements to Basis 'Bas' as long as the rank r |
---|
1448 | // of the matrix B increases: |
---|
1449 | for(int j = 1; j <= size(MON); j++){ |
---|
1450 | // translate MON[j] to v= -e_j |
---|
1451 | // the new row: e_j, i.e., all zeros (standard) and j-th entry is a 1: |
---|
1452 | if(r == 0){ // this means B is not defined! |
---|
1453 | matrix Btmp[1][size(MON)]; |
---|
1454 | Btmp[1,j] = 1; |
---|
1455 | } else { |
---|
1456 | matrix Btmp[nrows(B)+1][ncols(B)] = B[1..nrows(B), 1..ncols(B)]; |
---|
1457 | Btmp[nrows(B)+1, j] = 1; |
---|
1458 | } |
---|
1459 | |
---|
1460 | int rv = rank(Btmp); |
---|
1461 | |
---|
1462 | if(rv > r){ |
---|
1463 | if(r == 0){ |
---|
1464 | matrix B = Btmp; |
---|
1465 | } else{ |
---|
1466 | B = Btmp; |
---|
1467 | } |
---|
1468 | |
---|
1469 | r = rv; |
---|
1470 | |
---|
1471 | Bas[size(Bas)+1] = MON[j]; |
---|
1472 | } |
---|
1473 | |
---|
1474 | kill rv, Btmp; |
---|
1475 | } |
---|
1476 | |
---|
1477 | // remove columns of RL_w in B: |
---|
1478 | matrix B2[nrows(B) - size(RLw)][ncols(B)] = B[size(RLw)+1..nrows(B), 1..ncols(B)]; |
---|
1479 | |
---|
1480 | list LL; |
---|
1481 | LL[1] = Bas; |
---|
1482 | LL[2] = B2; |
---|
1483 | LL[3] = MON; |
---|
1484 | |
---|
1485 | return(LL); |
---|
1486 | } |
---|
1487 | example |
---|
1488 | { |
---|
1489 | echo=2; |
---|
1490 | |
---|
1491 | // example with torsion: ZZ + ZZ/5ZZ: |
---|
1492 | intmat Q[2][5] = |
---|
1493 | 1,1,1,1,1, |
---|
1494 | 2,3,1,4,0; |
---|
1495 | |
---|
1496 | list TOR = 5; |
---|
1497 | ring R = 0,T(1..5),dp; |
---|
1498 | |
---|
1499 | // attach degree matrix Q to R: |
---|
1500 | setBaseMultigrading(Q); |
---|
1501 | |
---|
1502 | ideal I = T(1)*T(2) + T(3)*T(4) + T(5)^2; |
---|
1503 | intvec w = 3,1; |
---|
1504 | |
---|
1505 | // as list of polynomials: |
---|
1506 | gradedCompRing(I, w, TOR); |
---|
1507 | kill R, w, Q, TOR; |
---|
1508 | |
---|
1509 | // 2nd example: no torsion |
---|
1510 | ring R = 0, T(1..4), dp; |
---|
1511 | |
---|
1512 | // Weights of variables |
---|
1513 | intmat Q[2][4] = |
---|
1514 | -2, 2, -1, 1, |
---|
1515 | 1, 1, 1, 1; |
---|
1516 | |
---|
1517 | // attach degree matrix Q to R: |
---|
1518 | setBaseMultigrading(Q); |
---|
1519 | |
---|
1520 | ideal I = T(1)*T(2) + T(3)*T(4); |
---|
1521 | intvec w = 0,2; |
---|
1522 | |
---|
1523 | gradedCompRing(I, w); |
---|
1524 | } |
---|
1525 | |
---|
1526 | |
---|
1527 | //////////////////////// |
---|
1528 | |
---|
1529 | // helper |
---|
1530 | // sort the monomials in the 3rd entry |
---|
1531 | // of each element of BL into a list |
---|
1532 | // according to the dp-ordering (descendingly) |
---|
1533 | static proc sortMons(list BL){ |
---|
1534 | list MON; |
---|
1535 | |
---|
1536 | // extract monomials |
---|
1537 | for(int i = 1; i <= size(BL); i++){ |
---|
1538 | list Tmp = BL[i][3]; |
---|
1539 | |
---|
1540 | for(int j = 1; j <= size(Tmp); j++){ |
---|
1541 | MON[size(MON)+1] = Tmp[j]; |
---|
1542 | } |
---|
1543 | |
---|
1544 | kill Tmp, j; |
---|
1545 | } |
---|
1546 | |
---|
1547 | // sort MON according to the degrevlex ordering (dp) |
---|
1548 | def R = basering; |
---|
1549 | |
---|
1550 | list l6; |
---|
1551 | for (int zz = 1; zz <= nvars(basering); zz++) |
---|
1552 | { |
---|
1553 | l6[zz] = "Y("+string(nvars(basering)+1-zz)+")"; |
---|
1554 | } |
---|
1555 | ring S = create_ring(ring_list(basering)[1], l6, "dp", "no_minpoly"); |
---|
1556 | setring S; //creates a ring of shape '0, Y(nvars(basering)..1), dp' |
---|
1557 | |
---|
1558 | map f = R, Y(1..nvars(basering)); |
---|
1559 | list fMON = f(MON); |
---|
1560 | list MONS; |
---|
1561 | |
---|
1562 | // to sort MONS: |
---|
1563 | // form Sum_i fMON[i] and look |
---|
1564 | // for the leading terms iteratively: |
---|
1565 | poly p = 0; |
---|
1566 | |
---|
1567 | for(i = 1; i <= size(fMON); i++){ |
---|
1568 | p = p + fMON[i]; |
---|
1569 | } |
---|
1570 | |
---|
1571 | // sort it descendingly: |
---|
1572 | int c = size(fMON); |
---|
1573 | |
---|
1574 | while(p != 0){ |
---|
1575 | MONS[c] = lead(p); |
---|
1576 | p = p - lead(p); |
---|
1577 | c--; |
---|
1578 | } |
---|
1579 | |
---|
1580 | // translate back to R: |
---|
1581 | setring R; |
---|
1582 | |
---|
1583 | ideal mi=maxideal(1); |
---|
1584 | map g = S, mi[nvars(basering)..1]; |
---|
1585 | MON = g(MONS); |
---|
1586 | |
---|
1587 | return(MON); |
---|
1588 | } |
---|
1589 | |
---|
1590 | //////////////////////// |
---|
1591 | |
---|
1592 | // helper: |
---|
1593 | // computes the nxn identity matrix: |
---|
1594 | static proc getIdMat(int n){ |
---|
1595 | intmat A[n][n]; |
---|
1596 | |
---|
1597 | for(int i = 1; i <= n; i++){ |
---|
1598 | A[i, 1..n] = 0:n; |
---|
1599 | A[i,i] = 1; |
---|
1600 | } |
---|
1601 | |
---|
1602 | return(A); |
---|
1603 | } |
---|
1604 | |
---|
1605 | ////////////////////////////////////// |
---|
1606 | |
---|
1607 | // helper |
---|
1608 | static proc insertIntvecIfNew(list L0, intvec v){ |
---|
1609 | list L = L0; |
---|
1610 | int takeit = 1; |
---|
1611 | |
---|
1612 | for(int j = 1; j <= size(L); j++){ |
---|
1613 | if(L[j] == v){ |
---|
1614 | takeit = 0; |
---|
1615 | break; |
---|
1616 | } |
---|
1617 | } |
---|
1618 | |
---|
1619 | if(takeit == 1){ |
---|
1620 | L[size(L) + 1] = v; |
---|
1621 | } |
---|
1622 | |
---|
1623 | return(L); |
---|
1624 | } |
---|
1625 | |
---|
1626 | ////////////////////////////////////// |
---|
1627 | |
---|
1628 | // helper: |
---|
1629 | // returns the set |
---|
1630 | // OmegaI = {w in K; I_w \subseteq I_<w} |
---|
1631 | // of unique minimal generator degrees. |
---|
1632 | // Elements of components are of shape |
---|
1633 | // [w, I_w, K[T]_w] |
---|
1634 | // We assume that if I = <f_1,...,f_s> we have |
---|
1635 | // OmegaI = {deg(f_1),...,deg(f_s)} |
---|
1636 | static proc getOmegaI(ideal I, list #){ |
---|
1637 | dbprint(printlevel, "Warning: we assume that OmegaI = {deg(f_1),...,deg(f_s)} if the input ideal is generated by f_1,...,f_s. "); |
---|
1638 | |
---|
1639 | list TOR; |
---|
1640 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
1641 | TOR = #[1]; |
---|
1642 | } |
---|
1643 | |
---|
1644 | // OmegaI = {w in K; I_w \subseteq I_<w} |
---|
1645 | list OmegaI; |
---|
1646 | |
---|
1647 | // go through each generator f of the ideal: |
---|
1648 | // compute a basis for I_deg(f) |
---|
1649 | for(int i = 1; i <= size(I); i++){ |
---|
1650 | poly f = I[i]; |
---|
1651 | intvec w = multiDegRedTor(f, TOR); |
---|
1652 | |
---|
1653 | OmegaI = insertIntvecIfNew(OmegaI, w); |
---|
1654 | kill f, w; |
---|
1655 | } |
---|
1656 | |
---|
1657 | return(OmegaI); |
---|
1658 | } |
---|
1659 | example |
---|
1660 | { |
---|
1661 | echo=2; |
---|
1662 | |
---|
1663 | intmat Q[1][5] = 3,3,2,2,1; |
---|
1664 | ring R = 0,T(1..5),dp; |
---|
1665 | |
---|
1666 | // attach degree matrix Q to R: |
---|
1667 | setBaseMultigrading(Q); |
---|
1668 | |
---|
1669 | ideal I = T(1)*T(2) + T(3)^2*T(4) + T(5)^6, |
---|
1670 | T(3)*T(4) + T(1)*T(5) + T(5)^4 + T(4)^2; |
---|
1671 | |
---|
1672 | // expect: OmegaI = {4,6} |
---|
1673 | getOmegaI(I); |
---|
1674 | } |
---|
1675 | |
---|
1676 | ///////////////////////////////////// |
---|
1677 | |
---|
1678 | // helper: |
---|
1679 | // write the variables of the basering in the rows of a matrix. |
---|
1680 | // n0 is the number of variables of the very first basering, |
---|
1681 | // i.e., the number of 'meaningful' rows of the resulting matrix. |
---|
1682 | static proc getAgeneral(int n0){ |
---|
1683 | int coldim = (nvars(basering) - 1) div n0; // there is the additional var 'Z' |
---|
1684 | matrix A[coldim][coldim]; |
---|
1685 | |
---|
1686 | // write the variables into A |
---|
1687 | int j; |
---|
1688 | int v = 1; |
---|
1689 | for(int i = 1; i <= n0; i++){ |
---|
1690 | for(j = 1; j <= coldim; j++){ |
---|
1691 | A[i,j] = var(v); |
---|
1692 | v++; |
---|
1693 | } |
---|
1694 | } |
---|
1695 | |
---|
1696 | return(A); |
---|
1697 | } |
---|
1698 | example |
---|
1699 | { |
---|
1700 | echo = 2; |
---|
1701 | ring R = 0,(Y(1..45),Z),dp; |
---|
1702 | matrix A = getAgeneral(5); |
---|
1703 | print(A); |
---|
1704 | } |
---|
1705 | |
---|
1706 | ///////////////////////////////////// |
---|
1707 | |
---|
1708 | // helper |
---|
1709 | // Applies the action given by A. |
---|
1710 | static proc applyA(matrix A, int n0, list AMONS){ |
---|
1711 | // Define the map |
---|
1712 | // T_i --> phi(T_i) = A * T_i: |
---|
1713 | // store A*T[i] in H[1,i]: |
---|
1714 | matrix H[1][n0]; |
---|
1715 | |
---|
1716 | int l; |
---|
1717 | for(int k = 1; k <= n0; k++){ |
---|
1718 | for(l = 1; l <= size(AMONS); l++){ |
---|
1719 | H[1,k] = H[1,k] + A[k,l] * AMONS[l]; |
---|
1720 | } |
---|
1721 | } |
---|
1722 | |
---|
1723 | return(H); |
---|
1724 | } |
---|
1725 | example |
---|
1726 | { |
---|
1727 | echo = 2; |
---|
1728 | ring R = 0,(Y(1..45),Z),dp; |
---|
1729 | matrix AR = getAgeneral(5); |
---|
1730 | |
---|
1731 | ring S = 0,(T(1..5),Y(1..45),Z),dp; |
---|
1732 | matrix A = imap(R, AR); |
---|
1733 | print(A); |
---|
1734 | list AMONS = T(1), T(2), T(3), T(4), T(5), T(3)*T(5), T(4)*T(5), T(5)^2, T(5)^3; |
---|
1735 | int n0 = 5; |
---|
1736 | applyA(A, n0, AMONS); |
---|
1737 | } |
---|
1738 | |
---|
1739 | ////////////////////////////////////// |
---|
1740 | |
---|
1741 | // helper |
---|
1742 | // returns A * (generator corresponding to u) |
---|
1743 | // The action of A is stored in H. |
---|
1744 | static proc applyAtoGenerator(matrix u, matrix H, list KTwS, int n0){ |
---|
1745 | poly g0 = 0; |
---|
1746 | |
---|
1747 | for(int m = 1; m <= nrows(u); m++){ |
---|
1748 | if(u[m,1] != 0){ |
---|
1749 | intvec v = leadexp(KTwS[m]); // only the first n0 entries are interesting |
---|
1750 | poly g = u[m,1]; |
---|
1751 | |
---|
1752 | for(int x = 1; x <= n0; x++){ |
---|
1753 | g = g * H[1,x]^v[x]; |
---|
1754 | } |
---|
1755 | |
---|
1756 | g0 = g0 + g; |
---|
1757 | |
---|
1758 | kill v, g, x; |
---|
1759 | } |
---|
1760 | } |
---|
1761 | |
---|
1762 | return(g0); |
---|
1763 | } |
---|
1764 | |
---|
1765 | |
---|
1766 | ///////////////////////////////////// |
---|
1767 | |
---|
1768 | // helper |
---|
1769 | static proc linEqs(list IBw, list MON){ |
---|
1770 | // compute linear equations for the vector space IBw: |
---|
1771 | // the elements of the kernel MM correspond to linear forms: |
---|
1772 | // e.g., [1,0,-1,0] means S[1] - S[3]: |
---|
1773 | matrix B[size(IBw)][size(MON)]; |
---|
1774 | |
---|
1775 | for(int j = 1; j <= size(IBw); j++){ |
---|
1776 | B[j, 1..size(MON)] = IBw[j]; |
---|
1777 | } |
---|
1778 | |
---|
1779 | dbprint(printlevel, "Graded component as vectors:"); |
---|
1780 | dbprint(printlevel, B); |
---|
1781 | |
---|
1782 | // in this case, this computes the kernel |
---|
1783 | // (since there are only constant entries in B): |
---|
1784 | matrix M = syz(B); |
---|
1785 | |
---|
1786 | dbprint(printlevel, "Kernel:"); |
---|
1787 | dbprint(printlevel, M); |
---|
1788 | |
---|
1789 | return(M); |
---|
1790 | } |
---|
1791 | |
---|
1792 | ////////////////////////////////////// |
---|
1793 | |
---|
1794 | // helper |
---|
1795 | // creates the ring 0,(T(1..n0), Y(1..n1),Z),(dp(n0+n1),dp(1)); |
---|
1796 | static proc getFullRingS(intmat Q, int n0, int n1){ |
---|
1797 | intmat QS[nrows(Q)][ncols(Q)+1]; |
---|
1798 | |
---|
1799 | for(int l = 1; l <= ncols(Q); l++){ |
---|
1800 | QS[1..nrows(Q),l] = Q[1..nrows(Q),l]; |
---|
1801 | } |
---|
1802 | |
---|
1803 | int ii; |
---|
1804 | list l4; |
---|
1805 | for (ii = 1; ii <= n0; ii++) |
---|
1806 | { |
---|
1807 | l4[ii] = "T("+string(ii)+")"; |
---|
1808 | } |
---|
1809 | for (ii = 1; ii <= n1; ii++) |
---|
1810 | { |
---|
1811 | l4[n0+ii] = "Y("+string(ii)+")"; |
---|
1812 | } |
---|
1813 | l4[size(l4)+1] = "Z"; |
---|
1814 | ring S = create_ring(ring_list(basering)[1], l4, "(dp("+string(n0+n1)+"),dp(1))", "no_minpoly"); |
---|
1815 | setring S; |
---|
1816 | setBaseMultigrading(QS); |
---|
1817 | |
---|
1818 | return(S); |
---|
1819 | } |
---|
1820 | |
---|
1821 | ///////////////////////////////////// |
---|
1822 | |
---|
1823 | // helper |
---|
1824 | static proc monomials2basis(poly g0, list BMONS, int n0){ |
---|
1825 | // replace in g0 e.g. T[1]^2 by e_1: |
---|
1826 | // store coeffs for e_1 in Res[1] etc: |
---|
1827 | list RES; |
---|
1828 | |
---|
1829 | for(int m = 1; m <= size(BMONS); m++){ |
---|
1830 | RES[m] = 0; |
---|
1831 | } |
---|
1832 | |
---|
1833 | for(m = 1; m <= size(BMONS); m++){ |
---|
1834 | poly u0 = moncoef2(g0, BMONS[m], n0); |
---|
1835 | |
---|
1836 | if(u0 != 0){ |
---|
1837 | RES[m] = RES[m] + u0; |
---|
1838 | } |
---|
1839 | |
---|
1840 | kill u0; |
---|
1841 | } |
---|
1842 | |
---|
1843 | return(RES); |
---|
1844 | } |
---|
1845 | |
---|
1846 | ///////////////////////////////////// |
---|
1847 | |
---|
1848 | // helper |
---|
1849 | static proc getEqs(matrix MS, list RES){ |
---|
1850 | ideal GL; |
---|
1851 | |
---|
1852 | for(int m = 1; m <= ncols(MS); m++){ |
---|
1853 | poly f0 = 0; |
---|
1854 | |
---|
1855 | for(int l = 1; l <= nrows(MS); l++){ |
---|
1856 | f0 = f0 + RES[l] * MS[l,m]; |
---|
1857 | } |
---|
1858 | |
---|
1859 | GL[size(GL)+1] = f0; |
---|
1860 | kill f0, l; |
---|
1861 | } |
---|
1862 | |
---|
1863 | return(GL); |
---|
1864 | } |
---|
1865 | |
---|
1866 | ////////////////////////////////////// |
---|
1867 | |
---|
1868 | proc stabilizer(ideal RL, list #) |
---|
1869 | "USAGE: stabilizer(RL, A, BB, AMON, n0): RL is an ideal, A a matrix (standing for a subgroup of GL(n)), BB is an intmat (standing for an automorphism of the grading group), AMON a list of monomials corresponding to the rows/columns of A, n0 an integer such that the first n0 variables of the basering are the T(i), optional: a list of elementary divisors if there is torsion. |
---|
1870 | ASSUME: the basering must be graded (see setBaseMultigrading()) and the cone over |
---|
1871 | the degrees of the variables must be pointed; there mustn't be 0-degrees. The vector |
---|
1872 | w must be an element of the cone over the degrees of the variables. Moreover, B must be such that it permutes the degrees of the variables and the degrees of the generators of RL. |
---|
1873 | PURPOSE: returns relations such that A*I = I. |
---|
1874 | RETURN: a ring. Also exports an ideal Jexported and a list stabExported. |
---|
1875 | EXAMPLE: shows an example." |
---|
1876 | { |
---|
1877 | def R = basering; |
---|
1878 | intmat Q = getVariableWeights(); |
---|
1879 | |
---|
1880 | // torsion: |
---|
1881 | list TOR; |
---|
1882 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
1883 | TOR = #[1]; |
---|
1884 | } |
---|
1885 | |
---|
1886 | //list Components = stabPreprocessor(RL); |
---|
1887 | //list OmegaI = getOmegaI(RL); |
---|
1888 | int nTvars = nvars(basering); |
---|
1889 | list DegL = getGenDegsDims(RL, TOR); |
---|
1890 | |
---|
1891 | // output of aut_K(S): |
---|
1892 | def RR = autKS(TOR); // = \KK[Y_{ij},Z] |
---|
1893 | setring RR; |
---|
1894 | string AMONSstring = MONexported; // needed for later, dependent on T_i |
---|
1895 | |
---|
1896 | list listAutKS = autKSexported; |
---|
1897 | ideal J; // will be the result |
---|
1898 | list stabExported; // also export this |
---|
1899 | |
---|
1900 | // ring with both T_i and the Y_{ij} |
---|
1901 | int nYvars = nvars(basering) - 1; // no of Y variables |
---|
1902 | def S = getFullRingS(Q, nTvars, nYvars); // = \KK[T_k,Y_{ij},Z] |
---|
1903 | setring S; |
---|
1904 | |
---|
1905 | // prepare for later: |
---|
1906 | // matrix A = imap(RR, Agen); |
---|
1907 | execute("list AMONS = " + AMONSstring); // define AMONS |
---|
1908 | |
---|
1909 | int i; |
---|
1910 | setring RR; // KK[Y_i] |
---|
1911 | for(int p = 1; p <= size(listAutKS); p++){ |
---|
1912 | |
---|
1913 | matrix A = listAutKS[p][1]; |
---|
1914 | intmat BB = listAutKS[p][2]; |
---|
1915 | ideal JBA = listAutKS[p][3]; |
---|
1916 | |
---|
1917 | if(!isDimInvar(BB, DegL, TOR)){ |
---|
1918 | dbprint(printlevel, "B was not diminvar."); |
---|
1919 | |
---|
1920 | kill A, BB, JBA; |
---|
1921 | } else { |
---|
1922 | |
---|
1923 | setring S; |
---|
1924 | ideal EQs; |
---|
1925 | |
---|
1926 | setring R; // = \KK[T_i] |
---|
1927 | |
---|
1928 | // go through each generator g of the ideal RL: |
---|
1929 | // compute a basis for RL_deg(g), linear equations |
---|
1930 | for(i = 1; i <= size(RL); i++){ |
---|
1931 | dbprint(printlevel, "computing stabilizer: starting loop with i = " + string(i)); |
---|
1932 | |
---|
1933 | // w = multidegree and reduce later entries |
---|
1934 | // if there is torsion: |
---|
1935 | intvec w = multiDegRedTor(RL[i], TOR); |
---|
1936 | |
---|
1937 | // compute a basis KTw of \KK[T_1,...,T_r]_w: |
---|
1938 | list KTw = wmonomials(w, 0, TOR); |
---|
1939 | |
---|
1940 | // compute a basis Iw of I_w: |
---|
1941 | // I_w will be a subspace of KTw: |
---|
1942 | list L = gradedCompIdeal(RL, w, 0, TOR); |
---|
1943 | list Iw = L[2]; |
---|
1944 | list MON = L[1]; |
---|
1945 | |
---|
1946 | // compute a basis IBw of I_(B*w) if B*w != w: |
---|
1947 | // I_(B*w) will be a subspace of KTBw: |
---|
1948 | intvec Bw = BB * w; |
---|
1949 | |
---|
1950 | if(Bw != w){ |
---|
1951 | list LB = gradedCompIdeal(RL, Bw, 0, TOR); |
---|
1952 | list IBw = LB[2]; |
---|
1953 | list BMON = LB[1]; |
---|
1954 | } else { |
---|
1955 | list LB = L; |
---|
1956 | list IBw = L[2]; |
---|
1957 | list BMON = L[1]; |
---|
1958 | } |
---|
1959 | |
---|
1960 | matrix M = linEqs(IBw, MON); // kernel computation |
---|
1961 | |
---|
1962 | // Define the map |
---|
1963 | // T_i --> phi(T_i) = A * T_i: |
---|
1964 | // store A*T[i] in H[1,i]: |
---|
1965 | setring S; |
---|
1966 | matrix A = imap(RR, A); |
---|
1967 | matrix H = applyA(A, nTvars, AMONS); |
---|
1968 | |
---|
1969 | dbprint(printlevel, "the map phi (stored in H[1,j] = A*T_j):"); |
---|
1970 | dbprint(printlevel, H); |
---|
1971 | |
---|
1972 | // for each basis vector u of Iw: |
---|
1973 | // form g0 := A * (current generator) = Sum A*u_i where u_i <> 0: |
---|
1974 | list IwS = imap(R, Iw); |
---|
1975 | list KTwS = imap(R, KTw); |
---|
1976 | |
---|
1977 | for(int k = 1; k <= size(IwS); k++){ |
---|
1978 | matrix u = IwS[k]; // ... x 1 matrix |
---|
1979 | |
---|
1980 | poly g0 = applyAtoGenerator(u, H, KTwS, nTvars); |
---|
1981 | |
---|
1982 | dbprint(printlevel, "Iw="); |
---|
1983 | dbprint(printlevel, IwS); |
---|
1984 | dbprint(printlevel, "current basis vector u of Iw:"); |
---|
1985 | dbprint(printlevel, u); |
---|
1986 | dbprint(printlevel, "g0 = A * (current generator) = "); |
---|
1987 | dbprint(printlevel, g0); |
---|
1988 | |
---|
1989 | // replace in g0 e.g. T[1]^2 by e_1: |
---|
1990 | // store coeffs for e_1 in Res[1] etc: |
---|
1991 | |
---|
1992 | // NOTE: I_w gets mapped to I_(B*w): |
---|
1993 | list BMON = imap(R, BMON); |
---|
1994 | list RES = monomials2basis(g0, BMON, nTvars); // AMONS = BMONS? |
---|
1995 | |
---|
1996 | dbprint(printlevel, "after replacing in g0 e.g. T[1]^2 by e_1: (where e_i <--> KTwS[i]), RES = "); |
---|
1997 | dbprint(printlevel, RES); |
---|
1998 | matrix MS = imap(R, M); |
---|
1999 | ideal GL = getEqs(MS, RES); |
---|
2000 | |
---|
2001 | if(size(EQs) == 0){ |
---|
2002 | EQs = GL[1..size(GL)]; |
---|
2003 | } else { |
---|
2004 | EQs = EQs, GL[1..size(GL)]; |
---|
2005 | } |
---|
2006 | |
---|
2007 | dbprint(printlevel, "applying the linear forms of the kernel MS = "); |
---|
2008 | dbprint(printlevel, MS); |
---|
2009 | dbprint(printlevel, "yields equations GL = "); |
---|
2010 | dbprint(printlevel, GL); |
---|
2011 | dbprint(printlevel, "Total equations so far = EQs = "); |
---|
2012 | dbprint(printlevel, EQs); |
---|
2013 | |
---|
2014 | kill MS, GL, u, BMON, g0, RES; |
---|
2015 | setring S; |
---|
2016 | } |
---|
2017 | |
---|
2018 | setring S; |
---|
2019 | kill IwS, KTwS, A, H; |
---|
2020 | |
---|
2021 | setring R; |
---|
2022 | kill M, w, KTw, L, Iw, MON, Bw, LB, IBw, BMON; |
---|
2023 | } |
---|
2024 | |
---|
2025 | //kill k; |
---|
2026 | |
---|
2027 | setring RR; |
---|
2028 | ideal JBAnew = JBA + imap(S, EQs); |
---|
2029 | stabExported[size(stabExported) + 1] = list(A, BB, JBAnew); |
---|
2030 | |
---|
2031 | if(size(J) == 0){ |
---|
2032 | J = JBAnew; |
---|
2033 | } else { |
---|
2034 | //J = J * JBAnew; |
---|
2035 | J = intersect(J, JBAnew); // in Singular, intersection is usually quicker |
---|
2036 | } |
---|
2037 | |
---|
2038 | |
---|
2039 | setring S; |
---|
2040 | kill EQs; |
---|
2041 | setring RR; |
---|
2042 | kill JBAnew, BB, A, JBA; |
---|
2043 | } |
---|
2044 | } |
---|
2045 | |
---|
2046 | // return the ring RR, export the ideal J: |
---|
2047 | setring RR; |
---|
2048 | ideal Jexported = J; |
---|
2049 | export(Jexported); |
---|
2050 | export(stabExported); |
---|
2051 | |
---|
2052 | return(RR); |
---|
2053 | } |
---|
2054 | example |
---|
2055 | { |
---|
2056 | echo=2; |
---|
2057 | |
---|
2058 | ////////////////// |
---|
2059 | // example: fano 15: |
---|
2060 | intmat Q[1][5] = 3,3,2,2,1; |
---|
2061 | ring R = 0,T(1..5),dp; |
---|
2062 | |
---|
2063 | // attach degree matrix Q to R: |
---|
2064 | setBaseMultigrading(Q); |
---|
2065 | |
---|
2066 | ideal I = T(1)*T(2) + T(3)^2*T(4) + T(5)^6; |
---|
2067 | def RR = stabilizer(I); |
---|
2068 | setring RR; |
---|
2069 | RR; |
---|
2070 | Jexported; |
---|
2071 | dim(std(Jexported)); |
---|
2072 | getVariableWeights(); |
---|
2073 | |
---|
2074 | kill RR, Q, R; |
---|
2075 | |
---|
2076 | ///////////// |
---|
2077 | // example 3.14 from the paper |
---|
2078 | intmat Q[3][5] = |
---|
2079 | 1,1,1,1,1, |
---|
2080 | 1,-1,0,0,1, |
---|
2081 | 1,1,1,0,0; |
---|
2082 | |
---|
2083 | list TOR = 2; |
---|
2084 | ring R = 0,T(1..5),dp; |
---|
2085 | |
---|
2086 | // attach degree matrix Q to R: |
---|
2087 | setBaseMultigrading(Q); |
---|
2088 | |
---|
2089 | ideal I = T(1)*T(2) + T(3)^2 + T(4)^2; |
---|
2090 | list TOR = 2; |
---|
2091 | def RR = stabilizer(I, TOR); |
---|
2092 | setring RR; |
---|
2093 | RR; |
---|
2094 | Jexported; |
---|
2095 | dim(std(Jexported)); |
---|
2096 | |
---|
2097 | stabExported; |
---|
2098 | kill RR, Q, R; |
---|
2099 | } |
---|
2100 | |
---|
2101 | |
---|
2102 | ////////////////////////////////// |
---|
2103 | |
---|
2104 | // helper: |
---|
2105 | // computes the multidegree of f and reduces its later entries |
---|
2106 | // mod TOR[i] if # is defined |
---|
2107 | static proc multiDegRedTor(poly f, list #){ |
---|
2108 | list TOR; |
---|
2109 | |
---|
2110 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
2111 | TOR = #[1]; |
---|
2112 | } |
---|
2113 | |
---|
2114 | intmat Q = getVariableWeights(); |
---|
2115 | intvec w = multiDeg(f); |
---|
2116 | |
---|
2117 | // reduce w if there is torsion: |
---|
2118 | return(reduceIntvec(w, TOR)); |
---|
2119 | } |
---|
2120 | |
---|
2121 | /////////////////////////////////// |
---|
2122 | |
---|
2123 | // helper: |
---|
2124 | // reduce later entries of w0 |
---|
2125 | // mod TOR[i] if TOR is non-empty: |
---|
2126 | static proc reduceIntvec(intvec w0, list TOR){ |
---|
2127 | intvec w = w0; |
---|
2128 | int pos; |
---|
2129 | |
---|
2130 | for(int l = 1; l <= size(TOR); l++){ |
---|
2131 | pos = size(w0) - size(TOR) + l; |
---|
2132 | w[pos] = w[pos] mod TOR[l]; |
---|
2133 | |
---|
2134 | // make positive: |
---|
2135 | if(w[pos] < 0){ |
---|
2136 | w[pos] = w[pos] + TOR[l]; |
---|
2137 | } |
---|
2138 | } |
---|
2139 | |
---|
2140 | return(w); |
---|
2141 | } |
---|
2142 | |
---|
2143 | ///////////////////////////////// |
---|
2144 | |
---|
2145 | static proc moncoef2(poly f, poly mon, int n0) |
---|
2146 | "USAGE: moncoef(f, mon): f is a polynomial, mon a monomial, n0 an integer. |
---|
2147 | PURPOSE: returns the (not necessarily constant) coefficient of mon in f. |
---|
2148 | For example the coefficient of T(1)^2 in T(1)^2*T(2) will be T(2). |
---|
2149 | The integer n0 means that the variables var(n0+1),... will be considered to be coefficients. |
---|
2150 | RETURN: a ring element. |
---|
2151 | EXAMPLE: shows an example." |
---|
2152 | { |
---|
2153 | poly p = f; |
---|
2154 | poly res = 0; |
---|
2155 | |
---|
2156 | while (p != 0) { |
---|
2157 | poly lp = lead(p) / mon; |
---|
2158 | |
---|
2159 | if(lp != 0){ |
---|
2160 | // test whether there is still some T(i)-factor in lp: |
---|
2161 | for(int i = 1; i <= n0; i++){ |
---|
2162 | lp = subst(lp, var(i), 0); |
---|
2163 | } |
---|
2164 | |
---|
2165 | if(lp != 0){ |
---|
2166 | //return(lp); |
---|
2167 | res = res + lp; |
---|
2168 | } |
---|
2169 | |
---|
2170 | kill i; |
---|
2171 | } |
---|
2172 | |
---|
2173 | p = p - lead(p); |
---|
2174 | |
---|
2175 | kill lp; |
---|
2176 | } |
---|
2177 | |
---|
2178 | return(res); |
---|
2179 | } |
---|
2180 | example |
---|
2181 | { |
---|
2182 | echo=2; |
---|
2183 | |
---|
2184 | ring R = (0,a),T(1..4),dp; |
---|
2185 | |
---|
2186 | poly f = T(1)^2*T(2)^3 + 7*a*T(3)^3 -8*T(3)^2 +7; |
---|
2187 | poly mon = T(3)^3; |
---|
2188 | poly mon = 1; |
---|
2189 | |
---|
2190 | moncoef2(f, mon,4); |
---|
2191 | |
---|
2192 | // example 2: |
---|
2193 | ring R = 0,(T(1..3), Y(1..12)),dp; |
---|
2194 | poly f = T(1)^2*Y(1)^2 + Y(9)*T(1)*Y(10)^2*T(1); |
---|
2195 | poly m = T(1)^2; |
---|
2196 | |
---|
2197 | // only the first three variables |
---|
2198 | // are considered for taking the coefficient: |
---|
2199 | moncoef2(f, m, 3); |
---|
2200 | } |
---|
2201 | |
---|
2202 | ////////////////////////////// |
---|
2203 | |
---|
2204 | // helper |
---|
2205 | static proc fixesDim(intmat B, list ABL, list TOR){ |
---|
2206 | // ABL is of shape |
---|
2207 | // [1]: |
---|
2208 | // [1]: |
---|
2209 | // 2 = dim of basering_w |
---|
2210 | // [2]: |
---|
2211 | // 1,0 = w |
---|
2212 | // [3]: = generators of basering_w |
---|
2213 | // [1]: |
---|
2214 | // T(2) |
---|
2215 | // [2]: |
---|
2216 | // T(1) |
---|
2217 | int j; |
---|
2218 | for(int i = 1; i <= size(ABL); i++){ |
---|
2219 | int dimw = ABL[i][1]; // the dimension |
---|
2220 | intvec w = ABL[i][2]; // the weight vector |
---|
2221 | intvec Bw = reduceIntvec(B * w, TOR); // Bw must appear among some ABL[j] |
---|
2222 | |
---|
2223 | for(j = 1; j <= size(ABL); j++){ |
---|
2224 | if(ABL[j][2] == Bw){ |
---|
2225 | if(ABL[j][1] != dimw){ |
---|
2226 | return(0); |
---|
2227 | } |
---|
2228 | |
---|
2229 | break; |
---|
2230 | } |
---|
2231 | } |
---|
2232 | |
---|
2233 | if(j > size(ABL)){ |
---|
2234 | ERROR("B * w not found among generator weights. This should not happen."); |
---|
2235 | } |
---|
2236 | |
---|
2237 | kill dimw, w, Bw; |
---|
2238 | } |
---|
2239 | |
---|
2240 | return(1); |
---|
2241 | } |
---|
2242 | |
---|
2243 | /////////////////////////////// |
---|
2244 | |
---|
2245 | proc autKS(list #) |
---|
2246 | "USAGE: autKS(TOR); TOR: optional list of elementary divisors in case of torsion. |
---|
2247 | ASSUME: the basering is multigraded having used the command setBaseMultigrading(Q) from 'multigrading.lib'. |
---|
2248 | PURPOSE: Compute the subgroup Aut_K(S) of GL(n) of graded automorphisms of the polynomial ring S (the basering). |
---|
2249 | RETURN: returns a ring S and exports an ideal ideal Iexported in the coordinate ring S = K[Y_ij] of GL(n) such that Aut_K(S) = V(I). |
---|
2250 | EXAMPLE: shows an example." |
---|
2251 | { |
---|
2252 | list TOR; |
---|
2253 | |
---|
2254 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
2255 | TOR = #[1]; |
---|
2256 | } |
---|
2257 | |
---|
2258 | def R = basering; |
---|
2259 | int n0 = nvars(basering); |
---|
2260 | |
---|
2261 | list BL = allMonomials(ideal(0), TOR); //list BL = allMonomials(RL, TOR); |
---|
2262 | list MONS = sortMons(BL); |
---|
2263 | |
---|
2264 | // needed later for the S-admissible equations: |
---|
2265 | // Note: the 2nd components of the elements of BL |
---|
2266 | // are the origs so use BL: |
---|
2267 | // form first the list of all [w, [M, N]] |
---|
2268 | // where M and N contain the monomials T^a, T^b |
---|
2269 | // such that deg(T^(a+b)) = w. |
---|
2270 | list sums = origsSumUp(BL, TOR); |
---|
2271 | list adMons = getSadmissibleMonomials(sums); |
---|
2272 | int nadMons = size(adMons); |
---|
2273 | |
---|
2274 | // create formal matrix out of this: |
---|
2275 | // phi(R_w) = R_w for all w in Orig and all automorphisms phi. |
---|
2276 | def S = createAutMatrix(BL); |
---|
2277 | setring S; |
---|
2278 | basering; |
---|
2279 | |
---|
2280 | matrix A = Aexported; |
---|
2281 | list AMON = MONSexported; |
---|
2282 | string MONexported = string(AMON); |
---|
2283 | int n = size(MONSexported); |
---|
2284 | |
---|
2285 | dbprint(printlevel, ".... matrix A, the monomials corresponding to rows/cols:..."); |
---|
2286 | dbprint(printlevel, A); |
---|
2287 | dbprint(printlevel, AMON); |
---|
2288 | |
---|
2289 | // will be needed later for grading: |
---|
2290 | intmat Qnew = getCharacteristicGrading(AMON, TOR); |
---|
2291 | |
---|
2292 | // compute originary weights: |
---|
2293 | // return the elements as permutation matrices |
---|
2294 | setring R; |
---|
2295 | intmat Q = getVariableWeights(); |
---|
2296 | list Autor0 = autGenWeights(Q, TOR); |
---|
2297 | |
---|
2298 | // redefine omegaSdiagEqs as the product |
---|
2299 | // omegaSdiagEqs * (B^* omegaSdiagEqs) |
---|
2300 | // where B in Autor0: |
---|
2301 | setring S; |
---|
2302 | ideal J; |
---|
2303 | |
---|
2304 | // convert the elements of Autor0 to permutation matrices: |
---|
2305 | list Autor; |
---|
2306 | list ABL = imap(R, BL); |
---|
2307 | |
---|
2308 | for(int i = 1; i <= size(Autor0); i++){ |
---|
2309 | if(fixesDim(Autor0[i], ABL, TOR)){ |
---|
2310 | Autor[size(Autor)+1] = matrix2compper(Autor0[i], ABL, AMON, Q, TOR); // of type matrix |
---|
2311 | } else { |
---|
2312 | dbprint(printlevel, "Did not fix dimension."); |
---|
2313 | dbprint(printlevel, string(Autor0[i])); |
---|
2314 | } |
---|
2315 | } |
---|
2316 | |
---|
2317 | dbprint(printlevel, ".... automorphisms of the originary weight vectors as permutation matrices: ...."); |
---|
2318 | dbprint(printlevel, Autor); |
---|
2319 | |
---|
2320 | // equations cutting out the S-admissible matrices: |
---|
2321 | // A*(T^(a + b)) - (A*(T^a))(A*(T^b)) = 0 |
---|
2322 | // where A*f means to substitute T_i by sum(A_(i*)T_j): |
---|
2323 | list MONS = imap(R, MONS); |
---|
2324 | list adMons; |
---|
2325 | |
---|
2326 | if(nadMons > 0){ |
---|
2327 | adMons = imap(R, adMons); |
---|
2328 | } |
---|
2329 | |
---|
2330 | // add determinant condition: |
---|
2331 | // we need one more variable for this: |
---|
2332 | def S2 = adjoinVar(n0); |
---|
2333 | setring S2; |
---|
2334 | poly deter; // we need to store it for the grading |
---|
2335 | |
---|
2336 | ideal J = 1; // init with 1 --> used for product later |
---|
2337 | list listAutKS; // elements are of shape list(BA, B, equations for BA) |
---|
2338 | setring S; |
---|
2339 | |
---|
2340 | for(i = 1; i <= size(Autor); i++){ |
---|
2341 | intmat BB = Autor0[i]; |
---|
2342 | intmat B = Autor[i]; |
---|
2343 | |
---|
2344 | // compute B*A and then admissible equations |
---|
2345 | matrix BA = B*A; |
---|
2346 | |
---|
2347 | dbprint(printlevel, "Starting computation of admissible equations:"); |
---|
2348 | |
---|
2349 | ideal Iadmiss = getAdmissibleEquations(adMons, BA, n0, MONS); |
---|
2350 | dbprint(printlevel, "Admissible equations: done."); |
---|
2351 | |
---|
2352 | BA = relabelEntries(BA, n0); |
---|
2353 | |
---|
2354 | dbprint(printlevel, "permuted matrix A:"); |
---|
2355 | dbprint(printlevel, BA); |
---|
2356 | |
---|
2357 | // add variable Y((i-1)*n+j) to Jexported |
---|
2358 | // if BA[i,j] = 0 where 1 <= i <= n0: |
---|
2359 | ideal omegaSdiagEqs = getMonomialsFor0Entries(BA, n0); |
---|
2360 | |
---|
2361 | // also store in listAutKS |
---|
2362 | // will be useful for stabilizer |
---|
2363 | setring S2; |
---|
2364 | matrix BA2 = imap(S, BA); |
---|
2365 | |
---|
2366 | dbprint(printlevel, "computing determinant..."); |
---|
2367 | deter = det(BA2); |
---|
2368 | poly detEq = deter * Z - 1; |
---|
2369 | dbprint(printlevel, "Determinant: done."); |
---|
2370 | |
---|
2371 | ideal JBA = imap(S, Iadmiss) + imap(S, omegaSdiagEqs) + ideal(detEq); |
---|
2372 | listAutKS[size(listAutKS) + 1] = list(BA2, BB, JBA, MONexported); // MONexported is an explanatory string |
---|
2373 | |
---|
2374 | dbprint(printlevel, "Computing intersection of ideals..."); |
---|
2375 | //J = J * JBA; |
---|
2376 | J = intersect(J, JBA); // in Singular, this is usually quicker |
---|
2377 | dbprint(printlevel, "Intersection of ideal: done."); |
---|
2378 | |
---|
2379 | kill JBA, detEq, BA2; |
---|
2380 | setring S; |
---|
2381 | kill BA, B, BB, Iadmiss, omegaSdiagEqs; |
---|
2382 | } |
---|
2383 | |
---|
2384 | setring S2; |
---|
2385 | |
---|
2386 | setGLnGrading(Qnew, deter); |
---|
2387 | ideal Iexported = J; //imap(S, J) + ideal(detEq); |
---|
2388 | list autKSexported = listAutKS; |
---|
2389 | |
---|
2390 | export(Iexported); |
---|
2391 | export(autKSexported); // list of elements of shape list(A, B) |
---|
2392 | export(MONexported); |
---|
2393 | return(S2); |
---|
2394 | } |
---|
2395 | example |
---|
2396 | { |
---|
2397 | echo=2; |
---|
2398 | |
---|
2399 | ////////////////// |
---|
2400 | // example: fano 15: |
---|
2401 | intmat Q[1][5] = 3,3,2,2,1; |
---|
2402 | ring R = 0,T(1..5),dp; |
---|
2403 | |
---|
2404 | // attach degree matrix Q to R: |
---|
2405 | setBaseMultigrading(Q); |
---|
2406 | |
---|
2407 | //ideal I = T(1)*T(2) + T(3)^2*T(4) + T(5)^6; |
---|
2408 | def S = autKS(); |
---|
2409 | setring S; |
---|
2410 | |
---|
2411 | dim(std(Iexported)); |
---|
2412 | basering; |
---|
2413 | |
---|
2414 | autKSexported; |
---|
2415 | getVariableWeights(); |
---|
2416 | |
---|
2417 | kill S, Q, R; |
---|
2418 | |
---|
2419 | ///////////// |
---|
2420 | // example 3.14 from the paper |
---|
2421 | intmat Q[3][5] = |
---|
2422 | 1,1,1,1,1, |
---|
2423 | 1,-1,0,0,1, |
---|
2424 | 1,1,1,0,0; |
---|
2425 | |
---|
2426 | list TOR = 2; |
---|
2427 | ring R = 0,T(1..5),dp; |
---|
2428 | |
---|
2429 | // attach degree matrix Q to R: |
---|
2430 | setBaseMultigrading(Q); |
---|
2431 | |
---|
2432 | //ideal I = T(1)*T(2) + T(3)^2 + T(4)^2; |
---|
2433 | def S = autKS(); |
---|
2434 | setring S; |
---|
2435 | |
---|
2436 | Iexported; |
---|
2437 | print(getVariableWeights()); |
---|
2438 | kill S, R, Q; |
---|
2439 | } |
---|
2440 | |
---|
2441 | //////////////////////////////////// |
---|
2442 | |
---|
2443 | // helper |
---|
2444 | static proc setGLnGrading(intmat Qnew, poly deter, list #){ |
---|
2445 | list TOR; |
---|
2446 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
2447 | TOR = #[1]; |
---|
2448 | } |
---|
2449 | |
---|
2450 | // set deg(Y_ij) = j-th col of Qnew |
---|
2451 | int k = nrows(Qnew); |
---|
2452 | int n = ncols(Qnew); |
---|
2453 | |
---|
2454 | // the last variable is the rabinovic trick variable Z |
---|
2455 | int nv = nvars(basering) - 1; |
---|
2456 | intmat Qlarge[nrows(Qnew)][nv + 1]; |
---|
2457 | |
---|
2458 | for(int m = 1; m <= nv; m++){ |
---|
2459 | int i = m div n; // 9 div 9 = 1 |
---|
2460 | int j = m mod n; |
---|
2461 | |
---|
2462 | if(j == 0){ // 9 --> i = 0, j = 9 |
---|
2463 | i = i - 1; |
---|
2464 | j = n; |
---|
2465 | } |
---|
2466 | |
---|
2467 | Qlarge[1..k, m] = Qnew[1..k, j]; |
---|
2468 | kill i, j; |
---|
2469 | } |
---|
2470 | |
---|
2471 | // grading is valid for all variables |
---|
2472 | // except Z: use it for deter: |
---|
2473 | setBaseMultigrading(Qlarge); |
---|
2474 | |
---|
2475 | // treat the last variable Z: |
---|
2476 | // deg(deter) = - deg(Z): |
---|
2477 | intvec degDet = reduceIntvec(- multiDeg(deter), TOR); |
---|
2478 | Qlarge[1..k, nvars(basering)] = degDet[1..size(degDet)]; |
---|
2479 | |
---|
2480 | // now set the final grading |
---|
2481 | setBaseMultigrading(Qlarge); |
---|
2482 | } |
---|
2483 | |
---|
2484 | //////////////////////////////////// |
---|
2485 | |
---|
2486 | // helper |
---|
2487 | static proc adjoinVar(int n0){ |
---|
2488 | string S2str = "("; |
---|
2489 | |
---|
2490 | for(int i = n0 + 1; i <= nvars(basering); i++){ |
---|
2491 | S2str = S2str + string(var(i)); |
---|
2492 | |
---|
2493 | if(i < nvars(basering)){ |
---|
2494 | S2str = S2str + ", "; |
---|
2495 | } else { |
---|
2496 | S2str = S2str + ", Z)"; |
---|
2497 | } |
---|
2498 | } |
---|
2499 | ring S2 = create_ring(ring_list(basering)[1], S2str, "(dp("+string(nvars(basering) - n0)+"), dp(1))", "no_minpoly"); |
---|
2500 | return(S2); |
---|
2501 | } |
---|
2502 | |
---|
2503 | |
---|
2504 | ///////////////////////////////////// |
---|
2505 | |
---|
2506 | // helper |
---|
2507 | static proc getMonomialsFor0Entries(matrix BA, int n0) { |
---|
2508 | ideal relB; |
---|
2509 | int j; |
---|
2510 | |
---|
2511 | for(int i = 1; i <= n0; i++){ |
---|
2512 | for( j = 1; j <= ncols(BA); j++){ |
---|
2513 | if(BA[i,j] == 0){ |
---|
2514 | relB[size(relB)+1] = Y((i-1)*ncols(BA)+j); |
---|
2515 | } |
---|
2516 | } |
---|
2517 | } |
---|
2518 | |
---|
2519 | return(relB); |
---|
2520 | } |
---|
2521 | |
---|
2522 | ///////////////////////////////////// |
---|
2523 | |
---|
2524 | // helper |
---|
2525 | static proc getAdmissibleEquations(list adMons, matrix A, int n0, list MONS){ |
---|
2526 | // S is the name of the current basering when calling this function, S --> S |
---|
2527 | ideal imageVars; // images of the T_i |
---|
2528 | int j; |
---|
2529 | |
---|
2530 | for(int i = 1; i <= n0; i++){ |
---|
2531 | imageVars[i] = 0; |
---|
2532 | |
---|
2533 | for(j = 1; j <= size(MONS); j++){ |
---|
2534 | imageVars[i] = imageVars[i] + A[i, j] * MONS[j]; |
---|
2535 | } |
---|
2536 | } |
---|
2537 | def S = basering; |
---|
2538 | map applyA = S, imageVars; |
---|
2539 | |
---|
2540 | // go through all possibilities: |
---|
2541 | int k, l; |
---|
2542 | ideal IadmissibleExported; |
---|
2543 | |
---|
2544 | for(i = 1; i <= size(adMons); i++){ |
---|
2545 | list AL = adMons[i][2][1]; |
---|
2546 | list BL = adMons[i][2][2]; |
---|
2547 | |
---|
2548 | for(k = 1; k <= size(AL); k++){ |
---|
2549 | for(l = 1; l <= size(BL); l++){ |
---|
2550 | poly Ta = AL[k]; |
---|
2551 | poly Tb = BL[l]; |
---|
2552 | poly Tab = Ta * Tb; |
---|
2553 | |
---|
2554 | IadmissibleExported[size(IadmissibleExported) + 1] = applyA(Tab) - applyA(Ta) * applyA(Tb); |
---|
2555 | |
---|
2556 | kill Tab, Ta, Tb; |
---|
2557 | } |
---|
2558 | } |
---|
2559 | |
---|
2560 | kill AL, BL; |
---|
2561 | } |
---|
2562 | |
---|
2563 | return(IadmissibleExported); |
---|
2564 | } |
---|
2565 | |
---|
2566 | /////////////////////////////////////// |
---|
2567 | |
---|
2568 | // Get generator degrees w1,...,ws --> store in GenDegs. |
---|
2569 | // Also compute dimensions of I_wi where wi in GenDegs: |
---|
2570 | static proc getGenDegsDims(ideal I, list TOR){ |
---|
2571 | list GenDegs; |
---|
2572 | list GenDims; |
---|
2573 | int i = 1; |
---|
2574 | int c = 1; |
---|
2575 | |
---|
2576 | while(c <= size(I)){ // there may be 0-entries in I |
---|
2577 | if(I[i] != 0){ |
---|
2578 | c++; |
---|
2579 | |
---|
2580 | GenDegs[i] = reduceIntvec(multiDeg(I[i]), TOR); |
---|
2581 | list L = gradedCompIdeal(I, GenDegs[i], 0, TOR); //gradedcompbasis(I, GenDegs[i], TOR); |
---|
2582 | GenDims[i] = size(L[1]); // this is dim(R_(deg f_i)), not dim(I_(deg f_i)) |
---|
2583 | |
---|
2584 | kill L; |
---|
2585 | } |
---|
2586 | |
---|
2587 | i++; |
---|
2588 | } |
---|
2589 | |
---|
2590 | list RES; |
---|
2591 | RES[1] = GenDegs; |
---|
2592 | RES[2] = GenDims; |
---|
2593 | |
---|
2594 | return(RES); |
---|
2595 | } |
---|
2596 | |
---|
2597 | /////////////////////////////////// |
---|
2598 | |
---|
2599 | // helper: |
---|
2600 | // checks a given matrix B in Autor0 |
---|
2601 | // if it satisfies dim(I_wi) = dim(I_(B*wi)) |
---|
2602 | // for all i where w1,...,ws |
---|
2603 | // are the generator degrees of I. |
---|
2604 | static proc isDimInvar(intmat B, list GenDegsDims, list TOR){ |
---|
2605 | list GenDegs = GenDegsDims[1]; |
---|
2606 | list GenDims = GenDegsDims[2]; |
---|
2607 | |
---|
2608 | // return true iff B is |
---|
2609 | // such that dim(I_wi) = dim(I_(B*wi)) |
---|
2610 | // holds for all i, where wi in GenDegs: |
---|
2611 | int j; |
---|
2612 | |
---|
2613 | for(j = 1; j <= size(GenDegs); j++){ |
---|
2614 | intvec wj = GenDegs[j]; |
---|
2615 | |
---|
2616 | //intvec Bwj = B * wj; |
---|
2617 | intvec Bwj = reduceIntvec(intvec(B * wj), TOR); |
---|
2618 | |
---|
2619 | // find out dim(I_Bwj): |
---|
2620 | int dimBwj; |
---|
2621 | for(int bb = 1; bb <= size(GenDegs); bb++){ |
---|
2622 | if(GenDegs[bb] == Bwj){ |
---|
2623 | dimBwj = GenDims[bb]; |
---|
2624 | break; |
---|
2625 | } |
---|
2626 | } |
---|
2627 | |
---|
2628 | // if dim(I_wi) != dim(I_(B*wi)), then B is not valid. |
---|
2629 | if(dimBwj != GenDims[j]){ |
---|
2630 | return(0); |
---|
2631 | } |
---|
2632 | |
---|
2633 | kill bb, wj, Bwj, dimBwj; |
---|
2634 | } |
---|
2635 | |
---|
2636 | return(1); |
---|
2637 | } |
---|
2638 | |
---|
2639 | |
---|
2640 | //////////////////////////////////////// |
---|
2641 | |
---|
2642 | static proc getSadmissibleMonomials(list sums){ |
---|
2643 | list Combs; |
---|
2644 | |
---|
2645 | // the elements of sums are of shape |
---|
2646 | // Comb = [c, [a,b]] where c = a + b. |
---|
2647 | // Form first the list of all [K, [M, N]] |
---|
2648 | // where M and N contain the monomials T^a, T^b |
---|
2649 | // such that deg(T^(a+b)) = deg(T^c) for each T^c in K. |
---|
2650 | for(int i= 1; i <= size(sums); i++){ |
---|
2651 | list K = wmonomials(sums[i][1], 0); // as monomials |
---|
2652 | list M = wmonomials(sums[i][2][1], 0); // as monomials |
---|
2653 | list N = wmonomials(sums[i][2][2], 0); // as monomials |
---|
2654 | |
---|
2655 | Combs[size(Combs) + 1] = list(K, list(M, N)); |
---|
2656 | |
---|
2657 | kill M, N, K; |
---|
2658 | } |
---|
2659 | |
---|
2660 | return(Combs); |
---|
2661 | } |
---|
2662 | |
---|
2663 | /////////////////////////////////////// |
---|
2664 | |
---|
2665 | // helper |
---|
2666 | static proc origsSumUp(list BL, list #) { |
---|
2667 | list TOR; |
---|
2668 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
2669 | TOR = #[1]; |
---|
2670 | } |
---|
2671 | |
---|
2672 | list RES; // elements of the form [c, [a,b]] such that c=a+b |
---|
2673 | int j, k; |
---|
2674 | |
---|
2675 | // note: the 2nd components of the entries |
---|
2676 | // of BL are the origs |
---|
2677 | for(int i = 1; i <= size(BL); i++){ |
---|
2678 | for(j = 1; j <= size(BL); j++){ |
---|
2679 | intvec sum = reduceIntvec(BL[i][2] + BL[j][2], TOR); |
---|
2680 | |
---|
2681 | // check if sum is not already in RES (no duplicates): |
---|
2682 | for(k = 1; k <= size(RES); k++){ |
---|
2683 | if(RES[k][1] == sum){ |
---|
2684 | break; |
---|
2685 | } |
---|
2686 | } |
---|
2687 | |
---|
2688 | if(k > size(RES)){ // then is it new |
---|
2689 | // check if sum is also in Origs |
---|
2690 | for(k = 1; k <= size(BL); k++){ |
---|
2691 | if(BL[k][2] == sum){ |
---|
2692 | RES[size(RES) + 1] = list(sum, list(BL[i][2], BL[j][2])); |
---|
2693 | } |
---|
2694 | } |
---|
2695 | } |
---|
2696 | |
---|
2697 | kill sum; |
---|
2698 | } |
---|
2699 | } |
---|
2700 | |
---|
2701 | return(RES); |
---|
2702 | } |
---|
2703 | example |
---|
2704 | { |
---|
2705 | echo = 2; |
---|
2706 | |
---|
2707 | // 1st example |
---|
2708 | list BL = |
---|
2709 | list(0, intvec(2,1), 0), |
---|
2710 | list(0, intvec(1,1), 0), |
---|
2711 | list(0, intvec(1,0), 0), |
---|
2712 | list(0, intvec(2,1), 0); |
---|
2713 | |
---|
2714 | list TOR = 2; |
---|
2715 | |
---|
2716 | origsSumUp(BL); |
---|
2717 | kill BL; |
---|
2718 | |
---|
2719 | // 2nd example |
---|
2720 | list BL = |
---|
2721 | list(0, intvec(3), 0), |
---|
2722 | list(0, intvec(3), 0), |
---|
2723 | list(0, intvec(2), 0), |
---|
2724 | list(0, intvec(2), 0), |
---|
2725 | list(0, intvec(1), 0); |
---|
2726 | |
---|
2727 | origsSumUp(BL); |
---|
2728 | } |
---|
2729 | |
---|
2730 | /////////////////////////////////////// |
---|
2731 | |
---|
2732 | // helper |
---|
2733 | // Given an automorphism K --> K as an invertible matrix A |
---|
2734 | // this function computes a permutation matrix. |
---|
2735 | // Assume: BL is of the form [[dim, w, [T[1]*T[2]^2, T[3]]],...] |
---|
2736 | // where the w are the originary weights. |
---|
2737 | // we also assume that Origs are ordered as occurring in Q |
---|
2738 | // returns a list of matrices (not intmats). |
---|
2739 | static proc matrix2compper(intmat B, list BL, list MON, intmat Q, list TOR){ |
---|
2740 | list Origs; |
---|
2741 | list MONweight; |
---|
2742 | |
---|
2743 | for(int i = 1; i <= size(BL); i++){ |
---|
2744 | Origs[i] = BL[i][2]; |
---|
2745 | } |
---|
2746 | |
---|
2747 | for(int j = 1; j <= size(MON); j++){ |
---|
2748 | MONweight[j] = reduceIntvec(multiDeg(MON[j]), TOR); |
---|
2749 | } |
---|
2750 | |
---|
2751 | // initialize the permutation matrix (to be returned) |
---|
2752 | // as zero-rows: |
---|
2753 | int d = size(MON); |
---|
2754 | intmat PER[d][d]; |
---|
2755 | intvec done = 0:size(MONweight); |
---|
2756 | |
---|
2757 | for(i = 1; i <= size(MONweight); i++){ |
---|
2758 | intvec w = MONweight[i]; |
---|
2759 | intmat mW[size(w)][1]; |
---|
2760 | mW[1..size(w),1] = w; |
---|
2761 | |
---|
2762 | intmat BmW = B * mW; // ...x1 matrix |
---|
2763 | intvec v = BmW[1..nrows(BmW), 1]; |
---|
2764 | intvec vred = reduceIntvec(v, TOR); |
---|
2765 | |
---|
2766 | for(j = 1; j <= size(MONweight); j++){ |
---|
2767 | intvec zred = reduceIntvec(MONweight[j], TOR); |
---|
2768 | |
---|
2769 | // test if MONweight[j] == v (possibly with Torsion, i.e. after reduction) |
---|
2770 | // and done[j] != 1: |
---|
2771 | if(vred == zred and done[j] != 1){ |
---|
2772 | PER[i,j] = 1; |
---|
2773 | done[j] = 1; |
---|
2774 | |
---|
2775 | kill zred; |
---|
2776 | break; |
---|
2777 | } |
---|
2778 | |
---|
2779 | kill zred; |
---|
2780 | } |
---|
2781 | |
---|
2782 | kill vred, v, w, mW, BmW; |
---|
2783 | } |
---|
2784 | |
---|
2785 | return(PER); |
---|
2786 | } |
---|
2787 | |
---|
2788 | //////////////////////////////////////// |
---|
2789 | |
---|
2790 | proc autXhat(ideal RL, intvec w, list #) |
---|
2791 | " |
---|
2792 | USAGE: autXhat(RL, w0, TOR): RL is an ideal, w an intvec, TOR a list of integers |
---|
2793 | ASSUME: the basering is multigraded, the elements of TOR stand for the torsion rows of the matrix getVariableWeights(), w is an ample class or the free part of such a class. |
---|
2794 | PURPOSE: compute an ideal J such that V(J) in some GL(n) is isomorphic to the H-equivariant automorphisms \widehat X --> \widehat X. |
---|
2795 | EXAMPLE: example autXhat; shows an example |
---|
2796 | " |
---|
2797 | { |
---|
2798 | list TOR; |
---|
2799 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
2800 | TOR = #[1]; |
---|
2801 | } |
---|
2802 | |
---|
2803 | def R = basering; |
---|
2804 | intmat Q = getVariableWeights(); |
---|
2805 | intmat Q0 = gradingFreePart(Q, TOR); |
---|
2806 | |
---|
2807 | int k0 = nrows(Q0); |
---|
2808 | bigintmat w0[1][k0] = w[1..k0]; |
---|
2809 | |
---|
2810 | if(size(w) < nrows(Q)){ |
---|
2811 | ERROR("size of w must equal the number of rows of the degree matrix."); |
---|
2812 | } |
---|
2813 | |
---|
2814 | cone c = gitCone(RL, Q0, w0); |
---|
2815 | def RR = stabilizer(RL, TOR); |
---|
2816 | setring RR; |
---|
2817 | |
---|
2818 | list RES; |
---|
2819 | for(int i = 1; i <= size(stabExported); i++){ |
---|
2820 | intmat B = stabExported[i][2]; |
---|
2821 | intvec wB = B * w; |
---|
2822 | |
---|
2823 | setring R; |
---|
2824 | bigintmat wB0[1][k0] = wB[1..k0]; |
---|
2825 | cone cB = gitCone(RL, Q0, wB0); |
---|
2826 | kill wB0; |
---|
2827 | setring RR; |
---|
2828 | |
---|
2829 | if(cB == c){ |
---|
2830 | RES[size(RES) + 1] = stabExported[i]; |
---|
2831 | } |
---|
2832 | |
---|
2833 | kill B, wB, cB; |
---|
2834 | } |
---|
2835 | |
---|
2836 | export(RES); |
---|
2837 | return(RR); |
---|
2838 | } |
---|
2839 | example |
---|
2840 | { |
---|
2841 | echo=2; |
---|
2842 | |
---|
2843 | intmat Q[3][5] = |
---|
2844 | 1,1,1,1,1, |
---|
2845 | 1,-1,0,0,1, |
---|
2846 | 1,1,1,0,0; |
---|
2847 | |
---|
2848 | list TOR = 2; |
---|
2849 | ring R = 0,T(1..5),dp; |
---|
2850 | |
---|
2851 | setBaseMultigrading(Q); |
---|
2852 | |
---|
2853 | ideal I = T(1)*T(2) + T(3)^2 + T(4)^2; |
---|
2854 | |
---|
2855 | intvec w0 = 2,1,0; |
---|
2856 | def RR = autXhat(I, w0, TOR); |
---|
2857 | setring RR; |
---|
2858 | |
---|
2859 | RES; |
---|
2860 | kill RR, Q, R; |
---|
2861 | } |
---|
2862 | |
---|
2863 | ////////////////////////////////////////// |
---|
2864 | |
---|
2865 | static proc gradingFreePart(intmat Q, list TOR){ |
---|
2866 | intmat Q0[nrows(Q) - size(TOR)][ncols(Q)]; |
---|
2867 | |
---|
2868 | for(int i = 1; i <= nrows(Q0); i++){ |
---|
2869 | Q0[i, 1..ncols(Q0)] = Q[i, 1..ncols(Q0)]; |
---|
2870 | } |
---|
2871 | |
---|
2872 | return(Q0); |
---|
2873 | } |
---|
2874 | |
---|
2875 | //////////////////////////////////////// |
---|
2876 | |
---|
2877 | // helper |
---|
2878 | // compute the image cone by applying A to |
---|
2879 | // the rays of c. We assume c is pointed. |
---|
2880 | static proc imageCone(cone c, intmat A){ |
---|
2881 | bigintmat R = rays(c); |
---|
2882 | cone cc = coneViaPoints(R * transpose(A)); |
---|
2883 | |
---|
2884 | return(cc); |
---|
2885 | } |
---|
2886 | example |
---|
2887 | { |
---|
2888 | echo=2; |
---|
2889 | |
---|
2890 | intmat M[2][2] = |
---|
2891 | 1,0, |
---|
2892 | 1,2; |
---|
2893 | |
---|
2894 | cone c = coneViaPoints(M); |
---|
2895 | rays(c); |
---|
2896 | |
---|
2897 | intmat A[2][2] = |
---|
2898 | 1,0, |
---|
2899 | 0,2; |
---|
2900 | |
---|
2901 | cone cc = imageCone(c, A); |
---|
2902 | rays(cc); |
---|
2903 | } |
---|
2904 | |
---|
2905 | ///////////////////////////////////// |
---|
2906 | |
---|
2907 | // helper. |
---|
2908 | // Assume: basering variables are |
---|
2909 | // T(1..n0), Y(1..), Z: |
---|
2910 | // Replaces in M each non-zero entry M[i,j] by V[i,j]: |
---|
2911 | // however: replace the lower entries, e.g., Y(11)^2, |
---|
2912 | // by Y(2)^2 if Y(11) --> Y(2). |
---|
2913 | static proc relabelEntries(matrix A, int n0){ |
---|
2914 | def R = basering; |
---|
2915 | int n = ncols(A); |
---|
2916 | |
---|
2917 | // define a map: A[i,j] --> V[i,j] |
---|
2918 | // for the case of M[i,j] being a single variable; |
---|
2919 | // this is true for the first n0 rows of A: |
---|
2920 | // Store in Img[i] the image of Y(i): |
---|
2921 | intvec v = 0:(nvars(R)); |
---|
2922 | list Img = v[1..size(v)]; // initialize with 0. |
---|
2923 | |
---|
2924 | int i, j, pos, k, m; |
---|
2925 | |
---|
2926 | for(i = 1; i <= n0; i++){ |
---|
2927 | for(j = 1; j <= ncols(A); j++){ |
---|
2928 | if(A[i,j] != 0){ // then it is some Y(k) --> replace it by V[i,j]: |
---|
2929 | // ordering of variables in R is |
---|
2930 | // T(1..n0), Y(1..n^2), Z: |
---|
2931 | pos = (i-1)*n + j; |
---|
2932 | |
---|
2933 | // A[i,j] should be Y(pos): |
---|
2934 | // find out k with A[i,j] = Y(k): |
---|
2935 | for(m = 1; m <= n*n; m++){ |
---|
2936 | if(A[i,j] == Y(m)){ |
---|
2937 | Img[m + n0] = Y(pos); //+n0 since vars T(1..n0) are at front |
---|
2938 | break; |
---|
2939 | } |
---|
2940 | } |
---|
2941 | } |
---|
2942 | } |
---|
2943 | } |
---|
2944 | |
---|
2945 | // define the map: |
---|
2946 | map f = R, Img[1..size(Img)]; |
---|
2947 | |
---|
2948 | return(f(A)); |
---|
2949 | } |
---|
2950 | |
---|
2951 | /////////////////////// |
---|
2952 | |
---|
2953 | // helper |
---|
2954 | // returns the grading by the characteristic |
---|
2955 | // quasitorus. |
---|
2956 | static proc getCharacteristicGrading(list BL, list #){ |
---|
2957 | list TOR; |
---|
2958 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
2959 | TOR = #[1]; |
---|
2960 | } |
---|
2961 | |
---|
2962 | intmat Q = getVariableWeights(); |
---|
2963 | intmat Qnew[nrows(Q)][size(BL)]; |
---|
2964 | |
---|
2965 | for(int i = 1; i <= size(BL); i++){ |
---|
2966 | intvec w = reduceIntvec(multiDeg(BL[i]), TOR); |
---|
2967 | |
---|
2968 | Qnew[1..nrows(Q), i] = w[1..size(w)]; |
---|
2969 | kill w; |
---|
2970 | } |
---|
2971 | |
---|
2972 | return(Qnew); |
---|
2973 | } |
---|
2974 | |
---|
2975 | /////////////////////////// |
---|
2976 | |
---|
2977 | // helper |
---|
2978 | // compute the linear hull over the columns of P. |
---|
2979 | static proc linearHull(intmat P){ |
---|
2980 | // the linear hull over the columns of P: |
---|
2981 | // equivalent: cone(P, -P): |
---|
2982 | intmat M[ncols(P) * 2][nrows(P)]; // rows are rays |
---|
2983 | |
---|
2984 | for(int i = 1; i <= ncols(P); i++){ |
---|
2985 | M[i, 1..nrows(P)] = P[1..nrows(P), i]; |
---|
2986 | M[i + ncols(P), 1..nrows(P)] = -P[1..nrows(P), i]; |
---|
2987 | } |
---|
2988 | cone V = coneViaPoints(M); |
---|
2989 | |
---|
2990 | return(V); |
---|
2991 | } |
---|
2992 | |
---|
2993 | //////////////////////////// |
---|
2994 | |
---|
2995 | // helper |
---|
2996 | // compute generators for the veronese subalgebra |
---|
2997 | static proc veron(intmat P){ |
---|
2998 | // linear hull of P: |
---|
2999 | intmat Pplusminus[2*nrows(P)][ncols(P)] = P, -P; |
---|
3000 | |
---|
3001 | cone V = coneViaPoints(Pplusminus); |
---|
3002 | cone posorth = coneViaPoints(getIdMat(ambientDimension(V))); |
---|
3003 | cone c = convexIntersection(V, posorth); |
---|
3004 | |
---|
3005 | intmat H = hilbertBas(c); |
---|
3006 | return(H); |
---|
3007 | } |
---|
3008 | example |
---|
3009 | { |
---|
3010 | echo = 2; |
---|
3011 | |
---|
3012 | intmat P[2][3] = |
---|
3013 | 1,0,-1, |
---|
3014 | 0,1,-1; |
---|
3015 | |
---|
3016 | veron(transpose(P)); |
---|
3017 | } |
---|
3018 | |
---|
3019 | ////////////////////////// |
---|
3020 | |
---|
3021 | proc autX(ideal RL, intvec w, list #) |
---|
3022 | " |
---|
3023 | USAGE: autX(RL, w, TOR); RL: ideal, w: intvec, TOR: optional list of integers. |
---|
3024 | PURPOSE: compute generators for the hopf algebra O(Aut(X)) |
---|
3025 | of the Mori dream space X given by Cox(X) := basering/RL and |
---|
3026 | the ample class w. |
---|
3027 | ASSUME: there is no torsion. |
---|
3028 | RETURN: a ring. Also exports an ideal Iexported. |
---|
3029 | EXAMPLE: example autX; shows an example |
---|
3030 | " |
---|
3031 | { |
---|
3032 | list TOR; |
---|
3033 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
3034 | TOR = #[1]; |
---|
3035 | } |
---|
3036 | |
---|
3037 | def RR = autXhat(RL, w, TOR); |
---|
3038 | setring RR; |
---|
3039 | intmat Q = getVariableWeights(); |
---|
3040 | |
---|
3041 | // ideal: |
---|
3042 | ideal J = RES[1][3]; |
---|
3043 | for(int i = 2; i <= size(RES); i++){ |
---|
3044 | J = J * RES[i][3]; |
---|
3045 | } |
---|
3046 | |
---|
3047 | intmat P = degMat2P(Q, TOR); |
---|
3048 | intmat V = veron(transpose(P)); |
---|
3049 | |
---|
3050 | // creat map: |
---|
3051 | ideal imMap; |
---|
3052 | for(int i = 1; i <= nrows(V); i++){ |
---|
3053 | intvec v = V[i, 1..ncols(V)]; |
---|
3054 | |
---|
3055 | imMap[size(imMap) + 1] = vec2monomial(v); |
---|
3056 | kill v; |
---|
3057 | } |
---|
3058 | |
---|
3059 | list l5; |
---|
3060 | for (int ii = 1; ii <= size(V); ii++) |
---|
3061 | { |
---|
3062 | l5[ii] = "T("+string(ii)+")"; |
---|
3063 | } |
---|
3064 | ring S = create_ring(ring_list(RR)[1], l5, "dp", "no_minpoly"); |
---|
3065 | setring S; |
---|
3066 | setring RR; |
---|
3067 | |
---|
3068 | // f: S --> RR, T_i --> T^mu |
---|
3069 | map f = S, imMap; |
---|
3070 | |
---|
3071 | setring S; |
---|
3072 | ideal Iexported = preimage(RR, f, J); |
---|
3073 | |
---|
3074 | export(Iexported); |
---|
3075 | return(S); |
---|
3076 | } |
---|
3077 | example |
---|
3078 | { |
---|
3079 | /////////////// |
---|
3080 | //// CAREFUL: the following examples seems to be unfeasible at the moment, see remark in the paper |
---|
3081 | |
---|
3082 | //echo=2; |
---|
3083 | /////////////// |
---|
3084 | //// PP2 |
---|
3085 | //intmat Q[1][4] = |
---|
3086 | // 1,1,1,1; |
---|
3087 | |
---|
3088 | //ring R = 0,T(1..ncols(Q)),dp; |
---|
3089 | |
---|
3090 | //// attach degree matrix Q to R: |
---|
3091 | //setBaseMultigrading(Q); |
---|
3092 | //ideal I = 0; |
---|
3093 | //intvec w0 = 1; |
---|
3094 | |
---|
3095 | //def RR = autX(I, w0); |
---|
3096 | //setring RR; |
---|
3097 | //Iexported; |
---|
3098 | |
---|
3099 | //basering; |
---|
3100 | //dim(std(Iexported)); |
---|
3101 | |
---|
3102 | //kill RR, Q, R; |
---|
3103 | |
---|
3104 | /////////////// |
---|
3105 | //// example 3.14 from the paper |
---|
3106 | //intmat Q[3][5] = |
---|
3107 | // 1,1,1,1,1, |
---|
3108 | // 1,-1,0,0,1, |
---|
3109 | // 1,1,1,0,0; |
---|
3110 | |
---|
3111 | //list TOR = 2; |
---|
3112 | //ring R = 0,T(1..5),dp; |
---|
3113 | |
---|
3114 | //// attach degree matrix Q to R: |
---|
3115 | //setBaseMultigrading(Q); |
---|
3116 | |
---|
3117 | //ideal I = T(1)*T(2) + T(3)^2 + T(4)^2; |
---|
3118 | //list TOR = 2; |
---|
3119 | |
---|
3120 | //intvec w0 = 2,1,0; |
---|
3121 | //def RR = autX(I, w0, TOR); |
---|
3122 | //setring RR; |
---|
3123 | |
---|
3124 | //kill RR, Q, R; |
---|
3125 | } |
---|
3126 | |
---|
3127 | //////////////////// |
---|
3128 | |
---|
3129 | // helper |
---|
3130 | // compute the hilbert basis of c. |
---|
3131 | static proc hilbertBas(cone c){ |
---|
3132 | intmat A = intmat(rays(c)); |
---|
3133 | intmat B = normaliz(A, 0); |
---|
3134 | |
---|
3135 | return(B); |
---|
3136 | } |
---|
3137 | example |
---|
3138 | { |
---|
3139 | echo = 2; |
---|
3140 | |
---|
3141 | intmat A[2][2] = |
---|
3142 | 1,0, |
---|
3143 | 1,3; |
---|
3144 | |
---|
3145 | intmat B = hilbertBas(coneViaPoints(A)); |
---|
3146 | print(B); |
---|
3147 | |
---|
3148 | // 2nd ex: |
---|
3149 | intmat C[3][4] = |
---|
3150 | 1, 0, 0, 0, |
---|
3151 | 1, 1, 0, 0, |
---|
3152 | 6, 3, 4, 2; |
---|
3153 | |
---|
3154 | intmat D = hilbertBas(coneViaPoints(C)); |
---|
3155 | print(D); |
---|
3156 | |
---|
3157 | } |
---|
3158 | |
---|
3159 | ///////////////////////// |
---|
3160 | |
---|
3161 | // helper |
---|
3162 | // compute the gale dual P of the degree matrix Q |
---|
3163 | static proc degMat2P(intmat Q, list TOR){ |
---|
3164 | int k = nrows(Q); |
---|
3165 | int nfree = nrows(Q) - size(TOR); |
---|
3166 | |
---|
3167 | int n = size(TOR); |
---|
3168 | if(n == 0){ |
---|
3169 | n = 1; // then L = (0,...0)^t |
---|
3170 | } |
---|
3171 | |
---|
3172 | intmat L[k][n]; |
---|
3173 | |
---|
3174 | for(int i = 1; i <= size(TOR); i++){ |
---|
3175 | L[i + nfree, i] = TOR[i]; |
---|
3176 | } |
---|
3177 | intmat U0 = preimageLattice(Q, L); |
---|
3178 | |
---|
3179 | return(transpose(U0)); |
---|
3180 | } |
---|
3181 | example |
---|
3182 | { |
---|
3183 | echo = 2; |
---|
3184 | |
---|
3185 | intmat Q[2][4] = |
---|
3186 | 1,1,1,1, |
---|
3187 | 1,0,1,0; |
---|
3188 | |
---|
3189 | list TOR = 2; |
---|
3190 | intmat P = degMat2P(Q, TOR); |
---|
3191 | print(P); |
---|
3192 | |
---|
3193 | } |
---|
3194 | |
---|
3195 | /////////////////////////// |
---|
3196 | |
---|
3197 | // helper |
---|
3198 | static proc compsAreTrivial(ideal RL, list TOR){ |
---|
3199 | // for each col w of Q |
---|
3200 | // we require I_w = 0: |
---|
3201 | intmat Q = getVariableWeights(); |
---|
3202 | |
---|
3203 | for(int i = 1; i <= ncols(Q); i++){ |
---|
3204 | intvec w = Q[1..nrows(Q), i]; |
---|
3205 | |
---|
3206 | // compute a basis Iw of I_w: |
---|
3207 | // I_w will be a subspace of KTw: |
---|
3208 | list Iw = gradedCompIdeal(RL, w, 0, TOR)[2]; |
---|
3209 | |
---|
3210 | if(size(Iw) > 0){ |
---|
3211 | return(0); |
---|
3212 | } |
---|
3213 | |
---|
3214 | kill w, Iw; |
---|
3215 | } |
---|
3216 | |
---|
3217 | return(1); |
---|
3218 | } |
---|
3219 | |
---|
3220 | /////////////////////////// |
---|
3221 | |
---|
3222 | proc autGradAlg(ideal I, list #) |
---|
3223 | " |
---|
3224 | USAGE: autGradAlg(I, TOR); I is an ideal, TOR is an optional list of integers representing the torsion part of the grading group. |
---|
3225 | ASSUME: minimally presented, degrees of the generators of I |
---|
3226 | are the minimal degrees, basering multigraded pointedly, I_w = 0 |
---|
3227 | for all w = deg(var(i)) |
---|
3228 | RETURN: a ring. Also exports an ideal Jexported and a list stabExported. |
---|
3229 | EXAMPLE: example autGradAlg; shows an example |
---|
3230 | " |
---|
3231 | { |
---|
3232 | list TOR; |
---|
3233 | if(size(#) > 0 and typeof(#) == "list" and typeof(#[1]) == "int"){ |
---|
3234 | TOR = #[1]; |
---|
3235 | } |
---|
3236 | |
---|
3237 | // for each col w of Q |
---|
3238 | // we require I_w = 0: |
---|
3239 | if(!compsAreTrivial(I, TOR)){ |
---|
3240 | ERROR("For some i, w = deg(var(i)) did not satisfy I_w = 0. We do not support this case at the moment. Sorry."); |
---|
3241 | } |
---|
3242 | |
---|
3243 | def RR = stabilizer(I, TOR); |
---|
3244 | setring RR; |
---|
3245 | |
---|
3246 | // Jexported is already being exported |
---|
3247 | return(RR); |
---|
3248 | } |
---|
3249 | example |
---|
3250 | { |
---|
3251 | echo = 2; |
---|
3252 | |
---|
3253 | intmat Q[1][3] = |
---|
3254 | 1,1,1; |
---|
3255 | |
---|
3256 | ring R = 0,T(1..3), dp; |
---|
3257 | setBaseMultigrading(Q); |
---|
3258 | |
---|
3259 | ideal I = 0; //T(1)*T(2) + T(3)*T(4); |
---|
3260 | def RR = autGradAlg(I); |
---|
3261 | setring RR; |
---|
3262 | |
---|
3263 | "resulting ideal:"; |
---|
3264 | Jexported; |
---|
3265 | |
---|
3266 | "dimension:"; |
---|
3267 | dim(std(Jexported)); |
---|
3268 | |
---|
3269 | "as a detailed list:"; |
---|
3270 | stabExported; |
---|
3271 | } |
---|
3272 | |
---|
3273 | ////////////////////////// |
---|
3274 | |
---|
3275 | static proc shrink(ideal I) |
---|
3276 | " |
---|
3277 | USAGE: shrink(I): I is an ideal |
---|
3278 | ASSUME: There are variable generators in I and I does not contain 0-generators. |
---|
3279 | PURPOSE: returns a new ring S obtained from the old one by removing all variables Y(i) |
---|
3280 | such that some generator of I is Y(i). Exports the image of the given ideal in S. |
---|
3281 | RETURN: a ring. Exports an ideal Ishrink. |
---|
3282 | EXAMPLE: example shrink; shows an example |
---|
3283 | " |
---|
3284 | { |
---|
3285 | // create smaller ring: |
---|
3286 | // we will map Y(i) --> 0 if i is in the list zeros: |
---|
3287 | def R = basering; |
---|
3288 | ideal take = var(1); // initialize as the list of all variables |
---|
3289 | |
---|
3290 | for(int i = 2; i <= nvars(R); i++){ |
---|
3291 | take = take, var(i); |
---|
3292 | } |
---|
3293 | |
---|
3294 | list zeros; |
---|
3295 | int j; |
---|
3296 | i = 1; |
---|
3297 | int c = 1; |
---|
3298 | |
---|
3299 | while(c <= size(I)){ // i may have 0-entries: |
---|
3300 | for(j = 1; j <= nvars(R); j++){ |
---|
3301 | if(I[i] == var(j)){ |
---|
3302 | zeros[size(zeros)+1] = j; |
---|
3303 | break; |
---|
3304 | } |
---|
3305 | } |
---|
3306 | |
---|
3307 | if(I[i] != 0){ |
---|
3308 | c++; |
---|
3309 | } |
---|
3310 | |
---|
3311 | i++; |
---|
3312 | } |
---|
3313 | |
---|
3314 | for(i = 1; i <= size(zeros); i++){ |
---|
3315 | take[zeros[i]] = 0; |
---|
3316 | } |
---|
3317 | |
---|
3318 | string s = "("; |
---|
3319 | int firstone = 1; |
---|
3320 | int n = 0; // no of variables in new ring |
---|
3321 | |
---|
3322 | for(i = 1; i <= nvars(R); i++){ |
---|
3323 | if(take[i] != 0){ |
---|
3324 | if(firstone){ |
---|
3325 | s = s + string(take[i]); |
---|
3326 | firstone = 0; |
---|
3327 | n++; |
---|
3328 | } else { |
---|
3329 | s = s + ", " + string(take[i]); |
---|
3330 | n++; |
---|
3331 | } |
---|
3332 | } |
---|
3333 | } |
---|
3334 | |
---|
3335 | s = s + ")"; |
---|
3336 | ring S = create_ring(ring_list(R)[1], s, "(dp(" + string(n-1) + "), dp(1))", "no_minpoly"); //charstr: e.g. 0,a |
---|
3337 | |
---|
3338 | // map the ideal to Small: |
---|
3339 | ideal Ishrink = imap(R, I); |
---|
3340 | |
---|
3341 | export Ishrink; |
---|
3342 | return(S); |
---|
3343 | } |
---|
3344 | example |
---|
3345 | { |
---|
3346 | echo=2; |
---|
3347 | |
---|
3348 | ring R = 0,(T(1..10),Y(1..3),Z),dp; |
---|
3349 | ideal I = |
---|
3350 | T(1)*T(3), |
---|
3351 | Y(1), |
---|
3352 | T(2)*Y(1), |
---|
3353 | Y(3), |
---|
3354 | Z*T(7)-7*Y(3), |
---|
3355 | T(8); |
---|
3356 | |
---|
3357 | def S = shrink(I); |
---|
3358 | setring S; |
---|
3359 | |
---|
3360 | Ishrink; |
---|
3361 | |
---|
3362 | } |
---|