[gmx-users] energy is coming "nan"

[Justin A. Lemkul]

sarbani chattopadhyay wrote:
> Hi ,
> I am trying to simulate a protein in vacuum for 10 ns.
> However at around 8.8 ns the log file shows the following
> Step Time Lambda
> 4411700 8823.40000 0.00000
> 
> Energies (kJ/mol)
> Bond Angle Proper Dih. Ryckaert-Bell. LJ-14
> 1.37934e+04 1.02343e+04 6.09706e+02 6.42224e+03 7.97160e+03
> Coulomb-14 LJ (SR) Coulomb (SR) Potential Kinetic En.
> 4.08238e+04 -1.23373e+04 -1.31884e+05 -6.43660e+04 2.34001e+04
> Total Energy Temperature Pressure (bar)
> -4.09659e+04 2.93347e+02 0.00000e+00
> 
> Large VCM(group rest): -0.15115, -0.85874, -0.88165, Temp-cm: 3.61048e+00
> Large VCM(group rest): -3.30356, -1.50152, -11.40962, Temp-cm: 3.36602e+02
> Group rest with mass 4.50563e+04, Ekrot 3.33396e+19 Det(I) = 4.98670e+95
> COM: 4.40566 4.40778 4.74707
> P: -148846.33496 -67653.09375 -514075.76562
> V: -3.30356 -1.50152 -11.40962
> J: 457896798773861765087232.00000 1805633659184017884315648.00000
> 86113127573044898824192.00000
> w: -0.00006 0.00002 0.00075
> Inertia tensor (3x3):
> Inertia tensor[ 0]={ 7.62652e+36, -2.04805e+36, -9.25856e+37}
> Inertia tensor[ 1]={-2.04805e+36, 5.49988e+35, 2.48632e+37}
> Inertia tensor[ 2]={-9.25856e+37, 2.48632e+37, 1.12399e+39}
> Large VCM(group rest): 0.16121, 0.70906, 46.02300, Temp-cm: 4.97490e+03
> Group rest with mass 4.50563e+04, Ekrot 8.32461e+18 Det(I) = -9.80133e+94
> COM: 4.35971 4.41628 4.51177
> P: 7263.49854 31947.81250 2073627.46875
> V: 0.16121 0.70906 46.02300
> J: 58359181927255633821696.00000 -2981584234146221780893696.00000
> 7570244118347216060416.00000
> w: 0.00002 -0.00001 -0.00026
> Inertia tensor (3x3):
> Inertia tensor[ 0]={ 1.16144e+37, -3.11896e+36, -1.40998e+38}
> Inertia tensor[ 1]={-3.11896e+36, 8.37572e+35, 3.78640e+37}
> Inertia tensor[ 2]={-1.40998e+38, 3.78640e+37, 1.71171e+39}
> Step Time Lambda
> 4411800 8823.60000 0.00000
> 
> Energies (kJ/mol)
> Bond Angle Proper Dih. Ryckaert-Bell. LJ-14
> nan 9.22750e+06 0.00000e+00 1.80795e+04 nan
> Coulomb-14 LJ (SR) Coulomb (SR) Potential Kinetic En.
> nan nan nan nan nan
> Total Energy Temperature Pressure (bar)
> nan nan 0.00000e+00
> 
> After that the energy, temperature are coming as "nan".
> 
> The ".mdp" file that I have used is this:
> integrator = md
> tinit = 0
> dt = 0.002
> nstcomm =1
> comm_mode = ANGULAR
> nsteps = 5000000
> nstxout = 100
> nstvout = 500
> nstfout = 0
> nstlog = 100
> nstenergy = 10
> nstxtcout = 100
> xtc-precision = 100
> energygrps = PROTEIN
> nstlist = 10
> ns_type = simple
> pbc = no
> rlist = 0.0
> domain-decomposition = no
> coulombtype = Cut-off
> rcoulomb-switch = 0
> rcoulomb = 0.0
> epsilon-r = 1
> ; Van der Waals =
> vdw-type = Cut-off
> rvdw-switch = 0
> rvdw = 0.0
> ; Berendsen temperature coupling is on in two groups
> Tcoupl = berendsen
> tc-grps = Protein
> tau_t = 0.1
> ref_t = 300
> ; Generate velocites is off at 300 K.
> gen_vel = yes
> gen_temp = 300
> gen_seed = 173529
> 
> What can be the reason for this?

Your system is blowing up.  Something about the model physics is failing and so 
everything is going haywire, and hence "nan" (not a number) values for energies, 
temperature, etc.

The force fields packaged with Gromacs were all designed for condensed-phase 
simulations.  There may be some precedence for using them in a vacuum 
simulation, so check the literature for recommended settings, shortcomings, etc.

-Justin

> Any suggestion regarding this will be very helpful.
> Thanks in advance
> Sarbani
> <http://sigads.rediff.com/RealMedia/ads/click_nx.ads/www.rediffmail.com/signatureline.htm@Middle?>
> 

-- 
========================================

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================