水下圓柱殼沖擊響應分析的有限元模型簡化方法

孫占忠 胡剛義 張新宇 黃國兵

中國艦船研究設計中心，湖北 武漢430064

水下圓柱殼沖擊響應分析的有限元模型簡化方法

孫占忠 胡剛義 張新宇 黃國兵

中國艦船研究設計中心，湖北 武漢430064

基于顯式非線性動力學分析程序ABAQUS，對水下圓柱殼受到水下爆炸沖擊作用下的非線性動態響應進行研究。詳細計算不同簡化模型受到水下爆炸沖擊波作用下的響應并與完整模型的結果進行對比。研究結果表明，用圓柱殼結構的部分模型代替完整模型，可以得到較為滿意的計算結果。

模型；簡化；水下爆炸；沖擊響應

1 引言

艦船的抗沖擊研究是振動領域比較熱門的一個課題，但可應用于工程實際的研究成果還較少。目前，研究的熱點集中在船體結構的沖擊響應、艦船內部設備的沖擊環境及對其進行有效的沖擊隔離研究等方面，在各大高校和研究院所都有所研究。數值模擬分析在艦船的抗沖擊研究中具有重要地位。為了取得良好的數值模擬結果，需要建立較準確的有限元模型。

通常對船體采用三種有限元模型：

1）全船船體梁模型；

2）全船三維模型；

3）船舶三維模型和梁模型相結合的整船混合有限元模型。

全船船體梁模型是將整船簡化為沿一條直線的梁單元，其優點是劃分的單元和節點數少，計算速度快，模型前處理工作量不大；缺點是不能反映內部結構的真實情況，無法分析船舶內部的薄弱環節。全船三維模型將船體離散成許多殼單元和梁單元的組集，其優點是能相對真實地反映船舶內部的情況，缺點是前處理工作量大，計算時間長。

文獻［1，2］對船體尾部混合有限元模型的建立提出了有效的方法，并在實際工程應用中證實，說明這些方法對船體的振動分析有較好的適用性。然而沖擊與振動不同，這些方法是否適用于沖擊響應分析有待于進一步研究。本文對圓柱殼采用三維結構與一維梁混合有限元模型和全部三維模型進行對比分析，探討混合模型對于抗沖擊模擬的可行性和準確性。本文主要運用ABAQUS等有限元程序，對圓柱殼體、外部水體進行建模，并模擬爆炸沖擊波，同時考慮流固耦合作用，計算沖擊響應。

2 有限元建模

2.1 三維圓筒整體模型

以鋼質圓筒為研究對象，其原始數據為：長L＝10 m；直徑R＝1 m；筒體厚度0.028 m；兩端蓋板厚度為0.05 m；E＝2.1×1011Pa；密度ρ＝7 850 kg／m3；泊松比ε＝0.3。用四邊形殼單元模擬。殼體網格劃分時，考慮到單元尺寸和時間步長對計算精度的影響，殼體網格尺寸取沖擊波長［3］（λ＝θc）的1／6。

殼體外部水體需要包圍整個圓筒，考慮到流固耦合的影響，流體模型范圍取約為殼體半徑的6倍時，附加質量與無限水域相比滿足精度要求，即兩端半球體半徑為3 m，中部環形水體長度與圓筒長度相同。為了達到更好的精度，采用六面體聲學單元模擬流體，對水體進行網格劃分，由于水體內部表面與圓筒外表面存在流固耦合作用，水體內表面單元不能太大，否則引起殼體表面應力應變極不均勻，其值取決于沖擊波長的大小。本模型取水體內部表面單元大小為λ／6，水體外表面單元大小為λ／3。

將圓柱殼體與水體模型組裝在一起，并定義結構模型外表面與水體內表面間的流固耦合面（ABAQUS中并不要求兩個耦合面間的節點重合），這樣就完成了整個系統建模過程［4］，系統模型如圖1所示，為表示方便隱去上部流體單元。

