搜档网
当前位置:搜档网 › Matlab求解线性方程组、非线性方程组.docx

Matlab求解线性方程组、非线性方程组.docx

Matlab求解线性方程组、非线性方程组.docx
Matlab求解线性方程组、非线性方程组.docx

求解线性方程组solve ,linsolve

例:

A=[5 0 4 2;1 -1 2 1;4 1 2 0;1 1 1 1]; %矩阵的行之间用分号隔开,元素之间用逗号或空格

B=[3;1;1;0]

X=zeros(4,1);% 建立一个4 元列向量

X=linsolve(A,B)

diff (fun , Var, n):对表达式fun中的变量Var求n阶导数。

例如:F=sym('u(x,y)*v(x,y)' ) ; %sym ()用来定义一个符号表达式diff(F); %matlab 区分大小写

pretty(ans) %pretty ():用习惯书写方式显示变量;ans 是答案表达式非线性方程求解

fsolVe(fun,x0,options)

其中fun 为待解方程或方程组的文件名;

x0 位求解方程的初始向量或矩阵;

option 为设置命令参数

建立文件fun.m :

function y=fun(x) y=[x(1)-0.5*sin(x(1))-0.3*cos(x(2)), ...

x(2) -

0.5*cos(x(1))+0.3*sin(x(2))]; >>clear;x0=[0.1,0.1];fsolVe(@fun,x0,optimset('fsolVe'))

注:

...为续行符

m 文件必须以function 为文件头,调用符为@;文件名必须与定义的函数名相同;fsolVe ()主要求解复杂非线性方程和方程组,求解过程是一个逼近过程。

Matlab 求解线性方程组

AX=B 或XA=B

在MATLAB 中,求解线性方程组时,主要采用前面章节介绍的除法运算符“和/ ” “”。如:

X=A?B表示求矩阵方程AX = B的解;

X= B/A表示矩阵方程XA=B的解。

对方程组X = A?B ,要求A和B用相同的行数,X和B有相同的列数,它的行数等于矩阵A 的列数,方程X= B/A 同理。

如果矩阵A不是方阵,其维数是m× n,则有:m = n 恰定方程,求解精确解;m>n 超定方程,寻求最小二乘解;

m

一.恰定方程组

恰定方程组由n 个未知数的n 个方程构成,方程有唯一的一组解,其一般形式可用矩阵,向量写成如下形式:

Ax=b 其中A 是方阵,b 是一个列向量;在线性代数教科书中,最常用的方程组解法有:

(1)利用cramer 公式来求解法;

(2)利用矩阵求逆解法,即x=A-1b ;

(3)利用gaussian 消去法;

(4)利用lu 法求解。

一般来说,对维数不高,条件数不大的矩阵,上面四种解法所得的结果差别不大。前三种解法的真正意义是在其理论上,而不是实际的数值计算。MATLAB 中,出于对算法稳定性的考虑,行列式及逆的计算大都在lu 分解的基础上进行。在MATLAB 中,求解这类方程组的命令十分简单,直接采用表达式:x=A?b 。在MATLAB 的指令解释器在确认变量 A 非奇异后,就对它进行lu 分解,并最终

给出解x;若矩阵A的条件数很大,MATLAB会提醒用户注意所得解的可靠性。如果矩阵A 是奇异的,则Ax=b 的解不存在,或者存在但不唯一;如果矩阵 A 接近奇异时,MATLAB将给出警告信息;如果发现A是奇异的,则计算结果为inf,并且给出警

告信息;如果矩阵A是病态矩阵,也会给出警告信息。注意:在求解方程时,尽量不要用inv(A)*b 命令,而应采用A\b 的解法。因为后者的计算速度比前者快、精度高,尤其当矩阵A的维数比较大时。另外,除法

命令的适用行较强,对于非方阵 A ,也能给出最小二乘解。

二.超定方程组

对于方程组Ax=b, A 为n×m 矩阵,如果 A 列满秩,且n>m 。则方程组没有精确解,此时称方程组为超定方程组。线性超定方程组经常遇到的问题是数据的曲线拟合。对于超定方程,在MATLAB 中,利用左除命令( x=A\b )来寻求它的最小二乘解;还可以用广义逆来求,即x=pinv(A), 所得的解不一定满足Ax=b,x 只是最小二乘意义上的解。左除的方法是建立在奇异值分解基础之上,由此获得的解最可靠;广义逆法是建立在对原超定方程直接进行householder 变换的基础上,其算法可靠性稍逊与奇异值求解,但速度较快;

【例7】

求解超定方程组

A=[2 -1 3;3 1 -5;4 -1 1;1 3 -13]

A=

2 -1 3

3 1 -5

4 -1 1

1 3 -13

