復合固體推進劑“脫濕”點的率相關特性研究

摘 要：復合固體推進劑顆粒和基體之間的界面“脫濕”是推進劑在加載過程中一種重要的損傷形式，“脫濕”點的強度和伸長率均與載荷相關。本文通過某推進劑在低溫、不同圍壓和應變率下加載時的應力-應變響應曲線特征，發現推進劑的“脫濕”點表現出明顯的率相關和圍壓相關性，圍壓載荷會抑制推進劑的“脫濕”，尤其是對于高應變率下具有“脫濕”點峰值的應力-應變曲線影響最大。針對推進劑“脫濕”點的率相關特征，借助分子動力學思想創建復合推進劑的細觀顆粒夾雜幾何模型，建立推進劑有限元計算細觀模型對推進劑界面在不同應變率下的“脫濕”損傷進行了仿真計算。計算結果揭示了推進劑在不同應變率下加載時的“脫濕”損傷機理，通過應變率對“脫濕”點的影響規律闡釋了推進劑在低溫、高應變率下伸長率驟降的主要原因，發現推進劑的破壞取決于損傷界面數量和界面損傷程度在損傷演化過程中的交互作用，是一個動態發展的過程。

關鍵詞：復合固體推進劑； “脫濕”點； 應變率； 細觀模型； 數值計算

中圖分類號： TJ760; V43

文獻標識碼： A

文章編號：1673-5048（2023）05-0086-06

DOI： 10.12132/ISSN.1673-5048.2023.0073

0 引" 言

復合固體推進劑是一種由氧化劑、金屬燃燒劑和黏合劑等多種組分以一定的配比混合而成的固體推進劑。其中，黏合劑充當基體即連續相，其力學響應依賴于時間、溫度和加載歷史等多種因素。氧化劑和其他固體填料為分散相，用來提高推進劑的剛度和強度，因此復合推進劑是一種結構不均勻的顆粒增強復合材料。推進劑的分散相表面通過吸附作用與連續相在生產過程中混合潤濕，而這種依靠吸附形成的黏合力并不牢靠，在受載過程中界面之間極易產生分離，這種現象被稱為“脫濕”。“脫濕”會大大削弱固體填料對于基體的強化作用，進而導致推進劑體積膨脹，其拉伸應力-應變曲線也表現出很強的非線性。在固體發動機點火增壓的過程中，極易由于“脫濕”造成藥柱承載能力下降，藥柱出現裂紋導致燃面大幅增加，使燃燒室壓力急升而引起發動機爆炸。

“脫濕”損傷問題的深入研究需要借助細觀試驗和細觀力學方法，研究推進劑內部各細觀組分的力學特性和受載過程中的相互作用，掃描電鏡［1-3］和微CT掃描［4-5］是目前較為常用的細觀試驗手段。此外，掃描電鏡還可以結合光學顯微鏡［6］和CT無損檢測［7］，分析推進劑內部的裂紋擴展和損傷情況。細觀試驗方法雖然簡單直觀，但是需要借助專業的試驗儀器和設備，推進劑從細觀尺度來看是一種典型的顆粒夾雜結構，通過細觀力學仿真方法，建立推進劑的細觀模型，針對不同的工況進行推進劑細觀損傷仿真也是一種行之有效的方法。顆粒夾雜模型建模過程中通常采用的顆粒填充算法有順序算法［8-10］和并發算法［11-13］兩種，分子動力學法和蒙特卡洛法又是應用較廣的兩種并發算法且分子動力學法在計算效率上有一定優勢，由于并發算法以隨機方式填充顆粒，因此更適用于推進劑的細觀建模。

內聚力模型常用來定義推進劑內部各組分之間的界面失效，內聚力模型采用統一的數學模型來描述裂紋從出現到擴展的整個過程，其數學模型是連續的且沒有奇異性的。雙線性內聚力模型［14］結構形式簡單，是目前研究人員應用最多的內聚力模型。隨后不少學者基于雙線性模型引申出更為復雜的內聚力模型，如Hou等引入了分段函數式內聚力模型，用來計算雙基推進劑 ［15］和復合推進劑［16］的裂紋擴展和界面“脫濕”損傷問題。

