《数值分析与算法》 第三讲-线性方程组的直接解法

非奇异矩阵, 线性方程组解的存在性与唯一性
Wenjian Yu 4
线性方程组与矩阵的基本概念

表示稀疏矩阵的Wilkinson图
Wenjian Yu
5
向量的范数

(正定条件)
(三角不等式)
Wenjian Yu
6
向量的范数

(曼哈顿范数)
(欧氏范数=内积范数)
(“最大”范数)
观看演示3.1
Wenjian Yu 7
美国TAMU, Prof. Davis, 2004, http://faculty.cse.tamu.edu/davis/
Wenjian Yu 61
Wenjian Yu
62
非线性方程组求解
Wenjian Yu
63
非线性方程组的求解

最大特征值<1
Wenjian Yu
牛顿法解非线性方程组

交互演示2.6

存储方式
稠密矩阵(二维数组) 、稀疏矩阵(压缩稀疏列) 以相同的方式支持一般的矩阵操作/运算; 通过whos
命令看两者区别

稠密矩阵的命令
生成矩阵:
ones(m,n), zeros(m,n), eye(m,n), rand(m,n) 矩阵维数等: size(A), diag(A)/diag(v), tril(A), triu(A) 载入/导出数据文件 (ascii或binary格式): load, save
43
应用实例

R4
V3
R5 R3
V1
R1
V2
R2
is1
1 1 R R 1 R 1 R
1

2
1 R4
2
1 R2 1 1 R2 R3
4
0
V1 0 V 0 0 2 1 1 i V s 1 3 R R
27
矩阵的LU分解

算法3.5 用高斯消去过程做LU分解
乘子的相反数
Wenjian Yu
28
矩阵的LU分解

Wenjian Yu
29
直接LU分解算法 (算法3.6)
u
l u
l
l
u
• 数学上等价算法3.5 • 算法描述不同，存 取、计算的顺序不同 • 对于稠密矩阵(二维 数组)计算效率一样， 但对于稀疏矩阵(特殊 数据结构)，效率可能 有较大差别

a1n a2 n a3n ann
Wenjian Yu
21
矩阵的LU分解
Wenjian Yu
22
高斯消去法

Wenjian Yu
23
矩阵的LU分解

高斯消去过程对矩阵A的变换 (3阶矩阵为例)
行倍加变换
等价于左乘消去矩阵:
乘子
两次倍加变换 单位下三角阵
向量有关的定理

(Hint: 利用范数，课后思考)
Wenjian Yu
8
矩阵的范数

Байду номын сангаас
(不满足交换律)
Wenjian Yu
9
矩阵的范数

后面只考虑 算子范数
0
1
思考: 如何证明？
Wenjian Yu
10
矩阵的范数

Wenjian Yu
11
问题的敏感性
Wenjian Yu
12
线性方程组求解问题的敏感性
压缩稀疏行 (CSR格式)
aa 1 2 3 4 5 6 7 8 9 10 11 col 1 3 1 2 4 1 3 5 2 4 5 prow 1 3 6 9 11 12 非零元按第1行, 第2行, …, 第n行的顺序存储 prow为各行在数组中的开始位置 如何得每行非零元数目?
Wenjian Yu 55
1
2 max lij
m ax a ii
1
(Matlab演示)
48
Wenjian Yu
带状矩阵解法与 稀疏矩阵简介
Wenjian Yu
49
带状线性方程组

Wenjian Yu
50
带状线性方程组

计算量: M.D.= 2n-2
Wenjian Yu
51
带状线性方程组

且至少有一个不等号严格成立
稀疏矩阵

广泛用于非结 构化稀疏阵
Wenjian Yu
56
稀疏矩阵

思考: 按稀疏矩阵存 储格式写相应程序
X X 0 0
Wenjian Yu
0 X 0 X X X X 0 X X 0 X 0 X X
57
X
X
X
Matlab中的相关技术
Wenjian Yu
58
Matlab中的矩阵
30
Wenjian Yu
LU分解的用途

也可求矩阵的逆
Wenjian Yu
31
选主元技术与稳定性
Wenjian Yu
32
主元为零的可能性

单位下三角阵
Wenjian Yu
33
主元为零的可能性

无LU分解 高斯消去中断
？
？ ？
无穷多LU分解
思考：这三条性质 可以互相推出吗？
34
矩阵奇异
矩阵不可LU分解

Wenjian Yu
13
矩阵的条件数

Wenjian Yu
14
矩阵的条件数

符合上界的意义
Wenjian Yu 15
矩阵的条件数

良态问题！
见pp. 68
Wenjian Yu
16
矩阵条件数的性质

0
1
(演示3.2)
Wenjian Yu 17
矩阵条件数的性质

p-范数下等号成立!
正交变换不改条件数!
数 值 分 析 (3)
Numerical Analysis
第三章 线性方程组的直接解法

线性方程组求解是一个基本的数学问题
新的挑战:
方程的规模大，例如变量数>106 困难: 存储、计算速度、舍入误差 是数值线性代数的基础 (我们分两章讨论)

