Tricubic interpolation

In the mathematical subfield numerical analysis, tricubic interpolation is a method for obtaining values at arbitrary points in 3D space of a function defined on a regular grid. The approach involves approximating the function locally by an expression of the form $$f(x,y,z)=\sum _{i=0}^{3}\sum _{j=0}^{3}\sum _{k=0}^{3}a_{ijk}x^{i}y^{j}z^{k}.$$

This form has 64 coefficients $a_{ijk}$; requiring the function to have a given value or given directional derivative at a point places one linear constraint on the 64 coefficients.

The term tricubic interpolation is used in more than one context; some experiments measure both the value of a function and its spatial derivatives, and it is desirable to interpolate preserving the values and the measured derivatives at the grid points. Those provide 32 constraints on the coefficients, and another 32 constraints can be provided by requiring smoothness of higher derivatives.

In other contexts, we can obtain the 64 coefficients by considering a 3×3×3 grid of small cubes surrounding the cube inside which we evaluate the function, and fitting the function at the 64 points on the corners of this grid.

The cubic interpolation article indicates that the method is equivalent to a sequential application of one-dimensional cubic interpolators. Let $\mathrm {CINT} _{x}(a_{-1},a_{0},a_{1},a_{2})$ be the value of a monovariable cubic polynomial (e.g. constrained by values, $a_{-1}$, $a_{0}$, $a_{1}$, $a_{2}$ from consecutive grid points) evaluated at $x$. In many useful cases, these cubic polynomials have the form $\mathrm {CINT} _{x}(u_{-1},u_{0},u_{1},u_{2})=\mathbf {v} _{x}\cdot \left(u_{-1},u_{0},u_{1},u_{2}\right)$ for some vector $\mathbf {v} _{x}$ which is a function of $x$ alone. The tricubic interpolator is equivalent to:

$${\begin{aligned}s(i,j,k)&{}={\text{The value at grid point }}(i,j,k)\\t(i,j,z)&{}=\mathrm {CINT} _{z}\left(s(i,j,-1),s(i,j,0),s(i,j,1),s(i,j,2)\right)\\u(i,y,z)&{}=\mathrm {CINT} _{y}\left(t(i,-1,z),t(i,0,z),t(i,1,z),t(i,2,z)\right)\\f(x,y,z)&{}=\mathrm {CINT} _{x}\left(u(-1,y,z),u(0,y,z),u(1,y,z),u(2,y,z)\right)\end{aligned}}$$ where $i,j,k\in \{-1,0,1,2\}$ and $x,y,z\in [0,1]$.

At first glance, it might seem more convenient to use the 21 calls to $\mathrm {CINT}$ described above instead of the ${64\times 64}$ matrix described in Lekien and Marsden. However, a proper implementation using a sparse format for the matrix (that is fairly sparse) makes the latter more efficient. This aspect is even much more pronounced when interpolation is needed at several locations inside the same cube. In this case, the ${64\times 64}$ matrix is used once to compute the interpolation coefficients for the entire cube. The coefficients are then stored and used for interpolation at any location inside the cube. In comparison, sequential use of one-dimensional integrators $\mathrm {CINT} _{x}$ performs extremely poorly for repeated interpolations because each computational step must be repeated for each new location.

See also

• Cubic interpolation
• Bicubic interpolation
• Trilinear interpolation

External links

• Java/C++ implementation of tricubic interpolation
• C++ implementation of tricubic interpolation
• Python implementation
• NumPy implementation
• [1] einspline library