Hello all,
I am trying to find out equilibrium lattice parameter of Cu by applying NPT thermodynamics. However, it keeps on changing in output. With each given input of lattice constant, it is coming out to be different How to find accurately the Equilibrium lattice constant?
Following is the program:
---------- Initialize Simulation ---------------------
clear
units metal
dimension 3
boundary p p p
atom_style atomic
atom_modify map array
---------- Create Atoms ---------------------
lattice fcc 3.631618387
region box block 0 1 0 1 0 1 units lattice
create_box 1 box
lattice fcc 3.631618387 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
create_atoms 1 box
replicate 10 10 10
---------- Define Interatomic Potential ---------------------
pair_style eam
pair_coeff * * Cu_u3.eam
neighbor 2.0 bin
neigh_modify delay 10 check yes
reset_timestep 0
compute new all temp
compute eng all pe/atom
compute eatoms all reduce sum c_eng
velocity all create 300 2223333 temp new
fix 1 all npt temp 300 300 0.2 iso 0.0 0.0 1
###fix 2 all box/relax iso 0.0 vmax 0.001
#------ output the lattice constant -----#
thermo 100
#thermo_style custom step temp lx ly lz
thermo_style custom step temp press pxx pyy pzz lx ly lz c_eatoms
run 200000
variable natoms equal “count(all)”
variable teng equal “c_eatoms”
variable a equal “lx/10”
variable ecoh equal “v_teng/v_natoms”
print “Lattice constant(Angstroms) = {a};"
print "No. of atoms = {natoms};”
print “Total energy (eV) = {teng};"
print "Cohesive energy (eV) = {ecoh};”
variable length equal “lx/10”
print “Lattice constant (Angstoms) = ${length};”
print “All done!”