基于MATLAB的迭代求解线性方程组（附完整代码与算法）前言一. 三个变换二. Sylvester方程三. Jacobi迭代法四. Gauss-Seidel迭代法五. SOR迭代法六. 两步迭代法

115 阅读 0 评论 76 点赞

我是靠谱客的博主平淡高跟鞋，最近开发中收集的这篇文章主要介绍基于MATLAB的迭代求解线性方程组（附完整代码与算法）前言一. 三个变换二. Sylvester方程三. Jacobi迭代法四. Gauss-Seidel迭代法五. SOR迭代法六. 两步迭代法，觉得挺不错的，现在分享给大家，希望可以做个参考。

概述

前言

在之前的文章中更新了线性方程组的基本解法，大型方程组的分解求法。本节将介绍线性方程组的迭代求解。

一. 三个变换

在线性方程组的迭代求解中，会用到系数矩阵A的上三角矩阵、对角矩阵和下三角矩阵。这三种变化在MATLAB中，可以直接利用函数实现。

1.1 上三角变换

triu(A,1)

1.2 对角变换

diag(A)

1.3 下三角变换

tril(A,-1)

针对此类矩阵的解释，可以参看这篇文章：

MATLAB:矩阵基础_唠嗑！的博客-CSDN博客

例题1

对以下矩阵A做此三种变换。

$a=begin{bmatrix}1&2&-2\1&1&1\2&2&1 end{}$

解：

MATLAB代码如下：

clc;clear;
A=[1 2 -3;1 1 1;2 2 1];
%上三角变换
up=triu(A,1)
%对角变换
diag=diag(A)
%下三角变换
down=tril(A,-1)

运行结果：

up =

0 2 -3
0 0 1
0 0 0

diag =

1
1
1

down =

0 0 0
1 0 0
2 2 0
分析：

二. Sylvester方程

Sylvester方程的形式如下：

$ax+xb=c$

此方程中，A是 $mtimes m$ 矩阵，B是 $ntimes n$ 矩阵，C和X都是 $mtimes n$ 矩阵。

在MATLAB中可以直接调用函数求解，如下：

X=sylvester(A,B,C)

例题2

求解Sylvester方程AX+XB=C的解X。其中A,B,C矩阵如下：

$a=begin{bmatrix}1&0&2&3\ 4&1&0&2\ 0&5&5&6\ 1&7&9&0 end{}$ ， $b=begin{bmatrix}0&-1\1&0 end{}$ ， $c=begin{bmatrix}1&0\2&0\0&3\1&1 end{}$

解：

MATLAB代码如下：

clc;clear;
A=[1 0 2 3;4 1 0 2;0 5 5 6;1 7 9 0];
B=[0 -1;1 0];
C=[1 0;2 0;0 3;1 1];
X=sylvester(A,B,C)

运行结果如下：

X =

0.4732 -0.3664
-0.4006 0.3531
0.3305 -0.1142
0.0774 0.3560

三. Jacobi迭代法

原方程组Ax=b中的矩阵A可以写成，如下：

$a=d-l-u$

上式子中-L和-U分别为A矩阵的严格下，上三角部分，当然不包含对角元素。其中矩阵D的表达如下：

$d=diag[a_{11},a_{22},cdots,a_{nn}]$

根据迭代的思想，x可得如下表达式：

$x^{(k+1)}=bx^{(k)}+f$

上式子中的B，f理解如下：

$b=d^{-1}(l+u)=i-d^{-1}a,;f=d^{-1}b$

为了更好理解这种迭代形式，将原线性方程组写全，如下：

如果系数矩阵非奇异，即 $a_{ii}neq0quad (i=1,2,cdots,n)$ ，那么可得：

运用迭代公式 $x^{(k+1)}=bx^{(k)}+gquad (k=0,1,2,cdots)$ ，最终的方程组如下：

在调用jacobi()函数之前，需要提前编码，我们来看下jacobi函数的内部结构，如下：

function y=jacobi(a,b,x0)
D=diag(diag(a));
U=-triu(a,1);
L=-tril(a,-1);
B=D(L+U);
f=Db;
y=B*x0+f;
n=1;
while(norm(y-x0))>=1.0e-6
x0=y;
y=B*x0+f;
n=n+1;
end
n
end

例题3

设x0=0，精度为 $10^{-6}$ ，利用Jacobi方法求解以下方程组。

$begin{cases}10x_1-x_2=9\-x_1+10x_2-2x_3=7\-2x_2+10x_3=6 end{}$

解：

MATLAB代码如下：


