朱俊杰，盧永剛，馮星奎，陶俊林

(1. 西南科技大學(xué) 土木工程與建筑學(xué)院，四川 綿陽(yáng) 621010；2. 中國(guó)工程物理研究院 總體工程研究所，四川 綿陽(yáng) 621000；3. 火箭軍駐中國(guó)工程物理研究院 型號(hào)辦事處，四川 綿陽(yáng) 621900)

殼單元和實(shí)體單元模擬爆炸荷載作用鋼板動(dòng)力響應(yīng)的差異性比較

朱俊杰1，盧永剛2，馮星奎3，陶俊林1

(1. 西南科技大學(xué) 土木工程與建筑學(xué)院，四川 綿陽(yáng) 621010；2. 中國(guó)工程物理研究院 總體工程研究所，四川 綿陽(yáng) 621000；3. 火箭軍駐中國(guó)工程物理研究院 型號(hào)辦事處，四川 綿陽(yáng) 621900)

以艦船艙室鋼板為研究對(duì)象，利用動(dòng)力分析有限元程序 Ls-dyan，分別采用殼單元和實(shí)體單元 2 種建模方式對(duì)爆炸荷載作用下鋼板的動(dòng)力響應(yīng)進(jìn)行數(shù)值模擬研究。對(duì)比鋼板的應(yīng)力分布云圖、中線上各點(diǎn)最大位移曲線圖，發(fā)現(xiàn) 2 種方法計(jì)算結(jié)果與試驗(yàn)值大致吻合，采用殼單元建模在保證計(jì)算結(jié)果準(zhǔn)確性的前提下可以大大提高計(jì)算效率。如何正確地選擇單元進(jìn)行分析是數(shù)值模擬很重要的一個(gè)環(huán)節(jié)，通過(guò)對(duì)殼單元和實(shí)體單元計(jì)算結(jié)果的初步對(duì)比，為艦船艙室數(shù)值模擬研究如何選取單元提供參考。

船用鋼板；爆炸荷載；數(shù)值模擬；殼單元

0 引 言

艦艇是由縱橫構(gòu)架和殼板結(jié)構(gòu)組成，其中殼板主要承受水平荷載。艦船艙室結(jié)構(gòu)承受的主要載荷之一是戰(zhàn)斗部爆炸產(chǎn)生的沖擊波載荷[1]，鋼板在受到具有傳播速度快、峰值荷載大、作用時(shí)間短等特征的爆炸沖擊波載荷作用時(shí)，其受力特性和動(dòng)態(tài)響應(yīng)顯著異于靜載情況下，且十分復(fù)雜，故而研究鋼板在爆炸荷載作用下的動(dòng)力響應(yīng)很有必要，是研究艙室結(jié)構(gòu)的基礎(chǔ)。目前國(guó)內(nèi)對(duì)艦船艙室方面的研究大多集中在艙室結(jié)構(gòu)的仿真分析上，朱建方等[2]對(duì)艦艇艙室在爆炸荷載下的動(dòng)態(tài)響應(yīng)過(guò)程進(jìn)行了數(shù)值模擬，得到了艦艇艙室在爆炸沖擊波作用下的響應(yīng)特性和破壞機(jī)理。鄢順偉等[3]仿真計(jì)算了不同壁厚艙室結(jié)構(gòu)在反艦導(dǎo)彈戰(zhàn)斗部?jī)?nèi)爆作用下的毀傷效應(yīng)。但是，對(duì)船用鋼板在爆炸荷載作用下的響應(yīng)研究較少。本文采用數(shù)值模擬的方法對(duì)爆炸載荷下四邊固支的鋼板進(jìn)行研究。在利用數(shù)值分析研究時(shí)，單元類型的選取對(duì)計(jì)算效率有很大的影響。鋼板屬于板結(jié)構(gòu)，實(shí)體單元與殼單元均可以對(duì)其進(jìn)行模擬，但是 2 種建模方法得到結(jié)果的差異值得探討。通過(guò)分析兩者的差別，可以更好地為鋼板使用何種單元在保證其結(jié)果準(zhǔn)確性的基礎(chǔ)上又能提高計(jì)算效率提供參考。為此，本文以文獻(xiàn)中試驗(yàn)為基礎(chǔ)[4–5]，在驗(yàn)證數(shù)值計(jì)算方法可靠的基礎(chǔ)上，分別采用殼單元與實(shí)體單元對(duì)四邊固支的方形板在爆炸荷載作用下的動(dòng)力響應(yīng)進(jìn)行數(shù)值模擬研究，考察使用不同的建模單元類型對(duì)計(jì)算結(jié)果的影響[6]，并對(duì)結(jié)果的差異性進(jìn)行了比較與分析。

