基于区间数的高速磨床主轴系统的不确定性优化研究＊

崔 中，文桂林，赵子衡，姜 潮

（湖南大学 汽车车身先进设计制造国家重点实验室，湖南 长沙 410082）

机械加工正在向着高精度、高效率的方向发展，高速磨削已成为现阶段精加工的重要手段.在精密磨削中，主轴系统的振动对工件的表面质量影响很大，而且对主轴轴承和砂轮的微刃都有损害.不平衡的砂轮在高速运转中会导致主轴系统的强迫振动.因此减小主轴系统的振动和优化主轴系统的性能成为国内外提高磨削质量的重要研究方向.在传统工程优化中，结构参数和作用载荷都被视为确定参数，而对于高速磨床这种工程实际问题，其材料特性、几何参数和作用载荷等并不是确定的，特别是实际工况下的载荷.为了得到有用且可靠的结构分析和设计，在结构设计和优化过程中考虑不确定性因素是很必要的.不确定性优化是传统优化方法的延伸，使得优化的结果更具有实际意义.

对于工程实际中的不确定性问题，近年来，大多采用概率方法处理不确定问题，把不确定变量作为随机变量.然而由于概率方法对于实际工程问题有其本身的局限性，因为有时很难提供充分的概率分布信息来描述这些随机变量，从而限制了基于概率的优化方法在工程中的应用[1－2].而基于区间方法的不确定性优化把不确定参数看作区间数，只需要不确定参数的上下界信息即可，这些信息在工程实际中是比较容易获得的，如此以来基于区间数不确定性优化在工程应用中的优势就显而易见了[3].对于主轴系统的优化问题，目标函数是通过建立柔性体动力学模型进行数值计算得到的，由于结构较大，计算费时，而且基于区间数的不确定优化是一个嵌套优化，如果调用动力学模型进行优化，就会造成优化的效率非常低，所以本文引入近似模型方法来代替主轴系统的动力学模型进行不确定优化，从而很大程度上提高了求解和寻优的效率.

1 基于区间数方法的不确定优化理论

式中：f为包含设计变量和不确定参数的目标函数；gi为第i个不确定约束函数；m为不确定约束的个数；X为n维设计向量；Ω为n维设计空间；p为q维不确定向量，上标L和R表示区间的下界和上界.

利用区间数理论中的区间数序关系对区间进行比较，则式（1）的不确定优化问题可以转化为式（2）的多目标优化问题：

一般的不确定性优化问题可以表述为式（1）：

式中：m（f（X，p））为目标函数区间的中点；w（f（X，p））为目标函数区间的半径.fL（X）和fR（X）分别由式（3）表达式得出：

对于每个确定的X，f（X，p）的取值为一个由式（3）确定的区间，如此以来，式（2）所表述的不确定两目标优化问题就可以转化为一个确定的两目标优化问题.

由于线性加权组合法相对简单并广泛地被应用于工程问题中，所以本文采用线性加权组合法来处理两目标优化问题，式（2）的两目标不确定优化问题可转化为如式（4）的单目标优化问题：

式中：α为权重系数，其取值范围为α∈[0，1]，α的选取依靠对实际问题充分的了解.β和γ为使m（f（X，p））＋β和w（f（X，p））＋γ非负的常数.φ和ψ为标准化因子，其表达式如下：

采用罚函数方法来处理不确定优化问题的约束条件，使其转化为单目标无约束优化问题.式（4）可写为：

式中：σ为惩罚因子，通常取一个较大的值，当不满足约束条件时将受到惩罚，使得搜索保持在可行域内.

2 高速磨床主轴系统柔性体模型

2.1 ADAMS柔性体理论

在多刚体动力学软件ADAMS中建立柔性体模型，可以在考虑到物体的弹性的基础上，对复杂系统进行有效的动力学仿真，提高了系统仿真的精度，同时可以显示出物体在运动受力时的弹性变形和应力情况.

