查看: 2433|回复: 4|关注: 0
求助matlab实现牛顿迭代法求解非线性方程中遇到的问题
<h1 style="color:# 麦片财富积分
新手, 积分 7, 距离下一级还需 43 积分
最近在研究weibull分布的参数估计，其中有一个环节为matlab编程利用最大似然函数求解估计参数，在此就不可避免的要面临求解非线性方程。
本人在网上搜索多个文献，其中被大家引用最多的一个例子就是&matlab实现牛顿迭代法求解非线性方程&。通过试着运行作者的matlab code，产生了如下一些疑问，请各位老师帮忙解答，不胜感谢。
syms x1 x2 x3 % 定义三个未知数
f1=3*x1-cos(x2*x3)-1/2;
f2=x1^2-81*(x2+0.1)^2+sin(x3)+1.06;
f3=exp(-x1*x2)+20*x3+(10*pi-3)/3; % 三个方程组成一个方程组并且f1=0;f2=0;f3=0
f=[f1 f2 f3];
df=[diff(f,x1);diff(f,x2);diff(f,x3)]; % 分别对x1 x2 x3 三个未知数求偏导
x0=[0.1 0.1 -0.1];
eps=0.00001;
N=20; % 设置初始值 x1 x2 x3 的初始值分别为0.1 0.1 -0.1，精确度0.00001 迭代运算次数设置为最多20次
for i=1:N; % 运算次数
& & f=subs(f,{'x1' 'x2' 'x3'},{x0(1) x0(2) x0(3)}); % 为subs函数的赋值运算 1x3的矩阵
& & df=subs(df,{'x1' 'x2' 'x3'},{x0(1) x0(2) x0(3)});% 为subs函数的赋值运算 3x3的矩阵
& & x=x0-f/ % 得到的一组新数值 1x3的矩阵 理解为 新的x0=[新x0(1) 新x0(2) 新x0(3)]
& &if norm(x-x0)&eps
& &end % 条件判断 得出的新值是否满足精确度的要求 可以求x矩阵元素平方和的开根
& &x0=x;% 当条件满足 输出求得的x0
运行上述程式出现的问题如下
当第一次循环初始值为x0=[0.1 0.1 -0.1]代入，得到第一次循环的计算结果 f/df得到一个结果 进而x得到一个结果
但是当第二次循环 f/df的值就已经保持不变 仅有x和x0在循环计算
当然这不是我们想要
请各位老师指点
PS 有没有更好的迭代的在matlab里的实现方法可以更好的缩短结算时间？
<h1 style="color:# 麦片财富积分
|此回复为最佳答案
lz方法没问题，但是程序里面有错误，按照正确的程序可以算出：
<h1 style="color:# 麦片财富积分
我也重新查看一下 问题出在在循环语句中 f df 均被替代 成以定义的 f df 所以保留了代入初始值所得的结果 无法随着循环 运算 可命名改一下 即可 fnew dfnew 再次感谢 您的解答 不过 我们的结果有些出入
<h1 style="color:# 麦片财富积分
楼组，采用的几参数的威布尔分布呀？
<h1 style="color:# 麦片财富积分
楼组，采用的几参数的威布尔分布呀？
采用威布尔分布2参数
站长推荐 /2
Powered by9484人阅读
课程实验（6）
题目：用牛顿法求方程x-cos(x)=0的实根（精确到1E-6）。
（1）要求用函数调用。
（2）进一步研究和弦截法作比较。
算法分析：
（1）&&&&&& 此题是利用牛顿方法解一元非线性方程的根。（牛顿法是把非线性方程局部线性化的一种方法，它在单根附近具有较高的收敛速度。）所以首先我们应先给出估计的根，先对方程x-cos(x)=0变形，令y1=x,y2=cos(x),则两函数图象的交点，就是方程x-cos(x)=0的根，这里利用Matlab作图估计根的&#20540;。
在Matlab命令行中输入，并运行：
&&x=-2:0.01:2;
&&y1=cos(x);
&&plot(x,y1,x,y2);
可得下图：
从图中可以很容易得到根x的初&#20540;可选0.6。
（2）利用牛顿法的迭代公式x1=x0-(f(x0)/f’(x0));这里首先取x0=0.6,代入迭代公式，然后判断x1与x0之差的绝对&#20540;是否小于精度，如果小于精度，则停止，即得出根x的&#20540;，如果不小于精度，则继续迭代，直到符合精度为止。
Matlab代码如下：
（1）funNewton.m函数文件为：
%此文件为被调函数文件。
function[y,dirv_y]=funNewton(x)
%y与dirv_y为函数返回的参数，x为调用时的传递参数。
y=x-cos(x);
dirv_y=1+sin(x);
%dirv_y为y的一阶导函数
（2）Newton.m脚本文件为：
%清除所有变量
Error=1e-6;
%根所要求的误差精度
formatlong
%格式化此行一下的变量。
%迭代次数20次
[y,dirv_y]=funNewton(x);
%调用函数文件
%xk可保留前一次x的值。
x=x-y/dirv_y;
%牛顿法的核心，即此迭代公式。
if(abs(xk-x)&=Error)
%判断当此x的值与前一次x的值即xk的差值，即误差是否小于题目给定的误差。
%输出根x的值。
在Matlab的命令行中输入（因为我的程序文件保存为：Newton.m脚本文件，funNewton.m函数文件）：
结果分析：
在Matlab的命令行中输入：
&&format long
&&fzero('x-cos(x)',0)
经比较可知：牛顿法的结果有一定的误差，但是牛顿方法由于在单根附近有良好的收敛性，所以与其他方法得出的结果相比误差较小。但是牛顿方法有个缺点：它只在根附近局部收敛。所以所我们在给定X的初&#20540;时尤为重要，如果给的初&#20540;离真&#20540;过远。那么用牛顿法可能永远也找不到此方程的解。对于这种情况，我们可以先利用工具软件（如：Matlab）画出草图，确定根的大致位置。牛顿法的每一步迭代都要计算一次导数&#20540;，而在计算机中，计算一次导数的近&#20284;&#20540;比计算函数&#20540;要麻烦的多。所以我们可以用弦截法，其迭代公式为：x2=x1-((x1-x0)/f(x1)-f(x0))*f(x1).
当f(x)在x*的某一个邻域内具有二阶连续导数，且f’(x)!=0,初&#2落在此邻域内，弦截
法是收敛的。
此文为本人原创！如需使用，请注明出处！
&&相关文章推荐
参考知识库
* 以上用户言论只代表其个人观点，不代表CSDN网站的观点或立场
访问：64265次
积分：1557
积分：1557
排名：千里之外
原创：34篇
评论：48条
(5)(2)(2)(1)(1)(3)(1)(2)(1)(2)(2)(1)(11)(1)(1)(3)分享给朋友：通用代码: <input id="link4" type="text" class="form_input form_input_s" value="" />复 制flash地址: 复 制html代码: <input type="text" class="form_input form_input_s" id="link3" value="" />复 制分享视频到站外获取收益&&手机扫码分享视频二维码2小时内有效MATLAB教学视频：非线性方程（组）在MATLAB中的求解方法下载至电脑扫码用手机看用或微信扫码在手机上继续观看二维码2小时内有效MATLAB教学视频：非线性方程（组）在MATLAB中的求解方法扫码用手机继续看用或微信扫码在手机上继续观看二维码2小时内有效，扫码后可分享给好友没有优酷APP？立即下载请根据您的设备选择下载版本
药品服务许可证(京)-经营- 节目制作经营许可证京字670号 请使用者仔细阅读优酷、、Copyright(C)2017 优酷
版权所有不良信息举报电话:君，已阅读到文档的结尾了呢~~
扫扫二维码，随身浏览文档
手机或平板扫扫即可继续访问
第3章 Matlab非线性方程求解与迭代法
举报该文档为侵权文档。
举报该文档含有违规或不良信息。
反馈该文档无法正常浏览。
举报该文档为重复文档。
推荐理由：
将文档分享至：
分享完整地址
文档地址：
粘贴到BBS或博客
flash地址：
支持嵌入FLASH地址的网站使用
html代码：
&embed src='/DocinViewer-4.swf' width='100%' height='600' type=application/x-shockwave-flash ALLOWFULLSCREEN='true' ALLOWSCRIPTACCESS='always'&&/embed&
450px*300px480px*400px650px*490px
支持嵌入HTML代码的网站使用
您的内容已经提交成功
您所提交的内容需要审核后才能发布，请您等待！
3秒自动关闭窗口MATLAB求解非线性方程（转）
最近有多人问如何用matlab解方程组的问题，其实在matlab中解方程组还是很方便的，例如，对于代数方程组Ax=b(A为系数矩阵，非奇异)的求解，MATLAB中有两种方法：
(1)x=inv(A)*b —
采用求逆运算解方程组；
(2)x=A — 采用左除运算解方程组。
2x1+3x2=13
&&A=[1,2;2,3];b=[8;13];
&&x=inv(A)*b
即二元一次方程组的解x1和x2分别是2和3。
对于同学问到的用matlab解多次的方程组，有符号解法，方法是：先解出符号解,然后用vpa(F,n)求出n位有效数字的数值解.具体步骤如下：
第一步:定义变量syms x y z ...;
第二步:求解[x,y,z,...]=solve('eqn1','eqn2',...,'eqnN','var1','var2',...'varN');
第三步:求出n位有效数字的数值解x=vpa(x,n);y=vpa(y,n);z=vpa(z,n);...。
如：解二（多）元二（高）次方程组：
x^2+3*y+1=0
y^2+4*x+1=0
解法如下：
&&[x,y]=solve('x^2+3*y+1=0','y^2+4*x+1=0');
&&x=vpa(x,4);
&&y=vpa(y,4);
1.635+3.029*i
1.635-3.029*i
1.834-3.301*i
1.834+3.301*i
二元二次方程组，共4个实数根；
还有的同学问，如何用matlab解高次方程组（非符号方程组）？举个例子好吗？
解答如下：
基本方法是：solve(s1,s2,…,sn,v1,v2,…,vn)，即求表达式s1,s2,…,sn组成的方程组，求解变量分别v1,v2,…,vn。
具体例子如下：
x^2 + x*y + y = 3
x^2 - 4*x + 3 = 0
[x,y] = solve('x^2 + x*y + y = 3','x^2 - 4*x + 3 = 0')
运行结果为
即x等于1和3；y等于1和-1.5
= solve('x^2 + x*y + y = 3','x^2 - 4*x + 3=
0','x','y')
结果一样，二元二方程都是4个实根。
通过这三个例子可以看出，用matlab解各类方程组都是可以的，方法也有多种，只是用到解方程组的函数，注意正确书写参数就可以了，非常方便。
from:http://bbs./pc/pccon.php?id=950&nid=14498&tid=0
2、变参数非线性方程组的求解
对于求解非线性方程组一般用fsolve命令就可以了，但是对于方程组中某一系数是变化的，该怎么求呢？
%定义方程组如下，其中k为变量
function F = myfun(x,k)
Pc0=0.23;W=0.18;
F=[Pc0+H*(1+1.5*(x(1)/W-1)-0.5*(x(1)/W-1)^3)-x(2);
x(1)-k*sqrt(x(2))];
Pc0=0.23;W=0.18;
x0 = [2*W; Pc0+2*H]; %
optimset('Display','off');
k=0:0.01:1; % 变量取值范围[0
for i=1:1:length(k)
x = fsolve(@(x) myfun(x,kk), x0,
options);%求解非线性方程组
x1(i)=x(1);
x2(i)=x(2);
plot(k,x1,'-b',k,x2,'-r');
xlabel('k')
legend('x1','x2')
from:/archiver/tid-836299.html
3、非线性方程数值求解
matlab里solve如何使用,是否有别的函数可以代替它.
matlab里我解y=9/17*exp(-1/2*t)*17^(1/2)*sin(1/2*17^(1/2)*t)=0这样的方程为什么只得到0这一个解,如何可以的到1/2*17^(1/2)*t=n*(pi)这样一族解??
在matlab里面solve命令主要是用来求解代数方程（即多项式）的解,但是也不是说其它方程一个也不能解，不过求解非代数方程的能力相当有限，通常只能给出很特殊的实数解。(该问题给出的方程就是典型的超越方程，非代数方程)
从计算机的编程实现角度讲，如今的任何算法都无法准确的给出任意非代数方程的所有解，但是我们有很多成熟的算法来实现求解在某点附近的解。matlab也不例外，它也只能给出任意非代数方程在某点附近的解，函数有两个：fzero和fsolve,具体用法请用help或doc命令查询吧。如果还是不行，你还可以将问题转化为非线性最优化问题，求解非线性最优化问题的最优解，可以用的命令有：fminbnd,
fminsearch, fmincon等等。
*非线性方程数值求解
*单变量非线性方程求解
在MATLAB中提供了一个fzero函数，可以用来求单变量非线性方程的根。该函数的调用格式为：
z=fzero('fname',x0,tol,trace)
其中fname是待求根的函数文件名，x0为搜索的起点。一个函数可能有多个根，但fzero函数只给出离x0最近的那个根。tol控制结果的相对精度，缺省时取tol=eps，trace
指定迭代信息是否在运算中显示，为1时显示，为0时不显示，缺省时取trace=0。
例 求f(x)=x-10x+2=0在x0=0.5附近的根。
步骤如下：
建立函数文件funx.m。
function fx=funx(x)
fx=x-10.^x+2;
(2) 调用fzero函数求根。
z=fzero('funx',0.5)
**非线性方程组的求解
对于非线性方程组F(X)=0，用fsolve函数求其数值解。fsolve函数的调用格式为：
X=fsolve('fun',X0,option)
其中X为返回的解，fun是用于定义需求解的非线性方程组的函数文件名，X0是求根过程的初值，option为最优化工具箱的选项设定。最优化工具箱提供了20多个选项，用户可以使用optimset命令将它们显示出来。如果想改变其中某个选项，则可以调用optimset()函数来完成。例如，Display选项决定函数调用时中间结果的显示方式，其中‘off’为不显示，‘iter’表示每步都显示，‘final’只显示最终结果。
optimset(‘Display’,‘off’)将设定Display选项为‘off’。
例 求下列非线性方程组在(0.5,0.5) 附近的数值解。
(1) 建立函数文件myfun.m。
q=myfun(p)
q(1)=x-0.6*sin(x)-0.3*cos(y);
q(2)=y-0.6*cos(x)+0.3*sin(y);
在给定的初值x0=0.5,y0=0.5下，调用fsolve函数求方程的根。
x=fsolve('myfun',[0.5,0.5]',optimset('Display','off'))
将求得的解代回原方程，可以检验结果是否正确，命令如下：
q=myfun(x)
可见得到了较高精度的结果。
from:.cn/s/blog_56ef652d0100ebew.html
4、fsolve函数解方程
[X,FVAL,EXITFLAG,OUTPUT,JACOB]=FSOLVE(FUN,X0,...)
returns the
Jacobian of FUN at X.
FUN can be specified using
x = fsolve(@myfun,[2 3
4],optimset('Display','iter'))
where myfun is a MATLAB function such
function F = myfun(x)
F = sin(x);
FUN can also be an anonymous
x = fsolve(@(x) sin(3*x),[1
4],optimset('Display','off'))
If FUN is parameterized, you can use
anonymous functions to capture the
problem-dependent parameters. Suppose
you want to solve the system of
nonlinear equations given in the
function myfun, which is parameterized
by its second argument c. Here myfun
is an M-file function such as
function F = myfun(x,c)
F = [ 2*x(1) - x(2) -
exp(c*x(1))
-x(1) + 2*x(2) -
exp(c*x(2))];
To solve the system of equations for
a specific value of c, first assign the
value to c. Then create a
one-argument anonymous function that captures
that value of c and calls myfun with
two arguments. Finally, pass this anonymous
function to FSOLVE:
c = -1; % define parameter
x = fsolve(@(x)
myfun(x,c),[-5;-5])
已投稿到：
以上网友发言只代表其个人观点，不代表新浪网的观点或立场。