{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# Start-to-Finish Example: Setting up Polytropic [TOV](https://en.wikipedia.org/wiki/Tolman%E2%80%93Oppenheimer%E2%80%93Volkoff_equation) Initial Data, in Curvilinear Coordinates\n",
"\n",
"## Authors: Zach Etienne, Phil Chang, and Leo Werneck\n",
"### Formatting improvements courtesy Brandon Clark\n",
"\n",
"## This module sets up initial data for a TOV star in *spherical, isotropic coordinates*, using the *Numerical* ADM Spherical to BSSN Curvilinear initial data module (numerical = BSSN $\\lambda^i$'s are computed using finite-difference derivatives instead of exact expressions).\n",
"\n",
"**Notebook Status:** Validated \n",
"\n",
"**Validation Notes:** This module has been validated to exhibit convergence to zero of the Hamiltonian constraint violation at the expected order to the exact solution (see [plots](#convergence) at the bottom). Note that convergence at the surface of the star will be lower-order due to the sharp drop to zero in $T^{\\mu\\nu}$.\n",
"\n",
"### NRPy+ Source Code for this module: \n",
"\n",
"* [TOV/TOV_Solver.py](../edit/TOV/TOV_Solver.py); ([**NRPy+ Tutorial module reviewing mathematical formulation and equations solved**](Tutorial-ADM_Initial_Data-TOV.ipynb)); ([**start-to-finish NRPy+ Tutorial module demonstrating that initial data satisfy Hamiltonian constraint**](Tutorial-Start_to_Finish-BSSNCurvilinear-TOV_initial_data.ipynb)): Tolman-Oppenheimer-Volkoff (TOV) initial data; defines all ADM variables and nonzero $T^{\\mu\\nu}$ components in Spherical basis.\n",
"* [BSSN/ADM_Initial_Data_Reader__BSSN_Converter.py](../edit/BSSN/ADM_Initial_Data_Reader__BSSN_Converter.py); ([**NRPy+ Tutorial documentation**](Tutorial-ADM_Initial_Data_Reader__BSSN_Converter.ipynb)): Registers the C function for our \"universal\" initial data reader/converter: `initial_data_reader__convert_ADM_Cartesian_to_BSSN()`.\n",
"* [CurviBoundaryConditions/CurviBoundaryConditions.py](../edit/CurviBoundaryConditions/CurviBoundaryConditions.py); ([**NRPy+ Tutorial documentation**](Tutorial-Start_to_Finish-Curvilinear_BCs.ipynb)): Registers the C function for our \"universal\" initial data reader/converter initial_data_reader__convert_ADM_Cartesian_to_BSSN().\n",
"* [BSSN/BSSN_constraints.py](../edit/BSSN/BSSN_constraints.py); [\\[**tutorial**\\]](Tutorial-BSSN_constraints.ipynb): Hamiltonian constraint in BSSN curvilinear basis/coordinates.\n",
"\n",
"## Introduction:\n",
"Here we use NRPy+ to set up initial data for a [simple polytrope TOV star](https://en.wikipedia.org/wiki/Tolman%E2%80%93Oppenheimer%E2%80%93Volkoff_equation).\n",
"\n",
"The entire algorithm is outlined as follows, with links to the relevant NRPy+ tutorial notebooks listed at each step:\n",
"\n",
"1. Allocate memory for gridfunctions.\n",
"1. Set gridfunction values to initial data \n",
" * [**NRPy+ tutorial on TOV initial data**](Tutorial-ADM_Initial_Data-TOV.ipynb)\n",
"1. Convert ADM initial data quantities to BSSN-with-a-reference-metric quantities\n",
" * [**NRPy+ tutorial on BSSN initial data reader/converter**](Tutorial-ADM_Initial_Data_Reader__BSSN_Converter.ipynb)\n",
"1. Apply boundary conditions to BSSN quantities (including $\\lambda^i$, which is computed via finite difference derivatives and thus only defined in grid interior)\n",
" * [**NRPy+ tutorial on boundary conditions**](Tutorial-Start_to_Finish-Curvilinear_BCs.ipynb)\n",
"1. Evaluate the Hamiltonian constraint violation\n",
" * [**NRPy+ tutorial on BSSN constraints**](Tutorial-BSSN_constraints.ipynb)\n",
"1. Repeat the above steps at two numerical resolutions to confirm the convergence of Hamiltonian constraint violation to zero."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Table of Contents\n",
"$$\\label{toc}$$\n",
"\n",
"This notebook is organized as follows\n",
"\n",
"1. [Step 1](#initializenrpy): Set core NRPy+ parameters for numerical grids and reference metric\n",
"1. [Step 2](#tov_solution): Solve TOV equations to generate initial data\n",
"1. [Step 3](#adm_id_spacetime): Register C functions for interpolating TOV initial data onto our numerical grids, converting TOV metric quantities to BSSN, and constructing $T^{\\mu\\nu}$\n",
"1. [Step 4](#ham_const_output): Register C functions for outputting the Hamiltonian constraint\n",
"\n",
"1. [Step 5](#cparams_rfm_and_domainsize) Output `free_parameters.h`, which sets numerical grid parameters\n",
"1. [Step 6](#mainc): Register the `main()` C code for generating the `TOV_Playground` executable\n",
"1. [Step 7](#nrpy_infra_codes): Register the core NRPy+ infrastructure C codes that set up needed data structures, C macros, etc.\n",
"1. [Step 8](#output_compile_and_run): Output all C codes, compile `TOV_Playground`, and run at two resolutions\n",
"1. [Step 9](#visualize_output): Visualize output\n",
"1. [Step 10](#convergence): Code validation test: confirm numerical errors (Hamiltonian constraint violation) converge to zero as expected\n",
"1. [Step 11](#latex_pdf_output): Output this notebook to $\\LaTeX$-formatted PDF file"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 1: Set core NRPy+ parameters for numerical grids and reference metric \\[Back to [top](#toc)\\]\n",
"$$\\label{initializenrpy}$$\n"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:30.577985Z",
"iopub.status.busy": "2022-12-11T12:47:30.577297Z",
"iopub.status.idle": "2022-12-11T12:47:31.521663Z",
"shell.execute_reply": "2022-12-11T12:47:31.520748Z"
}
},
"outputs": [],
"source": [
"# Step P1: Import needed NRPy+ core modules:\n",
"from outputC import outputC,add_to_Cfunction_dict # NRPy+: Core C code output module\n",
"import NRPy_param_funcs as par # NRPy+: Parameter interface\n",
"import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends\n",
"import finite_difference as fin # NRPy+: Finite difference C code generation module\n",
"import grid as gri # NRPy+: Functions having to do with numerical grids\n",
"import indexedexp as ixp # NRPy+: Symbolic indexed expression (e.g., tensors, vectors, etc.) support\n",
"import reference_metric as rfm # NRPy+: Reference metric support\n",
"import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface\n",
"import shutil, os, sys # Standard Python modules for multiplatform OS-level functions\n",
"import diagnostics_generic.process_2D_data as plot2D # NRPy+: analysis of output data\n",
"import diagnostics_generic.output_yz_or_xy_plane as planar_diags # NRPy+: C code for generating output data\n",
"\n",
"# Step P2: Create C code output directory:\n",
"Ccodesrootdir = os.path.join(\"TOVID_Ccodes/\")\n",
"# First remove C code output directory if it exists\n",
"# Courtesy https://stackoverflow.com/questions/303200/how-do-i-remove-delete-a-folder-that-is-not-empty\n",
"# !rm -r ScalarWaveCurvilinear_Playground_Ccodes\n",
"shutil.rmtree(Ccodesrootdir, ignore_errors=True)\n",
"# Then create a fresh directory\n",
"cmd.mkdir(Ccodesrootdir)\n",
"\n",
"# Step 1: Set the spatial dimension parameter\n",
"# to three this time, and then read\n",
"# the parameter as DIM.\n",
"par.set_parval_from_str(\"grid::DIM\",3)\n",
"DIM = par.parval_from_str(\"grid::DIM\")\n",
"\n",
"# Step 2: Set some core parameters, including CoordSystem MoL timestepping algorithm,\n",
"# FD order, floating point precision, and CFL factor:\n",
"# Choices are: Spherical, SinhSpherical, SinhSphericalv2, Cylindrical, SinhCylindrical,\n",
"# SymTP, SinhSymTP\n",
"CoordSystem = \"Spherical\"\n",
"par.set_parval_from_str(\"reference_metric::CoordSystem\", CoordSystem)\n",
"rfm.reference_metric()\n",
"\n",
"# Step 2.a: Set defaults for Coordinate system parameters.\n",
"# These are perhaps the most commonly adjusted parameters,\n",
"# so we enable modifications at this high level.\n",
"\n",
"# domain_size = 7.5 # SET BELOW BASED ON TOV STELLAR RADIUS\n",
"\n",
"# sinh_width sets the default value for:\n",
"# * SinhSpherical's params.SINHW\n",
"# * SinhCylindrical's params.SINHW{RHO,Z}\n",
"# * SinhSymTP's params.SINHWAA\n",
"sinh_width = 0.4 # If Sinh* coordinates chosen\n",
"\n",
"# sinhv2_const_dr sets the default value for:\n",
"# * SinhSphericalv2's params.const_dr\n",
"# * SinhCylindricalv2's params.const_d{rho,z}\n",
"sinhv2_const_dr = 0.05# If Sinh*v2 coordinates chosen\n",
"\n",
"# SymTP_bScale sets the default value for:\n",
"# * SinhSymTP's params.bScale\n",
"SymTP_bScale = 0.5 # If SymTP chosen\n",
"\n",
"# Step 2.b: Set the order of spatial finite difference derivatives;\n",
"# and the core data type.\n",
"FD_order = 4 # Finite difference order: even numbers only, starting with 2. 12 is generally unstable\n",
"REAL = \"double\" # Best to use double here.\n",
"\n",
"# Step 3: Polytropic EOS setup\n",
"# For EOS_type, choose either \"SimplePolytrope\" or \"PiecewisePolytrope\"\n",
"EOS_type = \"SimplePolytrope\"\n",
"# If \"PiecewisePolytrope\" is chosen as EOS_type, you\n",
"# must also choose the name of the EOS, which can\n",
"# be any of the following:\n",
"# 'PAL6', 'SLy', 'APR1', 'APR2', 'APR3', 'APR4',\n",
"# 'FPS', 'WFF1', 'WFF2', 'WFF3', 'BBB2', 'BPAL12',\n",
"# 'ENG', 'MPA1', 'MS1', 'MS2', 'MS1b', 'PS', 'GS1',\n",
"# 'GS2', 'BGN1H1', 'GNH3', 'H1', 'H2', 'H3', 'H4',\n",
"# 'H5', 'H6', 'H7', 'PCL2', 'ALF1', 'ALF2', 'ALF3',\n",
"# 'ALF4'\n",
"EOS_name = 'SLy' # <-- IGNORED IF EOS_type is not PiecewisePolytrope."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:31.527182Z",
"iopub.status.busy": "2022-12-11T12:47:31.526896Z",
"iopub.status.idle": "2022-12-11T12:47:31.531572Z",
"shell.execute_reply": "2022-12-11T12:47:31.530922Z"
}
},
"outputs": [],
"source": [
"# Step 4: Set the finite differencing order to FD_order (set above).\n",
"par.set_parval_from_str(\"finite_difference::FD_CENTDERIVS_ORDER\", FD_order)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 2: Solve TOV equations to generate initial data \\[Back to [top](#toc)\\]\n",
"$$\\label{tov_solution}$$\n",
"\n",
"As documented [in the TOV Initial Data NRPy+ Tutorial Module](Tutorial-TOV_Initial_Data.ipynb) ([older version here](Tutorial-GRMHD_UnitConversion.ipynb)), we will now set up TOV initial data, storing the densely-sampled result to file (***Courtesy Phil Chang***)."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"`scipy` provides the ODE integration routine used in our TOV solver, so we first make sure that it is installed:"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:31.535729Z",
"iopub.status.busy": "2022-12-11T12:47:31.535451Z",
"iopub.status.idle": "2022-12-11T12:47:34.092351Z",
"shell.execute_reply": "2022-12-11T12:47:34.090919Z"
},
"scrolled": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\r\n",
"\u001b[1m[\u001b[0m\u001b[34;49mnotice\u001b[0m\u001b[1;39;49m]\u001b[0m\u001b[39;49m A new release of pip is available: \u001b[0m\u001b[31;49m23.0.1\u001b[0m\u001b[39;49m -> \u001b[0m\u001b[32;49m23.1\u001b[0m\r\n",
"\u001b[1m[\u001b[0m\u001b[34;49mnotice\u001b[0m\u001b[1;39;49m]\u001b[0m\u001b[39;49m To update, run: \u001b[0m\u001b[32;49mpip install --upgrade pip\u001b[0m\r\n"
]
}
],
"source": [
"!pip install scipy > /dev/null"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next we call the [`TOV.TOV_Solver()` function](../edit/TOV/TOV_Solver.py) ([NRPy+ Tutorial module](Tutorial-ADM_Initial_Data-TOV.ipynb)) to set up the initial data, using the default parameters for initial data. This function outputs the solution to a file named `outputTOVpolytrope.txt`."
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:34.098442Z",
"iopub.status.busy": "2022-12-11T12:47:34.098121Z",
"iopub.status.idle": "2022-12-11T12:47:34.582698Z",
"shell.execute_reply": "2022-12-11T12:47:34.581967Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"1256 1256 1256 1256 1256 1256\n",
"Just generated a TOV star with\n",
"* M = 1.405030336771405e-01 ,\n",
"* R_Schw = 9.566044579232513e-01 ,\n",
"* R_iso = 8.100085557410308e-01 ,\n",
"* M/R_Schw = 1.468768334847266e-01 \n",
"\n"
]
}
],
"source": [
"##########################\n",
"# Polytropic EOS example #\n",
"##########################\n",
"import TOV.Polytropic_EOSs as ppeos\n",
"\n",
"if EOS_type == \"SimplePolytrope\":\n",
" # Set neos = 1 (single polytrope)\n",
" neos = 1\n",
"\n",
" # Set rho_poly_tab (not needed for a single polytrope)\n",
" rho_poly_tab = []\n",
"\n",
" # Set Gamma_poly_tab\n",
" Gamma_poly_tab = [2.0]\n",
"\n",
" # Set K_poly_tab0\n",
" K_poly_tab0 = 1. # ZACH NOTES: CHANGED FROM 100.\n",
"\n",
" # Set the eos quantities\n",
" eos = ppeos.set_up_EOS_parameters__complete_set_of_input_variables(neos,rho_poly_tab,Gamma_poly_tab,K_poly_tab0)\n",
" rho_baryon_central = 0.129285\n",
"elif EOS_type == \"PiecewisePolytrope\":\n",
" eos = ppeos.set_up_EOS_parameters__Read_et_al_input_variables(EOS_name)\n",
" rho_baryon_central=2.0\n",
"else:\n",
" print(\"\"\"Error: unknown EOS_type. Valid types are 'SimplePolytrope' and 'PiecewisePolytrope' \"\"\")\n",
" sys.exit(1)\n",
"\n",
"import TOV.TOV_Solver as TOV\n",
"M_TOV, R_Schw_TOV, R_iso_TOV = TOV.TOV_Solver(eos,\n",
" outfile=os.path.join(Ccodesrootdir, \"TOVdata.txt\"),\n",
" rho_baryon_central=rho_baryon_central,\n",
" return_M_RSchw_and_Riso = True,\n",
" verbose = True)\n",
"\n",
"# domain_size sets the default value for:\n",
"# * Spherical's params.RMAX\n",
"# * SinhSpherical*'s params.AMAX\n",
"# * Cartesians*'s -params.{x,y,z}min & .{x,y,z}max\n",
"# * Cylindrical's -params.ZMIN & .{Z,RHO}MAX\n",
"# * SinhCylindrical's params.AMPL{RHO,Z}\n",
"# * *SymTP's params.AMAX\n",
"domain_size = 2.0 * R_iso_TOV"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 3: Register C functions for interpolating TOV initial data onto our numerical grids, converting TOV metric quantities to BSSN, and constructing $T^{\\mu\\nu}$ \\[Back to [top](#toc)\\]\n",
"$$\\label{adm_id_spacetime}$$\n",
"\n",
"First, we register all the C functions needed to read the TOV data file, interpolate the TOV solution to our numerical grid, and convert the TOV quantities to ADM variables and $T^{\\mu\\nu}$."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:34.632361Z",
"iopub.status.busy": "2022-12-11T12:47:34.632071Z",
"iopub.status.idle": "2022-12-11T12:47:34.776634Z",
"shell.execute_reply": "2022-12-11T12:47:34.775945Z"
}
},
"outputs": [],
"source": [
"import TOV.TOV_Ccodegen_library as TOVCL # NRPy+: TOV C codegen library\n",
"TOVCL.add_to_Cfunction_dict_TOV_read_data_file_set_ID_persist(interp_stencil_size=12)\n",
"TOVCL.add_to_Cfunction_dict_TOV_ID_function()\n",
"TOVCL.add_to_Cfunction_dict_TOV_interpolate_1D()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next we register the C function for our \"universal\" initial data reader/converter `initial_data_reader__convert_ADM_spherical_to_BSSN()`, provided by [`BSSN.Initial_Data_Reader__BSSN_Converter`](../edit/BSSN.Initial_Data_Reader__BSSN_Converter.py) ([documentation](Tutorial-ADM_Initial_Data_Reader__BSSN_Converter.ipynb)).\n",
"\n",
"This function first calls `TOV_ID_function()` to read the initial data from file. Spacetime data are given in terms of ADM quantities, and fluid quantity data are given in terms of the stress-energy tensor. Then the function converts the ADM quantities to BSSN in our desired reference metric, and performs the basis transform on the stress-energy tensor as well."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:34.780636Z",
"iopub.status.busy": "2022-12-11T12:47:34.780486Z",
"iopub.status.idle": "2022-12-11T12:47:37.076931Z",
"shell.execute_reply": "2022-12-11T12:47:37.076250Z"
}
},
"outputs": [],
"source": [
"import BSSN.ADM_Initial_Data_Reader__BSSN_Converter as IDread\n",
"IDread.add_to_Cfunction_dict_initial_data_reader__convert_ADM_Sph_or_Cart_to_BSSN(input_Coord=\"Spherical\",\n",
" include_T4UU=True)\n",
"\n",
"IDread.register_NRPy_basic_defines(ID_persist_struct_contents_str=TOVCL.ID_persist_str(), include_T4UU=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 4: Register C functions for outputting the Hamiltonian constraint \\[Back to [top](#toc)\\]\n",
"$$\\label{ham_const_output}$$\n",
"\n",
"In our start-to-finish C code below, we will evaluate the Hamiltonian constraint violation of the TOV data on our numerical grid. The idea is, this quantity should drop to zero with increasing numerical resolution, at a rate consistent with our finite differencing operator's truncation error.\n",
"\n",
"The Hamiltonian constraint $\\mathcal{H}$ [as implemented in this corresponding NRPy+ tutorial notebook](Tutorial-BSSN_constraints.ipynb) depends on the Ricci tensor. As is done when we evolve our initial data forward in time, we first compute Ricci, then $\\mathcal{H}$, as this breaks up the large symbolic expression for $\\mathcal{H}$, decreasing the time needed to generate the C codes. Here we register functions for both computing the Ricci tensor and $\\mathcal{H}$:"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:37.081515Z",
"iopub.status.busy": "2022-12-11T12:47:37.081247Z",
"iopub.status.idle": "2022-12-11T12:47:49.399964Z",
"shell.execute_reply": "2022-12-11T12:47:49.399525Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Generating symbolic expressions for 3-Ricci tensor (Spherical coords)...\n",
"Finished generating symbolic expressions for 3-Ricci tensor (Spherical coords) in 0.3 seconds. Next up: C codegen...\n",
"Generating C code for 3-Ricci tensor (FD order=4) (Spherical coords)...\n",
"Finished generating C code for 3-Ricci tensor (FD order=4) (Spherical coords) in 7.4 seconds.\n",
"Generating symbolic expressions for BSSN constraints (Spherical coords)...\n",
"Finished generating symbolic expressions for BSSN constraints (Spherical coords) in 1.1 seconds. Next up: C codegen...\n",
"Generating C code for BSSN constraints (FD order=4) (Spherical coords)...\n",
"Finished generating C code for BSSN constraints (FD order=4) (Spherical coords) in 1.1 seconds.\n"
]
}
],
"source": [
"import BSSN.BSSN_Ccodegen_library as BCl\n",
"_ignore = BCl.add_Ricci_eval_to_Cfunction_dict(includes=[\"NRPy_basic_defines.h\"], rel_path_to_Cparams=os.path.join(\".\"),\n",
" enable_rfm_precompute=False, enable_golden_kernels=False, enable_SIMD=False,\n",
" enable_split_for_optimizations_doesnt_help=False, OMP_pragma_on=\"i2\")\n",
"\n",
"_ignore = BCl.add_BSSN_constraints_to_Cfunction_dict(includes=[\"NRPy_basic_defines.h\"],\n",
" rel_path_to_Cparams=os.path.join(\".\"), output_H_only=True,\n",
" enable_rfm_precompute=False, enable_SIMD=False,\n",
" enable_stress_energy_source_terms=True,\n",
" leave_Ricci_symbolic=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 5: Output `free_parameters.h`, which sets numerical grid parameters \\[Back to [top](#toc)\\]\n",
"$$\\label{cparams_rfm_and_domainsize}$$\n",
"\n",
"Here we output `free_parameters.h`, which sets the size of the numerical grid based on the star's radius, as prescribed in the variable `domain_size` (set above). Also, free parameters based on chosen reference metric (set above) overwrite defaults set in [`reference_metric.py`](reference_metric.py)."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:49.404156Z",
"iopub.status.busy": "2022-12-11T12:47:49.403880Z",
"iopub.status.idle": "2022-12-11T12:47:49.407673Z",
"shell.execute_reply": "2022-12-11T12:47:49.407279Z"
}
},
"outputs": [],
"source": [
"# Step 3.a.i: Set free_parameters.h\n",
"# Output to $Ccodesdir/free_parameters.h reference metric parameters based on generic\n",
"# domain_size,sinh_width,sinhv2_const_dr,SymTP_bScale,\n",
"# parameters set above.\n",
"\n",
"outstr = \"\"\n",
"outstr += rfm.out_default_free_parameters_for_rfm(\"returnstring\",\n",
" domain_size,sinh_width,sinhv2_const_dr,SymTP_bScale)\n",
"with open(os.path.join(Ccodesrootdir,\"free_parameters.h\"),\"w\") as file:\n",
" file.write(outstr.replace(\"params.\", \"griddata.params.\"))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 6: Register the `main()` C code for generating the `TOV_Playground` executable \\[Back to [top](#toc)\\]\n",
"$$\\label{mainc}$$\n",
"\n",
"As a core diagnostic, we will output grid data on the $yz$ plane. This diagnostic is constructed within the Python module [`diagnostics_generic.output_yz_or_xy_plane`](diagnostics_generic/output_yz_or_xy_plane.py)."
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:49.411420Z",
"iopub.status.busy": "2022-12-11T12:47:49.411246Z",
"iopub.status.idle": "2022-12-11T12:47:49.414923Z",
"shell.execute_reply": "2022-12-11T12:47:49.414533Z"
}
},
"outputs": [],
"source": [
"# T4UU00 = e^{-nu} rho, while alpha = e^{nu/2} -> T4UU00 = rho/alpha^2 -> rho = T4UU00*alpha^2\n",
"list_of_outputs = [\"y_n_gfs[IDX4ptS(CFGF,idx)]\",\n",
" \"log10(fabs(diagnostic_output_gfs[IDX4ptS(HGF,idx)]))\",\n",
" \"y_n_gfs[IDX4ptS(ALPHAGF,idx)]*y_n_gfs[IDX4ptS(ALPHAGF,idx)] * griddata->gridfuncs.auxevol_gfs[IDX4ptS(T4UU00GF,idx)]\",]\n",
"planar_diags.add_to_Cfunction_dict__plane_diagnostics(plane=\"yz\", include_ghosts=False,\n",
" list_of_outputs=list_of_outputs, num_sig_figs=4)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next we register the `main()` function for our `TOV_Playground` executable."
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:49.418588Z",
"iopub.status.busy": "2022-12-11T12:47:49.418493Z",
"iopub.status.idle": "2022-12-11T12:47:49.422823Z",
"shell.execute_reply": "2022-12-11T12:47:49.422443Z"
}
},
"outputs": [],
"source": [
"def add_to_Cfunction_dict_main__TOV_Playground():\n",
" includes = [\"NRPy_basic_defines.h\", \"NRPy_function_prototypes.h\", \"time.h\"]\n",
" desc = \"\"\"// main() function:\n",
"// Step 0: Read command-line input, set up grid structure, allocate memory for gridfunctions, set up coordinates\n",
"// Step 1: Set up initial data to an exact solution\n",
"// Step 2: Output data on xy plane to file.\n",
"// Step 3: Free all allocated memory\n",
"\"\"\"\n",
" c_type = \"int\"\n",
" name = \"main\"\n",
" params = \"int argc, const char *argv[]\"\n",
" body = r\"\"\" griddata_struct griddata;\n",
" set_Cparameters_to_default(&griddata.params);\n",
"\n",
" // Step 0.b: Read command-line input, error out if nonconformant\n",
" if((argc != 4) || atoi(argv[1]) < NGHOSTS || atoi(argv[2]) < NGHOSTS || atoi(argv[3]) < 2 /* FIXME; allow for axisymmetric sims */) {\n",
" fprintf(stderr,\"Error: Expected three command-line arguments: ./BrillLindquist_Playground Nx0 Nx1 Nx2,\\n\");\n",
" fprintf(stderr,\"where Nx[0,1,2] is the number of grid points in the 0, 1, and 2 directions.\\n\");\n",
" fprintf(stderr,\"Nx[] MUST BE larger than NGHOSTS (= %d)\\n\",NGHOSTS);\n",
" exit(1);\n",
" }\n",
" // Step 0.c: Check grid structure, first in space...\n",
" const int Nxx[3] = { atoi(argv[1]), atoi(argv[2]), atoi(argv[3]) };\n",
" if(Nxx[0]%2 != 0 || Nxx[1]%2 != 0 || Nxx[2]%2 != 0) {\n",
" fprintf(stderr,\"Error: Cannot guarantee a proper cell-centered grid if number of grid cells not set to even number.\\n\");\n",
" fprintf(stderr,\" For example, in case of angular directions, proper symmetry zones will not exist.\\n\");\n",
" exit(1);\n",
" }\n",
"\n",
"#include \"free_parameters.h\"\n",
"\n",
" // Step 0.d: Uniform coordinate grids are stored to *xx[3]\n",
" // Step 0.d.i: Set bcstruct\n",
" {\n",
" // Step 0.f: Call set_Nxx_dxx_invdx_params__and__xx(), which sets\n",
" // params Nxx,Nxx_plus_2NGHOSTS,dxx,invdx, and xx[] for the\n",
" // chosen (non-Eigen) CoordSystem.\n",
" int EigenCoord = 0;\n",
" set_Nxx_dxx_invdx_params__and__xx(EigenCoord, Nxx, &griddata.params, griddata.xx);\n",
" }\n",
"\n",
" // Step 0.j: Allocate memory for y_n_gfs gridfunctions\n",
" const int Nxx_plus_2NGHOSTS0 = griddata.params.Nxx_plus_2NGHOSTS0;\n",
" const int Nxx_plus_2NGHOSTS1 = griddata.params.Nxx_plus_2NGHOSTS1;\n",
" const int Nxx_plus_2NGHOSTS2 = griddata.params.Nxx_plus_2NGHOSTS2;\n",
" const int grid_size = Nxx_plus_2NGHOSTS0*Nxx_plus_2NGHOSTS1*Nxx_plus_2NGHOSTS2;\n",
" griddata.gridfuncs.y_n_gfs = (REAL *restrict)malloc(sizeof(REAL)*grid_size*NUM_EVOL_GFS);\n",
" griddata.gridfuncs.auxevol_gfs = (REAL *restrict)malloc(sizeof(REAL)*grid_size*NUM_AUXEVOL_GFS);\n",
"\n",
" // Step 0.l: Set up TOV initial data\n",
" ID_persist_struct ID_persist;\n",
" TOV_read_data_file_set_ID_persist(\"TOVdata.txt\", &ID_persist);\n",
"\n",
" initial_data_reader__convert_ADM_Spherical_to_BSSN(&griddata, &ID_persist, TOV_ID_function);\n",
"\n",
" REAL *restrict aux_gfs = (REAL *restrict)malloc(sizeof(REAL)*grid_size*NUM_AUX_GFS);\n",
"\n",
" // To simplify the expressions somewhat, we compute & store the Ricci tensor separately\n",
" // from the BSSN constraints.\n",
" Ricci_eval(&griddata.params, griddata.xx, griddata.gridfuncs.y_n_gfs, griddata.gridfuncs.auxevol_gfs);\n",
" BSSN_constraints(&griddata.params, griddata.xx, griddata.gridfuncs.y_n_gfs, griddata.gridfuncs.auxevol_gfs, aux_gfs);\n",
"\n",
" yz_plane_diagnostics(&griddata, griddata.gridfuncs.y_n_gfs, aux_gfs);\n",
"\n",
" // Step 4: Free all allocated memory\n",
" free(griddata.gridfuncs.y_n_gfs);\n",
" free(griddata.gridfuncs.auxevol_gfs);\n",
" free(aux_gfs);\n",
" for(int i=0;i<3;i++) free(griddata.xx[i]);\n",
"\n",
" return 0;\n",
"\"\"\"\n",
" add_to_Cfunction_dict(\n",
" includes=includes,\n",
" desc=desc,\n",
" c_type=c_type, name=name, params=params,\n",
" body=body,\n",
" rel_path_to_Cparams=os.path.join(\".\"), enableCparameters=False)"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:49.426501Z",
"iopub.status.busy": "2022-12-11T12:47:49.426408Z",
"iopub.status.idle": "2022-12-11T12:47:49.428929Z",
"shell.execute_reply": "2022-12-11T12:47:49.428494Z"
}
},
"outputs": [],
"source": [
"add_to_Cfunction_dict_main__TOV_Playground()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 7: Register the core NRPy+ infrastructure C codes that set up needed data structures, C macros, etc. \\[Back to [top](#toc)\\]\n",
"$$\\label{nrpy_infra_codes}$$"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:49.432593Z",
"iopub.status.busy": "2022-12-11T12:47:49.432501Z",
"iopub.status.idle": "2022-12-11T12:47:49.445481Z",
"shell.execute_reply": "2022-12-11T12:47:49.445015Z"
}
},
"outputs": [],
"source": [
"import outputC as outC\n",
"import finite_difference as fin\n",
"fin.register_C_functions_and_NRPy_basic_defines(NGHOSTS_account_for_onezone_upwind=True, enable_SIMD=False)\n",
"\n",
"outC.outC_NRPy_basic_defines_h_dict[\"MoL\"] = \"\"\"\n",
"typedef struct __MoL_gridfunctions_struct__ {\n",
" REAL *restrict y_n_gfs;\n",
" REAL *restrict auxevol_gfs;\n",
"} MoL_gridfunctions_struct;\n",
"\"\"\"\n",
"par.register_NRPy_basic_defines() # add `paramstruct params` to griddata struct.\n",
"list_of_extras_in_griddata_struct = [\"MoL_gridfunctions_struct gridfuncs;\"]\n",
"gri.register_C_functions_and_NRPy_basic_defines(list_of_extras_in_griddata_struct=list_of_extras_in_griddata_struct) # #define IDX3S(), etc.\n",
"\n",
"rfm.register_C_functions(use_unit_wavespeed_for_find_timestep=True)\n",
"rfm.register_NRPy_basic_defines()\n",
"\n",
"outC.outputC_register_C_functions_and_NRPy_basic_defines()\n",
"outC.NRPy_param_funcs_register_C_functions_and_NRPy_basic_defines(os.path.join(Ccodesrootdir))\n",
"\n",
"# Call this last: Set up NRPy_basic_defines.h and NRPy_function_prototypes.h.\n",
"outC.construct_NRPy_basic_defines_h(Ccodesrootdir, enable_SIMD=False)\n",
"outC.construct_NRPy_function_prototypes_h(Ccodesrootdir)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 8: Output all C codes, compile `TOV_Playground`, and run at two resolutions \\[Back to [top](#toc)\\]\n",
"$$\\label{output_compile_and_run}$$\n",
"\n",
"Finally, output all the C codes (and auto-generated `Makefile`) to `Ccodesrootdir/`, compile them, and run the executable at two numerical resolutions. We run at two resolutions to check that the Hamiltonian constraint violation converges to zero at a rate consistent with the order of our finite-difference operators."
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:49.449587Z",
"iopub.status.busy": "2022-12-11T12:47:49.449494Z",
"iopub.status.idle": "2022-12-11T12:47:51.921927Z",
"shell.execute_reply": "2022-12-11T12:47:51.920731Z"
},
"scrolled": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"(EXEC): Executing `make -j18`...\n",
"TOV_ID_function.c: In function ‘TOV_ID_function’:\n",
"TOV_ID_function.c:11:21: warning: variable ‘phi’ set but not used [-Wunused-but-set-variable]\n",
" 11 | REAL rbar, theta, phi;\n",
" | ^~~\n",
"(BENCH): Finished executing in 2.00 seconds.\n",
"Finished compilation.\n",
"(EXEC): Executing `taskset -c 1,3,5,7,9,11,13,15 ./TOV_Playground 48 12 2`...\n",
"(BENCH): Finished executing in 0.20 seconds.\n",
"(EXEC): Executing `taskset -c 1,3,5,7,9,11,13,15 ./TOV_Playground 96 12 2`...\n",
"(BENCH): Finished executing in 0.20 seconds.\n"
]
}
],
"source": [
"cmd.new_C_compile(Ccodesrootdir, \"TOV_Playground\", mkdir_Ccodesrootdir=False,\n",
" uses_free_parameters_h=True, compiler_opt_option=\"fast\") # fastdebug or debug also supported\n",
"\n",
"# Change to output directory\n",
"os.chdir(os.path.join(Ccodesrootdir))\n",
"# Clean up existing output files\n",
"cmd.delete_existing_files(\"out??.txt\")\n",
"\n",
"cmd.Execute(\"TOV_Playground\", \"48 12 2\", file_to_redirect_stdout=\"out48.txt\")\n",
"cmd.Execute(\"TOV_Playground\", \"96 12 2\", file_to_redirect_stdout=\"out96.txt\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 9: Visualize output \\[Back to [top](#toc)\\]\n",
"$$\\label{visualize_output}$$\n",
"\n",
"Next, we visualize the output on the $yz$ plane, analyzing the Hamiltonian constraint violation $\\mathcal{H}$ and the mass-energy density $\\rho$."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:51.925271Z",
"iopub.status.busy": "2022-12-11T12:47:51.924983Z",
"iopub.status.idle": "2022-12-11T12:47:52.526366Z",
"shell.execute_reply": "2022-12-11T12:47:52.525469Z"
}
},
"outputs": [],
"source": [
"import numpy as np\n",
"from scipy.interpolate import griddata\n",
"from pylab import savefig\n",
"import matplotlib.pyplot as plt\n",
"from IPython.display import Image"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:52.530849Z",
"iopub.status.busy": "2022-12-11T12:47:52.530669Z",
"iopub.status.idle": "2022-12-11T12:47:52.646455Z",
"shell.execute_reply": "2022-12-11T12:47:52.645547Z"
}
},
"outputs": [],
"source": [
"xy_extent=domain_size\n",
"# Data are in format x,y,z, CF,Ham,rho\n",
"output_grid_x, output_grid_y, output_grid_data_Ham = \\\n",
" plot2D.generate_uniform_2D_grid('out96.txt', 1,2,4, [-xy_extent,xy_extent], [-xy_extent,xy_extent])\n",
"# Ham data are in column 4 --------------------------^\n",
"output_grid_x, output_grid_y, output_grid_data_rho = \\\n",
" plot2D.generate_uniform_2D_grid('out96.txt', 1,2,5, [-xy_extent,xy_extent], [-xy_extent,xy_extent])\n",
"# rho data are in column 5 --------------------------^\n",
"\n",
"output_grid_data = []\n",
"output_grid_data += [output_grid_data_Ham]\n",
"output_grid_data += [output_grid_data_rho]"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"execution": {
"iopub.execute_input": "2022-12-11T12:47:52.650874Z",
"iopub.status.busy": "2022-12-11T12:47:52.650705Z",
"iopub.status.idle": "2022-12-11T12:47:53.873835Z",
"shell.execute_reply": "2022-12-11T12:47:53.873103Z"
}
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAABFAAAAH8CAYAAAD2YVsOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAADJY0lEQVR4nOzdeXhTVf4/8PdN2qYUaNkKZakURBFkXwd0BJ064OCC4yA4KosODiqOiqMjfpXFZcABFXFDnQHUEcUVHRdGBNFxRFEYFDdcRgdEWkClBQptk9zfH/xIzvnc5twmTdokvF/P08fc3JObk5ub5Hi5n/exbNu2QUREREREREREEXkaugNERERERERERMmOJ1CIiIiIiIiIiFzwBAoRERERERERkQueQCEiIiIiIiIicsETKERERERERERELngChYiIiIiIiIjIBU+gEBERERERERG5yGjoDhAREaWagwcPoqqqKmHbz8rKQnZ2dsK2T0RERKTi2KZ2eAKFiIgoCgcPHkSnjk1QsjOQsOcoKCjAN998kxYDDSIiIkpuHNvUHk+gEBERRaGqqgolOwP434Yi5DaNfyVs+d4gOvb/FlVVVSk/yCAiIqLkx7FN7fEEChERUQyaNLXQpKkV9+0GEf9tEhEREbnh2MYdQ2SJiIiIiIiIiFzwChQiIqIYBOwgAnZitktERERU3zi2cccrUIiIiIiIiIiIXPAKFCIiohgEYSOI+P8zTSK2SUREROSGYxt3vAKFiIiIiIiIiMgFr0AhIiKKQRBBJKKiNzFbJSIiIjLj2MYdr0AhIiIiIiIiInLBK1CIiIhiELBtBOz41/QmYptEREREbji2cccTKERERDFg0BoRERGlE45t3LGEh4iIiIiIiIjIBa9AISIiikEQNgL8VxoiIiJKExzbuOMVKERERERERERELngChYiIKAaH64QT8UdERERU35JtbHPfffehqKgI2dnZGDx4MNavXx+x7SeffIJzzjkHRUVFsCwLCxYsMG577ty5sCwLV111VVR94gkUIiIiIiIiIkoay5cvx7Rp0zBz5kxs3LgRvXv3xogRI7Bz584a21dUVKBz586YO3cuCgoKjNt+//338eCDD6JXr15R94snUIiIiGJweKq/RPwRERER1bdkGtvceeedmDx5MiZNmoTu3btj0aJFyMnJweLFi2tsP3DgQMybNw/jxo2Dz+eLuN19+/bh/PPPx8MPP4zmzZtH3S+eQCEiIiIiIiKihCovL9f+Kisra2xXVVWFDRs2oLi4OHSfx+NBcXEx1q1bV6c+XH755Rg1apS27WjwBAoREVEMggn8IyIiIqpviR7bFBYWIi8vL/Q3Z86cGvuxe/duBAIBtGnTRru/TZs2KCkpifn1Pfnkk9i4cWPE560NTmNMRERERERERAm1bds25ObmhpZNpTaJeO4rr7wSq1atQnZ2dszb4QkUIiKiGARgI5CAGXMSsU0iIiIiN4ke2+Tm5monUCJp1aoVvF4vSktLtftLS0tdA2Ij2bBhA3bu3Il+/fqF+xUI4K233sK9996LyspKeL1e1+2whIeIiIiIiIiIkkJWVhb69++P1atXh+4LBoNYvXo1hgwZEtM2f/GLX2Dz5s3YtGlT6G/AgAE4//zzsWnTplqdPAF4BQoREVFMAvahv0Rsl4iIiKi+JdPYZtq0aZgwYQIGDBiAQYMGYcGCBdi/fz8mTZoEABg/fjzat28fyjOpqqrCp59+Grq9fft2bNq0CU2aNEGXLl3QtGlT9OjRQ3uOxo0bo2XLlo77TXgChYiIiIiIiIiSxtixY7Fr1y7MmDEDJSUl6NOnD1auXBkKlt26dSs8nnBBzffff4++ffuGlufPn4/58+dj2LBhWLt2bdz6Zdl2DJMyExERHaHKy8uRl5eHTZ+2RtOm8a+E3bs3iD7dd6KsrKxWdcJEREREdcGxTe3xChQiIqIYBGEhACsh2yUiIiKqbxzbuGOILBERERERERGRC16BQkREFIOgfegvEdslIiIiqm8c27jjFShERERERERERC54BQoREVEMAgmqE07ENomIiIjccGzjjlegEBERERERERG54BUoREREMeC/0hAREVE64djGHa9AoaRzySWXoG3btsjNzUXPnj3xj3/8o6G7RERERFRnHOMQEaU2nkChpDNt2jR8++23KC8vx+LFi3HBBRfghx9+aOhuERFpgraVsD8iSk8c4xBRMuPYxh1PoFDSOe644+Dz+QAAlmWhqqoK27dvb+BeNaylS5fCsix8++23Kbn9+rRv3z54PB7cddddDd2VGhUVFWHWrFlJuz0iIkocjnEoWaXTWJAokXgChepsw4YNGDlyJHJzc9G0aVP88pe/xKZNm2psu3HjRpx55plo0aIFcnJy0KNHDyxcuNDR7rLLLkOjRo0wcOBAnHLKKejZs6exD4e/9LOzs2sciAwfPhw9evSI6fVR2OH9HOnv3Xffbegu4uOPP4Zt267HDFFdHa4TTsQfESWHZBrjWJaFt99+27Hetm0UFhbCsiycfvrpMb3OREqFsUMiydefnZ2Ndu3aYcSIEVi4cCH27t3b0F2s0TvvvINZs2Zhz549Dd0Vqkcc27hjiCzVycaNG3HiiSeisLAQM2fORDAYxP33349hw4Zh/fr16Nq1a6jta6+9hjPOOAN9+/bFTTfdhCZNmuDrr7/Gd99959ju/fffj3vuuQdr167Fxx9/DMuq3YeusrISc+fOxT333BO315gMLrzwQowbNy70r1YN7eabb0anTp0c93fp0qUBeqPbvHkzAPAEChER1UmyjXGys7OxbNkynHjiidr9b775Jr777rukGSNEksxjh/pw+PVXV1ejpKQEa9euxVVXXYU777wTL774Inr16tWg/ZNjzXfeeQezZ8/GxIkT0axZswbtG1Ey4QkUqpObbroJjRo1wrp169CyZUsAwAUXXIBjjz0WN9xwA5599lkAQHl5OcaPH49Ro0bhmWeegcfjfvGT1+vFL37xCyxYsADHHHMMfvWrX7k+pk+fPnj44Ycxffp0tGvXrm4vLgns378fjRs3htfrhdfrbejuhJx22mkYMGBA1I87/Hpqe38sNm/ejPz8fLRp0yYu26sPlZWVWLRoEVasWIGPPvoIwWAQa9asQd++fbV2Tz75JGbMmIExY8bgtttua6De0mEBeBBIwIWcgbhvkYhikWxjnF/96ld4+umnsXDhQmRkhIfwy5YtQ//+/bF79+4YX2n9iHXsEG/xHHNEQ77+6dOnY82aNTj99NNx5pln4rPPPkOjRo3qvV+HJdtYkxoGxzbuWMJDdfKvf/0LxcXFoYEFALRt2xbDhg3DSy+9hH379gE49ONeWlqK2267DR6PB/v370cwGKzVc/j9fnz11Ve1anvDDTcgEAhg7ty5xnYTJ05EUVGR4/5Zs2Y5/iXo8H1ffPEFLrjgAuTl5SE/Px833XQTbNvGtm3bcNZZZyE3NxcFBQW44447HNvdvn07LrroIrRp0wY+nw/HH388Fi9eXOPzfPrpp/jtb3+L5s2bh/6VKVJd6vbt23HxxRejXbt28Pl86NSpEy699FJUVVUBAP73v//hsssuQ9euXdGoUSO0bNkSY8aMqZf61kivx/Q6I/n888+xdevWWj3v5s2bHVefvPLKKxg+fDhyc3PRrFkznH/++REHmo888gh69eqFRo0aoX///li3bh3OPfdc9O7du3YvPEo//vgjZs2ahRNOOAFr1qzBBRdcgD179uB///uf1s7v9+MPf/gDAoEAT54QEdWDZBvjnHfeefjhhx+watWq0H1VVVV45pln8Nvf/tbRPpoxwN69e3HVVVehqKgIPp8PrVu3xqmnnoqNGzfWan08HB4ffPXVV6GrHvLy8jBp0iRUVFQ42td1bAUAa9euxYABA5CdnY2jjz4aDz74oGMs+MYbb8CyLDz//POOPixbtgyWZWHdunUxveZTTjkFN910E/73v//h73//e1SvLZp9Vpv3Tx1rzpo1C9deey0AoFOnTqHyo2+//Tau++ONN95A//790bhxYzRu3Bhjx45FWVlZrfcfUUPgCRSqk8rKyhrPlufk5KCqqgoff/wxAOD1119Hbm4utm/fjq5du6JJkybIzc3FpZdeioMHD4YeV1ZWhmXLlmHfvn3w+/14+umn8cYbb+Ckk06qVX86deqE8ePH4+GHH8b3338fnxf5/40dOxbBYBBz587F4MGDceutt2LBggU49dRT0b59e9x+++3o0qUL/vjHP+Ktt94KPa60tBQ/+9nP8Prrr2Pq1Km4++670aVLF1x88cVYsGCB43nGjBmDiooK/PnPf8bkyZMj9uf777/HoEGD8OSTT2Ls2LFYuHAhLrzwQrz55puhH833338f77zzDsaNG4eFCxdiypQpWL16NYYPH17jYKS2ysrKsHv3bu0v0iwCkV5PbV8nAHTr1g3jx4+vVd/kCZT58+dj1KhRaN68OebNm4crrrgCL7zwAs4991zHY+fMmYOJEyeiS5cuWLhwIX72s5/hzDPPxMaNGxNWEtSiRQvMmTMHAwYMCJ2o83g8jpNKa9aswa5du3DmmWcmpB8UPTtBKfV2GiXVE6WyZBvjFBUVYciQIXjiiSdC97366qsoKyvDuHHjHO2jGQNMmTIFDzzwAM455xzcf//9+OMf/4hGjRrhs88+q9X62qjt2OHcc8/F3r17MWfOHJx77rlYunQpZs+erbWJx9jqP//5D0aOHIkffvgBs2fPxsUXX4ybb74ZK1as0B47fPhwFBYW4vHHH3ds9/HHH8fRRx+NIUOG1Ho/SBdeeCGAQ2Vgsbw2wH2fRfv+/frXv8Z5550HALjrrrvw2GOP4bHHHkN+fn7c9sdf//pXFBcXo1u3bpg3bx5GjRqFp556CldddZXrYylxOLapBZuoDnr27Gkfe+yxtt/vD91XWVlpH3XUUTYA+5lnnrFt27Z79epl5+Tk2Dk5OfYVV1xhP/vss/YVV1xhA7DHjRsXemxZWZk9fPhwOy8vz87NzbX79etnP/vss679WLJkiQ3Afv/99+2vv/7azsjIsP/whz+E1g8bNsw+/vjjQ8sTJkywO3bs6NjOzJkzbfmxOHzfJZdcErrP7/fbHTp0sC3LsufOnRu6/6effrIbNWpkT5gwIXTfxRdfbLdt29bevXu3tt1x48bZeXl5dkVFhfY85513XsTX980334TuGz9+vO3xeOz333/f0T4YDNq2bYe2rVq3bp0NwH700UeN26/J4XY1/fl8Pq1tpNdjep2RALCHDRvm2u7777+3AdgPP/ywbdu2/dZbb9mWZdkzZszQ2t133302AHv9+vWh+zZs2GB7vV57+vTpWtuLL77YBmDPmTOn1v016dixoz1z5swa11VWVtpNmjSp8bVOnjzZBmCvWrWq1tujxCgrK7MB2K9t7mj/+9tOcf97bXNHG4BdVlbW0C+V6IiWjGOce++9127atGno933MmDH2ySefbNv2od+DUaNGhR5X2zGAbdt2Xl6effnll0fsg9v62vTfbexweHxw0UUXaY8/++yz7ZYtW2r3xWNsdcYZZ9g5OTn29u3bQ/d9+eWXdkZGhmMsOH36dNvn89l79uwJ3bdz5047IyPD9TdYff8iycvLs/v27RvVa1Nfn9s+q837J8eC8+bNizg2rMv+sG3b/uyzz+yMjAx74cKF2v0nnXSSnZ2dbVdXV7tug+KLY5va4xUoVCeXXXYZvvjiC1x88cX49NNP8fHHH2P8+PHYsWMHAODAgQMADk0tW1FRgfHjx2PhwoX49a9/jYULF+L3v/89nnzySXz55ZcAgNzcXLzxxhvYs2cPysrKsGHDBvz617+Oqk+dO3fGhRdeiIceeijUj3j43e9+F7rt9XoxYMAA2LaNiy++OHR/s2bN0LVrV/z3v/8FcCgZ/9lnn8UZZ5wB27a1f3UZMWIEysrKHJe/TpkyxbUvwWAQK1aswBlnnFFjPfHhS0/Vfzmrrq7GDz/8gC5duqBZs2Z1uuz2vvvuw6pVq7S/V199tca2kV5PbV7nYbZtY+3ata7tPvroIwDhANnbbrsN7du3x4wZM7R2J5xwAoBDpUGH/fnPf0bbtm0d/8r185//XNtmIr311lvYt28fdu3aheHDh4f+zjjjDDzxxBNo0qRJrf+lkoiI6iYZxzjnnnsuDhw4gJdeegl79+7FSy+9VGP5DhDdGKBZs2Z47733Il6967a+Nmo7dpDjg5///Of44YcfUF5eDiA+Y6tAIIDXX38do0eP1jLzunTpgtNOO83Rp/Hjx6OyshLPPPNM6L7ly5fD7/fjggsuiH5nCE2aNMHevXtjem01vT65z+Lx/qnquj9mzZqFXr16YerUqdr9J510Eg4ePIgff/wxLv0kSgSeQKE6mTJlCm644QYsW7YMxx9/PHr27Imvv/4a1113HYBDPwhA+Ef88OWAhx3+0Y+1djSSG2+8EX6/3zULJRpHHXWUtpyXl4fs7Gy0atXKcf9PP/0EANi1axf27NmDhx56CPn5+drfpEmTAAA7d+7UHl9TQr20a9culJeXu07NfODAAcyYMQOFhYXw+Xxo1aoV8vPzQ4O3WA0aNAjFxcXa38knn1xj20ivpzavM1qbN2+GZVk4/vjjUVFRgdWrV+O3v/2tIxRt//79AICsrCwAhy7TfvXVVzF27FhkZmZqbQ/XuB8+gfLAAw+gX79+yMzMxKxZsxx92LVrF0aNGoXGjRuja9euWL16da37//LLLwMAXnjhBaxduzb0d9ttt2Hfvn0oLi4O9ZkaXsD2JOyPiBpeMo5x8vPzUVxcjGXLluG5555DIBDAb37zmxrbRjMG+Mtf/oKPP/4YhYWFGDRoEGbNmhX6x6DarK+qqkJJSYn2FwjosZG1HTvI8Vbz5s0BIK5jq507d+LAgQM1zgBU033HHXccBg4cqJWtPP744/jZz34Wl1mE9u3bh6ZNm8b02gD3feb2/kWrLvvD7/fjlVdewW9+8xtH7uD+/fthWRZyc3Nj7hvVDcc27jgLD9XZbbfdhj/+8Y/45JNPkJeXh549e+KGG24AABx77LEAgHbt2uGTTz5xzIzSunVrAOEv+Hjp3LkzLrjgAjz00EO4/vrrHesjTRkof+xVNSWTR0ort20bAEIhchdccAEmTJhQY1s5bV08E9ivuOIKLFmyBFdddRWGDBmCvLw8WJaFcePG1Trgrq4ivZ5EJM1v3rwZnTp1QpMmTfDhhx/C7/dr00wedjiwr3v37gCAr7/+GhUVFTUGxf73v/9Fbm5uaHDStm1bzJo1C8uWLauxD5dffjkKCgqwa9cuvP766zj33HPx5ZdfokWLFq79f+WVV3D00Uc7Bh8rV64EAIwaNcp1G0REFD/JOMb57W9/i8mTJ6OkpASnnXZaxClmoxkDnHvuufj5z3+O559/Hq+99hrmzZuH22+/Hc899xxOO+001/XvvPOO42TIN998U2Ngv5tkHVuNHz8eV155Jb777jtUVlbi3Xffxb333lvn7X733XcoKytDly5dYnptgPs+c3v/YhHr/ti4cSP27t2LPn36ONZt2rQJvXv3RnZ2dkx9IqoPPIFCcSFTzV9//XV06NABxx13HACgf//+WLVqVShg7bDDlxLm5+fHvU833ngj/v73v+P222+vsb979uxx3C9nPqmr/Px8NG3aFIFAAMXFxXHdbm5ubijALpJnnnkGEyZM0GYGOnjwYI2vPR2oAbKHT5L5fD5Hu7/97W/o2LFj6Aqew5dhy6knKysr8cQTT2hX+owePRrAoZMd0r59+7BixQr897//RU5ODs4880z07NkTL7zwQuhfjiL5+uuv8cUXX+Cyyy5zrDt8AqU201xS/QnCQjABF3IGYcd9m4edeeaZ2LRpE3bu3InmzZujuLgYt99+e1pM+06UKMk2xjn77LPx+9//Hu+++y6WL18esV20Y4C2bdvisssuw2WXXYadO3eiX79+uO2220L/g21a37t3b212IAAoKCio+4utQTzGVq1bt0Z2dnaNMyBFmhVp3LhxmDZtGp544gkcOHAAmZmZGDt2bEzPr3rssccAACNGjEjYuBFwf3+lSP/YeFis+2PTpk0A4JhKeseOHXj77bdx44031u4FUUKk4timvqXPtTSUNJYvX473338fV111Veh/SA/PePK3v/1Na/vXv/4VGRkZGD58eNz7cfTRR+OCCy7Agw8+iJKSEse6srKyUGYGcOiLu6Yp2erC6/XinHPOwbPPPlvjyY5du3bFtF2Px4PRo0fjH//4Bz744APH+sP/4uD1ekO3D7vnnnuMV9o0pIqKCnz++eeOKYZrM41xIBDAZ599FjqBcvTRRyMjIwP//ve/tXZPPfUU3nrrLfzpT38KDQ46duwIAHj77be1tnPnzsWOHTtcS6UO+/LLL9GkSRN06NAhdF/Pnj3xySefuD728AmZESNGaPfv378f//73v9G7d2/+Ty7V2cknn4ynnnoKW7ZswbPPPouvv/464uX/ROSUDGOcJk2a4IEHHsCsWbNwxhlnRGxX2zFAIBBwlPS0bt0a7dq1Q2Vlpet6AKETsupfoq4iiMfYyuv1ori4GCtWrNByQb766quImW6tWrXCaaedhr///e94/PHHMXLkSEcZd7TWrFmDW265BZ06dcL555+fkHFjbd6/mhw+wRHphFus++Pw2PvNN98M3ef3+3HppZciLy8Pv//97123QdSQeAUK1clbb72Fm2++Gb/85S/RsmVLvPvuu1iyZAlGjhyJK6+8MtSub9++uOiii7B48WL4/X4MGzYMa9euxdNPP43p06cn7H8M/+///g+PPfYYtmzZguOPPz50/7hx4/CnP/0JZ599Nv7whz+goqICDzzwAI499tg6havWZO7cuXjjjTcwePBgTJ48Gd27d8ePP/6IjRs34vXXX485KOvPf/4zXnvtNQwbNgyXXHIJunXrhh07duDpp5/G22+/jWbNmuH000/HY489hry8PHTv3h3r1q3D66+/jpYtW9bpNb366qtaAOthQ4cORefOnWPe7vr163HyySdj5syZWr5It27dQsdMJF9++SUOHjwYOoHSuHFj/P73v8eiRYuQkZGBXr164b333sPixYtx4YUX4tJLLw09tlWrVhg5ciQeeuihUNvXX38dGzZsAFD7ANl9+/Y56nZzc3MjTvF8WCAQwJIlS9CmTRuceuqp2rovvvgCVVVVoRKi3bt313nARvERgIUA4j8tXyK2edjVV18dut2xY0dcf/31GD16NKqrqx35P0RHumQe40Qq71DVdgywd+9edOjQAb/5zW/Qu3dvNGnSBK+//jref/993HHHHa7rayueY4d4jK1mzZqF1157DSeccAIuvfRSBAIB3HvvvejRo0foKglp/PjxoZPOt9xyS1R9Pvz6/X4/SktLsWbNGqxatQodO3bEiy++GDrhFO9xY6zvX//+/QEcGkuPGzcOmZmZOOOMM7QrR2LZH5s3b0b37t1x2223Yf/+/Wjbti2efPJJvP/++3jmmWccpXBUv1JxbFPfeAKF6qR9+/bwer2YN28e9u7di06dOuHWW2/FtGnTkJGhH16LFi3CUUcdhSVLluD5559Hx44dcddddyV0vvcuXbrgggsuwCOPPKLd37JlSzz//POYNm0arrvuOnTq1Alz5szBl19+GfcTKG3atMH69etx880347nnnsP999+Pli1b4vjjj6+xvKi22rdvj/feew833XQTHn/8cZSXl6N9+/Y47bTTkJOTAwC4++674fV68fjjj+PgwYM44YQT8PrrrzuucoiWnNXmsCVLltTpBEpdbN68GQC0q0XmzZsHy7LwxBNP4K9//SuOPfZY3HfffTX+68YjjzyC3//+93j00UfRuHFjnHvuuZg/fz7OPvtsDBw4sFZ9aNKkSSjx/rDy8vJQ0GAk999/P0aNGoXx48c76rS7d++OUaNGhWYL+N3vfofTTz+9Vv2h1CaPJZ/PV2NJWqx+/PFHPP744xg6dChPnhDVINnHOG5qOwbIycnBZZddhtdeew3PPfccgsEgunTpgvvvvx+XXnopqqqqjOtrK55jh3iMrfr3749XX30Vf/zjH3HTTTehsLAQN998Mz777LMaT/QAwBlnnIHmzZsjGAzizDPPjKrPh19/VlYWWrRogZ49e2LBggWYNGkSmjZtGtfXpnJ7fyMZOHAgbrnlFixatAgrV65EMBjEN998o51AiWV/bN68GZdffjmOOuoo3HLLLaFyotdeew2nnHJK1K+PqL5Ztry2j4iIAABTp07F8uXLsWPHDsdgecqUKSgoKNCulNm3bx9atGiBb775Bu3btwdwqGRi/PjxoQyUoqIiTJw4scYZfGIR7+2Ru/LycuTl5eH5D49B46Y1B/fVxf69AZzd+0vH/fLKrFj96U9/wr333ouKigr87Gc/w0svvVTnq9KIiNLF6NGj8cknn4Smn1b5/X60a9cOZ5xxhqNk60gU7f7Ytm0bjjrqKCxbtswxaxU1rPoa25SVlaX8LEvMQCGiI15lZaWjTvytt97Cgw8+iClTpmgnT/x+Pw4ePIhAIKDdBg5dgXLWWWdh5syZOHDgAF566SV89NFHOOuss+r19VB62LZtG8rKykJ/06dPr7Hd9ddfD8uyjH/qv6Zee+21+M9//oPXXnsNXq8X48ePdxz/RERHgsMh8od9+eWXeOWVVyLm1qxYsQK7du3C+PHj66F3yS/a/XH4amG1rJ4o1bCEh4iOeP/85z9x88034ze/+Q2aNWuGjRs3YsmSJejbt29ousrDbr31VsyePTu0fNttt2HJkiWYOHEigEPlOBMmTEDLli3RoUMHLF++vFZTGFPqOZRUH/+a3sPbzM3NrdW/0lxzzTWh4y8S9fL4Vq1aoVWrVjj22GPRrVs3FBYW4t1338WQIUPq1G8iolTTuXNnTJw4EZ07d8b//vc/PPDAA8jKysJ1112ntXvvvffw0Ucf4ZZbbkHfvn0xbNiwBupxcoh1f3z00Ufwer3abFWUXBI9tkkHPIFCREe8Fi1aICsrC3PnzkVFRQWKiorwpz/9CdOnT3dkksyaNctYRpGfn1/jFMdEiZKfnx/zNKnBYBAAjDMxEBGlq5EjR+KJJ55ASUkJfD4fhgwZgj//+c845phjtHYPPPAA/v73v6NPnz5YunRpw3Q2icS6PzZv3oyjjz46rnleRPWNGShERERROFwn/PSHxyEnAXXCFXsDGNP787jXCb/33nt4//33ceKJJ6J58+b4+uuvcdNNN6G0tBSffPIJB7RERERHqFQd2zQEZqAQEREdAXJycvDcc8/hF7/4Bbp27YqLL74YvXr1wptvvsmTJ0RERES1wBIeIiKiGARsDwJ2/P8dIpCgC0N79uyJNWvWJGTbRERElPpSbWzTEHgCpRaCwSC+//57NG3aFJaVPgE4RETpyrZt7N27F+3atYPHk5iLLYPwIJiACzmDSJ9BBiUvjm2IiFILxzbJgSdQauH7779HYWFhQ3eDiIiitG3bNnTo0KGhu0GUdDi2ISJKTRzbNCyeQKmFpk2bAjh0sKZ66A0R0ZGgvLwchYWFoe/vRAjYFgJ2/P/lPhHbJJI4tiEiSi0c2yQHnkCphcOXtubm5nKQQUSUQliaQFQzjm2IiFITxzYNiydQiJLIqZ4xidlwNF+09RTytCr4dL08D1GiBOBBIAF1woE0qhMmIkrY2CYJcWxDqY5jG3ecxpiIiIiIiIiIyAWvQCEiIopB0PYgmICp/oJpNNUfERERpQ6ObdzxBApRgjkuXVXLaSz9C8rKiNNH0krQxWV2sA4P1b84T80YV+ttrwo8FfPzEhERUXwdSWU50Yhmv7Dchyg18QQKERFRDFgnTEREROmEYxt3zEAhIiIiIiIiInLBK1CIiIhiEAQQsOM/lWDshXJEREREsePYxl3KnUB56623MG/ePGzYsAE7duzA888/j9GjR0dsv3btWpx88smO+3fs2IGCgoIE9pTSyanec80NlMwRy6N/6ThyTUz5JJ7af2E1zBzw5ovWbENAlOUVdwRlW6WByEP5ZdZ5+vM4Hqut1BaZn0JEyY5jG2oIzDFpWHXZ/8xPIWo4KXcCZf/+/ejduzcuuugi/PrXv67147Zs2YLc3NzQcuvWrRPRPSIiOkIE4UEwAZWwidgmJTeObYiIKBlwbOMu5U6gnHbaaTjttNOiflzr1q3RrFmz+HeIiIiIqA44tiEiIkoN6XMqyEWfPn3Qtm1bnHrqqfj3v/9tbFtZWYny8nLtj4iISBWwPQn7I6oNjm2IiCieOLZxl3JXoESrbdu2WLRoEQYMGIDKykr89a9/xfDhw/Hee++hX79+NT5mzpw5mD17dj33lOqDVm8qM0QM2SSWVwR4yLaG7BJjVonH5cukQXJOYmfsrcxHcWSiqNkl+n6R2SpanorMQxHBV6dmjDP0SURaKc/D+mJyE4SFoPmoj3m7RCYc25CKWSZHnljfc45tyA3HNu7S/gRK165d0bVr19Dy0KFD8fXXX+Ouu+7CY489VuNjpk+fjmnTpoWWy8vLUVhYmPC+EhEREbnh2IaIiKhhpP0JlJoMGjQIb7/9dsT1Pp8PPp+vHntERESpJlGXpKbTZa5Ufzi2ISKiuuLYxt0ReQJl06ZNaNu2bUN3gxLAcUmjLIHxhGs/5HTDprIc1ymD1VKcaNoK5nKfOF76ZppKWZJlLtEwTDfsmPJYLZMK6s/peOXqYy3ZP/Ha1LYu5T7q9MiOqatFf3kZLBElE45t0hdLdChe3I4ljm2I3KXcCZR9+/bhq6++Ci1/88032LRpE1q0aIGjjjoK06dPx/bt2/Hoo48CABYsWIBOnTrh+OOPx8GDB/HXv/4Va9aswWuvvdZQL4GIiNJAAB4EEpDFnohtUnLj2IaIiJIBxzbuUu4EygcffICTTz45tHy4nnfChAlYunQpduzYga1bt4bWV1VV4ZprrsH27duRk5ODXr164fXXX9e2QURERNRQOLYhIiJKDSl3AmX48OHOS/8VS5cu1Zavu+46XHfddQnuFRERHWmCtoWgnYCk+gRsk5IbxzZERJQMOLZxl3InUIgkLavCo8+N68g5kdMRq20deSmx5Zo4t2N4rFsWSTS5Jwmb8jiavBTDVMUif8SZa6JkmXjFc4rHav+j4fK6LUNeii2+zLXpkQMB0T29D+pxtyrwlLEPRERE0WDuCTUE9bhjHgpRzXgChYiIKAbBBNUJB9OoTpiIiIhSB8c27tLnlRARERERERERJQhPoBAREcUgaHsS9kdERERU35JtbHPfffehqKgI2dnZGDx4MNavXx+x7SeffIJzzjkHRUVFsCwLCxYscLSZM2cOBg4ciKZNm6J169YYPXo0tmzZElWfWMJDSU/LOAEcuSGWmmsiMk7ilmsiMzlM2SVuuSWm563Ldo1PGZ98FFPIYY3U3BC37001A0U+j2XIT7H1XBPHZgPKerEfnDkshmwV2SclI+XUjHGird4nZqQQEZGKGSeU7NyOUWakUKItX74c06ZNw6JFizB48GAsWLAAI0aMwJYtW9C6dWtH+4qKCnTu3BljxozB1VdfXeM233zzTVx++eUYOHAg/H4/brjhBvzyl7/Ep59+isaNG9eqXzyBQkREFIMALAScp+Hisl0iIiKi+pZMY5s777wTkydPxqRJkwAAixYtwssvv4zFixfj+uuvd7QfOHAgBg4cCAA1rgeAlStXastLly5F69atsWHDBpx00km16hevEyYiIiIiIiKihCovL9f+Kisra2xXVVWFDRs2oLi4OHSfx+NBcXEx1q1bF7f+lJWVAQBatGhR68fwChRKCo7LBJXpiC059bCjhEdZ9shSm8glPK7TDavP41gXRRmOofTG0YdoSm3ka60Hrr0LinIaUxdtOTWxYetiCmG9REY8idiupb4fstwniumRLUdZkbI+ILYrXotW4hPUp0fmJbCpK1F5JcxAIUoPLNOhdGY6vjm2SV2JHtsUFhZq98+cOROzZs1ytN+9ezcCgQDatGmj3d+mTRt8/vnn8elTMIirrroKJ5xwAnr06FHrx/EEChERUQwCSEy5TcC9CREREVHcJXpss23bNuTm5obu9/l8cX+u2rr88svx8ccf4+23347qcTyBQkREREREREQJlZubq51AiaRVq1bwer0oLS3V7i8tLUVBQUGd+zF16lS89NJLeOutt9ChQ4eoHsvrhImIiGKQbFP9EREREdVFsoxtsrKy0L9/f6xevTrct2AQq1evxpAhQ2J+fbZtY+rUqXj++eexZs0adOrUKept8AoUahByamIrQxyKlmEKYVPOiVjnyBiReSrautinKtaexzR1suSWeWJaX4dpjRMmaHitLlMgOzJGtO3KjJHw8zimVpZ5KZYhL0Vkl2h7NKAXUjieR3lvLPkey/6qzyPeN/lZ4JTHRESpiZknRIfIzwIzUSgW06ZNw4QJEzBgwAAMGjQICxYswP79+0Oz8owfPx7t27fHnDlzABwKnv30009Dt7dv345NmzahSZMm6NKlC4BDZTvLli3DCy+8gKZNm6KkpAQAkJeXh0aNGtWqXzyBQkREFIOA7UEgAVeLJGKbRERERG6SaWwzduxY7Nq1CzNmzEBJSQn69OmDlStXhoJlt27dCo/yj5jff/89+vbtG1qeP38+5s+fj2HDhmHt2rUAgAceeAAAMHz4cO25lixZgokTJ9aqXzyBQkRERERERERJZerUqZg6dWqN6w6fFDmsqKjIecW44La+NngChYiIKAY2LAQTkFRvJ2CbRERERG44tnHHEyiUMI5aYE84f8SSWSRi2TLkmpjbig+nzDFRc07kurrkmpjaSp4o2prWuz22IXgMZ3Xdzviq62WOidzfSltHdko0eSnyGLDFY9Wmcp2y/035KICerWLbkdcBwC+zzgu39fu1dawhJiJqWMw5IYqe6XPDsQ2lGp5AISIiikEy1QkTERER1RXHNu7S55UQERERERERESUIr0ChOtEuyZNlC4mamtg03bDbVMQew9SzpumGHW2jKK1xW1bYjhKkKMp0GqKkJ5ogJtHWkmU6KlkKJJ/HOOVx5CmFHXtIlPvYQbW8yqW8R5maWJb3mPaKo+RIvG+2ul3xGZJTHqv7gZfA1r+gbSFox/9zl4htElHtsUyHqP6wvCe5cGzjjlegEBERERERERG54BUoREREMQjAg0AC/h0iEdskIiIicsOxjbv0eSVERERERERERAnCK1AoKsapiWV+h2m6YTmNsSnnxGW6YbUPjqmITTknbtkkpumGDcuuOSbRTE1sylpxe2x9iCYDRWSe2NFMeWzKT3FkisgpkJX9InNXZG6PenyY8lEONVbaBsRmI+enqBknNfZBfZjLlMe28nrkZ5N1w4nHOmGi9MDME6LkxLFN/ePYxh1PoBAREcUgCA+CCbiQMxHbJCIiInLDsY279HklREREREREREQJwitQiIiIYhCwLQQScElqIrZJRERE5IZjG3c8gUJGpswTQOSeyMwTmXOiZD041pnyUzzmtsZsFVP+iCNLxROxbVS5JtHkmLjko9imXJNkuH4s6N7kMEtmnrjknOjPY8hPiSYvRfZB5Jxoz2PKRwGAQDj3xNFzS2xX+dFwvKPyPVdXidcin8eC2gf92GfdMBFRzZh5QpSaOLahZMATKERERDFg0BoRERGlE45t3CXDv2ETERERERERESU1XoFCDtrlcaaSHQBWhnIIGUpgAMDKVNrKaV698rHKcjRTEzvWGcp0ZOmEqUzHbbuGEh5HGY7ptGUUUx4by3vcRHPq1FSmIyumDGU4dlQlO2K7pvIfuR1Z7qOW08jyHsd7o6w3lfcIcjuOKY8NpUEQ0xqrnzG72q+vM5Wo+fW2ppIeXvIaH7btQdCO/79D2AnYJtGRjmU7ROmHY5v449jGXfq8EiIiIiIiIiKiBOEVKERERDEIwELAGQscl+0SERER1TeObdzxChQiIiIiIiIiIhe8AoWMUxXLzAVLZpWo2SBubdVcE9M6AFaGOo2xOVtF25Yj00J/rJZz4rZdQ1tjrkk0UxO7ncI05ZwY81FcthsN5a2wDLElQA05J/pKfdmQrSIzT4z5KW55KWoeiWk6ZEDLPXHko8ipidXjIxDFlMcBfZXMkVFmJtZzgwDYgcg7TcsjAhzZKuor5TSA8RG0E5Mqb4jbIaJaYuYJ0ZGFY5v44NjGHa9AISIiIiIiIiJywStQiIiIYhBMUFJ9IrZJRERE5IZjG3fp80qIiIioViorK9GnTx9YloVNmzY1dHeIiIiIUgKvQDkCmTJPAD33xJh5Ami5J462IhNFex6PS2aEuux4zsjZJbZhnXO7LnkppqySaNoaMlFcSww9aluXxvUQbu0oX3SpZ7S0rBLR2CvbKps1ZZ4AWu6JIy8lKENRvBHbwtFW3d+GzBMAlpoxIt8bU1aJzMQJyuNDOYCCemCKfItNu990OMjHsW44NkFYCCbgg5eIbUrXXXcd2rVrhw8//DDhz0VUH5h5QkQqjm1ik8pjm/rCK1CIiIiSUHl5ufZXWVkZl+2++uqreO211zB//vy4bI+IiIjoSMETKERERDEI2FbC/gCgsLAQeXl5ob85c+bUuc+lpaWYPHkyHnvsMeTk5NR5e0RERJQ+Ej22SQcs4TlCaJexGUp2Di1HnppYlq6Y2xq265iKWDw2mqmJjduNXBpknIrYra03irIcUfqhbUu2jWZqYreSHhPTqVPD9MKuHKU3Sh/ldL22nKo4fNtylOzoy2qVi6PcxzGlcHi9qbxHtoUsrZHTAiv73zKU7AAwT3ks2yq1TG4zvqmPjaatZCrp4SWvkSU6aG3btm3Izc0N3e/z+eq0Xdu2MXHiREyZMgUDBgzAt99+W6ftETU0lu0QUW1xbFM7DJF1lz6vhIiIKI3k5uZqf5FOoFx//fWwLMv49/nnn+Oee+7B3r17MX369Hp+JURERETpgVegEBERxSAIC8EEXJIabdDaNddcg4kTJxrbdO7cGWvWrMG6descJ2IGDBiA888/H4888ki0XSUiIqI0kixjm2TGEyhEREQpLD8/H/n5+a7tFi5ciFtvvTW0/P3332PEiBFYvnw5Bg8enMguEhEREaUFnkBJU6apio2ZJ4BxCmHjVMWmzBPZNpqsEreplNVtiba2Y3rkyNsx5pzIqWe9seeaaCd1XXJYjFMXG1a5TnlsIrNK1KeUeSPOFno/jLkmoo/qennmW05VbNiuZcm2yrIhHwUQGSkyH8UwVbFjf4vX5pqRonUi3Ke6TFsc6xTH8rGcBjAyO0FT/dkJ+leao446Sltu0qQJAODoo49Ghw4dEvKcRPHEzBMiigeObSJLtbFNQ2AGChERERERERGRi5Q7gfLWW2/hjDPOQLt27WBZFlasWOH6mLVr16Jfv37w+Xzo0qULli5dmvB+EhFRegvaVsL+6kNRURFs20afPn3q5fkoMo5tiIgoGaT62KY+pNwJlP3796N379647777atX+m2++wahRo3DyySdj06ZNuOqqq/C73/0O//znPxPcUyIiIiJ3HNsQERGlhpTLQDnttNNw2mmn1br9okWL0KlTJ9xxxx0AgG7duuHtt9/GXXfdhREjRtT4mMrKSlRWVoaWy8vL69bpemDKPAEAyxM+62fMPAG0rBJjjol4Hve2ynqXPthKf539M+WliLObMn9E64NLW2W9I+PCU4dcE49hu8ZcE3GHp/Zncuty0leNFHGtXwyKPBLTY0VIh5ZlEow9L8W2ZB+UTBFTPgqgZaRYlp5bIvN0tEyUoDnjRH0WR5aKIVtFqrdMFFvJdxGhOKwbDgvaHgTt+P87RCK2ScmNY5uaMfOEiOoDxzZhHNu4S59XEsG6detQXFys3TdixAisW7cu4mPmzJmDvLy80F9hYWGiu0lERERUKxzbEBERNYy0P4FSUlKCNm3aaPe1adMG5eXlOHDgQI2PmT59OsrKykJ/27Ztq4+uEhFRCmGdMDUUjm2IiCgROLZxl3IlPPXB5/PB5/M1dDeiI8oALFneoU03bC610UpxDCU77m0NpTeyZMdU0lOXqYllmY5a7iOrQhx9sCK2leU/aimOLfe9fB512VQKBGjvq+v3ToK+l4wTF8uVYr9oFTNyumHxWFudQVjuQ1nuo5b4OMpwZB/VfSj6IMpl1PIaW1TSWAFR/oNopjxW+yCeE7FTHxvNBNNubbUyPrFD7TT6ASQ6kqTk2IaIiCjJpP0JlIKCApSWlmr3lZaWIjc3F40aNWqgXhERUaoLwkIwAWcuE7FNSi8c2xARUSJwbOMu7Ut4hgwZgtWrV2v3rVq1CkOGDGmgHhERERHFjmMbIiKihpFyJ1D27duHTZs2YdOmTQAOTeW3adMmbN26FcChGt/x48eH2k+ZMgX//e9/cd111+Hzzz/H/fffj6eeegpXX311Q3SfiIjSBOuEKV44tiEiomTAsY27lCvh+eCDD3DyySeHlqdNmwYAmDBhApYuXYodO3aEBhwA0KlTJ7z88su4+uqrcffdd6NDhw7461//GnGav1SiTrllZYi30pFr4jWsM2SimDJPHG2jmW7YbSplpa1jXe2nJpZ5JNrzypgNR25M5OmGHdvVcliiyDWRpzAdeSmGLxtH/yM3TRRLzrgrc03UtnKdDOJQdpRlmA4ZEHkpYopjOVWx9kQyH8XyRGxqOdpGMeWxODdt+smQuyHWn5d4TnFsykuRj1W/h460af8SNSBIp0EG1Q7HNmGcupiIGhrHNhzbmKTcCZThw4fDdvyfV9jSpUtrfMx//vOfBPaKiIiIKDYc2xAREaWGlDuBQkRElAz4rzRERESUTji2cZdyGShERERERERERPWNV6CkMjWfRGQ5WDI7w7IM62QeSXi7xswTsV1T5gkg8kdc8lLMuSbisRlqtorMHxGvVcknceajRM45ceSLODJQ1LwU2VZf1PrklmOirHe8FkdeCuqdjAWRgRmWckm6LTpozE+R+0y0tdT9L/JSHG2D6rEvc0xE1oqyLZljAo9oq2zL9iPiOrne7W3ScmNc2mqhLTLbyLBduc8c/TV8X9iG74sjDf+VhoiIiNIJxzbueAUKEREREREREZELXoFCREQUAxtAMOZ5k8zbJSIiIqpvHNu44wmUFHKq91xtWZ2a2K3UxjSNsbP0Ri0xiVwKBADIiDyNsXGqYlPJDmAsy3E+j6EsxzQ1sUtb83b1zWqlNi5ttfVuZTiG7QZl5YS6WVnG4ihVQUSO/mr7QWzX7XnUF+Qo75FPbNiu7K9aaiP2i5zWWCtXcUxNLLer3JZtDfvMEhfz2cFgxPWOcp/Im3Wf4lj9BveLDhvK1xyfecPzuPZBqZmS31GrAk8Zn4eIiNMWE1Eyk99RR9q0xuTEEyhEREQxYJ0wERERpROObdwxA4WIiIiIiIiIyAWvQCEiIooB/5WGiIiI0gnHNu54AiWJyZo7K0O8Xd7I+SPOaYLVDBQZZCGmQNayStzyUtSsErc+KG0NmSeHtqV8yAyZJ/J56zI1samtY1pdx/MY1hlyTYIZ5ql9MyrCd1Q11dtWtNE7Vd67MnS7TZsybV1TX6W2nOkJ52VkefTsDI/oxO4DTUK3t33fQluX+5FPW260U39s1r7wsj/bfMGbJ2DINZH5KcrTWI5cE9OUvKKpIzak9lMe688TOfMEAGxlvWNdvDJRvOZ4Li3XRExjDJmno84F7figyMcqWUxiM6wbJiKJmSdElMo4tiGeQCEiIooB/5WGiIiI0gnHNu6YgUJERERERERE5IInUIiIiGJw+F9pEvFHREREVN+SbWxz3333oaioCNnZ2Rg8eDDWr18fse0nn3yCc845B0VFRbAsCwsWLKjzNmvCEp5k5hHBBDKrRA10EOEOliOPRG0rtpMhnkfLQDHnmug5IZEzTw61VZZNmSdiW8EMw2sRfbAzDDkm4nlMOSaAnp3hWCfzUqzI64LiU+b7KRy8Ud5RX1lxyj5t+cXBi0K3O2Vka+u2+g9oywElgSIQx/8B8ypZH4GuYrsn64vtxP5/80DL0O2r/n2etq75uixtuenW8H6pbK7vcE+1nruhZZXISB+Za6Iuy1PGMhRFzQZxZKuIZS27RBzPIj/Fo7S13fJS/Oo6M1vJI7FkVomzcbitOCZlDosp18QOiv6rnxP5ZsjvsDRj2xbsBJzsSMQ2iYhSnvzNjheZ7UV0BEumsc3y5csxbdo0LFq0CIMHD8aCBQswYsQIbNmyBa1bt3a0r6ioQOfOnTFmzBhcffXVcdlmTXgFChEREREREREljTvvvBOTJ0/GpEmT0L17dyxatAg5OTlYvHhxje0HDhyIefPmYdy4cfD5fDW2iXabNeEJFCIiohgEYSXsj4iIiKi+JXpsU15erv1VVlbW2I+qqips2LABxcXFofs8Hg+Ki4uxbt26mF5bvLbJEp4kc6r33NBtS0wh7CjLMU1j7KhrUEttIpfLOJYdUwgb2ppKdgCtbMdUsgOIsh1ZamMo05FlOXKaYNPUxKa2jplcxXLAp5SUiNKPkmK9PuLek/4eul2U8ZO2zuOYqDbstQONteUVP5ykLb/x1bHhPmzTy30yKiJPC+xW+eFvFG4c7HBQWze407fa8rmt39eWj8ncFbr98vB7tHXB4XqfvvU3D92e+q/ztXVtV2ZqywGl+sdbFbm8B9A/CnKaaLm/TVMee2SZi/rtKbZribbq8Sy34yjpUS5xtMVB6vjfauWyY+MUxwDgMVyi7BEvQL3M0lGWE/mAcfZP3676/bYq8FTk/hBRWuHUxaRRyjstOSaVY19t/OoyYDGV+JjKdERpqi3bBvSBna2W+wZlvS8dadTvN05p7K6wsFBbnjlzJmbNmuVot3v3bgQCAbRp00a7v02bNvj8889jeu54bZMnUIiIiGLAqf6IiIgonSR6bLNt2zbk5uaG7o9UapPMeAKFiIiIiIiIiBIqNzdXO4ESSatWreD1elFaWqrdX1paioKCgpieO17bZAYKERFRDA4n1Sfij4iIiKi+JcvYJisrC/3798fq1atD9wWDQaxevRpDhgyJ6bXFa5u8AiXZaIENLtkDyrJx2mJAn5bUkZci23oirzO0dWSeOKb6jVzHapqq2HVqYmXZkWMi+6Bu121qYqVLgSx9XcYBvT62dHD49kfjFurrAlXa8l47/LF790Anbd1fPvqlttxkTTj3JPfban07R+m5IE0aKfkdhkiLQw3UlWKVo0xY2b87GmlrPn/7OG35jp3HasvlheGdXHXiXm3dH3us0pYHNvo2dPv1U+7W1rUo1o+PgU9eE7rdYa1ef1yZp7+xWkaKOH7FbtLza+QMx6KtIxNFJb5Z1UwUeaw7MlHUj594VkcmilcNs5FvpOFzLdrK7wRtrWwrs44MbZ1ZTJwqkogoLam5JlFk+Mm2MgNFG886xqR1+HdgNaNL/s4FXTJQ1GXHumDktsxLIYrKtGnTMGHCBAwYMACDBg3CggULsH//fkyaNAkAMH78eLRv3x5z5swBcCgk9tNPPw3d3r59OzZt2oQmTZqgS5cutdpmbfAEChERUQyYgUJERETpJJnGNmPHjsWuXbswY8YMlJSUoE+fPli5cmUoBHbr1q3wKP8I+P3336Nv376h5fnz52P+/PkYNmwY1q5dW6tt1gZPoBARERERERFRUpk6dSqmTp1a47rDJ0UOKyoqcs6iFeU2a4MnUIiIiGKQqLwSZqAQERFRQ+DYxh1PoDSwU73nasuWVhdqziqx1GVHzag3cltHVknk57Ed2SqGZbHKkYmiZD/IvBGZ2aJllRjWAXruiSnzBACU+BEEDZknALT9cLCZvvLuW+/Xlo/L3B+6/UW1/rH6y/dnaMsfP9ctdDt7t36WNKeJeG1KzElZZz3zRL42NYzC9TvKcHLW9FiPKN+V+3t/a5HvobT3vdVUW/fgS2dry3e3DG+r+zn6XOzXt3tVW147dl7o9sdnt9TW3TD/Im05R9nHQfExscRnyqNmdPjNO1GtcjbmoQDaZ0NuVe5DS2lsi+1atnwD1M+f/qbKt9hSz8jLz7Ws/1aWbfFdIqm14rb8jnI8T3ivye++VYGnjM9DRKnjVM+Yhu4CxZscV2bqYx0rI7zsyDURbfW8lMj5KIfaqr9zLrl80VB/92TmSTBoXvZHzjWxqvUfbjUDxfaLdaItM1JSn/zuWxV8uoF6QonEEyhEREQxsBNUJ5xO/0pDREREqYNjG3ecxpiIiIiIiIiIyAWvQKlnjsta5eXx6qWLpnIZ2VZexmiYxthtu7ZpGmM5halafuA6jbFy26UsR5262FSyI/vgKOERu1ct25HbydyvX6K544IDodv/HKKX7FSIDT9S1it0+68r9KmI877U+2A1C9+ubiwvR9UXtX2WYW4b6XFx5ZLL5CgrUsgZbOVrV8tgvnmwq7Zu7DH68u/Ofi10+1dNPtbWPTl9vrZ8+rpLQ7fbPKVPwyyPAa0sR77YQOTX5jrFsfpNK9fJN0u5lNiSpW6izsxWntmSNWiiLEf9XDteicdU0yXWyf56tA+2+bHq95C4ZDoVL3u14XyJ8douUSphyU6akGWYGeHSYStLLyOWJTxQS3gyIq87tBx5TGosHW+gEh5LTE1sqev9ouwmU/zIK2U7lizhyRDLVdXh2/5q6Cv5y5CKOLbRt5sueAKFiIgoBkFYsExnMuuwXSIiIqL6xrGNO5bwEBERERERERG54BUoREREMeBUf0RERJROOLZxxxMo9U3Wl4p6Tm1aVdO0xYBeC2rpFxNZMo8k0uNqeB5tWWaeyO1qUx675Zp4Iq8zPNaUeSLbmjJPACCYGV6WmSdtp3+lLS9q/1Lo9o/BLG3dee/9TlvOfzacreE7Sn9Of45e9adO7Sv7a8pAierKtwR9R0VdE2lHuA3UMJ9v+GZ1Y32V7ye98QuzikO3H/7NCdq6pwY/rC2/NOSB0O1bO56mrftyYXe9S4YMF0cmihZQo68Kyhkb1TJn+Z7LHxTlcwK/SFeRUx4rj5WfTcsr+6ssu2UmKduV3yW2rBVXp14X6xzTGmvTqYvPcRr9sBIRpQSRw+eROSdZytjHkYGiL2s5JyIfxc4QP3zKsjo2PNQJQ9aeabwKc/6bzGHTfhPF4MYKiN8ymYnij5yBYsllZRpjVIlckwx9WZ3+2a7S90NQPpZTHhM1GJ5AISIiikHQtrSTWPHcLhEREVF949jGHTNQiIiIiIiIiIhc8AoUIiKiGNh2gqb6S6e5/oiIiChlcGzjjidQ6oE2B7ioN5XZJfp892KdVz42cp6A47GG7co8EmMfxKJtbGvKVnHJS/Fakdc5slXCt02ZJwDgqQ5/envP/I+27sr8N7Tlfx8oCt2+8/5ztXWN/doiKvLDt72V+jqZ4aJmfziuZnNbVhivhEvUVXIuX36OGmP1tTrqjyM/Vu4zjyj93d86/Flo/FYTbd3Ed67Wlq+67JnQ7dntXtHWPXzdHm35zT8PCfch07wT9UwUc1st60a+bnHMausd2STys6tsxxb10MHIn2tb1nuL7ULNMnF8rkUui3ogyu8o+aZrOSwyT0lvq35vrgo+DSJKLtrYhlKGpWSVWD6fvi5Lz3tTc08c62TOibJsO9bJDJTw938ww5B5Ahiz9py5cZF/iy3H71HkdZZfBpuJLDstA0Xko1SLDJTq8IDREr+R6nsBALaaKSbaesSyXRkebNp+MSilpMWxTXrgCRQiIqIYMKmeiIiI0gnHNu6YgUJERERERERE5IJXoNQzxzSphjIXx7TFzo1FbiuX1emIHVOYRi69cVwO6SgrUtvqq0ylN851hsfKdYapih3lMsLImW+Gbp+f94G27omy/tryijtPCT9nnr4duVu0/pumIpbro2krNMiJXHmYGWb2dTCU7ACi2sOtrXIMWIaKEgC4f+45odvfTdPLtKa0fEdb9t0Yvgx29cwTtXV+X+TzzR5HSYzok9opjx15HaAf73KduIxYK/8xlPcc6pMdeaVjamJlvXxtxmmNxaXMjikmI5cdysR3W7yvyYj/SkNEyc7K1EtvrOxw2Y5lmrYYAHzhZVu0lSU8wSy1hEf8TogSnqCyXpblyHGcqaTbMYYy/LOwHCtoJTzytzVgXvYoJT6ealnC4xHL4dfuqYo8nTMAWN4oyuLV9Qf12nG7ugpEseLYxh2vQCEiIiIiIiIicsErUIiIiGIQtC3HlTPx2i4RERFRfePYxh2vQCEiIiIiIiIicsErUOqDWvcvpu405py45aWoy6YpjsWys4bUNI2xvsqRiaJNL2eeHtk4NbFjWuPwbVM9rGzrrdLrVPvfuEFbVnNPXtp3vLZOzTwBgIMtw8/jPaitqiGzRVlwyTHRTsC6nIxtkKmKo+Don2GaY0dTw2ON+Sguz+MRs/mp7+OLd56srWtyrf7Gjm+2PnT7u+nNtXVf3KofL/7syFMrylJrLSNFHr/ytWhTHhuyjABYSlCI/Gw6p0dWP1Qih8WUR+KamaS8Wse0y2JZmd3RmY8iQ1uS4AB3Yds1vH9x2i4RUSxMmScAYCm5JnIaY4icEztbzUARU+5myQwUr3Jb/z4PiGVbGdcFM8Xvp/yN1MaOenejyY0zZavJfBRH5olcrg4vW369U94qfdlTpWxcTNlsibbqb7z8jZT5Y9H8RjIThaLBsY07nkAhIiKKwaFBRiKC1uK+SSIiIiJXHNu4YwkPEREREREREZGLlDyBct9996GoqAjZ2dkYPHgw1q9fH7Ht0qVLYVmW9pednV2PvSUionR0eKq/RPzRkYdjGyIiamgc27hLuRKe5cuXY9q0aVi0aBEGDx6MBQsWYMSIEdiyZQtat25d42Nyc3OxZcuW0LIxdyQBLDWfxDSPu1x2yQSIKj/FMmUayCwCZVn2z1Cb6nhOR1ZJ5MwIRyaK2taxTn8atY619Fw90+KyVm9qyxsr24VuL134K21ddUv9edTcE2PmCaDX3brEQKjrXb9L4nWo1uWyuWj6YGjr1gX1obKtsazZ5TSw+j4eFO/xo/ecpi0XTvsxdHtam1XautEXHK0tt3qmUeQ+2IZjVua7iPgiLQvGkWMiHxte70hMl31QPteW/KwGDd8XjrayqfJZdWQkGb7f5LWcYtmSuU5EaSwVxzaUPKyM8HDelHlyaFlZL9apmScAEPSFM1GC2TIDReSc+JQMFJFrEvDpbdVsu4Aeu+LMFPOq6/S2ruMttakhS80K6Ks8ATHeFhko3mqlrV9fF8wUmSiZ4R9Nr8j0s0Qmikf5jfS4jPmj+ryreWl+v6EhEdVGyl2Bcuedd2Ly5MmYNGkSunfvjkWLFiEnJweLFy+O+BjLslBQUBD6a9OmTT32mIiI0pGdwL9EKSoqcly5MHfu3AQ+I9UGxzZERJQMUnFsU99S6gRKVVUVNmzYgOLi4tB9Ho8HxcXFWLduXcTH7du3Dx07dkRhYSHOOussfPLJJ8bnqaysRHl5ufZHRESUDm6++Wbs2LEj9HfFFVc0dJeOaBzbEBERpY6UKuHZvXs3AoGA419Z2rRpg88//7zGx3Tt2hWLFy9Gr169UFZWhvnz52Po0KH45JNP0KFDhxofM2fOHMyePTvmfp7qGaMta5dWGqYBdqx3lOEYSnpctquVwZjKhgDttJpp2mK53llBELmMwTGVsiyRUS/ZFOsC4tLQqtzw8jNDHtTW7Qnql6PefP8F4YVG2irjVMXGkh25PopLSqMqj3E5deu4PDWKx5o3bNhsNCVIblM2m57TdOmtKCkxlfTI9ziQrT+Renw8fOXd2rp//OwBbXnc2j+Gbvv2iMt7RTmKejw7p2iWn4VwA9tRWiPnYVTKcgyfTUB8fcid5BH9Vfpvmcp75LKjZEeW5YTX20H9jXP0X7ktv1NXBZ9GMkhUTW+i64SbNm2KgoKChD4H1V6qjm2oAXn0Wha1LMcSUxE7pipWynbsRvo6WaYT9IWXA2KdoyzHp5TlOKYx1rugjuOCooQnKMpcbEMJj9vYTG8smio/QbKEx1HSI6peAkoJj1rOAzjLdGxlt8kpmr1eWROrPqnog/xNV2+7lMRC+b21g3KdeLFUbzi2SV0pdQVKLIYMGYLx48ejT58+GDZsGJ577jnk5+fjwQcfjPiY6dOno6ysLPS3bdu2euwxERERHFcLVFZWxmW7c+fORcuWLdG3b1/MmzcPftbEpxyObYiIiBpGSl2B0qpVK3i9XpSWlmr3l5aW1vpf0zIzM9G3b1989dVXEdv4fD745Bl7IiIiVaKKev//NgsLC7W7Z86ciVmzZtVp03/4wx/Qr18/tGjRAu+88w6mT5+OHTt24M4776zTdil2HNsQEVHSSPDYJh2k1BUoWVlZ6N+/P1avXh26LxgMYvXq1RgyZEitthEIBLB582a0bds2Ud0kIiKqs23btmlXDEyfPr3Gdtdff70jGFb+HS4FmTZtGoYPH45evXphypQpuOOOO3DPPffE7eoWih7HNkRERKkjpa5AAQ4N/iZMmIABAwZg0KBBWLBgAfbv349JkyYBAMaPH4/27dtjzpw5AA6F5f3sZz9Dly5dsGfPHsybNw//+9//8Lvf/S5xnXSbytNEzQZxmW5Yyy7xGvJR5LJp2mLHdg3TFsttuUxNDG1qYhjbqs8j22bt1WtGL535HCK54F39fW5cpSyYckzksktb49TEdcg50SIvDPW7clmuixe36ZzrtF+imAJZm/JYHvpRZKJ4RO2y+kTy2JH5On++ZrFye4K2LpBlOJ7leyOn7zV8TkzTGsvPpnNaY0P+iGlaY/n9IL9r1P575M43vKkyeykgi9DVKd2TtG42QXXCh/dbbm4ucnNzXZtfc801mDhxorFN586da7x/8ODB8Pv9+Pbbb9G1a9eou0rxkRJjG2o44jvQI3NOsrJqvA3AOFWxI/MkW99uINur3Na/s/1iOaBkoPh9en+DMgNF+Y2UGSi2+D8T0zTGicpA8cifI7/Y/8rYIVilt5XTMmcoUxV7vaITjv4rL97ld8+jZpXJcYT4jYeSe+IRv7XBSvm7nUaXCVBsEjy2SQcpdwJl7Nix2LVrF2bMmIGSkhL06dMHK1euDIWvbd26VZtH/aeffsLkyZNRUlKC5s2bo3///njnnXfQvXv3hnoJREREcZOfn4/8/PyYHrtp0yZ4PB60bt06zr2iaHBsQ0RElBpS7gQKAEydOhVTp06tcd3atWu15bvuugt33XVXPfSKiIiOJLadmH+sS9Q/AK5btw7vvfceTj75ZDRt2hTr1q3D1VdfjQsuuADNmzdPzJNSrXFsQ0REDS3VxjYNISVPoBAREVF0fD4fnnzyScyaNQuVlZXo1KkTrr76akybNq2hu0ZERESUEngCJRFMmSey7l8uGx5rmXJNnI0jLzsyT/RF29BW5qWo5WyOdaacExnPIOtalfXBTL3x9wP0DZ/Q6NvQ7WfK+2rr2i7TZxwoP8rUP31Zy+9wqbPVyvrqknkic02UUlVHTW5Af3CjH8LLOTv0QEjvQTFNqXoaWJ4SFu9joFG4oLeiQC9kPtBK3zEyh0OtXUZdapcjR384HmbKRDHloQD6/pbHzj97Hq8tX5D7Yej29yfqL67VJlGPfDC8LI91Z2aLmlekb8cOys+1sl68b/L7Qv18Or5mHNs1fF8Ylh3PKZ5Gf2Lxwh2ZKGqBenLmndsJqhNOSO0xgH79+uHdd99NyLaJKHGsjMiZJwAANRNFZp749McGleWgT//fADXz5NBy+LvX30hmoFgRlwMiAyUgJn9SM1FkBkowU//lCCpddIzbxO+pHhwniO9V9fde/g57xJDJU23IQBFjVK/4Pyu1j85xTxS/bXKspvTZIzNQHG3DjS2RgWIF9BdvV4tQFzripNrYpiHwBAoREVEsbCsxoWhpNMggIiKiFMKxjavk/Gc9IiIiIiIiIqIkwitQEsAxJahcVtuapiqWl/ZFUbLjmEI4mkvylad1n3pWnXLV/Lq1ch+5zrGsPIUoXXl53HxteY9yLeiS507V1vm6iO6ql2VGM62uqWSnpsdqjUVTtepCXjYqptX1VIUbt/j0gLauOk+/7rWsKLxc3rGR/pxRzC4rr4BV34sMvQvI/Va/zjVrj37p54/dw/0Iiql9HVMXRnXcKW1dmkZTGqQe+3u66F+Pf3t2hLZ8yvjPQrffH3eHtm7UBj1TQpvG2LB/AWjTDcrvEkd5kva5FuU+8mOuVc/UvixHfjaNU2S7TeGulCSZSowOtVVKg5L0Xy0YtEZECeMJ135YYtpiGKYxtsU6O0tOVRxeDmTLEh5RpqOU7VSLEp5Att4lrYRHrJMlPAFf+EvOUcKTJX7LMtUSWLdpgGtfwgNtGmPxOydLdkRVizpW82aI3zJHubKhfN1lqmK9sawNUjYrpy0WZTlqmY7tKOHRl7X1QTEApyMCxzbueAUKEREREREREZELXoFCREQUCxs1JOXGabtERERE9Y1jG1e8AoWIiIiIiIiIyAWvQImDUz1jtGUrQ9+tlnEKYUP9oyNTxJCJ4ppromQaRDEtqaMPpimQ3XJC1IwLx9TJ+rI6LVzJL8R8csK870eGbjf7Uq/7rMyV0+wa+if7H00mh/ZA0VQuq1MTi5eW97UeguKtCr+e0sE5EfsHAAFlWkBHqWwd6oTVzIuAmDnxx27iWA/oy/kfhqdTDmTpnSg7WkytqD5U1hDLPpryaWRbNU/HLQtGPfT1maAdx5Z63N3W4UVt3Y5f6m9swevhFyczfRy5JupxJ6e1NnymHLkm8rMbjPx9Ib8TtOyVKPJSnNO0ixegLsr+yWm61RwWfZXjO3dV8Gk0BE71R+lMfs6oflmZGTXePrQsgkPU9aKtIwPFF/6BDfjk1MSRpyp2ZJ40EtMYK9FrjgyUbP1bXM1AsX0y80T8bijLngyR7SEyUTzK2EZ+j8rfEVvJPQn6xVix2rwcrFR+nxzTFhuyyxw/4nL6YS0sUG8qxwPB8Pvo8Yvjwy8a+wPKOnHs+PXxiqVsy65kBkp94tgmdfAKFCIiIiIiIiIiF7wChYiIKFZpVNNLRERExLGNGa9AISIiIiIiIiJywStQ4sGRCWA4LyUzAty2ZeKJnGkQzXZMmShynSxfs7Vslcjr5LbkOplj4fspXHd577C/a+sO2npAxqfPHhe67Wmhb8dxBtWUnRGn0jxT5olcbrX5gLbuh2564bC/cbhTfj0CBcFM/YnUDBFn5onLssqWNbnqk4g++PRlT5W+4V19wg0y9uvbla99V69GiEhmosT4XsnHOT4mavmxeM6DLfSdqh53+//wirbu3pP0Y3be0xeGbvsbi9pqx+ckcv9scXBpj3UceCLXxPC5dnxlGfJSjDlO4vhwMH0vye/GgPJBkR2M5nsygVgnTERx49F/dLQ8PZGthyw9x8JWck+CMvMkyyuWw9+nQZ/+XRMQy2oGiinzBAACyrJfZp7k6D8Oti+8bGXrg6QMkYGSmRXO6MjIEG1FzpbXE/mfzQMiz8sfDO8Hv1/fR9VV+j4MyAyUjHB7O0P8pkfze+UYOCsPE7+nnoDeNhhQ3sdqMWARmSgeNQOlWhxLMrdRWbarRVBfkJkoRwKObdzxChQiIiIiIiIiIhe8AoWIiCgWNhJTJ8zaYyIiImoIHNu44gmUeDCV7ADmsh3HpX7hZSua6YbdaFON1v5hztIPw/SnLmU5pvKZYIZ+R3nH8KF5TOYP2rp/VRytLWf/GP5EVueI/plKP9ymMTa1lQyXXcoSnlabD4Zu7+ynXwMrpwlWL4k1lewAgJ2hTN8np4mO0zTGjmmA5XR9YkpBW3lf5XssX3v+f8L7ZXdPMQdiNO+NYVpjx9to2K6jIkZcyZpdHm7w/oEibd3PGn2jLavHc/ZP+oa9QflESn8cXw9if6vvjeP7QbwX6mtD7GTJkWUo95F90i5tdtYN6Yvq92ZAfqh48SQRpRfL6424bMkSHtFWLeGxM0W5SZa+HFCWA2Kd31DS45yaWF9Wy3YCjUXJTrYotckJ/6Bm+aq1ddlZ+nKOsuzL0H+IfV592WMY2wTF2KYyEN5nlaLkpaJKL5E6KJarMsLLAY8o93H8wCr72DC+OrQcXu8s/9bfK3VM4hElRpachtk0JbYs6VGPO3mcsYSHCABPoBAREcXIQt1OR5m2S0RERFTfOLZxwxMoREREseBlrkRERJROOLZxxeugiYiIiIiIiIhc8AqUGJ3qGRNeEFPPObJAFI6cApNopiJ25I8Ylh3TkIptafkS5rYRHwc48yUil4E6zkpWnLIvdNsrVs7/+FRt2ddUqRk1TVsslus0m5aMrVDqWD16+S7yvtbv+KlreGrfoF5W65yq2Bd+omCGzBcRXVLyR2QWSV2mMVYXLVHca8mMGcOyY7psQd0vzb7S91nZ0WLKRkNOSDRXCBqnNZbbEc9TpRx38ph8YdAibVk9nhs909jYB/X1ODJnTDFDpqmIJbe2pu+LaETxWEdeirogv99EGbb6fbwq+HStn7PO+K80lGa0sQ3VK8srvpgzI09jbGfIDBRvjbcBZ85JMEvJNRG5a0GZw+ZDxLYBw1TFjsyTxvpveqNGVaHbTbIrtXW5Wfpy06xwPlpORpW2zufRfwwyPZEzOqqD+n6pVJYr/PqL25upB7yUZ/q05X3e8Os7IAaeAejjlYA6hhIDFseYSum+R84gLMaW2vsoMk88VfL4UN6bDMN02YB23FlVIktF9IESi2Ob5MUrUIiIiIiIiIiIXPAKFCIioljYVh0vZTNsl4iIiKi+cWzjilegEBEREREREVFSue+++1BUVITs7GwMHjwY69evN7Z/+umncdxxxyE7Oxs9e/bEK6+8oq3ft28fpk6dig4dOqBRo0bo3r07Fi1aFGFrNeMVKHFgGTJPDjUwZRHUIRPFmGsi8gRMbV0eG3E7gB7CIDMMjG31VZkVer3sy4MfCN2uEjWjWf9uqm9Wqalz5EmYsj+iaStXyYgRtW61Sl9pBfVlNTvDlHkCAMFMu8bbAGCLTBTtdKjHLQPFUIjo2InKTZmXI7JWHJ8Fwz4MiOepUmqB5T6T+zToVY9n0SdT7k3kVc62LseDugvlMRkYqDd+cXD4i3nyo1frbbMj7zP5GXJ8lSidcLSVr1Z9b0yZJ2Jbjvc0XnkpLrkmpm3JPtnBCA0TzLYdkUFx2y4RpTmZn+eVeXrKssitcC6Hf5yDmSK3IkP/vgwqy4FMfV0gS7RVokECYnwil21f+IvYm6MHeKiZJwCQ2yica9I8+4C2roVvv7bcLDO8volXz0fJ8erbzbQMGSgipK1CCXXZF9AzTvaIrJVMr55d5jGMoQ6IwUMgGP5frYBf/E6LnBOPsj5Qrbf16F2CV/k/OMd7LI4Bj3J8uB5L6nFnOiYBIBh5f1PqSqaxzfLlyzFt2jQsWrQIgwcPxoIFCzBixAhs2bIFrVu3drR/5513cN5552HOnDk4/fTTsWzZMowePRobN25Ejx49AADTpk3DmjVr8Pe//x1FRUV47bXXcNlll6Fdu3Y488wza9UvXoFCREREREREREnjzjvvxOTJkzFp0qTQlSI5OTlYvHhxje3vvvtujBw5Etdeey26deuGW265Bf369cO9994bavPOO+9gwoQJGD58OIqKinDJJZegd+/erle2qHgChYiIKBZ2Av+IiIiI6luCxzbl5eXaX2WlfmXZYVVVVdiwYQOKi4tD93k8HhQXF2PdunU1PmbdunVaewAYMWKE1n7o0KF48cUXsX37dti2jTfeeANffPEFfvnLX9Z6F7GEp765leyo6z0Nc35LvXzfdWriKKYF1rYr9oNa1gIAR2c2Cd1esb+Jtq7pNv2SwYrW4f0kL2OMZkpbI/k/NLI6QulSiy0HtXUlgxppywFldrxAlihNEWU5atmOo2RHTlWsLotDx3KU9NS+hMdWSzZEaQ3EFHxyq0HTSjElotUovK3dvfSVBev1S3x39wzv04C4otS9TicyO4pqFEuZm7DpDr2GZEu1flnhGTnlodtVTcWltX5xGbRaPgO5ztBfRF53qIG63Xoiv8M8yn6SZTduJT1ERGlEliJaolxCm9bYI8tyxNTEyrLtleUcoqQkU10H0TbyslxnixIeKzv8pZ3l0+e7lVMVq2U7+dn7tHWtfXu15RYZ4ZKePG+Fti7Ho283y1DCUyVLeILhsp2ygF5P3cir999rGDMFZclOQH9vDiplObZfvG9V8r2q+XZNy+r76BXvseMYUI4PSxw7lji2tONOHpNJUj5Lqa2wsFBbnjlzJmbNmuVot3v3bgQCAbRp00a7v02bNvj8889r3HZJSUmN7UtKSkLL99xzDy655BJ06NABGRkZ8Hg8ePjhh3HSSSfV+jXwBAoREVEsmFRPRERE6STBY5tt27YhNzc3dLfP54v0iIS455578O677+LFF19Ex44d8dZbb+Hyyy9Hu3btHFevRMITKERERERERESUULm5udoJlEhatWoFr9eL0tJS7f7S0lIUFBTU+JiCggJj+wMHDuCGG27A888/j1GjRgEAevXqhU2bNmH+/Pm1PoHCDBQiIqIYWHbi/oiIiIjqW7KMbbKystC/f3+sXr06dF8wGMTq1asxZMiQGh8zZMgQrT0ArFq1KtS+uroa1dXV8IjSNa/Xi2Cw9jVpvAIlAaxopvKMhmOK2Niexzm9sGm+3igyWxyPjbwcFLkVao4JAHxTHa6JfWrnyca2Wp6Kaw5Lzbej5ZzGOHxHdRPzx0rLNRH7wc6Qy8oTycwTkYliaRkoYp2jw4YOynnGDNPfuubeqLflY0WpsrpfvAf1tnKfqvvbkts1nRaOPEPzodWmL3f5WGW/yGPyqZ0DteWeHf4RbttGb9u4RHxhK6/NNaxEO/ble+xIRYm4zjkFchTUxzqmPI5mQ6anMGftxPpdSETUYOQUsYZl2yvzpORy+DtQTmEbFHkYaj5GUOTGyTGI+rscFJltdqb+BZ+hLGdn6RkiuVl6Vok6VbHMPGmTWa4tt8oILzcTGSi5Hj1zLlPOC6yoFi9uvx3OWsvx6/0zTYcMAAFl8FMtgtgqq/Xnqa4KL/vFPgvKaaOVqYpN008D+vhRvseyrdeQryiPLS2Lx+0Y9Ufe30TxMG3aNEyYMAEDBgzAoEGDsGDBAuzfvx+TJk0CAIwfPx7t27fHnDlzAABXXnklhg0bhjvuuAOjRo3Ck08+iQ8++AAPPfQQgENXvwwbNgzXXnstGjVqhI4dO+LNN9/Eo48+ijvvvLPW/eIJFCIiolgkasYcXoFCREREDSGJxjZjx47Frl27MGPGDJSUlKBPnz5YuXJlKCh269at2tUkQ4cOxbJly3DjjTfihhtuwDHHHIMVK1agR48eoTZPPvkkpk+fjvPPPx8//vgjOnbsiNtuuw1Tpkypdb94AoWIiIiIiIiIksrUqVMxderUGtetXbvWcd+YMWMwZsyYiNsrKCjAkiVL6tQnnkAhIiKKBWfhISIionTCsY0rnkCJVbxq7hO1HbkcTVywoUuO/BTD42QWhfa5EW3Le+v1p9VKg/e+KdLWNRE1ox5zqWrtmV6aI/NEX260O1zXWtYpU1sXyNLbBpVPnS1yTGyZc6LuQ68h8wSA5Q33wZLvt3gBxggUw0NtR0iIqKWVq5U3Xb42R611INxW7rOyIn2f5uwMvwH72+o1ubKLCfu+VrYbEMekPGYDHcLr5bGes1N/bVpOj3gf6/R1Ec1j1ecNuny3xCodckuS6DJXIkotjrw8mSGlZlU4citk9pfyWyvHXjIfw6uuQ8R1gBivZIovJpHnkZkVzsPIERkoTbP0rJJmmQdCt1tk7NfWqZknAFCQURZ+nMhAaWrpz5NlRQ7eqhI7Zq8d/u3NhHkgWS121IHM8GMrsvQBy/5qfflgVritP1Pfji3zapTsEvlemN4rx/FgOj7EOsuRp6Med6KtWx4ZpQeObVxxFh4iIiIiIiIiIhe8AoWIiCgW/FcaIiIiSicc27jiCZT65qitEKsb4LL2hihJk1dZtmlTpi2rU8R5vsvW1smSHePUxG7LMZL9zykJzz1X3qmRtk6WqqjXfTkvtZUlPOrUxHI7oixHWS+nLbY8tZ9P1jH5bTDy5Zy2acpj0UfHlM1yW+p+EfvMn6Mvq/u7oo2+v+vEUGZmOrbkMSmPWfV4lsc67FbR9TEC51TEkX+p5Gupr68A9futTr+jju/ReNXxERElkDbtu0uNpjb1rLmUUi3LMP22Hmpb8+0al9XfZTE+8WSIaYwzwt/Dvgx9etucjCptuYk3XMqaJ8pyWmbs05bVsp2WHr0Etql4rZmGC+urxa9OVrAyQkugGvqOqAj6tOUyb3jcIV+bfO3qfpH7LChLm7X9rffJtOx4j+UxoJbtuJaOmdoajllHDTdR+uIJFCIioljwX2mIiIgonXBs4yqqDJRbb70Vr7zyCkpLSxPVHyIiIqJ6w7ENERER1VZUV6DMmDEjdAl2QUEB+vXrh/79+4f+2759+4R0koiIKOlwqr+0wLENERHR/8exjauoTqAMHDgQO3bswKRJk9CqVSts3LgRzz33HP785z8jEAggPz8f/fr1wyuvvJKo/lKSifWzYAX067ia+vRaVK+S4ZFxQIZRyI3F1gdXdoTbNfAeCE+lZwVEBophOme5zpTZYsnME8fcyrbSNijaiqaO6YiVPsk3VdmWlodSw3ZlLouWy+EyzTUM+0Vmzqj7GxAZKI6plGt+jkRyHLMKeazvr308jYMdefcSUS1xbENJQWY7WYYsCsdyhNsw53c5szPEspq7Ibsn8jsylLGCz6vngPhEUFiON5wbkiNyTRpbeqaIOlWxzDxpYmVqy5mWCApRVNsyrCy83Spbnw65QvRB9lHtv3xt8rWr+0XuM7lPzbkm+rL2PkaT/+d6LGkDTxCRU1QnUN577z0sXboUN9xwAwYOHIg777wTRx99NCorK7Fp0yZs3LgR//nPfxLVVyIioqRh2dp5y7hul+oPxzZERESHcGzjLupTixMnTsQXX3yBrl27ol+/fpg+fToCgQAGDx6MSy+9FA899FAi+klERESUEBzbEBERUW3EdG1WkyZN8Je//AUffPABPv74Y3Tp0gWPPvpovPtGRESUvOwE/lG949iGiIiOeBzbuIp5GmO/34/Kykqcd955KC0txaRJk3D66aejRYsW8exfaki1GsFoupuglyYzLTJFDalX/ZTVV+ZJXdjhTjprUWVWSS3XyfUubU27xZR54tZWzUSRz+H+WtX+yw7LtlbEdbYjbCXJiZcWVPacPNatYIJ+UeRnN1Bjq+jbJkqqfY9SWuLYhupNHX7X5G+v+hvplq1mG/JSohmDeMSyV8lAk+vk716mFV7OsuQ6PUMkSxkwZoofK5l5YspAkTLt8PNkiUGpsw+R+y9fWzT7JRDVmElfNL2Pzuw95fioy3Aq1cZiRAkS1Yh17ty5OP/889GrVy/k5ORg6NChuP/++zFo0CA89NBDyMvLS1Q/Nffddx+KioqQnZ2NwYMHY/369cb2Tz/9NI477jhkZ2ejZ8+eDIIjIiIiABzbEBERUe1FdQXKDTfcgKKiIkyYMAHnnXcejj322ET1K6Lly5dj2rRpWLRoEQYPHowFCxZgxIgR2LJlC1q3bu1o/8477+C8887DnDlzcPrpp2PZsmUYPXo0Nm7ciB49etR7/4mIiCh5cGxDREREtRXVFSg///nP8cMPP2D27Nno378/TjjhBFxxxRVYsmQJPvzwQwQCib/u+84778TkyZMxadIkdO/eHYsWLUJOTg4WL15cY/u7774bI0eOxLXXXotu3brhlltuQb9+/XDvvfdGfI7KykqUl5drfxRfhxOeD/9leQLaXwBW6M8KQvtLyno6ywr/UeIl2/4Wx6Q8ZtXjWR7r8rNAqcOC87ssLn8N/cKOMBzbEBERHcKxjbuoTqC8+eabKCsrw5YtW/Dwww/jhBNOwGeffYZrrrkGffv2RZMmTTBo0KBE9RVVVVXYsGEDiouLQ/d5PB4UFxdj3bp1NT5m3bp1WnsAGDFiRMT2ADBnzhzk5eWF/goLC+PzAoiIKH3YVuL+qN5wbENERPT/cWzjKqYQ2WOOOQbHHHMMxo0bF7rvm2++wQcffID//Oc/ceuctHv3bgQCAbRp00a7v02bNvj8889rfExJSUmN7UtKSiI+z/Tp0zFt2rTQcnl5OQcaREREaYxjGyIiInIT1QmUGTNm4KyzzkL//v0d6zp16oROnTphzJgxcetcQ/H5fPD5fA3dDSIiSmaJKiNkKVe94tiGiIjo/+PYxlVUJ1C+++47nHbaacjKysIZZ5yBM888E7/4xS+QlZWVqP5pWrVqBa/Xi9LSUu3+0tJSFBQU1PiYgoKCqNofEcQUwjDN+hZN22i6ILbjkfMaK9ym5EsKShaHI8fCMd+gXbt1cr1LW3XJMd2weKxpWmPZNtJz1HiHc27FKNpGXpf02SAuUwiq03KbjvW4qqenIUp1HNtQg5C/91Fwzn4bvsPxE2OYKdf5uyyXI49BgmI5ELQirqsWg75qO7xcZct1+v+aVCk/qNWig9V27fOJZFt1W1XiR9vZh8j9l68tmv0S3ZhJXzS9j/IY0I6Puoyn6nDMEqWTqDJQFi9ejJKSEjzxxBNo2rQprrrqKrRq1QrnnHMOHn30Ufz444+J6icAICsrC/3798fq1atD9wWDQaxevRpDhgyp8TFDhgzR2gPAqlWrIrYnIiKqFRloHc8/qjcc2xAREf1/HNu4iuoECnAo2OznP/85/vKXv2DLli147733MHjwYDz44INo164dTjrpJMyfPx/bt29PRH8xbdo0PPzww3jkkUfw2Wef4dJLL8X+/fsxadIkAMD48eMxffr0UPsrr7wSK1euxB133IHPP/8cs2bNwgcffICpU6cmpH9ERESUWji2ISIiotqIKUR24cKF+M1vfoN27dqhW7du6NatG6677jrs2rUL//jHP/DCCy8AAP74xz/GtbMAMHbsWOzatQszZsxASUkJ+vTpg5UrV4bC1LZu3QqPJ3xeaOjQoVi2bBluvPFG3HDDDTjmmGOwYsUK9OjRI+59IyKiI0eipp5O+pK1NMWxDRERHek4tnEX0wmUq666CnfffTfWrl2rJbg3a9YMvXv3xkUXXRS3DtZk6tSpEf+VZe3atY77xowZk9gAOFsWnBou7ImmbQpQPwzRfC5sr173uftAE205oNSB+nNE1sdP9RSCYkW4XYOAL3I4jKkeWa5zlJeq+zcoOuGJnBNiB8Vx5dGfyJRzImnbEo+T23H0Maoa3cjrZKZIMMsQxiNfWj0cLnJ3ymPWo7x4eawHxWfB46/9JyluP0bJkJfi+G4kql8c21CDkt+BWl6azEeTyxFuo6ZcNmWd/B2Wy2psiOxeQP/t8itjhcqA/r8XlSInpCIQzheqCOrBxvttPXtor50Zup0VrNQ74anWFjNtPyKR+Sl7g+HlveI5y4PZen9FH9X+y9cmX7u6X+Q+k/tU3d9u7432Prpm2ai33Y4lbTAJInKK+f/ei4uLMWzYMGzbti10348//ohBgwbFpWNERERJjXXCaYdjGyIiOqJxbOMqpitQLMvCLbfcgtatW2PYsGF48803Q/9aYzOhmYiIiFIMxzZERETkJqYTKIfdcsstsCwrNNDIysqCZdVTeQUREVFDStS/qPD/1RsUxzZERHTE4tjGVUwnUNR/ibn55ptDA40nn3wybh07Uqn7tr6Ga7J20nR8W+Jf4ewYeykzLbZ930K/o6vStvCgvm57o4jbdbwWUx1oHXaw7H9F23B9bMZ+/Un92SKrJFjzbQCwRH2s9jxB8WJE3oitvniPzCbRO2x66Y73Xwn4cJTDyswTw3r52hz1vOp+EWXMGfv15f3twvtbvhd1Ek0duWkz8phVyGO9QMa5RC7hNpKfTXPb2J6jruL2r/isy66Tl19+GTfffDM++ugjZGdnY9iwYVixYkVDd6vBcWxDCad+BwZdQtDU33z5+y/aqt//lmhryjXRMk5qWtZ+l/Xf8KBf//H1+8M/ZpV+/X8vKvx6xsi+QPg3vCyQo63L8es5J5kQnVJU2XoGSpZjYKG21fur5p7sEX2Qy7KPav/la5OvXd0vcp+ZxkWu740xLyXy8WE8zuSyo61bcB/RkSGmEyi33XYbGjduHFqePXs2AOCMM86IT6+IiIiSXCom1T/77LOYPHky/vznP+OUU06B3+/Hxx9/nLgnTCEc2xAR0ZEuFcc29S2mEyjTp0933Dd79mxkZmZi/vz5de4UERHRka68vFxb9vl88Pl8EVq78/v9uPLKKzFv3jxcfPHFofu7d+8e8zbTCcc2RERE5KZOGSjSjTfeiBtvvDGem0x/cbuk3a12JZptxedx0UxTm/uR+J+Ck8M3hx39pbbqwzd6assBn1JiYpjdVnKU+xjqWuQ6+TwHWoUvy8zdqtdgVDfN1JbV2fDkdM6W2K7tVTopy3sc/VWmypOX98r+R1HDo01VLMuG5JR8hin6HJeqisuBPcpu81bpm8n7n36ZbvlR4a8u+V5EMUNznc6Gq8ewt1LfkDxmvcoTOY51UY5i+pzUqX405s91HL9bErGdhmRb0R1w0WwX0KbSBYCZM2di1qxZMW9248aN2L59OzweD/r27YuSkhL06dMH8+bNQ48ePerS47TGsQ0lguN32lFKEaz5NgArIB7rtyOvE8uegFpaqz+l/M1Rf5flb7ZdrZejVFeFf5crqvRxz95MfVrgPRnhH/lGXjEVseyU4qCtb7fCo5fLZsr6X7V/tv6/POp0yT/4m2jrdvtzteUf/Y215T3V4VLyvVViymPx2tX9ArHPrOrI4yCPSwmPRyvFMr/npuPDUZajHXeibD8dfrfJXYLHNukgridQiIiIjhgJDlrbtm0bcnPDA/m6XH0CAP/9738BALNmzcKdd96JoqIi3HHHHRg+fDi++OILtGjRwmULRERElNYYIusqnvGLREREFCe5ubnaX6QTKNdffz0syzL+ff755wj+/39Z/L//+z+cc8456N+/P5YsWQLLsvD000/X50sjIiIiSkm8AoWIiCgGyRK0ds0112DixInGNp07d8aOHTsA6JknPp8PnTt3xtatW6PtJhEREaWZZBnbJDOeQIlVQ9QByvpYTxR9ULMoPFHUNDryD2TGhbLeK4M25LaUPogazEY79RrMdhnhbf261Qfauv9+101bLuuk5GF46lBfJ/sbVSZK+I7MfaIG19LrYT1KzWswQ9apiocqr8ft3dbeKrkf5LFi+hYTL85Wc09kJodLBopaM+14bbKeV60FFpuV+9T2hvdpVCWVLjvR+OXuyPgJ39F0m17DLY/ZNt7wxX7yWJd1zqaplB206Qbldmr/uXZMgRzNLMHatJwNNj9ywzxvksjPz0d+fr5ru/79+8Pn82HLli048cQTAQDV1dX49ttv0bFjx0R3k4hUgYB5WcmisAL6l7LtyK2InHHhMeRjeMTvskdmcijLHpFNJjNQAsryQZEDUp6pXz2X6Q1ninhd/q+qWgk6qwjq29nj0ac8zjLkp1SJwDR1W3KaYpl5srOyqb6+Mry+vErvk3zt6n6xxD7zVEXe33LM5PEb3keXDBR9SmyZp2PIQHE7RomOUDyBQkREFIsUqxPOzc3FlClTMHPmTBQWFqJjx46YN28eAGDMmDGJeVIiIiJKHSk2tmkIPIFCRER0hJg3bx4yMjJw4YUX4sCBAxg8eDDWrFmD5s2bN3TXiIiIiJIeT6AQERHFIkF1won8V5rMzEzMnz8f8+fPT9yTEBERUWpKwbFNfeMJlASQmSJahaPMCNBLMo11ivDUz6RJWh6CzMNwtK35ds1tlRpdsd2sffprffNAy9DtYzJ3a+vKO+qHrS2zVwydsCK/NDPZVkbBKO/jj12ztXWtPtQLh3f3zAr3QfZdvsWG+BHH95D6gmTOjSWeJ4oMFO2JgvKFi836DfXTYp1X1P56D4Rvt9qs7zO5T9X97fbexMzli17N2ykv0mueO2f8qC2/eaBV6HbWfn2nBTNEbowhz0O+baa30dm2IXKbItfpO9tGXmfMaSIiSkG2+M6zZb6EP7xsydwKv8jSUpZlVoanWl/2KpFdAT2+Cx7DslwXrNR/u4IZ4R/mqgz9N3GfV/TX8OMVEGOQA5nhbZV5G2nrcrz6WCHTkIFSLTNQAuGx2L6AyFap1p9HzTwBgJ8OhtfvO6g/tqpSf+32wfDzeirlGAkRl03rAP19lO+x4xhQjg957MjfaVs57uQxKY9ZoiMVT6AQERHFgnXCRERElE44tnFVP5c0EBERERERERGlMF6BEgfykjZLluWYLj93TE1seiJZI2NHXif7pJR0OC6HNz1WvBZZBmCr9RKOqVHl06ht9XX+RvoLv+rf54Vu/2P4vdq6fafs15Zz3moSfk4IppcqSxxMD5VVLbKxsp+CWaLsSUwpnLVXeS9ECY+j1CZih8T0wgBsr7pd8VhHmUvtn0efflqUmzimLdYfqpbtyEtX1ZIdQOwXsc8c+1R5fa6lWIYyM+NrjaJcRh6T0lXvjAvdbiOOdW9V5No3R9mNsbzH5dR+0PB9YZrWWH5HRfM9VNv+1MSwraS5jJj/SkNEsQq6TWOsLPsjl/ccWlZKeKpluY+YOlcp7/DKshwxVXEwM/wD6xXlprb4Pwg7I/w8Aa8+CDlgZSGSoPgRrw6IUpus8GNzMvQO+jz6fsj0GEp4gvp2K5XlCr/ev71VetmwnKpYLds5cEB/bKBC3zFWZXi/eOU4SJb0KC9Pjg28hjIdy1G2Jcp01LIdt2NJPe5MxySlL45tXPEKFCIiIiIiIiIiF7wChYiIKAZWgpLqE5J+T0REROSCYxt3vAKFiIiIiIiIiMgFr0CJ0arg06Hbp2aM01capio2TnEMALZSp2jX4fxWvHIKHPkokae4dWae6MtWMHJbWRfXfF24pjQwXH/O63q9pi0/+MrZodvVjeV0vZGfx9HfOkxrrL5VQX0GO+zpot/RavPB0O1dffU6W0mdzs8S25WlqGo9sswQcRxKUUxjrL1vbtMWywwUZRpjmXkiZh9E8y2Vodu7e+r7xVFrrb6eOkxbbDwO5Trx2jP3hxvIYzIoOtX8nXC9tMwqcWStqM/jlsOiTTEt15myl1zaGr4vrGjyRxzblU+sNo1timNA/z4mothpYxvPmAbsyZHHDojsEjV/wq//uFoit8KqVqY8rpY5W/qyrWSZqFMPA0BAjDO86rhCZKs5M9yUDBTxu+yYmlhZDgT0/lVW6z/4+6vD40Ffhr4ffF592TQ9ssxaqQyEn6fSrz9nRZW+Iw6KZXWqYkfmyUGx/yuUDJSDMgMFEZcd68SYyaNkpHirxLFTbTg+HBkoYuAWUKcxjvybTYnHsU3y4gkUIiKiWDBojYiIiNIJxzauWMJDREREREREROSCV6AQERHFgEFrRERElE44tnHHEyjx4KjrF4WiQTVEQoZnyMcaLgoy5BTIXIKo4lNkroXyWEcsgSO/IbwsX4oVNORwiO16AvodTbeGazC/rW6hrRvY6Ftt+e5W4efxVMv+ij6ZomBk28irnLkbyrJj34vDQc33aL1RDwb5oZue/eFXMkb8OSKbJFPvVVAtOhZ9cOS7yONQaxx5WWagODJPqvTtZigvL2O/vuGWnx3Ulnf1ahR+Sllr7chwQWSm3BC399wUwyFqww8qx508Jr/1N9eWm24N7yh/Y/nmGPrkyGGJ/Plzfj9AtEXtxZqTFNXjxMEUdFk2PZaIKMXZAT2bwq5WfmBF9odVHXnZqtZ/QD0iH8OTqWRyiHUZGSIvRc3wc81WU9frK8UQD4Fg+H8/Doosteoq/X9NDmaFf3wzMvR9lOHR++/1RP4NCogxqT8Y7qPfr+8z2YeAyJWxD4bbW5WRM08AIEPJPfHqw54aMlDUXBP9tcj3Sl2W77FVLZeV40UcO7bIQFGPO3lMEtEhPIFCREQUqzT6FxUiIiIijm3MmIFCREREREREROSCV6AkQjRlOY7HKiUxblMeq5fzy8sWTdOSOkp2RPmPWhrk2I7ohdIHS17eaZiuVZaCyEtDK5uHL4+88u3ztHVrTrlbW+51zqeh218+0E1bV91EdDdDrbXR1zmmNdZWotYc5TLeGpsBAHb3bKQt532t1yCpl2ju6u3T1skpBANZyjrTtL9A7NMYi5IdOa2efF9bfRi+PjWYqXdCvnbtUmFZwlOHqYpVbtNna+VK4spVX5ne+Jgpn4VuNxY7Zuq/zteWC5TjWb1EF3CW31mG0jznlMd25HWmkj/H9OSIzPRdAhinGHZMTRxNiY/2uCQt2WFSPRHFS1CU8CilFZacalZOa6yUXXiqxA+oKMvxqtMYi3GE12soB4+qFFiMI8QgJKCU7dh+fZ0/MyiWw6/HkyFKU0R/1We1xI+inMbYVsqeg6IPECU7cmpoT6VSllMp9qGcqlgp28nQq7aRcUDvY8ZBpYRHjBU8clkr4QmIdYaSL5djSSvpCbKE54jEsY0rXoFCREREREREROSCV6AQERHFgEn1RERElE44tnHHK1CIiIiIiIiIKKncd999KCoqQnZ2NgYPHoz169cb2z/99NM47rjjkJ2djZ49e+KVV15xtPnss89w5plnIi8vD40bN8bAgQOxdevWWveJV6DEwarAU9ryL7POi9CyBjI/wGtYJ6f1VOtPo8gpMGWeOB4r+yCzVoJRZIqoWRrypcinqQ7f0faf+vyx7X6pZ4H8se0/Q7fP7aZnoDTaqdeiatPuuuwy7ZHmGWJF4a1LW+U9FrPqoexo/bV6lGns2qyv0NZVNcvSlss6hT/OAT1exJHnYRtql01niDMq9JW5W/Xa2cwyPcPlx+7hjgSz9OeU0wKrJdLOaZcj98ntuNOONbe2ynJAn1EaP3XQl9XjroVHPxfddqX+4oxZJTLeQ/2cuGW2aK/NJZvEMN2wMxPFND2y6ftCTkVsmhe69v8UYYu28ju3wbBOmNLYquDT2vKpnjEN1JMjkzadrFf/bUWGPny31OUMPQPFqhJ5HkomitcrvrMdeWnhO5zjhsgZeY4xXkBkoqgxG1VibCDGCraSnxYUmSfO/tY+3037rZX9qxZjR9FHj/J2ODJQ5NTEagbKQTviOrk+46CYtrhSTEddGR7YWSLzRE5VrC7b1fqxZFeJZflYqjfyO7fBJNHYZvny5Zg2bRoWLVqEwYMHY8GCBRgxYgS2bNmC1q1bO9q/8847OO+88zBnzhycfvrpWLZsGUaPHo2NGzeiR48eAICvv/4aJ554Ii6++GLMnj0bubm5+OSTT5Cdne3YXiS8AoWIiIiIiIiIksadd96JyZMnY9KkSejevTsWLVqEnJwcLF68uMb2d999N0aOHIlrr70W3bp1wy233IJ+/frh3nvvDbX5v//7P/zqV7/CX/7yF/Tt2xdHH300zjzzzBpPyETCEyhEREQxOFwnnIg/IiIiovqW6LFNeXm59ldZWVljP6qqqrBhwwYUFxeH7vN4PCguLsa6detqfMy6deu09gAwYsSIUPtgMIiXX34Zxx57LEaMGIHWrVtj8ODBWLFiRVT7iCdQiIiIiIiIiCihCgsLkZeXF/qbM2dOje12796NQCCANm3aaPe3adMGJSUlNT6mpKTE2H7nzp3Yt28f5s6di5EjR+K1117D2WefjV//+td48803a/0amIGSALJeX80IkFkDMn9EzUpwFJHatT/fZYk+2Gp2iWteirIdj7mp9trkc8oaWI9SJxqUr1u0VWptZVbG8cv+oC2/O25+6PaVo1/S1j31p9O05fJC5ZB3Kec1rTM+1CW/Qy3DlSW5jhJdb/iO3X1y9M0G9E41LgnXw+aUVGnrvAdE/bT6WK94UlHnHMgO77OKtnr+TPlR+leI7dXfrKBXXac/jeNwNuWcSOrHxC0nxK7lOrHc7Gu9Dnjs717VllsoNekDn7xWW9cqQ99whlLnLD8njmV10aWtnmuir3LsF1NeimHZ+JySYzuRM1GM2Us1rU9GSVQnTERpJhj+TZe5FJbMrfCGf2Atr/hxFcseNbPLMV6p/f8WWI4BjLJOBLzJHDaPP/IYL1glljPCbeW4QY4roslAUfskx6AeEQPikZkoyu73iP46MlAq1VwTOTYQywfCHfHKDJSDeqc8yrLMQHEcH1VKJ90yT4LizaIjT4LHNtu2bUNubm7obp/PF+EB8Rf8/1l9Z511Fq6++moAQJ8+ffDOO+9g0aJFGDZsWK22wxMoREREseAJFCIiIkonCR7b5ObmaidQImnVqhW8Xi9KS0u1+0tLS1FQUFDjYwoKCoztW7VqhYyMDHTv3l1r061bN7z99tu1fSUs4SEiIiIiIiKi5JCVlYX+/ftj9erVofuCwSBWr16NIUOG1PiYIUOGaO0BYNWqVaH2WVlZGDhwILZs2aK1+eKLL9CxY8da941XoNQHU1mOnBZOu+RdP7/lKPcxTEtqXDaU7AD6NMfG5wREPYosN5D1KMp2PebLO9XSIXVKYwBo/6be4S3nhKfKPSVH/0Dcc+FwbTnnjSZKf8w1I+qlom7VVOqWHCdtTU8j335DWUtAXKoq9+/+tuEGFW3EPMaQy7FxXD7rsgxDuVKsJTuAuEpXrpPTAps+foa2uy44oK36ec6X2vKnVS1Dtzus1S+Brc4R00aqfZBXxzr6r5bPiHWOz2Pkdcapih2vO5rvi8htXctu1O83x5THUUyBnCQSFfjKEFkiUtl+WZIhfmOUEh54xDoxzvSoy47f7Ch+mMUPvqVNC6w3lSUxgerIJTx2ht6H6EqBDf03jBUc/ZVjUkdJT/i2t0rfsCzpyVBKeNRyHqCGqYoPqiU8eic8laKEpzLcCeugeNJKsayU7WjlPHAeW0TJNLaZNm0aJkyYgAEDBmDQoEFYsGAB9u/fj0mTJgEAxo8fj/bt24dyVK688koMGzYMd9xxB0aNGoUnn3wSH3zwAR566KHQNq+99lqMHTsWJ510Ek4++WSsXLkS//jHP7B27dpa94snUIiIiIiIiIgoaYwdOxa7du3CjBkzUFJSgj59+mDlypWhoNitW7dquU5Dhw7FsmXLcOONN+KGG27AMcccgxUrVqBHjx6hNmeffTYWLVqEOXPm4A9/+AO6du2KZ599FieeeGKt+8UTKERERLFgBgoRERGlkyQb20ydOhVTp06tcV1NV42MGTMGY8aMMW7zoosuwkUXXRRbh8AMFCIiIiIiIiIiV7wCJRHk1J1a+EMU56zcsgfUzICAeSpabVtyalSPvmwbtmuJAjZbXZY1o7KtWm8q80dEW216OfFaqpvo+/CKv1weuv3s9HnauicH/lVbnvDO1YjEYygDdWSeGDJRopkdOZopj+WGHH2KYhrmqJhKit3KpeuQc6JtxjDdsCnHRK6XbWXtteqJAfqxky2KpG+YHz5z7cvTn1TWOatTTjunMRZPHIzcFoHI3wlyWmtn2/ATOfJSjNMju+Qrqd8X0Ux57Eb9wpDfqckiyf6VhojSlPguDYqpaD1qBoqYttgyZKJ43DJPtN8C8b8McupfZbzo8evPGRDjKzUnxDGNschAsaPKQEFkUWSgODNcRM6JmoEicvq8clrjqvATeSrlWCHyVMWOzBMxjbGWe2Katlgsy2OnTr/TlJ44tnHFK1CIiIiIiIiIiFyk1AmUH3/8Eeeffz5yc3PRrFkzXHzxxdi3b5/xMcOHD4dlWdrflClT6qnHRESUrg4n1Sfij44cHNsQEVGy4NjGXUqV8Jx//vnYsWMHVq1aherqakyaNAmXXHIJli1bZnzc5MmTcfPNN4eWc3JyEt1VIiIiIlcc2xAREaWOlDmB8tlnn2HlypV4//33MWDAAADAPffcg1/96leYP38+2rVrF/GxOTk5KCgoqK+uYlXgKW351Ixx4QVRayizCOyAklNgRRG84ZJToGYp2I4gCJlHEi4ydeSjmPovAzEcfVDWi6wEW9ThahkoYjfIWtSc3eHlEesu09a98LNF2vLVU54J3V644Bxtnb+x3gfvQURkzESRNbiG67xcT8ZaEW7XtK0oan+jEk2OiYlLH7Qz07JW2ZCBYqprlsuBbH1dxn79wX+46tnQ7RyPXm98+rpLteXWuyPncjj6oNZTO/oXORNFfj/I/WAZskqcOSfqbfkdoC9rfXLLQNH6ILcb+ftNbkd+t9jKY+V3atJgnTDVUUqNbYJPa8unesyzHFACBfWQDruyMrzgMV9cbhpWeEy/G47fWj2QJBgIP6+l/3zC4xfjK+X/PgIiA0XmnAQNGSjG3DjZ1JSlJjNPHJkoIrtEiRFx5KNUyUyU8E70Vuobtqr030yPloGiZ5VomSeAlnuivf8AbJmJoq4Xxw41HPmdmjQ4tnGVMiU869atQ7NmzUIDDAAoLi6Gx+PBe++9Z3zs448/jlatWqFHjx6YPn06KioqjO0rKytRXl6u/RERERHFE8c2REREqSVlrkApKSlB69attfsyMjLQokULlJSURHzcb3/7W3Ts2BHt2rXDRx99hD/96U/YsmULnnvuuYiPmTNnDmbPnh23vhMRUfpJVE1vOtUJkxnHNkRElEw4tnHX4CdQrr/+etx+++3GNp999lnM27/kkktCt3v27Im2bdviF7/4Bb7++mscffTRNT5m+vTpmDZtWmi5vLwchYWFMfdBI6f5dFzrr1x/aJwOWb/k3XHVonwej+GSfMO0xra81N8wrbEtPxmmaY0d6yDuUG7KSgnRWL28s83yRtq6mR3O0JbnFL4Qul16+evaumfvKNaWD7YIP49XvzrSwTGVnsJReaW+NvNmo7vazbSxeJXhuDF02PWLU72c1lSyI9saSnYAvWwn+wd9Q+dcox8Dwxr9N3R75ve/0ta1eUo/ttTjznG5r/xMqZ9VuU5cGqytN0xb7HisYdpi2QfX6YZNbQXbOOWxoVTIUUaURr+sdMQ74sY21KBsv1Izc9BlwKJwDL3Ed7hH++0S5SZ+/X8hgtXhH0VPtZzGWF+2lamKvZlyTCenMQ4vy7GWo2QnmmmMTeMI8XvqkcvK1MWWKOHxyrIcdRrjKlnC44+4bIkSHlSKqYmVshxbrhPHgHZ8EFGdNfgJlGuuuQYTJ040tuncuTMKCgqwc+dO7X6/348ff/wxqhrgwYMHAwC++uqriIMMn88Hn89X620SEdERiHXCFAHHNkRElJI4tnHV4CdQ8vPzkZ+f79puyJAh2LNnDzZs2ID+/fsDANasWYNgMBgaONTGpk2bAABt27aNqb9EREQAOMigiDi2ISKilMSxjauUCZHt1q0bRo4cicmTJ2P9+vX497//jalTp2LcuHGhlPrt27fjuOOOw/r16wEAX3/9NW655RZs2LAB3377LV588UWMHz8eJ510Enr16tWQL4eIiIiOcBzbEBERpZYGvwIlGo8//jimTp2KX/ziF/B4PDjnnHOwcOHC0Prq6mps2bIllESflZWF119/HQsWLMD+/ftRWFiIc845BzfeeGP9dlzNAXDkmIhMEXXBkZdiyEMIuISKKLWrllhnewzZCZY+Z5xpWmMrIItIxfk5tYZUzpwsZ1VT+iiniHNGooTbBrL0DX937zHa8l3X/CJ0+5rWq7V1Ta7T5y1+5C+nh25X5unb9ZjKSV2m1dN2v1tbw9M4an8b4MxuVIFQdZmaOIoMlKD4VvP9GG484bqXtHW/aqxnEDz849DQ7S8Xdtf7IKZaVI9LR720OGa1ablNmSdivXOdDNRR2pqmnwT0z7VhmnPHerld+V2jPlZux5BrIqctdmY+RZ4mOllYSEy8UH1FFlFySNmxDSUlu7rKvVGosfwtMGRnid8fyy++o5VMFEtknniqxJTHmeH1tsg8UdcBeu6J7TEPkqLJo9PHES7ZZIbfeE+12C/VIuekWs1AEYPHapmBouSeyFyTKrmsTGMsM0+iOQaIBI5t3KXUCZQWLVpg2bJlEdcXFRVpg/LCwkK8+eab9dE1IiIioqhxbENERJQ6UuoEChERUdJgnTARERGlE45tXKVMBgoRERERERERUUPhFSj1Qbn01pZ1ll5DW1mb6thuUGkrslQceQjKsinzBNBCOhyZJ7JtQGkrs1VkjoWWlSCfUma2KPtBrJN5KR5DEEcwU3/sprl9QrcfnbFXWze+2Xptuf3//RS6Peue8RGfA9Bra2UNrqN7SpeM+SiirfNJY1xXF1GePdZeuynzxKWtYx8a4jFk2/+74bHQ7e5Zpdq6R/cM0pbf/POQ8HbEZ1PmmnjUrBKZ4SOPb+V4lnXYjtprddltP6h12abadUBklUTOUjm03ch5Ka7ZJca2hkwUUw5LkrLsKHOAotguEVE8OPIw1O9s+VvgyLFS8vNkBpdf/+HzKMt2pv6/F3am/lhPRnigFMzQB01ekXOi5p7IvBRnjlzkwY/jN9GOvM7yG34/AXjU/Be/WwaKv8bbABwZKFrOiZqHItcBsCvDuSe23xTMRxQdjm3c8QoUIiIiIiIiIiIXvAKFiIgoFqwTJiIionTCsY0rnkCpB6uCT4dun+o9V19pmn5YXooo26qXOcrLLk3blZf6y9IbtcYgKC+XFNMaK23llfyyrMFWpzV2dE/2X10HF8qUxy6fTrWkZ82NJ2rrfpjZWFuemr82dPuhq+/W1p2/brK23Pp5X+h2eUd9H2Uc0PuglviYynsAUeIj3wrTS02SEh79Etnab0uWy/gb6S8od2v4ACn9tT593+NDHtaWm1rhy2DnlfxSW7dlzvF6l7KVY6k6cskOAFjKFbPGMhyYpyZGQF5WbJgaXH521WXTtMVi2dEHU0mPadpi17aG/su24gtE/d4kouSjjW08YxqwJxQNtdxDlk56xPeypSzbcp0/U1+uVqYxliU8GaImVlm2RAkPPGIaY69hXOwYvyIi5xjErvk2nNMWO357tRIeuV/Eb5tapiNKbexqvUxHLduRJTtBUdKDoKwdplTAsU164AkUIiKiWKXRv6gQERERcWxjxgwUIiIiIiIiIiIXvAKFiIgoBkyqJyIionTCsY07nkCpb3JKUEOegCXqTW1RF6rlGIi5UY3TGsvnlNMaa5kGLnkpSp8cNaPykyKfV+V4rNI9MTtbUB61apcCev9MmSjVOfr+/HBWX2359Au6hW4//7MHtXXPDl2kLb/Uo3fo9qMrTtHWNdqlP29lMzVnQ3RKXhMW41TFpjrgunD98jOsd5uK2Fbe18z9+rqKAn155Iw3Q7fPzv2Pti5bhO/86t3LQrcL/p6trfOLY0DNPXFknsjcHnVqYtFWZriY2pryU0xZKgD07wvHFOOGXBPTFMdiveu0xUFDW/kdZpjSPRWmLSYiSisiRyNYKcZ8AXUaY/EjKKfOzVAyUKr1gZqVIQZuagaKI/NEjHXV9R5zBkpUTONimXliyiOTmScy58RvykCR0xgrGSh+MUDkbyRR0uAJFCIiolgwqZ6IiIjSCcc2rngChYiIKAa8zJWIiIjSCcc27hgiS0RERERERETkgleg1DM5//ep3nP1BmqNqawZdWQcKG1lXais55S1q1pbkZei3LZlfakhE8WWoRYij8QynXqU69SyULEbHJkokbfqCAMxZaL4G+lt85/LCd2+4K1rtHV/uvIJbfm8vA2h2yMu3Kytm/f9SG350+eOC93O3KP3p7qJzK5Rbov8C9sra4Ej3E4gEb3jyOxQ+yhzWbL2620Ptgw36HjJF9q6P7V/VVtu560K3f6wqqW27oY7LtKWW+0Od9KvR6A4skrU3BNLHGemPBLXXBPleUyZJ45tycwTR1125LwiWS+t1XDL7xLTd4tbVonW1vhp1NuKLBX53ZgSeJkrEYAaxjaeMQ3UE6oTmdNXHf6ttcU40vIbck68Xm0dMuVATslAydDbykwUJEEGiiPPS/s9FeNrmWuiPNZ2yUBxbItSHsc2YrtpglegEBERERERERG54AkUIiKiGByuE07EXyKsXbsWlmXV+Pf+++8n5kmJiIgoZaTa2KYhsISnga0KPKUtn5oxLrwgy2XEY/UphM2X/dnKZY/O6U4jXxIpp25zlPSo27LEZZhyamXlSkXLce5OXvrvibxKtlS269IU8CtTCMtLVUVJjK1M75y9R9/y/X/SL0++8cRwf1eNma+tu7nDP7Tlg1e8Err9/oEibd1fPvqlttx4bePQ7dyt+qWe+9rpH19ZgqRyfGmpy44DSz44fDMorsr1VumNm2zXj8PyjuE+Vg7bq627quc/teV+2VtDt5uK+plm4pLegU9eG7rd7m39OX05onRFLYOSJUeGqYpdy3KU9Y4pj01TE8vPtSzT8SvTRjpKbSCWI5flOKZdNExPLh+rTUcsLyk2TE1svMwZ+tTF8ruPEm/o0KHYsWOHdt9NN92E1atXY8CAAQ3UK0pHLOlJQ+K3wK4Uy0o5iiVKeKwqMeZT18tyH7mslu3IMahVh38HVstIZWmqSxm8Vs7k9runlbmyRCfdpWTJDkWNJ1CIiIhikWJ1wllZWSgoKAgtV1dX44UXXsAVV1yh5VkRERHRESrFxjYNgSdQiIiIklB5ebm27PP54PP54rb9F198ET/88AMmTZoUt20SERERpTNmoBAREcXCTuAfgMLCQuTl5YX+5syZE9fu/+1vf8OIESPQoUOHuG6XiIiIUlSCxzbpgFegJButJlNeUi3Od6lTo4mWlqwLNbWVU8ap9ZumzBPHc4oMEdFfdRpjW+Y1iLa20lZOJ2s6al2nOFbzMGTJq6Nt+PV4K/XGgUz9tbb6T/j2bzf8UVtX8stqbfnenz8euj200X+1df/42QP68wwOP8+W6tbaumd26ZkF677uFLrt+V6frzejQkyPrObRyKgMcTj4m4Rfe7DjAW3dyUd/qS2f0/IDbblz5o+IpFo80bf+5qHbU9+6QFvX9jX9TW/lDfcp4NO3I3NZ1Ncns0nkNMaIoq2Wa2KYtvjQdtUMFLlOfBbUZUe2iswvMtRwO6YbjtxWfh5N0xjL+m7jFOm2PLhcU4pIsW3bNuTm5oaWI119cv311+P22283buuzzz7DcceFp1D/7rvv8M9//hNPPcUsGiKKAyXfw5Z5KdWirTqNsZyaWOanqGNLOV6VTKWI8jdRZcoBA5w5J+rvInNNiI44PIFCREQUg0Slyh/eZm5urnYCJZJrrrkGEydONLbp3LmztrxkyRK0bNkSZ555ZqzdJCIiojST6LFNOuAJFCIiohSWn5+P/Pz8Wre3bRtLlizB+PHjkZmZmcCeEREREaUXnkAhIiKKRYom1a9ZswbffPMNfve73yX2iYiIiCi1pOjYpj7xBEqSWRUI16OfmjFOXxmQuSHqgqj7dGQPKOtd5rfXBEV2hlitbSng0lbN3ZCZJyKBRM0yCWaILBWZiaKuFke0zESxldJaW2TMeGQuhLJdWy/Jhe0Vr/Wgku8iSnQLVun/wjvvqQtDt/cepXf4p6GV2vI9JywL3S5utFtb1739y9pyoJ2S2SKukws48nTCguKd8hi+4eR2M0XbozJytOWvlf1/5ntTtHU5a5poy7n/CzcuaC5qoMXhnFGtZIqI7sqMES0DxbDu0HLkrBJTzolcJz9jHr/yRDLzRGaiBNTtGjJPxLYswzrHsvzMy/6rnwX5XWLYrsxHscV21e83ajh/+9vfMHToUC0ThSiRVgWfDt0+1TOmAXtCSUHLSxHr/PrATfsVSdR066Z8FCIX6vcbHRl4AoWIiCgGlm3DSsDAOxHbVC1btsy9ERERER1xUnVsU584jTERERERERERkQtegZLMTGU40EtQ5Fk9eSm9Pq2xmI5NTCGnT6NqvlxSnV7OlvUQsq36nH65LnJJj2NqYlnSoy7I8h5ReqNdC+oRZ0K9sgxKeQ5Z4SDfGkNpkCzv8DcO999Xpq9rszJLW779ufGh27c00V/3gdb6cnmvcPlPYTt9+uBWjfZpy0GlzqgqqO+karG8tzI8dWppaZ62LvdDfVrVnFJ9x2TtDS+3zIk8VTWg7xc5bbSjTEcrtYm8DoBxamJjuY9rW3UaQ0PJDiBKbUT/RFv982cun9H7UPtyH8cUjXIaRmW7jmmL5XeNqdwn3actZp0wEVH9SaN/wSZKWhzbuOIJFCIiohhwqj8iIiJKJxzbuGMJDxERERERERGRC16BQkREFAte5kpERETphGMbVzyBksTklJ+nes/Vlo3TGMtlbdpSr2Gdfny7ThinPI9ximP5MNnWkIlimuIYEFMKy+l6ZR6J0lZmlcjSWlOuiWysZmfIaYydy8p0wyI7Qz5NMFNpW62vyynR90tOaXi6ZH+gjbZuZ6A1InG8NJmJo/S/QB5mLlPyBrKtiOtkBI2WOSPzRwzTDcsDrU5TE6vLUbR15JqY9osp8wTQs0sc2zFNTSy3Iz7X0UxjrK43TYcs2tpiu5y2mIgkOeUnpzUmolTCaYuJJ1CIiIhiwDphIiIiSicc27hjBgoRERERERERkQtegUJERBQL1gkTERFROuHYxhVPoKQQRyZKxrjwgsg/MOWRuOaaaBsyZ3/Ar2QeiKMpXpkolswqEddNqXkpjiwHrwwgUW7L669ENIzeVuRfeGSflGWZlSG7oD7UkT8i+mCId3FkwUT1xipP4cgQqf03nONyPLkt9XhxySrRMlDkcSZzTQx5KY79H4zc1piJ4mgrPmPqaplr4vichNcbM08ALWPEmHkC6J8/R4iPIdfELS9FayvX6Y+1lf3EzBMiihYzUYgomTHzhCSeQCEiIooB64SJiIgonXBs444ZKERERERERERELngFSioLKpfWy5ISOX2vcjm/nGpUlqOotSDycn1HlYhXqXuRZQGC+li3k5B6W49YJ0oI1HIfr+yh6JNtqIkxTXkcNJcRqSU+jvIeuVvUxzqmBZZdUqeJljUwom0DnA41leE42jpKTAzbktNEy+2qUwibpjiWzxvNNMamdWJbzrIc+VjDdMOmqYrdyn2012YutdHaioNSfieobW3H+yb7Lx57JGGdMBEREaUTjm1c8QoUIiIiIiIiIiIXvAKFiIgoRulU00tERETEsY0Zr0AhIiIiIiIioqRy3333oaioCNnZ2Rg8eDDWr19vbP/000/juOOOQ3Z2Nnr27IlXXnklYtspU6bAsiwsWLAgqj7xCpQUpk6rdar3XG2dcRpjrz5fr2vOSYTtuLU1qdMUxzITxVanchXnBGUmih3jlMfiVKzMWlEzUhzTLDuiVpRtybwUmYliiGyR5PPUC5cz1NpuiybXxG16ZK1t5GwSx7aimZrYkU0ip0tWs0pk/wzZJabME7nebRpjJbvEkXkis0lMbcU+NE9jrC8f0dP72bZzuuh4bZfoCKWNbTilMRE1AI5tkmNss3z5ckybNg2LFi3C4MGDsWDBAowYMQJbtmxB69atHe3feecdnHfeeZgzZw5OP/10LFu2DKNHj8bGjRvRo0cPre3zzz+Pd999F+3atYu6X7wChYiIKAaHp/pLxB8RERFRfUumsc2dd96JyZMnY9KkSejevTsWLVqEnJwcLF68uMb2d999N0aOHIlrr70W3bp1wy233IJ+/frh3nvv1dpt374dV1xxBR5//HFkZmZG3S+eQCEiIiIiIiKihCovL9f+Kisra2xXVVWFDRs2oLi4OHSfx+NBcXEx1q1bV+Nj1q1bp7UHgBEjRmjtg8EgLrzwQlx77bU4/vjjY3oNPIFCREQUCzuBf0RERET1LcFjm8LCQuTl5YX+5syZU2M3du/ejUAggDZt2mj3t2nTBiUlJTU+pqSkxLX97bffjoyMDPzhD3+o1e6oSUploNx22214+eWXsWnTJmRlZWHPnj2uj7FtGzNnzsTDDz+MPXv24IQTTsADDzyAY445JvEdrk8yw0BmPUDJLZBhGWJZi/5we9oo2mp9FEeeI+dEbSuzM7zi/y6U0BEbMtPCitTUkWNi2fpjbXW9I8hEdFi5Ls2SbWW0ipJ7YsvsDLFZdbuODBRP7UNP3PJTTKK65E4ed6bHylwT5X22ZD6N2E/qMWHKUnFsV+aYGPJTHDkm8jhUt2XKPBHrna9NtPWrmT6yfzKPRM1WiZx54mgr1hmzVWR/mc9BFFcc2xAR0ZFi27ZtyM3NDS37fL56e+4NGzbg7rvvxsaNG53/vxaFlLoCpaqqCmPGjMGll15a68f85S9/wcKFC7Fo0SK89957aNy4MUaMGIGDBw8msKdERJTurGDi/ujIwbENEREli0SPbXJzc7W/SCdQWrVqBa/Xi9LSUu3+0tJSFBQU1PiYgoICY/t//etf2LlzJ4466ihkZGQgIyMD//vf/3DNNdegqKio1vsopU6gzJ49G1dffTV69uxZq/a2bWPBggW48cYbcdZZZ6FXr1549NFH8f3332PFihURH1dZWemozyIiIiKKN45tiIiIdFlZWejfvz9Wr14dui8YDGL16tUYMmRIjY8ZMmSI1h4AVq1aFWp/4YUX4qOPPsKmTZtCf+3atcO1116Lf/7zn7XuW0qV8ETrm2++QUlJiRYmk5eXh8GDB2PdunUYN25cjY+bM2cOZs+eXV/djAs53Zac+s+GMnVxgqYtdm2rTp/sFyUEXll2oZblyCcS9yhHsbNkp/ZTHtuWLOkJL8vLvGw5PbK6XtaUBGXbKMp9LLUP+jo5Xa/sv/aUEddExzFlsGRY7XhsNGU58rHKotyusUzHMcWxoZwmmqmJTdsB9KmK3cp9tP4aynAAvWzHtW3kPsi2WtlOkNMWR5SovBJWSZEBxzZERPHFsY0iicY206ZNw4QJEzBgwAAMGjQICxYswP79+zFp0iQAwPjx49G+fftQjsqVV16JYcOG4Y477sCoUaPw5JNP4oMPPsBDDz0EAGjZsiVatmypPUdmZiYKCgrQtWvXWvcrpa5AidbhwJhowmcAYPr06SgrKwv9bdu2LaH9JCIiIqoNjm2IiOhIMHbsWMyfPx8zZsxAnz59sGnTJqxcuTL0+7d161bs2LEj1H7o0KFYtmwZHnroIfTu3RvPPPMMVqxYgR49esS1Xw1+Bcr111+P22+/3djms88+w3HHHVdPPToUZlOfgTZERJR6LDvKkOUotkupjWMbIiJKRck2tpk6dSqmTp1a47q1a9c67hszZgzGjKn91Yrffvtt1H1q8BMo11xzDSZOnGhs07lz55i2fTgwprS0FG3btg3dX1paij59+sS0TSIiIiITjm2IiIjSU4OfQMnPz0d+fn5Ctt2pUycUFBRg9erVoUFFeXk53nvvvajS7lORqW7Y7QRgvWSiqHkogDMHwvQchmmX5TrHlMdKHomIQIElC9rUHBZDPsqhp1W26zbdcBR5KZby6hzblX2KprjQVLxXlxlAZB6JcRrj2ueaGKdHNkxF7NiWzB8x5Zw4ths5q8Qy5ZjI53Vrq04hbMo8cWtrmKrYkXniaBteZl2wgW0nZlpnThWd8ji2SQxmohBRPHBsY8CxjauUykDZunUrNm3ahK1btyIQCITSc/ft2xdqc9xxx+H5558HcOh/aq+66irceuutePHFF7F582aMHz8e7dq1w+jRoxvoVRAREREdwrENERFR6mjwK1CiMWPGDDzyyCOh5b59+wIA3njjDQwfPhwAsGXLFpSVlYXaXHfdddi/fz8uueQS7NmzByeeeCJWrlyJ7Ozseu07ERGll2SrE6bUxLENERElC45t3KXUCZSlS5di6dKlxja2LCGwLNx88824+eabE9gzIiIiouhxbENERJQ6UuoECsVOrfWTNcOmE4KO/BExiFOXTHkohx4aOV/C8nhk4/Btj+ihV2+rZlzYXsN2AEB5HkfOhghFsZXMCFM+CqBnpDjaylwTU16KpOanWLKx6L/btlQB9ya14nI2Wc8fMTc25po4MkUiPAdgzDlxb2vIS5HHbCCKvBR1W2K7tnweNY/EkI8CiCwTQ+aJbGvKPAFYG1xrNtxDpWLdLhG5Mo1tiIhUHNvUEsc2rngChYiIKAa8zJWIiIjSCcc27lIqRJaIiIiIiIiIqCHwCpQjkNs0gFGV9Bge55x+OPL5OsdjTVNdyXVK2Y7rlMfqY+U8xnKqXLV8xlDec6ituh3DtMUQ5T6yv4byH7fXpj9HxFV14nr2OJr3zTBdsrO8yrAspyI2lenIMiJDWyuKtlFNTexSlqNty45chuN4rKFk59Cy2pYlO3HBqf6IkganOCYiFcc2MeLYxhWvQCEiIiIiIiIicsErUIiIiGLAOmEiIiJKJxzbuOMVKERERERERERELngFChnrht1OFhozUWROiNdraCymR/aE27pmf5imPPaIKYTVLJOgfE6xXTXLxGPIR5FtDZknhx6LyFzyU0zb1R5meAoA5lOnhmwSN1Fl15jauuWaaFk2MpvEkFUipyI25ZyYpjiW612mPNbySIJyCmFDXorMRzHksBgzT8Tzsi44TjjVH1HSYiYK0ZGFY5s44djGFa9AISIiIiIiIiJywStQiIiIYsA6YSIiIkonHNu44xUoREREREREREQueAUKOag1hLJm2HHy0FZyF0TGicw8UTMZHBkdtiHrwxb5EvJ5TLkVMndFzZsw5aMAWv6I7cieMGSVODJPRFs1mkLmmLjkp+jbjbjKXcC9Sa1EkZdizDGR3HJNTBkoppwTQ4aI43lkW0MeieP4MGWXOLJJZF6Kst6RpSIeq+alODJa9LasDU6AoO08NuO1XSKKK9PYhohSE8c2CcCxjStegUJERERERERE5IJXoBAREcWCSfVERESUTji2ccUTKGTkNg2gDXW64YBYp7PUkhj5RKJ8RitckbUqhimPHWU58nlUpvIesS3LbXph9bWZpjiu6bGGdY7yH5UswzFtN1FMZTiSqQzHbbvRlOWYlt36IEtk1Kamch+3qYntyNMNO8qI1LIcl/7Zhj7wslYioppximOi1MSxDSUDnkAhIiKKgYUEJdXHf5NERERErji2cccTKERERLGw7eiuwopmu0RERET1jWMbVwyRJSIiIiIiIiJywStQKCqmumHbNueP2EqWiXMaY5Froi3pGRGW13DeT/TBMbWrmmsSeSuHqBkpbvkiMealSK75KSbJnoFiyjFx2240uSambRkyTgCRc2LIMXHbriPnROWYxlhvqz3WlusiT4/MuuD6Z9kJusw1ff6RhiglMBOFKDlxbFP/OLZxxytQiIiIiIiIiIhc8AoUIiKiWHCqPyIiIkonHNu44hUoREREREREREQueAUK1YmpNvFU77n6HWrxm8yTUPNGAFhqLoRYJ09galkaMkNEPFbNSDHloxzqrrJe5ot4xLlHNXslirwU2daS0RmmbTVE5okbUyaKW15KrLkmcrvByFkltqOtIedEtnVkl5i2K/qgPNbRVuanBCN/TlYFngIlD8u29e+JOG6XiBqOcWzDfBSiuGLOSXLh2MYdr0AhIiIiIiIiInLBK1CIiIhiEYScJCx+2yUiIiKqbxzbuOIJFEoYU7mBLO+p05THaimLrV9U5ZgqVyvpkWVE+mO1R1riYi2PmFrZMpTwyHIf9aFuZTim9bJcKRqm7dblErtopiaOZn28ynLcnlN9rCytMZXpuE2PHMXUxCzTISJKXizvIYoey3QonbCEh4iIKAaH64QT8ZcoX3zxBc466yy0atUKubm5OPHEE/HGG28k7PmIiIgodaTi2Ka+8QQKERHREeL000+H3+/HmjVrsGHDBvTu3Runn346SkpKGrprREREREmPJ1CIiIhiYSfwLwF2796NL7/8Etdffz169eqFY445BnPnzkVFRQU+/vjjxDwpERERpY4UG9s0BGagUIOQOQ+OumGPklUiMyxEHomlZpeY8lHEesuxzvA8HpmroT825rwUyZh5YshSSVYyu0RluJQvYbkm8nFiu9rzGnJX5Hot46SG59GmJg6KzBPWBVME5eXl2rLP54PP54t5ey1btkTXrl3x6KOPol+/fvD5fHjwwQfRunVr9O/fv67dJTriye9zZqLQkYpjG0pnPIFCREQUC9uuW+iyabsACgsLtbtnzpyJWbNmxbxZy7Lw+uuvY/To0WjatCk8Hg9at26NlStXonnz5nXpMREREaWDBI9t0gFLeIiIiGJg2Yn7A4Bt27ahrKws9Dd9+vQa+3H99dfDsizj3+effw7btnH55ZejdevW+Ne//oX169dj9OjROOOMM7Bjx4563HNERESUjBI9tkkHvAKFiIgoCeXm5iI3N9e13TXXXIOJEyca23Tu3Blr1qzBSy+9hJ9++im03fvvvx+rVq3CI488guuvvz4e3SYiIiJKWzyBQknBVCt5qvdcbdny6uttv5I34dVXOnNOwqc/HfkoejSFljmi5azUxJSXIuMx1OeVeSkmnthDT4y5Ky4ceSTRkNklxidSXp/bc8Yr18TRh8htHfshED5gbJeMFpn5Q2kiSS5zzc/PR35+vmu7iooKAIBH5Cl5PB4E3T4bRBQ149iG+SiU4phzkqaSZGyTzFjCQ0REdAQYMmQImjdvjgkTJuDDDz/EF198gWuvvRbffPMNRo0a1dDdIyIiIkp6PIFCREQUAyuYuL9EaNWqFVauXIl9+/bhlFNOwYABA/D222/jhRdeQO/evRPzpERERJQyUm1s0xBYwkNJz638QS3xkYUqjovFtFIbvbWx3Ede3i6nFFbqdBzbEc+jl+24fJuoj3X74jGU6STsorloLsdzK+cxTU1cX2U5ptIguVmlhIclOpQqBgwYgH/+858N3Q2iI55b+QNLfKihsUSHqGY8gUJERBQL1gkTERFROuHYxhVLeIiIiIiIiIiIXPAKFCIioljYSEx9XPr8Iw0RERGlEo5tXPEECqU8NX/CUTMsc0EsW7mpr7PllMLKpWaOXJOAmPNYyUSJZnpkyfk86uNin4o4qumSTVxyQRwMuSfO/JEoth1FBoqel2Kebtg2ZaCI/rI2mIiIEkX9jWEeCtUXjm2I3PEEChERUQws24aVgJreRGyTiIiIyA3HNu6YgUJERERERERE5IJXoBAREcWCSfVERESUTji2ccUTKJRW3Go31TpiR5yHyB+xvN5wW7mhaPJSJDU/ReShOJ9H2ZbMUkk1bl+chiwTZ16KYVt1yDVRsQ6YiIiSQTRjGyITjm2I6i6lSnhuu+02DB06FDk5OWjWrFmtHjNx4kRYlqX9jRw5MrEdJSKi9GcDCCbgL33+kYZqgWMbIiJKGhzbuEqpK1CqqqowZswYDBkyBH/7299q/biRI0diyZIloWWfz5eI7hER0RGEQWsUDxzbEBFRsuDYxl1KnUCZPXs2AGDp0qVRPc7n86GgoCABPaJU8//au/+gKsv8/+OvAwlHLXAL5Eey/kgHakFRVIRxNkwMq02ZbVt/tIs6pZ+ctRExFR3DJtt1q7WlVot+rAPbZLVNm7rk0CLGNiVhQkzmiKNmYukhMzkIraCc8/3DrydPIOdw4Pz0+Zi5xzk3132f680Z4OV1X/d1uzp10fHjka+YzOXoccNXtu3FLxOHtwp5QadbbbrTg9twOn+Zxw0DCAxkG/RWn2Ub+A2yDeA9fnULj6sqKys1ePBgxcfHa/HixTpz5ky37dva2tTc3Gy3AQBgx6ofF1vr083bhcEfkG0AAH2ObONQwA+gTJ8+Xf/4xz9UUVGhp556Sv/973911113qaPj6itybtiwQeHh4bYtLi7Ogz0GAAC4OrINAADe4fUBlPz8/E4Lof10q6+vd/n8s2fP1owZM5SUlKTs7GyVlpbq008/VWVl5VWPWb16tcxms207ceKEy+8PAAhQbrlC46bHB8KjyDYAAL/kY9lm8+bNGjZsmIxGo1JTU7V3795u27/99ttKSEiQ0WhUUlKSdu7cafvahQsXtGrVKiUlJWngwIGKjY1VTk6OTp482aM+eX0NlOXLl2v+/PndthkxYkSfvd+IESMUERGhI0eOaOrUqV22CQ0NZTE22OnNvabTgn/7kz1XXCH86eOQu2H4ydoqvfovVg/e19F6JN0f2oNedvM+5R3/dLkPAOBpZBv4g15lG9ZP6TXWMQG699ZbbykvL09FRUVKTU1VYWGhsrKydOjQIQ0ePLhT+z179mjOnDnasGGDfvWrX2nr1q3Kzs5WbW2tEhMT9cMPP6i2tlaPPfaYxowZo7Nnz2rp0qWaMWOG9u3b53S/vD6AEhkZqcjISI+939dff60zZ84oJibGY+8JAAhAFknuWMvZ9TFL+AiyDQDAL/lQtnn22We1cOFCLViwQJJUVFSk9957T1u2bFF+fn6n9s8995ymT5+uFStWSJLWr1+v8vJybdq0SUVFRQoPD1d5ebndMZs2bdLEiRPV0NCgn//85071y+u38PREQ0OD6urq1NDQoI6ODtXV1amurk4tLS22NgkJCXr33XclSS0tLVqxYoU++eQTffXVV6qoqNDMmTM1cuRIZWVleasMAAAASWQbAMC146eLmbe1tXXZrr29XTU1NcrMzLTtCwoKUmZmpqqqqro8pqqqyq69JGVlZV21vSSZzWYZDAYNGjTI6Rq8PgOlJwoKClRSUmJ7PXbsWEnSBx98oIyMDEnSoUOHZDabJUnBwcH6/PPPVVJSoqamJsXGxurOO+/U+vXrmcYKAOgVg9UqgxvWK3HHOeG7yDYAAF/h7mzz0wXM161bp8cff7xT+++++04dHR2Kioqy2x8VFXXVNcRMJlOX7U0mU5ftz58/r1WrVmnOnDkKCwtzthT/GkApLi5WcXFxt22sV3zg/fv31/vvv+/mXgHd68n6Hd3dU2y19mI+Xac1T67+pIbOb9yL+wm6+QXMvb8AQLaBf+rJ3/Brab0Usg3QvRMnTtgNVnhr4P/ChQv67W9/K6vVqhdffLFHx/rVAAoAAD7DXU/MYQYKAADwBjdnm7CwMKdme0RERCg4OFiNjY12+xsbGxUdHd3lMdHR0U61vzx4cvz4ce3evbtHs08kP1sDBQAAAAAABK6QkBClpKSooqLCts9isaiiokJpaWldHpOWlmbXXpLKy8vt2l8ePDl8+LB27dqlm266qcd9YwYK4EOY+gn4EWagAIBDZBvAj/hQtsnLy9O8efM0fvx4TZw4UYWFhWptbbU9lScnJ0c333yzNmzYIElaunSpbr/9dm3cuFH33HOP3nzzTe3bt08vv/yypEuDJ7/5zW9UW1ur0tJSdXR02NZHufHGGxUSEuJUvxhAAQAAAAAAPmPWrFk6ffq0CgoKZDKZlJycrLKyMttCsQ0NDQoK+vGGmvT0dG3dulVr167VmjVrNGrUKG3btk2JiYmSpG+++UY7duyQJCUnJ9u915ULtzvCAAoAAK7woas0AAAAveZj2WbJkiVasmRJl1+rrKzstO/+++/X/fd3vXD1sGHD7BZldxUDKE64/I1ubm72ck8AAM64/Pu6L/5QXpVFUi8ejtXteQE3I9sAgH8h2/gGBlCccO7cOUmdn1sNAPBt586dU3h4uLe7Afgcsg0A+CeyjXcxgOKE2NhYnThxQjfccIMMBueG5JqbmxUXF9fpWdeBJNBrDPT6JGoMBIFen+RajVarVefOnVNsbKzb+mWwWmVww1Ugd5wT+ClXso0U+L9zAr0+KfBrDPT6pMCvMdDrk8g2/owBFCcEBQVpyJAhLh3r7LOu/Vmg1xjo9UnUGAgCvT6p5zVydQa4ut5kGynwf+cEen1S4NcY6PVJgV9joNcnkW38EQMoAAC4wscWWgMAAOgVso1DQY6bAAAAAAAAXNuYgeImoaGhWrdunUJDQ73dFbcJ9BoDvT6JGgNBoNcn+XCNFqtkcMMVFUvgXKVB4PHZn8c+Euj1SYFfY6DXJwV+jYFen+TDNZJtHDJY3focJAAAAktzc7PCw8OVeUuurgvu++BzsaNNu44Wymw2B/y93wAAwPvINs5jBgoAAK7gPmEAABBIyDYOsQYKAAAAAACAA8xAAQDAJW66SqPAuUoDAAD8CdnGEWagAAAAAAAAOMAMFAAAXMF9wgAAIJCQbRxiBkof+uMf/6j09HQNGDBAgwYNcuqY+fPny2Aw2G3Tp093b0d7wZUarVarCgoKFBMTo/79+yszM1OHDx92b0dd9P333+uBBx5QWFiYBg0apAcffFAtLS3dHpORkdHpM3z44Yc91GPHNm/erGHDhsloNCo1NVV79+7ttv3bb7+thIQEGY1GJSUlaefOnR7qqWt6Ul9xcXGnz8poNHqwtz334Ycf6t5771VsbKwMBoO2bdvm8JjKykqNGzdOoaGhGjlypIqLi93eT1f1tL7KyspOn6HBYJDJZPJMh4FrDNmma2Qb7wr0bCMFdr4h29gj2/gXBlD6UHt7u+6//34tXry4R8dNnz5dp06dsm1vvPGGm3rYe67U+PTTT+v5559XUVGRqqurNXDgQGVlZen8+fNu7KlrHnjgAR04cEDl5eUqLS3Vhx9+qEWLFjk8buHChXaf4dNPP+2B3jr21ltvKS8vT+vWrVNtba3GjBmjrKwsffvtt12237Nnj+bMmaMHH3xQn332mbKzs5Wdna0vvvjCwz13Tk/rk6SwsDC7z+r48eMe7HHPtba2asyYMdq8ebNT7Y8dO6Z77rlHU6ZMUV1dnXJzc/XQQw/p/fffd3NPXdPT+i47dOiQ3ec4ePBgN/WwGxar+zbAR5Btuka28Z5AzzZS4Ocbsk3XyDb+wWC1BtB8Gh9RXFys3NxcNTU1OWw7f/58NTU1OTXy6kucrdFqtSo2NlbLly/Xo48+Kkkym82KiopScXGxZs+e7YHeOufgwYO67bbb9Omnn2r8+PGSpLKyMt199936+uuvFRsb2+VxGRkZSk5OVmFhoQd765zU1FRNmDBBmzZtkiRZLBbFxcXpkUceUX5+fqf2s2bNUmtrq0pLS237Jk2apOTkZBUVFXms387qaX09+dn0RQaDQe+++66ys7Ov2mbVqlV677337ILh7Nmz1dTUpLKyMg/00nXO1FdZWakpU6bo7NmzTl8p7mvNzc0KDw9X5tAlui4otM/Pf9HSpl3HN8lsNissLKzPzw+4gmzzI7KNdwV6tpGurXxDtiHb+BtmoPiAyspKDR48WPHx8Vq8eLHOnDnj7S71mWPHjslkMikzM9O2Lzw8XKmpqaqqqvJizzqrqqrSoEGDbAFDkjIzMxUUFKTq6upuj3399dcVERGhxMRErV69Wj/88IO7u+tQe3u7ampq7L73QUFByszMvOr3vqqqyq69JGVlZfncZyW5Vp8ktbS0aOjQoYqLi9PMmTN14MABT3TXY/zpM+yN5ORkxcTEaNq0afr444+90wmrxX0b4OfINr6BbON/fxfJN53522foKrKNf2ARWS+bPn26fv3rX2v48OE6evSo1qxZo7vuuktVVVUKDg72dvd67fK9e1FRUXb7o6KifO6+PpPJ1Gmq3HXXXacbb7yx277OnTtXQ4cOVWxsrD7//HOtWrVKhw4d0r/+9S93d7lb3333nTo6Orr83tfX13d5jMlk8ovPSnKtvvj4eG3ZskWjR4+W2WzWX/7yF6Wnp+vAgQMaMmSIJ7rtdlf7DJubm/W///1P/fv391LP+kZMTIyKioo0fvx4tbW16dVXX1VGRoaqq6s1btw4z3aGhdaALpFtfAfZxr+yjUS+6QrZxoPINg4xgOJAfn6+nnrqqW7bHDx4UAkJCS6d/8ppnklJSRo9erRuueUWVVZWaurUqS6ds6fcXaO3OVufq668jzgpKUkxMTGaOnWqjh49qltuucXl86LvpaWlKS0tzfY6PT1dt956q1566SWtX7/eiz2Ds+Lj4xUfH297nZ6erqNHj+qvf/2rXnvtNS/2DPAfZJtLyDZXR7bxL+Qb/0a28S8MoDiwfPlyzZ8/v9s2I0aM6LP3GzFihCIiInTkyBGPhQx31hgdHS1JamxsVExMjG1/Y2OjkpOTXTpnTzlbX3R0dKfFuS5evKjvv//eVoczUlNTJUlHjhzxasiIiIhQcHCwGhsb7fY3NjZetZ7o6OgetfcmV+r7qX79+mns2LE6cuSIO7roFVf7DMPCwvz+Cs3VTJw4UR999JHn39hileSGKyoBtNAafBPZ5sd+uYJsQ7ZxJ/JNZ2QbDyLbOMQAigORkZGKjIz02Pt9/fXXOnPmjN0fZHdzZ43Dhw9XdHS0KioqbKGiublZ1dXVPV7R31XO1peWlqampibV1NQoJSVFkrR7925ZLBZbcHBGXV2dJHn0M+xKSEiIUlJSVFFRYVu4ymKxqKKiQkuWLOnymLS0NFVUVCg3N9e2r7y83O6qhq9wpb6f6ujo0P79+3X33Xe7saeelZaW1unxjL76GfaVuro6r/+8Af6EbNM7ZBvvCfRsI5FvukK2gS9hEdk+1NDQoLq6OjU0NKijo0N1dXWqq6tTS0uLrU1CQoLeffddSZcWe1qxYoU++eQTffXVV6qoqNDMmTM1cuRIZWVleauMbvW0RoPBoNzcXD355JPasWOH9u/fr5ycHMXGxna7GrU33HrrrZo+fboWLlyovXv36uOPP9aSJUs0e/Zs2yr133zzjRISErR3715J0tGjR7V+/XrV1NToq6++0o4dO5STk6Nf/vKXGj16tDfLkSTl5eXplVdeUUlJiQ4ePKjFixertbVVCxYskCTl5ORo9erVtvZLly5VWVmZNm7cqPr6ej3++OPat2+f03+wPa2n9T3xxBP6z3/+oy+//FK1tbX63e9+p+PHj+uhhx7yVgkOtbS02H7OpEuLF17+GZSk1atXKycnx9b+4Ycf1pdffqmVK1eqvr5eL7zwgv75z39q2bJl3ui+Qz2tr7CwUNu3b9eRI0f0xRdfKDc3V7t379Yf/vAHz3f+8n3C7tgAH0G2uYRsQ7bxpEDPN2Qbso0/YwZKHyooKFBJSYnt9dixYyVJH3zwgTIyMiRder632WyWJAUHB+vzzz9XSUmJmpqaFBsbqzvvvFPr169XaGjfPz6qL/S0RklauXKlWltbtWjRIjU1NWny5MkqKyuT0Wj0aN+d8frrr2vJkiWaOnWqgoKCdN999+n555+3ff3ChQs6dOiQbSX6kJAQ7dq1S4WFhWptbVVcXJzuu+8+rV271lsl2Jk1a5ZOnz6tgoICmUwmJScnq6yszLYQV0NDg4KCfhxHTU9P19atW7V27VqtWbNGo0aN0rZt25SYmOitErrV0/rOnj2rhQsXymQy6Wc/+5lSUlK0Z88e3Xbbbd4qwaF9+/ZpypQpttd5eXmSpHnz5qm4uFinTp2y/UGWLl0Zfe+997Rs2TI999xzGjJkiF599VWf/Y9LT+trb2/X8uXL9c0332jAgAEaPXq0du3aZXcOAH2HbJMhiWxDtvGsQM83ZBuyjT8zWK0BNBwEAICbNTc3Kzw8XJkx/6frgkL6/PwXLe3adeolmc1mhYWF9fn5AQAArkS2cR638AAAAAAAADjALTwAALjCXff0MjEUAAB4A9nGIWagAAAAAAAAOMAMFAAAXGGxSLK46bwAAAAeRrZxiBkoAAAAAAAADjADBQAAV3CfMAAACCRkG4eYgQIAAAAAAOAAM1AAAHAFV2kAAEAgIds4xAAKAACusFgluSEQWAInZAAAAD9CtnGIW3gAAAAAAAAcYAAFuIYtWLBAa9eutb2+/fbbZTAY9Kc//cmundVqVWpqqgwGg5544glPdxPwSVarxW2bu9TW1mratGkaNGiQbrrpJi1atEgtLS1uez8A8DSyDeA6f8w2nsYACnCN6ujoUGlpqWbMmCHpUpD47LPPNHToUO3fv9+ubUlJiU6ePClJGjdunMf7CqD3Tp48qczMTI0cOVLV1dUqKyvTgQMHNH/+fG93DQD6BNkGgLsxgAIEoCFDhuiFF16w27dnzx4NGDBAx48ft73u16+fJkyYIEk6fPiwzp07p3nz5tmFjHPnzmn16tW2/2SlpKR4pgjA11mtl+7p7evt/y+01tzcbLe1tbX1qrulpaXq16+fNm/erPj4eE2YMEFFRUV65513dOTIkb74jgCA25BtAA9wc7YJBAygAAEoNTVVn376qe211WpVbm6uli1bpqFDh0qSduzYoXvvvVcGg0GSVFNTowEDBmjOnDk6dOiQ2tvbJUnr16/X+PHjFRkZqejoaMXExHi+IOAaFBcXp/DwcNu2YcOGXp2vra1NISEhCgr68U9///79JUkfffRRr84NAO5GtgHgCxhAAQLQpEmT7ELGa6+9phMnTmj16tW2fdu3b7dNcZUurY0wevRoxcfHy2g0qr6+XocPH9aLL76oZ599VrW1tUxxBa50+VF/7tgknThxQmaz2bZd+fPrijvuuEMmk0nPPPOM2tvbdfbsWeXn50uSTp061etvBwC4E9kG8AA3Z5tAwAAKEIAmTZqkgwcPqqWlRa2trVqzZo2efPJJXX/99ZKkgwcP6uTJk5o6dartmMshwmAwaPTo0dq/f7+WLVumxYsXa9SoUaqpqWGKK+BBYWFhdltoaGiX7fLz82UwGLrd6uvr9Ytf/EIlJSXauHGjBgwYoOjoaA0fPlxRUVF2s1IAwBeRbQD4guu83QEAfS8lJUVBQUGqra3Vrl27FBkZqQULFti+vmPHDk2bNk1Go9G2r7a2VnPnzpUkJScnq7CwUCdOnNAbb7yh8+fPq76+nqs0wJUsFsnghlXle7hS/fLlyx0uBDtixAhJ0ty5czV37lw1NjZq4MCBMhgMevbZZ21fBwBfRbYBPMBHso0vYwAFCEADBgxQUlKS3nnnHb3yyivauXOn3RXm7du3a9GiRbbXX375pZqammwhYuzYsXrhhRe0ZcsW3XDDDaqurtbFixe5SgP4oMjISEVGRvbomKioKEnSli1bZDQaNW3aNHd0DQD6DNkGgC9gAAUIUJMmTdLf/vY3zZw5UxkZGbb93377rfbt26cdO3bY9tXU1CgkJESJiYmSpHnz5ik7O1s33XSTpEtXcCIjIxUXF+fRGgCfZrVKcsM9vW68T3jTpk1KT0/X9ddfr/Lycq1YsUJ//vOfNWjQILe9JwD0FbIN4GZ+mG08jQEUIECNGTNG/fr10zPPPGO3/9///rcmTpyoiIgI277a2lolJiaqX79+kqR+/fp1+vrYsWM903EAbrN3716tW7dOLS0tSkhI0EsvvaTf//733u4WADiFbAPA2wxWawANBwGwmTJlisaNG6eNGzfa7Z8xY4YmT56slStXeqlngH9rbm5WeHi47hgwW9cZQvr8/Bet7dr9w5sym80KCwvr8/MDgL8i2wDuQbZxHjNQgABisVh0+vRp/f3vf9fhw4e1ffv2Tm0mT56sOXPmeKF3AAAAPUO2AeBLGEABAsiHH36oO+64QwkJCXrnnXe6HOHl6gzQR7hPGADcjmwDeBDZxiEGUIAAkpGRIYslcB4TBvg0i1UyEDIAwJ3INoAHkW0cCnLcBAAAAAAA4NrGDBQAAFxhtUpyw1XRALpKAwAA/AjZxiFmoAAAAAAAADjADBQAAFxgtVhldcN9wtYAukoDAAD8B9nGMWagAAAAAAAAOMAMFAAAXGG1yD33CfO0CQAA4AVkG4eYgQIAAAAAAOAAAygAALjAarG6bQMAAPA0X8s2mzdv1rBhw2Q0GpWamqq9e/d22/7tt99WQkKCjEajkpKStHPnTvv6rFYVFBQoJiZG/fv3V2Zmpg4fPtyjPjGAAgAAAAAAfMZbb72lvLw8rVu3TrW1tRozZoyysrL07bffdtl+z549mjNnjh588EF99tlnys7OVnZ2tr744gtbm6efflrPP/+8ioqKVF1drYEDByorK0vnz593ul8GayAtiQsAgJs1NzcrPDxcGZqp6wz9+vz8F60XVKntMpvNCgsL6/PzAwAAXMkXs01qaqomTJigTZs2SZIsFovi4uL0yCOPKD8/v1P7WbNmqbW1VaWlpbZ9kyZNUnJysoqKimS1WhUbG6vly5fr0UcflSSZzWZFRUWpuLhYs2fPdqpfzEABAMAFF3VBF61u2HTB26UBAIBrkLuzTXNzs93W1tbWZT/a29tVU1OjzMxM276goCBlZmaqqqqqy2Oqqqrs2ktSVlaWrf2xY8dkMpns2oSHhys1NfWq5+wKT+EBAKAHQkJCFB0drY9MOx03dlF0dLRCQkLcdn4AAIDLPJFtrr/+esXFxdntW7dunR5//PFObb/77jt1dHQoKirKbn9UVJTq6+u7PL/JZOqyvclksn398r6rtXEGAygAAPSA0WjUsWPH1N7e7rb3CAkJkdFodNv5AQAALvNEtrFarTIYDHb7QkND3fZ+7sIACgAAPWQ0GhngAAAAAcOXsk1ERISCg4PV2Nhot7+xsVHR0dFdHhMdHd1t+8v/NjY2KiYmxq5NcnKy031jDRQAAAAAAOATQkJClJKSooqKCts+i8WiiooKpaWldXlMWlqaXXtJKi8vt7UfPny4oqOj7do0Nzerurr6qufsCjNQAAAAAACAz8jLy9O8efM0fvx4TZw4UYWFhWptbdWCBQskSTk5Obr55pu1YcMGSdLSpUt1++23a+PGjbrnnnv05ptvat++fXr55ZclSQaDQbm5uXryySc1atQoDR8+XI899phiY2OVnZ3tdL8YQAEAAAAAAD5j1qxZOn36tAoKCmQymZScnKyysjLbIrANDQ0KCvrxhpr09HRt3bpVa9eu1Zo1azRq1Cht27ZNiYmJtjYrV65Ua2urFi1apKamJk2ePFllZWU9unXJYLVarX1XJgAAAAAAQOBhDRQAAAAAAAAHGEABAAAAAABwgAEUAAAAAAAABxhAAQAAAAAAcIABFAAAAAAAAAcYQAEAAAAAAHCAARQAAAAAAAAHGEABAAAAAABwgAEUAAAAAAAABxhAAQAAAAAAcIABFAAAAAAAAAf+H06F1VzL6qORAAAAAElFTkSuQmCC\n",
"text/plain": [
"