摘要
针对广泛地应用在航空航天箭体结构中的螺栓法兰连接结构的动力学响应特征会显著受到螺栓预紧力大小影响的问题，在螺栓法兰连接结构不同工况下螺栓预紧力出现差异时，研究箭体连接结构的静动力学响应特征。通过对特定的筒壳螺栓法兰连接结构，基于有限元分析和敲击试验来提出一套考虑螺栓预紧力差异的箭体连接结构动力学响应分析方法。并由此研究螺栓松动时的响应特征，基于互补集合经验模态分解算法设计一种新的螺栓松动探测流程，可以量化螺栓发生松动时连接结构的中高频段频响幅值升高的特征，能无损快速地识别松动工况的发生。
Abstract
In order to address the problem that the dynamic response characteristics of the bolted flange connection structure, which is widely used in aerospace rocket structures, are significantly affected by the magnitude of bolt preload, the static and dynamic response characteristics of the rocket connection structure were studied taking into account the discrepancy of bolt preloads under different conditions. A set of methods for analyzing the dynamic response of rocket body connection structures considering the discrepancy of bolt preloads was proposed based on finite element analysis and knock tests for a specific barrel-shell bolted flange connection structure. By investigating the response characteristics of the connection structure with a loose bolt, a new bolt loosening detection procedure based on the complementary ensemble empirical mode decomposition algorithm is proposed which can quantify the characteristics of the increase of the frequency response amplitude in the middle and high frequency bands of the connection structure when the loosening happens, and can quickly identify the occurrence of loose conditions without damage.
箭体级间段对接面上多采用螺栓法兰连接，由于此处无法保持结构整体连续性，导致动力学响应在不同工况和装备条件下呈现非线性，影响结构整体安全性及可靠度，因此有必要深入研究螺栓法兰连接结构模型，分析其结构动力学响应。
建立螺栓法兰连接结构的精细有限元分析模型，可以有效实现这类结构的力学响应分析。比如，Kim等对比螺栓连接结构有限元模型和梁模型的差异，将螺栓连接结构精细有限元模型引入实际工程结构中并验证其有效性[1]；Jamia等通过螺栓连接结构精细有限元模型研究连接面微滑移力学行为，建立了等效预测模型及参数识别框架[2]；Beaudoin等为止口法兰结构建立机理分析模型和精细有限元模型，通过仿真分析和实验验证机理分析模型有效性，确认了双线性刚度模型适用范围[3]；Wu等通过精细有限元模型分析螺栓法兰连接结构在拉弯扭工况下的非线性表现行为[4]。但是，在动力学响应分析方面，精细有限元模型始终面临计算效率低的问题，通常需要建立满足精度要求的简化模型来解决。
在筒壳连接结构的简化动力学建模方面，国内研究者走到了前列。栾宇等通过理论分析及静力试验将拉压不同刚度模型引入螺栓法兰连接结构建模中，对航天器级间段的横纵耦合动力学响应特征实现了有效表征，由此解释了产生横纵耦合运动的机理[5]。芦旭等在此基础上，研究航天器对接面带剪力销时整体结构受弯扭剪耦合作用下的动力学响应，分析剪力销侧向倾角等对连接结构内力的影响[6-7]。蒋国庆等为箭体连接结构设计了参数化有限元模型，有效分析了不同几何参数、网格参数等对计算精度的影响[8-9]。田彤辉等对冲击载荷下某弹体螺栓法兰连接结构失效问题进行了仿真分析和破坏性试验[10-12]。潘嘉诚等分析航天器局部螺栓法兰连接结构对接面的分离效应，并进行近似理论建模，获得较高计算精度[13]。上述针对连接结构的简化建模与分析，在整体结构失效工况分析及微观滑移行为研究方面，都假设预紧力施加状态理想，尚未针对实际工况中不同螺栓位置预紧力差异对结构整体及局部力学特性的影响开展研究。
因此，本文基于某等效的箭体连接结构模型，采用精细有限元分析的形式，以预应力的形式施加螺栓预紧力；分析不同位置螺栓预紧力出现差异时，连接结构非线性静力与动力学响应特性；基于响应分析结果以模态分解算法设计一种信号处理流程，来识别螺栓预紧力下降即螺栓松动的情况。
1 分析方法
1.1 箭体连接结构分析方法
通常箭体连接结构的分析采用简化模型计算校核静刚度和动响应，但这往往是对整体结构的分析，忽略了预紧力变化的影响。本文综合了静力计算的特点和动响应的特点并考虑预紧力影响，针对箭体连接结构设计了如图1所示的分析流程。预紧力对于连接结构的刚度有着很大的影响，通过静力学分析获得箭体连接结构的刚度特点，进而在动力学方程式（1）中考虑K的变化；通过动响应分析得到连接结构加速度信号的时域和频域特征及预紧力对其的影响，在时域和频域上的特征有些并不特别直观，进而需要将这些特征以特定分析方法量化，该量化方法通过模态分解算法和离散积分来处理动力学响应中的加速度信号；量化后可以识别箭体连接结构的螺栓预紧力下降即螺栓松动的情况，可以此为依据调整连接结构。
(1)
式中：y表示系统位移，其关于时间的微分用上标点表示；M，C为系统的质量阵及阻尼阵；K（y）表示与y相关的系统刚度；F为系统外载荷。
图1分析方法流程图
Fig.1Flow chart of analysis methods
1.2 螺栓松动识别方法
单个螺栓预紧力减小对连接结构冲击响应有着显著影响，但实际工况响应信号会受到不同来源的干扰，不利于准确探测螺栓松动，因此基于螺栓松动响应特性，借助特定频域信号分析实现松动探测。主要思路是先提取时域信号即式（1）中的，再采用互补集合经验模态分解（complementary ensemble empirical mode decomposition，CEEMD）方法进行信号分解，分解后采用无限脉冲响应（infinite impulse response，IIR）递归滤波器提取部分高频特征信号，将这部分信号进行快速傅里叶变换（fast fourier transform，FFT）后，进行频域上的离散积分计算量化参数指标，实现松动状态识别。
1.2.1 CEEMD方法原理
经验模态分解（empirical mode decomposition，EMD）是一种针对非线性、非平稳信号的自适应时频分析方法，能够将不同时间尺度的信号分解为相对平滑的本征模态函数（intrinsic mode function，IMF）。通过EMD处理的信号被分解为n个基本分量fi（t）和1个残余项rn（t）之和，即
(2)
为了避免分解过程中IMF不连续造成的模态混叠，Huang等[14]通过在EMD方法中加入高斯白噪声提出改进经验模态分解（ensemble empirical mode decomposition，EEMD）方法，使信号在不同时间尺度上具有连续性，如式（3）所示。
(3)
其中，h（t）为信号序列，Iij为第i次加入白噪声分解后所得的第j个I，c为残余项，M为加入白噪声的次数。
虽然EEMD引入白噪声辅助分析改善了模态混叠问题，但噪声对原始信号有一定程度的破坏和干扰作用。因此，Yeh等[15]提出CEEMD方法，引入独立同分布且完全负相关的互补噪声，在重构信号时，基本消除了冗余噪声，原理为：
(4)
式中，S为原始信号序列，N为白噪声，M1、M2为正负互补的白噪声。
1.2.2 松动探测流程
设计松动探测方法流程如图2所示，对响应信号通过CEEMD进行本征模态分解；为了提取特征，在高频本征模态信号下进行IIR高通滤波；对分离后的信号进行FFT变换后，利用数据点在频域上进行离散积分；设计量化参数指标ΔSi如下（指标值大于零时判定出现松动）：
(5)
(6)
其中，C为经过CEEMD分解后的振动时域分解信号，Ci表示C的第i阶分解信号，本文采用IMF一阶信号即包含高频信号成分最多的信号（针对不同螺栓法兰连接结构可采用不同阶次的信号），下标after-IIR-FFT表示经过IIR滤波后进行FFT的结果；Si的右上标代表不同的螺栓状态下的取值；积分号表示在频域对FFT后的信号进行积分。为使设计量化参数指标具有更好的普适性，式（6）中ΔSi为无量纲比值。
图2探测方法流程图
Fig.2Flow chart of detection methods
2 基于有限元仿真的方法验证
为说明上述分析方法的有效性，本节针对某型号箭体级间段结构，建立其缩比尺寸等效模型，基于精细有限元分析开展方法验证。
考虑到真实级间段连接结构特征复杂（如图3所示），其主要承力结构为典型的含有螺栓连接的薄壁筒壳结构。为方便实现仿真分析与实验室实验，保留此主要承力结构特征后，对该结构适当简化并忽略螺纹影响（如图4所示）。根据相关箭体级间段参数背景，选择两个尺寸和材料属性一致的航空铝合金筒壳结构以法兰面处的六个螺栓连接，其长度为340 mm、壁厚为4 mm、直径为295 mm，弹性模量为70 GPa、泊松比为0.3、密度为2 700 kg/m3；法兰面螺栓孔径为8.4 mm，连续法兰厚为10 mm。高强钢螺栓弹性模量为210 GPa、泊松比为0.3、密度为7 850 kg/m3，头部直径为12 mm、螺栓杆直径为5.7 mm。利用Abaqus软件建立有限元模型，因螺纹简化可令螺栓螺母组成一体，以Bolt load形式施加预紧力，采用C3D8R类型的六面体网格。螺栓与法兰接触设置为面接触，接触摩擦系数设置为0.2，采用罚函数摩擦公式。
图4有限元模型
Fig.4FEM Model
2.1 理想预紧下的仿真分析
2.1.1 理想预紧下的静力分析
静力分析加载方式如图5所示，下端筒壳底部固定约束，先加载预紧力使螺栓紧固，再将筒壳顶面耦合在中心点处施加20 kN轴向载荷，分别进行压缩和拉伸，计算连接面上下两个测点间距变化，得到如图6所示的连接结构静力分析的载荷位移曲线。
图5静力计算模型
Fig.5Static calculation model
图6静力分析的载荷位移曲线
Fig.6Load displacement curve in static analysis
由图6可知，静力分析的载荷位移曲线呈现明显的非线性，结构在受拉和受压时的刚度差异较大。当螺栓预紧力从1 kN变化到5 kN时，若拉伸载荷较小则连接结构局部的刚度随预紧力增大而增大。而随着拉伸载荷增大到一定程度后，曲线斜率保持不变，继续增加拉伸载荷将引发连接面局部分离，而局部分离状态下螺栓预紧力对刚度影响很小。由于这时拉伸刚度来自螺栓自身，不同初始预紧力下保持不变，因此从图6中可以看到载荷增大到一定程度后呈现三条平行直线。而当施加压缩载荷时，由于主要为筒壁结构承压，预紧力的变化几乎不影响压缩刚度，因此各条曲线几乎完全重合。
图7冲击计算模型
Fig.7Impact calculation model
2.1.2 理想预紧下的冲击响应分析
航天器设计的静载强度一般偏向安全区间，且容易进行校核和实验比对，而动载荷是造成航天器连接结构失效的重要形式，在临界工况突变时会造成冲击载荷。因此，如图7所示本节对六螺栓均布的连接结构依前文施加5 kN预紧力后，采用隐式动力学分析，计算幅值为1 kN、脉宽为2 ms的结构顶端脉冲冲击载荷下的短时加速度响应（如图8所示）。
图8横向冲击下的加速度响应
Fig.8Acceleration response under transversal impact
横向冲击作用下，不仅产生横向加速度响应，也诱发了同量级的纵向加速度响应，这与栾宇等[5]提出的简化模型分析结果一致，即产生横纵耦合效应，其由连接面非线性刚度造成。对图8所示的时域响应信号进行快速傅里叶变换得到幅频曲线，如图9所示，可知响应主要集中在低频率段，其一阶弯曲模态频率为340 Hz。随着整体预紧力从5 kN降低到1 kN，此一阶频率几乎没有变化，但是其响应幅值降低；预紧力降低造成的主要频率变化在[500， 1 000] Hz频率区间段，降低预紧力使该频率段峰值向低频率区间偏移，且幅值增大。
图9横向冲击响应幅频曲线
Fig.9Transversal impact response amplitude and frequency curves
除上述频率特征之外，连接面测点处在冲击载荷下的加速度响应峰值也显著受到预紧力影响，如表1所示。在横向冲击载荷下，随着预紧力降低，横向加速度响应峰值显著增大，而纵向加速度响应峰值变化相较并不显著，这说明此时降低预紧力对结构的主要影响体现在横向响应上；由预紧力降低带来的连接面处横向约束作用下降，对螺栓和箭体局部连接结构的抗弯剪能力会造成负面影响。
表1横向冲击下的结构峰值响应
Tab.1 Peak response of structure under transversal impact
此外，当依照图5所示的轴向静载作用方式来施加纵向冲击激励时，如表2所示横向加速度响应峰值远小于纵向加速度响应峰值。这时，螺栓预紧力减小对纵向加速度响应峰值的影响仍不显著，虽然横向加速度响应峰值随着预紧力减小还是会增大，但是由于其量级较小，对整体结构响应的影响远没有受到来自横向冲击时的效果明显。
表2纵向冲击下的结构峰值响应
Tab.2 Peak response of structure under longitudinal impact
由此可知，预紧力变化对于加速度响应的影响主要体现在横向冲击工况时。因此，研究螺栓松动引发特定螺栓预紧力减小的问题时，围绕横向冲击工况下的连接结构响应开展分析为宜。
2.2 单个螺栓松动时的仿真分析
本节针对单个螺栓预紧力从5 kN下降到0.3 kN时，分析连接结构在静载和冲击载荷下的响应变化，以便为局部松动预测提供理论支撑，以避免极端载荷情况下因局部连接预紧力不足造成严重后果。
2.2.1 螺栓松动下的静力分析
单个螺栓松动静刚度如图10所示，单个螺栓松动时，施加2.1节中的轴向拉伸载荷会造成此处提前进入局部分离。对比图11中给出的Mises应力云图，能够看到紧固螺栓侧与松动螺栓侧的局部应力云图差异显著，螺栓预紧力降低会使此螺栓附近连接面部分发生局部翘曲，产生非平直的形变，即导致局部刚度损失。
图10单个螺栓松动静刚度
Fig.10Single bolt loose static stiffness
2.2.2 螺栓松动下的冲击响应分析
既然某处螺栓预紧力降低会造成局部刚度损失，那么对连接结构整体的冲击响应必然产生影响，因此本节采用隐式动力学计算流程，分析人为设置松动即预紧力大幅度降低后的结构动力学响应。采用与2.1节相同的冲击载荷与响应测点，对加速度响应进行短时傅里叶变换后得到如图12所示响应频谱，可以发现，在不同时刻下，当单个螺栓发生松动时，一阶频率变化不明显，但在[500，1000]Hz区间内会出现更多的频响峰值，这一特征预示螺栓松动的存在。
图11静载状态下紧固螺栓和松动螺栓应力云图比较
Fig.11Comparison of stress clouds near fastening bolt and loosening bolts under static load
横向冲击下考虑松动的加速度响应峰值如表3所示。施加横向冲击载荷时，随着单个螺栓松动程度的增加，测点处无论是横向响应峰值还是纵向响应峰值均有明显增加，比降低整体预紧力的响应增加幅度还要明显，可见局部螺栓松动会显著放大连接结构的动力学响应，对结构安全造成不利影响，提前探测这种松动现象十分必要。
图12不同时刻的加速度响应频谱
Fig.12Frequency spectrum of acceleration response for different time points
表3横向冲击下考虑松动的加速度响应峰值
Tab.3 Peak response of structure with loosening under transversal impact
3 松动探测
对与前文相同的横向冲击作用下的连接结构模型，分别在不含松动螺栓和含单个松动螺栓（指单个螺栓松动到0.3 kN即典型的完全松脱失效）两种情况下计算加速度时程响应来验证1.2节中方法的有效性，由前文分析取得本结构的IIR高通滤波为500 Hz以上频率。
图13和图14分别为无松动和单个螺栓松动时的IMF信号分解结果，可以发现：在含有松动螺栓的一阶本征模态信号即IMF第1阶（图13和图14中的IMF1）含有明显的高频成分，后续阶次中并没有明显的特征性成分，因此以IMF1分解信号作为设计流程中的指标有一定代表性。按照式（5）和式（6）进行量化参数指标计算，结果见表4。由表4可知，随着螺栓松动程度加深（即预紧力下降程度加大），振动时域信号中高频成分增多，Si值显著增大，ΔSi同时增大，该比值反映松动程度加深。
图13无松动IMF
Fig.13IMF without loosening
图14含单个螺栓松动的IMF
Fig.14IMF with single loose bolt
表4松动指标（有限元）
Tab.4 Loose indicators（FEM）
4 敲击实验验证
为进一步验证方法的有效性，如图15所示采用DH5920振动信号采集仪及加速度传感器进行力锤冲击实验，实验件尺寸与前文一致，材质为7系铝合金，采用M4高强钢螺栓连接，连接面处采用均布六螺栓分布装配。用压电式螺栓力传感器标定预紧力，六个螺栓初始预紧力均为5 kN，松动工况下测点处螺栓预紧力人为放松至3 kN及0.3 kN，分别模拟初期松动和几乎完全松动。在测点对侧的筒段顶端加载接近1 kN横向冲击力，采集横向加速度信号，进行前面设计的信号处理流程，进行CEEMD模态分解及量化参数计算，验证松动识别流程。
图15实验设备示意图
Fig.15Experiment equipment
实验中无松动和单个螺栓松动时的IMF信号分解结果分别如图16和图17所示，实验中的松动指标见表5。实验结果表明：与仿真结果一样，图16和图17所示的IMF1信号包含高频成分可用于松动检测。表4和表5给出的Si具体数值方面，仿真与实验之间存在误差，这是因为有限元模型较为理想，接触面均为平面，且只考虑理想接触面上线性摩擦关系，而实际实验工况不仅摩擦关系复杂，在连接接触面处也有复杂的局部翘曲和磨损情况存在。由于有限元采用隐式求解，其采样率不可避免地与实际设备有一定差异，这种采样的波动在高频段的能量影响更为明显，以上原因造成的高频能量差异导致了频域数值积分的误差，这也是表5中误差的来源。虽然有限元与实验结果的绝对数值有一定误差，但二者的定性趋势规律是一致且有效的，螺栓松动程度增大导致频域能量分布变化，使得量化指标同步增大，即Si和ΔSi越大，螺栓的松动程度越高，ΔSi超过1后螺栓的松动均已经相当严重，因此可基于该规律识别连接结构状态。另外，相比具有量纲的Si绝对数值，本文更关注无量纲量ΔSi，该指标结果与有限元结果相近且定性规律一致。由此可知，通过实验的指标趋势性，验证了本文方法的有效性，具体应用中可基于该流程来对实际连接结构进行标定后快速无损地识别螺栓预紧力差异状态。后续的工作将致力于开发与预紧力变化映射更为精准的、定量更为准确的、判断区间更小的识别方法。
图16无松动IMF（实验）
Fig.16IMF without loosening (experiments)
图17含单个螺栓松动的IMF（实验）
Fig.17IMF with single loose bolt (experiments)
表5松动指标（实验）
Tab.5 Loose indicators（experiments）
5 结论
本文基于螺栓法兰连接结构有限元模型分析了螺栓预紧力对静动力学特性的影响，并根据单个螺栓松动下整体及局部结构冲击响应特征，设计了松动检测流程，得到以下结论：
1）螺栓连接结构的局部松动会使得横向约束显著变低，且会出现更多的中高频频率偏移现象，该频段频响幅值升高；
2）实验验证了本文基于CEEMD的特定本征模态离散积分方法，能无损快速地识别松动工况的发生，具有一定的工程指导意义。
