线性方程组迭代法matlab

合集下载

MATLAB中对线性方程组求解一种方法

MATLAB中对线性方程组求解一种方法

实验报告MATLAB中对线性方程组求解一种方法一.前言MATLAB是一个功能强大的线性方程组求解工具,它特别适合求解大规模的线性方程组。

由于MATLAB是专为矩阵运算而设计的,它的运算效率较高,而且可以使用循环结构,因此用MATLAB可以求解需要千万次运算量的复杂方程组。

下面介绍的求解方法对于小规模方程组也是非常有效的。

二.实验原理迭代法就是用某种极限国采取逐步逼近线性方程组精确解的方法。

方法具有对计算机的存贮单元需求少,程序设计简单、原始系数矩阵在计算过程中不便等优点,是求解大型稀疏矩阵方程组的重要方法。

迭代法不是用有限步运算求精确解,而是通过迭代产生近似逼近精确解。

如Jacobi迭代法、Gauss-Serdel迭代法。

三.实验目的1.熟练MATLAB中运用迭代方法求解线性方程组的原理;2.熟练MATLAB中矩阵的运用和程序设计。

四.实验步骤1.Jacobi迭代法对于线性方程组Ax=b,如果A为非奇异方阵,即aii≠0(i=1,2,…,n),则可将A分解为A=D-L-U,其中D为对角阵,其元素为A的对角元素,L与U为A的下三角阵和上三角阵,于是Ax=b化为:x=D-1(L+U)x+D-1b与之对应的迭代公式为:x(k+1)=D-1(L+U)x(k)+D-1b这就是Jacobi迭代公式。

如果序列{x(k+1)}收敛于x,则x必是方程Ax=b的解。

Jacobi迭代法的MATLAB函数文件Jacobi.m如下:function[x,n]=jacobi(A,b,x0,eps,varargin)format longif nargin==3eps= 1.0e-6; %默认精度M = 200; %参数不足时默认后两个条件elseif nargin<3error('参数不足');returnelseif nargin==5;M=varargin;end[n,m]=size(A);nb=length(b);%当方程组行与列的维数不相等时,停止计算,并输出出错信息if n~=merror('矩阵A行数和列数必须相等!');return;end%当方程组与右端项的维数不匹配时,停止计算,并输出出错信息if n~=nberror('矩阵A的行数必须和b的长度相等!');return;endD =zeros(n,n);for i=1:nif A(i,i)==0error('A对角线元素为零!')return;endD(i,i)=A(i,i); %得到矩阵DendB=inv(D)*(D-A); %B为迭代矩阵g=inv(D)*b; %g为右端项pr=max(abs(eig(B))) %求迭代矩阵谱半径if pr>=1error('迭代矩阵谱半径大于1迭代法不收敛');return;endk=0;tol=1;while tol>=epsx = B*x0+gk = k+1; %迭代步数tol = norm(x-x0);%前后两步迭代结果的误差x0 = x;if(k>=200)disp('Warning: 迭代次数太多,可能不收敛!');return;endend2. 例题用Jacobi迭代法求解下列线性方程组。

MATLAB课程设计可视化(GUI)的线性方程组的Jacobi迭代解法。

MATLAB课程设计可视化(GUI)的线性方程组的Jacobi迭代解法。

华东交通大学课程设计(论文)任务书软件学院学院09 软件工程+电气专业 3 班一、课程设计(论文)题目可视化(GUI)的线性方程组的Jacobi迭代解法二、课程设计(论文)工作自 2011年6月27日起至2011 年 7月1 日止。

三、课程设计(论文) 地点: 电气学院机房四、课程设计(论文)内容要求:1.本课程设计的目的(1)熟练掌握MATLAB语言的基本知识和技能;(2)熟悉MA TLAB下的GUI程序设计;(3)熟悉多项式曲线拟合,MA TLAB的绘图功能;(4)培养分析、解决问题的能力;提高学生的科技论文写作能力。

2.课程设计的任务及要求1)基本要求:(1)利用matlab中的GUI设计窗口设计一个界面程序。

其中主界面包含控制背景颜色与图形坐标的菜单;(2)含有一个按钮控件,它的作用能够对一个文件的数据进行多项式曲线拟合;(3)文件名通过一个编辑控件由用户给定,给定文件内包含要拟合曲线的数据;(4)拟合好的多项式曲线能够在另一个坐标控件中显示;(5)拟合好的曲线与实际数据曲线用不同的颜色并加各种必要标注在坐标中显示。

2)创新要求:GUI界面使程序更加友好、美观和合理3)课程设计论文编写要求(1)要按照课程设计模板的规格书写课程设计论文(2)论文包括目录、正文、心得体会、参考文献等(3)课程设计论文用B5纸统一打印,装订按学校的统一要求完成4)答辩与评分标准:(1)完成原理分析:20分;(2)完成设计过程:40分;(3)完成调试:20分;(4)回答问题:20分;5)参考文献:(1)刘卫国.MATLAB程序设计与应用(第二版). 北京:高等教育出版社,2008.(2)刘志刚.电力电子学.北京:清华大学出版社、北京交通大学出版社,2004.(3)李传琦. 电力电子技术计算机仿真实验.电子工业出版社,2006.6)课程设计进度安排内容天数地点构思及收集资料2图书馆编程设计与调试1实验室撰写论文2图书馆、实验室学生签名:2011 年月日课程设计(论文)评审意见(1)完成原理分析(20分):优()、良()、中()、一般()、差();(2)设计分析(20分):优()、良()、中()、一般()、差();(3)完成调试(20分):优()、良()、中()、一般()、差();(4)翻译能力(20分):优()、良()、中()、一般()、差();(5)回答问题(20分):优()、良()、中()、一般()、差();(6)格式规范性及考勤是否降等级:是()、否()(7) 总评分数优()、良()、中()、一般()、差();评阅人:职称:讲师2011年月日Matlab应用课程设计目录课程设计任务书 (1)一、Matlab 软件简介 (4)1.1 MATLAB产生的历史背景 (4)1.2 MATLAB的语言特点和开发环境 (4)1.3 基本语法 (5)二、URI简介 (8)2.1特点 (8)2.2组成部分 (8)2.3实现方法 (10)三、设计题目 (10)四、设计内容 (10)4.1Jacobi迭代法基本原理 (10)4.2实验内容 (11)4.3实验结果 (19)五、课程设计心得 (21)六、参考文献 (21)一、Matlab 软件简介MATLAB是美国MathWorks公司生产的一个为科学和工程计算专门设计的交互式大型软件,是一个可以完成各种精确计算和数据处理的、可视化的、强大的计算工具。

matlab迭代法解线性方程组

