下载文档原格式
非线性方程组求解-Matlab-fsolve实例一：①建立文件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（）主要求解复杂非线性方程和方程组，求解过程是一个逼近过程。
实例二：①建立文件fun.mfunction F=myfun(x)F=[x(1)-3*x(2)-sin(x(1));2*x(1)+x(2)-cos(x(2))];②然后在命令窗口求解：>> x0=[0;0]; %设定求解初值>> options=optimset('Display','iter'); %设定优化条件>> [x,fv]=fsolve(@myfun,x0,options) %优化求解%MATLAB显示的优化过程Norm of First-order Trust-region Iteration Func-count f(x) step optimality radius0 3 1 2 11 6 0.000423308 0.5 0.0617 12 9 5.17424e-010 0.00751433 4.55e-005 1.253 12 9.99174e-022 1.15212e-005 9.46e-011 1.25 Optimization terminated: first-order optimality is less than options.TolFun.x =0.49660.0067fv =1.0e-010 *0.31610.0018实例三：求下列非线性方程组在(0.5,0.5) 附近的数值解。
⾮线性⽅程求解基于MATLAB的⾮线性⽅程的五种解法探讨摘要：本⽂利⽤matlab软件对⾮线性⽅程解法中的⼆分法、简单迭代法、⽜顿法、割线法以及Steffensen法的数值分析⽅法的算法原理及实现⽅法进⾏了探讨。
对f x x x=+-()2ln2的零点问题，分别运⽤以上五种不同的⽅法进⾏数值实验，⽐较⼏种解法的优缺点并进⾏初步分析评价。
关键词：⼆分法、简单迭代法、⽜顿法、割线法、Steffensen法1、引⾔在很多实际问题中，经常需要求⾮线性⽅程f(x) =0的根。
⽅程f(x) =0的根叫做函数f(x)的零点。
由连续函数的特性知：若f(x)在闭区间[a,b ]上连续,且()()0f a f b<.则f(x) =0在开区间(a,b)内⾄少有⼀个实根。
这时称[a,b]为⽅程f(x) =0的根的存在区间。
本⽂主要对⾮线性⽅程的数值解法进⾏分析，并介绍了⾮线性⽅程数值解法的五种⽅法。
并设=+-.f x x x()2ln2f x在[1,2]上的图形，如图1：. 显然，函数在[1,2]之间有⼀个零点。
⾸先画出()2、计算机配置操作系统Windows 7 旗舰版内存2GB处理器AMD 4核 A6-3400M APU 1.4GHz图.13、⼆分法⼆分法的基本思想是将⽅程根的区间平分为两个⼩区间,把有根的⼩区间再平分为两个更⼩的区间,进⼀步考察根在哪个更⼩的区间内。
如此继续下去,直到求出满⾜精度要求的近似值。
设函数()f x 在区间[a,b ]上连续,且f(a)·f(b) <0，则[a,b ]是⽅程f(x) =0的根的存在区间，设其内有⼀实根，记为x*。
取区间[a,b ]的中点()2k a b x +=并计算1()f x ,则必有下列三种情况之⼀成⽴： (1) 1()f x =0，x1就是⽅程的根x*;(2)()f a .1()f x <0,⽅程的根x*位于区间[a, 1x ]之中,此时令111,a a b x ==; (3)1()f x .()f b <0,⽅程的根x*位于区间[1x ,b ]之中,此时令11a x =,1b b =。
第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 xa=[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 的导函数[p,q]=polyder(P ,Q)：求P/Q 的导函数，导函数的分子存入p ，分母存入q 。
非线性方程组求解的牛顿迭代法用MATLAB实现首先，我们需要定义非线性方程组。
假设我们要求解方程组：```f1(x1,x2)=0f2(x1,x2)=0```其中，`x1`和`x2`是未知数，`f1`和`f2`是非线性函数。
我们可以将这个方程组表示为向量的形式：```F(x)=[f1(x1,x2);f2(x1,x2)]=[0;0]```其中，`F(x)`是一个列向量。
为了实现牛顿迭代法，我们需要计算方程组的雅可比矩阵。
雅可比矩阵是由方程组的偏导数组成的矩阵。
对于方程组中的每个函数，我们可以计算其对每个变量的偏导数，然后将这些偏导数组成一个矩阵。
在MATLAB中，我们可以使用`jacobi`函数来计算雅可比矩阵。
以下是一个示例函数的定义：```matlabfunction J = jacobi(x)x1=x(1);x2=x(2);J = [df1_dx1, df1_dx2; df2_dx1, df2_dx2];end```其中，`x`是一个包含未知数的向量，`df1_dx1`和`df1_dx2`是`f1`对`x1`和`x2`的偏导数，`df2_dx1`和`df2_dx2`是`f2`对`x1`和`x2`的偏导数。
下一步是实现牛顿迭代法。
牛顿迭代法的迭代公式为：```x(k+1)=x(k)-J(x(k))\F(x(k))```其中，`x(k)`是第`k`次迭代的近似解，`\`表示矩阵的求逆操作。
在MATLAB中，我们可以使用如下代码来实现牛顿迭代法：```matlabfunction x = newton_method(x_initial)max_iter = 100; % 最大迭代次数tol = 1e-6; % 收敛阈值x = x_initial; % 初始解for k = 1:max_iterF=[f1(x(1),x(2));f2(x(1),x(2))];%计算F(x)J = jacobi(x); % 计算雅可比矩阵 J(x)delta_x = J \ -F; % 计算增量 delta_xx = x + delta_x; % 更新 xif norm(delta_x) < tolbreak; % 达到收敛条件，停止迭代endendend```其中，`x_initial`是初始解的向量，`max_iter`是最大迭代次数，`tol`是收敛阈值。
matlab fsolve解方程组在数学和工程领域中，解方程组是一项重要的任务。
而在MATLAB 软件中，一个强大的求解方程组的函数是fsolve。
本文将详细介绍MATLAB中fsolve函数的用法和示例，以帮助读者更好地理解如何使用该函数来解决方程组问题。
一、fsolve函数的基本概念1. fsolve函数是MATLAB中用于求解非线性方程组的函数。
它可以通过数值方法找出一组方程组的近似解。
fsolve函数的基本语法如下：x = fsolve(fun,x0,options)其中，fun为一个函数句柄，表示待求解方程组的函数；x0为方程组的初值；options为可选参数，用于指定fsolve函数的求解选项。
2. 待求解方程的形式应为fun(x)=0，其中x为方程组的未知量。
在MATLAB中，方程组可以表示为一个函数句柄fun，该句柄接受一个向量x作为输入参数，并返回一个向量f(x)。
方程组的解即为fun(x)=0的解。
3. 在使用fsolve函数之前，需要先定义一个函数句柄fun，表示待求解方程组。
函数句柄可以通过匿名函数或函数文件的形式定义。
例如，可以使用匿名函数定义一个函数句柄fun：fun = @(x) [x(1)^2 + x(2)^2 - 1; x(1) - x(2)^2];二、fsolve函数的使用步骤1. 定义待求解方程组的函数句柄fun。
根据待求解方程组的数学模型，使用匿名函数或函数文件的方式定义一个函数句柄。
例如： fun = @(x) [x(1)^2 + x(2)^2 - 1; x(1) - x(2)^2];2. 定义待求解方程组的初值x0。
初值应与方程组的解接近，并且尽量选择不容易导致计算失败的初值。
3. 定义fsolve函数的求解选项options（可选）。
options参数用于控制fsolve函数的计算选项，例如迭代次数、容差等。
可以使用optimset函数创建一个options结构体变量，并通过修改该结构体的属性来控制求解选项。
matlab求解非线性方程组及极值默认分类2010-05-18 15:46:13 阅读1012 评论2 字号：大中小订阅一、概述：求函数零点和极值点：Matlab中三种表示函数的方法: 1. 定义一个m函数文件, 2.使用函数句柄; 3.定义inline函数, 其中第一个要掌握简单函数编写, 二, 三中掌握一个。
函数的'常规'使用有了函数了, 我们怎么用呢, 一种是直接利用函数来计算, 例如: sin(pi), 还有我们提到的mysqr(3)...另一种是函数画图, 例如Plottools中提到的ezplot, ezsurf... 但是这也太小儿科了, 有没有想过定义函数后, 利用它来: 求解零点(即解f(x)=0方程), 最优化(求最值/极值点), 求定积分, 常微分方程求解等. 当然这里由于篇幅有限(空间快满了)以及这个只是'基础教程'的缘故, 只提及一些皮毛知识, 掌握这些后, 如果需要你可以进一步学习.解f(x)=0已知函数求解函数值=0所表示的方程, Matlab中有两个函数可以做到, fzero和fsolve前者只能解一元方程, 后者可以解多元方程组, 不过基本使用形式上差不多:解=fzero(函数, 初值, options)解=fsolve(函数, 初值, options)关于解: fzero给出的是x单值的解, fsolve给出的是解x可能处于的区间, 当然, 这个区间很窄.关于'函数', 还记得前面提到的三种表示方法吧, 在这里都可以用, 记住就是: 如果直接使用函数名, 要用单引号将它括起来, 而函数句柄, inline函数可以直接使用.关于'初值': 电脑比较笨, 它寻找解的办法是尝试不同地x值, 摸索解在哪里, 所以我们一开始就要给它指明从哪里开始下手, 初值这里, 可以只给它一个值, 让它在这个值附近找解, 也可以给它一个区间(区间用[下限,上限]这种方式表示), 它会在这个区间内找解.fzero的一些局限, 如果你给定的初值是区间, 而恰好函数在区间端点处同号, fzero会出错, 而如果你只给一个初值, fezro又有可能'走错方向', 例如给初值2让它解mysqr这个函数方程就出错了, FT!寻找函数极值/最值Matlab中也有两个函数可以做到, 是: fminbnd: 寻找一元函数极小值; fminsearch: 寻找多元函数极小值(当然一元也行). 别问我怎么没有找极大值的Matlab函数, 你把原函数取负数, 寻找它的极小值不就行了. 相关语法:x=fminbnd(函数, 区间起始值, 区间终止值)x=fminsearch(函数, 自变量初值)相关说明: fminbnd中指定要查找极小值的自变量区间, 好像不指定也行, 不过那样的话, 如果函数有多个极小值就可能比较难以预料结果了.fminsearch中要给定一个初值, 这个初值可以是自变量向量(将自变量依次排在一起组成向量)的初值, 也可以是表示向量初值区间的一个矩阵.函数: 那三种形式都适用, 但是记住, 直接使用函数名称需要加单引号!cite from:/qq529312840/blog/item/3687e4c7e7e2d6d9d0006049.html二、实例+讲解（1）非线性方程数值求解：1 单变量非线性方程求解在MATLAB中提供了一个fzero函数，可以用来求单变量非线性方程的根。
matlab-fsolve函数求解多元⾮线性⽅程记录⼀下代码，⽅便下次套⽤模板options=optimset('MaxFunEvals',1e4,'MaxIter',1e4);[x,fval,exitflag] = fsolve(@(x) myfun1(x),[75;1.5],options)function f = myfun1(x)f=tan(x(1)*pi/180) - ( ( 1025*9.8*pi*x(2)/4-980 )/(0.625*4*(2-x(2))*24*24) );%有两个未知数x(1)和x(2),从参数⾥传进来endoptions理解成设定要求，精度范围，没有则默认，是多少问题不⼤。
[75;1.5]是x(1)和x(2)的初值，如果是同⼀个数不同初值则是[ 70 1;75 1.5 ]，在初值附近找最优解。
理解成：或许有多个最优解，如果初值不⼀样，最优解也不⼀样。
⾮线性⼏乎都是近似解。
⾄于初值怎么设置，结合问题分析，⽐如杆⼦靠墙的倾斜⾓度⼤约在60度以上，⽽不是⼗⼏⼆⼗度。
函数myfun1的求解情况是f=0。
fval表⽰误差，越⼩越好。
exitflag表⽰迭代退出条件，为1的时候最理想。
1 fsolve converged to a root.2 Change in X too small.3 Change in residual norm too small.4 Computed search direction too small.0 Too many function evaluations or iterations.-1 Stopped by output/plot function.-2 Converged to a point that is not a root.-3 Trust region radius too small (Trust-region-dogleg).最终求出来两个值，分别表⽰两个未知数x(1)和x(2)。
非线性方程组求解1.mulStablePoint用不动点迭代法求非线性方程组的一个根function [r,n]=mulStablePoint(F,x0,eps)%非线性方程组：f%初始解：a%解的精度：eps%求得的一组解：r%迭代步数：nif nargin==2eps=1.0e-6;endx0 = transpose(x0);n=1;tol=1;while tol>epsr= subs(F,findsym(F),x0); %迭代公式tol=norm(r-x0); %注意矩阵的误差求法，norm为矩阵的欧几里德范数n=n+1;x0=r;if(n>100000) %迭代步数控制disp('迭代步数太多，可能不收敛！');return;endend2.mulNewton用牛顿法法求非线性方程组的一个根function [r,n]=mulNewton(F,x0,eps)if nargin==2eps=1.0e-4;endx0 = transpose(x0);Fx = subs(F,findsym(F),x0);var = findsym(F);dF = Jacobian(F,var);dFx = subs(dF,findsym(dF),x0);r=x0-inv(dFx)*Fx;n=1;tol=1;while tol>epsx0=r;Fx = subs(F,findsym(F),x0);dFx = subs(dF,findsym(dF),x0);r=x0-inv(dFx)*Fx; %核心迭代公式tol=norm(r-x0);n=n+1;if(n>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endend3.mulDiscNewton用离散牛顿法法求非线性方程组的一个根function [r,m]=mulDiscNewton(F,x0,h,eps)format long;if nargin==3eps=1.0e-8;endn = length(x0);fx = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-fx)/h(i);endr=transpose(x0)-inv(J)*fx;m=1;tol=1;while tol>epsxs=r;fx = subs(F,findsym(F),xs);J = zeros(n,n);for i=1:nx1 = xs;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-fx)/h(i);endr=xs-inv(J)*fx; %核心迭代公式tol=norm(r-xs);m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endendformat short;4.mulMix用牛顿-雅可比迭代法求非线性方程组的一个根function [r,m]=mulMix(F,x0,h,l,eps)if nargin==4eps=1.0e-4;endn = length(x0);J = zeros(n,n);Fx = subs(F,findsym(F),x0);for i=1:nx1 = x0;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-Fx)/h(i);endD = diag(diag(J));C =D - J;inD = inv(D);H = inD*C;Hm = eye(n,n);for i=1:l-1Hm = Hm + power(H,i);enddr = Hm*inD*Fx;r = transpose(x0)-dr; m=1;tol=1;while tol>epsx0=r;Fx = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-Fx)/h(i);endD = diag(diag(J));C =D - J;inD = inv(D);H = inD*C;Hm = eye(n,n);for i=1:l-1Hm = Hm + power(H,i);enddr = Hm*inD*Fx;r = x0-dr; %核心迭代公式tol=norm(r-x0);m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endend5.mulNewtonSOR用牛顿-SOR迭代法求非线性方程组的一个根function [r,m]=mulNewtonSOR(F,x0,w,h,l,eps)if nargin==5eps=1.0e-4;endn = length(x0);J = zeros(n,n);Fx = subs(F,findsym(F),x0);for i=1:nx1 = x0;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-Fx)/h(i);endD = diag(diag(J));L = -tril(J-D);U = -triu(J-D);inD = inv(D-w*L);H = inD*(D - w*D+w*L);;Hm = eye(n,n);for i=1:l-1Hm = Hm + power(H,i);enddr = w*Hm*inD*Fx;r = transpose(x0)-dr;m=1;tol=1;while tol>epsx0=r;Fx = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-Fx)/h(i);endD = diag(diag(J));L = -tril(J-D);U = -triu(J-D);inD = inv(D-w*L);H = inD*(D - w*D+w*L);;Hm = eye(n,n);for i=1:l-1Hm = Hm + power(H,i);enddr = w*Hm*inD*Fx;r = x0-dr; %核心迭代公式tol=norm(r-x0);m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endend6.mulDNewton用牛顿下山法求非线性方程组的一个根function [r,m]=mulDNewton(F,x0,eps)%非线性方程组：F%初始解：x0%解的精度：eps%求得的一组解：r%迭代步数：nif nargin==2eps=1.0e-4;endx0 = transpose(x0);dF = Jacobian(F);m=1;tol=1;while tol>epsttol=1;w=1;Fx = subs(F,findsym(F),x0);dFx = subs(dF,findsym(dF),x0);F1=norm(Fx);while ttol>=0 %下面的循环是选取下山因子w的过程r=x0-w*inv(dFx)*Fx; %核心的迭代公式Fr = subs(F,findsym(F),r);ttol=norm(Fr)-F1;w=w/2;endtol=norm(r-x0);m=m+1;x0=r;if(m>100000) %迭代步数控制disp('迭代步数太多，可能不收敛！');return;endend7.mulGXF1用两点割线法的第一种形式求非线性方程组的一个根function [r,m]=mulGXF1(F,x0,x1,eps)format long;if nargin==3eps=1.0e-4;endx0 = transpose(x0);x1 = transpose(x1);n = length(x0);fx = subs(F,findsym(F),x0);fx1 = subs(F,findsym(F),x1);h = x0 - x1;J = zeros(n,n);for i=1:nxt = x1;xt(i) = x0(i);J(:,i) = (subs(F,findsym(F),xt)-fx1)/h(i);endr=x1-inv(J)*fx1;m=1;tol=1;while tol>epsx0 = x1;x1 = r;fx = subs(F,findsym(F),x0);fx1 = subs(F,findsym(F),x1);h = x0 - x1;J = zeros(n,n);for i=1:nxt = x1;xt(i) = x0(i);J(:,i) = (subs(F,findsym(F),xt)-fx1)/h(i);endr=x1-inv(J)*fx1;tol=norm(r-x1);m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endendformat short;8.mulGXF2用两点割线法的第二种形式求非线性方程组的一个根function [r,m]=mulGXF2(F,x0,x1,eps)format long;if nargin==3eps=1.0e-4;endx0 = transpose(x0);x1 = transpose(x1);n = length(x0);fx = subs(F,findsym(F),x0);fx1 = subs(F,findsym(F),x1);h = x0 - x1;J = zeros(n,n);xt = x1;xt(1) = x0(1);J(:,1) = (subs(F,findsym(F),xt)-subs(F,findsym(F),x1))/h(1);for i=2:nxt = x1;xt(1:i) = x0(1:i);xt_m = x1;xt_m(1:i-1) = x0(1:i-1);J(:,i) = (subs(F,findsym(F),xt)-subs(F,findsym(F),xt_m))/h(i);endr=x1-inv(J)*fx1;m=1;tol=1;while tol>epsx0 = x1;x1 = r;fx = subs(F,findsym(F),x0);fx1 = subs(F,findsym(F),x1);h = x0 - x1;J = zeros(n,n);xt = x1;xt(1) = x0(1);J(:,1) = (subs(F,findsym(F),xt)-subs(F,findsym(F),x1))/h(1);for i=2:nxt = x1;xt(1:i) = x0(1:i);xt_m = x1;xt_m(1:i-1) = x0(1:i-1);J(:,i) = (subs(F,findsym(F),xt)-subs(F,findsym(F),xt_m))/h(i);endr=x1-inv(J)*fx1;tol=norm(r-x1);m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endendformat short;9.mulVNewton用拟牛顿法求非线性方程组的一组解function [r,m]=mulVNewton(F,x0,A,eps)%方程组：F%方程组的初始解：x0% 初始A矩阵：A%解的精度：eps%求得的一组解：r%迭代步数：mif nargin==2A=eye(length(x0)); %A取为单位阵eps=1.0e-4;elseif nargin==3eps=1.0e-4;endendx0 = transpose(x0);Fx = subs(F, findsym(F),x0);r=x0-A\Fx;m=1;tol=1;while tol>epsx0=r;Fx = subs(F, findsym(F),x0);r=x0-A\Fx;y=r-x0;Fr = subs(F, findsym(F),r);z= Fr-Fx;A1=A+(z-A*y)*transpose(y)/norm(y); %调整A A=A1;m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endtol=norm(r-x0);end10.mulRank1用对称秩1算法求非线性方程组的一个根function [r,n]=mulRank1(F,x0,A,eps)if nargin==2l = length(x0);A=eye(l); %A取为单位阵eps=1.0e-4;elseif nargin==3eps=1.0e-4;endendfx = subs(F,findsym(F),x0);r=transpose(x0)-inv(A)*fx;n=1;tol=1;while tol>epsx0=r;fx = subs(F,findsym(F),x0);r=x0-inv(A)*fx;y=r-x0;fr = subs(F,findsym(F),r);z = fr-fx;A1=A+ fr *transpose(fr)/(transpose(fr)*y); %调整A A=A1;n=n+1;if(n>100000) %迭代步数控制disp('迭代步数太多，可能不收敛！');return;endtol=norm(r-x0);end11.mulDFP用D-F-P算法求非线性方程组的一组解function [r,n]=mulDFP(F,x0,A,eps)if nargin==2l = length(x0);B=eye(l); %A取为单位阵eps=1.0e-4;elseif nargin==3eps=1.0e-4;endendfx = subs(F,findsym(F),x0);r=transpose(x0)-B*fx;n=1;tol=1;while tol>epsx0=r;fx = subs(F,findsym(F),x0);r=x0-B*fx;y=r-x0;fr = subs(F,findsym(F),r);z = fr-fx;B1=B+ y*y'/(y'*z)-B*z*z'*B/(z'*B*z); %调整AB=B1;n=n+1;if(n>100000) %迭代步数控制disp('迭代步数太多，可能不收敛！');return;endtol=norm(r-x0);end12.mulBFS用B-F-S算法求非线性方程组的一个根function [r,n]=mulBFS(F,x0,B,eps)if nargin==2l = length(x0);B=eye(l); %B取为单位阵eps=1.0e-4;elseif nargin==3eps=1.0e-4;endendfx = subs(F,findsym(F),x0);r=transpose(x0)-B*fx;n=1;tol=1;while tol>epsx0=r;fx = subs(F,findsym(F),x0);r=x0-B*fx;y=r-x0;fr = subs(F,findsym(F),r);z = fr-fx;u = 1 + z'*B*z/(y'*z);B1= B+ (u*y*y'-B*z*y'-y*z'*B)/(y'*z); %调整B B=B1;n=n+1;if(n>100000) %迭代步数控制disp('迭代步数太多，可能不收敛！');return;endtol=norm(r-x0);end13.mulNumYT用数值延拓法求非线性方程组的一组解function [r,m]=mulNumYT(F,x0,h,N,eps)format long;if nargin==4eps=1.0e-8;endn = length(x0);fx0 = subs(F,findsym(F),x0);x0 = transpose(x0);J = zeros(n,n);for k=0:N-1fx = subs(F,findsym(F),x0);for i=1:nx1 = x0;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-fx)/h(i);endinJ = inv(J);r=x0-inJ*(fx-(1-k/N)*fx0);x0 = r;endm=1;tol=1;while tol>epsxs=r;fx = subs(F,findsym(F),xs);J = zeros(n,n);for i=1:nx1 = xs;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-fx)/h(i);endr=xs-inv(J)*fx; %核心迭代公式tol=norm(r-xs);m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endendformat short;14.DiffParam1用参数微分法中的欧拉法求非线性方程组的一组解function r=DiffParam1(F,x0,h,N)%非线性方程组：f%初始解：x0%数值微分增量步大小：h%雅可比迭代参量:l%解的精度：eps%求得的一组解：r%迭代步数：nx0 = transpose(x0);n = length(x0);ht = 1/N;Fx0 = subs(F,findsym(F),x0);for k=1:NFx = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)+h(i);J(:,i) = (subs(F,findsym(F),x1)-Fx)/h(i);endinJ = inv(J);r = x0 - ht*inJ*Fx0;x0 = r;end15.DiffParam2用参数微分法中的中点积分法求非线性方程组的一组解function r=DiffParam2(F,x0,h,N)%非线性方程组：f%初始解：x0%数值微分增量步大小：h%雅可比迭代参量:l%解的精度：eps%求得的一组解：r%迭代步数：nx0 = transpose(x0);n = length(x0);ht = 1/N;Fx0 = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nxt = x0;xt(i) = xt(i)+h(i);J(:,i) = (subs(F,findsym(F),xt)-Fx0)/h(i);endinJ = inv(J);x1 = x0 - ht*inJ*Fx0;for k=1:Nx2 = x1 + (x1-x0)/2;Fx2 = subs(F,findsym(F),x2);J = zeros(n,n);for i=1:nxt = x2;xt(i) = xt(i)+h(i);J(:,i) = (subs(F,findsym(F),xt)-Fx2)/h(i);endinJ = inv(J);r = x1 - ht*inJ*Fx0;x0 = x1;x1 = r;end16.mulFastDown用最速下降法求非线性方程组的一组解function [r,m]=mulFastDown(F,x0,h,eps)format long;if nargin==3eps=1.0e-8;endn = length(x0);x0 = transpose(x0);m=1;tol=1;while tol>epsfx = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)+h;J(:,i) = (subs(F,findsym(F),x1)-fx)/h;endlamda = fx/sum(diag(transpose(J)*J));r=x0-J*lamda; %核心迭代公式fr = subs(F,findsym(F),r);tol=dot(fr,fr);x0 = r;m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endendformat short;17.mulGSND用高斯牛顿法求非线性方程组的一组解function [r,m]=mulGSND(F,x0,h,eps)format long;if nargin==3eps=1.0e-8;endn = length(x0);x0 = transpose(x0);m=1;tol=1;while tol>epsfx = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)+h;J(:,i) = (subs(F,findsym(F),x1)-fx)/h;endDF = inv(transpose(J)*J)*transpose(J);r=x0-DF*fx; %核心迭代公式tol=norm(r-x0);x0 = r;m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endendformat short;18.mulConj用共轭梯度法求非线性方程组的一组解function [r,m]=mulConj(F,x0,h,eps)format long;if nargin==3eps=1.0e-6;endn = length(x0);x0 = transpose(x0);fx0 = subs(F,findsym(F),x0);p0 = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)*(1+h);p0(:,i) = -(subs(F,findsym(F),x1)-fx0)/h;endm=1;tol=1;while tol>epsfx = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)+h;J(:,i) = (subs(F,findsym(F),x1)-fx)/h;endlamda = fx/sum(diag(transpose(J)*J));r=x0+p0*lamda; %核心迭代公式fr = subs(F,findsym(F),r);Jnext = zeros(n,n);for i=1:nx1 = r;x1(i) = x1(i)+h;Jnext(:,i) = (subs(F,findsym(F),x1)-fr)/h;endabs1 = transpose(Jnext)*Jnext;abs2 = transpose(J)*J;v = abs1/abs2;if (abs(det(v)) < 1)p1 = -Jnext+p0*v;elsep1 = -Jnext;endtol=norm(r-x0);p0 = p1;x0 = r;m=m+1;if(m>100000) %迭代步数控制 disp('迭代步数太多，可能不收敛！');return;endendformat short;19.mulDamp用阻尼最小二乘法求非线性方程组的一组解function [r,m]=mulDamp(F,x0,h,u,v,eps)format long;if nargin==5eps=1.0e-6;endFI = transpose(F)*F/2;n = length(x0);x0 = transpose(x0);m=1;tol=1;while tol>epsj = 0;fx = subs(F,findsym(F),x0);J = zeros(n,n);for i=1:nx1 = x0;x1(i) = x1(i)+h;afx = subs(F,findsym(F),x1);J(:,i) = (afx-fx)/h;endFIx = subs(FI,findsym(FI),x0);for i=1:nx2 = x0;x2(i) = x2(i)+h;gradFI(i,1) = (subs(FI,findsym(FI),x2)-FIx)/h;ends=0;while s==0A = transpose(J)*J+u*eye(n,n);p = -A\gradFI;r = x0 + p;FIr = subs(FI,findsym(FI),r);if FIr<FIxif j == 0u = u/v;j = 1;elses=1;endelseu = u*v;j = 1;if norm(r-x0)<epss=1;endendendx0 = r;tol = norm(p);m=m+1;if(m>100000) %迭代步数控制disp('迭代步数太多，可能不收敛！');return;endendformat short;。
matlab中fsolve函数的用法在MATLAB中，fsolve函数用于求解非线性方程组。
它的基本语法如下：[x, fval, exitflag] = fsolve(fun, x0)其中，fun是一个函数句柄，表示要求解的非线性方程组。
x0是一个初始猜测值，表示方程组的解的初始估计。
x是求解得到的方程组的解向量。
fval是方程组在x处的函数值向量。
exitflag是一个整数，表示求解的终止条件。
以下是一些使用fsolve函数的例子，用于说明其用法：1.求解一元非线性方程：```matlabx0=1;[x, fval, exitflag] = fsolve(fun, x0);```在这个例子中，我们定义了一个匿名函数fun，表示非线性方程2*sin(x) - x = 0。
然后我们使用初始猜测值x0 = 1来调用fsolve函数求解方程。
得到的解存储在x变量中，函数值存储在fval变量中。
2.求解多元非线性方程组：```matlabx0=[1;2];[x, fval, exitflag] = fsolve(fun, x0);```在这个例子中，我们定义了一个匿名函数fun，表示非线性方程组{x1^2 + x2^2 - 1 = 0, x1 - x2 = 0}。
然后我们使用初始猜测向量x0 = [1; 2]来调用fsolve函数求解方程组。
得到的解存储在x变量中，函数值存储在fval变量中。
3.使用辅助函数求解方程：```matlabx0=1;[x, fval, exitflag] = fsolve(fun, x0);function y = myfun(x)y=x^2-2;end```在这个例子中，我们定义了一个辅助函数myfun，表示非线性方程x^2 - 2 = 0。
然后我们使用初始猜测值x0 = 1来调用fsolve函数求解方程。
得到的解存储在x变量中，函数值存储在fval变量中。
需要注意的是，fsolve函数对于非线性方程组的求解是基于数值方法的，所以有时候可能无法找到方程组的解，或者找到的解可能是局部最优解。
Matlab方程组解1. 引言方程组是数学中一个重要的概念，它描述了多个未知数之间的关系。
解方程组的过程在科学、工程和计算机科学等领域中有着广泛的应用。
Matlab作为一种高级数值计算环境，提供了丰富的工具和函数来解决方程组的求解问题。
本文将介绍如何使用Matlab解方程组，包括线性方程组和非线性方程组的求解方法。
2. 线性方程组的求解2.1 利用矩阵求解线性方程组可以表示为矩阵形式，例如：Ax = b，其中A是系数矩阵，x是未知数向量，b是常数向量。
在Matlab中，可以使用线性代数工具箱中的函数来求解线性方程组。
2.1.1 使用inv函数求解如果系数矩阵A是可逆的，可以使用inv函数求解线性方程组。
具体步骤如下： 1. 计算A的逆矩阵：A_inv = inv(A) 2. 计算解向量：x = A_inv * b2.1.2 使用linsolve函数求解linsolve函数可以直接求解线性方程组，无需计算逆矩阵。
具体步骤如下： 1. 调用linsolve函数：x = linsolve(A, b)2.2 利用高斯消元法求解高斯消元法是一种常用的求解线性方程组的方法，它通过矩阵的行变换将方程组转化为上三角矩阵，然后通过回代得到解。
在Matlab中，可以使用lu函数来进行高斯消元法求解。
2.2.1 使用lu函数求解lu函数可以将方程组的系数矩阵A分解为一个下三角矩阵L和一个上三角矩阵U，即A = LU。
具体步骤如下： 1. 调用lu函数：[L, U] = lu(A) 2. 解得方程组：x = U \ (L \ b)3. 非线性方程组的求解非线性方程组是指未知数与其函数之间存在非线性关系的方程组。
与线性方程组不同，非线性方程组的求解通常需要借助数值方法。
Matlab提供了多种函数和工具箱来解决非线性方程组的求解问题。
3.1 利用fsolve函数求解fsolve函数是Matlab中用于求解非线性方程组的函数，它通过迭代的方式逼近方程组的解。
Matlab⾮线性⽅程求解器fsolve总结（含实例）Matlab⾮线性⽅程求解器fsolve总结（含实例）fsolve是采⽤最⼩⼆乘法来求解⾮线性⽅程（组）。
它的⼀般求解⽅式为：X=FSOLVE(FUN,X0,OPTIONS)其中，fun是要求解的⾮线性⽅程，X0是变量初值，options由optimset函数产⽣的结构体，⽤于对优化参数的设置，可以省略（采⽤默认值）。
fsolve可以求解简单的⼀维⾮线性⽅程，如：x = fsolve(@myfun,[0.5 2 4],optimset('Display','iter')); %求解在初值分别为0.5，2和4时⽅程的解其中，函数myfun的定义为：function F = myfun(x)F = sin(x);fsolve还可以求解⼤型的⾮线性⽅程组，如x0 = [51.6;rand;unifrnd(-1,1);rand];h=optimset;h.MaxFunEvals=20000;h.MaxIter=5000;h.Display='off';[p,fval] = fsolve(@f,x0,options);此时，⽅程组可以写成矩阵形式：function F=f(x)f=[x(1)+x(2)*(1-exp(-(x(3)*(0)^x(4))))-51.61;x(1)+x(2)*(1-exp(-(x(3)*(9.78)^x(4))))-51.91;x(1)+x(2)*(1-exp(-(x(3)*(30.68)^x(4))))-53.27;x(1)+x(2)*(1-exp(-(x(3)*(59.7)^x(4))))-59.68;];以上是的⾮线性⽅程（组）都是简单的函数关系式，但是实际应⽤中的函数，例如在图像处理中的函数关系可能往往是⼀个隐函数，所以优化的⽬标⽅程可能是⼀个函数⽂件的结果，这种时候也是可以⽤来求解的。
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 comes 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%% alpha = 1.d-4, parameter to measure sufficient decrease %% maxarm = 10, maximum number of steplength reductions before % failure is reported%% set the debug parameter, 1 turns display on, otherwise off%debug=1;%% initialize it_hist, ierr, and set the iteration parameters%ierr = 0; maxit=40; maxdim=39;it_histx=zeros(maxit,3);maxarm=10;%if nargin == 4maxit=parms(1); maxdim=parms(2)-1;endrtol=tol(2); atol=tol(1); n = length(x); fnrm=1; itc=0; nbroy=0; %% evaluate f at the initial iterate% compute the stop tolerance%f0=feval(f,x);fc=f0;fnrm=norm(f0)/sqrt(n);it_hist(itc+1)=fnrm;it_histx(itc+1,1)=fnrm; it_histx(itc+1,2)=0;it_histx(itc+1,3)=0;fnrmo=1;stop_tol=atol + rtol*fnrm;outstat(itc+1, :) = [itc fnrm 0 0];%% terminate on entry?%if fnrm < stop_tolsol=x;returnend%% initialize the iteration history storage matrices%stp=zeros(n,maxdim);stp_nrm=zeros(maxdim,1);lam_rec=ones(maxdim,1);%% Set the initial step to -F, compute the step norm%lambda=1;stp(:,1) = -fc;stp_nrm(1)=stp(:,1)'*stp(:,1);%% main iteration loop%while(itc < maxit)%nbroy=nbroy+1;%% keep track of successive residual norms and% the iteration counter (itc)%fnrmo=fnrm; itc=itc+1;%% compute the new point, test for termination before% adding to iteration history%xold=x; lambda=1; iarm=0; lrat=.5; alpha=1.d-4;x = x + stp(:,nbroy);fc=feval(f,x);fnrm=norm(fc)/sqrt(n);ff0=fnrmo*fnrmo; ffc=fnrm*fnrm; lamc=lambda;%%% Line search, we assume that the Broyden direction is an% ineact Newton direction. If the line search fails to% find sufficient decrease after maxarm steplength reductions % brsola returns with failure.%% Three-point parabolic line search%while fnrm >= (1 - lambda*alpha)*fnrmo && iarm < maxarm% lambda=lambda*lrat;if iarm==0lambda=lambda*lrat;elselambda=parab3p(lamc, lamm, ff0, ffc, ffm);endlamm=lamc; ffm=ffc; lamc=lambda;x = xold + lambda*stp(:,nbroy);fc=feval(f,x);fnrm=norm(fc)/sqrt(n);ffc=fnrm*fnrm;iarm=iarm+1;end%% set error flag and return on failure of the line search%if iarm == maxarmdisp('Line search failure in brsola ')ierr=2;it_hist=it_histx(1:itc+1,:);sol=xold;return;end%% How many function evaluations did this iteration require?%it_histx(itc+1,1)=fnrm;it_histx(itc+1,2)=it_histx(itc,2)+iarm+1;if(itc == 1) it_histx(itc+1,2) = it_histx(itc+1,2)+1; end;it_histx(itc+1,3)=iarm;%% terminate?%if fnrm < stop_tolsol=x;rat=fnrm/fnrmo;outstat(itc+1, :) = [itc fnrm iarm rat];it_hist=it_histx(1:itc+1,:);% it_hist(itc+1)=fnrm;if debug==1disp(outstat(itc+1,:))endreturnend%%% modify the step and step norm if needed to reflect the line % search%lam_rec(nbroy)=lambda;if lambda ~= 1stp(:,nbroy)=lambda*stp(:,nbroy);stp_nrm(nbroy)=lambda*lambda*stp_nrm(nbroy);end%%% it_hist(itc+1)=fnrm;rat=fnrm/fnrmo;outstat(itc+1, :) = [itc fnrm iarm rat];if debug==1disp(outstat(itc+1,:))end%%% if there's room, compute the next search direction and step norm and% add to the iteration history%if nbroy < maxdim+1z=-fc;if nbroy > 1for kbr = 1:nbroy-1ztmp=stp(:,kbr+1)/lam_rec(kbr+1);ztmp=ztmp+(1 - 1/lam_rec(kbr))*stp(:,kbr);ztmp=ztmp*lam_rec(kbr);z=z+ztmp*((stp(:,kbr)'*z)/stp_nrm(kbr));endend%% store the new search direction and its norm%a2=-lam_rec(nbroy)/stp_nrm(nbroy);a1=1 - lam_rec(nbroy);zz=stp(:,nbroy)'*z;a3=a1*zz/stp_nrm(nbroy);a4=1+a2*zz;stp(:,nbroy+1)=(z-a3*stp(:,nbroy))/a4;stp_nrm(nbroy+1)=stp(:,nbroy+1)'*stp(:,nbroy+1);%%%else%% out of room, time to restart%stp(:,1)=-fc;stp_nrm(1)=stp(:,1)'*stp(:,1);nbroy=0;%%%end%% end whileend%% We're not supposed to be here, we've taken the maximum% number of iterations and not terminated.%sol=x;it_hist=it_histx(1:itc+1,:);ierr=1;if debug==1disp(' outstat')endfunction lambdap = parab3p(lambdac, lambdam, ff0, ffc, ffm)% Apply three-point safeguarded parabolic model for a line search. %% This code comes with no guarantee or warranty of any kind.%% function lambdap = parab3p(lambdac, lambdam, ff0, ffc, ffm)%% input:% lambdac = current steplength% lambdam = previous steplength% ff0 = value of \| F(x_c) \|^2% ffc = value of \| F(x_c + \lambdac d) \|^2% ffm = value of \| F(x_c + \lambdam d) \|^2%% output:% lambdap = new value of lambda given parabolic model%% internal parameters:% sigma0 = .1, sigma1=.5, safeguarding bounds for the linesearch%%% set internal parameters%sigma0=.1; sigma1=.5;%% compute coefficients of interpolation polynomial%% p(lambda) = ff0 + (c1 lambda + c2 lambda^2)/d1%% d1 = (lambdac - lambdam)*lambdac*lambdam < 0% so if c2 > 0 we have negative curvature and default to% lambdap = sigam1 * lambda%c2 = lambdam*(ffc-ff0)-lambdac*(ffm-ff0);if c2 >= 0lambdap = sigma1*lambdac; returnendc1=lambdac*lambdac*(ffm-ff0)-lambdam*lambdam*(ffc-ff0);lambdap=-c1*.5/c2;if (lambdap < sigma0*lambdac) lambdap=sigma0*lambdac; endif (lambdap > sigma1*lambdac) lambdap=sigma1*lambdac; end2.应用举例把以下代码复制在command 窗口中x=[1 2 3]’;f=@(x)[3*x(1)-cos(x(2)*x(3))-1/2;x(1)^2-81*(x(2)+0.1)^2+sin(x(3))+1.06;exp(-x(1)*x(2))+20*x(3)+(10*pi-3)/3;];tol=[3,-5];[sol, it_hist, ierr] = brsola(x,f,tol)说明：以上应用举例只是给出了上文中代码的一个应用实例，具体能否得到方程的满意数值解还需要进一步调节初始给的x和tol的值。
实验四非线性方程近似解一、按揭还贷㈠问题描述（1）小张夫妇以按揭方式贷款买了一套价值20万元的房子，首付5万元，每月还款1000元，15年还清。
问贷款利率是多少？（2）某人想贷款50万元购房，他咨询了两家银行，第一家银行开出的条件是每月还4500元，15年还清；第二家银行开出的条件是每年还45000元，20年还清。
从利率方面看，哪家银行较优惠？（简单假设年利率=月利率*12）㈡简要分析初看本题，一个简单的思路是每次测试一个利率值，以这个值为基础计算15年后所剩还款数量，通过结果判断应将利率值增大或减小，从而实现迭代。
这其实是一个二重迭代的过程，之所以这样是因为不容易一眼看出本题的非线性方程。
事实上，转换思路后，可以利用一个简单的方程描述整个迭代过程。
这样就将二重迭代转化为了一层迭代。
使得处理更加简便。
㈢方法与公式1、解题方法(1)二次迭代给定总的本金，从每一次还款中扣去这段时间中增加的利息，再将其还到本金，使本金总量逐渐减少。
代码：for i = 1:time*12less = (repay-left*interest);left = left - less;(2)方程描述虽然并不是所有本金都在还款的整个期间中产生了相应的利息，但是可以设想成这样，与此同时，还款从在相应的还款时间开始产生利息，这样可以得出，两者最终的“本息和”相等，即nA(1+q)n=P(1+q)n−ii=1其中A为总还款金额，q为了利率，P为每次还款金额。
2、解方程方法(1)牛顿法x k+1=x k−f(x k) f′(x k)(2)直接使用公式fzero()㈣结果与分析1、第一问：(1)二次迭代[i,q]=iterate(150000,1000,15,2,0,1,100,10^-6); 公式表意为：总贷款量=200000-50000=150000；每月还款100元；还款期限15年；还款方式为按月还款；迭代区间设定为[0,1]；最大迭代次数为100次；精度要求为10^-6；最终结果为：迭代次数：45；使用时间0.003030989435705s；利率为0.002081163889457。
