7.10. IMDReader — MDAnalysis.coordinates.IMD
This module provides support for reading molecular dynamics simulation data via the Interactive Molecular Dynamics (IMD) protocol v3. The IMD protocol allows two-way communicating molecular simulation data through a socket. Via IMD, a simulation engine sends data to a receiver (in this case, the IMDClient) and the receiver can send forces and specific control requests (such as pausing, resuming, or terminating the simulation) back to the simulation engine.
Note
This reader only supports IMDv3, which is implemented in GROMACS, LAMMPS, and NAMD at varying stages of development. See the imdclient simulation engine docs for more. While IMDv2 is widely available in simulation engines, it was designed primarily for visualization and gaps are allowed in the stream (i.e., an inconsistent number of integrator time steps between transmitted coordinate arrays is allowed)
The IMDReader
connects to a simulation via a socket and receives coordinate,
velocity, force, and energy data as the simulation progresses. This allows for real-time
monitoring and analysis of ongoing simulations. It uses the imdclient package
(dependency) to implement the IMDv3 protocol and manage the socket connection and data parsing.
See also
IMDReader
Technical details and parameter options for the reader class
- imdclient documentation
Complete documentation for the IMDClient package
- IMDClient GitHub repository
Source code and development resources
7.10.1. Usage Example
As an example of reading a stream, after configuring GROMACS to run a simulation with IMDv3 enabled (see the imdclient simulation engine docs for up-to-date resources on configuring each simulation engine), use the following commands:
gmx mdrun -v -nt 4 -imdwait -imdport 8889
The IMDReader
can then connect to the running simulation and stream data in real time:
import MDAnalysis as mda
u = mda.Universe("topol.tpr", "imd://localhost:8889", buffer_size=10*1024*1024)
print(" time [ position ] [ velocity ] [ force ] [ box ]")
sel = u.select_atoms("all") # Select all atoms; adjust selection as needed
for ts in u.trajectory:
print(f'{ts.time:8.3f} {sel[0].position} {sel[0].velocity} {sel[0].force} {u.dimensions[0:3]}')
Important
Jupyter Notebook Users: When using IMDReader in Jupyter notebooks, be aware that kernel restarts will not gracefully close active IMD connections. This can leave socket connections open, potentially preventing new connections to the same stream.
Always use try/except/finally
blocks to ensure proper cleanup:
import MDAnalysis as mda
try:
u = mda.Universe("topol.tpr", "imd://localhost:8889")
except Exception as e:
print(f"Error during connection: {e}")
else:
try:
# Your analysis code here
for ts in u.trajectory:
# Process each frame
pass
finally:
# Ensure connection is closed
u.trajectory.close()
Always explicitly call u.trajectory.close()
when finished with analysis to
ensure connection is closed properly.
7.10.2. Important Limitations
Warning
The IMDReader has some important limitations that are inherent in streaming data.
Since IMD streams data in real-time from a running simulation, it has fundamental constraints that differ from traditional trajectory readers:
No random access: Cannot jump to arbitrary frame numbers or seek backwards
Forward-only access: You can only move forward through frames as they arrive
No trajectory length: The total number of frames is unknown until the simulation ends
Single-use iteration: Cannot restart or rewind once the stream has been consumed
No independent copies: Cannot create separate reader instances for the same stream
No stream restart: Cannot reconnect or reopen once the connection is closed
No bulk operations: Cannot extract all data at once using timeseries methods
Limited multiprocessing: Cannot split reader across processes for parallel analysis
Single client connection: Only one reader can connect to an IMD stream at a time
No trajectory Writing: Complimentary IMD Writer class is not available for streaming data
See also
See StreamReaderBase
for technical details.
7.10.3. Multiple Client Connections
The ability to establish multiple simultaneous connections to the same IMD port is MD engine implementation dependent. Some simulation engines may allow multiple clients to connect concurrently, while others may reject or fail additional connection attempts.
See the imdclient simulation engine docs for further details.
Important
Even when multiple connections are supported by the simulation engine, each connection receives its own independent data stream. These streams may contain different data depending on the simulation engine’s configuration, so multiple connections should not be assumed to provide identical data streams.
7.10.4. Classes
- class MDAnalysis.coordinates.IMD.IMDReader(filename, n_atoms=None, buffer_size=10485760, **kwargs)[source]
Coordinate reader implementing the IMDv3 protocol for streaming simulation data.
This class handles the technical aspects of connecting to IMD-enabled simulation engines and processing the incoming data stream. For usage examples and protocol overview, see the module documentation above.
The reader manages socket connections, data buffering, and frame parsing according to the IMDv3 specification. It automatically handles different data packet types (coordinates, velocities, forces, energies, timing) and populates MDAnalysis timestep objects accordingly.
- Parameters:
filename (a string of the form "imd://host:port" where host is the hostname) – or IP address of the listening simulation engine’s IMD server and port is the port number.
n_atoms (int (optional)) – number of atoms in the system. defaults to number of atoms in the topology. Don’t set this unless you know what you’re doing.
buffer_size (int (optional) default=10*(1024**2)) – number of bytes of memory to allocate to the
IMDClient
’s internal buffer. Defaults to 10 megabytes. Larger buffers can improve performance for analyses with periodic heavy computation.**kwargs (dict (optional)) – keyword arguments passed to the constructed
IMDClient
Notes
The IMDReader provides access to additional simulation data through the timestep’s data attribute (ts.data). The following keys may be available depending on what the simulation engine transmits:
- dtfloat
Time step size in picoseconds (from the IMD_TIME packet of the IMDv3 protocol)
- stepint
Current simulation step number (from the IMD_TIME packet of the IMDv3 protocol)
- Energy termsfloat
Various energy components (e.g., ‘potential’, ‘kinetic’, ‘total’, etc.) from the IMD_ENERGIES packet of the IMDv3 protocol.
Note
For important limitations inherent to streaming data, see the module documentation above and
StreamReaderBase
for more technical details.Added in version 2.10.0.
- OtherWriter(filename, **kwargs)
Writer creation is not supported for streaming trajectories.
OtherWriter initialization requires frame-based parameters and trajectory indexing information. Streaming readers process data sequentially without meaningful frame indexing, making writer setup impossible.
- Parameters:
filename (str) – Output filename (ignored, as method is not supported)
**kwargs – Additional keyword arguments (ignored, as method is not supported)
- Raises:
RuntimeError – Always raised, as writer creation is not supported for streaming trajectories
- Writer(filename, **kwargs)
Writer creation is not supported for streaming trajectories.
Writer creation requires trajectory metadata that streaming readers cannot provide due to their sequential processing nature.
- Parameters:
filename (str) – Output filename (ignored, as method is not supported)
**kwargs – Additional keyword arguments (ignored, as method is not supported)
- Raises:
RuntimeError – Always raised, as writer creation is not supported for streaming trajectories
- add_auxiliary(aux_spec: str | Dict[str, str] | None = None, auxdata: str | AuxReader | None = None, format: str | None = None, **kwargs) None
Add auxiliary data to be read alongside trajectory.
Auxiliary data may be any data timeseries from the trajectory additional to that read in by the trajectory reader. auxdata can be an
AuxReader
instance, or the data itself as e.g. a filename; in the latter case an appropriateAuxReader
is guessed from the data/file format. An appropriate format may also be directly provided as a key word argument.On adding, the AuxReader is initially matched to the current timestep of the trajectory, and will be updated when the trajectory timestep changes (through a call to
next()
or jumping timesteps withtrajectory[i]
).The representative value(s) of the auxiliary data for each timestep (as calculated by the
AuxReader
) are stored in the current timestep in thets.aux
namespace under aux_spec; e.g. to add additional pull force data stored in pull-force.xvg:u = MDAnalysis.Universe(PDB, XTC) u.trajectory.add_auxiliary('pull', 'pull-force.xvg')
The representative value for the current timestep may then be accessed as
u.trajectory.ts.aux.pull
oru.trajectory.ts.aux['pull']
.The following applies to energy readers like the
EDRReader
.All data that is present in the (energy) file can be added by omitting aux_spec like so:
u.trajectory.add_auxiliary(auxdata="ener.edr")
aux_spec is expected to be a dictionary that maps the desired attribute name in the
ts.aux
namespace to the precise data to be added as identified by adata_selector
:term_dict = {"temp": "Temperature", "epot": "Potential"} u.trajectory.add_auxiliary(term_dict, "ener.edr")
Adding this data can be useful, for example, to filter trajectory frames based on non-coordinate data like the potential energy of each time step. Trajectory slicing allows working on a subset of frames:
selected_frames = np.array([ts.frame for ts in u.trajectory if ts.aux.epot < some_threshold]) subset = u.trajectory[selected_frames]
See also
Note
Auxiliary data is assumed to be time-ordered, with no duplicates. See the Auxiliary API.
- add_transformations(*transformations)
Add all transformations to be applied to the trajectory.
This function take as list of transformations as an argument. These transformations are functions that will be called by the Reader and given a
Timestep
object as argument, which will be transformed and returned to the Reader. The transformations can be part of thetransformations
module, or created by the user, and are stored as a list transformations. This list can only be modified once, and further calls of this function will raise an exception.u = MDAnalysis.Universe(topology, coordinates) workflow = [some_transform, another_transform, this_transform] u.trajectory.add_transformations(*workflow)
The transformations are applied in the order given in the list transformations, i.e., the first transformation is the first or innermost one to be applied to the
Timestep
. The example above would be equivalent tofor ts in u.trajectory: ts = this_transform(another_transform(some_transform(ts)))
- Parameters:
transform_list (list) – list of all the transformations that will be applied to the coordinates in the order given in the list
See also
- property aux_list
Lists the names of added auxiliary data.
- check_slice_indices(start, stop, step)
Check and validate slice indices for streaming trajectories.
Streaming trajectories have fundamental constraints that differ from traditional trajectory files:
No start/stop indices: Since streams process data continuously without knowing the total length,
start
andstop
must beNone
Step-only slicing: Only the
step
parameter is meaningful, controlling how many frames to skip during iterationForward-only:
step
must be positive (> 0) as streams cannot be processed backward in time
- Parameters:
- Returns:
(start, stop, step) with validated values
- Return type:
- Raises:
ValueError – If
start
orstop
are notNone
, or ifstep
is not a positive integer.
Examples
Valid streaming slices:
traj[:] # All frames (step=None, equivalent to step=1) traj[::2] # Every 2nd frame traj[::10] # Every 10th frame
Invalid streaming slices:
traj[5:] # Cannot specify start index traj[:100] # Cannot specify stop index traj[5:100:2] # Cannot specify start or stop indices traj[::-1] # Cannot go backwards (negative step)
See also
__getitem__
,StreamFrameIteratorSliced
Added in version 2.10.0.
- convert_forces_from_native(force, inplace=True)
Conversion of forces array force from native to base units
- Parameters:
force (array_like) – Forces to transform
inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input force is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Added in version 0.7.7.
- convert_forces_to_native(force, inplace=True)
Conversion of force array force from base to native units.
- Parameters:
force (array_like) – Forces to transform
inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input force is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Added in version 0.7.7.
- convert_pos_from_native(x, inplace=True)
Conversion of coordinate array x from native units to base units.
- Parameters:
x (array_like) – Positions to transform
inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input x is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
- convert_pos_to_native(x, inplace=True)
Conversion of coordinate array x from base units to native units.
- Parameters:
x (array_like) – Positions to transform
inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input x is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
- convert_time_from_native(t, inplace=True)
Convert time t from native units to base units.
- Parameters:
t (array_like) – Time values to transform
inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input t is modified in place and also returned (although note that scalar values t are passed by value in Python and hence an in-place modification has no effect on the caller.) In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
- convert_time_to_native(t, inplace=True)
Convert time t from base units to native units.
- Parameters:
t (array_like) – Time values to transform
inplace (bool, optional) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input t is modified in place and also returned. (Also note that scalar values t are passed by value in Python and hence an in-place modification has no effect on the caller.)
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
- convert_velocities_from_native(v, inplace=True)
Conversion of velocities array v from native to base units
- Parameters:
v (array_like) – Velocities to transform
inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input v is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Added in version 0.7.5.
- convert_velocities_to_native(v, inplace=True)
Conversion of coordinate array v from base to native units
- Parameters:
v (array_like) – Velocities to transform
inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input v is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Added in version 0.7.5.
- copy()
Reader copying is not supported for streaming trajectories.
Streaming readers maintain internal state and connection resources that cannot be duplicated. Each stream connection is unique and cannot be copied.
- Raises:
RuntimeError – Always raised, as copying is not supported for streaming trajectories
- property frame: int
Frame number of the current time step.
This is a simple short cut to
Timestep.frame
.
- get_aux_attribute(auxname, attrname)
Get the value of attrname from the auxiliary auxname
- Parameters:
See also
- get_aux_descriptions(auxnames=None)
Get descriptions to allow reloading the specified auxiliaries.
If no auxnames are provided, defaults to the full list of added auxiliaries.
Passing the resultant description to
add_auxiliary()
will allow recreation of the auxiliary. e.g., to duplicate all auxiliaries into a second trajectory:descriptions = trajectory_1.get_aux_descriptions() for aux in descriptions: trajectory_2.add_auxiliary(**aux)
- Returns:
List of dictionaries of the args/kwargs describing each auxiliary.
- Return type:
- iter_as_aux(auxname)
Iterate through timesteps for which there is at least one assigned step from the auxiliary auxname within the cutoff specified in auxname.
See also
- iter_auxiliary(auxname, start=None, stop=None, step=None, selected=None)
Iterate through the auxiliary auxname independently of the trajectory.
Will iterate over the specified steps of the auxiliary (defaults to all steps). Allows to access all values in an auxiliary, including those out of the time range of the trajectory, without having to also iterate through the trajectory.
After interation, the auxiliary will be repositioned at the current step.
- Parameters:
auxname (str) – Name of the auxiliary to iterate over.
(start (optional) – Options for iterating over a slice of the auxiliary.
stop (optional) – Options for iterating over a slice of the auxiliary.
step) (optional) – Options for iterating over a slice of the auxiliary.
selected (lst | ndarray, optional) – List of steps to iterate over.
- Yields:
AuxStep
object
See also
- property n_frames
Changes as stream is processed unlike other readers
- next()
Advance to the next timestep in the streaming trajectory.
Streaming readers process frames sequentially and cannot rewind once iteration completes. Use
for ts in trajectory
for iteration.- Returns:
The next timestep in the stream
- Return type:
- Raises:
StopIteration – When the stream ends or no more frames are available
- next_as_aux(auxname)
Move to the next timestep for which there is at least one step from the auxiliary auxname within the cutoff specified in auxname.
This allows progression through the trajectory without encountering
NaN
representative values (unless these are specifically part of the auxiliary data).If the auxiliary cutoff is not set, where auxiliary steps are less frequent (
auxiliary.dt > trajectory.dt
), this allows progression at the auxiliary pace (rounded to nearest timestep); while if the auxiliary steps are more frequent, this will work the same as callingnext()
.See the Auxiliary API.
See also
- classmethod parse_n_atoms(filename, **kwargs)
Read the coordinate file and deduce the number of atoms
- Returns:
n_atoms – the number of atoms in the coordinate file
- Return type:
- Raises:
NotImplementedError – when the number of atoms can’t be deduced
- rename_aux(auxname, new)
Change the name of the auxiliary auxname to new.
Provided there is not already an auxiliary named new, the auxiliary name will be changed in ts.aux namespace, the trajectory’s list of added auxiliaries, and in the auxiliary reader itself.
- Parameters:
- Raises:
ValueError – If the name new is already in use by an existing auxiliary.
- rewind()
Rewinding is not supported for streaming trajectories.
Streaming readers process data continuously from streams and cannot restart or go backward in the stream once consumed.
- Raises:
RuntimeError – Always raised, as rewinding is not supported for streaming trajectories
- set_aux_attribute(auxname, attrname, new)
Set the value of attrname in the auxiliary auxname.
- Parameters:
See also
- property time
Time of the current frame in MDAnalysis time units (typically ps).
This is either read straight from the Timestep, or calculated as time =
Timestep.frame
*Timestep.dt
- timeseries(**kwargs)
Timeseries extraction is not supported for streaming trajectories.
Streaming readers cannot randomly access frames or store bulk coordinate data in memory, which
timeseries()
requires. Use sequential frame iteration instead.- Parameters:
**kwargs – Any keyword arguments (ignored, as method is not supported)
- Raises:
RuntimeError – Always raised, as timeseries extraction is not supported for streaming trajectories
- property totaltime: float
Total length of the trajectory
The time is calculated as
(n_frames - 1) * dt
, i.e., we assume that the first frame no time as elapsed. Thus, a trajectory with two frames will be considered to have a length of a single time step dt and a “trajectory” with a single frame will be reported as length 0.
- property transformations
Returns the list of transformations
- units = {'length': None, 'time': None, 'velocity': None}
dict with units of of time and length (and velocity, force, … for formats that support it)