Simple simulation of TIP4P2005/cut water with NVE ensemble

Dear all,

Rencently, I changed my previous simulation of water from NPT ensemble to NVE ensemble. However, the calculation crashes quickly with TotalEnergy becoming NaN. I wonder if anyone has an idea since the previous simulation is correct and produces similar result to journal result.

Thanks in advance.

The input file is as the following. The initial water configuration is generated by moltemplate with fcc structure of 256 water molecules and 1g/cm3 density. Input file, log file and initial configuration are attached.

water.equil.data (45.7 KB)

in.water.equil (1.26 KB)

log.lammps (3.22 KB)

first off, this really has little to do with LAMMPS, but is about
understanding MD in general, and thus is mainly off-topic for this
list. find yourself an adviser/tutor that can instruct you on how to
do MD properly and how to assess the quality of an MD simulation. the
way this is typically done (or at least it used to be this way when i
learned MD), is that you attend a suitable (summer) school and then
join a group with the necessary expertise for several weeks/months so
you can learn from the group experts.

Dear all,

Rencently, I changed my previous simulation of water from NPT ensemble to
NVE ensemble. However, the calculation crashes quickly with TotalEnergy
becoming NaN. I wonder if anyone has an idea since the previous simulation
is correct and produces similar result to journal result.

"correct" or "not correct" is a matter of definition. i would consider
*both* simulations "bad". you cannot conserve energy (please note the
significant discrepancy between the desired temperature and the actual
average temperature, when you turn on the npt fix. with fix nve, you
have no longer the option to exchange excess energy with the
environment and hence you see your system's kinetic energy build until
the overall velocity is too large for the chosen time step and the
numerical integration become unstable. part of that is due to the
large time step, part it is due to using an unshifted cutoff for both
LJ and Coulomb. turning on long-range electrostatics, should improve
the situation significantly, also a shorter time step should help (but
not as much).

axel.