求解非线性方程组遇到的问题
这是一个比较复杂的四元非线性方程组m文件代码如下：
function =mf8(X)
rou=7;
kai=8.2;
niu=3;
you=0.1;
yeta=0.1;
n=0.2;
sigma=10;
theta=1.9;
a=0.3;
E=1;
Cj=1;
Ck=1;
Kj=1;
Kk=1;
Ij=you*Kj;
Ik=you*Kk;
MCh=((X(1)/a)^a)*(((1-a)/X(3))^(a-1));
MCf=((X(2)/a)^a)*(((1-a)/X(4))^(a-1));
ph=sigma*MCh/(sigma-1);
phx=sigma*MCh/((sigma-1)*E);
pfx=sigma*MCf/(sigma-1);
pf=E*pfx;
PH=ph;
PF=pf;
PHx=phx;
PFx=pfx;
P=(yeta*PH^(1-theta)+(1-yeta)*PF^(1-theta))^(1/(1-theta));
Px=(yeta*PHx^(1-theta)+(1-yeta)*PFx^(1-theta))^(1/(1-theta));
STj=niu*P*Cj^rou/X(3);
STk=niu*Px*Ck^rou/X(4);
Yh=(Kj^a)*(STj^(1-a));
Yf=(Kk^a)*(STk^(1-a));
yh=yeta*(ph/PH)^(-sigma)*(PH/P)^(-theta)*(n*Cj+n*Ij)/n;
yf=(1-yeta)*(pf/PF)^(-sigma)*(PF/P)^(-theta)*(n*Cj+n*Ij)/(1-n);
yhx=yeta*(phx/PHx)^(-sigma)*(PHx/Px)^(-theta)*((1-n)*Ck+(1-n)*Ik)/n;
yfx=(1-yeta)*(pfx/PFx)^(-sigma)*(PFx/Px)^(-theta)*((1-n)*Ck+(1-n)*Ik)/(1-n);
x1=Yh-yh-yhx;
x2=Yf-yf-yfx;
x3=yh*ph+yhx*phx-P*Cj-P*Ij;
x4=yf*pf+yfx*pfx-Px*Ck-Px*Ik;
end
%上述函数是一组非线性方程，X为自变量，包含4个元素：X(1),X(2),X(3),X(4)。我想用寻找到这样的X，使得最后四个方程中x1,x2,x3,x4等于0
%求解的代码：
options=optimset('Display','iter','MaxIter',200000, 'NonlEqnAlgorithm', 'gn','MaxFunEvals',10000000);
=fsolve(@mf8,,options);
%反复测试后用了上面这一组初值
返回的结果：
Directional
IterationFunc-count Residual Step-size derivative
0 5 22.1216
1 15 3.53439 4.47 -1.77
Iteration matrix ill-conditioned - Switching to LM method.
2 23 1.7281 0.235 -11.5 0.100634
3 31 0.0453446 0.806 -1.3 2.93274
4 39 2.82673e-006 0.906 -0.00088 1.44781
5 47 1.60298e-019 1.01 -1.34e-012 0.757026
Optimization terminated: directional derivative along
search direction less than TolFun and infinity-norm of
gradient less than 10*(TolFun+TolX).
X =
0.1304 2.9124 8.4168 8.9318
FVAL =
4.0037e-010
EXITFLAG =
1
OUTPUT =
iterations: 6
funcCount: 47
stepsize: 1
cgiterations: []
firstorderopt: []
algorithm: 'medium-scale: Levenberg-Marquardt, line-search'
message:
%根据返回结果看来，应该找到解了，因为EXITFLAG = 1
%但是把求得的解带入函数：
=mf8();
%得到的结果是：
x1 =
-3.8486e-004
x2 =
1.6040
x3 =
0.6598
x4 =
-0.1650
%x1,x2,x3,x4都不等于0。为什么会出现这样矛盾的情况呢？另外，如果把初值略作改变，比如用，也可以得到另一组解：，同样地，带入函数，并不能使函数值为0
求高手帮助，谢谢！ 为了看看我对fsolve的理解有没有问题，我用了下面这个简单的方程组做个试验：
m函数：
function =xxx(X)
y1=X(1)-X(2);
y2=X(1)^2-X(2);
用下面的代码求解
=fsolve(@xxx,,options)
结果返回了
>> =fsolve(@xxx,,options)
Directional
IterationFunc-count Residual Step-size derivative
0 3 0
Optimization terminated: directional derivative along
search direction less than TolFun and infinity-norm of
gradient less than 10*(TolFun+TolX).
X =
0.1000 0.1000
FVAL =
0
EXITFLAG =
1
OUTPUT =
iterations: 1
funcCount: 3
stepsize: 1
cgiterations: []
firstorderopt: []
algorithm: 'medium-scale: Gauss-Newton, line-search'
message:
很显然，两个正确的解应该是,。这么简单的方程，都出错了，fsolve应该不至于这么弱智吧。是我对fsolve的理解和使用有问题吗？ 建议你重新装一下matlab,我也用matlab解方程，之所以不为零，是因为迭代法判定结果收敛到所要求的精度，就不再进行迭代了，所以结果不为零，第二个问题是迭代法最还给定的初值要接近真实值这样的结果才具有可靠性。
个人建议，仅供参考 试下1stOpt，不需要初值。 谢谢楼上两位！
我来找个matlab新版本试试看，是不是能解决。
1stopt我之前也用过，解相对比较简单的方程组没有问题，比matlab厉害，像我这样的复杂一点的方程组好像就不行了，每一次都是不同的解，代入结果测试，也不为0。
而且我故意用一个无解的方程组测试，1stopt还能像模像样地算出“结果”来，每次也给出不同的解。
如果matlab新版本也解决不了这个问题，那我该怎么办呢。难道说可以认定这个方程组无解？
自己再顶顶 本帖最后由 ChaChing 于 2012-7-10 23:19 编辑
回复 2 # zijiu 的帖子
当然fsolve不至于这么弱智, 一定是使用错误造成!:@)
首先, 需承认个人时间有限&懒惰, 1F的内容太长不太想细看!
就仅试试2F的东东, 我想问题原因应该一个样
同样的, 需请LZ再看下fsolve的help, 我也是仅能如此 function ff=xxx(X)
ff=[X(1)-X(2);
X(1)^2-X(2)];=fsolve(@xxx,,options)
=fsolve(@xxx,,options) 回复 2 # zijiu 的帖子
LZ知道为何错误了吗?
请将xxx.m改为下式再试试function y1=xxx(X)
y1=X(1)-X(2);
y2=X(1)^2-X(2);是不是与LZ原来情况相同?
原因是仅y1被求解!
请LZ细看下fsolve的help, 裡头有例子:@) 看看问题
页:
[1]