范玉川¹, 赵春雨¹, 鲁艳², 张义民¹
1. 东北大学 机械工程与自动化学院, 辽宁 沈阳 110819;
2. 郑州信大先进技术研究院, 河南 郑州 450001
收稿日期：2018-04-19
基金项目：国家自然科学基金资助项目(51775094)。
作者简介：范玉川(1988-)，男，河南新乡人，东北大学博士研究生;
赵春雨(1963-), 男, 辽宁黑山人，东北大学教授，博士生导师;
张义民(1958-), 男, 吉林长春人，东北大学"****奖励计划"特聘教授，博士生导师。

摘要：推导出多自由度动力学方程的Wilson-θ数值算法显式表达形式, 进而提出了一种显式Wilson-θ的动载荷识别算法.该算法避免了Wilson-θ算法的隐式迭代形式的迭代误差, 在拥有显式算法特性的同时具备隐式算法的特性.当θ取合适的值时, 该算法是无条件稳定的.通过悬臂梁的算例和实验对算法的识别效果进行了验证, 并与传统的状态空间法的识别结果进行了对比.结果表明:该算法不仅能够对矩形载荷、谐波载荷和随机载荷进行准确地识别, 并且比状态空间法的识别精度更高.
关键词：Wilson-θ法显式表达载荷识别无条件稳定状态空间法
Research on Dynamic Load Identification Based on Explicit Wilson-θ Method
FAN Yu-chuan¹, ZHAO Chun-yu¹, LU Yan², ZHANG Yi-min¹
1. School of Mechanical & Automation, Northeastern University, Shenyang 110819, China;
2. Zhengzhou Xinda Institute of Advanced Technology, Zhengzhou 450001, China
Corresponding author: ZHAO Chun-yu, E-mail: chyzhao@mail.neu.edu.cn
Abstract: The explicit expression of Wilson-θ numerical algorithm for multi-dofs(degree of freedoms) dynamic equation is derived, as well as an explicit Wilson-θ dynamic load identification algorithm is proposed, avoiding the iteration error while keeping the characteristics of the implicit iteration algorithm. The algorithm is unconditionally stable when applying appropriate θ value. The recognition effect of the algorithm is verified by an example and an experiment of a cantilever beam, and the results were compared with those from the traditional state space method. The results show that the algorithm not only can accurately identify the rectangular load, the harmonic load and the random load, but also has higher recognition accuracy than state space method.
Key words: Wilson-θexplicit formularload identificationunconditionally stablestate-space method
在机械系统设计过程中, 动载荷是机械结构进行疲劳分析以及进行可靠性计算的基本依据.但是, 在动态系统中, 有些力是很难直接测量的, 特别是结构系统内部各部件之间的相互作用力, 在难以直接测量的情况下, 就需要通过逆动力学的分析技术来得到这些力, 因而开展动载荷识别技术的研究是很有必要的.
动载荷识别方法主要分为频域法^[1-2]和时域法^[3]两大类.由于时域法能够识别各种类型的载荷, 识别精度高, 并且其识别结果具有明确的物理意义, 因而, 时域法越来越受到专家学者们的青睐.Liu等^[4]推导了一种基于Newmark-β法的动载荷识别方法, 将传统的隐式Newmark-β算法转化为Ax=b方程解的显式形式, 与隐式的Newmark-β法具有相同的特点, 并与状态空间法进行了对比, 表明该算法具有明显的优势.Li等^[5]提出了一种基于二阶泰勒级数展开的时域动态结构载荷识别方法, 该算法将响应表示为一种泰勒级数的逼近形式, 推导出一系列公式, 并建立了系统响应、系统特性和输入激励相结合的显式离散方程.Liu等^[6]提出了一种新的时域动态Galerkin算法, 将形状函数作为加权函数, 建立了前向模型TDGM, 与传统的格林函数法相比, TDGM能有效地克服噪声的影响, 提高动态载荷识别的精度.在求逆运算中, 载荷识别结果通常对结构模型中的响应和误差测量中的噪声非常敏感^[7].张方等^[8-9]对复杂结构的载荷识别进行了研究, 推导了一种基于广义正交多项式特征技术的动载荷识别模型, 可以在一定精度范围内通过有限的测量点信息对无限未知量的分布动载荷进行识别.Allen等^[10]采用一种SWAT方法, 可以同时识别出冲击型载荷和稳态动态载荷.
本文将对Wilson-θ法进行变换, 推导其显式表达形式, 提出一种基于显式Wilson-θ法的动载荷识别算法, 并通过算例和实验验证了该算法的有效性.
1 动载荷识别算法线性阻尼结构中, 多自由度结构的动力学方程可以表示为

(1)

其中:M, C和K分别表示质量矩阵、阻尼矩阵以及刚度矩阵; P_t是作用在结构上的外加载荷; , 和x_t分别表示加速度响应、速度响应以及位移响应.本文假定阻尼为瑞利阻尼:

(2)

其中, α₁和α₂表示阻尼系数.
1.1 Wilson-θ法的显式表达Wilson-θ法假定在[t, t+θΔt](θ≥1)的时间间隔内, 加速度呈线性变化, 如图 1所示.
图 1(Fig. 1)

