计算方法第3章线性方程组数值解法

计算方法北京工业大学应用数理学院杨中华第三章线性方程组的数值解法

线性方程组是应用最为广泛的数学模型，很多复杂问题中都含有线性方程组子问题，因此讨论线性方程组问题的求解很有必要，本章将讨论线性方程组的数值解法。

线性方程组的一般形式：

记：

这里A称为系数矩阵，x称为解向量，b称为右端项。(3.1)(3.2)

求解线性方程组问题的数值方法可分为两类：直接解法和迭代解法。

直接解法是通过有限次初等运算，求得其解，虽然直接解法的推导过程都是无误差的，但是由于计算机的运算都是有舍入误差的，所求的解其实是一个有误差的近似解。

迭代解法则是从某个初始近似解出发，按照一个确定的迭代公式得到一个更好的近似解，反复迭代，直到求得一个满足精度要求的近似解。

本章首先讨论线性方程组的直接解法然后再介绍迭代解法。3.1消去法

1.顺序Gauss消去法

首先回顾一下线性代数中所讲的线性方程组消去法过程，然后归纳出消去法的数值算法，请看如下的例子：

例3.1求解线性方程组

解：求解线性方程组的第一阶段称为消元过程，其方法是第2个方程减去第1个方程的1/2倍，第3个方程减去第1个方程的2倍，得

第3个方程减去第2个方程的7/3倍，得消去法

这一过程就是消元过程，即把方程化为等价的上三角方程(对角线下变为0)。

第一个阶段完成后，进入第二个阶段，称为回代过程，其方法是：先由第3个方程解出

，将

代入第2个方程解出

，再将

和

代入第1个方程解出

也就解出所有的未知量。如下就是所求的解：

归纳以上求解方法，求解线性方程组包括两个过程，消元过程和回代过程。

首先给出回代过程的算法，回代过程其实是一个特殊形式的方程组的求解方法，就是一个上三角线性方程组的求解方法，如：消去法(3.3)

第1步：根据第n个方程解出

，得消去法假如已经求出

，代入第k个方程得：

第2步：根据已求出的

和第n-1个方程求

，得到回代过程就是对

实施这一公式，注意必须从后向前计算方可，所以此过程叫做回代过程。(3.4)

算法3.1

上三角线性方程组的回代算法 0)[初始化]设置上三角方程系数矩阵A，右端项向量b 1)[回代过程]对于

循环 1.1)计算 1.2)对于

循环

1.3)计算 2)[算法结束]消去法

再看Gauss消去法的消元过程，对于线性方程组

的消去法本质上可看作将其增广矩阵用初等行变换化为梯形矩阵的过程。

为清楚的表示每次消元前后系数矩阵和右端项的状态，通常以

表示系数矩阵和右端项，其元素分别记作

，其中上标(k)表示在第k次消元前的状态，其初始增广矩阵为：

首先对第1列消元，也就是对增广矩阵做初等行变换将元素

约化为0。(3.5)消去法

希望消元之后形如：

为此，对

，令消元因子

，将第i各行减去第1行的倍，也就完成了消元过程。之后的论述中，在不引起混淆的情况下，将(3.6*)式略去上标(k)记为：

第2列消元，对

，令

，将第i各行减去第2行的倍，则有(3.6*)(3.6)消去法

再对第3、4、…、n-1列消元，最后得到上三角线性方程组的增广矩阵为：

以上的消元过程可用一个统一公式表示为：

我们可以根据这一公式设计如下算法：(3.7)消去法

通常在求解线性方程的编程时，消元因子

存储于矩阵A对角线以下相应的位置上，因为此时这些位置的元素将被消元为0。

算法3.2

顺序Gauss消去法 0)[初始化]置系数矩阵A，右端项向量b； 1)[消元过程]对于

循环

对于

循环 1.1) 1.2) 1.3) 2)[回代过程]执行算法3.1，即 3)[算法结束]消去法

从前述的算法可知，消元过程能够顺利进行的前提是步骤1.1)中对角线元素

，关于顺序Gauss消去法有如下结论。

定理3.1顺序Gauss消去法能够顺利进行的充要条件是系数矩阵A的顺序主子式

，且

。特别有

证明略。

消去法 2.顺序Gauss消去法的计算量

