基于頻響函數(shù)提升小波總能量的桁架結構隨機模型修正

趙永鵬， 殷 紅， 彭珍瑞

(蘭州交通大學 機電工程學院，蘭州 730070)

1 引 言

隨著工程技術領域的不斷發(fā)展，有限元技術的應用更加廣泛，但是各種不確定性因素，如建模時的簡化及安裝誤差等對結果的影響很大，因此建立一個可靠有效的有限元模型就變得尤為重要。當前，通過模型修正來獲得準確的有限元模型(FEM)是主流的方法。模型修正方法分為確定性和不確定性兩種[1]，確定性方法在近年來有了廣泛發(fā)展。但是，由于各種不確定性的情況，出現(xiàn)的誤差較多。所以，在模型修正中充分考慮不確定性的影響，從而構建一種隨機性的模型修正方法勢在必行。

參數(shù)不確定性的模型修正方法指通過多次模態(tài)實驗得到相應數(shù)據(jù)，計算出具有統(tǒng)計特征的有限元模型，來表明實際結構的靜動力學行為、機理和演化規(guī)律等。方圣恩等[2]充分考慮了修正參數(shù)的不確定性，通過計算參數(shù)的統(tǒng)計特征量對一組鋼板進行了修正。Fang等[3]通過Hermite多項式的顯式多項式混沌展開，把不確定的參數(shù)和響應通過構造一個隨機響應面模型(RSM)對一組金屬板進行了修正。Zhai等[4]以模態(tài)參數(shù)為響應，基于改進的響應面模型和先進的蒙特卡洛(MC)方法，對一航空發(fā)動機定子系統(tǒng)(殼體)進行了修正。Bi等[5]構建了一個完全嵌入Bhattacharyya距離的近似貝葉斯模型的修正框架，對一質量彈簧系統(tǒng)進行了模型修正，證明了方法的有效性。Deng等[6]通過構建代理模型、設置校準參數(shù)以及將距離函數(shù)作為參數(shù)相關性的衡量標準等步驟進行了隨機模型修正。Khodaparast等[7]利用模態(tài)響應(如固有頻率和模態(tài)振型)中得到的可變性來估計隨機修正參數(shù)的第一和第二統(tǒng)計矩，通過一個三自由度的質量彈簧系統(tǒng)對所提方法進行了驗證。

通過結構模態(tài)參數(shù)進行修正將不可避免地在識別的過程中引入一些誤差，實際工程中由于實驗模態(tài)參數(shù)的識別而出現(xiàn)的誤差有可能會在某些特定情況下比建模精度低引起的誤差更大。而基于頻響函數(shù)(FRFs)的模型修正方法因為可以略去結構模態(tài)參數(shù)識別而得到廣泛應用，該方法減少了各種誤差的出現(xiàn)，適于對模態(tài)分布相對比較密的結構進行修正。李偉明等[8]利用頻響函數(shù)為輸出對一個二維桁架進行了修正。王巨濤等[9]利用加速度頻響函數(shù)曲線對應頻率點處的響應值之差構造目標函數(shù),對GARTEUR飛機模型進行了修正。但是頻響函數(shù)包含的數(shù)據(jù)信息量很大，怎樣選出好的頻率點數(shù)據(jù)以得到期望的修正效果，現(xiàn)在還沒有形成統(tǒng)一的標準。因此，可以考慮引入一個特征向量對頻響函數(shù)進行表征，以避免上述問題。張勇等[10]提取頻響函數(shù)的特征向量對一扭力梁模型進行了修正。上述研究大都應用于確定性的模型修正，但鮮見應用于不確定性模型修正。提升小波變換LWT(Lifting Wavelet Transform)屬于第二代小波變換，與傳統(tǒng)小波變換相比，提升小波變換不僅擁有多分辨率的特點，而且具有計算時間少，整個流程相對簡單以及變換后數(shù)據(jù)長度不變等優(yōu)勢[11]。宋萌萌等[12]選擇提升小波變換對得到的初始信號進行降低噪聲的處理，得到的結果滿足要求。袁海英等[13]利用提升小波變換對振動信號進行時頻特性分析和信息預處理，通過預測器和更新器的設計取代了小波基函數(shù)選取過程。提升小波變換不再使用傅立葉變換，也沒有利用函數(shù)的平移或伸縮，而是在時域中就可以完成原始信號不同頻帶上的分離。經過變換后各層分量的小波能量對信號的變化有很高的敏感性，且對噪聲污染、小波種類選擇及小波層數(shù)具有魯棒性。因此，選擇提升小波變換對頻響函數(shù)進行分解，提取總能量值來表征頻響函數(shù)可以最大程度上保證頻響函數(shù)的特征不發(fā)生變化，而且也減少了頻響函數(shù)信息過多引起的各種處理困難。

