， ， ， ， ，

(西北工業(yè)大學(xué) 機(jī)電學(xué)院，西安 710072)

1 引 言

隨著三維打印技術(shù)的發(fā)展，基于微滴噴射原理的三維粉末粘結(jié)打印技術(shù)已在很多領(lǐng)域得到應(yīng)用，特別是在骨支架制備方面技術(shù)已比較成熟[1,2]。該技術(shù)主要是通過打印噴頭噴射膠水，使膠水鋪展在粉末材料的表面，從而實現(xiàn)粉末材料的粘結(jié)成型。其結(jié)合了快速成型技術(shù)的特點，在計算機(jī)的輔助下可以實現(xiàn)任意結(jié)構(gòu)復(fù)雜零部件的制造，因此該技術(shù)在人工骨支架的生物制造領(lǐng)域有著極其重要的作用。但是，基體材料的選擇一直是阻礙該技術(shù)廣泛應(yīng)用的主要因素，因為骨支架植入體內(nèi)不僅會與人體組織直接接觸，還需要暫時起到支撐人體組織和器官的作用，所以材料必須具有一定的生物相容性及力學(xué)性能。目前，羥基磷灰石HA(Hydroxyapatite)因其具有優(yōu)良的生物相容性、骨傳導(dǎo)性以及與天然人骨相類似的化學(xué)成分，是制備人工骨最為理想的材料。因此，制備工藝中粘結(jié)劑的選擇成了人工骨支架制備的關(guān)鍵。當(dāng)前使用最多的是酸性粘結(jié)劑，因為酸性粘結(jié)劑能與HA界面發(fā)生一系列的化學(xué)反應(yīng)，使制備的骨支架具有較好的力學(xué)性能[3,4]。但是酸性粘結(jié)劑會給骨支架帶來污染，帶酸性的骨支架植入體內(nèi)后很容易引起周圍組織的炎癥和排異。解決這一問題最有效的方法是先采用聚合物粘結(jié)劑制備骨支架雛形，然后再使用酸性粘結(jié)劑或者其他粘結(jié)劑進(jìn)行固化處理。然而，該方法還存在很多的不確定因素，特別是內(nèi)部的微觀結(jié)構(gòu)和力學(xué)性能。在粉末粘結(jié)三維打印骨支架的制備過程中，液滴與粉末的相互作用決定了骨支架制件的微觀結(jié)構(gòu)和宏觀機(jī)械性能[5]，由于HA陶瓷材料的力學(xué)性能遠(yuǎn)高于聚合物粘結(jié)劑，所以聚合物粘結(jié)劑的力學(xué)性能對骨支架的機(jī)械強度起著決定性的作用，研究粘結(jié)劑的基本性能對骨支架的制備工藝具有十分重要的意義。而最為常用的三種聚合物粘結(jié)劑為聚乙烯毗咯烷酮PVP(polyvinylpyrrolidone)、聚丙烯酰胺PAM Polyacrylamide和聚乙烯醇PVA(polyvinyl alcohol)，主要是因為這三種粘結(jié)劑具有良好的生物惰性和生物相容性，并且不會給HA帶來污染。

粘結(jié)劑的微觀結(jié)構(gòu)及其與HA的相互作用是一種微觀現(xiàn)象，通過傳統(tǒng)實驗方法很難實現(xiàn)觀察和研究。而分子動力學(xué)MD(molecular dynamics)模擬方法是研究高分子與無機(jī)鹽材料界面相互作用的一種最有效的方法[6]。MD是一種基于原子力場的計算方法，能夠模擬不同物質(zhì)之間的力學(xué)及電學(xué)等性質(zhì)，計算空間距離從幾個埃到幾個納米，可以準(zhǔn)確解釋物質(zhì)之間的相互作用機(jī)理。目前該技術(shù)已經(jīng)在很多領(lǐng)域得到了成熟的運用[7,8]，其中圍繞HA展開的研究也有很多，Zhang等[9]利用MD模擬研究了在沒有耦合劑作用的條件下，人工骨組成材料羥基磷灰石與復(fù)合材料界面的相互作用，比較分析了聚合物在三個不同平面的約束能，對人工骨制造具有一定的指導(dǎo)作用。潘海華等[10]研究了人工骨材料羥基磷灰石HA晶面上的水分子行為，發(fā)現(xiàn)HA晶面附近的結(jié)構(gòu)化水層可阻止溶液離子自由進(jìn)入晶面，對HA顆粒在水溶液中的動力學(xué)穩(wěn)定性具有重要影響。李敏等[11]利用MD技術(shù)研究了HA及復(fù)合材料的制備與性能，特別是人工骨的制備工藝及彈性模量的研究對生產(chǎn)制造意義重大。

