[gmx-users] Regarding value of Elj, Eqq, Clj and Cqq values using LIE Approach
Justin Lemkul
jalemkul at vt.edu
Thu Jan 7 14:29:40 CET 2016
On 1/7/16 5:54 AM, shagun krishna wrote:
> Hii Justin..
>
> Since I am using PME option in my both simulation (Bound as well as alone),
> I think that it is necessary to use rerun option brfore parsing the values
> to g_lie. When I am running mdrun with new .mdp settings without PME option
> I am getting series of warnings and therefore program terminates.
>
> This is a look of my terminal:
>
> cbb at cbb-Precision-T1700:~/Desktop/md_ligand_RJC02836$ grompp -f
> md_SOL_rerun.mdp -c equil.gro -o md_sol_rerun.tpr
> :-) G R O M A C S (-:
>
> Gromacs Runs On Most of All Computer Systems
>
> :-) VERSION 4.6.5 (-:
>
> Contributions from Mark Abraham, Emile Apol, Rossen Apostolov,
> Herman J.C. Berendsen, Aldert van Buuren, Pär Bjelkmar,
> Rudi van Drunen, Anton Feenstra, Gerrit Groenhof, Christoph Junghans,
> Peter Kasson, Carsten Kutzner, 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-2012,2013, 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 Lesser General Public License
> as published by the Free Software Foundation; either version 2.1
> of the License, or (at your option) any later version.
>
> :-) grompp (-:
>
> Option Filename Type Description
> ------------------------------------------------------------
> -f md_SOL_rerun.mdp Input grompp input file with MD parameters
> -po mdout.mdp Output grompp input file with MD parameters
> -c equil.gro Input Structure file: gro g96 pdb tpr etc.
> -r conf.gro Input, Opt. Structure file: gro g96 pdb tpr etc.
> -rb conf.gro Input, Opt. Structure file: gro g96 pdb tpr etc.
> -n index.ndx Input, Opt. Index file
> -p topol.top Input Topology file
> -pp processed.top Output, Opt. Topology file
> -o md_sol_rerun.tpr Output Run input file: tpr tpb tpa
> -t traj.trr Input, Opt. Full precision trajectory: trr trj cpt
> -e ener.edr Input, Opt. Energy file
> -ref rotref.trr In/Out, Opt. Full precision trajectory: trr trj cpt
>
> Option Type Value Description
> ------------------------------------------------------
> -[no]h bool no Print help info and quit
> -[no]version bool no Print version info and quit
> -nice int 0 Set the nicelevel
> -[no]v bool no Be loud and noisy
> -time real -1 Take frame at or first after this time.
> -[no]rmvsbds bool yes Remove constant bonded interactions with virtual
> sites
> -maxwarn int 0 Number of allowed warnings during input
> processing. Not for normal use and may generate
> unstable systems
> -[no]zero bool no Set parameters for bonded interactions without
> defaults to zero instead of generating an error
> -[no]renum bool yes Renumber atomtypes and minimize number of
> atomtypes
>
> Ignoring obsolete mdp entry 'title'
>
> Back Off! I just backed up mdout.mdp to ./#mdout.mdp.9#
>
> WARNING 1 [file md_SOL_rerun.mdp, line 58]:
> Unknown left-hand '*coulombtype' in parameter file
>
Well, *coulombtype is not a keyword, but coulombtype is. Looks like you've
copied and pasted someone's attempt to emphasize a few lines.
Note that the settings below are for GROMOS96, and generally the point of
re-running without PME is to somehow still capture long-range electrostatics so
longer cutoffs are normally employed. So I don't think these settings are
appropriate for what you're trying to do.
>
>
> NOTE 1 [file md_SOL_rerun.mdp]:
> nstcomm < nstcalcenergy defeats the purpose of nstcalcenergy, setting
> nstcomm to nstcalcenergy
>
>
> WARNING 2 [file md_SOL_rerun.mdp]:
> You are generating velocities so I am assuming you are equilibrating a
> system. You are using Parrinello-Rahman pressure coupling, but this can
> be unstable for equilibration. If your system crashes, try equilibrating
> first with Berendsen pressure coupling. If you are not equilibrating the
> system, you can probably ignore this warning.
>
Harmless for the purposes of a re-run calculation.
-Justin
> Generated 168 of the 1653 non-bonded parameter combinations
> Excluding 3 bonded neighbours molecule type 'M_R'
> turning all bonds into constraints...
> Excluding 2 bonded neighbours molecule type 'SOL'
> turning all bonds into constraints...
> Velocities were taken from a Maxwell distribution at 300 K
>
> NOTE 2 [file topol.top]:
> The largest charge group contains 11 atoms.
> Since atoms only see each other when the centers of geometry of the charge
> groups they belong to are within the cut-off distance, too large charge
> groups can lead to serious cut-off artifacts.
> For efficiency and accuracy, charge group should consist of a few atoms.
> For all-atom force fields use: CH3, CH2, CH, NH2, NH, OH, CO2, CO, etc.
>
> Analysing residue names:
> There are: 1 Other residues
> There are: 859 Water residues
> Analysing residues not classified as Protein/DNA/RNA/Water and splitting
> into groups...
> Number of degrees of freedom in T-Coupling group M_R is 106.94
> Number of degrees of freedom in T-Coupling group SOL is 5151.06
>
> NOTE 3 [file md_SOL_rerun.mdp]:
> You are using a plain Coulomb cut-off, which might produce artifacts.
> You might want to consider using PME electrostatics.
>
>
> Largest charge group radii for Van der Waals: 0.370, 0.293 nm
> Largest charge group radii for Coulomb: 0.370, 0.307 nm
> This run will generate roughly 199 Mb of data
>
> There were 3 notes
>
> There were 2 warnings
>
> -------------------------------------------------------
> Program grompp, VERSION 4.6.5
> Source code file: /build/buildd/gromacs-4.6.5/src/kernel/grompp.c, line:
> 1910
>
> Fatal error:
> Too many warnings (2), grompp terminated.
> If you are sure all warnings are harmless, use the -maxwarn option.
> For more information and tips for troubleshooting, please check the GROMACS
> website at http://www.gromacs.org/Documentation/Errors
> -------------------------------------------------------
>
> "Garbage Collecting..." (GNU Emacs)
>
> And this is my new .mdp file
>
> title = only ligand MD simulation for rerun option
> constraints = all-bonds
> integrator = md
> dt = 0.002 ; ps !
> nsteps = 5000000 ; total 10000 ps.
> nstcomm = 1
> nstxout = 5000
> nstvout = 5000
> nstfout = 0
> nstlog = 5000
> nstenergy = 500
> nstxtcout = 500
> nstlist = 10
> ns_type = grid
> rlist = 0.9 ; short-range neighborlist cutoff
> (in nm)
> rcoulomb = 0.9 ; short-range electrostatic cutoff
> (in nm)
> vdwtype = cut-off
> rvdw = 1.4 ; gromos96 force field/ short-range
> van der Waals cutoff (in nm)
> ;coulombtype = PME ; Particle Mesh Ewald for long-range
> electrostatics
> ;fourierspacing = 0.12 ; grid spacing for FFT
> ;fourier_nx = 0
> ;fourier_ny = 0
> ;fourier_nz = 0
> ;pme_order = 4
> ;ewald_rtol = 1e-5
> *coulombtype = Reaction-Field-zero
> epsilon_rf = 0*
> optimize_fft = yes
> pbc = xyz
> ; modified Berendsen temperature coupling is on in two groups
> Tcoupl = V-rescale ; modified Berendsen thermostat
> tc-grps = M_R SOL
> tau_t = 0.1 0.1
> ref_t = 300 300
> ; Energy monitoring
> energygrps = M_R SOL
>
> ; Mode for center of mass motion removal
> comm-mode = Linear
> ; Groups for center of mass motion removal
> comm-grps = System
>
> ; Pressure coupling is now on
> Pcoupl = parrinello-rahman
> Pcoupltype = isotropic
> tau_p = 1.0
> compressibility = 4.5e-5
> ref_p = 1.0
>
> ; Generate velocites is on at 300 K.
> gen_vel = yes
> gen_temp = 300.0
> gen_seed = 173529
>
> Please suggest me from where I can take the new .mdp settings for
> performing mdrun with rerun option. Because I know that it is necessary to
> postprocess the simulation for g_lie calculation.
> I got this file from the archive itself.
>
> Thanks in advance...
>
> Best Regards...
>
>
> On Thu, Jan 7, 2016 at 12:16 PM, shagun krishna <krishna.shagun123 at gmail.com
>> wrote:
>
>> Hiii Justin,
>>
>> Thanks again for replying to my query. I just wanted to know your opinion
>> that while performing ligand alone simulation in water, is it correct to
>> just take the .mdp files from lysozyme tutorial and change accordingly
>> without changing any parameters like Bond parameter,Neighborsearching,
>> Electrostatics, Temperature coupling and pressure coupling. Or can we take
>> these setting (except md; because in md.mdp it is using sd integrator, how
>> ever in all other settings it is using md integrator) from the Solvation
>> free energy of ethanol tutorial by Sander Pronk
>>
>> I am using nvery expensive settings, I guess that this may be a probable
>> reason for getting higher values of Elj and Eqq. Kindly let me know that
>> from where should I select a the .mdp setting for simulation of ligand in
>> water only, without going through free energy codes.
>>
>> Thanks in advance.
>>
>> Regards,
>>
>> Shagun
>>
>> On Wed, Jan 6, 2016 at 11:27 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>>
>>>
>>>
>>> On 1/6/16 12:52 PM, shagun krishna wrote:
>>>
>>>> Thanks for your reply Justin. I know it well that your that tutorial is
>>>> nothing to do with g_lie. But being as a gromacs user I am not getting
>>>> this
>>>> point so I am asking this question.
>>>>
>>>
>>> I was pointing it out for clarity in the archive only; I know you
>>> understand this point, but when someone types in the right words to Google,
>>> they may get seriously confused by what they see :)
>>>
>>> For the protein I am taking all the default values which has been set by
>>>> g_lie program. And for the ligand we have to calculate it using g_energy
>>>> module. I have done with that part also. But f*or my ligand I am
>>>> obtaining
>>>> Elj (LJ-14)= 226.365 and Eqq (*Coulomb-14) *= -125.048 (both in
>>>> KJ/mol).**
>>>> When I am comparing them with their respective counterpart for the
>>>> protein
>>>> (Clj and Cqq value 0.181 and 0.5 respectively; the Elj and Eqq values are
>>>> set 0 for the protein), they seem to be very large. Is it natural or
>>>> there
>>>> is something wrong in my procedure. *
>>>>
>>>
>>> 1-4 interactions are intramolecular terms. These are not what you want.
>>>
>>> You need:
>>>
>>> (1) Protein-ligand interaction energy (from energygrps = Protein LIG) in
>>> the complex simulation. You pass this .edr file to gmx lie. Of course, if
>>> water plays any role in the interaction of the protein, or it is partially
>>> hydrated in the binding site, the result you get will not be right.
>>> (2) Ligand-water interaction energy (from energygrps = LIG SOL) from a
>>> simulation of the ligand in water. You pass these values to -Eqq and -Elj
>>> to gmx lie.
>>>
>>> As for the -Clj and -Cqq values, I still don't know. Hopefully someone
>>> who does will chime in with something useful. But again, the primary
>>> literature (and there are a number of LIE papers, I'm just not well-versed
>>> in their details) should tell you everything you need.
>>>
>>> -Justin
>>>
>>> *PS I have performed simple gromacs simulation of ligand in water (that
>>>> has
>>>> not included your free energy codes).*
>>>>
>>>> *Regards*
>>>>
>>>> *Shagun*
>>>>
>>>>
>>>> On Wed, Jan 6, 2016 at 9:42 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>>>>
>>>>
>>>>>
>>>>> On 1/6/16 1:20 AM, shagun krishna wrote:
>>>>>
>>>>> Dear Gmx-user,
>>>>>>
>>>>>> I am using LIE approach to calculate the binding energy of my protein
>>>>>> and
>>>>>> ligand. When I am running g_lie program for my protein I am getting the
>>>>>> value of DeltaG= 0. Can you please suggest me a way to get the values
>>>>>> of
>>>>>> Elj, Eqq, Clj and Cqq and also about the post-processing necessary to
>>>>>> run
>>>>>> a
>>>>>> LIE calculation. It is not given in tutorial. I have ran two separate
>>>>>> simulations: one for receptor-ligand complex and another for ligand
>>>>>> alone
>>>>>> in solvent. I am following Justin tutorial for free energy
>>>>>> calculation.
>>>>>>
>>>>>>
>>>>>> For the sake of clarity in the archive, I again reiterate that my free
>>>>> energy tutorial has nothing to do with LIE.
>>>>>
>>>>> The values of Elj and Eqq are what you get from simulating the ligand in
>>>>> water with appropriate energygrps.
>>>>>
>>>>> The values of Clj and Cqq are the parameters you have to figure out.
>>>>> What
>>>>> do the papers on LIE say? There are defaults given in gmx lie, where do
>>>>> they come from?
>>>>>
>>>>> -Justin
>>>>>
>>>>> --
>>>>> ==================================================
>>>>>
>>>>> Justin A. Lemkul, Ph.D.
>>>>> Ruth L. Kirschstein NRSA Postdoctoral Fellow
>>>>>
>>>>> Department of Pharmaceutical Sciences
>>>>> School of Pharmacy
>>>>> Health Sciences Facility II, Room 629
>>>>> University of Maryland, Baltimore
>>>>> 20 Penn St.
>>>>> Baltimore, MD 21201
>>>>>
>>>>> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
>>>>> http://mackerell.umaryland.edu/~jalemkul
>>>>>
>>>>> ==================================================
>>>>> --
>>>>> Gromacs Users mailing list
>>>>>
>>>>> * Please search the archive at
>>>>> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
>>>>> posting!
>>>>>
>>>>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>>>>
>>>>> * For (un)subscribe requests visit
>>>>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>>>>> send a mail to gmx-users-request at gromacs.org.
>>>>>
>>>>>
>>> --
>>> ==================================================
>>>
>>> Justin A. Lemkul, Ph.D.
>>> Ruth L. Kirschstein NRSA Postdoctoral Fellow
>>>
>>> Department of Pharmaceutical Sciences
>>> School of Pharmacy
>>> Health Sciences Facility II, Room 629
>>> University of Maryland, Baltimore
>>> 20 Penn St.
>>> Baltimore, MD 21201
>>>
>>> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
>>> http://mackerell.umaryland.edu/~jalemkul
>>>
>>> ==================================================
>>> --
>>> Gromacs Users mailing list
>>>
>>> * Please search the archive at
>>> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
>>> posting!
>>>
>>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>>
>>> * For (un)subscribe requests visit
>>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>>> send a mail to gmx-users-request at gromacs.org.
>>>
>>
>>
--
==================================================
Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow
Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201
jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul
==================================================
More information about the gromacs.org_gmx-users
mailing list