[gmx-users] Problem with incorrect GB-Polarization Energy Value
Justin Lemkul
jalemkul at vt.edu
Wed Aug 29 17:14:28 CEST 2012
On 8/29/12 11:11 AM, jesmin jahan wrote:
> Thanks Mark for your reply.
>
> For the time being, I admit your claim that I am comparing apple with orange.
> So, to investigate more, I run the simulation without any modification
> in parameter fields and force field I am using. My test data is CMV
> virus shell.
> I am using the following commands.
>
> pdb2gmx -f 1F15-full.pdb -ter -ignh -ff amber99 -water none
> grompp -f mdr.mdp -c conf.gro -p topol.top -o imd.tpr
> OMP_NUM_THREADS=12 mdrun -nt 16 -s imd.tpr
>
>
> The log file looks like this:
> :-) G R O M A C S (-:
>
> GROningen MAchine for Chemical Simulation
>
> :-) VERSION 4.6-dev-20120820-87e5bcf (-:
>
> Written by Emile Apol, Rossen Apostolov, Herman J.C. Berendsen,
> Aldert van Buuren, Pär Bjelkmar, Rudi van Drunen, Anton Feenstra,
> Gerrit Groenhof, Peter Kasson, Per Larsson, Pieter Meulenhoff,
> Teemu Murtola, Szilard Pall, Sander Pronk, Roland Schulz,
> Michael Shirts, Alfons Sijbers, Peter Tieleman,
>
> Berk Hess, David van der Spoel, and Erik Lindahl.
>
> Copyright (c) 1991-2000, University of Groningen, The Netherlands.
> Copyright (c) 2001-2010, The GROMACS development team at
> Uppsala University & The Royal Institute of Technology, Sweden.
> check out http://www.gromacs.org for more information.
>
> This program is free software; you can redistribute it and/or
> modify it under the terms of the GNU General Public License
> as published by the Free Software Foundation; either version 2
> of the License, or (at your option) any later version.
>
> :-) mdrun_mpi (-:
>
>
> ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
> B. Hess and C. Kutzner and D. van der Spoel and E. Lindahl
> GROMACS 4: Algorithms for highly efficient, load-balanced, and scalable
> molecular simulation
> J. Chem. Theory Comput. 4 (2008) pp. 435-447
> -------- -------- --- Thank You --- -------- --------
>
>
> ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
> D. van der Spoel, E. Lindahl, B. Hess, G. Groenhof, A. E. Mark and H. J. C.
> Berendsen
> GROMACS: Fast, Flexible and Free
> J. Comp. Chem. 26 (2005) pp. 1701-1719
> -------- -------- --- Thank You --- -------- --------
>
>
> ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
> E. Lindahl and B. Hess and D. van der Spoel
> GROMACS 3.0: A package for molecular simulation and trajectory analysis
> J. Mol. Mod. 7 (2001) pp. 306-317
> -------- -------- --- Thank You --- -------- --------
>
>
> ++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
> H. J. C. Berendsen, D. van der Spoel and R. van Drunen
> GROMACS: A message-passing parallel molecular dynamics implementation
> Comp. Phys. Comm. 91 (1995) pp. 43-56
> -------- -------- --- Thank You --- -------- --------
>
> Input Parameters:
> integrator = md
> nsteps = 0
> init-step = 0
> ns-type = Grid
> nstlist = 10
> ndelta = 2
> nstcomm = 10
> comm-mode = Linear
> nstlog = 1000
> nstxout = 0
> nstvout = 0
> nstfout = 0
> nstcalcenergy = 10
> nstenergy = 100
> nstxtcout = 0
> init-t = 0
> delta-t = 0.001
> xtcprec = 1000
> nkx = 0
> nky = 0
> nkz = 0
> pme-order = 4
> ewald-rtol = 1e-05
> ewald-geometry = 0
> epsilon-surface = 0
> optimize-fft = FALSE
> ePBC = no
> bPeriodicMols = FALSE
> bContinuation = FALSE
> bShakeSOR = FALSE
> etc = No
> bPrintNHChains = FALSE
> nsttcouple = -1
> epc = No
> epctype = Isotropic
> nstpcouple = -1
> tau-p = 1
> ref-p (3x3):
> ref-p[ 0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> ref-p[ 1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> ref-p[ 2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> compress (3x3):
> compress[ 0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> compress[ 1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> compress[ 2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> refcoord-scaling = No
> posres-com (3):
> posres-com[0]= 0.00000e+00
> posres-com[1]= 0.00000e+00
> posres-com[2]= 0.00000e+00
> posres-comB (3):
> posres-comB[0]= 0.00000e+00
> posres-comB[1]= 0.00000e+00
> posres-comB[2]= 0.00000e+00
> rlist = 1
> rlistlong = 1
> rtpi = 0.05
> coulombtype = Cut-off
> rcoulomb-switch = 0
> rcoulomb = 1
> vdwtype = Cut-off
> rvdw-switch = 0
> rvdw = 1
> epsilon-r = 1
> epsilon-rf = inf
> tabext = 1
> implicit-solvent = GBSA
> gb-algorithm = HCT
> gb-epsilon-solvent = 80
> nstgbradii = 1
> rgbradii = 1
> gb-saltconc = 0
> gb-obc-alpha = 1
> gb-obc-beta = 0.8
> gb-obc-gamma = 4.85
> gb-dielectric-offset = 0.009
> sa-algorithm = None
> sa-surface-tension = 2.25936
> DispCorr = No
> bSimTemp = FALSE
> free-energy = no
> nwall = 0
> wall-type = 9-3
> wall-atomtype[0] = -1
> wall-atomtype[1] = -1
> wall-density[0] = 0
> wall-density[1] = 0
> wall-ewald-zfac = 3
> pull = no
> rotation = FALSE
> disre = No
> disre-weighting = Conservative
> disre-mixed = FALSE
> dr-fc = 1000
> dr-tau = 0
> nstdisreout = 100
> orires-fc = 0
> orires-tau = 0
> nstorireout = 100
> dihre-fc = 0
> em-stepsize = 0.01
> em-tol = 10
> niter = 20
> fc-stepsize = 0
> nstcgsteep = 1000
> nbfgscorr = 10
> ConstAlg = Lincs
> shake-tol = 0.0001
> lincs-order = 4
> lincs-warnangle = 30
> lincs-iter = 1
> bd-fric = 0
> ld-seed = 1993
> cos-accel = 0
> deform (3x3):
> deform[ 0]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> deform[ 1]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> deform[ 2]={ 0.00000e+00, 0.00000e+00, 0.00000e+00}
> adress = FALSE
> userint1 = 0
> userint2 = 0
> userint3 = 0
> userint4 = 0
> userreal1 = 0
> userreal2 = 0
> userreal3 = 0
> userreal4 = 0
> grpopts:
> nrdf: 9534
> ref-t: 0
> tau-t: 0
> anneal: No
> ann-npoints: 0
> acc: 0 0 0
> nfreeze: N N N
> energygrp-flags[ 0]: 0
> efield-x:
> n = 0
> efield-xt:
> n = 0
> efield-y:
> n = 0
> efield-yt:
> n = 0
> efield-z:
> n = 0
> efield-zt:
> n = 0
> bQMMM = FALSE
> QMconstraints = 0
> QMMMscheme = 0
> scalefactor = 1
> qm-opts:
> ngQM = 0
>
> Initializing Domain Decomposition on 16 nodes
> Dynamic load balancing: auto
> Will sort the charge groups at every domain (re)decomposition
> Minimum cell size due to bonded interactions: 0.000 nm
> Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25
> Optimizing the DD grid for 16 cells with a minimum initial size of 0.000 nm
> Domain decomposition grid 4 x 4 x 1, separate PME nodes 0
> Domain decomposition nodeid 0, coordinates 0 0 0
>
> Detecting CPU-specific acceleration. Present hardware specification:
> Vendor: GenuineIntel
> Brand: Intel(R) Xeon(R) CPU X5680 @ 3.33GHz
> Family: 6 Model: 44 Stepping: 2
> Features: htt sse2 sse4.1 aes rdtscp
> Acceleration most likely to fit this hardware: SSE4.1
> Acceleration selected at Gromacs compile time: SSE4.1
>
> Table routines are used for coulomb: FALSE
> Table routines are used for vdw: FALSE
> Cut-off's: NS: 1 Coulomb: 1 LJ: 1
> System total charge: 6.000
> Configuring nonbonded kernels...
> Configuring standard C nonbonded kernels...
>
>
>
> Linking all bonded interactions to atoms
>
> The initial number of communication pulses is: X 2 Y 2
> The initial domain decomposition cell size is: X 0.79 nm Y 0.89 nm
>
> The maximum allowed distance for charge groups involved in interactions is:
> non-bonded interactions 1.000 nm
> (the following are initial values, they could change due to box deformation)
> two-body bonded interactions (-rdd) 1.000 nm
> multi-body bonded interactions (-rdd) 0.794 nm
>
> When dynamic load balancing gets turned on, these settings will change to:
> The maximum number of communication pulses is: X 2 Y 2
> The minimum size for domain decomposition cells is 0.500 nm
> The requested allowed shrink of DD cells (option -dds) is: 0.80
> The allowed shrink of domain decomposition cells is: X 0.63 Y 0.56
> The maximum allowed distance for charge groups involved in interactions is:
> non-bonded interactions 1.000 nm
> two-body bonded interactions (-rdd) 1.000 nm
> multi-body bonded interactions (-rdd) 0.500 nm
>
>
> Making 2D domain decomposition grid 4 x 4 x 1, home cell index 0 0 0
>
> Center of mass motion removal mode is Linear
> We have the following groups for center of mass motion removal:
> 0: rest
> There are: 3179 Atoms
> Charge group distribution at step 0: 84 180 252 196 237 210 255 157
> 254 197 266 176 186 104 224 201
> Grid: 4 x 4 x 4 cells
> Initial temperature: 0 K
>
> Started mdrun on node 0 Wed Aug 29 02:32:21 2012
>
> Step Time Lambda
> 0 0.00000 0.00000
>
> Energies (kJ/mol)
> GB Polarization LJ (SR) Coulomb (SR) Potential Kinetic En.
> -1.65116e+04 5.74908e+08 -2.37699e+05 5.74654e+08 6.36009e+11
> Total Energy Temperature Pressure (bar)
> 6.36584e+11 1.60465e+10 0.00000e+00
>
> <====== ############### ==>
> <==== A V E R A G E S ====>
> <== ############### ======>
>
> Statistics over 1 steps using 1 frames
>
> Energies (kJ/mol)
> GB Polarization LJ (SR) Coulomb (SR) Potential Kinetic En.
> -1.65116e+04 5.74908e+08 -2.37699e+05 5.74654e+08 6.36009e+11
> Total Energy Temperature Pressure (bar)
> 6.36584e+11 1.60465e+10 0.00000e+00
>
> Total Virial (kJ/mol)
> -1.13687e+09 1.14300e+07 -1.23884e+07
> 1.14273e+07 -1.15125e+09 -5.31658e+06
> -1.23830e+07 -5.31326e+06 -1.16512e+09
>
> Pressure (bar)
> 0.00000e+00 0.00000e+00 0.00000e+00
> 0.00000e+00 0.00000e+00 0.00000e+00
> 0.00000e+00 0.00000e+00 0.00000e+00
>
> Total Dipole (D)
> 1.35524e+03 -4.39059e+01 2.16985e+03
>
>
> M E G A - F L O P S A C C O U N T I N G
>
> RF=Reaction-Field FE=Free Energy SCFE=Soft-Core/Free Energy
> T=Tabulated W3=SPC/TIP3p W4=TIP4p (single or pairs)
> NF=No Forces
>
> Computing: M-Number M-Flops % Flops
> -----------------------------------------------------------------------------
> Generalized Born Coulomb 0.006162 0.296 0.2
> GB Coulomb + LJ 0.446368 27.228 19.8
> Outer nonbonded loop 0.015554 0.156 0.1
> Born radii (HCT/OBC) 0.452530 82.813 60.3
> Born force chain rule 0.452530 6.788 4.9
> NS-Pairs 0.940291 19.746 14.4
> Reset In Box 0.003179 0.010 0.0
> CG-CoM 0.006358 0.019 0.0
> Virial 0.003899 0.070 0.1
> Stop-CM 0.006358 0.064 0.0
> Calc-Ekin 0.006358 0.172 0.1
> -----------------------------------------------------------------------------
> Total 137.361 100.0
> -----------------------------------------------------------------------------
>
>
> D O M A I N D E C O M P O S I T I O N S T A T I S T I C S
>
> av. #atoms communicated per step for force: 2 x 7369.0
>
>
> R E A L C Y C L E A N D T I M E A C C O U N T I N G
>
> Computing: Nodes Number G-Cycles Seconds %
> -----------------------------------------------------------------------
> Domain decomp. 16 1 0.210 0.1 11.4
> Comm. coord. 16 1 0.006 0.0 0.3
> Neighbor search 16 1 0.118 0.1 6.4
> Force 16 1 1.319 0.8 71.4
> Wait + Comm. F 16 1 0.016 0.0 0.9
> Update 16 1 0.003 0.0 0.2
> Comm. energies 16 1 0.093 0.1 5.0
> Rest 16 0.082 0.1 4.4
> -----------------------------------------------------------------------
> Total 16 1.847 1.1 100.0
> -----------------------------------------------------------------------
>
> NOTE: 5 % of the run time was spent communicating energies,
> you might want to use the -gcom option of mdrun
>
>
> Parallel run - timing based on wallclock.
>
> NODE (s) Real (s) (%)
> Time: 0.036 0.036 100.0
> (Mnbf/s) (GFlops) (ns/day) (hour/ns)
> Performance: 12.702 3.856 2.425 9.896
> Finished mdrun on node 0 Wed Aug 29 02:32:21 2012
>
>
>
> The GB- energy value reported is half of that reported by Amber 11 and
> Octree based Molecular dynamic package.
>
> Although I guess the difference can be due to the difference in
> algorithms they are using, but there could be some other reason.
> If anyone knows what are the possible reasons behind this, please let
> me know. May be fixing them will give me same value for all different
> Molecular Dynamic Package.
>
I wouldn't trust the result you're getting here - the energy values and
temperature (10^10, yikes!) suggest there is something very wrong with the
starting configuration.
-Justin
--
========================================
Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
More information about the gromacs.org_gmx-users
mailing list