Knuth–Eve algorithm

From HandWiki
Short description: Algorithm for evaluating polynomials


In computer science, the Knuth–Eve algorithm is an algorithm for polynomial evaluation. It preprocesses the coefficients of the polynomial to reduce the number of multiplications required at runtime.

Ideas used in the algorithm were originally proposed by Donald Knuth in 1962. His procedure opportunistically exploits structure in the polynomial being evaluated.[1] In 1964, James Eve determined for which polynomials this structure exists, and gave a simple method of "preconditioning" polynomials (explained below) to endow them with that structure.[2][note 1]

Algorithm

Preliminaries

Consider an arbitrary polynomial p[x] of degree n. Assume that n3. Define m such that: if n is odd then n=2m+1, and if n is even then n=2m+2.[2]

Unless otherwise stated, all variables in this article represent either real numbers or univariate polynomials with real coefficients.[1][2] All operations in this article are done over .[2]

Again, the goal is to create an algorithm that returns p(x) given any x. The algorithm is allowed to depend on the polynomial p itself, since its coefficients are known in advance.[1]

Overview

Key idea

Using polynomial long division, we can write

p(x)=q(x)(x2α)+(βx+γ),

where x2α is the divisor. Picking a value for α fixes both the quotient q and the coefficients in the remainder β and γ. The key idea is to cleverly choose α such that β=0, so that

p(x)=q(x)(x2α)+γ.

[4] This way, no operations are needed to compute the remainder polynomial, since it's just a constant. We apply this procedure recursively to q, expressing

p(x)=((q(x)(x2αm)+γm))(x2α1)+γ1.

After m recursive calls, the quotient q is either a linear or a quadratic polynomial. In this base case, the polynomial can be evaluated with (say) Horner's method.[1][4][5]

"Preconditioning"

For arbitrary p, it may not be possible to force β=0 at every step of the recursion.[1] Consider the polynomials pe and po with coefficients taken from the even and odd terms of p respectively, so that

p(x)=pe(x2)+xpo(x2).

If every root of po is real, then it is possible to write p in the form given above. Each αi is a different root of po, counting multiple roots as distinct.[4] Furthermore, if at least n1 roots of p lie in one half of the complex plane, then every root of po is real.[2]

Ultimately, it may be necessary to "precondition" p by shifting it — by setting p(x)p(x+t) for some t — to endow it with the structure that most of its roots lie in one half of the complex plane. At runtime, this shift has to be "undone" by first setting xxt.[2]

Preprocessing step

The following algorithm is run once for a given polynomial p.[1][4] At this point, the values of x that p will be evaluated on are not known.[1]

  • Let r1,,rn be the complex roots of p, sorted in descending order by real part
  • Choose any tRe(r2)
  • Set pp(x+t)

  • Let pe and po be the polynomials such that p(x)=pe(x2)+xpo(x2)
  • Let α1,αm be the roots of po. All of its roots will be real.

  • Initialize qp
  • For i1,,m:
    • Divide q by x2αi to get quotient q[x] and remainder γi. The remainder will be a constant polynomial — a number.
    • Set qq

  • Output: The derived values t, α1,,αm, and γ1,,γm; as well as the base-case polynomial q

Better choice of t

While any tRe(r2) can work, it is possible to remove one addition during evaluation if t is also chosen such that two roots of p(x+t) are symmetric about the origin. In that case, α1 can be chosen such that the shifted polynomial has a factor of x2α1, so γ1=0. It is always possible to find such a t.[2]

One possible algorithm for choosing t is:[citation needed]

  • If r1:
    • If r2: t=12(r1+r2)
    • Else: t=Re(r2)
  • Else: t=Re(r1)

Evaluation step

The following algorithm evaluates p at some, now known, point x.[1][2][4][5]

  • Set xxt
  • Let s=x2. Compute this once so it can be reused.

  • Compute yq(x) using Horner's method
  • For im,,2,1:
    • Let yy(sαi)+γi
  • Output: y

Assuming t is chosen optimally, γ1=0. So, the final iteration of the loop can instead run

yy(sαi),

saving an addition.[2]

Analysis

In total, evaluation using the Knuth–Eve algorithm for a polynomial of degree n requires n additions and n/2+2 multiplications, assuming t is chosen optimally.[2]


The Knuth–Eve algorithm is not well-conditioned.[6]

Footnotes

  1. Article published under the name J. Eve, which is associated with the name James Eve by the ACM Digital Library.[3]

References

  1. 1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.7 Knuth, Donald (December 1962). "Evaluation of polynomials by computer". Communications of the ACM 5 (12): 595–599. doi:10.1145/355580.369074. https://dl.acm.org/doi/10.1145/355580.369074. Retrieved 25 July 2025. 
  2. 2.00 2.01 2.02 2.03 2.04 2.05 2.06 2.07 2.08 2.09 Eve, James (December 1964). "The evaluation of polynomials". Numerische Mathematik 6 (1): 17–21. doi:10.1007/BF01386049. https://link.springer.com/article/10.1007/BF01386049. Retrieved 25 July 2025. 
  3. "James Eve (Newcastle University)" (in en). ACM Digital Library. doi:10.1145/contrib-81100250587/abs. https://dl.acm.org/do/10.1145/contrib-81100250587/abs/. 
  4. 4.0 4.1 4.2 4.3 4.4 Overill, Richard (12 June 1997). "Data parallel evaluation of univariate polynomials by the Knuth-Eve algorithm". Parallel Computing 23 (13): 2115–2127. doi:10.1016/S0167-8191(97)00096-3. https://www.sciencedirect.com/science/article/pii/S0167819197000963. Retrieved 25 July 2025. 
  5. 5.0 5.1 Muller, Jean-Michel (17 November 2016). Elementary functions: Algorithms and implementation. Boston, MA: Birkhäuser Boston. pp. 82–84. doi:10.1007/978-1-4899-7983-4_5. ISBN 978-1-4899-7983-4. https://link.springer.com/chapter/10.1007/978-1-4899-7983-4_5. Retrieved 25 July 2025. 
  6. Mesztenyi, Charles (January 1967). "Stable evaluation of polynomials". Journal of Research of the National Bureau of Standards, Section B 71B (1): 11–17. doi:10.6028/jres.071B.003. https://nistdigitalarchives.contentdm.oclc.org/digital/collection/p16009coll6/id/87671. Retrieved 25 July 2025. 
Cite error: <ref> tag with name "erickson2003" defined in <references> is not used in prior text.