Mean square deviation
=====================

The ``msd`` tool computes **mean square deviation (MSD)** using the Einstein
relation, based on the displacement of atom groups over time.  
It supports:

- MSD in 1D, 2D, or full 3D  
- interactive or scripted atom-group selection  
- PBC-aware trajectory processing (``NoJump``)  
- customizable time windows and sampling intervals  
- output to an XVG file

This tool is implemented in ``msd.py``.


Overview
--------

Mean square deviation (MSD) provides quantitative insight into:

- diffusion behavior  
- chain mobility  
- local fluctuations  
- polymer dynamics inside and outside LLPS condensates  

The Einstein MSD is defined as:

.. math::

   \mathrm{MSD}(t) = \langle \lVert \mathbf{r}(t) - \mathbf{r}(0) \rVert^2 \rangle

DROPPS uses MDAnalysis’ ``EinsteinMSD`` class with optional reduced dimensions.

The workflow of ``dps msd``:

1. Load TPR + XTC (+ optional NDX)  
2. Choose an atom group (via ``--selection`` or interactively)  
3. Convert time selection (ns) → frame indices  
4. Apply ``NoJump`` to remove PBC-induced discontinuities  
5. Run the Einstein MSD calculation  
6. Compute lag times (ns)  
7. Write the MSD vs. lag time into an XVG file  


Usage
-----

.. code-block:: bash

   dps msd -s run.tpr -f run.xtc -sel 0 -o msd.xvg

Specify MSD dimensions:

.. code-block:: bash

   dps msd -s run.tpr -f run.xtc -sel 2 -t xy -o msd_xy.xvg


Arguments
---------

Required
--------

.. option:: -s, --run-input TPR

   TPR file containing topology, atom positions, and time-step metadata.

.. option:: -f, --input XTC

   Trajectory from which MSD is computed.

.. option:: -o, --output FILE

   XVG file where the MSD vs time curve will be written.

Optional
--------

Index file
~~~~~~~~~~

.. option:: -n, --index NDX
   Optional NDX file for defining custom atom groups.

Atom group selection
~~~~~~~~~~~~~~~~~~~~

.. option:: -sel, --selection INT

   Specify an index group for MSD calculation.  
   If not given, group selection is performed interactively.

MSD type (dimensions)
~~~~~~~~~~~~~~~~~~~~~

.. option:: -t, --msd-type {xyz, xy, yz, xz, x, y, z}

   Select which coordinate dimensions to include:

   - ``xyz``: full 3D MSD (default)  
   - ``xy``, ``yz``, ``xz``: 2D MSD  
   - ``x``, ``y``, ``z``: 1D MSD

Time selection
~~~~~~~~~~~~~~

.. option:: -b, --start-time INT  
.. option:: -e, --end-time INT  
.. option:: -dt, --delta-time INT

Times are specified in **nanoseconds**.

Converted to frame indices by:

.. code-block:: python

   start_frame, end_frame, interval_frame = trajectory.time2frame(b, e, dt)

Trajectory processing
~~~~~~~~~~~~~~~~~~~~~

NoJump transformation is always applied:

.. code-block:: python

   trajectory.Universe.trajectory.add_transformations(nojump.NoJump())

This ensures atomic trajectories are continuous (molecules are unwrapped)
before MSD is computed.


Program Workflow
----------------

(1) **Load trajectory**

.. code-block:: python

   trajectory = trajectory_class(run_input, index, input)

(2) **Determine frames**

Start, end, and interval frames are computed from ns time input.

(3) **Select atom group**

If ``--selection``:

.. code-block:: python

   selection = trajectory.getSelection(f"group {args.selection}")

Otherwise:

- all groups printed  
- user selects interactively  

(4) **Initialize MSD engine**

.. code-block:: python

   msd_run = MDAnalysis.analysis.msd.EinsteinMSD(selection, msd_type=args.msd_type)

(5) **Run MSD**

.. code-block:: python

   msd_run.run(start=start_frame, stop=end_frame, step=interval_frame)

Results stored in:

.. code-block:: python

   msd_timeseries = msd_run.results.timeseries

(6) **Compute lag times**

.. code-block:: python

   lagtimes = np.arange(msd_run.n_frames) * interval_frame * trajectory.time_step().value_in_unit(nanosecond)

(7) **Write output**

.. code-block:: python

   write_xvg(output_filename, lagtimes, msd_timeseries,
             title, xlabel="Lag time (ns)", ylabel="Mean square deviation")


Output Format
-------------

The XVG file contains:

- Column 1: lag time (ns)  
- Column 2: MSD (nm\*\*2)  

Labels automatically include the selected group ID and MSD type.

Example output structure:

::

   # Mean square deviation of group group0
   @    xaxis  label "Lag time (ns)"
   @    yaxis  label "Mean square deviation"
     0.0   0.000
     0.1   0.023
     0.2   0.041
     ...


Examples
--------

Compute full 3D MSD for group 0:

.. code-block:: bash

   dps msd -s run.tpr -f run.xtc -sel 0 -o msd.xvg

Compute planar (xy) MSD:

.. code-block:: bash

   dps msd -s run.tpr -f run.xtc -sel 1 -t xy -o msd_xy.xvg

Interactive selection:

.. code-block:: bash

   dps msd -s run.tpr -f run.xtc -o msd_mass.xvg


Error Messages
--------------

**“ERROR: No output file specified.”**  
You must specify ``--output``.

**“An exception occurred when trying to open trajectory file …”**  
Trajectory or topology cannot be read.

**Invalid selection**  
Errors may be raised by ``getSelection`` if the group does not exist.


Summary
-------

``dps msd`` provides a robust implementation of Einstein mean square deviation
for molecular trajectories. Features include:

- 1D / 2D / 3D MSD modes  
- PBC-corrected trajectories (``NoJump``)  
- flexible time-window selection  
- clean XVG output for further analysis  

It is an essential tool for studying polymer mobility, diffusive behavior, and
internal chain dynamics in LLPS simulations.