我用fmincon做一个多变量约束优化:
代码如下,NonLinF为目标函数,NonLinCon为约束方程
function FrictionDriving
x0=[(80/180)*pi,0.004,0.046,0.003,0.0015];
lb=[(80/180)*pi,0.004,0.046,0.002,0.001];
ub=[(120/180)*pi,0.007,0.06,0.005,0.003];
options=optimset('LargeScale','off');
[x,f]=fmincon(@NonLinF,x0,[0 0 0 -1 1;0 0 6 0 0],[0;0.371],[],[],lb,ub,@NonLinCon,options)
结果:
"Optimization terminated: no feasible solution found. Magnitude of search
direction less than 2*options.TolX but constraints are not satisfied."
x =
1.3963 0.0040 0.0460 0.0030 0.0015
f =
1.6086e+004
有几个问题,1. 结果提示中“Magnitude of search direction less than 2*options.TolX”是什么意思?
2. 为什么优化结果中,设计变量x的值即为初始值x0?
function [C,Ceq]=NonLinCon(x)
% 本函数约束函数
% T0为单齿作用产生的转矩
% x=0-2*pi 摩擦齿轮的旋转角度,w是设计变量
% w(1)-切入角
% w(2)-齿高H
% w(3)-基圆半径
% w(4)-齿距P
% w(5)-齿根宽S
b=0.133;%摩擦齿轮的厚度
Lmax=1080; %机器人VOD在爬坡时(20度)的最大负载,单位为Nm
deltaF=194000000; %齿根许用弯曲应力,单位为Pa
w=0:0.01:2*pi;
i=-4:1:10;
c1=433900; %轮胎径向刚度,单位为N/m
%c2=w3; %摩擦齿轮基圆半径
%c3=w2; %齿顶高度
c4=0.371; %轮胎外径
c5=0.35; %损耗因子,来自参考文献“基于粗糙接触理论的橡胶-金属摩擦副的摩擦分析”
c6=0.4997; %橡胶材料泊松比,来自参考文献“轮胎胶料有限元分析的实验基础及计算”
c7=5000000; %橡胶材料弹性模量,单位Pa,来自参考文献“轮胎胶料有限元分析的实验基础及计算”
c8=3/8; %公式M中的常系数
c9=3/4; %公式M中的常系数
c10=2; %公式M中的常系数
c11=0.5; %公式M中的常系数
e1=2; %公式M中的常指数
e2=1/3; %公式M中的常指数
e3=4/3; %公式M中的常指数
e4=2/3; %公式M中的常指数
e5=0.5; %公式M中的常指数
a3=x(4)/(x(2)+x(3));%齿距p对应的弧度
w0=pi/2-acos(((x(3)+c4)^e1+(x(3)+x(2))^e1-c4^e1)/(c10*(x(3)+c4)*(x(3)+x(2))));%单齿旋转初始角
w1=pi-w0;%单齿旋转结束角度
%a1=w1;%轮齿的齿侧面和齿顶面的夹角,也是摩擦轮胎的切入角
%y1=zeros(size(x));
%y2=zeros(size(x));
%y3=zeros(size(x));
z=zeros(size(w));
N=length(w);
M=length(i);
A=zeros(M,N);
for J=1:M
for k=1:N
if w(k)<w0+i(J)*a3&w(k)>=0;
z(k)=0;
elseif w(k)>=w0+i(J)*a3&w(k)<=w1+i(J)*a3;
a2=atan(((x(3)+x(2))*cos(w(k)-i(J)*a3))./((x(3)+c4)-(x(3)+x(2))*sin(w(k)-i(J)*a3)));%论文图5-10中的角度beta
c12=((c10*(c4+x(3))-(x(3)+x(2))*sin(w(k)-i(J)*a3)-c4*cos(a2)).^e1+((x(3)+x(2))*cos(w(k)-i(J)*a3)-c4*sin(a2)).^e1).^e5;%轮胎的径向变形量
p1=c1*c11*c12.*sin(c10*(x(1)-w(k)+i(J)*a3+a2));
F1=c1*c12.*sin(x(1)-w(k)+i(J)*a3+a2);%作用在齿侧面S1上的法向力
p2=c8*c5*((c9*((1-c6)^e1)/c7)^e2);
p3=((c1*c12.*sin(w(k)-i(J)*a3-a2)).^e3)/((x(3)+x(2))^e4);
p4=cos(pi/2-w(k)+i(J)*a3+a2);
p5=c4-c12;
z(k)=-((p1+p2*p3.*p4).*p5);
else w(k)<=(2*pi)&w(k)>w1+i(J)*a3;
z(k)=0;
end
end
B(J,:)=z;
end
maxF1=max(F1);
y=B;
sumT=zeros(size(w));
for J=1:M
sumT=sumT+y(J,:);
end
WsumT=sumT.*W(x);
indxx=find(w>w0&w<w1);
Tsum=WsumT(indxx);
maxT=max(Tsum);
minT=min(Tsum);
C(1)=Lmax-maxT; %最大驱动转矩maxT要大于最大负载转矩Lmax
C(2)=(6*maxF1*x(2))./(b*x(5).^2)-deltaF; %齿根弯曲强度要求
C(3)=-minT;
Ceq=[];
function f=NonLinF(x)
% 本函数求解多齿作用下的转矩优化目标函数
% T0为单齿作用产生的转矩
% x=0-2*pi 摩擦齿轮的旋转角度,w是设计变量
% w(1)-切入角
% w(2)-齿高H
% w(3)-基圆半径
% w(4)-齿距P
% w(5)-齿根宽S
w=0:0.01:2*pi;
i=-4:1:10;
c1=433900; %轮胎径向刚度,单位为N/m
%c2=w3; %摩擦齿轮基圆半径
%c3=w2; %齿顶高度
c4=0.371; %轮胎外径
c5=0.35; %损耗因子,来自参考文献“基于粗糙接触理论的橡胶-金属摩擦副的摩擦分析”
c6=0.4997; %橡胶材料泊松比,来自参考文献“轮胎胶料有限元分析的实验基础及计算”
c7=5000000; %橡胶材料弹性模量,单位Pa,来自参考文献“轮胎胶料有限元分析的实验基础及计算”
c8=3/8; %公式M中的常系数
c9=3/4; %公式M中的常系数
c10=2; %公式M中的常系数
c11=0.5; %公式M中的常系数
e1=2; %公式M中的常指数
e2=1/3; %公式M中的常指数
e3=4/3; %公式M中的常指数
e4=2/3; %公式M中的常指数
e5=0.5; %公式M中的常指数
a3=x(4)/(x(2)+x(3));%齿距p对应的弧度
w0=pi/2-acos(((x(3)+c4)^e1+(x(3)+x(2))^e1-c4^e1)/(c10*(x(3)+c4)*(x(3)+x(2))));%单齿旋转初始角
w1=pi-w0;%单齿旋转结束角度
%a1=w1;%轮齿的齿侧面和齿顶面的夹角,也是摩擦轮胎的切入角
%y1=zeros(size(x));
%y2=zeros(size(x));
%y3=zeros(size(x));
z=zeros(size(w));
N=length(w);
M=length(i);
A=zeros(M,N);
for J=1:M
for k=1:N
if w(k)<w0+i(J)*a3&w(k)>=0;
z(k)=0;
elseif w(k)>=w0+i(J)*a3&w(k)<=w1+i(J)*a3;
a2=atan(((x(3)+x(2))*cos(w(k)-i(J)*a3))./((x(3)+c4)-(x(3)+x(2))*sin(w(k)-i(J)*a3)));%论文图5-10中的beta
c12=((c10*(c4+x(3))-(x(3)+x(2))*sin(w(k)-i(J)*a3)-c4*cos(a2)).^e1+((x(3)+x(2))*cos(w(k)-i(J)*a3)-c4*sin(a2)).^e1).^e5;%轮胎的径向变形量
p1=c1*c11*c12.*sin(c10*(x(1)-w(k)+i(J)*a3+a2));
p2=c8*c5*((c9*((1-c6)^e1)/c7)^e2);
p3=((c1*c12.*sin(w(k)-i(J)*a3-a2)).^e3)/((x(3)+x(2))^e4);
p4=cos(pi/2-w(k)+i(J)*a3+a2);
p5=c4-c12;
z(k)=-((p1+p2*p3.*p4).*p5);
else w(k)<=(2*pi)&w(k)>w1+i(J)*a3;
z(k)=0;
end
end
B(J,:)=z;
end
y=B;
sumT=zeros(size(w));
for J=1:M
sumT=sumT+y(J,:);
end
WsumT=sumT.*W(w);
indxx=find(w>w0&w<w1);
Tsum=WsumT(indxx);
minT=min(Tsum);
maxT=max(Tsum);
f=maxT-minT;
function [y]=W(x)
% 本函数为一窗口函数
% T0为单齿作用产生的转矩
% x=0-2*pi 摩擦齿轮的旋转角度,w是设计变量
%c2=w3; %摩擦齿轮基圆半径
%c3=w2; %齿顶高度
w=0:0.01:2*pi;
c4=0.371; %轮胎外径
c10=2; %公式M中的常系数
e1=2; %公式M中的常指数
w0=pi/2-acos(((x(3)+c4)^e1+(x(3)+x(2))^e1-c4^e1)/(c10*(x(3)+c4)*(w(3)+w(2))));%单齿旋转初始角
w1=pi-w0;%单齿旋转结束角度
%a1=w1;%轮齿的齿侧面和齿顶面的夹角,也是摩擦轮胎的切入角
z=zeros(size(w));
N=length(w);
for k=1:N
if w(k)<w0;
z(k)=0;
elseif w(k)>=w0&&w(k)<=w1;
z(k)=1;
else w(k)>w1;
z(k)=0;
end
end
y=z;
vBulletin® v3.8.3,版权所有 ©2000-2025,Jelsoft Enterprises Ltd.