Monotone cubic interpolation

Last updated March 11, 2024

In the mathematical field of numerical analysis, monotone cubic interpolation is a variant of cubic interpolation that preserves monotonicity of the data set being interpolated.

Monotone cubic Hermite interpolation

Example showing non-monotone cubic interpolation (in red) and monotone cubic interpolation (in blue) of a monotone data set. MonotCubInt.png — Example showing non-monotone cubic interpolation (in red) and monotone cubic interpolation (in blue) of a monotone data set.

Monotone interpolation can be accomplished using cubic Hermite spline with the tangents $m_{i}$ modified to ensure the monotonicity of the resulting Hermite spline.

An algorithm is also available for monotone quintic Hermite interpolation.

Interpolant selection

There are several ways of selecting interpolating tangents for each data point. This section will outline the use of the Fritsch–Carlson method. Note that only one pass of the algorithm is required.

Let the data points be $(x_{k},y_{k})$ indexed in sorted order for $k=1,\,\dots \,n$ .

Compute the slopes of the secant lines between successive points:
$\delta _{k}={\frac {y_{k+1}-y_{k}}{x_{k+1}-x_{k}}}$
for $k=1,\,\dots \,n-1$ .
These assignments are provisional, and may be superseded in the remaining steps. Initialize the tangents at every interior data point as the average of the secants,
$m_{k}={\frac {\delta _{k-1}+\delta _{k}}{2}}$
for $k=2,\,\dots \,n-1$ .

For the endpoints, use one-sided differences:
$m_{1}=\delta _{1}\quad {\text{ and }}\quad m_{n}=\delta _{n-1}\,$ .
If $\delta _{k-1}$ and $\delta _{k}$ have opposite signs, set $m_{k}=0$ .
For $k=1,\,\dots \,n-1$ , where ever $\delta _{k}=0$ (where ever two successive $y_{k}=y_{k+1}$ are equal),
set $m_{k}=m_{k+1}=0,$ as the spline connecting these points must be flat to preserve monotonicity.
Ignore steps 4 and 5 for those $k\,$ .
Let
$\alpha _{k}=m_{k}/\delta _{k}\quad {\text{ and }}\quad \beta _{k}=m_{k+1}/\delta _{k}$ .
If either $\alpha _{k}$ or $\beta _{k}$ is negative, then the input data points are not strictly monotone, and $(x_{k},\,y_{k})$ is a local extremum. In such cases, piecewise monotone curves can still be generated by choosing $m_{k}=0\,$ if $\alpha _{k}<0$ or $m_{k+1}=0\,$ if $\beta _{k}<0$ , although strict monotonicity is not possible globally.
To prevent overshoot and ensure monotonicity, at least one of the following three conditions must be met:

(a) the function

$\phi _{k}=\alpha _{k}-{\frac {(2\alpha _{k}+\beta _{k}-3)^{2}}{3(\alpha _{k}+\beta _{k}-2)}}>0\,$ , or

(b)

\alpha _{k}+2\beta _{k}-3\leq 0\,

, or

(c)

2\alpha _{k}+\beta _{k}-3\leq 0\,

.

Only condition (a) is sufficient to ensure strict monotonicity:

\phi _{k}

must be positive.

One simple way to satisfy this constraint is to restrict the vector

(\alpha _{k},\,\beta _{k})

to a circle of radius 3. That is, if

\alpha _{k}^{2}+\beta _{k}^{2}>9\,

, then set

$\tau _{k}={\frac {3}{\sqrt {\alpha _{k}^{2}+\beta _{k}^{2}}}}\,$ ,

and rescale the tangents via

$m_{k}=\tau _{k}\,\alpha _{k}\,\delta _{k}\quad {\text{ and }}\quad m_{k+1}=\tau _{k}\,\beta _{k}\,\delta _{k}\,$ .

Alternatively it is sufficient to restrict

\alpha _{k}\leq 3

and

\beta _{k}\leq 3\,

. To accomplish this if

\alpha _{k}>3{\text{   or   }}\beta _{k}>3\,

, then set

m_{k}=3\,\delta _{k}\,

.

Cubic interpolation

After the preprocessing above, evaluation of the interpolated spline is equivalent to cubic Hermite spline, using the data $x_{k}$ , $y_{k}$ , and $m_{k}$ for $k=1,\,\dots \,n$ .

To evaluate at $x$ , find the index $k$ in the sequence where $x$ , lies between $x_{k}$ , and $x_{k+1}$ , that is: $x_{k}\leq x\leq x_{k+1}$ . Calculate

