面对资源日益紧缺的趋势，机械结构的轻量化设计显得尤为重要. 使用新型轻量化材料、开发先进的成形方法、合理的结构优化设计是实现轻量化的主要途径^[1].

叉车是工业运输的重要工具^[2]，门架是叉车主要的承载机构，通过不断增加立柱厚度来提高门架强度和刚度的现象在国内一些叉车企业普遍存在，造成材料的严重浪费^[3-4].

导重法是一种实用的机械结构优化设计方法，由陈树勋等^[5]提出. 相关学者将该方法运用到实际工程的优化设计中，取得了显著效果，大大提高了优化效率和设计质量. 李枝东等^[6]结合导重法与SIMP模型，成功求解了单工况拓扑优化问题；刘辛军等^[7]针对多工况优化问题，推出导重法的迭代准则，采用对偶法改进了拉格朗日乘子的求法；陈垂福等^[8]对导重法步长因子进行深入研究，提出2种变步长因子的控制策略，将密度补偿法引入导重法中，提高了求解效率；针对导重法在求解拓扑优化模型中需要多次迭代才能得到收敛解的问题，张昆鹏等^[9]提出基于向量Epsilon算法加速迭代序列收敛的方法，大大减少了迭代次数. 利用导重法进行结构优化的关键在于容重与导重的求解，陈树勋等^[10]提出将ANSYS作为求解器，利用ANSYS的梯度法，求解设计变量差分灵敏度的方法.

本文提出结合带有经典误差修正项的灵敏度求解半解析法ESA（exact semi-analytical technique） ^[11-12]与导重法的集成优化法，利用该方法进行某型号内燃叉车的轻量化设计，推导强度、刚度等多性态约束，建立单目标多性态约束优化模型，解决了导重法容重和导重难以求解的问题，提高了求解效率，实现了门架的轻量化. 采用有限元法对优化后的门架进行强度和刚度的校核，验证了轻量化设计结果的合理性.

一个完整的两级叉车门架系统包括货叉、货叉架、内门架、外门架、起升链条、起升液压缸及倾斜液压缸等结构，如图1所示.

内门架和外门架是门架系统的主要构件，采用C-L型组合方式，如图2所示，内门架以外门架为活动导轨，通过滚轮接触来实现伸缩运动.

门架立柱主体采用C型槽钢制成，材料为Q345，具体的材料属性见表1. 表中，E为弹性模量，ρ为密度，σ_s为屈服强度.

传统伸缩式门架安全系数较高，强度和刚度都有很大富余，存在很大的设计余量. 采用导重法，对门架进行轻量化设计. 门架的轻量化主要以内外门架立柱为主，机械性能包括强度、刚度，所以优化过程不是单一性态约束，而是多性态约束的优化. 将门架立柱截面尺寸作为优化设计变量，立柱质量作为优化目标量，强度和刚度作为约束条件，运用导重法展开轻量化设计.

门架立柱受载时可以简化为悬臂梁，上、下翼缘受到滚轮的压力作用，相当于受到一个集中力和一个力偶的作用. 产生的弯曲正应力远大于切应力，可以将弯曲正应力看作立柱的主应力，从而简化为平面应力状态，如图3所示.

门架立柱受载时危险截面应为翼缘与腹板连接处的根部，翼缘受力状态为单向应力状态，包含3项应力：由集中力造成的整体弯曲应力σ_P、由力偶造成的扭转应力σ_M、翼缘与腹板的局部弯曲应力σ_W.

门架立柱翼缘合成应力为

(1) $ \sigma = {\sigma _{\rm{P}}} + {\sigma _{\rm{M}}} + {\sigma _{\rm{W}}}. $

需要满足的强度条件为

(2) $ \sigma \leqslant {{{\sigma _{\rm{s}}}}}/{n}. $

式中：σ_s为材料的屈服极限，n为安全系数.

危险截面的当量应力σ_d应按照第四强度理论计算：

(3) $ {\sigma _{\rm{d}}} = \sqrt {{\sigma ^2} + \sigma _x^2 - {\sigma _x}\sigma } \leqslant {{{\sigma _{\rm{s}}}}}/{n}. $

