Papers
Topics
Authors
Recent
Search
2000 character limit reached

Piecewise Polynomial Collocation

Updated 7 July 2026
  • Piecewise polynomial collocation is a discretization technique that approximates solutions by enforcing governing equations at carefully chosen collocation points on subintervals.
  • It employs various bases such as B-splines, cubic Hermite, and Gauss/Radau/Lobatto schemes to address differential, integral, and nonlocal problems.
  • Adaptive segmentation and specialized schemes enhance convergence, stability, and computational efficiency across a range of applications.

Piecewise polynomial collocation is a strong-form discretization technique in which an unknown solution is approximated by polynomials or splines on subintervals or elements, and the governing equation is enforced at selected collocation points. In current usage, the term covers piecewise linear, quadratic, cubic Hermite, B-spline, multivariate spline, simplex, sparse-grid, and Gauss/Radau/Lobatto constructions, and it appears in second-order ordinary differential equations, Volterra integral equations of the first kind, nonlocal diffusion, functional equations, elliptic partial differential equations, and direct optimal control (Darbyshire, 2018, Tynda et al., 2021, Lai et al., 2021).

1. Mathematical definition and approximation spaces

A collocation method seeks an approximate solution yh(x)y_h(x) in a finite-dimensional space of sufficiently smooth functions and imposes that the residual of the governing equation vanishes at a selected set of points. For a second-order ODE initial value problem,

y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',

the residual is

R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),

and collocation enforces R(xj;yh)=0R(x_j;y_h)=0 at chosen nodes xjx_j together with the initial data (Darbyshire, 2018). On a partition

a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,

the approximation space is typically piecewise polynomial of degree kk on each subinterval, often with at least C1C^1 continuity when second derivatives are needed. In spline form one writes

yh(x)=∑j=1Mcj φj(x),y_h(x)=\sum_{j=1}^M c_j\,\varphi_j(x),

with basis functions φj\varphi_j such as B-splines (Darbyshire, 2018).

The same structural idea recurs in other settings. For direct optimal control, the time interval is partitioned as

y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',0

and the state, control, and later costate are represented by piecewise polynomials on each subinterval. The framework explicitly includes Gauss, Radau, and Lobatto schemes, and in the Hermite–Simpson example the state is represented by a cubic Hermite polynomial while the control is piecewise linear (Zheng et al., 7 Apr 2026). For weakly regular Volterra equations of the first kind, the approximation is a continuous local spline that is polynomial of degree y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',1 on each subinterval y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',2, with local nodes chosen from endpoints and mapped Legendre roots (Tynda et al., 2021). For a nonlocal functional equation on y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',3, the discrete space is the set of continuous, piecewise linear functions satisfying the boundary conditions, and collocation is imposed at the mesh nodes (Caballero et al., 2024).

A central feature of the subject is therefore not a specific basis, but the combination of a local polynomial representation with pointwise enforcement of an operator equation.

2. Major formulations across equations and applications

The literature uses the same collocational principle in markedly different operator classes.

Setting Piecewise polynomial structure Representative detail
Second-order ODE IVP Piecewise polynomial or spline y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',4 on y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',5 Van der Pol collocation with initial data (Darbyshire, 2018)
Optimal control Piecewise polynomial state, control, and costate on a time mesh Gauss, Radau, Lobatto; Hermite–Simpson uses cubic Hermite state and piecewise linear control (Zheng et al., 7 Apr 2026)
Nonlocal diffusion Piecewise quadratic polynomial collocation Standard and shifted-symmetric quadratic schemes (Chen et al., 2020)
Volterra equations Polynomial spline collocation with Gauss-type quadrature Degree y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',6 local splines on y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',7 (Tynda et al., 2021)
Functional equations Piecewise linear collocation Nonlocal operator evaluated at non-grid points (Caballero et al., 2024)
PDEs on triangulations Multivariate polynomial splines over triangles or tetrahedra y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',8 spline collocation for Poisson and non-divergence elliptic PDEs (Lai et al., 2021)
High-dimensional UQ Local polynomials on simplices or localized sparse-grid supports Kink-aware simplex and hp-sparse-grid collocation (Fuchs et al., 2019, Wilka et al., 2024)

For ODEs, a standard choice is to enforce the differential equation at interior Gauss–Legendre points mapped into each subinterval,

