復合圓柱殼沖擊壓縮數值模擬及穩定性研究*

謝富佩，徐 緋，曾 卓，周中玉，谷卓偉

（1. 西北工業大學航空學院計算力學與工程應用研究所，陜西 西安 710072；2. 中國航空發動機集團有限公司湖南動力機械研究所，湖南 株洲 412002；3. 中國工程物理研究院流體物理研究所，四川 綿陽 621999）

復合圓柱殼體結構具有高比剛度、高比吸能以及輕質等優點，在航空、航海、汽車工業及軍事等領域有廣泛應用。在強磁場物理、材料合成等領域，使用了一種利用炸藥爆轟驅動復合圓柱殼內部預置磁場以得到高能磁場的裝置，通常將此裝置稱為柱面內爆磁通量壓縮發生器(explosive magnetic flux compression generator)，簡稱MC-1 裝置[1]。與常見的纖維編織復合圓柱殼不同，由于考慮到復合圓柱殼通電后產生并壓縮磁場的要求，此類復合結構中銅線纖維是連續平行排列的；在炸藥爆轟產生的沖擊載荷作用下，MC-1 裝置中復合圓柱殼容易發生界面失穩和擾動增長而導致結構失穩[2]；另外，結構也可能在沖擊壓縮過程中發生較大的塑性變形而失效，而復合圓柱殼的結構響應及穩定性是影響MC-1 裝置性能的重要因素。目前MC-1 裝置在高能量密度物理實驗技術研究領域應用較為廣泛[3]，因此對復合圓柱殼體結構在沖擊壓縮下的結構響應模擬及穩定性研究有十分重要的意義。

目前國內外對于MC-1 裝置中的復合圓柱殼在沖擊壓縮下的數值模擬研究較少。俄羅斯實驗物理科學研究所于2015 年公布了幾種MC-1 裝置的設計參數，但未公開其參數設計細節；陸禹等[4]使用了引入解析形式狀態方程組的一維內爆磁流體力學編碼MC11D，對MC-1 裝置工作過程進行了一維條件下的數值模擬；張春波等[5]利用AUTODYN 軟件二次開發接口建立了內爆壓縮多層密繞螺線管復合圓柱殼三維模型，并實現了周期性邊界條件，使用光滑粒子流體動力學方法對其動態響應及不穩定性開展了數值模擬；趙士操等[6]使用SPH 方法構建了纖維織物的細觀尺度計算模型，使用AUTODYN 軟件對超高速碰撞下纖維增強復合材料模型進行了分析計算；劉軍等[7]對爆轟驅動內壁刻有正弦擾動的金屬鋼殼與內部硅橡膠界面產生不穩定性問題進行了數值模擬。由以上研究可以看出，對炸藥爆轟驅動作用下的復合圓柱殼動力學響應及穩定性數值模擬問題，由于涉及炸藥爆炸、高速沖擊等過程，傳統網格法面臨大變形、高度非均勻性、運動物質交界面及邊界處理等諸多問題，而粒子法又有計算耗時較長的問題，因此對復合圓柱殼細節模型的數值模擬常存在較大困難，導致現有的研究缺乏對復合圓柱殼細節模型的分析手段，且實驗對結構失穩現象未提出合理的失穩判據，復合圓柱殼穩定性的影響因素難以評估。

基于以上缺陷，本文以MC-1 裝置中的復合圓柱殼作為研究對象，使用光滑粒子動力學-有限元方法（smooth particle hydrodynamic-finite element method, SPH-FEM）耦合算法，構建復合圓柱殼二維細節爆轟模型并對其在炸藥爆轟沖擊下的結構響應進行分析計算，提出一種判斷圓柱殼結構在沖擊壓縮下的穩定性評價方式，分析影響MC-1 裝置性能的因素，以期為MC-1 裝置設計及相關實驗研究提供參考。

1 數值計算模型

1.1 模型介紹

