Changeset 89b9ee in git


Ignore:
Timestamp:
Jan 12, 2014, 8:47:25 PM (10 years ago)
Author:
Janko Boehm <boehm@…>
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
f84dfa3c651996cde0306f1b733e44e9aa69055b
Parents:
04498979b0db3e78f460ad02a00810c294c537f8
Message:
reesclos.lib rewritten to use normal.lib
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • Singular/LIB/reesclos.lib

    r044989 r89b9ee  
    66LIBRARY:     reesclos.lib   PROCEDURES TO COMPUTE THE INT. CLOSURE OF AN IDEAL
    77AUTHOR:      Tobias Hirsch, email: hirsch@math.tu-cottbus.de
     8             Janko Boehm, email: boehm@mathematik.uni-kl.de
     9             Magdaleen Marais, email: magdaleen@aims.ac.za
    810
    911OVERVIEW:
    1012 A library to compute the integral closure of an ideal I in a polynomial ring
    1113 R=k[x(1),...,x(n)] using the Rees Algebra R[It] of I. It computes the integral
    12  closure of R[It] (in the same manner as done in the library 'normal.lib'),
     14 closure of R[It],
    1315 which is a graded subalgebra of R[t]. The degree-k-component is the integral
    1416 closure of the k-th power of I.
     17
     18 In contrast to the previous version, the library uses 'normal.lib' to compute the
     19 integral closure of R[It]. This improves the performance considerably.
    1520
    1621PROCEDURES:
    … …  
    118123  setring R(1);                   // declaration of variables used later
    119124  ideal ker(1)=ker;               // in STEP 2
    120   poly p;
    121 
    122   // some auxiliary variables
    123 
    124   int i=1;         // counts the number of steps to reach the closure of R(1)
    125   int check=0;     // a 'boolean' variable for several checks
    126 
    127   /////// STEP 1: ///////////////////////////////////////////////////////////
    128   // construct R(i) step by step as done in normal.lib; 2 differences:     //
    129   //  - since the input is a prime ideal, no splitting will occur          //
    130   //  - the intermediate rings and nonzerodivisors for J are remembered    //
    131   ///////////////////////////////////////////////////////////////////////////
    132 
    133   if (dblvl>0)
    134   {
    135     "// STEP 1: Compute the integral closure of R[It]";
    136   }
    137 
    138   list RS;
    139 
    140   while (check==0)  // repeat until the closure is reached
    141   {
    142     // construction of HomJJ, J an ideal containing the non-normal locus
    143     // of R(i)/ker(i), as done in normalizationPrimes in normal.lib for
    144     // the special case that we are working with a prime ideal
    145 
    146     if (dblvl>0)
    147     {
    148       "";
    149       "// We are in step",i;
    150       pause();
    151     }
    152 
    153     if (homog(ker(i))==1)
    154     {
    155       list SM=mstd(ker(i));
    156     }
    157     else
    158     {
    159       list SM=groebner(ker(i)),ker(i);
    160     }
    161 
    162     if (dblvl>0)
    163     {
    164       "// Standard basis of the current ideal:";
    165       SM[1];
    166     }
    167 
    168     // In the first iteration, we have to compute the singular locus
    169     // "from scratch". In further iterations, we can fetch it from the
    170     // previous one but have to compute its radical.
    171 
    172     if (i==1)
    173     {
    174       ideal J=minor(jacob(SM[2]),nvars(basering)-dim(SM[1]));
    175       ideal sin=J+SM[2];
    176       if (homog(sin)==1)
    177       {
    178         list JM=mstd(sin);
    179       }
    180       else
    181       {
    182         list JM=groebner(sin),sin;
    183       }
    184 
    185       J=radical(JM[2]);
    186     }
    187     else
    188     {
    189       ideal J=radical(JM[2]);
    190     }
    191 
    192     if (dblvl>0)
    193     {
    194       "// Radical of the singular locus:";
    195       J;
    196     }
    197 
    198     ideal SL=simplify(reduce(JM[2],SM[1]),2);
    199     JM =J,J;
    200     poly nzd=SL[1];           // universal denominator for HomJJ
    201 
    202     if (dblvl>0)
    203     {
    204       "// The non-zerodivisor";
    205       nzd;
    206       pause();
    207     }
    208 
    209     list RR=SM[1],SM[2],JM[2],SL[1];
    210     RS=HomJJ(RR);
    211 
    212     if (RS[2]==1)             // we've reached the integral closure
    213     {
    214       if (dblvl>0)
    215       {
    216         "";
    217         "// We've reached the integral closure after",i,"iterations";
    218         pause();
    219       }
    220       check=1;
    221     }
    222     else                      // prepare the next iteration with new
    223     {                         // ring R(i+1)
    224       ideal MJ=JM[2];
    225 
    226       def R(i+1)=RS[1];       // the data of and some variable decla-
    227       setring R(i+1);         // rations in R(i+1) needed in STEP 2
    228       ideal ker(i+1)=endid;
    229       map phi=R(i),endphi;
    230       poly p;
    231 
    232       list JM=mstd(simplify(phi(MJ)+ker(i+1),4)); // fetch the singular locus
    233       i++;
    234     }
    235   }
    236 
    237   if (i==1)                     // R[It] (and thus I) was integrally
    238   {                             // closed ==> we're already done
    239         list result="closed";
    240         return(result);
    241   }
    242 
    243 
    244   /////// STEP 2: ////////////////////////////////////////////////////////
    245   // compute representations of the ring variables of R(i) as fractions //
    246   // of elements of R(1);                                               //
    247   ////////////////////////////////////////////////////////////////////////
    248 
    249   int length=nvars(R(i));  // the number of variables of the last ring
    250   int j,k,n;               // some counters
    251   string mapstr;           // will be used while constructing preimages
    252 
    253   list preimages;          // here the fractions are stored in the
    254                            // form var(j)=preimages[j]/preimages[length+1]
    255                            // ('=' means identification via the inclusion)
    256                            // the last entry corresponds to nzd in R(i)
    257 
    258   // For each variable (counter j) and for each intermediate ring (k):
    259   // find preimages of var(j)*endphi(nzd_k-1) in R(k-1).
    260   // Finally, do the same for nzd.
    261 
    262   if (dblvl>0)
    263   {
    264     "// STEP 2: Compute fractions representing the ring variables of";
    265     "           the last ring";
    266   }
    267 
    268   for (j=1;j<=length+1;j++)
    269   {
    270     setring R(i);
    271 
    272     if (j<=length)          // do it with for ring variables...
    273     {
    274       p=var(j);
    275     }
    276     else                   // ...and finally for nzd in R(i)
    277     {
    278       p=1;
    279     }
    280 
    281     // get back from R(i) to R(1) step by step
    282 
    283     for (k=i;k>1;k--)
    284     {
    285       // clear the fraction in the representation in R(i)
    286 
    287       p=p*phi(nzd);
    288 
    289       // compute the preimage of [p mod ker(k)] under phi in R(k-1):
    290       // as p is an element of im(phi), there is a polynomial h such that
    291       // h(vars(R(k-1)) is mapped to [p mod ker (k)], and h can be com-
    292       // puted as the normal form of a w.r.t <ker(k),Z(n)-phi(k)(n)> in R(k)[Z]
    293 
    294          // compute this normal form h...
    295 
    296       if (j==1)    // in the first iteration: construct S(k)=R(k)[Z], fetch
    297                    // endphi (the ideal defining phi) and ker(k) and construct
    298                    // the ideal from above
    299       {
    300         execute ("ring S(k) = "+charstr(R(k))+",("+varstr(R(k))+",Z(1.."+string(ncols(endphi))+")),(dp("+string(nvars(R(k)))+"),dp("+string(ncols(endphi))+"));");
    301         ideal endphi = imap(R(k),endphi);
    302         ideal J = imap(R(k),ker(k));
    303         for (n=1;n<=ncols(endphi);n++)
    304         {
    305           J=J+(Z(n)-endphi[n]);
    306         }
    307         J=groebner(J);
    308         poly h=NF(imap(R(k),p),J);
    309       }
    310       else
    311       {
    312         setring S(k);
    313         h=NF(imap(R(k),p),J);
    314       }
    315 
    316          // and compute h(vars(R(k-1)))
    317 
    318       setring R(k-1);
    319       if (j==1)  // in the first iteration: compute backmap:S(k)-->R(k-1)
    320       {
    321         mapstr="map backmap = S(k),";
    322         for (n=1;n<=nvars(R(k));n++)
    323         {
    324           mapstr=mapstr+"0,";
    325         }
    326         execute (mapstr+"maxideal(1);");
    327       }
    328 
    329       p=NF(backmap(h),std(ker(k-1)));
    330     }
    331 
    332     // when down to R(1), store the result in the list preimages
    333 
    334     preimages=insert(preimages,p,j-1);
    335     if (dblvl>0)
    336     {
    337       if (j<=length)
    338       {
    339         "numerator of variable ",j,":",p;
    340       }
    341       else
    342       {
    343         "";
    344         "and finally the 'universal' denominator:",p;
    345       }
    346     }
    347   }
    348 
    349   // at the end: go back to the original basering and construct gene-
    350   // rators of the closure of I
    351 
     125  list nor = normal(ker);
     126  list preimages=nor[2];
    352127  setring Kxt;
    353128  map psi=R(1),mapI;              // from ReesAlgebra: the map Rees->Kxt
    354 
    355   list images=psi(preimages);
    356 
    357   if (dblvl>-1)
    358   {
    359     pause();
    360     "";
    361     "// Get back to the original basering and construct the";
    362     "// generators of the closure of I";
    363     "";
    364     "// Back in k[x,t], the fractions, stored in the list images:";
    365     for (int j=1;j<=size(images);j++)
    366     {
    367       if (j<size(images))
    368       {
    369         "numerator",j,":",images[j];
    370       }
    371       else
    372       {
    373         "";
    374         "denominator:  ",images[j];
    375       }
    376     }
    377     pause();
    378   }
    379   return (images);
     129  ideal images=(psi(preimages))[1];
     130  images=images[size(images)]*ideal(psi)+images;
     131  list imagesl = images[1..size(images)];
     132  return(imagesl);
    380133}
    381134

  • Tst/Manual/normalI.res.gz.uu

    r044989 r89b9ee  
    1 begin 640 normalI.res.gz
    2 M'XL("`Y?<DX``VYO<FUA;$DN<F5S`'U0P4K$,!2\YRN&Q4,+M6O35L72',1+
    3 MBGAPO<FR5#>N@9A*DL7D[TUK:6^^RPS,O#?OO=W+`W\"4#`\\GMLG'6YDF^;
    4 M!I$=I)8N21LR(AB#'LQ7KWBNQ4]N7>_(;FZG<[L1PKZKP?[-6.22P4A]PG-[
    5 ME24^"VEV_%[5BD$>1:_`T<+3S(<J"_6JUPQ*6H<NRO,&"4]7_9J!-X2_%OO6
    6 MTX@T8J@B*?=MJ!?;#4/7X+_:;M'%*9`6[E-`:B=.)NXU7G0V`L,'.(F&.Q+-
    7 DASEOHG1*FF@,O?2A7&)OV?3*\5UGFQ1I<T%^`0%*@8-U`0``
     1begin 664 normalI.res.gz
     2M'XL("/'PTE(``VYO<FUA;$DN<F5S`'U034O$,!2\YU<,BX<60M=^J1B:@WA)
     3M$0^N-UF6ZL8U$%-)LMC\>]-N:6_F,L.;>>]-WN[U43P#R#F>Q`,VWOE,J_<-
     4M0V0'991/4D9&!.<PO?WNM,B,_,V<[SS9S>W%W&ZE=!^Z=Y<9BUQR6&5.>&FN
     5M:3+0D-+CSZI6'.HH.PV!!D-!AU#14*]ZS:&5\VBC/"=(1+KJ-QR"$?&6[YNA
     6MB%A$#%4DY;X)]6*[Y6@9_GO;+=HX!<K!?TDHX^7)QESCC\Y6HO^$(-%P3Z+Y
     7E,.^;Z+QRXI>M$ZW&<KDDN./35<?+G5V2I^R*_`'F.0@I@`$`````
    88`
    99end

  • Tst/Manual/normalI.stat

    r044989 r89b9ee  
    1 1 >> tst_memory_0 :: 1316104113:3132- exportiert :3-1-3:ix86-Linux:mamawutz:333020
    2 1 >> tst_memory_1 :: 1316104113:3132- exportiert :3-1-3:ix86-Linux:mamawutz:823296
    3 1 >> tst_memory_2 :: 1316104113:3132- exportiert :3-1-3:ix86-Linux:mamawutz:888832
    4 1 >> tst_timer_1 :: 1316104113:3132- exportiert :3-1-3:ix86-Linux:mamawutz:28
     11 >> tst_memory_0 :: 1389555953:4.0.0, 32 bit:4.0.0:i686-Linux:ubuntu:354340
     21 >> tst_memory_1 :: 1389555953:4.0.0, 32 bit:4.0.0:i686-Linux:ubuntu:2461696
     31 >> tst_memory_2 :: 1389555953:4.0.0, 32 bit:4.0.0:i686-Linux:ubuntu:2588236
     41 >> tst_timer_1 :: 1389555953:4.0.0, 32 bit:4.0.0:i686-Linux:ubuntu:7
Note: See TracChangeset for help on using the changeset viewer.