Download presentation
Presentation is loading. Please wait.
1
第五章 解线性方程组的直接方法 §5.1 引 言 线性方程组: (5.1) 1 结束
2
在自然科学和工程技术中,很多问题归结为解线性方程组
在自然科学和工程技术中,很多问题归结为解线性方程组.有的问题的数学模型中虽不直接表现为含线性方程组,但它的数值解法中将问题“离散化”或“线性化”为线性方程组.因此线性方程组的求解是数值分析课程中最基本的内容之一. 常记为矩阵形式 Ax=b (5.2) 此时A是一个n×n方阵,x 和 b是n维列向量. 根据线性代数知识若|A | ≠0,(5.2)的解存在且唯一. 2 结束
3
关于线性方程组的解法一般分为两大类,一类是直接法,即经过有限次的算术运算,可以求得(5. 1)的精确解(假定计算过程没有舍入误差)
关于线性方程组的解法一般分为两大类,一类是直接法,即经过有限次的算术运算,可以求得(5.1)的精确解(假定计算过程没有舍入误差).如线性代数课程中提到的克莱姆算法就是一种直接法.但该法对高阶方程组计算量太大,不是一种实用的算法〔1〕.实用的直接法中具有代表性的算法是高斯消元法,其它算法都是它的变形和应用. 另一类是迭代法,它将(5.1)变形为某种迭代公式,给出初 始解 x0 ,用迭代公式得到近似解的序列{xk},k=0,1,2, ,在一定的条件下 xk→x* (精确解). 3 结束
4
§5.2 高斯(Gauss)消元法 迭代法显然有一个收敛条件和收敛速度问题. 这两种解法都有广泛的应用,我们将分别讨论,本章介绍 直接法.
高斯消元法是一种古老的方法.我们在中学学过消元法,高斯消元法就是它的标准化的、适合在计算机上自动计算的一种方法。 4 结束
5
5.2.1 高斯消元法的基本思想 例1 解方程组 (5.3) (5.4) (5.5)
高斯消元法的基本思想 例1 解方程组 (5.3) (5.4) (5.5) 第一步,将(5.3)乘-2加到(5.4);(5.3)乘-1加到(5.5), 得到 (5.3) (5.6) (5.7) 5 结束
6
第二步,将(5.6)乘-2/3加到(5.7),得到 (5.3) (5.6) (5.8) 回代:解(5.8)得x3,将x3代入(5.6)得x2,将x2, x3代入(5.3)得x1,得到解 x*=(2,1,-1)T 容易看出第一步和第二步相当于增广矩阵[A:b]在作行变换,用ri表示增广阵[A:b]的第i行: 6 结束
7
7 结束
8
5.2.2 高斯消元法公式 由此看出上述过程是逐次消去未知数的系数,将Ax=b化为等价的三角形方程组,然后回代解之,这就是高斯消元法.
高斯消元法公式 记 Ax=b 为 A(1)x=b(1) , A(1) 和 b(1) 的元素记为 和 , i , j=1,2,,n .第一次消元,目的是消掉第二个方程到第n个方程中的 x1 项,得到 A(2)x=b(2),这个过程须假定 ≠0. 8 结束
9
在[A(1):b(1)]中,红方框中的元素是要化为0的部分;[A(2):b(2)]中,红方框中的元素全部已发生变化,故上标由(1)改(2),计算公式为:
9 结束
10
(i=2,3,,n) (i,j=2,3,,n) (i=2,3,,n) (i=2,3,,n) 第 k 次消元 (1≤k≤n-1) 设第k-1次消元已完成,且 ≠0,此时增广矩阵如下: 10 结束
11
本次消元的目的是对框内部分作类似第一次消元的处理,消掉第k+1个方程到第n个方程中的xk项,即把 到 化为零.计算公式如下:
11 结束
12
(i=k+1,,n) (i,j=k+1,,n) (i=k+1,,n) (i=k+1,,n) 只要 ≠0,(k=1,2,,n-1)消元过程就可以进行下去.当k=n-1时,消元过程完成,得: 12 结束
13
它的方阵部分A(n)是一个上三角形矩阵,它对应的方程组是一个上三角形方程组,只要 ≠0,就可以回代求解,公式为
(i=n-1,n-2,,1) 综合以上讨论,高斯消元法解线性方程的公式为: 13 结束
14
1消元 ①令 (i,j=1,2,…,n) ②对k=1到n-1,若 ≠0 ,进行: (i=k+1,k+2,,n) (i=k+1,k+2,,n) (5.9) (i,j=k+1,k+2,,n) (i=k+1,k+2,,n) 14 结束
15
2回代,若 ≠0 (i=n-1,n-2,,1) (5.10) 高斯消元法的条件 以上过程中,消元过程要求 ≠0 (i=1,2,,n-1),回代过程则进一步要求 ≠0,但就方程组Ax=b讲, 是否等于0是无法事先看出的. 15 结束
16
注意A的顺序主子式Di (i=1,2,,n)在消元过程中不变.这是因为消元所作的变换是“将某行的若干倍加到另一行”
上,据线性代数知识,此类变换不改变行列式的值.若高斯消元过程已进行了k-1步(此时当然应有 ≠0,i≤k-1),这时计算 A(k) 的顺序主子式: 16 结束
17
有递推公式 (i=2,3,,k) 显然, , 可知,消元过程能进行到底的充要 条件是Di≠0 ,(i=1,2,,n-1),若要回代过程也能完成,还应加 上Dn=|A|≠0,综合上述有: 定理5.1 高斯消元法消元过程能进行到底的充要条件是系数阵A的1到n-1阶顺序主子式不为零;Ax=b能用高斯消元法解的充要条件是A的各阶顺序主子式不为零. 17 结束
18
高斯消元法的计算量估计 消元过程的工作量,参看公式(5.9),k是消元次数,k=1,2,,n-1,第k步消元时,计算 lik(i=k+1,,n)需要 n-k次除法;计算 (i,j=k+1, ,n)需要 (n-k)2 次乘法及 (n-k)2次加减法;计算 需要 n-k次乘法及 n-k次减法,合计: 乘除法次数: 18 结束
19
加减法的次数: 回代过程的工作量,参见公式(5.10),求 xk需 n-k 次加减法, n-k次乘法和1次除法,合计为 乘除法次数: 加减法次数: 19 结束
20
总的运算次数为 乘除法 (当n较大时) 加减法 (当n较大时) 一般讲乘除法的运算比加减法占用机时多得多,往往只统计乘除法次数而称高斯消元法的运算量为 次. 20 结束
21
§5.3 高斯主元素消元法 在上节的算法中,消元时可能出现 =0的情况,高斯消元法将无法继续;即使 ≠0,但 <<1,此时用它作除数,也会导致其它元素数量级严重增加,带来舍入误差的扩散,使解严重失真. 例2 线性方程组 准确解是(1,1)T. 21 结束
22
现设我们的计算机为四位浮点数,方程组输入计算机后成为
(5.11) 用高斯消元法: l12=0.2×106 , r2=r2-l12×r1 回代:x2= ×10=1, x1=0, 解严重失真. 22 结束
23
若将 r1和 r2交换. 消元,l12=0.5×10-5,r2=r2-l12×r1 回代 , x1=0.1000×10 , x2=0.1000×10 ,得到准确解. 23 结束
24
5.3.1 列主元消元法 从上例中可以看出,对方程组作简单的行交换有时会显著 改善解的精度.在实际使用高斯消元法时,常结合使用 “选
改善解的精度.在实际使用高斯消元法时,常结合使用 “选 主元”技术以避免零主元或“小主元”出现,从而保证高斯消 元法能进行或保证解的数值稳定性. 列主元消元法 设已用列主元消元法完成Ax=b的第k-1次消元(1≤k≤n-1), 此时方程组Ax=b→A(k)x=b(k),有如下形式: 24 结束
25
(5.12) 进行第k次消元前,先进行①、②两个步骤: ①在方框内的一列内选出绝对值最大者,即 ,确定ik.若 =0,则必有方框内的元素 25 结束
26
全为零,此时易证|A|=|A(k)|=0,即方程组Ax=b无确定解,应给出信息退出计算.
②若 ≠0,则交换ik行和k行元素即 (k j n) 然后用高斯消元法进行消元. 这样从k=1做到k=n-1,就完成了消元过程,只要|A|≠0,列主元高斯消元法必可进行下去. 26 结束
27
5.3.2 全主元消元法 在(5.12)中,若每次选主元不局限在方框列中,而在整个主子矩阵 中选取,便称为全主元高斯消元法,
全主元消元法 在(5.12)中,若每次选主元不局限在方框列中,而在整个主子矩阵 中选取,便称为全主元高斯消元法, 此时增加的步骤为: ① ,确定ik , jk,若 =0,给出 |A|=0的信息,退出计算. 27 结束
28
②作如下行交换和列交换 行交换: (k j n) 列交换: (k i n) 值得注意的是,在全主元的消元法中,由于进行了列交换 ,x 各分量的顺序已被打乱.因此必须在每次列交换的同时, 让机器“记住”作了一次怎样的交换,在回代解出后将x各 分量换回原来相应的位置,这样增加了程序设计的复杂性. 此外作①步比较大小时,全主元消元法将耗用更多的机时. 28 结束
29
§5.4 高斯-若当(Gauss-Jordan)消元法
但全主元消元法比列主元消元法数值稳定性更好一些.实际应用中,这两种选主元技术都在使用. 选主元素的高斯消元法是一种实用的算法,它可以应用于任意的方程组Ax=b,只要|A|≠0. §5.4 高斯-若当(Gauss-Jordan)消元法 很多问题都需要计算方阵的逆阵.线性代数中有多种求逆方法,本节讨论求A-1的数值方法.A非奇异,求A-1 29 结束
30
5.4.1 高斯-若当消元法 设X=A-1,AX=I,将X和I列分块,得n个线性方程组: Axi=ei , i=1,2,,n.
解这n个线性方程组就求得A-1,原则上,所有解的方程组的方法都能求A-1.实际计算中使用较多的是高斯-若当消元法. 高斯-若当消元法 高斯-若当消元法是高斯消元法的一种变形.高斯消元法是消去对角元下方的元素.若同时消去对角元上方和下方的元素,而且将对角元化为1,就是高斯-若当消元法. 设高斯-若当消元法已完成k-1步,于是Ax=b化为等价方程组A(k)x=b(k),增广阵 30 结束
31
(5.13) 在第k步计算时,考虑将第k行上下的第k列元素都化为零,且akk化为1.对k=1,2,,n 1按列选主元,确定ik使 31 结束
32
2换行,交换增广阵第k行和第ik行 (j=k,k+1,,n), 3计算乘数 mik= - aik/akk(i=1,2,,n, i≠k) mkk=1/akk (5.14) 4消元 aij=aij+mikakj(i=1,2,,n,且i≠k,j=k+1,,n) bi=bi+mikbk (i=1,2,,n,且i≠k) (5.15) 32 结束
33
5主行计算 akj=akj×mkk (j=k,k+1,,n) bk=bk×mkk (5.16) 当k=n时, 显然xi=bi,i=1,2,,n 就是 Ax=b的解. 33 结束
34
高斯-若当消元法的消元过程比高斯消元法略复杂,但省去了回代过程,它的计算量约为n3/2,大于高斯消元法.也称为无回代的高斯消元法.
求方阵的逆 高斯-若当消元法解方程组并不比高斯消元法优越,但用于矩阵求逆是适宜的,实际上它是初等变换方法求逆的一种规范化算法: 例3 34 结束
35
35 结束
36
所以 将高斯-若当消元法算法略加改造便成了求逆的算法. AX=I 对k=1,2,,n 1按列选主元,确定ik 2换行,交换第k行和ik行 (j=k,k+1,,n) (j=1,2,,n) 3计算乘数 mik=-aik/akk , (i=1,2,, i≠k) mkk=1/akk, (5.16) 36 结束
37
4消元 aij=aij+mik akj,(i=1,2,…,n,且i≠k,j=k+1,…,n) xij=xij+mik xkj,(i=1,2,…,n,且i≠k,j=1,2,…,n) (5.18) 5主行计算 akj=akj×mkk (j=k,k+1,…,n) xkj=xkj×mkk (j=1,2,…,n) (5.19) A-1=(xij)n×n 应该注意到,在上述算法中都加进了选列主元的步骤,它可以保证在方阵A非奇异的情况下都能求得A-1.求逆算法中不宜使用全面选主元,否则会改变逆阵元素的排列,增加 程序设计的复杂性. 37 结束
38
§5.5 矩阵的 LU分解 Ax=b是线性方程组,A是n×n方阵,并设A的各阶顺序主
子式不为零。令 A(1)=A,当高斯消元法进行第一步后,相当于 用一个初等矩阵左乘A(1).不难看出,这个初等矩阵为 其中li1 (i=2,3,…,n) 由式(5.9)确定,即 38 结束
39
同样第k步消元有 进行n-1步后,得到A(n) ,记U=A(n),显然U的下三角部分全化为零元素,它是一个上三角阵。整个消元过程可表达如下: 已知U是上三角矩阵,下面讨论L的性态,首先指出 39 结束
40
其次指出 证明见习题。 40 结束
41
L相当是由各L-1k (k=1,2,…,n-1)的所有左下元素拼凑后加上对角元1而得.
L是下三角矩阵,且所有的对角元为1,故称为单位下三角阵.这样,A=LU称为A的LU分解,其中L是单位下三角阵,U是上三角阵. 定理5.2 矩阵 An×n,只要A的各阶顺序主子式非零,则A可以分解为一个单位下三角阵L和一个上三角阵U的乘积,即A=LU,且这种分解是唯一的. 当A进行LU分解后,Ax=b就容易解了. 即Ax=b等价于: 这样,Ax=b分解为解两个三角形方程组,三角形方程组是极易求解的.(5.24)即是高斯消元法的矩阵表述. 41 结束
42
5.5.2 直接LU分解 A的LU分解可以用高斯消元法完成,但也可以用矩阵乘法原理推出另一种方法,结果是完全一致的.设: 由矩阵乘法公式可推出直接LU分解算法,也称杜利特(Doolittle) 算法,现总结如下: 1.矩阵分解:A=LU 对k=1,2,…,n 42 结束
43
2.解 L y = b 3.解 U x = y 杜利特算法实际上就是高斯消元法的另一种形式.它的计算量与高斯消元法一样.但它不是逐次对A进行变换,而是一次性地算出L和U的元素.L和U的元素算出后,不必另辟存贮单元存放,可直接存放在A的对应元素的位置,节省存贮单元,因此也称为紧凑格式法. 43 结束
44
这样我们可以在A上直接修改计算,作如下的操作:
杜利特算法可用文字描述如下: 将A的元素划分为形如“┏”的n框, 如A的第一行和第一列元素为第一框,紧靠第一框内部的第二行和第二列元素为第二框,依此类推,ann一个元素为第n框.为便于描述我们称原A矩阵中元素为a元素,计算后每框中的每行元素为u元素(包含对角元素),每列元素为l元素(不包含对角元素). 这样我们可以在A上直接修改计算,作如下的操作: (1)第1框:u元素等于对应的a元素;l元素等于对应的a元素除以本框对角元素; (2)第2到n框:u元素等于对应的a元素减去已经算过的各框中同行同列元素的乘积;l元素除进行以上操作外,还要除以本框对角元素 44 结束
45
得到的矩阵由L和U的元素拼合而成,它的上三角部分(包含 对角元素)即是LU分解后的U矩阵;它的下三角部分(不包含对角
也可以用增广矩阵[A:b]进行以上操作,这时每框要多算一个 u 元 素,结果矩阵中的最后一列就是 Ly=b 的解,回代时只要解 Ux=y 就 行了. 当用手工计算小型线性方程组的解时,用此方法比较方便. 45 结束
46
5.5.3 方阵行列式求法 在实际问题中,有时会遇到求方阵的行列式.在线性代数中讲到的行列式定义算法和展开算法均不适用于阶数较高的行列式.而LU分解是计算行列式的十分方便和适用的算法. A=LU 即只要将U阵的对角元相乘就可得A的行列式. 为避免计算中断,还应加上选主元素过程,此时每做一次行交换(或列交换),行列式要改变一次符号,有: 46 结束
47
其中 p 是进行的行交换次数(对列主元消元法而言),或是行交换和列交换次数的总和(对全主元消元法而言)
其中 p 是进行的行交换次数(对列主元消元法而言),或是行交换和列交换次数的总和(对全主元消元法而言).若选不出非零的主元素,则必有 det A=0. 5.5.4 克劳特(Crout)分解 将LU分解换一个提法:要求L为一般下三角阵,U为单位上三角阵,就是克劳特分解. 实际上将A转置为AT,进行LU分解 AT=LU 则 A=UTLT 47 结束
48
UT显然是一般下三角阵,记 有 LT显然是单位上三角阵,记 这就是克劳特分解。 显然当A的各阶顺序主子式非零时,它是存在且唯一的.克劳特分解对应的解法称克劳特算法,也可用于解线性方程组。它的特点是在回代时不做除法.它在下文的追赶法中有应用.限于篇幅,公式从略。 48 结束
49
§5.6 平方根法 实际问题中Ax=b,A若是对称正定矩阵,则高斯消元法简化为 平方根法或改进的平方根法. 5.6.1 矩阵的LDU分解
§5.6 平方根法 实际问题中Ax=b,A若是对称正定矩阵,则高斯消元法简化为 平方根法或改进的平方根法. 矩阵的LDU分解 定理5.3 若A的各阶顺序主子式非零,则A可以分解为A=LDU,其中 L 是单位下三角阵, U是单位上三角阵, D是对角阵, 且这 种分解是唯一的. 证: 由条件有A=LU,由分解过程知 取 49 结束
50
5.6.2 对称正定矩阵的乔累斯基(Cholesky)分解
记 记 记 于是 5.6.2 对称正定矩阵的乔累斯基(Cholesky)分解 定理5.4 设 A 对称正定,则存在三角分解 A=LLT,L是非奇异下三角矩阵 ,且当限定 L 的对角元为正时,这种分解是唯一的. 50 结束
51
这是一个二次型,由A 对称正定知 di > 0 (i=1,2, … ,n).
证: 因为 A 对称正定, 所以 A 各阶顺序主子式为正. 所以有 A=LDU, AT=UTDLT=A=LDU. 由LDU分解的唯一性知: UT=L, U=LT 所以 A=LDLT. 下面证 D 的对角元为正数, 因为 |L|=1≠0,所以 LTyi=ei 必有非零解 yi ≠ 0, yiTAyi=yiTLDLTyi=(LTyi)TD(LTyi)=eiTDei=di 这是一个二次型,由A 对称正定知 di > 0 (i=1,2, … ,n). 由证明过程容易看出 是唯一的,形如 A=LLT 的分解称为 对称正定矩阵的乔累斯基分解。 51 结束
52
由矩阵乘法原理容易推出L的元素 lij 的算法:
5.6.3 平方根法和改进的平方根法 对称正定矩阵的A=LLT分解对应于解对称正定方程组Ax=b 的平方根法. 由矩阵乘法原理容易推出L的元素 lij 的算法: 对 j=1,2,…,n,计算 Ax=b可化为 52 结束
53
改进的平方根法计算量仍约为n3/6,但回避了开平方运算.
解法是: 这就是平方根法,它适用于系数阵是对称正定矩阵的方程组.它的运算量以乘除法计是n3/6左右,只是高斯消元法的一半,显然这是因为只计算L不算U的缘故. 平方根法不用考虑选主元,这也是它的优点.它的缺点是要计算n次开平方,为避免开平方运算,发展了平方根法的改进形式,它对应于A=LDLT分解.(方法略) 改进的平方根法计算量仍约为n3/6,但回避了开平方运算. 53 结束
54
§ 5.7 追赶法 在很多情况下,如三次样条插值,常微分方程的边值问题等都归结为求解系数矩阵为对角占优的三对角方程组 Ax=f,即:
其中|i-j|>1时,aij=0,且满足如下的对角占优条件: (1)|b1|>|c1|>0,|bn|>|an|>0 (2)|bi|≥|ai|+|ci|, aici≠0, i=2,3,…,n-1. 54 结束
55
用矩阵乘法比较和Crout 分解可推导总结算法步骤如下:
1.分解计算 2.解Ly=f 3.解Ux=y 55 结束
56
实际计算中Ax=f的阶数往往很高,应注意A的存贮技术. 已知数据只用4个一维数组就可存完
实际计算中Ax=f的阶数往往很高,应注意A的存贮技术.已知数据只用4个一维数组就可存完.即{ai},{bi},{ci},{fi}各占一个一维数组, { i}和{i}可存放在{bi},{ci}的位置,{yi} 和{xi}则可放在{fi}的位置,整个运算可在4个一维数组中运行.追赶法的计算量很小,只是5(n-1)次乘除法.追赶法的计算也不要选主元素. 56 结束
57
§5.8 向量和矩阵的范数 生一个问题,即如何判断向量 x 的“大小”,对矩阵也有类似的问 题. 本节介绍n维向量和 n×n 矩阵的范数.
在分析方程组的解的误差及下章中迭代法的收敛时,常产 生一个问题,即如何判断向量 x 的“大小”,对矩阵也有类似的问 题. 本节介绍n维向量和 n×n 矩阵的范数. 2.8.1 向量范数 定义2.1 x 和 y 是 Rn 中的任意向量 , 向量范数‖•‖是定义 在Rn上的实值函数,它满足: (1) ‖ x ‖≥0,并且,当且仅当x=0时, ‖ x ‖=0; (2) ‖k x ‖=|k| ‖ x ‖,k是一个实数; (3) ‖ x + y ‖≤ ‖ x ‖+ ‖ y ‖ 57 结束
58
常使用的向量范数有三种,设 x=(x1,x2,…,xn)T
容易看出,实数的绝对值,复数的模,三维向量的模都满足以上三 条, n 维向量的范数概念是它们的自然推广. 常使用的向量范数有三种,设 x=(x1,x2,…,xn)T 容易验证,它们都满足三个条件. 58 结束
59
例5 x=(1,0.5,0,-0.3)T, 求 解: 矩阵范数 从向量范数出发,可以定义矩阵的范数. 定义5.2 设A是n×n矩阵,定义 为矩阵A的范数. 59 结束
60
这样定义的范数有如下性质: 60 结束 (1) ‖A‖≥0,并且,当且仅当A是零矩阵时, ‖A‖ =0 (2) ‖kA‖= |k|× ‖A‖
(3) 两个同阶方阵A,B,有‖A+B‖≤‖A‖+‖B‖ (4) A是n×n矩阵,x是n维向量,有‖Ax‖≤‖A‖×‖x‖ (5) A,B都是n×n矩阵,有‖AB‖≤‖A‖×‖B‖ 矩阵范数最常用的有以下三种: (λ1 是ATA的最大特征值) 60 结束
61
它们分别与向量的三种范数对应,即用一种向量范数可定义相应的矩阵范数.
定理5.5 Rn 空间上的范数等价.即对任意给定的两种范数 ‖• ‖α ‖• ‖β 有下列关系: m×‖• ‖α ≤ ‖• ‖β ≤ M ×‖• ‖α 其中的m,M是正的常数, ‖• ‖α 表示向量或矩 阵的α范数. 从以上定理看出,当向量或矩阵的任一种范数趋于零时,其它各种范数也趋于零.因此讨论向量和矩阵序列的收敛性时,可不指明使用的何种范数;证明时,也只要就某一种范数证明就行了. 有了向量和矩阵范数的概念,就可以定义向量和矩阵序列的收敛. 61 结束
62
定义5.3 如果 称向量序列x (k)收敛于向量x。
定义5.4 如果 称方阵序列A (k)收敛于方阵A。 定理5.6 向量序列 x (k)收敛于向量 x的充要条件是: 定理5.7 方阵序列A (k)收敛于方阵A的充要条件是: 谱半径 定义5.5 设 n 阶方阵A的特征值为λj (j=1,2,…,n), 则 称为A的谱半径. 62 结束
63
定理5.8 方阵谱半径和方阵范数有如下关系: 证: 设λi是A的任一特征值,xi为对应的特征向量 Axi= λixi 两边取范数,用矩阵性质(4)有 |λi|×‖x i‖≤‖A‖×‖x i‖ 因为 x i ≠ 0,所以‖x i‖>0 所以|λi |≤‖A‖ i=1,2,…,n, 所以ρ(A)=max|λi |≤‖A‖. 定理5.9 设A是n×n阶矩阵,A的各次幂组成的矩阵序列 I,A,A2, …,Ak, … 收敛于零,即 的充要 条件是ρ(A)< 证明从略. 63 结束
64
例6 求‖A‖1, ‖A‖2 , ‖A‖∞ , ρ(A). 解: 显然‖A‖1 =4,‖A‖∞=4 64 结束
65
这里,我们指出,对于实对称矩阵 A , 有 65 结束
66
条件数及病态方程组 线性方程组 Ax=b 的解是由系数阵A及右端向量b决定的.由实际问题中得到的方程组中,A的元素和b的分量,总不可避免地带有误差,因此也必然对解向量x产生影响. 这就提出一个问题:当A有误差ΔA,b有误差Δb时,解向量x有多大误差?即当A和b有微小变化时,x的变化有多大? 若A和b的微小变化,也只导致x的微小变化,则称此问题是“良态”的;反之,若A和b的微小变化会导致x的很大变化,则称此问题为“病态”问题.在以下的讨论中,设A非奇异,b≠0,所以||x||≠0. 1.设Ax=b中仅b向量有误差Δb ,对应的解x发生误差Δx ,即: 66 结束
67
即x的相对误差小于等于b的相对误差的‖A‖×‖A-1‖倍.
注意到Ax=b,所以AΔx=Δb,若A非奇异,有 Δx=A-1Δb. 所以 ‖Δx‖≤‖A-1‖× ‖Δb‖ . 又因为‖b‖=‖Ax‖≤‖A‖×‖x‖, 所以 两式相除,有 即x的相对误差小于等于b的相对误差的‖A‖×‖A-1‖倍. 67 结束
68
3.由以上分析不难看出,当Δb和ΔA一定时,‖A‖×‖A-1‖的大 小,决定了x 的相对误差限. ‖A‖×‖A-1‖越大时,x 可能产生
2. A有误差ΔA ,b无误差, 此时有类似的结论。 3.由以上分析不难看出,当Δb和ΔA一定时,‖A‖×‖A-1‖的大 小,决定了x 的相对误差限. ‖A‖×‖A-1‖越大时,x 可能产生 的相对误差越大,即问题的“病态”程度越严重.同时, 我们看出 ,Ax=b 的“病态”程度,只与A的元素有关,而与b的分量是无关的. 为此,我们有: 定义5.6 若n×n方阵A非奇异,则称‖A‖×‖A-1‖为A的条件数, 记为 Cond(A)=‖A‖×‖A-1‖ 由于选用的范数不同,条件数也不同,在有必要时,可记为 Condp(A)=‖A‖p×‖A-1‖p (p=1,2, ∞ ). 68 结束
69
由于1=‖I‖=‖AA-1‖≤‖A‖×‖A-1‖=cond(A),可知
”.条件数越小,方程组的状态越好,条件数很大时,称方程组为病 态方程组.但多大的条件数才算病态则要视具体问题而定,病态 的说法只是相对而言. 条件数的计算是困难的,这首先在于要算A-1,而求A-1比解Ax=b的工作量还大,当A确实病态时,A-1也求不准确;其次要求范数,特别是求‖A‖2, ‖A-1‖2又十分困难,因此实际工作中一般不先去判断方程组的病态.但是必须明白,在解决实际问题的全过程中,发现结果有问题,同时数学模型中有线性方程组出现,则方程组的病态可能是出问题的环节之一. 69 结束
70
病态方程组无论选用什么方法去解,都不能根本解决原始误差的扩大,即使采用全主元消去法也不行
病态方程组无论选用什么方法去解,都不能根本解决原始误差的扩大,即使采用全主元消去法也不行.可以试用加大计算机字长,比如用双精度字长计算,或可使问题相对得到解决.如仍不行,则最好考虑修改数学模型,避开病态方程组. 如第七章中提到的拟合问题中出现的正规方程组,就往往呈现病态,此时解决问题的方法之一是避开正规方程组,采用正交多项式拟合的方法,尽管后者比前者在理论上和实际计算中都复杂得多. 70 结束
71
例6 n阶希尔伯特矩阵 当n较大时,是有名的病态矩阵.在函数逼近时,有时得到方程组Hnx=b,当n稍大,用任何方法都难以解出理想的结果.考查它的条件数表2-1可见Hn在n稍大时,即呈严重病态. 例7 求cond2(H2),cond1(H2)和cond∞(H2). 71 结束
72
72 结束
Similar presentations