[gmx-users] Free Energy Calculations: Error during minimization step.

Abhishek Acharya abhi117acharya at gmail.com
Tue Sep 12 15:36:15 CEST 2017


Hello Wes,

I ran some more test simulations. I think the error is because I added
counterions to neutralize my system (it has a +2 charge) prior to the free
energy. Without the counterions, the simulations are running fine. Though,
I don't know whether such a behavior is expected. Such a setup (with
counterions) is not something unheard of.

Nevertheless, I read elsewhere that for solvation free energy calculation
of charged system, the uniform background charge as implemented in PME for
non-neutral systems may be desirable. Hence, I plan to run the simulations
without counterions.

Sincerely,
Abhishek




Abhishek Acharya
Research Associate,
Department of Molecular Nutrition
CSIR-Central Food Technological Research Institute,
Mysuru-570020

On Tue, Sep 12, 2017 at 5:22 PM, Wes Barnett <w.barnett at columbia.edu> wrote:

> On Tue, Sep 12, 2017 at 12:17 AM, Abhishek Acharya <
> abhi117acharya at gmail.com
> > wrote:
>
> > Just to add to the above post, running a none to vdw-q tranformation
> using
> > the following parameters also results in the same error for simulation at
> > first lambda.
> >
>
>
> Do you get this error on a normal simulation as well? That is, if you run
> state 0 (with vdw and electrostatics on) without doing an free energy
> calculations, do you still get the error? Perhaps you are using too big a
> of a time step or your system is not equilibrated enough before using this
> time step. Or you may just have a bad starting configuration.
>
>
> >
> > ************************************************************
> > ************************************************************
> > ********************************************************
> > vdw_lambdas              = 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40
> > 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.95 1.00 1.00 1.00
> 1.00
> > 1.00 1.00 1.00 1.00 1.00 1.00 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.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.05 0.10
> 0.15
> > 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85
> 0.90
> > 0.95 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 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 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 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 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
> > sc-sigma                 = 0.3
> > couple-moltype           = Protein  ; name of moleculetype to decouple
> > couple-lambda0           = none
> > couple-lambda1           = vdw-q
> > couple-intramol          = no
> > nstdhdl                  = 10
> > ************************************************************
> > ************************************************************
> > ********************************************************
> >
> > Thanks in advance.
> >
> > Abhishek Acharya
> > Research Associate,
> > Department of Molecular Nutrition
> > CSIR-Central Food Technological Research Institute,
> > Mysuru-570020
> >
> > On Tue, Sep 12, 2017 at 9:42 AM, Abhishek Acharya <
> > abhi117acharya at gmail.com>
> > wrote:
> >
> > > Hello GROMACS users,
> > >
> > > As Wes had advised previously in this thread, the issue with the vdw
> > > tranformations was rectified by removing the charges from the topology
> > > file. Thereafter, I also tried running free energy calculations using
> the
> > > following free-energy parameters:
> > >
> > > ************************************************************
> > > ************************************************************
> > > ********************************************************
> > > 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 for each
> simulation
> > > coul_lambdas              = 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35
> 0.40
> > > 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85 0.90 0.95 1.00 1.00 1.00
> > 1.00
> > > 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00
> > 1.00
> > > 1.00 1.00
> > > vdw_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 0.05 0.10
> > 0.15
> > > 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.55 0.60 0.65 0.70 0.75 0.80 0.85
> > 0.90
> > > 0.95 1.00
> > > 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 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 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
> > > 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 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
> > > 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 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
> > > sc-sigma                 = 0.3
> > > couple-moltype           = Protein  ; name of moleculetype to decouple
> > > couple-lambda0           = vdw-q
> > > couple-lambda1           = none
> > > couple-intramol          = no
> > > nstdhdl                  = 10
> > >
> > > ************************************************************
> > > ************************************************************
> > > ********************************************************
> > >
> > > As you can see, I want to first annihilate the charges, followed by a
> vdw
> > > tranformation, using a total of 41 lambdas. Since both charges and vdw
> > are
> > > being interpolated, the topology files contain the normal charges for
> the
> > > transformed molecule. This setup is similar to what Wes had suggested
> in
> > a
> > > previous response to this thread.
> > >
> > > The issue is that I am facing the exactly the same error in
> minimization
> > > step at final lambda state where both vdw and electrostatics are turned
> > off.
> > >
> > > ------------------------------------------------------------
> > > ------------------------------------------------------------
> > > ------------------------------------------------------
> > > Program gmx mdrun, VERSION 5.1.2
> > > Source code file: /home/bp-lab/Downloads/gromacs-5.1.2/src/gromacs/
> > mdlib/constr.cpp,
> > > line: 555
> > >
> > > Fatal error:
> > >
> > > step 10: Water molecule starting at atom 120 can not be settled.
> > > Check for bad contacts and/or reduce the timestep if appropriate.
> > >
> > > For more information and tips for troubleshooting, please check the
> > GROMACS
> > > website at http://www.gromacs.org/Documentation/Errors
> > >
> > > ------------------------------------------------------------
> > > ------------------------------------------------------------
> > > -------------------------------------------------------
> > >
> > > Now, what I understood is that in case of just the vdw transformation,
> > > since we are not calculating any electrostatics, the charges should be
> > zero
> > > in topology file. Whereas, in case of just the charge transformation,
> the
> > > charges are required in the topology file as charges are interpolated
> in
> > > this step. And performing the simulation as such give no errors.
> > >
> > > I do not understand why the error is surfacing in case of a [vdw-q] to
> > > [none ] tranformation., especially now that the charges are gradually
> > > decoupled prior to decoupling the vdw interactions.  Any help will be
> > > greatly appreciated.
> > >
> > > Sincerely,
> > > Abhishek
> > >
> > >
> > > Abhishek Acharya
> > > Research Associate,
> > > Department of Molecular Nutrition
> > > CSIR-Central Food Technological Research Institute,
> > > Mysuru-570020
> > >
> > > On Sun, Sep 10, 2017 at 8:02 AM, Abhishek Acharya <
> > > abhi117acharya at gmail.com> wrote:
> > >
> > >>
> > >>
> > >>
> > >> You can do them separately, but the charges need to have been turned
> off
> > >> for the molecule you are transforming for the second step when turning
> > off
> > >> vdw. In other words, all charges must be 0.0 for the molecule in the
> > >> topology file for the vdw transformation. Is that what you've done?
> > >>
> > >> Justin's tutorial gives an explanation on why this needs to be done,
> and
> > >> in
> > >> his topology he had all charges at 0.0 in the methane topology file.
> > >>
> > >>
> > >> Okey, one thing that I presumed while performing the second
> > >> transformation is that it is sufficient to 'turn-off' charge
> > electrostatics
> > >> just by stating c-lambda0=none and c-lambda1 = vdw (or vice-versa). It
> > >> doesn't matter if the topology files have the charges or not.
> > Therefore, i
> > >> had included the charges in my topology file.
> > >>
> > >> Then this explains the problem in my case. Thank you so much Wes for
> > >> clarifying that out.
> > >>
> > >> Sincerely,
> > >> Abhishek Acharya
> > >>
> > >>
> > >>
> > >
> > --
> > 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.
> >
>
>
>
> --
> James "Wes" Barnett
> Postdoctoral Research Scientist
> Department of Chemical Engineering
> Kumar Research Group <http://www.columbia.edu/cu/kumargroup/>
> Columbia University
> w.barnett at columbia.edu
> http://wbarnett.us
> --
> 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