Measuring diffusion coefficient using fix smd

Dear Axel and lammps users,

I want to measure the diffusion coefficient of a K+ ion in a water box. My approach is to connect a spring to the ion and pull “the other side of the spring” with a constant velocity. I use fix smd in the following way:

fix pull K smd cvel 250.0 -0.00003 tether 10.0 NULL NULL 0.0

Actually, when I monitor the f_pull[7] which is the accumulated PMF, I do not get what I expect. I expect it to increase as the ion moves forward, but it fluctuates around 0. I will appreciate it if you give me a key to find the problem.

Thanks.

########## input script ##########

variable JobName string ion6
log ${JobName}.lammps
shell rm log.lammps

dimension 3
boundary p p p
units real
neighbor 2 bin
neigh_modify delay 2 every 1 check yes

atom_style full
bond_style harmonic
angle_style charmm

pair_style lj/charmm/coul/long 12 14
pair_modify mix arithmetic
kspace_style pppm 1.0e-6

read_data wt.data

set type 2 charge -0.830
set type 1 charge 0.415
group wt type 1 2
group K type 3

pair_coeff 1 1 0.0 0.0
pair_coeff 2 2 0.102 3.188
pair_coeff 1 2 0.0 1.7753
pair_coeff 3 3 0.0 0.0
pair_coeff 1 3 0.0 0.0
pair_coeff 2 3 0.20943 3.014

balance 1.0 shift xyz 50 1.0
minimize 1.0e-4 1.0e-6 100 1000
velocity all create 300.0 12345678 dist gaussian

compute Temp all temp/com

fix 1 all nvt temp 300.0 300.0 100.0
fix_modify 1 temp Temp
fix 2 wt shake 1e-6 100 0 b 1 a 1
fix 3 K recenter INIT INIT INIT units box

thermo 1000
thermo_style custom elapsed cpu pe ke c_Temp
timestep 1

restart 10000 {JobName}.restart1 {JobName}.restart2

run 10000

unfix 3
fix pull K smd cvel 250.0 -0.00003 tether 10.0 NULL NULL 0.0
variable t equal step
variable fx equal f_pull[1]
variable fy equal f_pull[2]
variable fz equal f_pull[3]
variable F equal f_pull[4]
variable r0 equal f_pull[5]
variable r equal f_pull[6]
variable PMF equal f_pull[7]
variable xk equal xcm(K,x)
fix PMF all print 1000 “t {xk} {fx} {F} {r0} {r} ${PMF}” file PMF_{JobName}.txt

run 10000000

Dear Axel and lammps users,

I want to measure the diffusion coefficient of a K+ ion in a water box. My approach is to connect a spring to the ion and pull “the other side of the spring” with a constant velocity.

what is wrong with computing the (self-)diffusion coefficient the traditional way via MSD or VACF?

axel.

Dear Axel,

My adviser says this is a faster way. I have not had any background about this before.

Dear Axel,

My adviser says this is a faster way. I have not had any background about this before.

then you should ask your adviser about how to handle your problems.

axel.

Dear Axel,

He has done this before via NAMD.
I do not know how it works in lammps. This is why I can not figure out the problem.

Dear Axel,

He has done this before via NAMD.
I do not know how it works in lammps. This is why I can not figure out the problem.

the COLVARS package in LAMMPS is the same software as used in NAMD.

axel.

then what is the problem with this fix smd?just 2 questions:

  1. how can I detect the coordinate of the moving side of the spring?
  2. how it calculates the accumulated PMF?

then what is the problem with this fix smd?

there is no problem with fix smd, it works just like advertised in the documentation.

just 2 questions:

  1. how can I detect the coordinate of the moving side of the spring?

that is prescribed by the starting point and the velocity you are using.

  1. how it calculates the accumulated PMF?

by numerically integrating the force in direction of the spring (essentially summing the product of the projected force times the distance).
if you need more details, you need to read the source.

but as i said, just do SMD with colvars and you can use the same input as in NAMD and you can have your adviser explain things to you.

axel.

Ok, thanks.but I think this part of the manual is so vague as I saw many people have asked questions in the mailing list about that.
maybe adding a picture and defining the parameters visually can help.

The purpose of the manual is to describe how to use the features. It is not to explain the methodology. For that you need to study the corresponding text books and publications.

Axel

I know, but many questions were about the parameters’ definitions.

For example, the manual’s descriptions give the depicture as pic1 (attachment).
But your explanation in mailing list (“consider that you have a …”) depicts it as pic2.

pic1.png

pic2.png

I know, but many questions were about the parameters’ definitions.

For example, the manual’s descriptions give the depicture as pic1 (attachment).
But your explanation in mailing list (“consider that you have a …”) depicts it as pic2.

neither picture is correct. but it is mostly irrelevant because R0 is practically always best chosen as 0.0 for fix smd.

R0 is a feature that was inherited from fix spring on which fix smd is based on.

axel.

neither picture is correct.

then maybe I am totally wrong about the way how this fix works.
I want to model the case as in the attached pic3.
Is that possible with this fix?

pic3.png

Hi all,

I have figured out my problem and changed the fix I was using. Now, I use fix spring between two particles, one of which is considered to have no effect on the system ( lj_epsilon=0 , charge=0). Then I pull the right-side particle with a constant velocity using fix move linear.

Just one question. As I want to calculate the work of the spring force on the K+ ion by the formula W=F_spring * ( X_now - X_last), is there a way to store the last position of the ion in each step ??

Thanks.