本文計算模型的結構參數參照俄羅斯實驗物理科學研究所核中心于2015 年公布的實驗參數，裝置結構如圖1(a)所示，復合圓柱殼外側包裹著一層環形炸藥，炸藥外側均勻布置48 個起爆點，在同一時刻起爆。裝置中的復合圓柱殼在制備時首先使用一定直徑的漆包銅導線按照規定的角度螺旋繞制成多匝多層的螺線管內筒體，將繞制的螺線管表面使用環氧混合劑固定并形成一層絕緣層，然后將從螺線管伸出的銅線覆蓋在絕緣層表面，覆蓋的銅線方向與復合圓柱殼軸線方向平行，最后在其表面使用環氧混合劑凝固。因此MC-1 裝置中復合圓柱殼由內向外分別為銅線層、內環氧層、銅線折返層及外環氧層，如圖1(b)所示。圖1(c)為公布的實驗裝置的結構參數示意圖，考慮到理想情況下銅線層及銅線折返層的厚度為銅線直徑的整數倍，因此計算時將銅線層的厚度調整為2 mm，進行模擬計算時裝置的參數如圖1(d)所示。

圖1 復合圓柱殼示意圖Fig. 1 Diagram of composite cylindrical shell

1.2 數值模型的建立

對在爆轟驅動作用下復合圓柱殼結構響應的模擬計算問題，由于復合圓柱殼三維模型較為復雜，且計算量較大，而復合圓柱殼結構沿軸線方向的截面形狀及尺寸無明顯變化，因此本文選用對復合圓柱殼具有代表性的二維截面模型進行計算。考慮到只使用有限元算法時裝置內圓柱殼結構可能會出現的網格畸變問題以及只使用SPH 算法可能導致的計算規模較大的問題，本文使用SPH-FEM 耦合算法[8]對爆轟沖擊下的復合圓柱殼進行數值模擬，即對裝置炸藥部分使用FEM 方法劃分網格，對本問題的主要研究對象復合圓柱殼部分使用SPH 方法建立模型，由于兩種方法都具有明顯的拉格朗日特性，因此只需要在軟件中設置接觸域即可實現數據傳遞，完成兩種算法的耦合。在建立模型時，由于復合圓柱殼結構較為復雜，直接在AUTODYN 軟件建模較為困難，本文使用導入外部文件[9]的方式建立模型。首先在CATIA 軟件中建立復合圓柱殼細節模型并導入到有限元軟件HyperMesh 中劃分三角形網格，在劃分網格時，網格尺寸即為后續計算的粒子間距，本文將粒子間距取為0.1，導出K 文件后編寫C 程序將模型節點格式轉化為AUTODYN 軟件可識別的SPH 粒子格式，將其存儲為Formatted 文件后導入AUTODYN軟件，進一步添加有限元炸藥，完成模型邊界和起爆點等的設置后即可完成建模。圖2 為最終在AUTODYN 中建立的15°銅線螺旋角下的復合圓柱殼數值模型。

圖2 15°螺旋角下的爆轟驅動復合圓柱殼模型Fig. 2 Detonation driven composite cylindrical shell model with a spiral angle of 15°

在數值模擬中，復合圓柱殼外側用于提供爆轟壓力的炸藥部件使用JWL 狀態方程[10]描述：

表1 Comp B 炸藥主要材料參數Table 1 Main material parameters of Comp B explosive

復合圓柱殼中銅線類型為高導無氧銅，其材料參數[11]如表2 所示，塑性部分使用Johnson-Cook 材料模型，材料的屈服應力可定義為：

表2 銅線主要材料參數Table 2 Main material parameters of copper wire

圓柱殼的環氧部分的材料模型則使用Shock 狀態方程描述:

式中：U為激波速度；up為粒子速度；c0和s為常數，參數取值參照AUTODYN 材料庫，c0取為2.73 km/s，s取為1.493。

數值計算模型中圓柱殼與炸藥之間的接觸定義為Lagrange/Lagrange 接觸，炸藥外側設置48 個均布起爆點，均在0 μs 時刻同時起爆。

2 結構響應與穩定性分析

2.1 數值模型結果驗證與分析

為了驗證建立的數值模型計算結構響應的準確性，在AUTODYN 中建立炸藥爆轟驅動金屬環氧復合套筒模型進行數值模擬，套筒的內徑為94 mm，外徑為100 mm，沿徑向由內向外分別為厚度為1.5 mm鋼套筒和1.5 mm 的環氧套筒，外部環形炸藥的厚度為55 mm。將計算的套筒內側半徑變化曲線與使用已經過驗證的Triangels-MHD 磁流體力學程序[12]計算的結果進行對比，如圖3 所示。

