{ "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", "; 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", "# IMD setup - full parameter specification\n", "fix ID group-ID imd trate version 3 unwrap time box coordinates velocities forces \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", "# 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 }