本文主要借助分子動力學(xué)方法對骨支架制備工藝中常用三種聚合物粘結(jié)劑聚乙烯毗咯烷酮PVP、聚丙烯酰胺PAM和聚乙烯醇PVA的結(jié)構(gòu)和性能進(jìn)行研究，從分子層面對粘結(jié)劑的相關(guān)性質(zhì)進(jìn)行了預(yù)估，為骨支架3D打印工藝制備工藝中膠水的選擇提供理論指導(dǎo)。

2 模型的建立與模擬方法

對于高分子聚合物體系的模擬仿真，合適聚合度的選取尤為重要，關(guān)系到模擬結(jié)果的準(zhǔn)確性。為了確定聚合物分子鏈的聚合度，本文根據(jù)三種粘結(jié)劑的分子結(jié)構(gòu)式構(gòu)建不同聚合度的分子鏈模型,如圖1所示。對構(gòu)建的不同重復(fù)單元數(shù)的三種粘結(jié)劑分子鏈Smart minimizer 方法進(jìn)行優(yōu)化，使其能量最小化，當(dāng)?shù)芰窟_(dá)到0.001 kcal/mol時停止計算。然后對所建模型進(jìn)行熱處理,(1) NVT系綜下分子動力學(xué)模擬，模擬溫度以50 K為步長從298 K逐漸升到598 K，然后從598 K又逐漸降到298 K，每步計算50 ps；(2) NVT系綜下分子動力學(xué)模擬，在298 K的溫度下對步驟(1)的最后結(jié)構(gòu)進(jìn)行200 ps的NVT分子動力學(xué)模擬，從而得到室溫下高分子鏈的平衡結(jié)構(gòu)。最后，對平衡結(jié)構(gòu)進(jìn)行數(shù)值分析，得到高分子鏈的內(nèi)聚能密度，并根據(jù)粘結(jié)劑的內(nèi)聚能密度變化來確定合適的聚合度選擇范圍。此外，聚合度的選擇還必須滿足每條粘結(jié)劑高分子鏈的分子量接近，這樣才能使其具有可比性，詳細(xì)的計算過程可以參考文獻(xiàn)[12]。

圖1 PVP，PAM和PVA粘結(jié)劑的分子結(jié)構(gòu)式

Fig.1 Molecular structural formula of PVP,PAM and PVA

本文PVP，PAM和PVA三種粘結(jié)劑分子鏈選取的聚合度分別為20，31和50，為了對三種粘結(jié)劑的基本性能進(jìn)行計算，分別建立了三種粘結(jié)劑的計算模型，模型的能量最小化過程同上面分子鏈的優(yōu)化過程一致。然后，對能量最小化的三種分子模型進(jìn)行分子動力學(xué)模擬，模擬過程如下。(1) NPT系綜模擬，模擬時間為200 ps，或者待密度達(dá)到平衡后停止計算。(2) NVT系綜模擬，模擬時間為50 ps，并對該過程的軌跡文件進(jìn)行數(shù)值分析。最終平衡后的分子模型如圖2所示。

