[gmx-users] solvation free energy

[Hannes Loeffler]

The couple-* parameters take already care of including the non-bonded
terms internal to your molecule to correctly describe the transfer of M
to vacuum.  That's the point of those parameters so that you would not
have to run an additional "correction" in vacuo.  See the discussion in
the manual (it's section 5.3 for Gromacs 4.6.x).  You would also need to
compute the electrostatic transformation.

On Thu, 16 Jul 2015 16:26:44 +0200
Daniele Veclani <danieleveclani at gmail.com> wrote:

> I'm trying to calculate the  solvation free energy of a molecule (M).
> I have done:
> M+water ---> dum+water
> 
> Now I have to do:
> M(vacuo) --> dum(vacuo)
> 
> In this case I have a problem, in fact I find a DG = 0.00000 and
> within the .xvg file there are only zeros.
> 
> Where is the problem?
> 
> I use gromacs 5.0.4 and my .mdp file (for lambda 00) is:
> integrator               = sd
> tinit                    = 0
> dt                       = 0.002
> nsteps                   = 2500000
> comm_mode                = angular
> nstcomm                  = 100
> nstxout                  = 500
> nstvout                  = 500
> nstfout                  = 0
> nstlog                   = 500
> nstenergy                = 500
> nstxout-compressed       = 0
> cutoff-scheme            = group
> nstlist                  = 0
> ns_type                  = simple
> pbc                      = no
> rlist                    = 0
> coulombtype              = cutoff
> rcoulomb                 = 0
> epsilon_r                = 1
> vdwtype                  = cutoff
> rvdw                     = 0
> DispCorr                  = no
> fourierspacing           = 0.12
> pme_order                = 6
> ewald_rtol               = 1e-06
> epsilon_surface          = 0
> tcoupl                   = berendsen
> tc_grps                  = system
> tau_t                    = 0.1
> ref_t                    = 300
> Pcoupl                   = no
> tau_p                    = 1.0
> compressibility          = 4.5e-05
> ref_p                    = 1.0
> free_energy              = yes
> init_lambda_state        = 0
> delta_lambda             = 0
> calc_lambda_neighbors    = 1
> vdw_lambdas              = 0.00 0.10 0.20 0.30 0.40 0.50 0.60 0.70
> 0.80 0.90 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
> bonded_lambdas           = 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
> mass_lambdas             = 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
> sc-alpha                 = 0.5
> sc-coul                  = no
> sc-power                 = 1.0
> sc-sigma                 = 0.3
> couple-moltype           = QUI
> couple-lambda0           = vdw
> 
> couple-lambda1           = none
> couple-intramol          = no
> nstdhdl                  = 10
> gen_vel                  = no
> constraints              = all-bonds  ; we only have C-H bonds here
> constraint-algorithm     = lincs
> continuation             = no
> lincs-order              = 12
> 
> Best regards
> Daniele Veclani.