}
cout<<"\n";
}
system("pause");
return0;
}
Ø程序运行结果如下:
三、注释最佳再现给定运动规律连杆机构优化设计
问题模型子程序
Ø作业要求
对给定的最佳再现给定运动规律连杆机构优化设计问题模型子程序(FORTRAN语言)进行必要的注释。
Ø注释结果
连杆机构问题函数子程序
目标函数==============
SUBROUTINEFFX(N,X,FX)/定义目标函数(子函数)FFX/
======================
DIMENSIONX(N)/定义一维数组X/
COMMON/ONE/I1,I2,I3,I4,NFX,I6/开辟公共区域ONE,其中有整型变量I1、I2、I3、I4、I6,实型变量NFX/
NFX=NFX+1
P0=ACOS(((1.0+X
(1))**2-X
(2)**2+25.0)/(10.0*(1.0+X
(1))))
Q0=ACOS(((1.0+X
(1))*+*2-X
(2)**2-25.0)/(10.0*X
(2)))
T=90.0*3.1415926/(180.0*30.0)将输入角30等份
FX=0.0
DO10K=0,30/DO语句实现循环功能,循环初值为0,终值为30/
PI=P0+K*T
QE=Q0+2.0*(PI-P0)**2/(3.0*3.1415926)
D=SQRT(26.0-10.0*COS(PI))
AL=ACOS((D*D+X
(2)*X
(2)-X
(1)*X
(1))/(2.0*D*X
(2)))
BT=ACOS((D*D+24.0)/(10.0*D))
IF(PI.GE.0.0.AND.PI.LT.3.1415926)THEN/IF条件语句:
若PI>0‘且’PI<π,则执行件/
QI=3.1415926-AL-BT/THEN块,条件满足选择此句/
ELSE
QI=3.1415926-AL+BT/ELSE块,条件不满足选择此句/
ENDIF/结束IF选择结构/
IF(K.NE.0.OR.k.NE.30)THEN
FX=FX+(QI-QE)**2*T
ELSE
FX=FX+(QI-QE)**2*T/2.0
ENDIF
10CONTINUE
END/程序停止运行/
不等约束=================
SUBROUTINEGGX(N,KG,X,GX)/计算不等式约束函数子程序GGX/
=========================
DIMENSIONX(N),GX(KG)/一维数组X、GX/
GX
(1)=1.0-X
(1)
GX
(2)=1.0-X
(2)
GX(3)=1.0-5.0
GX(4)=(1.0+X
(1))-(X
(2)+5.0)
GX(5)=(1.0+X
(2))-(X
(1)+5.0)
GX(6)=(1.0+5.0)-(X
(1)+X
(2))
GX(7)=-(1.4142*X
(1)*X
(2)-X
(1)**2-X
(2)**2)-16.0
GX(8)=-(X
(1)**2+X
(2)**2+1.4142*X
(1)*X
(2))+36.0/不等式约束条件GX
(1)、GX
(2)、GX(3)、GX(4)、GX(5)、GX(6)、GX(7)、GX(8)/
END/使程序停止运行/
等式约束=================
SUBROUTINEHHX(N,KH,X,HX)/计算等式约束函数子程序HHX/
=========================
DIMENSIONX(N),HX(KH)/定义一维数组X、HX/
X
(1)=X
(1)
END/使程序停止运行/
四、连杆机构问题+其他工程优化问题
1、连杆机构问题
设计一曲柄连杆摇杆机构,要求曲柄
从
转到
时,摇杆
的转角最佳再现已知的运动规律:
且
,
,
为极位角,其传动角允许在
范围内变化。
解:
1、设计变量的确定
已知
,根据余弦定理有:
因此,只有
为独立变量,则设计变量为
2、目标函数的建立
目标函数可根据已知的运动规律与机构实际运动规律之间的偏差最小为指标来建立,先假设在
和
间进行30等分,则可得目标函数的表达式:
式中:
3、约束条件的确定
曲柄存在的条件:
传动角要求:
4、优化计算结果
运用MATLAB编程求解:
●定义目标函数M文件
functionf=myobj(x)
f=0;%f为目标函数
Ao=acos(((1+x
(1))^2-x
(2)^2+25)/(10*(1+x
(1))));
%Ao--初始极位角,x
(1)--L2,x
(2)--L3
Bo=acos(((1+x
(1))^2-x
(2)^2-25)/(10*x
(2)));%Bo--初始位置角
forAi=Ao:
(pi/2)/30:
(Ao+pi/2)
%极位角Ai从Ao开始,分100等分,一直转到Ao+90°的角度
Bei=Bo+2*(Ai-Ao)^2/(3*pi);
%Bei为理想情况下曲柄每转过一个角度时,摇杆转过的角度
r=sqrt(26-10*cos(Ai));
m=acos((r^2+x
(2)^2-x
(1)^2)/(2*r*x
(2)));
n=acos((r^2+24)/(10*r));
if(Ai>=0&&Ai<=pi)
Bi=pi-m-n;
else
Bi=pi-m+n;
end
f=f+(Bi-Bei)^2;%Bi为曲柄每转过一个角度时,摇杆实际转过的角度
end
●定义非线性约束条件M文件
function[c,ceq]=mycon(x)
c=[acos((x
(1)^2+x
(2)^2-36)/(2*x
(1)*x
(2)))-135*(pi/180);
45*(pi/180)-acos((x
(1)^2+x
(2)^2-16)/(2*x
(1)*x
(2)))];
%c为非线性不等式约束(传动角要求)
ceq=[];%ceq为非线性等式约束(无)
●调用主程序
x0=[6,3];
A=[-10;0-1;-1-1;1-1;-11];
b=[-1;-1;-6;4;4];
[x,fval]=fmincon(@myobj,x0,A,b,[],[],[],[],@nonlcon)
程序运行结果如下:
通过MATLAB编程求解得到该问题的最优解为
2、其他工程问题:
运装计划:
某航空公司的运输机分前、后舱两部分装运客货,前舱容积为
,最大装载量为
,后舱容积为
,最大装载量为
。
今有客货两种,其单位体积及总重量如下表所示,表中还列出两种客货的运输单价。
装载时要求前后舱的载重量保持在
的比例,若货物一次装不完,剩下的货物可随其他客货第二次运出。
试为公司安排一个装货计划,使该次航班的收益最大。
货品
总重量/t
单位体积/(
)
运费/(元/t)
甲
20
20
200
乙
12
40
300
建立数学模型:
设a、b分别表示装载甲、乙货品的重量,下标1、2分别表示运输机的前、后舱。
根据题意,此问题的目标函数为:
按规定,运输机对装载重量和体积有要求,所以有下列约束条件:
装载货物的总重量不可能超过所需运输的货物总重量,且货物重量不能为负值,故:
同时,运输机要求前、后舱的载重量保持1:
1.5的比例:
转化为MATLAB的数学模型为:
运用MATLAB求解代码如下:
clear,clc
f=[-200;-200;-300;-300];
A=[1100;0011;1010;0101;200400;020040];
b=[20;12;10;15;160;320];
Aeq=[1.5-11.5-1];
beq=[0];
lb=zeros(4,1);
[x,fval]=linprog(f,A,b,Aeq,beq,lb)
求解结果如下:
结果显示:
当
时,可使该航班的收益最大,最大收益为4400元
五、课程实践心得体会