matlab迭代法解线性方程组

function x=ak(a,b)%a为系数矩阵,b为初始向量(默认为零向量)
%e为精度(默认为1e-6),N为最大迭代次数(默认为100),x为返回解向量
n=length(b);
N=100;
e=1e-6;
x0=zeros(n,1)
%生成一n*1阶零矩阵
x=x0;
x0=x+2*e;
k=0;
d=diag(diag(a));
%生成一个除对角线上元素不为零外其他元素皆为零的矩阵d,且d对角线上的元素为矩阵a 对角线上的元素
l=-tril(a,-1);
%生成一个下三角矩阵
u=-triu(a,1);
%生成一个上三角矩阵
while norm(x0-x,inf)>e & k<N %norm(x0-x,inf)为矩阵(x0-x)的无穷范数
k=k+1;
x0=x;
x=inv(d)*(l+u)*x+inv(d)*b;%雅可比迭代公式k
disp(x')
end
if k==N warning('已达最大迭代次数'); end
function X=BDD(f,x0,TOL)
%X用来存储迭代过程所有的根;
%f是符合不动点迭代要求的迭代方程;%x0设定的迭代初值;
%TOL允许的误差值;
x=feval(f,x0);
n=1;
X(:,n)=x;
while abs(x-x0)>TOL
x0=x;
x=feval(f,x0);
n=n+1;
X(:,n)=x;
end。

jacobi迭代法matlab

jacobi迭代法matlab

jacobi迭代法matlabJacobi迭代法是一种常用的线性方程组求解方法,它是一种迭代法,通过不断迭代来逼近线性方程组的解。

Jacobi迭代法的基本思想是将线性方程组的系数矩阵分解为一个对角矩阵和一个非对角矩阵的和,然后通过迭代求解对角矩阵和非对角矩阵的乘积,最终得到线性方程组的解。

Jacobi迭代法的具体步骤如下:1. 将线性方程组的系数矩阵A分解为一个对角矩阵D和一个非对角矩阵R的和,即A=D+R。

2. 将线性方程组的右端向量b分解为一个对角矩阵D和一个非对角矩阵N的乘积,即b=Dx。

3. 对于任意的初始解向量x0,计算下一次迭代的解向量x1,即x1=D^(-1)(b-Rx0)。

4. 重复步骤3,直到达到预定的精度或迭代次数。

Jacobi迭代法的优点是简单易懂,易于实现,收敛速度较快。

但是,它的缺点也很明显,即收敛速度较慢,需要进行大量的迭代才能达到较高的精度。

在Matlab中,可以使用以下代码实现Jacobi迭代法:function [x,k]=jacobi(A,b,x0,tol,maxit)% Jacobi迭代法求解线性方程组Ax=b% 输入:系数矩阵A,右端向量b,初始解向量x0,精度tol,最大迭代次数maxit% 输出:解向量x,迭代次数kn=length(b); % 系数矩阵A的阶数D=diag(diag(A)); % 对角矩阵DR=A-D; % 非对角矩阵Rx=x0; % 初始解向量for k=1:maxitx1=D\(b-R*x); % 计算下一次迭代的解向量if norm(x1-x)<tol % 判断是否达到精度要求break;endx=x1; % 更新解向量end输出结果可以使用以下代码实现:A=[4 -1 0; -1 4 -1; 0 -1 4]; % 系数矩阵b=[15; 10; 10]; % 右端向量x0=[0; 0; 0]; % 初始解向量tol=1e-6; % 精度要求maxit=1000; % 最大迭代次数[x,k]=jacobi(A,b,x0,tol,maxit); % Jacobi迭代法求解线性方程组fprintf('解向量x=[%f; %f; %f]\n',x(1),x(2),x(3)); % 输出解向量fprintf('迭代次数k=%d\n',k); % 输出迭代次数以上就是Jacobi迭代法的主要内容,通过Matlab实现Jacobi迭代法可以更好地理解其基本思想和具体步骤。

matlab 解线性方程组的迭代法

matlab 解线性方程组的迭代法
MATLAB的遗传算法与直接搜索工具箱(Genetic Algorithm and Direct Search Toolbox,简称GADS)是 MATLAB的一个优化工具箱。它有两种使用方式:一种是 通过命令行调用ga函数,另一种是通过图形界面调用。
小结
➢ 线性方程组求根方法的几何意义
➢ 线性方程组求根函数的理解与应用
设线性代数方程组为
展开为
若对角元素 逐一变量分离得方程组

此即为迭代公式
简单迭代解法的过程如下:
1 设定一组初值 2 第一次迭代:
得到
第k次迭代 第i个变量
3 第二次迭代: 得到
4 同样做法,得到第k+1次迭代:
迭代次数k的取值与精度要求有关,按下式判断:
若满足则停止迭代 为了便于编程,迭代公式可改写为:
matlab 解线性方程组的 迭代法
2020年4月22日星期三
第十讲 解线性方程组的迭代解法
内容提要
引言 简单迭代法 赛得尔迭代法 迭代解法的收敛性 MATLAB的线性方程组求解函数2 小结
1、引言
迭代解法的基本思想
根据给定方程组,设计出一个迭代公式,构造一 数组的序列 ,代入迭代公式,计算出 ,再代 入迭代公式,经过k次迭代运算后得到 ,若 收敛于某一极限数组xi,则xi就是方程组的近似解。
while(norm(x-x1)>eps) x1=x; x=(I-A)*x1+b; n = n + 1; if(n>=M) disp('Warning: 迭代次数太多,现
在退出!'); return;
end end
例:求解方程组
clear all; A =[ 1.0170 -0.0092 0.0095;

matlab jacobi迭代法代码

matlab jacobi迭代法代码

matlab jacobi迭代法代码Matlab是一种常用的数学软件,它具有强大的矩阵计算和绘图功能。

在数值计算中,迭代法是一种重要的求解方法。

本文将介绍如何使用Matlab实现Jacobi迭代法,并运用实例来说明其应用。

Jacobi迭代法是一种经典的迭代法,用于解线性方程组。

它的基本思想是通过迭代逐步逼近方程组的解。

具体而言,对于线性方程组Ax=b,Jacobi迭代法通过以下步骤进行计算:1. 将方程组表示为x=D^(-1)(L+U)x+b的形式,其中D为A的对角矩阵,L为A的严格下三角矩阵,U为A的严格上三角矩阵。

2. 初始化解向量x^(0)为一个初始猜测值,通常取零向量。

3. 根据迭代公式x^(k+1)=D^(-1)(b-(L+U)x^(k)),计算下一迭代解x^(k+1)。

4. 重复步骤3,直到解向量收敛于方程组的解。

下面是一个使用Matlab实现Jacobi迭代法的示例代码:```matlabfunction x = Jacobi(A, b, maxIter, tolerance)n = size(A, 1);x = zeros(n, 1);xPrev = x;iter = 0;while iter < maxIterfor i = 1:nsigma = A(i, 1:i-1) * xPrev(1:i-1) + A(i, i+1:n) * xPrev(i+1:n);x(i) = (b(i) - sigma) / A(i, i);endif norm(x - xPrev) < tolerancebreak;endxPrev = x;iter = iter + 1;endend```在上面的代码中,函数Jacobi接受四个参数:系数矩阵A,右侧常数向量b,最大迭代次数maxIter和收敛容限tolerance。

函数返回解向量x。

在迭代过程中,我们使用了一个for循环来更新解向量x的每个分量。

gauss-seidel迭代法matlab代码

gauss-seidel迭代法matlab代码

Gauss-Seidel迭代法是一种用于解线性方程组的数值方法,特别适用于稀疏矩阵。

以下是一个使用Matlab实现Gauss-Seidel迭代法的简单示例代码:```matlabfunction [x, iteration] = gaussSeidel(A, b, tol, maxIter)% 输入参数:% A:系数矩阵% b:右侧常数向量% tol:迭代收敛容差% maxIter:最大迭代次数n = length(b);x = zeros(n, 1); % 初始化解向量iteration = 0;while iteration < maxIterx_new = x; % 存储前一次迭代的解for i = 1:n% 计算新的解x(i) = (b(i) - A(i,1:i-1)*x(1:i-1) - A(i,i+1:n)*x_new(i+1:n)) / A(i,i);end% 计算迭代误差error = norm(x - x_new, inf);if error < tolreturn;enditeration = iteration + 1;endwarning('Gauss-Seidel迭代未收敛到指定容差。

');end```使用这个函数时,您需要提供系数矩阵A、右侧常数向量b、迭代收敛容差tol和最大迭代次数maxIter作为输入参数。

函数将返回解向量x和迭代次数iteration。

示例用法:```matlabA = [4, -1, 0; -1, 4, -1; 0, -1, 3];b = [12; -1; 0];tol = 1e-6;maxIter = 1000;[x, iteration] = gaussSeidel(A, b, tol, maxIter);fprintf('解向量x = \n');disp(x);fprintf('迭代次数= %d\n', iteration);```这将求解线性方程组Ax = b,并返回解向量x以及迭代次数。

matlab中jacobi迭代法

matlab中jacobi迭代法

一、简介Matlab中jacobi迭代法是一种用于求解线性方程组的迭代方法,适用于系数矩阵为对称、正定矩阵的情况。

该迭代方法通过将系数矩阵分解为对角矩阵、上三角矩阵和下三角矩阵的形式,然后通过迭代计算得到方程组的解。

在Matlab中,可以利用矩阵运算和迭代循环来实现jacobi迭代法。

二、 jacobi迭代法原理1. 基本思想jacobi迭代法的基本思想是将系数矩阵分解为对角矩阵D、上三角矩阵U和下三角矩阵L的形式,即A=D+L+U,其中D为系数矩阵A 的对角线元素组成的对角矩阵,L为系数矩阵A的下三角部分,U为系数矩阵A的上三角部分。

令x为方程组的解向量,b为方程组的右端向量,则方程组可表示为Ax=b。

根据方程组的性质,可将方程组表示为(D+L+U)x=b,然后利用迭代的方式逐步逼近方程组的解。

2. 迭代公式假设迭代到第k次,方程组可表示为(D+L+U)x=b,将其转化为迭代形式x(k+1)=(D+L)^(-1)(b-Ux(k)),利用迭代公式可以逐步计算出方程组的解。

3. 收敛条件对于jacobi迭代法,收敛条件为系数矩阵A为对角占优矩阵或正定矩阵。

如果满足这一条件,迭代计算会逐步收敛于方程组的解。

三、 Matlab中jacobi迭代法实现在Matlab中,可以利用矩阵运算和迭代循环来实现jacobi迭代法。

具体步骤如下:1. 对系数矩阵进行分解将系数矩阵A分解为对角矩阵D、上三角矩阵U和下三角矩阵L的形式。

2. 初始化迭代变量初始化迭代的初始值x0、迭代次数k、逐次逼近解向量x(k+1)。

3. 迭代计算利用迭代公式x(k+1)=(D+L)^(-1)(b-Ux(k))来逐步计算出方程组的解。

4. 判断收敛条件在迭代计算过程中,需要实时判断迭代计算是否满足收敛条件,如果满足则停止迭代计算,得到方程组的解。

四、实例分析假设有如下方程组:2x1 + x2 + 4x3 = 103x1 + 4x2 - x3 = 10x1 + 2x2 + 3x3 = 0可以利用jacobi迭代法来求解该方程组,在Matlab中可以通过编程实现迭代计算过程。

基于Matlab的解线性方程组的几种迭代法的实现及比较

基于Matlab的解线性方程组的几种迭代法的实现及比较

基于Matlab的解线性方程组的几种迭代法的实现及比较线性方程组的解法有很多种,其中一类常用的方法是迭代法。

迭代法根据一个初值逐步逼近方程组的解,在每一次迭代中利用现有的信息产生新的近似值,并不断地修正。

下面介绍基于Matlab的三种迭代法:雅可比迭代法、高斯-赛德尔迭代法和超松弛迭代法,并进行比较。

1. 雅可比迭代法雅可比迭代法是迭代法中最简单的一种方法。

对于线性方程组Ax=b,雅可比迭代法的迭代公式为:x_{i+1}(j)=1/a_{jj}(b_j-\\sum_{k=1,k\eq j}^n a_{jk}x_i(k))其中,i表示迭代次数,j表示未知数的下标,x_i表示第i次迭代的近似解,a_{jk}表示系数矩阵A的第j行第k列元素,b_j 表示方程组的常数项第j项。

在Matlab中,可以使用以下代码实现雅可比迭代:function [x,flag]=jacobi(A,b,X0,tol,kmax)n=length(b);x=X0;for k=1:kmaxfor i=1:nx(i)=(b(i)-A(i,:)*x+A(i,i)*x(i))/A(i,i);endif norm(A*x-b)<tolflag=1;returnendendflag=0;return其中,参数A为系数矩阵,b为常数项列向量,X0为初值列向量,tol为迭代误差容许值(默认为1e-6),kmax为最大迭代次数(默认为1000)。

函数返回值x为近似解列向量,flag表示是否满足容许误差要求。

2. 高斯-赛德尔迭代法高斯-赛德尔迭代法是雅可比迭代法的改进。

其基本思想是,每次迭代时,利用已经求出的新解中的信息来更新其他未知数的值。

迭代公式为:x_{i+1}(j)=(1/a_{jj})(b_j-\\sum_{k=1}^{j-1}a_{jk}x_{i+1}(k)-\\sum_{k=j+1}^n a_{jk}x_i(k))与雅可比迭代法相比,高斯-赛德尔迭代法的每一次迭代都利用了前面已求得的近似解,因此可以更快地收敛。

matlab中快速求解xa=b的方法

matlab中快速求解xa=b的方法

matlab中快速求解xa=b的方法在Matlab中,要快速求解线性方程组xa=b,可以使用以下几种方法:1. 直接求解法(\):直接使用斜杠操作符(\)可以求解线性方程组。

例如,对于方程组xa=b,可以直接使用x = A\b来解决,其中A是系数矩阵,b是常数向量。

这种方法使用了高效的LU分解算法,并且能够自动适应方程组的类型(如稀疏矩阵或密集矩阵),因此是一种快速求解线性方程组的常用方法。

2. QR分解法:QR分解是一种将矩阵分解为正交矩阵和上三角矩阵的方法。

在Matlab中,可以使用qr函数对系数矩阵进行QR分解,然后使用这个分解求解线性方程组。

具体而言,可以使用[q,r] = qr(A)将系数矩阵A分解为正交矩阵q和上三角矩阵r,然后使用x = r\(q'*b)求解方程组。

这种方法通常适用于方程组的系数矩阵具有较大的条件数或者方程组数目较多的情况。

3. Cholesky分解法:如果线性方程组的系数矩阵是对称正定的,那么可以使用Cholesky分解来求解方程组。

在Matlab中,可以使用chol函数对系数矩阵进行Cholesky分解,然后使用这个分解求解线性方程组。

具体而言,可以使用R = chol(A)将系数矩阵A分解为上三角矩阵R,然后使用x = R'\(R\b)求解方程组。

Cholesky分解法通常适用于系数矩阵具有良好的性质(如对称正定)的情况。

4. 迭代法:如果线性方程组的系数矩阵是稀疏的,那么可以使用迭代法来求解方程组。

迭代法的基本思想是通过迭代改进解的逼近值。

在Matlab中,可以使用pcg函数(预处理共轭梯度法)或者bicg函数(双共轭梯度法)来求解稀疏线性方程组。

这些函数需要提供一个预处理矩阵,用于加速迭代过程。

预处理矩阵可以根据具体问题进行选择,常见的预处理方法包括不完全LU分解(ilu)和代数多重网格(amg)等。

通过使用上述方法,可以在Matlab中快速求解线性方程组xa=b。

MATLAB代码解线性方程组的迭代法

MATLAB代码解线性方程组的迭代法

MATLAB代码解线性方程组的迭代法解线性方程组的迭代法1.rs里查森迭代法求线性方程组Ax=b的解function[x,n]=rs(A,b,x0,eps,M)if(nargin==3)eps=1.0e-6;%eps表示迭代精度M=10000;%M表示迭代步数的限制值elseif(nargin==4) M=10000;endI=eye(size(A));n=0;x=x0;tol=1;%迭代过程while(tol>eps)x=(I-A)*x0+b;n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!'); return;endend2.crs里查森参数迭代法求线性方程组Ax=b的解function[x,n]=crs(A,b,x0,w,eps,M)if(nargin==4)eps=1.0e-6;%eps表示迭代精度M=10000;%M表示迭代步数的限制值elseif(nargin==5)M=10000;endI=eye(size(A));n=0;x=x0;tol=1;%迭代过程while(tol>eps)x=(I-w*A)*x0+w*b;n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!'); return;endend3.grs里查森迭代法求线性方程组Ax=b的解function[x,n]=grs(A,b,x0,W,eps,M)if(nargin==4)eps=1.0e-6;%eps表示迭代精度M=10000;%M表示迭代步数的限制值elseif(nargin==5)M=10000;endI=eye(size(A));n=0;x=x0;tol=1;%前后两次迭代结果误差%迭代过程while(tol>eps)x=(I-W*A)*x0+W*b;%迭代公式n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!'); return;endend4.jacobi雅可比迭代法求线性方程组Ax=b的解function[x,n]=jacobi(A,b,x0,eps,varargin)if nargin==3eps=1.0e-6;M=200;elseif nargin<3errorreturnelseif nargin==5M=varargin{1};endD=diag(diag(A));%求A的对角矩阵L=-tril(A,-1);%求A的下三角阵U=-triu(A,1);%求A的上三角阵B=D\(L+U);f=D\b;x=B*x0+f;n=1;%迭代次数while norm(x-x0)>=epsx0=x;x=B*x0+f;n=n+1;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!');return;endend5.gauseidel高斯-赛德尔迭代法求线性方程组Ax=b的解function[x,n]=gauseidel(A,b,x0,eps,M)if nargin==3eps=1.0e-6;M=200;elseif nargin==4M=200;elseif nargin<3errorreturn;endD=diag(diag(A));%求A的对角矩阵L=-tril(A,-1);%求A的下三角阵U=-triu(A,1);%求A的上三角阵G=(D-L)\U;f=(D-L)\b;x=G*x0+f;n=1;%迭代次数while norm(x-x0)>=epsx0=x;x=G*x0+f;n=n+1;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!');return;endend6.SOR超松弛迭代法求线性方程组Ax=b的解function[x,n]=SOR(A,b,x0,w,eps,M)if nargin==4eps=1.0e-6;M=200;elseif nargin<4errorreturnelseif nargin==5M=200;endif(w<=0||w>=2)error;return;endD=diag(diag(A));%求A的对角矩阵L=-tril(A,-1);%求A的下三角阵U=-triu(A,1);%求A的上三角阵B=inv(D-L*w)*((1-w)*D+w*U);f=w*inv((D-L*w))*b;x=B*x0+f;n=1;%迭代次数while norm(x-x0)>=epsx0=x;x=B*x0+f;n=n+1;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!');return;endend7.SSOR对称逐次超松弛迭代法求线性方程组Ax=b的解function[x,n]=SSOR(A,b,x0,w,eps,M)if nargin==4eps=1.0e-6;M=200;elseif nargin<4errorreturnelseif nargin==5M=200;endif(w<=0||w>=2)error;return;endD=diag(diag(A));%求A的对角矩阵L=-tril(A,-1);%求A的下三角阵U=-triu(A,1);%求A的上三角阵B1=inv(D-L*w)*((1-w)*D+w*U);B2=inv(D-U*w)*((1-w)*D+w*L);f1=w*inv((D-L*w))*b;f2=w*inv((D-U*w))*b;x12=B1*x0+f1;x=B2*x12+f2;n=1;%迭代次数while norm(x-x0)>=epsx0=x;x12=B1*x0+f1;x=B2*x12+f2;n=n+1;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!');return;endend8.JOR雅可比超松弛迭代法求线性方程组Ax=b的解function[x,n]=JOR(A,b,x0,w,eps,M)if nargin==4eps=1.0e-6;M=10000;elseif nargin==5M=10000;endif(w<=0||w>=2)%收敛条件要求error;return;endD=diag(diag(A));%求A的对角矩阵B=w*inv(D);%迭代过程x=x0;n=0;%迭代次数tol=1;%迭代过程while tol>=epsx=x0-B*(A*x0-b);n=n+1;tol=norm(x-x0);x0=x;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!'); return;endend9.twostep两步迭代法求线性方程组Ax=b的解function[x,n]=twostep(A,b,x0,eps,varargin) if nargin==3eps=1.0e-6;M=200;elseif nargin<3errorreturnelseif nargin==5M=varargin{1};endD=diag(diag(A));%求A的对角矩阵L=-tril(A,-1);%求A的下三角阵U=-triu(A,1);%求A的上三角阵B1=(D-L)\U;B2=(D-U)\L;f1=(D-L)\b;f2=(D-U)\b;x12=B1*x0+f1;x=B2*x12+f2;n=1;%迭代次数while norm(x-x0)>=epsx0=x;x12=B1*x0+f1;x=B2*x12+f2;n=n+1;if(n>=M)disp('Warning:迭代次数太多,可能不收敛!');return;endend10.fastdown最速下降法求线性方程组Ax=b的解function[x,n]=fastdown(A,b,x0,eps)if(nargin==3)eps=1.0e-6;endx=x0;n=0;tol=1;while(tol>eps)%以下过程可参考算法流程r=b-A*x0;d=dot(r,r)/dot(A*r,r);x=x0+d*r;tol=norm(x-x0);x0=x;n=n+1;end11.conjgrad共轭梯度法求线性方程组Ax=b的解function[x,n]=conjgrad(A,b,x0)r1=b-A*x0;p=r1;n=0;for i=1:rank(A)%以下过程可参考算法流程if(dot(p,A*p)< 1.0e-50)%循环结束条件break;endalpha=dot(r1,r1)/dot(p,A*p);x=x0+alpha*p;r2=r1-alpha*A*p;if(r2< 1.0e-50)%循环结束条件break;endbelta=dot(r2,r2)/dot(r1,r1);p=r2+belta*p;n=n+1;end12.preconjgrad预处理共轭梯度法求线性方程组Ax=b的解function[x,n]=preconjgrad(A,b,x0,M,eps)if nargin==4eps=1.0e-6;endr1=b-A*x0;iM=inv(M);z1=iM*r1;p=z1;n=0;tol=1;while tol>=epsalpha=dot(r1,z1)/dot(p,A*p);x=x0+alpha*p;r2=r1-alpha*A*p;z2=iM*r2;belta=dot(r2,z2)/dot(r1,z1);p=z2+belta*p;n=n+1;tol=norm(x-x0);x0=x;%更新迭代值r1=r2;z1=z2;end13.BJ块雅克比迭代法求线性方程组Ax=b的解function[x,N]=BJ(A,b,x0,d,eps,M)if nargin==4eps=1.0e-6;M=10000;elseif nargin<4errorreturnelseif nargin==5M=10000;%参数的默认值endNS=size(A);n=NS(1,1);if(sum(d)~=n)disp('分块错误!');return;endbnum=length(d);bs=ones(bnum,1);for i=1:(bnum-1)bs(i+1,1)=sum(d(1:i))+1;%获得对角线上每个分块矩阵元素索引的起始值endDB=zeros(n,n);for i=1:bnumendb=bs(i,1)+d(i,1)-1;DB(bs(i,1):endb,bs(i,1):endb)=A(bs(i,1):endb,bs(i,1):endb);%求A的对角分块矩阵endfor i=1:bnumendb=bs(i,1)+d(i,1)-1;invDB(bs(i,1):endb,bs(i,1):endb)=inv(DB(bs(i,1):endb,bs(i,1):e ndb));%求A的对角分块矩阵的逆矩阵endN=0;tol=1;while tol>=epsx=invDB*(DB-A)*x0+invDB*b;%由于LB+DB=DB-AN=N+1;%迭代步数tol=norm(x-x0);%前后两步迭代结果的误差x0=x;if(N>=M)disp('Warning:迭代次数太多,可能不收敛!');return;endend14.BGS块高斯-赛德尔迭代法求线性方程组Ax=b的解function[x,N]=BGS(A,b,x0,d,eps,M)if nargin==4eps=1.0e-6;M=10000;elseif nargin<4errorreturnelseif nargin==5M=10000;endNS=size(A);n=NS(1,1);bnum=length(d);bs=ones(bnum,1);for i=1:(bnum-1)bs(i+1,1)=sum(d(1:i))+1;%获得对角线上每个分块矩阵元素索引的起始值endDB=zeros(n,n);for i=1:bnumendb=bs(i,1)+d(i,1)-1;DB(bs(i,1):endb,bs(i,1):endb)=A(bs(i,1):endb,bs(i,1):endb); %求A的对角分块矩阵endLB=-tril(A-DB);%求A的下三角分块阵UB=-triu(A-DB);%求A的上三角分块阵N=0;tol=1;while tol>=epsinvDL=inv(DB-LB);x=invDL*UB*x0+invDL*b;%块迭代公式N=N+1;tol=norm(x-x0);x0=x;if(N>=M)disp('Warning:迭代次数太多,可能不收敛!');return;end15.BSOR块逐次超松弛迭代法求线性方程组Ax=b的解function[x,N]=BSOR(A,b,x0,d,w,eps,M)if nargin==5eps=1.0e-6;M=10000;elseif nargin<5errorreturnelseif nargin==6M=10000;%参数默认值endNS=size(A);n=NS(1,1);bnum=length(d);bs=ones(bnum,1);for i=1:(bnum-1)bs(i+1,1)=sum(d(1:i))+1;%获得对角线上每个分块矩阵元素索引的起始值endDB=zeros(n,n);for i=1:bnumendb=bs(i,1)+d(i,1)-1;DB(bs(i,1):endb,bs(i,1):endb)=A(bs(i,1):endb,bs(i,1):endb); %求A的对角矩阵endLB=-tril(A-DB);%求A的下三角阵UB=-triu(A-DB);%求A的上三角阵N=0;iw=1-w;while tol>=epsinvDL=inv(DB-w*LB);x=invDL*(iw*DB+w*UB)*x0+w*invDL*b;%块迭代公式N=N+1;tol=norm(x-x0);x0=x;if(N>=M)disp('Warning:迭代次数太多,可能不收敛!'); return;endend。

