高溫環(huán)境下結(jié)構(gòu)動(dòng)力學(xué)模型修正方法研究

袁昭旭， 于開平

(哈爾濱工業(yè)大學(xué) 航天學(xué)院，哈爾濱 150001)

高溫環(huán)境下結(jié)構(gòu)動(dòng)力學(xué)模型修正方法研究

袁昭旭， 于開平

(哈爾濱工業(yè)大學(xué) 航天學(xué)院，哈爾濱 150001)

提出了一種新的基于貝葉斯參數(shù)估計(jì)的高溫環(huán)境動(dòng)力學(xué)模型修正方法，利用此方法對(duì)一復(fù)合材料層合板模型進(jìn)行了修正。該方法的主旨在于，利用高溫環(huán)境下的模態(tài)頻率試驗(yàn)數(shù)據(jù)修正有限元模型。該方法的特點(diǎn)在于，可以通過(guò)參數(shù)和響應(yīng)的離散度來(lái)描述其不確定性，從而更加合理地解決實(shí)際問題。在研究中，進(jìn)行了高溫環(huán)境下復(fù)合材料結(jié)構(gòu)的動(dòng)力學(xué)試驗(yàn)，通過(guò)試驗(yàn)得到其高溫環(huán)境下的模態(tài)頻率。以此為基礎(chǔ)，對(duì)復(fù)合材料有限元模型進(jìn)行修正。在修正前，對(duì)模型進(jìn)行了靈敏度分析，探討了參數(shù)與高溫環(huán)境動(dòng)力學(xué)響應(yīng)之間的關(guān)系。根據(jù)試驗(yàn)所得數(shù)據(jù)和對(duì)待修正參數(shù)的預(yù)估設(shè)置其對(duì)應(yīng)的離散度。利用自編譯程序進(jìn)行修正計(jì)算。通過(guò)修正，模型的溫度-固有頻率對(duì)應(yīng)關(guān)系得到了明顯的改善，證明了方法的有效性，可以適用于工程實(shí)際問題。

高溫環(huán)境結(jié)構(gòu)動(dòng)力學(xué)；模型修正；復(fù)合材料；靈敏度分析

近幾十年來(lái)，高超聲速飛行器的動(dòng)力學(xué)問題一直為研究者們所關(guān)注。文獻(xiàn)[1-7]高超聲速飛行所帶來(lái)的高溫環(huán)境對(duì)于結(jié)構(gòu)的動(dòng)力學(xué)特性有著較大的影響，故高溫環(huán)境結(jié)構(gòu)動(dòng)力學(xué)是一項(xiàng)重要的研究課題。本文提出了一種新的基于貝葉斯參數(shù)估計(jì)的模型修正方法，并將其應(yīng)用于復(fù)合材料層合板模型的修正工作中。

對(duì)于高溫環(huán)境下的動(dòng)力學(xué)問題，學(xué)者們進(jìn)行了很多基礎(chǔ)性的研究。Bailey[8]研究了當(dāng)存在熱應(yīng)力影響時(shí)不同邊界條件下板的振動(dòng)問題，他考慮了在非耦合情況下各向異性板的情況。Chang等[9-10]進(jìn)一步研究了耦合情況下考慮大變形非線性的熱彈性動(dòng)力學(xué)問題。Brown[11]通過(guò)實(shí)驗(yàn)和數(shù)值仿真的手段，研究了X-34 FASTRAC復(fù)合材料噴嘴的高溫動(dòng)力學(xué)特性。針對(duì)復(fù)合材料板、殼結(jié)構(gòu)的研究也有很多進(jìn)展[12-13]。

在反問題領(lǐng)域，學(xué)者們對(duì)高溫環(huán)境結(jié)構(gòu)的參數(shù)辨識(shí)也進(jìn)行了大量的研究。這些研究分為兩個(gè)主要方面，一方面是對(duì)熱物理參數(shù)的辨識(shí)[14-18]，一方面是對(duì)高溫環(huán)境下熱彈性參數(shù)的辨識(shí)[19-21]。其中，前者主要是針對(duì)研究熱傳導(dǎo)問題，并不是本文的重點(diǎn)。

熱彈性動(dòng)力學(xué)問題，其實(shí)質(zhì)是熱與彈性的耦合問題。但是在實(shí)際研究中，可以簡(jiǎn)化為熱對(duì)彈性的單向耦合，即所謂高溫結(jié)構(gòu)動(dòng)力學(xué)問題。高溫環(huán)境對(duì)于結(jié)構(gòu)動(dòng)力學(xué)特性的影響主要來(lái)自兩個(gè)方面：①溫度對(duì)材料屬性的影響，主要是材料的彈性模量受較高溫度的影響會(huì)減小，隨之改變的是結(jié)構(gòu)的整體剛度；②溫度梯度引起的結(jié)構(gòu)的熱應(yīng)力會(huì)造成一定的結(jié)構(gòu)附加剛度。結(jié)構(gòu)剛度受熱環(huán)境的影響往往會(huì)造成結(jié)構(gòu)固有頻率的下降，并在一定程度上縮小各階固有頻率之間的差距，進(jìn)而影響結(jié)構(gòu)的振動(dòng)特性[22-23]。

在實(shí)際工程中，常常使用有限元等數(shù)學(xué)模型來(lái)對(duì)系統(tǒng)的動(dòng)力學(xué)特征及相應(yīng)進(jìn)行預(yù)示。但是由于數(shù)學(xué)模型本身以及建模者的因素總會(huì)造成一定的誤差，這種誤差體現(xiàn)在數(shù)學(xué)模型與試驗(yàn)?zāi)Ｐ椭g的不吻合。為了提高有限元模型的預(yù)示精度，學(xué)者們提出了模型修正的概念。利用試驗(yàn)數(shù)據(jù)對(duì)模型進(jìn)行修改，從而使有限元模型與試驗(yàn)?zāi)Ｐ透游呛稀?0世紀(jì)80年代開始，模型修正的基礎(chǔ)理論開始逐步形成，并延伸出大量的方法。 Mottershead等[24-25]對(duì)這些方法進(jìn)行了系統(tǒng)的總結(jié)，在此不再?gòu)?fù)述。近些年來(lái)，模型修正問題的一個(gè)熱門方向是基于不確定性的模型修正方法的研究[26-29]。這一方面研究方法有很多，包括概率型的[30]、非概率型的[31-32]。其中基于貝葉斯參數(shù)估計(jì)的方法，由于可以將試驗(yàn)與仿真二者的不確定性同時(shí)考慮入模型修正流程，被廣泛使用于模型修正問題的研究[33-36]。