y′′(x)=f(x,y(x),y′(x)),y(a)=y0,y′(a)=y0′,y''(x) = f\bigl(x,y(x),y'(x)\bigr), \qquad y(a)=y_0,\quad y'(a)=y_0',9

and to determine the unknown polynomial coefficients from the collocation equations plus initial or boundary conditions (Darbyshire, 2018). In direct collocation for Bolza problems, the discrete nonlinear program is built from collocation equations, endpoint conditions, and a discretized cost, and the continuous-time state, control, and costate trajectories are then reconstructed as piecewise polynomials from the discrete primal–dual solution (Zheng et al., 7 Apr 2026).

For PDEs and nonlocal models, the same pointwise logic survives but the polynomial space changes. In nonlocal diffusion, the unknown is approximated by piecewise quadratic polynomials on a spatial grid, and the nonlocal operator is enforced at collocation points; in multivariate spline collocation, the unknown belongs to R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),0 spline spaces over triangulations or tetrahedralizations, and the PDE is enforced at domain points of higher degree (Chen et al., 2020, Lai et al., 2021). In high-dimensional uncertainty quantification, simplex stochastic collocation builds local polynomial interpolants on a Delaunay triangulation, while adaptive hp sparse-grid collocation uses localized hierarchical basis functions with varying order on a hierarchical knot tree (Fuchs et al., 2019, Wilka et al., 2024).

3. Bases, continuity, and algebraic structure

The basis determines both regularity and the algebraic structure of the collocation system. In spline collocation for second-order ODEs, one usually requires at least R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),1 continuity so that R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),2 exists piecewise and the equation can be enforced at collocation points (Darbyshire, 2018). In polynomial spline collocation for Volterra equations, continuity across subinterval endpoints is enforced by shared nodal values, while on each R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),3 the local polynomial is determined by Lagrange interpolation conditions at R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),4 nodes (Tynda et al., 2021). In the functional-equation setting, the reason for choosing a continuous piecewise linear space is explicit: the operator evaluates R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),5 at non-grid points R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),6 and R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),7, so the approximation must be defined for all R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),8, not only at nodes (Caballero et al., 2024).

Several papers show that the governing operator can force nontrivial modifications of the polynomial representation. For the KdV–Burgers equation,

