Re: [gmx-users] potential energy NaN and strange dependence on cut-offs

[David van der Spoel]

Manik Mayur wrote:


2009/9/8 Mark Abraham >


Jennifer Williams wrote:

Hi users,

I am running a very simple simulation of methane inside a pore
(v.much like a carbon nanotube but in my case the tube is
supposed to represent silica.) I keep this tube frozen.

I start with an energy minimisation-however this runs to
completion almost instantly and I keep get NaN for my potential
energy:

Steepest Descents converged to machine precision in 18 steps,
but did not reach the requested Fmax < 10.
Potential Energy  =nan
Maximum force =  6.5738518e+01 on atom 2133
Norm of force =  1.5461593e+00


This nan suggests some kind of severe atomic overlap. Reconsider
your coordinates and the box size implied by your coordinate file.

Mark


Otherwise the trajectory looks OK (methane moving around inside
the cylinder). If I go on to use the conf.gro file for an mdrun,
it runs to completion and generates what looks like a reasonable
trajectory, however the output again contains NaN i.e:

  Energies (kJ/mol)
   LJ (SR)   Coulomb (SR)  PotentialKinetic En.  
Total Energy
   nan0.0e+00nan3.36749e+01
   nan

 Conserved En.Temperature Pressure (bar)
   nan3.00010e+02nan

and calculating the Diffusion coefficient gives:
D[   CH4] 613.6682 (+/- 97.0563) 1e-5 cm^2/s

If I do the same calculation but reduce the cut-offs to 0.9. I get

  Energies (kJ/mol)
   LJ (SR)   Coulomb (SR)  PotentialKinetic En.  
Total Energy
   nan0.0e+00nan3.36750e+01
   nan

 Conserved En.Temperature Pressure (bar)
   nan3.00011e+02nan

D[   CH4] 237.8712 (+/- 53.5975) 1e-5 cm^2/s

And for a cut-off of 1.3nm I get

  Energies (kJ/mol)
   LJ (SR)   Coulomb (SR)  PotentialKinetic En.  
Total Energ

y
   nan0.0e+00nan3.36737e+01
   na

n
 Conserved En.Temperature Pressure (bar)
   nan2.9e+02nan


D[   CH4] 19.7953 (+/- 154.0168) 1e-5 cm^2/s


For this system, the cut-off shouldn?t need to be larger than
0.8 (I have plotted graphs of calculated V vs r) so it is
worrying that the diffusion coefficient is showing such
dependence on the cut-offs when they should all give the same
result.

Can anyone offer any insight into this? I?ve tried changing the
timestep making it both larger and smaller and many other
things. I?ve pasted the relevant parts of my files below:

I?m using gromacs 4.0.5 ?at the moment running in serial.

Thanks for any advice,

Top file

[ defaults ]
; nbfunccomb-rulegen-pairsfudgeLJfudgeQQ
12yes1.0   1.0
;
;
[ atomtypes ]
;   typemasschargeptype   c6c12
   OSM15.99940.00 A 0.2708   1.538176



Check your tpr file: are the C6 and C12 reasonable?
If you get NaN something is wrong.


;
; Include forcefield parameters
#include "CH4.itp"
;
;
[ moleculetype ]
;Namenrexcl
MCM3
[ atoms ]
;nrtyperesnrresidueatomcgnrcharge  
 mass

1   OSM 1   MCM OSM 1   0   15.9994
2   OSM 1   MCM OSM 2   0   15.9994
.etc
2127OSM 1   MCM OSM 21270   15.9994
2128OSM 1   MCM OSM 21280   15.9994


[ system ]
; Name
CH4 in MCM

[ molecules ]
; Compound#mols
MCM1
CH410

CH4.itp file

[ atomtypes ]
;   type  masschargeptype   c6c12
   CH416.043 0.00 A0.37321.24650457
;
[ moleculetype ]
; name  nrexcl
CH42

[ atoms ]
;   nr  typeresnr   residu  atomcgnrchargemass
1   CH4  1   CH4 CH4 10.00  16.043



.mdp file

;
;File 'mdout.mdp' was generated
;By user: jwillia4 (353773)
;On host: vlxhead2
;At date: Fri Jun 26 15:47:37 2009
;
; VARIOUS PREPROCESSING OPTIONS
; Preprocessor information: use cpp syntax.

Re: [gmx-users] potential energy NaN and strange dependence on cut-offs

[Manik Mayur]

2009/9/8 Mark Abraham 