通常我们在度量一个有限步终止的算法的计算量时，主要考虑整个计算过程中使用乘除法运算的总次数，在以后统计算法的计算量时，我们仅统计算法中乘除运算的次数。 (1)消元过程的计算量

(2)回代过程的计算量 (3)总的计算量

消去法

3.列主元Gauss消去法

前节所讲顺序Gauss消去法，有以下两个问题：在消元过程中当某个元素对角线为零时，算法无法进行，此时方程组并不一定无解；如果某个对角线元素的绝对值非常小，将会引起很大的误差。

因此顺序Gauss消去法不具实用性，但由此方法改进的列主元消去法是最为广泛使用的方法，以下我们通过一个实例体会一下顺序Gauss消去法的问题和稍作改进之后的效果。消去法

例3.2用顺序Gauss消去法求解线性方程组(运算过程中保留4位有效数字)

解：写出方程组的增广矩阵

消元过程：

，第2行减去第1行的倍，得到

回代过程求得计算解为：

，此方程的准确解为

，显然此计算解已经没有任何意义。消去法

如果将两个方程交换一下位置再求解，增广矩阵：

消元过程：

，第2行减去

乘以第1行，得到

回代过程求得解为：

此解与准确解完全相同。消去法

第1种解法产生非常大的误差的原因是

的绝对值太小，做为分母将产生太大的误差而导致最后的解面目全非；用第2种方法是通过交换方程的位置使做分母的

的绝对值尽可能的大，正如第一章误差分析所述，避免绝对值太小的数做分母，可提高计算结果的精度。这一过程称为选主元(所选的对角线元素称为主元)，一般的做法如下。

设进行第k步计算时增广矩阵如下：消去法

为了在第k次消元过程中使消元因子中分母的绝对值尽可能的大，在第k列中取：第k到第n个元素中绝对值最大者

在选定主元之后，考虑以下情况：如果

，则线性方程组的行列式的值为0，其解不定(无解或有无穷多组解，停止计算)；如果r=k，则

的绝对值最大，以

作为主元；如果

则交换矩阵的第r、k行，交换之后第k行对角线元素绝对值变为最大。交换增广矩阵的两行相当于交换两个方程的位置，显然不改变原方程的解。关于

的判断：取充分小的数

，当

时

增加选主元机制的Gauss消去法即为列主元Gauss消去法，列主元Gauss消去法的计算过程包括：选主元，交换主元行，消元。以下我们写出适合编程的Gauss算法。(3.8)消去法

算法3.3列主元Gauss消去法 0)[初始化]置系数矩阵A，右端项向量b，取充分小的 1)[消元过程]对于

循环 1.1)[选主元] 1.2)如果

，停止计算，计算失败 1.3)如果r>k交换

方程

1.4)[消元]对于

循环 1.4.1) 1.4.2) 1.4.3) 2)[回代过程]执行算法3.1，即计算 3)[算法结束]消去法 4.全主元Gauss消去法

全主元Gauss消去法与列主元Gauss消去法的区别在于取主元的方法是：

如果r>k交换增广矩阵的第r、k行，如果q>k交换增广矩阵的第q、k列(本质上是交换未知量

，计算结束时应该交换回来)。

显然全主元消去法取主元的方法更可靠，主元的绝对值会更大，也就更稳定，但计算量稍大。实际计算效果表明，列主元Gauss消去法已经很稳定，全主元法的改进甚微，是没必要的，所以通常使用列主元Gauss消去法。消去法

5.Gauss-Jordan消去法

列主元Gauss消去法在消元过程中是将对角线以下的元素约化为0，将系数矩阵化为上三角矩阵，通过回代得到方程的解；Gauss-Jordan消去法是将系数矩阵约化成单位矩阵，无需回代就直接得到了方程组的解，也称为无回代过程的消去法。

设线性方程组的增广矩阵为：

第1步计算，设

，第1行除以

，将

行减去新的第1行的

倍，则增广矩阵变为如上右端形式消去法

第2步计算，设

，第2行除以

，第

行减去新第2行的

倍，则增广矩阵变为如下右端形式消去法

设第k-1步计算结束时增广矩阵为如下形式

第k步计算，设

，第k行除以

，第

行减去新第k行的

倍，则增广矩阵变为如上右端形式。消去法

