Problem regarding TIP4P/SPC water model

LAMMPS LAMMPS Mailing List Mirror
1

Hello everyone,

Thank you Axel for helping me. I resolved those issues. However, when I am simulating with 100 molecules of water, the simulation works fine. I checked from 3, 10, 50, 100 molecules. But when I increased to 120 molecules, the output shows NAN in all thermodynamic output in every timestep extracted from thermo. I know I can have several issues here. But why it would work for 100 molecules but not with 120 molecules?The box dimensions are 22x22x22 angstrom. timestep is 1 fs, Tdamp=100*timestep. I have changed parameters in neigh_modify also. Here is my input file:

SPC water

echo screen
log debug.log

units real

neighbor 2.0 bin
neigh_modify delay 0 every 1 one 10000 check yes

boundary p p p

atom_style full
bond_style harmonic
angle_style harmonic
dihedral_style none
improper_style none

pair_style lj/cut/coul/cut 9.0 9.0

read_data water_diff.data

pair_modify mix arithmetic tail yes

group hydrogen type 1
group oxygen type 2
group water type 1 2

set group oxygen charge -0.82
set group hydrogen charge 0.41

pair_coeff 1 1 0.0 0.0
pair_coeff 1 2 0.0 0.0
pair_coeff 2 2 0.1553 3.166

bond_coeff 1 0.0 1.0
angle_coeff 1 0.0 109.47

velocity water create 300.0 12345 dist gaussian

fix 1 water shake 1.0e-8 100 0 b 1 a 1
fix 2 water npt temp 300.0 300.0 100.0 iso 1.01325 1.01325 50

Tdamp=100*timesteps

variable timesteps equal 1 #femtosecond
timestep ${timesteps}

thermo 50
thermo_style custom step temp press vol etotal ke pe evdwl ecoul elong

write_data data.spc_out pair ij
dump 1 all xyz 1 spc.xyz
dump 2 all image 100 dump.*.jpg type type

run 1000

Thank you Axel for helping me. I resolved those issues. However, when I am simulating with 100 molecules of water, the simulation works fine. I checked from 3, 10, 50, 100 molecules. But when I increased to 120 molecules, the output shows NAN in all thermodynamic output in every timestep extracted from thermo. I know I can have several issues here. But why it would work for 100 molecules but not with 120 molecules?The box dimensions are 22x22x22 angstrom. timestep is 1 fs, Tdamp=100*timestep. I have changed parameters in neigh_modify also. Here is my input file:

SPC water

echo screen
log debug.log

units real

neighbor 2.0 bin
neigh_modify delay 0 every 1 one 10000 check yes

boundary p p p

atom_style full
bond_style harmonic
angle_style harmonic
dihedral_style none
improper_style none

pair_style lj/cut/coul/cut 9.0 9.0

read_data water_diff.data

pair_modify mix arithmetic tail yes

group hydrogen type 1
group oxygen type 2
group water type 1 2

set group oxygen charge -0.82
set group hydrogen charge 0.41

pair_coeff 1 1 0.0 0.0
pair_coeff 1 2 0.0 0.0
pair_coeff 2 2 0.1553 3.166

bond_coeff 1 0.0 1.0
angle_coeff 1 0.0 109.47

velocity water create 300.0 12345 dist gaussian

fix 1 water shake 1.0e-8 100 0 b 1 a 1
fix 2 water npt temp 300.0 300.0 100.0 iso 1.01325 1.01325 50

Tdamp=100*timesteps

variable timesteps equal 1 #femtosecond
timestep ${timesteps}

thermo 50
thermo_style custom step temp press vol etotal ke pe evdwl ecoul elong

write_data data.spc_out pair ij
dump 1 all xyz 1 spc.xyz
dump 2 all image 100 dump.*.jpg type type

run 1000

Hello everyone,

Thank you Axel for replying to me. I changed the input file for SPC model to
a low density water box. Now, I am having numeric index out of bound errors.
I defined types of atoms, bonds, angles, impropers and dihedrals in the data
file, and the corresponding coefficients are in the input file. I don’t
understand where my atom IDs might be wrong because there’s only 9 atoms in
the box. Can someone please help me? I have copied both input file and data
file for you to see.