本文構建了一種以提升小波總能量為響應，通過篩選響應樣本反求出待修正參數(shù)的隨機模型修正方法。首先對加速度頻響函數(shù)進行提升小波變換，提取對應的提升小波總能量作為輸出，以待修正參數(shù)作為輸入，構造一個多項式響應面代理模型代替有限元模型；然后利用蒙特卡洛抽樣得到大量的響應樣本，并通過設置閾值對響應樣本進行優(yōu)化篩選；最后構造一個優(yōu)化反問題，利用布谷鳥優(yōu)化算法反求出輸入即待修正參數(shù)，計算修正后參數(shù)的統(tǒng)計特征。分別對平面桁架和三維鋼桁架結構模型進行修正，證明了所構建方法的有效性。

2 頻響函數(shù)提升小波總能量

2.1 加速度頻響函數(shù)

在具有阻尼的系統(tǒng)中，通過一個簡諧激勵作用所得的動力學方程為

(1)

在頻域內的輸入輸出關系可以表示為

A(ω)=H(ω)F(ω)

(2)

式中A(ω)為穩(wěn)態(tài)響應，H(ω)為加速度頻響函數(shù)，F(xiàn)(ω)為簡諧激勵，ω為激勵頻率。

相應的加速度頻響函數(shù)可以表示為

(3)

2.2 提升小波變換

對于得到的頻響函數(shù)，選擇提升小波變換對其能量值進行提取。因為常見的小波濾波器的輸出是浮點數(shù)，所以對變換后的值選擇壓縮時需要量化才能得出對應的整數(shù)。Swelden[14]選擇了根據(jù)提升方式改進小波特性得到小波基函數(shù)的新模式即提升小波變換。提升小波變換共分為三個過程: 分裂、預測和更新。

(1) 分裂。根據(jù)原始信號sj奇偶性的不同可以得到兩個沒有相關性的子集。所得子集的長度是初始信號的1/2。如式(4)所示：

sj=(ej - 1,oj - 1)

(4)

式中ej - 1={ej - 1,k=sj,2k}為偶序列，oj - 1={oj - 1,k=sj,2k + 1}為奇序列。

(2) 預測。根據(jù)第一步提取的偶序列ej - 1來預測奇序列oj - 1。如式(5)所示，用實際值oj - 1與預測值ej - 1的差值dj - 1來表示其接近程度，稱作細節(jié)系數(shù)，同時也屬于原始信號的高頻分量。

dj - 1=oj - 1-P(ej - 1)

(5)

式中P為預測函數(shù)，用相鄰數(shù)據(jù)的平均值來表示：

pk(ej - 1)=(ej - 1,k+ej - 1,k + 1)/2

(6)

(3) 更新。這一步的目的是要使第二步產生的某些特征(如均值等)與原始信號的特征保持不變。可以通過算子U來實現(xiàn)：

sj - 1=ej - 1+U(dj - 1)

(7)

式中sj - 1為原始信號的低頻部分，U為更新算子函數(shù)，可以表示為Uk(dj - 1)=dj - 1,k/2。

通過提升小波變換之后，可以將原始信號sj變?yōu)楦哳l分量dj - 1和低頻分量sj - 1；將低頻部分sj - 1繼續(xù)進行相同的步驟，又可以分解為二層高頻dj - 2和二層低頻sj - 2；以此類推，可以根據(jù)要求進行n層的提升小波分解，原始信號就可以表示為{sj - n,dj - n,…,dj - 1}。其中sj - n為信號的低頻部分，{dj - n,…,dj - 1}為信號的高頻部分。