所有分子動力學(xué)過程均是在COMPASS[13]力場下進(jìn)行。且在模擬中，采用Andersen恒溫器[14]，用Velocity Verlet積分法求解牛頓運動方程，原子初始速度按Boltzmann隨機(jī)分布方法確定。計算采用的時間步長為1 fs，模擬過程中每500步記錄一次體系的軌跡。模擬過程中Vander Waals和Coulomb分別按Atom based[15]和Ewald[16]方法計算，截斷半徑為9.5 ?，樣條寬度取1 ?，緩沖寬度取0.5 ?，溫度設(shè)置為298 K。

3 結(jié)果和討論

計算模型達(dá)到平衡是研究材料性質(zhì)的前提，關(guān)系到計算結(jié)果的可靠性。判斷體系是否達(dá)到平衡所依據(jù)的參數(shù)一般有能量、溫度、密度以及徑向分布函數(shù)等[17]。本文采用最常用的密度和能量來判斷體系的平衡，在體系密度和能量等圍繞一恒定值發(fā)生微小波動時，認(rèn)為體系達(dá)到了平衡。PAM優(yōu)化后體系的密度和能量的變化曲線如圖3所示?？梢钥闯?，經(jīng)過NPT分子動力學(xué)優(yōu)化后，模擬體系已達(dá)到了平衡。

圖2 三種粘結(jié)劑的平衡分子模型

Fig.2 Equilibrium models of three binders

內(nèi)聚能密度CED(Cohesive Energy Density)是表征高分子材料的一個重要參數(shù)，是評價分子間作用力大小的一個物理量，主要反映分子基團(tuán)間的相互作用。凡與物質(zhì)間相互作用有關(guān)的性質(zhì)和物性屬性與內(nèi)聚能密度都具有內(nèi)在聯(lián)系，如物質(zhì)的溶解性、粘度和相容性等，其表達(dá)式為

CED=Ecoh/v

(1)

式中Ec oh指內(nèi)聚能，v為物質(zhì)的量混合體積。一般來說，分子中所含基團(tuán)的極性越大，分子間的作用力就越大，則相應(yīng)的內(nèi)聚能密度就越大。

從圖4三種粘結(jié)劑內(nèi)聚能的計算結(jié)果可以看出,三種粘結(jié)劑內(nèi)聚能的大小關(guān)系為PAM>PVA>PVP，PAM的內(nèi)聚能密度最大，PVP的內(nèi)聚能密度最小。而CED的大小可以反映體系內(nèi)部分子間相互作用的強弱，即內(nèi)聚能密度越大，體系內(nèi)部分子間相互作用力越大。因此，三種高分子粘結(jié)劑

圖3 PAM無定形晶胞體系密度和能量隨模擬時間的變化曲線

Fig.3 Plot of the density and total potential energy,nonbond energy of PAM binder system as a function of MD simulation time

體系內(nèi)部相互作用力的大小關(guān)系為PAM>PVA>PVP，這與平衡結(jié)構(gòu)密度的計算結(jié)果一致。

力學(xué)性能是粘結(jié)劑材料的重要性能之一，因為在骨支架的粘結(jié)打印工藝中，粘結(jié)劑的力學(xué)性能直接關(guān)系到所制備骨支架的力學(xué)性能。本文通過對體系平衡結(jié)構(gòu)進(jìn)行靜態(tài)力學(xué)性能分析來計算粘結(jié)劑的力學(xué)性能，如彈性系數(shù)、模量和泊松比等。在原子水平計算中，靜態(tài)模型計算內(nèi)應(yīng)力張量σ[18]采用的是維里方程：

(2)

式中mi和Vi分別為原子的質(zhì)量和速度，V0為無形變時的體積。

彈性系數(shù)矩陣由數(shù)值法求解相應(yīng)應(yīng)力與應(yīng)變的一階導(dǎo)數(shù)而得到，采用最小二乘法擬合求得拉伸模量E和泊松比γ，進(jìn)而可求得其他有效的各向同性力學(xué)性能，如剪切模量G、體積模量K及柯西壓(C12-C44)等。對于各項同性材料，其模量可用拉梅參數(shù)表示為