完成n-1步计算之后，增广矩阵变为

显然最后一列就是方程组的解。

注意到在求解过程中，总是假设对角线上的主元

，若

算法将无法进行下去，如同列主元消去法，在消元过程中也可用选列主元的策略来避免这种情况发生，除非系数矩阵奇异。消去法

算法3.4求解线性方程组的Gauss-Jordan消去法 0)[初始化]置系数矩阵A，右端项向量b，取充分小的 1)[消元过程]对于

循环1.1)[选主元]

1.2)如果

停止计算，计算失败！1.3)如果r>k交换r、k方程

1.4)[规格化主元行]1.5)[消元]即对于

循环 1.5.1) 1.5.2) 2)[算法结束]消去法

6.Gauss-Jordan消去法求逆矩阵及算法设计

首先通过一个求逆过程来写出相应算法。

例3.3用Gauss-Jordan消去法求下矩阵的逆阵：

解：在矩阵A的右边添加单位矩阵进行初等行变换消去法

此时对应初始状态的单位矩阵就是所求的逆矩阵：

通常编程实现求逆矩阵时并不添加右边的单位矩阵，而是对原始矩阵直接进行变换，在完成计算时原矩阵就变成了逆矩阵。其实我们注意一下求逆的过程，有如下两个特征：

首先、从初始状态到最终状态，增广矩阵的每个状态都有n个列是单位列向量，这也就意味着每次变换，都会将原矩阵的一列变为单位列向量，而将原单位矩阵的一个列向量变为逆矩阵的一列，每次变换时可让这一对列向量占用共同的列，也就是可将右边新产生的一列直接放在左边将要变为单位列向量的位置上。消去法

其次、在求逆过程中如果没有行交换，则增广矩阵中单位列向量的消失与新增，在左右两部分都是按自然顺序递增的。从例3.3可以看出，因交换主元行，右边的单位列向量顺序随主元行而变，因此在完成计算后需要换回原位置，这正说明原矩阵的行交换对应逆矩阵的列交换这一性质。

以下过程是通常编程时的实际计算过程：

整个计算有n-1个消元过程，但是最后一次消元不须选主元，所以主元行标分别是{2,3}，按逆序分别交换列

，则得到正确的逆矩阵，此过程直观性太差，仅作为调试程序时核对之用。消去法

算法3.5

矩阵求逆 0)[初始化]输入矩阵A，并以向量p记录主元行下标(共n-1个

主元行下标)，取充分小的 1)[行变换]对

循环 1.1)[选主元] ， 1.2)如果

，矩阵不可逆，算法失败；否则令

1.3)[交换主元行]如果r>k交换r、k行 1.4)[规格化主元行] 1.5)[消元]对

循环 1.5.1) 1.5.2) 2)[交换列]对

循环

如果

，交换

列 3)[算法结束]消去法

以下是计算机求逆验证矩阵消去法(3.9)(3.10)3.2LU分解

1.消去法的矩阵运算表示形式

从矩阵计算的角度讲Gauss消去法是通过一系列初等行变换将矩阵约化成上三角矩阵；Guass-Jordan消去法是通过一系列初等行变换将矩阵约化成单位矩阵。矩阵的初等行变换本质上是矩阵左乘一个初等矩阵，因此消去法也可以用矩阵乘积形式进行描述。 Gauss消去法的矩阵形式：LU分解

Gauss消去法的矩阵形式：LU分解

Gauss-Jordan消去法的矩阵形式：LU分解

前述两式中的

就是我们消去法的消元因子，如果记：则Gauss消去法的整个过程可以表示为：注意到形如(3.11)式的单位下三角矩阵具有以下两个性质：(1)逆阵仍然是单位下三角阵；(2)任意两个单位下三角阵的乘积仍然是单位下三角阵。所以上式也可以表示为：(3.11)LU分解因此Gauss消去法的本质上是将系数矩阵A进行如下的矩阵分解：其中L是单位下三角矩阵，U是上三角矩阵，此种分解也称为杜里特尔(Doolittle)分解，如果L是下三角矩阵，U是单位上三角矩阵，则称为克洛特(Crout)分解，以下内容讨论的是杜里特尔分解方法。(3.12)LU分解2.LU分解

上述A=LU也可以写成：

