考慮熱－結(jié)構(gòu)耦合關(guān)系下空間結(jié)構(gòu)動(dòng)力響應(yīng)隨機(jī)性分析

閻 彬，陳建軍，馬洪波，靳紅玲

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

空間結(jié)構(gòu)的主要熱交換形式為熱輻射及結(jié)構(gòu)內(nèi)的導(dǎo)熱，而無(wú)對(duì)流，在其工作過(guò)程中受到太陽(yáng)照射的持續(xù)熱流。結(jié)構(gòu)在照射面和背陰面產(chǎn)生的截面溫差作用，將導(dǎo)致結(jié)構(gòu)的熱變形?？紤]結(jié)構(gòu)熱變形會(huì)改變光照入射角度，從而影響結(jié)構(gòu)表面的溫度分布，而溫度的分布差異又會(huì)使截面溫差產(chǎn)生影響，進(jìn)而影響結(jié)構(gòu)變形，從而產(chǎn)生了熱－結(jié)構(gòu)之間的耦合關(guān)系。另一方面，結(jié)構(gòu)溫度的變化使得結(jié)構(gòu)的各項(xiàng)參數(shù)都有不同程度的改變，體現(xiàn)出一定的隨機(jī)性，故有必要考慮參數(shù)隨機(jī)性的影響。綜合考慮上述因素的空間結(jié)構(gòu)的熱響應(yīng)情況，將使分析結(jié)果更為準(zhǔn)確，并可為空間結(jié)構(gòu)的熱疲勞可靠性預(yù)測(cè)奠定基礎(chǔ)。

目前航天器等空間結(jié)構(gòu)主要采用薄壁圓管構(gòu)件。文獻(xiàn)［1－2］給出了一種同時(shí)考慮截面平均溫度和截面溫差的薄壁圓管有限單元模型，可同時(shí)分析結(jié)構(gòu)熱軸力和熱彎矩，并應(yīng)用于哈勃望遠(yuǎn)鏡太陽(yáng)能帆板的算例中。文獻(xiàn)［3］在文獻(xiàn)［1］的基礎(chǔ)上考慮了熱變形與溫度的耦合關(guān)系，但基于在溫度變化范圍內(nèi)，結(jié)構(gòu)參數(shù)不變化的假定。實(shí)際中空間結(jié)構(gòu)由于溫變其參數(shù)必然發(fā)生改變，只是這些變化是微小的，但對(duì)于高精度的設(shè)備結(jié)構(gòu)，這種由于溫變而導(dǎo)致的結(jié)構(gòu)參數(shù)變化將是不容忽略的。文獻(xiàn)［4－5］分析了梁結(jié)構(gòu)在熱載荷作用下的動(dòng)力響應(yīng)問(wèn)題。文獻(xiàn)［7－9］先后對(duì)多種結(jié)構(gòu)的隨機(jī)響應(yīng)問(wèn)題進(jìn)行了研究，并獲得了隨機(jī)參數(shù)對(duì)響應(yīng)的影響。在目前諸多的研究成果中，空間結(jié)構(gòu)在熱載荷下的動(dòng)力響應(yīng)分析已有相應(yīng)文獻(xiàn)，但針對(duì)其隨機(jī)性分析的工作尚鮮見，故尋求一種簡(jiǎn)單有效的空間結(jié)構(gòu)在熱載荷下的隨機(jī)性分析方法是十分必要的。

本文基于以上研究成果，綜合考慮熱－結(jié)構(gòu)耦合效應(yīng)和參數(shù)的隨機(jī)性，對(duì)空間薄壁圓管的溫度場(chǎng)及動(dòng)力響應(yīng)進(jìn)行求解，利用基于靈敏度分析的矩法和隨機(jī)因子法推導(dǎo)出結(jié)構(gòu)溫度和動(dòng)力響應(yīng)的數(shù)字特征計(jì)算表達(dá)式，分析了熱－結(jié)構(gòu)耦合效應(yīng)對(duì)結(jié)構(gòu)動(dòng)力響應(yīng)的影響，獲得了各響應(yīng)量均值和均方差的時(shí)間歷程，并考察了各隨機(jī)參數(shù)對(duì)溫度、位移和應(yīng)力響應(yīng)分散性的影響，最后通過(guò)Monte-Carlo模擬法對(duì)計(jì)算結(jié)果進(jìn)行驗(yàn)證。