在ADAMS中是采用模态柔性体来表示物体弹性的，将弹性体离散为有限元模型，从而可以用有限的自由度来表示无限的自由度.基本方法是计算物体的模态，不同的模态是互相垂直的，从而构成一个线性的模态空间，在模态空间中用模态特征向量的线性组合来表示物体有限元节点的弹性位移，通过计算每一时刻物体的弹性位移来描述其变形运动.物体的变形可以在模态空间中通过模态的线性叠加关系得到，用如下表达式描述：

式中：d为各个节点的位移矢量；αi是模态参与因子；φ为物体的特征位移矢量.

本文采用的是利用ADAMS／FLEX模块在模型中引用模态中性文件的方法建立柔性体系统模型，即通过有限元分析软件对物体进行模态计算，生成模态中性文件——MNF文件，导入ADAMS中建立柔性体.

2.2 主轴系统的柔性体建模

某型高速磨床的主轴系统包括：砂轮架箱体，主轴，主轴静压轴承，砂轮等.

采用通用有限元软件MSC.Patran建立主轴和砂轮架箱体的有限元模型.主轴是一个多阶梯带有锥度的圆柱体，砂轮架的内部结构比较复杂，为了方便进行有限元建模和计算，两者的结构必须经过一定的简化，简化的原则为：不考虑各处的倒角，忽略空刀槽和储油槽；忽略很小的台阶和螺纹孔.

主轴是主轴系统的重要零件，使用六面体单元对其进行网格划分，因为六面体单元形状规则、单元质量好，计算精度高.而砂轮架箱体结构复杂，四面体单元的适应度较好，所以采用四面体单元对其进行网格划分.其有限元模型如图1和图2所示.

图1 主轴有限元模型Fig.1 Finite element model of spindle

图2 砂轮架有限元模型Fig.2 Finite element model of grinding carriage

由于主轴和砂轮架箱体的有限元模型要导入ADAMS中进行加载仿真，所以需要在其有限元模型中设置连接点.在主轴固定砂轮的轴端中心和前、后轴承中心以及安装电机的轴段中心位置分别设置一个连接点，用于在导入ADAMS后与刚性体的连接、创建运动副和成为载荷的受力点.同样在砂轮架箱体有限元模型的两个轴承座的中心处建立连接点.连接点采用Patran中的多点约束（multi－point constraint，MPC），具体如图2所示.

主轴的材料为65Mn，弹性模量E＝0.21TPa，泊松比μ＝0.3，密度为ρ＝7.81×103kg／m3.砂轮架的材料为 HT200，E＝0.12TPa，μ＝0.25，ρ＝7.4×103kg／m3.

采用MSC.Nastran对主轴和砂轮架箱体的有限元模型进行模态分析，并生成MNF文件.

2.3 ADAMS中主轴系统动力学建模

通过ADAMS／FLEX模块将生成的主轴和砂轮架箱体的模态中性文件导入ADAMS中，进行约束和加载.

由于在ADAMS中并不是所有的约束副和运动副都可以直接加在柔性体上的，需要通过哑物体过渡，在主轴和砂轮架箱体上的连接点处建立哑物体.采用轴套约束来模拟主轴静压轴承对主轴的支撑和受力关系.在主轴末端的连接点上建立圆柱副并施加旋转驱动，来模拟电机.

在磨削加工中，磨削力有3个分量：Ft为切向磨削力，Fn为法向磨削力，Fs为纵向进给方向的分力，如图3所示.对于纵向进给方向的力Fs来说，由于各个磨粒具有随机分布的正负倾角，使其各分力互相抵消，从而相对于其他两个磨削力很小，所以只考虑切向磨削力和法向磨削力[4].

图3 外圆磨削过程中的磨削力Fig.3 Grinding force in the cylindrical grinding process

切向磨削力Ft和法向磨削力Fn可以用式（7）和式（8）来计算：