對頻響函數(shù)進行提升小波變換后，即可得出各層的高頻部分和低頻部分。設Ej d為原始信號sj在尺度n上的高頻分量總能量：

(8)

相應的低頻分量能量Ej s也可以得到，信號的提升小波總能量為各個尺度下的小波能量的總和[15]。

Ez=E1d+E2d+…+Ej d+Ej s

(9)

3 響應面模型

對于一些大型復雜結構，由于其內部存在很多的零部件，以及很多邊界尺寸數(shù)據(jù)的不確定性，所以在有限元建模的過程中存在困難，從而建立的有限元模型變得十分復雜，在模型修正過程中可能會出現(xiàn)計算無法收斂或者所用時間很長等問題。選擇一個代理模型來代替這個復雜的有限元模型，可以方便計算，提高效率。多項式響應面代理模型具有建立相對簡單且預測精度高的優(yōu)點，因此本文選擇建立響應面模型來代替有限元模型。這個代理模型將一個物理系統(tǒng)的輸入?yún)?shù)x和輸出響應y通過數(shù)學多項式結合起來：

y=f(x1,x2,…，xn)+ε

(10)

式中f為輸入與輸出間的映射，ε為建模誤差，n為輸入?yún)?shù)的數(shù)量。

對于一般性的工程問題，大多采用二階多項式模型作為一種基本的形式，即

(11)

式中bo,bi,bi i和bi j為回歸系數(shù)，xi和xj為輸入?yún)?shù)。

響應面構建完成后，需對其進行精度檢驗，判斷是否滿足后續(xù)計算要求，常用的檢驗方法為計算真實值與預測值間的均方根誤差，即

(12)

因此,根據(jù)提取的輸出與輸入?yún)?shù)就可以構建二階響應面代理模型，然后再對其進行精度檢驗。

4 響應樣本獲取

4.1 蒙特卡洛抽樣

蒙特卡洛方法提供了一種可替代分析數(shù)學評估統(tǒng)計量在隨機樣本中行為的方法。根據(jù)樣本的概率分布，用計算機實現(xiàn)統(tǒng)計模擬或抽樣，屬于一種隨機抽樣方法。常見的數(shù)值計算方法由于需要多次迭代，導致最后誤差累積很大，影響結果，而蒙特卡洛抽樣雖然也有迭代，但是沒有像其他方法那樣頻繁迭代，產生的誤差一般可以忽略不計，對最后的抽樣結果也是基本沒有影響。因此，在已知響應概率分布的情況下，運用蒙特卡洛抽樣可以得到大量的響應樣本，認為這些樣本就是試驗值，然后用于隨機模型修正過程。

4.2 響應樣本篩選

通過蒙特卡洛抽樣得到大量響應樣本，由于是一種隨機的抽樣方法，所以可能會存在抽取的樣本波動范圍過大，以至于超出了實際試驗所能測得的范圍，樣本稱之為失真樣本。如果將這些數(shù)據(jù)也代入修正過程，就失去了實際的意義，可能導致修正的結果很不理想。因此，可以通過設定一個閾值來對抽樣數(shù)據(jù)進行一次優(yōu)化篩選，使得篩選后的數(shù)據(jù)與實際情況更加接近。本文通過計算抽樣響應值與真實響應值間的均方根誤差進行樣本的篩選，

(13)

為了保證所抽取樣本的真實性并且選擇有效的樣本，按照至少保留80%左右樣本的準則進行篩選。計算所有的RMSE值，根據(jù)該準則設定一個閾值，當RMSE值比閾值小時保留，反之則淘汰。

5 模型修正

根據(jù)上述理論，通過響應面模型代替有限元模型后，以蒙特卡洛抽樣后優(yōu)選的響應值與響應面模型預測響應值差值最小來構造目標函數(shù)。

(14)

建立目標函數(shù)之后，模型修正變成對目標函數(shù)求解的問題。由于布谷鳥優(yōu)化算法具有全局搜索能力強和參數(shù)少等優(yōu)勢，因此本文選擇該算法進行尋優(yōu)以得到修正后的參數(shù)值。

綜合頻響函數(shù)提升小波總能量的提取，響應面模型的建立，響應樣本的篩選，目標函數(shù)的求解等步驟，可得隨機模型修正流程如圖1所示。

