|||
原帖在http://www.52souji.net/how-to-calculate-equilibrium-lattice-constant-using-lammps.html
正文中我加了少许改动,还是感谢原作者的贡献,红色字体是我添加的。
平衡晶格常数(equilibrium lattice constant)对应的体系能量是最低的,因此只需要计算一些列不同晶格常数下体系的能量,那么体系能量最小时对应的晶格常数就是平衡晶格常数。
在初次计算时,建议将晶格常数增加的步长取得大一点,比如0.1A,这样通过结果数据就可以大致确定晶格常数的范围。
然后再将晶格常数的范围缩小,减小步长,比如0.01A,这样做就可以获得精度比较高的结果。下面以金刚石为例介绍具体如何计算平衡晶格常数。
输入文件: in.diamond
# This input script is used to calculate
# the lattice constant of diamond
# Powered by Xianbao Duan
# Email: xianbao.d@gmail.com
# Website: http://www.52souji.net/
units metal#单位是金属
boundary p p p#周期性边界条件
atom_style atomic
variable i loop 20#20次循环
variable x equal 3.5+0.01*$i#晶格常数
# build the model
lattice diamond $x
#lattice diamond 3.5
region box block 0 10 0 10 0 10
create_box 1 box
create_atoms 1 box
# specify the potential
pair_style tersoff
pair_coeff * * SiC.tersoff C #C 是使用原子
mass 1 12
variable n equal count(all)#原子的总数
variable P equal pe/$n
timestep 0.005
thermo 10
# minimize the total energy
min_style cg#算法
minimize 1.0e-12 1.0e-12 1000 1000#最小化条件
print "@ $x $P"
clear
next i
jump in.diamond
脚本很简单,只要一个简单的能量最小化就可以了,然后输出晶格常数和结合能log文件中。
然后从log文件中将晶格常数和结合能的数据提取出来。
grep ^@ log.lammps > lat.vs.Ecoh.step1
grep是一个linux命令(搜索命令),上面的意思就是将log.lammps文件中以@开始的行输出到 lat.vs.Ecoh.step1 文件中。
得到的文件中,每一行开头会有一个@,通过列编辑器一起删掉,就可以获得有如下内容的文件。
3.51 -7.35079243 3.52 -7.356562846 3.53 -7.361141466 3.54 -7.364556691 3.55 -7.366836396 3.56 -7.368007936 3.57 -7.368098161 3.58 -7.367133418 3.59 -7.365139565 3.6 -7.362141974 3.61 -7.358165545 3.62 -7.353234707 3.63 -7.347373432 3.64 -7.340605242 3.65 -7.332953213 3.66 -7.324439984 3.67 -7.315087768 3.68 -7.304918353 3.69 -7.293953116 3.7 -7.282213024
上面的数据实际上对应的就是结合能曲线,其能量最低点对应的就是平衡晶格常数。当然,能量最低点的能量就是晶体的结合能。
具体如何从上面的数据中获得最终结果,可以参考原作者之前写的一篇文章:MATLAB计算平衡晶格常数。
正是我的第一个运行的程序,很开心
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2024-9-25 15:01
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社