[gmx-users] Regarding value of Elj, Eqq, Clj and Cqq values using LIE Approach

shagun krishna krishna.shagun123 at gmail.com
Thu Jan 7 11:54:16 CET 2016


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



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.

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.
>>
>
>


More information about the gromacs.org_gmx-users mailing list