CN116230104B - 基于分子动力学和多变量分析的低共熔溶剂筛选方法 - Google Patents

Abstract

本发明属于木质纤维素生物质高值化利用领域，涉及一种基于分子动力学和多变量分析的低共熔溶剂筛选方法，包括以下步骤：构建半纤维素分子、乳酸分子、每一种氢键受体阳离子、氢键受体阴离子的结构模型，进行结构优化和能量优化，计算RESP电荷并导入拓扑文件中；借助分子动力学模拟，通过氢键数量差值和相互作用能来定量分析不同氢键受体的低共熔溶剂对半纤维素的溶解机理及效果；采用多变量分析确定分子动力学模拟结果与不同低共熔溶剂溶解半纤维素所有变量之间的相互关性，从而在不同季铵盐类氢键受体中筛选出溶解半纤维素效果最好的低共熔溶剂，该方法为不同氢键受体低共熔溶剂溶解半纤维素提供了理论基础。

Description

基于分子动力学和多变量分析的低共熔溶剂筛选方法

技术领域

本发明属于木质纤维素生物质高值化利用领域，具体涉及一种基于分子动力学和多变量分析的低共熔溶剂筛选方法。

背景技术

木质纤维素生物质是最丰富的自然资源之一，作为一种廉价的可再生能源，常年用于绿色能源生产，能够解决全球不断增长地对能源，特别是对石油基燃料的需求。木质纤维素生物质的主要成分是纤维素、半纤维素和木质素，其中半纤维素占15%-35%。半纤维素是一种复杂的支链杂聚物，主要结构是β-1,4连接的糖单元的直链骨架，其主链的长度和类型，侧链的分布和类型随木质纤维素的种类而变化，木聚糖是半纤维素中含量最丰富的一种多糖聚合物。半纤维素的含量丰富，仅次于纤维素，在食品、化妆品和造纸方面具有重要应用，此外，由于其独特的理化性质和功能，也被广泛应用于制药领域，具有广阔的应用前景。但其复杂的结构和多样化的性质，以及木质纤维素生物质中其他组分（纤维素和木质素）之间的相互作用，使得半纤维素的高值化利用在长期以来具有挑战性。

利用不同的预处理技术，能够打破木质纤维素生物质的顽固结构，将其分解成纤维素、半纤维素和木质素，但是传统的预处理方法具有耗时大、效率低等缺点，无法满足工业规模化的可持续性和实用性。此外，由于预处理技术中采用的离子液体预处理的成本高、有毒，严重阻碍了其在木质纤维素生物质中的应用，因此，需要寻求一种低成本、无毒的绿色溶剂处理技术。

低共熔溶剂作为一种新型的绿色溶剂在众多的生物质预处理方法中脱颖而出，目前报道的关于低共熔溶剂预处理的研究中，大多数的低共熔溶剂都是以氯化胆碱作为氢键受体合成的。但是，目前研究的由不同氢键受体组成的低共熔溶剂对半纤维素的溶解效果的影响主要集中在实验和表征分析上，缺乏相关的理论基础。因此，由不同氢键受体组成的低共熔溶剂溶解半纤维素的机理还在深入探究中。

发明内容

本发明为了克服不同氢键受体低共熔溶剂溶解半纤维素的实验方法带来的耗时耗力且成本高等诸多缺点，提出了一种基于分子动力学和多变量分析的低共熔溶剂筛选方法，采用分子动力学和多变量分析方法研究由不同氢键受体组成的低共熔溶剂对半纤维素的溶解机理以及溶解效果，筛选不同氢键受体和摩尔比的低共熔溶剂；该方法能够快速高效地从原子和分子尺度模拟计算出低共熔溶剂-半纤维素体系的详细信息，直观展现出不同氢键受体的性质对半纤维素作用效果的影响程度，弥补了实验表征在机理分析方面的不足，为不同氢键受体的低共熔溶剂在木质纤维素预处理领域的高值化利用提供了理论基础。

本发明的技术方案是：

一种基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，包括以下步骤：

（1）选择乳酸作为低共熔溶剂的氢键供体，季铵盐类作为低共熔溶剂的氢键受体，构建乳酸分子、每一种氢键受体阳离子、氢键受体阴离子的结构模型；并构建半纤维素分子的结构模型；

（2）将步骤（1）构建的所有结构模型进行结构优化和能量优化，得到分子静电势信息，计算相应的标准静电势（RESP）电荷；

（3）将步骤（2）优化后的结构模型通过CGenFF得到CHARMM36力场分子拓扑文件，并将计算得到的RESP电荷导入拓扑文件中；

其中，CHARMM General Force Field（CGenFF）程序是专门用于生成CGenFF拓扑文件以及CHARMM36力场参数的工具，通过分类确定原子类型、参数和电荷；