圖1 圓柱殼體及水體耦合模型

本文所用TNT炸藥重量W＝225 kg，爆炸點與入射點間距為R＝13.6 m，則入射點的壓力—時間函數如圖2。

2.2 兩種簡化方案模型

采用以下兩種方案。

方案1：取圓筒的前后1／3部分用梁模型簡化與圓筒三維模型聯結，聯結方式為梁節點和筒體圓周上的單元連接，聯結處采用運動耦合約束，且假設該圓周上單元的剛度很大（采用MPC剛性連接），從而保證聯結剖面處的變形一致，見圖3（a）。

方案2：采用局部圓筒結構，取總長度的1／3，見圖3（b），即在方案1的基礎上，拋棄前后的梁模型。

圖2 沖擊波壓力—時間函數

圖3 兩種方案有限元模型

2.3ABAQUS軟件的特殊梁模型

ABAQUS軟件提供一種劃分了網格的橫剖面用于程序的梁分析并能給梁單元施加載荷。對于某些結構的響應和梁的響應相似，然而這些梁剖面的幾何或這個剖面是由多種材料組成的原因，不能使用預先存在的梁單元庫模擬。在這些方案中，劃分網格的橫剖面來模擬梁剖面并產生梁剖面屬性適用于Timoshenko梁分析。產生的梁剖面屬性包括軸向、彎曲、扭轉和橫向剪切剛度、質量、轉動慣量和阻尼屬性，還包括這個剖面的質心和剪切中心。產生典型梁剖面屬性包含下列數據行，以本文使用的梁剖面為例：

需要劃分網格截面結構的典型應用是船體的振蕩分析，船體是由多種部件和多種材料組成的建筑。其他應用包括槳葉和機翼的翼型，分層的復合工字梁等。在ABAQUS中使用下列命令將簡化的梁單元和三維立體結構進行運動耦合連接，并提供附連水質量。

3 基本理論

在目前的數值分析中，有限元法和邊界元法均被采用。有限元法用來模擬殼體結構，邊界元法用來模擬其周圍的水介質。在結構分析中，采用VEC／DYNA3D程序［5］，使用不同的殼體單元類型，其中包括Hughes和Liu以及Belytschko［6］等人采用的殼單元和厚殼單元，幾乎得出了相同的結果。本文使用的鋼材的彈塑性本構方程未考慮應變強化假設。

水作為聲介質，用于距結構一定距離處爆炸物的傳播媒介。殼結構和聲學水介質之間的相互作用是通過雙重漸進近似法（DAA）［7］的近似技術來模擬的。DAA技術使用邊界元法來減少用于模擬聲學介質所需要的單元數量。在計算條件較差時，DAA法得到了廣泛的應用。但隨著計算機技術的發展和精確的無反射邊界條件的實現，采用經典的有限元法精確逼近求解流固耦合問題又重新具有了優勢。

下面來分析流固耦合作用。結構的運動方程可以用如下的有限元法給出：

式中，［Ms］、［Cs］和［Ks］分別表示結構的質量、阻尼和剛度矩陣；｛F｝表示力矢量；｛x｝是節點位移矢量，并具有二階時間導數。一部分力取決于流固耦合作用，一部分力取決于其他作用在結構上的外力。取決于流固耦合作用的力只發生在結構的濕表面上。因此，這個力矢量可以寫為：

式中，［G］是流固耦合傳遞矩陣；［Af］是與流體網格有關的對角線面積矩陣；Pi和Ps分別是入射和散射壓力，這里所說的散射壓力是波浪理論中繞射壓力和輻射壓力的總和；｛Fd｝是作用于干結構上的力矢量。在方程（2）中，當給定一個爆炸載荷時，散射壓力是未知的，入射壓力是已知的。

為了得到散射壓力值，下面方程是通過使用DAA理論和邊界元技術一起得到的：

