運用數值迭代的動載荷識別算法

徐 菁， 張 方， 姜金輝， 浦玉學， 蔣 祺

(南京航空航天大學機械結構力學及控制國家重點實驗室， 江蘇 南京 210016)

引 言

在動力學問題研究中，準確知道結構的動載荷是非常重要的。動載荷識別技術主要發展了頻域和時域兩類載荷識別方法，通過時域法識別的載荷是一個載荷時間歷程，可以直觀了解載荷隨時間的變化規律，更適用于工程實踐。張方等用廣義正交多項式特征技術建立載荷識別模型[1],解決了復雜結構的分布動態載荷識別問題。秦遠田采用矩量法基函數來擬合未知動載荷[2]，將載荷識別轉化成為在廣義正交域中求解基函數擬合系數的問題。Allen和Carne采用SWAT方法對沖擊載荷進行識別[3]，該方法可以同時識別穩態動態載荷與沖擊型載荷。Gunawan等提出采用正則化二次樣條函數來擬合沖擊載荷[4]，利用基于L曲線的TSVD方法求解動載荷。Chen和Lee利用模態疊加正交技術獲得梁的狀態方程[5]，通過在線基于自適應權系數的遞歸算法來識別作用于橋梁結構的動態載荷，該方法具有較好的跟蹤性能和誤差估計特性。Mao和Guo將激勵力描述為基函數與權系數之積的和[6]，為避免得到不穩定的解，在求解權系數的過程中，他們采用了基于GCV準則的Tikhonov正則化方法，在實驗室中通過該方法來識別無人機機翼的瞬態沖擊激勵。韓旭等將動態載荷表示為一系列的脈沖或者階躍函數的疊加[7]，并運用零相位濾波器、正則化技術和優化策略來實現載荷的穩定重構。高峰、李德武在Newmark法求解運動平衡方程的基礎上[8]，推導了求解振動荷載的公式，提出由結構對振動荷載的反應求振動荷載時程的有限元方法。姜金輝等根據譜分解的思想[9]，提出多點任意相關的隨機載荷識別方法，指出病態發生的原因和相應對策，提高了識別精度。Lin采用了廣義卡爾曼濾波器以及基于常權系數的遞歸最小二乘法[10]，對單自由度非線性系統進行載荷識別，該方法能夠較準確地估計正弦激勵和多種沖擊激勵，并具有較好的抗噪性能。其后，他采用線性化卡爾曼濾波器以及基于自適應權系數的遞歸最小二乘法來識別非線性系統的時域載荷[11]，并將兩種方法進行了對比，驗證了后者的識別精度更高。Lee提出采用智能模糊權系數以改進自適應權系數的缺點[12]，識別結果表明智能模糊權收斂性好，能夠有效降低測試誤差、模型誤差和儀器偏差的影響，且跟蹤能力好。Lee和Chen在假設回轉機械為剛體的基礎上[13]，提出采用智能模糊權方法估計動態載荷。

本文對Wilson-θ反分析法的推導過程做了理論分析，找出其算法發散的原因。針對該原因提出了一種新的數值迭代修正方法，從數學上證明了修正的可行性，并且運用仿真和實驗結果證明了修正算法的收斂。

1 Wilson-θ反分析法研究

Wilson-θ法是線性加速度法的擴展，在t至t+θΔt的時間區間內利用了線性加速度假設，并且計算得到系統的動力學方程。經嚴格數學證明，Wilson-θ法具有良好的數值穩定性，取θ>1.37即可保證結果無條件收斂，實際計算中一般取θ=1.4。文獻[14]運用反問題的有限元方法[8]，推導得到一種反分析法來進行載荷識別。

在t+θΔt時刻的系統擬靜力方程為

(1)

式(1)可分開寫為

(2)

(3)

假設在每一微小的時間段內，可視為線彈性問題，根據疊加原理，有

x(t+θΔt)=x′(t+θΔt)+x″(t+θΔt)

(4)

(5)

引入位移系數，得到如下方程

(6)

(7)

(8)

所以，式(2)的解滿足下列線性關系

(9)

由式(4)得

xk(t+θΔt)-(x″(t+θΔt))k=v

(10)

(11)

求出λt+θΔt后，根據式(5)，可以求出t+θΔt時刻的動載荷值qi(t+θΔt)，由于載荷qi(t)已知，故可得

qi(t+Δt)=(qi(t+θΔt)-qi(t))/θ+qi(t)

(12)

2 利用二分法迭代修正

