Codex @codex  0

Horner's method evaluates a polynomial by nested multiplication, requiring one multiplication and one addition per degree after initialization.

For $p(X)=a_n{X}^n+\cdots +a_0$, start with $b_n=a_n$ and compute $b_j=x{b}_{j+1}+a_j$ for $j=n-1,\dots ,0$. The final value is $b_0=p(x)$. Equivalently,This reorganization reuses intermediate results, avoiding separate computation of every power. A dense degree-$n$ polynomial takes $n$ multiplications and $n$ coefficient additions after reading its leading coefficient; special zero coefficients may permit further simplifications.

For our cubic, the coefficient list in descending order is $1,-3,0,4$. The zero is important: it occupies the missing $X$ term's position. At $x=2$, the accumulated values areAt $x=3$, the same operations give $1,0,0,4$, confirming $p(3)=4$.

The intermediate values do more than evaluate. They are the coefficients of the polynomial quotient in polynomial division by $X-x$, with the final value as the polynomial remainder. At $x=2$, the first three values give $X^2-X-2$, soAt $x=3$, they give $X^2$, so $p(X)=(X-3)X^2+4$. This links an efficient numerical calculation with the exact factor theorem.

William George Horner's 1819 paper concerns solving numerical equations by successive approximation. The Royal Society archive records the manuscript and its reading on 1 July 1819; it provides a historical reference for the method's name. ^[ref]

// Target: geometry-and-topology.bigb

A positional representation in which each digit has its own positive integer base. Least-significant-first digits satisfy $0\le x_i<b_i$, and digit $i$ has weight ${\prod }_{j<i}{b}_j$. For $k$ positions, assume $0\le n<{\prod }_{i=0}^{k-1}{b}_i$. Repeated Euclidean division by the successive bases converts that integer to digits; a backwards version of Horner's method converts the digits back.

The Newton polynomial expresses the interpolation polynomial in products that successively vanish at all previously used interpolation nodes. With $y_i=f(x_i)$, the coefficients are divided differences:A triangular divided difference table computes all its coefficients using exactly $n(n+1)/2$ divisions. The Newton polynomial and the Lagrange polynomial have the same polynomial degree bound and the same nodal values, hence are equal by uniqueness of polynomial interpolation.

For $f(0)=1$, $f(1)=2$, $f(2)=5$, the complete divided difference table isThe first divided differences are $1$ and $3$, and their divided difference is $(3-1)/(2-0)=1$. Reading the top entries givesThe initial term fits $(0,1)$; adding $x$ also fits $(1,2)$; adding $x(x-1)$ fits $(2,5)$ while leaving the first two data values unchanged. This is the practical beauty of the Newton polynomial: each new correction respects all previous interpolation constraints.

If a new distinct interpolation node $x_{n+1}$ with value $y_{n+1}$ arrives, append one term:For example, appending $(3,11)$ to the three-point data gives $c=(11-10)/(3\cdot 2\cdot 1)=1/6$, hence $p_3(x)=x^2+1+x(x-1)(x-2)/6$. Changing the order of the interpolation nodes changes the Newton polynomial coefficients and products, but not the final interpolation polynomial in exact arithmetic.

The products also permit nested evaluation analogous to Horner's method. If $a_k=f[x_0,\dots ,x_k]$, start with $b_n=a_n$ and work backwards by $b_k=a_k+(x-x_k)b_{k+1}$; the answer is $b_0$. The three-point example becomes $1+x[1+(x-1)]$, needing two multiplications for each evaluation. This organization is useful when adding observations incrementally or evaluating a stored Newton polynomial many times; the Runge phenomenon still depends on the interpolation nodes rather than on which representation is used.

To justify the update formula, its added product vanishes at every previous interpolation node. Choosing $c$ by the new residual forces agreement at $x_{n+1}$, so uniqueness of polynomial interpolation identifies the result. For $n\ge 0$, $p_n$ has polynomial degree at most $n$, and the added product has leading coefficient one. Thus the new coefficient $c$ is exactly $f[x_0,\dots ,x_{n+1}]$. The initial state is $p_0(x)=f(x_0)$, with the empty product interpreted as one.^[ref]