圖3 套筒內側半徑變化曲線對比Fig. 3 Comparison of radius change curves within the sleeve

根據半徑變化對比曲線可以很明顯地看出，除了曲線后半段外，使用兩種方法計算得到的套筒內壁位移響應曲線大部分階段都比較接近，而曲線后半段出現一定差異是因為TriAngels-MHD程序主要是基于磁流體力學方程建立的，在求解過程中考慮到了磁力的影響，在壓縮后期磁力會由于半徑減小而突增，使用AUTODYN 軟件求解時則暫未考慮到后期磁力影響。但是由于磁力影響只占套筒壓縮過程中的極小一部分，因此可以認為本文使用的二維數值計算模型能較好地模擬圓柱殼在炸藥沖擊壓縮下的結構響應。

而由于金屬環氧套筒的計算模型中未引入初始擾動，整個爆轟過程中各個起爆點的起爆時間相同，起爆點分布也較為密集，使得套筒外側的爆轟載荷沿環向分布較為均勻，因此在壓縮過程中套筒無較大的環向擾動，基本保持為一個圓環，如圖4(a)所示。同樣地，復合圓柱殼的變形也有此特點，以15°銅線螺旋角的復合圓柱殼模型為例，在爆轟作用下其不同時刻的變形如圖4(b)所示，從圖中可以看出，在壓縮過程中，復合圓柱殼的環氧層及銅線折返層有較大的環向變形，而決定裝置整體性能的銅線層則在壓縮前中期基本保持為環向相對光滑的狀態，直至壓縮后期才有輕微的環向擾動。

圖4 圓柱殼在不同時刻的變形Fig. 4 Deformation patterns of cylindrical shell at different times

2.2 穩定性判據

在MC-1 裝置的實際應用中，復合圓柱殼受到爆轟載荷后向內收縮是一種高壓高應變率條件下的大變形行為，在此階段中極易發生界面不穩定增長，而復合圓柱殼失穩是影響MC-1 裝置性能的重要因素，因此在數值模擬中需引入一種失穩判據，以對復合圓柱殼在沖擊壓縮過程中的穩定性進行評價。

首先以金屬環氧套筒為研究對象，圖5(a)為計算的金屬環氧套筒在沖擊壓縮過程中內壁粒子的速度變化曲線，曲線可大致分為3 段：(1) 炸藥爆轟波尚未到達套筒表面，結構尚無響應，速度為0；(2) 爆轟波到達套筒表面，套筒在爆轟波作用下開始向內收縮，內壁粒子速度急劇增加后逐漸趨于平穩；(3) 結構出現擾動增長，內壁粒子速度進一步增加并出現不穩定增長。根據套筒內壁粒子速度變化歷史曲線可大致地判斷出結構在第3 階段出現失穩。為了獲得更加精確的失穩時間和失穩時對應的壓縮率等數據，對第3 階段的速度曲線進行擬合并對擬合曲線進行求導，圖5(a)左上圖為對速度曲線第3 階段擬合曲線求一階導的導數曲線，導數曲線出現明顯轉折的時間即為套筒結構的失穩時間。經過曲線處理后，可判斷出金屬環氧套筒結構的失穩時間大致在15.20 μs 左右，此時套筒的內半徑為18.02 mm，整個結構壓縮率大致為61.66%。

圖5 金屬環氧套筒失穩時間判斷曲線Fig. 5 Judgment curves of instability time of metal epoxy sleeve

為了證明提出的失穩判據的合理性，沿金屬環氧套筒徑向由外向內等間距取點，并將其標記為Node 1～Node 5，提取其沖擊壓縮作用下的速度變化曲線，如圖5(b)所示，可以看到在炸藥爆轟載荷的作用下，金屬環氧套筒徑向各層粒子在壓縮前期變化規律基本一致，而壓縮至15 μs 左右時各層粒子的速度曲線開始出現較為明顯的差異，內層粒子的速度增加，而外層粒子速度減小，這種速度差異勢必會導致結構徑向的失穩，而此失穩時間也與前面提出的使用內壁粒子速度曲線變化判斷的失穩時間十分接近，由此可以證明提出的失穩判據較為合理。

