mdrun-openmm produces nans

[Dirk Colbry]

I have been having some trouble getting pre-compiled mdrun to produce the correct output (gromacs 4.0.5 and OpenMMPreview4 on sles11). The program runs successfully but prodces "nan"s (not a number) in the output instead of the correct numrical values.

I also have a standard version of gromacs installed on the same system and mdrun produces the correct values with the same input file.

I was able to successfully compile GromacsOpenMM from source on sles10 but I get the same behavior as the precompiled binaries.

I am not sure what the best direction to continue my debugging. Has anyone seen similar behavior? Does anyone have any suggestions as to what I should look into next?

Thank you,

- Dirk

-----------------------------------------------
Dr. Dirk Joel Luchini Colbry
Research Specialist
Institute for Cyber-Enabled Research (iCER)
Michigan State University
web: http://www.dirk.colbry.com

[Romica Sandu]

I have seen that mdrun-openmm produces nans when initial system was not sufficient minimized.

Romica

[Peter Eastman]

Look at the trajectory for your simulation. Does it immediately go to nan's after one time step? Or does it run correctly for a while first? Often you'll see a simulation blow up over a few time steps. A few atoms will jump out of place, then the errors will quickly propagate through the rest of the system. That's especially true for an insufficiently minimized starting structure. Also make sure you aren't relying on any features of Gromacs that OpenMM doesn't support (vsites, triclinic boxes, etc.).

Peter

[Lisa Lapidus]

Hi Peter,

Yes, the .gro file is nan from the very beginning. I tried increasing the number of steps for energy minimization from 20000 to 100000 and get the same output. We also tried using the villin example and did get numerical output, but the radius of gyration was always zero.

Lisa Lapidus

[Peter Eastman]

> Yes, the .gro file is nan from the very beginning.

Make sure you have nstxout set to 1. Otherwise, all we know is that it happens within the first reporting interval.

What do the forces look like in the first time step? Are they reasonable values, or huge, or nan?

Peter

[Dirk Colbry]

Peter and Romica,

Thank you for all of your help. We where able to get the minimization to work correctly by running pdb2gxm and mdrun using the standard build of gromacs (4.0.5) that was not compiled with openmm. We then did the main simulation with the versions compiled with openmm (Preview4) and everything seemed to run fine.

I am not sure why the standard minimization seemed to work while the openmm minimization did not.

I recently rebuilt the newest version of gromacs (4.0.7) using the newest version of openmm (1.0beta) and my preliminary test show that this build seems to be working much correctly. I am now looking at my build procedure to see if the problem is related to a specific build option or if the newest version of the program somehow fixes this problem.

In either case, we seem to be working so thank you again for your help.

- Dirk

[Peter Eastman]

Hi Dirk,

Gromacs-OpenMM does not currently support energy minimization, which is why you need to use standard Gromacs for it. This is discussed in the readme file, but it's worth emphasizing again: the *only* supported integrators are "md", "sd", and "bd". So if you specify one of the minimization integrators like "steep" or "cg", it will not work correctly!

Peter

[Kyle Hong]

Dear all,
I have the same problem. My mdrun-openmm produces nans, but mdrun is ok. I set nstxout to 1 as Peter saids to check what happens, but the x output is normal only 17 steps (my simulation is 50000 steps) ,then is nan, and the other ouputs (et energy and v)is still nan . I can compile all the exampls download from this website perfectly,why my system goes wrong?

My md.mdp file is:

integrator = sd
ld_seed = -1
nsteps = 50000
dt = 0.002
ns_type = Grid
nstlist = 10
rlist = 1.0
coulombtype = pme
rcoulomb = 1.0
vdw-type = cut-off
rvdw = 1.4
Tcoupl = no
tc_grps = protein sol CL-
tau_t = 0.1 0.1 0.1
ref_t = 310 310 310
Pcoupl = no
gen_vel = no
nstenergy = 10
nstlog = 10
nstxtcout = 1

Any help will be appreciate.
Kyle

[Peter Eastman]

So the values are normal for the first 17 steps, and then turn into nans? Look through the values just before that happens. (You'll want to output the forces and velocities at every step, not just the positions.) Can you see what is happening? For example, you might see that two atoms are coming very close, leading to a very large force. This often indicates insufficient energy minimization or too large a step size.

In particular, I notice that you're using a time step of 0.002 without any constraints, which is generally not stable. Try reducing the time step to 0.001, or setting constraints=hbonds.

Peter

[Romica Sandu]

You probably have in system more than one molecule.

Romica