1. 来源和背景

对于一个（主）三对角方程组，我们常用“追赶法”来进行求解. 而三对角方程组常常出现于微分方程的数值求解，例如热传导方程的边值问题

{y′′(x)=f(x,y,y′), a≤x≤by(a)=η1, y(b)=η2

\begin{cases} y^{''}(x)=f(x,y,y^{'}),\ a\leq x\leq b\\ y(a)=\eta_1,\ y(b)=\eta_2 \end{cases}
当 f(x,y,y′)f(x,y,y^{'})是一个线性函数时，对该边值问题的数值解转化为一个典型的三对角方程组求解.

“追赶法”目前比较可靠的来源是下面的文章：
Thomas, L.H., Elliptic Problems in Linear Differential Equations over a Network. Watson Science Computer Laboratory Report, 1949.
其中的一个依据是，在国外的文章和教材中，“追赶法”被称为“Thomas算法”.

2. 追赶法的基本原理

追赶法的基本原理是矩阵的LU分解，即将矩阵AA分解为

A=LU

A=LU
其中，LL为一个对角线上元素为11的下三角矩阵，UU为一个上三角矩阵. 容易验证，一个三对角矩阵作LU分解以后，得到一个下二对角矩阵与一个上二对角矩阵的乘积，即

A=⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢a11a21a12a22a32a23a33⋱a34⋱⋱an−1,n−2an−1,n−1an,n−1an−1,nan−1,n⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥

A=\left[\begin{array}{cccccc} a_{11} & a_{12} & & & & & \\ a_{21} & a_{22} & a_{23} & & & & \\ & a_{32} & a_{33} & a_{34} & & \\ & & & & & &\\ & & \ddots & \ddots & \ddots & & \\ & & & & & & \\ & & & & a_{n-1,n-2} & a_{n-1,n-1} & a_{n-1,n} \\ & & & & & a_{n,n-1} & a_{n-1,n} \\\end{array}\right]

L=⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢1ℓ211ℓ321⋱⋱ℓn−1,n−21ℓn,n−11⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥

L=\left[\begin{array}{cccccc} 1 & & & & & & \\ \ell_{21} & 1 & & & & & \\ & \ell_{32} & 1 & & & \\ & & & & & &\\ & & \ddots & \ddots & & & \\ & & & & & & \\ & & & & \ell_{n-1,n-2} & 1 & \\ & & & & & \ell_{n,n-1} & 1 \\\end{array}\right]

U=⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢u11u12u22u23u33u34⋱⋱un−1,n−1un−1,nun−1,n⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥

U=\left[\begin{array}{cccccc} u_{11} & u_{12} & & & & & \\ & u_{22} & u_{23} & & & & \\ & & u_{33} & u_{34} & & \\ & & & & & &\\ & & & \ddots & \ddots & & \\ & & & & & & \\ & & & & & u_{n-1,n-1} & u_{n-1,n} \\ & & & & & & u_{n-1,n} \\\end{array}\right]

三对角矩阵AA的LU分解计算过程如下：

for i = 2 to nA(i,i-1) = A(i,i-1)/A(i-1,i-1);A(i,i) = A(i,i) - A(i-1,i) * A(i,i-1);
end

在计算过程中，将下三角矩阵LL和上三角矩阵UU的值保存在原矩阵AA中. 计算结束以后，矩阵AA中的元素为

⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢u11ℓ21u12u22ℓ32u23u33⋱u34⋱⋱ℓn−1,n−2un−1,n−1ℓn,n−1un−1,nun−1,n⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥

\left[\begin{array}{cccccc} u_{11} & u_{12} & & & & & \\ \ell_{21} & u_{22} & u_{23} & & & & \\ & \ell_{32} & u_{33} & u_{34} & & \\ & & & & & &\\ & & \ddots & \ddots & \ddots & & \\ & & & & & & \\ & & & & \ell_{n-1,n-2} & u_{n-1,n-1} & u_{n-1,n} \\ & & & & & \ell_{n,n-1} & u_{n-1,n} \\\end{array}\right]

注：三对角矩阵AA做LU分解以后，严格上三角部分的元素没有发生变化，即上三角矩阵UU中的元素

ui,i+1=ai,i+1, i=1,2,…,n−1

u_{i,i+1}=a_{i,i+1},\ i=1,2,\dots,n-1

3. 追赶法求解三对角方程组

使用LU分解的求解线性方程组时，不需要存储下三角矩阵，而上三角矩阵将被用于回代求解.

3.1 “追”的过程：分解

对于nn阶的三对角方程组

Ax=b

Ax=b
我们先用LU分解得到

Ux=y

注：由Ax=LUx=bAx=LUx=b，得

Ux=L−1b

Ux=L^{-1}b
记y=L−1by=L^{-1}b，即得到方程组Ux=yUx=y.

计算过程如下：

for i = 2 to nA(i,i-1) = A(i,i-1)/A(i-1,i-1);A(i,i) = A(i,i) - A(i-1,i) * A(i,i-1);b(i) = b(i) - b(i-1) * A(i,i-1);
end

循环里面的前两行与LU分解完全相同，第三行负责对常数项做相应的变换. 在计算过程中，上三角矩阵UU的值保存在原矩阵AA中，变换后的常数y=L−1by=L^{-1}b保存在bb中.

3.2 “赶”的过程：回代

接着，我们用回代法求解上三角形方程组. 从三对角矩阵得到的上三角形方程组如下：

⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢u11u12u22u23u33u34⋱⋱un−1,n−1un−1,nun−1,n⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎛⎝⎜⎜⎜⎜⎜⎜⎜x1x2⋮xn−1xn⎞⎠⎟⎟⎟⎟⎟⎟⎟=⎛⎝⎜⎜⎜⎜⎜⎜⎜y1y2⋮yn−1yn⎞⎠⎟⎟⎟⎟⎟⎟⎟

\left[\begin{array}{cccccc} u_{11} & u_{12} & & & & & \\ & u_{22} & u_{23} & & & & \\ & & u_{33} & u_{34} & & \\ & & & & & &\\ & & & \ddots & \ddots & & \\ & & & & & & \\ & & & & & u_{n-1,n-1} & u_{n-1,n} \\ & & & & & & u_{n-1,n} \\\end{array}\right] \left(\begin{array}{c} x_1\\ x_2\\ \vdots\\ x_{n-1} \\ x_n \end{array}\right)= \left(\begin{array}{c} y_1\\ y_2\\ \vdots\\ y_{n-1} \\ y_n \end{array}\right)
注意在前面的计算过程中，我们将上三角矩阵UU保存在AA中，常数项yy保存在bb中. 因此，我们得到如下的回代过程：

x(n) = b(n) / A(i,i);
for i = n-1 to 1x(i) = (b(i) - A(i,i+1) * x(i+1)) / A(i,i);
end

4. 实用的程序代码

在三对角矩阵中，三对角线以外的元素均为00，为了提高存储的效率，我们只需存储三对角线上的元素即可. 因此，对于前面的矩阵AA，我们只存储三个向量：

d=[A(1,1),A(2,2),...,A(n,n)];
u=[A(1,2),A(2,3),...,A(n-1,n)];
l=[A(2,1),A(3,2),...,A(n,n-1)];

这三个向量分别为矩阵A<script type="math/tex" id="MathJax-Element-38">A</script>三条对角线上的元素. 假定常数向量为

b=[b(1),b(2),...,b(n)];

则实用的追赶法（亦称为“Thomas算法”）求解三对角方程组的过程如下：

% 追
for i = 2 to nl(i-1) = l(i-1)/d(i-1);d(i,i) = d(i,i) - u(i-1) * l(i-1);b(i) = b(i) - b(i-1) * l(i-1);
end
% 赶
x(n) = b(n) / d(i);
for i = n-1 to 1x(i) = (b(i) - u(i) * x(i+1)) / d(i);
end

追赶法求解三对角方程组相关推荐

用追赶法求方程组c语言程序,数值计算——追赶法求解三对角方程组（附代码）...
目录追赶法基础理论追赶法c++程序代码程序运行结果源码文件下载地址追赶法基础理论在数值计算中,对三次样条曲线插值和用差分方法求解常微分方程边值问题时,通常会遇到Ax=d三对角形式的方程组: ...
【数学建模之Python】12.追赶法求解三对角方程组
如果解决了你的问题,点个赞再走嘛٩(๑❛ᴗ❛๑)۶ 目录一.前言二.方程组的特点 1.类型 2.要求 3.方法的优势三.算法步骤四.例题+代码 1.用追赶法求解以下五阶方程组 2.代码一.前 ...
Python解线性方程组的直接法(6)————求解三对角方程组的追赶法
求解三对角方程组的追赶法 import numpy as npdef zuiganfa(A, d):n = A.shape[0]l = np.mat(np.zeros(n, dtype=float
追赶法求解三对角线性方程组的MATLAB程序
function [M] = after_method(a,b,c,g,tolerance) %% 追赶法求解三对角矩阵 % a为三对角矩阵左下对角线上的值 % b为三对角矩阵中间对角线上的值 % c ...
c语言程序用追赶法求解方程组,编写用追赶法解三对角线性方程组的程序,并解下列方程组（3页）-原创力文档...
计算方法与实习上机实验(二) 实验名称: 编写用追赶法解三对角线性方程组的程序,并解下列方程组: (1) (2)Ax=b,其中 A10×10=, b10×1= 程序代码: #include using ...
c语言程序用追赶法求解方程组,编写用追赶法解三对角线性方程组的程序,并解下列方程组...
计算方法与实习上机实验(二) 实验名称: 编写用追赶法解三对角线性方程组的程序,并解下列方程组: 2x1 x2 5, x 2x x 12, 123(1) x 2x x 11,34 2 x3 2x4 1 ...
Matlab实现——严格对角占优三对角方程组求解（高斯赛尔德Gauss-Seidel迭代、超松弛）
欢迎前往个人博客驽马点滴和视频空间哔哩哔哩-<挨踢日志> 严格对角占优三对角方程组求解对中等规模的n阶的(n<100)线性方程组,直接法的准确性和可靠性,所以常采用直接法对 ...
matlab trisys,Matlab实现——严格对角占优三对角方程组求解（高斯赛尔德Gauss-Seidel迭代、超松弛） | 学步园...
严格对角占优三对角方程组求解对中等规模的n阶的(n<100)线性方程组,直接法的准确性和可靠性,所以常采用直接法对于较高阶的方程组,特别是地于某些偏微分方程离散化后得到的大型稀疏方程组(系统 ...
追赶法求解方程组备忘
本篇内容为数值分析中,用追赶法求解方程组的方法,备忘如下: 1. 原理部分追赶法求解的矩阵格式一般如下: a1 c1 0 0 b2 a2 c2 0 0 b3 a3 c3 0 0 b4 a4 如果矩阵 ...

追赶法求解三对角方程组