{ "cells": [ { "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 executed:** Feb 06, 2020 with MDAnalysis 0.20.1\n", "\n", "**Last updated:** January 2020\n", "\n", "**Minimum version of MDAnalysis:** 0.20.1\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": {}, "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": {}, "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": {}, "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": {}, "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": {}, "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. Note that every chain we pass into it must have the same length." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plen = polymer.PersistenceLength(sorted_bb)\n", "plen.run()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The average bond length is found at `plen.lb`, the calculated persistence length at `plen.lp`, the measured autocorrelation at `plen.results` and the modelled decorrelation fit at `plen.fit`." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(179,)\n", "The persistence length is 6.917464580166461\n" ] } ], "source": [ "print(plen.results.shape)\n", "print('The persistence length is {}'.format(plen.lp))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`MDAnalysis.analysis.polymer.PersistenceLength` provides a convenience method to plot the results." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEGCAYAAACKB4k+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deXxU5b3H8c9vJgkhEJE1IEuCCwgoawShV0tREavitbUVRe2i5VLreq+tcrG1raUbXeyCtdT22iqKbV2LWrRUqnUPisoOAoGAEhbZIdv87h8zCSEkZEgms+X7fr3yypxznpz5MSTznXOe8zzH3B0RERGAQKILEBGR5KFQEBGRGgoFERGpoVAQEZEaCgUREamRkegCmqNLly5eUFCQ6DJERFLKokWLtrl71/q2pXQoFBQUUFRUlOgyRERSipkVN7RNp49ERKSGQkFERGooFEREpEZK9ymIiNSnoqKCkpISDh48mOhSEio7O5tevXqRmZkZ9c8oFEQk7ZSUlJCbm0tBQQFmluhyEsLd2b59OyUlJfTt2zfqn4vL6SMz+4OZlZrZkga2m5n90szWmNl7ZjY8qh0vWgQFBTBnTizLFZEUd/DgQTp37txqAwHAzOjcufMxHy3Fq0/hAWDCUbZfAJwS+ZoC/CbqPRcXw5QpCgYROUxrDoRqTXkN4hIK7v4SsOMoTS4B/uRhrwPHm1mPxvZb2q5j+MH+/TB9eixKFRFp1ZLl6qOewMZayyWRdUcwsylmVmRmRaXtO1EWjHSLbNjQ4kWKiEQrGAwydOhQTjvtNC6++GJ27twZ0/0/8MAD3HDDDQA8+eSTLFu2LCb7TZZQqO8Yp967/7j7bHcvdPdCN2Ndx0h29OnTguWJSFqbMyfcPxkIxKyfsm3btixevJglS5bQqVMnZs2a1ex9NiQdQ6EE6F1ruRewOZofXNk1H3JyYMaMFilMRNLcnDnhfsniYnBvkX7K0aNHs2nTpprlmTNncsYZZzB48GDuuusuAPbt28eFF17IkCFDOO2003j00UeB8HQ+27ZtA6CoqIixY8cetu9XX32Vp59+mq9//esMHTqUDz74oFm1JsslqU8DN5jZXGAUsMvdP2zsh8ydVScNhulTYPLkFi9SRNLQ9OnhfsnaqvspY/C+UlVVxYIFC7j22msBeP7551m9ejVvvvkm7s7EiRN56aWX2Lp1KyeccALPPPMMALt27Ypq/2PGjGHixIlcdNFFXHbZZc2uN16XpD4CvAb0N7MSM7vWzKaa2dRIk2eBtcAa4HfA9dHsNyszyMrLrlEgiEjTNdQf2cx+ygMHDjB06FA6d+7Mjh07OO+884BwKDz//PMMGzaM4cOHs2LFClavXs3pp5/OP/7xD26//XZefvllOnTo0Kznb6q4HCm4+xWNbHfga8e63+zMICu37GlyXSIi9OkTPmVU3/pmqO5T2LVrFxdddBGzZs3ipptuwt2ZNm0a//Vf/3XEzyxatIhnn32WadOmMX78eL71rW+RkZFBKBQCiMsI7WTpU2iS7MwgG3ccYF9ZZaJLEZFUNWNGuF+ythj2U3bo0IFf/vKX/OQnP6GiooLzzz+fP/zhD+zduxeATZs2UVpayubNm8nJyeGqq67itttu4+233wbCfQqLFi0C4LHHHqv3OXJzc9mzJzYfkFM7FDLC5a8u3ZvgSkQkZU2eDLNnQ34+mIW/z54d09PSw4YNY8iQIcydO5fx48dz5ZVXMnr0aE4//XQuu+wy9uzZw/vvv8/IkSMZOnQoM2bM4M477wTgrrvu4uabb+ass84iGAzWu/9JkyYxc+ZMhg0b1uyOZgufuUlNg4cO990T7ubHnx3M58/o3fgPiEirsHz5cgYMGJDoMpJCfa+FmS1y98L62qf0kUJWRoDszID6FUREYiSlQwGgX14uqxQKIiIxkRahsOIjhYKISCykfCj0z8tl654yduwrT3QpIiIpL+VDoV/3XACdQhIRiYGUD4X+eQoFEZFYSflQyDuuDcdlZ7BS/QoikkSqp86u/lq/fj1FRUXcdNNNACxcuJBXX301wVUeKVkmxGsyM+PU7sfpSEFEkkr1NBe1FRQUUFgYHh6wcOFC2rdvz5gxYxJRXoNS/kgBoF/39qz8aA+pPBBPRNLfwoULueiii1i/fj333XcfP//5zxk6dCgvv/xyokurkfJHChDuV9h9sJItu8vo3iE70eWISBL5zt+Wsmzz7pjuc+AJx3HXxYOO2qZ6llSAvn378sQTT9RsKygoYOrUqbRv357bbrstprU1V1qEQr9IZ/OKj3YrFEQkKdR3+igVpFUorNqyh7H9uyW4GhFJJo19opfDpUWfQsd2WXTLbcPKjzRbqoikhlhOdx1LaREKAP27aw4kEUkdF198MU888YQ6mltK/7xcHnqjmKqQEwxYossRkVau+iY6tY0dO5axY8cC0K9fP9577704V9W4tDlS6Nc9l4MVITbu2N94YxERqVfahEL1dBe6t4KISNOlTSicktceQNNdiAiABrPStNcgbUIhJyuDPp1ydKQgImRnZ7N9+/ZWHQzuzvbt28nOPraxW2nT0QyRu7DpSEGk1evVqxclJSVs3bo10aUkVHZ2Nr169Tqmn0mrUOjfvT0LV5ZSVllFm4xgossRkQTJzMykb9++iS4jJaXN6SOA/t2PozLkrNu2L9GliIikpPQKheorkHQKSUSkSdIqFPp2aUdGwDSyWUSkidIqFLIyApzYtZ2OFEREmiitQgHCVyDpslQRkaaJWyiY2QQzW2lma8zsjnq2dzCzv5nZu2a21My+1JTn6Z+Xy8YdB9hXVtn8okVEWpm4hIKZBYFZwAXAQOAKMxtYp9nXgGXuPgQYC/zUzLKO9bn6dw93Nq8u1TTaIiLHKl5HCiOBNe6+1t3LgbnAJXXaOJBrZga0B3YAx/xxvzoUNIhNROTYxSsUegIbay2XRNbV9mtgALAZeB+42d1DdXdkZlPMrMjMiuobrdi7Yw7ZmQH1K4iINEG8QqG+GxzUnZTkfGAxcAIwFPi1mR13xA+5z3b3Qncv7Nq16xE7DQQsPN2FQkFE5JjFKxRKgN61lnsRPiKo7UvA4x62BlgHnNqUJ+uXl6vLUkVEmiBeofAWcIqZ9Y10Hk8Cnq7TZgNwDoCZ5QH9gbVNebL+ebmU7inj433lzShZRKT1iUsouHslcAMwH1gO/Nndl5rZVDObGml2NzDGzN4HFgC3u/u2pjxfv+664Y6ISFPEbZZUd38WeLbOuvtqPd4MjI/Fc51afQXSlj2ceWLnWOxSRKRVSLsRzQDdctvQoW2m+hVERI5RWoaCmdFfVyCJiByztAwFgH7d27Pyoz2t+nZ8IiLHKm1DoX9eLrsPVrJld1miSxERSRlpGwr9IjfcWf7R7gRXIiKSOtI2FAb17EDA4J0NOxNdiohIykjbUGjfJoMBPY5jUfGORJciIpIy0jYUAArzO/LOhp1UVh0xr56IiNQjvUOhoBP7y6tY/qEuTRURiUaah0JHAIp0CklEJCppHQo9OrSl5/FtKVr/caJLERFJCWkdChA+Wigq3qFBbCIiUUj/UMjvyJbdZZR8fCDRpYiIJL20D4UR+Z0A9SuIiEQj7UOhf/dccttkqF9BRCQKaR8KwYAxLL+jQkFEJAppHwoQ7ldYVbqHXQcqEl2KiEhSS99QmDMHCgogEKBw2vW4w9sbdLQgInI06RkKc+bAlClQXAzuDF38MsFQFUXzXk50ZSIiSS1u92iOq+nTYf/+msWcijIGbfmAom0JrElEJAWk55HChg1HrCosWca7HfMpr9TkeCIiDUnPUOjT54hVhZuWczCzDUs370pAQSIiqSE9Q2HGDMjJOWxV4fZ1ACwqVmeziEhD0jMUJk+G2bMhPx/MID+fbvf8mD6dcjReQUTkKNKzoxnCwTB58mGrCh9dzEurt+LumFmCChMRSV7peaTQgBEFHdm2t5zi7fsbbywi0gq1qlA4oyA8Od5b6zU5nohIfVpVKJzctT3HZWeos1lEpAGtKhQCAWNEfkeKFAoiIvWKWyiY2QQzW2lma8zsjgbajDWzxWa21Mz+1RJ1FBZ0Yk3pXj7eV94SuxcRSWlxCQUzCwKzgAuAgcAVZjawTpvjgXuBie4+CPhcS9RSmN8R0HgFEZH6xOtIYSSwxt3Xuns5MBe4pE6bK4HH3X0DgLuXtkQhQ3ofT2bQdApJRKQe8QqFnsDGWsslkXW19QM6mtlCM1tkZtfUtyMzm2JmRWZWtHXr1mMuJDszyGk9O1CkK5BERI4Qr1Cob6SY11nOAEYAFwLnA980s35H/JD7bHcvdPfCrl27NqmYwvyOvLdpF2WVVU36eRGRdBWvUCgBetda7gVsrqfN3919n7tvA14ChrREMSPyO1FeGWLJJk2OJyJSW7xC4S3gFDPra2ZZwCTg6TptngLOMrMMM8sBRgHLW6KYwoJwZ/NbmgdJROQwcZn7yN0rzewGYD4QBP7g7kvNbGpk+33uvtzM/g68B4SA+919SUvU06V9G/p2aReeHO+TLfEMIiKpKW4T4rn7s8CzddbdV2d5JjAzHvWMyO/IguVbNDmeiEgtrWpEc21nFHTk4/0VfLB1X6JLERFJGq02FEbkhyfHW1SsS1NFRKq12lA4qWs7OuZk8uY6dTaLiFRrtaFgZow5uQv/WlVKVajukAkRkdap1YYCwPiBeWzbW87ijTpaEBGBVh4Knzq1G5lB4/mlWxJdiohIUmjVoXBcdiZnntiZ55eFL00VEWntWnUoAIwf1J112/bxwda9iS5FRCThWn0onDcgD4D5OoUkInLsoWBm7SI3zUkL3TtkM6RXB55fplAQEWk0FMwsYGZXmtkzZlYKrAA+jNwyc6aZndLyZbas8YO68+7GnXy062CiSxERSahojhReBE4CpgHd3b23u3cDzgJeB35oZle1YI0tbvzA8CmkF5braEFEWrdoJsQ7190r6q509x3AY8BjZpYZ88ri6ORu7enbpR0vLNvC1WfmJ7ocEZGEafRIoToQzOwea2A60fpCI+nNmQMFBRAIYH37Mp7tvPbBNnYfTL1/iohIrBxLR/Ne4GkzawdgZuPN7JWWKauFzZkDU6ZAcTG4Q3Ex5903g4oqZ+HKY7/vs4hIuoj6fgrufqeZXQksNLMyYB9wR4tV1pKmT4f9+w9bNWztu3Q5uJvnl37ExCEnJKgwEZHEijoUzOwc4CuEw6AHcK27r2ypwlrUhg1HrAp6iHNXvsa8Dp0oq6yiTUbaXHUrIhK1Yzl9NB34pruPBS4DHjWzcS1SVUvr06fe1eN3fsDeskpeX6t7LIhI6xR1KLj7OHf/d+Tx+8AFwPdaqrAWNWMG5OQcvi4nhzHXX0lOVpDnl36UmLpERBIsmsFrDV1x9CFwztHaJK3Jk2H2bMjPB7Pw99mzyb56Mp/s15UXlm0hpHssiEgrFNXgNTO70cwOO+diZlnAaDP7I/CFFqmuJU2eDOvXQygU/j55MgDjB+VRuqeMd0t2JrQ8EZFEiCYUJgBVwCNmttnMlpnZWmA1cAXwc3d/oAVrjKtx/fMIBkxzIYlIqxTN1UdfB/a7+yciI5e7AAfcPS0/SnfIyeTMEzvxwrIt3D7h1ESXIyISV9EcKVwN/AbCI5fd/UN332lm15nZtJYtLzHGD+zOmtK9useCiLQ60YTCAXffX8/6B4GUngivIedVT5CnU0gi0spEFQpm1qPuSncvAypjX1LinXB8W07reZwuTRWRVieaUPgp8JSZHTZ9qJl1A0ItUlUSGD+wO+9s3Enpbt1jQURaj2hmSf0LMAtYZGbzzOx7ZvZ94BXgJy1dYKKMH5SHO/xjeWmiSxERiZuoRjS7+x+BvsCfgUzgIHCFu89pwdoSqn9eLn065fDCMp1CEpHW41imudjj7n9y99vd/bvuXnQsT2RmE8xspZmtMbMGZ1c1szPMrMrMLjuW/ceamTF+YB6vrNnOrv26x4KItA7HMiFek5lZkPApqAuAgcAVZjawgXY/AubHo67GXDq8J+VVIf6yaGOiSxERiYu4hAIwEljj7mvdvRyYC1xST7sbCd/iM/En8ufMYdCYIYwoWcaDc18i9FDanikTEakRr1DoCdT+uF0SWVfDzHoClwL3HW1HZjbFzIrMrGjr1ha6S1qtO7Nd8/Y8inO78tKMWeH1IiJpLF6hUN8sqnWnIb0HuN3dq462I3ef7e6F7l7YtWvXmBV4mFp3Zrtg5at02fsxDw46L7xeRCSNxSsUSoDetZZ7AZvrtCkE5prZesI38bnXzP4zPuXVUevObFmhSq54bz7/PKmQjbvKElKOiEi8xCsU3gJOMbO+kSm3JwFP127g7n3dvcDdC4C/Ate7+5Nxqu9wde7MduXi5wi489DZlyekHBGReIlLKLh7JXAD4auKlgN/dvelZjbVzKbGo4ZjUufObD32bGf82rd49LRzOFhx1LNbIiIpLZqps2PC3Z8Fnq2zrt5OZXf/YjxqalDkhjtMnx4+ldSnD1efdxrPrQvwt3c387nC3kf/eRGRFBWv00epp86d2UZP+TyndGvPn14rxl236hSR9KRQiJKZcc3ofN7ftIvFG9Py/kIiIgqFY3Hp8F60b5PBg68VJ7oUEZEWoVA4Bu3bZPCZ4T2Z996HbNury1NFJP0oFI7RNaPzKa8K8ehbmg9JRNKPQuEYndwtlzEndebhNzZQWZW29xgSkVZKodAE14zOZ9POAyxYkfh5+0REYkmh0ATnDsijR4dsdTiLSNpRKDRBRjDA5FF9+Peabaw5fRQEAlBQoFlURSTlKRSa6PINb5FZVclDeUPBHYqLw9NtKxhEJIUpFJqo67f/lwtXvMxjp5/Dvszs8Mr9+zW9toikNIVCU23YwNVvP8OeNu2YO2T8YetFRFKVQqGp+vRh+OYVfGL9Yn41ZhK72rSrWS8ikqoUCk01YwaWk8P0f/6eXdnt+dWYSeHptmfMSHRlIiJNplBoqsmTYfZsBuaE+Pz7/+CPIy5m3a/uPzTttohIClIoNEdkeu3/efRHZLZtww8zT0l0RSIizaJQiIFuudlcP/Yk5i/dwutrtye6HBGRJlMoxMh1Z53ICR2y+d4zywiFdBMeEUlNCoUYyc4McvsFp7Jk024ef2dTossREWkShUIMXTz4BIb0Pp6Z81ewv7wyPLq5oEDTYIhIylAoxFAgYHzzwgFs2V3G7F8/GZ72orhY02CISMpQKMRYYUEnLjy9B7/dHOCjQNvDN2oaDBFJcgqFFnDHBadSRYCZZ1995EZNgyEiSUyh0AJ6d8rhS6tf5LHTz+X9vJMO36hpMEQkiSkUWsjXPjuSzvt38b1x11FzgaqmwRCRJKdQaCHHXTOZW0/O5I0+pzO/32jIz4fZszUNhogkNYVCC5r0tc/SL6893/3i3exctkqBICJJT6HQgjKCAWZeNoSte8v4nz+/q5HOIpL0FAotbEjv47nzwoEsWFHK7JfXHr5Rg9tEJMnELRTMbIKZrTSzNWZ2Rz3bJ5vZe5GvV81sSLxqa2nXjM7nwtN7MHP+St5ctyO8cs4cDW4TkaQTl1AwsyAwC7gAGAhcYWYD6zRbB3zS3QcDdwOz41FbPJgZP/zs6fTu2JYbH3mbbXvLwoPY9u8/vKEGt4lIgsXrSGEksMbd17p7OTAXuKR2A3d/1d0/jiy+DvSKU21xkZudyb2TR/Dx/gpumbuYqo0l9TfU4DYRSaB4hUJPYGOt5ZLIuoZcCzxX3wYzm2JmRWZWtHXr1hiW2PIGnnAc3504iH+v2cavJ3yl/kYa3CYiCRSvULB61tV7KY6ZfYpwKNxe33Z3n+3uhe5e2LVr1xiWGB+Xn9GbzwzryT2nX8gr/UYevlGD20QkweIVCiVA71rLvYDNdRuZ2WDgfuASd0/LW5iZGd+79DRO7pbLzZ+7ky39B4OZBreJSFKIVyi8BZxiZn3NLAuYBDxdu4GZ9QEeB65291VxqishcrIyuHfycPZZJjfeeh+VFZWwfr0CQUQSLi6h4O6VwA3AfGA58Gd3X2pmU81saqTZt4DOwL1mttjMiuJRW6KckpfL9z9zGm+u28HPXjhKBmosg4jEUUa8nsjdnwWerbPuvlqPrwOui1c9yeDSYb14c90O7l34Ab065nDlqDqdzNVjGaovXa0eywA6qhCRFqERzQl218WD+FT/rvzvE+/z4GvrD9+osQwiEmcKhQTLzgxy39UjOHdAHt98ain/98q6QxsbGrOgsQwi0kIUCkmgTUaQeycP5/xBeXznb8v43UuROZIaGrOgsQwi0kIUCkkiKyPAr68czoWn92DGs8u5d+Ga8JiFnJzDG2osg4i0IIVCEskMBvjFpKFMHHICP/77Sn7VfWR47EJ+fsNjGXR1kojEUNyuPpLoZAQD/PzyoWQEjJ++sIrKc87glnXrMKtnULiuThKRGNORQhIKBoyZnxvC50b04hcLVvOT51fiXs+sILo6SURiTEcKSSoYMH702cFkBAPMevEDtu8t566LB9E2K3ioka5OEpEY05FCEgsEjBn/eRrXjz2JuW9t5JJZ/2blR3sONdDVSSISYwqFJBcIGN+YcCp/+vJIduyrYOKv/81DrxeHTydFe3WSOqNFJEoKhRRxdr+uPHfzWYw6sTN3PrmE6+e8za5LPx/d1Um67aeIRMnq7cBMEYWFhV5UlNbz5h0hFHLu//dafvz3leQdl80vJg2lsKBTwz9QUBAOgrry88Mzs4pIq2Nmi9y9sL5tOlJIMYGAMeXsk3jsq2MIBozLZ7/OrxaspirUQLirM1pEjoFCIUUN6X08z9z0H1w0uAc/fWEVk2a/xrsbdx7ZUJ3RInIMFAopLDc7k3suH8pPPzeEtVv3ccmsV7jh4bcp3r7vUKNoOqPVES0iEepTSBN7Dlbwu5fW8ruX11EZCjF5VD43jjuZzu3bhN/kp08PnzLq0yccCNWd0XVHRUM4NHRrUJG0dbQ+BYVCmindfZB7Fqzm0bc20jYzyJSzT+S6s/qSk9XAOEV1RIu0OgqFVmhN6V5mzl/B/KVb6JrbhlvOPYXLRvSiTUbw8IaBQPhS1brMIBSKT7EiEle6+qgVOrlbe357dSGPfXU0+Z1ymP7EEs78/gK+/+xy1m2r1ecQTUe0+hxEWg0dKbQC7s6rH2xnzhvFPL90C5UhZ8xJnblyVB/Gv/siWVOP0qegPgeRtKPTR1KjdM9B/lJUwsNvbGDTzgN0aZ/F59ru5oo/fJ8+y985siNafQ4iaUehIEeoCjkvrd7Kw29sYMHyLThwRn4nzhnQjXMG5HFS13bhezioz0Ek7ahPQY4QDBif6t+N311TyCt3jOOWc/qxp6ySHzy3gnN/9i8+9ZOF3D1vGa8WnkNFIHjkDur2RajfQSQtKBSEHh3acvO5p/DczWfxyh3juPuSQeR3bseDrxVz5bhbGH7Tw9ww8Rs8NmgcGzrk4fUNftOke5IojX0gae72WO0jVbh7yn6NGDHCpeXsPVjhf1/yoX/9h4/5iJsf9vzb53n+7fN8xPSn/bo/vuX3vrjGX/9gm+8/8WT3cBwc/pWff2hnDz0UXjYLf3/ooQT9qySlNPZ789BD7jk5h//e5eQcatfc7bHaRzT/ljgCiryB99WEv7E350uhED9VVSFfummXP/jaer/10Xd87MwXa0LipNue9Iuv+Znfdc4Unzv4PC864VTflZUT/uV3j/6PRsKS6M3jqFq6zmh+b/Lzj/6BpLnbY7WPJPsbUChIi9i256C/sPQj/9HFN/rlk77vp97615qgyL99no+68UG/6v7X/Tv/+d/+yODxXtTzVN+Rneuh+v5oWpOjvZnG6lNnc7fHqs7mvBbRvNma1d+m+gNJc7fHah/RBkdLfxiIPMcIcFcoSIuJvEFUWsCLO+T5CyeN9Hv/4wq/9QeP+0W/fPmIsBhw61/8vC/P8i9edpff+cT7ft/CNf63Xz/q7ww727e07+iVBQXJ++m4uRp7M43Fp854nDKJps7q/TQ1WGLxZpssRwqN/Vui/T+JUZAnRSgAE4CVwBrgjnq2G/DLyPb3gOGN7VOhkESO8stYlV/gGzrk+T9OOsPvL5zo3xl3nX/l0ul+wZTf+OBvzz8sMPJvn+d9v/6Uj7jhIZ9w11N+1f2v+62PvuPff2aZ/+6lD/yJt0v8xRVb/J0NH/u6rXt9x94yr6wKNbvGqNs099NcY28gqfJGGE2dzQ2WVAnIeJzminGQHy0U4jJOwcyCwCrgPKAEeAu4wt2X1WrzaeBG4NPAKOAX7j7qaPvVOIUU0cio6N2nnMqmXWVsOq4bmzp0Y1u749mWczxbu57A1jPPZtueMrbuLaO8suFxEbnZGRxfVcbxW0rI3bOTnAwjZ+Cp5PQ7ibZZQXJWrSDnqcfJ2b+HnIqD5JQfJCfgtL3pBtp9+nxysoK0ffZv5Nz23+Ts2kGbqgqsTp1Rj+4+2qy0jY37iGawYGP7aO72WNXZWJvGniMWr3cstsfqOY72b4nz610IFLnbkY3jNHjNzEYD33b38yPL0wDc/Qe12vwWWOjuj0SWVwJj3f3DhvarUEghzXmjBPyhOey+8Ra2Bdqws20uu9q0Z2eHzuz6wrXsHDiEXYuXsOuVN9iZ2ZY9bXLYn5nNgay27OvanQPBLPYfKKOqvvEWR5FZVUFWZQWZhMjs0pmsjzaTVXaQzKrK8Lbq71kZZI39JJlBI7NkI1lvvEZm2UGCoRBBryIYDBIY9ymCAwcQmP1bgjt3EvQQgVCIoIcIhqoIdOhAxv9OI/DO2wTnPkKw7GDN9kBmBsEvfoHgJz5BIGAEb/gawdItBEIhAu4Yjrlj3brBA/+HfeGLWOmW8LrqbRDe/sjD2BVXwEdbDt/mjnXPwx77K2DYZy7FPvzw0Lbqtj16wLx52HPPYXd/Fztw4NA+sttg3/42dvFF4X0MHICFQjX7OOz/ddUqGDsWNm0Kr6LW9p494V//Cj9+6mn42U9h8+bwc992G0y8pN7/L6v3Le7YNLYPi+JJGmzxxBMw88dYyabwv/H2b8Cll4a3nTkaNpWEf772n0KvnvD6G9C7F7VfoprXywxKSiJte9X8HTdPgC4AAAgBSURBVB3xem/aBCNH1rT99L6PEx4KlwET3P26yPLVwCh3v6FWm3nAD93935HlBcDt7t7gu75CIU3E4lNQI9s9EKAskMGBzGz2Z2azP6vW97+/wP7yKg5c80X2RbaXZWRSEcikIphBeUYmFVO/SvnvH6AiGKQ8mElFMJPyYAYVgcj2kWdSURWifNUaykNOeTCTUCBAlQWoCgQJBYNU5bSnqrKSUGUVVYEAbhomJInR5UcXNRgKDUyyH3P1PXndNIqmDWY2BZgC0Ee3lEwPM2bUf2hde4BcY/eabmS79elDdnEx2VUVdDy459D2/HwYkBd+vGtVw8Hymd/Afz/e8PYHbws/Doxt/NRM5KjJN2ygKr+Aqu/eTejySVSGQoRCUOVOVcgJRb7Xfhz+DlXz5hGadS+hLVvw7t3xr16Pn38+TuTp58/HZ8/Gt5Ti3fPwr3wFP/e88AllHP6xAH7/e7y0FM/Lw7/0ZXzcuMjPe/gP758v4n/6E75tG961K3711fjZnzz0HHhkf9Ts99Cywyuv4Pf/Hi8vP/Q6tGkDX/4yfGJMePmVV/G//hW2bYcuneGyy2DMmCNevsY+u0bz0baxD8CN7iOKJ/FGGh21hDfewJ94AnbsgE6dwkcRo0bVbOPBB6G8/NAzZGXB1VfDyEibN8NtvLzi0D6zsuCqq2DUyMh+3sSfeop7jl5kXDqZRwPzay1PA6bVafNbwv0M1csrgR5H2686mtNIYx24ydBRF4sOxdYmVcZcpIJYXCgRQaLHKRA+IlkL9AWygHeBQXXaXAg8R/iI4Uzgzcb2q1BoRWJ1FUlLX32UZIOUROqT8FAI18CnCV+B9AEwPbJuKjA18tiAWZHt7wOFje1TodDKtPTlorGSLHWINOBooaCps0VEWhlNnS0iIlFRKIiISA2FgoiI1FAoiIhIDYWCiIjUUCiIiEgNhYKIiNRQKIiISA2FgoiI1FAoiIhIDYWCiIjUUCiIiEgNhYKIiNRI6VlSzWwP4ZvxJLMuwLZEF9EI1RgbqrH5kr0+SI8a8929a30b4nU7zpaysqHpX5OFmRWpxuZTjbGR7DUme32Q/jXq9JGIiNRQKIiISI1UD4XZiS4gCqoxNlRjbCR7jcleH6R5jSnd0SwiIrGV6kcKIiISQwoFERGpkbKhYGYTzGylma0xszsSXQ+Amf3BzErNbEmtdZ3M7AUzWx353jGB9fU2sxfNbLmZLTWzm5Owxmwze9PM3o3U+J1kq7FWrUEze8fM5iVjjWa23szeN7PFZlaUpDUeb2Z/NbMVkd/L0clUo5n1j7x+1V+7zeyWJKvx1sjfyhIzeyTyN9Tk+lIyFMwsCMwCLgAGAleY2cDEVgXAA8CEOuvuABa4+ynAgshyolQC/+PuA4Azga9FXrdkqrEMGOfuQ4ChwAQzOzPJaqx2M7C81nIy1vgpdx9a65r1ZKvxF8Df3f1UYAjh1zNpanT3lZHXbygwAtgPPJEsNZpZT+AmoNDdTwOCwKRm1efuKfcFjAbm11qeBkxLdF2RWgqAJbWWVwI9Io97EB5wl/A6I/U8BZyXrDUCOcDbwKhkqxHoFfljGwfMS8b/a2A90KXOuqSpETgOWEfkgpdkrLFOXeOBV5KpRqAnsBHoRHgw8rxInU2uLyWPFDj0QlQriaxLRnnu/iFA5Hu3BNcDgJkVAMOAN0iyGiOnZRYDpcAL7p50NQL3AN8AQrXWJVuNDjxvZovMbEpkXTLVeCKwFfi/yGm4+82sXZLVWNsk4JHI46So0d03AT8BNgAfArvc/fnm1JeqoWD1rNO1tVEys/bAY8At7r470fXU5e5VHj5c7wWMNLPTEl1TbWZ2EVDq7osSXUsjPuHuwwmfZv2amZ2d6ILqyACGA79x92HAPhJ/OqteZpYFTAT+kuhaaov0FVwC9AVOANqZ2VXN2WeqhkIJ0LvWci9gc4JqacwWM+sBEPlemshizCyTcCDMcffHI6uTqsZq7r4TWEi4nyaZavwEMNHM1gNzgXFm9hDJVSPuvjnyvZTwefCRJFeNJUBJ5EgQ4K+EQyKZaqx2AfC2u2+JLCdLjecC69x9q7tXAI8DY5pTX6qGwlvAKWbWN5Lgk4CnE1xTQ54GvhB5/AXC5/ETwswM+D2w3N1/VmtTMtXY1cyOjzxuS/iXfgVJVKO7T3P3Xu5eQPh375/ufhVJVKOZtTOz3OrHhM8zLyGJanT3j4CNZtY/suocYBlJVGMtV3Do1BEkT40bgDPNLCfy930O4c76pteX6M6bZnSwfBpYBXwATE90PZGaHiF8Xq+C8Kega4HOhDskV0e+d0pgff9B+DTbe8DiyNenk6zGwcA7kRqXAN+KrE+aGuvUO5ZDHc1JUyPh8/XvRr6WVv+NJFONkXqGAkWR/+8ngY5JWGMOsB3oUGtd0tQIfIfwB6clwINAm+bUp2kuRESkRqqePhIRkRagUBARkRoKBRERqaFQEBGRGgoFERGpoVAQEZEaCgUREamhUBCJITM7w8zei8xp3y4yz31Szd0kcjQavCYSY2b2PSAbaEt4bp8fJLgkkagpFERiLDIf11vAQWCMu1cluCSRqOn0kUjsdQLaA7mEjxhEUoaOFERizMyeJjyldl/Cd7+6IcEliUQtI9EFiKQTM7sGqHT3hyP3En/VzMa5+z8TXZtINHSkICIiNdSnICIiNRQKIiJSQ6EgIiI1FAoiIlJDoSAiIjUUCiIiUkOhICIiNf4fBzzy18lcVYQAAAAASUVORK5CYII=\n", "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": "Python (mdanalysis)", "language": "python", "name": "mdanalysis" }, "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.7.3" }, "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 } }, "nbformat": 4, "nbformat_minor": 2 }