[gmx-users] Fwd: Running free energy tutorial for ionic molecules

vivek sharma viveksharma.iitb at gmail.com
Sun Nov 9 13:21:24 CET 2014


Hi Dr Warren,
Sorry for my late reply as I was trying few variation in my simulation.
I got over the earlier reported problem by replacing the constraints from
"h-bonds" to "all-bonds" in my mdp file.
I dont have any idea if it will affect my results?.

Now I am trying to run the same simulation set with OPLS-AA parameters,
where the simulation are successful in water system but crashing in
1-octanol system with following warnings:

-------------------------------------
WARNING: Listed nonbonded interaction between particles 842 and 851
at distance 427.422 which is larger than the table limit 2.200 nm.

This is likely either a 1,4 interaction, or a listed interaction inside
a smaller molecule you are decoupling during a free energy calculation.
Since interactions at distances beyond the table cannot be computed,
they are skipped until they are inside the table limit again. You will
only see this message once, even if it occurs for several interactions.

IMPORTANT: This should not happen in a stable simulation, so there is
probably something wrong with your system. Only change the table-extension
distance in the mdp file if you are really sure that is the reason.
-------------------------------------

I am using the same set of steps as tried in earlier simulation the only
difference is I am using OPLS-AA FF parameters (generated using acpype
utility) instead of gromos parameters.  pasting below the mdp file content
used for the simulation.

-----------------------------------------------------------------------------------------------------------
; Run control
integrator               = sd       ; Langevin dynamics
tinit                    = 0
dt                       = 0.002
nsteps                   = 2500000   ; 5 ns
nstcomm                  = 100
; Output control
nstxout                  = 500
nstvout                  = 500
nstfout                  = 0
nstlog                   = 500
nstenergy                = 500
nstxout-compressed       = 0
; Neighborsearching and short-range nonbonded interactions
cutoff-scheme            = verlet
nstlist                  = 20
ns_type                  = grid
pbc                      = xyz
rlist                    = 1.2
; Electrostatics
coulombtype              = PME
rcoulomb                 = 1.2
; van der Waals
vdwtype                  = cutoff
vdw-modifier             = potential-switch
rvdw-switch              = 1.0
rvdw                     = 1.2
; Apply long range dispersion corrections for Energy and Pressure
DispCorr                  = EnerPres
; Spacing for the PME/PPPM FFT grid
fourierspacing           = 0.12
; EWALD/PME/PPPM parameters
pme_order                = 6
ewald_rtol               = 1e-06
epsilon_surface          = 0
; Temperature coupling
; tcoupl is implicitly handled by the sd integrator
tc_grps                  = system
tau_t                    = 1.0
ref_t                    = 300
; Pressure coupling is on for NPT
Pcoupl                   = Parrinello-Rahman
tau_p                    = 1.0
compressibility          = 4.5e-05
ref_p                    = 1.0
; Free energy control stuff
free_energy              = yes
init_lambda_state        = 0
delta_lambda             = 0
calc_lambda_neighbors    = 1        ; only immediate neighboring windows
; Vectors of lambda specified here
; Each combination is an index that is retrieved from init_lambda_state for
each simulation
; init_lambda_state        0    1    2    3    4    5    6    7    8    9
 10   11   12   13   14   15   16   17   18   19   20
vdw_lambdas              = 0.00 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80
0.90 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00
coul_lambdas             = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 1.00
; We are not transforming any bonded or restrained interactions
bonded_lambdas           = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
restraint_lambdas        = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
; Masses are not changing (particle identities are the same at lambda = 0
and lambda = 1)
mass_lambdas             = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
; Not doing simulated temperting here
temperature_lambdas      = 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
; Options for the decoupling
sc-alpha                 = 0.5
sc-coul                  = no       ; linear interpolation of Coulomb (none
in this case)
sc-power                 = 1.0
sc-sigma                 = 0.3
couple-moltype           = FEL  ; name of moleculetype to decouple
couple-lambda0           = none      ; only van der Waals interactions
couple-lambda1           = vdw-q     ; turn off everything, in this case
only vdW
couple-intramol          = no
nstdhdl                  = 10
; Do not generate velocities
gen_vel                  = no
; options for bonds
constraints              = all-bonds  ; we only have C-H bonds here
; Type of constraint algorithm
constraint-algorithm     = lincs
; Constrain the starting configuration
; since we are continuing from NPT
continuation             = yes
; Highest order in the expansion of the constraint coupling matrix
lincs-order              = 12
------------------------------------------------------------------------------------------------------------

