表面能的计算 & lammps输入脚本

大家好，我是木南

---

-------------------------------------------------------------------

1. 什么是表面能（Surface energy）？

A. 晶体表面的原子处于不平衡的状态，其能量要高于处于平衡态的晶体内部原子，高出的这部分能量就是表面能。

B. 将体相变成表面晶体所需要的能量。

C. 表面能是创造物质表面时，破坏分子间作用力所需消耗的能量。

1. 计算过程：

首先计算平衡态时的晶体总的势能E0，然后把晶体从中间切开，形成两个表面，把晶体移开足够远的距离，计算此时系统的总势能E1。

则该晶体的表面能E=（E1-E0）/2A,其中A为表面面积。

因为用到了，所以做一个总结，记录一下。

对原作者表示感谢！

参考资料：

1.lammps案例：Cu（111）表面能的 

  https://zhuanlan.zhihu.com/p/347255328

2.表面能的理解+lammps计算

  https://www.cnblogs.com/sysu/p/14279286.html

以下是Cu（100）、（110）和（111）面的表面能计算in文件：

1. ##----------------------基本参数设置----------------------
   units     metal
   boundary   p p p
   atom_style   atomic
   
   ##----------------------时间步长和近邻原子设置----------------------
   timestep 0.005
   neighbor 0.6 bin
   neigh_modify every 5 delay 0 check yes
   
   
   
   #计算Cu（100）面表面能，定义晶体取向
   lattice     fcc 3.615 origin 0.0 0.0 0.0 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1 
   
   
   #计算Cu（110）面表面能，定义晶体取向
   #lattice     fcc 3.615 origin 0.0 0.0 0.0 orient x 1 1 0 orient y 0 0 1 orient z 1 -1 0 
   
   
   #计算Cu（111）面表面能，定义晶体取向
   #lattice fcc 3.615 origin 0.0 0.0 0.0 orient x 1 1 1 orient y 1 -1 0 orient z 1 1 -2 
   
   
   ##----------------------box大小----------------------
   region    box block 0 20 0 20 0 40
   
   
   ##----------------------填充Cu原子----------------------
   create_box   1 box
   create_atoms 1 box
   
   
   ##--------势函数的选取-----------------------------------
   pair_style eam/alloy 
   pair_coeff * * Cu_mishin1.eam.alloy Cu 
   
   
   ##--------定义上下两个表面的边界区域-----------------------------------
   region boundary1 block INF INF INF INF 29.9 INF
   region boundary2 block INF INF INF INF INF 9.9
   group boundary1 region boundary1
   group boundary2 region boundary2
   group boundary union boundary1 boundary2
   
   
   #周期性边界条件+真空层=自由表面
   delete_atoms group boundary
   
   
   ##--------输出设置，每隔10步进行输出----------------
   thermo 10
   dump              1 all custom 1000 pos.lammpstrj id type x y z
   dump_modify       1 sort id
   
   
   ##--------采用最速下降算法进行能量最小化----------------
   min_style sd
   minimize 1.0e-30 1.0e-15 10000 10000
   
   
   #系统达到平衡后，输出能量E0
   variable E equal pe
   print "---------------------- E0=$E---------------------------"
   run 0
   
   
   #定义上下区域，目的是把晶体切为上下两部分
   region down block INF INF INF INF INF 19.94
   region up block INF INF INF INF 19.95 INF
   group up region up
   group down region down
   
   
   #将box在z方向上的尺寸扩大40
   change_box all z delta 0 40 units lattice
   
   #将上部分原子沿z轴移动40个晶格
   displace_atoms up move 0 0 40 units lattice
   
   
   #再次能量最小化，并输出新的能量E1
   minimize 1.0e-30 1.0e-20 10000 100000
   
   print "--------------E1=$E--------------------------------"

2. 

   以上是我们分享的一些经验或者文章的搬运，或有不足，欢迎大家指出！

    欢迎大家关注微信公众号： 原子与分子模拟