HTPB推進劑中增塑劑擴散系數計算①

李紅霞，強洪夫，李新其，王洪偉

(1.第二炮兵工程大學，西安 710025;2.中國航天科工集團公司六院四十六所，呼和浩特 010076)

HTPB推進劑中增塑劑擴散系數計算①

李紅霞1，強洪夫1，李新其1，王洪偉2

(1.第二炮兵工程大學，西安 710025;2.中國航天科工集團公司六院四十六所，呼和浩特 010076)

針對丁羥推進劑/襯層界面增塑劑的遷移問題，研究增塑劑的擴散特性。采用分子動力學方法模擬增塑劑在丁羥粘合劑體系中的運動，再通過愛因斯坦關系式求得擴散系數;采用對襯層加厚的粘接試件進行加速老化實驗方法，再通過費克第二定律計算得到增塑劑的擴散系數;考察不同環境溫度和不同增塑劑含量條件下癸二酸二辛酯的擴散系數的變化。分析認為，增塑劑在丁羥粘合體系高分子鏈段運動產生的間隙中，以“跳躍”方式發生空間位置遷移。結果表明，擴散系數模擬值和實驗值基本一致，數量級為10-12m2/s;溫度升高，增塑劑運動活性加強，有效活動空間增大，擴散系數增大;粘合體系與增塑劑的共容限量使得增塑劑含量大于3%時，模擬得到的擴散系數依次略有下降。分子動力學方法計算增塑劑擴散系數更具優勢。

丁羥推進劑;增塑劑;擴散系數;實驗;模擬

0 引言

目前，對于固體推進劑組分擴散遷移研究，多以現象分析和實驗為主。如徐思羽［1］分析了二茂鐵衍生物在固體推進劑中的遷移問題;尹華麗［2-3］用加速老化試驗研究了NEPE推進劑界面遷移組分含量的變化，還用浸泡增重法和氣相色譜儀考察了HTPB襯層、絕熱層對DOS和T27的吸收以及DOS、T27和GFP的遷移，他們的結果對后續研究工作具有重要參考價值。由于推進劑成分多、貯存時間長、反應復雜，實驗中簡化研究對象、設置極端條件，往往只能得到定性的結果，很少揭示其遷移機理。而分子模擬方法既能省時省力省錢達到實驗的目的，還能提供物質的微觀信息，為機理分析提供理論依據，現已成功用于推進劑領域［4-5］。筆者曾嘗試過推進劑中增塑劑擴散行為的分子模擬工作［6-7］，盡管當時缺少擴散系數的試驗驗證，但分子模擬方法在其他領域的成功應用，展示出了該方法的強大優勢。

擴散系數是物質的一種很重要的輸運系數，很難用實驗方法精確測量和確定，已有很多學者通過分子模擬手段來計算擴散系數或研究擴散過程，Kotelyanskii［8］等采用分子動力學模擬方法計算了水和鹽在反滲透膜中的擴散過程和擴散系數，對膜的截鹽效果做了微觀的解釋，Hofmann［9］等模擬了 H2、O2在無定形聚酰胺等幾種高分子聚合物膜中的擴散過程，并計算了擴散系數。擴散系數的模擬結果得到了實驗驗證［10-12］，如乙醇在聚硅氧烷(PDMS)中擴散系數的分子模擬結果4.6×10-6cm2/s與實驗結果4.5×10-6cm2/s接近;用分子動力學模擬得到的甲基丙烯酸甲酯寡聚物在高分子量的聚甲基丙烯酸甲酯母體中擴散系數的模擬結果，與實驗得到的數值在趨勢上吻合得很好。

增塑劑的擴散特性直接影響丁羥推進劑/襯層界面增塑劑的遷移速度，其遷移速度主要取決于增塑劑在推進劑的含量和擴散系數。本文主要開展HTPB推進劑中增塑劑擴散的分子模擬研究和擴散系數的實驗測定工作，旨在分析擴散機理，以期為界面增塑劑遷移問題的解決提供依據。

1 增塑劑擴散系數的分子模擬計算

