一、三对角方程
形如:
A
=
[
b
1
c
1
a
2
b
2
c
2
a
3
b
3
⋱
⋱
⋱
c
n
−
1
a
n
b
n
]
A= \begin{bmatrix} b_1&c_1& \\ a_2&b_2&c_2\\ &a_3&b_3 &\ddots\\ &&\ddots&\ddots&c_{n-1} \\ &&&a_n&b_n \end{bmatrix}
A=⎣⎢⎢⎢⎢⎡b1a2c1b2a3c2b3⋱⋱⋱ancn−1bn⎦⎥⎥⎥⎥⎤
的矩阵
A
A
A被称为三对角矩阵。
形如
A
x
=
f
\boldsymbol {Ax=f}
Ax=f
的方程被称为三对角方程。其中
x
,
f
x,f
x,f的形式如下:
x
=
[
x
1
x
2
⋮
x
n
]
,
f
=
[
f
1
f
2
⋮
f
n
]
x= \begin{bmatrix} x_1 \\ x_2 \\ \vdots \\ x_n \end{bmatrix}, f= \begin{bmatrix} f_1 \\ f_2 \\ \vdots \\ f_n \end{bmatrix}
x=⎣⎢⎢⎢⎡x1x2⋮xn⎦⎥⎥⎥⎤,f=⎣⎢⎢⎢⎡f1f2⋮fn⎦⎥⎥⎥⎤
二、求方程的解 x x x
1. Crout分解
假设三对角矩阵
A
A
A有分解,则:
A
=
L
U
\boldsymbol {A=LU}
A=LU
其中, L L L为下三角矩阵, U U U为单位上三角矩阵(对角线元素为1)。
假设
L
,
U
L,U
L,U的形式如下:
L
=
[
α
1
γ
2
α
2
γ
3
α
3
⋱
⋱
γ
n
α
n
]
,
U
=
[
1
β
1
1
β
2
⋱
⋱
1
β
n
]
L = \begin{bmatrix} \alpha_1&&& \\ \gamma_2&\alpha_2 \\ &\gamma_3&\alpha_3 \\ &&\ddots&\ddots \\ &&&\gamma_n&\alpha_n \end{bmatrix}, U = \begin{bmatrix} 1&\beta_1\\ &1&\beta_2 \\ &&\ddots&\ddots \\ &&&1&\beta_n \end{bmatrix}
L=⎣⎢⎢⎢⎢⎡α1γ2α2γ3α3⋱⋱γnαn⎦⎥⎥⎥⎥⎤,U=⎣⎢⎢⎡1β11β2⋱⋱1βn⎦⎥⎥⎤
根据矩阵的乘法,对比 A A A中的元素,可以求出 α i , β i , γ i \alpha_i,\beta_i,\gamma_i αi,βi,γi中的所有元素,从而得到矩阵 L , U L,U L,U。
这里直接给出结果:
{
α
1
=
b
1
,
初
值
β
i
=
c
i
α
i
,
i
=
1
,
2
,
3
,
.
.
.
,
n
−
1
α
i
+
1
=
b
i
+
1
−
γ
i
+
1
β
i
,
i
=
1
,
2
,
3
,
.
.
.
,
n
−
1
\begin{cases} \alpha_1=b_1,&初值 \\ \beta_i=\frac {c_i}{\alpha_i},&i=1,2,3,...,n-1 \\ \alpha_{i+1}=b_{i+1}-\gamma_{i+1}\beta_i,& i=1,2,3,...,n-1 \end{cases}
⎩⎪⎨⎪⎧α1=b1,βi=αici,αi+1=bi+1−γi+1βi,初值i=1,2,3,...,n−1i=1,2,3,...,n−1
2. 追赶法求解
将
A
=
L
U
A=LU
A=LU带入三对角方程可得到:
L
U
x
=
f
\boldsymbol {LUx=f}
LUx=f
令 U x = y \boldsymbol {Ux=y} Ux=y,则 L y = f \boldsymbol {Ly=f} Ly=f。追赶法的做法是先求解 y \boldsymbol y y,再通过 y \boldsymbol y y求解 x \boldsymbol x x。
同样,求解也是通过矩阵的乘法。结果如下:
y
1
=
f
1
α
1
,
y
i
+
1
=
f
i
+
1
−
γ
i
+
1
y
i
α
i
+
1
,
i
=
0
,
1
,
2
,
.
.
.
,
n
−
1
x
n
=
y
n
,
x
i
=
y
i
−
β
i
x
i
+
1
,
i
=
n
−
1
,
n
−
2
,
.
.
.
,
1
y_1=\frac {f_1}{\alpha_1}, y_{i+1}=\frac {f_{i+1}-\gamma_{i+1}y_i}{\alpha_{i+1}},i=0,1,2,...,n-1 \\ x_n=y_n,x_i=y_i-\beta_ix_{i+1},i=n-1,n-2,...,1
y1=α1f1,yi+1=αi+1fi+1−γi+1yi,i=0,1,2,...,n−1xn=yn,xi=yi−βixi+1,i=n−1,n−2,...,1