实现矩阵的LU分解算法既可以仿照消去法设计算法也可以直接从该式推导出计算公式，以下我们讨论直接LU分解方法。

首先用L的第1行乘以U的各列，则有

由此可以得到U的第1行：再让L的各行乘以U的第1列得到

，由此又可以得到L的第1列：假如我们已经得到了L的前k-1列和U的前k-1行：按如下方法确定U的第k行和L的第k列，将L的第k行与U的第

各列相乘，有得到U的第k行将L的第i(i>k)各行与U的第k列相乘，有得到L的第k列LU分解LU分解则有如下结果：按如此方法，从k=1执行到n-1则完成了矩阵的LU分解。

综上所述，我们得到直接LU分解的计算公式：(3.13)LU分解定理3.2设n阶矩阵A的顺序主子式均不为0，则A存在LU分解，且分解是唯一的。

事实上，如果分解成功，A的第k个顺序主子式的值就是

，因此结论是显然的。本定理表明LU分解可以正常进行的前提是顺序主子式不为0，表现在分解过程中就是所有

。

通常可以将下三角矩阵L和上三角矩阵U放在一个矩阵内，即：(3.14)

例3.4

求矩阵A的LU分解

解：设

先计利用(3.13)式算U的第1行

根据

的计算公式计算L的第1列

当前的分解为LU分解

再根据

的计算公式计算U的第2行

再根据

的计算公式计算L的第2列

当前的分解为

现在仅剩

，根据

的计算公式计算

LU分解

最终得到分解为：

将下三角矩阵L和上三角矩阵U放在一个矩阵内，如：LU分解 3.用LU分解求解线性方程组

考虑线性方程组

如果系数矩阵A存在LU分解A=LU，代入原方程组

令

，则可先求解

然后再求解

即可得到原线性方程组Ax=b的解。此处两个方程组，一个是下三角方程，一个是上三角方程。上三角方程

的求解可以使用算法3.1，以下考虑下三角方程的求解。LU分解

将下三角方程展开

求解上三角线性方程组是按

次序求出，所以称为回代过程，求解下三角线性方程组的次序正好相反是按

求出，所以称为顺代过程，明显的看出计算公式如下：LU分解

算法3.6求解线性方程组的直接LU分解算法 0)[初始化]输入系数矩阵A，右端项b 1)[LU分解]对于

循环 1.1)[计算U的第k行] 1.2)[计算L的第k列] 2)[解下三角方程] 3)[解上三角方程] 4)[算法结束]LU分解

例3.5用LU分解方法求解下列线性方程组

解：先对系数矩阵LU分解得到：

解下三角方称：

解上三角方程：LU分解

前述的LU分解相当于顺序Gauss消去法，如同顺序Gauss消去法不具实用性一样，算法3.6也不具实用性，应该使用更稳定的列主元策略。由线性代数知识可知，交换两行相当于系数矩阵左乘一个交换矩阵，将所有交换矩阵的乘积记为P，就得到如下定理：

定理3.3如果矩阵A非奇异(矩阵可逆)，则存在交换矩阵P，单位下三角阵L和上三角阵U，使得

用列主元实现矩阵的直接LU分解远不如对算法3.3稍作改进更容易，以下我们给出这一改进的算法。LU分解(3.15)

算法3.7

求解线性方程组的直接LU分解算法0)[初始化]设置系数矩阵A，向量p={1,2,…,n}存放主元行下标，取充分小的1)[分解过程]对于

循环1.1)[选主元]1.2)如果

，停止计算，算法失败！1.3)如果r>k交换第r、k行1.4)对于

循环1.4.1)1.4.2)2)[算法结束]LU分解

本算法的几点说明：1)向量p记录主元行的下标，必要时可以生成交换矩阵P，完成矩阵或方程的交换功能，在求解方程时根据主元行下标对右端项作相应的交换；2)交换主元行时，要对n个元素进行交换，这是因为此时矩阵L也要做相应的交换；3)每次循环都要对k行k列以后的元素进行更新，步骤1.4.1)是对

的最后计算结果，对

的计算则隐含在步骤1.4.2)中；4)完成算法之后，下三角矩阵L除主对角上的1之外其他元素置于矩阵A对应的位置上，上三对角矩阵U置于矩阵A的对应位置上。LU分解3.3三对角线性方程的追赶法

