第二章有限差分法初步

1、一、差商与微商一、差商与微商 第二章：第二章：有限差分法初步有限差分法初步1 有限差分法基本概念有限差分法基本概念（i）、有限差分的数学基础是用）、有限差分的数学基础是用差商差商代替代替微商微商。 有如下两种数学形式：有如下两种数学形式： （i）微商微商（导数）的定义（导数）的定义)(xT是连续函数，则它的导数为：是连续函数，则它的导数为：若若xTxxTxxTdxdTxx00lim)()(lim（2.1）式（式（2.1）右边）右边xT是有限的是有限的差商差商。 x与与T都不为零，都不为零， 而式（而式（2.1）左边）左边 dxdT是是 xT当当 x趋于零时极限情形下的差商，称之趋于零时极限情形

2、下的差商，称之微商微商。 在在x没有到达零之前，没有到达零之前， xT只是只是 dxdT的近似。的近似。 xT趋于趋于dxdT的过程认为是的过程认为是近似近似向向精确精确过渡，过渡， 用用xT代替代替dxdT就是就是精确精确向向近似近似过渡过渡。 两者的差值两者的差值dxdTxT表示表示差商差商代替代替微商微商的的偏差偏差。 （ii）偏差偏差-Taylor级数展开级数展开222! 2)()()(dxTdxdxdTxxTxxTnnndxTdnx!)(（2.2） 稍加整理后可写成：稍加整理后可写成：22! 2)()(dxTdxdxdTxxTxxTxTnnndxTdnx!)(1（2.3） xT可见可

3、见 与与 只能是近似相等。只能是近似相等。 dxdT偏差偏差为：为：)(0 x（iii）微商与差商的几何意义）微商与差商的几何意义 xx+xx-xT(x+x)T(x-x)T(x)T(x+x)-T(x-x)2xT(x+x)-T(x)xT(x)-T(x-x)xdT(x)dx图图2-1 差商与微商的比较差商与微商的比较图图2.1表示了差商与微商之间的关系。应当指表示了差商与微商之间的关系。应当指出，用不同方法得到的差商去代替微商，它们出，用不同方法得到的差商去代替微商，它们带来的偏差是不同的。带来的偏差是不同的。向向右右（前）差商：（前）差商：（iV）差商的几种表示）差商的几种表示 xxTxxTdx

4、dT)()(（2.4） 向向左左（后）差商：（后）差商：xxxTxTdxdT)()(（2.5） 中心差商中心差商，取向右差商与向左差商的平均值：，取向右差商与向左差商的平均值：xxxTxTxxTxxTdxdT)()()()(21xxxTxxT2)()(（2.6） 偏差分析：偏差分析： 将将Taylor级数写成：级数写成：)(2)()()()(2xTxxTxxTxxT （2.7） Taylor级数还可写成：级数还可写成：)()(! 3)(43xOxTx （2.8） )()(! 3)(43xOxTx )(! 2)()()()(2xTxxTxxTxxT 由式（由式（2.7）可得）可得)(! 2)()

