Error: Molecule topology/atom exceeds system topology/atom

LAMMPS LAMMPS Beginners
1

I want to simulate a mixture of one monoatomic and diatomic molecules. I am using this script.

atom_style full
units real
timestep 10

# set up the system
boundary p p p
region simbox block 0 50 0 50 0 50
create_box 3 simbox &
bond/types 1
# set up variables
variable natoms   equal 4000
variable molfrac  equal 0.1
variable nmet     equal ${natoms}*${molfrac}
variable nwater   equal ${natoms}-${nmet}
bond_style harmonic
bond_coeff 1 1.43 77200
molecule met molecule.dat

# force field
pair_style          hybrid sw lj/cut 9.5
pair_coeff          * * sw mW_real.sw NULL mW mW  
pair_coeff          1 2 lj/cut 3.7425 0.1876 9.5
pair_coeff          1 1 lj/cut 3.7750 0.2070 9.5 # CH3
pair_coeff          1 3 lj/cut 3.7425 0.1876 9.5 

# create water and methanol
create_atoms 0 random ${nmet} 12432513 simbox mol met 123324
create_atoms 3 random ${nwater} 12345678  simbox
mass 3 18.0350  
mass 1 15.0350  
mass 2 18.015

My molecule.dat is this:

# Methanol UA geometry
# header section:
2 atoms
1 bonds

# body section:
Coords

1 52.80  38.18  48.77
2 53.36  36.92  49.15

Types

1        1   
2        2  

Charges

1  0.2650
2 -0.2650

Bonds

1   1      1      2

I get this error:

ERROR: Molecule topology/atom exceeds system topology/atom (src/src/molecule.cpp:1857)
Last command: create_atoms 0 random ${nmet} 12432513 simbox mol met 123324

I am unable to see why is this error coming because my system should only have three atom types. Please let me know how to resolve this.

Thank you so much.

2

In addition to storing atom types (and bond types), you also need to store the definition of those bonds. In LAMMPS bond topology date is stored with individual atoms so that the domain decomposition parallelization is more efficient, especially for very large systems running across a very large number of processors (this is what the letters “LA” in the LAMMPS acronym stand for).

In most cases people start simulations of molecular systems from a data file and when using read_data, LAMMPS determines from the data file how much data needs to be reserved per atom for bonds, angles, dihedrals, impropers and exclusions (aka special neighbors). When using “create_box” this is not the case and you need to provide the numbers. You already provide the number of bond types, but you will also provide other properties. In your case this would be “extra/bond/per/atom” and “extra/special/per/atom”. These are missing and thus assumed to be 0 and this is what the molecule command complains about. For more details, please see the create_box documentation.