分子动力学模拟 心得 适合新手！！！

本文来源：“分子动力学”公众号

原文链接：https://mp.weixin.qq.com/s/fpiUACqfYyNq5y5ExuKiqQ

看完《分子模拟从算法到应用》那本书的第四章，不用全看完，但是至少要对分子动力学模拟过程有一个了解。试着按照书的过程做个Ar的NVE，其实Ar和离子晶体以及其它的任何材料的差别仅仅是势函数的问题，虽然由势函数带来了一些问题，但是这些都不是本质问题。从初始化的原子数，原子位置，初始速度，时间步长，初始温度等等这些初始化结束了以后，选择一个简单的积分算法，如6阶的Gear预测校正，不要控温控压，就是一个简单的NVE，不要考虑任何的提高效率的邻位算法，因为这个时候我们可以选择5×5×5的超原胞，总共的原子数也就500个，不需要考虑邻位算法。开始循环计算：预测----计算原子的力和能量—校正—输出能量。

这样最简单的NVE就编成了，总共也就1000多行，是个很小的程序。自己先试着感受一下。其实当这样的小程序完成以后，你会觉得分子动力学编程也很简单，那么接下来的复杂的分子动力学也不会是什么问题。

做完了这一些，你需要知道的是那些是和材料无关的东西，那么就尽量的分离，开始使用一个个的函数。例如，原子的位置是和具体的材料相关，但是初始速度却和材料一点关系都没有，同样的数值积分中的预测和校正也是和材料无关的，以后的控温和控压算法也是和材料无关的。当规模大了以后，邻位算法也是和材料无关的，像这些和材料无关的部分最好自己做成小函数，选择调用。以后换材料的时候程序也不会有太大的改动。

编程小技巧

1.选择用intel编译器，个人喜欢用10.1或者9.1的版本，打开优化选择，类似的/QaxS /QxS /Qipo /Qprec-div-等等，以后可以使用Openmp的并行计算（具体的可以参考intel编译器的帮助手册）

2.尽量的简化计算，例如2×a就要写成a+a，在计算机中，加减是一个数量的计算，乘是一个，除是一个。计算量逐渐增加，所以选择a+a来替代2×a 会减少一些计算量

3.选择数据来替代结构体，结构体看的比较方便，但是计算效率要低，尤其是在编译过程中的矢量化的时候，而数组则可以很好的矢量化，也更适合并行。结构体的计算是先找到结构体指针，然后再找里面的参数，当并行的时候，多线程同时找结构体指针，会很大的减低速度。

4.除法尽量用乘法来替代，有了程序以后，自己可以仔细的分析里面的计算消耗，这样可以更好的优化计算，当然最主要的消耗是在力和能量的计算中，可以选择离散-插值的方法，使得不同的势函数有相同的计算效率，也可以选择邻位算法来分块计算。

NVT程序

当简单的NVE做完以后，可以直接做NVT，这仅仅是增加一个简单的控温，是很简单的，早期编程的时候，建议选择Nose-Hoover控温，你会发现 NVT也是如此的简单。接下来，你可以考虑做做NPH，但是由于引入了压强算法，原来的原子的位置和速度等等和笛卡尔系有关的一切内容都要发生改变，这将让你重新写程序，可以说是一个完全不同的程序。不过，幸好我们已经有了前面的一些经验，NPH虽然复杂，但是并不是不可能完成的，注意的是和笛卡尔系相关的量要变化。早期的计算，可以选择Anderson控压算法。在完成了这些，可以说你已经把最主要的分子动力学程序都完成了。

分子动力学算法

那么，接下来就是复杂的分子动力学算法了，试着选择PR算法来代替Anderson压强算法，然后选择Metric-tensor来代替PR算法；选择 Nose-Poincare来代替Nose-Hoover算法，选择Generalized-Leap-Frog算法来替代预测校正算法，用Wolf来替代Ewalds算法，这些一步步的改进都会让你有很多新的发现。你会觉得你的程序达到了现在大部分软件包没有的功能。随着规模的增大，你可以选择 Verlet列表，原胞列表、结合法以及快速排序来实现邻位算法，如果再有了openmp并行计算的加入，你会发现你的程序已经有了一个质的突破。

END

当这些都完成了以后，接下来就可以改变势函数来模拟不同的材料了，而且你对分子动力学的理解也会上升到一个新的高度，此外，其余相关的结果也不会是问题了，例如，弹性模量的计算，MSD，径向分布函数等等。你可以直接的加入代码来实现你想要的过程。