{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Computing mass and charge density on each axis\n", "\n", "Here we compute the mass and charge density of water along the three cartesian axes of a fixed-volume unit cell (i.e. from a simulation in the NVT ensemble). \n", "\n", "**Last executed:** May 18, 2021 with MDAnalysis 1.1.1\n", "\n", "**Last updated:** February 2020\n", "\n", "**Minimum version of MDAnalysis:** 0.17.0\n", "\n", "**Packages required:**\n", " \n", "* MDAnalysis (Michaud-Agrawal *et al.*, 2011, Gowers *et al.*, 2016)\n", "* MDAnalysisTests\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:12.926693Z", "iopub.status.busy": "2021-05-19T06:03:12.926156Z", "iopub.status.idle": "2021-05-19T06:03:14.201615Z", "shell.execute_reply": "2021-05-19T06:03:14.202150Z" } }, "outputs": [], "source": [ "import MDAnalysis as mda\n", "from MDAnalysis.tests.datafiles import waterPSF, waterDCD\n", "from MDAnalysis.analysis import lineardensity as lin\n", "\n", "import pandas as pd\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading files" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The test files we are working with are a cube of water." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:14.208174Z", "iopub.status.busy": "2021-05-19T06:03:14.207100Z", "iopub.status.idle": "2021-05-19T06:03:14.365047Z", "shell.execute_reply": "2021-05-19T06:03:14.365692Z" } }, "outputs": [], "source": [ "u = mda.Universe(waterPSF, waterDCD)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`MDAnalysis.analysis.lineardensity.LinearDensity` ([API docs](https://docs.mdanalysis.org/stable/documentation_pages/analysis/lineardensity.html#MDAnalysis.analysis.lineardensity.LinearDensity)) will partition each of your axes into bins of user-specified `binsize` (in angstrom), and give the average mass density and average charge density of your atom group selection. \n", "\n", "This analysis is only suitable for a trajectory with a fixed box size. While passing a trajectory with a variable box size will not raise an error, `LinearDensity` will not account for changing dimensions. It will only evaluate the density of your atoms in the bins created from the trajectory frame when the class is first initialised.\n", "\n", "Below, we iterate through the trajectory to verify that its box dimensions remain constant." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:14.370749Z", "iopub.status.busy": "2021-05-19T06:03:14.369941Z", "iopub.status.idle": "2021-05-19T06:03:14.376456Z", "shell.execute_reply": "2021-05-19T06:03:14.376999Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n", "[50. 50. 50. 90. 90. 90.]\n" ] } ], "source": [ "for ts in u.trajectory:\n", " print(ts.dimensions)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can choose to compute the density of individual atoms, residues, segments, or fragments (groups of bonded atoms with no bonds to any atom outside the group). By default, the grouping is for `atoms`. " ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:14.381930Z", "iopub.status.busy": "2021-05-19T06:03:14.381082Z", "iopub.status.idle": "2021-05-19T06:03:14.399906Z", "shell.execute_reply": "2021-05-19T06:03:14.400301Z" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/Users/lily/anaconda3/envs/mda-user-guide/lib/python3.7/site-packages/MDAnalysis/analysis/lineardensity.py:111: DeprecationWarning: The structure of the `results` dictionary will change in MDAnalysis version 2.0.\n", " category=DeprecationWarning\n" ] } ], "source": [ "density = lin.LinearDensity(u.atoms,\n", " grouping='atoms').run()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The results of the analysis are in `density.results`." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:14.407350Z", "iopub.status.busy": "2021-05-19T06:03:14.406373Z", "iopub.status.idle": "2021-05-19T06:03:14.409688Z", "shell.execute_reply": "2021-05-19T06:03:14.410098Z" } }, "outputs": [ { "data": { "text/plain": [ "200" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "density.nbins" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:14.415895Z", "iopub.status.busy": "2021-05-19T06:03:14.415329Z", "iopub.status.idle": "2021-05-19T06:03:14.417659Z", "shell.execute_reply": "2021-05-19T06:03:14.418186Z" } }, "outputs": [ { "data": { "text/plain": [ "array([0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0.00053564, 0.00080345, 0.00876966, 0.03507863, 0.00107127,\n", " 0.00348163, 0.00241036, 0.02791588, 0.04277702, 0.01753932,\n", " 0.00160691, 0.00133909, 0.00026782, 0. , 0.00107127,\n", " 0.00107127, 0.00053564, 0. , 0.03400736, 0.01968186,\n", " 0.02339714, 0.01355621, 0.00026782, 0.00107127, 0.00107127,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ,\n", " 0. , 0. , 0. , 0. , 0. ])" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "density.results['x']['pos']" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:14.422349Z", "iopub.status.busy": "2021-05-19T06:03:14.421606Z", "iopub.status.idle": "2021-05-19T06:03:14.424303Z", "shell.execute_reply": "2021-05-19T06:03:14.424709Z" } }, "outputs": [ { "data": { "text/plain": [ "dict_keys(['dim', 'slice volume', 'pos', 'pos_std', 'char', 'char_std'])" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "density.results['x'].keys()" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:14.428405Z", "iopub.status.busy": "2021-05-19T06:03:14.427650Z", "iopub.status.idle": "2021-05-19T06:03:14.430358Z", "shell.execute_reply": "2021-05-19T06:03:14.430951Z" } }, "outputs": [ { "data": { "text/plain": [ "1" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "density.results['y']['dim']" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "execution": { "iopub.execute_input": "2021-05-19T06:03:14.447065Z", "iopub.status.busy": "2021-05-19T06:03:14.444367Z", "iopub.status.idle": "2021-05-19T06:03:14.550228Z", "shell.execute_reply": "2021-05-19T06:03:14.550713Z" } }, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAfhElEQVR4nO3df/Ac9X3f8ed79+4rgQEL8BcCEkI4UWLLdgxYxrgkTcrYroTdKGmbGJwED+1UZQot7iRNcTrTTjp1/+mM4yFlYIhNYhrHDG0SR0PVIUwMtdPUNuJHwLJQrWEwKMhIHluA/ZX0vbt994/dvdvb27vv3n1v7/v97vf1mNHo+93bO+0O4vX96P157+dj7o6IiNRXsNIXICIi1VLQi4jUnIJeRKTmFPQiIjWnoBcRqbnGSl9Akbe85S2+bdu2lb4MEZE148knn/yeu88XvbYqg37btm0cOHBgpS9DRGTNMLPvDHtNpRsRkZpT0IuI1JyCXkSk5hT0IiI1p6AXEak5Bb2ISM0p6EVEak5BL6vCmXaH/37gZbRstsj0KehlVfjq//se/+Z/PMuhY2+s9KWI1I6CXlaF0+0OAIudaIWvRKR+FPSyKrSSgO9EKt2ITJuCXlaFVicOeAW9yPQp6GVV0IhepDoKelkVWu046CN13YhMnYJeVoV2MpJva0QvMnUKelkV0m6bSEEvMnUKelkVWm2N6EWqoqCXVUGTsSLVKRX0ZrbLzA6b2REzu7PgdTOzu5LXnzWzq3Ovh2b2tJk9PK0Ll3ppRZqMFanKkkFvZiFwN7Ab2AHcZGY7cqftBrYnv/YC9+RevwM4tOyrldpS6UakOmVG9NcAR9z9BXdfBB4E9uTO2QM84LGvAZvM7BIAM9sCfBj47BSvW2qmpclYkcqUCfrNwMuZ748mx8qe8xngt4CRi5iY2V4zO2BmB06cOFHisqRO2pFq9CJVKRP0VnAs/39j4Tlm9hHguLs/udQf4u73uftOd985Pz9f4rKkThbbWgJBpCplgv4ocFnm+y3AKyXPuQ74BTN7kbjkc72Z/dHEVyu11e260WSsyNSVCfongO1mdoWZzQE3Avty5+wDbk66b64FXnP3Y+7+SXff4u7bkvd92d1/bZo3IPWg9kqR6jSWOsHd22Z2O/AIEAL3u/tBM7s1ef1eYD9wA3AEWABuqe6SpY60eqVIdZYMegB3308c5tlj92a+duC2JT7jceDxsa9Q1gWN6EWqoydjZVXotleqRi8ydQp6WRXaHT0wJVIVBb2sCosq3YhURkEvq4KejBWpjoJeVgWVbkSqo6CXVUGTsSLVUdDLqpDW6DWiF5k+Bb2sCqrRi1RHQS+rQltPxopURkEvq4JKNyLVUdDLqqDJWJHqKOhlVVDpRqQ6CnpZcVHk3ZKNgl5k+hT0suJaUW+XSQW9yPQp6GXFpWvRg3aYEqmCgl5WXLujEb1IlRT0MjO3/rcnefjZ/HbDvdZKUNCLVEFBLzPz2OHjPPPSyYHj2dKN2itFpk9BLzMTuRfW4LOlm3ZHQS8ybQp6mZl25BQN2FuZoNeIXmT6FPQyE+5xyBfV4Bfbma4b1ehFpk5BLzORBnjRiD07op9krRt35x/d89c8cvC7k1+gSI0p6GUm0tp8UdC3o+WVblod58nv/ICDr7w++QWK1JiCXmYizfJMpndlSzeTTMamPxyy/zIQkR4FvcxEOmov6rpJA3quEUw0ou8GfVtBL1JEQS8z0R3RjyjdbGgEE03Gpm/RiF6kmIJeZqJbox/RdXNWM5wo6NP3LKoHX6SQgl5motd1M/haOhLf2AwnWtTMVaMXGUlBLzORBn3hk7F9pZvJP7utoBcppKCXmUgD3osmY9PSzVxIp6gtZwm9Gr1KNyJFFPQyE9GIHaTS1Ss3Niar0acTvIsa0YsUUtDLTJSp0W9oBoWvL0V99CKjKehlJkZ13aQPSW1shn1PyZal9kqR0RT0MhOj1rpZzHTdTJDz3R8erbZq9CJFFPQyE72um8HXuu2VEz8wpRq9yCgKepmJNMCLum7aHScMjEYYTLR6ZfrZKt2IFFPQy0yko+6iEXurE9EMjTCYbPVK1ehFRisV9Ga2y8wOm9kRM7uz4HUzs7uS1581s6uT4xvN7Btm9jdmdtDMfmfaNyBrQ3uJGn0zCGgEyyvdqI9epNiSQW9mIXA3sBvYAdxkZjtyp+0Gtie/9gL3JMfPANe7+7uBK4FdZnbtlK5d1pB0wrRosrXViWg2AgKzZQa9RvQiRcqM6K8Bjrj7C+6+CDwI7Mmdswd4wGNfAzaZ2SXJ9z9MzmkmvzTsWodGdd20O94t3SxnUTMFvUixMkG/GXg58/3R5Fipc8wsNLNngOPAo+7+9aI/xMz2mtkBMztw4sSJstcva0TaR1+01s1iJ6IZBoQTlm68W6PXGEKkSJmgt4Jj+f+jhp7j7h13vxLYAlxjZu8s+kPc/T533+nuO+fn50tclqwlo5+M9SToi38QlP1sbTwiUqxM0B8FLst8vwV4Zdxz3P0k8Diwa+yrlDWvG/SFT8YmXTdJjb6oBXMU9dGLjFYm6J8AtpvZFWY2B9wI7Mudsw+4Oem+uRZ4zd2Pmdm8mW0CMLOzgA8Az0/x+mWNiEZsDt7KlG7icyb7bNXoRYo1ljrB3dtmdjvwCBAC97v7QTO7NXn9XmA/cANwBFgAbknefgnw+aRzJwAecveHp38bstqlGVy8eqXTSEo36TlhUFQNLJZ+ZOTjv1dkPVgy6AHcfT9xmGeP3Zv52oHbCt73LHDVMq9RaqD3ZOzga612xFxoBElAj/vQVPaHR6sTEQbh5BcqUkN6MlZmYqkdppphQCMJ+nGXQcj+YFCdXmSQgl5mojOiRp+WbgKLg37cFsvsQ1jqvBEZpKCXmYhGdN2kpZt0RF90zsjP9mzpRr30InkKepmJUX307SiiEQTdSdRxSzcd76/Ri0g/Bb3MRLdGP2SHqcYyJmNdQS8ykoJeZiIddRc9DNXxZD36YLIafTbbVboRGaSgl5kY1XWT9r5PPBmrEb3ISAp6mYnek7GDr3UiJzTr1ujH77pRe6XIKAp6mYlRa92kI/pu0I+91k3va7VXigxS0MtMjFqPPvJc0C+r60Y1epE8Bb3MxMium3REP2GNXl03IqMp6GUmel03Ba9FTrCMGn1HNXqRkRT0MhPRiK6bKHIayyjd9NXoFfQiAxT0MhPtETX6tHQTTDoZG6l0IzKKgl5morfWTcFr7gTLeGBKa92IjKagl5kYtXplJy3dTDgZq7VuREZT0MtMdHeYygW9uxM5fZOx469e2ftaffQigxT0MhNRpusm2w6Zjt6zffRjbzwSqXQjMoqCXmainQngbI6nI/xlTcZqhymRkRT0MhPZMI6GjOgn3Xgkv2esiPRT0MtMZMO46OvQeqtXjlu6yf4DQEEvMkhBLzORLcdkgzlttwyCySdjtdaNyGgKepmJTiaAs8HcTpK+ke2jn7BGv6ERsKiuG5EBCnqZic6wGn3ydZCdjJ2w62ZDI1DpRqSAgl5mIluO6f86/j20yR+YSk/f2AxpdSJu+YNv8OfP/O3yLlikRhT0MhP9I/re8WzpZrmrV25oBvxoscNjh0/w9Esnl3nFIvWhoJeZyIZ3X6tlwWTsJOvRm0EzDPjua6cBON3qLPOKRepDQS8z0RlSuuk9MMXEWwl2PF7Pfi4MOHbyFABnNCkr0qWgl5noH9Fnj8eBHAbBstajD81ohgGvvnEG0IheJEtBLzNR1GkDvcXOljUZG6WlG+u+V0Ev0qOgl5kYWrqJeqWbidsrk83Fm2Hvr7NKNyI9CnqZifaQydhe0AcTbzzSieJljucavb/OGtGL9CjoZSaiIe2V05iMjdwJkq6blEb0Ij0KepmJpRY1W97GI/FWhM3Qusc0ohfpUdDLTGT3ii3aeKQRBN3J2LE3HnHvdt2kTrc0ohdJKehlJjpDu27StW56k7Hjr0cPlvTRp1S6EelR0MtM9E3GZjI4rd2no/lGYGPX6N2dMMjV6FW6EekqFfRmtsvMDpvZETO7s+B1M7O7ktefNbOrk+OXmdljZnbIzA6a2R3TvgFZG6IhXTfpD4BGUl8PAhu7dNOJ4idjm434M87Z0OB0W0Evkloy6M0sBO4GdgM7gJvMbEfutN3A9uTXXuCe5Hgb+A13fztwLXBbwXtlHegkDzVBfq2b3mQsxCP78Sdj4/enI/rNm86i1fGx2zRF6qrMiP4a4Ii7v+Dui8CDwJ7cOXuABzz2NWCTmV3i7sfc/SkAd38DOARsnuL1yxoRuXeDuHArwSBTuhmzvB533dCt0W85/ywAzmhULwKUC/rNwMuZ748yGNZLnmNm24CrgK8X/SFmttfMDpjZgRMnTpS4LFlL2pF3g7iojz4d0QeBdde/KSvfddMNenXeiADlgt4KjuX/TTzyHDM7B/gT4BPu/nrRH+Lu97n7TnffOT8/X+KyZC2JIu/W4YuejE1fCyeYjO3W6NPSTRL0qtOLxBolzjkKXJb5fgvwStlzzKxJHPJfcPc/nfxSZS3rZEo3hWvdWCboxxyIu8f/EvjgjotpdSLmz90AqJdeJFVmRP8EsN3MrjCzOeBGYF/unH3AzUn3zbXAa+5+zMwM+BxwyN0/PdUrlzWlkyndFO0fm/bQhzZZ6SYw2HHpefzm3/8pNjZCQDV6kdSSI3p3b5vZ7cAjQAjc7+4HzezW5PV7gf3ADcARYAG4JXn7dcCvA8+Z2TPJsd929/3TvQ1Z7TqRsyFZdCxbmWl30idjJx/Rp6Wb1IZm/OdoRC8SK1O6IQnm/blj92a+duC2gvf9FcX1e1lnOpHTKOq6yU3GhhNNxtIX9N0RvR6aEgH0ZKzMSLa9sqiPPsyO6Mdsf0/bK1MbmnHQn9YyCCKAgl5mJK7RD3bddJ+MTYI+sMlWrwyzpZtGWrrRiF4EFPQyI50o23XTO56fjG0EAe0xSzfxU7eZ0k0znYzViF4EFPQyI51MH33R6pVh3wNT4322e6/0AxrRi+Qp6GUmsn30RevRh2FvCYRoogemet93R/QKehFAQS8zEkVk1rrpHS8a0U+y8Uhf103SXqnSjUhMQS8zEY/oC5ZA8P6um8aEa9309dEn7ZUq3YjEFPRSOXfvn4zNBn2nP+ibodFqj79McbZG3wyNwDSiF0kp6KVyaSVmrijocztMNcOA1kRdN73vzYwNjVAjepGEgl4ql1+hMluj7248kindtMd8YireSrD/AeyNzUBLIIgkFPRSuXQEX1S6aUfefVgqPac1Zn9lJ1ejh7jzRouaicQU9FK5dERfuEyxe3c0n54zbtBHEQNBv6GhEb1ISkEvlUvbJecagztMRVH/8gWNcNL2yv5jG5uq0YukFPRSuSi3nk1nqdLNmN0yUUGNfkMjUNeNSEJBL5Xr5Gr0nlu9Msi1RrbGHNHn16OHeAVLjehFYgp6qVyUK93k16PPjsYbQUB7zBp9upVgVjwZqxG9CCjoZQba3cnY9MnY3mv50XgzDMZur+wU1OjjyViN6EVAQS8z0O2jD4o3B2/kSjeL43bd5NajB43oRbIU9FK5bh99o+DJ2Kh/+YKJum4i+tajh2QyViN6EUBBLzPQ7aMv6LrpRFHfNoDNMKAT+Vi7TMVdN/3HNjYDbSUoklDQS+XyT8Zml5vveK+kkz1nnPVuirpuNjZCjehFEgp6qVx3Mrag6ybKbRqS1uvHmZCNCrpuNmhEL9KloJfKpcFevDl4lFtiOP4rOV7QFzwZ2wjpRD72cgoidaSgl8qlVZiwsOumdxx6LZjjdN4Udd1s0C5TIl0KeqlcOvnaCOINQfrWuslNpDbSEf2YNfrBrhvtMiWSUtBL5dKtAYPACAMbWOsmtOWVbjy3wxT0nsJV6UZEQS8zkGZtaIaZ9dXoo8gHtgGE8Uo3nWiwRt/t3hlzW0KROlLQS+XSydggiMM+/2Rsfq0bmGAydsiIftynbEXqSEEvlYsy+8Lma/SDa93EX49TcokKdphKO3wWNRkroqCX6mX3jA0CG1i9Mt1LFjIll7GCnoGum0k+R6SuFPRSuW7pxozArG89+vyIPg39cda7KarRq3Qj0qOgl8qlQR8WdN3kd4cadySe/tDI1+h7k7EKehEFvVQuDfagoEbf7gwuUwzQKjkZm/3XQpZG9CI9CnqpXJSt0ee6bvITqb2um3IBnX7UQOkmGdFrMlZEQS8zkNbbw6RGH/nw9spe6abciD5aqnQz5m5VInWkoJfKZcM4DIzsYL1TsDk4lK/RRz66dKOuGxEFvcxAJzOiN6O/68b7a/TjrnWT/eyspvroRboU9FK5fNfNQOmm8IGpsqWb+HcbVqPXiF6kXNCb2S4zO2xmR8zszoLXzczuSl5/1syuzrx2v5kdN7NvTvPCZe3IBn1gRif/ZOwy2iujzGdndbtuNKIXWTrozSwE7gZ2AzuAm8xsR+603cD25Nde4J7Ma38I7JrGxcralLZXht1livtH9H2lmzF3mBpWo9eTsSI9ZUb01wBH3P0Fd18EHgT25M7ZAzzgsa8Bm8zsEgB3/wrw/WletKwtUe7J2IH2yoIafdmA7izZdaOgFykT9JuBlzPfH02OjXvOSGa218wOmNmBEydOjPNWWeVG1ejz69HPjdkW6UP66DUZK9JTJuit4Fj+/8Iy54zk7ve5+0533zk/Pz/OW2WVSzM7DOL16PPtlX3LFKdr3ZQd0Q/pujEz5sKARfXRi5QK+qPAZZnvtwCvTHCOrFPpDlPxiL6/vTK/8Uhao2+VXNRsWI0e4glZlW5EygX9E8B2M7vCzOaAG4F9uXP2ATcn3TfXAq+5+7EpX6usUdkdpuKum1zpJhP0ZkYztDG6buLf8zV6iMs3Kt2IlAh6d28DtwOPAIeAh9z9oJndama3JqftB14AjgC/D/yL9P1m9kXg/wI/ZWZHzeyfTvkeZJXrPRlLsgRC/2v51shGEIyx1k06oh98rRlqRC8C0ChzkrvvJw7z7LF7M187cNuQ9960nAuUtS8N20YQxO2V+a0Ec2WXRmjlV6/04j56iEs3GtGL6MlYmYFTrQ5zjWCg68bdiXyw7DI3xkg8rfdbUY0+DPRkrAgKepmB04sdzmqGAEnXTRzO3S0G86Wb0Eo/MJWt/+epdCMSU9BL5RYWO5w9Fwd9aNbtfR9WdmkEAa2Si5qNqtGrdCMSU9BL5RZaHc5Kgz6zleCo3aHKj+iLn4wFku4d9dGLKOilcv2lm94ofGjpJijfXtl7MnbIZKxKNyIKeqleX+km6K11M6wHvhEGY+8wFRb8TW6GKt2IgIJeZuBUq8PGZESf7aNPNxcJB9aSt/IbjyzRdaPJWBEFvczAqcyIPsh23aSj8dxwvDFBe2VR140mY0ViCnqp3EKrzdlz8bN52fXooyGtkXGNfrz2yqIavdorRWIKeqncqcWoW7rJPjDVLd3k/hbGXTdjtlcOqdGr60ZEQS8zcGqx3Ve6SWv00ZDR+Dgj+mhIiybEPzDOqHQjoqCXark7C61MjT7TdZPW6BthUddN2RF9/HvhWjdjrIIpUmcKeqnUmXaEO5mum2wffZQcG1zrpl1yPfrOEqtXajJWREEvFTu12AHoWwKh92RsfM7AEgjjrEevjUdElqSgl0qdasVBn13ULK3ND38ytvwSCKNq9M3kXwZRyX8diNSVgl4qtZCM6Htr3QwugTC41s04I3qSzy0e0QNaBkHWPQW9VKpXukn76HvtlSNXrxxzc/CCAT1zSd+myjey3inopVL50k0QWLc2n4Z0UY2+bOnGR+ww1Uy6edRLL+udgl4qtbDYBnqlm8B64Tws6OfC8uvRd0ZOxsZ/pjpvZL1T0EulTrdGdd0Ur1Mzzp6xaY2+eDI2HdEr6GV9U9BLpbqTsX1dN8laNyNq9J3IuyP/UXpdN4OvpZOxejpW1jsFvVRqId9HH/SWQBhaummkk6glgn5U6UaTsSKAgl4qlpZuNs4VPRlbvA1g2ldfJqCH/bCAuI++7OeI1FljpS9A6q07ou/runH+1Ref5rVTLaDggakkoMt03qTVncL2yrSPXqUbWecU9FKphcUOc2HQDe/AjHbk7H/uWHc9m8G1bpIRfYnOm2G9+NAb0euBKVnvFPRSqdOtDhubvQphmNlhqntsyIi+TMll9Fo36qMXAdXopWILi73dpaC4O2ZwrZv4+zKlm5Hr0YfqoxcBBb1U7FQr6j4sBf0Tr5e8eePAMch23ZQZ0SefW7RMcUN99CKgoJeKnVpsd3vooTfy3tAI+Oh7LwN6bZCpRlC+vbJM141G9LLeqUYvlVpY7O0uBb1A3nrB2fyTn7mCi87dyJbzz+p7T2OMJ1rTGr2N6KPXZKysdxrRS6VOtTp9pZs0jy+/8GzO29jkY+/bOhDS6dIFZXaZGvZ0LYxXAhKpMwW9VOrUYqevdJOua7P1gjcNfc84DzqNrNGrdCMCKOilYqda/aWbtEZ/+YVnD33P+WfPAfDZr77Aj860R37+sM1LQCN6kZSCXiq1sNgp7LrZOiLo33HpeXxy99t49FuvcseDT4/8fB/RR6/16EVimoyVSsWlm8E++ssvGB70ZsY//7kf5wcLLX7/qy9wcmGRTWfPsbDY5rYvPMX3F1pcvXUT//4jO4ZuMA69ydhnj57kH9/z17Qi51ev2cqvJN0+IuuFRvRSmVYnSiZje3/Nfu4n5/lnP3sFl184vEaf2vXOH6MTOV9+/jgAjz1/gscOn6DdifiD//MiX3rmbzNPxg6+38xohsYjB1/lW8de5/jrp/m9x75davljkTpR0EslfvCjRX79c1+nEznvuPTN3eNvnT+Hf/fhHYUj8Lyf3vxmLjp3A39x8FUAHv3Wdzn/7CZfuu06rtq6if/08CF+sLCIWXF7JfQmZD/63sv4l9dv5+Xvn+Lwq29M4Q5F1o5SQW9mu8zssJkdMbM7C143M7sref1ZM7u67Hulftydf/3QMzz10kk+/Svv5oZ3XTLR5wSB8YEdF/OVb5/gh2fafPn541z/totphgGf+sV3cfJUiwefeLmwPp9KJ2R/9X2X84G3XwTAo8kPDpH1YsmgN7MQuBvYDewAbjKzHbnTdgPbk197gXvGeK/UzP/65nd5/PAJ/u2ut/EPr96yrM/60I6LWVjscNsXnuL1020+9I6LAdhx6Xnc8ne2sdiOBrYizHrTXIP3v/VCfuKic7jovI1cedkmHj2koJf1pcxk7DXAEXd/AcDMHgT2AN/KnLMHeMDj4ufXzGyTmV0CbCvx3qn5B7/3V92NLmTlvHLyFO+49Dw+/v7Ll/1Zf3f7PB9731b++OsvsaER8LPb39J97RMf/En+53PH+P6PFoe+/66bruLHkjV1AD6442L+yyOH+eCn//eyr01k2s4/e46Hbn3/1D+3TNBvBl7OfH8UeF+JczaXfC8AZraX+F8DbN26tcRlDfrx+TfpcfdV4B2Xnsft1/9Ed7nh5QgC4z//0rt477bzabW9byXMczY0+K8fu4onXvzB0Pe/5/Lz+77/5Z1bOHL8h5xpa0Agq895G5uVfG6ZoC/6d3G+bWHYOWXeGx90vw+4D2Dnzp0TtUV85sarJnmbrAG/dFVxCeg9l1/Aey6/oPTnXHTuRn73o1dO67JE1oQyQX8UyDYebwFeKXnOXIn3iohIhcr82/oJYLuZXWFmc8CNwL7cOfuAm5Pum2uB19z9WMn3iohIhZYc0bt728xuBx4BQuB+dz9oZrcmr98L7AduAI4AC8Ato95byZ2IiEghW41PCe7cudMPHDiw0pchIrJmmNmT7r6z6DU9GSsiUnMKehGRmlPQi4jUnIJeRKTmVuVkrJmdAL4z4dvfAnxvipezFuie62+93S/onsd1ubvPF72wKoN+OczswLCZ57rSPdffertf0D1Pk0o3IiI1p6AXEam5Ogb9fSt9AStA91x/6+1+Qfc8NbWr0YuISL86juhFRCRDQS8iUnO1Cfr1sAm5md1vZsfN7JuZYxeY2aNm9u3k9/NHfcZaY2aXmdljZnbIzA6a2R3J8dret5ltNLNvmNnfJPf8O8nx2t4zxHtMm9nTZvZw8n2t7xfAzF40s+fM7BkzO5Acm/p91yLo19Em5H8I7ModuxP4S3ffDvxl8n2dtIHfcPe3A9cCtyX/bet832eA69393cCVwK5kn4c63zPAHcChzPd1v9/U33P3KzP981O/71oEPZkNzN19EUg3Ia8Vd/8K8P3c4T3A55OvPw/84kwvqmLufszdn0q+foM4CDZT4/v22A+Tb5vJL6fG92xmW4APA5/NHK7t/S5h6vddl6Aftjn5enBxspsXye8XrfD1VMbMtgFXAV+n5vedlDGeAY4Dj7p73e/5M8BvAVHmWJ3vN+XAX5jZk2a2Nzk29fsus2fsWlB6E3JZm8zsHOBPgE+4++tmRf/J68PdO8CVZrYJ+DMze+dKX1NVzOwjwHF3f9LMfn6lr2fGrnP3V8zsIuBRM3u+ij+kLiP6MhuY19WrZnYJQPL78RW+nqkzsyZxyH/B3f80OVz7+wZw95PA48RzM3W95+uAXzCzF4nLrteb2R9R3/vtcvdXkt+PA39GXIae+n3XJejX8ybk+4CPJ19/HPjzFbyWqbN46P454JC7fzrzUm3v28zmk5E8ZnYW8AHgeWp6z+7+SXff4u7biP/f/bK7/xo1vd+Umb3JzM5NvwY+BHyTCu67Nk/GmtkNxHW+dBPyT63wJU2dmX0R+HnipUxfBf4D8CXgIWAr8BLwy+6en7Bds8zsZ4CvAs/Rq9/+NnGdvpb3bWY/TTwJFxIPxh5y9/9oZhdS03tOJaWb33T3j9T9fs3srcSjeIjL6H/s7p+q4r5rE/QiIlKsLqUbEREZQkEvIlJzCnoRkZpT0IuI1JyCXkSk5hT0IiI1p6AXEam5/w822hiTJhmRNQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.plot(np.linspace(0, 50, 200), density.results['x']['pos'])" ] }, { "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 (mda-user-guide)", "language": "python", "name": "mda-user-guide" }, "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.8.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 }