式中：σ为翼缘合成应力，σ_x为局部弯曲应力在x方向的分力.

由于σ_x总小于σ，根据式（1）、（3）可以看出，当量应力总小于合成应力，危险点总出现在翼缘的自由边.

门架的刚度是指门架在受载时的弹性变形，包括门架的转角变形及水平位移. 门架的刚度不足直接影响叉车的堆垛作业及稳定性，若水平位移量过大，则会使整车重心前移，容易造成叉车纵向失稳^[13].

叉车门架的刚度计算工况如下：门架处于直立状态，满载3 t的重物上升至3 m高的位置，起升重物的重心位于载荷中心.

从图4可以看出，门架总体顶端水平位移可以表示为

(4) $ f = {f_0} + {f_1} + {f_2}. $

式中：f为门架总体顶端水平位移，f₀为外门架顶端水平位移，f₁为由外门架转角引起的水平位移，f₂为内门架顶端水平位移.

根据弯矩图乘法^[14]，可以得到各个位移的表达式：

(5) $ \left. {\begin{aligned} & {{f_0} = [{{{M_0}}}/({{E{J_0}}})]\left( {{{H_0^2}}/{2} - {2}{H_0}{h_0}/{3} + {{h_0^2}}/{6} - {{h_1^2}}/{6}} \right)},\\ & {{f_1} = \theta \left( {{H_1} - {h_1}} \right)},\\ & {{f_2} = [{{{M_1}}}\;({{E{J_1}}})]\left( {{{H_1^2}}/{2} - {2}{H_1}{h_1}/{3} + {{h_1^2}}/{6} - {{h_2^2}}/{6}} \right)}. \end{aligned}} \right\} $

式中：M₀为外门架受到滚轮压力产生的弯矩，M₁为内门架受到的滚轮压力产生的弯矩，E为材料的弹性模量，J₀、J₁分别为外门架和内门架的惯性矩.

根据门架刚度要求可知，门架总体顶端水平位移应控制在一定范围内，即f≤[f]，其中 $\left[ f \right] = {H}/{{100}}$.

针对机械结构的优化设计，目前主要有数学规划法和准则法两大类. 数学规划法适用范围广且较严谨，依据设计点函数的局部性状来决定迭代寻优方向，需要多次迭代才能找到收敛解，优化效率低，难以满足复杂工程实际的需要. 结构优化准则法需要预先给定优化准则，因此迭代过程可以“定向”进行，具有目的明确、求解速度快和计算量小的优点，但是在数学上不够严谨，往往找不到优化设计的最优解.

导重法是在优化设计过程中利用“导重”在机械结构空间内将结构材料重量在构件间进行合理分配的新型结构优化法，它以“各设计变量对应的重量与导重成正比”为导重准则，使结构材料重量分配达到以导重准则衡量的最优结构状态. 采用导重法与灵敏度计算相结合的集成优化法来进行门架的轻量化设计，该方法同时具有数学规划法适用范围广及准则法求解速度快的优点.

以外门架为例建立优化模型，外门架立柱截面如图5所示.

以门架立柱上下翼缘厚度、腹板厚度、立柱宽度为优化变量，以立柱质量为优化目标，以强度和刚度为性态约束，采用导重法进行门架的优化设计.

叉车门架槽钢多约束质量优化问题的数学模型可以描述为

(6) $ \left. {\begin{array}{*{20}{c}} {{\rm{find}}\;{{X}} = \left[ {{x_1},{x_2},{x_3},{x_4}} \right]}^{\rm{T}},\\ {\min M\left( {{X}} \right)}.\\ {{\rm{s}}{\rm{.t}}.\;{g_1}\left( {{X}} \right) \leqslant 0}, {{g_2}\left( {{X}} \right) \leqslant 0},\\ {{x_{i{\rm{min}}}} \leqslant {x_i} \leqslant {x_{i{\rm{max}}}},i = 1,2,3,4}. \end{array}} \right\} $

式中：X为槽钢截面的尺寸变量；M（X）为目标函数，表征槽钢的质量；g₁（X）和g₂（X）为约束函数，分别表示强度、刚度两种约束；x_imin和x_imax分别为变量的上限和下限.

