//f defines a trimodal singularity for generic moduli ring R = 0,(x,y),ds; int a1,a2,a3=random(1,100),random(-100,1),random(1,100); poly f = (x^2-y^3)*(y+a1*x)*(y+a2*x)*(y+a3*x); ideal J = jacob(f); ideal I = f; // J:I, ideal of the closure of V(J) \ V(I) ideal Q = quotient(J,I); //the Hessian of f poly Hess = det(jacob(jacob(f))); //Hess is contained in Q iff NF is 0 reduce(Hess,std(Q));