实验五 线性代数方程组的数值解matlab代码
matlab求解方程组代码
matlab求解方程组代码
要在MATLAB中求解方程组，你可以使用`linsolve`函数或者反斯密特正交分解（QR分解）来求解线性方程组。
假设你有一个形如Ax = b的线性方程组，其中A是系数矩阵，x是未知向量，b是常数向量。
首先，使用`linsolve`函数可以直接求解线性方程组。
例如，如果你有一个3x3的系数矩阵A和一个3x1的常数向量b，你可以这样做：
matlab.
A = [1 2 3; 4 5 6; 7 8 10];
b = [3; 6; 10];
x = linsolve(A, b);
另一种方法是使用QR分解来求解方程组。
你可以使用MATLAB 中的`qr`函数来进行QR分解，然后使用得到的分解来求解方程组。
这是一个示例代码：
matlab.
A = [1 2 3; 4 5 6; 7 8 10];
b = [3; 6; 10];
[Q, R] = qr(A);
y = Q'b;
x = R\y;
以上是两种常见的方法，你可以根据具体情况选择合适的方法来求解你的线性方程组。
希望这些信息能帮助到你。
数值分析中求解线性方程组的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。
实验5_线性代数方程组的数值解法
实验5 线性代数方程组的数值解法化工系 毕啸天 2010011811【实验目的】1. 学会用MATLAB 软件数值求解线性代数方程组，对迭代法的收敛性和解的稳定性作初步分析；2. 通过实例学习用线性代数方程组解决简化的实际问题。
【实验内容】题目3已知方程组Ax=b ，其中A ，定义为⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡--------------=32/14/12/132/14/14/12/132/14/14/12/132/14/12/13A试通过迭代法求解此方程组，认识迭代法收敛的含义以及迭代初值和方程组系数矩阵性质对收敛速度的影响。
实验要求：（1）选取不同的初始向量x (0)和不同的方程组右端项向量b ，给定迭代误差要求，用雅可比迭代法和高斯-赛德尔迭代法计算，观测得到的迭代向量序列是否均收敛？若收敛，记录迭代次数，分析计算结果并得出你的结论；（2）取定右端向量b 和初始向量x (0)，将A 的主对角线元素成倍增长若干次，非主对角线元素不变，每次用雅可比迭代法计算，要求迭代误差满足 ，比较收敛速度，分析现象并得出你的结论。
3.1 模型分析选取初始向量x(0) =(1,1,…,1)T ，b=(1,1,…,1)T ，迭代要求为误差满足 ，编写雅各比、高斯-赛德尔迭代法的函数，迭代求解。
3.2 程序代码function x = Jacobi( x0,A,b,m ) D=diag(diag(A)); U=-triu(A,1); L=-tril(A,-1); B1=D\(L+U); f1=D\b; x(:,1)=x0;x(:,2)=B1*x(:,1)+f1; k=1;while norm((x(:,k+1)-x(:,k)),inf)>m x(:,k+2)=B1*x(:,k+1)+f1; k=k+1; endendfunction x = Gauss( x0,A,b,m )D=diag(diag(A));U=-triu(A,1);L=-tril(A,-1);B2=(D-L)\U;f2=(D-L)\b;x(:,1)=x0;x(:,2)=B2*x(:,1)+f2;k=1;while norm((x(:,k+1)-x(:,k)),inf)>mx(:,k+2)=B2*x(:,k+1)+f2;k=k+1;endendA1=3.*eye(20,20);A2=sparse(1:19,2:20,-1/2,20,20);A3=sparse(1:18,3:20,-1/4,20,20);AA=A1+A2+A3+A2'+A3';A=full(AA);b=ones(20,1); %输入自选右端项向量b x0=ones(20,1); %输入自选初始向量x0 m=1e-5;x1=Jacob(x0,A,b,m);x2=Gauss(x0,A,b,m);结果输出数据：3.3.1将x0各分量初值置为0【分析】从数据中可以看出，当迭代的初值变化了，达到相同精度所需要的迭代次数也变化了。
线性方程组的数值算法C语言实现(附代码)
线性方程组AX=B 的数值计算方法实验一、 实验描述：随着科学技术的发展，线性代数作为高等数学的一个重要组成部分，在科学实践中得到广泛的应用。
本实验的通过C 语言的算法设计以及编程，来实现高斯消元法、三角分解法和解线性方程组的迭代法（雅可比迭代法和高斯-赛德尔迭代法），对指定方程组进行求解。
二、 实验原理：1、高斯消去法：运用高斯消去法解方程组，通常会用到初等变换，以此来得到与原系数矩阵等价的系数矩阵，达到消元的目的。
初等变换有三种：(a)、(交换变换)对调方程组两行；(b)、用非零常数乘以方程组的某一行；(c)、将方程组的某一行乘以一个非零常数，再加到另一行。
通常利用(c)，即用一个方程乘以一个常数，再减去另一个方程来置换另一个方程。
在方程组的增广矩阵中用类似的变换，可以化简系数矩阵，求出其中一个解，然后利用回代法，就可以解出所有的解。
2、选主元：若在解方程组过程中，系数矩阵上的对角元素为零的话，会导致解出的结果不正确。
所以在解方程组过程中要避免此种情况的出现，这就需要选择行的判定条件。
经过行变换，使矩阵对角元素均不为零。
这个过程称为选主元。
选主元分平凡选主元和偏序选主元两种。
平凡选主元：如果()0p pp a ≠，不交换行；如果()0p pp a =，寻找第p 行下满足()0p pp a ≠的第一行，设行数为k ，然后交换第k 行和第p 行。
这样新主元就是非零主元。
偏序选主元：为了减小误差的传播，偏序选主元策略首先检查位于主对角线或主对角线下方第p 列的所有元素，确定行k ，它的元素绝对值最大。
然后如果k p >，则交换第k 行和第p 行。
通常用偏序选主元，可以减小计算误差。
3、三角分解法：由于求解上三角或下三角线性方程组很容易所以在解线性方程组时，可将系数矩阵分解为下三角矩阵和上三角矩阵。
其中下三角矩阵的主对角线为1，上三角矩阵的对角线元素非零。
有如下定理：如果非奇异矩阵A 可表示为下三角矩阵L 和上三角矩阵U 的乘积： A LU = (1) 则A 存在一个三角分解。
实验五(线性方程组的数值解法和非线性方程求解)
1大学数学实验 实验报告 | 2014/4/5一、 实验目的1、学习用Matlab 软件数值求解线性代数方程组，对迭代法的收敛性和解的稳定性作初步分析；2、通过实例学习用线性代数方程组解决简化问题。
二、 实验内容项目一：种群的繁殖与稳定收获：种群的数量因繁殖而增加，因自然死亡而减少，对于人工饲养的种群（比如家畜）而言，为了保证稳定的收获，各个年龄的种群数量应维持不变。
种群因雌性个体的繁殖而改变，为方便起见以下种群数量均指其中的雌性。
种群年龄记作k=1,2,…,n ,当年年龄k 的种群数量记作x k ，繁殖率记作b k （每个雌性个体1年的繁殖的数量），自然存活率记作s k （s k =1−d k ，d k 为1年的死亡率），收获量记作ℎk ，则来年年龄k 的种群数量x ̌k 应该为x ̌k =∑b k n k=1x k , x ̌k+1=s k x k −ℎk , (k=1,2,…,n -1)。
要求各个年龄的种群数量每年维持不变就是要求使得x ̌k =x k , (k=1,2,…,n -1).(1) 如果b k , s k 已知，给定收获量ℎk ，建立求各个年龄的稳定种群数量x k 的模型（用矩阵、向量表示）.(2) 设n =5,b 1=b 2=b 5=0,b 3=5,b 4=3,s 1=s 4=0.4,s 2=s 3=0.6,如要求ℎ1~ℎ5为500,400,200,100,100，求x 1~x 5.(3) 要使ℎ1~ℎ5均为500，如何达到？问题分析：该问题属于简单的种群数量增长模型，在一定的条件（存活率，繁殖率等）下为使各年龄阶段的种群数量保持不变，各个年龄段的种群数量将会满足一定的要求，只要找到种群数量与各个参量之间的关系，建立起种群数量恒定的方程就可以求解出各年龄阶段的种群数量。
模型建立：根据题目中的信息，令x ̌k =x k ，得到方程组如下：{x ̌1=∑b k nk=1x k =x 1x ̌k+1=s k x k −ℎk =x k+1整理得到：{−x 1∑b k nk=1x k =0−x k+1+s k x k =ℎk2 大学数学实验 实验报告 | 2014/4/52写成系数矩阵的形式如下：A =[b 1−1b 2b 3s 1−100s 2−1…b n−1b n0000⋮⋱⋮000000000⋯00−10s n−1−1]令h =[0, ℎ1,ℎ2,ℎ3,…,ℎn−2,ℎn−1]Tx =[x n , x n−1,…,x 1]T则方程组化为矩阵形式：Ax =h ，即为所求模型。
数学实验——线性代数方程组的数值解
数学实验——线性代数⽅程组的数值解实验5 线性代数⽅程组的数值解法分1 黄浩 43⼀、实验⽬的1.学会⽤MATLAB软件数值求解线性代数⽅程组，对迭代法的收敛性和解的稳定性作初步分析；2.通过实例学习⽤线性代数⽅程组解决简化的实际问题。
⼆、For personal use only in study and research; not forcommercial use三、四、实验内容1.《数学实验》第⼆版（问题1）问题叙述：通过求解线性⽅程组，理解条件数的意义和⽅程组性态对解的影响，其中是n阶范德蒙矩阵，即是n阶希尔伯特矩阵，b1，b2分别是的⾏和。
（1）编程构造（可直接⽤命令产⽣）和b1，b2；你能预先知道⽅程组和的解吗？令n=5，⽤左除命令求解（⽤预先知道的解可验证程序）。
（2）令n=5,7,9，…，计算和的条件数。
为观察他们是否病态，做以下试验：b1，b2不变，和的元素，分别加扰动后求解；和不变，b1，b2的分量b1(n),b2(n)分别加扰动后求解。
分析A与b的微⼩扰动对解的影响。
取10^-10,10^-8,10^-6。
（3）经扰动得到的解记做，计算误差，与⽤条件数估计的误差相⽐较。
模型转换及实验过程：（1）⼩题.由b1,b2为,的⾏和，可知⽅程组和的精确解均为n ⾏全1的列向量。
在n=5的情况下，⽤matlab编程（程序见四.1），构造,和b1，b2，使⽤⾼斯消去法得到的解x1,x2及其相对误差e1,e2（使⽤excel计算⽽得）为：由上表可见，当n=5时，所得的解都接近真值，误差在10^-12的量级左右。
（2）⼩题分别取n=5,7,9,11,13,15，计算和的条件数c1和c2，（程序见四.2），结果如下：由上表可见，⼆者的条件数都⽐较⼤，可能是病态的。
为证实和是否为病态，先保持b不变，对做扰动，得到该情况下的⾼斯消元解，（程序见四.3），结果如下：(为使结果清晰简洁，在此仅列出n=5,9,13的情况，n=7,11,15略去)=10^-10时：=10^-8时：=10^-6时：由上表可见：a)对于希尔伯特阵，随着阶数的增加，微⼩扰动对解带来的影响越来越⼤，到了n=9时，已经有了6倍误差的解，到了n=13时，甚⾄出现了22倍误差的解元素；⽽随着的增加，解的偏差似乎也有增加的趋势，但仅凭上述表格⽆法具体判断（在下⼀⼩题中具体叙述）。
用MATLAB做线性代数实验
2
0
, 2
5 3
, 3
1
3
, 4
1
4
, 5
1
2
。
3
6
0Hale Waihona Puke 73【程序如下】：
% (1)
A=[1 2 1 3;4 -1 5 6;1 -3 -4 7;1 2 1 1]' r=rank(A) [R,IP]=rref(A) % (2) A=[1 2 0 2 1;-2 -5 1 -1 1;0 -3 3 4 2;P3 6 0 -7 3] r=rank(A) [R,IP]=rref(A)
例如：
已知
A
1 3
2 4
,
B
1 1
2 0
，解矩阵方程
(1)
AX
B , (2) XA B 。
MATLAB 程序如下：
A=[1 2;3 4];
B=[1 2;-1 1];
X1=inv(A)*B % AX=B or
X1=A\B
X2=B*inv(A) % XA=B
X2=A/B
将 p(x) 分解为最简分式之和 q( x)
[p,q]=residue(a,b,r) 将简单分式之和合并为有理分式
例如，将有理分式
f
(x)
x2 x3 2x2 3x2
分解为最简分式之和的程序如下：
p=[1 2];
q=[1 2 3 2];
[a,b,r]=residue(p,q)
输出：a =
-0.2500 - 0.4725i
p=[1 -6 11 -6];
MATLAB解代数方程
例 用Jacobi迭代法求解下列线性方程组。设迭代初值为0， 迭代精度为10-6。 在命令中调用函数文件Jacobi.m，命令如下： A=[10,-1,0;-1,10,-2;0,-2,10]; b=[9,7,6]'; x= [x,n]=jacobi(A,b,[0,0,0]',1.0e-6) 0.9958
x= -66.5556 25.6667 -18.7778 26.5556
(2) QR分解 对矩阵X进行QR分解，就是把X分解为一个正交矩阵Q和一 个上三角矩阵R的乘积形式。QR分解只能对方阵进行。 MATLAB的函数qr可用于对矩阵进行QR分解，其调用格 式为： [Q,R]=qr(X)：产生一个正交矩阵Q和一个上三角矩阵R，使 之满足X=QR。 [Q,R,E]=qr(X)：产生一个一个正交矩阵Q、一个上三角矩阵 R以及一个置换矩阵E，使之满足XE=QR。 实现QR分解后，线性方程组Ax=b的解x=R\(Q\b)或 x=E(R\(Q\b))。
-1.7556 n= n1 =
4
1011
非线性方程数值求解
单变量非线性方程求解 在MATLAB中提供了一个fzero函数，可以用来求单变量 非线性方程的根。该函数的调用格式为： z=fzero('fname',x0,tol,trace) 其中fname是待求根的函数文件名，x0为搜索的起点。一个 函数可能有多个根，但fzero函数只给出离x0最近的那个 根。tol控制结果的相对精度，缺省时取tol=eps，trace• 指 定迭代信息是否在运算中显示，为1时显示，为0时不显示， 缺省时取trace=0。
例 用LU分解求解线性方程组。 命令如下： A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4]; b=[13,-9,6,0]'; [L,U]=lu(A); x=U\(L\b) 采用LU分解的第2种格式，命令如下： [L,U ,P]=lu(A); x=U\(L\P*b)
线性代数方程组的数值解法_百度文库
线性代数方程组的数值解法【实验目的】1. 学会用MATLAB 软件数值求解线性代数方程组，对迭代法的收敛性和解的稳定性作初步分析；2. 通过实例学习用线性代数方程组解决简化的实际问题。
【实验内容】【题目1】通过求解线性方程组A1x=b1和A2x=b2，理解条件数的意义和方程组的性态对解的影响。
其中A1是n阶范德蒙矩阵，即⎡1x0⎢1x1⎢A1=⎢⎢⎢⎣1xn-12x0x12 2xn-1n-1⎤ x0⎥ x1n-1⎥1,...,n-1 ，xk=1+0.1k，k=0,⎥ n-1⎥ xn-1⎥⎦A2是n阶希尔伯特矩阵，b1，b2分别是A1，A2的行和。
（1）编程构造A1（A2可直接用命令产生）和b1，b2；你能预先知道方程组A1x=和A2x=。
b2的解吗？令n=5，用左除命令求解（用预先知道的解可检验程序）b1（2）令n=5,7,9,…，计算A1，A2的条件数。
为观察它们是否病态，做以下试验：b1，b2不变，A1和A2的元素A1(n,n)，A2(n,n)分别加扰动ε后求解；A1和A2不变，b1，b2的分量b1(n)，分析A和b的微小扰动对解的影响。
b2(n)分别加扰动ε求解。
ε取10-1010，-8，10-6。
（3）经扰动得到的解记做x~，计算误差-x~x，与用条件数估计的误差相比较。
1.1构造A1，A2和b1，b2首先令n=5，构造出A1，A2和b1，b2。
首先运行以下程序，输出A1。
运行以下程序对A1,A2求行和：由于b1，b2分别是A1，A2的行和，所以可以预知x1=运行下列程序，用左除命令对b1，b2进行求解：得到以下结果： T。
x2=(1,1, ,1)1.2 计算条件数并观察是否为病态1.不加扰动，计算条件数。
运行以下程序：由此可知，A1，A2的条件数分别是3.574∗10, 4,766∗10。
2.b1，b2不变，A1(n,n)，A2(n,n)分别加扰动（1）n=5时设x11,x12,x13分别为A1添加扰动10−10，10−8,10−6后的解。
利用Matlab求解线性方程组
利⽤Matlab求解线性⽅程组利⽤⾼斯消元法编写了⼀个能够计算线性⽅程组，⽆解，有唯⼀解，⽆穷多解情况的matlab代码。
程序说明：变量n1表⽰系数矩阵或者增⼴矩阵的列数。
当增⼴矩阵的秩与系数矩阵的秩相等时（⽅程有唯⼀解时），n1表⽰系数矩阵的列数。
当⽅程组⽆解或者有⽆数多解时，n1表⽰增⼴矩阵的列数。
处理办法为：if sum(C)~=num1&&j==n1&&flag1==0%系数矩阵在消元过程中，若出现对⾓线及其⼀下元素均为0时，将n1变为增⼴矩阵的列数。
n1=n1+1;%在j等于系数矩阵的列时，n1增加1,变为增⼴矩阵的列。
flag1=1;%flag1保证if内的语句，只执⾏1次。
end当j执⾏到系数矩阵的列n1，且sum(C)~=num1（即系数消元过程中，出现了对⾓线及其⼀下元素均为0，如图1所⽰）时，将n1+1.图1function x=liner_equ_v2(A,b)%该函数⽤于求解线性⽅程组%输⼊参数,A:⽅程组的系数矩阵，b:⽅程组的常数向量（列向量）%输出参数,x:⽅程组的解%时间，2021.10.3%版权所有⼈，zsy%%使⽤实例% A=[1,1,-3,-1;% 3,-1,-3,4;% 1,5,-9,-8];% b=[1;4;0];B=[A,b];%增⼴矩阵[m,n]=size(B);num1=0;for i=1:mnum1=num1+i;endC=zeros(1,n);i=1;j=1;n1=n-1;%系数矩阵或增⼴矩阵的列数flag1=0;while j<=n1if B(i,j)~=0B(i,:)=B(i,:)/B(i,j);for k=i+1:mB(k,:)=B(k,:)-B(k,j)*B(i,:);endC(1,j)=i;if sum(C)~=num1&&j==n1&&flag1==0%系数矩阵在消元过程中，若出现对⾓线及其⼀下元素均为0时，将n1变为增⼴矩阵的列数。
matlab实验5 线性方程组的解法
马千里 热动71 970669实验5 线性方程组的解法实验目的1. 用MA TLAB 软件掌握线性方程组的解法，对迭代法的收敛性和解的稳定性作初步分析。
2. 通过实例练习用线性方程组求解实际问题。
实验内容预备：编写雅可比迭代和高斯—塞得尔迭代的程序 雅可比迭代的程序：function xxx=ykb(A,b) D=diag(diag(A)); U=-triu(A,1); L=-tril(A,-1); B1=inv(D)*(L+U); f=inv(D)*b;x(1:length(b),1)=0; xx=x+1;while norm(x-xx)>1e-6 xx=x;x=B1*x+f; end xxx=x;1. 用MA TLAB 软件的“lu ”(LU 分解)，“ \ ”，以及雅可比迭代和高斯—塞德尔迭代解方程组Ax=b (A 如下，b 任意，比较分析其结果包括迭代法收敛或不收敛的原因）。
a. A=[1,2,-1;1,1,1;2,2,1]; 设b=[1; 2; 1]在MA TLAB 下运行：》A=[1,2,-1;1,1,1;2,2,1];b=[1; 2; 1]; 用LU 分解方法： 》[L U p]=lu(A)L =1.0000 0 00.5000 1.0000 00.5000 0 1.0000》x=(L*U\p*b)’ x =-9 8 3 用除号“\”： 》x=(A\b)’ x =-9 8 3 用雅可比迭代法： 》x=ykb(A,b) x =高斯—塞得尔迭代的程序：function xxx=guass(A,b) D=diag(diag(A)); U=-triu(A,1); L=-tril(A,-1); B1=inv(D-L)*U; f=inv(D-L)*b;x(1:length(b),1)=0; xx=x+1;while norm(x-xx)>1e-8 xx=x;x=B1*x+f; end xxx=x;U = 2.0000 2.0000 1.0000 0 1.0000 -2.5000 0 0 0.5000 p = 0 0 1 1 0 0 0 1 0-9 8 3迭代公式中B1 = f=0 -2 2 1-1 0 -1 2-2 -2 0 1 因为ρ(B1)=0.126<1所以是收敛的用高斯—塞德尔迭代法：因为迭代公式中B1= f=0 -2 2 10 2 -3 10 0 2 3因为ρ(B1)=2>1所以迭代是不收敛的b.A为n阶Hilbert矩阵，n=3~10。
MATLAB求解方程解析解和数值解
辽宁工程技术大学上机实验报告用MATLAB求解质点振动方程振动是日常生活和工程技术中常见的一种运动形式。
利用常系数线性微分方程的理论来讨论有关自由振动和强迫振动的相关问题。
利用MATLAB数学软件大致可分四类情况：(1)无阻尼自由振动情况；(2)有阻尼自由振动；(3)无阻尼强迫振动；(4)有阻尼强迫振动求其数值解和解析解；MATLAB软件求解微分方程解析解的命令“ dsolve() ”求通解的命令格式：(’微分方程’自变量'注：微分方程在输入时，一阶导数y'应输入Dy，y'应输入D2y等，D应大写。
1, 无阻尼自由振动情况：常见的数学摆的无阻尼微小振动方程代码如下：>> t=0:pi/50:2*pi;>> y=2*s in( 3*t+2);>>plot(t,y,'b')2 ! ----------------- B------------------ , ------------------- p ----------------- , ------------------ , ------------------- , -----------------1.5 -\ ；1X ; '| - '■' ■I,°.5 ■-! I-0.5 --1 ■-\ : ■-1.5-2 ........................ ....................................................................................................................0 1 2 3 4 5 6 72, 有阻尼自由振动由无阻尼振动的通解可以看出，无阻尼振动是按照正弦规律运动的，摆动似乎可以无限期的进行下去，但事实上，空气从在阻力，在运动时，我们必须把空气阻力考虑在内，所以我们得到有阻尼摆动方程为：记u/m=2n,g/l=w A2,这里n,w是正常数，所以：y=dsolve('D2y+2* n*Dy+wA2*y=0','t'); ( 4.43)解得：y = C3*exp(-t*(n + ((n + w)*(n - w))A(1/2))) + C2*exp(-t*(n - ((n + w)*(n - w))A(1/2)))(1)小阻尼情形:n<w时，方程(4.43)的通解为：y=exp(-n*t)*(c1*cos(w1*t)+c2*si n(w1*t))和前面无阻尼的情形一样，可以把上式的通解改写为一下形式：y=A*exp(-n*t)*si n(w1*t+Q), (4.45)这里的A，Q为任意常数。
线性代数方程组的求解与Matlab编程
在数学领域中，线性代数方程组经常会需要使用到，这也是线性代数的重要内容。
基础理论知识，矩阵的构成，这些是学习线性代数求解的基本模块,矩阵的初等行（列）变换则是线性方程组求解的基础步骤。
这也延伸出了很多种线性代数方程组的求解办法，比如行阶梯法、矩阵运算法等一系列办法，根据具体线性代数题U的不同，我们在求解中也可以灵活运用不同的方法。
这些方法的不同优点可以更方便的帮助我们解决线性方程组求解中不同类型的题U,在具体求解过程中，矩阵和向量是基础也是重中之重，它使我们解决问题的关键。
我们生活中也有很多问题会用到线性代数方程组方面的知识，所以它与我们的生活也是密不可分的。
MATLAB作为一门科学与计算机方向的重要编程语言，是现如今我们可以在计算机上找到的最优秀的工具之一，它拥有强大的计算能力、图形编程能力、符号推算能力、文字处理能力、动态模拟仿真能力。
Matlab拥有的强大的编程讣算能力，让它可以根据本文中线性方程组的这些求解方法，分别给出在求解线性方程组中的Matlab应用，并且利用计算机对线性代数方程组进行求解。
［关键词］Mat lab线性方程组矩阵求解初等变换摘要........................................................... -1 一ABSTRACT ............................................................................................. 错误!未定义书签。
目录........................................................... -1 -1前言............................................................. - 2 -1.1国内外研究现状............................................. -2 -1.2课题的研究意义............................................. _ 3 -1.3本文所完成的工作............................................ -3 - 2预备知识......................................................... -4 -2.1线性方程组的定义............................................ -4 -2.2线性方程组有解判别定理...................................... -4 -2.3线性方程组解的结构.......................................... -5 -2. 3. 1齐次线性方程组的性质................................. -6 -2. 3. 2基础解系及其存在性................................... -6 -2. 3. 3 一般线性方程组的解的结构............................. -7 - 3线性方程组的求解方法............................................. -8 -3.1行阶梯法解线性方程.......................................... -8 -3. 1. 1初等行变换........................................... -8 -3.1.2行阶梯矩阵的生成（高斯消元过程）...................... - 8 -3.1. 3行阶梯法解欠定方程组................................ -10 -3.2用矩阵运算法解线性方程组 .................................. -10 -3.2.1矩阵运算规则及MATLAB实例........................... - 10 -3.2. 2 行列式.............................................. - 12 -3.2. 3矩阵的秩和矩阵求逆.................................. - 12 -3.2.4条件数一衡量奇异程度的量.............................. -12 -3.2.5用矩阵“除法”解线性方程.............................. -13 -4 Mat lab编程方法及实例 ........................................... -13 -4. 1 Matlab中的行阶梯法解线性方程.............................. -13 -4. 1. 1线性方程的Mat lab表示方法 ........................... - 13 -4.1.2基于Matlab的高斯消元及实例 ........................... - 15 -4. 1. 3 MATLAB中的行阶梯生成函数.......................... -17 -4. 1. 4解欠定方程组实例.................................... - 18 -4. 2 Matlab中的矩阵运算法解线性方程............................ - 20 -4. 2. 1矩阵运算规则的matlab实例............................. -20 -4. 2. 2初等变换乘子矩阵的生成及MATLAB实例................ -22 -4.2.3行列式的Mat lab实例................................... -26 -4. 2. 4矩阵求逆实例......................................... -27 -4. 2. 5矩阵“除法”解欠定方程实例........................... -28 -结论............................................................ 一30 -参考文献......................................................... 一31 一1.1国内外研究现状线性方程组求解这个数学领域，国内外数学学家对此的研究非常广泛，也为这个课题做出了许多贡献，在科技和学术日益进步的今天，数学已经与我们的生活息息相关，影响着社会的进程和科技的发展。
MATLAB线性方程组求解方法
——GW318 物联网实验室学术活动
MATLAB 线性方程组求解方法
1．线性方程组的问题 在工程计算中，一个重要的问题是线性方程组的求解。在矩阵表示法中，上述问题可以 表述为给定两个矩阵 A 和 B，是否存在惟一的解 X 使得 AX=B 或 XA=B。 例如，求解方程 3x=6 就可以将矩阵 A 和 B 看成是标量的一种情况，最后得出该方 程的 解为 x=6/3=2。 尽管在标准的数学中并没有矩阵除法的概念，但 MATLAB 7.0 采用了与解标量方程中 类 似的约定，用除号来表示求解线性方程的解。MATLAB 7.0 采用第 2 章介绍过的运算 符斜杠 ’/’和运算符反斜杠’\ ’来表示求线性方程的解，其具体含义如下： • X=A\B 表示求矩阵方程 AX=B 的解； • X=B/A 表示求矩阵方程 XA=B 的解。 对于 X=A\B，要求矩阵 A 和矩阵 B 有相同的行数，X 和 B 有相同的列数，X 的行 数等于 矩阵 A 的列数。X=B\A 行和列的性质则与之相反。 在实际情况中，形式 AX=B 的线性方程组比形式 XA=B 的线性方程组要常见得多。 因此 反斜杠’\’用得更多。本小节的内容也主要针对反斜杠’\’除法进行介绍。斜杠’/’ 除法的性质可以 由恒等变换式得到(B/A)’=(A’\B’)。 系数矩阵 A 不一定要求是方阵，矩阵 A 可以是 m×n 的矩阵，有如下 3 种情况： • m=n 恰定方程组，MATLAB7.0 会寻求精确解； • m>n 超定方程组，MATLAB7.0 会寻求最小二乘解； • m<n 欠定方程组，MATLAB7.0 会寻求基本解，该解最多有 m 个非零元素。 值得注意的是用 MATLAB 7.0 求解这种问题时，并不采用计算矩阵的逆的方法。针对 不 同的情况，MATLAB7.0 会采用不同的算法来解线性方程组。 2．线性方程组的一般解 线性方程组 AX=B 的一般解给出了满足 AX=B 的所有解。线性方程组的一般解可以通 过 下面的步骤得到。 • 解相应的齐次方程组 AX=0， 求得基础解。 可以使用函数 null()来得到基础解。 语 句 null(A) 返回齐次方程组 AX=0 的一个基础解，其他基础解与 null(A)是线性关系。 • 求非齐次线性方程组 AX=B，得到一个特殊解。 • 非齐次线性方程组 AX=B 的一般解等于基础解的线性组合加上特殊解。 在后面的章节将介绍求非齐次线性方程组 AX=B 特殊解的方法。 3．恰定方程组的求解
matlab 数值解方程
MATLAB 是一个强大的数值计算和数据分析工具，可以用来求解各种数学方程。
以下是一些基本的 MATLAB 函数，用于求解不同类型的方程。
1. **符号方程求解**：使用 `syms` 和 `solve` 函数。
```matlab
syms x
solution = solve(x^2 + 3*x - 4 = 0);
```
2. **数值方程求解**：使用 `fzero` 函数。
```matlab
f = @(x) x^2 - 3; % 定义函数
solution = fzero(f, 1); % 在区间 [1,2] 内寻找零点
```
3. **系统方程求解**：使用 `fsolve` 函数。
```matlab
fun = @(x) [x(1)^2 + x(2)^2 - 4; x(1) + x(2) - 2]; % 定义系统函数
solution = fsolve(fun, [1,1]); % 使用初始值 [1,1] 进行求解
```
4. **非线性方程求解**：使用 `fminbnd` 或 `fminunc` 函数。
```matlab
fun = @(x) x^3 - x - 1; % 定义函数
solution = fminbnd(fun, -10, 10); % 在区间 [-10,10] 内寻找最小值点
```
注意：以上所有的解决方案都需要先定义方程或系统，然后选择合适的初始值和参数范围进行求解。
MATLAB 的文档提供了更详细的例子和说明，可以帮助你理解如何使用这些函数。
线性代数方程组数值解法及MATLAB实现综述
线性代数方程组数值解法及MATLAB 实现综述廖淑芳 20122090 数计学院 12计算机科学与技术1班（职教本科） 一、分析课题随着科学技术的发展，提出了大量复杂的数值计算问题，在建立电子计算机成为数值计算的主要工具以后，它以数字计算机求解数学问题的理论和方法为研究对象。
其数值计算中线性代数方程的求解问题就广泛应用于各种工程技术方面。
因此在各种数据处理中，线性代数方程组的求解是最常见的问题之一。
关于线性代数方程组的数值解法一般分为两大类：直接法和迭代法。
直接法就是经过有限步算术运算，可求的线性方程组精确解的方法（若计算过程没有舍入误差），但实际犹如舍入误差的存在和影响，这种方法也只能求得近似解，这类方法是解低阶稠密矩阵方程组级某些大型稀疏矩阵方程组的有效方法。
直接法包括高斯消元法，矩阵三角分解法、追赶法、平方根法。
迭代法就是利用某种极限过程去逐步逼近线性方程组精确解的方法。
迭代法具有需要计算机的存储单元少，程序设计简单，原始系数矩阵在计算过程始终不变等优点，但存在收敛性级收敛速度问题。
迭代法是解大型稀疏矩阵方程组（尤其是微分方程离散后得到的大型方程组）的重要方法。
迭代法包括Jacobi 法SOR 法、SSOR 法等多种方法。
二、研究课题-线性代数方程组数值解法 一、 直接法 1、 Gauss 消元法通过一系列的加减消元运算，也就是代数中的加减消去法，以使A 对角线以下的元素化为零，将方程组化为上三角矩阵；然后，再逐一回代求解出x 向量。
1.1消元过程1. 高斯消元法(加减消元)：首先将A 化为上三角阵，再回代求解。
11121121222212n n n n nn n a a a b a a a b a a a b ⎛⎫ ⎪ ⎪ ⎪ ⎪⎝⎭ (1)(1)(1)(1)11121311(2)(2)(2)(2)222322(3)(3)(3)3333()()000000n n n n n nn n a a a a b a a a b a a b a b ⎛⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎝⎭ 步骤如下：第一步：1111,2,,i a i i n a -⨯+=第行第行11121121222212n n n n nnn a a a b a a a b a a a b ⎛⎫⎪ ⎪⎪⎪⎝⎭ 111211(2)(2)(2)2222(2)(2)(2)200n n n n n n a a a b a a b a a b ⎛⎫ ⎪⎪ ⎪ ⎪⎝⎭ 第二步：(2)2(2)222,3,,i a i i n a -⨯+=第行第行111211(2)(2)(2)2222(2)(2)(2)200n nn nnn a a a b a a b a a b ⎛⎫⎪ ⎪ ⎪ ⎪⎝⎭11121311(2)(2)(2)(2)222322(3)(3)(3)3333(3)(3)(3)300000n n n n nn n a a a a b a a a b a a b a a b ⎛⎫ ⎪⎪ ⎪ ⎪ ⎪ ⎪⎝⎭类似的做下去，我们有：第k 步：()()k ,1,,k ikk kka i i k n a -⨯+=+第行第行。
MatLab解线性方程组
MatLab解线性方程组当齐次线性方程AX=0,rank(A)=r<n时,该方程有无穷多个解,怎样用MATLAB求它的一个基本解呢?用matlab 中的命令x=null(A, r )即可.其中:r=rank(A)A=[ 1 1 1 1 -3 -1 11 0 0 0 1 1 0-2 0 0 -1 0 -1 -2]用matlab 求解程序为:A=[1 1 1 1 -3 -1 1;1 0 0 0 1 1 0;-2 0 0 -1 0 -1 -2];r=rank(A);y=null(A, r )得到解为:y=[ 0 -1 -1 0-1 2 1 11 0 0 00 2 1 -20 1 0 00 0 1 00 0 0 1]其列向量为Ay=0的一个基本解MatLab解线性方程组一文通！-------------------作者：liguoy（2005-2-3）写在阅读本文前的引子。
一：读者对线性代数与Matlab 要有基本的了解；二：文中的通用exp.m文件，你须把具体的A和b代进去。
一：基本概念1．N级行列式A：|A|等于所有取自不同行不同列的n个元素的积的代数和。
2．矩阵B：矩阵的概念是很直观的，可以说是一张表。
3．线性无关：一向量组(a ,a ,…. a )不线性相关，即没有不全为零的数k ,k ,……kn使得：k1* a +k2* a +…..+kn*an=04. 秩：向量组的极在线性无关组所含向量的个数称为这个向量组的秩。
5．矩阵B的秩：行秩，指矩阵的行向量组的秩；列秩类似。
记：R(B)6．一般线性方程组是指形式： (1)其中x1,x2,…….xn为n个未知数，s为方程个数。
记：A*X=b7．性方程组的增广矩阵：=8. A*X=0 (2)二：基本理论三种基本变换：1，用一非零的数乘某一方程；2，把一个方程的倍数加到另一个方程；3互换两个方程的位置。
以上称初等变换。
消元法（理论上分析解的情况，一切矩阵计算的基础）首先用初等变换化线性方程组为阶梯形方程组，把最后的一些恒等式”0=0”（如果出现的话）去掉，1：如果剩下的方程当中最后的一个等式是零等于一非零数，那么方程组无解；否则有解，在有解的情况下，2：如果阶梯形方程组中方程的个数r等于未知量的个数，那么方程组有唯一的解，3：如果阶梯形方程组中方程的个数r小于是未知量的个数，那么方程组就有无穷个解。
线性方程组求解Matlab程序
线性方程组求解1.直接法Gauss消元法：function x=DelGauss(a,b)% Gauss 消去法[n,m]=size(a);nb=length(b);det=1;%存储行列式值x=zeros(n,1);for k=1:n-1for i=k+1:nif a(k,k)==0returnendm=a(i,k)/a(k,k);for j=k+1:na(i,j)=a(i,j)-m*a(k,j);endb(i)=b(i)-m*b(k);enddet=det*a(k,k);enddet=det*a(n,n);for k=n:-1:1 %回代for j=k+1:nb(k)=b(k)-a(k,j)*x(j);endx(k)=b(k)/a(k,k);endExample:>> A=[1.0170 -0.0092 0.0095;-0.0092 0.9903 0.0136;0.0095 0.0136 0.9898]; >> b=[1 0 1]';>> x=DelGauss(A,b)x =0.9739-0.00471.0010列主元Gauss 消去法：function x=detGauss(a,b) % Gauss 列主元消去法[n,m]=size(a);nb=length(b);det=1;%存储行列式值x=zeros(n,1);for k=1:n-1amax=0;% 选主元for i=k:nif abs(a(i,k))>amaxamax=abs(a(i,k));r=i;endendif amax<1e-10return;endif r>k %交换两行for j=k:nz=a(k,j);a(k,j)=a(r,j);a(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;det=-det;endfor i=k+1:n %进行消元m=a(i,k)/a(k,k);for j=k+1:na(i,j)=a(i,j)-m*a(k,j);endb(i)=b(i)-m*b(k);enddet=det*a(k,k);enddet=det*a(n,n);for k=n:-1:1 %回代for j=k+1:nb(k)=b(k)-a(k,j)*x(j);endx(k)=b(k)/a(k,k);endExample:>> x=detGauss(A,b)x =0.9739-0.00471.0010Gauss-Jorda n 消去法：function x=GaussJacobi(a,b)% Gauss-Jacobi 消去法[n,m]=size(a);nb=length(b);x=zeros(n,1);for k=1:namax=0;% 选主元for i=k:nif abs(a(i,k))>amaxamax=abs(a(i,k));r=i;endendif amax<1e-10return;endif r>k %交换两行for j=k:nz=a(k,j);a(k,j)=a(r,j);a(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;end%进行消元b(k)=b(k)/a(k,k);for j=k+1:na(k,j)=a(k,j)/a(k,k);endfor i=1:nif i~=kfor j=k+1:na(i,j)=a(i,j)-a(i,k)*a(k,j);endb(i)=b(i)-a(i,k)*b(k);endendendfor i=1:nx(i)=b(i);endExample:>> x=GaussJacobi(A,b)x =0.9739-0.00471.0010LU 分解法：function [l,u]=lu(a) %LU 分解n=length(a);l=eye(n);u=zeros(n);for i=1:nu(1,i)=a(1,i);endfor i=2:nl(i,1)=a(i,1)/u(1,1);end for r=2:n%%%%for i=r:nuu=0;for k=1:r-1uu=uu+l(r,k)*u(k,i);endu(r,i)=a(r,i)-uu;end%%%%for i=r+1:nll=0;for k=1:r-1ll=ll+l(i,k)*u(k,r);endl(i,r)=(a(i,r)-ll)/u(r,r);end%%%% function x=lusolv(a,b)End%LU 分解求解线性方程组aX=bif length(a)~=length(b) error('Error in inputing!') return;endn=length(a);[l,u]=lu(a);y(1)=b(1);for i=2:nz=0;for k=1:i-1z=z+l(i,k)*y(k);endy(i)=b(i)-z;endx(n)=y(n)/u(n,n);for i=n-1:-1:1z=0;for k=i+1:nz=z+u(i,k)*x(k);endx(i)=(y(i)-z)/u(i,i);endExample:>> x=lusolv(A,b)x =0.9739 -0.0047 1.0010对称正定矩阵之Cholesky 分解法：function L=Cholesky(A)%对对称正定矩阵 A 进行Cholesky 分解n=length(A); L=zeros(n);for k=1:ndelta=A(k,k);for j=1:k-1delta=delta-L(k,j)A2;endif delta<1e-10return;L(k,k)=sqrt(delta);endfor i=k+1:nL(i,k)=A(i,k);for j=1:k-1L(i,k)=L(i,k)-L(i,j)*L(k,j);endL(i,k)=L(i,k)/L(k,k);endendfunction x=Chol_Solve(A,b)%利用对称正定矩阵之Cholesky 分解求解线性方程组n=length(b);l=Cholesky(A);x=ones(1,n);y=ones(1,n);for i=1:nz=0;for k=1:i-1z=z+l(i,k)*y(k);endy(i)=(b(i)-z)/l(i,i);for i=n:-1:1Ax=bendz=0;for k=i+1:nz=z+l(k,i)*x(k);endx(i)=(y(i)-z)/l(i,i);endExample:>> a=[9 -36 30 ;-36 192 -180;30 -180 180]; >> b=[1 1 1]';>> x=Chol_Solve(a,b)x =1.8333 1.0833 0.7833对称正定矩阵之LDL'分解法：function [L,D]=LDL_Factor(A)%对称正定矩阵 A 进行LDL' 分解n=length(A);L=eye(n);D=zeros(n);d=zeros(1,n);for k=1:nd(k)=A(k,k);for j=1:k-1d(k)=d(k)-L(k,j)*T(k,j);endif abs(d(k))<1e-10return;endfor i=k+1:nT(i,k)=A(i,k);for j=1:k-1T(i,k)=T(i,k)-T(i,j)*L(k,j);endL(i,k)=T(i,k)/d(k);endendD=diag(d);function x=LDL_Solve(A,b)Ax=b %利用对对称正定矩阵 A 进行LDL' 分解法求解线性方程组n=length(b);[l,d]=LDL_Factor(A);for i=2:nz=0;for k=1:i-1z=z+l(i,k)*y(k);endy(i)=b(i)-z;endx(n)=y(n)/d(n,n);for i=n-1:-1:1z=0;for k=i+1:nz=z+l(k,i)*x(k);endx(i)=y(i)/d(i,i)-z; endExample:>> x=LDL_Solve(a,b)1.8333 1.0833 0.78332.迭代法Richardson 迭代法：function [x,n]=richason(A,b,x0,eps,M) %Richardson 法求解线性方程组Ax=b %方程组系数矩阵:A%方程组之常数向量:b%迭代初始向量:X0%e 解的精度控制:eps%迭代步数控制：M%返回值线性方程组的解：x%返回值迭代步数：nif(nargin == 3)eps = 1.0e-6;M = 200;elseif(nargin == 4)M = 200;endI =eye(size(A));x1=x0;x=(I-A)*x0+b;n=1;while(norm(x-x1)>eps)x1=x;x=(I-A)*x1+b;n = n + 1;if(n>=M)disp('Warning: 迭代次数太多，现在退出！');return;endendExample:>> A=[1.0170 -0.0092 0.0095;-0.0092 0.9903 0.0136;0.0095 0.0136 0.9898]; >> b=[1 0 1]';x0=[0 0 0]';>> [x,n]=richason(A,b,x0)x =0.9739-0.00471.0010Jacobi 迭代法：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,-%求 A 的下三角阵1);0.9739U=-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;endendExample:>> [x,n]=Jacobi(A,b,x0)x =0.9739-0.00471.0010Gauss-Seidel 迭代法：function [x,n]=gauseidel(A,b,x0,eps,M) if nargin==3 eps= 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;endendExample:>> [x,n]=gauseidel(A,b,x0)x =0.97394超松驰迭代法：function [x,n]=SOR(A,b,x0,w,eps,M)if nargin==4-0.00471.0010D=diag(diag(A)); %求 A 的对角矩阵>> [x,n]=SOR(A,b,x0,1)x =eps= 1.0e-6;M = 200; elseif nargin<4errorreturnelseif nargin ==5M = 200;end if(w<=0 || w>=2)error;return;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;endendExample:0.9739-0.00471.0010 n =4对称逐次超松驰迭代法：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;endendExample:>> [x,n]=SSOR(A,b,x0,1)x =0.9739-0.00471.0010n =3两步迭代法：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;endendExample:>> [x,n]=twostep(A,b,x0)x =0.9739-0.00471.0010最速下降法：function [x,n]=fastdown(A,b,x0,eps)if(nargin == 3)eps = 1.0e-6;n =end r = b-A*x0;d = dot(r,r)/dot(A*r,r);x = x0+d*r;n=1;while(norm(x-x0)>eps)x0 = x;r = b-A*x0;d = dot(r,r)/dot(A*r,r);x = x0+d*r;n = n + 1;endExample: >> [x,n]=fastdown(A,b,x0) x =0.9739-0.00471.0010共轭梯度法：function [x,n]=conjgrad(A,b,x0)if(nargin == 3)eps = 1.0e-6;end r1 = b-A*x0;p1 = r1;d = dot(r1,r1)/dot(p1,A*p1);x = x0+d*p1;r2 = r1-d*A*p1;f = dot(r2,r2)/dot(r1,r1);p2 = r2+f*p1;n = 1;for(i=1:(rank(A)-1))x0 = x;p1 = p2;r1 = r2;d = dot(r1,r1)/dot(p1,A*p1);x = x0+d*p1;r2 = r1-d*A*p1;f = dot(r2,r2)/dot(r1,r1);p2 = r2+f*p1;n = n + 1;end d = dot(r2,r2)/dot(p2,A*p2);x = x+d*p2;n = n + 1;Example: >> [x,n]=conjgrad(A,b,x0)0.9739-0.00471.0010 n =4预处理的共轭梯度法：当AX=B 为病态方程组时，共轭梯度法收敛很慢。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
当n很大时Hilbert 矩阵呈病态
例：观察Hilbert矩阵的病态性 观察 矩阵的病态性
例. Hx=b, 其中 H=hilb(5), b=[1,…1]T
H=hilb(5), h=rats(H), b=ones(5,1); x=H\b; b(5)=1.1; x1=H\b; [x,x1], n1=cond(H), n2=rcond(H), x x1 1.0e+003 * 0.0050 0.0680 -0.1200 -1.3800 0.6300 6.3000 -1.1200 -9.9400 0.6300 5.0400 cond(H)=4.7661e+005
4 1 1 4 1 0 设 A= 1 4 O 0 O O 1 1 4n×n
b = [ ,2,L ] 1 n
T
n=500;b=[1:n]'; a1=sparse(1:n,1:n,4,n,n); a2=sparse(2:n,1:n-1,1,n,n); a=a1+a2+a2'; tic;x=a\b;t1=toc aa=full(a); tic;xx=aa\b;t2=toc y=sum(x) yy=sum(xx)
线性方程组数值解法的MATLAB实现
4． ilbert 矩 ： h=hilb(n) H ． 阵 输 出h为n阶Hilbert 矩 阵
L 1 1/ 2 L 1/ n 1/ 2 1/ 3 L 1/(n +1) L H= L L 1/ n 1/(n +1) L /(2n −1) 1
例. 解
10x1 + 3x2 + x3 =14
并对系数矩阵 2x1 −10x2 + 3x3 =−5 作LU分解 分解
x1 + 3x2 +10x3 =14
A=[10 3 1;2 -10 3;1 3 10], b=[14 -5 14]', x=A\b, [L1,U1]=lu(A); L1,U1, A1=L1*U1, [L2,U2,P]=lu(A); L2,U2,P, A2=L2*U2, A3=inv(P)*A2 [Q,R]=qr(A) Q*R [U,S,V]=svd(A) U*S*V’
其他相关的MATLAB函数
提取（产生） 1. 提取（产生）对角阵 v=diag(x) 输入向量 ，输出 是以 为对角元素的对角阵； 输入向量x，输出v是以 为对角元素的对角阵； 是以x为对角元素的对角阵 输入矩阵x，输出v是 的对角元素构成的向量 的对角元素构成的向量； 输入矩阵 ，输出 是x的对角元素构成的向量； trace(x) 返回 的迹，即对角元素的和 返回x的迹 的迹， 提取上（ 2. 提取上（下）三角阵 y=triu(x) v=tril(x) v=triu(x,1) v=tril(x,-1) flipud(x) fliplr(x) ) 输入矩阵x，输出 是 的上三角阵 的上三角阵； 输入矩阵 ，输出v是x的上三角阵； 输入矩阵x，输出v是 的下三角阵 的下三角阵； 输入矩阵 ，输出 是x的下三角阵； 同上，但对角元素为0， 同上，但对角元素为 ， 同上，但对角元素为0。 同上，但对角元素为 。 矩阵x上下翻转 矩阵 上下翻转 矩阵x 矩阵x左右翻转
例：
clear，clc A=[1 1/4 0;0 1/2 0;0 1/4 0]; x=[0 9 0.1 0]; [norm(x),norm(x,1),norm(x,inf)] [norm(A),norm(A,1),norm(A,inf)] [cond(A),cond(A,1),cond(A,inf)] rcond(A)
线性方程组数值解法的MATLAB实现
1. 求解 Ax = b 用左除： x = A\ b 。 用左除：
若A为可逆方阵，输出原方程的解x 若A为n×m矩阵（n>m），且ATA可逆，输出原 方程的最小二乘 矩阵 分解 矩阵LU分解
[x,y]=lu(A) 若A可逆且顺序主子式不为零， 输出 为单位下 可逆且顺序主子式不为零， 可逆且顺序主子式不为零 输出x为单位下 三角阵L， 为上三角阵 为上三角阵U， 可逆， 为一交 三角阵 ，y为上三角阵 ，使A=LU； 若A可逆，x为一交 ； 可逆 换阵与单位下三角阵之积. 换阵与单位下三角阵之积 [x,y,p]=lu(A) 若A可逆，输出 为单位下三角阵 ， y为上三 可逆， 为单位下三角阵L， 为上三 可逆 输出x为单位下三角阵 角阵U，p为一交换阵 ，使PA=LU. 为一交换阵P， 角阵 ， 为一交换阵 u =chol(A) 对正定对称矩阵 的Cholesky分解，输出 为上 对正定对称矩阵A的 分解， 分解 输出u为上 三角阵U， 三角阵 ，使A=UTU [Q,R]=qr(A) 正交三角分解，其中，Q为正交阵，R为上三角阵 正交三角分解，其中， 为正交阵 为正交阵， 为上三角阵 [U,S,V]=svd(A) 奇异值分解，即A=USV ’ ，U，V为正交阵，S 奇异值分解， 为正交阵， ， 为正交阵 为对角阵。 为对角阵。
t1, t2相差巨大，说明用稀疏矩阵计算的优点 相差巨大， 相差巨大 说明用稀疏矩阵计算的优点 用于简单地验证两种方法结果的一致） （y=yy 用于简单地验证两种方法结果的一致）
线性方程组数值解法的MATLAB实现
3. 范数 条件数 特征值 n=norm (x, p) 输入x为向量或矩阵，输出为x的p-范 数,p可取1，2，‘inf’，分别为1-范数， 2-范数，∞-范数，p缺省值为2. c=cond (x，p) 输入x为矩阵, 输出为x的p-条件数，p 的含义同norm r=rcond (x) 输入x为方阵, 输出为x条件数倒数 e=eig(x) 输入x为矩阵，输出x的全部特征值
例： A=[1 2 3;4 5 6;7 8 9]; t=trace(A),d=diag(A) l=tril(A),u=triu(A) b=flipud(A),c=fliplr(A)
MATLAB处理稀疏矩阵: 进行大规模计算的优点
实际应用中， 实际应用中，许多大型矩阵都含有很 元素，称为稀疏矩阵。 多0元素，称为稀疏矩阵。 a=sparse(r,c,v,m,n) 在第 行、第c列 在第r行 列 输入非0元素v，矩阵共m行 列 输出a为稀 输入非0元素 ，矩阵共 行n列，输出 为稀 疏矩阵，只给出(r,c)及v 疏矩阵，只给出 及 a=sparse(aa) 将满元素矩阵 转化为稀 将满元素矩阵aa转化为稀 疏矩阵a 疏矩阵 aa=full(a) 输入稀疏矩阵 ，输出aa为满 输入稀疏矩阵a，输出 为满 矩阵（包含零元素） 矩阵（包含零元素）
例：a=sparse(2,2:3,8,2,4), aa=full(a), a =(2,2) 8 输出 (2,3) 8 aa= 0 0 0 0 0 8 8 0
例：s=sparse([1 6 2],[2 4 6],[0.1 0.2 0.3],6 ,6) A=full(s)
例. 分别用稀疏矩阵和满矩阵求解Ax=b, 比较计算时间 分别用稀疏矩阵和满矩阵求解