图 1 Wilson-θ法模型Fig.1 Model of Wilson-θ

令δ为自t时刻开始的时间变量, 适用于0≤δ≤θΔt, 由线性加速度的假设可知, 在适用范围内的加速度为

(3)

积分后可得

(4)

(5)

若δ=θΔt, 由式(4)和式(5)可得t+θΔt瞬时的速度和位移：

(6)

(7)

联立式(6)和式(7), 可以用t+θΔt时刻的位移表示t+θΔt时刻的加速度和速度, 即

(8)

(9)

则t+θΔt时刻的动力方程可以表示为

(10)

式中:

(11)

将式(8)和式(9)以及式(11)代入式(10), 即得关于x_t+θΔt的求解方程为

(12)

式中:

(13)

(14)

求解方程式(12), 得x_t+θΔt.再把x_t+θΔt代入式(8)就可获得.在式(3)中取δ=Δt, 同时将式(8)代入, 可得

(15)

取δ=Δt, 将式(3)分别代入式(4)和式(5), 有

(16)

(17)

由式(12)~式(14)可得

(18)

将式(18)代入式(15), 已知恒等式:
可得

(19)

其中:
将式(19)代入式(16)可得

(20)

其中:
将式(19)代入式(17)可得

(21)

其中:
在此用和P_i表示第i时刻式(19)~式(21)中的位移、速度、加速度和激励, 并将其表示为矩阵形式:

(22)

在第i时刻的响应可以表示为

(23)

式中, 两个指数i和j分别表示相应矩阵的幂.式(23)是一种新型的Wilson-θ法的显式表达, 每一个时间步中位移、速度、加速度响应可被同时求解出来, 而通常的方法是利用每个时间步的迭代算法来计算, 显然这种显式的算法更有优势.
1.2 基于显式Wilson-θ法的动载荷识别令:
则式(23)可以改写为

(24)

令:

(25)

(26)

式(29)可以写成从时间段1到nt上的矩阵的卷积形式:

(27)

其中:
式(27)可以改写成载荷识别的表达式

(28)

对于一个给定的系统, H是常数, Y可以从系统测量的响应中得到, 考虑到式(28)可能存在不适定性, 可以通过Tikhonov正则化方法来确定F, 在阻尼最小二乘意义下对目标函数进行优化.

(29)

其中, λ是正则化参数, 它的值可以通过L曲线法来确定.基于显式Wilson-θ法的动载荷识别方法可以归纳为以下步骤:
1) 建立系统的有限元模型, 确定系统的质量矩阵M、刚度矩阵K和阻尼矩阵C;
2) 选取时间步Δt和Wilson-θ法的参数θ, 参数θ的值可以取为1.4, 从而保证算法的无条件稳定(在Wilson-θ法中, 只要θ的值大于1.37, 那么该算法就无条件稳定, 但是θ取的过大, 截断误差增大, 精度会下降).
3) 为递推式(22)计算矩阵A₁, A₂, A_d, A_v, A_a，B₁, B₂, B_d, B_v, B_a, C₁, C₂, C_d, C_v和C_a;
4) 计算式(25)中的矩阵H_k(k=0, …, nt-1)和式(28)中的组合矩阵H;
5) 通过实验或仿真计算获得系统的响应数据和x;
6) 利用式(29)中的Tikhonov正则化方法来确定F.
2 仿真算例及抗噪性能研究如图 2所示, 一矩形截面悬臂梁长0.64 m, 截面尺寸为:宽0.056 m、高0.008 m, 弹性模量E=200 GPa, 材料密度为7 840 kg/m³, 该悬臂梁被划分为18个单元, 共19个节点.
图 2(Fig. 2)

图 2 悬臂梁结构图Fig.2 Structure of cantilever beam

在第10个节点上施加载荷, 为了尽可能全面地验证本文载荷识别算法的准确性, 分别施加方波载荷、谐波载荷和随机载荷, 考虑到噪声的影响, 在响应数据中加入5%的随机噪声, 选取第12个节点的响应信息进行载荷识别, 时间间隔Δt=0.001 s.
选取第12个节点上的数据进行载荷识别计算, 本文算法与状态空间法进行对比, 如图 3所示, p₂为方波载荷、p₃为谐波载荷、p₄为随机载荷.两种方法对于这三种不同类型的载荷都能够比较准确的识别, 但是状态空间法的识别效果明显要差一些, 而基于显式Wilson-θ法的载荷识别方法识别得很好, 对实际载荷的还原度较高.
图 3(Fig. 3)

图 3 载荷识别结果Fig.3 Results of load identification (a)—p2；(b)—p3；(c)—p4.

为便于量化分析基于显式Wilson-θ法的动载荷识别算法的识别精度, 可引入式(30)计算误差:

(30)

式中:Error表示相对识别误差; F_id表示识别得到的外激励; F_real表示真实的外激励.运用不同的方法, 在不同的采样频率和采样时间下研究动载荷识别的误差, 结果列于表 1中.
表 1(Table 1)