CHARMM36为分子模拟常用力场；

（4）构建分子动力学初始结构盒子，将200~1000个氢键受体分子，500~2000个乳酸分子和20~100条半纤维素分子放入盒子中，其中，将半纤维素分子放在盒子中间位置，氢键受体分子和乳酸分子放在盒子四周位置，建立多个稳定的不同低共熔溶剂-半纤维素体系；

（5）运行不同低共熔溶剂-半纤维素体系，在298.15 K和1 bar压力下对低共熔溶剂-半纤维素体系采用最陡下降法模拟10000~15000步使系统能量达到最小，做限制性动力学，利用PME方法计算截止半径为1.0~1.8 nm的静电相互作用，在363.15~413.15K处进行步长为25000000~50000000步、持续时间为50~100 ns的长时间动力学模拟，得到每一帧下的模拟轨迹和结果；

其中，PME（Particle Mesh Ewald）方法是分子动力学模拟中采用的一种优化算法。

（6）将步骤（5）模拟产生的每一帧轨迹结果导入可视化软件中，将得到的轨迹文件中的半纤维素的羟基进行分组，在分子动力学仿真计算程序Gromacs中运行氢键命令和相互作用能命令，计算每类羟基分别与每类低共熔溶剂的乳酸分子、氢键受体阳离子、氯离子之间的氢键数量和相互作用能，以及半纤维素分子之间的氢键数量，通过分析半纤维素与每类低共熔溶剂相互作用后的半纤维素分子之间的氢键数量H₂和半纤维素没有与低共熔溶剂相互作用的半纤维素分子之间的氢键数量H₁的差值ΔH，对低共熔溶剂进行筛选；

每类羟基指的是半纤维素分子上不同位置的羟基；步骤（6）的目的是，将模拟产生的每一帧轨迹结果导入分子可视化软件Visual molecular dynamics (VMD)中，然后在分子动力学仿真计算程序Gromacs中运行氢键命令和相互作用能命令，进行氢键和相互作用能的分析，通过计算氢键数量和相互作用能的大小，探究不同氢键受体的低共熔溶剂对半纤维素的溶解机理以及对比溶解效果，进而进行低共熔溶剂的筛选。

（7）将步骤（6）中计算得到的氢键数量差值ΔH和相互作用能导入到多变量分析软件SIMCA-P中，在多变量分析软件SIMCA-P中，对由不同官能团或不同碳链长度的氢键受体组成的低共熔溶剂与半纤维素之间的相互作用的变量，采用主成分分析和偏最小二乘法进行多变量分析；其中，对变量的主成分分析得到主成分的解释方差为0.911~0.998，预测方差为0.932~0.968；偏最小二乘分析得到的模型能够解释95.8%~99.8%的因变量变化，预测能力为80.1%~99.5%；根据主成分分析和偏最小二乘分析的模型得到氢键受体的官能团、碳链长度与氢键数量差值ΔH之间存在强烈正相关，氢键数量差值ΔH与相互作用能之间存在强烈负相关，说明碳链长度越长其氢键数量差值ΔH越大，从而快速地筛选出对半纤维素的溶解效果最好的低共熔溶剂；最后结合偏最小二乘法得到的变量重要性分析，得出低共熔溶剂与半纤维素之间形成的氢键数量越多，说明低共熔溶剂溶解半纤维能力越强。

进一步的，所述步骤（1）中，低共熔溶剂的氢键供体与氢键受体的摩尔比为HBA：HBD =1：（5~9）；

所述低共熔溶剂的氢键受体包括盐酸胍、氯化胆碱、甜菜碱、乙胺氯、烯丙基三甲基氯化铵、四甲基氯化铵、四乙基氯化铵、四甲基溴化铵、苄基三甲基氯化铵、苄基三乙基氯化铵、苄基二甲基氯化铵中的任一种或几种。

进一步的，所述步骤（1）中半纤维素的结构模型包含6~8个五碳糖组成的聚合物结构、或者包含2~7个六碳糖组成的聚合物结构、或者包含6~8个五碳糖和2~7个六碳糖组成的聚合物结构。

进一步的，所述步骤（2）中，将步骤（1）构建的所有结构模型保存为gjf文件，在量子化学仿真软件Gaussian中基于密度泛函理论进行结构优化和能量优化，得到各自的分子静电势信息，结合量子化学波函数分析程序Multiwfn计算相应的RESP电荷。

进一步的，所述步骤（4）中，采用分子动力学初始结构构建程序Packmol构建10 nm×10 nm×16 nm的分子动力学初始结构盒子。