1 Shell163 單元介紹

Shell163 是一個(gè)具有 4 節(jié)點(diǎn)單元，有彎曲和膜特征，可施加平面和法向載荷。單元在每個(gè)節(jié)點(diǎn)上有 12個(gè)自由度，分別為在節(jié)點(diǎn) x，y 和 z 方向的平動(dòng)、加速度、速度和繞 x，y 和 z 軸的轉(zhuǎn)動(dòng)。Shell163 單元一共有 12 種算法，這里選第 10 種算法，對(duì)于翹曲情況可以得到正確結(jié)果。和實(shí)體單元一樣，積分點(diǎn)的個(gè)數(shù)直接影響著 CPU 的時(shí)間。典型地，對(duì)于彈性材料沿厚度方向需要 2 個(gè)積分點(diǎn)，而對(duì)于塑性材料則需要 3 個(gè)或更多的積分點(diǎn)[7]。

2 數(shù)值計(jì)算模型

2.1 算法選取

分析主要涉及到流體流動(dòng)、固體大變形以及流體與固體相互作用等問(wèn)題，Ls-dyna 程序中提供了 Lagrange，Euler 和 ALE（任意拉格郎日歐拉算法）3 種不同的算法，其中 ALE 算法兼具 Lagrange 算法和 Euler算法的優(yōu)點(diǎn)，可以描述流體與固體的相互耦合作用；同時(shí)在求解過(guò)程中網(wǎng)格還可以根據(jù)定義的參數(shù)適當(dāng)調(diào)整位置，克服網(wǎng)格嚴(yán)重畸變的問(wèn)題[8]。因此，本文中選取 Ls-dyna 程序提供的 ALE 算法建立流體與固體間的聯(lián)系，實(shí)現(xiàn)流固耦合動(dòng)態(tài)分析，模擬方形鋼板在爆炸荷載作用下的響應(yīng)和破壞。

2.2 材料本構(gòu)與狀態(tài)方程

結(jié)構(gòu)材料在爆炸沖擊載荷作用下會(huì)產(chǎn)生巨大的變形，與受靜載荷作用的情況相比，動(dòng)態(tài)屈服應(yīng)力比靜態(tài)屈服應(yīng)力有明顯提高，材料力學(xué)性能表現(xiàn)相對(duì)也更為復(fù)雜。固支方板采用低碳鋼，鋼板選用與應(yīng)變率相關(guān)的各向同性塑性隨動(dòng)硬化模型，材料參數(shù)見(jiàn)表 1。

表 1 鋼板材料參數(shù)Tab. 1 Steel plate material parameters

炸藥在空氣中爆炸時(shí)，會(huì)在極短的時(shí)間內(nèi)把化學(xué)能全部轉(zhuǎn)變?yōu)闊崮?，?dǎo)致爆轟產(chǎn)物向周圍迅速膨脹而做功，形成空氣沖擊波。文中采用 *HIGH_EXPLOSIVE_BURN 模型作為炸藥的材料模型，采用 Jones-Wilkins-Lee（JWL）狀態(tài)方程來(lái)描述，即

式中：A ，B，R1，R2和 ω 為 JWL 狀態(tài)方程參數(shù)；E為炸藥的內(nèi)能；V 為當(dāng)前相對(duì)體積。