matlab求线性方程组的解

matlab求线性方程组的解

matlab求线性方程组的解求解线性方程分为两种方法–直接法和迭代法常见的方法一共有8种直接法Gauss消去法Cholesky分解法迭代法Jacobi迭代法Gauss-Seidel迭代法超松弛迭代法共轭梯度法Bicg迭代法Bicgstab迭代法这里我从计算代码的角度来解释一下,代码按以下顺序给出。

把方程组直接带入已知条件,就可以得到答案。

适用条件Gauss消去法:求解中小规模线性方程(阶数不过1000),一般用于求系数矩阵稠密而且没有任何特殊结构的线性方程组Cholesky分解法:对称正定方程优先使用,系数矩阵A是n 阶对称正定矩阵Jacobi迭代法非奇异线性方程组,分量的计算顺序没有关系Gauss-Seidel迭代法与Jacobi迭代法相似,但计算的分量不能改变超松弛迭代法Jacobi迭代法和Gauss-Seidel迭代法的加速版,由Gauss-Seidel迭代法改进而来,速度较快共轭梯度法需要确定松弛参数w,只有系数矩阵具有较好的性质时才可以找到最佳松弛因子。

但好处是不用确定任何参数,他是对称正定线性方程组的方法也是求解大型稀疏线性方程组最热门的方法Bicg迭代法本质是用双共轭梯度求解线性方程组的方法,对求解的方程没有正定性要求Bicgstab迭代法本质是用稳定双共轭梯度求解线性方程组的方法,对求解的方程没有正定性要求Gauss消去法第一、二个函数ltri、utri是一定要掌握的,后面的几乎每个函数都要用到ltri简单来说,当Ly=bb,L(非奇异下三角矩阵)已知求yfunction y =ltri(L,b)n=size(b,1);y=zeros(n,1);for j =1:n-1y(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-y(j)*L(j+1:n,j); endy(n)=b(n)/L(n,n);utri简单来说,当Ux=yy,U(非奇异上三角矩阵)已知求xfunction x =utri(U,y)n=size(y,1);x=zeros(n,1);for j = n:-1:2x(j)=y(j)/U(j,j);y(1:j-1)=y(1:j-1)-x(j)*U(1:j-1,j);endx(1)=y(1)/U(1,1);gauss算法,计算时粘贴过去就好function[L,U]=gauss(A)n=size(A,1);for k =1:n-1A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k +1:n)-A(k+1:n,k)*A(k,k+1:n);endL=tril(A,-1)+eye(n);U=triu(A);使用例子已经知道一个线性方程组,这里我就不写出数学形式了,A是系数矩阵,直接把上面写好的函数复制过来在运算就可以。