圖1 模型修正流程

6 數(shù)值算例

6.1 平面桁架

如圖2所示平面桁架，該桁架結構共有14個節(jié)點和25個自由度。桿單元的彈性模量為210 GPa，材料密度為7850 kg/m3。將該平面桁架結構作為試驗模型，選擇彈性模量和密度為待修正參數(shù)。對參數(shù)進行偏離，其中密度減小10%，彈性模量增加10%作為相應的有限元模型列入表1。選擇第11節(jié)點y方向為激勵點位置，第3節(jié)點y方向為測點位置。

圖2 平面桁架模型

表1 試驗模型與有限元模型參數(shù)

選擇拉丁超立方法獲取樣本點。在待修正參數(shù)的±20%區(qū)間內抽取40個樣本點。對于抽樣所得40組樣本，選擇前20組為樣本集擬合響應面模型，后20組為測試集，用來驗證代理模型的精度。通過回歸擬合出響應面模型如圖3所示。

圖3 平面桁架的響應面模型

假設試驗值均服從正態(tài)分布。依據(jù)抽樣所得樣本點，可以得到其加速度頻響函數(shù)值，對其進行三層的提升小波變換，提取相應的提升小波總能量。根據(jù)蒙特卡洛方法完成抽樣，共抽取200個樣本。對于抽取的樣本通過設定閾值進行篩選，按照至少保留80%數(shù)據(jù)篩選過后，代入修正過程進行求解。

對構建的響應面模型進行精度檢驗后，所得均方根誤差RMSE值如圖4所示?？梢钥闯?，RMSE值最大不超過0.06。由此可得，所構建的響應面模型精度滿足計算要求，該模型有效。

圖4 響應面模型RMSE值

通過上述方法進行模型修正，修正結果列入 表2。由表2可知，修正誤差小于2.1%，達到了較高的修正精度。

為了進一步說明方法的有效性，修正前后的頻響函數(shù)曲線如圖5所示。可以看出，修正后的頻響函數(shù)曲線與真實值的曲線重合度高，證明用提升小波總能量來表征頻響函數(shù)作為結構響應進行模型修正是可行的。

圖5 修正前后頻響函數(shù)曲線

6.2 三維桁架

如圖6所示三維桁架。該桁架共有28個節(jié)點，66個單元，48個自由度。總長2.8 m，寬0.39 m，高0.27 m。桁架共四個支座固定(節(jié)點編號為1,8,9,16)，節(jié)點鉸接。桁架的材料參數(shù)彈性模量為190 GPa，密度為7800 kg/m3。尺寸參數(shù)上弦桿和下弦桿橫截面積為85.5 mm2，直腹桿橫截面積為141.0 mm2，斜腹桿橫截面積為45.0 mm2。結構的激勵位置和響應位置分別如圖6所示。

圖6 三維桁架模型

6.2.1 試驗設計與樣本選取

選擇如圖6所示三維桁架結構作為本算例的試驗模型，相應的有限元模型通過對待修正參數(shù)進行一定范圍的偏離得到。選擇材料參數(shù)的彈性模量和密度以及尺寸參數(shù)的上弦桿橫截面積共三個參數(shù)為待修正參數(shù)。對參數(shù)進行偏離，其中密度減小10%，彈性模量增加10%，上弦桿橫截面積減少10%，有限元模型參數(shù)列入表3。

表3 試驗模型與有限元模型參數(shù)

建立響應面模型，樣本點的選擇很重要。常用的試驗設計方法有中心復合法、全因子法和拉丁超立方法等。本文選擇拉丁超立方法獲取樣本點，該方法采樣比較均勻且樣本具有代表性，符合模型建立的要求。在三個待修正參數(shù)的±20%區(qū)間內利用拉丁超立方抽樣抽取40個樣本點。

假設試驗值均服從正態(tài)分布。根據(jù)蒙特卡洛方法完成抽樣，共抽取2000個樣本。對于抽取的樣本通過設定閾值進行篩選，按照至少保留80%數(shù)據(jù)(即1600個)的要求篩選后，代入修正過程進行求解。

6.2.2 提升小波總能量的提取

