區(qū)間參數(shù)梁結(jié)構(gòu)熱彈耦合效應(yīng)分析

云永琥， 陳建軍， 趙 寬， 閻 彬， 曹鴻鈞

(西安電子科技大學(xué) 電子裝備結(jié)構(gòu)設(shè)計(jì)教育部重點(diǎn)實(shí)驗(yàn)室，陜西 西安 710071)

近年來，結(jié)構(gòu)在突加熱載作用下的動(dòng)力響應(yīng)問題引起了人們的關(guān)注．例如在航天器結(jié)構(gòu)周期性進(jìn)出陰影區(qū)和日照區(qū)時(shí)，其中柔性結(jié)構(gòu)由于瞬態(tài)溫度梯度的巨變將有可能導(dǎo)致熱致振動(dòng)[1]．反之，結(jié)構(gòu)的變形也會(huì)引起其熱輻射邊界條件的改變，從而導(dǎo)致結(jié)構(gòu)溫度場(chǎng)的變化，顯然，此時(shí)結(jié)構(gòu)的溫度場(chǎng)與變形場(chǎng)是相互耦合的．但在計(jì)算精度要求不高的情況下，對(duì)于受到瞬態(tài)溫變作用的結(jié)構(gòu)動(dòng)力學(xué)分析，通?？梢院雎云錈釓楍詈闲?yīng)，在計(jì)算中只需將溫度梯度產(chǎn)生的熱彎矩當(dāng)作等效載荷引入到動(dòng)力學(xué)方程中即可．如文獻(xiàn)[2]用有限元方法對(duì)矩形板進(jìn)行了熱致振動(dòng)分析，并通過實(shí)驗(yàn)驗(yàn)證了該方法的有效性．文獻(xiàn)[3]對(duì)層疊板的熱誘發(fā)振動(dòng)進(jìn)行了有限元計(jì)算．文獻(xiàn)[4]對(duì)結(jié)構(gòu)的熱動(dòng)力響應(yīng)和熱致振動(dòng)等問題進(jìn)行了分析．以上文獻(xiàn)均未考慮結(jié)構(gòu)變形對(duì)溫度場(chǎng)的影響，而在計(jì)算精度要求較高的情況下，結(jié)構(gòu)中熱彈耦合項(xiàng)的作用將不能忽視．文獻(xiàn)[5]提出了耦合系數(shù)的概念，并對(duì)梁結(jié)構(gòu)進(jìn)行了熱彈耦合分析．在Boley工作的基礎(chǔ)上，國(guó)內(nèi)外學(xué)者相繼開展了一些研究．文獻(xiàn)[6]對(duì)功能梯度梁的熱彈耦合進(jìn)行了有限元分析．文獻(xiàn)[7]利用有限元方法分析了二維平面板考慮熱彈耦合的熱動(dòng)力響應(yīng)問題，給出了熱彈耦合項(xiàng)對(duì)響應(yīng)的影響．文獻(xiàn)[8]在考慮熱彈耦合效應(yīng)下，分析了Timoshenko梁同時(shí)受到熱載荷和力載荷作用時(shí)的動(dòng)力特性，并提出了求解耦合動(dòng)力學(xué)方程的方法．文獻(xiàn)[9]對(duì)輻射換熱條件下的空間薄壁桿件的熱-動(dòng)力學(xué)耦合問題進(jìn)行了有限元分析，得到了熱載荷對(duì)結(jié)構(gòu)響應(yīng)的影響．

目前，在結(jié)構(gòu)熱彈耦合效應(yīng)分析中，研究對(duì)象基本上為確定性參數(shù)模型．然而，由于結(jié)構(gòu)在制造過程中的各種不確定性因素以及誤差的影響，導(dǎo)致結(jié)構(gòu)的物性參數(shù)具有一定的不確定性，并且結(jié)構(gòu)所受的熱載荷和力載荷有時(shí)亦具有不確定性，此時(shí)結(jié)構(gòu)的動(dòng)力響應(yīng)也將呈現(xiàn)出不確定性．因此，分析結(jié)構(gòu)所受載荷及其物性參數(shù)的不確定性對(duì)動(dòng)力響應(yīng)的影響，其結(jié)果較確定性問題的解無疑更加符合客觀實(shí)際．文獻(xiàn)[10]將不確定結(jié)構(gòu)參數(shù)視為隨機(jī)變量，利用隨機(jī)因子法推導(dǎo)出結(jié)構(gòu)動(dòng)力響應(yīng)的均值和方差的計(jì)算表達(dá)式，對(duì)隨機(jī)參數(shù)彈性桿在瞬態(tài)溫度場(chǎng)下的響應(yīng)問題進(jìn)行了分析．文獻(xiàn)[11-12]利用區(qū)間分析方法對(duì)區(qū)間不確定性結(jié)構(gòu)的動(dòng)態(tài)特性和動(dòng)態(tài)響應(yīng)進(jìn)行了分析．