进一步的，所述步骤（5）中，在分子动力学仿真计算程序Gromacs中运行步骤（4）建立的低共熔溶剂-半纤维素体系；当能量达到最小后逐渐加热到363.15~413.15 K，在达到所需温度后，能量最小的低共熔溶剂-半纤维素体系在NVT体系（等温等压体系）合成中对半纤维素产生限制性动力学。

进一步的，所述步骤（6）中采用的分子可视化软件为Visual molecular dynamics(VMD)。

进一步的，所述步骤（7）中的变量包括不同氢键受体的物理化学性质和分子模拟结果；

其中，不同氢键受体的物理化学性质包括：熔点、分子量、分子极性指数、氨基数量、羟基数量、碳碳双键数量、碳链长度、氢键供体数量、氢键受体数量、氢键供体与氢键受体数量之和、可旋转化学键数量、重原子数量、分子极性表面积和复杂度；分子模拟结果包括：氢键数量差值、氢键寿命、相互作用能、径向分布函数和空间分布函数。

本发明的有益效果：

（1）本发明提供的筛选方法，借助分子动力学模拟，通过氢键数量差值和相互作用能来定量分析不同氢键受体的低共熔溶剂对半纤维素的溶解机理及效果；采用多变量分析确定分子动力学模拟结果与不同低共熔溶剂溶解半纤维素所有变量之间的相互关系，氢键受体的官能团与氢键数量差值之间存在强烈正相关，氢键数量差值与相互作用能之间存在强烈负相关，从而在不同季铵盐类氢键受体中筛选出溶解半纤维素效果最好的低共熔溶剂。

（2）本发明基于分子动力学和多变量分析准确地分析了不同氢键受体的低共熔溶剂对半纤维素的溶解机理，可以提供低共熔溶剂-半纤维素体系含有的变量之间详细的相互关系，克服了传统实验费时费力进行评价低共熔溶剂溶解半纤维素效果的缺点，能够快速地筛选出对半纤维素溶解效果最好的最优氢键受体低共熔溶剂。

（3）本发明的方法不仅能快速高效地筛选最优低共熔溶剂，还能研究不同低共熔溶剂作用于半纤维素的作用机理，为不同氢键受体低共熔溶剂溶解半纤维素提供了理论基础，可为开发新型不同氢键受体可溶半纤维素低共熔溶剂提供理论依据和指导。

附图说明

图1为实施例1和2提供的低共熔溶剂-半纤维素的氢键数量柱形图，包括氯离子、氢键供体阳离子、乳酸与半纤维素之间的氢键数量；

图2为实施例1和2提供的低共熔溶剂-半纤维素的相互作用能的示意图，包括氯离子、氢键供体阳离子、乳酸与半纤维素之间的静电相互作用能；

图3为实施例1提供的低共熔溶剂-半纤维素在进行分子模拟的相互作用的过程示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

为了进一步理解本发明，将结合附图和实施例对本发明作进一步的说明。

实施例1

本实施例提供了一种基于分子动力学和多变量分析的低共熔溶剂筛选方法，包括以下步骤：

（1）选择乳酸作为低共熔溶剂的氢键供体（HBD），含不同官能团的盐酸胍、氯化胆碱、甜菜碱作为低共熔溶剂的氢键受体（HBA），其摩尔比为HBA：HBD=1：5；在仿真计算软件Gaussian View中构建乳酸分子、盐酸胍阳离子、氯化胆碱阳离子、甜菜碱阳离子和氯离子的结构模型；

在仿真计算软件Gaussian View中构建半纤维素分子的结构模型；选择6个五碳糖组成的聚合物结构作为半纤维素的模型；

（2）将步骤（1）构建的半纤维素分子、乳酸分子、盐酸胍阳离子、氯化胆碱阳离子、甜菜碱阳离子和氯离子的结构模型文件保存为gjf文件，在仿真软件Gaussian中基于密度泛函理论进行结构优化和能量优化，得到各自的分子静电势信息，结合Multiwfn程序计算相应的Restrained Electrostatic Potential (RESP)电荷；

（3）将步骤（2）优化后的文件保存为mol2格式，上传至CGenFF；在CGenFF软件的More Info&Tools输入模式下生成基于CHARMM36力场的itp拓扑文件，并将计算出的RESP电荷复制导入到相应的itp拓扑文件中；

（4）采用Packmol程序构建10 nm×10 nm×16 nm的分子动力学初始结构盒子，分别将500个氢键受体分子，1800个乳酸分子和30条半纤维素分子放入盒子中；其中，将半纤维素分子放在盒子中间位置，氢键受体分子和乳酸分子放在盒子四周位置（如图3所示），建立多个稳定的不同低共熔溶剂-半纤维素体系；