本章内容
基本概念与问题敏感性
高斯消去法
Wenjian Yu
52
带状线性方程组

Wenjian Yu
53
稀疏矩阵

非零元 分布图
Wenjian Yu
54
稀疏矩阵

存储稀疏矩阵的数据结构
三元组
(COO格式)
aa 1 2 3 4 5 6 7 8 9 10 11 row 1 1 2 2 2 3 3 3 4 4 5 col 1 3 1 2 4 1 3 5 2 4 5 非零元在数组中顺序可以任意 有点冗余: 一些连续存储单元可能有相同的行(列)编号
k 1
j 1
End End
Wenjian Yu
47
Cholesky分解算法

T A L0U L0D 2D 2LT LL 0 1 1
增长因子 ρ
m ax U T (i, j ) m ax a ij



m ax l ij l jj
ij
m ax a
U D 2LT 1 T U LD 2
矩阵的LU分解 选主元技术与算法稳定性
直接解法，就是理论上 经过有限步计算能得到 准确解的方法。
对称正定矩阵与带状矩阵解法
稀疏矩阵与实用技术
(选学)
2
Wenjian Yu
基本概念
Wenjian Yu
3
线性方程组与矩阵的基本概念
复习: 线性空间, 线性相关(无关), 基, 维数, 内积；
Wenjian Yu
牛顿法解非线性方程组

解方程
解方程
Wenjian Yu
稀疏的对角阵 较稠密的带状方阵 上三角或下三角矩阵 对三角矩阵作行排列 形成的矩阵 对称矩阵, 且对角线元 素大于零
求解算法
右端项元素除以矩阵对角元 部分主元的带状矩阵LU分解 (LAPACK软件包) 回代法或前代法 重排序后用回代或前代法
尝试Cholesky分解算法(稠密矩阵: LAPACK, 稀 疏矩阵: CHOLMOD), 若稠密、且不正定使用 选主元的对称矩阵求解算法(LAPACK) 稠密上Hessenberg矩阵 消为上三角矩阵再回代求解 针对稀疏矩阵的直接解法(UMFPACK) 一般的稀疏方阵 部分主元LU分解(LAPACK) 一般的稠密方阵 通过矩阵的QR分解得到最小二乘解 不是方阵
18
Wenjian Yu
高斯消去法
Wenjian Yu
19
高斯消去法

行倍乘变换 行倍加变换
(此小节不要求)
Wenjian Yu
20
高斯消去法

算法3.1
a11 a21 a 31 解线性方程组的高斯消去过程 an1
a12 a22 a32 an 2
a13 a23 a33 an3
Wenjian Yu
24
矩阵的LU分解

第2步消去
行倍加变换
等价于左乘消去矩阵
一次倍加变换
Wenjian Yu
25
矩阵的LU分解
逆矩阵也是消 去矩阵
注: 矩阵的下标指 示了消去矩阵的 类型数(定义3.12)
定理 3.13
Wenjian Yu
26
矩阵的LU分解

(充要条件)
Wenjian Yu
1 R4
4 5
0
Wenjian Yu
44
对称矩阵的LU分解

希望快速求解对称 或对称正定矩阵
Wenjian Yu
45
对称矩阵的LU分解

Wenjian Yu
46
Cholesky分解算法

也叫平方根法. 下面按直接LU分解思想推导算法
l 11 0 L l 21 l 22 O M M O l l L n1 n2 0 l 11 M 0 0 M l nn 0 l 21 L l n 1 l 22 L l n 2 O O M L 0 l nn
行交换
消去
验证
Wenjian Yu
39
采用部分主元的LU分解算法 (算法3.7)
(Matlab命令: lu)
Wenjian Yu 40
部分主元技术与其他主元技术

(见演示3.5)
41
Wenjian Yu
算法的稳定性

(思考题)
Wenjian Yu
42
对称正定矩阵的 Cholesky分解
Wenjian Yu
Wenjian Yu
如何解决主元为零的问题

定理 3.16
Wenjian Yu
35
使用部分主元技术的LU分解

Wenjian Yu
36
使用部分主元技术的LU分解

Wenjian Yu
37
使用部分主元技术的LU分解

Wenjian Yu
38
使用部分主元技术的LU分解

行交换 消去
(见演示3.4)

稀疏矩阵的常用命令
生成与转换: sparse(X) or sparse(i,j,s,m,n), X = full(A)
非零元的信息:
nnz(A), spy(A), find(A)
Wenjian Yu
59
Matlab中的矩阵

Wenjian Yu
60
Matlab“\”的内部算法
矩阵A (按如下顺序)
=
算法3.10 (按从1到n的顺序逐列算出L的元素值)
For j =1 to n
a jj
2 a jj l j21 L l jj
a jj
2 a jk k 1
j 1
{对角元}
{当前列}
For i = j+1 to n
a ij (a ij a ika jk ) / a jj