1.1 模擬方法

分子模擬的總體趨勢是向著基于量子理論的方向發展，但流體的擴散性質研究一般不涉及到電子的轉移，采用經驗的位能函數就能達到精度要求，故不需采用費時的量子力學方法。相對于Monte Carlo模擬，分子動力學方法可更容易得到擴散性質，本文采用基于分子力學的分子動力學方法。

分子動力學模擬與真實實驗方法具有相似的過程，首先要選擇一個模擬體系，即樣品的準備，然后在一定的模擬條件(如，恒溫、恒壓等)下未解一系列運動方程，等待體系平衡后，開始測定物性性質及進行數據分析。分子動力學模擬是一種“計算實驗方法”，用來計算經典多體體系的平衡和傳遞性質。它通過求解體系中粒子的運動方程獲得粒子的運動速度和運動軌跡，并按照統計物理和統計熱力學原理得出系統相應的宏觀物理特性。在分子動力學模擬中，擴散性質通過時間平均獲得，即通過模擬中的均方位移MSD隨時間變化的Einstein關系［13］，求得擴散系數D。

式中r(t)為t時刻分子的坐標;r(0)為初始坐標;〈│r(t)-r(0)│2〉為均方位移。

擴散系數既是材料的屬性，又是溫度的函數，本文將用分子動力學方法計算不同環境溫度和增塑劑含量條件下的增塑劑癸二酸二辛酯(DOS)擴散系數。

1.2 模擬內容

整個模擬工作均在分子模擬軟件Materials Studio 4.3上進行，全部計算在高性能計算平臺IBM Blade上完成。用Visualizer模塊構建DOS、IPDI、HTPB分子模型，選用COMPASS力場，用分子力學法進行幾何結構優化。用Amorphous Cell模塊構造由10個DOS、5個IPDI、5個HTPB分子組成的混合物(共3 566個原子)模型，其體系密度 ρ=0.91 g/cm3，胞元模型的a、b、c都是3.29 nm。對混合物結構進行幾何優化，對優化過的結構進行弛豫。優化方法采用Smart Minimization方法，分別用atom-based和Ewald方法求范德華作用和靜電作用。模擬過程中位能采用球形截斷法，截斷半徑 cutoff取0.95 nm，spline width 取 0.1 nm，緩沖寬度buffer width取0.05 nm，時間步長為1 fs，截斷距離之外的分子間作用按平均密度近似方法進行校正。在真空條件下，采用Andersen控溫方法和Berendsen控壓方法，各分子起始速度按Maxwell-Boltzmann分布取樣，MD模擬積分方法為Velocity Verlet算法。

在Discover模塊下，選擇恒溫恒壓系綜(NPT)，在1標準大氣壓下，增塑劑含量為3%時，分別在273、298、310、323、348 K下進行MD模擬。達到平衡后，在相同條件下再進行1 000 ps的MD模擬，保存全運動軌跡，分析不同溫度下DOS在混合體系中的均方位移。根據增塑劑在推進劑中的含量為3%、6%、9%、15%，溫度為298 K，構建由DOS、IPDI、HTPB 不同分子數組成的混合物模型。用同樣方法分析不同含量的DOS在混合體系中的均方位移。圖1是溫度298 K、增塑劑3%含量下DOS的MSD-t曲線(略去其他MSD-t曲線)。

1.3 擴散系數的模擬結果

由均方位移隨時間的變化曲線，經擬合得到各直線斜率，計算得到 273、298、310、323、348 K 下 DOS 在粘接體系中的擴散系數分別為1.67×10-12、3.33×10-12、4.17 × 10-12、5.17 × 10-12、7.17 × 10-12m2/s，DOS在粘接體系中的擴散系數隨溫度的升高而逐漸增大。計算得到增塑劑在推進劑中的含量為3%、6%、9%、15%時的擴散系數分別為 3.33 ×10-12、3.00 ×10-12、2.33 × 10-12、2.00 × 10-12m2/s，隨增塑劑含量增大，擴散系數反而下降。

