不同細觀力學方法預測高聚物粘結炸藥有效模量的比較

王竟成， 羅景潤

(中國工程物理研究院 總體工程研究所， 四川 綿陽 621900)

不同細觀力學方法預測高聚物粘結炸藥有效模量的比較

王竟成， 羅景潤

(中國工程物理研究院 總體工程研究所， 四川 綿陽 621900)

高聚物粘結炸藥(PBX)不同于普通顆粒增強復合材料，其顆粒含量超過85%，組分彈性模量相差3～4個數量級，導致其有效模量的細觀力學理論預測出現很大偏差。結合有限元細觀模擬，對Mori-Tanaka法、自洽法、微分法3種細觀力學方法的預測結果進行了比較分析。結果表明：界限法上下界之間有量級上的差異；當顆粒含量小于10%，顆粒間相互作用較小，不同方法計算的有效模量差異不大；含量大于20%時，顆粒間相互作用增強，3種解析法預測的結果逐漸出現差異，而微分法與有限元結果比較接近；當顆粒含量為94.9%，微分法預測的PBX楊氏模量比實測值高3.7%，Mori-Tanaka法和自洽法結果都有量級上的偏差；對于顆粒含量高、組分性能反差大的復合材料，微分法較合理地計及了顆粒間的相互作用，能較準確地預測其有效模量。

兵器科學與技術； 細觀力學； 高聚物粘結炸藥； 有效模量； 有限元； 微分法

0 引言

高聚物粘結炸藥(PBX)廣泛運用于軍事領域，包含約90%的炸藥顆粒(奧克托今(HMX)、黑索今、梯恩梯等)，它們之間用高聚物粘接。細觀層次上可以將其劃歸為高密度顆粒填充復合材料，具有復雜的微細觀結構，分為基體相和顆粒相。PBX不僅是一個功能構件，在加工、運輸、儲存等過程中還要承載一定的外力載荷，其有效力學性能引起廣泛關注[1]，而有效模量[2-5]則是最基本的力學參量之一。由于PBX不同于普通顆粒增強復合材料，炸藥顆粒楊氏模量高于粘結劑基體3～4個數量級，泊松比相差近1倍，加之其顆粒含量極高，顆粒間相互作用強烈，這些都對基于細觀力學理論方法計算PBX的有效模量提出了巨大挑戰。經典細觀力學方法對PBX的適用性問題值得研究，而國內外也鮮有相關研究報道。

細觀力學利用多尺度的連續介質力學理論和方法，研究材料細觀結構和宏觀性能間的定量關系，預測材料的有效彈性模量、熱膨脹系數和強度等宏觀性能。細觀力學方法計算復合材料的有效模量已有半個多世紀的歷史，廣泛應用于各類復合材料：混凝土、金屬基復合材料、高聚物復合材料等，分為界限法和解析法。Voigt模型和Reuss模型分別考慮均勻應變和均勻應力邊界條件，根據最小勢能原理和最小余能原理得到有效剛度的V-R上下界；采用Hashin-Shtrikman (HS)變分原理適當地選擇均勻參考材料的剛度和極應力場，得到非均勻材料有效剛度的HS上下界。解析法主要以特征應變(Eshelby)問題的解答與等效夾雜理論為基礎，假設橢球形顆粒和基體均為各向同性材料，選擇不同的參考應變和參考剛度推導預測復合材料有效性能的公式，主要有稀疏法，Mori-Tanaka(MT)法，自洽(SC)法，微分(differential)法等[6]。

同時，經典的細觀力學方法在各個領域不斷地應用和發展，盧子興等[7]利用四相球模型對復合泡沫塑料的模量和屈服強度進行了理論預測；沈珉等[8]采用MT法與稀疏法結合內聚力模型對含非理想界面的顆粒增強復合材料的有效模量進行了理論估計；胡更開等[9]總結了細觀力學建立非均質材料有效性質的方法、主要結果和最新進展。對于顆粒含量極高、組分彈性模量相差幾個量級的PBX材料，其有效模量的細觀力學預測方法尚需研究。本文結合有限元法，研究了界限法和MT法、SC法、微分法3種經典細觀力學方法預測PBX有效模量的適用性。

1 細觀力學方法預測的對比