依據(jù)拉丁超立方抽樣所得樣本點，可以得到其加速度頻響函數(shù)值，如圖7所示。以試驗值的頻響函數(shù)曲線為例，對其進行三層的提升小波變換，可以求得各層高頻分量和低頻分量如圖8所示，然后提取相應的提升小波能量，將提取的各層能量求和，就可以得到提升小波總能量。

圖7 試驗模型頻響函數(shù)曲線

圖8 三層提升小波變換分量

6.2.3 響應面代理模型的構造及檢驗

對于抽樣所得40組樣本，選擇前20組為樣本集作為響應面模型的輸入x，后20組為測試集，用來驗證代理模型的精度。將20組輸入與其所對應的輸出(即頻響函數(shù)提升小波總能量)通過回歸擬合求出響應面模型的各個參數(shù)，得到代理模型表達式為

y=18.2007-0.0011x1-0.0091x2-215.08x3+

0.0013x1x3+0.043x2x3

(15)

將后20組樣本代入計算得到相應的輸出，然后對構建的響應面模型進行精度檢驗。所得均方根誤差RMSE值如圖9所示。均方根誤差越接近0表明響應面精度越高。可以看出，RMSE值不超過0.07，所構建的響應面模型精度滿足計算要求，該模型有效。

圖9 響應面模型RMSE值

6.2.4 模型修正

通過蒙特卡洛抽樣得到2000個響應樣本如 圖10 所示，通過多次預實驗發(fā)現(xiàn)設定閾值為0.092時滿足篩選準則要求。對所有樣本排序后進行篩選，經過篩選后共有1643個樣本符合要求。如圖11所示為篩選后的樣本。將這些樣本依次代入目標函數(shù)中，通過布谷鳥算法對目標函數(shù)求得最小值，可得到對應的1643組輸入即為修正后的參數(shù)。算出修正后各參數(shù)的均值，修正前后的參數(shù)均值列入表4。

圖10 蒙特卡洛抽取的響應值

圖11 優(yōu)化篩選的響應值

表4 修正前后參數(shù)均值

由表4可知，修正后參數(shù)均值的誤差相比于修正前有了明顯的減小，說明修正效果良好。修正前后的提升小波總能量值列入表5。由表5可知，提升小波總能量的誤差比修正前也有了大幅的減小，表明該方法完成隨機模型修正是可行的。采用修正后的參數(shù)均值可得其加速度頻響函數(shù)如 圖12 所示。圖中修正后的頻響函數(shù)曲線與真實值的曲線基本重合，相較于修正前的曲線，修正效果達到了預期的目的。

表5 修正前后提升小波總能量

圖12 修正前后頻響函數(shù)曲線

圖13和圖14分別給出修正前后頻響函數(shù)實部和虛部曲線，可以看出，修正值與真實值的實部虛部曲線基本重合，進一步說明修正結果良好。

圖13 頻響函數(shù)實部曲線

圖14 頻響函數(shù)虛部曲線

為了進一步驗證所提方法的有效性，以平面桁架為例，分別將頻響函數(shù)和提升小波總能量作為響應進行模型修正,結果列入表6，所有結果均是在同一臺計算機上運行20組?？梢钥闯?，將提升小波總能量作為響應進行模型修正所用時間更短；在修正精度方面，彈性模量修正精度有了很大提高，雖然材料密度修正精度略有降低，但是總體來說效率有了一定提高，計算量也有所減小，能夠達到隨機模型修正的目的。

表6 誤差與時間對比

7 結 論

(1) 本文構建了一種基于加速度頻響函數(shù)提升小波總能量的隨機模型修正方法，提取頻響函數(shù)的提升小波總能量作為響應，應用于隨機模型修正過程，既利用了頻響函數(shù)，又減小了計算量，提高了計算效率，修正的結果達到預期目標。

(2) 對于蒙特卡洛抽樣得到的響應樣本，通過設定閾值進行了一次優(yōu)化篩選，充分考慮了修正過程的隨機性，保證了樣本可靠性的同時也使后續(xù)的修正過程變得更加容易。

(3) 通過平面桁架和三維桁架模型對所提方法進行驗證，修正后參數(shù)均值的最大誤差小于 3.3%，對應點的頻響函數(shù)曲線吻合良好，達到了隨機模型修正的目的。