表 1 两种载荷识别方法的识别误差对比Table 1 Comparison of load identification errors with different load identification methods % 载荷 Wilson-θ法 状态空间法 p2 3.21 5.04 P3 0.61 3.84 P4 3.30 5.16

表 1 两种载荷识别方法的识别误差对比 Table 1 Comparison of load identification errors with different load identification methods

数据表明, 本文提出的基于显式Wilson-θ法的载荷识别方法的识别精度明显高于基于状态空间法的载荷识别方法的识别精度, 并且两种算法都显示谐波载荷的识别精度要远远好于方波载荷和随机载荷的识别精度.
图 4显示了随着噪声的增大载荷识别精度的变化情况, 两种方法的识别误差都是随着噪声的增大而增大, 但本文算法随噪声的增大识别误差变化较小, 表明本文算法抗噪性能更好、识别误差更小.
图 4(Fig. 4)

图 4 载荷识别误差随噪声变化情况Fig.4 Variation of load identification error with noise

3 实验验证实验用的悬臂梁模型参数与仿真悬臂梁模型参数一致, 选用YE6251振动力学实验系统.
激励载荷为谐波载荷, 激励频率设置为10 Hz, 算法中采用有限元法进行建模, 代入悬臂梁的相关参数, 选用测量得到的第11个节点上的响应数据进行载荷识别, 识别结果如图 5所示.由图可知, 该算法对实验数据的识别较准确, 识别的载荷曲线与实际的载荷曲线基本吻合.
图 5(Fig. 5)

图 5 实验识别结果Fig.5 Results of experiment

由于实际工程中响应数据的测量位置往往是有限的, 不能够随意选取, 因而在本实验中分别选取第4、第7、第10、第13和第16个节点的测量数据进行载荷的识别, 同样利用式(30)进行识别误差的计算, 计算结果分别为5.67%，5.79%，5.38%，5.86%，5.98%.
从结果可以看到, 选取不同的测量点, 识别得到的载荷与实际载荷的误差没有太大差别, 也就是说测量点的选取对识别结果影响不大; 该算法能够适应实际工程的应用环境.
4 结论1) 利用隐式Wilson-θ法, 推导出了它的一种显式表达形式, 进而提出了一种基于显式Wilson-θ法的载荷识别算法, 该算法既具有显式算法的优势, 同时具有隐式算法的特点.
2) 通过算例验证了本文算法在识别矩形波载荷、谐波载荷和随机载荷方面的识别精度都高于状态空间法的识别精度, 并且在算例中加入了5%的随机噪声.实验结果表明:算法不仅能够对矩形载荷、谐波载荷和随机载荷准确识别, 并且识别精度比状态空间法的识别精度更高.
参考文献

[1]	Jie H, Zhang X. An optimization method of load iden-tification in frequency domain[J].Noise & Vibration Control, 2009, 29(6): 34–36.
[2]	He Z C, Lin X Y, Li E. A novel method for load bounds identification for uncertain structures in frequency domain[J].International Journal of Computational Methods, 2017(3): 1850051.
[3]	Law S S, Chan T H T, Zhu Q X. Regularization in moving force identification[J].Journal of Engineering Mechanics ASCE, 2001, 127(2): 136–148.DOI:10.1061/(ASCE)0733-9399(2001)127:2(136)
[4]	Liu K, Law S S, Zhu X Q, et al. Explicit form of an imply-CIT method for inverse force identification[J].Journal of Sound and Vibration, 2014, 33: 730–744.
[5]	Li X W, Deng Z M. Identification of dynamic loads based on second-order taylor-series expansion method[J].Shock and Vibration, 2016(2016): 1–9.
[6]	Liu J, Meng X, Jiang C, et al. Time-domain Galerkin method for dynamic load identification[J].International Journal for Mumerical Methods in Engineering, 2016, 105: 620–640.DOI:10.1002/nme.v105.8
[7]	Qian B, Zhang X, Wang C, et al. Sparse regularization for force identification using dictionaries[J].Journal of Sound and Vibration, 2016, 368: 71–86.DOI:10.1016/j.jsv.2016.01.030
[8]	张方, 秦远田, 邓吉宏. 复杂分布动载荷识别技术研究[J].振动工程学报, 2006, 19(1): 81–85. ( Zhang Fang, Qin Yuan-tian, Deng Ji-hong. Research of identification technology of dynamic load distributed on the structure[J].Journal of Vibration Engineering, 2006, 19(1): 81–85.DOI:10.3969/j.issn.1004-4523.2006.01.014)
[9]	徐菁, 张方, 姜金辉, 等. 运用数值迭代的动载荷识别算法[J].振动工程学报, 2014, 27(15): 702–707. ( Xu Jing, Zhang Fang, Jiang Jin-hui, et al. An algorithm of dynamic load identification based on numerical iteration,[J].Journal of Vibration Engineering, 2014, 27(15): 702–707.)
[10]	Allen M S, Carne T G. Delayed multi-step inverse structural filter for robust force identification[J].Mechanical Systems and Signal Processing, 2008, 22(5): 1036–1054.DOI:10.1016/j.ymssp.2007.11.022