王 春，校金友，梁騰飛，文立華，張承雙，程 文，候 曉

(1.西北工業(yè)大學(xué) 航天學(xué)院，西安 710072；2.中國航天科技集團公司四院四十三所，西安 710025；3.中國航天科技集團公司四院，西安 710025)

0 引言

近年來，復(fù)合材料因其優(yōu)良的熱學(xué)與力學(xué)性能(高熱穩(wěn)定性、高比模量與高比強度等)在航空航天領(lǐng)域獲得了越來越廣泛的應(yīng)用。其中，聚對苯撐苯并二噁唑(Poly-p-phenylene-benzobisoxazole，PBO)具有優(yōu)異的力學(xué)性能與較強耐高溫性能，被視為新一代先進樹脂基復(fù)合材料的增強體[1]。但PBO纖維表面非常光滑且呈化學(xué)惰性，嚴(yán)重影響了其與環(huán)氧樹脂基體的界面結(jié)合性能，而界面作為基體與增強體的過渡部分，很大程度上影響了復(fù)合材料的性能[2-3]。

由于界面在復(fù)合材料性能中的關(guān)鍵作用，界面性能的研究受到科學(xué)界與工程界的廣泛關(guān)注。但界面區(qū)域尺度小、形態(tài)復(fù)雜，界面實驗較困難，實驗結(jié)果離散度大，事實上，界面的精確實驗面臨巨大挑戰(zhàn)。雖然許多實驗被用來測試界面的剪切模量與剪切強度，但對界面的法向力學(xué)性能的測試仍缺乏可靠實驗[4-5]。面對實驗無法解決的問題，許多理論及計算模擬方法被用來研究復(fù)合材料界面。例如，基于從頭算起的量子化學(xué)、密度泛函理論(DFT)和基于經(jīng)典力學(xué)的分子動力學(xué)(MD)等。其中，分子動力學(xué)基于原子層面建模，在納秒和納米的時空尺度復(fù)現(xiàn)真實的物理化學(xué)過程，從而克服了許多實驗的局限性，在材料領(lǐng)域得到廣泛應(yīng)用[6-7]，特別是對樹脂基體熱力學(xué)性質(zhì)的研究[8]。

目前，分子動力學(xué)已經(jīng)被廣泛用于碳納米管、石墨烯和納米二氧化硅等納米粒子聚合物復(fù)合材料的研究中，且取得了重要進展。Joel等[9]建立了炭纖維/環(huán)氧樹脂界面分子模型，運用分子動力學(xué)方法預(yù)測了界面的模量、強度等力學(xué)性能，并研究了炭纖維表面孔隙對力學(xué)性能的影響，但其模型中界面區(qū)域是人為定性確定的，并沒有定量地計算界面的位置和厚度，所以該方法預(yù)測的界面力學(xué)性能缺乏準(zhǔn)確性。Varshney等[10]將固化后的基體與纖維進行拼接，并弛豫獲得一個由纖維到基體的過渡部分，即界面模型，并運用該界面模型計算了炭纖維/環(huán)氧樹脂界面的熱阻系數(shù)，深入研究了不同纖維表面形態(tài)對熱阻系數(shù)的影響。但該界面建模方法是一個靜態(tài)過程，未考慮纖維與基體分子的相互作用對固化過程以及最終所形成的界面結(jié)構(gòu)的影響，這必將影響對界面性能的預(yù)測。因此，為了獲得真實的界面結(jié)構(gòu)，必須動態(tài)地模擬樹脂基體在纖維表面的固化過程。總的來說，運用分子動力學(xué)對纖維增強復(fù)合材料界面的研究還不是很完善，具體到PBO/環(huán)氧樹脂界面，相關(guān)的研究仍很匱乏。利用分子動力學(xué)研究纖維增強復(fù)合材料性能的報道較少，對PBO纖維增強復(fù)合材料的研究，就作者對目前的文獻查找，還未見相關(guān)公開的報道。

本文將采用分子動力學(xué)方法，運用動態(tài)建模方式，建立PBO/環(huán)氧樹脂界面，并研究其性能。工作主要包含兩部分：首先，基于通用的原子力場參數(shù)及X射線實驗結(jié)果建立PBO纖維的單胞模型，運用分子動力學(xué)研究其性能；其次，動態(tài)地模擬PBO纖維與TGDDM-DDS環(huán)氧樹脂基體界面的形成過程，并由體系原子濃度剖視圖確定界面區(qū)域。通過虛擬單軸拉伸模擬，研究界面的力學(xué)性能。本文分子動力學(xué)模擬均采用全原子力場CVFF(consistent valence force field)來描述原子間作用勢。

