////////////////////////////////////////////////////////////////////////////// version="version normal.lib 4.0.1.1 Dec_2014 "; // $Id$ category="Commutative Algebra"; info=" LIBRARY: normal.lib Normalization of Affine Rings AUTHORS: G.-M. Greuel, greuel@mathematik.uni-kl.de, @* S. Laplagne, slaplagn@dm.uba.ar, @* G. Pfister, pfister@mathematik.uni-kl.de @* Peter Chini, chini@rhrk.uni-kl.de (normalConductor) PROCEDURES: normal(I,[...]); normalization of an affine ring normalP(I,[...]); normalization of an affine ring in positive characteristic normalC(I,[...]); normalization of an affine ring through a chain of rings HomJJ(L); presentation of End_R(J) as affine ring, J an ideal genus(I); computes the geometric genus of a projective curve primeClosure(L); integral closure of R/p, p a prime ideal closureFrac(L); writes a poly in integral closure as element of Quot(R/p) iMult(L); intersection multiplicity of the ideals of the list L deltaLoc(f,S); sum of delta invariants at conjugated singular points locAtZero(I); checks whether the zero set of I is located at 0 norTest(I,nor); checks the output of normal, normalP, normalC getSmallest(J); computes the polynomial of smallest degree of J getOneVar(J, vari); computes a polynomial of J in the variable vari changeDenominator(U1, c1, c2, I); computes ideal U2 such that 1/c1*U1=1/c2*U2 normalConductor(ideal); computation of the conductor as ideal in the basering SEE ALSO: locnormal_lib;modnormal_lib "; LIB "general.lib"; LIB "poly.lib"; LIB "sing.lib"; LIB "primdec.lib"; LIB "elim.lib"; LIB "presolve.lib"; LIB "inout.lib"; LIB "ring.lib"; LIB "hnoether.lib"; LIB "reesclos.lib"; LIB "algebra.lib"; LIB "curveInv.lib"; /////////////////////////////////////////////////////////////////////////////// proc normal(ideal id, list #) "USAGE: normal(id [,choose]); id = radical ideal, choose = list of options. @* Optional parameters in list choose (can be entered in any order):@* Decomposition:@* - \"equidim\" -> computes first an equidimensional decomposition of the input ideal, and then the normalization of each component (default).@* - \"prim\" -> computes first the minimal associated primes of the input ideal, and then the normalization of each prime. (When the input ideal is not prime and the minimal associated primes are easy to compute, this method is usually faster than \"equidim\".)@* - \"noDeco\" -> no preliminary decomposition is done. If the ideal is not equidimensional radical, output might be wrong.@* - \"isPrim\" -> assumes that the ideal is prime. If this assumption does not hold, the output might be wrong.@* - \"noFac\" -> factorization is avoided in the computation of the minimal associated primes; Other:@* - \"useRing\" -> uses the original ring ordering.@* If this option is set and if the ring ordering is not global, normal will change to a global ordering only for computing radicals and prime or equidimensional decompositions.@* If this option is not set, normal changes to dp ordering and performs all computations with respect to this ordering.@* - \"withDelta\" (or \"wd\") -> returns also the delta invariants.@* If the optional parameter choose is not given or empty, only \"equidim\" but no other option is used.@* - list(\"inputJ\", ideal inputJ) -> takes as initial test ideal the ideal inputJ. This option is only for use in other procedures. Using this option, the result might not be the normalization.@* (Option only valid for global algorithm.)@* - list(\"inputC\", ideal inputC) -> takes as initial conductor the ideal inputC. This option is only for use in other procedures. Using this option, the result might not be the normalization.@* (Option only valid for global algorithm.)@* Options used for computing integral basis (over rings of two variables):@* - \"var1\" -> uses a polynomial in the first variable as universal denominator.@* - \"var2\" -> uses a polynomial in the second variable as universal denominator.@* If the optional parameter choose is not given or empty, only \"equidim\" but no other option is used.@* ASSUME: The ideal must be radical, for non-radical ideals the output may be wrong (id=radical(id); makes id radical). However, when using the \"prim\" option the minimal associated primes of id are computed first and hence normal computes the normalization of the radical of id.@* NOTE: \"isPrim\" should only be used if id is known to be prime. RETURN: a list, say nor, of size 2 (resp. 3 with option \"withDelta\"). Let R denote the basering and id the input ideal. * nor[1] is a list of r rings, where r is the number of associated primes P_i with option \"prim\" (resp. >= no of equidimenensional components P_i with option \"equidim\").@* Each ring Ri := nor[1][i], i=1..r, contains two ideals with given names @code{norid} and @code{normap} such that: @* - Ri/norid is the normalization of the i-th component, i.e. the integral closure of R/P_i in its field of fractions (as affine ring); - @code{normap} gives the normalization map from R/id to Ri/norid for each i.@* - the direct sum of the rings Ri/norid, i=1,..r, is the normalization of R/id as affine algebra; @* * nor[2] is a list of size r with information on the normalization of the i-th component as module over the basering R:@* nor[2][i] is an ideal, say U, in R such that the integral closure of basering/P_i is generated as module over R by 1/c * U, with c the last element U[size(U)] of U.@* * nor[3] (if option \"withDelta\" is set) is a list of an intvec of size r, the delta invariants of the r components, and an integer, the total delta invariant of basering/id (-1 means infinite, and 0 that R/P_i resp. R/id is normal). THEORY: We use here a general algorithm described in [G.-M.Greuel, S.Laplagne, F.Seelisch: Normalization of Rings (2009)].@* The procedure computes the R-module structure, the algebra structure and the delta invariant of the normalization of R/id:@* The normalization of R/id is the integral closure of R/id in its total ring of fractions. It is a finitely generated R-module and nor[2] computes R-module generators of it. More precisely: If U:=nor[2][i] and c:=U[size(U)], then c is a non-zero divisor and U/c is an R-module in the total ring of fractions, the integral closure of R/P_i. Since U[size(U)]/c is equal to 1, R/P_i resp. R/id is contained in the integral closure.@* The normalization is also an affine algebra over the ground field and nor[1] presents it as such. For geometric considerations nor[1] is relevant since the variety of the ideal norid in Ri is the normalization of the variety of the ideal P_i in R.@* The delta invariant of a reduced ring A is dim_K(normalization(A)/A). For A=K[x1,...,xn]/id we call this number also the delta invariant of id. nor[3] returns the delta invariants of the components P_i and of id. NOTE: To use the i-th ring type e.g.: @code{def R=nor[1][i]; setring R;}. @* Increasing/decreasing printlevel displays more/less comments (default: printlevel=0). @* Implementation works also for local rings. @* Not implemented for quotient rings. @* If the input ideal id is weighted homogeneous a weighted ordering may be used together with the useRing-option (qhweight(id); computes weights). KEYWORDS: normalization; integral closure; delta invariant. SEE ALSO: normalC, normalP. EXAMPLE: example normal; shows an example " { ASSUME(0, not isQuotientRing(basering) ); intvec opt = option(get); // Save current options int i,j; int decomp; // Preliminary decomposition: // 0 -> no decomposition (id is assumed to be prime) // 1 -> no decomposition // (id is assumed to be equidimensional radical) // 2 -> equidimensional decomposition // 3 -> minimal associated primes int noFac, useRing, withDelta; int dbg = printlevel - voice + 2; int nvar = nvars(basering); int chara = char(basering); int denomOption; // Method for choosing the conductor ideal inputJ = 0; // Test ideal given in the input (if any). ideal inputC = 0; // Conductor ideal given in the input (if any). list result, resultNew; list keepresult; list ringStruc; ideal U; poly c; int sp; // Number of components. // Default methods: noFac = 0; // Use facSTD when computing minimal associated primes decomp = 2; // Equidimensional decomposition useRing = 0; // Change first to dp ordering, and perform all // computations there. withDelta = 0; // Do not compute the delta invariant. denomOption = 0; // The default universal denominator is the smallest // degree polynomial. //--------------------------- define the method --------------------------- for ( i=1; i <= size(#); i++ ) { if ( typeof(#[i]) == "string" ) { //--------------------------- chosen methods ----------------------- if ( (#[i]=="isprim") or (#[i]=="isPrim") ) {decomp = 0;} if ( (#[i]=="nodeco") or (#[i]=="noDeco") ) {decomp = 1;} if (#[i]=="prim") {decomp = 3;} if (#[i]=="equidim") {decomp = 2;} if ( (#[i]=="nofac") or (#[i]=="noFac") ) {noFac=1;} if ( ((#[i]=="useRing") or (#[i]=="usering")) and (ordstr(basering) != "dp("+string(nvars(basering))+"),C")) {useRing = 1;} if ( (#[i]=="withDelta") or (#[i]=="wd") or (#[i]=="withdelta")) { if((decomp == 0) or (decomp == 3)) { withDelta = 1; } else { decomp = 3; withDelta = 1; //Note: the delta invariants cannot be computed with an equidimensional //decomposition, hence we compute first the minimal primes } } if (#[i]=="var1") {denomOption = 1;} if (#[i]=="var2") {denomOption = 2;} } if(typeof(#[i]) == "list"){ if(size(#[i]) == 2){ if (#[i][1]=="inputJ"){ if(typeof(#[i][2]) == "ideal"){ inputJ = #[i][2]; } } } if (#[i][1]=="inputC"){ if(size(#[i]) == 2){ if(typeof(#[i][2]) == "ideal"){ inputC = #[i][2]; } } } } } kill #; //------------------------ change ring if required ------------------------ // If the ordering is not global, we change to dp ordering for computing the // min ass primes. // If the ordering is global, but not dp, and useRing = 0, we also change to // dp ordering. int isGlobal = attrib(basering,"global");// Checks if the original ring has // global ordering. def origR = basering; // origR is the original ring // R is the ring where computations will be done if((useRing == 1) and (isGlobal == 1)) { def globR = basering; } else { // We change to dp ordering. list rl = ringlist(origR); list origOrd = rl[3]; list newOrd = list("dp", intvec(1:nvars(origR))), list("C", 0); rl[3] = newOrd; def globR = ring(rl); setring globR; ideal id = fetch(origR, id); } //------------------------ trivial checkings ------------------------ id = groebner(id); if((size(id) == 0) or (id[1] == 1)) { // The original ring R/I was normal. Nothing to do. // We define anyway a new ring, equal to R, to be able to return it. setring origR; list lR = ringlist(origR); def ROut = ring(lR); setring ROut; ideal norid = fetch(origR, id); ideal normap = maxideal(1); export norid; export normap; setring origR; if(withDelta) { result = list(list(ROut), list(ideal(1)), list(intvec(0), 0)); } else { result = list(list(ROut), list(ideal(1))); } sp = 1; // number of rings in the output option(set, opt); normalOutputText(dbg, withDelta, sp); return(result); } //------------------------ preliminary decomposition----------------------- list prim; if(decomp == 2) { dbprint(dbg, "// Computing the equidimensional decomposition..."); prim = equidim(id); } if((decomp == 0) or (decomp == 1)) { prim = id; } if(decomp == 3) { dbprint(dbg, "// Computing the minimal associated primes..."); if( noFac ) { prim = minAssGTZ(id,1); } else { prim = minAssGTZ(id); } } sp = size(prim); if(dbg>=1) { prim; ""; "// number of components is", sp; ""; } //----------------- back to the original ring if required ------------------ // if ring was not global and useRing is on, we go back to the original ring if((useRing == 1) and (isGlobal != 1)) { setring origR; def R = basering; list prim = fetch(globR, prim); } else { def R = basering; ideal inputJ = fetch(origR, inputJ); ideal inputC = fetch(origR, inputC); if(useRing == 0) { ideal U; poly c; } } // ---------------- normalization of the components------------------------- // calls normalM to compute the normalization of each component. list norComp; // The normalization of each component. int delt; int deltI = 0; int totalComps = 0; setring origR; def newROrigOrd; list newRListO; setring R; def newR; list newRList; for(i=1; i<=size(prim); i++) { if(dbg>=2){pause();} if(dbg>=1) { "// start computation of component",i; " --------------------------------"; } if(groebner(prim[i])[1] != 1) { if(dbg>=2) { "We compute the normalization in the ring"; basering; } printlevel = printlevel + 1; norComp = normalM(prim[i], decomp, withDelta, denomOption, inputJ, inputC); printlevel = printlevel - 1; for(j = 1; j <= size(norComp); j++) { newR = norComp[j][3]; if(!defined(savebasering)) { def savebasering;} savebasering=basering; setring newR; // must be in a compatible ring to newR // as ringlist may produce ring-dep. stuff if(!defined(newRList)) { list newRList;} newRList = ringlist(newR); setring savebasering; U = norComp[j][1]; c = norComp[j][2]; if(withDelta) { delt = norComp[j][4]; if((delt >= 0) and (deltI >= 0)) { deltI = deltI + delt; } else { deltI = -1; } } // -- incorporate result for this component to the list of results --- if(useRing == 0) { // We go back to the original ring. setring origR; U = fetch(R, U); c = fetch(R, c); newRListO = imap(newR, newRList); // We change the ordering in the new ring. if(nvars(newR) > nvars(origR)) { newRListO[3]=insert(origOrd, newRListO[3][1]); } else { newRListO[3] = origOrd; } newROrigOrd = ring(newRListO); setring newROrigOrd; ideal norid = imap(newR, norid); ideal normap = imap(newR, normap); export norid; export normap; setring origR; totalComps++; result[totalComps] = list(U, c, newROrigOrd); if(withDelta) { result[totalComps] = insert(result[totalComps], delt, 3); } setring R; } else { setring R; totalComps++; result[totalComps] = norComp[j]; } } } } // -------------------------- delta computation ---------------------------- if(withDelta == 1) { // Intersection multiplicities of list prim, sp=size(prim). if ( dbg >= 1 ) { "// Sum of delta for all components: ", deltI; } if(size(prim) > 1) { dbprint(dbg, "// Computing the sum of the intersection multiplicities of the components..."); int mul = iMult(prim); if ( mul < 0 ) { deltI = -1; } else { deltI = deltI + mul; } if ( dbg >= 1 ) { "// Intersection multiplicity is : ", mul; } } } // -------------------------- prepare output ------------------------------ setring origR; list RL; // List of rings list MG; // Module generators intvec DV; // Vector of delta's of each component for(i = 1; i <= size(result); i++) { RL[i] = result[i][3]; MG[i] = lineUpLast(result[i][1], result[i][2]); if(withDelta) { DV[i] = result[i][4]; } } if(withDelta) { resultNew = list(RL, MG, list(DV, deltI)); } else { resultNew = list(RL, MG); } sp = size(RL); //RL = list of rings option(set, opt); normalOutputText(dbg, withDelta, sp); return(resultNew); } example { "EXAMPLE:"; printlevel = printlevel+1; echo = 2; ring s = 0,(x,y),dp; ideal i = (x2-y3)*(x2+y2)*x; list nor = normal(i, "withDelta", "prim"); nor; // 2 branches have delta = 1, and 1 branch has delta = 0 // the total delta invariant is 13 def R2 = nor[1][2]; setring R2; norid; normap; echo = 0; printlevel = printlevel-1; pause(" hit return to continue"); echo=2; ring r = 2,(x,y,z),dp; ideal i = z3-xy4; list nor = normal(i, "withDelta", "prim"); nor; // the delta invariant is infinite // xy2z/z2 and xy3/z2 generate the integral closure of r/i as r/i-module // in its quotient field Quot(r/i) // the normalization as affine algebra over the ground field: def R = nor[1][1]; setring R; norid; normap; } /////////////////////////////////////////////////////////////////////////////// // Prints the output text in proc normal. // static proc normalOutputText(int dbg, int withDelta, int sp) // int dbg: printlevel // int withDelta: output contains information about the delta invariant // int sp: number of output rings. { if ( dbg >= 0 ) { ""; if(!withDelta) { "// 'normal' created a list, say nor, of two elements."; } else { "// 'normal' created a list, say nor, of three elements."; } "// To see the list type"; " nor;"; ""; "// * nor[1] is a list of", sp, "ring(s)."; "// To access the i-th ring nor[1][i], give it a name, say Ri, and type"; " def R1 = nor[1][1]; setring R1; norid; normap;"; "// For the other rings type first (if R is the name of your base ring)"; " setring R;"; "// and then continue as for R1."; "// Ri/norid is the affine algebra of the normalization of R/P_i where"; "// P_i is the i-th component of a decomposition of the input ideal id"; "// and normap the normalization map from R to Ri/norid."; ""; "// * nor[2] is a list of", sp, "ideal(s). Let ci be the last generator"; "// of the ideal nor[2][i]. Then the integral closure of R/P_i is"; "// generated as R-submodule of the total ring of fractions by"; "// 1/ci * nor[2][i]."; if(withDelta) { ""; "// * nor[3] is a list of an intvec of size", sp, "the delta invariants "; "// of the components, and an integer, the total delta invariant "; "// of R/id (-1 means infinite, and 0 that R/P_i resp. R/id is normal)."; } } } /////////////////////////////////////////////////////////////////////////////// proc HomJJ (list Li) "USAGE: HomJJ (Li); Li = list: ideal SBid, ideal id, ideal J, poly p ASSUME: R = P/id, P = basering, a polynomial ring, id an ideal of P, @* SBid = standard basis of id, @* J = ideal of P containing the polynomial p, @* p = nonzero divisor of R COMPUTE: Endomorphism ring End_R(J)=Hom_R(J,J) with its ring structure as affine ring, together with the map R --> Hom_R(J,J) of affine rings, where R is the quotient ring of P modulo the standard basis SBid. RETURN: a list l of three objects @format l[1] : a polynomial ring, containing two ideals, 'endid' and 'endphi' such that l[1]/endid = Hom_R(J,J) and endphi describes the canonical map R -> Hom_R(J,J) l[2] : an integer which is 1 if phi is an isomorphism, 0 if not l[3] : an integer, = dim_K(Hom_R(J,J)/R) (the contribution to delta) if the dimension is finite, -1 otherwise @end format NOTE: printlevel >=1: display comments (default: printlevel=0) EXAMPLE: example HomJJ; shows an example " { ASSUME(0, not isQuotientRing(basering) ); //---------- initialisation --------------------------------------------------- int isIso,isPr,isHy,isCo,isRe,isEq,oSAZ,ii,jj,q,y; intvec rw,rw1; list L; y = printlevel-voice+2; // y=printlevel (default: y=0) def P = basering; ideal SBid, id, J = Li[1], Li[2], Li[3]; poly p = Li[4]; int noRed = 0; if(size(Li) > 4) { if(Li[5] == 1) { noRed = 1; } } attrib(SBid,"isSB",1); int homo = homog(Li[2]); //is 1 if id is homogeneous, 0 if not //---- set attributes for special cases where algorithm can be simplified ----- if( homo==1 ) { rw = ringweights(P); } if( typeof(attrib(id,"isPrim"))=="int" ) { if(attrib(id,"isPrim")==1) { isPr=1; } } if( typeof(attrib(id,"onlySingularAtZero"))=="int" ) { if(attrib(id,"onlySingularAtZero")==1){oSAZ=1; } } if( typeof(attrib(id,"isIsolatedSingularity"))=="int" ) { if(attrib(id,"isIsolatedSingularity")==1) { isIso=1; } } if( typeof(attrib(id,"isCohenMacaulay"))=="int" ) { if(attrib(id,"isCohenMacaulay")==1) { isCo=1; } } if( typeof(attrib(id,"isRegInCodim2"))=="int" ) { if(attrib(id,"isRegInCodim2")==1) { isRe=1; } } if( typeof(attrib(id,"isEquidimensional"))=="int" ) { if(attrib(id,"isEquidimensional")==1) { isEq=1; } } //-------------------------- go to quotient ring ------------------------------ qring R = SBid; ideal id = fetch(P,id); ideal J = fetch(P,J); poly p = fetch(P,p); ideal f,rf,f2; module syzf; //---------- computation of p*Hom(J,J) as R-ideal ----------------------------- if ( y>=1 ) { "// compute p*Hom(J,J) = p*J:J"; "// the ideal J:";J; } f = quotient(p*J,J); //### (neu GMG 4.10.08) divide by the greatest common divisor: poly gg = gcd( f[1],p ); for(ii=2; ii <=ncols(f); ii++) { gg=gcd(gg,f[ii]); } for(ii=1; ii<=ncols(f); ii++) { f[ii]=f[ii]/gg; } p = p/gg; if ( y>=1 ) { "// the non-zerodivisor p:"; p; "// the module p*Hom(J,J) = p*J:J :"; f; ""; } f2 = std(p); //---------- Test: Hom(J,J) == R ?, if yes, go home --------------------------- //rf = interred(reduce(f,f2)); //### interred hier weggelassen, unten zugefuegt rf = reduce(f,f2); //represents p*Hom(J,J)/p*R = Hom(J,J)/R if ( size(rf) == 0 ) { if ( homog(f) && find(ordstr(basering),"s")==0 ) { ring newR1 = char(P),(X(1..nvars(P))),(a(rw),dp); } else { ring newR1 = char(P),(X(1..nvars(P))),dp; } ideal endphi = maxideal(1); ideal endid = fetch(P,id); endid = simplify(endid,2); L = substpart(endid,endphi,homo,rw); //## hier substpart def lastRing = L[1]; setring lastRing; attrib(endid,"onlySingularAtZero",oSAZ); attrib(endid,"isCohenMacaulay",isCo); attrib(endid,"isPrim",isPr); attrib(endid,"isIsolatedSingularity",isIso); attrib(endid,"isRegInCodim2",isRe); attrib(endid,"isEqudimensional",isEq); attrib(endid,"isHypersurface",0); attrib(endid,"isCompleteIntersection",0); attrib(endid,"isRadical",0); L=lastRing; L = insert(L,1,1); dbprint(y,"// case R = Hom(J,J)"); if(y>=1) { "// R=Hom(J,J)"; lastRing; "// the new ideal"; endid; " "; "// the old ring"; P; "// the old ideal"; setring P; id; " "; setring lastRing; "// the map to the new ring"; endphi; " "; pause(); ""; } setring P; L[3]=0; return(L); } if(y>=1) { "// R is not equal to Hom(J,J), we have to try again"; pause(); ""; } //---------- Hom(J,J) != R: create new ring and map from old ring ------------- // the ring newR1/SBid+syzf will be isomorphic to Hom(J,J) as R-module // f2=p (i.e. ideal generated by p) //f = mstd(f)[2]; //### geaendert GMG 04.10.08 //ideal ann = quotient(f2,f); //### f durch rf ersetzt rf = mstd(rf)[2]; //rf = NF(f,p), hence = ideal ann = quotient(f2,rf); //p:f = p:rf //------------- compute the contribution to delta ---------- //delt=dim_K(Hom(JJ)/R (or -1 if infinite) int delt=vdim(std(modulo(f,ideal(p)))); f = p,rf; // generates pJ:J mod(p), i.e. p*Hom(J,J)/p*R as R-module q = size(f); syzf = syz(f); if ( homo==1 ) { rw1 = rw,0; for ( ii=2; ii<=q; ii++ ) { rw = rw, deg(f[ii])-deg(f[1]); rw1 = rw1, deg(f[ii])-deg(f[1]); } ring newR1 = char(R),(X(1..nvars(R)),T(1..q)),(a(rw1),dp); } else { ring newR1 = char(R),(X(1..nvars(R)),T(1..q)),dp; } //map psi1 = P,maxideal(1); //### psi1 durch fetch ersetzt //ideal SBid = psi1(SBid); ideal SBid = fetch(P,SBid); attrib(SBid,"isSB",1); qring newR = std(SBid); //map psi = R,ideal(X(1..nvars(R))); //### psi durch fetch ersetzt //ideal id = psi(id); //ideal f = psi(f); //module syzf = psi(syzf); ideal id = fetch(R,id); ideal f = fetch(R,f); module syzf = fetch(R,syzf); ideal pf,Lin,Quad,Q; matrix T,A; list L1; //---------- computation of Hom(J,J) as affine ring --------------------------- // determine kernel of: R[T1,...,Tq] -> J:J >-> R[1/p]=R[t]/(t*p-1), // Ti -> fi/p -> t*fi (p=f1=f[1]), to get ring structure. This is of course // the same as the kernel of R[T1,...,Tq] -> pJ:J >-> R, Ti -> fi. // It is a fact, that the kernel is generated by the linear and the quadratic // relations // f=p,rf, rf=reduce(f,p), generates pJ:J mod(p), // i.e. p*Hom(J,J)/p*R as R-module pf = f[1]*f; T = matrix(ideal(T(1..q)),1,q); Lin = ideal(T*syzf); if(y>=1) { "// the ring structure of Hom(J,J) as R-algebra"; "// the linear relations:"; Lin; } poly ff; for (ii=2; ii<=q; ii++ ) { for ( jj=2; jj<=ii; jj++ ) { ff = NF(f[ii]*f[jj],std(0)); //this makes lift much faster A = lift(pf,ff); //ff lin. comb. of elts of pf mod I Quad = Quad, ideal(T(jj)*T(ii) - T*A); //quadratic relations } } if(y>=1) { "// the quadratic relations"; Quad; pause(); newline; } Q = Lin,Quad; Q = subst(Q,T(1),1); //Q = mstd(Q)[2]; //### sehr aufwendig, daher weggelassen (GMG) //### ev das neue interred //mstd dient nur zum verkleinern, die SB-Eigenschaft geht spaeter verloren //da in neuen Ring abgebildet und mit id vereinigt //---------- reduce number of variables by substitution, if possible ---------- if (homo==1) { ring newRing = char(R),(X(1..nvars(R)),T(2..q)),(a(rw),dp); } else { ring newRing = char(R),(X(1..nvars(R)),T(2..q)),dp; } ideal endid = imap(newR,id),imap(newR,Q); //hier wird Q weiterverwendet, die SB-Eigenschaft wird nicht verwendet. endid = simplify(endid,2); ideal endphi = ideal(X(1..nvars(R))); if(noRed == 0) { L = substpart(endid,endphi,homo,rw); def lastRing=L[1]; setring lastRing; //return(lastRing); } else { list RL = ringlist(newRing); def lastRing = ring(RL); setring lastRing; ideal endid = fetch(newRing, endid); ideal endphi = fetch(newRing, endphi); export(endid); export(endphi); //def lastRing = newRing; //setring R; //return(newR); } // L = substpart(endid,endphi,homo,rw); // def lastRing=L[1]; // setring lastRing; attrib(endid,"onlySingularAtZero",0); map sigma=R,endphi; ideal an=sigma(ann); export(an); //noetig? //ideal te=an,endid; //if(isIso && (size(reduce(te,std(maxideal(1))))==0)) //#### ok??? // { // attrib(endid,"onlySingularAtZero",oSAZ); // } //kill te; attrib(endid,"isCohenMacaulay",isCo); //#### ok??? attrib(endid,"isPrim",isPr); attrib(endid,"isIsolatedSingularity",isIso); attrib(endid,"isRegInCodim2",isRe); attrib(endid,"isEquidimensional",isEq); attrib(endid,"isHypersurface",0); attrib(endid,"isCompleteIntersection",0); attrib(endid,"isRadical",0); if(y>=1) { "// the new ring after reduction of the number of variables"; lastRing; "// the new ideal"; endid; ""; "// the old ring"; P; "// the old ideal"; setring P; id; " "; setring lastRing; "// the map to the new ring"; endphi; " "; pause(); ""; } L = lastRing; L = insert(L,0,1); L[3] = delt; setring(P); return(L); } example {"EXAMPLE:"; echo = 2; ring r = 0,(x,y),wp(2,3); ideal id = y^2-x^3; ideal J = x,y; poly p = x; list Li = std(id),id,J,p; list L = HomJJ(Li); def end = L[1]; // defines ring L[1], containing ideals endid, endphi setring end; // makes end the basering end; endid; // end/endid is isomorphic to End(r/id) as ring map psi = r,endphi;// defines the canonical map r/id -> End(r/id) psi; L[3]; // contribution to delta } /////////////////////////////////////////////////////////////////////////////// //compute intersection multiplicities as needed for delta(I) in //normalizationPrimes and normalP: proc iMult (list prim) "USAGE: iMult(L); L a list of ideals RETURN: int, the intersection multiplicity of the ideals of L; if iMult(L) is infinite, -1 is returned. THEORY: If r=size(L)=2 then iMult(L) = vdim(std(L[1]+L[2])) and in general iMult(L) = sum{ iMult(L[j],Lj) | j=1..r-1 } with Lj the intersection of L[j+1],...,L[r]. If I is the intersection of all ideals in L then we have delta(I) = delta(L[1])+...+delta(L[r]) + iMult(L) where delta(I) = vdim (normalisation(R/I)/(R/I)), R the basering. EXAMPLE: example iMult; shows an example " { ASSUME(0, not isQuotientRing(basering) ); int i,mul,mu; int sp = size(prim); int y = printlevel-voice+2; if ( sp > 1 ) { ideal I(sp-1) = prim[sp]; mu = vdim(std(I(sp-1)+prim[sp-1])); mul = mu; if ( y>=1 ) { "// intersection multiplicity of component",sp,"with",sp-1,":"; mu; } if ( mu >= 0 ) { for (i=sp-2; i>=1 ; i--) { ideal I(i) = intersect(I(i+1),prim[i+1]); mu = vdim(std(I(i)+prim[i])); if ( mu < 0 ) { break; } mul = mul + mu; if ( y>=1 ) { "// intersection multiplicity of components",sp,"...",i+1,"with",i; mu; } } } } return(mul); } example { "EXAMPLE:"; echo = 2; ring s = 23,(x,y),dp; list L = (x-y),(x3+y2); iMult(L); L = (x-y),(x3+y2),(x3-y4); iMult(L); } /////////////////////////////////////////////////////////////////////////////// //check if I has a singularity only at zero, as needed in normalizationPrimes proc locAtZero (ideal I) "USAGE: locAtZero(I); I = ideal RETURN: int, 1 if I has only one point which is located at zero, 0 otherwise ASSUME: I is given as a standard bases in the basering NOTE: only useful in affine rings, in local rings vdim does the check EXAMPLE: example locAtZero; shows an example " { ASSUME(0, not isQuotientRing(basering) ); int ii,jj, caz; //caz: conzentrated at zero int dbp = printlevel-voice+2; int nva = nvars(basering); int vdi = vdim(I); if ( vdi < 0 ) { if (dbp >=1) { "// non-isolated singularitiy";""; } return(caz); } //Now the ideal is 0-dim //First an easy test //If I is homogenous and not constant it is concentrated at 0 if( homog(I)==1 && size(jet(I,0))==0) { caz=1; if (dbp >=1) { "// isolated singularity and homogeneous";""; } return(caz); } //Now the general case with I 0-dim. Choose an appropriate power pot, //and check each variable x whether x^pot is in I. int mi1 = mindeg1(lead(I)); int pot = vdi; if ( (mi1+(mi1==1))^2 < vdi ) { pot = (mi1+(mi1==1))^2; //### alternativ: pot = vdi lassen } while ( 1 ) { caz = 1; for ( ii=1; ii<= nva; ii++ ) { if ( NF(var(ii)^pot,I) != 0 ) { caz = 0; break; } } if ( caz == 1 || pot >= vdi ) { if (dbp >=1) { "// mindeg, exponent, vdim used in 'locAtZero':", mi1,pot,vdi; ""; } return(caz); } else { if ( pot^2 < vdi ) { pot = pot^2; } else { pot = vdi; } } } } example { "EXAMPLE:"; echo = 2; ring r = 0,(x,y,z),dp; poly f = z5+y4+x3+xyz; ideal i = jacob(f),f; i=std(i); locAtZero(i); i= std(i*ideal(x-1,y,z)); locAtZero(i); } /////////////////////////////////////////////////////////////////////////////// //The next procedure normalizationPrimes computes the normalization of an //irreducible or an equidimensional ideal i. //- If i is irreducuble, then the returned list, say nor, has size 2 //with nor[1] the normalization ring and nor[2] the delta invariant. //- If i is equidimensional, than the "splitting tools" can create a //decomposition of i and nor can have more than 1 ring. static proc normalizationPrimes(ideal i,ideal ihp,int delt,intvec delti,list #) "USAGE: normalizationPrimes(i,ihp,delt[,si]); i = equidimensional ideal, ihp = map (partial normalization), delt = partial delta-invariant, si = ideal s.t. V(si) contains singular locus (optional) RETURN: a list of rings, say nor, and an integer, the delta-invariant at the end of the list. each ring nor[j], j = 1..size(nor)-1, contains two ideals with given names norid and normap such that - the direct sum of the rings nor[j]/norid is the normalization of basering/i; - normap gives the normalization map from basering/id to nor[j]/norid (for each j) nor[size(nor)] = dim_K(normalisation(P/i) / (P/i)) is the delta-invariant, where P is the basering. EXAMPLE: example normalizationPrimes; shows an example " { ASSUME(1, not isQuotientRing(basering) ); //Note: this procedure calls itself as long as the test for //normality, i.e if R==Hom(J,J), is negative. int printlev = printlevel; //store printlevel in order to reset it later int y = printlevel-voice+2; // y=printlevel (default: y=0) if(y>=1) { ""; "// START a normalization loop with the ideal"; i; ""; "// in the ring:"; basering; ""; pause(); ""; } def BAS=basering; list result,keepresult1,keepresult2,JM,gnirlist; ideal J,SB,MB; int depth,lauf,prdim,osaz; int ti=timer; gnirlist = ringlist(BAS); //----------- the trivial case of a zero ideal as input, RETURN ------------ if(size(i)==0) { if(y>=1) { "// the ideal was the zero-ideal"; } def newR7 = ring(gnirlist); setring newR7; ideal norid=ideal(0); ideal normap=fetch(BAS,ihp); export norid; export normap; result=newR7; result[size(result)+1]=list(delt,delti); setring BAS; return(result); } //--------------- General NOTATION, compute SB of input ----------------- // SM is a list, the result of mstd(i) // SM[1] = SB of input ideal i, // SM[2] = (minimal) generators for i. // We work with SM and will copy the attributes from i to SM[2] // JM will be a list, either JM[1]=maxideal(1),JM[2]=maxideal(1) // in case i has onlySingularAtZero, or JM = mstd(si) where si = #[1], // or JM = mstd(J) where J is the ideal of the singular locus // JM[2] must be (made) radical if(y>=1) { "// SB-computation of the ideal"; } list SM = mstd(i); //Now the work starts int dimSM = dim(SM[1]); //dimension of variety to normalize if(y>=1) { "// the dimension is:"; dimSM; } //----------------- the general case, set attributes ---------------- //Note: onlySingularAtZero is NOT preserved under the ring extension //basering --> Hom(J,J) (in contrast to isIsolatedSingularity), //therefore we reset it: attrib(i,"onlySingularAtZero",0); if(attrib(i,"isPrim")==1) { attrib(SM[2],"isPrim",1); } else { attrib(SM[2],"isPrim",0); } if(attrib(i,"isIsolatedSingularity")==1) { attrib(SM[2],"isIsolatedSingularity",1); } else { attrib(SM[2],"isIsolatedSingularity",0); } if(attrib(i,"isCohenMacaulay")==1) { attrib(SM[2],"isCohenMacaulay",1); } else { attrib(SM[2],"isCohenMacaulay",0); } if(attrib(i,"isRegInCodim2")==1) { attrib(SM[2],"isRegInCodim2",1); } else { attrib(SM[2],"isRegInCodim2",0); } if(attrib(i,"isEquidimensional")==1) { attrib(SM[2],"isEquidimensional",1); } else { attrib(SM[2],"isEquidimensional",0); } if(attrib(i,"isCompleteIntersection")==1) { attrib(SM[2],"isCompleteIntersection",1); } else { attrib(SM[2],"isCompleteIntersection",0); } if(attrib(i,"isHypersurface")==1) { attrib(SM[2],"isHypersurface",1); } else { attrib(SM[2],"isHypersurface",0); } if(attrib(i,"onlySingularAtZero")==1) { attrib(SM[2],"onlySingularAtZero",1); } else { attrib(SM[2],"onlySingularAtZero",0); } //------- an easy and cheap test for onlySingularAtZero --------- if( (attrib(SM[2],"isIsolatedSingularity")==1) && (homog(SM[2])==1) ) { attrib(SM[2],"onlySingularAtZero",1); } //-------------------- Trivial cases, in each case RETURN ------------------ // input ideal is the ideal of a partial normalization // ------------ Trivial case: input ideal contains a unit --------------- if( dimSM == -1) { ""; " // A unit ideal was found."; " // Stop with partial result computed so far";""; MB=SM[2]; intvec rw; list LL=substpart(MB,ihp,0,rw); def newR6=LL[1]; setring newR6; ideal norid=endid; ideal normap=endphi; kill endid,endphi; export norid; export normap; result=newR6; result[size(result)+1]=list(delt,delti); setring BAS; return(result); } // --- Trivial case: input ideal is zero-dimensional and homog --- if( (dim(SM[1])==0) && (homog(SM[2])==1) ) { if(y>=1) { "// the ideal was zero-dimensional and homogeneous"; } MB=maxideal(1); intvec rw; list LL=substpart(MB,ihp,0,rw); def newR5=LL[1]; setring newR5; ideal norid=endid; ideal normap=endphi; kill endid,endphi; export norid; export normap; result=newR5; result[size(result)+1]=list(delt,delti); setring BAS; return(result); } // --- Trivial case: input ideal defines a line --- //the one-dimensional, homogeneous case and degree 1 case if( (dim(SM[1])==1) && (maxdeg1(SM[2])==1) && (homog(SM[2])==1) ) { if(y>=1) { "// the ideal defines a line"; } MB=SM[2]; intvec rw; list LL=substpart(MB,ihp,0,rw); def newR4=LL[1]; setring newR4; ideal norid=endid; ideal normap=endphi; kill endid,endphi; export norid; export normap; result=newR4; result[size(result)+1]=list(delt,delti); setring BAS; return(result); } //---------------------- The non-trivial cases start ------------------- //the higher dimensional case //we test first hypersurface, CohenMacaulay and complete intersection if( ((size(SM[2])+dim(SM[1])) == nvars(basering)) ) { //the test for complete intersection attrib(SM[2],"isCohenMacaulay",1); attrib(SM[2],"isCompleteIntersection",1); attrib(SM[2],"isEquidimensional",1); if(y>=1) { "// the ideal is a complete intersection"; } } if( size(SM[2]) == 1 ) { attrib(SM[2],"isHypersurface",1); if(y>=1) { "// the ideal is a hypersurface"; } } //------------------- compute the singular locus ------------------- // Computation if singular locus is critical // Notation: J ideal of singular locus or (if given) containing it // JM = mstd(J) or maxideal(1),maxideal(1) // JM[1] SB of singular locus, JM[2] minbasis, dimJ = dim(JM[1]) // SM[1] SB of the input ideal i, SM[2] minbasis // Computation if singular locus is critical, because it determines the // size of the ring Hom_R(J,J). We only need a test ideal contained in J. //----------------------- onlySingularAtZero ------------------------- if( attrib(SM[2],"onlySingularAtZero") ) { JM = maxideal(1),maxideal(1); attrib(JM[1],"isSB",1); attrib(JM[2],"isRadical",1); if( dim(SM[1]) >=2 ) { attrib(SM[2],"isRegInCodim2",1); } } //-------------------- not onlySingularAtZero ------------------------- if( attrib(SM[2],"onlySingularAtZero") == 0 ) { //--- the case where an ideal #[1] is given: if( size(#)>0 ) { J = #[1],SM[2]; JM = mstd(J); if( typeof(attrib(#[1],"isRadical"))!="int" ) { attrib(JM[2],"isRadical",0); } } //--- the case where an ideal #[1] is not given: if( (size(#)==0) ) { if(y >=1 ) { "// singular locus will be computed"; } J = SM[1],minor(jacob(SM[2]),nvars(basering)-dim(SM[1]),SM[1]); if( y >=1 ) { "// SB of singular locus will be computed"; } JM = mstd(J); } int dimJ = dim(JM[1]); attrib(JM[1],"isSB",1); if( y>=1 ) { "// the dimension of the singular locus is"; dimJ ; ""; } if(dim(JM[1]) <= dim(SM[1])-2) { attrib(SM[2],"isRegInCodim2",1); } //------------------ the smooth case, RETURN ------------------- if( dimJ == -1 ) { if(y>=1) { "// the ideal is smooth"; } MB=SM[2]; intvec rw; list LL=substpart(MB,ihp,0,rw); def newR3=LL[1]; setring newR3; ideal norid=endid; ideal normap=endphi; kill endid,endphi; export norid; export normap; result=newR3; result[size(result)+1]=list(delt,delti); setring BAS; return(result); } //------- extra check for onlySingularAtZero, relatively cheap ---------- //it uses the procedure 'locAtZero' from for testing //if an ideal is concentrated at 0 if(y>=1) { "// extra test for onlySingularAtZero:"; } if ( locAtZero(JM[1]) ) { attrib(SM[2],"onlySingularAtZero",1); JM = maxideal(1),maxideal(1); attrib(JM[1],"isSB",1); attrib(JM[2],"isRadical",1); } else { attrib(SM[2],"onlySingularAtZero",0); } } //displaying the attributes: if(y>=2) { "// the attributes of the ideal are:"; "// isCohenMacaulay:", attrib(SM[2],"isCohenMacaulay"); "// isCompleteIntersection:", attrib(SM[2],"isCompleteIntersection"); "// isHypersurface:", attrib(SM[2],"isHypersurface"); "// isEquidimensional:", attrib(SM[2],"isEquidimensional"); "// isPrim:", attrib(SM[2],"isPrim"); "// isRegInCodim2:", attrib(SM[2],"isRegInCodim2"); "// isIsolatedSingularity:", attrib(SM[2],"isIsolatedSingularity"); "// onlySingularAtZero:", attrib(SM[2],"onlySingularAtZero"); "// isRad:", attrib(SM[2],"isRad");""; } //------------- case: CohenMacaulay in codim 2, RETURN --------------- if( (attrib(SM[2],"isRegInCodim2")==1) && (attrib(SM[2],"isCohenMacaulay")==1) ) { if(y>=1) { "// the ideal was CohenMacaulay and regular in codim 2, hence normal"; } MB=SM[2]; intvec rw; list LL=substpart(MB,ihp,0,rw); def newR6=LL[1]; setring newR6; ideal norid=endid; ideal normap=endphi; kill endid,endphi; export norid; export normap; result=newR6; result[size(result)+1]=list(delt,delti); setring BAS; return(result); } //---------- case: isolated singularity only at 0, RETURN ------------ // In this case things are easier, we can use the maximal ideal as radical // of the singular locus; // JM mstd of ideal of singular locus, SM mstd of input ideal if( attrib(SM[2],"onlySingularAtZero") ) { //------ check variables for being a non zero-divizor ------ // SL = ideal of vars not contained in ideal SM[1]: attrib(SM[2],"isIsolatedSingularity",1); ideal SL = simplify(reduce(maxideal(1),SM[1]),2); ideal Ann = quotient(SM[2],SL[1]); ideal qAnn = simplify(reduce(Ann,SM[1]),2); //NOTE: qAnn=0 if and only if first var (=SL[1]) not in SM is a nzd of R/SM //------------- We found a non-zerodivisor of R/SM ----------------------- // here the enlarging of the ring via Hom_R(J,J) starts if( size(qAnn)==0 ) { if(y>=1) { ""; "// the ideal rad(J):"; maxideal(1); ""; } // ------------- test for normality, compute Hom_R(J,J) ------------- // Note: // HomJJ (ideal SBid, ideal id, ideal J, poly p) with // SBid = SB of id, J = radical ideal of basering P with: // nonNormal(R) is in V(J), J contains the nonzero divisor p // of R = P/id (J = test ideal) // returns a list l of three objects // l[1] : a polynomial ring, containing two ideals, 'endid' and 'endphi' // s.t. l[1]/endid = Hom_R(J,J) and endphi= map R -> Hom_R(J,J) // l[2] : an integer which is 1 if phi is an isomorphism, 0 if not // l[3] : an integer, = dim_K(Hom_R(J,J)/R) if finite, -1 otherwise list RR; RR = SM[1],SM[2],maxideal(1),SL[1]; RR = HomJJ(RR,y); // --------------------- non-normal case ------------------ //RR[2]==0 means that the test for normality is negative if( RR[2]==0 ) { def newR=RR[1]; setring newR; map psi=BAS,endphi; list JM = psi(JM); //### ideal J = JM[2]; if ( delt>=0 && RR[3]>=0 ) { delt = delt+RR[3]; } else { delt = -1; } delti[size(delti)]=delt; // ---------- recursive call of normalizationPrimes ----------- //normalizationPrimes(ideal i,ideal ihp,int delt,intvec delti,list #) //ihp = (partial) normalisation map from basering //#[1] ideal s.t. V(#[1]) contains singular locus of i (test ideal) if ( y>=1 ) { "// case: onlySingularAtZero, non-zerodivisor found"; "// contribution of delta in ringextension R -> Hom_R(J,J):"; delt; } //intvec atr=getAttrib(endid); //"//### case: isolated singularity only at 0, recursive"; //"size endid:", size(endid), size(string(endid)); //"interred:"; //endid = interred(endid); //endid = setAttrib(endid,atr); //"size endid:", size(endid), size(string(endid)); printlevel=printlevel+1; list tluser = normalizationPrimes(endid,psi(ihp),delt,delti); //list tluser = // normalizationPrimes(endid,psi(ihp),delt,delti,J); //#### ??? improvement: give also the old ideal of sing locus??? printlevel = printlev; //reset printlevel setring BAS; return(tluser); } // ------------------ the normal case, RETURN ----------------- // Now RR[2] must be 1, hence the test for normality was positive MB=SM[2]; def newR7 = ring(gnirlist); setring newR7; ideal norid=fetch(BAS,MB); ideal normap=fetch(BAS,ihp); if ( delt>=0 && RR[3]>=0 ) { delt = delt+RR[3]; } else { delt = -1; } delti[size(delti)]=delt; intvec atr = getAttrib(norid); //"//### case: isolated singularity only at 0, final"; //"size norid:", size(norid), size(string(norid)); //"interred:"; //norid = interred(norid); //norid = setAttrib(norid,atr); //"size norid:", size(norid), size(string(norid)); export norid; export normap; result=newR7; result[size(result)+1]=list(delt,delti); setring BAS; return(result); } //------ zerodivisor of R/SM was found, gives a splitting ------------ //Now the case where qAnn!=0, i.e. SL[1] is a zero divisor of R/SM //and we have found a splitting: id and id1 //id = Ann defines components of R/SM in the complement of V(SL[1]) //id1 defines components of R/SM in the complement of V(id) else { ideal id = Ann; attrib(id,"isCohenMacaulay",0); attrib(id,"isPrim",0); attrib(id,"isIsolatedSingularity",1); attrib(id,"isRegInCodim2",0); attrib(id,"isHypersurface",0); attrib(id,"isCompleteIntersection",0); attrib(id,"isEquidimensional",0); attrib(id,"onlySingularAtZero",1); ideal id1 = quotient(SM[2],Ann); attrib(id1,"isCohenMacaulay",0); attrib(id1,"isPrim",0); attrib(id1,"isIsolatedSingularity",1); attrib(id1,"isRegInCodim2",0); attrib(id1,"isHypersurface",0); attrib(id1,"isCompleteIntersection",0); attrib(id1,"isEquidimensional",0); attrib(id1,"onlySingularAtZero",1); // ---------- recursive call of normalizationPrimes ----------- if ( y>=1 ) { "// case: onlySingularAtZero, zerodivisor found, splitting:"; "// total delta before splitting:", delt; "// splitting in two components:"; } printlevel = printlevel+1; //to see comments in normalizationPrimes keepresult1 = normalizationPrimes(id,ihp,0,0); //1st split factor keepresult2 = normalizationPrimes(id1,ihp,0,0); //2nd split factor printlevel = printlev; //reset printlevel int delt1 = keepresult1[size(keepresult1)][1]; int delt2 = keepresult2[size(keepresult2)][1]; intvec delti1 = keepresult1[size(keepresult1)][2]; intvec delti2 = keepresult2[size(keepresult2)][2]; if( delt>=0 && delt1>=0 && delt2>=0 ) { ideal idid1=id,id1; int mul = vdim(std(idid1)); if ( mul>=0 ) { delt = delt+mul+delt1+delt2; } else { delt = -1; } } if ( y>=1 ) { "// delta of first component:", delt1; "// delta of second componenet:", delt2; "// intersection multiplicity of both components:", mul; "// total delta after splitting:", delt; } else { delt = -1; } for(lauf=1;lauf<=size(keepresult2)-1;lauf++) { keepresult1=insert(keepresult1,keepresult2[lauf]); } keepresult1[size(keepresult1)]=list(delt,delti); return(keepresult1); } } // Case "onlySingularAtZero" has finished and returned result //-------------- General case, not onlySingularAtZero, RETURN --------------- //test for non-normality, i.e. if Hom(I,I)<>R //we can use Hom(I,I) to continue //------ check variables for being a non zero-divizor ------ // SL = ideal of vars not contained in ideal SM[1]: ideal SL = simplify(reduce(JM[2],SM[1]),2); ideal Ann = quotient(SM[2],SL[1]); ideal qAnn = simplify(reduce(Ann,SM[1]),2); //NOTE: qAnn=0 <==> first var (=SL[1]) not contained in SM is a nzd of R/SM //------------- We found a non-zerodivisor of R/SM ----------------------- //SM = mstd of ideal of variety, JM = mstd of ideal of singular locus if( size(qAnn)==0 ) { list RR; list RS; // ----------------- Computation of the radical ----------------- if(y>=1) { "// radical computation of singular locus"; } J = radical(JM[2]); //the radical of singular locus JM = mstd(J); if(y>=1) { "// radical is equal to:";""; JM[2]; ""; } // ------------ choose non-zerodivisor of smaller degree ---------- //### evtl. fuer SL[1] anderen Nichtnullteiler aus J waehlen ? if( deg(SL[1]) > deg(J[1]) ) { Ann=quotient(SM[2],J[1]); qAnn=simplify(reduce(Ann,SM[1]),2); if(size(qAnn)==0) { SL[1]=J[1]; } } // --------------- computation of Hom(rad(J),rad(J)) -------------- RR=SM[1],SM[2],JM[2],SL[1]; if(y>=1) { "// compute Hom(rad(J),rad(J))"; } RS=HomJJ(RR,y); //most important subprocedure // ------------------ the normal case, RETURN ----------------- // RS[2]==1 means that the test for normality was positive if(RS[2]==1) { def lastR=RS[1]; setring lastR; map psi1=BAS,endphi; ideal norid=endid; ideal normap=psi1(ihp); kill endid,endphi; intvec atr=getAttrib(norid); //"//### general case: not isolated singularity only at 0, final"; //"size norid:", size(norid), size(string(norid)); //"interred:"; //norid = interred(norid); //norid = setAttrib(norid,atr); //"size norid:", size(norid), size(string(norid)); export norid; export normap; result=lastR; if ( y>=1 ) { "// case: not onlySingularAtZero, last ring Hom_R(J,J) computed"; "// delta before last ring:", delt; } if ( delt>=0 && RS[3]>=0 ) { delt = delt+RS[3]; } else { delt = -1; } // delti = delti,delt; delti[size(delti)]=delt; if ( y>=1 ) { "// delta of last ring:", delt; } result[size(result)+1]=list(delt,delti); setring BAS; return(result); } // ----- the non-normal case, recursive call of normalizationPrimes ------- // RS=HomJJ(RR,y) was computed above, RS[1] contains endid and endphi // RS[1] = new ring Hom_R(J,J), RS[2]= 0 or 1, RS[2]=contribution to delta // now RS[2]must be 0, i.e. the test for normality was negative int n = nvars(basering); ideal MJ = JM[2]; def newR=RS[1]; setring newR; map psi=BAS,endphi; if ( y>=1 ) { "// case: not onlySingularAtZero, compute new ring = Hom_R(J,J)"; "// delta of old ring:", delt; } if ( delt>=0 && RS[3]>=0 ) { delt = delt+RS[3]; } else { delt = -1; } if ( y>=1 ) { "// delta of new ring:", delt; } delti[size(delti)]=delt; intvec atr=getAttrib(endid); //"//### general case: not isolated singularity only at 0, recursive"; //"size endid:", size(endid), size(string(endid)); //"interred:"; //endid = interred(endid); //endid = setAttrib(endid,atr); //"size endid:", size(endid), size(string(endid)); printlevel = printlevel+1; list tluser= normalizationPrimes(endid,psi(ihp),delt,delti,psi(MJ)); printlevel = printlev; //reset printlevel setring BAS; return(tluser); } //---- A whole singular component was found, RETURN ----- if( Ann == 1) { "// Input appeared not to be a radical ideal!"; "// A (everywhere singular) component with ideal"; "// equal to its Jacobian ideal was found"; "// Procedure will stop with partial result computed so far";""; MB=SM[2]; intvec rw; list LL=substpart(MB,ihp,0,rw); def newR6=LL[1]; setring newR6; ideal norid=endid; ideal normap=endphi; kill endid,endphi; export norid; export normap; result=newR6; result[size(result)+1]=lst(delt,delti); setring BAS; return(result); } //------ zerodivisor of R/SM was found, gives a splitting ------------ //Now the case where qAnn!=0, i.e. SL[1] is a zero divisor of R/SM //and we have found a splitting: new1 and new2 //id = Ann defines components of R/SM in the complement of V(SL[1]) //id1 defines components of R/SM in the complement of V(id) else { if(y>=1) { "// zero-divisor found"; } int equi = attrib(SM[2],"isEquidimensional"); int oSAZ = attrib(SM[2],"onlySingularAtZero"); int isIs = attrib(SM[2],"isIsolatedSingularity"); ideal new1 = Ann; ideal new2 = quotient(SM[2],Ann); //ideal new2=SL[1],SM[2]; def newR1 = ring(gnirlist); setring newR1; ideal vid = fetch(BAS,new1); ideal ihp = fetch(BAS,ihp); attrib(vid,"isCohenMacaulay",0); attrib(vid,"isPrim",0); attrib(vid,"isIsolatedSingularity",isIs); attrib(vid,"isRegInCodim2",0); attrib(vid,"onlySingularAtZero",oSAZ); attrib(vid,"isEquidimensional",equi); attrib(vid,"isHypersurface",0); attrib(vid,"isCompleteIntersection",0); // ---------- recursive call of normalizationPrimes ----------- if ( y>=1 ) { "// total delta before splitting:", delt; "// splitting in two components:"; } printlevel = printlevel+1; keepresult1 = normalizationPrimes(vid,ihp,0,0); //1st split factor list delta1 = keepresult1[size(keepresult1)]; setring BAS; def newR2 = ring(gnirlist); setring newR2; ideal vid = fetch(BAS,new2); ideal ihp = fetch(BAS,ihp); attrib(vid,"isCohenMacaulay",0); attrib(vid,"isPrim",0); attrib(vid,"isIsolatedSingularity",isIs); attrib(vid,"isRegInCodim2",0); attrib(vid,"isEquidimensional",equi); attrib(vid,"isHypersurface",0); attrib(vid,"isCompleteIntersection",0); attrib(vid,"onlySingularAtZero",oSAZ); keepresult2 = normalizationPrimes(vid,ihp,0,0); list delta2 = keepresult2[size(keepresult2)]; //2nd split factor printlevel = printlev; //reset printlevel setring BAS; //compute intersection multiplicity of both components: new1 = new1,new2; int mul=vdim(std(new1)); // ----- normalizationPrimes finished, add up results, RETURN -------- for(lauf=1;lauf<=size(keepresult2)-1;lauf++) { keepresult1 = insert(keepresult1,keepresult2[lauf]); } if ( delt >=0 && delta1[1] >=0 && delta2[1] >=0 && mul >=0 ) { delt = delt+mul+delta1[1]+delta2[1]; } else { delt = -1; } delti = -2; if ( y>=1 ) { "// zero divisor produced a splitting into two components"; "// delta of first component:", delta1; "// delta of second componenet:", delta2; "// intersection multiplicity of both components:", mul; "// total delta after splitting:", delt; } keepresult1[size(keepresult1)] = list(delt,delti); return(keepresult1); } } example { "EXAMPLE:";echo = 2; // Huneke ring qr=31991,(a,b,c,d,e),dp; ideal i= 5abcde-a5-b5-c5-d5-e5, ab3c+bc3d+a3be+cd3e+ade3, a2bc2+b2cd2+a2d2e+ab2e2+c2de2, abc5-b4c2d-2a2b2cde+ac3d2e-a4de2+bcd2e3+abe5, ab2c4-b5cd-a2b3de+2abc2d2e+ad4e2-a2bce3-cde5, a3b2cd-bc2d4+ab2c3e-b5de-d6e+3abcd2e2-a2be4-de6, a4b2c-abc2d3-ab5e-b3c2de-ad5e+2a2bcde2+cd2e4, b6c+bc6+a2b4e-3ab2c2de+c4d2e-a3cde2-abd3e2+bce5; list pr=normalizationPrimes(i); def r1 = pr[1]; setring r1; norid; normap; } /////////////////////////////////////////////////////////////////////////////// static proc substpart(ideal endid, ideal endphi, int homo, intvec rw) "//Repeated application of elimpart to endid, until no variables can be //directy substituded. homo=1 if input is homogeneous, rw contains //original weights, endphi (partial) normalization map"; //NOTE concerning iteration of maps: Let phi: x->f(y,z), y->g(x,z) then //phi: x+y+z->f(y,z)+g(x,z)+z, phi(phi):x+y+z->f(g(x,z),z)+g(f(y,z),z)+z //and so on: none of the x or y will be eliminated //Now subst: first x and then y: x+y+z->f(g(x,z),z)+g(x,z)+z eliminates y //further subst replaces x by y, makes no sense (objects more compicated). //Subst first y and then x eliminates x //In our situation we have triangular form: x->f(y,z), y->g(z). //phi: x+y+z->f(y,z)+g(z)+z, phi(phi):x+y+z->f(g(z),z)+g(z)+z eliminates x,y //subst x,y: x+y+z->f(g(z),z)+g(z)+z, eliminates x,y //subst y,x: x+y+z->f(y,z)+g(z)+z eliminates only x //HENCE: substitute vars depending on most other vars first //However, if the sytem xi-fi is reduced then xi does not appear in any of the //fj and hence the order does'nt matter when substitutinp xi by fi { ASSUME(1, not isQuotientRing(basering) ); def newRing = basering; int ii,jj; map phi = newRing,maxideal(1); //identity map list Le = elimpart(endid); //this proc and the next loop try to substitute as many variables as //possible indices of substituted variables int q = size(Le[2]); //q vars, stored in Le[2], have been substitutet intvec rw1 = 0; //will become indices of substituted variables rw1[nvars(basering)] = 0; rw1 = rw1+1; //rw1=1,..,1 (as many 1 as nvars(basering)) while( size(Le[2]) != 0 ) { endid = Le[1]; if ( defined(ps) ) { kill ps; } map ps = newRing,Le[5]; phi = ps(phi); for(ii=1;ii<=size(Le[2]);ii++) { phi=phi(phi); } //eingefuegt wegen x2-y2z2+z3 for( ii=1; ii<=size(rw1); ii++ ) { if( Le[4][ii]==0 ) //ii = index of var which was substituted { rw1[ii]=0; //substituted vars have entry 0 in rw1 } } Le=elimpart(endid); //repeated application of elimpart q = q + size(Le[2]); } endphi = phi(endphi); //---------- return ----------------------------------------------------------- // first the trivial case, where all variable have been eliminated if( nvars(newRing) == q ) { ring lastRing = char(basering),T(1),dp; ideal endid = T(1); ideal endphi = T(1); for(ii=2; ii<=q; ii++ ) { endphi[ii] = 0; } export(endid,endphi); list L = lastRing; setring newRing; return(L); } // in the homogeneous case put weights for the remaining vars correctly, i.e. // delete from rw those weights for which the corresponding entry of rw1 is 0 if (homo==1 && nvars(newRing)-q >1 && size(endid) >0 ) { jj=1; for( ii=2; ii1){ERROR("no hypersurface");} ideal J=std(slocus(I)); if(dim(J)<=0){return(0);} poly h; d=1; while((d)&&(i10){ERROR("delta not found, please inform the authors")}; h=randomLast(100)[n]; d=dim(std(J+ideal(h))); } I=I,h-1; if(char(R)<=19) { nor=normalP(I); } else { nor=normal(I); } return(nor[2][2]); } proc genus(ideal I,list #) "USAGE: genus(I) or genus(I,