本节考虑如下形式的线性方程组

称上述方程组为三对角方程组，其系数矩阵称为三对角矩阵。其矩阵形式为：

其中系数矩阵A是三对角矩阵。在一些实际问题中，例如偏微分方程的差分方程、三次样条函数、数值微分等，都会遇到此类线性方程组，以下我们通过LU分解，导出求解三对角方程更为简便有效的算法，称其为追赶法。(3.16)(3.17)

对三对角矩阵进行求解LU分解：

右端两矩阵相乘我们可以得到：

从而得到以下三对角阵LU分解的递推公式：追赶法(3.18)

在得到三对角阵的LU分解之后，先解下三角方程：

由此得到y的计算公式：

再求解上三角方程：

由此得到x的计算公式：追赶法

算法3.8求解三对角方程的追赶法 0)[初始化]输入数据

，其中A是三对角阵，令

1)[矩阵分解] 2)[解下三角方程] 3)[解上三角方程] 4)[算法结束]追赶法

例3.6

求解三对角方程(计算过程保留小数点以后3位数字)

解：首先对三对角矩阵进行LU分解，结果为：追赶法

先解下三角矩阵

再解上三角矩阵

追赶法

有时还会遇到如下形式称为广义三对角方程组：

对(3.19)式也可以设计出O(n)计算量的算法，仿照Gauss消去法，可用如下方法求解。

首先消元第一列：显然仅需要消元

，在完成对应的两个方程的消元之后会增加

两个位置的非零元素，我们仍然用

表示，成为如下形式：追赶法(3.19)

一般情况，对于

，每列仅两个位置

的元素

应该消元，同时在

的位置上会新增两个非零元素，我们仍记为本次被消元的元素符号

，当完成了

的消元过程时是如下形式

追赶法(3.19)(3.20)

显然只要再完成(3.20)式对

的消元就可以将其变为如下形式

剩下的问题就是一个简单的回代算法。

算法3.8A求解广义三对角方程算法 0)[初始化]输入数据

a,b,c,d 1)对于

执行追赶法(3.21)

2)对于

做 3)对于

4)[回代过程] 5)[算法结束]

显然容易算出该算法的计算量约为11n

次乘除法运算。追赶法3.4平方根法

1.平方根法

本节我们讨论对称正定矩阵的一个更为有效的分解方法—平方根法或Cholesky分解，首先回顾线性代数中对称正定矩阵的有关结论。

定义3.1设矩阵

满足

，则称矩阵A是对称矩阵，如果对任意非零向量

均有

则称A是正定矩阵，如果有

则称A是半正定矩阵。

定理3.5

设A是对称正定矩阵，则有以下结论成立：所有特征值为正；所有顺序主子式为正；如果对称正定矩阵可逆则逆矩阵也是对称正定矩阵。

定理3.6设

，且A的所有顺序主子式均大于0，则存在唯一的单位下三角阵L和对角阵U，使得

证明：因为

的顺序主子式不为零，所以

存在唯一的LU分解A=LU，注意到上三角矩阵U可以分解为：

将其带入

则有

，由对称性得

再根据LU分解的唯一性，则有

得证。平方根法(3.22)注意到如果矩阵A对称正定的，则有

，令将上式可以改写成如下形式：

则有如下结论。

定理3.7

(Cholesky分解)设A是正定对称矩阵，则存在唯一的非奇异下三角阵L，使得

其中L的对角元素均为正数，称为平方根分解或Cholesky(乔莱斯基)分解。平方根法(3.23)

至此我们仅给出了Cholesky分解的存在性条件，以下讨论分解算法。如同LU分解，也考虑直接的分解方法，从下式推导出计算公式：

首先注意到

假设已得到L的第1列到第

列，现要计算第k列，则有平方根法(3.25)

例3.7求对称矩阵的Cholesky分解，并求其行列式的值

解：

最后得到平方根法 2．改进平方根法

平方根算法是稳定的，缺点是算法中有开平方运算。一次开平方运算的CPU时间远远大于乘除运算(开平方运算一般是迭代法实现的)，因此通常我们应该尽量避免使用平方根运算，改进平方根法就是无开方运算的平方根法。

其实，所谓无开方运算的平方根法就是基于

