How to use compute TI with soft-core potentials

LAMMPS
1

Hi all,

Since the soft core potential already have lambda which scales the interaction, how exactly should I write the compute ti command? I am attatching the input file that I am using here:

# -- Initialise system --
units real
atom_style full

# -- read data --
read_data system.data

# -- set charges --
set  type  1  charge  -0.8476
set  type  2  charge  0.4238
set  type  3  charge  -0.683
set  type  4  charge  0.418
set  type  5  charge  0.04
set  type  6  charge  0.145

# -- set variables for TI --
variable  lj   equal  0.3
variable  lj1   equal  1
variable  dlj  equal  1
variable  ks   equal  1
variable  ks1   equal  0
variable  dks  equal  1

# -- non bonded interactions --
pair_style lj/cut/coul/long/soft 1  0.5 0.0 10.0 10.0
pair_coeff  1  1  0.1554     3.16557 1.0  
pair_coeff  2  2  0.0        0.0     1.0  
pair_coeff  3  3  0.17       3.12    1.0  
pair_coeff  4  4  0.0        0.0     1.0  
pair_coeff  5  5  0.03       2.5     1.0  
pair_coeff  6  6  0.066      3.5     1.0  
# mixed rules
pair_coeff  1  2  0          0       1.0   
pair_coeff  1  3  0.162536   3.1427  v_lj  
pair_coeff  1  4  0          0       v_lj  
pair_coeff  1  5  0.0682788  2.81317 v_lj  
pair_coeff  1  6  0.101274   3.32859 v_lj  
pair_coeff  2  3  0          0       v_lj  
pair_coeff  2  4  0          0       v_lj  
pair_coeff  2  5  0          0       v_lj  
pair_coeff  2  6  0          0       v_lj  
pair_coeff  3  4  0          0       1.0   
pair_coeff  3  5  0.0714143  2.79285 1.0   
pair_coeff  3  6  0.105925   3.30454 1.0   
pair_coeff  4  5  0          0       1.0   
pair_coeff  4  6  0          0       1.0   
pair_coeff  5  6  0.0444972  2.95804 1.0   


# -- kspace --
kspace_style pppm 0.0001

# -- special bonds --
special_bonds lj/coul 0.0 0.0 0.5

# -- bonded interactions --
bond_style harmonic
bond_coeff 1 553.0 0.945
bond_coeff 2 320.0 1.41
bond_coeff 3 340.0 1.09
bond_coeff 4 600.0 1.0

# -- angular potential --
angle_style harmonic
angle_coeff 1 55.0 108.5
angle_coeff 2 33.0 107.8
angle_coeff 3 35.0 109.5
angle_coeff 4 75.0 109.47

# -- dihedral potential ( only methanol ) --
dihedral_style opls
dihedral_coeff 1 0.0 0.0 0.352 0.0


# -- write ff data --
write_data ff.data pair ij

# -- water shake --
group spc type 1:2
fix myshake spc shake 0.0001 10 100 b 4 a 4

# -- TI --
compute TI all ti lj/cut/coul/long 3*6 v_lj1 v_dlj kspace * v_ks1 v_dks

# -- simulation protocol --
timestep 0.5
minimize 0.0 0.0 1000 10000
velocity all create 300.0 5463576
fix relax1 all npt temp 300.0 300.0 50.0 iso 1.0 1.0 100.0
run 20000
unfix relax1

fix relax2 all npt temp 300.0 300.0 50.0 iso 1.0 1.0 100.0
thermo_style custom step temp press density c_TI
thermo 10
log cg.log
run 20000
write_data prod.data nocoeff

Does my TI command make sense? My resoning here is that since the lambda is already incorporated in the soft-core potential, I dont have to scale again, so I have the scaling factor in compute ti as 1. Is that right? If not, please let me know how to do this.

Thank you and regards,
Ved

2

Sorry, but this is not the place to learn how to do TI. For that you have to consult the corresponding text books and publications. In general, free energy methods can be quite complex and errors may be subtle and difficult to detect. So you are well advised to make some simple tests and try to reproduce published results and follow protocols described there. If you need in-person tutoring instead or in addition, you would have to find a suitable tutor or collaborator. This is essentially true for any “advanced” method included in LAMMPS except for the few cases where the developer of such method is active here, but that applies to a very, very small subset.