Changeset ea3c28a in git


Ignore:
Timestamp:
Jan 8, 2007, 1:59:38 PM (17 years ago)
Author:
Gerhard Pfister <pfister@…>
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
d1b006571ed1b6586799448402366ab672020d06
Parents:
f183dd53eec7359ae570b06c0f66b28fb8a64b3b
Message:
~2 bugs in charakteristik p gefixed


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

Legend:

Unmodified
Added
Removed
  • Singular/LIB/noether.lib

    rf183dd rea3c28a  
    1 // AH/GP last modified:  12.05.2006
     1// AH last modified:  01.07.2007
    22//////////////////////////////////////////////////////////////////////////////
    3 version = "$Id: noether.lib,v 1.3 2006-12-11 18:57:42 Singular Exp $";
     3version = "$Id: noether.lib,v 1.4 2007-01-08 12:59:38 pfister Exp $";
    44category="Commutative Algebra";
    55info="
    6 LIBRARY: noether.lib   Noether normalization of an ideal (not nessecary homogeneous)
     6LIBRARY: noether.lib   Noether normalization of an ideal (not nessecary
     7                       homogeneous)
    78AUTHORS: A. Hashemi,  Amir.Hashemi@lip6.fr
    89
    910
    1011OVERVIEW:
    11  A library for computing the Noether normalization of an ideal that DOES NOT require the computation of the dimension of the ideal.
    12 It checks whether an ideal is in Noether position.  A modular version of these algorithms is also provided.
    13 The procedures are based on a paper of Amir Hashemi 'Efficient Algorithms for Computing Noether Normalization'
    14 Submitted to: Special Issue of Mathematics in Computer Science on Symbolic and Numeric Computation.
    15 
    16 This library computes also Castelnuovo-Mumford regularity and satiety of an ideal.  A modular version of these algorithms is also provided.
    17 The procedures are based on a paper of Amir Hashemi 'Computation of  Castelnuovo-Mumford regularity and satiety'
     12A library for computing the Noether normalization of an ideal that DOES NOT
     13require the computation of the dimension of the ideal.
     14It checks whether an ideal is in Noether position.  A modular version of
     15these algorithms is also provided.
     16The procedures are based on a paper of Amir Hashemi 'Efficient Algorithms for
     17Computing Noether Normalization'
     18Submitted to: Special Issue of Mathematics in Computer Science on Symbolic
     19and Numeric Computation.
     20
     21This library computes also Castelnuovo-Mumford regularity and satiety of an
     22ideal.  A modular version of these algorithms is also provided.
     23The procedures are based on a paper of Amir Hashemi 'Computation of 
     24Castelnuovo-Mumford regularity and satiety'
    1825Submitted to: ISSAC 2007.
    1926
     
    2128PROCEDURES:
    2229 NP_test(id);  checks whether monomial ideal id is in Noether position
    23  MNP_test(id);  checks whether ideal id is in Noether position by modular methods
    24  NP(id);  Noether normalization of ideal id
    25  MNP(id);  Noether normalization of ideal id by modular methods
    26  nsatiety(id);  Satiety of ideal id
     30 MNP_test(id); the same as above using modular methods
     31 NP(id);       Noether normalization of ideal id
     32 MNP(id);      Noether normalization of ideal id by modular methods
     33 nsatiety(id); Satiety of ideal id
    2734 msatiety(id)  Satiety of ideal id by modular methods
    28  regCM(id);  Castelnuovo-Mumford regularity of ideal id
    29  mregCM(id);  Castelnuovo-Mumford regularity of ideal id by modular methods
     35 regCM(id);    Castelnuovo-Mumford regularity of ideal id
     36 mregCM(id);   Castelnuovo-Mumford regularity of ideal id by modular methods
    3037";
    3138LIB "elim.lib";
     
    3643
    3744///////////////////////////////////////////////////////////////////////////////
     45
    3846proc NP_test (ideal I)
    3947"
    40 USAGE:  NP-test (I); I monomial ideal
    41 RETURN:  A list which first element is 1 if i is in Noether position
    42   0  otherwise. The second element of this list is the list of variable which
    43   its first part is the variable such that a power of this varaibles belong to the initial of i.
    44   It return also the dimension of i if i is in Noether position
    45 ASSUME:  i is a nonzero monomial ideal.
     48USAGE:  NP-test (I); I monomial ideal
     49RETURN: A list which first element is 1 if i is in Noether position
     50        0  otherwise. The second element of this list is the list of variable which
     51        its first part is the variable such that a power of this varaibles belong to the initial of i.
     52        It return also the dimension of i if i is in Noether position
     53ASSUME: i is a nonzero monomial ideal.
    4654"
    4755{
     
    5462   if (I[1]==1)
    5563   {
    56   print("The ideal is 1");return(1);
     64        print("The ideal is 1");return(1);
    5765   }
    5866   for ( ii = 1; ii <= n+1; ii++ )
    5967     {
    60   L[ii]=0;
     68        L[ii]=0;
    6169     }
    6270   for ( ii = 1; ii <= size(I); ii++ )
    6371   {
    64     Y=findvars(I[ii],1)[1];
    65   l=rvar(Y[1][1]);
    66   if (size(Y[1])==1)
    67   {
    68     L[l]=1;
    69     P1=insert(P1,Y[1][1]);
    70   }
    71   if (L[l]==0)
    72   {
    73     L[l]=-1;
    74   }
     72        Y=findvars(I[ii],1)[1];
     73        l=rvar(Y[1][1]);
     74        if (size(Y[1])==1)     
     75        {
     76                L[l]=1;
     77                P1=insert(P1,Y[1][1]);
     78        }
     79        if (L[l]==0)
     80        {
     81                L[l]=-1;
     82        }
    7583    }
    7684   t=size(P1);
     
    7987         for ( jj = 1; jj <= n+1; jj++ )
    8088         {
    81          P3=insert(P3,varstr(jj));
    82          }
     89               P3=insert(P3,varstr(jj));
     90         }       
    8391   }
    8492   else
     
    8795         for ( jj = 1; jj <= size(P2[1]); jj++ )
    8896         {
    89          P3=insert(P3,P2[1][jj]);
     97               P3=insert(P3,P2[1][jj]);
    9098         }
    9199   }
    92100   if (L[n+1]==-1)
    93101        {
    94      return(list(0,P1+P3));
    95   }
     102                 return(list(0,P1+P3));
     103        }
    96104   for ( ii = 1; ii <= n; ii++ )
    97105     {
    98   if (L[ii]==-1)
     106        if (L[ii]==-1)
    99107        {
    100      return(list(0,P1+P3));
    101   }
    102   if (L[ii]==0 and L[ii+1]==1)
     108                 return(list(0,P1+P3));
     109        }
     110        if (L[ii]==0 and L[ii+1]==1)
    103111        {
    104      return(list(0,P1+P3));
    105   }
     112                 return(list(0,P1+P3));
     113        }
    106114    }
    107115    d=n+1-sum(L);
     
    112120//////////////////////////////////////////
    113121proc MNP_test (ideal i)
    114 "USAGE:  MNP-test (i); i an ideal
    115 RETURN:  1 if i is in Noether position 0  otherwise.
    116 NOTE:  This test is a probabilistic test, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31).
     122"USAGE: MNP-test (i); i an ideal
     123RETURN: 1 if i is in Noether position 0  otherwise.
     124NOTE:   This test is a probabilistic test, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31).     
    117125"
    118126{
     
    145153"
    146154USAGE:   NP (i); i ideal
    147 RETURN:  A linear map phi such that  phi(i) is in Noether position
     155RETURN:  A linear map phi such that  phi(i) is in Noether position 
    148156"
    149157{
     
    168176   if ( #[1]== 1 )
    169177     {
    170   return ("The ideal is in Noether position and the time of this computation is:",rtimer-time,"/10 sec.");
     178        return ("The ideal is in Noether position and the time of this computation is:",rtimer-time,"/10 sec.");
    171179     }
    172180   else
    173181     {
    174     L=maxideal(1);
    175     chcoord=maxideal(1);
    176      for ( ii = 1; ii<=n+1; ii++ )
    177     {
    178        chcoord[rvar(#[2][ii])]=L[ii];
    179     }
    180     phi=r1,chcoord;
     182                L=maxideal(1);
     183                chcoord=maxideal(1);
     184                for ( ii = 1; ii<=n+1; ii++ )
     185                {
     186                        chcoord[rvar(#[2][ii])]=L[ii];
     187                }
     188                phi=r1,chcoord;
    181189                sbi=phi(sbi);
    182190                if ( NP_test(sbi)[1] == 1 )
    183               {
    184            setring r0;
    185            chcoord=fetch(r1,chcoord);
    186     return (chcoord,"and the time of this computation is:",rtimer-time,"/10 sec.");
     191                {
     192                 setring r0;
     193                 chcoord=fetch(r1,chcoord);
     194                return (chcoord,"and the time of this computation is:",rtimer-time,"/10 sec.");
    187195                }
    188196     }
     
    190198   {
    191199     nl=nl+1;
    192      I=i;
     200     I=i;       
    193201     L=maxideal(1);
    194202     for ( ii = n; ii>=0; ii-- )
     
    201209       for ( jj = 1; jj<=ii+1; jj++ )
    202210       {
    203            P=P+ran[1,jj]*m[jj];
     211                P=P+ran[1,jj]*m[jj];
    204212       }
    205213       chcoord[ii+1]=P;
    206   L[ii+1]=P;
     214        L[ii+1]=P;
    207215       P=0;
    208216       phi=r1,chcoord;
    209217       I=phi(I);
    210218       if ( NP_test(sort(lead(std(I)))[1])[1] == 1 )
    211          {
    212     K=x(ii..n);
    213            setring r0;
    214            K=fetch(r1,K);
    215                  ideal L=fetch(r1,L);
    216     return (L,"and the time of this computation is:",rtimer-time,"/10 sec.");
     219           {
     220                K=x(ii..n);
     221                 setring r0;
     222                 K=fetch(r1,K);
     223                 ideal L=fetch(r1,L); 
     224                return (L,"and the time of this computation is:",rtimer-time,"/10 sec.");
    217225           }
    218226
     
    226234///////////////////////////////////////////////////////////////////////////////
    227235proc MNP (ideal i)
    228 "USAGE:  MNP (i); i ideal
    229 RETURN:  A linear map phi such that  phi(i) is in Noether position
    230 NOTE:  It uses the procedure  MNP_test to test Noether position.
     236"USAGE: MNP (i); i ideal
     237RETURN: A linear map phi such that  phi(i) is in Noether position
     238NOTE:   It uses the procedure  MNP_test to test Noether position.       
    231239"
    232240{
     
    246254   time=rtimer;
    247255   system("--ticks-per-sec",10);
    248     #=MNP_test(i);
     256   #=MNP_test(i);
    249257   if ( #[1]== 1 )
    250258     {
    251   return ("The ideal is in Noether position and the time of this computation is:",rtimer-time,"/10 sec.");
     259        return ("The ideal is in Noether position and the time of this computation is:",rtimer-time,"/10 sec.");
    252260     }
    253261   else
    254262     {
    255     L=maxideal(1);
    256     chcoord=maxideal(1);
    257      for ( ii = 1; ii<=n+1; ii++ )
    258     {
    259        chcoord[rvar(#[2][ii])]=L[ii];
    260     }
    261     phi=r1,chcoord;
     263                L=maxideal(1);
     264                chcoord=maxideal(1);
     265                for ( ii = 1; ii<=n+1; ii++ )
     266                {
     267                        chcoord[rvar(#[2][ii])]=L[ii];
     268                }
     269                phi=r1,chcoord;
    262270                I=phi(i);
    263271                if ( MNP_test(I)[1] == 1 )
    264               {
    265            setring r0;
    266            chcoord=fetch(r1,chcoord);
    267     return (chcoord,"and the time of this computation is:",rtimer-time,"/10 sec.");
     272                {
     273                 setring r0;
     274                 chcoord=fetch(r1,chcoord);
     275                return (chcoord,"and the time of this computation is:",rtimer-time,"/10 sec.");
    268276                }
    269277     }
     
    271279   {
    272280     nl=nl+1;
    273      I=i;
     281     I=i;       
    274282     L=maxideal(1);
    275283     for ( ii = n; ii>=0; ii-- )
     
    282290       for ( jj = 1; jj<=ii+1; jj++ )
    283291       {
    284            P=P+ran[1,jj]*m[jj];
     292                P=P+ran[1,jj]*m[jj];
    285293       }
    286294       chcoord[ii+1]=P;
     
    290298       I=phi(I);
    291299       if ( MNP_test(I)[1] == 1 )
    292          {
    293     K=x(ii..n);
    294            setring r0;
    295            K=fetch(r1,K);
    296                  ideal L=fetch(r1,L);
    297     return (L,"and the time of this computation is:",rtimer-time,"/10 sec.");
     300           {
     301                K=x(ii..n);
     302                 setring r0;
     303                 K=fetch(r1,K);
     304                 ideal L=fetch(r1,L); 
     305                return (L,"and the time of this computation is:",rtimer-time,"/10 sec.");
    298306           }
    299307
     
    309317proc Test (ideal i)
    310318"
    311 USAGE:   Test (i); i a monomial ideal,
     319USAGE:   Test (i); i a monomial ideal, 
    312320RETURN:  1 if the last variable is in generic position for i and 0 otherwise.
    313 THEORY:  The last variable is in generic position if the quotient of the ideal
    314   with respect to this variable is equal to the quotient of the ideal with respect to the maximal ideal.
     321THEORY:  The last variable is in generic position if the quotient of the ideal 
     322        with respect to this variable is equal to the quotient of the ideal with respect to the maximal ideal.
    315323"
    316324{
     
    329337   if (size(reduce(I,i)) <> 0)
    330338   {
    331   ret=0;
     339        ret=0;
    332340   }
    333341return(ret);
     
    338346proc nsatiety (ideal i)
    339347"
    340 USAGE:   nsatiety (i); i ideal,
     348USAGE:   nsatiety (i); i ideal, 
    341349RETURN:  an integer, the satiety of i.
    342350         (returns -1 if i is not homogeneous)
     
    344352THEORY:  The satiety, or saturation index, of a homogeneous ideal i is the
    345353         least integer s such that, for all d>=s, the degree d part of the
    346          ideals i and isatiety=satiety(i,maxideal(1))[1] coincide.
     354         ideals i and isat=sat(i,maxideal(1))[1] coincide.
    347355"
    348356{
     
    365373   {
    366374       dbprint(2,"The ideal is not homogeneous, and time for this test is: " + string(rtimer-time) + "/100sec.");
    367   return ();
     375        return ();
    368376   }
    369377   I=simplify(lead(sbi),1);
     
    372380   if (size(K) == 0)
    373381   {
    374   dbprint(2,"satiety(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");
    375   return();
     382        dbprint(2,"sat(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");
     383        return();
    376384   }
    377385   if (Test(I) == 1 )
    378386   {
    379   dbprint(2,"satiety(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
    380   return();
     387        dbprint(2,"sat(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
     388        return();
    381389   }
    382390   while ( nl < 5 )
    383391   {
    384      nl=nl+1;
     392     nl=nl+1;   
    385393     chcoord=select1(maxideal(1),1,(n));
    386394     ran=random(100,1,n);
     
    390398     for ( jj = 1; jj<=n+1; jj++ )
    391399       {
    392            P=P+ran[1,jj]*m[jj];
     400                P=P+ran[1,jj]*m[jj];
    393401       }
    394402     chcoord[n+1]=P;
     
    398406     I=simplify(lead(L),1);
    399407     attrib(I,"isSB",1);
    400      K=select(I,n+1);
     408     K=select(I,n+1); 
    401409     if (size(K) == 0)
    402410     {
    403   dbprint(2,"satiety(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");
    404   return();
     411        dbprint(2,"sat(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");
     412        return();
    405413     }
    406414     if (Test(I) == 1 )
    407415     {
    408   dbprint(2,"satiety(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
    409   return();
     416        dbprint(2,"sat(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
     417        return();
    410418     }
    411419   }
     
    416424proc msatiety (ideal i)
    417425"
    418 USAGE:   msatiety (i); i ideal,
     426USAGE:   msatiety (i); i ideal, 
    419427RETURN:  an integer, the satiety of i.
    420428         (returns -1 if i is not homogeneous)
     
    422430THEORY:  The satiety, or saturation index, of a homogeneous ideal i is the
    423431         least integer s such that, for all d>=s, the degree d part of the
    424          ideals i and isatiety=satiety(i,maxideal(1))[1] coincide.
    425 NOTE:  This is a probabilistic procedure, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31).
     432         ideals i and isat=sat(i,maxideal(1))[1] coincide.
     433NOTE:    This is a probabilistic procedure, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31).
    426434"
    427435{
    428436//--------------------------- initialisation ---------------------------------
    429        "// WARNING:
    430 // The procedure is probabilistic and  it computes the initial of the ideal modulo the prime number 2147483647";
     437       "// WARNING: The characteristic of base field must be zero.
     438// The procedure is probabilistic and  it computes the
     439//initial ideals modulo the prime number 2147483647.";
    431440   int e,ii,jj,h,d,time,lastv,nl,ret,s1,d1,siz,j,si,u,k,p;
    432441   intvec v1;
     
    449458       "// WARNING: The ideal is not homogeneous.";
    450459       dbprint(2,"Time for this test is: " + string(rtimer-time) + "/100sec.");
    451   return ();
     460        return ();
    452461   }
    453462   option(redSB);
     
    467476   if (size(K) == 0)
    468477   {
    469   dbprint(2,"msatiety(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");
    470   return();
     478        dbprint(2,"msat(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");
     479        return();
    471480   }
    472481   if (Test(I) == 1 )
    473482   {
    474   dbprint(2,"msatiety(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
    475   return();
     483        dbprint(2,"msat(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
     484        return();
    476485   }
    477486   while ( nl < 30 )
    478487   {
    479      nl=nl+1;
     488     nl=nl+1;   
    480489     chcoord=select1(maxideal(1),1,(n));
    481490     ran=random(100,1,n);
     
    485494     for ( jj = 1; jj<=n+1; jj++ )
    486495       {
    487            P=P+ran[1,jj]*m[jj];
     496                P=P+ran[1,jj]*m[jj];
    488497       }
    489498     chcoord[n+1]=P;
     
    502511     lsbi1=lead(sbi);
    503512     attrib(lsbi1,"isSB",1);
    504      K=select(lsbi1,n+1);
     513     K=select(lsbi1,n+1); 
    505514     if (size(K) == 0)
    506515     {
    507   dbprint(2,"msatiety(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");
    508   return();
     516        dbprint(2,"msat(i)=0 and the time of this computation: " + string(rtimer-time) + "/100sec.");
     517        return();
    509518     }
    510519     if (Test(lsbi1) == 1 )
    511520     {
    512   dbprint(2,"msatiety(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
    513   return();
     521        dbprint(2,"msat(i)=" + string(maxdeg1(K)) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
     522        return();
    514523     }
    515524   }
     
    518527//////////////////////////////////////////////////////////////////////////////
    519528//
    520 proc regCM (ideal i)
    521 "
    522 USAGE:  regCM (i); i ideal
    523 RETURN:  the Castelnuovo-Mumford regularity of i.
     529proc reg (ideal i)
     530"
     531USAGE:  reg (i); i ideal
     532RETURN: the Castelnuovo-Mumford regularity of i.
    524533         (returns -1 if i is not homogeneous)
    525534ASSUME:  i is a homogeneous ideal.
     
    527536{
    528537//--------------------------- initialisation ---------------------------------
    529    int e,ii,jj,H,h,d,time,lastv,nv,nl;
    530    int lastind,ch,nesttest,NPtest,nl,N,acc;
    531    intmat ran;
     538   int e,ii,jj,H,h,d,time,nl;
    532539   def r0 = basering;
    533540   int n = nvars(r0)-1;
    534541   string s = "ring r1 = ",charstr(r0),",x(0..n),dp;";
    535542   execute(s);
    536    ideal i,sbi,I,J,K,chcoord,m,L;
     543   ideal i,sbi,I,J,K,L;
    537544   list #;
    538545   poly P;
     
    540547   i = fetch(r0,i);
    541548   time=rtimer;
     549   system("--ticks-per-sec",100);
    542550   sbi=std(i);
    543551//----- Check ideal homogeneous
     
    550558   attrib(I,"isSB",1);
    551559   d=dim(I);
     560   if (char(r1) > 0 and d == 0)
     561   {
     562        def r2=changechar("0",r1);
     563        setring r2;
     564        ideal sbi,I,i,K,T;
     565        map phi;
     566        I = fetch(r1,I);
     567        i=I;
     568        attrib(I,"isSB",1);
     569   }
     570   else
     571   {
     572        def r2=changechar(charstr(r1),r1);
     573        setring r2;
     574        ideal sbi,I,i,K,T,ic,Ic;
     575        map phi;
     576        I = imap(r1,I);
     577        Ic=I;
     578        attrib(I,"isSB",1);
     579        i = imap(r1,i);
     580        ic=i;
     581   }
    552582   K=select(I,n+1);
    553583   if (size(K) == 0)
    554584   {
    555   h=0;
    556    }
    557    if (Test(I) == 1)
    558    {
    559   h=maxdeg1(K);
     585        h=0;
    560586   }
    561587   else
    562588   {
    563   while ( nl < 5 )
    564   {
    565     nl=nl+1;
    566          chcoord=select1(maxideal(1),1,(n));
    567          ran=random(100,1,n);
    568          ran=intmat(ran,1,n+1);
    569          ran[1,n+1]=1;
    570          m=select1(maxideal(1),1,(n+1));
    571          for ( jj = 1; jj<=n+1; jj++ )
    572            {
    573              P=P+ran[1,jj]*m[jj];
    574            }
    575          chcoord[n+1]=P;
    576          P=0;
    577          phi=r1,chcoord;
    578     i=phi(i);
    579          I=simplify(lead(std(i)),1);
    580          attrib(I,"isSB",1);
    581          K=select(I,n+1);
    582          if (size(K) == 0)
    583          {
    584       h=0;break
    585          }
    586          if (Test(I) == 1 )
    587          {
    588       h=maxdeg1(K);break;
    589     }
    590   }
    591     }
     589        if (Test(I) == 1)
     590        {
     591                h=maxdeg1(K);
     592        }
     593        else
     594        {       
     595                while ( nl < 30 )
     596                {
     597                        nl=nl+1;
     598                        phi=r2,randomLast(100);
     599                        T=phi(i);
     600                        I=simplify(lead(std(T)),1);
     601                        attrib(I,"isSB",1);
     602                        K=select(I,n+1);
     603                        if (size(K) == 0)
     604                        {
     605                                h=0;break
     606                        }
     607                        if (Test(I) == 1 )
     608                        {
     609                                h=maxdeg1(K);break;
     610                        }
     611                }
     612                i=T;
     613        }
     614   }
    592615   for ( ii = n; ii>=n-d+1; ii-- )
    593616   {
    594   i=subst(i,x(ii),0);
    595   s = "ring mr = ",charstr(r0),",x(0..ii-1),dp;";
    596   execute(s);
    597   ideal i,sbi,I,J,K,chcoord,m,L;
    598      poly P;
    599      map phi;
    600   i=imap(r1,i);
    601   sbi=std(i);
    602         I=simplify(lead(sbi),1);
    603   attrib(I,"isSB",1);
    604        K=select(I,ii);
    605        if (size(K) == 0)
    606        {
    607     H=0;
    608        }
    609   if (Test(I) == 1)
    610   {
    611     H=maxdeg1(K);
    612      }
    613      else
    614      {
    615     while ( nl < 5 )
    616     {
    617       nl=nl+1;
    618            chcoord=select1(maxideal(1),1,(ii-1));
    619            ran=random(100,1,ii-1);
    620            ran=intmat(ran,1,ii);
    621            ran[1,ii]=1;
    622            m=select1(maxideal(1),1,(ii));
    623            for ( jj = 1; jj<=ii; jj++ )
    624              {
    625                P=P+ran[1,jj]*m[jj];
    626              }
    627            chcoord[ii]=P;
    628            P=0;
    629            phi=mr,chcoord;
    630            i=phi(i);
    631            I=simplify(lead(std(i)),1);
    632            attrib(I,"isSB",1);
    633            K=select(I,ii);
    634            if (size(K) == 0)
    635            {
    636         H=0;break;
    637            }
    638            if (Test(I) == 1 )
    639            {
    640         H=maxdeg1(K);break;
    641       }
    642     }
    643     setring r1;
    644     i=imap(mr,i);
    645     kill mr;
    646       }
    647   if (H > h)
    648   {
    649     h=H;
    650   }
    651    }
    652 dbprint(2,"regCM(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
    653 return();
    654 }
     617        i=subst(i,x(ii),0);
     618        s = "ring mr = ",charstr(r1),",x(0..ii-1),dp;";
     619        execute(s);
     620        ideal i,sbi,I,J,K,L,T;
     621        poly P;
     622        map phi;
     623        i=imap(r2,i);
     624        I=simplify(lead(std(i)),1);
     625        attrib(I,"isSB",1);
     626        K=select(I,ii);
     627        if (size(K) == 0)
     628        {
     629                H=0;
     630        }
     631        else
     632        {
     633                if (Test(I) == 1)
     634                {
     635                        H=maxdeg1(K);
     636                }
     637                else
     638                {
     639                        while ( nl < 30 )
     640                        {
     641                                nl=nl+1;       
     642                                phi=mr,randomLast(100);
     643                                T=phi(i);
     644                                I=simplify(lead(std(T)),1);
     645                                attrib(I,"isSB",1);
     646                                K=select(I,ii);
     647                                if (size(K) == 0)
     648                                {
     649                                        H=0;break;
     650                                }
     651                                if (Test(I) == 1 )
     652                                {
     653                                        H=maxdeg1(K);break;
     654                                }
     655                        }
     656                        setring r2;
     657                        i=imap(mr,T);
     658                        kill mr;
     659                }
     660        }
     661        if (H > h)
     662        {
     663                h=H;
     664        }
     665   }
     666   if (nl < 30)
     667   {   
     668        dbprint(2,"reg(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + " sec./100");
     669        return();
     670   }
     671   else
     672   {
     673        I=Ic;
     674        attrib(I,"isSB",1);
     675        i=ic;
     676        K=subst(select(I,n+1),x(n),1);
     677        K=K*maxideal(maxdeg1(I));
     678        if (size(reduce(K,I)) <> 0)
     679        {
     680                nl=0;
     681                while ( nl < 30 )
     682                {
     683                        nl=nl+1;
     684                        phi=r1,randomLast(100);
     685                        sbi=phi(i);
     686                        I=simplify(lead(std(sbi)),1);
     687                        attrib(I,"isSB",1);
     688                        K=subst(select(I,n+1),x(n),1);
     689                        K=K*maxideal(maxdeg1(I));
     690                        if (size(reduce(K,I)) == 0)
     691                        {
     692                                break;
     693                        }
     694                }
     695        }       
     696        h=maxdeg1(simplify(reduce(quotient(I,maxideal(1)),I),2))+1;
     697        for ( ii = n; ii> n-d+1; ii-- )
     698        {
     699                sbi=subst(sbi,x(ii),0);
     700                s = "ring mr = ",charstr(r0),",x(0..ii-1),dp;";
     701                execute(s);
     702                ideal sbi,I,L,K,T;
     703                map phi;
     704                sbi=imap(r1,sbi);
     705                I=simplify(lead(std(sbi)),1);
     706                attrib(I,"isSB",1);
     707                K=subst(select(I,ii),x(ii-1),1);
     708                K=K*maxideal(maxdeg1(I));
     709                if (size(reduce(K,I)) <> 0)
     710                {
     711                        nl=0;
     712                        while ( nl < 30 )
     713                        {
     714                                nl=nl+1;
     715                                L=randomLast(100);
     716                                phi=mr,L;
     717                                T=phi(sbi);
     718                                I=simplify(lead(std(T)),1);
     719                                attrib(I,"isSB",1);
     720                                K=subst(select(I,ii),x(ii-1),1);
     721                                K=K*maxideal(maxdeg1(I));
     722                                if (size(reduce(K,I)) == 0)
     723                                {
     724                                        sbi=T;
     725                                        break;
     726                                }
     727                        }       
     728                }
     729                H=maxdeg1(simplify(reduce(quotient(I,maxideal(1)),I),2))+1;
     730                if (H > h)
     731                {
     732                        h=H;
     733                }
     734                setring r1;
     735                sbi=fetch(mr,sbi);
     736                kill mr;
     737        }
     738   sbi=subst(sbi,x(n-d+1),0);
     739   s = "ring mr = ",charstr(r0),",x(0..n-d),dp;";       
     740   execute(s);
     741   ideal sbi,I,L,K,T;
     742   map phi;
     743   sbi=imap(r1,sbi);
     744   I=simplify(lead(std(sbi)),1);
     745   attrib(I,"isSB",1);
     746   H=maxdeg1(simplify(reduce(quotient(I,maxideal(1)),I),2))+1;
     747   if (H > h)
     748   {
     749        h=H;
     750   }
     751   dbprint(2,"reg(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + " sec./100");
     752   return();
     753   }
     754}
     755
    655756//////////////////////////////////////////////////////////////////////////////
    656757//
    657 proc mregCM (ideal i)
    658 "
    659 USAGE:  mregCM (i); i ideal
    660 RETURN:  an integer, the Castelnuovo-Mumford regularity of i.
     758proc mreg (ideal i)
     759"
     760USAGE:   mreg (i); i ideal
     761RETURN:  an integer, the Castelnuovo-Mumford regularity of i.
    661762         (returns -1 if i is not homogeneous)
    662 ASSUME:  i is a homogeneous ideal.
    663 NOTE:  This is a probabilistic procedure, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31).
     763ASSUME:  i is a homogeneous ideal and the characteristic of base field is zero..
     764NOTE:    This is a probabilistic procedure, and it computes the initial of the ideal modulo the prime number 2147483647 (the biggest prime less than 2^31).
    664765"
    665766{
    666767//--------------------------- initialisation ---------------------------------
    667        "// WARNING:
    668 // The procedure is probabilistic and  it computes the initial of the ideal modulo the prime number 2147483647";
    669    int e,ii,jj,H,h,d,time,lastv,sat,firstind,nv,nl,p;
    670    int lastind,ch,nesttest,NPtest,nl,N,acc;
    671    intmat ran;
     768       "// WARNING: The characteristic of base field musr be zero.
     769// This procedure is probabilistic and  it computes the initial
     770//ideals modulo the prime number 2147483647";
     771   int e,ii,jj,H,h,d,time,p,nl;
    672772   def r0 = basering;
    673773   int n = nvars(r0)-1;
    674774   string s = "ring r1 = ",charstr(r0),",x(0..n),dp;";
    675775   execute(s);
    676    ideal i,sbi,I,J,K,chcoord,m,L,lsbi1,lsbi2;
     776   ideal i,sbi,I,J,K,L,lsbi1,lsbi2;
    677777   list #;
    678778   poly P;
     
    680780   i = fetch(r0,i);
    681781   time=rtimer;
     782   system("--ticks-per-sec",100);
    682783//----- Check ideal homogeneous
    683784   if ( homog(i) == 0 )
     
    702803   I=lsbi1;
    703804   d=dim(I);
    704    K=select(I,n+1);
     805   K=select(I,n+1); 
    705806   if (size(K) == 0)
    706807   {
    707   h=0;
    708    }
    709    if (Test(I) == 1)
    710    {
    711   h=maxdeg1(K);
     808        h=0;
    712809   }
    713810   else
    714811   {
    715   while ( nl < 5 )
    716   {
    717     nl=nl+1;
    718          chcoord=select1(maxideal(1),1,(n));
    719          ran=random(100,1,n);
    720          ran=intmat(ran,1,n+1);
    721          ran[1,n+1]=1;
    722          m=select1(maxideal(1),1,(n+1));
    723          for ( jj = 1; jj<=n+1; jj++ )
    724            {
    725              P=P+ran[1,jj]*m[jj];
    726            }
    727          chcoord[n+1]=P;
    728          P=0;
    729          phi=r1,chcoord;
    730     i=phi(i);
    731        #=ringlist(r1);
    732        #[1]=p;
    733        def oro=ring(#);
    734        setring oro;
    735        ideal sbi,lsbi;
    736        sbi=fetch(r1,i);
    737        lsbi=lead(std(sbi));
    738        setring r1;
    739        lsbi1=fetch(oro,lsbi);
    740        lsbi1=simplify(lsbi1,1);
    741        attrib(lsbi1,"isSB",1);
    742        kill oro;
    743     I=lsbi1;
    744          K=select(I,n+1);
    745          if (size(K) == 0)
    746          {
    747       h=0;break
    748          }
    749          if (Test(I) == 1 )
    750          {
    751       h=maxdeg1(K);break;
    752     }
    753   }
    754     }
     812        if (Test(I) == 1)
     813        {
     814                h=maxdeg1(K);
     815        }
     816        else
     817        {
     818                while ( nl < 30 )
     819                {
     820                        nl=nl+1;       
     821                        phi=r1,randomLast(100);
     822                        sbi=phi(i);
     823                        #=ringlist(r1);
     824                        #[1]=p;
     825                        def oro=ring(#);
     826                        setring oro;
     827                        ideal sbi,lsbi;
     828                        sbi=fetch(r1,sbi);
     829                        lsbi=lead(std(sbi));
     830                        setring r1;
     831                        lsbi1=fetch(oro,lsbi);
     832                        lsbi1=simplify(lsbi1,1);
     833                        attrib(lsbi1,"isSB",1);
     834                        kill oro;
     835                        I=lsbi1;
     836                        K=select(I,n+1);
     837                        if (size(K) == 0)
     838                        {
     839                                h=0;break
     840                        }
     841                        if (Test(I) == 1 )
     842                        {
     843                                h=maxdeg1(K);break;
     844                        }
     845                }
     846                i=sbi;
     847        }
     848   }
    755849   for ( ii = n; ii>=n-d+1; ii-- )
    756850   {
    757   i=subst(i,x(ii),0);
    758   s = "ring mr = ",charstr(r0),",x(0..ii-1),dp;";
    759   execute(s);
    760   ideal i,sbi,I,J,K,chcoord,m,L,lsbi1;
    761      poly P;
     851        i=subst(i,x(ii),0);
     852        s = "ring mr = ","0",",x(0..ii-1),dp;";
     853        execute(s);
     854        ideal i,sbi,I,J,K,L,lsbi1;
     855        poly P;
    762856        list #;
    763      map phi;
    764   i=imap(r1,i);
    765      #=ringlist(mr);
    766      #[1]=p;
    767      def oro=ring(#);
    768      setring oro;
    769      ideal sbi,lsbi;
    770      sbi=fetch(mr,i);
    771      lsbi=lead(std(sbi));
    772      setring mr;
    773      lsbi1=fetch(oro,lsbi);
    774      lsbi1=simplify(lsbi1,1);
    775      attrib(lsbi1,"isSB",1);
    776      kill oro;
    777   I=lsbi1;
    778        K=select(I,ii);
    779        if (size(K) == 0)
    780        {
    781     H=0;
    782        }
    783   if (Test(I) == 1)
    784   {
    785     H=maxdeg1(K);
    786      }
    787      else
    788      {
    789     while ( nl < 5 )
    790     {
    791       nl=nl+1;
    792            chcoord=select1(maxideal(1),1,(ii-1));
    793            ran=random(100,1,ii-1);
    794            ran=intmat(ran,1,ii);
    795            ran[1,ii]=1;
    796            m=select1(maxideal(1),1,(ii));
    797            for ( jj = 1; jj<=ii; jj++ )
    798              {
    799                P=P+ran[1,jj]*m[jj];
    800              }
    801            chcoord[ii]=P;
    802            P=0;
    803            phi=mr,chcoord;
    804            i=phi(i);
    805          #=ringlist(mr);
    806          #[1]=p;
    807          def oro=ring(#);
    808          setring oro;
    809          ideal sbi,lsbi;
    810          sbi=fetch(mr,i);
    811          lsbi=lead(std(sbi));
    812          setring mr;
    813          lsbi1=fetch(oro,lsbi);
    814          lsbi1=simplify(lsbi1,1);
    815          attrib(lsbi1,"isSB",1);
    816          kill oro;
    817       I=lsbi1;
    818            K=select(I,ii);
    819            if (size(K) == 0)
    820            {
    821         H=0;break;
    822            }
    823            if (Test(I) == 1 )
    824            {
    825         H=maxdeg1(K);break;
    826       }
    827     }
    828     setring r1;
    829     i=imap(mr,i);
    830     kill mr;
    831       }
    832   if (H > h)
    833   {
    834     h=H;
    835   }
    836   }
    837 dbprint(2,"mregCM(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + "/100sec.");
     857        map phi;
     858        i=imap(r1,i);
     859        #=ringlist(mr);
     860        #[1]=p;
     861        def oro=ring(#);
     862        setring oro;
     863        ideal sbi,lsbi;
     864        sbi=fetch(mr,i);
     865        lsbi=lead(std(sbi));
     866        setring mr;
     867        lsbi1=fetch(oro,lsbi);
     868        lsbi1=simplify(lsbi1,1);
     869        attrib(lsbi1,"isSB",1);
     870        kill oro;
     871        I=lsbi1;
     872        K=select(I,ii);
     873        if (size(K) == 0)
     874        {
     875                H=0;
     876        }
     877        else
     878        {
     879                if (Test(I) == 1)
     880                {
     881                        H=maxdeg1(K);
     882                }
     883                else
     884                {
     885                        nl=0;
     886                        while ( nl < 30 )
     887                        {
     888                                nl=nl+1;       
     889                                phi=mr,randomLast(100);
     890                                sbi=phi(i);
     891                                #=ringlist(mr);
     892                                #[1]=p;
     893                                def oro=ring(#);
     894                                setring oro;
     895                                ideal sbi,lsbi;
     896                                sbi=fetch(mr,sbi);
     897                                lsbi=lead(std(sbi));
     898                                setring mr;
     899                                lsbi1=fetch(oro,lsbi);
     900                                lsbi1=simplify(lsbi1,1);
     901                                kill oro;
     902                                I=lsbi1;
     903                                attrib(I,"isSB",1);
     904                                K=select(I,ii);
     905                                if (size(K) == 0)
     906                                {
     907                                        H=0;break;
     908                                }
     909                                if (Test(I) == 1 )
     910                                {
     911                                        H=maxdeg1(K);break;
     912                                }
     913                        }
     914                        setring r1;
     915                        i=imap(mr,sbi);
     916                        kill mr;
     917                }
     918        }
     919        if (H > h)
     920        {
     921                h=H;
     922        }
     923   }
     924dbprint(2,"mreg(i)=" + string(h) + " and the time of this computation: " + string(rtimer-time) + "sec./100");
    838925return();
    839926}
     927
     928
     929
     930
     931
     932
     933//////////////////////////////////////////////////////////////
    840934example
    841935{ "EXAMPLE:"; echo = 2;
     
    9171011}
    9181012
    919 
Note: See TracChangeset for help on using the changeset viewer.