陳先智, 周新元, 曾耀祥, 張亞輝

(1. 大連理工大學(xué) 工業(yè)裝備結(jié)構(gòu)分析國家重點實驗室, 遼寧 大連 116024;2. 北京宇航系統(tǒng)工程研究所, 北京 100076)

0 引 言

隨著航空航天技術(shù)的快速發(fā)展,對運載火箭的優(yōu)化設(shè)計也提出了更高的要求,準確描述其自由運行過程中的外激勵,是結(jié)構(gòu)優(yōu)化設(shè)計的一個重要方面．火箭在運行過程中,外激勵通過傳感器直接測量往往比較困難,在此背景下,發(fā)展通過響應(yīng)信息間接獲取外部激勵的荷載重構(gòu)方法有著重要意義．

隨著荷載重構(gòu)方法的發(fā)展,其不適定性引起了許多學(xué)者的關(guān)注,而正則化理論被廣泛應(yīng)用于解決各種不適定性問題[1-3]．彭凡等[4]使用Green函數(shù)構(gòu)建荷載和響應(yīng)之間的關(guān)系,引入奇異值分解(SVD)改善矩陣病態(tài),進而重構(gòu)了自由結(jié)構(gòu)的沖擊荷載．周玙等[5]綜合利用結(jié)構(gòu)的位移模態(tài)和應(yīng)變模態(tài)參數(shù)實現(xiàn)了應(yīng)變到位移的轉(zhuǎn)換,避免了空間積分,反求出位移響應(yīng)后,結(jié)合擴展Tikhonov正則化實現(xiàn)了基于應(yīng)變響應(yīng)的自由體荷載識別．Liu等[6]組合Tikhonov正則化和截斷廣義奇異值分解(TGSVD),綜合兩者的優(yōu)點,成功識別了動力荷載．Qiao等[7]采用三次B樣條函數(shù)結(jié)合截斷奇異值分解法(TSVD),解決了沖擊荷載識別問題．張超東和黎劍安[8]利用Kalman濾波算法獲得了增秩狀態(tài)向量的最小方差無偏估計,實現(xiàn)了狀態(tài)和荷載向量的同時識別．

由于l2正則化存在過擬合的現(xiàn)象[9],因此近幾年有學(xué)者通過稀疏正則化改善這一現(xiàn)象．稀疏正則化所對應(yīng)的問題為l0正則化,在數(shù)學(xué)上是難以求解的,一般凸釋放為l1正則化求解．Qiao等[10]使用字典建立識別方程,并采用可分離近似法(SpaRSA)重構(gòu)荷載．由于字典質(zhì)量對識別結(jié)果有影響,Zhang等[11]使用雙稀疏字典學(xué)習(xí)改善字典的質(zhì)量,進而改善了識別結(jié)果．Liu和Li[12]使用盲源分離(BSS)算法結(jié)合正交基追蹤算法(OMP),等效識別了時空耦合的分布荷載．Samagassi等[13]使用縱向應(yīng)變響應(yīng)在低噪聲情況下結(jié)合小波相關(guān)矢量機(WRVM)稀疏重構(gòu)了沖擊荷載．Yan和Sun[14]通過脈沖響應(yīng)函數(shù)建立荷載與響應(yīng)的關(guān)系,利用沖擊荷載的稀疏特性,使用非負Bayes學(xué)習(xí)算法識別了沖擊荷載．

針對大維度荷載重構(gòu)問題時,其計算成本會變得非常大,因此一般常用的手段為分段擬合．對于線性時不變系統(tǒng)而言,結(jié)構(gòu)響應(yīng)由初始條件引起的自由振動和荷載引起的強迫振動線性疊加而成,而各分段的初始條件難以通過測量得到,因此分段擬合最重要的是如何合理表達未知初始條件．Prawin和Rama Mohan Rao[15]表達初始條件的連續(xù)遞推公式受噪聲影響較大,且有誤差累積和不適定現(xiàn)象．張開華[16]假設(shè)初始條件為上一段重構(gòu)荷載作用下結(jié)構(gòu)末狀態(tài)的響應(yīng),但是僅選取末狀態(tài)響應(yīng)作為初始條件,會導(dǎo)致重構(gòu)精度嚴重依賴上一步重構(gòu)結(jié)果．考慮到對于一個多自由度系統(tǒng)而言,在外力作用時,通常只會激起一些較低的振型,大部分高階振型被激起的分量很小,一般可忽略不計．Pan等[17]用低階振型表示初始條件,使用快速軟閾值迭代算法(FIST)結(jié)合Bayes信息準則(BIC)成功識別了未知初始條件下的荷載,但是識別效率較低．在此基礎(chǔ)上,Pan等[18]采用矩陣正則化構(gòu)建長時間荷載改善了識別效率,為進一步提高效率,利用各分段可以單獨重構(gòu)的特點,使用并行算法同時重構(gòu)所有分段荷載[19],但是忽略了各分段之間的聯(lián)系,不能充分利用上一段的重構(gòu)結(jié)果．為了能高效處理含有未知初始條件與大維度荷載重構(gòu)的問題,在現(xiàn)有文獻基礎(chǔ)上,本文提出利用線性系統(tǒng)的疊加原理,將荷載使用三角基函數(shù)擬合,分別求解了所有荷載基函數(shù)作用下的強迫振動響應(yīng)和低階振型,作為初始位移、初始速度的自由振動響應(yīng),直接使用稀疏Bayes學(xué)習(xí)(SBL)算法的快速算法[20-21]對未知分解系數(shù)進行求解．提出了改進分段擬合的思想,將上一分段求解得到的荷載作用在結(jié)構(gòu)上,產(chǎn)生的末狀態(tài)響應(yīng)作為可能初始條件,并使用低階振型做補充,通過不斷更新的系數(shù)矩陣,高效完成了大維度荷載重構(gòu)問題．

