weixin_44873868的博客

本案例将通过模拟退火方式获得无定形二氧化硅结构，即令晶态二氧化硅升温至融化以破坏晶态结构，如何降温成固态，最后做结构优化。注意不能对熔融状态下的二氧化硅直接做结构优化，否则会收敛到能量很高，无化学意义的极小点结构。首先进行高温动力学模拟，CP2K输入文件中任务类型设为分子动力学，理论方法选择GFN1-XTB，温度设为3000，采用CSVR热浴，计算10 ps。接着进行降温动力学的模拟，CP2K输入文件温度设为800，采用CSVR热浴，模拟时间计算6 ps。图5 模拟产生的非晶二氧化硅的径向分布函数。

GROMACS（GROningen MAchine for Chemical Simulations）是当前广泛应用的分子动力学模拟软件，具有高效的计算能力和灵活的配置选项，支持多种不同类型的分子动力学模拟，并广泛应用于纳米物理学、化学反应动力学等领域。如图4所示，随着初始水球的逐渐蒸发，最大的簇中所含的分子数逐渐减小。对常温平衡后的纳米水球进行450 K高温下的蒸发模拟（图三），可以看到，随着模拟的进行，大的水球逐渐被分裂为一个个小的水团簇，说明水的蒸发过程正在进行。四、纳米水滴蒸发的模拟应用。

例如，在自由基反应中，自由基的生成通常依赖于键解离能的大小，键解离能较小的键更容易断裂，从而生成自由基。在苯甲酸酯（如methyl benzoate、phenyl benzoate）中，产生的自由基（如芳香酰自由基）同样可以通过苯环的π电子进行一定的稳定化，但相较于苯甲酸，酯的碳-氧键解离后的自由基可能稳定性较差，因此BDE稍低。溶剂的多样性：有机分子可能在多种溶剂中表现出不同的反应性，因此，对于不同溶剂环境下的计算，需要采取不同的溶剂模型和方法，而这些模型的选择和调优本身也需要大量的经验和试探。

（3）非平衡分子动力学方法。由于流体具有黏性，粒子或分子的运动不会立即跟上外力的变化，而是会有一定的滞后。(3)计算黏度：通过测量流体系统在周期扰动下的应力响应和外力的关系，结合相关的理论公式来计算黏度。在周期扰动法中，通常会测量系统在多个周期内的响应，并对数据进行统计分析，以获得准确的黏度值。(1)施加周期扰动：在流体系统中引入一个周期性变化的外力场，比如在模拟体系的某个方向上按照一定的频率和振幅施加外力。周期扰动法计算黏度的原理主要基于对流体系统施加周期性的外力扰动，然后通过分析系统的响应来计算黏度。

然后点击Task右边的More..Quality选择Fine，然后切换到option选项卡勾选Use delocalized internals。具体做法：首先双击打开SiC_beta，然后点击Modules | CASTEP | Calculation 选择几何优化任务。然后回到CASTEP calculation对话框的Properties选项卡，将System type勾选为Crystal。关键词：Materials Studio，DFT，castep，介电常数。优化的时候记得把晶胞设置成P1。

差分电荷密度（Differential Charge Density）、d带中心理论（d-Band Center）与晶体轨道哈密顿布居（Crystal Orbital Hamiltonian Population, COHP）作为量子尺度的重要分析工具，可直观揭示化学键形成机制、电荷转移路径及催化活性位点特性，为材料理性设计提供理论基石。1. 差分电荷密度：通过计算吸附/解吸附过程的电荷的转移，可以直观看出材料分界面处的电荷的得失（如图2、图3、图4所示）。图3 金属离子在MXene上的扩散路径。

我们进一步统计正十六烷的溶剂可及表面积（SASA）的变化（图5），发现正十六烷从一开始的混乱无序分布自发组装成膜，SASA逐渐减小，到2ns时油分子膜已趋于平衡，说明该体系已自发完成了油水相分离过程。经过能量最小化和2ns的平衡模拟后，我们可以观察到，在无外界干预的条件下（常温常压），油水分子的相互作用（疏水作用）驱动体系自发发生相分离，逐渐形成油滴与水相的分层结构。本案例基于GROMACS，建立油水混合体系的分子动力学模型，模拟其在常温下的自发相分离过程，并考察油分子的溶剂可及表面积的变化。

