第三章 解线性方程组的直接法
/* Direct Method for Solving Linear Systems */
求解 bxA
§ 1 高斯消元法 /* Gaussian Elimination */
高斯消元法:
思路首先将 A化为上三角阵 /* upper-triangular matrix */,
再回代求解 /* backward substitution */。
=
§ 1 Gaussian Elimination – The Method
消元 记
,)( )1()1( nnijaAA
)1(
)1(
1
)1(
...
nb
b
bb
Step 1,设,计算因子0)1(
11?a )...,,2(/ )1(11)1( 11 niaam ii
将增广矩阵 /* augmented matrix */ 第 i 行? mi1? 第 1行,
得到
)1(1)1(1)1(12)1(11,,,baaa n
)2(A )2(b? ).,,,,2,(
)1(
11
)1()2(
)1(
11
)1()2(
nji
bmbb
amaa
iii
jiijij
其中
Step k,设,计算因子且计算
0)(?kkka )...,,1(/ )()( nkiaam kkkkikik
).,,,,1,(
)()()1(
)()()1(
nkji
bmbb
amaa
k
kik
k
i
k
i
k
kjik
k
ij
k
ij
共进行? 步n? 1
)(
)2(
2
)1(
1
2
1
)(
)2(
2
)2(
22
)1(
1
)1(
12
)1(
11
..
.
..
.
..
....
...
...
n
nn
n
nn
n
n
b
b
b
x
x
x
a
aa
aaa
回代
)()( / nnnnnn abx?
0)(?nnnaWhat if?No unique solution exists.
)1.,,,,1()( 1
)()(
ni
a
xab
x i
ii
n
ij
j
i
ij
i
i
i
0)(?iiiahat if?
Then we must find the
smallest integer k? i with
,and interchange
the k-th row with the i-th
row,
0)(?ikia
What if we can’t
find such k?No unique
solution exists.
定理 若 A的所有 顺序主子式 /* determinant of leading
principal submatrices */ 均不为 0,则高斯消元无需换行即可进行到底,得到唯一解。
iii
i
i
aa
aa
A
.,,
.,,.,,.,,
.,,
)d e t (
1
111
注,事实上,只要 A 非奇异,即 A?1 存在,则可通过逐次消元及行交换,将方程组化为三角形方程组,求出唯一解。
§ 1 Gaussian Elimination – The Method
选主元消去法 /* Pivoting Strategies */
例,单精度解方程组
2
110
21
21
9
xx
xx
/* 精确解为 和 */...1000...00.1
101
1
91
x
8个
...8999...99.02 12 xx
8个用 Gaussian Elimination计算:
9112121 10/ aam
999
2122 10101010.,,0.011ma
8个
9212 1012mb
99
9
10100
1110
0,1 12 xx
小主元 /* Small pivot
element */ 可能导致计算失败。
§ 1 Gaussian Elimination – Pivoting Strategies
全主元消去法 /* Complete Pivoting */
每一步选绝对值最大的元素为主元素,保证 。1||?ikm
Step k,① 选取 ;0||m a x||
, ijnjikji aa kk
② If ik? k then 交换第 k 行与第 ik 行 ;
If jk? k then 交换第 k 列与第 jk 列 ;
③ 消元注,列交换改变了 xi 的顺序,须记录 交换次序,解完后再换回来。
列主元消去法 /* Partial Pivoting,or maximal column pivoting */
省去换列的步骤,每次仅选一列中最大的元。
0||m a x||, iknikki aa k
§ 1 Gaussian Elimination – Pivoting Strategies
例:
211
1110 9
1,1 12 xx
110
211
1110
211
9
注,列主元法没有全主元法稳定。
0,1 12 xx
例:
211 10101
99
99
99
10100
10101
注意:这两个方程组在数学上 严格等价 。
标度化列主元消去法 /* Scaled Partial Pivoting */
对每一行计算 。为省时间,si 只在初始时计算一次。以后每一步考虑子列 中 最大的 aik 为主元。
||m a x1 ijnji as
nk
kk
a
a
..,i
ik
s
a
注,稳定性介于列主元法和全主元法之间。
§ 1 Gaussian Elimination – Pivoting Strategies
§ 1 Gaussian Elimination – Pivoting Strategies
实际应用中直接调用 Gauss
Elimination 解 3阶线性方程组的结果:
结合全主元消去后的结果:
高斯 -若当 消去法 /* Gauss-Jordan Method */
与 Gaussian Elimination 的主要区别:
每步不计算 mik,而是先将当前主元 akk(k) 变为 1;
把 akk(k) 所在列的上、下元素全消为 0;
bAxIbxA 1
Hey! Isn’t it better than
Gaussian Elimination?
What makes you
say so?
Obviously we no longer
need the backward
substitution!
You’d better wait till we go
through the next section
to draw your
conclusion…
§ 1 Gaussian Elimination – Gauss-Jordan Method
运算量 /* Amount of Computation */
§ 1 Gaussian Elimination – Amount of Computation
由于计算机中乘除 /* multiplications / divisions */ 运算的时间远远超过加减 /* additions / subtractions */ 运算的时间,故估计某种算法的运算量时,往往只估计 乘除的次数,而且通常以乘除次数 的 最高次幂 为运算量的 数量级 。
Gaussian Elimination:
Step k,设,计算因子且计算
0)(?kkka )...,,1(/ )()( nkiaam kkkkikik
).,,,,1,(
)()()1(
)()()1(
nkji
bmbb
amaa
k
kik
k
i
k
i
k
kjik
k
ij
k
ij
共进行 n? 1步
)(
)2(
2
)1(
1
2
1
)(
)2(
2
)2(
22
)1(
1
)1(
12
)1(
11
............
...
...
n
nn
n
nn
n
n
b
b
b
x
x
x
a
aa
aaa
)()( / nnnnnn abx?
)1.,,,,1()( 1
)()(
ni
a
xab
x i
ii
n
ij
j
i
ij
i
i
i
(n? k) 次(n )2 次次
(n? k) (n? k + 2) 次
nnn
knkn
n
k
6
5
23
)2)((
23
1
1
消元乘除次数:
1 次(n? i +1) 次
22)1(1
21
1
nninn
i
回代乘除次数:
Gaussian Elimination 的总乘除次数为
,运算量为 级。nnn
3
1
3 2
3
3
3n
§ 1 Gaussian Elimination – Amount of Computation
Complete Pivoting:
比 Gaussian Elimination多出 比较,保证稳定,但费时。
3
3nO
Partial Pivoting:
比 Gaussian Elimination只多出 比较,略省时,但不保证稳定。
3
2nO
Scaled Partial Pivoting:
比 Gaussian Elimination多出 除法和 比较,比列主元法稳定。但若逐次计算 si(k),则比全主元法还慢。
)( 2nO 3
2nO
Gauss-Jordan Method:
运算量约为 。故通常只用于求逆矩阵,而不用于解方程组。求逆矩阵即 。
23nO
1|| AIIA HW,p.42 #1
p.43 #6
§ 1 Gaussian Elimination – Amount of Computation
Lab 05,Matrix Inversion
Use Gauss-Jordan Method with Partial Pivoting to find the
inverse of a given matrix.
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?n real numbers
which are the entries of the matrix,The numbers are separated by spaces
or new lines.
Output
Each entry of the inverse matrix is to be printed as in the C fprintf:
fprintf(outfile,"%16.8e?",a ); /* here?represents a space*/
nnnnn aaaaaaa,,,.,,.,,.,,122111211
§ 1 Gaussian Elimination – Amount of Computation
The entries of an n?n matrix are to be printed in the following format:
If the matrix is not invertable,print the message,The?matrix?is?
singular.\n”,
The outputs of two test cases must be seperated by a blank line.
Sample Input
3
–38?5
2?–7?4
1?9?–6
2
0?1
0?2
–1
nnnn
n
n
aaa
aaa
aaa
...
..
.
..
.
..
.
..
.
...
...
21
22221
11211
Sample Output (?represents a space)
2.55319149e-0023.95744681e-0012.85106383e-001?
6.80851064e-0025.53191489e-0029.36170213e-002?
1.06382979e-0011.48936170e-0012.12765957e-002?
The?matrix?is?singular.
/* Direct Method for Solving Linear Systems */
求解 bxA
§ 1 高斯消元法 /* Gaussian Elimination */
高斯消元法:
思路首先将 A化为上三角阵 /* upper-triangular matrix */,
再回代求解 /* backward substitution */。
=
§ 1 Gaussian Elimination – The Method
消元 记
,)( )1()1( nnijaAA
)1(
)1(
1
)1(
...
nb
b
bb
Step 1,设,计算因子0)1(
11?a )...,,2(/ )1(11)1( 11 niaam ii
将增广矩阵 /* augmented matrix */ 第 i 行? mi1? 第 1行,
得到
)1(1)1(1)1(12)1(11,,,baaa n
)2(A )2(b? ).,,,,2,(
)1(
11
)1()2(
)1(
11
)1()2(
nji
bmbb
amaa
iii
jiijij
其中
Step k,设,计算因子且计算
0)(?kkka )...,,1(/ )()( nkiaam kkkkikik
).,,,,1,(
)()()1(
)()()1(
nkji
bmbb
amaa
k
kik
k
i
k
i
k
kjik
k
ij
k
ij
共进行? 步n? 1
)(
)2(
2
)1(
1
2
1
)(
)2(
2
)2(
22
)1(
1
)1(
12
)1(
11
..
.
..
.
..
....
...
...
n
nn
n
nn
n
n
b
b
b
x
x
x
a
aa
aaa
回代
)()( / nnnnnn abx?
0)(?nnnaWhat if?No unique solution exists.
)1.,,,,1()( 1
)()(
ni
a
xab
x i
ii
n
ij
j
i
ij
i
i
i
0)(?iiiahat if?
Then we must find the
smallest integer k? i with
,and interchange
the k-th row with the i-th
row,
0)(?ikia
What if we can’t
find such k?No unique
solution exists.
定理 若 A的所有 顺序主子式 /* determinant of leading
principal submatrices */ 均不为 0,则高斯消元无需换行即可进行到底,得到唯一解。
iii
i
i
aa
aa
A
.,,
.,,.,,.,,
.,,
)d e t (
1
111
注,事实上,只要 A 非奇异,即 A?1 存在,则可通过逐次消元及行交换,将方程组化为三角形方程组,求出唯一解。
§ 1 Gaussian Elimination – The Method
选主元消去法 /* Pivoting Strategies */
例,单精度解方程组
2
110
21
21
9
xx
xx
/* 精确解为 和 */...1000...00.1
101
1
91
x
8个
...8999...99.02 12 xx
8个用 Gaussian Elimination计算:
9112121 10/ aam
999
2122 10101010.,,0.011ma
8个
9212 1012mb
99
9
10100
1110
0,1 12 xx
小主元 /* Small pivot
element */ 可能导致计算失败。
§ 1 Gaussian Elimination – Pivoting Strategies
全主元消去法 /* Complete Pivoting */
每一步选绝对值最大的元素为主元素,保证 。1||?ikm
Step k,① 选取 ;0||m a x||
, ijnjikji aa kk
② If ik? k then 交换第 k 行与第 ik 行 ;
If jk? k then 交换第 k 列与第 jk 列 ;
③ 消元注,列交换改变了 xi 的顺序,须记录 交换次序,解完后再换回来。
列主元消去法 /* Partial Pivoting,or maximal column pivoting */
省去换列的步骤,每次仅选一列中最大的元。
0||m a x||, iknikki aa k
§ 1 Gaussian Elimination – Pivoting Strategies
例:
211
1110 9
1,1 12 xx
110
211
1110
211
9
注,列主元法没有全主元法稳定。
0,1 12 xx
例:
211 10101
99
99
99
10100
10101
注意:这两个方程组在数学上 严格等价 。
标度化列主元消去法 /* Scaled Partial Pivoting */
对每一行计算 。为省时间,si 只在初始时计算一次。以后每一步考虑子列 中 最大的 aik 为主元。
||m a x1 ijnji as
nk
kk
a
a
..,i
ik
s
a
注,稳定性介于列主元法和全主元法之间。
§ 1 Gaussian Elimination – Pivoting Strategies
§ 1 Gaussian Elimination – Pivoting Strategies
实际应用中直接调用 Gauss
Elimination 解 3阶线性方程组的结果:
结合全主元消去后的结果:
高斯 -若当 消去法 /* Gauss-Jordan Method */
与 Gaussian Elimination 的主要区别:
每步不计算 mik,而是先将当前主元 akk(k) 变为 1;
把 akk(k) 所在列的上、下元素全消为 0;
bAxIbxA 1
Hey! Isn’t it better than
Gaussian Elimination?
What makes you
say so?
Obviously we no longer
need the backward
substitution!
You’d better wait till we go
through the next section
to draw your
conclusion…
§ 1 Gaussian Elimination – Gauss-Jordan Method
运算量 /* Amount of Computation */
§ 1 Gaussian Elimination – Amount of Computation
由于计算机中乘除 /* multiplications / divisions */ 运算的时间远远超过加减 /* additions / subtractions */ 运算的时间,故估计某种算法的运算量时,往往只估计 乘除的次数,而且通常以乘除次数 的 最高次幂 为运算量的 数量级 。
Gaussian Elimination:
Step k,设,计算因子且计算
0)(?kkka )...,,1(/ )()( nkiaam kkkkikik
).,,,,1,(
)()()1(
)()()1(
nkji
bmbb
amaa
k
kik
k
i
k
i
k
kjik
k
ij
k
ij
共进行 n? 1步
)(
)2(
2
)1(
1
2
1
)(
)2(
2
)2(
22
)1(
1
)1(
12
)1(
11
............
...
...
n
nn
n
nn
n
n
b
b
b
x
x
x
a
aa
aaa
)()( / nnnnnn abx?
)1.,,,,1()( 1
)()(
ni
a
xab
x i
ii
n
ij
j
i
ij
i
i
i
(n? k) 次(n )2 次次
(n? k) (n? k + 2) 次
nnn
knkn
n
k
6
5
23
)2)((
23
1
1
消元乘除次数:
1 次(n? i +1) 次
22)1(1
21
1
nninn
i
回代乘除次数:
Gaussian Elimination 的总乘除次数为
,运算量为 级。nnn
3
1
3 2
3
3
3n
§ 1 Gaussian Elimination – Amount of Computation
Complete Pivoting:
比 Gaussian Elimination多出 比较,保证稳定,但费时。
3
3nO
Partial Pivoting:
比 Gaussian Elimination只多出 比较,略省时,但不保证稳定。
3
2nO
Scaled Partial Pivoting:
比 Gaussian Elimination多出 除法和 比较,比列主元法稳定。但若逐次计算 si(k),则比全主元法还慢。
)( 2nO 3
2nO
Gauss-Jordan Method:
运算量约为 。故通常只用于求逆矩阵,而不用于解方程组。求逆矩阵即 。
23nO
1|| AIIA HW,p.42 #1
p.43 #6
§ 1 Gaussian Elimination – Amount of Computation
Lab 05,Matrix Inversion
Use Gauss-Jordan Method with Partial Pivoting to find the
inverse of a given matrix.
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?n real numbers
which are the entries of the matrix,The numbers are separated by spaces
or new lines.
Output
Each entry of the inverse matrix is to be printed as in the C fprintf:
fprintf(outfile,"%16.8e?",a ); /* here?represents a space*/
nnnnn aaaaaaa,,,.,,.,,.,,122111211
§ 1 Gaussian Elimination – Amount of Computation
The entries of an n?n matrix are to be printed in the following format:
If the matrix is not invertable,print the message,The?matrix?is?
singular.\n”,
The outputs of two test cases must be seperated by a blank line.
Sample Input
3
–38?5
2?–7?4
1?9?–6
2
0?1
0?2
–1
nnnn
n
n
aaa
aaa
aaa
...
..
.
..
.
..
.
..
.
...
...
21
22221
11211
Sample Output (?represents a space)
2.55319149e-0023.95744681e-0012.85106383e-001?
6.80851064e-0025.53191489e-0029.36170213e-002?
1.06382979e-0011.48936170e-0012.12765957e-002?
The?matrix?is?singular.