筆者針對(duì)Euler-Bernouli梁結(jié)構(gòu)，考慮其物性參數(shù)、溫度和外力載荷均為區(qū)間變量，建立了其熱彈耦合動(dòng)力學(xué)區(qū)間有限元方程，提出了結(jié)構(gòu)熱彈耦合動(dòng)力響應(yīng)范圍的區(qū)間計(jì)算方法．并通過算例證明了該方法的可行性與有效性．考察了結(jié)構(gòu)的區(qū)間物性參數(shù)對(duì)梁動(dòng)力響應(yīng)的影響，以及結(jié)構(gòu)受力變形對(duì)熱彈耦合效應(yīng)的影響．

圖1 懸臂梁結(jié)構(gòu)圖

1 熱彈耦合懸臂梁有限元方程的建立與求解

以圖1矩形截面懸臂梁為分析對(duì)象，環(huán)境溫度為T0，在t=0的初始時(shí)刻，梁的上端面同時(shí)受到階躍均布力f和熱流q的共同作用，下端面絕熱，忽略固定端面的熱傳導(dǎo)和自由端截面以及前后表面的換熱．由于熱流在梁上端面均勻分布，沿著梁軸向任意截面上的溫度處處相等，因此熱流僅沿梁厚度方向傳導(dǎo)．鑒于此，求解該梁熱彈耦合動(dòng)力響應(yīng)將需分別構(gòu)建以下兩種有限元模型．

1.1 梁動(dòng)力學(xué)有限元模型

圖2 熱分析模型

通過對(duì)各個(gè)單元進(jìn)行組集，可得梁的動(dòng)力學(xué)有限元方程為[13]

(1)

其中，u為結(jié)構(gòu)的位移向量，M和K分別為結(jié)構(gòu)的總質(zhì)量矩陣和總剛度矩陣，F(xiàn)B為總力載荷向量，F(xiàn)T為總溫度載荷向量．

1.2 梁熱分析有限元模型

梁熱分析模型如圖2所示．將梁橫截面沿厚度方向離散分為m個(gè)單元和m+1個(gè)節(jié)點(diǎn)．由自由能密度和熵密度的計(jì)算，并利用最小勢(shì)能原理，導(dǎo)出一維耦合熱傳導(dǎo)有限元方程為[7]

其中，T為待求的節(jié)點(diǎn)瞬態(tài)溫度向量；C為熱容矩陣；KT為熱傳導(dǎo)剛度矩陣；P為節(jié)點(diǎn)熱載荷列陣；N為單元節(jié)點(diǎn)溫度的形函數(shù)；B為單元節(jié)點(diǎn)位移應(yīng)變矩陣；H為溫度與彈性變形的耦合矩陣項(xiàng)，它表明溫度場(chǎng)不僅與熱源、熱力學(xué)物性參數(shù)及換熱邊界條件有關(guān)，還受到彈性變形應(yīng)變率的影響，其在一定程度上改變物體內(nèi)部的熱量傳遞；k為熱傳導(dǎo)系數(shù)；ρ為質(zhì)量密度；c為比熱容；μ為泊松比；q為熱流量；T0為結(jié)構(gòu)初始溫度．

1.3 熱彈耦合下動(dòng)力響應(yīng)求解

對(duì)于文中熱彈耦合動(dòng)力學(xué)有限元方程求解的問題，采用聯(lián)立求解式(1)和式(2)，并相互交替迭代的計(jì)算方法，其中，運(yùn)動(dòng)方程式(1)則由Newmark法求解，熱傳導(dǎo)方程式(2)由Galerkin迭代格式進(jìn)行求解．