（5）对步骤（4）建立的体系模型添加边界周期性条件，在分子动力学仿真计算程序Gromacs中运行不同的低共熔溶剂-半纤维素稳定体系，在298.15 K和1 bar压力下对低共熔溶剂-半纤维素体系采用最陡下降法模拟12000步使系统能量达到最小，能量最小后逐渐加热到373.15 K，在达到所需温度后，能量最小的体系在NVT体系合成中对半纤维素产生限制性动力学，然后利用PME方法计算截止半径为1.8 nm的静电相互作用，在373.15 K处进行步长为40,000,000步、持续时间为80 ns的长时间动力学模拟，得到每一帧下的模拟轨迹和结果；

（6）将步骤（5）模拟产生的每一帧轨迹结果导入可视化软件Visual moleculardynamics (VMD)中，将得到的tpr轨迹文件中的半纤维素的羟基进行分组，在软件Gromacs中运行氢键命令和相互作用能命令，计算半纤维素分子上不同位置的羟基分别与每类低共熔溶剂的乳酸分子、盐酸胍阳离子、氯化胆碱阳离子、甜菜碱阳离子、氯离子之间的氢键数量和相互作用能，以及半纤维素之间的氢键数量，通过分析半纤维素与盐酸胍、氯化胆碱或甜菜碱之间的氢键数量H₂和半纤维素分子之间的氢键数量H₁的差值ΔH得到如下结论：

含氨基的盐酸胍/乳酸-半纤维素体系的氢键数量和半纤维素分子之间的氢键数量的差值最大，说明盐酸胍/乳酸低共熔溶剂与半纤维素之间形成了强大的氢键作用，更大程度的破坏了半纤维素之间的氢键，能够更好地溶解半纤维素，通过计算相互作用能得到盐酸胍/乳酸-半纤维素体系的总相互作用能E_total的绝对值最高（如表1所示），此外，含氨基的氢键受体阳离子与半纤维素间主要是静电相互作用，说明由含氨基的氢键受体合成的低共熔溶剂以及含碳链越长的氢键受体合成的低共熔溶剂对半纤维素具有更好的溶解效果，进而达到对低共熔溶剂的筛选；

在本实施例对盐酸胍/乳酸-半纤维素体系的氢键分析中，盐酸胍阳离子与半纤维素羟基间形成了强相互作用。如图1所示，盐酸胍/乳酸与半纤维素形成的氢键总数最多，说明以盐酸胍为氢键受体合成的低共熔溶剂对半纤维素具有更好的溶解效果，其次是氯化胆碱/乳酸低共熔溶剂；

如图2所示，在盐酸胍/乳酸-半纤维素体系的相互作用能中，氯离子与半纤维素之间的静电相互作用（Eelec）比范德华相互作用（EvdW）更负，说明静电相互作用是氯离子-半纤维素相互作用的控制参数，且盐酸胍阳离子与半纤维素的相互作用以静电相互作用为主；

表1 半纤维素与不同氢键受体低共熔溶剂的总相互作用能

（7）将步骤（6）中计算得到的氢键数量差值ΔH和相互作用能导入到多变量分析软件SIMCA-P中，将含不同官能团的盐酸胍、氯化胆碱、甜菜碱的性质以及模拟结果作为变量在SIMCA-P软件中采用主成分分析和偏最小二乘法作多变量分析，过程中考虑的变量包括不同氢键受体的物理化学性质和分子模拟结果；不同氢键受体的物理化学性质包括：熔点、分子量、分子极性指数、氨基数量、羟基数量、碳碳双键数量、碳链长度、氢键供体数量、氢键受体数量、氢键供体与氢键受体数量之和、可旋转化学键数量、重原子数量、分子极性表面积和复杂度；分子模拟结果包括：氢键数量差值、氢键寿命、相互作用能、径向分布函数和空间分布函数；

对所有变量进行主成分分析，得到主成分的解释方差和预测方差分别为0.995和0.942，说明主成分分析模型的结果足以可靠地解释分子动力学模拟结果与不同低共熔溶剂溶解半纤维素之间的相关性；将不同氢键受体的物理化学性质设为自变量，分子模拟的氢键数量差值或静电相互作用能变量设为因变量作偏最小二乘分析，得到的模型能解释95.6%的因变量变化，预测能力为90.5%，说明该模型具有较好的评价能力，从而证明分子模拟的氢键数量差值ΔH或静电相互作用能结果能够准确地解释低共熔溶剂与半纤维素之间的作用机理，以及溶解半纤维素的效果；结合变量的重要性分析，发现低共熔溶剂与半纤维素之间形成的氢键数量越多，说明低共熔溶剂溶解半纤维能力越强；

根据主成分分析和偏最小二乘分析模型，得到氢键受体的官能团与氢键数量差值之间存在强烈正相关，氢键数量差值ΔH与相互作用能之间存在强烈负相关，说明碳链长度越长其氢键数量差值ΔH越大，从而快速地筛选出烷基链长的氢键受体，例如可快速筛选出以盐酸胍作为氢键受体合成的低共熔溶剂，对半纤维素的溶解效果最好。