(3)

通常材料的硬度與其拉伸模量E和剪切模量G密切相關(guān)，其值越大，材料的剛性和硬度越高。材料斷裂強度的大小由體模量K度量，K值越大，表明材料的抗斷裂能力越大?？挛鲏和ǔＳ脕眍A(yù)估體系的延展性，若為負(fù)值則材料顯脆性，若是正值則表明材料延展性較好，其值越高表明體系的延展性越好[19]。

由表2可知，三種粘結(jié)劑的力學(xué)性能存在一定的差異，且三種粘結(jié)劑模量和柯西壓值的大小關(guān)系表現(xiàn)一致，都為PVA>PVP>PAM，也就是說相對于PVP和PAM，PVA的力學(xué)性能較好。該計算結(jié)果與粘結(jié)劑的內(nèi)聚能結(jié)果不一致，說明粘結(jié)劑的力學(xué)性能與內(nèi)聚能密度和質(zhì)量密度不存在特定的關(guān)系。此外，三種粘結(jié)劑的泊松比都在塑料的泊松比范圍之內(nèi)，表明三種粘結(jié)劑具有一定的塑性。同時，柯西壓值都為正值，說明三種粘結(jié)劑都具有一定的延展性，且PVA的延展性最好，PAM最差。

圖4 三種粘結(jié)劑的內(nèi)聚能柱狀圖

Fig.4 CED Histograms of three polymer binders

為了進(jìn)一步探索粘結(jié)劑在生物陶瓷表面的粘結(jié)行為，建立了三種粘結(jié)劑在羥基磷灰石HA(hydroxyapatite)界面的相互作用模型，所選HA界面是表面原子密度最大的HA(110)界面[5]。所建計算模型的平衡構(gòu)型如圖5所示，該模型的分子動力學(xué)計算模擬過程與粘結(jié)劑體系的模擬細(xì)節(jié)類似。

表2 粘結(jié)劑在298 K溫度下的力學(xué)性能

(單位：GPa)

Tab.2 Mechanical properties of binders at 298 K(unit:GPa)

SystemPVPPAMPVAC119.535.7715.99C228.546.0814.23C337.354.5614.95C441.802.392.23C552.912.763.64C662.392.124.02C124.504.277.89C133.873.557.16C233.163.5810.25C15-1.140.84-0.27C25-0.650.32-0.62C350.17-0.480.12C46-0.20-0.190.13Tensile modulus6.814.289.22Bulk modulus5.273.588.09Shear modulus2.651.653.41Poisson ratio0.280.290.31C12-C442.701.885.66

高分子粘結(jié)劑與生物陶瓷間的粘結(jié)力大小可以通過材料間界面結(jié)合能的大小來反映，界面結(jié)合能越大，表示粘結(jié)劑與生物陶瓷材料間的粘結(jié)力就越大。因此，本文通過計算粘結(jié)劑與HA(110)界面間的結(jié)合能來評價三種粘結(jié)劑在生物陶瓷表面的粘結(jié)特性。其表達(dá)式如下，

Ebind=-Einter=-(Etotal-Epolymer-EHA)

(4)

式中Etotal，Epolymer和EHA(110)分別為高聚物與HA相互作用體系、高聚物和HA(110)的能量。通過對PVP/HA(110)和PAM/HA(110)和PVA/HA(110)的相互作用平衡體系進(jìn)行MD模擬計算，得到的結(jié)合能列入表3。