對(duì)熱傳導(dǎo)方程式(2)，取θ=2/3的無條件穩(wěn)定迭代格式[14]如下:

(3)

而動(dòng)力學(xué)方程式(1)則利用Newmark法將轉(zhuǎn)化為求解如下的擬靜力方程:

其中，γ是根據(jù)積分的精度和穩(wěn)定性要求給定的可調(diào)參數(shù)．當(dāng)γ=0.25時(shí)，積分格式無條件穩(wěn)定[14]．

2 區(qū)間參數(shù)熱彈耦合響應(yīng)分析

考慮到結(jié)構(gòu)制造誤差和外界環(huán)境等多種不確定性因素的影響，將結(jié)構(gòu)的物性參數(shù)ρ、c、k、α、μ和E等均視為區(qū)間參數(shù)，同時(shí)，將結(jié)構(gòu)所受外力f、溫度載荷q以及結(jié)構(gòu)初始溫度T0亦視為區(qū)間參數(shù)，將它們統(tǒng)一以區(qū)間參數(shù)向量形式表示為β= (β1，β2，…，βm)T，β中既有區(qū)間結(jié)構(gòu)參數(shù)也有區(qū)間載荷參數(shù)，其所在的范圍為

(7)

將式(7)代入方程(1)和(2)中，則得描述區(qū)間參數(shù)梁結(jié)構(gòu)的動(dòng)力學(xué)和瞬態(tài)溫度場(chǎng)方程如下:

由式(8)和式(9)可見，當(dāng)結(jié)構(gòu)參數(shù)和載荷均具有區(qū)間不確定性時(shí)，梁結(jié)構(gòu)的時(shí)變位移響應(yīng)u(β，t)和瞬態(tài)溫度場(chǎng)T(β，t)將分別是兩個(gè)區(qū)間參數(shù)時(shí)變函數(shù)的集合，即

u(β，t)和T(β，t)的上下界可表示為

以下將給出求解區(qū)間參數(shù)結(jié)構(gòu)熱彈耦合動(dòng)力方程的響應(yīng)上下界的計(jì)算方法．

3 熱彈耦合區(qū)間結(jié)構(gòu)有限元方程的求解方法

由于區(qū)間有限元方程的參數(shù)是區(qū)間變量，所有區(qū)間變量在各自的區(qū)間范圍內(nèi)的取值和變量的分布類型未知，對(duì)此問題為能夠近似有效求解，可假設(shè)各變量在其區(qū)間范圍內(nèi)均服從具有最大熵的矩形分布[15]．此假設(shè)的理由在于，將變量取到區(qū)間兩端點(diǎn)的概率密度等同于取到其中值點(diǎn)的概率密度，從而將得到最為保守的計(jì)算結(jié)果，這對(duì)于結(jié)構(gòu)的可靠性預(yù)測(cè)和設(shè)計(jì)結(jié)果是最為安全的．基于此假設(shè)，則可利用蒙特卡羅方法對(duì)每一區(qū)間變量在給定的區(qū)間內(nèi)生成均勻分布的隨機(jī)樣本，進(jìn)而再按照確定性有限元模型進(jìn)行計(jì)算，最后得到結(jié)構(gòu)響應(yīng)的區(qū)間范圍．

(1) 確定結(jié)構(gòu)中獨(dú)立的各區(qū)間變量及抽樣次數(shù)l．

(3) 對(duì)于當(dāng)前時(shí)間步長(zhǎng)，隨機(jī)抽取1～N之間的正整數(shù)并提取出其所對(duì)應(yīng)子區(qū)間的兩個(gè)端點(diǎn)值，記當(dāng)前抽樣次數(shù)i=1．(令f1和f2用來存儲(chǔ)當(dāng)前抽樣所計(jì)算出來的最大和最小值，fmax和fmin用來存儲(chǔ)當(dāng)前時(shí)間步下的最大和最小值)．

(4) 將各個(gè)參數(shù)子區(qū)間的端點(diǎn)值進(jìn)行組合，代入式(3)和式(4)，計(jì)算出函數(shù)值并比較．令當(dāng)前計(jì)算結(jié)果中f1為最大值，f2為最小值，同時(shí)令fmax=f1，fmin=f2．

(5) 當(dāng)il時(shí)，則轉(zhuǎn)入步驟(6)．

