source: git/Singular/LIB/sheafcoh.lib @ a890be

spielwiese
Last change on this file since a890be was a890be, checked in by Oleksandr Motsak <motsak@…>, 14 years ago
minor fix of behaviour git-svn-id: file:///usr/local/Singular/svn/trunk@12764 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 37.6 KB
Line 
1///////////////////////////////////////////////////////////////////////////////
2version="$Id$";
3category="Commutative Algebra";
4info="
5LIBRARY:  sheafcoh.lib   Procedures for Computing Sheaf Cohomology
6AUTHORS:  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
11PROCEDURES:
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
21AUXILIARY PROCEDURES:
22 displayCohom(B,l,h,n);  display intmat as Betti diagram (with zero rows)
23
24KEYWORDS: sheaf cohomology
25";
26
27///////////////////////////////////////////////////////////////////////////////
28LIB "matrix.lib";
29LIB "nctools.lib";
30LIB "homolog.lib";
31
32
33///////////////////////////////////////////////////////////////////////////////
34static proc jacobM(matrix M)
35{
36   int n=nvars(basering);
37   matrix B=transpose(diff(M,var(1)));
38   int i;
39   for(i=2;i<=n;i++)
40   {
41     B=concat(B,transpose(diff(M,var(i))));
42   }
43   return(transpose(B));
44}
45
46
47///////////////////////////////////////////////////////////////////////////////
48/**
49  let M = { w_1, ..., w_k }, k = size(M) == ncols(M), n = nvars(currRing).
50  assuming that nrows(M) <= m*n;
51  computes transpose(M) * transpose( var(1) I_m | ... | var(n) I_m ) :== transpose(module{f_1, ... f_k}),
52  where f_i = \sum_{j=1}^{m} (w_i, v_j) gen(j),  (w_i, v_j) is a `scalar` multiplication.
53  that is, if w_i = (a^1_1, ... a^1_m) | (a^2_1, ..., a^2_m) | ... | (a^n_1, ..., a^n_m) then
54
55  (a^1_1, ... a^1_m) | (a^2_1, ..., a^2_m) | ... | (a^n_1, ..., a^n_m)
56*  var_1  ... var_1  |  var_2  ...  var_2  | ... |  var_n  ...  var(n)
57*  gen_1  ... gen_m  |  gen_1  ...  gen_m  | ... |  gen_1  ...  gen_m
58+ =>
59  f_i =
60
61   a^1_1 * var(1) * gen(1) + ... + a^1_m * var(1) * gen(m) +
62   a^2_1 * var(2) * gen(1) + ... + a^2_m * var(2) * gen(m) +
63                             ...
64   a^n_1 * var(n) * gen(1) + ... + a^n_m * var(n) * gen(m);
65
66   NOTE: for every f_i we run only ONCE along w_i saving partial sums into a temporary array of polys of size m
67
68*/
69static proc TensorModuleMult(int m, module M)
70{
71  return( system("tensorModuleMult", m, M) ); // trick!
72
73  int n = nvars(basering);
74  int k = ncols(M);
75
76  int g, cc, vv;
77
78  poly h;
79
80  module Temp; // = {f_1, ..., f_k }
81
82  intvec exp;
83  vector pTempSum, w;
84
85  for( int i = k; i > 0; i-- ) // for every w \in M
86  {
87    pTempSum[m] = 0;
88    w = M[i];
89
90    while(w != 0) // for each term of w...
91    {
92      exp = leadexp(w);
93      g = exp[n+1]; // module component!
94      h = w[g];
95
96      w = w - h * gen(g);
97
98      cc = g % m;
99
100      if( cc == 0)
101      {
102        cc = m;
103      }
104
105      vv = 1 + (g - cc) / m;
106
107      pTempSum = pTempSum + h * var(vv) * gen(cc);
108    }
109
110    Temp[i] = pTempSum;
111  }
112
113  Temp = transpose(Temp);
114
115  return(Temp);
116}
117
118///////////////////////////////////////////////////////////////////////////////
119static proc max(int i,int j)
120{
121  if(i>j){return(i);}
122  return(j);
123}
124
125///////////////////////////////////////////////////////////////////////////////
126static proc min(int i,int j)
127{
128  if(i>j){return(j);}
129  return(i);
130}
131
132
133///////////////////////////////////////////////////////////////////////////////
134proc truncate(module phi, int d)
135"USAGE:   truncate(M,d);  M module, d int
136ASSUME:  @code{M} is graded, and it comes assigned with an admissible degree
137         vector as an attribute
138RETURN:  module
139NOTE:    Output is a presentation matrix for the truncation of coker(M)
140         at degree d.
141EXAMPLE: example truncate; shows an example
142KEYWORDS: truncated module
143"
144{
145  if ( typeof(attrib(phi,"isHomog"))=="string" ) {
146    if (size(phi)==0) {
147      // assign weights 0 to generators of R^n (n=nrows(phi))
148      intvec v;
149      v[nrows(phi)]=0;
150      attrib(phi,"isHomog",v);
151    }
152    else {
153      ERROR("No admissible degree vector assigned");
154    }
155  }
156  else {
157   intvec v=attrib(phi,"isHomog");
158  }
159  int i,m,dummy;
160  int s = nrows(phi);
161  module L; // TOO BIG!!!
162  for (i=1; i<=s; i++) {
163    if (d>v[i]) {
164      L = L+maxideal(d-v[i])*gen(i);
165    }
166    else {
167      L = L+gen(i);
168    }
169  }
170  L = modulo(L,phi);
171  L = minbase(prune(L));
172  if (size(L)==0) {return(L);}
173
174  // it only remains to set the degrees for L:
175  // ------------------------------------------
176  m = v[1];
177  for(i=2; i<=size(v); i++) {  if(v[i]<m) { m = v[i]; } }
178  dummy = homog(L);
179  intvec vv = attrib(L,"isHomog");
180  if (d>m) { vv = vv+d; }
181  else     { vv = vv+m; }
182  attrib(L,"isHomog",vv);
183  return(L);
184}
185example
186{"EXAMPLE:";
187   echo = 2;
188   ring R=0,(x,y,z),dp;
189   module M=maxideal(3);
190   homog(M);
191   // compute presentation matrix for truncated module (R/<x,y,z>^3)_(>=2)
192   module M2=truncate(M,2);
193   print(M2);
194   dimGradedPart(M2,1);
195   dimGradedPart(M2,2);
196   // this should coincide with:
197   dimGradedPart(M,2);
198   // shift grading by 1:
199   intvec v=1;
200   attrib(M,"isHomog",v);
201   M2=truncate(M,2);
202   print(M2);
203   dimGradedPart(M2,3);
204}
205
206///////////////////////////////////////////////////////////////////////////////
207
208proc truncateFast(module M, int d)
209"USAGE:   truncateFast(M,d);  M module, d int
210ASSUME:  @code{M} is graded, and it comes assigned with an admissible degree
211         vector as an attribute 'isHomog'
212RETURN:  module
213NOTE:    Output is a presentation matrix for the truncation of coker(M)
214         at d.
215         Fast + experimental version. M shoud be a SB!
216DISPLAY: If @code{printlevel}>=1, step-by step timings will be printed.
217         If @code{printlevel}>=2 we add progress debug messages
218         if @code{printlevel}>=3, even all intermediate results...
219EXAMPLE: example truncateFast; shows an example
220KEYWORDS: truncated module
221"
222{
223//  int PL = printlevel + 1;
224  int PL = printlevel - voice + 2;
225
226  dbprint(PL-1, "// truncateFast(M: "+ string(nrows(M)) + " x " + string(ncols(M)) +", " + string(d) + "):");
227  dbprint(PL-2, M);
228
229  intvec save = option(get);
230  if( PL >= 2 )
231  {
232    option(prot);
233    option(mem);
234  }
235
236  int tTruncateBegin=timer;
237
238  if (attrib(M,"isSB")!=1)
239  {
240    ERROR("M must be a standard basis!");
241  };
242
243  dbprint(PL-1, "// M is a SB! ");
244
245  if ( typeof(attrib(M,"isHomog"))=="string" ) {
246    if (size(M)==0) {
247      // assign weights 0 to generators of R^n (n=nrows(M))
248      intvec v;
249      v[nrows(M)]=0;
250      attrib(M,"isHomog",v);
251    }
252    else {
253      ERROR("No admissible degree vector assigned");
254    }
255  }
256  else {
257   intvec v=attrib(M,"isHomog");
258  }
259
260  dbprint(PL-1, "// weighting(M): ["+ string(v) + "]");
261
262  int i,m,dummy;
263  int s = nrows(M);
264
265   int tKBaseBegin = timer;
266    module L = kbase(M, d); // TODO: check whether this is always correct!?!
267
268
269    dbprint(PL-1, "// L = kbase(M,d): "+string(nrows(L)) + " x " + string(ncols(L)) +"");
270    dbprint(PL-2, L);
271    dbprint(PL-1, "// weighting(L): ["+ string(attrib(L, "isHomog")) + "]");
272
273   int tModuloBegin = timer;
274    L = modulo(L,M);
275
276    dbprint(PL-1, "// L = modulo(L,M): "+string(nrows(L)) + " x " + string(ncols(L)) +"");
277    dbprint(PL-2, L);
278    dbprint(PL-1, "// weighting(L): ["+ string(attrib(L, "isHomog")) + "]");
279
280   int tPruneBegin = timer;
281    L = prune(L);
282
283    dbprint(PL-1, "// L = prune(L): "+string(nrows(L)) + " x " + string(ncols(L)) +"");
284    dbprint(PL-2, L);
285    dbprint(PL-1, "// weighting(L): ["+ string(attrib(L, "isHomog")) + "]");
286
287   int tPruneEnd = timer;
288    L = minbase(L);
289   int tMinBaseEnd = timer;
290
291    dbprint(PL-1, "// L = minbase(L): "+string(nrows(L)) + " x " + string(ncols(L)) +"");
292    dbprint(PL-2, L);
293    dbprint(PL-1, "// weighting(L): ["+ string(attrib(L, "isHomog")) + "]");
294
295
296
297
298  if (size(L)!=0)
299  {
300
301    // it only remains to set the degrees for L:
302    // ------------------------------------------
303    m = v[1];
304    for(i=2; i<=size(v); i++) {  if(v[i]<m) { m = v[i]; } }
305    dummy = homog(L);
306    intvec vv = attrib(L,"isHomog");
307    if (d>m) { vv = vv+d; }
308    else     { vv = vv+m; }
309    attrib(L,"isHomog",vv);
310  }
311
312  int tTruncateEnd=timer;
313
314  dbprint(PL-1, "// corrected weighting(L): ["+ string(attrib(L, "isHomog")) + "]");
315
316
317  if(PL > 0)
318  {
319  "
320  -------------- TIMINGS --------------
321  Trunc        Time: ", tTruncateEnd - tTruncateBegin, "
322    :: Before .Time: ", tKBaseBegin - tTruncateBegin, "
323    :: kBase   Time: ", tModuloBegin - tKBaseBegin, "
324    :: Modulo  Time: ", tPruneBegin - tModuloBegin, "
325    :: Prune   Time: ", tPruneEnd - tPruneBegin, "
326    :: Minbase Time: ", tMinBaseEnd - tPruneEnd, "
327    :: After  .Time: ", tTruncateEnd - tMinBaseEnd;
328  }
329
330  option(set, save);
331
332  return(L);
333}
334example
335{"EXAMPLE:";
336   echo = 2;
337   ring R=0,(x,y,z,u,v),dp;
338   module M=maxideal(3);
339   homog(M);
340   // compute presentation matrix for truncated module (R/<x,y,z,u>^3)_(>=2)
341   int t=timer;
342   module M2t=truncate(M,2);
343   t = timer - t;
344   "// Simple truncate: ", t;
345   t=timer;
346   module M2=truncateFast(std(M),2);
347   t = timer - t;
348   "// Fast truncate: ", t;
349   print(M2);
350   "// Check: M2t == M2?: ", size(NF(M2, std(M2t))) + size(NF(M2t, std(M2)));
351
352   dimGradedPart(M2,1);
353   dimGradedPart(M2,2);
354   // this should coincide with:
355   dimGradedPart(M,2);
356   // shift grading by 1:
357   intvec v=1;
358   attrib(M,"isHomog",v);
359   t=timer;
360   M2t=truncate(M,2);
361   t = timer - t;
362   "// Simple truncate: ", t;
363   t=timer;
364   M2=truncateFast(std(M),2);
365   t = timer - t;
366   "// Fast truncate: ", t;
367   print(M2);
368   "// Check: M2t == M2?: ", size(NF(M2, std(M2t))) + size(NF(M2t, std(M2))); //?
369   dimGradedPart(M2,3);
370}
371
372///////////////////////////////////////////////////////////////////////////////
373
374
375
376
377
378proc dimGradedPart(module phi, int d)
379"USAGE:   dimGradedPart(M,d);  M module, d int
380ASSUME:  @code{M} is graded, and it comes assigned with an admissible degree
381         vector as an attribute
382RETURN:  int
383NOTE:    Output is the vector space dimension of the graded part of degree d
384         of coker(M).
385EXAMPLE: example dimGradedPart; shows an example
386KEYWORDS: graded module, graded piece
387"
388{
389  if ( typeof(attrib(phi,"isHomog"))=="string" ) {
390    if (size(phi)==0) {
391      // assign weights 0 to generators of R^n (n=nrows(phi))
392      intvec v;
393      v[nrows(phi)]=0;
394    }
395    else { ERROR("No admissible degree vector assigned"); }
396  }
397  else {
398    intvec v=attrib(phi,"isHomog");
399  }
400  int s = nrows(phi);
401  int i,m,dummy;
402  module L,LL;
403  for (i=1; i<=s; i++) {
404    if (d>v[i]) {
405      L = L+maxideal(d-v[i])*gen(i);
406      LL = LL+maxideal(d+1-v[i])*gen(i);
407    }
408    else {
409      L = L+gen(i);
410      if (d==v[i]) {
411        LL = LL+maxideal(1)*gen(i);
412      }
413      else {
414        LL = LL+gen(i);
415      }
416    }
417  }
418  LL=LL,phi;
419  L = modulo(L,LL);
420  L = std(prune(L));
421  if (size(L)==0) {return(0);}
422  return(vdim(L));
423}
424example
425{"EXAMPLE:";
426   echo = 2;
427   ring R=0,(x,y,z),dp;
428   module M=maxideal(3);
429   // assign compatible weight vector (here: 0)
430   homog(M);
431   // compute dimension of graded pieces of R/<x,y,z>^3 :
432   dimGradedPart(M,0);
433   dimGradedPart(M,1);
434   dimGradedPart(M,2);
435   dimGradedPart(M,3);
436   // shift grading:
437   attrib(M,"isHomog",intvec(2));
438   dimGradedPart(M,2);
439}
440
441///////////////////////////////////////////////////////////////////////////////
442
443proc CM_regularity (module M)
444"USAGE:   CM_regularity(M);    M module
445ASSUME:   @code{M} is graded, and it comes assigned with an admissible degree
446         vector as an attribute
447RETURN:  integer, the Castelnuovo-Mumford regularity of coker(M)
448NOTE:    procedure calls mres
449EXAMPLE: example CM_regularity; shows an example
450KEYWORDS: Castelnuovo-Mumford regularity
451"
452{
453  if ( typeof(attrib(M,"isHomog"))=="string" ) {
454    if (size(M)==0) {
455      // assign weights 0 to generators of R^n (n=nrows(M))
456      intvec v;
457      v[nrows(M)]=0;
458      attrib(M,"isHomog",v);
459    }
460    else {
461      ERROR("No admissible degree vector assigned");
462    }
463  }
464
465  if( attrib(CM_regularity,"Algorithm") == "minres_res" )
466  {
467    def L = minres( res(M,0) ); // let's try it out!
468  } else
469  {
470    def L = mres(M,0);
471  }
472
473  intmat BeL = betti(L);
474  int r = nrows(module(matrix(BeL)));  // last non-zero row
475  if (typeof(attrib(BeL,"rowShift"))!="string") {
476    int shift = attrib(BeL,"rowShift");
477  }
478  return(r+shift-1);
479}
480example
481{"EXAMPLE:";
482   echo = 2;
483   ring R=0,(x,y,z,u),dp;
484   resolution T1=mres(maxideal(1),0);
485   module M=T1[3];
486   intvec v=2,2,2,2,2,2;
487   attrib(M,"isHomog",v);
488   CM_regularity(M);
489}
490
491///////////////////////////////////////////////////////////////////////////////
492proc sheafCohBGG(module M,int l,int h)
493"USAGE:   sheafCohBGG(M,l,h);    M module, l,h int
494ASSUME:  @code{M} is graded, and it comes assigned with an admissible degree
495         vector as an attribute, @code{h>=l}, and the basering has @code{n+1}
496         variables.
497RETURN:  intmat, cohomology of twists of the coherent sheaf F on P^n
498         associated to coker(M). The range of twists is determined by @code{l},
499         @code{h}.
500DISPLAY: The intmat is displayed in a diagram of the following form: @*
501  @format
502                l            l+1                      h
503  ----------------------------------------------------------
504      n:     h^n(F(l))    h^n(F(l+1))   ......    h^n(F(h))
505           ...............................................
506      1:     h^1(F(l))    h^1(F(l+1))   ......    h^1(F(h))
507      0:     h^0(F(l))    h^0(F(l+1))   ......    h^0(F(h))
508  ----------------------------------------------------------
509    chi:     chi(F(l))    chi(F(l+1))   ......    chi(F(h))
510  @end format
511         A @code{'-'} in the diagram refers to a zero entry; a @code{'*'}
512         refers to a negative entry (= dimension not yet determined).
513         refers to a not computed dimension. @*
514NOTE:    This procedure is based on the Bernstein-Gel'fand-Gel'fand
515         correspondence and on Tate resolution ( see [Eisenbud, Floystad,
516         Schreyer: Sheaf cohomology and free resolutions over exterior
517         algebras, Trans AMS 355 (2003)] ).@*
518         @code{sheafCohBGG(M,l,h)} does not compute all values in the above
519         table. To determine all values of @code{h^i(F(d))}, @code{d=l..h},
520         use @code{sheafCohBGG(M,l-n,h+n)}.
521SEE ALSO: sheafCoh, dimH
522EXAMPLE: example sheafCohBGG; shows an example
523"
524{
525  int i,j,k,row,col;
526  if( typeof(attrib(M,"isHomog"))!="intvec" )
527  {
528     if (size(M)==0) { attrib(M,"isHomog",0); }
529     else { ERROR("No admissible degree vector assigned"); }
530  }
531  int n=nvars(basering)-1;
532  int ell=l+n;
533  def R=basering;
534  int reg = CM_regularity(M);
535  int bound=max(reg+1,h-1);
536  module MT=truncate(M,bound);
537  int m=nrows(MT);
538  MT=transpose(jacobM(MT));
539  MT=syz(MT);
540  matrix ML[n+1][1]=maxideal(1);
541  matrix S=transpose(outer(ML,unitmat(m)));
542  matrix SS=transpose(S*MT);
543  //--- to the exterior algebra
544  def AR = Exterior();
545  setring AR;
546  option(redSB);
547  option(redTail);
548  module EM=imap(R,SS);
549  intvec w;
550  //--- here we are with our matrix
551  int bound1=max(1,bound-ell+1);
552  for (i=1; i<=nrows(EM); i++)
553  {
554     w[i]=-bound-1;
555  }
556  attrib(EM,"isHomog",w);
557  resolution RE=mres(EM,bound1);
558  intmat Betti=betti(RE);
559  k=ncols(Betti);
560  row=nrows(Betti);
561  int shift=attrib(Betti,"rowShift")+(k+ell-1);
562  intmat newBetti[n+1][h-l+1];
563  for (j=1; j<=row; j++)
564  {
565    for (i=l; i<=h; i++)
566    {
567      if ((k+1-j-i+ell-shift>0) and (j+i-ell+shift>=1))
568      {
569        newBetti[n+2-shift-j,i-l+1]=Betti[j,k+1-j-i+ell-shift];
570      }
571      else { newBetti[n+2-shift-j,i-l+1]=-1; }
572    }
573  }
574  for (j=2; j<=n+1; j++)
575  {
576    for (i=1; i<j; i++)
577    {
578      newBetti[j,i]=-1;
579    }
580  }
581  int d=k-h+ell-1;
582  for (j=1; j<=n; j++)
583  {
584    for (i=h-l+1; i>=k+j; i--)
585    {
586      newBetti[j,i]=-1;
587    }
588  }
589  displayCohom(newBetti,l,h,n);
590  setring R;
591  return(newBetti);
592  option(noredSB);
593  option(noredTail);
594}
595example
596{"EXAMPLE:";
597   echo = 2;
598   // cohomology of structure sheaf on P^4:
599   //-------------------------------------------
600   ring r=0,x(1..5),dp;
601   module M=0;
602   def A=sheafCohBGG(M,-9,4);
603   // cohomology of cotangential bundle on P^3:
604   //-------------------------------------------
605   ring R=0,(x,y,z,u),dp;
606   resolution T1=mres(maxideal(1),0);
607   module M=T1[3];
608   intvec v=2,2,2,2,2,2;
609   attrib(M,"isHomog",v);
610   def B=sheafCohBGG(M,-8,4);
611}
612
613
614///////////////////////////////////////////////////////////////////////////////
615static proc showResult( def R, int l, int h )
616{
617  int PL = 1; // printlevel - voice + 2;
618// int PL = printlevel + 1;
619
620  intmat Betti;
621  if(typeof(R)=="resolution")
622  {
623    Betti = betti(R);
624  } else
625  {
626    if(typeof(R)!="intmat")
627    {
628      ERROR("Wrong input!!!");
629    };
630
631    Betti = R;
632  }
633
634  int n=nvars(basering)-1;
635  int ell = l + n;
636
637  int k     = ncols(Betti);
638  int row   = nrows(Betti);
639  int shift = attrib(Betti,"rowShift") + (k + ell - 1);
640
641  intmat newBetti[ n + 1 ][ h - l + 1 ];
642
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  int iWTH = h-l+1;
659  for (j=2; j<=n+1; j++) {
660    for (i=1; i<min(j, iWTH); i++) {
661      newBetti[j,i]=-1;
662    }
663  }
664  int d = k - h + ell - 1;
665  for (j=1; j<=n; j++) {
666    for (i=iWTH; i>=k+j; i--) {
667      newBetti[j,i]=-1;
668    }
669  }
670
671  if( PL > 0 )
672  {
673    "Cohomology table:";
674    displayCohom(newBetti, l, h, n);
675  }
676
677  return(newBetti);
678}
679///////////////////////////////////////////////////////////////////////////////
680
681
682
683///////////////////////////////////////////////////////////////////////////////
684proc sheafCohBGG2(module M,int l,int h)
685"USAGE:   sheafCohBGG2(M,l,h);    M module, l,h int
686ASSUME:  @code{M} is graded, and it comes assigned with an admissible degree
687         vector as an attribute, @code{h>=l}, and the basering has @code{n+1}
688         variables.
689RETURN:  intmat, cohomology of twists of the coherent sheaf F on P^n
690         associated to coker(M). The range of twists is determined by @code{l},
691         @code{h}.
692DISPLAY: The intmat is displayed in a diagram of the following form: @*
693  @format
694                l            l+1                      h
695  ----------------------------------------------------------
696      n:     h^n(F(l))    h^n(F(l+1))   ......    h^n(F(h))
697           ...............................................
698      1:     h^1(F(l))    h^1(F(l+1))   ......    h^1(F(h))
699      0:     h^0(F(l))    h^0(F(l+1))   ......    h^0(F(h))
700  ----------------------------------------------------------
701    chi:     chi(F(l))    chi(F(l+1))   ......    chi(F(h))
702  @end format
703         A @code{'-'} in the diagram refers to a zero entry; a @code{'*'}
704         refers to a negative entry (= dimension not yet determined).
705         refers to a not computed dimension. @*
706         If @code{printlevel}>=1, step-by step timings will be printed.
707         If @code{printlevel}>=2 we add progress debug messages
708         if @code{printlevel}>=3, even all intermediate results...
709NOTE:    This procedure is based on the Bernstein-Gel'fand-Gel'fand
710         correspondence and on Tate resolution ( see [Eisenbud, Floystad,
711         Schreyer: Sheaf cohomology and free resolutions over exterior
712         algebras, Trans AMS 355 (2003)] ).@*
713         @code{sheafCohBGG(M,l,h)} does not compute all values in the above
714         table. To determine all values of @code{h^i(F(d))}, @code{d=l..h},
715         use @code{sheafCohBGG(M,l-n,h+n)}.
716         Experimental version. Should require less memory.
717SEE ALSO: sheafCohBGG
718EXAMPLE: example sheafCohBGG2; shows an example
719"
720{
721  int PL = printlevel - voice + 2;
722//  int PL = printlevel;
723
724  dbprint(PL-1, "// sheafCohBGG2(M: "+ string(nrows(M)) + " x " + string(ncols(M)) +", " + string(l) + ", " + string(h) + "):");
725  dbprint(PL-2, M);
726
727  intvec save = option(get);
728
729  if( PL >= 2 )
730  {
731    option(prot);
732    option(mem);
733  }
734
735  def isCoker = attrib(M, "isCoker");
736  if( typeof(isCoker) == "int" )
737  {
738    if( isCoker > 0 )
739    {
740      dbprint(PL-1, "We are going to assume that M is given by coker matrix (that is, M is not a submodule presentation!)");
741    }
742  }
743
744  int i,j,k,row,col;
745
746  if( typeof(attrib(M,"isHomog"))!="intvec" )
747  {
748     if (size(M)==0) { attrib(M,"isHomog",0); }
749     else { ERROR("No admissible degree vector assigned"); }
750  }
751
752  dbprint(PL-1, "// weighting(M): ["+ string(attrib(M, "isHomog")) + "]");
753
754  option(redSB); option(redTail);
755
756  def R=basering;
757
758  int n = nvars(R) - 1;
759  int ell = l + n;
760
761
762/////////////////////////////////////////////////////////////////////////////
763// computations
764
765  int tBegin=timer;
766  int reg   = CM_regularity(M);
767  int tCMEnd = timer;
768
769  dbprint(PL-1, "// CM_reg(M): "+ string(reg));
770
771  int bound = max(reg + 1, h - 1);
772
773  dbprint(PL-1, "// bound: "+ string(bound));
774
775  ///////////////////////////////////////////////////////////////
776  int tSTDBegin=timer;
777  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!)
778  int tSTDEnd = timer;
779
780  dbprint(PL-1, "// M = std(M: "+string(nrows(M)) + " x " + string(ncols(M))  + ")");
781  dbprint(PL-2, M);
782  dbprint(PL-1, "// weighting(M): ["+ string(attrib(M, "isHomog")) + "]");
783
784
785  printlevel = printlevel + 1;
786  int tTruncateBegin=timer;
787  module MT = truncateFast(M, bound);
788  int tTruncateEnd=timer;
789  printlevel = printlevel - 1;
790
791  dbprint(PL-1, "// MT = truncateFast(M: "+string(nrows(MT)) + " x " + string(ncols(MT)) +", " + string(bound) + ")");
792  dbprint(PL-2, MT);
793  dbprint(PL-1, "// weighting(MT): ["+ string(attrib(MT, "isHomog")) + "]");
794
795  int m=nrows(MT);
796
797  ///////////////////////////////////////////////////////////////
798  int tTransposeJacobBegin=timer;
799  MT = jacob(MT); // ! :(
800  int tTransposeJacobEnd=timer;
801
802  dbprint(PL-1, "// MT = jacob(MT: "+string(nrows(MT)) + " x " + string(ncols(MT))  + ")");
803  dbprint(PL-2, MT);
804  dbprint(PL-1, "// weighting(MT): ["+ string(attrib(MT, "isHomog")) + "]");
805
806  int tSyzBegin=timer;
807  MT = syz(MT);
808  int tSyzEnd=timer;
809
810  dbprint(PL-1, "// MT = syz(MT: "+string(nrows(MT)) + " x " + string(ncols(MT))  + ")");
811  dbprint(PL-2, MT);
812  dbprint(PL-1, "// weighting(MT): ["+ string(attrib(MT, "isHomog")) + "]");
813
814  int tMatrixOppBegin=timer;
815  module SS = TensorModuleMult(m, MT);
816  int tMatrixOppEnd=timer;
817
818  dbprint(PL-1, "// SS = TensorModuleMult("+ string(m)+ ", MT: "+string(nrows(MT)) + " x " + string(ncols(MT))  + ")");
819  dbprint(PL-2, SS);
820  dbprint(PL-1, "// weighting(SS): ["+ string(attrib(SS, "isHomog")) + "]");
821
822  //--- to the exterior algebra
823  def AR = Exterior(); setring AR;
824
825  dbprint(PL-1, "// Test: var(1) * var(1): "+ string(var(1) * var(1)));
826
827  int maxbound = max(1, bound - ell + 1);
828//  int maxbound = max(1, bound - l + 1); // As In M2!!!
829
830  dbprint(PL-1, "// maxbound: "+ string(maxbound));
831
832  //--- here we are with our matrix
833  module EM=imap(R,SS);
834  intvec w;
835  for (i=1; i<=nrows(EM); i++)
836  {
837     w[i]=-bound-1;
838  }
839
840  attrib(EM,"isHomog",w);
841
842  ///////////////////////////////////////////////////////////////
843
844  dbprint(PL-1, "// EM: "+string(nrows(EM)) + " x " + string(ncols(EM))  + ")");
845  dbprint(PL-2, EM);
846  dbprint(PL-1, "// weighting(EM): ["+ string(attrib(EM, "isHomog")) + "]");
847
848  int tResulutionBegin=timer;
849  resolution RE = nres(EM, maxbound);
850  int tMinResBegin=timer;
851  RE = minres(RE);
852  int tBettiBegin=timer;
853  intmat Betti = betti(RE); // betti(RE, 1);?
854  int tResulutionEnd=timer;
855
856  int tEnd = tResulutionEnd;
857
858  if( PL > 0 )
859  {
860    "
861        ----      RESULTS  ----
862        Tate Resolution (Length: ", size(RE), "):
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}
900example
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
975proc dimH(int i,module M,int d)
976"USAGE:   dimH(i,M,d);    M module, i,d int
977ASSUME:  @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.
980RETURN:  int,  vector space dimension of @math{H^i(F(d))} for F the coherent
981         sheaf on P^n associated to coker(M).
982NOTE:    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)].
985SEE ALSO: sheafCoh, sheafCohBGG
986EXAMPLE: 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}
1026example
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
1043proc sheafCoh(module M,int l,int h,list #)
1044"USAGE:   sheafCoh(M,l,h);    M module, l,h int
1045ASSUME:  @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.
1048RETURN:  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}.
1051DISPLAY: 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.
1063NOTE:    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.
1069SEE ALSO: dimH, sheafCohBGG
1070EXAMPLE: 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}
1112example
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///////////////////////////////////////////////////////////////////////////////
1133proc displayCohom (intmat data, int l, int h, int n)
1134"USAGE:   displayCohom(data,l,h,n);  data intmat, l,h,n int
1135ASSUME:  @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.
1138RETURN:  none
1139NOTE:    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
1220proc 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}
1277example
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
1320proc 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
1331proc 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
1341proc getIdealSheaf(ideal I)
1342{
1343  resolution FI = mres(I,2); // Syz + grading...
1344  module M = FI[2];
1345  attrib(M, "isCoker", 1);
1346  attrib(M);
1347  return(M);
1348};
1349
1350
1351
1352
1353
1354/*
1355Examples:
1356---------
1357 LIB "sheafcoh.lib";
1358
1359 ring S = 32003, x(0..4), dp;
1360 module MI=maxideal(1);
1361 attrib(MI,"isHomog",intvec(-1));
1362 resolution kos = nres(MI,0);
1363 print(betti(kos),"betti");
1364 LIB "random.lib";
1365 matrix alpha0 = random(32002,10,3);
1366 module pres = module(alpha0)+kos[3];
1367 attrib(pres,"isHomog",intvec(1,1,1,1,1,1,1,1,1,1));
1368 resolution fcokernel = mres(pres,0);
1369 print(betti(fcokernel),"betti");
1370 module dir = transpose(pres);
1371 attrib(dir,"isHomog",intvec(-1,-1,-1,-2,-2,-2,
1372                             -2,-2,-2,-2,-2,-2,-2));
1373 resolution fdir = mres(dir,2);
1374 print(betti(fdir),"betti");
1375 ideal I = groebner(flatten(fdir[2]));
1376 resolution FI = mres(I,0);
1377 print(betti(FI),"betti");
1378 module F=FI[2];
1379 int t=timer;
1380 def A1=sheafCoh(F,-8,8);
1381 timer-t;
1382 t=timer;
1383 def A2=sheafCohBGG(F,-8,8);
1384 timer-t;
1385
1386 LIB "sheafcoh.lib";
1387 LIB "random.lib";
1388 ring S = 32003, x(0..4), dp;
1389 resolution kos = nres(maxideal(1),0);
1390 betti(kos);
1391 matrix kos5 = kos[5];
1392 matrix tphi = transpose(dsum(kos5,kos5));
1393 matrix kos3 = kos[3];
1394 matrix psi = dsum(kos3,kos3);
1395 matrix beta1 = random(32002,20,2);
1396 matrix tbeta1tilde = transpose(psi*beta1);
1397 matrix tbeta0 = lift(tphi,tbeta1tilde);
1398 matrix kos4 = kos[4];
1399 matrix tkos4pluskos4 = transpose(dsum(kos4,kos4));
1400 matrix tgammamin1 = random(32002,20,1);
1401 matrix tgamma0 = tkos4pluskos4*tgammamin1;
1402 matrix talpha0 = concat(tbeta0,tgamma0);
1403 matrix zero[20][1];
1404 matrix tpsi = transpose(psi);
1405 matrix tpresg = concat(tpsi,zero);
1406 matrix pres = module(transpose(talpha0))
1407                    + module(transpose(tpresg));
1408 module dir = transpose(pres);
1409 dir = prune(dir);
1410 homog(dir);
1411 intvec deg_dir = attrib(dir,"isHomog");
1412 attrib(dir,"isHomog",deg_dir-2);        // set degrees
1413 resolution fdir = mres(prune(dir),2);
1414 print(betti(fdir),"betti");
1415 ideal I = groebner(flatten(fdir[2]));
1416 resolution FI = mres(I,0);
1417
1418 module F=FI[2];
1419 def A1=sheafCoh(F,-5,7);
1420 def A2=sheafCohBGG(F,-5,7);
1421
1422*/
Note: See TracBrowser for help on using the repository browser.