1 空間結(jié)構(gòu)熱－動(dòng)力學(xué)有限元模型的建立

1.1 考慮熱－結(jié)構(gòu)耦合關(guān)系的熱傳導(dǎo)有限元方程

以圖1和圖2所示的空間結(jié)構(gòu)最常見的薄壁圓管為分析對(duì)象，該結(jié)構(gòu)具有壁薄且軸向尺寸長(zhǎng)的特點(diǎn)。據(jù)此特點(diǎn)和空間環(huán)境的熱傳導(dǎo)機(jī)理，在建立熱傳導(dǎo)微分方程時(shí)基于三點(diǎn)假定:① 輻射換熱條件只考慮結(jié)構(gòu)外表面對(duì)太空的散熱，而忽略結(jié)構(gòu)內(nèi)壁間的換熱;②結(jié)構(gòu)截面平均溫度遠(yuǎn)大于截面擾動(dòng)溫度;③ 忽略沿壁厚方向的導(dǎo)熱，只考慮結(jié)構(gòu)周向?qū)?。考慮熱－結(jié)構(gòu)的耦合關(guān)系，即計(jì)及結(jié)構(gòu)變形對(duì)表面接收熱流量的影響。由圖1幾何關(guān)系可知，θz為結(jié)構(gòu)變形后表面法線與豎直方向的夾角，其大小等于入射點(diǎn)處的截面轉(zhuǎn)角θ=dw/(ldξ)，其中l(wèi)為管長(zhǎng)，w為入射點(diǎn)處橫向位移，ξ=x/l為軸向x方向的無(wú)量綱化坐標(biāo)。轉(zhuǎn)角變化會(huì)導(dǎo)致入射熱流量的改變，考慮圖1位置處dw/dξ為負(fù)值，故S=S0cos［β－dw/(ldξ)］表示結(jié)構(gòu)表面的凈輻射熱流量，S0為太陽(yáng)輻射熱流密度，β為光照入射角。據(jù)文獻(xiàn)［1－2］和以上假定，考慮由太陽(yáng)輻射的熱量輸入，結(jié)構(gòu)向外的輻射散熱，結(jié)構(gòu)在軸向和周向的導(dǎo)熱，以及熱－結(jié)構(gòu)的耦合關(guān)系，利用熱輻射的Stefan-Boltzmann定律建立熱平衡微分方程為:

式中:R，h分別為圓管的半徑和管壁厚，ρ為質(zhì)量密度，c為比熱容，ε為管外表面散射系數(shù)，T為管內(nèi)的瞬態(tài)溫度場(chǎng)，kφ、kξ分別為周向和軸向?qū)嵯禂?shù)，αs為管外表面吸收系數(shù)，φ為管截面的周向角度，σ為Stefan-Boltzmann常數(shù)，δ為跟蹤功能函數(shù)，反映了太陽(yáng)光照射角度與照射表面間的關(guān)系，它取決于β、θz角，即有:

為了便于考察由管截面溫差導(dǎo)致的熱致振動(dòng)，將管截面溫度T分解為兩部分之和:

其中:T0(ξ，t)為管截面內(nèi)平均溫度，T1(ξ，t)為管截面內(nèi)最大擾動(dòng)溫度，兩者都為軸向坐標(biāo)和時(shí)間的函數(shù)。將式(1)中的溫度非線性項(xiàng)T4在截面平均溫度處一階泰勒展開，即有:

將式(2)、(3)代入式(1)中，整理可得兩個(gè)獨(dú)立的分別由截面平均溫度和截面最大溫差為變量的微分方程，對(duì)此微分方程利用伽遼金方法沿軸向積分，可得到薄壁管的熱傳導(dǎo)有限元方程。圓管的有限元離散采用兩節(jié)點(diǎn)Lagrange單元，每個(gè)節(jié)點(diǎn)有兩個(gè)溫度自由度，即截面平均溫度和截面最大溫差，其有限元方程式為:

其中各矩陣和向量表達(dá)式如下［1－3］:

