Mathematical foundations of Cubic Spline Interpolation.
While global polynomial interpolation attempts to fit a single high-degree polynomial through all data points—often incurring Runge's wild boundary oscillations—spline interpolation adopts a piecewise paradigm. Originating from the physical flexible wooden splines used by naval draftsmen to draw smooth ship hulls, mathematical cubic splines connect adjacent data pairs using cubic polynomials joined with C² continuity (matching values, first derivatives, and second derivatives). The resulting system of continuity equations is tridiagonal and strictly diagonally dominant, enabling direct solution in linear O(N) time via the Thomas algorithm. In this laboratory, we construct natural and clamped cubic splines, examine the tridiagonal curvature system, and confirm the optimal O(h⁴) error convergence rate.
Viva voce questions & analytical theory
What is the formal definition of a cubic spline interpolant?
A cubic spline S(x) for data (x_i, y_i) is a function such that on each sub-interval [x_i, x_{i+1}], S(x) is a cubic polynomial, S(x_i) = y_i, and S(x) is twice continuously differentiable across the entire interval (S ∈ C²[a, b]).
Why is the tridiagonal matrix for cubic splines guaranteed to be nonsingular?
For interior knots, the diagonal entry is 2(h_{i-1} + h_i) while the sum of the off-diagonal entries is h_{i-1} + h_i. Thus the matrix is strictly diagonally dominant by a factor of 2, guaranteeing existence, uniqueness, and numerical stability by the Levy-Desplanques theorem.
What is Holladay's minimum curvature theorem for natural cubic splines?
Holladay's theorem states that among all interpolating functions g ∈ C²[a, b] satisfying g(x_i) = y_i, the natural cubic spline S(x) minimizes the total bending energy functional J(g) = ∫_a^b [g''(x)]² dx.
What is the difference between natural and clamped boundary conditions?
Natural boundary conditions set the second derivative to zero at the endpoints (S''(a) = S''(b) = 0), representing a free flexible beam. Clamped boundary conditions set the first derivative to specified values (S'(a) = f'(a), S'(b) = f'(b)).
How does the asymptotic error convergence rate differ between natural and clamped splines?
Clamped splines achieve optimal O(h⁴) uniform error convergence throughout the entire domain: ||f - S||_∞ ≤ (5/384) h⁴ ||f^{(4)}||_∞. Natural splines degrade to O(h²) near the endpoints unless f''(a) = f''(b) = 0.
What is the computational complexity of solving a cubic spline with N knots?
Because the continuity matrix is tridiagonal and strictly diagonally dominant, the Thomas algorithm solves it in O(N) operations and O(N) storage, without requiring pivoting.
Why are cubic splines preferred over high-degree global polynomials in computer graphics and CAD?
Splines completely eliminate Runge's phenomenon, provide localized control (adjusting one knot affects primarily nearby intervals), minimize bending curvature, and yield smooth C² visual appearance.
Algorithm precautions and stability rules
- Ensure knot coordinates x_0 < x_1 < ... < x_N are strictly monotonically increasing before assembling the tridiagonal matrix.
- Verify diagonal dominance when setting custom boundary conditions; clamped conditions preserve strict diagonal dominance.
- In the Thomas algorithm, check for division by near-zero pivot elements, though strict diagonal dominance theoretically guarantees non-zero pivots.
- When querying an interpolant at point x, use binary search (std::upper_bound) in O(log N) time to locate the knot interval [x_i, x_{i+1}].
- Store the moments M_i and compute local cubic coefficients a_i, b_i, c_i, d_i explicitly to allow fast repeated point evaluations.
- Ensure floating point precision buffer remains 64-bit IEEE 754 to prevent error accumulation in the backward substitution phase.
Theoretical limitations
- Natural cubic splines exhibit boundary layer degradation to O(h²) error near endpoints unless the underlying function naturally satisfies f''(a) = f''(b) = 0.
- Standard cubic splines are not shape-preserving: they can introduce spurious overshoot or undershoot in regions with sharp gradients, requiring tension splines or monotonic Hermite splines.
- Splines do not easily handle parametric curves with self-intersections without parameterizing by chord length.