实施例2

本实施例提供了一种基于分子动力学和多变量分析的低共熔溶剂筛选方法，包括以下步骤：

（1）选择乳酸作为低共熔溶剂的氢键供体（HBD），烯丙基三甲基氯化铵、四甲基氯化铵、乙胺氯、四甲基溴化铵、四乙基氯化铵作为低共熔溶剂的氢键受体（HBA），其摩尔比为HBA：HBD=1：9；在仿真计算软件Gaussian View中构建乳酸分子、烯丙基三甲基铵根阳离子、四甲基铵根阳离子、乙铵根阳离子、溴离子、四乙基铵根阳离子和氯离子的结构模型；

在仿真计算软件Gaussian View中构建半纤维素分子的结构模型；选择8个五碳糖和4个六碳糖组成的聚合物结构作为半纤维素的模型；

（2）将步骤（1）构建的半纤维素模型分子、乳酸分子、烯丙基三甲基铵根阳离子、四甲基铵根阳离子、乙铵根阳离子、溴离子、四乙基铵根阳离子和氯离子的结构模型文件保存为gjf文件，在仿真软件Gaussian中基于密度泛函理论进行结构优化和能量优化，得到各自的分子静电势信息，结合Multiwfn程序计算相应的Restrained Electrostatic Potential(RESP)电荷；

（3）将步骤（2）优化后的文件保存为mol2格式，上传至CGenFF；在CGenFF软件的More Info&Tools输入模式下生成基于CHARMM36力场的itp拓扑文件，并将计算出的RESP电荷复制导入到相应的itp拓扑文件中；

（4）采用Packmol程序构建10 nm×10 nm×16 nm的分子动力学初始结构盒子，分别将200个氢键受体分子，500个乳酸分子和20条半纤维素分子放入盒子中；其中，将半纤维素模型分子放在盒子中间位置，氢键受体分子和乳酸分子放在盒子四周位置，建立多个稳定的不同低共熔溶剂-半纤维素体系；

（5）对步骤（4）建立的体系模型添加边界周期性条件，在分子动力学仿真计算程序Gromacs中运行不同低共熔溶剂-半纤维素体系稳定体系，在298.15 K和1 bar压力下对低共熔溶剂-半纤维素体系采用最陡下降法模拟10000步使系统能量达到最小，能量最小后逐渐加热到383.15 K，在达到所需温度后，能量最小的体系在NVT体系合成中对半纤维素产生限制性动力学，然后利用PME方法计算截止半径为1.0 nm的静电相互作用，在383.15 K处进行步长为25,000,000步、持续时间为50 ns的长时间动力学模拟，得到每一帧下的模拟轨迹和结果；

（6）将步骤（5）模拟产生的每一帧轨迹结果导入可视化软件Visual moleculardynamics (VMD)中，将得到的tpr轨迹文件中的半纤维素的羟基进行分组，在软件Gromacs中运行氢键命令和相互作用能命令，计算每类羟基分别与每类低共熔溶剂的乳酸分子、烯丙基三甲基铵根阳离子、四甲基铵根阳离子、乙铵根阳离子、溴离子、四乙基铵根阳离子、氯离子之间的氢键数量和相互作用能，以及半纤维素之间的氢键数量，通过分析半纤维素与烯丙基三甲基氯化铵、四甲基氯化铵、乙胺氯、四甲基溴化铵或四乙基氯化铵之间的氢键数量H₂和半纤维素分子之间的氢键数量H₁的差值ΔH得到如下结论：

碳链长度越长的烯丙基三甲基氯化铵/乳酸-半纤维素体系的氢键数量和半纤维素分子之间的氢键数量的差值最大，碳链长度越长的烯丙基三甲基氯化铵/乳酸-半纤维素体系，低共熔溶剂与半纤维素的相互作用能E_total的绝对值越大（如表2所示）。

在烯丙基三甲基氯化铵/乳酸-半纤维素体系的氢键分析中，氯离子与半纤维素羟基间形成了强氢键相互作用，如图1所示，氯离子与半纤维素形成的氢键总数最多，说明烯丙基三甲基氯化铵/乳酸中的氯离子更倾向于与半纤维素间形成强氢键作用；其次是四甲基氯化铵/乳酸低共熔溶剂；

如图2所示，在烯丙基三甲基氯化铵/乳酸-半纤维素体系的相互作用能中，氯离子与半纤维素之间的静电相互作用（Eelec）比范德华相互作用（EvdW）更负，说明静电相互作用是氯离子-半纤维素相互作用的控制参数，且烯丙基三甲基铵根阳离子与半纤维素的相互作用以范德华相互作用为主。

