[gmx-users] Re: Free energy calculation question
Fabian Casteblanco
fabian.casteblanco at gmail.com
Wed Aug 31 20:26:08 CEST 2011
Hello Justin,
You mean that only for vdW decoupling, you would need to use soft-core
potentials? I had soft -core potentials on for decoupling the
electrostatic interactions (see below). What would I use in its
place? Thanks again for your help!
;Production MD
;------------------------------------------------------------------------------------
include =-I/mphase/users2/fabian/CGenff
title =CGenFF Lov/Eth Solution MD run
;PARAMETERS - describing what to do, when to stop and what to save
;------------------------------------------------------------------------------------
;Run parameters
integrator =sd ;leap-frog integrator
ld_seed =-1
nsteps =500000 ;2*500000=1000ps, 1 ns
dt =0.002 ;2 fs
nstcomm = 100 ;*** - frequency for center of mass
motion removal
;Output control
nstxout =1000 ;save coordinates every 2ps
nstvout =1000 ;save velocities every 2 ps
nstenergy =1000 ;save energies every 2 ps
nstlog =1000 ;update log file every 2 ps
nstxtcout =1000 ;xtc compressed trajectory output every 2 ps
xtc-precision =1000 ;*** - precision to write to xtc trajectory
;Bond Parameters
continuation =yes ;Restarting after NPT
constraint_algorithm =lincs ;holonomic constraints
constraints =all-bonds ;all bonds (even heavy atom-H bonds) constr;ained
lincs_iter =1 ;accuracy of LINCS
lincs_order =12 ;also related to accuracy
;Neighborhood searching
ns_type =grid ;search neighboring grid cells
nstlist =5 ;#steps. 5*0.002 ps = 5* 2 fs = 10 fs - Frequency to
update the neighbor list (and the long-range forces, when ;using
twin-range cut-off’s). When this is 0, the neighbor list is made only
once.
rlist =1.1 ;short-range neighborlist cutoff (in nm)
rcoulomb =1.1 ;short-range electrostatic cutoff (in nm)
pbc =xyz ; 3-D PBC
;Electrostatics
coulombtype =PME ;Particle Mesh Ewald for long-range electrostat;ics
pme_order =4 ;cubic interpolation
fourierspacing =0.16 ;grid spacing for FFT
; van der Waals
vdwtype =Shift ;Van der Waals for CHARMM
rvdw_switch =0.8
rvdw =1.0 ;Short-range Van der Waals cut-off
;Dispersion correction
DispCorr =EnerPres ;account for cut-off vdW scheme
;Temperature coupling is on
tcoupl =V-rescale ;modified Berendsen thermostat
tc-grps =SYSTEM ;two coupling groups - more accurate
tau_t =0.1 ;time constant, in ps
ref_t =298 ;reference temperature, on for each group, in K
;Pressure coupling is on
pcoupl =Parrinello-Rahman ;Pressure coupling on in NPT
pcoupltype =isotropic ;uniform scaling of box vect;ors
tau_p =2.0 ;time constant, in ps
ref_p =1.0 ;reference pressure, in bar
compressibility =4.5e-5 ;isothermal compr of H2O, ba;r^(-1)
; Free energy control stuff
free_energy = yes ;*** - Indicates we are doing a free
energy calculation, and that interpolation between the A and B states
of the ;chosen molecule (defined below) will occur.
init_lambda = 0.0 ;*** - Value of λ
delta_lambda = 0 ;*** - The value of λ can be incremented
by some amount per timestep (i.e., δλ/δt) in a technique called "slow
;growth." This method can have significant errors associated
with it, and thus we will make no time-dependent ;changes to our
λ values.
foreign_lambda = 0.05 ;*** - Additional values of λ for
which ΔH will be written to dhdl.xvg (with frequency nstdhdl). The
;configurations generated in the trajectory at λ = init_lambda
will have ΔH calculated for these same ;configurations at all
values of λ = foreign_lambda.
sc-alpha = 0.5 ;*** - the soft-core parameter, a value
of 0 results in linear interpolation of the LJ and Coulomb
;interactions
sc-power = 1.0 ;*** - the power for lambda in the
soft-core function, only the values 1 and 2 are supported
sc-sigma = 0.3 ;*** - the soft-core sigma for
particles which have a C6 or C12 parameter equal to zero or a sigma
;smaller than sc_sigma
couple-moltype = LOV ;*** - name of moleculetype to decouple
couple-lambda0 = vdw-q ;*** - all interactions are on at lambda=0
couple-lambda1 = vdw ;*** - only vdw interactions are on
at lambda=1
couple-intramol = no ;*** - Do not decouple intramolecular
interactions. That is, the λ factor is applied to only solute-solvent
;nonbonded interactions and not solute-solute nonbonded
interactions.
nstdhdl = 100 ;*** - The frequency with which ∂H/∂λ
and ΔH are written to dhdl.xvg. A value of 100 would probably suffice,
since ;the resulting values will be highly correlated and the
files will get very large. You may wish to increase this ;value
to 100 for your own work.
;Velocity generation
gen_vel =no ;Velocity generation is off
;END
On Wed, Aug 31, 2011 at 11:24 AM, Fabian Casteblanco
<fabian.casteblanco at gmail.com> wrote:
> Hello Justin,
>
> I'm calculating the free energy of a drug in an alcohol solvent. I
> have a question referring to your free energy tutorial. You mentioned
> that decoupling of electrostatic interactions is linear and decoupling
> of vdW can vary. Is this true for your case of methanol in water or
> for all cases? When I ran my system of drug in ethanol solvent, I got
> a non linear dG for both electrostatic and vdW. Also, is there no
> need to find dG of cav ( the free energy required to form the solute
> cavity within the solvent) ? I have attached some graphs.
>
> Thanks for your help. Your tutorial was extremely useful.
>
> --
> Best regards,
>
> Fabian F. Casteblanco
> Rutgers University --
> Chemical Engineering PhD Student
> C: +908 917 0723
> E: fabian.casteblanco at gmail.com
>
--
Best regards,
Fabian F. Casteblanco
Rutgers University --
Chemical Engineering PhD Student
C: +908 917 0723
E: fabian.casteblanco at gmail.com
More information about the gromacs.org_gmx-users
mailing list