圖1 增塑劑DOS的MSD-t曲線(298 K、含量3%)Fig.1 Curve of MSD-t of plasticizer DOS(298 K，content 3%)

2 增塑劑擴散系數的加速老化實驗測定

2.1 實驗內容

固體火箭研制定型，需要評定其貯存壽命，在沒有充分的現場試驗結果時，通過加速貯存試驗預測固體火箭使用壽命是一個重要的途徑。本文采用對粘接試件進行加速老化試驗方法，測定襯層中增塑劑的含量。按照QJ 2038.1A—2004矩形件制作方法制作丁羥推進劑/襯層/絕熱層/殼體模擬粘接試件，增塑劑含量為3%、襯層配方與模擬部分中丁羥粘合系統一致(HTPB-IPDI)，襯層加厚到10 mm，丁羥推進劑加厚到200 mm。將粘接試件分別在25、50、70℃溫度下貯存50 d后，將襯層取下由內向外切成厚度約1 mm左右的切片，記錄襯層切片到界面的距離，用內標法測定切片中增塑劑的濃度。采用同樣方法測定推進劑中增塑劑含量為2%、3%、4%的粘接試件在70℃溫度下貯存50 d后襯層切片中增塑劑的濃度。

2.2 計算擴散系數

根據費克第二定律計算增塑劑擴散系數。在丁羥推進劑/襯層界面附近區域，假設增塑劑擴散符合一維情況，擴散系數D在整個擴散過程中不變，且聚合物體積不發生變化;由于推進劑相對襯層而言質量大得多，推進劑作為增塑劑擴散遷移的擴散源，可認為推進劑中增塑劑濃度始終不變;初始時襯層中增塑劑含量為零。可采用以下邊界條件來求解費克第二定律的表達式［14］。

邊界條件:

費克第二定律的表達式:

式中t為擴散時間，s;x為擴散距離，m;C為擴散劑濃度，mol/m3;Q為推進劑中增塑劑的表面濃度，mol/m2。

求解得

又因為D、Q均為常數，所以當擴散時間t一定時，式(3)可簡化為

作lnC-x2直線，由直線斜率可換算得到增塑劑表觀擴散系數D［15-16］。

2.3 擴散系數的實驗結果

根據式(4)對實驗數據進行處理，得到lnC-x2直線，由斜率可計算出增塑劑含量為3%，25、50、70℃溫度下增塑劑擴散系數分別為 9.08×10-12、9.68×10-12、10.18 ×10-12m2/s。隨溫度升高，擴散系數依次升高。計算70℃溫度下增塑劑含量為2%、3%、4%時，擴散系數分別為 9.68 ×10-12、10.18 ×10-12、9.37×10-12m2/s，擴散系數的變化趨勢不明顯，增塑劑含量為3%時，擴散系數相比而言最大。表1是增塑劑含量為3%、70℃，老化50 d后，襯層不同位置測得的增塑劑濃度和擴散系數。略去其他試驗數據。

3 結果與討論

表2是通過分子動力學模擬及加速老化實驗方法得到不同溫度下增塑劑DOS的擴散系數，DOS擴散系數模擬值和實驗值數量級為10-12m2/s。隨溫度升高，DOS擴散系數均逐漸增大。增塑劑與推進劑粘接體系相容性好［17］，由親合性拉力引起擴散運動，增塑劑不參與固化系統的反應，在高分子鏈段運動產生的間隙中(即自由體積孔穴之間)以“跳躍”方式發生空間位置遷移［18］。這在用分子動力學方法模擬DOS分子在HTPB-IPDI體系中的微觀運動過程中得到直觀演示。隨溫度升高，DOS擴散系數均逐漸增大的原因，經分析主要有2點:一是對于粘接體系而言，溫度升高，運動單元熱運動能量隨之升高，聚合物鏈段間的活動性增加，分子間距離增加，聚合物中的自由體積膨脹、溶脹度升高，使得增塑劑DOS進行擴散的有效空間增大;二是溫度升高，增塑劑分子的活動性增加，加上粘接體系的粘度下降，使得DOS的擴散行為更活躍。

