28-2 Splines

A pratical method for interpolating a set of points with a curve is to use cubic splines. We are given a set \(\\{(x_i, y_i): i = 0, 1, \ldots, n\\}\) of \(n + 1\) point-value pairs, where \(x_0 < x_1 < \cdots < x_n\). We wish to fit a piecewise-cubic curve (spline) \(f(x)\) to the points. That is, the curve \(f(x)\) is made up of \(n\) cubic polynomials \(f_i(x) = a_i + b_ix + c_ix^2 + d_ix^3\) for \(i = 0, 1, \ldots, n - 1\), where if \(x\) falls in the range \(x_i \le x \le x_{i + 1}\), then the value of the curve is given by \(f(x) = f_i(x - x_i)\). The points \(x_i\) at which the cubic polynomials are "pasted" together are called knots. For simplicity, we shall assume that \(x_i = i\) for \(i = 0, 1, \ldots, n\).

To ensure continuity of \(f(x)\), we require that

\[ \begin{aligned} f(x_i) & = f_i(0) = y_i, \\\\ f(x_{i + 1}) & = f_i(1) = y_{i + 1} \end{aligned} \]

for \(i = 0, 1, \ldots, n - 1\). To ensure that \(f(x)\) is sufficiently smooth, we also insist that the first derivative be continuous at each knot:

\[f'(x_{i + 1}) = f'\_i(1) = f'\_{i + 1}(0)\]

for \(i = 0, 1, \ldots, n - 2\).

a. Suppose that for \(i = 0, 1, \ldots, n\), we are given not only the point-value pairs \(\\{(x_i, y_i)\\}\) but also the first derivatives \(D_i = f'(x_i)\) at each knot. Express each coefficient \(a_i\), \(b_i\), \(c_i\) and \(d_i\) in terms of the values \(y_i\), \(y_{i + 1}\), \(D_i\), and \(D_{i + 1}\). (Remember that \(x_i = i\).) How quickly can we compute the \(4n\) coefficients from the point-value pairs and first derivatives?

The question remains of how to choose the first derivatives of \(f(x)\) at the knots. One method is to require the second derivatives to be continuous at the knots:

\[f''(x_{i + 1}) = f''\_i(1) = f''\_{i + 1}(0)\]

for \(i = 0, 1, \ldots, n - 2\). At the first and last knots, we assume that \(f''(x_0) = f''\_0(0) = 0\) and \(f''(x_n) = f''_{n - 1}(1) = 0\); these assumptions make \(f(x)\) a natural cubic spline.

b. Use the continuity constraints on the second derivative to show that for \(i = 1, 2, \ldots, n - 1\),

\[D_{i - 1} + 4D_i + D_{i + 1} = 3(y_{i + 1} - y_{i - 1}). \tag{23.21}\]

c. Show that

\[ \begin{aligned} 2D_0 + D_1 & = 3(y_1 - y_0), & \text{(28.22)} \\\\ D_{n - 1} + 2D_n & = 3(y_n - y_{n - 1}). & \text{(28.23)} \end{aligned} \]

d. Rewrite equations \(\text{(28.21)}\)–\(\text{(28.23)}\) as a matrix equation involving the vector \(D = \langle D_0, D_1, \ldots, D_n \rangle\) or unknowns. What attributes does the matrix in your equation have?

e. Argue that a natural cubic spline can interpolate a set of \(n + 1\) point-value pairs in \(O(n)\) time (see Problem 28-1).

f. Show how to determine a natural cubic spline that interpolates a set of \(n + 1\) points \((x_i, y_i)\) satisfying \(x_0 < x_1 < \cdots < x_n\), even when \(x_i\) is not necessarily equal to \(i\). What matrix equation must your method solve, and how quickly does your algorithm run?

a. We have \(a_i = f_i(0) = y_i\) and \(b_i = f_i'(0) = f'(x_i) = D_i\). Since \(f_i(1) = a_i + b_i + c_i + d_i\) and \(f_i'(1) = b_i + 2c_i + 3d_i\), we have \(d_i = D_{i + 1} - 2y_{i + 1} + 2y_i + D_i\) which implies \(c_i = 3y_{i + 1} - 3y_i - D_{i + 1} - 2D_i\). Since each coefficient can be computed in constant time from the known values, we can compute the \(4n\) coefficients in linear time.

b. By the continuity constraints, we have \(f_i''(1) = f_{i + 1}''(0)\) which implies that \(2c_i + 6d_i = 2c_{i + 1}\), or \(c_i + 3d_i = c_{i + 1}\). Using our equations from above, this is equivalent to

Rearranging gives the desired equation

c. The condition on the left endpoint tells us that \(f_0''(0) = 0\), which implies \(2c_0 = 0\). By part (a), this means \(3(y_1 − y_0) = 2D_0 + D_1\). The condition on the right endpoint tells us that \(f_{n - 1}''(1) = 0\), which implies \(c_{n - 1} + 3d_{n - 1} = 0\). By part (a), this means \(3(y_n - y_{n - 1}) = D_{n - 1} + 2D_n\).

d. The matrix equation has the form \(AD = Y\), where \(A\) is symmetric and tridiagonal. It looks like this:

e. Since the matrix is symmetric and tridiagonal, Problem 28-1 (e) tells us that we can solve the equation in \(O(n)\) time by performing an LUP decomposition. By part (a), once we know each \(D_i\) we can compute each \(f_i\) in \(O(n)\) time.

f. For the general case of solving the nonuniform natural cubic spline problem, we require that \(f(x_{i + 1}) = f_i(x_{i + 1} − x_i) = y_{i + 1}\), \(f'(x_{i + 1}) = f_i'(x_{i + 1} - x_i) = f_{i + 1}'(0)\) and \(f''(x_{i + 1}) = f_i''(x_{i + 1} - x_i) = f_{i + 1}''(0)\). We can still solve for each of \(a_i\), \(b_i\), \(c_i\) and \(d_i\) in terms of \(y_i\), \(y_{i + 1}\), \(D_i\) and \(D_{i + 1}\), so we still get a tridiagonal matrix equation. The solution will be slightly messier, but ultimately it is solved just like the simpler case, in \(O(n)\) time.

本页面的全部内容在 小熊老师 - 莆田青少年编程俱乐部 0594codes.cn 协议之条款下提供，附加条款亦可能应用