(给学生)第五章线性方程组的迭代解法

⎡10 −1 −2⎤ ⎡ x1 ⎤ ⎡7.2⎤
⎢⎢−1
10
−2⎥⎥
⎢ ⎢
x2
⎥ ⎥
=
⎢⎢8.3⎥⎥
⎢⎣−1 −1 5 ⎥⎦ ⎢⎣x3 ⎥⎦ ⎢⎣4.2⎥⎦
试用 Jacobi 迭代法求解，使得 x(k+1) − x(k) < 10−3 。 ∞
解：程序见下页，结果： k= 7 x=
1.099811159 1.199811158 1.29977665
0⎥
gJ = D−1b = ⎢⎢1.5⎥⎥ ⎢⎣ 2 ⎥⎦
⎢⎣ 5 5 ⎦⎥
计算方法第五章
9
原方程改写为
⎡ x1 ⎤ ⎡ 0 0.2 0.1⎤ ⎡ x1 ⎤ ⎡0.3⎤
⎢ ⎢
x2
⎥ ⎥
=
Hale Waihona Puke ⎢⎢0.200.1⎥⎥
⎢ ⎢
x2
⎥ ⎥
+
⎢⎢1.5
⎥ ⎥
⎢⎣ x3 ⎥⎦ ⎢⎣0.2 0.4 0 ⎥⎦ ⎢⎣ x3 ⎥⎦ ⎢⎣ 2 ⎥⎦
第五章 方程组的迭代解法
5.1 引言 5.2 线性方程组的迭代解法 5.3 非线性方程组的迭代解法
计算方法第五章
1
5.1 引言
• 直接法是通过有限步运算后得到线性方程组的 解，解线性方程组还有另一种解法，称为迭代 法，它的基本思想是将线性方程组 Ax=b 化 为 x = Bx + g
• 再由此构造向量序列{x (k)}:
end if norm(x-x0,inf)<ep, break; end； x0=x;k=k+1; end if k==N,warning('已达到迭代次数上限');end disp(['k= ',num2str(k)]), x
x2(
k
+1)
=
−1 a22
(a21 x1(k+1)
+
a23 x3(k )
+"+
a2n xn(k )
−
b2 )
⎪ ⎪
#
⎪ ⎪
xn(
k
+1)
⎩
=
−1 ann
(an1 x1(k+1)
+"+
a x (k+1) n,n−1 n−1
−
bn )
计算方法第五章
14
就是
∑ ∑ x(k+1) i
=
1 aii
(bi
−
i −1
for i=1:n
x(i)=(b(i)-A(i,[1:i-1,i+1:n])*x0([1:i-1,i+1:n]))/A(i,i);
end
if norm(x-x0,inf)<ep,break;end
x0=x;k=k+1;
end if k==N,Warning('已达到迭代次数上限');end disp(['k= ',num2str(k)]),x
end
if k==N,Warning('已达到迭代次数上限');end
disp(['k= ',num2str(k)]),x
计算方法第五章
11
例2 试用Jacobi迭代法解线性方程组
⎡ 0.76
⎢ ⎢
−0.0
⎢−0.14
⎢⎣−0.16
−0.01 0.88 −0.03 0.06
−0.14 −0.03 1.01 −0.12
+
0.4 x2(k+1)
+
2
x 取初值 (0) = (0, 0, 0)T，迭代到第6次得
x(6) = (1.000, 2.000, 3.000)T
显然，Seidel迭代法收敛速度快。
计算方法第五章
16
用 Gauss-Seidel 迭代法解线性方程组程序
%maseidel.m function x=maseidel(A,b,x0,ep,N) %用途：用 Gauss-Seidel 迭代法解线性方程组 Ax=b %格式：x=maseidel(A,b,x0,ep,N) A 为系数矩阵，b 为右端向量， % x0 为初始向量（默认零向量），ep 为精度（默认 1e-6），N 为迭代次数上限 %（默认 500 次），x 返回近似解向量 n=length(b);if nargin<5,N=500;end；if nargin<4,ep=1e-6;end if nargin<3,x0=zeros(n,1);end；x=zeros(n,1); k=0; while k<N
计算方法第五章
⎧ ⎪⎪ ⎨ ⎪
x (k+1) 1
x (k+1) 2
= = #
(b1 (b2
− (a12 x2(k ) − (a21 x1(k )
+ " + a1n xn(k) )) + a23 x3(k ) + " +
/ a11 a1n xn(k ) ))
/
a22
8 ⎪⎩ xn(k+1)
=
(bn
−
(an1 x1(k )
BJ
=
I
−
D−1 A =
⎢ ⎢
⎢⎣
1
⎥ ⎥
−
⎢ ⎢
1⎥⎦ ⎢
⎢
1 10
⎥ ⎥
⎢⎢−2
1
⎥ ⎥
⎢⎣ −1
10 −2
−1⎥⎥ 5 ⎥⎦
⎢⎣
5 ⎥⎦
⎡ ⎢
0
⎢
2 1⎤
10
10
⎥ ⎥
⎡0
0.2 0.1⎤
⎡0.3⎤
=
⎢2 ⎢10
⎢ ⎢
1
0 2
1⎥ 10 ⎥
⎥
=
⎢⎢0.2 ⎢⎣0.2
0 0.4
0.1⎥⎥ , 0 ⎥⎦
−0.16⎤ ⎡ x1 ⎤ ⎡0.68⎤
0.06
⎥ ⎥
⎢ ⎢
x2
⎥ ⎥
=
⎢⎢1.18⎥⎥
−0.12⎥
0.72
⎥ ⎦
⎢ ⎢ ⎣
x3 x4
⎥ ⎥ ⎦
⎢0.12⎥ ⎢⎣0.74⎥⎦
解： 程序见下页，结果为：
K=13 x=
1.2715 1.2844 0.4858 1.2843
计算方法第五章
12
2
例 3：设方程组为
clear;format long;clc;
%用途：用 Jacobi 迭代法解线性方程组 Ax=b,A 为系数矩阵，b 为右端向量，
% x0 为初始向量，ep 为精度，N 为最大次数，x 是近似解向量
A=[10 -2 -1;-2 10 -1;-1 -2 5];b=[3 15 10]';N=500;ep=1e-6;
%用途：用 Jacobi 迭代法解线性方程组 Ax=b
%格式：x=majacobi(A,b,x0,ep,N) A 为系数矩阵，b 为右端向量，
% x0 为初始向量（默认零向量），ep 为精度（默认 1e-6），N 为
%（默认 500 次），x 返回近似解向量
x=zeros(n,1); k=0;
while k<N
计算方法第五章
3
⎧ ⎪ ⎪
x1
=
−1 a11
(a12 x2
+"+
a1n xn
−
b1 )
⎪
⎪ ⎨
x2
=
−1 a22
(a21 x1
+
a23 x3
+"+
a1n xn
−
b2 )
⎪ ⎪
#
⎪ ⎪ xn ⎩
=
−1 ann
(an1 x1
+"+
a x n,n−1 n−1
− bn )
于是有 Jacobi 迭代公式
计算方法第五章
x(11) = (1.0000, 2.0000, 3.0000)T
计算方法第五章
10
例 3、试用 Jacobi 迭代法解线性方程组
⎡10 −2 −1⎤ ⎡ x1 ⎤ ⎡ 3 ⎤
⎢⎢−2
10
−1⎥⎥
⎢ ⎢
x2
⎥ ⎥
=
⎢⎢15⎥⎥
⎢⎣−1 −2 5 ⎥⎦ ⎢⎣ x3 ⎥⎦ ⎢⎣10⎥⎦
-----------------------------------------------------------------------
⎟
⎜
⎟
⎜
0
⎟
⎜⎝ −an1
" −an,n−1 0 ⎟⎠
⎛0 ⎜ ⎜
−a12 0
" %
−a1n #
⎞ ⎟ ⎟
U =⎜
%
⎟
⎜ ⎜
0
−an
−1,
n
⎟ ⎟
⎜⎝ 0
0 ⎟⎠
计算方法第五章
6
1
方程组可写为
(D− L−U)x = b
即
即
Dx = (L + U )x + b
x = D−1(L + U ) x + D−1b
a x a x a x b ⎧
+
+"+
=
⎪ 11 1
12 2
1n n
1
a x a x a x b ⎪
⎨
21
⎪
+
1
+"+