对于式（6）表示的数学模型，采用拉格朗日乘数法求解，构造拉格朗日函数：

(7) $ L = M\left( {{X}} \right) + {\lambda _1}{g_1}\left( {{X}} \right) + {\lambda _2}{g_2}\left( {{X}} \right). $

式中：λ₁、λ₂分别为g₁（X）和g₂（X）所对应的拉格朗日乘子.

针对优化参数X取得最优值X^*时的情况，当x_i=x_imin时，优化参数的下限起约束作用；当x_i=x_imax时，优化参数的上限起约束作用. 基于Kuhn-Tucker条件^[15]，在最优值X^*处必须满足以下条件：

(8) $\left. \begin{aligned} & \frac{{\partial L}}{{\partial {x_i}}} = \left\{\begin{array}{l}\dfrac{{\partial M}}{{\partial {x_i}}} + {\lambda _1}\dfrac{{\partial {g_1}}}{{\partial {x_i}}} + {\lambda _2}\dfrac{{\partial {g_2}}}{{\partial {x_i}}}\leqslant 0,\; x_i^* = {x_{i{\rm{max}}}};\\ \dfrac{{\partial M}}{{\partial {x_i}}} + {\lambda _1}\dfrac{{\partial {g_1}}}{{\partial {x_i}}} + {\lambda _2}\dfrac{{\partial {g_2}}}{{\partial {x_i}}}= 0,\; {x_{i{\rm{min}}}} \leqslant x_i^* \leqslant {x_{i{\rm{max}}}};\\ \dfrac{{\partial M}}{{\partial {x_i}}} + {\lambda _1}\dfrac{{\partial {g_1}}}{{\partial {x_i}}} + {\lambda _2}\dfrac{{\partial {g_2}}}{{\partial {x_i}}}\geqslant 0,\; x_i^* = {x_{i{\rm{min}}}};\end{array}\right.\\& {\lambda _j} \geqslant 0,j = 1,2;\\ & {\lambda _j}{g_j}\left( {{X}} \right) = 0,\;{{\lambda _j} = 0,\;{g_j}\left( {{X}} \right) < 0};\\ & {\lambda _j}{g_j}\left( {{X}} \right) = 0,\;{{\lambda _j} \geqslant 0,\;{g_j}\left( {{X}} \right) = 0}; \\ & {g_j}\left( {{X}} \right),j = 1,2;\\ & {x_{i{\rm{min}}}} \leqslant {x_i^*} \leqslant {x_{i{\rm{max}}}},i = 1,2,3,4. \end{aligned} \!\!\! \right\} $

根据 ${{\partial L}}/{{\partial {x_i}}} = {{\partial M}}/{{\partial {x_i}}} + {\lambda _1}{{\partial {g_1}}}/{{\partial {x_i}}} + {\lambda _2}{{\partial {g_2}}}/{{\partial {x_i}}} = 0$，可得

(9) $ \frac{{ - {\lambda _1}\dfrac{{\partial {g_1}}}{{\partial {x_i}}} - {\lambda _2}\dfrac{{\partial {g_2}}}{{\partial {x_i}}}}}{{\dfrac{{\partial M}}{{\partial {x_i}}}}} = \dfrac{{ - {\lambda _1}{x_i}\dfrac{{\partial {g_1}}}{{\partial {x_i}}} - {\lambda _2}{x_i}\dfrac{{\partial {g_2}}}{{\partial {x_i}}}}}{{{x_i}\dfrac{{\partial M}}{{\partial {x_i}}}}} = 1. $

根据导重法构造迭代式的准则，令

(10) $ {H_i} = {{\partial M}}/{{\partial {x_i}}}, $

(11) $ {W_i} = {x_i}{H_i}, $

(12) $ G_i^j = - {x_i}{{\partial {g_j}}}/{{\partial {x_i}}}, $

(13) $ {G^j} = \mathop \sum \nolimits_{i = 1}^4 G_i^j. $

将H_i看作变量x_i的容重，W_i看作变量x_i的等效质量，G^j_i看作变量在约束条件j下的导重，G^j看作约束条件j下的总导重.

将式（10）、（12）代入式（9），可以得到导重法的迭代准则：