對(duì)于高溫環(huán)境下的模型修正方法，目前的研究還很少[41]。Humbert等[37]提出了基于熱彈性場(chǎng)的修正方法，這種方法強(qiáng)調(diào)了熱應(yīng)力對(duì)結(jié)構(gòu)動(dòng)特性的影響，但是由于熱彈性場(chǎng)本身的測(cè)試有一定難度，故工程應(yīng)用價(jià)值有限。Cheng等[38-39]提出了高溫環(huán)境分級(jí)模型修正方法，將溫度場(chǎng)參數(shù)和結(jié)構(gòu)物理參數(shù)分為兩個(gè)級(jí)別進(jìn)行修正。Sun等[40]利用模型修正的方法，對(duì)結(jié)構(gòu)熱參數(shù)進(jìn)行了辨識(shí)，在他的研究中考慮了材料熱物理參數(shù)隨溫度變化的特性。在這些研究中，都沒有將參數(shù)和試驗(yàn)響應(yīng)的不確定性考慮進(jìn)去，而本文提出的基于貝葉斯理論的修正方法彌補(bǔ)了這項(xiàng)空白。同時(shí)，現(xiàn)有的方法都是針對(duì)傳統(tǒng)各向同性材料的，本文研究了在復(fù)合材料結(jié)構(gòu)的高溫環(huán)境模型修正方法，其響應(yīng)和參數(shù)特性更加復(fù)雜。

1 高溫環(huán)境結(jié)構(gòu)動(dòng)力學(xué)基礎(chǔ)

從原理上來(lái)講，高溫環(huán)境對(duì)結(jié)構(gòu)振動(dòng)特性的影響歸根結(jié)底是對(duì)結(jié)構(gòu)整體剛度矩陣的影響。考慮結(jié)構(gòu)整體剛度受溫度影響主要有兩方面：一方面在高溫環(huán)境下，材料彈性模量E將隨溫度的改變而發(fā)生較大幅度的變化，進(jìn)而改變了結(jié)構(gòu)本身的整體剛度。在高溫環(huán)境下，結(jié)構(gòu)的剛度矩陣與常溫剛度矩陣具有相同的形式：

(1)

式中：B為幾何常數(shù)矩陣，由結(jié)構(gòu)及其約束形式確定；D是結(jié)構(gòu)的彈性常數(shù)矩陣，只與結(jié)構(gòu)材料彈性模量E以及泊松比μ有關(guān)。當(dāng)溫度變化時(shí)彈性模量隨之變化，則矩陣D也將隨之改變。

另外，結(jié)構(gòu)在高溫環(huán)境下不可避免的會(huì)在內(nèi)部產(chǎn)生溫度梯度，這將會(huì)引起結(jié)構(gòu)內(nèi)部一定的熱應(yīng)力，進(jìn)而產(chǎn)生附加的初始應(yīng)力剛度矩陣改變結(jié)構(gòu)整體剛度分布。以下簡(jiǎn)單介紹附加剛度矩陣的形成。

彈性體受溫度載荷將產(chǎn)生初應(yīng)變?chǔ)?，在彈性應(yīng)力的作用下產(chǎn)生彈性應(yīng)變D-1σ，彈性體的總應(yīng)變?chǔ)攀沁@兩者的和，即

ε=D-1σ+ε0

(2)

式中，{σ}是考慮變溫影響的彈性應(yīng)力，也就是熱應(yīng)力。

σ=D(ε-ε0)

(3)

綜上，結(jié)構(gòu)的熱剛度矩陣為：

K=KT+Kσ

(4)

該式即為結(jié)構(gòu)熱振動(dòng)問題的有限元計(jì)算中所用到得結(jié)構(gòu)受溫度及溫度梯度影響存在預(yù)應(yīng)力情況下的總剛度矩陣。

從數(shù)學(xué)方法上，傳統(tǒng)的模型修正方法多為基于靈敏度的方法。即通過(guò)微分方法或有限差分方法找到待修正參數(shù)與目標(biāo)響應(yīng)之間的微分關(guān)系，然后通過(guò)建立修正方程求解待修正參數(shù)的改變量。基于靈敏度的方法關(guān)鍵是要尋找修正參數(shù)與動(dòng)力學(xué)響應(yīng)之間的關(guān)系。即修正方程：

ΔR=SΔP

(5)

式中：S為靈敏度矩陣；ΔP為修正參數(shù)改變量向量；ΔR為有限元與試驗(yàn)響應(yīng)之間的殘差向量。對(duì)于高溫環(huán)境模型修正，修正參數(shù)可以選擇為材料剛度隨溫度改變率、材料熱膨脹系數(shù)等熱物理參數(shù)，也可以選擇為溫度場(chǎng)分布、應(yīng)變場(chǎng)等熱載荷變量。對(duì)于修正目標(biāo)響應(yīng)，可以選擇為熱模態(tài)頻率及對(duì)應(yīng)的模態(tài)振型。需要注意的是，在求解靈敏度方程過(guò)程中，剛度陣對(duì)于熱物理參數(shù)的微分關(guān)系應(yīng)該分為兩部分：

(6)

式中：γ=(ET-E0)/ΔT為材料剛度隨溫度變化率，α=(LT-L0)/ΔT為材料熱膨脹系數(shù)。這兩個(gè)參數(shù)為材料本身的物理屬性，不會(huì)隨外界條件所改變，可以作為修正參數(shù)。