clc;clear;
a=[10 -1 0;-1 10 -2;0 -2 10];
b=[9;7;6];
jacobi(a,b,[0;0;0])
%注意初始值也为一个矩阵
%对jacobi函数的定义
function y=jacobi(a,b,x0)
D=diag(diag(a));
U=-triu(a,1);
L=-tril(a,-1);
B=D(L+U);
f=Db;
y=B*x0+f;
n=1;
while(norm(y-x0))>=1.0e-6
x0=y;
y=B*x0+f;
n=n+1;
end
n
end

运行结果：

n =

ans =

0.9958
0.9579
0.7916

四. Gauss-Seidel迭代法

借助Jacobi的方法，迭代方程也可以如下形式：

$x^{(k+1)}=gx^{(k)}+f$

同样地，-L、-U分别为A的严格下、上三角部分，且不包含对角线元素。G，f,D的定义如下：

$g=(d-l)^{-1}u,;f=(d-l)^{-1}b,;d=diag[a_{11},a_{22},cdots,a_{nn}]$

为了深入理解此种迭代方式，可列出矩阵中的元素。

在迭代计算过程中，需要同时保留两个近似解向量 $x^{(k)}$ 和 $x^{(k+1)}$ 。由此可以把迭代公式改写为如下：

备注：由于显示的原因，上式子中的“L"代表省略号

由此可以在MATLAB中编码seidel()函数，如下：

function y=seidel(a,b,x0)
D=diag(diag(a));
U=-triu(a,1);
L=-tril(a,-1);
G=(D-L)U;
f=(D-L)b;
y=G*x0+f;
n=1;
while norm(y-x0)>=1.0e-6
x0=y;
y=G*x0+f;
n=n+1;
end
n
end

例题4

设x0=0，精度为 $10^{-6}$ ，利用Gauss-Seidel方法求解以下方程组。

$begin{cases}10x_1-x_2=9\-x_1+10x_2-2x_3=7\-2x_2+10x_3=6 end{}$

解：

MATLAB代码如下：

clc;clear;
a=[10 -1 0;-1 10 -2;0 -2 10];
b=[9;7;6];
seidel(a,b,[0;0;0])
function y=seidel(a,b,x0)
D=diag(diag(a));
U=-triu(a,1);
L=-tril(a,-1);
G=(D-L)U;
f=(D-L)b;
y=G*x0+f;
n=1;
while norm(y-x0)>=1.0e-6
x0=y;
y=G*x0+f;
n=n+1;
end
n
end

运行结果：

n =

ans =

0.9958
0.9579
0.7916

分析：此例子可以直接求解。也可能会出现使用Seidel迭代时，结果不收敛的情况。

五. SOR迭代法

在某些情况下，Jacobi法和Gauss-Seidel法收敛的过程较慢，为了进一步改进，可以引入一种新的迭代法：逐次超松弛迭代法（Succesise Over Relaxation)，简记为SQR法。

迭代的公式如下：

$x^{(k+1)}=(d-wl)^{-1}((1-w)d+wu)x^{(k)}+w(d-wl)^{-1}b$

上式子中，w的最佳值一般在[1,2)之间，具体视情况而定。

松弛法的主题思想是将 $delta x$ 乘上一个参数因子w来作为修正项，从而得到新的近似解。迭代公式可见如下：

$x^{(k+1)}=x^{(k)}+wdelta x$

上式子中的 $delta x$ 的计算，可见如下：

由此可得迭代的化简过程，如下：

按照上述式子计算Ax=b的近似解序列的方法，就被称之为松弛法，其中w为松弛因子。当w<1时称为低松弛；当w>1时称为超松弛；当w=1时，其实就是Gauss-Seidel迭代。

根据以上过程构造sor函数，MATLAB代码如下：

function y=sor(a,b,w,x0)
D=diag(diag(a)); %形成对角矩阵，其余元素为0
U=-triu(a,1);
L=-tril(a,-1);
M=(D-w*L)((1-w)*D+w*U);
f=(D-w*L)b*w;
y=M*x0+f;
n=1;
while norm(y-x0)>=1.0e-6
x0=y;
y=M*x0+f;
n=n+1;
end
n
end

例题5

设x0=0，精度为 $10^{-6}$ ，w=1.103，利用SOR方法求解以下方程组。

$begin{cases}10x_1-x_2=9\-x_1+10x_2-2x_3=7\-2x_2+10x_3=6 end{}$

解：

MATLAB代码如下：

clc;clear;
a=[10 -1 0;-1 10 -2;0 -2 10];
b=[9;7;6];
sor(a,b,1.103,[0;0;0])
%函数部分
function y=sor(a,b,w,x0)
D=diag(diag(a)); %形成对角矩阵，其余元素为0
U=-triu(a,1);
L=-tril(a,-1);
M=(D-w*L)((1-w)*D+w*U);
f=(D-w*L)b*w;
y=M*x0+f;
n=1;
while norm(y-x0)>=1.0e-6
x0=y;
y=M*x0+f;
n=n+1;
end
n
end