5、()(xTxxTxxTxxT )( xO 由式（由式（2.8）可得）可得)(! 2)()()(xTxxTxxxTxT （2.9） )( xO （2.10） （2.9）+（2.10），得到），得到)(! 3)()(2)()(2xTxxTxxxTxxT 2)( xO （2.11） 比较式（比较式（2.9）、（）、（2.10）、（）、（2.11）可看到，用不）可看到，用不同的差商形式去代替微商，所带来的偏差是不同的。同的差商形式去代替微商，所带来的偏差是不同的。这些这些偏差偏差都是截去了都是截去了Taylor级数展开式中的高阶项而级数展开式中的高阶项而引起的，常称引起的，常称“截断误差截断误差”。

6、用向右差商与向左差商代替微商，其截断误差为与用向右差商与向左差商代替微商，其截断误差为与 x同量级的小量同量级的小量 ；)( xO 2)( x同量级的小量；同量级的小量； 中心差商的截断误差小于向右差商或向左差商。中心差商的截断误差小于向右差商或向左差商。 而用中心差商代替微商，其截断误差是与而用中心差商代替微商，其截断误差是与 讨论：讨论：上述一阶差商一般仍是上述一阶差商一般仍是x的函数，对它们还可以的函数，对它们还可以求差商。这种求差商。这种一阶差商的差商一阶差商的差商称为称为二阶差商二阶差商，它是二阶微商的近似，常用向右差商的向左差它是二阶微商的近似，常用向右差商的向左差商来近似二阶微商

7、，即：商来近似二阶微商，即：xxxxTxTxxTxxTdxTd)()()()(222)()()(2)(xxxTxTxxT（V）二阶差商二阶差商根据式（根据式（2.7）+（2.8），可得），可得)()()()(2)(2xTxxxTxTxxT 2)( xO （2.12） 由式（由式（2.12）知，二阶差商的截断误差也为与）知，二阶差商的截断误差也为与 同阶的小量。同阶的小量。 2)( xn结论：结论： 由于用差商代替微商必然带来截断误差，相应地用由于用差商代替微商必然带来截断误差，相应地用差分方程代替微分方程也必然带来截断误差。这是差分方程代替微分方程也必然带来截断误差。这是有限差分法固有的。因此

8、，在应用有限差分法进行有限差分法固有的。因此，在应用有限差分法进行数值解时，必须对差分的构成及其对方程造成的误数值解时，必须对差分的构成及其对方程造成的误差引起注意。差引起注意。 二、从微分形式出发的差分格式二、从微分形式出发的差分格式n图图2.2给出了一个简单边界值问题。给出了一个简单边界值问题。 (i, j)(i+1, j)(i-1, j)(i, j-1)(i, j+1).qTwThxyL1L2图图2.2 矩形区域离散化矩形区域离散化问题是求图问题是求图2.2所示的边值问题的解，其数学表达所示的边值问题的解，其数学表达如下，方程：如下，方程：02222 kqyTxT（2.13） n边界条件

9、：边界条件：0 x20Ly )(TThxTk0y10Lx qyTk 1Lx 20Ly 0 xT2Ly 10Lx wTT （2.14） （2.15） （2.15） （2.16） 式（式（2.13）、（）、（2.14）、（）、（2.15）、（）、（2.16） 所示所示定解问题定解问题解法。解法。 在问题的提法已经明白之后，差分在问题的提法已经明白之后，差分格式的构成可通过以下几步来实现：格式的构成可通过以下几步来实现： （i） 区域离散法区域离散法下面分别予以说明下面分别予以说明: （vi） 构成差分格式构成差分格式（ii） 建立区域内差分方程建立区域内差分方程（iii） 边界条件的差分形式边界条

10、件的差分形式 1区域离散化区域离散化 所谓离散化，就是把几何上连续的区域用一系列所谓离散化，就是把几何上连续的区域用一系列网格线把它划分开。一般说来，网格形式应视几何网格线把它划分开。一般说来，网格形式应视几何区域的不同而不同，对于矩形区域而言，用矩形的区域的不同而不同，对于矩形区域而言，用矩形的网格，如图网格，如图2.2，用五条水平网线与五条垂直网线，用五条水平网线与五条垂直网线把矩形区域离散掉。网线与网线的交点称之为把矩形区域离散掉。网线与网线的交点称之为“节节点点”，节点与节点的距离称之为，节点与节点的距离称之为步长步长，x方向的步方向的步长表示为长表示为 ，y方向的步长表示为方向的步长

11、表示为 。xy 节点编号：节点编号：为便于计算，需对节点逐个编号。为便于计算，需对节点逐个编号。常用（常用（i，j） 表示节点位置，其中，表示节点位置，其中，i、j是是与网线相对应的正整数。与网线相对应的正整数。i，j 的排列：的排列：可有不同的方式。可有不同的方式。习惯上，与习惯上，与x、y 轴相一致，轴相一致，i 由左而右逐个增长，由左而右逐个增长，j 由下而上逐个增长。由下而上逐个增长。但也有，考虑到与矩阵的格式相一致，但也有，考虑到与矩阵的格式相一致， i 表示行数，表示行数，由上而下逐个增长，由上而下逐个增长，j 表示列数，由左而右逐个增长。表示列数，由左而右逐个增长。这种从上到下，

12、从左到右的编排与一般书写这种从上到下，从左到右的编排与一般书写习惯也是一致的，因此，在计算机上算题也常被采用。习惯也是一致的，因此，在计算机上算题也常被采用。在本章中，大都采用与坐标相一致的编排方法。在本章中，大都采用与坐标相一致的编排方法。n在区域内的节点称在区域内的节点称“内节点内节点”，在边界上的，在边界上的节点称节点称“边界节点边界节点”。图。图2.2所示边界是规所示边界是规则的则的,则节点或在区域内，或正好落在边界上。则节点或在区域内，或正好落在边界上。xy 步长步长 或或 可以是不变的常量，即可以是不变的常量，即等步长等步长，也可以，也可以 在区域内的不同处是不同的，即在区域内的不

13、同处是不同的，即变步长变步长。如果区域。如果区域 内各处的温度梯度变化很大，则在温度变化剧烈处，内各处的温度梯度变化很大，则在温度变化剧烈处， 网格布得密些；在温度变化不剧烈处，网格布得疏网格布得密些；在温度变化不剧烈处，网格布得疏 些。至于网格布置多少，步长取多大为宜，要根据些。至于网格布置多少，步长取多大为宜，要根据 具体问题，兼顾到计算的精确度与计算的工作量等具体问题，兼顾到计算的精确度与计算的工作量等 因素而定。因素而定。步长步长:n从物理方面对区域离散化可作这样的理解，即从物理方面对区域离散化可作这样的理解，即认为区域内离散的每个节点，都集中着它周围认为区域内离散的每个节点，都集中着

14、它周围区域（尺度为步长）的热容，或者说，区域内区域（尺度为步长）的热容，或者说，区域内连续分布的热容都被分别地集中到离散的节点连续分布的热容都被分别地集中到离散的节点上去了。这样，节点的温度代表着它周围区域上去了。这样，节点的温度代表着它周围区域的某种平均温度。一系列离散的节点温度值代的某种平均温度。一系列离散的节点温度值代表着连续区域内的温度分布。表着连续区域内的温度分布。 jiT,区域离散化物理理解区域离散化物理理解:节点（节点（i，j）处的温度表示成）处的温度表示成 。 2差分方程代替微分方程差分方程代替微分方程 在上节我们已对有限差分法的数学基础作了简要在上节我们已对有限差分法的数学基

15、础作了简要 的介绍，说明了如何用差商代替微商，以及由此的介绍，说明了如何用差商代替微商，以及由此带来的误差。这里介绍用差商代替微商的办法来带来的误差。这里介绍用差商代替微商的办法来处理导热方程（处理导热方程（2.13），得到相应的），得到相应的差分方程差分方程。 方程（方程（2.13） 02222 kqyTxT对区域内各个点都成立的，当然对任意一个内节对区域内各个点都成立的，当然对任意一个内节点（点（i，j）也成立。）也成立。 或者说，在（或者说，在（i，j）处存在二阶偏微商）处存在二阶偏微商 与与 ，这些二阶偏微商所对应的差商可表示成：，这些二阶偏微商所对应的差商可表示成：jixT,22ji

16、yT,222,1,122)(2xTT,TxTjijiji,ji)(2xO i=2，3，4；j=2，3，4（2.17） 21,122)(2yTT,TyTjijiji,ji)(2yO （2.18） i=2，3，4；j=2，3，4其中，其中， 与与 表示相应的表示相应的二阶差商二阶差商与与二阶偏微商二阶偏微商的差别为的差别为 与与 的数量级。的数量级。)(2xO )(2yO 2)( x2)( y将式（将式（2.17）与（）与（2.18）代入方程（）代入方程（2.13），得），得21,1,2, 1, 1)(2)(2yTTTxTTTjijijijijiji0)()(22 yxOkq（2.19） 式（式（

17、2.19）中去掉）中去掉 项，得到项，得到)()(22yxO21,1,2, 1, 1)(2)(2yTTTxTTTjijijijijiji0 kq（2.20） i=2，3，4；j=2，3，4式（式（2.20）被称为对应于方程（）被称为对应于方程（2.13）的的差分方程差分方程。方程（。方程（2.20）被改写成：）被改写成：1,2, 12,22)(1)(1)(2)(2jijijiTxTxTyxkqTyTyjiji 1,21,2)(1)(1（2.21） 若若 ， ，则式（，则式（2.21）又被）又被改写成：改写成： yx0q 041,1, 1, 1,jijijijijiTTTTT或或)(411,1,

18、 1, 1,jijijijijiTTTTT （2.22） 物理意义：物理意义：一点（一点（i, j）处的温度是它周围）处的温度是它周围4点温度的平均值。点温度的平均值。)()(22yxO)()(22yxO由于差分方程（由于差分方程（2.20）是从式（）是从式（2.19）中）中去掉项去掉项 得来的，称去掉得来的，称去掉的项的项 为差分方程（为差分方程（2.19）的）的截断误差。当截断误差。当 与与 趋于零时，差分方程趋于零时，差分方程的截断误差也趋于零，即差分方程逼近微的截断误差也趋于零，即差分方程逼近微分方程。我们称这种逼近的差分方程与相分方程。我们称这种逼近的差分方程与相应的微分方程为应的微

19、分方程为“相容相容”。 xy3边界条件的差分形式边界条件的差分形式对流换热边界条件：对流换热边界条件：20, 0)(LyxTThxTk（2.14） -用用T 对对 x 的向前差商代替式（的向前差商代替式（2.14）中的）中的 T 对对 x 的一阶偏微商，使式（的一阶偏微商，使式（2.14）变）变 成为如下差分形式成为如下差分形式 ：这里介绍用这里介绍用差商代替微商差商代替微商的办法把定解问题中的办法把定解问题中的各种的各种边界条件边界条件表示成表示成差分差分的形式的形式.)(, 1TThxTTkjijijiTkxhTTkxhjiji, 1,1或或 （2.23） i=1；j=2，3，4热流边界条

20、件：热流边界条件：10, 0LxyqyTk （2.15） 用用T 对对 y 的的向前差商向前差商代替式（代替式（2.15）中）中T 对对 y的一阶偏微商，使式（的一阶偏微商，使式（2.15）变成为如下差分）变成为如下差分形式形式: qyTTkjiji ,1,kyqTTjiji 1,i=1，2，3，4，5；j=1 （2.24） 或或绝热边界条件：绝热边界条件： 2Ly0Lx0 xT，变成为：变成为：0TTj1iji ,i=5；j=2，3，4 （2.25） （2.16） 给定温度边界条件：给定温度边界条件：wjiTT ,i=1，2，3，4，5；j=5 （2.26） 至此，我们对至此，我们对全部节点

21、全部节点，包括，包括内节点内节点与与边界节点边界节点，都用差分形式代替了原来的函数形式。对于内节都用差分形式代替了原来的函数形式。对于内节点上差分形式，我们通称差分方程，因为内节点点上差分形式，我们通称差分方程，因为内节点上温度都是未知的。对于边界节点的差分形式，上温度都是未知的。对于边界节点的差分形式，在边界节点的温度为未知量时，它是差分方程。在边界节点的温度为未知量时，它是差分方程。而对于边界节点为给定的温度时，得到的就不是而对于边界节点为给定的温度时，得到的就不是差分方程了。但在实际应用中，人们往往习惯地差分方程了。但在实际应用中，人们往往习惯地把由内节点与边界节点建立起来的差分形式，都

22、把由内节点与边界节点建立起来的差分形式，都统称为统称为差分方程差分方程。笼统地讲，。笼统地讲，一个节点对应一个一个节点对应一个差分方程。差分方程。在边界节点的处理方面还有几点需要在边界节点的处理方面还有几点需要强调强调： （i） 每一边界节点只应属于一种边界条件。如图每一边界节点只应属于一种边界条件。如图2.2 中，中，i=1，j=5的节点只属于边界条件式（的节点只属于边界条件式（2.16）。）。 （iii）若边界节点不正好落在区域的边界）若边界节点不正好落在区域的边界 上，则需对它们进行特殊处理。上，则需对它们进行特殊处理。 （ii）对应不同边界的差分方程式（）对应不同边界的差分方程式（2.

23、23）、）、2.24）、）、 （2.25）都是用一阶向前差商代替一阶微商得到的，）都是用一阶向前差商代替一阶微商得到的， 也即它们的截断误差为也即它们的截断误差为)( xO 或或)( yO 量级，与内节点差分方程的截断误差相比，量级，与内节点差分方程的截断误差相比，低了一个量级。这一点也是从微分形式出发低了一个量级。这一点也是从微分形式出发建立差分格式的建立差分格式的弱点弱点。 4差分格式的构成差分格式的构成由于式（由于式（2.13）、（）、（2.14）、（）、（2.15）、）、（2.16）、（）、（2.17）所表示的方程式与边界条件）所表示的方程式与边界条件都是线性的，由此而得到的内节点与边

24、界节点的差都是线性的，由此而得到的内节点与边界节点的差分方程也都是分方程也都是线性代数方程线性代数方程。由全部节点的差分方。由全部节点的差分方程构成一个程构成一个线性代数方程组线性代数方程组。在这个方程组里，。在这个方程组里，方方程式的个数等于节点的个数程式的个数等于节点的个数。针对前面讨论的例子，。针对前面讨论的例子，我们可以看到，这个有我们可以看到，这个有25个节点所构成的方程组，个节点所构成的方程组，只需要只需要5个式子，即式（个式子，即式（2.21）、（）、（2.23）、）、（2.24）、（）、（2.25）、（）、（2.26）就可表示，这里）就可表示，这里每个式子都表示几个节点方程。每

25、个式子都表示几个节点方程。 n也就是说，在组成方程组时，不必把每个节点也就是说，在组成方程组时，不必把每个节点方程都写出来，而只要写出几个规格化了的方方程都写出来，而只要写出几个规格化了的方程就可以了。因此，人们把规格化了的，由内程就可以了。因此，人们把规格化了的，由内节点与边界节点全部差分方程所构成的线性代节点与边界节点全部差分方程所构成的线性代数方程组，称之为数方程组，称之为“差分格式差分格式”。 一般地说，差分格式被写成如下的形式：一般地说，差分格式被写成如下的形式： nnnnnnnnnnbTTTbTTTbTTT22112222212111212111（2.27） ija其中，其中，n是

26、节点数，也即方程个数，每个方是节点数，也即方程个数，每个方程对应一个节点。程对应一个节点。 和和bi（i=1，2，n；j=1，2，n）都是常数）都是常数 ，ija方程组（方程组（2.27）可被进一步写成矩阵的形式）可被进一步写成矩阵的形式 : IIBTA（2.28） nnnnnnnITTTTA21212222111211nIbbbB21其中其中n如果在方程组中去掉其中已知温度节点的那如果在方程组中去掉其中已知温度节点的那些方程，由此构成的线性代数方程组中，方些方程，由此构成的线性代数方程组中，方程的个数等于温度未知的节点个数，也即方程的个数等于温度未知的节点个数，也即方程组的未知数。这样的方程

27、组也是差分格式。程组的未知数。这样的方程组也是差分格式。也可写成式（也可写成式（2.28）。）。n综上所述，用有限差分法对式（综上所述，用有限差分法对式（2.13）、）、（2.14）、（）、（2.15）、（）、（2.16）、（）、（2.17）所）所组成的边值问题的数值处理，组成的边值问题的数值处理，最终归结成求解线最终归结成求解线性代数方程组（性代数方程组（2.28），方程组的解即各节点的，方程组的解即各节点的温度。如果整个区域的节点足够多，那么，离散温度。如果整个区域的节点足够多，那么，离散节点的温度分布就近似代替了区域内的连续温度节点的温度分布就近似代替了区域内的连续温度分布。分布。n为便

28、于讨论各种差分格式的优缺点，最好为便于讨论各种差分格式的优缺点，最好把方程组（把方程组（2.28）中系数矩阵具体地写出）中系数矩阵具体地写出来。但当我们着手书写由式（来。但当我们着手书写由式（2.21）、）、（2.23）、（）、（2.24）、（）、（2.25）、）、（2.26）组成的代数方程组时，发现它所）组成的代数方程组时，发现它所占的版面太大，造成印刷的困难。占的版面太大，造成印刷的困难。n所以，为便于书写，采用图所以，为便于书写，采用图2.3所示的网格所示的网格 (2, 3).qTwThxyL1L2(1, 3)(3, 3)(2, 2)(3, 2)(1, 2)(2, 1)(3, 1)(1,

29、 1)并假定并假定 ckyqbkxhakxqyx ,)(, 12n得到由全部节点组成的线性代数方程组为：得到由全部节点组成的线性代数方程组为：cTTcTTcTTTTaTTTTTbTTTbTTTTTTwww1 , 32, 31 , 22, 21 , 12, 12, 32, 21 , 22, 32, 22, 13 , 22, 22, 13 , 33 , 23 , 104)1 (（2.29）n表示成矩阵形式：表示成矩阵形式：1 , 31 , 21 , 12, 32, 22, 13 , 33 , 231111111111141111111TTTTTTTTTbIcccabTTTTwww0（2.30） n

30、式（式（2.29）或（）或（2.30）构成差分格式。若在方程）构成差分格式。若在方程组中去掉已知温度节点所对应的方程，即第组中去掉已知温度节点所对应的方程，即第1至第至第3个方程，则式（个方程，则式（4.2.24）被改写成：）被改写成： cccTabTTTTTTTbw0111111111141111 , 31 , 21 , 12, 32, 22, 1（2.31） n将式（将式（2.31）与式（）与式（4.2.28）进行对照，）进行对照，即可得到矩阵即可得到矩阵 、 、 的各个元素。的各个元素。 这里特别提醒读者注意，在式（这里特别提醒读者注意，在式（2.31）与）与（2.28）中温度）中温度T

31、下角码的对应关系，在讨下角码的对应关系，在讨论二维稳定导热问题时，人们习惯用二个论二维稳定导热问题时，人们习惯用二个角码（角码（i，j）来表示节点位置及对应的曾）来表示节点位置及对应的曾度度 。IA TIBjiT, （二）解线性代数方程组的直接法（二）解线性代数方程组的直接法n这里介绍计算机上常用的解线性代数方程这里介绍计算机上常用的解线性代数方程组的直接法。组的直接法。n大家知道，线性方程组大家知道，线性方程组 BTA（2.32） 中只要矩阵的行列式中只要矩阵的行列式 ，方程组，方程组（2.32）就有唯一解，其表达式为）就有唯一解，其表达式为:0detAnjAATjj, 2 , 1det/d

32、et（2.33） 其中其中 为用右端向量为用右端向量 替换行列式替换行列式 的第的第j列而得的，这一公式为著名的列而得的，这一公式为著名的Carmer法则。法则。 jAdet BAdetn显然，按显然，按Cramer法则求解方程组（法则求解方程组（2.32）需要）需要计算计算n+1个个n阶行列式。每个阶行列式。每个n阶行列式按直接展阶行列式按直接展开办法来算，需作（开办法来算，需作（n-1）n！次乘法和！次乘法和n!次加次加法运算。当法运算。当n=30时，共约需完成时，共约需完成 次乘法和加次乘法和加法运算，这是一个十分惊人的数字，即使在一台法运算，这是一个十分惊人的数字，即使在一台每秒作一亿

33、次运算的计算机上完成这一计算也是每秒作一亿次运算的计算机上完成这一计算也是不可能的。所以，尽管这种办法也是一种直接法，不可能的。所以，尽管这种办法也是一种直接法，并且理论上可行，但实际上是无法进行求解的。并且理论上可行，但实际上是无法进行求解的。即使采用其它办法来计算行列式，按即使采用其它办法来计算行列式，按Cramer法法则求解的工作量也比通常的直接法大得多。因而，则求解的工作量也比通常的直接法大得多。因而，Cramer法则对于数值计算来说是没有什么用处法则对于数值计算来说是没有什么用处的，仅在一些特殊场合才有用。的，仅在一些特殊场合才有用。3310n另外，矩阵求逆也是人们熟悉的求解线性方程

34、另外，矩阵求逆也是人们熟悉的求解线性方程的一种方法。的一种方法。1BAT(2.34) 但求逆矩阵但求逆矩阵 时，要计算时，要计算 个个( n -1)阶行列阶行列式，和一个式，和一个n阶行列式，它的计算工作量也是阶行列式，它的计算工作量也是相当可观的，对于相当可观的，对于n较大的情况也没有什么现较大的情况也没有什么现实意义。实意义。2n1An现在计算机上常用的直接解法大多数是以系数矩现在计算机上常用的直接解法大多数是以系数矩阵的三角形化为基础的。就是说说，先对方程组阵的三角形化为基础的。就是说说，先对方程组进行变换，使其化为等价的（即具有相同解的）进行变换，使其化为等价的（即具有相同解的）三角形

35、方程组。由于三角形方程组的求解十分容三角形方程组。由于三角形方程组的求解十分容易，原方程的求解问题即告解决。下面我们简要易，原方程的求解问题即告解决。下面我们简要地讨论这一类型的方法。为讨论方便起见，我们地讨论这一类型的方法。为讨论方便起见，我们首先叙述三角形方程的解法，然后叙述将原方程首先叙述三角形方程的解法，然后叙述将原方程化为等价三角形方程组的方法。由于计算机的字化为等价三角形方程组的方法。由于计算机的字长是有限的，每次运算之后还要对结果进行舍入，长是有限的，每次运算之后还要对结果进行舍入，所以，虽然理论上直接法在有限步内可以得到精所以，虽然理论上直接法在有限步内可以得到精确解，但计算机

36、上实际得到的只是近似解。确解，但计算机上实际得到的只是近似解。n1三角形方程组的解法三角形方程组的解法 所谓三解形方程组是指下面两种形式的方程组所谓三解形方程组是指下面两种形式的方程组：nnnnnnnbTlTlTlTlbTlTlTlbTlTlbTl332211333323213122221211111(2.35) nnnnnnnnnnnnnnndTudTuTudTuTuTudTuTuTuTu1, 111, 12232222211313212111或或(2.36) 方程组方程组(2.35)叫作叫作下三角形方程组下三角形方程组，方程组，方程组(2.36)叫作叫作上三角形方程组上三角形方程组。 n若

37、用矩阵符号可分别写为：若用矩阵符号可分别写为：,DTUBTL其中其中L为方程组为方程组(2.35)的系数所构成的下三的系数所构成的下三角形矩阵，其元素满足关系：角形矩阵，其元素满足关系：filij 0U为方程组为方程组(2.36)的系数所构成的上三角的系数所构成的上三角形矩阵，其元素满足关系：形矩阵，其元素满足关系： fiuij 0n三角形方程组的求解是很简单的。方程组三角形方程组的求解是很简单的。方程组(2.35)的计算公式可归结为：的计算公式可归结为：nilTlTlTlbTlbTiiiiiiiii, 3, 2/ ),(/1122111111(2.37) 这个计算过程通常也只作这个计算过程通

38、常也只作前推前推过程。过程。n对于方程组对于方程组(2.36)，其计算公式可归结为：，其计算公式可归结为：1 , 2, 1/ ),(/22,11,nniuTuTuTudTudTiinniiiiiiiiinnnn(2.38) 这个计算过程通常也叫作这个计算过程通常也叫作回代方程回代方程。由以上分析可以看到，只要把方程组化成了等由以上分析可以看到，只要把方程组化成了等价的三角形方程组，求解就容易了。价的三角形方程组，求解就容易了。 2Gauss消去法消去法 Gauss消去法（简称消去法）的提出已有相当长消去法（简称消去法）的提出已有相当长的时间了，是一种古老的方法。然而，近年来在的时间了，是一种古

39、老的方法。然而，近年来在计算机上求解线性代数方程组的实践表明，它仍计算机上求解线性代数方程组的实践表明，它仍是直接法中最常用的一种方法，也是最有效的方是直接法中最常用的一种方法，也是最有效的方法之一。其基本思想是，用逐次消去一个未知数法之一。其基本思想是，用逐次消去一个未知数的办法把原来的方程组化为等价的（具有相同解）的办法把原来的方程组化为等价的（具有相同解）三角形方程组。这样，求解就很容易了。三角形方程组。这样，求解就很容易了。n假定把要求的假定把要求的n阶线性的方程组阶线性的方程组(2.32)改写成如改写成如下形式：下形式： ) 1 () 1 (3) 1 (32) 1 (21) 1 (1

40、) 1 (2) 1 (23) 1 (232) 1 (221) 1 (21) 1 (1) 1 (13) 1 (132) 1 (121) 1 (11nnnnnnnnnnnbTaTaTaTabTaTaTaTabTaTaTaTa(2.39) 用矩阵符号记为用矩阵符号记为)1()1(BTAn其中其中 为为 方阵，方阵， 为为 向量，它们向量，它们分别为：分别为：分别从原方程组的第二个方程减去第一个方分别从原方程组的第二个方程减去第一个方程乘以程乘以 ，第三个方程减去第一个方程，第三个方程减去第一个方程乘以乘以 ，如此等等，即可消去后面，如此等等，即可消去后面n-1个方程中的未知量个方程中的未知量 。 )