將上述失穩判據推廣到復合圓柱殼結構，可以發現相同的規律。以15°銅線螺旋角的復合圓柱殼為例，圖6(a)為復合圓柱殼結構內壁粒子的速度變化曲線，也有明顯的三階段特征，圖6(a)左上圖為將復合圓柱殼內壁粒子速度變化曲線的第3 階段曲線擬合求導后的曲線，可以看到曲線在24.76 μs 時開始有突增的趨勢，而此時刻與圖6(b)所示的復合圓柱殼沿徑向各層粒子速度開始出現明顯差異的時間相近，圖6(b)中Node 1～Node 11 指沿復合套筒徑向由外向內等距標記的11 個粒子。綜合以上分析，可以證明本文提出的判斷失穩時間的判據對復合圓柱殼結構的穩定性分析有一定的參考價值。在后面的研究中，將繼續使用根據結構內壁粒子速度變化歷史判斷失穩時刻的判據對不同條件下的復合圓柱殼結構的穩定性進行分析討論。

圖6 復合圓柱套筒失穩時間判斷曲線Fig. 6 Judgment curves of instability time of composite cylindrical shell

3 影響復合圓柱殼穩定性因素分析

3.1 復合圓柱殼制備缺陷對穩定性的影響

由復合圓柱殼的制備過程可知，銅線折返層的銅線是復合圓柱殼內層的銅線層延伸出來的，因此在理想情況下復合圓柱殼的銅線層和銅線折返層銅線數量應相同。在繞制復合圓柱殼結構的內層銅線層時，銅線層橢圓截面的長軸長b計算公式為

式中：D為銅線直徑，θ 為銅線的螺旋角度，同樣以15°螺旋角的復合圓柱殼模型為例，內層銅線層銅線截面的長軸長約為0.965 9 mm，則內層每層銅線的數量ni可表示為

式中：Ri為每層銅線的銅線中心距圓柱殼中心的距離，l為銅線層數，按下式計算

式中：h為模型中銅線層的厚度，對本文中計算的銅線直徑為0.25 mm 的模型，結合式(5)和式(6)估算出銅線層銅線的數量約為3668 根。

而理想條件下，將模型的銅線折返層全部排滿大致需要的銅線數量為5623 根，顯然在實際制備復合圓柱殼時并不能將銅線折返層完全由銅線排滿，在銅線折返層極有可能存在一定的缺陷。因此建立含有缺陷的復合圓柱殼模型來研究復合圓柱殼在制備過程中的缺陷對其穩定性的影響，具體做法為在復合圓柱殼折返層建模時隨機刪去若干銅線，本文中含缺陷模型隨機在銅線折返層外側刪去100 根銅線，缺陷部分細節如圖7(a)所示，模型共有4 處這樣的缺陷。根據計算結果可以發現，由于銅線折返層缺陷部分強度較小，在炸藥爆轟壓縮復合圓柱殼的過程中會導致強度較弱的部分銅線折返層及內環氧層侵入銅線層，并且隨著爆轟波的運動，侵入部分的深度增加，影響銅線層的變形形態，進而會影響整個結構的穩定性，不同時刻含缺陷模型的復合圓柱殼局部形態如圖7(b)所示，圖7(c)為含缺陷模型復合圓柱殼在24.58 μs 的形態，可以看到其結構內壁環向不再保持為圓形，圖7(d)為失穩時間判斷曲線，含缺陷結構的失穩時間為22.63 μs，裝置的穩定性明顯降低。以上對含缺陷結構的復合圓柱殼的計算表明：復合圓柱殼結構折返層的制備缺陷會影響結構穩定性。

圖7 含缺陷復合圓柱殼變形圖及失穩時刻判斷曲線Fig. 7 Deformation diagrams of defective composite cylindrical shell and the judgement curve of instability time

3.2 銅線螺旋角對穩定性的影響