300 K 时的循环应力-应变曲线如图2所示，当高温合金受到循环加载时，最大应力随循环次数的增加而增大，即首先发生应力循环硬化，这主要是由于初始缺陷的积累，如位错、堆垛层错等。图3和图4分别为不同应变下Cu单晶的循环载荷、剪切应变、位错分析、公共邻域分析的可视化图，通过ovito可视化后，可以发现循环载荷下Cu单晶存在明显的应力集中现象，同时发生均匀相变，在Cu单晶内部可以发现存在少量的bcc以及Other原子，这对Cu单晶的变形和力学性能有显著的影响。(c)循环载荷下位错分析；(c)循环载荷下位错分析；

利用 amorphous cell 模块分别构建包含缓蚀剂分子的无定形组织结构，采用(NVT)进行分子动力学模拟，平衡后统计体系密度的平均值，并把该值作为计算体系中缓蚀剂膜的密度。腐蚀介质粒子在缓蚀剂膜中的扩散行为的模拟通过forcite模块的正则系综(NVT)来实现,模拟温度为 298 K，温度采用 Andersen方法控制，各分子起始速度由Maxwell-Boltzmann分布随机产生，运用 velocityverlet 算法叫求解牛顿运动方程.在众多的防腐蚀措施之中，添加缓蚀剂是广泛采用的方法之一。

pKa 值反映了分子在水溶液中的酸性或碱性强度，通常通过实验测定，但实验方法常常受到溶剂效应、温度、离子强度等因素的影响，且对于复杂分子的测定具有较大的挑战。胺类分子的 pKa 值受到多种因素的影响，包括分子的结构（如取代基、环结构、电子效应等）、溶剂的极性、温度和分子的构象变化等。酸碱反应的平衡 计算 pKa 值要求对酸碱平衡进行全面分析，通常涉及到化学反应的自由能变化，这不仅仅是计算单一反应物和产物的能量差异，还要考虑过渡态的能量和解离过程中的所有中间态。对这些现象的精确建模仍然具有挑战性。

在x方向再进行同样的操作，固定干酪根z方向壁面，在体系的x方向左右两侧设置LJ势能壁，固定左侧壁面，使右壁面向左缓慢移动，对干酪根进行横向压缩。因此，本文采用分子动力学模拟方法，研究体相及纳米孔隙中CO2/原油的混相机理，以及纳米孔隙中CO2驱替原油的作用规律：另外，根据现场应用时CO2的注采工艺特点，分别开展了CO2加压驱替以及降压抽提原油时的混相行为研究：最后，根据气驱的基本特点，提出并分析了CO2/N2段塞驱提高致密油采收率的方法与微观机理。而在干酪根孔隙内，仅有驱替速度为2 m/s时满足条件。

这就导致计算思路产生难以逾越的困境，针对于次，我们给出了在Forcite的模拟中施加电场的解决方案，通过简单的参数设置，用户只需输入电场强度、方向等基本信息，脚本即可自动完成电场设置，并能够进行运算。电场与力场的结合：在科研计算中，电场和力场往往需要结合使用。类比于电场，在具体的研究中，我们还需要研究下压作用力或者上拉作用力对体系的影响，有特别是在做钻头方面研究的同学，需要研究不同下压力下的作用情况。在进行电池等领域的模拟的时候，我们往往需要体系处于电场当中，然而，ms的设置面板没有进行电场设置的选择。

首先，通过DFT计算可以优化稀土化合物的几何结构，获得最低能量构型，从而帮助确定金属中心的配位数和几何构型，分析化合物的稳定性。其次，Gaussian能够精确计算稀土化合物的电子结构，揭示电子轨道分布和f电子的贡献，为理解其光学性质（如吸收光谱、发光特性）和磁性（如磁矩、磁各向异性）提供重要信息。Gaussian采用DFT方法时，常常需要选择合适的交换-关联泛函（如LDA、GGA或混合泛函），但是传统的DFT方法可能无法精确描述f电子的局部化特性或强电子相关效应，这可能导致优化结果与实验观察之间存在差异。