(6) 結(jié)束循環(huán)抽樣，輸出當(dāng)前時(shí)間步的最大與最小的響應(yīng)fmax和fmin．

(7) 返回步驟(3)，對(duì)下一時(shí)間迭代步進(jìn)行新的抽樣計(jì)算．

當(dāng)在給定時(shí)間域上逐步計(jì)算完成之后，則得到整個(gè)時(shí)間域上動(dòng)力響應(yīng)的最大和最小包絡(luò)線．

4 算 例

如圖1所示的矩形截面懸臂深梁，長(zhǎng)L=500 mm，寬b=100 mm，高h(yuǎn)=50 mm，材料為鋁，梁結(jié)構(gòu)的區(qū)間參數(shù)分別為:ρ= [2 560，2 730] kg/m3，k= [190，210] W/(m·K)，α= [20×10-6，23×10-6]/℃，μ= [0.31，0.35]，q= [1.55×106，1.71×106] W/m2，c= [825，915] J/(kg·℃)，E= [70×109，77×109] Pa，T0= [21，22]℃，f= [9，11] N/m2，其中梁上端面熱流q作用的時(shí)間為 10 s．

求解懸臂梁熱彈耦合動(dòng)力響應(yīng)分別采用兩種有限元模型進(jìn)行計(jì)算．圖1為梁的動(dòng)力學(xué)分析模型，沿其長(zhǎng)度方向被離散為4個(gè)單元、5個(gè)節(jié)點(diǎn)，利用結(jié)構(gòu)動(dòng)力學(xué)方程計(jì)算得出節(jié)點(diǎn)的動(dòng)力響應(yīng)．圖2為梁的熱分析模型，沿梁橫截面的厚度方向亦被離散為4個(gè)單元、5個(gè)節(jié)點(diǎn)，利用熱傳導(dǎo)有限元方程計(jì)算得出截面溫差，從而計(jì)算出梁截面的熱彎矩．

圖3為梁自由端橫截面中點(diǎn)溫度場(chǎng)均值分別在耦合與非耦合(即耦合項(xiàng)H為0時(shí)，結(jié)構(gòu)變形對(duì)溫度場(chǎng)未產(chǎn)生影響)情況下的時(shí)間歷程計(jì)算結(jié)果．從圖3可以看出，耦合效應(yīng)對(duì)橫截面中點(diǎn)溫度分布的影響，這是由于耦合與結(jié)構(gòu)變形有關(guān)．由式(2)可知，由于耦合項(xiàng)的作用，溫度場(chǎng)產(chǎn)生小幅度的波動(dòng)．

圖3 梁自由端中點(diǎn)溫度的時(shí)間歷程圖4 梁自由端節(jié)點(diǎn)位移的時(shí)間歷程

圖4為梁自由端節(jié)點(diǎn)位移響應(yīng)均值隨時(shí)間變化歷程的計(jì)算結(jié)果．在不考慮熱彈耦合效應(yīng)時(shí)，位移響應(yīng)在各個(gè)周期的振動(dòng)狀態(tài)相同．而在考慮熱彈耦合效應(yīng)時(shí)，位移振幅隨時(shí)間逐漸減弱，這種現(xiàn)象說明了耦合對(duì)振動(dòng)產(chǎn)生的抑制作用．由阻尼的效應(yīng)可知，這種減小的變化實(shí)際上相當(dāng)于阻尼的作用．

圖5給出了僅當(dāng)E增大(取E=210×109Pa)以及E和α兩者都增大(取E=210×109Pa 和α=60×10-6/℃) 的情況下，自由端節(jié)點(diǎn)位移均值響應(yīng)的時(shí)間歷程．由圖5可見，E的增大使耦合效果明顯，振幅衰減迅速．而當(dāng)E和α都增大時(shí)，耦合效應(yīng)不僅影響了振幅，還使得平衡位置發(fā)生了改變．說明了E和α對(duì)耦合的影響較大．該計(jì)算結(jié)果也可以通過熱傳導(dǎo)方程式(2)中的熱彈耦合項(xiàng)的表達(dá)式進(jìn)行分析．當(dāng)參數(shù)E和α增大時(shí)，熱彈耦合項(xiàng)就越大，其所表現(xiàn)的阻尼作用就越明顯．