公式，将其展开

以

作为中间矩阵，则有

，直接对

按矩阵乘法规则推导出计算公式为：平方根法

如果求解线性方程组Ax=b，以下只要分别求解下三角方程

和下三角方程

，即用如下公式计算：

在编程实现改进平方根算法时，对角阵D可以存放在矩阵L的对角线上，中间矩阵T可以借用矩阵L的上三角部分，改进平方根算法如下：平方根法(3.26)(3.27)算法3.9

改进平方根矩阵分解算法0)[初始化]输入系数矩阵和右端项

1)[分解]对于

循环1.1)[计算L]对

1.2)[计算D]2)[算法结束]以上算法中没有出现中间矩阵T，矩阵T隐含在矩阵L中。平方根法在编程时，如果希望将下三角矩阵L和中间矩阵T存放在原矩阵A内，只需要将算法3.9中的所有

改为

即可，其它不须作任何修改；如果利用改进平方根算法求解线性方程组，可以在算法3.9的基础上按照(3.27)公式继续计算则可以得到方程组的解，即：平方根法例3.8

用改进平方根算法分解以下矩阵解：按算法3.9步骤1)可求得如下结果：即：平方根法3.5迭代法

以上所述Guass消去法、LU分解求解线性方程组的方法均属于直接方法，还有一类求解线性方程组的方法属于迭代法，迭代法在求解大型稀疏问题、病态方程组等特殊情况时是一个非常有效的方法。

直接法可以用有限次初等运算求出其解，比如Gauss消去法可以用大约

次乘除法运算和

次加减法运算求出方程组的解；而迭代法思想，先给出一个初始近似解，根据一个迭代公式得到一个更好的近似解，通过反复执行迭代过程逐步逼近到精确解，当近似解达到指定的精度时结束计算。

本课程仅介绍Jacobi迭代法和Gauss-Seidel迭代法。

如同求解非线性方程的迭代法，求解线性方程组的迭代法也是先把所给的线性方程组Ax=b化成同解的线性方程组

这里称B为迭代矩阵，给定初始向量

并按迭代格式

进行迭代计算。

如果按上述迭代格式所得到的向量序列

收敛于某一向量x*，则x*就是方程组Ax=b的解，并称此迭代法收敛；否则称为不收敛或发散。迭代法(3.28)(3.29)

1.Jacobi迭代法

设有线性方程组

写成其紧凑形式为

假设系数矩阵非奇异且主对角线元素

，将方程组的非对角线项移至等式的右端并两端除以对角线元素，则方程组的等价形式为：

按此式建立迭代格式：

则称为Jacobi迭代格式，此迭代方法称为Jacobi迭代法。迭代法(3.31)(3.30)

为便于今后的收敛性分析，现将迭代公式改写成矩阵形式。记

则有

，方程组

改写成：

也就有

相应矩阵形式的迭代公式为：

其中

称为Jacobi迭代矩阵，

。迭代法(3.32)

在非线性方程的迭代法中，迭代过程的终止准则是两相邻迭代近似解之差的绝对值小于给定精度

，而对于线性方程组的迭代法，其终止准则是两次相邻迭代近似解之间的某种范数小于给定精度，即

例3.9用Jacobi迭代法求解方程组

解：等价形式的方程组为：迭代法

对应的迭代格式为：

取初始点

，终止精度为

，迭代过程如下表：迭代法迭代序号00.0000.0000.00010.444-0.1251.00020.417-0.8610.90930.918-0.7971.12740.851-1.0590.96651.043-0.9411.05960.954-1.0480.97171.035-0.9701.02780.977-1.0260.981

经21次迭代，收敛到线性方程组的精确解

，收敛速度较慢。迭代法迭代序号91.019-0.9831.014100.987-1.0130.990111.010-0.9901.008120.993-1.0070.994131.005-0.9951.004140.996-1.0040.997151.003-0.9971.002160.998-1.0020.998171.002-0.9981.001180.999-1.0010.999191.001-0.9991.001200.999-1.0010.999211.000-1.0001.000

算法3.10求解线性方程组的Jacobi迭代法 0)[初始化]初始解

，最大迭代次数N，精度

1)[迭代]当k<N时循环执行下列各步骤 1.1)对于

循环 1.2)如果

