Changeset 08fa62 in git


Ignore:
Timestamp:
Oct 14, 2013, 6:15:36 PM (11 years ago)
Author:
Hans Schoenemann <hannes@…>
Branches:
(u'fieker-DuVal', '117eb8c30fc9e991c4decca4832b1d19036c4c65')(u'spielwiese', 'b21a664aa22dc6e196223af8a74ad4885e83547c')
Children:
76d26c42ef4754edaa37b1777b4d2831f6f7e349
Parents:
ac665c6c03e4e4354b586d1fd770d8d5669253cf
Message:
fix: examples from derham.lib
File:
1 edited

Legend:

Unmodified
Added
Removed
  • Singular/LIB/derham.lib

    rac665c r08fa62  
    11///////////////////////////////////////////////////////////////////////////////
    2 version="";
    3 category="Noncommutative";
    4 info="
    5 LIBRARY:  derham.lib      Computation of deRham cohomology
    6 AUTHORS:  Cornelia Rottner, rottner@mathematik.uni-kl.de
    7 OVERVIEW:
    8 PROCEDURES:
    9 
    10 ";
    11 
    12 
    13 LIB "nctools.lib";
    14 LIB "matrix.lib";
    15 LIB "qhmoduli.lib";
    16 LIB "general.lib";
    17 LIB "dmod.lib";
    18 LIB "bfun.lib";
    19 LIB "dmodapp.lib";
    20 LIB "poly.lib";
    21 
    22 /////////////////////////////////////////////////////////////////////////////
    23 
    24 static proc divdr(matrix m, matrix n)
    25 {
    26   m=transpose(m);
    27   n=transpose(n);
    28   matrix con=concat(m,n);
    29   matrix s=syz(con);
    30   s=submat(s,1..ncols(m),1..ncols(s));
    31   s=transpose(compress(s));
    32   return(s);
    33 }
    34 
    35 /////////////////////////////////////////////////////////////////////////////
    36 
    37 
    38 static proc matrixlift(matrix M, matrix N)
    39 {
    40   // option(noreturnSB);
    41   matrix l=transpose(lift(transpose(M),transpose(N)));
    42   return(l);
    43 }
    44 
    45 ///////////////////////////////////////////////////////////////////////////////
    46 
    47 proc shortexactpieces(list #)
    48 {
    49   matrix  Bnew= divdr(#[2],#[3]);
    50   matrix Bold=Bnew;
    51   matrix Z=divdr(Bnew,#[1]);
    52   list bzh; list zcb;
    53   bzh=list(list(),list(),Z,unitmat(ncols(Z)),Z);
    54   zcb=(Z, Bnew, #[1], unitmat(ncols(#[1])), Bnew);
    55   list sep;
    56   sep[1]=(list(bzh,zcb));
    57   int i;
    58   list out;
    59   for (i=3; i<=size(#)-2; i=i+2)
    60     {
    61       out=bzhzcb(Bold, #[i-1] , #[i], #[i+1], #[i+2]);
    62       sep[size(sep)+1]=out[1];
    63       Bold=out[2];
    64     }
    65   bzh=(divdr(#[size(#)-2], #[size(#)-1]),#[size(#)-2], #[size(#)-1],unitmat(ncols(#[size(#)-1])),transpose(concat(transpose(#[size(#)-2]),transpose(#[size(#)-1]))));
    66   zcb=(#[size(#)-1], unitmat(ncols(#[size(#)-1])), #[size(#)-1],list(),list());
    67   sep[size(sep)+1]=list(bzh,zcb);
    68   return(sep);
    69 }
    70 
    71 ////////////////////////////////////////////////////////////////////////////////////////
    72 
    73 static proc bzhzcb (matrix Bold, matrix f0, matrix C1, matrix f1,matrix C2)
    74 {
    75   matrix Bnew=divdr(f1,C2);
    76   matrix Z= divdr(Bnew,C1);
    77   matrix lift1= matrixlift(Bnew,f0);
    78   list bzh=(Bold, lift1, Z, unitmat(ncols(Z)), transpose(concat(transpose(lift1),transpose(Z))));
    79   list zcb=(Z, Bnew, C1, unitmat(ncols(C1)),Bnew);
    80   list out=(list(bzh, zcb), Bnew);
    81   return(out);
    82 }
    83 
    84 //////////////////////////////////////////////////////////////////////////////////////
    85 
    86 proc VdstrictGB (matrix M, int d ,list #);
    87 "USAGE:VdstrictGB(M,d[,v]); M a matrix, d an integer, v an optional intvec(shift vector)
    88 RETURN:matrix M; the rows of M foem a Vd-strict Groebner basis for imM
    89 ASSUME:1<=d<=nvars(basering)/2; size(v)=ncols(M)
    90 "
    91 {
    92   if (M==matrix(0,nrows(M),ncols(M)))
    93     {
    94       return (matrix(0,1,ncols(M)));
    95     }
    96   def W =basering;
    97   int ncM=ncols(M);
    98   list Data=ringlist(W);
    99   Data[2]=list("nhv")+Data[2];
    100   Data[3][3]=Data[3][1];
    101   Data[3][1]=Data[3][2];
    102   Data[3][2]=list("dp",intvec(1));
    103   matrix re[size(Data[2])][size(Data[2])]=UpOneMatrix(size(Data[2]));
    104   Data[5]=re;
    105   int k; int l;
    106   Data[6]=transpose(concat(matrix(0,1,1),transpose(concat(matrix(0,1,1),Data[6]))));
    107   def Whom=ring(Data);
    108   setring Whom;
    109   matrix Mnew=imap(W,M);
    110   intvec v;
    111   if (size(#)!=0)
    112     {
    113       v=#[1];
    114     }
    115   if (size(v) < ncM)
    116     {
    117       v=v,0:(ncM-size(v));
    118     }
    119   Mnew=homogenize(Mnew, d, v);
    120   Mnew=transpose(Mnew);
    121   Mnew=std(Mnew);
    122   Mnew=subst(Mnew,nhv,1);
    123   Mnew=transpose(Mnew);
    124   setring W;
    125   M=imap(Whom,Mnew);
    126   return(M);
    127 }
    128 
    129 ////////////////////////////////////////////////////////////////////////////////////
    130 
    131 static proc Vdnormalform(matrix F, matrix M, int d, intvec v)
    132 {
    133   def W =basering;
    134   int c=ncols(M);
    135   F=submat(F,intvec(1..nrows(F)),intvec(1..c));
    136   list Data=ringlist(W);
    137   Data[2]=list("nhv")+Data[2];
    138   Data[3][3]=Data[3][1];
    139   Data[3][1]=Data[3][2];
    140   Data[3][2]=list("dp",intvec(1));
    141   matrix re[size(Data[2])][size(Data[2])]=UpOneMatrix(size(Data[2]));
    142   Data[5]=re;
    143   int k;
    144   int l;
    145   matrix rep[size(Data[2])][size(Data[2])];
    146   for (l=size(Data[2])-1;l>=1; l--)
    147     {
    148       for (k=l-1; k>=1;k--)
    149         {
    150           rep[k+1,l+1]=Data[6][k,l];
    151         }
    152     }
    153   Data[6]=rep;
    154   def Whom=ring(Data);
    155   setring Whom;
    156   matrix Mnew=imap(W,M);
    157   Mnew=(homogenize(Mnew, d, v));//doppelte Berechung unnötig->muss noch geändert werden!!!
    158   matrix Fnew=imap(W,F);
    159   matrix Fb;
    160   for (l=1; l<=nrows(Fnew); l++)
    161     {
    162       Fb=homogenize(submat(Fnew,l,intvec(1..ncols(Fnew))),d,v);
    163       Fb=transpose(reduce(transpose(Fb),std(transpose(Mnew))));// doppelte Berechnung unnötig, unterdrückt aber Fehler meldung
    164       for (k=1; k<=ncols(Fnew);k++)
    165         {
    166           Fnew[l,k]=Fb[1,k];
    167         }
    168     }
    169   Fnew=subst(Fnew,nhv,1);
    170   setring W;
    171   F=imap(Whom,Fnew);
    172   return(F);
    173 }
    174 
    175 
    176 ///////////////////////////////////////////////////////////////////////////////
    177 
    178 static proc homogenize (matrix M, int d, intvec v)
    179 {
    180   int l; poly f; int s; int i; intvec vnm;int kmin; list findmin;
    181   int n=(nvars(basering)-1) div 2;
    182   list rempoly;
    183   list remk;
    184   list rem1;
    185   list rem2;
    186   for (int k=1; k<=nrows(M); k++) //man könnte auch paralell immer weiter homogenisieren, d.h. immer ein enues Minimum finden und das dann machen
    187     {
    188       for (l=1; l<=ncols (M); l++)
    189         {
    190           f=M[k,l];
    191           s=size(f);
    192           for (i=1; i<=s; i++)
    193             {
    194               vnm=leadexp(f);
    195               vnm=vnm[n+2..n+d+1]-vnm[2..d+1];
    196               kmin=sum(vnm)+v[l];
    197               rem1[size(rem1)+1]=lead(f);
    198               rem2[size(rem2)+1]=kmin;
    199               findmin=insert(findmin,kmin);
    200               f=f-lead(f);
    201             }
    202           rempoly[l]=rem1;
    203           remk[l]=rem2;
    204           rem1=list();
    205           rem2=list();
    206         }
    207       if (size(findmin)!=0)
    208         {
    209           kmin=Min(findmin);
    210         }
    211       for (l=1; l<=ncols(M); l++)
    212         {
    213           if (M[k,l]!=0)
    214             {
    215               M[k,l]=0;
    216               for (i=1; i<=size(rempoly[l]);i++)
    217                 {
    218                   M[k,l]=M[k,l]+nhv^(remk[l][i]-kmin)*rempoly[l][i];
    219                 }
    220             }
    221         }
    222       rempoly=list();
    223       remk=list();
    224       findmin=list();
    225     }
    226   return(M);
    227 }
    228 
    229 
    230 //////////////////////////////////////////////////////////////////////////////////////
    231 
    232 static proc soldr (matrix M, matrix N)
    233 {
    234   int n=nrows(M);
    235   int q=ncols(M);
    236   matrix S=concat(transpose(M),transpose(N));
    237   def W=basering;
    238   list Data=ringlist(W);
    239   list Save=Data[3];
    240   Data[3]=list(list("c",0),list("dp",intvec(1..nvars(W))));
    241   def Wmod=ring(Data);
    242   setring Wmod;
    243   matrix Smod=imap(W,S);
    244   matrix E[q][1];
    245   matrix Smod2;
    246   matrix Smodnew;
    247   option(returnSB);
    248   int i; int j;
    249   for (i=1;i<=q;i++)
    250     {
    251       E[i,1]=1;
    252       Smod2=concat(E,Smod);
    253       print (Smod2);
    254       Smod2=syz(Smod2);
    255       E[i,1]=0;
    256       for (j=1;j<=ncols(Smod2);j++)
    257         {
    258           if (Smod2[1,j]==1)
    259             {
    260               Smodnew=concat(Smodnew,(-1)*(submat(Smod2,intvec(2..n+1),j)));
    261               break;
    262             }
    263         }
    264     }
    265   Smodnew=transpose(submat(Smodnew,intvec(1..n),intvec(2..q+1)));
    266   setring W;
    267   matrix  Snew=imap(Wmod,Smodnew);
    268   return (Snew);
    269 }
    270 
    271 
    272 /////////////////////////////////////////////////////////////////////////////
    273 
    274 proc toVdstrictsequence (list C,int n, intvec v)
    275 
    276 {
    277   matrix J_C=VdstrictGB(C[5],n,list(v));
    278   matrix J_A=C[1];
    279   matrix f_CB=C[4];
    280   matrix f_ACB=transpose(concat(transpose(C[2]),transpose(f_CB)));
    281   matrix J_AC=divdr(f_ACB,C[3]);
    282   matrix P=matrixlift(J_AC * prodr(ncols(C[1]),ncols(C[5])) ,J_C);
    283   list storePi;
    284   matrix Pi[1][ncols(J_AC)];
    285   int i;int j;
    286   for (i=1; i<=nrows(J_C); i++)
    287     {
    288       for (j=1; j<=nrows(J_AC);j++)
    289         {
    290           Pi=Pi+P[i,j]*submat(J_AC,j,intvec(1..ncols(J_AC)));
    291         }
    292       storePi[i]=Pi;
    293       Pi=0;
    294     }
    295   intvec m_a;
    296   list findMin;
    297   int comMin;
    298   for (i=1; i<=ncols(J_A); i++)
    299     {
    300       for (j=1; j<=size(storePi);j++)
    301         {
    302           if (storePi[j][1,i]!=0)
    303             {
    304               comMin=Vddeg(storePi[j]*prodr(ncols(J_A),ncols(C[5])),n,v)-Vddeg(storePi[j][1,i],n,intvec(0));
    305               findMin[size(findMin)+1]=comMin;
    306             }
    307         }
    308       if (size(findMin)!=0)
    309         {
    310           m_a[i]=Min(findMin);
    311           findMin=list();
    312         }
    313       else
    314         {
    315           m_a[i]=0;
    316         }
    317     }
    318   matrix zero[ncols(J_A)][ncols(J_C)];
    319 
    320   matrix g_AB=concat(unitmat(ncols(J_A)),zero);
    321   matrix g_BC= transpose(concat(transpose(zero),transpose(unitmat(ncols(J_C)))));
    322   intvec m_b=m_a,v;
    323 
    324   J_A=VdstrictGB(J_A,n,m_a);
    325   J_AC=transpose(storePi[1]);
    326   for (i=2; i<= size(storePi); i++)
    327     {
    328       J_AC=concat(J_AC, transpose(storePi[i]));
    329     }
    330   J_AC=transpose(concat(transpose(matrix(J_A,nrows(J_A),nrows(J_AC))),J_AC));
    331 
    332   list Vdstrict=(list(J_A),list(g_AB),list(J_AC),list(g_BC),list(J_C),list(m_a),list(m_b),list(v));
    333   return (Vdstrict);
    334 }
    335 
    336 
    337 /////////////////////////////////////////////////////////////////////////
    338 
    339 static proc prodr (int k, int l)
    340 {
    341   if (k==0)
    342     {
    343       matrix P=unitmat(l);
    344       return (P);
    345     }
    346   matrix O[l][k];
    347   matrix P=transpose(concat(O,unitmat(l)));
    348   return (P);
    349 }
    350 
    351 /////////////////////////////////////////////////////////////////////////
    352 
    353 proc Vddeg(matrix M, int d, intvec v, list #)//Aternative: in WHom Leadmonom ausrechnen!!
    354   "USAGE: Vddeg(M,d,v); M 1xr-matrix, d int, v intvec of size r
    355 RETURN:int; the Vd-degree of M
    356 "
    357 {
    358   int i;int j;
    359   int n=nvars(basering) div 2;
    360   intvec  e;
    361   int etoint;
    362   list findmax;
    363   int c=ncols(M);
    364   poly l;
    365   list positionpoly;
    366   list positionVd;
    367   for (i=1; i<=c; i++)
    368     {
    369       positionpoly[i]=list();
    370       positionVd[i]=list();
    371       while (M[1,i]!=0)
    372         {
    373           l=lead(M[1,i]);
    374           positionpoly[i][size(positionpoly[i])+1]=l;
    375           e=leadexp(l);
    376           e=-e[1..d]+e[n+1..n+d];
    377           e=sum(e)+v[i];
    378           etoint=e[1];
    379           positionVd[i][size(positionVd[i])+1]=etoint;
    380           findmax[size(findmax)+1]=etoint;
    381           M[1,i]=M[1,i]-l;
    382         }
    383     }
    384   if (size(findmax)!=0)
    385     {
    386       int maxVd=Max(findmax);
    387       if (size(#)==0)
    388         {
    389           return (maxVd);
    390         }
    391     }
    392   else // M is 0-modul
    393     {
    394       return(int(0));
    395     }
    396   l=0;
    397   for (i=c; i>=1; i--)
    398     {
    399       for (j=1; j<=size(positionVd[i]); j++)
    400         {
    401           if (positionVd[i][j]==maxVd)
    402             {
    403               l=l+positionpoly[i][j];
    404             }
    405         }
    406       if (l!=0)
    407         {
    408           return (list(l,i));
    409         }
    410     }
    411 
    412 }
    413 
    414 
    415 ///////////////////////////////////////////////////////////////////////////////
    416 
    417 proc toVdstrictsequences (list L,int d, intvec v)
    418   "USAGE: toVdstrictsequences(L,d,v); L list, d int, v intvec, L contains two lists of short exact sequences(D,f_DA,A,f_AF,F) and (A,f_AB,B,f_BC,C), v is a shift vector on the range of C
    419 RETURN: list of two lists; each lists contains Vd-strict exact sequences with corresponding shift vectors
    420 "
    421 {
    422   matrix J_F=L[1][5];
    423   matrix J_D=L[1][1];
    424   matrix f_FA=L[1][4];
    425   matrix f_DFA=transpose(concat(transpose(L[1][2]),transpose(f_FA)));
    426   matrix J_DF=divdr(f_DFA,L[1][3]);
    427   matrix J_C=L[2][5];
    428   matrix f_CB=L[2][4];
    429   matrix f_DFCB=transpose(concat(transpose(f_DFA*L[2][2]),transpose(f_CB)));
    430   matrix J_DFC=divdr(f_DFCB,L[2][3]);
    431   matrix P=matrixlift(J_DFC*prodr(ncols(J_DF),ncols(L[2][5])),J_C);
    432   list storePi;
    433   matrix Pi[1][ncols(J_DFC)];
    434   int i; int j;
    435   for (i=1; i<=nrows(J_C); i++)
    436     {
    437 
    438       for (j=1; j<=nrows(J_DFC);j++)
    439         {
    440           Pi=Pi+P[i,j]*submat(J_DFC,j,intvec(1..ncols(J_DFC)));
    441         }
    442       storePi[i]=Pi;
    443       Pi=0;
    444     }
    445   intvec m_a;
    446   list findMin;
    447   list noMin;
    448   int comMin;
    449   for (i=1; i<=ncols(J_DF); i++)
    450     {
    451       for (j=1; j<=size(storePi);j++)
    452         {
    453           if (storePi[j][1,i]!=0)
    454             {
    455               comMin=Vddeg(storePi[j]*prodr(ncols(J_DF),ncols(J_C)),d,v)-Vddeg(storePi[j][1,i],d,intvec(0));
    456               findMin[size(findMin)+1]=comMin;
    457             }
    458         }
    459       if (size(findMin)!=0)
    460         {
    461           m_a[i]=Min(findMin);
    462           findMin=list();
    463           noMin[i]=0;
    464         }
    465       else
    466         {
    467           noMin[i]=1;
    468         }
    469     }
    470   if (size(m_a) < ncols(J_DF))
    471     {
    472       m_a[ncols(J_DF)]=0;
    473     }
    474   intvec m_f=m_a[ncols(J_D)+1..size(m_a)];
    475   J_F=VdstrictGB(J_F,d,m_f);
    476   P=matrixlift(J_DF * prodr(ncols(L[1][1]),ncols(L[1][5])) ,J_F);// selbe Prinzip wie oben--> evtl auslagern
    477   list storePinew;
    478   matrix Pidf[1][ncols(J_DF)];
    479   for (i=1; i<=nrows(J_F); i++)
    480     {
    481       for (j=1; j<=nrows(J_DF);j++)
    482         {
    483           Pidf=Pidf+P[i,j]*submat(J_DF,j,intvec(1..ncols(J_DF)));
    484         }
    485       storePinew[i]=Pidf;
    486       Pidf=0;
    487     }
    488   intvec m_d;
    489   for (i=1; i<=ncols(J_D); i++)
    490     {
    491       for (j=1; j<=size(storePinew);j++)
    492         {
    493           if (storePinew[j][1,i]!=0)
    494             {
    495               comMin=Vddeg(storePinew[j]*prodr(ncols(J_D),ncols(L[1][5])),d,m_f)-Vddeg(storePinew[j][1,i],d,intvec(0));
    496               findMin[size(findMin)+1]=comMin;
    497             }
    498         }
    499       if (size(findMin)!=0)
    500         {
    501           if (noMin[i]==0)
    502             {
    503               m_d[i]=Min(insert(findMin,m_a[i]));
    504               m_a[i]=m_d[i];
    505             }
    506           else
    507             {
    508               m_d[i]=Min(findMin);
    509               m_a[i]=m_d[i];
    510             }
    511         }
    512       else
    513         {
    514           m_d[i]=m_a[i];
    515         }
    516       findMin=list();
    517     }
    518   J_D=VdstrictGB(J_D,d,m_d);
    519   J_DF=transpose(storePinew[1]);
    520   for (i=2; i<=nrows(J_F); i++)
    521     {
    522       J_DF=concat(J_DF,transpose(storePinew[i]));
    523     }
    524   J_DF=transpose(concat(transpose(matrix(J_D,nrows(J_D),nrows(J_DF))),J_DF));
    525   J_DFC=transpose(storePi[1]);
    526   for (i=2; i<=nrows(J_C); i++)
    527     {
    528       J_DFC=concat(J_DFC,transpose(storePi[i]));
    529     }
    530   J_DFC=transpose(concat(transpose(matrix(J_DF,nrows(J_DF),nrows(J_DFC))),J_DFC));
    531   intvec m_b=m_a,v;
    532   matrix zero[ncols(J_D)][ncols(J_F)];
    533   matrix g_DA=concat(unitmat(ncols(J_D)),zero);
    534   matrix g_AF=transpose(concat(transpose(zero),unitmat(ncols(J_F))));
    535   matrix zero1[ncols(J_DF)][ncols(J_C)];
    536   matrix g_AB=concat(unitmat(ncols(J_DF)),zero1);
    537   matrix g_BC=transpose(concat(transpose(zero1),unitmat(ncols(J_C))));
    538   list out=(list(list(J_D),list(g_DA),list(J_DF),list(g_AF),list(J_F),list(m_d),list(m_a),list(m_f)),list(list(J_DF),list(g_AB),list(J_DFC),list(g_BC),list(J_C),list(m_a),list(m_b),list(v)));
    539   return(out),
    540     }
    541 
    542 ///////////////////////////////////////////////////////////////////////////////////////////
    543 
    544 proc shortexactpiecestoVdstrict(list C, int d,list #)
    545 
    546 {
    547 
    548   int s =size(C);
    549   if (size(#)==0)
    550     {
    551       intvec v=0:ncols(C[s][1][5]);
    552     }
    553   else
    554     {
    555       intvec v=#[1];
    556     }
    557   list out;
    558   out[s]=list(toVdstrictsequence(C[s][1],d,v));
    559   out[s][2]=list(list(out[s][1][3][1]),list(unitmat(ncols(out[s][1][3][1]))),list(out[s][1][3][1]),list(list()),list(list()));
    560   out[s][2][6]=list(out[s][1][7][1]);
    561   out[s][2][7]=list(out[s][2][6][1]);
    562   out[s][2][8]=list(list());
    563   int i;
    564   for (i=s-1; i>=2; i--)
    565     {
    566       C[i][2][5]=out[i+1][1][1][1];
    567       out[i]=toVdstrictsequences(C[i],d,out[i+1][1][6][1]);
    568     }
    569   out[1]=list(list());
    570   out[1][2]=toVdstrictsequence(C[1][2],d,out[2][1][6][1]);
    571   out[1][1][3]=list(out[1][2][1][1]);
    572   out[1][1][5]=list(out[1][2][1][1]);
    573   out[1][1][4]=list(unitmat(ncols(out[1][1][3][1])));
    574   out[1][1][7]=list(out[1][2][6][1]);
    575   out[1][1][8]=list(out[1][2][6][1]);
    576   out[1][1][1]=list(list());
    577   out[1][1][2]=list(list());
    578   out[1][1][6]=list(list());
    579   list Hi;
    580   for (i=1; i<=size(out); i++)
    581     {
    582       Hi[i]=list(out[i][1][5][1],out[i][1][8][1]);
    583     }
    584   list outall;
    585   outall[1]=out;
    586   print (out);
    587   outall[2]=Hi;
    588   return(outall);
    589 
    590 
    591 }
    592 
    593 ///////////////////////////////////////////////////////////////////////////////////////////
    594 
    595 proc toVdstrict2x3complex(list L, int d, list #)
    596 {
    597   matrix rem; int i; int j;
    598   list J_A=list(list());
    599   list J_B=list(list());
    600   list J_C=list(list());
    601   list g_AB=list(list());
    602   list g_BC=list(list());
    603   list n_a=list(list());
    604   list n_b=list(list());
    605   list n_c=list(list());
    606   intvec n_b1;
    607   if (size(L[5])!=0)
    608     {
    609       intvec n_c1;
    610       for (i=1; i<=nrows(L[5]); i++)
    611         {
    612           rem=submat(L[5],i,intvec(1..ncols(L[5])));
    613           n_c1[i]=Vddeg(rem,d, L[8]);
    614         }
    615       n_c[1]=n_c1;
    616       J_C[1]=transpose(syz(transpose(L[5])));
    617       if (J_C[1]!=matrix(0,nrows(J_C[1]),ncols(J_C[1])))
    618         {
    619           J_C[1]=VdstrictGB(J_C[1],d,n_c1);
    620           if (size(#[2])!=0)
    621             {
    622               n_a[1]=#[2];
    623               n_b1=n_a[1],n_c[1];
    624               n_b[1]=n_b1;
    625               matrix zero[nrows(L[1])][nrows(L[5])];
    626               g_AB=concat(unitmat(nrows(L[1])),matrix(0,nrows(L[1]),nrows(L[5])));
    627               if (size(#[1])!=0)
    628                 {
    629                   J_A=#[1];
    630                   J_B=transpose(matrix(syz(transpose(L[3]))));
    631                   matrix P=matrixlift(J_B[1] * prodr(nrows(L[1]),nrows(L[5])) ,J_C[1]);
    632 
    633                   matrix Pi[1][ncols(J_B[1])];
    634                   matrix Picombined;
    635                   for (i=1; i<=nrows(J_C[1]); i++)
    636                     {
    637                       for (j=1; j<=nrows(J_B[1]);j++)
    638                         {
    639                           Pi=Pi+P[i,j]*submat(J_B[1],j,intvec(1..ncols(J_B[1])));
    640 
    641                         }
    642                       if (i==1)
    643                         {
    644                           Picombined=transpose(Pi);
    645                         }
    646                       else
    647                         {
    648                           Picombined=concat(Picombined,transpose(Pi));
    649                         }
    650                       Pi=0;
    651                     }
    652                   Picombined=transpose(Picombined);
    653                   Picombined=concat(Vdnormalform(Picombined,J_A[1],d,n_a[1]),submat(Picombined,intvec(1..nrows(Picombined)),intvec((ncols(J_A[1])+1)..ncols(Picombined))));
    654                   J_B[1]=transpose(concat(transpose(matrix(J_A[1],nrows(J_A[1]),ncols(J_B[1]))),transpose(Picombined)));
    655                   g_BC=transpose(concat(transpose(zero),unitmat(nrows(L[5]))));
    656                 }
    657               else
    658                 {
    659                   J_B[1]=concat(matrix(0,nrows(J_C[1]),nrows(L[3])-nrows(L[5])),J_C[1]);
    660                   g_BC=transpose(concat(transpose(zero),unitmat(nrows(L[5]))));
    661                 }
    662             }
    663           else
    664             {
    665               n_b=n_c[1];
    666               J_B[1]=J_C[1];
    667               g_BC=unitmat(ncols(J_C[1]));
    668             }
    669         }
    670       else
    671         {
    672           J_C=list(list());
    673           if (size(#[2])!=0)
    674             {
    675               matrix zero[nrows(L[1])][nrows(L[5])];
    676               g_BC=transpose(concat(transpose(zero),unitmat(nrows(L[5]))));
    677               n_a[1]=#[2];
    678               n_b1=n_a[1],n_c[1];
    679               n_b[1]=n_b1;
    680               g_AB=concat(unitmat(nrows(L[1])),matrix(0,nrows(L[1]),nrows(L[5])));;
    681 
    682               if (size(#[1])!=0)
    683                 {
    684                   J_A=#[1];
    685                   J_B=concat(J_A[1],matrix(0,nrows(J_A[1]),nrows(L[3])-nrows(L[1])));
    686                 }
    687             }
    688           else
    689             {
    690               n_b=n_c[1];
    691               g_BC=unitmat(ncols(L[5]));
    692             }
    693 
    694         }
    695     }
    696   else
    697     {
    698       if (size(#[2])!=0)
    699         {
    700           n_a[1]=#[2];
    701           n_b=n_a[1];
    702           g_AB=unitmat(size(n_b[1]));
    703           if (size(#[1])!=0)
    704             {
    705               J_A=#[1];
    706               J_B[1]=J_A[1];
    707             }
    708         }
    709     }
    710   list out=(J_A[1],g_AB[1],J_B[1],g_BC[1],J_C[1],n_a[1],n_b[1],n_c[1]);
    711   return (out);
    712 }
    713 
    714 
    715 //////////////////////////////////////////////////////////////////////////
    716 
    717 proc Vdstrictdoublecompexes(list L, int d)
    718 {
    719   int i; int k; int c; int j;
    720   intvec n_b;
    721   matrix rem;
    722   matrix J_B;
    723   list store;
    724   int t=size(L)+nvars(basering) div 2-2;
    725   for (k=1; k<=(size(L)+nvars(basering) div 2-3); k++)//
    726     {
    727       L[1][1][1][k+1]=list();
    728       L[1][1][2][k+1]=list();
    729       L[1][1][6][k+1]=list();
    730       if (size(L[1][1][3][k])!=0)
    731         {
    732           for (i=1; i<=nrows(L[1][1][3][k]); i++)
    733             {
    734               rem=submat(L[1][1][3][k],i,(1..ncols(L[1][1][3][k])));
    735               n_b[i]=Vddeg(rem,d,L[1][1][7][k]);
    736             }
    737           J_B=transpose(syz(transpose(L[1][1][3][k])));
    738           L[1][1][7][k+1]=n_b;
    739           L[1][1][8][k+1]=n_b;
    740           L[1][1][4][k+1]=unitmat(nrows(L[1][1][3][k]));
    741           if (J_B!=matrix(0,nrows(J_B),ncols(J_B)))
    742             {
    743               J_B=VdstrictGB(J_B,d,n_b);
    744               L[1][1][3][k+1]=J_B;
    745               L[1][1][5][k+1]=J_B;
    746             }
    747           else
    748             {
    749               L[1][1][3][k+1]=list();
    750               L[1][1][5][k+1]=list();
    751             }
    752           n_b=0;
    753         }
    754       else
    755         {
    756           L[1][1][3][k+1]=list();
    757           L[1][1][5][k+1]=list();
    758           L[1][1][7][k+1]=list();
    759           L[1][1][8][k+1]=list();
    760           L[1][1][4][k+1]=list();
    761         }
    762       for (i=1; i<size(L); i++)
    763         {
    764           store=toVdstrict2x3complex(list(L[i][2][1][k],L[i][2][2][k],L[i][2][3][k],L[i][2][4][k],L[i][2][5][k],L[i][2][6][k],L[i][2][7][k],L[i][2][8][k]),d,L[i][1][3][k+1],L[i][1][7][k+1]);
    765           for (j=1; j<=8; j++)
    766             {
    767               L[i][2][j][k+1]=store[j];
    768             }
    769 
    770           store=toVdstrict2x3complex(list(L[i+1][1][1][k],L[i+1][1][2][k],L[i+1][1][3][k],L[i+1][1][4][k],L[i+1][1][5][k],L[i+1][1][6][k],L[i+1][1][7][k],L[i+1][1][8][k]),d,L[i][2][5][k+1],L[i][2][8][k+1]);
    771 
    772           for (j=1; j<=8; j++)
    773             {
    774               L[i+1][1][j][k+1]=store[j];
    775             }
    776         }
    777       if (size(L[size(L)][1][7][k+1])!=0)
    778         {
    779           L[size(L)][2][4][k+1]=list();
    780           L[size(L)][2][5][k+1]=list();
    781           L[size(L)][2][6][k+1]=L[size(L)][1][7][k+1];
    782           L[size(L)][2][7][k+1]=L[size(L)][1][7][k+1];
    783           L[size(L)][2][8][k+1]=list();
    784           L[size(L)][2][2][k+1]=unitmat(size(L[size(L)][1][7][k+1]));
    785 
    786           if (size(L[size(L)][1][3][k+1])!=0)
    787             {
    788               L[size(L)][2][1][k+1]=L[size(L)][1][3][k+1];
    789               L[size(L)][2][3][k+1]=L[size(L)][1][3][k+1];
    790             }
    791           else
    792             {
    793               L[size(L)][2][1][k+1]=list();
    794               L[size(L)][2][3][k+1]=list();
    795             }
    796         }
    797       else
    798         {
    799           for (j=1; j<=8; j++)
    800             {
    801               L[size(L)][2][j][k+1]=list();
    802             }
    803         }
    804     }
    805 
    806 
    807   k=t;
    808   intvec n_c;
    809   intvec vn_b;
    810   list N_b;
    811   int n;
    812   for (i=1; i<=size(L); i++)
    813     {
    814       for (n=1; n<=2; n++)
    815         {
    816           if (i==1 and n==1)
    817             {
    818               L[i][n][6][k+1]=list();
    819             }
    820           else
    821             {
    822               if (n==1)
    823                 {
    824                   L[i][1][6][k+1]=L[i-1][2][8][k+1];
    825                 }
    826               else
    827                 {
    828                   L[i][2][6][k+1]=L[i][1][7][k+1];
    829                 }
    830             }
    831           N_b[1]=L[i][n][6][k+1];
    832           if (size(L[i][n][5][k])!=0)
    833             {
    834               for (j=1; j<=nrows(L[i][n][5][k]); j++)
    835                 {
    836                   rem=submat(L[i][n][5][k],j,(1..ncols(L[i][n][5][k])));
    837                   n_c[j]=Vddeg(rem,d,L[i][n][8][k]);
    838                 }
    839               L[i][n][8][k+1]=n_c;
    840             }
    841           else
    842             {
    843               L[i][n][8][k+1]=list();
    844             }
    845           N_b[2]=L[i][n][8][k+1];
    846           n_c=0;
    847           if (size(N_b[1])!=0)
    848             {
    849               vn_b=N_b[1];
    850               if (size(N_b[2])!=0)
    851                 {
    852                   vn_b=vn_b,N_b[2];
    853                 }
    854               L[i][n][7][k+1]=vn_b;
    855             }
    856           else
    857             {
    858               if (size(N_b[2])!=0)
    859                 {
    860                   L[i][n][7][k+1]=N_b[2];
    861                 }
    862               else
    863                 {
    864                   L[i][n][7][k+1]=list();
    865                 }
    866             }
    867 
    868         }
    869     }
    870   return(L);
    871 }
    872 
    873 ////////////////////////////////////////////////////////////////////////////
    874 
    875 proc assemblingdoublecomplexes(list L)
    876 {
    877   list out;
    878   int i; int j;int k;int l; int oldj; int newj;
    879   for (i=1; i<=size(L); i++)
    880     {
    881       out[i]=list(list());
    882       out[i][1][1]=ncols(L[i][2][3][1]);
    883       if (size(L[i][2][5][1])!=0)
    884         {
    885           out[i][1][4]=prodr(ncols(L[i][2][3][1])-ncols(L[i][2][5][1]),ncols(L[i][2][5][1]));
    886         }
    887       else
    888         {
    889           out[i][1][4]=matrix(0,ncols(L[i][2][3][1]),1);
    890         }
    891 
    892       oldj=newj;
    893       for (j=1; j<=size(L[i][2][3]);j++)
    894         {
    895           out[i][j][2]=L[i][2][7][j];
    896           if (size(L[i][2][3][j])==0)
    897             {
    898               newj =j;
    899               break;
    900             }
    901           out[i][j+1]=list();
    902           out[i][j+1][1]=nrows(L[i][2][3][j]);
    903           out[i][j+1][3]=L[i][2][3][j];
    904           if (size(L[i][2][5][j])!=0)
    905             {
    906               out[i][j+1][4]=(-1)^j*prodr(nrows(L[i][2][3][j])-nrows(L[i][2][5][j]),nrows(L[i][2][5][j]));
    907             }
    908           else
    909             {
    910               out[i][j+1][4]=matrix(0,nrows(L[i][2][3][j]),1);
    911             }
    912           if(j==size(L[i][2][3]))
    913             {
    914               out[i][j+1][2]=L[i][2][7][j+1];
    915               newj=j+1;
    916             }
    917         }
    918       if (i>1)
    919         {
    920           for (k=1; k<=Min(list(oldj,newj)); k++)
    921             {
    922               out[i-1][k][4]=matrix(out[i-1][k][4],nrows(out[i-1][k][4]),out[i][k][1]);
    923             }
    924           for (k=newj+1; k<=oldj; k++)
    925             {
    926               out[i-1][k]=delete(out[i-1][k],4);
    927             }
    928         }
    929     }
    930   return (out);
    931 }
    932 
    933 //////////////////////////////////////////////////////////////////////////////
    934 
    935 proc totalcomplex(list L);
    936 {
    937   list out;intvec rem1;intvec v; list remsize; int emp;
    938   int i; int j; int c; int d; matrix M; int k; int l;
    939   int n=nvars(basering) div 2;
    940   list K;
    941   for (i=1; i<=n; i++)
    942     {
    943       K[i]=list();
    944     }
    945   L=K+L;
    946   for (i=1; i<=size(L); i++)
    947     {
    948       emp=0;
    949       if (size(L[i])!=0)
    950         {
    951           out[3*i-2]=L[i][1][1];
    952           v=L[i][1][1];
    953           rem1=L[i][1][2];
    954           emp=1;
    955         }
    956       else
    957         {
    958           out[3*i-2]=0;
    959           v=0;
    960         }
    961 
    962       for (j=i+1; j<=size(L); j++)
    963         {
    964           if (size(L[j])>=j-i+1)
    965             {
    966               out[3*i-2]=out[3*i-2]+L[j][j-i+1][1];
    967               if (emp==0)
    968                 {
    969                   rem1=L[j][j-i+1][2];
    970                   emp=1;
    971                 }
    972               else
    973                 {
    974                   rem1=rem1,L[j][j-i+1][2];
    975                 }
    976               v[size(v)+1]=L[j][j-i+1][1];
    977             }
    978           else
    979             {
    980               v[size(v)+1]=0;
    981             }
    982         }
    983       out[3*i-1]=rem1;
    984       v[size(v)+1]=0;
    985       remsize[i]=v;
    986     }
    987   int o1;
    988   int o2;
    989   for (i=1; i<=size(L)-1; i++)
    990     {
    991       o1=1;
    992       o2=1;
    993       if (size(out[3*i-2])!=0)
    994         {
    995           o1=out[3*i-2];
    996         }
    997       if (size(out[3*i+1])!=0)
    998         {
    999           o2=out[3*i+1];
    1000         }
    1001       M=matrix(0,o1,o2);
    1002       if (size(L[i])!=0)
    1003         {
    1004           if (size(L[i][1][4])!=0)
    1005             {
    1006               M=matrix(L[i][1][4],o1,o2);
    1007             }
    1008         }
    1009       c=remsize[i][1];
    1010       // d=remsize[i+1][1];
    1011       for (j=i+1; j<=size(L); j++)
    1012         {
    1013           if (remsize[i][j-i+1]!=0)
    1014             {
    1015               for (k=c+1; k<=c+remsize[i][j-i+1]; k++)
    1016                 {
    1017                   for (l=d+1; l<=d+remsize[i+1][j-i];l++)
    1018                     {
    1019                       M[k,l]=L[j][j-i+1][3][(k-c),(l-d)];
    1020                     }
    1021                 }
    1022               d=d+remsize[i+1][j-i];
    1023               if (remsize[i+1][j-i+1]!=0)
    1024                 {
    1025                   for (k=c+1; k<=c+remsize[i][j-i+1]; k++)
    1026                     {
    1027                       for (l=d+1; l<=d+remsize[i+1][j-i+1];l++)
    1028                         {
    1029                           M[k,l]=L[j][j-i+1][4][k-c,l-d];
    1030                         }
    1031                     }
    1032                   c=c+remsize[i][j-i+1];
    1033                 }
    1034             }
    1035           else
    1036             {
    1037               d=d+remsize[i+1][j-i];
    1038             }
    1039         }
    1040       out[3*i]=M;
    1041       d=0; c=0;
    1042     }
    1043   out[3*size(L)]=matrix(0,out[3*size(L)-2],1);
    1044   return (out);
    1045 }
    1046 
    1047 /////////////////////////////////////////////////////////////////////////////////////
    1048 
    1049 proc toVdstrictfreecomplex(list L,list #)
    1050 {
    1051   def B=basering;
    1052   int n=nvars(B) div 2+2;
    1053   int d=nvars(B) div 2;
    1054   intvec v;
    1055   list out;list outall;
    1056   int i;int j;
    1057   matrix mem;
    1058   int k;
    1059   if (size(#)!=0)
    1060     {
    1061       for (i=1; i<=size(#); i++)
    1062         {
    1063           if (typeof(#[i])==intvec)
    1064             {
    1065               v=#[i];
    1066             }
    1067           if (typeof(#[i])==int)
    1068             {
    1069               d=#[i];
    1070             }
    1071         }
    1072     }
    1073   if (size(L)==2)
    1074     {
    1075       v=(0:ncols(L[1]));
    1076       out[3*n-1]=v;
    1077       out[3*n-2]=ncols(L[1]);
    1078       out[3*n]=L[2];
    1079       out[3*n-3]=VdstrictGB(L[1],d,v);
    1080       for (i=n-1; i>=1; i--)
    1081         {
    1082           out[3*i-2]=nrows(out[3*i]);
    1083           v=0;
    1084           for (j=1; j<=out[3*i-2]; j++)
    1085             {
    1086               mem=submat(out[3*i],j,intvec(1..ncols(out[3*i])));
    1087               v[j]=Vddeg(mem,d, out[3*i+2]);
    1088             }
    1089           out[3*i-1]=v;
    1090           if (i!=1)
    1091             {
    1092               out[3*i-3]=transpose(syz(transpose(out[3*i])));
    1093               if (out[3*i-3]!=matrix(0,nrows(out[3*i-3]),ncols(out[3*i-3])))
    1094                 {
    1095                   out[3*i-3]=VdstrictGB(out[3*i-3],d,out[3*i-1]);
    1096                 }
    1097               else
    1098                 {
    1099                   out[3*i-3]=matrix(0,1,ncols(out[3*i-3]));
    1100                   out[3*i-4]=intvec(0);
    1101                   out[3*i-5]=int(0);
    1102                   for (j=i-2; j>=1; j--)
    1103                     {
    1104                       out[3*j]=matrix(0,1,1);
    1105                       out[3*j-1]=intvec(0);
    1106                       out[3*j-2]=int(0);
    1107                     }
    1108                   break;
    1109                 }
    1110             }
    1111         }
    1112       outall[1]=out;
    1113       outall[2]=list(list(out[3*n-3],out[3*n-1]));
    1114       return(outall);
    1115     }
    1116   out=shortexactpieces(L);
    1117   list rem;
    1118   if (v!=intvec(0:size(v)))
    1119     {
    1120       rem=shortexactpiecestoVdstrict(out,d,v);
    1121     }
    1122   else
    1123     {
    1124       rem=shortexactpiecestoVdstrict(out,d);
    1125     }
    1126   out=Vdstrictdoublecompexes(rem[1],d);
    1127   out=assemblingdoublecomplexes(out);
    1128   out=totalcomplex(out);
    1129   outall[1]=out;
    1130   outall[2]=rem[2];
    1131   return (outall);
    1132 }
    1133 
    1134 ////////////////////////////////////////////////////////////////////////////////
    1135 
    1136 proc derhamcohomology(list L)
    1137 {
    1138   def R=basering;
    1139   int n=nvars(R);int le=2*size(L)+n-1;
    1140   def W=makeWeyl(n);
    1141   setring W;
    1142   list man=ringlist(W);
    1143   if (n==1)
    1144     {
    1145       man[2][1]="x(1)";
    1146       man[2][2]="D(1)";
    1147       def Wi=ring(man);
    1148       setring Wi;
    1149       kill W;
    1150       def W=Wi;
    1151       setring W;
    1152       list man=ringlist(W);
    1153     }
    1154   man[2][size(man[2])+1]="s";;
    1155   man[3][3]=man[3][2];
    1156   man[3][2]=list("dp",intvec(1));
    1157   matrix N=UpOneMatrix(size(man[2]));
    1158   man[5]=N;
    1159   matrix M[1][1];
    1160   man[6]=transpose(concat(transpose(concat(man[6],M)),M));
    1161   def Ws=ring(man); setring R;  int r=size(L); int i;  int j;int k; int l; int count;  list Fi; list subsets; list maxnum; list bernsteinpolys; list annideals; list minint; list diffmaps;
    1162   for (i=1; i<=r; i++)
    1163     {
    1164       Fi[i]=list(); bernsteinpolys[i]=list(); annideals[i]=list(); subsets[i]=list();
    1165       maxnum[i]=list();
    1166       Fi[1][i]=L[i];
    1167       maxnum[1][i]=i;
    1168       subsets[1][i]=intvec(i);
    1169     }
    1170   intvec v;
    1171   for (i=2; i<=r; i++)
    1172     {
    1173       count=1;
    1174       for (j=1; j<=size(Fi[i-1]);j++)
    1175         {
    1176           for (k=maxnum[i-1][j]+1; k<=r; k++)
    1177             {
    1178               maxnum[i][count]=k;
    1179               v=subsets[i-1][j],k;
    1180               subsets[i][count]=v;
    1181               Fi[i][count]=lcm(Fi[i-1][j],L[k]);/////////
    1182               count=count+1;
    1183             }
    1184         }
    1185     }
    1186   for (i=1; i<=r; i++)
    1187     {
    1188       for (j=1; j<=size(Fi[i]); j++)
    1189         {
    1190           bernsteinpolys[i][j]=bfct(Fi[i][j])[1];
    1191           for (k=1; k<=ncols(bernsteinpolys[i][j]); k++)
    1192             {
    1193               if (isInt(number(bernsteinpolys[i][j][k]))==1)
    1194                 {
    1195                   minint[size(minint)+1]=int(bernsteinpolys[i][j][k]);
    1196                 }
    1197             }
    1198           def D=Sannfs(Fi[i][j]);
    1199           setring Ws;
    1200           annideals[i][j]=fetch(D,LD);
    1201           kill D;
    1202           setring R;
    1203         }
    1204     }
    1205   int m=Min(minint);
    1206   list zw;
    1207   for (i=1; i<r; i++)
    1208     {
    1209       diffmaps[i]=matrix(0,size(subsets[i]),size(subsets[i+1]));
    1210       for (j=1; j<=size(subsets[i]); j++)
    1211         {
    1212           for (k=1; k<=size(subsets[i+1]); k++)
    1213             {
    1214               zw=mysubset(subsets[i][j],subsets[i+1][k]);
    1215               diffmaps[i][j,k]=zw[2]*(L[zw[1]]/gcd(L[zw[1]],Fi[i][j]))^(-m);
    1216             }
    1217         }
    1218     }
    1219   diffmaps[r]=matrix(0,1,1);
    1220   setring Ws;
    1221   for (i=1; i<=r; i++)
    1222     {
    1223       for (j=1; j<=size(annideals[i]); j++)
    1224         {
    1225           annideals[i][j]=subst(annideals[i][j],s,m);
    1226         }
    1227     }
    1228   setring W;
    1229   list annideals=imap(Ws,annideals);
    1230   list diffmaps=fetch(R,diffmaps);
    1231   list fortoVdstrict;
    1232   ideal IFourier=var(n+1);
    1233   for (i=2;i<=n;i++)
    1234     {
    1235       IFourier=IFourier,var(n+i);
    1236     }
    1237   for (i=1; i<=n;i++)
    1238     {
    1239       IFourier=IFourier,-var(i);
    1240     }
    1241   map cFourier=W,IFourier;
    1242   matrix sup;
    1243   for (i=1; i<=r; i++)
    1244     {
    1245       sup=matrix(annideals[i][1]);
    1246       fortoVdstrict[2*i-1]=transpose(cFourier(sup));
    1247       for (j=2; j<=size(annideals[i]); j++)
    1248         {
    1249           sup=matrix(annideals[i][j]);
    1250           fortoVdstrict[2*i-1]=dsum(fortoVdstrict[2*i-1],transpose(cFourier(sup)));
    1251         }
    1252       sup=diffmaps[i];
    1253       fortoVdstrict[2*i]=cFourier(sup);
    1254     }
    1255   list rem=toVdstrictfreecomplex(fortoVdstrict);
    1256   list newcomplex=rem[1];
    1257   list minmaxk=globalbfun(rem[2]);
    1258   if (size(minmaxk)==0)
    1259     {
    1260       return (0);
    1261     }
    1262   list truncatedcomplex; list shorten; list  upto;
    1263   for (i=1; i<=size(newcomplex) div 3; i++)
    1264     {
    1265       shorten[3*i-1]=list();
    1266       for (j=1; j<=size(newcomplex[3*i-1]); j++)
    1267         {
    1268           shorten[3*i-1][j]=list(minmaxk[1]-newcomplex[3*i-1][j]+1,minmaxk[2]-newcomplex[3*i-1][j]+1);
    1269           upto[size(upto)+1]=shorten[3*i-1][j][2];
    1270           if (shorten[3*i-1][j][2]<=0)
    1271             {
    1272               shorten[3*i-1][j]=list();
    1273             }
    1274           else
    1275             {
    1276               if (shorten[3*i-1][j][1]<=0)
    1277                 {
    1278                   shorten[3*i-1][j][1]=1;
    1279                 }
    1280             }
    1281         }
    1282     }
    1283   int iupto=Max(upto);
    1284   if (iupto<=0)
    1285     {
    1286       /////die Kohomologie ist dann überall 0, muss noch entsprechend ausgegeben werden
    1287     }
    1288   list allpolys;
    1289   allpolys[1]=list(1);
    1290   list minvar;
    1291   minvar[1]=list(1);
    1292   for (i=1; i<=iupto-1; i++)
    1293     {
    1294       allpolys[i+1]=list();
    1295       minvar[i+1]=list();
    1296       for (k=1; k<=size(allpolys[i]); k++)
    1297         {
    1298           for (j=minvar[i][k]; j<=nvars(W) div 2; j++)
    1299             {
    1300               allpolys[i+1][size(allpolys[i+1])+1]=allpolys[i][k]*D(j);
    1301               minvar[i+1][size(minvar[i+1])+1]=j;
    1302             }
    1303         }
    1304     }
    1305   list keepformatrix;list sizetruncom;int stc;list fortrun;
    1306   for (i=1; i<=size(newcomplex) div 3; i++)
    1307     {
    1308       truncatedcomplex[2*i-1]=list();
    1309       sizetruncom[2*i-1]=list();
    1310       sizetruncom[2*i]=list();
    1311       truncatedcomplex[2*i]=newcomplex[3*i];
    1312       v=0;count=0;
    1313       sizetruncom[2*i][1]=0;
    1314       for (j=1; j<=newcomplex[3*i-2]; j++)
    1315         {
    1316           if (size(shorten[3*i-1][j])!=0)
    1317             {
    1318               fortrun=sublist(allpolys,shorten[3*i-1][j][1],shorten[3*i-1][j][2]);
    1319               truncatedcomplex[2*i-1][size(truncatedcomplex[2*i-1])+1]=fortrun[1];
    1320               count=count+fortrun[2];
    1321               sizetruncom[2*i-1][size(sizetruncom[2*i-1])+1]=list(int(shorten[3*i-1][j][1])-1,int(shorten[3*i-1][j][2])-1);
    1322               sizetruncom[2*i][size(sizetruncom[2*i])+1]=count;
    1323               if (v!=0)
    1324                 {
    1325                   v[size(v)+1]=j;
    1326                 }
    1327               else
    1328                 {
    1329                   v[1]=j;
    1330                 }
    1331             }
    1332         }
    1333 
    1334       if (v!=0)
    1335         {
    1336           truncatedcomplex[2*i]=submat(truncatedcomplex[2*i],v,1..ncols(truncatedcomplex[2*i]));
    1337           if (i!=1)
    1338             {
    1339               truncatedcomplex[2*(i-1)]=submat(truncatedcomplex[2*(i-1)],1..nrows(truncatedcomplex[2*(i-1)]),v);
    1340             }
    1341         }
    1342       else
    1343         {
    1344           truncatedcomplex[2*i]=matrix(0,1,ncols(truncatedcomplex[2*i]));
    1345           if (i!=1)
    1346             {
    1347               truncatedcomplex[2*(i-1)]=matrix(0,nrows(truncatedcomplex[2*(i-1)]),1);
    1348             }
    1349         }
    1350     }
    1351   int b;int d;poly form;poly lform; poly nform;int ideg;int kplus; int lplus;
    1352   for (i=1; i<size(truncatedcomplex) div 2; i++)
    1353     {
    1354       M=matrix(0,max(1,sizetruncom[2*i][size(sizetruncom[2*i])]),sizetruncom[2*i+2][size(sizetruncom[2*i+2])]);
    1355       for (k=1; k<=size(truncatedcomplex[2*i-1]);k++)
    1356         {
    1357           for (l=1; l<=size(truncatedcomplex[2*(i+1)-1]); l++)
    1358             {
    1359               if (size(sizetruncom[2*i])!=1)//?
    1360                 {
    1361                   for (j=1; j<=size(truncatedcomplex[2*i-1][k]); j++)
    1362                     {
    1363                       for (b=1; b<=size(truncatedcomplex[2*i-1][k][j]); b++)
    1364                         {
    1365                           form=truncatedcomplex[2*i-1][k][j][b][1]*truncatedcomplex[2*i][k,l];
    1366                           while (form!=0)
    1367                             {
    1368                               lform=lead(form);
    1369                               v=leadexp(lform);
    1370                               v=v[1..n];
    1371                               if (v==(0:n))
    1372                                 {
    1373                                   ideg=deg(lform)-sizetruncom[2*(i+1)-1][l][1];
    1374                                   if (ideg>=0)
    1375                                     {
    1376                                       for (d=1; d<=size(truncatedcomplex[2*(i+1)-1][l][ideg+1]);d++)
    1377                                         {
    1378                                           if (leadmonom(lform)==truncatedcomplex[2*(i+1)-1][l][ideg+1][d][1])
    1379                                             {
    1380                                               M[sizetruncom[2*i][k]+truncatedcomplex[2*i-1][k][j][b][2],sizetruncom[2*(i+1)][l]+truncatedcomplex[2*(i+1)-1][l][ideg+1][d][2]]=leadcoef(lform);
    1381                                               break;
    1382                                             }
    1383                                         }
    1384                                     }
    1385                                 }
    1386                               form=form-lform;
    1387                             }
    1388                         }
    1389                     }
    1390                 }
    1391             }
    1392         }
    1393       truncatedcomplex[2*i]=M;
    1394       truncatedcomplex[2*i-1]=sizetruncom[2*i][size(sizetruncom[2*i])];
    1395     }
    1396   truncatedcomplex[2*i-1]=sizetruncom[2*i][size(sizetruncom[2*i])];
    1397   if (truncatedcomplex[2*i-1]!=0)
    1398     {
    1399       truncatedcomplex[2*i]=matrix(0,truncatedcomplex[2*i-1],1);
    1400     }
    1401   setring R;
    1402   list truncatedcomplex=imap(W,truncatedcomplex);
    1403   list derhamhom=findhomology(truncatedcomplex,le);
    1404   return (derhamhom);
    1405 }
    1406 
    1407 ///////////////////////////////////
    1408 static proc sublist(list L, int m, int n)
    1409 {
    1410   list out;
    1411   int i; int j;
    1412   int count;
    1413   for (i=m; i<=n; i++)
    1414     {
    1415       out[size(out)+1]=list();
    1416       for (j=1; j<=size(L[i]); j++)
    1417         {
    1418           count=count+1;
    1419           out[size(out)][j]=list(L[i][j],count);
    1420         }
    1421     }
    1422   list o=list(out,count);
    1423   return(o);
    1424 }
    1425 
    1426 //////////////////////////////////////////////////////////////////////////
    1427 static proc mysubset(intvec L, intvec M)
    1428 {
    1429   int i;
    1430   int j=1;
    1431   list position=(M[size(M)],(-1)^(size(L)));
    1432   for (i=1; i<=size(L); i++)
    1433     {
    1434       if (L[i]!=M[j])
    1435         {
    1436           if (L[i]!=M[j+1] or j!=i)
    1437             {
    1438               return (L[i],0);
    1439             }
    1440           else
    1441             {
    1442               position=(M[i],(-1)^(i-1));
    1443               j=j+i;
    1444             }
    1445         }
    1446       j=j+1;
    1447     }
    1448   return (position);
    1449 }
    1450 
    1451 
    1452 
    1453 ////////////////////////////////////////////////////////////////////////////
    1454 
    1455 proc globalbfun(list L)
    1456 {
    1457   int i; int j;
    1458   def W=basering;
    1459   int n=nvars(W) div 2;
    1460   list G0;
    1461   ideal I;
    1462   for (j=1; j<=size(L); j++)
    1463     {
    1464       G0[j]=list();
    1465       for (i=1; i<=ncols(L[j][1]); i++)
    1466         {
    1467           G0[j][i]=I;
    1468         }
    1469     }
    1470   list out;
    1471   for (j=1; j<=size(L); j++)
    1472     {
    1473       for (i=1; i<=nrows(L[j][1]); i++)
    1474         {
    1475           out=Vddeg(submat(L[j][1],i,(1..ncols(L[j][1]))),n,L[j][2],1);
    1476           G0[j][out[2]][size(G0[j][out[2]])+1]=(out[1]);
    1477         }
    1478     }
    1479   list Data=ringlist(W);
    1480   for (i=1; i<=n; i++)
    1481     {
    1482       Data[2][2*n+i]=Data[2][i];
    1483       Data[2][3*n+i]=Data[2][n+i];
    1484       Data[2][i]="v("+string(i)+")";
    1485       Data[2][n+i]="w("+string(i)+")";
    1486     }
    1487   Data[3][1][1]="M";
    1488   intvec mord=(0:16*n^2);
    1489   mord[1..2*n]=(1:2*n);
    1490   mord[6*n+1..8*n]=(1:2*n);
    1491   for (i=0; i<=2*n-2; i++)
    1492     {
    1493       mord[(3+i)*4*n-i]=-1;
    1494       mord[(2*n+2+i)*4*n-2*n-i]=-1;
    1495     }
    1496   Data[3][1][2]=mord;//ordering mh?????????
    1497   matrix Ones=UpOneMatrix(4*n);
    1498   Data[5]=Ones;
    1499   matrix con[2*n][2*n];
    1500   Data[6]=transpose(concat(con,transpose(concat(con,Data[6]))));
    1501 
    1502   def Wuv=ring(Data);
    1503   setring Wuv;
    1504   list G0=imap(W,G0); list G3; poly lterm;intvec lexp;
    1505   list G1;  list G2; intvec e; intvec f; int  kapp; int k; int l; poly h; ideal I;
    1506   for (l=1; l<=size(G0); l++)
    1507     {
    1508       G1[l]=list();  G2[l]=list(); G3[l]=list();
    1509       for (i=1; i<=size(G0[l]); i++)
    1510         {
    1511           for (j=1; j<=ncols(G0[l][i]);j++)
    1512             {
    1513               G0[l][i][j]=mhom(G0[l][i][j]);
    1514             }
    1515           for (j=1; j<=nvars(Wuv) div 4; j++)
    1516             {
    1517               G0[l][i][size(G0[l][i])+1]=1-v(j)*w(j);
    1518             }
    1519           G1[l][i]=std(G0[l][i]);
    1520           G2[l][i]=I;
    1521           G3[l][i]=list();
    1522           for (j=1; j<=ncols(G1[l][i]); j++)
    1523             {
    1524               e=leadexp(G1[l][i][j]);
    1525               f=e[1..2*n];
    1526               if (f==intvec(0:(2*n)))
    1527                 {
    1528                   for (k=1; k<=n; k++)
    1529                     {
    1530                       kapp=-e[2*n+k]+e[3*n+k];
    1531                       if (kapp>0)
    1532                         {
    1533                           G1[l][i][j]=(x(k)^kapp)*G1[l][i][j];
    1534                         }
    1535                       if (kapp<0)
    1536                         {
    1537                           G1[l][i][j]=(D(k)^(-kapp))*G1[l][i][j];
    1538                         }
    1539                     }
    1540                   G2[l][i][size(G2[l][i])+1]=G1[l][i][j];
    1541                   G3[l][i][size(G3[l][i])+1]=list();
    1542                   while (G1[l][i][j]!=0)
    1543                     {
    1544                       lterm=lead(G1[l][i][j]);
    1545                       G1[l][i][j]=G1[l][i][j]-lterm;
    1546                       lexp=leadexp(lterm);
    1547                       lexp=lexp[2*n+1..3*n];
    1548                       G3[l][i][size(G3[l][i])][size(G3[l][i][size(G3[l][i])])+1]=list(lexp,leadcoef(lterm));
    1549                     }
    1550 
    1551                 }
    1552             }
    1553         }
    1554     }
    1555   ring r=0,(s(1..n)),dp;
    1556   ideal I;
    1557   map G3forr=Wuv,I;
    1558   list G3=G3forr(G3);
    1559   poly fs;
    1560   poly gs;
    1561   int a;
    1562   list G4;
    1563   for (l=1; l<=size(G3); l++)
    1564     {
    1565       G4[l]=list();
    1566       for (i=1; i<=size(G3[l]);i++)
    1567         {
    1568           G4[l][i]=I;
    1569           for (j=1; j<=size(G3[l][i]); j++)
    1570             {
    1571               fs=0;
    1572               for (k=1; k<=size(G3[l][i][j]); k++)
    1573                 {
    1574                   gs=1;
    1575                   for (a=1; a<=n; a++)
    1576                     {
    1577                       if (G3[l][i][j][k][1][a]!=0)
    1578                         {
    1579                           gs=gs*permutevar(list(G3[l][i][j][k][1][a]),a);
    1580                         }
    1581                     }
    1582                   gs=gs*G3[l][i][j][k][2];
    1583                   fs=fs+gs;
    1584                 }
    1585               G4[l][i]=G4[l][i],fs;
    1586             }
    1587         }
    1588     }
    1589   if (n==1)
    1590     {
    1591       ring rnew=0,t,dp;
    1592     }
    1593   else
    1594     {
    1595       ring rnew=0,(t,s(2..n)),dp;
    1596     }
    1597   ideal Iformap;
    1598   Iformap[1]=t;
    1599    poly forel=1;
    1600    for (i=2; i<=n; i++)
    1601      {
    1602        Iformap[1]=Iformap[1]-s(i);
    1603        Iformap[i]=s(i);
    1604        forel=forel*s(i);
    1605      }
    1606    map rtornew=r,Iformap;
    1607    list G4=rtornew(G4);
    1608    list getintvecs=fetch(W,L);
    1609    ideal J;
    1610    option(redSB);
    1611    for (l=1; l<=size(G4); l++)
    1612      {
    1613        J=1;
    1614        for (i=1; i<=size(G4[l]); i++)
    1615          {
    1616            G4[l][i]=eliminate(G4[l][i],forel);
    1617            G4[l][i]=subst(G4[l][i],t,t-getintvecs[l][2][i]);
    1618            J=intersect(J,G4[l][i]);
    1619          }
    1620        G4[l]=poly(std(J)[1]);
    1621      }
    1622    list minmax=minmaxintroot(G4);//besser factorize nehmen
    1623    // Fall: keine Nullstelle muss noch weiter beruecksichtigt werden
    1624   return(minmax);
    1625 }
    1626 
    1627 
    1628 
    1629 
    1630 //////////////////////////////////////////////////////////////////////////
    1631 
    1632 proc minmaxintroot(list L);
    1633 {
    1634   int i; int j; int k; int l; int sa; int s; number d; poly f; poly rest; list a0; list possk; list alldiv; intvec e;
    1635   possk[1]=list();
    1636   for (i=1; i<=size(L); i++)
    1637     {
    1638       d=1;
    1639       f=L[i];
    1640       while (f!=0)
    1641         {
    1642           rest=lead(f);
    1643           d=d*denominator(leadcoef(rest));
    1644           f=f-rest;
    1645         }
    1646       e=leadexp(rest);
    1647       if (e[1]!=0)
    1648         {
    1649           rest=rest/(t^(e[1]));
    1650           possk[1][size(possk[1])+1]=i;
    1651         }
    1652       a0[i]=int(absValue(d*rest));
    1653     }
    1654   int m=Max(a0);
    1655   for (i=2; i<=m+1; i++)
    1656     {
    1657       possk[i]=list();
    1658     }
    1659   list allprimefac;
    1660   for (i=1; i<=size(L); i++)
    1661     {
    1662       allprimefac=primefactors(a0[i]);
    1663       alldiv=1;
    1664       possk[2][size(possk[2])+1]=i;
    1665 
    1666       for (j=1; j<=size(allprimefac[1]); j++)
    1667         {
    1668           s=size(alldiv);
    1669           for (k=1; k<=s; k++)
    1670             {
    1671               for (l=1; l<=allprimefac[2][j]; l++)
    1672                 {
    1673                   alldiv[size(alldiv)+1]=alldiv[k]*allprimefac[1][j]^l;
    1674                   possk[alldiv[size(alldiv)]+1][size(possk[alldiv[size(alldiv)]+1])+1]=i;
    1675                 }
    1676             }
    1677         }
    1678     }
    1679   int mink;
    1680   int maxk;
    1681   int indi;
    1682   for (i=m+1; i>=1; i--)
    1683     {
    1684       if (size(possk[i])!=0)
    1685         {
    1686           for (j=1; j<=size(possk[i]); j++)
    1687             {
    1688               if (subst(L[possk[i][j]],t,(i-1))==0)
    1689                 {
    1690                   maxk=i-1;
    1691                   indi=1;
    1692                   break;
    1693                 }
    1694             }
    1695         }
    1696       if (maxk!=0)
    1697         {
    1698           break;
    1699         }
    1700     }
    1701   int indi2;
    1702   for (i=m+1; i>=1; i--)
    1703     {
    1704       if (size(possk[i])!=0)
    1705         {
    1706           for (j=1; j<=size(possk[i]); j++)
    1707             {
    1708               if (subst(L[possk[i][j]],t,-(i-1))==0)
    1709                 {
    1710                   mink=-i+1;
    1711                   indi2=1;
    1712                   break;
    1713                 }
    1714             }
    1715         }
    1716       if (mink!=0)
    1717         {
    1718           break;
    1719         }
    1720     }
    1721   list mima=mink,maxk;
    1722   if (indi==0)
    1723     {
    1724       if (indi2==0)
    1725         {
    1726           mima=list();//es gibt keine ganzzahlige NS
    1727         }
    1728       else
    1729         {
    1730           mima[2]=mima[1];
    1731         }
    1732     }
    1733   else
    1734     {
    1735       if (indi2==0)
    1736         {
    1737           mima[1]=mima[2];
    1738         }
    1739     }
    1740   return (mima);
    1741 }
    1742 
    1743 ///////////////////////////////////////////////////////
    1744 
    1745 
    1746 proc findhomology(list L, int le)
    1747 {
    1748   int li;
    1749   matrix M; matrix N;
    1750   matrix N1;
    1751   matrix lift1;
    1752   list out;
    1753   int i;
    1754   option (redSB);
    1755   for (i=2; i<=size(L); i=i+2)
    1756     {
    1757       if (L[i-1]==0)
    1758         {
    1759           li=0;
    1760           out[i div 2]=0;
    1761         }
    1762       else
    1763         {
    1764 
    1765           if (li==0)
    1766             {
    1767 
    1768 
    1769               li=L[i-1];
    1770               N1=transpose(syz(transpose(L[i])));
    1771               out[i div 2]=matrix(transpose(syz(transpose(N1))));
    1772               out[i div 2]=transpose(matrix(std(transpose(out[i div 2]))));
    1773 
    1774             }
    1775 
    1776           else
    1777             {
    1778 
    1779 
    1780               li=L[i-1];
    1781               N1=transpose(syz(transpose(L[i])));
    1782               N=transpose(syz(transpose(N1)));
    1783               lift1=matrixlift(N1,L[i-2]);
    1784               out[i div 2]=transpose(concat(transpose(lift1),transpose(N)));
    1785               out[i div 2]=transpose(matrix(std(transpose(out[i div 2]))));
    1786             }
    1787         }
    1788       if (out[i div 2]!=matrix(0,1,ncols(out[i div 2])))
    1789         {
    1790           out[i div 2]=ncols(out[i div 2])-nrows(out[i div 2]);
    1791         }
    1792       else
    1793         {
    1794           out[i div 2]=ncols(out[i div 2]);
    1795         }
    1796     }
    1797   if (size(out)>le)
    1798     {
    1799       out=delete(out,1);
    1800     }
    1801   return(out);
    1802 }
    1803 
    1804 
    1805 
    1806 
    1807 /////////////////////////////////////////////////////////////////////
    1808 
    1809 static proc mhom(poly f)
    1810 {
    1811   poly g;
    1812   poly l;
    1813   poly add;
    1814   intvec e;
    1815   list minint;
    1816   list remf;
    1817   int i;
    1818   int j;
    1819   int n=nvars(basering) div 4;
    1820   if (f==0)
    1821     {
    1822       return(f);
    1823     }
    1824   while (f!=0)
    1825     {
    1826       l=lead(f);
    1827       e=leadexp(l);
    1828       remf[size(remf)+1]=list();
    1829       remf[size(remf)][1]=l;
    1830       for (i=1; i<=n; i++)
    1831         {
    1832           remf[size(remf)][i+1]=-e[2*n+i]+e[3*n+i];
    1833           if (size(minint)<i)
    1834             {
    1835               minint[i]=list();
    1836             }
    1837           minint[i][size(minint[i])+1]=-e[2*n+i]+e[3*n+i];
    1838         }
    1839       f=f-l;
    1840     }
    1841   for (i=1; i<=n; i++)
    1842     {
    1843       minint[i]=Min(minint[i]);
    1844     }
    1845   for (i=1; i<=size(remf); i++)
    1846     {
    1847       add=remf[i][1];
    1848       for (j=1; j<=n; j++)
    1849         {
    1850           add=v(j)^(remf[i][j+1]-minint[j])*add;
    1851         }
    1852       g=g+add;
    1853     }
    1854   return (g);
    1855 }
    1856 
    1857 
    1858 
    1859 
    1860 //////////////////////////////////////////////////////////////////////////
    1861 
    1862 static proc permutevar(list L,int n)
    1863 {
    1864   if (typeof(L[1])=="intvec")
    1865     {
    1866       intvec v=L[1];
    1867     }
    1868   else
    1869     {
    1870       intvec v=(1:L[1]),(0:L[1]);
    1871     }
    1872   int i;int k; int indi=0;
    1873   int j;
    1874   int s=size(v);
    1875   poly e;
    1876   intvec fore;
    1877   for (i=2; i<=size(v); i=i+2)
    1878     {
    1879       if (v[i]!=0)
    1880         {
    1881           j=i+1;
    1882           while (v[j]!=0)
    1883             {
    1884               j=j+1;
    1885             }
    1886           v[i]=0;
    1887           v[j]=1;
    1888           fore=0;
    1889           indi=0;
    1890           for (k=1; k<=size(v); k++)
    1891             {
    1892               if (k!=i and k!=j)
    1893                 {
    1894                   if (indi==0)
    1895                     {
    1896                       indi=1;
    1897                       fore[1]=v[k];
    1898                     }
    1899                   else
    1900                     {
    1901                       fore[size(fore)+1]=v[k];
    1902                     }
    1903                 }
    1904             }
    1905           e=e-(j-i)*permutevar(list(fore),n);
    1906         }
    1907     }
    1908   e=e+s(n)^(size(v) div 2);
    1909   return (e);
    1910 }
    1911 
    1912 ///////////////////////////////////////////////////////////////////////////////
    1913 static proc max(int i,int j)
    1914 {
    1915   if(i>j){return(i);}
    1916   return(j);
    1917 }
    1918 
    1919 ////////////////////////////////////////////////////////////////////////////////////
    1920 version="$Id$";
     2version="version derham.lib 4.0.0.0 Jun_2013 ";
    19213category="Noncommutative";
    19224info="
     
    195436LIB "qhmoduli.lib";
    195537LIB "general.lib";
    1956 LIB "dmod.lib";
     38//LIB "dmod.lib";
    195739LIB "bfun.lib";
    195840LIB "dmodapp.lib";
     
    2643725             annihilators we computed for L[i]  ->this will ensure  we can compute
    2644726             the maps of the MV complex*/
    2645           minroot=minIntRoot(findminintroot);
     727          minroot=Derham::minIntRoot(findminintroot);
    2646728          for (j=1; j<=i; j++)
    2647729            {
     
    2729811                   s-parametric annihilators we computed for L[i]  ->this will
    2730812                   ensure  we can compute the maps of the MV complex*/
    2731               minroot=minIntRoot(findminintroot);
     813              minroot=Derham::minIntRoot(findminintroot);
    2732814              for (j=1; j<=i; j++)
    2733815                {
     
    29691051    }
    29701052  def B=basering;
    2971   int n=nvars(B) div 2 + 1;//+1 müsste stimmen! bitte kontrollieren!
     1053  int n=nvars(B) div 2 + 1;//+1 muesste stimmen! bitte kontrollieren!
    29721054  int d=nvars(B) div 2;
    29731055  int r=size(L) div 2;
     
    32661348  matrix zerom;
    32671349  list rofA;
    3268   for (i=3; i<=n+3; i++)////////////////////////////////////////////////////////////////////////////n und si müssen noch definiert werden
     1350  for (i=3; i<=n+3; i++)////////////////////////////////////////////////////////////////////////////n und si muessen noch definiert werden
    32691351    {
    32701352      if (size(RES)>=i)
     
    57963878                component of submat(L[j][1],i,(1..ncols(L[j][1])))*/
    57973879              i1=(1..ncols(L[j][1]));
    5798               out=VdDegTilde(submat(L[j][1],i,i1),n,intvec(0:size(L[j][2])),1);//hier könnte es evtl noch einen Fehler geben!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3880              out=VdDegTilde(submat(L[j][1],i,i1),n,intvec(0:size(L[j][2])),1);//hier koennte es evtl noch einen Fehler geben!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    57993881              f=L[j][1][i,out[2]];
    58003882              if (out[2]==1)
     
    61484230   for (i=1; i<=size(G4); i++)
    61494231     {
    6150        minmax[i]=minIntRoot(G4[i],1);
     4232       minmax[i]=Derham::minIntRoot(G4[i],1);
    61514233       if (size(minmax[i])!=0)
    61524234         {
     
    61964278          else
    61974279            {
    6198               out[i div 2]=0;////achtung geändert??????????????????????????????????????????????????!!!!!!!!!!!!!!!!!!!!!!!!!war mal out[i-1]
     4280              out[i div 2]=0;////achtung geaendert??????????????????????????????????????????????????!!!!!!!!!!!!!!!!!!!!!!!!!war mal out[i-1]
    61994281              im=L[i-1];
    62004282            }
     
    65344616              and submat(Fnew,l,intvec(1..ncols(Fnew)))!=zeromat)
    65354617            {
    6536               //print("Reduzierung des V_d-Grades nötig");
     4618              //print("Reduzierung des V_d-Grades noetig");
    65374619              /*We need to reduce the V_d-degree. First we homogenize the
    65384620                lth row of Fnew*/
Note: See TracChangeset for help on using the changeset viewer.