§ 2 三角分解法 /* Matrix Factorization */
高斯消元法的矩阵形式 /* Matrix Form of G.E,*/：
Step 1,)0(/
111111 aaam ii
记 L1 =
1
..
.
1
1
1
21
nm
m
,则
][ )1()1(1 bAL?
)1(1)1(1)1(11,.,baa n
)2(A )2(b?
Step n? 1:
)(
)2(
2
)1(
1
)(
)2(
2
)2(
22
)1(
1
)1(
12
)1(
11
121
..
.
..
....
...
...
...
n
n
n
nn
n
n
nn
b
b
b
a
aa
aaa
bALLL
其中
Lk =
1
..
.
1
1
,
,1
kn
kk
m
m
§ 2 Matrix Factorization – Matrix Form of G.E.
1kL
1
..
.
1
1
,
,1
kn
kk
m
m
1 11211,.,nLLL
1
1
1
jim,
记为 L
单位下三角阵
/* unitary lower-triangular matrix */
记 U =
)(
)2(
2
)2(
22
)1(
1
)1(
12
)1(
11
..
..,,
.,,
.,,
n
nn
n
n
a
aa
aaa
LUA
A 的 LU 分解
/* LU factorization */
Hey hasn’t GE given me
enough headache? Why
do I have to know its
matrix form!
When you have to
solve the system for
different with a fixed
A.
b?
Could you be more
specific,please?
Factorize A first,then
for every you nly have to
solve two simple triangular
systems and
.
b?
byL
yxU
§ 2 Matrix Factorization – Matrix Form of G.E.
定理 若 A的所有顺序主子式 /* determinant of leading
principal submatrices */ 均不为 0，则 A 的 LU 分解唯一（其中 L 为 单位 下三角阵）。
证明,由 § 1中定理可知,LU 分解存在。下面证明唯一性。
若不唯一，则可设 A = L1U1 = L2U2,推出
121 UU 211122211 LLUULL
Upper-triangular
Lower-triangular
With diagonal
entries 1
I
注,L 为一般下三角阵而 U 为 单位 上三角阵的分解称为
Crout 分解 。
实际上只要考虑 A* 的 LU 分解，即,则即是 A 的 Crout 分解。
ULA ~~*?
*~*~ LUA?
§ 2 Matrix Factorization – Doolittle
道立特分解法 /* Doolittle Factorization */：
—— LU 分解的紧凑格式 /* compact form */
反复计算,
很浪费哦 ……
通过比较法直接导出 L 和 U 的计算公式。思路
nn
n
nnnn
n
u
uu
l
l
aa
aa
..
.
..
....
...
1...
......
1
1
...
..
.
..
.
..
.
..
.
..,111
1
21
1
111
),m i n (
1
ji
k
jkki ul
jia
§ 2 Matrix Factorization – Doolittle
),m i n( 1 jik jkkiji ula
固定 i,
对 j = i,i+1,…,n 有
ijkj
i
k
ikij uula
1
1
lii = 1
kj
i
k
ikijij ulau?
1
1
a
将 i,j 对换，对 j = i,i+1,…,n 有
iijiki
i
k
jkji ulula
1
1
ii
i
k
kijkjiji uulal /)(
1
1
b
Algorithm,Doolittle Factorization
Step 1,u1j = a1j; lj1 = aj1 / u11; ( j = 1,…,n )
Step 2,compute and for i = 2,…,n?1;
Step 3,
a b
11nk knnknnnn ulau
一般采用 列主元法增强稳定性。但注意也必须做相应的行交换。
b?
§ 2 Matrix Factorization – Choleski
平方根法 /* Choleski’s Method */：
——对称 /* symmetric */ 正定 /* positive definite */
矩阵的分解法定义 一个矩阵 A = ( aij )n?n 称为 对称阵,如果 aij = aji 。
定义 一个矩阵 A 称为 正定阵,如果 对任意非零向量 都成立。
0?xAx T
x?
回顾,对称正定阵的几个重要性质
A?1 亦对称正定，且 aii > 0?
若不然，则 0xA 存在非零解，即
0?xAx T 存在非零解。
1111 )()(, AAIAAIAA TTT
对任意,存在,使得,
即 。
0x 0y xyA
xAy 1 011 yAyyAAAyxAx TTT
其中0 xAxa Tii Tx )0...1...0(
第 i 位
A 的顺序主子阵 /* leading principal submatrices */ Ak 亦对称正定对称性显然。对任意 有
,其中 。
kk Rx 0
0 xAxxAx TkkTk nk Rxx
0
0
A 的特征值 /* eigen value */?i > 0
设对应特征值? 的非零特征向量为,则 。 20 xxxxAx TTx?
A 的全部顺序主子式 det ( Ak ) > 0
因为?
n
i
iA
1
)de t (?
§ 2 Matrix Factorization – Choleski
将 对称 正定阵 A 做 LU 分解
U = uij =
u11 u
ij / uii
1 1
1
u22
unn
记为
UD~
A 对称 TUL ~? 即 TLDLA?
记 D1/2 =
11u
22u
nnu
Why is uii > 0?Since det(Ak) >
0
2/1~ LDL?则 仍是下三角阵
TLLA ~~?
nnRL
定理 设矩阵 A对称正定，则存在非奇异下三角阵使得 。若限定 L 对角元为正，则分解唯一。TLLA?
注,对于对称正定阵 A,从 可知对任意 k? i
有 。 即 L 的元素不会增大，误差可控，不需选主元。
ik ikii la 1 2
iiik al?||
§ 2 Matrix Factorization – Choleski
Algorithm,Choleski’s Method
To factor the symmetric positive definite n?n matrix A into LLT,where L
is lower triangular.
Input,the dimension n; entries aij for 1? i,j? n of A.
Output,the entries lij for 1? j? i and 1? i? n of L.
Step 1 Set ;
Step 2 For j = 2,…,n,set ;
Step 3 For i = 2,…,n?1,do steps 4 and 5
Step 4 Set ;
Step 5 For j = i+1,…,n,set ;
Step 6 Set ;
Step 7 Output ( lij for j = 1,…,i and i = 1,…,n );
STOP,
1111 al?
1111 / lal jj?
11 2ik ikiiii lal
iiik ikjkjiji lllal 11
11 2nk nknnnn lal
因为 A 对称，所以只需存半个 A，即其中
nnn aaaaannA,..,,...,,,,2/)1( 1222111
jiiAa ij 2/)1(
运算量为 O( 3/6),比普通 LU
分解少一半，但有 n 次开方。用 A = LDLT
分解，可省开方时间 (p.50-51)。
HW,p.54
#2,#5,#6
§ 2 Matrix Factorization – Tridiagonal System
追赶法解 三对角 方程组
/* Crout Reduction for Tridiagonal Linear System */
nnnn
nnn
f
f
f
x
x
x
ba
cba
cba
cb
2
1
2
1
111
222
11
Step 1,对 A 作 Crout 分解
1
1
1
1
2
1
n
nn
A
直接比较等式两边的元素，可得到计算公式 （ p.52） 。
Step 2,追 ——即解,fyL,
1
11?fy? ),.,,,2()( 1 niyrfy
i
iiii
Step 3,赶 ——即解,yxU )1,.,,,1(,1 nixyxyx iiiinn?
与 G.E.类似，一旦
i = 0 则算法中断，故并非任何三对角阵都可以用此方法分解。
§ 2 Matrix Factorization – Tridiagonal System
定理 若 A 为 对角占优 /* diagonally dominant */ 的三对角阵，且满足,则追赶法可解以 A 为系数矩阵的方程组。
0,0,0||||,0|||| 11 iinn caabcb
Hey,what does diagonally
dominant meanIt means that the diagonal entries of the matrix are very LARGE.
Well,how large is
LARGE?
They satisfy the following inequality:
ij ijii
aa ||||
注,
如果 A 是 严格 对角占优阵，则不要求三对角线上的所有元素非零。
根据不等式可知：分解过程中，矩阵元素不会过分增大，算法保证稳定。
运算量为 O(6n)。
||||||||||,1|| 1 iiiiiiii abbab
§ 2 Matrix Factorization – Tridiagonal System
Lab 06,Crout Reduction for Tridiagonal Linear Systems
Apply Crout Reduction to solve a given n× n tridiagonal linear
system
Input
There are several sets of inputs,For each set:
The 1st line contains an integer 100? n? 0 which is the size of a
matrix,n =?1 signals the end of file.
The 2nd line contains n?1 real numbers,
The 3rd line contains n real numbers,
The 4th line contains n?1 real numbers,
The 5th line contains n real numbers,
The numbers are separated by spaces.
naaa,..32
nnnn
nnn
f
f
f
x
x
x
ba
cba
cba
cb
2
1
2
1
111
222
11
nbbb,..21
121,.,?nccc
nfff,..21
§ 2 Matrix Factorization – Tridiagonal System
Output
Each entry of the solution is to be printed as in the C fprintf:
fprintf(outfile,"%16.8e\n",x );
If the method fails to give a solution,print the message
,The?Crout?method?failed.\n”,/* here?represents a space*/
The outputs of two test cases must be seperated by a blank line.
Sample Input
(?represents a space)
5
–1?–1?–1?–1
4?4?4?4?4
–1?–1?–1?–1
100?200?200?200?100
3
–1?–1
0?3?3
–0.5?–1
8?8?7
–1
Sample Output
(?represents a space)
4.61538462e+001
8.46153846e+001
9.23076923e+001
8.46153846e+001
4.61538462e+001
The Crout method failed.
高斯消元法的矩阵形式 /* Matrix Form of G.E,*/：
Step 1,)0(/
111111 aaam ii
记 L1 =
1
..
.
1
1
1
21
nm
m
,则
][ )1()1(1 bAL?
)1(1)1(1)1(11,.,baa n
)2(A )2(b?
Step n? 1:
)(
)2(
2
)1(
1
)(
)2(
2
)2(
22
)1(
1
)1(
12
)1(
11
121
..
.
..
....
...
...
...
n
n
n
nn
n
n
nn
b
b
b
a
aa
aaa
bALLL
其中
Lk =
1
..
.
1
1
,
,1
kn
kk
m
m
§ 2 Matrix Factorization – Matrix Form of G.E.
1kL
1
..
.
1
1
,
,1
kn
kk
m
m
1 11211,.,nLLL
1
1
1
jim,
记为 L
单位下三角阵
/* unitary lower-triangular matrix */
记 U =
)(
)2(
2
)2(
22
)1(
1
)1(
12
)1(
11
..
..,,
.,,
.,,
n
nn
n
n
a
aa
aaa
LUA
A 的 LU 分解
/* LU factorization */
Hey hasn’t GE given me
enough headache? Why
do I have to know its
matrix form!
When you have to
solve the system for
different with a fixed
A.
b?
Could you be more
specific,please?
Factorize A first,then
for every you nly have to
solve two simple triangular
systems and
.
b?
byL
yxU
§ 2 Matrix Factorization – Matrix Form of G.E.
定理 若 A的所有顺序主子式 /* determinant of leading
principal submatrices */ 均不为 0，则 A 的 LU 分解唯一（其中 L 为 单位 下三角阵）。
证明,由 § 1中定理可知,LU 分解存在。下面证明唯一性。
若不唯一，则可设 A = L1U1 = L2U2,推出
121 UU 211122211 LLUULL
Upper-triangular
Lower-triangular
With diagonal
entries 1
I
注,L 为一般下三角阵而 U 为 单位 上三角阵的分解称为
Crout 分解 。
实际上只要考虑 A* 的 LU 分解，即,则即是 A 的 Crout 分解。
ULA ~~*?
*~*~ LUA?
§ 2 Matrix Factorization – Doolittle
道立特分解法 /* Doolittle Factorization */：
—— LU 分解的紧凑格式 /* compact form */
反复计算,
很浪费哦 ……
通过比较法直接导出 L 和 U 的计算公式。思路
nn
n
nnnn
n
u
uu
l
l
aa
aa
..
.
..
....
...
1...
......
1
1
...
..
.
..
.
..
.
..
.
..,111
1
21
1
111
),m i n (
1
ji
k
jkki ul
jia
§ 2 Matrix Factorization – Doolittle
),m i n( 1 jik jkkiji ula
固定 i,
对 j = i,i+1,…,n 有
ijkj
i
k
ikij uula
1
1
lii = 1
kj
i
k
ikijij ulau?
1
1
a
将 i,j 对换，对 j = i,i+1,…,n 有
iijiki
i
k
jkji ulula
1
1
ii
i
k
kijkjiji uulal /)(
1
1
b
Algorithm,Doolittle Factorization
Step 1,u1j = a1j; lj1 = aj1 / u11; ( j = 1,…,n )
Step 2,compute and for i = 2,…,n?1;
Step 3,
a b
11nk knnknnnn ulau
一般采用 列主元法增强稳定性。但注意也必须做相应的行交换。
b?
§ 2 Matrix Factorization – Choleski
平方根法 /* Choleski’s Method */：
——对称 /* symmetric */ 正定 /* positive definite */
矩阵的分解法定义 一个矩阵 A = ( aij )n?n 称为 对称阵,如果 aij = aji 。
定义 一个矩阵 A 称为 正定阵,如果 对任意非零向量 都成立。
0?xAx T
x?
回顾,对称正定阵的几个重要性质
A?1 亦对称正定，且 aii > 0?
若不然，则 0xA 存在非零解，即
0?xAx T 存在非零解。
1111 )()(, AAIAAIAA TTT
对任意,存在,使得,
即 。
0x 0y xyA
xAy 1 011 yAyyAAAyxAx TTT
其中0 xAxa Tii Tx )0...1...0(
第 i 位
A 的顺序主子阵 /* leading principal submatrices */ Ak 亦对称正定对称性显然。对任意 有
,其中 。
kk Rx 0
0 xAxxAx TkkTk nk Rxx
0
0
A 的特征值 /* eigen value */?i > 0
设对应特征值? 的非零特征向量为,则 。 20 xxxxAx TTx?
A 的全部顺序主子式 det ( Ak ) > 0
因为?
n
i
iA
1
)de t (?
§ 2 Matrix Factorization – Choleski
将 对称 正定阵 A 做 LU 分解
U = uij =
u11 u
ij / uii
1 1
1
u22
unn
记为
UD~
A 对称 TUL ~? 即 TLDLA?
记 D1/2 =
11u
22u
nnu
Why is uii > 0?Since det(Ak) >
0
2/1~ LDL?则 仍是下三角阵
TLLA ~~?
nnRL
定理 设矩阵 A对称正定，则存在非奇异下三角阵使得 。若限定 L 对角元为正，则分解唯一。TLLA?
注,对于对称正定阵 A,从 可知对任意 k? i
有 。 即 L 的元素不会增大，误差可控，不需选主元。
ik ikii la 1 2
iiik al?||
§ 2 Matrix Factorization – Choleski
Algorithm,Choleski’s Method
To factor the symmetric positive definite n?n matrix A into LLT,where L
is lower triangular.
Input,the dimension n; entries aij for 1? i,j? n of A.
Output,the entries lij for 1? j? i and 1? i? n of L.
Step 1 Set ;
Step 2 For j = 2,…,n,set ;
Step 3 For i = 2,…,n?1,do steps 4 and 5
Step 4 Set ;
Step 5 For j = i+1,…,n,set ;
Step 6 Set ;
Step 7 Output ( lij for j = 1,…,i and i = 1,…,n );
STOP,
1111 al?
1111 / lal jj?
11 2ik ikiiii lal
iiik ikjkjiji lllal 11
11 2nk nknnnn lal
因为 A 对称，所以只需存半个 A，即其中
nnn aaaaannA,..,,...,,,,2/)1( 1222111
jiiAa ij 2/)1(
运算量为 O( 3/6),比普通 LU
分解少一半，但有 n 次开方。用 A = LDLT
分解，可省开方时间 (p.50-51)。
HW,p.54
#2,#5,#6
§ 2 Matrix Factorization – Tridiagonal System
追赶法解 三对角 方程组
/* Crout Reduction for Tridiagonal Linear System */
nnnn
nnn
f
f
f
x
x
x
ba
cba
cba
cb
2
1
2
1
111
222
11
Step 1,对 A 作 Crout 分解
1
1
1
1
2
1
n
nn
A
直接比较等式两边的元素，可得到计算公式 （ p.52） 。
Step 2,追 ——即解,fyL,
1
11?fy? ),.,,,2()( 1 niyrfy
i
iiii
Step 3,赶 ——即解,yxU )1,.,,,1(,1 nixyxyx iiiinn?
与 G.E.类似，一旦
i = 0 则算法中断，故并非任何三对角阵都可以用此方法分解。
§ 2 Matrix Factorization – Tridiagonal System
定理 若 A 为 对角占优 /* diagonally dominant */ 的三对角阵，且满足,则追赶法可解以 A 为系数矩阵的方程组。
0,0,0||||,0|||| 11 iinn caabcb
Hey,what does diagonally
dominant meanIt means that the diagonal entries of the matrix are very LARGE.
Well,how large is
LARGE?
They satisfy the following inequality:
ij ijii
aa ||||
注,
如果 A 是 严格 对角占优阵，则不要求三对角线上的所有元素非零。
根据不等式可知：分解过程中，矩阵元素不会过分增大，算法保证稳定。
运算量为 O(6n)。
||||||||||,1|| 1 iiiiiiii abbab
§ 2 Matrix Factorization – Tridiagonal System
Lab 06,Crout Reduction for Tridiagonal Linear Systems
Apply Crout Reduction to solve a given n× n tridiagonal linear
system
Input
There are several sets of inputs,For each set:
The 1st line contains an integer 100? n? 0 which is the size of a
matrix,n =?1 signals the end of file.
The 2nd line contains n?1 real numbers,
The 3rd line contains n real numbers,
The 4th line contains n?1 real numbers,
The 5th line contains n real numbers,
The numbers are separated by spaces.
naaa,..32
nnnn
nnn
f
f
f
x
x
x
ba
cba
cba
cb
2
1
2
1
111
222
11
nbbb,..21
121,.,?nccc
nfff,..21
§ 2 Matrix Factorization – Tridiagonal System
Output
Each entry of the solution is to be printed as in the C fprintf:
fprintf(outfile,"%16.8e\n",x );
If the method fails to give a solution,print the message
,The?Crout?method?failed.\n”,/* here?represents a space*/
The outputs of two test cases must be seperated by a blank line.
Sample Input
(?represents a space)
5
–1?–1?–1?–1
4?4?4?4?4
–1?–1?–1?–1
100?200?200?200?100
3
–1?–1
0?3?3
–0.5?–1
8?8?7
–1
Sample Output
(?represents a space)
4.61538462e+001
8.46153846e+001
9.23076923e+001
8.46153846e+001
4.61538462e+001
The Crout method failed.
