1.看完《分子模擬從算法到應用》那本書的第四章，不用全看完，但是至少要對分子動力學模擬過程有一個了解。

2.試著按照書的過程做個Ar的NVE，其實Ar和離子晶體以及其它的任何材料的差別僅僅是勢函數的問題，雖然由勢函數帶來了一些問題，但是這些都不是本質問題。
3.從初始化的原子數，原子位置，初始速度，時間步長，初始溫度等等這些初始化結束了以后，選擇一個簡單的積分算法，如6階的Gear預測校正，不要控溫控壓，就是一個簡單的NVE，不要考慮任何的提高效率的鄰位算法，因為這個時候我們可以選擇5×5×5的超原胞，總共的原子數也就500個，不需要考慮鄰位算法。
4.開始循環計算：預測----計算原子的力和能量---校正
5.輸出能量

這樣最簡單的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.除法盡量用乘法來替代，有了程序以后，自己可以仔細的分析里面的計算消耗，這樣可以更好的優化計算，當然最主要的消耗是在力和能量的計算中，可以選擇離散-插值的方法，使得不同的勢函數有相同的計算效率，也可以選擇鄰位算法來分塊計算。

當簡單的NVE做完以后，可以直接做NVT，這僅僅是增加一個簡單的控溫，是很簡單的，早期編程的時候，建議選擇Nose-Hoover控溫，你會發現 NVT也是如此的簡單。接下來，你可以考慮做做NPH，但是由于引入了壓強算法，原來的原子的位置和速度等等和笛卡爾系有關的一切內容都要發生改變，這將讓你重新寫程序，可以說是一個完全不同的程序。不過，幸好我們已經有了前面的一些經驗，NPH雖然復雜，但是并不是不可能完成的，注意的是和笛卡爾系相關的量要變化。早期的計算，可以選擇Anderson控壓算法。在完成了這些，可以說你已經把最主要的分子動力學程序都完成了。

那么，接下來就是復雜的分子動力學算法了，試著選擇PR算法來代替Anderson壓強算法，然后選擇Metric-tensor來代替PR算法；選擇 Nose-Poincare來代替Nose-Hoover算法，選擇Generalized-Leap-Frog算法來替代預測校正算法，用Wolf來替代Ewalds算法，這些一步步的改進都會讓你有很多新的發現。你會覺得你的程序達到了現在大部分軟件包沒有的功能。隨著規模的增大，你可以選擇 Verlet列表，原胞列表、結合法以及快速排序來實現鄰位算法，如果再有了openmp并行計算的加入，你會發現你的程序已經有了一個質的突破。

當這些都完成了以后，接下來就可以改變勢函數來模擬不同的材料了，而且你對分子動力學的理解也會上升到一個新的高度，此外，其余相關的結果也不會是問題了，例如，彈性模量的計算，MSD，徑向分布函數等等。你可以直接的加入代碼來實現你想要的過程。

很多好的文章請關注微信公眾號： 互動派教育

http://flac3d.cn/hdp/grom/zwt.html

http://www.hdpaii.com/