復合圓柱殼是非密實、非均勻的結構，銅線的螺旋角度可能是影響其穩定性的重要因素[13]?？紤]到在實際應用中復合圓柱殼中銅線螺旋角通常在20°以內，故建立螺旋角度分別為0°和10°的復合圓柱殼模型，不同銅線螺旋角的圓柱殼細節如圖8 所示。除螺旋角度外其余參數及條件均與前面的15°螺旋角模型相同，建模時不同角度模型的粒子數也相近。比較不同銅線螺旋角度下復合圓柱殼的結構響應，使用本文第2 節提到的速度判斷失穩判據對模型的失穩時間進行判斷，結合計算的復合圓柱殼位移變化結果，不同螺旋角模型的失穩時間及失穩時對應的壓縮率比較如表3 所示。

表3 不同螺旋角復合圓柱殼失穩時間及壓縮率對比Table 3 Comparison of instability time and compression ratio of composite cylindrical shell with different spiral angles

圖8 不同螺旋角復合圓柱殼細節模型Fig. 8 Detailed composite cylindrical shell models at different spiral angle

實際上，在建立幾何模型時，可以推算出復合圓柱殼銅線部分的截面積S與銅線的螺旋角度無關的結論：

因此，當銅線直徑相同時，不同螺旋角度銅線在復合圓柱殼結構中的截面積相同，建立的計算模型中粒子數目相近，因此可得到螺旋角度對復合圓柱殼影響不大的結論。

3.3 銅線直徑對穩定性的影響

由式(7) 可知復合圓柱殼銅線部分的截面積與銅線直徑有關，另外當銅線直徑發生變化時，相應地，復合圓柱殼各層的厚度也有可能發生改變，因此本節討論銅線直徑對復合圓柱殼的影響，主要討論0.25 和0.50 mm 兩種直徑，銅線螺旋角度均為15°。當銅線直徑為0.50 mm時，可以估算出將2 mm 的銅線層排滿需917 根銅線，這些銅線在銅線折返層只需排列一層，圖9為0.50 mm 銅線直徑下復合圓柱殼模型的結構參數及細節示意圖。使用前文提到的失穩判據，可以判斷對0.50 mm 銅線直徑的復合圓柱殼模型的失穩時間大致為23.62 μs，如圖10 所示，此時對應的壓縮率為66.41%，比0.25 mm 直徑的模型失穩時對應的壓縮率小。另外當選用0.50 mm直徑的銅線時，由于在復合圓柱殼只有一層銅線折返層，制備時極易出現折返層缺陷而導致裝置穩定性變差，圖11 為0.50 mm 銅線直徑下復合圓柱殼銅線折返層僅有一根銅線缺失時在23.20 μs 時的變形形態，可看到僅缺失一根銅線就會導致復合圓柱殼內壁在壓縮過程中有較大的環向擾動。

圖9 0.50 mm 銅線直徑下復合圓柱殼示意圖Fig. 9 Schematic diagrams of composite cylindrical shell under 0.50 mm copper wire diameter

圖10 0.50 mm 直徑下失穩時間判斷曲線Fig. 10 Judgment curve of 0.50 mm model instability time

圖11 t=23.20 μs 時含缺陷復合圓柱殼 變形Fig. 11 Deformation of defective cylindrical shell at t=23.20μs

4 結 論

本文利用AUTODYN 軟件建立了復合圓柱殼二維細節模型，基于SPH-FEM 耦合算法，對在炸藥爆轟沖擊壓縮下的復合圓柱殼的結構響應過程進行了數值模擬，并根據圓柱殼內壁粒子速度歷史提出了一種判斷失穩時間的方法。通過對金屬環氧套筒爆轟模型的計算，驗證了數值模型計算的準確性，另外由復合圓柱殼及金屬環氧套筒沿徑向不同位置的速度出現明顯差異的時間與使用提出的失穩判據判斷的失穩時間一致，證明了提出的失穩判據的合理性。

通過對不同數值模型的計算對比，可以發現：在復合圓柱殼制備過程中存在的折返層缺陷對其結構響應及穩定性影響較大，因此在制備過程中需要盡可能提高工藝，保證銅線折返層的完整性及均勻性；通過對0°、10°和15°等3 種螺旋角度復合圓柱殼模型的計算，可以發現復合圓柱殼的螺旋角度對其穩定性影響不大；而由于銅線直徑會直接影響復合圓柱殼各層材料的厚度，因此對穩定性影響較大，是進行裝置設計及實驗需重點考慮的參數。