1 重構(gòu)方程的建立

對于自由度為n的結(jié)構(gòu)體而言,其運動微分方程為

(1)

假設(shè)阻尼為比例阻尼,令X=Φu,將其代入式(1),并前乘ΦT,根據(jù)振型的正交性,假設(shè)模態(tài)截斷數(shù)為q,式(1)可解耦為

(2)

對于線性系統(tǒng)而言,含初始條件的結(jié)構(gòu)響應(yīng)可以分解為無初始條件強迫振動響應(yīng)和有初始條件自由振動響應(yīng)．因此可以先建立無初始條件的響應(yīng)-荷載關(guān)系,由Duhamel積分可知,式(2)的模態(tài)加速度響應(yīng)[22]為

(3)

(4)

(5)

假設(shè)結(jié)構(gòu)上作用有Nf個荷載,測點有Nm個．由于是線性系統(tǒng),那么Nf個荷載作用下的結(jié)構(gòu)響應(yīng)可以分解為每個荷載單獨作用下的響應(yīng)之和．令fk為第k個荷載(k=1,2,…,Nf)作用在自由度lk的荷載,則其相應(yīng)的模態(tài)荷載pik可表示為

pik=φilkfk,

(6)

式中φilk表示第i階模態(tài)中自由度為lk處的值,fk=[fk(t1)fk(t2) …fk(tNt)]T．

將式(6)代入式(4)得到第i階模態(tài)加速度為

(7)

式中Hilk=φilkhi．當不考慮噪聲影響的情況下, 由振型分解法可得當fk單獨作用時, 第j個測點的加速度響應(yīng)為

(8)

在荷載重構(gòu)領(lǐng)域中,由于通過函數(shù)擬合可以改善求解時矩陣的性態(tài),故一般將荷載函數(shù)通過一系列的正交基函數(shù)近似展開,因此fk(t)可表達為

(9)

(10)

式中odd表示奇數(shù),even表示偶數(shù),Δω表示頻率間隔,nb和Δω可根據(jù)響應(yīng)的頻譜和所關(guān)心頻率綜合估算．對式(9)在時域離散,并代入式(8),得

(11)

稱ajk=GmjlkBk為元胞,表示當全體荷載基函數(shù)分別作用在加載點k時,測點j響應(yīng)的集合,式中yk=[y1ky2k…ynbk]T,

(12)

根據(jù)疊加原理,測點響應(yīng)與未知系數(shù)的關(guān)系為

(13)

(14)

測量響應(yīng)一般含有噪聲,因此式(13)可以改寫為

(15)

至此建立了已知測點響應(yīng)和未知系數(shù)之間的關(guān)系．由于未知噪聲的存在,直接通過對系數(shù)矩陣A求逆,解未知系數(shù)y,往往會導(dǎo)致較大的誤差,因此本文使用SBL算法求解未知系數(shù)y,并用式(9)對荷載進行重構(gòu)．

2 稀疏Bayes學(xué)習(xí)算法推導(dǎo)

在Bayes模型中,誤差被概率性的建模為零均值Gauss分布,方差為σ2[20]．本文假設(shè)噪聲項服從N(0,σ2I),I表示單位陣,方差σ2從數(shù)據(jù)中估計得出,則的似然函數(shù)為

(16)

式中N=Nt×Nm．

為了得到稀疏解,假設(shè)存在超參γi,y中元素yi的均值為0,方差為1/γi．那么y的先驗表示為

(17)

式中M=Nf×nb,γ=diag(γ1,γ2,…,γM)．

(18)

由Bayes公式,可以得到y(tǒng)的后驗:

(19)

式中

Σy=(σ-2ATA+γ)-1,

(20)

(21)

