有限體積法的彈性結構動力學隨機分析

陳衛東，陳浩，于艷春

(哈爾濱工程大學航天與建筑工程學院，黑龍江哈爾濱150001)

艦船、潛艇等武器平臺在戰斗中不可避免地遭到離其數米甚至百米位置的爆炸沖擊，這種爆炸通常不會使其產生嚴重的破損，但可能引起劇烈的振動，從而破壞武器上的機械設備和電子儀器，使其失去戰斗力［1］.因此，研究遠場爆炸載荷作用下結構的動力學隨機響應問題具有重要的現實意義.

結構在隨機載荷作用下的隨機響應一直受到研究者的關注并不斷取得進展.目前，求解隨機響應問題的方法有Monte-Carlo法、隨機有限元法和虛擬激勵法等.Monte-Carlo法具有普遍的適用性，但其巨大的計算量限制了其在大型工程中的應用;隨機有限元法雖然在靜力學隨機分析中得到廣泛應用，但在動力學分析中發展緩慢［2］;而基于虛擬激勵法的動力學隨機響應分析主要應用在抗震、抗風等工程中［3-4］，對于爆炸等沖擊載荷作用下隨機動響應分析涉及較少.所以，文章將有限體積法與泰勒展開相結合，在考慮爆炸載荷的隨機性時，求解了遠場水下爆炸載荷作用下彈性結構動響應的統計特性，并對結構動響應的變異系數做了相關討論.由于采用顯式求解算法，同步求解結構動響應及其對基本隨機變量的偏導，不需求解大型線性方程組，因而在計算爆炸沖擊這一類問題時具有較高的效率.

1 有限體積法

有限體積法由于其思想簡潔，用其導出的離散方程，物理上表示的是控制體積的通量平衡，方程中的各項有明確的物理意義，所以近來有不少學者將其應用于固體力學［5-6］.

圖1 控制體與其輔助網格示意Fig.1 Relation between control volume and auxiliary mesh

圖1是非結構化網格示意圖.對于圖1(a)中平面網格:虛線網格為輔助網格(三角形)，連接輔助網格的重心及其各邊中點(實心圓)就構成了控制體的邊界，實線圍成的閉合多邊形即為控制體.顯然，虛心圓即為輔助網格的節點，又為控制體的中心.對于三維問題，輔助網格為四面體，圖1(b)給出了輔助網格與其內部控制體邊界面的關系:每個四面體都被其內部控制體邊界面分成體積相等的4部分.

在不考慮能量變化時，有限體積法的基本控制方程為動量守恒方程:

若將加速度、速度等變量定義在控制體中心，且輔助網格中應力、應變為常量，那么將式(1)在空間網格上離散，最終可得到求解三維彈性動力學問題的方程為［7］

式中:m為控制體中心連接的四面體的個數.引入初始條件后，可采用文獻［8］中的顯式中心差分法求解式(2):

常數an1、an2和an3與控制體中心在四面體中的節點編號有關，若其在四面體節點編號中為1，則

對其他的編號有同樣形式的表達式，需要注意的是當前需求系數的節點編號與行列式中節點編號應符合右手定則:行列式中節點按逆時針排列，方向指向當前系數的節點編號.

2 結構動響應的統計特性

在不考慮結構尺寸的隨機性時，設有一組基本隨機變量X=(x1x2… xn)T，將式(2)對X求偏導得

式中:i、j=1，2，…n.若 t時刻各物理量已知，則可利用顯式算法得到t+Δt時刻的物理量對基本隨機變量的偏導:

式中:D和B分別為彈性矩陣和幾何矩陣，具體可參見文獻［9］中有關常應力四面體元的表達式.

若要求解等效應力的統計特性，則還需要知道等效應力對基本隨機變量的偏導.若skl為偏應力的分量，那么

式中:Cov(xi，xj)為基本隨機變量 xi和 xj的協方差，對于位移和應力的統計特征也是同樣的方法.

3 固支板的隨機動響應分析

為了驗證前述理論的可行性，利用FORTRAN語言編寫了計算程序，并計算了一個固支板遭受水下遠場爆炸載荷的沖擊問題.

3.1 爆炸載荷模型

水中沖擊波過后壓力隨時間變化關系呈指數衰減規律:

式中:Pm為沖擊波超壓峰值，θ為衰減時間常數.計算時，采用文獻［10］中的經驗公式:

式中:w為TNT當量，kg;R為結構表面到爆源的距離，m.

對于遠場爆炸，近似認為沖擊波是平面壓力波，并垂直作用于結構表面［11］.若忽略結構變形的影響，將沖擊波看作是作用在剛體上，則作用于結構表面的壓力p(t)與沖擊波P(t)之間的關系為［12］

當把w和R作為基本隨機變量時，該載荷模型就是一個隨機爆炸載荷.

3.2 結構模型

計算模型見圖2.材料密度ρ=7 800 kg/m3，彈性模量E=210 GPa，泊松比為0.3.將炸藥的裝藥量w、板面到爆源的距離R看作相互獨立的基本隨機變量.由于結構在水中可能遭受各種炸藥的爆炸沖擊，為簡化問題，將不同炸藥都看成TNT，并對TNT的裝藥量和爆距取較大的變異系數來代表不同炸藥性能的變異.各基本隨機變量的取值見表1.

圖2 固支板模型Fig.2 Model of a clamped plane

表1 基本隨機變量統計特征Table 1 Statistical characterization of basic random variables

3.3 結構的隨機動響應分析

