[gmx-developers] Free energy calculations in 3.3_rc1

[David Evans]

Hi,

I just tested the (single processor) 3.3_rc1 version on one
of my systems, involving perturbing a ligand in water.
I compared the results to those from identical input
files run with a CVS version checked out on 23/11/2004.
I have been experimenting with this CVS version for a while,
since I believed it to include free-energy PME loops (see previous postings
on this list).

I found significantly (100 kJ) different values of dg/dl between the two 
versions,
although when you integrate across lambda the sums differ by only 1-2 kJ
(which is still significant according to my error analysis, but not by 
much).
Does anyone knows why this is, and whether the 3.3_rc1 version
is likely to be more accurate? My .mdp file is below, and I can provide
any more details anyone needs.

Thanks,

Dave

cpp                 =  /lib/cpp
;define              =  -DFLEX_SPC
integrator          =  md
nsteps              =  1000000
dt                  =  0.002
                                                                                                                                               

comm_mode           = linear
nstcomm             = 1000
comm_grps           = system
                                                                                                                                               

nstxout             = 500000
nstvout             = 500000
nstlog              = 100
nstenergy           = 100
nstxtcout           = 2000
                                                                                                                                               

tcoupl              = nose-hoover
tc_grps             = system
tau_t               = 1
ref_t               = 300
gen_vel             = yes
gen_temp            = 300
                                                                                                                                               

pcoupl              = berendsen
pcoupltype          = isotropic
tau_p               = 1
ref_p               = 1
compressibility     = 4.46e-5
                                                                                                                                               

constraints         = hbonds
constraint_algorithm = lincs
                                                                                                                                               

coulombtype         = PME
rvdw                = 0.9
rlist               = 0.9
rcoulomb            = 0.9
fourierspacing      = 0.12
pme_order           = 4
ewald_rtol          = 1e-5
                                                                                                                                               

vdwtype         = cut-off
pbc                 = xyz
dispcorr            = enerpres
                                                                                                                                               

free_energy  = yes
init_lambda = 0.1
sc_alpha  = 1.0

--
David Evans
School of Pharmacy
University of London