Labfans是一个针对大学生、工程师和科研工作者的技术社区。 论坛首页 | 联系我们(Contact Us)
MATLAB爱好者论坛-LabFans.com
返回   MATLAB爱好者论坛-LabFans.com > 工程数学软件 > MATLAB论坛
MATLAB论坛 一切MATLAB相关问题在此讨论。
 
 
主题工具 显示模式
旧 2010-04-08, 23:17   #1
ncgcake
初级会员
 
注册日期: 2009-07-22
帖子: 1
声望力: 0
ncgcake 正向着好的方向发展
默认 求助:求最小值,程序老出错

目标函数:122.5221*sqrt(0.9025+x1^2)*x2 约束:(1.2732e+007*sqrt(0.9025+x1^2)/(x1*x2)-4.2e+008)<=0,(25908*(0.0025^2+x2^2)*x1*x2-1.2732*(0.9025+x1^2)^1.5)+1/(x2-0.005)+1/(x1)>=0
clc
m=zeros(1,100);a=zeros(1,100);b=zeros(1,100);f0=zeros(1,100);%a b为最优点坐标,f0为最优点函数值,f1 f2最优点梯度。
syms x1 x2 e; %e为罚因子。
m(1)=1;c=0.2;a(1)=0.75;b(1)=0.08; %c为递增系数,赋初值。
f=122.5221*sqrt(0.9025+x1^2)*x2+e*(1/(42*(x1*x2)-1.2732*sqrt(0.9025+x1^2))+1/(25908*(0.0025^2+x2^2)*x1*x2-1.2732*(0.9025+x1^2)^1.5)+1/(x2-0.005)+1/(x1));
%x1为高度h,x2为管直径D。
fx1=diff(f,'x1');fx2=diff(f,'x2');fx1x1=diff(fx1,'x1');fx1x2=diff(fx1,'x2');fx2x1=diff(fx2,'x1');fx2x2=diff(fx2,'x2');%求偏导、海森元素。
for k=1:5 %外点法e迭代循环.
x1=a(k);x2=b(k);e=m(k);
for n=1:1000 %牛顿法求最优值。
f1=subs(fx1); %求解梯度值和海森矩阵
f2=subs(fx2);
f11=subs(fx1x1);
f12=subs(fx1x2);
f21=subs(fx2x1);
f22=subs(fx2x2);
if(double(sqrt(f1^2+f2^2))<0.01) %最优值收敛条件
a(k+1)=double(x1);b(k+1)=double(x2);f0(k+1)=double(subs(f));
break;
else %如未最优值收敛条件,对x1,x2进行迭代,直到寻找到合适的迭代点
X=[x1 x2]-(([f1 f2]*[f1 f2]')/([f1 f2]*[f11 f12;f21 f22]*[f1 f2]'))*[f1 f2]';
x1=X(1,1);x2=X(2,1);
end
end
if(double(sqrt((a(k+1)-a(k))^2+(b(k+1)-b(k))^2))<=0.01)&&(double(abs((f0(k+1)-f0(k))/f0(k)))<=0.01) %罚因子迭代收敛条件
h=a(k+1) %输出最优点h,D值,罚因子迭代次数,最优值
D=b(k+1)
times=k
zhiliang=f0(k+1)
break;
else
m(k+1)=c*m(k);

h=a(k+1) %输出最优点h,D值,罚因子迭代次数,最优值
D=b(k+1)
times=k
zhiliang=f0(k+1)
end
end
有什么错误,请详细指出来,谢谢了

此帖于 2010-04-08 23:20 被 ncgcake 编辑。
ncgcake 当前离线   回复时引用此帖
 


发帖规则
不可以发表新主题
不可以发表回复
不可以上传附件
不可以编辑自己的帖子

启用 BB 代码
论坛启用 表情符号
论坛启用 [IMG] 代码
论坛禁用 HTML 代码



所有时间均为北京时间。现在的时间是 09:04


Powered by vBulletin
版权所有 ©2000 - 2025,Jelsoft Enterprises Ltd.