(14) $ {x_i} = \frac{{{\lambda _1}G_i^1 + {\lambda _2}G_i^2}}{{{H_i}}};\;i = 1,2,3,4. $

在具体迭代求解时，将其写成如下形式：

(15) $ x_i^{\left( {k + 1} \right)} = {\left[ {{{{\big(\lambda _1}G_i^1 + {\lambda _2}G_i^2\big)}}/{{{H_i}}}} \right]^{\left( k \right)}};\;i = 1,2,3,4. $

引入步长因子α来保证迭代的收敛性，考虑优化参数的上、下限：

(16) $ x_i^{\left( {k + 1} \right)} = \left\{ {\begin{array}{*{20}{l}} {{x_{i{\rm{max}}}},{x_i} \geqslant {x_{i{\rm{max}}}}};\\ {\alpha {{\left( {\dfrac{{{\lambda _1}G_i^1 + {\lambda _2}G_i^2}}{{{H_i}}}} \right)}^{\left( k \right)}} + \left( {1 - \alpha } \right)x_i^{\left( k \right)},}\\ \qquad{{x_{i{\rm{min}}}} \leqslant {x_i} \leqslant {x_{i{\rm{max}}}}};\\ {{x_{i{\rm{min}}}},\;{x_i} \leqslant {x_{i{\rm{min}}}}}; \end{array}} \right. $

式（16）即为槽钢多约束质量优化问题的导重法迭代式.

根据约束条件，推出求解2个Kuhn-Tucker乘子λ₁、λ₂的线性不等式方程组为

(17) $ \left[ \!\!\!{\begin{array}{*{20}{c}} {\displaystyle \sum \limits_{i = 1}^4 (G_i^1G_i^1/{W_i})}&{\displaystyle \sum \limits_{i = 1}^4 (G_i^1G_i^2/{W_i})}\\ {\displaystyle \sum \limits_{i = 1}^4 (G_i^2G_i^1/{W_i})}&{\displaystyle \sum \limits_{i = 1}^4 (G_i^2G_i^2/{W_i})} \end{array}}\!\!\! \right] \times \left[ \!\!\!{\begin{array}{*{20}{c}} {{\lambda _1}}\\ {{\lambda _2}} \end{array}}\!\!\! \right] \geqslant \left[ \!\!\!{\begin{array}{*{20}{c}} {{G^1} + \dfrac{{{g_1}}}{\alpha }}\\ {{G^2} + \dfrac{{{g_2}}}{\alpha }} \!\!\! \end{array}} \right]. $

将式（17）化为线性互补问题，采用Lamker算法^[16]，可以求得λ₁、λ₂.

根据导重法的迭代式（16），迭代求解的关键是H_i与G^j_i的求解，即目标函数与性态约束对设计变量的偏导数.

灵敏度是设计变量对于目标函数的影响程度，数学意义是目标函数对设计变量的偏导数的求解. 结构性能参数的设计变量u_j对设计变量x_i的灵敏度可以定义为

(18) $ {\rm{Sen}}\left( {{{{u_j}}}/{{{x_i}}}} \right) = {{\partial {u_j}}}/{{\partial {x_i}}}. $

对结构进行灵敏度分析，可以求得结构静力学特性对所有优化参数的差分灵敏度，即为 ${{\partial M}}/{{\partial {x_i}}}$、 ${{\partial {g_j}}}/{{\partial {x_i}}}$的近似值.

在结构灵敏度求解方法上，目前有全局差分法、解析法与半解析法^[17]. Hansen等^[18-20]提出具有高收敛性的灵敏度计算的新解析法；为了提高灵敏度求解的精度，Keulen等^[21]提出改进的灵敏度计算的半解析法，该方法通过分解节点位移，创建虚载荷来获得高精度；Bletzinger等^[22]引入计算矫正因子来提高半解析法的求解精度；Voorhees等^[23]提出复变函数灵敏度求解法，使得计算步长独立，不会产生精度问题.

从结构单元入手，采用半解析法得到的静力学问题位移灵敏度为

(19) $ \frac{{\partial {{U}}}}{{\partial {{X}}}} = {{{K}}^{ - 1}}\left( {\frac{{\partial {{F}}}}{{\partial {{X}}}} - \frac{{\partial {{K}}}}{{\partial {{X}}}}{{U}}} \right). $

应力灵敏度为

(20) $ \frac{{\partial {\rm{\sigma }}}}{{\partial {{X}}}} = \frac{{\partial {{B}}}}{{\partial {{X}}}}{{U}} + {{B}}\frac{{\partial {{U}}}}{{\partial {{X}}}}. $

式中：K为结构总体刚度矩阵，B为单元应变矩阵，U为节点位移矩阵，σ为结构应力，F为载荷向量.

利用差分法，将位移矩阵、刚度矩阵及应变矩阵对设计变量的偏导数转化为差分形式：

(21) $ \frac{{\partial {{U}}}}{{\partial {{X}}}} = \frac{{ {{{U}}\left( {{{X}} + \Delta {{X}}} \right) - {{U}}\left( {{X}} \right)} }}{{\Delta {{X}}}}, $

(22) $ \frac{{\partial {{K}}}}{{\partial {{X}}}} = \frac{{ {{{K}}\left( {{{X}} + \Delta {{X}}} \right) - {{K}}\left( {{X}} \right)} }}{{\Delta {{X}}}}, $

(23) $ \frac{{\partial {{B}}}}{{\partial {{X}}}} = \frac{{ {{{B}}\left( {{{X}} + \Delta {{X}}} \right) - {{B}}\left( {{X}} \right)}}}{{\Delta {{X}}}}. $

因为F与设计变量无关，结合式（19）~（23）可以求出位移灵敏度、应力灵敏度.

为了消除刚体转动位移带来的影响，以位移灵敏度为例，引入单元误差修正项E_e. 带有修正项的位移灵敏度方程为

(24) $ \frac{{\partial {{U}}}}{{\partial {{X}}}} = {{{K}}^{ - 1}}\mathop \sum \nolimits_e \left[ {\frac{{\partial {{F}}}}{{\partial {{X}}}} - \left( {\frac{{\partial {{{K}}_e}}}{{\partial {{X}}}} + {{{E}}_e}} \right){{{U}}_e}} \right]. $

式中：E_e为单元误差修正项，K_e为单元刚度矩阵，U_e为单元位移矩阵. 将单元位移累加于总体，考虑摄动前的刚度阵. 对于静力问题，误差修正表达式^[24]为

(25) $ \left. {\begin{array}{*{20}{c}}\displaystyle {e = \mathop \sum \nolimits_e {\mu _e} = \mathop \sum \nolimits_e {\beta _i}{{\varPhi }}_{\rm{r}}^j = \mathop \sum \nolimits_e {a_{ij}}{{\varPhi }}_{\rm{r}}^{j{\rm{T}}}{{{U}}_e}{{\varPhi }}_{\rm{r}}^j},\\ {{\beta _i} = {a_{ij}}{c_j} = {a_{ij}}{{\varPhi }}_{\rm{r}}^{j{\rm{T}}}{{{U}}_e}},\\ {{c_j} = {{\varPhi }}_{\rm{r}}^{j{\rm{T}}}{{{U}}_e};\; {i,j = 1,2, \cdots, n}}. \end{array}} \right\} $