b = [3 0 3 -6]';

rank(A) ans=

3

x1=A\b x1=

1.0000

2.0000

1.0000

x2=pinv(A)*b

x2=

1.0000

2.0000

1.0000

A*x1-b

ans=

1.0e-014

-0.0888

-0.0888

-0.1332

可见x1 并不是方程Ax=b 的精确解,用x2=pinv(A)*b 所得的解与x1 相同。三.欠定方程组欠定方程组未知量个数多于方程个数,但理论上有无穷个解。MATLAB 将寻求一个基本解,其中最多只能有m 个非零元素。特解由列主元qr 分解求得。【例8】

解欠定方程组

A = [1-2 1 1;1 -2 1 -1;1 -2 1 5]

A=

1 -

2 1 1

1 -

2 1 -1

1 -

2 1 -1

1 -

2 1 5

b=[1 -1 5] '

x1=A\b

Warning:Rank deficient,rank=2 tol=4.6151e-015 x1= 0

-0.0000

1.0000

x2=pinv(A)*b

x2=

-0.0000

0.0000

1.0000

四.方程组的非负最小二乘解

在某些条件下,所求的线性方程组的解出现负数是没有意义的。虽然方程组可以得到精确解,但却不能取负值解。在这种情况下,其非负最小二乘解比方程的精确解更有意义。在MATLAB 中,求非负最小二乘解常用函数nnls ,其调用格式为:( 1)X=nnls(A,b) 返回方程Ax=b 的最小二乘解,方程的求解过程被限制在x 的条件下;

(2)X=nnls(A,b,TOL)指定误差ToL来求解,ToL的默认值为

TOL=max(size(A))*norm(A,1)*eps ,矩阵的-1 范数越大,求解的误差越大;

(3)[X,W]=nnls(A,b)当x(i)=0 时,w(i)<0 ;当下x(i)>O 时,w(i)O,同时返回一个双向量w。

【例9】求方程组的非负最小二乘解

A=[3.4336 -0.5238 0.6710

-0.5238 3.2833 -0.7302

0.6710 -0.7302 4.0261];

b=[-1.000 1.5000 2.5000];

[X,W]=nnls(A,b)

X= 0

0.6563

0.6998

W=

-3.6820

-0.0000

-0.0000 x1=A\b x1= -0.3569

0.5744

0.7846

A*X-b ans=

1.1258

0.1437 -0.1616

A*x1-b ans= 1.0e-0.15 -0.2220 0.4441

推荐-Broyden方法求解非线性方程组的Matlab实现 精品

Broyden方法求解非线性方程组的Matlab实现 注:matlab代码来自网络,仅供学习参考。 1.把以下代码复制在一个.m文件上 function [sol, it_hist, ierr] = brsola(x,f,tol, parms) % Broyden's Method solver, globally convergent % solver for f(x) = 0, Armijo rule, one vector storage % % This code es with no guarantee or warranty of any kind. % % function [sol, it_hist, ierr] = brsola(x,f,tol,parms) % % inputs: % initial iterate = x % function = f % tol = [atol, rtol] relative/absolute % error tolerances for the nonlinear iteration % parms = [maxit, maxdim] % maxit = maxmium number of nonlinear iterations % default = 40 % maxdim = maximum number of Broyden iterations % before restart, so maxdim-1 vectors are % stored % default = 40 % % output: % sol = solution % it_hist(maxit,3) = scaled l2 norms of nonlinear residuals % for the iteration, number function evaluations, % and number of steplength reductions % ierr = 0 upon successful termination % ierr = 1 if after maxit iterations % the termination criterion is not satsified. % ierr = 2 failure in the line search. The iteration % is terminated if too many steplength reductions % are taken. % % % internal parameter: % debug = turns on/off iteration statistics display as % the iteration progresses

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; end I=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; end end 2.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; end I=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; end end 3.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; end I=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; end end 4.jacobi雅可比迭代法求线性方程组Ax=b的解 function[x,n]=jacobi(A,b,x0,eps,varargin) if nargin==3 eps=1.0e-6; M=200; elseif nargin<3 error return elseif nargin==5 M=varargin{1}; end D=diag(diag(A));%求A的对角矩阵 L=-tril(A,-1);%求A的下三角阵

MatLab求解线性方程组

MatLab解线性方程组一文通 当齐次线性方程AX=0,rank(A)=r

MATLAB解线性方程组的直接方法