2 基于貝葉斯參數(shù)估計(jì)的高溫環(huán)境模型修正方法

對(duì)于修正方程的求解，學(xué)術(shù)界進(jìn)行了大量的研究。比較適合本研究的是貝葉斯參數(shù)估計(jì)方法，現(xiàn)介紹其理論基礎(chǔ)：

對(duì)于一個(gè)模型修正問題，可以寫為如下修正方程的形式：

Re=Ra+S(Pu-Po)

(7)

或簡(jiǎn)寫為：

ΔR=SΔP

(8)

式中：Re為參考系統(tǒng)響應(yīng)向量(試驗(yàn)數(shù)據(jù))；Ra在參數(shù)向量{Po}情況下計(jì)算得到的響應(yīng)向量；Pu為待求的修正后的參數(shù)向量；S為靈敏度矩陣。

這是基于靈敏度的模型修正方法的一般公式，其中對(duì)于靈敏度矩陣的求解一般是通過(guò)待求參數(shù)和響應(yīng)之間的微分關(guān)系獲得的。當(dāng)待求參數(shù)和關(guān)心響應(yīng)之間不易獲得理論微分關(guān)系時(shí)，可以使用有限差分法獲得對(duì)應(yīng)靈敏度方程。

這個(gè)方程在一般情況下是超定的，對(duì)于傳統(tǒng)的參數(shù)估計(jì)理論一般使用廣義逆(最小二乘)方法求解。但是當(dāng)需要考慮參數(shù)的不確定時(shí)，則可以使用貝葉斯參數(shù)估計(jì)方法。

根據(jù)貝葉斯參數(shù)估計(jì)理論，需要對(duì)待估計(jì)參數(shù)和響應(yīng)同時(shí)進(jìn)行加權(quán)處理，對(duì)應(yīng)的加權(quán)誤差可以寫成如下形式：

E=ΔRTCRΔR+ΔPTCPΔP

(9)

這個(gè)誤差可以通過(guò)如下算法得到最小值：

Pu=Po+G(-ΔR)

(10)

其中增益矩陣

G=(CP+STCRS)-1STCR

(11)

這個(gè)公式在響應(yīng)個(gè)數(shù)多于參數(shù)個(gè)數(shù)時(shí)是有效的，當(dāng)參數(shù)個(gè)數(shù)多于方程個(gè)數(shù)時(shí)，可以使用如下公式：

(12)

這個(gè)方程可以適用于參數(shù)個(gè)數(shù)多于方程個(gè)數(shù)的欠定方程的求解。

在上述方程中，權(quán)系數(shù)矩陣CP和CR分別表征了參數(shù)和響應(yīng)的置信程度。在貝葉斯估計(jì)理論中稱之為離散度(Scatter Value)。一般來(lái)講，離散度可以通過(guò)統(tǒng)計(jì)量進(jìn)行定義，這里使用定義公式：

(13)

式中：μ為平均值；σ為標(biāo)準(zhǔn)差。對(duì)于離散度的選擇可以按照如下原則：對(duì)于你認(rèn)為置信程度低的參數(shù)設(shè)置較高的離散度，這樣在方程求解過(guò)程中此參數(shù)可以發(fā)生較大改變；對(duì)于你認(rèn)為置信程度高的參數(shù)設(shè)置較低的離散度，這樣在方程求解過(guò)程中此參數(shù)不會(huì)發(fā)生較大改變。對(duì)于響應(yīng)也是同樣的道理。

貝葉斯方法的優(yōu)點(diǎn)在于，可以同時(shí)引入待求參數(shù)和響應(yīng)的不確定性。在此基礎(chǔ)上，對(duì)式(10)進(jìn)行反復(fù)迭代計(jì)算，以各階計(jì)算固有頻率與試驗(yàn)結(jié)果誤差的平均值作為收斂標(biāo)準(zhǔn)，最終求解得到合理的修正結(jié)果。算法的整體流程如圖1所示。

圖1 修正算法流程

修正具體算法流程如下：

步驟1 初始化收斂系數(shù)ε>0，初始參數(shù)向量Po，離散度矩陣CP、CR；

步驟2 調(diào)用Nastran，計(jì)算響應(yīng)殘差向量ΔR；

步驟3 計(jì)算響應(yīng)對(duì)參數(shù)的靈敏度矩陣S；

步驟4 根據(jù)式(10)計(jì)算修正方程，得到修正后參數(shù)向量Pu；

步驟7 輸出Pu。

對(duì)于修正算法流程，做如下幾點(diǎn)說(shuō)明：

(1) 在模型修正流程開始之前，先要進(jìn)行靈敏度分析，通過(guò)靈敏度分析的結(jié)果來(lái)選擇合適的修正參數(shù)；

(2) 離散度的設(shè)置需要根據(jù)實(shí)際工程情況，其中參數(shù)的離散度是根據(jù)對(duì)參數(shù)的認(rèn)識(shí)進(jìn)行設(shè)定，試驗(yàn)數(shù)據(jù)的離散度是通過(guò)試驗(yàn)數(shù)據(jù)通過(guò)統(tǒng)計(jì)算法計(jì)算而得；

(3) 對(duì)于步驟3中靈敏度方程的求解，可以根據(jù)相關(guān)理論公式計(jì)算，或者利用有限差分方法進(jìn)行計(jì)算；

(4) 對(duì)于步驟4的修正方程求解，可以根據(jù)參數(shù)與響應(yīng)之間的個(gè)數(shù)關(guān)系來(lái)決定使用式(11)或式(12)進(jìn)行參數(shù)估計(jì)。

對(duì)于一個(gè)模型修正算法，除了算法本身的性能之外，如何將其實(shí)現(xiàn)并應(yīng)用于工程實(shí)際也是一項(xiàng)重要挑戰(zhàn)。本文通過(guò)商用有限元模型修正軟件FEMTools的二次開發(fā)功能，將所提出的算法嵌入到FEMTools修正流程中。