式中，［Mf］是流體質量矩陣；ρ和c是水密度和水中聲速；｛u˙s｝是結構表面法線方向的散射波流體質點速度矢量。剩下的變量如上所述。流固耦合作用通過方程（2）聯立方程（1）和（3）。散射波流體質點速度與結構速度有關，如方程（4）所示：

式中，T是轉置矩陣的標志；i和s分別代表入射波和散射波。

方程（1）～（4）是完整的流固耦合表達式。方程（3）是精確的。在數值求解過程中引入一個增量，用交互求解的方法使方程求解時具有更好的穩定性［8］。基于以上數值方法理論，國外的科研機構和公司開發了一些通用計算軟件。本文的計算采用基于聲固耦合算法的ABAQUS軟件。

4 結果對比分析

4.1 方案1結果分析

三維長圓筒整體模型與用特殊梁簡化的混合模型進行比較分析。選取圓筒中心截面處4個典型節點加速度、位移曲線進行對比，如圖4、圖5所示。并對典型節點的沖擊譜速度VS進行對比分析（表1）。

圖4 典型節點加速度曲線對比

圖5 典型節點位移曲線對比

在圖4和圖5中，虛線表示混合有限元模型的結果曲線，實線表示三維長圓筒模型的結果曲線。從圖4中看出，沖擊波經過一定的時間傳播到結構上，在前2 ms沖擊波還未傳遞到梁單元上時，結果吻合較好。傳遞到梁單元之后簡化模型有較寬幅的振蕩，使得結果產生較大的偏差。從圖5的位移曲線能清晰地看出這種現象。

導致兩者結果產生較大偏差的原因有：

1）將部分圓筒簡化為梁單元，雖然其剛度質量等剖面屬性相等，梁單元與三維殼單元相比，梁單元只能施加入射壓力，而殼單元與流體介質之間耦合產生的散射壓力無法計算，因此忽略了殼體振動對入射壓力的緩沖作用，使得結構速度｛x˙｝增大。由式（4）解得的散射波流體質點速度｛us｝增大，由此導致式（3）中散射壓力的增大；所以結構的運動方程中的載荷矢量增大，導致響應增大；

2）應力波在結構中的傳播與結構的形狀、彈性模量、密度等有關，將部分圓筒簡化為梁單元以后，與應力波相關的這些參數發生很大變化，所以應力波的變化是導致結果規律不一致的另一可能原因。

從表1中看出譜速度的結果規律完全不一致，并有較大誤差。

表1 方案1的譜速度（單位：m／s）

4.2 方案2結果分析

下面將三維長圓筒整體模型與1／3圓筒局部模型進行比較分析。選取圓筒中心截面處4個典型節點加速度、位移曲線進行對比。如圖6、圖7所示。并對典型節點的沖擊譜速度進行對比分析（表2）。

在圖6和圖7中，虛線表示1／3圓筒局部模型的結果曲線，實線表示三維長圓筒模型的結果曲線。從圖6中看出，沖擊波經過一定的時間傳播到結構上，在整個時間段上結果吻合較好，偏差較小，沒有寬幅的振蕩。從圖7的位移曲線能清晰看出局部圓筒模型與長圓筒模型的位移變化規律一致性較好。由于位移曲線反映的是圓柱筒整體位移，局部模型受到的總載荷比長圓筒模型小，所以其整體位移較后者小。從迎爆面和背爆面曲線看出，位移曲線不完全相同，說明存在局部位移，但由于圓筒剛度較大，局部位移較小。由于在沖擊環境的考核中主要以加速度為考核對象，所以位移的影響可以忽略。

圖6 典型節點加速度曲線對比

圖7 典型節點位移曲線對比

導致兩者結果產生偏差的原因如下：

1）從結構運動方程式（1）中看出，簡化的局部圓筒模型和長圓筒整體模型相比，結構質量矩陣和剛度矩陣不同；