1 PBO纖維分子動力學(xué)建模

PBO分子結(jié)構(gòu)式如圖1所示。由于PBO纖維結(jié)構(gòu)較為復(fù)雜，運用分子動力學(xué)對PBO纖維的研究在文獻中報道很少。本文首先運用CVFF原子力場參數(shù)及X射線實驗結(jié)果建立PBO分子模型，并將分子動力學(xué)模擬結(jié)果與實驗結(jié)果進行對比，驗證分子模型的準(zhǔn)確性。

1.1 PBO分子模型

PBO分子鏈骨架是由苯環(huán)和苯并二噁唑環(huán)耦合而成。X射線實驗表明，苯環(huán)和苯并二噁唑環(huán)以180°鍵角連接形成共平面[11]，故PBO單鏈骨架上的所有原子近似處于同一平面。根據(jù)CVFF力場，可確定PBO分子中各原子的類型，如圖1所示，其中氫原子類型均為h。根據(jù)原子類型，從CVFF力場中可獲取PBO分子中的鍵長參數(shù)如表1所列，鍵角除C—O—C為109°，其余鍵角均為120°。根據(jù)以上參數(shù)，運用商業(yè)軟件Material Stu-dio[12]，即可建立近似的PBO纖維單鏈分子結(jié)構(gòu)。

圖1 PBO分子結(jié)構(gòu)式與其原子類型Fig.1 Chemical structures and atom types of PBO fiber

表1 PBO纖維分子鍵長參數(shù)Table1 Bond length of PBO fiber

PBO纖維由多條單鏈分子按照特定的晶格結(jié)構(gòu)堆疊而成，由X射線實驗獲取的纖維晶胞的晶格常數(shù)依次為a=11.2 ?，b=3.57 ?，c=11.02 ?(沿纖維方向)，γ=100.6°。每個纖維晶胞中包含兩條相同的PBO單鏈，兩條鏈在纖維方向的位錯距離為d1=2.63 ?，在a方向上間距為d2=5.58 ?[11]。由此可建立PBO纖維單胞模型如圖2所示。

(a)晶胞 (b)分子鏈

1.2 PBO分子模型驗證與性能預(yù)測

本文首先使用PBO單胞模型堆疊建立尺寸為47×49×55 ?3的模擬盒子，三個方向均采用周期性邊界條件。使用Material Studio中的Forcite模塊，對PBO纖維分子結(jié)構(gòu)進行幾何優(yōu)化，以得到CVFF力場下較為準(zhǔn)確的PBO纖維模型。然后，運用該模塊在300 K與1 atm下，算得PBO纖維的密度與模量，如表2所示。表2中，PBO批量生產(chǎn)數(shù)據(jù)由相關(guān)研究院所提供。

表2 PBO性能對比Table2 Properties comparison of the PBO

運用構(gòu)建的分子模型預(yù)測PBO纖維的宏觀性能，并與實驗對比，可驗證該模型的準(zhǔn)確性。如表2所示，在PBO密度方面，分子動力學(xué)模擬(MD)的計算結(jié)果與批量生產(chǎn)和實驗室獲得的PBO纖維的測量結(jié)果均吻合較好(相對誤差為5%)。在PBO的模量方面，分子動力學(xué)模擬的計算結(jié)果與實驗室獲得的PBO纖維的測量結(jié)果接近(相對誤差為10.6%)，而遠(yuǎn)大于批量生產(chǎn)的PBO纖維的模量(相對誤差為88%)。其主要原因是MD模擬中建立的是理想的無任何缺陷的PBO模型，而由于工藝技術(shù)的限制，批量生產(chǎn)的PBO難免存在缺陷，如納米級孔隙、雜質(zhì)以及纖維取向不完全一致等，這些缺陷致使纖維的力學(xué)性能有所下降。通過以上MD計算結(jié)果與實驗結(jié)果的對比，驗證了本文建立的PBO纖維模型的準(zhǔn)確性。