圖5 E和α變化時(shí)梁自由端節(jié)點(diǎn)位移的時(shí)間歷程圖6 梁自由端節(jié)點(diǎn)位移響應(yīng)區(qū)間的時(shí)間歷程

圖6給出了當(dāng)梁各參數(shù)取區(qū)間變量時(shí)，利用區(qū)間有限元模型計(jì)算獲得的自由端節(jié)點(diǎn)位移響應(yīng)的區(qū)間時(shí)間歷程．與確定性模型的計(jì)算方法相比，采用區(qū)間的計(jì)算方法可以得到更多的計(jì)算結(jié)果信息，不僅給出了結(jié)構(gòu)位移響應(yīng)的中值，亦可獲得位移響應(yīng)的區(qū)間范圍．

以上算例均考慮了階躍均布力f作用時(shí)的響應(yīng)情況．為了考察該作用力對(duì)熱彈耦合效應(yīng)的影響，對(duì)無作用力和有作用力兩種情況分別計(jì)算懸臂梁自由端位移響應(yīng)均值，其計(jì)算結(jié)果見表1．從表1可見，當(dāng)梁承受作用力時(shí)，變形比無承受作用力時(shí)的大，此時(shí)耦合效應(yīng)的影響也更為明顯．因此，當(dāng)結(jié)構(gòu)同時(shí)受到力和熱載荷共同作用時(shí)，則需考慮耦合效應(yīng)才能得到符合實(shí)際情況的計(jì)算結(jié)果．

表1 有無作用力時(shí)對(duì)懸臂梁自由端撓度的影響

5 結(jié) 論

研究了結(jié)構(gòu)參數(shù)、力和熱載荷同時(shí)具有區(qū)間不確定性時(shí)，梁結(jié)構(gòu)在熱彈耦合情況下的動(dòng)力響應(yīng)求解問題．通過懸臂梁算例，獲得了以下結(jié)論:

(1) 在已知結(jié)構(gòu)參數(shù)和載荷關(guān)于不確定性信息較少時(shí)，利用區(qū)間分析模型，并通過文中改進(jìn)的蒙特卡羅數(shù)值仿真方法可以得到結(jié)構(gòu)響應(yīng)的區(qū)間范圍．該方法無需不確定參數(shù)的概率信息，為解決區(qū)間參數(shù)梁結(jié)構(gòu)的熱彈耦合計(jì)算問題提供了一種途徑．

(2) 在不考慮熱與結(jié)構(gòu)變形耦合的情況下，溫度不受變形的影響．由于耦合效應(yīng)，則結(jié)構(gòu)變形對(duì)其溫度場(chǎng)的分布產(chǎn)生影響，表現(xiàn)為溫度的明顯波動(dòng)．

(3) 熱彈耦合效應(yīng)對(duì)結(jié)構(gòu)位移振動(dòng)有抑制作用，使結(jié)構(gòu)的振幅隨時(shí)間不斷減小，并趨于穩(wěn)態(tài)的平衡位置．

(4) 結(jié)構(gòu)的彈性模量和熱膨脹系數(shù)與耦合效應(yīng)密切相關(guān)，兩參數(shù)越大，熱彈耦合效應(yīng)越明顯．

(5) 在有力的作用時(shí)，不考慮熱彈耦合效應(yīng)對(duì)結(jié)構(gòu)位移的影響將可能造成一定的計(jì)算誤差．因此，在高精度要求的結(jié)構(gòu)分析和設(shè)計(jì)中，考慮熱彈耦合效應(yīng)的影響是十分必要的．

[1] Thornton E A, Kim Y A. Thermally Induced Bending Vibrations of a Flexible Rolled-up Solar Array[J]. Journal of Spacecraft and Rockets, 1993, 30(4): 438-448.

[2] Lee Y S, Jeon B H, Kang H W. Free Vibration Characteristics of Thermally Loaded Rectangular Plates[J]. Key Engineering Materials, 2011, 478(81): 81-86.

[3] 趙壽根, 王靜濤, 黎康, 等. 考慮輻射散熱疊層板熱誘發(fā)振動(dòng)的有限元分析[J]. 力學(xué)學(xué)報(bào), 2010, 42(5): 978-982.

