{ "cells": [ { "cell_type": "markdown", "id": "29896764-6d44-4af7-b98e-1259562e0097", "metadata": {}, "source": [ "# Euler's Elastica" ] }, { "cell_type": "markdown", "id": "c3eafa01-bf1d-4269-b259-0bf27836ad6e", "metadata": {}, "source": [ "The **Euler elastica problem** is a classical problem in mechanics describing the equilibrium shapes of an inextensible, slender elastic rod subject to bending. Originally studied by Euler in the 18th century, it remains a fundamental model in continuum mechanics, geometric variational calculus, and modern applications ranging from soft robotics to biological filaments.\n", "\n", "Let $\\mathcal{I} = [0,L] \\subset \\mathbb{R}$ denote a one-dimensional interval representing the material coordinates of the longitudinal axis of the rod in its undeformed, straight configuration, where $L$ is the total length of the rod. The configuration of the rod is described by the angle field $u(\\sigma)$, defined as the angle between the tangent to the rod axis and the horizontal direction, evaluated at the dimensionless arc-length parameter $\\sigma := s/L$.\n", "\n", "We consider a rod with constant cross-section, clamped at the left end and subjected to a vertical end load at the right end. Under these assumptions, the boundary-value problem associated with the dimensionless continuous Euler's elastica equation reads\n", "$$\\frac{\\mathrm{d}^2 u}{\\mathrm{d}\\sigma^2} + f \\cos u = 0 \\quad \\text{in } (0,1), \\quad u(0) = 0, \\quad \\frac{\\mathrm{d}u}{\\mathrm{d}\\sigma}(1) = 0,$$\n", "where $f = PL^2/B$ is the **dimensionless load parameter**. Here, $P < 0$ denotes the **vertical component of the load** applied at the right end, and $B$ is the **bending stiffness** of the rod.\n", "\n", "This problem also admits the following variational formulation: find $u$ minimizing the energy functional\n", "$$\\min_{u} \\;\n", "\\frac{1}{2} \\int_0^1 \\left( \\frac{\\mathrm{d}u}{\\mathrm{d}\\sigma} \\right)^2 \\, \\mathrm{d}\\sigma\n", "\\;-\\;\n", "\\int_0^1 f \\sin u \\, \\mathrm{d}\\sigma,\n", "\\qquad\n", "u(0) = 0.$$\n", "\n", "In this notebook, we aim to solve this problem using the framework of DEC. Two complementary discretization strategies are considered. First, we formulate a **variational DEC approach**, in which the continuous energy functional is discretized directly using discrete differential forms and Hodge star operators, and the equilibrium configuration is obtained by minimizing the resulting discrete energy. Second, we investigate a **non-variational DEC formulation**, where the Euler–Lagrange equation is discretized directly at the level of the governing differential equation." ] }, { "cell_type": "markdown", "id": "7e68e97b-57d7-435f-b8a5-58a0c7204466", "metadata": {}, "source": [ "## Imports" ] }, { "cell_type": "code", "execution_count": 1, "id": "eb2aba03-c5ac-4612-bd66-197604fc930a", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import dctkit as dt\n", "import jax.numpy as jnp\n", "from jax import jit, grad\n", "from dctkit.math.opt import optctrl\n", "from scipy import sparse\n", "import os\n", "from dctkit.dec import cochain as C\n", "from dctkit.mesh import util\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "code", "execution_count": 2, "id": "66aefb84-5f95-4914-8ba5-6c7bd1961d02", "metadata": {}, "outputs": [], "source": [ "dt.config()" ] }, { "cell_type": "markdown", "id": "d270545e-38e5-4d97-bee3-8090ddfe96a1", "metadata": {}, "source": [ "## Complex generation and util" ] }, { "cell_type": "code", "execution_count": 3, "id": "1895bce5-1f16-445e-9050-acfa64779bb2", "metadata": {}, "outputs": [], "source": [ "data = \"data/xy_F_35.txt\"\n", "F = -35.\n", "np.random.seed(42)\n", "\n", "# load true data\n", "data = np.genfromtxt(data)\n", "\n", "# sampling factor for true data\n", "sampling = 10\n", "\n", "num_elements = 10\n", "\n", "L = 1\n", "h = L/(num_elements)\n", "\n", "B = 7.854\n", "\n", "# define dimensionless load\n", "A = F*L**2/B" ] }, { "cell_type": "code", "execution_count": 4, "id": "70e30c3c-7cf7-4441-95e7-20bcbf82ef6b", "metadata": {}, "outputs": [], "source": [ "def compute_transform_u_xy(num_nodes):\n", " # bidiagonal matrix to transform u in (x,y)\n", " diag = [1]*(num_nodes)\n", " upper_diag = [-1]*(num_nodes-1)\n", " upper_diag[0] = 0\n", " diags = [diag, upper_diag]\n", " transform = sparse.diags(diags, [0, -1]).toarray()\n", " transform[1, 0] = -1\n", " return transform\n", "\n", "def compute_true_solution(data, sampling, num_elements):\n", " x_true = data[:, 1][::sampling]\n", " y_true = data[:, 2][::sampling]\n", " u_true = np.empty(num_elements, dtype=dt.float_dtype)\n", " for i in range(num_elements):\n", " u_true[i] = np.arctan(\n", " (y_true[i+1]-y_true[i])/(x_true[i+1]-x_true[i]))\n", "\n", " return u_true, x_true, y_true\n", "\n", "def reconstruct_xy(u, h, num_nodes):\n", " transform = compute_transform_u_xy(num_nodes)\n", " cos_u = h*jnp.cos(u)\n", " sin_u = h*jnp.sin(u)\n", " b_x = jnp.insert(cos_u, 0, 0)\n", " b_y = jnp.insert(sin_u, 0, 0)\n", " x = jnp.linalg.solve(transform, b_x)\n", " y = jnp.linalg.solve(transform, b_y)\n", " return x, y" ] }, { "cell_type": "code", "execution_count": 5, "id": "e6be0253-62b0-4ab6-9166-dd600f0a0833", "metadata": {}, "outputs": [], "source": [ "# compute true solution\n", "u_true, x_true, y_true = compute_true_solution(data, sampling, num_elements)" ] }, { "cell_type": "code", "execution_count": 6, "id": "b3365a34-ea63-41f9-9215-34fd0d5b04b6", "metadata": {}, "outputs": [], "source": [ "num_nodes = num_elements + 1\n", "mesh, _ = util.generate_line_mesh(num_nodes, 1.)\n", "S = util.build_complex_from_mesh(mesh)\n", "S.get_hodge_star()\n", "\n", "# define internal cochain to compute the energy only in the interior points\n", "int_vec = np.ones(S.num_nodes, dtype=dt.float_dtype)\n", "int_vec[0] = 0\n", "int_vec[-1] = 0\n", "int_coch = C.CochainP0(complex=S, coeffs=int_vec)\n", "# define the unit cochain on the dual nodes\n", "ones_coch = C.CochainD0(complex=S, coeffs=np.ones(num_elements,dtype=dt.float_dtype))\n", "\n", "# initial guess for the angles (EXCEPT PRESCRIBED ANGLE AT LEFT END)\n", "u_0 = np.zeros(num_elements-1, dtype=dt.float_dtype)\n", "\n", "# cochain to zero residual on elements where BC is prescribed\n", "mask = np.ones(num_elements, dtype=dt.float_dtype)\n", "mask[0] = 0." ] }, { "cell_type": "markdown", "id": "b1d6d867-1fc1-4c33-b56c-4626c87e188c", "metadata": {}, "source": [ "## Variational Formulation" ] }, { "cell_type": "markdown", "id": "b8d6dd2c-9ee5-43c9-95e4-bba1060effdb", "metadata": {}, "source": [ "The discrete elastica energy can be formulated as\n", "$$\\mathcal{E}_{\\text{el}}(u) := \\frac{1}{2} \\langle k, \\star d^\\star u \\rangle - \\langle f\\mathbb{1}, \\sin u \\rangle,$$\n", "where $\\mathbb{1}$ is the dual one $0$-cochain, $k := \\mathbb{1}_{\\text{int}} \\odot \\star d^\\star u$ is the primal $0$-cochain that encodes the (discrete) curvatures on the nodes and $\\mathbb{1}_{\\text{int}}$ is the primal $0$-cochain that is identically $1$ on the interior and $0$ on the boundary nodes." ] }, { "cell_type": "code", "execution_count": 7, "id": "c7d74298-6d9e-46ef-bcd9-093fca5d09bc", "metadata": {}, "outputs": [], "source": [ "def obj(u):\n", " \"\"\"Computes the total potential energy.\n", "\n", " Args:\n", " u: current configuration angles.\n", "\n", " Returns:\n", " value of the energy.\n", " \"\"\"\n", " # apply bc at left end\n", " u = jnp.insert(u, 0, u_true[0])\n", " u_coch = C.CochainD0(complex=S, coeffs=u)\n", "\n", " # curvature at internal nodes\n", " curvature = C.cochain_mul(int_coch, C.star(C.coboundary(u_coch)))\n", "\n", " # bending moment\n", " moment = curvature\n", "\n", " # potential of the applied load\n", " A_coch = C.scalar_mul(ones_coch, A)\n", " load = C.inner(C.sin(u_coch), A_coch)\n", "\n", " energy = 0.5*C.inner(moment, curvature) - load\n", "\n", " return energy" ] }, { "cell_type": "code", "execution_count": 8, "id": "290b4d1f-ed0e-48f5-9eb9-01275d430ba9", "metadata": {}, "outputs": [], "source": [ "prb = optctrl.OptimizationProblem(\n", " dim=num_elements-1, state_dim=num_elements-1, objfun=obj)\n", "\n", "prb.set_obj_args({})\n", "sol = prb.solve(x0=u_0)" ] }, { "cell_type": "code", "execution_count": 9, "id": "8d550ca9-57a5-438c-855f-f34f0ab62fbd", "metadata": {}, "outputs": [], "source": [ "u = jnp.insert(sol, 0, u_true[0])\n", "\n", "x, y = reconstruct_xy(u, h, num_nodes)" ] }, { "cell_type": "code", "execution_count": 11, "id": "a020f353-556e-4f83-8bc2-e91ee7d39fb8", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA/YAAAHVCAYAAABbv2DlAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjgsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvwVt1zgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAmw9JREFUeJzs3QdYFEcbB/A/vQkoKoK9g70LYu81GktsUaPGkhg1xsQWY9QUjcY0NWqsSdQYTWLX2HvDir33BlhBRZBy3/OOOT5AQPrd3v1/z7Pe3d7eMTt33u67M/OOhU6n04GIiIiIiIiINMnS0AUgIiIiIiIiorRjYE9ERERERESkYQzsiYiIiIiIiDSMgT0RERERERGRhjGwJyIiIiIiItIwBvZEREREREREGsbAnoiIiIiIiEjDGNgTERERERERaRgDeyIiIiIiIiINY2BPZGR+/fVXWFhYxC729vbw8PBA/fr1MXHiRAQHB6fr/bdu3YqqVavCyclJvf/KlSthCnbs2KH2R24NpXDhwujZs6fB/j4REWWsuMfj5BZDHnuMiRwD5VgYl9TPuHHjUvU+d+7cUa8JCAjItPOsa9euGfz4HxYWpvaT3x/KCNYZ8i5ElOEWLFgAb29vREZGqmB+z549mDRpEqZMmYKlS5eiUaNGqX5PnU6Hjh07omTJkli9erUK7r28vPjpERERJWL//v3xHn/55ZfYvn07tm3bFm996dKlWX/J1GH+/PlTHdiPHz9eBcwVK1Y02bqVwF72U9SrV8/QxSGNY2BPZKTKli2rWtb12rdvj48++gi1atVCu3btcPHiReTJkyfVB8qHDx+ibdu2aNiwYYaU8/nz56pXgVz91vKB1dHR0dDFICIiI+Pr6xvvce7cuWFpafnKeq0fVzLzWP66uiKijMGu+EQaUrBgQXz33Xd48uQJfvnll3jPHT58GK1bt4abm5s6OFeqVAnLli2LfV66eumvmI8YMUIdvON2l5MeARLsOzs7q5MRPz8/rFu3LtHua5s2bULv3r3VCY5sGxERoa40y8UIuTIvr3VwcFDvLz0PhLxX5cqV1fblypXDhg0bXtk/uVjRtWtXuLu7w87ODqVKlcLPP//8ynbnzp1Ds2bN1HvlypUL7733nqqTlJB6kH04evQoOnTogBw5cqBYsWLqufDwcIwaNQpFihSBra0t8uXLhw8++ACPHz+O9x7Si2L48OFqiISUQS62HDx4MEV/n4iITIv++Ldr1y51/JPjghwjk+uGnljX7cDAQPTv318dq+UYJMciac2NiopKUTnkbw0cOFCdH0jPPDmOSk+CP//8M8XHciG9AmvUqKF69WXLlg1NmzbFsWPHXvl78j7S609/vP7999+TLFfCOrh9+zb69euHAgUKqH3NmzevOiYHBQWpbunVqlVT2/Xq1St2uEPc93jdOY/egQMHULNmTbWN/A05xssxPCXk85H9P336tDo/kvqQupI6lgs3r3Pjxg1069Yt3jmNnMPFxMSo52UogLyfkM9Zv58c0kdpxRZ7Io1p0aIFrKys1AmEnnQLlEDXx8cHs2bNgqurqzqQd+rUSR185CDRp08fVKhQQbX2Dxo0SAXQcqARO3fuROPGjVG+fHnMmzdPrZ8xYwbeeOMNLFmyRL1PXHIi0LJlSyxcuBDPnj2DjY1N7EmJHIQl6JUTk2nTpqltb968ib///huffvqpKtsXX3yBN998E1euXFEHWnHmzBl1QqS/eCFB88aNGzF48GDcv38fY8eOVdvJQb9u3brqb0oZpdfC4sWL1YE2NaQeOnfurC4KyD7IMAUpk+QgkAN/7dq1ceLECfV35WKFLPr66tu3rzqB+eSTT1S9nTp1Sr1fSi8uEBGRabl7964K4uT4N2HCBNWqnxpy/Kxevbp63eeff64uOMtx56uvvlIBoP4i+evIMDs5J5DjrASicpzs0qULrK2tVeD8umO5lP2zzz5Tx3K5ffHiBb799lt1TJQL2PohBxLUyzZt2rRRx+yQkBAVeMvFgdftuwT1ErhLgC3nBXLu8eDBA3XMf/TokWoEkP3Vl0HKKPSNEyk559GfV0hALhdRpLxy8ULq448//kjx5yJllPMuueAycuRI7Nu3T30m169fx5o1a5J83b1799Q5jdSfDN+QMqxdu1adN1y+fFmVw9PTUzVyyL68++676jxN6IN9olTTEZFRWbBggU7+ax46dCjJbfLkyaMrVapU7GNvb29dpUqVdJGRkfG2a9Wqlc7T01MXHR2tHl+9elW997fffhtvO19fX527u7vuyZMnseuioqJ0ZcuW1eXPn18XExMTr2w9evR4pUx169ZVzx0+fDh23YMHD3RWVlY6BwcH3e3bt2PXBwQEqG2nTp0au65p06bqb4WEhMR734EDB+rs7e11Dx8+VI9HjBihs7CwUO8RV+PGjdV7bt++XZecsWPHqu0+//zzeOs3bNig1k+ePDne+qVLl6r1s2fPVo/Pnj2rHn/00Ufxtlu8eLFa/8477yT794mISLvkN97JySnR49/WrVtf2V7Wy3EnoUKFCsU7XvTv31+XLVs23fXr1+NtN2XKFPUep0+ffm3ZZDs53gYGBsY7lss5QvHixWPXJXUsv3Hjhs7a2lo3aNCgeOvl3MDDw0PXsWNH9VjOKfLmzaurXLly7PmBuHbtms7GxkbtW3J10Lt3b7XdmTNnktwXOQeS10lZE0rpOU+nTp2SrA95bzknSo58PrLdTz/9FG/9119/rdbv2bMnyc9z5MiRaht/f/94r33//ffVOcz58+fV43v37iX5HSFKLXbFJ9Kgl8fJly5duqS6pr/99tvqsXTZ0y9ylVlaEc6fP5/ke8lVen9/f3UlX7qc6UmvgO7du+PWrVuvvF7G+ydGrj5XqVIl9rF0kZMuaJL4Rt8yL6Q7mpAr3vou8NJSLmP/5Yp6wn2Q56U7nf5KfZkyZVTvg7ikB0JqJNwHfSKkhF3g3nrrLdXqIeXT/32hr289SUooLSJERGR+ZFhXgwYN0vx6ac2V2W/kWBn3GNi8efPYnnUiOjo63vP6bt160kIdN/+OHMulJVvOFeR4ntxxUFrM5T179OgR729IN3bpKafP3C7nBJKzR467ccfkFypUSLVSv86///6r9lV/LpAaqTnnkeN1UvWRGgmP9/rzDf35QGLknEJ6N0gvjLjkHEPO4RImXyTKCAzsiTRGAnHpsqYPlKVrupDuXdKNLu4yYMAA9Zx0ZU+KdHuTg4wE5Qnp/4b8vbgS21YfyCckY+cSrpd1QgJ2/fvLQVm67ifcBzlQx90H2Va66SeU2LrkJNwHeV8JzBN2gZOTFnlvfR3obxP+PXltzpw5U1UGIiIyDUkdF1NKjuXStTvhMVAuZMc9BkqgGvd5/Vh+veSOj687luvPJ6SbfMJyyLj7uMfh1/2t5Eg39dRmyU9YxpSc82TE+UJix/ak6jMueS4151VEGYHNS0QaI0no5Iq9floUSR4nZFy4jPNOTHJT2kkrg4yHk6vcCckV+bh/Qy+js+ZKGfQ9BCRZXWIkiZCQA6yMRUwosXXJSbgP8r5ycUFOOOIG93LRQ95bn8hHf4CXdZJcT09eywM1EZF5Suq4KLlZ9Enp4kp4vJDjrIw1//rrrxN9H31AKInx4uZzSXh8Tu74mDBATVhm/XtJThxpfU9K3ONgUn8rOXKMTdh7IKVSc86TEecL+mN73LpLqj7jkudSc15FlBEY2BNpiGRYlavUkihGErnoD2AlSpTA8ePHVdKb1JJu5pKAZvny5ZgyZYrKZi+ke9+iRYvUVXXJrpuZpPu9dMuTrLtyYqNv0U+MbDd58mS1v3G746cmGU5ipBVE3lf2WaYV1Pvnn39ULwn99ID6CyqSsC/usAPJxpvSzMVERGQeJGmaJGKNS7phP336NN66Vq1aYf369SppnlzsTkpyF+qFDBuTVm1993NpCJDWdnnf17WSS/Z7aaGW5G5JDbnTl0FaoyW57tChQ2MvEMjwOkkuF3foXWJkeIEk7JMu80ntjz5ZrUzDl/Bvp/ScR84XJJlgYvWRGnK8l0S+Cc83kpt3Xs4ZJk6cqGbgkWSAepJ4V+pLypbcfhKlBQN7IiMlmdb148aCg4Oxe/dulSVWWrZXrFgRr1VZruDLgVIOyjJ+S1qSZb76s2fPqoPKX3/9lezfkoOPZHeXA41cOJDAWjK2ShnkwJ0Vc9T/9NNPato4ybz7/vvvq5MhaZWQ8XTSPVE/Hm3IkCGYP3++ypIrmWn1WfFlzF16yP5L/clUgKGhoWp6HH1WfJlGR3oTCBkTKJmPf/zxR9X1r1GjRqqe5KKIi4tLhtQFERGZBjl2jBkzRmW6l3Hqkql9+vTp6gJ9XJLFfvPmzWqMugSREsDKcDXJiC8Bv2R/T0n3dWkJlrH+8jf1WfHl+JhwyrvEyHFXyjF69Gg1a41ka5eLDBIYS0Z8eT+Zlk16+Ummd8niLrlxZKYYmRZWsuKnpJu7/A0ZZ1+nTh2VFV+mwJXXS4Z4uVDg7e2tLkRIQ4Mc3+W4KzmA5IKBLCk955GM+hLYS31I/UsjgkyhKxfrU0rOhyTrv1yIkZ57+qz48vflnCUp0kAgQbycq8j+Sg8I6XEpn4ec4+gbTGSKYXlu1apV6mKADF2UzzDudMREKZbqdHtElKn02Wr1i62trcpYL1l3J0yYoAsODk70dcePH1cZa2VbyTYrGWwbNGigmzVrVuw2SWXFF7t371bbS7ZfySIrmfLXrFmT4oz9Ur4yZcq8sl4yxbZs2fKV9fI+H3zwQbx1Uj7JlpsvXz61D7lz59b5+fnpvvrqq3jbSSZdyYIv2fLd3Nx07777rm7VqlWpyoovmWgTev78ucq6L2WWvy/ZdSWD7aNHj+JtFxERofv4449VXUsZpK7279//SlZcIiIyj6z4iR3/9MeL4cOH6woUKKCOrbKtzOqS2PFCjkuDBw/WFSlSRB2D5PhWpUoV3ejRo3VPnz59bdn0x9UZM2boihUrpt5DMsDLrC2pmX1n5cqVuvr16+tcXFx0dnZ2qqwdOnTQbdmyJd52c+fO1ZUoUUKdp5QsWVI3f/58tU+vy4ovbt68qY73cq4i5ZQs+3IOExQUFLvNkiVLVPnl+YTvkZJzHrF37151jJb9kG2GDRumZrlJaVZ8+axPnDihq1evnvr85DOR84KEn0din6fMcNC1a1ddzpw5VRm9vLzU+Zc+a7+e1Ktk+ZcycnYdSg8L+SfllwGIiIiIiMjYSO86yVMjPQIo/aQ3gOQbSDhsgshYMSs+ERERERERkYYxsCciIiIiIiLSMHbFJyIiIiIiItIwttgTERERERERaRgDeyIiIiIiIiINY2BPREREREREpGHWhi6AFsTExODOnTtwdnZWU4kQEREZmsxW++TJE+TNmxeWlrxOn1481hMRkZaP9QzsU0CC+gIFCmTU50NERJRhbt68ifz587NG04nHeiIi0vKxnoF9CkhLvb5CXVxc0vWhREZGYtOmTWjSpAlsbGxg6ri/poufrWkzp89Xq/saGhqqLjrrj1GUecd6rX5HshrriXXF7xX//2lBpIZ+01NzrGdgnwL67vdyoM+IwN7R0VG9j7F/kTIC99d08bM1beb0+Wp9XzlELPOP9Vr/jmQV1hPrit8r/v/TgkgN/qan5FjPQXlEREREREREGsbAnoiIiIiIiEjDGNgTERERERERaRjH2BMRmbHo6Gg11iwhWWdtbY3w8HC1jSkz1n2VcX9WVlaGLgYRERFpAAN7IiIznRc1MDAQjx8/TvJ5Dw8PlSHc1JOzGfO+Zs+eXZXN2MpFRERExoWBPRGRGdIH9e7u7iozbMLAMSYmBk+fPkW2bNlgaWnao7aMcV/lYkNYWBiCg4PVY09PT0MXiYiIiIwYA3siIjMj3c31QX3OnDmTDHZfvHgBe3t7owl2M4ux7quDg4O6leBePit2yyciIqKkGM8ZDBERZQn9mHppqSfjpv+MEsuDQERERKTHwJ6IyExx3Lbx42dEREREKcHAnoiIiIiIiEjDNBnYz5gxA0WKFFHjIatUqYLdu3cnu/3OnTvVdrJ90aJFMWvWrCwrKxEREREREZm+6Bgd9l9+gFUBt9WtPM4qmgvsly5diiFDhmD06NE4duwYateujebNm+PGjRuJbn/16lW0aNFCbSfbf/rppxg8eDD++eefLC87ERGlT8+ePVX3dFlknvc8efKgcePGmD9/vkqCl1K//vqrmkqOiIiIKCNsOHUXtSZtQ5c5B/DhnwHqVh7L+qygucD++++/x7vvvos+ffqgVKlS+PHHH1GgQAHMnDkz0e2ldb5gwYJqO9leXte7d29MmTIFBnH4cPxbIiJKlWbNmuHu3bu4du0a/v33X9SvXx8ffvghWrVqhaioKNYmERERZSkJ3t9fdBR3Q8LjrQ8MCVfrsyK411RgL9MRHTlyBE2aNIm3Xh7v27cv0dfs37//le2bNm2Kw4cPZ32W4REjcKJLP9x8Cjxo1Ra6ESOy9u8TEZkAOzs7eHh4IF++fKhcubLqibVq1SoV5EtLvP4icLly5eDk5KQu/g4YMEDNVS927NiBXr16ISQkRLX8yzRy33zzjXpu0aJFqFq1KpydndXf6Nq1a+xc8kREREQJSXf78WvOILFO9/p18nxmd8vX1Dz29+/fV/MvS9fLuORxYGBgoq+R9YltL6068n6enp6vvCYiIkIteqGhoepWLgSk+WKAtNBPm4avOn6NYyetMaXvPNhFvYDnl/8ib57s8HC1R15Xe3jKkl3uO8DT1Q6Otpr6iF6hry9zmarJnPbXnPbV1PZX9kGn06mu60l1X5fn9bev7eLu7w9cuACULAn4+GRGkeOVK7Ey1atXDxUqVFDDrKRXlgTs0lOrcOHCakjWwIEDMWzYMPz888/w9fXFDz/8gLFjx+Ls2bOx7ylLeHg4xo8fDy8vLxXQf/zxx3jnnXewbt06GILsp5RLPrOE89ibwneRiIhI6w5efRjbUh8dFoKwC/vgXLF57PNyRiXPy3Y1iuXMtHJYm8L0P3LSk9yUQIltn9h6vYkTJ6oTu4Q2bdqUvnmflyxB5DlLuDzVITTSAhHWtrj2LAbXrjxM8iWO1jrksAVy2L28zS63dkB225e3rraAVdK7bjQ2b94Mc2JO+2tO+2oq+2ttba1ao6UFW3pCJefJkyfJPm8/bhzsf/op9nH4hx8ifNw4ZBYJZuXCrP6Ca1ySHPXMmTPqOWmR18uZMydGjhypgnT5fRe2trbqNu5vuuxrhw4dYh/nypULX3/9NRo2bIg7d+4gW7ZsyGry+Tx//hy7du16ZZhBWFhYlpeHiIiI4gt+8jKojwoJQtCyzxET/gyOJf1g5eia6HaZRVOBvZxkSYtFwtZ5aVVJ2CqvJyeviW0vJ7ZyspeYUaNGYejQobGP5SRRunJKl34XF5e0t9g3bIjGDg7YPH8+6vbphwdWDrg753fc9SiEO4/DcTc0HHflNiQct0Oe41lENMKiLBAWBdwOSzx6t7QA3J3tkDe7Azxd9K39/7X8/7fkcLQx2FzIchIugZAkt5JEV6bOnPbXnPbV1PZXWqVv3rypAlWZLSQxcgFUAl3pkp7k74e/PyzjBPVCgnzbTp0yreVe6l5+vxP7LZb1coyQ57Zv366CeGmRl99wCYplv+V56Z4v+y37JdvG3deAgAB1Yff48eN4+PBhbM+Ax48fI2/evMhqUmYHBwfUqVPnlc8qsYsbRERElLXcne3x4t41BC/7HBbWtvDoNvmVoF6/XWbSVGAvLSwybZ2cXLdt2zZ2vTxu06ZNoq+pUaMG1qxZ80rLu4yhTOrkXMZvypKQbJ/mE/oaNYBBg1R3fOH07AmyD+6FYm/USfIloeGRKtC/8/g57oQ8f3kb57EkY4iM1iEwNEItSbG3sVSBv3Tvz5tdgn0H5JPHcjHgv27/Drbxu3hmtHTVnQaZ0/6a076ayv7KkCYJai0tLdWSGH1Aq98uUZcuJbraUtbLb14m0GfET6xM586dU1OhykULSaT33nvv4auvvoKbmxv27NmjEq/Kvsfdb7nV76u0gEtiPrmIK2Ptc+fOrWZckbwscmEgyXrIRPI39TMAJPzeaf17SEREZAqqF3HDs53zYOWYHe5vjYdVthzxnpfmERl2LdtlJk0F9kJa0rt3764CcwnaZ8+erU685ARO39p++/Zt/P777+qxrJ8+fbp6Xd++fVUyvXnz5mHJkiVZX/hJk4A335QuA8DWra898XWxt4GLhw28PJwTfT4mRof7TyNw+/Fz1covAb+6L8G/uhAQrp4Pj4zBlXvP1JIUadVXwb+6AGD///sS+Gd3UFeYrKR7ABGRnoypT836TLRt2zacPHkSH330kUqOKoH4d999FxuML1u27JULxRLkJ7wwILlXJJGe9NIS8l5EREREiZG8bNIg/Mu83zBi5VlY2jnFS6Knj57GvlE602MpzQX2nTp1woMHD/DFF1+o6Y7Kli2L9evXo1ChQup5WRd3TntpvZHn5WRPkiZJV8qpU6eiffv2htmBqlWB9etf3qaTpaUF3F3s1VIpiW3CI6NVy74+0L/7X2v/bf39x8/x7EU0HoVFquX0ncS7dsoX0cPFPrbFPzboj3Pf1SFBl399Qq20Dl8gIuMm3e2HDwcmT/7/OpntI5MT6MlBVIZYSWAeFBSEDRs2qG730krfo0cPFeBLYD9t2jS88cYb2Lt3r5r6NC5Jqic5BrZu3aqy58v2MjWqBPzyOrkofOrUKXz55ZeZui9ERESkTQsWLFDnH3Ke0bluOWTPmUtlv4875Z201EtQ36zsqwnbYe6BvZBpi2RJjH6qo7jq1q2Lo0ePwhzZ21ihcC4ntSRGxpaGhkf9181fgv7/uvr/1/IvPQCCQsMRFaNT92UBHiX6Xo62VmpMvwT6+S6chOeerSj09B7CPx8EjB0LTJiQyXtLRAbpidSuXZZlxRcSyMuMJjKmPkeOHCobvlywlez10kJfsWJFNd3dpEmTVC8uGZ8uB14J+vX8/PxU8K6/WDxixAhMmDBBHUNk+jx5P5lKb8qUKWjdunWm7xMRERFpg06nw7fffqvOHfr166eG/AkJ3huX9lDZ7yVRnvR4lu73WdXrWZOBPWUcaWGXlnZZSnkm3rIucy7ee6Lv8v/qWH+5APDg2QuEvYjG5XvP1AILT6B2N/V6q8M6bL5sizf+3omGrfzgbM9xoUQmRYL5LAjohQTeiV3ATUh6ackSlwzjimvmzJlqkTH2+kR0Xbp0UUtiM6kQERGReYuJiVHT50oDwpgxY1TC3bg9liWIz8wp7ZLDwJ5eS3XDd7VXCxA/GUTcLv+qlV8y+q/firtLV+KOS24cKlAGV9zyY2vRath6+ClsA7agXsncaFneE41K5YGTHb+CRERERERk/E6cOKHyt8mwvYEDB8KYMKqiDOvyXzR3NrXApzDw4Z9q/QsHB8ybtwRPZ6/Evw074crTaGw6E6QWO2tLNPB2V0G+3Dra8utIRERERETG5fnz5ypJngz3u3TpUmySXWOS9XP3kPkk1PovE2ReR2CIX15sHd0U/35YGx/UL4bCOR0RERWDf08FYuAfx1D5y834YPFR/HvyLp6/iJ+pmoiIiIiIyBAePHiA+vXrY/To0eqxMQb1gk2klPkJtcS4cWr8iYzjl+WTJl4qA/+6k3ex7sRd3HgY9vL+ybsqCV/DUnnQspwn6nnlVr0BiIiIiIiIstLNmzfRtGlTNR2udME3ZgzsKXNb7itXfjm9XwIS5JfN56qW4U29cOp2KNaevKOC/FuPnmPN8TtqcbK1QqPSL4P8OiUZ5BMRERERUeY7c+aMCuqtrKywZ88elJTZf4wYA3syOAnyy+V3VcvIZt44fisE6068DPJl+r1VAXfU4mxnjcYS5Jf3RK0SuWBnzZZ8IiIiIiLKeDL1bfbs2bFx40bkzZsXxo6BPRldkF+xQHa1jGpeCgG3HmPt8btYf/IuAkPDsfzYbbU421ujSWkPtCrviZrFc8HWmukiiIiIiIgofR4/fqwCegnsw8LC1H0tYGBPRsvS0gKVC+ZQy2ctS+HojUdYe+JlkB/8JAL/HL2lFlcHGzQtIy35eeFXLCdsrBjkExERERFR6ixatAiDBg3Cvn37UKpUKdja2kIrGNiTZoL8qoXd1PJ5q9I4fF2C/DtYfzIQ959GYNnhW2rJ7miDZmU8VHf9GkVzwppBPhERERERvcb333+Pjz/+GL1790aJEiWgNWzaJE0G+dWLuOGLNmXh/2lDLOnri26+BZHTyRaPwyLx56Gb6D7vIKpP2IpPV5zEvkv3ERUdY+hiE5GRGDdunJqHNr1+/fVXzXTPIyIiosTpdDqMGDFCBfUjR47E3LlzYW2tvfZvBvakaVaWFqhRLCe+erOcCvL/6OODLtULIoejDR4+e4E//G+g61x/+E7cis9WnsT+yw8QHaMzdLGJKB2Cg4PRv39/FCxYEHZ2dvDw8FBZa/fv359p9Vq4cGH8+OOP8dZ16tQJF/RTelKGefToEbp37w5XV1e1yH0Z75iUyMhIdUJWrlw5ODk5qQRHPXr0wJ07d/ipEBHRa929e1ddrJcW+4kTJ6qcX1qkvUsRREmQbvd+xXOp5cs2ZbD/ygOVWX/Daemu/wKLDtxQS25nO7QoK93186JqoRyqBwARaUf79u1VMPfbb7+haNGiCAoKwtatW/Hw4cMsLYeDg4NaKGN17doVt27dwoYNG9Tjfv36qeB+zZo1iW4viY2OHj2KMWPGoEKFCurCwJAhQ9C6dWscPnyYHw8RESXq+fPnCA8PVxeEz58/r/leeGyxJ5MN8muXyI1v2pfHodGN8Fvv6uhYNT9c7K1x70kEftt/HR1/2Y8a32zFuNWncfjaQ8SwJZ/I6EnLrcwlO2nSJNSvXx+FChVC9erVMWrUKLRs2VJtc+PGDbRp0wbZsmWDi4sLOnbsqIL/pDRo0EC9Pq4333wTPXv2VPfr1auH69ev46OPPlJX8fVX8hPrij9z5kwUK1ZMJdvx8vLCwoUL4z0vr5Uufm3btoWjo6Maw7d69eoMqx+tO3v2rAropY5q1Kihljlz5mDt2rXqpCsx0qq/efNm9TlLnfv6+mLatGk4cuSI+i4QEREl9PTpUzRv3lxdOBZaD+oFA3syeZIlv27J3JjcoQIOf9YYC3pWQ/vK+dWUeUGhEfh13zV0mLUfNSdtw5drz6js+zLWhsicyHc+7EVUvOX5i+hX1mXGkpr/bxKsy7Jy5UpEREQkuh8SlEvr/c6dO1XAd/nyZdVtPq2WL1+O/Pnz44svvlDd9WRJzIoVK/Dhhx+qMXqnTp1SwwV69eqF7du3x9tu/PjxKgg9ceIEWrRogbfffjvLexsYKxlOIYG6j49P7DoJ1GWdZChOqZCQEHURxRRO1IiIKGPdvn0bn376Kc6dO/fKhX0tY1d8Misy3319b3e1RESVxZ6L91V3/U1ngnA3JBzz9lxVS77sDmhRzgOtyudF+fyumh1rQ5RSzyOjUfrzjQapsDNfNIWjbcoOR5LMRlrK+/bti1mzZqFy5cqoW7cuOnfujPLly2PLli0qYL569SoKFCigXiOt5mXKlMGhQ4dQrVq1VJfPzc0NVlZWcHZ2VuP5kzJlyhTVyj9gwAD1eOjQoThw4IBaL70L9GSbLl26qPsTJkxQrcsHDx5Es2bNYO4CAwPh7u7+ynpZJ8+lhHSrlORH0qVfemwkRS4Mxb04FBoaqm5lmIcscekfJ1xP8bGeUo51xbrKaPxOpYz0/pIefjKMSy7+S34WY/5tT03ZGNiT2bKztkLDUnnUEh4ZjV0X7mHdybvYciYItx8/x5zdV9WSP4eDmj6vVbm8KJvP5f9Bvr8/IImzSpYE4rQuEVHmj7GXg/Lu3btVC6903Z48ebLqvi3BmQT0+qBelC5dWrXcSjfvtAT2KSXvL+PB46pZsyZ++umneOvkAoSeJHuTCwaSENDUZyKQngrJkQsvIrELqdITIyUXWOUESC7yxMTEYMaMGcluKwmSEivTpk2b1DCJxMhJIL0e6ynlWFesq4zG71TyZPibHCPkGHDz5k21GDO5AJFSDOyJANjbWKFJGQ+1SJC/4/zLIH/r2SDcevQcv+y8opaCbo4qyG+5fRnKTB6L2NPM4cOBSZNYl6RZDjZWquVcTw56T0KfwNnFGZaWlpn+t1PL3t4ejRs3Vsvnn3+OPn36YOzYsaqVPLWBoexfwuEAab16n/BvJPZ3bWxsXnmN1LcpGzhwoAq4XzfzgPS2SCwfwr1795AnT55kXy+fmQxxkN4a27ZtS7a1Xkj3S/m+6OkvCjVp0uSV18p7y8myfN8Sfn7EekoLfqdYVxmN36nkSVJWGVYn4+pluJYM79LCb7q+N1lKMLAnSiTIb1bWQy0yxnj7+WDVXX/ruSDceBiGmTsuY6ZFFRTp+wt6H1qJzic2wWbyZKBdO6ByZdYnaZIEl3G7w0ugGWVrpdZldmCfEaRVXsbdy60kTJMr8PpW+zNnzqiDeKlSpRJ9ba5cueIFk9HR0WqMfNzu85IMT9YnR95fEvvJVGt6cuKQ1N81J1LHsryOJMuTz0qGJkhSROHv76/W+fn5vTaov3jxosppkDNnztf+LZkqUZaE5CQvqRO95J4j1lNa8DvFuspo/E69asmSJWoY3Lp169CoUaPY/CtaqKvUlM/4z9aIDMjB1gotynni57cr4+iYxpjetRKau7yAXWQErrrlw5imH6BZ7+nYXLw6dOc5nzVRZnvw4IHKYr9o0aLYsfR//fWX6oovmfDlgC1d3SUhnUyBJgGiBNoyDr9q1aqJvqcE8NL9Wg74kkhHxsgnnDddWpN37dqlEu7cv38/0fcZNmyYGv8vY/8lwJT5cCXx3ieffJIpdWGK5CKI5BqQHAqSn0AWud+qVSuV8V7P29tbJSsUUVFR6NChg5rabvHixeoCjIzHl+XFixcG3BsiIjK0adOmqXMCSaIr5wKmjIE9UQpJy6Uk05vZMC+OTO+G8ZtnwS0sBJdzFkDf9p+jS3AenLqd8u4yRJR6khFfMqb/8MMPqFOnDsqWLavmL5fgb/r06arngbTc58iRQz0vgb7Mdb906dIk37N3796qm7hczZeDfpEiReK11gvJiH/t2jU1lV3u3LkTfR/Jxi/j6b/99luVrO+XX37BggUL1HR5lHISnEsyI+kSL4tcqEk4baAkP5JWfH33ShkzKbcVK1aEp6dn7JKaTPpERGQ6dDqdOj8YPHiwGnYlF96NvXU+vdgVnyi1fHyQbcggvDN5Mtqe2oaZvm9hXo32OHA/Em1nHUDVXJao+Pg5CuU27R8PIkOQrtOS8EaWpBQsWBCrVq1KNpGbLHpyoP/uu+/UfOlJDTuQKdeOHz8eb51cCNDPda/3/vvvqyUpiU3tl7B3gLmTWQikR0Zy4taj9KbgFKVERJQw6dzatWtVjz7pUWcOGNgTpYUkymvXDi4XLmBEyZLo5lUe3208j+XHbuPwfUs0/mkvetcsggH1i8HFngE+EREREVFmCw8PV0PmJFGeDOdKLJeKqWJXfKK0kinuundXtzLv/fedKmLFe74o7hKDF1ExmLXzMup9uwO/7buGyGjTznhNRERERGRIISEhKk9LixYtVL4VcwrqBQN7ogwk89wPLB2DX7pVQrHcTnj47AXGrj6NJj/swsbTgewuSkRERESUwQIDA1WeHBk2N3PmTFhZpX4qXa1jYE+UwWTK6gZeubFxSB189WZZ5HSyxdX7z9B/4RF0+uUAAm5yPC0RERERUUa4fPkyatasiXv37mH37t3qvjliYE+USaytLNHNtxB2DKuHgfWLw87aEgevPcSbP+/F4CXHcPNhGOueiIiIiCgdLly4AAcHBzUTisyWY64Y2BNlMmd7G3zS1EsF+O0r51ct+quP30HD73ZiwvqzCAmL5GdABhETw9wPxo6fERERUeJOnTqljpPNmzdHQEAAChUqZNZVxaz4RFnE09UB33WsgF41C6uAft/lB5i96wqWHb6JwQ1KqNZ9W2tea6PMZ2trq6Z1u3PnjpqTXR7L/O9xyYHyxYsXKrtsUlPAmQpj3FeZvk3KJN0KpUzyGREREdFLf//9N95++21Mnz4dffv2hbU1w1rWAFEWK5vPFYv7+GDH+XsqwL8Y/BRfrD2D3/dfw4hm3mhW1uOVIIsoI0mgWKRIEdy9e1cF90kFls+fP1dd20z9+2jM++ro6IiCBQsazQUHIiIiQ5s1axYGDBiAzp0745133jF0cYwGA3siA5Dgob63O2qXyIVlh2/h+80XcO1BGN5ffBRVCuXA6JalULlgDn42lGmkBVgCxqioKDUlTEKRkZHYtWsX6tSpAxsbG5P+JIx1XyWjr7RAGNvFBiIiIkNdiP/iiy8wbtw4DB48GD/88AMvfMfBwJ7IwAn2uvoUROuKeVW3/Nm7LuPI9UdoN2MfWpb3xIim3iiY05GfEWUKCRglkE0smJWgUoJ+e3t7owp2M4M57SsREZGWA/vTp09jwoQJGDlyJC98J8DAnsgIZLOzxtDGJdG1ekF8v/k8/jpyC+tO3MWm04F4p0ZhDGxQHNkdOcaWiIiIiMxLREQEzp49i4oVK+LPP/9kK30SOGiPyIh4uNpjcocKWDeotuqmHxmtw9w9V1H32x2Yu/sKIqJe7TJNRERERGSKQkND0aJFCzRp0gTPnj1jUJ8MBvZERqh0XhcsfNcHv/WuDq88zgh5Homv1p1F4+93qZZ86YpERERERGSqgoKCUL9+fRw5ckRlwXdycjJ0kYwau+ITGbG6JXOjVvFc+PvITXy36QJuPAzDB38cRaWC2fFZy1KoUsjN0EUkIiIiIspQV69ejW2llwS35cuXZw2/BlvsiYyclaUFOlUriO2f1MOQRiXgYGOFYzceo/3M/Riw+Aiu3X9m6CISEREREWXouPrcuXNj7969DOpTiIE9kUY42VljSKOS2DmsHjpXKwBLC2D9yUA0/mEnxq85jUfPXhi6iEREREREaebv74+nT5/C29tbBfVFihRhbaYQA3sijXF3scc37cvj3w/roJ5XbpVgb8Hea6jz7XY1XV54JBPsEREREZG2rFixAnXr1sWkSZNip+WllGNgT6RRXh7O+LVXdSx8tzq8PZzxJDwKE9afQ6Pvd2L18TtMsEdEREREmjB37lx06NABrVu3xmeffWbo4miSpgL7R48eoXv37nB1dVWL3H/8+HGyr1m+fDmaNm2KXLlyqas+AQEBWVZeoqxQu0RurBtcG992KI88Lna49eg5Bi85hjdn7MPBqw/5IRARERGR0ZowYQL69u2L999/H0uWLIGdnZ2hi6RJmgrsu3btqgLzDRs2qEXuS3CfHMmkWLNmTXzzzTdZVk4iQyTYe6tqAez4pD4+blwSjrZWOH7zMTr+sh/9Fx7GlXtP+aEQERERkdGRaZzHjx+PadOmwcrKytDF0SzNTHd39uxZFcwfOHAAPj4+at2cOXNQo0YNnD9/Hl5eXom+Th/4X7t2LUvLS2QIDrZWGNSwBDpVL4Aft1zEnwdvYOPpIGw9G4xuvoUwuGEJuDnZ8sMhIiIiIoN58eIFtm3bhmbNmmH06NH8JMwpsN+/f7/qfq8P6oWvr69at2/fviQD+7ROryCLXmhoqLqNjIxUS3roX5/e99EK7q9h5LC3wvhW3uhWPT8mb7yAHRfu49d91/D3kVt4v24RvONbEHY26bsiys/WtJnT56vVfdVaeYmIiIRkvW/fvr2an/7SpUvIly8fK8acAvvAwEC4u7u/sl7WyXMZaeLEiao7SEKbNm2Co6NjhvyNzZs3w5xwfw2nbU6gTGkLrLpuiVvPovDtpouYu+MCWhWMQeVcOjVtXnrwszVt5vT5am1fw8LCDF0EIiKiVLl37x5atmyJc+fOYf369QzqTSmwHzduXKJBdFyHDh1KcsoDGZOR0VMhjBo1CkOHDo3XYl+gQAE0adIELi4u6W5hkZPHxo0bw8bGBqaO+2scWgAYHKPD6hN38d3miwgMjcDCS1YIuPQEI/6ZAp/bZ15uOGQI8Jr/j3r8bE2bOX2+Wt1XfW8yIiIiLbh586Y61kpC9J07d6JSpUqGLpJJMXhgP3DgQHTu3DnZbQoXLowTJ04gKCgo0as+efLkydAySSbGxLIxyglfRp30ZeR7aQH31zi8Va0Q3qiYH/P2XMXMrRdwMsoZ3dqPR9dj/2Ls1l9gN3Ei0KYNEGfIy+vwszVt5vT5am1ftVRWIiKibNmywdvbG1OmTEHx4sVZIaYW2Ms0dLK8jiTJCwkJwcGDB1G9enW1zt/fX63z8/PLgpISmQZ7Gyt8UL84Ol07gB+X7MPiSs3xR6XmOJ2nGGasnIh8Fy6kKrAnIiIiIkqK5EPz9PREkSJFsHLlSlaUuU93V6pUKZU1UeY4lMz4ssj9Vq1axUucJ1eBVqxYEfv44cOHalq8M2dedjWWDPryOKPH5RNpTa7SJfDV5plY8Nc4uD5/guN5S6JVzx+xJ3thQxeNiIiIiEzA2rVr0bBhQ3z55ZeGLorJ00xgLxYvXoxy5cqpse6ylC9fHgsXLoy3jQTu0oqvt3r1ajV+Q5I0COn2L49nzZqV5eUnMirSKj98OOpdPYq1vw1B2cBLeOToih77QvHz9kuIidEZuoREREREpFG//vor3nzzTTRv3hwzZswwdHFMnsG74qeGm5sbFi1alOw2kkwvrp49e6qFiBIxaRLQrh0KXLiAv4uVwNg7jlh6+Ca+3Xgex248xncdK8DVgeN4iYiIiCjlvv/+e3z88ceqh/XMmTNhZZW+aZbJxFrsiSiTWu67d4e9ny8mdSiPb9qVg621JbacDULr6Xtw9i4zbxMRERFRyklyvDFjxuCXX35hUJ9FGNgTUTydqxfE3+/VQL7sDrj+IAxtZ+zFimO3WEtERERElOz0sfPnz1c9qFu3bo0vvvgiw6clp6QxsCeiV5TPnx1rB9VCnZK5ER4Zg4+WHsfnq07hRVQMa4uIiIiI4nn27JkaT//ee+/h+PHjrB0DYGBPRInK4WSLBT2rYXDDEurx7/uvo9Ps/bgb8pw1RkRERETKgwcP0KhRI+zcuRPr1q1DxYoVWTMGwMCeiJJkZWmBoY1LYn7PqnCxt1YJ9VpN3YN9l++z1oiIiIjM3L1791C7dm1cunQJ27dvR+PGjQ1dJLPFwJ6IXquBdx6sHVQbpTxd8ODZC3Sb64/Zu68iwSQURERERGRGZNYymad+z549qFatmqGLY9YY2BNRihTM6YgVA/zQvnJ+yBT33266iPkXLPEkPIo1SERERGRG/P39Vdd7mcZu2rRp8PLyMnSRzB4DeyJKMXsbK0x5qzy+blsWNlYWOPHQEu1nHcCFoCesRSIiIiIzsGHDBjRo0ACTJ082dFEoDgb2RJQqMm3J2z6FsKRPdWS31eHqgzC0mb4XqwJusyaJiIiITNjixYvxxhtvqO73f//9t6GLQ3EwsCeiNKmQ3xXDykfDr5gbnkdG48M/AzBu9WlOiUdERERkgubMmYNu3bqhe/fuWL58ORwcHAxdJIqDgT0RpVk2G2B+jyr4oH4x9fjXfdfQZc4BBIWGs1aJiIiITIifnx/GjRuHefPmwdra2tDFoQQY2BNRuqfEG9bUG3N6VIWznTWOXH+EllP34MCVB6xZIiIiIg2LiorCN998g7CwMJQpUwZjx45VwzLJ+DCwJ6IM0bh0HqwZVAveHs64/zQCb8/1x5xdV6DjnHhEREREmvP8+XO0b98eY8aMwb59+wxdHHoNBvZElGEK53LCigE10bZSPkTH6PD1+rP44I+jeBrBKfGIiIiItOLRo0do0qQJtmzZgtWrV6NRo0aGLhK9BgN7IspQDrZW+L5jBXzRpoyaEm/9yUC0mb4Hl4I5JR4RERGRFlrq69atizNnzmDr1q1o3ry5oYtEKcDAnogynIy96lGjMJb2rwEPF3tcvvcMrafvxdoTd1jbREREREZMst337t0be/bsga+vr6GLQynEwJ6IMk3lgjmwdnAt1CiaE2EvojHwj2P4cu0ZREbHsNaJiIiIjMihQ4cwf/58dX/IkCEoVaqUoYtEqcDAnogyVa5sdlj4bnW8V/fllHjz9lzF23P8EfyEU+IRERERGYPNmzejfv36aio7yYRP2sPAnogynbWVJUY298asblWQzc4aB689VFPiHbr2kLVPREREZEB//vknWrZsiTp16mDTpk2co16jGNgTUZZpVtYDqwfWRMk82XDvSQS6zD6A+Xuucko8IiIiIgP466+/0LVrV3Tu3BmrVq2Ck5MTPweNYmBPRFmqaO5sakq81hXyIipGhy/WnsGgJcfwjFPiEREREWU4mYJ4/+UHWBVwW93KY72GDRti8uTJ+PXXX2FjY8Pa1zBrQxeAiMyPk501fupcEZUKZsfX685i7Ym7OB/4BDO7VUFx92yGLh4RERGRSdhw6i7GrzmDuyH/z23k4WyDQlfX4NvPh6FgwYL45JNPDFpGyhhssScig02J16tmEfzZzxfuzna4GPxUzXf/78m7/ESIiIiIMiCof3/R0XhBvS7qBU7+Pg7L5s/A7H82sY5NCAN7IjKoqoXd1JR4PkXc8OxFNN5ffBQT1p9FFKfEIyIiIkoT6W4vLfX/73QPxEQ8Q9Cyz/H8ylG4t/sMW8OLxOuWT9rGwJ6IDM7d2R6L+/igX52i6vHsXVfQbZ6/SrBHRERERKlz8OrD+C31uhgELf0ckcFX4d7pKzgUr66el+3INDCwJyKjmRLv0xalMOPtynCytcKBKw/RatpuHLnOAw4RERFRagQ/+X9QLywsLOHq1wl53p4E+/ylktyOtIuBPREZlRblPLFqYC2VRC8oNAKdfjmA3/Zd45R4RERERKnoDSleBF3Gox0L1HmUY/HqsM1dONHtSPsY2BOR0ZGgfuUHNdGynKeaEm/s6tMYsjQAYS+iAH9/YOHCl7dERERE9IrqRdzgeP8cAv8YifAbJ6B78Tze8xYAPF3t1XZkGhjYE5FRymZnjeldK+GzlqVgZWmBVQF30HbMclxt3g7o0QPw9QVGjDB0MYkogz169Ajdu3eHq6urWuT+48ePU/z6/v37q1k3fvzxR342RGS2Viz/Bxd//xR2eb3h0XkCLO0c4wX1YuwbpdU5FpkGBvZEZLTk5LxP7aL4o48PcttZ4LyFE1q/8wP2Fyj3coPJk9lyT2RiunbtioCAAGzYsEEtcl+C+5RYuXIl/P39kTdv3kwvJxGRsdq+fTs6duyItzq0x7LlK5A3d454z3u42mNmt8poVtbTYGWkjGedCe9JRJShfIrmxLpCD/HB1js4VKAMencYi1//GgufW6eBCxcAHx/WOJEJOHv2rArmDxw4AJ///l/PmTMHNWrUwPnz5+Hl5ZXka2/fvo2BAwdi48aNaNmyZRaWmojIuNSuXVv9dvbq1QuWlpZoUaGgyn4vifJkTL10v2dLvelhiz0RaYJ7mRJYuPQz1LlyBM9t7dHrrXE4nK8UULKkoYtGRBlk//79qvu9PqgXvr6+at2+ffuSfF1MTIxq1R82bBjKlCnDz4OIzE50dDTmzZuHPXv2wNraGu+++64K6oUE8TWK5USbivnULYN608QWeyLSBh8f2H/8EWZ//zX6tB+DPYUr4Z1uE/G7R0lUMXTZiChDBAYGwt3d/ZX1sk6eS8qkSZPUiezgwYNT/LciIiLUohcaGqpuIyMj1RKX/nHC9RQf6ynlWFesq4wkv2U9evTAunXrVI+lWrVqZej7m5pIDf2mp6aMDOyJSDsmTYJ9u3aYc+4Cet+zwf77QM/5B7Gwjw8qFshu6NIRURLGjRuH8ePHJ1s/hw4dis2tkZBM05TYenHkyBH89NNPOHr0aJLbJGbixImJlmnTpk1wdPx/kqm4Nm/enOL3N2esJ9YVv1dZ5/nz55gwYQLOnTuH4cOHw8PDA+vXr8/CEmjXZg38poeFhaV4Wwb2RKQtPj5w8PHBvBdR6LXgEPyvPkT3ef5Y3McH5fMzuCcyRjL2vXPnzsluU7hwYZw4cQJBQUGvPHfv3j3kyZMn0dft3r0bwcHBKFiwYLwuqR9//LHKjH/t2rVEXzdq1CgMHTo0Xot9gQIF0KRJE7i4uLzSYiIngI0bN4aNjc1r99dcsZ5YV/xeZb22bdvixo0bWLt2rQry+TtlWr9V+t5kKcHAnog0ydHWGvN7VkPPBQdx6NojdJvrjz/6+qJsPldDF42IEsiVK5daXkeS5IWEhODgwYOoXr26WidZ7mWdn59foq+RsfWNGjWKt65p06ZqvSSOSoqdnZ1aEpKTvKRO9JJ7jlhPacHvFOsqvaT3keQZKV26tGqp53cq5bRQV6kpH5PnEZFmOdlZY0Gv6qhSKAdCw6Pw9lx/nL4TYuhiEVEalSpVCs2aNUPfvn1VZnxZ5H6rVq3iZcT39vbGihUr1P2cOXOibNmy8RY5EZLuqMll0Sci0qrjx4+jQ4cOqpt2uXLlUKFCBUMXiYwAA3si0rRsdtb4tVc1NcY+5Hmkark/F5jybktEZFwWL16sTlSlS7ws5cuXx8KFC+NtI1PfSSs+EZG52bVrF+rUqYMrV67g2bNnhi4OGRF2xScizXO2t8Hv71ZH97n+OH4rBG/P8ceSfr4omcfZ0EUjolRyc3PDokWLkt1GkuklJ6lx9UREWrZq1Sp06tQJNWvWVL2WEuYDIfOmqRb7R48eqTFzMp+tLHL/8ePHySZGGDFihLry7+TkhLx586qpIO7cuZOl5SaizOciwX1vH5TN54IHz16g65wDuBT8hFVPREREmnf27Fm0a9cOrVu3VmPpGdSTpgP7rl27IiAgABs2bFCL3JfgPiky7kSmvxkzZoy6Xb58OS5cuKD+QxCR6XF1tMGid31Q2tMF95++QJc5/rh876mhi0VERESU7hwk0mK/ZMmSRBN/Ellr6SqVBPOSSMfHx0etmzNnjsqgK2PtEkuQI636CecnnDZtmsq0K9NCxJ0ah4hMQ3ZHWzX1XZc5B3Au8Am6zD6Apf1roEguJ0MXjYiIiCjFJNu9TMspQX3//v1VIlEizbfY79+/XwXq+qBe+Pr6qnX79u1L8ftIsh0LCwtkz875rolMVQ6nl8G9Vx5nBD+JUMH99QdMMENERETa8OLFC3Tr1g1Tp059bV4RIk212AcGBsLd3f2V9bJOnkuJ8PBwjBw5UnXpT25cSkREhFr0QkNDY8fsy5Ie+ten9320gvtruoz9s3Wxs8RvPSuj2/zDuHTvGTrPPoDF71ZFgRyOJrm/Gc2c9ler+6q18hIRUco8ffoU7du3x44dO7Bs2TI1tR2R0Qf248aNw/jx45Pd5tChQ+pWWtoTkitYia1P7ASoc+fOqkvLjBkzkt124sSJiZZp06ZNcHRMW1CQUMIhAqaO+2u6jP2zfacgMP2pFe6GhKPDz7sxsHQ0ctqb7v5mNHPaX63tq+SRISIi0/PJJ5+o3sr//vsvGjRoYOjikEYYPLAfOHCgCriTU7hwYZw4cQJBQUGvPHfv3j3kyZPntUF9x44dcfXqVWzbtu21WSRHjRqlxrPEbbEvUKCAmk83vRkopSxy8ti4cWPY2NjA1HF/TZeWPtsGDSPQbd4hXH0QhvnXnLG4d1Xkze5gsvubEcxpf7W6r/reZEREZBr0DZZff/013nvvPVSsWNHQRSINMXhgnytXLrW8jiTJk/HxBw8eVMnvhL+/v1rn5+f32qD+4sWL2L59O3LmzPnavyWZJhPLNiknfBl10peR76UF3F/TpYXPNp+bDZb0q4FOs/fj+oMwdF9wBEv7+8LTNXXBvVb2NyOZ0/5qbV+1VFYiIkre6dOn0adPH9X1XhoUUxKzEGkyeZ5kg2zWrBn69u2rMuPLIvclO2TcjPje3t5YsWKFuh8VFaXGpBw+fBiLFy9GdHS0Go8viySkICLz4eFqjyV9fVHQzRE3HoaphHpBoeGGLhYRERGZOUkEXrt2bTx//hzW1gZvdyWN0kxgLyQ4L1eunOoSL0v58uWxcOHCeNvI1HfSii9u3bqF1atXq1vpyuLp6Rm7pCaTPhGZBul+v6SfL/LncMC1By+D+2AG90RERGQga9euRaNGjVSMs3PnThWnEKWFpi4Jubm5YdGiRcluE3c6CBmbz+khiCiufBLc9/VVWfKv3H+m5rv/s18N5HZ+dfgNERERUWaRXGGdOnVC06ZNsWTJEtjbpyO7L5k9TbXYExFlhAJujiq493S1x+V7z9B1zgHcf/r/KS6JiIiIMpPM1JU7d241pd1ff/3FoJ7SjYE9EZmlgjlfBvd5XOxwMfgpus31x8NnzL1BREREmRvQDxs2DAMGDFA9i6tVq8Zx9ZQhGNgTkdkqnMtJBffuznY4F/gEb8/1xyMG90RERJQJZLauXr164bvvvkOZMmXU1HZEGYWBPRGZtaK5s+GPvr7Ilc0OZ++Gots8f4SERRq6WERERGRCnj17hjfffFONpf/jjz8waNAgQxeJTAwDeyIye8Xds2FJXx/kdLLF6Tuh6D7fHyHPGdwTERFRxpg6darKer9u3Tp07tyZ1UoZjoE9ERGAEnmcVcu9m5MtTtwKQY/5BxEazuCeiIiI0i4qKkrdfvLJJzh48CAaN27M6qRMwcCeiOg/Xh7OWNzHB9kdbXD85mO8M/8gnjC4JyIiojQ4e/YsypYti/3798PGxgalS5dmPVKmYWBPRBRHKU8XFdy7Otjg2I3H6LXgEJ5GvLzaTkRERJQS/v7+qFWrlgroCxYsyEqjTMfAnogogTJ5XVVw72JvjcPXH6H3gkMIe8HgnoiIiF5vw4YNaNCgAUqVKoVdu3YhX758rDbKdAzsiYgSUTafKxa+6wNnO2scvPYQ/RYdw4toVhURERElLSIiAv3790fDhg2xadMm5MiRg9VFWYKBPRFREioUyI7f362ObHbW8L/6CHPOWyI8ktE9ERERverFixews7PDjh07sHz5cjg6OrKaKMswsCciSkalgjnwW+9qcLK1woUQS7y3OIDBPREREcXS6XQYNWoUGjVqhMjISBQpUgTW1tasIcpSDOyJiF6jSiE3zO1RGbaWOuy9/AD9Fx5BRBRb7omIiMydTGfXp08ffPPNN3jzzTdVsjwiQ2BgT0SUAlUL5UB/72g42Fhi54V7eH/RUQb3REREZiwsLAzt2rXD77//joULF2Lo0KGGLhKZMQb2REQpVNwVmN2tMuxtLLHtXDA+WHwML6JiWH9ERERmaM2aNdi6dStWr16Nbt26Gbo4ZOYY2BMRpYJvUemWXw121pbYcjYIg5YcRWQ0g3siIiJzaqkXnTp1wrlz59C8eXNDF4mIgT0RUWrVKpELc3pUha21JTaeDsLgJccQuf8AsHAh4O/PCiUiIjJRFy5cQJkyZbBo0SL1uECBAoYuEpHCFnsiojSoUzI3fulWBbZWlvj3VCA++mYFot7pCfj6AiNGsE6JiIhMzOHDh1GzZk04ODigbt26hi4OUTwM7ImI0qi+tztmVssGm+hIrC1VB581/QA6eWLyZLbcExERmZDNmzejfv36KF68OHbv3s2WejI6DOyJiNKh4ZNrmL5qEixjovFnhab4o0Kzl09cuMB6JSIiMpF56sePH4/atWtjy5YtyJkzp6GLRPQKBvZEROlRsiSaXjyAYbt+Vw/HNe6PI3m91XoiIiLStsePH8PCwkJlwF+1ahWcnJwMXSSiRDGwJyJKDx8fYPhwvOf/D1qc24NIKxu83/1rBJeqwHolIiLScCv9mDFjUKFCBRXc58iRAzY2NoYuFlGSGNgTEaXXpEmwOHAA375ZGiWdrRAMO7y/+CjnuCciItKg6OhovPfee/jqq68wcOBAZM+e3dBFInotBvZERBnBxwdOPbvjl/614WxvjSPXH+GLtadZt0RERBoSHh6Ot956C/PmzcOCBQswbNgwQxeJKEUY2BMRZaAiuZzwU+eKsLAAFh24gWWHbrJ+iYiINCIgIADbtm3DihUr0LNnT0MXhyjFGNgTEWWwBt558FGjl8nzPlt5CgE3H7OOiYiIjNiDBw9UF3xfX19cu3YNb7zxhqGLRJQqDOyJiDLBwPrF0bh0HryIjsF7C4/g3pMI1jMREZERunz5MqpXr46xY8eqxxxTT1rEwJ6IKDN+XC0t8H3HCiia2wmBoeH44I+jiIyOYV0TEREZkWPHjsHPzw/W1tbo27evoYtDlGYM7ImIMomzvQ1md6+KbHbWOHj1ISasP8u6JiIiMhLbt29H3bp1UahQIezZs0fdEmkVA3siokxU3D0bvuv4ck77BXuvYfnRW6xvIiIiI7B48WI1pl6S5eXOndvQxSFKFwb2RESZrGkZDwxqUFzdH7X8JE7dDmGdExERGcitWy8vss+cORNr165FtmzZ+FmQ5jGwJyLKAkMalUR9r9yIiIpB/4VH8PDZC9Y7ERFRFtLpdBg/fjy8vLxw5coV2NjYwNbWlp8BmQQG9kREWcDK0gI/dq6EwjkdcfvxcwxachRRTKZHRESUJWQqu4EDB2LcuHH49NNPUaRIEdY8mRQG9kREWcTVwQa/dK8KR1sr7L30AJM3nmfdExERZbKIiAh07doVs2bNwuzZszF69GhYWFiw3smkMLAnIspCXh7O+LbDy2R6s3ddwZrjd1j/REREmSgwMBD+/v74+++/OaUdmSxrQxeAiMjctCzviZO3i2HWzssY/vcJlTm/lKeLoYtFRERkUoKDg2FnZ6emsTt//ry6T2Sq2GJPRGQAw5p6oXaJXHgeGa2S6T0OYzI9IiKijHL16lXUrFkT/fv3V48Z1JOpY2BPRGSgZHpTO1dC/hwOuPEwDB/+GYDoGB0/CyIionQ6ceIE/Pz81P0JEyawPsksMLAnIjKQHE62+KV7FdjbWGLnhXv4fjOT6REREaXHrl27UKdOHXh6emLPnj0oWrQoK5TMAgN7IiIDKpPXFZPal1f3f95+GRtO3eXnQURE9BrSy23/5QdYFXBb3ep7vR05cgRVqlTBjh07kCdPHtYjmQ1NBfaPHj1C9+7d4erqqha5//jx42RfI3NVent7w8nJCTly5ECjRo1UVkwiImPRpmI+vFvr5Xy6Hy87jotBTwxdJCIiIqMlF8FrTdqGLnMOqKFsclv5o7lq/UcffYSNGzfCxYVJacm8aCqwl/knAwICsGHDBrXIfQnuk1OyZElMnz4dJ0+eVN1xChcujCZNmuDevXtZVm4iotcZ1dwbvkXd8OxFNPotPILQ8EhWGhERUQISvL+/6CjuhoSrxzqdDiH7l+HE1H7o9c1i9by1NSf+IvOjmcD+7NmzKpifO3cuatSooZY5c+Zg7dq1avqK5C4GSCu9jK8pU6YMvv/+e4SGhqqkGkRExsLayhI/d62MvK72uHr/GT76MwAxTKZHREQUS7rbj19zBvpUszpdDB5tnYPHu36Ha82usMvnrZ5nMloyR5oJ7Pfv36+63/v4+MSu8/X1Vev27duXovd48eIFZs+erV5ToUKFTCwtEVHq5cxmh1+6V4WttSW2ngvG1G0XWY1ERET/OXj14f9b6qMjcX/td3hyZA3cmgxA9lpdAQsL9bxsR2RuNNNPJTAwEO7u7q+sl3XyXHKkVb9z584ICwtTGTI3b96MXLlyJbl9RESEWvSkhV9ERkaqJT30r0/v+2gF99d08bPNHN55HPHFG6UwcsVp/LjlIrzdndCw1Ku/fZnNnD5fre6r1spLRJRewU9eBvVCFxWJ6MdByNVmBJy8ayW5HZG5MHhgL8ntxo8fn+w2hw4dUrcWFhavPCfjahJbH1f9+vXVePz79++r7vsdO3ZUCfQSu1AgJk6cmGiZNm3aBEdHR2QEubhgTri/poufbcZzAFDbwxK7Ay0xZOkxDC0XjTyy0gDM6fPV2r7KxWoiInPi7myP6LAQxESEwSaHJ/J0mwwLC8tEtyMyNwYP7AcOHKha05MjCe9kTHxQUNArz0kSvNdNZSEZ8YsXL64W6b5fokQJzJs3D6NGjUp0e1k/dOjQeC32BQoUUEn30pthU1pY5OSxcePGsLGxganj/poufraZq3F0DHosOIzD1x9j6W1X/N3fB9nssu4n25w+X63uq743mamRGXAGDx6M1atXq8etW7fGtGnTkD179tfm4hkxYgR27tyJmJgYlVdn2bJlKFiwYBaVnIgym4fVU9xfMhIxNo7w6D7llaBemvo8XO1RvYgbPwwyOwYP7KVLfHLd4vUkWV5ISAgOHjyI6tWrq3XS6i7r/Pz8UvU3pZU/blf7hOzs7NSSkJzwZdRJX0a+lxZwf00XP9vMqldgRrcqeGPaHly+90x1zZ/5dhVYWibfQynjy2E+v1Va21ctlTU1JOntrVu3VMJc0a9fPzUDzpo1a5J8zeXLl1GrVi28++67qsed5NKRQN/enq12RKbi9OnTaNWqFVztLGHdYigsLSxik+gJ/dFx7BulYZXFx0oiY6CZ5HmlSpVCs2bN0LdvXxw4cEAtcl/+g3t5ecVuJ3PWr1ixQt1/9uwZPv30U7Xt9evXcfToUfTp00edMLz11lsG3BsioteTroQzu1WBrZUlNp4Owsydl1ltZNLSOgPO6NGj0aJFC0yePBmVKlVSM+G0bNkyySF3RKQt586dQ4MGDZAzZ04cO+SPuYPfUC3zccnjmd0qo1lZT4OVk8iQNBPYi8WLF6NcuXKqS7ws5cuXx8KFC+NtIwd+acUXVlZW6oegffv2aj57uQggXfd3796tuugRERm7ygVzYHybl79XUzadx47zwYYuEpFRzYAj3e7XrVunjvNNmzZVwby8fuXKlfykiEzE8+fPUblyZezatUslwpbgfc+IBljS1xc/da6obuUxg3oyZwbvip8abm5uWLRo0Wu72etJF7zly5dnQcmIiDJPl+oFceJWCJYcvIHBS45hzaBaKJTTiVVOJictM+AEBwfj6dOn+Oabb/DVV19h0qRJqtW/Xbt22L59O+rWrZvuGXC0OnNCVmM9sa4ymgTycqFOeuIMGzYMtra28f4fVi0oua9e5r+KiY5CTDTMFv//mWZdpaaMmgrsiYjM1bjWpXEuMBTHbjxGv9+PYPkAPzhlYTI9ImOdAUda7EWbNm3w0UcfqfsVK1ZULfyzZs1KMrBPyww4Wps5wVBYT6yrjCBDa3/77Td88sknKofGli1bMuR9TR3//5lWXaVmBhyeFRIRaYCdtRVmdauCVtP24HzQEwz/5wSmd6n02uk+iYxBZs6AIwl4ra2tUbp06Vdy8+zZsyfJv5eaGXC0OnNCVmM9sa4yglyskxxZEtSPHDkSn332mQrq+f+P///M8bcqNBUz4DCwJyLSiDwu9pjxdmV0mX0A607cRYX8ruhXp5ihi0Vk0BlwpGtutWrVXkmud+HCBRQqVChDZ8DR2swJhsJ6Yl2lVVRUFPr376+C+qlTp2LQoEGxXZH5veL/v4xmo4Hf9NSUT1PJ84iIzF21wm5qKh/xzb/nsOfifUMXicigM+AIGXu7dOlSlUH/0qVLmD59upoeb8CAAfx0iDTE0vJlaLJkyRIV1BNRKv7/sLKIiLSlm28hvFUlP2J0wMAlR3HzYcrHXxGZ2gw4om3btmo8vUx3J6+V6fL++ecfNS6XiIzfw4cPsXfvXhXY//rrr68dukNEr2JXfCIijZFx9V++WVaNtZds+f0XHsE/7/vBwdbK0EUjyvIZcPR69+6tFiLSllu3bqmpKp89e6aG0MjwGiJKPbbYExFpkL3Ny2R6OZ1sceZuKEYtP5FosENERGSszp49q/JnSFC/ceNGBvVE6cDAnohIo/Jmd8D0rpVhZWmBlQF3sOCPnYB0Wfb3N3TRyATI1HHSJfbo0aPx5nsnIsoIMsWlDJdxdXVV3fDj5tEgotRjYE9EpGE1iuXEpy1KqftfHw/FgdGTAV9fYMQIQxeNNK5evXq4efOmmu+9UqVKauy6jHudMGEC1q5da+jiEZHGZc+eHXXq1MGuXbuQL18+QxeHSPMY2BMRaVxv6yC8eXo7oi2tMLDNCDxwcAEmT2bLPaVLmzZtMGbMGPz11184c+aMmoLu448/VvPJb926lbVLRGkiM1bI3NwlSpRQs1vkyJGDNUmUARjYExFpnMXFi5i4YTpK3ruO+0458Hnj914+ceGCoYtGJjKvtCSzk6zzkrlaEtT98MMPhi4WEWnQjz/+iNatW6uZK4goYzGwJyLSupIl4RAVge/WfQ+rmGisK1UH671qqvVE6dWlSxfs2bNHzcbw999/q275Fy9eZMUSUYpJctdRo0ap3B3Dhw9Xt0SUsTjdHRGR1vn4AMOHo9zkyXj/wF+Y7tcZY978GD5lKiKnoctGmidzxp84cSL2sSTT69u3L3bs2GHQchGRdoJ6+c2YN28evvvuOwwdOtTQRSIySalqsZckOkREZIQmTQIOHMCgd5vAy8UKD2CLz1edNnSpyARky5YNly9fjn1cuXJl1SWfiCglpLdPqVKl8PvvvzOoJzKWFntvb2/1H3LkyJFwcnLKvFIREVHq+fjAzscH390OQZuf92LdybtoceIuWpb3ZG1Smv3yyy9488030bx5c3VyLvNOFyxYkDVKRMl6/PgxNmzYoGbTkMSbRGRELfabN2/Gpk2bVBbLBQsWZF6piIgozcrmc8UH9Yqp+2NWncL9p5yDnNJOprk7fPgwqlSpguvXr6N48eJYtmwZq5SIknTnzh01ld3AgQPZw4fIGAN7Pz8/+Pv745tvvsHnn3+uEuhwjB0RkfEZ2KAEvD2c8fDZC3y+6pShi0Ma1q1bN0RERKBTp06oVq0acufODUdHR0MXi4iM1IULF1TM8OjRI+zevRtubm6GLhKRWUhTVvwePXqo/7RvvPEGWrZsibZt2+LSpUsZXzoiIkoTW2tLTHmrAqwtLbD+ZCDWnrjDmqQ0kcR5Li4uai77YcOGqa61Q4YMYW0S0StOnTqFmjVrqot/+/btU8N3iMjIp7uTDJdNmjRBv379sHr1apQtW1aNn3ny5EnGlpCIiNLcJX9A/eLqviTSY5d8SgsbGxt1zP/111/x6aefqjH30gpHRJRQgQIF0KZNG/UbIfeJyEgD+1mzZuHdd99F+fLl4erqikaNGmHv3r344IMPMGPGDAQEBKB06dJqLB4RERnewPrFUcrTRXXJH7PylArQiFKjf//+qgu+zGEvPfTEs2fPWIlEFOuff/5Rs2dIfDB37lzkzMnJVomMOrD/+uuvERoainfeeUeNrQ8JCcHBgwcxdepU9O7dG1u3bsX777+Pnj17Zl6JiYgolV3yy6su+f+eki75d1l7lKyEU9lJz7wtW7aoLvkyI44MvfPx8WEtEpEyffp0vPXWW5g9ezZrhEgr092lZB57adEfM2ZMespEREQZqExeV3xQvzh+2npRJdLzLZoTuZ3tWMeUqFy5ciF//vyoUKFCvEVmxBGSFf+3335j7RGZOekBNnbsWHz55ZdqOuyJEycaukhEZi3NY+yT4u7ujm3btmX02xIRUTp88F+X/EdhkeyST8mSJHmTJ09WQ+sOHTqEAQMGqARYzs7ObKknolgSzEtQP2nSJEyZMgWWlhkeVhBRZrXYp4SFhQXq1q2b0W9LREQZ0CW/zfS92HA6EGtO3EXrCnlZp/QKb29vtXTu3Dm2VU4y4Q8aNAgNGzZkjRGRIkm0Je9Wr169WCNERoCX1oiIzKhL/sAG+iz5pxD8JNzQRSINkAv2zZs3x6JFi3DnDqdNJDJnkl9LevTExMSo3wUG9UTGg4E9EZGZdckv7emCx2GR+GwFs+TTq+SEPTG+vr4qcS4RmafAwEDUq1dPjaW/cuWKoYtDRJndFZ+IiIyXjZV0ya+A1tP3YNOZIKw+fgdtKuYzdLHIiGTLlg1ly5ZFxYoVVdI8ufXy8lKz4Dx9+tTQxSMiA5Cp7KTrfXh4OHbt2qWSaBKRcWFgT0RkZkrndcGgBiXww5YLGLv6NGoUywl3Z3tDF4uMxPLly3H8+HG1/Pzzz7h48aJqxZcu+ZIoi4jMy/Xr11GzZk01R/3evXtRuHBhQxeJiBLBwJ6IyAwNqF8Mm84E4vSdUIxecQqzu1dRgRtRs2bN1KInLXTSWpczZ054eHiwgojMjEx/2a9fP5VAM3fu3IYuDhElgWPsiYjMuEu+jZUFNp8JwqoAJkWjxNnb26NMmTIM6onMsPfOli1bYGVlhS+++IJBPZGRY2BPRGSmZF576ZIvpEt+cCiz5BMREfDLL7+gQ4cOWLp0KauDSCMY2BMRmbH36xVD2XwuCHkeiU+ZJZ+IyKzpdDrVOv/ee+9h4MCBKsAnIm1gYE9EZMbidsnfcjYIKwNuG7pIRERkIJIgc+zYsfjqq6/w008/wdKSoQKRVvB/KxGRmfP2cMHg/7rkj1t9hl3yiYjMVKdOnTB//nyMHj2aCVWJNIaBPRER4b14XfJPqu6YRERk+p48eYKPPvoIT58+hZeXF3r16mXoIhFRGjCwJyKiBF3yg7Hq+F3WChGRiQsODkb9+vVVK/358+cNXRwiSgcG9kREFNsl/8OGL7vkf7nuHEJesGKIiEzV1atXUatWLdy6dQs7d+5ElSpVDF0kIkoHBvZERBTrvbrFUC6fK0LDo7D0iiW75BMRmaDHjx+roD4mJgb79u1DxYoVDV0kIkonBvZERBTLOk6X/NOPLLEygF3yiYhMTfbs2VUG/L1796Jo0aKGLg4RZQAG9kREFI+XhzMG1S+m7n+1/hyCQsNZQ0REJmD16tWYPn26ut+7d2/kyZPH0EUiInMM7B89eoTu3bvD1dVVLXJfuhKlVP/+/dXUHT/++GOmlpOISOv61iqMAk461SV/1HJmySci0jpJkNe2bVvs2rWLw6yITJCmAvuuXbsiICAAGzZsUIvcl+A+JVauXAl/f3/kzZs308tJRGQKXfK7FY9WXfK3nQvGP0dvG7pIRESUBjJ96TfffIN3330X/fr1w5IlSzhHPZEJ0kxgf/bsWRXMz507FzVq1FDLnDlzsHbt2tdOz3H79m0MHDgQixcvho2NTZaVmYhIyzwcgQ8bFFf3x685jcAQdsknItKaWbNmYdSoURg3bhxmzJgBKysrQxeJiDKBNTRi//79qvu9j49P7DpfX1+1TrJ5enl5Jfo6yfYprfrDhg1DmTJlUvS3IiIi1KIXGhqqbiMjI9WSHvrXp/d9tIL7a7r42ZrH59vDJy82nQ3CiVuhGPH3cczpXsnkWnq0+l3WWnmJKPNEx+hw8OpDBD8Jh7uzPaoXcYOV5cvf6i5duqhkeXJLRKZLM4F9YGAg3N3dX1kv6+S5pEyaNAnW1tYYPHhwiv/WxIkTMX78+FfWb9q0CY6OjsgImzdvhjnh/pouframbfvWrWiZEzh92wo7L97H2N82wNddB1Okte9yWFiYoYtAREZgw6m7GL/mDO7G6VXlbq+Dy7HfMefHSShSpAiDeiIzYPDAXroFJRZEx3Xo0CF1m1grkYwbSqr16MiRI/jpp59w9OjRVLUwSXeloUOHxmuxL1CgAJo0aQIXFxekt4VFTh4bN25sFsMCuL+mi58tzOrzjXS/iimbL2LNLTu839YPnq72MBVa/S7re5MRkXkH9e8vOoq4l1ujw0Jw/PfxiHxwE3/t6IzhRYoYsIREZDaBvYx979y5c7LbFC5cGCdOnEBQUNArz927dy/JqTp2796N4OBgFCxYMHZddHQ0Pv74Y5UZ/9q1a4m+zs7OTi0JyQlfRp30ZeR7aQH313TxszWPz/e9esWx+dw9HL/5GGNWn8WvvaqZXJd8rX2XtVRWIsqc7vfSUh83qI8KDUbQ0s8RE/4EHl0mYkWgKz6O0cV2yyci02XwwD5XrlxqeR1JlhcSEoKDBw+ievXqap1kuZd1fn5+ib5GxtY3atQo3rqmTZuq9b169cqgPSAiMo8s+d+9VR4tpu7Bzgv38NfhW+hYrYChi0VEZLZkTH3c7ve66CgE/fkZdDHR8Hh7Mmzc8qnnZbsaxXIatKxElPk0kxW/VKlSaNasGfr27YsDBw6oRe63atUqXuI8b29vrFixQt3PmTMnypYtG2+RFg4PD48kk+0REVHiirs7Y2jjkur+l2vP4M7j56wqIiIDkUR5cVlYWcOtyQB4dPtWBfVJbUdEpkkzgb2Q6erKlSunxrrLUr58eSxcuDDeNjL1nbTiExFRxutbuygqFsiOJxFRGLn8pMpzQkREWU+y34uwy4fwYOPP6vfYoXBFWGdzS3Q7IjJtBu+Knxpubm5YtGhRstu87iQzqXH1RET0ejJOc8pbFdBi6m7sunAPyw7fRKdq/89jQkREWUOmtLO+vAv3/pkCh+LVgegowPr/uTdkVL2H68up74jI9GmqxZ6IiAyvuHs2fPxfl/yv1p7Fne37AOk95e9v6KIREZmN77+bgst/T0a2co3g/uYoWCQI6sXYN0ozcR6RmWBgT0REqdandlFUKvhfl/yZW6Dr0QPw9QVGjGBtEhFlsuXLl2P48OEYPXo0li1aAM8cTvGel5b6md0qo1lZT34WRGZCU13xiYjIeLrkf1sSaHHlBXYVrYJ/yjZAh1PbgMmTgXbtAB8fQxeRiMhktWnTBmvWrFFJpEWTMp4q+70kypMx9dL9nlPcEZkXttgTEVGaFL97BUP2/KHuf1O3F57YOrx84sIF1igRUQYLCwtD+/btsW3bNlhZWcUG9UKCeJnSrk3FfOqWQT2R+WFgT0REaVOyJN49vBJFHt7G/Ww5MM2vc+x6IiLKOA8fPkSjRo2wceNGREVFsWqJ6BUM7ImIKG18fGD38VCM2TpHPVxQtTWujBjHbvhERBno1q1bqF27Ni5cuKBa62XKZyKihBjYExFR2k2ahAZ/TEf9bC8QaWWDL72aszaJiDKITOPcuXNnPH36FHv37kX16tVZt0SUKAb2RESUPj4+GNO/MWysLLD9/D1sOxfEGiUiyoCg3sLCAnPnzsW+ffvg5eXFOiWiJDGwJyKidCuaOxt61yyi7n+59iwioqJZq0REaSRj6WVM/bNnz+Dt7Y18+fKxLokoWQzsiYgoQwxsUBy5stnh6v1nWLD3GmuViCgN/vjjD5Xx3tHRUbXYExGlBAN7IiLKEM72NhjZ3Fvdn7b1IoJDw1mzRESp8NNPP+Htt99Gt27dsGLFChXcExGlBAN7IiLKMO0q5UOFAtnx7EU0vtlwjjVLRJRCBw8exJAhQzB8+HDMnz8f1tbWrDsiSjEG9kRElGEsLS0wvnUZdX/50ds4euMRa5dS5dGjR+jevTtcXV3VIvcfP36c7GskY/jAgQORP39+ODg4oFSpUpg5cyZrnjQhJiZG3UrG+/3792PSpEnsgk9EqcbAnoiIMlTFAtnRoUp+dX/86tOIidGxhinFunbtioCAAGzYsEEtcl+C++R89NFHattFixbh7Nmz6vGgQYOwatUq1jwZtefPn6Ndu3aYPn26euzr62voIhGRRjGwJyKiDDe8mRey2Vnj+K0Q/H30FmuYUkSCcgnQZXqvGjVqqGXOnDlYu3Ytzp8/n+TrpJXznXfeQb169VC4cGH069cPFSpUwOHDh1nzZLSkp0nLli2xadMmFC1a1NDFISKNY2BPREQZzt3ZHoMbFlf3J284h9DwSNYyvZYE6NL93sfHJ3adtGDKOpnHOym1atXC6tWrcfv2bTX39/bt23HhwgU0bdqUtU5G6c6dOxg9ejROnz6NrVu3okWLFoYuEhFpHLNyEBFRpujpVwR/HryJK/efqSz5o1uWZk1TsgIDA+Hu7v7KelknzyVl6tSp6Nu3rxpjLwnHLC0tVau/BPxJiYiIUIteaGiouo2MjFRLXPrHCddTfKynlJMEedJiv3nzZpQvX57frWTwe5UyrCfTrKvUlJGBPRERZQpba0uMeaM0ei04pOa171StIIq7Z2Ntm6Fx48Zh/PjxyW5z6NAhdZvYvN3SCp/cfN4S2B84cEC12hcqVAi7du3CgAED4OnpiUaNGiX6mokTJyZaJukWndQUYxKE0euxnpIWHR0NKysrtG7dGo0bN8atW7fUQvxeZRT+/zOtugoLC0vxtgzsiYgo09T3ckdDb3dsPReML9aewW+9qjHbsxmSjPWdO3dOdhsZG3/ixAkEBQW98ty9e/eQJ0+eJJOPffrpp2rObxmvLKQFVJLuTZkyJcnAftSoURg6dGi8FvsCBQqgSZMmcHFxeaXFRE4AJRCzsbFJ0T6bI9ZT8qTL/Ycffoj169eri078TvF7xf9/hhGpod90fW+ylGBgT0REmeqzVqWx6+I97LpwD1vPBqNR6cQDNDJduXLlUsvrSLK8kJAQNZ+3TP0l/P391To/P79EX6PvOi/d7+OSVlH9NGKJsbOzU0tCcpKX1Ilecs8R6yk5y5YtQ7du3dCwYUN4eHjEfo/4nUo51hXryRy/UzapKB+T5xERUaYqkssJ79Z6mfH5y3VnEBEVzRqnRMn8882aNVPj5aVrvSxyv1WrVvDy8ordztvbW7XQC2ldr1u3LoYNG4YdO3bg6tWr+PXXX/H777+jbdu2rGkyuJ9//ln1WOnYsaMaLuLk5GToIhGRCWJgT0REmW5gg+Jwd7bD9QdhmLfnKmuckrR48WKUK1dOdYmXRbrVL1y4MN42MvWdtOLr/fnnn6hWrRrefvttlC5dGt988w2+/vprvPfee6xpMqgbN27gk08+wZAhQ9TFJmNvHSQi7WJXfCIiynQyp/3I5t4Yuuw4pm+7hPaV8yOPiz1rnl7h5uaGRYsWJVszkkwvLunavGDBAtYmGVWSPPmeFixYUOWOKF68OPOLEFGmYos9ERFliTcr5kOlgtkR9iIa3/x7jrVORCYpPDwcnTp1wgcffKAelyhRgkE9EWU6BvZERJQlLC0tMO6NMpBZy1Ycu40j1x+y5onIpMgQkebNm2PdunWxszQQEWUFBvZERJRlKhTIjreq5Ff3x60+g5iY+F2qiYi0KjAwEPXq1cOxY8ewadMmNVc9EVFWYWBPRERZalhTbzjbWePk7RD8deQma5+ITML06dMRFBSE3bt3o3bt2oYuDhGZGQb2RESUpXI72+HDRiXU/ckbziPkeSQ/ASLSrLCwMHU7btw4HD58WM3qQESU1RjYExFRlutRozCK5XbCg2cvMHXrRX4CRKRJO3bsQJEiReDv7w9ra2vkzZvX0EUiIjPFwJ6IiLKcrbUlPn+jjLr/275ruBT8hJ8CEWnK8uXL0bRpU5QvXx6lS5c2dHGIyMwxsCciIoOoWzI3GpXKg6gYHcavOfPK3ORERMZq9uzZeOutt9C2bVusXbsWzs7Ohi4SEZk5BvZERGQwY1qVgq2VJXZfvI/NZ4L4SRCR0Xvy5Am++OILDBgwAH/88Qfs7OwMXSQiIgb2RERkOIVyOqFP7SLq/lfrziI8MpofBxEZpZiYGBXUS+v80aNHMXXqVFhaso2MiIwDf42IiMigPqhfHHlc7HDjYRjm7bnKT4OIjE5ERAS6du2Kli1bqgDf3d0dFhYWhi4WEVEsBvZERGRQTnbWGNW8lLr/8/ZLCAwJ5ydCREZDWulbtWqFlStXYsiQIWylJyKjxMCeiIgMrk3FvKhSKAfCXkTjm3/PGro4RETKvXv30KBBAxw8eBAbNmxAu3btWDNEZJQY2BMRkcFJl9Zxb5SB9GxdGXAHh689NHSRiIhUMH/r1i3s3LkT9erVY40QkdFiYE9EREahXH5XdKpaQN0ft+Y0omM4/R0RZT75rdl/+QFWBdxWt/JYWupF9+7dcfbsWVSsWJEfBREZNWtDF4CIiEjvk6ZeWHfyLk7dDsWywzfRpXpBVg4RZZoNp+5i/JozuBsnt4fTw4u48ec4zJk9C507d0b27Nn5CRCR0dNUi/2jR4/UlVNXV1e1yP3Hjx8n+5qePXuqLp5xF19f3ywrMxERpVyubHYY0qikuv/txvMIeR7J6iOiTAvq3190NF5QH3bRH2cXjEBUjkKwKVyZNU9EmqGpwF6mGQkICFDjnWSR+xLcv06zZs1w9+7d2GX9+vVZUl4iIkq9HjUKobh7Njx89gI/brnAKiSiDCfd7aWlPu6An6cnNuHeiq/hULQq3N8ahynbb3JIEBFphmYCexnfJMH83LlzUaNGDbXMmTMHa9euxfnz55N9rZ2dHTw8PGIXNze3LCs3ERGljo2VJca+UVrd/33/dVwMesIqJKIMdfDqw3gt9bqYaDw9sRnZKjRFrjYjYGFtq56X7YiItEAzgf3+/ftV93sfH5/YddKlXtbt27cv2dfu2LED7u7uKFmyJPr27Yvg4OAsKDEREaVV7RK50bh0nv+3qumYSI+IMk7wk5dBvU4Xg6gn92FhaQX3jl/CrckAdT/hdkRExk4zyfMCAwNVcJ6QrJPnktK8eXO89dZbKFSoEK5evYoxY8ao+UiPHDmiWvITExERoRa90NBQdRsZGamW9NC/Pr3voxXcX9PFz9a0GcPnO7JpCey8cA97Lt3HvyfuoHHpV48BprKvaaG18hIZE3dne+iiI3F//Y+IuHUGefvMgqWtfaLbERFpgcED+3HjxmH8+PHJbnPo0CF1K4nvEpJWnMTW63Xq1Cn2ftmyZVG1alUV5K9btw7t2rVL9DUTJ05MtEybNm2Co6MjMsLmzZthTri/poufrWkz9OdbN48lNt+2xJjlxxB2JRo2lqa7r6kVFhZm6CIQaVbp3LYIWfU1wi4HIFerj2FpE7+xR84sPVztUb0Ih28SkTYYPLAfOHCgmkokOYULF8aJEycQFBT0ynMyz2iePHlS/Pc8PT1VYH/x4sUktxk1ahSGDh0ar8W+QIECaNKkCVxcXJDeFhY5eWzcuDFsbGxg6ri/poufrWkzls+3bkQUmk7di6DQCNwOccYAq9tA8eJA1aomt6+ppe9NRkSpc//+fbRs2RIvbp9FnrfGwaFwxXhJ9PTNRZLrw8oy6cYjIiJjYvDAPleuXGp5HUmWFxISgoMHD6J69epqnb+/v1rn5+eX4r/34MED3Lx5UwX4SZEu+ol105cTvow66cvI99IC7q/p4mdr2gz9+Wa3scGnLUrhwz8DMOvEI3Sc+xE8nzwAhg8HJk0yqX1NLS2VlciYSNJlmSVp966duGeX95V57KWlXoL6ZmWTPlckIjI2Bg/sU6pUqVJq2jpJfvfLL7+odf369UOrVq3g5eUVu523t7fqSt+2bVs8ffpUdfVv3769CuSvXbuGTz/9VF1IkOeJiMj4tQ6/iYW3TuNw/jKYWK8Xpq6ZAkyeDMhwqjgJVYmIkiO5lgoWLIiaNWuqnpv6RpzGpT1U9ntJlCdj6qX7PVvqiUhrNJMVXyxevBjlypVTXeJlKV++PBYuXPjKVVhpxRdWVlY4efIk2rRpozLiv/POO+pWMuw7OzsbaC+IiCg1LC5exLjNv8BCF4PVpevhVJ5iL5+4wDnuiShlZAalKlWqqMYfEbdnpgTxNYrlRJuK+dQtg3oi0iLNtNgLmX9+0aJFyW4Td0okBwcHbNy4MQtKRkREmaZkSZQNvoI2Z3ZiZZn6mFT3HSxc9rlaT0T0OpIwWWZIqlatmsrtRERkijTVYk9ERGZIutsPH46huxfBOjoKu4tUxr4RE9gNn4he6/fff1c9N5s2bYoNGzYge/bsrDUiMkkM7ImIyPhNmoSCG1eha+6olw8L1InXQ4uIKDE7d+5Ez5498ddff6menEREpoqBPRERaYOPDwa93wqOtlY4fvMxNp5+dQpUIiK56Hfq1ClVEbNnz8acOXNgba2p0adERKnGwJ6IiDQjt7Md3q1VRN2fsuk8oqJjDF0kIjIikZGR6NWrF3x8fBAYGKgSKVtYcC56IjJ9DOyJiEhT+tYpiuyONrgU/BTLj902dHGIyEiEhYWp6YxlFiVppffw8DB0kYiIsgwDeyIi0hQXext8UK+4uv/j5gsIj4w2dJGIyMAePnyIxo0bY8eOHVi7di26du1q6CIREWUpBvZERKQ53WsUgqerPe6EhGPRgeuGLg4RGVhoaCiePn2Kbdu2qQz4RETmhoE9ERFpjr2NFT5q9HIe+5+3X0JoeKShi0REBnDhwgU8fvwYhQsXxrFjx1C9enV+DkRklhjYExGRJrWrnA/FcjvhUVgk5u66YujiEFEW8/f3h5+fHz755BP12NKSp7VEZL74C0hERJpkbWWJYU291P25e67i3pMIQxeJiLLIxo0b0aBBA3h5eWHy5MmsdyIyewzsiYhIs5qW8UCF/K4IexGN6dsuGro4RJQFlixZglatWqF+/frYvHkz3NzcWO9EZPYY2BMRkWbJ/NQjmnmr+38cvIEbD8IMXSQiymT3799Ht27dsGLFCjg6OrK+iYjYFZ+IiLTOr3gu1C6RC5HROvyw5YKhi0NEmUCn06mp7MSgQYMwf/582NjYsK6JiP7DFnsiItK84U1fttqvDLiNs3dDDV0cIspAUVFR6Nu3r+p6f/LkydjeOkRE9H8M7ImISPPK5XdFy/Ke0OmAKRvPG7o4RJRBnj9/jg4dOuDXX3/Fb7/9hnLlyrFuiYgSwcCeiIhMwseNS8LK0gJbzwXj0LWHhi4OEaVTSEgImjZtik2bNmH16tXo0aMH65SIKAkM7ImIyCQUzZ0NHasWUPcn/XtOjcklIu2ysrKCk5MTtm7dihYtWhi6OERERo2BPRERmYwPG5aAnbUlDl9/hG3ngg1dHCJKg4sXL+LcuXPIli0b/v33X9SoUYP1SET0GgzsiYjIZHi42qNXzSLq/uQN5xEdw1Z7Ii05cuQIatasiQ8//NDQRSEi0hQG9kREZFLer1sMLvbWOB/0BKuP3zZ0cYgohaTLfb169VC0aFEsXryY9UZElAoM7ImIyKS4OtrgvXrF1P3vNl3Ai6gYQxeJiF7j77//VuPopbVeAvxcuXKxzoiIUoGBPRERmZxefkXg7myHW4+eY8nBG4YuDhG9hpubG7p27aqy30vCPCIiSh0G9kREZHIcbK0wuGEJdX/atot4FhFl6CIRmTXJd7H/8gOsCritbuWxzFyxbNkyREdHo0GDBliwYAFsbW0NXVQiIk2yNnQBiIiIMkOnagUwd/cVXHsQhnl7rsYG+kSUtTacuovxa87gbkh47DoPZxvkPLEY65f9jo0bN6JJkyb8WIiI0oEt9kREZJJsrCzxcRMvdX/2rit4+OyFoYtEZJZB/fuLjsYL6nVRL3Dy9/FY/9ciDPniOwb1REQZgIE9ERGZrJblPFEmrwueRkRhxvZLhi4OkVmR7vbSUh930smYyHAE/TUOz68chnvb0ThgXYHTUhIRZQAG9kREZLIsLS0wvJm3uv/7geu4/fi5oYtEZDYOXn0Yr6VeWFjbwjZXQbh3/AIOJXzU87IdERGlDwN7IiIyaXVK5IJvUTc17d1PWy4YujhEZiP4yf+D+sjHgXh+9RgsLCzh1vg92Bcom+h2RESUNgzsiYjIpFlY/L/V/u8jt3Ax6Imhi0RkFtyd7dXti6ArCFz0CR7tmA9dTHSS2xERUdoxsCciIpNXuWAONC2TBzE6YMqm84YuDpFZqF7EDY4PziHwj5Gwds6FPB2/hIWlVezzFgA8Xe3VdkRElD4M7ImIyCx80sQLlhbAxtNBOHbjkaGLQ0n4+uuv4efnB0dHR2TPnj1F9STzoY8bNw558+aFg4MD6tWrh9OnT7OODWzjhn9x6ffRsPMsCY/OE2DllD1eUC/GvlEaVvIfk4iI0oWBPRERmYUSeZzRvnJ+dX/ShnMqGCTj8+LFC7z11lt4//33U/yayZMn4/vvv8f06dNx6NAheHh4oHHjxnjyhMMuDMnb2xv9+vbBshUrkdc9fqu8h6s9ZnarjGZlPQ1WPiIiU2Jt6AIQERFllSGNS2JVwB0cuPIQuy/eR52SuVn5Rmb8+PHq9tdff03R9nKB5scff8To0aPRrl07te63335Dnjx58Mcff6B///6ZWl569fOYO3cuunbtiqJFi6qLLaJFhYIq+70kypMx9dL9ni31REQZhy32RERkNvJld0D3GoXU/ckbzyFGBt2Tpl29ehWBgYFo0qRJ7Do7OzvUrVsX+/btM2jZzE1MTAzmzJmDAQMGYPXq1fGekyC+RrGcaFMxn7plUE9ElLHYYk9ERGblg/rFsfTQTZy6HYp1J++iWWm22muZBPVCWujjksfXr19P8nURERFq0QsNDVW3kZGRaolL/zjheopfn7169cKGDRswbdo0dO7cmfWVDH6nUo51xXoy5+9UZCrKyMCeiIjMipuTLfrVKYrvN1/Ad5vOo6FXTkMXyeRJYjt9F/ukyNj4qlWrpmtaw4RdwhOui2vixImJlmnTpk0qcV9iNm/enObymXpL/RdffKESFg4bNgwFChTA+vXrDV0sTeB3inXF7xT//yUnLCwMKcXAnoiIzM67tYrgt33XcO1BGP4+ehuuhi6QiRs4cKBqwU1O4cKF0/TekihP33Lv6fn/RGzBwcGvtOLHNWrUKAwdOjRei70EpNKl38XF5ZUWEwnAJCGfjY1Nmspp6m7evImSJUuqlnvW0+vxO5VyrCvWkzl/p0L/602WEgzsiYjI7DjZWWNQg+IYt+YMpm+/gmGlDF0i05YrVy61ZIYiRYqo4F5O0ipVqhSbWX/nzp2YNGlSkq+TcfiyJCQneUmd6CX3nDm6du2a6uHQr18/DB48WJ0sS0s96ynlWFesq4zG75Rp1VVqysfkeUREZJa6+BRE/hwOCH4SgV2BnEfbWNy4cQMBAQHqNjo6Wt2X5enTp/GmUVuxYoW6L93thwwZggkTJqh1p06dQs+ePVV3esnMTpnjxIkT8PPzU1MNPnv2jNVMRGRgbLEnIiKzZGdthaGNS2LosuPYctsSIc8jkcvIr9ybg88//1xNV6enb4Xfvn076tWrp+6fP38eISEhsdsMHz4cz58/V9nYHz16BB8fH9WS7OzsbIA9MH27d+/GG2+8oXpLSLI8JycnQxeJiMjsaarFXg7W3bt3h6urq1rk/uPHj1/7urNnz6J169bqNXKQ9/X1VS0BRERk3mTqLa882fA82gKzd181dHHov/nrJfFdwkUf1At5LK3yetJqLwn67t69i/DwcNUNv2zZsqzPTArqJQ+BXHDZsWNHsnkMiIgo62gqsJcuddIdT64OyyL3JbhPzuXLl1GrVi3VbU8OQMePH8eYMWNgb2+fZeUmIiLjJHNpD21cQt3/bf8NBO7YByxcCPj7G7poREapYsWKKungv//+qxpMiIjIOGimK760ukswf+DAAdXFTsyZMwc1atRQXfK8vLwSfd3o0aPRokULNQZMr2jRollWbiIiMm71S+ZCEWcdrj6JwdSJf2DCpp9fPjF8OJBM8jUicyE9JH766Sc0b95cnW99/fXXhi4SERFpNbDfv3+/ujKsD+qFdKmXdfv27Us0sJd5VdetW6fG3jVt2hTHjh1T48Fkips333wzyb8lU7XIknCaAcn2Kkt66F+f3vfRCu6v6eJna9rM6fONiorCGwWjMfW0NZZWaIKeJ/9Fkcd3gWnTADlWpGNu9cxkDp8NGZ6cS3388cf48ccf1ZCHpBpSiIjIsDQT2Mv8tO7u7q+sl3XyXGJkDlvJovvNN9/gq6++UtPeSKt/u3btVBKeunXrJvq6iRMnYvz48a+sl0Q8kmU3I8i0POaE+2u6+NmaNnP5fIu5AKWzx+DMYyuMGvojepSIeflEcDCwfj2MUVhYmKGLQCZOpg3s3bs3/vjjD/z8888qOSERERkngwf2kuwmsSA6rkOHDqlbuVKcWPewxNbrrzKLNm3a4KOPPoodGyYt/LNmzUoysJcWfRk/FrfFvkCBAipZjIuLC9LbwiInyo0bNzb6eRMzAvfXdPGzNW3m9Pnq9/WL+WPRod2XOHoP+OqHT1D80W1g61ajbbHX9yYjyizdunXDypUr8eeff6Jjx46saCIiI2bwwH7gwIHo3LlzstsULlxYzZcaFBT0ynP37t1LMiNrrly5YG1tjdKlS8dbX6pUKezZsyfJv2dnZ6eWhOTkNqNOcDPyvbSA+2u6+NmaNnP6fCt0aIamZ/dho5cfZlRpi2mlLYEaNWCszOVzIcPp27cv+vXrh0aNGvFjICIycgYP7CX4luV1JEmezFl78OBBVK9eXa3z9/dX6/z8/BJ9ja2tLapVq6aS68V14cIFFCpUKIP2gIiITML48Riy/SA2bnuMtaXrYuCQOuBoYjI3Mh2wdLuXYYnSY4eIiLRBM9PdSSt7s2bN1NVjyYwvi9xv1apVvEQuMq3dihUrYh8PGzYMS5cuVRn0L126hOnTp2PNmjUcJ0ZERK8ea5rURItyHtAB+GnrBdYQmaToGB32X36AVQG31a08FmfOnEHNmjWxbNmyRHtJEhGR8TJ4i31qLF68GIMHD1Zj3UXr1q1VoB6XtM5LK75e27Zt1Xh6ufIsr5WLAP/884+a256IiCihDxuWxL+nArH+ZCDO3AlF6bzpy61CZEw2nLqL8WvO4G5IeOw6T1d7dCwQhi8Hv4P8+fOrRMOenp4GLScREZlwYO/m5oZFixYlu40k00tIMrrKQkRE9DpeHs5oVT4v1hy/gx+3XMDsHsaZPI8oLUH9+4uOqh4pcd24dB4ffz4UZcpXxK6tG5A9e3ZWLhGRxmimKz4REVFW+bBhcciEK5vOBOHU7f/3AiPSKuluLy31rzZ/ANa5CiBH3Xfg2OZzOLu4GqB0RESUXgzsiYiIEiju7ow2FfKq+9JqT6R1B68+jNf9XoQeXIHnV47AwsISzlVbIyhMp7YjIiLtYWBPRESUiMENS8DSAthyNhjHbz5mHZGmBT8Jjzds8dH2+Xi0fR4i7pxPcjsiItIOBvZERESJKJo7G9pWyq/u/8BWe9I4d2d7dauLicaD9T8h9OBy5GjYF9lrdU10OyIi0hYG9kREREkY3LA4rCwtsOP8PRy5/oj1RJpVvYibyn7/aMsveHZmO3K98QlcqraJfd7iv+z4sh0REWkPA3siIqIkFMrphPaV86n7HGtPWiYXqMa+UVoF83k6jIVT6Xrxgnohz8t2RESkPQzsiYiIkjGoQQlYW1pg98X7OHSNicVIe27fvo1u3brBr6AT5g5+A0Uq+sV73sPVHjO7VUazspy7nohIqzQ1jz0REVFWK+DmiLeqFsCSgzfww+YL+KOvLz8E0ozz58+jSZMmKmFecHAwmpUtjsalPVT2e0mUJ2Pqpfs9W+qJiLSNLfZERESvMbBBcdhYWWDf5QfYf/kB64s04dChQ6hVqxayZcuGffv2oXjx4mq9BPE1iuVEm4r51C2DeiIi7WNgT0RE9Br5sjugU7UCsRnypfWTyJgFBQWhQYMGKFGiBHbv3o38+V/O8EBERKaJgT0REVEKfFC/OGytLFUXZrbak7HLkycPfvvtN2zZsgVubsx0T0Rk6hjYExERpYCnqwO6+hRU97/fzFZ7Mk5Tp07Fd999p+63a9cOjo6Ohi4SERFlAQb2REREKfR+vWKws7bE4euPVJZ8ImMhw0NGjx6NDz/8UHXDJyIi88LAnoiIKIXyuNijm28hdZ+t9mQsoqKi0K9fP0yYMAHffvstJk+ebOgiERFRFmNgT0RElArv1S0GextLBNx8jB3n77HuyOAmTpyIBQsWqDH1n3zyiaGLQ0REBsDAnoiIKBVyO9uhR43C6j4z5JMxkO73mzZtQo8ePQxdFCIiMhAG9kRERKnUv05RONpa4cStEGw9G8z6o0wVHaNTMzGsCritbuXx3bt30aRJE1y8eBEuLi5qajsiIjJf1oYuABERkdbkzGaHd/wKY+aOy2qsfcNS7rCwsDB0scgEbTh1F+PXnMHdkPDYddkj7yNo6eewQjRevHhh0PIREZFxYIs9ERFRGvSrXRROtlY4czcUG08zCzllTlD//qKj8YL6iMBLODnrQ9x/HoMJC1agTJkyrHoiImJgT0RElBY5nGzRq2YRdf/HLRcQE6NjRVKGke720lIf91sVExmBe3+Ph7WrBzzfnoSZh0PVdkRERGyxJyIiSqM+tYvA2c4a5wKfYMPpQNYjZZiDVx/Ga6mXeeotbeyQu+1o5On8FSwdXdXzsh0REREDeyIiojTK7miL3rVettr/sPkCW08pwwQ/+X9Q/+TYejxY9z10uhjY5fOGpa1DotsREZH5YmBPRESUDhLYu9hb42LwU6w7eZd1SRnC3dletdI/3rMYDzfNgKWDS5LbERERMbAnIiJKB1cHG/SpXVTd/2kLW+0pY1Qp6IrwHb8gZO8SZK/bEzka9IGFxf9P22QOBk9Xe1Qv4sYqJyIiBvZERETp1atmYRXgX773DGuO32GFUrotXrQQ9w+vR87mg5Hdt0O86RT198a+URpWlpxmkYiIGNgTERGlm7O9DfrV+a/VfutFREXHsFYpTWJiXn533nnnHezfvx+LJo+Eh2v87vbyeGa3ymhW1pO1TEREivXLGyIiIkqPd/wKY96eq7h6/xlWBtxBhyr5WaGUKkFBQWjTpg2+/PJLNG7cGNWrV1frG5f2UNnvJVGejKmX7vdsqSciorg4xp6IiCgDZLOzRv//Wu2nbr2ISLbaUypcuXIFNWvWxI0bN5AnT554z0kQX6NYTrSpmE/dMqgnIqKEGNgTERFlkO41CiFXNlvceBiG5UdvsV4pRQICAuDn5wdLS0vs27cP5cuXZ80REVGqMLAnIiLKII621nivbjF1f9q2S3gRxbH29H/RMTrsv/wAqwJuq1t5LGPqe/bsifz582PPnj0oXLgwq4yIiFKNY+yJiIgy0Ns+hfDLriu49eg5/j5yC119CrJ+CRtO3cX4NWdwNyQ8tjY8stlg3JvlsHz5cuTOnRvOzs6sKSIiShO22BMREWUgB1srDKj3stV++raLiIiKZv2aOQnq3190NF5Q/yRgA47NGIz+8/fhQpgDg3oiIkoXBvZEREQZrEv1gsjjYoc7IeFYdugm69eMSXd7aanX/fdYp9Ph8b4/8XDjdNh5loCFja16XrYjIiJKKwb2REREGczexgof1C+u7k/ffgnhkWy1N1cyTZ2+pV6ni8GjLbMQsnsRXGu9jRyN+gMWlup52Y6IiCitGNgTERFlgk7VCsDT1R5BoRH48+AN1rGZkrnn9cJvnMSTY//CrekHyF6zCywsLBLdjoiIKLUY2BMREWUCO2srDGzwstX+5x2X2Wpvptnv7z+JgC46Uq13KFQBefvMhHPF5q+8xt3Z3gAlJSIiU8Gs+ERERJnkrSoFMGP7Zdx+/ByLDlxHn9pFWddmlv0+OiwEwX+Ng1PZBnCp8gZs3PLF217a7D1c7VG9iJuBSkxERKaALfZERESZxNbaEoP+a7WftfMywvYeABYuBPz9WedmkP0+KiQIgYuHI+rJPdjnL/PK9vqO+GPfKA0ry/93yyciIkotBvZERESZqH2V/Cjg5oD7T19g4dBvgR49AF9fYMQI1rsJZ79/ce8aAhcNA2Ki4fH2t7DNUxQJY3dpqZ/ZrTKalfU0RJGJiMiEaCqwf/ToEbp37w5XV1e1yP3Hjx8n+xpJTJPY8u2332ZZuYmIyHzZWFlicGErdf8Xn/Z4ZvPfWOrJk9lyb6LZ78Xj3Ytg6ZhdBfU2OV4G7jKj3ZiWpfBT54pY0tcXe0Y0YFBPRETmN8a+a9euuHXrFjZs2KAe9+vXTwX3a9asSfI1d+/ejff433//xbvvvov27dtnenmJiIhE2/AbmH3/CareOoMXVjZwivwvALxwAfDxYSWZgIRZ7XO1/EjdWto5xV/vbIc2FeOPsyciIjKbwP7s2bMqoD9w4AB8/jsJmjNnDmrUqIHz58/Dy8sr0dd5eHjEe7xq1SrUr18fRYsygREREWUNa6+SWNerFmxjouI/UbIkPwITkTCrfcKAPqntiIiIzKor/v79+1X3e31QL3x9fdW6ffv2peg9goKCsG7dOtViT0RElGV8fGD7ydD462SMPVvrTYZktfd0tY9NiJeQrJfnmf2eiIjMusU+MDAQ7u7ur6yXdfJcSvz2229wdnZGu3btkt0uIiJCLXqhoaHqNjIyUi3poX99et9HK7i/poufrWkzp883y/b1q6+AN98ELl0CihcHqlaVP5rmtzOHz0ZLJKu9ZLeXrPgSxOuT6AlmvyciIpMP7MeNG4fx48cnu82hQ4fUrSS9S0in0yW6PjHz58/H22+/DXv75LvBTZw4MdEybdq0CY6OjsgImzdvhjnh/pouframzZw+3yzbVxcXIDgYWL8+XW8TFhaWYUWijCHZ7SXLfdx57PXZ7yXoZ/Z7IiIy2cB+4MCB6Ny5c7LbFC5cGCdOnFBd6RO6d+8e8uTJ89q/s3v3bjUWf+nSpa/ddtSoURg6dGi8FvsCBQqgSZMmcJETsnSQFhY5eWzcuDFsbGxg6ri/poufrWkzp89Xq/uq701GxkWC98alPVSWfEmoJ2Pqpfs956knIiKTDuxz5cqllteRJHkhISE4ePAgqlevrtb5+/urdX5+fq99/bx581ClShVUqFDhtdva2dmpJSE54cuok76MfC8t4P6aLn62ps2cPl+t7auWympuJIivUSynoYtBRERmRDPJ80qVKoVmzZqhb9++KjO+LHK/VatW8TLie3t7Y8WKFa+0avz111/o06ePAUpOREREKfX111+rC/Yy9C179uwp6nExYsQIlCtXDk5OTsibNy969OiBO3fusNKJiMhsaCawF4sXL1YHbukSL0v58uWxcOHCeNtId3tpxY/rzz//VGPxu3TpksUlJiIiotR48eIF3nrrLbz//vspzjVw9OhRjBkzRt0uX74cFy5cQOvWrVnxRERkNgzeFT813NzcsGjRomS3kQA+oX79+qmFiIiIjJs+ee2vv/6aou1l2tuEiQ+nTZumhu3duHEDBQsWzJRyEhERGRNNtdgTERERvY703JMZc1LSlZ+IiMgUaKrFnoiIiCg54eHhGDlyJLp27ZrsTDYRERFqSTjLgIzZlyUu/eOE6yk+1lPKsa5YVxmN3ynTrKvUlJGBPREREWWqcePGxXaxT8qhQ4dQtWrVdJ8AyRS6MTExmDFjRrLbTpw4MdEybdq0SSXuS0zCLv+UONZTyrGuWFcZjd8p06orySOTUgzsiYiIKFMNHDhQBdzJKVy4cLqD+o4dO+Lq1avYtm1bsq31YtSoURg6dGi8FvsCBQqo5LwJXyvvLSeAjRs35jSDr/kMWE8p/76yrlhXGYnfKdOsK31vspRgYE9ERESZKleuXGrJLPqg/uLFi9i+fTty5nz9HPJ2dnZqSUhO8pI60UvuOWI9pQW/U6yrjMbvlGnVVWrKx+R5REREZDQkk31AQIC6jY6OVvdlefr0aew23t7eWLFihbofFRWFDh064PDhw2paXHlNYGCgWmTqPCIiInPAFvsU0E+hl5quEMm1KshYCXkvY79ClBG4v6aLn61pM6fPV6v7qj8mJTbNq5Z9/vnn+O2332IfV6pUSd1KS3y9evXU/fPnz6vM9+LWrVtYvXq1ul+xYsV47xX3Nek51mv1O5LVWE+sK36v+P9PCyI19JuemmO9hc7UzggygZw0yLg7IiIiY3Pz5k3kz5/f0MXQPB7riYhIy8d6BvYpINl179y5A2dnZzUvbnrok/PIh/O6xD6mgPtruvjZmjZz+ny1uq9yXf7JkyfImzcvLC05si4zj/Va/Y5kNdYT64rfK/7/04JQDf2mp+ZYz674KSCVmNGtIfIlMvYvUkbi/pouframzZw+Xy3uq6urq6GLYFbHei1+RwyB9cS64veK//+0wEUjv+kpPdbzEj8RERERERGRhjGwJyIiIiIiItIwBvZZTObMHTt2bKJz55oi7q/p4mdr2szp8zWnfaW04XeE9ZTR+J1iXfE7ZTh2JnrcZ/I8IiIiIiIiIg1jiz0RERERERGRhjGwJyIiIiIiItIwBvZEREREREREGsbAnoiIiIiIiEjDGNhnghkzZqBIkSKwt7dHlSpVsHv37mS337lzp9pOti9atChmzZoFU93fu3fvomvXrvDy8oKlpSWGDBkCU93X5cuXo3HjxsidOzdcXFxQo0YNbNy4Eaa6v3v27EHNmjWRM2dOODg4wNvbGz/88ANM+f+u3t69e2FtbY2KFSvCVPd3x44dsLCweGU5d+4cTPGzjYiIwOjRo1GoUCGVNbdYsWKYP39+lpWXDOvrr7+Gn58fHB0dkT179tduHxkZiREjRqBcuXJwcnJC3rx50aNHD9y5cwemLrV1JXQ6HcaNG6fqSY4X9erVw+nTp2HKHj16hO7du8PV1VUtcv/x48fJvubp06cYOHAg8ufPr+qpVKlSmDlzJkxdWupKnD17Fq1bt1avcXZ2hq+vL27cuAFTlta60uvfv786lv/4448wZY9SWU+a/U3XUYb6888/dTY2Nro5c+bozpw5o/vwww91Tk5OuuvXrye6/ZUrV3SOjo5qO9leXiev//vvv01yf69evaobPHiw7rffftNVrFhRba8Vqd1XeX7SpEm6gwcP6i5cuKAbNWqUev3Ro0d1pri/sl9//PGH7tSpU+pzXrhwofpu//LLLzpT3F+9x48f64oWLapr0qSJrkKFCjqtSO3+bt++XSeHjPPnz+vu3r0bu0RFRelM8bNt3bq1zsfHR7d582b1ffb399ft3bs3S8tNhvP555/rvv/+e93QoUN1rq6ur91efgcaNWqkW7p0qe7cuXO6/fv3q+9PlSpVdKYutXUlvvnmG52zs7Pun3/+0Z08eVLXqVMnnaenpy40NFRnqpo1a6YrW7asbt++fWqR+61atUr2NX369NEVK1ZM/f7K75AcT62srHQrV67UmbK01NWlS5d0bm5uumHDhqnzkcuXL+vWrl2rCwoK0pmytNSV3ooVK9R5S968eXU//PCDzpQ1S2U9afU3nYF9Bqtevbruvffei7fO29tbN3LkyES3Hz58uHo+rv79++t8fX11pri/cdWtW1dTgX169lWvdOnSuvHjx+vMZX/btm2r69atm86U91dOSD/77DPd2LFjNRXYp3Z/9YH9o0ePdFqT2n39999/VYDy4MGDLCohGasFCxakOFhNSC7qyv+Z110cNLe6iomJ0Xl4eKjgXi88PFy9dtasWTpTJBcU5btw4MCB2HUSKMg6CRqSUqZMGd0XX3wRb13lypXVMcdUpbWu5FislfMNQ9eVuHXrli5fvnyqMaZQoUImHdifSUc9ae03nV3xM9CLFy9w5MgRNGnSJN56ebxv375EX7N///5Xtm/atCkOHz6suoGY2v5qVUbsa0xMDJ48eQI3NzeYw/4eO3ZMbVu3bl2Y6v4uWLAAly9fxtixY6El6fl8K1WqBE9PTzRs2BDbt2+HKe7r6tWrUbVqVUyePBn58uVDyZIl8cknn+D58+dZVGoyBSEhIaqLa0q7p5uLq1evIjAwMN7/SRnuIscKUzt3iHuuJ91/fXx8YtdJN3FZl9w+16pVS/0e3b59Ww1fkN/cCxcuqPNEU5WWupLzq3Xr1qnfaqkbd3d39fqVK1fClKX1eyX1JV3Rhw0bhjJlysDU7U9jPWnxN52BfQa6f/8+oqOjkSdPnnjr5bEcxBIj6xPbPioqSr2fqe2vVmXEvn733Xd49uwZOnbsCFPeXxkLKCdpEhh98MEH6NOnD0xxfy9evIiRI0di8eLFany9lqRlfyWYnz17Nv755x+VP0LyZEhwv2vXLpjavl65ckXljDh16hRWrFihxh7+/fff6vtMlBLh4eHq90FyykiOFfo//f87czh30JP9kmAzIVmX3D5PnToVpUuXVsdVW1tbNGvWTOULkYDfVKWlroKDg1U+gm+++UbV0aZNm9C2bVu0a9dO5bEyVWn9Xk2aNEmdtwwePBjmIDCN9aTF33QG9plArubEJVdZE6573faJrTeV/dWytO7rkiVLVKKgpUuXJvrjYkr7K0nJpMeJJIGUgEj23dT2VwJF+XEfP368aiHQqtR8vhLI9+3bF5UrV1aJIOXksmXLlpgyZQpMbV+lNUOek4s21atXR4sWLfD999/j119/Zau9hslvcGIJIOMu8tuVXtLbrnPnzup7JP9PtCgr6soUzh1SU0+J7dvr9lkC+wMHDqhWe+l5JA0EAwYMwJYtW6A1mVlX8n9NtGnTBh999JFKZCtBWKtWrTSXkDqz60q+Rz/99JM6nmnt/1tW///T4m+6tpqZjFyuXLlgZWX1ytUfuZKY8Mq0noeHR6Lby5U0yS5uavurVenZVwnm3333Xfz1119o1KgRTH1/JfO4kEyiQUFB6oe3S5cuMKX9lSEVcrCQ4QaSsVjID74cJOT/rrQWNGjQAKb+f1e6si1atAjGLC37Kr0TpAu+dNPTk2zU8vneunULJUqUyPRyU8aT/6tycpacwoULp+tvyAmg9MqS7ubbtm0z6pYdQ9WVnPcI+T8p/9e0fO6Q0no6ceKEOh4mdO/evST3WYb+fPrpp6rXkFxEFeXLl0dAQIC6oKqV84msqCv5nZdjr/RuiEt+t6X3ldZkZl1J44v8XytYsGC8xoqPP/5YNcZcu3YNWpGZ9aTV33QG9hlIuknJNEqbN29WXYD05LFcRUyMtHytWbMm3joJCqQbs42NDUxtf7UqrfsqrdW9e/dWt/oDsxZk1GcrgZBMG2Zq+ys/7CdPnoy3Tq7iyo++dNnWX9ww9c9XLmzEPTE3lX2VaRvlQpx07cyWLZtaJ+NaZYpO6RJL2iQn/7JkFv0JoAzTkbHQxn5x3lB1Jb+PEtzL/0HJ2aHPhSFdpqWLsCnWk5zryfjcgwcPql5Awt/fX62TqQKT+j7JIr87ccmFSn0LtZZkZl3J73y1atVw/vz5eOvld1umLNWazKwrGVuf8KKQ5CWQ9b169YKWZGY9afY33dDZ+0yNflqlefPmqSyMQ4YMUdMqXbt2TT0vWZi7d+/+ynR3H330kdpeXqfF6e5Sur/i2LFjapEpI7p27arunz59Wmdq+ypTv1lbW+t+/vnneNODyRQaWpDa/Z0+fbpu9erVamo/WebPn69zcXHRjR49Wmeq3+W4tJYVP7X7KxlzZWoc+Wwli648L4cQma7K1Pb1yZMnuvz58+s6dOigfpt27typK1GihJp6isyDZD2WY5PMYpItW7bY45Z8N/S8vLx0y5cvV/cjIyPVFInyvQkICIj3mx8REaEzZamtKyEZ8SULvqyT6e66dOliFtPdlS9fXmXjlqVcuXKvTLeVsJ5k9iDJjC+zksj5osw8YG9vr5sxY4bOlKWlruS+/M7Pnj1bd/HiRd20adPU1IC7d+/WmbK01FVCpp4VPy31pNXfdAb2mUACOflPYmtrq6YlkZNCvXfeeUf9UMe1Y8cOXaVKldT2hQsX1s2cOVNnyvsrwUDCRV5vavsq9xPbV9lOK1Kzv1OnTlUnIHKhSgJ6+U7LyUd0dLTOVL/LWg7sU7u/kyZNUvMpy0lljhw5dLVq1dKtW7dOZ6qf7dmzZ9Uctg4ODurALnN0h4WFGaDkZAjynUjs91sCLD15LIGWkDnGE9s+4WtMUWrrSj/lnfxmyrR3dnZ2ujp16qgA35TJ9Jlvv/22ztnZWS1yP+H0oQnrSYKInj17qnnG5bdXAo/vvvtO1Z8pS0tdCbl4W7x4cVVXcjxeuXKlztSlta7MLbB/kMp60upvuoX8Y+heA0RERERERESUNsyKT0RERERERKRhDOyJiIiIiIiINIyBPREREREREZGGMbAnIiIiIiIi0jAG9kREREREREQaxsCeiIiIiIiISMMY2BMRERERERFpGAN7IiIiIiIiIg1jYE9ERERERESkYQzsiYiIiIgoSRs2bICDgwOioqJi1509exYWFha4f/8+a47ICDCwJyIiIiKiJAUEBKBMmTKwtraOty5fvnzIlSsXa47ICDCwJ6JMt2TJEtjb2+P27dux6/r06YPy5csjJCSEnwAREZERO378OCpWrBhv3bFjx1ChQgV1/8qVK1izZo2BSkdEgoE9EWW6zp07w8vLCxMnTlSPx48fj40bN+Lff/+Fq6srPwEiIiIjJq3z+iA+sXVyPD937lyir42Ojs6SMhKZOwb2RJTpZAze119/jblz52LChAn46aef1Hg96cJHRERExuv58+e4ePFivBb7mJgYHD16VAX2O3fuxGeffYY5c+agUqVKavvmzZtj+PDhqFOnDn7//XfVQ+/Ro0fqtXv37sU777yj7p8/fx4tWrRAlSpVUK9ePY7XJ0oHBvZElCVatWqF0qVLq9b6FStWqLF6REREZNwuX76sWt2l552e9Lp78OCBCuzr1q2LsmXLYuvWrap7viTZO3XqlLp4v2vXLnTv3h1Pnz5Fjhw51GtPnDihzgEiIiLwwQcfYPbs2Thy5Ag6dOigGgCIKG0Y2BNRlpCTAOmmJycHefLkYa0TERFpQM6cOVXPu4MHD6rHBw4cwMCBA1UAX6JECbXu1q1bKFCggLovuXNk+w8//DC2Vb5kyZKx76cP7FeuXIkzZ86oC//SG+Dnn3+GjY2NQfaRyBQwsCeiTCfd9d566y388ssvaNq0KcaMGcNaJyIi0gBPT098+eWX6NGjBwoWLIgZM2aoY7oE51ZWViqojzu0Tlrr/fz84j2WFn29w4cPq8cnT57Ed999p8bqyyLT53388cdZvn9EpuL/c1YQEWWCa9euoWXLlhg5cqTqjifd8atVq6a63cmYOiIiIjJuo0ePVktirl69irx588YL5MuVKxf7+OHDh6p1X0jXfGnBlwsEHh4eqjdfly5d1HMS6Md9HRGlDlvsif7Xzh3bMAhDURR1toFNaGAAOtZgIKZgDfZgA0d2lyJNpAg96ZySCldwZfvzN+1j3gboLMtS9n3vz1rMz/P89QcBAMjRdt/bcL0W5e3K3XVdH4E+TVOfrbOuaznPs4zj2I/qb9tW7vsuwzD0u/rHcTy6Dkj3qrXWp18CAAAA+I0dewAAAAgm7AEAACCYsAcAAIBgwh4AAACCCXsAAAAIJuwBAAAgmLAHAACAYMIeAAAAggl7AAAACCbsAQAAIJiwBwAAgGDCHgAAAEquN5Tuwl9iqkdnAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.figure(1, figsize=(12, 5))\n", "_, ax = plt.subplots(1, 2, num=1)\n", "\n", "# plot (x,y) true and (x,y) predicted\n", "ax[0].scatter(x,y, s=10, c=\"r\", label=\"Data\")\n", "ax[0].plot(x_true, y_true, label=\"Solution\")\n", "ax[0].legend()\n", "ax[0].grid()\n", "ax[0].set_xlabel(r\"$x$\")\n", "ax[0].set_ylabel(r\"$y$\")\n", "ax[0].set_title(\"Deformed rod\")\n", "\n", "# u_true vs u_predicted plot\n", "u_min = min(u_true.min(), u.min())\n", "u_max = max(u_true.max(), u.max())\n", "ax[1].scatter(u_true, u)\n", "ax[1].grid()\n", "ax[1].plot([u_min, u_max], [u_min, u_max], \"k--\", linewidth=1)\n", "ax[1].set_xlabel(r\"$u_{true}$\")\n", "ax[1].set_ylabel(r\"$u_{sol}$\")\n", "ax[1].set_title(\"True-predicted plot\")\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "1c75f436-76c3-4dec-baf5-a60bc8edc9be", "metadata": {}, "source": [ "## Non-variational formulation" ] }, { "cell_type": "markdown", "id": "314f4414-e7f0-4c32-99a7-fbd48383a65f", "metadata": {}, "source": [ "Computing the first variation of the discrete elastica energy, we obtain the discrete Euler's Elastica equation:\n", "$$\\delta \\star k - f \\cos(u) = 0,$$\n", "which is equivalent, using the definition of the codifferential, to\n", "$$-\\star d k - f \\cos(u) = 0.$$" ] }, { "cell_type": "code", "execution_count": 12, "id": "e3580f9f-fcd2-4f3e-b8d2-f20acce5f7a9", "metadata": {}, "outputs": [], "source": [ "def obj(x):\n", " # apply Dirichlet BC at left end\n", " u = jnp.insert(x, 0, u_true[0])\n", " u_coch = C.CochainD0(S, u)\n", "\n", " # dimensionless curvature at primal nodes (primal 0-cochain)\n", " du = C.coboundary(u_coch)\n", " curv = C.cochain_mul(int_coch, C.star(du))\n", "\n", " load = C.scalar_mul(C.cos(u_coch), A)\n", "\n", " # dimensionless bending moment\n", " moment = curv\n", "\n", " residual = C.sub(C.codifferential(C.star(moment)), load)\n", "\n", " return jnp.linalg.norm(residual.coeffs[1:])" ] }, { "cell_type": "code", "execution_count": 13, "id": "56c51252-73bf-4cdd-976c-4632f47949da", "metadata": {}, "outputs": [], "source": [ "prb = optctrl.OptimizationProblem(\n", " dim=num_elements-1, state_dim=num_elements-1, objfun=obj)\n", "\n", "prb.set_obj_args({})\n", "sol = prb.solve(x0=u_0)" ] }, { "cell_type": "code", "execution_count": 14, "id": "360bdada-720a-4a02-ae65-09a3de5dd707", "metadata": {}, "outputs": [], "source": [ "u = jnp.insert(sol, 0, u_true[0])\n", "\n", "x, y = reconstruct_xy(u, h, num_nodes)" ] }, { "cell_type": "code", "execution_count": 15, "id": "a8b02125-91db-404a-ab64-8e800c933c22", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA/YAAAHVCAYAAABbv2DlAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjgsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvwVt1zgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAmw9JREFUeJzs3QdYFEcbB/A/vQkoKoK9g70LYu81GktsUaPGkhg1xsQWY9QUjcY0NWqsSdQYTWLX2HvDir33BlhBRZBy3/OOOT5AQPrd3v1/z7Pe3d7eMTt33u67M/OOhU6n04GIiIiIiIiINMnS0AUgIiIiIiIiorRjYE9ERERERESkYQzsiYiIiIiIiDSMgT0RERERERGRhjGwJyIiIiIiItIwBvZEREREREREGsbAnoiIiIiIiEjDGNgTERERERERaRgDeyIiIiIiIiINY2BPZGR+/fVXWFhYxC729vbw8PBA/fr1MXHiRAQHB6fr/bdu3YqqVavCyclJvf/KlSthCnbs2KH2R24NpXDhwujZs6fB/j4REWWsuMfj5BZDHnuMiRwD5VgYl9TPuHHjUvU+d+7cUa8JCAjItPOsa9euGfz4HxYWpvaT3x/KCNYZ8i5ElOEWLFgAb29vREZGqmB+z549mDRpEqZMmYKlS5eiUaNGqX5PnU6Hjh07omTJkli9erUK7r28vPjpERERJWL//v3xHn/55ZfYvn07tm3bFm996dKlWX/J1GH+/PlTHdiPHz9eBcwVK1Y02bqVwF72U9SrV8/QxSGNY2BPZKTKli2rWtb12rdvj48++gi1atVCu3btcPHiReTJkyfVB8qHDx+ibdu2aNiwYYaU8/nz56pXgVz91vKB1dHR0dDFICIiI+Pr6xvvce7cuWFpafnKeq0fVzLzWP66uiKijMGu+EQaUrBgQXz33Xd48uQJfvnll3jPHT58GK1bt4abm5s6OFeqVAnLli2LfV66eumvmI8YMUIdvON2l5MeARLsOzs7q5MRPz8/rFu3LtHua5s2bULv3r3VCY5sGxERoa40y8UIuTIvr3VwcFDvLz0PhLxX5cqV1fblypXDhg0bXtk/uVjRtWtXuLu7w87ODqVKlcLPP//8ynbnzp1Ds2bN1HvlypUL7733nqqTlJB6kH04evQoOnTogBw5cqBYsWLqufDwcIwaNQpFihSBra0t8uXLhw8++ACPHz+O9x7Si2L48OFqiISUQS62HDx4MEV/n4iITIv++Ldr1y51/JPjghwjk+uGnljX7cDAQPTv318dq+UYJMciac2NiopKUTnkbw0cOFCdH0jPPDmOSk+CP//8M8XHciG9AmvUqKF69WXLlg1NmzbFsWPHXvl78j7S609/vP7999+TLFfCOrh9+zb69euHAgUKqH3NmzevOiYHBQWpbunVqlVT2/Xq1St2uEPc93jdOY/egQMHULNmTbWN/A05xssxPCXk85H9P336tDo/kvqQupI6lgs3r3Pjxg1069Yt3jmNnMPFxMSo52UogLyfkM9Zv58c0kdpxRZ7Io1p0aIFrKys1AmEnnQLlEDXx8cHs2bNgqurqzqQd+rUSR185CDRp08fVKhQQbX2Dxo0SAXQcqARO3fuROPGjVG+fHnMmzdPrZ8xYwbeeOMNLFmyRL1PXHIi0LJlSyxcuBDPnj2DjY1N7EmJHIQl6JUTk2nTpqltb968ib///huffvqpKtsXX3yBN998E1euXFEHWnHmzBl1QqS/eCFB88aNGzF48GDcv38fY8eOVdvJQb9u3brqb0oZpdfC4sWL1YE2NaQeOnfurC4KyD7IMAUpk+QgkAN/7dq1ceLECfV35WKFLPr66tu3rzqB+eSTT1S9nTp1Sr1fSi8uEBGRabl7964K4uT4N2HCBNWqnxpy/Kxevbp63eeff64uOMtx56uvvlIBoP4i+evIMDs5J5DjrASicpzs0qULrK2tVeD8umO5lP2zzz5Tx3K5ffHiBb799lt1TJQL2PohBxLUyzZt2rRRx+yQkBAVeMvFgdftuwT1ErhLgC3nBXLu8eDBA3XMf/TokWoEkP3Vl0HKKPSNEyk559GfV0hALhdRpLxy8ULq448//kjx5yJllPMuueAycuRI7Nu3T30m169fx5o1a5J83b1799Q5jdSfDN+QMqxdu1adN1y+fFmVw9PTUzVyyL68++676jxN6IN9olTTEZFRWbBggU7+ax46dCjJbfLkyaMrVapU7GNvb29dpUqVdJGRkfG2a9Wqlc7T01MXHR2tHl+9elW997fffhtvO19fX527u7vuyZMnseuioqJ0ZcuW1eXPn18XExMTr2w9evR4pUx169ZVzx0+fDh23YMHD3RWVlY6BwcH3e3bt2PXBwQEqG2nTp0au65p06bqb4WEhMR734EDB+rs7e11Dx8+VI9HjBihs7CwUO8RV+PGjdV7bt++XZecsWPHqu0+//zzeOs3bNig1k+ePDne+qVLl6r1s2fPVo/Pnj2rHn/00Ufxtlu8eLFa/8477yT794mISLvkN97JySnR49/WrVtf2V7Wy3EnoUKFCsU7XvTv31+XLVs23fXr1+NtN2XKFPUep0+ffm3ZZDs53gYGBsY7lss5QvHixWPXJXUsv3Hjhs7a2lo3aNCgeOvl3MDDw0PXsWNH9VjOKfLmzaurXLly7PmBuHbtms7GxkbtW3J10Lt3b7XdmTNnktwXOQeS10lZE0rpOU+nTp2SrA95bzknSo58PrLdTz/9FG/9119/rdbv2bMnyc9z5MiRaht/f/94r33//ffVOcz58+fV43v37iX5HSFKLXbFJ9Kgl8fJly5duqS6pr/99tvqsXTZ0y9ylVlaEc6fP5/ke8lVen9/f3UlX7qc6UmvgO7du+PWrVuvvF7G+ydGrj5XqVIl9rF0kZMuaJL4Rt8yL6Q7mpAr3vou8NJSLmP/5Yp6wn2Q56U7nf5KfZkyZVTvg7ikB0JqJNwHfSKkhF3g3nrrLdXqIeXT/32hr289SUooLSJERGR+ZFhXgwYN0vx6ac2V2W/kWBn3GNi8efPYnnUiOjo63vP6bt160kIdN/+OHMulJVvOFeR4ntxxUFrM5T179OgR729IN3bpKafP3C7nBJKzR467ccfkFypUSLVSv86///6r9lV/LpAaqTnnkeN1UvWRGgmP9/rzDf35QGLknEJ6N0gvjLjkHEPO4RImXyTKCAzsiTRGAnHpsqYPlKVrupDuXdKNLu4yYMAA9Zx0ZU+KdHuTg4wE5Qnp/4b8vbgS21YfyCckY+cSrpd1QgJ2/fvLQVm67ifcBzlQx90H2Va66SeU2LrkJNwHeV8JzBN2gZOTFnlvfR3obxP+PXltzpw5U1UGIiIyDUkdF1NKjuXStTvhMVAuZMc9BkqgGvd5/Vh+veSOj687luvPJ6SbfMJyyLj7uMfh1/2t5Eg39dRmyU9YxpSc82TE+UJix/ak6jMueS4151VEGYHNS0QaI0no5Iq9floUSR4nZFy4jPNOTHJT2kkrg4yHk6vcCckV+bh/Qy+js+ZKGfQ9BCRZXWIkiZCQA6yMRUwosXXJSbgP8r5ycUFOOOIG93LRQ95bn8hHf4CXdZJcT09eywM1EZF5Suq4KLlZ9Enp4kp4vJDjrIw1//rrrxN9H31AKInx4uZzSXh8Tu74mDBATVhm/XtJThxpfU9K3ONgUn8rOXKMTdh7IKVSc86TEecL+mN73LpLqj7jkudSc15FlBEY2BNpiGRYlavUkihGErnoD2AlSpTA8ePHVdKb1JJu5pKAZvny5ZgyZYrKZi+ke9+iRYvUVXXJrpuZpPu9dMuTrLtyYqNv0U+MbDd58mS1v3G746cmGU5ipBVE3lf2WaYV1Pvnn39ULwn99ID6CyqSsC/usAPJxpvSzMVERGQeJGmaJGKNS7phP336NN66Vq1aYf369SppnlzsTkpyF+qFDBuTVm1993NpCJDWdnnf17WSS/Z7aaGW5G5JDbnTl0FaoyW57tChQ2MvEMjwOkkuF3foXWJkeIEk7JMu80ntjz5ZrUzDl/Bvp/ScR84XJJlgYvWRGnK8l0S+Cc83kpt3Xs4ZJk6cqGbgkWSAepJ4V+pLypbcfhKlBQN7IiMlmdb148aCg4Oxe/dulSVWWrZXrFgRr1VZruDLgVIOyjJ+S1qSZb76s2fPqoPKX3/9lezfkoOPZHeXA41cOJDAWjK2ShnkwJ0Vc9T/9NNPato4ybz7/vvvq5MhaZWQ8XTSPVE/Hm3IkCGYP3++ypIrmWn1WfFlzF16yP5L/clUgKGhoWp6HH1WfJlGR3oTCBkTKJmPf/zxR9X1r1GjRqqe5KKIi4tLhtQFERGZBjl2jBkzRmW6l3Hqkql9+vTp6gJ9XJLFfvPmzWqMugSREsDKcDXJiC8Bv2R/T0n3dWkJlrH+8jf1WfHl+JhwyrvEyHFXyjF69Gg1a41ka5eLDBIYS0Z8eT+Zlk16+Ummd8niLrlxZKYYmRZWsuKnpJu7/A0ZZ1+nTh2VFV+mwJXXS4Z4uVDg7e2tLkRIQ4Mc3+W4KzmA5IKBLCk955GM+hLYS31I/UsjgkyhKxfrU0rOhyTrv1yIkZ57+qz48vflnCUp0kAgQbycq8j+Sg8I6XEpn4ec4+gbTGSKYXlu1apV6mKADF2UzzDudMREKZbqdHtElKn02Wr1i62trcpYL1l3J0yYoAsODk70dcePH1cZa2VbyTYrGWwbNGigmzVrVuw2SWXFF7t371bbS7ZfySIrmfLXrFmT4oz9Ur4yZcq8sl4yxbZs2fKV9fI+H3zwQbx1Uj7JlpsvXz61D7lz59b5+fnpvvrqq3jbSSZdyYIv2fLd3Nx07777rm7VqlWpyoovmWgTev78ucq6L2WWvy/ZdSWD7aNHj+JtFxERofv4449VXUsZpK7279//SlZcIiIyj6z4iR3/9MeL4cOH6woUKKCOrbKtzOqS2PFCjkuDBw/WFSlSRB2D5PhWpUoV3ejRo3VPnz59bdn0x9UZM2boihUrpt5DMsDLrC2pmX1n5cqVuvr16+tcXFx0dnZ2qqwdOnTQbdmyJd52c+fO1ZUoUUKdp5QsWVI3f/58tU+vy4ovbt68qY73cq4i5ZQs+3IOExQUFLvNkiVLVPnl+YTvkZJzHrF37151jJb9kG2GDRumZrlJaVZ8+axPnDihq1evnvr85DOR84KEn0din6fMcNC1a1ddzpw5VRm9vLzU+Zc+a7+e1Ktk+ZcycnYdSg8L+SfllwGIiIiIiMjYSO86yVMjPQIo/aQ3gOQbSDhsgshYMSs+ERERERERkYYxsCciIiIiIiLSMHbFJyIiIiIiItIwttgTERERERERaRgDeyIiIiIiIiINY2BPREREREREpGHWhi6AFsTExODOnTtwdnZWU4kQEREZmsxW++TJE+TNmxeWlrxOn1481hMRkZaP9QzsU0CC+gIFCmTU50NERJRhbt68ifz587NG04nHeiIi0vKxnoF9CkhLvb5CXVxc0vWhREZGYtOmTWjSpAlsbGxg6ri/poufrWkzp89Xq/saGhqqLjrrj1GUecd6rX5HshrriXXF7xX//2lBpIZ+01NzrGdgnwL67vdyoM+IwN7R0VG9j7F/kTIC99d08bM1beb0+Wp9XzlELPOP9Vr/jmQV1hPrit8r/v/TgkgN/qan5FjPQXlEREREREREGsbAnoiIiIiIiEjDGNgTERERERERaRjH2BMRmbHo6Gg11iwhWWdtbY3w8HC1jSkz1n2VcX9WVlaGLgYRERFpAAN7IiIznRc1MDAQjx8/TvJ5Dw8PlSHc1JOzGfO+Zs+eXZXN2MpFRERExoWBPRGRGdIH9e7u7iozbMLAMSYmBk+fPkW2bNlgaWnao7aMcV/lYkNYWBiCg4PVY09PT0MXiYiIiIwYA3siIjMj3c31QX3OnDmTDHZfvHgBe3t7owl2M4ux7quDg4O6leBePit2yyciIqKkGM8ZDBERZQn9mHppqSfjpv+MEsuDQERERKTHwJ6IyExx3Lbx42dEREREKcHAnoiIiIiIiEjDNBnYz5gxA0WKFFHjIatUqYLdu3cnu/3OnTvVdrJ90aJFMWvWrCwrKxEREREREZm+6Bgd9l9+gFUBt9WtPM4qmgvsly5diiFDhmD06NE4duwYateujebNm+PGjRuJbn/16lW0aNFCbSfbf/rppxg8eDD++eefLC87ERGlT8+ePVX3dFlknvc8efKgcePGmD9/vkqCl1K//vqrmkqOiIiIKCNsOHUXtSZtQ5c5B/DhnwHqVh7L+qygucD++++/x7vvvos+ffqgVKlS+PHHH1GgQAHMnDkz0e2ldb5gwYJqO9leXte7d29MmTIFBnH4cPxbIiJKlWbNmuHu3bu4du0a/v33X9SvXx8ffvghWrVqhaioKNYmERERZSkJ3t9fdBR3Q8LjrQ8MCVfrsyK411RgL9MRHTlyBE2aNIm3Xh7v27cv0dfs37//le2bNm2Kw4cPZ32W4REjcKJLP9x8Cjxo1Ra6ESOy9u8TEZkAOzs7eHh4IF++fKhcubLqibVq1SoV5EtLvP4icLly5eDk5KQu/g4YMEDNVS927NiBXr16ISQkRLX8yzRy33zzjXpu0aJFqFq1KpydndXf6Nq1a+xc8kREREQJSXf78WvOILFO9/p18nxmd8vX1Dz29+/fV/MvS9fLuORxYGBgoq+R9YltL6068n6enp6vvCYiIkIteqGhoepWLgSk+WKAtNBPm4avOn6NYyetMaXvPNhFvYDnl/8ib57s8HC1R15Xe3jKkl3uO8DT1Q6Otpr6iF6hry9zmarJnPbXnPbV1PZX9kGn06mu60l1X5fn9bev7eLu7w9cuACULAn4+GRGkeOVK7Ey1atXDxUqVFDDrKRXlgTs0lOrcOHCakjWwIEDMWzYMPz888/w9fXFDz/8gLFjx+Ls2bOx7ylLeHg4xo8fDy8vLxXQf/zxx3jnnXewbt06GILsp5RLPrOE89ibwneRiIhI6w5efRjbUh8dFoKwC/vgXLF57PNyRiXPy3Y1iuXMtHJYm8L0P3LSk9yUQIltn9h6vYkTJ6oTu4Q2bdqUvnmflyxB5DlLuDzVITTSAhHWtrj2LAbXrjxM8iWO1jrksAVy2L28zS63dkB225e3rraAVdK7bjQ2b94Mc2JO+2tO+2oq+2ttba1ao6UFW3pCJefJkyfJPm8/bhzsf/op9nH4hx8ifNw4ZBYJZuXCrP6Ca1ySHPXMmTPqOWmR18uZMydGjhypgnT5fRe2trbqNu5vuuxrhw4dYh/nypULX3/9NRo2bIg7d+4gW7ZsyGry+Tx//hy7du16ZZhBWFhYlpeHiIiI4gt+8jKojwoJQtCyzxET/gyOJf1g5eia6HaZRVOBvZxkSYtFwtZ5aVVJ2CqvJyeviW0vJ7ZyspeYUaNGYejQobGP5SRRunJKl34XF5e0t9g3bIjGDg7YPH8+6vbphwdWDrg753fc9SiEO4/DcTc0HHflNiQct0Oe41lENMKiLBAWBdwOSzx6t7QA3J3tkDe7Azxd9K39/7X8/7fkcLQx2FzIchIugZAkt5JEV6bOnPbXnPbV1PZXWqVv3rypAlWZLSQxcgFUAl3pkp7k74e/PyzjBPVCgnzbTp0yreVe6l5+vxP7LZb1coyQ57Zv366CeGmRl99wCYplv+V56Z4v+y37JdvG3deAgAB1Yff48eN4+PBhbM+Ax48fI2/evMhqUmYHBwfUqVPnlc8qsYsbRERElLXcne3x4t41BC/7HBbWtvDoNvmVoF6/XWbSVGAvLSwybZ2cXLdt2zZ2vTxu06ZNoq+pUaMG1qxZ80rLu4yhTOrkXMZvypKQbJ/mE/oaNYBBg1R3fOH07AmyD+6FYm/USfIloeGRKtC/8/g57oQ8f3kb57EkY4iM1iEwNEItSbG3sVSBv3Tvz5tdgn0H5JPHcjHgv27/Drbxu3hmtHTVnQaZ0/6a076ayv7KkCYJai0tLdWSGH1Aq98uUZcuJbraUtbLb14m0GfET6xM586dU1OhykULSaT33nvv4auvvoKbmxv27NmjEq/Kvsfdb7nV76u0gEtiPrmIK2Ptc+fOrWZckbwscmEgyXrIRPI39TMAJPzeaf17SEREZAqqF3HDs53zYOWYHe5vjYdVthzxnpfmERl2LdtlJk0F9kJa0rt3764CcwnaZ8+erU685ARO39p++/Zt/P777+qxrJ8+fbp6Xd++fVUyvXnz5mHJkiVZX/hJk4A335QuA8DWra898XWxt4GLhw28PJwTfT4mRof7TyNw+/Fz1covAb+6L8G/uhAQrp4Pj4zBlXvP1JIUadVXwb+6AGD///sS+Gd3UFeYrKR7ABGRnoypT836TLRt2zacPHkSH330kUqOKoH4d999FxuML1u27JULxRLkJ7wwILlXJJGe9NIS8l5EREREiZG8bNIg/Mu83zBi5VlY2jnFS6Knj57GvlE602MpzQX2nTp1woMHD/DFF1+o6Y7Kli2L9evXo1ChQup5WRd3TntpvZHn5WRPkiZJV8qpU6eiffv2htmBqlWB9etf3qaTpaUF3F3s1VIpiW3CI6NVy74+0L/7X2v/bf39x8/x7EU0HoVFquX0ncS7dsoX0cPFPrbFPzboj3Pf1SFBl399Qq20Dl8gIuMm3e2HDwcmT/7/OpntI5MT6MlBVIZYSWAeFBSEDRs2qG730krfo0cPFeBLYD9t2jS88cYb2Lt3r5r6NC5Jqic5BrZu3aqy58v2MjWqBPzyOrkofOrUKXz55ZeZui9ERESkTQsWLFDnH3Ke0bluOWTPmUtlv4875Z201EtQ36zsqwnbYe6BvZBpi2RJjH6qo7jq1q2Lo0ePwhzZ21ihcC4ntSRGxpaGhkf9181fgv7/uvr/1/IvPQCCQsMRFaNT92UBHiX6Xo62VmpMvwT6+S6chOeerSj09B7CPx8EjB0LTJiQyXtLRAbpidSuXZZlxRcSyMuMJjKmPkeOHCobvlywlez10kJfsWJFNd3dpEmTVC8uGZ8uB14J+vX8/PxU8K6/WDxixAhMmDBBHUNk+jx5P5lKb8qUKWjdunWm7xMRERFpg06nw7fffqvOHfr166eG/AkJ3huX9lDZ7yVRnvR4lu73WdXrWZOBPWUcaWGXlnZZSnkm3rIucy7ee6Lv8v/qWH+5APDg2QuEvYjG5XvP1AILT6B2N/V6q8M6bL5sizf+3omGrfzgbM9xoUQmRYL5LAjohQTeiV3ATUh6ackSlwzjimvmzJlqkTH2+kR0Xbp0UUtiM6kQERGReYuJiVHT50oDwpgxY1TC3bg9liWIz8wp7ZLDwJ5eS3XDd7VXCxA/GUTcLv+qlV8y+q/firtLV+KOS24cKlAGV9zyY2vRath6+ClsA7agXsncaFneE41K5YGTHb+CRERERERk/E6cOKHyt8mwvYEDB8KYMKqiDOvyXzR3NrXApzDw4Z9q/QsHB8ybtwRPZ6/Evw074crTaGw6E6QWO2tLNPB2V0G+3Dra8utIRERERETG5fnz5ypJngz3u3TpUmySXWOS9XP3kPkk1PovE2ReR2CIX15sHd0U/35YGx/UL4bCOR0RERWDf08FYuAfx1D5y834YPFR/HvyLp6/iJ+pmoiIiIiIyBAePHiA+vXrY/To0eqxMQb1gk2klPkJtcS4cWr8iYzjl+WTJl4qA/+6k3ex7sRd3HgY9vL+ybsqCV/DUnnQspwn6nnlVr0BiIiIiIiIstLNmzfRtGlTNR2udME3ZgzsKXNb7itXfjm9XwIS5JfN56qW4U29cOp2KNaevKOC/FuPnmPN8TtqcbK1QqPSL4P8OiUZ5BMRERERUeY7c+aMCuqtrKywZ88elJTZf4wYA3syOAnyy+V3VcvIZt44fisE6068DPJl+r1VAXfU4mxnjcYS5Jf3RK0SuWBnzZZ8IiIiIiLKeDL1bfbs2bFx40bkzZsXxo6BPRldkF+xQHa1jGpeCgG3HmPt8btYf/IuAkPDsfzYbbU421ujSWkPtCrviZrFc8HWmukiiIiIiIgofR4/fqwCegnsw8LC1H0tYGBPRsvS0gKVC+ZQy2ctS+HojUdYe+JlkB/8JAL/HL2lFlcHGzQtIy35eeFXLCdsrBjkExERERFR6ixatAiDBg3Cvn37UKpUKdja2kIrGNiTZoL8qoXd1PJ5q9I4fF2C/DtYfzIQ959GYNnhW2rJ7miDZmU8VHf9GkVzwppBPhERERERvcb333+Pjz/+GL1790aJEiWgNWzaJE0G+dWLuOGLNmXh/2lDLOnri26+BZHTyRaPwyLx56Gb6D7vIKpP2IpPV5zEvkv3ERUdY+hiE5GRGDdunJqHNr1+/fVXzXTPIyIiosTpdDqMGDFCBfUjR47E3LlzYW2tvfZvBvakaVaWFqhRLCe+erOcCvL/6OODLtULIoejDR4+e4E//G+g61x/+E7cis9WnsT+yw8QHaMzdLGJKB2Cg4PRv39/FCxYEHZ2dvDw8FBZa/fv359p9Vq4cGH8+OOP8dZ16tQJF/RTelKGefToEbp37w5XV1e1yH0Z75iUyMhIdUJWrlw5ODk5qQRHPXr0wJ07d/ipEBHRa929e1ddrJcW+4kTJ6qcX1qkvUsRREmQbvd+xXOp5cs2ZbD/ygOVWX/Daemu/wKLDtxQS25nO7QoK93186JqoRyqBwARaUf79u1VMPfbb7+haNGiCAoKwtatW/Hw4cMsLYeDg4NaKGN17doVt27dwoYNG9Tjfv36qeB+zZo1iW4viY2OHj2KMWPGoEKFCurCwJAhQ9C6dWscPnyYHw8RESXq+fPnCA8PVxeEz58/r/leeGyxJ5MN8muXyI1v2pfHodGN8Fvv6uhYNT9c7K1x70kEftt/HR1/2Y8a32zFuNWncfjaQ8SwJZ/I6EnLrcwlO2nSJNSvXx+FChVC9erVMWrUKLRs2VJtc+PGDbRp0wbZsmWDi4sLOnbsqIL/pDRo0EC9Pq4333wTPXv2VPfr1auH69ev46OPPlJX8fVX8hPrij9z5kwUK1ZMJdvx8vLCwoUL4z0vr5Uufm3btoWjo6Maw7d69eoMqx+tO3v2rAropY5q1Kihljlz5mDt2rXqpCsx0qq/efNm9TlLnfv6+mLatGk4cuSI+i4QEREl9PTpUzRv3lxdOBZaD+oFA3syeZIlv27J3JjcoQIOf9YYC3pWQ/vK+dWUeUGhEfh13zV0mLUfNSdtw5drz6js+zLWhsicyHc+7EVUvOX5i+hX1mXGkpr/bxKsy7Jy5UpEREQkuh8SlEvr/c6dO1XAd/nyZdVtPq2WL1+O/Pnz44svvlDd9WRJzIoVK/Dhhx+qMXqnTp1SwwV69eqF7du3x9tu/PjxKgg9ceIEWrRogbfffjvLexsYKxlOIYG6j49P7DoJ1GWdZChOqZCQEHURxRRO1IiIKGPdvn0bn376Kc6dO/fKhX0tY1d8Misy3319b3e1RESVxZ6L91V3/U1ngnA3JBzz9lxVS77sDmhRzgOtyudF+fyumh1rQ5RSzyOjUfrzjQapsDNfNIWjbcoOR5LMRlrK+/bti1mzZqFy5cqoW7cuOnfujPLly2PLli0qYL569SoKFCigXiOt5mXKlMGhQ4dQrVq1VJfPzc0NVlZWcHZ2VuP5kzJlyhTVyj9gwAD1eOjQoThw4IBaL70L9GSbLl26qPsTJkxQrcsHDx5Es2bNYO4CAwPh7u7+ynpZJ8+lhHSrlORH0qVfemwkRS4Mxb04FBoaqm5lmIcscekfJ1xP8bGeUo51xbrKaPxOpYz0/pIefjKMSy7+S34WY/5tT03ZGNiT2bKztkLDUnnUEh4ZjV0X7mHdybvYciYItx8/x5zdV9WSP4eDmj6vVbm8KJvP5f9Bvr8/IImzSpYE4rQuEVHmj7GXg/Lu3btVC6903Z48ebLqvi3BmQT0+qBelC5dWrXcSjfvtAT2KSXvL+PB46pZsyZ++umneOvkAoSeJHuTCwaSENDUZyKQngrJkQsvIrELqdITIyUXWOUESC7yxMTEYMaMGcluKwmSEivTpk2b1DCJxMhJIL0e6ynlWFesq4zG71TyZPibHCPkGHDz5k21GDO5AJFSDOyJANjbWKFJGQ+1SJC/4/zLIH/r2SDcevQcv+y8opaCbo4qyG+5fRnKTB6L2NPM4cOBSZNYl6RZDjZWquVcTw56T0KfwNnFGZaWlpn+t1PL3t4ejRs3Vsvnn3+OPn36YOzYsaqVPLWBoexfwuEAab16n/BvJPZ3bWxsXnmN1LcpGzhwoAq4XzfzgPS2SCwfwr1795AnT55kXy+fmQxxkN4a27ZtS7a1Xkj3S/m+6OkvCjVp0uSV18p7y8myfN8Sfn7EekoLfqdYVxmN36nkSVJWGVYn4+pluJYM79LCb7q+N1lKMLAnSiTIb1bWQy0yxnj7+WDVXX/ruSDceBiGmTsuY6ZFFRTp+wt6H1qJzic2wWbyZKBdO6ByZdYnaZIEl3G7w0ugGWVrpdZldmCfEaRVXsbdy60kTJMr8PpW+zNnzqiDeKlSpRJ9ba5cueIFk9HR0WqMfNzu85IMT9YnR95fEvvJVGt6cuKQ1N81J1LHsryOJMuTz0qGJkhSROHv76/W+fn5vTaov3jxosppkDNnztf+LZkqUZaE5CQvqRO95J4j1lNa8DvFuspo/E69asmSJWoY3Lp169CoUaPY/CtaqKvUlM/4z9aIDMjB1gotynni57cr4+iYxpjetRKau7yAXWQErrrlw5imH6BZ7+nYXLw6dOc5nzVRZnvw4IHKYr9o0aLYsfR//fWX6oovmfDlgC1d3SUhnUyBJgGiBNoyDr9q1aqJvqcE8NL9Wg74kkhHxsgnnDddWpN37dqlEu7cv38/0fcZNmyYGv8vY/8lwJT5cCXx3ieffJIpdWGK5CKI5BqQHAqSn0AWud+qVSuV8V7P29tbJSsUUVFR6NChg5rabvHixeoCjIzHl+XFixcG3BsiIjK0adOmqXMCSaIr5wKmjIE9UQpJy6Uk05vZMC+OTO+G8ZtnwS0sBJdzFkDf9p+jS3AenLqd8u4yRJR6khFfMqb/8MMPqFOnDsqWLavmL5fgb/r06arngbTc58iRQz0vgb7Mdb906dIk37N3796qm7hczZeDfpEiReK11gvJiH/t2jU1lV3u3LkTfR/Jxi/j6b/99luVrO+XX37BggUL1HR5lHISnEsyI+kSL4tcqEk4baAkP5JWfH33ShkzKbcVK1aEp6dn7JKaTPpERGQ6dDqdOj8YPHiwGnYlF96NvXU+vdgVnyi1fHyQbcggvDN5Mtqe2oaZvm9hXo32OHA/Em1nHUDVXJao+Pg5CuU27R8PIkOQrtOS8EaWpBQsWBCrVq1KNpGbLHpyoP/uu+/UfOlJDTuQKdeOHz8eb51cCNDPda/3/vvvqyUpiU3tl7B3gLmTWQikR0Zy4taj9KbgFKVERJQw6dzatWtVjz7pUWcOGNgTpYUkymvXDi4XLmBEyZLo5lUe3208j+XHbuPwfUs0/mkvetcsggH1i8HFngE+EREREVFmCw8PV0PmJFGeDOdKLJeKqWJXfKK0kinuundXtzLv/fedKmLFe74o7hKDF1ExmLXzMup9uwO/7buGyGjTznhNRERERGRIISEhKk9LixYtVL4VcwrqBQN7ogwk89wPLB2DX7pVQrHcTnj47AXGrj6NJj/swsbTgewuSkRERESUwQIDA1WeHBk2N3PmTFhZpX4qXa1jYE+UwWTK6gZeubFxSB189WZZ5HSyxdX7z9B/4RF0+uUAAm5yPC0RERERUUa4fPkyatasiXv37mH37t3qvjliYE+USaytLNHNtxB2DKuHgfWLw87aEgevPcSbP+/F4CXHcPNhGOueiIiIiCgdLly4AAcHBzUTisyWY64Y2BNlMmd7G3zS1EsF+O0r51ct+quP30HD73ZiwvqzCAmL5GdABhETw9wPxo6fERERUeJOnTqljpPNmzdHQEAAChUqZNZVxaz4RFnE09UB33WsgF41C6uAft/lB5i96wqWHb6JwQ1KqNZ9W2tea6PMZ2trq6Z1u3PnjpqTXR7L/O9xyYHyxYsXKrtsUlPAmQpj3FeZvk3KJN0KpUzyGREREdFLf//9N95++21Mnz4dffv2hbU1w1rWAFEWK5vPFYv7+GDH+XsqwL8Y/BRfrD2D3/dfw4hm3mhW1uOVIIsoI0mgWKRIEdy9e1cF90kFls+fP1dd20z9+2jM++ro6IiCBQsazQUHIiIiQ5s1axYGDBiAzp0745133jF0cYwGA3siA5Dgob63O2qXyIVlh2/h+80XcO1BGN5ffBRVCuXA6JalULlgDn42lGmkBVgCxqioKDUlTEKRkZHYtWsX6tSpAxsbG5P+JIx1XyWjr7RAGNvFBiIiIkNdiP/iiy8wbtw4DB48GD/88AMvfMfBwJ7IwAn2uvoUROuKeVW3/Nm7LuPI9UdoN2MfWpb3xIim3iiY05GfEWUKCRglkE0smJWgUoJ+e3t7owp2M4M57SsREZGWA/vTp09jwoQJGDlyJC98J8DAnsgIZLOzxtDGJdG1ekF8v/k8/jpyC+tO3MWm04F4p0ZhDGxQHNkdOcaWiIiIiMxLREQEzp49i4oVK+LPP/9kK30SOGiPyIh4uNpjcocKWDeotuqmHxmtw9w9V1H32x2Yu/sKIqJe7TJNRERERGSKQkND0aJFCzRp0gTPnj1jUJ8MBvZERqh0XhcsfNcHv/WuDq88zgh5Homv1p1F4+93qZZ86YpERERERGSqgoKCUL9+fRw5ckRlwXdycjJ0kYwau+ITGbG6JXOjVvFc+PvITXy36QJuPAzDB38cRaWC2fFZy1KoUsjN0EUkIiIiIspQV69ejW2llwS35cuXZw2/BlvsiYyclaUFOlUriO2f1MOQRiXgYGOFYzceo/3M/Riw+Aiu3X9m6CISEREREWXouPrcuXNj7969DOpTiIE9kUY42VljSKOS2DmsHjpXKwBLC2D9yUA0/mEnxq85jUfPXhi6iEREREREaebv74+nT5/C29tbBfVFihRhbaYQA3sijXF3scc37cvj3w/roJ5XbpVgb8Hea6jz7XY1XV54JBPsEREREZG2rFixAnXr1sWkSZNip+WllGNgT6RRXh7O+LVXdSx8tzq8PZzxJDwKE9afQ6Pvd2L18TtMsEdEREREmjB37lx06NABrVu3xmeffWbo4miSpgL7R48eoXv37nB1dVWL3H/8+HGyr1m+fDmaNm2KXLlyqas+AQEBWVZeoqxQu0RurBtcG992KI88Lna49eg5Bi85hjdn7MPBqw/5IRARERGR0ZowYQL69u2L999/H0uWLIGdnZ2hi6RJmgrsu3btqgLzDRs2qEXuS3CfHMmkWLNmTXzzzTdZVk4iQyTYe6tqAez4pD4+blwSjrZWOH7zMTr+sh/9Fx7GlXtP+aEQERERkdGRaZzHjx+PadOmwcrKytDF0SzNTHd39uxZFcwfOHAAPj4+at2cOXNQo0YNnD9/Hl5eXom+Th/4X7t2LUvLS2QIDrZWGNSwBDpVL4Aft1zEnwdvYOPpIGw9G4xuvoUwuGEJuDnZ8sMhIiIiIoN58eIFtm3bhmbNmmH06NH8JMwpsN+/f7/qfq8P6oWvr69at2/fviQD+7ROryCLXmhoqLqNjIxUS3roX5/e99EK7q9h5LC3wvhW3uhWPT8mb7yAHRfu49d91/D3kVt4v24RvONbEHY26bsiys/WtJnT56vVfdVaeYmIiIRkvW/fvr2an/7SpUvIly8fK8acAvvAwEC4u7u/sl7WyXMZaeLEiao7SEKbNm2Co6NjhvyNzZs3w5xwfw2nbU6gTGkLrLpuiVvPovDtpouYu+MCWhWMQeVcOjVtXnrwszVt5vT5am1fw8LCDF0EIiKiVLl37x5atmyJc+fOYf369QzqTSmwHzduXKJBdFyHDh1KcsoDGZOR0VMhjBo1CkOHDo3XYl+gQAE0adIELi4u6W5hkZPHxo0bw8bGBqaO+2scWgAYHKPD6hN38d3miwgMjcDCS1YIuPQEI/6ZAp/bZ15uOGQI8Jr/j3r8bE2bOX2+Wt1XfW8yIiIiLbh586Y61kpC9J07d6JSpUqGLpJJMXhgP3DgQHTu3DnZbQoXLowTJ04gKCgo0as+efLkydAySSbGxLIxyglfRp30ZeR7aQH31zi8Va0Q3qiYH/P2XMXMrRdwMsoZ3dqPR9dj/2Ls1l9gN3Ei0KYNEGfIy+vwszVt5vT5am1ftVRWIiKibNmywdvbG1OmTEHx4sVZIaYW2Ms0dLK8jiTJCwkJwcGDB1G9enW1zt/fX63z8/PLgpISmQZ7Gyt8UL84Ol07gB+X7MPiSs3xR6XmOJ2nGGasnIh8Fy6kKrAnIiIiIkqK5EPz9PREkSJFsHLlSlaUuU93V6pUKZU1UeY4lMz4ssj9Vq1axUucJ1eBVqxYEfv44cOHalq8M2dedjWWDPryOKPH5RNpTa7SJfDV5plY8Nc4uD5/guN5S6JVzx+xJ3thQxeNiIiIiEzA2rVr0bBhQ3z55ZeGLorJ00xgLxYvXoxy5cqpse6ylC9fHgsXLoy3jQTu0oqvt3r1ajV+Q5I0COn2L49nzZqV5eUnMirSKj98OOpdPYq1vw1B2cBLeOToih77QvHz9kuIidEZuoREREREpFG//vor3nzzTTRv3hwzZswwdHFMnsG74qeGm5sbFi1alOw2kkwvrp49e6qFiBIxaRLQrh0KXLiAv4uVwNg7jlh6+Ca+3Xgex248xncdK8DVgeN4iYiIiCjlvv/+e3z88ceqh/XMmTNhZZW+aZbJxFrsiSiTWu67d4e9ny8mdSiPb9qVg621JbacDULr6Xtw9i4zbxMRERFRyklyvDFjxuCXX35hUJ9FGNgTUTydqxfE3+/VQL7sDrj+IAxtZ+zFimO3WEtERERElOz0sfPnz1c9qFu3bo0vvvgiw6clp6QxsCeiV5TPnx1rB9VCnZK5ER4Zg4+WHsfnq07hRVQMa4uIiIiI4nn27JkaT//ee+/h+PHjrB0DYGBPRInK4WSLBT2rYXDDEurx7/uvo9Ps/bgb8pw1RkRERETKgwcP0KhRI+zcuRPr1q1DxYoVWTMGwMCeiJJkZWmBoY1LYn7PqnCxt1YJ9VpN3YN9l++z1oiIiIjM3L1791C7dm1cunQJ27dvR+PGjQ1dJLPFwJ6IXquBdx6sHVQbpTxd8ODZC3Sb64/Zu68iwSQURERERGRGZNYymad+z549qFatmqGLY9YY2BNRihTM6YgVA/zQvnJ+yBT33266iPkXLPEkPIo1SERERGRG/P39Vdd7mcZu2rRp8PLyMnSRzB4DeyJKMXsbK0x5qzy+blsWNlYWOPHQEu1nHcCFoCesRSIiIiIzsGHDBjRo0ACTJ082dFEoDgb2RJQqMm3J2z6FsKRPdWS31eHqgzC0mb4XqwJusyaJiIiITNjixYvxxhtvqO73f//9t6GLQ3EwsCeiNKmQ3xXDykfDr5gbnkdG48M/AzBu9WlOiUdERERkgubMmYNu3bqhe/fuWL58ORwcHAxdJIqDgT0RpVk2G2B+jyr4oH4x9fjXfdfQZc4BBIWGs1aJiIiITIifnx/GjRuHefPmwdra2tDFoQQY2BNRuqfEG9bUG3N6VIWznTWOXH+EllP34MCVB6xZIiIiIg2LiorCN998g7CwMJQpUwZjx45VwzLJ+DCwJ6IM0bh0HqwZVAveHs64/zQCb8/1x5xdV6DjnHhEREREmvP8+XO0b98eY8aMwb59+wxdHHoNBvZElGEK53LCigE10bZSPkTH6PD1+rP44I+jeBrBKfGIiIiItOLRo0do0qQJtmzZgtWrV6NRo0aGLhK9BgN7IspQDrZW+L5jBXzRpoyaEm/9yUC0mb4Hl4I5JR4RERGRFlrq69atizNnzmDr1q1o3ry5oYtEKcDAnogynIy96lGjMJb2rwEPF3tcvvcMrafvxdoTd1jbREREREZMst337t0be/bsga+vr6GLQynEwJ6IMk3lgjmwdnAt1CiaE2EvojHwj2P4cu0ZREbHsNaJiIiIjMihQ4cwf/58dX/IkCEoVaqUoYtEqcDAnogyVa5sdlj4bnW8V/fllHjz9lzF23P8EfyEU+IRERERGYPNmzejfv36aio7yYRP2sPAnogynbWVJUY298asblWQzc4aB689VFPiHbr2kLVPREREZEB//vknWrZsiTp16mDTpk2co16jGNgTUZZpVtYDqwfWRMk82XDvSQS6zD6A+Xuucko8IiIiIgP466+/0LVrV3Tu3BmrVq2Ck5MTPweNYmBPRFmqaO5sakq81hXyIipGhy/WnsGgJcfwjFPiEREREWU4mYJ4/+UHWBVwW93KY72GDRti8uTJ+PXXX2FjY8Pa1zBrQxeAiMyPk501fupcEZUKZsfX685i7Ym7OB/4BDO7VUFx92yGLh4RERGRSdhw6i7GrzmDuyH/z23k4WyDQlfX4NvPh6FgwYL45JNPDFpGyhhssScig02J16tmEfzZzxfuzna4GPxUzXf/78m7/ESIiIiIMiCof3/R0XhBvS7qBU7+Pg7L5s/A7H82sY5NCAN7IjKoqoXd1JR4PkXc8OxFNN5ffBQT1p9FFKfEIyIiIkoT6W4vLfX/73QPxEQ8Q9Cyz/H8ylG4t/sMW8OLxOuWT9rGwJ6IDM7d2R6L+/igX52i6vHsXVfQbZ6/SrBHRERERKlz8OrD+C31uhgELf0ckcFX4d7pKzgUr66el+3INDCwJyKjmRLv0xalMOPtynCytcKBKw/RatpuHLnOAw4RERFRagQ/+X9QLywsLOHq1wl53p4E+/ylktyOtIuBPREZlRblPLFqYC2VRC8oNAKdfjmA3/Zd45R4RERERKnoDSleBF3Gox0L1HmUY/HqsM1dONHtSPsY2BOR0ZGgfuUHNdGynKeaEm/s6tMYsjQAYS+iAH9/YOHCl7dERERE9IrqRdzgeP8cAv8YifAbJ6B78Tze8xYAPF3t1XZkGhjYE5FRymZnjeldK+GzlqVgZWmBVQF30HbMclxt3g7o0QPw9QVGjDB0MYkogz169Ajdu3eHq6urWuT+48ePU/z6/v37q1k3fvzxR342RGS2Viz/Bxd//xR2eb3h0XkCLO0c4wX1YuwbpdU5FpkGBvZEZLTk5LxP7aL4o48PcttZ4LyFE1q/8wP2Fyj3coPJk9lyT2RiunbtioCAAGzYsEEtcl+C+5RYuXIl/P39kTdv3kwvJxGRsdq+fTs6duyItzq0x7LlK5A3d454z3u42mNmt8poVtbTYGWkjGedCe9JRJShfIrmxLpCD/HB1js4VKAMencYi1//GgufW6eBCxcAHx/WOJEJOHv2rArmDxw4AJ///l/PmTMHNWrUwPnz5+Hl5ZXka2/fvo2BAwdi48aNaNmyZRaWmojIuNSuXVv9dvbq1QuWlpZoUaGgyn4vifJkTL10v2dLvelhiz0RaYJ7mRJYuPQz1LlyBM9t7dHrrXE4nK8UULKkoYtGRBlk//79qvu9PqgXvr6+at2+ffuSfF1MTIxq1R82bBjKlCnDz4OIzE50dDTmzZuHPXv2wNraGu+++64K6oUE8TWK5USbivnULYN608QWeyLSBh8f2H/8EWZ//zX6tB+DPYUr4Z1uE/G7R0lUMXTZiChDBAYGwt3d/ZX1sk6eS8qkSZPUiezgwYNT/LciIiLUohcaGqpuIyMj1RKX/nHC9RQf6ynlWFesq4wkv2U9evTAunXrVI+lWrVqZej7m5pIDf2mp6aMDOyJSDsmTYJ9u3aYc+4Cet+zwf77QM/5B7Gwjw8qFshu6NIRURLGjRuH8ePHJ1s/hw4dis2tkZBM05TYenHkyBH89NNPOHr0aJLbJGbixImJlmnTpk1wdPx/kqm4Nm/enOL3N2esJ9YVv1dZ5/nz55gwYQLOnTuH4cOHw8PDA+vXr8/CEmjXZg38poeFhaV4Wwb2RKQtPj5w8PHBvBdR6LXgEPyvPkT3ef5Y3McH5fMzuCcyRjL2vXPnzsluU7hwYZw4cQJBQUGvPHfv3j3kyZMn0dft3r0bwcHBKFiwYLwuqR9//LHKjH/t2rVEXzdq1CgMHTo0Xot9gQIF0KRJE7i4uLzSYiIngI0bN4aNjc1r99dcsZ5YV/xeZb22bdvixo0bWLt2rQry+TtlWr9V+t5kKcHAnog0ydHWGvN7VkPPBQdx6NojdJvrjz/6+qJsPldDF42IEsiVK5daXkeS5IWEhODgwYOoXr26WidZ7mWdn59foq+RsfWNGjWKt65p06ZqvSSOSoqdnZ1aEpKTvKRO9JJ7jlhPacHvFOsqvaT3keQZKV26tGqp53cq5bRQV6kpH5PnEZFmOdlZY0Gv6qhSKAdCw6Pw9lx/nL4TYuhiEVEalSpVCs2aNUPfvn1VZnxZ5H6rVq3iZcT39vbGihUr1P2cOXOibNmy8RY5EZLuqMll0Sci0qrjx4+jQ4cOqpt2uXLlUKFCBUMXiYwAA3si0rRsdtb4tVc1NcY+5Hmkark/F5jybktEZFwWL16sTlSlS7ws5cuXx8KFC+NtI1PfSSs+EZG52bVrF+rUqYMrV67g2bNnhi4OGRF2xScizXO2t8Hv71ZH97n+OH4rBG/P8ceSfr4omcfZ0EUjolRyc3PDokWLkt1GkuklJ6lx9UREWrZq1Sp06tQJNWvWVL2WEuYDIfOmqRb7R48eqTFzMp+tLHL/8ePHySZGGDFihLry7+TkhLx586qpIO7cuZOl5SaizOciwX1vH5TN54IHz16g65wDuBT8hFVPREREmnf27Fm0a9cOrVu3VmPpGdSTpgP7rl27IiAgABs2bFCL3JfgPiky7kSmvxkzZoy6Xb58OS5cuKD+QxCR6XF1tMGid31Q2tMF95++QJc5/rh876mhi0VERESU7hwk0mK/ZMmSRBN/Ellr6SqVBPOSSMfHx0etmzNnjsqgK2PtEkuQI636CecnnDZtmsq0K9NCxJ0ah4hMQ3ZHWzX1XZc5B3Au8Am6zD6Apf1roEguJ0MXjYiIiCjFJNu9TMspQX3//v1VIlEizbfY79+/XwXq+qBe+Pr6qnX79u1L8ftIsh0LCwtkz875rolMVQ6nl8G9Vx5nBD+JUMH99QdMMENERETa8OLFC3Tr1g1Tp059bV4RIk212AcGBsLd3f2V9bJOnkuJ8PBwjBw5UnXpT25cSkREhFr0QkNDY8fsy5Ie+ten9320gvtruoz9s3Wxs8RvPSuj2/zDuHTvGTrPPoDF71ZFgRyOJrm/Gc2c9ler+6q18hIRUco8ffoU7du3x44dO7Bs2TI1tR2R0Qf248aNw/jx45Pd5tChQ+pWWtoTkitYia1P7ASoc+fOqkvLjBkzkt124sSJiZZp06ZNcHRMW1CQUMIhAqaO+2u6jP2zfacgMP2pFe6GhKPDz7sxsHQ0ctqb7v5mNHPaX63tq+SRISIi0/PJJ5+o3sr//vsvGjRoYOjikEYYPLAfOHCgCriTU7hwYZw4cQJBQUGvPHfv3j3kyZPntUF9x44dcfXqVWzbtu21WSRHjRqlxrPEbbEvUKCAmk83vRkopSxy8ti4cWPY2NjA1HF/TZeWPtsGDSPQbd4hXH0QhvnXnLG4d1Xkze5gsvubEcxpf7W6r/reZEREZBr0DZZff/013nvvPVSsWNHQRSINMXhgnytXLrW8jiTJk/HxBw8eVMnvhL+/v1rn5+f32qD+4sWL2L59O3LmzPnavyWZJhPLNiknfBl10peR76UF3F/TpYXPNp+bDZb0q4FOs/fj+oMwdF9wBEv7+8LTNXXBvVb2NyOZ0/5qbV+1VFYiIkre6dOn0adPH9X1XhoUUxKzEGkyeZ5kg2zWrBn69u2rMuPLIvclO2TcjPje3t5YsWKFuh8VFaXGpBw+fBiLFy9GdHS0Go8viySkICLz4eFqjyV9fVHQzRE3HoaphHpBoeGGLhYRERGZOUkEXrt2bTx//hzW1gZvdyWN0kxgLyQ4L1eunOoSL0v58uWxcOHCeNvI1HfSii9u3bqF1atXq1vpyuLp6Rm7pCaTPhGZBul+v6SfL/LncMC1By+D+2AG90RERGQga9euRaNGjVSMs3PnThWnEKWFpi4Jubm5YdGiRcluE3c6CBmbz+khiCiufBLc9/VVWfKv3H+m5rv/s18N5HZ+dfgNERERUWaRXGGdOnVC06ZNsWTJEtjbpyO7L5k9TbXYExFlhAJujiq493S1x+V7z9B1zgHcf/r/KS6JiIiIMpPM1JU7d241pd1ff/3FoJ7SjYE9EZmlgjlfBvd5XOxwMfgpus31x8NnzL1BREREmRvQDxs2DAMGDFA9i6tVq8Zx9ZQhGNgTkdkqnMtJBffuznY4F/gEb8/1xyMG90RERJQJZLauXr164bvvvkOZMmXU1HZEGYWBPRGZtaK5s+GPvr7Ilc0OZ++Gots8f4SERRq6WERERGRCnj17hjfffFONpf/jjz8waNAgQxeJTAwDeyIye8Xds2FJXx/kdLLF6Tuh6D7fHyHPGdwTERFRxpg6darKer9u3Tp07tyZ1UoZjoE9ERGAEnmcVcu9m5MtTtwKQY/5BxEazuCeiIiI0i4qKkrdfvLJJzh48CAaN27M6qRMwcCeiOg/Xh7OWNzHB9kdbXD85mO8M/8gnjC4JyIiojQ4e/YsypYti/3798PGxgalS5dmPVKmYWBPRBRHKU8XFdy7Otjg2I3H6LXgEJ5GvLzaTkRERJQS/v7+qFWrlgroCxYsyEqjTMfAnogogTJ5XVVw72JvjcPXH6H3gkMIe8HgnoiIiF5vw4YNaNCgAUqVKoVdu3YhX758rDbKdAzsiYgSUTafKxa+6wNnO2scvPYQ/RYdw4toVhURERElLSIiAv3790fDhg2xadMm5MiRg9VFWYKBPRFREioUyI7f362ObHbW8L/6CHPOWyI8ktE9ERERverFixews7PDjh07sHz5cjg6OrKaKMswsCciSkalgjnwW+9qcLK1woUQS7y3OIDBPREREcXS6XQYNWoUGjVqhMjISBQpUgTW1tasIcpSDOyJiF6jSiE3zO1RGbaWOuy9/AD9Fx5BRBRb7omIiMydTGfXp08ffPPNN3jzzTdVsjwiQ2BgT0SUAlUL5UB/72g42Fhi54V7eH/RUQb3REREZiwsLAzt2rXD77//joULF2Lo0KGGLhKZMQb2REQpVNwVmN2tMuxtLLHtXDA+WHwML6JiWH9ERERmaM2aNdi6dStWr16Nbt26Gbo4ZOYY2BMRpYJvUemWXw121pbYcjYIg5YcRWQ0g3siIiJzaqkXnTp1wrlz59C8eXNDF4mIgT0RUWrVKpELc3pUha21JTaeDsLgJccQuf8AsHAh4O/PCiUiIjJRFy5cQJkyZbBo0SL1uECBAoYuEpHCFnsiojSoUzI3fulWBbZWlvj3VCA++mYFot7pCfj6AiNGsE6JiIhMzOHDh1GzZk04ODigbt26hi4OUTwM7ImI0qi+tztmVssGm+hIrC1VB581/QA6eWLyZLbcExERmZDNmzejfv36KF68OHbv3s2WejI6DOyJiNKh4ZNrmL5qEixjovFnhab4o0Kzl09cuMB6JSIiMpF56sePH4/atWtjy5YtyJkzp6GLRPQKBvZEROlRsiSaXjyAYbt+Vw/HNe6PI3m91XoiIiLStsePH8PCwkJlwF+1ahWcnJwMXSSiRDGwJyJKDx8fYPhwvOf/D1qc24NIKxu83/1rBJeqwHolIiLScCv9mDFjUKFCBRXc58iRAzY2NoYuFlGSGNgTEaXXpEmwOHAA375ZGiWdrRAMO7y/+CjnuCciItKg6OhovPfee/jqq68wcOBAZM+e3dBFInotBvZERBnBxwdOPbvjl/614WxvjSPXH+GLtadZt0RERBoSHh6Ot956C/PmzcOCBQswbNgwQxeJKEUY2BMRZaAiuZzwU+eKsLAAFh24gWWHbrJ+iYiINCIgIADbtm3DihUr0LNnT0MXhyjFGNgTEWWwBt558FGjl8nzPlt5CgE3H7OOiYiIjNiDBw9UF3xfX19cu3YNb7zxhqGLRJQqDOyJiDLBwPrF0bh0HryIjsF7C4/g3pMI1jMREZERunz5MqpXr46xY8eqxxxTT1rEwJ6IKDN+XC0t8H3HCiia2wmBoeH44I+jiIyOYV0TEREZkWPHjsHPzw/W1tbo27evoYtDlGYM7ImIMomzvQ1md6+KbHbWOHj1ISasP8u6JiIiMhLbt29H3bp1UahQIezZs0fdEmkVA3siokxU3D0bvuv4ck77BXuvYfnRW6xvIiIiI7B48WI1pl6S5eXOndvQxSFKFwb2RESZrGkZDwxqUFzdH7X8JE7dDmGdExERGcitWy8vss+cORNr165FtmzZ+FmQ5jGwJyLKAkMalUR9r9yIiIpB/4VH8PDZC9Y7ERFRFtLpdBg/fjy8vLxw5coV2NjYwNbWlp8BmQQG9kREWcDK0gI/dq6EwjkdcfvxcwxachRRTKZHRESUJWQqu4EDB2LcuHH49NNPUaRIEdY8mRQG9kREWcTVwQa/dK8KR1sr7L30AJM3nmfdExERZbKIiAh07doVs2bNwuzZszF69GhYWFiw3smkMLAnIspCXh7O+LbDy2R6s3ddwZrjd1j/REREmSgwMBD+/v74+++/OaUdmSxrQxeAiMjctCzviZO3i2HWzssY/vcJlTm/lKeLoYtFRERkUoKDg2FnZ6emsTt//ry6T2Sq2GJPRGQAw5p6oXaJXHgeGa2S6T0OYzI9IiKijHL16lXUrFkT/fv3V48Z1JOpY2BPRGSgZHpTO1dC/hwOuPEwDB/+GYDoGB0/CyIionQ6ceIE/Pz81P0JEyawPsksMLAnIjKQHE62+KV7FdjbWGLnhXv4fjOT6REREaXHrl27UKdOHXh6emLPnj0oWrQoK5TMAgN7IiIDKpPXFZPal1f3f95+GRtO3eXnQURE9BrSy23/5QdYFXBb3ep7vR05cgRVqlTBjh07kCdPHtYjmQ1NBfaPHj1C9+7d4erqqha5//jx42RfI3NVent7w8nJCTly5ECjRo1UVkwiImPRpmI+vFvr5Xy6Hy87jotBTwxdJCIiIqMlF8FrTdqGLnMOqKFsclv5o7lq/UcffYSNGzfCxYVJacm8aCqwl/knAwICsGHDBrXIfQnuk1OyZElMnz4dJ0+eVN1xChcujCZNmuDevXtZVm4iotcZ1dwbvkXd8OxFNPotPILQ8EhWGhERUQISvL+/6CjuhoSrxzqdDiH7l+HE1H7o9c1i9by1NSf+IvOjmcD+7NmzKpifO3cuatSooZY5c+Zg7dq1avqK5C4GSCu9jK8pU6YMvv/+e4SGhqqkGkRExsLayhI/d62MvK72uHr/GT76MwAxTKZHREQUS7rbj19zBvpUszpdDB5tnYPHu36Ha82usMvnrZ5nMloyR5oJ7Pfv36+63/v4+MSu8/X1Vev27duXovd48eIFZs+erV5ToUKFTCwtEVHq5cxmh1+6V4WttSW2ngvG1G0XWY1ERET/OXj14f9b6qMjcX/td3hyZA3cmgxA9lpdAQsL9bxsR2RuNNNPJTAwEO7u7q+sl3XyXHKkVb9z584ICwtTGTI3b96MXLlyJbl9RESEWvSkhV9ERkaqJT30r0/v+2gF99d08bPNHN55HPHFG6UwcsVp/LjlIrzdndCw1Ku/fZnNnD5fre6r1spLRJRewU9eBvVCFxWJ6MdByNVmBJy8ayW5HZG5MHhgL8ntxo8fn+w2hw4dUrcWFhavPCfjahJbH1f9+vXVePz79++r7vsdO3ZUCfQSu1AgJk6cmGiZNm3aBEdHR2QEubhgTri/poufbcZzAFDbwxK7Ay0xZOkxDC0XjTyy0gDM6fPV2r7KxWoiInPi7myP6LAQxESEwSaHJ/J0mwwLC8tEtyMyNwYP7AcOHKha05MjCe9kTHxQUNArz0kSvNdNZSEZ8YsXL64W6b5fokQJzJs3D6NGjUp0e1k/dOjQeC32BQoUUEn30pthU1pY5OSxcePGsLGxganj/poufraZq3F0DHosOIzD1x9j6W1X/N3fB9nssu4n25w+X63uq743mamRGXAGDx6M1atXq8etW7fGtGnTkD179tfm4hkxYgR27tyJmJgYlVdn2bJlKFiwYBaVnIgym4fVU9xfMhIxNo7w6D7llaBemvo8XO1RvYgbPwwyOwYP7KVLfHLd4vUkWV5ISAgOHjyI6tWrq3XS6i7r/Pz8UvU3pZU/blf7hOzs7NSSkJzwZdRJX0a+lxZwf00XP9vMqldgRrcqeGPaHly+90x1zZ/5dhVYWibfQynjy2E+v1Va21ctlTU1JOntrVu3VMJc0a9fPzUDzpo1a5J8zeXLl1GrVi28++67qsed5NKRQN/enq12RKbi9OnTaNWqFVztLGHdYigsLSxik+gJ/dFx7BulYZXFx0oiY6CZ5HmlSpVCs2bN0LdvXxw4cEAtcl/+g3t5ecVuJ3PWr1ixQt1/9uwZPv30U7Xt9evXcfToUfTp00edMLz11lsG3BsioteTroQzu1WBrZUlNp4Owsydl1ltZNLSOgPO6NGj0aJFC0yePBmVKlVSM+G0bNkyySF3RKQt586dQ4MGDZAzZ04cO+SPuYPfUC3zccnjmd0qo1lZT4OVk8iQNBPYi8WLF6NcuXKqS7ws5cuXx8KFC+NtIwd+acUXVlZW6oegffv2aj57uQggXfd3796tuugRERm7ygVzYHybl79XUzadx47zwYYuEpFRzYAj3e7XrVunjvNNmzZVwby8fuXKlfykiEzE8+fPUblyZezatUslwpbgfc+IBljS1xc/da6obuUxg3oyZwbvip8abm5uWLRo0Wu72etJF7zly5dnQcmIiDJPl+oFceJWCJYcvIHBS45hzaBaKJTTiVVOJictM+AEBwfj6dOn+Oabb/DVV19h0qRJqtW/Xbt22L59O+rWrZvuGXC0OnNCVmM9sa4ymgTycqFOeuIMGzYMtra28f4fVi0oua9e5r+KiY5CTDTMFv//mWZdpaaMmgrsiYjM1bjWpXEuMBTHbjxGv9+PYPkAPzhlYTI9ImOdAUda7EWbNm3w0UcfqfsVK1ZULfyzZs1KMrBPyww4Wps5wVBYT6yrjCBDa3/77Td88sknKofGli1bMuR9TR3//5lWXaVmBhyeFRIRaYCdtRVmdauCVtP24HzQEwz/5wSmd6n02uk+iYxBZs6AIwl4ra2tUbp06Vdy8+zZsyfJv5eaGXC0OnNCVmM9sa4yglyskxxZEtSPHDkSn332mQrq+f+P///M8bcqNBUz4DCwJyLSiDwu9pjxdmV0mX0A607cRYX8ruhXp5ihi0Vk0BlwpGtutWrVXkmud+HCBRQqVChDZ8DR2swJhsJ6Yl2lVVRUFPr376+C+qlTp2LQoEGxXZH5veL/v4xmo4Hf9NSUT1PJ84iIzF21wm5qKh/xzb/nsOfifUMXicigM+AIGXu7dOlSlUH/0qVLmD59upoeb8CAAfx0iDTE0vJlaLJkyRIV1BNRKv7/sLKIiLSlm28hvFUlP2J0wMAlR3HzYcrHXxGZ2gw4om3btmo8vUx3J6+V6fL++ecfNS6XiIzfw4cPsXfvXhXY//rrr68dukNEr2JXfCIijZFx9V++WVaNtZds+f0XHsE/7/vBwdbK0EUjyvIZcPR69+6tFiLSllu3bqmpKp89e6aG0MjwGiJKPbbYExFpkL3Ny2R6OZ1sceZuKEYtP5FosENERGSszp49q/JnSFC/ceNGBvVE6cDAnohIo/Jmd8D0rpVhZWmBlQF3sOCPnYB0Wfb3N3TRyATI1HHSJfbo0aPx5nsnIsoIMsWlDJdxdXVV3fDj5tEgotRjYE9EpGE1iuXEpy1KqftfHw/FgdGTAV9fYMQIQxeNNK5evXq4efOmmu+9UqVKauy6jHudMGEC1q5da+jiEZHGZc+eHXXq1MGuXbuQL18+QxeHSPMY2BMRaVxv6yC8eXo7oi2tMLDNCDxwcAEmT2bLPaVLmzZtMGbMGPz11184c+aMmoLu448/VvPJb926lbVLRGkiM1bI3NwlSpRQs1vkyJGDNUmUARjYExFpnMXFi5i4YTpK3ruO+0458Hnj914+ceGCoYtGJjKvtCSzk6zzkrlaEtT98MMPhi4WEWnQjz/+iNatW6uZK4goYzGwJyLSupIl4RAVge/WfQ+rmGisK1UH671qqvVE6dWlSxfs2bNHzcbw999/q275Fy9eZMUSUYpJctdRo0ap3B3Dhw9Xt0SUsTjdHRGR1vn4AMOHo9zkyXj/wF+Y7tcZY978GD5lKiKnoctGmidzxp84cSL2sSTT69u3L3bs2GHQchGRdoJ6+c2YN28evvvuOwwdOtTQRSIySalqsZckOkREZIQmTQIOHMCgd5vAy8UKD2CLz1edNnSpyARky5YNly9fjn1cuXJl1SWfiCglpLdPqVKl8PvvvzOoJzKWFntvb2/1H3LkyJFwcnLKvFIREVHq+fjAzscH390OQZuf92LdybtoceIuWpb3ZG1Smv3yyy9488030bx5c3VyLvNOFyxYkDVKRMl6/PgxNmzYoGbTkMSbRGRELfabN2/Gpk2bVBbLBQsWZF6piIgozcrmc8UH9Yqp+2NWncL9p5yDnNJOprk7fPgwqlSpguvXr6N48eJYtmwZq5SIknTnzh01ld3AgQPZw4fIGAN7Pz8/+Pv745tvvsHnn3+uEuhwjB0RkfEZ2KAEvD2c8fDZC3y+6pShi0Ma1q1bN0RERKBTp06oVq0acufODUdHR0MXi4iM1IULF1TM8OjRI+zevRtubm6GLhKRWUhTVvwePXqo/7RvvPEGWrZsibZt2+LSpUsZXzoiIkoTW2tLTHmrAqwtLbD+ZCDWnrjDmqQ0kcR5Li4uai77YcOGqa61Q4YMYW0S0StOnTqFmjVrqot/+/btU8N3iMjIp7uTDJdNmjRBv379sHr1apQtW1aNn3ny5EnGlpCIiNLcJX9A/eLqviTSY5d8SgsbGxt1zP/111/x6aefqjH30gpHRJRQgQIF0KZNG/UbIfeJyEgD+1mzZuHdd99F+fLl4erqikaNGmHv3r344IMPMGPGDAQEBKB06dJqLB4RERnewPrFUcrTRXXJH7PylArQiFKjf//+qgu+zGEvPfTEs2fPWIlEFOuff/5Rs2dIfDB37lzkzMnJVomMOrD/+uuvERoainfeeUeNrQ8JCcHBgwcxdepU9O7dG1u3bsX777+Pnj17Zl6JiYgolV3yy6su+f+eki75d1l7lKyEU9lJz7wtW7aoLvkyI44MvfPx8WEtEpEyffp0vPXWW5g9ezZrhEgr092lZB57adEfM2ZMespEREQZqExeV3xQvzh+2npRJdLzLZoTuZ3tWMeUqFy5ciF//vyoUKFCvEVmxBGSFf+3335j7RGZOekBNnbsWHz55ZdqOuyJEycaukhEZi3NY+yT4u7ujm3btmX02xIRUTp88F+X/EdhkeyST8mSJHmTJ09WQ+sOHTqEAQMGqARYzs7ObKknolgSzEtQP2nSJEyZMgWWlhkeVhBRZrXYp4SFhQXq1q2b0W9LREQZ0CW/zfS92HA6EGtO3EXrCnlZp/QKb29vtXTu3Dm2VU4y4Q8aNAgNGzZkjRGRIkm0Je9Wr169WCNERoCX1oiIzKhL/sAG+iz5pxD8JNzQRSINkAv2zZs3x6JFi3DnDqdNJDJnkl9LevTExMSo3wUG9UTGg4E9EZGZdckv7emCx2GR+GwFs+TTq+SEPTG+vr4qcS4RmafAwEDUq1dPjaW/cuWKoYtDRJndFZ+IiIyXjZV0ya+A1tP3YNOZIKw+fgdtKuYzdLHIiGTLlg1ly5ZFxYoVVdI8ufXy8lKz4Dx9+tTQxSMiA5Cp7KTrfXh4OHbt2qWSaBKRcWFgT0RkZkrndcGgBiXww5YLGLv6NGoUywl3Z3tDF4uMxPLly3H8+HG1/Pzzz7h48aJqxZcu+ZIoi4jMy/Xr11GzZk01R/3evXtRuHBhQxeJiBLBwJ6IyAwNqF8Mm84E4vSdUIxecQqzu1dRgRtRs2bN1KInLXTSWpczZ054eHiwgojMjEx/2a9fP5VAM3fu3IYuDhElgWPsiYjMuEu+jZUFNp8JwqoAJkWjxNnb26NMmTIM6onMsPfOli1bYGVlhS+++IJBPZGRY2BPRGSmZF576ZIvpEt+cCiz5BMREfDLL7+gQ4cOWLp0KauDSCMY2BMRmbH36xVD2XwuCHkeiU+ZJZ+IyKzpdDrVOv/ee+9h4MCBKsAnIm1gYE9EZMbidsnfcjYIKwNuG7pIRERkIJIgc+zYsfjqq6/w008/wdKSoQKRVvB/KxGRmfP2cMHg/7rkj1t9hl3yiYjMVKdOnTB//nyMHj2aCVWJNIaBPRER4b14XfJPqu6YRERk+p48eYKPPvoIT58+hZeXF3r16mXoIhFRGjCwJyKiBF3yg7Hq+F3WChGRiQsODkb9+vVVK/358+cNXRwiSgcG9kREFNsl/8OGL7vkf7nuHEJesGKIiEzV1atXUatWLdy6dQs7d+5ElSpVDF0kIkoHBvZERBTrvbrFUC6fK0LDo7D0iiW75BMRmaDHjx+roD4mJgb79u1DxYoVDV0kIkonBvZERBTLOk6X/NOPLLEygF3yiYhMTfbs2VUG/L1796Jo0aKGLg4RZQAG9kREFI+XhzMG1S+m7n+1/hyCQsNZQ0REJmD16tWYPn26ut+7d2/kyZPH0EUiInMM7B89eoTu3bvD1dVVLXJfuhKlVP/+/dXUHT/++GOmlpOISOv61iqMAk461SV/1HJmySci0jpJkNe2bVvs2rWLw6yITJCmAvuuXbsiICAAGzZsUIvcl+A+JVauXAl/f3/kzZs308tJRGQKXfK7FY9WXfK3nQvGP0dvG7pIRESUBjJ96TfffIN3330X/fr1w5IlSzhHPZEJ0kxgf/bsWRXMz507FzVq1FDLnDlzsHbt2tdOz3H79m0MHDgQixcvho2NTZaVmYhIyzwcgQ8bFFf3x685jcAQdsknItKaWbNmYdSoURg3bhxmzJgBKysrQxeJiDKBNTRi//79qvu9j49P7DpfX1+1TrJ5enl5Jfo6yfYprfrDhg1DmTJlUvS3IiIi1KIXGhqqbiMjI9WSHvrXp/d9tIL7a7r42ZrH59vDJy82nQ3CiVuhGPH3cczpXsnkWnq0+l3WWnmJKPNEx+hw8OpDBD8Jh7uzPaoXcYOV5cvf6i5duqhkeXJLRKZLM4F9YGAg3N3dX1kv6+S5pEyaNAnW1tYYPHhwiv/WxIkTMX78+FfWb9q0CY6OjsgImzdvhjnh/pouframbfvWrWiZEzh92wo7L97H2N82wNddB1Okte9yWFiYoYtAREZgw6m7GL/mDO7G6VXlbq+Dy7HfMefHSShSpAiDeiIzYPDAXroFJRZEx3Xo0CF1m1grkYwbSqr16MiRI/jpp59w9OjRVLUwSXeloUOHxmuxL1CgAJo0aQIXFxekt4VFTh4bN25sFsMCuL+mi58tzOrzjXS/iimbL2LNLTu839YPnq72MBVa/S7re5MRkXkH9e8vOoq4l1ujw0Jw/PfxiHxwE3/t6IzhRYoYsIREZDaBvYx979y5c7LbFC5cGCdOnEBQUNArz927dy/JqTp2796N4OBgFCxYMHZddHQ0Pv74Y5UZ/9q1a4m+zs7OTi0JyQlfRp30ZeR7aQH313TxszWPz/e9esWx+dw9HL/5GGNWn8WvvaqZXJd8rX2XtVRWIsqc7vfSUh83qI8KDUbQ0s8RE/4EHl0mYkWgKz6O0cV2yyci02XwwD5XrlxqeR1JlhcSEoKDBw+ievXqap1kuZd1fn5+ib5GxtY3atQo3rqmTZuq9b169cqgPSAiMo8s+d+9VR4tpu7Bzgv38NfhW+hYrYChi0VEZLZkTH3c7ve66CgE/fkZdDHR8Hh7Mmzc8qnnZbsaxXIatKxElPk0kxW/VKlSaNasGfr27YsDBw6oRe63atUqXuI8b29vrFixQt3PmTMnypYtG2+RFg4PD48kk+0REVHiirs7Y2jjkur+l2vP4M7j56wqIiIDkUR5cVlYWcOtyQB4dPtWBfVJbUdEpkkzgb2Q6erKlSunxrrLUr58eSxcuDDeNjL1nbTiExFRxutbuygqFsiOJxFRGLn8pMpzQkREWU+y34uwy4fwYOPP6vfYoXBFWGdzS3Q7IjJtBu+Knxpubm5YtGhRstu87iQzqXH1RET0ejJOc8pbFdBi6m7sunAPyw7fRKdq/89jQkREWUOmtLO+vAv3/pkCh+LVgegowPr/uTdkVL2H68up74jI9GmqxZ6IiAyvuHs2fPxfl/yv1p7Fne37AOk95e9v6KIREZmN77+bgst/T0a2co3g/uYoWCQI6sXYN0ozcR6RmWBgT0REqdandlFUKvhfl/yZW6Dr0QPw9QVGjGBtEhFlsuXLl2P48OEYPXo0li1aAM8cTvGel5b6md0qo1lZT34WRGZCU13xiYjIeLrkf1sSaHHlBXYVrYJ/yjZAh1PbgMmTgXbtAB8fQxeRiMhktWnTBmvWrFFJpEWTMp4q+70kypMx9dL9nlPcEZkXttgTEVGaFL97BUP2/KHuf1O3F57YOrx84sIF1igRUQYLCwtD+/btsW3bNlhZWcUG9UKCeJnSrk3FfOqWQT2R+WFgT0REaVOyJN49vBJFHt7G/Ww5MM2vc+x6IiLKOA8fPkSjRo2wceNGREVFsWqJ6BUM7ImIKG18fGD38VCM2TpHPVxQtTWujBjHbvhERBno1q1bqF27Ni5cuKBa62XKZyKihBjYExFR2k2ahAZ/TEf9bC8QaWWDL72aszaJiDKITOPcuXNnPH36FHv37kX16tVZt0SUKAb2RESUPj4+GNO/MWysLLD9/D1sOxfEGiUiyoCg3sLCAnPnzsW+ffvg5eXFOiWiJDGwJyKidCuaOxt61yyi7n+59iwioqJZq0REaSRj6WVM/bNnz+Dt7Y18+fKxLokoWQzsiYgoQwxsUBy5stnh6v1nWLD3GmuViCgN/vjjD5Xx3tHRUbXYExGlBAN7IiLKEM72NhjZ3Fvdn7b1IoJDw1mzRESp8NNPP+Htt99Gt27dsGLFChXcExGlBAN7IiLKMO0q5UOFAtnx7EU0vtlwjjVLRJRCBw8exJAhQzB8+HDMnz8f1tbWrDsiSjEG9kRElGEsLS0wvnUZdX/50ds4euMRa5dS5dGjR+jevTtcXV3VIvcfP36c7GskY/jAgQORP39+ODg4oFSpUpg5cyZrnjQhJiZG3UrG+/3792PSpEnsgk9EqcbAnoiIMlTFAtnRoUp+dX/86tOIidGxhinFunbtioCAAGzYsEEtcl+C++R89NFHattFixbh7Nmz6vGgQYOwatUq1jwZtefPn6Ndu3aYPn26euzr62voIhGRRjGwJyKiDDe8mRey2Vnj+K0Q/H30FmuYUkSCcgnQZXqvGjVqqGXOnDlYu3Ytzp8/n+TrpJXznXfeQb169VC4cGH069cPFSpUwOHDh1nzZLSkp0nLli2xadMmFC1a1NDFISKNY2BPREQZzt3ZHoMbFlf3J284h9DwSNYyvZYE6NL93sfHJ3adtGDKOpnHOym1atXC6tWrcfv2bTX39/bt23HhwgU0bdqUtU5G6c6dOxg9ejROnz6NrVu3okWLFoYuEhFpHLNyEBFRpujpVwR/HryJK/efqSz5o1uWZk1TsgIDA+Hu7v7KelknzyVl6tSp6Nu3rxpjLwnHLC0tVau/BPxJiYiIUIteaGiouo2MjFRLXPrHCddTfKynlJMEedJiv3nzZpQvX57frWTwe5UyrCfTrKvUlJGBPRERZQpba0uMeaM0ei04pOa171StIIq7Z2Ntm6Fx48Zh/PjxyW5z6NAhdZvYvN3SCp/cfN4S2B84cEC12hcqVAi7du3CgAED4OnpiUaNGiX6mokTJyZaJukWndQUYxKE0euxnpIWHR0NKysrtG7dGo0bN8atW7fUQvxeZRT+/zOtugoLC0vxtgzsiYgo09T3ckdDb3dsPReML9aewW+9qjHbsxmSjPWdO3dOdhsZG3/ixAkEBQW98ty9e/eQJ0+eJJOPffrpp2rObxmvLKQFVJLuTZkyJcnAftSoURg6dGi8FvsCBQqgSZMmcHFxeaXFRE4AJRCzsbFJ0T6bI9ZT8qTL/Ycffoj169eri078TvF7xf9/hhGpod90fW+ylGBgT0REmeqzVqWx6+I97LpwD1vPBqNR6cQDNDJduXLlUsvrSLK8kJAQNZ+3TP0l/P391To/P79EX6PvOi/d7+OSVlH9NGKJsbOzU0tCcpKX1Ilecs8R6yk5y5YtQ7du3dCwYUN4eHjEfo/4nUo51hXryRy/UzapKB+T5xERUaYqkssJ79Z6mfH5y3VnEBEVzRqnRMn8882aNVPj5aVrvSxyv1WrVvDy8ordztvbW7XQC2ldr1u3LoYNG4YdO3bg6tWr+PXXX/H777+jbdu2rGkyuJ9//ln1WOnYsaMaLuLk5GToIhGRCWJgT0REmW5gg+Jwd7bD9QdhmLfnKmuckrR48WKUK1dOdYmXRbrVL1y4MN42MvWdtOLr/fnnn6hWrRrefvttlC5dGt988w2+/vprvPfee6xpMqgbN27gk08+wZAhQ9TFJmNvHSQi7WJXfCIiynQyp/3I5t4Yuuw4pm+7hPaV8yOPiz1rnl7h5uaGRYsWJVszkkwvLunavGDBAtYmGVWSPPmeFixYUOWOKF68OPOLEFGmYos9ERFliTcr5kOlgtkR9iIa3/x7jrVORCYpPDwcnTp1wgcffKAelyhRgkE9EWU6BvZERJQlLC0tMO6NMpBZy1Ycu40j1x+y5onIpMgQkebNm2PdunWxszQQEWUFBvZERJRlKhTIjreq5Ff3x60+g5iY+F2qiYi0KjAwEPXq1cOxY8ewadMmNVc9EVFWYWBPRERZalhTbzjbWePk7RD8deQma5+ITML06dMRFBSE3bt3o3bt2oYuDhGZGQb2RESUpXI72+HDRiXU/ckbziPkeSQ/ASLSrLCwMHU7btw4HD58WM3qQESU1RjYExFRlutRozCK5XbCg2cvMHXrRX4CRKRJO3bsQJEiReDv7w9ra2vkzZvX0EUiIjPFwJ6IiLKcrbUlPn+jjLr/275ruBT8hJ8CEWnK8uXL0bRpU5QvXx6lS5c2dHGIyMwxsCciIoOoWzI3GpXKg6gYHcavOfPK3ORERMZq9uzZeOutt9C2bVusXbsWzs7Ohi4SEZk5BvZERGQwY1qVgq2VJXZfvI/NZ4L4SRCR0Xvy5Am++OILDBgwAH/88Qfs7OwMXSQiIgb2RERkOIVyOqFP7SLq/lfrziI8MpofBxEZpZiYGBXUS+v80aNHMXXqVFhaso2MiIwDf42IiMigPqhfHHlc7HDjYRjm7bnKT4OIjE5ERAS6du2Kli1bqgDf3d0dFhYWhi4WEVEsBvZERGRQTnbWGNW8lLr/8/ZLCAwJ5ydCREZDWulbtWqFlStXYsiQIWylJyKjxMCeiIgMrk3FvKhSKAfCXkTjm3/PGro4RETKvXv30KBBAxw8eBAbNmxAu3btWDNEZJQY2BMRkcFJl9Zxb5SB9GxdGXAHh689NHSRiIhUMH/r1i3s3LkT9erVY40QkdFiYE9EREahXH5XdKpaQN0ft+Y0omM4/R0RZT75rdl/+QFWBdxWt/JYWupF9+7dcfbsWVSsWJEfBREZNWtDF4CIiEjvk6ZeWHfyLk7dDsWywzfRpXpBVg4RZZoNp+5i/JozuBsnt4fTw4u48ec4zJk9C507d0b27Nn5CRCR0dNUi/2jR4/UlVNXV1e1yP3Hjx8n+5qePXuqLp5xF19f3ywrMxERpVyubHYY0qikuv/txvMIeR7J6iOiTAvq3190NF5QH3bRH2cXjEBUjkKwKVyZNU9EmqGpwF6mGQkICFDjnWSR+xLcv06zZs1w9+7d2GX9+vVZUl4iIkq9HjUKobh7Njx89gI/brnAKiSiDCfd7aWlPu6An6cnNuHeiq/hULQq3N8ahynbb3JIEBFphmYCexnfJMH83LlzUaNGDbXMmTMHa9euxfnz55N9rZ2dHTw8PGIXNze3LCs3ERGljo2VJca+UVrd/33/dVwMesIqJKIMdfDqw3gt9bqYaDw9sRnZKjRFrjYjYGFtq56X7YiItEAzgf3+/ftV93sfH5/YddKlXtbt27cv2dfu2LED7u7uKFmyJPr27Yvg4OAsKDEREaVV7RK50bh0nv+3qumYSI+IMk7wk5dBvU4Xg6gn92FhaQX3jl/CrckAdT/hdkRExk4zyfMCAwNVcJ6QrJPnktK8eXO89dZbKFSoEK5evYoxY8ao+UiPHDmiWvITExERoRa90NBQdRsZGamW9NC/Pr3voxXcX9PFz9a0GcPnO7JpCey8cA97Lt3HvyfuoHHpV48BprKvaaG18hIZE3dne+iiI3F//Y+IuHUGefvMgqWtfaLbERFpgcED+3HjxmH8+PHJbnPo0CF1K4nvEpJWnMTW63Xq1Cn2ftmyZVG1alUV5K9btw7t2rVL9DUTJ05MtEybNm2Co6MjMsLmzZthTri/poufrWkz9OdbN48lNt+2xJjlxxB2JRo2lqa7r6kVFhZm6CIQaVbp3LYIWfU1wi4HIFerj2FpE7+xR84sPVztUb0Ih28SkTYYPLAfOHCgmkokOYULF8aJEycQFBT0ynMyz2iePHlS/Pc8PT1VYH/x4sUktxk1ahSGDh0ar8W+QIECaNKkCVxcXJDeFhY5eWzcuDFsbGxg6ri/poufrWkzls+3bkQUmk7di6DQCNwOccYAq9tA8eJA1aomt6+ppe9NRkSpc//+fbRs2RIvbp9FnrfGwaFwxXhJ9PTNRZLrw8oy6cYjIiJjYvDAPleuXGp5HUmWFxISgoMHD6J69epqnb+/v1rn5+eX4r/34MED3Lx5UwX4SZEu+ol105cTvow66cvI99IC7q/p4mdr2gz9+Wa3scGnLUrhwz8DMOvEI3Sc+xE8nzwAhg8HJk0yqX1NLS2VlciYSNJlmSVp966duGeX95V57KWlXoL6ZmWTPlckIjI2Bg/sU6pUqVJq2jpJfvfLL7+odf369UOrVq3g5eUVu523t7fqSt+2bVs8ffpUdfVv3769CuSvXbuGTz/9VF1IkOeJiMj4tQ6/iYW3TuNw/jKYWK8Xpq6ZAkyeDMhwqjgJVYmIkiO5lgoWLIiaNWuqnpv6RpzGpT1U9ntJlCdj6qX7PVvqiUhrNJMVXyxevBjlypVTXeJlKV++PBYuXPjKVVhpxRdWVlY4efIk2rRpozLiv/POO+pWMuw7OzsbaC+IiCg1LC5exLjNv8BCF4PVpevhVJ5iL5+4wDnuiShlZAalKlWqqMYfEbdnpgTxNYrlRJuK+dQtg3oi0iLNtNgLmX9+0aJFyW4Td0okBwcHbNy4MQtKRkREmaZkSZQNvoI2Z3ZiZZn6mFT3HSxc9rlaT0T0OpIwWWZIqlatmsrtRERkijTVYk9ERGZIutsPH46huxfBOjoKu4tUxr4RE9gNn4he6/fff1c9N5s2bYoNGzYge/bsrDUiMkkM7ImIyPhNmoSCG1eha+6olw8L1InXQ4uIKDE7d+5Ez5498ddff6menEREpoqBPRERaYOPDwa93wqOtlY4fvMxNp5+dQpUIiK56Hfq1ClVEbNnz8acOXNgba2p0adERKnGwJ6IiDQjt7Md3q1VRN2fsuk8oqJjDF0kIjIikZGR6NWrF3x8fBAYGKgSKVtYcC56IjJ9DOyJiEhT+tYpiuyONrgU/BTLj902dHGIyEiEhYWp6YxlFiVppffw8DB0kYiIsgwDeyIi0hQXext8UK+4uv/j5gsIj4w2dJGIyMAePnyIxo0bY8eOHVi7di26du1q6CIREWUpBvZERKQ53WsUgqerPe6EhGPRgeuGLg4RGVhoaCiePn2Kbdu2qQz4RETmhoE9ERFpjr2NFT5q9HIe+5+3X0JoeKShi0REBnDhwgU8fvwYhQsXxrFjx1C9enV+DkRklhjYExGRJrWrnA/FcjvhUVgk5u66YujiEFEW8/f3h5+fHz755BP12NKSp7VEZL74C0hERJpkbWWJYU291P25e67i3pMIQxeJiLLIxo0b0aBBA3h5eWHy5MmsdyIyewzsiYhIs5qW8UCF/K4IexGN6dsuGro4RJQFlixZglatWqF+/frYvHkz3NzcWO9EZPYY2BMRkWbJ/NQjmnmr+38cvIEbD8IMXSQiymT3799Ht27dsGLFCjg6OrK+iYjYFZ+IiLTOr3gu1C6RC5HROvyw5YKhi0NEmUCn06mp7MSgQYMwf/582NjYsK6JiP7DFnsiItK84U1fttqvDLiNs3dDDV0cIspAUVFR6Nu3r+p6f/LkydjeOkRE9H8M7ImISPPK5XdFy/Ke0OmAKRvPG7o4RJRBnj9/jg4dOuDXX3/Fb7/9hnLlyrFuiYgSwcCeiIhMwseNS8LK0gJbzwXj0LWHhi4OEaVTSEgImjZtik2bNmH16tXo0aMH65SIKAkM7ImIyCQUzZ0NHasWUPcn/XtOjcklIu2ysrKCk5MTtm7dihYtWhi6OERERo2BPRERmYwPG5aAnbUlDl9/hG3ngg1dHCJKg4sXL+LcuXPIli0b/v33X9SoUYP1SET0GgzsiYjIZHi42qNXzSLq/uQN5xEdw1Z7Ii05cuQIatasiQ8//NDQRSEi0hQG9kREZFLer1sMLvbWOB/0BKuP3zZ0cYgohaTLfb169VC0aFEsXryY9UZElAoM7ImIyKS4OtrgvXrF1P3vNl3Ai6gYQxeJiF7j77//VuPopbVeAvxcuXKxzoiIUoGBPRERmZxefkXg7myHW4+eY8nBG4YuDhG9hpubG7p27aqy30vCPCIiSh0G9kREZHIcbK0wuGEJdX/atot4FhFl6CIRmTXJd7H/8gOsCritbuWxzFyxbNkyREdHo0GDBliwYAFsbW0NXVQiIk2yNnQBiIiIMkOnagUwd/cVXHsQhnl7rsYG+kSUtTacuovxa87gbkh47DoPZxvkPLEY65f9jo0bN6JJkyb8WIiI0oEt9kREZJJsrCzxcRMvdX/2rit4+OyFoYtEZJZB/fuLjsYL6nVRL3Dy9/FY/9ciDPniOwb1REQZgIE9ERGZrJblPFEmrwueRkRhxvZLhi4OkVmR7vbSUh930smYyHAE/TUOz68chnvb0ThgXYHTUhIRZQAG9kREZLIsLS0wvJm3uv/7geu4/fi5oYtEZDYOXn0Yr6VeWFjbwjZXQbh3/AIOJXzU87IdERGlDwN7IiIyaXVK5IJvUTc17d1PWy4YujhEZiP4yf+D+sjHgXh+9RgsLCzh1vg92Bcom+h2RESUNgzsiYjIpFlY/L/V/u8jt3Ax6Imhi0RkFtyd7dXti6ArCFz0CR7tmA9dTHSS2xERUdoxsCciIpNXuWAONC2TBzE6YMqm84YuDpFZqF7EDY4PziHwj5Gwds6FPB2/hIWlVezzFgA8Xe3VdkRElD4M7ImIyCx80sQLlhbAxtNBOHbjkaGLQ0n4+uuv4efnB0dHR2TPnj1F9STzoY8bNw558+aFg4MD6tWrh9OnT7OODWzjhn9x6ffRsPMsCY/OE2DllD1eUC/GvlEaVvIfk4iI0oWBPRERmYUSeZzRvnJ+dX/ShnMqGCTj8+LFC7z11lt4//33U/yayZMn4/vvv8f06dNx6NAheHh4oHHjxnjyhMMuDMnb2xv9+vbBshUrkdc9fqu8h6s9ZnarjGZlPQ1WPiIiU2Jt6AIQERFllSGNS2JVwB0cuPIQuy/eR52SuVn5Rmb8+PHq9tdff03R9nKB5scff8To0aPRrl07te63335Dnjx58Mcff6B///6ZWl569fOYO3cuunbtiqJFi6qLLaJFhYIq+70kypMx9dL9ni31REQZhy32RERkNvJld0D3GoXU/ckbzyFGBt2Tpl29ehWBgYFo0qRJ7Do7OzvUrVsX+/btM2jZzE1MTAzmzJmDAQMGYPXq1fGekyC+RrGcaFMxn7plUE9ElLHYYk9ERGblg/rFsfTQTZy6HYp1J++iWWm22muZBPVCWujjksfXr19P8nURERFq0QsNDVW3kZGRaolL/zjheopfn7169cKGDRswbdo0dO7cmfWVDH6nUo51xXoy5+9UZCrKyMCeiIjMipuTLfrVKYrvN1/Ad5vOo6FXTkMXyeRJYjt9F/ukyNj4qlWrpmtaw4RdwhOui2vixImJlmnTpk0qcV9iNm/enObymXpL/RdffKESFg4bNgwFChTA+vXrDV0sTeB3inXF7xT//yUnLCwMKcXAnoiIzM67tYrgt33XcO1BGP4+ehuuhi6QiRs4cKBqwU1O4cKF0/TekihP33Lv6fn/RGzBwcGvtOLHNWrUKAwdOjRei70EpNKl38XF5ZUWEwnAJCGfjY1Nmspp6m7evImSJUuqlnvW0+vxO5VyrCvWkzl/p0L/602WEgzsiYjI7DjZWWNQg+IYt+YMpm+/gmGlDF0i05YrVy61ZIYiRYqo4F5O0ipVqhSbWX/nzp2YNGlSkq+TcfiyJCQneUmd6CX3nDm6du2a6uHQr18/DB48WJ0sS0s96ynlWFesq4zG75Rp1VVqysfkeUREZJa6+BRE/hwOCH4SgV2BnEfbWNy4cQMBAQHqNjo6Wt2X5enTp/GmUVuxYoW6L93thwwZggkTJqh1p06dQs+ePVV3esnMTpnjxIkT8PPzU1MNPnv2jNVMRGRgbLEnIiKzZGdthaGNS2LosuPYctsSIc8jkcvIr9ybg88//1xNV6enb4Xfvn076tWrp+6fP38eISEhsdsMHz4cz58/V9nYHz16BB8fH9WS7OzsbIA9MH27d+/GG2+8oXpLSLI8JycnQxeJiMjsaarFXg7W3bt3h6urq1rk/uPHj1/7urNnz6J169bqNXKQ9/X1VS0BRERk3mTqLa882fA82gKzd181dHHov/nrJfFdwkUf1At5LK3yetJqLwn67t69i/DwcNUNv2zZsqzPTArqJQ+BXHDZsWNHsnkMiIgo62gqsJcuddIdT64OyyL3JbhPzuXLl1GrVi3VbU8OQMePH8eYMWNgb2+fZeUmIiLjJHNpD21cQt3/bf8NBO7YByxcCPj7G7poREapYsWKKungv//+qxpMiIjIOGimK760ukswf+DAAdXFTsyZMwc1atRQXfK8vLwSfd3o0aPRokULNQZMr2jRollWbiIiMm71S+ZCEWcdrj6JwdSJf2DCpp9fPjF8OJBM8jUicyE9JH766Sc0b95cnW99/fXXhi4SERFpNbDfv3+/ujKsD+qFdKmXdfv27Us0sJd5VdetW6fG3jVt2hTHjh1T48Fkips333wzyb8lU7XIknCaAcn2Kkt66F+f3vfRCu6v6eJna9rM6fONiorCGwWjMfW0NZZWaIKeJ/9Fkcd3gWnTADlWpGNu9cxkDp8NGZ6cS3388cf48ccf1ZCHpBpSiIjIsDQT2Mv8tO7u7q+sl3XyXGJkDlvJovvNN9/gq6++UtPeSKt/u3btVBKeunXrJvq6iRMnYvz48a+sl0Q8kmU3I8i0POaE+2u6+NmaNnP5fIu5AKWzx+DMYyuMGvojepSIeflEcDCwfj2MUVhYmKGLQCZOpg3s3bs3/vjjD/z8888qOSERERkngwf2kuwmsSA6rkOHDqlbuVKcWPewxNbrrzKLNm3a4KOPPoodGyYt/LNmzUoysJcWfRk/FrfFvkCBAipZjIuLC9LbwiInyo0bNzb6eRMzAvfXdPGzNW3m9Pnq9/WL+WPRod2XOHoP+OqHT1D80W1g61ajbbHX9yYjyizdunXDypUr8eeff6Jjx46saCIiI2bwwH7gwIHo3LlzstsULlxYzZcaFBT0ynP37t1LMiNrrly5YG1tjdKlS8dbX6pUKezZsyfJv2dnZ6eWhOTkNqNOcDPyvbSA+2u6+NmaNnP6fCt0aIamZ/dho5cfZlRpi2mlLYEaNWCszOVzIcPp27cv+vXrh0aNGvFjICIycgYP7CX4luV1JEmezFl78OBBVK9eXa3z9/dX6/z8/BJ9ja2tLapVq6aS68V14cIFFCpUKIP2gIiITML48Riy/SA2bnuMtaXrYuCQOuBoYjI3Mh2wdLuXYYnSY4eIiLRBM9PdSSt7s2bN1NVjyYwvi9xv1apVvEQuMq3dihUrYh8PGzYMS5cuVRn0L126hOnTp2PNmjUcJ0ZERK8ea5rURItyHtAB+GnrBdYQmaToGB32X36AVQG31a08FmfOnEHNmjWxbNmyRHtJEhGR8TJ4i31qLF68GIMHD1Zj3UXr1q1VoB6XtM5LK75e27Zt1Xh6ufIsr5WLAP/884+a256IiCihDxuWxL+nArH+ZCDO3AlF6bzpy61CZEw2nLqL8WvO4G5IeOw6T1d7dCwQhi8Hv4P8+fOrRMOenp4GLScREZlwYO/m5oZFixYlu40k00tIMrrKQkRE9DpeHs5oVT4v1hy/gx+3XMDsHsaZPI8oLUH9+4uOqh4pcd24dB4ffz4UZcpXxK6tG5A9e3ZWLhGRxmimKz4REVFW+bBhcciEK5vOBOHU7f/3AiPSKuluLy31rzZ/ANa5CiBH3Xfg2OZzOLu4GqB0RESUXgzsiYiIEiju7ow2FfKq+9JqT6R1B68+jNf9XoQeXIHnV47AwsISzlVbIyhMp7YjIiLtYWBPRESUiMENS8DSAthyNhjHbz5mHZGmBT8Jjzds8dH2+Xi0fR4i7pxPcjsiItIOBvZERESJKJo7G9pWyq/u/8BWe9I4d2d7dauLicaD9T8h9OBy5GjYF9lrdU10OyIi0hYG9kREREkY3LA4rCwtsOP8PRy5/oj1RJpVvYibyn7/aMsveHZmO3K98QlcqraJfd7iv+z4sh0REWkPA3siIqIkFMrphPaV86n7HGtPWiYXqMa+UVoF83k6jIVT6Xrxgnohz8t2RESkPQzsiYiIkjGoQQlYW1pg98X7OHSNicVIe27fvo1u3brBr6AT5g5+A0Uq+sV73sPVHjO7VUazspy7nohIqzQ1jz0REVFWK+DmiLeqFsCSgzfww+YL+KOvLz8E0ozz58+jSZMmKmFecHAwmpUtjsalPVT2e0mUJ2Pqpfs9W+qJiLSNLfZERESvMbBBcdhYWWDf5QfYf/kB64s04dChQ6hVqxayZcuGffv2oXjx4mq9BPE1iuVEm4r51C2DeiIi7WNgT0RE9Br5sjugU7UCsRnypfWTyJgFBQWhQYMGKFGiBHbv3o38+V/O8EBERKaJgT0REVEKfFC/OGytLFUXZrbak7HLkycPfvvtN2zZsgVubsx0T0Rk6hjYExERpYCnqwO6+hRU97/fzFZ7Mk5Tp07Fd999p+63a9cOjo6Ohi4SERFlAQb2REREKfR+vWKws7bE4euPVJZ8ImMhw0NGjx6NDz/8UHXDJyIi88LAnoiIKIXyuNijm28hdZ+t9mQsoqKi0K9fP0yYMAHffvstJk+ebOgiERFRFmNgT0RElArv1S0GextLBNx8jB3n77HuyOAmTpyIBQsWqDH1n3zyiaGLQ0REBsDAnoiIKBVyO9uhR43C6j4z5JMxkO73mzZtQo8ePQxdFCIiMhAG9kRERKnUv05RONpa4cStEGw9G8z6o0wVHaNTMzGsCritbuXx3bt30aRJE1y8eBEuLi5qajsiIjJf1oYuABERkdbkzGaHd/wKY+aOy2qsfcNS7rCwsDB0scgEbTh1F+PXnMHdkPDYddkj7yNo6eewQjRevHhh0PIREZFxYIs9ERFRGvSrXRROtlY4czcUG08zCzllTlD//qKj8YL6iMBLODnrQ9x/HoMJC1agTJkyrHoiImJgT0RElBY5nGzRq2YRdf/HLRcQE6NjRVKGke720lIf91sVExmBe3+Ph7WrBzzfnoSZh0PVdkRERGyxJyIiSqM+tYvA2c4a5wKfYMPpQNYjZZiDVx/Ga6mXeeotbeyQu+1o5On8FSwdXdXzsh0REREDeyIiojTK7miL3rVettr/sPkCW08pwwQ/+X9Q/+TYejxY9z10uhjY5fOGpa1DotsREZH5YmBPRESUDhLYu9hb42LwU6w7eZd1SRnC3dletdI/3rMYDzfNgKWDS5LbERERMbAnIiJKB1cHG/SpXVTd/2kLW+0pY1Qp6IrwHb8gZO8SZK/bEzka9IGFxf9P22QOBk9Xe1Qv4sYqJyIiBvZERETp1atmYRXgX773DGuO32GFUrotXrQQ9w+vR87mg5Hdt0O86RT198a+URpWlpxmkYiIGNgTERGlm7O9DfrV+a/VfutFREXHsFYpTWJiXn533nnnHezfvx+LJo+Eh2v87vbyeGa3ymhW1pO1TEREivXLGyIiIkqPd/wKY96eq7h6/xlWBtxBhyr5WaGUKkFBQWjTpg2+/PJLNG7cGNWrV1frG5f2UNnvJVGejKmX7vdsqSciorg4xp6IiCgDZLOzRv//Wu2nbr2ISLbaUypcuXIFNWvWxI0bN5AnT554z0kQX6NYTrSpmE/dMqgnIqKEGNgTERFlkO41CiFXNlvceBiG5UdvsV4pRQICAuDn5wdLS0vs27cP5cuXZ80REVGqMLAnIiLKII621nivbjF1f9q2S3gRxbH29H/RMTrsv/wAqwJuq1t5LGPqe/bsifz582PPnj0oXLgwq4yIiFKNY+yJiIgy0Ns+hfDLriu49eg5/j5yC119CrJ+CRtO3cX4NWdwNyQ8tjY8stlg3JvlsHz5cuTOnRvOzs6sKSIiShO22BMREWUgB1srDKj3stV++raLiIiKZv2aOQnq3190NF5Q/yRgA47NGIz+8/fhQpgDg3oiIkoXBvZEREQZrEv1gsjjYoc7IeFYdugm69eMSXd7aanX/fdYp9Ph8b4/8XDjdNh5loCFja16XrYjIiJKKwb2REREGczexgof1C+u7k/ffgnhkWy1N1cyTZ2+pV6ni8GjLbMQsnsRXGu9jRyN+gMWlup52Y6IiCitGNgTERFlgk7VCsDT1R5BoRH48+AN1rGZkrnn9cJvnMSTY//CrekHyF6zCywsLBLdjoiIKLUY2BMREWUCO2srDGzwstX+5x2X2Wpvptnv7z+JgC46Uq13KFQBefvMhHPF5q+8xt3Z3gAlJSIiU8Gs+ERERJnkrSoFMGP7Zdx+/ByLDlxHn9pFWddmlv0+OiwEwX+Ng1PZBnCp8gZs3PLF217a7D1c7VG9iJuBSkxERKaALfZERESZxNbaEoP+a7WftfMywvYeABYuBPz9WedmkP0+KiQIgYuHI+rJPdjnL/PK9vqO+GPfKA0ry/93yyciIkotBvZERESZqH2V/Cjg5oD7T19g4dBvgR49AF9fYMQI1rsJZ79/ce8aAhcNA2Ki4fH2t7DNUxQJY3dpqZ/ZrTKalfU0RJGJiMiEaCqwf/ToEbp37w5XV1e1yP3Hjx8n+xpJTJPY8u2332ZZuYmIyHzZWFlicGErdf8Xn/Z4ZvPfWOrJk9lyb6LZ78Xj3Ytg6ZhdBfU2OV4G7jKj3ZiWpfBT54pY0tcXe0Y0YFBPRETmN8a+a9euuHXrFjZs2KAe9+vXTwX3a9asSfI1d+/ejff433//xbvvvov27dtnenmJiIhE2/AbmH3/CareOoMXVjZwivwvALxwAfDxYSWZgIRZ7XO1/EjdWto5xV/vbIc2FeOPsyciIjKbwP7s2bMqoD9w4AB8/jsJmjNnDmrUqIHz58/Dy8sr0dd5eHjEe7xq1SrUr18fRYsygREREWUNa6+SWNerFmxjouI/UbIkPwITkTCrfcKAPqntiIiIzKor/v79+1X3e31QL3x9fdW6ffv2peg9goKCsG7dOtViT0RElGV8fGD7ydD462SMPVvrTYZktfd0tY9NiJeQrJfnmf2eiIjMusU+MDAQ7u7ur6yXdfJcSvz2229wdnZGu3btkt0uIiJCLXqhoaHqNjIyUi3poX99et9HK7i/poufrWkzp883y/b1q6+AN98ELl0CihcHqlaVP5rmtzOHz0ZLJKu9ZLeXrPgSxOuT6AlmvyciIpMP7MeNG4fx48cnu82hQ4fUrSS9S0in0yW6PjHz58/H22+/DXv75LvBTZw4MdEybdq0CY6OjsgImzdvhjnh/pouframzZw+3yzbVxcXIDgYWL8+XW8TFhaWYUWijCHZ7SXLfdx57PXZ7yXoZ/Z7IiIy2cB+4MCB6Ny5c7LbFC5cGCdOnFBd6RO6d+8e8uTJ89q/s3v3bjUWf+nSpa/ddtSoURg6dGi8FvsCBQqgSZMmcJETsnSQFhY5eWzcuDFsbGxg6ri/poufrWkzp89Xq/uq701GxkWC98alPVSWfEmoJ2Pqpfs956knIiKTDuxz5cqllteRJHkhISE4ePAgqlevrtb5+/urdX5+fq99/bx581ClShVUqFDhtdva2dmpJSE54cuok76MfC8t4P6aLn62ps2cPl+t7auWympuJIivUSynoYtBRERmRDPJ80qVKoVmzZqhb9++KjO+LHK/VatW8TLie3t7Y8WKFa+0avz111/o06ePAUpOREREKfX111+rC/Yy9C179uwp6nExYsQIlCtXDk5OTsibNy969OiBO3fusNKJiMhsaCawF4sXL1YHbukSL0v58uWxcOHCeNtId3tpxY/rzz//VGPxu3TpksUlJiIiotR48eIF3nrrLbz//vspzjVw9OhRjBkzRt0uX74cFy5cQOvWrVnxRERkNgzeFT813NzcsGjRomS3kQA+oX79+qmFiIiIjJs+ee2vv/6aou1l2tuEiQ+nTZumhu3duHEDBQsWzJRyEhERGRNNtdgTERERvY703JMZc1LSlZ+IiMgUaKrFnoiIiCg54eHhGDlyJLp27ZrsTDYRERFqSTjLgIzZlyUu/eOE6yk+1lPKsa5YVxmN3ynTrKvUlJGBPREREWWqcePGxXaxT8qhQ4dQtWrVdJ8AyRS6MTExmDFjRrLbTpw4MdEybdq0SSXuS0zCLv+UONZTyrGuWFcZjd8p06orySOTUgzsiYiIKFMNHDhQBdzJKVy4cLqD+o4dO+Lq1avYtm1bsq31YtSoURg6dGi8FvsCBQqo5LwJXyvvLSeAjRs35jSDr/kMWE8p/76yrlhXGYnfKdOsK31vspRgYE9ERESZKleuXGrJLPqg/uLFi9i+fTty5nz9HPJ2dnZqSUhO8pI60UvuOWI9pQW/U6yrjMbvlGnVVWrKx+R5REREZDQkk31AQIC6jY6OVvdlefr0aew23t7eWLFihbofFRWFDh064PDhw2paXHlNYGCgWmTqPCIiInPAFvsU0E+hl5quEMm1KshYCXkvY79ClBG4v6aLn61pM6fPV6v7qj8mJTbNq5Z9/vnn+O2332IfV6pUSd1KS3y9evXU/fPnz6vM9+LWrVtYvXq1ul+xYsV47xX3Nek51mv1O5LVWE+sK36v+P9PCyI19JuemmO9hc7UzggygZw0yLg7IiIiY3Pz5k3kz5/f0MXQPB7riYhIy8d6BvYpINl179y5A2dnZzUvbnrok/PIh/O6xD6mgPtruvjZmjZz+ny1uq9yXf7JkyfImzcvLC05si4zj/Va/Y5kNdYT64rfK/7/04JQDf2mp+ZYz674KSCVmNGtIfIlMvYvUkbi/pouframzZw+Xy3uq6urq6GLYFbHei1+RwyB9cS64veK//+0wEUjv+kpPdbzEj8RERERERGRhjGwJyIiIiIiItIwBvZZTObMHTt2bKJz55oi7q/p4mdr2szp8zWnfaW04XeE9ZTR+J1iXfE7ZTh2JnrcZ/I8IiIiIiIiIg1jiz0RERERERGRhjGwJyIiIiIiItIwBvZEREREREREGsbAnoiIiIiIiEjDGNhnghkzZqBIkSKwt7dHlSpVsHv37mS337lzp9pOti9atChmzZoFU93fu3fvomvXrvDy8oKlpSWGDBkCU93X5cuXo3HjxsidOzdcXFxQo0YNbNy4Eaa6v3v27EHNmjWRM2dOODg4wNvbGz/88ANM+f+u3t69e2FtbY2KFSvCVPd3x44dsLCweGU5d+4cTPGzjYiIwOjRo1GoUCGVNbdYsWKYP39+lpWXDOvrr7+Gn58fHB0dkT179tduHxkZiREjRqBcuXJwcnJC3rx50aNHD9y5cwemLrV1JXQ6HcaNG6fqSY4X9erVw+nTp2HKHj16hO7du8PV1VUtcv/x48fJvubp06cYOHAg8ufPr+qpVKlSmDlzJkxdWupKnD17Fq1bt1avcXZ2hq+vL27cuAFTlta60uvfv786lv/4448wZY9SWU+a/U3XUYb6888/dTY2Nro5c+bozpw5o/vwww91Tk5OuuvXrye6/ZUrV3SOjo5qO9leXiev//vvv01yf69evaobPHiw7rffftNVrFhRba8Vqd1XeX7SpEm6gwcP6i5cuKAbNWqUev3Ro0d1pri/sl9//PGH7tSpU+pzXrhwofpu//LLLzpT3F+9x48f64oWLapr0qSJrkKFCjqtSO3+bt++XSeHjPPnz+vu3r0bu0RFRelM8bNt3bq1zsfHR7d582b1ffb399ft3bs3S8tNhvP555/rvv/+e93QoUN1rq6ur91efgcaNWqkW7p0qe7cuXO6/fv3q+9PlSpVdKYutXUlvvnmG52zs7Pun3/+0Z08eVLXqVMnnaenpy40NFRnqpo1a6YrW7asbt++fWqR+61atUr2NX369NEVK1ZM/f7K75AcT62srHQrV67UmbK01NWlS5d0bm5uumHDhqnzkcuXL+vWrl2rCwoK0pmytNSV3ooVK9R5S968eXU//PCDzpQ1S2U9afU3nYF9Bqtevbruvffei7fO29tbN3LkyES3Hz58uHo+rv79++t8fX11pri/cdWtW1dTgX169lWvdOnSuvHjx+vMZX/btm2r69atm86U91dOSD/77DPd2LFjNRXYp3Z/9YH9o0ePdFqT2n39999/VYDy4MGDLCohGasFCxakOFhNSC7qyv+Z110cNLe6iomJ0Xl4eKjgXi88PFy9dtasWTpTJBcU5btw4MCB2HUSKMg6CRqSUqZMGd0XX3wRb13lypXVMcdUpbWu5FislfMNQ9eVuHXrli5fvnyqMaZQoUImHdifSUc9ae03nV3xM9CLFy9w5MgRNGnSJN56ebxv375EX7N///5Xtm/atCkOHz6suoGY2v5qVUbsa0xMDJ48eQI3NzeYw/4eO3ZMbVu3bl2Y6v4uWLAAly9fxtixY6El6fl8K1WqBE9PTzRs2BDbt2+HKe7r6tWrUbVqVUyePBn58uVDyZIl8cknn+D58+dZVGoyBSEhIaqLa0q7p5uLq1evIjAwMN7/SRnuIscKUzt3iHuuJ91/fXx8YtdJN3FZl9w+16pVS/0e3b59Ww1fkN/cCxcuqPNEU5WWupLzq3Xr1qnfaqkbd3d39fqVK1fClKX1eyX1JV3Rhw0bhjJlysDU7U9jPWnxN52BfQa6f/8+oqOjkSdPnnjr5bEcxBIj6xPbPioqSr2fqe2vVmXEvn733Xd49uwZOnbsCFPeXxkLKCdpEhh98MEH6NOnD0xxfy9evIiRI0di8eLFany9lqRlfyWYnz17Nv755x+VP0LyZEhwv2vXLpjavl65ckXljDh16hRWrFihxh7+/fff6vtMlBLh4eHq90FyykiOFfo//f87czh30JP9kmAzIVmX3D5PnToVpUuXVsdVW1tbNGvWTOULkYDfVKWlroKDg1U+gm+++UbV0aZNm9C2bVu0a9dO5bEyVWn9Xk2aNEmdtwwePBjmIDCN9aTF33QG9plArubEJVdZE6573faJrTeV/dWytO7rkiVLVKKgpUuXJvrjYkr7K0nJpMeJJIGUgEj23dT2VwJF+XEfP368aiHQqtR8vhLI9+3bF5UrV1aJIOXksmXLlpgyZQpMbV+lNUOek4s21atXR4sWLfD999/j119/Zau9hslvcGIJIOMu8tuVXtLbrnPnzup7JP9PtCgr6soUzh1SU0+J7dvr9lkC+wMHDqhWe+l5JA0EAwYMwJYtW6A1mVlX8n9NtGnTBh999JFKZCtBWKtWrTSXkDqz60q+Rz/99JM6nmnt/1tW///T4m+6tpqZjFyuXLlgZWX1ytUfuZKY8Mq0noeHR6Lby5U0yS5uavurVenZVwnm3333Xfz1119o1KgRTH1/JfO4kEyiQUFB6oe3S5cuMKX9lSEVcrCQ4QaSsVjID74cJOT/rrQWNGjQAKb+f1e6si1atAjGLC37Kr0TpAu+dNPTk2zU8vneunULJUqUyPRyU8aT/6tycpacwoULp+tvyAmg9MqS7ubbtm0z6pYdQ9WVnPcI+T8p/9e0fO6Q0no6ceKEOh4mdO/evST3WYb+fPrpp6rXkFxEFeXLl0dAQIC6oKqV84msqCv5nZdjr/RuiEt+t6X3ldZkZl1J44v8XytYsGC8xoqPP/5YNcZcu3YNWpGZ9aTV33QG9hlIuknJNEqbN29WXYD05LFcRUyMtHytWbMm3joJCqQbs42NDUxtf7UqrfsqrdW9e/dWt/oDsxZk1GcrgZBMG2Zq+ys/7CdPnoy3Tq7iyo++dNnWX9ww9c9XLmzEPTE3lX2VaRvlQpx07cyWLZtaJ+NaZYpO6RJL2iQn/7JkFv0JoAzTkbHQxn5x3lB1Jb+PEtzL/0HJ2aHPhSFdpqWLsCnWk5zryfjcgwcPql5Awt/fX62TqQKT+j7JIr87ccmFSn0LtZZkZl3J73y1atVw/vz5eOvld1umLNWazKwrGVuf8KKQ5CWQ9b169YKWZGY9afY33dDZ+0yNflqlefPmqSyMQ4YMUdMqXbt2TT0vWZi7d+/+ynR3H330kdpeXqfF6e5Sur/i2LFjapEpI7p27arunz59Wmdq+ypTv1lbW+t+/vnneNODyRQaWpDa/Z0+fbpu9erVamo/WebPn69zcXHRjR49Wmeq3+W4tJYVP7X7KxlzZWoc+Wwli648L4cQma7K1Pb1yZMnuvz58+s6dOigfpt27typK1GihJp6isyDZD2WY5PMYpItW7bY45Z8N/S8vLx0y5cvV/cjIyPVFInyvQkICIj3mx8REaEzZamtKyEZ8SULvqyT6e66dOliFtPdlS9fXmXjlqVcuXKvTLeVsJ5k9iDJjC+zksj5osw8YG9vr5sxY4bOlKWlruS+/M7Pnj1bd/HiRd20adPU1IC7d+/WmbK01FVCpp4VPy31pNXfdAb2mUACOflPYmtrq6YlkZNCvXfeeUf9UMe1Y8cOXaVKldT2hQsX1s2cOVNnyvsrwUDCRV5vavsq9xPbV9lOK1Kzv1OnTlUnIHKhSgJ6+U7LyUd0dLTOVL/LWg7sU7u/kyZNUvMpy0lljhw5dLVq1dKtW7dOZ6qf7dmzZ9Uctg4ODurALnN0h4WFGaDkZAjynUjs91sCLD15LIGWkDnGE9s+4WtMUWrrSj/lnfxmyrR3dnZ2ujp16qgA35TJ9Jlvv/22ztnZWS1yP+H0oQnrSYKInj17qnnG5bdXAo/vvvtO1Z8pS0tdCbl4W7x4cVVXcjxeuXKlztSlta7MLbB/kMp60upvuoX8Y+heA0RERERERESUNsyKT0RERERERKRhDOyJiIiIiIiINIyBPREREREREZGGMbAnIiIiIiIi0jAG9kREREREREQaxsCeiIiIiIiISMMY2BMRERERERFpGAN7IiIiIiIiIg1jYE9ERERERESkYQzsiYiIiIgoSRs2bICDgwOioqJi1509exYWFha4f/8+a47ICDCwJyIiIiKiJAUEBKBMmTKwtraOty5fvnzIlSsXa47ICDCwJ6JMt2TJEtjb2+P27dux6/r06YPy5csjJCSEnwAREZERO378OCpWrBhv3bFjx1ChQgV1/8qVK1izZo2BSkdEgoE9EWW6zp07w8vLCxMnTlSPx48fj40bN+Lff/+Fq6srPwEiIiIjJq3z+iA+sXVyPD937lyir42Ojs6SMhKZOwb2RJTpZAze119/jblz52LChAn46aef1Hg96cJHRERExuv58+e4ePFivBb7mJgYHD16VAX2O3fuxGeffYY5c+agUqVKavvmzZtj+PDhqFOnDn7//XfVQ+/Ro0fqtXv37sU777yj7p8/fx4tWrRAlSpVUK9ePY7XJ0oHBvZElCVatWqF0qVLq9b6FStWqLF6REREZNwuX76sWt2l552e9Lp78OCBCuzr1q2LsmXLYuvWrap7viTZO3XqlLp4v2vXLnTv3h1Pnz5Fjhw51GtPnDihzgEiIiLwwQcfYPbs2Thy5Ag6dOigGgCIKG0Y2BNRlpCTAOmmJycHefLkYa0TERFpQM6cOVXPu4MHD6rHBw4cwMCBA1UAX6JECbXu1q1bKFCggLovuXNk+w8//DC2Vb5kyZKx76cP7FeuXIkzZ86oC//SG+Dnn3+GjY2NQfaRyBQwsCeiTCfd9d566y388ssvaNq0KcaMGcNaJyIi0gBPT098+eWX6NGjBwoWLIgZM2aoY7oE51ZWViqojzu0Tlrr/fz84j2WFn29w4cPq8cnT57Ed999p8bqyyLT53388cdZvn9EpuL/c1YQEWWCa9euoWXLlhg5cqTqjifd8atVq6a63cmYOiIiIjJuo0ePVktirl69irx588YL5MuVKxf7+OHDh6p1X0jXfGnBlwsEHh4eqjdfly5d1HMS6Md9HRGlDlvsif7Xzh3bMAhDURR1toFNaGAAOtZgIKZgDfZgA0d2lyJNpAg96ZySCldwZfvzN+1j3gboLMtS9n3vz1rMz/P89QcBAMjRdt/bcL0W5e3K3XVdH4E+TVOfrbOuaznPs4zj2I/qb9tW7vsuwzD0u/rHcTy6Dkj3qrXWp18CAAAA+I0dewAAAAgm7AEAACCYsAcAAIBgwh4AAACCCXsAAAAIJuwBAAAgmLAHAACAYMIeAAAAggl7AAAACCbsAQAAIJiwBwAAgGDCHgAAAEquN5Tuwl9iqkdnAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.figure(1, figsize=(12, 5))\n", "_, ax = plt.subplots(1, 2, num=1)\n", "\n", "# plot (x,y) true and (x,y) predicted\n", "ax[0].scatter(x,y, s=10, c=\"r\", label=\"Data\")\n", "ax[0].plot(x_true, y_true, label=\"Solution\")\n", "ax[0].legend()\n", "ax[0].grid()\n", "ax[0].set_xlabel(r\"$x$\")\n", "ax[0].set_ylabel(r\"$y$\")\n", "ax[0].set_title(\"Deformed rod\")\n", "\n", "# u_true vs u_predicted plot\n", "u_min = min(u_true.min(), u.min())\n", "u_max = max(u_true.max(), u.max())\n", "ax[1].scatter(u_true, u)\n", "ax[1].plot([u_min, u_max], [u_min, u_max], \"k--\", linewidth=1)\n", "ax[1].grid()\n", "ax[1].set_xlabel(r\"$u_{true}$\")\n", "ax[1].set_ylabel(r\"$u_{sol}$\")\n", "ax[1].set_title(\"True-predicted plot\")\n", "plt.show()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.11" } }, "nbformat": 4, "nbformat_minor": 5 }