Hello,
I have been trying to test a few free energy methods to calculate the solvation free energy of a LJ solute (B) in the same LJ solvent (A), a very simple model system.
These methods include free energy perturbation (FEP), Widom insertion, and thermodynamic integration (TI) methods.
The results I’ve got are as following:
FEP: -3.61 (0.05)
Widom: -3.78 (0.2)
TI: -4.22 (0.02)
I can’t see why TI method doesn’t give the right answer. Here is my input file:
begin
units lj
dimension 3
boundary p p p
atom_style atomic
read_restart run_b.rest
pair_style hybrid lj/cut 2.5 lj/cut/soft 1.0 0.5 2.5
mass 1 1.0000 # A
mass 2 1.0000 # B
pair_coeff 1 1 lj/cut 1.0 1.0 # A A
pair_coeff 1 2 lj/cut/soft 1.0 1.0 ${Lambda} # A B
pair_coeff 2 2 lj/cut 1.0 1.0 # B B
restart 10000 run_1.rest run_2.rest
thermo 10000
timestep 0.01
run_style verlet
variable lambda equal 1
variable dlambda equal 1
compute TI all ti lj/cut/soft 1*2 v_lambda v_dlambda
fix NPT all npt temp 0.8 0.8 1.0 iso 0.1 0.1 10.0
run 300000
reset_timestep 0
fix thermodynamic all ave/time 1 100 100 c_TI file Thermodynamic.dat
run 1000000
end
Please note that I am doing separate independent simulations in a parallel way by having ${Lambda} from 0 to 1.0 with 0.025 intervals (41 simulations).
I’m using lammps 16Mar18.
Any thoughts will be appreciated.