Math 441

Interpolation: Cubic Splines Implementation

Setup

Given: \(x_i\), \(y_i\), \(i = 0, 1, \dots, n\).

  • \(h_i = x_i - x_{i-1}\), \(i = 1, 2, \dots, n\).
  • \(s_i = \frac{y_i - y_{i-1}}{h_i}\), \(i = 1, 2, \dots, n\).
  • \(\delta_i = h_i + h_{i+1}\), \(i = 1, 2, \dots, n-1\).
  • \(B_i = 6\left(s_{i+1} - s_i\right)\), \(i = 1, 2, \dots, n-1\).

These lead to a tridiagonal system of equations.

Gauss Elimination

\[ \begin{bmatrix} 1 & 0 & 0 & & & & & \\ h_1 & 2\delta_1 & h_2 & & & & & \\ &h_2 & 2\delta_2 & h_3 & & & & \\ & & \ddots & \ddots & \ddots & & & \\ & & & h_{i} & 2\delta_i & h_{i+1} & & \\ & & & & \ddots & \ddots & \ddots & \\ & & & & & h_{n-1} & 2\delta_{n-1} & h_{n} \\ & & & & & 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} z_0\\z_1\\z_2\\\vdots\\z_i\\\vdots\\z_{n-1}\\z_n \end{bmatrix} = \begin{bmatrix} 0\\B_1\\B_2\\\vdots\\B_i\\\vdots\\B_{n-1}\\0 \end{bmatrix} \]

  • \(h'_1 = 0\), \(h'_i = \frac{h_i}{2\delta_{i-1} - h'_{i-1}h_{i-1}}\), \(i = 2, 3, \dots, n\)
  • \(\beta_0 = 0\), \(\beta_k = \frac{B_k - h_k\beta_{k-1}}{2\delta_k - h'_{k}h_{k}}\), \(k = 1, 2, \dots, n-1\), \(\beta_n=0\).

Back-substitution

\[ \begin{bmatrix} 1 & h'_1 & & & & & & \\ & 1 & h'_2 & & & & & \\ & & 1 & h'_3 & & & & \\ & & & \ddots & \ddots & & & \\ & & & & 1 & h'_{i+1} & & \\ & & & & & \ddots & \ddots & \\ & & & & & & 1 & h'_{n} \\ & & & & & & & 1 \end{bmatrix} \begin{bmatrix} z_0\\z_1\\z_2\\\vdots\\z_i\\\vdots\\z_{n-1}\\z_n \end{bmatrix} = \begin{bmatrix} \beta_0\\\beta_1\\\beta_2\\\vdots\\\beta_i\\\vdots\\\beta_{n-1}\\\beta_n \end{bmatrix} \]

  • \(z_n = 0\), \(z_k = \beta_k - h'_{k+1}z_{k+1}\), \(k = n-1, n-2, \dots, 0\).

Finally

  • \(a_i = y_{i-1}\), \(y = 1, 2, \dots, n\).
  • \(b_i = s_i - \frac{z_i + 2z_{i-1}}{6}h_i\), \(i = 1, 2,\dots, n\)
  • \(c_i = \frac{z_{i-1}}{2}\), \(i = 1, 2,\dots, n\)
  • \(d_i = \frac{z_i - z_{i-1}}{6h_i}\), \(i = 1, 2,\dots, n\)

\(h_i\)s and \(\delta\)s

  • \(h_i = x_i - x_{i-1}\), \(i = 1, 2, \dots, n\).
  • \(\delta_i = h_i + h_{i+1}\), \(i = 1, 2, \dots, n-1\).

  • What are the “inputs”?

  • What are the indices of \(h_i\)s and \(\delta\)s?

    • Do they start with 0 or 1?
    • What do they end with?

  • Come up with an example of input data and output data.

\(s_i\)s and \(B_i\)s

  • \(s_i = \frac{y_i - y_{i-1}}{h_i}\), \(i = 1, 2, \dots, n\).
  • \(B_i = 6\left(s_{i+1} - s_i\right)\), \(i = 1, 2, \dots, n-1\).

  • What are the “inputs”?

  • What are the indices of \(s_i\)s and \(B_i\)s?

    • Do they start with 0 or 1?
    • What do they end with?

  • Come up with an example of input data and output data.

\(h'_i\)s and \(\beta\)s

\[ \begin{bmatrix} 1 & 0 & 0 & & & & & \\ h_1 & 2\delta_1 & h_2 & & & & & \\ &h_2 & 2\delta_2 & h_3 & & & & \\ & & \ddots & \ddots & \ddots & & & \\ & & & h_{i} & 2\delta_i & h_{i+1} & & \\ & & & & \ddots & \ddots & \ddots & \\ & & & & & h_{n-1} & 2\delta_{n-1} & h_{n} \\ & & & & & 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} z_0\\z_1\\z_2\\\vdots\\z_i\\\vdots\\z_{n-1}\\z_n \end{bmatrix} = \begin{bmatrix} 0\\B_1\\B_2\\\vdots\\B_i\\\vdots\\B_{n-1}\\0 \end{bmatrix} \]

  • \(h'_1 = 0\), \(h'_i = \frac{h_i}{2\delta_{i-1} - h'_{i-1}h_{i-1}}\), \(i = 2, 3, \dots, n\)
  • \(\beta_0 = 0\), \(\beta_k = \frac{B_k - h_k\beta_{k-1}}{2\delta_k - h'_{k}h_{k}}\), \(k = 1, 2, \dots, n-1\), \(\beta_n=0\).

Again, what are the inputs, what are the indices, give an example…

\(z\)s

\[ \begin{bmatrix} 1 & h'_1 & & & & & & \\ & 1 & h'_2 & & & & & \\ & & 1 & h'_3 & & & & \\ & & & \ddots & \ddots & & & \\ & & & & 1 & h'_{i+1} & & \\ & & & & & \ddots & \ddots & \\ & & & & & & 1 & h'_{n} \\ & & & & & & & 1 \end{bmatrix} \begin{bmatrix} z_0\\z_1\\z_2\\\vdots\\z_i\\\vdots\\z_{n-1}\\z_n \end{bmatrix} = \begin{bmatrix} \beta_0\\\beta_1\\\beta_2\\\vdots\\\beta_i\\\vdots\\\beta_{n-1}\\\beta_n \end{bmatrix} \]

  • \(z_n = 0\), \(z_k = \beta_k - h'_{k+1}z_{k+1}\), \(k = n-1, n-2, \dots, 0\).

Again, what are the inputs, what are the indices, give an example…

\(h'\)s, \(\beta_i\)s and \(z_i\)s again

Purpose of these is to solve the tridiagonal system of equations:

\[ \begin{bmatrix} 1 & 0 & 0 & & & & & \\ h_1 & 2\delta_1 & h_2 & & & & & \\ &h_2 & 2\delta_2 & h_3 & & & & \\ & & \ddots & \ddots & \ddots & & & \\ & & & h_{i} & 2\delta_i & h_{i+1} & & \\ & & & & \ddots & \ddots & \ddots & \\ & & & & & h_{n-1} & 2\delta_{n-1} & h_{n} \\ & & & & & 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} z_0\\z_1\\z_2\\\vdots\\z_i\\\vdots\\z_{n-1}\\z_n \end{bmatrix} = \begin{bmatrix} 0\\B_1\\B_2\\\vdots\\B_i\\\vdots\\B_{n-1}\\0 \end{bmatrix} \]

Put it all together

  • What are the inputs?

  • What are the outputs?

  • How do we evaluate the output?