1 | /////////////////////////////////////////////////////////////////////////// |
---|
2 | version="version sheafcoh.lib 4.1.1.0 Dec_2017 "; // $Id$ |
---|
3 | category="Commutative Algebra"; |
---|
4 | info=" |
---|
5 | LIBRARY: sheafcoh.lib Procedures for Computing Sheaf Cohomology |
---|
6 | AUTHORS: Wolfram Decker, decker@mathematik.uni-kl.de |
---|
7 | Christoph Lossen, lossen@mathematik.uni-kl.de |
---|
8 | Gerhard Pfister, pfister@mathematik.uni-kl.de |
---|
9 | Oleksandr Motsak, U@D, where U={motsak}, D={mathematik.uni-kl.de} |
---|
10 | |
---|
11 | PROCEDURES: |
---|
12 | truncate(phi,d); truncation of coker(phi) at d |
---|
13 | truncateFast(phi,d); truncation of coker(phi) at d (fast+ experimental) |
---|
14 | CM_regularity(M); Castelnuovo-Mumford regularity of coker(M) |
---|
15 | sheafCohBGG(M,l,h); cohomology of sheaf associated to coker(M) |
---|
16 | sheafCohBGG2(M,l,h); cohomology of sheaf associated to coker(M), experimental version |
---|
17 | sheafCoh(M,l,h); cohomology of sheaf associated to coker(M) |
---|
18 | dimH(i,M,d); compute h^i(F(d)), F sheaf associated to coker(M) |
---|
19 | dimGradedPart() |
---|
20 | |
---|
21 | displayCohom(B,l,h,n); display intmat as Betti diagram (with zero rows) |
---|
22 | |
---|
23 | KEYWORDS: sheaf cohomology |
---|
24 | "; |
---|
25 | |
---|
26 | /////////////////////////////////////////////////////////////////////////////// |
---|
27 | LIB "matrix.lib"; |
---|
28 | LIB "nctools.lib"; |
---|
29 | LIB "homolog.lib"; |
---|
30 | |
---|
31 | |
---|
32 | /////////////////////////////////////////////////////////////////////////////// |
---|
33 | static proc jacobM(matrix M) |
---|
34 | { |
---|
35 | int n=nvars(basering); |
---|
36 | matrix B=transpose(diff(M,var(1))); |
---|
37 | int i; |
---|
38 | for(i=2;i<=n;i++) |
---|
39 | { |
---|
40 | B=concat(B,transpose(diff(M,var(i)))); |
---|
41 | } |
---|
42 | return(transpose(B)); |
---|
43 | } |
---|
44 | |
---|
45 | |
---|
46 | /////////////////////////////////////////////////////////////////////////////// |
---|
47 | /** |
---|
48 | let M = { w_1, ..., w_k }, k = size(M) == ncols(M), n = nvars(currRing). |
---|
49 | assuming that nrows(M) <= m*n; |
---|
50 | computes transpose(M) * transpose( var(1) I_m | ... | var(n) I_m ) :== transpose(module{f_1, ... f_k}), |
---|
51 | where f_i = \sum_{j=1}^{m} (w_i, v_j) gen(j), (w_i, v_j) is a `scalar` multiplication. |
---|
52 | that is, if w_i = (a^1_1, ... a^1_m) | (a^2_1, ..., a^2_m) | ... | (a^n_1, ..., a^n_m) then |
---|
53 | |
---|
54 | (a^1_1, ... a^1_m) | (a^2_1, ..., a^2_m) | ... | (a^n_1, ..., a^n_m) |
---|
55 | * var_1 ... var_1 | var_2 ... var_2 | ... | var_n ... var(n) |
---|
56 | * gen_1 ... gen_m | gen_1 ... gen_m | ... | gen_1 ... gen_m |
---|
57 | + => |
---|
58 | f_i = |
---|
59 | |
---|
60 | a^1_1 * var(1) * gen(1) + ... + a^1_m * var(1) * gen(m) + |
---|
61 | a^2_1 * var(2) * gen(1) + ... + a^2_m * var(2) * gen(m) + |
---|
62 | ... |
---|
63 | a^n_1 * var(n) * gen(1) + ... + a^n_m * var(n) * gen(m); |
---|
64 | |
---|
65 | NOTE: for every f_i we run only ONCE along w_i saving partial sums into a temporary array of polys of size m |
---|
66 | |
---|
67 | */ |
---|
68 | static proc TensorModuleMult(int m, module M) |
---|
69 | { |
---|
70 | return( system("tensorModuleMult", m, M) ); // trick! |
---|
71 | |
---|
72 | int n = nvars(basering); |
---|
73 | int k = ncols(M); |
---|
74 | |
---|
75 | int g, cc, vv; |
---|
76 | |
---|
77 | poly h; |
---|
78 | |
---|
79 | module Temp; // = {f_1, ..., f_k } |
---|
80 | |
---|
81 | intvec exp; |
---|
82 | vector pTempSum, w; |
---|
83 | |
---|
84 | for( int i = k; i > 0; i-- ) // for every w \in M |
---|
85 | { |
---|
86 | pTempSum[m] = 0; |
---|
87 | w = M[i]; |
---|
88 | |
---|
89 | while(w != 0) // for each term of w... |
---|
90 | { |
---|
91 | exp = leadexp(w); |
---|
92 | g = exp[n+1]; // module component! |
---|
93 | h = w[g]; |
---|
94 | |
---|
95 | w = w - h * gen(g); |
---|
96 | |
---|
97 | cc = g % m; |
---|
98 | |
---|
99 | if( cc == 0) |
---|
100 | { |
---|
101 | cc = m; |
---|
102 | } |
---|
103 | |
---|
104 | vv = 1 + (g - cc) / m; |
---|
105 | |
---|
106 | pTempSum = pTempSum + h * var(vv) * gen(cc); |
---|
107 | } |
---|
108 | |
---|
109 | Temp[i] = pTempSum; |
---|
110 | } |
---|
111 | |
---|
112 | Temp = transpose(Temp); |
---|
113 | |
---|
114 | return(Temp); |
---|
115 | } |
---|
116 | |
---|
117 | /////////////////////////////////////////////////////////////////////////////// |
---|
118 | proc truncate(module phi, int d) |
---|
119 | "USAGE: truncate(M,d); M module, d int |
---|
120 | ASSUME: @code{M} is graded, and it comes assigned with an admissible degree |
---|
121 | vector as an attribute |
---|
122 | RETURN: module |
---|
123 | NOTE: Output is a presentation matrix for the truncation of coker(M) |
---|
124 | at degree d. |
---|
125 | EXAMPLE: example truncate; shows an example |
---|
126 | KEYWORDS: truncated module |
---|
127 | " |
---|
128 | { |
---|
129 | if ( typeof(attrib(phi,"isHomog"))=="string" ) { |
---|
130 | if (size(phi)==0) { |
---|
131 | // assign weights 0 to generators of R^n (n=nrows(phi)) |
---|
132 | intvec v; |
---|
133 | v[nrows(phi)]=0; |
---|
134 | attrib(phi,"isHomog",v); |
---|
135 | } |
---|
136 | else { |
---|
137 | ERROR("No admissible degree vector assigned"); |
---|
138 | } |
---|
139 | } |
---|
140 | else { |
---|
141 | intvec v=attrib(phi,"isHomog"); |
---|
142 | } |
---|
143 | int i,m,dummy; |
---|
144 | int s = nrows(phi); |
---|
145 | module L; // TOO BIG!!! |
---|
146 | for (i=1; i<=s; i++) { |
---|
147 | if (d>v[i]) { |
---|
148 | L = L+maxideal(d-v[i])*gen(i); |
---|
149 | } |
---|
150 | else { |
---|
151 | L = L+gen(i); |
---|
152 | } |
---|
153 | } |
---|
154 | L = modulo(L,phi); |
---|
155 | L = minbase(prune(L)); |
---|
156 | if (size(L)==0) {return(L);} |
---|
157 | |
---|
158 | // it only remains to set the degrees for L: |
---|
159 | // ------------------------------------------ |
---|
160 | m = v[1]; |
---|
161 | for(i=2; i<=size(v); i++) { if(v[i]<m) { m = v[i]; } } |
---|
162 | dummy = homog(L); |
---|
163 | intvec vv = attrib(L,"isHomog"); |
---|
164 | if (d>m) { vv = vv+d; } |
---|
165 | else { vv = vv+m; } |
---|
166 | attrib(L,"isHomog",vv); |
---|
167 | return(L); |
---|
168 | } |
---|
169 | example |
---|
170 | {"EXAMPLE:"; |
---|
171 | echo = 2; |
---|
172 | ring R=0,(x,y,z),dp; |
---|
173 | module M=maxideal(3); |
---|
174 | homog(M); |
---|
175 | // compute presentation matrix for truncated module (R/<x,y,z>^3)_(>=2) |
---|
176 | module M2=truncate(M,2); |
---|
177 | print(M2); |
---|
178 | dimGradedPart(M2,1); |
---|
179 | dimGradedPart(M2,2); |
---|
180 | // this should coincide with: |
---|
181 | dimGradedPart(M,2); |
---|
182 | // shift grading by 1: |
---|
183 | intvec v=1; |
---|
184 | attrib(M,"isHomog",v); |
---|
185 | M2=truncate(M,2); |
---|
186 | print(M2); |
---|
187 | dimGradedPart(M2,3); |
---|
188 | } |
---|
189 | |
---|
190 | /////////////////////////////////////////////////////////////////////////////// |
---|
191 | |
---|
192 | proc truncateFast(module M, int d) |
---|
193 | "USAGE: truncateFast(M,d); M module, d int |
---|
194 | ASSUME: @code{M} is graded, and it comes assigned with an admissible degree |
---|
195 | vector as an attribute 'isHomog' |
---|
196 | RETURN: module |
---|
197 | NOTE: Output is a presentation matrix for the truncation of coker(M) |
---|
198 | at d. |
---|
199 | Fast + experimental version. M shoud be a SB! |
---|
200 | DISPLAY: If @code{printlevel}>=1, step-by step timings will be printed. |
---|
201 | If @code{printlevel}>=2 we add progress debug messages |
---|
202 | if @code{printlevel}>=3, even all intermediate results... |
---|
203 | EXAMPLE: example truncateFast; shows an example |
---|
204 | KEYWORDS: truncated module |
---|
205 | " |
---|
206 | { |
---|
207 | // int PL = printlevel + 1; |
---|
208 | int PL = printlevel - voice + 2; |
---|
209 | |
---|
210 | dbprint(PL-1, "// truncateFast(M: "+ string(nrows(M)) + " x " + string(ncols(M)) +", " + string(d) + "):"); |
---|
211 | dbprint(PL-2, M); |
---|
212 | |
---|
213 | intvec save = option(get); |
---|
214 | if( PL >= 2 ) |
---|
215 | { |
---|
216 | option(prot); |
---|
217 | option(mem); |
---|
218 | } |
---|
219 | |
---|
220 | int tTruncateBegin=timer; |
---|
221 | |
---|
222 | if (attrib(M,"isSB")!=1) |
---|
223 | { |
---|
224 | ERROR("M must be a standard basis!"); |
---|
225 | } |
---|
226 | |
---|
227 | dbprint(PL-1, "// M is a SB! "); |
---|
228 | |
---|
229 | if ( typeof(attrib(M,"isHomog"))=="string" ) { |
---|
230 | if (size(M)==0) { |
---|
231 | // assign weights 0 to generators of R^n (n=nrows(M)) |
---|
232 | intvec v; |
---|
233 | v[nrows(M)]=0; |
---|
234 | attrib(M,"isHomog",v); |
---|
235 | } |
---|
236 | else { |
---|
237 | ERROR("No admissible degree vector assigned"); |
---|
238 | } |
---|
239 | } |
---|
240 | else { |
---|
241 | intvec v=attrib(M,"isHomog"); |
---|
242 | } |
---|
243 | |
---|
244 | dbprint(PL-1, "// weighting(M): ["+ string(v) + "]"); |
---|
245 | |
---|
246 | int i,m,dummy; |
---|
247 | int s = nrows(M); |
---|
248 | |
---|
249 | int tKBaseBegin = timer; |
---|
250 | module L = kbase(M, d); // TODO: check whether this is always correct!?! |
---|
251 | |
---|
252 | |
---|
253 | dbprint(PL-1, "// L = kbase(M,d): "+string(nrows(L)) + " x " + string(ncols(L)) +""); |
---|
254 | dbprint(PL-2, L); |
---|
255 | dbprint(PL-1, "// weighting(L): ["+ string(attrib(L, "isHomog")) + "]"); |
---|
256 | |
---|
257 | int tModuloBegin = timer; |
---|
258 | L = modulo(L,M); |
---|
259 | |
---|
260 | dbprint(PL-1, "// L = modulo(L,M): "+string(nrows(L)) + " x " + string(ncols(L)) +""); |
---|
261 | dbprint(PL-2, L); |
---|
262 | dbprint(PL-1, "// weighting(L): ["+ string(attrib(L, "isHomog")) + "]"); |
---|
263 | |
---|
264 | int tPruneBegin = timer; |
---|
265 | L = prune(L); |
---|
266 | |
---|
267 | dbprint(PL-1, "// L = prune(L): "+string(nrows(L)) + " x " + string(ncols(L)) +""); |
---|
268 | dbprint(PL-2, L); |
---|
269 | dbprint(PL-1, "// weighting(L): ["+ string(attrib(L, "isHomog")) + "]"); |
---|
270 | |
---|
271 | int tPruneEnd = timer; |
---|
272 | L = minbase(L); |
---|
273 | int tMinBaseEnd = timer; |
---|
274 | |
---|
275 | dbprint(PL-1, "// L = minbase(L): "+string(nrows(L)) + " x " + string(ncols(L)) +""); |
---|
276 | dbprint(PL-2, L); |
---|
277 | dbprint(PL-1, "// weighting(L): ["+ string(attrib(L, "isHomog")) + "]"); |
---|
278 | |
---|
279 | |
---|
280 | |
---|
281 | |
---|
282 | if (size(L)!=0) |
---|
283 | { |
---|
284 | |
---|
285 | // it only remains to set the degrees for L: |
---|
286 | // ------------------------------------------ |
---|
287 | m = v[1]; |
---|
288 | for(i=2; i<=size(v); i++) { if(v[i]<m) { m = v[i]; } } |
---|
289 | dummy = homog(L); |
---|
290 | intvec vv = attrib(L,"isHomog"); |
---|
291 | if (d>m) { vv = vv+d; } |
---|
292 | else { vv = vv+m; } |
---|
293 | attrib(L,"isHomog",vv); |
---|
294 | } |
---|
295 | |
---|
296 | int tTruncateEnd=timer; |
---|
297 | |
---|
298 | dbprint(PL-1, "// corrected weighting(L): ["+ string(attrib(L, "isHomog")) + "]"); |
---|
299 | |
---|
300 | |
---|
301 | if(PL > 0) |
---|
302 | { |
---|
303 | " |
---|
304 | -------------- TIMINGS -------------- |
---|
305 | Trunc Time: ", tTruncateEnd - tTruncateBegin, " |
---|
306 | :: Before .Time: ", tKBaseBegin - tTruncateBegin, " |
---|
307 | :: kBase Time: ", tModuloBegin - tKBaseBegin, " |
---|
308 | :: Modulo Time: ", tPruneBegin - tModuloBegin, " |
---|
309 | :: Prune Time: ", tPruneEnd - tPruneBegin, " |
---|
310 | :: Minbase Time: ", tMinBaseEnd - tPruneEnd, " |
---|
311 | :: After .Time: ", tTruncateEnd - tMinBaseEnd; |
---|
312 | } |
---|
313 | |
---|
314 | option(set, save); |
---|
315 | |
---|
316 | return(L); |
---|
317 | } |
---|
318 | example |
---|
319 | {"EXAMPLE:"; |
---|
320 | echo = 2; |
---|
321 | ring R=0,(x,y,z,u,v),dp; |
---|
322 | module M=maxideal(3); |
---|
323 | homog(M); |
---|
324 | // compute presentation matrix for truncated module (R/<x,y,z,u>^3)_(>=2) |
---|
325 | int t=timer; |
---|
326 | module M2t=truncate(M,2); |
---|
327 | t = timer - t; |
---|
328 | "// Simple truncate: ", t; |
---|
329 | t=timer; |
---|
330 | module M2=truncateFast(std(M),2); |
---|
331 | t = timer - t; |
---|
332 | "// Fast truncate: ", t; |
---|
333 | print(M2); |
---|
334 | "// Check: M2t == M2?: ", size(NF(M2, std(M2t))) + size(NF(M2t, std(M2))); |
---|
335 | |
---|
336 | dimGradedPart(M2,1); |
---|
337 | dimGradedPart(M2,2); |
---|
338 | // this should coincide with: |
---|
339 | dimGradedPart(M,2); |
---|
340 | // shift grading by 1: |
---|
341 | intvec v=1; |
---|
342 | attrib(M,"isHomog",v); |
---|
343 | t=timer; |
---|
344 | M2t=truncate(M,2); |
---|
345 | t = timer - t; |
---|
346 | "// Simple truncate: ", t; |
---|
347 | t=timer; |
---|
348 | M2=truncateFast(std(M),2); |
---|
349 | t = timer - t; |
---|
350 | "// Fast truncate: ", t; |
---|
351 | print(M2); |
---|
352 | "// Check: M2t == M2?: ", size(NF(M2, std(M2t))) + size(NF(M2t, std(M2))); //? |
---|
353 | dimGradedPart(M2,3); |
---|
354 | } |
---|
355 | |
---|
356 | /////////////////////////////////////////////////////////////////////////////// |
---|
357 | |
---|
358 | |
---|
359 | |
---|
360 | |
---|
361 | |
---|
362 | proc dimGradedPart(module phi, int d) |
---|
363 | "USAGE: dimGradedPart(M,d); M module, d int |
---|
364 | ASSUME: @code{M} is graded, and it comes assigned with an admissible degree |
---|
365 | vector as an attribute |
---|
366 | RETURN: int |
---|
367 | NOTE: Output is the vector space dimension of the graded part of degree d |
---|
368 | of coker(M). |
---|
369 | EXAMPLE: example dimGradedPart; shows an example |
---|
370 | KEYWORDS: graded module, graded piece |
---|
371 | " |
---|
372 | { |
---|
373 | if ( typeof(attrib(phi,"isHomog"))=="string" ) { |
---|
374 | if (size(phi)==0) { |
---|
375 | // assign weights 0 to generators of R^n (n=nrows(phi)) |
---|
376 | intvec v; |
---|
377 | v[nrows(phi)]=0; |
---|
378 | } |
---|
379 | else { ERROR("No admissible degree vector assigned"); } |
---|
380 | } |
---|
381 | else { |
---|
382 | intvec v=attrib(phi,"isHomog"); |
---|
383 | } |
---|
384 | int s = nrows(phi); |
---|
385 | int i,m,dummy; |
---|
386 | module L,LL; |
---|
387 | for (i=1; i<=s; i++) { |
---|
388 | if (d>v[i]) { |
---|
389 | L = L+maxideal(d-v[i])*gen(i); |
---|
390 | LL = LL+maxideal(d+1-v[i])*gen(i); |
---|
391 | } |
---|
392 | else { |
---|
393 | L = L+gen(i); |
---|
394 | if (d==v[i]) { |
---|
395 | LL = LL+maxideal(1)*gen(i); |
---|
396 | } |
---|
397 | else { |
---|
398 | LL = LL+gen(i); |
---|
399 | } |
---|
400 | } |
---|
401 | } |
---|
402 | LL=LL,phi; |
---|
403 | L = modulo(L,LL); |
---|
404 | L = std(prune(L)); |
---|
405 | if (size(L)==0) {return(0);} |
---|
406 | return(vdim(L)); |
---|
407 | } |
---|
408 | example |
---|
409 | {"EXAMPLE:"; |
---|
410 | echo = 2; |
---|
411 | ring R=0,(x,y,z),dp; |
---|
412 | module M=maxideal(3); |
---|
413 | // assign compatible weight vector (here: 0) |
---|
414 | homog(M); |
---|
415 | // compute dimension of graded pieces of R/<x,y,z>^3 : |
---|
416 | dimGradedPart(M,0); |
---|
417 | dimGradedPart(M,1); |
---|
418 | dimGradedPart(M,2); |
---|
419 | dimGradedPart(M,3); |
---|
420 | // shift grading: |
---|
421 | attrib(M,"isHomog",intvec(2)); |
---|
422 | dimGradedPart(M,2); |
---|
423 | } |
---|
424 | |
---|
425 | /////////////////////////////////////////////////////////////////////////////// |
---|
426 | |
---|
427 | proc CM_regularity (module M) |
---|
428 | "USAGE: CM_regularity(M); M module |
---|
429 | ASSUME: @code{M} is graded, and it comes assigned with an admissible degree |
---|
430 | vector as an attribute |
---|
431 | RETURN: integer, the Castelnuovo-Mumford regularity of coker(M) |
---|
432 | NOTE: procedure calls mres |
---|
433 | EXAMPLE: example CM_regularity; shows an example |
---|
434 | KEYWORDS: Castelnuovo-Mumford regularity |
---|
435 | " |
---|
436 | { |
---|
437 | if ( typeof(attrib(M,"isHomog"))=="string" ) { |
---|
438 | if (size(M)==0) { |
---|
439 | // assign weights 0 to generators of R^n (n=nrows(M)) |
---|
440 | intvec v; |
---|
441 | v[nrows(M)]=0; |
---|
442 | attrib(M,"isHomog",v); |
---|
443 | } |
---|
444 | else { |
---|
445 | ERROR("No admissible degree vector assigned"); |
---|
446 | } |
---|
447 | } |
---|
448 | |
---|
449 | if( attrib(CM_regularity,"Algorithm") == "minres_res" ) |
---|
450 | { |
---|
451 | def L = minres( res(M,0) ); // let's try it out! |
---|
452 | } else |
---|
453 | { |
---|
454 | def L = mres(M,0); |
---|
455 | } |
---|
456 | |
---|
457 | intmat BeL = betti(L); |
---|
458 | int r = nrows(module(matrix(BeL))); // last non-zero row |
---|
459 | if (typeof(attrib(BeL,"rowShift"))!="string") { |
---|
460 | int shift = attrib(BeL,"rowShift"); |
---|
461 | } |
---|
462 | return(r+shift-1); |
---|
463 | } |
---|
464 | example |
---|
465 | {"EXAMPLE:"; |
---|
466 | echo = 2; |
---|
467 | ring R=0,(x,y,z,u),dp; |
---|
468 | resolution T1=mres(maxideal(1),0); |
---|
469 | module M=T1[3]; |
---|
470 | intvec v=2,2,2,2,2,2; |
---|
471 | attrib(M,"isHomog",v); |
---|
472 | CM_regularity(M); |
---|
473 | } |
---|
474 | |
---|
475 | /////////////////////////////////////////////////////////////////////////////// |
---|
476 | proc sheafCohBGG(module M,int l,int h) |
---|
477 | "USAGE: sheafCohBGG(M,l,h); M module, l,h int |
---|
478 | ASSUME: @code{M} is graded, and it comes assigned with an admissible degree |
---|
479 | vector as an attribute, @code{h>=l}, and the basering has @code{n+1} |
---|
480 | variables. |
---|
481 | RETURN: intmat, cohomology of twists of the coherent sheaf F on P^n |
---|
482 | associated to coker(M). The range of twists is determined by @code{l}, |
---|
483 | @code{h}. |
---|
484 | DISPLAY: The intmat is displayed in a diagram of the following form: @* |
---|
485 | @format |
---|
486 | l l+1 h |
---|
487 | ---------------------------------------------------------- |
---|
488 | n: h^n(F(l)) h^n(F(l+1)) ...... h^n(F(h)) |
---|
489 | ............................................... |
---|
490 | 1: h^1(F(l)) h^1(F(l+1)) ...... h^1(F(h)) |
---|
491 | 0: h^0(F(l)) h^0(F(l+1)) ...... h^0(F(h)) |
---|
492 | ---------------------------------------------------------- |
---|
493 | chi: chi(F(l)) chi(F(l+1)) ...... chi(F(h)) |
---|
494 | @end format |
---|
495 | A @code{'-'} in the diagram refers to a zero entry; a @code{'*'} |
---|
496 | refers to a negative entry (= dimension not yet determined). |
---|
497 | refers to a not computed dimension. @* |
---|
498 | NOTE: This procedure is based on the Bernstein-Gel'fand-Gel'fand |
---|
499 | correspondence and on Tate resolution ( see [Eisenbud, Floystad, |
---|
500 | Schreyer: Sheaf cohomology and free resolutions over exterior |
---|
501 | algebras, Trans AMS 355 (2003)] ).@* |
---|
502 | @code{sheafCohBGG(M,l,h)} does not compute all values in the above |
---|
503 | table. To determine all values of @code{h^i(F(d))}, @code{d=l..h}, |
---|
504 | use @code{sheafCohBGG(M,l-n,h+n)}. |
---|
505 | SEE ALSO: sheafCoh, dimH |
---|
506 | EXAMPLE: example sheafCohBGG; shows an example |
---|
507 | " |
---|
508 | { |
---|
509 | int i,j,k,row,col; |
---|
510 | if( typeof(attrib(M,"isHomog"))!="intvec" ) |
---|
511 | { |
---|
512 | if (size(M)==0) { attrib(M,"isHomog",0); } |
---|
513 | else { ERROR("No admissible degree vector assigned"); } |
---|
514 | } |
---|
515 | int n=nvars(basering)-1; |
---|
516 | int ell=l+n; |
---|
517 | def R=basering; |
---|
518 | int reg = CM_regularity(M); |
---|
519 | int bound=max(reg+1,h-1); |
---|
520 | module MT=truncate(M,bound); |
---|
521 | int m=nrows(MT); |
---|
522 | MT=transpose(jacobM(MT)); |
---|
523 | MT=syz(MT); |
---|
524 | matrix ML[n+1][1]=maxideal(1); |
---|
525 | matrix S=transpose(outer(ML,unitmat(m))); |
---|
526 | matrix SS=transpose(S*MT); |
---|
527 | //--- to the exterior algebra |
---|
528 | def AR = Exterior(); |
---|
529 | setring AR; |
---|
530 | intvec saveopt=option(get); |
---|
531 | option(redSB); |
---|
532 | option(redTail); |
---|
533 | module EM=imap(R,SS); |
---|
534 | intvec w; |
---|
535 | //--- here we are with our matrix |
---|
536 | int bound1=max(1,bound-ell+1); |
---|
537 | for (i=1; i<=nrows(EM); i++) |
---|
538 | { |
---|
539 | w[i]=-bound-1; |
---|
540 | } |
---|
541 | attrib(EM,"isHomog",w); |
---|
542 | resolution RE=mres(EM,bound1); |
---|
543 | intmat Betti=betti(RE); |
---|
544 | k=ncols(Betti); |
---|
545 | row=nrows(Betti); |
---|
546 | int shift=attrib(Betti,"rowShift")+(k+ell-1); |
---|
547 | intmat newBetti[n+1][h-l+1]; |
---|
548 | for (j=1; j<=row; j++) |
---|
549 | { |
---|
550 | for (i=l; i<=h; i++) |
---|
551 | { |
---|
552 | if ((k+1-j-i+ell-shift>0) and (j+i-ell+shift>=1)) |
---|
553 | { |
---|
554 | newBetti[n+2-shift-j,i-l+1]=Betti[j,k+1-j-i+ell-shift]; |
---|
555 | } |
---|
556 | else { newBetti[n+2-shift-j,i-l+1]=-1; } |
---|
557 | } |
---|
558 | } |
---|
559 | for (j=2; j<=n+1; j++) |
---|
560 | { |
---|
561 | for (i=1; i<j; i++) |
---|
562 | { |
---|
563 | newBetti[j,i]=-1; |
---|
564 | } |
---|
565 | } |
---|
566 | int d=k-h+ell-1; |
---|
567 | for (j=1; j<=n; j++) |
---|
568 | { |
---|
569 | for (i=h-l+1; i>=k+j; i--) |
---|
570 | { |
---|
571 | newBetti[j,i]=-1; |
---|
572 | } |
---|
573 | } |
---|
574 | displayCohom(newBetti,l,h,n); |
---|
575 | option(set,saveopt); |
---|
576 | setring R; |
---|
577 | return(newBetti); |
---|
578 | } |
---|
579 | example |
---|
580 | {"EXAMPLE:"; |
---|
581 | echo = 2; |
---|
582 | // cohomology of structure sheaf on P^4: |
---|
583 | //------------------------------------------- |
---|
584 | ring r=0,x(1..5),dp; |
---|
585 | module M=0; |
---|
586 | def A=sheafCohBGG(M,-9,4); |
---|
587 | // cohomology of cotangential bundle on P^3: |
---|
588 | //------------------------------------------- |
---|
589 | ring R=0,(x,y,z,u),dp; |
---|
590 | resolution T1=mres(maxideal(1),0); |
---|
591 | module M=T1[3]; |
---|
592 | intvec v=2,2,2,2,2,2; |
---|
593 | attrib(M,"isHomog",v); |
---|
594 | def B=sheafCohBGG(M,-8,4); |
---|
595 | } |
---|
596 | |
---|
597 | |
---|
598 | /////////////////////////////////////////////////////////////////////////////// |
---|
599 | static proc showResult( def R, int l, int h ) |
---|
600 | { |
---|
601 | int PL = 1; // printlevel - voice + 2; |
---|
602 | // int PL = printlevel + 1; |
---|
603 | |
---|
604 | intmat Betti; |
---|
605 | if(typeof(R)=="resolution") |
---|
606 | { |
---|
607 | Betti = betti(R); |
---|
608 | } else |
---|
609 | { |
---|
610 | if(typeof(R)!="intmat") |
---|
611 | { |
---|
612 | ERROR("Wrong input!!!"); |
---|
613 | } |
---|
614 | |
---|
615 | Betti = R; |
---|
616 | } |
---|
617 | |
---|
618 | int n=nvars(basering)-1; |
---|
619 | int ell = l + n; |
---|
620 | |
---|
621 | int k = ncols(Betti); |
---|
622 | int row = nrows(Betti); |
---|
623 | int shift = attrib(Betti,"rowShift") + (k + ell - 1); |
---|
624 | |
---|
625 | int iWTH = h-l+1; |
---|
626 | |
---|
627 | int d = k - h + ell - 1; |
---|
628 | |
---|
629 | if( PL > 1 ) |
---|
630 | { |
---|
631 | "// l: ", l; |
---|
632 | "// h: ", h; |
---|
633 | "// n: ", n; |
---|
634 | "// ell: ", ell; |
---|
635 | "// k: ", k; |
---|
636 | "// row: ", row; |
---|
637 | "// shift: ", shift; |
---|
638 | "// iWTH: ", iWTH; |
---|
639 | "// d: ", d; |
---|
640 | } |
---|
641 | |
---|
642 | intmat newBetti[ n + 1 ][ iWTH ]; |
---|
643 | int i, j; |
---|
644 | |
---|
645 | for (j=1; j<=row; j++) { |
---|
646 | for (i=l; i<=h; i++) { |
---|
647 | if( (n+2-shift-j)>0 ) { |
---|
648 | |
---|
649 | if ( (k+1-j-i+ell-shift>0) and (j+i-ell+shift>=1)) { |
---|
650 | newBetti[n+2-shift-j,i-l+1]=Betti[j,k+1-j-i+ell-shift]; |
---|
651 | } |
---|
652 | else { newBetti[n+2-shift-j,i-l+1]=-1; } |
---|
653 | |
---|
654 | } |
---|
655 | } |
---|
656 | } |
---|
657 | |
---|
658 | for (j=2; j<=n+1; j++) { |
---|
659 | for (i=1; i<min(j, iWTH); i++) { |
---|
660 | newBetti[j,i]= -1; |
---|
661 | } |
---|
662 | } |
---|
663 | |
---|
664 | for (j=1; j<=n; j++) { |
---|
665 | for (i=iWTH; i>=k+j; i--) { |
---|
666 | newBetti[j,i]=0; // -1; |
---|
667 | } |
---|
668 | } |
---|
669 | |
---|
670 | if( PL > 0 ) |
---|
671 | { |
---|
672 | "Cohomology table:"; |
---|
673 | displayCohom(newBetti, l, h, n); |
---|
674 | } |
---|
675 | |
---|
676 | return(newBetti); |
---|
677 | } |
---|
678 | /////////////////////////////////////////////////////////////////////////////// |
---|
679 | |
---|
680 | |
---|
681 | |
---|
682 | /////////////////////////////////////////////////////////////////////////////// |
---|
683 | proc sheafCohBGG2(module M,int l,int h) |
---|
684 | "USAGE: sheafCohBGG2(M,l,h); M module, l,h int |
---|
685 | ASSUME: @code{M} is graded, and it comes assigned with an admissible degree |
---|
686 | vector as an attribute, @code{h>=l}, and the basering has @code{n+1} |
---|
687 | variables. |
---|
688 | RETURN: intmat, cohomology of twists of the coherent sheaf F on P^n |
---|
689 | associated to coker(M). The range of twists is determined by @code{l}, |
---|
690 | @code{h}. |
---|
691 | DISPLAY: The intmat is displayed in a diagram of the following form: @* |
---|
692 | @format |
---|
693 | l l+1 h |
---|
694 | ---------------------------------------------------------- |
---|
695 | n: h^n(F(l)) h^n(F(l+1)) ...... h^n(F(h)) |
---|
696 | ............................................... |
---|
697 | 1: h^1(F(l)) h^1(F(l+1)) ...... h^1(F(h)) |
---|
698 | 0: h^0(F(l)) h^0(F(l+1)) ...... h^0(F(h)) |
---|
699 | ---------------------------------------------------------- |
---|
700 | chi: chi(F(l)) chi(F(l+1)) ...... chi(F(h)) |
---|
701 | @end format |
---|
702 | A @code{'-'} in the diagram refers to a zero entry; a @code{'*'} |
---|
703 | refers to a negative entry (= dimension not yet determined). |
---|
704 | refers to a not computed dimension. @* |
---|
705 | If @code{printlevel}>=1, step-by step timings will be printed. |
---|
706 | If @code{printlevel}>=2 we add progress debug messages |
---|
707 | if @code{printlevel}>=3, even all intermediate results... |
---|
708 | NOTE: This procedure is based on the Bernstein-Gel'fand-Gel'fand |
---|
709 | correspondence and on Tate resolution ( see [Eisenbud, Floystad, |
---|
710 | Schreyer: Sheaf cohomology and free resolutions over exterior |
---|
711 | algebras, Trans AMS 355 (2003)] ).@* |
---|
712 | @code{sheafCohBGG(M,l,h)} does not compute all values in the above |
---|
713 | table. To determine all values of @code{h^i(F(d))}, @code{d=l..h}, |
---|
714 | use @code{sheafCohBGG(M,l-n,h+n)}. |
---|
715 | Experimental version. Should require less memory. |
---|
716 | SEE ALSO: sheafCohBGG |
---|
717 | EXAMPLE: example sheafCohBGG2; shows an example |
---|
718 | " |
---|
719 | { |
---|
720 | int PL = printlevel - voice + 2; |
---|
721 | // int PL = printlevel; |
---|
722 | |
---|
723 | dbprint(PL-1, "// sheafCohBGG2(M: "+ string(nrows(M)) + " x " + string(ncols(M)) +", " + string(l) + ", " + string(h) + "):"); |
---|
724 | dbprint(PL-2, M); |
---|
725 | |
---|
726 | intvec save = option(get); |
---|
727 | |
---|
728 | if( PL >= 2 ) |
---|
729 | { |
---|
730 | option(prot); |
---|
731 | option(mem); |
---|
732 | } |
---|
733 | |
---|
734 | def isCoker = attrib(M, "isCoker"); |
---|
735 | if( typeof(isCoker) == "int" ) |
---|
736 | { |
---|
737 | if( isCoker > 0 ) |
---|
738 | { |
---|
739 | dbprint(PL-1, "We are going to assume that M is given by coker matrix (that is, M is not a submodule presentation!)"); |
---|
740 | } |
---|
741 | } |
---|
742 | |
---|
743 | int i,j,k,row,col; |
---|
744 | |
---|
745 | if( typeof(attrib(M,"isHomog"))!="intvec" ) |
---|
746 | { |
---|
747 | if (size(M)==0) { attrib(M,"isHomog",0); } |
---|
748 | else { ERROR("No admissible degree vector assigned"); } |
---|
749 | } |
---|
750 | |
---|
751 | dbprint(PL-1, "// weighting(M): ["+ string(attrib(M, "isHomog")) + "]"); |
---|
752 | |
---|
753 | option(redSB); option(redTail); |
---|
754 | |
---|
755 | def R=basering; |
---|
756 | |
---|
757 | int n = nvars(R) - 1; |
---|
758 | int ell = l + n; |
---|
759 | |
---|
760 | |
---|
761 | ///////////////////////////////////////////////////////////////////////////// |
---|
762 | // computations |
---|
763 | |
---|
764 | int tBegin=timer; |
---|
765 | int reg = CM_regularity(M); |
---|
766 | int tCMEnd = timer; |
---|
767 | |
---|
768 | dbprint(PL-1, "// CM_reg(M): "+ string(reg)); |
---|
769 | |
---|
770 | int bound = max(reg + 1, h - 1); |
---|
771 | |
---|
772 | dbprint(PL-1, "// bound: "+ string(bound)); |
---|
773 | |
---|
774 | /////////////////////////////////////////////////////////////// |
---|
775 | int tSTDBegin=timer; |
---|
776 | M = std(M); // for kbase! // NOTE: this should be after CM_regularity, since otherwise CM_regularity computes JUST TOOOOOOO LONG sometimes (see Reg_Hard examples!) |
---|
777 | int tSTDEnd = timer; |
---|
778 | |
---|
779 | dbprint(PL-1, "// M = std(M: "+string(nrows(M)) + " x " + string(ncols(M)) + ")"); |
---|
780 | dbprint(PL-2, M); |
---|
781 | dbprint(PL-1, "// weighting(M): ["+ string(attrib(M, "isHomog")) + "]"); |
---|
782 | |
---|
783 | |
---|
784 | printlevel = printlevel + 1; |
---|
785 | int tTruncateBegin=timer; |
---|
786 | module MT = truncateFast(M, bound); |
---|
787 | int tTruncateEnd=timer; |
---|
788 | printlevel = printlevel - 1; |
---|
789 | |
---|
790 | dbprint(PL-1, "// MT = truncateFast(M: "+string(nrows(MT)) + " x " + string(ncols(MT)) +", " + string(bound) + ")"); |
---|
791 | dbprint(PL-2, MT); |
---|
792 | dbprint(PL-1, "// weighting(MT): ["+ string(attrib(MT, "isHomog")) + "]"); |
---|
793 | |
---|
794 | int m=nrows(MT); |
---|
795 | |
---|
796 | /////////////////////////////////////////////////////////////// |
---|
797 | int tTransposeJacobBegin=timer; |
---|
798 | MT = jacob(MT); // ! :( |
---|
799 | int tTransposeJacobEnd=timer; |
---|
800 | |
---|
801 | dbprint(PL-1, "// MT = jacob(MT: "+string(nrows(MT)) + " x " + string(ncols(MT)) + ")"); |
---|
802 | dbprint(PL-2, MT); |
---|
803 | dbprint(PL-1, "// weighting(MT): ["+ string(attrib(MT, "isHomog")) + "]"); |
---|
804 | |
---|
805 | int tSyzBegin=timer; |
---|
806 | MT = syz(MT); |
---|
807 | int tSyzEnd=timer; |
---|
808 | |
---|
809 | dbprint(PL-1, "// MT = syz(MT: "+string(nrows(MT)) + " x " + string(ncols(MT)) + ")"); |
---|
810 | dbprint(PL-2, MT); |
---|
811 | dbprint(PL-1, "// weighting(MT): ["+ string(attrib(MT, "isHomog")) + "]"); |
---|
812 | |
---|
813 | int tMatrixOppBegin=timer; |
---|
814 | module SS = TensorModuleMult(m, MT); |
---|
815 | int tMatrixOppEnd=timer; |
---|
816 | |
---|
817 | dbprint(PL-1, "// SS = TensorModuleMult("+ string(m)+ ", MT: "+string(nrows(MT)) + " x " + string(ncols(MT)) + ")"); |
---|
818 | dbprint(PL-2, SS); |
---|
819 | dbprint(PL-1, "// weighting(SS): ["+ string(attrib(SS, "isHomog")) + "]"); |
---|
820 | |
---|
821 | //--- to the exterior algebra |
---|
822 | def AR = Exterior(); setring AR; |
---|
823 | |
---|
824 | dbprint(PL-1, "// Test: var(1) * var(1): "+ string(var(1) * var(1))); |
---|
825 | |
---|
826 | int maxbound = max(1, bound - ell + 1); |
---|
827 | // int maxbound = max(1, bound - l + 1); // As In M2!!! |
---|
828 | |
---|
829 | dbprint(PL-1, "// maxbound: "+ string(maxbound)); |
---|
830 | |
---|
831 | //--- here we are with our matrix |
---|
832 | module EM=imap(R,SS); |
---|
833 | intvec w; |
---|
834 | for (i=1; i<=nrows(EM); i++) |
---|
835 | { |
---|
836 | w[i]=-bound-1; |
---|
837 | } |
---|
838 | |
---|
839 | attrib(EM,"isHomog",w); |
---|
840 | |
---|
841 | /////////////////////////////////////////////////////////////// |
---|
842 | |
---|
843 | dbprint(PL-1, "// EM: "+string(nrows(EM)) + " x " + string(ncols(EM)) + ")"); |
---|
844 | dbprint(PL-2, EM); |
---|
845 | dbprint(PL-1, "// weighting(EM): ["+ string(attrib(EM, "isHomog")) + "]"); |
---|
846 | |
---|
847 | int tResulutionBegin=timer; |
---|
848 | resolution RE = nres(EM, maxbound); // TODO: Plural computes one too many syzygies...?! |
---|
849 | int tMinResBegin=timer; |
---|
850 | RE = minres(RE); |
---|
851 | int tBettiBegin=timer; |
---|
852 | intmat Betti = betti(RE); // betti(RE, 1);? |
---|
853 | int tResulutionEnd=timer; |
---|
854 | |
---|
855 | int tEnd = tResulutionEnd; |
---|
856 | |
---|
857 | if( PL > 0 ) |
---|
858 | { |
---|
859 | // list L = RE; // TODO: size(L/RE) is wrong! |
---|
860 | " |
---|
861 | ---- RESULTS ---- |
---|
862 | Tate Resolution: |
---|
863 | "; |
---|
864 | RE; |
---|
865 | "Betti numbers for Tate resolution (diagonal cohomology table):"; |
---|
866 | print(Betti, "betti"); // Diagonal form! |
---|
867 | } |
---|
868 | |
---|
869 | // printlevel = printlevel + 1; |
---|
870 | Betti = showResult(Betti, l, h ); // Show usual form of cohomology table |
---|
871 | // printlevel = printlevel - 1; |
---|
872 | |
---|
873 | if(PL > 0) |
---|
874 | { |
---|
875 | " |
---|
876 | ---- TIMINGS ------- |
---|
877 | Trunc Time: ", tTruncateEnd - tTruncateBegin, " |
---|
878 | Reg Time: ", tCMEnd - tBegin, " |
---|
879 | kStd Time: ", tSTDEnd - tSTDBegin, " |
---|
880 | Jacob Time: ", tTransposeJacobEnd - tTransposeJacobBegin, " |
---|
881 | Syz Time: ", tSyzEnd - tSyzBegin, " |
---|
882 | Mat Time: ", tMatrixOppEnd - tMatrixOppBegin, " |
---|
883 | ------------------------------ |
---|
884 | Res Time: ", tResulutionEnd - tResulutionBegin, " |
---|
885 | :: NRes Time: ", tMinResBegin - tResulutionBegin, " |
---|
886 | :: MinRes .Time: ", tBettiBegin - tMinResBegin, " |
---|
887 | :: Betti .Time: ", tResulutionEnd - tBettiBegin, " |
---|
888 | --------------------------------------------------------- |
---|
889 | Total Time: ", tEnd - tBegin, " |
---|
890 | --------------------------------------------------------- |
---|
891 | "; |
---|
892 | } |
---|
893 | |
---|
894 | setring R; |
---|
895 | |
---|
896 | option(set, save); |
---|
897 | |
---|
898 | return(Betti); |
---|
899 | } |
---|
900 | example |
---|
901 | {"EXAMPLE:"; |
---|
902 | echo = 2; |
---|
903 | int pl = printlevel; |
---|
904 | int l,h, t; |
---|
905 | |
---|
906 | //------------------------------------------- |
---|
907 | // cohomology of structure sheaf on P^4: |
---|
908 | //------------------------------------------- |
---|
909 | ring r=32001,x(1..5),dp; |
---|
910 | |
---|
911 | module M= getStructureSheaf(); // OO_P^4 |
---|
912 | |
---|
913 | l = -12; h = 12; // range of twists: l..h |
---|
914 | |
---|
915 | printlevel = 0; |
---|
916 | ////////////////////////////////////////////// |
---|
917 | t = timer; |
---|
918 | |
---|
919 | def A = sheafCoh(M, l, h); // global Ext method: |
---|
920 | |
---|
921 | "Time: ", timer - t; |
---|
922 | ////////////////////////////////////////////// |
---|
923 | t = timer; |
---|
924 | |
---|
925 | A = sheafCohBGG(M, l, h); // BGG method (without optimization): |
---|
926 | |
---|
927 | "Time: ", timer - t; |
---|
928 | ////////////////////////////////////////////// |
---|
929 | t = timer; |
---|
930 | |
---|
931 | A = sheafCohBGG2(M, l, h); // BGG method (with optimization) |
---|
932 | |
---|
933 | "Time: ", timer - t; |
---|
934 | ////////////////////////////////////////////// |
---|
935 | printlevel = pl; |
---|
936 | |
---|
937 | kill A, r; |
---|
938 | |
---|
939 | //------------------------------------------- |
---|
940 | // cohomology of cotangential bundle on P^3: |
---|
941 | //------------------------------------------- |
---|
942 | ring R=32001,(x,y,z,u),dp; |
---|
943 | |
---|
944 | module M = getCotangentialBundle(); |
---|
945 | |
---|
946 | l = -12; h = 11; // range of twists: l..h |
---|
947 | |
---|
948 | |
---|
949 | ////////////////////////////////////////////// |
---|
950 | printlevel = 0; |
---|
951 | t = timer; |
---|
952 | |
---|
953 | def B = sheafCoh(M, l, h); // global Ext method: |
---|
954 | |
---|
955 | "Time: ", timer - t; |
---|
956 | ////////////////////////////////////////////// |
---|
957 | t = timer; |
---|
958 | |
---|
959 | B = sheafCohBGG(M, l, h); // BGG method (without optimization): |
---|
960 | |
---|
961 | "Time: ", timer - t; |
---|
962 | ////////////////////////////////////////////// |
---|
963 | t = timer; |
---|
964 | |
---|
965 | B = sheafCohBGG2(M, l, h); // BGG method (with optimization) |
---|
966 | |
---|
967 | "Time: ", timer - t; |
---|
968 | ////////////////////////////////////////////// |
---|
969 | printlevel = pl; |
---|
970 | } |
---|
971 | |
---|
972 | |
---|
973 | /////////////////////////////////////////////////////////////////////////////// |
---|
974 | |
---|
975 | proc dimH(int i,module M,int d) |
---|
976 | "USAGE: dimH(i,M,d); M module, i,d int |
---|
977 | ASSUME: @code{M} is graded, and it comes assigned with an admissible degree |
---|
978 | vector as an attribute, @code{h>=l}, and the basering @code{S} has |
---|
979 | @code{n+1} variables. |
---|
980 | RETURN: int, vector space dimension of @math{H^i(F(d))} for F the coherent |
---|
981 | sheaf on P^n associated to coker(M). |
---|
982 | NOTE: The procedure is based on local duality as described in [Eisenbud: |
---|
983 | Computing cohomology. In Vasconcelos: Computational methods in |
---|
984 | commutative algebra and algebraic geometry. Springer (1998)]. |
---|
985 | SEE ALSO: sheafCoh, sheafCohBGG |
---|
986 | EXAMPLE: example dimH; shows an example |
---|
987 | " |
---|
988 | { |
---|
989 | if( typeof(attrib(M,"isHomog"))=="string" ) |
---|
990 | { |
---|
991 | if (size(M)==0) |
---|
992 | { |
---|
993 | // assign weights 0 to generators of R^n (n=nrows(M)) |
---|
994 | intvec v; |
---|
995 | v[nrows(M)]=0; |
---|
996 | attrib(M,"isHomog",v); |
---|
997 | } |
---|
998 | else |
---|
999 | { |
---|
1000 | ERROR("No admissible degree vector assigned"); |
---|
1001 | } |
---|
1002 | } |
---|
1003 | int Result; |
---|
1004 | int n=nvars(basering)-1; |
---|
1005 | if ((i>0) and (i<=n)) { |
---|
1006 | list L=Ext_R(n-i,M,1)[2]; |
---|
1007 | def N=L[1]; |
---|
1008 | return(dimGradedPart(N,-n-1-d)); |
---|
1009 | } |
---|
1010 | else |
---|
1011 | { |
---|
1012 | if (i==0) |
---|
1013 | { |
---|
1014 | list L=Ext_R(intvec(n+1,n+2),M,1)[2]; |
---|
1015 | def N0=L[2]; |
---|
1016 | def N1=L[1]; |
---|
1017 | Result=dimGradedPart(M,d) - dimGradedPart(N0,-n-1-d) |
---|
1018 | - dimGradedPart(N1,-n-1-d); |
---|
1019 | return(Result); |
---|
1020 | } |
---|
1021 | else { |
---|
1022 | return(0); |
---|
1023 | } |
---|
1024 | } |
---|
1025 | } |
---|
1026 | example |
---|
1027 | {"EXAMPLE:"; |
---|
1028 | echo = 2; |
---|
1029 | ring R=0,(x,y,z,u),dp; |
---|
1030 | resolution T1=mres(maxideal(1),0); |
---|
1031 | module M=T1[3]; |
---|
1032 | intvec v=2,2,2,2,2,2; |
---|
1033 | attrib(M,"isHomog",v); |
---|
1034 | dimH(0,M,2); |
---|
1035 | dimH(1,M,0); |
---|
1036 | dimH(2,M,1); |
---|
1037 | dimH(3,M,-5); |
---|
1038 | } |
---|
1039 | |
---|
1040 | |
---|
1041 | /////////////////////////////////////////////////////////////////////////////// |
---|
1042 | |
---|
1043 | proc sheafCoh(module M,int l,int h,list #) |
---|
1044 | "USAGE: sheafCoh(M,l,h); M module, l,h int |
---|
1045 | ASSUME: @code{M} is graded, and it comes assigned with an admissible degree |
---|
1046 | vector as an attribute, @code{h>=l}. The basering @code{S} has |
---|
1047 | @code{n+1} variables. |
---|
1048 | RETURN: intmat, cohomology of twists of the coherent sheaf F on P^n |
---|
1049 | associated to coker(M). The range of twists is determined by @code{l}, |
---|
1050 | @code{h}. |
---|
1051 | DISPLAY: The intmat is displayed in a diagram of the following form: @* |
---|
1052 | @format |
---|
1053 | l l+1 h |
---|
1054 | ---------------------------------------------------------- |
---|
1055 | n: h^n(F(l)) h^n(F(l+1)) ...... h^n(F(h)) |
---|
1056 | ............................................... |
---|
1057 | 1: h^1(F(l)) h^1(F(l+1)) ...... h^1(F(h)) |
---|
1058 | 0: h^0(F(l)) h^0(F(l+1)) ...... h^0(F(h)) |
---|
1059 | ---------------------------------------------------------- |
---|
1060 | chi: chi(F(l)) chi(F(l+1)) ...... chi(F(h)) |
---|
1061 | @end format |
---|
1062 | A @code{'-'} in the diagram refers to a zero entry. |
---|
1063 | NOTE: The procedure is based on local duality as described in [Eisenbud: |
---|
1064 | Computing cohomology. In Vasconcelos: Computational methods in |
---|
1065 | commutative algebra and algebraic geometry. Springer (1998)].@* |
---|
1066 | By default, the procedure uses @code{mres} to compute the Ext |
---|
1067 | modules. If called with the additional parameter @code{\"sres\"}, |
---|
1068 | the @code{sres} command is used instead. |
---|
1069 | SEE ALSO: dimH, sheafCohBGG |
---|
1070 | EXAMPLE: example sheafCoh; shows an example |
---|
1071 | " |
---|
1072 | { |
---|
1073 | int use_sres; |
---|
1074 | if( typeof(attrib(M,"isHomog"))!="intvec" ) |
---|
1075 | { |
---|
1076 | if (size(M)==0) { attrib(M,"isHomog",0); } |
---|
1077 | else { ERROR("No admissible degree vector assigned"); } |
---|
1078 | } |
---|
1079 | if (size(#)>0) |
---|
1080 | { |
---|
1081 | if (#[1]=="sres") { use_sres=1; } |
---|
1082 | } |
---|
1083 | int i,j; |
---|
1084 | module N,N0,N1; |
---|
1085 | int n=nvars(basering)-1; |
---|
1086 | intvec v=0..n+1; |
---|
1087 | int col=h-l+1; |
---|
1088 | intmat newBetti[n+1][col]; |
---|
1089 | if (use_sres) { list L=Ext_R(v,M,1,"sres")[2]; } |
---|
1090 | else { list L=Ext_R(v,M,1)[2]; } |
---|
1091 | for (i=l; i<=h; i++) |
---|
1092 | { |
---|
1093 | N0=L[n+2]; |
---|
1094 | N1=L[n+1]; |
---|
1095 | newBetti[n+1,i-l+1]=dimGradedPart(M,i) - dimGradedPart(N0,-n-1-i) |
---|
1096 | - dimGradedPart(N0,-n-1-i); |
---|
1097 | } |
---|
1098 | for (j=1; j<=n; j++) |
---|
1099 | { |
---|
1100 | N=L[j]; |
---|
1101 | attrib(N,"isSB",1); |
---|
1102 | if (dim(N)>=0) { |
---|
1103 | for (i=l; i<=h; i++) |
---|
1104 | { |
---|
1105 | newBetti[j,i-l+1]=dimGradedPart(N,-n-1-i); |
---|
1106 | } |
---|
1107 | } |
---|
1108 | } |
---|
1109 | displayCohom(newBetti,l,h,n); |
---|
1110 | return(newBetti); |
---|
1111 | } |
---|
1112 | example |
---|
1113 | {"EXAMPLE:"; |
---|
1114 | echo = 2; |
---|
1115 | // |
---|
1116 | // cohomology of structure sheaf on P^4: |
---|
1117 | //------------------------------------------- |
---|
1118 | ring r=0,x(1..5),dp; |
---|
1119 | module M=0; |
---|
1120 | def A=sheafCoh(0,-7,2); |
---|
1121 | // |
---|
1122 | // cohomology of cotangential bundle on P^3: |
---|
1123 | //------------------------------------------- |
---|
1124 | ring R=0,(x,y,z,u),dp; |
---|
1125 | resolution T1=mres(maxideal(1),0); |
---|
1126 | module M=T1[3]; |
---|
1127 | intvec v=2,2,2,2,2,2; |
---|
1128 | attrib(M,"isHomog",v); |
---|
1129 | def B=sheafCoh(M,-6,2); |
---|
1130 | } |
---|
1131 | |
---|
1132 | /////////////////////////////////////////////////////////////////////////////// |
---|
1133 | proc displayCohom (intmat data, int l, int h, int n) |
---|
1134 | "USAGE: displayCohom(data,l,h,n); data intmat, l,h,n int |
---|
1135 | ASSUME: @code{h>=l}, @code{data} is the return value of |
---|
1136 | @code{sheafCoh(M,l,h)} or of @code{sheafCohBGG(M,l,h)}, and the |
---|
1137 | basering has @code{n+1} variables. |
---|
1138 | RETURN: none |
---|
1139 | NOTE: The intmat is displayed in a diagram of the following form: @* |
---|
1140 | @format |
---|
1141 | l l+1 h |
---|
1142 | ---------------------------------------------------------- |
---|
1143 | n: h^n(F(l)) h^n(F(l+1)) ...... h^n(F(h)) |
---|
1144 | ............................................... |
---|
1145 | 1: h^1(F(l)) h^1(F(l+1)) ...... h^1(F(h)) |
---|
1146 | 0: h^0(F(l)) h^0(F(l+1)) ...... h^0(F(h)) |
---|
1147 | ---------------------------------------------------------- |
---|
1148 | chi: chi(F(l)) chi(F(l+1)) ...... chi(F(h)) |
---|
1149 | @end format |
---|
1150 | where @code{F} refers to the associated sheaf of @code{M} on P^n.@* |
---|
1151 | A @code{'-'} in the diagram refers to a zero entry, a @code{'*'} |
---|
1152 | refers to a negative entry (= dimension not yet determined). |
---|
1153 | " |
---|
1154 | { |
---|
1155 | int i,j,k,dat,maxL; |
---|
1156 | intvec notSumCol; |
---|
1157 | notSumCol[h-l+1]=0; |
---|
1158 | string s; |
---|
1159 | maxL=4; |
---|
1160 | for (i=1;i<=nrows(data);i++) |
---|
1161 | { |
---|
1162 | for (j=1;j<=ncols(data);j++) |
---|
1163 | { |
---|
1164 | if (size(string(data[i,j]))>=maxL-1) |
---|
1165 | { |
---|
1166 | maxL=size(string(data[i,j]))+2; |
---|
1167 | } |
---|
1168 | } |
---|
1169 | } |
---|
1170 | string Row=" "; |
---|
1171 | string Row1="----"; |
---|
1172 | for (i=l; i<=h; i++) { |
---|
1173 | for (j=1; j<=maxL-size(string(i)); j++) |
---|
1174 | { |
---|
1175 | Row=Row+" "; |
---|
1176 | } |
---|
1177 | Row=Row+string(i); |
---|
1178 | for (j=1; j<=maxL; j++) { Row1 = Row1+"-"; } |
---|
1179 | } |
---|
1180 | print(Row); |
---|
1181 | print(Row1); |
---|
1182 | for (j=1; j<=n+1; j++) |
---|
1183 | { |
---|
1184 | s = string(n+1-j); |
---|
1185 | Row = ""; |
---|
1186 | for(k=1; k<4-size(s); k++) { Row = Row+" "; } |
---|
1187 | Row = Row + s+":"; |
---|
1188 | for (i=0; i<=h-l; i++) |
---|
1189 | { |
---|
1190 | dat = data[j,i+1]; |
---|
1191 | if (dat>0) { s = string(dat); } |
---|
1192 | else |
---|
1193 | { |
---|
1194 | if (dat==0) { s="-"; } |
---|
1195 | else { s="*"; notSumCol[i+1]=1; } |
---|
1196 | } |
---|
1197 | for(k=1; k<=maxL-size(s); k++) { Row = Row+" "; } |
---|
1198 | Row = Row + s; |
---|
1199 | } |
---|
1200 | print(Row); |
---|
1201 | } |
---|
1202 | print(Row1); |
---|
1203 | Row="chi:"; |
---|
1204 | for (i=0; i<=h-l; i++) |
---|
1205 | { |
---|
1206 | dat = 0; |
---|
1207 | if (notSumCol[i+1]==0) |
---|
1208 | { |
---|
1209 | for (j=0; j<=n; j++) { dat = dat + (-1)^j * data[n+1-j,i+1]; } |
---|
1210 | s = string(dat); |
---|
1211 | } |
---|
1212 | else { s="*"; } |
---|
1213 | for (k=1; k<=maxL-size(s); k++) { Row = Row+" "; } |
---|
1214 | Row = Row + s; |
---|
1215 | } |
---|
1216 | print(Row); |
---|
1217 | } |
---|
1218 | /////////////////////////////////////////////////////////////////////////////// |
---|
1219 | |
---|
1220 | proc getStructureSheaf(list #) |
---|
1221 | { |
---|
1222 | |
---|
1223 | if( size(#) == 0 ) |
---|
1224 | { |
---|
1225 | module M = 0; |
---|
1226 | intvec v = 0; |
---|
1227 | attrib(M,"isHomog",v); |
---|
1228 | // homog(M); |
---|
1229 | |
---|
1230 | attrib(M, "isCoker", 1); |
---|
1231 | |
---|
1232 | // attrib(M); |
---|
1233 | return(M); |
---|
1234 | } |
---|
1235 | |
---|
1236 | if( typeof(#[1]) == "ideal") |
---|
1237 | { |
---|
1238 | ideal I = #[1]; |
---|
1239 | |
---|
1240 | if( size(#) == 2 ) |
---|
1241 | { |
---|
1242 | if( typeof(#[2]) == "int" ) |
---|
1243 | { |
---|
1244 | if( #[2] != 0 ) |
---|
1245 | { |
---|
1246 | qring @@@@QQ = std(I); |
---|
1247 | |
---|
1248 | module M = getStructureSheaf(); |
---|
1249 | |
---|
1250 | export M; |
---|
1251 | |
---|
1252 | // keepring @@@@QQ; // This is a bad idea... :(? |
---|
1253 | return (@@@@QQ); |
---|
1254 | } |
---|
1255 | } |
---|
1256 | } |
---|
1257 | |
---|
1258 | /* |
---|
1259 | // This seems to be wrong!!! |
---|
1260 | module M = I * gen(1); |
---|
1261 | homog(M); |
---|
1262 | |
---|
1263 | M = modulo(gen(1), module(I * gen(1))); // basering^1 / I |
---|
1264 | |
---|
1265 | homog(M); |
---|
1266 | |
---|
1267 | attrib(M, "isCoker", 1); |
---|
1268 | |
---|
1269 | attrib(M); |
---|
1270 | return(M); |
---|
1271 | */ |
---|
1272 | } |
---|
1273 | |
---|
1274 | ERROR("Wrong argument"); |
---|
1275 | |
---|
1276 | } |
---|
1277 | example |
---|
1278 | {"EXAMPLE:"; |
---|
1279 | echo = 2; int pl = printlevel; |
---|
1280 | printlevel = voice; |
---|
1281 | |
---|
1282 | |
---|
1283 | //////////////////////////////////////////////////////////////////////////////// |
---|
1284 | ring r; |
---|
1285 | module M = getStructureSheaf(); |
---|
1286 | "Basering: "; |
---|
1287 | basering; |
---|
1288 | "Module: ", string(M), ", grading is given by weights: ", attrib(M, "isHomog"); |
---|
1289 | |
---|
1290 | def A=sheafCohBGG2(M,-9,9); |
---|
1291 | print(A); |
---|
1292 | |
---|
1293 | //////////////////////////////////////////////////////////////////////////////// |
---|
1294 | setring r; |
---|
1295 | module M = getStructureSheaf(ideal(var(1)), 0); |
---|
1296 | |
---|
1297 | "Basering: "; |
---|
1298 | basering; |
---|
1299 | "Module: ", string(M), ", grading is given by weights: ", attrib(M, "isHomog"); |
---|
1300 | |
---|
1301 | def A=sheafCohBGG2(M,-9,9); |
---|
1302 | print(A); |
---|
1303 | |
---|
1304 | //////////////////////////////////////////////////////////////////////////////// |
---|
1305 | setring r; |
---|
1306 | def Q = getStructureSheaf(ideal(var(1)), 1); // returns a new ring! |
---|
1307 | setring Q; // M was exported in the new ring! |
---|
1308 | |
---|
1309 | "Basering: "; |
---|
1310 | basering; |
---|
1311 | "Module: ", string(M), ", grading is given by weights: ", attrib(M, "isHomog"); |
---|
1312 | |
---|
1313 | def A=sheafCohBGG2(M,-9,9); |
---|
1314 | print(A); |
---|
1315 | |
---|
1316 | printlevel = pl; |
---|
1317 | } |
---|
1318 | |
---|
1319 | |
---|
1320 | proc getCotangentialBundle() |
---|
1321 | { |
---|
1322 | resolution T1=mres(maxideal(1),3); |
---|
1323 | module M=T1[3]; |
---|
1324 | // attrib(M,"isHomog"); |
---|
1325 | // homog(M); |
---|
1326 | attrib(M, "isCoker", 1); |
---|
1327 | // attrib(M); |
---|
1328 | return (M); |
---|
1329 | } |
---|
1330 | |
---|
1331 | proc getIdealSheafPullback(ideal I, ideal pi) |
---|
1332 | { |
---|
1333 | def save = basering; |
---|
1334 | map P = save, pi; |
---|
1335 | return( P(I) ); |
---|
1336 | } |
---|
1337 | |
---|
1338 | // TODO: set attributes! |
---|
1339 | |
---|
1340 | |
---|
1341 | proc getIdealSheaf(ideal I) |
---|
1342 | { |
---|
1343 | int i = homog(I); |
---|
1344 | resolution FI = mres(I,2); // Syz + grading... |
---|
1345 | module M = FI[2]; |
---|
1346 | attrib(M, "isCoker", 1); |
---|
1347 | // attrib(M); |
---|
1348 | return(M); |
---|
1349 | } |
---|
1350 | |
---|
1351 | |
---|
1352 | |
---|
1353 | |
---|
1354 | |
---|
1355 | /* |
---|
1356 | Examples: |
---|
1357 | --------- |
---|
1358 | LIB "sheafcoh.lib"; |
---|
1359 | |
---|
1360 | ring S = 32003, x(0..4), dp; |
---|
1361 | module MI=maxideal(1); |
---|
1362 | attrib(MI,"isHomog",intvec(-1)); |
---|
1363 | resolution kos = nres(MI,0); |
---|
1364 | print(betti(kos),"betti"); |
---|
1365 | LIB "random.lib"; |
---|
1366 | matrix alpha0 = random(32002,10,3); |
---|
1367 | module pres = module(alpha0)+kos[3]; |
---|
1368 | attrib(pres,"isHomog",intvec(1,1,1,1,1,1,1,1,1,1)); |
---|
1369 | resolution fcokernel = mres(pres,0); |
---|
1370 | print(betti(fcokernel),"betti"); |
---|
1371 | module dir = transpose(pres); |
---|
1372 | attrib(dir,"isHomog",intvec(-1,-1,-1,-2,-2,-2, |
---|
1373 | -2,-2,-2,-2,-2,-2,-2)); |
---|
1374 | resolution fdir = mres(dir,2); |
---|
1375 | print(betti(fdir),"betti"); |
---|
1376 | ideal I = groebner(flatten(fdir[2])); |
---|
1377 | resolution FI = mres(I,0); |
---|
1378 | print(betti(FI),"betti"); |
---|
1379 | module F=FI[2]; |
---|
1380 | int t=timer; |
---|
1381 | def A1=sheafCoh(F,-8,8); |
---|
1382 | timer-t; |
---|
1383 | t=timer; |
---|
1384 | def A2=sheafCohBGG(F,-8,8); |
---|
1385 | timer-t; |
---|
1386 | |
---|
1387 | LIB "sheafcoh.lib"; |
---|
1388 | LIB "random.lib"; |
---|
1389 | ring S = 32003, x(0..4), dp; |
---|
1390 | resolution kos = nres(maxideal(1),0); |
---|
1391 | betti(kos); |
---|
1392 | matrix kos5 = kos[5]; |
---|
1393 | matrix tphi = transpose(dsum(kos5,kos5)); |
---|
1394 | matrix kos3 = kos[3]; |
---|
1395 | matrix psi = dsum(kos3,kos3); |
---|
1396 | matrix beta1 = random(32002,20,2); |
---|
1397 | matrix tbeta1tilde = transpose(psi*beta1); |
---|
1398 | matrix tbeta0 = lift(tphi,tbeta1tilde); |
---|
1399 | matrix kos4 = kos[4]; |
---|
1400 | matrix tkos4pluskos4 = transpose(dsum(kos4,kos4)); |
---|
1401 | matrix tgammamin1 = random(32002,20,1); |
---|
1402 | matrix tgamma0 = tkos4pluskos4*tgammamin1; |
---|
1403 | matrix talpha0 = concat(tbeta0,tgamma0); |
---|
1404 | matrix zero[20][1]; |
---|
1405 | matrix tpsi = transpose(psi); |
---|
1406 | matrix tpresg = concat(tpsi,zero); |
---|
1407 | matrix pres = module(transpose(talpha0)) |
---|
1408 | + module(transpose(tpresg)); |
---|
1409 | module dir = transpose(pres); |
---|
1410 | dir = prune(dir); |
---|
1411 | homog(dir); |
---|
1412 | intvec deg_dir = attrib(dir,"isHomog"); |
---|
1413 | attrib(dir,"isHomog",deg_dir-2); // set degrees |
---|
1414 | resolution fdir = mres(prune(dir),2); |
---|
1415 | print(betti(fdir),"betti"); |
---|
1416 | ideal I = groebner(flatten(fdir[2])); |
---|
1417 | resolution FI = mres(I,0); |
---|
1418 | |
---|
1419 | module F=FI[2]; |
---|
1420 | def A1=sheafCoh(F,-5,7); |
---|
1421 | def A2=sheafCohBGG(F,-5,7); |
---|
1422 | |
---|
1423 | */ |
---|