對(duì)空氣采用 *MAT_NULL 材料模型以及線性多項(xiàng)式 *EOS_LINEAR_POLYNOMI 狀態(tài)方程加以描述，即

式中：C0～C6為常量；e 為初始體積內(nèi)能；μ = 1/V0–1；V0為相對(duì)體積；C1～C4均為 0。

2.3 有限元模型

1） 空氣網(wǎng)格的驗(yàn)證

在數(shù)值模擬中，網(wǎng)格的劃分至關(guān)重要，網(wǎng)格的大小不僅關(guān)系著數(shù)值模擬計(jì)算的精度，也關(guān)系著數(shù)值模擬計(jì)算所用的時(shí)間。一般而言，網(wǎng)格越小則計(jì)算的結(jié)果越接近于真實(shí)情況，但隨著網(wǎng)格的減少計(jì)算時(shí)間也急劇增大，所以選擇合理的網(wǎng)格尺寸很有必要。網(wǎng)格驗(yàn)證幾何模型尺寸為空氣 120 mm × 120 mm × 60 mm，炸藥為 27 mm × 27 mm × 27 mm。利用對(duì)稱性取 1/4 模型進(jìn)行構(gòu)建，模型如圖 2 所示。

表 2 炸藥與狀態(tài)方程參數(shù)Tab. 2 Explosive and equation of state parameters

表 3 多項(xiàng)式狀態(tài)方程參數(shù)Tab. 3 Polynomial equation of state parameters

空氣與炸藥均采用 8 節(jié)點(diǎn)實(shí)體單元，TNT 當(dāng)量為 20 g，炸藥中心位置起爆。在對(duì)稱面上施加對(duì)稱約束，四周（除對(duì)稱面）施加固定邊界條件。為了模擬出無(wú)限空氣域，空氣模型邊界施加透射邊界。工況如表 4 所示。

表 4 空氣網(wǎng)格驗(yàn)證工況表Tab. 4 Air grid authentication working table

上述工況在距離炸藥中心 6 cm 和 8 cm 位置處的超壓如圖 3～圖 4 所示。

由圖 5～圖6 可以看出，隨著網(wǎng)格的減少，壓力的峰值越來(lái)越大，當(dāng)網(wǎng)格取 0.3 cm 和 0.2 cm 時(shí)壓力峰值的差值約為 5%，考慮到計(jì)算精度和時(shí)間的因素空氣網(wǎng)格取 0.3 cm。

2） 炸藥網(wǎng)格的驗(yàn)證

根據(jù)上面網(wǎng)格確定炸藥網(wǎng)格的驗(yàn)證工況如表 5 所示。

表 5 炸藥網(wǎng)格驗(yàn)證工況表Tab. 5 Explosive grid authentication working table

由圖 5～圖 6 可知，隨著炸藥網(wǎng)格的減少，壓力峰值有所增加，當(dāng)網(wǎng)格取為 0.2 cm 和 0.1 cm 時(shí)的壓力峰值差值約為 2%，炸藥網(wǎng)格選取 0.2 cm 較為合適。

3） 2 種單元計(jì)算模型

采用 Ansys/Ls-dyna 建立有限元模型，對(duì)于固支鋼板在爆炸載荷的作用過(guò)程進(jìn)行了數(shù)值模擬研究。固支鋼板尺寸為 63 mm × 63 mm，厚度為 1.6 mm，建立 1/4模型。方板尺寸如圖 7 所示。

數(shù)值模型由炸藥、空氣、鋼板 3 部分組成，采用殼單元和實(shí)體單元分別建立板模型，2 種模擬方式除了方板的單元類型不同，其他參數(shù)均一致。空氣單元尺寸為 3 mm，炸藥單元尺寸為 2 mm，板單元為 0.4 mm，有限元模型如圖 8～圖9 所示。

2.4 失效準(zhǔn)則