Zhao Shougen, Wang Jingtao, Li Kang, et al. Finite Element Method Analysis of Thermally Induced Vibration of Laminated Plates Considering Radiation[J]. Chinese Journal of Theoretical and Applied Mechanics, 2010, 42(5): 978-982.

[4] 范緒箕. 高速飛行器熱結(jié)構(gòu)分析與應(yīng)用[M]. 北京: 國(guó)防工業(yè)出版社, 2008.

[5] Boly B A. Thermally Induced Vibrations of Beams[J]. Journal of the Aeroautical Science, 1956, 23(2): 179-181.

[6] Babaei M H, Abbasi M, Eslami M R. Coupled Thermoelasticity of Functionally Graded Beams[J]. Journal of Thermal Stresses, 2008, 31(8): 680-697.

[7] 李智勇, 劉錦陽, 洪嘉振. 作平面運(yùn)動(dòng)的二維平面板的熱耦合動(dòng)力學(xué)問題[J]. 動(dòng)力學(xué)與控制學(xué)報(bào), 2006, 4(2): 114-121.

Li Zhiyong, Liu Jinyang, Hong Jiazhen. Coupled Thermoelastic Dynamics of a Two-dimensional Plate Undergoing Planar Motion[J]. Journal of Dynamics and Control, 2006, 4(2): 114-121.

[8] Manoach E, Ribeiro P. Coupled, Thermoelastic, Large Amplitude Vibrations of Timoshenko Beams[J]. International Journal of Mechanical Sciences, 2004, 46(11): 1589-1606.

[9] 薛明德, 向志海. 大型空間結(jié)構(gòu)的熱-動(dòng)力學(xué)耦合問題及其有限元分析[J]. 固體力學(xué)學(xué)報(bào), 2011, 32(s1): 318-328.

Xue Mingde, Xiang Zhihai. Thermal-dynamic Coupling Problem of Large Space Structures and Its FEM Analysis[J]. Chinese Journal of Solid Mechanics, 2011, 32(s1): 318-328.

[10] 閻彬, 陳建軍, 馬洪波. 隨機(jī)彈性桿在隨機(jī)瞬態(tài)溫度場(chǎng)下的熱響應(yīng)分析[J]. 電子科技大學(xué)學(xué)報(bào), 2012, 41(4): 631-636.

Yan Bin, Chen Jianjun, Ma Hongbo. Thermal Response Analysis of Stochastic Pole Structures under Random Transient Temperature Field[J]. Journal of University of Electronic Science and Technology of China, 2012, 41(4): 631-636.

[11] 王敏娟, 陳建軍, 林立廣, 等. 區(qū)間參數(shù)智能梁結(jié)構(gòu)閉環(huán)系統(tǒng)動(dòng)態(tài)特性分析[J]. 電子科技大學(xué)學(xué)報(bào), 2011, 40(1): 152-156.

Wang Minjuan, Chen Jianjun, Lin Liguang, et al. Dynamic Characteristic Analysis of Closed Loop Systemsfor the Intelligent Beam with Interval Parameters[J]. Journal of University of Electronic Science and Technology of China, 2011, 40(1): 152-156.

[12] 梁震濤, 陳建軍, 朱增青, 等. 不確定結(jié)構(gòu)動(dòng)力區(qū)間分析方法研究[J]. 應(yīng)用力學(xué)學(xué)報(bào), 2008, 25(1): 46-50.

Liang Zhentao, Chen Jianjun, Zhu Zengqing, et al. Dynamic Interval Analysis for Uncertain Structures[J]. Chinese Journal of Applied Mechanics, 2008, 25(1): 46-50.

[13] Marakala N, Kadoli R. Thermally Induced Vibration of a Simply Supported Beam Using Finite Element Method[J]. International Journal of Engineering Science, 2010, 2(12): 7874-7879.

[14] 王勖成. 有限單元法[M]. 北京: 清華大學(xué)出版社, 2003.

[15] 易平. 對(duì)區(qū)間不確定性問題的可靠性度量的探討[J]. 計(jì)算力學(xué)學(xué)報(bào), 2006, 23(2): 152-156.

Yi Ping. Discussions on Reliability Measure for Problems with Bounded-but-unknown Uncertainties[J]. Chinese Journal of Computational Mechanics, 2006, 23(2): 152-156.