表2 半纤维素与不同氢键受体低共熔溶剂的总相互作用能

（7）将步骤（6）中计算得到的氢键数量差值ΔH和相互作用能导入到多变量分析软件SIMCA-P中，将不同烷基链长度的烯丙基三甲基氯化铵、四甲基氯化铵、乙胺氯、四甲基溴化铵、四乙基氯化铵的性质以及模拟结果作为变量在SIMCA-P软件中采用主成分分析和偏最小二乘法作多变量分析，过程中考虑的变量包括不同氢键受体的物理化学性质和分子模拟结果；不同氢键受体的物理化学性质包括：熔点、分子量、分子极性指数、氨基数量、羟基数量、碳碳双键数量、碳链长度、氢键供体数量、氢键受体数量、氢键供体与氢键受体数量之和、可旋转化学键数量、重原子数量、分子极性表面积和复杂度；分子模拟结果包括：氢键数量差值、氢键寿命、相互作用能、径向分布函数和空间分布函数；

对所有变量进行主成分分析，得到主成分的解释方差和预测方差分别为0.945和0.932，说明主成分分析模型的结果足以可靠地解释分子动力学模拟结果与不同低共熔溶剂溶解半纤维素之间的相关性；将不同氢键受体的物理化学性质设为自变量，分子模拟的氢键数量差值或静电相互作用能变量设为因变量作偏最小二乘分析，得到的模型能解释96.3%的因变量变化，预测能力为94.5%，说明该模型具有较好的评价能力，从而证明分子模拟的氢键数量差值ΔH或静电相互作用能结果能够准确地解释低共熔溶剂与半纤维素之间的作用机理，以及溶解半纤维素的效果；结合变量的重要性分析，发现低共熔溶剂与半纤维素之间形成的氢键数量越多，说明低共熔溶剂溶解半纤维能力越强；

根据主成分分析和偏最小二乘分析模型，得到氢键受体的碳链长度与氢键数量差值之间存在强烈正相关，氢键数量差值ΔH与相互作用能之间存在强烈负相关，说明碳链长度越长其氢键数量差值ΔH越大，从而快速地筛选出烷基链长的氢键受体，如烯丙基三甲基氯化铵作为氢键受体对半纤维素的溶解效果最好。

实施例3

本实施例提供了一种基于分子动力学和多变量分析的低共熔溶剂筛选方法，包括以下步骤：

（1）选择乳酸作为低共熔溶剂的氢键供体（HBD），苄基三甲基氯化铵、苄基三乙基氯化铵、苄基二甲基氯化铵作为低共熔溶剂的氢键受体（HBA），其摩尔比为HBA：HBD=1：7；在仿真计算软件Gaussian View中构建乳酸分子、苄基三甲基铵根阳离子、苄基三乙基铵根阳离子、苄基二甲基铵根阳离子和氯离子的结构模型；

在仿真计算软件Gaussian View中构建半纤维素分子的结构模型；选择8个五碳糖组成的聚合物结构作为半纤维素的模型；

（2）将步骤（1）构建的半纤维素分子、乳酸分子、苄基三甲基铵根阳离子、苄基三乙基铵根阳离子、苄基二甲基铵根阳离子和氯离子的结构模型文件保存为gjf文件，在仿真软件Gaussian中基于密度泛函理论进行结构优化和能量优化，得到各自的分子静电势信息，结合Multiwfn程序计算相应的Restrained Electrostatic Potential (RESP)电荷；

（3）将步骤（2）优化后的文件保存为mol2格式，上传至CGenFF；在CGenFF软件的More Info&Tools输入模式下生成基于CHARMM36力场的itp拓扑文件，并将计算出的RESP电荷复制导入到相应的itp拓扑文件中；

（4）采用Packmol程序构建10 nm×10 nm×16 nm的分子动力学初始结构盒子，分别将1000个氢键受体分子，2000个乳酸分子和100条半纤维素分子放入盒子中；其中，将半纤维素分子放在盒子中间位置，氢键受体分子和乳酸分子放在盒子四周位置，建立多个稳定的不同低共熔溶剂-半纤维素体系；

（5）对步骤（4）建立的体系模型添加边界周期性条件，在分子动力学仿真计算程序Gromacs中运行不同低共熔溶剂-半纤维素稳定体系，在298.15 K和1 bar压力下对低共熔溶剂-半纤维素体系采用最陡下降法模拟15000步使系统能量达到最小，能量最小后逐渐加热到383.15 K，在达到所需温度后，能量最小的体系在NVT体系合成中对半纤维素产生限制性动力学，然后利用PME方法计算截止半径为1.4 nm的静电相互作用，在413.15 K处进行步长为50000000步、持续时间为100 ns的长时间动力学模拟，得到每一帧下的模拟轨迹和结果；

