Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions doc/source/examples/other/README.rst
Original file line number Diff line number Diff line change
Expand Up @@ -9,3 +9,4 @@ Other
:maxdepth: 1

parmed_sim
streaming_imd
325 changes: 325 additions & 0 deletions doc/source/examples/other/streaming_imd.ipynb
Original file line number Diff line number Diff line change
@@ -0,0 +1,325 @@
{
Copy link

@ljwoods2 ljwoods2 Oct 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line #2.    # !pip install imdclient>=0.2.2

Don't we want people to install with Conda for consistency?

Also, will people be running this notebook, or reading it as a guide? I think this cell can potentially be removed, seems like a bit of overkill


Reply via ReviewNB

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think we need package installation here. imdclient will always be installed when MDA is installed. imdclient will need to be installed separately. However, I'd mention this fact and then link to installation instructions instead of repeating instructions that may change over time.

In general, always link to authoritative instructions; you're the expert on the matter so you should know what the authoritative sources are.

Copy link

@ljwoods2 ljwoods2 Oct 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why is IMD-coords essential for GROMACS? Will the stream fail out without it? Or do you just mean its essential for this example?

Maybe it's worth specifying exactly what's necessary for this example: you need to have coordinates turned on regardless of simulation engine, you need to have freq set to 1, and you need to turn off coordinate unwrapping (so that there are consistent configs across engines)

Also, should this link the docs for each of the simulation engines?


Reply via ReviewNB

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Link to

Copy link

@ljwoods2 ljwoods2 Oct 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line #42.    # simple_imd_analysis()

Why is this commented out? Will github run it automatically as a part of docstring tests if it's uncommented?


Reply via ReviewNB

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The notebook needs to be primarily a readable document, don't expect people to edit it. Thus, it should show _working code_, so avoid "un-commenting something so that it works".

If there are issues with it being run, comment on the PR and then let's look for a solution.

Copy link

@ljwoods2 ljwoods2 Oct 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line #41.                potential_energy = -1000 + np.random.normal(0, 50)

Not all simulation engines implement energy packet, so maybe remove this?


Reply via ReviewNB

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That's true but for the ones that do, we should show how one can get the information.

Copy link

@ljwoods2 ljwoods2 Oct 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line #51.                    ax1.clear()

Just double checking, did you test this method of live updating plot locally? Looks a bit different than what we did in workshop

https://github.com/MDAnalysis/imd-workshop-2024/blob/main/activity/graph_utils.py


Reply via ReviewNB

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Include links to the authoritative docs — imdclient and IMDReader in the main docs.

The imdclient package should be installed as a dependency of MDAnalysis (in 2.10.0...IIRC(


Reply via ReviewNB

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For all codes:

1) Do we want to send IMD data every integrator step in the example? It's ok to do so but write some text above that says explicitly that we are doing it in the example but in production simulations this may be every 5 or every 10 steps.

For GROMACS:

2) Let's model best practice and encourage everyone to send time and box information so that this is always available in MDAnalysis:

IMD-time         = Yes         ; Send time information
...
IMD-box          = Yes         ; Send box dimensions (recommended)

Reply via ReviewNB

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This section appears to be a duplicate of what's in the imd.rst. Can you just omit all of it and just link?

It's basically impossible to keep documentation in sync so the best you can do is to only have it in a single place and link — DRY, don't repeat yourself.

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Explain what analysis you're performing. Link to https://docs.mdanalysis.org/stable/documentation_pages/core/groups.html#MDAnalysis.core.groups.AtomGroup.radius_of_gyration (or somewhere in the User Guide) for background.

Also explain that **as an example** you only process the first 10 frames but normally you'd read all frames.

You should also discuss

  • skipping (trajectory[::5])
  • limitations (no start/stop)

Reply via ReviewNB

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't like the structure of the code in that the contents of the try block are the whole code. If you want to protect code then make a function out of the inner block. However, I would just not use the try/except and just state that the example requires a running simulation.

(Additionally, I don't like catching all exceptions and then just printing an error — this is bad practice because errors can slip through.)


Reply via ReviewNB

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line #21.            frame_count = 0

In the text, explain that one of the limitations is that one cannot use trajectory[:10] to make clear why you're keeping track of frames. (You can link to StreamReaderBase for more details on limitations).

By the way, couldn't you use ts.frame? That would at least not be quite as clunky.


Reply via ReviewNB

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line #15.        try:

Same comment as above: kind of ugly to have so much code inside try/except.

I'd rewrite the code block inside the try block as a function and then keep the existing try/except structure to demonstrate the overall control structure.


Reply via ReviewNB

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line #41.                potential_energy = -1000 + np.random.normal(0, 50)

Why do we fake it? Instead we should get it from the energies packet. You just need to say above that one needs to enable sending of energies and then pull the correct data out of the ts.data. This would be a good example for showing how to do it, don't squander the opportunity.


Reply via ReviewNB

Copy link
Member

@orbeckst orbeckst Oct 15, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line #96.    # advanced_live_streaming()

No commented out code, see above.


Reply via ReviewNB

"cells": [
{
"cell_type": "markdown",
"id": "d090b35a",
"metadata": {},
"source": [
"# Real-time Streaming Analysis with IMDv3\n",
"\n",
"This tutorial demonstrates how to use MDAnalysis for real-time streaming analysis of molecular dynamics simulations using the Interactive Molecular Dynamics (IMD) v3 protocol. You'll learn how to connect to running simulations and perform live analysis as the simulation progresses.\n",
"\n",
"**Streaming** involves processing data in real-time as it is generated, rather than storing it for later analysis. In molecular dynamics, this means sending simulation data to a client on-the-fly while the simulation is running, without writing large trajectory files to disk.\n",
"\n",
"This is achieved through a TCP/IP socket connection between the simulation engine and receiving client, transmitting coordinates, velocities, forces, energies, and timing information using the IMDv3 protocol.\n",
"\n",
"## What it covers\n",
"\n",
"- How to set up streaming connections to MD engines\n",
"- Real-time monitoring\n",
"- Live analysis workflows\n",
"\n",
"## Prerequisites\n",
"\n",
"Before starting, you'll need:\n",
"- MDAnalysis with IMD support\n",
"- The `imdclient` package (≥ 0.2.2)\n",
"- A running MD simulation with IMD enabled (examples are engine agnostic for the most part)"
]
},
{
"cell_type": "markdown",
"id": "e2168297",
"metadata": {},
"source": [
"## Installation and Setup\n",
"\n",
"The IMDReader requires the `imdclient` package. Let's check if everything is properly installed:"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "6024bb10",
"metadata": {},
"outputs": [],
"source": [
"# Install required packages (uncomment if needed)\n",
"# !pip install imdclient>=0.2.2\n",
"\n",
"import warnings\n",
"warnings.filterwarnings('ignore')\n",
"\n",
"import MDAnalysis as mda\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"from datetime import datetime\n",
"import time\n",
"\n",
"# Check if IMD support is available\n",
"try:\n",
" from MDAnalysis.coordinates.IMD import IMDReader, HAS_IMDCLIENT\n",
" print(f\"IMD support available: {HAS_IMDCLIENT}\")\n",
" if HAS_IMDCLIENT:\n",
" import imdclient\n",
" print(f\"imdclient version: {imdclient.__version__}\")\n",
" print(\"✅ Ready for streaming analysis!\")\n",
" else:\n",
" print(\"❌ IMD support not available\")\n",
"except ImportError as e:\n",
" print(f\"❌ IMD support not available: {e}\")\n",
" print(\"Please install imdclient: pip install imdclient>=0.2.2\")"
]
},
{
"cell_type": "markdown",
"id": "e698ca2d",
"metadata": {},
"source": [
"## Setting Up a Simulation with IMD\n",
"\n",
"Before we can demonstrate streaming analysis, we need a simulation running with IMD enabled. Here are configuration examples for different MD engines:\n",
"\n",
"### GROMACS Setup\n",
"\n",
"Add these comprehensive IMD settings to your `.mdp` file:\n",
"```code\n",
Copy link

Copilot AI Oct 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The code block language identifier should be a valid language (e.g., 'bash', 'text', or 'ini') instead of 'code' for proper syntax highlighting in the documentation.

Copilot uses AI. Check for mistakes.

"; IMD settings for v3 protocol\n",
"IMD-group = System ; Group to stream (typically System)\n",
"IMD-version = 3 ; Use IMDv3 protocol (required for MDAnalysis)\n",
"IMD-nst = 1 ; Frequency of data transmission (every step)\n",
"IMD-time = No ; Send time information\n",
"IMD-coords = Yes ; Send atomic coordinates (essential)\n",
"IMD-vels = No ; Send velocities (optional)\n",
"IMD-forces = No ; Send forces (optional)\n",
"IMD-box = No ; Send box dimensions (optional)\n",
"IMD-unwrap = No ; Unwrap coordinates across PBC\n",
"IMD-energies = No ; Send energy information (optional)\n",
"```\n",
"\n",
"Run the simulation:\n",
"```bash\n",
"gmx mdrun -v -nt 4 -imdwait -imdport 8889\n",
"```\n",
"\n",
"### LAMMPS Setup\n",
"\n",
"Use the comprehensive IMD fix in your input script:\n",
"```code\n",
Copy link

Copilot AI Oct 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The code block language identifier should be a valid language (e.g., 'bash', 'text', or 'lammps') instead of 'code' for proper syntax highlighting in the documentation.

Suggested change
"```code\n",
"```lammps\n",

Copilot uses AI. Check for mistakes.

"# IMD setup - full parameter specification\n",
"fix ID group-ID imd <port> trate <frequency> version 3 unwrap <on/off> time <on/off> box <on/off> coordinates <on/off> velocities <on/off> forces <on/off>\n",
"\n",
"# Example with specific values:\n",
"fix imd all imd 8889 trate 1 version 3 unwrap on time on box on coordinates on velocities on forces on\n",
"```\n",
"\n",
"Run your LAMMPS simulation as usual.\n",
"\n",
"### NAMD Setup\n",
"\n",
"Add comprehensive IMD configuration to your NAMD configuration file:\n",
"```code\n",
Copy link

Copilot AI Oct 10, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The code block language identifier should be a valid language (e.g., 'bash', 'text', or 'tcl') instead of 'code' for proper syntax highlighting in the documentation.

Copilot uses AI. Check for mistakes.

"# IMD Settings\n",
"IMDon yes\n",
"IMDport 8889 # Must match client port\n",
"IMDwait yes # Wait for client connection\n",
"IMDfreq 1 # Frequency of sending data\n",
"\n",
"# Data transmission settings\n",
"IMDsendPositions yes # Send atomic coordinates\n",
"IMDsendEnergies yes # Send energy information\n",
"IMDsendTime yes # Send timing data\n",
"IMDsendBoxDimensions yes # Send simulation box info\n",
"IMDsendVelocities yes # Send atomic velocities\n",
"IMDsendForces yes # Send atomic forces\n",
"IMDwrapPositions no # Don't wrap coordinates\n",
"```\n",
"\n",
"Run your NAMD simulation as usual."
]
},
{
"cell_type": "markdown",
"id": "47c27d17",
"metadata": {},
"source": [
"## Example 1: Simple IMD Analysis\n",
"\n",
"This example shows the basics of connecting to a live simulation and performing simple analysis."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d9181c5c",
"metadata": {},
"outputs": [],
"source": [
"# Simple Example: Basic IMD Connection and Analysis\n",
"\n",
"import MDAnalysis as mda\n",
"import numpy as np\n",
"\n",
"def simple_imd_analysis():\n",
" \"\"\"\n",
" Simple example of connecting to a live simulation and analyzing it.\n",
" \n",
" This assumes you have GROMACS running with IMD enabled on localhost:8889\n",
" \"\"\"\n",
" try:\n",
" # Connect to the live simulation\n",
" u = mda.Universe(\"topol.tpr\", \"imd://localhost:8889\")\n",
" \n",
" print(f\"Connected to simulation: {u.atoms.n_atoms} atoms\")\n",
" \n",
" # Basic analysis on streaming frames\n",
" protein = u.select_atoms(\"protein\")\n",
" \n",
" frame_count = 0\n",
" for ts in u.trajectory:\n",
" # Calculate radius of gyration\n",
" rg = protein.radius_of_gyration()\n",
" \n",
" print(f\"Frame {frame_count}: Time = {ts.time:.1f} ps, Rg = {rg:.2f} Å\")\n",
" \n",
" frame_count += 1\n",
" \n",
" # Stop after 10 frames for this simple example\n",
" if frame_count >= 10:\n",
" break\n",
" \n",
" print(\"Simple analysis completed!\")\n",
" \n",
" except Exception as e:\n",
" print(f\"Connection failed: {e}\")\n",
" print(\"Make sure your simulation is running with IMD enabled\")\n",
" print(\"Example: gmx mdrun -s topol.tpr -imdport 8889 -imdwait\")\n",
"\n",
"# To run this example, uncomment the line below:\n",
"# simple_imd_analysis()"
]
},
{
"cell_type": "markdown",
"id": "43179ebe",
"metadata": {},
"source": [
"## Example 2: Advanced Live Streaming with Real-time Plotting\n",
"\n",
"This example demonstrates continuous monitoring with live plots that update as the simulation runs."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "30bafb63",
"metadata": {},
"outputs": [],
"source": [
"# Advanced Example: Live Streaming Analysis with Real-time Plotting\n",
"\n",
"import MDAnalysis as mda\n",
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"from collections import deque\n",
"import time\n",
"\n",
"def advanced_live_streaming():\n",
" \"\"\"\n",
" Advanced example with live plotting and continuous monitoring.\n",
" \n",
" This connects to your simulation and creates real-time plots.\n",
" \"\"\"\n",
" try:\n",
" # Connect to the simulation\n",
" u = mda.Universe(\"system.tpr\", \"imd://localhost:8889\")\n",
" \n",
" print(f\"Starting live analysis of {u.atoms.n_atoms} atoms\")\n",
" \n",
" # Setup for real-time plotting\n",
" plt.ion() # Interactive mode\n",
" fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(10, 8))\n",
" \n",
" # Data storage for plotting\n",
" times = deque(maxlen=100) # Keep last 100 points\n",
" rg_values = deque(maxlen=100)\n",
" energies = deque(maxlen=100)\n",
" \n",
" protein = u.select_atoms(\"protein\")\n",
" \n",
" frame_count = 0\n",
" for ts in u.trajectory:\n",
" current_time = ts.time\n",
" \n",
" # Calculate properties\n",
" rg = protein.radius_of_gyration()\n",
" \n",
" # Simulate potential energy (real simulations would have this)\n",
" # In real IMD, you might get this from the simulation engine\n",
" potential_energy = -1000 + np.random.normal(0, 50)\n",
" \n",
" # Store data\n",
" times.append(current_time)\n",
" rg_values.append(rg)\n",
" energies.append(potential_energy)\n",
" \n",
" # Update plots every 5 frames\n",
" if frame_count % 5 == 0:\n",
" # Clear and replot\n",
" ax1.clear()\n",
" ax2.clear()\n",
" \n",
" # Plot radius of gyration\n",
" ax1.plot(list(times), list(rg_values), 'b-', linewidth=2)\n",
" ax1.set_xlabel('Time (ps)')\n",
" ax1.set_ylabel('Radius of Gyration (Å)')\n",
" ax1.set_title('Real-time Rg Evolution')\n",
" ax1.grid(True)\n",
" \n",
" # Plot energy\n",
" ax2.plot(list(times), list(energies), 'r-', linewidth=2)\n",
" ax2.set_xlabel('Time (ps)')\n",
" ax2.set_ylabel('Potential Energy (kJ/mol)')\n",
" ax2.set_title('Real-time Energy Evolution')\n",
" ax2.grid(True)\n",
" \n",
" plt.tight_layout()\n",
" plt.draw()\n",
" plt.pause(0.01) # Small pause to update display\n",
" \n",
" # Print status\n",
" print(f\"Frame {frame_count}: Time = {current_time:.1f} ps, \"\n",
" f\"Rg = {rg:.2f} Å, Energy = {potential_energy:.1f} kJ/mol\")\n",
" \n",
" frame_count += 1\n",
" \n",
" # Run for 50 frames in this example\n",
" if frame_count >= 50:\n",
" break\n",
" \n",
" print(\"Live streaming analysis completed!\")\n",
" plt.ioff() # Turn off interactive mode\n",
" plt.show()\n",
" \n",
" except KeyboardInterrupt:\n",
" print(\"\\nAnalysis stopped by user\")\n",
" plt.ioff()\n",
" \n",
" except Exception as e:\n",
" print(f\"Error during live analysis: {e}\")\n",
" print(\"Make sure your simulation is running with IMD enabled\")\n",
" plt.ioff()\n",
"\n",
"# To run this advanced example, uncomment the line below:\n",
"# advanced_live_streaming()"
]
}
],
"metadata": {
"language_info": {
"name": "python"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Loading