三角分解解线性方程组的公式47页
合集下载
研究生数值分析(8)直接三角分解法
§2 直接三角分解法
三角分解法也是直接法，基本思想是： 将系数矩阵A分解为两个三角形矩阵L和U的 乘积A=LU ，将方程组AX=b的求解问题归结为 两个三角形方程组 LY=b与UX=Y的求解问题。
即：先由LY=b求出Y ，然后由UX=Y求出X ， 从而获得AX=b的解。
(1) A为一般稠密（零元素占很小比例）矩阵 的杜利特尔(Doolittlr)和克劳特(Crout)分解法；
1 0 0
3
7
2
2 2
1 12
7 6
1
0
0
16
3
等价的三角方程组为
12
3 3 2
3 7
2
x1
x2
15
15
2
16
x3
16
an1
bn1
cn1
xn1
fn1
an bn xn fn
其中方程组AX=f 的系数矩阵A的元素满足条件：
b1 c1 0
bi ai ci
bn
an
0
（1） 且 aici 0
(i 2,3, , n 1)
由以上推导过程知，方程组AX=f 有唯一解
由（3）式可得计算 i 的递推公式
1i
c1 ci
/ b1 /(bi
ai i1)
(i 2,3,
（4）
三角分解法也是直接法，基本思想是： 将系数矩阵A分解为两个三角形矩阵L和U的 乘积A=LU ，将方程组AX=b的求解问题归结为 两个三角形方程组 LY=b与UX=Y的求解问题。
即：先由LY=b求出Y ，然后由UX=Y求出X ， 从而获得AX=b的解。
(1) A为一般稠密（零元素占很小比例）矩阵 的杜利特尔(Doolittlr)和克劳特(Crout)分解法；
1 0 0
3
7
2
2 2
1 12
7 6
1
0
0
16
3
等价的三角方程组为
12
3 3 2
3 7
2
x1
x2
15
15
2
16
x3
16
an1
bn1
cn1
xn1
fn1
an bn xn fn
其中方程组AX=f 的系数矩阵A的元素满足条件：
b1 c1 0
bi ai ci
bn
an
0
（1） 且 aici 0
(i 2,3, , n 1)
由以上推导过程知，方程组AX=f 有唯一解
由（3）式可得计算 i 的递推公式
1i
c1 ci
/ b1 /(bi
ai i1)
(i 2,3,
（4）
三角分解
l
ik
uk j
ai j
min( i , j ) k 1
l
ik
uk j
固定 i : i 1 对 j = i, i+1, …, n 有 aij l ik ukj uij
k 1
lii = 1
uij aij l ik ukj
k 1
i 1
a
i 1
将 i ，j 对换，对 j = i, i+1, …, n 有 一般采用列主元 ii a ji l jk uki l ji u k 1 法增强稳定性。但注意 i 1 b 也必须做相应的 l ji (a ji l jk uki ) / uii b 行交换。 k 1 算法:道立特分解法
得到原方程组的解, 求解方程组计算公式:
yi bi lik yk , i 1,2,, n.
k 1 i 1
xi ( yi uik xk ) uii , i n,,1.
k i 1
n
说明: L和U的存放;计算∑aibi; 运算量n3/3; Doolittle分解法;求解方程组优点. 注： L 为一般下三角阵而 U 为单位上三角阵的分解称为 Crout 分解。 ~~ 实际上只要考虑 A* 的 LU 分解，即A* LU ，则 保存L,U;对解A不变而仅b变 ~ ~ A U * L * 即是 A 的 Crout 分解。 化的方程组很方便
列主元素高斯消去法相当于先进行一系列行交换后再对 PAX Pb 应用顺序高斯消去法.
定理8(列主元三角分解) 若A为非奇异矩阵, 则存在排列 矩阵P使得 PA LU 其中L为单位下三角矩阵,U为上三角矩阵. 说明: L, U, Ip的存贮.
用直接三角分解法解线性方程组
三角阵。等式左边是单位下三角阵，右边是上三角阵，要使等式
成则立L ， L只1 ,能U等于U单1，位即矩此阵三I。角于分是解L唯1一L1。
UU
1 1
I,
1 2 1
1 2 1
例7 解：
设 A 3 7
1
，试将A进行三角分解。
1 1 3
由高斯消去法得到
m21
3 1
3，m31
1 1
1
m 32
L1
1 1
0 0 1 2
例：
求
0 1 2
0 1 0
3 0 1
103的PLU分 解 。
解：用1，2， ，n的排列表示n阶置换阵P，其中排列的第i个元素
j，表示P的i行非零元素位于j列。则分解过程如下：
1 0 0 1 2
3 1 1 0 1
3 1 1 0
2 0
43
1 2
0 1 0
3 0 1
0 1 3
Ux j y j
Ly
j
bj
n1
k
n(n 1)
n2
n 次乘法
k 1
2
22
Ux j y j n k n(n 1) n2 n 次乘除法
k 1
2
22
即共需n 2 次乘除法运算。
n 2 次 乘 除 法
三角分解法的存放元素的方法：
以A (a ij )33 为例，
a11 A a21
1 mk1,k 1
k，
Lk1
1 mk1,k 1
k
mnk
1
mn,k
1
A ( L11 L21
L1 n1
)U
LU，1
a (1) 11
解线性方程组的矩阵三角分解法
为了节省存储空间，通常用 A 的绝对下三角部分来存放 L (对角线元素无需存储)，用 A 的上三角部分来存放 U
6
PLU 分解
矩阵的 PLU 分解
PA LU
k -1
ukj akj lki uij i 1
k 1
lik aik lijujk ukk
j 1
l21 M
l22
O
l22 O
lnM2
a21 M
a22
O
aM2n
ln1 L ln,n1 lnn
lnn an1 an2 L ann
计算公式
n
j 1
aij likl jk likl jk l jjlij
n
j 1
aij likl jk likl jk l jjlij
k 1
k 1
10
平方根法
Ax b A 对称正定
算法 ：(解对称正定线性方程组的平方根法 )
计算 A 的 Cholesky 分解 解方程：Ly = b 和 LTx = y
y1 b1 l11 ,
yi
bi
i1
lik yk
lii
,
k 1
xn yn lnn ,
xi
yi
n
lki
xk
lii
ki1
i = 2, 3, …, n i = n-1, …, 2, 1
11
改进的 Cholesky 分解
改进的 Cholesky 分解
6
PLU 分解
矩阵的 PLU 分解
PA LU
k -1
ukj akj lki uij i 1
k 1
lik aik lijujk ukk
j 1
l21 M
l22
O
l22 O
lnM2
a21 M
a22
O
aM2n
ln1 L ln,n1 lnn
lnn an1 an2 L ann
计算公式
n
j 1
aij likl jk likl jk l jjlij
n
j 1
aij likl jk likl jk l jjlij
k 1
k 1
10
平方根法
Ax b A 对称正定
算法 ：(解对称正定线性方程组的平方根法 )
计算 A 的 Cholesky 分解 解方程：Ly = b 和 LTx = y
y1 b1 l11 ,
yi
bi
i1
lik yk
lii
,
k 1
xn yn lnn ,
xi
yi
n
lki
xk
lii
ki1
i = 2, 3, …, n i = n-1, …, 2, 1
11
改进的 Cholesky 分解
改进的 Cholesky 分解
线性代数方程组的解法
2 3 2 n O( n ) 3
mult a(i , j ) a( j, j ); for k j 1 : n a(i , k ) a(i , k ) mult * a( j , k ); end b(i ) b(i ) mult * b( j ); end
end
LU分解
求A的LU分解（L是下三角矩阵，U是上三角矩阵）
1 1 1 1 3 4 3 4
LU分解
性质1 设向量
, xn ) 且 xk 0 T 则存在唯一的下三角阵 Lk I lk ek ,满足 x ( x1 , x2 ,
T
Lk x ( x1 ,
第三章 线性方程组的直接解法
/*Direct Method for Solving Linear Systems*/
求解 A x b, A R
Cramer法则:
nn
det( A) 0
Di xi D
i 1, 2,
,n
所需乘除法的运算量大约为(n+1)!+n
n=20时，每秒1亿次运算速度的计算机要算30多万年！
Gauss消去法的消元过程算法
for for
j 1: n 1
i j 1: n
2 3 2 n O( n ) 3
mult a(i , j ) a( j, j ); for k j 1 : n a(i , k ) a(i , k ) mult * a( j , k ); end b(i ) b(i ) mult * b( j ); end
方程组可化为下面两个易求解的三角方程组
Ly b Ux y
二、 高斯消去法
三角分解法
. . . . l21 1 . . = . . . . ... . . . . . an1 ... ann ln1 ... 1 ... . . . . . . unn
ai j =
min( i , j ) k =1
∑l
ik
uk j
ai j =
min( i , j )
一般采用列主元 对换, 将 i ,j 对换,对 j = i, i+1, …, n 有 一般采用列主元 ii a ji = ∑ l jk uki + l ji u k= k =1 法增强稳定性. 法增强稳定性.但注意 v i 1 b 也必须做相应的 l ji = ( a ji ∑ l jk uki ) / uii b 行交换. 行交换. k =1
=I
Upper-triangular
Lower-triangular With diagonal entries 1
注: L 为一般下三角阵而 U 为单位上三角阵的分解称为 单位上三角阵的分解称为 Crout 分解. 分解. ~~ 分解, 实际上只要考虑 A* 的 LU 分解,即A* = L U ,则 ~ ~ A= U * L* 即是 A 的 Crout 分解. 分解. =
(
)
n1 Step 6 Set l = 运算量为2 O(n3/6), 比普通 ann ∑ k =1 lnk ; , 比普通LU nn Step 7 Output ( lij for j = 1, …, i and i = 1, …, n );A = LDLT 分解少一半, 次开方. 分解少一半,但有 n 次开方.用
mn1
v A b
( 2) (2)
(1 (1 ( a11) a12) ... a11) n Step n 1: (2 ( v a22) ... a22) n Ln1Ln2 ... L1 A b = ... . . . (n ann)
ai j =
min( i , j ) k =1
∑l
ik
uk j
ai j =
min( i , j )
一般采用列主元 对换, 将 i ,j 对换,对 j = i, i+1, …, n 有 一般采用列主元 ii a ji = ∑ l jk uki + l ji u k= k =1 法增强稳定性. 法增强稳定性.但注意 v i 1 b 也必须做相应的 l ji = ( a ji ∑ l jk uki ) / uii b 行交换. 行交换. k =1
=I
Upper-triangular
Lower-triangular With diagonal entries 1
注: L 为一般下三角阵而 U 为单位上三角阵的分解称为 单位上三角阵的分解称为 Crout 分解. 分解. ~~ 分解, 实际上只要考虑 A* 的 LU 分解,即A* = L U ,则 ~ ~ A= U * L* 即是 A 的 Crout 分解. 分解. =
(
)
n1 Step 6 Set l = 运算量为2 O(n3/6), 比普通 ann ∑ k =1 lnk ; , 比普通LU nn Step 7 Output ( lij for j = 1, …, i and i = 1, …, n );A = LDLT 分解少一半, 次开方. 分解少一半,但有 n 次开方.用
mn1
v A b
( 2) (2)
(1 (1 ( a11) a12) ... a11) n Step n 1: (2 ( v a22) ... a22) n Ln1Ln2 ... L1 A b = ... . . . (n ann)
三角分解解线性方程组报告
三角分解解线性方程组报告三角分解是一种常用的解决线性方程组的方法，其基本思想是将一个复杂的线性方程组化简为几个简单的三角形式方程，从而简化计算过程，提高求解效率。
在实际应用中，三角分解可以用于求解线性方程组、求解矩阵的逆等问题。
本报告将介绍三角分解的原理、步骤和应用，并通过一个具体的实例来说明三角分解方法的有效性。
一、三角分解的原理三角分解的基本原理是将线性方程组的系数矩阵通过一系列的初等行变换，转化为一个上三角形矩阵和一个下三角形矩阵的乘积形式，即LU 分解。
其中，L为下三角形矩阵，U为上三角形矩阵。
通过LU分解，可以将原始的线性方程组转化为两个三角形矩阵的乘积形式，从而简化求解步骤。
二、三角分解的步骤1.对系数矩阵进行初等行变换，将其转化为一个上三角形矩阵。
2.对每个行变换操作，记录一个对应的下三角形矩阵的乘积，以便最后求解时使用。
3.将转化后的上三角形矩阵和对应的下三角形矩阵相乘，得到原始系数矩阵。
三、三角分解的应用1.求解线性方程组：通过三角分解可以将复杂的线性方程组转化为简单的三角形式方程组，从而提高求解效率。
通过求解上三角形矩阵和下三角形矩阵，可以得到原始线性方程组的解。
2.求解矩阵的逆：通过三角分解可以将原始的矩阵转化为一个上三角形矩阵和一个下三角形矩阵的乘积形式。
通过求解上三角形矩阵和下三角形矩阵，可以得到原始矩阵的逆矩阵。
这在计算机图形学、信号处理等领域中具有广泛的应用。
四、三角分解实例假设有一个线性方程组：$$\begin{aligned}2x+3y+4z&=8\\3x+4y+5z&=9\\4x+5y+6z&=10\\\end{aligned}$$首先，将系数矩阵进行初等行变换，将其转化为上三角形矩阵：$$\begin{bmatrix}2&3&4\\3&4&5\\4&5&6\\\end{bmatrix}\xrightarrow{R_2 - \dfrac{3}{2}R_1}\begin{bmatrix}2&3&4\\-\dfrac{1}{2} & -\dfrac{1}{2} & -\dfrac{1}{2} \\ 4&5&6\\\end{bmatrix}\xrightarrow{R_3 - 2R_1}\begin{bmatrix}2&3&4\\-\dfrac{1}{2} & -\dfrac{1}{2} & -\dfrac{1}{2} \\ 0&-1&-2\\\end{bmatrix}$$得到上三角形矩阵U为：$$\begin{bmatrix}2&3&4\\0 & -\dfrac{1}{2} & -\dfrac{1}{2} \\0 & 0 & -\dfrac{8}{7} \\$$同时，记录每个行变换对应的下三角形矩阵，即：$$L_1=\begin{bmatrix}1&0&0\\\dfrac{3}{2} & 1 & 0 \\2&0&0\\\end{bmatrix},L_2=\begin{bmatrix}1&0&0\\0&1&0\\-2&3&1\\\end{bmatrix}$$将上三角形矩阵U和下三角形矩阵L相乘，得到原始的系数矩阵：$$2&3&4\\3&4&5\\4&5&6\\\end{bmatrix}= L_2 \cdot L_1 \cdot U$$通过求解上三角形矩阵U和下三角形矩阵L，可以得到原始线性方程组的解。
解线性方程组的三角分解法
⎛1 0 ⎛1 3 2 1 ⎞ ⎜3 1 ⎜ ⎟ ⎜ 3 13 12 9 ⎜ ⎟ = ⎜2 3 ⎜ 2 12 29 15 ⎟ ⎜ 2 ⎜ ⎟ ⎜ 1 9 15 34 3 ⎝ ⎠ ⎜1 ⎜ 2 ⎝ ⎛1 ⎜ ⎜3 =⎜ ⎜2 ⎜ ⎜ ⎜1 ⎝ 0 1 3 2 3 2 0⎞ ⎟ 0⎟⎛1 ⎜ ⎟⎜0 1 0⎟⎜ 0 ⎟⎜ 1 ⎟ 0 1⎟⎝ 4 ⎠ 0 0
⎛1⎞ 0 0 0 ⎞ ⎛ y1 ⎞ ⎛ 1 ⎞ ⎛ y1 ⎞ ⎜ 5 ⎟ ⎜ ⎟ ⎜ ⎟ ⎟⎜ ⎟ ⎜ ⎟ y 4 0 0 ⎟ ⎜ y2 ⎟ ⎜ 10 ⎟ 2 ⎟。 = ，得到 ⎜ 2 ⎟ = ⎜ ⎜ ⎜ ⎟ ⎜ ⎟ ⎟ ⎜ ⎟ y3 0 16 0 y3 24 3⎟ ⎜ ⎟ ⎜ ⎟ ⎟⎜ ⎟ ⎜ ⎟ 0 0 25 ⎠ ⎝ y4 ⎠ ⎝ 50 ⎠ ⎝ y4 ⎠ ⎜ 2 ⎟ ⎜2⎟ ⎝ ⎠
⎛ 1 0 0 0⎞ ⎜ ⎟ 6 ⎞ ⎜ 2 1 0 0⎟⎛3 ⎟ห้องสมุดไป่ตู้⎟⎜0 6 ⎟ ⎜3 ⎜ ⎟⎜ = 1 8⎟ ⎜ 0 1 0⎟⎜0 ⎟ ⎜3 ⎟⎜ 12 ⎠ ⎜ ⎟⎝0 4 1 0 1 ⎜ ⎟ ⎝3 ⎠
6 3 6⎞ ⎟ 1 0 2⎟ 0 2 6⎟ ⎟ 0 0 2⎠
解方程组
⎛ 1 0 0 0⎞ ⎜ ⎟ ⎜ 2 1 0 0 ⎟ ⎛ y1 ⎞ ⎛ 9 ⎞ ⎜3 ⎟⎜ y ⎟ ⎜ 7 ⎟ ⎜1 ⎟⎜ 2 ⎟ = ⎜ ⎟ 0 1 0 ⎟ ⎜ y3 ⎟ ⎜ 1 ⎟ ⎜ ⎜3 ⎟⎜ ⎟ ⎜ ⎟ ⎜4 ⎟ ⎝ y4 ⎠ ⎝13 ⎠ 1 0 1⎟ ⎜ ⎝3 ⎠
⎛3 ⎜ ⎜2 ⎜1 ⎜ ⎝4
6 3 5 2 2 3
6 6 8
1 0 0 0 1 0 0 0 1
9 4 12 0 0 0
2 3 1 0 2 6 − 3 4 1 0 4 − 3 1 0 2 −
2019一基本的三角分解法LU分解.ppt
A的第r行元素主对角线以右元 素arj ( j r,, n)为
r
arj lrkukj k 1
j r,,n r 1,2,, n
同样
可知A的第r列元素主对角线以下元 素 air (i r 1,, n)为
r
air lik ukr k 1
i r 1,, n r 1,2,, n 1
设 A LU
即 a1 c1
b2 a2 c2
p1 b2 p2
1 q1
1 q2
b3
bn1 an1 cn1
pn1
1 qn1
bn
an
bn
pn
解Ly b,得
y1 b1
j 1
yr br lrj y j
r1
y1 y2 y3 y4 T 10 20 17 /11 16T
解Ux y,得
x1 x2 x3 x4 T 1 2 3 4T
xn
yn unn
n
yr urj x j
xr
jr1
L的第r列 ------(4)
称上述(1) ~ (4)式所表示的分解过程为LU分解
对于线性方程组
Ax b
系数矩阵非奇异,经过LU分解后
A LU
线性方程组可化为下面两个三角形方程组
Ly b
Ux y
y为中间未知量向量
1
l21
1
L
l31
l32
2_2 三角分解法
( 1) ( 1) ( 1) a11 ... a1 b n 1
(2) A b
( 2)
(1) (1) ) a11 a12 ... a1(1 n Step n 1: ( 2) ( 2) a ... a 22 2n Ln1 Ln 2 ... L1 A b ... . . . ( n) a nn
T
3 4 LU . 24 x (1, 2, 3) .
T
y (14, 10, 72)T ,
Ux (14, 10, 72)
§2 Matrix Factorization – Tridiagonal System
追赶法解三对角方程组
/* Crout Reduction for Tridiagonal Linear System */
思 通过比较法直接导出L 和 U 的计算公式。 路 反复计算, a1n 1 a11 ... …… u11 ... u1n 很浪费哦
l . . . . . 21 1 . . . . . ... . . . . . an1 ... ann l n1 ... 1 ... . . . . . . unn
L L ... L
1
1 1
1 2
1 n1
1
记为
mi , j
L
1
1) 12 ( 2) 22
Hey hasn’t GE given me 单位下三角阵 (1) Factorize A first, then When you have to ... /* a1unitary lower-triangular n headache? enough Why matrix */ Could you be more for (every you only have to solve the system for ) b know its aI 2 do have to ... 2n solve two simple triangular specific, please? different with a fixed b matrix form??! . systems and L ... . A .y b . Ux y .
三角分解解线性方程组的公式
1 u11 u1n A l21 ln1 1 unn
因为 k u11u22 ukk 0 k 1,2,n ukk 0
1 u11 l u22 21 A ln1 ln 2 1
• 追赶法解三对角方程组
b1 a 2 c1 b2 a3 c2 b3 c3 ai bi ci an an 1 bn 1 x1 d1 x d 2 2 x3 d 3 xi d i cn 1 xn 1 d n 1 bn xn d n
分解公式
k 1 m 1
2 2 akk lkm lkm lkm lkk
k 1
2 lkk akk lkm m 1
第一步 ： a l 2 l a 11 11 11 11
ai1 l11 li1 li1 ai1 / l11
ik
aik lim lkm
其中L为单位下三角阵，D是对角元全为正的对角 阵且这种分解式唯一;
2A LLT
其中L为下三角阵，当限定L的对角元为正时的分 解式唯一（Cholesky分解）.
August 5, 2013 yfnie@ 5
平方根法（Cholesky分解）
定理证明
证明： 因为 A对称正定，故其顺序主式 k 0k 1,2,n ， 从而A 有唯一的Doolittle分解
August 5, 2013
yfnie@
8
平方根法（Cholesky分解）：
l11 l11 l21 ln1 l l22 l22 A 21 ln1 lnn lnn
因为 k u11u22 ukk 0 k 1,2,n ukk 0
1 u11 l u22 21 A ln1 ln 2 1
• 追赶法解三对角方程组
b1 a 2 c1 b2 a3 c2 b3 c3 ai bi ci an an 1 bn 1 x1 d1 x d 2 2 x3 d 3 xi d i cn 1 xn 1 d n 1 bn xn d n
分解公式
k 1 m 1
2 2 akk lkm lkm lkm lkk
k 1
2 lkk akk lkm m 1
第一步 ： a l 2 l a 11 11 11 11
ai1 l11 li1 li1 ai1 / l11
ik
aik lim lkm
其中L为单位下三角阵，D是对角元全为正的对角 阵且这种分解式唯一;
2A LLT
其中L为下三角阵，当限定L的对角元为正时的分 解式唯一（Cholesky分解）.
August 5, 2013 yfnie@ 5
平方根法（Cholesky分解）
定理证明
证明： 因为 A对称正定，故其顺序主式 k 0k 1,2,n ， 从而A 有唯一的Doolittle分解
August 5, 2013
yfnie@
8
平方根法（Cholesky分解）：
l11 l11 l21 ln1 l l22 l22 A 21 ln1 lnn lnn
方程组直接三角分解法
如果U的第1至k-1列和L的第1至k-1列已经算出，则由
ak1 lk1 , k 2,3, u11
k r 1
, n.
（4.2.2）
a l u ,j k , k 1 , , n , kj kr rj
可得U的第k行元素
ukj =akj 同理，由
k
k 1 r 1
l kr u rj
用向后回代的方法即可求得x。设x=（x1 ，x2， · · · xn) T, y=（y1, y2, · · · yn) T,b= （b1 ，b2， · · · bn) T, 则有计算公式
y b 1 1 i 1 （4.2.5） y b l y 1 , 2 ,..., n i i ir r ,i r 1
(4.2.8)
利用(4.2.7)和(4.2.7)可得
u1 b 1 n li ai / ui1, i 2,3,... u b l c , i 2,3,... n i i i1 i
(4.2.9)
由此可求得L和U的所有元素.。解原方程组Ax=b可分为两步Ly=d 和Ux=y,计算公式为
由于方车程组的右端参与了消元计算，所以Ly=Pb的解为y=b（3）= （20，14/3，216/39） T 。解Ux=y得x=（1，2，3） T
4.2.2
三对角方程组的追赶法
b a1 A c1 b c2 a n1 bn1 an c n1 bn
(k ) akk
uk 1,k 1 l k , k 1 l n , k 1
(k ) ank
u1n u2 n uk 1,n (k ) akn (k ) ann
研究生数值分(8)直接三角分解法
u1i (i 1, 2, , n) 与L的第1列元素 li1 (i 2,3, , n)
(b) 对k+2,3,…,n 按计算公式(3),(4)依次
计算U的第k行元素 uki (i k, k 1, , n) 与L的第
k列元素 lik (i k 1, , n; k n)
20 求解三角形方程组LY=b，即按计算公
(i k, k 1, , n) （3）
k 1
lik (aik liju jk ) / ukk j 1
(i k 1, , n; k n) （4）
在我们利用杜利特尔矩阵分解解线性方程 组AX=b时,只要实现矩阵分解A=LU,依次解三角 形方程组LY=b与UX=Y即可。
计算公式：
y1
yk
对那些明确是1或是0的元素不再求。 由矩阵乘法规则与相等条件,
利用 aij 在上述计算过程中，
导出计算 lij 或 uij 的公式。
例如
第一步计算由 ai1 li1u11 得
u1i a1i (i 1,2, ,n)
第二步计算由 a1i u1i 得 li1 ai1 / u11 (i 2,3, ,n)
, n 1)
因此有 1 c1 / a1且0 1 1 由 a2 b2 a21 b2 a2 1 b2 a2 c2 0 有 2 c2 / a2且0 2 1
一般地，用归纳法可以证明
ai ci 0 （即0 i 1） (i 1, 2, , n 1)
因此我们从关系式（2）解出待定系数为
5 3 2, 2 3 5
3
2
3
4
b 7
1
0
2、用追赶法求方程组的解
4 1 0 0 x1 3
1
4
1
(b) 对k+2,3,…,n 按计算公式(3),(4)依次
计算U的第k行元素 uki (i k, k 1, , n) 与L的第
k列元素 lik (i k 1, , n; k n)
20 求解三角形方程组LY=b，即按计算公
(i k, k 1, , n) （3）
k 1
lik (aik liju jk ) / ukk j 1
(i k 1, , n; k n) （4）
在我们利用杜利特尔矩阵分解解线性方程 组AX=b时,只要实现矩阵分解A=LU,依次解三角 形方程组LY=b与UX=Y即可。
计算公式：
y1
yk
对那些明确是1或是0的元素不再求。 由矩阵乘法规则与相等条件,
利用 aij 在上述计算过程中，
导出计算 lij 或 uij 的公式。
例如
第一步计算由 ai1 li1u11 得
u1i a1i (i 1,2, ,n)
第二步计算由 a1i u1i 得 li1 ai1 / u11 (i 2,3, ,n)
, n 1)
因此有 1 c1 / a1且0 1 1 由 a2 b2 a21 b2 a2 1 b2 a2 c2 0 有 2 c2 / a2且0 2 1
一般地，用归纳法可以证明
ai ci 0 （即0 i 1） (i 1, 2, , n 1)
因此我们从关系式（2）解出待定系数为
5 3 2, 2 3 5
3
2
3
4
b 7
1
0
2、用追赶法求方程组的解
4 1 0 0 x1 3
1
4
1
数值分析ch05三角分解法
17
计算 1 2 n1 及 y1 y2 yn 的过程,称为追的过程，计算方程组的解 xn xn1 x1 的过程称为赶的过程，因此上述
方法称为追赶法。
计算量：5n-4
18
追赶法Matlab程序
function pursue[A,b] [n,m]=size(A); a=zeros(1,n); a(2:n)=diag(A,-1); c=diag(A,1); b=diag(A); if b(1)==0 error('主对角元素不能为0'); return; end alpha(1)=b(1); beta(1)=c(1)/b(1);
j 1 j i
n
15
追赶法
对角占优的三对角矩阵的 LU 分解
b1 c1 1 1 1 a2 a 1 2 2 c n 1 n an n an bn 1
第五章
解线性方程组的直接方法 —— 矩阵三角分解法
1
直接三角分解法
Matrix Factorization method
直接三角分解的基本思想 如果将线性方程组Ax=b的系数矩阵A分解成两个 三角矩阵L和U，即A=LU（也称为LU分解）， 则：
Ax b
LUx b
Ly b Ux y
14
定理：若A为对角占优三对角阵，且满足
| b1 | | c1 | 0, | bn | | an | 0, ai 0 , ci 0
| bi || ai | | ci | , ai ci 0 i 2 , , n 1
则方程组有唯一的LU分解。
| aii | | aij |
计算 1 2 n1 及 y1 y2 yn 的过程,称为追的过程，计算方程组的解 xn xn1 x1 的过程称为赶的过程，因此上述
方法称为追赶法。
计算量：5n-4
18
追赶法Matlab程序
function pursue[A,b] [n,m]=size(A); a=zeros(1,n); a(2:n)=diag(A,-1); c=diag(A,1); b=diag(A); if b(1)==0 error('主对角元素不能为0'); return; end alpha(1)=b(1); beta(1)=c(1)/b(1);
j 1 j i
n
15
追赶法
对角占优的三对角矩阵的 LU 分解
b1 c1 1 1 1 a2 a 1 2 2 c n 1 n an n an bn 1
第五章
解线性方程组的直接方法 —— 矩阵三角分解法
1
直接三角分解法
Matrix Factorization method
直接三角分解的基本思想 如果将线性方程组Ax=b的系数矩阵A分解成两个 三角矩阵L和U，即A=LU（也称为LU分解）， 则：
Ax b
LUx b
Ly b Ux y
14
定理：若A为对角占优三对角阵，且满足
| b1 | | c1 | 0, | bn | | an | 0, ai 0 , ci 0
| bi || ai | | ci | , ai ci 0 i 2 , , n 1
则方程组有唯一的LU分解。
| aii | | aij |
第一节三角形方程组和三角分解-PPT
L1A 0
3
6
0 6 11
然后再计算Gauss变换 L2 使得L2 (L1A) 的第2列的最 后一个元素为0，再由命题1和4得
上页 下页 返回 结束
1 0 0 L2 0 1 0
0 2 1
且 由命题2
1 4 7
L2 (L1A) 0 3 6
0 0 1
1
1
L11 2 1
,
L21
0
1
,
3 0 1
命题1 对于一个给定的向量 x (x1, x2, 若 xk 0 ，取
lk (0, , 0, lk1,k , , lnk )T ,
其中：
lik
xi xk
,i
k
1,
, n.
则Gauss变换：Lk In lkekT 使
Lk x (x1, , xk , 0, , 0)T .
证明：由于 Lk x (x1, , xk , xk1 xklk1,k ,
第一节三角形方程组和三角分解
如何利用电子计算机来快速、有效地求解线 性方程组得问题就是数值线性代数研究得核心问 题,而且也就是目前仍在继续研究得重大课题之一。 这就是因为各种各样得科学与工程问题往往最终 都要归结为一个线性方程组得求解问题。例如结 构分析、网络分析、大地测量、数据分析、最优 化及非线性方程组与微分方程组数值解等,都常遇 到线性方程组得求解问题。
(1)求解中小规模线性方程组(即阶数不要太 高,例如不超过1000)最常用得方法;
(2)一般用于系数矩阵稠密(即矩阵得绝大多 数元素都就是非零得)而又没有任何特殊结构 得线性方程组。
如若系数矩阵具有某种特殊形式,则为了尽可 能地减少计算量与存储量,需采用其她专门得 方法来求解。
05三角分解法
教学重点及难点： 杜里特尔分解法的计算过程
教 学 基 本 内 容
教学方法、 手段 及时间设计
复习： 直接法： 在不考虑舍入误差的情况下， 通过有限步四则运算可以求 得准确解的方法。
a 11 a 21 A a n1
a 12 a 22 an2
a 1n x1 b1 a 2n b x 2 ，b 2 ，x a nn xn b n
2 10 44 190 2 1 8 i4 1 18 1 190 1 2 2 3 7 3 6 6 6 4 12 24 24
思考 3：A 矩阵发 生变化后，b 向量 如何变化？（b 向 量可以看做 A 矩 阵的增广矩阵的 最后一列，此时公 式中需发生变化）
区别：先行后列，先列后行 杜里特尔分解和克洛特分解法都可以选主元 下面给出按列选主元 Corut 分解方法解 n 阶线性方程组的步骤 Ax=b，其中 A a ij n n , b a1, n 1 , , a n , n 1 输入：方程组的阶数 n；增广矩阵[A,b] 输出：方程组的解 x1,…,xn 或系数矩阵奇异的信息 Step1：对 i=1,2,…,n 做 Step2~5 Step2：对 j=i,…,n 求 a ji a ji
22三角分解法221杜里特尔分解法高斯消去法的第一步等价于用单位下三角阵????????????????????1001011131211?????nllll左乘方程组axb即l1axl1b记为a1xb1继续下去直到第n1步即得到????????????????2222211112111121nnnnnnaaaaaaalll?????????????????????2211121nnnbbbblll??即an1xbn1注意到an1为上三角阵将其记为r于是rallln??121?则rlllan111211??????容易验证此式即高斯消去法的消去结果3?????????????????????1111321323121111211??????nnnnlllllllll令111211?????nllll?则l是单位下三角矩阵且alr高斯消去法其实就是将矩阵a分解成了单位下三角矩阵l和上三角矩阵r??????????????????????????????????????nnnnnnnnnnnnrrrrrrllllraaaaaaaaaa???????????????222112112121212222111211111如果我们能通过某种方法直接将a矩阵进行lr分解就可以得到一种新的分解法这就是我们要讲的杜里特尔分解法
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
1
8/30/2019
6
平方根法（Cholesky分解）
续1
AT ALT D R T D R T L T R T ( D T ) L ( D ) R
由Doolittle分解的唯一性有
R T L DLT DR
（D可逆）
L R
9
平方根法（Cholesky分解）
k1
aikk1lim lk m
l11 l21 ln1
l22
ln1
lnn
lnn
第一步 ： a11l121l11 a11
ai1 l11li1li1 ai1/l11
i2,3 n
设L前k-1列元素已求出，则 第k步
n
k1
ak k lk m lk m lk2mlk2k
续2
L LD
这时 L 为一般的下三角矩阵，故 ALLT，若 L 的对角 元全为正时，由Doolittle分解的唯一性及上述分解 的推理过程，可以得到Cholesky分解的唯一性。
8/30/2019
8
平方根法（Cholesky分解）： 分解公式
l11
Al21 l22
8/30/2019
5
平方根法（Cholesky分解） 定理证明
证明：因为 A对称正定，故其顺序主式 k0 k 1 ,2 , n，
1
u11 u1n
Al21
ln1 1
m1
m1
k1
lkk akk lk2m m1
i k n
a ik lim l km m 1
k 1
n
lim lkm lik lkk
lim l km
m 1
m k 1
k 1
lim l km lik l kk m 1
8/30/2019
进行直接三角分解即可。
Reamrk：上述直接三角分解法所对应的是Gauss 顺序消去法，二者的乘除运算次数是相当的。实 际中对阶数较高的线性方程组，应采用选主元的 三角分解法求解，以保证计算结果的可靠性。
8/30/2019
2
3.3 列主元直接三角分解法
设对 Ab的分解已完成k-1步，即L 的前k-1列,U之
定义：设A为n阶n2 对称正定矩阵，L是非奇异下 三角矩阵，称 ALLT为矩阵A的Cholesky分解。
定理: n阶n2对称正定矩阵A一定存在分解: 1ALDTL
其中L为单位下三角阵，D是对角元全为正的对角 阵且这种分解式唯一;
2ALT L
其中L为下三角阵，当限定L的对角元为正时的分 解式唯一（Cholesky分解）.
•三角分解解线性方程组的公式
Axb的求解过程为:
L
y
b
U x y
可推导求解单位下三角形方程组
Ly
b的递归公式为
:
k1
y1 b1 yk bk lkmym k2,3,,n m1
求解上三角形方程组Uxy的递归公式为:
xn yn/unn
的前k-1行已求出：
u11 u12 u1,k 1
u1k u1,n1
l21
u22 u2,k 1
u2k u2,n1
l l u u u k 1,1 k 1,2
k 1,k 1
k 1,k
k 1,n1
lk1
lk 2 lk ,k 1
akk ak ,n1
RT LT
即 ALDTL,由Doolittle分解的唯一性，及 UDR 的分解过程可知该分解式的唯一性。
8/30/2019
7
平方根法（Cholesky分解）
改写
DDD
u11
u1
unn
unn
则 A L D D L T L D L D T L L T
liksik/sk,k ik1, ,n,
之后按原k公 行u式 中 的计 其算 它 uk,k1 第 元 , ,uk素 ,n1
若 sikk maxsik,则需将第k行与第ik行完全交换，这样 满足前面情形，按第一种情形实施计算。
3.4 平方根法（Cholesky分解）
1. 正定矩阵的Cholesky分解
yk
n
ukmxm
xk
mk1
ukk
kn 1 ,n 2 , ,1
8/30/2019
1
直接三角分解法
续
k1
ukj akj lkmum,j m1
k1
yk bk lkmym m1
发算式现取y的计 代求算 。方事y k程实的组上规的，律过这与程只计可要算用把u kj三的b 角做规分为律解A相的的同第紧,因n凑＋此格1计列
ln1
ln2 ln,k 1
ank an,n1
第k步计算：先选主元，再计算列，最后计算行
8/30/2019
3
u
kj
a kj
k 1
lkm u mj
m 1
aik
k 1
lim u mk
l
ik
m 1
unn
因为 k u 1 u 2 1 u 2 k 0 k k 1 , 2 , n u k 0 k
1 Al21
ln1
u11
ln2 1
u22
1 unn
u12 u11 1
u1nu11 u2nu22 LDR
u kk
jk,k1 ,n1 参选主元量
ik1,n
sik aik k1 limumk
m1
i k, k 1,n.
若 skkmsai k,x则无需 skk 交 akk 换 k 1lkm u ， m ku 且 kk
m 1
8/30/2019
6
平方根法（Cholesky分解）
续1
AT ALT D R T D R T L T R T ( D T ) L ( D ) R
由Doolittle分解的唯一性有
R T L DLT DR
（D可逆）
L R
9
平方根法（Cholesky分解）
k1
aikk1lim lk m
l11 l21 ln1
l22
ln1
lnn
lnn
第一步 ： a11l121l11 a11
ai1 l11li1li1 ai1/l11
i2,3 n
设L前k-1列元素已求出，则 第k步
n
k1
ak k lk m lk m lk2mlk2k
续2
L LD
这时 L 为一般的下三角矩阵，故 ALLT，若 L 的对角 元全为正时，由Doolittle分解的唯一性及上述分解 的推理过程，可以得到Cholesky分解的唯一性。
8/30/2019
8
平方根法（Cholesky分解）： 分解公式
l11
Al21 l22
8/30/2019
5
平方根法（Cholesky分解） 定理证明
证明：因为 A对称正定，故其顺序主式 k0 k 1 ,2 , n，
1
u11 u1n
Al21
ln1 1
m1
m1
k1
lkk akk lk2m m1
i k n
a ik lim l km m 1
k 1
n
lim lkm lik lkk
lim l km
m 1
m k 1
k 1
lim l km lik l kk m 1
8/30/2019
进行直接三角分解即可。
Reamrk：上述直接三角分解法所对应的是Gauss 顺序消去法，二者的乘除运算次数是相当的。实 际中对阶数较高的线性方程组，应采用选主元的 三角分解法求解，以保证计算结果的可靠性。
8/30/2019
2
3.3 列主元直接三角分解法
设对 Ab的分解已完成k-1步，即L 的前k-1列,U之
定义：设A为n阶n2 对称正定矩阵，L是非奇异下 三角矩阵，称 ALLT为矩阵A的Cholesky分解。
定理: n阶n2对称正定矩阵A一定存在分解: 1ALDTL
其中L为单位下三角阵，D是对角元全为正的对角 阵且这种分解式唯一;
2ALT L
其中L为下三角阵，当限定L的对角元为正时的分 解式唯一（Cholesky分解）.
•三角分解解线性方程组的公式
Axb的求解过程为:
L
y
b
U x y
可推导求解单位下三角形方程组
Ly
b的递归公式为
:
k1
y1 b1 yk bk lkmym k2,3,,n m1
求解上三角形方程组Uxy的递归公式为:
xn yn/unn
的前k-1行已求出：
u11 u12 u1,k 1
u1k u1,n1
l21
u22 u2,k 1
u2k u2,n1
l l u u u k 1,1 k 1,2
k 1,k 1
k 1,k
k 1,n1
lk1
lk 2 lk ,k 1
akk ak ,n1
RT LT
即 ALDTL,由Doolittle分解的唯一性，及 UDR 的分解过程可知该分解式的唯一性。
8/30/2019
7
平方根法（Cholesky分解）
改写
DDD
u11
u1
unn
unn
则 A L D D L T L D L D T L L T
liksik/sk,k ik1, ,n,
之后按原k公 行u式 中 的计 其算 它 uk,k1 第 元 , ,uk素 ,n1
若 sikk maxsik,则需将第k行与第ik行完全交换，这样 满足前面情形，按第一种情形实施计算。
3.4 平方根法（Cholesky分解）
1. 正定矩阵的Cholesky分解
yk
n
ukmxm
xk
mk1
ukk
kn 1 ,n 2 , ,1
8/30/2019
1
直接三角分解法
续
k1
ukj akj lkmum,j m1
k1
yk bk lkmym m1
发算式现取y的计 代求算 。方事y k程实的组上规的，律过这与程只计可要算用把u kj三的b 角做规分为律解A相的的同第紧,因n凑＋此格1计列
ln1
ln2 ln,k 1
ank an,n1
第k步计算：先选主元，再计算列，最后计算行
8/30/2019
3
u
kj
a kj
k 1
lkm u mj
m 1
aik
k 1
lim u mk
l
ik
m 1
unn
因为 k u 1 u 2 1 u 2 k 0 k k 1 , 2 , n u k 0 k
1 Al21
ln1
u11
ln2 1
u22
1 unn
u12 u11 1
u1nu11 u2nu22 LDR
u kk
jk,k1 ,n1 参选主元量
ik1,n
sik aik k1 limumk
m1
i k, k 1,n.
若 skkmsai k,x则无需 skk 交 akk 换 k 1lkm u ， m ku 且 kk
m 1