其中，對(duì)于FEMTools的自定義算法實(shí)現(xiàn)，有如下幾個(gè)重點(diǎn)問題：

(1) 由于待修正參數(shù)中，有一些并不在FEMTools的標(biāo)準(zhǔn)庫(kù)中，需要自定義修正參數(shù)，這可以通過(guò)parameter text path 'modal.bdf' line 1254 column 17 length 8 label 'TSM_E11' script 'bulk_format_8.bas'命令來(lái)實(shí)現(xiàn)；

(2) 對(duì)于自定義參數(shù)的實(shí)現(xiàn)，實(shí)質(zhì)上是讀取bdf文件中的對(duì)應(yīng)行，但是需要注意的是，bdf中由于參數(shù)位數(shù)只能為8字節(jié)，需要自行編寫一個(gè)腳本來(lái)處理卡片的數(shù)據(jù)，使其盡量不丟失精度；

(3) 由于熱模態(tài)計(jì)算使用了Nastran的SOL106模塊，需要編寫自定義腳本，令FEMTools讀取其計(jì)算結(jié)果op2文件，可以使用FEMTools的Ft_Import "mode", "nastran.bin", sOP2File命令來(lái)實(shí)現(xiàn)op2文件的讀取。

(4) 具體實(shí)現(xiàn)方法，請(qǐng)參考FEMTools官方手冊(cè)二次開發(fā)部分。

3 復(fù)合材料層合板高溫環(huán)境振動(dòng)試驗(yàn)

本研究的修正對(duì)象為碳纖維/雙馬樹脂基復(fù)合材料層合板，采用試件由16層纖維鋪設(shè)而成，鋪層順序[02/±45/02/90/0]s。尺寸為400 mm×300 mm，厚度2.4 mm。這是一種典型的層合板結(jié)構(gòu)，被廣泛應(yīng)用于飛行器結(jié)構(gòu)。

如圖2，3所示，試驗(yàn)邊界條件為單邊固支。利用專用夾具將試驗(yàn)件固定在工作臺(tái)面上，一面采用石英燈陣進(jìn)行加熱，另一面采用激振器進(jìn)行激勵(lì)。傳感器與激振器均處于燈陣加熱陰影中，這樣可以保護(hù)設(shè)備。由于本實(shí)驗(yàn)關(guān)心高溫環(huán)境下熱應(yīng)力對(duì)結(jié)構(gòu)的影響，為了防止激振桿帶來(lái)的連接剛度，采用隨機(jī)脈沖載荷進(jìn)行激勵(lì)。隨機(jī)脈沖激勵(lì)是指，在傳統(tǒng)沖擊載荷的基礎(chǔ)上，利用隨機(jī)信號(hào)發(fā)生器創(chuàng)造連續(xù)的隨機(jī)脈沖激勵(lì)信號(hào)。當(dāng)信號(hào)的功率譜能夠完整覆蓋所關(guān)心頻段時(shí)，即可通過(guò)這種信號(hào)進(jìn)行基于頻響方法的模態(tài)參數(shù)辨識(shí)。其優(yōu)點(diǎn)在于，相對(duì)傳統(tǒng)脈沖激勵(lì)，能夠獲得更加寬的激勵(lì)頻帶，同時(shí)在安裝上不需要跟試件連接，減小了連接剛度的影響。根據(jù)板型結(jié)構(gòu)的模態(tài)振型特點(diǎn)，對(duì)于測(cè)試前三階模態(tài)，在中線位置均勻布置三個(gè)傳感器。所采用的傳感器為PCB高溫環(huán)境傳感器，適用溫度環(huán)境-60 ℃～240 ℃。

為了適應(yīng)高溫環(huán)境下的操作，采用耐高溫膠對(duì)傳感器進(jìn)行粘貼安裝。采用的耐高溫膠溫度適應(yīng)范圍為-60 ℃～300 ℃，可以滿足本試驗(yàn)的需求。采用德國(guó)M+P振動(dòng)測(cè)試儀進(jìn)行試驗(yàn)數(shù)據(jù)的采集、處理。同時(shí)，利用M+P數(shù)字信號(hào)生成功能產(chǎn)生隨機(jī)信號(hào)，經(jīng)過(guò)功率放大器(江蘇聯(lián)能)放大后由激振器輸出隨機(jī)脈沖激勵(lì)。通過(guò)石英燈陣的控制程序，輸入設(shè)置溫度對(duì)結(jié)構(gòu)進(jìn)行加熱。同時(shí)打開激振器對(duì)結(jié)構(gòu)進(jìn)行激勵(lì)，采集相關(guān)振動(dòng)信號(hào)，對(duì)結(jié)構(gòu)振動(dòng)特性進(jìn)行分析。

圖2 試驗(yàn)概況

圖3 試驗(yàn)概況

試驗(yàn)溫度采用梯度增加的形式，在每個(gè)梯度段采集相關(guān)振動(dòng)信號(hào)，溫度曲線，如圖4所示。

圖4 試驗(yàn)溫度曲線

本試驗(yàn)每個(gè)溫度工況下至少進(jìn)行五次數(shù)據(jù)采集，其中溫度由常溫最高升至180℃。利用M+P的基于頻響函數(shù)的模態(tài)參數(shù)辨識(shí)方法進(jìn)行模態(tài)參數(shù)辨識(shí)。如圖5，其中三個(gè)顏色的頻響曲線分別代表三個(gè)傳感器得到的對(duì)應(yīng)頻響函數(shù)。

將識(shí)別出的各階模態(tài)頻率取平均數(shù)繪制曲線，如圖6所示。

圖5 在M+P中利用試驗(yàn)得到的頻響函數(shù)進(jìn)行模態(tài)辨識(shí)

圖6 模態(tài)頻率變化曲線

4 復(fù)合材料層合板高溫環(huán)境結(jié)構(gòu)動(dòng)力學(xué)建模

利用MSC Patran/Nastran建立目標(biāo)結(jié)構(gòu)的有限元模型，如圖7所示。

圖7 有限元模型