由表3可知，三種粘結(jié)劑與HA(110)的結(jié)合能大小關(guān)系是PAM>PVA>PVP。這表明相比PVA和PVP兩種粘結(jié)劑，PAM與HA(110)表面的粘結(jié)作用力更強，這與三種粘結(jié)劑的內(nèi)聚能大小關(guān)系是一致的。產(chǎn)生這種現(xiàn)象的主要原因是，在相同質(zhì)量的三種粘結(jié)劑中含有不同數(shù)量和類型的極性官能團(tuán)。其中，一條PVP分子鏈含有20個羰基(-C=O)，PVA分子鏈中含有50個羥基(-OH)，但是在PAM分子鏈中卻含有31個羰基(-C=O)和31個氨基(-NH2)。因此，在三種粘結(jié)劑中，PAM的極性是最強的，這也是PAM/HA(110)界面結(jié)合能最大的根本原因所在，同時也說明高分子粘結(jié)劑的化學(xué)結(jié)構(gòu)和官能團(tuán)是影響粘結(jié)劑粘結(jié)性能的關(guān)鍵所在。

表3 三種相互作用模型的結(jié)合能(單位：kJ/mol)

Tab.3 Binding energies of three interaction models(unit:kJ/mol)

SystemEtotalEploymerEHAEbindPVP/HA(110)39806733.5736773.0339917760.62147800.08PAM/HA(110)39706546.8958660.4639917652.18269765.65PVA/HA(110)39744011.1349372.8839917777.98223139.73

圖5 三種粘結(jié)劑與HA(110)相互作用模型

Fig.5 Interaction models of three binders on HA(110) interface

4 結(jié) 論

本文借助分子動力學(xué)方法，對骨支架粉末粘結(jié)工藝中三種常用聚合物粘結(jié)劑PVP，PAM和PVA進(jìn)行建模，并且對其平衡結(jié)構(gòu)進(jìn)行了計算分析，得到了三種粘結(jié)劑的基本性質(zhì)。結(jié)果表明，粘結(jié)劑的內(nèi)聚能大小關(guān)系為PAM>PVA>PVP，這與粘結(jié)劑體系的密度大小關(guān)系一致。靜態(tài)力學(xué)分析結(jié)果表明，三種粘結(jié)劑體系的機(jī)械強度大小關(guān)系為PVA>PVP>PAM，其與體系內(nèi)聚能和密度的大小關(guān)系并不一致，說明內(nèi)聚能與力學(xué)性能之間并不存在特定的關(guān)系。此外，界面結(jié)合能的計算結(jié)果表明，三種粘結(jié)劑與HA 的結(jié)合能大小關(guān)系是PAM>PVA>PVP，與粘結(jié)劑內(nèi)聚能的大關(guān)系一致，產(chǎn)生這一現(xiàn)象的根本原因是相同質(zhì)量高分子粘結(jié)劑中所含官能團(tuán)的數(shù)量和類型不一樣。因此，高分子粘結(jié)劑的化學(xué)結(jié)構(gòu)和官能團(tuán)是影響粘結(jié)劑粘結(jié)性能的關(guān)鍵所在。

:

[1] Chiu C,Ferreira J,Luo T J,et al.Direct scaffolding of biomimetic hydroxyapatite -gelatin nanocomposites using aminosilane cross-linker for bone regeneration [J].JournalofMaterialsScience-MaterialsinMe-dicine,2012,23(9):2115-2126.

[2] Wang Y E,Li X P,Wei Q H,et al.Study on the mechanical properties of three -dimensional directly binding hydroxyapatite powder [J].CellBiochemBiophys,2015,72(1):289-295.

[3] Bergmann C,Lindner M,Zhang W,et al.3D printing of bone substitute implants using calcium phosphate and bioactive glasses [J].JournaloftheEuropeanCeramicSociety,2010,30(12):2563-2567.

[4] Vorndran E,Klarner M,Klammert U,et al.3D powder printing of b-tricalcium phosphate ceramics using different strategies [J].AdvancedEngineeringMaterials,2008,10:B67-B71.

[5] Wei Q H,Wang Y N,Li X P,et al.Study the bonding mechanism of binders on hydroxyapatite surface and mechanical properties for 3DP fabrication bone scaffolds [J].JournaloftheMechanicalBehaviorofBiomedicalMaterters,2016,57:190-200.