本研究以某高固體含量復合推進劑為對象，通過推進劑在低溫下以不同圍壓和應變率加載時的應力-應變響應特征，找到“脫濕”點的率相關變化特征。借助分子動力學思想建立推進劑細觀顆粒夾雜模型，并通過有限元軟件二次開發生成推進劑細觀有限元計算模型，在顆粒和基體界面上引入內聚力模型，最終形成推進劑的細觀有限元損傷計算模型，基于該模型對推進劑在不同應變率下的損傷過程進行仿真計算，分析了推進劑的率相關“脫濕”損傷特征，探究了推進劑“脫濕”點與應變率的關系，揭示了推進劑的率相關“脫濕”損傷機理，闡釋了推進劑在低溫、高應變率下伸長率下降的主要原因。

1 試驗方案和試驗結果

1.1 試驗方案

復合推進劑由于基體的粘彈特性，具有很強的率相關性、溫度相關性和壓力相關性。在拉伸試驗中應變率、溫度和環境壓力是影響拉伸曲線形狀特征最重要的3個個因素。常常將低溫點火（低溫+高應變率+圍壓）作為戰術導彈發動機考核的重要試驗之一。

針對戰術導彈的低溫點火極端考核工況，本研究選擇在-55 ℃下開展推進劑在不同應變率和圍壓下的單軸定速拉伸試驗。選取5個拉伸速率100 mm/min、500 mm/min、1 000 mm/min、2 100 mm/min、4 200 mm/min進行單軸拉伸試驗，試驗采用標準啞鈴形試件（長度120 mm，標距70 mm），5個拉伸速率根據試件的標距折算的應變率分別為2.38 %/s、11.9 %/s、23.8 %/s、50 %/s和100 %/s，選取4個圍壓0.1 MPa（常壓）、3 MPa、6 MPa和10 MPa，開展推進劑在低溫、不同圍壓和應變率下的力學性能試驗研究。試驗最大拉伸速率是6 000 mm/min，加載過程中溫度變化范圍從低溫-80 ℃到高溫+180 ℃，該試驗還可加載圍壓，最大圍壓可達15 MPa，加載過程中兩個安裝在對稱夾頭上的試件同時拉伸得到兩組試驗數據，如果試驗數據一致性不好，再補充2～4個平行試樣，整個試驗系統如圖1所示。

1.2 試驗結果

在常壓、低溫-55 ℃和不同圍壓下，推進劑在不同應變率下單軸定速拉伸的應力-應變曲線如圖2所示。在常壓下，如圖2（a）所示，隨著應變率的增加，推進劑的最大伸長率呈下降趨勢，尤其是當應變率為100 %/s時，推進劑的最大伸長率下降到10%以內。從試驗曲線看，當應變率超過50 %/s時，應力-應變曲線上在“脫濕”點附近出現峰值，推進劑在“脫濕”點附近強度達到最大，因此應力-應變曲線上并未出現推進劑在較低應變率下拉伸時“脫濕”后的線性強化段，而是在經過“脫濕”點峰值后很快破壞。當圍壓增加到3 MPa后，如圖2（b）所示，當應變率超過50 %/s時，“脫濕”點峰值消失，說明圍壓在一定程度上抑制了推進劑的“脫濕”。從總體來看，隨著應變率的增加，推進劑的最大伸長率仍呈下降趨勢。當圍壓增加到6 MPa時，如圖2（c）所示，推進劑在不同應變率下的拉伸曲線均看不出明顯的“脫濕”點，在高應變率下推進劑的拉伸曲線變得更加平滑，由于圍壓抑制了推進劑的“脫濕”，高應變率下推進劑的最大伸長率顯著增加，且推進劑的最大伸長率并未隨應變率的增加單調下降，應變率23.8 %/s下推進劑最大伸長率出現極小值。隨著圍壓的繼續增大，當圍壓達到10 MPa時，如圖2（d）所示，與6 MPa圍壓下推進劑的拉伸曲線類似，推進劑在不同應變率下的拉伸曲線均沒有明顯“脫濕”點，且最大伸長率的極小值出現在23.8%應變率下，但是應變率100 %/s下推進劑的最大伸長率超過了50 %/s下的最大伸長率。

