Create_atoms for H2 and O2 pairs in REAXFF

To learn REAXFF and REACTION, im making a minimal example called in.hydrogen_combustion. id like to create H2 and O2 molecules but since reaxff is a bond-order potential there is no molecular topology. how can i create pairs of hydrogen 0.74A apart and pairs of oxygen 1.21A apart in a 2:1 ratio using create_atoms ?

if i create random H and O atoms, then even after a few timesteps i still have mostly H and O species, only a few H2 and O2.

#  Timestep    No_Moles    No_Specs           H          H2          HO         H2O         H3O         HO2           O          O2
       2500        2712           8        1627         118         117           8           1           1         828          12

Another thing i dont understand from the literature and examples ive seen is why use a thermostat to hold a constant temperature which doesnt model the chain reaction increase in energy explosive nature of this exothermic reaction ??

This is what i have so far:

units real
dimension 3
boundary p p p
atom_style charge
newton on

region box block 0 99 0 99 0 99
create_box 2 box

mass 1 1.008
mass 2 15.999

create_atoms 1 random 2000 12345 NULL overlap 0.31 maxtry 100
create_atoms 2 random 1000 12345 NULL overlap 0.66 maxtry 100
velocity all create 2500 12345

pair_style      reaxff lmp_control
pair_coeff      * * ffield.reax.cho H O

neighbor        2 bin
neigh_modify    every 10 delay 0 check no

variable dt equal 0.1
timestep        ${dt}

fix             1 all nvt temp 2500 2500 $(100.0*dt)
fix             2 all qeq/reax 1 0.0 10.0 1e-6 param.qeq
fix             3 all reaxff/species 10 10 100 species.out

thermo 100
dump 1 all movie 100 hydrogen_combustion.mkv type type size 800 800
dump_modify 1 acolor * white/red/green/blue/aqua/magenta

run           $(500*1000/dt)

Molecule files?

  1. Create isolated atoms with create_atoms
  2. Group the existing atoms
  3. replicate the box in the x-direction
  4. Use group ... subtract to group newly created H and O atoms
  5. Now displacing the newly created atoms leftward by one box length minus your desired “bond” length gives you the desired dimers. (You can displace different types by different amounts using type-based grouping.)
  6. Use change_box to set the box back to its original size.

Hey, I’m just saying how I would do it. That doesn’t make it a good idea. :smiling_imp:

actually the solution was molecule files as suggested by @akohlmey. sorry i didnt have time earlier today to post my solution. i do appreciate your effort @srtee however with such a brute force :hammer: !

at first i thought i couldnt use molecule command because REAXFF doesnt have molecular topology, ie. no bonds.

my trick is to use molecule files but with NO bonds as a workaround, so here it is for future googlers:

# H2.txt molecule file, with no bonds as workaround for REAXFF example

2 atoms

Coords

1    0.0000    0.0000    0.0000
2    0.7414    0.0000    0.0000

Types

1        1
2        1
# O2 molecule file, with no bonds as workaround for REAXFF example

2 atoms

Coords

1    0.0000    0.0000    0.0000
2    1.2075    0.0000    0.0000

Types

1        2
2        2
mass 1 1.008
mass 2 15.999

molecule H2 H2.txt
molecule O2 O2.txt

create_atoms 0 random 2000 12345 NULL overlap 1.0 maxtry 100 mol H2 12345
create_atoms 0 random 1000 12345 NULL overlap 1.0 maxtry 100 mol O2 12345