//////////////////////////////////////////////////////////////////////// version="$Id$"; category="Noncommutative"; info=" LIBRARY: dmodvar.lib Algebraic D-modules for varieties AUTHORS: Daniel Andres, daniel.andres@math.rwth-aachen.de @* Viktor Levandovskyy, levandov@math.rwth-aachen.de @* Jorge Martin-Morales, jorge@unizar.es Support: DFG Graduiertenkolleg 1632 'Experimentelle und konstruktive Algebra' OVERVIEW: Let K be a field of characteristic 0. Given a polynomial ring R = K[x_1,...,x_n] and polynomials f_1,...,f_r in R, define F = f_1*...*f_r and F^s = f_1^s_1*...*f_r^s_r for symbolic discrete (that is shiftable) variables s_1,..., s_r. The module R[1/F]*F^s has the structure of a D-module, where D = D(R) tensored with S over K, where @* - D(R) is the n-th Weyl algebra K @* - S is the universal enveloping algebra of gl_r, generated by s_i = s_{ii}. @* One is interested in the following data: @* - the left ideal Ann F^s in D, usually denoted by LD in the output @* - global Bernstein polynomial in one variable s = s_1+...+s_r, denoted by bs, @* - its minimal integer root s0, the list of all roots of bs, which are known to be negative rational numbers, with their multiplicities, which is denoted by BS @* - an r-tuple of operators in D, denoted by PS, such that the functional equality sum(k=1 to k=r) P_k*f_k*F^s = bs*F^s holds in R[1/F]*F^s. References: (BMS06) Budur, Mustata, Saito: Bernstein-Sato polynomials of arbitrary varieties (2006). @* (ALM09) Andres, Levandovskyy, Martin-Morales: Principal Intersection and Bernstein-Sato Polynomial of an Affine Variety (2009). PROCEDURES: bfctVarIn(F[,L]); computes the roots of the Bernstein-Sato polynomial b(s) of the variety V(F) using initial ideal approach bfctVarAnn(F[,L]); computes the roots of the Bernstein-Sato polynomial b(s) of the variety V(F) using Sannfs approach SannfsVar(F[,O,e]); computes the annihilator of F^s in the ring D makeMalgrange(F[,ORD]); creates the Malgrange ideal, associated with F = F[1],..,F[P] SEE ALSO: bfun_lib, dmod_lib, dmodapp_lib, gmssing_lib KEYWORDS: D-module; D-module structure; Bernstein-Sato polynomial for variety; global Bernstein-Sato polynomial for variety; Weyl algebra; parametric annihilator for variety; Budur-Mustata-Saito approach; initial ideal approach "; /* // Static procs: // coDim(I); compute the codimension of the leading ideal of I // dmodvarAssumeViolation() // ORDstr2list (ORD, NN) // smallGenCoDim(I,k) */ /* CHANGELOG 11.10.10 by DA: - reformated help strings - simplified code - add and use of safeVarName - renamed makeIF to makeMalgrange */ LIB "bfun.lib"; // for pIntersect LIB "dmodapp.lib"; // for isCommutative etc. /////////////////////////////////////////////////////////////////////////////// // testing for consistency of the library: proc testdmodvarlib () { example makeMalgrange; example bfctVarIn; example bfctVarAnn; example SannfsVar; } // example coDim; /////////////////////////////////////////////////////////////////////////////// static proc dmodvarAssumeViolation() { // char K = 0, no qring if ( (size(ideal(basering)) >0) || (char(basering) >0) ) { ERROR("Basering is inappropriate: characteristic>0 or qring present"); } return(); } static proc safeVarName (string s, string cv) // assumes 's' to be a valid variable name // returns valid var name string @@..@s { string S; if (cv == "v") { S = "," + "," + varstr(basering) + ","; } if (cv == "c") { S = "," + "," + charstr(basering) + ","; } if (cv == "cv") { S = "," + charstr(basering) + "," + varstr(basering) + ","; } s = "," + s + ","; while (find(S,s) <> 0) { s[1] = "@"; s = "," + s; } s = s[2..size(s)-1]; return(s) } // da: in smallGenCoDim(), rewritten using mstd business static proc coDim (ideal I) "USAGE: coDim (I); I an ideal RETURN: int PURPOSE: computes the codimension of the ideal generated by the leading monomials of the given generators of the ideal. This is also the codimension of the ideal if it is represented by a standard basis. NOTE: The codimension of an ideal I means the number of variables minus the Krull dimension of the basering modulo I. EXAMPLE: example coDim; shows examples " { int n = nvars(basering); int d = dim(I); // to insert: check whether I is in GB return(n-d); } example { "EXAMPLE:"; echo = 2; ring R = 0,(x,y,z),Dp; ideal I = x^2+y^3, z; coDim(std(I)); } static proc ORDstr2list (string ORD, int NN) { /* convert an ordering defined in NN variables in the */ /* string form into the same ordering in the list form */ string st; st = "ring @Z = 0,z(1.." + string(NN) + "),"; st = st + ORD + ";"; execute(st); kill st; list L = ringlist(@Z)[3]; kill @Z; return(L); } proc SannfsVar (ideal F, list #) "USAGE: SannfsVar(F [,ORD,eng]); F an ideal, ORD an optional string, eng an optional int RETURN: ring (Weyl algebra tensored with U(gl_P)), containing an ideal LD PURPOSE: compute the D-module structure of D*f^s where f = F[1]*...*F[P] and D is the Weyl algebra D tensored with K=U(gl_P), according to the generalized algorithm by Briancon and Maisonobe for affine varieties ASSUME: The basering is commutative and over a field of characteristic 0. NOTE: Activate the output ring D with the @code{setring} command. In the ring D, the ideal LD is the needed D-module structure. @* The value of ORD must be an elimination ordering in D for Dt written in the string form, otherwise the result may have no meaning. By default ORD = '(a(1..(P)..1),a(1..(P+P^2)..1),dp)'. @* If eng<>0, @code{std} is used for Groebner basis computations, otherwise, and by default @code{slimgb} is used. DISPLAY: If printlevel=1, progress debug messages will be printed, @* if printlevel>=2, all the debug messages will be printed. EXAMPLE: example SannfsVar; shows examples " { dmodvarAssumeViolation(); if (!isCommutative()) { ERROR("Basering must be commutative"); } def save = basering; int N = nvars(basering); int P = ncols(F); //ncols better than size, since F[i] could be zero // P is needed for default ORD int i,j,k,l; // st = "(a(1..(P)..1),a(1..(P+P^2)..1),dp)"; string st = "(a(" + string(1:P); st = st + "),a(" + string(1:(P+P^2)); st = st + "),dp)"; // default values string ORD = st; int eng = 0; if ( size(#)>0 ) { if ( typeof(#[1]) == "string" ) { ORD = string(#[1]); // second arg if (size(#)>1) { // exists 2nd arg if ( typeof(#[2]) == "int" ) { // the case: given ORD, given engine eng = int(#[2]); } } } else { if ( typeof(#[1]) == "int" ) { // the case: default ORD, engine given eng = int(#[1]); // ORD = "(a(1..(P)..1),a(1..(P+P^2)..1),dp)"; //is already set } else { // incorr. 1st arg ORD = string(st); } } } // size(#)=0, i.e. there is no elimination ordering and no engine given // eng = 0; ORD = "(a(1..(P)..1),a(1..(P^2)..1),dp)"; //are already set int ppl = printlevel-voice+2; // returns a list with a ring and an ideal LD in it // save, N, P and the indices are already defined int Nnew = 2*N+P+P^2; list RL = ringlist(basering); list L; L[1] = RL[1]; //char L[4] = RL[4]; //char, minpoly // check whether vars have admissible names list Name = RL[2]; list RName; // (i,j) <--> (i-1)*p+j for (i=1; i<=P; i++) { RName[i] = safeVarName("Dt("+string(i)+")","cv"); for (j=1; j<=P; j++) { RName[P+(i-1)*P+j] = safeVarName("s("+string(i)+")("+string(j)+")","cv"); } } // now, create the names for new vars list DName; for(i=1; i<=N; i++) { DName[i] = safeVarName("D"+Name[i],"cv"); //concat } list NName = RName + Name + DName; L[2] = NName; // Name, Dname will be used further kill NName; //block ord (a(1..(P)..1),a(1..(P+P^2)..1),dp); //export Nnew; L[3] = ORDstr2list(ORD,Nnew); // we are done with the list def @R@ = ring(L); setring @R@; matrix @D[Nnew][Nnew]; // kronecker(i,j) equals (i==j) // (i,j) <--> (i-1)*p+j for (i=1; i<=P; i++) { for (j=1; j<=P; j++) { for (k=1; k<=P; k++) { //[sij,Dtk] = djk*Dti // @D[k,P+(i-1)*P+j] = (j==k)*Dt(i); @D[k,P+(i-1)*P+j] = (j==k)*var(i); for (l=1; l<=P; l++) { if ( (i-k)*P < l-j ) { //[sij,skl] = djk*sil - dil*skj // @D[P+(i-1)*P+j,P+(k-1)*P+l] = -(j==k)*s(i)(l) + (i==l)*s(k)(j); @D[P+(i-1)*P+j,P+(k-1)*P+l] = -(j==k)*var(i*P+l) + (i==l)*var(k*P+j); } } } } } for (i=1; i<=N; i++) { //[Dx,x]=1 @D[P+P^2+i,P+P^2+N+i] = 1; } def @R = nc_algebra(1,@D); setring @R; //@R@ will be used further dbprint(ppl,"// -1-1- the ring @R(_Dt,_s,_x,_Dx) is ready"); dbprint(ppl-1, @R); // create the ideal I // (i,j) <--> (i-1)*p+j ideal F = imap(save,F); ideal I; for (i=1; i<=P; i++) { for (j=1; j<=P; j++) { // I[(i-1)*P+j] = Dt(i)*F[j] + s(i)(j); I[(i-1)*P+j] = var(i)*F[j] + var(i*P+j); } } poly p,q; for (i=1; i<=N; i++) { p=0; for (j=1; j<=P; j++) { // q = Dt(j); q = var(j); q = q*diff(F[j],var(P+P^2+i)); p = p + q; } I = I, p + var(P+P^2+N+i); } // -------- the ideal I is ready ---------- dbprint(ppl,"// -1-2- starting the elimination of Dt(i) in @R"); dbprint(ppl-1, I); ideal J = engine(I,eng); ideal K = nselect(J,1..P); kill I,J; dbprint(ppl,"// -1-3- all Dt(i) are eliminated"); dbprint(ppl-1, K); //K is without Dt(i) // ----------- the ring @R2(_s,_x,_Dx) ------------ //come back to the ring save, recover L and remove all Dt(i) //L[1],L[4] do not change setring save; list Lord, tmp; // variables tmp = L[2]; Lord = tmp[P+1..Nnew]; L[2] = Lord; // ordering // st = "(a(1..(P^2)..1),dp)"; st = "(a(" + string(1:P^2); st = st + "),dp)"; tmp = ORDstr2list(st,Nnew-P); L[3] = tmp; def @R2@ = ring(L); kill L; // we are done with the list setring @R2@; matrix tmpM,LordM; // non-commutative relations intvec iv = P+1..Nnew; tmpM = imap(@R@,@D); kill @R@; LordM = submat(tmpM,iv,iv); matrix @D2 = LordM; def @R2 = nc_algebra(1,@D2); setring @R2; kill @R2@; dbprint(ppl,"// -2-1- the ring @R(_s,_x,_Dx) is ready"); dbprint(ppl-1, @R2); ideal K = imap(@R,K); kill @R; dbprint(ppl,"// -2-2- starting cosmetic Groebner basis computation"); dbprint(ppl-1, K); K = engine(K,eng); dbprint(ppl,"// -2-3- the cosmetic Groebner basis has been computed"); dbprint(ppl-1,K); ideal LD = K; attrib(LD,"isSB",1); export LD; return(@R2); } example { "EXAMPLE:"; echo = 2; ring R = 0,(x,y),Dp; ideal F = x^3, y^5; //ORD = "(a(1,1),a(1,1,1,1,1,1),dp)"; //eng = 0; def A = SannfsVar(F); setring A; A; LD; } proc bfctVarAnn (ideal F, list #) "USAGE: bfctVarAnn(F[,gid,eng]); F an ideal, gid,eng optional ints RETURN: list of an ideal and an intvec PURPOSE: computes the roots of the Bernstein-Sato polynomial and their multiplicities for an affine algebraic variety defined by F = F[1],..,F[r]. ASSUME: The basering is commutative and over a field in char 0. NOTE: In the output list, the ideal contains all the roots and the intvec their multiplicities. @* If gid<>0, the ideal is used as given. Otherwise, and by default, a heuristically better suited generating set is used. @* If eng<>0, @code{std} is used for GB computations, otherwise, and by default, @code{slimgb} is used. @* Computational remark: The time of computation can be very different depending on the chosen generators of F, although the result is always the same. @* Further note that in this proc, the annihilator of f^s in D[s] is computed and then a system of linear equations is solved by linear reductions in order to find the minimal polynomial of S = s(1)(1) + ... + s(P)(P). The resulted is shifted by 1-codim(Var(F)) following (BMS06). DISPLAY: If printlevel=1, progress debug messages will be printed, @* if printlevel=2, all the debug messages will be printed. EXAMPLE: example bfctVarAnn; shows examples " { dmodvarAssumeViolation(); if (!isCommutative()) { ERROR("Basering must be commutative"); } int gid = 0; // default int eng = 0; // default if (size(#)>0) { if (typeof(#[1])=="int" || typeof(#[1])=="number") { gid = int(#[1]); } if (size(#)>1) { if (typeof(#[2])=="int" || typeof(#[2])=="number") { eng = int(#[2]); } } } def save = basering; int ppl = printlevel - voice + 2; printlevel = printlevel+1; list L = smallGenCoDim(F,gid); F = L[1]; int cd = L[2]; kill L; def @R2 = SannfsVar(F,eng); printlevel = printlevel-1; int sF = size(F); // no 0 in F setring @R2; // we are in D[s] and LD is a std of SannfsVar(F) ideal F = imap(save,F); ideal GF = std(F); ideal J = NF(LD,GF); J = J, F; dbprint(ppl,"// -3-1- starting Groebner basis of ann F^s + F "); dbprint(ppl-1,J); ideal K = engine(J,eng); dbprint(ppl,"// -3-2- finished Groebner basis of ann F^s + F "); dbprint(ppl-1,K); poly S; int i; for (i=1; i<=sF; i++) { // S = S + s(i)(i); S = S + var((i-1)*sF+i); } dbprint(ppl,"// -4-1- computing the minimal polynomial of S"); dbprint(ppl-1,"S = "+string(S)); vector M = pIntersect(S,K); dbprint(ppl,"// -4-2- the minimal polynomial has been computed"); ring @R3 = 0,s,dp; vector M = imap(@R2,M); poly p = vec2poly(M); dbprint(ppl-1,p); dbprint(ppl,"// -5-1- codimension of the variety"); dbprint(ppl-1,cd); dbprint(ppl,"// -5-2- shifting BS(s)=minpoly(s-codim+1)"); p = subst(p,var(1),var(1)-cd+1); dbprint(ppl-1,p); dbprint(ppl,"// -5-3- factorization of the minimal polynomial"); list BS = bFactor(p); setring save; list BS = imap(@R3,BS); kill @R2,@R3; return(BS); } example { "EXAMPLE:"; echo = 2; ring R = 0,(x,y,z),Dp; ideal F = x^2+y^3, z; bfctVarAnn(F); } proc makeMalgrange (ideal F, list #) "USAGE: makeMalgrange(F [,ORD]); F an ideal, ORD an optional string RETURN: ring (Weyl algebra) containing an ideal IF PURPOSE: create the ideal by Malgrange associated with F = F[1],...,F[P]. NOTE: Activate the output ring with the @code{setring} command. In this ring, the ideal IF is the ideal by Malgrange corresponding to F. @* The value of ORD must be an arbitrary ordering in K<_t,_x,_Dt,_Dx> written in the string form. By default ORD = 'dp'. DISPLAY: If printlevel=1, progress debug messages will be printed, @* if printlevel>=2, all the debug messages will be printed. EXAMPLE: example makeMalgrange; shows examples " { string ORD = "dp"; if ( size(#)>0 ) { if ( typeof(#[1]) == "string" ) { ORD = string(#[1]); } } int ppl = printlevel-voice+2; def save = basering; int N = nvars(save); int P = ncols(F); int Nnew = 2*P+2*N; int i,j; string st; list RL = ringlist(save); list L,Lord; list tmp; intvec iv; L[1] = RL[1]; L[4] = RL[4]; //check whether vars have admissible names list Name = RL[2]; list TName, DTName; for (i=1; i<=P; i++) { TName[i] = safeVarName("t("+string(i)+")","cv"); DTName[i] = safeVarName("Dt("+string(i)+")","cv"); } //now, create the names for new vars list DName; for (i=1; i<=N; i++) { DName[i] = safeVarName("D"+Name[i],"cv"); //concat } list NName = TName + Name + DTName + DName; L[2] = NName; // Name, Dname will be used further kill NName, TName, Name, DTName, DName; // ORD already set, default ord dp; L[3] = ORDstr2list(ORD,Nnew); // we are done with the list def @R@ = ring(L); setring @R@; def @R = Weyl(); setring @R; kill @R@; //dbprint(ppl,"// -1-1- the ring @R(_t,_x,_Dt,_Dx) is ready"); //dbprint(ppl-1, @R); // create the ideal I ideal F = imap(save,F); ideal I; for (j=1; j<=P; j++) { // I[j] = t(j) - F[j]; I[j] = var(j) - F[j]; } poly p,q; for (i=1; i<=N; i++) { p=0; for (j=1; j<=P; j++) { // q = Dt(j); q = var(P+N+j); q = diff(F[j],var(P+i))*q; p = p + q; } I = I, p + var(2*P+N+i); } // -------- the ideal I is ready ---------- ideal IF = I; export IF; return(@R); } example { "EXAMPLE:"; echo = 2; ring R = 0,(x,y,z),Dp; ideal I = x^2+y^3, z; def W = makeMalgrange(I); setring W; W; IF; } proc bfctVarIn (ideal I, list #) "USAGE: bfctVarIn(I [,a,b,c]); I an ideal, a,b,c optional ints RETURN: list of ideal and intvec PURPOSE: computes the roots of the Bernstein-Sato polynomial and their multiplicities for an affine algebraic variety defined by I. ASSUME: The basering is commutative and over a field of characteristic 0. @* Varnames of the basering do not include t(1),...,t(r) and Dt(1),...,Dt(r), where r is the number of entries of the input ideal. NOTE: In the output list, say L, @* - L[1] of type ideal contains all the rational roots of a b-function, @* - L[2] of type intvec contains the multiplicities of above roots, @* - optional L[3] of type string is the part of b-function without rational roots. @* Note, that a b-function of degree 0 is encoded via L[1][1]=0, L[2]=0 and L[3] is 1 (for nonzero constant) or 0 (for zero b-function). @* If a<>0, the ideal is used as given. Otherwise, and by default, a heuristically better suited generating set is used to reduce computation time. @* If b<>0, @code{std} is used for GB computations in characteristic 0, otherwise, and by default, @code{slimgb} is used. @* If c<>0, a matrix ordering is used for GB computations, otherwise, and by default, a block ordering is used. @* Further note, that in this proc, the initial ideal of the multivariate Malgrange ideal defined by I is computed and then a system of linear equations is solved by linear reductions following the ideas by Noro. The result is shifted by 1-codim(Var(F)) following (BMS06). DISPLAY: If printlevel=1, progress debug messages will be printed, @* if printlevel>=2, all the debug messages will be printed. EXAMPLE: example bfctVarIn; shows examples " { dmodvarAssumeViolation(); if (!isCommutative()) { ERROR("Basering must be commutative"); } int ppl = printlevel - voice + 2; int idealasgiven = 0; // default int whicheng = 0; // default int whichord = 0; // default if (size(#)>0) { if (typeof(#[1])=="int" || typeof(#[1])=="number") { idealasgiven = int(#[1]); } if (size(#)>1) { if (typeof(#[2])=="int" || typeof(#[2])=="number") { whicheng = int(#[2]); } if (size(#)>2) { if (typeof(#[3])=="int" || typeof(#[3])=="number") { whichord = int(#[3]); } } } } def save = basering; int i; int n = nvars(basering); // step 0: get small generating set I = simplify(I,2); list L = smallGenCoDim(I,idealasgiven); I = L[1]; int c = L[2]; kill L; // step 1: setting up the multivariate Malgrange ideal int r = size(I); def D = makeMalgrange(I); setring D; dbprint(ppl-1,"// Computing in " + string(n+r) + "-th Weyl algebra:", D); dbprint(ppl-1,"// The Malgrange ideal: ", IF); // step 2: compute the b-function of the Malgrange ideal w.r.t. approriate weights intvec w = 1:r; w[r+n] = 0; dbprint(ppl,"// Computing the b-function of the Malgrange ideal..."); list L = bfctIdeal(IF,w,whicheng,whichord); dbprint(ppl,"// ... done."); dbprint(ppl-1,"// The b-function: ",L); // step 3: shift the result ring S = 0,s,dp; list L = imap(D,L); kill D; if (size(L)==2) { ideal B = L[1]; ideal BB; int nB = ncols(B); for (i=nB; i>0; i--) { BB[i] = -B[nB-i+1]+c-r-1; } L[1] = BB; } else // should never get here: BS poly has non-rational roots { string str = L[3]; L = delete(L,3); str = "poly @b = (" + str + ")*(" + string(fl2poly(L,"s")) + ");"; execute(str); poly b = subst(@b,s,-s+c-r-1); L = bFactor(b); } setring save; list L = imap(S,L); return(L); } example { "EXAMPLE:"; echo = 2; ring R = 0,(x,y,z),dp; ideal F = x^2+y^3, z; list L = bfctVarIn(F); L; } static proc smallGenCoDim (ideal I, int Iasgiven) { // call from K[x], returns list L // L[1]=I or L[1]=smaller generating set of I // L[2]=codimension(I) int ppl = printlevel - voice + 2; int n = nvars(basering); int c; if (attrib(I,"isSB") == 1) { c = n - dim(I); if (!Iasgiven) { list L = mstd(I); } } else { def save = basering; list RL = ringlist(save); list @ord; @ord[1] = list("dp", intvec(1:n)); @ord[2] = list("C", intvec(0)); RL[3] = @ord; kill @ord; if (size(RL)>4) // commutative G-algebra present { RL = RL[1..4]; } def R = ring(RL); kill RL; setring R; ideal I = imap(save,I); if (!Iasgiven) { list L = mstd(I); c = n - dim(L[1]); setring save; list L = imap(R,L); } else { I = std(I); c = n - dim(I); setring save; } kill R; } if (!Iasgiven) { if (size(L[2]) < size(I)) { I = L[2]; dbprint(ppl,"// Found smaller generating set of the given variety: ", I); } else { dbprint(ppl,"// Have not found smaller generating set of the given variety."); } } dbprint(ppl-1,"// The codim of the given variety is " + string(c) + "."); if (!defined(L)) { list L; } L[1] = I; L[2] = c; return(L); } /* // Some more examples static proc TXcups() { "EXAMPLE:"; echo = 2; //TX tangent space of X=V(x^2+y^3) ring R = 0,(x0,x1,y0,y1),Dp; ideal F = x0^2+y0^3, 2*x0*x1+3*y0^2*y1; printlevel = 0; //ORD = "(a(1,1),a(1,1,1,1,1,1),dp)"; //eng = 0; def A = SannfsVar(F); setring A; LD; } static proc ex47() { ring r7 = 0,(x0,x1,y0,y1),dp; ideal I = x0^2+y0^3, 2*x0*x1+3*y0^2*y1; bfctVarIn(I); // second ex - too big ideal J = x0^4+y0^5, 4*x0^3*x1+5*y0^4*y1; bfctVarIn(J); } static proc ex48() { ring r8 = 0,(x1,x2,x3),dp; ideal I = x1^3-x2*x3, x2^2-x1*x3, x3^2-x1^2*x2; bfctVarIn(I); } static proc ex49 () { ring r9 = 0,(z1,z2,z3,z4),dp; ideal I = z3^2-z2*z4, z2^2*z3-z1*z4, z2^3-z1*z3; bfctVarIn(I); } static proc ex410() { LIB "toric.lib"; ring r = 0,(z(1..7)),dp; intmat A[3][7]; A = 6,4,2,0,3,1,0,0,1,2,3,0,1,0,0,0,0,0,1,1,2; ideal I = toric_ideal(A,"pt"); I = std(I); //ideal I = z(6)^2-z(3)*z(7), z(5)*z(6)-z(2)*z(7), z(5)^2-z(1)*z(7), // z(4)*z(5)-z(3)*z(6), z(3)*z(5)-z(2)*z(6), z(2)*z(5)-z(1)*z(6), // z(3)^2-z(2)*z(4), z(2)*z(3)-z(1)*z(4), z(2)^2-z(1)*z(3); bfctVarIn(I,1); // no result yet } */