{ "cells": [ { "cell_type": "markdown", "id": "814578a2-71ad-4d3d-9004-221a0f66fe87", "metadata": {}, "source": [ "# RBFE calculations of a Protein-Membrane System\n", "\n", "This tutorial walks you through the process of setting up a relative binding free energy (RBFE) simulation campaign for a membrane protein system using OpenFE.\n", "\n", "RBFE calculations for membrane proteins require additional preparation compared to soluble protein systems. For soluble proteins, solvation and box setup are handled automatically within the OpenFE protocol. In contrast, membrane protein systems must be provided as fully built, solvated, and pre-equilibrated systems, including correctly defined box vectors.\n", "You can prepare these systems using tools such as `packmol-memgen`, `CHARMM-GUI`, or `Maestro`.\n" ] }, { "cell_type": "markdown", "id": "1a1c5af8-354d-49e3-a24a-91055b885add", "metadata": {}, "source": [ "
Note: This tutorial focuses on importing inputs prepared with Maestro. However, the overall procedure is applicable to systems generated with other preparation tools.\n", "
" ] }, { "cell_type": "markdown", "id": "b3ed2d3d-bad9-41fd-bf8a-10dadb3be6df", "metadata": {}, "source": [ "
Warning: The membrane must be aligned with the coordinate axes to ensure correct pressure coupling with a membrane barostat during the simulation.\n", "
" ] }, { "cell_type": "code", "execution_count": 1, "id": "5e8a1b0e-0c58-4743-abb2-c7b0a7646116", "metadata": { "scrolled": true }, "outputs": [], "source": [ "import numpy as np\n", "import openmm\n", "import openfe\n", "from openfe import ProteinMembraneComponent, ProteinComponent\n", "from openff.units import unit as offunit\n", "from rdkit import Chem\n", "\n", "import logging\n", "logger = logging.getLogger()\n", "logger.setLevel(logging.ERROR)" ] }, { "cell_type": "markdown", "id": "62ad3b9e-ee68-4735-ba77-8cf2aa585bbc", "metadata": {}, "source": [ "## (Optional) Combining System Components into a Single PDB File\n", "\n", "Maestro exports solvated protein-membrane systems in a non-standard PDB format that OpenMM cannot directly read. To simplify the workflow, it is often easier to save different parts of the system separately (e.g. the protein, lipids, and water) and then combine them into a single PDB file.\n", "\n", "Using OpenMM’s Modeller, the components can be merged as follows:" ] }, { "cell_type": "code", "execution_count": 2, "id": "e6270a97-1073-471e-9ebf-8744050b0c34", "metadata": { "scrolled": true }, "outputs": [], "source": [ "# Specifying input and output filenames\n", "pdb_lipids = 'a2a/a2a_lipids.pdb'\n", "pdb_protein = 'a2a/a2a_protein.pdb'\n", "pdb_water = 'a2a/a2a_water.pdb'\n", "pdb_complex = 'a2a/a2a_complex.pdb' # specify the output filename" ] }, { "cell_type": "markdown", "id": "7d8d621b-5179-4fdb-9e52-2d0f84e7fc9a", "metadata": {}, "source": [ "#### Step 1: Load each component\n", "\n", "We start by loading the protein, lipids, and water as `ProteinComponent` objects from their respective PDB files:" ] }, { "cell_type": "code", "execution_count": 3, "id": "3aec23b3-0120-4ec4-87c2-dd390f293406", "metadata": { "scrolled": true }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/Users/atravitz/micromamba/envs/openfe_env/lib/python3.13/site-packages/gufe/vendor/pdb_file/pdbstructure.py:465: UserWarning: WARNING: two consecutive residues with same number (ATOM 3205 N GLU A 219 9.932 6.197 -38.597 1.00 0.00 N , ATOM 3220 HH33 ACE A 219 11.901 4.558 -37.817 0.00 0.00 H )\n", " warnings.warn(\n" ] } ], "source": [ "protein = ProteinComponent.from_pdb_file(pdb_protein)\n", "lipids = ProteinComponent.from_pdb_file(pdb_lipids)\n", "water = ProteinComponent.from_pdb_file(pdb_water)" ] }, { "cell_type": "markdown", "id": "080472b2-3a5d-4cff-af6e-7f6cfc4d6e61", "metadata": {}, "source": [ "#### Step 2: Merge components using OpenMM’s Modeller and save the merged system to a single PDB file\n", "\n", "We create a Modeller with the protein first, then add lipids and water. This ensures all atoms are included in a single system. Finally, we write the combined system to a new PDB file that is ready for simulation." ] }, { "cell_type": "code", "execution_count": 4, "id": "80298646-6013-421e-9a7c-3dba1fd145d3", "metadata": {}, "outputs": [], "source": [ "# Create modeller using the protein inputs\n", "modeller = openmm.app.Modeller(protein.to_openmm_topology(), protein.to_openmm_positions())\n", "# Add the lipids\n", "modeller.add(lipids.to_openmm_topology(), lipids.to_openmm_positions())\n", "# Add the solvent\n", "modeller.add(water.to_openmm_topology(), water.to_openmm_positions())\n", "# Write out the new topology and positions\n", "mergedTopology = modeller.topology\n", "mergedPositions = modeller.positions\n", "with open(pdb_complex, \"w\") as f:\n", " openmm.app.pdbfile.PDBFile.writeFile(modeller.topology, modeller.positions, f)" ] }, { "cell_type": "markdown", "id": "47f2268f-826a-4b99-96fe-606b45a3b0a6", "metadata": {}, "source": [ "## Using the openfe CLI with membrane systems\n", "For basic membrane support, you may use the CLI command [openfe plan-rbfe-network](https://docs.openfree.energy/en/latest/reference/cli/plan_rbfe_network.html) with the ``--protein-membrane`` argument in place of the usual ``--protein`` argument.\n", "You will still need to have a fully built, solvated, and pre-equilibrated system that includes box vectors in a PDB or PDBx/mmCIF file when using the CLI.\n", "\n", "For the full functionality available in the openfe Python API, continue this tutorial." ] }, { "cell_type": "markdown", "id": "44aaba48-9664-4887-a6ec-f2946bda94b9", "metadata": {}, "source": [ "## Loading the system as a `ProteinMembraneComponent` with Periodic Box Vectors\n", "\n", "The **`ProteinMembraneComponent`** requires periodic box vectors to define the simulation box. There are several ways to provide these vectors:\n", "\n", "1. **`CRYST` record in the PDB file** \n", " If your PDB file includes a CRYST record (common in files exported from Maestro or other modeling tools), OpenMM can automatically read the box vectors from it.\n", "\n", "2. **Manually specifying box vectors** \n", " You can provide the box vectors explicitly in OpenMM format (see below).\n", "\n", "3. **Inferring from atomic positions** \n", " Box vectors can be estimated from the atomic coordinates in the PDB file.\n", " Caveat: This approach can be inaccurate if the PDB comes from a previous simulation and some atoms are positioned in periodic images. \n" ] }, { "cell_type": "markdown", "id": "52d7bb89-4b6e-473c-8be8-9cd555756c45", "metadata": {}, "source": [ "
Maestro-specific note: When merging multiple components exported from Maestro, CRYST records are not currently preserved in the combined output PDB file. If one of the original input files contains valid CRYST information (for example, a2a_water.pdb in this case), you can load that file directly as a ProteinMembraneComponent and extract the box vectors from it, as shown below.\n", "
" ] }, { "cell_type": "markdown", "id": "328fbbf9-72c7-4ed8-a304-c7623f0dd2c4", "metadata": {}, "source": [ "#### Option 1: `CRYST` record in the PDB file\n", "\n", "Box vectors are automatically recognized when stored in the `CRYST` record. In our case, only the `a2a_water.pdb` file had a `CRYST` record. " ] }, { "cell_type": "code", "execution_count": 5, "id": "e5adf727-158b-4567-bb66-d5d2909b729a", "metadata": {}, "outputs": [], "source": [ "a2a_cryst = \"a2a/a2a_water.pdb\"" ] }, { "cell_type": "code", "execution_count": 6, "id": "4e8d7607-8a47-4840-8000-d0f55a0248a1", "metadata": {}, "outputs": [], "source": [ "system_with_box = openfe.ProteinMembraneComponent.from_pdb_file(a2a_cryst)" ] }, { "cell_type": "code", "execution_count": 7, "id": "3ab4214a-070a-4bbe-b051-0451ec904ee6", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
Magnitude
[[6.9587 0.0 0.0] [0.0 5.9164 0.0] [0.0 0.0 9.2692]]
Unitsnanometer
" ], "text/latex": [ "$\\begin{pmatrix} & & \\\\ \n", " & & \\\\ \n", " & & \\end{pmatrix}\\ \\mathrm{nanometer}$" ], "text/plain": [ "" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "box_vector = system_with_box.box_vectors\n", "box_vector" ] }, { "cell_type": "markdown", "id": "c0d92085-1dfc-4307-911a-e2ce59116564", "metadata": {}, "source": [ "#### Option 2: Manually provide box vectors\n", "\n", "Some membrane setup tools (e.g. CHARMM-GUI or Maestro) store the periodic box vectors in the outputs of the membrane-building step or a prior equilibration simulation. In these cases, you can manually supply the box vectors to OpenFE.\n", "\n", "Box vectors must be provided as a NumPy array with OpenFF units. The vectors must be specified in their reduced form, as required by OpenMM.\n", "\n", "For details on reduced box vectors and periodic boundary conditions, see the OpenMM documentation:\n", "https://docs.openmm.org/latest/userguide/theory/05_other_features.html#periodic-boundary-conditions" ] }, { "cell_type": "code", "execution_count": 8, "id": "e6c47bf9-aee6-4f87-bca4-16ae8c0aad91", "metadata": {}, "outputs": [], "source": [ "box_vector = np.array([\n", " [6.9587, 0.0, 0.0],\n", " [0.0, 5.9164, 0.0],\n", " [0.0, 0.0, 9.2692]\n", "]) * offunit.nanometer" ] }, { "cell_type": "code", "execution_count": 9, "id": "911ee888-2563-4221-aecc-9a7b48f65e7e", "metadata": {}, "outputs": [], "source": [ "membrane_protein_user_box = openfe.ProteinMembraneComponent.from_pdb_file(pdb_complex, box_vectors=box_vector)" ] }, { "cell_type": "code", "execution_count": 10, "id": "fb648237-c527-4658-abb1-e7a6d0ed0df2", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
Magnitude
[[6.9587 0.0 0.0] [0.0 5.9164 0.0] [0.0 0.0 9.2692]]
Unitsnanometer
" ], "text/latex": [ "$\\begin{pmatrix} & & \\\\ \n", " & & \\\\ \n", " & & \\end{pmatrix}\\ \\mathrm{nanometer}$" ], "text/plain": [ "" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "membrane_protein_user_box.box_vectors" ] }, { "cell_type": "markdown", "id": "aca20b2a-5ecd-440f-8b98-7d3f92e096d2", "metadata": {}, "source": [ "#### Option 3: Infer box vectors from atomic positions\n", "\n", "As a fallback, box vectors can be estimated from the atomic coordinates in the PDB file. An orthorhombic box is constructed by taking the minimum and maximum coordinates along each axis and expanding the resulting bounding box by an optional padding." ] }, { "cell_type": "markdown", "id": "68e1135c-7640-4db8-b7eb-7a2666e64bd1", "metadata": {}, "source": [ "
Caveat: When merging multiple components exported from Maestro, CRYST records are not currently preserved in the combined output PDB file. If one of the original input files contains valid CRYST information (for example, a2a_water.pdb in this case), you can load that file directly as a ProteinMembraneComponent and extract the box vectors from it, as shown below.\n", "
" ] }, { "cell_type": "markdown", "id": "8532d0af-0d63-423e-b7ad-c4b7f10a4b6d", "metadata": {}, "source": [ "This approach is convenient when no box information is available, but should generally be used as a last resort." ] }, { "cell_type": "code", "execution_count": 11, "id": "10652b44-f348-449a-a37e-36ce49dc8b51", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/Users/atravitz/micromamba/envs/openfe_env/lib/python3.13/site-packages/gufe/components/solvatedpdbcomponent.py:230: UserWarning: Box vectors were inferred from the atomic coordinates.\n", "Note: This heuristic assumes that the coordinates reflect the true periodic unit cell. It may produce incorrect box dimensions for structures that were post-processed (e.g., unwrapped after MD).\n", "Inferred box vectors: [[9.445200000000002 0.0 0.0] [0.0 8.196100000000001 0.0] [0.0 0.0 9.8092]] nanometer\n", "Applied padding: 0.2 nanometer\n", " warnings.warn(\n" ] } ], "source": [ "membrane_protein_infer_box = openfe.ProteinMembraneComponent.from_pdb_file(pdb_complex, infer_box_vectors=True)" ] }, { "cell_type": "markdown", "id": "5826d248-4e03-45ed-ad57-b577fcfdba35", "metadata": {}, "source": [ "This box is much bigger than the other boxes since PBC is not explicitly handled in `infer_box_vector`." ] }, { "cell_type": "code", "execution_count": 12, "id": "389247f6-e411-4b7e-bfb6-a84a990bc2ed", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
Magnitude
[[9.445200000000002 0.0 0.0] [0.0 8.196100000000001 0.0] [0.0 0.0 9.8092]]
Unitsnanometer
" ], "text/latex": [ "$\\begin{pmatrix} & & \\\\ \n", " & & \\\\ \n", " & & \\end{pmatrix}\\ \\mathrm{nanometer}$" ], "text/plain": [ "" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "membrane_protein_infer_box.box_vectors" ] }, { "cell_type": "code", "execution_count": 13, "id": "1f93ac80-1718-4b84-951b-c44f629fe93d", "metadata": {}, "outputs": [ { "data": { "text/html": [ "519.1071333278928 gram/liter" ], "text/latex": [ "$519.1071333278928\\ \\frac{\\mathrm{gram}}{\\mathrm{liter}}$" ], "text/plain": [ "" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# The density is quite low for this case since the box is too big\n", "membrane_protein_infer_box.density" ] }, { "cell_type": "markdown", "id": "27178c69-b664-42ab-b3ee-1c83f8e8fe6a", "metadata": {}, "source": [ "## Validating Your `ProteinMembraneComponent`\n", "\n", "After adding periodic box vectors, it’s a good practice to validate your `ProteinMembraneComponent`. This helps catch common issues early, specifically:\n", "\n", "- **Number of waters:** \n", " Ensures that the system contains at least the specified `min_waters`. If there are too few, it may indicate missing solvent. \n", " **Default:** 50\n", "\n", "- **System density:** \n", " Checks that the density is above `min_density`. A very low density often signals incorrect box vectors or missing solvent. \n", " **Default:** 700 g/L\n", "\n", "If any of these checks fail, the `validate` method will raise a `ValueError`, indicating that the system may be missing solvent or have incorrect box vectors." ] }, { "cell_type": "code", "execution_count": 14, "id": "386937a7-4bb2-4880-8175-d601407ffed9", "metadata": {}, "outputs": [], "source": [ "# Our created `ProteinMembraneComponent` should not give any errors\n", "membrane_protein_user_box.validate()" ] }, { "cell_type": "code", "execution_count": 15, "id": "327e343d-4ee5-49e2-b4cd-b0a8618fc648", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "ProteinMembraneComponent validation failed:\n", "- Estimated system density is very low.\n", " Density: 519.107 gram / liter (expected ≥ 0.7 gram / milliliter). This usually indicates missing solvent or incorrect box vectors.\n", "This usually indicates missing solvent or incorrect box vectors.\n" ] } ], "source": [ "# The `ProteinMembraneComponent` with inferred box vectors should give an error\n", "from gufe.components.errors import ComponentValidationError\n", "try:\n", " membrane_protein_infer_box.validate()\n", "except ComponentValidationError as e:\n", " print(e)" ] }, { "cell_type": "code", "execution_count": 16, "id": "f8f3ead4-5fc2-41f7-b4bd-886fb21e5176", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "ProteinMembraneComponent validation failed:\n", "- Estimated system density is very low.\n", " Density: 141.802 gram / liter (expected ≥ 0.7 gram / milliliter). This usually indicates missing solvent or incorrect box vectors.\n", "- Only 0 water molecules detected (expected ≥ 50).\n", "This usually indicates missing solvent or incorrect box vectors.\n" ] } ], "source": [ "# The pure protein, without lipids and waters, should give errors\n", "only_protein = openfe.ProteinMembraneComponent.from_pdb_file(pdb_protein, box_vectors=box_vector)\n", "try:\n", " only_protein.validate()\n", "except ComponentValidationError as e:\n", " print(e)\n", "\n" ] }, { "cell_type": "markdown", "id": "6e4578df-6207-4fbd-baeb-7513e12adea3", "metadata": {}, "source": [ "## Loading the ligands and creating a `LigandNetwork`\n", "\n", "Next, we load the small molecules for which we want to calculate free energies. In this example, the ligands are stored in an SDF file containing multiple molecules. We can load them using RDKit’s SDMolSupplier and pass them to OpenFE:" ] }, { "cell_type": "code", "execution_count": 17, "id": "78a27762-f5eb-41ed-93af-3d0f17823171", "metadata": { "scrolled": true }, "outputs": [], "source": [ "ligands_sdf = Chem.SDMolSupplier(\"a2a/ligands_am1bcc.sdf\", removeHs=False)\n", "\n", "# Now pass these to form a list of Molecules\n", "ligands = [openfe.SmallMoleculeComponent(sdf) for sdf in ligands_sdf]" ] }, { "cell_type": "markdown", "id": "69ba8591-027d-4494-9147-956b44a2464c", "metadata": {}, "source": [ "The `LigandNetwork` is then created using three main components:\n", "\n", "- **Atom Mapper** – Proposes potential atom mappings between ligand pairs. \n", "- **Scorer** – Evaluates the quality of each proposed mapping. \n", "- **Network Planner** – Generates the network itself." ] }, { "cell_type": "code", "execution_count": 18, "id": "d2e64a3e-87f6-44d6-b2d9-92451082051a", "metadata": {}, "outputs": [], "source": [ "mapper = openfe.setup.KartografAtomMapper()\n", "scorer = openfe.lomap_scorers.default_lomap_score\n", "network_planner = openfe.ligand_network_planning.generate_lomap_network" ] }, { "cell_type": "code", "execution_count": 19, "id": "7dc15ad6-5b3a-4c10-aa9f-e369414a6c37", "metadata": { "scrolled": true }, "outputs": [], "source": [ "ligand_network = network_planner(\n", " ligands=ligands,\n", " mappers=[mapper],\n", " scorer=scorer\n", ")" ] }, { "cell_type": "code", "execution_count": 20, "id": "e621fb10-7630-44bf-9530-ec527e1f7c54", "metadata": {}, "outputs": [], "source": [ "# get the first edge; it automatically displays in a Jupyter notebook\n", "mapping = next(iter(ligand_network.edges))" ] }, { "cell_type": "code", "execution_count": 21, "id": "3dbcf3f3-76dd-4f07-a6ab-2855d2ff6e56", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAEsCAIAAACQX1rBAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydZ1wT2ffwT0IHRVSKiLKCKC42lLWi2EAB0bUsKgoWYLFjWw2uBTtBXQV7sKyIa0GxoIISRKWsDRtNUIo0UXoLNck8L+6z888vIIIkmUDu94WfzJ1yz2BOzi2n0AiCAAwGg8FgZBU61QJgMBgMBkMl2BBiMBgMRqbBhhCDwWAwMg02hBgMBoORabAhxGAwGIxMgw0hBoPBYGQabAgxGAwGI9NgQ4jBYDAYmQYbQgwGg8HINNgQYjAYDEamwYYQg8FgMDINNoQYDAaDkWmwIcRgMBiMTIMNIQaDwWBkGmwIMRgMBiPTYEOIwWAwGJkGG0IMBoPByDTYEGIwGAxGpsGGEIPBYDAyDTaEGAwGg5FpsCHEYDAYjEyDDSEGg8FgZBpsCDEYDAYj02BDKL1kZWUlJycDAEEQ8fHxX758oVoiDEa2+Pz5c0JCAvqclJSUk5NDrTwYMUEjCIJqGTCN8PbtWzabnZOTY21tXVRUpKSkdOnSpTNnznTt2pVq0TAYmSApKen69eu1tbUDBgxQUlKqra29deuWt7d3r169qBYNI2LwjFBKMTU13bhxo52dXW5urqOjo66urpaWlpqaGtVyYTCygomJyfbt22fOnJmTkzNr1qzevXt36dJFQ0ODarkwogcbQumlvLz86tWrDg4OAKChoSEvL19SUkK1UBiMDFFdXe3n5+fs7AwAnTp1UlZWLigooFoojOjBhlBKSUlJcXBwWLBgQUlJSVBQUE1NTXV1dVVVFdVyYTCyQmZmpr29/ezZsysrK+/cuVNaWlpfX19ZWUm1XBjRg/cIpZSkpKTExEQAMDY2NjAwiIiI6NOnj4mJCdVyYTCyQmpq6ps3bwDAwMCgf//+4eHhPXv2NDU1pVoujOjBhhCDwWAwMg1eGsVgMBiMTIMNobTz+vXr9+/fC7ZERkZmZ2dTJQ8GI2vEx8fHxcUJtjx9+jQ9PZ0qeTAiBxtCacfLy+v06dOCLcuWLQsLC6NKHgxG1jh69Ojhw4cFWzZu3BgUFESVPBiRgw0hBoPBYGQabAgxGAwGI9PIUy0A5vvExsbu37+fPCwqKqJQGAxGBklISBDUwdzcXAqFwYgcbAjbABwOJy8vjzzkcrkUCoPByCBVVVWCOlhXV0ehMBiRgw1hG2DcuHGHDh0iDx88eEChMBiMDDJ8+HBBf5mXL19SKAxG5OA9QgwGg8HINNgQYjAYDEamwUuj0o6pqamWlpZgy5gxY3r06EGVPBiMrNG/f38ejyfYMmLECAMDA6rkwYgcnGu0LcHlcmNiYuh0+tixY6mWBYORRXg83osXLzgcjqWlJdWyYEQGXhptS1y7dm38+PG7du2iWhAMRkYJDw8fPXo0g8GgWhCMKMEzwrZEcXGxjo4OnU4vLCzs2LEj1eJgMDJHbW2tpqYmh8PJzs7W09OjWhyMaMAzwrZEly5dRowYUVdXFx4eTrUsGIwsoqSkNGHCBIIg7t+/T7UsGJGBDWEbw9bWFgBCQkKoFgSDkVGwDrY/sCFsY0ydOhUA7t27J/1r2qmpqcXFxYItiYmJHA4Hfc7MzFy/fv2iRYsSExPPnj3r4eExc+bMt2/fUiEpBtMCpk6dSqPRwsLCamtrqZblO2RkZBQUFAi2JCcnl5eXo8/h4eGbNm1auHBhSUkJACQkJAwbNqy6upoCQakGG8I2xuDBg/X19fPy8qTfZvz2229Xr14VbBkyZEhsbCz6rKWl5eXltWrVqnv37rm4uDCZTE1NzUGDBlEhKQbTAnr27DlgwIDKysro6GiqZfkOLi4ufn5+gi3jxo2LiIhAn3/55Zf9+/ePGzfu+fPnfD7/6NGjI0aMkP4RtjjAhrDtYW1tDQD37t2jWpBWoaqqqqSkdOPGDRsbGwCIjY0dPnw4nY6/kJg2AFqYaeuroxoaGjweLzIycuTIkceOHVu0aJGioiLVQlED/t1pe7SPLQoej+fl5WVlZTVw4EAAOH/+vKOjI9VCYTDNAulgWx+MVlRUeHh4bNmyRUNDo6CgIDg4ODo6uq3/sPwYOLNM28PS0lJZWfn58+f5+fna2tpUi9MUr169Elwd5fP55OeIiIi0tLSysrLy8vIpU6YMGzZMRUWFChkxmBYzevTorl27pqSkpKamGhkZUS1OU8THxwvqoOC+5tWrV6uqqs6dOzdz5szdu3cDgL+//7Rp0yiQkmqwIWx7qKmpWVhYhIWFhYWFSfks6tWrV2gfHiFoCK2srKysrMjDRYsWSVQyDKYVyMnJWVlZXblyJSQkxN3dnWpxmiIhIaG+vp48rKmpIT+7uroKXSyzaogNYZvE1tY2LCwsJCREyg2hm5vb8uXLyUOZ3YHAtD9sbW3bhCF0cHDYsmULeaijo0OhMFIL3iNsk6Dli/v37+MivRgMJdjY2MjJyT1+/LiyspJqWTCtBRvCNomhoWHfvn1LSkqePXtGtSwYjCyiqak5bNiw2tpaMhoB03bBhrCtIrW+o3w+/9q1awBgYGDQuXNnwVP9+vVTVVWlSC4MRsRIrQ4SBBEYGEgQhL6+vqampuCpPn364DTFjUBg2iZsNhsABg0aRLUg/wOfz//9998B4M8//6RaFgxGvLx69QoAevTowefzqZblf9i4cSMALF++nGpB2gx4RthWsbCw6NixY1xcXGZmJtWy/B+bNm06ffq0ioqKoEcoBtMuGTJkiJ6eXk5OTnx8PNWy/B87duw4cOCAgoICmrBimgM2hG0VRUVFVBpUerLgb9my5eDBg4qKitevXx8/fjzV4mAw4oVGo0lbmicfH5+dO3fKyckFBATY2dlRLU6bARvCNoxUbVHs3bt33759cnJyFy9exENRjIwgVTp49uzZ9evX02g0Fos1d+5cqsVpS+DCvG2YvLw8PT09VVXVwsJCZWVlCiU5evSou7s7Goc6ODhQKAkGI0kqKys1NTW5XO7Xr1+7du1KoSQXLlxYsmQJQRAnTpxYtmwZhZK0RfCMsA2jq6tramrK4XCePHlCoRh///33mjVraDTayZMnsRXEyBQdOnQYM2YMj8cLCwujUIwbN264uLjw+Xxvb29sBX+AtmoIS0pKPD09V65cmZqaGh0dvXnz5v3798vg7JbylZnr16///vvvBEEcPHgQ+YtiZAQOh7Nr166VK1fGxcW9e/du48aN3t7eMljNjnIdfPDgwfz587lc7p49e5C/KKaltNWl0eLiYnl5+bi4uPj4eDqdPnLkyMGDB1MtFAU8ffp09OjRhoaGaWlpku/91q1b9vb2XC7Xy8vLw8ND8gJgKKS8vJzL5ebl5d2+fdvIyEhHR2fs2LEyWEgrJSWlX79+Xbp0yc/Pl5OTk3DvDx8+tLOzq6mpWbdu3aFDhyTce7uhrX5ru3TpkpycvH379kmTJhkbG7PZbAcHB6FazLLAiBEjtLW109PTU1JSJNw1m82eN28el8v19PTEVlAGUVdXz83NXbt2rY2NzU8//fTixQtHR0epCuaRDMbGxkZGRsXFxS9evJBw10+fPp0xY0ZNTc2qVauwFWwVlEYx/jiVlZX19fUZGRlbt27lcDgEQZw9ezYyMpJquShgwYIFALB7925JdhodHa2mpgYAa9askWS/GOmhpqampqamtLR05cqVSAdv3rx569YtquWigNWrVwPA+vXrJdnpmzdvUOamxYsX83g8SXbd/mir1Se+fv3KYrFqampWrFgRFhb27NkzRUVFJycnquWigM6dO3fq1Gnbtm1///23paWlubn5+PHj9fX1xdfjs2fPbGxsOByOs7Pz4cOHxdcRRpopLi728fGpqalxcXGJiYlhs9ny8vLbtm2jWi4K0NbW7tix4+HDh2/dumVubj5mzBgrKysDAwPx9RgfH29paVlSUvLbb7+dOXNGBlekRUtb3SOUEerr64OCgkaNGvXTTz/dvHkzJSXF1dVVMHngrVu35syZU19fr6SkJFhys2/fvhYWFhYWFuPGjROtUXz37t3EiROLi4udnJzOnz+PNRAj40RERNjZ2VVXVysrKwtW+zM0NBw7duz48eMtLCwMDQ1F2OPHjx8tLCy+fPny66+/Xrt2TUFBQYQPl1GonpJimuLFixeenp5sNjs/P//OnTvJyckrV64kz7LZbBQ+uHXr1vr6+tjYWB8fH3t7ew0NDcH/Yl1dXXt7exaLlZCQ0MqkiMnJyaie2axZs+rr61v9fhhM2+bp06coh7WLi4ugDnbp0kVQB7t162Zvb+/j4xMbG9tKHczMzPzpp58AwMrKqrq6WlQvIuNgQyjtBAYGstls9Dk+Pt7T0xN9jomJQbt07u7uQrdwudyEhAQWi2Vvby8U5Nsahfz48WP37t0BYMqUKTU1Na1+MwymbfP27Vu0S7dw4UKhXToej0fqoFD9B21tbTs7OyaTGRsb29K9vZycHDS5NDc3r6ysFOnbyDRt3hBu2bJl27Ztgi1z58719/enSh6RQxrCnJycVatWVVRUEATx7NkzNA5dsmTJd+1ZWloai8VycnLq0aOHoEJqaWk1XyGzsrJ69eoFAJMmTcLjUIwg3t7eQj5Tv//++9GjR6mSRzK0aHWE1EGhfYqOHTtaWloymcyoqKja2tqmH5Kfn//zzz8DwIgRI8rLy0X3Kpi2bwiXLFni6uoq2DJ8+PDDhw9TJY9oSUlJWbBggaur64sXL/r06bNz585Lly69e/cOLbw4OjoiA8bn848ePdqcEWJaWpq/v7+bmxtaXWmmQn758sXY2BgARo0ahSwxBkOyfv362bNnC7ZMmTJl+/btVMkjAT5+/KirqwsA1tbWgqsjXl5e48eP9/T0fPjwIfKkbQipg0LeNGpqapaWlmgrpOGKS0lJyZAhQwBg0KBBRUVFYnw3mQQbQoLD4bBYrJMnT3K5XIIguFzuoUOHSkpKRCyoiEhJSenWrRsAzJgxgxyHbtiwAe0ZtOhRubm5gYGBbm5uJiYm31LI6urqgoKC/v37A4CpqWlxcbEY3gnTthGJIQwJCdm7d29ubm5UVBSTyWQymTExMSIVU2Q0sToybtw4Uo/k5eXNzMzc3d0DAwO/pTiCOkij0ch7VVVVzc3NGQwGm82uqqoqKysbNmwYAPTt2/fLly8SeUvZop0YQp4ALTWE//77b1JS0qlTp1BN58OHD8+ZMyctLU1sIv84qampaJdu8uTJ5JgROawrKirevXv3h5/86dOnCxcuuLi49O3bV9AoqqqqamtrA8DAgQMLCwtF9B6YdgUyhII62FJDWF5efu3ataysrEWLFlVVVRUXF588eTIiIkJ8Mv8wTa+OlJeXs9lsBoNhbm4u6MwpJydnYmLi5uYWGBj4LT3Kzc29dOnSsmXLhAamSkpKaA3WyMgoNzdX/K8oi7QHQwgN+IGl0WPHjkVERGRkZGzbts3Ly4s0hJs2bZo5c6bglXZ2dlu2bBGN9C0hOzsbraVMnDixqqoKNaIwPjk5uatXr4qqoy9fvgQHBzMYDDMzMzqdrqen16VLFzwOxXyL9evXN9TBH1gazc7OJvcaly1bRrYzmcxx48YJXunk5LRixYrWSf0j5OfnIyvV6OpIQUHBp0+fCIIoLi4ODw/Py8tjs9menp6WlpZKSkqCfxxDQ0M3Nzd/f//s7OxvdSSog926dVNXV0cPx4iDthpQL4iLiwuLxSIPR48e3dInREdHFxUVTZgwYdasWd27d4+MjOzYsePKlSvhv/ASwYsbtkiA/Px8KyurjIyMkSNH3rp1S0VFBQCOHz++bt06Op3u7+8/Z86clj7z8ePHhw8f9vT0lJOTu3z5ck1NjZyc3F9//aWjozNt2rRp06YBQE5OTt++fWtqang8nujfCtNemD179tWrV8nDqVOntvQJhYWFTCZz586dAPDo0aOxY8eSp6REB8vKymxsbJKSkgYOHBgeHo78RUmysrJYLJaiouL27dsPHTpkbm6+Zs2aq1evourZVVVVr1+/jomJCQ8Pj46OTk9P9/Pz8/PzAwBDQ0MUgz9lyhRy515LS4vUwaKiImNj4/LycsEgRYxoaQ/R0DQaTU6Alt4eFBS0atUqJSWlmJiYGzduHDt2zNra2sbGRhyi/hiFhYWTJk1KTk4ePHhwSEgI8hf19/d3d3en0WgnTpxAWdZayvjx42fMmFFXVzd48GAmk2liYtLQmvbo0WPixIkEQTx48EAEb4Jpv7RGB3Nzc8eOHdulS5fg4GAAuHLlir29vRhk/HE4HI6dnd2rV6/69OkTFhbWsPSgvr4+GjrTaLTdu3dbW1sL7fmNGTMG7flVVFTExsYymUw7O7tOnTqlp6cHBAQsXbq0V69e3bt3nzNnjp+fX2JiInlv165d0c/RvXv3JPKuskh7mBG2Eltb24kTJwIAmmYBwNq1axUVFckLKioqkpKSyMPKykpJiofGoQkJCcbGxg8ePEDj0KCgIFdXVz6ff+DAgaVLl7a+F4Ig3rx54+bm1vCUra3tvXv3QkJCGl2FxmBaj46Ozr///gsAyIgePnxYKFtKVVWVoA6WlZWh4aBkqK6utrOzi46O1tfXZ7PZyFutaQ4ePOjq6troKeREY2ZmxmAw6uvrX758GRkZGRkZGRMTk5eXd+3atWvXrgGAvr5+ZGQkmiPa2tpevHgxJCSk0VVoTOvBhhBUVFRIE4gQqvb+7NmzX3/9lTz8/Pmzubm5ZGSrqqqaNm1abGyskZHRo0eP0J55cHCwg4MDl8vdu3fvH3/8IZKO7t27JzgJXr169b179169etW5c2c7O7uVK1eGhYXV1dUJjg8wGFEhLy8vuNKoqqoqdEFiYqKgDn79+tXR0VEystXV1f3222+PHz/W09N79OiRUNyR4GVRUVHp6enZ2dk7duz46aefmpN9UEFBYfTo0aNHj/bw8ODxeMnJyWj5NCIioqysjAz8tba2lpeXj4yMLCsr69SpkyhfD4OQ9KakqMnJycnJyRFs+fDhgwj9Gzdu3DhjxgzBlqlTp/7555/JyclMJlOsTlxVVVUTJkwAgJ49e2ZkZKBGMq1a6x12Xr16tW7dug0bNpSVlV28eFEwph71e+XKFXSIwiek04sPQzl5eXmZmZmCLWlpaV+/fhXV8728vCwsLARbHB0dly9fnpOTs3PnTrG6kHC5XLRIq6WllZSU1MSVtbW1bDabzWZnZmaiDw8fPmxRX6mpqSiENzMzMzY2Vui9xowZAwBBQUE/8BaY79LmDaEgL168GDly5IIFC0T4zG8ZQjQVo9PplpaW/v7+3wqe/WHq6urs7OwAQEdHJzk5GTXGxMR06NABAFavXi3a7oQ4cOAAACxcuBAdorLXf/zxh1g7xbQDEhISLCwsbG1tRfjMbxnCffv2IR00NzdnsVgiTznG4/Hmz58PABoaGq9fvxbtw4W4ffu2lZVVZmZmamrq4cOHd+3adeHCBcEL0Mu6uLiIVQyZpT04y5B07dr12bNnISEhXC5X3H3Z2NjY29srKCiEh4cvWrSoZ8+e7u7ur1+/FsnDeTyek5PT3bt3tbS0IiIiUNzSixcvrK2tKysrlyxZ4uvrK5KOvgXy+gsNDeXz+eQh3qvHfBcdHZ2YmJiHDx9KYCt9/Pjxjo6OyM1t6dKlenp6bm5uT58+FcnDCYJYvnz5pUuX1NXVw8LCUFYX8TF9+nQUjN+7d+8ZM2ZUVVX16dNH8AKkgyEhIQSuFyQG2pUhNDQ07Nu3b0lJybNnz8Td18SJEwMDA798+cJisczNzYuLi48ePWpmZta/f39vb++vX7/+8JMJgli6dOnVq1c7deoUGhqK4pbi4uJsbGwqKirs7e1Pnz4t6JAmDn7++efevXsXFBS8fPkSAMzNzTU0NN6/f5+WlibWfjFtHU1NzWHDhtXW1kZERIi7r1GjRgUEBHz58sXf39/S0rK8vPz06dOjR4/u16/fjh07MjMzW/PwjRs3+vn5qaqq3rlzB2V1kRgqKioGBgYfP34UbBw0aJC+vn5eXt6bN28kKYysQPGMVNSsW7cOADZv3iyqB/J4PJR6jYTL5TZMUZ2YmMhgMFASFgCQk5OztLQMDAysq6trUXd8Pn/ZsmUAoKamFhUVhRobTasmbpAvOBkWjXZKjh07JpneMW2XXbt2AcDSpUtF9cBm6mBycrKnpyfpzNKabQsPDw8AUFRURDMwCZCXl7dixYq7d+8mJCRcv37dz8/P19dX6Brk1L17927JiCRTtDdDyGazAWDQoEGU9M7lctlsNloyRdrYpUsXNze3N2/eNPMJaDdORUXl0aNHqCUzMxNlrLeyspJk8aOQkBAAMDMzQ4fnz58HABsbG4kJgGmjvHr1CgB69OjRysJ7PwaPx2Oz2U5OTqQruIaGhpubGzms/C4oqF9BQSE4OFisogqSkZGBXGwIgoiOjo6IiGj417t9+zYAjBw5UmJSyQ7tzRDW1taiACNq0xEVFRWxWCxTU1Ny5m1iYsJkMgsKCpq4688//0Tj0Hv37qEWMq2a5MuPVVdXq6qq0mg05Bmbn59Pp9OVlZVF7haEaWfw+XyUEffdu3cUilFSUoK2LYR0sGl3Vh8fH7Sic/nyZYmJ2kwqKyuVlZXpdLoIPXIxiPZmCAmCmDlzJgCcOnWKakEIgiASEhIYDAZZmVNJScne3j44OLjhCudff/2FNBDl/iYI4uvXr6j82MiRIykpP4b258+dO4cO0U7JnTt3JC8Jpm3h7OwMAF5eXlQLQhAt2bY4e/YsjUaj0WinT5+mRNTvMnnyZAAQcijFtJ52aAhPnz4NANOnT6dakP+jpqYmMDDQzs5OXv7/ZzDo3r07g8H48OEDuuDIkSNoV+Off/5BLWT5scGDB1NVfuzEiRMAQFbY2bFjBwAsX76cEmEwbYigoCAAGDNmDNWC/B/f3ba4cOECnU5HOQupFbUJkLv4vHnzqBakvdEODeHnz59pNJqampqYCqmfPXt2zpw5gi0eHh7NDG/Pysras2cP6RhNo9HGjBnj4uKCxqEsFgtdVlZW9ssvvwCAsbExhWUfsrKyAEBdXR3F+SIP0p49e1IlD6atUFFRoaSkJCcnJ6bSXdevX7e2thZs8fLyWrVqVXPuzcvLO3DggGCpo+HDh7u5uaFBqre3tzgEFhXIbbtz584Sc5qTEdqhISQIAs2l7t+/L46HHzx4cNSoUYItDYsDf5fY2Fg3Nze0nYn+PXr0KDrF4XAsLCxAOsqPDRgwAABQjgw+n4+qcsfHx1MrFUb6mTRpEgBcunRJHA8/e/Zsv379BFsaFgf+LrGxse7u7mjbAnnW7Ny5U6RiigUUVRwZGUm1IO2KdhVHSEIGn1ItyDcxMzNjsVh5eXnHjx/ncDgKCgpkSut//vkHJdt9+PAhcjqgEFtbW/jvL0mj0aZMmQLS/YfFSAmC3xzpxMzMzNfXNzs729/fn8vl0mg0waTz5eXldXV1FIr3LaT/D9sWaZ+GEH1X7t69K6bnEwRRLwDxv7keCgsLr1+/TsbzZmRkfCvbhZqa2ooVK4YPH15fX//w4UPU+Pvvvx88eDA8PBxFTVCLUE4ZnGIG00zIwaj4KlkK6iBKgURSWlp648YNMiY9Nzf3yZMnjT5EWVl54cKFVlZWBEHcv38fNa5atUpLSys8PFxMkrcGrIPioH0awhEjRmhra6enp6ekpIjj+c+ePVMUAMXYkdy6dat3795r1qwBAC6Xy2Qyz50718TTGn6zN2zYYGRkJA7JW8ro0aM1NDSSk5PR5sTkyZMVFBRiYmJQee6QkBAmkwkAFRUVe/bs8fDwwLVDMQhjY2MjI6Pi4uIXL16I4/nJycmCOojCHkhu3bplYGCA0msQBLF79+4W6WDXrl3r6uqkc9ZlYWHRqVOn+Pj4VqbOwQjSPg0hnU5HfsaiHTcVFRVVVFQAQMM9QsHLXF1dw8PDBw8eDAC+vr7Ozs5Nlyol1zoI6csiKC8vj/6SoaGhAKCurm5ubs7j8dBguVu3bhkZGQCwd+/eefPmMZlMoQpWGFkGFfYSrTkpLy8vKSkBgIZ7hIKXLV68+N9//0XRR+fOnbO3t2/6m4kS3IeFhdXX18N/dvHOnTsilFxUKCgooPqp5PwV03rapyEEUa+kp6SkrFmz5qeffjp58mRzrp87d25GRsanT5+CgoLevXv3/v37/Pz8b108ZMgQPT29nJyc+Ph41FJeXo6WeqqrqwsKCkTyCj/Mt1ZHaTTa0KFDUdbT9+/fX7p0aeHChUVFRRSKipEqRLuIl5qa6uHh8dNPP6HSKN9l5syZBQUFeXl5fn5+aWlp79+/z87O/tbF+vr6JiYm5eXlMTExAPDLL7/o6OhkZWUJVgOWHvDqqMhpt4ZQsJTlDz+Ex+PdvHlz4sSJ/fr1O3LkSHV1dXOyTh88ePDly5eVlZU6OjoXL160tLTU09NropwmjUaztraG/77ZwcHBU6dOzcnJSUtL2717965du6gdmdrY2NDp9MePH1dVVcF/IwyyMAWiS5cuW7ZsmTx5MgqxwGAAYPz48R06dHj79m1ubu4PP4TP54eGhtrY2PTt29fb27usrKw5OnjkyJGnT5+WlJR06tTp8uXLSAe7du3axC2C1oVOpwuqJIWUl5fv378fJVcDgJcvX65fv97W1pZGo4WHh1dXV1MrXruh3RrCzp07jxw5UtAJpUWUl5f7+voaGRnNmjXr0aNHHTt2dHNzi4+PZ7FY3713+fLlhoaGFy9eVFFRMTQ0NDQ0ZDKZSkpKTdwiOH+dPn06UkK00aivr9/0veJGS0vrl19+qampQfUETExMyMIUvr6+aWlply5dcnFxWbNmzdOnT1H5UAwGAJSUlCZMmCDohNIiKioq/Pz8Bg4caGtre//+fUVFRScnp3fv3l29evW797q5uRkaGv7zzz+qqqqkDjYsfC+I0BqSlDhncrncCRMmoIkpl8s9c+ZMeXm5rq6uqalpdRGVda0AACAASURBVHX1tzyAMC1GQmEaVPBjpSxTUlLc3d3V1NTQ36d3795MJhP5hiCKiorIevGIz58/f/78uTWiCgUg79mzB5X8zsjI2LNnD+VZzVAa4mXLlqFDocIUGMy3OHXqFADMnDmzRXelpqYyGIzOnTsjHezevbunp6dgbH5JSUlqaqrgLV++fMnOzm5RLyhVRWFhYVpaWlpaWmFhoYaGBgCgJ5eWliooKMjLy5eUlLTosSInIyNj3759BEHs27fv9evXbm5uBEFs3boVxF+gW3Zoz4bw3bt3AKCrq9ucLPgoab2dnR1Z6s/c3DwwMFCo/ov4IAOQS0pKNmzYEBUVlZSU9OLFi7t37yI1oBChnDJChSkwmG+RnZ1No9E6dOjQzMIpUVFR9vb2pHMZ0kFxZFG5ePHi0KFD6+vr3717FxgY+Mcff1y5ckWo1hiqlEvm/qUKZAhrampsbW0ZDMbgwYOfPXuGIrIMDAyola3d0J4NIUEQqDjZq1evmrimrKyMxWL169cP6V6HDh3c3NwSEhIkJiQCJd12dHR88+YNi8VisVhcLvfSpUv+/v4tLWoocsh6AiinjFBhCgymCQYOHAgA4eHhTVxTUVHBYrH69++PdFBJSQmtgopVsJUrV5ImdtmyZbW1tX///TcA2NraokZvb28AWLx4sVjFaJry8vJt27YtXryY/EXas2cPQRA8Hg+lEX///j2F4rUb2rkhXLp0KQDs2rWr0bMfPnxwd3fv0KEDUj+0kUBVhuvk5GQA6NKli8TmoC0ChYgwmUx0KFSYAoP5FqjI7bp16xo9m5aWxmAwunTpgnRQV1fX09Oz6WplooI0hCkpKQcPHiQI4uvXr3Q6XUVFBdUaS0hIAABtbe2GRYClAUdHRwD466+/qBakPdDODWGjpSwpXwX9FiiI/t9//6VWjEa5du0aAFhYWKBDocIUGMy3iIyMBABjY2OhdqFVUDMzM39/f4nlki4sLHR0dExLSyMIYsOGDeRGIKo1dvfuXXTYq1cvAHjx4oVkpGoRly9fBoCJEydSLUh7oJ0bQqFSluXl5SwWC4XZAoCysrKTk5P0pJB2d3cHgGYWspAwZWVlioqKcnJyaMYsVJgCg/kWXC4XxS18/PiRIIjq6mp/f3+0Xgr/Veh8+vSphKV68+ZNYGAgqkH/8uVLsl2o1tiKFSsAwNPTU8LifQvBuWlxcbG8vLyCgkJpaSmFIrUP2rkhJAgC5Yk+cOAAg8FAXmFok5nCVdBvgbzMhwwZQrUgjePs7Lx27dq8vDx0iHZ0IiIiqJUKI/3MmzcP7VB4enqSwXzdunXz9PTMz8+nWrr/ASWEI/3CUL7iYcOGUSsVorS0dNSoUdevXydbULSSYAvmx2jnhjAtLe23334TDBeZNGnSrVu3pHPRv6ampkOHDjQaLScnh2pZvs/UqVNVVVXV1dXt7OyYTGZsbKx0/lUx1JKVlYU2mAV3Iq5evSqdFfV4PF63bt0AADmnVFVVqaqq0ul0cvxHFeXl5SNGjACAAQMGkH86R0dHZWVlZWVlS0tLJpMZFRVFuWNdG6UdGkIulxsVFcVgMMzMzJDiKSgo0Ol0R0fHuLg4qqX7DtOnTweA06dPUy3Id7hy5YqcnBxZ7BvRpUuX6dOn//XXXy9fvqR8wxVDIVwuNzY21tPT08zMDNk/OTk5Op0+c+ZM6dwCF2Tx4sUgUKEXRdb//fffFIpUVVU1fvx4ANDX1//06RNqDAkJQbsVgjqorq4+depUb2/vp0+fSudQQzppP4YwLy/v7Nmzs2fPVldXJ78WGhoac+bMQeVk20Qpyx8LQJYwt2/fRiZw7969aWlp/v7+KJGHoEKqqalZWlp6enqy2exmxpBh2joFBQUBAQHz5s0jY+EBoGPHjrNmzUJfj9u3b1Mt4/cJDAwEgHHjxqHDY8eOAYC9vT1V8tTW1iInbT09PeTdQxBEVFQUSvqxZs2a3NzcwMBANzc3ExMTQR1UVVU1NzdnMBhsNru6upoq+dsEbdsQ8ni82NhYJpNpbm5Op/9fujhDQ0N3d3c2m41cOVBmeg8PD6rl/T4tDUCWPGw2GyXyX79+vdApQYUk18GEFLKqqqrp51dWVgousfL5/IqKCtG/BkZ0JCQkMJlMS0tLeXl5QR10c3MLDg5G3+Rdu3YBwNKlS6kW9vuUlZUJ5pRB1Y7U1dUpWXXkcrkozF9LSyspKQk1Pn36tGPHjgDg7OwslC0kLy8vMDDQ3d2dnIsj5OXlzczMGAxGcHDwd51rOByO0IpORUVFc9KStF3apCGsrKwMDg52c3MTLOCuoqJiaWnp4+NDLh2QoJpBAwcOpETaljJo0CAAYLPZVAvSCDExMSjs8ru5nX5YIRUUFB4/fkwefvjwAQAkE1iGaT4cDofNZru7u/fo0UPwP9fc3JzJZCYmJgpd/+rVKwDo0aNHm/g9ReuQV69eRYdopvXo0SMJi8Hj8ebPn49Wtl6/fo0a3759iybcCxcubHpX/uvXr8HBwWiTSHCegHTQ3d09MDBQMHkkiY6Ozo0bN8jDwsJCAEhJSRHt20kVFBvCvLy8y5cvC7YkJyeHhoY2enFaWpqPj4+lpaWioiL5n9qrVy83N7fAwMAm5g11dXWo+ENDGymFbN68Gb4dgEwhr1+/Rm63ixcvbtHPWX5+fqMKKScnRyok6cGLDaGEKS4uDggIEPxJTU9Pv3XrVqMXp6WlsVgsOzs7wUTw2traTk5OgYGBZWVl3+qFz+fr6ekBgLjzxYiE/fv3A8CiRYvQ4caNGwFg48aNkpSBz+e7ubmhySgZyJicnKyjowMAs2bNatEWYHl5OZvNZjAY5ubmgr+fcnJyJiYm6CeUzOaKDaGkYbPZqqqqgi0+Pj6CzspVVVVo4ImSpZH/eWjgGRsb28yOZs2aBQAnT54UpfTiISoqCgD69u1LtSD/Q1xcHHJ8/+2331rjCFNcXBwcHLxhw4Zhw4YJrqTJyckNHTo0NzcXG0IJg+ZqgkvxFy5c6N27N3lYX1+PvM/IAFxyEOPp6RkbG9vMUZGLiwsAUJ44tzkkJiai1Ug0Pnj06BEAmJiYSFKGDRs2oG0FUh0+fvyI1sCmTJnSmq2TsrKye/fueXh4jB49WtDfjU6nDxo06MOHD9gQSppvGcK8vLxjx45ZW1sL1pXW0dFZsmTJtWvXfiCA9MyZMwAwbdo00ckuLsgA5A8fPlAty//nw4cPyKd8xowZInRFq6ysZLPZnp6elpaWysrKnTp14nK5CgoKp06dev4f169fx4ZQrHzLEBYWFp46dWr69OlkJRYA0NTUXLBgwaVLl34gBjcoKAgAxowZI1LxxUXv3r0B4Pnz5wRB1NfXCxamkABoWUhRUTEkJAS1ZGVloTQ3kyZNEqHnC4fDiYqKQlu8KioqCgoKHA5HR0dn//79pA6y2WxsCMXLtwwh+tMjTExMGAxGVFTUd8PUHj16tHnz5r179xIE8eDBAwaDQSbD/Pz5M41GU1FR+a6zhjTg4OAAAD4+PlQLQhAEkZmZiabjVlZW4nPh4XA4b9++JQhCQUHByMjI9D/QLETQEMbFxZEJiBMTE2NiYnD8Ymv4liFE7YI6yGazv+sw8uzZsz///HPr1q08Hi86Onrjxo3Hjx9Hp4RqjUk5qNYYmVNGqDCFWEGORQoKCijxDUEQX758MTY2BoBRo0aJz3esuroa7UTq6OgYGBiQOojSAAkawpSUFLQaV1BQwGaz2Wx2m9j6bQLqDaGCgoKHAJaWlsOGDautrXV0dDx//jxKjdZM3rx5QxDEunXrMjMz582bRxDEmjVrUNUxgiCGDh0KAN/agJQqAgICAGDy5MlUC0JkZ2cbGBgAgLm5eWVlpQR6bHpp9MmTJ+fPn1+2bFl0dPTZs2d9fHxCQkKwIWwNyOBt3LiR1MFp06b17t2bz+cvXrz49OnTLaox8vbtWx6Pt2vXrjdv3syfP5/H423fvp1c2yBrjYnnVUQJqk3/yy+/oEOhwhTiw9fXF608k84TBQUFKIuTqalpo74tIqfppdHXr1/7+Phs27bt6tWrFy9ePH78eNOlRdoEkqhQz+VyP3/+jD6Xl5eXlpYKXaAiAFqzVlRUDAgIWLRoESo10kxMTU15PN7nz581NDSQe3Hfvn0/ffqEzkpJyenmYGtrKycn9+TJk4qKCgrFyM/Pnzx5ckZGxogRI0JDQwWXyKjCwsJi0aJFY8aMKSoqCgsL4/P5RUVFgg44mEbh8Xi5ubnoc2VlZVFRkdAFDXWQRqP9/fffrq6ugr7Z32Xw4ME0Gu3jx4+Ghoby8vJ0Ov3nn39OTU1FZ9uQDk6YMEFVVfXVq1fot8vW1pZOpz969Kiqqkp8nZ47d27t2rU0Gu3UqVMoL11ZWZm1tXViYuLAgQPDw8MFAzSpYsiQIWvWrJk0aVJhYeHXr1+zsrLIn/c2jLgtbVlZ2YYNG3bs2HHo0KEnT554eXk5ODigRTCiGc4yLaK6uprBYKB1M0dHR4IgPDw8yLrVbauU5ahRowBA0HkvNjY2ICAgKytLMgKUlJQMGTIEAAYPHizJpKzfdZbJzc11c3PjcrkzZswoKyvbtm1b832mZJOampo//vjD29t7+/btr1+/3rlzp7Oz85MnT9DZ7zrLtIi6urrt27cjL8cFCxYQBMFkMsm89lJea0wIoVpjQoUpRM6FCxfodDqNRiMXkysrK1E20T59+nz+/FlM/Tbku84ypaWlLi4uFRUVdXV1PB7P29tb+hMGNY3Yh9Lq6urbtm3r3bu3qqqqhYWFg4NDz549Bf13RYinp2dVVdXdu3ezs7OnTJnCYDA6d+5MxjkNHz5cW1s7IyMDaaOU03Ds7O/v7+TkpK+v37179zlz5vj5+aFoLXH0Xl5ePnny5Ddv3hgbGz948IAsFycBlJSUBGd4dDpdUVGRjER89uyZo6Pjr7/+mp2d3b9///Dw8KysLJQ5CPMtlJSUPD09jYyMOnbsOGTIEBcXl549ewq6oYmQ/fv3f/36NSIi4uPHj/b29n/88Ud1dfWAAQPQWWNjYyMjo+LiYpTbWsoR0sF169YdO3aMTNwoWm7evImi4728vFDJi+rqajs7u+joaH19fTabLckvOdrKJQ9pNJqioiKplcnJyTNmzJg6dWpOTs6dO3fCw8M/fPiAYmPaMBIwtgUFBcePH0ehC/Hx8Rs3bkSbeYSoZ4TfxcnJCQBQEU5pwN/ff8uWLTt27KiqqvL29t65cyd56vXr1wDQvXt3chc6ICDAzs6OLKCB0NPTmz9//qlTp8isE62Hw+FYWFgAgJGRkbTVoE9OTkab8xkZGTwe7/HjxxKbIrdpSktLT58+jb75KSkpW7duJTMOinZG+F1QrbGtW7eK6fkt5dq1a1u3bmUwGHV1dYcOHdq0aRN5iswpI+5aYw8ePEBxmeQvQKNp1aSET58+IR18//59VVVVREREmygS0DRiN4QFBQX3799/9+7d0qVLY2JiUlJSjh49GhgYiM7Gx8eTUauI0NBQ8RXkk7ZSlijA39HRsbKysqKiQjABFRmATC4jI3g8XkJCAovFsre319LSEjSK2trarS8EUVtba21tDQA9e/bMyMhoxcthpIWKiorbt2+npKQ4OTnFxsbGxcUFBASwWCx0NiMjY+7cuYJRMZGRkWvWrBGTMNJWawzp4OrVqz9//lxaWurm5iZ4VgK1xh4+fIhm52vXrkUtgmnVGubowYgDScwI7969y2KxysrKOBzOmTNngoKCqPK1lcJSlqWlpWgrhSAIoUyMc+fOBQArK6uoqKhvxS2gTB9OTk6CCQcAoGPHjmRlluaPZ+vq6uzs7JBNff/+fWveCyNVsNnsU6dOFRQU1NbWnj9//vLly1S52kphrbGqqqq5c+eiP4iQIUSrlKNHj0ZuMiLv+unTpyhn4cqVK1GLYFq1V69eibxHTKO0yVyjrWHs2LEgNaUsP3/+vHbtWtIThDSE9fX1Pj4+Kioq5JxPQUGBTNH5rURWrSwEweVykaOapqYmGaiHwYgcqao1VlhYuHbtWtIqk4aQz+f7+/t36tSJ3CAXTJOL8nG3kjdv3iAv0EWLFiEz3GhaNYwEkDlD6OXlBQDOzs5UC0IQBLFs2bLVq1czGIz8/PydO3c6ODjcu3cvNjYWuWsCgLW1tZub24ABAwTzVisqKpqbm2/evDk0NLS8vLzRJ3/69Mnf39/FxaVPnz6CRlFVVXXixIk7d+4UckLj8/murq4A0KlTJ+yEiRErUlVrbOPGjcuXL2cwGJmZmfv37583b97169cTEhLMzc2RykyYMGHp0qWmpqZCeatHjBixcePGO3fu/JhRjI+PRwmkZs+eTa5L//HHH/C/adUwkkHmDGFcXBwAdOvWTQpTIXA4HAaDgfy1DAwMBGP/BXPmCqYHbDRnrhBfvnwhC0GQyiy4/8fn85cvX440sE1UbcS0aaS51lhdXR2TyUSuK926dfP39ydPVVRUkBkBBf3e6XQ6qYP5+fnN7Cg9Pd3AwGD69Olksp6GadUwEkPmDCFBEGg7TdrmPffu3UOCycvLu7u7f/36dffu3Xv27EFnk5KSxo8fj9ZPBBVSsAgA/FcEzt/fn4yeFCI/Pz8oKEjIHWnTpk0AoKKiIlanAAyGRDprjUVGRqKUfjQazcnJqaio6O+//16yZAlBECUlJSj5DhpuVlZWkik6hQJRSB38bq2b7Oxsciiwe/du+N+0ahhJIouGcOnSpQCwa9cuCfd7+vTp+fPnC7YwGAxPT8+8vDwU1wEApqamaG+gvr6+qqoK7RryeDwXFxdnZ+eGYciCOXMbKqSTkxOLxWpaIbdt24bGoeKLFMZghPDw8AAqao0FBgba2dkJtuzdu3ft2rUlJSVubm5oA2LAgAFkeHhpaSnSQU9Pz6SkpEZTrdbV1aHy4HZ2durq6o3qYNMhEGRatTaRfK5dIouGMDg4GABGjBgh4X69vb2FUu8vXLhwwoQJaKtAVVWVyWQKmTqkhMeOHYuIiFixYkXT+Tjq6+tJhUT1F0l0dXXt7e1ZLJaQF8zhw4eRBpI1SDEYCUBVrTE/P7/+/fsLtqxZs2bUqFEolaOysrKnp6fQgi3SwQULFuzbt2/u3LlxcXFNPB/poI+Pj729vVDIr6AOCu7LnDt3jkaj0Wg0Pz8/kb4rpgXIoiGsqqpSUVGh0+ktyujdeoQMIVnbCACmTp3a6LwNKSFyqBkyZEjzl03q6upiYmK8vLxsbGyERqn6+vpOTk6nT5/29PREOxwXL14UyQtiMM2ErDX28eNHSfYrZAjT09PJuCMLC4tGQ4aQDq5bt+7Tp08hISEXLlxoZl/19fUvXrw4ePDgtGnThHKEdu/e3cHB4cSJE97e3kJp1TCUIIuGkCCIKVOmAEDzv9MiARlCLpdbXV29b98+tL2nqqoquCEvyNatWx0cHMLCwtChr6/vj8V+cblcMgYf/fogOnfuTKPRpMSLHSNroFgdX19fSXaKDCGXy62pqTl06BCK4VNUVGSxWI16z507d87BweHUqVNZWVnr16/38PD44RosZMgvmfQRANCs8cCBA617LUxrkVFDeOTIEQCYO3euJDv19vaGBggl1hE3PB4vLi7uyJEjkyZNotPp+vr6kuwdgyGhpNaYn59fQx389ddfJSkDn89PTEw8efKknZ2dnJxc165dhWxwVVXVnj17kEfbmTNndu/ejUer4kZG69eg/Cn379+vr6+XWKfV1dXdu3dHAQx9+vR5+PDhwoULxZT7+FvQ6fSBAweuXr06NDS0Y8eOWVlZZJkqDEaSUFJrrK6uTlNTE0UoGRoahoaGrlmzRsj1WtzQaDQTE5Nly5YFBwd369atqKgIxXSRKCkprV+/vqCgAAASEhKWL1/u7OwsSQllEBk1hAYGBv369SsrK0O1mZqDt7f3xo0bkXfJsmXLvL29X7582fwer127dvDgwc+fP8vJyTEYjPj4+IkTJ/6I6CJCQUHB0tISAEJDQykUAyOzdOnSZfjw4bW1tREREc285fjx4wwGY//+/QBw6tSpP//8Ezm+NZO7d+9u3769sLCQRqO5u7vHxcWhtLpUQaPRbGxsAAAVASah0+kqKiro88SJE69everg4MDlcikQUWaQUUMILa8RunDhwgMHDiC/anl5eVdXV1Sf7LtkZGTY2NjMmTOnsrJSXV393bt3ZMQutbShKqmYdklLv4GzZs1CA9DPnz+np6evX78eZWv7Lnl5eQsXLpw2bVpxcbGysvKzZ898fX2lodD0d/8C5ubmK1as0NPTKysrk6BcModsGcK9e/d6eHjMnz+/srISVTkRGos1ga6ublpaWvfu3Xk8Xu/evVks1qFDh9CpT58+CZX8TklJqays5HK5vr6+gwYNun//voaGxqxZswYNGoQidqUBW1tbGo328OFDsRbdxmAEOXTokIeHh5OTU15eHtJBFMDanHt1dXUzMzN1dHQyMjI+fvx49uxZlJMMALKzs9FaIklqaipKQHjhwoUBAwYEBASoqanNmjXL0NBQTDUFfwArKyslJaVnz54J/YAcOXJETk7u0aNHd+7cYTAY/fv3F3Rzw4geivcoJQ6Px1uxYgVBEHV1dSjY7rsJIBBv375lMBg1NTVcLrekpKS2ttbV1RWdMjc3379/v+DFGhoa169fJ1OGLliwIDc3Nzk5OT09XfCynJwcsRb8i4yMXLVqFUEQ6enpK1euZDAYxcXFghegXwSc0gkjYVxdXfl8/rdqjX2LhISETZs2VVdXp6enb926FRVqQKfs7Ow2b94seLGBgcH58+fJlKEzZszIyspKTk4WCtjIy8ujtqQl2qH4559/KJQBI1szQgC4ffs2GocqKChMmjQJAG7cuPHdu/h8/uHDhzt16uTj41NWVvbXX38xGIyNGzc2cYuKikq/fv0MDAzu379/8eLFmzdvjhw5MjIyUvAaPT297t27t+6FmkJLS6umpgYAPnz4YGNjw2QyheKZ8OooRvI8fvx4/PjxKIoc7dIFBQU158ZDhw516tTJ19dXT0+vW7duDAZj1apVTVwvLy8/ePBgXV1df3//mzdvRkVFmZmZPXjwQPCabt269ezZszWv00qwDkoFVFtiSePs7ExG47m6umpoaAjmzCUrIrWIRmeE9+7dKy4u5nA4BEEUFxdramoCwO3bt1v/Ci0ChQM/f/58586dCxcuFKq1++zZMwAwMDCQsFQYWWbp0qVkjczNmzd37NgRBFJ0/tj8rNEZ4cWLF8vLy1HZsqqqKhQ7f/78+da/gghJSUkBgC5dujSdNwojVuQptcKSJi0tDcXPAUBQUND58+e5XK6CgkJSUlJSUpKfnx+NRuvfv/+4cePGjh07btw4MvPLd4mPj79+/Tp5iKIyyOnXzp07CwsLJ06c2My9fZFjZmY2fPjwx48fR0ZG9urVi2wfNmyYtrZ2RkbG+/fvpWfzEtOO+fr169ChQ1H1hrCwsEOHDtXW1ioqKqanp/v5+aE4v379+llYWFhYWIwbN04w/LxpkpOTBXWQw+EAALKyAODt7Z2ZmTlkyBAyr6+U0LdvXyMjo9TU1OfPn48ePZpqcWQU2TKEvXv37t27NwCEhYUtWLCAy+Xu2LFj06ZNr169iomJCQ8Pj4mJSUhISEhIOH78OAAYGhqam5uPGTPGysrKwMCgiSe/fv26vLycPKyrqyM/p6SknDhxQk5ODoVeSJKrV6+mpaUdPnx44MCBYWFh1dXVO3fuFLyATqdPmTIlICDg3r172BBiJICOjg6qPfvvv//Onj27trZ29erVhw4devfuXXh4eHR0dHR0dHJycnJyMjKKurq6Y8aMsbS0NDc379+/fxNPTkxMvHjxInlYWVlJfs7NzT148CAA+Pr6CpYVlBJsbW2PHDkSEhKCDSFlUD0lpYCHDx+iMPa1a9cKnRLMW93MnLnEt5dG0We0B7Bs2TKxvtQPc+XKFQCYMGEC1YJgZIjXr18j/Vq8eLGQNnG5XDJvNVkdHtGtWzd7e3sfH5/Y2Fihu761NIo+Ozo6gsQzSTUftG1pampKtSCyi8wZwqdPn6IEgy4uLk3X5hVUSCHfZSGFbMIQhoeHA4C6uvrHjx8vXLgQExMj3tdrOaWlpQoKCgoKCj9WaBuDaSlxcXFIoX777bemN8YE0+SiXXYSbW1tOzs7JpMZGxvL4/GaMISxsbF0Ol1ZWfn9+/f//POPFBbdrKmp6dixI41G+1YZUYy4kS1D+ObNG7Rvt3DhQuQyw+fzjx8//t1K2TweDymkUM5cANDS0urSpYutrS1SSHQ9MoRcLnfAgAEAcPDgwVWrVr148cLFxYVaX+1GsbCwAIBr165RLQim/UMWXZkxY0Z9fT1qZLFYFRUV372XzFutr68vqIMdO3bU1NScMGFCVFQUWTIQGUI+nz9mzBgA2LJly/bt28PDw9euXRsfHy/GN/whfv31VwDAlZioQoYMYXx8PBqHzpo1i9RAd3d3pJMtelRaWpq/v7+bmxtZw4VUSEtLSyaTaWxszGazjx49CgC9e/euqan5559/5s+f7+rq2mhtT2phMpkAgCpxYzDiIzMzE6mMlZUVOfpE5cDGjRvX9AqNEKQOCm3ed+jQwdLS0tPT09TUNCgo6J9//gEAHR2dsrKykJCQGTNmODs7/3AFCfHBYrF+4IcIIypkxRB+/PhRV1cXAKytrUkN3Lx5MwAoKiq2JqL8w4cPZ86ccXJyEjKKHTp0QDuRt27dIghi+fLlOTk569evRwXopYr4+Hi01vRjZZ4wmOaQnZ2NjJa5uTlpisjS0FeuXPnhJ2dkZJw/f37JkiVGRkaCOqimpqaqqgoA586dIwhi7dq1nz598vT0ZLPZonkl0ZGbm0uj0dTU1L67OoURBzJhCLOyslDMgKWlZXV1NWrcvXs3ACgoKIgwti83LY1dtQAAIABJREFUNzcwMNDd3d3MzIxGo/Xo0aNPnz7o1Jo1a86fP7948eKkpCRRdSdCkBV/+fIl1YJg2idfv35FbskjRoxAmc8Igjhz5gwKqz9z5oyoOsrLyxPUQV1d3Z49e6IR3pYtW/z8/JYuXfr8+XNRdSdCBg8eDABk/VGMJGn/hvDLly/GxsYAMHr0aHIfwtfXF41DL126JKZ+k5OT5eXl5eTkioqKiP+KkKHPUsiyZcsAYOfOnVQLgmmHlJSUoHSDgwcPJlXA398fFWc/efKkmPrNyclBcYRkGsX379/n5+eLqbtWghaoGrqyYyRAOzeE+fn5JiYmADBkyBDSK/LcuXNoHCrurenx48cDwNWrV8Xai0i4c+cOAAwfPpxqQTDtjbKyMlSnxdjY+MuXL6gxKChIXl4eAITcrUXOrFmzAEB8tlaEREdHI5cCqgWRRdqzISwtLR06dCgADBw4sLCwEDUGBASgiNq//vpL3AIcOHAAeaiKu6PWU1VVpaKiQqfTyZ8qDKb1cDgc5JNsZGRE5pcPDQ1FZcj27NkjbgHOnj0LAHZ2duLuqPVwuVzkzffhwweqZZE52q0hrKysRG7Tffr0ycvLQ403b95E41AmkymOTlFpi759+yJfgMTERBRf0SacUFD6Y2nLxIhpu9TW1qIvVc+ePckkt2w2GzmRrVu3Tkz9btq0ycjICA3p8vLyaDSaiopKVVWVmLoTIQ4ODgBw+PBhqgWROdqnIayqqkLLkvr6+uT2wIMHD9A4dMeOHeLresSIEQAQHByMDlFGN+ncnBcCBXvMmTOHakEw7YG6ujo7Ozvkjfz+/XvUGBMTg9JZoOpgYgKVfSeHdGhZKDQ0VHw9igqUIs7KyopqQWSOdmgIa2trUVazHj16kPX/oqKiUEHqNWvWiLV3lMyTTKi2cuVKANi+fbtYOxUJ6enpANCpUycpjHTEtC24XO68efMAQFNTMyEhATUKplUT6xqJ0JBu27Zt4ja9oqKoqEhOTk5RUZF0rMVIhvZmCLlcrr29PRqHkoEKT58+Rc5jzs7OLQra/QFiY2PRWhDqCJUZMzMzE2unogI5uD9+/JhqQTBtGD6f7+LiggZVsbGxqLH5adVaj9CQrm3VGkN5t2/evEm1ILJFuzKEqGI1AGhoaLx+/Ro1vn37Viitmljh8/mo1m5cXBxBENXV1aqqqjQaTayV6EXFhg0bAGDTpk1UC4Jpq/D5/OXLlwOAqqpqZGQkamw0rZpYERzS8Xg8bW1tAJDOEF4h9uzZAwC///471YLIFlJXkeSHIQhi+fLlly5dUldXDwsLQ3FLKSkpU6ZMKSkpmTVr1tmzZyVQgYVGo02ZMgX+KzmtrKyMCjsI1cWWTqZOnQoA9+7do1oQTFvFw8Pj5MmTKioqd+/eHTt2LABkZWVZWVl9+fLFysrqypUryFtN3AiWfUe1xqCNVIFHOnj37l2CIKiWRZag2BCLjvv37wOAmppadHQ0aiHTqk2ZMkWSiYtQddCxY8eiwxMnTgDA7NmzJSbAD1NXV9epUycAECpkj8E0h+fPn9PpdEVFRdIzJScnp2FaNQkQEREBAP3790eHbajWGJ/P19PTA4A3b95QLYsMIdWG8MSJE5s2bWro0O/p6Tl79mzStlVXV1taWmZnZx85ciQ8PBw1kmnVJk2aRKZVkwxlZWWKiopycnIoeDErKwsA1NXVa2trJSlGiyD/mLNnzwaAEydOUCsPRkrw9/ffsGFDw4D0gwcPTps2TbBkhLW1dWJi4vnz52/cuIFa8vPzG6ZVkwxCQ7ri4mJ5eXkFBYXS0lJJitEiSB10dXUFgL1791Irj0whvYawurra2dk5Nja24UBy5syZNBpt165d6BCVok5OTiYvINOqjRo1qjnlXUTOxIkTAeDy5cvoEFXWlsJCaARB8Hg8BoOxefPmrVu3EgRx7tw5AJg6dSrVcmGkggULFrx69aqsrEyo3dnZmUajCcYCKioq/vvvv+Rho2nVJInQkA6t016/fl3ykjSHbdu2bdu2DeVXu3nzJgCMHj2aaqFkCOndIywsLPz06dOnT5+WLFnS8Ky9vb23t3dqamrDU1wud/LkySkpKcOGDbt//z6KW5IwglsUIN17b3Q6fe/evXPmzOFwOAAwdepUOp0eERFRVVUlSTGePn06cuRISfaI+S4cDic1NTUjI2PhwoUNz86YMcPPz+/t27cNTxEEYWdn9+bNGxMTEzabLVRoXjIIKZ006yAA7NixY968eXV1dQBgZWWlpKT0/PnzwsJCScqQkJAwYsSIiooKSXYqJUivIdTW1tbT05sxY0ajZwcMGODo6Iii9ISQl5ffs2fPsGHDQkND1dXVxSxm4yCtCw0N5fF40MAuShscDuf+/fvI2VVbW3vw4MHV1dX79u0rKCiQjAB1dXVnzpwZOHCgZLrDNBM1NTUDA4OZM2cqKSmhb7IgBgYG7u7ubm5uDU/RaLQ9e/aYmpqy2WwtLS1Jyfs/2NjY0Gg0ckiHdPDevXt8Pp8SeZqmvr4+JCRER0cHANTU1EaPHs3j8Xbt2vXlyxfJCMDn848ePTp69Gjp/PuIHaqnpE0RERGxefPmhrWBZs6cuWvXroKCgq5du165cqXh0ihBEJRnNUM5ZZ4+fUoQRH19PQolTk1NpVaqhtTV1UVGRnI4HHt7e4IgSktLe/Xq1aNHD/T1MDQ0dHNz8/f3JxP0iIPdu3e/fft26dKl4usC82M8e/bMw8ODDIQgcXZ2Xr9+fUVFRc+ePdHyo9DSKCEFOmhmZgYAZLVRVGuMDG2UHvh8fkRERHV19axZswiC4HA4P//8M3KZQTro5OTEYrHI9CDiwNfXNzo6etOmTdK8jSo+pHdGCAATJkzYt2/fL7/80uhZTU1NLy+v9evXI0MohAQiJZoG5XlCs0B5efnJkycDQGhoKLVSNURBQaGoqOjAgQO7du2qqKiYMmXKp0+f+Hy+hYWFqqpqenq6n5/fokWLevXq1bdvX1dX14CAgMzMTNHKUF5efvny5efPnyPXX4z0MGLECC8vL7TB1pAOHTocPnx4y5YtjS7iUa6DQsuhSCWlcHWURqPV1tZ6e3tv27atpqZm+vTp79+/r6urGz9+fIcOHdLT0wMCApYuXWpoaGhoaLh48eK///47LS1NtDIUFRXduXPnyZMnUrtwJV6otsQ/ApoREgTB4/FGjhy5du1aaDAjpBxk84YOHYoO/f39AcDa2ppaqZqgYYLW+vr62NhYHx8fe3t7lJSARFdX197e3sfHJzY2VlTJei5evCiS52AkAJoRos/W1tYrVqxoOCOkHKGcMtJfa6xhglYul0vqoNBWa7du3USug0FBQW0iO7nIaduGkCCI169fq6ioSKEhrKmp6dChA41Gy8nJIQgiPz+fTqcrKytzOByqRWsEMkGrnp5eWlpawwu4XG5CQgKLxbK3t0e5skh0dHTs7OyYTGZsbGwTq2F//fUXuUhFEERKSgoO0mi7CBrCjx8/dujQgU6nS5sh5PF4aNcN5ZSR8lpjjSZoFYTH45E6KLTzqqWl1RwdPH78eFBQEHmYmZmJK10g2rwhJAhi9erVUmgICYKYNm0aAJw5cwYdDh8+HADu3LlDrVQNIRO0amlpNTMNVVpaGovFcnJy6tmzp6BCqqurW1paMpnMqKgooeTdqqqqqqqqZKh+UFDQTz/9JOpXwUgIQUNIEISnpycASJshJAgC+bseOHAAHaKyUP7+/tRK1ZBGE7Q2DamDaO+TpGPHjqQOCsUu9+jRQ1FRMTExER0+evRIXV1d9C/TBmmThjAmJob8vyQIoqysLDAwUArztZ88eRIA0AY40aAwhZTwrQStVlZWZDXjpklLS/P393dzc0MZDEg6dOhgaWnp6enJZrNra2tVVVVHjRr166+/oruwIWzTvHz58u3bt+RhdXV1YGBgM78wkkQop4x01hr7VoLWiRMnNjNHMamDJiYmgjqopqZG6mBNTU2PHj1GjRo1btw4tJSKDSFJmzSEgoSHh0+fPv3UqVNUC9IIWVlZNBqtQ4cOKGfEy5cvAaBnz55Uy/V/8Pl8Nzc3NJN78eIFakxISNDU1ASAP//8s6UPzM3NDQwMbKiQWVlZqqqq9+/f79q16+3btwlsCNsR//777+zZs729vakWpBFKS0sVFBQUFBRKSkoIaa01tmnTJgBQUVEhc26kp6cjz+0fGDcjHXR3dzczM6PRaKQOxsXF9ejR48aNG/r6+ihdFzaEJG3eEEp5KUsUG4cSv5GFKeLj46mW6/+Dyk2oqqo+efIEtXz8+BEJ2foErVlZWQEBAa6urpMnTyYIQlVV9fnz58ePH9fX16+srMSGsN0g5bXGLCwsAODatWvosF+/fgBAfuEpB5VLVFRURIm2CYLIyckxNDQEUSRozc3NvXz58vLly1FgYo8ePUJDQy9duqSlpVVYWIgNIUmbN4RSXsqSwWAAALmbgrLkMJlMaqVCbN68GWkg6cMi1gStyBByuVwzMzMGg4ENYbtBymuNMZlMAFiyZAk6lKpaY4cPHwYAOTm5q1evopb8/Hy0miKOBK3IEBIEMXHiRDc3N2wISdq8ISSku5RlZGQkABgbG6PDa9euAYCFhQV5QWVlpVjjZL/Frl27AEBBQSE4OBi1iDtBKzKEBEE8f/5cRUVl//792BC2G5DT/9mzZ6kWpBHi4+MBQFtbG7lTPnz4EAAGDBhAXlBdXU1Jpotjx44BAJ1OJwOHSkpKhg4dCgCDBg0SR4JW0hAmJSUpKysfPnwYG0JEezCE0lzKksvlomCDjx8/EgKFKdC3/OvXr/PmzTt48KCvr68kpfL19UXjUDIteEFBAcoMbmpqWlxcLI5OSUNIEMTSpUs1NTWxIWw3SHmtMbTOgXJUCRWmqKiomDNnztGjR3fu3ClJkc6fP0+n02k0GunfUFZWNmzYMADo27evmAI8SENIEISHh4empiY2hIj2YAjfvHkDALq6uqKKKhUtKDaINHUTJkwAgCtXrhAEERUV5ePjw+VyJenGdvbsWRqNRqPRTp8+jVpKS0tRMqqBAweKz/FP0BAWFxdraWlhQ9hukPJaY8gnc8eOHehQsDBFYmKip6cnn8+fNm2axOS5fv06KlBMxnVwOJxx48YBQO/evcW3wixoCDkcTq9evbAhREh1irVmMnjwYD09vby8vHfv3lEtSyM0UYli1KhRHA7Hy8uLy+VKRpiAgIDff/8dAI4dO4bKnnE4HDs7u1evXvXp0+fBgwdCwfIiZP369ahOMgB07tw5ICBg2bJlYuoLI2F69uw5YMCA8vLy6OhoqmVphCZ08Oeff+7UqROaDkpGmODgYAcHBy6Xu3fv3j/++AMA6urqfvvttydPnvTo0YPNZiNvNXGwYsUK5IYDAKqqqv7+/igtF6Y9zAgJ6S5lWVBQICcnp6SkhDbekpKSAEBTUxPtWMTHxwcGBh46dEgCkty4cQONQ0lvnYZp1SRDVFTUH3/88fnzZ4n1iBE3KAZgw4YNVAvSCFVVVaqqqmROmby8PBqNpqKigtI8JSYmhoSEkPNFscJms5WVlQFgy5YtqKWurg5l3iDTqkmGFy9ebN68GW3Z/L/2zjOuiayLwzcJSWgCioCIsqCgYgddim3VRUVF1oYKotjrqqurKzYUK/a+ClbsDUUsVBUXbAgI6CsICUjvJfTUeT/c3TGGFiDJDHKfHx+YyTBzApycW875H8QPEghJ3srSyclp9erV+fn58BA2pnj79i2fz79x40ZwcLACbAgKCmIymQAAfC+Ey+XCoXF9smqyJTc3Ny4uDn4vobmD+AEICwsDAJiZmRFtSN0sXrx4+fLlGRkZ8BBvTCESiW7fvv348WMFbKy8evUKtkddtWoVPCMuq6aAqqqCggJcNENCc6eN84MEwoqKCiaTSaPRCgoKiLalcWxtbVVUVLS0tGSumVsfz549g+NQ2AIbwzC4MQkA0NHREZfpkROhoaFUKnXo0KHwUEJzB/EDQOZeY7VxdHRUVlZu166dNBKdMuHdu3ft2rUDAMybNw/6u0gkgktZmpqatZvNyZyoqCgajda7d294KKG508b5QQIhhmFjxowBAJC/g8HVq1epVCqDwRBfoNbT03N0dDx58mR8fLzMHfL169dwHLpy5Up4RigUzp49GwCgpaUVHR0t28fVSXl5ORypwGQcCc0dxI8BHFqdPHmSaEMawd/fn06nw20CHG1t7cmTJx89ejQ6OlogEMj2iXFxcbB3xPTp0+HNRSLRihUrwPeyanJFYqQiobnTxvlxAiEsTXV2dibakIZ48OABdL99+/bhmrmGhobiDimumdtyIagPHz7ADkqurq4wxIpEoqVLl4LvZdUUwK+//goAuHHjBjwU19xp7ZSXl+Nz+hZKgbRqLl++DAAYP3480YY0RGhoKFwdWbt2LS7RaWxsLO6D4jK5LR+rffnypVOnTgCAyZMn8/l8eBJKbYjLqimA6dOnAwBOnToFDyU0d1o1FRUV+BSiGT744wRCFosFAOjQoQP+r0Y28F267du3S7wkJ4f8+PEjzAKdNm0a/muBiWqqqqphYWEtf1PSc+TIEQDA7Nmz4aGbmxsQ09whLYWFhdeuXRNfu2az2eItRPbs2bNnz56ZM2fyeLzZs2e7u7sfOnSICEuJh+S9xjCx1ZHff/9d4qX6ZHJVVVWHDh26cePGkJCQZsgtpaWlwZHumDFjcBd2d3cHANDpdAX3orl06RIAYMKECfBw//79cKlWkTY0g7KysmvXrol/AGZkZNy7dw8/PHLkyJ49eyZPniwUCl1dXT08PNzd3Zv0CBkEwvLycvxDlsPhtPyGzcbU1BQAEBERQaAN9REREaGmpgYAWLNmTcNXijukuGauuENK0zwzOTkZlis4ODjgM8vasmoKIzExEY5U4NKQhOYOaXn37h0AQHyt7OLFi7XNnjVrVnJy8qZNmzAMI/myhFwhba8xrK7VkfrIycmpU7daSUkJCgT6+/uXlpY2+sSMjAw4tBVXDcVl1WAxsSLJy8ujUql4uqyE5g5pSUpKAgCI53/4+fnp6OhIXLZw4cL09HQoUz5lypQmPaKlgfDMmTMHDhxwcHDg8/lr1qzZt2/f0qVLW3jPZrNmzRrQrJ4J8ubNmzdwn3zBggVNyotptkPiqqG2trb4MHbXrl1wHAr7PygeExMT8F/XOgnNHdIiTSCMiopyd3cXiUQ7d+7csWNHW05A2LFjBwBg+fLlRBsiSXx8PPx/w3fppCQvL8/f33/jxo2DBg2iUqkSPrh69eo7d+7UKcaUl5dnZmYGALC2tsZVQ6GsGoVCISplGorX4ALf4po7pEWaQJiQkACXlw4ePLht2zZra+smPUIGM8LExMT58+cXFxevWLECw7D169fn5OS0/LbNIDg4GAAwYMCA2i/FxcXBPfzc3FxPT095aGnWR2xsLNwnnzt3bktGXvn5+XU6JI1Gk3DIzMxMWKExZMgQ/J3ismr4Lp3iWb16NRCroJLQ3CEnjQZCf39/Nze3oqIioVBYWFj46dOn1atXE2EpKWig11hiYiJcNC4qKjp8+DBeTaQAkpKS4C7db7/91pJ997KyspCQkI0bNw4dOpROp4v7YO/evZcsWYI3ZSwpKTE3N4cfR7hqKC6rdubMGdm8saYDWyjDz2qsluYOOWk0ED579mzNmjUFBQUCgaCwsJDNZi9ZsqRJj5BBILx///7SpUszMjI2bNiAYZiHhwdR+dNcLhdOvNLS0iReCgwMnDt3LoZhERERy5cvV1ioTkxM1NPTAwBMnTpVhpuXRUVFDx8+XLdu3eDBg2k0mrhDmpubwydaWlriK9UXL16Esmre3t6ysqEZBAYGAgDMzc3h4ZUrVwAA48aNI9CkRoGB0M3NbdN/TJw4UTwQHj9+3NPT09PTs6ys7MyZM97e3jJv3FGb8vLyt2/fws/0goICOcnDNoMGeo2FhIS4uLhgGPbu3buNGzd+/vxZMSalpaXBHu5jxoyR4Z+Gw+E8efJk48aNNjY24kGRSqX269cPdhPs06cP/vHt6+sLE+UOHDggKxuaQWRkJADA0NAQHj569Ah+VhBoUqPAQPjHH3/gPjh9+nTxQOjl5QV9MC8v7/Lly6dOnWrqVKelgTAoKKisrGzevHlZWVmOjo4FBQUwa6CFt202v/32GwCgzo97fIzg7u6umEAow95+DVBRURESErJ9+3ZbW1uYjGNsbNy5c2dcNRQWbAAADh8+LCcbpKSmpkZdXZ1CoWRmZmK1NHfICQyE27dv3/EfkydPJnZrk8/n79ix49atW0uXLo2NjZ03b97Vq1cJtEeCBnqN4fsmR48eVUwglGFvvwaorKwMDw/39PS0tbWFKalGRkba2tq4amhAQAD0TcLVr4RCIRwof/r0CauluUNOYCD866+/cB90cnKqvUfYEloaCBMTE/fv3w8zgJOTkw8fPqywgV6deHt7wwWQ2i8pOBDKtbdffVRWVsIlUGNjY/wk7MdG7DgUx8HBAQCA631bW1sDAIjas5QGKZNlFEx5efn58+fhqvLLly9JFQhr9xrDUXAglGtvv/qorq6+evUqAEBbWxv/tzlz5gyVSoW5VITj6uoKANi/fz88HD9+PAAA9qwnJ1Imy7SElopu9+zZ86+//oIdFQxpYJGhnl5YQOn1C9UxkUAkauHNm8HEiRMpFEpoaGhNTY34+Zs3byYlJd25c+fVq1cxMTFeXl5yNSMvL2/MmDFfv361sbHx8/ODg0T5ERAQsHnzZicnJz6fv3Llyo4dO6ampsL/HgDAxo0b3759u2HDBrnaICUS8scShwgpodFoBgYGmZmZEuf5mWllD24Vnz1aeu18deQrTKggMXdxxo4dy2AwXr16VVxcLH7+wYMHSUlJV65ciY6OfvXq1aVLl3g8nvzMKC0ttbOz+/z5c//+/Z8+fQo3TeRHRETEX3/9NWvWrJKSEhcXFxMTk6KiIrhjCgBYtmzZ69ev9+7dK1cbpAQKKyIf/A6ZhFNuSnLGnN9ShpqlDDVjWRixLU1Thvf9avtzxQtFqGhKMHDgQABAUFAQfiYyMvL27dvNngWWl5d7eXnBHe9Hjx7t2bOnYSE3BfT2kwAmQ586dQqWBkLVGMUIeTeVjIwMcU2Z6OhoAECXLl3I2UILI+WMsKam5uzZsxEREXPmzCkvLz937tyuXbuKvyRkLZ7FtunFHtqbZWHEtjRJGdEvddTA8qcENKwW7zUGiY+Pv3btGlwSbwbV1dUXL16EK42hoaG7d+9u+Fbivf0Usw8CffDGjRuwQ/iqVasAANu2bVPAo5sKh8OBwjpQUyYlJQUAoKmpSeCWVsO0ghkhAKA6+m3W3N/YsR8qKitF1dUAAEzAF1RWJOdkZ21eU3zueMsf0SRqD3DOnz8/c+ZMfX397t27z50719vbOzU1VfobxsTEfPr0CSYEvn79etGiRbAar044HI6dnd3//ve/fv36hYaGwroleaOqqsrlct+8eQOdn8xDvC5duvTr16+ioiI8PBwAYG5u3rlz58zMTFjSRCqePXt27949bW3tGTNmiKfpGhsbw98wUTCZzF9++aWwsPD06dNcLldTU7ObmkrK3Ck1se9Ty8pTSjmpXH5KZTW7qJhTXFywZ0vhQQ8FW1j7P/D69esuLi5dunTp3LnzjBkzvL29//e//0l/w9jY2OTk5PT09IyMjAcPHqxYsaIBH6yqqnJwcHj//n337t1fvHgB80XljaqqqlAoDAoKgm0F4W8AdnoiGxoaGkOHDhUIBDDN3tjY2MzMjMPhvH79mmjTJImMjLx06VK7du1mzJgBN1khBgYGkydPluWTWhhI+Xk5KcP6siyM9Oi0k8Y6LAsj+PVhgCEAILC3QcqQ3gqeF7569QoA0K1bN/zMuXPnxo4dCxUlcIyNjV1dXS9evChNjuuhQ4e+fPni4+MDJcFgAmptKisrhw8fDgAwNTVVZI+h/Pz8NWvW4GPkoqIiJSUlBoOhsH2RJgE/wtauXQsPFyxYAADYt28fsVZJAAUQaDTamzdviLalEYScktRRA1gWxiwLIyUKRVuJps9Qgl/HjXVYFkYpQ83KHt5WpEmw15iOjg5eL3T9+vWJEydqaGiI+2CXLl1cXFy8vb2laT904cKFN2/ePHr06M6dO1j9PsjlcuGmV5cuXVJSUmT4phqmtLR03bp1+IeJRF4Y2Thw4AAAwNXVFR5CtSmY9k8eYmNj27dvT6FQAgMD5f2sls4Ii08dwHjcBi4Q1VQV7N0CRMIWPkh6rKystLW1U1JS8E2yRYsWBQUFlZaWRkVFHTt2zNHRsUOHDqmpqT4+PgsWLDAxMdHX158xY8bx48ehAnV9d9bU1CwpKQEA1HlNdXW1vb19eHi4oaFhSEgI3oRWATx69IhKpZ48eTImJgYA0KFDBysrKx6PFxoaqjAbpEdisCzeJZUkxMbGTpo0qbKy0sXFBUqlkJni86dFNTUA/Ps/ebqbTnjfLvBrYns1AICourrw8J6G/VS2mJmZde/evaCgAN8kc3Z2fvz4cXFx8adPn7y8vBwdHTt27JiZmXnt2rUlS5aYmZnp6elNmjRp//790dHRovrTCzQ1NUtLS0E9Psjn86dPnx4QEKCrqxsSEiIhWChXAgMDRSLRuXPnYGtiJpMJdRVgyRDZwLcJ4a+ahGtISUlJdnZ2JSUlv/32G5Qpli8tiaIiPp9t0wtOAeubEbIsjFKG962OUZy4M/aftqyDg8O7d+/qrN4TCoW4Q+ro6Ij/QnR1dSU6s7x//37ixInr16/ncDizZ89ev369j4+PxA0V3NuvUfbs2QMAWLRoEdGG1AGuKZOUlIRhWFlZGYPBwBtTEA4uQSIu0EpmUn/pj/udEoVyq0cn/BD/ShneryIsRJFWLVmyBAAwevToiIgILpdb5zVoSlUmAAAgAElEQVTS685//Phx+vTpK1asKCkpmTt3rpubW+0eFwru7dcoJO81BqtK3r17h2EYj8eDjSlSU1OJtgvD6hFolSstCoRc1peU4X3xQLhIT/Pvbrrw67CRDh4I2VY9Sq4pTk+Iw+H07NkT1soAANTU1BrVzMUdEhbe1umQuDPX3lJWcG8/aYiNjQUA6OvrkzMJxcnJCQBw7NgxeDh69GgAwM2bN4m1CsOwpKSk2gKtZIafl4MPRmEgXNZJ0/Onjp4/dTwqNjBlD+5WdEpx9TNVVVUWFhYdO3aEfiSNTG6TdOdr/2kU3NtPGkjeaww2gcIbAMDJw+nTpwk1CsPqEWiVNy0KhNVx0Sn/jUb16LQ+qoxRmirwa7iGCh4IWRZGxV7HZGVxw1RUVAwZMgQA0LVr19mzZ/fs2VPcqVRUVEaNGrV9+/bnz5/X55Bfvnw5d+5cnaPU8ePH79u3T0IaSvG9/aQE2o83pCYVsNBq7Nix8PDw4cMAACg7gmEYl8uFOhFfv35VpFW4BIm4QCvJ4aWy2P8NRmEgtG6nPL692vj2ao7a6uKTwvw9WxRjUk1Nzbhx4wAAenp6Li4uEtrxTCZz2LBhW7ZsCQoKqk9IISUl5fLly/Pnz4dKgThqamq2trY7d+6U2HhTfG8/KSFzrzG4GTF48GB4ePHiRQDAxIkT8QsOHDjg6emZmJioSKvqFGhVAC0KhLzMdPZQs0aXRtlDe3PuXZeVxQ3A5XLt7OxgFMTn+Lm5uY1q5tbXmrLORhC4WgT2fW8/uMhAHuDa1K5du4g2pA6Kioqgpgz8X09ISIArWrBKISUlRfG1H7gEibhAK/kRckrZVqaNLo2yLU1LLipisC8QCBwdHeHqCF4yz+Fw6pPoxH0QF+SUoE7deaiKggN7+zGZTAVkVTQJaBg5e41BTRkKhQI/0CQaU5SUlCi+9qNOgVbF0LKsUZEodbR5Y4Hwp5ShZryvct824/F49vb2cBxa3yhGes3c2mRlZd24cQMXjIYQ1dtPGvz8/AAANjY2RBtSNzY2NgAAPz8/eAgbU8AUzYSEhJUrV65cuVJhWX+4BImlpSU5U20bIM1hRON7hMP6VMfJfblCKBQ6OzvD1ZH6liLKy8txRUAGg4H7IJVKxX2wvjrdvLy8u3fvrl+/XnzBH6pIK763nzSQvNcYTGu4ePEiPBw8eDAA4MmTJxiGZWVlLV68eNWqVRJjDvnB4XCgAT179lS83ltLyyeKL5yGRfT1BsLBxhlzJ8vE1gYQCAQzZ86Eswop/3LiDileoQIA6Nat25IlS3x8fDIyMhq4w+bNmwFBvf2koaKiQllZmUql5uXlEW1LHcCeULjuHSxA3rp1K37B169fa3cwlgclJSUWFhYAgP79+yt4HCoTOL43Uob1bigQDjJKmzxa3maIRCK4CKGhoREZKVVyHCwnFZforO2DtQX0xSGwt580CAQC2HmGnL3GTp8+DQCYPn06PIRDipUrV+IXlJWVrVq1SgGWVFZWjhgxAgBgYmIivuSmMFoaCEVcbtrkUWxLEws15gUTPdz34gYYdmEoPetjwLYyKXsirbbFzZs34UZRcnLy6tWr161bJ81mqUgkWrhwIQBAU1MzKiqqGe+itmauuEPOmTPHy8tLYr9q9+7dgNDeftIwduxYAACphChxYKVH586d4eg+ICAAAGBhYQFfWr9+/YIFCxSQ9cDhcGCBhMIkSGSPUJDpMolt1YNlYdRNmf6gl37tGWHpjYtS3uzhw4fTpk3DMCwzM3P16tVr1qyRcnDw559/wtWRly9fNuNN8Hi8qKgoT09Pe3t7iXJD3Acl8rFPnz4Ne6oQ1dtPGuAA/cSJE0QbUgdpaWlw4AIzAaGIEmxMkZSU9McffyxduhTqSMuVOre0FIwMJNYEhflpU0an/LdZ+O3L0jRlSJ8Um17sIWaVr6XyjaKiIijLu2bNmry8vEePHjUqBSsSiZYtWwY30sPDw1v+dqqrq8PCwjw8PH799VdVVVVxhzQxMVmwYIGPjw/sPkpsbz9pOHbsGADAycmJaEPqQCQSHT58ODY2Fh7iBcgvXrzAMIzP5yuga3ZVVRXUATE0NGx45kFyhJzSDKcJ7Fo+yP7ZJMWmV8rwvmwr0/IQqdYtiouLoQ/u2LEjMTExJiZGmv1aKJIgq9URLpcbERGxZ8+ecePGSWiEGhkZzZ0798KFCwcPHiS8t580kLzX2IkTJ96/fw8Ho3hjCj8/P5FIJBQKm9TBuHlIs6WlAGSjNSricosvnE4dbZ4y1CxleF/28L5sm1557uv4+bk1n+JSRg5kW5lWPAuQ5lbQCf39/RctWrRy5cpGeyZALWkVFRX4ASpb+Hw+XoMvLpampaVFoVDwtXXSwmazAQDt27cnfz1cQUGBjo4OzHStXcopD7hcLqwjVrAEiZwQCfgl1y+mjv2ZPcQsZXi/lOF92dY9cjau5GdlcFlfUsdasgZ3L/O/K82toA+GhYW5urquX7/ezc2t4et37twJV0f8/f1l8E6+RyAQ4CW/sL4T90EAwNGjR2X+RNlSUFBApVJJ3msMwuFwDAwMYOJ07VJOeSBe+qmwzcg6kU0gxOFlptd8jueyvogE3z55axI+po42Zw3uXvbIt9E74I1ahEKhj49PwwPMLVu2wHHo48ePW2h5owgEgvfv3x8+fHjYsGFUKtXExETeT5QJPXr0AADIZK4sP/BsMT09PQlFng4dOjg4OBw+fPj9+/cyHJ/iyY26urrENg6TOfzc7JrP8dzkRJHY5xcvlfXVzpo1yLj0ZuPddsR98NGjR9evN5TyDdt+0Wg0BZSBCoXC2NjY48eP29ra0mi0Tp06yfuJMsHKygoAII9RggzB5SG1tbVhV2EcDQ2NiRMnenp6vn79WoZBseVbWjJExoGwPrgpyV/HWbEGd+c8aGhP+/r1605OTidOnIiPj9+0adOBAwcaqAc/cuQI9ECoPagwWkUrS5y1a9cCAEjSCK1Oagu04rXVsJ4BB5aRSdRWN4+5c+dCnyeDBIli4GdlpDn8whpkXHL1XAOX3b9/38nJydPTk8Viubm57dmzp4HPvgsXLsBdOrzBpGLAxYm+fPmiyOc2DzhjxocXJKS2QKt45Zi4D4prI7Sw1vaPP/6AU08yyPkqKBBiGMbLTE+bNJw1yLj0+oWW3+3kyZMAACqV2vBwVU6Qv5UlTkhICACgf//+RBtSN1VVVf82szQ0rLN8vs5STiCdWAmGYZWVlXfu3BHP2UtISHj58uWLFy8MDAxIIkGiMPi52WmTR7IsjApP1NE+vqlcuXIF7tIRIkcC6zTIvzSKYVhUVBQgca8xHo8H22Xr6urWqX5eZykn+K8Ue+PGjf7+/qWlpfXdn8/n37lzR1xyi81mh4aGRkZGGhgYyGNLqxkoLhBiGMbPyUr77ReWhVGx9/GW3Ofy5cvQA728vGRlW5M4deoUAMDR0ZGQpzcJLpcL0w0UrNIiDU0VaJXGISW0EdLT0wEAffv2xac1e/fuhYo2DYTPHxhBUUH6TDuWhVHh8Ra1+7h//76SkhIAwNNTBjG1GVy/fh0AYGtrS8jTm4RIJOrcuTMAID4+nmhbJBEIBFDvUEtL68OHD41en5eX16g+iUQT1vLycgDATz/9hJcAnD179ueff8bI5IMKDYQYhgkK89Mdx7ZkTHr37l0ajQYAOHTokGxtkx7YylJDQ6NVyFFOmTIFAHD27FmiDfmOFgq05ufn1+mQEmIlMBD27t0bz7rCA2GbRcgpzZw7mWVhlL93K9asXKSgoCBYeuvh4SFz86SE5L3GJCBnrzGRSLR48WL4UdaM1RFxfRJxbQRxfZKCggIYCHv37o1v0OCBkDwoOhBiGCYoKsxwmsCyMCrwdMeauFbg5+cHFWH27t0rJ/OkBArikVBQpjbnzp0DADg4OBBtyDdkK9BaXFzs7+//559/WlpawmkKRF1dHQ5Znjx5oq6uDgskUCDEMExYXpY5fxrLwihv6x/ieW3S8OzZM1hr+8cff8jJPCkZOnQoAOD+/fvEmiEN9+7dAwAMHz6caEO+IS7Q2rzST3HKysqePn3q5uY2ZMgQcdEuJSWl7OxsAIC/v7+Kigoc8qJA+C/CMk7mvKksC6P83ZukH5MGBwdDD1S8CF5tyNnKsk6ys7MpFIqamhpJhKTlKtBaUVGBCwZNmDABzghLS0tnz54NhwIoEEKEVZVZy11YFka5m1aJpK6uef36NWxwLa4/QhRk7jUmAYfDIVWvMew/HVQGgxEQIFVhm/SI65PY2NjAGWFqaurvv/8+YsQIkUiEAuE3hFWVWUudWRZGuZvXSDMmhR3DAQCrV69WgHmN8vz5cwBAnz59iDZEKmBxAkkkieUh0Irv1WdmZkLVYAgeCHNycjQ1NR8+fIgCIY6Iy81Zt5hlYZS9ZoGI23gW7ocPH2A1raurqwLkDhqF5L3GJCBPrzFMPgKtuA9mZ2eLr1fjgbCkpERPT+/SpUsoEH6HqLoqe+VcloVRztpFonpad0Levn0LMz7mz59Pkn96srWybJitW7eSZAwhc4FWHo/n7u4+Y8YMDMP27t176NChqVOn4tv1eCDEMOz48eMmJiY7duxAgRBHxOPl/rWcZWGUtcRJKDaAqM3Hjx9J2K+YzL3GJDh06BAAYM6cOUQbInuBVpFItHv3bnt7ewzDzpw54+HhMWPGDFyoGQ+EGIZduXKlU6dOBw8eJFsg/JZloHgoyiqdjp5XG21X+TI0d/0yjFtT52Xx8fETJkwoLy93cXE5f/68eLoggdDpdFtbWwDA06dPibalcaZMmbJs2TKYNUMge/bs2bt3L51Ov3v3LixBaTl0Ot3DwwPOVEQikZaWloGBAUynkmDlypXq6uqw7xoCQqHT9fadbDdpenXUm5zV80SVFXVexmKxxo4dW1RU5ODgcPPmTfGNWGKBjQ9bhQ86ODgsWrQIVn0QyIULF9atWwdT7qEOasuhUChbtmwxMDAAAAiFwnbt2hkaGtb5T+Li4tKzZ09Y/EYqiAyEAAAKnd5p/6l29lOrXr3IWTVPVFUpcUFSUtLYsWOLi4unTJly6dIl8fxAwoECXa3CCZ89ezZlypSRI0fCw8TERG9vb/h9SEgInDbJm5MnT27dupVGo/n4+MC6JZlDo9HodHpxcTGHw6nzVW9v78zMTHk8uhVDpeluP6A5y7Xmw/vsZc5CTonE6xkZGWPGjMnJybG1tb19+7Z4KgThwPIb2GOW5AQGBk6aNAmqSwMAUlNTT5w4Ab8PCwuDWV3y5sqVK7BDyN9//w1VXWQOHIwWFBSUlpbWfpVCofz99985OTnyeHSLIHpKimEYhgkFeR5/sSyMMudPF5Z/W1xmsViw/mbs2LEtVBKRBxKtLMkMnU7X0NDA+3rfunXL1NQUwzA+nw/XnPX19R0dHb28vD59+iSPxeeLFy9CCRJvb2+Z3zw9Pd3JySk/P9/FxSUjI8PDwwNvU15YWOjo6Cj+B9q7dy85mxUTjEhUeGQ3y8Iow2mCoPhbu4nMzEzYJp6c/YpJ3mtMHF1dXSaTietKBwUFaWtrw+/hp1ynTp0cHR2PHTsWFRUlDx/09fWFs7RGBZybQWZmppOTU25u7oIFC1JSUo4fP4635amurnZ0dMzPz8cvPnHiRKMCtgqGHIEQwzCRqODQTpaFUYazvaCkGMOwjIwMY2NjAMDo0aPJU3cpgXgrSzJDp9Otra1nzpwJD/FAmJuba29vDzc7cTp37uzk5HTmzJlmlPfVydWrV+FU/vDhwzK5oQQ3b9708vJ6+fJldnb2iRMnGsgJ+vjx49GjR0myzUxCSi6dYVkYpU39lZ+Xg4n1KzY3N5dQKiAPZO41Jo6urq61tTW+P40HwrKyssmTJ3fs2FHcB/X09BwdHU+ePBkfHy+TvKSAgABY+rl79+6W36029+/f9/LyCg4OLiwsPHXqFN5wuzaJiYlHjhwh28SGNIEQwzAMK/r7MMvCKH2abc7nT7169QIAWFtbk3AcilO7lSU5odPpgYGBmpqaUJ0cD4Q4bDbby8trzpw5Enq7Ojo6LWwE8eDBAzgOJUM1McytaHnl4g8MjIVf7YcVJXyC/Yr79etHnqT/2pC515g4urq6fn5++vr6MHFUfEYIwX0Q/pfiiDeC4DaYVFgfISEhsPBs7dq1snkzLaB///4AgJCQEKIN+Q5yBUIMw4q9jn2x6BY6qJeqcrtBgwY1IGFHBsRbWZIZOp0OG8t17969qqqqdiAUB9e8hg1ZWuKQuASJYnrNNwrcIEFLow1TevMSa5Bx5GCTDioavXr1IvmqY2vpNaarq/vs2bNLly516tSppKSkdiAUB/dBuCqG0wzdebz08/fff5fdu2k+bm5uJAnJ4pAuEApFoul/3u0y9y5zxs33iVlEm9MIeCtLWa0iygkYCPl8fr9+/dzd3WEgFAgE48ePd3d3Dw0NrW+bMykp6fz583Pnzq0dFO3s7Pbu3RsREVGfzhxe+rlmzRp5vrkm8PDhQ7jMQLQhpEYkwv7ecMR52k6lmbcD3xLWK1V6WkWvMRgIRSLRsGHDVq1ahQfCqVOnbt68OTAwsL6lr9TUVB8fnwULFpiYmIj7oKqq6q+//urh4REWFlZfUMRLP+fNm0eG0k8Mw8LDwwEAPXr0INqQ7yBXIBSJsKUXPgLnJ703vKQ6P+20IvRjBtmFBF1dXeW0/yxDYCDEMCw8PFxVVfXAgQOmpqbR0dG4UzWgmYuTnZ1dp+Z1nV1t37x5A9NwFixYQJ49uVaUW0EgG28mAucnpuvC6HOetl8c/I5F0t1BHNjQh8y9xrD/AiGGYfHx8crKykePHtXW1k5NTcX9SFyis7616Pp05+vs5xcfHw9LP6dPn66AXvNSgrfQEu8JQzjkCoTrrycA5yfDPF5X1Aiuv8pSmvNUd1lIbBqHaLsa4vbt2wCAkSNHEm1IQ+CBEMOwefPmdezY0dTUtKyszM/Pb+3atYMGDRIvvKPRaIMHD163bt3Dhw+LiorqvGF2dvbNmzdXrFhRn/z/77//DkgjQSIOrDy7cuUK0YaQlO33koDzk4GbwosreI9i8pRdAzQXBb36UvfYiCSQvNcYBA+EGIatW7euY8eO2traVVVVjx8/3rBhg7W1tXjhHZVKHThw4OrVq319fcXzLcXJy8u7d+/e6tWrhwwZUmec27ZtGwBgypQpZFs0hl3pjx9vUQ8i2UKiQLjpViJwfmLt/qqs+t8/m19ULtM1QGtR8Jtk8o5JS0tL6XS6kpISadPqsO8DYWFhoba2tsQeYXl5OS7RCXf1cLp167ZkyRIfH5/09PTady4rK2vfvv2SJUvwMxcuXLC1tRUKhRcuXCDPOBQHdlSfNWsW0YaQkaMBKcD5Sb+N/xSW/7sHHBCbrzIvUG1+YMjHAmJtawAy9xrDEQ+EZWVlBgYGEnuE4hKdMLdF3AfnzJnj5eVV3xvU0dER7wrn6+s7aNAgDMMuXbrUvPwauXL16lUAAKkEnsgSCHc9SAbOTwa4hRdXfLfh9ORDvsq8QPUFgc8+kTdv7ZdffgEA3Llzh2hD6mXz5s1ZWd82XB8/fnzw4MH6LhZ3SBUVlTodEheWg6XrKioquHDo33//bWVlJc930yJaS26F4jn/Ip0y+4npurCcku82nP5JLNJYGMR0DfCLyiXKtkYhZ68xcXbt2sVisfDD58+fN9DEis/nR0VFeXp62tvba2pqivugeMkvfr2SkpKysjJetHDz5s0GsuEIp6CggEajMZlM8lQEkCIQHg9MBc5Pev4Zlltax5Zv2OeidguCVOcFBsWTdEy6f/9+uB1NtCH/IhKJgoKCZCKnW11d/fLly507d9ra2sLMF5zu3bvPnz8fikRs2LChd+/ecOxJ8kCIYVjPnj0BAHjRPQLDMJ9/Mqmznxquev61oI6a3ffsUu0lIYy5AXff5SjeNmkgYa+xsLCwu3fvtvw+PB7v1atX+/btGz9+PJz44hgaGs6ZM6egoEBJSWnDhg2GhoYwtJA8EGIYZm1tDQDAi+4Jh3jFsksvM/+4+tlET/XFVms9TWbtC34x6/D0r5/pSpRJh6IeROUq3sJGgTpPT58+FYlERNsCAABFRUUUCiU4OLjlt1JWVh4xYsS2bdtCQkI4HE5UVNSxY8ccHR07dOjAZrOfPHkCp4zLli2jUChQVpj8tCJtPMVw/33uwnPx+lrMF1utfuqoUvuCwd00QzdbaqoqzTr54fI/ZNSomzhxIoVCefbsWU1N3ZLFCqa8vJzH48nEB+l0+pAhQ9zc3J4+fVpSUvLp0ycvLy9HR0dtbe309HRfX18NDQ0AgJOTU5cuXTw8PFr+RAVAOm08YuPw1fBM6uynXVc9S81vRDsmOrW049IQmsvTK+GZirGtSRgZGQEAIiMjiTbkX7hc7qpVq+R3f4FAEB0d/ejRI7g0mpaW9uLFC1VVVTabTf4ZYWhoKACgX79+RBtCCgJi85muATrLQj5nNrJOlZBVbrDyGWX2k5NBZNyKI1WvMcjSpUvld3OhUBgXF+fr64thmJKS0ocPH2JiYphMZlxcHPlnhDBlvUuXLiRJKSdyRugXlTffO15HgxGyycpIp45xqDgWRpr/bLPW02TO94q/+JJ0Y1I4ySDRAEfO0Gg0CwsLe3t7/MzIkSOnTJmybt06Aq2SkhEjRmhqan78+DEtLY1oWwgm9FPhlKPRKnRa0EZLMwP1hi/u1Vk9YruNsY7q6iv/O/I0teGLFU9bm+hTqdT+/ftPnToVP2Nubr548WLYd57kmJubGxgYZGZmfvz4kWhbACCw+0Twx8JZpz5oqdKfbbbqqa/W+A8AYGag/nyLVef2zEXn4o8HfpWzgU1jwoQJSkpKeXl5RBsCAAAikWjLli2vX79++fKlIp976NChsLAw2LKYzNDp9F9//RUAEBAQQLQtRPI6qWTK0WiGEjV4k6W5kYY0P2Kko/Jiq5WJntqf1xPcbn2Rt4VNYuLEibBClGhD/mXHjh3v3r0LDAxU5EN3797NZrN9fX0V+dBmQKFQYCMOskweCJmHPvtUqDIvUHNRUFRKkxXUvhZUmax9AZyf7LxPonrMiooKiTr0+vRWfjDwpVF4eOLECQqFQvKlUQzDzp8/DwCYNGkS0YYQRkwqR2tRsOq8wJcJdZeKNkBuaU2/jf8A5ycbb5JId4bL5UroJLQRH8T+WxqF31+9epVCoZB8aRTDMBithw0bRrQhGEbI0uhbVunkI9FUCni8/udBxpqN/8D3/NRRJdzdpm+Xdu73kkgyJq2pqVFXV//69St+JiYmhsFgcLlc4owihhUrVgwcOJBoKxpnwoQJFAolNDS0urqaaFsI4GNG+Zh9kVU84b0/LEb06tDUH9fTZL7cZm3ZXWv/I/aKS59EGCYPI5uKlpaW+Dpbamoqg8HIzSVjep1ccXFxGTVqFNFWNM7YsWOZTOabN2+KioqItkXhS6Nx6WUTDrwXiLDHGwYP69m+eTfppMV8tsVqgKHG/kfsVT7/I4cbtlHU1dXZbDZsqAYAoNFoz58/J//KjL6+vrm5OSwOIdoWRZOcWznWM7Ksmn9vjcX4ATrNu0l7NXqwm+WQHu3PhKYvu0CWWNhmSUpKgg2zIA8ePHj27BmB9kiDurr6sGHDhEKhTHJrW4hCA+GXnMpxnu8rucK7q81Hmmm35Fa6GoywbVbWJlqngtOWk2ZM2gahUqndunUTV4fS0tIyMDAg0CQpaWu5FZD0ouox+yILynhXlg+YZKHbkltpqioFu1na9u147kXG7NNxfCHyQcIwNjZmMBj4oYaGRteuXQm0R0rI44OKC4SsvKpRu98WlvOurRgw0bxFHgjRUqWHbLIa3Ufb61n6nL/jBET7YVJSUux/JCUlEWsMolFgJdPjx4+JNkRxZBbXjNz1Lr2o+uyCvrNsOrf8hmpM2uP1g38bpHfrTfa0Y9E1fILraNlsNu6DCQkJxBqDaBS8AlsoFBJrCQVTyFwqo6hmxK436YU1V1cMcB4iAw/EqeIJpxyJDv5YOMNa/9qKgXQaBQDAzqsKSygqLOepMGg/d9O0MtGiiim1y5yamhoVFZVu3brhCoE1NTUpKSk1NTUSup1tgcjIyLdv386aNUtXVwbDHfkhEon09fXz8/MTEhJgF+gfm/wy3sjdbxOzK07P67vc1rDxH5AankDkfDrWNzLXboDO/T8sVBg0AEB6UfWLz0W5pVwmnWphpDm0R3saVY4+CABQVVXV1dXF9Y/4fH5ycnJOTk6nTp3k+lwSEhcXFxYWNm3aNIk+2yTE1NSUxWK9fv3axsaGQDOUGr+kKSTnVoYlFBdXfBeB8jjcMfvepRVWn13QV7ZREACgyqA9Wj941snYO29zKmoE+516rbj4v+ivHCoFVHJFTDqVTqNoqCidXdDXXhbT0Aa4d+8eLOkFAMTExAwaNEiujyMn+fn5Xl5eHh4e69evv3LlCtHmNASVSh03btzVq1eDg4N/pED4taD6xeei/LLvIlBpFd9uf2RCVoXnrJ6yjYIAAIYS9fYq84XeH33CM+32vz+7sO/aqwn/JBbRqJQqrohBpzJoFCademJun1k2+rJ9tAQXLlyAhTEAgNTU1G7dusn1ceSkoqLiwIEDhw4dWrdu3c2bN4k2pxEmTJhw4sSJ4ODgHyQQfsmpnO8VF5dWRqFQK7kCZTqVrkTVUFE66GS22y/5S07lodlmS0bL2AMhDCXq7dXmc/6Ovf0259n/igTl+UKldkCJCQCo4QlrAFZemDXzePXmKT23/GbS6N0QLSE2Nnb06NFdunShUolX72uUDRs2rFq16ocZsqQVVi/0/vg6uZhGpVbWCBl0KlOJwqRTDzr1Oh2S9uFrmcc0042Tusvj0TQq5dLS/u1UaKeC0/ptDKdQMEFZAQAiAEANADUAAKbGonPxkdG+A40AABDMSURBVCmlR2abycMABE5CQoKNjY2+vr66ujqPxxPfOyQhK1asmDVrlqWlJbFmyObTKvxL8eCtEe/ikquqqiq5AgBADV9UXsXPSv+60DumpFKwd2bPPycYy+RZdUKnUY64mKkwqVy+SBi4HmS9+/aaoAY8XFBVkLbvYcqjmHz52YAAAKipqVVUVAAASCK72gDh4eHr16//+eef8Zjt4+OzfPlyYq1qNjFfOf03hb9MKKouzqkozMQqcrglWWUFmQXFnIXeH9+ncP6y7+Y+1VR+BlAoYNf0nioMqlCECYQAPF4CAtaAoHX/fmW+reQKvZ+lX/4nQ342IAAAampq5eXlAAA+n0+n04k2pyHi4uLmzZtnY2OD90N98ODB7NmzFW+JDGaE2SU19gejKmqEIGgdGLwMdB367wv8KuC/sNr+LJVq1L9ruwbvIQPc7yU1vOFayRUsufAxy3y0zPcLYato8YbRtc+0EaysrP7+++/CwkLxZG5yUlJSEhMTI34mOzu7lWZYlFTyx3pGllcLMAyAR4sBQw1Q/3Nti8XCn0bQaZR+8vfB/Y/Z36UcjNgKdPuKX1DJFf5xNWGWjYEyXfYLBhQKBfkgAMDMzGz//v0HDhwwMDAg+dsvLy+PjIwUP5Ofn0+I6JoMAuGm20k1vIZCUCVXuPTip/QTo+SXsSIQYjff5PAEjcxCqriiiC8lzaggbhgmkymRc2Rubq6YLCSyoaSkdO3aNQ6Ho6WlRbQtbYjdfqxKrlgN0S/uoON3K5B8IbbqymdHK32mHCIQzrnnGY0mjmIYFhhXMHmwnsyfXllZKX5oZGTUNn2QQqH4+PiUlJS0b9/MQu02SEu9gi/E7r7N4TVWulBWzX+TXNrCZzVAcm4lTTzKlrBB9vt/v3K+jfpreIJIthzNQAAAKBQKioKKBMPAxZeZDQ9GAQAiEQj5VCg/M/I43Eru9zbkxoKvL8DXFyA9Aj9XUSN8yyqRnxkICIqCTaKlM8LE7Aq6EqWa/99xMQtQ/1uVFnxrDMYXYO9TSof2kNffpqxaQBVPzs6OBsXsf7/Hvg1ReUKMUy2Qkw2IVkdhYaG6+reWC3w+n9jUteaRXVojuRaSHQ3KcwAAgMYAhsPguQqu4B2bI7/c6bJqAZ1GqeGLncqLB6WpAACgpIybIcJAUQVPTjYgWh0ikUjcBwUCQY8ePRRvRksDYa0IFAWKWf9+LxaBuAJRWZUcI5CeJpMv/lnQZwb4acS/3wuqwZ3p8FsVBq1TXb1/ZciHDx98fX2trKwmTZp069atuLg4V1fXHyk7/0dCW1s7OTkZPzx27FhrVFwrqxIoSZTo5X8EpV8BAICu+i0CibCicjmK3+pqMLkS66ID5krsEQIAlGjUrtqN9FxrIfHx8Xfv3u3Xr9+MGTMePnz45s2bmTNn4qVNCFJBpVIzMr7lT12+fPnSpUuKN6NpgdDb2zsjI6Oqqmrnzp1bt27lcrm9zG34ArExZt9Z35JleJXg3gz4rSqD1klLjhHISEdFlUmramx1iEoBv5jJeINQgry8PA8Pj+nTp48cObJjx45btmxZtGjRrVu35PpQRPOgUCjiK0i4HgKZOX78eF5enpKS0s6dO93c3DAM62psyhV8Xzc9cJ7EHiEAgE6jGMozAmmqKnXRVk7Jr2r4MhUGZVTvFskrNkpWVtaOHTtcXV1//fVXBoPh4eHh7OxMfv3bNou4D6qqqhJiQ9P2CJ2dnbds2ZKbm6uqqrpz585Tp069fvaEqdT4TTAMyDxFRYK1441VGbQGLqBQKD311fp2kW/uHGyyRafT27VrZ2trW1JSoqcn+7wARJtl4cKFW7duTU1N/fz5s4aGxv79+9+9eqkvxTqHMp0m71Hgxknd1JQb8kEAgE475rAe8jVj/PjxVCpVIBBoamqOHz++srKyQwf5PhHR2mlaIFRXV9+7d6+enh6FQmnXrp2vr+/48ePX2BmrMhv676dSQL+u6r06N9L/uoX8OcG4i7YynUYFKu1hNT3+fKCqDahKqgyqz/IBcrUBAMDn893d3Tdv3gwAyM7OPnTo0M6dO+X9UETbQV1dfevWrcbGxvn5+fr6+gAAGo223t644QhEoQB9Laa1iXwTKBaO7NpLX50BR8aqHQFNsohNlUm7tmKAvFP6hULh7t27V61apaSkVFRUtGPHjj179sj3kYhWTtMCIYvF2rlzZ0lJSXV19Y0bNzgcjqur61+TuulrMZVoFKCiDWhiEYhCAao6gEJTZdAuLZV7BGLSqWFbrX/qqKLqcAIYWH17QYmpNO2auu5PfusGyXs6CABwdnYuLCwMCgpKTk6eOHFihw4d7ty5I++HIprB8OHDJbpjz50799SpU0TZIyUsFuvIkSPJycmmpqbx8fFCoZDP5y8d/ZOpntq/EUhN51vC2n+o0GnXVg6UdwSiUSnBbpZmBupqTBqYdA5o9xR/SY1Ju7Z8oI2p3LMZly9fzmKxIiIi2Gy2nZ2dhobG7du35f1QRDMYMGBAeHi4+JnJkycTIs3YNNFtX1/fDx8+mJqa2tra7t+/38DAQFdXd/78+VklNb/septTyq36Pn+aTqMqM6j+fw5qYdMl6anhiw48SjkWmMoXYHQlilCEcQWiyYP1Djj1kusGCQKhGK5evZqYmDhw4EBHR8ebN2/Gx8dPmzZt8ODBRRW8UbvfpeRXSdQwKNEoTDr1xgpzh0EK0kDnC7HjgakHHqdUcYUMJapQhHH5onEDdA459zLtpKYYGxCIJiGz7hM1fNF+f/bxoK98oYhOowpFGE8omvZzp/1OvQzaKzoHAcNAcm5lVkmNliq9t4G6XIuIEQiSwBdixwJSDz5JqeIK6TSqCMO4fNFEc92Dzr266RKQg8DOq8oorm6nrNTbQF2lwf17BIJYZNyGCcMAK68yq6SmgxqjV+f/1moIIjc318vLq6qqatu2bbGxsXfu3Dlx4gSB9iBIS2JiYmJiop2dnbKy8osXL6qrq+3s7FqFbnidpORXZRTVaKoq9eqsLg8xM+lhsVjXrl0TCoXu7u7+/v6RkZHm5uazZs0i0CQEOWGz2XFxcePGjVNTU4uIiCgpKbGzs1OYVqqMnYRCAaad1Eaaafc3bEdsFAQAcDicv/76y8zM7NWrV9ra2tXV1cTagyAnJSUljx8/ZjAYmzZtCg4ODg8Pz8nJIaSYSVZ001X9xazDwJ80iI2CAIDq6upt27apqqomJCR07dp1//79d+/eJdYkBAnhcrl+fn7t27dfs2ZNVFTUgwcPuFzu8ePHFWZAax3zSkPPnj2fP39++fJlKysrMzMzXOAcgRCnffv269evHzx4MJ/P79Gjx6tXr2JiYqytrYm260egX79+t2/ffvHiRc+ePS0tLcvKyogqFEOQGSaT+eeff1paWopEIkNDw48fP4aHhw8fPlxhBsi4MS+p4PP5EydOZDAY/v7+c+fOJdocBHmpqamB9en//PPP0qVLeTze48eP+/TpQ7RdrR4+n+/s7CwQCMLDw83NzXft2nXgwAGijUKQEYFAsG3btk2bNkVHRzs5Oeno6Pj7+1tZWTX+k7LgR54Rvn792s3N7e7du6NGjbp16xabzT527BjRRiFIR15e3pgxY3r16pWUlESlUuPj4wsLC1uFygz5efr06ZYtW169ejVo0KBZs2bB9S60SYGQgMPh2NnZGRgYsNlsKpX6+fPngoICRbYUlnGyDALR6igvL09KSgIAqKmpwXDI4/H69pUUyUQgEHKiqqoKtgJVVlbu06dPamoqh8MZMGCAwvopokCIQCAQiDbNj7w0ikAgEAhEo6BAiEAgEIg2DQqECAQCgWjToECIQCAQiDYNCoQIBAKBaNOgQIhAIBCINg0KhAgEAoFo06BAiEAgEIg2DQqECAQCgWjToECIQCAQiDYNCoQIBAKBaNOgQIhAIBCINg0KhAgEAoFo06BAiEAgEIg2DQqECAQCgWjToECIQCAQiDYNCoQIBAKBaNOgQIhAIBCINg0KhAgEAoFo0ygRbQCiyfD5/IcPH2IY1q9fv169ehFtDgKBQLRu0Iyw9XHw4ME5c+bMmDHDz8+PaFsQCASi1YMCYSsjKSlp165dW7ZsIdoQBAKB+EFAgbA1gWHY8uXLTU1NV69eXftVHo8XFhZ27do1Pz+/2NjY8PBwxVuIQLQRysvLQ0NDCwoKiDYEIQPQHmFrwsvLKywsLCIigk6nS7z0/v37GTNmZGZmdu/evbi4uKCgwMDAIDMzkxA7EYgfD19fXxaLpa2tvWjRIgAAi8UaM2bMgwcPJk+eLH5ZYmLiw4cPXV1dO3XqhJ8UiUQHDx4cMmTI8OHDFW03QgrQjLDVkJ2dvWnTpuXLl9vY2Ei8VFZWZm9vr6WllZqampiYmJ+fv3DhQkKMRCB+SGJjY2fNmrVv377Dhw83fOXHjx/d3NwkxqBCodDNzS04OFieNiKaD5oRthpWrlyprKy8e/fu2i/du3cvPz//3r17Xbp0gWcoFIpirUMgfliEQuHixYttbGz09fXj4+OJNgche9CMsHVw584dPz+/48ePa2lp1X718+fPFArl559/VrxhCMQPz5EjR+Lj48+ePdvA+PL9+/f37t3j8XiKNAwhK1AgbB1s27YNALBs2bIOHTp06NChc+fOAICdO3cOHToUAKCqqkqhUKhU9NdEIGTM169fPTw8Nm7c2Lt37/quiY6OHjly5KdPnxgMhiJtQ8gKtDTaOnB3d6+srMQP+Xz+77//bmdn5+zsDAAwMTERiUQxMTHW1tbwgqSkJGIMRSB+IDAMW7Jkib6+/ubNm+u7Jjs7+7fffpswYcL27dsVaRtChqBA2DqYPXu2+GF1dfXvv/9uaWk5ffp0AIC9vb22tvby5cvPnj2rq6t7+vTpf/75x8DAgCBjEYgfhIsXL4aEhISEhCgrK9d5QXV19eTJkzt16uTj44M25lsvaDHtR6BDhw5Pnz6l0WhDhw7t0aNHXFycvb090UYhEK2b3NzcDRs2uLi42Nra1ndNQkLChw8fLC0tVVVVFWkbQragGWGrREVFBcMw8TOWlpZRUVE1NTUYhqmoqMybN692rSECgZAeb2/vkpKSa9euXbt2Tfw8hUJ5//49jUYDAFhYWOzatWvTpk22trZTp04FAKirqwMAysvLxX+koqICAIBckrSgQPhDgS/gfPnypXv37sQag0C0ambPnj1kyBDxM3v27GGz2RcvXuzRowebzYYnN27c+Pr164ULF5qbmxsbG/fp0wcAEBAQMGrUKPwHAwICAAANpNsgiAUFwh+BqqoqOzu7KVOmDBw4kMFg+Pv7v3371tvbm2i7EIhWTPfu3SVGk+fPn8/NzZVYKaVQKBcuXDA3N585c2ZERIShoaGjo+Phw4ezsrJ++eUXBoMRHR196dKlgQMHOjg4KPYdIKQFBcIfAYFA4ODgcPv27b1799bU1Jiamp46dWrx4sVE24VA/FDo6up27doVfs9kMrt166ampgYA0NHRuX79+oIFC7y8vFatWnX16lVzc3NfX9+goCAej2dkZLRkyZJt27ah4grSQpHYakIgEAgEok2BskYRCAQC0aZBgRCBQCAQbRoUCBEIBALRpkGBEIFAIBBtGhQIEQgEAtGmQYEQgUAgEG0aFAgRCAQC0aZBgRCBQCAQbRoUCBEIBALRpkGBEIFAIBBtGhQIEQgEAtGmQYEQgUAgEG0aFAgRCAQC0aZBgRCBQCAQbZr/A/pfalWeKZBLAAACTnpUWHRyZGtpdFBLTCByZGtpdCAyMDI1LjA5LjUAAHice79v7T0GIBAAYiYGCFABYjUgbmBkd9AA0szMMJoNQrPAaHRxXOqoRbMzQM1nUADSLFCKESzKyEhvmpuBMYOJgTWBkSmBiTmDiYlTgZklgZlfgYU1g4lVSIGNPYOJjTuBnSOBXYSBgzODiUOUgZMrgYs7g4lLjIFbmIGHN4OJRzCBR5yBly+BV4KBjz+DiU+SgV8gQUAwg0lAikFQmkFIhkFIlkFYjkFYnkFYgUFEkUFEiUFEmcGJGegENgZWFmYmRjY2bi5ODnY2HkEBfj5e8V2MQMcxwOJy13sjh9MTQw6AOHdiVtt/PGYLZi/XzrZ3sT28H8RODvxln8f5aB+I/TczykHRwAcs3uhb6hBuwg9W/+tHzIFJH16C1RjNljrQUatgD2LPTC/db72kyw7E7rjJub9/+VSw3uAVJ/fbGAqA9ZaKOh04dOYTWHzyD0P7feYVDiC2+N91+7xuGYHZ2jULbL8nrwWbGdV1xL5lTj2Y/XWRisNUF1awmoL9nA5PchLB7DT3lQ5G/1XB5s+zbDpw9bAdmL2tzvTARncWsJpNRgL7N7eXgc3xfjh5f1dkDlgNy9UeO3aDVWA1Qf7392+6HgJmbw2esd/zLC9YvcKnBofSz3fAbAtZHwed1MlgNZo7tzkcEG4Em+NVdsqh2vkg2F+RL5YckD7zA8w+KWp/4NyEBrCapozFB5YKNEHYlzsORJ64BDbz4qxd+9cba4HNtEoKOPA/PwHMFgMAi5KsdvIFRqUAAAMzelRYdE1PTCByZGtpdCAyMDI1LjA5LjUAAHicfVZbbuRGDPz3KXQBC803+3NtL7JBsGMgcXKFRf5zf6TYM5a0CBF5BpCoaj6qSI71x9NW1+9vv/39z3Zc/PYE+/ifz5xz+0vGGE/ft7rZXr7+8utte/348vJpeX3/8/bxxya+SeIM/n7Gfvl4//5poe228R7TMnl7ll2YBBHGPtZ1HmUAaRfJCNqeeZ9u7NEAZXuFdVK62PZMu7HJ6IAKIO3TzMfcnseuShrZAA2hZVfioAWMGXXkv0CHR9lTOM0rR5awGA0w4BHFKkqYvmKbT9IGmXAJT0oWzHjvOEHWAGcBUbelUwE5RmTnkcYdaZIaXgwRXOKugVJBQWGY6lwVgcvRhUfM1yppCIPPOpTTM9oE7gr5NILS4At5Du6Aes9UlCJstYlaSIe05RIVi3O1ialTC/QlugEQFXLO1GxjB4CMylnVcTOSBXQ0wFxAMrxH86HdGH3SASdEt51zTl/dYVNidEheCuk+QJHoavYcMbp6mO79kWJps/J0sN5pie75tmR34SXLdEjZtTtLIdEWZgISIapPxrh1UAW0es2GC4oTggBtnnb3GWATQWvOzXuX/og+0WgFQBb9aHAACZa4RhiiBjtlN22cAMpOOTkFJ3w4jRY4AbQ9aBAashQQce8kwufb5qjCp65utzEzukYSqnoMbEpRA+hUCe7WhyyJSm1sDVvxMZfRLjl5eB1uWG6rWdSkTaA0AgAiJXqQdlcPaeM/RMKQpYGN3Q3VdSLJEkmwDmxOq64fDuUb5Nfb20/7/v4L8PJ+ezt/ARRfPvd8Pfq5zRUmOXe2lunczIRHuizgAti5Zpfzw1nWy9P3rECYy8tmVJzJy/qrZzpyo0out/MAYuOZ6bKxtNxeFxPBwKeLSpHK92XTaIW6hEkcgoWPsku3is7XHaEr/sEF08Ny5I+Rp7LwiSk6kSQdFHGxDQufFiuMXWjDSFJZ+LSsnGGOy4xplXKJvnIOYC9DQ1UcHxgMB1VxfFAqfLfItdWpypWDMdGH5WBM7GE5GJPKGZTIUbvEw3JUKvmwHHVVt157s54//7PB/dO/xpSz/SZaXy4AAANjelRYdFNNSUxFUyByZGtpdCAyMDI1LjA5LjUAAHicdVXBahwxDP2VQqGkMDGWZEmWc+ylp5B76SGE0hbapoQc+/GVnN0dz6xnIRNWo2c9Pb+Vvnz++gRPN18+f/14ft4+YTzo/P1T/98f8ffnKvJEPe4v7tcsPJ8H/nz37yYnVM5CCyeCSrLc5STGYLZQqlrLcnebExVQ5QWTWmH1ECQtRcBzCFk6KqvDoS6QiONVJIHlzEtOzGTgR0Fif6keMaulYkSM2bO9SClQtMZRBlWIFz+AkQk6BaZaVCIG4u97zLgK4BI95DNVFkItUUEo58jDVIA10gQdyj1UybHeUBapkVRCBycESYpo6ey91RpNC1sEKHmr3rHLVd++F0btzFl6ex7KoUewrCZV+8HZ9aTFa1rNEeEETNSTrJBS52PVdY1kUKlvSVkYekgKUyfEyoXjJDwFgEMeLyuGpqEmUVXtxSTI3cWVsXd4GzcVHP1WCb31OIU0LpMTVrMOYqMeUMhvpZ2nX644VlyDoJAt2vLeAV1Tb15NJVtUcpFL6JUrxqU5HM3FcK9gGMMzgOONS4iI/RSohnUpSbJ74OPy+Pr8++Hl+W/L6eHx5fXn469PPx5fvn9L+cP7IndAnFGHj+cF5P759Zt302CHuj3BAExx/RRPvMAAGx7ASMoWhgOMGh1w1A3H6nkrx9LKvFi2LUciz1yrlcbzaoV4rMaet1bjJgetSdm1JkMxbjotVmRTrKLnXVDU6pyhKW/vrA4Mmx3IoVvtXQ5bYdZgb5ATDsmNtb3qwSK1wd4iZxnrhqRE4soSGuw9clIEfP5sFIHBI9pg75Gz/j7odjRHl+QGe5ucZHGWG54Wmesd+E9g7pMZcHAKOYG9VY6Bg1eIG+zNcupRxEZHUyRecC7w3C4TRQe/gN/EgWMmkg6OAWs4nymTBnGwDOWGc8/MgINpCBrOTTMDDq4hH0gHtrl2Nw62Aac6t821vXFwDXrBvWvOBdF2kxMH13iHB/OFKm1wPfOCc57zAVMRbOQZeRcUN9xb5lTNCDfDQiJzZVkb7j1z8prv1e1owsEy6JN7bplrHI17yBodLCL1RbuVkwbH+LeDTTRZDjTuIml0sIyu1gONc8ZxB+voevfRuI6844N9dL39xoWk//4DTIrJzZ7U44oAAAJXelRYdHJka2l0UEtMMSByZGtpdCAyMDI1LjA5LjUAAHice79v7T0GIBAAYiYGCFABYjUgbmBkd9AA0szMMJoNQrPAaHRxXOqoRbMzQM1nUADSLFCKESzKyEgpzQk2Do3CpZibgTGDiYE1gZEpgYk5g4mJU4GZJYGZX4GFNYOJVUiBjT2DiY07gZ0jgV2EgYMzg4lDlIGTK4GLO4OJS4yBW5iBhzeDiUcwgUecgZcvgVeCgY8/g4lPkoFfIEFAMINJQIpBUJpBSIZBSJZBWI5BWJ5BWIFBRJFBRIlBRJnBiRnoBDYGVhZmJkY2Nm4uTg52Nh5BAX4+XvFdjEDHMcAictd7I4fTE0MOgDh3YlbbfzxmC2Yv1862d7E9vB/ETg78ZZ/H+WgfiP03M8pB0cAHLN7oW+oQbsIPVv/rR8yBSR9egtUYzZY60FGrYA9iz0wv3W+9pMsOxO64ybm/f/lUsN7gFSf32xgKgPWWijodOHTmE1h88g9D+33mFQ4gdsF+TocnOYlg9tdFKg5TXVjB7KiuI/Ytc+rB5mvXLLD9nrwWzBb/u26f1y0jsJo095UORv9VwebPs2w6cPWwHZi9rc70wEZ3FrCaTUYC+ze3l4H1ej+cvL8rMgeshuVqjx27wSqwGgtZHwed1MlgtsKnBofSz3fA6rcGz9jveZYXzA7yv79/0/UQsBrNndscDgg3gs3xKjvlUO18EOyvpozFB5YKNIHFI18sOSB95gdY/KSo/YFzExrA4k2XOw5EnrgENvPirF371xtrgc20Sgo48D8/AcwWAwDUYqxzURTi1QAAAzd6VFh0TU9MMSByZGtpdCAyMDI1LjA5LjUAAHicfVZbbhsxDPzPKXQBC+Jb+myctCmKOECb9A5Ff3t/dCg7uxtU6NoGdukRH8Mhbf39dlfy+v7w7defsl38cAd7+897jFF+Smvt7rnkTbl//PL1Us6vn+7fLeeXt8vrjyJepOMMXh+xn15fnt8tVC6FawzrnctJqjAJIrTa5rUfZQCpivQIKieuw409FkApZ1gHdRcrJ6rGJm0FVACpDjNvo5xaVSWNvgAaQktV4qAJjBF55F+gw6PULtzNM0eWsGgLYMAjilWUMHzGNh+kC2SHS3hSsmDG944TZAvgSCDqtu6UQI4WfeWR2hVp0jU8GSK4xN0CSgkFhWGqY1YELtsqPGKes6QmDD7zUB/eY5nAtUM+jNBp8IU8G6+A2SGuZCxIWtBJRgtWQJtAUKTquGkdJ5ZAnz03MY2MOEbXvgwdM0dQKM6pO1NH3gtgv7IpShE2pawWS+RA061yH8OnOmxItFWWPDuktYEi0Sn23qKtnDJd9dHFuo0s3cH6qpdQz9NsuwvPtgxHK1dyZ0kkZGEmqBlN9cEYtxVUAU2tWXNBcULga5mnASiV+uAuRas3p7aaDHYAUTrnXIL4YKe+BMYtywFBZm7Idj1C3K/IgDhQcO4Y83U5A0irQY0gyOyAiPuqRXg/FYcnHzrVbm30WClJqHwuJ6vU3LCGZlvVZAnlG9QkWYTXoRK82jQiCU1hYMHY9IoRjuU+zB4BgCZ1aJCqq4csndqVKAi4G9iobqhuRZRklyAMNhvDcjSbo/ML5OPl4cO+v/4C3L9cHvZfAMWH9z2fj75vc4VJ9p2tado3M+GRDgs4Abav2el8c9bzy933yECYy8NmVJzph/WXz7TlRplcL/sBxMYz02Fjabo9biaCgXcXmSKl78Oq0Qx1CNNxCBbeys6+ZXQ+7gid8TcumG6WLX+MPKWFd0zSiSRpo4iTbVh4t1hi7EAbZpLSwrtl5gxzHIZMs5RD9JlzAHsYGsrieMNgOCiL441S4atFjlKnLFc2xkRvlo0xsZtlY0wyZ1AiW+0SN8tWqfSbZasr1XrUZj6//7PB/d1fzRG0UIFMjdIAAAOAelRYdFNNSUxFUzEgcmRraXQgMjAyNS4wOS41AAB4nHVVwWocSQz9lcBCsKFdlKRSqVRzDIScQu4hB2PCZmF3vRgf8/H71DPurm53D3YPVf1e6Ul6U/r+5ccTPd19//Lj/u358MTxkLf1p/l7fsT/v7Hz+R5/n++fZN7B1tf1Pb+dRHh++H2XE5vmKpMmoSZ1uuRUXcl9ktSslenCiZSFsFZmtukiiZpzm0qquVIA1LjoxCm3AF5KyuxaJ0rGldp0oaSixaac3FtpHDuuCvJDTqVQsRaBnVoVnR4AZ42D8FalFauxR9XjLOy5tko8hfY8SwxcFbYSEarkHDhOhSAMW5VB1XmrCbihtNYWoBL5QxClWqrFUZSKWotsqnpsSGLkguxzbdd1UZQhlCvOvW5lCRDIzWuz+eCMOsqEmN5y7GiiXBUpx6uCV/OWisw8L2IyS/SmNoPIapsFqSnKi1c8s1AKjfIgbHV2i2qKNLM5WA1x6ImhSi0wwqER3RRG6nGKGCJcNHFzn0nqMm8Y5as+iIIVKrgVNQgJ2dvc+kKMmiJ5c6vZo23ogUSboUJLGCIUes6K0kPotUQ5SSGzqKt5ifhAWQl8SFTEu58eX5//+fby/F/P6dvjy+tfj39/+vX48ufPlD/+UeqFxDLb8mmABePr8+tPJNVpR3q4sSi78fABcKERdz6hSSlbGg806XIi0XWUSMCtGksvx8Gyl61GBnKNVroeRysyFgQV0SGa9nqSGpq6Ta0OwbTbYTBQxtTMgFtY0tuxQvyURpIDtyrsflIOk41AMiAXmnfa++PG44xg21YPFmmd9ha5iVTNW5E0WISp094jt4pQrWOzOYALzzrtPfJWf82ySS+Qa7zcaW+TW1m0CTq3BmyBXHuAn8CxT46Ig1OkdNpb5Zw4eEW0094stxzjzh1zrIFciFgd+0XYNyUN3GpNFGrvmCtNvY00qwFced75zDHsLj62ggfHSO58cqscEAfPCJYn98oBcTCN4Jo5cw0u8LGiAVxThNJj05Dt7wgePBPH7D2zCBXfupQHz2B1crtIoz1vsAx0Hl8vyM5HnRTAhaad94a5hfPdPV0DucpEuseWwRCX7XXGg2W4dTkeQu95Mo4h73LiGEPjR6ElkAsRqxPDHMwGGUdR7XIyi95NBxmvGfBOptH70SfjNELGJ+Po/fAb55H9/h/t18bLkoGoiwAAAABJRU5ErkJggg==", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "mapping" ] }, { "cell_type": "markdown", "id": "674f34f1-27e5-433c-96cd-14ece7d12eaf", "metadata": {}, "source": [ "## Creating the AlchemicalNetwork\n", "\n", "The **`AlchemicalNetwork`** contains all the information needed to run a complete RBFE campaign. Each `Transformation` represents one alchemical mutation between a pair of ligands in a specific environment (either **solvent** or **protein-membrane complex**).\n", "\n", "In this example, we loop over all ligand mappings, and for each mapping, we loop over the legs (`solvent` and `complex`) to create the corresponding `Transformation` objects.\n", "\n", "The default settings apply the **AMBER lipid17** forcefield to the lipids. Additional lipid forcefields from OpenMMForcefields, such as AMBER lipid21, can also be added (note: not yet tested with this Protocol).\n", "\n", "In this example, we use **`_adaptive_settings`** in `RelativeHybridTopologyProtocol` to generate protocol settings based on the system composition:\n", "\n", "- **Complex leg (protein-membrane)** \n", " - Uses a `MonteCarloMembraneBarostat` for appropriate membrane pressure coupling. \n", "\n", "- **Solvent leg** \n", " - Uses a standard `MonteCarloBarostat` for bulk solvent pressure regulation." ] }, { "cell_type": "code", "execution_count": 22, "id": "92790417-ad0f-41f3-958a-ce3a64286c86", "metadata": { "scrolled": true }, "outputs": [], "source": [ "from openfe.protocols.openmm_rfe import RelativeHybridTopologyProtocol\n", "transformations = []\n", "for mapping in ligand_network.edges:\n", " for leg in ['solvent', 'complex']:\n", " if leg == 'solvent':\n", " sysA_dict = {'ligand': mapping.componentA,\n", " 'solvent': openfe.SolventComponent()}\n", " sysB_dict = {'ligand': mapping.componentB,\n", " 'solvent': openfe.SolventComponent()}\n", "\n", " if leg == 'complex':\n", " sysA_dict = {'ligand': mapping.componentA,\n", " 'protein': membrane_protein_user_box}\n", " sysB_dict = {'ligand': mapping.componentB,\n", " 'protein': membrane_protein_user_box}\n", "\n", " # we don't have to name objects, but it can make things (like filenames) more convenient\n", " system_a = openfe.ChemicalSystem(sysA_dict, name=f\"{mapping.componentA.name}_{leg}\")\n", " system_b = openfe.ChemicalSystem(sysB_dict, name=f\"{mapping.componentB.name}_{leg}\")\n", "\n", " prefix = \"rbfe_\" \n", " protocol_settings = RelativeHybridTopologyProtocol.default_settings()\n", " protocol_settings.protocol_repeats = 1\n", "\n", " protocol = RelativeHybridTopologyProtocol(\n", " settings=RelativeHybridTopologyProtocol._adaptive_settings(\n", " stateA=system_a,\n", " stateB=system_b,\n", " mapping=mapping,\n", " initial_settings=protocol_settings,\n", " ),\n", " )\n", " \n", "\n", " transformation = openfe.Transformation(\n", " stateA=system_a,\n", " stateB=system_b,\n", " mapping=mapping,\n", " protocol=protocol, # use protocol created above\n", " name=f\"{prefix}{system_a.name}_{system_b.name}\"\n", " )\n", " transformations.append(transformation)\n", "\n", "network = openfe.AlchemicalNetwork(transformations)" ] }, { "cell_type": "code", "execution_count": 23, "id": "2086b75d-e03b-4135-bb62-1b315fad50e1", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "{'barostat': 'MonteCarloBarostat',\n", " 'barostat_frequency': {'unit': 'timestep', 'val': 25.0},\n", " 'constraint_tolerance': 1e-06,\n", " 'langevin_collision_rate': {'unit': '1 / picosecond', 'val': 1.0},\n", " 'n_restart_attempts': 20,\n", " 'reassign_velocities': False,\n", " 'remove_com': False,\n", " 'surface_tension': {'unit': 'bar * nanometer', 'val': 0},\n", " 'timestep': {'unit': 'femtosecond', 'val': 4.0}}\n" ] } ], "source": [ "# A solvent leg transformation uses the `MonteCarloBarostat`\n", "transformations[-2].protocol.settings.integrator_settings" ] }, { "cell_type": "code", "execution_count": 24, "id": "ee608914-b917-411f-a726-14000d426156", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "{'barostat': 'MonteCarloMembraneBarostat',\n", " 'barostat_frequency': {'unit': 'timestep', 'val': 25.0},\n", " 'constraint_tolerance': 1e-06,\n", " 'langevin_collision_rate': {'unit': '1 / picosecond', 'val': 1.0},\n", " 'n_restart_attempts': 20,\n", " 'reassign_velocities': False,\n", " 'remove_com': False,\n", " 'surface_tension': {'unit': 'bar * nanometer', 'val': 0},\n", " 'timestep': {'unit': 'femtosecond', 'val': 4.0}}\n" ] } ], "source": [ "# A complex leg transformation uses the `MonteCarloMembraneBarostat`\n", "transformations[-1].protocol.settings.integrator_settings" ] }, { "cell_type": "markdown", "id": "e504e3b3-801f-4501-9207-225f8e462675", "metadata": {}, "source": [ "## Writing the `Transformation`s to disk\n", "\n", "We'll write out each transformation to disk, so that they can be run independently using the `openfe quickrun` command:" ] }, { "cell_type": "code", "execution_count": 25, "id": "4e6cbea7-c357-432f-a115-a9f3b9c4e453", "metadata": {}, "outputs": [], "source": [ "import pathlib\n", "# first we create the directory\n", "transformation_dir = pathlib.Path(\"transformations\")\n", "transformation_dir.mkdir(exist_ok=True)\n", "\n", "# then we write out each transformation\n", "for transformation in network.edges:\n", " transformation.to_json(transformation_dir / f\"{transformation.name}.json\")" ] }, { "cell_type": "markdown", "id": "6f0085a5-d5a9-43b9-bb38-43deca490d2d", "metadata": {}, "source": [ "Each of these individual `.json` files contains a `Transformation`, which contains all the information to run the calculation. These could be farmed out as individual jobs on a HPC cluster. For details on running them, please follow from the section on running simulations in the CLI tutorial (https://docs.openfree.energy/en/latest/tutorials/rbfe_cli_tutorial.html)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.13" } }, "nbformat": 4, "nbformat_minor": 5 }