matlab计算方程组

matlab计算方程组

matlab计算方程组Matlab作为一款试用范围广泛的科学计算软件,其计算方程组的能力也是非常强大的。

在Matlab中,可以通过多种方式计算方程组,比如使用直接法、迭代法、线性方程组求解器等等。

下面将分步骤阐述使用Matlab计算方程组的方法。

一、使用直接法求解直接法是一种将系数矩阵直接求逆再与常数向量相乘的方法,通常在方程组的规模较小时使用。

下面是使用Matlab求解线性方程组的示例代码:```matlab% 定义系数矩阵和常数向量A = [1 2 3; 4 5 6; 7 8 9];b = [3; 6; 9];% 求解方程组x = A\b;disp(x);```这段代码首先定义了一个3x3的系数矩阵A和一个3x1的常数向量b,然后使用反斜线符号来求解方程组。

该符号将A的逆矩阵乘上b,得到解向量x。

二、使用迭代法求解当方程组的规模较大时,直接法的计算量可能会非常大,在这种情况下可以使用迭代法来求解方程组。

迭代法的主要思想是通过反复迭代求解来逼近方程组的解。

常见的迭代法有Jacobi迭代法、Gauss-Seidel迭代法等。

以Jacobi迭代法为例,下面是使用Matlab求解线性方程组的示例代码:```matlab% 定义系数矩阵和常数向量A = [1 2 3; 4 5 6; 7 8 9];b = [3; 6; 9];% 定义Jacobi迭代法函数function [x, k] = jacobi(A, b, x0, tol, max_iter)D = diag(diag(A));L = -tril(A, -1);U = -triu(A, 1);x = x0;for k = 1:max_iterx = inv(D)*(b + L*x + U*x);if norm(A*x - b) < tolreturnendendend% 求解方程组x0 = [0; 0; 0];tol = 1e-6;max_iter = 1000;[x, k] = jacobi(A, b, x0, tol, max_iter);disp(x);```这段代码首先定义了一个3x3的系数矩阵A和一个3x1的常数向量b,然后定义了一个Jacobi迭代法的函数来求解方程组。