在这章中我们要学习线性方程组的直接法,特别是适合用数学软件在计算机上求解的方法. 3.1 方程组的逆矩阵解法及其MATLAB 程序 3.1.3 线性方程组有解的判定条件及其MATLAB 程序 判定线性方程组A n m ?b X =是否有解的MATLAB 程序 function [RA,RB,n]=jiepb(A,b) B=[A b];n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0, disp('请注意:因为RA~=RB ,所以此方程组无解.') return end if RA==RB if RA==n disp('请注意:因为RA=RB=n ,所以此方程组有唯一解.') else disp('请注意:因为RA=RB> A=[2 3 -1 5;3 1 2 -7;4 1 -3 6;1 -2 4 -7]; b=[ 0; 0; 0; 0]; [RA,RB,n]=jiepb(A,b) 运行后输出结果为 请注意:因为RA=RB=n ,所以此方程组有唯一解. RA = 4,RB =4,n =4 在MATLAB 工作窗口输入 >>X=A\b, 运行后输出结果为 X =(0 0 0 0)’. (2) 在MATLAB 工作窗口输入程序 >> A=[3 4 -5 7;2 -3 3 -2;4 11 -13 16;7 -2 1 3];b=[ 0; 0; 0; 0]; [RA,RB,n]=jiepb(A,b)

Matlab求解线性方程组非线性方程组

求解线性方程组 solve,linsolve 例: A=[5 0 4 2;1 -1 2 1;4 1 2 0;1 1 1 1]; %矩阵的行之间用分号隔开,元素之间用逗号或空格 B=[3;1;1;0] X=zeros(4,1);%建立一个4元列向量 X=linsolve(A,B) diff(fun,var,n):对表达式fun中的变量var求n阶导数。 例如:F=sym('u(x,y)*v(x,y)'); %sym()用来定义一个符号表达式 diff(F); %matlab区分大小写 pretty(ans) %pretty():用习惯书写方式显示变量;ans是答案表达式 非线性方程求解 fsolve(fun,x0,options) 为待解方程或方程组的文件名;fun其中 x0位求解方程的初始向量或矩阵; option为设置命令参数 建立文件fun.m: function y=fun(x) y=[x(1)-0.5*sin(x(1))-0.3*cos(x(2)), ... x(2) - 0.5*cos(x(1))+0.3*sin(x(2))]; >>clear;x0=[0.1,0.1];fsolve(@fun,x0,optimset('fsolve')) 注: ...为续行符 m文件必须以function为文件头,调用符为@;文件名必须与定义的函数名相同;fsolve()主要求解复杂非线性方程和方程组,求解过程是一个逼近过程。Matlab求解线性方程组 AX=B或XA=B 在MATLAB中,求解线性方程组时,主要采用前面章节介绍的除法运算符“/”和“\”。如: X=A\B表示求矩阵方程AX=B的解; 的解。XA=B表示矩阵方程B/A=X. 对方程组X=A\B,要求A和B用相同的行数,X和B有相同的列数,它的行数等于矩阵A的列数,方程X=B/A同理。 如果矩阵A不是方阵,其维数是m×n,则有: m=n 恰定方程,求解精确解; m>n 超定方程,寻求最小二乘解; m

matlab程序设计实践-牛顿法解非线性方程

中南大学MATLAB程序设计实践学长有爱奉献,下载填上信息即可上交,没有下载券的自行百度。所需m文件照本文档做即可,即新建(FILE)→脚本(NEW-Sscript)→复制本文档代码→运行(会跳出保存界面,文件名默认不要修改,保存)→结果。第一题需要把数据文本文档和m文件放在一起。全部测试无误,放心使用。本文档针对做牛顿法求非线性函数题目的同学,当然第一题都一样,所有人都可以用。←记得删掉这段话 班级: ? 学号: 姓名:

一、《MATLAB程序设计实践》Matlab基础 表示多晶体材料织构的三维取向分布函数(f=f(φ1,φ,φ2))是一个非常复杂的函数,难以精确的用解析函数表达,通常采用离散 空间函数值来表示取向分布函数,是三维取向分布函数的一个实例。 由于数据量非常大,不便于分析,需要借助图形来分析。请你编写一 个matlab程序画出如下的几种图形来分析其取向分布特征: (1)用Slice函数给出其整体分布特征; " ~ (2)用pcolor或contour函数分别给出(φ2=0, 5, 10, 15, 20, 25, 30, 35 … 90)切面上f分布情况(需要用到subplot函数);

(3) 用plot函数给出沿α取向线(φ1=0~90,φ=45,φ2=0)的f分布情况。 (

备注:数据格式说明 解: (1)( (2)将文件内的数据按照要求读取到矩阵f(phi1,phi,phi2)中,代码如 下: fid=fopen(''); for i=1:18 tline=fgetl(fid); end phi1=1;phi=1;phi2=1;line=0; f=zeros(19,19,19); [ while ~feof(fid) tline=fgetl(fid); data=str2num(tline); line=line+1;数据说明部分,与 作图无关此方向表示f随着 φ1从0,5,10,15, 20 …到90的变化而 变化 此方向表示f随着φ 从0,5,10,15, 20 … 到90的变化而变化 表示以下数据为φ2=0的数据,即f(φ1,φ,0)

利用MATLAB求线性方程组

《MATLAB语言》课成论文 利用MATLAB求线性方程组 姓名:郭亚兰 学号:12010245331 专业:通信工程 班级:2010级通信工程一班 指导老师:汤全武 学院:物电学院 完成日期:2011年12月17日

利用MATLAB求解线性方程组 (郭亚兰 12010245331 2010 级通信一班) 【摘要】在高等数学及线性代数中涉及许多的数值问题,未知数的求解,微积分,不定积分,线性方程组的求解等对其手工求解都是比较复杂,而MATLAB语言正是处理线性方程组的求解的很好工具。线性代数是数学的一个分支,它的研究对象是向量,向量空间(或称线性空间),线性变换和有限维的线性方程组。因而,线性代数被广泛地应用于抽象代数和泛函分析中;由于科学研究中的非线性模型通常可以被近似为线性模型,使得线性代数被广泛地应用于自然科学和社会科学中。线性代数是数学的一个分支,它的研究对象是向量,向量空间(或称线性空间),线性变换和有限维的线性方程组。因而,线性代数被广泛地应用于抽象代数和泛函分析中;由于科学研究中的非线性模型通常可以被近似为线性模型,使得线性代数被广泛地应用于自然科学和社会科学中。线性代数是讨论矩阵理论、与矩阵结合的有限维向量空间及其线性变换理论的一门学科。 【关键字】线性代数MATLAB语言秩矩阵解 一、基本概念 1、N级行列式A:A等于所有取自不同性不同列的n个元素的积的代数和。 2、矩阵B:矩阵的概念是很直观的,可以说是一张表。 3、线性无关:一向量组(a1,a2,…,an)不线性相关,既没有不全为零的数 k1,k2,………kn使得:k1*a1+k2*a2+………+kn*an=0 4、秩:向量组的极在线性无关组所含向量的个数成为这个向量组的秩。 5、矩阵B的秩:行秩,指矩阵的行向量组的秩;列秩类似。记:R(B)

线性方程组求解matlab实现

3.1 方程组的逆矩阵解法及其MATLAB 程序 3.1.3 线性方程组有解的判定条件及其MATLAB 程序 判定线性方程组A n m ?b X =是否有解的MATLAB 程序 function [RA,RB,n]=jiepb(A,b) B=[A b];n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0, disp('请注意:因为RA~=RB ,所以此方程组无解.') return end if RA==RB if RA==n disp('请注意:因为RA=RB=n ,所以此方程组有唯一解.') else disp('请注意:因为RA=RB> A=[2 3 -1 5;3 1 2 -7;4 1 -3 6;1 -2 4 -7]; b=[ 0; 0; 0; 0]; [RA,RB,n]=jiepb(A,b) 运行后输出结果为 请注意:因为RA=RB=n ,所以此方程组有唯一解. RA = 4,RB =4,n =4 在MATLAB 工作窗口输入 >>X=A\b, 运行后输出结果为 X =(0 0 0 0)’. (2) 在MATLAB 工作窗口输入程序 >> A=[3 4 -5 7;2 -3 3 -2;4 11 -13 16;7 -2 1 3];b=[ 0; 0; 0; 0]; [RA,RB,n]=jiepb(A,b) 运行后输出结果 请注意:因为RA=RB> A=[4 2 -1;3 -1 2;11 3 0]; b=[2;10;8]; [RA,RB,n]=jiepb(A,B) 运行后输出结果 请注意:因为RA~=RB ,所以此方程组无解. RA =2,RB =3,n =3 (4)在MATLAB 工作窗口输入程序

遗传算法解非线性方程组的Matlab程序

遗传算法解非线性方程组的Matlab程序 程序用MATLAB语言编写。之所以选择MATLB,是因为它简单,但又功能强大。写1行MATLAB程序,相当于写10行C++程序。在编写算法阶段,最好用MATLAB语言,算法验证以后,要进入工程阶段,再把它翻译成C++语言。 本程序的算法很简单,只具有示意性,不能用于实战。 非线性方程组的实例在函数(2)nonLinearSumError1(x)中,你可以用这个实例做样子构造你自己待解的非线性方程组。 %注意:标准遗传算法的一个重要概念是,染色体是可能解的2进制顺序号,由这个序号在可能解的集合(解空间)中找到可能解 %程序的流程如下: %程序初始化,随机生成一组可能解(第一批染色体) %1: 由可能解的序号寻找解本身(关键步骤) %2:把解代入非线性方程计算误差,如果误差符合要求,停止计算 %3:选择最好解对应的最优染色体 %4:保留每次迭代产生的最好的染色体,以防最好染色体丢失 %5: 把保留的最好的染色体holdBestChromosome加入到染色体群中 %6: 为每一条染色体(即可能解的序号)定义一个概率(关键步骤) %7:按照概率筛选染色体(关键步骤) %8:染色体杂交(关键步骤) %9:变异 %10:到1 %这是遗传算法的主程序,它需要调用的函数如下。 %由染色体(可能解的2进制)顺序号找到可能解: %(1)x=chromosome_x(fatherChromosomeGroup,oneDimensionSet,solutionSum); %把解代入非线性方程组计算误差函数:(2)functionError=nonLinearSumError1(x); %判定程是否得解函数:(3)[solution,isTrue]=isSolution(x,funtionError,solutionSumError); %选择最优染色体函数: %(4)[bestChromosome,leastFunctionError]=best_worstChromosome(fatherChromosomeGroup,functionError); %误差比较函数:从两个染色体中,选出误差较小的染色体 %(5)[holdBestChromosome,holdLeastFunctionError]... % =compareBestChromosome(holdBestChromosome,holdLeastFunctionError,... % bestChromosome,leastFuntionError) %为染色体定义概率函数,好的染色体概率高,坏染色体概率低 %(6)p=chromosomeProbability(functionError); %按概率选择染色体函数: %(7)slecteChromosomeGroup=selecteChromome(fatherChromosomeGroup,p); %父代染色体杂交产生子代染色体函数 %(8)sonChrmosomeGroup=crossChromosome(slecteChromosomeGroup,2); %防止染色体超出解空间的函数 %(9)chromosomeGroup=checkSequence(chromosomeGroup,solutionSum) %变异函数 %(10)fatherChromosomeGroup=varianceCh(sonChromosomeGroup,0.8,solutionN); %通过实验有如下结果: %1。染色体应当多一些

实验一用matlab求解线性方程组

实验1.1 用matlab 求解线性方程组 第一节 线性方程组的求解 一、齐次方程组的求解 rref (A ) %将矩阵A 化为阶梯形的最简式 null (A ) %求满足AX =0的解空间的一组基,即齐次线性方程组的基 础解系 【例1】 求下列齐次线性方程组的一个基础解系,并写出通解: 我们可以通过两种方法来解: 解法1: >> A=[1 -1 1 -1;1 -1 -1 1;1 -1 -2 2]; >> rref(A) 执行后可得结果: ans= 1 -1 0 0 0 0 -1 1 0 0 0 0 由最简行阶梯型矩阵,得化简后的方程 ??? ??=+--=+--=-+-0 22004321 43214321x x x x x x x x x x x x

取x2,x4为自由未知量,扩充方程组为 即 提取自由未知量系数形成的列向量为基础解系,记 所以齐次方程组的通解为 解法2: clear A=[1 -1 1 -1;1 -1 -1 1;1 -1 -2 2]; B=null(A, 'r') % help null 看看加个‘r’是什么作用, 若去掉r ,是什么结果? 执行后可得结果: B= 1 0 1 0 0 1 0 1 ?? ?=-=-0 04321x x x x ?????? ?====4 4432221x x x x x x x x ??? ??? ??????+????????????=????? ???????1100001142 4321x x x x x x , 00111????? ? ??????=ε, 11002????? ???????=ε2 211εεk k x +=

Matlab线性方程组求解(Gauss消去法)

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-1 for i=k+1:n if a(k,k)==0 return end m=a(i,k)/a(k,k); for j=k+1:n a(i,j)=a(i,j)-m*a(k,j); end b(i)=b(i)-m*b(k); end det=det*a(k,k); %计算行列式 end det=det*a(n,n); for k=n:-1:1 %回代求解 for j=k+1:n b(k)=b(k)-a(k,j)*x(j); end x(k)=b(k)/a(k,k);

end Example: >> 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.0047 1.0010 2. 列主元Gauss消去法: function x=detGauss(a,b) % Gauss列主元消去法 [n,m]=size(a); nb=length(b); det=1; %存储行列式值 x=zeros(n,1); for k=1:n-1 amax=0; %选主元 for i=k:n if abs(a(i,k))>amax amax=abs(a(i,k));r=i; end end if amax<1e-10 return; end if r>k %交换两行 for j=k:n

MATLAB应用 求解非线性方程

第7章 求解非线性方程 7.1 多项式运算在MATLAB 中的实现 一、多项式的表达 n 次多项式表达为:n a +??++=x a x a x a p(x )1-n 1-n 1n 0,是n+1项之和 在MATLAB 中,n 次多项式可以用n 次多项式系数构成的长度为n+1的行向量表示 [a0, a1,……an-1,an] 二、多项式的加减运算 设 有 两 个 多 项 式 n a +??++=x a x a x a p1(x )1-n 1-n 1n 0和 m b +??++=x b x b x b p2(x )1-m 1-m 1m 0。它们的加减运算实际上就是它们的对应系 数的加减运算。当它们的次数相同时,可以直接对多项式的系数向量进行加减运算。当它们的次数不同时,应该把次数低的多项式无高次项部分用0系数表示。 例2 计算()()1635223-+++-x x x x a=[1, -2, 5, 3]; b=[0, 0, 6, -1]; c=a+b 例3 设()6572532345++-+-=x x x x x x f ,()3532-+=x x x g ,求f(x)+g(x) f=[3, -5, 2, -7, 5, 6]; g=[3, 5, -3]; g1=[0, 0, 0, g];%为了和f 的次数找齐 f+g1, f-g1 三、多项式的乘法运算 conv(p1,p2) 例4 在上例中,求f(x)*g(x) f=[3, -5, 2, -7, 5, 6]; g=[3, 5, -3]; conv(f, g) 四、多项式的除法运算 [Q, r]=deconv(p1, p2) 表示p1除以p2,给出商式Q(x),余式r(x)。Q,和r 仍为多项式系数向量 例4 在上例中,求f(x)/g(x) f=[3, -5, 2, -7, 5, 6]; g=[3, 5, -3]; [Q, r]=deconv(f, g) 五、多项式的导函数 p=polyder(P):求多项式P 的导函数 p=polyder(P,Q):求P·Q 的导函数

3-7变量非线性方程组的逆Broyden解法matlab程序

3-7变量非线性方程组的逆Broyden解法 matlab程序 function n_broyden clear all %清内存 clc %清屏 format long i=input('请输入非线性方程组个数(3-7)i='); switch i case 3 syms x y z x0=input('请输入初值(三维列向量[x;y;z])='); eps=input('请输入允许的误差精度eps='); f1=input('请输入f1(x,y,z)='); f2=input('请输入f2(x,y,z)='); f3=input('请输入f3(x,y,z)='); F=[f1;f2;f3]; J=jacobian(F,[x,y,z]); %求jacobi矩阵 Fx0=subs(F,{x,y,z},x0); Jx0=subs(J,{x,y,z},x0); H=inv(Jx0); x1=x0-H*Fx0; k=0; while norm(x1-x0)>eps %用2范数来做循环条件 p=x1-x0; q=subs(F,{x,y,z},x1)-subs(F,{x,y,z},x0); H=H-((H*q-p)*p'*H)*(p'*H*q)^-1; x0=x1; Fx0=subs(F,{x,y,z},x0); x1=x1-H*Fx0; k=k+1; end x1 k case 4 syms a b c d x0=input('请输入初值(列向量[a;b;c;d])=');

eps=input('请输入允许的误差精度eps='); f1=input('请输入f1(a,b,c,d)='); f2=input('请输入f2(a,b,c,d)='); f3=input('请输入f3(a,b,c,d)='); f4=input('请输入f4(a,b,c,d)='); F=[f1;f2;f3;f4]; J=jacobian(F,[a,b,c,d]); %求jacobi矩阵 Fx0=subs(F,{a,b,c,d},x0); Jx0=subs(J,{a,b,c,d},x0); H=inv(Jx0); x1=x0-H*Fx0; k=0; while norm(x1-x0)>eps %用2范数来做循环条件 p=x1-x0; q=subs(F,{a,b,c,d},x1)-subs(F,{a,b,c,d},x0); H=H-((H*q-p)*p'*H)*(p'*H*q)^-1; x0=x1; Fx0=subs(F,{a,b,c,d},x0); x1=x1-H*Fx0; k=k+1; end x1 k case 5 syms a b c d e x0=input('请输入初值(列向量[a;b;c;d;e])='); eps=input('请输入允许的误差精度eps='); f1=input('请输入f1(a,b,c,d,e)='); f2=input('请输入f2(a,b,c,d,e)='); f3=input('请输入f3(a,b,c,d,e)='); f4=input('请输入f4(a,b,c,d,e)='); f5=input('请输入f5(a,b,c,d,e)='); F=[f1;f2;f3;f4;f5]; J=jacobian(F,[a,b,c,d,e]); %求jacobi矩阵 Fx0=subs(F,{a,b,c,d,e},x0); Jx0=subs(J,{a,b,c,d,e},x0); H=inv(Jx0); x1=x0-H*Fx0;

matlab程序设计实践-牛顿法解非线性方程

中南大学 MATLAB程序设计实践学长有爱奉献,下载填上信息即可上交,没有下载券 的自行百度。所需m文件照本文档做即可,即新建(FILE)→脚本(NEW-Sscript)→复制本文档代码→运行(会跳出 保存界面,文件名默认不要修改,保存)→结果。第 一题需要把数据文本文档和m文件放在一起。全部测 试无误,放心使用。本文档针对做牛顿法求非线性函 数题目的同学,当然第一题都一样,所有人都可以用。 ←记得删掉这段话 班级: 学号: 姓名: 一、《MATLAB程序设计实践》Matlab基础

表示多晶体材料织构的三维取向分布函数(f=f(φ1,φ,φ2))是一个非常复杂的函数,难以精确的用解析函数表达,通常采用离散空间函数值来表示取向分布函数,是三维取向分布函数的一个实例。由于数据量非常大,不便于分析,需要借助图形来分析。请你编写一个matlab程序画出如下的几种图形来分析其取向分布特征:(1)用Slice函数给出其整体分布特征; (2)用pcolor或contour函数分别给出(φ2=0, 5, 10, 15, 20, 25, 30, 35 … 90)切面上f分布情况(需要用到subplot函数);

(3) 用plot函数给出沿α取向线(φ1=0~90,φ=45,φ2=0)的f分布情况。

备注:数据格式说明 解: (1)将文件内的数据按照要求读取到矩阵f(phi1,phi,phi2)中,代码如下: fid=fopen(''); for i=1:18 tline=fgetl(fid); end phi1=1;phi=1;phi2=1;line=0; f=zeros(19,19,19); while ~feof(fid) tline=fgetl(fid); data=str2num(tline); line=line+1; if mod(line,20)==1 phi2=(data/5)+1; phi=1; 数据说明部分,与作图无关 此方向表示f 随着φ1从0,5,10,15, 20 …到90的变化而变化 此方向表示f 随着φ从0,5,10,15, 20 …到90的变化而变化 表示以下数据为φ2=0的数据,即f (φ1,φ,0)

Matlab求解线性方程组、非线性方程组

Matlab求解线性方程组、非线性方程组 姓名:罗宝晶学号:15 专业:材料学院高分子系 第一部分数值计算 Matlab求解线性方程组AX=B或XA=B 在MATLAB中,求解线性方程组时,主要采用除法运算符“/”和“\”。如:X=A\B表示求矩阵方程AX=B的解; X=B/A表示矩阵方程XA=B的解。 对方程组X=A\B,要求A和B用相同的行数,X和B有相同的列数,它的行数等于矩阵A的列数,方程X=B/A同理。 如果矩阵A不是方阵,其维数是m×n,则有: m=n 恰定方程,求解精确解; m>n 超定方程,寻求最小二乘解; mm。则方程组没有精确解,此时称方程组为超定方程组。线性超定方程组经常遇到的问题是数据的曲线拟合。对于超定方程,在MATLAB中,利用左除命令(x=A\b)来寻求它的最小二乘解;还可以用广义逆来求,即x=pinv(A),所得的解不一定满足Ax=b,x只是最小二乘意义上的解。左除的方法是建立在奇异值分解基础之上,由此获得的解最可靠;广义逆法是建立在对原超定方程直接进行householder变换的基础上,其算法可靠性稍逊与奇异值求解,但速度较快;

线性方程组求解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-1 for i=k+1:n if a(k,k)==0 return end m=a(i,k)/a(k,k); for j=k+1:n a(i,j)=a(i,j)-m*a(k,j); end b(i)=b(i)-m*b(k); end det=det*a(k,k); end

det=det*a(n,n); for k=n:-1:1 %回代 for j=k+1:n b(k)=b(k)-a(k,j)*x(j); end x(k)=b(k)/a(k,k); end Example: >> 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.0047 1.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-1 amax=0;% 选主元 for i=k:n if abs(a(i,k))>amax amax=abs(a(i,k));r=i; end end if amax<1e-10 return; end if r>k %交换两行 for j=k:n z=a(k,j);a(k,j)=a(r,j);a(r,j)=z; end z=b(k);b(k)=b(r);b(r)=z;det=-det; end

基于Matlab的牛顿迭代法解非线性方程组

基于Matlab 实现牛顿迭代法解非线性方程组 已知非线性方程组如下 2211221212 10801080x x x x x x x ?-++=??+-+=?? 给定初值0(0,0)T x =,要求求解精度达到0.00001 首先建立函数F(x),方程组编程如下,将F.m 保存到工作路径中: function f=F(x) f(1)=x(1)^2-10*x(1)+x(2)^2+8; f(2)=x(1)*x(2)^2+x(1)-10*x(2)+8; f=[f(1) f(2)]; 建立函数DF(x),用于求方程组的Jacobi 矩阵,将DF.m 保存到工作路径中: function df=DF(x) df=[2*x(1)-10,2*x(2);x(2)^2+1,2*x(1)*x(2)-10]; 编程牛顿迭代法解非线性方程组,将newton.m 保存到工作路径中: clear; clc x=[0,0]'; f=F(x); df=DF(x); fprintf('%d %.7f %.7f\n',0,x(1),x(2)); N=4; for i=1:N y=df\f'; x=x-y; f=F(x); df=DF(x); fprintf('%d %.7f %.7f\n',i,x(1),x(2)); if norm(y)<0.0000001 break ; else end end

运行结果如下: 0 0.0000000 0.0000000 1 0.8000000 0.8800000 2 0.9917872 0.9917117 3 0.9999752 0.9999685 4 1.0000000 1.0000000

MATLAB 非线性方程(组)求根

实用数值方法(Matlab) 综述报告题目:非线性方程(组)求根问题 小组成员

许多数学和物理问题归结为解函数方程f(x)=0。方程f(x)=0的解称为方程的根。对于非线性方程,在某个范围内往往不止一个根,而且根的分布情况可能很复杂,面对这种情况,通常先将考察的范围花费为若干个子段,然后判断哪些子段内有根,然后再在有根子段内找出满足精度要求的近似根。为此适当选取有根子段内某一点作为根的初始值近似,然后运用迭代方法使之足部精确化。这就是方程求根的迭代法。下面介绍书上的几种方法: 1、二分法 (1)方法概要: 假定函数f(x)在[a,b]上连续,且f(a)f(b)=0,则方程f(x)=0在[a,b]内一定有实根。取其中 将其二分,判断所求的根在的左侧还是右侧,得到一个新的有根区间 点 [],长度为[a,b]的一半。对新的有根区间继续实行上述二分手段,直至二分k次后有根区间[]长度 可见,如果二分过程无限继续下去,这些有限根区间最终必收敛于一点,该点就是所求的根。在实际计算过程中不可能完成这个无限过程,允许有一定的误差,则二分k+1次后 只要有根区间[]的长度小于,那么结果关于允许误差就能“准确”地满足方程f(x)=0。 (2)计算框图:

2、开方法 对于给定,求开方值 为此,可以运用校正技术设计从预报值生成校正值的迭代公式。自然希望校正值 能更好满足所给方程: 这是个关于校正量的近似关系式,如果从中删去二次项,即可化归为一次方程 解之有 从而关于校正值有如下开方公式 上述演绎过程表明,开方法的设计思想是逐步线性化,即将二次方程的求解画归为一次方程求解过程的重复。开方公式规定了预报值与校正值之间的一种函数关系,这里 为开方法的迭代函数。 3、Newton法 (1)方法概要

Matlab求解线性方程组、非线性方程组.docx

求解线性方程组solve ,linsolve 例: A=[5 0 4 2;1 -1 2 1;4 1 2 0;1 1 1 1]; %矩阵的行之间用分号隔开,元素之间用逗号或空格 B=[3;1;1;0] X=zeros(4,1);% 建立一个4 元列向量 X=linsolve(A,B) diff (fun , Var, n):对表达式fun中的变量Var求n阶导数。 例如:F=sym('u(x,y)*v(x,y)' ) ; %sym ()用来定义一个符号表达式diff(F); %matlab 区分大小写 pretty(ans) %pretty ():用习惯书写方式显示变量;ans 是答案表达式非线性方程求解 fsolVe(fun,x0,options) 其中fun 为待解方程或方程组的文件名; x0 位求解方程的初始向量或矩阵; option 为设置命令参数 建立文件fun.m : function y=fun(x) y=[x(1)-0.5*sin(x(1))-0.3*cos(x(2)), ... x(2) - 0.5*cos(x(1))+0.3*sin(x(2))]; >>clear;x0=[0.1,0.1];fsolVe(@fun,x0,optimset('fsolVe')) 注: ...为续行符 m 文件必须以function 为文件头,调用符为@;文件名必须与定义的函数名相同;fsolVe ()主要求解复杂非线性方程和方程组,求解过程是一个逼近过程。 Matlab 求解线性方程组 AX=B 或XA=B

在MATLAB 中,求解线性方程组时,主要采用前面章节介绍的除法运算符“和/ ” “”。如: X=A?B表示求矩阵方程AX = B的解; X= B/A表示矩阵方程XA=B的解。 对方程组X = A?B ,要求A和B用相同的行数,X和B有相同的列数,它的行数等于矩阵A 的列数,方程X= B/A 同理。 如果矩阵A不是方阵,其维数是m× n,则有:m = n 恰定方程,求解精确解;m>n 超定方程,寻求最小二乘解; m

相关主题