RESP电荷方法，即Kollman等人在其开创性论文《A well-behaved electrostatic potential based method using charge restraints for deriving atomic charges: the RESP model》中提出的技术，是当前最适合用于柔性小分子模拟的电荷模型。通过本教程，您将掌握从Gaussian优化结构到AmberTools计算RESP电荷及力场参数的整个流程，为您的分子模拟提供坚实的基础。

无论您是研究新分子的性质，还是优化现有分子的结构，Gaussian16都能够为您的工作提供精确的计算结果和深刻的洞见。跟随我们的教程，您将能够顺利安装Gaussian16，并熟练掌握其结构优化功能，为您的科研工作奠定坚实的基础。使用tar -jxvf G16-A03-AVX2.tbz命令解压Gaussian16的压缩包，在本教程中解压路径为/root/autodl-tmp/student1，涉及路径的命令请根据实际情况修改。本教程以JZ4为例，从Pubchem上获取JZ4的3D结构SDF文件。

但是在生成力场和电荷的时候，根据分子特征还要进行许多调整得到适合这个分子的具体参数，例如在计算电荷时使用什么基组和泛函进行计算等等，这些选择的合适与否会直接影响到计算精度和结果。)可以产生opls力场的作用参数。此时会生成一个amb2gmx的文件夹，文件夹中包含了模拟需要的所有文件例如分子结构，力场参数和mdp文件，该参数基于GAFF力场，准确性得到了许多相关工作的验证。然后通过使用gaussian计算拟合得到分子的RESP电荷，最终将原子电荷写入之前生成的力场文件中得到最终可用来执行模拟的力场文件。

相比之下，图(f)中的裂纹在相同应变下的扩展宽度d2较小，这表明裂纹的扩展相对较慢，或者裂纹的扩展受到了某种程度的抑制。由于纳米晶体金属及合金材料具有优越的物理、化学、力学特性，越来越受到人们的重 视，但是材料的缺陷严重影响着人们的安全，所以研究裂纹的扩展机制成为一项重要的研究课题。图(g)和图(h)分别为图(a)和(b)在拉伸载荷作用下的原子应力集中示意图，值得注意的是，尽管预裂纹在图(a)和图(b)中的建立方式有所不同，但两者在拉伸载荷下的原子应力集中现象都呈现出相似的模式。

该脚本具备高度的可定制性，您可以根据具体的研究需求和计算条件，灵活设置脚本参数，实现个性化的计算任务配置。然而，我们的计算往往是相同的参数进行大量的计算，面多几十甚至几百个结构的计算量时，这种手动寻找能量的方式就难以使用。，因此，使用我们所提供的脚本，可以批量的提交计算任务，获得计算结果。特别是需要获取结果的能量等信息的时候，要进入每一个结果文件进行查找，当你面对少量的模型时，能量提取还能hold住。在一般的计算中，我们要实现对结构的计算需要依次打开每一个结构文件，进行参数设置和计算，并获得计算结果文件。

利用同样的方法我们可以查看每个构型在不同入射频率下的的比旋光度，下表列出了所有计算和实验得到的比旋光度，可以看出在这个例子中，S构型计算得到的比旋光度与实验测得的更为接近，因此我们认为实验测量的分子为S构型，这也与事实一致。使用第一步优化得到的结构作为初始结构，计算S和R型甲基环氧乙烷的比旋光度。计算S和R型甲基环氧乙烷在气相条件下，355nm和633nm入射光下比旋光度的输入文件分别如下。以甲基环氧乙烷作为例子，通过计算两种手性分子构型的比旋光度，并与实验值对比，以确定分子的绝对构型。

3、大量分子体系的内聚能密度和溶解性参数对于估算自扩散系数的均方位移和速度相关函数研究润滑剂及润滑能力的限制性剪切模拟在学习表中观察并绘制轨迹数据，按任意性质排序，如按能量排序，找到最低能量构型。DMoL3 可以模拟气相、溶液、表面及固体等过程及性质，应用于化学、材料、化工、固体物理等许多领域，研究分子结构、分子反应、均相催化、多相催化等。CASTEP：先进的量子力学程序，广泛应用于陶瓷、半导体、金属等多种材料的研究，包括晶体材料的性质、表面和表面重构的性质等。