表3是通過分子動力學模擬及加速老化實驗方法得到不同增塑劑含量下DOS的擴散系數，增塑劑含量為2%、3%、4%時，實驗得到的擴散系數變化不明顯，含量為3%時的擴散系數相比而言最大;當推進劑中增塑劑的含量為3%、6%、9%、15%(已超過推進劑中的真實含量2% ～4%)時，模擬得到的擴散系數依次略有下降。

表1 襯層不同位置的增塑劑濃度和擴散系數Table 1 Plasticizer concentration in the different position of liner and diffusion coefficient

表2 不同溫度下的增塑劑擴散系數(含量3%)Table 2 Plasticizer diffusion coefficient on the different temperatures(content 3%)

表3 不同增塑劑含量下的擴散系數Table 3 Plasticizer diffusion coefficient at different contents

共容增塑限量主要取決于粘合劑相的網絡結構，在固體火箭發動機推進劑中通常可加入2%～4%的增塑劑，推進劑中粘合劑相的含量一般不超過10%，那么增塑劑約占其粘接體系的20% ～35%(質量比)，這大概就是粘接體系的共容增塑限量［1］。而增塑劑在襯層中其含量極低，甚至為零，界面處濃度差引起的一種拉曳力會使增塑劑從推進劑遷移到襯層中。襯層中粘合體系約占40%，根據這個共容增塑限量值，可知襯層最多可吸收約8% ～14%的增塑劑，這與試驗現象(增塑劑遷移的相對含量反而減小，詳見文獻［19］)吻合。

在交聯聚合物的溶脹過程中，存在著2種相反趨勢的平衡過程:增塑劑力圖滲入高聚物內部，使高聚物體積膨脹，引起分子網交聯點之間分子鏈的伸展，降低高聚物的構象熵值;構象熵值降低又必然引起分子網的彈性收縮力，力圖使分子網收縮。當粘接體系中DOS含量超過共容增塑限量后，DOS的增加對溶脹度的影響很小，幾乎保持不變;過量的DOS占據大量空間，不僅束縛分子網的伸展，而且促使分子網收縮，使其自由體積減少，降低DOS的擴散活動。所以，擴散系數有所下降。

4 結論

(1)增塑劑擴散系數試驗值和模擬值的數量級基本一致，都為10-12m2/s。

(2)隨溫度升高，DOS擴散系數均逐漸增大。

(3)增塑劑含量為2%、3%、4%時，試驗得到的擴散系數變化不明顯;當推進劑中增塑劑的含量為3%、6%、9%、15%，模擬得到的擴散系數依次略有下降。

(4)雖然試驗過程中人為因素干擾較大、模擬計算的體系較小，且對實物作了必要簡化，但從試驗結果和模擬結果的對比分析可知，分子模擬可得到與實驗相一致的結果，增塑劑擴散系數的試驗結果可驗證其分子模擬結果，分子動力學方法計算增塑劑擴散系數更具優勢。

［1］徐思羽，唐大森.二茂鐵衍生物在固體推進劑應用中的遷移問題［J］.推進技術，1983，28(3):40-45.

［2］尹華麗，李東峰，王玉，等.組分遷移對NEPE推進劑界面粘接性能的影響［J］.固體火箭技術，2005，28(2):126-129.

［3］尹華麗，李東峰，張綱要，等.丁羥推進劑粘接體系中的組分遷移［J］.固體火箭技術，2007，30(5):416-419.

［4］姚維尚，李倩，譚惠民.NEPE推進劑粘合劑性能的分子模擬研究［J］.含能材料，2007，15(6):650-655.

［5］楊月誠，焦東明，強洪夫，等.HTPB推進劑組分溶度參數的分子模擬研究［J］.含能材料，2008，16(2):191-195.

［6］李紅霞，強洪夫，等.丁羥推進劑粘接體系中增塑劑遷移的分子模擬［J］.火炸藥學報，2008，31(5):74-78.

［7］李紅霞，強洪夫，王廣，武文明.基于MD方法的增塑劑擴散行為的模擬研究［J］.含能材料，2009，17(1):36-41.