從試驗數據來看，圍壓會抑制推進劑的“脫濕”，推進劑在圍壓下進行單軸拉伸時“脫濕”點不明顯［17］，尤其是對于高應變率下具有“脫濕”點峰值的應力-應變曲線影響最大，圍壓載荷的增加不僅會導致推進劑的“脫濕”點峰值消失，而且會使推進劑的最大伸長率顯著增加，且推進劑的“脫濕”點表現出明顯的率相關性，本研究將針對推進劑的“脫濕”點率相關性展開。

1.3 “脫濕”點伸長率

典型的復合推進劑以較低拉伸速率，單軸拉伸的應力-應變曲線，如圖3（a）所示。曲線特征主要包括線彈性階段、“脫濕”損傷段、線性強化段和拉伸斷裂段等4個階段。用K1表示線彈性階段的斜率，用K2表示線性強化段的斜率，工業部門通常會將K2和K1的交點所形成的鈍角的角平分線與應力-應變曲線的交點C定義為“脫濕”點。對于“脫濕”點附近存在峰值的應力-應變曲線，如圖3（b）所示。由于在應力-應變曲線上出現了除抗拉強度以外的應力極值點，如果通過角平分線法來確定“脫濕”點，“脫濕”點的位置很可能并非出現在曲線上斜率變化最大的位置，如圖3（b）中的點C1所示。為了便于比較“脫濕”點位置的率相關性，針對出現明顯“脫濕”點峰值的應力-應變曲線，將該曲線上出現的第一個極值點C定義為“脫濕”點，如圖3（b）所示。

采用如圖3所示的方法確定“脫濕”點的位置，通過該推進劑在低溫、常壓下的率相關應力-應變曲線（圖2（a））提取不同應變率下推進劑的“脫濕”點伸長率，在應力-應變曲線上將“脫濕”點位置通過紅色圓點高亮標注如圖4（a）所示，“脫濕”點伸長率變化趨勢如圖4（b）所示。整體來看，隨著應變率的增大，“脫濕”點伸長率隨應變率的增加呈先減小后增大的趨勢，在某個應變率下推進劑進行單軸拉伸時最容易出現“脫濕”。

2 “脫濕”細觀模型和計算結果

推進劑是一種典型的粘彈性材料，溫度和時間（應變率）對推進劑力學性質的影響可以近似等效。因此，本文主要針對推進劑“脫濕”點的率相關性展開研究。首先，基于分子動力學思想，獲取顆粒粒徑信息和代表性體積單元的尺寸數據等，并借助有限元軟件的二次開發功能，創建顆粒和基體的原始模型； 其次， 自動生成顆粒夾雜有限元計算幾何模型； 再次，通過顆粒和基體的材料參數定義、接觸對定義、損傷設置、載荷和邊界條件設置和網格模型創建等步驟，生成顆粒夾雜有限元計算模型； 最后，基于創建的有限元模型進行有限元分析和結果后處理，分析推進劑的損傷演化情況，闡釋“脫濕”點特征與應變率的相關性。

2.1 模型建立過程

本研究所采用的推進劑AP顆粒質量分數約為70%，由于AP的平均粒徑遠大于Al粒子，在建模過程中僅將AP粒子建立為顆粒模型，將Al粒子和其他組分如丁羥膠等均歸為基體材料。采用的推進劑級配如圖5（a）所示，創建的有限元細觀損傷計算模型如圖5（b）所示，該代表性體積單元的邊長為1.1 mm。

AP顆粒被認為是完全彈性的，彈性模量取32 450 MPa，泊松比取0.143 3?？紤]基體的粘彈性，松弛模量采用Prony級數的形式表示：

