这次来实现三对角线方程组的追赶法，追赶法的本质还是高斯消元法，而且是没选主元的高斯消元法，只是因为Ax=b中系数矩阵A非常特殊，所以就可以采用相对特殊的方法来解方程组。同样，按照常规的步骤，先分析什么是追赶法，再给出追赶法的数学步骤，最后用C++实现这种算法。
（一）追赶法的功能和步骤
明确好目的，正所谓磨刀不误砍柴工，做一件事情事先规划好，那重要性真的是不言而喻。在一些实际问题中，对角占优的三对角线方程组很常见，如热传导方程，形式如下图1：
应用追赶法能解这种三对角线方程组还有很严苛的条件，那就是对角占优，主对角线元素的绝对值要最大，大到什么程度呢？大到绝对值比旁边两条次对角线的值的绝对值之和还要大，这样才能用追赶法来接三对角线方程组，苛刻的条件如下图2：
所以，追赶法就是用来解这类线性方程组的。
追赶法的步骤：追赶法的本质就是没有选主元的高斯消元法，当然系数矩阵***A***都这么特殊了，再用原始的高斯消元法就实在是浪费空间和时间了，特殊的特性就应该用起来，就像做数学证明题时一样，要把事物本身的性质尽可能都用上。实现追赶法最有趣的不是数学方法，而是编程实现，因为之前的系数矩阵***A***是低阶的稠密矩阵，所以编程采用的数据结构就是二维的矩阵，或者说是二重指针，但是实现追赶法，没必要用到二重指针，只要一维数组就够用了，因为就三对角线上有元素，所以一维数组的大小是**n+2(n-1)***，即**3n-1***就够了，然后最重要的就是将系数矩阵***A**的（i,j）***元素同一维数组***B[]***的位置***k***的元素相对应，将这层关系理清楚了，那么追赶法就很容易编程实现了，二者的对应关系如下图3。
所以，数据结构确定，且对应关系理清之后，就又是没选主元的高斯消元法了。数学语言的步骤如下图4。
（二）代码实现
笔者分别将追赶法写成了函数和类两种方式，仔细介绍函数double trde(int n, double A, double b)***，函数就是一个处理器，或是黑匣子，或说是类似映射一样的东西，输入参数，计算处理后得到我们想要的输出，这里参数就是系数矩阵的阶数***n***，以及储存有系数矩阵元素的一维数组***A***（当然这里以指针的形式传输），以及常数向量***b***，输出就是得到解向量，同样是指针形式。代码如下：
double* trde(int n, double* A, double* b)
{
int k, j, m; m = 3 * n - 2;
double s;
for (k = 0; k < n - 1; k++)
{
j = 3 * k; s = A[j];
if (fabs(3) + 1.0 == 1.0)
{
cout << "主对角元素是0，方程无解！！" << endl;
return b;
}
A[j + 1] = A[j + 1] / s;//一下这四行是归一化和消元，对角线下的元素没必要进行操作
b[k] = b[k] / s;//因为后面的回代没用到这些数据
A[j + 3] = A[j + 3] - A[j + 2] * A[j + 1];
b[k + 1] = b[k + 1] - A[j + 2] * b[k];
}
s = A[3 * n - 3];
if (fabs(s) + 1.0 == 1.0)
{
cout << "消元后，主对角元素是0，无解！！" << endl;
return b;
}
b[n - 1] = b[n - 1] / s;
for (k = n - 2; k >= 0; k--)//回代求解方程，相比之下很简单了，累积求和都没必要的
b[k] = b[k] - b[k + 1] * A[3 * k + 1];
return b;
}
主函数的调用：
int main()
{
int n = 5;
double A[13]{ 2,-1,-1,2,-1,-1,2,-1 ,-1,2,-1 ,-1,2 };
double b[5]{ 1,0,0,0,0 };
double *d;
d = trde(n, A, b);
for (int i = 0; i < n; i++)
cout << d[i] << " ";
//这里直接将解显示在黑框框里了
system("pause");
return 0;
}
运行结果，原线性方程组，以及输出的解，如图5,6所示。
最后把类的形式给出来，和全选主元高斯消元法的形式是一样的，input（），输入txt文件名，然后进行操作，a_trde（），输出结果保存在txt文件里，自己命名，接下来是完整的追赶法的类的源码。
//追赶求解三对角方程组
#include <iostream>
#include <fstream>
#include <cmath>
using namespace std;
class trde
{
private:
int n;
double *b, *d;
public:
trde(int nn)
{
n = nn;
b = new double[3 * n - 2]; //动态分配内存空间，一维数组，存放系数矩阵A的元素
d = new double[n];//一维数组存放常数向量b的元素
}
void input(); //从文件读入存放三对角矩阵的向量B以及常数向量D
void a_trde(); //执行追赶法
void output(); //输出结果到文件并显示
~trde()
{
delete[] b;
delete[] d;
}
};
void trde::input() //从文件读入存放三对角矩阵的向量B以及常数向量D
{
int i;
char str1[20];
cout << "\n输入文件名: ";
cin >> str1;
ifstream fin(str1);
if (!fin)
{
cout << "\n不能打开这个文件 " << str1 << endl; exit(1);
}
for (i = 0; i<3 * n - 2; i++) fin >> b[i]; //读入三对角矩阵A中的非0元素
for (i = 0; i<n; i++) fin >> d[i]; //读入常数向量d
fin.close();
}
void trde::a_trde() //执行追赶法，核心的函数
{
int k, j, m;
double s;
m = 3 * n - 2;//第一步，对于k从0~n-2，执行归一和消元：
for (k = 0; k <= n - 2; k++)
{
j = 3 * k; s = b[j];
if (fabs(s) + 1.0 == 1.0)
{
cout << "\n程序工作失败！无解. " << endl;
return;
}
b[j + 1] = b[j + 1] / s;//系数矩阵的归一化
d[k] = d[k] / s;//常数向量的归一化
b[j + 3] = b[j + 3] - b[j + 2] * b[j + 1];//系数矩阵的消元
d[k + 1] = d[k + 1] - b[j + 2] * d[k];//消元过程常数向量的变化
}//第二步：进行回代，倒序进行求解
s = b[3 * n - 3];
if (fabs(s) + 1.0 == 1.0)
{
cout << "\n程序工作失败！无解. " << endl;
return;
}
d[n - 1] = d[n - 1] / s;
for (k = n - 2; k >= 0; k--)
d[k] = d[k] - b[3 * k + 1] * d[k + 1];
}
void trde::output()//结果输出，自己命名txt文件，解向量保存在文件里
{
int i;
char str[20];
cout << "解保存在txt文件里，输入文件名！" << endl;
cin >> str;
ofstream fout(str);
if (!fout)
{
cout << "不能打开这个文件！！" << str << endl;
exit(1);
}
for (i = 0; i < n; i++)
{
fout << d[i] << " ";
cout << d[i] << " ";
}
fout << endl; cout << endl;
fout.close();
}
//接下来是主函数的调用
/*
int main() //主函数
{
trde c(5);
c.input(); //从文件读入存放三对角矩阵的向量B以及常数向量D
c.a_trde(); //执行追赶法
c.output(); //输出结果到文件并显示
system("pause");
return 0;
}
*/
上周六一看数值分析老师的进度，吓一跳，原来自己已经落后这么多了（捂脸哭），所以又学习编程实现解三对角线方程组的追赶法，再补充一篇，自己本没有时间观念，但是也算体会到了动漫团队持续一周更新一集动漫的艰辛，若是让我每一周都学习并用C++实现一个算法，那我肯定坚持不住了。佩服那些坚持周更的作者们，真是不容易，正好前阵子看了鸣人的原画师作者黄成希的一席演讲，真的很不容易啊，那样的劳动强度，为我们展示了如此精彩的动漫热血视频，人家用十年时间成为了最优秀的原画师，想想我十年能够成为最优秀的某某某，还真是羞愧难当啊，只能去追求美好生活了，哈哈~