一种基于materials studio的聚乳酸复合体系玻璃化转变温度预测方法
技术领域
本发明属于高分子材料模拟
技术领域
，具体涉及一种基于materials studio的聚乳酸复合体系玻璃化转变温度预测方法。背景技术
随着人们环保意识的不断增强，能循环利用的可降解材料研究和应用在环境治理方面越来越受到关注。利用天然植物纤维资源制备可降解复合材料已成为国内外科学家研究的热点。
聚乳酸（PLA）作为人工合成的热塑性聚酯材料，在应用价值方面是生物可降解材料中最具潜力的一种，但聚乳酸本身的耐热性较差，针对该方面的不足，需要改善其耐热性，使其能够应用到工作温度较广的工况环境中，扩大其应用范围。当前有关聚乳酸性能的改善研究多数都是试验研究，耗时耗力且成本较高。随着计算机技术的高速发展，借助计算机模拟对开发新材料的实验方案辅助研究的方法已被广泛应用。
分子动力学模拟有着从分子尺度模拟材料性质的独特优势，还能克服传统试验对样品表面难以处理的困难，已成为材料领域的重要研究手段。分子动力学模拟能从分子尺度解释物质性质，还弥补了实验的短处，能够基于宏观层面的研究成果对现有和未出现的材料进行性能预测，得到实验难以获得的结果，从而指导实验。
发明内容
本发明的目的在于提供一种基于materials studio的聚乳酸复合体系玻璃化转变温度预测方法，该方法有利于便捷、准确地预测聚乳酸复合体系的玻璃化转变温度，降低实验成本，缩短研发周期。
为实现上述目的，本发明采用的技术方案是：一种基于materials studio的聚乳酸复合体系玻璃化转变温度预测方法，包括以下步骤：
S1）基于蒙特卡洛方法搭建聚乳酸体系无定型晶胞模型；
S2）对所构建的无定型晶胞模型进行几何优化和退火模拟，得到能量最低构型；
S3）对得到的能量最低构型进行动力学计算，得到聚乳酸体系的降温模拟初始平衡构型；
S4）对平衡构型进行降温模拟，得到各个降温温阶下的密度和自由体积分数，线性拟合密度-温度关系得到玻璃化转变温度预测值，依据自由体积分数分析聚乳酸体系的玻璃化转变温度变化趋势。
进一步地，所述步骤S1中，在构建无定型晶胞模型时，预先设置分子链的聚合度、数量、温度以及晶胞模型的初始参数；为了消除非键合碳原子的影响，对构建的碳纳米管两端进行加氢处理；几何优化时，预先设置优化算法，相关力、能量和位移的最大收敛容差，以及几何优化过程最大迭代次数。
进一步地，所述步骤S2中，搭建的无定型晶胞模型并不处于能量最低状态，需要进行结构优化和退火模拟，搜寻全局条件下的能量最低状态；退火模拟时，预先设置退火过程中的温度最高值和最低值、循环次数、模拟时长参数。
进一步地，所述步骤S3中，对能量最低构型进行NVT系综下的动力学计算充分驰豫结构和NPT系综下的动力学计算最优结构，使模拟体系接近真实条件下的材料。
进一步地，所述步骤S4中，在降温模拟过程中，预先设置模拟体系的初始降温温度和最终降温截止温度、每级降温的温差、每个温阶在NVT和NPT系综下动力学计算的时长。
进一步地，每级温阶下的模拟体系密度为NPT动力学计算最后与设定温度相近的n个帧文件的平均值；自由体积分数由平衡后的构型生成Connoly表面划分自由体积和占有体积计算得到。
进一步地，得到每级温阶下的NPT计算后的模拟体系密度和自由体积分数；每一级降温后的构型作为下一级的初始构型进行计算；在得到所有降温温阶下的密度后，得到模拟体系密度-温度关系，通过分段线性拟合得到模拟体系的玻璃化转变温度预测值。
进一步地，基于以下公式计算聚乳酸复合体系每级温阶的自由体积分数：
其中，FFV是自由体积分数，V free 是自由体积，V occupied 是占有体积。
与现有技术相比，本发明具有以下有益效果：提供了一种基于materials studio的聚乳酸复合体系玻璃化转变温度预测方法，能够构建出精确的复合材料模型，预测出材料的玻璃化转变温度，通过将性能模拟结果参数与相关实验文献对比，验证模型的构建与性能预测方法的正确性；基于计算机模拟的优势，可以降低实验成本，缩短研发周期。同时，本发明通过对模拟设置合适的力场参数和优化方案，使整个模拟体系接近于真实状态下的物质体系。此外，本发明可通过改变分子链的长度、碳纳米管的状态、模拟温度条件等参数得到不同聚乳酸模拟体系的玻璃化转变温度预测，可以更直观表现出对聚乳酸体系玻璃化转变温度产生影响的因素。
附图说明
图1是本发明实施例的方法实现流程图。
图2是本发明实施例中聚乳酸单体分子和优化后的聚乳酸分子链和碳纳米管模型图。
图3是本发明实施例中动力学平衡后聚乳酸体系模型图。
图4是本发明实施例中降温模拟流程图。
图5是本发明实施例中聚乳酸体系计算自由体积分数的Connoly表面模型。
图6是本发明实施例中聚乳酸模拟体系的密度-温度曲线。
图7是本发明实施例中聚乳酸体系的自由体积分数-温度曲线。
具体实施方式
下面结合附图及具体实施例对本发明作进一步的详细说明。
如图1所示，本实施例提供了一种基于materials studio的聚乳酸复合体系玻璃化转变温度预测方法，包括以下步骤：
S1）基于蒙特卡洛方法搭建聚乳酸体系无定型晶胞模型。
在构建无定型晶胞模型时，需要预先设置分子链的聚合度、数量、温度以及晶胞模型的初始参数。为了消除非键合碳原子的影响，对构建的碳纳米管两端进行加氢处理。几何优化时，需要预先设置优化算法，相关力、能量和位移的最大收敛容差，以及几何优化过程最大迭代次数。
所述步骤S1具体包括以下步骤：
首先构建聚乳酸单体分子模型如图2(a)，将建好的单体分子模型初步优化后通过build polymers构造出聚合度为50的聚乳酸分子链，之后将构建好的聚乳酸分子链进行几何结构优化，优化后的结构如图2(b)。
借助materials studio中的Build Nanostructure模块构建出（6,6）型单壁碳纳米管，键长设置为1.42Å，通过改变碳纳米管的轴向长度参数，构建出不同长度碳纳米管模型。本实施例采用碳纳米管长度为24.6Å，如图2(c)。为防止非键合碳原子的影响，通过Adjust Hydrogen对碳纳米管模型进行加氢处理。
借助materials studio平台的Amorphous Cell模块按照蒙特卡洛方法加载20条PLA分子链和碳纳米管构建无定型晶胞模型；在搭建含有碳纳米管的晶胞模型时，需将碳纳米管的质心与晶胞体心重合，碳纳米管的轴向与x正向方向重合；立方晶胞边长参数设置为50Å。
结构优化参数设置为：优化算法选择Smart，几何优化收敛精度为能量收敛容差2×10-5kcal/mol，力的收敛容差是1×10-3kcal/mol/Å，不施加外部压强，在优化原子坐标的同时，晶格参数也会得到优化；力场选用 COMPASSⅡ力场，电荷由力场分配，库伦作用与范德华作用分别选用Ewald算法与Atom based算法，由于库伦势是长程势，所以选择Ewald算法合适。
本实施例共构建三种聚乳酸体系，分别是：纯聚乳酸体系（PLA），未官能化碳纳米管/聚乳酸体系（p-CNTs/PLA），氨基化碳纳米管/聚乳酸体系（CNTs-NH2/PLA）。氨基化碳纳米管是以碳纳米管的碳原子数的6%的数量为接枝位点接枝-NH2。
S2）对所构建的无定型晶胞模型进行几何优化和退火模拟，得到能量最低构型。
初步建立的无定型晶胞模型并不处于能量最低状态，需要达到能量最低状态，得到能量最低构型才能进行动力学计算。而要达到能量最低状态，需要借助几何结构优化和退火模拟，通过结构优化和退火模拟，搜寻全局条件下的能量最低状态。
体系几何结构优化参数设置为：力场选择COMPASS II力场，能量和力的收敛容差分别为2×10-5kcal/mol和1×10-3kcal/mol/Å，几何优化最大迭代次数设置为20000。
在对模拟几何优化初步处理之后，为了防止模拟体系陷入局部能量的最低点，需要借助退火模拟Anneal模块搜寻全局能量最低点状态，这一步通过施加外界温度变化来实现。
退火模拟时，需要预先设置退火过程中的温度最高值和最低值、循环次数、模拟时长等参数。
退火模拟的具体参数设置为：预设置COMPASS II力场，计算库伦相互作用和范德华相互作用分别用Ewald求和方法和Atom based求和方法；在NVT系综下设置退火循环次数为10，初始温度为200K，最高温度为400K，循环设置20个温阶，每个温阶的动力学计算步数为1000，步长为1fs，温控选择Andersen方法，总退火时长为400ps；退火模拟是允许模拟体系在升温与降温过程中克服能量壁垒进入其他能量低状态区域。
S3）对得到的能量最低构型进行动力学计算，得到聚乳酸体系的降温模拟初始平衡构型。
在经过几何结构优化和退火模拟后，为了使模拟体系接近真实情况下的材料体系，需要进行动力学计算。对能量最低的聚乳酸复合体系模型的动力学计算分为两个过程：NVT系综下驰豫结构和NPT系综下最优结构。
其中，模拟体系NVT系综下驰豫结构的动力学计算参数设置为：预设置COMPASS II力场，计算库伦相互作用和范德华相互作用分别用Ewald求和方法和Atom based求和方法；温度设置298K，步长为1fs，总时长为200ps，每计算5000步输出一帧结构文件；温控选择Andersen方法。
NPT动力学计算的参数设置为：预设置COMPASS II力场，计算库伦相互作用和范德华相互作用分别用Ewald求和方法和Atom based求和方法；温度设置298K，步长为1fs，压强为1.01×10-4GPa，总时长为400ps，每计算5000步输出一帧结构文件；温控压控分别选择Andersen方法和Berendsen方法；最后100ps计算过程温度波动不超过10K，能量波动不超过总能量值的3%，可认为结构达到平衡状态。
动力学计算完成后得到的纯聚乳酸体系（PLA）如图3（a），未官能化碳纳米管/聚乳酸体系（p-CNTs/PLA）如图3（b），氨基化碳纳米管/聚乳酸体系（CNTs-NH2/PLA）如图3（c）。
S4）对平衡构型进行降温模拟，得到各个降温温阶下的密度和自由体积分数，线性拟合密度-温度关系得到玻璃化转变温度预测值，依据自由体积分数分析聚乳酸体系的玻璃化转变温度变化趋势。
在降温模拟过程中，需要预先设置模拟体系的初始降温温度和最终降温截止温度、每级降温的温差、每个温阶在NVT和NPT系综下动力学计算的时长。
每级温阶下的模拟体系密度为NPT动力学计算最后与设定温度相近的n个帧文件的平均值；自由体积分数由平衡后的构型生成Connoly表面划分自由体积和占有体积计算得到。
得到每级温阶下的NPT计算后的模拟体系密度和自由体积分数；每一级降温后的构型作为下一级的初始构型进行计算；在得到所有降温温阶下的密度后，得到模拟体系密度-温度关系，通过分段线性拟合得到模拟体系的玻璃化转变温度预测值。
具体的，借助Dynamics模块对模拟体系进行降温模拟，模拟过程如图4所示。
首先采用动力学计算平衡后的构型作为降温模拟初始的模型，将初始模拟升温到420K；借助MS平台的Forcite模块中动力学计算模块，在NVT系综下将初始模型升温到420K，以此获取模型在420K温度状态下的参数，之后依次进行50ps的NVT动力学计算和200ps的NPT动力学计算。
然后将模型通过逐步降温冷却到280K，降温冷却速度为10K/250ps，即先进行50ps的NVT系综下的动力学计算，再保持温度恒定进行200ps下的NPT下的动力学计算，总共模拟15个温度，从420K到280K以10K为间隔将模型逐步降温，获取每个温度阶段下NPT平衡之后的密度值。
为了减少密度取值的误差，每个降温温阶的密度值取NPT平衡最后与设定温度最近的5个帧文件密度的平均值。
降温模拟过程的参数设置为：选用COMPASS II力场，压强设置为1个大气压强，控温方法为Andersen，控压方法为Berendsen，每一步都从上一步降温计算完成的构型开始计算。
模拟体系降温的各个温阶的自由体积分数通过对模拟体系构建Connoly表面，计算模拟体系中自由体积和占用体积，自由体积分数的计算公式如下：
其中，FFV是自由体积分数，V free 是自由体积，V occupied 是占有体积。Connoly表面模型如图5所示。
降温模拟得到各个温阶下的密度值，通过分段拟合密度-温度曲线，两条拟合线的交点的横坐标即为模拟体系的Tg值。
本实施例经过降温模拟后，提取得到模拟体系密度-温度曲线；纯聚乳酸体系（PLA）如图6（a），未官能化碳纳米管/聚乳酸体系（p-CNTs/PLA）如图6（b），氨基化碳纳米管/聚乳酸体系（CNTs-NH2/PLA）如图6（c）；聚乳酸体系的自由体积分数-温度关系如图7。
本发明提供的基于materials studio的聚乳酸复合体系玻璃化转变温度预测方法，能够构建出精确的复合材料模型，通过将性能模拟结果参数与相关实验文献对比，验证模型的构建与性能预测方法的正确。基于计算机模拟的优势，可以降低实验成本，缩短研发周期。本发明可通过改变分子链的长度和碳纳米管的状态，预测待研究聚乳酸体系的玻璃化转变温度，借助模型直观的表现出引起聚乳酸体系玻璃化转变温度的微观因素。本发明还可以进一步借助模型从界面状态分析界面相互作用，对材料的性质进行深一步的探索。
以上详细描述本发明较佳的具体实例的细节和操作内容。本领域的技术人员会理解，但本发明并不局限于以上所述的实例。对本领域的技术人员来说，可以很熟练的对参数细节以及操作步骤顺序进行各种修改和变化，在不脱离本发明构思的情况下，在本技术领域中的技术人员依据理论修改和相当的实践可以得到的技术方案，均应该在本权利要求书所确定的保护范围之中。