[6] Wang Y N,Wei Q H,Pan F L,et al.Molecular dynamics simulations of the mechanical properties of HA/NBCA under Additive Manufacturing[J].BiomedicalMaterialsandEngineering,2014,24:825-833.

[7] 干 湧,江 山,蘇昱臻,等.基于分子動力學(xué)與耗散粒子動力學(xué)串行耦合的面心立方金屬粗粒化模型[J].計算力學(xué)學(xué)報,2016,33(4):621-628.(GAN Yong,JIANG Shan，SU Yu-zhen.A coarse -grained model for fcc metals based on hierarchical coupling between molecular dynamics and isothermal dissipative particle dynamics [J].ChineseJournalofComputationalMechanics,2016,33(4):621-628.(in Chinese))

[8] 田 霞,崔俊芝.基于分子動力學(xué)模擬的多晶結(jié)構(gòu)微觀熱-力耦合行為的計算[J].計算力學(xué)學(xué)報,2012,29(1):95-98.(TIAN Xia,CUI Jun-zhi.Computations of the thermo -mechanical behaviors of polycrystalline structures at micro -scale based on molecular dynamics [J].ChineseJournalofComputationalMecha-nics,2012,29(1):95-98.(in Chinese))

[9] Zhang H P，Lu X，Leng Y，et al.Molecular dynamics simulations on the interaction between polymers and hydroxyapatite with and without coupling agents[J].ActaBiomaterialia,2009,5:1169-1181.

[10] 潘海華,陶錦輝,吳 韜,等.羥基磷灰石界面水行為的分子模擬 [J].無機(jī)化學(xué)學(xué)報，2006,22(8):1392-1400.(PAN Hai-hua,TAO Jing-hui,WU Tao,et al.Molecular simulation of water behavior on Hydroxya-patite crystal faces [J].ChineseJournalofInorganicChemistry,2006,22(8):1392-1400.(in Chinese))

[11] 李 敏.骨組織工程用多孔羥基磷灰石及其復(fù)合材料的制備與性能[M].湖南：湖南科技大學(xué)，2008.(LI Min.PreparationandPropertiesofHydroxyapatiteandItsCompositesinBoneTissueEngineering[M].Hunan:University of Science and Technology of Hunan,2008.(in Chinese))

[12] Wei Q H,Wang Y E,Chai W H,et al.Effects of composition ratio on the properties of poly (vinyl alco -hol)/poly (acrylic acid) blend membrane:a molecular dynamics simulation study [J].MaterialsandDesign,2016,89:848-855.

[13] Sun H.COMPASS:an ab initio force -field optimized for condensed-phase application overview with details on alkane and benzene compounds [J].TheJournalofPhysicalChemistryB,1998,102(38):7338-7364.

[14] Andersen H C.Molecular dynamics simulations at constant pressure and/or temperature [J].JournalofChemicalPhysics,1980,72(4):2374-2383.

[15] Tosi M P.Cohesion of ionic solids in the Born model [J].SolidStatePhysics,1964,16:1-120.

[16] Ewald P P.Evaluation of optical and electrostatic lattice potentials [J].AnnalsofPhysic,1921,369(3):253-287.

[17] Dong S J,Yan J T,Xu N,et al.Molecular dynamics simulation on surface modification of carbon black with polyvinyl alcohol [J].SurfaceScience,2011,605(9-10)：868-874.

[18] Swenson R J.Comments on virial theorems for bounded systems [J].AmericanJournalofPhysics,1983,51(10):940-942.

[19] 魏慶華,汪焰恩,楊明明,等.HA/PLA復(fù)合材料界面相互作用及其力學(xué)性能的MD模擬[J].功能材料,2013,44(21):3089-3093.(WEI Qing-hua,WANG Yan-en,YANG Ming-ming,et al.Molecular dynamics simulation of mechanical properties and surface inte -raction for HA/PLA [J].ChineseFunctionalMaterials,2013,44(21):3089-3093.(in Chinese))