采用 *MAT_ADD_EROSION 對(duì)鋼板材料失效準(zhǔn)則進(jìn)行定義，塑性失效應(yīng)變 εmax= 0.30，考慮到在沖擊荷載作用下薄壁材料容易出現(xiàn)斷裂、絕熱剪切等破壞情況，在計(jì)算模型中增加最大剪應(yīng)變失效準(zhǔn)則，剪切失效應(yīng)變 τmax= 0.20[9]。利用單元?jiǎng)h除技術(shù)來(lái)模擬鋼板的破壞情況，即當(dāng)單元變形引起的塑性應(yīng)變或剪切應(yīng)變超過(guò)臨界值時(shí)，就認(rèn)為該單元破壞。

3 殼單元與實(shí)體單元模擬結(jié)果比較

3.1 方形鋼板的動(dòng)態(tài)響應(yīng)

1）t = 0.15 ms 時(shí)刻，方形鋼板 Von Mises 應(yīng)力云圖對(duì)比。

從圖 10 可知，無(wú)論是采用殼單元還是實(shí)體單元，方形板中的應(yīng)力分布大致相同，最大值均出現(xiàn)在炸藥最上方中心位置處。同一時(shí)刻，實(shí)體單元中最大應(yīng)力值為 503.6 MPa，殼單元中最大應(yīng)力值為 504.7 MPa。

從圖 11 可知，方形鋼板的變形大致呈現(xiàn)為中心處的外鼓，呈“金字塔”狀，開(kāi)始先發(fā)生中心處的翹曲，隨著時(shí)間的增加，翹曲越來(lái)越嚴(yán)重，而后慢慢下降，2 種方法模擬的鋼板變形發(fā)展趨勢(shì)與試驗(yàn)結(jié)果變形趨勢(shì)相同。

2）鋼板中線上各點(diǎn)處最大位移值對(duì)比

圖 12給出了在 20 g 炸藥量產(chǎn)生的爆炸荷載作用下鋼板中線上各點(diǎn)的最大位移曲線。由圖中可以看出，2種方法得到的各點(diǎn)處最大位移曲線基本一致，與試驗(yàn)值較為吻合，各點(diǎn)處最大位移差值均在 10% 以內(nèi)。最大的位移發(fā)生在鋼板中心位置處，2 種模擬方法得到的結(jié)果與試驗(yàn)值差值均在 5% 以內(nèi)，其他位置處的差值較小。

3.2 兩種方法計(jì)算時(shí)間對(duì)比

2 種模型的炸藥和空氣域單元數(shù)目相同，而使用殼單元建模使得鋼板單元數(shù)目減少 50%。由圖 13 可以看出使用殼單元建模計(jì)算時(shí)間降低約 50%，顯著提高計(jì)算效率。尤其對(duì)于單元數(shù)目比較多的模型，殼單元的優(yōu)勢(shì)更為突出。

無(wú)論是從鋼板的動(dòng)力響應(yīng)歷程還是鋼板的變形特征層面上的對(duì)比驗(yàn)證，其結(jié)構(gòu)均表明：本次數(shù)值模擬的建模技術(shù)和參數(shù)選取合理、可靠，其數(shù)值模擬的結(jié)果可以艦船艙室內(nèi)爆選取單元提供依據(jù)。

4 結(jié) 語(yǔ)

對(duì)方形鋼板在爆炸荷載作用下的響應(yīng)，分別采用殼單元和實(shí)體單元進(jìn)行模擬。對(duì)比 2 種方法中板的應(yīng)力分布和中線處各點(diǎn)的最大位移值，得到的結(jié)果大致相同，差值在 10% 以內(nèi)。結(jié)果表明使用實(shí)體單元建模和使用殼單元建模對(duì)于等厚度方板的計(jì)算結(jié)果影響不大。而采用殼單元建模相比采用實(shí)體單元建?？梢燥@著提高計(jì)算效率。

數(shù)值分析方法在爆炸荷載作用下艦船艙室的動(dòng)力響應(yīng)研究中有著廣泛的應(yīng)用，但如何正確地選擇單元進(jìn)行分析是數(shù)值模擬中很重要的一個(gè)環(huán)節(jié)。本文通過(guò)對(duì)殼單元和實(shí)體單元模擬方形鋼板計(jì)算結(jié)果的初步對(duì)比，為艦船艙室模型中如何選擇單元提供參考。

