Some interesting issues encountered when simulating the adsorption of atoms on iron-carbon surfaces using the ReaxFF force field in LAMMPS with the GCMC method

Dear all
I have modified my file( fix_gcmc.cpp) according to Gcmc region triclinical box error - #2 by akohlmey
The follows are my input file

boundary     p  p  p 
units        real 
atom_style   charge 
read_data    1.data
mass 1 55.845  # Fe
mass 2 12.0107  # C
mass 3 1.00794  # H
mass 4 15.9994  # O
mass 5 12.0107  # C
timestep       0.1
pair_style    reax/c   NULL  
pair_coeff    * * ffield.reax.Fe_C_O Fe C H O C          
neighbor      2.0 bin  
thermo_style      custom step atoms vol temp press etotal pe evdwl ecoul 
thermo            100
velocity all create 500.0 4928459
fix qeq all qeq/reax 1 0.0 10.0 1e-6 reax/c #
region          bottom block INF INF INF INF INF 6   units box#Fe5C2 fix
group           bottom region bottom
velocity        bottom set 0 0 0  
fix             bottom bottom setforce 0 0 0 
region          mobile1 block INF INF INF INF 6.1 13 #Fe5C2 unfix
group           mobile1 region mobile1 #gcmc 
region          mobile2 prism 0.0 18.2908 0.0 16.3580969934 6.1 20 -4.5727105522 0.0 0.0 side in #Fe5C2 unfix  
group           mobile2 region mobile2 #gcmc
group           FeC type 1:2 # 1=Fe, 2=C
group           H type 3   # 3=H
group           O type 4   # 4=O
group           C type 5   # 5=C
group           gas type 3:5 #gas

compute          Htemp H temp
compute_modify   Htemp dynamic/dof yes
fix              nvtH H nvt temp 500.0 500.0 10.0
fix_modify       nvtH temp Htemp 
fix             g3 gas gcmc 10 10 10 3 7641 500 -5.46 0.1 pressure 28  region mobile2 #group H

compute          Otemp O temp
compute_modify   Otemp dynamic/dof yes
fix              nvtO O nvt temp 500.0 500.0 10.0
fix_modify       nvtO temp Otemp 
fix             g4 gas gcmc 10 10 10 4 7641 500 -5.46 0.1 pressure 28  region mobile2 #group O

I encountered two interesting situations. The first was that if I did not annotate “fix qeq all qeq/reax 1 0.0 10.0 1e-6 reax/c ”, it would report an error (ERROR: Fix gcmc region extends outside simulation box (…/fix_gcmc.cpp:145)) when the first ‘gcmc’ simulation was conducted. If I annotate “fix qeq all qeq/reax 1 0.0 10.0 1e-6 reax/c ”, it would report an error (ERROR: Fix gcmc region extends outside simulation box (…/fix_gcmc.cpp:145)) when the second ‘gcmc’ simulation was conducted. Request help to answer the questions below
1.Where did I go wrong?
2.Why does balancing charges lead to exceeding the simulated range?
3.The parameters of gcmc are exactly the same twice. The first time it works, but the second time it doesn’t. What’s the matter?
4.Is it not possible to perform two GCMC simulations simultaneously in one simulation?

Last time I checked it was necessary to use the more memory-robust Kokkos version of ReaxFF with GCMC since the regular ReaxFF crashed, are you doing that, or maybe it works now? However that seems unrelated to the issues above.

It seems that I have solved this problem. By chance, I replaced the previous version of Lammps, which was “7Aug2019,” with “2Aug2023,” and since then, this kind of error has not continued to occur.
I apologize for any inconvenience caused.

1 Like