其中Σy和μy分別表示y的方差和均值．計算出y的均值μy后,此時μy即為所求基底系數(shù),結(jié)合式(9)即可重構(gòu)荷載．為快速確定未知的參數(shù)σ2和超參γ,使用文獻[21]中的快速算法．

3 未知初始條件的處理

采用較小的時間間隔能提高荷載重構(gòu)的精度,但這也會導(dǎo)致更大的預(yù)處理成本．為了提高預(yù)處理的效率往往會將測量響應(yīng)分成多個沿時間軸移動的分段,荷載按順序在每個分段單元重構(gòu),分段示意圖如圖1所示．

圖1 分段重疊擬合Fig. 1 Piecewise overlap fitting

為合理表達每個分段的初始條件,本文結(jié)合SBL算法的稀疏特性,提出改進分段擬合思想,其選用末狀態(tài)響應(yīng)作為可能初始條件,由于可能的初始條件與真實的初始條件有差異,因此引入低階振型作為初始條件的補充．

(22)

式中djs表示初始位移向量僅為Φs時,引起的第j個測點自由振動的加速度響應(yīng);vjs表示初始速度向量僅為Φs時,引起的第j個測點自由振動的加速度響應(yīng);Q為選取參與組合初始條件模態(tài)的數(shù)量,本文選取前50階彈性模態(tài)．

將初始條件使用上一分段的末狀態(tài)響應(yīng)和補充初始條件表示,則第j個測點的加速度響應(yīng)可以近似表示為

(23)

(24)

圖2 所提方法荷載重構(gòu)流程圖Fig. 2 The flowchart of the proposed method for load reconstruction

4 無約束結(jié)構(gòu)的處理

火箭在自由飛行的過程中是一種無約束結(jié)構(gòu),這種結(jié)構(gòu)存在剛體運動,而傳感器測量的響應(yīng)為結(jié)構(gòu)的彈性響應(yīng),因此需要剔除結(jié)構(gòu)運動中的剛體運動成分．于是將結(jié)構(gòu)響應(yīng)分解為剛體響應(yīng)和彈性響應(yīng):

(25a)

(25b)

(25c)

(26)

式中

(27)

5 數(shù)值算例驗證

本文采用文獻[23]中的模擬無約束火箭模型,模型的直徑為3 m,壁厚為0.005 m,長50 m,101個節(jié)點分成100個單元,相鄰節(jié)點間隔為0.5 m,彈性模量E=7×1010N/m2,剪切模量G=2.65×1010N/m2,各節(jié)點加集中質(zhì)量1 600 kg,繞縱軸的轉(zhuǎn)動慣量為100 kg·m2,繞橫軸的轉(zhuǎn)動慣量為50 kg·m2,基頻為1.05 Hz,阻尼比為0.01,采樣頻率為1 000 Hz．響應(yīng)計算使用前50階彈性模態(tài)．本文選取401個基函數(shù),采用等間隔頻率Δω=0.25 Hz．

采用整體相對誤差(R)定量地衡量荷載重構(gòu)精度,其表達式如下:

(28)

式中,Rj表示第j個重構(gòu)荷載的整體相對誤差,fTj表示第j個真實荷載向量,fIj表示重構(gòu)的第j個荷載向量,‖·‖1表示向量的1-范數(shù)．

為了模擬真實情況,對加速度響應(yīng)數(shù)據(jù)施加Gauss白噪聲,公式如下所示:

(29)

5.1 所提方法有效性驗證

為了驗證本文所提方法的有效性, 假設(shè)有3個荷載作用在結(jié)構(gòu)上, 荷載作用在距離頂端0.5 m,1.5 m,2.5 m處,荷載作用方向為橫向,荷載的時間歷程函數(shù)為

(30)

模型的初始位移和初始速度如圖3、圖4所示．

圖3 初始位移Fig. 3 Initial displacements

6個加速度傳感器安裝在距離頂端1 m,2 m,3 m,4 m,5 m,5.5 m處,使用5%,10%和15%三個噪聲等級驗證方法的抗噪性．

在不同噪聲等級下,使用本文所提的方法三點加載荷載重構(gòu)結(jié)果如圖5所示．由圖中可以看出,本文所提方法重構(gòu)的荷載都基本和原荷載重合,驗證了本文方法對于不同噪聲水平的魯棒性和可行性．

為了說明所提方法的優(yōu)勢,將本文所提方法與文獻[24]中加權(quán)l(xiāng)1正則化方法(weightedl1-norm regularization method)在計算精度和效率上進行比較．由于兩種方法采用相同的基函數(shù),所以預(yù)處理即計算荷載基函數(shù)結(jié)構(gòu)響應(yīng)的時間成本相同,因此在兩種方法的計算成本對比上,僅列出計算式(26)所需的時間．

圖4 初始速度Fig. 4 Initial velocities