式中:C為熱容矩陣，KT為熱傳導(dǎo)的等效剛度矩陣，RT為與輻射相關(guān)的矩陣項(xiàng)，N為兩節(jié)點(diǎn)Lagrange插值函數(shù)，B=dN/dξ為應(yīng)變矩陣，m為單元數(shù)，Q0為變形前等效的熱載荷向量，Q1由于包含dw/dξ項(xiàng)，為結(jié)構(gòu)變形與溫度場(chǎng)的耦合項(xiàng)。諸式中上角標(biāo)為“0”者是以截面平均溫度建立的有限元列式中的相關(guān)矩陣及向量，而為“1”者則對(duì)應(yīng)于截面擾動(dòng)溫度部分。

對(duì)式(4)的求解分兩部進(jìn)行，由于R0T項(xiàng)的存在，使式(4a)成為非線性瞬態(tài)溫度場(chǎng)問(wèn)題。本文采用時(shí)間差分并結(jié)合牛頓內(nèi)迭代的數(shù)值解法。首先對(duì)式(4a)利用Galerkin格式的時(shí)間差分將其展開為［10］:

式中:Δt為時(shí)間步長(zhǎng)，當(dāng)已知t時(shí)刻平均溫度T0t，且結(jié)構(gòu)在t時(shí)刻的動(dòng)力響應(yīng)已得出，代入式(5)成為關(guān)于的非線性方程，利用牛頓內(nèi)迭代法即可近似求解

對(duì)式(4b)，由于平均溫度T0在t和t+Δt時(shí)刻都已知，故R1T即為已知，同樣按照式(5)的形式展開，得到以T1t+Δt為待求變量的線性方程，直接利用時(shí)間差分法即可求解。

1.2 溫度荷載下空間結(jié)構(gòu)動(dòng)力學(xué)有限元分析

將空間薄壁圓管結(jié)構(gòu)作為Euler-Bernoulli梁，并離散為m個(gè)單元、m+1個(gè)節(jié)點(diǎn)。對(duì)該梁考慮軸向位移的動(dòng)力學(xué)有限元方程為［6］:

式中，a=［u1w1θ1… um+1wm+1θm+1］T為結(jié)構(gòu)節(jié)點(diǎn)位移向量，每一節(jié)點(diǎn)具有:軸向位移u、橫向位移w和截面轉(zhuǎn)角θ三個(gè)自由度;單元軸向位移由兩節(jié)點(diǎn)Lagrange插值形函數(shù)N插值表示，單元橫向位移則利用Hermite插值函數(shù)Nw;M，K和D分別為結(jié)構(gòu)的質(zhì)量矩陣、剛度矩陣和阻尼矩陣，在此采用 Rayleigh阻尼［2，6］，即取 D= α2M+ β2K，其中 α2、β2為 Rayleigh 阻尼系數(shù)，Pf、PT分別為結(jié)構(gòu)節(jié)點(diǎn)力載荷向量和由溫度引起的熱載荷向量，M、K和Pf的表達(dá)式可由文獻(xiàn)［6］得知，PT的表達(dá)式如下:

式中:E為彈性模量，α為線膨脹系數(shù)，T0為結(jié)構(gòu)初始溫度，F(xiàn)T為軸向熱載荷，MT為熱彎矩。

由上式可看出，軸向熱載荷主要由截面平均溫度產(chǎn)生，而熱彎矩僅由截面擾動(dòng)溫度引起，兩者分別產(chǎn)生的軸向位移和橫向位移(包括截面轉(zhuǎn)角)可以分開求解。

式(4)和式(6)共同構(gòu)成了空間結(jié)構(gòu)的熱－動(dòng)力學(xué)耦合有限元方程，對(duì)其求解須采用時(shí)間積分法在每個(gè)時(shí)間步處將溫度場(chǎng)和動(dòng)力響應(yīng)交替迭代進(jìn)行，這是一種在時(shí)間域內(nèi)的近似數(shù)值計(jì)算方法，其中關(guān)于瞬態(tài)溫度場(chǎng)的計(jì)算見1.1節(jié)，而對(duì)動(dòng)力響應(yīng)的求解則利用了Newmark積分法，具體求解步驟為:

(1)給定初始值，即t0時(shí)刻溫度分布T0(t0)、T1(t0)，位移a(t0)，速度(t0)和加速度(t0)均為已知，并確定時(shí)間步長(zhǎng)Δt;

(2)以Δt為時(shí)間步長(zhǎng)，將T0(t0)、a(t0)代入 式(4a)，由于非線性項(xiàng)R0T，利用時(shí)間差分法即式(5)格式和牛頓迭代法計(jì)算得到T0(t0+Δt);

