matlab trisys,Matlab实现——严格对角占优三对角方程组求解(高斯赛尔德Gauss-Seidel迭代、超松弛) | 学步园...
严格对角占优三对角方程组求解
对中等规模的n阶的(n<100)线性方程组,直接法的准确性和可靠性,所以常采用直接法
对于较高阶的方程组,特别是地于某些偏微分方程离散化后得到的大型稀疏方程组(系统矩
阵绝大多数为零元素),由于直接解法的计算代价较高,使得迭代法更具有竞争力。
于是设计以下的2种
——(1)
其系数矩阵是对角的,且元素满足严格对角占优:
1)追赶法:
利用方程组(1)的特点,应用Gauss消元法求解时,每步只需消一个元素。其消元过程为:
——(2)
得到同解方程组(仍然严格对角占优)为:
——(3)
回代过程(对角占优,不必选主元)为:
——(4)
用追赶法解方程组(2)仅需(5n-4)次乘除过程,(3n-3)次加减过程,算法时间复杂度O(n)。
程序:
function X=trisys(A,D,C,B)
%Input- A is the subdiagonal of the coefficient matrix
% - D is the main diagonal of the coefficient matrix
% - C is the superdiagonal of the coefficient matrix
% - B is the constant vector of the linear system
%Output - X is the solution vector
N=length(B);
X=zeros(N,1);
for k=2:N
mult=A(k-1)/D(k-1);
D(k)=D(k)-mult*C(k-1);
B(k)=B(k)-mult*B(k-1);
end
X(N)=B(N)/D(N);
for k= N-1:-1:1
X(k)=(B(k)-C(k)*X(k+1))/D(k);
end
这个方法的精度很高,和系统的内置函数linsolve(H,B)的求解结果一致
2)迭代法(采用改进的Gauss-Seidel迭代)(这个方法是看了超松弛迭代后,得出的类似方法):
原理介绍:
——(1)
它的Gauss-Seidel迭代方法我们已经很熟悉了:
【1】给一个初始列向量:
【2】利用迭代公式:
经过一定的迭代次数以后,就能得到近似解:
现在对上述的Gauss-Seidel迭代进行加速
得到迭代方法:
(注意:当ω=1时,就是我们所熟悉的Gauss-Seidel迭代)
其中:
ω是迭代加速的相关系数——松弛因子
上述方法可解释为第k+1次迭代近似解的各分量依次为用Gauss-Seidel方法求得的第k+1次迭代近似值和第次近似值的加权平均值。适当选取收敛因子ω(事实上叫做松弛因子),可望该方法比Gauss-Seidel迭代法收敛得更快。
根据以上的原理分析,作出程序如下:
function X=acc(A,D,C,B,P,delta, max1,w)
%Input- A is the subdiagonal of the coefficient matrix
% - D is the main diagonal of the coefficient matrix
% - C is the superdiagonal of the coefficient matrix
% - B is the constant vector of the linear system
% - P is an N x 1 matrix; the initial guess
% - w is the convergence multiplicate
% - delta is the tolerance for P
% - max1 is the maximum number of iterations
% Output - X is an N x 1 matrix: the gauss-seidel approximation
% to the solution of AX = B
N = length(B);
L=P; %L is a mediut
for k=1:max1 %max1th iteration
X=L; %initial the X=[x1;x2;…;xN]=L=[d01;d02;…;d0N]
% the kth iteration of valuing the X
for j=1:N
if j==1
X(1)=(1-w)*X(1)+w*(B(1)-C(1)*X(2))/D(1);
elseif j==N
X(N)=(1-w)*X(N)+w*(B(N)-A(N-1)*X(N-1))/D(N);
else
%X contains the kth approximations
X(j)=(1-w)*X(j)+w*(B(j)-A(j-1)*X(j-1)-C(j)*X(j+1))/D(j);
end
end
err=abs(norm(X-L)); %get the error
L=X;
relerr=err/(norm(X)+eps);
if (err
break
end
end
分析误差:
这时候我们看到w=0.2时误差是e=0.01550147497154
经过类似的试验可以知道:
在w=0.98-w=1之间的时候存在最优松弛因子
我们看到迭代次数的增加,带来误差的显著减小,且迭代次数max1=20的时候精度达到1.0e-11
可见,该方法的求解精度还是令人满意的。
在求解该问题的过程中,对于求解方程组的方法选择是一个很重要的因素,注意到这个系数矩阵是50阶严格对角占优三对角稀疏矩阵,查询了相关知识后,我个人认为,50阶的严格对角占优三对角稀疏矩阵,完全可以用高斯消去法,这是因为高斯消去后的上三角(或者下三角)仍然是严格对角占优,而对于这个稀疏矩阵,迭代法是一个非常不错的选择,而我采取的迭代法受限制的就是这个松弛因子w,注意到0
matlab trisys,Matlab实现——严格对角占优三对角方程组求解(高斯赛尔德Gauss-Seidel迭代、超松弛) | 学步园...相关推荐
- Matlab实现——严格对角占优三对角方程组求解(高斯赛尔德Gauss-Seidel迭代、超松弛)
欢迎前往个人博客 驽马点滴 和视频空间 哔哩哔哩-<挨踢日志> 严格对角占优三对角方程组求解 对中等规模的n阶的(n<100)线性方程组,直接法的准确性和可靠性,所以常采用直接法 对 ...
- 数值计算大作业:Jacobi与Gauss -Seidel迭代求解线性方程组(Matlab实现)
作为研究生的入门课,数值计算的大作业算是所有研究生开学的重要编程作业. 我把Jacobi与Gauss -Seidel迭代求解线性方程组的数值计算作业在MATLAB中编程实现.具体的程序详细标注后放在文 ...
- 高斯赛尔德、牛顿拉尔逊matlab潮流计算
个人电气博文目录传送门 :学好电气全靠它,个人电气博文目录(持续更新中-) 本文讲解复杂电力系统的潮流计算.具体文档程序见下载链接: 文章目录 一.复杂电力系统潮流计算 二.程序 三.下载链接 一.复 ...
- 追赶法求解三对角方程组
1. 来源和背景 对于一个(主)三对角方程组,我们常用"追赶法"来进行求解. 而三对角方程组常常出现于微分方程的数值求解,例如热传导方程的边值问题 {y′′(x)=f(x,y,y′ ...
- 10种基于MATLAB的方程组求解方法
线性方程组的求解包括直接法和迭代法,其中迭代法包括传统的高斯消元法,最速下降法,牛顿法,雅克比迭代法,共轭梯度法,以及智能启发式算法求解法和神经网络学习算法,传统算法可以相互组合改进,智能仿生启发式算 ...
- 随机变量序列 收敛 检验 matlab,利用MATLAB验证样条插值的收敛性
问题提出 多项式插值是不收敛的,即插值的节点多,效果不一定就好.对样条函数插值又如何呢?理论上证明样条插值的收敛性是比较困难的,但通过本实验可以验证这一理论结果.. 实验内容 请按一定的规则分别选择等 ...
- c++引用matlab类,matlab调用C++函数浅谈(一)
由于在下才疏学浅,在网上看各高手指南时亦觉云里雾里,遂决定一切说明从最基础说起,一是方便自己(记性奇差),二是方便似我的小白.以下部分是我从各网站论坛等摘抄.重组.改写过的,以求更加详实明朗,由于参考 ...
- 磁盘驱动读取系统MATLAB仿真,matlab读写..doc
matlab读写. MATLAB二进制数据文件的读写 (2011-06-04 19:44:27) 转载▼ 标签: easleyhux matlab 二进制读写 杂谈分类: MATLAB 所谓二进制格式 ...
- 【 MATLAB 】MATLAB 实现模拟信号采样后的重建(三)一阶保持(FOH)内插
上篇博文采用了零阶保持(ZOH)的方式进行了重构:[ MATLAB ]MATLAB 实现模拟信号采样后的重建(二)零阶保持(ZOH) 这篇博文我们使用一阶保持(FOH)内插来重建信号,采用的案例依然是 ...
最新文章
- Debian8.8解决双系统访问windows磁盘时,有时能成功挂载,有时不能成功挂载的情况...
- 概率分布函数的关系框架
- Distributed Systems笔记-Cryptographic Protocols
- websocket多人聊天php,php-notes/基于websocket实现多人聊天室.md at master · dd-code-site/php-notes · GitHub...
- pku 2387 Til the Cows Come Home
- eclipse 中使用Git
- Redis 安装与简单示例 第一篇
- 怎么维护 SQLite
- 《Android第一行代码(第二版)》源码
- Unity3D客户端实时同步技术
- [译]区块链技术全解析入门版
- 如何使用Git Pull覆盖本地文件
- 平安夜html5源码,平安夜
- 战争迷雾效果 第1章_要探索,不要地图全开!
- oracle如何按时间索引降序,Oracle 降序索引
- VMware Fusion安装cloudera manager
- 什么是session
- 动画中的“表现手法”和“表现方式”分别是什么?今天艺点动画大家庭来给你唠嗑唠嗑
- 中国大公司罹患“大企业病” 国企成亏损重灾区
- 【bzoj3698】【XWW的难题】【有上下界的网络流】
热门文章
- SSD内部逻辑结构——Chip与Package的误解
- utrack调试 艾肯icon_Icon Utrack 设置调试以及常用软件设置[超详细]
- 自适应端在线玩俄罗斯方块源码开源版源码
- 港大全面封杀ChatGPT!AI肆虐高校,作弊丑闻震惊全世界
- 堡垒机、跳板机JumpServer的搭建使用
- 滴答顺风车怎么抢90%以上的订单_顺风车,又来了!!!
- 小程序全时全域自动化箱况检测+信息识别,中集飞瞳全球领先集装箱管理方案,AI自动化箱信息识别+箱况检测+地点报备,智慧港航智能化
- 30个极致实用的谷歌浏览器插件,让你开发事半功倍
- w7设置双显示器_Win7系统如何设置两个显示器即多屏幕模式
- epass 1000ND 使用心得