2 界面的分子動力學(xué)建模

本文界面的動態(tài)建模過程主要分兩步：(1)將固化前的環(huán)氧樹脂與固化劑分子在纖維表面充分混合；(2)通過分子動力學(xué)方法動態(tài)地模擬其固化反應(yīng)過程。此過程充分考慮了纖維與基體分子的相互作用對固化過程的影響，動態(tài)模擬獲得的界面結(jié)構(gòu)更真實合理。

2.1 基體模型

樹脂基體由4,4-二氨基二本甲烷環(huán)氧樹脂(TGDDM)與3,3-二氨基二苯砜固化劑(33DDS)固化而成，其化學(xué)結(jié)構(gòu)式如圖3所示。固化劑與樹脂基體的交聯(lián)反應(yīng)(圖4)主要分兩個過程：(1)環(huán)氧基開環(huán)與固化劑上的伯胺基團加成反應(yīng)生成仲胺；(2)環(huán)氧基團開環(huán)與仲胺基團反應(yīng)生成叔胺，氧原子與氫原子形成羥基。

圖3 環(huán)氧樹脂與固化劑的化學(xué)結(jié)構(gòu)Fig.3 Chemical structures of epoxy monomers and curing agent

圖4 固化反應(yīng)過程Fig.4 Curing reaction path

在MD模擬中，整個固化過程采用以截斷距離(環(huán)氧基團與氨基之間的距離)為基準(zhǔn)，分步成鍵的方法實現(xiàn)[13]。首先給定一個初始判定距離r，當(dāng)反應(yīng)基團之間距離小于此判定距離時，則固化反應(yīng)發(fā)生；待系統(tǒng)中所有此類反應(yīng)完成后，則逐步增大判定距離r，繼續(xù)發(fā)生反應(yīng)；當(dāng)判定距離r達(dá)到給定的判定距離的最大值或達(dá)到目標(biāo)交聯(lián)度時，反應(yīng)結(jié)束。根據(jù)多次模擬測試，在充分交聯(lián)且盡量減少模擬時間的條件下，確定本次模擬的判定距離為3～10 ?，判定距離增量步長為1 ?。對于每個特定的判定距離，交聯(lián)鍵生成后，首先使用共軛梯度法對系統(tǒng)進行能量最小化，然后在NPT系綜下進行10 ps的模擬，使系統(tǒng)弛豫到穩(wěn)定狀態(tài)，以最大限度地保證固化過程的合理性。

整個交聯(lián)過程使用perl腳本在Material Studio實現(xiàn)，表3給出樹脂基體的各組分信息及摩爾比。

表3 基體組分Table3 Resin components

2.2 界面模型

界面是從纖維到基體的過渡部分，界面模型的建立過程如圖5所示。將固化劑分子與環(huán)氧樹脂分子在PBO纖維右側(cè)按摩爾比隨機混合，模型的初始尺寸為47×49×123 ?3，三個方向均為周期性邊界條件。首先，運用共軛梯度法對初始模型進行能量最小化；然后，在NPT系綜下以1 fs為時間步長，弛豫1 ns。溫度和壓強分別為實驗中的固化溫度447 K與固化壓強1 atm。通過弛豫過程，纖維將與固化前的環(huán)氧體系充分作用，系統(tǒng)勢能逐漸趨于一個穩(wěn)定值，表明未固化的基體混合物在纖維表面附近達(dá)到平衡狀態(tài)。

圖5 界面建模過程Fig.5 Process of interphase model establishment

系統(tǒng)達(dá)到平衡狀態(tài)后，開始進行交聯(lián)反應(yīng)。交聯(lián)反應(yīng)過程中參數(shù)控制與混合弛豫過程相同。在反應(yīng)過程中，當(dāng)交聯(lián)度依次達(dá)到50%、60%、70%和80%時，將所得模型分別保存下來，以供下一步分析。本文中，交聯(lián)度的定義為實際參與反應(yīng)的環(huán)氧基團數(shù)目與體系中總的環(huán)氧基團數(shù)目的比值。當(dāng)交聯(lián)過程完成后，再次在NPT系綜下進行弛豫(300 K和1 atm)，弛豫時間為2 ns，最終得到常溫常壓下的界面模型。界面模型建立完成后，由體系相對原子濃度剖視圖，可得原子密度分布和界面區(qū)域厚度。