(3)將T1(t0)、T0(t0)、T0(t0+Δt)、a(t0)代入式(4a)，式(4b)為線性方程，直接利用時(shí)間差分求得T1(t0+Δt);

(4)利用已求得的T0(t0+Δt)和T1(t0+Δt)代入式(7)，得到溫度載荷PT，結(jié)合a(t0)(t0)和(t0)對(duì)(6)式利用Newmark法進(jìn)行動(dòng)力響應(yīng)計(jì)算，獲得a(t0+Δt)、(t0+Δt)和·a·(t0+ Δt)。至此 t0+ Δt時(shí)刻的溫度場(chǎng)及動(dòng)力響應(yīng)均已求出，繼續(xù)重復(fù)(2)～(4)步即可得到各個(gè)時(shí)間步的溫度場(chǎng)及動(dòng)力響應(yīng);

2 結(jié)構(gòu)響應(yīng)隨機(jī)性分析

在建立了上述確定性分析模型后，現(xiàn)將模型中的物性參數(shù):質(zhì)量密度ρ，彈性模量E，比熱容c，熱傳導(dǎo)系數(shù)k，太陽(yáng)光入射角度β，結(jié)構(gòu)外表面散射系數(shù)ε，熱吸收系數(shù)αs，熱膨脹系數(shù)α等均視為服從正態(tài)分布的隨機(jī)參數(shù)，則該模型成為含有隨機(jī)參數(shù)的不確定性模型，需分別對(duì)溫度場(chǎng)和動(dòng)力響應(yīng)的隨機(jī)性進(jìn)行分析。鑒于各自方程的不同特點(diǎn)和已有的研究成果，對(duì)兩者隨機(jī)性分析采用了不同的方法。

2.1 溫度場(chǎng)隨機(jī)性分析

將上述所有隨機(jī)參數(shù)代入式(4)中，則式(4)成為含隨機(jī)參數(shù)的有限元方程，對(duì)于其隨機(jī)性分析，是以1.1節(jié)溫度場(chǎng)求解方法為基礎(chǔ)，利用求解隨機(jī)變量函數(shù)數(shù)字特征的矩法進(jìn)行分析，首先計(jì)算出溫度對(duì)各物性參數(shù)的靈敏度，再根據(jù)矩法，即可獲得溫度場(chǎng)的均值和方差。

同前述瞬態(tài)溫度場(chǎng)的求解順序，先對(duì)平均溫度的靈敏度進(jìn)行分析。將結(jié)構(gòu)物性參數(shù)記為ψ，從式(5)出發(fā)，兩邊對(duì)ψ求導(dǎo)，整理得:

式(8)即為平均溫度靈敏度求解的遞推迭代格式，其中:(·)'表示 ?(·)/?ψ，對(duì)變量 ψ 的靈敏度;為將中平均溫度視為常數(shù)后，對(duì)ψ求偏導(dǎo)所得部分。式中各個(gè)矩陣和向量對(duì)ψ的偏導(dǎo)數(shù)由以下差分格式求得:

其中:X代表式中各相關(guān)矩陣和向量。

(1)對(duì)于各時(shí)不變矩陣及向量如K0T，C，Q0，可由式(9)計(jì)算出對(duì)ψ的靈敏度，由于它們與時(shí)間無(wú)關(guān)，故為常值，在迭代過(guò)程中只需計(jì)算一次。由于Q1和 R0T包含 dw/dξ和T0，故為時(shí)變向量和矩陣，其對(duì)ψ的偏導(dǎo)仍用差分法計(jì)算，只是在每個(gè)時(shí)間步都需重新計(jì)算一次。

(2)根據(jù)前述求解熱－動(dòng)力學(xué)耦合有限元方程的方法，已知t時(shí)刻的溫度分布和動(dòng)力響應(yīng)，計(jì)算一步得到，代 入中 則為 已 知，仍可利用差分法計(jì)算出，此時(shí)為僅包含待求量'的矩陣，將前兩步所得數(shù)據(jù)代入式(8)則可計(jì)算出

在得到溫度場(chǎng)T0和T1對(duì)各隨機(jī)變量的靈敏度后，則可利用矩法［8］求得溫度場(chǎng)的數(shù)字特征，推導(dǎo)得其均值和方差的計(jì)算表達(dá)式為:

