Algebra`Horner`
This package applies Horner's rule to rearrange polynomials in Horner form. This is useful for efficient and stable numerical evaluation. Any polynomial can be rewritten in Horner, or nested, form.
Assume that can be calculated using only multiplications for integer . For a polynomial of degree , the Horner form requires multiplications and additions. The expanded form, however, requires multiplications, which is already more than twice as expensive for a polynomial of degree 10. Thus, one advantage of Horner form is that the work involved in exponentiation is distributed across addition and multiplication, resulting in savings of some basic arithmetic operations. Another advantage is that Horner form is more stable to evaluate numerically when compared with the expanded form. The reason for this is that each sum and product involve quantities which vary on a more evenly distributed scale.
Factoring polynomials in Horner or nested form.
The function Horner can be used for both univariate and multivariate polynomials. It also works efficiently on both sparse and dense polynomials and fractional exponents.
This loads the package.
In[1]:= <<Algebra`Horner`
By using the Variables function, Horner puts the polynomial into Horner form with respect to the variables identified using Variables.
In[2]:= Horner[11 x^3 4 x^2 + 7 x + 2]
Out[2]=
Here the coefficients of the polynomial are symbolic, so it is important to specify the variable .
In[3]:= Horner[a x^3 + b x^2 + c x + d, x]
Out[3]=
The exponents of a polynomial do not need to be integers, nor do they need to be of equal increments. The nested polynomial is scaled by the smallest exponent.
In[4]:= Horner[x^(1/3) + x + x^(3/2)]
Out[4]=
Exponents must be integers or rational numbers for the expression to be considered a valid polynomial.
In[5]:= Horner[1 + x^a, x]
Out[5]=
This constructs a bivariate polynomial.
In[6]:= bipoly = Sum[i j x^i y^j, {i, 0, 2}, {j, 0, 2}]
Out[6]=
Here the polynomial is put in Horner form first with respect to and then with respect to .
In[7]:= Horner[bipoly, {x, y}]
Out[7]=
The factorization that you obtain is dependent upon the ordering of the variables which you specify. In this example, the polynomial is put in Horner form first with respect to and then with respect to .
In[8]:= hornerbipoly = Horner[bipoly, {y, x}]
Out[8]=
Polynomials are often evaluated using machineprecision floatingpoint arithmetic. This is an inexact process because there are only finitely many numbers available in the representation. Therefore, different methods of evaluation will give different answers. The following example illustrates how numerical stability of a polynomial can be improved using Horner form. The example also illustrates how the time required for evaluation can be significantly reduced. Timings will be different for different computers, but the relative timings should be similar.
Define the derivative of a Legendre polynomial.
In[9]:= dpoly[x_] = D[LegendreP[50, x], x];
This plots the polynomial and gives the time required for evaluation. Plot uses machineprecision numbers to evaluate the polynomial by default. The plot is not very smooth due to the accumulation of representation errors associated with machineprecision computations.
In[10]:= Plot[Evaluate[dpoly[x]], {x, 1, 1.01}, PlotRange > All, PlotPoints > 2000]//Timing
Out[10]=
Here the polynomial is put in Horner form before plotting. The time required to rearrange the polynomial is usually small compared to the time required to render the plot, and the time required to evaluate the polynomial is significantly reduced. The plot reveals that the Horner form of the polynomial is much more stable to evaluate using machineprecision numbers.
In[11]:= Plot[Evaluate[Horner[dpoly[x]]], {x, 1, 1.01}, PlotRange > All, PlotPoints > 2000]//Timing
Out[11]=
Highprecision arithmetic can be used to overcome the problem of cancellation. However, highprecision arithmetic is expensive when compared with machineprecision arithmetic. It is preferable to seek alternative strategies, such as structural transformations, whenever possible. For example, the Horner form can be used to reduce the amount of extra precision required to obtain reliable solutions.
Factoring rational functions in Horner or nested form.
Several types of expressions are detected by Horner, and appropriately nested forms are produced. The following example illustrates how to improve the efficiency of evaluating a rational polynomial approximant to a function.
This loads the package Calculus`Pade`.
In[12]:= <<Calculus`Pade`
Here is the rational Padé approximation to Exp[Log[x]x].
In[13]:= approx = Pade[Exp[Log[x]x], {x, 0, 3, 2}]
Out[13]=
This puts the Padé approximation in Horner form. Both the numerator and the denominator are nested, and the variable to use in each is determined automatically.
In[14]:= Horner[approx]
Out[14]=
Expanding the terms illustrates that the Horner form is a nested restructuring of the original expression.
In[15]:= Expand[hornerbipoly  bipoly] == 0
Out[15]=
We can verify that the nested form is indeed a valid polynomial in the specified variables.
In[16]:= PolynomialQ[hornerbipoly, {x, y}]
Out[16]=
Nesting polynomials in Horner form is a general technique, but more efficient schemes for evaluating certain polynomials sometimes exist . Additional considerations also arise if polynomials are not dense or computations are to be performed in parallel .
References
The Art of Computer Programming Volume 2: Seminumerical Algorithms, D. E. Knuth, Second edition, Addison Wesley, London, 1981.
"Optimal Code for Serial and Parallel Computation", R. J. Fateman, Communications of the ACM, 12 (12), pp. 694695, 1969.
