[gmx-users] strange vdw energy from rerun with GBSA model.
Justin A. Lemkul
jalemkul at vt.edu
Thu Aug 11 15:03:51 CEST 2011
Da-Wei Li wrote:
> Dear Justin
>
> An implicit water simulaiton with this short cutoff is problematic but I
> think it is fine for rerun. I want to exactly repeat the original
> energies in the explicit water MD.
>
> The manu say "neighbor list searching will be performed for every frame"
> with rerun option. So that I don't think this is the cause.
>
Right, forgot about that. I still think the cutoffs are a problem, though.
-Justin
> best,
>
> dawei
>
> On Thu, Aug 11, 2011 at 8:47 AM, Justin A. Lemkul <jalemkul at vt.edu
> <mailto:jalemkul at vt.edu>> wrote:
>
>
>
> Da-Wei Li wrote:
>
> Dear Mark and others
>
> I did more tests and thought that it might come from numerical
> error. The reasons are
>
> 1. If I use .trr file instead of the low precision xtc file,
> things become better, ie, I get much less snapshots that has
> high energy.
>
> 2. I supplied -pforce in my mdrun -rerun and found that the high
> vdw energy was usually caused by one pair of atoms, whose
> distance was just very near the clash zone, so that small error
> on the coordinates would cause large energy error. The force is
> always around 10000.
>
> 3. Actually bond length and bond angle energies are also
> affected. I can fully reproduce these two energies if I use .trr
> file in my rerun but will get several tens of kj/mol error if I
> use .xtc file, for a protein of size of 100 AA.
>
>
> Now the question I still have is whether numerical error can be
> so large? The xtc file has a precision of 0.001 nm. Anyway, I
> will test more by using double precision Gromacs and define
> energy groups so that I can compare energy of protein directly
> between original MD and rerun.
>
>
>
> To Mark only
>
> Thanks. Here it is my script for
>
> rerun: mdrun -v -s pbsa.tpr -rerun coor.xtc -e rerun
> superpose: trjconv -s em.tpr -f coor.xtc -o nojump.xtc -pbc
> nojump (em.tpr is generated for energy minimization, protein is
> in the middle of the box)
>
> rerun .mdp file:
>
> ******************************__****************************
>
> ; Run parameters
> integrator = md ; leap-frog integrator
> nsteps = 50000000 ; 100 ns
> dt = 0.002 ; 2 fs
> ; Output control
> nstxout = 500000 ; save coordinates every 1000 ps
> nstvout = 500000 ; save velocities every 1000 ps
> nstxtcout = 5000 ; xtc compressed trajectory output
> every 1 ps
> nstenergy = 5000 ; save energies every 1 ps
> nstlog = 5000 ; update log file every 1 ps
> xtc_grps = Protein ; save protein part only
> ; Bond parameters
> continuation = yes ; Restarting after NPT
> constraint_algorithm = lincs ; holonomic constraints
> constraints = hbonds ; all bonds (even heavy atom-H bonds)
> constrained
> lincs_iter = 1 ; accuracy of LINCS
> lincs_order = 4 ; also related to accuracy
> ; Neighborsearching
> ns_type = grid ; search neighboring grid cels
> nstlist = 10 ; 20 fs
> rlist = 0.8 ; short-range neighborlist cutoff (in nm)
> rcoulomb = 0.8 ; short-range electrostatic cutoff (in nm)
> rvdw = 1.0 ; short-range van der Waals cutoff (in nm)
> ; Electrostatics
> coulombtype = cut-off ; Particle Mesh Ewald for long-range
> electrostatics
> pme_order = 4 ; cubic interpolation
> fourierspacing = 0.12 ; grid spacing for FFT
> ; Temperature coupling is on
> tcoupl = no ; modified Berendsen thermostat
> tc-grps = System ; two coupling groups - more accurate
> tau_t = 0.1 ; time constant, in ps
> ref_t = 300 ; reference temperature, one for each
> group, in K
> ; Pressure coupling is on
> pcoupl = no ; Pressure coupling on in NPT
> pcoupltype = isotropic ; uniform scaling of box vectors
> tau_p = 2.0 ; time constant, in ps
> ref_p = 1.0 ; reference pressure, in bar
> compressibility = 4.5e-5 ; isothermal compressibility of
> water, bar^-1
> ; Periodic boundary conditions
> pbc = no ; 3-D PBC
> ; Dispersion correction
> ;DispCorr = EnerPres ; account for cut-off vdW scheme
> DispCorr = no
> ; Velocity generation
> gen_vel = no ; Velocity generation is off
>
>
>
> ; IMPLICIT SOLVENT ALGORITHM
> implicit_solvent = GBSA
> gb_algorithm = OBC
> nstgbradii = 1
> rgbradii = 0.8
> gb_epsilon_solvent = 80
> gb_saltconc = 0
> gb_obc_alpha = 1
> gb_obc_beta = 0.8
> gb_obc_gamma = 4.85
> gb_dielectric_offset = 0.009
> sa_algorithm = Ace-approximation
> sa_surface_tension = 2.25936
>
> ******************************__******************************__***************************
>
> Thanks all.
>
>
> Using cutoffs this small may be the source of your problem. Proper
> implicit solvent calculations require longer cutoffs than would
> normally be used in explicit solvent MD. Try with longer (2.0 nm)
> or infinite cutoffs and a fixed neighbor list (nstlist = 0) and see
> if that smooths out the problem. What's likely happening now is
> that you've got interactions moving very quickly in and out of the
> very short cutoff, causing spikes in energy in between neighbor list
> updates.
>
> -Justin
>
> --
> ==============================__==========
>
> Justin A. Lemkul
> Ph.D. Candidate
> ICTAS Doctoral Scholar
> MILES-IGERT Trainee
> Department of Biochemistry
> Virginia Tech
> Blacksburg, VA
> jalemkul[at]vt.edu <http://vt.edu> | (540) 231-9080
> <tel:%28540%29%20231-9080>
> http://www.bevanlab.biochem.__vt.edu/Pages/Personal/justin
> <http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin>
>
> ==============================__==========
>
> --
> gmx-users mailing list gmx-users at gromacs.org
> <mailto:gmx-users at gromacs.org>
> http://lists.gromacs.org/__mailman/listinfo/gmx-users
> <http://lists.gromacs.org/mailman/listinfo/gmx-users>
> Please search the archive at
> http://www.gromacs.org/__Support/Mailing_Lists/Search
> <http://www.gromacs.org/Support/Mailing_Lists/Search> before posting!
> Please don't post (un)subscribe requests to the list. Use the www
> interface or send it to gmx-users-request at gromacs.org
> <mailto:gmx-users-request at gromacs.org>.
> Can't post? Read http://www.gromacs.org/__Support/Mailing_Lists
> <http://www.gromacs.org/Support/Mailing_Lists>
>
>
--
========================================
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
========================================
More information about the gromacs.org_gmx-users
mailing list