41、 1 (11) 1 (21/aa)1(11)1(31/aa1T)1(Ann)1(B1n) 1 () 1 (2) 1 (1) 1 () 1 () 1 (3) 1 (2) 1 (1) 1 (2) 1 (23) 1 (22) 1 (21) 1 (1) 1 (13) 1 (12) 1 (11) 1 (,nnnnnnnnbbbBaaaaaaaaaaaaAn这时方程组这时方程组(2.39)就变为如下等价方程组：就变为如下等价方程组： )2()2(3)2(3)2(2)2(3)2(33)2(33)2(32)2(2)2(23)2(23)2(22) 1 (1) 1 (13) 1 (13) 1 (12) 1 (11

42、_0_0_0_0nnnnnnnnnnnnbTaTaabTaTaabTaTaabTaTaaan表示成矩阵形式为：表示成矩阵形式为：)2()2(BTA其中其中)2()2(3)2(2) 1 (1)2()2()2(3)2(2)2(3)2(33)2(23)2(2)2(23)2(22) 1 (1) 1 (13) 1 (12) 1 (11)2(b,000nnnnnnnnbbbBaaaaaaaaaaaaaAn若令若令 ，则系数，则系数 和和 的计的计 算公式应为：算公式应为：) 1 (11) 1 (11/aamaii)2(ija)2(ibnjibmbbamaaiiijiijij, 3, 2,) 1 (11)