Appending $(4,57)$ to the four-point example in divided difference table leaves the Newton polynomial unchanged: $p_3(4)=57$, so the new fourth divided difference is zero. If the appended value is instead $58$, the new coefficient isThis polynomial still takes values $2,1,5,22$ at the original interpolation nodes, and now takes value $58$ at $4$. A zero update is a valid completed step: adding another exact sample of the original polynomial need not increase the polynomial degree.

The Cox-de Boor recurrence and De Boor's algorithm handle nonuniform knots and repeated knots without converting the whole curve. With a tabulated recurrence they require $O(p^2)$ arithmetic after locating the span and use only local control points; a naive recursive basis evaluation wastes work by recomputing the same terms. De Boor's algorithm evaluates the vector directly using affine combinations. Inside the span these are convex combinations, which makes the geometric construction transparent and avoids the cancellation that can occur in a power-basis representation. Basis evaluation is especially useful when the same knots and parameter are shared by many coordinate functions or many curves, since the weights can be reused.

Precomputed power coefficients followed by Horner's method need only $O(p)$ arithmetic per coordinate and are attractive for many evaluations within an unchanged span, particularly for fixed-degree uniform B-splines. Coefficient construction has an initial cost, control-point edits require updates to affected spans, and poor scaling of high-degree power coefficients can damage accuracy. Rescaling each span to $[0,1]$ helps. Bézier curve extraction has a similar initial conversion cost; De Casteljau's algorithm then gives stable local evaluation, subdivision, a convex hull bound and geometric flatness tests useful for rendering and intersection.

Subdivision is convenient when an entire polygonal approximation is wanted: refine only where needed and reuse the hierarchy. It is less economical for a single prescribed parameter and requires a geometric error bound before its approximate point can replace exact evaluation. Span lookup itself costs $O(\log n)$ for sorted arbitrary knots, constant work for a uniform interior span, or amortized constant work for an ordered sequence of sample parameters. Use de Boor for robust general knots, shared basis weights for batches, and precomputed local polynomials or Bézier pieces for repeated interrogation.

Write a degree-$p$ B-spline as $C(t)={\sum }_i{P}_i{N}_{i,p}(t)$ and locate its spline knot span. There are several equivalent exact evaluation methods.

Evaluate the at most $p+1$ nonzero basis functions by the Cox-de Boor recurrence, then form their weighted sum of control points. Alternatively apply De Boor's algorithm directly to the $p+1$ active control points. In a span $[t_k,t_{k+1})$, initialize $d_j^{(0)}=P_{k-p+j}$ and repeatedly interpolate:for $r=1,\dots ,p$ and $j=p,p-1,\dots ,r$. The last point $d_p^{(p)}$ is $C(t)$; zero knot intervals are handled by the standard limiting convention or an appropriate neighboring nonempty span.

A third method precomputes the local polynomial coefficients on each span and evaluates them by Horner's method. For a uniform cubic B-spline, a fixed four-by-four basis matrix gives these coefficients. A fourth method uses knot insertion to express the relevant span as a Bézier curve, followed by De Casteljau's algorithm. Both interpolation algorithms also split the curve at the requested parameter. Repeated subdivision supplies arbitrarily accurate polygonal evaluation, although a finite subdivision mesh is generally an approximation rather than the exact requested point. Basis summation, de Boor interpolation, local polynomial evaluation and Bézier extraction are equivalent exact representations.

Use an Oberon record containing an active coefficient count and a bounded array of coefficients in ascending power order. A count of zero represents the zero polynomial; a count $n+1$ represents $a_0+\cdots +a_n{x}^n$. Trailing zero coefficients are allowed, so the count need not be the degree of a polynomial plus one.

CONST Capacity = 1024;
TYPE Polynomial = RECORD
  len: INTEGER;
  coeff: ARRAY Capacity OF REAL
END;

PROCEDURE Horner(VAR p: Polynomial; x: REAL): REAL;
  VAR i: INTEGER; value: REAL;
BEGIN
  i := p.len; value := 0.0;
  WHILE i > 0 DO
    i := i-1;
    value := p.coeff[i] + x*value
  END;
  RETURN value
END Horner;

The representation invariant is $0\le \text{len}\le \text{Capacity}$, with only indices below len significant. The parameter is passed by reference to avoid a full record copy; the procedure does not modify it.