其中，使用四節(jié)點(diǎn)殼單元(QUAD4)作為板結(jié)構(gòu)的單元，使用質(zhì)量點(diǎn)單元來(lái)對(duì)傳感器質(zhì)量進(jìn)行建模。邊界條件方面，在對(duì)應(yīng)約束邊界處添加彈簧單元作為約束條件，從而滿足在邊界上的熱應(yīng)力情況。

材料屬性上，首先使用二維各向異性材料建立單層板的材料屬性，然后使用基于層合板理論的Laminate模塊建立層合板材料，其中鋪層形式按照TSM材料設(shè)計(jì)角度進(jìn)行鋪層，如圖8所示。

熱環(huán)境的建模方面，對(duì)板面整體施加均一的溫度場(chǎng)。需要注意的是，由于是單邊加熱，所以板的兩面溫度并不一致。所以在建模過(guò)程中，需要考慮版中的溫度梯度分布問題。故建立溫度場(chǎng)時(shí)，采用了Element Variable類型建立對(duì)應(yīng)的溫度場(chǎng)。

模型建立完成后，使用MSC Nastran的SOL106非線性靜力計(jì)算模塊進(jìn)行計(jì)算，此計(jì)算模塊會(huì)首先計(jì)算其熱應(yīng)力場(chǎng)，然后形成對(duì)應(yīng)的幾何微分剛度陣添加到系統(tǒng)剛度陣中，最后計(jì)算高溫環(huán)境下的結(jié)構(gòu)模態(tài)頻率及振型。得到的熱應(yīng)力場(chǎng)及前四階振型，如圖9所示。

5 高溫環(huán)境結(jié)構(gòu)參數(shù)靈敏度分析

先對(duì)常溫環(huán)境下的參數(shù)靈敏度進(jìn)行計(jì)算，得到結(jié)果如圖11所示。

其中，parameter代表不同物理參數(shù)，response代表不同階次的模態(tài)頻率。這里的參數(shù)編號(hào)對(duì)應(yīng)，見表1。

圖10 前四階模態(tài)振型(左上：一階，右上：二階，左下：三階，右下：四階)

圖11 常溫環(huán)境靈敏度分析結(jié)果

編號(hào)名稱說(shuō)明1TSM_E11纖維拉伸模量E112TSM_E22纖維拉伸模量E223TSM_S12纖維剪切模量S124TSM_S23纖維剪切模量S235TSM_S13纖維剪切模量S136TSM_EP11熱膨脹系數(shù)EP117TSM_EP22熱膨脹系數(shù)EP22

可以看到，在常溫環(huán)境下熱膨脹系數(shù)對(duì)固有頻率的影響基本為0。而影響最大的為E11即沿著約束方向的拉伸模量。而對(duì)于高溫環(huán)境，靈敏度分析的結(jié)果與傳統(tǒng)的靈敏度分析有所區(qū)別。

圖12 高溫環(huán)境靈敏度分析

可以看到，在高溫環(huán)境下，各參數(shù)的影響發(fā)生了巨大的變化。本來(lái)靈敏度數(shù)值為正的彈性模量等參數(shù)變?yōu)樨?fù)數(shù)，而熱膨脹系數(shù)的影響逐漸加強(qiáng)，這是由于高溫動(dòng)力學(xué)問題本身的非線性特征所帶來(lái)的。

隨著溫度的升高，各個(gè)參數(shù)呈現(xiàn)出非常復(fù)雜的變化規(guī)律，如圖13所示。

圖13 靈敏度隨溫度變化曲線(一階固有頻率)

根據(jù)本節(jié)的分析，所選擇的各項(xiàng)參數(shù)均有著不同的物理意義，且對(duì)模型計(jì)算響應(yīng)都有較大影響。其中，各向拉伸模量TSM_E11、TSM_E22和剪切模量TSM_S12 對(duì)于常溫環(huán)境固有頻率的影響較大。另一方面，隨著溫度的變化，這些彈性參數(shù)也會(huì)發(fā)生變化，這也將帶來(lái)誤差。而對(duì)于高溫環(huán)境模態(tài)計(jì)算，熱膨脹系數(shù)決定了熱應(yīng)力的計(jì)算結(jié)果，其對(duì)于高溫環(huán)境模態(tài)計(jì)算結(jié)果的影響是非常重要的。而之前影響不大的TSM_S23、TSM_E13也體現(xiàn)出了影響。所以，在后續(xù)的模型修正工作中，將繼續(xù)沿用本節(jié)靈敏度分析所使用的參數(shù)作為修正參數(shù)。

6 模型修正

在進(jìn)行高溫環(huán)境模型修正之前，首先對(duì)模型進(jìn)行驗(yàn)證，利用試驗(yàn)得到的前三階固有頻率與有限元模型預(yù)示結(jié)果進(jìn)行對(duì)比，情況如表2所示。

表2 常溫環(huán)境下試驗(yàn)結(jié)果與有限元結(jié)果對(duì)比

可以看到結(jié)果對(duì)應(yīng)良好，說(shuō)明結(jié)構(gòu)的常溫物理參數(shù)以及邊界條件建模準(zhǔn)確。

然后，先根據(jù)之前的離散度計(jì)算公式(13)，計(jì)算各階固有頻率測(cè)試結(jié)果離散度，并代入修正參數(shù)。從實(shí)驗(yàn)結(jié)果離散度的設(shè)置可以看出，二階模態(tài)試驗(yàn)結(jié)果不確定性較低，故其在修正計(jì)算中的占比也較低。

另一方面，按照對(duì)各個(gè)參數(shù)的預(yù)估，設(shè)置修正參數(shù)的離散度。其中，拉伸模量E11在測(cè)量上本身比較準(zhǔn)確，不確定性較小，故離散度設(shè)置較低。剪切模量S23與S13對(duì)響應(yīng)的影響較小，為保證計(jì)算效果設(shè)置較低，其他參數(shù)離散度設(shè)置相對(duì)較高，使其在修正計(jì)算中的占比更高。離散度的設(shè)置相當(dāng)于對(duì)參數(shù)和響應(yīng)的加權(quán)處理，這種處理描述了各個(gè)未知量的不確定性，從而更加合理地描述實(shí)際物理情況。離散度設(shè)置情況如表3所示。