随后，我们将详细阐述模拟的具体步骤，包括结构准备（如蛋白质结构预测、小分子结构优化等）、相互作用模拟（如分子对接、分子动力学模拟等）以及结果分析（如相互作用能计算、轨迹分析等）。只有在力场的.rtp文件中存在构建块的条目时，拓扑才能自动组装。通过本教程的学习，您将能够掌握从蛋白质与小分子配体的结构准备、相互作用模拟到结果分析的全流程，从而能够自主进行相关的模拟研究。本教程使用的力场为amber14sb.ff，选择默认的水模型TIP3P，然后为封端选择“NH3+”和“COO-”，获得力场文件及完整坐标文件。

总之，我们的模拟拉伸-断裂过程的Perl脚本为材料科学研究提供了强有力的工具。它不仅能够提高研究效率，节省实验成本，还能够揭示材料性能背后的微观机制，为新材料的设计与开发提供有力支持。这款脚本能够高效模拟材料的拉伸-断裂过程，帮助科学家们在计算机中重现实验场景，从而更深入地理解材料在受力过程中的行为。该脚本通过强大的算法和精确的物理模型，能够模拟材料在不同条件下的拉伸过程，包括应力-应变关系、微观裂纹的萌生与扩展等关键参数。这一过程涉及复杂的力学行为和材料内部微观结构的变化，对于理解材料的性能至关重要。

由于NICS值为各向同性化学屏蔽值的负值，因此苯的NICS(0)和NICS(1)分别为-7.5538和-10.5301，这也表明苯分子具有芳香性，当比较不同的分子时，NICS值越负，表明具有越强的芳香性。其中第1行的Bq即为NICS(0)虚原子对应的位置，7.5538为NICS(0)对应虚原子位置处的各向同性化学位移屏蔽值。第14行为NICS(1)设置的虚原子，通过复制第13行NICS(0)虚原子而来，仅是修改了坐标z的值。设置NICS(0)和NICS(1)的计算位点。计算分子NICS值的基本流程。

Panagiotopoulos 等 人提出的 Gibbs 系综 Monte Carlo 方法（Gibbs ensemble Monte Carlo，GEMC）是近十几年发展起来的应用最为广泛的计算流体相平衡的方法，也是目前模拟相平衡的主流方法之一。GEMC 方法可同时在两个彼此相对独立但在热力学上相关的盒子中进行 MC 模拟，模拟时需要满足相平衡条件（压力、温度和化学势相等），且模拟过程中温度 T、总体积 V 和两个盒子中的总粒子数 N 保持不变。2）构建异丁烷的PDB文件，描述异丁烷残基的拓扑文件；

图8 优化后模型及MSD参数（a：PPI）（b：PPI-PPTA（PPI-2））图2 两个模型的AC参数（a：PPI）（b：PPI-PPTA）（c：模型参数）图10 计算平衡后的能量和温度数据（a：PPI）（b：PPI-PPTA）图6 计算平衡后的能量和温度数据（a：PPI）（b：PPI-PPTA）图8 计算平衡后的能量和密度数据（a：PPI）（b：PPI-PPTA）图9 MSD结果（a：PPI）（b：PPI-PPTA（PPI-2））图4 模型的收敛以及能量数据（a：PPI）（b：PPI-PPTA）

而在具体的计算中，不需要将配体分子和受体分子之间的相互作用按照热力学循环每一个步骤进行单独计算，只需要慢慢的打开和关闭两者之间的相互作用，通过转换，就可以得到最后的自由能。Lamda为0时是该体系的第一个模拟，在该状态下，相互作用是全部打开的，而随后慢慢打开Coul作用，之后继续打开vdw作用，这相当于上述循环右侧的耦合步骤，逐步增加配体-环境或者受体-配体之间的作用。如下图所示，直接计算受体-配体之间的结合自由能是比较困难的，但是通过炼金术方法，可以将自由能的计算转化在热力学循环中。

盐湖中高浓度的镁，钠，钾等离子是萃取锂资源的主要影响因素。文献通过Materials Studio软件模拟计算，采用分子动力学中的径向分布函数（RDF）研究TBP对不同金属离子(Li+ , Na+ , K+ , Mg2+)的内在相互作用强度。利用Materials Studio (MS)中的Amorphous cell模块构建混合有机相模型，模型包含10个LiBPh4分子、10个KBPh4分子、10个Mg(BPh4)2分子、40个TBP分子、40个CH2ClBr分子和60个H2O分子。结构模型如图2所示。

