Changeset 34a9eb1 in git


Ignore:
Timestamp:
Jul 31, 2001, 7:37:43 PM (22 years ago)
Author:
Mathias Schulze <mschulze@…>
Branches:
(u'spielwiese', '91fdef05f09f54b8d58d92a472e9c4a43aa4656f')
Children:
1418c4ecec6f1dda76f77dca9a541a4f29998c7e
Parents:
fa01b75313d60b772595052a5fd57fb8c2ec49e2
Message:
*mschulze: procedure to compute spectral pairs


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

Legend:

Unmodified
Added
Removed

  • Singular/LIB/gaussman.lib

    rfa01b7 r34a9eb1  
    11///////////////////////////////////////////////////////////////////////////////
    2 version="$Id: gaussman.lib,v 1.45 2001-06-21 14:59:10 mschulze Exp $";
     2version="$Id: gaussman.lib,v 1.46 2001-07-31 17:37:43 mschulze Exp $";
    33category="Singularities";
    44
    … …  
    1515 gmsnf(p,K[,Kmax]);         Gauss-Manin system normal form
    1616 gmscoeffs(p,K[,Kmax]);     Gauss-Manin system basis representation
    17  monodromy(f[,opt]);        monodromy matrix, spectrum of monodromy of f
    18  vfiltration(f[,opt]);      V-filtration on H''/H', singularity spectrum of f
    19  spectrum(f);               singularity spectrum of f
     17 monodromy(f[,opt]);        monodromy matrix or spectrum of monodromy of f
     18 vfilt(f[,opt]);            V-filtration on H''/H' or spectrum of f
    2019 endfilt(poly f,list V);    endomorphism filtration of filtration V
    21  spprint(list S);           print spectrum S
    22  spadd(list S1,list S2);    sum of spectra S1 and S2
    23  spsub(list S1,list S2);    difference of spectra S1 and S2
    24  spmul(list S,int k);       product of spectrum S and integer k
    25  spmul(list S,intvec k);    linear combination of spectra S with coefficients k
    26  spissemicont(list S[,opt]);        test spectrum S for semicontinuity
    27  spsemicont(list S0,list S[,opt]);  relative semicontinuity of spectra S0 and S
    28  spmilnor(list S);          milnor number of spectrum S
    29  spgeomgenus(list S);       geometrical genus of spectrum S
    30  spgamma(list S);           gamma invariant of spectrum S
     20 spectrum(f);               spectrum of f
     21 sppairs(f[,opt]);          spectral pairs or spectrum of f
     22 spgen(a);                  generate spectrum defined by a
     23 sppgen(a,w);               generate spectral pairs defined by a and w
     24 spprint(list Sp);          print spectrum or spectral pairs Sp
     25 spadd(list Sp1,list Sp2);  sum of spectra Sp1 and Sp2
     26 spsub(list Sp1,list Sp2);  difference of spectra Sp1 and Sp2
     27 spmul(list Sp,int k);      product of spectrum Sp and integer k
     28 spmul(list Sp,intvec k);   linear combination of spectra Sp with coeffs k
     29 spissemicont(list Sp[,opt]);         test spectrum Sp for semicontinuity
     30 spsemicont(list Sp0,list Sp[,opt]);  semicontinuity of spectra Sp0 and Sp
     31 spmilnor(list Sp);         milnor number of spectrum Sp
     32 spgeomgenus(list Sp);      geometrical genus of spectrum Sp
     33 spgamma(list Sp);          gamma invariant of spectrum Sp
    3134
    3235SEE ALSO: mondromy_lib, spectrum_lib
    3336
    34 KEYWORDS: singularities; Gauss-Manin connection; monodromy; spectrum;
    35           Brieskorn lattice; Hodge filtration; V-filtration
     37KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice;
     38          monodromy; spectrum; spectral pairs;
     39          Hodge filtration; V-filtration; weight filtration
    3640";
    3741
    … …  
    4448  def R=basering;
    4549
    46   string s=ordstr(R);
    47   int j=find(s,",C");
     50  string os=ordstr(R);
     51  int j=find(os,",C");
    4852  if(j==0)
    4953  {
    50     j=find(s,"C,");
     54    j=find(os,"C,");
    5155  }
    5256  if(j==0)
    5357  {
    54     j=find(s,",c");
     58    j=find(os,",c");
    5559  }
    5660  if(j==0)
    5761  {
    58     j=find(s,"c,");
     62    j=find(os,"c,");
    5963  }
    6064  if(j>0)
    6165  {
    62     s[j..j+1]="  ";
    63   }
    64 
    65   execute("ring S="+charstr(R)+",(gmspoly,"+varstr(R)+"),(c,dp,"+s+");");
     66    os[j..j+1]="  ";
     67  }
     68
     69  execute("ring S="+charstr(R)+",(gmspoly,"+varstr(R)+"),(c,dp,"+os+");");
    6670
    6771  ideal I=homog(imap(R,I),gmspoly);
    … …  
    103107proc gmsring(poly t,string s)
    104108"USAGE:    gmsring(f,s); poly f, string s;
    105 ASSUME:   basering has characteristic 0 and local degree ordering,
    106           f has isolated singularity at 0
     109ASSUME:   basering with characteristic 0 and local degree ordering,
     110          f with isolated singularity at 0
    107111RETURN:
    108112@format
    … …  
    120124@end format
    121125NOTE:     do not modify gms variables if you want to use gms procedures
    122 KEYWORDS: singularities; Gauss-Manin connection
     126KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice
    123127EXAMPLE:  example gms; shows examples
    124128"
    … …  
    163167  }
    164168
     169  string os=ordstr(R);
     170  int j=find(os,",C");
     171  if(j==0)
     172  {
     173    j=find(os,"C,");
     174  }
     175  if(j==0)
     176  {
     177    j=find(os,",c");
     178  }
     179  if(j==0)
     180  {
     181    j=find(os,"c,");
     182  }
     183  if(j>0)
     184  {
     185    os[j..j+1]="  ";
     186  }
     187
    165188  execute("ring G="+string(charstr(R))+",("+s+","+varstr(R)+"),(ws("+
    166     string(deg(highcorner(g))+2*gmsmaxweight)+"),"+ordstr(R)+");");
     189    string(deg(highcorner(g))+2*gmsmaxweight)+"),"+os+",c);");
    167190
    168191  poly gmspoly=imap(R,t);
    … …  
    202225NOTE:     by setting p=l[2] the computation can be continued up to order
    203226          at most Kmax, by default Kmax=infinity
    204 KEYWORDS: singularities; Gauss-Manin connection
     227KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice
    205228EXAMPLE:  example gmsnf; shows examples
    206229"
    … …  
    305328NOTE:     by setting p=l[2] the computation can be continued up to order
    306329          at most Kmax, by default Kmax=infinity
    307 KEYWORDS: singularities; Gauss-Manin connection
     330KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice
    308331EXAMPLE:  example gmscoeffs; shows examples
    309332"
    … …  
    387410proc monodromy(poly f,list #)
    388411"USAGE:    monodromy(f[,opt]); poly f, int opt
    389 ASSUME:   basering has characteristic 0 and local degree ordering,
    390           f has isolated singularity at 0
     412ASSUME:   basering with characteristic 0 and local degree ordering,
     413          f with isolated singularity at 0
    391414RETURN:
    392415@format
    393 if opt==0:
    394   matrix M: exp(-2*pi*i*M) is a monodromy matrix of f
    395 if opt==1:
    396   ideal e: exp(-2*pi*i*e) are the eigenvalues of the monodromy of f
     416if opt=0:
     417list l:
     418  ideal l[1]: exp(-2*pi*i*l[1]) are the eigenvalues of the monodromy
     419if opt=1:
     420list l: Jordan data jordan(M) of a monodromy matrix exp(-2*pi*i*M)
    397421default: opt=1
    398422@end format
    399423SEE ALSO: mondromy_lib
    400 KEYWORDS: singularities; Gauss-Manin connection; monodromy
     424KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice; monodromy
    401425EXAMPLE:  example monodromy; shows examples
    402426"
    … …  
    412436
    413437  def R=basering;
     438  int n=nvars(R)-1;
    414439  def G=gmsring(f,"s");
    415440  setring G;
    416 
    417   int n=nvars(R)-1;
    418441  int mu,modm=ncols(gmsbasis),maxorddif(gmsbasis);
    419   ideal w=gmspoly*gmsbasis;
     442  ideal r=gmspoly*gmsbasis;
    420443  list l;
    421444  matrix U=freemodule(mu);
    422   matrix A0[mu][mu],A,C;
     445  matrix A0[mu][mu],C;
    423446  module H,dH=freemodule(mu),freemodule(mu);
    424447  module H0;
    425   int sdH=1;
    426448  int k=-1;
    427 
    428   while(sdH>0)
     449  while(size(reduce(H,std(H0*s)))>0)
    429450  {
    430451    k++;
    … …  
    433454    if(opt==0)
    434455    {
    435       l=gmscoeffs(w,k,mu);
     456      l=gmscoeffs(r,k,mu);
    436457    }
    437458    else
    438459    {
    439       l=gmscoeffs(w,k,mu+n);
    440     }
    441     C,w=l[1..2];
     460      l=gmscoeffs(r,k,mu+n);
     461    }
     462    C,r=l[1..2];
    442463    A0=A0+C;
    443464
    … …  
    448469
    449470    dbprint(printlevel-voice+2,"//gaussman::monodromy: test dH==0");
    450     sdH=size(reduce(H,std(H0*s)));
    451   }
    452 
    453   A0=A0-s^k;
     471  }
     472  A0=A0-k*s;
    454473
    455474  dbprint(printlevel-voice+2,
    456475    "//gaussman::monodromy: compute basis of saturation");
    457   H=minbase(H0);
     476  H=std(H0);
    458477  int modH=maxorddif(H)/deg(s);
    459478  dbprint(printlevel-voice+2,"//gaussman::monodromy: k="+string(modH+1));
    … …  
    461480  if(opt==0)
    462481  {
    463     l=gmscoeffs(w,modH+1,modH+1);
     482    l=gmscoeffs(r,modH+1,modH+1);
    464483  }
    465484  else
    466485  {
    467     l=gmscoeffs(w,modH+1,modH+1+n);
    468   }
    469   C,w=l[1..2];
     486    l=gmscoeffs(r,modH+1,modH+1+n);
     487  }
     488  C,r=l[1..2];
    470489  A0=A0+C;
    471 
    472490  dbprint(printlevel-voice+2,
    473491    "//gaussman::monodromy: compute A on saturation");
    474492  l=division(H*s,A0*H+s^2*diff(matrix(H),s));
    475   A=jet(l[1],l[2],0);
     493  matrix A=jet(l[1],l[2],0);
    476494
    477495  dbprint(printlevel-voice+2,
    … …  
    488506    setring(R);
    489507    ideal e=imap(G,e);
    490     return(e);
     508    return(list(e));
    491509  }
    492510
    493511  int mide=maxintdif(e);
    494 
    495512  if(mide>0)
    496513  {
    … …  
    498515      "//gaussman::monodromy: k="+string(modH+1+mide));
    499516    dbprint(printlevel-voice+2,"//gaussman::monodromy: compute C");
    500     l=gmscoeffs(w,modH+1+mide,modH+1+mide);
    501     C,w=l[1..2];
     517    l=gmscoeffs(r,modH+1+mide,modH+1+mide);
     518    C,r=l[1..2];
    502519    A0=A0+C;
     520    l=division(H*s,A0*H+s^2*diff(matrix(H),s));
     521    A=jet(l[1],l[2],mide);
    503522
    504523    intvec ide;
    … …  
    510529      {
    511530        k=int(e[j]-e[i]);
    512         if(k>ide[i])
    513         {
    514           ide[i]=k;
    515         }
    516         if(-k>ide[j])
    517         {
    518           ide[j]=-k;
     531        if(k>ide[j])
     532        {
     533          ide[j]=k;
     534        }
     535        if(-k>ide[i])
     536        {
     537          ide[i]=-k;
    519538        }
    520539      }
    … …  
    528547      }
    529548    }
    530   }
    531   while(mide>0)
    532   {
    533     dbprint(printlevel-voice+2,"//gaussman::monodromy: mide="+string(mide));
    534 
    535     U=0;
    536     A0=jet(A,0);
    537     for(i=ncols(e);i>=1;i--)
    538     {
    539       U=syz(power(A0-e[i],n+1))+U;
    540     }
    541     A=division(U,A*U)[1];
    542 
    543     for(i=mu;i>=1;i--)
    544     {
    545       for(j=mu;j>=1;j--)
    546       {
    547         if(ide[i]==0&&ide[j]>0)
    548         {
    549           A[i,j]=A[i,j]*s;
    550         }
    551         else
    552         {
    553           if(ide[i]>0&&ide[j]==0)
     549
     550    while(mide>0)
     551    {
     552      dbprint(printlevel-voice+2,"//gaussman::monodromy: mide="+string(mide));
     553
     554      A0=jet(A,0);
     555      U=0;
     556      for(i=ncols(e);i>=1;i--)
     557      {
     558        U=syz(power(A0-e[i],n+1))+U;
     559      }
     560      A=lift(U,A*U);
     561
     562      for(i=mu;i>=1;i--)
     563      {
     564        for(j=mu;j>=1;j--)
     565        {
     566          if(ide[i]==0&&ide[j]>0)
    554567          {
    555568            A[i,j]=A[i,j]/s;
    556569          }
    557         }
    558       }
    559     }
    560     for(i=mu;i>=1;i--)
    561     {
    562       if(ide[i]>0)
    563       {
    564         A[i,i]=A[i,i]+1;
    565         e[i]=e[i]+1;
    566         ide[i]=ide[i]-1;
    567       }
    568     }
    569     mide--;
    570   }
    571   A=jet(A,0);
     570          else
     571          {
     572            if(ide[i]>0&&ide[j]==0)
     573            {
     574              A[i,j]=A[i,j]*s;
     575            }
     576          }
     577        }
     578      }
     579      for(i=mu;i>=1;i--)
     580      {
     581        if(ide[i]>0)
     582        {
     583          A[i,i]=A[i,i]-1;
     584          ide[i]=ide[i]-1;
     585        }
     586      }
     587      mide--;
     588    }
     589
     590    A=jet(A,0);
     591  }
    572592
    573593  setring(R);
    574594  matrix A=imap(G,A);
    575   return(A);
     595  return(jordan(A));
    576596}
    577597example
    … …  
    583603///////////////////////////////////////////////////////////////////////////////
    584604
    585 proc vfiltration(poly f,list #)
    586 "USAGE:    vfiltration(f[,opt]); poly f, int opt
    587 ASSUME:   basering has characteristic 0 and local degree ordering,
    588           f has isolated singularity at 0
     605proc vfilt(poly f,list #)
     606"USAGE:    vfilt(f[,opt]); poly f, int opt
     607ASSUME:   basering with characteristic 0 and local degree ordering,
     608          f with isolated singularity at 0
    589609RETURN:
    590610@format
    591611list V: V-filtration of f on H''/H'
    592 if opt==0 or opt==1:
     612if opt=0 or opt=1:
    593613  intvec V[1]: numerators of spectral numbers
    594614  intvec V[2]: denominators of spectral numbers
    595615  intvec V[3]:
    596616    int V[3][i]: multiplicity of spectral number V[1][i]/V[2][i]
    597 if opt==1:
     617if opt=1:
    598618  list V[4]:
    599619    module V[4][i]: vector space basis of V[1][i]/V[2][i]-th graded part
    … …  
    604624NOTE:     H' and H'' denote the Brieskorn lattices
    605625SEE ALSO: spectrum_lib
    606 KEYWORDS: singularities; Gauss-Manin connection;
    607           Brieskorn lattice; Hodge filtration; V-filtration; spectrum
    608 EXAMPLE:  example vfiltration; shows examples
     626KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice;
     627          Hodge filtration; V-filtration; spectrum
     628EXAMPLE:  example vfilt; shows examples
    609629"
    610630{
    … …  
    624644  int n=nvars(R)-1;
    625645  int mu,modm=ncols(gmsbasis),maxorddif(gmsbasis);
    626   ideal w=gmspoly*gmsbasis;
     646  ideal r=gmspoly*gmsbasis;
    627647  list l;
    628648  matrix U=freemodule(mu);
    … …  
    637657  {
    638658    k++;
    639     dbprint(printlevel-voice+2,"//gaussman::vfiltration: k="+string(k));
    640     dbprint(printlevel-voice+2,"//gaussman::vfiltration: compute C");
    641     l=gmscoeffs(w,k);
    642     C,w=l[1..2];
     659    dbprint(printlevel-voice+2,"//gaussman::vfilt: k="+string(k));
     660    dbprint(printlevel-voice+2,"//gaussman::vfilt: compute C");
     661    l=gmscoeffs(r,k);
     662    C,r=l[1..2];
    643663    A=A+C;
    644664
    645665    H0=H;
    646     dbprint(printlevel-voice+2,"//gaussman::vfiltration: compute dH");
     666    dbprint(printlevel-voice+2,"//gaussman::vfilt: compute dH");
    647667    dH=jet(module(A*dH+s^2*diff(matrix(dH),s)),k);
    648668    H=H*s+dH;
    649669
    650     dbprint(printlevel-voice+2,"//gaussman::vfiltration: test dH==0");
     670    dbprint(printlevel-voice+2,"//gaussman::vfilt: test dH==0");
    651671    sdH=size(reduce(H,std(H0*s)));
    652672  }
    653 
    654   A=A-s^k;
    655 
    656   dbprint(printlevel-voice+2,
    657     "//gaussman::vfiltration: compute basis of saturation");
     673  A=A-k*s;
     674
     675  dbprint(printlevel-voice+2,"//gaussman::vfilt: compute basis of saturation");
    658676  H=minbase(H0);
    659677  int modH=maxorddif(H)/deg(s);
    660678  dbprint(printlevel-voice+2,"//gaussman::monodromy: k="+string(N+modH));
    661679  dbprint(printlevel-voice+2,"//gaussman::monodromy: compute C");
    662   l=gmscoeffs(w,N+modH,N+modH);
    663   C,w=l[1..2];
     680  l=gmscoeffs(r,N+modH,N+modH);
     681  C,r=l[1..2];
    664682  A=A+C;
    665683
    666   dbprint(printlevel-voice+2,
    667     "//gaussman::vfiltration: transform H0 to saturation");
     684  dbprint(printlevel-voice+2,"//gaussman::vfilt: transform H0 to saturation");
    668685  l=division(H,freemodule(mu)*s^k);
    669686  H0=jet(l[1],l[2],N-1);
    670687
    671688  dbprint(printlevel-voice+2,
    672     "//gaussman::vfiltration: compute H0 as vector space V0");
     689    "//gaussman::vfilt: compute H0 as vector space V0");
    673690  dbprint(printlevel-voice+2,
    674     "//gaussman::vfiltration: compute H1 as vector space V1");
     691    "//gaussman::vfilt: compute H1 as vector space V1");
    675692  poly p;
    676693  int i0,j0,i1,j1;
    … …  
    699716  }
    700717
    701   dbprint(printlevel-voice+2,
    702     "//gaussman::vfiltration: compute A on saturation");
     718  dbprint(printlevel-voice+2,"//gaussman::vfilt: compute A on saturation");
    703719  l=division(H*s,A*H+s^2*diff(matrix(H),s));
    704720  A=jet(l[1],l[2],N-1);
    705721
    706   dbprint(printlevel-voice+2,"//gaussman::vfiltration: compute matrix M of A");
     722  dbprint(printlevel-voice+2,"//gaussman::vfilt: compute matrix M of A");
    707723  matrix M[mu*N][mu*N];
    708724  for(i0=mu;i0>=1;i0--)
    … …  
    730746  }
    731747
    732   dbprint(printlevel-voice+2,
    733     "//gaussman::vfiltration: compute eigenvalues eA of A");
     748  dbprint(printlevel-voice+2,"//gaussman::vfilt: compute eigenvalues eA of A");
    734749  ideal eA=eigenval(jet(A,0))[1];
    735   dbprint(printlevel-voice+2,"//gaussman::vfiltration: eA="+string(eA));
    736 
    737   dbprint(printlevel-voice+2,
    738     "//gaussman::vfiltration: compute eigenvalues eM of M");
     750  dbprint(printlevel-voice+2,"//gaussman::vfilt: eA="+string(eA));
     751
     752  dbprint(printlevel-voice+2,"//gaussman::vfilt: compute eigenvalues eM of M");
    739753  ideal eM;
    740754  k=0;
    … …  
    762776    }
    763777  }
    764   dbprint(printlevel-voice+2,"//gaussman::vfiltration: eM="+string(eM));
     778  dbprint(printlevel-voice+2,"//gaussman::vfilt: eM="+string(eM));
    765779
    766780  dbprint(printlevel-voice+2,
    767     "//gaussman::vfiltration: compute V-filtration on H0/H1");
     781    "//gaussman::vfilt: compute V-filtration on H0/H1");
    768782  ideal a;
    769783  k=0;
    … …  
    776790    {
    777791      dbprint(printlevel-voice+2,
    778         "//gaussman::vfiltration: compute V["+string(i)+"]");
     792        "//gaussman::vfilt: compute V["+string(i)+"]");
    779793      V1=module(V1)+syz(power(M-eM[i],n+1));
    780794      V[i]=interred(intersect(V1,V0));
    … …  
    789803
    790804    dbprint(printlevel-voice+2,
    791       "//gaussman::vfiltration: symmetry index found");
     805      "//gaussman::vfilt: symmetry index found");
    792806    int j=k;
    793807
    … …  
    795809    {
    796810      dbprint(printlevel-voice+2,
    797         "//gaussman::vfiltration: compute V["+string(i)+"]");
     811        "//gaussman::vfilt: compute V["+string(i)+"]");
    798812      V1=module(V1)+syz(power(M-eM[i],n+1));
    799813      V[i]=interred(intersect(V1,V0));
    … …  
    807821    }
    808822
    809     dbprint(printlevel-voice+2,"//gaussman::vfiltration: apply symmetry");
     823    dbprint(printlevel-voice+2,"//gaussman::vfilt: apply symmetry");
    810824    while(j>=1)
    811825    {
    … …  
    828842    {
    829843      dbprint(printlevel-voice+2,
    830         "//gaussman::vfiltration: compute V["+string(i)+"]");
     844        "//gaussman::vfilt: compute V["+string(i)+"]");
    831845      V1=module(V1)+syz(power(M-eM[i],n+1));
    832846      V[i]=interred(intersect(V1,V0));
    … …  
    839853          a[k]=eM[i]-1;
    840854          dbprint(printlevel-voice+2,
    841             "//gaussman::vfiltration: transform to V0");
     855            "//gaussman::vfilt: transform to V0");
    842856          v[k]=matrix(freemodule(ncols(V[i])),mu,mu*N)*division(V0,V[i])[1];
    843857        }
    … …  
    860874          v[k]=v[j];
    861875          dbprint(printlevel-voice+2,
    862             "//gaussman::vfiltration: transform to V0");
     876            "//gaussman::vfilt: transform to V0");
    863877          v[j]=matrix(freemodule(ncols(V[i])),mu,mu*N)*division(V0,V[i])[1];
    864878          j--;
    … …  
    867881    }
    868882
    869     dbprint(printlevel-voice+2,
    870       "//gaussman::vfiltration: compute graded parts");
     883    dbprint(printlevel-voice+2,"//gaussman::vfilt: compute graded parts");
    871884    for(k=1;k<size(v);k++)
    872885    {
    … …  
    888901  ring R=0,(x,y),ds;
    889902  poly f=x5+x2y2+y5;
    890   vfiltration(f);
    891 }
    892 ///////////////////////////////////////////////////////////////////////////////
    893 
    894 proc spectrum(poly f)
    895 "USAGE:    spectrum(f); poly f
    896 ASSUME:   basering has characteristic 0 and local degree ordering,
    897           f has isolated singularity at 0
    898 RETURN:
    899 @format
    900 list S: singularity spectrum of f
    901   ideal S[1]: spectral numbers in increasing order
    902   intvec S[2]:
    903     int S[2][i]: multiplicity of spectral number S[1][i]
    904 @end format
    905 SEE ALSO: spectrum_lib
    906 KEYWORDS: singularities; Gauss-Manin connection; spectrum
    907 EXAMPLE:  example spectrum; shows examples
    908 "
    909 {
    910   return(vfiltration(f,0));
    911 }
    912 example
    913 { "EXAMPLE:"; echo=2;
    914   ring R=0,(x,y),ds;
    915   poly f=x5+x2y2+y5;
    916   spprint(spectrum(f));
     903  vfilt(f);
    917904}
    918905///////////////////////////////////////////////////////////////////////////////
    … …  
    920907proc endfilt(poly f,list V)
    921908"USAGE:   endfilt(f,V); poly f, list V
    922 ASSUME:  basering has characteristic 0 and local degree ordering,
    923          f has isolated singularity at 0
     909ASSUME:  basering with characteristic 0 and local degree ordering,
     910         f with isolated singularity at 0
    924911RETURN:
    925912@format
    … …  
    934921@end format
    935922SEE ALSO: spectrum_lib
    936 KEYWORDS: singularities; Gauss-Manin connection; spectrum;
    937           Brieskorn lattice; Hodge filtration; V-filtration
     923KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice; spectrum;
     924          Hodge filtration; V-filtration
    938925EXAMPLE: example endfilt; shows examples
    939926"
    … …  
    11001087  ring R=0,(x,y),ds;
    11011088  poly f=x5+x2y2+y5;
    1102   endfilt(f,vfiltration(f));
    1103 }
    1104 ///////////////////////////////////////////////////////////////////////////////
    1105 
    1106 proc spprint(list S)
    1107 "USAGE:   spprint(S); list S
    1108 RETURN:  string: spectrum S
    1109 EXAMPLE: example spprint; shows examples
     1089  endfilt(f,vfilt(f));
     1090}
     1091///////////////////////////////////////////////////////////////////////////////
     1092
     1093proc spectrum(poly f)
     1094"USAGE:    spectrum(f); poly f
     1095ASSUME:   basering with characteristic 0 and local degree ordering,
     1096          f with isolated singularity at 0
     1097RETURN:
     1098@format
     1099list Sp: spectrum of f
     1100  ideal Sp[1]: spectral numbers in increasing order
     1101  intvec Sp[2]:
     1102    int Sp[2][i]: multiplicity of spectral number Sp[1][i]
     1103@end format
     1104SEE ALSO: spectrum_lib
     1105KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice; spectrum
     1106EXAMPLE:  example spnumbers; shows examples
    11101107"
    11111108{
    1112   string s;
    1113   for(int i=1;i<size(S[2]);i++)
    1114   {
    1115     s=s+"("+string(S[1][i])+","+string(S[2][i])+"),";
    1116   }
    1117   s=s+"("+string(S[1][i])+","+string(S[2][i])+")";
    1118   return(s);
     1109  return(sppairs(f,0));
    11191110}
    11201111example
    11211112{ "EXAMPLE:"; echo=2;
    11221113  ring R=0,(x,y),ds;
    1123   list S=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1124   spprint(S);
    1125 }
    1126 ///////////////////////////////////////////////////////////////////////////////
    1127 
    1128 proc spadd(list S1,list S2)
    1129 "USAGE:   spadd(S1,S2); list S1,S2
    1130 RETURN:  list: sum of spectra S1 and S2
     1114  poly f=x5+x2y2+y5;
     1115  spprint(spectrum(f));
     1116}
     1117///////////////////////////////////////////////////////////////////////////////
     1118
     1119proc sppairs(poly f,list #)
     1120"USAGE:    sppairs(f[,opt]); poly f, int opt
     1121ASSUME:   basering with characteristic 0 and local degree ordering,
     1122          f with isolated singularity at 0
     1123RETURN:
     1124@format
     1125if opt=0
     1126list Sp: spectrum of f
     1127  ideal Sp[1]: spectral numbers in increasing order
     1128  intvec Sp[2]: corresponding multiplicities of spectral numbers
     1129if opt=1:
     1130list Spp: spectral pairs of f
     1131  ideal Spp[1]: spectral numbers in increasing order
     1132  intvec Spp[2]: corresponding nilpotency indices in increasing order
     1133  intvec Spp[3]: corresponding multiplicities of spectral pairs
     1134default: opt=1
     1135@end format
     1136SEE ALSO: spectrum_lib
     1137KEYWORDS: singularities; Gauss-Manin connection; Brieskorn lattice;
     1138          spectrum; spectral pairs
     1139EXAMPLE:  example sppairs; shows examples
     1140"
     1141{
     1142  int opt=1;
     1143  if(size(#)>0)
     1144  {
     1145    if(typeof(#[1])=="int")
     1146    {
     1147      opt=#[1];
     1148    }
     1149  }
     1150
     1151  def R=basering;
     1152  int n=nvars(R)-1;
     1153  def G=gmsring(f,"s");
     1154  setring(G);
     1155  int mu,modm=ncols(gmsbasis),maxorddif(gmsbasis);
     1156  ideal r=gmspoly*gmsbasis;
     1157  list l;
     1158  matrix U=freemodule(mu);
     1159  matrix A0[mu][mu],C;
     1160  module H0;
     1161  module H,dH=freemodule(mu),freemodule(mu);
     1162  int k=-1;
     1163  while(size(reduce(H,std(H0*s)))>0)
     1164  {
     1165    k++;
     1166    l=gmscoeffs(r,k,mu+n);
     1167    C,r=l[1..2];
     1168    A0=A0+C;
     1169    H0=H;
     1170    dH=jet(module(A0*dH+s^2*diff(matrix(dH),s)),k);
     1171    H=H*s+dH;
     1172  }
     1173  A0=A0-k*s;
     1174
     1175  H=std(H0);
     1176  int modH=maxorddif(H)/deg(s);
     1177  l=division(H,freemodule(mu)*s^k);
     1178  H0=l[1];
     1179  l=gmscoeffs(r,modH+1,modH+1+n);
     1180  C,r=l[1..2];
     1181  A0=A0+C;
     1182  l=division(H*s,A0*H+s^2*diff(matrix(H),s));
     1183  matrix A=jet(l[1],l[2],0);
     1184
     1185  int i,j;
     1186  matrix V;
     1187  if(!opt)
     1188  {
     1189    U=jordanbasis(A);
     1190    V=lift(U,freemodule(mu));
     1191    A=V*A*U;
     1192    H0=std(V*H0);
     1193    ideal a;
     1194    for(i=1;i<=mu;i++)
     1195    {
     1196      j=leadexp(H0[i])[nvars(basering)+1];
     1197      a[i]=A[j,j]+deg(lead(H0[i]))/deg(s)-1;
     1198    }
     1199    setring(R);
     1200    return(spgen(imap(G,a)));
     1201  }
     1202
     1203  l=eigenval(A);
     1204  def e,b=l[1..2];
     1205  int mide=maxintdif(e);
     1206  if(mide>0)
     1207  {
     1208    l=gmscoeffs(r,modH+1+mide,modH+1+mide);
     1209    C,r=l[1..2];
     1210    A0=A0+C;
     1211    l=division(H*s,A0*H+s^2*diff(matrix(H),s));
     1212    A=jet(l[1],l[2],mide);
     1213
     1214    intvec ide;
     1215    ide[mu]=0;
     1216    for(i=ncols(e);i>=1;i--)
     1217    {
     1218      for(j=i-1;j>=1;j--)
     1219      {
     1220        k=int(e[j]-e[i]);
     1221        if(k>ide[j])
     1222        {
     1223          ide[j]=k;
     1224        }
     1225        if(-k>ide[i])
     1226        {
     1227          ide[i]=-k;
     1228        }
     1229      }
     1230    }
     1231    for(j,k=ncols(e),mu;j>=1;j--)
     1232    {
     1233      for(i=b[j];i>=1;i--)
     1234      {
     1235        ide[k]=ide[j];
     1236        k--;
     1237      }
     1238    }
     1239
     1240    while(mide>0)
     1241    {
     1242      A0=jet(A,0);
     1243      U=0;
     1244      for(i=ncols(e);i>=1;i--)
     1245      {
     1246        U=syz(power(A0-e[i],n+1))+U;
     1247      }
     1248      V=lift(U,freemodule(mu));
     1249      A=V*A*U;
     1250      H0=V*H0;
     1251
     1252      for(i=mu;i>=1;i--)
     1253      {
     1254        for(j=mu;j>=1;j--)
     1255        {
     1256          if(ide[i]==0&&ide[j]>0)
     1257          {
     1258            A[i,j]=A[i,j]/s;
     1259          }
     1260          else
     1261          {
     1262            if(ide[i]>0&&ide[j]==0)
     1263            {
     1264              A[i,j]=A[i,j]*s;
     1265            }
     1266          }
     1267        }
     1268      }
     1269      H0=transpose(H0);
     1270      for(i=mu;i>=1;i--)
     1271      {
     1272        if(ide[i]>0)
     1273        {
     1274          A[i,i]=A[i,i]-1;
     1275          ide[i]=ide[i]-1;
     1276          H0[i]=H0[i]*s;
     1277        }
     1278      }
     1279      H0=transpose(H0);
     1280      mide--;
     1281    }
     1282
     1283    A=jet(A,0);
     1284  }
     1285
     1286  U=jordanbasis(A);
     1287  V=lift(U,freemodule(mu));
     1288  A0=V*A*U;
     1289
     1290  intvec w;
     1291  w[mu]=1;
     1292  j=1;
     1293  for(i=mu-1;i>=1;i--)
     1294  {
     1295    j++;
     1296    if(A0[i,i+1]==0)
     1297    {
     1298      j=1;
     1299    }
     1300    w[i]=j;
     1301  }
     1302
     1303  vector u;
     1304  for(i=1;i<ncols(A);i++)
     1305  {
     1306    j=i+1;
     1307    while(j<ncols(A)&&A[i,i]==A[j,j])
     1308    {
     1309      if(w[i]<w[j])
     1310      {
     1311        k=w[i];
     1312        w[i]=w[j];
     1313        w[i]=k;
     1314        u=U[i];
     1315        U[i]=U[j];
     1316        U[j]=u;
     1317      }
     1318      j++;
     1319    }
     1320  }
     1321
     1322  V=lift(U,freemodule(mu));
     1323  A=V*A*U;
     1324  H0=std(V*H0);
     1325  ideal a;
     1326  for(i=1;i<=mu;i++)
     1327  {
     1328    j=leadexp(H0[i])[nvars(basering)+1];
     1329    a[i]=A[j,j]+deg(lead(H0[i]))/deg(s)-1;
     1330  }
     1331  setring(R);
     1332  return(sppgen(imap(G,a),w));
     1333}
     1334example
     1335{ "EXAMPLE:"; echo=2;
     1336  ring R=0,(x,y),ds;
     1337  poly f=x5+x2y2+y5;
     1338  spprint(sppairs(f));
     1339}
     1340///////////////////////////////////////////////////////////////////////////////
     1341
     1342proc spgen(ideal a)
     1343"USAGE:   spgen(a); ideal a
     1344RETURN:
     1345@format
     1346list Sp: numbers in a with multiplicities
     1347  ideal Sp[1]: numbers in a in increasing order
     1348  intvec Sp[2]: corresponding multiplicities
     1349@end format
     1350EXAMPLE: example spgen; shows examples
     1351"
     1352{
     1353  ideal a0=jet(a,0);
     1354  int i,j;
     1355  number n;
     1356  for(i=1;i<=ncols(a0);i++)
     1357  {
     1358    for(j=i+1;j<=ncols(a0);j++)
     1359    {
     1360      if(number(a0[i])>number(a0[j]))
     1361      {
     1362        n=a0[i];
     1363        a0[i]=a0[j];
     1364        a0[j]=n;
     1365      }
     1366    }
     1367  }
     1368  j=1;
     1369  a=a0[1];
     1370  intvec m=1;
     1371  for(i=2;i<=ncols(a0);i++)
     1372  {
     1373    if(a0[i]==a[j])
     1374    {
     1375      m[j]=m[j]+1;
     1376    }
     1377    else
     1378    {
     1379      j++;
     1380      a[j]=a0[i];
     1381      m[j]=1;
     1382    }
     1383  }
     1384  return(list(a,m));
     1385}
     1386example
     1387{ "EXAMPLE:"; echo=2;
     1388  ring R=0,(x,y),ds;
     1389  ideal a=-1/2,-3/10,-3/10,-1/10,-1/10,0,1/10,1/10,3/10,3/10,1/2;
     1390  spprint(spgen(a));
     1391}
     1392///////////////////////////////////////////////////////////////////////////////
     1393
     1394proc sppgen(ideal a,intvec w)
     1395"USAGE:   sppgen(a); ideal a
     1396RETURN:
     1397@format
     1398list Spp: pairs in a and w with multiplicities
     1399  ideal Spp[1]: numbers in a
     1400  intvec Spp[2]: corresponding integers in w
     1401  intvec Spp[3]: corresponding multiplicities
     1402@end format
     1403EXAMPLE: example sppgen; shows examples
     1404"
     1405{
     1406  ideal a0=jet(a,0);
     1407  intvec w0=w;
     1408  int i,j,k;
     1409  number n;
     1410  for(i=1;i<=ncols(a0);i++)
     1411  {
     1412    for(j=i+1;j<=ncols(a0);j++)
     1413    {
     1414      if(number(a0[i])>number(a0[j])||a0[i]==a0[j]&&w0[i]>w0[j])
     1415      {
     1416        n=a0[i];
     1417        a0[i]=a0[j];
     1418        a0[j]=n;
     1419        k=w0[i];
     1420        w0[i]=w0[j];
     1421        w0[j]=k;
     1422      }
     1423    }
     1424  }
     1425  j=1;
     1426  a=a0[1];
     1427  w=w0[1];
     1428  intvec m=1;
     1429  for(i=2;i<=ncols(a0);i++)
     1430  {
     1431    if(a0[i]==a[j]&&w0[i]==w[j])
     1432    {
     1433      m[j]=m[j]+1;
     1434    }
     1435    else
     1436    {
     1437      j++;
     1438      a[j]=a0[i];
     1439      w[j]=w0[i];
     1440      m[j]=1;
     1441    }
     1442  }
     1443  return(list(a,w,m));
     1444}
     1445example
     1446{ "EXAMPLE:"; echo=2;
     1447  ring R=0,(x,y),ds;
     1448  ideal a=-1/2,-3/10,-3/10,-1/10,-1/10,0,1/10,1/10,3/10,3/10,1/2;
     1449  intvec w=2,1,1,1,1,1,1,1,1,1,1;
     1450  spprint(sppgen(a,w));
     1451}
     1452///////////////////////////////////////////////////////////////////////////////
     1453
     1454proc spprint(list Sp)
     1455"USAGE:   spprint(Sp); list Sp
     1456RETURN:  string: spectrum or spectral pairs Sp
     1457EXAMPLE: example spprint; shows examples
     1458"
     1459{
     1460  string s;
     1461  if(size(Sp)==2)
     1462  {
     1463    for(int i=1;i<size(Sp[2]);i++)
     1464    {
     1465      s=s+"("+string(Sp[1][i])+","+string(Sp[2][i])+"),";
     1466    }
     1467    s=s+"("+string(Sp[1][i])+","+string(Sp[2][i])+")";
     1468  }
     1469  else
     1470  {
     1471    for(int i=1;i<size(Sp[3]);i++)
     1472    {
     1473      s=s+"(("+string(Sp[1][i])+","+string(Sp[2][i])+"),"+string(Sp[3][i])+"),";
     1474    }
     1475    s=s+"(("+string(Sp[1][i])+","+string(Sp[2][i])+"),"+string(Sp[3][i])+")";
     1476  }
     1477  return(s);
     1478}
     1479example
     1480{ "EXAMPLE:"; echo=2;
     1481  ring R=0,(x,y),ds;
     1482  list Sp=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1483  spprint(Sp);
     1484}
     1485///////////////////////////////////////////////////////////////////////////////
     1486
     1487proc spadd(list Sp1,list Sp2)
     1488"USAGE:   spadd(Sp1,Sp2); list Sp1,Sp2
     1489RETURN:  list: sum of spectra Sp1 and Sp2
    11311490EXAMPLE: example spadd; shows examples
    11321491"
    … …  
    11351494  intvec m;
    11361495  int i,i1,i2=1,1,1;
    1137   while(i1<=size(S1[2])||i2<=size(S2[2]))
    1138   {
    1139     if(i1<=size(S1[2]))
    1140     {
    1141       if(i2<=size(S2[2]))
    1142       {
    1143         if(number(S1[1][i1])<number(S2[1][i2]))
    1144         {
    1145           s[i]=S1[1][i1];
    1146           m[i]=S1[2][i1];
     1496  while(i1<=size(Sp1[2])||i2<=size(Sp2[2]))
     1497  {
     1498    if(i1<=size(Sp1[2]))
     1499    {
     1500      if(i2<=size(Sp2[2]))
     1501      {
     1502        if(number(Sp1[1][i1])<number(Sp2[1][i2]))
     1503        {
     1504          s[i]=Sp1[1][i1];
     1505          m[i]=Sp1[2][i1];
    11471506          i++;
    11481507          i1++;
    … …  
    11501509        else
    11511510        {
    1152           if(number(S1[1][i1])>number(S2[1][i2]))
     1511          if(number(Sp1[1][i1])>number(Sp2[1][i2]))
    11531512          {
    1154             s[i]=S2[1][i2];
    1155             m[i]=S2[2][i2];
     1513            s[i]=Sp2[1][i2];
     1514            m[i]=Sp2[2][i2];
    11561515            i++;
    11571516            i2++;
    … …  
    11591518          else
    11601519          {
    1161             if(S1[2][i1]+S2[2][i2]!=0)
     1520            if(Sp1[2][i1]+Sp2[2][i2]!=0)
    11621521            {
    1163               s[i]=S1[1][i1];
    1164               m[i]=S1[2][i1]+S2[2][i2];
     1522              s[i]=Sp1[1][i1];
     1523              m[i]=Sp1[2][i1]+Sp2[2][i2];
    11651524              i++;
    11661525            }
    … …  
    11721531      else
    11731532      {
    1174         s[i]=S1[1][i1];
    1175         m[i]=S1[2][i1];
     1533        s[i]=Sp1[1][i1];
     1534        m[i]=Sp1[2][i1];
    11761535        i++;
    11771536        i1++;
    … …  
    11801539    else
    11811540    {
    1182       s[i]=S2[1][i2];
    1183       m[i]=S2[2][i2];
     1541      s[i]=Sp2[1][i2];
     1542      m[i]=Sp2[2][i2];
    11841543      i++;
    11851544      i2++;
    … …  
    11911550{ "EXAMPLE:"; echo=2;
    11921551  ring R=0,(x,y),ds;
    1193   list S1=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1194   spprint(S1);
    1195   list S2=list(ideal(-1/6,1/6),intvec(1,1));
    1196   spprint(S2);
    1197   spprint(spadd(S1,S2));
    1198 }
    1199 ///////////////////////////////////////////////////////////////////////////////
    1200 
    1201 proc spsub(list S1,list S2)
    1202 "USAGE:   spsub(S1,S2); list S1,S2
    1203 RETURN:  list: difference of spectra S1 and S2
     1552  list Sp1=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1553  spprint(Sp1);
     1554  list Sp2=list(ideal(-1/6,1/6),intvec(1,1));
     1555  spprint(Sp2);
     1556  spprint(spadd(Sp1,Sp2));
     1557}
     1558///////////////////////////////////////////////////////////////////////////////
     1559
     1560proc spsub(list Sp1,list Sp2)
     1561"USAGE:   spsub(Sp1,Sp2); list Sp1,Sp2
     1562RETURN:  list: difference of spectra Sp1 and Sp2
    12041563EXAMPLE: example spsub; shows examples
    12051564"
    12061565{
    1207   return(spadd(S1,spmul(S2,-1)));
     1566  return(spadd(Sp1,spmul(Sp2,-1)));
    12081567}
    12091568example
    12101569{ "EXAMPLE:"; echo=2;
    12111570  ring R=0,(x,y),ds;
    1212   list S1=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1213   spprint(S1);
    1214   list S2=list(ideal(-1/6,1/6),intvec(1,1));
    1215   spprint(S2);
    1216   spprint(spsub(S1,S2));
     1571  list Sp1=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1572  spprint(Sp1);
     1573  list Sp2=list(ideal(-1/6,1/6),intvec(1,1));
     1574  spprint(Sp2);
     1575  spprint(spsub(Sp1,Sp2));
    12171576}
    12181577///////////////////////////////////////////////////////////////////////////////
    … …  
    12211580"USAGE:
    12221581@format
    1223 1) spmul(S,k); list S, int k
    1224 2) spmul(S,k); list S, intvec k
     15821) spmul(Sp,k); list Sp, int k
     15832) spmul(Sp,k); list Sp, intvec k
    12251584@end format
    12261585RETURN:
    12271586@format
    1228 1) list: product of spectrum S and integer k
    1229 2) list: linear combination of spectra S with coefficients k
     15871) list: product of spectrum Sp and integer k
     15882) list: linear combination of spectra Sp with coefficients k
    12301589@end format
    12311590EXAMPLE: example spmul; shows examples
    … …  
    12421601      if(typeof(#[2])=="intvec")
    12431602      {
    1244         list S0=list(ideal(),intvec(0));
     1603        list Sp0=list(ideal(),intvec(0));
    12451604        for(int i=size(#[2]);i>=1;i--)
    12461605        {
    1247           S0=spadd(S0,spmul(#[1][i],#[2][i]));
    1248         }
    1249         return(S0);
     1606          Sp0=spadd(Sp0,spmul(#[1][i],#[2][i]));
     1607        }
     1608        return(Sp0);
    12501609      }
    12511610    }
    … …  
    12561615{ "EXAMPLE:"; echo=2;
    12571616  ring R=0,(x,y),ds;
    1258   list S=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1259   spprint(S);
    1260   spprint(spmul(S,2));
    1261   list S1=list(ideal(-1/6,1/6),intvec(1,1));
    1262   spprint(S1);
    1263   list S2=list(ideal(-1/3,0,1/3),intvec(1,2,1));
    1264   spprint(S2);
    1265   spprint(spmul(list(S1,S2),intvec(1,2)));
    1266 }
    1267 ///////////////////////////////////////////////////////////////////////////////
    1268 
    1269 proc spissemicont(list S,list #)
    1270 "USAGE:   spissemicont(S[,opt]); list S, int opt
     1617  list Sp=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1618  spprint(Sp);
     1619  spprint(spmul(Sp,2));
     1620  list Sp1=list(ideal(-1/6,1/6),intvec(1,1));
     1621  spprint(Sp1);
     1622  list Sp2=list(ideal(-1/3,0,1/3),intvec(1,2,1));
     1623  spprint(Sp2);
     1624  spprint(spmul(list(Sp1,Sp2),intvec(1,2)));
     1625}
     1626///////////////////////////////////////////////////////////////////////////////
     1627
     1628proc spissemicont(list Sp,list #)
     1629"USAGE:   spissemicont(Sp[,opt]); list Sp, int opt
    12711630RETURN:
    12721631@format
    12731632int k=
    1274 if opt==0:
    1275   1, if sum of spectrum S over all intervals [a,a+1) is positive
    1276   0, if sum of spectrum S over some interval [a,a+1) is negative
    1277 if opt==1:
    1278   1, if sum of spectrum S over all intervals [a,a+1) and (a,a+1) is positive
    1279   0, if sum of spectrum S over some interval [a,a+1) or (a,a+1) is negative
     1633if opt=0:
     1634  1, if sum of spectrum Sp over all intervals [a,a+1) is positive
     1635  0, if sum of spectrum Sp over some interval [a,a+1) is negative
     1636if opt=1:
     1637  1, if sum of spectrum Sp over all intervals [a,a+1) and (a,a+1) is positive
     1638  0, if sum of spectrum Sp over some interval [a,a+1) or (a,a+1) is negative
    12801639default: opt=0
    12811640@end format
    … …  
    12921651  }
    12931652  int i,j,k=1,1,0;
    1294   while(j<=size(S[2]))
    1295   {
    1296     while(j+1<=size(S[2])&&S[1][j]<S[1][i]+1)
    1297     {
    1298       k=k+S[2][j];
     1653  while(j<=size(Sp[2]))
     1654  {
     1655    while(j+1<=size(Sp[2])&&Sp[1][j]<Sp[1][i]+1)
     1656    {
     1657      k=k+Sp[2][j];
    12991658      j++;
    13001659    }
    1301     if(j==size(S[2])&&S[1][j]<S[1][i]+1)
    1302     {
    1303       k=k+S[2][j];
     1660    if(j==size(Sp[2])&&Sp[1][j]<Sp[1][i]+1)
     1661    {
     1662      k=k+Sp[2][j];
    13041663      j++;
    13051664    }
    … …  
    13081667      return(0);
    13091668    }
    1310     k=k-S[2][i];
     1669    k=k-Sp[2][i];
    13111670    if(k<0&&opt==1)
    13121671    {
    … …  
    13201679{ "EXAMPLE:"; echo=2;
    13211680  ring R=0,(x,y),ds;
    1322   list S1=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1323   spprint(S1);
    1324   list S2=list(ideal(-1/6,1/6),intvec(1,1));
    1325   spprint(S2);
    1326   spissemicont(spsub(S1,spmul(S2,5)));
    1327   spissemicont(spsub(S1,spmul(S2,5)),1);
    1328   spissemicont(spsub(S1,spmul(S2,6)));
    1329 }
    1330 ///////////////////////////////////////////////////////////////////////////////
    1331 
    1332 proc spsemicont(list S0,list S,list #)
    1333 "USAGE:   spsemicont(S,k[,opt]); list S0, list S, int opt
     1681  list Sp1=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1682  spprint(Sp1);
     1683  list Sp2=list(ideal(-1/6,1/6),intvec(1,1));
     1684  spprint(Sp2);
     1685  spissemicont(spsub(Sp1,spmul(Sp2,5)));
     1686  spissemicont(spsub(Sp1,spmul(Sp2,5)),1);
     1687  spissemicont(spsub(Sp1,spmul(Sp2,6)));
     1688}
     1689///////////////////////////////////////////////////////////////////////////////
     1690
     1691proc spsemicont(list Sp0,list Sp,list #)
     1692"USAGE:   spsemicont(Sp,k[,opt]); list Sp0, list Sp, int opt
    13341693RETURN:  list of intvecs l:
    1335          spissemicont(sub(S0,spmul(S,k)),opt)==1 iff k<=l[i] for some i
    1336 NOTE:    if the spectra S occur with multiplicities k in a deformation
    1337          of the [quasihomogeneous] spectrum S0 then
    1338          spissemicont(sub(S0,spmul(S,k))[,1])==1
     1694         spissemicont(sub(Sp0,spmul(Sp,k)),opt)==1 iff k<=l[i] for some i
     1695NOTE:    if the spectra Sp occur with multiplicities k in a deformation
     1696         of the [quasihomogeneous] spectrum Sp0 then
     1697         spissemicont(sub(Sp0,spmul(Sp,k))[,1])==1
    13391698EXAMPLE: example spsemicont; shows examples
    13401699"
    … …  
    13421701  list l,l0;
    13431702  int i,j,k;
    1344   while(spissemicont(S0,#))
    1345   {
    1346     if(size(S)>1)
    1347     {
    1348       l0=spsemicont(S0,list(S[1..size(S)-1]));
     1703  while(spissemicont(Sp0,#))
     1704  {
     1705    if(size(Sp)>1)
     1706    {
     1707      l0=spsemicont(Sp0,list(Sp[1..size(Sp)-1]));
    13491708      for(i=1;i<=size(l0);i++)
    13501709      {
    … …  
    13581717          if(l[j]==l0[i])
    13591718          {
    1360             l[j][size(S)]=k;
     1719            l[j][size(Sp)]=k;
    13611720          }
    13621721          else
    13631722          {
    1364             l0[i][size(S)]=k;
     1723            l0[i][size(Sp)]=k;
    13651724            l=l+list(l0[i]);
    13661725          }
    … …  
    13721731      }
    13731732    }
    1374     S0=spsub(S0,S[size(S)]);
     1733    Sp0=spsub(Sp0,Sp[size(Sp)]);
    13751734    k++;
    13761735  }
    1377   if(size(S)>1)
     1736  if(size(Sp)>1)
    13781737  {
    13791738    return(l);
    … …  
    13871746{ "EXAMPLE:"; echo=2;
    13881747  ring R=0,(x,y),ds;
    1389   list S0=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1390   spprint(S0);
    1391   list S1=list(ideal(-1/6,1/6),intvec(1,1));
    1392   spprint(S1);
    1393   list S2=list(ideal(-1/3,0,1/3),intvec(1,2,1));
    1394   spprint(S2);
    1395   list S=S1,S2;
    1396   list l=spsemicont(S0,S);
     1748  list Sp0=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1749  spprint(Sp0);
     1750  list Sp1=list(ideal(-1/6,1/6),intvec(1,1));
     1751  spprint(Sp1);
     1752  list Sp2=list(ideal(-1/3,0,1/3),intvec(1,2,1));
     1753  spprint(Sp2);
     1754  list Sp=Sp1,Sp2;
     1755  list l=spsemicont(Sp0,Sp);
    13971756  l;
    1398   spissemicont(spsub(S0,spmul(S,l[1])));
    1399   spissemicont(spsub(S0,spmul(S,l[1]-1)));
    1400   spissemicont(spsub(S0,spmul(S,l[1]+1)));
    1401 }
    1402 ///////////////////////////////////////////////////////////////////////////////
    1403 
    1404 proc spmilnor(list S)
    1405 "USAGE:   spmilnor(S); list S
    1406 RETURN:  int: Milnor number of spectrum S
     1757  spissemicont(spsub(Sp0,spmul(Sp,l[1])));
     1758  spissemicont(spsub(Sp0,spmul(Sp,l[1]-1)));
     1759  spissemicont(spsub(Sp0,spmul(Sp,l[1]+1)));
     1760}
     1761///////////////////////////////////////////////////////////////////////////////
     1762
     1763proc spmilnor(list Sp)
     1764"USAGE:   spmilnor(Sp); list Sp
     1765RETURN:  int: Milnor number of spectrum Sp
    14071766EXAMPLE: example spmilnor; shows examples
    14081767"
    14091768{
    1410   return(sum(S[2]));
     1769  return(sum(Sp[2]));
    14111770}
    14121771example
    14131772{ "EXAMPLE:"; echo=2;
    14141773  ring R=0,(x,y),ds;
    1415   list S=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1416   spprint(S);
    1417   spmilnor(S);
    1418 }
    1419 ///////////////////////////////////////////////////////////////////////////////
    1420 
    1421 proc spgeomgenus(list S)
    1422 "USAGE:   spgeomgenus(S); list S
    1423 RETURN:  int: geometrical genus of spectrum S
     1774  list Sp=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1775  spprint(Sp);
     1776  spmilnor(Sp);
     1777}
     1778///////////////////////////////////////////////////////////////////////////////
     1779
     1780proc spgeomgenus(list Sp)
     1781"USAGE:   spgeomgenus(Sp); list Sp
     1782RETURN:  int: geometrical genus of spectrum Sp
    14241783EXAMPLE: example spgeomgenus; shows examples
    14251784"
    … …  
    14271786  int g=0;
    14281787  int i=1;
    1429   while(i+1<=size(S[2])&&number(S[1][i])<=number(0))
    1430   {
    1431     g=g+S[2][i];
     1788  while(i+1<=size(Sp[2])&&number(Sp[1][i])<=number(0))
     1789  {
     1790    g=g+Sp[2][i];
    14321791    i++;
    14331792  }
    1434   if(i==size(S[2])&&number(S[1][i])<=number(0))
    1435   {
    1436     g=g+S[2][i];
     1793  if(i==size(Sp[2])&&number(Sp[1][i])<=number(0))
     1794  {
     1795    g=g+Sp[2][i];
    14371796  }
    14381797  return(g);
    … …  
    14411800{ "EXAMPLE:"; echo=2;
    14421801  ring R=0,(x,y),ds;
    1443   list S=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1444   spprint(S);
    1445   spgeomgenus(S);
    1446 }
    1447 ///////////////////////////////////////////////////////////////////////////////
    1448 
    1449 proc spgamma(list S)
    1450 "USAGE:   spgamma(S); list S
    1451 RETURN:  number: gamma invariant of spectrum S
     1802  list Sp=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1803  spprint(Sp);
     1804  spgeomgenus(Sp);
     1805}
     1806///////////////////////////////////////////////////////////////////////////////
     1807
     1808proc spgamma(list Sp)
     1809"USAGE:   spgamma(Sp); list Sp
     1810RETURN:  number: gamma invariant of spectrum Sp
    14521811EXAMPLE: example spgamma; shows examples
    14531812"
    … …  
    14551814  int i,j;
    14561815  number g=0;
    1457   for(i=1;i<=ncols(S[1]);i++)
    1458   {
    1459     for(j=1;j<=S[2][i];j++)
    1460     {
    1461       g=g+(number(S[1][i])-number(nvars(basering)-2)/2)^2;
    1462     }
    1463   }
    1464   g=-g/4+sum(S[2])*number(S[1][ncols(S[1])]-S[1][1])/48;
     1816  for(i=1;i<=ncols(Sp[1]);i++)
     1817  {
     1818    for(j=1;j<=Sp[2][i];j++)
     1819    {
     1820      g=g+(number(Sp[1][i])-number(nvars(basering)-2)/2)^2;
     1821    }
     1822  }
     1823  g=-g/4+sum(Sp[2])*number(Sp[1][ncols(Sp[1])]-Sp[1][1])/48;
    14651824  return(g);
    14661825}
    … …  
    14681827{ "EXAMPLE:"; echo=2;
    14691828  ring R=0,(x,y),ds;
    1470   list S=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
    1471   spprint(S);
    1472   spgamma(S);
    1473 }
    1474 ///////////////////////////////////////////////////////////////////////////////
     1829  list Sp=list(ideal(-1/2,-3/10,-1/10,0,1/10,3/10,1/2),intvec(1,2,2,1,2,2,1));
     1830  spprint(Sp);
     1831  spgamma(Sp);
     1832}
     1833///////////////////////////////////////////////////////////////////////////////
Note: See TracChangeset for help on using the changeset viewer.