圖5 本文所提方法不同噪聲水平下荷載重構(gòu)結(jié)果Fig. 5 Reconstruction results with the proposed method under different noise levels

由表1可知,本文所提方法與加權(quán)l(xiāng)1正則化方法在重構(gòu)精度上相差無幾,但是計算的平均時長僅為0.68 s,而加權(quán)l(xiāng)1正則化方法的平均時間則需要617.02 s．以15%的噪聲為例,荷載重構(gòu)的結(jié)果如圖6所示,從圖中可以看出,本文所提方法和加權(quán)l(xiāng)1正則化方法均能比較準確地重構(gòu)出真實荷載,僅在幅值處有略微的差別．

由上可以得出,本文所提方法和加權(quán)l(xiāng)1正則化方法都能夠處理未知初始條件下的荷載重構(gòu)問題,且在幾乎相同求解精度的情況下,本文所提方法通常能夠更快地獲取重構(gòu)結(jié)果．

5.2 改進分段擬合驗證

為驗證本文選取初始條件表達方式的合理性,式(24)中的初始條件項僅由低階振型表達[18]和本文所提改進分段擬合做比較．至于僅考慮將上一分段末狀態(tài)作為初始條件的情況[16],由于受上一步計算結(jié)果影響,結(jié)合本文方法后,重構(gòu)效果非常不理想,因此本文便不再與其進行比較．

假設(shè)荷載加載點距離頂端0 m,2 m,4 m,加速度傳感器距離頂端0.5 m,1.5 m,2.5 m,3.5 m,4.5 m和5.5 m,加速度響應(yīng)數(shù)據(jù)施加15%的噪聲,測量的加速度響應(yīng)時間為4 s．對響應(yīng)曲線進行分段,本文所提改進分段擬合思想重疊長度為0,每段時間長為2 s．根據(jù)測點和加載點位置調(diào)用相應(yīng)的元胞,組成系數(shù)矩陣．荷載函數(shù)表達式為

(31)

僅使用低階振型表達初始條件的重構(gòu)結(jié)果和本文所提改進分段擬合的重構(gòu)結(jié)果如圖7、圖8所示．為了清楚對比兩種不同初始條件表達的荷載重構(gòu)效果,將圖7中重構(gòu)荷載f2第2秒的銜接點放大,如圖9(a)所示,同樣將圖8中荷載f2第2秒的銜接點放大,如圖9(b)所示,并結(jié)合表2可以明顯看出,僅使用低階振型表達初始條件的重構(gòu)結(jié)果并不理想,而本文所提改進分段擬合的重構(gòu)結(jié)果基本與真實荷載重合,并且兩者在計算效率上也幾乎相同．由此可以得出,本文所提改進的分段擬合是合理、有效的．

表2 3種重構(gòu)方式識別結(jié)果對比

對比不分段和分段的整體相對誤差,分段之后精度比不分段的精度略好,這是因為選擇了合適的分段長度并且荷載比較復(fù)雜．采用分段擬合重構(gòu)荷載所需時間約占不分段重構(gòu)荷載的1/2,這是由于不分段意味著預(yù)處理所需計算響應(yīng)的時間維度更大,相應(yīng)的時間成本也更大．由此,本文所提的改進分段擬合方法,在采用合適分段長度的情況下,能夠高效準確重構(gòu)大維度的復(fù)雜荷載．

由上可知,本文提出分段擬合的思想能明顯提升計算效率,對于復(fù)雜的荷載分段擬合也能夠提高反演精度．針對大維度荷載重構(gòu)問題,本文認為分段長度不宜太長,太長會導(dǎo)致預(yù)處理的成本過大,當荷載比較復(fù)雜時還會導(dǎo)致重構(gòu)精度降低．

圖8 改進分段擬合Fig. 8 Modified piecewise fitting time histories

(a) 僅使用低階振型表達初始條件(a) Only low-order modes expressed initial conditions

(b) 改進分段擬合(b) Modified piecewise fitting time histories圖9 銜接點放大圖Fig. 9 The articulation point enlargement

6 結(jié) 論

本文針對含有未知初始條件無約束結(jié)構(gòu)的荷載重構(gòu)問題,提出基于SBL算法荷載識別方法結(jié)合改進分段擬合提高了整體的計算效率,并在識別過程中利用了上一分段的識別結(jié)果,通過合理的初始條件表達方式,改善了識別的精度．數(shù)值結(jié)果也表明,本文方法較現(xiàn)在流行方法能夠更快速地重構(gòu)荷載,且通過不同等級噪聲下的重構(gòu)結(jié)果,驗證了本文方法具有較強的魯棒性;和僅使用低階振型表達初始條件相比,采用的改進分段擬合更加合理;改進的分段擬合手段在采用合適的分段長度的情況下,和不分段的識別結(jié)果相比,也能夠以較高的精度重構(gòu)荷載．