求教：非线性方程组迭代解法的编程问题
方程组如附件，X1=VA，X2=6A按照书上的例子，先编写bro：
function =bro(x0,eps)
%布鲁顿迭代法
if nargin==2, eps=1.0e-6
elseif nargin<2,error, return; end
A=eye(length(x0)); x1=x0-fbro(x0)/A; j=1;
while (norm(x1-x0)>=eps) & (j<=1000)
x0=x1; x1=x0-fbro(x0)/A;
p=x1-x0; q=fbro(x1)-fbro(x0);
A=A+(q-p*A)'*p/norm(p); j=j+1;
end
y=x1; n=j;
%% 再编写fbro
function y=fbro(x)
d=5; ve=3; ro=0.4;
y(1)=(x(1)*normcdf((log(x(1)/d)+0.02+0.5*x(2)*x(2))/x(2))-exp(-0.02)*d*normcdf((log(x(1)/d)+0.02-0.5*x(2)*x(2))/x(2)))-ve;
y(2)=(x(1)*x(2)*normcdf((log(x(1)/d)+0.02+0.5*x(2)*x(2))/x(2)))/ve-ro;
y=;
计算不出来结果，请大侠指点一二。谢谢
[ 本帖最后由 ChaChing 于 2010-3-1 08:26 编辑 ] 再命令窗口输入>>=bro(,eps)
算出结果。。。。。
y =
1.0e+003 *
0.4800 -9.1507
n =
1001
可是X2应该是大于零的才对。。。。。。
=bro(,eps)
y =
8.3036 0.0006
n =
1001
参数不同，答案也不同。。。。。很多时候参数输入，就报错
上面的也不是答案，是我J值设定的太小了，1000
方程里含正态分布，可以用牛顿法解麽？？
[ 本帖最后由 xinyuxf 于 2007-5-8 09:28 编辑 ] 请高手指点，用fsolve可以解决麽？
如果这个问题解决了，再变d/ve/ro为变量，读取一个2000*3的矩阵，求得x1/x2的2000个答案。
我用循环语句控制，i=i+1，可以求解，问题是我不会把答案输入，每次都是拷贝屏幕上的答案到TXT里整理半天，好傻哦，请高手指点如何保存答案？
*************************
我就这么整理，删除其他文字，只保留答案，再贴进EXCEL里面。。。。
va =
1.0e+007 *
Columns 1 through 5
0 0.000015191744873 0.000001959011473 0.000004423352757 0.000007518352076
Columns 6 through 10
0.000004005238744 0.000004556613173 0.000015579341581 0.000006500117041 0.000002330693525
Columns 11 through 15
0.000003743739380 0.000022860635561 0.000093524485905 0.000003339738823 0.000002639275663
Columns 16 through 20
0.000010126396430 0.000008367559064 0.000012712029811 0.000004602640151 0.000002300267011
Columns 21 through 25 用牛顿、布鲁顿迭代，都出不来结果
换成fsolve,JIU :victory:
就出来了，太好了
******************
>> x0=;
>> options=optimset('Display','iter');
>> =fsolve(@myfun,x0,options)
Norm of First-order Trust-region
IterationFunc-count f(x) step optimality radius
0 3 8.83889 0.721 1
1 6 5.56837 1 1.36 1
2 9 2.13659 2.5 0.363 2.5
3 10 2.13659 6.25 0.363 6.25
4 13 1.41694 1.5625 0.499 1.56
5 14 1.41694 3.90625 0.499 3.91
6 17 0.905391 0.976562 0.225 0.977
7 18 0.905391 2.44141 0.225 2.44
8 21 0.633925 0.610352 0.493 0.61
9 24 0.398121 1.52588 0.949 1.53
10 27 0.000471548 0.546743 0.0569 1.53
11 30 1.98947e-009 0.00828717 8.48e-005 1.53
12 33 7.81081e-020 3.3332e-005 6.61e-010 1.53
Optimization terminated: first-order optimality is less than options.TolFun.
x =
7.9008
0.1520
fv =
1.0e-009 *
0.1175
-0.2536
**************************
下面再接着研究如何带入2000*3的变量矩阵，呵呵
还有导出答案 function F=myfun(x)
d=5; ve=3; ro=0.4;
F=[((x(1)*normcdf((log(x(1)/d)+0.02+0.5*x(2)*x(2))/x(2))-exp(-0.02)*d*normcdf((log(x(1)/d)+0.02-0.5*x(2)*x(2))/x(2)))-ve);((x(1)*x(2)*normcdf((log(x(1)/d)+0.02+0.5*x(2)*x(2))/x(2)))/ve-ro)];
[ 本帖最后由 ChaChing 于 2010-3-1 08:32 编辑 ] 搞定了第二步，可以带入变量了
>> mymain
Optimization terminated: first-order optimality is less than options.TolFun.
Optimization terminated: first-order optimality is less than options.TolFun.(repeat)
ans =
39.4965 13.7142 27.6819 30.1427146.0034 15.6690 29.6886 43.6006 30.3374 54.3623 48.7404
>>
[ 本帖最后由 ChaChing 于 2010-3-1 08:27 编辑 ] function =mymain
i=1;
x1=[];
x2=[];
dd=[];
ld=;
ro=;
ve=;
std=;
while i<=11
x0=;
x=fsolve(@myfun2,x0,'[]',ro(i),ve(i),ld(i));
x1(i)=x(1);
x2(i)=x(2);
dd(i)=(x(1)-std(i))/(x(1)*x(2));
i=i+1;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function F=myfun2(x,ro,ve,ld)
F=[((x(1)*normcdf((log(x(1)/ld)+0.02+0.5*x(2)*x(2))/x(2))-exp(-0.02)*ld*normcdf((log(x(1)/ld)+0.02-0.5*x(2)*x(2))/x(2)))-ve);((x(1)*x(2)*normcdf((log(x(1)/ld)+0.02+0.5*x(2)*x(2))/x(2)))/ve-ro)]; 现在已经解决了参数的传递问题，最后一步，就是参数的导入，以及方程解的导出
这个该用什么命令呢？
为什么10楼的程序，会出来9楼的结果？？？只显示了X1，没有显示X2......
修改了options
options=optimset('Display','off');
x=fsolve(@myfun2,x0,options,ro(i),ve(i),ld(i));
结果就没有前面的信息了，只有X1的解，可是X2呢？？
[ 本帖最后由 ChaChing 于 2010-3-1 13:36 编辑 ] >> mymain
ans =
39.4965 13.7142 27.6819 30.1427146.0034 15.6690 29.6886 43.6006 30.3374 54.3623 48.7404
>> =mymain
%%%%%%%%%%%%%%%%%%%
hahahahah，自己全部搞定！～～～～～～～～～太佩服自己了，接触MATLAB20个小时，就解决了方案。
>> =mymain
x1 =
39.4965 13.7142 27.6819 30.1427146.0034 15.6690 29.6886 43.6006 30.3374 54.3623 48.7404
x2 =
0.0790 0.0866 0.0877 0.0263 0.0266 0.0709 0.0965 0.0428 0.0343 0.0184 0.0227
dd =
2.1214 2.3574 3.5911 1.7227 2.4851 2.6968 2.3246 2.8107 8.9279 1.8828 1.8188
:victory: :victory: :victory:
[ 本帖最后由 xinyuxf 于 2007-5-8 09:23 编辑 ] 自导自演! 5.8就要交导师, 5.1长假期间，请教了周围的朋友，没有人会MATLAB
昨天被逼急了，上网找，搜索到了这个论坛，20点发贴，也没人帮忙
只好自己翻这个论坛的贴，足足翻了60页，几乎是一年的全部帖子了，
发现自己一开始用布鲁顿是错误的，改用fsolve，就OK了
不过MATLAB的基本的命令都不会，就是看着别人的程序，找葫芦画瓢
一个人，就这么摸索着，在漫漫黑夜里。。真的是黑夜啊，我昨夜通宵
终于自己解决了，难道不高兴么？ 楼上的说什么自编自演。呵呵。。
10楼的真狭隘。不过我一点也不生气。论文完成了，真是爽啊！
[ 本帖最后由 ChaChing 于 2010-3-1 08:31 编辑 ] 评分分数: 振动币 +14 个
操作理由: 虽然你的自说自话让人觉得...
＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
我自説自话，是因为半夜三更的，没人来帮助我，我只好自己记录下来摸索的痕迹，万一不成功，白天有高手看见了可以指点一下。
例如，我现在还不知道如何导入大量的数据，只能事先手工把变量先输入M文件里。 例如，解出来的X1/X2/DD，我预想的是竖列的，结果答案里是横排的，拷贝回EXCEL出错，超过255列了，我不知道转置的命令。只好一批批的200个200个拷贝进EXCEL，呵呵，自己汗一个。。。。 你摸索问题的精神很可嘉，如果你有问题或者找到了答案想要分享，都可以发在论坛上，这也正是论坛所提倡的。
但是，我认为在论坛上发帖、回帖，内容应该都是你所要表达的思想的关键、要点，应该是比较简要、易懂的，而不是把论坛当作你解决问题这一过程的详细的记录本。你也提到了，自己看了60页的网页，肯定很累很繁琐，你可以换个角度思考，试想如果大家都像你这样发帖子，像是个演草本，把所有东西都罗列上去，其他人浏览帖子也会不舒服的。
我想其他人说的话，没有什么恶意，只是可能表达方式不对。尤其是xjzuo，都已经给你加分，是对你努力的肯定，对你发帖存在的问题提出意见也是应该的，他是本版的版主，有责任指正在本版发帖会员的不当。 用1stOpt求解：
代码：
Constant d=5,ve=3,ro=0.4;
Function
(x1*normcdf((log(x1/d)+0.02+0.5*x2*x2)/x2)-exp(-0.02)*d*normcdf((log(x1/d)+0.02-0.5*x2*x2)/x2))-ve;
(x1*x2*normcdf((log(x1/d)+0.02+0.5*x2*x2)/x2))/ve-ro;
结果：
x1: 8.00023237371198
x2: 0.161444269349853
页:
[1]
2