[1]陳攀, 劉志忠. 艙室內(nèi)爆沖擊波載荷特性及影響因素[J]. 艦船科學(xué)技術(shù), 2016, 38(2): 43–48. CHEN Pan, LIU Zhi-zhong. Research on loading of explosive and influencing factors inside closed cabin [J]. Ship Science and Technology, 2016, 38(2): 43–48.

[2]朱建方, 王偉力, 曾亮. 艦艇艙室內(nèi)爆毀傷的建模與仿真分析[J]. 系統(tǒng)仿真學(xué)報(bào), 2009, 21(22): 7066–7068.

[3]杜志鵬, 李曉彬, 夏利娟, 等. 反艦導(dǎo)彈攻擊艦船舷側(cè)防護(hù)結(jié)構(gòu)過(guò)程數(shù)值仿真[J]. 哈爾濱工程大學(xué)學(xué)報(bào), 2006. 27(4): 484–487.

[4]LONGDON G S, CHUNG Ki Y S, NURICK G N. Experimental and numercial studiesion the response of quadrangular stiffened plates. Part I: subjected to uniform blast load [J]. In-ternational Journal of Impact Engineering. 2005, (31): 55–83.

[5]LONGDON G S, CHUNG Ki Y S, NURICK G N. Experimental and numercial studiesion the response of quadrangular stiffened plates. Part II: localist blast loading [J]. International Journal of Impact Engineering. 2005, (31): 85–11.

[6]徐新光. 梁?jiǎn)卧蛯?shí)體單元模擬框架結(jié)構(gòu)的差異性比較[J].四川建材, 2014, 40(5): 46–48.

[7]ANSYS工程結(jié)構(gòu)數(shù)值分析[M]. 北京. 人民交通出版社, 2012: 11–13.

[8]ANSYS/LSDYNA在爆炸與沖擊領(lǐng)域內(nèi)的工程應(yīng)用[M]. 北京:中國(guó)建筑工業(yè)出版社, 2011.

[9]李帆. 船艦室內(nèi)爆炸破壞的數(shù)值模擬和試驗(yàn)研究[D]. 北京: 國(guó)防科技大學(xué), 2012.

Comparision of differentiation with the steel plates' dynamic response subjected to blast loads using shell element and solid element modeling methods

ZHU Jun-jie1, LU Yong-gang2, FENG Xing-kui3, TAO Jun-lin1
(1. School of Civil Engineering and Architecture, Southwest University of Science and Technology, Mianyang 621010, China; 2. Institute of Structural Mechanics,China Academy of Engineering Physics, Mianyang 621900, China; 3. The Office of Rocket Force, China Academy of Eugineering Physicies, Mianyang 621900, China)

Take the steel plates of warship cabins as research objects, numerical simulation study on the response of the quadrangular steel plates subjected to blast loads is presented by dynamic analysis finite element procedure Ls-dyan3d, using modeling methods of shell element and solid element respectively. Comparing pressure distribution cloud picture and each point of midcourt line miximum displacement curve graph, it can be found that transformation trend of steel plates is similar for different thickness while each point of midcourt line miximum displacement is discriminating. However results from two methods are alike for equal thickness. Meanwhile, it can improve calculative effciency. How to choose which kind of element types correctly to analysis is an important link. Through the preliminary contrast study of calculative results of shell element and solid element, it can provide reference for choosing element type in warship cabins numerical simulation.

marine steel plate；blast loads；numerical simulation；shell element

TJ03

A

1672–7619(2017)03–0018–05

10.3404/j.issn.1672–7619.2017.03.004

2016–09–18；

2016–10–27

國(guó)家自然科學(xué)基金資助項(xiàng)目(11502258)

朱俊杰(1991–)，女，碩士研究生，研究方向?yàn)榻Y(jié)構(gòu)沖擊與爆炸。

陶俊林(1972–)，男，博士，教授，E-mail:junlintao@126.com。