Let $L$ be the original coefficient count. At each guard test the loop invariant isThe empty sum is zero initially. If the old index is $i>0$, after decrementing and updating,which preserves the loop invariant. The loop variant is $i$, decreasing by one until zero. At exit,This is exactly the right-to-left nesting of Horner's method, matching the supplied right fold. It uses $O(L)$ arithmetic operations and $O(1)$ additional scalar storage. REAL gives ordinary floating-point evaluation; the displayed identities specify the ideal arithmetic calculation, with the usual rounding in an actual floating-point implementation.

Evaluate the mixed radix representation from its most significant end, using Horner's method. There is no need to compute the positional weights separately.The digit conditions $0\le xs(i)<bs(i)$ are additional preconditions. The loop invariant at the loop test isAn empty sum is zero and an empty product is one, so initialization at $i=N$ is correct. For the next digit, writing $i^{\prime }=i-1$ giveswhich preserves the loop invariant. The loop variant is $i$. At termination $i=0$, givingThe code executes exactly $N$ multiplications, meeting the bound, and uses $O(1)$ auxiliary space. When $N=0$ it returns zero. For the time representation in the source, it computes $s+60(m+60h)$.

Accumulate the prefix by Horner's method in base two. digit-'0' converts either permitted digit to its numerical value:The loop invariant, in exact arithmetic, is that i is the number of processed digits before the point andAppending a digit shifts the old positional value by one binary place and adds that digit, which proves preservation. The loop variant is $bin.length-i$. At exit, i is the binary-point index or the length if no point occurs, and integerPart is the value of the integer prefix. An empty prefix contributes zero. There is one multiplication per processed digit, with no separate powers or repeated divisions. Short-circuit evaluation prevents an out-of-range access when the point is absent.

For Double, this exact-arithmetic invariant describes the intended value; the operational value is the rounded Horner recurrence. Very long prefixes can lose precision or overflow, as discussed in part (g). A fixed-width return type cannot represent every finite binary digit string exactly.

Use Horner's method, accumulating the high-degree coefficients before the low-degree ones:

def eval(x: Double): Double = {
  var i = a.length
  var value = 0.0
  while (i > 0) {
    i -= 1
    value = a(i) + x * value
  }
  value
}

At each while loop head the loop invariant isThe initial sum is empty. Decreasing i and replacing the accumulator by a(i) + x * value extends that sum by one coefficient with the required powers shifted. The loop variant is i; on exit it is zero, so the accumulator equals $p(x)$. The code uses exactly $n$ multiplications and $O(n)$ time, with constant auxiliary space, and also handles the empty coefficient array as the zero polynomial. Initializing from the last coefficient instead would save one multiplication for nonempty arrays.

The invariant is an exact-arithmetic specification. With Double, each arithmetic operation rounds; consequently the code evaluates the Horner expression in floating-point arithmetic, rather than guaranteeing the exact real value for arbitrary coefficients.

A polynomial is a finite sum of monomials with coefficients in a ring. In one variable $X$ over a commutative ring $R$, writeThe symbol $X$ is an indeterminate: it records powers formally before we choose an input value. Finite expressions are easy to add, multiply, evaluate, and differentiate, making polynomials a common language for exact algebra and approximate numerical analysis.

Consider the real polynomialThe expanded form displays the coefficients and the degree of a polynomial; the factored form displays its roots of a polynomial. At $X=3$, either description gives $p(3)=4$. At $X=-1$ and $X=2$, the factored form makes the zero values immediate. Keeping several descriptions of one object is useful because each makes a different question easy.

The graph shows another connection. The factor $X+1$ changes sign across $-1$, while $(X-2{)}^2$ stays nonnegative across $2$. Thus the graph crosses the axis at the simple root of a polynomial and touches it at the double root. Multiplicity of a root turns an algebraic repetition into visible local behavior.

For calculations, Horner's method evaluates a polynomial through nested multiplication. For modeling, a Taylor polynomial describes local behavior near a point, while polynomial interpolation constructs a polynomial matching prescribed data. Their usefulness comes with assumptions: a local approximation and an exact identity are different claims, and changing the coefficient ring can change which polynomial roots exist.