> Jennifer Williams wrote:
>
>> Hi users,
>>
>> I am running a very simple simulation of methane inside a pore (v.much
>> like a carbon nanotube but in my case the tube is supposed to represent
>> silica.) I keep this tube frozen.
>>
>> I start with an energy minimisation-however this runs to completion almost
>> instantly and I keep get NaN for my potential energy:
>>
>> Steepest Descents converged to machine precision in 18 steps,
>> but did not reach the requested Fmax < 10.
>> Potential Energy  =nan
>> Maximum force =  6.5738518e+01 on atom 2133
>> Norm of force =  1.5461593e+00
>>
>
> This nan suggests some kind of severe atomic overlap. Reconsider your
> coordinates and the box size implied by your coordinate file.
>
> Mark
>
>
>  Otherwise the trajectory looks OK (methane moving around inside the
>> cylinder). If I go on to use the conf.gro file for an mdrun, it runs to
>> completion and generates what looks like a reasonable trajectory, however
>> the output again contains NaN i.e:
>>
>>   Energies (kJ/mol)
>>LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energy
>>nan0.0e+00nan3.36749e+01nan
>>  Conserved En.Temperature Pressure (bar)
>>nan3.00010e+02nan
>>
>> and calculating the Diffusion coefficient gives:
>> D[   CH4] 613.6682 (+/- 97.0563) 1e-5 cm^2/s
>>
>> If I do the same calculation but reduce the cut-offs to 0.9. I get
>>
>>   Energies (kJ/mol)
>>LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energy
>>nan0.0e+00nan3.36750e+01nan
>>  Conserved En.Temperature Pressure (bar)
>>nan3.00011e+02nan
>>
>> D[   CH4] 237.8712 (+/- 53.5975) 1e-5 cm^2/s
>>
>> And for a cut-off of 1.3nm I get
>>
>>   Energies (kJ/mol)
>>LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energ
>> y
>>nan0.0e+00nan3.36737e+01na
>> n
>>  Conserved En.Temperature Pressure (bar)
>>nan2.9e+02nan
>>
>>
>> D[   CH4] 19.7953 (+/- 154.0168) 1e-5 cm^2/s
>>
>>
>> For this system, the cut-off shouldn?t need to be larger than 0.8 (I have
>> plotted graphs of calculated V vs r) so it is worrying that the diffusion
>> coefficient is showing such dependence on the cut-offs when they should all
>> give the same result.
>>
>> Can anyone offer any insight into this? I?ve tried changing the timestep
>> making it both larger and smaller and many other things. I?ve pasted the
>> relevant parts of my files below:
>>
>> I?m using gromacs 4.0.5 ?at the moment running in serial.
>>
>> Thanks for any advice,
>>
>> Top file
>>
>> [ defaults ]
>> ; nbfunccomb-rulegen-pairsfudgeLJfudgeQQ
>> 12yes1.0   1.0
>> ;
>> ;
>> [ atomtypes ]
>> ;   typemasschargeptype   c6c12
>>OSM15.99940.00 A 0.2708   1.538176
>>
>> ;
>> ; Include forcefield parameters
>> #include "CH4.itp"
>> ;
>> ;
>> [ moleculetype ]
>> ;Namenrexcl
>> MCM3
>> [ atoms ]
>> ;nrtyperesnrresidueatomcgnrcharge
>>  mass
>> 1   OSM 1   MCM OSM 1   0   15.9994
>> 2   OSM 1   MCM OSM 2   0   15.9994
>> ..etc
>> 2127OSM 1   MCM OSM 21270   15.9994
>> 2128OSM 1   MCM OSM 21280   15.9994
>>
>>
>> [ system ]
>> ; Name
>> CH4 in MCM
>>
>> [ molecules ]
>> ; Compound#mols
>> MCM1
>> CH410
>>
>> CH4.itp file
>>
>> [ atomtypes ]
>> ;   type  masschargeptype   c6c12
>>CH416.043 0.00 A0.37321.24650457
>> ;
>> [ moleculetype ]
>> ; name  nrexcl
>> CH42
>>
>> [ atoms ]
>> ;   nr  typeresnr   residu  atomcgnrchargemass
>> 1   CH4  1   CH4 CH4 10.00  16.043
>>
>>
>>
>> .mdp file
>>
>> ;
>> ;File 'mdout.mdp' was generated
>> ;By user: jwillia4 (353773)
>> ;On host: vlxhead2
>> ;At date: Fri Jun 26 15:47:37 2009
>> ;
>> ; VARIOUS PREPROCESSING OPTIONS
>> ; Preprocessor information: use cpp syntax.
>> ; e.g.: -I/home/joe/doe -I/home/mary/hoe
>> include  = -I../top
>> ; e.g.: -DI_Want_Cookies -DMe_Too
>> define   =
>>
>> ; RUN CONTROL PARAMETERS
>> integrator   = steep
>> ; Start time and timestep in ps
>> tinit= 0
>> dt   = 0.0001
>>
>
dt = 0.1 fs are you sure??