E（t）=E0－∑ni=1Ei（1－e－tτi）

Prony級數的參數可由推進劑基體的應力松弛試驗數據擬合獲得，如圖6所示。

擬合得到推進劑的松弛模量數據如表1所示。 該基體的初始模量E0=6.629 4 MPa，泊松比μ=0.499 5。

本文采用雙線性內聚力模型定義顆粒和基體界面的損傷，如圖7所示。表征損傷的3個參數分別為初始剛度K0、破壞距離δf和臨界應力σ。通過參數反演得到的界面損傷參數如表2所示。

2.2 “脫濕”損傷計算結果及分析

為了比較不同應變率下推進劑的損傷演化過程，計算中所采用的應變率分別為： 2.5 %/s、12 %/s、50 %/s和100 %/s。以應變率為100 %/s的工況為例，推進劑代表性體積單元拉伸過程中的應力-應變曲線如圖8所示，從“脫濕”點對應的Mises應力云圖可以看出，在“脫濕”點附近，在一些大顆粒和基體的界面上已經出現了明顯的顆?！懊摑瘛?，然后觀察界面的接觸狀態，讀取接觸的界面上從面節點與主面的距離。設置一個大于0的容差可以觀察接觸界面的分離情況，本次提取的界面是從面節點與主面距離大于1 μm的顆粒與基體的接觸界面，由圖8可以看出，在“脫濕”點附近，大部分的顆粒都沿著拉伸方向出現了顆粒與界面的“脫濕”分離。

提取2.5 %/s、12 %/s、50 %/s和100 %/s等4個應變率下單軸拉伸時“脫濕”點對應的損傷界面進行進一步的對比研究。如圖9所示，在各應變率下均提取從面節點與主面距離超過1 μm的界面。比較4個不同應變率下單軸拉伸時的界面損傷情況，從界面損傷程度來看， 2.5 %/s應變率下從面節點與主面分離位移更大，在同樣的應變下“脫濕”損傷程度更高。從損傷區域的范圍來看，100 %/s應變率下界面“脫濕”損傷區域的數量明顯大于在2.5 %/s應變率下的損傷區域數量。由于計算過程中并未定義基體損傷，如果再考慮高應變率下基體的損傷，在應變率為100 %/s時推進劑內部損傷區域的數量將遠大于2.5 %/s應變率時損傷區域的數量。整體來看，隨著應變率的增加，損傷界面的數量逐漸增加，但界面的損傷程度逐漸下降。高應變率下應力-應變曲線上出現的“脫濕”點峰值與損傷界面數量急劇增加有關。

通過分析推進劑在不同應變率下拉伸時的損傷界面特征，發現推進劑的“脫濕”與損傷界面數量和界面損傷程度均相關。在低應變率下推進劑的顆粒“脫濕”是從界面產生微孔洞開始，在拉伸過程中損傷界面逐步擴大，“脫濕”區域的數量并不多且“脫濕”以大顆粒周圍界面的局部損傷演化為主。以較高應變率拉伸時，推進劑內部界面的損傷演化顯著加快，將會出現“脫濕”點逐漸提前的情況。隨著應變率的繼續增大，界面的損傷區域逐漸放大，在更多的顆粒和基體界面上出現了損傷，因此也有更多的界面去分攤載荷，在一定程度上抑制了推進劑局部的“脫濕”損傷演化，進而又出現了“脫濕”點相對滯后的情況。因此，增大應變率并不意味著一定會造成推進劑“脫濕”點靠前，而是在拉伸過程中一直伴隨著損傷程度加劇和損傷界面數量增加的博弈，推進劑的破壞取決于損傷界面數量和界面損傷程度在損傷演化過程中的交互作用，是一個動態發展的過程。此外，由于圍壓會抑制推進劑的“脫濕”，對高應變率下推進劑的最大伸長率影響也最大，在常壓、高應變率下出現的“脫濕”點峰值隨著圍壓的增加逐漸消失，應力-應變曲線上沒有明顯的“脫濕”點，且“脫濕”點伸長率最小的推進劑在較高圍壓下最大伸長率也最低。