运行结果：

n =

ans =

0.9958
0.9579
0.7916

六. 两步迭代法

当线性方程系数矩阵为对称正定时，还可以使用一种特殊的迭代法来解决：两步迭代法。迭代的公式，如下：

首先可得：

$(d-l)x^{(k+frac{1}{2})}=ux^{(k)}+b$

$(d-u)x^{(k+1)}=lx^{(k+1)}=lx^{(k+frac{1}{2})}+b$

进一步推导可得：

$x^{(k+frac{1}{2})}=(d-l)^{-1}ux^{(k)}+(d-l)^{-1}b$

$x^{(k+1)}=(d-u)^{-1}lx^{(k+frac{1}{2})}+(d-u)^{-1}b$

根据以上计算过程，可以编写twostp()函数，MATLAB代码如下：

function y=twostp(a,b,x0)
D=diag(diag(a));
U=-triu(a,1);
L=-tril(a,-1);
G1=(D-L)U;f1=(D-L)b;
G2=(D-U)L;f2=(D-U)b;
y=G1*x0+f1;
y=G2*y+f2;
n=1;
while norm(y-x0)>=1.0e-6
x0=y;
y=G1*x0+f1;
y=G2*y+f2;
n=n+1;
end
n
end

例题6

利用两步法，求解如下方程组：

$begin{pmatrix}10&-1&2&0\ -1&11&-1&3\ 2&-1&10&-1\ 0&3&-1&8 end{}begin{pmatrix}x_1\x_2\x_3\x_4 end{}=begin{pmatrix}6\25\-11\15 end{}$

解：

MATLAB代码如下：

clc;clear;
a=[10 -1 2 0;-1 11 -1 3;2 -1 10 3;0 3 -1 8];
b=[6;25;-11;15];
twostp(a,b,[0;0;0;0])
function y=twostp(a,b,x0)
D=diag(diag(a));
U=-triu(a,1);
L=-tril(a,-1);
G1=(D-L)U;f1=(D-L)b;
G2=(D-U)L;f2=(D-U)b;
y=G1*x0+f1;
y=G2*y+f2;
n=1;
while norm(y-x0)>=1.0e-6
x0=y;
y=G1*x0+f1;
y=G2*y+f2;
n=n+1;
end
n
end

运行结果：

n =

ans =

1.0791
1.9824
-1.4044
0.9560

最后

以上就是平淡高跟鞋为你收集整理的基于MATLAB的迭代求解线性方程组（附完整代码与算法）前言一. 三个变换二. Sylvester方程三. Jacobi迭代法四. Gauss-Seidel迭代法五. SOR迭代法六. 两步迭代法的全部内容，希望文章能够帮你解决基于MATLAB的迭代求解线性方程组（附完整代码与算法）前言一. 三个变换二. Sylvester方程三. Jacobi迭代法四. Gauss-Seidel迭代法五. SOR迭代法六. 两步迭代法所遇到的程序开发问题。

如果觉得靠谱客网站的内容还不错，欢迎将靠谱客网站推荐给程序员好友。

本图文内容来源于网友提供，作为学习参考使用，或来自网络收集整理，版权属于原作者所有。

本文分类：MATLAB
浏览次数：115 次浏览
发布日期：2023-03-11 12:58:02
本文链接：https://www.kaopuke.com/article/k-p-k_14_ujo_6_fy_13_z_10_1.html

基于MATLAB的迭代求解线性方程组（附完整代码与算法）前言一. 三个变换二. Sylvester方程三. Jacobi迭代法四. Gauss-Seidel迭代法五. SOR迭代法六. 两步迭代法

概述

前言

一. 三个变换

1.1 上三角变换

1.2 对角变换

1.3 下三角变换

例题1

二. Sylvester方程

例题2

三. Jacobi迭代法

例题3

四. Gauss-Seidel迭代法

例题4

五. SOR迭代法

例题5

六. 两步迭代法

最后

评论列表共有 0 条评论

发表评论取消回复

基于MATLAB的迭代求解线性方程组（附完整代码与算法）前言一. 三个变换二. Sylvester方程三. Jacobi迭代法四. Gauss-Seidel迭代法五. SOR迭代法六. 两步迭代法

概述

前言

一. 三个变换

1.1 上三角变换

1.2 对角变换

1.3 下三角变换

例题1

二. Sylvester方程

例题2

三. Jacobi迭代法

例题3

四. Gauss-Seidel迭代法

例题4

五. SOR迭代法

例题5

六. 两步迭代法

最后

相关文章

评论列表共有 0 条评论

发表评论 取消回复

发表评论取消回复