43、1 ()2() 1 (11) 1 ()2(类似地，分别从上述等价方程组的第三个方程减类似地，分别从上述等价方程组的第三个方程减去第二个方程乘以去第二个方程乘以 ，第四个方程减去第，第四个方程减去第二个方程乘以二个方程乘以 ，如此等等，即可进一步，如此等等，即可进一步消去后面消去后面n-2个方程中的未知量个方程中的未知量T2，而将方程，而将方程(2.39)变为如下等价形式：变为如下等价形式： )()(/222232aa)()(/222242aa)3()3(3)3(3)3(44)3(43)3(43)3(33)3(33)3(33)2(2)2(23)2(232)2(22)1(1)1(13)1(132)

44、1(121)1(11nnnnnnnnnnnbTaTabTaTabTaTabTaTaTabTaTaTaTa表示成矩阵形式为：表示成矩阵形式为：)3()3(BTAn其中其中)3()3(4)3(3)2(2) 1 (1)3()3()3(3)3(4)3(43)3(3)3(33)2(2)2(23)2(22) 1 (1) 1 (13) 1 (12) 1 (11)3(,000000nnnnnnnnbbbbbBaaaaaaaaaaaaaA若令若令 ，则系数，则系数 和和 的计算的计算公式应为：公式应为：)2(22)2(22/aamii)3(ija)3(ibnjibmbbamaaiiijiijij, 4, 3,)

