Changeset 2815e8 in git for Singular/LIB/multigrading.lib


Ignore:
Timestamp:
Mar 17, 2011, 2:22:01 PM (13 years ago)
Author:
Hans Schoenemann <hannes@…>
Branches:
(u'fieker-DuVal', '117eb8c30fc9e991c4decca4832b1d19036c4c65')(u'spielwiese', 'd25190065115c859833252500a64cfb7b11e3a50')
Children:
94312b5b33172e44faad48ff7a9912f1cd6aa31f
Parents:
fec3bde6a8e1670c13a5b3176f479ac2b3da7399
Message:
fix some errors of multigrading.lib

git-svn-id: file:///usr/local/Singular/svn/trunk@14007 2c84dea3-7e68-4137-9b89-c4e89433aadc
File:
1 edited

Legend:

Unmodified
Added
Removed
  • Singular/LIB/multigrading.lib

    rfec3bde r2815e8  
    1616OVERVIEW: This library allows one to virtually add multigradings to Singular.
    1717For more see http://code.google.com/p/convex-singular/wiki/Multigrading
    18 For theoretical references see: 
     18For theoretical references see:
    1919E. Miller, B. Sturmfels: 'Combinatorial Commutative Algebra' and
    2020M. Kreuzer, L. Robbiano: 'Computational Commutative Algebra'.
     
    3232createGroup(S,L);          create a group generated by S, with relations L
    3333createQuotientGroup(L);    create a group generated by the unit matrix whith relations L
    34 createTorsionFreeGroup(S); create a group generated by S which is torsionfree 
     34createTorsionFreeGroup(S); create a group generated by S which is torsionfree
    3535printGroup(G);             print a group
    3636
     
    5454intersectLattices(A,B);   compute an integral basis for the intersection of the lattices A and B.
    5555isIntegralSurjective(P);  test whether the map P of lattices is surjective.
    56 isPrimitiveSublattice(A); test whether A generates a primitive sublattice. 
     56isPrimitiveSublattice(A); test whether A generates a primitive sublattice.
    5757intInverse(A);            compute the integral inverse matrix of the intmat A
    5858intAdjoint(A,i,j);        delete row i and column j of the intmat A.
     
    7070isPositive();             test whether the current multigrading is positive
    7171isZeroElement(p);         test whether p has zero multidegree
    72 areZeroElements(M);       test whether an integer matrix M considered as a collection of columns has zero multidegree   
     72areZeroElements(M);       test whether an integer matrix M considered as a collection of columns has zero multidegree
    7373isHomogeneous(a);         test whether 'a' is multigraded-homogeneous
    7474
     
    7979mDegModulo(I,J);          compute the multigraded 'modulo' module of I and J
    8080mDegResolution(M,l[,m]);  compute the multigraded resolution of M
    81 mDegTensor(m,n);          compute the tensor product of multigraded modules m,n 
     81mDegTensor(m,n);          compute the tensor product of multigraded modules m,n
    8282mDegTor(i,m,n);           compute the Tor_i(m,n) for multigraded modules m,n
    8383
     
    128128proc createGradedRingHomomorphism(def src, ideal Im, def A)
    129129"USAGE: createGradedRingHomomorphism(R, f, A); ring R, ideal f, group homomorphism A
    130 PURPOSE: create a multigraded group ring homomorphism defined by 
     130PURPOSE: create a multigraded group ring homomorphism defined by
    131131a ring map from R to the current ring, given by generators images f
    132132and a group homomorphism A between grading groups
     
    137137  string isGRH = "isGRH";
    138138
    139   if( !isGradedRingHomomorphism(def src, ideal Im, def A) )
     139  if( !isGradedRingHomomorphism(src, Im, A) )
    140140  {
    141141    ERROR("Input data is wrong");
    142142  }
    143  
     143
    144144  list h;
    145145  h[3] = A;
    146  
     146
    147147//  map f = src, Im;
    148148  h[2] = Im; // f?
    149149  h[1] = src;
    150  
     150
    151151  attrib(h, isGRH, (1==1)); // mark it "a graded ring homomorphism"
    152152
     
    165165proc isGradedRingHomomorphism(def src, ideal Im, def A)
    166166"USAGE: isGradedRingHomomorphism(R, f, A); ring R, ideal f, group homomorphism A
    167 PURPOSE: test a multigraded group ring homomorphism defined by 
     167PURPOSE: test a multigraded group ring homomorphism defined by
    168168a ring map from R to the current ring, given by generators images f
    169169and a group homomorphism A between grading groups
     
    176176  intmat result_degs = mDeg(Im);
    177177  print(result_degs);
    178    
     178
    179179  setring src;
    180    
     180
    181181  intmat input_degs = mDeg(maxideal(1));
    182182  print(input_degs);
     
    184184  def image_degs = A * input_degs;
    185185  print( image_degs );
    186  
     186
    187187  def df = image_degs - result_degs;
    188188  print(df);
    189  
     189
    190190  setring dst;
    191    
     191
    192192  return (areZeroElements( df ));
    193193}
     
    195195{
    196196  "EXAMPLE:"; echo=2;
    197  
     197
    198198  ring r = 0, (x, y, z), dp;
    199199  intmat S1[3][3] =
     
    203203  intmat L1[3][1] =
    204204    0,
    205     0, 
     205    0,
    206206    0;
    207    
     207
    208208  def G1 = createGroup(S1, L1); // (S1 + L1)/L1
    209209  printGroup(G1);
    210    
     210
    211211  setBaseMultigrading(S1, L1); // to change...
    212212
     
    221221  def G2 = createGroup(S2, L2);
    222222  printGroup(G2);
    223    
     223
    224224  setBaseMultigrading(S2, L2); // to change...
    225225
     
    229229      1, 0, 0,
    230230      3, 2, -6;
    231    
     231
    232232  // graded ring homomorphism is given by (compatible):
    233233  print(F);
     
    238238
    239239  print(h);
    240  
     240
    241241  // not a homo..
    242242  intmat B[nrows(L2)][nrows(L1)] =
     
    247247  isGradedRingHomomorphism(r, ideal(F), B);
    248248  def hh = createGradedRingHomomorphism(r, ideal(F), B);
    249  
     249
    250250  if( defined(hh) ) { ERROR("That input was not valid"); }
    251251}
     
    254254proc createQuotientGroup(intmat L)
    255255"
    256 L - relations 
    257 TODO: bad name 
     256L - relations
     257TODO: bad name
    258258"
    259259{
     
    267267"
    268268S - generators
    269 TODO: bad name 
     269TODO: bad name
    270270"
    271271{
     
    279279proc createGroup(intmat S, intmat L)
    280280"USAGE: createGroup(S, L); S, L are integer matrices
    281 PURPOSE: create the group of the form (S+L)/L, i.e. 
     281PURPOSE: create the group of the form (S+L)/L, i.e.
    282282S specifies generators, L specifies relations.
    283283RETURN: group
     
    289289  string attrGroupSNF = "smith";
    290290
    291    
     291
    292292/*
    293293  if( size(#) > 0 )
     
    305305      } else { ERROR("Wrong optional argument: 2"); }
    306306    }
    307   } 
    308   if( !defined(S) ) 
     307  }
     308  if( !defined(S) )
    309309  {}
    310 */   
     310*/
    311311
    312312  if( nrows(L) != nrows(S) )
     
    314314    ERROR("Incompatible matrices!");
    315315  }
    316  
     316
    317317  def H = attrib(L, attrGroupHNF);
    318318  if( !defined(H) || typeof(H) != "intmat")
     
    320320    attrib(L, attrGroupHNF, hermiteNormalForm(L));
    321321  } else { kill H; }
    322  
     322
    323323  def HH = attrib(L, attrGroupSNF);
    324324  if( !defined(HH) || typeof(HH) != "intmat")
     
    337337{
    338338  "EXAMPLE:"; echo=2;
    339    
     339
    340340  intmat S[3][3] =
    341341    1, 0, 0,
     
    345345  intmat L[3][2] =
    346346    1, 1,
    347     1, 3, 
     347    1, 3,
    348348    1, 5;
    349    
     349
    350350  def G = createGroup(S, L); // (S+L)/L
    351351
    352352  printGroup(G);
    353    
     353
    354354  kill S, L, G;
    355355
     
    366366
    367367  printGroup(G);
    368    
     368
    369369  kill S, L, G;
    370    
     370
    371371  // ----------- extreme case ------------ //
    372372  intmat S[1][3] =
     
    396396  "Relations: ";
    397397  print(G[2]);
    398  
     398
    399399//  attrib(G[2]);
    400400}
     
    402402{
    403403  "EXAMPLE:"; echo=2;
    404    
     404
    405405}
    406406
     
    418418{
    419419  "EXAMPLE:"; echo=2;
    420    
     420
    421421}
    422422
     
    426426{
    427427  string isGroup = "isGroup";
    428  
     428
    429429  // valid?
    430430  if( typeof(G) != "list" ){ return(0); }
     
    435435
    436436//  if( !defined(a) ) { return(0); }
    437 //  if( typeof(a) != "int" ) { return(0); } 
     437//  if( typeof(a) != "int" ) { return(0); }
    438438//  if( !a ){ return(0); }
    439439
    440    
     440
    441441  if( size(G) != 2 ){ return(0); }
    442442  if( typeof(G[1]) != "intmat" ){ return(0); }
    443443  if( typeof(G[2]) != "intmat" ){ return(0); }
    444444  if( nrows(G[1]) != nrows(G[2]) ){ return(0); }
    445    
     445
    446446  return(1==1);
    447447}
     
    460460  string attrMgrad   = "mgrad";
    461461  string attrGradingGroup = "gradingGroup";
    462  
     462
    463463  if( size(#) > 0 )
    464464  {
     
    466466    {
    467467      def L = createGroup(M, #[1]);
    468     } 
     468    }
    469469
    470470    if( isGroup(#[1]) )
     
    477477      }
    478478    }
    479   } 
     479  }
    480480  else
    481481  {
     
    484484
    485485  if( !defined(L) ){ ERROR("Wrong arguments: no group given?"); }
    486    
    487 
    488  
    489   attrib(basering, attrMgrad, M); 
    490   attrib(basering, attrGradingGroup, L); 
     486
     487
     488
     489  attrib(basering, attrMgrad, M);
     490  attrib(basering, attrGradingGroup, L);
    491491
    492492  ideal Q = ideal(basering);
     
    495495    "Warning: your quotient ideal is not homogenous (multigrading was set anyway)!";
    496496  }
    497  
    498    
     497
     498
    499499
    500500
     
    503503{
    504504  "EXAMPLE:"; echo=2;
    505    
     505
    506506  ring R = 0, (x, y, z), dp;
    507507
     
    515515  intmat L[3][2] =
    516516    1, 1,
    517     1, 3, 
     517    1, 3,
    518518    1, 5;
    519    
     519
    520520  // attaches M & L to R (==basering):
    521521  setBaseMultigrading(M, L); // Grading: Z^3/L
     
    523523  // Weights are accessible via "getVariableWeights()":
    524524  getVariableWeights();
    525  
    526   // Test all possible usages: 
     525
     526  // Test all possible usages:
    527527  (getVariableWeights() == M) && (getVariableWeights(R) == M) && (getVariableWeights(basering) == M);
    528528
    529529  // Grading group is accessible via "getLattice()":
    530530  getLattice();
    531  
    532   // Test all possible usages: 
     531
     532  // Test all possible usages:
    533533  (getLattice() == L) && (getLattice(R) == L) && (getLattice(basering) == L);
    534534
     
    536536  getLattice("hermite");
    537537
    538   // Test all possible usages: 
     538  // Test all possible usages:
    539539  intmat H = hermiteNormalForm(L);
    540540  (getLattice("hermite") == H) && (getLattice(R, "hermite") == H) && (getLattice(basering, "hermite") == H);
     
    553553    0,
    554554    2;
    555    
     555
    556556  // attaches M & L to R (==basering):
    557557  setBaseMultigrading(M, L); // Grading: Z + (Z/2Z)
     
    573573  intmat L[1][1] =
    574574    0;
    575    
     575
    576576  // attaches M & L to R (==basering):
    577577  setBaseMultigrading(M); // Grading: Z^3
     
    614614  def M = attrib(R, attrMgrad);
    615615  if( typeof(M) == "intmat"){ return (M); }
    616   ERROR( "Sorry no multigrading matrix!" ); 
     616  ERROR( "Sorry no multigrading matrix!" );
    617617}
    618618example
    619619{
    620620  "EXAMPLE:"; echo=2;
    621    
     621
    622622  ring R = 0, (x, y, z), dp;
    623623
     
    631631  intmat L[3][2] =
    632632    1, 1,
    633     1, 3, 
     633    1, 3,
    634634    1, 5;
    635    
     635
    636636  // attaches M & L to R (==basering):
    637637  setBaseMultigrading(M, L); // Grading: Z^3/L
     
    653653    0,
    654654    2;
    655    
     655
    656656  // attaches M & L to R (==basering):
    657657  setBaseMultigrading(M, L); // Grading: Z + (Z/2Z)
     
    671671  intmat L[1][1] =
    672672    0;
    673    
     673
    674674  // attaches M & L to R (==basering):
    675675  setBaseMultigrading(M); // Grading: Z^3
     
    697697      def R = #[i];
    698698      i++;
    699     } 
    700   } 
     699    }
     700  }
    701701
    702702  if( !defined(R) )
     
    706706
    707707  def G = attrib(R, attrGradingGroup);
    708  
     708
    709709  if( !isGroup(G) )
    710710  {
    711711    ERROR("Sorry no grading group!");
    712   }   
    713 
    714   return(G); 
     712  }
     713
     714  return(G);
    715715}
    716716example
    717717{
    718718  "EXAMPLE:"; echo=2;
    719    
     719
    720720  ring R = 0, (x, y, z), dp;
    721721
     
    729729  intmat L[3][2] =
    730730    1, 1,
    731     1, 3, 
     731    1, 3,
    732732    1, 5;
    733    
     733
    734734  // attaches M & L to R (==basering):
    735735  setBaseMultigrading(M, L); // Grading: Z^3/L
     
    738738
    739739  printGroup( G );
    740  
     740
    741741  G[1] == M; G[2] == L;
    742742
     
    754754    0,
    755755    2;
    756    
     756
    757757  // attaches M & L to R (==basering):
    758758  setBaseMultigrading(M, L); // Grading: Z + (Z/2Z)
     
    761761
    762762  printGroup( G );
    763  
     763
    764764  G[1] == M; G[2] == L;
    765765
     
    774774  intmat L[1][1] =
    775775    0;
    776    
     776
    777777  // attaches M & L to R (==basering):
    778778  setBaseMultigrading(M); // Grading: Z^3
     
    781781
    782782  printGroup( G );
    783  
     783
    784784  G[1] == M; G[2] == L;
    785785
     
    792792"USAGE: getLattice([R[,opt]])
    793793PURPOSE: get associated grading group matrix, i.e. generators (cols) of the grading group
    794 RETURN: intmat, the grading group matrix, or 
     794RETURN: intmat, the grading group matrix, or
    795795its hermite normal form if an optional argument (\"hermiteNormalForm\") is given or
    796796smith normal form if an optional argument (\"smith\") is given
     
    804804    {
    805805      i++;
    806     } 
    807   } 
     806    }
     807  }
    808808
    809809  string attrGradingGroupHNF = "hermite";
     
    811811
    812812  def G = getGradingGroup(#);
    813  
    814 //  printGroup(G); 
     813
     814//  printGroup(G);
    815815
    816816
     
    824824      def M = attrib(T, attrGradingGroupHNF);
    825825      if( (!defined(M)) or (typeof(M) != "intmat") )
    826       { 
    827         M = hermiteNormalForm(T);
     826      {
     827        M = hermiteNormalForm(T);
    828828      }
    829829      return (M);
     
    834834      def M = attrib(T, attrGradingGroupSNF);
    835835      if( (!defined(M)) or (typeof(M) != "intmat") )
    836       { 
    837         M = smithNormalForm(T);
     836      {
     837        M = smithNormalForm(T);
    838838      }
    839839      return (M);
     
    841841  }
    842842
    843   return(T); 
     843  return(T);
    844844}
    845845example
    846846{
    847847  "EXAMPLE:"; echo=2;
    848    
     848
    849849  ring R = 0, (x, y, z), dp;
    850850
     
    858858  intmat L[3][2] =
    859859    1, 1,
    860     1, 3, 
     860    1, 3,
    861861    1, 5;
    862    
     862
    863863  // attaches M & L to R (==basering):
    864864  setBaseMultigrading(M, L); // Grading: Z^3/L
     
    883883    0,
    884884    2;
    885    
     885
    886886  // attaches M & L to R (==basering):
    887887  setBaseMultigrading(M, L); // Grading: Z + (Z/2Z)
     
    904904  intmat L[1][1] =
    905905    0;
    906    
     906
    907907  // attaches M & L to R (==basering):
    908908  setBaseMultigrading(M); // Grading: Z^3
     
    920920"
    921921{
    922   if( typeof(m) == "ideal" ) 
     922  if( typeof(m) == "ideal" )
    923923  {
    924924    return (m[i]);
     
    928928  {
    929929    def v = getModuleGrading(m);
    930    
     930
    931931    return ( setModuleGrading(m[i],v) );
    932932  }
     
    948948
    949949  def V = attrib(m, attrModuleGrading);
    950  
     950
    951951  if( typeof(V) != "intmat" )
    952952  {
     
    957957      return (VV);
    958958    }
    959      
     959
    960960    ERROR("Sorry: vector or module need module-grading-matrix! See 'getModuleGrading'.");
    961961  }
     
    970970    ERROR("Sorry wrong width of V: " + string(ncols(V)));
    971971  }
    972    
     972
    973973  return (V);
    974974}
     
    976976{
    977977  "EXAMPLE:"; echo=2;
    978    
     978
    979979   ring R = 0, (x,y), dp;
    980980   intmat M[2][2]=
     
    984984     1,  2,  3,  4, 0,
    985985     0, 10, 20, 30, 1;
    986    
     986
    987987   setBaseMultigrading(M, T);
    988    
     988
    989989   ideal I = x, y, xy^5;
    990990   isHomogeneous(I);
    991  
     991
    992992   intmat V = mDeg(I); print(V);
    993993
    994994   module S = syz(I); print(S);
    995    
     995
    996996   S = setModuleGrading(S, V);
    997997
    998998   getModuleGrading(S) == V;
    999    
     999
    10001000   vector v = getGradedGenerator(S, 1);
    10011001   getModuleGrading(v) == V;
    10021002   isHomogeneous(v);
    1003    print( mDeg(v) );   
    1004    
     1003   print( mDeg(v) );
     1004
    10051005   isHomogeneous(S);
    10061006   print( mDeg(S) );
     
    10221022  if(nrows(G) != nrows(R)){ ERROR("Incompatible gradings.");}
    10231023  if(ncols(G) < nrows(m)){ ERROR("Multigrading does not fit to module.");}
    1024  
     1024
    10251025  attrib(m, attrModuleGrading, G);
    10261026  return(m);
     
    10291029{
    10301030  "EXAMPLE:"; echo=2;
    1031    
     1031
    10321032   ring R = 0, (x,y), dp;
    10331033   intmat M[2][2]=
     
    10371037     1,  2,  3,  4, 0,
    10381038     0, 10, 20, 30, 1;
    1039    
     1039
    10401040   setBaseMultigrading(M, T);
    1041    
     1041
    10421042   ideal I = x, y, xy^5;
    10431043   intmat V = mDeg(I);
    1044    
     1044
    10451045   // V == M; modulo T
    10461046   print(V);
    10471047
    10481048   module S = syz(I);
    1049    
     1049
    10501050   S = setModuleGrading(S, V);
    10511051   getModuleGrading(S) == V;
    10521052
    10531053   print(S);
    1054    
     1054
    10551055   vector v = getGradedGenerator(S, 1);
    10561056   getModuleGrading(v) == V;
    10571057
    1058    print( mDeg(v) );   
     1058   print( mDeg(v) );
    10591059
    10601060   isHomogeneous(S);
     
    10891089    for(i = 1; i<=mr; i++)
    10901090    {
    1091       result[((i-1)*nr+1)..(i*nr),((i-1)*nc+1)..(i*nc)] = N[1..nr,1..nc]; 
     1091      result[((i-1)*nr+1)..(i*nr),((i-1)*nc+1)..(i*nc)] = N[1..nr,1..nc];
    10921092    }
    10931093  }
     
    11571157  else
    11581158  {
    1159    
     1159
    11601160    matrix fd[nrows(m)][0];
    11611161    matrix fd2[nrows(l[i+1])][0];
     
    11691169    freedim2 = setModuleGrading(freedim2,getModuleGrading(l[i+1]));
    11701170    freedim3 = setModuleGrading(freedim3, getModuleGrading(l[i]));
    1171    
     1171
    11721172    module mimag = mDegTensor(freedim3, m);
    11731173    //"mimag ok.";
     
    12411241"USAGE: gisGoupHomomorphism(L1,L2,A); L1 and L2 are groups, A is an integer matrix
    12421242PURPOSE: checks whether A defines a group homomorphism phi: L1 --> L2
    1243 RETURN: int, 1 if A defines the homomorphism and 0 otherwise 
     1243RETURN: int, 1 if A defines the homomorphism and 0 otherwise
    12441244EXAMPLE: example isGroupHomomorphism; shows an example
    12451245"
    12461246{
    1247   // TODO: L1, L2 
     1247  // TODO: L1, L2
    12481248  if( (ncols(A) != nrows(L1)) or (nrows(A) != nrows(L2)) )
    12491249  {
     
    12521252
    12531253  intmat im = A * L1;
    1254  
    1255   return  (areZeroElements(im, L2)); 
     1254
     1255  return  (areZeroElements(im, L2));
    12561256}
    12571257example
     
    12641264     0,
    12651265     2;
    1266      
     1266
    12671267   intmat L2[3][2]=
    12681268     0, 0,
     
    12701270     0, 3;
    12711271
    1272   intmat A[3][4] = 
     1272  intmat A[3][4] =
    12731273     1, 2, 3, 0,
    12741274     7, 0, 0, 0,
     
    12781278  isGroupHomomorphism(L1, L2, A);
    12791279
    1280   intmat B[3][4] = 
     1280  intmat B[3][4] =
    12811281     1, 2, 3, 0,
    12821282     7, 0, 0, 0,
    12831283     1, 2, 0, 2;
    1284   print( B ); 
     1284  print( B );
    12851285
    12861286  isGroupHomomorphism(L1, L2, B); // Not a homomorphism!
     
    13081308
    13091309    if(H[j, i]!=0)
    1310     { 
     1310    {
    13111311      d=d*H[j, i];
    13121312    }
     
    13141314
    13151315  if( (d*d)==1 )
    1316   { 
     1316  {
    13171317    return(1==1);
    13181318  }
     
    13301330     1, 2, 3, 4, 0,
    13311331     0,10,20,30, 1;
    1332    
     1332
    13331333   setBaseMultigrading(M,T);
    1334    
     1334
    13351335   // Is the resulting group  free?
    13361336   isTorsionFree();
     
    13401340
    13411341   ring R=0,(x,y,z),dp;
    1342    intmat A[3][3] = 
     1342   intmat A[3][3] =
    13431343     1,0,0,
    13441344     0,1,0,
     
    14301430      }
    14311431      result= result+ list (buf);
    1432      
     1432
    14331433    }
    14341434    return(result);
    1435   } 
     1435  }
    14361436  else
    14371437  {
     
    15221522    }
    15231523  }
    1524   //"here is the size ",size(#); 
     1524  //"here is the size ",size(#);
    15251525  if(size(#) == 0){
    15261526    return(D);
     
    15331533"EXAMPLE: "; echo=2;
    15341534
    1535 intmat A[5][7] = 
     1535intmat A[5][7] =
    153615361,0,1,0,-2,9,-71,
    153715370,-24,248,-32,-96,448,-3496,
     
    15501550
    15511551/******************************************************/
    1552 proc hermiteNormalForm(intmat A, list #) 
     1552proc hermiteNormalForm(intmat A, list #)
    15531553"USAGE: hermiteNormalForm( A );
    1554 PURPOSE: Computes the (lower triangular) Hermite Normal Form 
     1554PURPOSE: Computes the (lower triangular) Hermite Normal Form
    15551555           of the matrix A by column operations.
    15561556RETURN: intmat, the Hermite Normal Form of A
     
    15581558"
    15591559{
    1560  
     1560
    15611561  int row, column, i, j;
    15621562  int rr = nrows(A);
     
    15741574              if(A[row, j]!=0)
    15751575                {
    1576                   transform = unitMatrix(cc);
    1577                   transform[j,j] = 0;
    1578                   transform[column, column] = 0;
    1579                   transform[column,j] = 1;
    1580                   transform[j,column] = 1;
    1581                   q = q*transform;
     1576                  transform = unitMatrix(cc);
     1577                  transform[j,j] = 0;
     1578                  transform[column, column] = 0;
     1579                  transform[column,j] = 1;
     1580                  transform[j,column] = 1;
     1581                  q = q*transform;
    15821582                  A = A*transform;
    15831583                  break;
     
    16061606              q = q*transform;
    16071607              A = A*transform;
    1608               // A;
     1608              // A;
    16091609            }
    16101610        }
     
    16121612        {
    16131613          transform = unitMatrix(cc);
    1614           transform[column,column] = -1;
    1615           q = q*transform;
     1614          transform[column,column] = -1;
     1615          q = q*transform;
    16161616          A = A*transform;
    16171617        }
     
    16231623            transform[column,j]=transform[column,j]+1;}
    16241624          q = q*transform;
    1625           A = A*transform; 
     1625          A = A*transform;
    16261626        }
    16271627      }
     
    16371637  "EXAMPLE:"; echo=2;
    16381638
    1639    intmat M[2][5] = 
     1639   intmat M[2][5] =
    16401640     1, 2, 3, 4, 0,
    16411641     0,10,20,30, 1;
     
    16441644   print(hermiteNormalForm(M));
    16451645
    1646    intmat T[3][4] = 
     1646   intmat T[3][4] =
    16471647     3,3,3,3,
    16481648     2,1,3,0,
     
    16521652   print(hermiteNormalForm(T));
    16531653
    1654    intmat A[4][5] = 
     1654   intmat A[4][5] =
    16551655     1,2,3,2,2,
    16561656     1,2,3,4,0,
     
    16691669
    16701670  intvec v;
    1671  
     1671
    16721672  for( ; i > 0; i-- )
    16731673  {
    16741674    v = m[1..r, i];
    16751675    if( !isZeroElement(v, #) )
    1676     { 
     1676    {
    16771677      return (0);
    16781678    }
     
    17091709 intmat m[5][2]=mDeg(a)-mDeg(b),mDeg(b)-mDeg(c),mDeg(c)-mDeg(d),mDeg(d)-mDeg(e),mDeg(e)-mDeg(f);
    17101710 m=transpose(m);
    1711  areZeroElementes(m);
    1712  areZeroElementes(m,tt);
     1711 areZeroElements(m);
     1712 areZeroElements(m,tt);
    17131713}
    17141714
     
    17361736      }
    17371737    }
    1738    
     1738
    17391739  }
    17401740  if( !defined(H) )
     
    18011801  isZeroElement(v1);
    18021802  isZeroElement(v1, tt);
    1803  
     1803
    18041804  intvec v2 = mDeg(a) - mDeg(c);
    18051805  v2;
    18061806  isZeroElement(v2);
    18071807  isZeroElement(v2, tt);
    1808  
     1808
    18091809  intvec v3 = mDeg(e) - mDeg(f);
    18101810  v3;
    18111811  isZeroElement(v3);
    18121812  isZeroElement(v3, tt);
    1813  
     1813
    18141814  intvec v4 = mDeg(c) - mDeg(d);
    18151815  v4;
     
    18221822proc defineHomogeneous(poly f, list #)
    18231823"USAGE: defineHomogeneous(f[, G]); polynomial f, integer matrix G
    1824 PURPOSE: Yields a matrix which has to be appended to the grading group matrix to make the 
     1824PURPOSE: Yields a matrix which has to be appended to the grading group matrix to make the
    18251825polynomial f homogeneous in the grading by grad.
    18261826EXAMPLE: example defineHomogeneous; shows an example
    18271827"
    18281828{
    1829   if( size(#) > 0 ) 
    1830   {
    1831     if( typeof(#[1]) == "intmat" ) 
     1829  if( size(#) > 0 )
     1830  {
     1831    if( typeof(#[1]) == "intmat" )
    18321832    {
    18331833      intmat grad = #[1];
     
    18591859
    18601860  ring r =0,(x,y,z),dp;
    1861   intmat grad[2][3] = 
     1861  intmat grad[2][3] =
    18621862    1,0,1,
    18631863    0,1,1;
     
    18701870  M;
    18711871  defineHomogeneous(f, grad) == M;
    1872  
     1872
    18731873  isHomogeneous(f);
    18741874  setBaseMultigrading(grad, M);
     
    19061906    return(s);
    19071907  }
    1908   return(1==0); 
     1908  return(1==0);
    19091909}
    19101910example
     
    19441944//  listvar();
    19451945  def pre = preimage(f);
    1946  
     1946
    19471947//  "pre: ";  pre;
    19481948
     
    20352035
    20362036  ring r = 0,(x,y,z),dp;
    2037  
    2038  
     2037
     2038
    20392039
    20402040  // Setting degrees for preimage ring.;
    2041   intmat grad[3][3] = 
     2041  intmat grad[3][3] =
    20422042    1,0,0,
    20432043    0,1,0,
     
    20452045
    20462046  setBaseMultigrading(grad);
    2047  
     2047
    20482048  // grading on r:
    20492049  getVariableWeights();
     
    20602060
    20612061  listvar();
    2062  
     2062
    20632063  map f = r,-a2b6+b5+a3b+a2+ab,-a2b7-3a2b5+b4+a,a6-b6-b3+a2; f;
    20642064
    2065  
     2065
    20662066  // TODO: Unfortunately this is not a very spectacular example...:
    20672067  // Pushing forward f:
     
    20712071  getVariableWeights();
    20722072  getLattice();
    2073  
     2073
    20742074
    20752075  // only for the purpose of this example
     
    20822082proc equalMDeg(intvec exp1, intvec exp2, list #)
    20832083"USAGE: equalMDeg(exp1, exp2[, V]); intvec exp1, exp2, intmat V
    2084 PURPOSE: Tests if the exponent vectors of two monomials (given by exp1 and exp2) 
     2084PURPOSE: Tests if the exponent vectors of two monomials (given by exp1 and exp2)
    20852085represent the same multidegree.
    2086 NOTE: the integer matrix V encodes multidegrees of module components, 
     2086NOTE: the integer matrix V encodes multidegrees of module components,
    20872087if module component is present in exp1 and exp2
    20882088EXAMPLE: example equalMDeg; shows an example
     
    20942094  }
    20952095
    2096   if( exp1 == exp2) 
     2096  if( exp1 == exp2)
    20972097  {
    20982098    return (1==1);
     
    22282228    def g = groebner(a); // !!!!
    22292229
    2230     def b, aa; int j; 
     2230    def b, aa; int j;
    22312231    for( int i = ncols(a); i > 0; i-- )
    22322232    {
     
    22432243    }
    22442244    return(1==1);
    2245   } 
     2245  }
    22462246}
    22472247example
     
    22522252
    22532253  //Grading and Torsion matrices:
    2254   intmat M[3][3] = 
     2254  intmat M[3][3] =
    22552255    1,0,0,
    22562256    0,1,0,
     
    22782278  /////////////////////////////////////////////////////////
    22792279  // new Torsion matrix:
    2280   intmat T[3][4] = 
     2280  intmat T[3][4] =
    22812281    3,3,3,3,
    22822282    2,1,3,0,
    22832283    1,2,0,3;
    2284  
     2284
    22852285  setBaseMultigrading(M,T);
    22862286
     
    22892289  mDegPartition(f);
    22902290
    2291   // --------------------- 
     2291  // ---------------------
    22922292  g;
    22932293  isHomogeneous(g);
     
    22982298  ring R = 0, (x,y,z), dp;
    22992299
    2300   intmat A[2][3] = 
     2300  intmat A[2][3] =
    23012301    0,0,1,
    23022302    3,2,1;
    2303   intmat T[2][1] = 
    2304     -1, 
     2303  intmat T[2][1] =
     2304    -1,
    23052305     4;
    23062306  setBaseMultigrading(A, T);
     
    23142314  mDegPartition(x3 -y2z + x2 -y3 + z + 1);
    23152315
    2316  
     2316
    23172317  module N = gen(1) + (x+y) * gen(2), z*gen(3);
    23182318
    23192319  intmat V[2][3] = 0; // 1, 2, 3,  4, 5, 6; //  column-wise weights of components!!??
    2320  
     2320
    23212321  vector v1, v2;
    2322  
     2322
    23232323  v1 = setModuleGrading(N[1], V); v1;
    23242324  mDegPartition(v1);
     
    23332333  print( mDeg(N) );
    23342334
    2335   /////////////////////////////////////// 
    2336 
    2337   V = 
    2338     1, 2, 3, 
     2335  ///////////////////////////////////////
     2336
     2337  V =
     2338    1, 2, 3,
    23392339    4, 5, 6;
    23402340
     
    23512351  print( mDeg(N) );
    23522352
    2353   /////////////////////////////////////// 
    2354 
    2355   V = 
    2356     0, 0, 0, 
     2353  ///////////////////////////////////////
     2354
     2355  V =
     2356    0, 0, 0,
    23572357    4, 1, 0;
    23582358
     
    23942394  {
    23952395    intmat V = getModuleGrading(A);
    2396  
     2396
    23972397    if( nrows(V) != r )
    23982398    {
     
    24002400    }
    24012401  }
    2402                            
     2402
    24032403  if( A == 0 )
    24042404  {
     
    24192419    A = A - lead(A);
    24202420    while( size(A) > 0 )
    2421     { 
     2421    {
    24222422      v = leadexp(A); //  v;
    24232423      m = max( m, M * v, r ); // ????
     
    24402440    A = A - lead(A);
    24412441    while( size(A) > 0 )
    2442     { 
     2442    {
    24432443      v = leadexp(A); //  v;
    24442444
     
    24652465      {
    24662466        G[j, i] = d[j];
    2467       }     
     2467      }
    24682468    }
    24692469    return(G);
     
    24772477    for( i = ncols(A); i > 0; i-- )
    24782478    {
    2479       v = getGradedGenerator(A, i); 
    2480 
    2481       // G[1..r, i] 
     2479      v = getGradedGenerator(A, i);
     2480
     2481      // G[1..r, i]
    24822482      d = mDeg(v);
    24832483
     
    24852485      {
    24862486        G[j, i] = d[j];
    2487       }     
     2487      }
    24882488
    24892489    }
     
    24912491    return(G);
    24922492  }
    2493  
     2493
    24942494}
    24952495example
     
    25152515
    25162516  setBaseMultigrading(A, Ta);
    2517  
     2517
    25182518  mDeg( x*x*y );
    2519  
     2519
    25202520  mDeg( y*y*y*x );
    2521  
     2521
    25222522  mDeg( x*y + x + 1 );
    25232523
     
    25292529
    25302530//  "ideal:";
    2531  
     2531
    25322532  ideal I = y*x*x, x*y*y*y;
    25332533  print( mDeg(I) );
     
    25372537
    25382538//  "vectors:";
    2539  
     2539
    25402540  intmat B[2][2] = 0, 1, 1, 0;
    25412541  print(B);
    2542  
     2542
    25432543  mDeg( setModuleGrading(y*y*y*gen(2), B ));
    25442544  mDeg( setModuleGrading(x*x*gen(1), B ));
    25452545
    2546  
     2546
    25472547  vector V = x*gen(1) + y*gen(2);
    25482548  V = setModuleGrading(V, B);
     
    25512551  vector v1 = setModuleGrading([0, 0, 0], B);
    25522552  print( mDeg( v1 ) );
    2553  
     2553
    25542554  vector v2 = setModuleGrading([0], B);
    25552555  print( mDeg( v2 ) );
    25562556
    25572557//  "module:";
    2558  
     2558
    25592559  module D = x*gen(1), y*gen(2);
    25602560  D;
    25612561  D = setModuleGrading(D, B);
    25622562  print( mDeg( D ) );
    2563  
     2563
    25642564
    25652565  module DD = [0, 0],[0, 0, 0];
     
    25852585{ // TODO: What about an ideal or module???
    25862586
    2587   if( typeof(p) == "poly" ) 
    2588   {
    2589     ideal I;     
     2587  if( typeof(p) == "poly" )
     2588  {
     2589    ideal I;
    25902590    poly mp, t, tt;
    25912591  }
     
    25942594    if(  typeof(p) == "vector" )
    25952595    {
    2596       module I;     
     2596      module I;
    25972597      vector mp, t, tt;
    25982598    }
     
    26052605  if( typeof(p) == "vector" )
    26062606  {
    2607     intmat V = getModuleGrading(p); 
     2607    intmat V = getModuleGrading(p);
    26082608  }
    26092609  else
     
    26122612  }
    26132613
    2614   if( size(p) > 1) 
     2614  if( size(p) > 1)
    26152615  {
    26162616    intvec m;
     
    26192619    {
    26202620      m = leadexp(p);
    2621       mp = lead(p); 
     2621      mp = lead(p);
    26222622      p = p - lead(p);
    26232623      tt = p; t = 0;
    26242624
    26252625      while( size(tt) > 0 )
    2626       { 
     2626      {
    26272627        // TODO: we do not cache matrices (M,T,H,V), which remain the same :(
    26282628        // TODO: we need some low-level procedure with all these arguments...!
    2629         if( equalMDeg( leadexp(tt), m, V  ) ) 
     2629        if( equalMDeg( leadexp(tt), m, V  ) )
    26302630        {
    26312631          mp = mp + lead(tt); // "mp", mp;
     
    26742674
    26752675  mDegPartition(f);
    2676  
     2676
    26772677  vector v = xy*gen(1)-x3y2*gen(2)+x4y*gen(3);
    26782678  intmat B[2][3]=1,-1,-2,0,0,1;
    26792679  v = setModuleGrading(v,B);
    26802680  getModuleGrading(v);
    2681  
     2681
    26822682  mDegPartition(v, B);
    26832683}
     
    26892689{
    26902690  intmat A[n][n];
    2691  
     2691
    26922692  for( int i = n; i > 0; i-- )
    26932693  {
     
    27042704"
    27052705USAGE: finestMDeg(r); ring r
    2706 RETURN: ring, r endowed with the finest multigrading 
     2706RETURN: ring, r endowed with the finest multigrading
    27072707TODO: not yet...
    27082708"
     
    27332733
    27342734
    2735   if( n > 0) 
    2736   {
    2737 
    2738     intmat L[N][n]; 
     2735  if( n > 0)
     2736  {
     2737
     2738    intmat L[N][n];
    27392739    //  list L;
    27402740    int j = n;
     
    27442744      p = I[i];
    27452745
    2746       if( size(p) > 1 ) 
     2746      if( size(p) > 1 )
    27472747      {
    27482748        intvec m0 = leadexp(p);
     
    27592759
    27602760    print(L);
    2761     setBaseMultigrading(A, L);     
    2762   } 
     2761    setBaseMultigrading(A, L);
     2762  }
    27632763  else
    27642764  {
     
    27972797
    27982798  if( size(#) > 0 and typeof(#[1]) == "intmat" )
    2799   { 
     2799  {
    28002800    attrib(F, "P", #[1]);
    28012801  }
     
    28832883     ERROR("Sorry: cannot find 'hilbert' command from 4ti2. Please install 4ti2!");
    28842884   }
    2885    
     2885
    28862886//--------------------------------------------------------------------------
    28872887// Initialization and Sanity Checks
     
    29392939   j=system("sh","hilbert -q -n sing4ti2"); ////////// be quiet + no loggin!!!
    29402940
    2941    j=system("sh", "awk \'BEGIN{ORS=\",\";}{print $0;}\' sing4ti2.hil " + 
    2942                 "| sed s/[\\\ \\\t\\\v\\\f]/,/g " + 
    2943                 "| sed s/,+/,/g|sed s/,,/,/g " + 
    2944                 "| sed s/,,/,/g " + 
     2941   j=system("sh", "awk \'BEGIN{ORS=\",\";}{print $0;}\' sing4ti2.hil " +
     2942                "| sed s/[\\\ \\\t\\\v\\\f]/,/g " +
     2943                "| sed s/,+/,/g|sed s/,,/,/g " +
     2944                "| sed s/,,/,/g " +
    29452945                "> sing4ti2.converted" );
    29462946
     
    29622962   string ergstr = "intvec erglist = " + s + "0;";
    29632963   execute(ergstr);
    2964  
     2964
    29652965   //   print(erglist);
    2966  
     2966
    29672967   int Rnc = erglist[1];
    29682968   int Rnr = erglist[2];
    2969    
     2969
    29702970   intmat R[Rnr][Rnc];
    29712971
     
    30353035
    30363036/******************************************************/
    3037 proc mDegBasis(intvec d) 
     3037proc mDegBasis(intvec d)
    30383038"
    30393039USAGE: multidegree d
     
    30793079
    30803080    intmat AA[nr][nc + 2 * n];
    3081     AA[1..nr, 1.. nc] = A[1..nr, 1.. nc]; 
    3082     AA[1..nr, nc + (1.. n)] = T[1..nr, 1.. n]; 
    3083     AA[1..nr, nc + n + (1.. n)] = -T[1..nr, 1.. n]; 
     3081    AA[1..nr, 1.. nc] = A[1..nr, 1.. nc];
     3082    AA[1..nr, nc + (1.. n)] = T[1..nr, 1.. n];
     3083    AA[1..nr, nc + n + (1.. n)] = -T[1..nr, 1.. n];
    30843084
    30853085
    30863086    //      print ( AA );
    30873087
    3088     intmat K = leftKernelZ(( AA ) ); // 
     3088    intmat K = leftKernelZ(( AA ) ); //
    30893089
    30903090    //      print(K);
     
    30953095    //      "!";
    30963096
    3097     intmat B = hilbert4ti2intmat(transpose(KK), 1); 
     3097    intmat B = hilbert4ti2intmat(transpose(KK), 1);
    30983098
    30993099    //      "!";      print(B);
     
    31063106
    31073107
    3108   int i; 
     3108  int i;
    31093109  int nnr = nrows(B);
    31103110  int nnc = ncols(B);
     
    31653165  intvec v1=4,0;
    31663166  intvec v2=4,4;
    3167  
     3167
    31683168  intmat g3[1][2]=1,1;
    31693169  setBaseMultigrading(g3);
     
    31713171  v3;
    31723172  mDegBasis(v3);
    3173  
     3173
    31743174  setBaseMultigrading(g1,l);
    31753175  mDegBasis(v1);
    31763176  setBaseMultigrading(g2);
    31773177  mDegBasis(v2);
    3178  
     3178
    31793179  intmat M[2][2] = 1, -1, -1, 1;
    31803180  intvec d = -2, 2;
     
    32423242"
    32433243{
    3244   if( isHomogeneous(I, "checkGens") == 0) 
    3245   { 
    3246     ERROR ("Sorry: inhomogeneous input!"); 
    3247   } 
     3244  if( isHomogeneous(I, "checkGens") == 0)
     3245  {
     3246    ERROR ("Sorry: inhomogeneous input!");
     3247  }
    32483248  module S = syz(I);
    32493249  S = setModuleGrading(S, mDeg(I));
     
    32533253{
    32543254  "EXAMPLE:"; echo=2;
    3255  
     3255
    32563256
    32573257  ring r = 0,(x,y,z,w),dp;
     
    32613261  setBaseMultigrading(MM);
    32623262  module M = ideal(  xw-yz, x2z-y3, xz2-y2w, yw2-z3);
    3263  
    3264  
     3263
     3264
    32653265  intmat v[2][nrows(M)]=
    32663266    1,
    32673267    0;
    3268  
     3268
    32693269  M = setModuleGrading(M, v);
    32703270
     
    32853285PURPOSE: computes the multigraded 'modulo' module of I and J
    32863286RETURNS: module, see 'modulo' command
    3287 NOTE: I and J should have the same multigrading, and their 
     3287NOTE: I and J should have the same multigrading, and their
    32883288generators must be multigraded homogeneous
    32893289EXAMPLE: example mDegModulo; shows an example
     
    32913291{
    32923292  if( (isHomogeneous(I, "checkGens") == 0) or (isHomogeneous(J, "checkGens") == 0) )
    3293   { 
    3294     ERROR ("Sorry: inhomogeneous input!"); 
    3295   } 
     3293  {
     3294    ERROR ("Sorry: inhomogeneous input!");
     3295  }
    32963296  module K = modulo(I, J);
    32973297  K = setModuleGrading(K, mDeg(I));
     
    33123312
    33133313  "Multidegrees: "; print(mDeg(h1));
    3314    
     3314
    33153315  // Let's compute modulo(h1, h2):
    33163316  def K = mDegModulo(h1, h2); K;
     
    33263326"USAGE: mDegGroebner(I); I is a poly/vector/ideal/module
    33273327PURPOSE: computes the multigraded standard/groebner basis of I
    3328 NOTE: I must be multigraded homogeneous 
     3328NOTE: I must be multigraded homogeneous
    33293329RETURNS: ideal/module, the computed basis
    33303330EXAMPLE: example mDegGroebner; shows an example
    33313331"
    33323332{
    3333   if( isHomogeneous(I) == 0) 
    3334   { 
    3335     ERROR ("Sorry: inhomogeneous input!"); 
    3336   } 
     3333  if( isHomogeneous(I) == 0)
     3334  {
     3335    ERROR ("Sorry: inhomogeneous input!");
     3336  }
    33373337
    33383338  def S = groebner(I);
    3339  
     3339
    33403340  if( typeof(I) == "module" or typeof(I) == "vector" )
    33413341  {
    3342     S = setModuleGrading(S, getModuleGrading(I));     
     3342    S = setModuleGrading(S, getModuleGrading(I));
    33433343  }
    33443344
     
    33573357  setBaseMultigrading(MM);
    33583358
    3359  
     3359
    33603360  module M = ideal(  xw-yz, x2z-y3, xz2-y2w, yw2-z3);
    3361  
    3362  
     3361
     3362
    33633363  intmat v[2][nrows(M)]=
    33643364    1,
    33653365    0;
    3366  
     3366
    33673367  M = setModuleGrading(M, v);
    33683368
     
    33973397proc mDegResolution(def I, int ll, list #)
    33983398"USAGE: mDegResolution(I,l,[f]); I is poly/vector/ideal/module; l,f are integers
    3399 PURPOSE: computes the multigraded resolution of I of the length l, 
    3400 or the whole resolution if l is zero. Returns minimal resolution if an optional 
     3399PURPOSE: computes the multigraded resolution of I of the length l,
     3400or the whole resolution if l is zero. Returns minimal resolution if an optional
    34013401argument 1 is supplied
    34023402NOTE: input must have multigraded-homogeneous generators.
    3403 The returned list is truncated beginning with the first zero differential.   
     3403The returned list is truncated beginning with the first zero differential.
    34043404RETURNS: list, the computed resolution
    34053405EXAMPLE: example mDegResolution; shows an example
    34063406"
    34073407{
    3408   if( isHomogeneous(I, "checkGens") == 0) 
    3409   { 
    3410     ERROR ("Sorry: inhomogeneous input!"); 
    3411   } 
     3408  if( isHomogeneous(I, "checkGens") == 0)
     3409  {
     3410    ERROR ("Sorry: inhomogeneous input!");
     3411  }
    34123412
    34133413  def R = res(I, ll, #); list L = R; int l = size(L);
     
    34183418  }
    34193419
    3420   int i; 
     3420  int i;
    34213421  for( i = 2; i <= l; i++ )
    34223422  {
     
    34293429    }
    34303430  }
    3431  
     3431
    34323432  return (L);
    34333433
    3434  
     3434
    34353435}
    34363436example
    34373437{
    34383438  "EXAMPLE:"; echo=2;
    3439  
     3439
    34403440  ring r = 0,(x,y,z,w),dp;
    34413441
     
    34463446  setBaseMultigrading(M);
    34473447
    3448  
     3448
    34493449  module m= ideal(  xw-yz, x2z-y3, xz2-y2w, yw2-z3);
    3450  
     3450
    34513451  isHomogeneous(ideal(  xw-yz, x2z-y3, xz2-y2w, yw2-z3), "checkGens");
    3452  
     3452
    34533453  ideal A = xw-yz, x2z-y3, xz2-y2w, yw2-z3;
    34543454
    34553455  int j;
    3456  
     3456
    34573457  for(j=1; j<=ncols(A); j++)
    34583458  {
    34593459    mDegPartition(A[j]);
    34603460  }
    3461  
     3461
    34623462  intmat v[2][1]=
    34633463    1,
    34643464    0;
    3465  
     3465
    34663466  m = setModuleGrading(m, v);
    34673467
     
    34903490
    34913491  /////////////////////////////////////////////////////////////////////////////
    3492  
     3492
    34933493  L = mDegResolution(maxideal(1), 0, 1);
    34943494
     
    35003500    "Multigrading: "; print(mDeg(L[j]));
    35013501  }
    3502  
     3502
    35033503  kill v;
    3504  
     3504
    35053505
    35063506  def h = hilbertSeries(m);
     
    35093509  numerator1;
    35103510  factorize(numerator1);
    3511  
     3511
    35123512  denominator1;
    35133513  factorize(denominator1);
     
    35233523proc hilbertSeries(def I)
    35243524"USAGE: hilbertSeries(I); I is poly/vector/ideal/module
    3525 PURPOSE: computes the multigraded Hilbert Series of M 
    3526 NOTE: input must have multigraded-homogeneous generators. 
     3525PURPOSE: computes the multigraded Hilbert Series of M
     3526NOTE: input must have multigraded-homogeneous generators.
    35273527Multigrading should be positive.
    3528 RETURNS: a ring in variables t_(i), s_(i), with polynomials 
    3529 numerator1 and denominator1 and muturally prime numerator2 
     3528RETURNS: a ring in variables t_(i), s_(i), with polynomials
     3529numerator1 and denominator1 and muturally prime numerator2
    35303530and denominator2, quotients of which give the series.
    35313531EXAMPLE: example hilbertSeries; shows an example
    35323532"
    35333533{
    3534    
     3534
    35353535  if( !isFreeRepresented() )
    35363536  {
     
    35383538    //ERROR("SORRY: ONLY TORSION-FREE CASE (POSITIVE GRADING)");
    35393539  }
    3540    
     3540
    35413541  int i, j, k, v;
    35423542
    35433543  intmat M = getVariableWeights();
    3544  
     3544
    35453545  int cc = ncols(M);
    35463546  int n = nrows(M);
     
    35543554
    35553555  int l = size(RES);
    3556  
     3556
    35573557  list L; L[l + 1] = 0;
    35583558
     
    35613561    intmat zeros[n][1];
    35623562    L[1] = zeros;
    3563   } 
     3563  }
    35643564  else
    35653565  {
     
    35713571    L[j + 1] = mDeg(RES[j]);
    35723572  }
    3573  
     3573
    35743574  l++;
    35753575
    35763576  ring R = 0,(t_(1..n),s_(1..n)),dp;
    3577  
    3578   ideal units; 
     3577
     3578  ideal units;
    35793579  for( i=n; i>=1; i--)
    35803580  {
    35813581    units[i] = (var(i) * var(n + i) - 1);
    35823582  }
    3583  
     3583
    35843584  qring Q = std(units);
    3585  
     3585
    35863586  // TODO: should not it be a quotient ring depending on Torsion???
    35873587  // I am not sure about what to do in the torsion case, but since
     
    35923592  poly monom, summand, numerator;
    35933593  poly denominator = 1;
    3594  
     3594
    35953595  for( i = 1; i <= cc; i++)
    35963596  {
     
    36033603      {
    36043604        monom = monom * (var(k)^(v));
    3605       } 
     3605      }
    36063606      else
    36073607      {
     
    36093609      }
    36103610    }
    3611    
     3611
    36123612    if( monom == 1)
    36133613    {
     
    36173617    denominator = denominator * (1 - monom);
    36183618  }
    3619  
    3620   for( j = 1; j<= l; j++) 
     3619
     3620  for( j = 1; j<= l; j++)
    36213621  {
    36223622    summand = 0;
     
    36323632        {
    36333633          monom = monom * (var(k)^v);
    3634         } 
     3634        }
    36353635        else
    36363636        {
     
    36423642    numerator = numerator - (-1)^j * summand;
    36433643  }
    3644  
     3644
    36453645  if( denominator == 0 )
    36463646  {
    36473647    ERROR("Multigrading not positive.");
    3648   } 
    3649  
     3648  }
     3649
    36503650  poly denominator1 = denominator;
    36513651  poly numerator1 = numerator;
     
    36813681  "The s_(i)-variables are defined to be the inverse of the t_(i)-variables.";
    36823682  " ------------ ";
    3683  
     3683
    36843684  return(Q);
    36853685}
     
    36873687{
    36883688  "EXAMPLE:"; echo=2;
    3689  
     3689
    36903690  ring r = 0,(x,y,z,w),dp;
    36913691  intmat g[2][4]=
     
    36933693    0,1,3,4;
    36943694  setBaseMultigrading(g);
    3695  
     3695
    36963696  module M = ideal(xw-yz, x2z-y3, xz2-y2w, yw2-z3);
    36973697  intmat V[2][1]=
     
    37053705  factorize(numerator2);
    37063706  factorize(denominator2);
    3707  
     3707
    37083708  kill g, h; setring r;
    37093709
     
    37113711    1,2,3,4,
    37123712    0,0,5,8;
    3713  
     3713
    37143714  setBaseMultigrading(g);
    3715  
     3715
    37163716  ideal I = x^2, y, z^3;
    37173717  I = std(I);
     
    37283728  mDeg(I);
    37293729  def h = hilbertSeries(I); setring h;
    3730  
     3730
    37313731  factorize(numerator2);
    37323732  factorize(denominator2);
     
    37353735  ////////////////////////////////////////////////
    37363736  ring R = 0,(x,y,z),dp;
    3737   intmat W[2][3] = 
     3737  intmat W[2][3] =
    37383738     1,1, 1,
    37393739     0,0,-1;
    37403740  setBaseMultigrading(W);
    37413741  ideal I = x3y,yz2,y2z,z4;
    3742  
     3742
    37433743  def h = hilbertSeries(I); setring h;
    3744  
     3744
    37453745  factorize(numerator2);
    37463746  factorize(denominator2);
     
    37493749  ////////////////////////////////////////////////
    37503750  ring R = 0,(x,y,z,a,b,c),dp;
    3751   intmat W[2][6] = 
     3751  intmat W[2][6] =
    37523752     1,1, 1,1,1,1,
    37533753     0,0,-1,0,0,0;
    37543754  setBaseMultigrading(W);
    37553755  ideal I = x3y,yz2,y2z,z4;
    3756  
     3756
    37573757  def h = hilbertSeries(I); setring h;
    3758  
     3758
    37593759  factorize(numerator2);
    37603760  factorize(denominator2);
    3761  
     3761
    37623762  kill R, W, h;
    37633763  ////////////////////////////////////////////////
    37643764  // This is example 5.3.9. from Robbianos book.
    3765  
     3765
    37663766  ring R = 0,(x,y,z,w),dp;
    3767   intmat W[1][4] = 
     3767  intmat W[1][4] =
    37683768     1,1, 1,1;
    37693769  setBaseMultigrading(W);
     
    37713771
    37723772  hilb(std(I));
    3773  
     3773
    37743774  def h = hilbertSeries(I); setring h;
    3775  
     3775
    37763776  numerator1;
    37773777  denominator1;
     
    37793779  factorize(numerator2);
    37803780  factorize(denominator2);
    3781  
     3781
    37823782
    37833783  kill h;
     
    37883788
    37893789  hilb(std(I2));
    3790  
     3790
    37913791  def h = hilbertSeries(I2); setring h;
    37923792
     
    37983798  ////////////////////////////////////////////////
    37993799  setring R;
    3800  
     3800
    38013801  W = 2,2,2,2;
    3802  
     3802
    38033803  setBaseMultigrading(W);
    38043804
     
    38103810
    38113811  kill w;
    3812  
     3812
    38133813
    38143814  def h = hilbertSeries(I2); setring h;
    38153815
    3816  
     3816
    38173817  numerator1; denominator1;
    38183818  kill h;
    38193819
    3820  
     3820
    38213821  kill R, W;
    38223822
     
    38323832
    38333833  hilb(std(I));
    3834  
     3834
    38353835  def h = hilbertSeries(I); setring h;
    38363836
     
    38483848
    38493849  numerator1; denominator1;
    3850  
    3851   kill h; 
    3852   ////////////////////////////////////////////////
    3853   setring R;
    3854 
    3855   I = x^5; I;
    3856 
    3857   hilb(std(I));
    3858   hilb(std(I), 1);
    3859 
    3860   def h = hilbertSeries(I); setring h;
    3861 
    3862   numerator1; denominator1;
    3863  
    3864  
    3865   kill h; 
    3866   ////////////////////////////////////////////////
    3867   setring R;
    3868 
    3869   I = x^10; I;
    3870 
    3871   hilb(std(I));
    3872 
    3873   def h = hilbertSeries(I); setring h;
    3874 
    3875   numerator1; denominator1;
    38763850
    38773851  kill h;
     
    38793853  setring R;
    38803854
    3881   module M = 1;
    3882 
    3883   M = setModuleGrading(M, W);
    3884 
    3885  
    3886   hilb(std(M));
    3887 
    3888   def h = hilbertSeries(M); setring h;
     3855  I = x^5; I;
     3856
     3857  hilb(std(I));
     3858  hilb(std(I), 1);
     3859
     3860  def h = hilbertSeries(I); setring h;
    38893861
    38903862  numerator1; denominator1;
     3863
    38913864
    38923865  kill h;
     
    38943867  setring R;
    38953868
    3896   kill M; module M = x^5*gen(1);
    3897 //  intmat V[1][3] = 0; // TODO: this would lead to a wrong result!!!?
    3898   intmat V[1][1] = 0; // all gen(i) of degree 0!
    3899 
    3900   M = setModuleGrading(M, V);
    3901 
    3902   hilb(std(M));
    3903 
    3904   def h = hilbertSeries(M); setring h;
     3869  I = x^10; I;
     3870
     3871  hilb(std(I));
     3872
     3873  def h = hilbertSeries(I); setring h;
    39053874
    39063875  numerator1; denominator1;
     
    39083877  kill h;
    39093878  ////////////////////////////////////////////////
    3910   setring R;
    3911 
    3912   module N = x^5*gen(3);
    3913 
    3914   kill V;
    3915  
    3916   intmat V[1][3] = 0; // all gen(i) of degree 0!
    3917 
    3918   N = setModuleGrading(N, V);
    3919      
    3920   hilb(std(N));
    3921 
    3922   def h = hilbertSeries(N); setring h;
     3879  setring R;
     3880
     3881  module M = 1;
     3882
     3883  M = setModuleGrading(M, W);
     3884
     3885
     3886  hilb(std(M));
     3887
     3888  def h = hilbertSeries(M); setring h;
    39233889
    39243890  numerator1; denominator1;
     
    39263892  kill h;
    39273893  ////////////////////////////////////////////////
    3928   setring R;
    3929 
    3930  
     3894  setring R;
     3895
     3896  kill M; module M = x^5*gen(1);
     3897//  intmat V[1][3] = 0; // TODO: this would lead to a wrong result!!!?
     3898  intmat V[1][1] = 0; // all gen(i) of degree 0!
     3899
     3900  M = setModuleGrading(M, V);
     3901
     3902  hilb(std(M));
     3903
     3904  def h = hilbertSeries(M); setring h;
     3905
     3906  numerator1; denominator1;
     3907
     3908  kill h;
     3909  ////////////////////////////////////////////////
     3910  setring R;
     3911
     3912  module N = x^5*gen(3);
     3913
     3914  kill V;
     3915
     3916  intmat V[1][3] = 0; // all gen(i) of degree 0!
     3917
     3918  N = setModuleGrading(N, V);
     3919
     3920  hilb(std(N));
     3921
     3922  def h = hilbertSeries(N); setring h;
     3923
     3924  numerator1; denominator1;
     3925
     3926  kill h;
     3927  ////////////////////////////////////////////////
     3928  setring R;
     3929
     3930
    39313931  module S = M + N;
    3932  
     3932
    39333933  S = setModuleGrading(S, V);
    39343934
     
    39523952"
    39533953{
    3954   if( 2*size(v) != nvars(h) ) 
     3954  if( 2*size(v) != nvars(h) )
    39553955  {
    39563956    ERROR("Wrong input/size!");
     
    39763976    V[i] = var(i) - k;
    39773977  }
    3978    
     3978
    39793979  V = groebner(V);
    3980    
    3981   n = NF(n, V); 
    3982   d = NF(d, V); 
     3980
     3981  n = NF(n, V);
     3982  d = NF(d, V);
    39833983
    39843984  n;
     
    39893989    ERROR("Sorry: denominator is zero!");
    39903990  }
    3991  
     3991
    39923992  if( n == 0 )
    39933993  {
     
    39963996
    39973997  poly g = gcd(n, d);
    3998  
     3998
    39993999  if( g != leadcoef(g) )
    40004000  {
     
    40054005  n;
    40064006  d;
    4007    
    4008    
     4007
     4008
    40094009  for( i = N; i > 0; i -- )
    40104010  {
     
    40124012    n;
    40134013    d;
    4014      
     4014
    40154015    k = v[i];
    40164016    k;
    4017      
     4017
    40184018    n = subst(n, var(i), k);
    40194019    d = subst(d, var(i), k);
    4020    
     4020
    40214021    if( k != 0 )
    40224022    {
     
    40344034    ERROR("Sorry: denominator is zero!");
    40354035  }
    4036  
     4036
    40374037  if( n == 0 )
    40384038  {
     
    40414041
    40424042  poly g = gcd(n, d);
    4043  
     4043
    40444044  if( g != leadcoef(g) )
    40454045  {
     
    40504050  n;
    40514051  d;
    4052    
     4052
    40534053  if( n != leadcoef(n) || d != leadcoef(d) )
    40544054  {
     
    40804080ideal I = mDegBasis(0);
    40814081ideal J = attrib(I,"ZeroPart");
    4082 /* 
     4082/*
    40834083I am not quite sure what this ZeroPart is anymore. I thought it
    40844084should contain all monomials of degree 0, but then apparently 1 should
     
    40944094  setBaseMultigrading(A);
    40954095  isPositive();
    4096  
     4096
    40974097  intmat A[1][2]=1,1;
    40984098  setBaseMultigrading(A);
     
    41274127
    41284128  example mDegResolution;
    4129  
    4130   "// ******************* example hilbertSeries ************************//"; 
     4129
     4130  "// ******************* example hilbertSeries ************************//";
    41314131  example hilbertSeries;
    41324132
     
    41424142  "d: ";
    41434143  print(md);
    4144  
     4144
    41454145  "M: ";
    41464146  module LL = M; // + L for d+1
     
    41494149
    41504150
    4151   intmat V = getModuleGrading(M); 
     4151  intmat V = getModuleGrading(M);
    41524152  intvec vi;
    4153   int s = nrows(M); 
     4153  int s = nrows(M);
    41544154  int r = nrows(V);
    41554155  int i;
    41564156  module L; def B;
    4157   for (i=s; i>0; i--) 
     4157  for (i=s; i>0; i--)
    41584158  {
    41594159    "comp: ", i;
     
    41794179  print(L);
    41804180  print(mDeg(L));
    4181    
     4181
    41824182  V = getModuleGrading(L);
    41834183
     
    41904190    }
    41914191  }
    4192  
     4192
    41934193  L = simplify(L, 2);
    41944194  L = setModuleGrading(L, V);
    41954195  print(L);
    41964196  print(mDeg(L));
    4197  
     4197
    41984198  return(L);
    41994199}
     
    42044204  // TODO!
    42054205  ring r = 32003, (x,y), dp;
    4206  
    4207   intmat M[2][2] = 
    4208     1, 0, 
     4206
     4207  intmat M[2][2] =
     4208    1, 0,
    42094209    0, 1;
    42104210
    42114211  setBaseMultigrading(M);
    42124212
    4213   intmat V[2][1] = 
    4214     0, 
     4213  intmat V[2][1] =
     4214    0,
    42154215    0;
    4216  
     4216
    42174217  "X:";
    42184218  module h1 = x;
     
    42314231
    42324232/******************************************************/
    4233 /* Some functions on lattices. 
    4234 TODO Tuebingen: - add functionality (see wiki) and 
     4233/* Some functions on lattices.
     4234TODO Tuebingen: - add functionality (see wiki) and
    42354235- adjust them to work for groups as well.*/
    42364236/******************************************************/
     
    42414241proc imageLattice(intmat Q, intmat L)
    42424242"USAGE: imageLattice(Q,L); Q and L are of type intmat
    4243 PURPOSE: compute an integral basis for the image of the 
     4243PURPOSE: compute an integral basis for the image of the
    42444244lattice L under the homomorphism of lattices Q.
    42454245RETURN: intmat
     
    42554255{
    42564256  "EXAMPLE:"; echo=2;
    4257    
     4257
    42584258  intmat Q[2][3] =
    42594259    1,2,3,
     
    42754275"
    42764276USAGE: intRank(A); intmat A
    4277 PURPOSE: compute the rank of the integral matrix A 
     4277PURPOSE: compute the rank of the integral matrix A
    42784278by computing a hermite normalform.
    42794279RETURNS: int
     
    42914291  {
    42924292    iszero = 1;
    4293    
     4293
    42944294    for ( i = 1; i <= nrows(B); i++ )
    42954295    {
    4296       if ( B[i,j] != 0 ) 
     4296      if ( B[i,j] != 0 )
    42974297      {
    42984298        iszero = 0;
     
    43004300      }
    43014301    }
    4302    
     4302
    43034303    if ( iszero == 1 )
    43044304    {
     
    43134313  {
    43144314    iszero = 1;
    4315    
     4315
    43164316    for ( j = 1; j <= ncols(B); j++ )
    43174317    {
    4318       if ( B[i,j] != 0 ) 
     4318      if ( B[i,j] != 0 )
    43194319      {
    43204320        iszero = 0;
     
    43224322      }
    43234323    }
    4324    
     4324
    43254325    if ( iszero == 1 )
    43264326    {
     
    43314331  int r = nrows(B) - nzerorows;
    43324332
    4333   if ( ncols(B) - nzerocols < r ) 
     4333  if ( ncols(B) - nzerocols < r )
    43344334  {
    43354335    r = ncols(B) - nzerocols;
    43364336  }
    4337  
     4337
    43384338  return(r);
    43394339}
     
    43594359proc isSublattice(intmat L, intmat S)
    43604360"USAGE: isSublattice(L, S); L, S are of tpye intmat
    4361 PURPOSE: checks whether the lattice created by L is a 
     4361PURPOSE: checks whether the lattice created by L is a
    43624362sublattice of the lattice created by S.
    4363 The procedure checks whether each generator of L is 
     4363The procedure checks whether each generator of L is
    43644364contained in S.
    43654365RETURN: 0 if false, 1 if true
     
    43694369  int a,b,g,i,j,k;
    43704370  intmat Ker;
    4371    
     4371
    43724372  // check whether each column v of L is contained in
    43734373  // the lattice generated by S
    43744374  for ( i = 1; i <= ncols(L); i++ )
    43754375  {
    4376    
     4376
    43774377    // v is the i-th column of L
    43784378    intvec v;
     
    43984398    }
    43994399
    4400    
     4400
    44014401    // check gcd
    44024402    Ker = kernelLattice(B);
     
    44104410
    44114411    g = R[1];
    4412    
     4412
    44134413    for ( j = 2; j <= size(R); j++ )
    44144414    {
     
    44164416    }
    44174417
    4418     if ( g != 1 ) 
     4418    if ( g != 1 )
    44194419    {
    44204420      return(0);
     
    44304430{
    44314431  "EXAMPLE:"; echo=2;
    4432    
     4432
    44334433  //ring R = 0,(x,y),dp;
    44344434  intmat S2[3][2]=
     
    44574457proc latticeBasis(intmat B)
    44584458"USAGE: latticeBasis(B); intmat B
    4459 PURPOSE: compute an integral basis for the lattice defined by 
     4459PURPOSE: compute an integral basis for the lattice defined by
    44604460the columns of B.
    44614461RETURNS: intmat
     
    44644464{
    44654465  int n = ncols(B);
    4466   int r = intRank(B); 
    4467  
    4468   if ( r == 0 ) 
     4466  int r = intRank(B);
     4467
     4468  if ( r == 0 )
    44694469  {
    44704470    intmat H[nrows(B)][1];
     
    44734473    for ( j = 1; j <= nrows(B); j++ )
    44744474    {
    4475       H[j,1] = 0;   
     4475      H[j,1] = 0;
    44764476    }
    44774477  }
     
    44804480    intmat H = hermiteNormalForm(B);;
    44814481
    4482     if (r < n) 
     4482    if (r < n)
    44834483    {
    44844484      // delete columns r+1 to n
    44854485      // should be identical with the function
    4486       // H = submat(H,1..nrows(H),1..r); 
     4486      // H = submat(H,1..nrows(H),1..r);
    44874487      // for matrices
    44884488      intmat Hdel[nrows(H)][r];
    44894489      int k;
    44904490      int m;
    4491      
     4491
    44924492      for ( k = 1; k <= nrows(H); k++ )
    44934493      {
     
    44984498      }
    44994499
    4500       H = Hdel;     
    4501     }
    4502   }
    4503  
    4504   return(H); 
    4505 } 
     4500      H = Hdel;
     4501    }
     4502  }
     4503
     4504  return(H);
     4505}
    45064506example
    45074507{
    45084508  "EXAMPLE:"; echo=2;
    4509  
     4509
    45104510  intmat L[3][3] =
    45114511    1,4,8,
     
    45384538  else
    45394539  {
    4540     if ( r == n ) 
     4540    if ( r == n )
    45414541    {
    45424542      intmat U[n][1];  // now all entries are zero.
     
    45494549      // delete columns 1 to r
    45504550      // should be identical with the function
    4551       // U = submat(U,1..nrows(U),r+1..); 
     4551      // U = submat(U,1..nrows(U),r+1..);
    45524552      // for matrices
    45534553      intmat Udel[nrows(U)][ncols(U) - r];
    45544554      int k;
    45554555      int m;
    4556      
     4556
    45574557      for ( k = 1; k <= nrows(U); k++ )
    45584558      {
     
    45634563      }
    45644564
    4565       U = Udel; 
     4565      U = Udel;
    45664566
    45674567    }
     
    45744574  "EXAMPLE"; echo = 2;
    45754575
    4576   intmat LL[3][4] = 
     4576  intmat LL[3][4] =
    45774577    1,4,7,10,
    45784578    2,5,8,11,
     
    46344634  int k;
    46354635  int m;
    4636      
     4636
    46374637  for ( k = 1; k <= nrows(Del); k++ )
    46384638  {
     
    46424642    }
    46434643  }
    4644  
     4644
    46454645  L = latticeBasis(Del);
    46464646
    4647   return(L); 
     4647  return(L);
    46484648
    46494649}
     
    46524652  "EXAMPLE"; echo = 2;
    46534653
    4654   intmat P[2][3] = 
     4654  intmat P[2][3] =
    46554655    2,6,10,
    46564656    4,8,12;
     
    46714671proc isPrimitiveSublattice(intmat A);
    46724672"USAGE: isPrimitiveSublattice(A); intmat A
    4673 PURPOSE: check whether the given set of integral vectors in ZZ^m, 
    4674 i.e. the columns of A, generate a primitive sublattice in ZZ^m 
    4675 (a direct summand of ZZ^m). 
     4673PURPOSE: check whether the given set of integral vectors in ZZ^m,
     4674i.e. the columns of A, generate a primitive sublattice in ZZ^m
     4675(a direct summand of ZZ^m).
    46764676RETURNS: int, where 0 is false and 1 is true.
    46774677EXAMPLE: example isPrimitiveSublattice; shows an example
     
    46804680  intmat B = smithNormalForm(A);
    46814681  int r = intRank(B);
    4682  
    4683   if ( r == 0 ) 
     4682
     4683  if ( r == 0 )
    46844684  {
    46854685    return(1);
     
    46964696{
    46974697  "EXAMPLE"; echo = 2;
    4698  
     4698
    46994699  intmat A[3][2] =
    47004700    1,4,
     
    47044704  // should be 0
    47054705  int b = isPrimitiveSublattice(A);
    4706  
     4706
    47074707  kill A,b;
    47084708}
     
    47114711proc isIntegralSurjective(intmat P);
    47124712"USAGE: isIntegralSurjective(P); intmat P
    4713 PURPOSE: test whether the given linear map P of lattices is 
     4713PURPOSE: test whether the given linear map P of lattices is
    47144714surjective.
    47154715RETURNS: int, where 0 is false and 1 is true.
     
    47184718{
    47194719  int r = intRank(P);
    4720  
     4720
    47214721  if ( r < nrows(P) )
    47224722  {
     
    47344734{
    47354735  "EXAMPLE"; echo = 2;
    4736  
     4736
    47374737  intmat A[3][2] =
    47384738    1,3,5,
    47394739    2,4,6;
    4740    
     4740
    47414741  // should be 0
    47424742  int b = isIntegralSurjective(A);
    47434743  print(b);
    4744  
     4744
    47454745  kill A,b;
    47464746}
     
    47494749proc projectLattice(intmat B)
    47504750"USAGE: projectLattice(B); intmat B
    4751 PURPOSE: A set of vectors in ZZ^m is given as the columns of B. 
    4752 Then this function provides a linear map ZZ^m --> ZZ^n 
     4751PURPOSE: A set of vectors in ZZ^m is given as the columns of B.
     4752Then this function provides a linear map ZZ^m --> ZZ^n
    47534753having the primitive span of B its kernel.
    47544754RETURNS: intmat
     
    47654765  else
    47664766  {
    4767     if ( r == n ) 
     4767    if ( r == n )
    47684768    {
    47694769      intmat U[n][1]; // U now is the n-dim zero-vector
     
    47724772    {
    47734773      // we want a matrix with column operations so we transpose
    4774       list L = hermiteNormalForm(B, "transform"); //hermite(transpose(B), "transform"); 
    4775       intmat U = transpose(L[2]); 
     4774      list L = hermiteNormalForm(B, "transform"); //hermite(transpose(B), "transform");
     4775      intmat U = transpose(L[2]);
    47764776
    47774777      // delete rows 1 to r
     
    47794779      int k;
    47804780      int m;
    4781      
     4781
    47824782      for ( k = 1; k <= nrows(U) - r ; k++ )
    47834783      {
     
    47884788      }
    47894789
    4790       U = Udel; 
    4791      
    4792     }
    4793   }
    4794  
     4790      U = Udel;
     4791
     4792    }
     4793  }
     4794
    47954795  return(U);
    47964796}
     
    47984798{
    47994799  "EXAMPLE"; echo = 2;
    4800  
    4801   intmat B[4][2] = 
     4800
     4801  intmat B[4][2] =
    48024802    1,5,
    48034803    2,6,
    48044804    3,7,
    48054805    4,8;
    4806  
     4806
    48074807  // should result in a (2x4)-matrix with columns
    48084808  // [-1, 2], [2, −3], [-1, 0] and [0, 1].
    48094809  intmat U = projectLattice(B);
    4810  
     4810
    48114811}
    48124812
     
    48144814proc intersectLattices(intmat A, intmat B)
    48154815"USAGE: intersectLattices(A, B); intmat A, intmat B
    4816 PURPOSE: compute an integral basis for the intersection of the 
     4816PURPOSE: compute an integral basis for the intersection of the
    48174817lattices A and B.
    48184818RETURNS: intmat
     
    48454845  // delete all rows in K from ncols(A)+1 onwards
    48464846  intmat Bas[ncols(A)][ncols(K)];
    4847  
     4847
    48484848  for ( i = 1; i <= nrows(Bas); i++ )
    48494849  {
    4850     for ( j = 1; j <= ncols(Bas); j++ ) 
     4850    for ( j = 1; j <= ncols(Bas); j++ )
    48514851    {
    48524852      Bas[i,j] = K[i,j];
     
    48594859  int r = intRank(Cut);
    48604860
    4861   if ( r == 0 ) 
     4861  if ( r == 0 )
    48624862  {
    48634863    intmat Cutdel[nrows(Cut)][1]; // is now the zero-vector
     
    48724872    for ( i = 1; i <= nrows(Cutdel); i++ )
    48734873    {
    4874       for ( j = 1; j <= r; j++ ) 
     4874      for ( j = 1; j <= r; j++ )
    48754875      {
    48764876        Cutdel[i,j] = Cut[i,j];
     
    48864886{
    48874887  "EXAMPLE"; echo = 2;
    4888  
    4889   intmat A[3][2] = 
     4888
     4889  intmat A[3][2] =
    48904890    1,4,
    48914891    2,5,
    48924892    3,6;
    48934893
    4894   intmat B[3][2] = 
     4894  intmat B[3][2] =
    48954895    6,9,
    48964896    7,10,
    48974897    8,11;
    4898  
     4898
    48994899  // should result in a (2x4)-matrix with columns
    49004900  // [3, 3, 3], [0, 3, 6]
    49014901  intmat U = intersectLattices(A,B);
    4902  
     4902
    49034903}
    49044904
    49054905proc intInverse(intmat A);
    49064906"USAGE: intInverse(A); intmat A
    4907 PURPOSE: compute the integral inverse of the intmat A. 
     4907PURPOSE: compute the integral inverse of the intmat A.
    49084908If det(A) is neither 1 nor -1 an error is returned.
    49094909RETURNS: intmat
     
    49124912{
    49134913  int d = det(A);
    4914  
     4914
    49154915  if ( d * d != 1 ) // is d = 1 or -1? Else: error
    49164916  {
    49174917    ERROR("determinant of the given intmat has to be 1 or -1.");
    49184918  }
    4919  
     4919
    49204920  int c;
    49214921  int i,j;
     
    49304930      Ad = intAdjoint(A,i,j);
    49314931      s = 1;
    4932    
     4932
    49334933      if ( ((i + j) % 2) > 0 )
    49344934      {
     
    50055005{
    50065006  "EXAMPLE"; echo = 2;
    5007  
     5007
    50085008  intmat A[2][3] =
    50095009    1,3,5,
     
    50275027  int m = nrows(P);
    50285028  int n = ncols(P);
    5029  
     5029
    50305030  if ( m == n )
    50315031  {
    5032     intmat U = intInverse(P); 
     5032    intmat U = intInverse(P);
    50335033  }
    50345034  else
    50355035  {
    50365036    intmat U = (hermiteNormalForm(P, "transform"))[2];
    5037    
     5037
    50385038    // delete columns m+1 to n
    50395039    intmat Udel[nrows(U)][ncols(U) - (n - m)];
    50405040    int k;
    50415041    int z;
    5042      
     5042
    50435043    for ( k = 1; k <= nrows(U); k++ )
    50445044    {
     
    50485048      }
    50495049    }
    5050    
    5051     U = Udel; 
     5050
     5051    U = Udel;
    50525052  }
    50535053
     
    50625062    2,4,5,7;
    50635063
    5064   // should be a matrix with two columns 
     5064  // should be a matrix with two columns
    50655065  // for example: [−2, 1, 0, 0], [3, −3, 0, 1]
    50665066  intmat U = integralSection(P);
     
    50695069  print(P * U);
    50705070
    5071   kill U; 
     5071  kill U;
    50725072}
    50735073
     
    50875087  intmat L2 = H[2];
    50885088
    5089   // check whether G,H are subgroups of a common group, i.e. whether L1 and L2 span the same lattice 
     5089  // check whether G,H are subgroups of a common group, i.e. whether L1 and L2 span the same lattice
    50905090  if ( !isSublattice(L1,L2) || !isSublattice(L2,L1))
    50915091  {
     
    51095109{
    51105110  "EXAMPLE"; echo = 2;
    5111  
     5111
    51125112  intmat S1[2][2] =
    51135113    1,0,
     
    51175117    0;
    51185118
    5119   intmat S2[2][1] = 
     5119  intmat S2[2][1] =
    51205120    1,
    51215121    0;
     
    51315131
    51325132  kill G,H,N,S1,L1,S2,L2;
    5133    
     5133
    51345134}
    51355135
     
    51525152
    51535153  // concatinate matrices to get S
    5154   intmat A = concatintmat(S1,OS1); 
    5155   intmat B = concatintmat(OS2,S2); 
     5154  intmat A = concatintmat(S1,OS1);
     5155  intmat B = concatintmat(OS2,S2);
    51565156  intmat At = transpose(A);
    51575157  intmat Bt = transpose(B);
     
    51605160
    51615161  // concatinate matrices to get L
    5162   intmat C = concatintmat(L1,OL1); 
    5163   intmat D = concatintmat(OL2,L2); 
     5162  intmat C = concatintmat(L1,OL1);
     5163  intmat D = concatintmat(OL2,L2);
    51645164  intmat Ct = transpose(C);
    51655165  intmat Dt = transpose(D);
     
    51755175{
    51765176  "EXAMPLE"; echo = 2;
    5177  
     5177
    51785178  intmat S1[2][2] =
    51795179    1,0,
     
    51835183    0;
    51845184
    5185   intmat S2[2][2] = 
     5185  intmat S2[2][2] =
    51865186    1,0,
    51875187    0,2;
     
    51975197
    51985198  kill G,H,N,S1,L1,S2,L2;
    5199    
     5199
    52005200}
    52015201
     
    52135213  int r = intRank(V);
    52145214
    5215  
     5215
    52165216  if ( r == 0 )
    52175217  {
     
    52215221  {
    52225222    list L = smithNormalForm(V, "transform"); // L = [A,S,B] where S is the smith-NF and S = A*S*B
    5223     intmat P = intInverse(L[1]); 
     5223    intmat P = intInverse(L[1]);
    52245224
    52255225    print(L);
    5226    
    5227     if ( r < m ) 
     5226
     5227    if ( r < m )
    52285228    {
    52295229      // delete columns r+1 to m in P:
     
    52485248{
    52495249  "EXAMPLE"; echo = 2;
    5250  
     5250
    52515251  intmat V[2][4] =
    52525252    1,4,
Note: See TracChangeset for help on using the changeset viewer.