>  nsteps   = 10
>> ; For exact run continuation or redoing part of a run
>> ; Part index is updated automatically on checkpointing (keeps files
>> separate)
>> simulation_part  = 1
>> init_step= 0
>> ; mode for center of mass motion removal
>> comm-mode

Re: [gmx-users] potential energy NaN and strange dependence on cut-offs

[Mark Abraham]

Jennifer Williams wrote:

Hi users,

I am running a very simple simulation of methane inside a pore (v.much 
like a carbon nanotube but in my case the tube is supposed to represent 
silica.) I keep this tube frozen.


I start with an energy minimisation-however this runs to completion 
almost instantly and I keep get NaN for my potential energy:


Steepest Descents converged to machine precision in 18 steps,
but did not reach the requested Fmax < 10.
Potential Energy  =nan
Maximum force =  6.5738518e+01 on atom 2133
Norm of force =  1.5461593e+00


This nan suggests some kind of severe atomic overlap. Reconsider your 
coordinates and the box size implied by your coordinate file.


Mark

Otherwise the trajectory looks OK (methane moving around inside the 
cylinder). If I go on to use the conf.gro file for an mdrun, it runs to 
completion and generates what looks like a reasonable trajectory, 
however the output again contains NaN i.e:


   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energy
nan0.0e+00nan3.36749e+01nan
  Conserved En.Temperature Pressure (bar)
nan3.00010e+02nan

and calculating the Diffusion coefficient gives:
D[   CH4] 613.6682 (+/- 97.0563) 1e-5 cm^2/s

If I do the same calculation but reduce the cut-offs to 0.9. I get

   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energy
nan0.0e+00nan3.36750e+01nan
  Conserved En.Temperature Pressure (bar)
nan3.00011e+02nan

D[   CH4] 237.8712 (+/- 53.5975) 1e-5 cm^2/s

And for a cut-off of 1.3nm I get

   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energ
y
nan0.0e+00nan3.36737e+01na
n
  Conserved En.Temperature Pressure (bar)
nan2.9e+02nan


D[   CH4] 19.7953 (+/- 154.0168) 1e-5 cm^2/s


For this system, the cut-off shouldn?t need to be larger than 0.8 (I 
have plotted graphs of calculated V vs r) so it is worrying that the 
diffusion coefficient is showing such dependence on the cut-offs when 
they should all give the same result.


Can anyone offer any insight into this? I?ve tried changing the timestep 
making it both larger and smaller and many other things. I?ve pasted the 
relevant parts of my files below:


I?m using gromacs 4.0.5 ?at the moment running in serial.

Thanks for any advice,

Top file

[ defaults ]
; nbfunccomb-rulegen-pairsfudgeLJfudgeQQ
12yes1.0   1.0
;
;
[ atomtypes ]
;   typemasschargeptype   c6c12
OSM15.99940.00 A 0.2708   1.538176

;
; Include forcefield parameters
#include "CH4.itp"
;
;
[ moleculetype ]
;Namenrexcl
MCM3
[ atoms ]
;nrtyperesnrresidueatomcgnrcharge
mass

1   OSM 1   MCM OSM 1   0   15.9994
2   OSM 1   MCM OSM 2   0   15.9994
..etc
2127OSM 1   MCM OSM 21270   15.9994
2128OSM 1   MCM OSM 21280   15.9994


[ system ]
; Name
CH4 in MCM

[ molecules ]
; Compound#mols
MCM1
CH410

CH4.itp file

[ atomtypes ]
;   type  masschargeptype   c6c12
CH416.043 0.00 A0.37321.24650457
;
[ moleculetype ]
; name  nrexcl
CH42

[ atoms ]
;   nr  typeresnr   residu  atomcgnrchargemass
1   CH4  1   CH4 CH4 10.00  16.043



.mdp file

;
;File 'mdout.mdp' was generated
;By user: jwillia4 (353773)
;On host: vlxhead2
;At date: Fri Jun 26 15:47:37 2009
;
; VARIOUS PREPROCESSING OPTIONS
; Preprocessor information: use cpp syntax.
; e.g.: -I/home/joe/doe -I/home/mary/hoe
include  = -I../top
; e.g.: -DI_Want_Cookies -DMe_Too
define   =

; RUN CONTROL PARAMETERS
integrator   = steep
; Start time and timestep in ps
tinit= 0
dt   = 0.0001
nsteps   = 10
; For exact run continuation or redoing part of a run
; Part index is updated automatically on checkpointing (keeps files 
separate)

simulation_part  = 1
init_step= 0
; mode for center of mass motion removal
comm-mode= linear
; number of steps for center of mass motion removal
nstcomm  = 1
; group(s) for center of mass motion removal
comm-grps=

; LANGEVIN DYNAMICS OPTIONS
; Friction coefficient (amu/ps) and random seed
bd-fric  = 0
ld-seed  = 1993