look closely at the output to the screen and note the line in the
input where the error occurs.
it has nothing to do with atom indices, but bond (and angle) type.
your data file has only one type (1), but you are setting parameters
for type 2. that clearly means the type index is out of bounds.

axel.

SPC water box (3 molecules)

SPC water

echo screen
log debug.log

units real

neighbor 2.0 bin
#neigh_modify delay 0 every 1 one 10 check yes
neigh_modify delay 5

boundary p p p

atom_style full
bond_style harmonic
angle_style harmonic
dihedral_style none
improper_style none

pair_style lj/cut/coul/cut 9.0 9.0

read_data data.lowdensity
pair_modify mix arithmetic tail yes

group hydrogen type 1
group oxygen type 2
group water type 1 2

set group oxygen charge -0.82
set group hydrogen charge 0.41

pair_coeff 1 1 0.0 0.0
pair_coeff 1 2 0.0 0.0
pair_coeff 2 2 0.1553 3.166

bond_coeff 2 0.0 1.0
angle_coeff 2 0.0 109.47

velocity water create 300.0 12345 dist gaussian

fix 1 water shake 1.0e-8 100 0 b 1 a 1
fix 2 water npt temp 300.0 300.0 100.0 iso 1.01325 1.01325 500

timestep 1.0
thermo 50
thermo_style custom step temp press vol etotal ke pe evdwl ecoul elong

dump 1 all xyz 1 spc.xyz
run 10000

DATA file:
LAMMPS data file. CGCMM style. atom_style full generated by VMD/TopoTools
v1.5 on Tue Dec 01 16:55:25 CST 2015
9 atoms
6 bonds
6 angles
0 dihedrals
0 impropers
2 atom types
1 bond types
1 angle types
0 dihedral types
0 improper types

[…]

2

Hello everyone,
Thank you Axel for helping me. I resolved those issues. However, when I am simulating with 100 molecules of water, the simulation works fine. I checked from 3, 10, 50, 100 molecules. But when I increased to 120 molecules, the output shows NAN in all thermodynamic output in every timestep extracted from thermo. I know I can have several issues here. But why it would work for 100 molecules but not with 120 molecules?The box dimensions are 22x22x22 angstrom.

You may just have created a starting configuration with very high potential energy. That is the most common reason for NaN thermo data.

There is not much of a point to use fix npt for mostly empty simulation cells.

Axel

timestep is 1 fs, Tdamp=100*timestep. I have changed parameters in neigh_modify also. Here is my input file:

3

Hello everyone,

I resolved the NaN error by applying timestep 0.5 femtosecond. now the simulation is stable. Thank you Axel for helping me so far.

Regards,
Baig

Hello everyone,
Thank you Axel for helping me. I resolved those issues. However, when I am simulating with 100 molecules of water, the simulation works fine. I checked from 3, 10, 50, 100 molecules. But when I increased to 120 molecules, the output shows NAN in all thermodynamic output in every timestep extracted from thermo. I know I can have several issues here. But why it would work for 100 molecules but not with 120 molecules?The box dimensions are 22x22x22 angstrom.

You may just have created a starting configuration with very high potential energy. That is the most common reason for NaN thermo data.
There is not much of a point to use fix npt for mostly empty simulation cells.
Axel

timestep is 1 fs, Tdamp=100*timestep. I have changed parameters in neigh_modify also. Here is my input file:

4

Hello everyone,

I resolved the NaN error by applying timestep 0.5 femtosecond. now the
simulation is stable. Thank you Axel for helping me so far.

that confirms my hunch that you have a bad (i.e. unrelaxed) initial geometry.
you may consider adjusting the settings for packmol to space the
molecules a bit farther apart, or run a minimization before doing MD.
also keep in mind that nose-hoover thermostats are not designed to
effectively equilibrate configurations that are far away from
equilibrium. using fix nve in combination with a dissipative
thermostat with strong damping is often a better choice for the
initial part of equilibration. with fix shake on water models in a
relaxed configuration, you should be able to run a stable MD with up
2fs for the timestep at room temperature.

axel.