［8］Kotelyanskii M J，Wagner N J，Paulaitis M E.Molecular dynamics simulation study of the mechanisms of water diffusion in dydrated，amorphous polyamide［J］.Computational and Theoretical Polymer Science，1999，9:301-306.

［9］Hofmann D，Fritz L，Ulbrich J，et al.Molecular simulation of small molecule diffusion and solution in dense amorphous polysiloxanes and polyimides［J］.Computational and Theoretical Polymer Science，2000，10:419-436.

［10］Guo Xiang-dan，Teng Jia-wei，Xie Zai-ku.Investigation of water adsorbed on HZSM-5:experiment and molecular simulation［J］.Computer and Applied Chemistry，2006，6:486-490.

［11］Liu Qing-lin，Huang Yu.Transport behavior of oxygen and nitrogen through organasilicon-containing polystyrenes by molecular simulation［J］.J.Phys.Chem.B，2006，100:17375-17382.

［12］于坤千，李澤生，李志儒，等.寡聚物在高分子母體中的擴散——分子動力學模擬研究［J］.高等學校化學學報，2002，23(7):1327-1330.

［13］Allen M P，Tildesley D J.Computer simulation of liquids［M］.Oxford:Clarendon Press，1987.

［14］王國梅，萬發榮.材料物理［M］.武漢:武漢理工大學出版社，2004:31.

［15］Gonzalez S，Gil M A，Hernandez J O，et al.Resistance to corrosion of galvanized steel covered with an epoxy-polyamide primer coating［J］.Progress in Organic Coatings，2001，41:3167-3170.

［16］Glavchev I，Petrova K，Devedjiev I.Investigation of solvent interaction in epoxy coatings［J］.Polymer Testing，2000，19:111-114.

［17］張景春.固體推進劑化學及工藝學［M］.長沙:國防科技大學出版社，1987:36.

［18］潘福生.高分子-無機雜化膜擴散特性的分子動力學模擬［D］.天津:天津大學，2006.

［19］李紅霞.丁羥推進劑/襯層粘接界面失效的分子模擬與實驗研究［D］.西安:第二炮兵工程學院，2009:92-93.

Measurement of diffusion coefficient of plasticizer in HTPB propellant

LI Hong-xia1，QIANG Hong-fu1，LI Xin-qi1，WANG Hong-wei2
(1.The Second Artillery Engineering University，Xi＇an 710025，China;2.The 46th Institute of the Sixth Academy of CASIC，Huhhot 010076，China)

Aiming at the plasticizer migration in the interface between HTPB propellant and liner，the diffusibility of plasticizer was studied.The diffusion coefficients were obtained by two methods.Firstly，plasticizer diffusivity in the bond system of HTPB propellant was simulated by the molecular dynamics method and the diffusion coefficient was obtained by Einstein equation.Secondly，taking the bonding interface of propellant and thickening liner as the object，the accelerated migration test of plasticizer was performed and the diffusion coefficient was obtained by the second law(Fick).The diffusion coefficients of dioctyl sebacate were discussed under different temperatures and at different contents of plasticizer.It was analyzed that plasticizer migrated by jump on the space situation which were produced when macromolecule chain moved in the HTPB bond system.The results indicated that the magnitude order of DOS diffusion coefficient was 10-12m2/s.The diffusion coefficients gradually increased with increasing of environmental temperature because of more intense activity and larger active space of plasticizer.Yet the diffusion coefficients decreased a little with the increasing of DOS content when content was more than 3%.It is more advantageous to measure diffusion coefficient of plasticizer by molecule dynamics.

HTPB propellant;plasticizer;diffusion coefficient;experiment;simulation

V512

A

1006-2793(2012)03-0387-04

2011-01-23;

2012-03-06。

國家863項目(2010AAJ119);國家973項目(973-61338)。

李紅霞(1979—)，女，博士，研究方向為納米材料應用與含能材料模擬。E-mail:lhxlihongxia@163.com

(編輯:劉紅利)