R(x;yh):=yh′′(x)−f(x,yh(x),yh′(x)),R(x;y_h):=y_h''(x)-f\bigl(x,y_h(x),y_h'(x)\bigr),9

the extended cubic B-spline basis is only R(xj;yh)=0R(x_j;y_h)=00, so the nonexistence of third order derivatives of the cubic B-splines forces a reduction of the order of the term R(xj;yh)=0R(x_j;y_h)=01 and a coupled system with R(xj;yh)=0R(x_j;y_h)=02 (Hepson et al., 2017). In the LGL integral formulation for optimal control, the right-hand side vector field is approximated by a Lagrange polynomial at Legendre–Gauss–Lobatto points, and the state is recovered through an integration matrix R(xj;yh)=0R(x_j;y_h)=03,

R(xj;yh)=0R(x_j;y_h)=04

while a derivative-like form is obtained by multiplying by the inverse of an appropriate full-rank block of the integration matrix (Abadia-Doyle et al., 16 Jun 2025).

The algebraic consequences can be decisive. Standard piecewise quadratic collocation for the nonlocal diffusion model produces a nonsymmetric indefinite system and does not satisfy the discrete maximum principle, whereas the shifted-symmetric piecewise quadratic polynomial collocation produces a symmetric positive definite system and satisfies the discrete maximum principle (Chen et al., 2020). For multivariate spline collocation, rather than constructing a basis for R(xj;yh)=0R(x_j;y_h)=05 directly, the method works in the discontinuous space R(xj;yh)=0R(x_j;y_h)=06 and enforces R(xj;yh)=0R(x_j;y_h)=07 and R(xj;yh)=0R(x_j;y_h)=08 continuity by linear smoothness matrices R(xj;yh)=0R(x_j;y_h)=09 and xjx_j0 (Lai et al., 2021).

4. Convergence, stability, and certification

The convergence theory of piecewise polynomial collocation is highly problem dependent, but several recurring patterns are explicit in the literature. For polynomial spline collocation applied to weakly regular Volterra equations, the error estimate is tied to spline approximation: xjx_j1 where xjx_j2 is the number of local interpolation nodes per element and the local polynomials have degree xjx_j3 (Tynda et al., 2021). For the functional equation

xjx_j4

the collocation error is bounded by the error of projecting the Lipschitz solution onto the piecewise linear polynomial space,

xjx_j5

and, when xjx_j6, the method is of the second order measured in the supremum norm: xjx_j7 (Caballero et al., 2024).

For nonlocal diffusion with horizon xjx_j8, the piecewise quadratic collocation analysis is more delicate because of horizon-grid interaction. The global error is

xjx_j9

if a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,0 is not set as a grid point, but it recovers

a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,1

when a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,2 is set as a grid point. The shifted-symmetric scheme is also proved asymptotically compatible, with global error

a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,3

(Chen et al., 2020). For multivariate spline collocation of Poisson and non-divergence elliptic PDEs, the analysis proceeds through PDE-based norm equivalences. On domains with uniformly positive reach, the method yields estimates of the form

a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,4

with a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,5 the discrete collocation residual tolerance (Lai et al., 2021).

In optimal control, the emphasis is not only convergence but certification. Starting from a discrete KKT point, one reconstructs piecewise polynomial state, control, and costate trajectories, evaluates dynamics, stationarity, and boundary residuals, and obtains a computable lower bound for the continuous second variation. The resulting a posteriori test is expressed in terms of the discrete reduced curvature a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,6, the residual-based measure a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,7, and computable constants such as a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,8, a=x0<x1<⋯<xN=b,a=x_0<x_1<\cdots<x_N=b,9, kk0, and kk1 (Zheng et al., 7 Apr 2026). This makes piecewise polynomial reconstruction not merely a discretization device but a bridge between solver output and continuous second-order theory.

5. Segmentation, adaptivity, and non-smooth regimes

A recurrent motive in the literature is that a single global collocation problem may be unnecessary or even impractical. For the Van der Pol equation,

kk2

the collocation method in its original form is described as impractical for generating these approximations due to the numerical costs of producing such an approximation, particularly as the stiffness of the equation increases with parameter kk3. The proposed segmented collocation splits the interval into segments and solves a sequence of smaller collocation problems, using the final values of one segment as initial conditions for the next; it is shown theoretically and numerically to be capable of generating equivalent approximations for much superior costs (Darbyshire, 2018).

For piecewise smooth functions with kinks, adaptivity changes both the mesh and the local polynomial degree. Adaptive hp-polynomial based sparse grid collocation constructs localised hierarchical basis functions with varying order on a hierarchical multivariate knot tree, uses hierarchical surplus as an error indicator to detect the non-smooth region, and adaptively refines collocation points there. The local polynomial degrees are optionally selected by a greedy approach or a kink detection procedure (Wilka et al., 2024). Simplex stochastic collocation uses a Delaunay triangulation, builds local polynomial interpolants on simplices, and—when each sample can be assigned to its smooth region—approximates the function separately on each smooth region and recovers a global order of convergence of kk4, where kk5 is the degree of the employed polynomials and kk6 the dimension of the parameter space (Fuchs et al., 2019).

A nearby development replaces deterministic local fits by randomized discrete least-squares reconstructions on each mesh element. These are computable projection operators onto piecewise polynomial spaces, defined via sampling and discrete least-squares polynomial approximations, and they exhibit almost optimal approximation properties in kk7 and kk8 (Storn, 4 Feb 2026). This suggests a broader viewpoint in which collocation, local regression, and projection onto piecewise polynomial spaces form a continuum rather than disjoint techniques.

6. Applications, limitations, and recurrent misconceptions

Piecewise polynomial collocation is used across a wide range of models: the Van der Pol equation with varying stiffness parameter kk9; Bolza optimal control with reconstructed state, control, and costate trajectories; linear nonlocal diffusion with volume constraints; weakly regular Volterra equations of the first kind; a nonlocal functional equation arising in behavioral sciences; the KdV–Burgers equation; Poisson’s equation; and second-order elliptic PDEs in non-divergence form (Darbyshire, 2018, Zheng et al., 7 Apr 2026, Chen et al., 2020, Tynda et al., 2021, Caballero et al., 2024, Hepson et al., 2017, Lai et al., 2021). The same broad label therefore covers spline collocation, direct collocation transcriptions, nonlocal strong-form schemes, and high-dimensional interpolation methods.

Several misconceptions are repeatedly corrected by these studies. First, high order is not by itself a guarantee of monotonicity or favorable linear algebra: standard piecewise quadratic collocation for nonlocal diffusion yields a nonsymmetric indefinite system and fails the discrete maximum principle, whereas a shifted-symmetric modification restores symmetry, positive definiteness, and the discrete maximum principle (Chen et al., 2020). Second, global formulations are not automatically preferable: for stiff second-order IVPs, segmented collocation can preserve the approximation while sharply reducing cost (Darbyshire, 2018). Third, additional support points are not always essential: in the integral form of Legendre–Gauss–Lobatto collocation for optimal control, a second integral form can be constructed by including an additional noncollocated support point, but such a point is superfluous and has no impact on the solution to the nonlinear programming problem (Abadia-Doyle et al., 16 Jun 2025).

The literature therefore presents piecewise polynomial collocation less as a single algorithm than as a design principle: choose a local polynomial space compatible with the operator, enforce the governing equations at carefully selected points, and adapt the segmentation, continuity, quadrature, or algebraic form to the regularity, stiffness, or non-smoothness of the problem at hand.

Topic to Video (Beta)

No one has generated a video about this topic yet.

Whiteboard

No one has generated a whiteboard explanation for this topic yet.

Follow Topic

Get notified by email when new papers are published related to Piecewise Polynomial Collocation.