From: Daniel Lichtblau Newsgroups: sci.math.symbolic Subject: Re: Q:implicitation of curves Date: Mon, 21 Jul 1997 12:27:21 -0500 Andreas Unterkircher wrote: > > Say I have plane curve c=(f_1(t),f_2(t)) and I want to compute the > implicit form p(x,y) of this curve. This could be done with Groebner > Basis techniques. Unfortunately this approach often is not useful due to > the complexity of this algorithm. My question is: are there any other > techniques to compute the implicit form of c (maybe some special cases). > I am especially interested in the case where c is a Bezier curve or some > spline. Any hints are welcome ! > > Andreas > > -- > Andreas Unterkircher Institut fuer Umformtechnik > ETH Zentrum / CLA, Tannenstr. 3,CH-8092 Zuerich, Switzerland > Phone: ++43 1 632 6621 www.ifu.ethz.ch The reply below is based on e-mail back and forth between AU and myself. A first attempt at an answer is that it depends on the degree of the data and type of coefficients. For low degree, Groebner techniques should suffice. Alternatively one can form a resultant. Below I show a simple example with Mathematica, running on a Pentium-Pro under Linux. In[1]:= f[t] = Array[a,4,0] . t^Range[0,3]; In[2]:= g[t] = Array[b,4,0] . t^Range[0,3]; In[3]:= Timing[imp1 = GroebnerBasis[{x-f[t], y-g[t]}, {x,y}, t, MonomialOrder->EliminationOrder][[1]];] Out[3]= {0.16 Second, Null} In[4]:= Timing[imp2 = Resultant[x-f[t], y-g[t], t];] (* gives -imp1 *) Out[4]= {0.14 Second, Null} For much higher degree input, it is not likely that these will be viable using symbolic coefficients. In our development version of Mathematica one can find a Groebner basis with specified bignum precision coefficients, so for inexact data this might be a reasonable approach. In[6]:= SeedRandom[1111]; In[7]:= f[t] = Table[Random[Real, 10, 150], {10}] . t^Range[0,9]; In[8]:= g[t] = Table[Random[Real, 10, 150], {10}] . t^Range[0,9]; In[9]:= Timing[imp = GroebnerBasis[{x-f[t], y-g[t]}, {x,y}, t, MonomialOrder->EliminationOrder, CoefficientDomain->InexactNumbers[150]][[1]];] Out[9]= {3.75 Second, Null} Alternative methods for implicitization might include interpolating the polynomial using some a priori degree bounds. One can check the text "Effective Polynomial Computation" by Zippel for more information on this. I think this reference may also provide information about the resultant-based approach. Based on this preliminary response, and an invitation to send an example of interest, AU sent me the following with a note to the effect that Mathematica 3.0 could not handle it. polys = {x - (110*t^2 - 495*t^3 + 1320*t^4 - 2772*t^5 + 5082*t^6 - 7590*t^7 + 8085*t^8 - 5555*t^9 + 2189*t^10 - 374*t^11), y - (22*t - 110*t^2 + 330*t^3 - 1848*t^5 + 3696*t^6 - 3300*t^7 + 1650*t^8 - 550*t^9 + 88*t^10 + 22*t^11)}; The development version of Mathematica will handle this provided one starts at 250 or so digits of precision (200 did not suffice). While I have not tried an exact computation, AU says it runs over 3 hours with no result. So clearly the inexact method is of interest here. Results below are abbreviated for purposes of, well, brevity. In[1]:= InputForm[Timing[poly250 = First[GroebnerBasis[ polys, {x,y}, t, MonomialOrder->EliminationOrder, CoefficientDomain->InexactNumbers[250]]]]] Out[1]//InputForm= {20.689999999999998*Second, 1.97552413410743382899217474906499999999999\ 99999999999580685`46.7148*^30*x^3 - 6.7273185232570481684919210047500000000000000000000000198397`47.46\ 78*^29*x^4 + 1.753245020195242250921029593500000000000000000000000003732\ 38`47.7702*^29*x^5 - 1.954692034344520501241951172400000000000000000000000\ 000215795`48.9375*^28*x^6 + 1.13993498927601163681082731799999999999999999\ 99999999999999887325`52.0201*^27*x^7 - .... The result contains decimals that appear to end generally in long sequences of 0's or 9's. After playing a bit with Rationalize, it seemed that the common denominator is 1024. In fact we do get a polynomial with integer coefficients below. In[7]:= newpoly250 = Rationalize[Expand[poly250*1024]] // InputForm Out[7]//InputForm= 2022936713326012240887986943042560*x^3 - 688877416781521732453572710886400*x^4 + 179532290067992806494313430374400*x^5 - 20016046431687889932717580005376*x^6 + ... This is in fact the correct exact implicit polynomial in x and y for the curve defined parametrically by polys. First evidence that this is so comes when one repeats the computation at 500 digits of precision, rationalizes, and gets the same polynomial. More importantly, one can prove that this is the implicit polynomial using the following facts. (i) The polynomial is irreducible (Factor[] will show this) (ii) Polynomial reduction with respect to a Groebner basis will reduce it to zero. As it happens, the input is already a Groebner basis with respect to the lexicographic term order in the variables {x,y,t} (this is the case for ideals defined by parametric equations). A simple reduction with respect to polys, in this term order, gives zero. In[14]:= PolynomialReduce[newpoly250, {x - (110*t^2 - 495*t^3 + 1320*t^4 - 2772*t^5 + 5082*t^6 - 7590*t^7 + 8085*t^8 - 5555*t^9 + 2189*t^10 - 374*t^11), y - (22*t - 110*t^2 + 330*t^3 - 1848*t^5 + 3696*t^6 - 3300*t^7 + 1650*t^8 - 550*t^9 + 88*t^10 + 22*t^11)}, {x,y,t}] [[2]] Out[14]= 0 If the coefficients of the inexact polynomial did not lend themselves to easy rationalization, one might still use its "skeleton," that is, the power-products in x and y that show up with non-zero coefficients. A bit of linear algebra, described next, will allow one to determine the exact coefficients. This is pretty much the method of basis conversion, and the idea goes back to Bruno Buchberger's '85 work "Groebner bases: an algorithmic method in polynomial ideal theory." One multiplies each power-product by a distinct new indeterminate to form a new polynomial poly. One then calls PolynomialReduce on the whole thing (with respect to any Groebner bases; polys will suffice), and gets a "reduced" polynomial in the original variables, with coefficients in terms of the new indeterminates. One sets this to zero, which means each coefficient must be zero; this gives a system of equations in the indeterminates. Solving this system and back-substituting the solutions into poly gives the exact polynomial in x and y. These last paragraphs were a bit of a digression on how one can use inexact methods to arrive at an exact solution. In practice the approximate solutions might suffice for one's work, and indeed the input data might well be approximate in the setting described (Bezier or spline curves). Daniel Lichtblau Wolfram Research danl@wolfram.com