风洞模型支撑机构有限元模型修正与确认
2021-10-19麻越垠陈万华
麻越垠，郝 鹏，马 斌，陈万华
(1.中国空气动力研究与发展中心，绵阳 621000；2.大连理工大学工程力学系工业装备结构分析国家重点实验室，大连 116024)
0 引言
在工程设计领域，随着数值模拟技术的快速发展和计算机运算能力的不断提升，基于高可信度数值模拟的设计方法已经开始发挥越来越重要的作用。有限单元法是工程分析应用最为广泛的数值方法，有限元模拟分析目前已经成为与理论分析、试验研究鼎足而立的重要科学研究手段之一。为减小有限元模型参数误差，使有限元模型更可靠的反应实际结构特性，需应用理论数据或试验数据对有限元模型实施修正，减小数值计算结果和试验测试结果之间的偏差。结构有限元修正问题可转化为参数优化问题,即对目标函数求极值。使用优化方法最主要的一个障碍就是所采用的数值模拟耗时太长,实际工程问题往往需要多轮迭代优化，若每一次迭代都调用一次数值模拟，将非常耗时甚至不可行。近似建模技术是通过数学模型的方法逼近一组输入变量与输出变量的方法，避免优化设计中频繁调用仿真程序。20世纪80年代，波音777的GE90涡扇发动机设计项目中，不同专业的工程师通过手工“试算-评估-校正”方法进行发动机优化设计，耗时巨大且优化效果不明显。为解决此类多学科优化问题，计算机辅助优化平台iSIGHT由此诞生，具有自编程集成接口，融合多种优化探索策略，设计仿真流程集成自动化，有效提升优化效率。吴思远等提出基于iSIGHT-MSC Patran/Nastran联合仿真与特征灵敏度优化相结合的快速建模方法，完成弧形闸门动力模型修正。吴建云等利用ADAMS和iSIGHT软件联合仿真，采用优化拉丁方试验设计缩小修正参数范围，采用自适应模拟退火算法进行模型参数修正，修正后的仿真模型更符合实际试验样机。李杰等基于iSIGHT联合仿真的方案，采用最优拉丁超立方试验设计对重型汽车平顺性进行了灵敏度分析。李永洲等利用iSIGHT软件对反正切马赫数分布可控的轴对称基准流场设计参数进行灵敏度分析。
本文根据参数型有限元模型修正一般流程，联合ABAQUS、MATLAB搭建iSIGHT结构有限元模型修正平台，开展某模型支撑机构动力学模型修正，以试验设计变量对响应的贡献量为依据，筛选修正变量，构建前4阶频率差和振型相关性为目标函数，使用近似建模得到修正变量和目标函数的响应面模型，采用多目标优化方法非劣排序遗传算法NSGA-II完成风洞模型支撑机构有限元模型多目标修正；并对修正后的有限元模型开展动力学响应确认，计算锤击激励下的结构模态动力学响应，与锤击试验响应作对比，结果显示，修正后的模型满足工程需求。
1 多目标有限元模型修正
模型修正方法有矩阵型法和参数型法两大类。矩阵型法是最早发展一类模型修正方法，起初是根据飞机地面共振试验数据修改结构柔度矩阵，用直接求解的方法，计算效率很高；缺点是修正目标的物理意义不明确，操作困难。参数型法是以工程师建模时所使用的设计参数为修正对象，包括物理参数、几何参数、边界条件和连接刚度等，结果易于解释，便于工程判断，与CAE软件配合，可以实现参数优化设计，是目前模型修正主要发展方向。图1为参数型有限元模型修正一般流程，由有限元初始模型得到结构动力学计算结果，与动力学试验结果开展相关性分析；对潜在的待修正参数进行灵敏度分析，选择修正参数；以结构动力学特性构建目标函数，开展优化设计，也可以通过试验设计等手段构建修正变量的样本空间，建立目标函数的近似模型，根据近似模型直接开展数值优化，最终输出修正后的参数。采用修正后的有限元模型，来研究其在其他设计工况下的动力学行为。
图1 参数型有限元模型修正流程Fig.1 Procedure of parametric FEA model updating
工程问题中的优化设计存在数值模拟耗费时间过长的问题，近似模型方法是通过数学模型的方法逼近一组输入变量与输出变量的方法，通过显式函数来表达结构参数和响应值之间的隐式关系，从而获得结构的近似模型(也称之为代理模型)；在近似模型基础上进行的迭代计算，完全属于数学的运算过程，可以明显提高计算效率，减少数值迭代计算的次数。响应面模型是一种常用的近似建模方法，假设待修改的结构参数x
(i
=1,2,…,n
)与有限元模型计算所得的结构性能指标y
(位移、应力、频率、振型等)的关系为y
=f
(x
)+ε
(1)
式中，f
(x
)是x
区域上的一个曲面的响应函数；ε
为随机误差;实际问题中，f
是隐式的，通常采用多项式或其他响应面模型近似描述f
，2阶多项式是最常见的响应面模型，其形式为
(2)
式中，待定系数a
表示变量x
的线性效应；当i
≠j
时，a
表示x
和x
的交互作用效应；当i
=j
时，a
表示变量x
的二次效应。有限元动力学模型修正的优化目标函数通常同时包括频率信息和振型相关性信息，使用固有频率和振型相关性修正目标函数为
(3)
(4)
式中，f
,为第i
阶试验频率和计算频率的误差百分比，MAC
为第i
阶试验振型和计算振型的相关值，f
,和{},为第i
阶试验模态频率和振型，f
,和{},为第i
阶计算态频率和振型。2 模态试验
模型支撑机构是风洞试验中固定和调节模型姿态的机构，其动态特性决定风洞试验数据精度。图2为模型支撑机构三维模型，机构主要由框架、弯刀支板、支杆、直线导轨系统、圆弧导轨系统、丝杠系统和侧向压轮系统组成。弯刀支板和支杆是整个机构实现所有功能的结构载体，丝杠系统在模型支撑机构运动中起传动作用，圆弧导轨系统起导向作用，侧向压轮系统包括固定压轮组件和弹性压轮组件，固定压轮组件通过螺钉调节弯刀支板的位置，另一侧的弹性压轮组件通过碟簧施加作用在弯刀支板上的顶紧力。结构整体安装在框架上，框架与地面采用地脚螺栓固定。
图2 风洞模型支撑机构三维模型Fig.2 The 3D wind tunnel model of support mechanism
通过试验模态确定结构的动力学特性，包括模态频率、模态振型、模态阻尼比等参数，为结构有限元模型修正提供依据。图3为模态试验对象和试验测点布置结果，共有55个三向测点，规定竖直方向为Z
向，弯刀支板法向为Y
向，X
向根据右手坐标系确定。
图3 模态试验测点布置Fig.3 Layout of test points for modal test
采用移动力锤法、多个参考点进行测试，测得前4阶(40 Hz以内)的试验模态频率分别为8.6，21.2，23.5，39.2 Hz。对应的模态振型如图4所示，1阶振型为Y
向摆动，2阶振型为绕Z
轴扭转，3阶振型为X
向摆动，4阶振型为支杆Y
向摆动。
图4 模态试验结果Fig.4 Results of modal test
使用模态置信准则MAC评估前4阶模态振型的相关性，如图5所示。各试验模态阵型互相关值均较小，最大0.025，远小于工程允许值0.2，说明试验振型正交性较好，满足修正需要。
图5 试验模态振型相关性Fig.5 Correlation of test mode shapes
3 模态分析与匹配
采用ABAQUS开展模型支撑机构动力学仿真分析，对结构设计三维模型做出如下简化：
1)底板上的螺栓孔去除，设为固定连接；
2)轴承连接简化为铰链连接；
3)丝杠系统、直线导轨系统和圆弧导轨系统均只考虑其功能自由度，丝杠系统、直线导轨系统简化为直线连接，圆弧导轨系统简化为转动连接；
4)侧向压轮系统的固定压轮组件与弯刀支板简化为平面滑动，弹性压轮组件与弯刀支板简化为考虑弹性的轴向连接；
5)所有部件均为钢制结构，采用统一材料属性，不考虑材料差异性。
求解40 Hz以下的模态，模态分析结果如图6所示，40 Hz以内共有4阶模态信息，计算模态频率分别为11.5，24.7，31.0和36.4 Hz。
图6 初始模态分析结果Fig.6 Result of modal analysis
将初始模型计算模态与试验模态进行匹配，表1为试验模态和计算模态匹结果，前4阶振型相关性均在0.8以上，达到良好水平；频率差分别为34.12%、16.36%、31.86%和7.08%，本项目要求的误差率不大于10%，且振型相关性(MAC)不小于0.8，因此，有限元动力学模型不满足项目要求，需要开展模型修正。
表1 试验模态和计算模态匹配
4 有限元模型修正与确认
由表1可知，初始计算模态与试验模态振型相关性较好，说明模型整体边界设置相对合理，修正变量应为内部部件的材料属性和连接刚度。初步选择前、后框架和弯刀支板的材料属性以及弹性压轮组件中的碟簧刚度，材料属性包括密度和弹性模量。
在分析时忽略地板、轴承等小零件和动力部件，如电机、丝杠螺母机构、导轨滑块机构等，导致结构整体质量减小，表现为按照标准钢材选择的材料密度可能偏小；同时，对螺栓、弹簧、轴承的简化，也导致部件刚度发生改变，表现为材料弹性模量的变化。根据工程经验，设置修正变量水平变化范围为基准值±20%。
表2为修正变量的代表符号、初始值和水平范围，共7个变量，根据计算简化的不同，各个变量对结果的影响度也不同，因此对变量进行灵敏度分析，筛选对修正目标影响较大的变量，减少修正变量个数，提升优化效率。
表2 初选修正变量 mm
采用表1中试验模态和计算模态前4阶频率差f
,(i
=1,2,3,4)和MAC(i
=1,2,3,4)作为响应函数，联合ABAQUS、MATLAB搭建iSIGHT风洞结构有限元模型修正平台，如图7所示，主要包括试验模态导入模块、有限元参数化模型模块、相关性分析模块、试验设计模块、近似模型模块和优化模块。
图7 基于近似建模的有限元模型修正平台Fig.7 FEA model updating platform
选择适合高维设计空间探索的优化的拉丁超立方法开展试验设计，样本点为40个，将修正变量归一化至[-1,1]区间，用最小二乘法拟合得到修正变量线性项，反映修正变量对响应的贡献量，如图8所示。变量S
对各响应贡献率均不超过5%，D
和E
的贡献率在25%左右，D
和E
的贡献率在15%左右，D
和E
的贡献率在30%左右。因此，选择D
,E
,D
和E
为修正变量，即修正目标为前框架、弯刀支板的材料密度和弹性模量。
图8 各修正变量对响应的贡献率Fig.8 Contributions of the variables to the responses
使用优化的拉丁超立方法设计D
,E
,D
和E
的样本空间，共有100个样本点。采用4阶响应面拟合频率差、振型相关性与修正变量的多项式。一般采用相关指数R
作为衡量近似模型精度的指标，R
计算公式为
(5)
n
个变量，m
个目标函数的多目标优化问题可表达为
(6)
式中，x
=(x
,x
,…,x
)∈⊂R
为变量，y
=(f
(x
),f
(x
),…f
(x
))∈⊂R
为目标函数，和分别为n
和m
维目标空间，g
(x
)为q
个不等式约束，h
(x
)为p
个等式约束。对于给定的两个点x
,x
∈(为可行解集合)，仅当：条件1：f
(x
)≤f
(x
)(i
=1,2,…,m
)条件2：至少一个目标函数f
(x
)<f
(x
)条件1和条件2同时满足，则认为x
支配x
，表达为x
≻x
，其中为x
非支配解，为x
支配解，非支配解集即为全局最优解集，记为Pareto解集，其在目标空间的映射，表示为Pareto前沿。求解Pareto解集有两种方法:1)选择适当的权重系数，将多目标优化转化为单目标优化。更改权重系数，再次优化。
2)基于种群的搜索，同时评估一系列设计方案，最常用的方法是非劣排序遗传算法NSGA-II。
NSGA-II使用遗传算法进行搜索，具有全局寻优的特点，目标函数的驱动来自可行解之间的支配关系和相互距离，NSGA-II能进化出高质量的Pareto解集。
使用多目标优化算法NSGA-II对拟合响应式进行多目标优化，各变量求解结果如表3所示。相对于初始值，D
和D
修正后增大，表示修正后的整体质量增大，因为初始模型中忽略一部分零部件，所以密度修正方向与预测方向一致；E
修正后减小，表示修正后结构整体刚度变弱，因为螺栓连接简化连接面绑定，会导致局部刚度增大，所以弹性模量修正方向与预测方向一致。
表3 修正变量优化值 mm
表4列出各个目标函数优化值和提升系数，f
,代表第i
阶计算和试验模态频率差，MAC表示第i
阶模态振型相关性。可以看出,NSGA-II修正结果满足要求，频率差均在10%以内，振型相关性比初始值有些许增大，均大于0.8。
表4 频率差和阵型相关性优化结果
由表4可知，修正后的前4阶模态与试验频率差分别减小了81.36%，61.31%，87.41%和75.71%。图9为前4阶频率差的Pareto前沿，理想的优化点靠近0，优化后的推荐点与预测一致。
图9 前4阶频率的Pareto前沿Fig.9 Pareto fronts of the first four mode frequencies
使用模态动力学计算修正后的模型在冲击载荷激励下的结构响应，各模态阻尼采用对应的试验模态阻尼，计算结果与锤击试验响应做比较。锤击试验如图10所示，锤击点位于前框架上部顶端，锤击方向为Y
向，监测点位于后框架中上部，采用加速度传感器监测为Y
向响应。
图10 锤击试验Fig.10 Impacting test
图11和图12分别是监测点Y
向的响应时域和频域曲线，修正模型在振动响应量级和响应频率幅值均与试验接近，模型可靠性高，可用作下一步复杂载荷作用下的结构动力学响应预测。
图11 时域响应曲线Fig.11 Time domain response
图12 频域响应曲线Fig.12 Frequency domain response
5 结论
1)根据试验设计变量对响应的贡献量，筛选修正变量；采用近似建模方法构建修正变量和响应的响应面模型，开展多个多目标优化方法的对比。结果显示：使用NSGA-II修正后的模型支撑机构有限元模型修正结果满足工程需求，修正后频率差均在10%以内，振型相关性比初始值有些许增大，均大于0.8。
2) 基于修正后的有限元模型，开展结构在锤击载荷作用下的模态动力学计算，与锤击试验结果作比较，两者振动加速度响应在时域量级和频域幅值匹配较好，模型可靠性高，可用作下一步复杂载荷作用下的结构动力学响应预测。
3) 搭建的结构有限元模型修正平台可实现优化流程自动化，提升设计分析效率，使设计人员专注于设计结果对工程问题的解决。