表3 參數(shù)及響應(yīng)離散度設(shè)置

應(yīng)用本研究設(shè)計(jì)的方法，利用FEMtools軟件的二次開發(fā)功能編譯自定義求解序列。分別對(duì)各個(gè)溫度下的模型參數(shù)進(jìn)行修正，以140℃為例，依然選擇靈敏度分析時(shí)選擇的七個(gè)參數(shù)作為修正參數(shù)。修正計(jì)算收斂準(zhǔn)則采用前三階計(jì)算固有頻率與試驗(yàn)值之間誤差的平均值作為收斂目標(biāo)函數(shù)，當(dāng)目標(biāo)函數(shù)小于設(shè)定的標(biāo)準(zhǔn)時(shí)，迭代計(jì)算停止，修正完成。收斂曲線如圖14所示。

圖14 修正計(jì)算收斂曲線

從收斂曲線的情況來(lái)看，雖然整體上達(dá)到了所設(shè)定的誤差目標(biāo)，計(jì)算結(jié)果也比較合理，但是收斂情況并不是非常理想的單一方向下降。這是由于問題的非線性導(dǎo)致的，想要解決這一問題還需要進(jìn)一步研究適用于非線性問題的算法。目前的方法能夠解決現(xiàn)有問題，但還有改進(jìn)的空間。修正前后結(jié)果的對(duì)比情況，見表4和表5。

表4 修正前后結(jié)果對(duì)比情況(140 ℃)

表5 修正前后參數(shù)對(duì)比(140 ℃)

可以看到，經(jīng)過(guò)修正模態(tài)頻率誤差有了明顯下降。同時(shí)，參數(shù)改變主要集中在彈性模量和熱膨脹系數(shù)上。參數(shù)改變量合理，符合預(yù)期。

按照這種方式，對(duì)各個(gè)溫度下的物理參數(shù)均進(jìn)行修正，修正后各階頻率誤差對(duì)比情況，如表6所示。

表6 修正前后誤差對(duì)比

修正前后對(duì)應(yīng)一階、二階固有頻率的溫度變化曲線，如圖15～17所示。

圖15 溫度-固有頻率對(duì)比(一階固有頻率)

可以看到，經(jīng)過(guò)修正各個(gè)溫度下固有頻率的對(duì)應(yīng)情況均有所改善，誤差最大控制在2.68%。修正結(jié)果合理，修正有效。

7 結(jié) 論

本文提出了一種新的基于貝葉斯參數(shù)估計(jì)的高溫環(huán)境結(jié)構(gòu)動(dòng)力學(xué)模型修正方法。通過(guò)貝葉斯參數(shù)估計(jì)，可以將測(cè)試和參數(shù)的不確定性考慮入修正計(jì)算，從而更加準(zhǔn)確地描述實(shí)際問題。利用本方法，對(duì)一復(fù)合材料層合板結(jié)構(gòu)進(jìn)行了模型修正研究，通過(guò)修正改善了在高溫環(huán)境下模型的動(dòng)力學(xué)特性預(yù)示結(jié)果。

圖16 溫度-固有頻率對(duì)比(二階固有頻率)

圖17 溫度-固有頻率對(duì)比(三階固有頻率)

本方法中一個(gè)重要的特征是需要設(shè)定響應(yīng)和參數(shù)的離散度。其中，參數(shù)的離散度可以通過(guò)對(duì)參數(shù)的預(yù)估來(lái)設(shè)置。而對(duì)于響應(yīng)的離散度，可以按照試驗(yàn)情況進(jìn)行設(shè)置，或者按照離散度公式進(jìn)行計(jì)算而得到。

在本研究中，進(jìn)行了高溫環(huán)境下復(fù)合材料結(jié)構(gòu)的動(dòng)力學(xué)試驗(yàn)。得到的結(jié)果基本合理，溫度與固有頻率之間呈現(xiàn)線性對(duì)應(yīng)關(guān)系。但是，這類試驗(yàn)有很多需要注意的地方，包括試件的夾持與安裝需要符合對(duì)應(yīng)的工藝，高溫環(huán)境下部分模態(tài)會(huì)出現(xiàn)混疊現(xiàn)象等。

相對(duì)于傳統(tǒng)各項(xiàng)同性材料，復(fù)合材料在高溫環(huán)境下有著更加復(fù)雜的特性。從靈敏度分析可以看到，參數(shù)靈敏度和溫度之間呈現(xiàn)非線性的對(duì)應(yīng)關(guān)系。這代表參數(shù)和高溫動(dòng)態(tài)響應(yīng)之間也呈現(xiàn)一種非線性的關(guān)系。這對(duì)參數(shù)的設(shè)置和迭代計(jì)算的過(guò)程提出了更高的要求。在后續(xù)的研究中，可以考慮通過(guò)智能優(yōu)化算法來(lái)提升修正運(yùn)算的計(jì)算性能。

[1] THORNTON E A. Thermal structures-four decades of progress[J]. J. Aircr, 1992, 29(3):485-498.

[2] THORNTON E A. Thermal structures for aerospace applications[J]. AIAA, Reston, VA, 1996(2):153:161.

[3] 楊亞政,楊嘉陵,方岱寧. 高超聲速飛行器熱防護(hù)材料與結(jié)構(gòu)的研究進(jìn)展[J]. 應(yīng)用數(shù)學(xué)和力學(xué),2008(1):47-56.

YANG Yazheng, YANG Jialing, FANG Daining. Research progress on the thermal protection materials and structures in hypersonic vehicles[J]. Applied Mathematics and Mechanics, 2008(1):47-56.