2.2 動(dòng)力響應(yīng)隨機(jī)性分析

對(duì)于動(dòng)力響應(yīng)隨機(jī)性分析，與溫度場(chǎng)隨機(jī)性分析相似，同樣從方程的求解格式出發(fā)，利用隨機(jī)因子法推導(dǎo)得到響應(yīng)的數(shù)字特征表達(dá)式。首先給出利用Newmark法求解動(dòng)力響應(yīng)的計(jì)算公式，再將各隨機(jī)參數(shù)以隨機(jī)因子的表示形式代入，利用求解隨機(jī)變量函數(shù)數(shù)字特征的代數(shù)綜合法，導(dǎo)出位移響應(yīng)的均值和方差。

考慮到空間環(huán)境下熱載荷為主要載荷，故在本文中僅考察熱載荷下結(jié)構(gòu)的動(dòng)力響應(yīng)，重點(diǎn)分析在熱彎矩作用下結(jié)構(gòu)的橫向位移及截面轉(zhuǎn)角，此時(shí)結(jié)構(gòu)節(jié)點(diǎn)位移向量為 w=(0，w1，θ1，…，0，wm+1，θm+1)T，相應(yīng)的式(6)中各矩陣向量的形式都有所改變，只是將各矩陣向量中對(duì)應(yīng)軸向位移部分都取為0即可。利用Newmark法對(duì)式(6)求解，即將其轉(zhuǎn)化為如下擬靜力迭代方程求解:

式中:α1和β1為精度和穩(wěn)定性要求的可調(diào)參數(shù)，取α1=0.25，β1=0.5［10］時(shí)格式為無(wú)條件穩(wěn)定。

利用隨機(jī)因子法［7－8］，將所有隨機(jī)參數(shù)均以隨機(jī)因子形式表達(dá)，即令:E=，其余隨機(jī)參數(shù)表示類似，其中:(·)*項(xiàng)為隨機(jī)因子，其均值為1，變異系數(shù)為 ν(·)=均方差/均值;)項(xiàng)為隨機(jī)參數(shù)的均值即名義值。

圖3 懸臂薄壁圓管梁Fig.3 Thin-walled cantilever beam of circular cross-section

利用代數(shù)綜合法對(duì)式(12)兩邊分別求均值和方差，可得位移響應(yīng)在 t+Δt時(shí)刻的均值 μwt+Δt和方差分別為:

在得到位移響應(yīng)后，不考慮由軸向位移導(dǎo)致的軸向應(yīng)變，根據(jù)單元節(jié)點(diǎn)位移與單元應(yīng)力之間的關(guān)系，可得結(jié)構(gòu)軸向正應(yīng)力為:

將式(15)中的隨機(jī)參數(shù)以隨機(jī)因子形式表示，則單元軸向正應(yīng)力Seeξ可表示為:

對(duì)式(16)利用代數(shù)綜合法，可得單元應(yīng)力的均值和方差為:

式中，μwe，σ2we分別為第e個(gè)單元的節(jié)點(diǎn)位移均值和方差。

3 算例

以哈勃太空望遠(yuǎn)鏡上的太陽(yáng)能電池帆板中的懸臂薄壁圓管梁為例，見圖3。總長(zhǎng)5.91 m，將其劃分為50個(gè)單元51個(gè)節(jié)點(diǎn)，則l=0.118 2 m。該結(jié)構(gòu)在t=0 s時(shí)刻受到來(lái)自太陽(yáng)的熱流S0=1 350 W/m2，假設(shè)結(jié)構(gòu)表面初始溫度為290 K，σ =5.67 ×10－8W·m－2·K－4，給定Rayleigh阻尼系數(shù)α2=β2=0.002，結(jié)構(gòu)其它參數(shù)的名義值(均值)見表1［1］。由結(jié)構(gòu)固有頻率的計(jì)算結(jié)果，選取Δt為0.1 s，分析結(jié)構(gòu)的瞬態(tài)隨機(jī)溫度場(chǎng)分布及動(dòng)力響應(yīng)，并考察各隨機(jī)參數(shù)對(duì)響應(yīng)分散性的影響。

表1 結(jié)構(gòu)各參數(shù)的名義值Tab.1 The nominal value of material properties

