Open In App

Natural Cubic spline

Improve
Improve
Like Article
Like
Save
Share
Report

Cubic Spline:

The cubic spline is a spline that uses the third-degree polynomial which satisfied the given m control points. To derive the solutions for the cubic spline, we assume the second derivation 0 at endpoints, which in turn provides a boundary condition that adds two equations to m-2 equations to make them solvable. The system of equations for the Cubic spline for 1-dimension can be given by:

S_i (x) = a_i + b_i *x + c_i * x^{2} + d_i x^{3}

We take a set of points [xi, yi] for i = 0, 1, …, n for the function y = f(x). The cubic spline interpolation is a piecewise continuous curve, passing through each of the values in the table. 

  • Following are the conditions for the spline of degree K=3:
    • The domain of s is in intervals of [a, b].
    • S, S’, S” are all continuous function on [a,b].

S(x) =\begin{bmatrix} S_0 (x), x \epsilon [x_0, x_1] \\ S_1 (x), x \epsilon [x_1, x_2] \\ ... \\ ... \\ ... \\ S_{n-1} (x), x \epsilon [x_1, x_2] \end{bmatrix}

Here Si(x) is the cubic polynomial that will be used on the subinterval [xi, xi+1]

Since, there are n intervals and 4 coefficients of each equation, for that we require a 4n parameters to solve the spline, we can get 2n equation from the fact that each of the cubic spline equation should satisfy the value at both ends:

S_i(x_i) = y_i \\ S_i(x_{i+1}) = y_{i+1}

The above cubic spline equations should not only continuous and differentiable but also have defined first and second derivatives that are also continuous on control points.

S_{i-1}^{'}(x_i) = S_{i}^{'}(x_i)

and

S_{i-1}^{''}(x_i) = S_{i}^{''}(x_i)

For 1, 2, 3…n-1 provides the 2n -2 equation constraints. So, we need 2 more equations to solve the above cubic spline. For that, we will use some natural boundary conditions.

Natural Cubic Spline:

In Natural cubic spline, we assume that the second derivative of the spline at boundary points is 0:

S^{''}(x_0) = 0 \\ S^{''}(x_n) = 0    

Now, since the S(x) is a third-order polynomial we know that S”(x) is a linear spline which interpolates. Hence, first, we construct S”(x) then integrate it twice to obtain S(x). 

Now, let’s assume t_i = x_i for i= 0, 1,…n, and z_i = S^{''}(x_i), i=0,1,2..n     and from the natural boundary condition z_0= z_n =0   . Twice differentiating a cubic spline gives a linear spline which can be written as:

S_i^{''}(x) = z_i \frac{x- t_{i+1}}{t_i - t_{i+1}} + z_{i+1} \frac{x - t_i}{t_{i+1} - t_i}

where, 

h_i = t_{i+1} - t_i ; t=0,1, 2,...,n

Now, the equation becomes:

S''(x) = z_{i+1} \frac{x- t_i}{h_i} + z_{i} \frac{t_{i+1} - x}{h_i}

Integrating this equation two times to obtain the cubic spline:

S(x) =\frac{z_{i+1}}{6h_i}(x - t_i)^{3} + \frac{z_{i}}{6h_i}( t_{i+1} - x)^{3} + C_i(x_i -t)  + D_i(t_{i+1} -x)

where, 

C_i = \frac{y_{i+1}}{h_i} - \frac{z_{i+1} * h_i}{6} \\ D_i = \frac{y_{i}}{h_i} - \frac{z_{i} * h_i}{6}

Now, to check for the continuity of derivative at t_i i.e S^{‘}_{i}(t_i) = S^{‘}_{i-1}(t_{i}). We first need to find derivatives and put that condition:

S'(x) = \frac{z_{i+1}}{2h_i}(x - t_i)^{2} - \frac{z_{i}}{2h_i}( t_{i+1} - x)^{2} + \frac{1}{h_i}(y_{i+1}-y_i)- \frac{h_i}{6}(z_{i+1} -z_i)