2.3 界面厚度

確定界面區(qū)域與厚度是進行虛擬單軸拉伸的前提。由體系沿界面法向的相對原子濃度剖視圖，可確定界面位置及厚度。相對原子濃度是關(guān)于距離的函數(shù)，在界面的法向方向上將體系切成多個子域，然后計算每個子域的相對原子濃度，其定義為

式中Nslab和Vslab分別為各個子域包含的原子數(shù)目和子域體積；Ntot和Vtot分別為整個基體包含的原子數(shù)目和體積。

本文根據(jù)交聯(lián)完成后，體系弛豫過程的最后500 ps的原子軌跡統(tǒng)計得到相對原子濃度剖視圖，結(jié)果如圖6所示，橫坐標(biāo)為距離纖維表面的距離，縱坐標(biāo)為相對原子濃度。由圖6可見，原子濃度在纖維附近出現(xiàn)一個波峰，而在遠(yuǎn)離纖維后趨于平穩(wěn)。波峰的出現(xiàn)是基體在纖維表面的范德瓦爾長程吸引力的作用下自適應(yīng)形成的中間相，其原子濃度比平穩(wěn)區(qū)高20%～30%，這段波動區(qū)域即界面區(qū)域。本文以相對原子濃度曲線回歸到1為標(biāo)準(zhǔn)測量界面區(qū)域厚度。結(jié)果表明，不同交聯(lián)度下界面的厚度基本不變，都在5～6 ?，與炭纖維、硅顆粒等增強體增強復(fù)合材料的界面厚度相比較小[14-17]。

圖6 界面在不同交聯(lián)密度下的原子濃度剖視圖Fig.6 Profile of the relative atomic concentration of the interphase with different crosslink density

3 界面力學(xué)性能

3.1 虛擬單軸拉伸模擬

界面建模完成后，將Material Studio中建立的模型導(dǎo)入到分子動力學(xué)模擬工具Lammps[18]中進行虛擬單軸拉伸。拉伸方式如圖7所示，除界面區(qū)域允許產(chǎn)生應(yīng)變外，其余部分在拉伸過程中均固定不動。此設(shè)置可排除基體和纖維的影響，單獨測量界面的模量。虛擬拉伸過程通過階梯式改變系統(tǒng)在拉伸方向上的尺寸實現(xiàn)，即由一系列NVT過程組成。在每個NVT過程中，保持系統(tǒng)的原子數(shù)、體積和溫度不變，而在不同的NVT過程之間，系統(tǒng)會在拉伸方向上發(fā)生變形，所有原子在拉伸方向上的坐標(biāo)同時發(fā)生等比例調(diào)整。在兩個側(cè)向方向采用周期性邊界條件，且長度固定。虛擬拉伸中界面溫度控制為300 K，時間步長1 fs。

圖7 界面拉伸模型Fig.7 Schematic illustration of the process of uniaxial tensile test

由于分子動力學(xué)模擬在時間尺度上的限制，本文的模擬拉伸速度為5×108s-1。此虛擬拉伸速度雖然遠(yuǎn)大于實際的拉伸速度，但此前研究表明，分子動力學(xué)方法預(yù)測的材料楊氏模量值受拉伸速度的影響很小[14]。拉伸應(yīng)力由維里應(yīng)力公式求得：

式中Ω為研究對象的體積，在本次模擬中即界面區(qū)域體積；mk與uk分別為k號原子質(zhì)量與速度；rkl為k號與l號原子的距離；fkl為l號原子對k號原子的作用力。

本文重點研究界面的橫向模量，即沿界面法向方向(z方向)的楊氏模量，故虛擬單軸拉伸的方向沿界面的法向方向。界面處的本構(gòu)方程可寫為

式中σx、σy和σz分別為界面軸向、周向與法向應(yīng)力；εx、εy和εz分別為相應(yīng)方向上的應(yīng)變。

由于拉伸在NVT系綜下進行，系統(tǒng)尺寸在兩個側(cè)向方向保持不變，除法向方向應(yīng)變εz不為0外，其余方向應(yīng)變均為0。因此，由法向虛擬單軸拉伸得到的應(yīng)力-應(yīng)變曲線可解出與界面法向相關(guān)的剛度系數(shù)Kzz、Kxz和Kyz。

3.2 結(jié)果與討論