三种迭代法matlab程序 数值分析

三种迭代法matlab程序 数值分析


end
• end
• err=abs(norm(X'-P));
• P=X';
• if(err<delta)

break
• end
• end
• X=X';
• err,k
function X=jacobi(A,b,P,delta,max1) %A是n维非奇异阵。%b是n维向量。%P是初值。%delta是误差界。 %max1是给定的迭代最高次数。%X为所求的方程组AX=b的近似解。 N=length(b); for k=1:max1 for j=1:N
X(j)=(b(j)-A(j,[1:j-1,j+1:N])*P([1:j-1,j+1:N]))/A(j,j); end err=abs(norm(X'-P)); P=X'; if(err<delta)
break end end X=X';k,erction X=gseid(A,b,P,delta,max1)
• %A是n维非奇异阵。%b是n维向量。%P是初值。
• %delta是误差界。%max1是给定的迭代最高次数。%X为所求的方程组AX=b的近似解。
• N=length(b);
• for k=1:max1
• for j=1:N

if j==1

X(1)=(b(1)-A(1,2:N)*P(2:N))/A(1,1);

elseif j==N

X(N)=(b(N)-A(N,1:N-1)*(X(1:N-1))')/A(N,N);

else

X(j)=(b(j)-A(j,1:j-1)*X(1:j-1)-A(j,j+1:N)*P(j+1:N))/A(j,j);

数值分析中求解线性方程组的MATLAB程序(6种)

数值分析中求解线性方程组的MATLAB程序(6种)

数值分析中求解线性方程组的MATLAB程序(6种)1.回溯法(系数矩阵为上三角)function X=uptrbk(A,B)%求解方程组,首先化为上三角,再调用函数求解[N,N]=size(A);X=zeros(N,1);C=zeros(1,N+1);Aug=[A B];for p=1:N-1[Y,j]=max(abs(Aug(p:N,p)));C=Aug(p,:);Aug(p,:)=Aug(j+p-1,:);Aug(j+p-1,:)=C;if Aug(p,p)==0'A was singular.No unique solution.'break;endfor k=p+1:Nm=Aug(k,p)/Aug(p,p);Aug(k,p:N+1)=Aug(k,p:N+1)-m*Aug(p,p:N+1);endendD=Aug;X=backsub(Aug(1:N,1:N),Aug(1:N,N+1));2.系数矩阵为下三角function x=matrix_down(A,b)%求解系数矩阵是下三角的方程组n=length(b);x=zeros(n,1);x(1)=b(1)/A(1,1);for k=2:1:nx(k)=(b(k)-A(k,1:k-1)*x(1:k-1))/A(k,k);end3.普通系数矩阵(先化为上三角,在用回溯法)function X=uptrbk(A,B)%求解方程组,首先化为上三角,再调用函数求解[N,N]=size(A);X=zeros(N,1);C=zeros(1,N+1);Aug=[A B];for p=1:N-1[Y,j]=max(abs(Aug(p:N,p)));C=Aug(p,:);Aug(p,:)=Aug(j+p-1,:);Aug(j+p-1,:)=C;if Aug(p,p)==0'A was singular.No unique solution.'break;endfor k=p+1:Nm=Aug(k,p)/Aug(p,p);Aug(k,p:N+1)=Aug(k,p:N+1)-m*Aug(p,p:N+1);endendD=Aug;X=backsub(Aug(1:N,1:N),Aug(1:N,N+1));4.三角分解法function [X,L,U]=LU_matrix(A,B)%A是非奇异矩阵%AX=B化为LUX=B,L为下三角,U为上三角%程序中并没有真正解出L和U,全部存放在A中[N,N]=size(A);X=zeros(N,1);Y=zeros(N,1);C=zeros(1,N);R=1:N;for p=1:N-1[max1,j]=max(abs(A(p:N,p)));C=A(p,:);A(p,:)=A(j+p-1,:);A(j+p-1,:)=C;d=R(p);R(p)=R(j+p-1);R(j+p-1)=d;if A(p,p)==0'A is singular.No unique solution'break;endfor k=p+1:Nmult=A(k,p)/A(p,p);A(k,p)=mult;A(k,p+1:N)=A(k,p+1:N)-mult*A(p,p+1:N);endendY(1)=B(R(1));for k=2:NY(k)=B(R(k))-A(k,1:k-1)*Y(1:k-1);endX(N)=Y(N)/A(N,N);for k=N-1:-1:1X(k)=(Y(k)-A(k,k+1:N)*X(k+1:N))/A(k,k);endL=tril(A,-1)+eye(N)U=triu(A)5.雅克比迭代法function X=jacobi(A,B,P,delta,max1);%雅克比迭代求解方程组N=length(B);for k=1:max1for j=1:NX(j)=(B(j)-A(j,[1:j-1,j+1:N])*P([1:j-1,j+1:N]))/A(j,j);enderr=abs(norm(X'-P));relerr=err/(norm(X)+eps);P=X';if (err<delta)|(relerr<delta)breakendendX=X';k6.盖斯迭代法function X=gseid(A,B,P,delta,max1);%盖斯算法,求解赋初值的微分方程N=length(B);for k=1:max1for j=1:Nif j==1X(1)=(B(1)-A(1,2:N)*P(2:N))/A(1,1);elseif j==NX(N)=(B(N)-A(N,1:N-1)*(X(1:N-1))')/A(N,N);elseX(j)=(B(j)-A(j,1:j-1)*X(1:j-1)-A(j,j+1:N)*P(j+1:N))/A(j,j);endenderr=abs(norm(X'-P));relerr=err/(norm(X)+eps);P=X';if (err<delta)|(relerr<delta)break;endendX=X';k。

jacobi迭代法matlab编程例题

jacobi迭代法matlab编程例题

jacobi迭代法matlab编程例题Jacobi迭代法是一种常用的数值方法,用于求解线性方程组。

它的基本思想是通过迭代逐步逼近方程组的解。

在使用Jacobi迭代法解决线性方程组时,首先需要将方程组的系数矩阵A进行分解,即将A分解为D、L和U三个矩阵的和,其中D是A的对角线矩阵,L是A的下三角矩阵,U是A的上三角矩阵。

然后可以得到迭代公式:X(k+1) = D^(-1) * (B - (L+U) * X(k))其中,X(k)表示第k次迭代的解向量,X(k+1)表示第k+1次迭代的解向量,B是方程组的常数项向量。

下面我们通过一个具体的例子来展示如何使用Matlab编程实现Jacobi迭代法。

假设有如下线性方程组:2x + y + z = 9x + 3y - z = 43x - y + 2z = 8首先,我们可以将这个方程组转换为矩阵形式:Ax = B其中,A = [2 1 1;1 3 -1;3 -1 2]B = [9; 4; 8]然后,我们需要将矩阵A进行分解:D = diag(diag(A)) = [2 0 0;0 3 0;0 0 2]L = -tril(A) + D = [0 0 0;-1 0 0;-3 1 0]U = -triu(A) + D = [0 -1 -1;0 0 1;0 0 0]接下来,我们需要设置迭代的初始解向量X(0),这里可以选择一个任意的初始值。

假设我们选择X(0) = [0; 0; 0]。

然后,我们可以通过迭代公式来逐步逼近方程组的解。

根据公式,我们可以得到如下的迭代过程:X(1) = D^(-1) * (B - (L+U) * X(0))X(2) = D^(-1) * (B - (L+U) * X(1))X(3) = D^(-1) * (B - (L+U) * X(2))...直到满足停止条件,通常可以选择迭代次数或解的相对误差作为停止条件。

在Matlab中,我们可以使用for循环来实现迭代过程,具体代码如下:A = [2 1 1; 1 3 -1; 3 -1 2];B = [9; 4; 8];X = [0; 0; 0]; % 初始解向量D = diag(diag(A));L = -tril(A) + D;U = -triu(A) + D;for k = 1:100 % 设置最大迭代次数为100X = inv(D) * (B - (L+U) * X);% 判断停止条件% 可以根据需要设置不同的停止条件,比如迭代次数或解的相对误差if norm(A*X - B) < 1e-6break;endenddisp('解向量:');disp(X);在上述代码中,我们设置了最大迭代次数为100,并使用了解的相对误差作为停止条件。

matlabsor迭代法求解方程

matlabsor迭代法求解方程

一、概述Matlab是一种常用的数学软件,其sor迭代法是解决线性方程组的常用方法之一。

本文将介绍Matlab中利用sor迭代法求解方程的方法以及其应用。

二、sor迭代法简介1. sor迭代法是一种求解线性方程组的数值方法,其基本思想是利用矩阵的分解和迭代逼近的方式求解方程组。

2. sor迭代法是Jacobi迭代法和Gauss-Seidel迭代法的一种改进,能够加快收敛速度。

3. 对于给定的线性方程组Ax=b,sor迭代法的迭代公式为:x(k+1) = (1-w)x(k) + (w/D)(b-Lx(k+1)-Ux(k))其中,w为松弛因子,D为对角矩阵,L为下三角矩阵,U为上三角矩阵。

4. sor迭代法的收敛条件是矩阵A严格对角占优。

三、Matlab中sor迭代法的实现1. 在Matlab中,可以利用sor函数实现sor迭代法求解方程。

2. sor函数的调用格式为:[x,flag,relres,iter,resvec] =sor(A,b,w,tol,maxit)其中,A为系数矩阵,b为右端向量,w为松弛因子,tol为容许误差,maxit为最大迭代次数。

3. sor函数返回求解的近似解x,求解的标志flag,残差relres,迭代次数iter和残差向量resvec。

四、sor迭代法的应用示例1. 示例一:求解5x5线性方程组给定线性方程组Ax=b,其中:A = [4 -1 0 0 0; -1 4 -1 0 0; 0 -1 4 -1 0; 0 0 -1 4 -1; 0 0 0 -1 4] b = [10; 10; 10; 10; 10]利用Matlab的sor函数求解该线性方程组,设置松弛因子w=1.2,容许误差tol=1e-6,最大迭代次数maxit=100。

调用sor函数,得到近似解x,收敛标志flag,残差relres,迭代次数iter和残差向量resvec。

2. 示例二:求解100x100线性方程组给定100x100线性方程组Ax=b,利用sor迭代法求解并比较其结果与精确解的误差。

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

第三章 大作业
考虑线性方程组
12345615310002213131002201310010013101130013122150001322x x x x x x ⎛⎫⎛⎫- ⎪ ⎪ ⎪ ⎪⎛⎫ ⎪ ⎪-- ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪-- ⎪ ⎪= ⎪ ⎪ ⎪-- ⎪ ⎪ ⎪ ⎪ ⎪ ⎪-- ⎪ ⎪ ⎪ ⎪⎝⎭ ⎪ ⎪ ⎪ ⎪-⎝⎭⎝⎭
(1)分别用Jacobi 方法和G-S 方法求解上述方程组,并比较它们的收敛快慢;
(2)用SOR 方法求解,试分析松弛因子ω的选取对方法收敛的影响,并找到最佳松弛因子
J 法:
function [X]=jf(A,B,X)
m=rank(A);
X=reshape(X,length(X),1);
B=reshape(B,length(B),1);
k=1;
%
if m~=size(A,1) %求行数
disp('A 不可逆');
return
end
%
Bj=zeros(m,m);
%
while k<m
g=max_(A(k:m,k)')+k-1; %求列主元行数
R=A(g,:);A(g,:)=A(k,:);A(k,:)=R; %交换第k 行和主元行
M=B(g);B(g)=B(k);B(k)=M; %交换B 中元素
k=k+1;
end
%
k=1;
while k<=m
Bj(k,1:k-1)=-A(k,1:k-1)/A(k,k); Bj(k,1+k:m)=-A(k,k+1:m)/A(k,k); f(k)=B(k)/A(k,k);
k=k+1;
end
f=f';
%
k=0;
t=norm((X-f),inf);
while t>0.00000001 %回代求解
X1=Bj*X+f;
t=norm((X1-X),inf);
X=X1;
k=k+1;
end;
X=X';
disp('迭代次数')
k
GS法:
function [X]=gs1(A,B,X)
m=rank(A);
X=reshape(X,length(X),1);
B=reshape(B,length(B),1);
k=1;
%
if m~=size(A,1) %求行数
disp('A不可逆');
return
end
%
%
while k<=m
g=max_(A(k:m,k)')+k-1; %求列主元行数
R=A(g,:);A(g,:)=A(k,:);A(k,:)=R; %交换第k行和主元行M=B(g);B(g)=B(k);B(k)=M; %交换B中元素
k=k+1;
end
%
k=1;
while k<=m
G(k:m,k)=A(k:m,k); %G=D-L
U(k,k:m)=-A(k,k:m);
U(k,k)=0;
k=k+1;
end
Bgs=G\U;
fgs=G\B;
%
k=0;
t=norm((X-fgs),inf);
while t>0.00000001 %回代求解
X1=Bgs*X+fgs;
t=norm((X1-X),inf);
X=X1;
k=k+1;
end;
X=X';
disp('迭代次数')
k
SOR法
function [X]=sor1(A,B,X,w)
m=rank(A);
X=reshape(X,length(X),1);
B=reshape(B,length(B),1);
k=1;
%
if m~=size(A,1) %求行数
disp('A不可逆');
return
end
%
D=zeros(m,m);
%
while k<=m
g=max_(A(k:m,k)')+k-1; %求列主元行数
R=A(g,:);A(g,:)=A(k,:);A(k,:)=R; %交换第k行和主元行M=B(g);B(g)=B(k);B(k)=M; %交换B中元素
k=k+1;
end
%
k=1;
while k<=m
D(k,k)=A(k,k);
L(k:m,k)=-A(k:m,k);
U(k,k:m)=-A(k,k:m);
k=k+1;
end
L=L+D;
U=U+D;
Bw=(D-w*L)\((1-w)*D+w*U);
fw=w*((D-w*L)\B);
%
k=0;
t=norm((X-fw),inf);
while t>0.00000001 %回代求解X1=Bw*X+fw;
t=norm((X1-X),inf);
X=X1;
k=k+1;
end;
X=X';
disp('迭代次数')
k
最佳
function w=opt(A)
%opt
disp('best')
m=rank(A);
k=1;
%
if m~=size(A,1) %求行数
disp('A不可逆');
return
end
%
Bj=zeros(m,m);
%
while k<m
g=max_(A(k:m,k)')+k-1; %求列主元行数
R=A(g,:);A(g,:)=A(k,:);A(k,:)=R; %交换第k行和主元行k=k+1;
end
%
k=1;
while k<=m
Bj(k,1:k-1)=-A(k,1:k-1)/A(k,k);
Bj(k,1+k:m)=-A(k,k+1:m)/A(k,k);
k=k+1;
end
%
[x y]=eig(Bj);
r=max(sqrt(max(abs(y))));
w=2/(1+sqrt(1-r^2));。

相关文档
最新文档