; ENERGY MINIMIZATION OPTIONS
; Force tolerance and initial step-size
emtol=
emstep   = 0

[gmx-users] potential energy NaN and strange dependence on cut-offs

[Jennifer Williams]

Hi users,

I am running a very simple simulation of methane inside a pore (v.much  
like a carbon nanotube but in my case the tube is supposed to  
represent silica.) I keep this tube frozen.


I start with an energy minimisation-however this runs to completion  
almost instantly and I keep get NaN for my potential energy:


Steepest Descents converged to machine precision in 18 steps,
but did not reach the requested Fmax < 10.
Potential Energy  =nan
Maximum force =  6.5738518e+01 on atom 2133
Norm of force =  1.5461593e+00

Otherwise the trajectory looks OK (methane moving around inside the  
cylinder). If I go on to use the conf.gro file for an mdrun, it runs  
to completion and generates what looks like a reasonable trajectory,  
however the output again contains NaN i.e:


   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energy
nan0.0e+00nan3.36749e+01nan
  Conserved En.Temperature Pressure (bar)
nan3.00010e+02nan

and calculating the Diffusion coefficient gives:
D[   CH4] 613.6682 (+/- 97.0563) 1e-5 cm^2/s

If I do the same calculation but reduce the cut-offs to 0.9. I get

   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energy
nan0.0e+00nan3.36750e+01nan
  Conserved En.Temperature Pressure (bar)
nan3.00011e+02nan

D[   CH4] 237.8712 (+/- 53.5975) 1e-5 cm^2/s

And for a cut-off of 1.3nm I get

   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energ
y
nan0.0e+00nan3.36737e+01na
n
  Conserved En.Temperature Pressure (bar)
nan2.9e+02nan


D[   CH4] 19.7953 (+/- 154.0168) 1e-5 cm^2/s


For this system, the cut-off shouldn?t need to be larger than 0.8 (I  
have plotted graphs of calculated V vs r) so it is worrying that the  
diffusion coefficient is showing such dependence on the cut-offs when  
they should all give the same result.


Can anyone offer any insight into this? I?ve tried changing the  
timestep making it both larger and smaller and many other things. I?ve  
pasted the relevant parts of my files below:


I?m using gromacs 4.0.5 ?at the moment running in serial.

Thanks for any advice,

Top file

[ defaults ]
; nbfunccomb-rule   gen-pairs   fudgeLJ fudgeQQ
1   2   yes 1.01.0
;
;
[ atomtypes ]
;   typemasschargeptype   c6c12
OSM15.99940.00 A 0.2708   1.538176

;
; Include forcefield parameters
#include "CH4.itp"
;
;
[ moleculetype ]
;   Namenrexcl
MCM 3
[ atoms ]
;   nr  typeresnr   residue atomcgnrcharge  mass
1   OSM 1   MCM OSM 1   0   15.9994
2   OSM 1   MCM OSM 2   0   15.9994
..etc
2127OSM 1   MCM OSM 21270   15.9994
2128OSM 1   MCM OSM 21280   15.9994


[ system ]
; Name
CH4 in MCM

[ molecules ]
; Compound#mols
MCM1
CH410

CH4.itp file

[ atomtypes ]
;   type  masschargeptype   c6c12
CH416.043 0.00 A0.37321.24650457
;
[ moleculetype ]
; name  nrexcl
CH42

[ atoms ]
;   nr  typeresnr   residu  atomcgnrcharge  mass
1   CH4  1   CH4 CH4 10.00  16.043



.mdp file

;
;   File 'mdout.mdp' was generated
;   By user: jwillia4 (353773)
;   On host: vlxhead2
;   At date: Fri Jun 26 15:47:37 2009
;
; VARIOUS PREPROCESSING OPTIONS
; Preprocessor information: use cpp syntax.
; e.g.: -I/home/joe/doe -I/home/mary/hoe
include  = -I../top
; e.g.: -DI_Want_Cookies -DMe_Too
define   =

; RUN CONTROL PARAMETERS
integrator   = steep
; Start time and timestep in ps
tinit= 0
dt   = 0.0001
nsteps   = 10
; For exact run continuation or redoing part of a run
; Part index is updated automatically on checkpointing (keeps files separate)
simulation_part  = 1
init_step= 0
; mode for center of mass motion removal
comm-mode= linear
; number of steps for center of mass motion removal
nstcomm  = 1
; group(s) for center of mass motion removal
comm-grps=

; LANGEVIN DYNAMICS OPTIONS
; Friction coefficient (amu/ps) and random seed
bd-fric  = 0
ld-seed  = 1993

; ENERGY MINIMIZATION OPTIONS
; Force tolerance and initial step-size
emtol=
emstep   = 0.001
; Max number of iterations in relax_shells
niter=
; Step size (ps^2) for minimization of flexibl