（6）将步骤（5）模拟产生的每一帧轨迹结果导入可视化软件Visual moleculardynamics (VMD)中，将得到的tpr轨迹文件中的半纤维素的羟基进行分组，在软件Gromacs中运行氢键命令和相互作用能命令，计算每类羟基分别与每类低共熔溶剂的乳酸分子、苄基三甲基铵根阳离子、苄基三乙基铵根阳离子、苄基二甲基铵根阳离子、氯离子之间的氢键数量和相互作用能，以及半纤维素之间的氢键数量，通过分析半纤维素与苄基三甲基氯化铵、苄基三乙基氯化铵或苄基二甲基氯化铵之间的氢键数量H₂和半纤维素分子之间的氢键数量H₁的差值ΔH得到如下结论：

碳链长度越长的苄基三乙基氯化铵/乳酸-半纤维素体系的氢键数量和半纤维素分子之间的氢键数量的差值最大，碳链长度越长的苄基三乙基氯化铵/乳酸-半纤维素体系，低共熔溶剂与半纤维素的相互作用能Etotal的绝对值越大；

在苄基三乙基氯化铵/乳酸-半纤维素体系的氢键分析中，氯离子与半纤维素羟基间形成了强氢键相互作用，且氯离子与半纤维素形成的氢键总数最多，说明苄基三乙基氯化铵/乳酸中的氯离子更倾向于与半纤维素间形成强氢键作用；

（7）将步骤（6）中计算得到的氢键数量差值ΔH和相互作用能导入到多变量分析软件SIMCA-P中，将不同烷基链长度的苄基三甲基氯化铵、苄基三乙基氯化铵、苄基二甲基氯化铵的性质以及模拟结果作为变量在SIMCA-P软件中采用主成分分析和偏最小二乘法作多变量分析，过程中考虑的变量包括不同氢键受体的物理化学性质和分子模拟结果；不同氢键受体的物理化学性质包括：熔点、分子量、分子极性指数、氨基数量、羟基数量、碳碳双键数量、碳链长度、氢键供体数量、氢键受体数量、氢键供体与氢键受体数量之和、可旋转化学键数量、重原子数量、分子极性表面积和复杂度；分子模拟结果包括：氢键数量差值、氢键寿命、相互作用能、径向分布函数和空间分布函数；

对所有变量进行主成分分析，得到主成分的解释方差和预测方差分别为0.945和0.932，说明主成分分析模型的结果足以可靠地解释分子动力学模拟结果与不同低共熔溶剂溶解半纤维素之间的相关性；将不同氢键受体的物理化学性质设为自变量，分子模拟的氢键数量差值或静电相互作用能变量设为因变量作偏最小二乘分析，得到的模型能解释96.3%的因变量变化，预测能力为94.5%，说明该模型具有较好的评价能力，从而证明分子模拟的氢键数量差值ΔH或静电相互作用能结果能够准确地解释低共熔溶剂与半纤维素之间的作用机理，以及溶解半纤维素的效果；结合变量的重要性分析，发现低共熔溶剂与半纤维素之间形成的氢键数量越多，说明低共熔溶剂溶解半纤维能力越强；

根据主成分分析和偏最小二乘分析模型，得到氢键受体的碳链长度与氢键数量差值之间存在强烈正相关，氢键数量差值ΔH与相互作用能之间存在强烈负相关，说明碳链长度越长其氢键数量差值ΔH越大，从而快速地筛选出烷基链长的氢键受体，如苄基三乙基氯化铵作为氢键受体对半纤维素的溶解效果最好。

上述说明仅为本发明的优选实施例，并非是对本发明的限制，尽管参照前述实施例对本发明进行了详细的说明，对于本领域技术人员来说，其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改型等，均应包含在本发明的保护范围之内。

Claims (8)

1.一种基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，包括以下步骤：

(1)选择乳酸作为低共熔溶剂的氢键供体，季铵盐类作为低共熔溶剂的氢键受体，构建乳酸分子、每一种氢键受体阳离子、氢键受体阴离子的结构模型；并构建半纤维素分子的结构模型；

其中，所述低共熔溶剂的氢键受体包括盐酸胍、氯化胆碱、甜菜碱、乙胺氯、烯丙基三甲基氯化铵、四甲基氯化铵、四乙基氯化铵、四甲基溴化铵、苄基三甲基氯化铵、苄基三乙基氯化铵、苄基二甲基氯化铵中的任一种或几种；

(2)将步骤(1)构建的所有结构模型进行结构优化和能量优化，得到分子静电势信息，计算相应的RESP电荷；

(3)将步骤(2)优化后的结构模型通过CGenFF得到CHARMM36力场分子拓扑文件，并将计算得到的RESP电荷导入拓扑文件中；