[4] 孟松鶴, 杜善義, 韓杰才. 熱防護(hù)系統(tǒng)及材料的研究進(jìn)展[C]//復(fù)合材料——基礎(chǔ), 創(chuàng)新, 高效: 第十四屆全國(guó)復(fù)合材料學(xué)術(shù)會(huì)議論文集 (上). 2006.

[5] 蘇芳, 孟憲紅. 三種典型熱防護(hù)系統(tǒng)發(fā)展概況[J]. 飛航導(dǎo)彈, 2006(10): 57-60.

SU Fang, MENG Xianhong. Research progress on three classical thermal protection systems[J]. Winged Missiles Journal, 2006(10): 57-60.

[6] 郭朝邦,李文杰. 高超聲速飛行器結(jié)構(gòu)材料與熱防護(hù)系統(tǒng)[J]. 工藝與材料, 2010(4):88-94.

GUO Chaobang, LI Wenjie. A research on structural materials and thermal protection systems of hypersonic vehicles[J]. Processing and Materials, 2010(4):88-94.

[7] 夏德順.重復(fù)運(yùn)載器金屬熱防護(hù)系統(tǒng)的評(píng)述[J].導(dǎo)彈與航天運(yùn)載技術(shù)，2002(2):21-26.

XIA Deshun. A review on thermal protection systems for reusable vehicles[J]. Missiles and Astronautical Technology, 2002(2):21-26.

[8] BAILEY C D. Vibration of thermally stressed plates with various boundary conditions[J]. Am. Inst. Aeronaut. Astronaut, 1973,11(1):14-19.

[9] CHANG W P, SHOOU-CHIAN J. Nonlinear free vibration of heated orthotropic rectangular plates[J]. Int. J. Solids Struct., 1986, 22(3):267-281.

[10] CHANG W P, WAN S M. Thermomechanically coupled non-linear vibration of plates[J]. Int. J. Non-Linear Mecha. 1986, 21(5):375-389.

[11] BROWN A M. Temperature-dependent modal test/analysis correlation of X-34 FASTRAC composite rocket nozzle[J]. J. Propuls. Power, 2002, 18(2):284-288.

[12] 劉芹, 任建亭, 姜節(jié)勝, 等. 復(fù)合材料層合板非線性熱振動(dòng)分析[J]. 動(dòng)力學(xué)與控制學(xué)報(bào), 2005, 3(1): 78-83.

LIU Qin, REN Jianting, JIANG Jiesheng, et al. Nonlinear thermal vibration characteristic analysis composite laminated plates[J]. Journal of Dynamics and Control, 2005, 3(1): 78-83.

[13] 劉芹, 任建亭, 姜節(jié)勝, 等. 復(fù)合材料薄壁圓柱殼熱振動(dòng)特性分析[J]. 機(jī)械強(qiáng)度, 2006, 28(5): 643-648.

LIU Qin, REN Jianting, JIANG Jiesheng, et al. Nonlinear thermal vibration characteristic analysis of composite thin-cylindrical shells[J]. Journal of Mechanical Strength, 2006, 28(5): 643-648.

[14] DANTAS L B, ORLANDE H R B. A function estimation approach for determining temperature-dependent thermo-physical properties[J]. Inverse Probl. Eng., 1996, 3(4): 261-279.

[15] YANG C. Determination of the temperature dependent thermophysical properties from temperature responses measured at medium’s boundaries[J]. Int.J. Heat Mass Transf, 2000, 43(7): 1261-1270.

[16] LIU C. An efficient simultaneous estimation of temperature-dependent thermophysical properties[M]. CMC-Tech science press, 2006, 77.

[17] CZéL B, GRF G. Simultaneous identification of temperature-dependent thermal properties via enhanced genetic algorithm[J]. Int. J. Thermophys, 2012, 33(6): 1023-1041.

[19] MAHNKEN R. An inverse finite-element algorithm for parameter identification of thermoelastic damage models[J]. Int. J. Numer. Methods Eng., 2000, 48(7): 1015-1036.

[20] MERAGHNI F, CHEMISKY Y, PIOTROWSKI B, et al. Parameter identification of a thermodynamic model for superelastic shape memory alloys using analytical calculation of the sensitivity matrix, Eur[J]. J. Mech. -A/Solids, 2014, 45: 226-237.

[21] POTTIER T, TOUSSAINT F, LOUCHE H, et al. Inelastic heat fraction estimation from two successive mechanical and thermal analyses and full-field measurements[J]. Eur. J. Mechan.-A/Solids, 2013, 38: 1-11.

[22] 邢晶晶, 張能輝, 王明祿. 黏彈性薄板的熱振動(dòng)[J]. 振動(dòng)與沖擊, 2005, 24(5): 23-25.

XING Jingjing, ZHANG Nenghui, WANG Minglu. Thermal vibration of thin viscoelastic plates[J]. Journal of Vibration and Shock, 2015, 24(5): 23-25.

[23] 劉芹, 任建亭, 郭運(yùn)強(qiáng), 等. 埋入形狀記憶合金絲的復(fù)合材料圓柱殼壁板的熱振動(dòng)特性分析[J]. 振動(dòng)與沖擊, 2007, 26(7): 18-21.

LIU Qin, REN Jianting, GUO Yunqiang, et al. Thermal vibration characteristics analysis of a composite thin-cylindrical panel embedded with shape memory alloy wires[J]. Journal of Vibration And Shock, 2007, 26(7): 18-21.

[24] MOTTERSHEAD J E, FRISWELL E. Model updating in structural dynamics: a survey[J]. Journal of Sound and Vibration, 1993, 167(2): 347-375.

[25] MICHAEL F, MOTTERSHEAD J E. Finite element model updating in structural dynamics[M]. Springer, 1995，38.

[26] 姜東,費(fèi)慶國(guó),吳邵慶. 基于攝動(dòng)法的不確定性有限元模型修正方法研究[J]. 計(jì)算力學(xué)學(xué)報(bào),2014(4):431-437.

JIANG Dong, FEI Qingguo, WU Shaoqing. Perturbation based uncertain finite element model updating method[J]. Journal of Computational Mechanics, 2014(4):431-437.

