{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "
\n", "\n", "

Interoperability between scientific computing codes with Python

\n", "

James Kermode

\n", "Warwick Centre for Predictive Modelling / School of Engineering
\n", "University of Warwick\n", "
\n", "
\n", "**Centre for Scientific Computing seminar, University of Warwick - 16 Nov 2015**\n", "\n", "\n", "\n", "
\n", "\n", "
\n", "
\n", "
\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Introduction\n", "\n", "- Interfacing codes allows existing tools to be combined\n", "- Produce something that is more than the sum of the constituent parts \n", "- This is a general feature of modern scientific computing, with many well-documented packages and libraries available\n", "- Python has emerged as the de facto standard “glue” language\n", "- Codes that have a Python interface can be combined in complex ways\n", "\n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Motivation\n", "\n", "- My examples are from atomistic materials modelling and electronic structure, but approach is general\n", "- All the activities I'm interested in require *robust*, *automated* coupling of two or more codes\n", "- For example, my current projects include:\n", " - developing and applying multiscale methods\n", " - generating interatomic potentials\n", " - uncertainty quantification\n", "\n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Case study - concurrent coupling of density functional theory and interatomic potentials\n", "\n", "\n", "\n", "*J.R. Kermode, L. Ben-Bashat, F. Atrash, J.J. Cilliers, D. Sherman and A. De Vita, [Macroscopic scattering of cracks initiated at single impurity atoms.](http://www.nature.com/ncomms/2013/130912/ncomms3441/full/ncomms3441.html) Nat. Commun. 4, 2441 (2013).*\n", "\n", " - File-based communication good enough for many projects that combine codes\n", " - Here, efficient simlulations require a programmatic interface between QM and MM codes\n", " - Keeping previous solution in memory, wavefunction extrapolation, etc.\n", " - In general, deep interfaces approach bring much wider benefits" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Benefits of Scripting Interfaces\n", "\n", "**Primary:**\n", " - Automated preparation of input files\n", " - Analysis and post-processing\n", " - Batch processing\n", " \n", "**Secondary:**\n", " - Expand access to advanced feautures to less experienced programmers\n", " - Simplify too-level programs\n", " - Unit and regression testing framework\n", " \n", "**Longer term benefits:**\n", " - Encourages good software engineering in main code - modularity, well defined APIs\n", " - Speed up development of new algorithms by using an appropriate mixture of high- and low-level languages\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Python scripting for interoperability\n", "\n", "- Python has emerged as *de facto* \"glue\" language for atomistic simulations\n", " - [numpy](http://www.numpy.org/)/[scipy](http://scipy.org/) ecosystem\n", " - [Matplotlib](http://matplotlib.org/) plotting/interactive graphics\n", " - [Jupyter](https://jupyter.org/)/IPython notebooks encourage reproducible research\n", " - [Anaconda](https://jupyter.org/) distribution and package management system\n", " \n", "![scipy](scipy-stack.png)\n", "\n", "http://www.scipy.org" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Atomic Simulation Environment (ASE)\n", "\n", "- Within atomistic modelling, standard for scripting interfaces is ASE\n", "- Wide range of calculators, flexible (but not too flexible) data model\n", "- Can use many codes as drop-in replacements\n", "- Collection of parsers aids validation and verification (cf. DFT $\\Delta$-code project)\n", "- Coupling $N$ codes requires maintaining $N$ parsers/interfaces to be maintained, not $N^2$ input-output converters\n", "- High-level functionality can be coded generically, or imported from other packages (e.g. `spglib`, `phonopy`) using minimal ASE-compatible API\n", "\n", "https://wiki.fysik.dtu.dk/ase/" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## QUIP - **QU**antum mechanics and **I**nteratomic **P**otentials\n", "\n", "- General-purpose library to manipulate atomic configurations, grew up in parallel with ASE \n", "- Includes interatomic potentials, tight binding, external codes such as Castep\n", "- Developed at Cambridge, NRL and King's College London over $\\sim10$ years\n", "- Started off as a pure Fortran 95 code, but recently been doing more and more in Python\n", "- `QUIP` has a full, deep Python interface, `quippy`, auto-generated from Fortran code using `f90wrap`, giving access to all public subroutines, derived types (classes) and data\n", "\n", "![quippy](http://libatoms.github.io/QUIP/_static/hybrid.png)\n", "\n", "http://libatoms.github.io/QUIP" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Interactive demonstations\n", "\n", "- Using live IPython notebook linked to local Python kernel - could also be remote parallel cluster\n", "- Static view of notebook can be rendered as `.html` or `.pdf`, or as executable Python script" ] }, { "cell_type": "code", "execution_count": 57, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Populating the interactive namespace from numpy and matplotlib\n" ] } ], "source": [ "%pylab inline\n", "import numpy as np\n", "from chemview import enable_notebook\n", "from matscipy.visualise import view\n", "enable_notebook()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Example: vacancy formation energy in silicon\n", "\n", "- Typical example workflow, albeit simplified\n", "- Couples several tools\n", " - Generate structure with ASE lattice tools\n", " - Stillinger-Weber potential implementation from QUIP\n", " - Relaxations with generic LBFGS optimiser from ASE" ] }, { "cell_type": "code", "execution_count": 58, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "LBFGSLineSearch: 0 13:20:46 -927.656224 0.0402\n", "LBFGSLineSearch: 1 13:20:47 -927.656744 0.0233\n", "LBFGSLineSearch: 2 13:20:47 -927.657016 0.0109\n", "LBFGSLineSearch: 3 13:20:47 -927.657090 0.0094\n", "LBFGSLineSearch: 4 13:20:47 -927.657155 0.0055\n", "LBFGSLineSearch: 5 13:20:48 -927.657201 0.0031\n", "LBFGSLineSearch: 6 13:20:48 -927.657214 0.0017\n", "LBFGSLineSearch: 7 13:20:48 -927.657219 0.0018\n", "LBFGSLineSearch: 8 13:20:48 -927.657223 0.0012\n", "LBFGSLineSearch: 9 13:20:49 -927.657225 0.0011\n", "LBFGSLineSearch: 10 13:20:49 -927.657226 0.0004\n", "SW vacancy formation energy 4.33384049255 eV\n" ] } ], "source": [ "from ase.lattice import bulk\n", "from ase.optimize import LBFGSLineSearch\n", "from quippy.potential import Potential\n", "\n", "si = bulk('Si', a=5.44, cubic=True)\n", "sw_pot = Potential('IP SW') # call into Fortran code\n", "si.set_calculator(sw_pot)\n", "e_bulk_per_atom = si.get_potential_energy()/len(si)\n", "\n", "vac = si * (3, 3, 3)\n", "del vac[len(vac)/2]\n", "vac.set_calculator(sw_pot)\n", "p0 = vac.get_positions()\n", "opt = LBFGSLineSearch(vac)\n", "opt.run(fmax=1e-3)\n", "p1 = vac.get_positions()\n", "u = p1 - p0\n", "\n", "e_vac = vac.get_potential_energy()\n", "print 'SW vacancy formation energy', e_vac - e_bulk_per_atom*len(vac), 'eV'" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Inline visualisation" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "- Chemview IPython plugin uses WebGL for fast rendering of molecules\n", "- Didn't have an ASE interface, but was quick and easy to add one (~15 lines of code)" ] }, { "cell_type": "code", "execution_count": 59, "metadata": { "collapsed": false }, "outputs": [], "source": [ "view(vac, np.sqrt(u**2).sum(axis=1), bonds=False)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## DFT example - persistent connection, checkpointing\n", "\n", "- Doing MD or minimisation within a DFT code typically much faster than repeated single-point calls\n", "- `SocketCalculator` from `matscipy` package keeps code running, feeding it new configurations via POSIX sockets (local or remote)\n", "- Use high-level algorithms to move the atoms in complex ways, but still take advantage of internal optimisations like wavefunction extrapolation\n", "- Checkpointing after each force or energy call allows seamless continuation of complex workflows\n", "- Same Python script can be submitted as jobscript on cluster and run on laptop" ] }, { "cell_type": "code", "execution_count": 60, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "# {__init__}: No nsw key in vasp_args, setting nsw=1000000\n", "# {server_activate}: AtomsServer running on hudson.local 127.0.0.1:51013 with njobs=1\n", "FIRE: 0 13:23:03 -35.191978 0.2026\n", "FIRE: 1 13:23:03 -35.193557 0.1925\n", "FIRE: 2 13:23:03 -35.196437 0.1734\n", "FIRE: 3 13:23:04 -35.200040 0.1458\n", "FIRE: 4 13:23:04 -35.203660 0.1102\n", "FIRE: 5 13:23:04 -35.206586 0.0694\n", "FIRE: 6 13:23:04 -35.208236 0.0252\n", "FIRE: 7 13:23:04 -35.208278 0.0204\n", "FIRE: 8 13:23:04 -35.208282 0.0209\n", "FIRE: 9 13:23:04 -35.208290 0.0203\n", "FIRE: 10 13:23:05 -35.208303 0.0193\n", "FIRE: 11 13:23:05 -35.208319 0.0181\n", "FIRE: 12 13:23:05 -35.208337 0.0167\n", "FIRE: 13 13:23:05 -35.208357 0.0161\n", "FIRE: 14 13:23:05 -35.208377 0.0136\n", "FIRE: 15 13:23:05 -35.208400 0.0126\n", "FIRE: 16 13:23:05 -35.208422 0.0099\n", "FIRE: 17 13:23:05 -35.208443 0.0073\n", "FIRE: 18 13:23:06 -35.208459 0.0037\n", "FIRE: 19 13:23:06 -35.208465 0.0001\n", "VASP vacancy formation energy 2.09038054338 eV\n" ] } ], "source": [ "import distutils.spawn as spawn\n", "from matscipy.socketcalc import SocketCalculator, VaspClient\n", "from matscipy.checkpoint import CheckpointCalculator\n", "from ase.lattice import bulk\n", "from ase.optimize import FIRE\n", "\n", "mpirun = spawn.find_executable('mpirun')\n", "vasp = spawn.find_executable('vasp')\n", "\n", "vasp_client = VaspClient(client_id=0, npj=2, ppn=12,\n", " exe=vasp, mpirun=mpirun, parmode='mpi', \n", " lwave=False, lcharg=False, ibrion=13,\n", " xc='PBE', kpts=[2,2,2])\n", "vasp = SocketCalculator(vasp_client)\n", "chk_vasp = CheckpointCalculator(vasp, 'vasp_checkpoint.db')\n", "\n", "si = bulk('Si', a=5.44, cubic=True)\n", "si.set_calculator(chk_vasp)\n", "e_bulk_per_atom = si.get_potential_energy()/len(si)\n", "\n", "vac3 = si.copy()\n", "del vac3[0]\n", "vac3.set_calculator(chk_vasp)\n", "\n", "opt = FIRE(vac3)\n", "opt.run(fmax=1e-3)\n", "e_vac3 = vac3.get_potential_energy()\n", "print 'VASP vacancy formation energy', e_vac3 - e_bulk_per_atom*len(vac3), 'eV'" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## File-based interfaces vs. Native interfaces\n", "\n", "- ASE provides file-based interfaces to a number of codes - useful for high throughput\n", "- However, file-based interfaces to electronic structure codes can be slow and/or incomplete and parsers are hard to keep up to date and robust\n", "- Standardised output (e.g. chemical markup language, XML) and robust parsers are part of solution\n", "- [NoMaD Centre of Excellence](http://nomad-coe.eu/) will produce parsers for top ~40 atomistic codes\n", "- Alternative: native interfaces provide a much deeper wrapping, exposing full public API of code to script writers\n", " - e.g. GPAW, LAMMPSlib, new Castep Python interface\n", "\n", "![nomad](http://nomad-coe.eu/uploads/images/nomad_logo_v2.png)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Fortran/Python interfacing\n", "\n", "- Despite recent increases in use of high-level languages, there's still lots of high-quality Fortran code around\n", "- Adding \"deep\" Python interfaces to existing codes is another approach\n", "- Future proofing: anything accessible from Python is available from other high-level languages too (e.g. Julia)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "\n", "\n", "\n", "\n", "# `f90wrap` adds derived type support to `f2py`\n", "\n", "- There are many automatic interface generators for C++ codes (e.g. SWIG or Boost.Python), but not many support modern Fortran\n", " - [f2py](https://sysbio.ioc.ee/projects/f2py2e/) scans Fortran 77/90/95 codes, generates Python interfaces\n", " - Great for individual routines or simple codes: portable, compiler independent, good array support\n", " - No support for modern Fortran features: derived types, overloaded interfaces\n", "- Number of follow up projects, none had all features we needed\n", "- `f90wrap` addresses this by generating an additional layer of wrappers\n", " - Supports derived types, module data, efficient array access, Python 2.6+ and 3.x\n", "\n", "https://github.com/jameskermode/f90wrap" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Example: wrapping the `bader` code\n", "\n", "- Widely used code for post-processing charge densities to construct Bader volumes\n", "- Python interface would allow code to be used as part of workflows without needing file format converters, I/O, etc.\n", "- Downloaded [source](http://theory.cm.utexas.edu/henkelman/code/bader/), used `f90wrap` to *automatically* generate a deep Python interface with very little manual work\n", "\n", "Generation and compilation of wrappers:\n", "\n", " f90wrap -v -k kind_map -I init.py -m bader \\\n", " kind_mod.f90 matrix_mod.f90 \\\n", " ions_mod.f90 options_mod.f90 charge_mod.f90 \\\n", " chgcar_mod.f90 cube_mod.f90 io_mod.f90 \\\n", " bader_mod.f90 voronoi_mod.f90 multipole_mod.f90\n", "\n", " f2py-f90wrap -c -m _bader f90wrap_*.f90 -L. -lbader" ] }, { "cell_type": "code", "execution_count": 61, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", " ___ ___ ___ _ _ _ \n", " | | |_ | | | | \n", " | | | | | . | | | | \n", " |__ | _|___|_____| 0.12.0.13366\n", " |___|_| \n", "\n", "User: jameskermode@hudson.local\n", "Date: Mon Nov 16 13:29:50 2015\n", "Arch: x86_64\n", "Pid: 46591\n", "Python: 2.7.10\n", "gpaw: //anaconda/lib/python2.7/site-packages/gpaw\n", "_gpaw: //anaconda/lib/python2.7/site-packages/_gpaw.so\n", "ase: //anaconda/lib/python2.7/site-packages/ase (version 3.10.0)\n", "numpy: //anaconda/lib/python2.7/site-packages/numpy (version 1.10.1)\n", "scipy: //anaconda/lib/python2.7/site-packages/scipy (version 0.16.0)\n", "units: Angstrom and eV\n", "cores: 1\n", "\n", "Memory estimate\n", "---------------\n", "Process memory now: 480.42 MiB\n", "Calculator 67.31 MiB\n", " Density 27.12 MiB\n", " Arrays 4.35 MiB\n", " Localized functions 21.09 MiB\n", " Mixer 1.67 MiB\n", " Hamiltonian 7.66 MiB\n", " Arrays 2.85 MiB\n", " XC 0.00 MiB\n", " Poisson 3.24 MiB\n", " vbar 1.57 MiB\n", " Wavefunctions 32.53 MiB\n", " Arrays psit_nG 9.38 MiB\n", " Eigensolver 11.17 MiB\n", " Projections 0.04 MiB\n", " Projectors 2.55 MiB\n", " Overlap op 9.39 MiB\n", " .------------. \n", " /| | \n", " / | Si | \n", " * | | \n", " | | Si | \n", " |Si| | \n", " | | Si | \n", " | .SiSi--------. \n", " | / / \n", " |/ Si / \n", " *------------* \n", "\n", "Unit Cell:\n", " Periodic X Y Z Points Spacing\n", " --------------------------------------------------------------------\n", " 1. axis: yes 5.440000 0.000000 0.000000 28 0.1943\n", " 2. axis: yes 0.000000 5.440000 0.000000 28 0.1943\n", " 3. axis: yes 0.000000 0.000000 5.440000 28 0.1943\n", "\n", "Si-setup:\n", " name : Silicon\n", " id : ee77bee481871cc2cb65ac61239ccafa\n", " Z : 14\n", " valence: 4\n", " core : 10\n", " charge : 0.0\n", " file : /usr/local/share/gpaw-setups/Si.PBE.gz\n", " cutoffs: 1.06(comp), 1.86(filt), 2.06(core), lmax=2\n", " valence states:\n", " energy radius\n", " 3s(2.00) -10.812 1.058\n", " 3p(2.00) -4.081 1.058\n", " *s 16.399 1.058\n", " *p 23.131 1.058\n", " *d 0.000 1.058\n", "\n", "Using partial waves for Si as LCAO basis\n", "\n", "Using the PBE Exchange-Correlation Functional.\n", "Spin-Paired Calculation\n", "Total Charge: 0.000000\n", "Fermi Temperature: 0.100000\n", "Wave functions: Uniform real-space grid\n", "Kinetic energy operator: 6*3+1=19 point O(h^6) finite-difference Laplacian\n", "Eigensolver: Davidson(niter=1, smin=None, normalize=True)\n", "XC and Coulomb potentials evaluated on a 56*56*56 grid\n", "Interpolation: tri-quintic (5. degree polynomial)\n", "Poisson solver: Jacobi solver with 4 multi-grid levels\n", " Coarsest grid: 7 x 7 x 7 points\n", " Stencil: 6*3+1=19 point O(h^6) finite-difference Laplacian\n", " Tolerance: 2.000000e-10\n", " Max iterations: 1000\n", "Reference Energy: -55204.429291\n", "\n", "Total number of cores used: 1\n", "MatrixOperator buffer_size: default value or \n", " see value of nblock in input file\n", "Diagonalizer layout: Serial LAPACK\n", "Orthonormalizer layout: Serial LAPACK\n", "\n", "Symmetries present (total): 24\n", "\n", " ( 1 0 0) ( 1 0 0) ( 1 0 0) ( 1 0 0) ( 0 1 0) ( 0 1 0)\n", " ( 0 1 0) ( 0 0 1) ( 0 0 -1) ( 0 -1 0) ( 1 0 0) ( 0 0 1)\n", " ( 0 0 1) ( 0 1 0) ( 0 -1 0) ( 0 0 -1) ( 0 0 1) ( 1 0 0)\n", "\n", " ( 0 1 0) ( 0 1 0) ( 0 0 1) ( 0 0 1) ( 0 0 1) ( 0 0 1)\n", " ( 0 0 -1) (-1 0 0) ( 1 0 0) ( 0 1 0) ( 0 -1 0) (-1 0 0)\n", " (-1 0 0) ( 0 0 -1) ( 0 1 0) ( 1 0 0) (-1 0 0) ( 0 -1 0)\n", "\n", " ( 0 0 -1) ( 0 0 -1) ( 0 0 -1) ( 0 0 -1) ( 0 -1 0) ( 0 -1 0)\n", " ( 1 0 0) ( 0 1 0) ( 0 -1 0) (-1 0 0) ( 1 0 0) ( 0 0 1)\n", " ( 0 -1 0) (-1 0 0) ( 1 0 0) ( 0 1 0) ( 0 0 -1) (-1 0 0)\n", "\n", " ( 0 -1 0) ( 0 -1 0) (-1 0 0) (-1 0 0) (-1 0 0) (-1 0 0)\n", " ( 0 0 -1) (-1 0 0) ( 0 1 0) ( 0 0 1) ( 0 0 -1) ( 0 -1 0)\n", " ( 1 0 0) ( 0 0 1) ( 0 0 -1) ( 0 -1 0) ( 0 1 0) ( 0 0 1)\n", "\n", "8 k-points: 2 x 2 x 2 Monkhorst-Pack grid\n", "1 k-point in the Irreducible Part of the Brillouin Zone\n", "\n", " k-points in crystal coordinates weights\n", " 0: 0.25000000 0.25000000 0.25000000 8/8\n", "\n", "Mixer Type: Mixer\n", "Linear Mixing Parameter: 0.05\n", "Mixing with 5 Old Densities\n", "Damping of Long Wave Oscillations: 50\n", "\n", "Convergence Criteria:\n", " Total Energy Change: 0.0005 eV / electron\n", " Integral of Absolute Density Change: 0.0001 electrons\n", " Integral of Absolute Eigenstate Change: 4e-08 eV^2\n", "Number of Atoms: 7\n", "Number of Atomic Orbitals: 28\n", "Number of Bands in Calculation: 28\n", "Bands to Converge: Occupied States Only\n", "Number of Valence Electrons: 28\n", "Memory usage: 480.42 MiB\n", "==========================================\n", "Timing: incl. excl.\n", "==========================================\n", "Read: 0.003 0.001 0.0% |\n", " Band energies: 0.000 0.000 0.0% |\n", " Density: 0.001 0.001 0.0% |\n", " Hamiltonian: 0.001 0.001 0.0% |\n", " Projections: 0.000 0.000 0.0% |\n", " Redistribute: 0.000 0.000 0.0% |\n", "Set symmetry: 0.348 0.348 0.0% |\n", "Other: 4510.593 4510.593 100.0% |---------------------------------------|\n", "==========================================\n", "Total: 4510.945 100.0%\n", "==========================================\n", "date: Mon Nov 16 13:29:50 2015\n" ] } ], "source": [ "from gpaw import restart\n", "si, gpaw = restart('si-vac.gpw')\n", "rho = gpaw.get_pseudo_density()" ] }, { "cell_type": "code", "execution_count": 62, "metadata": { "collapsed": false, "scrolled": false, "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 62, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPYAAAD7CAYAAABZjGkWAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsvV3MLU121/df1b33c855x2MGG5P4hRkHG4xtMLaZsSKD\nXhsh8IgkxIneWEmUoEQiQrlCRMkFN0hc5DpRckEIClYGUBLzIowE5MOJZF4cJOIYEwhjEA4ZDx4b\nG2w8X+ecZ+/dtXJRtapXr171sZ+Pc4bJWY9KVd1Pd3V9/Wqt+ujexMx4I2/kjXx5SXjdCXgjb+SN\nPLy8AfuNvJEvQ3kD9ht5I1+G8gbsN/JGvgzlDdhv5I18GcobsN/IG/kylPm+ERDRm/WyN/JGXqMw\nM9lz9wY7yZ8yx5QdjO+d867RLlSOgwl7vnZ/AsAfcNLhhVn5tXBU56yL6v/RnNNuMeE/C+Bf69xf\nC9vruXJ83374RwD8DlxXR7U6a4UDgMn4Em7VfwDwJwH8/kYa0QkD7bpvuV69if+nAfzbaLcjNHyR\nfxOevDHF38gb+TKURwB7ZxVUzt0lLturjv5PX+Nd17unFibU762lsRXfNemslU3vmQ8l1+Z/NN86\n7tb/W+e8dI6kxbv3rvfc5RkPJw8Etq3kUcis74W949Z9XpoIwEcH7x09Z9Ni3TVpE/ctg/G0yqYG\nXO2ZI6Kv/7rBuHtpraVnJB5gH6++7jsa99ylnlv14JVlr94IwLdW4mqlbVweAey7JKRVKK1G2ypo\nG/9HO/e20tF7Vqvh1e7z8uSBfS3QNg0j6azVlfe/XzcQn713tB4J/vyIl5fWM77Dub5WVr32UDvX\nK8se6N/qnKvdd7080hh7NHG1yh5pfK3CaMVl4xutzJG4vYbnpbEWV61hjwI9ArIXRy19I/fUgGul\nrxfHXerw2mfA8b1ztXz10t+6t1Wvtf9dJw80K1578EjiRirCO2/vvbbRYMDnHNa+l34tci01wvdx\nkobWM4BtWr00yvmRWfJeo/TOtxrnQ5RDLV4vfi8vPYBt3ffElner/mt1Y+vCq8sxeSCwtdR6tZH7\nRhqMPV8LjzS8mm8BGalgC4sHoldx9wG81Wi8tNr/tTosmGt7YHr5qoXvC3PvGb10Xgt1rQ20AG3V\nv81L7/z18kiTZ/Z/3rWtOK6BvFZhvXhG/da5Vvpr6Wtd00u3lwcvbaPl58HgnW/lvZXuWnjk/l45\n2XNeXLXn1PI06rfSbq+rpa8Vz/3lkU3x3j29zIxUlHdtK46Rymtp6ZZGtP/XcZDzP++aGhyte1pD\nBi+dNU1g77X/89Jjr6k1bnvNSB5bQFwD+Eg6W/Xfs9q88m7Vv77O3m//13pWXR7BFG+Jl+BWpXqu\nNrl0zbVwwtZvVa6320hEjKBowjpN9theo//XGpv38nGt6IZbi/eh66q2a611HDrHI8+3+bL+yPyK\niNSpllo91RioyTXQr/KKNbbNZKvgWxVc24LYirO3pdDLh7eNz9s+aIGXuHVle+ds3FYY+9GSPqc7\nDzbnvE5HS6/xthr9SF31tv726tC6a86PdDC1fIrYLaTknLMduj6v680e2/sD6nXV+l9dhsAmok8B\n+BzShuYzM39n547O+Wug9sI9qEe1OirHFmQbloK2QOsKsM/yKqZ2Xp5hG4PkzVoD9txoo9Hp8MzN\nGgC9zre3H3wU5muvadW/17HbvOmyALb1ChMGtjBbqKV91Cw4fb8+rtVT6397GdXYDOB7mPmX/H/3\nQLbH10JdO9dqOK3eugZ2rde24RbULZBtWVwDtae5dXpCxRfRxxZkew7YloU9N9r5jobvAu7IdbW2\nYPNlw8C+br1zLc0t4nXCQB3oGkv2urZcY4q3nti41GsUI66mhUcaSy/OXjqseeyNc2pQt8bE2nkV\npLWyPmehBurmnvWtFWFFznkTQ7XjkY64B/J94B6FutbR9/JYA1j/v2aG67ZQ09Ze51vrkK3YjsOX\nazT2/0pEC4A/zsx/YvvvVoO5BmxbEfeBugV4z+lse0VhK99q8GvFA1rOt86NwA2sefLOsTn2xP5v\nxJIaBdw7Juf8iJsqcbTq2MtjDWxdvxZmq4GBtrbW52ydtOpiTEbB/m3M/HNE9KsA/DAR/V1m/qv+\npS2gxe9BXWswI1DXeukezNZMq0lLU7fG2L04vXNWc0s6vRn3lkZomXq1jqwHeqsD7oH8EJp66vy/\nV/89ab0n3TPDR4HW4dYw6noZApuZfy77/5iI/jyA7wSgwP5BdfW3APjNOWxh1uFewd9HG4wA3DpX\nLQnl16Ae1QpevPrYQm01d2vmtdVYrtHMret61pRXN9do7hbwHtT2XKuOe3mVTtSD2wLnaW5gHGgd\n1vVTg/tvZ9eWLthE9AzAxMyfJ6K3APxuAH90e5X9ikNNS+vwKNDXQn2NxvZ69hGRCtRQC3wW7sUc\n2zKwZdY6tvfb5bTohL1GdV+pWUojwBL2EF5z3Lt2RGP3xKtfb2mzptV1PEAbaEmfDQPbTkPOfSvW\nN8MA4H9wczCisX81gD9PRHL9n2Hm/2V7Sa1Rej3kNUDfpdF4mrsWrrmeaJC9Hl07ahwv6tgLe8eE\nfedhO5Ra2Jsg03myUrt2tINtad1rIfY+kzQK9l3r13baFuZgjiUsZRdNnCIj1pZn5j/gchcz/78A\nvq19VU+7WL8F9H0bjAe3B/tDV7zunUfBtuBqgC3M9lg77pzTpqXk7y6TfLZz7tVRq75GtXTtu2cj\nYEtaH6p+PcBJHZPxtbRgtsd6nG4n5R5+uash5IR7QLe0aw/wnhbw7n0ssO2xFH4PbAuzB3dLe3vP\nIXNtTWt7Y/qa1Drta8G+i9bu+S2wdVptG+zJyKQZOX4NahFvotMeB+PDCfc75QcCOzjnalq7BrXn\nP0TDqUFeq/RrKh7w4R4BW/stuL3jGsyjWlvnw+ZLxANa5KHqpmeG3xVsr25hzrXEs8Rs/QpsXifa\nk5621ufs88bkkfaKj2rth4TamxWtgV7T0tdUfM2/BmzP9L5mrN1zNa2t0+yFdRl4Ya98rwV7ZGw9\noQ14D2w0/Ja0oK5NkvZWRKxYuGuz6aMrGlt5TWAHJzwKtdeAemuaXsXr9NwXbHtuFGwPbpsuT1u3\nALf3s7lXzrXyIuLB0Kura2BuAe5BfRewvfT3pAa1HNespRGNKvUhsArEXljDrO/pyyseY4vfgnsk\n3GpINci9XUlohEVGJpksKNeA7QE9oq09uHUaapqkNoyw+Ripv5rmHgF4BGoL+DVg23xcKxZsHQ6o\nd669OD14W3Dba8fy80hj7F6jGGkcd9HWgxsXyMKs02aFtpfthIH9L6w4cWTHNc2sIZSK1L8U0jvW\nJryelNHXSzinu+r3OueRervvWnUDcup16rroG3XT7LPtJFk+xxr0Wty187U5jRbQuh4l/Momz67V\n2C2IeyBrjdAae1FuAAYeCzWZ9DfYruc9FzTrcwIxAWw7PnNducYCOGELrz2WaxcT1laDtSAeQmOP\ndLitTncQcqrcJ2CTer4F2NZrTTxGNvUoHXeGvFS3bqv23CYRFZFO2EtIbfwt/+vn6xVqbA/sEYh7\nvXNl4wKp+zXMZPxacr3sSHjXGBy4WcAOAOeK2oGvAKWYj2vwtmD2znlAvw6Nfe06tIHXhqU+pX6l\nbonq9afF1l9Nie6Kg7K2VpDrTpukPvWDWomwY+XWTHkLdF9e03JXD2qruUehlkahniHh0gCw+rtw\nw7dhwJSv9O7IUOdzjFT5rM03M+4t2lrg03DXQLcwa7+msVn5aPitghjV2CNWVQtyqTepR9NRk/jK\nlbp0oGrVXXORQ9XrBmqkeqOY64/zv63Wvka0uS3HXvhLzhTXrqa1W0DXGoo9p3px2IrXPbuC2uvx\na8BXe/zcq5eJVFrDETltvF5HWZuTHrt5mlqHBXQLeEBdY9twbWJI/F5BtOpvtI56mlvXmaOdg4Ha\nNi9TLZuwC7A9p+sy31TAVcJB1Z8B/87Sm0R7bRq7R0arUbRg9uCuja8tyLoBOA3Baxi1MFCHW2Au\nDBEQGQi5kcRCvYI7H5cGMqK1PaBbGrtlivc0tq3DkfrrrUU3wiTtqdI5B1WPuzBt66kWbvVvrsv1\nWjpmiU86ad05j0gP/t6EWlu+RE3xEVebRQ3r8zZQq+NNMq6E3OvxYc7FHG+xejPcQGoA+gb2WpcH\ntNXgwYQXc/6hNHatY27V433Woqf8GPVs3SHr40DqsbRNkk6yTb4GFqiADNVJex0elCYnU48t8Hrg\nM9raeiSOR92g0lJ7rYbRg7zWUHQDMWnZgEp+Emy4BTjQHpcJ1Po8aO1sWf1jo0XkYguxB7mFWY7l\nAa3xdQ9unbBa7zaisUfArtUf7cten9OP0pB7YHt114NaF1mgNbwTgbpmGmjxgLTaQXfYqIRfmykO\nXA/3iAavbT5R9+00MG0fUetXPLA90EVqWluzFMmf2N5oA+wb1aacCMCCdVJGT7wF9T9dDrZlWh/Y\nt259rqapJdwrQAUwNUxxquw76NVBo/o31+uiBPb1Vyt/KQLtS7F79Wel1IdE1oJRl73mSE+kvTZT\n3NPY1vd6/Z6qbEEcsJprKu5NlLSv8JaiaSVF/qelZYpvAIcDNvaNSpvuZQ1cPXjDIQFsdqexgFN7\niJ0Rr5kd2JZp8XW4ZupI3QikBuLNGnS+R4+je9H3moXXJ9lw1fSGX386HMw5QurAmRXktpwY616G\nWmcq8HuaxKuDtrwCsG3ievZupaHsNLPR0hJ/uYX2t3vR1a4ZAVukqbWxh7rm79jLYEujYHVOoC9r\n5BHrWnhPHaHjew3K+o2OWcOrXTmnr6HVwUTnwd2qy5a21nIN1N4xqWIqBlA+uYFbdczE2056B7QH\nt/QiTkfRkVcMdqdBVLtlfTxho6m93l7M8Z7l3oK7lbyWeGC3jr2wNvWKH7KfGw/rROYIeBOBcqic\ng/m/lZrqq2kVlaay1qzCQYNsoPaaRcuq6tVdS1pA1zrl4JzTz5EJU1LlyDlvBW6oYtYP13B75WoB\n78srmjzT51ogj0CuarO2pOXBO7KN3DaOWrJq2Qb2jaMHuXsNqYYW0qwsczb5YspjzOZdVBMqBWyd\nEB32IK6Faxn16hMmbKAVoMtyVFjDej36vlZWC2zPqmqNVqQu9CrjSAfPujxYwR3UeS8RrQ7Tlm9f\nXvOs+AjEHom5sEoDynHXoqmtstTgtklswW2zXQO7ZZ5XTXbaaguidTwHTnAHBT4B27XUnm/DcuzV\npw03fL0LTAMsrnS++n/oV38LcB32sqWzNwq1njDTTqY2dJy26ARoALsdat2xgGTGM8O92bq9vKKd\nZ14vfx/IpeGYRuTdUlt18cA2G9dcyC3YNlxrJPqc1gKeZohYoZYwZXBJxnG5NUV1nmvAtkCunfMy\nWMu0DUt9KLjDiI993XnnPKBbYNe0dQ9q6VS9l+gWJ07LHtNaJBuovYToRifX1DR3Xx4Z7BrkLZA9\nqmq2co6ztq2wBnho+PYRNbhttrTYhtIy8dzNZeQ0JjLP4a0FVx7cktb/W2NsTzoNTNdHKVeqdKa1\n85VzHtD6Oa0pBA9sz7cmuAe1xDOZeEmVT0mL7vWhEhpN2DYwT3P35UtAY7cArxCmtYG337vW87fg\n7pnk2m9lD6hrgtrYzUK+8R3II4CF/I1pWpMMy7UauyO9znNkj8pdwNb10xqR1OpEh2swk4prwpbN\nwqVArG+QU/K/4N2oMmFhfu1g6+OWehuF3DgycVi4a4q+BvMo3FaTt7JYM/esqecB24O+2xGoNGi5\n1gK34rWnWrXfF+ZW2KsHr16shtZhWw+eCS4wa6ih4gnK12HNo0C823paY8CDXF/nzX/48og/fH9X\nqB1T2wvbt7FaULc0tzcG77laNuRcb26kO76GD/QI1NpMhAlXh9pkT6jzlcNa1QJ3g9oDuFUfrSbk\ngazPeXXgge1BLSJxWrDFYtAdLMHALQU24tg57ssDgW1lFGp7fc/p6+U01aGuaW7PzejDLRM3vSTW\nGlYL7hbs18Ld0lZQx1pYJ16VrxavyrzwXaDugTwKtc7fXcpfWz010RCLk85AIC/lIUBjhbuqkW0E\nNbj78hq+Umqvq2nrhiN1fzlNvonWM8NnJ9zTJiOdLNA2A0cbWA/2lsauaaymmd7Q0nLstUd9/lVC\n7aWpBbU3OWYnyux4WpeTB7ULOKGsWOze576Lu2bc9Cgauwd1LeGDUEscrWhq5nfLHJ+V32uUtaSK\nGQbUTeGa1vDOeXCPauyea6UTqqi9cM+Nauwe1F59elDbNKJRDlKGtYlH/Qwrth7lXuvkf5tyoXy/\nBbRWgFZL24y35TVp7Gt6qsazpMB6cNfMcQ2zB3Ztos3rh2yDqGlF28Bq4WvH4Nr3rIMe2K0hdgts\nr5Mbhbo3hq6FazDrsrd512ENtP6iFME3w23d6fRruHfNN7fNXY852jteb4KLPPLkWcuXcI3KDuxe\nH9Hq6T0t7ZniMwB+DvzC+8Cn3wOefxr4wIeBb3gX+DXvAIdne7gD9vVRE88cHzXPRzW216At4Gj4\ntny98rZ5tv41GrsFtvW1VeQJAbh9DvzD94G/9x7w2U8DH8z197XvAPOzVVtLfAK1hLW0hlO1jkiu\nKWVGa3hTiNdoa5hwWx5p8qwmFmrxR6F3bvUurTWIlgYXqP/a9wO/9GPAi88A8bw+86c+AXzgbeCf\n+xjwe3+w3QBtVqy0tIkHsj3fWwKrxa/PoeHbsrV+rYxtWY9C7dVRqwNpleuf/X7gZ38M+Jypv0+q\n+vs9P+g/rzU/0Rpb1ywXj8udtC70eBmTBwLbGzew43v3XDcpsJMRi6YFuTax4/ME9Rc/tX9OPAOf\ny+fpOXB81m6QOn0w50Ymd1qmem2W3JriHtBiarbg7hlbPW1dM7e9DSYtc7sFiZf203Pg534M+OVP\nYSe6/vA8ae5aXXn10UpbDW4d573Fmg1teQSwa1Br/wHEqxSvsGsV4DW6n38/aeqWfOEzwM++D/yG\nj/sNuGU96YYjfsvU86Cuwd0C+zFM8Zo2bVlFNYh7YNTK0Kb7H7wPfHag/n7ufeAjH2/H2zO3r9HU\nMOGujFRQW4bAJqIJwP8J4GeY+V/pJ6oH9d0TvE/cgOs1QNHa//C9rfnmSTwDP/ke8Js+vm+0o+vc\nOsvsnGuZ6a1xuDcr7h17z7+rxvbKt+dqdSJh75k67HWEDOD/Hqy/v/ce8PUf9+tEl5OdHGsNE0br\n/d7ysBr7DwL4JICvqD9MpKexH8gE1zIKdGscGJAmykbk858GDtgCfRe4rdQarAdny/XufyiN7TXu\nnuZudQqeBq2Vk+c+e2X9ec+wZSx1KuFevh8UYkmQ9sekCzYR/RoAvwfAfwrgPxpPCNDX2DZ8hdRA\nGWl4tXH2Bz489uxf+WHgiD3Q4lqWQq3h1jSH51oTa7XlLQ04On5PY7e0VQ1s27n2Oj6vWdiOz1oz\nv2qw/j704QS2jbtlFV3bOT2KxrY9c11GNPZ/BuA/AfDB9gO1kDpX09jefYNizbOe80woD+5veDfN\nfrfMuXAAfuu7q8bWy2Ua7NZYsgeOp7n1cWtirQX2Y2jsHtgjZeAB0EqjByAD+G3vAj/6CWDp1N+3\nvbtt+V4H2LM0RjsqkXsBfj0nTbCJ6F8G8AvM/BNE9D1jD7Yg23OtxF6ZgdHG5lWE1dYBaZ36A2+r\n2VNHPvQ28E3vrBq7tiW11oG0GoCEbSOG8Xtr4PeBelRjj4Dd8r34JFzr0GpWi7hvfwf4qreBX/gU\nqvKht4Hf+M62E7WdhZjeMm8xMqRo1emdxLLgVV5dehr7uwD8XiL6PQCeAPggEX2CmX/f9rI/rsK/\nFcBHc9iDu1VjnoybH1cB7mmX47O0zgmk2VOtucMhNYp/4WPAW8/2u9Vm3uxcW38gktUxl2eW340r\nX0bB+o1G/SWU/G9fY3F75vsasG3Ya5jXdKSbMqaqtubdM9JJ3qQ5n5M0RgCR1h8rzV+b4eNT4Btz\n/f3iZ7aae8r195GPAc+eAWcdl3F2A8tIfm0Z1cpwJ6327VXc/5FdW4jZi9C5kOi7AfzHdlaciHj7\noFo3LL6esaq9atVwFLZAedtDD0ga9eC4o/JrYTxPSyI/+R7wuU8DX/XhZH5/0zsKal7hluMD8uez\nGRQYNGU/qHOBM8AMggpTDoMzyLyWHG/PgQHaQMx7sAG/Xdg21FICNbhrgOtzBYgVav1DlKx+rmeF\nO4dBKUv5U8speyocCRwpw03gxZz7wnPgb74PvP/ngF/I9fexd4Ff/93AlKE+AzihHj51wufKvWcA\nF6SOwfoS5qhOLOafEq5dI72PyLeD2f5S4PXr2JVeoKVt7TOlBVoNrs8H53+dx3qWi5WexpHGeHyW\n1qm/5eN7+BXE2zCDxA+MEBgUImhihBAL1CHEDLGArRzy/1g1cWYTxgoyJ4W//YVPU2StcK0cRzRQ\nK7xxZF5ooqyEKT+OTBhgDmCmdMy0C8cYCshxkXBAjAReCPTBp+Dv+l7gOz+eQaMVuAu2mtma27uO\nqZK/nnT1Zc+sqplfega0LsNgM/NfAfBX6onU4VbOWxnxYPee0Xh0Lc+1xteaNW+9CTYDODBwYFAJ\nI4E9RVCICCEiTNnPYQoRAQyi5AeKJVx8cdkeTUp9PZ+A1nArTd4yt3tQy7GtvpZ56fgaZAat4ay9\nS04EclAGPrkIAnNA5FQakcN6DhnsJSDGAMp+jAxaAjiGYthtTX0xFXLdeSZ3b/xcKxNddkPN1mu0\n18D9gGC3pfcg3Vo8gO05mZaslZQp1Z62bo11alDbyTULdTbF6YACOA4MOiRNHaYlw7xg0n5YECgi\nIDsyfmrCIPAKuYK6AC9afQe2yvw1QNtztbKqlan6nwbYAs60yUk5jtBg55Lgih8nLEsAxQlYAhAn\ncGTQgvyV31CeS2LvMGPze1u2XkfH0iPaWsrSbZdeJdwF7rY8ItgtzT2iqWutkbD/dUPntlqyaiZ4\nazlsBzXSeFoBjaOCWzTztGCaFkzhkv3kAkVMWBAQMVFEwIIJEYGWDdyBpZnrJm+gBkOPwcnL9zVA\ne+XVOudosgJxxY+q25JwpOznXC485ZKYsHD2M9xLnEBxxrIwsExp4iyPs7FQoptok0SuQW2tsxGt\nXdPW4lfL+FqoW3C35RVpbH2d1sTXamp92Oo4HLGNseb0bG5tDq+MrRXU2dEh5rH1gjBdEszTBfOU\n4J7DRZprcrSs4eykiYcCdTbbEdWYezXRU5YU4LbIPRmtsloZNs4zUExrAMXkhtLSUUHMmzBtS4PX\ncMg+LRGIDCwRvDA4AnEBaKE0Q05hm0fGVlsvqL9x1hqm2bza8qgBvYN6909sGRhZt2zLKwDbAtiC\nuqa5K/Fr82rkFltBLahr4+wMNW1m2zPUR9HYSzLFpwUhXDBNlwz3GXO4YC5NVYVpDa9QRwReAZdz\nBWoNuArfCdqHEhnXFtN6H46bHBIibXKcSoYnLLlELsWfEDCD4gwsDCycNTUQFkJcImhJioF1giJA\nkfPyGK0TzCNwt9qMJ24b7FXIXbR2Wx4IbD39LnBqH+ZcdML62JsVt5nS8aIN+Eh5eJUnlWxX3CYA\nU1rOwhxBc0ywzxF0iMX8njPMh3DGQYVTk704LjXfpNFTE1/D+zF4ze3EOXWFrVOLohpRO3VkILa5\nDaY01tKZMOOCBSEuIEorCAicfyKc0k+ZBQIzp0m0GMAxAksALwBN6Zohk7tmdo+0seJsqclxb4N/\nbceRPm7LI2hsC7X4ZK7zukIJ6999rnST5SdmgfITOPpncUY7vp5m3429Wa1Tr2NpmiJoSuFDOK+O\nxJ1wwBkHnEtTlfABF8zq/JTH3El7r2HR5ttxtx2HZ82t8rWGt+fvKmnVNJvZzvnN7PZuhoAKzMsu\nh0HB7JWUhA+YKOIsKwphXUFgADxN4CkgThM4hPyz4QEsG4RqUG8yqfxa2/H4LO2P8zW8hptg27Xq\n2nu58uC2PILGBvbAas0t56O5TkR/fMpzGWhkMzxSvlzB3Sp0D2gP8IaZntam19lvCusMeJguBeoj\nnXCgE450xgEnHJH81DSluZ4L8GvzXU11af7ruQXeyLRMuCECLGPuXLJqUq0GdqkBVVW1PqD8PHe+\nk3MEnH+rSs8MrKlaw9ucbXN3weSWzKGEL0W3i9ZOmjunhwhxYsRpAgWAAyNOAAdCzMflVz3vop29\ntrRrZ6zgzg4a7laPMOK+ZExxHfaAth1BK1PSqhzTu0DuOM+aGTHPLdTizPp0mC7FBE9gZ6jDCUdS\nrgC+h13Oy6hyzs3cjjS3I9K9I0AtiamwM2u+acus7CranN5IMpYow503oDCp+al2Cpdil2gbZfVX\n2+aoSuaSsT+kpcGsrSlwSTgTwJGwyK6/CYgTgECIIaRrWya4Jy3LzyrRTdtipa0V3PcCW5sNbXmF\npjhQ7yLvADYI5cfgo4rD3l4DugZ3bQSwMeHyLrJJLWlNS54kS6b3MWvsGzrhSLfJxy1ucMpAJ6iP\nJXwqoK+jSh1enR6VruPvNazmoxPMAjVvgdfiaXK9UbH8W2CG+KKx13Pr9N/khi/KrPZmGtbSOGPO\nUJ9wLHo+bebhdUUrpMSmTiYkoAOwqJ/u5RCBELIWV3XqAW3bhcdkD26rsVlf1Noq2jLDX6vGBuqm\nOCOVZnSu0/d2xtica0d+c5gp9Y5oaGzb4XmVZ6UC9ro9lDGFmNenL0lbT5cypj4K1Bno5N9mkNO5\nffikxt5ikK6QH0rjFuM17sJpEwvSOrgAHdcwgMoYXGlsKRYLt2hGInBYYY4QjU2IuxRux9I6d56f\nSuSYhysX3CrDfV0CZPXSTBoGREq71DCl5W1MBEwEnkLR4LBQ2xdTrFioW3Bbc1w0t0DNLVg9yD3A\niyZryiucFfdKzvvfwBgbOWr5zeGYx03R3OpVwl2hVnBTgBpjp1nweVrKktYBZxwpaeMbOiWo6SWe\nFMgF8C3wAv1BmeqHzUhTANcrvdvR6oSl7DUvYEdG0HCrPJext+RbVVkBWZVJWmFMUMeIDHeCirM5\nXhs7r+FGtTQiAAAgAElEQVRtrqxLUN/gdjOZuKixNa9NADIckJ1pUwEaIYBDQAgRMXfGm6Wt0TG2\nDo/CzVhN8S7UdzHJ2/IKTXHxtbb2ep6OKV5+TJy245byA/DObRbqmhmuk2CNDmOGpxnwtEVUb0I5\nTHr8nNwNbvGEbvEELwvAErbntNZeId+OxWeFih2hFlOVefUL3Mm3MBctrqqvtGcH8gR1hpvzhhPx\nQQXqSwF5O0ugc7QdR6fzArXtsNIkoSRWDQco5Ocmq4AF6ikgTAEhTAiBsegxttXWHtQ63NLcriLh\nrSvt9z4wf0mZ4tbX9NXuqZkqMra216tulWkPs1dGI85LtsBN6QUO2f89UdpVNtOlzOeuY+jbDbzF\n8S1u9DFe4oZvUxPnEw6c41HhA58SHpxRYYVMDoeYtHOI7Po1oHXbttpaA55mmKnqL5hwIbXBJIfF\nP9MBR5xxoiPOdMYRB5zonM+fcCibddTIXN6Gg2hrvYMtZKhnXHABU9LUkQIiTfnturia7i3Tu2Sy\n0gxrCmPTtrQZnrV1Mcmt2X2XcfZrAxvw1R5jO8b2pDXGlgybOIvmRrtXrZWhBz2wL7+cHSJsXrXU\n68zb1dZsjiu4E9AJ4qf8shyn8As84VvM8YxDvCR/UeF4xhwvCWBeMMXsc0wuRgTmAvYG5ojiNtra\n5lN3aDmvG5srd27yjnlI+0AQZFIqID0kEIiWXFbp/4EYC63zEyEw5iliDgsOYcE5nHEOB8z6xRiK\nCmgGyksj2w0uqeTPmDHnGolYyux5vl991KJqdtuO/xqFunFGUw81wl4D1QlsywOBrX8XxcJswyOi\ntbaFXF0jC6oUtjCOVEprXsIz0RXU8mGEtfFljZ3NR5n0OagZ8BucNpr5Kb/EU35R/Cfix1vMywXz\n5VL8yRyHGDHFZesvESEypgKxAnrBBu6SJTbZ08dKWxPUCAgJ4DReJSDPN5T5TNnEM2W41ZyETDaG\nmTFNjGmOuEwL5nnGMl9wmWZc5nMCm7KWzjCW0ReSub9q6WwFYM5/FyyYEbPGX/LVZb2bGGU6vWRo\nsN14zNl7NkM8a4p7k2M9TdPSPHV5RFPcs2Nr42p7rwVah/M1HHPryxNoxP3KGe0Ma+Nv6Ucgs7Jc\nlp62GvuszPF1nC1uhfoFnsaXeMYv8DS+WME+XzCdl+xSeC7HC8ISlbPHGmyAFgYEaJXXTRdZaScC\ndAnnwPrZpwQ4TwRWcGOKoIlAU0ScgTABcUowxykksA8R82HBcliwHC5YDhMuh1SCk3yMIg95KKzP\nF/N7QdjMLhxyicsAYMGS3pYjLh1E0xS/Rju3AK9BXTXDRxqn11jb8kgauwV2T3NrkDXgGm7JYIba\nvoPs9aSy6V9r65rWbkANWk1x/Q61wK1N8WNx+zH2U36ZoI4v8Iyf41nMcC8vEc4R02lBuDX+KSKc\nFoQLI1wi6BIRLlz8cGHQwuW9ZCjIobS2LitbG3rGYnOsNDaytqbA4EBp33zIcOe99DRH0MwIMyPO\nETwHxDkgzoTpGBGPC5abCfF4wRKn9DomJsQQMHH+4gy4jNw4IG8q1G9/pdI+YFarBhcsuGCiebPf\nbdXWJsN26FVTDiM7PUvbY9UuBeoaxCNLXK8NbN1aWkCPmuQWZAM3i/mdM0pqnD2itVvlpcutCjeX\ncZuY4nqXsx1jH70xdtbaz/hFgnt5jrfiCzzJYIfbiPCSQbcR4eXq6JZB54hwZtCZQZfsl2OlpRcA\nAnOBvV4DejJNr1/bSTXavBXF+eWKtE7MAeBDAjscCHxg8IHAc0Q8EPgQwDcL4pOAeAmIi/pSCoWk\n0SfZ+50l0PpZJLU+vk5Tbtf6L+VVWDNOz53y/gthjTZjlcOotSeRsoqYW5Nm15rjbXlkU9xzrckz\nSbDV0jB+hlpmgqzGtmMer2KWyjW2cnifBcrmeMifMpoQ8xh7ccfYN8WlmfCksV9kE/w5ni0v8FZ8\ngbeW53i6vEyAnhLU9JxBL4w7MXBGuubEoDMAHd6ALXnl9dgUt9fOaXtJOUFAhprTGHpCgVrWiPmY\nYMaRwEcksI8ESPgpgS+UvlUmO8ZCSCb9QnlJi9eZeKL8QdZ1ScvuJT85L9HIhyvK+FpMFKtfrKa+\ndpztwb1pQJ4p7pkB2q8l4pWOsWumeE1be3AztguLEa621plkWxvYw3lNxVR7Xp093mvsnXHoj7Gf\nIC9vsTHF43O8tST39PIiQwvgJQMvAPoigOdcfNwCOAF0C+BkjwFk09tae6KxrxFXuxeNzev3xfQX\nSW4S3LhB+gDkTQK8HF8oKS/Oi1cBaYeYfJmmtI08A57XyxeWsbWU7i1OOOCAY9m3JnVw0Rq7wI26\n0ei1mZFxd0spuOPBmqZu9RY1k7Iuj2SKj5jh3uhO398wxYsd6JHJDtR0nZVjO0Yz6EypUjPjkEm0\nLeBac28m0fgWx5jdcsLN5YSjuPMZeAngeXZfBPCFrc8ZYL4FcAvwKflRwM5thGtt575iP0ARoDR2\nApoyyHQD4MgIN1jPS7u2g/kc3yXOOM8nnHnGDWZcKDueceZZvUwjQG/filOfcMD6NTXA3xCPNtQj\nS876K7EAxhpazQz3ZuVem8ZmE+6BLdoZleMW3LoANOByTDnvtK+shcbnJpqQ51nx8pHB7ftLevPk\nCrjab8XntHR1XjCdI8IpZvMbCcwvOO7zKHDzCYin1ddhPmEzV6M3OvEDgk1qayZpwANAFyCcAToD\n4QSEY7JAwq0BW5e3lC0BlD8xFY4RIS6YOFtBQfaprfvpt0DrF2Ea76h7MFt+upPUvIV6bRxORD3A\ne9faxPXlNYyxe1AD9aWuWsHY32TJMySbXWi0vaQ3b+GNu/IQbR2bZnNcuTQ7HvNrl1vTvGwLjWlN\nejpfyqw33TJwm81qDfPnVfgLAH8xQ3wGllN2Z3V8xgbojc/Zv6fIctd2s8rqT+fsMtTTAZhuAToA\nfFRjf88qIgBPgHDktD7P666+Q0wLWrIZdS7mt31JRL39LZ1vDeremLrFYWRsXvYAY7spxWtYI8ta\n3hjbmhNtecRZcaCuvfX/LeQShwVan7NmuNLeTHu3udQeo13GdtyU4S7auviiJey8rdbYeUtoTBp7\nuixp+UpmvPOYegO1dhn4eE4AX7LT4Uv+IL7eyShg8wOCTbQFnBTo8wGYT9nP34Wj7MJBlbm1KmWC\nbsm75pC+4jrnz0wd4hlHTOatN/1hCv19OEbgqKAW6JRcA/aGRwHacWB7sRcB9o3vS3KMXTPFgTrY\nNc2tr7NAyzlrhhttXWbOsdXWotxrHWlLY+ts5jekBGr7sUH9btN2+cua4kpjv+A0phawLdyfw6q1\nLxnuSwL5fAbOl9UViLVT5+4rsmxUADfHhzm5eECaEMvfXw/y6eYa1NkRr8ObKSyY5gXzcsGBz7hw\nKJ+Zsho7QL0Bhu1W1DLO7o2pe1raAr5733rUpq9B/SU1xrZqwJredkb8IcbYLcDVODuGNKETkVrd\niKa2Y2sP8txQVnN8+8EDrbHtwswh7wGfLgum2wh6mcCm57xOkmkT/HMAC9ifB+IFWJYM9QU4XYDT\nsvrS1lylMlahTZGf49JQBxVeJiDOGeopAT2pY2GApKp1lefjQOnnkMK8YD5cMB8nHGLCdjXFL2XI\ns9XY9ptwvP+FFA9wD+4qc6zgzje7H1MYMQV6jdImrC+P/D62rTWgDzVw3RhbAy6TZ6Ktw9rCpTXe\nZU+A6SRlE1PRBOzBvWpsmR0/4owj55c5ZPKsaOyYNLadMBOos+PPA7wAcQEuC3BeEsy3C3B7SX7E\nCnE0/oNpbKyfDrN+nBTUEzBlF2XmPKdBWkmZEbcfspgXTIcF082CeTljiYQDKE+erd9f2b7tbT4T\npb7DnuDD+vC7QL2Z9LMae0Rr/DM1xvZM8VJt8AGXa+8yxva0tOnR0q+7KbDZ7xtanWRVYyvzjrcf\n95fJM19jpxndOZ7Tyx3nyzrGthrbM8Wz1uYlaewlZrBjAvpFXMHeGB58TZPoi+aQaP+KM+cZ8jAB\nU0hQH2RX2oTtqhOZmwXsKSIcIqbjgvl8wbLkd76BYv3MRmvrb7FvtTVvn1nT2D2luvEZ27F2D+R/\nJsfY1jzwgBafnOtg/t8bY1stLaDnc5zjYmD9eikA8PqllV5v3AG8fEdMNaDtNzn1BNp2rD3zJb1y\nuSwIl7zcdctpfP0cYOXbML8AzgroUwRulXsZ/Wagfa+mWmKbkfd1IR3GgvxGV3IC9xSAOQA0A3RM\n424cVfigzt0wwonT9tolYl6W9HEHkFnmkvG1/t6bwKzG1sUc5zGYm46Nxo4G7hrIvZ6jBfZr0di6\nuWhtXfNF7A40ObZg63CrcEL2gXVWPCa4pTBqII+a5wVuLoAjQ55yb5fAtGme3qWe87vUYYn5ZQ4u\n20gha9J5CYsvQFy27gUn9zK7W057WvJ+leoUgW4Ste5VlVQqRiesodYu2LAqH0Q1y7IA4ZJdXuum\nU/bzbrq07z1BPWWoJ05ae95MlunNKGZcLYmuaeYa3N5wbdc+RFtHxwyvuRrMXkJqUNuu2ZfX8Gkk\nuUYS6G0vJewhX5QfGk5NnJXlLpU+r6KqlQe3rGnTWGQMh82GlTKri/2OtInTr1kI2OHM69bQW4Bv\ngajWqJdLdtn8fsFJMwvYArcsg+9GDua4ZU/ZGvV8CzTsObFoSHV8uYqZ06va0wJMl+zyujflHXSU\nQadzemMtLBGBF0x5X/n2c0l6siyudYGkTUVLk9fbtYBuQi1Om+IL1pc8WspnRJO3oNa1WpdXaIpb\nuDXQFm6rrS3cNW2tdAZHbGfFswkmSr2msWvHjl2rNz1sTXKtre372vL1E/X+9CW9vEEn3kCt16Yv\nebLsIkBrsLHV2jlJVb8GpBx7nYE+7nUMGm4gl1dYq2ResrsAhzPK3ni6TRtacLtCTUtEWAhTJDAv\nYGCzy2yC3j+gdpqxpIf3dddzPY29ALuJs7ITaERb92CuAW7trrp0wSaiJ0g/eC9b+P8CM//h7VUe\n2DVtXQNaa2hWYW/cbWEO6zn5LDEi1q+YYm1pPY3dqgMLtRrHrTnfz5BvTHFc8meMEtx0iaBz2lKK\nrLFke+iSoT7JGnVMztPU2hT3RDcFr1S1sxrfswAkHk/KeX1xrk4m4LAAx0saZuCCAnYoGlteQ40I\nF0KIhCnK6Icrprg1w1eNvVpY3NbWtc7eM8OrUNc087Vwt7R2X7pgM/NLIvodzPyciGYAP0pEv52Z\nf3RbayIWZAt57r67cNechtrT1iG1Hvk/R2y+ddWrxJa21h2ncduft93CPW3gTt8oC3aMfVJg36Ls\nJpPNJ7dLXquOalydTe8CeLp9B5x33Boj17JpdYWnN0rN8/a56bPFSD92ma1WzlCHvP1UXmbBKcOd\nPygRllSNUwZ7u7S1hXufIAVhT2uPDM/EZxU3GPWJs2s0h9Ee7rkHNMWZ+XkOHpGmn3/JXGHCnra+\n1hSvOQ9qDbfMjBsSmdfb76KxN2WbtUFJ1XbCZgu1HWdfMDErU1x9NCHTGZXGlg0ot4tyWM3vAnY+\nPqkStL4teQ2zbBvSGttrUlJL1kTXvo6I87GsPkYAcUHS1HkSbToD8zl1aJiR5hvOjHBm8MKYlpgT\nkVJgvzO+34yyumrv5PEzqrEXVQLVTx/1tMW1mluf68sQ2EQUAPwNAF8P4I8x8ye3V9iH1bQ1mes1\n0Lp/7xmLHtTiUzbBZSeEFL6K+q7lvtPcog109NvJMzsrPvOCwBFTZISFy6eNNho7gy1j7FOG++WS\nxtaimYtT2vpkSsorTW3faJj1uV4T88zzzbFhqmzUytVPGuoTsByAKNtP8xib5LNPC4HLDq+gPqaw\n3ZCymxGXPDO2E569eh7R2DYy9m70tPeo5u6VfltGNXYE8G1E9JUA/mci+h5m/pH1CqtttWb2wp62\n1npE/BbYtRoR3SOFbZqg10tf4yqterumrTV52jBRGiAvCPmzwCEivfCgxpo4A3wG4gXr1tElj685\nK7P1UlyUk+ZiTWsoX4tXsnYwpK+zzammCCUdF5RvMWxmQ+aY3GXJM/1ilksZXJA+8XThYoYH0dic\nXu4omprsUteaOtKb4wVuW4d3qn9rodrI7jpZ5v2/Zj+15apZcWb+LBH9JQAfBfAj639+WF319QB+\nnXc3toCLr6G+dgKtRp1HoEmKDrcq2QPZU1XmAXpzRPnJHSigBer8meCi0lSbsHu8yzZR59EaQrtx\nxDtuudE2rmsJyrediFvEKm/SB+/mnnKZkIAdCYj5m+SmE10f7DR6z6Tw6rbTebcj9nqMnt9qbB7c\nDODvA/h/agkqMjIr/tUALsz8y0T0FMDvAvBHt1f9Tn2HybCnKzyYrVmu9YQGXK5rFdYAfT1gR2G2\nnQSwXwbTu6CYN9/9lkZM2YlFJ5OteleslzstnrlttaV3zh7rpjYyWtFij6tQA7tXSzXUqdOj9Sur\n0hmKA5fPUxGwrpl7usMD2OvMh4FuPWC0W6wpJE9ja98qTq1UVxnR2P88gP82j7MDgD/FzP+bn0kJ\nezBbkXss3BZkC3itYHrqtZGM0V6701do0V+7le2n6Vc69C91rH5q1LyjZ/MKppO08jzUoZ5MuOdr\nxamP9ZSl+F5xbrSoKeaS/pwnCNwa6sjgSLmzY9CSNLVoboq8vmGmY7efPrLNsgX4NfVfjXQETO8h\nNchb7botI8tdfxvAd3Su6sWCfVVLU9FhrbFrgI+YNJ4ucZJkO91rytXrsJ0HiLYOMbW7sPm5HV6/\n/61osdraQh3N40VqcNec/b8GW0Otu10Zw1uopQjs2NwWUakxlbddtSlTPGltBiIhREaMWL+KIq6l\nR3oddA3olmk0FHmvfXr/swm5OkFFHmFLqT7nlThjD7SEvXlcD/BWz9gjUQV7vfio5jZStLW4/PO1\n5dcuC9hicsKdPN2Z4ny9xp6QKtkLe3BrsPUiYlThnka2/9+xIx0Vq05Ma+1cNhpuyXz5wUH5UmyQ\nSUonUTXp1blXyG7HXcthC1IP6tr9tQT25RHBrl1n51elKQR1TU1Te2PsUdVaSfII0FdGLbKup2L9\nnWoN9gZk3sCtv1VmN0y1Hm8nyjTQLcBrYOsZbk8bSxF6s+m1Yi5wK7fR1AvSsCRSMc2pWDJJc6df\nCskRpcLej7GvdVfXcyuilmZu/b+VIG2rteWBwP4hAN8M4CNIe1i06NIWsTDrc8CYpq7BXdPajtTq\no1fJXmddjZsL1HZcLWPrPeTJ2cZvc2elZ4YL3F88AD/8EeBHvhn4/FcCH/ws8Ls/CXzfTwNffd4u\nnV2w7yxsFrU2l//XNHYJG6jtC1IyRLFlAik/IMHN6zzGbuJsVxeOqynUZj17/2hF2DvXMxlswvry\nQGD/BID/C8AHAXwtgO9X//NWQFvn79rF1ki0cd/zkWzuqUS4WeoqTmnwnFxSySaB2WSlPKryPG8t\numaO/4fvAv/0bSB+MJ9A+o7DD/wW4BOfA77qM8APvedvObXr2VLq8smLmpm+LZktH3oSzeOCduWv\nLKDsV39V0JbZtfVca0b7nFRu6gH6oI1uIw8ENnKCfzk/+IRVczPW5uAZabVm4BUaHL92rhWfE613\nqffYOzzagrH7/1hdNeMVkBn7JSyB+/khQ/0hJ7IJWD4E/CKALxyArzivz7Di2Ufeengrz718afih\nwlfF16qnVpO6qj5GYKxFWmtUtWeMXJskdK+4Wj4P4KfNOVua1yXSj0fHV3tGo8utlXOv37gSwl0j\nVArmrkDreEe0tbgf/kjW1A1ZPgj8+Y9stfxcCXuz7OT4Or13zaONozbed+UanXCPuh6P9Jr237MK\nfHkEsCMA2UreysA1vVXPv7I2vKT0QO7dY88bZ3/F9S4N3UoLag/uH/lmFPO7KhPwl795BbgHtLe7\nzTPfgfvlWXeG2jy/Ou5rNXa3OY1q6lb77yW2FY8vD2iKa/msCjP8WQ3PPK/JKNytntFLj/OvEY3d\nqaNmY5PGeA+t7QGtk1bboPL5rxyL/5e/MjUMb0kL2Mxxbda6W6b4naDOZaUmvksaSF3T7M97sF4N\n94gSabXPa+DvJbIujwS2bUE2IT2oWxkYgdxe76lkqt82qrEHxZrNBWrnmrtIS2truD/42TRQ6smH\nPttuGDJZpn0N972BNqLH3Loj9IY5blifG9HY1fh67bIWaasd9+R6qIFHMcUD0tLXSI9mw61rR3q6\nawqT/bIaOTeQ/N34UCXJAt6KpyaeudsaX08Aftcn4e8D1bIA//onfdPbmufePvQa3HeZQAPgWjZX\nWzlev27DXUXcass2wpp2rimZ0ed5cfnyCGB/BdJ6tk7MSKZatTVCWa0DGI13IKprkuzEW9Zbzb13\n1WoW6t44+/t+Gpg+145z+hzw+356C3FrwsybPKtNnF07DrbXb+YobCfZiGdzfG3dsv1H7QHXKJ7a\n/TVpdRy+PCDYAckEfxvA4eGivUqu7cr//yUfOqd16umfYq+5l3T+az4D/Mrz3eJ/CLP7jdTkurb9\nQGPsb8e688xCvet3G/+riXeP53uGX+sZ5B/uBsWVaK5MPgNgyofqXr4iKi26325tXtS7yP7ce2md\n+oc+AvyP3wx89ivTmPpf/STw7/408DXn7ccb9C60nrPbheyU6XDTNNMf+nyJS5dpJ57Ncatuq81M\nj/JrD2hFWGtE17TTayabHwzs39tIUO98z1Abhdq7tlOQ3Qp1zg32GzsjTN0n52i0XivxW1cD2u7/\n/sAZ+Pd+Cvj9P7U3tT2gW6B7HYlkZ8xorOevRIQVZBsekladDY8XRuDWD+spoNr9reeNw/1Is+KS\nCO94ZMTVUo8jfqurNvGPgOwl6Q4qVhq4AK2r6K5aW8er4bZAyyuXl3yP9yxpRhrgGujeZ+PEWajv\nK5sRJvnnAfR5GdXY1TZRA4uMa51rJdAT/Uztt+WRwL5399i5vgZzq1C99Dj/ulZjO1E2p0zIwE3j\nVazjr2nqGtwB6wsdPbFQW81dM8M12IT9+vrVQgNQ1+pE/x/w69P7f5fBlgJpRdoCvceDB3dbHlFj\nWxkuucH7al3vACajCr2lsVvnHFd+cQjXVE9bRs1xPUN+cWPaxumBPQJ06/UG4I5wA6X8gLUcS7le\nG/c1GtuGqxGOQFtLSC9uNuHXqrFFepoVuK5596C2z2h05y0o7f9HQXfEqwLR0nfR1jbeFtQe2J7G\ntospFuaaBm9pbP1xBvuMu+TRS+cw1CPgtprW1TKqsW1ianFpuMfkFWjsGgm9krOZv6Z2rkhaSxP3\nQB989KYROjDfFfBRqO2usF6cNfPb0976WfpNL/HvNXEG7Mr6qvh6AI9o7OEH3UVzjwIu/xsH/IHA\nHlF/rfBdXeutYa9nrBTItY/tRAeUX4/a3MdEqxkZAKIMtfw/Z4ecdx4JqD7PM8c11LqkWvfoez2I\na2FvjN9614RUgLDm36tabXrzxlHxmw3dawJ3cW7EtTbWa7PROV+LU+S6wdsjge31VjZsW6/3fpDd\noGj3V3nnrtjQ6JV57zWl0f5CGiMIkbILBA60SSpNGegpfb9L79OkDH0Q3yTLige3l2Sr4WW7qN73\nXTO5R0xx+Q0WW0Sb4s15km/flm/g5m1sHACeKPvpf1xcPl/g3k6suc2xxk+r+VXr2ftHqz16/5Oa\n8l7hAdaFQ5HgnKvLKwDbg7wGtW0CrY2StcJr1U4jia3kjADuCFPW3JQbIye4KTfkELjAjUCbn83Q\nDV5/ancEbtEHGm77f89Jn6LXwEfg9ibOrOzST1u32fsaAEy6AyREBXYMpKyfXDgwP8k3Wn+tZtOt\n516ktXZq7SULuz7WMI/D/YrBHilVD9IW1K177XMryav12CPa2xGZ/QZyVZH82A82GpuLduYEt/nI\nt274clvt8VpTy/8E6ou5zo7D7X5vvUZtw61dZx7cXnddOirJX1idzv9Wa9OqsSmDHZQ5PtDR7hJU\na4aeXtjFTSbstdUa5HYGwmpvYAXYAq6vqcsjgm1LZKTbHNXGrYIb6mp99u8CdfVRoq1zg2QUjR0y\nzJy1NE0EDqsprht6afzAzhy3JW9NcYFbNxX9uqWe6BLANdge3J5fg9qDu9Rg7qx21Sym+ASgwLz6\nUXyBm2jzzJ20ml0Nbu/46sh7Fqanqdncr4Ee6bVWeUVgexS0wGwB7XUENXcF4K0euwV1RURrp9+F\nzr+/mcfaCEAIAAcGQoY6N2RMvMmChruVKwu1Tp7VBwK03cDigW0B9tauR8zxHS9knIY6myd7qBPY\nMYglJGPrHImtEHIS0IN5tO6bkbbaaQ1mfU5qza4xvNYxtm31HtgtM6WloWvwehq/Q+FoT96Dupb1\nXF/l9yoojbHXhpoc5XG2aKmN1lYNX8+51borD259rJuIOL0UVgO75qxZ3wJ6U5Nk4HY1NoFzmSQz\nnMykmSpXrKC7dWHDPZgD/ALe5czL4Qjc1gSvmeIW6maCNvLIYNcgrxWG2wxQL6ja/TW4K8ntPXIU\n8ixrg6N8jHV2PBCQzW7WDVgg17PiGm6sftg/ctfPkznWTaSWVRmTe9B6ANfCLbhLURuovTG22Ous\nNHXMmnv9gR8DtNf4a02uVve1enYjFr8Fsm2/ehuPN5Fmr4M515cHAjuY4xbcLRB72tm6qeKrGrGL\npK0OttZ7D1S6/cn1DcwIq6OQwebkTwyeONXEITmakwtTcvMEzAE4EHAhZzIr3760k1jSRdg2Ja3d\npelosFsz3SK6iW6+uEKbrOGIlI85AFMAppxHmtRFM8DZxWldJkxDmVSG60/dB6w/qqvTSXlSjdbE\nPoQekR5WlwK3Iqm52nKX1dg6rKdG2/JIYMMkwIJeM1muBbkCN+k1IwU0/MvdKL1PcO56fK1OU/Rp\n9js1uNT4JkRMWLK70IyJIhAYmCJoiqBDAB8icGTgBqALMF2A6QwczsAyA3Fe0xSQbrfJkjbR6lZb\n7dh1y0EAACAASURBVFuHtU7xdpPB+Dp8A+CGkv8k+5vwBBwn4DgD8wzMByAck8MNwEcARwLPATwT\n4hSwhICFCAsCFsxYcrkK4HEDOG3TRmo87jFo67rXRjiXJNPazjjfzDI1qbtG+5U4XYq6BG0X6sEN\nc40vX0Jj7Gu1tl0fMos2xbYT+GiNfgRuawBsIE/xrVVhtfXa2JbiJlxyg0yD5giasmaeI/hI4COB\njoxwBsIpNfiooKY89g7A+hleKV2tip1a8WrE0+5aL2idUoMazvENAU/E14Dn8wL1cQYOB2A6ANMR\noAw2jgQ+EOIhQb2CLeUoHWUq27gBep0t33S6I5baSNuQLXWb3kJD7QFtodb78lrjbcCH+5WBbTV2\nTVfo/91Ha9fUa2uAmt1IxdUe7TWS7Czgq8YORVsL3KAICgtiSGDHAyMcVtVG59TQ5xMQD3lSbQJC\nQDJbGdvva4uY+q5VfyMbxdWmdVpxiojGfkIKaBU+TsBhAg5zcnPW1iSq/QbgAyVtPQcs04QYAhaa\ncDHlqYc5vEmFKFbVsY80rVGNvTMFlNZ2gdZQ12D27CQ4Ydtl7+UVjLGtX+syPcPwDoCTmOJqQTiY\nyr0P3DmJ666nrKlJtLVsRvFN8QUziJaUtDy+5JnBx5gGoBnscEqa7CDj7TmNR6eA9be/tLBvInu+\nrQ3bDVvT2h7Xalb8wmcG+amAHdLxQeYNshk+ZVOciimerBc+JG0dQ8ASpgz2ZEzxyZjgynpSW07X\noRPazWukbUikUWvrXKGsIdVAa+A3vQP2k2kij2iKE9GvBfAJAF+TY/yvmfm/2F5lwQb2zUeHa1D3\nznnEWTPcmOJ2veiuUNv+SDrrnC2vYVlz/IIZF0wpOYERJ0acY9JOh2yK36xQ4wDQYYVaJtHK1zrV\nrBcTyof1PV0AE27ZVMAe5Na91onpXVzY+nMApjm7w+qK1j4mjR1nWsfYlDS2NcWtxl6hVs1/xEDs\nTNvsrpGKF7gppAqgCelXFbWm9iCv2UMW2sczxc8A/hAz/00i+gCAHyeiH2bmn1wv0U2i1ZeLr0u3\nVdotTV0BWg9GC9DkV2LL9wAXi2vTktXYmrZAC9RWY4cMdpgYcQ6Ih5imikVj36ZGLjPkYU5QxwlY\nMthArt64Ql187PWBbTZejUm4Zt1bsK2NJU6PpwXup2EFW2bCg1gjGuo8ecYHgVu0tmjsWZnh02aM\nrc3xVWOvmttN+F3gZmSgs4uU2xuvlbLbMWAf1AJbj7GBRzHFmfkfAfhHOfwFIvpJpN/KVWBbjd2C\nulaqPah79rG5Ru/L1HA/gBmuzbr1lcK1QenJs1VjJ6gvNCEQI4SIaQrgKc388pHAeSqZ8pgzHHID\nz18b5Aw3lm3nUixD3sNsN5OMjJFFWtc2DJmksbEdW2u4Q540xJwsEsrrYHRE7txyecxbU/xCU9ba\n82bizNXY2L82220+tlNvtQ0peAF6U1o1ba1rQ8Osw3YNQuQRZ8WJ6OuQvjX817f/eRVgX6G9y2x4\nhrpMnFTgHgJ565IZLg1o1dR2fC0TZhfMOGPGBYeisZeJsRwi5mNIUD8BcALoBNA5hXEC+KTCZ2BZ\ngCUCUfkxApy3jEUkze3C3W8TXZGdY7uay+dvQlrSEv8Y8oRZSONreoKU16c5/FS5Z+l/fEOIR0I8\nBCzzhGWasIQZZxxyOc5lrK0BL3Whd6cBqZ5UhzysRzy/DIEs3CJS6t73azytbSHX8cA5fkCwsxn+\nHoA/yMxfMP91wmTC+lwNZn3s2UVet1pxu/E11XvjmtldcwI06Q0o240oSx5Pr1AfcMYBJxxSfzMx\npkNEPEbEmwX8lMB690lt8zUBtABTdoclg513q9CioOYVZu3fVwgKbuMTZainBPOcXdHSE9KPxXzQ\n+B/I7i2AnxH4aUC8CViOEy6HGZdpxpkOOONYytICrjcB7Ze/TAY803wYcmMqFbjlAfpGqUgLuQZa\nQz6pROoZUn38QGAT0QHAnwPwp5n5h/ZX/EUV/kYAvxF7mPXxNZramfnu2lAe2JUovfOefakahn75\nYN0+as3ubHpnLZ0a4hFnOiIEYJojpjliPi6ITwLiElaQ5ZWs7MrQTdrPAoQLMF+AeElQ8wWgCxAy\n2MUhA63cfcXuYSdz7jgBN7JOLbP5ea6AFNgkUGuwP5DBfkKpXI4TlnlKYIfUMQrcFxzKhORa5sY8\np9U0TwUomVCu1gxrc7XF/FGZJq1pa41LA+4BrXvwGtQ/CeDvdutoZFacAPw3AD7JzP+5f9X32btQ\nB7unqVslOgA0pm1ha5XSMwQ852rslK9UDav5beEWjX3OGvuUXaCIeVowHxYsxwnxEsAxvSQCwNXU\nBJTl0pB3pQnYfEmmu8BeII4G8kcAu7xPrY5ntfnkkHeXhTwRiBnAVwBktbX4bwF4BsSntGrsecZ5\nOuAc2hp7C7VZpag1yR7ctWbGBAReJ84Qse50bJniduLMg3tT2thC/psAfIs6/gtuHY1o7N8G4N8B\n8LeI6CfyuT/MzP/T9uE63CrBWinK8Sh5Le0tvbPR2D3zuwW0Ftl1tpkFF1Ncr1uvGltr7TlEXKaI\nZV4QjxfEqDZXSD06qx/y4gidgXBG2XIqUE9nYDkrqHM8BfDW7Nk1IjBnoCFgZ38+rOvT8yFra5kk\nOwAkIGunNfZTAj8JiDcTlmPAcpiKKX6i1DlejMa+YLvFVK9tF8tKytdrijWovVGgWFTSrmKunM04\n2dPSNcDtvaqgXc39AKY4M/+oeZoj9t81oG0JtjR3b4xdAd2qj6KtsxvR0h2NraEWM5xJj699jX3G\nESccMYcFh2nB5TBhiRMiZ7NR0qfrzRbZBNApQYw8yRby8XwC4jkDnGFGVMcSvqfoN7E2foZ7yjP6\n0xFlu6isU+OAFWKtrZXGZpk8yxp7mediiieNLaZ4W2NHUtpaOvpNRrBvgiM6ZMqNoLQJwnYroNDf\nAvxLYIzdl7uC3bN/PMB7E2mTSkLuUWvbCe2EegNmDbUo1+1GlOT0GPuyG2Mf0gRQuOA8XXA8zFg4\nT/YESq9CdcCmKYEss+ThBEwnIJ4APiRfw4xFhR8QbJ2eEs6+7CKTFztCXsITh7ewBdu64zojHo8T\nLlljX8KMk2OK6zG2aO3y1tcG7nbdDgNe+KS1vUjkhNyjWrg14KJ17zLGHpNHegmkV3o9uEfH2BX1\naytRT555ILe0tQa9ZFVra7vEpTekzJtZ8VPW2Ae64DKdceHcEENI72XLa4teMVqNfQuEW6SlsFuA\nb9O9LJ8bNd8vYv2No/tKLrsCtRniyNZQUptO1hc8sI6lNcxvrX7aJx7SPvFZaWxSHaQyxfWs+Kq1\nt0uQaME9qmM2jrZtY2MNtDS1nL/rGJuca/byCBr7LmDrc54Jfs1kGq1JGr1lCOYsZRNImuzaLnFN\nGxP8spk0O27APoQLDtMFJ1zy8RmH6YB5WtIL1hF5PzgXExeBU42dALoFcAvQiYHbrMVvAZzY/Roh\ntcButRNyzuWyZAN0CR8p7SATkG8ovYp5A0TR2G8R+FkOP8vhvFXtPB9wmrQ74hRucEvZ4WYD91l1\noOs4W2vuPL72fqZzRGNby054FWtQx1EKqNdb1ODWLHmVw865vbwGsGsZrWlt23pUK6Jdie7BrFWS\n93h7r9eo81oXc9YITLBbRmVcfcIBRxxxwg1OOOMWN8lwpJh3nzECM4jXbRSRA+jAoBsGLfl/zCBK\nH2WgA4NODDpz3sySjnHO4TN8mHNngWWfJS+bJbveyQn5W+hUylJDLvveZc93CefXMfkp5bXq7J7Q\n+r+Z8GJ6iufTMzwPz/A8PMWL8AQv6Ale0hO8xBPc4ia7YzHNt0tfameazF+wwF3J1DUmea3Juu0m\n/5O1BrcbTTTck7rXM8Ff6dtdkzm+C9jXQq1KU3/IzoPU09h2TF2rHLF8lEvrw1QaTWpAW229bkg5\n4rCBekl6hBiBIigkcAFZJAkIBwYtESFmW4AiaGKEOSIcI+jMW3eRcASds6YXsCNnn1ZLQFeV2rRi\nDJNtderzAeD83e8ENK0f9Z+QPpAgL7YcVmAlHJ+E4mSSLB5Dfpsrgf3F8DRD/QzP6WkCewP1TbaA\nZEJtVnMaebzN0wq1/lVEC/c1TbRlKG4svRzYfYzB09TWl5m52rjptWpsHb4G7BrUKrybsRG4Bx7R\nMsc9qLWURqE1dtIK69tbsxoHHnHCGQcc1aLMskJNXJaJZHfUggnTYUGI+RUHipimJUO9YDpF0CUi\nXBjhEkHFjwiXALpEUIaaI4My0JxNe1p8kAXwMoKzMOtLp7WNygcGob/7nd+jFj+9pbX6yzEtZcVj\nwKL9Q3r3+uX0BM8z2Anqp3hJT43GTpaQXtcus+Ss3tXmgJjrqqzjtzR2y4LWzdD+IJrXfuT3m66C\nW1hqTYa8MrBrGluH5XjUzqlArf8vmlrGOTDjnRFt3dLaWqzG5lVjbz595Ow2u1VQT4gp+gw1QlYk\nMcG90ISJc0xhwRQWTPOC6bBgOl8wXRaES0RYIsJC2UfacZbPJYAz1DFvoogMLMkv7Y7XGoI+50FN\nynAsGprWTwQryONE2YXVzQFxmhAnwnKYsRymtD6d/XJunrLp/RQvBGwxxbGCfSoaW7R2MsfL7DjL\n8ldejtQae6Sp3kVjF6cO2APaNKgN1JYlT177cpcOXwO2t4RVMcU9qEejbVYO4MINlJnwpLX9baQy\nxp5xxPqyYZ6vJaQYwtpJMKWOYqGQTfYL5nDBPF8wXy6YLwvmJWC6XDDFiBCX1V9oXaOOXDQzZY1N\nck6vb6vsMfujtu3S3lqd+vPJ5ZPA5RPBwBLWN7LSZ43SzH8ME5YQcJnnMtN9mbObZlzmCcs040V4\ngpfZ/H5R/DTGvsWNgnsLtTbHiynOqynOmqUa4CNA63bkKZJs1W1+bVEq2zW7NdQjk2OvzRQX39q0\nI6a4R2JFc5etT46mvsbqrwGuRaBWpvhq6k3F/NNj7Blp1nvCDdZXFLj0Q6tlL5NwAZcQcKAL5nDG\nYT7jEC+Y4xmH5YwlEuZISaMzIcYFU7YgprwhnJgy1EDIoKelUwJLOGdPA83A5gUlPYEsgJd/a6Cz\nhpbvfscARKL0YYT8OaOFpvIFlIUmXELaInoJadPJJRxwDjMuedvoS0rwvqSneIE8aUZPUliZ4rdl\ns8pxB/Vuz7gdZ6s63bSVljk+uqpiTXGBe1OiHsxyrgauMaca8oimeM1vEddbZLb/k7GLMsdHx0oj\nUFfM8aJl1XKXmOEBS4Z6wVmZ3/L+F4G3cziktqNywAUTDiHrep5x4RMOHLCAcMxaJ4LWd68ZYEkU\nGOCIIGPrmGAVwGWjiobaAm3bjgCu4dbaunwaOP9QXgzIHx0MuOSPD5YPJORPG53pgDPlsTEdcaLD\nukZNR7zEjRpPZxNcYHdMcT3GXnehTVg4gEVjQ0Hd0tgjTXTE2tOm+OaBNW2tnQe2Hn+/VlO85ntd\no1di3uKh0dS9L+pfY3J7zhsOifaTMXYMWHgC8YQQJ1w4aY0TDijbJUi/kqBMeTIdA004Y8YRBxxw\nxinPqa9z66dicGrDX/sTFoTICW7OvjkGFMxKa3sae2OOZ3+FePvNbwFczzd4vqwW6G22st6f5iRu\n8BI3uMXW9LbhE25w4tUcL7vReM5bdRPckQM4hrYprpvuqMVnlUe5z2vvZRznNKqRhNnety2v6PPD\n2h81wytus5+RVmff/vcKvga3FdvBFqCRTdrkYgxY4gRaFtAygZYZWHidzUZarxaghRx/U8s6o77d\nX6Xf5JbRu94Ks+wcgVNx5PXvQOkHCigiL6+hjF4kjxpuF2q5lLD+0CDp/dhij8CkZt6l9FLJoTi9\nTi0A3xr3km9wy3l2nLM5znlmfJmxxAz3MiEuATFSnkTE/pMyLQ2uwy29ZBVu+dihMcXFXGq285rG\n9sK+vEKN/UBga8B7UHsm+c5kUr4Wq6UV3ByBuASEhRCXgGWZ8hoxAwvjTIe8pLX9WM9qetuXRdYX\nRY445dCl+If8X3lPbIKertt+4DhIZ8JIM+/Im1uAsmauNTNJXk0xbCbL9HkBuzisb7lRskaiSmFU\nKZNjyUnNl6WsUx5H2/CtQM15VpyPOPEhgc0HLMuMyzJhWab05lwMpTNuQj0Cd61dWS4jUruMtI61\nyzh7FGwP5j7UwGv9JZCaGV4zwe2MeI5n8waXE/VDwl3AXjU2xQDEKUOd6i6EDJO88RNkeEeq4U8G\nagH7uNlqIUBrt26Y3O6Q1i8tFpjB+Y3CDLsBG8DueFMEqlxK0zK/S10AL/nTH3HUS4IprJcEL05u\nzzhulrJOeYJMzq2AH3HiG5z5WKA+xxnLIi6BHYuVha3GvtY0b0Gtw2LqRF6VDwesX8y4L9hfcqZ4\nb+AyMu04pTjs2/4t09tWjOdEWkMfY47HJQDLBF4y1Etazw2cTd8yvNqOqfXs+UGtectIep3bvZhw\n8revm+w/zFSABoACeMoc5arY9WO8LwZ7kTbHN9/szia4bIzdvhCzd9uZgX14HXuvH1bYhkVLZxeP\nCe54wCUeVjNcmeK8BGAxcNdebfY6fa8NNU1xrEqnfMHUgtobY3vX10z1rbzmH+W7qwP24+uKOd4D\nviUWaEZZL+YYEGM2v+MEXiib5snkBbjUwboFlRBp+5LIATPOWVeLyb2OTNfRqZrrhf5eixdOxcOq\n1BXQ5hMqmyJg/x+2GZUvsgrg+XrZ+rL97bJ9eDvm3o/B1wHIOiG2GYfzIY2r+YhTPKx+dgnm5FI9\nGVO8Ns72xFMCNaUhTsxwKbyN8pAVpJ7G3pW6OvfawAbqpvh9gM5Oqx3vK6StvqPVG2upamtkbU0I\nMSAuDF4AknMLbWeXgfIrkTLOPJiJstUo1TPeeuy8HbWuX9DWHztWb4cTSoMq2aJ10qwlcmtPWNWB\nHY/rr5dEJ7z91vo+l+vsgh57r2ENuGjqczzgvBxwiXMCeslj6wK3aGvyTXFgn/EW1E2tnUuRSbUt\nwmYpYhhsT1P33719BFNcjmtauzYwqb3w4W1SkaiUlvY+plAzu2tmuMjtc+AfvA/8nfeAz34a+OoP\nA9/1LvDt7wCHZ9mko/IBQloIWAIoxhTOkibLQtlPvprhF8yYMeNS1rz1ZlRtYq/j6OrXvHZuUyUa\n8J6FcoVsGNiZ7O0Ubn9MwY7F9Rh8u+lErxxceFbj6gPOy4xLhjt+4SXij/414C/+EPgzPwN8za8F\n3nkX+KbvBugDe63tKUgy4VEzvLCZQS6Wpe5V7wK2DnuKdCuvUGNrU7xmF4/SaJ8zeFsNZGAtuz/7\n/cDP/Rjwuc+k7wyJ/O+fAL7qbeA3fAz8R/570ELp4wUXAk8BuDAoEIgIS0yvWGay0waJicCc9kun\n8aQdO68Gqf5KtjdKTdnIo1oWA1gtq10jI+r5ykhlvA2BOfcq/TH4Ogfhr9DnV2mWGZesnVN4nTA7\n/wf/PvATPw78o58FLqr+/tKfAb76beDXfwz4Qz/YnkS71jSvOulRranU7RWwhdgLt+URNbb1PdvX\n2jg1uGuPpXbBWvE6P+1unwM/+2PAZz+1v3c5A7+Qz3/hBTg8ybsFqTyLARAn0DlM6be5JmCZkGam\nJwKHgIgFkSZELJgwYxGjlJR2Zj06tWPoPMOtwhqlar0P2diVsrPiXVPKQcbdYj9sAfdnBkL5Rrts\nz92Z6bJtN6YZ7yWvV5cx9RdfJqh/5qf3aVvOwM9/KuXv+XOkT6Fi8zGK7oz5XcAvJ1uF32KCnfN9\neQVg1yD3AK5B7WTIfpiuJ94Eo3U/8z7w+c+04/nFzwB/633gX/xet6yZV6gpTFjKCEKAXxARE8yk\nTG1ax9V63BzAIN6OVJMSyCDnyTDru3lulUlPGiDb8Dq5ZgCXr7silJykj0BuZwpKSbBa1GO1uBen\nDLbMfM/gOGH50b+WNHVLfvEzwN9+H/jmj/ehbpnqcM7VyogkUBvE63CZfcMWan3cl1cM9gjIHtQ2\neu+cua1ya7MX/nvvbc1vT5Yz8FffAz76vXnfgVLXQJpcE40dCDRRmkyZQtLW04RFRpa0fkhB+1RM\nbPMhXfnSipjf7AAtvxKg8wuqTxC1NHytDL2+euOvYKdH05qEzddd12+S6Q8Pyr55/QWUjR8nLDEg\nLqvPMYD/4l/Ymt+eLGfgr78HfOPH91B7S2GjZnpDB63/13D3uKhB/drBrmnqa+HuaO7a47XUGrXu\nmSOAz33azd1O/vGngXNuvGCQPJDTG1QcgJhN7/TOcgBPESFMCCGCKGngS/bTRxei1mM5pyvMxamf\n0yxmtwBem2+plUFPo9c6yRbYKsxqqCRaXMKbXJVjMePl4wgCdJ6fED8vYaWtomuYf/ZnnEw48k8+\nvWprrbVbIPdMcCteUy1wb05grP17mqsur2FWfMQED9hnopIhe0kt362xtbiv+HDlZiO/8sPAZX24\nvFxQ3ncOIf0ypkAdIuKUfmGTZFcaZB+32immtqAKuG6YsWpmC7ZnjXhhWy42PAJ1K2zATo7WXyeV\n7qsSLm/PZcDLO9VYVyRk9x9HAi+UZr9H5EMf3pvhvXG2lVpHONRsvUbbUnx2nN2XV7TzzIZrcLfW\nqQYfWzsGfLh1JX79u8BPfqJtjocD8NF/I4OdH5Tv57z0xRMnMzykHWgcAigAy6QAJgC0HgNYw5yT\nL2Y2zDnT6KjkYf+/XX49kO25VnVKuGZMKYhBSNujnf+tZrkUodLo5d1pFMiRX1lllo0mlDcK0br5\n5Le/C/zlP5PM7ZqEA/Dt7+4+zzxkfl+jtbv8jQJtoZYxeFseGewa5Ne4Kx4/clutsUcAX/sO8IG3\ngc99qn7/h94GvuGdZIpbUz6ifDIIAeDAaT9NAKDDGWpJK+ljSX5Jp2pFRVur59rOSf9Idskn7Rul\njV+f96rN+rWqopx/23+rRQ/WcQHQ36sruc3p3r5qSeXXTVC29q5hfNP3pCWtn/8UqvKht4Gvf2er\nsT247Y60FtC19uZ1ik1pAW21d1tek8a+BuIrIB/tbT0wIoDjM+BXfyzF9QWzjh0OwK94G/jIx4Dw\nDDireyekhrVg/Ryv2gnHqlHvDBKvWFwAVdimuwZ5S9u0/LuC7UC8W8UMJg4v3laaTUe6Oaa3gG/4\nWDr3i5/Zau5wSFDr+rtk15tEGy1TKyOavclArefss0B8z59fJCIGftxJkE2E+Hp3mbfTrLX7LKu8\nWZ2azSUz0q9pXOOO5hjPgZ97H/j77wGf/zTwKz4M/JZ3gd/wDvDWs7HNcV7DlnCtaMSvmdFy7DVq\nG74P3F6arjG6avmulYFXJl6H5nVsXhk8fw78nffT7Pc/+XSaE/m2d4Ff9w4wPUswC9jat+GTCtvj\ni/nf2fxPLAIdlmO2C+h2Js87r491RX0ruPxM6yoPBPbf8P4Dv1XoLaQeyANge0Drcx681wA/K382\n566BumaO6iLyxIPaWho1rVUD2zMta/6oxq5Nh9iNVPbctUabbaI1C0XCLV4ELnEW8BawLZBrYHuA\nXw22vaYP9gOZ4lZa6qgVvkI8E6gFQQ0KPStqG54n9h7nG4tV09SL1+tXPRhbYLfMxlocaPhAvfp6\nefPA9kAfLR8453r5qzEh4Yvjap3AiFluz9k03k93Khkzw4FXuqVUfA/0e0jL3PSA9uAQqMWvlR87\n90zGHx1D6zitX0u/TXuv4dUgrz1X/BGN7VkiFuzeh/9acXjPlXCtk7Jge+XkaW0P8FoZe22pBveD\nyN14eQSN3WsVvf9dKS2wR7Sdhdr5bavNc6ym9n7TvNef1RqkHNcmxKyp6eWl1Rk8tMaugdn6yGxt\nqGJ9+yx9bMvOhltlo7W2BrylsXW5tsq05R5MHkhjE9GfBPAvAfgFZv7NYw8bgfueYgusp62tbyvd\nNqpWvB7QWlPBxNPTODXtWtMSLY2yOPHZsOSr5veqr7KM5YLtfXzWmuvkhK2V07N4dLhWdp7WtmFP\na3vxjcLtpbMr9+dnRGP/AID/EsAn+gnxzj0y5CMaugW1B3gtfrluZDnHisQ72gG1GpXVJBZ0L982\njIav01vT2K3JQW94Uhuy1Ez1lskuUoOlVue6/FrOdgA9uEesonvLddx0wWbmv0pEX3f3BGnpJcqW\nQqX78xrkXWAWkLW558Vn4+g17F62Wqb2qKvBLWC34K6Vofj30die6e1BXTPRp0a8XsdZK++WxdZa\nVbLm+Wh99BSMLuNdYmvq/c4q/6HG2PahYi9Z397TUrPa7nJKqeSd9rfVINHjaA9onTR7bzBxeNrE\ni8+rm5omeQioRzRKLV0tjW3PeRNmOlyDuve9Smuqe+FaunR+alrU08berLlnrvdmyz3QAezeuGu2\n+1ovYe9tyyOAXYPa+h7INtzJAGMtNAt4T0t7ZjfU/ZOJKyjfgmzP9eDxzLprgB7V2C3t0UpfC2oJ\n19yoxm5pb89Et4DX0mebl8dJr2P0NpWMAL0rawu0ll4l9QBvywOB/V+p8EcBfAzjUPdaYENj2z3Q\nI9BooK35LX5Q92qgrbb2wNbx2SxI+BrtbIEe1ditIrR59o5tp9cC2567RmOPTK55kOtn2rDkpcZG\nq6P0tLY3Q37NOrdAvkmcDo9CzQD+OoAfQ08eCOw/oMIjGtuKZEbf33Dy9UcJS8FF2t5uYRZXW9KS\npGiYbdhC7GmtXv80CvBdoF7M84B9Gmx+N8eM3YcsPMhbcPU09gjcPbBbddDTEa3y1tB6gLfu08pE\nv2mnldGw8qpp6e9EUpwifwyejCx3/XcAvhvAVxHRPwTwR5j5B+p32G5/RGOjcr6lamgtKKatdvU0\ntD72kqydbUga7J62GslKr1GNhmu+pKHl2/y7qpz2UEt+W/4I1L0lsdoEm1cHXp2MdqxeRzs6a+5p\n6Wh8D/Bdgu7i2jIyK/5vdWPZPUgPNEfghgoLRa1r6P9r73xCJ7mKOP6t+WnQaCQHwahZjIiCkTFn\n9AAAEM9JREFUiGKIihgU/7OKRE/+ARE8eFIURBE9eRNP5uAx8RCEKChKogdFk0OCsHF1l8TdGBQ3\nIYkaBY0aAmrS5WH69VRXV9Wr7t/M9ORnFzRd73VP9+tX9al6/bpnZt1hxfHK5gZ9iHVZZm7rEsTh\nw3u8clle1rJGUtr4UwDOQF3OI69LX6d58UIls7CR2jA9m7E9uGu6dQsky9pdtB5BzYhnzK2gzIau\nz91l7UzUj6JSTvYweWbVabiB+GLFwsDmN5oZ65/kEA6oM3cBTU6a6SYzNpNlEmzppNa9tLe2QJbr\nsaBmht86Y2vJ+0RerAAJjIM688xbb9P20WADCioM3SwC0pst9+zjZe9u8qzWmDFZeksZOyfeiSy4\nAbuROmrJMbDaNvgyi7BmGZYXoOVaDsPlaVdiXU5pOY81I6uzFTB0lGjIV83YHIDO9udSFE8h3SNZ\n7aMnHI8APE1BNqb6yyyyLrJNGVFZoxa5ztgoG0g9oDudN7ob7a1sXYPclx2DrfexwqnubRL6StQ1\n6NMkUmd3ehrGBw102V/CXE5jvcvsgY1Aj7KB5TiDsgC20411mBW0Xtum66Kba6usdAa6ic1GRUEW\n23vzGEZQKDAX3bo1ssCOLrfGVRRsdSYfQA1hFw/KaByfBTyWHYMtYbb2lwB7mbvQpcFeqcEA9Q8Z\n3UsXh5C6N0nmDfWk6DoPbK1HEzE9gHlT7nTebBtkgwhya5tVti6sdg8i1iwXoP8n8GVh9CY+VwCI\n+kBLwCOoLbC9y4oSpTe89m6NLBt2Jij20SepRgWncbOC7cFs7V/L3JIoBTYDmz8Tp/V9tzUMh1gz\n+vfSGuwI6ghsLTWgo7qekwmoWeqN2tZ+yIXb0+HotWGJnmhQevdH7619mraMtrxqtzUkgKR+P2ug\nZX0N8JpYvGTgrmXyXtcrqHtwRyeIQO9FjVD2NBTXzmJlaT301hckpRyrPUaZIS/HK8+z5cdqMGeh\nJnVMSyJ4M5mchVNwY5ebRmyTjlMuOsoCuvFWnY6Kuq6ylIDbrNZQr1atzutt5c/qmpX4GPUhrWVm\nD3JLtN0sfvTwvMafyaQIsBLw0CHGZOrZwa5lbgmzzijy4qyMUByk3a/Ty3Cc0d3XWYf0gNZwe5B7\n/Vrqa5nZc6JeW2XEL/AKkNGIurERXzYWTl0EttVJusNWLdzt0hRbtbaj1XpB0+rtqEv+Fnk0QebZ\nroDtuZ7MGV5X1UZVXl3vmXXbn3IUxYzhc8mpkMeyx8mzsp+ettRZRANddKBnrY5jEp0ntjXtDtqQ\nK6VHUHv+azEh9dooKxoG9rpDQF1A7q2fFuWnBeD6oJ7TwFnr/o4ytTHMYXmDvAL4CN3vJlPrA90h\nV6Ke4sPXsrYHttYz3GRHWN2a1UhZHWwQfDOO4UFelz2ALR0FiIfh3kU4jsbSSuJcJRYwb+7dSjMk\n0EW3skK0eEzIy6nZMFwcg3bOoWZvBj+Olx0W6IsoegbsSjrlBhu4eQ237KgSlHV7an0fAS2h751H\nrTPs6HXEnvwt985+MA46FuooIMeyJbBlFJEUSai1XqjSabQcw9P1WmX23sv2tO70bjeVwS3gG8RO\nVSRKeD17sQ82tEMAm/ux6OG291qU90qUbpQFcw1srUc3v0eb87GcrTzarFnXiY41/x6ZNqcfwEzD\nZniXIU6TgtuNiwbM0q4hyNFUulX2grIvO5wV12sNelkXmoouAZf7WWJtM84x/HVWATL1myeHchpo\nnbEh9B4fBWYWtuLNtrJjzyF0ZBgLtixbaUbqXUOtxhsXDVVXm9kq2dpycBnMpV7eZBFksmoDUdvE\ndi3nUbo7PIqBtgZ5mpdBomS7roNYrzXEmdfVIsCtQB3LDjI2YJOgQdZA62PoY40VQvcPcOVRmD5e\nAwwmbHRcsPQiHuCdMxSwWQFcyp6Hjfn2R/RQ1XOKCGrRf9WhuAe4BbIsS4gtvbVd9xfFKsp2jzWL\nfbHZz2q21vXluqDrrCzrix2B4a1TxnaZxQO8Lnsaiusxr5etc43un0tLex75LLUYXwe6zhEU3Kis\nB8eROvtLAwyeOcvnm/0/poL/alrtZeZsxB+bscu6lrHFULynS4itrN3CXR6HgTZ6B/dqnalXbXtK\ncJZQy6+dWnDXLl/qEmZIXQXnweuimRfOp0CtA7AtexyKA/1obwGtJ9Ak6DrTRG0ps7MF7pXaRQAt\ndRn9o+GcPp0s9JxDw90IB2gE5JYBvYjvAW45TBTxp2ZsIM7WxZ4SaD1Ki7K2WJgAkuviH4Tu+Xhv\nUc2EUy4SAd6zpdZVUJY2HWTssVDrrB8FZl/2PBSXkGeytVevLagsw4zuby15BVCDzT/jlf2oXQm4\n9eEjfdC3IqJr40uAdVZ2H4PUXlb2nEVn7Cjie6nKuxeRdd5zJ5mx9bqSqU24Vy3UpY3qHpxWbZE3\nmZrkNajmF9G5yNW5H6hlPaQ9hZ2P9ZWwbNaOZY+z4jWgrYkz+TlPrItctacuDiBHAu0x5fe5o1PV\nTq/bwELpvZxggSwMzRmQPbA9p/EcpTTUS1nywi24PaCtjD0GZlHmFUBCL5Nx5c01hgAe2DwnV8E7\nI5YLWYCXwmDYbQTsSZk6A7WV6Iayx6F4WUtjQ+lWWadJL7zqfRj9Ry1lUW3S3+eWkvSLYRMMw/ey\ntDfpZU2AjSl7bzVZjuFBbQVjbyhey9gWzI3a3yuXR2ICahypfi6KuNUivcMx4PZGNu5kWcYuRY+y\ncWQ7zYAte8rYVgaPsrXcZqXU9Ph4WNcV9W0CMHAEt/9q54mMLAH2HldNgVrWZR0jm7H1ugYnB2UB\nrluWAbl8/kg00XNuVvbV1yTrI1EgD/SMDaYG3gzwewM7ytje/sWgNaAzw3Xr+HopztNgOLy0hpu1\n40d6Bsja8+gxzhMtlqPItlpQ636wMrbuW6vP5WKNnnRdZh8ZJHQAKXVe+7XuiYZZlzP9HtliDNQa\n6L1m7OhEFuiyowrgNaAjmCNjWU5RPqOX2rHkMaN1FuxIH5Olx0At+9GCWgdkq18smAnDvvaA1BnZ\nKkPUed+19cCO2p+1rwW3lbG9kZFnq7FQe3DHsucvgRSRMHtrBLrXBmvRDicz9jbAtuAeC7YF9XGA\nzoIt2631KMtZfZvJ2tlsDaNOw58BO5oAjCQC28raU7L4GKh1tt4b2NJZvHtq6zMa5qkz55ZYzlfO\n5YGdnS3zDF/KY2e1azOo24TacwyrzntM4MGcydo683rZW6/l+TJge1BvA+wajMeF2srSGu5Y9vQL\nKh7gGupstq4BLp2AxL7yXHopwGfEip6yPBXsSLfqIkeKHE9eR1Zk30igPT3KdLV7aw31ylhnwI7m\nUiLJBKUoq0Y2yUJtwZ1Nanv/oYVSJx5PdMAVXTtQlK2t4bmEmdA3dpSxveBjSc3oUx5RTb23thyp\n5jC1a6v1g87Onq5B0EBG99pHRv1UsKeMynSGtKD2dMsuUdnaZmXpXLYGZvk+NtDP0LKsszYSZWtS\nx7odKIs8twV2xqkjQ48Fu5alx75XXMvWuYhf7wMvU1tgN2qfTLa2oJ4CthfIa9dnQTWmLhNgI9tF\n56nLTJNnRWS21nUyu+hy7VGGB+wqqJ8KtjaChm7qo6vM57IOsm2wi5S+KoDJulIvt3mZJ1tmo6wz\nvQzu0SgtkikwW9vHAu75k3WOWPb4N7qeWHDrKC8Bb9DP2pauI7Y8lxe9dTCIrjUydJRZpw6ns1m5\nBnXGKYrNavtIO4yZG5HBWW7T5VKXFQ2+F7S1r0XHqwEc7ZOFemoAiWVLYEuZkr013BrqMc+5szB7\nWd2aF7Dam4F76jPnLLC1qJ+J9tb16TqvT/Q6CzerOgtyvY/VRs/XPKi9oG1d35TFs0EE9RjQZ83Y\nQD5bS9Fwe1mhlh28mW8rUmujS7AjAKIoXcvaNbinZuqMk5S1dU2Rw+ssLkGW/V7qgBhufY5a5q6J\ntoMHdamT57LaE8GUAW8M1LUAMSvYWrIwa9EOMuU5t4Q6k7nHZGtp3NoQbExmzmbqMYBHTqKvTesW\n0FL37KQnP2twbwPq0ja5RFDLWy/dluiYVl3Uxxk7ZSC3jluXHWZsq77URfc5x3WascNwGGXddq3X\nDLoNwI+bqXXbdFlfm1zLfhkDt9RrdpoCtefs1hyNB7W2s6VrkKNzR8Hds2EUfDMjgbrMlLGlI9Q+\nr51GGiJymmgYLo0MQ4+yWiZiZyGeCrcVQKxsEDmKvia9tiZBtVi2yYyqdObJQu35jL62GtRW1vbs\nH+le32YCcyb4etvqUgWbiE4DuBHrZwo3MfPXh3vlTtaXCG6ZoWVZZ23AdqDMMFxnJrkGYqevdbwH\nbWTs7D4WzBmnkcu/AVwCcAHA4wCuBPBqAC8DcBlsqCO4rf4GcvfZ2cwdDZX147TSHutWywJ72/bP\n2mwK1FsAm4iOAHwTwLsBPArgl0R0GzPfnzp6T0qDvOEuktutDobSrUm7swDeOLKt0ToyhAVaJkLf\nB+A1leNa1+xtt/a5FcAjAP6JPkDnAbwAwEsAfFjUa7j/AOAVqr90f8sJKDh6dA3WBJbctyQFWYao\n+xWAN8AftXlBSwavsWDD2e595jyA1yU+HzHgS+2h3psA/J6ZH2Tm/wL4DoAPVo96kHJ27gYk5MKO\nj/8frKF+HMOs2LT1j7b7eXJpN03bqvx67gYk5N6dHr0G9ksBPCzKj7R1izwj5RLWmTqSfwF4aA9t\nWWSXUgM7N6Bf5BkiF1CfpGoAXNxDWxbZpRCzzy4RvRnAV5n5dFv+MoBGTqAR0QL/IovMKMzD/7Cq\ngf0sAA8AeBeAPwK4B8DHpk2eLbLIIvuScFacmZ8ios8A+AnWj7tuXqBeZJHDlzBjL7LIIs9MyX6H\nzRQiOk1EvyWi3xHRl7bVqG0JEX2LiB4jovvmbosnRHSKiO4kogtE9Bsi+uzcbdJCRM8hojNEdJ6I\nLhLR1+ZukydEdERE54jo9rnbYgkRPUhE97ZtvGdn55masduXVx6AeHkFB3b/TURvBfAEgFuY+bVz\nt8cSIroKwFXMfJ6Ino/12xUfOqR+BAAiupyZn2znXe4G8AVmvnvudmkhos8DuA7AFcx8w9zt0UJE\nlwBcx8x/2+V5jpOxD/7lFWa+C8Df525HJMz8Z2Y+3+pPALgf69e/DkqY+clWvQzr+ZadOuYUIaKr\nAbwfwE2w3389FNl5244D9vLyypaFiK4BcC2AM/O2ZChEtCKi8wAeA3AnMx/iw+5vAPgi6g/r5xQG\n8DMiOktEn9rVSY4D9jLrtkVph+HfA/C5NnMflDBzw8yvB3A1gLcR0dtnblJPiOgDAP7CzOdw2Nn6\nema+FsD7AHy6vV3cuhwH7EcBnBLlU1hn7UVGChE9G8D3AXybmX84d3siYeZ/APgx1t+yOCR5C4Ab\n2nvYWwG8k4humblNA2HmP7XrvwL4Ada3tFuX44B9FsAriegaIroMwEcA3LadZv3/CBERgJsBXGTm\nG+dujyVE9EIiurLVnwvgPQDOzduqvjDzV5j5FDO/HMBHAdzBzJ+Yu11SiOhyIrqi1Z8H4L1Yf6Vv\n6zIZbGZ+CkB5eeUigO8e4EzurQB+AeBVRPQwEX1y7jYZcj2AjwN4R/sI5Fz7HfhDkhcDuKO9xz4D\n4HZm/vnMbarJId4qvgjAXaIff8TMP93FiZYXVBZZ5ATKsV5QWWSRRQ5TFrAXWeQEygL2IoucQFnA\nXmSREygL2IsscgJlAXuRRU6gLGAvssgJlAXsRRY5gfI/5/+sVSGEhPcAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "atom = 5\n", "plot(si.positions[:, 0], si.positions[:, 1], 'k.', ms=20)\n", "plot(si.positions[atom, 0], si.positions[5, 1], 'g.', ms=20)\n", "imshow(rho[:,:,0], extent=[0, si.cell[0,0], 0, si.cell[1,1]])" ] }, { "cell_type": "code", "execution_count": 63, "metadata": { "collapsed": false, "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "import bader\n", "bdr = bader.bader(si, rho)" ] }, { "cell_type": "code", "execution_count": 66, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "24" ] }, "execution_count": 66, "metadata": {}, "output_type": "execute_result" } ], "source": [ "bdr.nvols" ] }, { "cell_type": "code", "execution_count": 67, "metadata": { "collapsed": false, "scrolled": true, "slideshow": { "slide_type": "-" } }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 67, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPYAAAD7CAYAAABZjGkWAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsvV2sbUt2HvSNOddae5/bv073tZ2AaQd8Y2NdKyEXBSk3\nUjoSCYpDnJe8ICEkHniOggArL6A88GCkCAQPSEghwoB4MUqcyCJgwG3iTkAJmJ+mO8KOjR2itLvd\nTrfdfc/Za605Bw81Rs1RY46qWXP/nevWGVu1q2atmjXr76sxatQfMTPe0Bt6Q99eNLzuBLyhN/SG\nHp/eAPsNvaFvQ3oD7Df0hr4N6Q2w39Ab+jakN8B+Q2/o25DeAPsNvaFvQzo8NAIiejNf9obe0Gsk\nZibv92BgJ/oLgZ//1urbG2Gow95y2/B/FcCfanynRrzhrvn539kZAJjdb/8NgH9h4/2aOwpfS+dD\n6HMAPivu3jraU2fUMMAiZPp3rJ+tax8W7h3vjmhvG+ips58C8MMb8dTSYOnPhr5vRPE39Ia+DemR\nOPYWbfWIve/5Hvo50kJYesvIbf1q7/XEF6XRx996d086GCR+1lY3y/drdpnGWnn2csg9nLP1e286\n9n7ztyc9EbB7CmuPaNSqnEgki8J8f0ca1M8DoQeMW7/1hHun492t+GodzvJMYAyYMYAxYMrPaqdf\nyNgj5vz8vY00WNu79+bnIfQDnXE/tP63ytp30Nb9TvD+fTrumB4J2A+tiJ7x90PiAlJl91RqK95a\nhfeCfYtDfx+WsZhtCK13e6hMpwJ4xBSaxXXAhAETGMAob36vSdtDJadW5xyNi+HC1Nzf755769qG\nf+rO/fd0hLl/nT+DKH4f0fcx49lbqepuibVbFYsgXK2xAOueOfLjyu8tDhKncwBjxIwDrtkcccUo\n7sV3xhWjiOjaHbTyXOPaPWCt/V575z7i/J4OaQvc3m8vCFtc/uGKzicG9n1Ecu/XquitRhPFswfo\ne8TqPRWyB9Der5bGre+vRfEDJhxxyeaEM4644IIjzjhhAAM4AkjiOOV82rT5tNaeW+nqAfNWpxDF\nuyddvW2gV2yGCRe5W9/YSsM2PQGwHyqePUaceyvPP29x2V5u3NsxAA8HeRRv3OEoxx4x4YgrTjjj\nhDNucIcTzhhl3J3eTN3AiBnXldTg3b30WGDeSsdeTq7uLeBtSXS99HRc+zWNsVuVsMWt98a35z2b\nhvuAtgZ6oN5jR9/fq4Cz30XgX8ZPojRLInji1Le4ww1e4RavRGWWQs6iZrviANpMt7p7JSr/bitc\n7fdeqazmX3t/74zDXomtNdT60IviW9TTIbQqqadS7ys61sDZCrO3t7Zxtnpv/+293y3DENioxy6Z\nW7/AK9ziZf7ajAGTjLxVa36/jreVLrVbAN8CeyvuWrr21H/ruafeWh32nk6jv419SBao1Lj1nt44\nCrMVL6HdaPzvW+/2cKvWO1HYKD0Ifm/FX4bVOetRJrJUkZa4dxpvH3AV8M95aqxMTU9H2lOvrfy1\n3t3j3lOHUXp6y/Y+bdT798bTpmfm2D2Fou49jT9qRHsKfasAI05lp6W2evdaD97zLf+bTW8U1vt5\nLbm+H5XLlkHwTs3uBXgLfL2A3NsGLPXUfS1cNAPhyxpo13/0rN/birdOXcAmov8XwG8CmABcmPkP\nbLyxFWPFVncLxN5uNdLa91ruVnp7x9q1jiD6vQbOVtqi91t+NZExAq3/LfLHht3r7ulget/dC+w+\ngJTUAp0NQ1h3CDUdyNZzqw7r1MuxGcBnmfk34p97C6kHbA8BdW+cUZp7eu5WGF+ZNSDvjXsL5Pb9\nHnDX4ukB+WOCuvV+L5B7fo++4/1rVKuXFpdthd3i5DZdPZ1InfaI4vfp4iqv9fSgPY2ht7FEcbbS\nElFUkerfEstR8fO0B8DebZ+3wO3T1QJJT909Bqhr332In/9e9NxLEXBrfh6U1t+H6a2vVh3G1Ks8\nYwD/PRH9HSL61/pe6W0Ykd9jgLpW8T59Ww27Fm5w/v75Pg2yp+HX8r9VRj3l48smCtuKb6+7Bsit\nfPfEVauf1js934vKpKd+WumPwlm75lenXo79PjP/QyJ6G8BPE9HfZea/EQetfThKpHX3FMre32vP\nURr6CiyF2zvPSEh7r/03WiJWT8+ubgS/bY3Dt8Ae+dXKcatMe8AS/dbbGWwBJrJ7qVdBaePu0aHs\nqa+Ic7epC9jM/A/F/ioR/WUAfwCAAfZfN6G/D+tdSntAvRVuD8CjeHq//5jESByjpWjbGme3wO2f\ne8Gt1As6oK/ce+oJQfjHBrxPF7COu5daYPNhvHvr9y3R2/72iwD+n83UbgKbiN4CMDLzbxHRRwD8\nMQB/vgz1x2tvV9z63KoAH+Y+7t5v7a10W9D2OfKfA/9anBFXb4Xp6b1rjWkLDC0QPdTtT0G5D3C3\nwqJhPwZpnT6kDiJlawRqG887KBmnZaoL9XDs7wLwl4lIw/+XzPzftV9pAVr9toC31UB6gd1r19La\noqjHtRXuObUHujcRV/fv+Phqmvhaeme09QG1svbv4AHu+4B8D+gR2N7dSxHQahrsXuqRtuxvcL9v\n0yawmfmXAfy+dqgaeGtha6De0zB6wrXi9+m8D6hb2u894I3iroG59a0tspr7mhmc/Rjg7QHjnjA9\n78I8w/n1UKt+ayAH6iBtgXnrORL72/SMu7tqoNoL6tq7W42pZdf8ahRVuPpbQPswEdBrpod7w7n9\nt2ruPeXlAf9QcD8GgG0cNY7vyyNy16hVnhG16qMF0JrSraYj6UlLoife3bVVwD2NqububTCt+LfS\nuEVRDxuJ4xHQ93yjRxRH4Be5I2BHovlecb3m7gFy67f7gL1WtzW/GkXl521gDUD/jb3c2vrdj3M/\n05ln0fNTgHrPe1HaeivdF7T6+Qr3oGyZuRF+SxSv+UUNLgIsOv2Gxu9bz70ABhYu3NO5bLUlVJ5b\nFAHaun0518Jb6gFzzT8Cd5uemWOr30NA7eO5TyOtpW9P5UdkCz8CpQWuP2RoCMJH8dckgC2ge7cF\nai+A0Hj2v20BufVbzVjAb8Xty+Eh1AKuf64BFojfs23Gu6MOxccT02sYY0fPe0C91fB6K70F8F6q\n9bj6W2vsPBi/+CSxOJ5WmJpf1BhrHWqrnAfUD+y/D7hbYR6Daz+EPIh8ebamMVvg9mG2wO3D2rTV\n6TWOsfV5q2G1wt0H1N7uyUeLIkCpf+0bNW5dM/a9Wofhv8LGTvuvCTMA3YfNAlM2iy1TORGGDOF0\njJIefpj2badYdYe2HqKE/NwHwPtOc7WUZdEK6fuCvAWiHslqa02Cj68GaI2vH9TAs3DsWuNu9bJ7\nQd3bcweg7imnZhhT6Gz9jEdxtZI2it6xtwf3dtjlQKNZzg1Ptj6nY5DOcqjChANmAXcCjsJnlAMP\nT7hixgWMOxCACQfoeeN6culyHMN9Oew9OgOKOmz72FG5PUPWol5rY2zr1rUCfePh+IO1tL/WMfYW\nN+wBbuu3PaB2AKcgbVv1v9k+pKJ9fdvfNJL8e6Td3gJ2XydAAlQ9sHDMp6GkI4ZvcMYNLjjhikOG\nJxtgEwYwDphxxATGFYwzIP4TjnKw0ixnj6f30tf3ArtXyjLGgnarqT2UVhgjgD339Fy2F9Q1Tr31\n2za9JlH8qUAd9eDaAKIefcfzJnlwK5B5+XZUd0VvEAEcWF/e1zIKahbReTm0MJ0WrocLX3DEVSA6\nY8xJSxx7ROoYjulsDSRQE0YwrgLoFEPqBgA2XcpjABul23NnW5/2Nwv6HlrVh3uOftdvcAS4sHdv\nUE30jn7rz9drEMU7uOqDQO2BDdMIqJGcoNB8kiIq6jCqVPFj+UYhwdnGUePK/mNtUCvHVWAf8xHD\nCdAJ4BcB+yTAtxybsHBsBjAJpC85zkvm3wrq0Ywo9wC7BXBT6FR5VoBt1VNP3UWcOQ+lXDysnXgN\n3Hupl1tHw4CYfpuI4vdpKPI+RW73iZZfK0vq78FNAmSm5XdrA43GEXFqwjagl/cHlKL0CVfc4CLj\n6rvgWp/EsRO0BoFq+j5hypw6id4T7qQTSB3BiBksnJxQarPVvi+wUemUjduCXIHuqqOLbPGr7d8t\nBCtbf7UX9tBecLfpt4EoHjWCTmAX4lotGbROVpRU6x8xZSXFmu/9fWdrG8fqI72A9uPvJBInuC2i\ntAJbzw5fFGt61d6MaIydgJ3E+hlTDrlw6tR9THKZAIV1cV/R3JVfUTcezM5utXtfB0ADzFpv3Kg/\n579JLWB+qEVxuAREyNiqja2K7wR2mBwq/WrJ8e/WshRVaNGL+7BBeVQbRo/oTS7sDCtKL8A+4xZn\nvMAdEvj12lx7pa4FNjBgBkuc+gZgxe+DXOG3jKz3nSyz1TE72xVbWXeVDrr2HNZbxT+SvMAufA/o\ntPPdmgr70IrikV/teQvQNszgnnsaRQBkII6+1wbWnav6Fb06leE8g43qiIG+OVIbzs6DJ+gt4FRF\nGuOYlWGtTCwj7XW4ZF9kTjvNby/LVsi41qZ31VgA0K262+qcEfht1YN/b/U7mXdbILcfUFD76U5r\n+wT63z4Uonj0e40d1uLoQVnwComjp1FE7iip/pNh725s72b3ztazT0SWAgTIHHFxBZE1o5hl9NxO\nvPvuym+smAawQ612pd69d+15q/6ipGs2I7x4OxS/g7jtD2H9+Q+2bJ/gKBHb9Iy3bdZqqtb97jFR\n1NT3mdpzK8meQq6L+lA4akT+N/3oSusatK48VlfjgeYB3kq0+374HMVtjQW0vOeBTBU30K6LVnNp\nJTnKqu9Uve3B3QI4E0ARuPuAGCfSg3sVeZVe0xi79lsUTweYc8Oxz0Hwnp4+csP5tyhqGLiHW0Gd\nG5HkiRml0kYaFGcPrIHsuWor4ZZqmfVxR8DGAtgC3MbffmJPv19z36eOrF+NQ3tw298LNYcE8s9F\nxqLERD2FpShRbXrGMbb1qyGmhirrF4QpoquAes8CpyipvY0G2Afi1rOO0S2gNY0MA2oLbgtsBaB/\n9ondIp/pKE4LcgtmY/upqwLoWNfD1nPL3aKIOxPieqrt9VBQr3RhtEhZRYdrPwT33OLOEai3M/ka\n9mP31I712whXm9Kyz71bkB8L2EAfoHv87Dw3OX9YUNvGEI2rrfEJ3aII2DWOPbrXPIhbbpTuXj/v\n3qLWUEhtNa2Vv14HBhcOQH0RS8uvF+h1egblWQ3UtXd7kYc1B/CvAA+abWkCPSLf66vt/bZWiPrf\nLKhhfrfgLn6ojX1bovgeqsVrt3UKReD1U1N7OtzGArXCr9X2a3VUABL1urGk4PbxE9BeoebrzPrV\nEu9BX6cnVJ61wFv7veZXY6firkUfRREtjLJ2iyP0ZAVogzoar60aDxlOwK5RmY+xe86eLTHZc+z7\nUMStrTtIVrMzlcAdM2GbgK7VReS3JVFpfP5wGwtkz1B9nleLkNTtn20EtXd6pKtEr2mteK1GWihy\nfn5cpn6tnv++OwprSaolsQZo+9uWyeCn4HdOClhmkACfxI9Y903rPPOyc7rs2Sy1Gsy6J9PFpIPE\nrbvIdMMJiFKypT6SnfxCHch9DCp2BDKfvV5ge3Hc2nZxoL6Xi4qQd4ApuLWTLkBbExNqoK/1IGt6\nDUcj9cYVIaYRrw++BdbIRBwjet5KHrAP2MCuDVzEjIFnMQzK7mRucYcT7nDEWTZ8TBh5BnU0iB5K\nxyzMONCEIy7F2jMCMBNhpgEzDWAa3PNOYAPteolsSzWAh52oiad1qI3FWIRF61+Q96wBeauetnH1\nTKJ45NdijS0E0fqnVjT2OeLYLS7u4+ntc6LGtMUdBvfb7GzDrRXAI0+hucEdbviME844yh7sdNCC\nSVjRdmqNicJHXV8+Coc+4QImgp7SMtG4MqDEuXkvgGv1gcCuZaVWHzVw+/qIDryZzXstBrvCag3c\nkTuKtI+eSStu/Vq9TYszb/RSLUDXnoeKuxfYvVzbumsNy/vFq0Uzxx55wpGvGPmKA19xFPvAFxz5\n4rZpXjF4jh1yE/2hUtYZ2CJ+0xXHXDBJLB9oxkQHXOmAi9hELMAeMBOje0wd/YbA9u4oj3vL34LZ\nK8g8oHczXg/UGqBr7752Ubzn91qN1Tg41hVK4qhV/BaQvdu+02pwPj0+m1tioHW3/LwoDsY4zxhx\nxWlO+6tPsxygMC8HKhyQwD5iwsAzqg1ixb3rj8mLMWDCQRaaJOlaOfgFl+GIC06gQTqTIfF4GmYA\nY7tstwBdA3Ot7Hu4NWNdNNqh2mX5flxtAb0CNWG1Em2FyS1A7+fUSk+8pBSoc+8ad+7k6NHrUSMY\nKv69XLvV+HqTrdTiGtFzZdw9MGMYJhx4wmG44IbvcBrOSfwezss+a7Z22prZFFNr5PI2EMvRChOI\n0i6wkSZMuGLCiDNNGGhODVvG1yONmIRzdwO69hykKaQauHvMFqBt2tTthZ2axF2Q/6EVgfptZ/4Z\n14rXnms1FYE/ChKylJir94C6F9xo2BFFgPKAtm4Vw4EC4Kr5HnnGyNckdvMZN3yHW36FG77DgKRU\nGzAnxRrmRRSvATnyr+SHKHFskBxpSLJXmxJnHij9lsbUAyYacaUJRPP9OXSLF7TysyWC10C9BWib\nJnZum0ACljPSbMK3OLPvKfZx7WdUnvVw4oglAusarvzcyw08eD2ocQd88AXgW58Hrl8Fjm8DH3sf\n+Oi7wHBTb3hb2QTqDa3mtgC3HJsnGVuXwH7BrwCW3dOcXspTYxr5Hq4dVAVhxkBpzzfTBMBMb4FA\nAupZQD3SAaMZa98b0B19PQBgugO+9QXg658HLlJ/H38f+IjUXw3MhBLUfrlotNsyArclBlZHYlUB\n67l1lc1v0jOI4g+h3prceL2HQyiof+3HgLtfAK5fA+ze5W/8DHD8FPDWO8BnfnQb1A/hJnDz1EA5\nTy3KsmxkrH3EBcf5shF3kI7oudEvEyAgZeOx/HilEVcccB0OhX0Y0qGKTEM5zw1CnudWLndfqeiX\nfgz44BeAs6u/f2Tq73t+NB5faydqwUzGz4MZ7vleGLwfR96iJwR2NOB4AEhbtMXw/W/VirkTUH8l\niGACLl8BPgCAu5hz93JtpcoYkESMzvPUMPPUmHE7v8INXuHEZxz5gpFF662r1DTOWgfi3dFzlI+o\n4wq4KRFjJJnnpgsmGjHTAFDi5klsl3ltmHlumHnu2ndaNN0JqDvrrzWWrnHnqthdlsEj43Q3PTHH\nfmJwtzhkD8i9SP7yC8KpG3T5WhLzPvlevbPYSptSBDhwnqMeYOaqRQl2wDUpyfgORz4nrTdPaZ7a\nLz+tcesHiuKFOwD4QIyB0jnmRxox0zmDmmhOqj2d4xa3vjtFe7OjdETp/s0vCKdu0OVraZj1sffW\nY+mIW0dMwEs2Nf8ojc9EXcAmohHA3wHw/zHzn9z3iSfk1AiiroGsBnTr/63PoxDfQprS2O1T763f\nt6s1t8TzRoXryrJjnq4y89SQeWqWeWq+4FDj2C27Ix1huqO8uYZPSPPZB7pmsTuJ70nRdsVxmeeG\njL2BNbeupaOWh9/orL9vfB74+Hvlzq0at/ag9W0pAvmTA3w7sl6O/WcAfBHAx7Y/VOtanxDc+tkW\naGvitzXXr/Z96/LVbe25TQOwK/sETlpvJOWYAviIBOg0P32RBSlXHFiXjM4xd65x6z2NrSWNrJ45\nT4ExXRPYhYsfccWZrmmem0RTTzIzTjtbv++szjvrT0VxO/uw1UZqAI9E89dIm8Amon8cwA8D+HcB\n/OvbUdqce38K3Pek+77eAuPx7b44bt4uNelAfb7cf7eVLo2S0x1bB0w48AUnnPMS0RPOZvnoNYvo\niWNLXLUxdQ+oW1WzBWixB6QxNkjmuSktXpnomrTkco9YSrBMiWHAgEM/h4vyc9tZf6e3Y823B3jL\n9AxfHh3okQY0ph6O/e8D+DcBfHx/IiJwA2XO74HQTY5RMVthPvZ+0n43xbkRePv9eN7bnAoUArwG\nbvKPnFZxKZfGWdZ/v8IN7rJSLc9Ti1ItFMWt+5G04qv8WT9K6de57DzPTXKKOQ15DlyPNJ4wpOkw\nOcfcqRzW5PUI6v7u94Ff+xmAN+rv0++vQRyJ37X24/Lb5OpRk99N+19qApuI/kUAX2Hmnyeiz/Z9\nmJx/DbhPIJ5vRReBzJqPvpumRC6RVlXo5lPAJ99FvuwqAnfrWz6dgVu3Qh54wkmBza/wAi9xizs0\n56lrnHkPqKO0tZ5dHojStQIDuXnu9CMIcs+XaMKvOOBCjU0qtfx4gH/63VQ/rzbq7xPvrvPSw6V7\nuXYN3A+ifWOoLY79BwH8CBH9MIBbAB8noh9n5n+lDPaTxv39AH4A/Vy5Jrrb33d0AL2VEr0z3AAv\n3kl+FzcPijE1io+9Axxv1mAeuP0t97wsmOPldyCNRUXLPcqa7xPSIhTl3Pk1XdHUoyirARxBGJOe\nkFodVGhTkfekDUja8CsOeQnsAFmdpifCsIncCCMhqBnAcAI++Q7wdQB3Xys5N43AydSfXa67t520\nwu0Gcku8jn77EoC/uxkr8Wq5WyUg0R8G8G94rTgRMfAXa29VnqOBacT+Kn5bmzmi03paJwWtDgBx\nK5du3k7i9yffNaDm8nuDSwsYNIhdZIdNqfDqvrkBM17gA7yFl3gLH+AtfIAXnJ5f4AO8wKuHab5b\nbchTDdw1icM+N+yXuF1ySDmXeIkXeIm3ZG/3kqbURGnxModPsLhzGAZwuQN+/QvAl/8m8FLr7w8C\nn/whgATUaan74vZmqri3wupCFrUjd9Ez1QJ6Y3+z9K+CmVc1tXceu9IsWqK49dPu0dpb79j3Hp7S\nMJk+6vEmzVP/jvdi8CuoI2AP6cNpzjaF1aWUpIAXLVdiZupO9kBzntbS+6yL3Vq4rjm1dT8WqGtl\nE/lvgbvwo5yPyIx0BfEgyVLOvWw4ZdnQvdyVQBno2e90A3z3e8B3vlcCzWJD6w3m2ac3ysejUYtT\nA2WF9oybSuoGNjP/LICf7Qgpdg/A/e8w4bbeCbx78e/FY7jnlhi9ArL6LWCnDOYEbn3OfpJYaaKF\nXzrA4JrF8NzoeQF4zqxK8Vv13it+16hHLI+eqZSm9ccDDhnEIy4YccEBx5xvkkJVRmRvF0sjEAIx\nGW7NoOwmYHZtx+ZRNeG2rWyJ1ltl4KmnM129EAE5qsg+CfuRVp71ABPOz7t7OXTw2z0Y+opqYN7a\n9aVc2toZ1IHBYqfPlH4D5jUnY8u1l7XgVqRf++HhgI7KqOO3RXKmwE+BLfd00zGB3ORX7/5kzqWT\n4jDPPEsFzJ6DA8u5atIwIsWYP8jC5+Oh4M6Z9s97KsWL47Vwa3qiJaU1AOtvEYi3xHb7CgPRdk37\nae9W8pVU66WrIPbGgdpwbCIGDbLWezAAptw8V6bk2FfhZtZcJYls8idNvwes9wG0p06pyGoT7HMh\niVgxXJRoqjUPjXBqHgjzPCxNRAyzb0sNRhABHai3iR5qDneiHyPu3OLMz8qxa+TlZKAEMYwfgnCt\nqBvg3qIN8bHKtVcAN6CW8Tdl7XhaF51WXIl7kLln21RJTxGdDbCngoupfZQVZyqgrrMSVPpjAHmL\nwmpYyRNQUfpIC8e2ebSieO7mOOj6ZFw9EMQ9CMcm1y6k/dkzvSNQ+/pG5bma11VGO39sidpeNO8X\nw4EnFcX9b70lUuPm0fuuwmCCRxRV0haoW6J5yL0ZGGYngs9pbpfmBHIB8gDjFn8PalWkKZc75jG2\nJq3s6R86InkMWnc55fmoVxwzlz7igjOVnHvCCGbCjDmDXLs9BgOs4EZStPEM4iGNtefg+xYjypkj\ngANtMPeK5CsMbgGyxq3vz7WfANgeWZ5j16gmgnuu70u10nHUJJytz6vbVnA4vl7G02Q4Nw2MYUgn\njIzDnHY5kbWVY8vJJgXIE7e+xUvcIm3NvKFXuMFZDiY8G2CzaWe88gvL4iko+KAF1qIWXPwuOOMk\ngrcK4FoaaTGtbulUGUZLSOQdHjDxgJnHbGNOklIxVTYQmBkYDCdvgTkSvWv8xLvDMm4BM+LIUaQ1\n06YnEMVrII8A7gGsZFWX+o5drS/vrnpFqkszNpgvtxatKp+zm4ySDGZMnTS+kxs/LuPIoWimpXvE\nJDO6r/LMrj7f4hVOOAOZ0y+Z8s82X0/JxdlXpwEyG7TYZ80xy1ozq1s4yDnlnIFOhVtX0V/z0Yly\nTQEdcGUGhkMag8slBSqWFzvHInCvMmbsWnvxOCuwxwEOW+K1tnF9nhvvbNMTi+JRImqitg/jM+07\niYhTM1arlno7ui3xvWgUbvrKuAfZmnjEBSc640gXOQI4rSA74FIA25vEsV9ljn2bTboE4IS7DGIq\n4KI8MQI1t/NoS9EIRc0mRIujaLdkubTX/cPYKUTK84wDpMzMBQSxIchZMbjgiJHSG2RjpAE0YBla\nC49I4I6kvg5qMdkVQ/Ult8VZar3Hll2nD5EoHsVnAR2AenUvknktKi//2xZVxlYLI3Dj6CFtfjjg\nmgBNZ9zgjBPd4YbuskhdAnkN7BPu0qaPbM7ZXXJsD3ALJZONnPd1pntE97iolBOuw5UpK1MHUGaU\ny/VAss9cdq9ZYE8BuM844Y5SqYw8g+iEJAkQrhiT5E0E0JCFOBJQV28iqVENoz225dwhuKs9Q8XA\nuev0jKJ49NyKo5drix23sNj2fpHAUBtzG5Gc3NTWoBwb6UigG5xxQ5brpvHy4Jqud+v2zJN0BMtz\navwWMpGxmczHklm2yvt5Vn5dxW0t8lz0S4xzM4X6i3Zk18ypb3DGGXe5RCbX/cmZMrjDjQjjSSnJ\njHRwIkZcMadTWIphU5pAC+9268l0i6E2GbD5gaMAjFL83gPwNj2jKF4bT3vS3xvj7IxtUYqA1tH3\ndHhbZVRTsOgSUZSgTieEyAITSruybukV3sJL3OKl2K8ygJfx9vI8ikiqGmM1h2xfZWSq02QRdJZ5\nbRXN9TlnrZL3Vu2wC6hFz3JEisoKNjUKcjvJt4ypU1ldZAeb5lZLZMoyzGJPGHHEBYOcaQ4pCd0y\nM2AGD0looWcfAAAgAElEQVTBNmdgM3igpX2gYbeoB9CrcBbcrRda68RrH6jTaxLFPQrtOwroFsd2\nnNszcVTsvs5uTQ7cJKBGBvWcp7Qoc+xrAjbucIuX+Ai+hY9Q2uqgzXYsmuvi1lsrveLtIM23mCYr\nIFQq0RZQC+fmMjuZKhycg8BWmucC1Mto3s45l1NVyT1ixlVUX6dgmcrSxQ2FW0vpA0rr5bUEJhpF\nlklvzMKZB4JcBMi5vthzbV+/LaqC1xkftuvFWmT3a8ivURS3IPULeGuAdnHlHpjLVrhV2N4/oo2e\nXQ7dXIniyrGTCH2HW7zCC7zCW/QBPoJv4i18sNzUISBemuzVgX1yzTw9r+bAHXwgQFaoLfu1EYK4\ndSKR/0mLWXdTqRjOjOz2QPbuSTZrRqCdQjMUz4Q5nVvOiy79jBMOdE0SkNwRxrlusPT/kRS2Rb6t\nRBhbtTHz4wrgrUbpI6t9sE2vQSseUcSJayAPfiZa+wFxeUQtNUpmJHQUNofj63TIgIrUup86Kb5u\nBegHfwVPFiSXJp6PH165pxLUHHBuyZM9iGG5R3udrR4q2rUCxZwFXvh5QNMyH80YMGDAiBmz3iCC\nARNZBZldYLqcZGrnvdOBiMc8YDnjiBNOOIlyciLpKGgEETARkgKtJ/P3Yaa1NqYA5+jlmvgNtO9V\n3qYnAnZmpUIU+G8lUDPux9mWHKARBLW2v+0hqpyepOWssTyWS0LtVNZYmMShFdC2+R54abbhvdcK\nblbwcjLgNHDJzzBidwXUAdfO5KvOZFk5HiGBO03xkYi6RmEFWRNCvBQ1zSAS5ZldLy9r6AdizDSL\nEkzSnjfL8FInhKx/OGWp6JC7QwbhjGPSVJAOao4aIwhjoWdd5b3V8bdAbcMX7hagW7+3PrhNTzTG\n7knAlhLNsl8PcHmftYU6xchWRSDw2yq3iFsjUhVF89Nr4TIvsZCzwkcB9oGndPzwPOdjiGk2QJ+5\nAHV2C3Azt4YBszF7DwKNisDecb0cIpHGvP5QiQRw+TwlgGOgPIwZBk5i85BAzYPXqSN3Jlo9y1Dn\nmEFtgX3AjcyGH3ERkCTJQS4GzOMpkzkPyJbxYeDcOT52cfsXWhz7/qAGngTYEbfmwB09e3+72ix6\nBsrrSs13ozKyNzu0AB713Joskzu/S0vXUtEK1CW4F1BfC2CrPcwMmhk0z2IL0NVfwTwDxLMAVv1S\nwhdAO4CborsXaRmQApxljlg/SsUxUTQooPWSPgINBB4YNAzggTEMc9JiD3PuHArtvnxXk32QYU4C\n98EMTlLmDphwJ4o0AHIv94hJOWFtbL0FaM9wa8xiBepWBD1mFfkmPYNW/D7+XjOOii2/64n0Njlb\nPavtJ/aUWQHuaLrJLzhZT2dlTi1gzmeD44rDPC0AnsTMi41JuPPMgILbuos2wcUzdeRxS44CFlBn\npZSI5HpiTDpzOCmtWNfQE4EHCTumbZfDmKaheCQwEqiZKJftKkHCvlUUvzFcWrsAXfQyQE88TTVy\nRToFdQXoVpnUOn/r12o74bh6r4k+tk1PNMZWippJj1+NY3vbuLlAXBlkT0fI7v0wWwklBB1nelDz\nCtxrXe81i94HnvLh/yNPBtRIQJ4YdF3cUDDPMn6eAZqRNkKoXwTojXZhRxvN5lOAGgbQEI7MhivL\nHLLOJQ8Aj5TAzRBQm26SFjF9+R6Z5JPM6S/LeoBybnwgFb8HGbOnFeiD1Fk1b1vtpCU5+zgKjyjy\nKDIfLnrPxlmnJ57uUuoRwa1fC9CezbouuKfQW+Xps+CTDQuA9OOyD4mxjLFZNNilyYozVkXa1dzo\nkYCOGQnACmgF9RXAlQXEkEP07LO4ezq0dbZCv7AJhWNsfWZw3psOsNgYAQycOPiBisMHk1M2bgyL\nNFSsGGZkjfsBEyZcMreG6VTHrHqzoE6bcoZCaeUo4htq19qJ/72rx2gpziIZfz+ogUcD9mzcrX6/\nNd62II1EcCueN6gXzFFZ9ZSfNmAJoLykBPdaHNeNDosG3Ivi6X6ucZ4WDjxxArUYugC48nIqZmTr\nibsBwHtE8Si7kSc7QIOQ+2IaAYycQGrc0MMgmYppXhF9Upwz8oGQmEm2XibQq52OT1q6Ht0VdxAY\nk9RIWvByzptE0yKVWqawbgc1fPUw1VXELRDXIkPD3aZn4Ngtbh25t36f3XNQOCtv2lcheTEHFwBO\n7UEBesWB9NSP5XAE3TOdzLJabLUwkicMPGGc0zSWtTGhADQuxr4gHZdtgSzu7F9phOzzeV+KQG3N\naLi0gnk0/nbsvxR67hySBnvCMBDGmTAPCdAzCCOniwYOWVV5QUmMC455jf0d9DaVtCnniAsmsBSJ\nyAZyIgsDS1uRJHXpuXZzkx6DDXebnuEww5af989dtjxHbiskRibojvOYMwD4CuR6N7VZAMIlB87L\nPHk5zkcXRJ5wzmvClx1ZqXGl5aCz3LM1G+03J5FagXkxRgF9Nu55ATFbt9qV/D0msIvbbr1oPgun\nFsmDBnFPSABvcjdIB5C0/NCz4jBjGOQUGlMXdghkZye0s01r0GW+m9P89pVlfM5D0sbrAIrTarai\nnQBxm7HXKa2oxjFqv29xnei5Ta9BeRa9a8G4xbG3zLzEpxW09aoZtuuCjmjuebmn2q/hXtwnnOUq\nnpe4wR2OOGM5I9zsWWIznWUUZZhQAtsZviJxZwPq7BaTQYwlj8XzQykPRxY3sICdhgXQ2Z5kzD1i\nrQfw6WJ5Z8SiiKTBrLwbsRy74JfYolBYHrGs2Z8xYuYBFxwxsZmE5GSYxyINTWxVuXZP+9ziyH6B\nSuRu0xOeUtoaT7d+8/H4MFsFp5oW85pOh+mC5lqHKs968Xw6b2xarqyFjIOxHOgfAftW7tjSvdQn\n2ZVlOXZaeLJMaeWxcwTsMwoObkE8B7bmqbqx6KHkgJzd8tsgCrNBgD1YRZpKJa14maVlGi48zBiY\nMIDc4h87I7EoL0dMOPIFV1mhpgt2Zww44Jpqi9PEWWoaEguPpWTj8ajcWol9YduXerjxVpiWu07P\nKIrfx9/2YFZL3lMwhHzOVa1M7XPm2oyB0/3UB55w5GUf9LKF0h4NfC2eT7gUovgpv6fj7XlZTTYn\nmyZp8H5c7UEt4rgF9cxiG7cHNpt8PhquLddGCfJhkCWlAm4dOyvAQ8ZjhTZJcN4WOyTpZhgU2J5r\nL0tPF3Fc1+snUKc7z1IiD5hwlgkznRqcoWI46k1q5R/0oHHABxgfXx89gyhuf2eUNYiGf42r9xaG\nvKfLTqNVaZXyWk7WuhRrkm/0fuoC0KX7iIsRxV9lYKd91OY8EJb13iqCX7kEtBfBLcBlHM0W2AxM\nYiuQC9u7H0CrvRRU2qwacrHJPStHzB1E+A1eNtgMsscas+ztmjHlKUW7lLfk2El7PmLGHcAkU+2M\nO7lP/AxOhzTwgJFHTKwFq0O5HmPaZbjfuhe0tfeiMNv0xBcGWPIg3xLHOXB3jrN1AlRbjo2qtQGE\nZXzNnEXxI19w4jvc8qIMS9NTlwXY7qYOPaPsRjqEhctPRhRnUZ7Nyxjbc2sLagNsbX+zMdMsQ3QD\naJu1xwK1kgVywbUBjGb6S1epzWQAbtJGNkKnhEvcegbNJBtiAF6J4lYct8BOSk5VjgHLIhZVl4GF\nU0s3QDj0Mc0VM7gPqB8C/G16xm2bRTUG5MEduaM1oLUCsYAWBwU9cRAlqSiOKXPsG9zhBb+UbZeG\nU1uA56WOd3mrpmrEF449ZVF8ATWnhSdXXo+tvSh+QW5HCm4F9JUXjl00BQ/yB5IHs1eOl+BMZrDA\n5Upk0PcYeqwzjWmdfN79CBRieHmkg7h50YrbpaZ6iMVyNiphwoArH9LVxTK9mTtG20ZWzUwrAa6d\nPRWw7fM2PbEorhQBfSs8N9w9vVvwfVXaDDqPyrnn1h1SIyelyxGXZKvRuWl3gL89skj912Nxu//I\niI2ycQOOY3NebYY8p5214Ve5OHIW2wB75vSaL4lWs2h1tUq1dwhYzlTzcSmHFntSkEPG2RNAV8ha\ncjFXlPPfsnRWpShGJHovANeU2uW9I6YcMnX2hAOPOPABRx6TAs2YAx8wMwtuSfKf3ltuG/Gl8xDA\n1txb4dr0DKL4Ho14r7geFVawtdOfYkrpXb/fWbdHDqIwu8ErnPgORz4LkKccrsgix5oBLM0oN7Ll\nEARextez2bgh67zV5KkrHfbx4lbMXxXQKNeq+Opv9fURaLeaDptwjtmCsCx+I/eBLL5LvtJmFnme\nkOe+MWHZCDOT2Zo6i4Rvy5SLMtfcFx2Nw8aAJJHpDrsjLqINT+PwK6dpsZkH8Jzmt5nNPLeNM1fO\nY3FmNPx6GeOzXvGj/rZZAPWm5P13aMZZ3tdthPmAacj0C6dVX4E58BU38xknvsOJZQ5alC1kP90o\nhqWhrY1ybLvNslj7LSaao5ZNXAuosXZPLnk195LOuMR92BbXrsVZuCUC1WUOCmoDbgW4XSa77DVn\nqf1yaisGN/LHyufFf4DOelwErKmdDMy48jK3vZrnLpoaF3E+zm6uLVD3gfs1iuJ7lGc+TG/BAMUc\nNkjwbuapdY12dsvyQz6LGJ5+S1MjcxE1FZ9aGtnCMdYrogoRXBuHcu6J12DWOWpexs/KqVV6t2tb\nlFvWSnHp4kq3fW7VRBQmFOfZjLdd/0qSH5qBQQzpyjkV00UMT5IKSWeIQux2+8JiEKvbmCQ1LR05\nyzRXahczLlY8n4+5rpjtPHfjA7vBW2vXaLjbtAlsIrpFuvD+BsAJwE8y859rv1UDsH+O3D2dgP6+\nUWh6CIOd6prTqR1aqac8jj7n8bSCXZVjI09JDDTR63litUP4S5CXY0Hl2IOumbZiuKz55snMTTtj\nQT1zya09sCOqgTvi2L75er/uj5mXBk6d1sAohh+DrE6jidP91yqKg3NHkYbj62WkKMraVBSvO54B\nLIuQrjhyqstiURKfcOEDLnySb58w87DaDpsd3RpxuHD+ndpv3r1Nm8Bm5ldE9EeY+QMiOgD4OSL6\nQ8z8c2VCPbUAugVuH48vlC3tOIQl8FIBslCFaBlb2+msGz7jxGcc+bwW0WWMDXC1bIsD+U1lxAcx\nJM4DxiKK28FywLF1jjpzZgtoB+4thVg0PrbUAvNW0yo6Cgau8jGWDzOnzV65Gq0IPiz5x8wgpizN\ngNPPc+4k16J4nB8FH4yRKTGeJJ3Kwa+YeMSBJwx8Ss1nBpgJEw+YeMRqntu2t01gay9WA/MW+Htq\nIFGXKM7MH4jzhKSv/A0Xwr2xpRTzYXrH2RpHi2sb0BvxO0cjY+xlnjoB+oZf4VbG1XmTBps5Zxjl\nmdsA4I/xWTg2sgi+FsUhHIBX40srjs8G0AWwsQZ0tFozGvParrUG8qgpRYD3biBJEgexR/lxlA+N\n+oK0b7KKQgW3ATwXp8AwBpCRgLgoZ5uSUHGWf0vKszTWT0o0VZalxSpTPnIqrUgbcVWFJ3PWli+Z\n3wJvDbQ1cLdA/YjAJqIBwP8G4J8C8B8z8xfbb7QADPdcA3eLezsAh4VRT44dYx/5ihs+45bvcMsv\nccNn6DRYnhLLNm+Wc6mZXcaAK02ugntOcRcytTRyXf+t3NqOrZ2urfCrKbRqCq+ItvhFjX/o8xUC\nYvNDjkccWWk2peWmbDJFM4NnJK6dFVMp5dElCTBlHeXPKj5Jlu2NzGk4IFNZOq01Kqh5wDSPmPgg\nCtQgE5uA9KDvfadVC9vUy7FnAL+PiD4B4L8los8y8+dMiOCtXhBHwI/iLNmunoxFUmAkf8vmXwV+\naUae83z0UURxnbM+ibtZvhVaAF0+ZxBn91ycKKpHG62murg0eYyNtYmajRVNGW0we67e0+3CPdsi\nmk14ddtOZZa8zyzjbUae2vMmKdCktlXCsZ0mWVDbvPB6elLszOH9yh0xMxOubLb48AFHVabhIEc2\nMFgSqbu7OawJ33ha4EfwW+Tepl1acWb+BhH9FIB/FsDnll/+RxPqd4vxzSkCL0y4VoLL35NIC6R1\nQwtsU4VZcHtgjxgx5+WhJyzTWZQnjFEvxz1lzTGnJKAQLRcO4uIVk5mV+9RWaflvFt/f8AuSsfLT\nd21HAue2VOsj2XiuhqkSm4q/meuSAbKNrEdp4BMR1C8xzHTYFRPOeWkRANmjxwJfztBWmJcfaT17\n/57G9UsAfrmS0YV6tOKfBnBl5q8T0QsAfxTAny9D/RH3lgdrDeAc2Gi8A4HtjBGczVC4gTWwF/cA\nNmu4z7IN081Trz/b21GmNBZhWfxK7rCI+siNyXb0bhgPj/0oSVug7QG2BatvbnPwu83lFq5WTTmK\nxJaH6dTs2DavdhNHMZ6O1Dm1hLhvKqVlqbqC/IJjZiWJqaTFwcAEEgPh4iTl08hYWLK1xES2Mk6l\nnwky2sexfyeA/0zG2QOA/5yZ/4ft1zyA1R39FqlwomabMq3L9tNC/xlHzHLe1dwBbBTbMNNZ3jJP\nraqnrc6zA+QFtzaNU2/rKFp3o96jYX2tdPyzNbq0p8dWBZy3YZ4p8G8xz1W6TeD22o405aVzxwuo\nEfQmjYqp9S6RH8u93XxNoGZzCipP+SycCYSLbCRd1sP5falbAG6FqSV4m3qmu/4vAL9/I5R73uLQ\nNf9awpcwyw6dtKD/hCuOYidwM0pAlwAn6E0ScoggLnLQXTlPHXa6PbZ/gQ13xgL0UjFnzIxq3LUm\nYEspMkPFD5XfWjWipau2B3fOY0AV3DY7tjRjKVNf0gEkUDNAIgm1FAhRHdZ+k+c0cyHrzPmS/ZYr\nAEdcxGgb09PMy/bdC+BamK3GWKcnXivuAb4Fbks17s0g5GvZZDfVBSdchRPPWJ2gVwCbKscdLZvu\naxUeJqdB+cZLlnQ7MHu7xU0yCHq+ixjUNfeWKG7J9j/AGtxdSazlMwJ1DlNOe5EBtd5R0E2t+pX4\n7VHGSTRPbW7CiAMOGHEAYQRwBGPGBAbhAJK9Y/FH94C8lcBteoZtmzWARuBuxbn4Lxvp9RCEC25x\nxo1cxxYDO7kJJHPS9iZLXhagREm4d8eZOgpVmFEA6pVxClUP6FazAOpj51LBWOfiCmyvYbeTi+qe\nzbe2iqPZrBv5z2U1S1r0WQs2N7FgmsvXn/er4CfPczOy+G2vgTjjAMIRjCNmsJysRjLQa/RUfaUR\nJyp8rtMjAfuvAvhBAJ8BcDQf71Gg2ecakXHpsTd6zthyuskL3OGUZ1Ar599K4S/isZn/bGnFt1pt\n5KcN0MWXxfJWG2AUyqUoaTXy4PaA/tYR+PxngL/9g8C3PgF89BvAP/dF4LO/Anz8sp6uirBgwW3x\nFTFONmGyMR1XtQNlw7WhYvkytCGYsbbNeES1Djmob5KtdIsyttwBMOAAYAJjlrUDhBGDcPAwt+5j\nlcxWG59/d5seCdg/D+D/APBxAL8LwJ9ewCKJSUVia8DWAjm7Lijqxah64+JyhvdZjC7X0Lns2bkH\nrKins/ThmmVc+ZFNToO6yytg+2JbUY1bW3D/O38a+K1/DJg/jryC5AMAf+33Aj/1m8An/gHw7/3E\nmjNbbk0oQU3Bt7spArdr31ZDniWdvUqzltv6sWmrfmoiB2NMIFxBuGIQezFTlnG0W6ytPIjMnL8S\nJ7pM0BUxPeIYewbwdSQwnzHK6VTlPVY2UR7MrWZpgT1nsXs5/XNCPu5mL/UiaU8DceQbf9f3Oihi\nVr7rtOaDo4D6O4LIxuT/DSSO/pFLubvdjqn1ORqb22/bdN43X/ctmyr11N3mN5eZmTQdZvn31QDb\ng5rRB+6Iq0c28FuVFD6B8uy3QPhljPjd7lzuGfmaFQD1Pj7iNcj2CMYpK8z09E+dWayx2c7WUQu6\nE8iWonHf3obe+w1lZDXz+c8Ip27Q/HHgc58B/uQvxmDW77UMsK7Vx6AC8Huo9dI9wJ1apZ5brmvP\nZqQFU6OUUwvUdiKxBWygDfRnBfYM4EsY8U8UR/Omnk1WeN2TWwOEESyHwKdproVjL2L/veg+vXcP\nt3f19JSNfasU//YPwi3gDmgE/tYPAj/yi6WCrMWpe8H9GHkv4lyJ5Pek3Z11Oi11lOlVHWkPuCJd\nawC0Qb3FrXu4dpueZLqL8Q0jplzygX4nnCXbUdX3gXsA8mKUdNjgcgd1CewW+31seJXkcxTSvnoK\n4/eAVr9I5hmQFGU99M1PxGK4lZ+ib9n0+bw/lZTSXZ094+suSuvLSJRrurhZQT1Cj0O87xi7BuoP\nAbAHfCyfEnnCxZzYecah2H9U69/rwCbocoDS7Btjd7aG+4hwDSoaY+23B1KLe370G0lRtkUf+8by\nnuXUkcKs9j2bnkelh0g/Hh+RdNX8Lst3ZxF8CEknTrKsmUQsVzBHoN4zvv5QAXvAAb8Ho9Fc63nc\n6ehee9RdD6hTnPobCde2nGiA5SStkZgF9Aa4ox5+R9muGvhWQ3pIR1H5zZYRIU1p/bXfi7Y4PgHv\nf3E9T+214+Tc9ptRp7KXVrXoy+qxeoxucJc/DAJa5dwKaB1v9wEbKLXgj8e1nwDYH8MRv0tWhi1n\nct/iFV7gJY5ZQb+fW5ece4E4Zb9gKqtJlRbyQFDvpcfkaq2x7md/JU1phVpxoeE3gT/6KzGnjsBt\nO5DHLpoecEdSwm5qiea8DqBbhXU2BrJMZQG0KnLvqxF/ONfei4SNqD4B4HeCcFz9ujqLCuhwx5le\nKrOcJ38S1H0b0UcvaZ56+EdYH442Jf9P/oO0SOUN9ZE/vLIu6/eAuYWH1yKK/zPQlWcHADPuMGHE\nBUczxwwAhMMjcGwvgpOx9/XeldBeE9QI+hj02JJlrd+fAfzYTwDfPAL/02eS9vubsvLsD30R+Od/\nBfjEpc5bZmdHvz0mreIL6iUa4eymqPCLei9lB3uQdFqhYUXxrfnrx+DY2/RIwP6R7GK5VmWCXjZ7\nylyVQXlh/X5gA8CARXk2m301KhptZTxC7I6gNQ1RQKuhW6vxdCSn9o2tJhHtof7oBfgTv5jmqn3H\nGIHaN08P5hqf0Wzdp+PqAfWjULe2rwS3zlpP0qLLc909kCNQt8bWHxpgL6R7UhXYC6jTlP6YM+RZ\nYawsW4+tgSMmjOA83QXoDRE+7hp1towt7dROytVSGdY/Rnvd4gN2Ua2do7aTkDVQt5qeNTYfjymN\nAGVkuzl1jZcEccffTeBmmVzVZaUX2KWlAxgRl36oVjyy6/Qk010K7OUwvwXsD5/HZkyyn4bldG09\nIqk8uaLGWp9QphbqKv4HiPgRl2yBbTZhCXVw63MN1D1NFChrwKbxMUu+1UmGtMXxuxOnI+lBbkQf\ncJZZ7LM8LwtUekutBW5U/Nr0pBw7PS/ce1lSGoFZ7bYZMWPCFbNcQZn2VwMjZulGdvfjZRJ6/Wu/\n+SRoGF6sJ2ng6GsuXrsN1MFdE7t7miQ5t6bvoXkvOrXHKsid8ZRSabot/Q4HnHHAOQ82H7pA5f7c\nGnhCjp2SMCBdUn7Iagbl4At5zurBXJ71kW5PvCCNt3Ulmh4tVzt8d4fo/chj4agangLcLdF4D7jt\nONuDurdJRqB+rHzei+4zpGq8w2Z8vewtPOAVjriThdN9yjMt9Zrsg8rzNj0Jx4aMP5bGsux5bnNr\ntSMzSIJnAHeZUx/AmDBjkjMsdlOvWFYDd6circll7tn6W9Vv/Vqg9uBWvy1xvMZfbLruy6mLfDz2\nyKmn7ja/aU8FOOCMI+5wwivc4BWOuOZSfIhmHA17m56AYy8zy9wsxYhTW38PapI4Z1knvjYX6Hnj\nI+qHLAxFR5N1InAXu20pV5qVH/+oVyuTuH2oUHO+gyIAb4E60nbURPBe0RzYh8elDlxixI/NM0f1\n0vrang5ZvwFpGXKh99Kek7+5LV12QpwM5z6JXNoC9ha3joAcdd11eqIzz1qI0OeaOB4BeyG7EzaN\ntHUZfuoAptxsPbDVprw/nEj3iicblPxWyYmy0Mqe9bMIFlTnC0Cl7dSElHxxPLBcGo9S1I3IisO2\nCdXAzCinu2rArc1r+2FArfaL7Lk8rgMsZVUAHNIF09JRhp2wp14eQwDTcmTmMletrYRwxkEPr84b\nh9M+rwH1zHDgbskytUHNs4ritcQ1usXqs/dfg3xRXBxxhoI6TYBdQmCXZ56NmDCS2StOE0CQeXAu\nk+DtrV7f54GQ73HTBkkKbmPnuN1OiwyCSqlE4rf1V9uOoaOmp+CGCd8CdmtUWCmJepft818DOJAB\nnQvFRLJKg6+/WoJQuhmEmQZMJC2EypMFzhgF3AcBtp444E8prdWWT5CvNR9mvxz0TBzb+tVAXUOP\n9Uv+cm05rtk37a2ZcTDHD0fcOincDnLyygFXHHMDYQw09/U7veUroNYHbZALqHkN7kqx9DAlpZqq\nxYrnNs6Ij9TG51t+XVTLp3MrN7blpj0dI3FWYOcQZqt+Jf5JprPyiQK0HMilRw+fMcgZ42m5lN7Q\ntkTkgRrVYo171zh11J2v6Yk5tv+tF9Rhvw7AiuLAwqnTzPYk89kxsJN7AHDEiBMuuWFoV8E+SVHS\nt2zzQm6YkKpgqb4MZsrgLjhWYw1PC+ARuNTPA9oIBSvj10TZ0WBk10aInsKa9R7u5EWWQFo3GcQZ\n1NrZNz4auaPfSPNFYBKpkOzJejqWHvOClIvMYS8Lpy1Io24zArylCOARuNv0RBzbfzwq0S1k1EGu\nK3+AZX3uBF3MCtRuAQEGjEBx3goBGEhjCNhjT1IDKkRhuZubiXJDJHvMph3XO+PH1rUk1ETyuWLr\n+9EaPw/QrcWPtQ7F0irdtqgln3F/vhQAUzC2Dr8SUFSHNT/h2JNw6Qsd5RzcG9zRKd2VjfJ6nyS5\n2EjvK4q3bGBdsjE9iVZ8+7mGlFWtrvwZVHAMvVJ52RBCzqcEeDr1gpaQlPj/RNIt1PqjDjAXYdk+\nLJwlN1IkQK+4tQN0joXqpaZkAWZBrL9FuKk1vwi0LT/v9sXhbQX0KmAGMcBSdSyg1gisEk2/uyqQ\nqJBMY+MAACAASURBVL6igguaXB5jyzz1WTYfp+ms22KzptrLsKQmMtfadksU9+CO4o3pCUXxvX4e\n8L4gFltvRSbpNdvHL5TmoCI3lpnIBOphzR59kqMkBtmyIriK3+qvv1lAs3BtCrKeuVlQWjXyTcD3\n91slHQEZDXf0m6dausnltUwMLYGkPLNobr+11enWOhD/m+1UkJRnVzrkM4DSiQK3hbS43LXpD1io\nteEI8Jai2vI1upXhZx9j2+ctlLQKJk1HcAheD+5SFJ91soyuyehkhTyTTIOl9mPB5ua5HbHLqz5n\nbbia3HlwbrCUrpzIxk8FEQGDGAs2zal3+1KP0hwJd74Z9YA1+p4bJq+N5IWciaqRB8qitw5lmLQN\nIH9pnU4yOhSTwFyGS4+pwyWVDC50TMbOVxu/pP2eUA5wvNYC7pmdX0Q10Ttyt+lDMsaOCsI/W1OZ\nEwnfcUkgJMWI9sbDESPNohFnTDRgkPnsQfwHWua5lwu51p8A4BqbdgUG2Aq/PLZm0EDSB3G+24A4\nfW4YkkjKBBwU6FiS4c9L8NTqVr1t3RGPaAmNluxkY2FosQcCxgGgIeURY3IvOs8USMHGqgmHDse0\nc7dlDGOXlDm7drQ0YKYBMxHmIbkZyT7TCa8ojafPdMR1GDHNQzEckIZkSq6nLfqVBD6FvlRr4N6m\n1ySKb4G6A6hVwLuBavDpmXSeMs1DqiJrJsJEhzTHbQ3SPPdg73A1UZcNy49Fl6GDLngYMjIJNAA8\nJI6tFxXTWIK6yA4FYK7Uda1vr4Fb3VaQ3AJ1lIQQzHDAHoRzO3Brr8A572SUZiQLRhTUCvSlDlbp\niTpfIkwkJ4PbuWoxZzriTDc40wlXkgUoNGImJwXYUssgtwt0ffuM5BtN7RbY+7k18KEVxVvgduCt\nhc/KqhLcqZEkjn2hA2gQ1keUgX3AZRHV8/ypcu04l1YEtKbkMLqiSXgNzeCBgEG4tuXYtsFLtrPI\nyiZbXn7uJA9mb/voWqD2ZIF8QALwASW3VlAP0pGtZycFzAOJOF6W4WJbcNdySq5aUid7JZmnpmOa\npyaZ3lKxO3PsA+Z5ANMQ4IyQ7/Zl3y3WmM9s0lfjxhGgX7soXgOzfa6B2oeJwvnffcE5FldopBXA\nY+LUEE5Ao1TsFcdhxIkuaSxFBCLGSOU8typx1iugYBqbBbcVG9P7s4r9CmoVSUXqH3SzmgG1fm7q\nUJJ6757xtwew5zM1nmOfi5X6DtQW2GTBbQDOI5XghoriFtx2aLPuSHOKaF0OM5aVZRnIMECmBewX\nmfJSjl0wClXPk4C7AG6tffaOn59YFCei7wHw4wC+U2L9T5j5P3ShojcbfjVQ18BcC9NjckayVxpf\njZgI0niWMfeV0h3IMw24oSR+D5R2jx1oKJLPRXb8mLoEOAOFKJ4WqEhPL+I4j4TBieXa6HV1mkin\nkgBXtLz2rrltLfhn31dEgK51DoQSxAcsIvhBbAU0RcaI4ZxZewvcthNtJKowZihGB5zphPNwwlnG\n1Vc6SPtIJo2xU5somHIuVNr6oDMIIonw8rSi+AXAn2Xm/52IPgrgfyWin2bmL5UZqVGr+fQA/B7G\nT/oWoIb0vKkHnjBgoDEDeKQr5iGFT0ozxkiTKFdcUitZXZQ4JWdZ7HlRCgmoFchZFJcJegU1y5h0\n9t/mRQKsabGtFn2LfLt1WVv5RQBXbn2wbloAPmg16QyAdmR+U95gFGcZ1MtmDFu+KZ9OcUYVRRrp\n6sXEsXUByiu6xR3dYKLDolgj2QBCRnk2wPTqUvBkK78H2LWuUSPZ0wGsaRPYzPxlAF8W9zeJ6EtI\nd+V+qf7WVrNogbr2+30BTi46yj0vEWEakkiepPYE4hSMMYjybMK4VKzLQizulpzEgnsmAtEAlvE1\nE6c2YsfY08LBYBVoAm7PPnU9egTmyB1R1My2wtZqsRgy08KpM8c2apKiQzNuHmVsbcCdOldfpras\ny5SsFGhikihOuIqe5QLl1rd4SS8wSY+zdL6pR82bT5hND1j05vEHdwE7UlM+sVaciL4X6azh/8X9\n4kPWYgjsFsgfC+j+Va20cRXlTMMyv02TmOVywSvpYfAlWW4x5P1nM5Jgf8AkK9kPMpcOUForPhAw\nElhBfQDIb+VVpMwAHYBxFrDOggPRxbDa8m8PsPdQtVZI85+mskaxB2MPljMfsCjOlMUrm8+dmoLa\nziwMhc7Cj6ttuDSnMSDfqk4HGTsnTn0lM0dNR1xEDK82TTUZf1Fb32qP9+HIT6QVFzH8JwD8GWb+\npvu19lbl+TFB2wCzAXK46JpQKtjlnTzmxogzjhh0NbD09gcc8w2iCfSX5RkXWYuuEev5rHrqW9Ks\n6jQXjwSeCDgQeJZK9xudPTF0ChyDhGUGDuLnObh9fixg5+IyRa3PozE6ni7K+2jMSYzx44OUh3Dt\neRiWeWfZHlkCvAS5DrGWrb2nbC75wqkb8UtKsjydtWo7iPVhnoOy9etlNnvEbO0Q+qgL2ER0BPBf\nA/gvmPmvrEP8tHH/kwC+zyXIJ/A5QB0BfB2s+H3Q8RcJsA8YacaZZtnHkRpMmga7LACni1znm9Yp\nLWvRdS/aLDz7ihmj9DMMogE0MPhACXx6gYrddQEBq0lrXrzCwqU5ce18P5i2tQrAH4MiQKtH1n6r\n8o9EC64AF0CTAzSOyFybszg+JHCb8W5p1ppy5dS6YuycL5pKl03lKyIpXcZsL2LuaicO09D3sma8\nF9h26stPg+kHPPB/EcDf26yfHq04AfiLAL7IzP9BHOqP+bcazxGQ/fMDgW971GKqyxnPrYVYGtCV\n7LnolEW6dC7lEfbub925e8LZfG7OQuARh7wZn8AYiBPHVs3vwSQiAGDur8WRp8TMSTtkwQ0DbCOH\nP5oo7ovXVJOdo86AtssPjhVQG3CzcuxRAF0F9aJMKxcCjQLuA85IK8le0Qu8wi3OdJLloTq9dRSO\nrfOKJqO2rfitcVoQCuiiknqAvVcxxgDeEaP002HIHo79PoB/GcD/SUQ/L35/jpn/+hKkliAKfn8E\n4HYVGspWVwu28nfz3FJ59hKEA452C740n+S+CNtNRz9McqnBFVdcZEP+ACLGPDCGgcGjbiQwVJmj\n0jUQBajFZC5uRHN937ofjWp9KJVTWRbceWWdF8Uttz4iz5ctCrTEudOOK937XB5XVE4zljfRnPMm\njlu8xAs59eQoWzIPuV641k6qvEhAXcxj+4BbEQL9ANeeY7sie7TiP4dlnVwn1QDt7fsAfCuM/E7G\nthXREWUxzw0zzw1ZblqerVG4j7jDYEZ4qWmlRqT3MxIxZmLMw4yBCawLwV2/VBSdNQJsTAnINCeR\nXJVu7LRl+fkpgO2L3YDYzk1bYFc5tXJrs3Qtj7MLrj0G015DwbEnqMb7IMC+LYBdHHkk89UJ3Ijb\nRcGxaekt0xhtwVx+aDCcFRYikXwv6Et6oiWlrUyou9kVNtyt5aT6eyOaWqVlUCcFWbqCQOa5MaZm\nRNJkSAE9FcAeccUJYx5T65kbF5xlD5kqfhgzzWn8iDlNqUxAIQZWipkIy8aiEaApgZoF5DrmtsYD\n/cHkgOzLtljzrivKdFWdAvuAhVsrp1YzUBpji1Y8bdKgAtRegWbH2cu9Woe8nzodD3yLD+gtTJB5\naukE0vRnPJ25yU/Uv1itVFt5ZteRtzTjA5qa0w5wP8MmkAjQPe6oUKIOpNYr2mBUj2712WWeG8Jd\n8zZOaToK4pEWYCu4k6iu4+oLbnDOJ1kq156Jk1QwzMtUDVECZSQb+bRPSOK4njIxJxBkwPPasIzD\nH4UUwK2qGktAF2vBIzAfjG2muZa14vE4O1pmaodNC7Bv8ZISx54onVab56kl4XmeugfQGi7USG61\nSe/v3T3ieJue8ASVKOH3BXX0e4tzb0RV+1xOrqjMwpVLnBvTzKnXH0EyP51QuZyRZY9NPBQNcqIB\nA4bEjTBg4gEDDZij+SppQ9B92wrgaXGTcUfAHlrAbrWTWvuSFXEtYCMANut89ZHAOm9t56/HtGjI\nbqvV+eWLHFOkHeWlKOejdKX+rO+TgPuIC59kP/W4zuPe5pjDtQDYamgtYANtjl775kJPwLFbgFZ3\nlOHWbxGQG0mJPrcF8OjzESkHlADa288g6FUE5d3Ji7JHR+dEvJyz4Opp2e0FYLQroHXZKS8LWCYA\nMydQz+KODv8G8uKVCNytZhJifjB2UKbpSCPK2y/VDV1FdxBQy2JylnWmGvYiGzPOlOaZk31K+6Nx\nyisGFMjLCoLk/gBv4SVu8Qo30skeMndv6hr2ALoVtoi/9qIdV/eAtn98DTz5GLsF8lYptbjxBsCj\naGu/tfqUKGt+Xlh+VFDPhbGKHnMuNR0MqFl2mHH5wYFBYyn+pzXVDMrA5qxEKxRnc+o0JGEpxgao\n7RRYkNX1D0oCbh6CZ5m4zmeW6UESukRUJrr5ALPyTtyUgH1HJ9zRTTJyiKDOPxdabQG4db9Emta6\n4xvpCI648gHlud+OouapdqtttNpbGFAB7bXoQAn2iAotXZOemGPb51aJeL8ozIaJCnVvj7uVFSU5\nKpOBpNHmBObljpFSe7ssbRStONIVw+U3Uh54pgReNqPGwZiRE2dmAs0MVjv7WS146hxUoUZAs910\nN5tBwghOkptyNvQ4I+TxMfIWzOxnp7Osmwjn4Sgngt7gJb3AK1lYoocJWp3FspLAzEvwDV7iBc64\nwZlPuKjikuNNIZm2OnvFneq2aszDFmguUXLuHo4dcenXqjzbAvoDQezfrYHaR936VJTUiAquHW1I\nWCt4lvXKo9w8ar5l6okBzAMtKiFKC1loZtAwJy4+i2ER51nBjYRoAXUSBAgs4FYFGtmWHbTyasM3\nZVPcpQUgHzaofoRl84Rx511adm7auJMGnHCmdMndHaWTQV/Si2SQjAd1aUZZLppWmGknoNrzVQZr\nnLYFcOvn21ARBwXgjgDdaoD7RHClZwC2f34IqAMR3XG8Zr+xNVMWJc+SF8WlJdtpFg9uBbUVx/On\nqIxa4xkGEe5pxjAwhnkGzYSBZwH02hTTWcKtE8ATuAFgJfW7fIVUafT24P4F0AbIMEB2ZiYBMpnz\nxiiBmmnIy0Bf0W06GZTewgd4gQ/wFj7AW+a+LDvRmI+lzEqzs4rtbEVxU0a1/PY0Qzd6Kt4FkPfS\nEkxPGL0Q2Z72g/uZzjyrAX0LWR6RKMPadYy1z0XRI/BrARpYDTpZRXF4bu3H2BbUCehFtCTxSVwz\nJQAPpItXOHFxnsGcVrQRy50T0m4GBbKmbk5xJW6OcJkprfK0QWSCm3Ir7tLSZ4na7spawC5lJLMB\nOoec14LTsFxLSzd4iVsB9UfwTXwE38JHRLQuVxKU9iKWF9OMObFBnUbtwbq3GECE3RW4gwIM7ZoI\n3k/PxLEj/6i0aixVw0eKswqoWxVRiWaz4wRcb2/BPSAtDrWLJqxGfBlna/wFoyXDsUliYbXt2FC4\nOWf+gxmMwW/8VxFQbgQkcCGKF02ng2MXEryUVQIqFq7t8jFbNwZzGKHRP9BQ2hjy3PMrswz0W/QW\nvoWP4Jv4aF45pkBe3PosnSnroqD0m56Rtsmt/XML5CGgrTtqTH51mbXZ2fejZ75ts8Uih8pvERIJ\nq/PMfPRbn4P7rZV0jmzh2mKUGzDLBpK8glyna25wh0vRFA80yTq1hdfMGJbugKz6bcCcu4lF8B8W\n2GQIUQawaMh50ZSTiuU+by1y4LaX4WXxG5Zjlym0Sz3tJg0/Haj+ed4Zp7zFclnrfZtF7gkjJh4d\nsEfhztK58oCZRzAPMnTqyHfUybcA3QQ3VbSS2q49wO0Y/P7g/hBe8WPDtLpLG7UAfKtHbfW0nkIw\nL/ayc4qyZpw5LTTRky8vwnVGviKdcJo4uupxdUmFurXRpyWss0DZNv+DmSGPRvbmPBFpSHmqTNJO\nRntG+m8vx8ZS1mye7QSdT1Wkf6gZWzJ2DL1w5IPcn6XDG3mX1S47XWYzUrESV2usbfIeAtf7r9qU\nFE4Wxe2L6vaq9YhTRxx8uzd+RlG8F8QbYI5KdU9fcR9yDYFF40yzKIOYMPCAmYGBl73Ad7jBSJMK\ngJg57Ro7FTwp7RWbZR/3ARcnzOsGxRHxqmgUsKEMrwRqEMwimEVDjhxmu1w8o8macCyS5gLqtd7B\nzxzYXPlcJmAfirnqBdxjwaVV5J55zMBmDMKd5TlPTWoHHGQq8ouol5OvMGgftgAdVYZn+dsN+UMi\nikcmGm+LuF4cokDtwsaGu0VR7y6gBsvlfjqPraDmIXFrvmLECYNcNsAA9BwVubcRJ5wx4c7AMoWz\nO4yt6B2vjC4BvYz6UYJazulqdbm1Ilj55fLz4IakQFNhz/1eKxk5BDg5QB+M6H2AXQ8wC8CVU6vo\nzbPj2loB1vgM7hiSNHmN2hnchHK/NlcC9jbKPnqNoviWicICC7hRhtt63Selxcmjyl8ZAs+cGjcP\nAM+YGQAPmFjPq05TUUmRtOw2uuIVLjjgVhQ9Csl1s/cjU8ulPXQcqHMRcdHhtzi09W6184IJSYQl\nPiy4bZezpHqGB/hi6xr7iwH1AnAdRw8Z1AryBeBOFJ+BUl0f1HGrQFptqNZ0Id/UwifGWpGmzCpS\npLXoQy+K3wfgJj4CirF1LfqWf4sqvbmu7CKQWeU1CNCRTl6Rxdks529d+YAzHXEjWzhvrZY2J6kE\ndOyOgbwGtZO3M9NoNAr7U6N88lco9vfgjtzRPuoF2F4MP2QxfBHH7Zh6NO4hnR1ngA3oONuIFr0A\nb/ElVOximCweqwPJfeAeYPv36/QaRfGa3eLWGwXQs+VubxYiUdwAWjdj6VTSLBybhiPAAmqMuFA6\nHOmMSxYlFdgEFl3uVXaBlyNSC3o/jkb0TFg3VuJm9ptlUC0W9wKV/rbbWVJqRXOVN8pxeKk4s1x7\nmVPI43IL8NlwawF3AXAdRt2HtsBcGEK5nTP6ZmucbX/38fSB+xmu+Ont8ux0V2187dDpC3OrN0Xw\ne1RG0x3wW18AfuPzwPmrwO3bwHe9D3z6XWC4SRswCFnjWizfnBnX4QCeCdMwYMQhQzkdwHDNY0nl\ntMtEzxrY2+L20vMUWbHto5LNh1DRbwSRlwCPObkV05cBBqHcgmk59jJnvdKuz1KS84D5cga+/EXg\nV/8W8K1fB158Gvju94HveBeg2z5uHTXhXr5D+sAOn7VIPcjZ+cel26JnvG2z2r2536Lwe76z83Ul\nLatf+jHgg18Azl9Dca/lr/0McPMp4JPvAL//30o9vwE4zPLfeSaARvBMmIkx0IjrkBacHAbdub2M\nFIfMv1MkfiztZ4KXrLCp4w6O/FQUtLOCo0tD12C1GXjN9StZ533Oe6pVNJfucVaN96L5Vvf8s38B\n+NrfAz74GtKxNEK/+jng9lPAJ94BfuhHS1C3lGoRbTVT+yPBcW8NHAG6nyNv0TPcttn7uwf5znhr\n/YMlX76+Mqc7AfVXgncn4NVXgK8DuN4BdCPfsuISLX6EZfHGMAAzMA/APLMconfCnUzYkLyrYuhY\ngH7Z+KlHGSugy2yGwnE7/w+hxoc4CLQwLDK5sotTliW3utosbb98kQE+8SGH5ZnyopNsLucE6m9V\n6u/lV1JCrncAbpaERSCPMlTRuVTLhq3HVk/h2lAoivsP1OkJL773CYjk4UgWrvl5rwf0arWK++YX\nhFM36O5rwK9/Afiu91K6qKwABoFAsk95AIakLR/kEq6JkLXjKnYrX7b+5bKMZZQ5YypL0nODPeDd\nE3ZPcZuwCdAlwHXaL1rjnYC9bNO8M7u0Jt2lJYBO4jeB5zRM4y//34lTt+jua8DXvgD8jvfWFzNE\nXPvRyigCqffrEcX7EvSaxti1dzplnBaotxh+rfdlAN/4PIJr5d37E/Dlvwl853vSSm0FJE7KxTLC\nAZBTSfU443QnyBFn7QxIR9hpjbQ//CcB+oqZ06q1nDVeMpCzGgJ9YwzSK3b2htH5bRuAludlZd6h\nyKk9o0wPIDyziuTp5BRdjFKuLBPl2K/+z6X4HWZ3Ar7yeeA73luPtXtE8x5uXWOy3AJ3jWv3fqCk\nZxxje78trr3jWw8dkmglXr7aF/7VV92BBQpqA2iC2fXE2XDm2Hp2WjpVRUG9nISd7Ku4GWewRQcs\nuLGu663nqAw87R1h1Z7LZKd7qM3WyrNZi6cnnuTzzFi7NuXYYwBqpML+4Nc3Eiyk9ZfEhzWorV3z\n66GIIYfgbr1gyQO/Tq/h4vsWkP1AuWO8vWdIrlSrqMPbfe/fvF0CWxvWkB4YyU0gsHBpJgIRF0fc\nzplTiwgup57e4A5XPssqtRGMdKpmKhEO059F8hag94B7D6gr1bk6zlfsBOLl/qy80YP0KKNDWuRj\nteHybLk1GGadPoDbT28kWujm7TqHrnHqPeSZ72YcXgSvfbj1W0nPoDxridceyN6vJ/57JAmIx1Wf\neB/4+s+gLY6PwKffX3p7/y0WsdPWkdmrrJswGAMmmnDFYdmmKSqkJHa+MlsNU4R5oqg6DoxBv8pv\ni6JhnafaKKnwp8XtmkC+nYOW446SJjw953lqXkpFT4SdMMKuImMZajAjTWn9/c9tiOOm/jzHbpXZ\nfQDepC2QtsbWr41j+4/3jJ0fEdB7OLetvLfeBY6fAi6BVlXp5lPAJ95NY+siK1yI4cvJIrQMcylh\nfIacVy7Fn9xpDJ52dCUOlU/WhFwZxAn4Ydoju9fty6RWfltVat2RwCXPFxzTfdS4SSek4K2kBZdz\nv/MiE+jKMRR+1THxd7ybprRebtTfx9/d5th7xO+HDgWrEdXA/Gwc29JecG4hcWep1Sqp1ivbd4Yb\n4IVceHZx89gYU6P42DsA3SzLexlyVI4AvbatXPzsiSkRIHSGN939NQuYZXKIzeKVnBde9lyrW/ME\n567l3T/3VGHUP2cRPAFxOYBBO7nkr2L3IoKLoky2x+QZ+Zb0EQETt8DH30nuOzePHdWfmt62sqdN\nddPjiN6ennlJae9v96QeMVLD1cz3/CiAO+CDLyQt+eWrwOlt4FPvJ059vClFcAW3groQwStuGL/g\nebZLUfmY57ABgHkQgMsJKyx8nWcMnNRxTc7dA+5WGq1fhVMzlqOOdANMcic7gfp2OZOM5OgiUNyW\nt4YVFlw/9KNpnvprXwC++vmkKLt5O4nfH9OVg+gH9n2A3CMZ7SY/cG/TM4ni1q8mw+0gCx6f1/tU\nkjfDDfCx94CPv1de/0ooOTW531qgjrJbyT4zGaXahIFPANJppBMGjDwlg8m4AeIpbvSRjcpzRM2x\n9Npm0t1sy4V31n3GCWdRlF3omI8zYhri5LQ6owjcuEnz1DqlZcEbAboX4P77W+G8+94gt6Duo2fS\ninu/R+DWvdw5em/LRHeR+1NjPZgj0RuVZ58e95hvDeEDBtF2pyuFEkAOfMWRrxg56Y3BlG7/4YBb\nt7h37dlSqzOq5I0JclupXlN7wJXS3PWVxnwvdV4qSiNmGvMavM00tADTW7/RdFcv9+6hR+PUlvpF\n8w/JktJHoEhKiSqlh7Pr73o4fHQMtLr9dTc+nKWoGII6Us488AEXAMyya0zAfsQZE19w4nRO9sBJ\nFB/1Zr7HAvcWqCt20vgPuJLuSz+muWuSzR0O7BMOmERs3+QJW+lXe3Zu/c0C2Np7QL2HU78meuIl\npcCTg9pz7gjE3r8FcAtmb9v4BvN+JIZb8oCIGoHxS6L4gCuPKUlMy7wuT0ljznorHstYe1of/bMl\nCraee/QBFe6toviVDvmM8DMl8fuOTivRfMKYx+PNzkPdLXHX5zcCX0sMr4G8Fi+C37b8HkyPMMYm\nov8UwJ8A8BVm/qGHJeCRQe4B2gNk/1v0HIHZN2Q2fi0FmactYMuUTjpZkzKnHlgAjElAjSx+jzzh\nwNckys61eJ3tv+/z5dMcPVfyPQ/LiroLyYYXusErusUd3SRFGSWuns8ZFwUbLLh7gB1RT/167u3d\nW1y59i3v/5qoh2P/JQD/EYAf3xf1E4PakuWa9tmHiSqbUFasvZvcultj50hRtpVeny7jn7YiAuB0\nf2cCsU5nsXBmFb+vOHBaQ10c/+Pjjp6j9HiqicYtUZzIieKyGIXSXLUqycqji5cpsVW8e9UzUf5r\nwO4BNRr+/ruvGdBKm8Bm5r9BRN/79EkJvy72PTuFPWD23FrdLaWYt3uTWWtw2b2sf2b3OzHjyteV\nufABBz5g5APCeW7pFHIStzi3zZt7zqBEPE+dxs4HEcW9OWJ1AGVUxrVO875lbO0WiCOl2hbAHwTk\np+kJnnCM7dlmxEZr70TvwvxmidZe9wEz0AY0YZ20vdykxalrDSfwS2PuNAY/8xEDT/mSPuZ0tQDl\nW0TMGaDC5bvS18hLvpJHFsJynrNOz+n621ucKU1nXemAmcbMoatAbg1t7lvWkb1lauCubfNsGQDr\ngxbgAng/n+j99MTKsy1wByAtUAi0W5j8IxfOgtk+2yRYgAPxVFbLwIX3WWil+Z6ABkOuBRowzWmb\n5zDLQUsy0zXNw3InhpnrTkmb62nr4dgalGQMjRGTna8WW28+OQ9yQMIg89SDOWG2VqZ+5V7LHZEv\nf4+fHjB6EPesTgsZL3diswbyasSb9EjA/ivG/QNilCJO3QJti7Ob31p9hOXOPYD2bka5iqwlKgL1\nxuYblHXvBLT9beakXLvyiIGPedytBykeeIK9MTrd3SXz3C3q5NgzKAM6T13J1sqsCacjLvMxzV3P\nMk+t6+tbYnitgwXW5bzFwSNQW/9aGW+BuhaH9/dpCR+2uLb3/5KYNj0SsP+Uce8Fq/5u3+eKO4gu\n/8zI67UtiG2BbwFa3a0GFgEbzt9+M8rmHiD7OVhGmudmmecWLq6g1nnuo+yMAvS8cjPPHaWp8HCZ\nIx9CprPy4RDHctEJHfP4WsfaqgHv4tJbonqQpnuXf6seerTn/r38Ha53LKFHC8yW/mmUjPMn/Se2\neAAAENRJREFUfaQA+qa7/isAfxjAp4jo7wP4t5n5L9XfqHFji7aaeN7JrVs/Rxw6ega2FWU+Wb5R\n1RpXjUtYvx1A9ka585VHATVlUF94wo1sewRIFGZJcK4ed8j5n0sshUWuhxDmE0/ypYNpO2YC8rgy\nMxHCO9b8rEMP0IG4HqI6sM8ejNZdM9HWzlbHEH1/RVtcu+e9OvVoxf+lrpjCBNTArW6gzaF9J+Ge\nfV/hX/XPFshKXvy28T5UDGzZHQCODeV5bCt++3luzYIeGzhD5rnDdlFrLHFnqkc5qcCvV97eIc1V\n52tx3d3XVeVZC+g1wLfsWtZqDDICqNZPa2oMFbv22ypAT6Kid7bpkURx+9EaN94Tl+fsLS4feNUK\n1nNsJSt+E8rP14CNxnNUd73cIeIEc2nzLOeVzjLPPSOJ2zNAmJEOVJQ93HIY4sTL3u6HEkN2n9Fo\nzijTvdW3YAzp1iMQkG3ZqmqX4EZ2L+CBdZ14KatM9Pr5Ptx7K1zte81E7WL1XWGeQCsecVkPSv29\n5a6BW2y5ZC6PqwEUx5ZEY62IY+vvXlmGhtvSXm7dC+Lsz2HYZZ6bXBwDDjjggLEwelvlBSPsZTtL\nFpYLCfx1PPqrui9Id39feLlj6yzmguO6LIoq5aCKqazqXrCjYfs0RH5bzHFu2E2wmw+y/OPCI7C9\n35a7Tc80j+3BHv3eAjKwbg3Gz4LbRg8XHKiPpe3y0Sg5NWDD/O6/fV+uoA0hbDhc+gVcgMFIt1YS\nLhgx4oABRxAmMMqbPJebSBZ/e/NleSVg8n+FI+5wwBkjrhjybnEGYxErbJn4IQCVYeRcuNXps7aj\njQBdc9eKxsd93w64BujoG1UAtyJogb6PnkEUt34tAEfvtDKmXNuB2ydlD5fe4gR7JNk9HGGrsSiY\nW27zIoMxgTFhwAUDBowgHMA4YQIv89zZpIJL8JwEnoRJRG49yF9NOoDwKBw6gXsFbC00bdzKlW1y\nyYRRgBfPiOtlj0RVoz3ABup1FfnnuOSBez7gEwX3TpTQOj2DKO79LLj///bOJVSSs4rjv3PvnUei\nAUXRqBlIFgkGFAwREYPim1EkulNBBBeuFAVRRFfuxJVZuExcBCEqipLoQtFkkSAkRmbwMcmQQAbi\nRGMgMRom87i3j4uu6vvV6fO9qqu7bq71h6bO96jzna7z/c/3qOpqHNmeY+vQzV+QO1AdqrBT8pDM\n3ihtO4qadAlKSO3Jrp9N52g7jCV4oGyGMgN2EbbYakg9a2hLc9f5CjvADkp7C0wa0s+34GAX2KW9\nrbWz+AfM/dcFb3OlGbH3cIg9f3ujIbe5qNZ/TXYvQvf1kU2nfBMju9W9kGNKPeVeWezcNNYwYkOa\nrLkpt9VToL+jzrLbNBMjs7dqCNX06TSe/yiUF98lIHBL6KWy4Eg7YgtbzVR81pB9D2EX2GGLo8Cs\nIWLz3yVNOJiP6TOU3WY6Pyfwkf1XGS1e7r/TEN+bigvdZwuci96SOiR3agnkyesidpiX42HnExSE\nfosqTdWBZSNtno81rbFDdliCe6TO6fFgAoY2edKmpVOcJbPtILZz5UwNzc1F/lSdxVG7eSGp25Ms\nuZtP+5PiXdq/Dpq/23Qb4QrCsaaWLibcwnaH2DP2gD2EK2wtNsYuLV7ov7P458u95o1se7SBYm//\nYi6Ru7mY7QZah9zBtbbHWF6uTgoese2xhNAd2fGZ26j3yRkb0+ljg8+Kx8juybUwJIb9DpUyp7Rz\npDqLZ3IqqlNxtArsT700VLz/mY/Y0P5/pwBbzbR8qxmN2+21dgU+H3H3iT0fsWmIvc1FdrjU3K9u\n/1Cvu+U2b3lpd3Lp/6hthJQgYEmZD3IjuEXKb6ljibw4aoSjKTLniB6T81jjVNyWpUZsW8fqlsQx\nZU6kjmdCrKnQLK+Pxtou8ZPV4eW5ilL3XfZ3xeeT7PCx0naDrCX1drMCb/fLW2K30nyEv8w2lzjC\nyxzhIscI/6ieztEh9oLg5iH8zneUZdEug+wxNyW3iPkudf1LfRf7F5Zkg6mOEGsw28gCa9oVt+mw\nToxFVo/VFzt69cMq2p0GWrNLSJ1qzvNHe7SR3Au6RWTOfZYJ3lUbfu/5ja/wNtWsCQYsUrPFBtys\nmWbvNedcYTu4MDaAK8tM8x7CT8yTQ8NjZLREbqf6MWLnfFdD7sVR/fqdk7zIDunX3NhO0o/ca9wV\nB3/Ys0aV1OkLQ/KWZB2yZkb0MM+a68GSuc3TToXluktKaokNy886erL3Br9Z5LMXyQ8vYEwu+XF7\n7BOY3fqo44/Ar12BZcc5xSXB1ZI4yy9bkPqk6lgdnpzHmqfinhGpYdHWGcAsG9XDfAV3LR42Xxpz\nvM6y0BFGd0deUtD3QXJ7XkyOEdl+FJ/8KUJbv9YSOnzoILjo1n9LdQJYn3llS3pM3pKvvHNTvuzz\nseeSkNMYYSpOpO7Q5G7OV1h6iGXJHKWzMxsGhBJYXWrlMOKbXrLoFJ7zev9SJKIvRexYmUdub5Eb\nykL/kVq7+hY+awOwNMEY9m+Vwf5fGAfn5/zoBWMv7ea3/vMixbpJnefFBqfiYTpGaMu6IUbu1vkB\nub2mcr2gyAwT4aHreHdXuz1hDGJ7pI0RPywrIXbfkTrMb3yyCM7tdZN9Ui/Owe8uNd0nGqBtuSW1\nd1zFTzZdR2rY6FQ8NSLXlPVB2wmcETs2obDlWVhSB85YemghF61tWS3BbV2bTk2zS8jdBkE7Bcek\n+4zU9hOQG+iM2NAl+CJoZ4K0Ra57uX0jRbohgm+qLI+RdsVjckm9Gjj6Y0+mrQR1xW6VUiIPTWwr\n54jrkdvWt4S25F51pPbqLS5kQOCQ8IGaBfEK/VxFbDVyjpxDk3qjxA6RI7mSJm2K+CVfyjpUnLKg\nTm20LkLKyRC7/1zeCVYldgmZUwEgNQ2H9AvjYoQtDQbB62U1PBeT1yaGDOJqjmF+bTAmIcfKyjvj\nhn7d1eanCIxJ92GU0n2vsEVMb2w+3he1RPQeOIF858gFkFh5zQ54mB9+7BQ6lGcmL0X09jvaOt5P\n8tprEiNr25bNWwWWbKEc80WJ3+ghWxvi2MCIbfNjHaIts4anAoCF7VRh27GRwtpSC8+WEmcOPWp7\nZRSW27LU97Ck9WRMXswnsZ/aEaTtaB073waKEEP5NnatavxUS2rPljy5N7R5FiN4DaFtwEiR216E\nGJmtbbXIOT6USz9DjdIltqRG5FTgKb02uQCbIrwHS3YvyETW5dV+ThGpNoDGAmZOTpWnsYFfd+XS\nVi5dY8fInXNijNA1RI+1G8op59aM1il9XhkZ2U73S4NEeG4OKVLH6nujt7e2zr3wPfxOubyU/V46\n5uNcOlfHk1PleYz0BhVbZg1OkbvmXM/OVIcrXdtbB9syj8heWe0npr+G3GE6FWByI7b1ZcmU3MtL\njd4hZol6HtHbfDJ6PcT8W0q8VUmdasfa5GPEN6jE8mOEtsEjRe6wXdvhYu2Uktrq9+S+RK45l4ic\n6yS52QLO+V6dUlhShySM1S8ZnWNEDjdPUwNKqe32nBy5w7xcYM3JsXPzGOENKpj8FEFLR+scKT1S\nr7K2Du0I5RJC9CF7LblL5FRgKbGpBJZQMUJb8qZG5/Cc1OaZN4jU+Dzn39wxd+1zZbbcS8cxwhtU\nbHmYtufk0jXkDtv2RpBa5JxaSmRYz22tmKzEyVszNc9Nw3OksoT0ECOvN323dZRlctfA+rWVc2mb\nX+ub0rbT2MAbVCA/8rbp1FQ9TPe1JXT6KgQvcXYJKVetE5NjwcbLs+ekysLrFvNtTm715K51S+rc\nqJ0KDkP4uJTgsWONnNNl68Sx5s2zsDy13imJ+G29ko7ljfKenal0CrmLX0NYq7e0TkwuJXXs+5TU\ngdXIHVuehSjZ+4iVrUJoqydFulydEtme57VTR2rY2OaZzc+RO6xnZSrKMfkpPUNE8xSJvfxY3Rzh\nScieHX1stnlh3VZedeRu86zdFikCp540XAUlhPbycoG1T3mqTR8b+icQL78vufuM1qk2bb0S5CJ5\nX1J7bdSQurSTpDpXTE+M1EOQm6AsbMOWeXm522TW5hrEiFTif1tvc6SGja2xwZ8a1ZIbJy9F7tR0\nrs+6vTZilxAylVfbGUo7SY19FrnrWyPXTrNL8iyB23Q4M4PyYOClYz63dfoQ3Gu7jtSwkbeUenVy\no2jftXWsrm3XYuj1dXscguiezho5ZZP9Xl55TO/Q5O5DWi8vtrQK0338HZNLSV3i01g61l4ca77d\nFSI1YteQ2+pMlYf5uam4pzNELnrXEJtEWUpPqtyTvSOmrkf0XFslM6PadXYJufuS38L2sZy+2Pl9\n+wAVae+YR5bYInISuAPYBu5U1e8Va+8gNCrVMcK6JRtgfafiMf2xi1cauXPkjcmpvFQHyOlfrqNc\n5jLneJmzwIsc5xhv4HUc43jGhlDf0OQmUZY61/abFKlrSJ/KLyFfitQeUUtJXUbuJLFFZBv4AfBh\n4DzwRxG5V1UfK9JeDI/cUOa0VJ1Q71ngrU59InrDci8dc04o1xAc4AngRqddT6e1IeZ4m/9TLnIe\n5b+0D8ZcAJ7nac5xlNdyFTfztkj7AE8BNzAMub3v5021c3qszseAm+lei9wMcKg+4KU9nAVuypxf\nRmIPsd/HtXgX8KSqnlPVK8CPgU/2bm1UPD62AQV4cs36LwPPoLzI8u+WlYtc4t9cYK/5900f59Zm\n3XA4O7YBBXhirdpzxH4L8HSQ/nuTN+EViXPAf5I1LnKZFzJ1Jhx85Ijdfy4w4QDiDP4bRvahwHM8\ntxFrJqwPohrnroi8G/iOqp5s0t8CZuEGmohM5J8wYUSoLv9PVY7YO8wXLB8CngEeAT47/ObZhAkT\nhkRyV1xVd0Xky8BvmN/uumsi9YQJBx/JEXvChAmvTOQ2z5IQkZMi8riIPCEi3xzKqKEgIj8UkWdF\n5C9j2xKDiJwQkQdE5G8i8lcR+crYNlmIyHEReVhETovIGRH57tg2xSAi2yJySkTuG9sWDyJyTkT+\n3Nj4yNra6TtiNw+vnCV4eIUDtv4WkfcCLwF3q+rbx7bHg4hcC1yrqqdF5NXAn4BPHaTrCCAiV6vq\nhWbf5SHg66r60Nh2WYjI14BbgWtU9fax7bEQkaeAW1X1+XW2s8qIfeAfXlHVB4EXxrYjBVX9p6qe\nbuSXmD829eZxrVqGql5oxKPM91vW2jH7QESuAz4O3En+2dExsXbbViH29PDKwBCR64FbgIfHtWQZ\nIrIlIqeBZ4EHVPXM2DY5+D7wDXI368eFAr8TkUdF5IvramQVYk+7bgOimYb/DPhqM3IfKKjqTFXf\nAVwHvE9E3j+ySR2IyCeAf6nqKQ72aH2bqt4CfAz4UrNcHByrEPs8cCJIn2A+ak+ohIgcAX4O/EhV\nfzm2PSmo6ovAr4F3jm2LwXuA25s17D3AB0Xk7pFtWoKq/qM5Pgf8gvmSdnCsQuxHgRtF5HoROQp8\nGrh3GLP+fyAiAtwFnFHVO8a2x4OIvF5EXtPIVwEfAU6Na1UXqvptVT2hqjcAnwHuV9XPj21XCBG5\nWkSuaeRXAR8F1nLHpjexVXUXaB9eOQP85ADu5N4D/AG4SUSeFpEvjG2Tg9uAzwEfaG6BnGp+A3+Q\n8Cbg/maN/TBwn6r+fmSbcjiIS8U3Ag8G1/FXqvrbdTQ0PaAyYcIhxEoPqEyYMOFgYiL2hAmHEBOx\nJ0w4hJiIPWHCIcRE7AkTDiEmYk+YcAgxEXvChEOIidgTJhxC/A+SnGev7o1dugAAAABJRU5ErkJg\ngg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# collect Bader volumes associated with atom #5\n", "mask = np.zeros_like(rho, dtype=bool)\n", "for v in (bdr.nnion == atom+1).nonzero()[0]:\n", " mask[bdr.volnum == v+1] = True\n", " \n", "plot(si.positions[:, 0], si.positions[:, 1], 'k.', ms=20)\n", "plot(si.positions[atom, 0], si.positions[5, 1], 'g.', ms=20)\n", "imshow(rho[:,:,0], extent=[0, si.cell[0,0], 0, si.cell[1,1]])\n", "imshow(mask[:,:,0], extent=[0, si.cell[0,0], 0, si.cell[1,1]], alpha=.6)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Wrapping Castep with f90wrap\n", "\n", "- `f90wrap` can now wrap large and complex codes like the Castep electronic structure code, providing deep access to internal data on-the-fly \n", "- Summer internship by Greg Corbett at STFC in 2014 produced proof-of-principle implementation\n", "- Results described in [RAL technical report](https://epubs.stfc.ac.uk/work/18048381)\n", "- Since then I've tidied it up a little and added a minimal high-level ASE-compatibility layer, but there's plenty more still to be done\n", "\n", "![Castep](http://www.castep.org/files/CASTEP_Logo_mini-01.png)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Current Features\n", "\n", "- Implemented in Castep development release:\n", " make python\n", "- Restricted set of source files currently wrapped, can be easily expanded.\n", " Utility: constants.F90 algor.F90 \n", " comms.serial.F90\n", " io.F90 trace.F90 \n", " license.F90 buildinfo.f90\n", " Fundamental: parameters.F90 cell.F90 basis.F90 \n", " ion.F90 density.F90 wave.F90\n", " Functional: model.F90 electronic.F90 firstd.f90\n", " \n", "- Already **far** too much to wrap by hand!\n", " - 35 kLOC Fortran and 55 kLOC Python auto-generated \n", " - 23 derived types \n", " - ~2600 subroutines/functions" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## What is wrapped?\n", "- Module-level variables: `current_cell`, etc.\n", "- Fortran derived types visible as Python classes: e.g. `Unit_Cell`\n", "- Arrays (including arrays of derived types) - no copying necessary to access/modify data in numerical arrays e.g. `current_cell%ionic_positions`\n", "- Documentation strings extracted from source code\n", "- Dynamic introspection of data and objects\n", "- Error catching: `io_abort()` raises `RuntimeError` exception, allowing post mortem debugging\n", "- Minimal ASE-compatible high level interface" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Test drive" ] }, { "cell_type": "code", "execution_count": 68, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import castep" ] }, { "cell_type": "code", "execution_count": 69, "metadata": { "collapsed": false }, "outputs": [], "source": [ "#castep.\n", "#castep.cell.Unit_Cell.\n", "castep.model.model_wave_read?" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Single point calculation" ] }, { "cell_type": "code", "execution_count": 70, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Energy -401.419912223 eV\n", "Forces (eV/A):\n", "[[-0.12335761 -0.12340173 -0.12332035]\n", " [ 0.12335492 0.12337595 0.12332078]\n", " [-0.123409 -0.1233628 -0.12327315]\n", " [ 0.12330122 0.1233129 0.12334265]\n", " [-0.12330156 -0.1234069 -0.12336426]\n", " [ 0.12335953 0.12338503 0.12331035]\n", " [-0.12334791 -0.12329187 -0.12338976]\n", " [ 0.1234004 0.12338941 0.12337375]]\n" ] } ], "source": [ "from ase.lattice.cubic import Diamond\n", "atoms = Diamond('Si')\n", "calc = castep.calculator.CastepCalculator(atoms=atoms)\n", "atoms.set_calculator(calc)\n", "e = atoms.get_potential_energy()\n", "f = atoms.get_forces()\n", "print 'Energy', e, 'eV'\n", "print 'Forces (eV/A):'\n", "print f" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Interactive introspection" ] }, { "cell_type": "code", "execution_count": 74, "metadata": { "collapsed": false }, "outputs": [], "source": [ "#calc.model.eigenvalues\n", "#calc.model.wvfn.coeffs\n", "#calc.model.cell.ionic_positions.T\n", "#calc.model.wvfn.\n", "#calc.parameters.cut_off_energy" ] }, { "cell_type": "code", "execution_count": 75, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "(-0.5, 0.5)" ] }, "execution_count": 75, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfMAAAFrCAYAAADFOmBlAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHdNJREFUeJzt3XuQnXWd5/H3t28JuSdccmkigQIh4IqACYyDTIvDGLyh\n7CoybqS8IGOZ0Z2pUkBrNK5VO4M1rq6DMhllXHbK64pocEUFTCtlAUmmCIkQLpGEIQkJwSRASEi6\nO7/9o09C0+lOX87p8zy/7ver6lSe85znPL/vydOnP/37PbdIKSFJkvLVUHQBkiSpOoa5JEmZM8wl\nScqcYS5JUuYMc0mSMmeYS5KUuaaiCzgkIjxHTpI0pqSUohbrKVXPPKU04ONHD/2Iy39w+aCWHe7j\nuecSEyeO3PpH2+Pzn/984TX4cPuN1YfbL99HLZUqzAejqaGJzoOdI9rGpEmwbx90jmwzkiTVhGHe\nh4YGmDIFnn9+RJuRJKkmsgzzjq6OEW9n6lTYvXvEmxkV2traii5BVXD75c3tJ8gwzJsbm0e8Zw7d\nYf7ccyPezKjgL5O8uf3y5vYTZBjm9RhmB5g2zTCXJOXBMO+HPXNJUi4M8364z1ySlIssw7zj4Mgf\nAOcwuyQpF9mFeXODB8BJktRTdmHuPnNJkl7JMO/HtGnuM5ck5SHLMK/XRWPsmUuScpBdmHvRGEmS\nXim7MHeYXZKkVzLM+2HPXJKUC8O8H4a5JCkXWYa5F42RJOll2YV5vS4aM348dHXB/v0j3pQkSVXJ\nLszrNcweAdOnexCcJKn8sgvzhmjgYDrIwXRwxNuaPh127hzxZiRJqkp2YR4RdeudT58Ou3aNeDOS\nJFUluzCH+g21z5hhz1ySVH5Zhnm9DoKzZy5JykGWYV7PnrlhLkkqO8P8KDwATpKUg2zDvB53TnOY\nXZKUgyzDvF53TnOYXZKUg6rDPCIWRcQjEfF4RFx7lOUWRERnRFxebZsOs0uS9LKqwjwiGoEbgUXA\nmcCVETG/n+VuAH4BRDVtggfASZLUU7U984XAhpTSppRSB/B94LI+lvtr4EfAjirbA+yZS5LUU7Vh\n3go81eP55sq8wyKile6Av6kyK1XZZt3unOYBcJKkHDRV+f7BBPNXgetSSikigqMMsy9duvTwdFtb\nG21tbX0uV++LxqTUfeMVSZKGq729nfb29hFZd6Q0/I5yRFwALE0pLao8vx44mFK6occyT/BygB8H\n7AWuTikt77WuNNhaLvjWBXx10Ve54MQLhl37YE2YADt2wMSJI96UJGkMiQhSSjXpKlbbM18NnBYR\n84CtwBXAlT0XSCmdcmg6Ir4N3N47yIeqXvvM4eXeuWEuSSqrqvaZp5Q6gSXAL4GHgR+klNZHxDUR\ncU0tCuxLvS4aA95sRZJUftX2zEkp3QHc0Wvesn6W/WC17UH9LhoDHgQnSSq/LK8AV89hds81lySV\nnWE+AM81lySVnWE+AIfZJUlll22Y1+OiMeABcJKk8ssyzOt10RiwZy5JKr8sw9xhdkmSXmaYD8Bh\ndklS2WUb5vW6aIw9c0lS2WUb5vXqmR97LPzxj3VpSpKkYckyzOt5ANxxxxnmkqRyyzLM69kznzoV\nXnwRDhyoS3OSJA2ZYT6ACIfaJUnllm2Y1+uiMdA91P7ss3VrTpKkIckyzOt51zQwzCVJ5ZZlmNdz\nmB0Mc0lSuRnmg2CYS5LKzDAfBMNcklRm2YZ5va4AB4a5JKncsgzzel40BuD44w1zSVJ5ZRnmnpom\nSdLLsgzz5sZmw1ySpIo8w7yh2X3mkiRVZBnmLY0t9swlSarIMsybG5s50FW/O59MmAAHD8LevXVr\nUpKkQcszzOs8zB5h71ySVF5Zhnm9h9nBMJcklVeWYV7vYXYwzCVJ5ZVnmNd5mB0Mc0lSeWUZ5i2N\nLfbMJUmqyDLM633RGDDMJUnllWWYtzS2FDLMvmNHXZuUJGlQsgzz5ob6HwB3wgnwzDN1bVKSpEHJ\nM8wLGGafORO2b69rk5IkDUqWYV7EMLthLkkqqyzDvIhh9lmzDHNJUjnlGeYFDLNPmQIHDnh9dklS\n+WQZ5kWcZx7hULskqZyyDPMirgAHhrkkqZzyDPPKMHtKqa7tut9cklRGWYZ5QzTQGI10Huysa7v2\nzCVJZZRlmIPnmkuSdEi2YV7UuebbttW1SUmSBpRtmHuuuSRJ3fINc4fZJUkCMg7zIs41N8wlSWWU\nbZgXca65YS5JKqNsw7ylsaXuw+xTp8L+/bBvX12blSTpqLIN8+bG+h8A5yVdJUlllG+YF3RJ11mz\n4Omn696sJEn9yjbMixhmB2htha1b696sJEn9yjbMixhmh+4w37Kl7s1KktSvfMO8oGF2w1ySVDbZ\nhnkR55kDzJljmEuSyiXbMC/iCnBgz1ySVD7ZhnkRN1oBw1ySVD7ZhnkRN1qBl49mT6nuTUuS1Kd8\nw7ygYfbJk7svHvP883VvWpKkPmUb5i0NxQyzg0PtkqRyqTrMI2JRRDwSEY9HxLV9vP7+iHgwItZG\nxO8i4rXVtgnFnWcOhrkkqVyqCvOIaARuBBYBZwJXRsT8Xos9AVyUUnot8EXgX6pp85DmhmKG2cEw\nlySVS7U984XAhpTSppRSB/B94LKeC6SU7k0pPVd5ej9wYpVtAsWdZw5e0lWSVC7Vhnkr8FSP55sr\n8/rzYeDnVbYJVA6AK2ifuReOkSSVSVOV7x/0CVoR8SbgQ8Cf9rfM0qVLD0+3tbXR1tbW7/qK7pnf\nfXchTUuSMtXe3k57e/uIrLvaMN8CzO3xfC7dvfNXqBz09k1gUUppV38r6xnmAxnXOI7n9j838IIj\nYO5ceOqpgZeTJOmQ3p3UL3zhCzVbd7XD7KuB0yJiXkS0AFcAy3suEBGvAn4M/NeU0oYq2ztsXNM4\n9nfur9XqhuSkk+DJJwtpWpKkI1TVM08pdUbEEuCXQCNwc0ppfURcU3l9GfA5YDpwU0QAdKSUFlZX\ndnfPfH9XMWF+/PGwdy/s2QOTJhVSgiRJh1U7zE5K6Q7gjl7zlvWY/gjwkWrb6a3InnkEvOpV3b3z\ns84qpARJkg7L9gpwRfbMAebNc6hdklQO+YZ5U7Fh7n5zSVJZZBvmLY0thQ2zQ3eYb9pUWPOSJB2W\nbZiPaxxX2HnmYM9cklQe+YZ5wcPs7jOXJJVFvmHeWNzR7OAwuySpPPIN84J75rNnw86d8NJLhZUg\nSRKQc5gX3DNvbOy+RruXdZUkFS3fMC+4Zw7d+80dapckFS3fMC+4Zw5wyinwhz8UWoIkSRmHeQl6\n5qedBhtqdusYSZKGJ98wL0HP/NRT4fHHCy1BkqSMw7wkPXPDXJJUtGzD/NDlXFNKhdVw6qmwcSN0\ndRVWgiRJ+YZ5U0MTEUFXKi5JJ0yAY4+FzZsLK0GSpHzDHNxvLkkS5B7m7jeXJCnzMC9Bz9zT0yRJ\nRcs7zO2ZS5KUeZiXoGd+6qnw2GOFliBJGuPyDvOS9MyffBL2F1uGJGkMyzvMS9AzHz+++4Yr9s4l\nSUXJO8xL0DMHeM1r4Pe/L7oKSdJYlXeYl6BnDnDWWYa5JKk4WYd5S2NLaXrmDz1UdBWSpLEq6zAf\n11SOnrnD7JKkIuUd5o3jONB1oOgyOPVU2LIFXnyx6EokSWNR3mFekgPgmprg9NNh/fqiK5EkjUVZ\nh/n4xvG81PlS0WUA3QfBrVtXdBWSpLEo6zA/pvkY9nXsK7oMAM45B9asKboKSdJYlHeYNx3Dvs5y\nhPnrXw+rVxddhSRpLMo7zEvUMz/3XHjwQejsLLoSSdJYk3eYl6hnPmUKtLZ6EJwkqf7yDvMS9czB\noXZJUjHyDvMS9czBMJckFSPrMJ/QPKFUYb5gAaxcWXQVkqSxJuswP6b5GPZ27C26jMNe//rufeZ7\n9hRdiSRpLMk7zJvKtc98/Pju883vu6/oSiRJY0neYd5crn3mAG98I9xzT9FVSJLGkrzDvGQ9c+gO\n89/+tugqJEljSd5hXsKe+RveAKtWwYHib+YmSRoj8g7zEvbMp06F+fPh3nuLrkSSNFbkHeYl7JkD\nXHop3HFH0VVIksaKvMO8hD1zMMwlSfWVdZiX7aIxhyxcCFu2dD8kSRppWYf5+Kbx7O/cT0qp6FJe\nobERLrkEfv7zoiuRJI0FTUUXUI2IYFzTOF7qfIljmo+p2XqffuFpblp9Eyu3rKTjYAdnHHsGl51x\nGW8++c00NjQOah3vfjfcfDNcfXXNypIkqU9Z98yh9jdbWf7ocl77z69l175dLFm4hE+/4dOcNO0k\nPnP3ZzjzG2fynbXf4WA6OOB63vY2uP9+eOaZmpUmSVKfoixD1BGRhlNL6/9sZeVHVtI6pbXqGu78\nw50svm0xt195OwtaF7zitZQSKzat4Lq7rmN803iWvX0Z84+ff9T1/eVfdl9E5mMfq7o0SdIoExGk\nlKIW67JnXrHjxR0svm0xP3zPD48Icuj+T7/45Iu598P38t6z3stF//sivvS7Lx21l/6+98F3vlN1\naZIkHVX+YV6jO6dde9e1vP8/vZ+LTrroqMs1NjSyZOESVl+9mp899jMu+bdL2Pz85j6XvfRS2LgR\n1q2rujxJkvqVf5jX4FzzDTs3sPzR5Xzuzz436PecNO0kVly1govnXcx5/3Ietz586xHLNDfDRz8K\nX/96VeVJknRU+Yd5Da4C9/f3/D1LFi5h6vipQ3pfY0Mjn73osyx/33Kuu/s6PvjTD/L8/udfscxH\nPwo/+AHs2lVViZIk9Sv7MJ/QPKGqYfbdL+3m1vW3smThkmGv4/wTz+eBax6gpaGFs//5bO558uV7\noM6eDZdfDl/5yrBXL0nSUWUf5pNbJrPnwJ5hv/+7677LW059C8dNOK6qOia1TGLZO5bxtUVf44of\nXcH1d13Pga7uW6f93d91D7U/+2xVTUiS1Kfsw3xSy6Sqwvzba77Nh173oZrV847T38Gav1rD+mfX\nc86yc/jNpt8wb173ke2f/3zNmpEk6bBREeYv7H9hWO/dtHsTm3Zv4s9P+fOa1nTCxBO47Yrb+OKb\nvsji2xaz+LbFfPy6rdx2G/zudzVtSpKk/MO8mmH25Y8u5+2vfvugL9E6FBHB5fMv5+GPP8ycSXO4\n8LuvYcFnP8Xia55l9+6aNydJGsOqDvOIWBQRj0TE4xFxbT/LfK3y+oMRcU61bfY0qWUSLxwYXs/8\np4/+lMtOv6yW5RxhUsskbrjkBtZ9bB0nztvL0//5dF57/RLWbVs/ou1KksaOqsI8IhqBG4FFwJnA\nlRExv9cybwVOTSmdBnwUuKmaNnubPG54PfNd+3axassqLjnlklqW06/WKa18/W1f55FPrmVc1wwW\n3PgmLvrXNr6x6hts27OtLjVIkkananvmC4ENKaVNKaUO4PtA767uO4FbAFJK9wPTImJmle0eNtwD\n4O564i4ufNWFTGyZWKtSBuWkGa08dON/522PP8m2n/w3fvHw75j/9fks+OYCPvWrT3H7o7fz1HNP\nle62rpKk8qr2FqitwFM9nm8Gzh/EMicC26tsGxj+MPuKTSt488lvrkUJQ9bSAv/3e+O48cZ38YW/\nfReLr9rPBRfezx86f8M/rfwn1j2zjr0deznz+DM5aepJzJk8hzmT53DCxBOY3DKZSS2TDj+aG5tp\njEYaGxppiAYao5GgkYNdjXR2QGcndHRCV+Xf3n8jDO2PhsEtm46yXPS6pcBg2z/aOoezviGts8C2\na/F/PlKC6u8PEb1/IIa7niHW0lezOX+ePtdRo1rop5ahrL5Mn2fIPytD/tmq1f/70FQb5oP9DdL7\n0/X5vqVLlx6ebmtro62tbcAVD/cAuF9v/DVXn1vczcYbGuATn4D3vAf+8R/H8TeXX8SsWRdx/vmw\n6AyYMOeP7Gx8mD++sJlnd27h4Re38seX1rDnwB72du5hX9ce9qc9dKYOug52cTB1cTAd5CBdJLqg\noYvDP4bB4X/7NKSb9gxu2b6+AP1l3eC/LINcbgQ+T6FtD3qd9fwlUoM/HqJWf4DUYj2j6/PU7I+7\nmnymnP9vh7j8AHWm/zhA+o+OYdZydNWG+RZgbo/nc+nueR9tmRMr847QM8wHazinpm19YSs79u7g\n7FlnD7m9Wps9G778ZbjhBvj3f4dVq2DDBti+6lhefPGNRHT/BXz8ZDhtGkyrPKZO7X5MmXLkvxMn\ndv+xIEkqr1r24qsN89XAaRExD9gKXAFc2WuZ5cAS4PsRcQGwO6VUkyF2GN4BcCs2rqBtXhsNUZ7E\na2qC88/vfkiSNBRVhXlKqTMilgC/BBqBm1NK6yPimsrry1JKP4+It0bEBuBF4INVV93DcPaZ//bJ\n3/JnJ/1ZLcuQJKkw1fbMSSndAdzRa96yXs+HfxeTAQznaPaVW1fy4XM/PEIVSZJUX+UZZx6myS2T\nh7TPfG/HXh599lHOnln8/nJJkmoh+zAf3zSezoOddHR1DLwwsGbbGs48/kzGNY0b4cokSaqP7MM8\nIoY01L5qyyoWti4c4aokSaqf7MMchrbffOXWlSyYs2CEK5IkqX5GTZgP9oj2VVtWsaDVMJckjR6j\nIsynjZ/G7pcGvq/orn27eHrP08w/bv6Ay0qSlItREebTj5nOrn27Blxu9dbVnDv73BG5f7kkSUUZ\nHWE+fjq7Xho4zFdtXcXCOR78JkkaXUZPmA+iZ75yy0r3l0uSRp3REebHDL5n7pHskqTRZlSE+Yxj\nZgzYM9/y/BYOdB1g3rR59SlKkqQ6GRVhPph95od65UXdOF6SpJEyOsJ8EMPsXvlNkjRajY4wH8QB\ncF75TZI0Wo2OMD9mOjv37ez39ZQSq7eu9kh2SdKoNCrC/PgJx7Nj745+X9+wcwNTx03lhIkn1LEq\nSZLqY3SE+cTj2blvJ10Hu/p83fPLJUmj2agI86aGJqaPn95v79zzyyVJo9moCHOAmZNmsn3P9j5f\nW7llpUeyS5JGrVET5rMmzWL7i0eGeUdXB2u3r+W82ecVUJUkSSNv1IT5zIl998zXbl/LydNPZvK4\nyQVUJUnSyBtVYb5tz7Yj5t+3+T4uaL2ggIokSaqPURPmrVNa2fz85iPm37v5Xv5k7p8UUJEkSfUx\nasL8lOmnsHH3xiPm37f5Pi440Z65JGn0GlVh/sSuJ14xb8eLO3h277OccdwZBVUlSdLIGzVhfvK0\nk9m4eyMppcPz7tt8HwtbF9IQo+ZjSpJ0hFGTcpPHTWZC84RXnJ5298a7aZvXVlxRkiTVwagJc4DT\njz2dh5556PDzO5+4k0tOuaTAiiRJGnmjKswXti5k1dZVAGx+fjPb9mzj3NnnFlyVJEkja1SF+fmt\n53P/lvsB+PH6H3PpqZfS2NBYcFWSJI2sURXmb5j7Bu558h4OdB3glgdv4QNnf6DokiRJGnHR8+jv\nIkVEqkUtF99yMRNbJvLEridY+1dr7ZlLkkopIkgpRS3WNap65gDL3r6MyS2T+eF/+aFBLkkaE0Zd\nz1ySpBzYM5ckSYcZ5pIkZc4wlyQpc4a5JEmZM8wlScqcYS5JUuYMc0mSMmeYS5KUOcNckqTMGeaS\nJGXOMJckKXOGuSRJmTPMJUnKnGEuSVLmDHNJkjJnmEuSlDnDXJKkzBnmkiRlzjCXJClzhrkkSZkz\nzCVJypxhLklS5gxzSZIyZ5hLkpS5qsI8ImZExJ0R8VhE/CoipvWxzNyIWBERD0XE7yPiE9W0KUmS\nXqnanvl1wJ0ppVcDd1ee99YB/E1K6SzgAuDjETG/ynYlSVJFtWH+TuCWyvQtwLt6L5BS2pZSWlOZ\n3gOsB+ZU2a4kSaqoNsxnppS2V6a3AzOPtnBEzAPOAe6vsl1JklTRNNACEXEnMKuPlz7b80lKKUVE\nOsp6JgE/Aj5Z6aEfYenSpYen29raaGtrG6g8SZKy0N7eTnt7+4isO1LqN38HfnPEI0BbSmlbRMwG\nVqSUzuhjuWbgZ8AdKaWv9rOuVE0tkiTlJCJIKUUt1lXtMPty4KrK9FXAT3ovEBEB3Aw83F+QS5Kk\n4au2Zz4D+CHwKmAT8N6U0u6ImAN8M6X0toi4EPgtsBY41Nj1KaVf9FqXPXNJ0phRy555VWFeS4a5\nJGksKdMwuyRJKphhLklS5gxzSZIyZ5hLkpQ5w1ySpMwZ5pIkZc4wlyQpc4a5JEmZM8wlScqcYS5J\nUuYMc0mSMmeYS5KUOcNckqTMGeaSJGXOMJckKXOGuSRJmTPMJUnKnGEuSVLmDHNJkjJnmEuSlDnD\nXJKkzBnmkiRlzjCXJClzhrkkSZkzzCVJypxhLklS5gxzSZIyZ5hLkpQ5w1ySpMwZ5pIkZc4wlyQp\nc4a5JEmZM8wlScqcYS5JUuYMc0mSMmeYS5KUOcNckqTMGeaSJGXOMJckKXOGuSRJmTPMJUnKnGEu\nSVLmDHNJkjJnmEuSlDnDXJKkzBnmkiRlzjCXJClzhrkkSZkzzCVJypxhLklS5gxzSZIyZ5hLkpQ5\nw1ySpMwZ5pIkZc4wlyQpc4a5JEmZM8wlScqcYS5JUuaGHeYRMSMi7oyIxyLiVxEx7SjLNkbEAxFx\n+3DbkyRJfaumZ34dcGdK6dXA3ZXn/fkk8DCQqmhPkiT1oZowfydwS2X6FuBdfS0UEScCbwW+BUQV\n7UmSpD5UE+YzU0rbK9PbgZn9LPcV4FPAwSrakiRJ/Wg62osRcScwq4+XPtvzSUopRcQRQ+gR8Xbg\nmZTSAxHRNlAxS5cuPTzd1tZGW9uAb5EkKQvt7e20t7ePyLojpeHtxo6IR4C2lNK2iJgNrEgpndFr\nmf8BLAY6gfHAFODWlNIH+lhfGm4tkiTlJiJIKdVk93M1w+zLgasq01cBP+m9QErpMymluSmlk4H3\nAb/uK8glSdLwVRPm/wBcEhGPARdXnhMRcyLi//XzHrvekiTV2LCH2WvNYXZJ0lhSlmF2SZJUAoa5\nJEmZM8wlScqcYS5JUuYMc0mSMmeYS5KUOcNckqTMGeaSJGXOMJckKXOGuSRJmTPMJUnKnGEuSVLm\nDHNJkjJnmEuSlDnDXJKkzBnmkiRlzjCXJClzhrkkSZkzzCVJypxhLklS5gxzSZIyZ5hLkpQ5w1yS\npMwZ5pIkZc4wlyQpc4a5qtbe3l50CaqC2y9vbj+BYa4a8JdJ3tx+eXP7CQxzSZKyZ5hLkpS5SCkV\nXQMAEVGOQiRJqpOUUtRiPaUJc0mSNDwOs0uSlDnDXJKkzBUe5hGxKCIeiYjHI+LaoutR3yJiU0Ss\njYgHImJlZd6MiLgzIh6LiF9FxLQey19f2aaPRMRfFFf52BMR/xoR2yNiXY95Q95WEXFeRKyrvPa/\n6v05xqp+tt/SiNhc+f49EBGX9njN7VcSETE3IlZExEMR8fuI+ERl/sh//1JKhT2ARmADMA9oBtYA\n84usyUe/22ojMKPXvC8Bn65MXwv8Q2X6zMq2bK5s2w1AQ9GfYaw8gDcC5wDrhrmtDh1LsxJYWJn+\nObCo6M82Fh79bL/PA3/bx7JuvxI9gFnA6yrTk4BHgfn1+P4V3TNfCGxIKW1KKXUA3wcuK7gm9a/3\nUZfvBG6pTN8CvKsyfRnwvZRSR0ppE90/oAvrUqFIKd0D7Oo1eyjb6vyImA1MTimtrCz3f3q8RyOo\nn+0HR37/wO1XKimlbSmlNZXpPcB6oJU6fP+KDvNW4KkezzdX5ql8EnBXRKyOiKsr82amlLZXprcD\nMyvTc+jeloe4XYs31G3Ve/4W3IZF++uIeDAibu4xTOv2K6mImEf3CMv91OH7V3SYe15cPv40pXQO\ncCnw8Yh4Y88XU/dY0NG2p9u6JAaxrVQ+NwEnA68Dnga+XGw5OpqImATcCnwypfRCz9dG6vtXdJhv\nAeb2eD6XV/41opJIKT1d+XcHcBvdw+bbI2IWQGVY6JnK4r2364mVeSrOULbV5sr8E3vNdxsWJKX0\nTKoAvsXLu63cfiUTEc10B/m/pZR+Upk94t+/osN8NXBaRMyLiBbgCmB5wTWpl4iYEBGTK9MTgb8A\n1tG9ra6qLHYVcOgHdznwvohoiYiTgdPoPphDxRnStkopbQOej4jzIyKAxT3eozqrBMAh76b7+wdu\nv1Kp/F/fDDycUvpqj5dG/vtXgqP/LqX7iL8NwPVF1+Ojz210Mt1HXK4Bfn9oOwEzgLuAx4BfAdN6\nvOczlW36CPCWoj/DWHoA3wO2AgfoPiblg8PZVsB5dIfGBuBrRX+usfLoY/t9iO4DoNYCD1Z+qc90\n+5XvAVwIHKz8rnyg8lhUj++fl3OVJClzRQ+zS5KkKhnmkiRlzjCXJClzhrkkSZkzzCVJypxhLklS\n5gxzSZIyZ5hLkpS5/w/rS5hxx/t1FgAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "figsize(8,6)\n", "plot(castep.ion.get_array_core_radial_charge())\n", "plot(castep.ion.get_array_atomic_radial_charge())\n", "ylim(-0.5,0.5)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Visualise charge density isosurfaces on-the-fly" ] }, { "cell_type": "code", "execution_count": 76, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# grid points, in Angstrom\n", "real_grid = (castep.basis.get_array_r_real_grid()*\n", " castep.io.io_atomic_to_unit(1.0, 'ang'))\n", "resolution = [castep.basis.get_ngx(), \n", " castep.basis.get_ngy(), \n", " castep.basis.get_ngz()]\n", "origin = np.array([real_grid[i, :].min() for i in range(3)])\n", "extent = np.array([real_grid[i, :].max() for i in range(3)]) - origin\n", "\n", "# charge density resulting from SCF\n", "den = calc.model.den.real_charge.copy()\n", "den3 = (den.reshape(resolution, order='F') / \n", " castep.basis.get_total_fine_grid_points())\n", "\n", "# visualise system with isosurface of charge density at 0.002\n", "viewer = view(atoms)\n", "viewer.add_isosurface_grid_data(den3, origin, extent, resolution,\n", " isolevel=0.002, color=0x0000ff,\n", " style='solid')\n", "viewer" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Postprocessing/steering of running calculations\n", "\n", "- Connect Castep and Bader codes without writing any explicit interface or converter" ] }, { "cell_type": "code", "execution_count": 77, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "
IonChargeVolume
0.000.0378.00
1.000.06178.35
2.000.05141.46
3.000.07193.28
4.000.04106.92
5.000.06176.24
6.000.0376.36
7.000.04129.82
" ], "text/plain": [ "[['Ion', 'Charge', 'Volume'],\n", " ['0.00', '0.03', '78.00'],\n", " ['1.00', '0.06', '178.35'],\n", " ['2.00', '0.05', '141.46'],\n", " ['3.00', '0.07', '193.28'],\n", " ['4.00', '0.04', '106.92'],\n", " ['5.00', '0.06', '176.24'],\n", " ['6.00', '0.03', '76.36'],\n", " ['7.00', '0.04', '129.82']]" ] }, "execution_count": 77, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from display import ListTable\n", "from bader import bader\n", "bdr = bader(atoms, den3)\n", "\n", "rows = ListTable()\n", "rows.append(['{0}'.format(hd) for hd in ['Ion', 'Charge', 'Volume']])\n", "for i, (chg, vol) in enumerate(zip(bdr.ionchg, bdr.ionvol)):\n", " rows.append(['{0:.2f}'.format(d) for d in [i, chg, vol] ])\n", "rows" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "So far this is just analysis/post-processing, but could easily go beyond this and steer calculations based on results of e.g. Bader analysis." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Updating data inside a running Castep instance\n", "\n", "We can move the ions and continue the calculation without having to restart electronic minimisation from scratch (or do any I/O of `.check` files etc.). Here's how the core of the electronic minimisation is coded in Python, almost entirely calling auto-generated routines:\n", "\n", "```python\n", " new_cell = atoms_to_cell(atoms, kpts=self.kpts)\n", " castep.cell.copy(new_cell, self.current_cell)\n", " self.model.wvfn.have_beta_phi = False\n", " castep.wave.wave_sorthonormalise(self.model.wvfn)\n", " self.model.total_energy, self.model.converged = \\\n", " electronic_minimisation(self.model.wvfn, \n", " self.model.den, \n", " self.model.occ, \n", " self.model.eigenvalues, \n", " self.model.fermi_energy)\n", " self.results['energy'] = io_atomic_to_unit(self.model.total_energy, 'eV')\n", "```" ] }, { "cell_type": "code", "execution_count": 78, "metadata": { "collapsed": true }, "outputs": [], "source": [ "castep.wave.wave_orthogonalise?" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Example - geometry optimisation\n", "\n", "Use embedded Castep efficiently as a standard ASE calculator, giving access to all of the existing high-level algorithms: geometry optimisation, NEB, basin hopping, etc.\n", " - Compared to file-based interface, save overhead of restarting Castep for each call\n", " - Reuse electronic model from one ionic configuration to the next\n", " - Wavefunction and charge density extrapolation possible just as in MD" ] }, { "cell_type": "code", "execution_count": 80, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "LBFGS: 0 13:42:06 -401.372736 0.0903\n" ] } ], "source": [ "from ase.optimize import LBFGS\n", "atoms.rattle(0.01)\n", "opt = LBFGS(atoms)\n", "opt.run(fmax=0.1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Developing and testing new high-level algorithms\n", "\n", "Having a Python interface makes it quick to try out new high-level algorithms. \n", "\n", "- e.g. I'm working on a general-purpose preconditioner for geometry optimisation with Christoph Ortner (Warwick), let's try that with Castep\n", "- This was implemented in a general purpose Python code by Warwick summer student John Woolley, just plug in Castep and off we go!" ] }, { "cell_type": "code", "execution_count": 81, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "No preconditioner\n", "LBFGS: 0 13:43:33 -401.438808 0.2319\n", "LBFGS: 1 13:43:35 -401.442860 0.2484\n", "LBFGS: 2 13:43:38 -401.412385 0.0564\n", "LBFGS: 3 13:43:41 -401.411553 0.0407\n", "LBFGS: 4 13:43:42 -401.411726 0.0431\n", "LBFGS: 5 13:43:45 -401.410910 0.0318\n", "LBFGS: 6 13:43:47 -401.410539 0.0263\n", "LBFGS: 7 13:43:50 -401.410235 0.0194\n", "LBFGS: 8 13:43:52 -401.410229 0.0171\n", "LBFGS: 9 13:43:53 -401.410238 0.0178\n", "LBFGS: 10 13:43:55 -401.410182 0.0163\n", "LBFGS: 11 13:43:57 -401.410165 0.0156\n", "LBFGS: 12 13:43:59 -401.410104 0.0124\n", "LBFGS: 13 13:44:01 -401.410060 0.0096\n", "Exp preconditioner\n", "LBFGS: 0 13:44:15 -401.438808 0.2317\n", "LBFGS: 1 13:44:20 -401.412285 0.0498\n", "LBFGS: 2 13:44:23 -401.410115 0.0169\n", "LBFGS: 3 13:44:25 -401.410006 0.0030\n" ] } ], "source": [ "from ase.lattice import bulk\n", "import castep\n", "import preconpy.lbfgs as lbfgs\n", "import preconpy.precon as precon\n", "from preconpy.utils import LoggingCalculator\n", "\n", "atoms = bulk('Si', cubic=True)\n", "s = atoms.get_scaled_positions()\n", "s[:, 0] *= 0.98\n", "atoms.set_scaled_positions(s)\n", "initial_atoms = atoms\n", "log_calc = LoggingCalculator(None)\n", "\n", "for precon, label in zip([None, precon.Exp(A=3, use_pyamg=False)],\n", " ['No preconditioner', 'Exp preconditioner']):\n", " print label\n", " atoms = initial_atoms.copy()\n", " calc = castep.calculator.CastepCalculator(atoms=atoms)\n", " log_calc.calculator = calc\n", " log_calc.label = label\n", " \n", " atoms.set_calculator(log_calc) \n", " opt = lbfgs.LBFGS(atoms, \n", " precon=precon, \n", " use_line_search=False)\n", " opt.run(fmax=1e-2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Conclusions and Outlook\n", "\n", "- Scripting interfaces can be very useful for automating calculations, or connecting components in new ways\n", "- Can give lecacy C/Fortran code a new lease of life\n", "- Provides interactive environment for testing, debugging, development and visualisation\n", "- Appropriate mix of high- and low-level languages maximses overall efficiency\n", "- **CSC MSc [project](https://www2.warwick.ac.uk/fac/sci/csc/teaching/taughtdegrees/msc/projects/csc-msc-project-jameskermode-pythoncastep.pdf) available on extending Castep/Python interface**" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Links and References\n", "- `QUIP` developed with Gábor Csányi, Noam Bernstein, et al.\n", " - Code https://github.com/libAtoms/QUIP\n", " - Documentation http://libatoms.github.io/QUIP\n", "- `matscipy` https://github.com/libAtoms/matscipy, developed with Lars Pastewka, KIT\n", "- `f90wrap` https://github.com/jameskermode/f90wrap\n", "- `chemview` https://github.com/gabrielelanaro/chemview/ by Gabriele Lanaro\n", "- [RAL technical report](https://epubs.stfc.ac.uk/work/18048381) on Castep/Python interface:\n", " G Corbett, J Kermode, D Jochym and K Refson\n" ] } ], "metadata": { "kernelspec": { "display_name": "Python 2", "language": "python", "name": "python2" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 2 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython2", "version": "2.7.10" } }, "nbformat": 4, "nbformat_minor": 0 }