----- visualize your results -----
First of all, if you haven't done this yet, be sure to look at the
trajectory of the water molecules. (Don't "wrap" the molecule
coordinates to lie within the periodic boundaries. Let them diffuse
outside the boundaries.)
Are they moving?
Instructions for visualizing a trajectory using VMD can be found here:
https://github.com/jewettaij/moltemplate/blob/master/examples/all_atom/force_field_OPLSAA/butane/README_visualize.txt
Thank you for your response, I am using a periodic boundary along the axial direction of the CNT and it is completely filled with water molecules to get density of 998kg/m3, and the NVT thermostat is applied with no net flow.
Great
Also as far my knowledge "compute MSD" performs calculations on all the mobile molecules and gives the net MSD.
No. The "compute msd" command in LAMMPS does not know which atoms
are mobile and will not select them for you aromatically. "Compute
msd" performs the average only on the atoms belonging to the group-ID
that you specify. The group-ID is the 2nd argument following the
compute command:
https://lammps.sandia.gov/doc/compute_msd.html
It's important that you are clear about what atoms are in this group.
(If you defined a group of atoms with the group command which you
named "mobile" https://lammps.sandia.gov/doc/group.html ... then
"compute msd" will consider only those atoms belonging to the group
with that name. But you must be clear about which atoms belong to
this group. Don't use the "all" group. Your group should exclude all
the atoms in the nanotube, and the water outside the nanotube, if
there are any.)
--- writing your own script ----
In my previous email I suggested writing your own script to read
the LAMMPS trajectory file ("dump" file). This is easy enough to do
that it might be worth it. When you write you own code, you know
exactly what it does. It's less work initially to try using a 3rd
party tool like "compute msd", but it's hard to know what to do if
something doesn't behave the way you expect. If you are getting
strange results using "compute msd", it is relatively easy to extract
the coordinates from a LAMMPS dump file (ie. trajectory file). You
can do this using a combination of awk and pizza.py (or "dump2data.py
-raw < traj.lammpstrj"), or writing your own code to interpret the
trajectory file.
https://pizza.sandia.gov/
https://github.com/jewettaij/moltemplate/blob/master/doc/utils/docs_dump2data.txt
(You will need to discard the atoms belonging to the carbon nanotube.
Perhaps you could track just one oxygen atom, for example.)
Then you can write your own simple for-loop to calculate the MSD.
Also, when you create the dump file, be sure to include the image
flags and the molecule-ID numbers. This will make it easier to keep
track of molecules as they pass through the periodic boundaries. You
can do this using a "dump" command which includes columns for "mol,
ix, iy, iz", for example:
dump mydump1 all custom 200 traj.lammpstrj id mol type x y z ix iy iz
I don't want to discourage you from getting to the bottom of this
problem with "compute msd". If there is a bug in "compute msd", then
the developers should know about it. But be sure to double-check that
you are supplying it with a group of atoms which makes sense before
writing a bug report.
I am getting a curve that initially increases but then the fluctuations kicks in and it becomes pretty unstable.
In the past something similar would happen when I have calculate
MSD. I wrote my own code to calculate the MSD-vs-time, and I was only
tracking a single particle. I would get huge fluctuations in MSD for
time separations on the order of the length of the simulation. In
that case, it was unavoidable because there are only a few pairs of
time points in the simulation which are separated by long times (and
these were correlated with each other), compared to the many pairs of
time points which are separated by short times (which tended to be
uncorrelated). As a result, when I graphed MSD-vs-time, at moderately
short times times (larger than the velocity auto-correlation time, but
less than 1/10th of the simulation duration) the MSD would increase
linearly with time. Time points larger than say 1/10th of my
simulations were justifiably untrustworthy and were discarded.
However, this is probably less of an issue when calculating MSD using
the "compute msd" command in LAMMPS, because MSD is calculated using
the trajectories of many particles, not just one.
Finally, if you haven't tried this yet, try using the "-com" argument
with "compute msd".
Cheers
Andrew