計算中炸藥顆粒和粘結劑均視為各向同性的彈性體，炸藥顆粒選用HMX，粘結劑選用聚酯型聚氨酯(Estane)。材料參數見表1，粘結劑基體和炸藥顆粒的彈性性能有巨大的差異，泊松比ν相差1倍，楊氏模量E更是相差4個數量級。通過(1)式可進一步計算出體積模量K和剪切模量G.

. (1) 表1 PBX組分彈性參量[10]Tab.1 Elastic properties of PBX ingredients

細觀力學方法分為上下界限法和解析法。由于HS變分原理的條件比最小勢能和最小余能原理強，(3)式HS上下界相較于(2)式V-R上下界給出了較精細的范圍。

(2)

(3)

式中：c表示顆粒的含量；K0與G0為粘結劑的體積模量與剪切模量；K1與G1為炸藥顆粒的體積模量與剪切模量。

圖1 HS上下界和V-R上下界的對比Fig.1 Comparison between HS bounds and V-R bounds

圖1給出了不同含量下這兩種界限法預測的體積模量K和剪切模量G. 從圖1可以看出，HS上下界在V-R上下界之間，因此在后續的對比中，僅參考HS上下界。由于組分的材料參數相差巨大，無論是基于串并聯模型的V-R界限法，還是通過選取不同組分的剛度作為參考剛度的HS界限法，上下界之間的差別都非常大，并且差異隨顆粒含量的增加不斷增大。因此，上下界限法不能用來預測PBX材料的有效模量，不過它們可以用來界定其他方法的合理性。

經典的解析法中，稀疏法沒有考慮夾雜之間的相互作用，僅適用于夾雜含量很低時有效模量的預測，其他細觀力學方法則通過不同程度地考慮夾雜間的相互作用，提高顆粒增強復合材料有效模量的估計值。(4)式MT法通過將遠處應力(應變)改為基體平均應力(應變)，考慮了顆粒之間的相互作用，通過將局部集中張量轉換為整體集中張量的方式提高有效模量。

(4)

(5)式SC法假定把夾雜物嵌于彈性模量未知的等效介質之中，通過提高參考剛度的方式提高有效模量的估計。

(5)

(6)式微分法沿用SC法的思想，將高顆粒含量復合材料的形成劃分為諸多稀疏顆粒投入等效介質的過程，一定程度上弱化了SC法中有效模量的提高。

(6)

MT法給出了復合材料有效剛度簡單明了的顯示形式，而SC法的隱式方程組和微分法的非線性微分方程組都需要借助數值方法才能求解。此外，還有基于三相球模型的廣義自洽法，由于其求解剪切模量過于冗長復雜的方程組，這里不加以列舉。

利用(7)式可進一步計算出楊氏模量和泊松比為

(7)

圖2顯示了不同細觀力學方法計算的有效模量，MT法與微分法的趨勢基本一致，接近HS的下限值，可以證明MT法和HS下限的結果完全相同。SC法在顆粒含量c=50%附近脫離下邊界，開始向上邊界靠近。當含量為90%時，SC法計算的K值是MT法和微分法預測值的11倍左右，G值更是高達上百倍，這主要是由組分間彈性性能的巨大差異導致。將圖2中低含量部分放大后(見圖3)發現：盡管不同方法體積模量K的預測值比較一致，但剪切模量G從c=20%附近即分散開來，且差異逐漸增大；其中，SC法預測值最大，MT法最小。

圖2 不同細觀力學方法計算的有效模量Fig.2 Effective moduli estimated by different methods

圖3 不同方法計算值的分離Fig.3 Separation of estimates among different methods

2 有限元細觀模型及計算方法

有限元是解決細觀力學問題主要的數值方法，可以更精確地描述材料的細觀幾何結構，計算更復雜的力學行為。Ananda等[10]采用內聚力模型研究了PBX多種顆粒形貌和粒徑分布的細觀熱力響應；戴開達等[11]通過建立炸藥顆粒圓形隨機模型和六邊形模型，研究了粘結劑和顆粒性質對PBX有效彈性模量的影響；Biswajit等[12]采用循環單元法(RCM)和傳統有限元法對PBX9501的有效模量進行了分析計算。