的生长基元为[MgO6]正八面体，其通过Mg-O键以共顶点的方式紧密相连，沿化学键作用力较强的[010]方向无限连接形成长链。采用线性拟合计算得lgtind与f(lgs)的关系为y=3.828 25-5.199 65x，决定系数R2=0.983 49，拟合结果与试验结果良好吻合，重镁水溶液浓度增大，诱导期时间缩短。

式中：k为与工件材料相关的系数；b为磨削宽度；λ为有效磨粒间隔；γ为圆锥半顶角；dw为工件直径；ds为砂轮直径；ω1为工件转速；ω2为砂轮转速；Δ为工件每转磨削深度；ε为指数[5]，一般取0.2～0.5.在该磨床中dw＝100mm，ds＝500mm，ω1＝150 r／min，ω2＝360rad／s，Δ＝0.1mm，采用 CBN 砂轮，磨削宽度近似为砂轮的宽度，b＝25mm，ε＝0.5，γ＝60°，λ＝0.5mm.工件材料选用淬火45钢，k＝165kg／mm.将这些参数代入式（7）和式（8）中得，Fn＝910N，Ft＝412N.

考虑到磨削的实际工况，将法向磨削力Fn视为周期载荷，加载在砂轮的质心上，方向为Y轴的负方向.由于切向磨削力距离主轴250mm，因此在主轴上加一个扭矩T＝412×250＝103 000N·mm.

在ADAMS中用STEP函数来对力和力矩进行加载，使其在1～1.2s之间平滑过渡，如图4所示.

图4 磨削力和力矩Fig.4 Grinding force and moment

对该模型进行动力学仿真计算，主轴在运动过程中应力变化情况都可以在后处理动画中看到，如图5所示.从结果中可以看出，在运动过程中，主轴在靠近电机连接处的应力较大，最大应力为10.4 MPa.由于将砂轮视为刚体，所以砂轮质心的位移与边缘接触点的位移是一致的，从而砂轮质心的位移情况可以在一定程度上表征磨削时主轴的受力振动对工件表面质量的影响.观察砂轮架质心在X，Y和Z3个方向上的位移曲线，如图6所示，可以看到在Y方向上（主轴的径向）砂轮质心在一开始运转时就有17μm的位移，之后在周期力的作用下又有最大幅值为3.3μm的跳动.本文以表征主轴系统刚度的砂轮质心最大的跳动幅值为响应，建立目标函数.

图5 砂轮架主轴系统运动过程应力图Fig.5 Grinding stress figure of spindle－carriage system

图6 砂轮质心在X，Y和Z3个方向的位移变化图Fig.6 Displacement variation figure of grinding wheel’s barycenter

3 主轴系统近似模型的建立

由于对柔性体主轴系统模型进行数值计算的时间较长，如此以来基于该模型进行优化计算，特别是进行嵌套优化计算，每次计算响应值都要调用费时的柔性体模型，使计算效率低下，从而阻碍了优化算法在主轴系统优化设计中的应用.近似模型方法是目前较为流行的解决大规模模型优化问题的有效途径，用其代替真实模型进行优化，可以很大程度地提高求解工程优化问题的效率.

3.1 试验设计

为了建立有效的近似模型，需要借助实验设计方法，以获取足够的、合适的样本点.采样点的选择不当会导致近似模型的精度低、甚至错误，并且会提高构建近似模型的成本.

本文采用优化拉丁超立方试验设计方法，该方法的采样点在采样空间中分布均匀，可以在相对较少的采样个数的情况下，获得充分的模型信息.该方法的基本思想是将每个设计参数的设计空间均匀地划分为N×N的方阵，然后在其中随机生成不同行不同列的N个采样点.

由于基于区间的不确定性优化的目标函数是优化设计变量和不确定参数的函数，所以在试验设计中，采样空间应该为由设计空间和不确定参数空间组成的混合空间.

