{ "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 updated:** December 2022 with MDAnalysis 2.4.0-dev0\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": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.498830621Z", "start_time": "2023-06-09T12:02:38.566620541Z" }, "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": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.721579406Z", "start_time": "2023-06-09T12:02:39.539285035Z" }, "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": [ { "name": "stderr", "output_type": "stream", "text": [ "/Users/lily/micromamba/envs/mda-user-guide-dev/lib/python3.10/site-packages/MDAnalysis/coordinates/DCD.py:165: DeprecationWarning: DCDReader currently makes independent timesteps by copying self.ts while other readers update self.ts inplace. This behavior will be changed in 3.0 to be the same as other readers. Read more at https://github.com/MDAnalysis/mdanalysis/issues/3889 to learn if this change in behavior might affect you.\n", " warnings.warn(\"DCDReader currently makes independent timesteps\"\n" ] } ], "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": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.729113515Z", "start_time": "2023-06-09T12:02:39.720089223Z" }, "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": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.792023484Z", "start_time": "2023-06-09T12:02:39.729952761Z" }, "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": [], "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": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.809039922Z", "start_time": "2023-06-09T12:02:39.763907644Z" }, "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": [ "np.int64(200)" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "density.nbins" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.809831995Z", "start_time": "2023-06-09T12:02:39.807626584Z" }, "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" }, "tags": [ "nbval-ignore-output" ] }, "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.00053562, 0.00080344, 0.00876945, 0.03507781, 0.00107125,\n", " 0.00348155, 0.00241031, 0.02791523, 0.04277601, 0.0175389 ,\n", " 0.00160687, 0.00133906, 0.00026781, 0. , 0.00107125,\n", " 0.00107125, 0.00053562, 0. , 0.03400656, 0.0196814 ,\n", " 0.02339659, 0.0135559 , 0.00026781, 0.00107125, 0.00107125,\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": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "density.results['x']['mass_density']" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.810506364Z", "start_time": "2023-06-09T12:02:39.808001286Z" }, "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" }, "scrolled": false, "tags": [ "nbval-ignore-output" ] }, "outputs": [ { "data": { "text/plain": [ "KeysView({'dim': 0, 'slice_volume': np.float64(625.0), 'mass_density': 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.00053562, 0.00080344, 0.00876945, 0.03507781, 0.00107125,\n", " 0.00348155, 0.00241031, 0.02791523, 0.04277601, 0.0175389 ,\n", " 0.00160687, 0.00133906, 0.00026781, 0. , 0.00107125,\n", " 0.00107125, 0.00053562, 0. , 0.03400656, 0.0196814 ,\n", " 0.02339659, 0.0135559 , 0.00026781, 0.00107125, 0.00107125,\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. ]), 'mass_density_stddev': 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.00107125, 0.00122727, 0.01688797, 0.01691979, 0.00177646,\n", " 0.00241031, 0.00279604, 0.02179554, 0.02689655, 0.02096112,\n", " 0.001312 , 0.00133906, 0.00080344, 0. , 0.001312 ,\n", " 0.001312 , 0.00107125, 0. , 0.01700328, 0.03402765,\n", " 0.02131476, 0.01957657, 0.00080344, 0.001312 , 0.001312 ,\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. ]), 'charge_density': 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.00022158, 0.00033237, -0.00033237, -0.00132949, 0.00044316,\n", " 0.00144029, 0.00099712, -0.00033237, -0.00210503, -0.00066475,\n", " 0.00066475, 0.00055396, 0.00011079, 0. , 0.00044316,\n", " 0.00044316, 0.00022158, 0. , -0.00177266, 0.00022158,\n", " -0.00022158, -0.00033237, 0.00011079, 0.00044316, 0.00044316,\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. ]), 'charge_density_stddev': 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.00044316, 0.00050771, 0.00099712, 0.00108553, 0.00073491,\n", " 0.00099712, 0.00115669, 0.00111344, 0.00144029, 0.00112985,\n", " 0.00054276, 0.00055396, 0.00033237, 0. , 0.00054276,\n", " 0.00054276, 0.00044316, 0. , 0.00088633, 0.0018406 ,\n", " 0.00129204, 0.00111344, 0.00033237, 0.00054276, 0.00054276,\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." ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/html": [ "limit_output extension: Maximum message size of 10000 exceeded with 12668 characters" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "density.results['x'].keys()" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.810902504Z", "start_time": "2023-06-09T12:02:39.808816049Z" }, "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": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "density.results['y']['dim']" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "ExecuteTime": { "end_time": "2023-06-09T12:02:39.944318607Z", "start_time": "2023-06-09T12:02:39.809409947Z" }, "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": 9, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiwAAAGdCAYAAAAxCSikAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAAPB1JREFUeJzt3XuUVOWd7/9PXfqCCu0F001rgw1qwOAl6c4kzaTFjNoOZnIbskJyEpMTxLMYklHocZ0lklka4gTHsFg9jFyWinGcmZ/yBzrJ7xwmoZMowdCZCDaRYUhGJwhEukOaSWi80Jeqff7o3rt2Ve269q56dtHv11qslupdVbtrBfjk+32e7xOyLMsSAABAgIVN3wAAAEAuBBYAABB4BBYAABB4BBYAABB4BBYAABB4BBYAABB4BBYAABB4BBYAABB4UdM34Jd4PK4TJ05o6tSpCoVCpm8HAADkwbIsnTlzRo2NjQqHM9dRzpnAcuLECTU1NZm+DQAAUITjx4/r8ssvz/j9cyawTJ06VdLYDzxt2jTDdwMAAPIxODiopqYm59/xTM6ZwGK3gaZNm0ZgAQCgwuRazsGiWwAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFqDCnB2J6fGf/Fqvnzxj+lYAoGwILECFeeGXJ/U3Ow9r/Q/+0/StAEDZEFiACjN4dkSS9NbQqOE7AYDyIbAAFWYkZkmSYnHL8J0AQPkQWIAKMxqLS5JiFoEFwORBYAEqzOh4ZSVOhQXAJEJgASqM0xKiwgJgEiGwABXGbglRYQEwmRBYgAozEqfCAmDyIbAAFSZRYTF8IwBQRgQWoMI4i26psACYRAgsQIUZsbc1s4YFwCRCYAEqzCi7hABMQgQWoMKMxtklBGDyIbAAFYY5LAAmIwILUGHYJQRgMiKwABXGmcNCSwjAJEJgASoMhx8CmIwILECFsXcJsegWwGRSVGDZvHmzmpubVVtbq5aWFu3Zsyfr9bt371ZLS4tqa2s1e/Zsbd26NeO1zz77rEKhkD71qU8Vc2vAOY/R/AAmo4IDy/bt27Vy5UqtWbNGvb29am9v16JFi3Ts2DHP648cOaLbb79d7e3t6u3t1f3336+7775bO3bsSLv26NGjuvfee9Xe3l74TwJMEqMMjgMwCRUcWDZs2KA777xTy5Yt07x589TV1aWmpiZt2bLF8/qtW7dq5syZ6urq0rx587Rs2TItXbpU69evT7ouFovpC1/4gr7xjW9o9uzZxf00wCRASwjAZFRQYBkeHtb+/fvV0dGR9HhHR4f27t3r+Zyenp6062+77Tbt27dPIyMjzmNr167VpZdeqjvvvDOvexkaGtLg4GDSL2AyGImz6BbA5FNQYBkYGFAsFlN9fX3S4/X19erv7/d8Tn9/v+f1o6OjGhgYkCT99Kc/1bZt2/T444/nfS/r1q1TXV2d86upqamQHwWoWE6FhbwCYBIpatFtKBRK+r1lWWmP5brefvzMmTP64he/qMcff1zTp0/P+x5Wr16t06dPO7+OHz9ewE8AVK6RGKP5AUw+0UIunj59uiKRSFo15eTJk2lVFFtDQ4Pn9dFoVJdccokOHTqkN954Qx//+Med78fHS97RaFS/+tWvNGfOnLTXrampUU1NTSG3D5wTRtklBGASKqjCUl1drZaWFnV3dyc93t3drQULFng+p62tLe36Xbt2qbW1VVVVVZo7d64OHjyoAwcOOL8+8YlP6KMf/agOHDhAqwdIYe8SsqxEtRIAznUFVVgkqbOzU3fccYdaW1vV1tamxx57TMeOHdPy5csljbVq3nzzTT399NOSpOXLl+vRRx9VZ2en7rrrLvX09Gjbtm165plnJEm1tbWaP39+0ntceOGFkpT2OIDE4YfS2NbmaCRzOxYAzhUFB5YlS5bo1KlTWrt2rfr6+jR//nzt3LlTs2bNkiT19fUlzWRpbm7Wzp07tWrVKm3atEmNjY3auHGjFi9e7N9PAUwio65TD2OWVfgfYgCoQCHrHKkpDw4Oqq6uTqdPn9a0adNM3w5QMi3f7Napt4clSYfX/qmmVEcM3xEAFC/ff785SwioMPYuIYmFtwAmDwILUGFG48lrWABgMiCwABXGHViYxQJgsiCwABVmlJYQgEmIwAJUkHjcShrJT4UFwGRBYAEqyIhrS7NEhQXA5EFgASrIaCw5oFBgATBZEFiACpIWWEgsACYJAgtQQdJaQgQWAJMEgQWoIKkVFtawAJgsCCxABXFPuZVoCQGYPAgsQAUZjVNhATA5EViACjIaYw0LgMmJwAJUkJG0XUKGbgQAyozAAlSQUQbHAZikCCxABUmtsNASAjBZEFiACpK6hiVOhQXAJEFgASpI2i4hKiwAJgkCC1BBmMMCYLIisAAVhEm3ACYrAgtQQVJ3CVFgATBZEFiACpI+h4XEAmByILAAFSRtDguBBcAkQWABKkjaHBbWsACYJAgsQAVJXXRLSwjAZEFgASoIo/kBTFYEFqCCMJofwGRFYAEqCKP5AUxWBBaggqSP5jd0IwBQZgQWoIIwmh/AZEVgAQJieDR3uYTR/AAmKwILEACbXnhd133jBzr4m9NZrxthcByASYrAAgTA/qO/19mRuP79RPbAkjaHhQoLgEmCwAIEgF0pyRVA0nYJUWEBMEkQWIAAsINKrvwxkrpLiLwCYJIgsAABYFdYLCosAOCJwAIEgNMSyhFA2CUEYLIisAABYLeEcrV40lpCVFgATBIEFiAA7NxBSwgAvBFYgADId5eQffhhODT+PJ9bQr3Hfq8Hv3dIg2dHfH1dAJgoAgsQAPnuEhodHxxXE42MXe9zhWXLi/+lp/a+oR8fPunr6wLARBFYgADIfw7L2Pdrqsb+6PpdYTk7fjzAO8MxX18XACaKwAIEQGJbc/br7MMPa6LjgcXn05rtNTSjcY6BBhAsBBYgAOygkmvXz+j4952WkM8VFvv1RphIByBgCCxAAMSsfFtCqRUWf4OF/Xqpu5EAwDQCCxAA8Xh+i25HUtew+BxY7JcbIbAACBgCCxAAdoUl1xyWWIlbQhYtIQABRWABAiDvOSzx5JaQ/2tYxr6y6BZA0BBYgABILLrNfp2zrblEu4TsAJR6ZhEAmEZgAQKg0NOaSzU4LrGGhcACIFgILEAA5LtLyD78sFSD4+wAREsIQNAQWIAAyHeXkF1hqY6Ek57n2304i24JLACChcACBED+c1hKXGGhJQQgoAgsQADE8x3Nn3L4od9zWJzR/FRYAAQMgQUIgLizSyjPCkvJtjWPt4R8DkIAMFEEFiAA8pnDYllW2llCpZp0S4UFQNAQWIAASKxhyXzNqOubidH8/t6Hs0uINSwAAobAAgRAPI85LO4QUeqW0DAVFgABQ2ABAiCexy6hEddslNK3hKiwAAgWAgtgmGVZTlDI2hJyhYjqEldYGBwHIGgILIBh7pCSbRCcvRA2FJKqIqHx5/q9rXnsK3NYAAQNgQUwzN3Wyd4SGvteVTiscCiU9lw/UGEBEFQEFsAwd0jJ3hIaCxHRSEiR8HiFxedcEWOXEICAIrAAhuVdYRkPEdFwKFFhKdFofnYJAQgaAgtgmDukZMsfdpumKhJ2KiylG81PhQVAsBBYAMPcbZ1sFRY7RIy1hHJfX9S9cJYQgIAisACGuds62SomI/YalpIuuh1/L84SAhAwBBbAsOQ1LJmvs0fzV0VCpQsscSosAIKpqMCyefNmNTc3q7a2Vi0tLdqzZ0/W63fv3q2WlhbV1tZq9uzZ2rp1a9L3n3vuObW2turCCy/U+eefrxtuuEH/+I//WMytARUneQ1LHhUW1xqW0rWEqLAACJaCA8v27du1cuVKrVmzRr29vWpvb9eiRYt07Ngxz+uPHDmi22+/Xe3t7ert7dX999+vu+++Wzt27HCuufjii7VmzRr19PTo1Vdf1Ve+8hV95Stf0Q9+8IPifzKgQiRva85jDUu4hBUWdgkBCKiCA8uGDRt05513atmyZZo3b566urrU1NSkLVu2eF6/detWzZw5U11dXZo3b56WLVumpUuXav369c41N910kz796U9r3rx5mjNnju655x5dd911eumll4r/yYAKkX9LKH2XkN9LTRKD46iwAAiWggLL8PCw9u/fr46OjqTHOzo6tHfvXs/n9PT0pF1/2223ad++fRoZGUm73rIs/ehHP9KvfvUr3XjjjRnvZWhoSIODg0m/gEqU7y6hEY9dQv5va068brb2FACUW0GBZWBgQLFYTPX19UmP19fXq7+/3/M5/f39ntePjo5qYGDAeez06dO64IILVF1drY997GP6+7//e916660Z72XdunWqq6tzfjU1NRXyowCBESuwJVSO0fwS5wkBCJaiFt2Gxv+ytFmWlfZYrutTH586daoOHDigl19+WX/zN3+jzs5Ovfjiixlfc/Xq1Tp9+rTz6/jx40X8JIB5SS2hLEtH7JZQ0mh+n6sg7vDEeUIAgiRayMXTp09XJBJJq6acPHkyrYpia2ho8Lw+Go3qkksucR4Lh8O68sorJUk33HCDDh8+rHXr1ummm27yfN2amhrV1NQUcvtAIFl5VlgSLaFEhcXPwGJZVtKkXSosAIKkoApLdXW1Wlpa1N3dnfR4d3e3FixY4Pmctra2tOt37dql1tZWVVVVZXwvy7I0NDRUyO0BFSmW72j+8Z07VeGQazS/f/eR+t7MYgEQJAVVWCSps7NTd9xxh1pbW9XW1qbHHntMx44d0/LlyyWNtWrefPNNPf3005Kk5cuX69FHH1VnZ6fuuusu9fT0aNu2bXrmmWec11y3bp1aW1s1Z84cDQ8Pa+fOnXr66acz7jwCziV5H34Ydy+69b/CkvpaVFgABEnBgWXJkiU6deqU1q5dq76+Ps2fP187d+7UrFmzJEl9fX1JM1mam5u1c+dOrVq1Sps2bVJjY6M2btyoxYsXO9e8/fbbWrFihX7zm99oypQpmjt3rv7pn/5JS5Ys8eFHBIIt/7OEEoPjSrHoNvWlRqiwAAiQggOLJK1YsUIrVqzw/N5TTz2V9tjChQv1yiuvZHy9hx56SA899FAxtwJUvKSzhLK2hOxdQq4Ki6+BJfm1mMUCIEg4SwgwLO/R/HHXaH67wlLClhBrWAAECYEFMCye5xoWp8ISCSlcgsFx6S0hKiwAgoPAAhiW9xwWew1LuDTbmtNbQlRYAAQHgQUwLN9Jt167hPyssFgp+YRFtwCChMACGOYuZOQ1hyVpcFz2dS8F3QfbmgEEGIEFMCyetEsoj0m3rl1CY8/3/z6kxJoZAAgCAgtgWN6HH3rsEpL8awulhqUR1rAACBACC2CYe5dQ9pZQosISdv3J9WvhbfpofiosAIKDwAIYlvdo/lj6ottczykEc1gABBmBBTAsnmdLKBa3Dz9MLLode9yvwJL8+2ECC4AAIbAAhrmDQrZlI3aHJpy66NanXJE65p+WEIAgIbAAhuXbErIDRSSk5EW3pVrDwqJbAAFCYAEMy78lNB5YwiGFw2XYJUSFBUCAEFgAw5IrLJmvs8NMaLy64pzYzKJbAJMAgQUwLJa0rTlLS8hKVFikRFvIrwpL6nuP+jj2HwAmisACGJbcEsp8ndMSGg8qfp/YzC4hAEFGYAEMS9ollG0Ni2uXkJQILqVrCVFhARAcBBbAMHeFJFu1xNklNP6nNuxzSyh1UxBrWAAECYEFMMxd2chWLLGDSdhpCZW2wjLCGhYAAUJgAQzLdw5LLHXRbdiusPhzH+wSAhBkBBbAsMIHx41XWPxuCaW8DHNYAAQJgQUwzEpadJv5OrvC4iy6DdvPKVFLiAoLgAAhsACGxZLWsGSbwzL21a6s+L1LKG0OCxUWAAFCYAEMK3qXULjELSHOEgIQIAQWwLB4nqP5U3cJ+T2aPzX4UGEBECQEFsCw1EMHM7WFMo/m9+c+0nYJUWEBECAEFsCw1KpKpipL+mh+v88SSv49u4QABAmBBTAsnhI4MrV40nYJlXg0P7uEAAQJgQUwLLUllCmAJBbdlqbCkvoyrGEBECQEFsCwtApLhsKGU2EJJc9hSQ08ft0HFRYAQUJgAQxLrZBkrrCMfU1ddJsaNIqVvuiWCguA4CCwAIalL7rNsIYlZdFtqMSj+TlLCECQEFgAw1IDSqb8YV83nlN8n8OSvuiWCguA4CCwAIalVkhMzWFJfV/WsAAIEgILYFj6LqEM16XtEvJ+frHs97UrOKxhARAkBBbAsNRFs5nWpGQazZ/twMRC2K9fPb79iAoLgCAhsACGpa4dydwSGvvqVFh8X3Q79jo10bG/FpjDAiBICCyAYamFjHxH80dKNJq/OhqRxFlCAIKFwAIYlr5LKNdo/rHfl2o0v11hYZcQgCAhsACG5T84LtNofn/uw76NRGChwgIgOAgsgGFpu4RyjOZ3WkL2GhafKyzVrGEBEEAEFsCw1EW2XhUWy7KcNSahlDUsvo3mj6e0hFjDAiBACCyAYfm0hNyXlPq0ZrvCYln+vTYATBSBBTAsn11C7uCQaAnZ15emJSSNrWM5dOK0fjt41pf3AIBiEVgAw/KZw+K+xt4l5HeFxX5fe3CcJP3m9+/qE4/+VF9+8ue+vAcAFIvAAhiW3hLKfk3a4DifR/PXjM9hkaTXT76lWNzSiT+868t7AECxCCyAYakVFq+KiTuUhFN2Cfm26Hb8PaJ2r0nS786MtYKG2eIMwDACC2BYPoPj4l4VFp/nsLgPV6waDy2/OzMkSRoajft2ZhEAFIPAAhiWWlHxygWei259Pq3ZfplIKKTo+EKZk+OBxbI4vRmAWQQWwLDUcSdeFZakllA4uSXkV+XDft9QKOS0hezAIo1VWQDAFAILYFjapFvPltDYV7sdJJVuDks4JFVF7ApLYjvzMIEFgEEEFsCwfHYJ2SHGlVdKNpo/HEqsYfntoLvCEvPlfQCgGAQWwLB8RvPbocbeIST5P5rfvo9wWM4allNvJQILFRYAJhFYAMPSDz/MPDjOuyXk032Mv07IVWFx3wprWACYRGABDCtkNH/EXWGx57D43BKKhEKKRtL/aqDCAsAkAgtgWGpFJdto/nAJF91arnUyUfdimXGsYQFgEoEFMCx9l5DHNR67hPxfdDv2dawllP5XAy0hACYRWADD8pl0673odvx6n0fzh11zWNwILABMIrAAhqUGDq+KSWLRbeKxcsxhcRsaIbAAMIfAAhiWGlC81rBkW3Tr+xwW11lCbhyACMAkAgtgWNpofo9c4Lno1u/Tml1tJ3sOi9vQCItuAZhDYAEMs6sn1eNtGM/R/Fb6GhanJeTTmYTJLSEqLACChcACGGaHkajHsDab9y6h5Of7dR+ZKywEFgDmEFgAw5zAEs58+nJil1DisZKN5g/Jc5cQFRYAJhFYAMPsMGLvzMm+S6j0pzWHQiGnPSVJ1dGx/6bCAsAkAgtgmLMDKJytJeQxh8Xn0fyxDHNYGqbVSpKGYyy6BWAOgQUwzA4odoXFsyVUhgqL5Zr14j5LyA4sVFgAmFRUYNm8ebOam5tVW1urlpYW7dmzJ+v1u3fvVktLi2prazV79mxt3bo16fuPP/642tvbddFFF+miiy7SLbfcop///OfF3BpQcewKSeKEZI+WUDw9sCTmsPh0H+7Tml3v01BnV1gILADMKTiwbN++XStXrtSaNWvU29ur9vZ2LVq0SMeOHfO8/siRI7r99tvV3t6u3t5e3X///br77ru1Y8cO55oXX3xRn//85/XCCy+op6dHM2fOVEdHh958883ifzKgQtgVEruq4TWHxXs0v89zWJJaQq4KSx0VFgDmFRxYNmzYoDvvvFPLli3TvHnz1NXVpaamJm3ZssXz+q1bt2rmzJnq6urSvHnztGzZMi1dulTr1693rvnnf/5nrVixQjfccIPmzp2rxx9/XPF4XD/60Y+K/8mACpG6S8h7DsvY13IsunXvEoqEQ5p+QbUkKiwAzCoosAwPD2v//v3q6OhIeryjo0N79+71fE5PT0/a9bfddpv27dunkZERz+e88847GhkZ0cUXX5zxXoaGhjQ4OJj0C6hEiQpLtsDisa25VKP5XbuELjqvWjXRiCRpaJRFtwDMKSiwDAwMKBaLqb6+Punx+vp69ff3ez6nv7/f8/rR0VENDAx4Pue+++7TZZddpltuuSXjvaxbt051dXXOr6ampkJ+FCAQLMtyKhv2sLa8dwmV6LTmUChxLxedV6Wa8W3Nw5zWDMCgohbdhkLJQ6Usy0p7LNf1Xo9L0iOPPKJnnnlGzz33nGprazO+5urVq3X69Gnn1/Hjxwv5EYBAcGeNfEbzRzzOEvKvwiLnPexqz0XnVyfmsBBYABgULeTi6dOnKxKJpFVTTp48mVZFsTU0NHheH41GdckllyQ9vn79en3rW9/SD3/4Q1133XVZ76WmpkY1NTWF3D4QOO5wkn00v8cuoRIuuk20hKpcLSECCwBzCqqwVFdXq6WlRd3d3UmPd3d3a8GCBZ7PaWtrS7t+165dam1tVVVVlfPYt7/9bX3zm9/U97//fbW2thZyW0DFci+YjWabw+LREvK7wuIezX/j1ZdqbsNUffr9l1FhARAIBbeEOjs79cQTT+jJJ5/U4cOHtWrVKh07dkzLly+XNNaq+dKXvuRcv3z5ch09elSdnZ06fPiwnnzySW3btk333nuvc80jjzyir3/963ryySd1xRVXqL+/X/39/Xrrrbd8+BGB4EqqsGSpmGQfze/TvbjmsLy3Yaq+v/JG/en8GaxhARAIBbWEJGnJkiU6deqU1q5dq76+Ps2fP187d+7UrFmzJEl9fX1JM1mam5u1c+dOrVq1Sps2bVJjY6M2btyoxYsXO9ds3rxZw8PD+sxnPpP0Xg888IAefPDBIn80IPiSKizhzIPg7FDiNZrfqyJT1L1Y6VUcyXWWELuEABhUcGCRpBUrVmjFihWe33vqqafSHlu4cKFeeeWVjK/3xhtvFHMbQMVzD4nLbzR/4rHxjTwlGc3vRoUFQBBwlhBgkLslVPxofv9Pa3ZjDQuAICCwAAa5w0YkyxyWuEe7ppS7hNzsXUJUWACYRGABDIrHEztznEFw+e4SCpemwhJOGY9UwxoWAAFAYAEMci90DTuLaNOv89olZLeEvA5LLIaVscKSWMPi1wJfACgUgQUwyKmchEPO2hGvRbSeu4R8PvzQfp3UAdT2Gpa4JY369F4AUCgCC2CQXbCIhEJOKyb7aP7EY/6P5k+v4kiJNSwS61gAmENgAQxyj9y3A4i50fxjXzPNYZHYKQTAHAILYFDMci+6zTwILttpzaUYze8WCYecoXZUWACYQmABDIonrWEZf6zQ05p9rrB4naLOtFsAphFYAIOcCbah/FpCJuawSEy7BWAegQUwyN6SHA67Ft16Hn449tWzwuLXotu4d0tIYtotAPMILIBBcc8KS+aWkDtM2IPj/Npp7Cy69Ugs9k4hAgsAUwgsgEHu3T+hfFpCrjAR9XkOS7aWEGtYAJhGYAEMsts5oTxH80cyDI7zYwJtptH8EmtYAJhHYAEMinvMYcl3NH/U9d9+VFkyjeaXWMMCwDwCC2CQs5g25G4J5TeHJeoae+vHyPy4q9qTigoLANMILIBB7rUpdsHEq1qSq8LiR2DxCkW2ahbdAjCMwAIY5L1LKP06r9H8SYElNvEg4Zxr5LlLiAoLALMILIBBXhUW79H8Y1+9BsdJ0kjMv5ZQ9jks7BICYAaBBTDIfZZQYq5Kfqc1h0IhVUX829qcbTQ/FRYAphFYAIMsK32XkFf2yDQjxa6yjPjQEspnND9rWACYQmABDHK3epzR/HnuEpKkqvDYH2E/Ft3ab+s9h2Vs0S0VFgCmEFgAg2Iec1i8zxJKX3QrSRGnJTTxIJF9lxBrWACYRWABDHLvEip0NL8kRccrLP4uumUNC4DgIbAABtlBJBRSjpbQ2NdIaktovMIy6ktgGfsa9vhboTrCGhYAZhFYAIPcrR673ZN9NH/y4/ZzRn1oCWUbzV9TRYUFgFkEFsAgd2ApdDS/JFVF/Ft0m3UOCxUWAIYRWACD8t0llGnRrT3t1s+WkNccFkbzAzCNwAIY5HVas9dIFa/R/O7f+9EScu4l6xwWdgkBMIPAAhiUNOk2y2j+TDt4nJZQiXcJVbNLCIBhBBbAIPfalGxrWOIeZwlJUtTeJeTraP707zHpFoBpBBbAIPdo/ki2OSwZdgkl1rD4OJrfY9UtFRYAphFYAIOSTmse/9NYyC4hZ3BcmUbzs4YFgCkEFsAge+lJJJRYdJt9Dot3S8iX0fz5rGHxoZIDAMUgsAAGuXcJhZxdQlkqLBm2NZdrNP/QCIEFgBkEFsAgu6qRezS/95bjqE+7hCzLytESosICwCwCC2CQ+/DDolpCYX9aQu739K6wjK9hocICwBACC2CQ1+C4ghbdRvw5rdn9nqxhARBEBBbAIGc0fzjXaP6xr6kVlqpw5nUvhXA/PeTxt4LdEorFLV+2UANAoQgsgEExj5aQV7Ek12nNIxNsCeVbYZGosgAwg8ACGBR3Wj1y5rB4jea3KyipBxP6tejWHViynSUksY4FgBkEFsAg93TZ7KP5vXcJVfk0mj+pJeSxSygaCTstKyosAEwgsAAGebWEvLo7sQy7hCI+jebP1RKS2CkEwCwCC2CQe5dQJOsuobGvGU9rnmCFxXJlEK85LFJiHQvj+QGYQGABDPLaJVTIHJZEhaW025oltjYDMIvAAhiUGIefWFAbyzbpNuVPrL2teXSCu4Tc75khryTey4djAACgUAQWwKDkSbfJjyVdV6bBce7glKoqar8XFRYA5UdgAQxyH2poH2zo1RLKtejWr9H8mdpBUmK9jB8HLQJAoQgsgEF5V1gynKTsbGv2rcKSObAkToamwgKg/AgsgEFJFZasc1jkXOcWDfuzS8h+epa84iy6neh6GQAoBoEFMMguViQdfphtDkvaGhZ/Ft1mWiOT9F7jYWl4lJYQgPIjsAAGuRe75nVac8qfWLvC4uei20wSM1+osAAoPwILYJB77UgowxqWuKvdk6nC4tdpzaktJ7fEolsCC4DyI7AABsVck24TFZaUa9wHE6atYfFnIWw+i27tBb4jtIQAGEBgAQxyT7DNdFqzu3qStujWp9OarTxaQs7MF1pCAAwgsAAGxVyLXSMZKizuFlFaSyjsc0soS4Wl2g4sowQWAOVHYAEMch9qGMqw6NYdRjK2hCa6S2j8PTNNuZVcM18mGI4AoBgEFsAgy2kJJdoxqdWSeNJJyt6nNU+0wpKo9GS+xm4JcfghABMILIBBMddiVzuMpO5qjmdZdBtxFt1OdA2L9+u7Vfm0XgYAikFgAQzy3iWU0hJy/T41TziD48q5S4gKCwADCCyAQV67hDLNYQl5nKTsV0son9H8HH4IwCQCC2BQLJ5Y7JprDkvqDiHJ1RLyadFt1tH8VFgAGERgAQyyw0kkaQ2L9y4hrym0VeHyzWGpZtItAIMILIBB8Xj+u4SyVVgmutU4luE0aDdaQgBMIrAABsWSzhLK0RLyqrCUcdEtLSEAJhFYAIPiSbuEEo+720LZZqT4NZo/n9Oaq533IrAAKD8CC2CQU2EJh5IqKO4qSzxLhSXqU0vIzkdZKyw+zXwBgGIQWACD4h6j+aXkrc1ZA4szLr8Mo/mjLLoFYE5RgWXz5s1qbm5WbW2tWlpatGfPnqzX7969Wy0tLaqtrdXs2bO1devWpO8fOnRIixcv1hVXXKFQKKSurq5ibguoOHHXlmV3Hol7tIS8wkTU3iXk2+GHma+xdyQRWACYUHBg2b59u1auXKk1a9aot7dX7e3tWrRokY4dO+Z5/ZEjR3T77bervb1dvb29uv/++3X33Xdrx44dzjXvvPOOZs+erYcfflgNDQ3F/zRAhUlsWU5ux7gLJtl2CdltGsua2PA491qaTKqiHH4IwJyCA8uGDRt05513atmyZZo3b566urrU1NSkLVu2eF6/detWzZw5U11dXZo3b56WLVumpUuXav369c41H/zgB/Xtb39bn/vc51RTU1P8TwNUmFhShcW7JZRtl5DdEpImVvnIpyVkV3OGR6mwACi/ggLL8PCw9u/fr46OjqTHOzo6tHfvXs/n9PT0pF1/2223ad++fRoZGSnwdhOGhoY0ODiY9AuoNO7KRihHSyjs8afVno3ivq6o+8inJRTxp/0EAMUoKLAMDAwoFoupvr4+6fH6+nr19/d7Pqe/v9/z+tHRUQ0MDBR4uwnr1q1TXV2d86upqano1wJMSZzhk8cuoSyD46SJbW3OZw5LdZQ5LADMKWrRbWrZ2LKsrKVkr+u9Hi/E6tWrdfr0aefX8ePHi34twBSv05qlDHNYsmxrliZ2nlA+o/lpCQEwKVrIxdOnT1ckEkmrppw8eTKtimJraGjwvD4ajeqSSy4p8HYTampqWO+Cipd5l5Drv+OZKyx2ZSYWt3xpCWXd1kxLCIBBBVVYqqur1dLSou7u7qTHu7u7tWDBAs/ntLW1pV2/a9cutba2qqqqqsDbBc4t7vUp7rDgDh/ZFt1K7oFuxVc+YllCka2K0fwADCq4JdTZ2aknnnhCTz75pA4fPqxVq1bp2LFjWr58uaSxVs2XvvQl5/rly5fr6NGj6uzs1OHDh/Xkk09q27Ztuvfee51rhoeHdeDAAR04cEDDw8N68803deDAAb3++us+/IhAcKUOhbMzifdofu8wUeXDeH5nDUuWvxH8eB8AKFZBLSFJWrJkiU6dOqW1a9eqr69P8+fP186dOzVr1ixJUl9fX9JMlubmZu3cuVOrVq3Spk2b1NjYqI0bN2rx4sXONSdOnND73/9+5/fr16/X+vXrtXDhQr344osT+PGAYLMLKXZlIxwKKW5ZSS0hO7tkqrD4cWJzXqP5xyssw1RYABhQcGCRpBUrVmjFihWe33vqqafSHlu4cKFeeeWVjK93xRVXJP0/SmCySJ1iGw6HpLjlva05Q5ao8mE8fz5zWDj8EIBJnCUEGJQ6YdYOJV6D47x2CUmu8fwTagkp6f093ydij+bn/1wAKD8CC2BQLGXGit2ScRccs+0SkvxpCcVzrJORWHQLwCwCC2BQ6hRbOzB47RLKVGFxWkI+jObPOjguwuGHAMwhsACGWJalofEhbDXRiCQ54/m91rBkqrBEfZiPUkhLKD7BgxYBoBgEFsCQIdfE2CnVY4HFrnB4jubPMYel1KP5q3w6aBEAikFgAQw5OxJz/rs2OvZH0Q4lyXNYxr5mXHRrry3xYzR/HnNYJAILgPIjsACGnB0Z+0c/Gg457ZbELqHEdYlFt96vExlPGTEfdgnlM5pfYngcgPIjsACGvDteYamtijiPhZyWUCIQ5GoJVYUnPocl1zRd+/3tb1NhAVBuBBbAkLMegcXOJF67hDJVP5yWkA9rWDJVcWx2lWWERbcAyozAAhiSCCyJP4bFzGGxB8dNZOdOPqP5pUQ1Z2SUCguA8iKwAIbYa1iSKyzpLaFYPMcuIR8GuuUzml+SqqL2FmoCC4DyIrAAhnhWWMb/M3k0v/290lVY8pnDIiVaQsOjtIQAlBeBBTDEDixTPCssiety7RKy57BMZF1JPnNYJH8W+AJAMQgsgCFnR70W3XrMYcl1+KEfo/njueewSImWELuEAJQbgQUwxF7DYo/llxKj+ZN2CeVYdGu3afxpCWWvsNjVHFpCAMqNwAIY8u5w5l1Cni2hDBUW+3E/tjXnbAlFWHQLwAwCC2CI3RKa4jGHxUoaHDf+vRKe1uyM5s93DgstIQBlRmABDMm+rTlxXczK3hKKOAthSzuaX0qEo4lUcwCgGAQWwJBsg+OSRvPHs1c/7G3NE2nT5NsSilJhAWAIgQUwxHNbs+ccluy7hBItoeKrHk4VJ8ffCNX2GhYqLADKjMACGGIHlpock25zjeaPOBWWMozmHw9Hw1RYAJQZgQUwxGsNi3NasysP5BrN78eiWzsU5VrDEqXCAsAQAgtgyLuea1jGvhbSEor6UGHJdzR/NWtYABhCYAEM8VrDEsk6mj/XpNvSz2Hx46BFACgGgQUwZCjLtuaCRvM7ZwmVcw4LLSEA5UVgAQxJnCWU+GMYclpCievsYkbmCsvER/PnCkW2KiosAAwhsACGOKP5o+kVlpirwmLl2HJsV1gm1hJKfv9MnNH8BBYAZUZgAQxxKizV6XNYklpC8RwtIR+qHvm2hOwFvsO0hACUGYEFMMTZ1uxRYfHcJZSpJTSeMiZ0WvN41sk5mj868S3UAFAMAgtgyNlspzW78kDOXULjVY+RCW1rzm+XENuaAZhCYAEMcU5rrk4/rTm5wjL+vRwtoQkNjstzDosf4QgAikFgAQyIxS1na3CullCiwuL9Wv4Mjss+Tddmt4RGRqmwACgvAgtggD00Tsowmj9pW3P2MOFPhSW/0fxVPoQjACgGgQUw4F1XYKmJTmw0vz0bZUKLbvNsCXH4IQBTCCyAAc5JzdFwUhAJFzGa3z6teSLTZ/MfzT/+XrSEAJQZgQUwwOukZinR9nHPYYnnqrDYg+PKMJrf3iVESwhAuRFYAAPOepzULLlG88fTdwllrrDYgaX0c1g4/BCAKQQWwIBEYEmusCRG8yceizuTbr1fK+qMy5/4WUI5dwkxhwWAIQQWwAC7JTQlLbCMffUczZ+h+lHlwy4hTmsGEHQEFsAAZ9FthgqL1y6hTNUPX1pC40/Nua3Zh3AEAMUgsAAG2Nuaa6Opa1gK3yVU5cNC2Hx3CdnvxeGHAMqNwAIYkGkNy3geKGgOi3344UTWleQ9mp8KCwBDCCyAAWdHM61hsbc1Jx7L9/DDiZ3WzOGHAIItavoGgMloKOO25pDz/bue3qcr33NBzjUsiaqHDy2hHCWWKItuARhCYAEMeHc407bmsa+/+M1p7f7P3+lHh3+rGXVTxr6XK7BMYHBcPO9dQsxhAWAGLSHAgLOj2eewHP/9O5LG1pa8+Yd3JeVuCcWt5IFzhUisYWEOC4BgIrAABmQazW9XON78/btpz8k8OC4RMkaKrLIUOodlIu0nACgGgQUwINNofrvtM+RxuGCm6kfUlTKKXXhb6BwWTmsGUG4EFsCAd3OM5veScdGtq/RS7GLYXNN0bX7MfAGAYhBYAAOGcozmt9nbiMe+l7vCUux8FMvZiZT9OjuwxOJW0etlAKAYBBbAgIwtoZRQctv8Bue/M1VYwuGQE3RK3RLyY70MABSDwAIYkGmXUGpgWNLa5Px3pl1Ckms+StGBpbDBcRKzWACUF4EFMMCew1ITzdwSuvC8Kn149sWqm1IlSTqvJvlat6rxJ8aKDBF5j+b3of0EAMVgcBxggL2teUp16llCiUBw2YVTFI2Etel/fECvnzyj2dPPz/h69vMmvq05e2KJhEMKhcaODmCnEIByIrAABjgtoQynNUtS44VjE24/ctV0feSq6Vlfb6LzUfLdJRQKhVQVDms4FqclBKCsaAkBBpzNMZpfGquw5Gui4/nzHc0vJWax0BICUE4EFsAA+7TmbHNYCgos4YlVWOzToXMdfihJVVHG8wMoPwILYIC9rTnbHJbLLgpmhcUOR7SEAJQTgQUoM8uyco7mlwqtsIw9r/fYH4q6p3znsEhSNSc2AzCAwAKU2XAs7gSEmiwtocYCAsvt186QJD30fw/riT2/Lvie8p3DIrlmvlBhAVBGBBagzOwtzVLmllBNNKzpF1Tn/Zqrbrlayz7SLGkstHzvFycKuid7zH624XS2KiosAAwgsABlNjTeDgqHEv/42+wKx2UXTsmrPeM8LxzSmo/N0/9ccIUk6bu9bxZ0T4mWUO5rJ7qFGgCKwRwWoMzsCkttVSQtlNi/L6Qd5H7uZ1ou11N739DPfn1KI7G4Ey4kaeOPXtOmF15X3LJUHQnrW39+rT55w2WSCmsJ2a/57ydOa+X2Azr97rAk6ZZ59dryxZaC7xsA8kGFBSizd0e8Z7BI0rWX1ak6EtbCqy8t6rWvmTFNF59frbeHY0kLcONxS9/56RENjY4NfHt7OKZv/p//0NtDo2Pfd7Y1534Puyr05EtHNPDWkEZilkZilv713/t1/L/fKeq+ASAXAgtQZqfeHpKUvn5Fkv6o+WId/EaH7rpxdlGvHQ6HtGDOJZKkl177nfP4oROD+v07I7qgJqo9//ujuuKS8zTw1rCe2HNEUv6j+aXEotuTZ8Z+jq1f/IDeP/NCSdKe1waKum8AyIXAApTRL/sHdc+zByRJV9Vf4HlN6oGIhWofH+O/5/VEeNjz+lh4+fDsi9V08Xm697b3SpIe+8l/aeCtoYLmsLhPbJ5RV6uOaxqcitBLr/8u09MAYEIILECZHO4b1Ge39uh3Z4Y0t2GqHvnMdSV5n49cNRYefnH8Dzr97ogk6afj4eUjV46Fmdvnz9B1l9fp7eGYHv3x63mfJSQlhtRJ0iduaFQ4HHJC0k9fP+W8FgD4qajAsnnzZjU3N6u2tlYtLS3as2dP1ut3796tlpYW1dbWavbs2dq6dWvaNTt27NA111yjmpoaXXPNNXr++eeLuTUgkCzL0gPfPaTBs6P6wMwLtf1/tek9U2tL8l6XXThFs6efr7gl9fzXKZ0dienlN34vKRFmwuGQ7vvTuZKkf/63o3pn2N65lP+iW0n69PvHFu1ef/mFmloT1el3R3ToxGlffx4AkIoILNu3b9fKlSu1Zs0a9fb2qr29XYsWLdKxY8c8rz9y5Ihuv/12tbe3q7e3V/fff7/uvvtu7dixw7mmp6dHS5Ys0R133KFf/OIXuuOOO/TZz35W//Zv/1b8TwYEyI9/eVI/f+O/VRMNa9MXPqC686pK+n726c67//Okfn7kvzU8GteMulrNufR855oFV07XjVdfqpGYpdECKiz2otu5DVM1t2GapLF1LR8eXzvDOhYApVBwYNmwYYPuvPNOLVu2TPPmzVNXV5eampq0ZcsWz+u3bt2qmTNnqqurS/PmzdOyZcu0dOlSrV+/3rmmq6tLt956q1avXq25c+dq9erVuvnmm9XV1VX0DwYERSxu6W+//0tJ0lf+uFkz6grfslyoj859jyTp2ZeP61s7D0saawelbqP+3+NrWWz5zGGx7/+zrU1Jj9ttoZcILABKoKA5LMPDw9q/f7/uu+++pMc7Ojq0d+9ez+f09PSoo6Mj6bHbbrtN27Zt08jIiKqqqtTT06NVq1alXZMtsAwNDWloaMj5/eDgYCE/St62vXREv/k9WzVRvJNnhvSfv31LdVOq9BcL55TlPW+6+lJ9uW2W/qHnqH7Zf0ZSouriNv+yOn3yhkZ998DYZNx8Tmu+++ar9JGrpuumlK3X9vqY/Ud/r2/8/4cm+iMACKClf9yspovPM/LeBQWWgYEBxWIx1dfXJz1eX1+v/v5+z+f09/d7Xj86OqqBgQHNmDEj4zWZXlOS1q1bp2984xuF3H5R/u+rJ/RKkQfKAW4rbppT8laQLRQK6cFPvE9151Vr449eUyQc0h9fmR5YJOmvbn2vdh7sU9ySzq/OvUOpbkqVPvre96Q93jz9fDVdPEXH//tdfeenb0z0RwAQQB+/vrEyAosttaxsWVbWMeJe16c+Xuhrrl69Wp2dnc7vBwcH1dTUlPH6Yi1uuVxt4715oFgXn1+jL7fNKut7hkIhdd56teY3TlNVJKzpF9R4XjfzkvP0/931YZ05O6ILz8v//CKv99v8P1r0g0P9ssROIeBcVD+tNJsF8lFQYJk+fboikUha5ePkyZNpFRJbQ0OD5/XRaFSXXHJJ1msyvaYk1dTUqKbG+y9gP33hQ+X9RwbwW8f7GnJe88ErLvblva69vE7XXl7ny2sBgFtBi26rq6vV0tKi7u7upMe7u7u1YMECz+e0tbWlXb9r1y61traqqqoq6zWZXhMAAEwuBbeEOjs7dccdd6i1tVVtbW167LHHdOzYMS1fvlzSWKvmzTff1NNPPy1JWr58uR599FF1dnbqrrvuUk9Pj7Zt26ZnnnnGec177rlHN954o/72b/9Wn/zkJ/Xd735XP/zhD/XSSy/59GMCAIBKVnBgWbJkiU6dOqW1a9eqr69P8+fP186dOzVr1ljrpK+vL2kmS3Nzs3bu3KlVq1Zp06ZNamxs1MaNG7V48WLnmgULFujZZ5/V17/+df31X/+15syZo+3bt+tDH/qQDz8iAACodCHLXgFb4QYHB1VXV6fTp09r2rRppm8HAADkId9/vzlLCAAABB6BBQAABB6BBQAABB6BBQAABB6BBQAABB6BBQAABB6BBQAABB6BBQAABB6BBQAABF7Bo/mDyh7YOzg4aPhOAABAvux/t3MN3j9nAsuZM2ckSU1NTYbvBAAAFOrMmTOqq6vL+P1z5iyheDyuEydOaOrUqQqFQr697uDgoJqamnT8+HHOKCohPufy4bMuDz7n8uBzLo9Sfs6WZenMmTNqbGxUOJx5pco5U2EJh8O6/PLLS/b606ZN4w9DGfA5lw+fdXnwOZcHn3N5lOpzzlZZsbHoFgAABB6BBQAABB6BJYeamho98MADqqmpMX0r5zQ+5/Lhsy4PPufy4HMujyB8zufMolsAAHDuosICAAACj8ACAAACj8ACAAACj8ACAAACj8CSw+bNm9Xc3Kza2lq1tLRoz549pm+pov3kJz/Rxz/+cTU2NioUCulf/uVfkr5vWZYefPBBNTY2asqUKbrpppt06NAhMzdbwdatW6cPfvCDmjp1qt7znvfoU5/6lH71q18lXcNnPXFbtmzRdddd5wzTamtr07/+67863+czLo1169YpFApp5cqVzmN81hP34IMPKhQKJf1qaGhwvm/6MyawZLF9+3atXLlSa9asUW9vr9rb27Vo0SIdO3bM9K1VrLffflvXX3+9Hn30Uc/vP/LII9qwYYMeffRRvfzyy2poaNCtt97qnBWF/OzevVtf/epX9bOf/Uzd3d0aHR1VR0eH3n77becaPuuJu/zyy/Xwww9r37592rdvn/7kT/5En/zkJ52/xPmM/ffyyy/rscce03XXXZf0OJ+1P973vvepr6/P+XXw4EHne8Y/YwsZ/dEf/ZG1fPnypMfmzp1r3XfffYbu6NwiyXr++eed38fjcauhocF6+OGHncfOnj1r1dXVWVu3bjVwh+eOkydPWpKs3bt3W5bFZ11KF110kfXEE0/wGZfAmTNnrKuuusrq7u62Fi5caN1zzz2WZfG/Z7888MAD1vXXX+/5vSB8xlRYMhgeHtb+/fvV0dGR9HhHR4f27t1r6K7ObUeOHFF/f3/SZ15TU6OFCxfymU/Q6dOnJUkXX3yxJD7rUojFYnr22Wf19ttvq62tjc+4BL761a/qYx/7mG655Zakx/ms/fPaa6+psbFRzc3N+tznPqdf//rXkoLxGZ8zhx/6bWBgQLFYTPX19UmP19fXq7+/39Bdndvsz9XrMz969KiJWzonWJalzs5OfeQjH9H8+fMl8Vn76eDBg2pra9PZs2d1wQUX6Pnnn9c111zj/CXOZ+yPZ599Vq+88opefvnltO/xv2d/fOhDH9LTTz+tq6++Wr/97W/10EMPacGCBTp06FAgPmMCSw6hUCjp95ZlpT0Gf/GZ++trX/uaXn31Vb300ktp3+Oznrj3vve9OnDggP7whz9ox44d+vKXv6zdu3c73+cznrjjx4/rnnvu0a5du1RbW5vxOj7riVm0aJHz39dee63a2to0Z84c/cM//IM+/OEPSzL7GdMSymD69OmKRCJp1ZSTJ0+mJUz4w16Nzmfun7/8y7/U9773Pb3wwgu6/PLLncf5rP1TXV2tK6+8Uq2trVq3bp2uv/56/d3f/R2fsY/279+vkydPqqWlRdFoVNFoVLt379bGjRsVjUadz5PP2l/nn3++rr32Wr322muB+N8zgSWD6upqtbS0qLu7O+nx7u5uLViwwNBdnduam5vV0NCQ9JkPDw9r9+7dfOYFsixLX/va1/Tcc8/pxz/+sZqbm5O+z2ddOpZlaWhoiM/YRzfffLMOHjyoAwcOOL9aW1v1hS98QQcOHNDs2bP5rEtgaGhIhw8f1owZM4Lxv+eyLO2tUM8++6xVVVVlbdu2zfqP//gPa+XKldb5559vvfHGG6ZvrWKdOXPG6u3ttXp7ey1J1oYNG6ze3l7r6NGjlmVZ1sMPP2zV1dVZzz33nHXw4EHr85//vDVjxgxrcHDQ8J1Xlr/4i7+w6urqrBdffNHq6+tzfr3zzjvONXzWE7d69WrrJz/5iXXkyBHr1Vdfte6//34rHA5bu3btsiyLz7iU3LuELIvP2g9/9Vd/Zb344ovWr3/9a+tnP/uZ9Wd/9mfW1KlTnX/zTH/GBJYcNm3aZM2aNcuqrq62PvCBDzjbQlGcF154wZKU9uvLX/6yZVljW+ceeOABq6GhwaqpqbFuvPFG6+DBg2ZvugJ5fcaSrO985zvONXzWE7d06VLn74dLL73Uuvnmm52wYll8xqWUGlj4rCduyZIl1owZM6yqqiqrsbHR+vM//3Pr0KFDzvdNf8Yhy7Ks8tRyAAAAisMaFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHj/DzERACGncU4FAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.plot(np.linspace(0, 50, 200), density.results['x']['mass_density'])" ] }, { "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 3 (ipykernel)", "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.10.16" }, "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 }