對于單點輸入模型，可以采用二分法迭代對每一步的識別力進行修正。運用二分法迭代，必須先確定含根區間。由Wilson-θ法可知

(13)

(14)

(15)

于是

(16)

(17)

(18)

化簡得

(19)

對gi(f)求導

(20)

首先由Wilson-θ反分析法，得到識別力的時間序列q(t)。對于時間點tm，設其對應的載荷為q(tm)，令

a1=-rq(tm)，b1=rq(tm)

式中r為區間放大倍數，且滿足g(a1)g(b1)<0，則[a1,b1]為初始含根區間(假設q(tm)>0)。

一般地，如果已計算得到含根區間[ak,bk](k=1,2,…)，則令ffk+1=(ak+bk)/2，若

終值f*=(ak+bk)/2，設真實值為f0，對于給定精度|f*-f0|<ε，實際計算是所采用的終止原則為

bk-ak<2ε=10-p

(21)

取ε=10-p/2，其中p為正整數，則有

(22)

(23)

這里的[]為取整符號。

可以看到，這種步步迭代算法增加的計算量主要取決于每一次K*的大小，這又與初始含根區間的邊界值密切相關，而結果的精確性主要是由ε來決定的。

3 仿真算例及抗噪性能

仿真計算模型為一根長度為0.68 m，截面尺寸為0.04 m×0.008 m的兩端簡支矩形截面梁，彈性模量70 GPa，材料密度2 700 kg/m3，劃分為20個有限元單元。利用實驗測得系統固有頻率和模態阻尼比，調整建立的理論有限元模型，如表1所示。在第10自由度上加載f=10sin(4πt)N的正弦力，已知數據為第8自由度上的加速度響應。

設激勵時間為t，計算次數為i，計算結果表明Wilson-θ反分析法不是一種成熟穩定的方法，不同的時間步長對計算結果有很大的影響。從圖1(a)，(c)，(e)可以看出，無論時間步長的取何值，當計算次數達到一定數量，結果都會出現發散。而如圖1(b)，(d)，(f)所示，修正計算后，識別結果收斂到了真實力。

圖1 不同時間步長下識別的載荷

表1 模型參數

在實際工程問題中，噪聲的影響是必須要考慮的，一個真正成熟穩定的算法應該具有不錯的抗噪性能。在上述算例中，對于已知的加速度響應數據，加入5%的隨機噪聲(Δt=0.01 s)。圖2對比了二分法迭代修正前后對于噪聲的不同表現。可以看到修正后的識別力更能反映真實情況，而且發散的趨勢明顯減小。

圖2 加噪后的識別載荷(5%噪聲水平)

下面考察對于實際工程問題中更常遇到的幾種常見加載形式下的抗噪能力，圖3顯示了沖擊載荷、方波、三角波和鋸齒波的識別結果(Δt=0.001 s)。可見，由于噪聲的影響，加載函數平直段的識別效果較差(如沖擊載荷與方波載荷)，曲線段也不如沒有噪聲時光滑，但是基本能夠反映出所加載荷的情況，是可以接受的識別結果。以組合三角函數波的長時間識別為例，加入5%的隨機噪聲，圖4放大顯示了尾段2 s的識別情況，可以看到在半分鐘內，雖然些微不收斂，但該修正算法依然可以保持良好的穩定性(Δt=0.001 s)。

圖3 不同形式加載下的識別結果(5%噪聲水平)

圖4 長時間激勵識別結果(5%噪聲水平)

4 實驗驗證

實驗模型與上述仿真模型一致，主要儀器設備有：Agilent35670動態信號分析儀，HEW 系列高性能電磁激振器，ICP型加速度傳感器以及力傳感器等。

在第6個結點上加載頻率為16 Hz的正弦力，采樣率2 048 Hz。測得第13和17結點上的加速度響應，應用上述修正算法分別得到兩條識別力曲線，與實際所加載荷比較，結果如圖5所示。

圖5 實驗數據得到的識別結果

為進一步分析實驗結果數據，對于其中較為光滑的最后一個周期作分析，設每一個采樣點的誤差值為ei(i=1,…,N。其中，N為一周期內采樣點數)，這些數組成了一個隨機序列

(24)

式中fi為實際所加載荷，pi為計算識別力，fpeak為加載的峰值。該隨機序列的統計數據如表2所示。

表2 誤差統計

其中

(25)

(26)

不同于式(24)，這兩個誤差計算公式主要反映了每一個時間點的實時誤差，式(25)是常用的誤差計算公式，式(26)是一種考慮權重的誤差計算方法，在接近零點時由于fi很小，其實時誤差可能會很大，但是這個點的誤差參考價值并不大，加權計算可以有效避免這個問題。