It will be really useful to have expert comments and an insight into the
simulation in order to get a solution.

Thanks and regards,
Vivek Sharma

On 16 October 2014 04:35, Dallas Warren <Dallas.Warren at monash.edu> wrote:

> Which bonds are rotating too far?
>
> At what lambda values and which decoupling process is this occurring?
>
> Catch ya,
>
> Dr. Dallas Warren
> Drug Delivery, Disposition and Dynamics
> Monash Institute of Pharmaceutical Sciences, Monash University
> 381 Royal Parade, Parkville VIC 3052
> dallas.warren at monash.edu
> +61 3 9903 9304
> ---------------------------------
> When the only tool you own is a hammer, every problem begins to resemble a
> nail.
>
> > -----Original Message-----
> > From: gromacs.org_gmx-users-bounces at maillist.sys.kth.se
> > [mailto:gromacs.org_gmx-users-bounces at maillist.sys.kth.se] On Behalf Of
> > vivek sharma
> > Sent: Wednesday, 15 October 2014 8:07 PM
> > To: Discussion list for GROMACS users
> > Subject: [gmx-users] Fwd: Running free energy tutorial for ionic
> > molecules
> >
> > Dear Users,
> > I am resending this mail, as I am still stuck with the error. Now I
> > have
> > upgraded my Gromacs from version-4.5.5. to version-5.0.2, I am able to
> > run
> > the simulation using Gromacs-5.0.2 for sytems with 1-octanol as
> > solvent,
> > whereas the system with water as solvent is still giving errors with
> > LINCS
> > warnings.
> > I assume if the runs are going fine with water system then my topology
> > and
> > MDP (adapted from tutorial) files are in place, Please correct me if
> > this
> > is a wrong assumption.
> > It will be really helpfull if anybody can help me in understanding the
> > system behaviour or possible reasons for the error.
> >
> >
> > regards,
> > Vivek
> >
> > ---------- Forwarded message ----------
> > From: vivek sharma <viveksharma.iitb at gmail.com>
> > Date: 10 October 2014 09:01
> > Subject: Running free energy tutorial for ionic molecules
> > To: Discussion list for GROMACS users <gmx-users at gromacs.org>
> >
> >
> > Hi users and experts,
> >
> > I am trying to run a free energy tutorial by Dr Lemkul provided at:
> > http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-
> > tutorials/free_energy/
> > My system is an ionic pyridinium molecule unlike the methane molecvule
> > used
> > in the tutorial, following are the steps I am following for my runs:
> > 1. I downloaded the molecule structure and topology from ATB site:
> > http://compbio.biosci.uq.edu.au/atb/index.py?tab=existing_tab&nocache=1
> > 53
> > 2. I set up the system using octanol water box.
> > 3. I am using the MDP parameters as suggested in the advanced
> > application
> > section: decoupling vdw and coulombic interaction step by step.
> > 4. WHile running the simulation it goes fine till npt equilibration
> > step,
> > but it fails during the production run with the LINCS warnings, like
> > some
> > angles are rotating more than 30 degree.
> >
> > I am unable to figure out if I can apply this protocol on ionic
> > molecules
> > or I am doing something wrong in setting up the system.
> >
> > Can somebody help me in understanding the possible reason for error?
> >
> > It will be really helpful to get more insights into this type of
> > simulation.
> >
> > Thanks and regards,
> > Vivek
> > --
> > 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.
> --
> 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