From a6851a6a782ac0600bc129eeb786f2eff4e564e8 Mon Sep 17 00:00:00 2001
From: Armin <129873177+Lucky-Armin@users.noreply.github.com>
Date: Tue, 11 Jun 2024 16:33:19 +0200
Subject: [PATCH] Added the missing notebook, because of .gitignore

---
 examples/shake_topology/gen_shake.ipynb | 151 ++++++++++++++++++++++++
 1 file changed, 151 insertions(+)
 create mode 100644 examples/shake_topology/gen_shake.ipynb

diff --git a/examples/shake_topology/gen_shake.ipynb b/examples/shake_topology/gen_shake.ipynb
new file mode 100644
index 0000000..ec61629
--- /dev/null
+++ b/examples/shake_topology/gen_shake.ipynb
@@ -0,0 +1,151 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "In this example a short constraint free NVT simulation of the molecule perylene under vacuum conditions will be used to generate shake constraints. By that the size of the simulation step can be increased up to 2 fs.\n",
+    "\n",
+    "First, the PQAnalysis package is imported and the currently loaded version will be checked."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "'1.1.2'"
+      ]
+     },
+     "execution_count": 5,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "import PQAnalysis\n",
+    "PQAnalysis._version.version"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Then we read in the trajectory from a `.xyz` file generated by a constraint free `PQ` simulation.\n",
+    "\n",
+    "Note: In this example notebook a trajectory consisting of just 500 frames is used due to convenience. For determining shake constraints a larger number of frames is highly recommended."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "application/vnd.jupyter.widget-view+json": {
+       "model_id": "6772016017a94c9999aec38127c7a23e",
+       "version_major": 2,
+       "version_minor": 0
+      },
+      "text/plain": [
+       "  0%|          | 0/500 [00:00<?, ?it/s]"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "data": {
+      "text/plain": [
+       "Trajectory with 500 frames"
+      ]
+     },
+     "execution_count": 6,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "from PQAnalysis.io import TrajectoryReader\n",
+    "\n",
+    "input_file = \"./perylene-md-01.xyz\"\n",
+    "reader = TrajectoryReader(input_file)\n",
+    "\n",
+    "traj = reader.read()\n",
+    "traj"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "For reasons of symmetry the four H-C bonds closer to the molecule center (green) should attain the same length. Equivalently, the outer eight H-C bonds (blue) should also be of equal length. According to this premise, hydrogen atoms are selected from the read trajectory and split into two groups."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 7,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from PQAnalysis.topology import Selection\n",
+    "import numpy as np\n",
+    "\n",
+    "selection = Selection(\"H\")\n",
+    "h_indices = selection.select(traj.topology)\n",
+    "\n",
+    "inner_hydrogen = np.array(h_indices)[[0, 5, 9, 10]]\n",
+    "outer_hydrogen = np.setdiff1d(h_indices, inner_hydrogen)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Lastly, the shake constraints are calculated with averaged bond distances for both H-C groups and printed to an output file `shake.top`. This file can directly be used as input for subsequent simulations.\n",
+    "\n",
+    "Note: To enable constraint dynamics the shake keyword needs to be turned on as described in the documentation."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 8,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from PQAnalysis.topology.shake_topology import ShakeTopologyGenerator\n",
+    "\n",
+    "shake_generator = ShakeTopologyGenerator(\"H\")\n",
+    "shake_generator.generate_topology(traj)\n",
+    "shake_generator.average_equivalents([inner_hydrogen, outer_hydrogen], comments=[\"inner\", \"outer\"])\n",
+    "shake_generator.write_topology(\"shake.top\", mode=\"w\")\n",
+    "\n"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "pq",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.12.3"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}