有限元細觀模擬中，通常以代表性體積單元(RVE)為基礎，將整個宏觀結構視為由許多這樣的體積單元按周期排布組合而成，通過施加均勻位移邊界條件，提取RVE的形變和反力等信息，計算出RVE的有效彈性模量等。

考慮到真實炸藥顆粒的平均尺寸在100 μm左右，選取RVE尺寸為1 000 μm×1 000 μm的方型區域[13-14]。由于圓形顆粒的填充度難以達到80%，基于Voronoi法構建不規則多邊形顆粒模型對PBX的有效模量進行計算。Voronoi模型(見圖4(a))可以使顆粒含量達到95%以上，并且與PBX真實細觀形貌(見圖4(b))有較好的相似性。通過改變粘結劑層的厚度，能方便地構建不同含量下的Voronoi模型。

圖4 Voronoi 模型和PBX真實細觀形貌Fig.4 Voronoi model and PBX microstructure

值得注意的是，有限元模型的構建是為了預測不同顆粒含量下有效模量的變化趨勢，僅作參考使用。Voronoi模型采用的是不規則多邊形顆粒，而經典的細觀力學方法假設的是圓球形夾雜，因此，它們的預測值必然有差別。低含量下顆粒形貌對復合材料的有效性能影響不大，但隨著顆粒含量的增加，顆粒形貌的影響便會凸顯。除了借助有限元方法來模擬材料的力學行為，相關計及顆粒形貌的細觀力學方法更值得不斷發展完善。

有限元計算采用文獻[11]中的方法，左邊界和下邊界分別約束住x向和y向位移，為保證邊界節點的位移與邊界線平行，在模型的上邊界和右邊界進行耦合約束；同時，在上邊界節點施加y向位移來模擬單軸壓縮試驗。同樣采用表1中的材料參數，使用平面應變假設和plane183單元進行有限元計算。根據平面應變假設和邊界條件：

(8)

可得

(9)

(10)

式中：εx、εy分別為x方向、y方向上的應變分量；σx、σy、σz分別為x方向、y方向、z方向上的應力分量；Ux、Uy分別為x方向、y方向上的位移分量；lx和ly分別為x方向和y方向RVE的邊長。

根據文獻[15]，平均應力既可以使用體積分得到，也可以使用面積分得到

(11)

式中：F代表作用反力；N為RVE邊界上節點的個數；V代表整個計算域；t為面力；x為節點間距。對于當前二維問題，有

(12)

式中：Fr為有限元后處理中提取的上邊界的節點反力之和；l=1 000 μm.

考慮(9)式、(10)式、(12)式有

(13)

利用(1)式進一步求解出K和G.

3 結果及分析

表2列出了c<50%時有限元模型與細觀力學方法的計算結果，有限元結果在圖3中用三角形表示。當顆粒含量小于10%，顆粒之間的相互作用很小，各細觀力學方法計算結果都具有很高的精度，而此時顆粒形貌的影響可以忽略，有限元結果與解析法計算結果吻合，體積模量計算誤差小于1%，剪切模量計算誤差小于3%，證明了有限元計算結果的可靠性。由于組分材料性能的巨大差異，當顆粒含量大于20%，各細觀方法計算的剪切模量出現明顯差異，且隨著含量的增加，它們與有限元結果的偏差越來越大。從圖3不難發現，微分法與有限元結果最接近，SC法與MT法的計算結果都有較大偏差。

表2 c<50%時有限元結果與細觀力學結果Tab.2 FE results and micromechanical results for c<50% MPa

表3列出了c>50%時有限元模型與細觀力學方法的計算結果，圖5直觀地給出了有效模量隨顆粒含量的變化趨勢。需要指出的是，由于SC法預測的有效模量太大，已不能在圖5中顯示，這也說明SC法不能用來預測PBX材料的有效模量。高含量下，微分法與MT法預測的體積模量與有限元結果很接近，但剪切模量僅有微分法與有限元結果較符合，SC法太大，MT法太小，與低含量的情況一致。隨著顆粒含量的增加，微分法預測值同有限元結果的差異變大，這主要由顆粒形貌的不同引起。微分法基于圓球形夾雜推導，而有限元結果則由多邊形顆粒模型計算。