在高速磨削的实际工况中，磨削力不可能是一个确定量，它受到诸如磨削表面粗糙度、砂轮接触面上动态磨刃数、砂轮耐用度、磨削比能等因素的直接影响.本文以轴承的径向刚度k1和轴向刚度k2为设计变量，法向磨削力Fn和磨削力矩T为不确定参数，在其组成的采样空间内选取样本点.

3.2 径向基近似模型

近似模型的基本原理是，通过数理统计和试验设计的方法，建立设计变量和响应值之间的函数关系，用以代替复杂的真实模型.本文选用径向基函数近似模型，因其在考虑模型精度和鲁棒性的同时，相比其他近似模型要可靠[6－7].

径向基函数近似模型是以径向函数为基函数，通过线性加权插值构造出来的.径向基函数是以待测点与样本点之间的欧式距离为变量的函数.本文选用Gaussian函数作为径向基函数，其表达式如式（9）所示，则径向基函数近似模型的解析表达式如式（10）所示.

式中：δi为权系数；Li为待测点到样本点的欧氏距离，Li＝‖x－xi‖；α是给定的大于零的常数；n为样本点个数.当已知样本点 Xj（j＝1，...，n）和其响应值yj后，利用插值条件p（Xj）＝yj，计算出径向基函数近似模型的权系数δi，从而构建出显式表达的径向基近似模型[8].

将试验设计中得到的由径向刚度k1、轴向刚度k2、法向磨削力Fn和磨削力矩T组成的样本点，代入到柔性体主轴系统模型中进行计算，得出其相应的砂轮质心最大跳动幅值，根据上述方法，建立以主轴径向刚度k1、轴向刚度k2、法向磨削力Fn和磨削力矩T为自变量，砂轮质心最大跳动幅值为目标值的主轴系统近似模型F（k1，k2，Fn，T）.

4 主轴系统的不确定性优化研究

本文主轴系统不确定优化研究的目的是使砂轮质心在工作状态下最大跳动值最小，从而一定程度上提高主轴的回转精度和工件表面的磨削精度.以上述径向基函数近似模型F（k1，k2，Fn，T）为目标函数，约束条件为轴承的径向刚度和轴向刚度的取值范围，k1∈ [8×104N／mm，1×106N／mm]，k2∈[1×104N／mm，2×105N／mm].

结合本问题实际情况，在上述区间不确定性优化方法中，将区间中点和区间半径视为同等重要，则式（6）中权重系数α＝0.5.常系数β＝γ＝10.主轴系统的不确定优化问题的目标函数可用表达式（11）描述：

其中：

采用隔代遗传算法（intergeneration projection genetic algorithm，IPGA）对上述目标函数进行优化.隔代遗传算法是在微型遗传算法的基础上加入隔代映射算子，隔代映射算子是用来沿着连续两代中最好个体的方向上，寻找更好的个体，以代替当前代中最差的个体[9].隔代遗传算法不仅继承了微型遗传算法种群规模小、基因多样性和全局优化的特点，而且提高了收敛速率.

首先，外层遗传算法在主轴径向刚度k1和轴向刚度k2组成的设计空间内寻优，对于每个所取的设计向量进入内层遗传算法，在不确定磨削力Fn和T组成的不确定参数空间内搜索，通过计算径向基近似模型F（k1，k2，Fn，T）确定目标函数响应的上下界，进而得到目标函数响应区间的中点和半径.把内层优化结果反馈给外层优化算法，以帮助外层算法继续寻优，直到满足停止准则输出最后的设计变量作为优化结果.优化流程如图7所示.

最后的优化结果为，在不确定参数作用下，当k1＝771 573.8N／mm，k2＝43 802.5N／mm 时，砂轮质心最大跳动值区间为[1.61μm，2.49μm]，不确定区间中心为2.05μm，区间半径为0.44μm.以区间中点来说，相比优化前的砂轮质心最大跳动值降低了37.9%，即使最坏情况为2.49μm也比优化前降低了24.5%，结果令人满意.