点击选择Forcite，Calculation，弹出Forcite Calculation ，选择Energy 中Forcitefield中的COMPASSII，点击旁边的More，取消默认的两个勾选，单击面板中的两个Calculate，确保模型中所有的原子具有正确的力场，可以点击一个原子，在信息栏里查看力场是否正确建立。选择xtd轨迹文件，点击Modules，Forcite，Analysis，选择Radial distribution function，Sets中分别选择O和H，计算O和H之间的RDF。

使用 vaspkit-681命令生成高密度的KPOINTS，然后进行静态计算 （注意只有使用这项功能生成KPOINTS计算的结果才能继续使用vaspkit命令计算下一步，使用M-S方法自动生成K点的计算结果无法进行下一步）。这一步一般需要非常长的时间去进行计算，因此扩胞的大小可以稍微比二阶力常数矩阵的计算时的大小小一点。其中Phonopy用于计算声子谱及二阶力常数，Thirdorder用于计算三阶力常数，ShengBTE用于结合前面两者的结果计算晶格热导率。准备好计算的材料对应的POSCAR。

例如两个蛋白之间、蛋白和纳米材料之间、蛋白和配体之间的相互作用能的准确计算能够帮助研究人员更深层次的理解相关科学问题。这两个能够给出两个物体之间的整体的相互作用能、相互作用能的组成部分以及物体内每一部分例如每个氨基酸对相互作用能的贡献数值，从而在不同层面解析物体间的相互作用。此项工作的计算流程并不难，但是难点在于理解计算脚本中的每一项内容，从而根据自己的需要去得到合适准确的计算脚本，以及理解MMPBSA计算的物理含义。根据之前的定义，更改计算脚本中相应的文件和结构名字，并根据文件内容更改相应软件路径。

最后计算糊精抑制剂在煤样、石英、高岭石表面相互作用的分子动力学计算，采用Materials Studio 2019软件中Forcite模块对糊精与三种矿物相互作用模型进行搭建并计算。同时，将搭建好的糊精模型用Dmol3进行优化计算，计算参数选择GGA广义梯度近似和PBE泛函方法，在计算精度为Fine的基础上选择DNP 4.4基组，体系总能量收敛值取 1.0×10-5Ha，最大内应力为0.002 Ha/Å，最大位移为0.005Å，选用全电子核处理方式，使用TS色散修正方法进行修正。图3煤样的XRD分析。

由上图可以得出ChCl-2MA低共熔溶剂与氧化锌相互作用共分为三部分， ChCl-2MA中的Cl和C=O中的O与氧化锌中的Zn之间的距离在3~3.5内并且介于两种原子的共价键半径之和之间，表明Cl、O与Zn以化学作用的形式相互结合。H次之，O-Zn最弱。Forcite模块是分子动力学计算的主要模块，研究范围广，可以对多种周期性体系进行计算分析，在矿物分选领域中主要是计算分析药剂与矿物相互作用，在不同计算参数条件下可以实现药剂与矿物相互作用模型的预测与分析，从而得到表面相互作用机理。

由此，模拟技术与实验、理论三者的结合是现在以及未来被认可和推广的研究手段。. 在分子动力学的模拟研究中，一款开源、自由、免费的软件GROMACS得到了广泛的应用。. 它可以用于几百万个粒子体系的分子动力学模拟研究，特别是生物体系，比如磷脂双分子层生物膜、蛋白质、药物分子等。针对gromacs软件的常用模块进行教学，包括蛋白与配体模拟分析，离子液体，小分子与细胞膜相互作用，同时介绍gromacs中的各分析模块使用功能。2.2气相分子，分子团簇与纳米液滴的模拟。1.7.分子两端的距离和原了距离的分析。

本次模拟采用LJ约化单位，初始晶体模型为六方最密堆积结构，晶格参数为0.7，沿x（100）方向为20倍晶格长度，y方向（010）为20倍晶格长度。在进行流动模拟时，为顶端的固定区域设置沿x方向的初速度为5.0，其他方向速度为0。这样的滞后体现在空间和时间尺度上。这里也可以清楚的看的靠近顶端移动平板的原子在相同时间内有着更长的移动距离。这个流动是由作用在流体上的粘性力和与平板平行的外部压力推动的。图2.2展示了在初始1000步，50000步和100000步时流体原子沿着y方向的x方向速度（vx）的分布情况。