為減小系統(tǒng)初始不確定性對結(jié)果的影響[19]和應(yīng)力波動，本文在NPT系綜下，分別把系統(tǒng)弛豫1、1.2、1.4 ns得到三個子體系，然后分別對其進行拉伸再取系綜平均。由于界面實驗十分復(fù)雜，當(dāng)前并沒有可靠的實驗方案，本文采用純基體作為對照組來驗證建模方法的準(zhǔn)確性，并研究界面的性能相對于純基體的變化。圖8為在300 K下不同交聯(lián)度界面的拉伸曲線，圖9為同等溫度下60%交聯(lián)度純基體的拉伸曲線。界面與純基體的剛度系數(shù)均通過應(yīng)力-應(yīng)變曲線中εz<2%的線性段擬合得出。圖8、圖9中出現(xiàn)的應(yīng)力-應(yīng)變曲線的初始應(yīng)力不為0現(xiàn)象是由于MD系統(tǒng)的應(yīng)力漲落產(chǎn)生的，其漲落幅度由模擬系統(tǒng)的尺寸決定。本文模擬中的漲落幅度遠(yuǎn)小于拉伸過程中產(chǎn)生的真實應(yīng)力，故其對所測量的彈性模量影響很小。

(a)50%交聯(lián)度

(b)60%交聯(lián)度

(c)70%交聯(lián)度

(d)80%交聯(lián)度

根據(jù)圖9的拉伸曲線，可得到60%交聯(lián)度純基體的剛度系數(shù)為Kzz=8.356 GPa、Kxz=4.338 GPa、Kyz=4.087 GPa。由剛度系數(shù)求得60%交聯(lián)度純基體的楊氏模量為5.535 GPa，略大于Okabe[20]等的實驗值4.5 GPa(23%)。這主要是因為分子動力學(xué)建立的模型為理想模型，不考慮任何缺陷和孔隙，故由分子動力學(xué)獲得的預(yù)測值一般會偏大[21]，由以上純基體結(jié)果可驗證本文建模方法的準(zhǔn)確性。

圖9 60%交聯(lián)度純基體拉伸曲線Fig.9 Stress-stain curves of pure matrix with 60% convers

根據(jù)圖8的拉伸曲線，得到不同交聯(lián)度下界面的力學(xué)參數(shù)，如表4所示，其中E為界面的橫向模量。由表4可見，隨著交聯(lián)度的提升，界面力學(xué)性能有所增強，但幅度有限。80%交聯(lián)度與50%交聯(lián)度相比，彈性模量僅增加1.27 GPa。對比基體與界面的結(jié)果發(fā)現(xiàn)，同在60%交聯(lián)度下，界面的彈性模量較基體高出4.465 GPa。以上數(shù)據(jù)表明，交聯(lián)度對界面的模量影響較小，界面模量的增強主要是PBO纖維與基體之間的范德瓦爾力產(chǎn)生的吸附作用造成的。由于該吸附作用，纖維和基體緊緊吸附在一起，從而在纖維表面形成了一層致密的界面層，在單軸拉伸時，會產(chǎn)生更強的抵制外力變形的效果，從而界面比基體具有更高的彈性模量。類似現(xiàn)象在文獻[22-23]中也有報道。

表4 不同交聯(lián)度下界面力學(xué)參數(shù)Table4 Mechanics parameters of interphase in different cross-linking conversion GPa

4 結(jié)論

(1)本文通過建立PBO纖維分子模型，預(yù)測了其極限模量為527.7 GPa，略大于實驗室條件下獲取的最大模量，而遠(yuǎn)大于批量生產(chǎn)條件下的模量，表明PBO生產(chǎn)工藝還有較大的提升空間。

(2)本文建立了四種不同交聯(lián)密度的界面模型，由原子濃度剖視圖發(fā)現(xiàn)，不同交聯(lián)密度下界面厚度區(qū)別不大，均為5～6 ?，小于炭纖維、硅顆粒等增強體復(fù)合材料界面厚度。

(3)由界面拉伸曲線發(fā)現(xiàn)，隨著交聯(lián)度的提高，界面的力學(xué)性能得到增強，且交聯(lián)度對界面力學(xué)性能影響較小。界面力學(xué)性能的增強主要是由于纖維與界面之間的分子作用力造成的。