式中：a_ij、β_i、c_j均为计算系数，Φ^j_r为单元刚体转动位移矩阵.

刚体平动位移不影响转动位移，当节点产生平动位移时，可以表示成

(26) $ \left. {\begin{array}{*{20}{c}} {{{\varPhi }}_{\rm{t}}^1 = {{\left[ {1\;0\;0\;0\;0\;0\;1\;0\;0\;0\;0\;0} \right]}^{\rm{T}}}},\\ {{{\varPhi }}_{\rm{t}}^2 = {{\left[ {0\;1\;0\;0\;0\;0\;0\;1\;0\;0\;0\;0} \right]}^{\rm{T}}}},\\ {{{\varPhi }}_{\rm{t}}^3 = {{\left[ {0\;0\;1\;0\;0\;0\;0\;0\;1\;0\;0\;0} \right]}^{\rm{T}}}}. \end{array}} \right\} $

式中：Φ¹_t、Φ²_t、Φ³_t分别为刚体沿x、y、z方向的平动位移矩阵.

根据K_eΦ_r=0可以看出，刚体转动位移是刚度矩阵特征值为0时所对应的特征向量，可以通过求解刚度矩阵特征向量得到刚体的转动位移：

(27) $ \left. {\begin{array}{*{20}{c}} {{{\varPhi }}_{\rm{r}}^1 = {{\left[ {0\;0\;0\;1\;0\;0\;0\;0\;0\;1\;0\;0} \right]}^{\rm{T}}}},\\ {{{\varPhi }}_{\rm{r}}^2 = {{\left[ {0\;0\;l\;0\;1\;0\;0\;0\;0\;0\;1\;0} \right]}^{\rm{T}}}},\\ {{{\varPhi }}_{\rm{r}}^3 = {{\left[ {0\; - l\;0\;0\;0\;1\;0\;0\;0\;0\;0\;1} \right]}^{\rm{T}}}}. \end{array}} \right\} $