45、 2(22) 2() 3 () 2(22) 2() 3 (n上述的消去步骤还可以进行下去。如此继续之，上述的消去步骤还可以进行下去。如此继续之，重复上述步骤重复上述步骤(n -1)次以后，我们即可得到如次以后，我们即可得到如下等价三角形方程组：下等价三角形方程组：)()()3(3)3(33)3(33)2(2)2(23)2(232)2(22) 1 (1) 1 (13) 1 (132) 1 (121) 1 (11nnnnnnnnnnnnbTabTaTabTaTaTabTaTaTaTa(2.40) 表示为矩阵形式：表示为矩阵形式：)()(nnBTA其中其中 为如下上三角形矩阵，为如下上三角形矩阵，

46、为为 向量；向量；)(nA)(nB1n)() 3 (3) 2(2) 1 (1)()() 3 (3) 3 (33) 2(2) 2(23) 2(22) 1 (1) 1 (13) 1 (13) 1 (11)(,nnnnnnnnnnnbbbbBaaaaaaaaaaA0(2.41) n三角形方程组三角形方程组 很容易用前述的回很容易用前述的回代过程代过程(2.38)求解，这样就完成了消去法求解求解，这样就完成了消去法求解n阶线性代数方程组的过程。从原来方程组阶线性代数方程组的过程。从原来方程组(2.39)得出等价三角形方程组得出等价三角形方程组(2.40)的过程称之为消去的过程称之为消去过程。采用前面的

47、记号，我们可将消去过程的计过程。采用前面的记号，我们可将消去过程的计算公式归结为对于算公式归结为对于 ，递推地，递推地计算如下各量：计算如下各量：)()(nnBTA1n21k,nijkjabaabbnjknikaaaaanjkibbaakijkkkkkkikkikikkjkkkkikkijkijkikikijkij1,1011,) 1()()()()() 1()()()()() 1()() 1()() 1(2.42) 用用Gauss消去法解线性代数方程组时，为了能求到消去法解线性代数方程组时，为了能求到最后的解（尽管已经具备了最后的解（尽管已经具备了 的条件），并的条件），并使解尽可能的精确，