表3 c>50%時有限元結果與細觀力學結果Tab.3 FE results and micromechanical results for c>50% MPa

圖5 細觀力學預測值與有限元結果的對比(c>50%)Fig.5 Comparison between micromechanical estimate and FE result (c>50%)

顆粒含量c=94.9%時，文獻[16]給出了室溫、準靜態條件下PBX的楊氏模量為956.9 MPa±85.2 MPa. 利用(7)式分別計算c=94.9%時各種方法的楊氏模量，MT法為36.4 MPa，SC法為22.7 GPa，微分法為992.6 MPa. 由此可見，在顆粒含量很高時，3種細觀力學方法的計算結果已有量級上的差異，這主要是因為炸藥顆粒和基體粘結劑的彈性性能存在巨大差異。與實驗數據相比，MT法低估了顆粒間的相互作用導致預測值太小，SC法則高估了顆粒間的相互作用導致預測值太大，微分法最接近(誤差為3.7%)。以上結果表明，在預測顆粒含量高、組分性能反差大的復合材料有效模量時，微分法這種無窮次均勻化的方法較合理地計及了顆粒間的相互作用，其有效模量預測值比較準確。微分法最初用于研究懸浮液的性能，后來進一步發展用來預測復合材料的有效性能。其無窮次投放顆粒的過程理論上可以使顆粒含量無限趨近100%，因此，在計算高顆粒含量PBX材料的有效模量時仍較準確。

另外，在微分法基礎上進一步引入顆粒形貌、尺寸分布等特征參量，能更精確地估測PBX的有效模量。

4 結論

1) PBX組分泊松比相差近2倍，楊氏模量相差4個數量級，界限法上下界之間出現量級上的差異，界限法不能有效預測PBX的有效模量。

2) 隨著顆粒含量的增加，不同細觀力學解析法預測的有效模量差異變大。當顆粒含量大于20%，微分法與有限元結果比較接近。

3) 微分法計算的PBX(c=94.9%)楊氏模量比實測值高3.7%，MT法和SC法都有量級上的偏差。因此，對于顆粒含量高、組分性能反差大的復合材料，微分法較合理地計及了顆粒間的相互作用，能較準確地預測其有效模量。

References)

[1] 羅景潤. PBX的損傷、斷裂及本構關系研究[D]. 綿陽：中國工程物理研究院，2001. LUO Jing-run. Study on damage, fracture and constitutive relation of PBX[D]. Mianyang: China Academy of Engineering Physics, 2001. (in Chinese)

[2] Biswajit B, Daniel O A. On predicting the effective elastic properties of polymer bonded explosives using the recursive cell methods[J]. International Journal of Solid and Structure, 2004, 41(2): 481-509.

[3] Biswajit B. Effective elastic moduli of PBX from finite element simulations[EB/OL].[2016-09-29].arXiv:cond-mat/0510367.

[4] 敬仕明，李明，龍新平. 基于改進的Hahin-Shtrikman方法預測PBX的有效模量[J]. 含能材料，2009, 17(6): 664-667. JING Shi-ming, LI Ming, LONG Xin-ping. On predicting effective elastic modulus of PBX by modified Hashin-Shtrikman model[J]. Chinese Journal of Energetic Materials, 2009, 17(6):664-667. (in Chinese)

[5] 敬仕明，李明，龍新平. PBX有效彈性性能研究進展[J]. 含能材料，2009, 17(1): 119-124. JING Shi-ming, LI Ming, LONG Xin-ping. Progress in predicting the effective elastic properties of PBX[J]. Chinese Journal of Energetic Materials, 2009, 17(1):119-124. (in Chinese)

[6] 黃克智，黃永剛.固體本構關系[M].北京：清華大學出版社，1999: 120-172. HUANG Ke-zhi, HUANG Yong-gang. Solid constitutive relations [M]. Beijing: Tsinghua University Press, 1999: 120-172. (in Chinese)

[7] 盧子興，鄒波. 復合泡沫塑料模量和屈服強度的理論預測[J].復合材料學報，2014, 31(4): 998-1005. LU Zi-xing, ZOU Bo. Theoretical prediction for modulus and yield strength of syntactic foams[J]. Acta Composite Material Sinica, 2014, 31(4): 998-1005. (in Chinese)

