{ "cells": [ { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Determining the persistence length of a polymer\n", "\n", "Here we determine the persistence length of a polymer.\n", "\n", "**Last updated:** December 2022 with MDAnalysis 2.4.0-dev0\n", "\n", "**Minimum version of MDAnalysis:** 1.0.0\n", "\n", "**Packages required:**\n", " \n", "* MDAnalysis (Michaud-Agrawal *et al.*, 2011, Gowers *et al.*, 2016)\n", "* MDAnalysisTests" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:17:36.116304Z", "iopub.status.busy": "2021-05-19T06:17:36.115254Z", "iopub.status.idle": "2021-05-19T06:17:37.222491Z", "shell.execute_reply": "2021-05-19T06:17:37.222862Z" } }, "outputs": [], "source": [ "import MDAnalysis as mda\n", "from MDAnalysis.tests.datafiles import TRZ_psf, TRZ\n", "from MDAnalysis.analysis import polymer\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading files" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The test files we will be working with here feature a pure polymer melt of a polyamide." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:17:37.227353Z", "iopub.status.busy": "2021-05-19T06:17:37.226573Z", "iopub.status.idle": "2021-05-19T06:17:37.381502Z", "shell.execute_reply": "2021-05-19T06:17:37.380891Z" } }, "outputs": [], "source": [ "u = mda.Universe(TRZ_psf, TRZ)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Choosing the chains and backbone atoms" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can define the chains of polyamide to be the common definition of a molecule: where each atom is bonded to at least one other in the group, and not bonded to any atom outside the group. MDAnalysis provides these as `fragments`." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:17:37.425154Z", "iopub.status.busy": "2021-05-19T06:17:37.424424Z", "iopub.status.idle": "2021-05-19T06:17:37.535938Z", "shell.execute_reply": "2021-05-19T06:17:37.536374Z" } }, "outputs": [], "source": [ "chains = u.atoms.fragments" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We then want to select only the backbone atoms for each chain, i.e. only the carbons and nitrogens." ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:17:37.540838Z", "iopub.status.busy": "2021-05-19T06:17:37.540263Z", "iopub.status.idle": "2021-05-19T06:17:37.545444Z", "shell.execute_reply": "2021-05-19T06:17:37.545942Z" } }, "outputs": [], "source": [ "backbones = [ch.select_atoms('not name O* H*') for ch in chains]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This should give us AtomGroups where the spatial arrangement is linear. However, the atoms are in index order. We can use `sort_backbone` to arrange our atom groups into their linear arrangement order." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:17:37.553562Z", "iopub.status.busy": "2021-05-19T06:17:37.552437Z", "iopub.status.idle": "2021-05-19T06:17:40.714834Z", "shell.execute_reply": "2021-05-19T06:17:40.715215Z" } }, "outputs": [], "source": [ "sorted_bb = [polymer.sort_backbone(bb) for bb in backbones]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Calculating the persistence length" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The persistence length is the length at which two points on the polymer chain become decorrelated. This is determined by first measuring the autocorrelation $C(n)$ of two bond vectors $(\\mathbf{a}_i, \\mathbf{a}_{i + n})$ separated by $n$ bonds, where\n", "\n", "$$ C(n) = \\langle \\cos\\theta_{i, i+n} \\rangle = \\langle \\mathbf{a_i} \\cdot \\mathbf{a_{i+n}} \\rangle $$\n", "\n", "An exponential decay is then fitted to this, which yields the persistence length $l_P$ from the average bond length $\\bar{l_B}$.\n", "\n", "$$ C(n) \\approx \\exp\\left( - \\frac{n \\bar{l_B}}{l_P} \\right)$$\n", "\n", "We set up our `PersistenceLength` class ([API docs](https://docs.mdanalysis.org/stable/documentation_pages/analysis/polymer.html)). Note that every chain we pass into it must have the same length." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:17:40.718809Z", "iopub.status.busy": "2021-05-19T06:17:40.718167Z", "iopub.status.idle": "2021-05-19T06:17:40.853094Z", "shell.execute_reply": "2021-05-19T06:17:40.853886Z" } }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plen = polymer.PersistenceLength(sorted_bb)\n", "plen.run()" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "The average bond length is found at `plen.results.lb`, the calculated persistence length at `plen.results.lp`, the measured autocorrelation at `plen.results` and the modelled decorrelation fit at `plen.results.fit`." ] }, { "cell_type": "code", "execution_count": 20, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:17:40.859365Z", "iopub.status.busy": "2021-05-19T06:17:40.858478Z", "iopub.status.idle": "2021-05-19T06:17:40.861400Z", "shell.execute_reply": "2021-05-19T06:17:40.862437Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(179,)\n", "The persistence length is 6.917\n" ] } ], "source": [ "print(plen.results.fit.shape)\n", "print('The persistence length is {:.3f}'.format(plen.results.lp))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`MDAnalysis.analysis.polymer.PersistenceLength` provides a convenience method to plot the results." ] }, { "cell_type": "code", "execution_count": 21, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:17:40.887817Z", "iopub.status.busy": "2021-05-19T06:17:40.885343Z", "iopub.status.idle": "2021-05-19T06:17:41.304662Z", "shell.execute_reply": "2021-05-19T06:17:41.305599Z" } }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 21, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEGCAYAAACKB4k+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAAsTAAALEwEAmpwYAAAn10lEQVR4nO3deXxU5b3H8c9vJgkhENkJyJJEBGSRoEQQuFqLSnG/Wq0L1drW8rJWq96rVS/tbaultdW2dsFSqr22imKtWi1q0aJW64JGRWUXhEBAWWWHbPO7f8wkhJCQkMye7/v1yitzznly5seQzHfOec7zHHN3REREAAKJLkBERJKHQkFERGopFEREpJZCQUREaikURESkVkaiC2iN7t27e0FBQaLLEBFJKe+8885md+/R0LaUDoWCggJKSkoSXYaISEoxs9LGtun0kYiI1FIoiIhILYWCiIjUSuk+BRGRhlRWVlJWVsa+ffsSXUpCZWdn07dvXzIzM5v9MwoFEUk7ZWVl5ObmUlBQgJklupyEcHe2bNlCWVkZhYWFzf65uJw+MrM/mtlGM1vYyHYzs1+b2Qoz+8DMjm/Wjt95BwoKYNasaJYrIilu3759dOvWrc0GAoCZ0a1bt8M+WopXn8IDwKRDbD8DGBj5mgL8rtl7Li2FKVMUDCJygLYcCDVa8hrEJRTc/RVg6yGanAf82cPeBDqbWe+m9ruxQ5fwgz17YOrUaJQqItKmJcvVR32AtXWWyyLrDmJmU8ysxMxKNnbsSnkw0i2yZk3MixQRaa5gMMjIkSMZPnw455xzDtu2bYvq/gsKCti8eTPbtm3j3nvvjdp+kyUUGjrGafDuP+4+092L3b3YzVjVJZId/fvHsDwRSWuzZoX7JwOBqPVTtm/fngULFrBw4UK6du3K9OnTW73PhqRrKJQB/eos9wXWN+cHl/XIh5wcmDYtJoWJSJqbNSvcL1laCu4x6accO3Ys69atA2DlypVMmjSJUaNGcdJJJ7F06VIAHnvsMYYPH05RUREnn3wyAA888ADXXntt7X7OPvtsXn755QP2feutt7Jy5UpGjhzJzTff3Opak+WS1KeBa81sNjAG2O7unzT1Q+bO8gEjYOoUmDw55kWKSBqaOjXcL1lXTT9lFN5XqqurmTdvHl//+tcBmDJlCjNmzGDgwIHMnz+fa665hhdffJHbb7+duXPn0qdPn8M61XTnnXeycOFCFixY0OpaIU6hYGaPAKcA3c2sDPg+kAng7jOAZ4EzgRXAHuCrzdlvVmaQZRdeAZNPiEXZItIWNNYf2cp+yr179zJy5EhWr17NqFGjOP3009m1axevv/46F110UW278vJyAMaPH8+VV17Jl770JS644IJWPXdrxCUU3P3SJrY78K3D3W92ZpBlG3a2uC4REfr3D58yamh9K9T0KWzfvp2zzz6b6dOnc+WVV9K5c+cGP9XPmDGD+fPn88wzzzBy5EgWLFhARkYGoVCotk08RmgnS59Ci2RnBlm7dS+7y6sSXYqIpKpp08L9knVFsZ+yU6dO/PrXv+buu++mffv2FBYW8thjjwHhUcfvv/8+EO5rGDNmDLfffjvdu3dn7dq1FBQUsGDBAkKhEGvXruWtt946aP+5ubns3Bm9D8epHQoZ4fI/2rgrwZWISMqaPBlmzoT8fDALf585M6r9lMcddxxFRUXMnj2bWbNmcf/991NUVMSwYcN46qmnALj55ps59thjGT58OCeffDJFRUWMHz+ewsJCjj32WG666SaOP/7gyR66devG+PHjGT58eFQ6mi185iY1jRh5vO+YdAc/++IIvnRCv6Z/QETahCVLljBkyJBEl5EUGnotzOwddy9uqH1KHylkZQTIzgyoX0FEJEpSOhQABuXlslyhICISFWkRCks/VSiIiERDyofC4LxcNu0sZ+vuikSXIiKS8lI+FAb1ygXQKSQRkShI+VAYnKdQEBGJlpQPhbwj2nFEdgbL1K8gIkmkZursmq/Vq1czbtw4AFavXs3DDz+c4AobliwT4rWYmXFMryN0pCAiSaVmmou6Xn/9dWB/KFx22WUJqOzQUv5IAWBQr44s+3QnqTwQT0TSX8eOHYHwdNevvvoqI0eO5Je//GWCqzpQyh8pQLhfYce+KjbsKKdXp+xElyMiSeSHf1/E4vU7orrPoUcewffPGXbINjWzpAIUFhby5JNP1m678847ufvuu5kzZ05U64qGtAiFQZHO5qWf7lAoiEhSaOj0USpIq1BYvmEnpwzumeBqRCSZNPWJXg6UFn0KXTpk0TO3Hcs+1WypIpL8oj3ddTSlRSgADO6lOZBEJDWMGDGCjIwMioqK1NEcK4PzcnlofinVIScYsESXIyJt3K5dB5+5qFmXmZnJvHnz4l1Ss6TNkcKgXrnsqwyxduuephuLiEiD0iYUaqa70L0VRERaLm1CYWBeeFCIprsQEUCDWWnZa5A2oZCTlUH/rjk6UhARsrOz2bJlS5sOBndny5YtZGcf3tittOlohshd2HSkINLm9e3bl7KyMjZt2pToUhIqOzubvn37HtbPpFUoDO7VkZeXbaS8qpp2GcFElyMiCZKZmUlhYWGiy0hJaXP6CGBwryOoCjmrNu9OdCkiIikpvUKh5goknUISEWmRtAqFwu4dyAiYRjaLiLRQWoVCVkaAo3p00JGCiEgLpVUoQPgKJF2WKiLSMnELBTObZGbLzGyFmd3awPZOZvZ3M3vfzBaZ2Vdb8jyD83JZu3Uvu8urWl+0iEgbE5dQMLMgMB04AxgKXGpmQ+s1+xaw2N2LgFOAn5tZ1uE+1+Be4c7mjzZqGm0RkcMVryOF0cAKd//Y3SuA2cB59do4kGtmBnQEtgKH/XG/JhQ0iE1E5PDFKxT6AGvrLJdF1tX1W2AIsB74ELje3UP1d2RmU8ysxMxKGhqt2K9LDtmZAfUriIi0QLxCoaEbHNSflOQLwALgSGAk8FszO+KgH3Kf6e7F7l7co0ePg3YaCFh4uguFgojIYYtXKJQB/eos9yV8RFDXV4EnPGwFsAo4piVPNigvV5elioi0QLxC4W1goJkVRjqPLwGertdmDXAqgJnlAYOBj1vyZIPzctm4s5zPdle0omQRkbYnLqHg7lXAtcBcYAnwF3dfZGZXm9nVkWZ3AOPM7ENgHnCLu29uyfMN6qUb7oiItETcZkl192eBZ+utm1Hn8XpgYjSe65iaK5A27OTEo7pFY5ciIm1C2o1oBuiZ245O7TPVryAicpjSMhTMjMG6AklE5LClZSgADOrVkWWf7mzTt+MTETlcaRsKg/Ny2bGvig07yhNdiohIykjbUBgUueHOkk93JLgSEZHUkbahMKxPJwIG763ZluhSRERSRtqGQsd2GQzpfQTvlG5NdCkiIikjbUMBoDi/C++t2UZV9UHz6omISAPSOxQKurKnopoln+jSVBGR5kjzUOgCQIlOIYmINEtah0LvTu3p07k9Jas/S3QpIiIpIa1DAcJHCyWlWzWITUSkGdI/FPK7sGFHOWWf7U10KSIiSS/tQ2FUfldA/QoiIs2R9qEwuFcuue0y1K8gItIMaR8KwYBxXH4XhYKISDOkfShAuF9h+cadbN9bmehSRESSWvqGwqxZUFAAgQDFt12DO7y7RkcLIiKHkp6hMGsWTJkCpaXgzsgFrxIMVVMy59VEVyYiktTido/muJo6FfbsqV3MqSxn2IaVlGxOYE0iIikgPY8U1qw5aFVx2WLe75JPRZUmxxMRaUx6hkL//getKl63hH2Z7Vi0fnsCChIRSQ3pGQrTpkFOzgGriresAuCdUnU2i4g0Jj1DYfJkmDkT8vPBDPLz6XnPz+jfNUfjFUREDiE9O5ohHAyTJx+wqvjRBbzy0SbcHTNLUGEiIskrPY8UGjGqoAubd1VQumVP041FRNqgNhUKJxSEJ8d7e7UmxxMRaUibCoWje3TkiOwMdTaLiDSiTYVCIGCMyu9CiUJBRKRBcQsFM5tkZsvMbIWZ3dpIm1PMbIGZLTKzf8WijuKCrqzYuIvPdlfEYvciIiktLqFgZkFgOnAGMBS41MyG1mvTGbgXONfdhwEXxaKW4vwugMYriIg0JF5HCqOBFe7+sbtXALOB8+q1uQx4wt3XALj7xlgUUtSvM5lB0ykkEZEGxCsU+gBr6yyXRdbVNQjoYmYvm9k7ZnZFQzsysylmVmJmJZs2bTrsQrIzgwzv04kSXYEkInKQeIVCQyPFvN5yBjAKOAv4AvA9Mxt00A+5z3T3Yncv7tGjR4uKKc7vwgfrtlNeVd2inxcRSVfxCoUyoF+d5b7A+gba/MPdd7v7ZuAVoCgWxYzK70pFVYiF6zQ5nohIXfEKhbeBgWZWaGZZwCXA0/XaPAWcZGYZZpYDjAGWxKKY4oJwZ/PbmgdJROQAcZn7yN2rzOxaYC4QBP7o7ovM7OrI9hnuvsTM/gF8AISA+9x9YSzq6d6xHYXdO4Qnx/tcLJ5BRCQ1xW1CPHd/Fni23roZ9ZbvAu6KRz2j8rswb8kGTY4nIlJHmxrRXNcJBV34bE8lKzftTnQpIiJJo82Gwqj88OR475Tq0lQRkRptNhQG9OhAl5xM3lqlzmYRkRptNhTMjHFHd+dfyzdSHao/ZEJEpG1qs6EAMHFoHpt3VbBgrY4WRESgjYfC54/pSWbQeH7RhkSXIiKSFNp0KByRncmJR3Xj+cXhS1NFRNq6Nh0KABOH9WLV5t2s3LQr0aWIiCRcmw+F04fkATBXp5BERA4/FMysQ+SmOWmhV6dsivp24vnFCgURkSZDwcwCZnaZmT1jZhuBpcAnkVtm3mVmA2NfZmxNHNaL99du49Pt+xJdiohIQjXnSOElYABwG9DL3fu5e0/gJOBN4E4z+3IMa4y5iUPDp5BeWKKjBRFp25ozId5p7l5Zf6W7bwUeBx43s8yoVxZHR/fsSGH3DryweAOXn5if6HJERBKmySOFmkAws3uskelEGwqNpDdrFhQUQCCAFRYykS28sXIzO/al3j9FRCRaDqejeRfwtJl1ADCziWb2WmzKirFZs2DKFCgtBXcoLeX0GdOorHZeXnb4930WEUkXzb6fgrt/18wuA142s3JgN3BrzCqLpalTYc+eA1Yd9/H7dN+3g+cXfcq5RUcmqDARkcRqdiiY2anANwiHQW/g6+6+LFaFxdSaNQetCnqI05a9wZxOXSmvqqZdRtpcdSsi0myHc/poKvA9dz8FuBB41MwmxKSqWOvfv8HVE7etZFd5FW9+rHssiEjb1OxQcPcJ7v7vyOMPgTOAH8WqsJiaNg1ycg5cl5PDuGsuIycryPOLPk1MXSIiCdacwWuNXXH0CXDqodokrcmTYeZMyM8Hs/D3mTPJvnwynxvUgxcWbyCkeyyISBvUrMFrZnadmR1wzsXMsoCxZvYn4CsxqS6WJk+G1ashFAp/nzwZgInD8ti4s5z3y7YlsjoRkYRoTihMAqqBR8xsvZktNrOPgY+AS4FfuvsDMawxriYMziMYMM2FJCJtUnOuProZ2OPu4yMjl7sDe919W0wrS5BOOZmceFRXXli8gVsmHZPockRE4qo5RwqXA7+D8Mhld//E3beZ2VVmdltsy0uMiUN7sWLjLt1jQUTanOaEwl5339PA+geBlJ4IrzGn10yQp1NIItLGNCsUzKx3/ZXuXg5URb+kxDuyc3uG9zlCl6aKSJvTnFD4OfCUmR0wfaiZ9QRCMakqCUwc2ov31m5j4w7dY0FE2o7mzJL6GDAdeMfM5pjZj8zsx8BrwN2xLjBRJg7Lwx3+uWRjoksREYmbZo1odvc/AYXAX4BMYB9wqbvPimFtCTU4L5f+XXN4YbFOIYlI23E401zsdPc/u/st7n67u5cczhOZ2SQzW2ZmK8ys0dlVzewEM6s2swsPZ//RZmZMHJrHayu2sH2P7rEgIm3D4UyI12JmFiR8CuoMYChwqZkNbaTdT4G58airKecf34eK6hCPvbM20aWIiMRFXEIBGA2scPeP3b0CmA2c10C76wjf4jPxJ/JnzWLYuCJGlS3mwdmvEHoobc+UiYjUilco9AHqftwui6yrZWZ9gPOBGYfakZlNMbMSMyvZtClGd0mrc2e2K96dQ2luD16ZNj28XkQkjcUrFBqaRbX+NKT3ALe4e/WhduTuM9292N2Le/ToEa36DlTnzmxnLHud7rs+48Fhp4fXi4iksXiFQhnQr85yX2B9vTbFwGwzW034Jj73mtl/xqW6+urcmS0rVMWlH8zlxQHFrN1enpByRETiJV6h8DYw0MwKI1NuXwI8XbeBuxe6e4G7FwB/Ba5x97/Fqb4D1bsz22ULniPgzkMnX5yQckRE4iUuoeDuVcC1hK8qWgL8xd0XmdnVZnZ1PGo4LPXuzNZ75xYmfvw2jw4/lX2Vhzy7JSKS0pozdXZUuPuzwLP11jXYqezuV8ajpkZFbrjD1KnhU0n9+3P56cN5blWAv7+/nouK+x3650VEUlS8Th+lnnp3Zhs75UsM7NmRP79Rirtu1Ski6Umh0ExmxhVj8/lw3XYWrN2W6HJERGJCoXAYzj++Lx3bZfDgG6WJLkVEJCYUCoehY7sMLji+D3M++ITNu3R5qoikH4XCYbpibD4V1SEefVvzIYlI+lEoHKaje+YybkA3Hp6/hqrqtL3HkIi0UQqFFrhibD7rtu1l3tLEz9snIhJNCoUWOG1IHr07ZavDWUTSjkKhBTKCASaP6c+/V2xmxbFjIBCAggLNoioiKU+h0EIXr3mbzOoqHsobCe5QWhqeblvBICIpTKHQQj1+8D+ctfRVHj/2VHZnZodX7tmj6bVFJKUpFFpqzRouf/cZdrbrwOyiiQesFxFJVQqFlurfn+PXL2X86gX8ZtwlbG/XoXa9iEiqUii01LRpWE4OU1+8n+3ZHfnNuEvC021Pm5boykREWkyh0FKTJ8PMmQzNCfGlD//Jn0adw6rf3Ld/2m0RkRSkUGiNyPTa//3oT8ls3447MwcmuiIRkVZRKERBz9xsrjllAHMXbeDNj7ckuhwRkRZTKETJVScdxZGdsvnRM4sJhXQTHhFJTQqFKMnODHLLGcewcN0OnnhvXaLLERFpEYVCFJ0z4kiK+nXmrrlL2VNRFR7dXFCgaTBEJGUoFKIoEDC+d9YQNuwoZ+Zv/xae9qK0VNNgiEjKUChEWXFBV846tje/Xx/g00D7AzdqGgwRSXIKhRi49YxjqCbAXSdffvBGTYMhIklMoRAD/brm8NWPXuLxY0/jw7wBB27UNBgiksQUCjHyrS+Optue7fxowlXUXqCqaTBEJMkpFGLkiCsmc+PRmczvfyxzB42F/HyYOVPTYIhIUlMoxNAl3/oig/I6cvuVd7Bt8XIFgogkPYVCDGUEA9x1YRGbdpXz3395XyOdRSTpKRRirKhfZ7571lDmLd3IzFc/PnCjBreJSJKJWyiY2SQzW2ZmK8zs1ga2TzazDyJfr5tZUbxqi7UrxuZz1rG9uWvuMt5atTW8ctYsDW4TkaQTl1AwsyAwHTgDGApcamZD6zVbBXzO3UcAdwAz41FbPJgZd37xWPp1ac91j7zL5l3l4UFse/Yc2FCD20QkweJ1pDAaWOHuH7t7BTAbOK9uA3d/3d0/iyy+CfSNU21xkZudyb2TR/HZnkpumL2A6rVlDTfU4DYRSaB4hUIfYG2d5bLIusZ8HXiuoQ1mNsXMSsysZNOmTVEsMfaGHnkEt587jH+v2MxvJ32j4UYa3CYiCRSvULAG1jV4KY6ZfZ5wKNzS0HZ3n+nuxe5e3KNHjyiWGB8Xn9CPC47rwz3HnsVrg0YfuFGD20QkweIVCmVAvzrLfYH19RuZ2QjgPuA8d0/LW5iZGT86fzhH98zl+ou+y4bBI8BMg9tEJCnEKxTeBgaaWaGZZQGXAE/XbWBm/YEngMvdfXmc6kqInKwM7p18PLstk+tunEFVZRWsXq1AEJGEi0souHsVcC0wF1gC/MXdF5nZ1WZ2daTZ/wLdgHvNbIGZlcSjtkQZmJfLjy8YzlurtvKLFw6RgRrLICJxlBGvJ3L3Z4Fn662bUefxVcBV8aonGZx/XF/eWrWVe19eSd8uOVw2pl4nc81YhppLV2vGMoCOKkQkJjSiOcG+f84wPj+4B//z5Ic8+MbqAzdqLIOIxJlCIcGyM4PMuHwUpw3J43tPLeL/Xlu1f2NjYxY0lkFEYkShkATaZQS5d/LxfGFYHj/8+2L+8EpkjqTGxixoLIOIxIhCIUlkZQT47WXHc9axvZn27BLufXlFeMxCTs6BDTWWQURiSKGQRDKDAX51yUjOLTqSn/1jGb/pNTo8diE/v/GxDLo6SUSiKG5XH0nzZAQD/PLikWQEjJ+/sJyqU0/ghlWrMGtgULiuThKRKNORQhIKBoy7LiriolF9+dW8j7j7+WW4NzAriK5OEpEo05FCkgoGjJ9+cQQZwQDTX1rJll0VfP+cYbTPCu5vpKuTRCTKdKSQxAIBY9p/DueaUwYw++21nDf93yz7dOf+Bro6SUSiTKGQ5AIB4zuTjuHPXxvN1t2VnPvbf/PQm6Xh00nNvTpJndEi0kwKhRRx8qAePHf9SYw5qhvf/dtCrpn1LtvP/1Lzrk7SbT9FpJmswQ7MFFFcXOwlJWk9b95BQiHnvn9/zM/+sYy8I7L51SUjKS7o2vgPFBSEg6C+/PzwzKwi0uaY2TvuXtzQNh0ppJhAwJhy8gAe/+Y4ggHj4plv8pt5H1EdaiTc1RktIodBoZCiivp15plv/wdnj+jNz19YziUz3+D9tdsObqjOaBE5DAqFFJabnck9F4/k5xcV8fGm3Zw3/TWuffhdSrfs3t+oOZ3R6ogWkQj1KaSJnfsq+cMrH/OHV1dRFQoxeUw+1004mm4d24Xf5KdODZ8y6t8/HAg1ndH1R0VDODR0a1CRtHWoPgWFQprZuGMf98z7iEffXkv7zCBTTj6Kq04qJCerkXGK6ogWaXMUCm3Qio27uGvuUuYu2kCP3HbccNpALhzVl3YZwQMbBgLhS1XrM4NQKD7Fikhc6eqjNujonh35/eXFPP7NseR3zWHqkws58cfz+PGzS1i1uU6fQ3M6otXnINJm6EihDXB3Xl+5hVnzS3l+0QaqQs64Ad24bEx/Jr7/EllXH6JPQX0OImlHp4+k1sad+3ispIyH569h3ba9dO+YxUXtd3DpH39M/yXvHdwRrT4HkbSjUJCDVIecVz7axMPz1zBvyQYcOCG/K6cO6cmpQ/IY0KND+B4O6nMQSTvqU5CDBAPG5wf35A9XFPParRO44dRB7Cyv4ifPLeW0X/yLz9/9MnfMWczrxadSGQgevIP6fRHqdxBJCwoFoXen9lx/2kCeu/4kXrt1AnecN4z8bh148I1SLptwA8d/+2GuPfc7PD5sAms65eENDX7TpHuSKE19IGnt9mjtI1W4e8p+jRo1yiV2du2r9H8s/MRvvvNxH3X9w55/yxzPv2WOj5r6tF/1p7f93pdW+JsrN/ueo452D8fBgV/5+ft39tBD4WWz8PeHHkrQv0pSSlO/Nw895J6Tc+DvXU7O/nat3R6tfTTn3xJHQIk38r6a8Df21nwpFOKnujrki9Zt9wffWO03Pvqen3LXS7UhMeCmv/k5V/zCv3/qFJ894nQvOfIY356VE/7ld2/+H42EJdGbxyHFus7m/N7k5x/6A0lrt0drH0n2N6BQkJjYvHOfv7DoU//pOdf5xZf82I+58a+1QZF/yxwfc92D/uX73vQf/ud/+SMjJnpJn2N8a3auhxr6o2lLDvVmGq1Pna3dHq06W/NaNOfN1qzhNjUfSFq7PVr7aG5wxPrDQOQ5RoG7QkFiJvIGUWUBL+2U5y8MGO33/selfuNPnvCzf/3qQWEx5MbH/PSvTfcrL/y+f/fJD33Gyyv877991N877mTf0LGLVxUUJO+n49Zq6s00Gp8643HKpDl11uynpcESjTfbZDlSaOrf0tz/kygFeVKEAjAJWAasAG5tYLsBv45s/wA4vql9KhSSyCF+GavzC3xNpzz/54AT/L7ic/2HE67yb5w/1c+Y8jsf8YO5BwRG/i1zvPDmp3zUtQ/5pO8/5V++702/8dH3/MfPLPY/vLLSn3y3zF9ausHfW/OZr9q0y7fuKveq6lCra2x2m9Z+mmvqDSRV3gibU2drgyVVAjIep7miHOSHCoW4jFMwsyCwHDgdKAPeBi5198V12pwJXAecCYwBfuXuYw61X41TSBFNjIreMfAY1m0vZ90RPVnXqSebO3Rmc05nNvU4kk0nnszmneVs2lVORVXj4yJyszPoXF1O5w1l5O7cRk6GkTP0GHIGDaB9VpCc5UvJeeoJcvbsJKdyHzkV+8gJOO2/fS0dzvwCOVlB2j/7d3Ju+i9ytm+lXXUlVq/OZo/uPtSstE2N+2jOYMGm9tHa7dGqs6k2TT1HNF7vaGyP1nMc6t8S59e7GChxt4Mbx2nwmpmNBX7g7l+ILN8G4O4/qdPm98DL7v5IZHkZcIq7f9LYfhUKKaQ1b5SAPzSLHdfdwOZAO7a1z2V7u45s69SN7V/5OtuGFrF9wUK2vzafbZnt2dkuhz2Z2ezNas/uHr3YG8xiz95yqhsab3EImdWVZFVVkkmIzO7dyPp0PVnl+8isrgpvq/melUHWKZ8jM2hklq0la/4bZJbvIxgKEfRqgsEggQmfJzh0CIGZvye4bRtBDxEIhQh6iGComkCnTmT8z20E3nuX4OxHCJbvq90eyMwgeOVXCI4fTyBgBK/9FsGNGwiEQgTcMRxzx3r2hAf+D/vKldjGDeF1NdsgvP2Rh7FLL4VPNxy4zR3rlYc9/lfAsAvOxz75ZP+2mra9e8OcOdhzz2F33I7t3bt/H9ntsB/8ADvn7PA+hg7BQqHafRzw/7p8OZxyCqxbF15Fne19+sC//hV+/NTT8Iufw/r14ee+6SY497wG/7+swbe4w9PUPqwZT9JoiyefhLt+hpWtC/8bb/kOnH9+eNuJY2FdWfjn6/4p9O0Db86Hfn2p+xLVvl5mUFYWadu39u/ooNd73ToYPbq27Zm7P0t4KFwITHL3qyLLlwNj3P3aOm3mAHe6+78jy/OAW9y90Xd9hUKaiManoCa2eyBAeSCDvZnZ7MnMZk9Wne//eIE9FdXsveJKdke2l2dkUhnIpDKYQUVGJpVXf5OK+x+gMhikIphJZTCTimAGlYHI9tEnUlkdomL5CipCTkUwk1AgQLUFqA4ECQWDVOd0pLqqilBVNdWBAG4aJiSJ0f2nZzcaCo1Msh91DT15/TRqThvMbAowBaC/bimZHqZNa/jQuu4AuabuNd3Eduvfn+zSUrKrK+myb+f+7fn5MCQv/Hj78saD5YLfwX890fj2B28KPw6c0vSpmchRk69ZQ3V+AdW330Ho4kuoCoUIhaDaneqQE4p8r/s4/B2q58whNP1eQhs24L164d+8Bv/CF3AiTz93Lj5zJr5hI94rD//GN/DTTg+fUMbhn/Pg/vvxjRvxvDz8q1/DJ0yI/LyH//BefAn/85/xzZvxHj3wyy/HT/7c/ufAI/ujdr/7lx1eew2/7368omL/69CuHXztazB+XHj5tdfxv/4VNm+B7t3gwgth3LiDXr6mPrs256NtUx+Am9xHM57Em2h0yBLmz8effBK2boWuXcNHEWPG1G7jwQehomL/M2RlweWXw+hIm7fCbbyicv8+s7Lgy1+GMaMj+3kLf+op7jl0kXHpZB4LzK2zfBtwW702vyfcz1CzvAzofaj9qqM5jTTVgZsMHXXR6FBsa1JlzEUqiMaFEhEkepwC4SOSj4FCIAt4HxhWr81ZwHOEjxhOBN5qar8KhTYkWleRxPrqoyQbpCTSkISHQrgGziR8BdJKYGpk3dXA1ZHHBkyPbP8QKG5qnwqFNibWl4tGS7LUIdKIQ4WCps4WEWljNHW2iIg0i0JBRERqKRRERKSWQkFERGopFEREpJZCQUREaikURESklkJBRERqKRRERKSWQkFERGopFEREpJZCQUREaikURESkVkrPkmpmOwnfjCeZdQc2J7qIJqjG6FCNrZfs9UF61Jjv7j0a2hCv23HGyrLGpn9NFmZWohpbTzVGR7LXmOz1QfrXqNNHIiJSS6EgIiK1Uj0UZia6gGZQjdGhGqMj2WtM9vogzWtM6Y5mERGJrlQ/UhARkShSKIiISK2UDQUzm2Rmy8xshZndmuh6AMzsj2a20cwW1lnX1cxeMLOPIt+7JLC+fmb2kpktMbNFZnZ9EtaYbWZvmdn7kRp/mGw11qk1aGbvmdmcZKzRzFab2YdmtsDMSpK0xs5m9lczWxr5vRybTDWa2eDI61fztcPMbkiyGm+M/K0sNLNHIn9DLa4vJUPBzILAdOAMYChwqZkNTWxVADwATKq37lZgnrsPBOZFlhOlCvhvdx8CnAh8K/K6JVON5cAEdy8CRgKTzOxEkqvGGtcDS+osJ2ONn3f3kXWuWU+2Gn8F/MPdjwGKCL+eSVOjuy+LvH4jgVHAHuDJZKnRzPoA3waK3X04EAQuaVV97p5yX8BYYG6d5duA2xJdV6SWAmBhneVlQO/I496EB9wlvM5IPU8BpydrjUAO8C4wJtlqBPpG/tgmAHOS8f8aWA10r7cuaWoEjgBWEbngJRlrrFfXROC1ZKoR6AOsBboSHow8J1Jni+tLySMF9r8QNcoi65JRnrt/AhD53jPB9QBgZgXAccB8kqzGyGmZBcBG4AV3T7oagXuA7wChOuuSrUYHnjezd8xsSmRdMtV4FLAJ+L/Iabj7zKxDktVY1yXAI5HHSVGju68D7gbWAJ8A2939+dbUl6qhYA2s07W1zWRmHYHHgRvcfUei66nP3as9fLjeFxhtZsMTXNIBzOxsYKO7v5PoWpow3t2PJ3ya9VtmdnKiC6onAzge+J27HwfsJvGnsxpkZlnAucBjia6lrkhfwXlAIXAk0MHMvtyafaZqKJQB/eos9wXWJ6iWpmwws94Ake8bE1mMmWUSDoRZ7v5EZHVS1VjD3bcBLxPup0mmGscD55rZamA2MMHMHiK5asTd10e+byR8Hnw0yVVjGVAWORIE+CvhkEimGmucAbzr7hsiy8lS42nAKnff5O6VwBPAuNbUl6qh8DYw0MwKIwl+CfB0gmtqzNPAVyKPv0L4PH5CmJkB9wNL3P0XdTYlU409zKxz5HF7wr/0S0miGt39Nnfv6+4FhH/3XnT3L5NENZpZBzPLrXlM+DzzQpKoRnf/FFhrZoMjq04FFpNENdZxKftPHUHy1LgGONHMciJ/36cS7qxveX2J7rxpRQfLmcByYCUwNdH1RGp6hPB5vUrCn4K+DnQj3CH5UeR71wTW9x+ET7N9ACyIfJ2ZZDWOAN6L1LgQ+N/I+qSpsV69p7C/ozlpaiR8vv79yNeimr+RZKoxUs9IoCTy//03oEsS1pgDbAE61VmXNDUCPyT8wWkh8CDQrjX1aZoLERGplaqnj0REJAYUCiIiUkuhICIitRQKIiJSS6EgIiK1FAoiIlJLoSAiIrUUCiJRZGYnmNkHkTntO0TmuU+quZtEDkWD10SizMx+BGQD7QnP7fOTBJck0mwKBZEoi8zH9TawDxjn7tUJLkmk2XT6SCT6ugIdgVzCRwwiKUNHCiJRZmZPE55Su5Dw3a+uTXBJIs2WkegCRNKJmV0BVLn7w5F7ib9uZhPc/cVE1ybSHDpSEBGRWupTEBGRWgoFERGppVAQEZFaCgUREamlUBARkVoKBRERqaVQEBGRWv8PCAbeKAA5Cj4AAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plen.plot()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "[1] Richard J. Gowers, Max Linke, Jonathan Barnoud, Tyler J. E. Reddy, Manuel N. Melo, Sean L. Seyler, Jan Domański, David L. Dotson, Sébastien Buchoux, Ian M. Kenney, and Oliver Beckstein.\n", "MDAnalysis: A Python Package for the Rapid Analysis of Molecular Dynamics Simulations.\n", "Proceedings of the 15th Python in Science Conference, pages 98–105, 2016.\n", "00152.\n", "URL: https://conference.scipy.org/proceedings/scipy2016/oliver_beckstein.html, doi:10.25080/Majora-629e541a-00e.\n", "\n", "[2] Naveen Michaud-Agrawal, Elizabeth J. Denning, Thomas B. Woolf, and Oliver Beckstein.\n", "MDAnalysis: A toolkit for the analysis of molecular dynamics simulations.\n", "Journal of Computational Chemistry, 32(10):2319–2327, July 2011.\n", "00778.\n", "URL: http://doi.wiley.com/10.1002/jcc.21787, doi:10.1002/jcc.21787." ] } ], "metadata": { "kernelspec": { "display_name": "guide", "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.9.15 | packaged by conda-forge | (main, Nov 22 2022, 15:55:03) \n[GCC 10.4.0]" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": false, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false }, "vscode": { "interpreter": { "hash": "5edc5d8d8cbc0935a054a8e44024f729bc376180aae27775d15f2ff38c68f892" } } }, "nbformat": 4, "nbformat_minor": 2 }