lammps計算的應力有兩種：

一是體系整體的應力狀態(tài)，通過在thermo_style custom里加上pxx pyy pzz pxy pxz pyz字段可將給定時間步（由thermo N命令所指定）的體系應力值輸出，再求時間平均即可（實際上求出的是壓強張量，即負的應力值）；該應力的計算用的是統(tǒng)計力學里的Virial定理（參見<<Computer simulation of liquid>> by Allen & Tildesley），所算出來的應力與宏觀應力是一致的（強調(diào)一下，用于平衡態(tài)）；

compute 1 all pressure thermo_temp virial

二是單個原子的應力，也就是樓主所討論的，通過compute stress/atom命令所計算出來的六個值；這個應力的計算則是通過對上述Virial定理所定義的體系應力按原子分解，即將公式中的按原子求和的算符拿掉（同時應將體系的體積換為單個原子的體積），不過，正如樓主所指出的，由于單個原子的體積計算太麻煩，lammps在計算時干脆去掉了體積項，這就是為什么用compute stress/atom命令所算出來的“應力”具有能量的單位的原因。

可以看出，在lammps里，如果要計算體系中某個區(qū)域（由region定義，可以是整個模擬盒）所圍成的“塊”的應力，只需將該區(qū)域里的所有原子的單原子應力值加起來，再除以這個區(qū)域的體積即可，無須進行單個原子體積的計算。

對原子應力進行Voronoi體積加權平均即可得到系統(tǒng)瞬時應力；系統(tǒng)瞬時應力的系綜平均值為宏觀測量的系統(tǒng)應力值。

compute s all stress/atom
compute p all reduce sum c_s[1] c_s[2] c_s[3]
variable Press equal -(c_p[1]+c_p[2]+c_p[3])/(3*vol)
thermo_style custom step temp etotal press v_Press

如果要求X方向的應力,就把所有原子用compute stress/atom得到的sigma(xx)加起來,去除掉"所有原子佔有的體積",就可以得到系統(tǒng)宏觀的sigma(xx).得到的值與別人發(fā)表的計算paper相差無幾(與試驗值相差大約10%)。

文章轉(zhuǎn)載自：http://blog.sina.com.cn/s/blog_c0468c8f0101tfea.html#opennewwindow

本文轉(zhuǎn)載于網(wǎng)絡，版權歸原作者，若有不妥，請聯(lián)系刪除！