Practice of histogram of vacancy formation energy and migration energy

Dear lammps developers and users;
I recently encountered a problem. I read the simulation of vacancy formation energy and migration energy in “influence of compositional complexity on species diffusion behavior in high entropy solid solution alloys” in this article, and carefully read the histogram of vacancy formation energy of alloys in this article. I wrote an in file for calculation of vacancy formation energy of each element as follows:

    ##定义变量##
variable        a equal 2.0 ##定义原子模型截断半径
##模拟环境初始化##
units          metal
atom_style     atomic
boundary       p p p
neighbor       ${a} bin
neigh_modify   every 1 delay 0 check yes
##建立单晶原子模型##
lattice               fcc 3.56
region                box block 0 30 0 30  0 30 units box
create_box            3 box
create_atoms          1 box
####将Ni原子进行替换##
set                          type 1 type/subset 2 818 666 ##替换成Cr
set                          type 1 type/subset 3 818 666 ##替换成Co

##定义相互作用势##
pair_style       eam/alloy
pair_coeff       * * NiCoCr.lammps.eam Ni Cr Co


##对原子进行分组
group          Ni-1 type 1
group          Cr-2 type 3
group          Co-3 type 3
##计算各元素的总势能
compute        1 Ni-1 pe/atom
compute        pe1 Ni-1 reduce sum c_1
compute        2 Cr-2 pe/atom
compute        pe2 Cr-2 reduce sum c_2
compute        3 Co-3 pe/atom
compute        pe3 Co-3 reduce sum c_3

#输出原子热力学信息和位置
thermo           100
thermo_style     custom step lx ly lz c_pe1 c_pe2 c_pe3  temp vol dt time


###结构优化得到E_perfect能量
reset_timestep         0
minimize         1.0e-12 1.0e-12 10000 10000

###记录第一次优化后的能量和总原子数目
variable         E1 equal c_pe1
variable         E2 equal c_pe2
variable         E3 equal c_pe3         
variable         E1_perfect equal ${E1}
variable         E2_perfect equal ${E2}
variable         E3_perfect equal ${E3}
variable         N equal count(all)
variable         N0 equal ${N}
##对原子模型进行删除形成缺陷
region          2 sphere 10.0 10.0 10.0 0.1 side in
delete_atoms    region 2

write_data      NiCrCo.data

##对删除原子之后进行分组
group          Ni-4 type 1
group          Cr-5 type 3
group          Co-6 type 3
##计算各元素的总势能
compute        4 Ni-4 pe/atom
compute        pe4 Ni-4 reduce sum c_4
compute        5 Cr-5 pe/atom
compute        pe5 Cr-5 reduce sum c_5
compute        6 Co-6 pe/atom
compute        pe6 Co-6 reduce sum c_6

##结构优化得到E_defect的能iang
reset_timestep   0
thermo           100
thermo_style     custom step lx ly lz c_pe4 c_pe5 c_pe6 temp vol dt time
minimize        1.0e-12 1.0e-12 10000 10000

###记录第二次优化的能量值

variable        E4_defect equal c_pe4
variable        E5_defect equal c_pe5
variable        E6_defect equal c_pe6
           
###通过公式分别计算元素Ni、Cr和Co的空位形成能
variable       Ev_Ni equal ((1/${N0})*${E4_defect}-((${N0}-1)/${N0})*${E1_perfect})
variable       Ev_Cr equal ((1/${N0})*${E5_defect}-((${N0}-1)/${N0})*${E2_perfect})
variable       Ev_Co equal ((1/${N0})*${E6_defect}-((${N0}-1)/${N0})*${E3_perfect})

###导出直方图
fix            8 all ave/histo/weight  100 5 1000 -5 5 100 c_Ev_Ni mode vector file Ni.pe
fix            9 all ave/histo/weight  100 5 1000 -5 5 100 c_Ev_Cri  mode vector file Cr.pe
fix            10 all ave/histo/weight  100 5 1000 -5 5 100 c_Ev_Co mode vector file Co.pe

timestep   0.001
run            0

I can’t define the calculation formula of the ordinate of this figure in the in file.The figure is as follows:
Uploading: image.png…
I wonder if any of you know the method of this drawing or the specific calculation formula for changing the ordinate of the drawing. I hope you can get your help!
xue bao shuai

Hi
You have two times
group …. type 3
It could generate an issue in your simulation.
Best
Pascal

Hi,
First of all, thank you very much for your timely answers! I did not upload the vacancy formation energy and migration energy histogram in time.This figure is attached. Do you know how the ordinate in this figure is implemented in lammps?
This figure is from the article “influence of composite complexity on species diffusion behavior in High-enterprise solid solution alloys”.I wonder if you can give me some other suggestions. Thank you

Hi,
First of all, thank you very much for your timely answers! I did not upload the vacancy formation energy and migration energy histogram in time.
This figure is from the article “influence of composite complexity on species diffusion behavior in High-enterprise solid solution alloys”.I wonder if you can give me some other suggestions. Thank you