Ignore:
Timestamp:
Sep 24, 2010, 4:52:16 PM (14 years ago)
Author:
Stefan Steidel <steidel@…>
Branches:
(u'spielwiese', '17f1d200f27c5bd38f5dfc6e8a0879242279d1d8')
Children:
f0e080284bc1a4ba286ba66df791c785e808feb8
Parents:
1492bbce4e131a22de468dd41b9ced71d0fe3f93
Message:
new libs

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

Legend:

Unmodified
Added
Removed

  • Singular/LIB/surfacesignature.lib

    r1492bb r2e7410  
    22category="Singularities";
    33info="
    4 LIBRARY:  surfacesignature.lib        signature of surface singularity
    5 AUTHORS:  Gerhard Pfister,            pfister@mathematik.uni-kl.de
    6           Muhammad Ahsan Banyamin,    ahsanbanyamin@gmail.com
    7 
    8 OVERVIEW:
    9  A library for computing the signature of irreducible surface singularity.
    10  It contains also a procedure for computing the newton pairs of a surface
    11  singularity.
    12 
    13 THEORY: The signature of a surface singularity is defined in
    14         [Durfee,A.:The Signature of Smoothings of Complex Surface
    15         Singularities, Math. Ann.,232,85-98 (1978)].
    16         The algorithm we use has been proposed in
    17         [Nemethi,A.:The signature of f(x,y)+z^n, Proceedings of Singularity
    18         Conference (C.T.C Wall's 60th birthday meeting),Liverpool 1996,
    19         London Math.Soc. LN 263(1999),131-149].
     4LIBRARY:  signaturesurfsing.lib       signature of surface singularity
     5
     6AUTHORS:  Gerhard Pfister             pfister@mathematik.uni-kl.de
     7@*        Muhammad Ahsan Banyamin     ahsanbanyamin@gmail.com
     8@*        Stefan Steidel              steidel@mathematik.uni-kl.de
     9
     10OVERVIEW:
     11
     12  A library for computing the signature of irreducible surface singularity.
     13  The signature of a surface singularity is defined in [Durfee, A.: The
     14  Signature of Smoothings of Complex Surface Singularities, Math. Ann., 232,
     15  85-98 (1978)]. The algorithm we use has been proposed in [Nemethi, A.: The
     16  signature of f(x,y)+z^N, Proceedings of Singularity Conference (C.T.C. Wall's
     17  60th birthday meeting), Liverpool 1996, London Math.Soc. LN 263(1999),
     18  131-149].
    2019
    2120PROCEDURES:
    22  brieskornSign(a,b,c);   signature of Brieskorn singularity x^a+y^b+z^c=0
    23  newtonpairs(f);         newton pairs of surface singularity
    24  signature(N,f);         signature of singularity z^+f(x,y)=0,f irreducible
     21 brieskornSign(a1,a2,a3);  signature of Brieskorn singularity x^a1+y^a2+z^a3
     22 signature(N,f);           signature of singularity z^N+f(x,y)=0, f irreducible
    2523";
    2624
    2725LIB "hnoether.lib";
    28 ///////////////////////////////////////////////////////////////////////////////
    29 proc signature(int N,poly f)
     26LIB "alexpoly.lib";
     27LIB "gmssing.lib";
     28
     29///////////////////////////////////////////////////////////////////////////////
     30//------- sigma(z^N + f) in terms of Puiseux pairs of f for f irreducible ----
     31
     32static proc exponentSequence(poly f)
     33//=== computes the sequence a_1,...,a_s of exponents as described in [Nemethi]
     34//=== using the Puiseux pairs (m_1, n_1),...,(m_s, n_s) of f:
     35//===  - a_1 = m_1,
     36//===  - a_i = m_i - n_i * (m_[i-1] - n_[i-1] * a_[i-1]).
     37//===
     38//=== Return: list of two intvecs:
     39//===         1st entry: A = (a_1,...,a_s)
     40//===         2nd entry: N = (n_1,...,n_s)
     41{
     42   def R = basering;
     43   ring S = 0,(x,y),dp;
     44   poly f = fetch(R,f);
     45   list puiseuxPairs = invariants(f);
     46   setring R;
     47   
     48   intvec M = puiseuxPairs[1][3];
     49   intvec N = puiseuxPairs[1][4];
     50   
     51   int i;
     52   int a = M[1];
     53   intvec A = a;
     54   for(i = 2; i <= size(M); i++)
     55   {
     56      a = M[i] - N[i] * (M[i-1] - N[i-1] * a);
     57      A[size(A)+1] = a;
     58   }
     59 
     60   return(list(A,N));
     61}
     62example
     63{ "EXAMPLE:"; echo = 2;
     64   ring r = 0,(x,y),dp;
     65   exponentSequence(y4+2x3y2+x6+x5y);
     66}
     67
     68///////////////////////////////////////////////////////////////////////////////
     69
     70proc brieskornSign(a1,a2,a3)
     71"USAGE:  brieskornSign(a1,a2,a3); a1,a2,a3 = integers
     72RETURN:  signature of Brieskorn singularity x^a1+y^a2+z^a3
     73EXAMPLE: example brieskornSign; shows an example
    3074"
    31 USAGE:   signature(N,f); N=integer, f=irreducible poly in 2 variables
    32 RETURN:  signature of surface singularity defined by z^N+f=0
    33 EXAMPLE: example signature; shows an example
     75{
     76   int a_temp, t, k1, k2, k3, s_t, sigma;
     77   number s;
     78   
     79   if(a1 > a2) { a_temp = a1; a1 = a2; a2 = a_temp; }
     80   if(a2 > a3) { a_temp = a2; a2 = a3; a3 = a_temp; }
     81   if(a1 > a2) { a_temp = a1; a1 = a2; a2 = a_temp; }
     82   
     83   for(t = 0; t <= 2; t++)
     84   {
     85      s_t = 0;
     86      for(k1 = 1; k1 <= a1-1; k1++)
     87      {
     88         for(k2 = 1; k2 <= a2-1; k2++)
     89         {
     90            for(k3 = 1; k3 <= a3-1; k3++)
     91            {
     92               s = number(k1)/a1 + number(k2)/a2 + number(k3)/a3;
     93               if(t < s)
     94               {
     95                  if(s < t+1)
     96                  {
     97                     s_t = s_t + 1;
     98                  }
     99                  else
     100                  {
     101                     break;
     102                  }
     103               }
     104            }
     105            if(k3 == 1) { break; }
     106         }
     107         if(k2 == 1) { break; }
     108      }
     109      sigma = sigma + (-1)^t * s_t;
     110   }
     111   return(sigma);
     112}
     113example
     114{ "EXAMPLE:"; echo = 2;
     115   ring R = 0,x,dp;
     116   brieskornSign(11,3,5);
     117}
     118
     119///////////////////////////////////////////////////////////////////////////////
     120
     121static proc signatureP(int N,poly f)
     122"USAGE:  signatureP(N,f); N = integer, f = irreducible poly in 2 variables
     123RETURN:  signature of surface singularity defined by z^N + f(x,y) = 0
     124EXAMPLE: example signatureP; shows an example
    34125"
    35126{
    36    def R=basering;
    37    ring S=0,(x,y),dp;
    38    poly f = fetch(R,f);
    39    list L=newtonpairs(f);
    40    setring R;
    41    if(size(L)>2)
    42    {
    43      return(Generalcase(N,L));
    44    }
    45    if(size(L)==2)
    46    {
    47      return(signtwopairs(N,L));
    48    }
    49    return(brieskornSign(L[1][2],L[1][1],N));
    50  }
    51  example
     127   int i, d, prod, sigma;
     128   list L = exponentSequence(f);
     129   int s = size(L[2]);
     130   
     131   if(s == 1)
     132   {
     133      return(brieskornSign(L[1][1], L[2][1], N));
     134   }
     135   
     136   prod = 1;
     137   sigma = brieskornSign(L[1][s], L[2][s], N);
     138   for(i = s - 1; i >= 1; i--)
     139   {
     140      prod = prod * L[2][i+1];
     141      d = gcd(N, prod);
     142      sigma = sigma + d * brieskornSign(L[1][i], L[2][i], N/d);
     143   }
     144   
     145   return(sigma);
     146}
     147example
    52148{ "EXAMPLE:"; echo = 2;
    53149   ring r = 0,(x,y),dp;
    54150   int N  = 3;
    55    poly f= x15-21x14+8x13y-6x13-16x12y+20x11y2-x12+8x11y-36x10y2
    56       +24x9y3+4x9y2-16x8y3+26x7y4-6x6y4+8x5y5+4x3y6-y8;
    57    signature(N,f);
    58 }
    59 
    60 ///////////////////////////////////////////////////////////////////////////
    61 proc newtonpairs(poly f)
    62 "USAGE:   newtonpairs(f); f= irreducible bivariate poly
    63  RETURN:  newton pairs of curve singularity f(x,y)=0
    64  EXAMPLE: example newtonpairs; shows an example
    65 "
    66 {
    67  def R=basering;
    68  ring S=0,(x,y),dp;
    69  poly f = fetch(R,f);
    70  list M=invariants(f);
    71  setring R;
    72  list L=M[1][3];
    73  list K=M[1][4];
    74  int i,j;
    75  intvec v;
    76  list N1,N2,P,T;
    77  N1[1]=M[1][4][1];
    78  N2[1]=M[1][3][1];
    79  for(i=2;i<=size(M[1][3]);i++)
    80  {
    81   N1[i]=M[1][4][i];
    82   N2[i]=M[1][3][i]-N1[i]*M[1][3][i-1];
    83  }
    84   P[1]=N1;
    85   P[2]=N2;
    86   for(j=1;j<=size(P[1]);j++)
    87   {
    88      v=P[1][j],P[2][j];
    89      T[j]=v;
    90    }
    91   return(T);
    92 }
    93 example
    94 { "EXAMPLE:"; echo = 2;
    95    ring r=0,(x,y),dp;
    96    newtonpairs(y4+2x3y2+x6+x5y);
    97 }
    98 
    99 ///////////////////////////////////////////////////////////////////////////
    100 proc brieskornSign(a,b,c)
    101 "USAGE:   brieskornSign(a,b,c);a,b,c=integers
    102  RETURN:  signature of Brieskorn singularity x^a+y^b+z^c
    103  EXAMPLE: example brieskornSign; shows an example
    104 "
    105 {
    106    int i,j,k,d;
    107    for(i=1;i<=a-1;i++)
    108    {
    109        for(j=1;j<=b-1;j++)
    110        {
    111           for(k=1;k<=c-1;k++)
    112          {
    113             d=d+signa(number(i)/a+number(j)/b+number(k)/c);
    114          }
    115        }
    116     }
    117     return(d);
    118  }
    119 example
    120 { "EXAMPLE:"; echo = 2;
    121    ring R=0,x,dp;
    122    brieskornSign(11,3,5);
    123 }
    124 ///////////////////////////////////////////////////////////////////////////
    125 static proc signsin(number n)
    126 "USAGE:   signsin(n); n=number
    127  RETURN:  the sign of sin(n) (sin(n) = the value of the sine of n)
    128  EXAMPLE: example signSin; shows an example
    129 "
    130 {
    131    def r=basering;
    132    int a;
    133    int b=10;
    134    ring s=(real,10),x,dp;
    135    number m=imap(r,n);
    136    number pi = number_pi(11);
    137    number t=m*pi;
    138    poly f=sin(b);
    139    poly h=subst(f,x,t);
    140    if(h>0){a=1;}
    141    if(h<0){a=-1;}
    142    setring r;
    143    return(a);
    144 }
    145 example
    146 { "EXAMPLE:"; echo = 2;
    147    ring r=0,x,dp;
    148    signsin(11/3);
    149 }
    150 ///////////////////////////////////////////////////////////////////////////
    151 static proc split1(number n)
    152 "USAGE:   split1(n); n=number
    153  RETURN:  integral and fractional parts of number n
    154  EXAMPLE: example split1; shows an example
    155  "
    156 {
    157    number a,b;
    158    int r;
    159    a=numerator(n);
    160    b=denominator(n);
    161    int z=int(number(a));
    162    int y=int(number(b));
    163    r=z mod y;
    164    int q=(z-r) div y;
    165    number n1=q;
    166    number n2=n-n1;
    167    list l=n1,n2;
    168    return(l);
    169 }
    170 example
    171 { "EXAMPLE:"; echo = 2;
    172      ring r=0,x,dp;
    173      split1(11/3);
    174 }
    175 ///////////////////////////////////////////////////////////////////////////
    176 static proc sin( int n)
    177 "USAGE:   sin(n); n=integer
    178  RETURN:  approximate value of the sine of n by using maclaurin series
    179  EXAMPLE:  example sin; shows an example
    180  "
    181 {
    182    def R=basering;
    183    ring S=0,x,dp;
    184    int i;
    185    poly f;
    186    int z;
    187    for(i=1;i<=n;i=i+2)
    188    {
    189      f=f+(-1)^z*x^i/factorial(i) ;
    190      z++;
    191    }
    192    setring R;
    193    map phi=S,var(1);
    194    poly f=phi(f);
    195    return(f);
    196 }
    197 example
    198 { "EXAMPLE:"; echo = 2;
    199    ring r=0,x,dp;
    200    sin(10);
    201 }
    202 ///////////////////////////////////////////////////////////////////////////
    203 static proc cos(int n)
    204 "USAGE:   cos(n); n=integer
    205  RETURN:  approximate value of the cosine of n by using maclaurin series
    206  EXAMPLE: example cos; shows an example
    207  "
    208 {
    209    def R=basering;
    210    ring S=0,x,dp;
    211    poly f;
    212    int i;
    213    int z;
    214    for(i=0;i<=n;i=i+2)
    215    {
    216      // f=f+(-1)^z*x^i/factorial(i,0);
    217      f=f+(-1)^z*x^i/factorial(i) ;
    218      z++;
    219    }
    220    setring R;
    221    map phi=S,var(1);
    222    poly f=phi(f);
    223    return(f);
    224  }
    225 example
    226 { "EXAMPLE:"; echo = 2;
    227    ring r=0,x,dp;
    228    cos(10);
    229 }
    230 ///////////////////////////////////////////////////////////////////////////
    231 static proc signcos(number n)
    232 "USAGE:   signcos(n); n=number
    233  RETURN:  the sign of cosin(n) (cosin(n) = the value of the cosine of n)
    234  EXAMPLE: example signcos; shows an example
    235 "
    236 {
    237    def r=basering;
    238    int a;
    239    int b=10;
    240    ring s=(real,10),x,dp;
    241    number m=imap(r,n);
    242    number pi = number_pi(11);
    243    number t=m*pi;
    244    poly f=cos(b);
    245    poly h=subst(f,x,t);
    246    if(h>0){a=1;}
    247    if(h<0){a=-1;}
    248    setring r;
    249    return(a);
    250 }
    251 example
    252 { "EXAMPLE:"; echo = 2;
    253    ring r=0,x,dp;
    254    signcos(11/3);
    255 }
    256 ///////////////////////////////////////////////////////////////////////////
    257 static proc  signa( number u)
    258 "USAGE:   signa(u); u=number
    259  RETURN:  the signa of a number
    260  EXAMPLE: example signa; shows an example
    261 "
    262 {
    263    list l=split1(u);
    264    int z;
    265    if( l[1] mod 2==0 )
    266    { z=signsin(l[2]); }
    267    else
    268    { z=signcos(l[1]);}
    269    return(z);
    270 }
    271 example
    272 { "EXAMPLE:"; echo = 2;
    273    ring r=0,x,dp;
    274    signa(11/3);
    275 }
    276 ///////////////////////////////////////////////////////////////////////////
    277 static proc prods(list L)
    278 "USAGE:   product(L); L=list of intvec
    279  RETURN:  product of first components of Newton pairs in L
    280  EXAMPLE: example product; shows an example
    281 "
    282 {
    283    int a=L[2][1];
    284    int i;
    285    for(i=1;i<=size(L)-2;i++)
    286    {
    287       a=a*L[i+2][1];
    288    }
    289    return(a);
    290 }
    291 example
    292 { "EXAMPLE:"; echo = 2;
    293    list L=intvec(2,3),intvec(2,1);
    294    prods(L);
    295 }
    296 ///////////////////////////////////////////////////////////////////////////
    297 static proc Generalcase(int N, list L)
    298 "USAGE:   Generalcase(N,f);N=integer,list L of intvec
    299  RETURN:  signature of surface singularity with Newton pairs in L
    300  ASSUME:  number of newton pairs greater than 2
    301  EXAMPLE: example Generalcase; shows an example
    302 "
    303 {
    304    int i,j,k,n,m,t,p;
    305    int a=L[1][2];
    306    int b=prods(L);
    307    int d=gcd(N,b);
    308    int q=d*brieskornSign(a,L[1][1],N/d);
    309    list A;
    310    list B;
    311    list S;
    312    for(i=1;i<=size(L)-1;i++)
    313    {
    314      a=L[i+1][2]+L[i+1][1]*L[i][1]*a;
    315      A[i]=a;
    316      S[i]=L[i+1][1];
    317    }
    318    for(m=1;m<=size(L)-2;m++)
    319    {
    320      B[m]=L[m+2][1];
    321    }
    322    list C;
    323    int c=B[size(B)];
    324    C[size(B)]=c;
    325    for(j=1;j<=size(B)-1;j++)
    326    {
    327       c=c*B[size(B)-j];
    328       C[size(B)-j]=c;
    329    }
    330    list D;
    331    D[size(L)-1]=1;
    332    for(k=1;k<=size(L)-2;k++)
    333    {
    334       D[k]=gcd(N,C[k]);
    335    }
    336    for(n=1;n<=size(L)-1;n++)
    337    {
    338      t=t+D[n]*brieskornSign(A[n],S[n],N/D[n]);
    339    }
    340    return(q+t);
    341 }
    342 example
    343 { "EXAMPLE:"; echo = 2;
    344    int N=2;
    345    list L=intvec(2,3),intvec(2,1),intvec(2,1);
    346    Generalcase(N,L);
    347 }
    348 ///////////////////////////////////////////////////////////////////////////
    349 static proc signtwopairs(int N,list L)
    350 "USAGE:   signtwopairs(N,f);N=integer,L=list of intvec
    351  RETURN:  signature of surface singularity with Newton pairs in L
    352  ASSUME:  number of newton pairs equal to 2
    353  EXAMPLE: example signtwopairs; shows an example
    354 "
    355 {
    356   int a1,a2,b,d1,q,t;
    357   a1=L[1][2];
    358   b=prods(L);
    359   d1=gcd(N,b);
    360   q=d1*brieskornSign(a1,L[1][1],N/d1);
    361   a2=L[2][2]+L[2][1]*L[1][1]*a1;
    362   t=brieskornSign(a2,L[2][1],N);
    363   return(q+t);
    364 }
    365 example
    366 { "EXAMPLE:"; echo = 2;
    367    int N=2;
    368    list L=intvec(2,3),intvec(2,1);
    369    signtwopairs(N,L);
    370 }
    371 ///////////////////////////////////////////////////////////////////////////
    372 static proc DedekindSum(number b, number c, int a)
     151   poly f = x15-21x14+8x13y-6x13-16x12y+20x11y2-x12+8x11y-36x10y2
     152            +24x9y3+4x9y2-16x8y3+26x7y4-6x6y4+8x5y5+4x3y6-y8;
     153   signatureP(N,f);
     154}
     155
     156///////////////////////////////////////////////////////////////////////////////
     157//------- sigma(z^N + f) in terms of the imbedded resolution graph of f -------
     158
     159static proc dedekindSum(number b, number c, int a)
    373160{
    374161   number s,d,e;
    … …  
    385172   return(s);
    386173}
    387 ///////////////////////////////////////////////////////////////////////////
     174
     175///////////////////////////////////////////////////////////////////////////////
     176
     177static proc isRupture(intvec v)
     178//=== decides whether the exceptional divisor given by the row v in the
     179//=== incidence matrix of the resolution graph intersects at least 3 other divisors
     180{
     181   int i,j;
     182   for(i=1;i<=size(v);i++)
     183   {
     184       if(v[i]<0){return(0);}
     185       if(v[i]!=0){j++;}
     186   }
     187   return(j>=4);
     188}
     189
     190///////////////////////////////////////////////////////////////////////////////
     191
     192static proc sumExcepDiv(intmat N, list M, int K, int n)
     193//=== computes part of the formulae for eta(g,K), g defining an
     194//=== isolated curve singularity
     195//=== N the incidence matrix of the resolution graph of g
     196//=== M list of total multiplicities
     197//=== n = nrows(N)
     198{
     199   int i,j,m,d;
     200   for(i=1;i<=n;i++)
     201   {
     202      if(N[i,i]>0)
     203      {
     204         m=gcd(K,M[i]);
     205         for(j=1;j<=n;j++)
     206         {
     207            if((i!=j)&&(N[i,j]!=0))
     208            {
     209               if(m==1){break;}
     210               m=gcd(m,M[j]);
     211            }
     212         }
     213         d=d+m-1;
     214      }
     215   }
     216   return(d);
     217}
     218
     219///////////////////////////////////////////////////////////////////////////////
     220
     221static proc sumEdges(intmat N, list M, int K, int n)
     222//=== computes part of the formulae for eta(g,K), g defining an
     223//=== isolated curve singularity
     224//=== N the incidence matrix of the resolution graph of g
     225//=== M list of total multiplicities
     226//=== n = nrows(N)
     227{
     228   int i,j,d;
     229   for(i=1;i<=n-1;i++)
     230   {
     231      for(j=i+1;j<=n;j++)
     232      {
     233         if(N[i,j]==1)
     234         {
     235            d=d+gcd(K,gcd(M[i],M[j]))-1;
     236         }
     237      }
     238   }
     239   return(d);
     240}
     241
     242///////////////////////////////////////////////////////////////////////////////
     243
     244static proc etaRes(list L, int K)
     245//=== L total multiplicities
     246//=== eta-invariant in terms of the imbedded resolution graph of f
     247{
     248   int i,j,d;
     249   intvec v;
     250   number e;
     251   intmat N = L[1];         // incidence matrix of the resolution graph
     252   int n = ncols(L[1]);     // number of vertices in the resolution graph
     253   int a = ncols(L[2]);     // number of branches
     254   list M;                  // total multiplicities
     255   for(i=1;i<=n;i++)
     256   {
     257      d=L[2][i,1];
     258      for(j=2;j<=a;j++)
     259      {
     260         d=d+L[2][i,j];
     261      }
     262      if(d==0){d=1;}
     263      M[i]=d;
     264   }
     265   for(i=1;i<=n;i++)
     266   {
     267      v=N[i,1..n];
     268      if(isRupture(v))    // the divisor intersects more then two others
     269      {
     270         for(j=1;j<=n;j++)
     271         {
     272            if((i!=j)&&(v[j]!=0))
     273            {
     274               e=e+dedekindSum(M[j],K,M[i]);
     275            }
     276         }
     277      }   
     278   }
     279   if(a==1)
     280   {
     281      //the irreducible case
     282      return(4*e);
     283   }
     284   return(a-1+4*e+sumEdges(N,M,K,n)-sumExcepDiv(N,M,K,n));
     285}
     286
     287///////////////////////////////////////////////////////////////////////////////
     288//------------ sigma(z^N + f) in terms of the spectral pairs of f -------------
     289
     290static proc fracPart(number n)
     291//=== computes the fractional part n2 of n
     292//=== i.e. n2 is not in Z but n-n2 is in Z
     293{
     294   number a,b;
     295   int r;
     296   a = numerator(n);
     297   b = denominator(n);
     298   int z = int(number(a));
     299   int y = int(number(b));
     300   r = z mod y;
     301   int q = (z-r) div y;
     302   number n1 = q;
     303   number n2 = n-n1;
     304   return(n2);
     305}
     306
     307///////////////////////////////////////////////////////////////////////////////
     308
     309static proc etaSpec(list L, int N)
     310//=== L spectral numbers
     311//=== eta-invariant in terms of the spectral pairs of f
     312{
     313   int i;
     314   number e, h;
     315   
     316   int n = ncols(L[1]);
     317   
     318   if((n mod 2) == 0)               
     319   // 0 is not a spectral number, thus f is irreducible
     320   {
     321      for(i = n/2+1; i <= n; i++)
     322      {
     323         e = e + (1 - 2 * fracPart(N * number(L[1][i]))) * L[3][i];
     324      }
     325      return(2*e);
     326   }
     327   else
     328   // 0 is a spectral number, thus f is reducible
     329   {
     330      // sum of Hodge numbers in eta function
     331      for(i = 1; i <= n; i++)                   
     332      {
     333         if((L[2][i] == 2)&&((denominator(leadcoef(N*L[1][i]))==1)||(denominator(leadcoef(N*L[1][i]))==-1)))
     334         {
     335            h = h + L[3][i];
     336         }
     337      }
     338     
     339      // summand coming from spectral number 0 in eta function
     340      h = h + L[3][(n+1)/2];
     341     
     342      // sum coming from non-zero spectral numbers in eta function
     343      for(i = (n+3)/2; i <= n; i++)
     344      {
     345         if(!((denominator(leadcoef(N*L[1][i]))==1)||(denominator(leadcoef(N*L[1][i]))==-1)))
     346         {
     347            e = e + (1 - 2 * fracPart(N * number(L[1][i]))) * L[3][i];
     348         }
     349      }   
     350      return(h + 2*e);
     351   }
     352}
     353
     354///////////////////////////////////////////////////////////////////////////////
     355//---------------- Consolidation of the three former variants -----------------
     356
     357proc signature(int N, poly f, list #)
     358"USAGE:  signature(N,f); N = integer, f = reduced poly in 2 variables, # empty or 1,2,3
     359@*       - if # is empty or #[1] = 2 then resolution of singularities is used
     360@*       - if #[1] = 1 then f has to be analytically irreducible and Puiseux expansions are used
     361@*       - if #[1] = 3 then spectral pairs are used
     362RETURN:  signature of surface singularity defined by z^N + f(x,y) = 0
     363EXAMPLE: example signature; shows an example
     364"
     365{
     366   if(size(#) == 0)
     367   {
     368      list L = totalmultiplicities(f);
     369      return(etaRes(L,N) - N*etaRes(L,1));
     370   }
     371   
     372   if(#[1] == 1)
     373   {
     374      return(signatureP(N,f));
     375   }
     376   
     377   if(#[1] == 2)
     378   {
     379      list L = totalmultiplicities(f);
     380      return(etaRes(L,N) - N*etaRes(L,1));
     381   }
     382   
     383   if(#[1] == 3)
     384   {
     385      def R = basering;
     386      def Rds = changeord("ds");
     387      setring Rds;
     388      poly f = imap(R,f);
     389      list L = sppairs(f);
     390      setring R;
     391      list L = imap(Rds,L);
     392      return(etaSpec(L,N) - N*etaSpec(L,1));
     393   }
     394}
     395example
     396{ "EXAMPLE:"; echo = 2;
     397   ring r = 0,(x,y),dp;
     398   int N  = 3;
     399   poly f = x15-21x14+8x13y-6x13-16x12y+20x11y2-x12+8x11y-36x10y2
     400            +24x9y3+4x9y2-16x8y3+26x7y4-6x6y4+8x5y5+4x3y6-y8;
     401   signature(N,f,1);
     402   signature(N,f,2);
     403}
     404
     405///////////////////////////////////////////////////////////////////////////////
    388406
    389407/*
    390408Further examples
    391409
    392    ring r = 0,(x,y),dp;
    393    int N;
    394    poly f;
    395 
    396    N  = 5;
    397    poly f= x15-21x14+8x13y-6x13-16x12y+20x11y2-x12+8x11y-36x10y2    //3 characteristic pairs
    398       +24x9y3+4x9y2-16x8y3+26x7y4-6x6y4+8x5y5+4x3y6-y8;
    399 
    400    N=6;
    401    f= y4+2x3y2+x6+x5y;                                             //2 characteristic pairs
    402 
    403    N=7;
    404    f=x5+y11;                                                       //1 characteristc pair
     410ring r = 0,(x,y),dp;
     411int N;
     412poly f,g,g1,g2,g3;
     413
     414
     415// irreducible polynomials
     416
     417N = 5;
     418f = x15-21x14+8x13y-6x13-16x12y+20x11y2-x12+8x11y-36x10y2   
     419    +24x9y3+4x9y2-16x8y3+26x7y4-6x6y4+8x5y5+4x3y6-y8;     
     420g = f^3 + x17y17;                                               
     421
     422N = 6;
     423f = y4+2x3y2+x6+x5y;                                             
     424g1 = f^2 + x5y5;
     425g2 = f^3 + x11y11;
     426g3 = f^3 + x17y17;
     427
     428N = 7;
     429f = x5+y11;
     430g1 = f^3 + x11y11;
     431g2 = f^3 + x17y17;
     432
     433N = 6;     
     434// k0 = 30, k1 = 35, k2 = 71
     435f = x71+6x65+15x59-630x52y6+20x53+6230x46y6+910x39y12+15x47
     436    -7530x40y6+14955x33y12-285x26y18+6x41+1230x34y6+4680x27y12
     437    +1830x20y18+30x13y24+x35-5x28y6+10x21y12-10x14y18+5x7y24-y30;
     438
     439// k0 = 16, k1 = 24, k2 = 28, k3 = 30, k4 = 31
     440f = x31-781x30+16x29y-3010x29-2464x28y+104x27y2-2805x28-7024x27y
     441    -5352x26y2+368x25y3+366x27-7136x26y-984x25y2-8000x24y3
     442    +836x23y4+34x26-320x25y-6464x24y2+6560x23y3-8812x22y4+1392x21y5
     443    -12x25+256x24y-1296x23y2-1536x22y3+4416x21y4-8864x20y5+1752x19y6
     444    -x24+16x23y-88x22y2-16x21y3-404x20y4+3056x19y5-6872x18y6+1648x17y7
     445    +8x21y2-96x20y3+524x19y4-1472x18y5+3464x17y6-3808x16y7+1290x15y8
     446    -28x18y4+240x17y5-976x16y6+2208x15y7-2494x14y8+816x13y9+56x15y6
     447    -320x14y7+844x13y8-1216x12y9+440x11y10-70x12y8+240x11y9-344x10y10
     448    +240x9y11+56x9y10-96x8y11+52x7y12-28x6y12+16x5y13+8x3y14-y16;
     449
     450
     451// reducible polynomials
     452
     453N = 12;
     454f = ((y2-x3)^2 - 4x5y - x7)*(x2-y3);
     455
     456f = 2x3y3-2y5+x4-xy2;
     457
     458f = -x3y3+x6y+xy6-x4y4;
     459*/
Note: See TracChangeset for help on using the changeset viewer.