，则停止计算，

为近似解 1.3)令k=k+1，继续执行下一次循环 2)[算法结束]迭代法

2.Guass-Seidel迭代法

前述的Jacobi迭代的优点是算法简单，缺点是收敛速度较慢。我们稍加观察就会注意到，如果迭代是收敛的，则有理由认为一个新的近似解应该优于前一个近似解，但在迭代过程中并不是充分利用了最新的计算结果。比如在计算

时仍然在使用

，其实此时

已经计算出来，而近似值

应该比

更好，如果用新的结果随时参与新的计算，结果会更好些。Guass-Seidel迭代法就是基于此思想而设计的，因此我们将Jacobi迭代公式(3.31)改写为：

称其为Guass-Seidel迭代格式，相应的迭代法称为Guass-Seidel迭代法。(3.33)迭代法

以下我们给出Gauss-Seidel迭代格式的矩阵形式，沿用前述的矩阵记号，(3.33)式可以写成等价的形式：

迭代格式的矩阵形式为：

这里

称为Gauss-Seidel迭代矩阵，迭代法(3.34)

例3.10请用Gauss-Seidel迭代法求解例3.9中的线性方程组，初始点和终止精度仍按原例所取。

解：首先写出Gauss-Seidel迭代格式为：

如同例3.9，取初始点

，终止精度为

，迭代过程如下表：迭代次数10.444-0.2360.94020.497-0.8371.09730.881-1.0311.04341.016-1.0311.00451.020-1.0080.99661.006-0.9990.99871.000-0.9991.00080.999-1.0001.000迭代法

算法3.11

求解线性方程组的Gauss-Seidel迭代法 0)[初始化]初始解

，最大迭代次数N，精度 1)[迭代]当k<N时循环执行下列各步骤 1.1)对于

循环 1.2)如果

，则停止计算，

为近似解 1.3)令k=k+1，继续下一次循环 2)[算法结束]迭代法

3.超松弛(SOR)迭代

解线性方程超松驰迭代法是目前解大型方程组的一种最常用方法，其本质是对Gauss-Seidel迭代法的一种加速方法，其思想是根据两次相邻迭代的近似解构造一个新的更好的近似解。在Gauss-Seidel迭代过程中，对Gauss-Seidel迭代法进行如下修正：

该迭代法称为逐次超松弛(SOR)迭代法，其中系数

称为松弛因子，当

时称为低松弛迭代，取为

就是Gauss-Seidel迭代，当

称为超松弛迭代。迭代法

例3.11用超松弛迭代法求解方程组。

解：取初始点

，精度均取

。

松弛因子

的迭代12次的如下结果迭代法迭代次数01.0001.0001.00015.2503.813-5.047…………83.0083.993-5.00293.0053.996-5.001103.0033.997-5.001113.0023.998-5.000123.0013.999-5.000

松弛因子

的迭代7次的如下结果迭代法迭代次数01.0001.0001.00016.3133.520-6.65022.6223.959-4.60033.1334.010-5.09742.9574.007-4.97353.0044.003-5.00662.9964.001-4.99873.0004.000-5.000

算法3.12

求解线性方程组的逐次超松弛(SOR)迭代法 0)[初始化]初始解

，松弛因子

，最大迭代次

数N，精度 1)[迭代]当k<N时循环执行下列各步骤 1.1)对于

循环 1.2)如果

，则停止计算，

为近似解 1.3)令k=k+1，继续执行下一次循环 2)[算法结束]迭代法

4.迭代法的收敛性

我们使用任何迭代算法，都要清楚该算法必须具有收敛性、了解收敛速度如何?使用发散的迭代算法或者是收敛速度太慢的迭代算法是没有实际意义的，因此我们要研究一下前述的Jacobi迭代法和Gauss-Seidel迭代法的收敛性、收敛速度问题。

定义3.2设

是

上的任一向量范数。若存在

使

则称序列

收敛于x*。

无论在任何一种范数下，序列

收敛于x*，根据向量范数的等价性，必可推出在其它范数下也收敛于x*，即收敛性与范数的选择无关。迭代法

定理3.8

(迭代法收敛的充分条件)对于迭代形式的线性方程组

，如果

则对于任意的初始向量

及任意的右端向量f，解此线性程组的迭代格式