圖4～圖9為當(dāng)隨機(jī)參數(shù)變異系數(shù)均為0.01時(shí)，結(jié)構(gòu)的溫度和動(dòng)力響應(yīng)的均值及均方差隨時(shí)間的歷程圖。

圖4和圖5為自由端截面平均溫度和截面最大溫差均值的時(shí)間歷程圖，從中可見，截面溫差的變化率大于平均溫度，平均溫度在800 s后趨于穩(wěn)定值400 K，而截面溫差在100 s后基本穩(wěn)定達(dá)到10 K。由于溫度場(chǎng)的穩(wěn)態(tài)值由初始太陽(yáng)光入射角、熱流密度、散射和吸收系數(shù)決定［4－5］，故轉(zhuǎn)角變化導(dǎo)致熱流突變，進(jìn)而產(chǎn)生小幅度變化的截面溫差，突變的截面溫差又會(huì)產(chǎn)生變化的熱載荷，從而導(dǎo)致了不穩(wěn)定的顫振。為了看清截面溫差突變效果，給出圖5的局部放大圖，可見有突變現(xiàn)象存在，但是幅度不大，這是因?yàn)榻Y(jié)構(gòu)振動(dòng)的幅度在初始時(shí)段較小，截面轉(zhuǎn)角亦很小，故耦合作用不明顯造成的。

圖6為自由端橫向位移均值的時(shí)間歷程，可見振動(dòng)有增強(qiáng)且伴有小幅度的突變，表現(xiàn)為顫振，顫振的幅度由阻尼和光照入射角決定，其中以阻尼為最主要因素，算例選取較小的阻尼系數(shù)為了使耦合顫振效果更為明顯。另見自由端的準(zhǔn)靜態(tài)位移由0迅速增大并趨于穩(wěn)定，這是由截面溫差的變化規(guī)律決定的。

圖7為最大應(yīng)力固定端處應(yīng)力均值隨時(shí)間變化歷程，可見其變化規(guī)律與圖4平均溫度相似，隨著平均溫度的升高不斷增大，且伴隨著由熱振動(dòng)造成的應(yīng)力擾動(dòng)，擾動(dòng)隨著顫振而增強(qiáng)，這有可能造成結(jié)構(gòu)的破壞或功能失效，且長(zhǎng)期作用亦可能引發(fā)結(jié)構(gòu)的熱疲勞破壞。

由圖8和圖9響應(yīng)的均方差分析結(jié)果可見，截面最大溫差分散性與隨機(jī)參數(shù)的分散性基本一致，位移響應(yīng)分散性隨著時(shí)間的增加逐漸增大且顯出振動(dòng)特性，這是由熱－結(jié)構(gòu)耦合作用產(chǎn)生的熱顫振引起的。

為了考察各隨機(jī)參數(shù)分散性對(duì)響應(yīng)結(jié)果的影響，表2列出了所有隨機(jī)參數(shù)變異系數(shù)取不同值時(shí)自由端橫向位移及固定端應(yīng)力響應(yīng)的均值和均方差在100 s時(shí)刻的結(jié)果。由于E和ρ呈正相關(guān)，故兩者的變異系數(shù)相同。

表2 動(dòng)力響應(yīng)在100 s時(shí)的均值和均方差Tab.2 The mean value and standard deviation of dynamic responses at the one hundred seconds

從表2中的計(jì)算結(jié)果可見:參數(shù)的隨機(jī)性對(duì)響應(yīng)分散性(均方差)影響明顯，當(dāng)所有隨機(jī)參數(shù)變異系數(shù)均增大10倍時(shí)，響應(yīng)亦以近10倍關(guān)系增大。針對(duì)各參數(shù)隨機(jī)性對(duì)響應(yīng)均值和均方差影響程度評(píng)估上，E和ρ的隨機(jī)性對(duì)位移和應(yīng)力響應(yīng)均值均有一定影響。E，ρ和α的隨機(jī)性對(duì)位移和應(yīng)力響應(yīng)分散性影響很大。比較其它隨機(jī)參數(shù)，αs，k對(duì)位移和應(yīng)力響應(yīng)分散性影響較大，ε和c則影響很小。為了驗(yàn)證本文結(jié)果的正確性，這里利用Monte-Carlo方法檢驗(yàn)。由于MC方法需要大量的抽樣模擬，對(duì)于本文瞬態(tài)隨機(jī)性問(wèn)題，計(jì)算量巨大，為此，本文對(duì)所有隨機(jī)模型均給出了3 000次MC法的數(shù)值模擬結(jié)果，對(duì)比可見，本文方法的結(jié)果與MC法的結(jié)果相差不大，但本文方法可直接一次性的得到各時(shí)間步處溫度場(chǎng)及動(dòng)力響應(yīng)的數(shù)字特征，具有較高的計(jì)算效率。

