基于MATLAB的曲柄摇杆机构优化设计说明

课程作业曲柄摇杆优化设计姓名：XX学号：XXXXX班级：XXXXX XX大学机械与动力学院目录1摘要2问题研究2.1问题重述2.2问题分析3数学模型的建立3.1设计变量的确定3.2目标函数的建立3.3约束条件的确定3.4标准数学模型4使用MATLAB编程求解4.1调用功能函数4.2首先编写目标函数M 文件4.3编写非线性约束函数M 文件4.4编写非线性约束函数M 文件confun.m4.5运行结果5结果分析6结论推广7过程反思8个人小结9参考文献1摘要: 为分析机构能够满足给定的运动规律和运动空间的要求,运用Matlab 22.10(32πψψ+=式中0ϕ和0ψ得小于45=≥][min γγ1可空间，可以适当预选机架杆的长度，现取l 4 =5。

2.2 问题分析设计时，可在给定最大和最小传动角的前提下，当曲柄从0ϕ转到090ϕ︒+时，要求摇杆的输出角最优地实现一个给定的运动规律()f ϕ。

这里假设要求：()()20023E f φϕφϕϕπ==+- （1）图1 曲柄摇杆机构简图对于这样的设计问题，可以取机构的期望输出角()E f φϕ=和实际输出角()F φϕ=的平方误差之和作为目标函数，使得它的值达到最小。

在图 1 所示的曲柄摇杆机构中，1l 、2l 、3l 、 4l 分别是曲柄AB 、连杆BC 、摇杆CD 和机架AD 的长度。

这里规定0ϕ为摇杆在右极限位置0φ时的曲柄起始位置角，它们由1l 、2l 、3l 和4l 确定。

3数学模型的建立 3.1 设计变量的确定决定机构尺寸的各杆长度1l 、2l 、3l 和4l ，以及当摇杆按已知运动规律开始运行时，曲柄所处的位置角0ϕ应列为设计变量，所有设计变量有:[][]1234512340TTx x x x x x l l l l ϕ== （2）考虑到机构的杆长按比例变化时，不会改变其运动规律，通常设定曲柄长度1l =1.0，在这里可给定4l =5.0，其他杆长则按比例取为1l 的倍数。

若取曲柄的初始位置角为极位角，则ϕ及相应的摇杆l 位置角φ均为杆长的函数，其关系式为:()()()()2222212432301242125arccos 2101l l l l l l l l l l ϕ⎡⎤⎡⎤++-+-+==⎢⎥⎢⎥++⎢⎥⎢⎥⎣⎦⎣⎦ （3）()()222221243230343125arccos 210l l l l l l l l l φ⎡⎤⎡⎤+--+--==⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦（4）因此，只有2l 、3l 为独立变量，则设计变量为[][]1223T Tx x x l l ==。

3.2 目标函数的建立目标函数可根据已知-的运动规律与机构实际运动规律之间的偏差最小为指标来建立，即：()()21mi n mE i i i f x φφ==-→∑ （5）式中，Ei φ-期望输出角； m-输出角的等分数；i φ-实际输出角，由图 1 可知：图2 曲柄摇杆机构的运动学关系()()02i i i i i i i παβϕπφπαβπϕπ--≤≤⎧⎪=⎨-+≤≤⎪⎩ (6) 式中，222222322132arccos arccos 22i i i i i r l l r x x rl r x α⎛⎫⎛⎫+-+-== ⎪ ⎪⎝⎭⎝⎭ (7)222241424arccos arccos 210i i i i i r l l r rl r β⎛⎫⎛⎫+-+== ⎪ ⎪⎝⎭⎝⎭(8)i r == (9)3.3 约束条件 曲柄存在条件:12131423;,l l l l l l l l ≤≤+≤+ ()()24133412,l l l l l l l l ≤-+≤-+ 曲柄与机架共线位置时的传动角（连杆BC 和摇杆CD 之间的夹角）: 最小传动角min min 45r BCD ︒=∠≥ 最大传动角max max 135r BCD ︒=∠≤ 由上面的分析可以算出：()222222234112min231216arccos 4522l l l l x x r l l x x ︒⎡⎤+--⎡⎤+-⎢⎥==≥⎢⎥⎢⎥⎣⎦⎣⎦ （10）()222222234112max231236arccos 13522l l l l x x r l l x x ︒⎡⎤+-+⎡⎤+-⎢⎥==≤⎢⎥⎢⎥⎣⎦⎣⎦ （11）3.4 标准数学模型通过上面的分析后，将输入角分成 30 等分（m=30），经过转化为标准形式得到曲柄摇杆机构优化设计标准数学模型为:()()21min mEi i i f x φφ==-→∑[][]2312TTx l l x x ==()()()()()()()112231241252122612122271212101060..40401.41436036 1.4140g x x g x x g x x x s t g x x x g x x x g x x x x x g x x x x x =-≤⎧⎪=-≤⎪⎪=--≤⎪=--≤⎨⎪=--≤⎪=+--≤⎪⎪=---≤⎩ (12) 机械优化设计中的问题，大多数属于约束优化问题，此为非线性约束优化问题，运用 MATLAB 优化工具箱的命令函数 fmincon 来处理有约束的非线性多元函数最小化优化问题。

4使用MATLAB 编程求解4.1 本问题属于一般非线性规划问题，其标准型为：min ()f x,,()0..()0,AX b Aeq X beq C X s t Ceq X vlb X vub ≤∙=≤⎧⎨=≤≤⎩ （13）调用MATLAB 软件优化工具箱中非线性规划求解函数fmincon 来求解。

其命令的基本格式为： [函数] fmincon [格式]x = fmincon(fun,x0,A,b) x = fmincon(fun,x0,A,b,Aeq,beq) x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub) x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon) x = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options) [x,fval] = fmincon(…) [x,fval,exitflag] = fmincon(…) [x,fval,exitflag,output] = fmincon(…) [x,fval,exitflag,output,lambda] = fmincon(…) [x,fval,exitflag,output,lambda,grad] = fmincon(…) [x,fval,exitflag,out put,lambda,grad,hessian] = fmincon(…) [说明]fun 是目标函数options 设置优化选项参数fval 返回目标函数在最优解x点的函数值exitflag 返回算法的终止标志output 返回优化算法信息的一个数据结构grad 返回目标函数在最优解x点的梯度hessian 返回目标函数在最游解x点Hessian矩阵值编写程序求解4.2首先编写目标函数M 文件fun1.mfunction f=fun1(x)s=30;qb=1;jj=5;fx=0;ci0=acos(((qb+x(1))^2-x(2)^2+jj^2)/(2*(qb+x(1))*jj));%曲柄初始角fa0=acos(((qb+x(1))^2-x(2)^2-jj^2)/(2*x(2)*jj));%摇杆初始角for i=1:sci=ci0+(pi*i)/(2*s);fai(i)=fa0+(2*(ci-ci0)^2)/(3*pi);ri=sqrt(qb^2+jj^2-2*qb*jj*cos(ci));alfi=acos(((ri^2+x(2)^2)-x(1)^2)/(2*ri*x(2)));bati=acos((ri^2+jj^2-qb^2)/(2*ri*jj));if ci>0 && ci<=pipsi(i)=pi-alfi-bati;elseif ci>pi && ci<=2*pipsi(i)=pi-alfi+bati;endfx=fx+(fai(i)-psi(i))^2;endf=fx;i=1:1:30;plot(i,fai(i),i,psi(i),'--'); %画曲线图legend('期望曲线','实际曲线'); %标注曲线图对应名称4.3编写非线性约束函数M 文件confun.mfunction [c,ceq]=confun(x)qb=1;jj=5;m=45*pi/180;n=135*pi/180;c(1)=x(1)^2+x(2)^2-2*x(1)*x(2)*cos(m)-(jj-qb)^2;%重合时最小传动角的非线性约束条件c(2)=-x(1)^2-x(2)^2+2*x(1)*x(2)*cos(n)+(jj+qb)^2;%共线时最小传动角的非线性约束条件ceq=[];4.4在MATLAB 命令窗口调用优化程序x0=[6;6];lb=[1;1];ub=[];a=[-1 0;0 -1;-1 -1;1 -1; -1 1];b=[-1;-1;-6;4;4];options=optimset('LargeScale','off','display','iter');[x,fval,exitflag]=fmincon(@fun1,x0,a,b,[],[],lb,ub,@confun,options);4.5运行结果x =[4.1285 2.3226]fval =0.0076图3 输出角期望曲线与在MATLAB结果下的实际曲线对比图图4 传动角与曲柄输入角变化关系图5结果分析通过Matlab工具箱的优化求解，我们得到了最终的曲柄摇杆机构的最优杆长条件，即L2=4.1285，L3=2.3226。

从运行结果上面来看，得到的数据还是比较理想的，在输出角期望曲线与在MATLAB结果下的实际曲线对比图（图3）中，我们可以清楚地看到，期望曲线与实际曲线的拟合程度比较好。

在传动角6结论推广由于在本问题当中，曲柄长度L1和机架长度L4是预先取的L1=1，L4=5，我们通过对L2和L3的优化设计，最终得到了L2=4.1285，L3=2.3226，如果把1看作是单位长度，那么我们最终求解出来的其实是曲柄摇杆机构符合已知运动轨迹的杆长比例。

只要曲柄摇杆机构的四杆长度按照这个比例，即L1：L2：L3：L4=1:4.1285:2.3226:5,那么我们得到的曲柄摇杆机构的运动轨迹都是比较理想的。