Carlos,
I have attached a complete input script and data file for a system that
fails with a lost atoms error. As far as testing with other codes, I have
used an old in-house MC code without issues for this system at higher
concentrations, where the lammps GCMC fix frequently runs into lost atom
problems. Most of the low concentration systems complete without issue.
Thanks,
Chris
variable t equal 300.00
variable temp string "300"
variable seed equal 31798.000000
variable phi string "q"
variable numH equal 200
variable mu equal -2.4000e+00
variable bulkH equal 3.0000e-05
variable e equal -0.01
units metal
boundary p f p
atom_style atomic
### set up simulation domain
read_data phi_\{phi\}\.data
change\_box all x scale 1\.004347451179 y scale 1\.004347451179 z scale
1\.004347451179
variable rxlo equal xlo\+0\.0001
variable rxhi equal xhi\-0\.0001
variable rzlo equal zlo\+0\.0001
variable rzhi equal zhi\-0\.0001
region active block {rxlo} \{rxhi\} 47\.0000 107\.0000 {rzlo}
${rzhi} side in
group hydrogen type 4
group nickel type 1 2 3
# EAM potential -- Ni
pair_style eam/alloy
pair_coeff * * /Users/cjobrie/lammps-pot/NiH_jea.eam.alloy Ni Ni Ni H
mass * 58.710
mass 4 1.0080
#re-neighbor
neighbor 1.5 bin
neigh_modify every 5 delay 0
#THERMALLY EQUILIBRATE
thermo 0
log logqequ.lammps
velocity all create t {seed} mom yes rot yes dist gaussian
fix 7 all nvt temp $t t 10\.0
timestep 0\.001
run 100000
reset\_timestep 0
\#GCMC SETTINGS
\# fID gID gcmc Nmdsteps Xchange Moves attype seed T mu
maxtrans options
fix gcmc hydrogen gcmc 100 {numH} 0 4 \{seed\} {t}.0 \{mu\}
0\.01 region active full\_energy
compute\_modify thermo\_temp dynamic yes
variable nh equal count\(hydrogen\)
fix avenh all ave/time 50 100 5000 v\_nh ave one
thermo 100000
thermo\_style custom step temp pe f\_avenh
thermo\_modify format float %24\.16g
restart 100000 restarteq\.{phi}.*
timestep 0.001
#re-neighbor
neighbor 1.5 bin
neigh_modify every 1 delay 0
run 3000000
#GCMC AVERAGING RUN
unfix avenh
variable cex equal "(count(hydrogen) -
count(nickel,active)*v_bulkH) / (lx*lz)"
unfix avenh
variable cex equal "(count(hydrogen) -
count(nickel,active)*v_bulkH) / (lx*lz)"
variable nh equal count(hydrogen)
variable tau equal (pxx+pzz)/2
compute PE all pe/atom
compute mype all pe
fix avetau all ave/time 2 250 1000 v_tau start 1000001 ave
running
fix avetemp all ave/time 2 250 1000 c_thermo_temp start
1000001 ave one
fix avecex all ave/time 2 250 1000 v_cex start 1000001 ave
running
fix avepe all ave/time 2 250 1000 c_mype start 1000001 ave
running
log logqave.lammps
dump 4 all custom 1000 dump.phi_\{phi\}\.\* id type x y z c\_PE
dump\_modify 4 sort id pad 7
restart 20000 restart\_{phi}
thermo 1000
thermo_style custom step f_avetemp pxx pzz f_avetau f_avepe v_nh v_cex
f_avecex
thermo_modify format float %24.16g
timestep 0.001
run 200000
phi_q.data (276 KB)