利用前述的隨機有限體積法(SFVM)計算了固支板在遭受爆炸載荷時動響應的統計特性:圖3是板中心(0.04，0.01，0.04)處 y 向位移的均值和方差，圖4 是點(0.04，0.01，0.00)處等效應力的均值和方差.

為驗證計算結果的正確性，基于有限體積法進行了1 000組Monte-Carlo模擬.圖3、4結果對比表明:SFVM與Monte-Carlo結果符合較好，該方法能準確的計算彈性結構動響應的統計特征.

圖5(a)是板中心(0.04，0.01，0.04)處 y 向位移和點(0.02，0.01，0.02)處各向位移的變異系數絕對值;圖5(b)是板中心點和點(0.04，0.01，0.00)處等效應力的變異系數.從圖5可知，無論是節點位移變異系數的絕對值還是等效應力的變異系數，其隨時間始終近似等于0.289 1.雖然當響應量的方差接近零時，變異系數有很大的波動，但這種波動只會使響應量的均值在接近于零時產生影響，并且當均值在零值附近時，其數量級與均值極值的數量級相差較大，所以可以忽略零值附近變異系數的波動，近似認為響應量的變異系數是不變的，因而可以推論:在彈性體內，當動載荷的隨機性確定時，結構內各點動響應的變異系數絕對值是相同的，且不隨時間發生變化.

圖3 板中心的位移統計特征Fig.3 Statistical characterization of central of the plane

圖4 等效應力的統計特征Fig.4 Statistical characterization of equivalent stress

圖5 響應量的變異系數Fig.5 Coefficient of variation for responses

為了證明以上推論的正確性，根據響應量的方差和變異系數 0.289 1，預測了點(0.04，0.01，0.02)處 y 向位移和點(0.02，0.01，0.00)處等效應力的均值，并與計算結果做了對比(圖6)，結果表明該推論是可行的.利用該推論，只需計算響應量對基本隨機變量的一階偏導，因而提高了計算效率，為SFVM應用于大型工程提供了基礎.

圖6 響應量的均值預測Fig.6 Predicted mean value of responses

4 結束語

將有限體積法與泰勒展開結合，能夠對任意隨機載荷作用下彈性結構的動響應進行隨機分析.而且在相同計算模型下，結構內各點動響應的變異系數絕對值近似相同，且不隨時間發生變化.這樣只要求出初始一段時間內響應量的變異系數和方差，就可以對響應量的均值進行預測，這使SFVM在具有較高精度的同時兼顧有較高的效率，為大型工程的動力學隨機分析提供了一種途徑.此外，由于對響應變異系數的討論只限于文章中所給形式的隨機爆炸載荷，所以文章中的推論不一定適用于其他形式的載荷，作者會在后續工作中繼續討論載荷形式和約束條件等對響應變異系數的影響.

［1］姚熊亮，郭君，許維軍.船舶結構遠場爆炸沖擊動響應的數值試驗方法［J］.中國造船，2006，47(2):24-34.

YAO Xiongliang，GUO Jun，XU Weijun.Far field numerical experimental method on the explosion impact dynamic responses of ship［J］.Ship Building of China，2006，47(2):24-34.

［2］安偉光，蔡蔭林，陳衛東.隨機結構系統可靠性分析與優化設計［M］.哈爾濱:哈爾濱工程大學出版社，2005:65-93.

AN Weiguang，CAI Yinlin，CHEN Weidong.Reliability analysis and optimal design of stochastic structural system［M］.Harbin:Harbin Engineering University Press，2005:65-93.

［3］LIN J H，ZHANG W S，WILLIAMS F W.Pseudo-excitation algorithm for non-stationary random seismic responses［J］.Eng Struct，1994，16(4):270-276.

［4］慕文品.受演變隨機激勵結構響應的擴展精細積分方法［J］.振動與沖擊，2009，28(7):131-134.

MU Wenpin.An extended precise integration method for response of a structure subjected to evolutionary random exciation［J］.Journal of Vibration and Shock，2009，28(7):131-134.

［5］XIA Guohua，LIN Chinglong.An unstructured finite volume approach forstructuraldynamicsin responsetofluid motions［J］.Computers ＆ Structures，2008，86(7/8):684-701.

［6］LV X，ZHAO Y，HUANG X Y，XIA G H，SU X H.A matrix-free implicit unstructured multigrid finite volume method for simulating structural dynamics and fluid-structure interaction［J］.Journal of Computational Physics，2007，225:120-144.

［7］CHEN Weidong，CHEN Hao，ZHANG Wenping，et al.A finite volume method for 3-D elastodynamics［C］//Proceedings of the Third International Conference on Modeling and Simulation，VOL2-modelling and simulation in Engineering.Wuxi，China，2010:88-91.

［8］張雄，王天舒.計算動力學［M］.北京:清華大學出版社，2007:266-268.

［9］趙均海，汪夢甫.彈性力學及有限元［M］.2版.武漢:武漢理工大學出版社，2008:167-169.

［10］惲壽榕，趙衡陽.爆炸力學［M］.北京:國防工業出版社，2005:233-240.

［11］何建，肖玉鳳，陳振勇，等.空爆載荷作用下固支矩形鋼板的塑性極限變形［J］.哈爾濱工業大學學報，2007，39(2):310-313.

HE Jian，XIAO Yufeng，CHEN Zhenyong，et al.Plastic limited deformation analysis of the clamped rectangular steel plate subjected to air non-contact explosions［J］.Journal of Harbin Institute of Technology，2007，39(2):310-313.

［12］LIANG C C，TAI Y S.Shock responses of a surface ship subjected to noncontact underwater explosions［J］.Ocean Engineering，2006，33:748-772.