/////////////////////////////////////////////////////////////////////////////// version="$Id$"; category="Singularities"; info=" LIBRARY: gmssing.lib Gauss-Manin System of Isolated Singularities AUTHOR: Mathias Schulze, mschulze at mathematik.uni-kl.de OVERVIEW: A library for computing invariants related to the Gauss-Manin system of an isolated hypersurface singularity. REFERENCES: [Sch01] M. Schulze: Algorithms for the Gauss-Manin connection. J. Symb. Comp. 32,5 (2001), 549-564. [Sch02] M. Schulze: The differential structure of the Brieskorn lattice. In: A.M. Cohen et al.: Mathematical Software - ICMS 2002. World Scientific (2002). [Sch03] M. Schulze: Monodromy of Hypersurface Singularities. Acta Appl. Math. 75 (2003), 3-13. [Sch04] M. Schulze: A normal form algorithm for the Brieskorn lattice. J. Symb. Comp. 38,4 (2004), 1207-1225. PROCEDURES: gmsring(t,s); Gauss-Manin system of t with variable s gmsnf(p,K); Gauss-Manin normal form of p gmscoeffs(p,K); Gauss-Manin basis representation of p bernstein(t); Bernstein-Sato polynomial of t monodromy(t); Jordan data of complex monodromy of t spectrum(t); singularity spectrum of t sppairs(t); spectral pairs of t vfilt(t); V-filtration of t on Brieskorn lattice vwfilt(t); weighted V-filtration of t on Brieskorn lattice tmatrix(t); matrix of t w.r.t. good basis of Brieskorn lattice endvfilt(V); endomorphism V-filtration on Jacobian algebra sppnf(a,w[,m]); spectral pairs normal form of (a,w[,m]) sppprint(spp); print spectral pairs spp spadd(sp1,sp2); sum of spectra sp1 and sp2 spsub(sp1,sp2); difference of spectra sp1 and sp2 spmul(sp0,k); linear combination of spectra sp spissemicont(sp[,opt]); semicontinuity test of spectrum sp spsemicont(sp0,sp[,opt]); semicontinuous combinations of spectra sp0 in sp spmilnor(sp); Milnor number of spectrum sp spgeomgenus(sp); geometrical genus of spectrum sp spgamma(sp); gamma invariant of spectrum sp SEE ALSO: mondromy_lib, spectrum_lib, gmspoly_lib, dmod_lib KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; mixed Hodge structure; V-filtration; weight filtration; Bernstein-Sato polynomial; monodromy; spectrum; spectral pairs; good basis "; LIB "linalg.lib"; /////////////////////////////////////////////////////////////////////////////// static proc stdtrans(ideal I) { def @R=basering; string os=ordstr(@R); int j=find(os,",C"); if(j==0) { j=find(os,"C,"); } if(j==0) { j=find(os,",c"); } if(j==0) { j=find(os,"c,"); } if(j>0) { os[j..j+1]=" "; } execute("ring @S="+charstr(@R)+",(gmspoly,"+varstr(@R)+"),(c,dp,"+os+");"); ideal I=homog(imap(@R,I),gmspoly); module M=transpose(transpose(I)+freemodule(ncols(I))); M=std(M); setring(@R); execute("map h=@S,1,"+varstr(@R)+";"); module M=h(M); for(int i=ncols(M);i>=1;i--) { for(j=ncols(M);j>=1;j--) { if(M[i][1]==0) { M[i]=0; } if(i!=j&&M[j][1]!=0) { if(lead(M[i][1])/lead(M[j][1])!=0) { M[i]=0; } } } } M=transpose(simplify(M,2)); I=M[1]; attrib(I,"isSB",1); M=M[2..ncols(M)]; module U=transpose(M); return(list(I,U)); } /////////////////////////////////////////////////////////////////////////////// proc gmsring(poly t,string s) "USAGE: gmsring(t,s); poly t, string s ASSUME: characteristic 0; local degree ordering; isolated critical point 0 of t RETURN: @format ring G; Gauss-Manin system of t with variable s poly gmspoly=t; ideal gmsjacob; Jacobian ideal of t ideal gmsstd; standard basis of Jacobian ideal matrix gmsmatrix; matrix(gmsjacob)*gmsmatrix==matrix(gmsstd) ideal gmsbasis; monomial vector space basis of Jacobian algebra int gmsmaxdeg; maximal weight of variables @end format NOTE: gmsbasis is a C[[s]]-basis of H'' and [t,s]=s^2 KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice EXAMPLE: example gmsring; shows examples " { def @R=basering; if(charstr(@R)!="0") { ERROR("characteristic 0 expected"); } for(int i=nvars(@R);i>=1;i--) { if(var(i)>1) { ERROR("local ordering expected"); } } ideal dt=jacob(t); list l=stdtrans(dt); ideal g=l[1]; if(vdim(g)<=0) { if(vdim(g)==0) { ERROR("singularity at 0 expected"); } else { ERROR("isolated critical point 0 expected"); } } matrix B=l[2]; ideal m=kbase(g); int gmsmaxdeg; for(i=nvars(@R);i>=1;i--) { if(deg(var(i))>gmsmaxdeg) { gmsmaxdeg=deg(var(i)); } } string os=ordstr(@R); int j=find(os,",C"); if(j==0) { j=find(os,"C,"); } if(j==0) { j=find(os,",c"); } if(j==0) { j=find(os,"c,"); } if(j>0) { os[j..j+1]=" "; } execute("ring G="+string(charstr(@R))+",("+s+","+varstr(@R)+"),(ws("+ string(deg(highcorner(g))+2*gmsmaxdeg)+"),"+os+",c);"); poly gmspoly=imap(@R,t); ideal gmsjacob=imap(@R,dt); ideal gmsstd=imap(@R,g); matrix gmsmatrix=imap(@R,B); ideal gmsbasis=imap(@R,m); attrib(gmsstd,"isSB",1); export gmspoly,gmsjacob,gmsstd,gmsmatrix,gmsbasis,gmsmaxdeg; exportto(Top, gmsmaxdeg); return(G); } example { "EXAMPLE:"; echo=2; ring @R=0,(x,y),ds; poly t=x5+x2y2+y5; def G=gmsring(t,"s"); setring(G); gmspoly; print(gmsjacob); print(gmsstd); print(gmsmatrix); print(gmsbasis); gmsmaxdeg; } /////////////////////////////////////////////////////////////////////////////// proc gmsnf(ideal p,int K) "USAGE: gmsnf(p,K); poly p, int K ASSUME: basering returned by gmsring RETURN: list nf; ideal nf[1]; projection of p to C[[s]] mod s^(K+1) ideal nf[2]; p==nf[1]+nf[2] NOTE: computation can be continued by setting p=nf[2] KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice EXAMPLE: example gmsnf; shows examples " { if(system("with","gms")) { return(system("gmsnf",p,gmsstd,gmsmatrix,(K+1)*deg(var(1))-2*gmsmaxdeg,K)); } intvec v=1; v[nvars(basering)]=0; int k; ideal r,q; r[ncols(p)]=0; q[ncols(p)]=0; poly s; int i,j; for(k=ncols(p);k>=1;k--) { while(p[k]!=0&°(lead(p[k]),v)<=K) { i=1; s=lead(p[k])/lead(gmsstd[i]); while(i(K+1)*deg(var(1))-2*gmsmaxdeg&& deg(lead(p[k]),v)<=K) { q[k]=q[k]+lead(p[k]); p[k]=p[k]-lead(p[k]); } } q[k]=q[k]+p[k]; } return(list(r,q)); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; def G=gmsring(t,"s"); setring(G); list l0=gmsnf(gmspoly,0); print(l0[1]); list l1=gmsnf(gmspoly,1); print(l1[1]); list l=gmsnf(l0[2],1); print(l[1]); } /////////////////////////////////////////////////////////////////////////////// proc gmscoeffs(ideal p,int K) "USAGE: gmscoeffs(p,K); poly p, int K ASSUME: basering constructed by gmsring RETURN: @format list l; matrix l[1]; C[[s]]-basis representation of p mod s^(K+1) ideal l[2]; p==matrix(gmsbasis)*l[1]+l[2] @end format NOTE: computation can be continued by setting p=l[2] KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice EXAMPLE: example gmscoeffs; shows examples " { list l=gmsnf(p,K); ideal r,q=l[1..2]; poly v=1; for(int i=2;i<=nvars(basering);i++) { v=v*var(i); } matrix C=coeffs(r,gmsbasis,v); return(list(C,q)); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; def G=gmsring(t,"s"); setring(G); list l0=gmscoeffs(gmspoly,0); print(l0[1]); list l1=gmscoeffs(gmspoly,1); print(l1[1]); list l=gmscoeffs(l0[2],1); print(l[1]); } /////////////////////////////////////////////////////////////////////////////// static proc mindegree(matrix A) { int d=-1; int i,j; for(i=nrows(A);i>=1;i--) { for(j=ncols(A);j>=1;j--) { if(d==-1||(ord(A[i,j])-1)) { d=ord(A[i,j]); } } } return(d); } /////////////////////////////////////////////////////////////////////////////// static proc maxdegree(matrix A) { int d=-1; int i,j; for(i=nrows(A);i>=1;i--) { for(j=ncols(A);j>=1;j--) { if(deg(A[i,j])>d) { d=deg(A[i,j]); } } } return(d); } /////////////////////////////////////////////////////////////////////////////// static proc saturate() { int mu=ncols(gmsbasis); ideal r=gmspoly*gmsbasis; matrix A0[mu][mu],C; module H0; module H,H1=freemodule(mu),freemodule(mu); int k=-1; list l; dbprint(printlevel-voice+2,"// compute saturation of H''"); while(size(reduce(H,std(H0*var(1))))>0) { dbprint(printlevel-voice+2,"// compute matrix A of t"); k++; dbprint(printlevel-voice+2,"// k="+string(k)); l=gmscoeffs(r,k); C,r=l[1..2]; A0=A0+C; dbprint(printlevel-voice+2,"// compute saturation step"); H0=H; H1=jet(module(A0*H1+var(1)^2*diff(matrix(H1),var(1))),k); H=H*var(1)+H1; } A0=A0-k*var(1); dbprint(printlevel-voice+2,"// compute basis of saturation of H''"); H=std(H0); dbprint(printlevel-voice+2,"// transform H'' to saturation of H''"); H0=division(freemodule(mu)*var(1)^k,H,k*deg(var(1)))[1]; return(A0,r,H,H0,k); } /////////////////////////////////////////////////////////////////////////////// static proc tjet(matrix A0,ideal r,module H,int K0,int K) { dbprint(printlevel-voice+2,"// compute matrix A of t"); dbprint(printlevel-voice+2,"// k="+string(K0+K+1)); list l=gmscoeffs(r,K0+K+1); matrix C; C,r=l[1..2]; A0=A0+C; dbprint(printlevel-voice+2,"// transform A to saturation of H''"); matrix A=division(A0*H+var(1)^2*diff(matrix(H),var(1)),H, (K+1)*deg(var(1)))[1]/var(1); return(A,A0,r); } /////////////////////////////////////////////////////////////////////////////// static proc eigenval(matrix A0,ideal r,module H,int K0) { dbprint(printlevel-voice+2, "// compute eigenvalues e with multiplicities m of A1"); matrix A; A,A0,r=tjet(A0,r,H,K0,0); list l=eigenvals(A); def e,m=l[1..2]; dbprint(printlevel-voice+2,"// e="+string(e)); dbprint(printlevel-voice+2,"// m="+string(m)); return(e,m,A0,r); } /////////////////////////////////////////////////////////////////////////////// static proc transform(matrix A,matrix A0,ideal r,module H,module H0,ideal e, intvec m,int K0,int K,int opt) { int mu=ncols(gmsbasis); int i,j,k; intvec d; d[ncols(e)]=0; if(opt) { dbprint(printlevel-voice+2, "// compute rounded maximal differences d of e"); for(i=1;i<=ncols(e);i++) { d[i]=int(e[ncols(e)]-e[i]); } } else { dbprint(printlevel-voice+2, "// compute maximal integer differences d of e"); for(i=1;id[i]) { d[i]=k; } if(-k>d[j]) { d[j]=-k; } } } } } dbprint(printlevel-voice+2,"// d="+string(d)); for(i,k=1,0;i<=size(d);i++) { if(k0) { int i0,j0,i1,j1; list l; while(k>0) { dbprint(printlevel-voice+2,"// transform to separate eigenvalues"); U=0; for(i=1;i<=ncols(e);i++) { U=U+syz(power(jet(A,0)-e[i],m[i])); } V=inverse(U); A=V*A*U; H=H*U; H0=V*H0; dbprint(printlevel-voice+2,"// transform to reduce maximum of d by 1"); for(i0,i=1,1;i0<=ncols(e);i0++) { for(i1=1;i1<=m[i0];i1,i=i1+1,i+1) { for(j0,j=1,1;j0<=ncols(e);j0++) { for(j1=1;j1<=m[j0];j1,j=j1+1,j+1) { if(d[i0]==0&&d[j0]>=1) { A[i,j]=A[i,j]*var(1); } if(d[i0]>=1&&d[j0]==0) { A[i,j]=A[i,j]/var(1); } } } } } H0=transpose(H0); for(i0,i=1,1;i0<=ncols(e);i0++) { if(d[i0]>=1) { for(i1=1;i1<=m[i0];i1,i=i1+1,i+1) { H[i]=H[i]*var(1); } d[i0]=d[i0]-1; } else { for(i1=1;i1<=m[i0];i1,i=i1+1,i+1) { A[i,i]=A[i,i]-1; H0[i]=H0[i]*var(1); } e[i0]=e[i0]-1; } } H0=transpose(H0); l=sppnf(list(e,d,m)); e,d,m=l[1..3]; k--; K0++; } A=jet(A,K); } dbprint(printlevel-voice+2,"// transform to separate eigenvalues"); U=0; for(i=1;i<=ncols(e);i++) { U=U+syz(power(jet(A,0)-e[i],m[i])); } V=inverse(U); A=V*A*U; H=H*U; H0=V*H0; return(A,A0,r,H,H0,e,m,K0); } /////////////////////////////////////////////////////////////////////////////// proc bernstein(poly t) "USAGE: bernstein(t); poly t ASSUME: characteristic 0; local degree ordering; isolated critical point 0 of t RETURN: @format list bs; Bernstein-Sato polynomial b(s) of t ideal bs[1]; number bs[1][i]; i-th root of b(s) intvec bs[2]; int bs[2][i]; multiplicity of i-th root of b(s) @end format KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; Bernstein-Sato polynomial EXAMPLE: example bernstein; shows examples " { def @R=basering; int n=nvars(@R)-1; def @G=gmsring(t,"s"); setring(@G); matrix A; module U0; ideal e; intvec m; def A0,r,H,H0,K0=saturate(); A,A0,r=tjet(A0,r,H,K0,0); list l=minipoly(A); e,m=l[1..2]; e=-e; l=spnf(spadd(list(e,m),list(ideal(-1),intvec(1)))); setring(@R); list l=imap(@G,l); kill @G,gmsmaxdeg; return(l); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; bernstein(t); } /////////////////////////////////////////////////////////////////////////////// proc monodromy(poly t) "USAGE: monodromy(t); poly t ASSUME: characteristic 0; local degree ordering; isolated critical point 0 of t RETURN: @format list l; Jordan data jordan(M) of monodromy matrix exp(-2*pi*i*M) ideal l[1]; number l[1][i]; eigenvalue of i-th Jordan block of M intvec l[2]; int l[2][i]; size of i-th Jordan block of M intvec l[3]; int l[3][i]; multiplicity of i-th Jordan block of M @end format SEE ALSO: mondromy_lib, linalg_lib KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; monodromy EXAMPLE: example monodromy; shows examples " { def @R=basering; int n=nvars(@R)-1; def @G=gmsring(t,"s"); setring(@G); matrix A; module U0; ideal e; intvec m; def A0,r,H,H0,K0=saturate(); e,m,A0,r=eigenval(A0,r,H,K0); A,A0,r,H,H0,e,m,K0=transform(A,A0,r,H,H0,e,m,K0,0,0); list l=jordan(A,e,m); setring(@R); list l=imap(@G,l); kill @G,gmsmaxdeg; return(l); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; monodromy(t); } /////////////////////////////////////////////////////////////////////////////// proc spectrum(poly t) "USAGE: spectrum(t); poly t ASSUME: characteristic 0; local degree ordering; isolated critical point 0 of t RETURN: @format list sp; singularity spectrum of t ideal sp[1]; number sp[1][i]; i-th spectral number intvec sp[2]; int sp[2][i]; multiplicity of i-th spectral number @end format SEE ALSO: spectrum_lib KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; mixed Hodge structure; V-filtration; spectrum EXAMPLE: example spectrum; shows examples " { list l=vwfilt(t); return(spnf(list(l[1],l[3]))); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; spprint(spectrum(t)); } /////////////////////////////////////////////////////////////////////////////// proc sppairs(poly t) "USAGE: sppairs(t); poly t ASSUME: characteristic 0; local degree ordering; isolated critical point 0 of t RETURN: @format list spp; spectral pairs of t ideal spp[1]; number spp[1][i]; V-filtration index of i-th spectral pair intvec spp[2]; int spp[2][i]; weight filtration index of i-th spectral pair intvec spp[3]; int spp[3][i]; multiplicity of i-th spectral pair @end format SEE ALSO: spectrum_lib KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; mixed Hodge structure; V-filtration; weight filtration; spectrum; spectral pairs EXAMPLE: example sppairs; shows examples " { list l=vwfilt(t); return(list(l[1],l[2],l[3])); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; sppprint(sppairs(t)); } /////////////////////////////////////////////////////////////////////////////// proc vfilt(poly t) "USAGE: vfilt(t); poly t ASSUME: characteristic 0; local degree ordering; isolated critical point 0 of t RETURN: @format list v; V-filtration on H''/s*H'' ideal v[1]; number v[1][i]; V-filtration index of i-th spectral number intvec v[2]; int v[2][i]; multiplicity of i-th spectral number list v[3]; module v[3][i]; vector space of i-th graded part in terms of v[4] ideal v[4]; monomial vector space basis of H''/s*H'' ideal v[5]; standard basis of Jacobian ideal @end format SEE ALSO: spectrum_lib KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; mixed Hodge structure; V-filtration; spectrum EXAMPLE: example vfilt; shows examples " { list l=vwfilt(t); return(spnf(list(l[1],l[3],l[4]))+list(l[5],l[6])); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; vfilt(t); } /////////////////////////////////////////////////////////////////////////////// proc vwfilt(poly t) "USAGE: vwfilt(t); poly t ASSUME: characteristic 0; local degree ordering; isolated critical point 0 of t RETURN: @format list vw; weighted V-filtration on H''/s*H'' ideal vw[1]; number vw[1][i]; V-filtration index of i-th spectral pair intvec vw[2]; int vw[2][i]; weight filtration index of i-th spectral pair intvec vw[3]; int vw[3][i]; multiplicity of i-th spectral pair list vw[4]; module vw[4][i]; vector space of i-th graded part in terms of vw[5] ideal vw[5]; monomial vector space basis of H''/s*H'' ideal vw[6]; standard basis of Jacobian ideal @end format SEE ALSO: spectrum_lib KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; mixed Hodge structure; V-filtration; weight filtration; spectrum; spectral pairs EXAMPLE: example vwfilt; shows examples " { def @R=basering; int n=nvars(@R)-1; def @G=gmsring(t,"s"); setring(@G); int mu=ncols(gmsbasis); matrix A; ideal e; intvec m; def A0,r,H,H0,K0=saturate(); e,m,A0,r=eigenval(A0,r,H,K0); A,A0,r,H,H0,e,m,K0=transform(A,A0,r,H,H0,e,m,K0,0,1); dbprint(printlevel-voice+2,"// compute weight filtration basis"); list l=jordanbasis(A,e,m); def U,v=l[1..2]; kill l; vector u0; int v0; int i,j,k,l; for(k,l=1,1;l<=ncols(e);k,l=k+m[l],l+1) { for(i=k+m[l]-1;i>=k+1;i--) { for(j=i-1;j>=k;j--) { if(v[i]>v[j]) { v0=v[i];v[i]=v[j];v[j]=v0; u0=U[i];U[i]=U[j];U[j]=u0; } } } } dbprint(printlevel-voice+2,"// transform to weight filtration basis"); matrix V=inverse(U); A=V*A*U; dbprint(printlevel-voice+2,"// compute standard basis of H''"); H=H*U; H0=std(V*H0); dbprint(printlevel-voice+2,"// compute spectral pairs"); ideal a; intvec w; for(i=1;i<=mu;i++) { j=leadexp(H0[i])[nvars(basering)+1]; a[i]=A[j,j]+ord(H0[i]) div deg(var(1))-1; w[i]=v[j]+n; } H=H*H0; H=simplify(jet(H/var(1)^(mindegree(H) div deg(var(1))),0),1); kill l; list l=sppnf(list(a,w,H))+list(gmsbasis,gmsstd); setring(@R); list l=imap(@G,l); kill @G,gmsmaxdeg; attrib(l[5],"isSB",1); return(l); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; vwfilt(t); } /////////////////////////////////////////////////////////////////////////////// static proc fsplit(ideal e0,intvec m0,matrix A,module H,module H0) { int mu=ncols(gmsbasis); dbprint(printlevel-voice+2,"// compute standard basis of H''"); H0=std(H0); H0=simplify(H0,1); dbprint(printlevel-voice+2,"// compute Hodge filtration"); int i,j,k; ideal e; intvec m; e[mu]=0; for(i=1;i<=ncols(e0);i++) { for(j=m0[i];j>=1;j--) { k++; e[k]=e0[i]; m[k]=i; } } number n,n0; vector v,v0; list F; for(i=ncols(e0);i>=1;i--) { F[i]=module(matrix(0,mu,1)); } for(i=mu;i>=1;i--) { v=H0[i]; v0=lead(v); n0=leadcoef(e[leadexp(v0)[nvars(basering)+1]])+leadexp(v0)[1]; v=v-lead(v); while(v!=0) { n=leadcoef(e[leadexp(v)[nvars(basering)+1]])+leadexp(v)[1]; if(n==n0) { v0=v0+lead(v); v=v-lead(v); } else { v=0; } } j=m[leadexp(v0)[nvars(basering)+1]]; F[j]=F[j]+v0; } dbprint(printlevel-voice+2,"// compute splitting of Hodge filtration"); matrix A0=jet(A,0); module U,U0,U1,U2; matrix N; s=0; for(i=size(F);i>=1;i--) { N=A0-e0[i]; U0=0; while(size(F[i])>0) { U1=jet(F[i],0); k=0; while(size(U1)>0) { for(j=ncols(U1);j>=1;j--) { if(size(reduce(U1[j],std(U0)))>0) { U0=U0+U1[j]; } } U1=N*U1; k++; } F[i]=module(F[i]/var(1)); } U=U0+U; } dbprint(printlevel-voice+2,"// transform to Hodge splitting basis"); H=H*U; H0=lift(U,H0); A=lift(U,A*U); return(e,A,H,H0); } /////////////////////////////////////////////////////////////////////////////// static proc glift(ideal e,matrix A,module H,module H0,int K) { poly s=var(1); int mu=ncols(gmsbasis); dbprint(printlevel-voice+2,"// compute standard basis of H''"); H0=std(H0); H0=simplify(H0,1); int i,j,k; ideal v; for(i=mu;i>=1;i--) { v[i]=e[leadexp(H0[i])[nvars(basering)+1]]+leadexp(H0[i])[1]; } dbprint(printlevel-voice+2, "// compute matrix A0 of t w.r.t. good basis H0 of H''"); number c; matrix h0[mu][1]; matrix m[mu][1]; matrix a0[mu][1]; matrix A0[mu][mu]; module M=H0; module N=jet(s*A*matrix(H0)+s^2*diff(matrix(H0),s),K+1); while(size(N)>0) { j=mu; for(k=mu-1;k>=1;k--) { if(N[k]>N[j]) { j=k; } } i=mu; while(leadexp(M[i])[nvars(basering)+1]!=leadexp(N[j])[nvars(basering)+1]) { i--; } k=leadexp(N[j])[1]-leadexp(M[i])[1]; if(k==0||i==j) { dbprint(printlevel-voice+3,"// compute A0["+string(i)+","+string(j)+"]"); c=leadcoef(N[j])/leadcoef(M[i]); A0[i,j]=A0[i,j]+c*s^k; N[j]=jet(N[j]-c*s^k*M[i],K+1); } else { dbprint(printlevel-voice+3, "// reduce H0["+string(j)+"] with H0["+string(i)+"]"); c=leadcoef(N[j])/leadcoef(M[i])/(1-k-leadcoef(v[i])+leadcoef(v[j])); H0[j]=H0[j]+c*s^(k-1)*H0[i]; M[j]=M[j]+c*s^(k-1)*M[i]; h0=c*s^(k-1)*H0[i]; N[j]=N[j]+jet(s*A*h0+s^2*diff(h0,s),K+1)[1]; m=M[i]; a0=transpose(A0)[j]; N=N-jet(c*s^(k-1)*m*transpose(a0),K+1); } } H0=H*H0; H0=H0/var(1)^(mindegree(H0) div deg(var(1))); return(A0,H0); } /////////////////////////////////////////////////////////////////////////////// proc tmatrix(poly t) "USAGE: tmatrix(t); poly t ASSUME: characteristic 0; local degree ordering; isolated critical point 0 of t RETURN: @format list l=A0,A1,T,M; matrix A0,A1; t=A0+s*A1+s^2*(d/ds) on H'' w.r.t. C[[s]]-basis M*T module T; C-basis of C^mu ideal M; monomial C-basis of H''/sH'' @end format KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; mixed Hodge structure; V-filtration; weight filtration; monodromy; spectrum; spectral pairs; good basis EXAMPLE: example tmatrix; shows examples " { def @R=basering; int n=nvars(@R)-1; def @G=gmsring(t,"s"); setring(@G); int mu=ncols(gmsbasis); matrix A; module U0; ideal e; intvec m; def A0,r,H,H0,K0=saturate(); e,m,A0,r=eigenval(A0,r,H,K0); A,A0,r,H,H0,e,m,K0=transform(A,A0,r,H,H0,e,m,K0,K0+int(e[ncols(e)]-e[1]),1); A,H0=glift(fsplit(e,m,A,H,H0),K0); A0=jet(A,0); A=jet(A/var(1),0); list l=A0,A,H0,gmsbasis; setring(@R); list l=imap(@G,l); kill @G,gmsmaxdeg; return(l); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; list l=tmatrix(t); print(l[1]); print(l[2]); print(l[3]); print(l[4]); } /////////////////////////////////////////////////////////////////////////////// proc endvfilt(list v) "USAGE: endvfilt(v); list v ASSUME: v returned by vfilt RETURN: @format list ev; V-filtration on Jacobian algebra ideal ev[1]; number ev[1][i]; i-th V-filtration index intvec ev[2]; int ev[2][i]; i-th multiplicity list ev[3]; module ev[3][i]; vector space of i-th graded part in terms of ev[4] ideal ev[4]; monomial vector space basis of Jacobian algebra ideal ev[5]; standard basis of Jacobian ideal @end format KEYWORDS: singularities; Gauss-Manin system; Brieskorn lattice; mixed Hodge structure; V-filtration; endomorphism filtration EXAMPLE: example endvfilt; shows examples " { def a,d,V,m,g=v[1..5]; attrib(g,"isSB",1); int mu=ncols(m); module V0=V[1]; for(int i=2;i<=size(V);i++) { V0=V0,V[i]; } dbprint(printlevel-voice+2,"// compute multiplication in Jacobian algebra"); list M; module U=freemodule(ncols(m)); for(i=ncols(m);i>=1;i--) { M[i]=division(coeffs(reduce(m[i]*m,g,U),m)*V0,V0)[1]; } int j,k,i0,j0,i1,j1; number b0=number(a[1]-a[ncols(a)]); number b1,b2; matrix M0; module L; list v0=freemodule(ncols(m)); ideal a0=b0; list l; while(b0=1;j--) { for(i=ncols(a);i>=1;i--) { b2=number(a[i]-a[j]); if(b2>b0&&b2=2;k--) { l=l+list(ideal()); } dbprint(printlevel-voice+2,"// collect conditions for EV["+string(b0)+"]"); j=ncols(a); j0=mu; while(j>=1) { i0=1; i=1; while(i1) { l[k]=l[k],M0[1..i0-1,j0-d[j]+1..j0]; } } j0=j0-d[j]; j--; } dbprint(printlevel-voice+2,"// compose condition matrix"); L=transpose(module(l[1])); for(k=2;k<=ncols(m);k++) { L=L,transpose(module(l[k])); } dbprint(printlevel-voice+2,"// compute kernel of condition matrix"); v0=v0+list(syz(L)); a0=a0,b0; } dbprint(printlevel-voice+2,"// compute graded parts"); option(redSB); for(i=1;i0) { v1=v1+list(v0[i]); d1=d1,size(v0[i]); a1=a1,a0[i]; } i++; } return(list(a1,d1,v1,m,g)); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; poly t=x5+x2y2+y5; endvfilt(vfilt(t)); } /////////////////////////////////////////////////////////////////////////////// proc sppnf(list sp) "USAGE: sppnf(list(a,w[,m])); ideal a, intvec w, intvec m ASSUME: ncols(a)==size(w)==size(m) RETURN: order (a[i][,w[i]]) with multiplicity m[i] lexicographically EXAMPLE: example sppnf; shows examples " { ideal a=sp[1]; intvec w=sp[2]; int n=ncols(a); intvec m; list V; module v; int i,j; for(i=3;i<=size(sp);i++) { if(typeof(sp[i])=="intvec") { m=sp[i]; } if(typeof(sp[i])=="module") { v=sp[i]; for(j=n;j>=1;j--) { V[j]=module(v[j]); } } if(typeof(sp[i])=="list") { V=sp[i]; } } if(m==0) { for(i=n;i>=1;i--) { m[i]=1; } } int k; ideal a0; intvec w0,m0; list V0; number a1; int w1,m1; for(i=n;i>=1;i--) { if(m[i]!=0) { for(j=i-1;j>=1;j--) { if(m[j]!=0) { if(number(a[i])>number(a[j])|| (number(a[i])==number(a[j])&&w[i]0) { v=V[i]; V[i]=V[j]; V[j]=v; } } if(number(a[i])==number(a[j])&&w[i]==w[j]) { m[i]=m[i]+m[j]; m[j]=0; if(size(V)>0) { V[i]=V[i]+V[j]; } } } } k++; a0[k]=a[i]; w0[k]=w[i]; m0[k]=m[i]; if(size(V)>0) { V0[k]=V[i]; } } } if(size(V0)>0) { n=size(V0); module U=std(V0[n]); for(i=n-1;i>=1;i--) { V0[i]=simplify(reduce(V0[i],U),1); if(i>=2) { U=std(U+V0[i]); } } } if(k>0) { sp=a0,w0,m0; if(size(V0)>0) { sp[4]=V0; } } return(sp); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; list sp=list(ideal(-1/2,-3/10,-3/10,-1/10,-1/10,0,1/10,1/10,3/10,3/10,1/2), intvec(2,1,1,1,1,1,1,1,1,1,0)); sppprint(sppnf(sp)); } /////////////////////////////////////////////////////////////////////////////// proc sppprint(list spp) "USAGE: sppprint(spp); list spp RETURN: string s; spectral pairs spp EXAMPLE: example sppprint; shows examples " { string s; for(int i=1;inumber(sp2[1][i2])) { s[i]=sp2[1][i2]; m[i]=sp2[2][i2]; i++; i2++; } else { if(sp1[2][i1]+sp2[2][i2]!=0) { s[i]=sp1[1][i1]; m[i]=sp1[2][i1]+sp2[2][i2]; i++; } i1++; i2++; } } } else { s[i]=sp1[1][i1]; m[i]=sp1[2][i1]; i++; i1++; } } else { s[i]=sp2[1][i2]; m[i]=sp2[2][i2]; i++; i2++; } } return(list(s,m)); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; 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)); spprint(sp1); list sp2=list(ideal(-1/6,1/6),intvec(1,1)); spprint(sp2); spprint(spadd(sp1,sp2)); } /////////////////////////////////////////////////////////////////////////////// proc spsub(list sp1,list sp2) "USAGE: spsub(sp1,sp2); list sp1, list sp2 RETURN: list sp; difference of spectra sp1 and sp2 EXAMPLE: example spsub; shows examples " { return(spadd(sp1,spmul(sp2,-1))); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; 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)); spprint(sp1); list sp2=list(ideal(-1/6,1/6),intvec(1,1)); spprint(sp2); spprint(spsub(sp1,sp2)); } /////////////////////////////////////////////////////////////////////////////// proc spmul(list #) "USAGE: spmul(sp0,k); list sp0, int[vec] k RETURN: list sp; linear combination of spectra sp0 with coefficients k EXAMPLE: example spmul; shows examples " { if(size(#)==2) { if(typeof(#[1])=="list") { if(typeof(#[2])=="int") { return(list(#[1][1],#[1][2]*#[2])); } if(typeof(#[2])=="intvec") { list sp0=list(ideal(),intvec(0)); for(int i=size(#[2]);i>=1;i--) { sp0=spadd(sp0,spmul(#[1][i],#[2][i])); } return(sp0); } } } return(list(ideal(),intvec(0))); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; 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)); spprint(sp); spprint(spmul(sp,2)); list sp1=list(ideal(-1/6,1/6),intvec(1,1)); spprint(sp1); list sp2=list(ideal(-1/3,0,1/3),intvec(1,2,1)); spprint(sp2); spprint(spmul(list(sp1,sp2),intvec(1,2))); } /////////////////////////////////////////////////////////////////////////////// proc spissemicont(list sp,list #) "USAGE: spissemicont(sp[,1]); list sp, int opt RETURN: @format int k= 1; if sum of sp is positive on all intervals [a,a+1) [and (a,a+1)] 0; if sum of sp is negative on some interval [a,a+1) [or (a,a+1)] @end format EXAMPLE: example spissemicont; shows examples " { int opt=0; if(size(#)>0) { if(typeof(#[1])=="int") { opt=1; } } int i,j,k; i=1; while(i<=size(sp[2])-1) { j=i+1; k=0; while(j+1<=size(sp[2])&&number(sp[1][j])<=number(sp[1][i])+1) { if(opt==0||number(sp[1][j])1) { l0=spsemicont(sp0,list(sp[1..size(sp)-1])); for(i=1;i<=size(l0);i++) { if(size(l)>0) { j=1; while(j1) { return(l); } else { return(list(intvec(k-1))); } } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; 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)); spprint(sp0); list sp1=list(ideal(-1/6,1/6),intvec(1,1)); spprint(sp1); list sp2=list(ideal(-1/3,0,1/3),intvec(1,2,1)); spprint(sp2); list sp=sp1,sp2; list l=spsemicont(sp0,sp); l; spissemicont(spsub(sp0,spmul(sp,l[1]))); spissemicont(spsub(sp0,spmul(sp,l[1]-1))); spissemicont(spsub(sp0,spmul(sp,l[1]+1))); } /////////////////////////////////////////////////////////////////////////////// proc spmilnor(list sp) "USAGE: spmilnor(sp); list sp RETURN: int mu; Milnor number of spectrum sp EXAMPLE: example spmilnor; shows examples " { return(sum(sp[2])); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; 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)); spprint(sp); spmilnor(sp); } /////////////////////////////////////////////////////////////////////////////// proc spgeomgenus(list sp) "USAGE: spgeomgenus(sp); list sp RETURN: int g; geometrical genus of spectrum sp EXAMPLE: example spgeomgenus; shows examples " { int g=0; int i=1; while(i+1<=size(sp[2])&&number(sp[1][i])<=number(0)) { g=g+sp[2][i]; i++; } if(i==size(sp[2])&&number(sp[1][i])<=number(0)) { g=g+sp[2][i]; } return(g); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; 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)); spprint(sp); spgeomgenus(sp); } /////////////////////////////////////////////////////////////////////////////// proc spgamma(list sp) "USAGE: spgamma(sp); list sp RETURN: number gamma; gamma invariant of spectrum sp EXAMPLE: example spgamma; shows examples " { int i,j; number g=0; for(i=1;i<=ncols(sp[1]);i++) { for(j=1;j<=sp[2][i];j++) { g=g+(number(sp[1][i])-number(nvars(basering)-2)/2)^2; } } g=-g/4+sum(sp[2])*number(sp[1][ncols(sp[1])]-sp[1][1])/48; return(g); } example { "EXAMPLE:"; echo=2; ring R=0,(x,y),ds; 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)); spprint(sp); spgamma(sp); } ///////////////////////////////////////////////////////////////////////////////