[8] 沈珉，郝培. 顆粒增強復合材料非理想界面剛度和有效模量的理論估計[J]. 復合材料學報，2016, 33(1): 189-197. SHEN Min, HAO Pei. Theoretical estimation of imperfection interfacial stiffness and effective modulus in particle reinforced composites[J]. Acta Composite Material Sinica, 2016, 33(1): 189-197. (in Chinese)

[9] 胡更開，鄭泉水，黃筑平. 復合材料有效彈性性質分析方法[J].力學進展，2001,31(3): 361-393. HU Geng-kai, ZHEN Quan-shui, HUANG Zhu-ping. Micromechanics methods for effective elastic properties of composite materials[J]. Advances in Mechanics, 2001, 31(3): 361-393. (in Chinese)

[10] Barua A, Zhou M. A Lagrangian framework for analyzing microstructural level response of polymer-bonded explosives[J]. Modeling and Simulation in Materials Science and Engineering, 2011, 19(5): 1-24.

[11] 戴開達，劉龑龍，陳鵬萬. PBX炸藥有效彈性模量的有限元模擬[J].北京理工大學學報，2012, 32(11): 1154-1158. DAI Kai-da, LIU Yan-long, CHEN Peng-wan. Finite element simulation on effective elastic modulus of PBX explosives[J]. Transactions of Beijing Institute of Technology, 2012, 32(11): 1154-1158. (in Chinese)

[12] Biswajit B, Daniel O A. Micromechanics based determination of effective elastic properties of polymer bonded explosives[J]. Physica B, 2003, 338(1/2/3/4): 8-15.

[13] Lemaitre J. Formulation and identification of damage kinetic constitutive equations[J]. Continuum Damage Mechanics Theory and Applications, 1987, 295: 37-89.

[14] Ostoja S M. Random field models of heterogeneous materials [J].Solids and Structures, 1998, 35(19): 2429-2455.

[15] Nguyen V D, Bechet E, Geuzaine C, et al. Imposing periodic boundary condition on arbitrary meshes by polynomial interpolation[J]. Computational Material Science, 2012, 55(5): 390-406.

[16] Gray G T, Idar D J, et al. High- and low-strain rate compression properties of several energetic material composites as a function of strain rate and temperature [C]∥11th international Detonation Symposium. Snowmass, Colorado, US: Los Alamos National Laboratory, 1998.

Comparison of Effective Moduli of Polymer Bonded ExplosivePredicted by Different Micromechanical Methods

WANG Jing-cheng, LUO Jing-run

(Institute of Systems Engineering, China Academy of Engineering Physics, Mianyang 621900, Sichuan, China)

Polymer bonded explosives (PBXs) are different from conventional particle reinforced composites. PBXs possess over 85% of particle volume fraction, and the elastic modulus of particle can be three or four orders of magnitude higher than that of binder. Large deviations are observed on the effective moduli predicted by micromechanical methods. The predicted results of several micromechanical methods for the effective moduli of PBX are analyzed based on finite element method. Results indicate that both Voigt-Reuss bounds and Hashin-Shtrikman bounds give rather large estimated regions; when particle fraction is less than 10%, the differences among the effective moduli estimated by the different methods are not obvious since the particle interaction is negligible; when particle fraction exceeds 20%, the estimated results of three analytic methods are different from each other gradually due to the enhancement in particle interactions, and the estimated results of differential method agree with those of finite element method. For the particle fraction of 94.9%, Young’s modulus calculated by differential method is only 3.7% higher than the experimental value of PBX, while a huge error occurs using Mori-Tanaka method or self-consistent method. For the composites that possess high particle fraction and strong contrast in their ingredient properties, the differential method provides better estimates of effective modulus owing to its relatively favorable evaluation on particle interactions.

ordnance science and technology; micromechanics; PBX; effective modulus; finite element; differential method

2016-10-19

國家自然科學基金項目(11472257)

王竟成(1991—)，男，博士研究生。E-mail：wangjczzy@sina.cn

羅景潤(1966—)，男，研究員，博士生導師。E-mail：jrluo2000@vip.sina.com

TB332; O343.1

A

1000-1093(2017)06-1106-07

10.3969/j.issn.1000-1093.2017.06.009