where,

b_i = \frac{1}{h_i}(y_{i+1}-y_i)

Putting the above equation for the continuity and solving it gives the following equation:

6(b_i -b_{i-1}) =  h_{i-1}z_{i-1} + 2(h_{i-1} + h_i)z_i + h_{i}z_{i+1}

Take, v_i = 6(b_i – b_{i-1}) and the above equation can be written as form of matrix:

\begin{bmatrix} 1&  0 &  0&  .&  .&  .&  .&  .&  .&  .&  .&  0& \\ h_0&  2(h_0 + h_1)&  h_1&  &  &  &  &  &  &  &  &  .& \\ .&  h_1 &  2(h_1 + h_2)&  h_2&  &  &  &  &  &  &  &  .& \\ .&  .&  .&  .&  &  &  &  &  &  &  &  .& \\ .&  &  &  .&  .&  .&  &  &  &  &  &  .& \\ .&  &  &  &  .&  .&  .&  &  &  &  &  .& \\ .&  &  &  &  &  .&  .&  .&  &  &  &  .& \\ .&  &  &  &  &  &  .&  .&  .&  &  &  .& \\ .&  &  &  &  &  &  & & .&  h_{n-2}&  2(h_{n-1} + h_{n-2})&  h_{n-1}&   \\ 0&  .&  .&  .&  .&  .&  .&  .&  .&  0&  0&  1& \\ \end{bmatrix} \begin{bmatrix} z_0 \\ z_1\\ z_2\\ .\\ .\\ .\\ .\\ .\\ z_{n-1}\\ z_{n} \end{bmatrix} = \begin{bmatrix} 0\\ v_1 \\ .\\ .\\ .\\ .\\ .\\ .\\ v_{n-1}\\ 0 \end{bmatrix}

Implementation

  • In this implementation, we will be performing the spline interpolation for function f(x) = 1/x for points b/w 2-10 with cubic spline that satisfied natural boundary condition.

Python3

#imports
import matplotlib.pyplot as plt
import numpy as np
from scipy.interpolate import CubicSpline, interp1d
plt.rcParams['figure.figsize'] =(12,8)
  
x = np.arange(2,10)
y = 1/(x)
# apply cubic spline interpolation
cs = CubicSpline(x, y, extrapolate=True)
# apply natural cubic spline interpolation
ns = CubicSpline(x, y,bc_type='natural', extrapolate=True)
  
# Apply Linear interpolation
linear_int = interp1d(x,y)
  
xs = np.arange(2, 9, 0.1)
ys = linear_int(xs)
  
# plot linear interpolation
plt.plot(x, y,'o', label='data')
plt.plot(xs,ys,  label="interpolation", color='green')
plt.legend(loc='upper right', ncol=2)
plt.title('Linear Interpolation')
plt.show()
  
# define a new xs
xs = np.arange(1,15)
#plot cubic spline and natural cubic spline
plt.plot(x, y, 'o', label='data')
plt.plot(xs, 1/(xs), label='true')
plt.plot(xs, cs(xs), label="S")
plt.plot(xs, ns(xs), label="NS")
plt.plot(xs, ns(xs,2), label="NS''")
plt.plot(xs, ns(xs,1), label="NS'")
  
plt.legend(loc='upper right', ncol=2)
plt.title('Cubic/Natural Cubic Spline Interpolation')
plt.show()
  
# check for boundary condition
print("Value of double differentiation at boundary conditions are %s and %s"
      %(ns(2,2),ns(10,2)))

                    
Output:
Value of double differentiation at boundary conditions are -1.1102230246251565e-16 and -0.00450262550915248

Linear Interpolation

Cubic/Natural Spline Interpolation

References:



Last Updated : 18 Jul, 2021
Like Article
Save Article
Previous
Next
Share your thoughts in the comments
Similar Reads