2）從外力公式（2）中可以看出，濕表面積減小，接受入射波壓力的面積減小，產生的散射壓力也不同，由此導致結構運動方程中的載荷矢量發生變化。由于沖擊波屬于高頻脈沖，沖擊波本身衰減較快，且在迎爆面產生較大變形消耗大部分能量，因此認為其局部影響較大。所以簡化后的模型對結果的影響有限。

從表2中看出典型節點處的譜速度值規律完全一致，并且誤差在工程允許的范圍內。誤差的大小與節點所處位置的殼體運動有關，在迎爆面所受到的載荷最大，所以其誤差也最大。

下面從設計沖擊譜的角度，再次證實方案2的準確性（圖8、圖9）。取誤差最大的迎爆面處節點的設計沖擊譜進行比較，得出譜位移誤差為8%，譜加速度誤差為11.5%。

表2 方案2的譜速度（單位：m／s）

圖8 三維圓筒整體模型迎爆面節點沖擊譜

4.3 不同爆距的結果分析

為了充分驗證上述結論，采用與上文中相同的模型、相同的藥量，分別對爆距為10 m和26.2 m的工況對3個模型進行了計算，得出的結論與上述一致。由于篇幅的原因這里不再詳細列出加速度、位移曲線結果，僅列出譜速度的對比結果。

圖9 局部圓筒迎爆面節點沖擊譜

表3 爆距為10 m時，譜速度對比結果（單位：m／s）

表4 爆距為26.2 m時，譜速度對比結果（單位：m／s）

5 結語

從以上結果可以得出結論，考慮流固耦合作用下，大型水下圓柱殼結構遭受水下爆炸沖擊時，簡化的局部模型引起的質量和剛度的差別，對計算結果的影響在工程允許的范圍內，為了提高計算效率，建立部分結構模型進行沖擊響應計算是可行的。使用特殊梁剖面簡化的混合模型由于應力波和載荷加載的影響，無法正確模擬，有待于進一步研究。

［1］夏利娟，吳衛國，翁長儉，等.混合有限元模型的集成方法［J］.上海交通大學學報，2001，35（4）：544－547.

［2］吳梵，郭日修.艦船尾部振動分析整船混合有限元模型的集成方法研究［J］.海軍工程學院學報，1999（4）：.

［3］汪俊.加筋圓柱殼在水下爆炸和靜水壓力聯合作用下的動響應研究［D］.CSSRC.2005.

［4］ABAQUS 6.8 Documentation［C］，2006.

［5］STILLMAN D W，HALLQUIST J O.LS-DYNA3D User’s Manual［S］.LSTC Report，1990.

［6］BELYTSCHKO T.，LIN J.I.，TSAY C S.Explicit algorithms for the nonlinear dynamics of shells［M］.Elsevier，1984.

［7］GEERS T L.Boundary element analysis for transient response.In Computational Methods for Transient Analysis［M］.Elsevier，New York，1983.

［8］PARK K C.Partitioned transient analysis procedures for coupled field problems：stability analysis［J］.Journal of Applied Mechanics，1980，47，（4）：370－376.

Method of Simplifying FEM Model of Underwater Cylinder in the Shock Response Analysis

Sun Zhan－zhong Hu Gang－yi Zhang Xin－yu Huang Guo－bing
China Ship Development and Design Center，Wuhan 430064，China

Based on the explicit finite element code ABAQUS，non－linear dynamic response of an immersed cylinder subjected to underwater explosion was analyzed.The response of different simplified models subjected to underwater explosion and the comparisons with the results of the whole model were comprehensively discussed.The results show that satisfactory simulation can be obtained by establishing a part of the model of the cylinder instead of the whole model.

model；simplifying；underwater explosion；impact response

U663.1

：A

：1673－3185（2009）01－22－07

2008-12-20

“十一五”預先研究課題

孫占忠（1983－），男，碩士研究生。研究方向：水下結構沖擊響應分析。E-mail：a191221468＠163.com胡剛義（1966－），男，研究員，碩士生導師。研究方向：艦船結構動態響應分析