CO2能显著降低原油的黏度和界面张力，并可通过溶胀作用进一步提高原油采收率，因此，被认为是提高低渗透油藏原油采收率最有效的方法之一。因此,本文通过MS软件进行分子动力学模拟，分析沥青质沉积对CO2驱油的效果。使用力场为文献中常用的COMPASS力场，对体系分别进行200ps的NVT和200PS的NPT（可进行不同温度以及压力下的模拟）分子动力学模拟。对最终输出的体系进行500PS的NVE分子动力学体系。通过MS的建模功能，构建出，树脂，沥青质和烷烃等模型，沥青质-Fe模型，以及在体系内填充CO2后的模型。

在此我们还需要准备需要的力场文件，然后利用我们提供的脚本来进行力场的建立。如果此时使用gmx pdb2gmx -f AA.pdb -o AA.gro命令来建立力场时，会发现gromacs无法识别小分子并产生相关的力场。其中要注意的是生成力场是一定要对配体把缺失的H原子补上，以生成准确的全原子力场。本文介绍一种建立力场的方式以及在蛋白质-小分子模拟过程中需要注意的重要问题。这一较为精确的力场产生工具。但是在使用该工具是仍会遇到各种问题，例如脚本内容的更换，原子的检查以及成键信息或者位置限定的多种方式。

本次模拟采用LJ约化单位，晶体为面心立方结构，晶格参数为0.6，沿x（100）、y（010）方向为10倍晶格长度，z（001）方向为20倍晶格长度。采用compute chunk/atom将模型沿着z方向分成20块，设置底端0-1倍晶格长度为热端，中间10-11为冷端。模拟温度为1.35。热端设置为1.7，冷端设置为1.0，弛豫到热源和冷源的温差达到相对稳定的状态。在最后一1000步的温度梯度和温差波动如图1.2所示，可以明显看出热端和冷端的温度差，热流保持着相对稳定的状体。1.3结果整理与分析。

通过图3可以看出氯离子与氧、碳原子发生强烈共振,最高占据态内是由Cl-3s、Cl-3p、O-2s、O-2p、C2-2s、C2-2p、C4-2s、C4-2p、C5-2s、C5-2p贡献，其中Cl-3p、O-2p、C2-2p、C4-2p、C5-2p占主导地位，说明形成的p-p杂化轨道是阴阳离子间电子转移的主要通道。2.将氯化胆碱与丙二酸放入同一面板中应用Dmol3进行相互作用计算得到氯化胆碱-丙二酸低共熔溶剂分子构型（ChCl-MA）（如图4），得到稳定构型后可以进行径向分布函数性质分析（如图5）。

自由能形貌图一般通过两个描述体系特征的量来进行绘制，例如RMSD和Rg，也有文献中用主成分分析PC1和PC2绘制。对rmsd.xvg和rg.xvg文件进行处理，删除所有注释行（以#或者@开头的行和空行）以及rg.xvg中的XYZ三个方向上的Rg数据。利用主成分分析绘制自由能形貌图，一般来说就是用主成分1和2取代上文中的RMSD和Rg，后处理都是一样的。2. 合并rmsd.xvg和rg.xvg文件，把同一时间下的RMSD和Rg数据写在同一行，组合文件中包含三列：时间、RMSD、Rg，可使用如下命令实现。

我们可以通过试验数据比较，进行随机掺杂、指定位置掺杂等。可以计算掺杂原子的形成能垒，掺杂后对周围原子的势能和结构分布畸变、力学性能、热导率等，研究掺杂原子对结构的改变而引起的力学、热学和辐照性能的影响。同时通过TEM、XRD等实验手段观测到的团簇或者具有特定的分布函数，也可以通过编程现实。2.3在完成掺杂后可以采用displace_atoms等方式计算掺杂原子的移动能垒、也可以进一步计算其辐照、力学、热学性能。2.2. 通过自己的需要，通过编程或者lammps中的create_atom进行掺杂。