求解非线性方程组

huangwei_vicky

请大侠们帮忙看看，这是我的mfile

function y=myfun(x)

a=0.020;
b=0.024;
r=0.001;
aa=35;
L=0.300;
n=16;
S=L/sind(aa);

R=b+2*r;
At=pi*(b^2-a^2);
Ah=pi*r^2;
Eh=206000000000;
Gh=79000000000;
Et=7840000;
v=0.47;
I=pi*r^4/4;
J=I*2;
F=500;

%12个方程组解12个未知数
y(1)=F-x(1)-n*x(5);
y(2)=2*pi*b*L*x(2)-n*S*x(6);
y(3)=x(3)-x(7);
y(4)=x(4)-x(8);
y(5)=x(3)+b*v/Et/At*x(1)-b/Et*((a^2+b^2)/(b^2-a^2)-v)*x(2);
y(6)=x(4)-L/Et/At*x(1)-2*b^2*v*L/Et/(b^2-a^2)*x(2);

y(7)=x(6)-Eh*I*(cosd(x(10))^2/x(12)-cosd(aa)^2/R)*sind(x(10))^2*cosd(x(10))^2/x(12)^2+Gh*J*(sind(x(10))*cosd(x(10))/x(12)-sind(aa)*cosd(aa)/R)*sind(x(10))*cosd(x(10))^3/x(12)^2-x(9)*cosd(x(10))^2/x(12);

y(8)=x(5)-x(9)*sind(x(10))+Eh*I*(cosd(x(10))^2/x(12)-cosd(aa)^2/R)*sind(x(10))*cosd(x(10))^2/x(12)-Gh*J*(sind(x(10))*cosd(x(10))/x(12)-sind(aa)*cosd(aa)/R)*cosd(x(10))^3/x(12);

y(9)=x(12)-R*(1+x(9)/Eh/Ah)*sind(x(10))/sind(aa)*tand(aa)/(1+x(11)/2/pi)/tand(x(10));
y(10)=x(8)-L*((1+x(9)/Eh/Ah)*sind(x(10))/sind(aa)-1);
y(11)=x(11)-S*((1+x(9)/Eh/Ah)*cosd(x(10))/x(12)-cosd(aa)/R);
y(12)=x(12)-R;



调用函数：
x=fsolve(@myfun,[300; 10000; 0.005; 0.02; 500; 1000; 0.005; 0.02; 50; 35; 0; 0.026],optimset('Display','iter'))

算到92步停止，结果为
x =

  1.0e+006 *

   -0.0068
    3.8905
    0.0000
   -0.0000
    0.0005
    0.0212
    0.0000
   -0.0000
    0.0009
    0.0000
    0.0000
    0.0000

这里的1.0e+006 *是不是说明结果局部收敛就停止了？
除了调整初值还有没有别的方法可以解决这个问题呢？

小妹在此先谢过啦！！！

messenger

没什么好方法吧，12个方程，如果漫无目的地找理想解比较困难

huangwei_vicky

Solver stopped prematurely.

fsolve stopped because it exceeded the function evaluation limit,
options.MaxFunEvals = 1200 (the default value).

这是提示:dizzy:

shamohu

下面是用1stopt找的一组结果：

x1        191.131322330841
x2        74601.0177559414
x3        0.000661733604602548
x4        0.0219836958896179
x5        19.3042923542118
x6        403.280461304068
x7        0.000661733601993253
x8        0.0219581085629029
x9        22.0798887160174
x10        37.987438338195
x11        -0.177035915006289
x12        0.0252883935623459

messenger

这么重要的提示在问问题时就应该写出来。如果是MaxFunEvals的问题，试试用opiset改一下。

Solver stopped prematurely.

fsolve stopped because it exceeded the function evaluation limit,
options.MaxFunEvals = 1200 (the default value).

这是提示:dizzy:
huangwei_vicky 发表于 2011-6-27 23:40

huangwei_vicky

4# shamohu

感谢shamohu版主，原来是我式子有点小问题，哈哈，仍然非常感谢

huangwei_vicky

5# messenger

真正的问题是我式子的问题，不然在我给定的初值附近应该在超限之前就能得出答案的，不过也试着改了MaxFunEvals，学会了以后说不定也会用上

lin2009

式子有什么样问题，表达式有误？
建议在“致谢”帖子中贴出正确的式子。

huangwei_vicky

y(5)=x(3)+b*v/Et/At*x(1)+b/Et*((a^2+b^2)/(b^2-a^2)-v)*x(2);  %符号错了
y(12)=x(12)-R-x(7); %少了一项
修改过后16步就能求解出结果

qianxunzn

请问 如果我的方程在部分区域没有实数解，这样的关系曲线可以用matlab作出来不啊