数值代数上机习题:线性方程组的古典迭代解法

合集下载
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

x ( k ) BGS x ( k 1) g GS ,
其分量形式为:
k 1,2,3,...
(k ) 1 n ( k 1) b1 ), ( a1 j x j x1 a j 2 11 n ( k ) 1 i 1 (k ) ( k 1) bi ), x ( a x aij x j i ij j a j j i 1 1 ii (k ) 1 n 1 ( k 1) bn ), xn ( anj x j a j 1 nn
4.SOR 迭代中最佳松弛因子的选取 由于 SOR 迭代法的效果和其松弛因子 w 的选取有关, 所以有必要选取合适的松弛因子。 《数 值线性代数》 (徐树方,北京大学出版社)中的相关结论,当选择最佳松弛因子
wopt
2 1 1 ( B) 2
时,SOR 方法的迭代速度最快。
算法使用 matlab 实现, 程序中采用无穷范数 || . || , 每种迭代法最大允许的迭代次数为 100000 次。
b(1) ah 2 y (0) * b(n 1) ah 2 y (1) * ( h)
3.迭代终止条件 首先确定要求的精度 tol ,我们希望当
|| xk x* || tol
则停止迭代。 由《数值线性代数》 (徐树方,北京大学出版社) P105 定理 4.2.3,我们知道,对于迭代格 式 xk Bxk 1 g ,若 || B || q 1 且 || I || 1 ,则迭代序列 {xk } 的第 k 次近似解和精确解之 间有估计式 || xk x* ||
1 , n

yi 1 2 yi yi 1 yi 1 yi a h2 h
简化为:
( h) yi 1 (2 h) yi yi 1 ah 2
从而离散后得到的线性方程组的系数矩阵为:
h (2 h) (2 h) h A (2 h) h (2 h)
% ---------------------------- % % filename : jacobi.m % abstract : Jacobi迭代法解线性方程组 Ax = b % 调用方式: % [x, k, e] = jacobi(A,b,y) % 返回值: % x : 迭代解 % k : 迭代次数 % e : 与精确解的误差 % ---------------------------- % % Reference : 数值线性代数,徐树方,北京大学出版社 % ---------------------------- % function [x, k, e] = jacobi(A,b,y) N = 100000; tol = 1e-4; q D L U B g = = = = = = 1 - 1/(100*100); diag(diag(A)); triu(A)-A; tril(A)-A; D\(L+U); D\b; % y为精确解
然后在矩阵乘法时对下标处理一下即可。但是考虑到三种迭代方法的一般性,且本题中
n 100 并不是很大,所以实验中并没有采用紧缩存储,而是采用了直接存储。
2.边值条件的处理 由于差分得到的方程组的第一行和最后一行中分别出现了边值 y(0)与 y(1)作为常数项, 因此 要在常向量的第一项和最后一项作一些修改:
283 199 111 108
五.结果分析
1.迭代法的收敛速度比较
从表格数据可以看出,当 不是很小时,GS 迭代的速度要快于 Jacobi 迭代的速度,而 SOR 迭代(选取最佳松弛因子)的速度要远远快于前两者。当 减小时,三种迭代方法的收敛速 度均增快,而当 0.0001 时,求得 SOR 迭代的最佳松弛因子接近 1,基本退化为 GS 迭 代,这是因为当 很小时,方程系数矩阵的谱半径 ( A) 也很小,所以求出的最佳松弛因子
四.计算结果
题目要求 a
1 , n 100, 下面是 1,0.1,0.01,0.0001 时三种方法计算所得近似解需要的迭 2
代次数以及与精确解的误差。 1.Jacobi 迭代

迭代次数 k
与精确解的误差 || xk x* || 3.100e-4 8.826e-3 6.606e-2 4.950e-3
它的分量形式为:
xi
(k )
(1 ) xi
( k 1)

n 1 i 1 (k ) (k ) ( aij x j aij x j bi ), aii j 1 j i 1
其中 称为松弛因子,当 1 时称为超松弛;当 1 时叫低松弛; 1 时就是 G-S 迭 代。
wopt
2 1 1 ( B ) 2 接近于 1。
2.误差比较 从表格可以看出, 三种迭代方法所得近似解的误差基本一样。 当 1 时, 误差在 1e-4 左右, 符合定理 4.2.3 的结果。但当 减小时,误差却出现先增大再减小的情况,这可能是由于前
面确定迭代终止时的近似 q 1 h 不够准确所致。
3.超松弛迭代法(SOR 迭代法) :
i 2,..., n 1
SOR迭代法是对G-S迭代法的加权平均改造,即 x
(k )
(1 ) x ( k 1) x ,
(k )
x
(k )
为G-S迭代解,即
x ( k ) (1 ) x ( k 1) ( D 1 Lx ( k ) D 1Ux ( k 1) D 1b),
2
参考文献 [1] 徐树方, 《数值线性代数》 ,北京大学出版社
附录:matlab 代码
% ---------------------------- % % filename : test.m % abstract : 用于测试的m脚本文件,输入系数eps, 显示三种迭代方法得到的近似解,迭代次数及与精确解的误差 % % ---------------------------- % % ---------------------------- % a = 0.5; n = 100; eps = input('请输入系数eps:'); [A b] = genmatrix(eps, a, n); y = exactsolution(eps, a, n); % Jacobi迭代 [x k e] = jacobi(A,b,y); disp('% Jacobi迭代 '); disp('近似解:'); disp(x); fprintf('误差:%.6f;迭代次数:%d\n\n',e,k); % GS迭代 % 输入eps % 生成系数矩阵和常向量 % 获得精确解
a13 a1n a23 a2 n 0 , an 1,n 0
则原方程组可改写为:
x BJ x g J ,
其中
(2.2)
BJ D 1 ( L U ), g J D 1b,
给定初始向量:
x0 ( x1 , x2 ,..., xn ) T ,
数值代数上机习题一:线性方程组的古典迭代解法
一.题目要求
考虑两点边值问题
d 2 y dy a ,0 a 1 2 dx dx y (0) 0, y (1) 1.
容易知道它的精确解为
y
1 a (1 e ) ax 1 / 1 e
x
为了把微分方程离散,把[0, 1]区间 n 等分,令 h 得到差分方程为:
q || xk 1 xk || 。 1 q
4
由 题 目 要 求 知 我 们 需 要 有 || xk x* || 10
,而由上面的迭代估计,只要
q 1 q 即可。而本题中 q 可近似取为 || xk 1 xk || 10 4 ,即 || xk 1 xk || 10 4 1 q q q 1 h 2 ,因此最后令迭代终止条件为 || xk 1 xk || 10 4 1 q q

1 0 .1 0.01 0.0001
22274 8754 662 122
2.GS 迭代

迭代次数 k
与精确解的误差 || xk x* || 2.900e-4 8.821e-3 6.606e-2 4.950e-3
1 0 .1 0.01 0.0001
( 0)
(0)
(0)
三.算法实现
1.方程的表示及存储
由于本题中线性方程组的系数矩阵为三对角矩阵,所以可以采用紧缩方法存储,即
0 A
(2 h) h (2 h) h (2 h) h (2 h) 0
D diag (a11 , a22 ,..., ann ),
0 a 21 L a31 an1 , 0
(2.1)
0 a32 0
an 2 an ,n 1
0 a12 0 U
要求: (1)对 1, a
1 , n 100, 分别用 Jacobi, G-S 和 SOR 迭代方法求线性方程组的解,要求 2
有 4 位有效数字,然后比较与精确解的误差。 (2)对 0.1, 0.01, 0.0001, 考虑同样的问题。
二.算法原理
1.Jacobi 迭代法(J-迭代法) : 对于非奇异线性方程组 Ax b ,令 A D 1
2.Gauss-Seidel 迭代法(G-S 迭代法) : 类似于 J-迭代法,给定初值:
x0 ( x1 , x2 ,..., xn ) T ,
( 0)
(0)
(0)

BGS ( D L) 1U , g GS ( D L) 1 b,
则得到 G-S 迭代公式:
从上面的介绍可以看出, 三种迭代法的本质区别在于迭代公式的不同, 因此算法的基本框架 如下:
1)输入初值 x0 ( x1 , x2 ,..., xn ) T , 精度要求 tol,最大允许迭代次数 N,令 k := 0,转 2); 2)由迭代公式计算出 x ( k 1) , 转 3); 3)若 || x ( k ) x ( k 1) || tol ,则停止,输出 x ( k 1) 作为方程的解;若否,转 4) 4)若 k = N,则停止,输出迭代失败信息;若否,令 k := k+1,转 2)
[x k e] = gs(A,b,y); disp('% GS迭代 '); disp('近似解:'); disp(x); fprintf('误差:%.6f;迭代次数:%d\n\n',e,k); % SOR迭代 [x k e] = sor(A,b,y); disp('% SOR迭代 '); disp('近似解:'); disp(x); fprintf('误差:%.6f;迭代次数:%d\n\n',e,k);
10788 3855 344 109
3.SOR 迭代

迭代次数 k
误差 || xk x* || 3.000e-4 8.829e-3 6.606e-2 4.950e-3
松 弛 因 子w 1.9384 1.8921 1.4985 1.0021
1 0 .1 0.01 0.0001
由(2.2)可以构造迭代公式:
( 0)
(0)
(0)
xk BJ xk 1 g ,
其分量形式为:
k 1,2,3,...
(k ) 1 n ( k 1) b1 ), ( a1 j x j x1 a11 j 2 n ( k ) 1 i 1 ( k 1) ( k 1) aij x j bi ), xi ( aij x j aii j 1 j i 1 (k ) 1 n 1 ( k 1) bn ), xn ( anj x j ann j 1
相关文档
最新文档