從圖5來看，兩條識別曲線都與實際所加載荷波形基本吻合，第1測點的識別力曲線要優于第2測點，雖然有個別點誤差較大，但總體并沒有發散的趨勢。造成實驗出現誤差的原因有多種，計算所用的有限元模型與實際模型并不完全一致，實驗儀器和設備的不穩定，以及布點位置不準確等都會造成計算值的偏差和波動。從實驗可以看到，該修正算法計算結果基本反映了真實情況，誤差在可以接受的范圍內。

5 結 論

在以往一步到位直接得到載荷的算法中，難免會遇到誤差累積的問題，本文提出了一種載荷識別的新思路，即在每一個時間步長內部迭代修正。典型載荷的仿真和實驗結果表明，該修正算法可以有效的減小由于累積誤差導致的發散，得到收斂的識別結果，并且有著理想的抗噪性能，具備一定的工程實用價值。

參考文獻：

[1] 張方, 秦遠田, 鄧吉宏. 復雜動態載荷識別技術研究[J].振動工程學報, 2006, 19 (1): 81—85.Zhang Fang, Qin Yuantian, Deng Jihong. Research of identification technology of dynamic load distributed on the structure[J]. Journal of Vibration Engineering,2006,19 (1):81—85.

[2] 秦遠田. 動載荷識別應用技術研究[D]. 南京: 南京航空航天大學, 2007.Qin Yuantian. Study on dynamic load identification applications [D]. Nanjing: Nanjing University of Aeronautics and Astronautics, 2007.

[3] Matthew S Allen, Thomas G Carne． Delayed, multi-step inverse structural filter for robust force identification [J].Mechanical Systems and Signal Processing, 2008, 22:1 036—1 054.

[4] Fergyanto E Gunawan, Hiroomi Homma， Yuichi Morisawa. Impact force estimation by quadratic spline approximation [J]. Journal of Solid Mechanics and Materials Engineering, 2008, 2(8):1 092—1 103.

[5] Chen Tsung-Chien, Lee Ming-Hui. Determination of moving tank and missile impact forces on a bridge structure [J]. Defence Science Journal, 2008, 58(6): 752—761.

[6] Mao Yu-ming,Guo Xing-lin. Experiment study on dynamic force identification by a parameter estimation method [A]. Proceedings of the SEM Annual Conference [C]. Albuquerque,2009:1—4.

[7] 韓旭, 劉杰, 李偉杰, 等. 時域內多源動態載荷的一種計算反求技術[J]. 力學學報, 2009, 41(4): 595—602.Han Xu, Liu Jie, Li Weijie, et al. A computational inverse technique for reconstruction of multisource loads in time domain[J]. Chinese Journal of Theoretical and Applied Mechanics, 2009, 41(4): 595—602.

[8] 高峰, 李德武. 運用振動反分析方法計算振動荷載[J]． 工程力學, 1999, 16(4):91—96.Gao Feng, Li Dewu. Evalution of dynamic loading by using the inverse analysis method in vibration[J]. Engineering Mechanics, 1999, 16(4):91—96.

[9] 姜金輝, 陳國平, 張方. 多點平穩隨機載荷識別方法研究[J].振動工程學報, 2009, 22 (2): 162—167.Jiang Jinhui, Chen Guoping, Zhang Fang. Identification method of multi-point stationary random Load [J]. Journal of Vibration Engineering, 2009, 22 (2): 162—167.

[10] Lin Dong-Cherng. Input estimation for nonlinear systems [J]. Inverse Problems in Science and Engineering, 2010, 18(5): 673—689.

[11] Lin Dong-Cherng. Adaptive weighting input estimation for nonlinear systems[J].International Journal of Systems Science, 2012,43(1): 31—40.

[12] Ming-Hui Lee. Intelligent fuzzy weighted input estimation method for the input force on the plate structure [J]. Structural Engineering and Mechanics, 2010, 34,(1): 1—14.

[13] Lee M H, Chen T C. Intelligent fuzzy weighted input estimation method for the forces generated by an operating rotating machine [J]. Measurement, 2011, 44(5):917—926.

[14] 趙鳳遙. 水電站廠房結構及水力機械動力反分析[D]. 大連: 大連理工大學, 2006.Zhao Fengyao. The dynamic identification of hydropower house structure and hydraulic machinery[D]. Dalian: Dalian University of Technology, 2006.