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

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

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

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

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

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

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

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

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

圖2 熱分析模型

通過對各個單元進行組集，可得梁的動力學(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個單元和m+1個節(jié)點．由自由能密度和熵密度的計算，并利用最小勢能原理，導(dǎo)出一維耦合熱傳導(dǎo)有限元方程為[7]

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

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

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

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

(3)

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

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

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

考慮到結(jié)構(gòu)制造誤差和外界環(huán)境等多種不確定性因素的影響，將結(jié)構(gòu)的物性參數(shù)ρ、c、k、α、μ和E等均視為區(qū)間參數(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)的動力學(xué)和瞬態(tài)溫度場方程如下:

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

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

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

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

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

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

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

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

(5) 當il時，則轉(zhuǎn)入步驟(6)．

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

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

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

4 算 例

如圖1所示的矩形截面懸臂深梁，長L=500 mm，寬b=100 mm，高h=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作用的時間為 10 s．

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

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

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

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

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

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

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

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

表1 有無作用力時對懸臂梁自由端撓度的影響

5 結(jié) 論

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

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

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

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

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

(5) 在有力的作用時，不考慮熱彈耦合效應(yīng)對結(jié)構(gòu)位移的影響將可能造成一定的計算誤差．因此，在高精度要求的結(jié)構(gòu)分析和設(shè)計中，考慮熱彈耦合效應(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ā)振動的有限元分析[J]. 力學(xué)學(xué)報, 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]. 北京: 國防工業(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] 李智勇, 劉錦陽, 洪嘉振. 作平面運動的二維平面板的熱耦合動力學(xué)問題[J]. 動力學(xué)與控制學(xué)報, 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)的熱-動力學(xué)耦合問題及其有限元分析[J]. 固體力學(xué)學(xué)報, 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] 閻彬, 陳建軍, 馬洪波. 隨機彈性桿在隨機瞬態(tài)溫度場下的熱響應(yīng)分析[J]. 電子科技大學(xué)學(xué)報, 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)動態(tài)特性分析[J]. 電子科技大學(xué)學(xué)報, 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)動力區(qū)間分析方法研究[J]. 應(yīng)用力學(xué)學(xué)報, 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] 易平. 對區(qū)間不確定性問題的可靠性度量的探討[J]. 計算力學(xué)學(xué)報, 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.