3 結" 論

（1） 推進劑以較為寬泛的應變率進行單軸拉伸時存在某個應變率下推進劑的“脫濕”點伸長率最小，且在一定圍壓載荷下，“脫濕”點伸長率最小的推進劑，其最大伸長率也最小。

（2） 圍壓對高應變率下推進劑的最大伸長率影響最大，在常壓、高應變率下出現的“脫濕”點峰值隨著圍壓的增加逐漸消失，且隨著圍壓載荷的增加，推進劑的最大伸長率顯著增加，在較高圍壓下，“脫濕”點伸長率較大的推進劑，其可以達到的最大伸長率也較大。

（3） 推進劑在拉伸過程中一直伴隨著損傷程度加劇和損傷界面數量增加的博弈，推進劑的破壞取決于損傷界面數量和界面損傷程度在損傷演化過程中的交互作用，是一個動態發展的過程。

參考文獻：

［1］ 楊秋秋， 蔡如琳， 徐勝良， 等. 原位拉伸掃描電鏡法研究GAP推進劑的損傷行為［J］. 火炸藥學報， 2019， 42（5）： 511-515.

Yang Qiuqiu， Cai Rulin， Xu Shengliang， et al. Damage Behavior of GAP Solid Propellant by In-Situ Tensile SEM Method［J］. Chinese Journal of Explosives amp; Propellants， 2019， 42（5）： 511-515.（in Chinese）

［2］ 周水平， 唐根， 龐愛民， 等. GAP復合固體推進劑細觀結構演變特性［J］. 固體火箭技術， 2018， 41（5）： 580-585.

Zhou Shuiping， Tang Gen， Pang Aimin， et al. Evolution of Meso-Structure of GAP Composite Solid Propellants［J］. Journal of Solid Rocket Technology， 2018， 41（5）： 580-585.（in Chinese）

［3］ 伍鵬， 李高春， 王鑫. HTPB推進劑三點彎曲過程試驗與數值模擬［J］. 含能材料， 2020， 28（6）： 514-521.

Wu Peng， Li Gaochun， Wang Xin. Experimental and Numerical Simulation on the Damage Process of HTPB Propellant at the Crack Tip［J］. Chinese Journal of Energetic Materials， 2020， 28（6）： 514-521.（in Chinese）

［4］ Collins B， Maggi F， Matous K， et al. Using Tomography to Characterize Heterogeneous Propellants［C］∥46th AIAA Aerospace Sciences Meeting and Exhibit， 2008.

［5］ 劉新國， 劉佩進， 強洪夫， 等. 基于微CT技術的丁羥推進劑脫濕定量表征方法研究［J］. 推進技術， 2019， 40（5）： 1162-1168.

Liu Xinguo， Liu Peijin， Qiang Hongfu， et al. Quantificational Method of Dewetting Damage of HTPB Propellant Based on Micro-CT Detection［J］. Journal of Propulsion Technology， 2019， 40（5）： 1162-1168.（in Chinese）

［6］ 汪文強， 鄭健， 許進升， 等. AP/Al/CMDB推進劑的斷裂特性實驗研究［J］. 推進技術， 2015， 36（11）： 1728-1733.

Wang Wenqiang， Zheng Jian， Xu Jinsheng， et al. Research on Fracture Mechanism of AP/Al/CMDB Propellant［J］. Journal of Propulsion Technology， 2015， 36（11）： 1728-1733.（in Chinese）

［7］ 陳向東， 常新龍， 賴建偉， 等. 不同狀態HTPB推進劑動態沖擊損傷研究［J］. 推進技術， 2020， 41（7）： 1649-1659.

Chen Xiangdong， Chang Xinlong， Lai Jianwei， et al. Research on Dynamic Impact Damage of HTPB Propellant under Different States［J］. Journal of Propulsion Technology， 2020， 41（7）： 1649-1659.（in Chinese）

