{ "cells": [ { "attachments": {}, "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": { "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" }, "ExecuteTime": { "end_time": "2023-06-09T12:02:39.498830621Z", "start_time": "2023-06-09T12:02:38.566620541Z" } }, "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" }, "ExecuteTime": { "end_time": "2023-06-09T12:02:39.721579406Z", "start_time": "2023-06-09T12:02:39.539285035Z" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/home/pbarletta/mambaforge/envs/guide/lib/python3.9/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": { "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" }, "ExecuteTime": { "end_time": "2023-06-09T12:02:39.729113515Z", "start_time": "2023-06-09T12:02:39.720089223Z" } }, "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" }, "ExecuteTime": { "end_time": "2023-06-09T12:02:39.792023484Z", "start_time": "2023-06-09T12:02:39.729952761Z" } }, "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": { "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" }, "ExecuteTime": { "end_time": "2023-06-09T12:02:39.809039922Z", "start_time": "2023-06-09T12:02:39.763907644Z" } }, "outputs": [ { "data": { "text/plain": "200" }, "execution_count": 5, "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" }, "tags": [ "nbval-ignore-output" ], "ExecuteTime": { "end_time": "2023-06-09T12:02:39.809831995Z", "start_time": "2023-06-09T12:02:39.807626584Z" } }, "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": { "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" }, "tags": [ "nbval-ignore-output" ], "ExecuteTime": { "end_time": "2023-06-09T12:02:39.810506364Z", "start_time": "2023-06-09T12:02:39.808001286Z" } }, "outputs": [ { "data": { "text/plain": "KeysView({'dim': 0, 'slice_volume': 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. , 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. ]), 'hist_bin_edges': array([ 0. , 0.25, 0.5 , 0.75, 1. , 1.25, 1.5 , 1.75, 2. ,\n 2.25, 2.5 , 2.75, 3. , 3.25, 3.5 , 3.75, 4. , 4.25,\n 4.5 , 4.75, 5. , 5.25, 5.5 , 5.75, 6. , 6.25, 6.5 ,\n 6.75, 7. , 7.25, 7.5 , 7.75, 8. , 8.25, 8.5 , 8.75,\n 9. , 9.25, 9.5 , 9.75, 10. , 10.25, 10.5 , 10.75, 11. ,\n 11.25, 11.5 , 11.75, 12. , 12.25, 12.5 , 12.75, 13. , 13.25,\n 13.5 , 13.75, 14. , 14.25, 14.5 , 14.75, 15. , 15.25, 15.5 ,\n 15.75, 16. , 16.25, 16.5 , 16.75, 17. , 17.25, 17.5 , 17.75,\n 18. , 18.25, 18.5 , 18.75, 19. , 19.25, 19.5 , 19.75, 20. ,\n 20.25, 20.5 , 20.75, 21. , 21.25, 21.5 , 21.75, 22. , 22.25,\n 22.5 , 22.75, 23. , 23.25, 23.5 , 23.75, 24. , 24.25, 24.5 ,\n 24.75, 25. , 25.25, 25.5 , 25.75, 26. , 26.25, 26.5 , 26.75,\n 27. , 27.25, 27.5 , 27.75, 28. , 28.25, 28.5 , 28.75, 29. ,\n 29.25, 29.5 , 29.75, 30. , 30.25, 30.5 , 30.75, 31. , 31.25,\n 31.5 , 31.75, 32. , 32.25, 32.5 , 32.75, 33. , 33.25, 33.5 ,\n 33.75, 34. , 34.25, 34.5 , 34.75, 35. , 35.25, 35.5 , 35.75,\n 36. , 36.25, 36.5 , 36.75, 37. , 37.25, 37.5 , 37.75, 38. ,\n 38.25, 38.5 , 38.75, 39. , 39.25, 39.5 , 39.75, 40. , 40.25,\n 40.5 , 40.75, 41. , 41.25, 41.5 , 41.75, 42. , 42.25, 42.5 ,\n 42.75, 43. , 43.25, 43.5 , 43.75, 44. , 44.25, 44.5 , 44.75,\n 45. , 45.25, 45.5 , 45.75, 46. , 46.25, 46.5 , 46.75, 47. ,\n 47.25, 47.5 , 47.75, 48. , 48.25, 48.5 , 48.75, 49. , 49.25,\n 49.5 , 49.75, 50. ], dtype=float32)})" }, "execution_count": 7, "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" }, "ExecuteTime": { "end_time": "2023-06-09T12:02:39.810902504Z", "start_time": "2023-06-09T12:02:39.808816049Z" } }, "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": { "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" }, "ExecuteTime": { "end_time": "2023-06-09T12:02:39.944318607Z", "start_time": "2023-06-09T12:02:39.809409947Z" } }, "outputs": [ { "data": { "text/plain": "[]" }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": "
", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiwAAAGdCAYAAAAxCSikAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA8HUlEQVR4nO3de5RU5Z3v/09d+oIK7QXTTWuDDWrA4CXpziTNpMWM2g5mchuyQnISkxPEsxiSUehxnSWSWRriBMewWD2MXJaKcZyZn/IHOsnvHCahkyjB0JkINpFhSEYnCES6Q5pJaLzQl6p9/ujeu3ZV7br2rnp20e/XWqyW6l1Vu2sF+OT7fZ7vE7IsyxIAAECAhU3fAAAAQC4EFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHgEFgAAEHhR0zfgl3g8rhMnTmjq1KkKhUKmbwcAAOTBsiydOXNGjY2NCocz11HOmcBy4sQJNTU1mb4NAABQhOPHj+vyyy/P+P1zJrBMnTpV0tgPPG3aNMN3AwAA8jE4OKimpibn3/FMzpnAYreBpk2bRmABAKDC5FrOwaJbAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWoMKcHYnp8Z/8Wq+fPGP6VgCgbAgsQIV54Zcn9Tc7D2v9D/7T9K0AQNkQWIAKM3h2RJL01tCo4TsBgPIhsAAVZiRmSZJiccvwnQBA+RBYgAozGotLkmIWgQXA5EFgASrM6HhlJU6FBcAkQmABKozTEqLCAmASIbAAFcZuCVFhATCZEFiACjMSp8ICYPIhsAAVJlFhMXwjAFBGBBagwjiLbqmwAJhECCxAhRmxtzWzhgXAJEJgASrMKLuEAExCBBagwozG2SUEYPIhsAAVhjksACYjAgtQYdglBGAyIrAAFcaZw0JLCMAkQmABKgyHHwKYjAgsQIWxdwmx6BbAZFJUYNm8ebOam5tVW1urlpYW7dmzJ+v1u3fvVktLi2prazV79mxt3bo147XPPvusQqGQPvWpTxVza8A5j9H8ACajggPL9u3btXLlSq1Zs0a9vb1qb2/XokWLdOzYMc/rjxw5ottvv13t7e3q7e3V/fffr7vvvls7duxIu/bo0aO699571d7eXvhPAkwSowyOAzAJFRxYNmzYoDvvvFPLli3TvHnz1NXVpaamJm3ZssXz+q1bt2rmzJnq6urSvHnztGzZMi1dulTr169Pui4Wi+kLX/iCvvGNb2j27NnF/TTAJEBLCMBkVFBgGR4e1v79+9XR0ZH0eEdHh/bu3ev5nJ6enrTrb7vtNu3bt08jIyPOY2vXrtWll16qO++8M697GRoa0uDgYNIvYDIYibPoFsDkU1BgGRgYUCwWU319fdLj9fX16u/v93xOf3+/5/Wjo6MaGBiQJP30pz/Vtm3b9Pjjj+d9L+vWrVNdXZ3zq6mpqZAfBahYToWFvAJgEilq0W0oFEr6vWVZaY/lut5+/MyZM/riF7+oxx9/XNOnT8/7HlavXq3Tp087v44fP17ATwBUrpEYo/kBTD7RQi6ePn26IpFIWjXl5MmTaVUUW0NDg+f10WhUl1xyiQ4dOqQ33nhDH//4x53vx8dL3tFoVL/61a80Z86ctNetqalRTU1NIbcPnBNG2SUEYBIqqMJSXV2tlpYWdXd3Jz3e3d2tBQsWeD6nra0t7fpdu3aptbVVVVVVmjt3rg4ePKgDBw44vz7xiU/oox/9qA4cOECrB0hh7xKyrES1EgDOdQVVWCSps7NTd9xxh1pbW9XW1qbHHntMx44d0/LlyyWNtWrefPNNPf3005Kk5cuX69FHH1VnZ6fuuusu9fT0aNu2bXrmmWckSbW1tZo/f37Se1x44YWSlPY4gMThh9LY1uZoJHM7FgDOFQUHliVLlujUqVNau3at+vr6NH/+fO3cuVOzZs2SJPX19SXNZGlubtbOnTu1atUqbdq0SY2Njdq4caMWL17s308BTCKjrlMPY5ZV+B9iAKhAIescqSkPDg6qrq5Op0+f1rRp00zfDlAyLd/s1qm3hyVJh9f+qaZURwzfEQAUL99/vzlLCKgw9i4hiYW3ACYPAgtQYUbjyWtYAGAyILAAFcYdWJjFAmCyILAAFWaUlhCASYjAAlSQeNxKGslPhQXAZEFgASrIiGtLs0SFBcDkQWABKshoLDmgUGABMFkQWIAKkhZYSCwAJgkCC1BB0lpCBBYAkwSBBaggqRUW1rAAmCwILEAFcU+5lWgJAZg8CCxABRmNU2EBMDkRWIAKMhpjDQuAyYnAAlSQkbRdQoZuBADKjMACVJBRBscBmKQILEAFSa2w0BICMFkQWIAKkrqGJU6FBcAkQWABKkjaLiEqLAAmCQILUEGYwwJgsiKwABWESbcAJisCC1BBUncJUWABMFkQWIAKkj6HhcQCYHIgsAAVJG0OC4EFwCRBYAEqSNocFtawAJgkCCxABUlddEtLCMBkQWABKgij+QFMVgQWoIIwmh/AZEVgASoIo/kBTFYEFqCCpI/mN3QjAFBmBBaggjCaH8BkRWABAmJ4NHe5hNH8ACYrAgsQAJteeF3XfeMHOvib01mvG2FwHIBJisACBMD+o7/X2ZG4/v1E9sCSNoeFCguASYLAAgSAXSnJFUDSdglRYQEwSRBYgACwg0qu/DGSukuIvAJgkiCwAAFgV1gsKiwA4InAAgSA0xLKEUDYJQRgsiKwAAFgt4RytXjSWkJUWABMEgQWIADs3EFLCAC8EViAAMh3l5B9+GE4NP48n1tCvcd+rwe/d0iDZ0d8fV0AmCgCCxAA+e4SGh0fHFcTjYxd73OFZcuL/6Wn9r6hHx8+6evrAsBEEViAAMh/DsvY92uqxv7o+l1hOTt+PMA7wzFfXxcAJorAAgRAYltz9uvsww9rouOBxefTmu01NKNxjoEGECwEFiAA7KCSa9fP6Pj3nZaQzxUW+/VGmEgHIGAILEAAxKx8W0KpFRZ/g4X9eqm7kQDANAILEADxeH6LbkdS17D4HFjslxshsAAIGAILEAB2hSXXHJZYiVtCFi0hAAFFYAECIO85LPHklpD/a1jGvrLoFkDQEFiAAEgsus1+nbOtuUS7hOwAlHpmEQCYRmABAqDQ05pLNTgusYaFwAIgWAgsQADku0vIPvywVIPj7ABESwhA0BBYgADId5eQXWGpjoSTnufbfTiLbgksAIKFwAIEQP5zWEpcYaElBCCgCCxAAMTzHc2fcvih33NYnNH8VFgABAyBBQiAuLNLKM8KS8m2NY+3hHwOQgAwUQQWIADymcNiWVbaWUKlmnRLhQVA0BBYgABIrGHJfM2o65uJ0fz+3oezS4g1LAAChsACBEA8jzks7hBR6pbQMBUWAAFDYAECIJ7HLqER12yU0reEqLAACBYCC2CYZVlOUMjaEnKFiOoSV1gYHAcgaAgsgGHukJJtEJy9EDYUkqoiofHn+r2teewrc1gABA2BBTDM3dbJ3hIa+15VOKxwKJT2XD9QYQEQVAQWwDB3SMneEhoLEdFISJHweIXF51wRY5cQgIAisACG5V1hGQ8R0XAoUWEp0Wh+dgkBCBoCC2CYO6Rkyx92m6YqEnYqLKUbzU+FBUCwEFgAw9xtnWwVFjtEjLWEcl9f1L1wlhCAgCKwAIa52zrZKiYj9hqWki66HX8vzhICEDAEFsCw5DUsma+zR/NXRUKlCyxxKiwAgqmowLJ582Y1NzertrZWLS0t2rNnT9brd+/erZaWFtXW1mr27NnaunVr0vefe+45tba26sILL9T555+vG264Qf/4j/9YzK0BFSd5DUseFRbXGpbStYSosAAIloIDy/bt27Vy5UqtWbNGvb29am9v16JFi3Ts2DHP648cOaLbb79d7e3t6u3t1f3336+7775bO3bscK65+OKLtWbNGvX09OjVV1/VV77yFX3lK1/RD37wg+J/MqBCJG9rzmMNS7iEFRZ2CQEIqIIDy4YNG3TnnXdq2bJlmjdvnrq6utTU1KQtW7Z4Xr9161bNnDlTXV1dmjdvnpYtW6alS5dq/fr1zjU33XSTPv3pT2vevHmaM2eO7rnnHl133XV66aWXiv/JgAqRf0sofZeQ30tNEoPjqLAACJaCAsvw8LD279+vjo6OpMc7Ojq0d+9ez+f09PSkXX/bbbdp3759GhkZSbvesiz96Ec/0q9+9SvdeOONGe9laGhIg4ODSb+ASpTvLqERj11C/m9rTrxutvYUAJRbQYFlYGBAsVhM9fX1SY/X19erv7/f8zn9/f2e14+OjmpgYMB57PTp07rgggtUXV2tj33sY/r7v/973XrrrRnvZd26daqrq3N+NTU1FfKjAIERK7AlVI7R/BLnCQEIlqIW3YbG/7K0WZaV9liu61Mfnzp1qg4cOKCXX35Zf/M3f6POzk69+OKLGV9z9erVOn36tPPr+PHjRfwkgHlJLaEsS0fsllDSaH6fqyDu8MR5QgCCJFrIxdOnT1ckEkmrppw8eTKtimJraGjwvD4ajeqSSy5xHguHw7ryyislSTfccIMOHz6sdevW6aabbvJ83ZqaGtXU1BRy+0AgWXlWWBItoUSFxc/AYllW0qRdKiwAgqSgCkt1dbVaWlrU3d2d9Hh3d7cWLFjg+Zy2tra063ft2qXW1lZVVVVlfC/LsjQ0NFTI7QEVKZbvaP7xnTtV4ZBrNL9/95H63sxiARAkBVVYJKmzs1N33HGHWltb1dbWpscee0zHjh3T8uXLJY21at588009/fTTkqTly5fr0UcfVWdnp+666y719PRo27ZteuaZZ5zXXLdunVpbWzVnzhwNDw9r586devrppzPuPALOJXkffhh3L7r1v8KS+lpUWAAEScGBZcmSJTp16pTWrl2rvr4+zZ8/Xzt37tSsWbMkSX19fUkzWZqbm7Vz506tWrVKmzZtUmNjozZu3KjFixc717z99ttasWKFfvOb32jKlCmaO3eu/umf/klLlizx4UcEgi3/s4QSg+NKseg29aVGqLAACJCCA4skrVixQitWrPD83lNPPZX22MKFC/XKK69kfL2HHnpIDz30UDG3AlS8pLOEsraE7F1CrgqLr4El+bWYxQIgSDhLCDAs79H8cddofrvCUsKWEGtYAAQJgQUwLJ7nGhanwhIJKVyCwXHpLSEqLACCg8ACGJb3HBZ7DUu4NNua01tCVFgABAeBBTAs30m3XruE/KywWCn5hEW3AIKEwAIY5i5k5DWHJWlwXPZ1LwXdB9uaAQQYgQUwLJ60SyiPSbeuXUJjz/f/PqTEmhkACAICC2BY3ocfeuwSkvxrC6WGpRHWsAAIEAILYJh7l1D2llCiwhJ2/cn1a+Ft+mh+KiwAgoPAAhiW92j+WPqi21zPKQRzWAAEGYEFMCyeZ0soFrcPP0wsuh173K/Akvz7YQILgAAhsACGuYNCtmUjdocmnLro1qdckTrmn5YQgCAhsACG5dsSsgNFJKTkRbelWsPColsAAUJgAQzLvyU0HljCIYXDZdglRIUFQIAQWADDkissma+zw0xovLrinNjMolsAkwCBBTAslrStOUtLyEpUWKREW8ivCkvqe4/6OPYfACaKwAIYltwSynyd0xIaDyp+n9jMLiEAQUZgAQxL2iWUbQ2La5eQlAgupWsJUWEBEBwEFsAwd4UkW7XE2SU0/qc27HNLKHVTEGtYAAQJgQUwzF3ZyFYssYNJ2GkJlbbCMsIaFgABQmABDMt3DkssddFt2K6w+HMf7BICEGQEFsCwwgfHjVdY/G4JpbwMc1gABAmBBTDMSlp0m/k6u8LiLLoN288pUUuICguAACGwAIbFktawZJvDMvbVrqz4vUsobQ4LFRYAAUJgAQwrepdQuMQtIc4SAhAgBBbAsHieo/lTdwn5PZo/NfhQYQEQJAQWwLDUQwcztYUyj+b35z7SdglRYQEQIAQWwLDUqkqmKkv6aH6/zxJK/j27hAAECYEFMCyeEjgytXjSdgmVeDQ/u4QABAmBBTAstSWUKYAkFt2WpsKS+jKsYQEQJAQWwLC0CkuGwoZTYQklz2FJDTx+3QcVFgBBQmABDEutkGSusIx9TV10mxo0ipW+6JYKC4DgILAAhqUvus2whiVl0W2oxKP5OUsIQJAQWADDUgNKpvxhXzeeU3yfw5K+6JYKC4DgILAAhqVWSEzNYUl9X9awAAgSAgtgWPouoQzXpe0S8n5+sez3tSs4rGEBECQEFsCw1EWzmdakZBrNn+3AxELYr189vv2ICguAICGwAIalrh3J3BIa++pUWHxfdDv2OjXRsb8WmMMCIEgILIBhqYWMfEfzR0o0mr86GpHEWUIAgoXAAhiWvkso12j+sd+XajS/XWFhlxCAICGwAIblPzgu02h+f+7Dvo1EYKHCAiA4CCyAYWm7hHKM5ndaQvYaFp8rLNWsYQEQQAQWwLDURbZeFRbLspw1JqGUNSy+jeaPp7SEWMMCIEAILIBh+bSE3JeU+rRmu8JiWf69NgBMFIEFMCyfXULu4JBoCdnXl6YlJI2tYzl04rR+O3jWl/cAgGIRWADD8pnD4r7G3iXkd4XFfl97cJwk/eb37+oTj/5UX37y5768BwAUi8ACGJbeEsp+TdrgOJ9H89eMz2GRpNdPvqVY3NKJP7zry3sAQLEILIBhqRUWr4qJO5SEU3YJ+bbodvw9onavSdLvzoy1gobZ4gzAMAILYFg+g+PiXhUWn+ewuA9XrBoPLb87MyRJGhqN+3ZmEQAUg8ACGJZaUfHKBZ6Lbn0+rdl+mUgopOj4QpmT44HFsji9GYBZBBbAsNRxJ14VlqSWUDi5JeRX5cN+31Ao5LSF7MAijVVZAMAUAgtgWNqkW8+W0NhXux0klW4OSzgkVUXsCktiO/MwgQWAQQQWwLB8dgnZIcaVV0o2mj8cSqxh+e2gu8IS8+V9AKAYBBbAsHxG89uhxt4hJPk/mt++j3BYzhqWU28lAgsVFgAmEVgAw9IPP8w8OM67JeTTfYy/TshVYXHfCmtYAJhEYAEMK2Q0f8RdYbHnsPjcEoqEQopG0v9qoMICwCQCC2BYakUl22j+cAkX3VqudTJR92KZcaxhAWASgQUwLH2XkMc1HruE/F90O/Z1rCWU/lcDLSEAJhFYAMPymXTrveh2/HqfR/OHXXNY3AgsAEwisACGpQYOr4pJYtFt4rFyzGFxGxohsAAwh8ACGJYaULzWsGRbdOv7HBbXWUJuHIAIwCQCC2BY2mh+j1zguejW79OaXW0new6L29AIi24BmENgAQyzqyfV420Yz9H8VvoaFqcl5NOZhMktISosAIKFwAIYZoeRqMewNpv3LqHk5/t1H5krLAQWAOYQWADDnMASznz6cmKXUOKxko3mD8lzlxAVFgAmEVgAw+wwYu/Myb5LqPSnNYdCIac9JUnV0bH/psICwCQCC2CYswMonK0l5DGHxefR/LEMc1gaptVKkoZjLLoFYA6BBTDMDih2hcWzJVSGCovlmvXiPkvIDixUWACYVFRg2bx5s5qbm1VbW6uWlhbt2bMn6/W7d+9WS0uLamtrNXv2bG3dujXp+48//rja29t10UUX6aKLLtItt9yin//858XcGlBx7ApJ4oRkj5ZQPD2wJOaw+HQf7tOaXe/TUGdXWAgsAMwpOLBs375dK1eu1Jo1a9Tb26v29nYtWrRIx44d87z+yJEjuv3229Xe3q7e3l7df//9uvvuu7Vjxw7nmhdffFGf//zn9cILL6inp0czZ85UR0eH3nzzzeJ/MqBC2BUSu6rhNYfFezS/z3NYklpCrgpLHRUWAOYVHFg2bNigO++8U8uWLdO8efPU1dWlpqYmbdmyxfP6rVu3aubMmerq6tK8efO0bNkyLV26VOvXr3eu+ed//metWLFCN9xwg+bOnavHH39c8XhcP/rRj4r/yYAKkbpLyHsOy9jXciy6de8SioRDmn5BtSQqLADMKiiwDA8Pa//+/ero6Eh6vKOjQ3v37vV8Tk9PT9r1t912m/bt26eRkRHP57zzzjsaGRnRxRdfnPFehoaGNDg4mPQLqESJCku2wOKxrblUo/ldu4QuOq9aNdGIJGlolEW3AMwpKLAMDAwoFoupvr4+6fH6+nr19/d7Pqe/v9/z+tHRUQ0MDHg+57777tNll12mW265JeO9rFu3TnV1dc6vpqamQn4UIBAsy3IqG/awtrx3CZXotOZQKHEvF51XpZrxbc3DnNYMwKCiFt2GQslDpSzLSnss1/Vej0vSI488omeeeUbPPfecamtrM77m6tWrdfr0aefX8ePHC/kRgEBwZ418RvNHPM4S8q/CIuc97GrPRedXJ+awEFgAGBQt5OLp06crEomkVVNOnjyZVkWxNTQ0eF4fjUZ1ySWXJD2+fv16fetb39IPf/hDXXfddVnvpaamRjU1NYXcPhA47nCSfTS/xy6hEi66TbSEqlwtIQILAHMKqrBUV1erpaVF3d3dSY93d3drwYIFns9pa2tLu37Xrl1qbW1VVVWV89i3v/1tffOb39T3v/99tba2FnJbQMVyL5iNZpvD4tES8rvC4h7Nf+PVl2puw1R9+v2XUWEBEAgFt4Q6Ozv1xBNP6Mknn9Thw4e1atUqHTt2TMuXL5c01qr50pe+5Fy/fPlyHT16VJ2dnTp8+LCefPJJbdu2Tffee69zzSOPPKKvf/3revLJJ3XFFVeov79f/f39euutt3z4EYHgSqqwZKmYZB/N79O9uOawvLdhqr6/8kb96fwZrGEBEAgFtYQkacmSJTp16pTWrl2rvr4+zZ8/Xzt37tSsWbMkSX19fUkzWZqbm7Vz506tWrVKmzZtUmNjozZu3KjFixc712zevFnDw8P6zGc+k/ReDzzwgB588MEifzQg+JIqLOHMg+DsUOI1mt+rIlPUvVjpVRzJdZYQu4QAGFRwYJGkFStWaMWKFZ7fe+qpp9IeW7hwoV555ZWMr/fGG28UcxtAxXMPictvNH/isfGNPCUZze9GhQVAEHCWEGCQuyVU/Gh+/09rdmMNC4AgILAABrnDRiTLHJa4R7umlLuE3OxdQlRYAJhEYAEMiscTO3OcQXD57hIKl6bCEk4Zj1TDGhYAAUBgAQxyL3QNO4to06/z2iVkt4S8DksshpWxwpJYw+LXAl8AKBSBBTDIqZyEQ87aEa9FtJ67hHw+/NB+ndQB1PYalrgljfr0XgBQKAILYJBdsIiEQk4rJvto/sRj/o/mT6/iSIk1LBLrWACYQ2ABDHKP3LcDiLnR/GNfM81hkdgpBMAcAgtgUMxyL7rNPAgu22nNpRjN7xYJh5yhdlRYAJhCYAEMiietYRl/rNDTmn2usHidos60WwCmEVgAg5wJtqH8WkIm5rBITLsFYB6BBTDI3pIcDrsW3Xoefjj21bPC4tei27h3S0hi2i0A8wgsgEFxzwpL5paQO0zYg+P82mnsLLr1SCz2TiECCwBTCCyAQe7dP6F8WkKuMBH1eQ5LtpYQa1gAmEZgAQyy2zmhPEfzRzIMjvNjAm2m0fwSa1gAmEdgAQyKe8xhyXc0f9T1335UWTKN5pdYwwLAPAILYJCzmDbkbgnlN4cl6hp768fI/Lir2pOKCgsA0wgsgEHutSl2wcSrWpKrwuJHYPEKRbZqFt0CMIzAAhjkvUso/Tqv0fxJgSU28SDhnGvkuUuICgsAswgsgEFeFRbv0fxjX70Gx0nSSMy/llD2OSzsEgJgBoEFMMh9llBirkp+pzWHQiFVRfzb2pxtND8VFgCmEVgAgywrfZeQV/bINCPFrrKM+NASymc0P2tYAJhCYAEMcrd6nNH8ee4SkqSq8NgfYT8W3dpv6z2HZWzRLRUWAKYQWACDYh5zWLzPEkpfdCtJEaclNPEgkX2XEGtYAJhFYAEMcu8SKnQ0vyRFxyss/i66ZQ0LgOAhsAAG2UEkFFKOltDY10hqS2i8wjLqS2AZ+xr2+FuhOsIaFgBmEVgAg9ytHrvdk300f/Lj9nNGfWgJZRvNX1NFhQWAWQQWwCB3YCl0NL8kVUX8W3SbdQ4LFRYAhhFYAIPy3SWUadGtPe3Wz5aQ1xwWRvMDMI3AAhjkdVqz10gVr9H87t/70RJy7iXrHBZ2CQEwg8ACGJQ06TbLaP5MO3icllCJdwlVs0sIgGEEFsAg99qUbGtY4h5nCUlS1N4l5Oto/vTvMekWgGkEFsAg92j+SLY5LBl2CSXWsPg4mt9j1S0VFgCmEVgAg5JOax7/01jILiFncFyZRvOzhgWAKQQWwCB76UkklFh0m30Oi3dLyJfR/PmsYfGhkgMAxSCwAAa5dwmFnF1CWSosGbY1l2s0/9AIgQWAGQQWwCC7qpF7NL/3luOoT7uELMvK0RKiwgLALAILYJD78MOiWkJhf1pC7vf0rrCMr2GhwgLAEAILYJDX4LiCFt1G/Dmt2f2erGEBEEQEFsAgZzR/ONdo/rGvqRWWqnDmdS+FcD895PG3gt0SisUtX7ZQA0ChCCyAQTGPlpBXsSTXac0jE2wJ5VthkaiyADCDwAIYFHdaPXLmsHiN5rcrKKkHE/q16NYdWLKdJSSxjgWAGQQWwCD3dNnso/m9dwlV+TSaP6kl5LFLKBoJOy0rKiwATCCwAAZ5tYS8ujuxDLuEIj6N5s/VEpLYKQTALAILYJB7l1Ak6y6hsa8ZT2ueYIXFcmUQrzksUmIdC+P5AZhAYAEM8tolVMgclkSFpbTbmiW2NgMwi8ACGJQYh59YUBvLNuk25U+sva15dIK7hNzvmSGvJN7Lh2MAAKBQBBbAoORJt8mPJV1XpsFx7uCUqipqvxcVFgDlR2ABDHIfamgfbOjVEsq16Nav0fyZ2kFSYr2MHwctAkChCCyAQXlXWDKcpOxsa/atwpI5sCROhqbCAqD8CCyAQUkVlqxzWORc5xYN+7NLyH56lrziLLqd6HoZACgGgQUwyC5WJB1+mG0OS9oaFn8W3WZaI5P0XuNhaXiUlhCA8iOwAAa5F7vmdVpzyp9Yu8Li56LbTBIzX6iwACg/AgtgkHvtSCjDGpa4q92TqcLi12nNqS0nt8SiWwILgPIjsAAGxVyTbhMVlpRr3AcTpq1h8WchbD6Lbu0FviO0hAAYQGABDHJPsM10WrO7epK26Nan05qtPFpCzswXWkIADCCwAAbFXItdIxkqLO4WUVpLKOxzSyhLhaXaDiyjBBYA5UdgAQxyH2oYyrDo1h1GMraEJrpLaPw9M025lVwzXyYYjgCgGAQWwCDLaQkl2jGp1ZJ40knK3qc1T7TCkqj0ZL7Gbglx+CEAEwgsgEEx12JXO4yk7mqOZ1l0G3EW3U50DYv367tV+bReBgCKQWABDPLeJZTSEnL9PjVPOIPjyrlLiAoLAAMILIBBXruEMs1hCXmcpOxXSyif0fwcfgjAJAILYFAsnljsmmsOS+oOIcnVEvJp0W3W0fxUWAAYRGABDLLDSSRpDYv3LiGvKbRV4fLNYalm0i0AgwgsgEHxeP67hLJVWCa61TiW4TRoN1pCAEwisAAGxZLOEsrREvKqsJRx0S0tIQAmEVgAg+JJu4QSj7vbQtlmpPg1mj+f05qrnfcisAAoPwILYJBTYQmHkioo7ipLPEuFJepTS8jOR1krLD7NfAGAYhBYAIPiHqP5peStzVkDizMuvwyj+aMsugVgTlGBZfPmzWpublZtba1aWlq0Z8+erNfv3r1bLS0tqq2t1ezZs7V169ak7x86dEiLFy/WFVdcoVAopK6urmJuC6g4cdeWZXceiXu0hLzCRNTeJeTb4YeZr7F3JBFYAJhQcGDZvn27Vq5cqTVr1qi3t1ft7e1atGiRjh075nn9kSNHdPvtt6u9vV29vb26//77dffdd2vHjh3ONe+8845mz56thx9+WA0NDcX/NECFSWxZTm7HuAsm2XYJ2W0ay5rY8Dj3WppMqqIcfgjAnIIDy4YNG3TnnXdq2bJlmjdvnrq6utTU1KQtW7Z4Xr9161bNnDlTXV1dmjdvnpYtW6alS5dq/fr1zjUf/OAH9e1vf1uf+9znVFNTU/xPA1SYWFKFxbsllG2XkN0SkiZW+cinJWRXc4ZHqbAAKL+CAsvw8LD279+vjo6OpMc7Ojq0d+9ez+f09PSkXX/bbbdp3759GhkZKfB2E4aGhjQ4OJj0C6g07spGKEdLKOzxp9WejeK+rqj7yKclFPGn/QQAxSgosAwMDCgWi6m+vj7p8fr6evX393s+p7+/3/P60dFRDQwMFHi7CevWrVNdXZ3zq6mpqejXAkxJnOGTxy6hLIPjpIltbc5nDkt1lDksAMwpatFtatnYsqyspWSv670eL8Tq1at1+vRp59fx48eLfi3AFK/TmqUMc1iybGuWJnaeUD6j+WkJATApWsjF06dPVyQSSaumnDx5Mq2KYmtoaPC8PhqN6pJLLinwdhNqampY74KKl3mXkOu/45krLHZlJha3fGkJZd3WTEsIgEEFVViqq6vV0tKi7u7upMe7u7u1YMECz+e0tbWlXb9r1y61traqqqqqwNsFzi3u9SnusOAOH9kW3UrugW7FVz5iWUKRrYrR/AAMKrgl1NnZqSeeeEJPPvmkDh8+rFWrVunYsWNavny5pLFWzZe+9CXn+uXLl+vo0aPq7OzU4cOH9eSTT2rbtm269957nWuGh4d14MABHThwQMPDw3rzzTd14MABvf766z78iEBwpQ6FszOJ92h+7zBR5cN4fmcNS5a/Efx4HwAoVkEtIUlasmSJTp06pbVr16qvr0/z58/Xzp07NWvWLElSX19f0kyW5uZm7dy5U6tWrdKmTZvU2NiojRs3avHixc41J06c0Pvf/37n9+vXr9f69eu1cOFCvfjiixP48YBgswspdmUjHAopbllJLSE7u2SqsPhxYnNeo/nHKyzDVFgAGFBwYJGkFStWaMWKFZ7fe+qpp9IeW7hwoV555ZWMr3fFFVck/T9KYLJInWIbDoekuOW9rTlDlqjyYTx/PnNYOPwQgEmcJQQYlDph1g4lXoPjvHYJSa7x/BNqCSnp/T3fJ2KP5uf/XAAoPwILYFAsZcaK3ZJxFxyz7RKS/GkJxXOsk5FYdAvALAILYFDqFFs7MHjtEspUYXFaQj6M5s86OC7C4YcAzCGwAIZYlqWh8SFsNdGIJDnj+b3WsGSqsER9mI9SSEsoPsGDFgGgGAQWwJAh18TYKdVjgcWucHiO5s8xh6XUo/mrfDpoEQCKQWABDDk7EnP+uzY69kfRDiXJc1jGvmZcdGuvLfFjNH8ec1gkAguA8iOwAIacHRn7Rz8aDjntlsQuocR1iUW33q8TGU8ZMR92CeUzml9ieByA8iOwAIa8O15hqa2KOI+FnJZQIhDkaglVhSc+hyXXNF37/e1vU2EBUG4EFsCQsx6Bxc4kXruEMlU/nJaQD2tYMlVxbHaVZYRFtwDKjMACGJIILIk/hsXMYbEHx01k504+o/mlRDVnZJQKC4DyIrAAhthrWJIrLOktoVg8xy4hHwa65TOaX5KqovYWagILgPIisACGeFZYxv8zeTS//b3SVVjymcMiJVpCw6O0hACUF4EFMMQOLFM8KyyJ63LtErLnsExkXUk+c1gkfxb4AkAxCCyAIWdHvRbdesxhyXX4oR+j+eO557BIiZYQu4QAlBuBBTDEXsNij+WXEqP5k3YJ5Vh0a7dp/GkJZa+w2NUcWkIAyo3AAhjy7nDmXUKeLaEMFRb7cT+2NedsCUVYdAvADAILYIjdEpriMYfFShocN/69Ep7W7Izmz3cOCy0hAGVGYAEMyb6tOXFdzMreEoo4C2FLO5pfSoSjiVRzAKAYBBbAkGyD45JG88ezVz/sbc0TadPk2xKKUmEBYAiBBTDEc1uz5xyW7LuEEi2h4qseThUnx98I1fYaFiosAMqMwAIYYgeWmhyTbnON5o84FZYyjOYfD0fDVFgAlBmBBTDEaw2Lc1qzKw/kGs3vx6JbOxTlWsMSpcICwBACC2DIu55rWMa+FtISivpQYcl3NH81a1gAGEJgAQzxWsMSyTqaP9ek29LPYfHjoEUAKAaBBTBkKMu25oJG8ztnCZVzDgstIQDlRWABDEmcJZT4YxhyWkKJ6+xiRuYKy8RH8+cKRbYqKiwADCGwAIY4o/mj6RWWmKvCYuXYcmxXWCbWEkp+/0yc0fwEFgBlRmABDHEqLNXpc1iSWkLxHC0hH6oe+baE7AW+w7SEAJQZgQUwxNnW7FFh8dwllKklNJ4yJnRa83jWyTmaPzrxLdQAUAwCC2DI2WynNbvyQM5dQuNVj5EJbWvOb5cQ25oBmEJgAQxxTmuuTj+tObnCMv69HC2hCQ2Oy3MOix/hCACKQWABDIjFLWdrcK6WUKLC4v1a/gyOyz5N12a3hEZGqbAAKC8CC2CAPTROyjCaP2lbc/Yw4U+FJb/R/FU+hCMAKAaBBTDgXVdgqYlObDS/PRtlQotu82wJcfghAFMILIABzknN0XBSEAkXMZrfPq15ItNn8x/NP/5etIQAlBmBBTDA66RmKdH2cc9hieeqsNiD48owmt/eJURLCEC5EVgAA856nNQsuUbzx9N3CWWusNiBpfRzWDj8EIApBBbAgERgSa6wJEbzJx6LO5NuvV8r6ozLn/hZQjl3CTGHBYAhBBbAALslNCUtsIx99RzNn6H6UeXDLiFOawYQdAQWwABn0W2GCovXLqFM1Q9fWkLjT825rdmHcAQAxSCwAAbY25pro6lrWArfJVTlw0LYfHcJ2e/F4YcAyo3AAhiQaQ3LeB4oaA6LffjhRNaV5D2anwoLAEMILIABZ0czrWGxtzUnHsv38MOJndbM4YcAgi1q+gaAyWgo47bmkPP9u57epyvfc0HONSyJqocPLaEcJZYoi24BGEJgAQx4dzjTtuaxr7/4zWnt/s/f6UeHf6sZdVPGvpcrsExgcFw8711CzGEBYAYtIcCAs6PZ57Ac//07ksbWlrz5h3cl5W4Jxa3kgXOFSKxhYQ4LgGAisAAGZBrNb1c43vz9u2nPyTw4LhEyRoqsshQ6h2Ui7ScAKAaBBTAg02h+u+0z5HG4YKbqR9SVMopdeFvoHBZOawZQbgQWwIB3c4zm95Jx0a2r9FLsYthc03Rtfsx8AYBiEFgAA4ZyjOa32duIx76Xu8JS7HwUy9mJlP06O7DE4lbR62UAoBgEFsCAjC2hlFBy2/wG578zVVjC4ZATdErdEvJjvQwAFIPAAhiQaZdQamBY0trk/HemXUKSaz5K0YGlsMFxErNYAJQXgQUwwJ7DUhPN3BK68LwqfXj2xaqbUiVJOq8m+Vq3qvEnxooMEXmP5veh/QQAxWBwHGCAva15SnXqWUKJQHDZhVMUjYS16X98QK+fPKPZ08/P+Hr28ya+rTl7YomEQwqFxo4OYKcQgHIisAAGOC2hDKc1S1LjhWMTbj9y1XR95KrpWV9vovNR8t0lFAqFVBUOazgWpyUEoKxoCQEGnM0xml8aq7Dka6Lj+fMdzS8lZrHQEgJQTgQWwAD7tOZsc1gKCizhiVVY7NOhcx1+KElVUcbzAyg/AgtggL2tOdsclssuCmaFxQ5HtIQAlBOBBSgzy7JyjuaXCq2wjD2v99gfirqnfOewSFI1JzYDMIDAApTZcCzuBISaLC2hxgICy+3XzpAkPfR/D+uJPb8u+J7yncMiuWa+UGEBUEYEFqDM7C3NUuaWUE00rOkXVOf9mqtuuVrLPtIsaSy0fO8XJwq6J3vMfrbhdLYqKiwADCCwAGU2NN4OCocS//jb7ArHZRdOyas94zwvHNKaj83T/1xwhSTpu71vFnRPiZZQ7msnuoUaAIrBHBagzOwKS21VJC2U2L8vpB3kfu5nWi7XU3vf0M9+fUojsbgTLiRp449e06YXXlfcslQdCetbf36tPnnDZZIKawnZr/nvJ05r5fYDOv3usCTplnn12vLFloLvGwDyQYUFKLN3R7xnsEjStZfVqToS1sKrLy3qta+ZMU0Xn1+tt4djSQtw43FL3/npEQ2Njg18e3s4pm/+n//Q20OjY993tjXnfg+7KvTkS0c08NaQRmKWRmKW/vXf+3X8v98p6r4BIBcCC1Bmp94ekpS+fkWS/qj5Yh38RofuunF2Ua8dDoe0YM4lkqSXXvud8/ihE4P6/TsjuqAmqj3/+6O64pLzNPDWsJ7Yc0RS/qP5pcSi25Nnxn6OrV/8gN4/80JJ0p7XBoq6bwDIhcAClNEv+wd1z7MHJElX1V/geU3qgYiFah8f47/n9UR42PP6WHj58OyL1XTxebr3tvdKkh77yX9p4K2hguawuE9snlFXq45rGpyK0Euv/y7T0wBgQggsQJkc7hvUZ7f26HdnhjS3Yaoe+cx1JXmfj1w1Fh5+cfwPOv3uiCTpp+Ph5SNXjoWZ2+fP0HWX1+nt4Zge/fHreZ8lJCWG1EnSJ25oVDgcckLST18/5bwWAPipqMCyefNmNTc3q7a2Vi0tLdqzZ0/W63fv3q2WlhbV1tZq9uzZ2rp1a9o1O3bs0DXXXKOamhpdc801ev7554u5NSCQLMvSA989pMGzo/rAzAu1/X+16T1Ta0vyXpddOEWzp5+vuCX1/NcpnR2J6eU3fi8pEWbC4ZDu+9O5kqR//rejemfY3rmU/6JbSfr0+8cW7V5/+YWaWhPV6XdHdOjEaV9/HgCQiggs27dv18qVK7VmzRr19vaqvb1dixYt0rFjxzyvP3LkiG6//Xa1t7ert7dX999/v+6++27t2LHDuaanp0dLlizRHXfcoV/84he644479NnPflb/9m//VvxPBgTIj395Uj9/479VEw1r0xc+oLrzqkr6fvbpzrv/86R+fuS/NTwa14y6Ws259HznmgVXTteNV1+qkZil0QIqLPai27kNUzW3YZqksXUtHx5fO8M6FgClUHBg2bBhg+68804tW7ZM8+bNU1dXl5qamrRlyxbP67du3aqZM2eqq6tL8+bN07Jly7R06VKtX7/euaarq0u33nqrVq9erblz52r16tW6+eab1dXVVfQPBgRFLG7pb7//S0nSV/64WTPqCt+yXKiPzn2PJOnZl4/rWzsPSxprB6Vuo/7f42tZbPnMYbHv/7OtTUmP222hlwgsAEqgoDksw8PD2r9/v+67776kxzs6OrR3717P5/T09KijoyPpsdtuu03btm3TyMiIqqqq1NPTo1WrVqVdky2wDA0NaWhoyPn94OBgIT9K3ra9dES/+T1bNVG8k2eG9J+/fUt1U6r0FwvnlOU9b7r6Un25bZb+oeeoftl/RlKi6uI2/7I6ffKGRn33wNhk3HxOa7775qv0kaum66aUrdf2+pj9R3+vb/z/hyb6IwAIoKV/3Kymi88z8t4FBZaBgQHFYjHV19cnPV5fX6/+/n7P5/T393tePzo6qoGBAc2YMSPjNZleU5LWrVunb3zjG4XcflH+76sn9EqRB8oBbitumlPyVpAtFArpwU+8T3XnVWvjj15TJBzSH1+ZHlgk6a9ufa92HuxT3JLOr869Q6luSpU++t73pD3ePP18NV08Rcf/+11956dvTPRHABBAH7++sTICiy21rGxZVtYx4l7Xpz5e6GuuXr1anZ2dzu8HBwfV1NSU8fpiLW65XG3jvXmgWBefX6Mvt80q63uGQiF13nq15jdOU1UkrOkX1HheN/OS8/T/3fVhnTk7ogvPy//8Iq/32/w/WvSDQ/2yxE4h4FxUP600mwXyUVBgmT59uiKRSFrl4+TJk2kVEltDQ4Pn9dFoVJdccknWazK9piTV1NSopsb7L2A/feFD5f1HBvBbx/sacl7zwSsu9uW9rr28TtdeXufLawGAW0GLbqurq9XS0qLu7u6kx7u7u7VgwQLP57S1taVdv2vXLrW2tqqqqirrNZleEwAATC4Ft4Q6Ozt1xx13qLW1VW1tbXrsscd07NgxLV++XNJYq+bNN9/U008/LUlavny5Hn30UXV2duquu+5ST0+Ptm3bpmeeecZ5zXvuuUc33nij/vZv/1af/OQn9d3vflc//OEP9dJLL/n0YwIAgEpWcGBZsmSJTp06pbVr16qvr0/z58/Xzp07NWvWWOukr68vaSZLc3Ozdu7cqVWrVmnTpk1qbGzUxo0btXjxYueaBQsW6Nlnn9XXv/51/fVf/7XmzJmj7du360Mf+pAPPyIAAKh0IcteAVvhBgcHVVdXp9OnT2vatGmmbwcAAOQh33+/OUsIAAAEHoEFAAAEHoEFAAAEHoEFAAAEHoEFAAAEHoEFAAAEHoEFAAAEHoEFAAAEHoEFAAAEXsGj+YPKHtg7ODho+E4AAEC+7H+3cw3eP2cCy5kzZyRJTU1Nhu8EAAAU6syZM6qrq8v4/XPmLKF4PK4TJ05o6tSpCoVCvr3u4OCgmpqadPz4cc4oKiE+5/Lhsy4PPufy4HMuj1J+zpZl6cyZM2psbFQ4nHmlyjlTYQmHw7r88stL9vrTpk3jD0MZ8DmXD591efA5lwefc3mU6nPOVlmxsegWAAAEHoEFAAAEHoElh5qaGj3wwAOqqakxfSvnND7n8uGzLg8+5/Lgcy6PIHzO58yiWwAAcO6iwgIAAAKPwAIAAAKPwAIAAAKPwAIAAAKPwJLD5s2b1dzcrNraWrW0tGjPnj2mb6mi/eQnP9HHP/5xNTY2KhQK6V/+5V+Svm9Zlh588EE1NjZqypQpuummm3To0CEzN1vB1q1bpw9+8IOaOnWq3vOe9+hTn/qUfvWrXyVdw2c9cVu2bNF1113nDNNqa2vTv/7rvzrf5zMujXXr1ikUCmnlypXOY3zWE/fggw8qFAol/WpoaHC+b/ozJrBksX37dq1cuVJr1qxRb2+v2tvbtWjRIh07dsz0rVWst99+W9dff70effRRz+8/8sgj2rBhgx599FG9/PLLamho0K233uqcFYX87N69W1/96lf1s5/9TN3d3RodHVVHR4fefvtt5xo+64m7/PLL9fDDD2vfvn3at2+f/uRP/kSf/OQnnb/E+Yz99/LLL+uxxx7Tddddl/Q4n7U/3ve+96mvr8/5dfDgQed7xj9jCxn90R/9kbV8+fKkx+bOnWvdd999hu7o3CLJev75553fx+Nxq6GhwXr44Yedx86ePWvV1dVZW7duNXCH546TJ09akqzdu3dblsVnXUoXXXSR9cQTT/AZl8CZM2esq666yuru7rYWLlxo3XPPPZZl8b9nvzzwwAPW9ddf7/m9IHzGVFgyGB4e1v79+9XR0ZH0eEdHh/bu3Wvors5tR44cUX9/f9JnXlNTo4ULF/KZT9Dp06clSRdffLEkPutSiMVievbZZ/X222+rra2Nz7gEvvrVr+pjH/uYbrnllqTH+az989prr6mxsVHNzc363Oc+p1//+teSgvEZnzOHH/ptYGBAsVhM9fX1SY/X19erv7/f0F2d2+zP1eszP3r0qIlbOidYlqXOzk595CMf0fz58yXxWfvp4MGDamtr09mzZ3XBBRfo+eef1zXXXOP8Jc5n7I9nn31Wr7zyil5++eW07/G/Z3986EMf0tNPP62rr75av/3tb/XQQw9pwYIFOnToUCA+YwJLDqFQKOn3lmWlPQZ/8Zn762tf+5peffVVvfTSS2nf47OeuPe+9706cOCA/vCHP2jHjh368pe/rN27dzvf5zOeuOPHj+uee+7Rrl27VFtbm/E6PuuJWbRokfPf1157rdra2jRnzhz9wz/8gz784Q9LMvsZ0xLKYPr06YpEImnVlJMnT6YlTPjDXo3OZ+6fv/zLv9T3vvc9vfDCC7r88sudx/ms/VNdXa0rr7xSra2tWrduna6//nr93d/9HZ+xj/bv36+TJ0+qpaVF0WhU0WhUu3fv1saNGxWNRp3Pk8/aX+eff76uvfZavfbaa4H43zOBJYPq6mq1tLSou7s76fHu7m4tWLDA0F2d25qbm9XQ0JD0mQ8PD2v37t185gWyLEtf+9rX9Nxzz+nHP/6xmpubk77PZ106lmVpaGiIz9hHN998sw4ePKgDBw44v1pbW/WFL3xBBw4c0OzZs/msS2BoaEiHDx/WjBkzgvG/57Is7a1Qzz77rFVVVWVt27bN+o//+A9r5cqV1vnnn2+98cYbpm+tYp05c8bq7e21ent7LUnWhg0brN7eXuvo0aOWZVnWww8/bNXV1VnPPfecdfDgQevzn/+8NWPGDGtwcNDwnVeWv/iLv7Dq6uqsF1980err63N+vfPOO841fNYTt3r1ausnP/mJdeTIEevVV1+17r//fiscDlu7du2yLIvPuJTcu4Qsi8/aD3/1V39lvfjii9avf/1r62c/+5n1Z3/2Z9bUqVOdf/NMf8YElhw2bdpkzZo1y6qurrY+8IEPONtCUZwXXnjBkpT268tf/rJlWWNb5x544AGroaHBqqmpsW688Ubr4MGDZm+6Anl9xpKs73znO841fNYTt3TpUufvh0svvdS6+eabnbBiWXzGpZQaWPisJ27JkiXWjBkzrKqqKquxsdH68z//c+vQoUPO901/xiHLsqzy1HIAAACKwxoWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeAQWAAAQeP8PMREAIadxTgUAAAAASUVORK5CYII=\n" }, "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": "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" }, "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 }