4 結(jié)論

本文對(duì)綜合考慮熱－結(jié)構(gòu)耦合關(guān)系和參數(shù)具有隨機(jī)性的空間薄壁圓管梁在瞬態(tài)突加熱流作用下的動(dòng)力響應(yīng)問(wèn)題進(jìn)行了研究，理論分析與計(jì)算結(jié)果表明:

(1)熱－結(jié)構(gòu)耦合將有可能導(dǎo)致結(jié)構(gòu)不穩(wěn)定的熱顫振，且使動(dòng)力響應(yīng)的幅值增大。

(2)結(jié)構(gòu)參數(shù)的隨機(jī)性對(duì)響應(yīng)分散性影響不可忽略，其中彈性模量E，質(zhì)量密度ρ，線膨脹系數(shù)α，熱吸收系數(shù)αs和熱傳導(dǎo)系數(shù)k對(duì)響應(yīng)分散性影響較大，而散射系數(shù)ε和比熱容c則相對(duì)較小。

(3)通過(guò)與Monte-Carlo數(shù)值模擬法的結(jié)果對(duì)比表明，本文的求解方法可行且有效。

［1］ Ding Y，Xue M，Kim J K.Thermo-structural analysis of space structures using fourier tube elements［J］.Comput Mech，2005，36:289－297.

［2］ Xue M D，Duan J，Xiang Z H.Thermally-induced bendingtorsion coupling vibration of large scale space structures［J］.Comput Mech，2007，40:707－723.

［3］ 程樂(lè)錦，薛明德.大型空間結(jié)構(gòu)熱－動(dòng)力學(xué)耦合有限元分析［J］.清華大學(xué)學(xué)報(bào)(自然科學(xué)版)，2004，44(5):681－684.CHENG Le-jin，XUE Ming-de.Coupled thermal-dynamic FEM analysis of large scale structures［J］.J Tsinghua Univ(Sci＆Tech)，2004，44(5):681－684.

［4］ Ribeiro P，Manoach E.The effect of temperature on the large amplitude vibrations of curved beams［J］.Journal of Sound and Vibration，2005，285:1093－1107.

［5］ Hassan M A，Kanno I，Kotera H.Compound two-dimensional thermo-elastic and thermodynamic analysis for c-axis-oriented epitaxial lead titanate thin films［J］.Vacuum，2006，81:461－467.

［6］ 王勖成，邵 敏.有限單元法基本原理和數(shù)值方法［M］.北京:清華大學(xué)出版社，1997.

［7］ 王小兵，陳建軍，梁震濤，等.熱機(jī)電耦合隨機(jī)參數(shù)智能板結(jié)構(gòu)的數(shù)值求解［J］.西安電子科技大學(xué)學(xué)報(bào)，2008，35(4):592－599.WANG Xiao-bing，CHEN Jian-jun，LIANG Zhen-tao，et al.Randomicity analysis of the piezothermoelasticity intelligent thin plate by appling the numerical method［J］.Journal of Xidian University，2008，35(4):592－599.

［8］ 魏永祥，陳建軍，王敏娟.隨機(jī)參數(shù)齒輪－轉(zhuǎn)子系統(tǒng)的扭轉(zhuǎn)振動(dòng)分析［J］.航空動(dòng)力學(xué)報(bào)，2010，25(11):2637－2642.WEIYong-xiang， CHEN Jian-jun， WANG Min-juan.Dynamic response of torsional vibration of gear-rotor with random parameters［J］.Journal of Aerospace Power，2010，25(11):2637－2642.

［9］ Salam M R A E，Alla A M A，Hosham H A.A numerical solution of magneto-thermoelastic problem in nonhomogeneous isotropic cylinder［J］.Applied Mathematical Modelling，2007，31:1662－1670.

［10］ 吳永禮.計(jì)算固體力學(xué)方法［M］.北京:科學(xué)出版社，2005.