［8］ Widom B. Random Sequential Addition of Hard Spheres to a Vo-lume［J］. The Journal of Chemical Physics， 1966， 44（10）： 3888-3894.

［9］ Visscher W M， Bolsterli M. Random Packing of Equal and Unequal Spheres in Two and Three Dimensions［J］. Nature， 1972， 239（5374）： 504-507.

［10］ Sherwood J D. Packing of Spheroids in Three-Dimensional Space by Random Sequential Addition［J］. Journal of Physics A： Mathematical and General， 1997， 30（24）： L839-L843.

［11］ Lubachevsky B D， Stillinger F H， Pinson E N. Disks vs. Spheres： Contrasting Properties of Random Packings［J］. Journal of Statistical Physics， 1991， 64（3/4）： 501-524.

［12］ Lubachevsky B D， Stillinger F H. Geometric Properties of Random Disk Packings［J］. Journal of Statistical Physics， 1990， 60（5/6）： 561-583.

［13］ Speedy R J. Random Jammed Packings of Hard Discs and Spheres［J］. Journal of Physics： Condensed Matter， 1998， 10（19）： 4185-4194.

［14］ Hillerborg A， Modéer M， Petersson P E. Analysis of Crack Formation and Crack Growth in Concrete by Means of Fracture Mechanics and Finite Elements［J］. Cement and Concrete Research， 1976， 6（6）： 773-781.

［15］ Hou Y F， Xu J S， Zhou C S， et al. Microstructural Simulations of Debonding， Nucleation， and Crack Propagation in an HMX-MDB Propellant［J］. Materials amp; Design， 2021， 207： 109854.

［16］ Hou Y F， Xu J S， Zhou C S， et al. Microstructural Simulation of Initial Particle-Matrix Interfacial Defect in HTPB Solid Propellants［J］. Journal of Physics： Conference Series， 2021， 2011（1）： 012008.

［17］ 張亮， 申志彬， 虞跨海. 圍壓下固體推進劑的破壞機理分析［J］. 國防科技大學學報， 2022， 44（6）： 6-13.

Zhang Liang， Shen Zhibin， Yu Kuahai. Failure Mechanism Ana-lysis of Solid Propellant under Confining Pressure［J］. Journal of National University of Defense Technology， 2022， 44（6）： 6-13.（in Chinese）

Study on Rate-Dependent Characteristics of “Dewetting”

Point of Composite Solid Propellant

Shen Xin1，Zhang Liang2*， Dong Meng1， Chen Jing3

（1. The First Military Representative Office of Air Force Equipment Department in Luoyang， Luoyang 471009， China;

2. Henan University of Science and Technology， Luoyang 471003， China; 3. Unit 93160 of PLA， Beijing 100076， China）

Abstract： “Dewetting” can happen between solid particles and matrix of solid propellant during loading process. The strength and elongation at “dewetting” point are directly related to loads. The study based on the stress-strain response curve characteristics of a propellant under low temperature， different superimposed pressures and strain rates， the “dewetting” point of the propellant has shown obvious rate and superimposed pressure correlation， and the superimposed pressure load would inhibit “dewetting”， especially can influence the stress-strain curve with peak value of “dewetting” point at high strain rate. In order to obtain the rate-related characteristics of the “dewetting” point of the propellant， the mesoscale geometric particle inclusion model of the composite propellant is obtained based on the idea of molecular dynamics method， the mesoscale finite element model is created to calculate the “dewetting” damage of solid propellant under different strain rates. Finally， the damage mechanism of propellant “dewetting” under different strain rates is revealed， and the main reason of abruptly decrease of elongation of propellant at low temperature and high strain rate is explained by the influence rule of strain rate on “dewetting” point. It is found that the propellant failure depends on the interaction of the number of damage interfaces and the degree of interface damage in the process of damage evolution， which is a dynamic process.

Key words： composite solid propellant； “dewetting” point； strain rate； mesoscale model； numerical simulation