(4)构建分子动力学初始结构盒子，将200～1000个氢键受体分子，500～2000个乳酸分子和20～100条半纤维素分子放入盒子中，其中，将半纤维素分子放在盒子中间位置，氢键受体分子和乳酸分子放在盒子四周位置，建立多个稳定的不同低共熔溶剂-半纤维素体系；

(5)运行不同低共熔溶剂-半纤维素体系，在298.15K和1bar压力下对低共熔溶剂-半纤维素体系采用最陡下降法模拟10000～15000步使系统能量达到最小，做限制性动力学，利用PME方法计算截止半径为1.0～1.8nm的静电相互作用，在363.15～413.15K处进行步长为25000000～50000000步、持续时间为50～100ns的长时间动力学模拟，得到每一帧下的模拟轨迹和结果；

(6)将步骤(5)模拟产生的每一帧轨迹结果导入可视化软件中，将得到的轨迹文件中的半纤维素的羟基进行分组，在分子动力学仿真计算程序Gromacs中运行氢键命令和相互作用能命令，计算每类羟基分别与每类低共熔溶剂的乳酸分子、氢键受体阳离子、氯离子之间的氢键数量和相互作用能，以及半纤维素分子之间的氢键数量，通过分析半纤维素与每类低共熔溶剂相互作用后的半纤维素分子之间的氢键数量H₂和半纤维素没有与低共熔溶剂相互作用的半纤维素分子之间的氢键数量H₁的差值ΔH，对低共熔溶剂进行筛选；

(7)将步骤(6)中计算得到的氢键数量差值ΔH和相互作用能导入到多变量分析软件SIMCA-P中，在多变量分析软件SIMCA-P中，对由不同官能团或不同碳链长度的氢键受体组成的低共熔溶剂与半纤维素之间的相互作用的变量，采用主成分分析和偏最小二乘法进行多变量分析；其中，对变量的主成分分析得到主成分的解释方差为0.911～0.998，预测方差为0.932～0.968；偏最小二乘分析得到的模型能够解释95.8％～99.8％的因变量变化，预测能力为80.1％～99.5％；根据主成分分析和偏最小二乘分析的模型得到氢键受体的官能团、碳链长度与氢键数量差值ΔH之间存在强烈正相关，氢键数量差值ΔH与相互作用能之间存在强烈负相关，说明碳链长度越长其氢键数量差值ΔH越大，从而快速地筛选出对半纤维素的溶解效果最好的低共熔溶剂；最后结合偏最小二乘法得到的变量重要性分析，得出低共熔溶剂与半纤维素之间形成的氢键数量越多，说明低共熔溶剂溶解半纤维能力越强。

2.根据权利要求1所述的基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，所述步骤(1)中，低共熔溶剂的氢键供体与氢键受体的摩尔比为HBA：HBD＝1：(5～9)。

3.根据权利要求1所述的基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，所述步骤(1)中半纤维素的结构模型包含6～8个五碳糖组成的聚合物结构、或者包含2～7个六碳糖组成的聚合物结构、或者包含6～8个五碳糖和2～7个六碳糖组成的聚合物结构。

4.根据权利要求1所述的基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，所述步骤(2)中，将步骤(1)构建的所有结构模型保存为gjf文件，在量子化学仿真软件Gaussian中基于密度泛函理论进行结构优化和能量优化，得到各自的分子静电势信息，结合量子化学波函数分析程序Multiwfn计算相应的RESP电荷。

5.根据权利要求1所述的基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，所述步骤(4)中，采用分子动力学初始结构构建程序Packmol构建10nm×10nm×16nm的分子动力学初始结构盒子。

6.根据权利要求1所述的基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，所述步骤(5)中，在分子动力学仿真计算程序Gromacs中运行步骤(4)建立的低共熔溶剂-半纤维素体系；当能量达到最小后逐渐加热到363.15～413.15K，在达到所需温度后，能量最小的低共熔溶剂-半纤维素体系在NVT体系合成中对半纤维素产生限制性动力学。

7.根据权利要求1所述的基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，所述步骤(6)中采用的分子可视化软件为VMD软件。

8.根据权利要求1所述的基于分子动力学和多变量分析的低共熔溶剂筛选方法，其特征在于，所述步骤(7)中的变量包括不同氢键受体的物理化学性质和分子模拟结果；

其中，不同氢键受体的物理化学性质包括：熔点、分子量、分子极性指数、氨基数量、羟基数量、碳碳双键数量、碳链长度、氢键供体数量、氢键受体数量、氢键供体与氢键受体数量之和、可旋转化学键数量、重原子数量、分子极性表面积和复杂度；分子模拟结果包括：氢键数量差值、氢键寿命、相互作用能、径向分布函数和空间分布函数。