是收敛的，且

证明略。

本定理说明迭代法的收敛性与右端项无关，也与所取初始点无关，但初始点的选取与收敛的快慢有关。(3.36)(3.37)迭代法

例3.12设线性方程组的系数矩阵为

试判断使用Jacobi迭代法和Gauss-Seidel迭代法求解该方程组的收敛性。

解：分解系数矩阵为

对于Jacobi迭代法，有迭代法

对迭代矩阵

取1-范数，则有

根据定理3.8可知Jacobi迭代法收敛。对于Gauss-Seidel迭代法，有

对迭代矩阵

取∞-范数，则有

同样是根据定理3.8可得Gauss-Seidel迭代法收敛。迭代法

定理3.9

(迭代法收敛的充要条件)对于迭代形式的线性方程组

，迭代法收敛的充要条件是谱半径小于1，即

证明：如果方程组的解为x*，则应有

显然

的充要条件是

，根据定理1.7可知

，极限成立的充要条件是

，证毕。

与定理3.8相比，定理3.9的条件不容易验证，以下我们再给出一个容易验证的收敛条件，首先给出如下定义。迭代法

定义3.3如果矩阵

满足

则称矩阵A为对角占优矩阵，如果严格不等式成立，则称A是严格对角占优的。

定理3.10

在求解线性方程组Ax=b的迭代法中：

如果矩阵A严格对角占优，则Jaocbi迭代法和Gauss-Seidel迭代法均是收敛的，当松弛因子满足

时，SOR迭代法是收敛的；

如果矩阵A是对称正定的，则Gauss-Seidel迭代法是收敛的，当松弛因子

时，SOR迭代法是收敛的。迭代法3.6线性方程组的误差分析

在本章的前几节，我们给出了一些线性方程组的求解方法，一个还没有涉及到而且很重要的问题是，这些方法所得到的解可靠吗？如果是近似解其近似程度如何？在何种情况下所求的解是可靠的？何种方法求出的解更为可靠？这就是本节所要研究的问题。

1.病态线性方程组

在求解出方程组Ax=b的计算解x*后，直观上会认为可用如下方法度量其计算解的准确性，计算其解的残向量：

如果

的值很小，可能认为x*比较精确。但有时并非如此，看下述线性方程组：

此方程组的精确解是

，但将

带入方程组计算残量的∞-范数值为

此时残向量确实已经很小，但是

无论如何也不可以看作是

的近似解。线性方程组的误差分析(3.38)

其实可以验证

是下列方程组的精确解

而这两个线性方程组的差别仅仅是元素

相差

，这种差别可以当做

有一个微小的扰动。

定义3.4

如果线性方程组Ax=b的系数矩阵或右端项有一个微小的扰动将导致扰动前后两方程组的解有一个很大的变化，我们称此线性方程组是病态的，其系数矩阵A称为病态矩阵；反之称线性方程组是良态的。线性方程组的误差分析

众所周知，计算机求解线性方程组的每步计算都会产生舍入误差，而当前的舍入误差又会传播到之后的求解过程中去，导致误差不断的积累，而这些计算过程中的误差及误差积累势必影响到最终计算解误差的大小。目前对这类误差的分析采用一种称之为向后误差分析的方法，其思想是在误差分析过程中，把舍入误差导致最终解的误差看作是初始数据有一个扰动，把有舍入误差的计算过程看作是精确的运算。

之后的讨论就是基于如此思想的误差分析方法，为了度量方程组的病态或良性程度，首先给出关于矩阵条件数的概念。线性方程组的误差分析2.矩阵的条件数

定义3.5

设

是非奇异矩阵，对于某一矩阵范数，称

为矩阵A求解线性方程组的条件数。

特别如果取矩阵的2-范数，则可记作：

如果矩阵A是对称矩阵，

分别是A的最大、最小特征值，则有线性方程组的误差分析(3.39)

考虑之前讨论的病态方程组，系数矩阵为

而

的两个特征值分别为

则其系数矩阵A在2-范数下的条件数为

显然该条件数相比较于系数矩阵的元素已经非常大了。

通常我们使用条件数来度量线性方程组的病态程度，其中的矩阵范数可取为2-范数、1-范数、∞-范数中的任一范数，条件数愈大方程组的病态程度就愈