48、应注意如下两点：使解尽可能的精确，应注意如下两点：0Adet)(kiiaiia)(kiia)( kiia（i）系数矩阵中对角线上的元素）系数矩阵中对角线上的元素 都不应都不应 为零，因为在消元的过程中，不断用为零，因为在消元的过程中，不断用 作为除数，倘若有一个作为除数，倘若有一个 为零，为零， 计算就无法进行下去。计算就无法进行下去。（ii）在系数矩阵）在系数矩阵A每一行的元素中，每一行的元素中， 的绝对值最好比同一行的其他元素都的绝对值最好比同一行的其他元素都 大，这大，这 样在作除法运算时，引起的舍样在作除法运算时，引起的舍 入误差就比较小。入误差就比较小。n对照上节中差分格式，不难看到

49、，由有限对照上节中差分格式，不难看到，由有限差分法得到的系数矩阵是能够满足以上两差分法得到的系数矩阵是能够满足以上两个条件的。因为每个节点方程（第个条件的。因为每个节点方程（第i个方程）个方程）都代表着一个单元（第都代表着一个单元（第i个单元）与周围单个单元）与周围单元或外界环境的热量交换关系。在这些热元或外界环境的热量交换关系。在这些热量交换中，无疑都与该单元的温度量交换中，无疑都与该单元的温度(Ti)有关，有关，Ti的系数的系数aii当然不能为零。当然不能为零。 n而且在第而且在第i个方程中，个方程中，Ti的地位比其周围单的地位比其周围单元温度更为突出，表现在系数上，它的绝元温度更为突出，

50、表现在系数上，它的绝对值总是最大的。对于给定温度的节点方对值总是最大的。对于给定温度的节点方程而言，这种性质更明显。因为方程中除程而言，这种性质更明显。因为方程中除了该节点温度以外再也没有别的节点温度，了该节点温度以外再也没有别的节点温度，它的系数当然也就最大了。它的系数当然也就最大了。n综上所述，由于对稳定导热问题用有限综上所述，由于对稳定导热问题用有限差分法得到的代数方程具有上述性质，差分法得到的代数方程具有上述性质，因此在求解方程组时，可大胆放心使用因此在求解方程组时，可大胆放心使用Gauss消去法。（对于更一般的方程组，消去法。（对于更一般的方程组，目前更多采用主元素消去法，而不用目前

51、更多采用主元素消去法，而不用Gauss消去法。）消去法。）（三）解线性代数方程组的迭代法（三）解线性代数方程组的迭代法n前面介绍的解线性代数方程组的直接法对于前面介绍的解线性代数方程组的直接法对于阶数不是很高的问题是非常有效的，这种场阶数不是很高的问题是非常有效的，这种场合一般不使用下面介绍的迭代法。然而对于合一般不使用下面介绍的迭代法。然而对于阶数很高的稀疏矩阵，尽管提出了很多特殊阶数很高的稀疏矩阵，尽管提出了很多特殊的直接法来处理它们，在运算量和存储量的的直接法来处理它们，在运算量和存储量的节省方面也取得了很大的进展，但仍然难于节省方面也取得了很大的进展，但仍然难于克服存储需要量大的缺点，

52、特别在不具备大克服存储需要量大的缺点，特别在不具备大型计算机的条件下，采用下面介绍的迭代法型计算机的条件下，采用下面介绍的迭代法更为合适。更为合适。 n迭代法的优点：迭代法的优点：由于不需要存储系数矩阵的零由于不需要存储系数矩阵的零元素，所以占用的存储单元少。同时程序也比元素，所以占用的存储单元少。同时程序也比较简单，对于稳定导热用有限差分法所得到的较简单，对于稳定导热用有限差分法所得到的方程组，求解收敛较快，因此广泛地被采用。方程组，求解收敛较快，因此广泛地被采用。迭代法的缺点是：迭代法的缺点是：它所得到的是一种近似解，它所得到的是一种近似解，在运算过程中需要进行多次迭代才能达到收敛在运算过

53、程中需要进行多次迭代才能达到收敛指标的要求，而迭代次数事先是不知道的，这指标的要求，而迭代次数事先是不知道的，这样，往往要耗费较多的时间。因此，一般地讲，样，往往要耗费较多的时间。因此，一般地讲，直接法与迭代法各有优缺点。直接法与迭代法各有优缺点。 n迭代法所要讨论的问题，仍是如下线性代数方程组：迭代法所要讨论的问题，仍是如下线性代数方程组：nnnnnnnnnnbTaTaTabTaTaTabTaTaTa22112222212111212111(2.43) 简写成：简写成：nibTaijijni, 2, 11(2.44) n迭代法的基本思想是，构造一个由迭代法的基本思想是，构造一个由 组成的向量

54、序列，使其收敛于某个极限向量组成的向量序列，使其收敛于某个极限向量 ，并且，并且 就是方程组就是方程组(2.43)的精确解。根据构造向量序列的方法不的精确解。根据构造向量序列的方法不同，常用的有简单迭代法，同，常用的有简单迭代法，GaussSeidel迭迭代法与超松弛迭代法，下面分别予以介绍。代法与超松弛迭代法，下面分别予以介绍。n21TTT,*2*1,nTTT*2*1,nTTTn1简单迭代法简单迭代法 最简单的迭代法称为简单迭代法，也称同步迭代最简单的迭代法称为简单迭代法，也称同步迭代法，法，Jakobi迭代法。迭代的最终目的是求解方程迭代法。迭代的最终目的是求解方程组组(2.43)中的中的

55、 。 前面已经说到，用有限差分法（包括有限元法）前面已经说到，用有限差分法（包括有限元法）得到的代数方程组能保证系数矩阵得到的代数方程组能保证系数矩阵A对角线上对角线上元素不为零，即元素不为零，即 n21TTT,n21i0aii, n则可将式则可将式(2.43)改写成：改写成： )()()(1n1n11nn1nnnnn232312121222nn121211111TaTabaTTaTaTabaTTaTabaT(2.45) 其中任一方程均可写成：其中任一方程均可写成：niTabaTjijnijjiiii, 2, 1)(11(2.46) n任意给定各节点上的温度值任意给定各节点上的温度值 作为作为

56、解的第零次近似，把它们代入式解的第零次近似，把它们代入式(2.46)的右端，的右端，由此算得的由此算得的),()(n21iT0iniTabaTjijnijjiiii, 2, 1)()0(11) 1 (作为解的第一次近似，把第一次近似得到的解再作为解的第一次近似，把第一次近似得到的解再代入式代入式(2.46)的右端，得到解的第二次近似。一的右端，得到解的第二次近似。一般地讲，在已得到解的第般地讲，在已得到解的第k次近似次近似 后，代入式后，代入式(2.46)右端，得右端，得)(kiTniTabaTkjijnijjiiiki, 2, 1)()(11) 1( 为解的第为解的第k+1次的似。这样得到的

57、序列次的似。这样得到的序列 为方程组的近为方程组的近似解。只要方程组似解。只要方程组(2.43)存在唯一解，则不存在唯一解，则不论零次近似如何选取，当论零次近似如何选取，当 时，此序时，此序列列 必然收敛，且收敛于方程组的必然收敛，且收敛于方程组的解解 。实际计算中。实际计算中k不可能取不可能取 ，但可以说，当但可以说，当k充分大时，序列充分大时，序列 已足够精确地接近方程组的。已足够精确地接近方程组的。 ),(,)()()(210kTTTknk2k1kn21TTT,*2*1,nTTTn21TTT,n通常，对充分大的通常，对充分大的k，其相邻两次迭代解，其相邻两次迭代解 之间的偏差小于之间的偏差小于预先给定的适当小量预先给定的适当小量 （ 大于零），即大于零），即满足满足),(,)()(n21kTTk21k1niTTkiki, 2, 1|)() 1(就结束迭代过程，而取就结束迭代过程，而取 作作为方程组为方程组(2.43)的近似解。的近似解。 ), 2, 1()(niT