将刚体平动位移和转动位移向量正交单位化，代入误差修正式（25），可得误差修正项.

将求得的λ₁、λ₂、 ${{\partial M}}/{{\partial {x_i}}}$、 ${{\partial {g_j}}}/{{\partial {x_i}}}$代入式（16），即可完成1次迭代计算. 经过多次迭代求解，最终得到目标函数以及约束函数的收敛结果. 优化后设计变量的值见表2，将此时设计变量X的值作为外门架轻量化设计的最优解X^*.

在外门架立柱优化结束后，利用同样的方法对内门架进行优化，再对货叉、货叉架、内外门架横梁等一些应力较小的构件进行轻量化处理，优化前、后各部分质量对比见表3. 表中，m₀、m₁分别为优化前、后的质量，p为优化比例.

以优化后尺寸为设计参数，建立门架的三维模型，利用有限元法进行强度和刚度校核. 分析工况为门架直立起升至3 m高的位置；采用solid186实体单元划分网格，网格大小设置为10 mm，货叉采用六面体网格，内外门架以及货叉架采用四面体网格，其余部分采用自动网格划分；滚轮与槽钢接触设置为摩擦接触，摩擦系数取0.2；载荷为叉车额定载荷3 t，换算成均布载荷为129 360 Pa，施加在货叉上表面. 求解得到门架静力学分析结果，如图6所示.

如图6（a）所示为门架在满载工况下上升至3 m高位置的应力云图. 可以看出，门架在最危险工况下的最大应力发生在货叉架下横梁与背板焊接处，最大应力为364.94 MPa. 焊接处的失效形式以焊缝断裂为主，主要考虑抗拉强度. Q345抗拉强度为490~620 MPa，最大应力在安全范围之内.

内外门架经过轻量化设计，其安全性能需要校核，根据有限元分析结果可以看出，门架在优化后，内门架最大应力发生在滚轮与翼缘接触处，最大应力为294.59 MPa，小于材料的屈服极限345 MPa. 滚轮与翼缘的接触满足赫兹接触的条件，可以通过对滚轮母线的修形来降低接触应力^[25].

外门架最大应力发生在上横梁与立柱的焊接处，最大应力达到252.04 MPa，小于材料的屈服极限345 MPa. 有限元分析结果表明，门架在经过轻量化设计后，其强度条件满足设计要求.

为了校核门架优化后刚度是否满足设计要求，通过有限元计算门架顶端在水平方向的位移云图，如图6（b）所示. 在内门架顶端产生最大水平位移，最大位移为22.741 mm，即f=22.741 mm，门架许用位移 $\left[ f \right] ={H}/{{100}} = {{3\;000}}/{{100}} = 30\;{\rm{mm}}$. f<[f]，所以门架的刚度满足设计要求.

（1）将灵敏度求解的带有经典误差修正项的半解析法（exact semi-analytical technique，ESA）与导重法结合，运用到叉车门架的轻量化设计上效果显著，解决了导重法容重与导重难以求解的问题，提高了优化求解的效率.

（2）以叉车门架质量为优化目标，以门架强度和刚度为性态约束的轻量化设计切实可行，经过优化，门架质量降低了116.04 kg，比优化前质量降低了18.21%. 经过有限元分析校核可知，优化后的门架强度、刚度均满足要求.