第五章解线方程组的直接方法
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
第五章解线性方程组的直接方法
⏹预备知识
⏹消元法
⏹矩阵分解法
⏹追赶法
⏹误差分析
线性代数是数值计算方法的基础,学习它对数值计算方法其它内容的学习会有很大的帮助。
无论是插值公式的建立,还是微分方程的离散格式的构造,其基本思想都是转化为代数问题来处理,即归结为解线性方程组。
MATLAB的强大功能是建立在矩阵和向量运算基础上的,线性代数的学习也可以大大提高对MATLAB的掌握程度。
线性方程组的基本解法:
直接解法:经过有限步算术运算,在不考虑舍入误差的情况下求得方程组的精确解;
迭代解法:用某种极限过程逐步逼近方程组的精确解。
5.1 预备知识: 矩阵和向量及线性方程组的解
方阵:m=n 的矩阵;
零矩阵:所有元素都为0的矩阵。
在MATLAB中零矩阵由zeros 命令定义。
如A=zeros(m,n)定义一个m×n 零矩阵,n×n 零矩阵可以用命令A=zeros(n)定义。
单位矩阵:所有对角元为1而其余元素均为0的方阵。
单位矩阵记为I。
在MATLAB 中单位矩阵由eye命令定义。
如A=eye(n)定义一个n阶单位矩阵。
元素都是1的矩阵:在MATLAB中元素都是1的矩阵由ones命令定义。
如
A=ones(m,n)定义一个m×n阶的元素都是1的矩阵。
矩阵的加法和减法:行列数相同的矩阵之间才可以进行加法和减法。
矩阵的乘法:若A的行数和B的列数相等,则它们可以相乘C=AB。
其中C的第i 行第j列元素等于A的第i行和B的第j列对应元素乘积之和。
逆矩阵:若两个方阵A和B满足:AB=I且BA=I,则称A和B互为逆矩阵。
在MATLAB 中M的逆矩阵由inv(M) 命令计算。
对于任一非奇异矩阵都可用inv命令计算其逆矩阵。
若MATLAB拒绝计算一个方阵的逆矩阵,则此矩阵一定是奇异的。
一个奇异矩阵的行列式是0(或者至少有一行(列)可以用其它行(列)通过多次加法和减法表示)。
行列式:方阵A的行列式是一个标量值,用det(A)或|A|表示。
在MATLAB中矩阵A
的行列式由det(A)命令计算。
如d =det (A )计算一个n 阶方阵A 的行列式。
特征值和特征向量:在MATLAB 中用eig 函数计算矩阵A 的特征值和特征向量。
如D =eig(A)计算方阵A 的特征值向量,而[V,D]=eig(A)则计算V 为特征向量,D 为特征值组成的对角矩阵。
矩阵运算的性质:设A ,B ,C ,…是矩阵,p,q,r,…是标量,而且对应矩阵运算有定义,则
A(BC)=(AB)C IA =AI =A
A(B +C)=AB +AC (A +B)C =AC +BC P(A +B)=pA+pB p(AB)=(pA)B=A(pB) A+B=B+A A+0=0+A=A A+(B+C)=(A+B)+C p(qA)=(pq)A (AB)T =B T A T (AB)-1=B -1A -1 det(AB)=det(A)det(B)
正定矩阵及其性质:……。
给定下列n 元线性方程组
⎪⎪
⎩⎪⎪⎨
⎧=+++=+++=+++n
n nn n n n n n n b x a x a x a b x a x a x a b x a x a x a ΛΛΛΛΛ22112
22221211
1212111 矩阵形式 b AX = 线性方程组解唯一的条件:0||≠A 。
定理5.1.1(线性方程组的初等运算) 下面三种变换可以将一个线性方程组变成另一个等价线性方程组:
交换变换:对调方程组的任意两个方程; 比例变换:用非零常数乘某一个方程;
替换变换:将某一个方程乘一个常数,再加到另一个方程上去。
定理5.1.2(矩阵的初等行变换)对一个线性方程组的增广矩阵进行如下的变换可得到一个等价的线性方程组的增广矩阵:
交换行变换:对调任意两行; 比例行变换:用非零常数乘某一行;
替换行变换:将某一行乘一个常数,再加到另一行上去。
初等矩阵及其性质:……。
5.2 消去法
直接法:假设计算过程中不产生舍入误差,经过有限次运算可求得方程组的精确解。
思路:将线性方程组变形成等价的三角方程组。
最基本的直接法就是消去法。
5.5.1.高斯消去法
思路:先逐次消去变量,将方程组化解成同解的上三角方程组,此过程称为消去过程;然后按方程组相反顺序求解上三角方程组,得到原方程组的解,此过程称为回代过程。
一般地,对于n阶线性方程组:
⎪⎪
⎩⎪⎪⎨
⎧=+++=+++=+++n
n nn n n n n n n b x a x a x a b x a x a x a b x a x a x a ΛΛΛΛΛ22112
222212********* 矩阵形式 b AX = (5.5.1) 消元过程:
第一步:用第1个方程乘上乘数-
1
,11,a a i 分别加到第i(i=2,3,…,n)个方程上去,从而消
去第i 个方程的首项,方程组(5.3.1)变为下列形式:
a 1,1x 1+ a 1,2x 2 + … + a 1,n x n =
b 1
)1(2,2a x 2 + … + )
1(,2n a x n =)1(2b (5.5.2)
……
)1(2,n a x 2 + … + )
1(,n
n a x n =)1(n b 其中 )
1(1,i a =
1
,11,a a i , j=2,3,…,n;
)1(,j i a = a i,j –)1(1,i a * a 1,j ;
)1(i b = b i –)1(1,i a * b 1 , i=2,3, … , n, j=2,.3,…,n.
若经过k-1 步后,从方程组(5.5.1)得到下列形式的方程组:
)1(11a x 1+)1(12a x 2 + … + )1(1n a x n = )
1(1b
)2(2,2a x 2 + … +)
2(,2n a x n = )2(2b
…… (5.5.3)
)1(1,1---k k k a x k-1+ … +)
1(,1--k n
k a x n = )1(-k k b ………
)1(1,--k k n a x k-1+ … +)
1(,-k n
n a x n =)1(-k n b 则对第k 列进行消元可将变为如下形式:
)1(11a x 1+)1(12a x 2 + … + )1(1n a x n = )
1(1b
)2(2,2a x 2 + … +)2(,2n a x n = )2(2b
………… (5.5.4)
)(,k k k a x k + … +)
(,k n k a x n = )(k k b
………
)(,k k n a x k + … +)
(,k n n a x n =)(k n b
其中: )
(,k k
i a
=
)1(1
,1)1(1
,-----k k k k k i a
a , i=k,…,n,
)(,k j i a = )1(,-k j i a –)(,k k i a *)
1(,-k j k a ;
)(k i b = )1(-k i b –)(,k k i a *)1(-k k b , i=k, … , n, j=k,…,n.
经过(n-1)步消元后,最后从(5.5.1)可得到一个上三角方程组:
)1(11a x 1+)1(12a x 2 + … + )1(1n a x n = )1(1b
)2(2,2a x 2 + … +)
2(,2n
a x n = )2(2
b …………………… (5.5.5)
)1(1,1---n n n a x n-1 +)
1(,1--n n
n a x n = )1(1--n n b )
(,n n n a x n = )(n n b
回代过程:
逐步回代得到解:
x n = )(,)(n n
n n n
a b
x i = ()(i i
b -
∑
+=n
i j 1
)(,i j i a x j ) / )(,i i i a , i=n-1,n-2, …, 1. (5.5.6)
上述方法称为高斯消去法。
高斯消去法进行的条件:
由以上过程有:0det )()2(22)1(11≠=n nn a a a A Λ,其中)1(-k kk a 称为约化主元素。
定理 5.5.1 设n n ij a A ⨯=}{非奇异,且各阶顺序主子式0≠∆k (n k ,,2,1Λ=),则
0)
(≠k kk a ,从而高斯顺序消元法可以进行,且得到方程组(5.5.1)的解为(5.5.6).
高斯消元法
约化主元素0)
(,≠k k k a (k=1,2,…,n )
⇔ A 的顺序主子式0≠k D (k=1,2,…,n)
定义5.3.1 称n 阶方阵A =(a i,j )n 是对角占优的,若其主对角元素的绝对值大于同行其它元素绝对值之和
∑
≠=n
k
j j ,1|a k,j | < |a k,k |, k=1,2,…,n
定理5.5.2 设方程组(5.5.1)是对角占优的,则在高斯消去法的消元过程中约
化主元素0)
1(≠-k kk
a (k=1,2,…,n-1). 练习1 在MATLAB 上用高斯消去法求解方程组:
-0.04 x 1 + 0.04 x 2 + 0.12 x 3 = 3, 0.56 x 1 - 1.56 x 2 + 0.32 x 3 = 1,
0.24 x 1 + 1.24 x 2 - 0.28 x 3 = 0,
gauss.m(自编程作业)
>>A=[-0.04 0.04 0.12;0.56 –1.56 0.32;0.24 1.24 –0.28]; >>B=[3;1;0]; >>X=gauss(A,B)
高斯消去法的计算时间复杂性:
高斯消去法需要次31(n 3+3n 2-n)乘除法和61
(2 n 3+3n 2-5n)次减法,即它的计算时间
复杂性为O(n 3)。
高斯消去法的缺点:
考察高斯消去法的消去过程,我们可以看出为使消元过程顺利完成,
)(,k k k a (k=1,2,…,n )必须全不为0。
但即使在消元过程中主元素)
(,k k k a (k=1,2,…,n )全不为
0,也不能保证能够避免在计算过程中产生巨大的误差。
因为在计算过程中,若主元素的绝对值过小,则舍入误差的影响会严重地损失精度。
提出主元法是为了控制舍入误差。
5.5.5.选主元素消元法(以避免a (k)k,k = 0 或过小的情形) 例2 考察方程组: 10-5 x +y =1 x + y =2 精确解为 x=100000/99999, y=99998/99999.
若使用四位浮点十进制运算进行消元 10-5 x +y =1;
-105 y =-105. % 105-1≈105 利用回代过程得
x=0, y=1 若用消去法求解方程组:
x + y =5. 10-5 x +y =1;
则可得 x=1 , y=1。
显然第一个结果严重失真,而第二个结果则比较精确。
究其原因,就是由于消元过程中主元素的绝对值过小。
选主元法分为列选主元法和全选主元法。
列选主元法:
思路:列选主元法就是在待消元方程的所在列中选取主元,经方程的行交换,置主元于主对角线位置后进行消元的方法。
主元:绝对值最大的元素。
列选主元法交换的原则:在第k 列中将主元所在的方程与第k 个方程进行交换,使主元位于第k 个对角线元素位置上。
具体做法:从第k 列位于第k 行以下的(n-k)个元素中选出绝对值最大的)
(,k k p a : | )(,k k p a | =n
j k ≤≤max {| )
(,k k j a |},若p>k ,则对调第k 和第p 个方程(增广矩阵的对应行)
,从而使)(,k k p a 成为新的主元素,然后进行第k+1步消元。
(以例1为例说明)
定理5.5.3设方程组(5.5.1)的系数矩阵对角占优,则在高斯消去法的消元过程
中主元素)(,k k k a (k=1,2,…,n )全都是部分选主元过程中的主元素。
练习2 在MATLAB 中用列选主元高斯消去法求解方程组: -0.04 x 1 + 0.04 x 2 + 0.12 x 3 = 3, 0.56 x 1 - 1.56 x 2 + 0.32 x 3 = 1,
0.24 x 1 + 1.24 x 2 - 0.28 x 3 = 0,
全选主元法:如果不按列选主元,而是在全体待选系数中选取,则得全选主元法。
(以上例说明)
主元消元法较一般消元法好处:
1.在消去过程中减少了舍入误差,在回代过程中采用数值较大的主元作分母,可以减少除法运算的误差; 2.具有很好的数值稳定性;
3.全选主元法的精度优于列选主元法,但排序时间较列选主元法长,既要进行行交换,又要进行列交换,还需要进行未知量序号的恢复,导致运算较大。
实际应用多的是列主元消元法. 5.3 直接三角分解法 5.3.1.Doolittle 分解法
一、LU 分解
由于三角方程组容易求解,因此我们现介绍将给定矩阵A 分解成单位下三角矩阵L 和上三角阵U 的乘积的方法。
⎥⎥⎥⎥⎦⎤⎢
⎢⎢⎢⎣⎡=
1112121ΛO M M n n l l l L , ⎥⎥⎥⎥⎦
⎤
⎢⎢⎢⎢⎣⎡=nn n n u u u u u u U M O ΛΛ22111211 在Gauss 顺序消元法的消去过程中有单位下三角初等矩阵s P P P ,,,21Λ,使
U s s =-A P P P 11Λ, 使方程组化为:
⎪⎪
⎩⎪
⎪⎨
⎧==++=+++n n nn n n n n y x u y x u x u y x u x u x u
22122211212111ΛΛΛΛ 即: UX=Y 其中U 为上三角矩阵。
于是: U A s 1
1
21
1---=P P P Λ,令1
1
21
1---=s L P P P Λ,则L 为单位下三角矩阵,且有A=LU 。
(称将A 进行LU 分解)
由此,方程组AX=B 化为:LUX=B,于是方程组可以化为以下两个三角方程组求解: LY=B; UX=Y.
定理5.3.1 设A 为一个矩阵,且A 的前n-1个顺序主子式都不等于零,则存在一个单位下三角阵L 和上三角阵U 使得A =LU 。
由
⎪⎪⎪⎪
⎪⎭
⎫
⎝⎛⋅⎪⎪⎪⎪⎪⎭⎫
⎝⎛=⎪⎪⎪⎪⎪⎭⎫
⎝⎛nn n n n n nn n n n n u u u u u u l l l a a a
a a a a a a ΛM M M ΛΛ
ΛM M M Λ
ΛΛ
M M M ΛΛ000101001222112112121212222111211 可计算出ij ij u l 和.
注:1. 矩阵A 的L 与U 的求法: (1)可用选定系数法;
(2) MATLAB不选主元LU分解:[L1,U1]==lu(A);
MATLAB列选主元LU分解:[L1,U1,P]==lu(A)
5.若上述LU分解中, 取L是单位下三角阵,而U是上三角阵,称此分解为杜利特尔(Doolittle)分解;取U是单位上三角阵,而L是下三角阵,那么称此分解为克劳特(Courant)分解.
二、用直接三角分解法解线性方程组
若A=LU,AX=b ⇔LUX=b 则先由 LY=b求解Y;再由UX=Y求解X。
(不选主元)
若PA=LU,AX=b ⇔ LUX=Pb 则先由 LY=Pb求解Y;再由UX=Y求解X。
(列选主元)
例3用LU分解法求解
x
1 +
2 x
3
= 5,
x
2 +x
4
= 13,
x 1 + 2x
2
+4x
3
+3x
4
=17,
x
2
+3x
4
=7,
练习3在MATLAB中用三角分解法求解方程组:
-0.04 x
1 + 0.04 x
2
+ 0.12 x
3
= 3,
0.56 x
1 - 1.56 x
2
+ 0.32 x
3
= 1,
0.24 x
1 + 1.24 x
2
- 0.28 x
3
= 0,
tridiv.m(自编程作业)
>>A=[-0.04 0.04 0.12;0.56 –1.56 0.32;0.24 1.24 –0.28];
>>B=[3;1;0];
>>[L1,U1]==lu(A) % MATLAB不选主元LU分解
>>[L1,U1,P]==lu(A) % MATLAB列选主元LU分解
>>[L,U]=tridiv(A)
>>X=forsub(L,B)
>>X=backsub(U,X) (与练习1、2的结果比较)
3.LU分解法和Gauss消元法解线性方程组是等价的,它们具有基本相同的计算量和精度,但LU分解法具有更好的设计灵活性。
4.LU 分解法和Gauss 消元法具有相同的计算复杂性,即它们的运算量相当.
5.3.2 追赶法(求解三对角方程组)
对一些具有特殊形状或性质的系数矩阵的线性方程组,我们可用特殊的方法求解。
定义3.3.1下列形状的方程组称为三对角方程组:
d Ax =
其中:⎥⎥⎥
⎥⎥
⎥⎦
⎤
⎢⎢⎢⎢
⎢⎢⎣⎡=-n n
n b a c b a c b a c b A 133
222
1
1
000O O
O ΛΛ,⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=n x x x x M 2
1,⎥
⎥⎥⎥⎦⎤
⎢⎢⎢⎢⎣⎡=n d d d d M 21
即:
b 1 x 1+
c 1 x 2 =
d 1 a 2 x 1+ b 2 x 2 + c 2 x 3 = d 2
…………………………………… (5.3.1) a n-1 x n-2 + b n-1 x n-1 + c n-1 x n = d n-1 a n x n-1 + b n x n = d n
其系数矩阵A 称为三对角阵,这是一种特殊的稀疏矩阵(Sparse matrix )。
定理5.3.1若(3.3.1)的系数矩阵对角占优,即 | b 1 |>| c 1| ,
| b i | > | a i | + | c i | , i=2,3, …, n-1 | b n | > | a n | 则它非奇异。
追赶法:
(一) 消元过程:
用高斯消去法将(3.3.1)变为下列形式的单位上三角方程组:
x 1+ u 1 x 2 = y 1 + x 2 + u 2 x 3 = y 2
…………………………………… (5.3.2) x n-1 + u n-1 x n = y n-1 x n = y n
其中 u 1 = c 1/ b 1, y 1 = d 1 / b 1;
u i = c i / (b i – u i-1 a i ), i=2,3,…,n-1; y i = (d i – y i-1 a i ) / (b i - u i-1 a i ) , i=2,3,…,n (二) 回代过程
x n = y n ,
x i = y i – u i x i+1, i=n-1,n-2, (1)
练习4 在MATLAB 上用追赶法求解三对角方程组: -4 x 1 + x 2 = -27, x 1 -4 x 2 + x 3 = -15,
x 2 -4 x 3 +x 4 = -15,
x 3 -4x 4 = -15。
(自编程作业)
5.3.3 平方根法(求解对称正定方程组)
系数矩阵必须为对称矩阵才可用此法。
定理3.3.1(对称阵的三角分解定理) 设A 为n 阶,对称阵,且A 的所有顺序主子式均不为零,则A 可惟一分解为A =LDL T 。
其中L 为单位下三角阵,D 为对角阵.
若A 为对称正定矩阵,那么D 的对角元素全为正.
定理5.3.2(对称正定矩阵的三角分解或Cholesky 分解) 若A 为n 阶对称正定矩阵,则存在一个实的非奇异下三角L 使A =LL T 。
若限定L 的对角元素为正时,这种分解是唯一的。
(1) 平方根法(A=LL T )
AX =B 即LL T X=B LY =B , L T X =Y 。
设A =(a i,j )n , L=(l i,j )n ,则平方根法的计算公式为: l i,j = (a i,j - ∑-=1
1j k l i,k l j,k )/ l j,j ;j=1 ,2,…,i-1
l i,i =
∑-=-1
1
2,,i k k i j i l a ; i=1 ,2,…,n
求解下三角方程组LY =b 的公式为
y i =(b
i
- ∑-
=
1
1
i
k
l
i,k
y
k
) / l
i,i
, i=1,2,…,n;
求解上三角方程组L T X=Y的公式为
x i =(y
i
- ∑
+
=
n
i
k1
l
k,i
x
k
) / l
i,i
, i=n,…,2,1;
练习5 在MATLAB中用平方根法求解方程组:
4 x
1- 2 x
2
+ 4 x
3
= 10,
-2x
1+17x
2
+ 10x
3
= 3,
-4 x
1 + 10x
2
+9x
3
=-7,
chle.doc(自编程作业)
>>A=[4 –2 4;-217 10;-410 9];
>>L1=chol(A) %MATLAB平方根分解
>>L=chle(A)
(2) 改进平方根法(A=LDL T) (此时A为对称阵但不一定正定)
AX=B LY=B, L T X=D-1Y
这种分解的公式为
l i,j =(a
i,j
- ∑-
=
1
1
j
k
d
k
l
i,k
l
j,k
)/ d
j
; j=1 ,2,…,i-1
d
i = a
i,i
–∑-
=
1
1
i
k
d
k
(l2
i,k
);i=1 ,2,…,n
求解下三角方程组LZ=B的公式为
z i =b
i
- ∑-
=
1
1
i
k
l
i,k
z
k
; i=1,2,…,n;
求解对角方程组DY=Z的公式为
y i = z
i
/ d
i
;i=1,2,…,n;
求解上三角方程组L T X=Z的公式为
x i = y
i
- ∑
+
=
n
i
k1
∑n
k=i+1
l
k,i
y
k
; i=n,…,2,1;
例4 分别用平方根法和改进平方根法求解
4 x 1 - x 2 + x 3 = 6, -x 1 +4.25x 2 + 5.75x 3 =-0.5,
x 1 + 5.75x 2 +3.5x 3 =1.25,
练习6 在MATLAB 中用改进的平方根法求解方程组: 4 x 1 - 2 x 2 + 4 x 3 = 10, -2x 1 + 17x 2 + 10x 3 = 3,
-4 x 1 + 10x 2 +9x 3 =-7,
advchle.doc (自编程作业) >>A=[4 –2 4;-217 10;-410 9];
>>[L,D]=advchle(A) %改进平方根法解方程组 5.4.误差分析
求解AX=b 时A 有误差A δ或b 有误差b δ时都会引起解的误差。
解方程组AX=b 时,总假定系数矩阵A 和常数项b 是准确的,但实际是有误差的,这些误差对方程组Ax=b 的解x 产生多大的影响?矩阵的条件数将给出一种粗略的衡量尺度。
一、条件数
定义 若矩阵A 非奇异,称p
p
p A A A Cond 1)(-=为A 的条件数。
其中p ⋅代表矩阵的
某种范数。
条件数的性质:
(1)1)(≥A cond ;)()(1-=A cond A cond ; (2)0, ),()(≠∈=αααR A cond A cond ;
(3)设n λλ,1为A 的按模最大和最小的特征值,则:
n
A cond λλ1
)(≥
若A 是对称阵,则有:n
A cond λλ1
2)(=。
MATLAB :cond(A) 二、数据扰动的误差分析
对于方程组Ax=b (1)方程右端数据扰动
若常数b 有小扰动b δ,解x 的扰动x δ可以表示为:
b
b
A cond x
x
δδ)
(≤
因为b b x x A δδ+=+)(,b Ax =,于是b x A δδ=,b A x δδ1-=,有:
b A
x δδ1
-≤,x A b ≤,从而有
b
b
A cond b
b
A
A x
x
δδδ)
(1
=≤-,可见条件数
)(A cond 表示常数项b 的相对误差在解中的放大率。
(2)系数矩阵扰动
若系数矩阵A 有小的扰动A δ,这时方程组的解也有扰动x δ,于是:
b x x A A =++))((δδ
扰动A δ对x δ的影响可以表示为:
A
A
A cond A
A
A cond x
x
δδδ)
(1)(-≤
当A
A
A cond δ)
(很小时,1)
(1≈-A
A
A cond δ,条件数)(A cond 表示矩阵A 的相对误
差在解中的放大率。
因此,Cond(A)较大(Cond(A)>>1)的矩阵称为“坏矩阵”或“病态矩阵”,对条件数大的矩阵,小的误差可能会引起解的失真。
一般地,若A 的按模最大特征值与按模最小特征值之比值较大时,矩阵就会呈病态。
特别地,当detA 很小时,A 总是病态的。
(3)方程右端和系数矩阵都有扰动时,有:
b b x x A A δδδ+=++))((
若11<-A A δ,则有误差估计:
⎪
⎪⎭⎫ ⎝⎛+-≤b b A A A
A A cond A cond x x
δδδδ)(1)
( 说明:
⏹ 为了计算一个病态问题的解,计算的精度必须非常高. 对于轻微的病态问题,用双精度既可. 对于较严重的病态问题则需要更高的精度,如四倍精度. ⏹ 对于病态矩阵,逆矩阵和行列式的计算都会变得不精确. ⏹ 具备下列特征的问题可认为是病态的: (1) )det()det(1-A A 的计算值与1有较大偏差 (2) inv(inv(A))的计算值明显不等于A (3) A*inv(A)的计算与单位矩阵有较大偏差
(4) cond(A*inv(A))与1的偏差可作为A*inv(A)与单位矩阵的偏差的精确度量 ⏹ 线性方程组解的误差也取决于计算环境的精度
内容总结
线性方程组的直接解法(AX=b,A非奇异)
1、消元法:
顺序高斯消元法(要求A的前n-1个顺序主子式非零)
全选主元消元法、
列选主元消元法、
2、矩阵分解法:LU分解(不选主元时A=LU;选主元时PA=LU)
3、追赶法:针对三对角方程组,要求A对角占优
4、平方根法与改进平方根法:A=L
1L
1
T或A=LDL T,要求A对称正定
注:
1.在MATLAB中,可用下列方式求解线性方程组:
X=A\B (自动根据系数矩阵的特点选择算法)
X=inv(A)*B
X=(A)∧(-1)*B
5.直接法是一种计算量小而精度高的方法,选主元的算法有很好的数值稳定性。
从计算简单出发多选用列主元法。
线性方程组的病态程度是其本身的固有特性,因此即使用数值稳定的方法求解,也难以克服严重病态导致的解的失真。
在病态不十分严重时,用双精度求解可减轻病态的影响。
3.在实际应用中如何选择算法是一个重要问题,往往从三个方面考虑:①解的精度高;
②计算量小;③所需计算机内存小。
但这些条件相互间是矛盾而不能兼顾的,因此实际计算时应根据问题的特点和要求及所用计算机的性能来选择算法。
一般说,系数矩阵为中、小型满矩阵,用直接法较好;当系数矩阵为大型、稀疏矩阵时,有效的解法是下章将要讨论的迭代法。
4.Gauss顺序消元法直接三角分解算法数值不稳定,选主元消元法和列选三角分解算法是数值稳定的,平方根法和追赶法都是数值稳定的算法.
本章作业题
1、编写一个函数M 文件用高斯消去法的消去过程求矩阵A 的行列式。
⎥⎥
⎥⎥⎦
⎤
⎢
⎢⎢
⎢⎣⎡----=1115
142023123211A 2、编写一个函数M 文件用高斯消去法求矩阵B 的逆矩阵。
⎥⎥
⎥⎦
⎤
⎢⎢⎢⎣⎡=110121013B 3、求三阶Hilbert 矩阵3H 的条件数)(3H Cond ,并说明方程b X H =3是否“病态
⎥⎥⎥⎥⎥⎥⎦
⎤
⎢⎢⎢⎢
⎢⎢⎣⎡=514
13
1413121
3121
13H 4、编程作业(教程中)。