我们可以通过调整液体静压轴承一些参数，如供油压力、轴承长度、轴承间隙等[10]，实现优化出来的轴承刚度参数.

图7 主轴系统优化流程图Fig.7 Optimizaiton flowchart of spindle system

5 结 论

1）将基于区间数的不确定性优化方法引入到高速磨床主轴系统的优化中，在考虑实际工况中载荷不确定性的基础上，对主轴系统的结构参数进行优化.优化结果显示，当主轴系统刚度k1＝771 573.8 N／mm，k2＝43 802.5N／mm时，表征磨削质量的砂轮质心径向最大跳动值区间为[1.61μm，2.49μm].主轴系统的振动性能有了较大的提高，而且使优化结果更有实际意义.

2）建立高速磨床主轴系统的柔性体动力学模型，提高了磨削仿真的精度，可以更实际地反映出磨削时主轴系统的受力和振动情况.

3）根据高速磨床自身模型大、计算费时等特点，建立结构参数与目标值之间的径向基函数近似模型，用于代替动力学模型进行优化，从而大大提高了求解高速磨床这类复杂有限元模型优化问题的效率.

[1] CHO G M.Log－barrier method for two－stage quadratic stochastic programming[J].Applied Mathematics and Computation，2005，164（1）：45－69.

[2] SENGUPTA A，PAL T K，CHAKRABORTY D.Interpretation of inequality constraints involving interval coefficients and a solution to interval linear programming[J].Fuzzy Sets and Systems，2001，119：129－138.

[3] JIANG C，HAN X，GUAN F J，etal.An uncertain structural optimization method based on nonlinear interval number programming and interval analysis method[J].Engineering Structures，2007，Doi：10.1016／j.engstruct，2007.01.020.

[4] 李伯民，赵波.现代磨削技术[M].北京：机械工业出版社，2004.LI Bo－min，ZHAO Bo.Modern grinding technology[M].Beijing：China Machine Press，2004.（In Chinese）

[5] 王龙山，李国发.磨削过程模型的建立及其计算机仿真[J].中国机械工程，2002，13（1）：1－5.WANG Long－shan，LI Guo－fa.Modelling and computer simulation for grinding process[J].China Mechanical Engineering，2002，13（1）：1－5.（In Chinese）

[6] JIN R，CHEN W，SIMPSON T W.Comparative studies of metamodeling techniques under multiple modeling criteria[J].Journal of Structural and Multidisciplinary Optimization，2001，23（1）：1－13.

[7] 穆雪峰，姚卫星，余雄庆，等.多学科设计优化中常用代理模型的研究[J].计算力学学报，2005，22（5）：608－612.MU Xue－feng，YAO Wei－xing，YU Xiong－qing，etal.A survey of surrogate models used in MDO[J].Chinese Journal of Computational Mechanics，2005，22（5）：608－612.（In Chinese）

[8] 文桂林，崔中，彭克立.基于近似模型的高速磨床零部件结构优化设计研究[J].中国机械工程，2009，20（8）：906－910.WEN Gui－lin，CUI Zhong，PENG Ke－li.Structural optimization for high speed grinder’s components based on the approximate model[J].Chinese Journal of Computational Mechanics，2009，20（8）：906－910.（In Chinese）

[9] 刘桂萍，韩旭，姜潮.并行隔代映射遗传算法及其在材料参数反演中的应用[J].中国机械工程，2007，18（4）：387－390.LIU Gui－ping，HAN Xu，JIANG Chao.Parallel intergeneration projection genetic algorithm and its application to determination of material parameters[J].China Mechanical Engineering，2007，18（4）：387－390.（In Chinese）

[10] 丁叙生.以主轴系统刚度与总功耗之比值为目标函数的液体静压轴承优化设计方法[J].南昌航空工业学院学报，2000，14（3）：26－29.DING Xu－sheng.The optimum design of hydrostatic bearing[J].Journal of Nanchang Ininstitude of Aeronautical Technology，2000，14（3）：26－29.（In Chinese）