[27] 方圣恩,林友勤,夏樟華. 考慮結(jié)構(gòu)參數(shù)不確定性的隨機(jī)模型修正方法[J]. 振動(dòng).測(cè)試與診斷,2014(5):832-837.

FANG Shengen, LIN Youqin, XIA Zhanghua. Stochastic model updating method for structural parameters considering uncertainty[J]. Journal of Vibration，Measurement & Diagnosis, 2014(5):832-837.

[28] 陳學(xué)前,杜強(qiáng),馮加權(quán). 基于結(jié)構(gòu)參數(shù)不確定性的有限元模型修正[J]. 機(jī)械科學(xué)與技術(shù),2007(8):965-968.

CHEN Xueqian, DU Qiang, FENG Jiaquan. Modification of Finite Element Model (FEM) considering the uncertainty of structural parameters[J]. Mechanical Science And Technology For Aerospace Engineering, 2007(8):965-968.

[29] 房長(zhǎng)宇,張耀庭. 基于參數(shù)不確定性的預(yù)應(yīng)力混凝土梁模型修正[J]. 華中科技大學(xué)學(xué)報(bào)(自然科學(xué)版),2011,11:87-91.

FANG Changyu, ZHANG Yaoting. Model updating of prestressed concrete beams using parameters uncertainty[J]. Journal of Huazhong University of Science And Technology Nature Science, 2011,11:87-91.

[30] 張秋虎. 考慮結(jié)構(gòu)參數(shù)不確定性的有限元模型修正方法研究[D].合肥：合肥工業(yè)大學(xué),2014.

[31] 姜東,費(fèi)慶國(guó),吳邵慶. 基于區(qū)間分析的不確定性結(jié)構(gòu)動(dòng)力學(xué)模型修正方法[J]. 振動(dòng)工程學(xué)報(bào),2015(3):352-358.

JIANG Dong, FEI Qingguo, WU Shaoqing. Updating of structural dynamics model with uncertainty based on interval analysis[J]. Journal of Vibration Engineering, 2015(3):352-358.

[32] 何成,陳國(guó)平,何歡. 徑向基模型的不確定性模型區(qū)間修正與確認(rèn)[J]. 機(jī)械工程學(xué)報(bào),2013,11:128-134.

HE Cheng, CHEN Guoping, HE Huan. Interval model updating and validation with uncertainty based on the radial basis function[J]. Journal of Mechanical Engineering, 2013,11:128-134.

[33] 張建新. 基于貝葉斯方法的有限元模型修正研究[D].重慶：重慶大學(xué),2014.

[34] 韓芳,鐘冬望,汪君. 基于貝葉斯法的復(fù)雜有限元模型修正研究[J]. 振動(dòng)與沖擊,2012，31(1):39-43.

HAN Fang, ZHONG Dongwang, WANG Jun. Complicated finite element model updating based on Bayesian method[J]. Journal of Vibration and Shock, 2012，31(1):39-43.

[35] 張明亮. 基于貝葉斯理論的材料非線性橋梁結(jié)構(gòu)模型修正與損傷識(shí)別[D].長(zhǎng)春：吉林大學(xué),2011.

[36] 萬(wàn)華平,任偉新,黃天立. 基于貝葉斯推理的隨機(jī)模型修正方法[J]. 中國(guó)公路學(xué)報(bào),2016(4):67-76.

WAN Huaping, REN Weixin, HUANG Tianli. Stochastic model updating approach by using bayesian inference[J]. China Journal of Highway and Transport, 2016(4):67-76.

[37] HUMBERT L, THOUVEREZ F, JéZéQUEL L. Finite element dynamic model updating using modal thermoelastic fields[J]. Journal of Sound and Vibration, 1999, 228(2): 397-420.

[38] HE C, CHEN G, HE H, et al. Model updating of a dynamic system in a high-temperature environment based on a hierarchical method[J]. Finite Elements in Analysis and Design, 2013, 77: 59-68.

[39] 何成. 高溫環(huán)境下結(jié)構(gòu)動(dòng)力學(xué)建模關(guān)鍵技術(shù)研究[D].南京：南京航空航天大學(xué),2014.

[40] SUN Kaipeng, ZHAO Yonghui, HU Haiyan. Identification of temperature-dependent thermal-structural properties via finite element model updating and selection[J]. Mechanical Systems and Signal Processing, 2015, 52: 147-161.

[41] 張保強(qiáng). 熱結(jié)構(gòu)不確定性動(dòng)力學(xué)仿真及模型確認(rèn)方法研究[D].南京：南京航空航天大學(xué),2012.

Structural dynamic model updating under high temperature environment

YUAN Zhaoxu, YU Kaiping

(School of Aeronautics, Harbin Institute of Technology, Harbin 150001, China)

A new Bayesian parametric estimation-based model updating method for structures under high temperature environment (HTE) was proposed here. It was used to update the dynamic model of an actual composite laminate structure. The main idea of this method was to update the FE model of a structure with its modal frequencies test data under HTE. The feature of this method was using the dispersion of parameters and responses to describe their uncertainty and then to solve actual problems more reasonably. The dynamic tests of the composite structure were conducted under HTE to obtain its modal frequencies under HTE. The test data were used to update the FE model of the composite structure. Before updating, a sensitivity analysis was performed for the model to explore relations among its parameters and the structure’s dynamic responses. Test data and pre-estimations of parameters to be updated were used to set up their dispersions. Finally, the updating computation was done with a self-programmed algorithm. After updating, the relations among temperature and modal frequencies of the model were improved obviously to verify the effectiveness of the proposed method and its applicability for actual engineering problems.

structural dynamics under high temperature environment; model updating; composite material; sensitivity analysis

國(guó)家自然科學(xué)基金(11172078)

2016-03-08 修改稿收到日期：2016-06-28

袁昭旭 男，博士生，1988年生

于開平 男，博士，教授，1968年生

O327

A

10.13465/j.cnki.jvs.2017.15.026