前回からの続きです。
こちらを先に読んでからご覧ください。
スプライン補間について(その1)
3.導出計算
式の形と求める対象
求めるものを最初に確認しておきます。
いま $n+1$個すべてのデータ点を通るスプライン曲線を求めようとしているので、欲しいのは
S_0(x), S_1(x), \cdots, S_{n-1}(x)
という $n$個の多項式です。
いま各多項式はすべて3次の多項式としているので、結局求めるべきなのは $4 \times n$ 個の係数ということになります。
各多項式を以下のように設定します。
\begin{eqnarray}
S_0(x) &=& a_0 + b_0(x - x_0) + c_0(x - x_0)^2 + d_0(x - x_0)^3 \\
S_1(x) &=& a_1 + b_1(x - x_1) + c_1(x - x_1)^2 + d_1(x - x_1)^3 \\
S_2(x) &=& a_2 + b_2(x - x_2) + c_2(x - x_2)^2 + d_2(x - x_2)^3 \\
&\vdots& \\
S_j(x) &=& a_j + b_j(x - x_j) + c_j(x - x_j)^2 + d_j(x - x_j)^3 \\
&\vdots& \\
S_{n-1}(x) &=& a_{n-1} + b_{n-1}(x - x_{n-1}) + c_{n-1}(x - x_{n-1})^2 + d_{n-1}(x - x_{n-1})^3 \\
\end{eqnarray}
変数 $x$ から各区間の左側のデータ $x_j$ を引いているのは、計算をより簡単にするためです。
この形にせず、各多項式を $S_j(x)=a_j+b_jx+c_jx^2+d_jx^3$ として求めるとかなり大変になります。
興味のある方はぜひ挑戦してみてください。
これらの多項式の中で、
・$x$には求めたい点の$x$座標を任意に入れる。
・各$x_j$は既にデータとして持っている。
のでこれらは既知。
つまり結局求めなければならないのは
\begin{eqnarray}
&&a_0, b_0, c_0, d_0, \\
&&a_1, b_1, c_1, d_1, \\
&&\qquad \vdots \\
&&a_j, b_j, c_j, d_j, \\
&&\qquad \vdots \\
&&a_{n-1}, b_{n-1}, c_{n-1}, d_{n-1}
\end{eqnarray}
の各多項式の項の係数。
これが $4n$ 個ある、ということです。
未知変数が $4n$個なので、必要になる方程式の数も $4n$個になります。
では、先ほどまでに要請した4つの条件を使って式を立てていきましょう。
条件から導かれる式
条件1から
1つ目の条件はこれでした。
\begin{eqnarray}
S_0(x_1) &=& S_1(x_1) = y_1 \\
S_1(x_2) &=& S_2(x_2) = y_2 \\
S_2(x_3) &=& S_3(x_3) = y_3 \\
&\vdots& \\
S_j(x_{j+1}) &=& S_{j+1}(x_{j+1}) = y_{j+1}\\
\end{eqnarray}
これを$(x_j - x_j)=0$に注意して具体的に書くと、
\scriptsize{
\begin{eqnarray}
\begin{matrix}
S_0(x_0) & = & a_0 & + & 0 & + & 0 & + & 0 & = & y_0 \\
S_0(x_1) & = & a_0 & + & b_0(x_1 - x_0) & + & c_0(x_1 - x_0)^2 & + & d_0(x_1 - x_0)^3 & = & y_1 \\
S_1(x_1) & = & a_1 & + & 0 & + & 0 & + & 0 & = & y_1 \\
S_1(x_2) & = & a_1 & + & b_1(x_2 - x_1) & + & c_1(x_2 - x_1)^2 & + & d_1(x_2 - x_1)^3 & = & y_2 \\
& \vdots & \\
S_j(x_j) & = & a_j & + & 0 & + & 0 & + & 0 & = & y_j \\
S_j(x_{j+1}) & = & a_j & + & b_j(x_{j+1} - x_j) & + & c_j(x_{j+1} - x_j)^2 & + & d_j(x_{j+1} - x_j)^3 & = & y_{j+1} \\
& \vdots & \\
S_{n-1}(x_{n-1}) & = & a_{n-1} & + & 0 & + & 0 & + & 0 & = & y_{n-1} \\
S_{n-1}(x_{n}) & = & a_{n-1} & + & b_{n-1}(x_{n} - x_{n-1}) & + & c_{n-1}(x_{n} - x_{n-1})^2 & + & d_{n-1}(x_{n} - x_{n-1})^3 & = & y_{n} \\
\end{matrix}
\end{eqnarray}
}
となります。
$S_j(x)$1つにつき2つの式が出来ているので、これで $2n$個の式が得られました。
式が煩雑になるので、以下
$\Delta x_j = (x_{j+1} - x_j)$
と置き換えることにします。
【例】
$S_j(x_{j+1}) = a_j + b_j(x_{j+1} - x_j) + c_j(x_{j+1} - x_j)^2 + d_j(x_{j+1} - x_j)^3 = y_{j+1}$
↓ ↓ ↓
$S_j(x_{j+1}) = a_j + b_j {\Delta x_j} + c_j{\Delta x_j}^2 + d_j{\Delta x_j}^3 = y_{j+1} $
条件2から
2つ目の条件はこれでした。
\begin{eqnarray}
S_0'(x_1) &=& S_1'(x_1) \\
S_1'(x_2) &=& S_2'(x_2) \\
S_2'(x_3) &=& S_3'(x_3) \\
&\vdots& \\
S_j'(x_{j+1}) &=& S_{j+1}'(x_{j+1})
\end{eqnarray}
$S_j'(x) = b_j + 2c_j(x - x_j) + 3d_j(x - x_j)^2$ なので、これを具体的に書けば、
\begin{eqnarray}
\begin{matrix}
S_0'(x_1) & = & b_0 & + & 2c_0{\Delta x_0} & + & 3d_0{\Delta x_0}^2 & = & b_1 & = & S_1'(x_1) \\
S_1'(x_2) & = & b_1 & + & 2c_1{\Delta x_1} & + & 3d_1{\Delta x_1}^2 & = & b_2 & = & S_2'(x_2) \\
& & & & & \vdots \\
S_j'(x_{j+1}) & = & b_j & + & 2c_j{\Delta x_j} & + & 3d_j{\Delta x_j}^2 & = & b_{j+1} & = & S_{j+1}'(x_{j+1}) \\
& & & & & \vdots \\
S_{n-2}'(x_{n-1}) & = & b_{n-2} & + & 2c_{n-2}{\Delta x_{n-2}} & + & 3d_{n-2}{\Delta x_{n-2}}^2 & = & b_{n-1} & = & S_{n-1}'(x_{n-1}) \\
\end{matrix}
\end{eqnarray}
となります。
これらの式は $n-1$個あるので、さっきのものと併せて $3n-1$個になりました。
条件3から
3つ目の条件はこれでした。
\begin{eqnarray}
S_0''(x_1) &=& S_1''(x_1) \\
S_1''(x_2) &=& S_2''(x_2) \\
S_2''(x_3) &=& S_3''(x_3) \\
&\vdots& \\
S_j''(x_{j+1}) &=& S_{j+1}''(x_{j+1})
\end{eqnarray}
$S_j''(x) = 2c_j + 6d_j(x - x_j)$ なので、これを具体的に書けば、
\begin{eqnarray}
\begin{matrix}
S_0''(x_1) & = & 2c_0 & + & 6d_0{\Delta x_0} & = & 2c_1 & = & S_1''(x_1) \\
S_1''(x_2) & = & 2c_1 & + & 6d_1{\Delta x_1} & = & 2c_2 & = & S_2''(x_2) \\
&&&\vdots \\
S_j''(x_{j+1}) & = & 2c_j & + & 6d_j{\Delta x_j} & = & 2c_{j+1} & = & S_{j+1}''(x_{j+1}) \\
&&&\vdots \\
S_{n-2}''(x_{n-1}) & = & 2c_{n-2} & + & 6d_{n-2}{\Delta x_{n-2}} & = & 2c_{n-1} & = & S_{n-1}''(x_{n-1}) \\
\end{matrix}
\end{eqnarray}
となります。
これらの式は2で割れるので、あとで両辺を2で割ったものを使います。
これらの式も $n-1$個あるので、いままでと併せて $4n-2$個になりました。
あと2つです。
条件4から
4つ目の条件はこれでした。
\begin{eqnarray}
S_0''(x_0) &=& 0 \\
S_{n-1}''(x_n) &=& 0
\end{eqnarray}
これを具体的に書くと、
\begin{eqnarray}
\begin{matrix}
S_0''(x_1) & = & 2c_0 & & & = & 0 \\
S_{n-1}''(x_n) & = & 2c_{n-1} & + & 6d_{n-1}{\Delta x_{n-1}} & = & 0 \\
\end{matrix}
\end{eqnarray}
です。
これで残りの2個が出たので、全部で $4n$個の連立する方程式が用意できました。
求解方法
要請した4つの条件から必要な$4n$個の式が用意できました。
ではこれらの式を連立させて$4n$個の係数の値を求めていきます。
連立方程式を解く方法は、
・式変形と代入をして変数を一つづつ消していく方法
・行列計算を用いる方法
があります。
式変形と代入を繰り返す方法はかなり大変なので、行列計算を用いる方法を使うことにします。
いま求めたいのは係数です。
なので、係数を列ベクトルにすると、
A =
\begin{bmatrix}
a_0 \\ b_0 \\ c_0 \\ d_0 \\ a_1 \\ b_1 \\ \vdots \\ c_{n-1} \\ d_{n-1}
\end{bmatrix}
これに合わせてこれまでの条件から出た式を適用するためには、
\begin{bmatrix}
1 & 0 & 0 & 0 & 0 & 0 & & \cdots & & & \\
1 & \Delta x_0 & {\Delta x_0}^2 & {\Delta x_0}^3 & 0 & \cdots & \\
0 & 0 & 0 & 0 & 1 & 0 & \cdots & \\
0 & 0 & 0 & 0 & 1 & \Delta x_1 & {\Delta x_1}^2 & {\Delta x_1}^3 & 0 & \cdots & \\
&&& \vdots &&&
\end{bmatrix}
と係数の順番に合わせて、それぞれの条件式を不要な項に0をかけたりして合う形に調整した行列を作ってやればいいです。
具体的には、以下のようにすべての係数を足し合わせた形を考えます。
\Box a_0 + \Box b_0 + \Box c_0 + \Box d_0 + \Box a_1 + \Box b_1 + \cdots + \Box c_{n-1} + \Box d_{n-1} = \Box
そして、この$\Box$のところに、各条件式の係数にかかっている値をいれてやります。
例えば、条件4の2つ目
c_{n-1} + 3d_{n-1}{\Delta x_{n-1}} = 0
なら、
[0] \cdot a_0 + [0] \cdot b_0 + \cdots + [0] \cdot b_{n-1} + [1] \cdot c_{n-1} + [3{\Delta x_{n-1}}] \cdot d_{n-1} = [0]
といった形です。
各式において、出てこない係数のところには0を掛けて消してしまうわけです。
そして、これを係数ベクトルとその係数にかかっていた値との行列に分けてやると、先ほどの行列が得られます。
残念ながらちゃんと書くと大きくなりすぎてここに表示できません。
ちょっと理解しずらいと思うので、並べ方を理解してもらう為に3区間の場合を例として書いておきます。
【例】 $(x_0, y_0), (x_1, y_1), (x_2, y_2), (x_3, y_3)$ の4点を結ぶスプライン曲線を求める場合
区間多項式は $S_0(x), S_1(x), S_2(x)$ の3つ。
よって、これらに条件を適用した式からできる行列は、
\scriptstyle{
X=
\begin{bmatrix}
1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
1 & \Delta x_0 & {\Delta x_0}^2 & {\Delta x_0}^3 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 &1 & \Delta x_1 & {\Delta x_1}^2 & {\Delta x_1}^3 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & \Delta x_2 & {\Delta x_2}^2 & {\Delta x_2}^3 \\
0 & 1 & 2{\Delta x_0} & 3{\Delta x_0}^2 & 0 & -1 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 1 & 3{\Delta x_0} & 0 & 0 & -1 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 1 & 2{\Delta x_1} & 3{\Delta x_1}^2 & 0 & -1 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 1 & 3{\Delta x_1} & 0 & 0 & -1 & 0 \\
0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 3{\Delta x_2} \\
\end{bmatrix}
}
という行列になる。
例えば2行目は
\small{1 \cdot a_0 + {\Delta x_0} \cdot b_0 + {\Delta x_0}^2 \cdot c_0 + {\Delta x_0}^3 \cdot d_0 + 0 \cdot a_1 + 0 \cdot b_1 + \cdots + 0 \cdot c_2 + 0 \cdot d_2 = y_1
}
7行目は $b_1$ を移項して
\small{0 \cdot a_0 + 1 \cdot b_0 + 2{\Delta x_0} \cdot c_0 + 3{\Delta x_0}^2 \cdot d_0 + 0 \cdot a_1 + (-1) \cdot b_1 + \cdots + 0 \cdot c_2 + 0 \cdot d_2 = 0
}
という各条件式から作られている。
尚、これに対応する $y$ のベクトルは
Y =
\begin{bmatrix}
y_0 \\ y_1 \\ y_1 \\ y_2 \\ y_2 \\ y_3 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0
\end{bmatrix}
になる。
つまり、3区間の場合、条件式をまとめた行列の式は以下になる。
\scriptstyle{
\begin{bmatrix}
1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
1 & \Delta x_0 & {\Delta x_0}^2 & {\Delta x_0}^3 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 &1 & \Delta x_1 & {\Delta x_1}^2 & {\Delta x_1}^3 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & \Delta x_2 & {\Delta x_2}^2 & {\Delta x_2}^3 \\
0 & 1 & 2{\Delta x_0} & 3{\Delta x_0}^2 & 0 & -1 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 1 & 3{\Delta x_0} & 0 & 0 & -1 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 1 & 2{\Delta x_1} & 3{\Delta x_1}^2 & 0 & -1 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 1 & 3{\Delta x_1} & 0 & 0 & -1 & 0 \\
0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 3{\Delta x_2} \\
\end{bmatrix}
\begin{bmatrix}
a_0 \\ b_0 \\ c_0 \\ d_0 \\ a_1 \\ b_1 \\ c_1 \\ d_1 \\ a_2 \\ b_2 \\ c_2 \\ d_2
\end{bmatrix}
=
\begin{bmatrix}
y_0 \\ y_1 \\ y_1 \\ y_2 \\ y_2 \\ y_3 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0
\end{bmatrix}
}
行列とベクトルを大文字で表して、
XA=Y
こうしてできた行列を $X$とし、対応する $y$ のベクトルを $Y$とすると、
XA = Y
となります。
問題はこの行列$X$が正則かどうかですが、$\Delta x_j$ が0でない限りどの行も一次従属にならないことがいろいろ計算すると分かります。
(長いので割愛します、すみません)
なので、この$X$は正則です。
よって、この行列$X$には逆行列が存在して、
A = X^{-1}Y
と、$X^{-1}$を求めることで計算できます。
これで全部の係数を求めることができました。
あとは、この係数をそれぞれの$S_j(x)$に入れて、間の値を計算して補間するだけです。
A = X^{-1}Y = \begin{bmatrix} a_0 & b_0 & \cdots & a_j & b_j & c_j & d_j & \cdots & c_{n-1} & d_{n-1} \end{bmatrix}^{t}
と$A$が解かれるので、該当場所の係数を抜き出して各$S_j(x)$の係数に入れて利用する。
S_j(x) = a_j + b_j(x - x_j) + c_j(x - x_j)^2 + d_j(x - x_j)^3
実際の計算はプログラム等を使って行列の数値計算で解けばいいですね。
Excelでも計算できます。
別解
もう一つの求め方も書いておきます。
各係数を移項と代入で求めていきます。
ただ、結局は複雑な連立方程式を解かなければいけないので、私はあまりお勧めしません。
まず、条件1から以下が分かります。
a_j = y_j
多項式の $x$ を $(x - x_j)$ に変えたことにより、$x=x_j$を代入したときに定数項 $a_j$だけが残ります。
この時の $y$の値は $y_j$ですから、このことが分かります。
次に条件3から、
c_j + 3d_j{\Delta x_j} = c_{j+1} \\
\therefore d_j = \frac{c_{j+1} - c_j}{3{\Delta x_j}}
が出ます。
さらに条件1から、この$d_j$を使って、
a_j + b_j{\Delta x_j} + c_j{\Delta x_j}^2 + d_j{\Delta x_j}^3 = y_{j+1} \\
y_j + b_j{\Delta x_j} + c_j{\Delta x_j}^2 + d_j{\Delta x_j}^3 = y_{j+1} \\
b_j{\Delta x_j} = y_{j+1} - y_j - c_j{\Delta x_j}^2 - d_j{\Delta x_j}^3 \\
b_j{\Delta x_j} = y_{j+1} - y_j - c_j{\Delta x_j}^2 - \frac{c_{j+1} - c_j}{3{\Delta x_j}}{\Delta x_j}^3 \\
b_j{\Delta x_j} = y_{j+1} - y_j - c_j{\Delta x_j}^2 - \frac{c_{j+1} - c_j}{3}{\Delta x_j}^2 \\
b_j{\Delta x_j} = y_{j+1} - y_j - \frac{c_{j+1} + 2c_j}{3}{\Delta x_j}^2 \\
\\
\therefore b_j = \frac{y_{j+1} - y_j}{\Delta x_j} - \frac{{\Delta x_j}(c_{j+1} + 2c_j)}{3}
ここまでで、$a_j, b_j, d_j$ の3つを $c_j$ を使って表すことが出来ました。
\begin{eqnarray}
a_j &=& y_j \\
b_j &=& \frac{y_{j+1} - y_j}{\Delta x_j} - \frac{{\Delta x_j}(c_{j+1} + 2c_j)}{3} \\
d_j &=& \frac{c_{j+1} - c_j}{3{\Delta x_j}}
\end{eqnarray}
残りの $c_j$ですが、残念ながらこのように綺麗になりません。
条件2から、
\small{
\begin{eqnarray}
b_j + 2c_j{\Delta x_j} + 3d_j{\Delta x_j}^2 & = & b_{j+1} \\
\frac{y_{j+1} - y_j}{\Delta x_j} - \frac{{\Delta x_j}(c_{j+1} + 2c_j)}{3} + 2c_j{\Delta x_j} + 3 \cdot \frac{c_{j+1} - c_j}{3{\Delta x_j}}{\Delta x_j}^2 &=& \frac{y_{j+2} - y_{j+1}}{\Delta x_{j+1}} - \frac{{\Delta x_{j+1}}(c_{j+2} + 2c_{j+1})}{3}
\end{eqnarray}
}
これを頑張って整理すると、
{\Delta x_j}c_j + 2(\Delta x_{j+1} + \Delta x_j)c_{j+1} + {\Delta x_{j+1}}c_{j+2} = \frac{3(y_{j+2} - y_{j+1})}{\Delta x_{j+1}} - \frac{3(y_{j+1} - y_{j})}{\Delta x_{j}}
という、連続する $c_j, c_{j+1}, c_{j+2}$3項の逐次式が得られます。
漸化式っぽいですが、項が進むと$\Delta x_j$も進んで係数が変わるので漸化式ではありません。
これを解く為には、項を一つずつ進めた複数の式を並べて
\small{
\begin{bmatrix}
1 & 0 & 0 & 0 & \cdots & 0 & 0 & 0 \\
{\Delta x_0} & 2(\Delta x_{1} + \Delta x_0) & {\Delta x_1} & 0 & \cdots & 0 & 0 & 0 \\
0 & {\Delta x_1} & 2(\Delta x_{2} + \Delta x_1) & {\Delta x_2} & \cdots & 0 & 0 & 0 \\
\vdots & \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \vdots \\
0 & 0 & 0 & 0 & \cdots & {\Delta x_{n-2}} & 2(\Delta x_{n-1} + \Delta x_{n-2}) & {\Delta x_{n-1}} \\
0 & 0 & 0 & 0 & \cdots & 0 & 0 & 1 \\
\end{bmatrix}
\begin{bmatrix}
c_0 \\ c_1 \\ c_2 \\ \vdots \\ c_{n-2} \\ c_{n-1} \\ 0
\end{bmatrix}
=
\begin{bmatrix}
0 \\ \frac{3(y_{2} - y_{1})}{\Delta x_{1}} - \frac{3(y_{1} - y_{0})}{\Delta x_{0}} \\
\frac{3(y_{3} - y_{2})}{\Delta x_{2}} - \frac{3(y_{2} - y_{1})}{\Delta x_{1}} \\
\vdots \\
\frac{3(y_{n} - y_{n-1})}{\Delta x_{n-1}} - \frac{3(y_{n-1} - y_{n-2})}{\Delta x_{n-2}} \\
0
\end{bmatrix}
}
という行列の式をつくって解くという方法がありますが、これなら最初から先ほどの行列での解法を用いた方が早いです。
この行列が「三重対角行列」という特別な形なので、『TDMA法(Thomas’algorithm)』という方法で数値計算がしやすい、というのが他でこの方法の紹介が多い理由かとは思いますが、このTDMA法を使える状況でない限りこちらの方法を使うメリットはないかと思います。
あらかじめ適用する区間数を決めてしまってこの行列のサイズを決定し、掃き出し法(ガウスの消去法)などを用いて各$c_j$の具体的な式の形を求めておく、ということも考えられます。
しかし、式がかなり複雑になりとても大変です。
計算が好きな方は挑戦してみるといいかもしれません。
尚、最初の2区間、3区間において$c_1, c_2$を求めておいて、その値を初期値として用いて逐次計算で$c_j$を求めていく、という方法を思いつくかもしれません。
しかし残念ながらこの方法では正しい係数が求まっていきません。
これらは連立方程式なので、前からだけ、後ろからだけでは値が決まらないからです。
実際にそのような計算を行ってみると、値が発散していってしまいます。
したがって、逐次手法としてこの手法をとることはできません。
2区間(連続3点)でスプライン曲線を求めた場合の$c$の式を記載しておきます。
興味のある方は、上のような逐次計算を行うとどう値が発散していくか確認してみてください。
\begin{eqnarray}
c_0 &=& 0 \\
c_1 &=& \frac{1}{2(\Delta x_1 + \Delta x_0)} \left \{ \frac{3(y_2 - y_1)}{\Delta x_1} - \frac{3(y_1 - y_0)}{\Delta x_0} \right \}
\end{eqnarray}
$a_j,b_j,d_j$は既に示してある通りです。
式だけ取っていかれる方がいるかもしれないので、再度書いておきます。
{\Delta x_j}c_j + 2(\Delta x_{j+1} + \Delta x_j)c_{j+1} + {\Delta x_{j+1}}c_{j+2} = \frac{3(y_{j+2} - y_{j+1})}{\Delta x_{j+1}} - \frac{3(y_{j+1} - y_{j})}{\Delta x_{j}}
と
c_1 = \frac{1}{2(\Delta x_1 + \Delta x_0)} \left \{ \frac{3(y_2 - y_1)}{\Delta x_1} - \frac{3(y_1 - y_0)}{\Delta x_0} \right \}
は、同時使用できません!!!
導出計算のまとめ
解法を2つ紹介しました。
・条件式をまとめて行列の式にしてしまい、逆行列を使って解く
・各係数の具体的な形を代入と式変形を繰り返して求めて解く
2つ目は$c_j$を求めるのが大変。
行列計算ができるのであれば、1つ目の方法がおすすめ。
つづく
その3に続きます。