\Delta =x_{k+1}-x_{k}\quad {\text{ and }}\quad t={\frac {x-x_{k}}{\Delta }}

then the interpolated value is

f_{\text{interpolated}}(x)=y_{k}\cdot h_{00}(t)+\Delta \cdot m_{k}\cdot h_{10}(t)+y_{k+1}\cdot h_{01}(t)+\Delta \cdot m_{k+1}\cdot h_{11}(t)

where $h_{ii}$ are the basis functions for the cubic Hermite spline.

Example implementation

The following JavaScript implementation takes a data set and produces a monotone cubic spline interpolant function:

/* * Monotone cubic spline interpolation * Usage example listed at bottom; this is a fully-functional package. For * example, this can be executed either at sites like * https://www.programiz.com/javascript/online-compiler/ * or using nodeJS. */functionDEBUG(s){/* Uncomment the following to enable verbose output of the solver: *///console.log(s);}varj=0;varcreateInterpolant=function(xs,ys){vari,length=xs.length;// Deal with length issuesif(length!=ys.length){throw'Need an equal count of xs and ys.';}if(length===0){returnfunction(x){return0;};}if(length===1){// Impl: Precomputing the result prevents problems if ys is mutated later and allows garbage collection of ys// Impl: Unary plus properly converts values to numbersvarresult=+ys[0];returnfunction(x){returnresult;};}// Rearrange xs and ys so that xs is sortedvarindexes=[];for(i=0;i<length;i++){indexes.push(i);}indexes.sort(function(a,b){returnxs[a]<xs[b]?-1:1;});varoldXs=xs,oldYs=ys;// Impl: Creating new arrays also prevents problems if the input arrays are mutated laterxs=[];ys=[];// Impl: Unary plus properly converts values to numbersfor(i=0;i<length;i++){xs.push(+oldXs[indexes[i]]);ys.push(+oldYs[indexes[i]]);}DEBUG("debug: xs = [ "+xs+" ]")DEBUG("debug: ys = [ "+ys+" ]")// Get consecutive differences and slopesvardys=[],dxs=[],ms=[];for(i=0;i<length-1;i++){vardx=xs[i+1]-xs[i],dy=ys[i+1]-ys[i];dxs.push(dx);dys.push(dy);ms.push(dy/dx);}// Get degree-1 coefficientsvarc1s=[ms[0]];for(i=0;i<dxs.length-1;i++){varm=ms[i],mNext=ms[i+1];if(m*mNext<=0){c1s.push(0);}else{vardx_=dxs[i],dxNext=dxs[i+1],common=dx_+dxNext;c1s.push(3*common/((common+dxNext)/m+(common+dx_)/mNext));}}c1s.push(ms[ms.length-1]);DEBUG("debug: dxs = [ "+dxs+" ]")DEBUG("debug: ms = [ "+ms+" ]")DEBUG("debug: c1s.length = "+c1s.length)DEBUG("debug: c1s = [ "+c1s+" ]")// Get degree-2 and degree-3 coefficientsvarc2s=[],c3s=[];for(i=0;i<c1s.length-1;i++){varc1=c1s[i];varm_=ms[i];varinvDx=1/dxs[i];varcommon_=c1+c1s[i+1]-m_-m_;DEBUG("debug: "+i+". c1 = "+c1);DEBUG("debug: "+i+". m_ = "+m_);DEBUG("debug: "+i+". invDx = "+invDx);DEBUG("debug: "+i+". common_ = "+common_);c2s.push((m_-c1-common_)*invDx);c3s.push(common_*invDx*invDx);}DEBUG("debug: c2s = [ "+c2s+" ]")DEBUG("debug: c3s = [ "+c3s+" ]")// Return interpolant functionreturnfunction(x){// The rightmost point in the dataset should give an exact resultvari=xs.length-1;//if (x == xs[i]) { return ys[i]; }// Search for the interval x is in, returning the corresponding y if x is one of the original xsvarlow=0,mid,high=c3s.length-1,rval,dval;while(low<=high){mid=Math.floor(0.5*(low+high));varxHere=xs[mid];if(xHere<x){low=mid+1;}elseif(xHere>x){high=mid-1;}else{j++;i=mid;vardiff=x-xs[i];rval=ys[i]+diff*(c1s[i]+diff*(c2s[i]+diff*c3s[i]));dval=c1s[i]+diff*(2*c2s[i]+diff*3*c3s[i]);DEBUG("debug: "+j+". x = "+x+". i = "+i+", diff = "+diff+", rval = "+rval+", dval = "+dval);return[rval,dval];}}i=Math.max(0,high);// Interpolatevardiff=x-xs[i];j++;rval=ys[i]+diff*(c1s[i]+diff*(c2s[i]+diff*c3s[i]));dval=c1s[i]+diff*(2*c2s[i]+diff*3*c3s[i]);DEBUG("debug: "+j+". x = "+x+". i = "+i+", diff = "+diff+", rval = "+rval+", dval = "+dval);return[rval,dval];};};/*   Usage example below will approximate x^2 for 0 <= x <= 4.   Command line usage example (requires installation of nodejs):   node monotone-cubic-spline.js*/varX=[0,1,2,3,4];varF=[0,1,4,9,16];varf=createInterpolant(X,F);varN=X.length;console.log("# BLOCK 0 :: Data for monotone-cubic-spline.js");console.log("X"+"\t"+"F");for(vari=0;i<N;i+=1){console.log(F[i]+'\t'+X[i]);}console.log(" ");console.log(" ");console.log("# BLOCK 1 :: Interpolated data for monotone-cubic-spline.js");console.log("      x       "+"\t\t"+"     P(x)      "+"\t\t"+"    dP(x)/dx     ");varmessage='';varM=25;for(vari=0;i<=M;i+=1){varx=X[0]+(X[N-1]-X[0])*i/M;varrvals=f(x);varP=rvals[0];varD=rvals[1];message+=x.toPrecision(15)+'\t'+P.toPrecision(15)+'\t'+D.toPrecision(15)+'\n';}console.log(message);

Related Research Articles

In mathematical physics and mathematics, the Pauli matrices are a set of three $2 \times 2$ complex matrices that are Hermitian, involutory and unitary. Usually indicated by the Greek letter sigma, they are occasionally denoted by tau when used in connection with isospin symmetries.

The stress–energy tensor, sometimes called the stress–energy–momentum tensor or the energy–momentum tensor, is a tensor physical quantity that describes the density and flux of energy and momentum in spacetime, generalizing the stress tensor of Newtonian physics. It is an attribute of matter, radiation, and non-gravitational force fields. This density and flux of energy and momentum are the sources of the gravitational field in the Einstein field equations of general relativity, just as mass density is the source of such a field in Newtonian gravity.

In mathematics, particularly in linear algebra, tensor analysis, and differential geometry, the Levi-Civita symbol or Levi-Civita epsilon represents a collection of numbers; defined from the sign of a permutation of the natural numbers $1, 2, ..., n$ , for some positive integer $n$ . It is named after the Italian mathematician and physicist Tullio Levi-Civita. Other names include the permutation symbol, antisymmetric symbol, or alternating symbol, which refer to its antisymmetric property and definition in terms of permutations.

In probability theory and statistics, the gamma distribution is a versatile two-parameter family of continuous probability distributions. The exponential distribution, Erlang distribution, and chi-squared distribution are special cases of the gamma distribution. There are two equivalent parameterizations in common use:

With a shape parameter $k$ and a scale parameter $θ$
With a shape parameter $and an inverse scale parameter, called a rate parameter.$

<span class="mw-page-title-main">De Rham cohomology</span> Cohomology with real coefficients computed using differential forms

In mathematics, de Rham cohomology is a tool belonging both to algebraic topology and to differential topology, capable of expressing basic topological information about smooth manifolds in a form particularly adapted to computation and the concrete representation of cohomology classes. It is a cohomology theory based on the existence of differential forms with prescribed properties.

<span class="mw-page-title-main">Quartic function</span> Polynomial function of degree four

In algebra, a quartic function is a function of the form

In mathematics, the Hodge star operator or Hodge star is a linear map defined on the exterior algebra of a finite-dimensional oriented vector space endowed with a nondegenerate symmetric bilinear form. Applying the operator to an element of the algebra produces the Hodge dual of the element. This map was introduced by W. V. D. Hodge.

In the mathematical field of numerical analysis, spline interpolation is a form of interpolation where the interpolant is a special type of piecewise polynomial called a spline. That is, instead of fitting a single, high-degree polynomial to all of the values at once, spline interpolation fits low-degree polynomials to small subsets of the values, for example, fitting nine cubic polynomials between each of the pairs of ten points, instead of fitting a single degree-nine polynomial to all of them. Spline interpolation is often preferred over polynomial interpolation because the interpolation error can be made small even when using low-degree polynomials for the spline. Spline interpolation also avoids the problem of Runge's phenomenon, in which oscillation can occur between points when interpolating using high-degree polynomials.

In mathematics, bicubic interpolation is an extension of cubic spline interpolation for interpolating data points on a two-dimensional regular grid. The interpolated surface is smoother than corresponding surfaces obtained by bilinear interpolation or nearest-neighbor interpolation. Bicubic interpolation can be accomplished using either Lagrange polynomials, cubic splines, or cubic convolution algorithm.

In general relativity, the Gibbons–Hawking–York boundary term is a term that needs to be added to the Einstein–Hilbert action when the underlying spacetime manifold has a boundary.

In physics, Maxwell's equations in curved spacetime govern the dynamics of the electromagnetic field in curved spacetime or where one uses an arbitrary coordinate system. These equations can be viewed as a generalization of the vacuum Maxwell's equations which are normally formulated in the local coordinates of flat spacetime. But because general relativity dictates that the presence of electromagnetic fields induce curvature in spacetime, Maxwell's equations in flat spacetime should be viewed as a convenient approximation.

In applied mathematics, polyharmonic splines are used for function approximation and data interpolation. They are very useful for interpolating and fitting scattered data in many dimensions. Special cases include thin plate splines and natural cubic splines in one dimension.

In many-body theory, the term Green's function is sometimes used interchangeably with correlation function, but refers specifically to correlators of field operators or creation and annihilation operators.

In mathematics, the Schur orthogonality relations, which were proven by Issai Schur through Schur's lemma, express a central fact about representations of finite groups. They admit a generalization to the case of compact groups in general, and in particular compact Lie groups, such as the rotation group SO(3).

In numerical optimization, the nonlinear conjugate gradient method generalizes the conjugate gradient method to nonlinear optimization. For a quadratic function

Smoothing splines are function estimates, $, obtained from a set of noisy observations of the target, in order to balance a measure of goodness of fit of to with a derivative based measure of the smoothness of . They provide a means for smoothing noisy data. The most familiar example is the cubic smoothing spline, but there are many other possibilities, including for the case where is a vector quantity.$

A synchronous frame is a reference frame in which the time coordinate defines proper time for all co-moving observers. It is built by choosing some constant time hypersurface as an origin, such that has in every point a normal along the time line and a light cone with an apex in that point can be constructed; all interval elements on this hypersurface are space-like. A family of geodesics normal to this hypersurface are drawn and defined as the time coordinates with a beginning at the hypersurface. In terms of metric-tensor components $, a synchronous frame is defined such that$

In mathematics, a homogeneous distribution is a distribution S on Euclidean space Rⁿ or Rⁿ \ {0} that is homogeneous in the sense that, roughly speaking,

In mathematics, Ricci calculus constitutes the rules of index notation and manipulation for tensors and tensor fields on a differentiable manifold, with or without a metric tensor or connection. It is also the modern name for what used to be called the absolute differential calculus, developed by Gregorio Ricci-Curbastro in 1887–1896, and subsequently popularized in a paper written with his pupil Tullio Levi-Civita in 1900. Jan Arnoldus Schouten developed the modern notation and formalism for this mathematical framework, and made contributions to the theory, during its applications to general relativity and differential geometry in the early twentieth century.

The Einstein–Hilbert action for general relativity was first formulated purely in terms of the space-time metric. To take the metric and affine connection as independent variables in the action principle was first considered by Palatini. It is called a first order formulation as the variables to vary over involve only up to first derivatives in the action and so doesn't overcomplicate the Euler–Lagrange equations with higher derivative terms. The tetradic Palatini action is another first-order formulation of the Einstein–Hilbert action in terms of a different pair of independent variables, known as frame fields and the spin connection. The use of frame fields and spin connections are essential in the formulation of a generally covariant fermionic action which couples fermions to gravity when added to the tetradic Palatini action.

References

Fritsch, F. N.; Carlson, R. E. (1980). "Monotone Piecewise Cubic Interpolation". SIAM Journal on Numerical Analysis . SIAM. 17 (2): 238–246. doi:10.1137/0717021.
Dougherty, R.L.; Edelman, A.; Hyman, J.M. (April 1989). "Positivity-, monotonicity-, or convexity-preserving cubic and quintic Hermite interpolation". Mathematics of Computation . 52 (186): 471–494. doi: 10.2307/2008477 .

External links

GPLv2 licensed C++ implementation: MonotCubicInterpolator.cpp MonotCubicInterpolator.hpp

This page is based on this Wikipedia article
Text is available under the CC BY-SA 4.0 license; additional terms may apply.
Images, videos and audio are available under their respective licenses.