{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "90700428",
   "metadata": {},
   "source": [
    "***circle_packing.ipynb***\n",
    "\n",
    "Place $N$ disjoint circles in the unit square to maximize the sum of their radii.\n",
    "\n",
    "&copy; Copyright 2025 Fair Isaac Corporation\n",
    "\n",
    "Licensed under the Apache License, Version 2.0 (the \"License\"); you may not use this file except in compliance with the License. You may obtain a copy of the License at http://www.apache.org/licenses/LICENSE-2.0.\n",
    " \n",
    "Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an \"AS IS\" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License.\n",
    "\n",
    "This example uses FICO&reg; Xpress software. By running it, you agree to the Community License terms of the [Xpress Shrinkwrap License Agreement](https://www.fico.com/en/shrinkwrap-license-agreement-fico-xpress-optimization-suite-on-premises) with respect to the FICO&reg; Xpress software. See the [licensing options](https://www.fico.com/en/fico-xpress-trial-and-licensing-options) overview for additional details and information about obtaining a paid license."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "52204f86",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Install the xpress package\n",
    "%pip install -q xpress"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "68fa499f",
   "metadata": {},
   "source": [
    "## Problem description and formulation"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "79ee9221",
   "metadata": {},
   "source": [
    "In this example, we aim at packing $N$ circles inside the unit square to maximize the sum of their radii. The circles must not overlap with each other and must be fully contained inside the unit square. Let $CIRCLES$ be the set $\\{1,...,N\\}$.\n",
    "\n",
    "The continuous decision variables $x_i$ and $y_i$ represent the vector of $(x,y)$ coordinates of the center point of each circle $i \\in CIRCLES$, and variables $r_i$ the corresponding radius.\n",
    "\n",
    "The goal is to maximize the area inside the unit square occupied by circles, that is, the sum of the radii of all circles:\n",
    "\n",
    "$$\\max \\sum_{i \\in CIRCLES} r_{i}$$\n",
    "\n",
    "The circles must not overlap, that is, the Euclidean distance between the center of any two circles must not be shorter than the sum of their radii. This can be represented by the following set of quadratic constraints:\n",
    "\n",
    "$$ (x_i - x_j)^2 + (y_i - y_j)^2  \\geq (r_i + r_j)^2, \\qquad \\forall i \\in CIRCLES, \\forall j \\in i+1,...,N$$\n",
    "\n",
    "Moreover, we need each circle to be fully contained within the unit square. Therefore, both the $x_i$ and $y_i$ coordinates of each circle must not be lower than its radius (similarly, they must not be greater than (1 - $r_i$)). This can be represented by the following linking constraints:\n",
    "\n",
    "$$  \n",
    "\\begin{array}{llll}\n",
    "& \\qquad  x_i \\geq r_i, \\qquad \\forall i \\in CIRCLES \\\\\n",
    "& \\qquad  x_i \\leq 1 - r_i, \\qquad \\forall i \\in CIRCLES \\\\\n",
    "& \\qquad  y_i \\geq r_i, \\qquad \\forall i \\in CIRCLES \\\\\n",
    "& \\qquad  y_i \\leq 1 - r_i, \\qquad \\forall i \\in CIRCLES \\\\\n",
    "\\end{array}\n",
    "$$\n",
    "\n",
    "Additionally, the radius of each circle has a trivial upper bound equal to 0.5, since the diameter of any circle can never be greater than 1:\n",
    "\n",
    "$$0 \\leq r_i \\leq 0.5, \\qquad \\forall i \\in 1..N $$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f7e9ef6f",
   "metadata": {},
   "source": [
    "## Model implementation"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "39c67d63",
   "metadata": {},
   "source": [
    "The previous model is implemented below. Note that the upper bound is defined on variable creation with [p.addVariables](https://www.fico.com/fico-xpress-optimization/docs/latest/solver/optimizer/python/HTML/problem.addVariables.html), which creates a *NumPy* array of $N$ variables when an integer is passed as the first argument."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "ff48ce99",
   "metadata": {},
   "outputs": [],
   "source": [
    "import xpress as xp\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "N = 6   # Number of circles to place.\n",
    "CIRCLES = range(N)\n",
    "\n",
    "# Create a problem instance.\n",
    "p = xp.problem()\n",
    "\n",
    "# Decision variables.\n",
    "x = p.addVariables(N, name=\"x\")         # x-coordinate of center points.\n",
    "y = p.addVariables(N, name=\"y\")         # y-coordinate of center points.\n",
    "r = p.addVariables(N, name=\"r\", ub=0.5) # Radii of circles, upper bound is 0.5 to fit in unit square.\n",
    "\n",
    "# Non-overlap constraints.\n",
    "p.addConstraint((x[i] - x[j])**2 + (y[i] - y[j])**2 >= (r[i] + r[j])**2 for i in CIRCLES for j in range(i + 1, N))\n",
    "\n",
    "# Ensure that each circle is contained in the unit square.\n",
    "p.addConstraint(x[i] >= r[i] for i in CIRCLES)\n",
    "p.addConstraint(x[i] <= 1 - r[i] for i in CIRCLES)\n",
    "p.addConstraint(y[i] >= r[i] for i in CIRCLES)\n",
    "p.addConstraint(y[i] <= 1 - r[i] for i in CIRCLES)\n",
    "\n",
    "# Objective function.\n",
    "p.setObjective(xp.Sum(r), sense=xp.maximize)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c72a2c05",
   "metadata": {},
   "source": [
    "Before triggering the optimization, a number of optimizer controls are set. For exact arithmetic, <!--we set the tolerance for matrix coefficients to zero via [MATRIXTOL](https://www.fico.com/fico-xpress-optimization/docs/latest/solver/optimizer/HTML/MATRIXTOL.html), and --> the feasibility tolerance is lowered to a very small value with [FEASTOL](https://www.fico.com/fico-xpress-optimization/docs/latest/solver/optimizer/HTML/FEASTOL.html). \n",
    "\n",
    "Moreover, we set a time limit of 3 seconds before defining the global solver (default) to solve the problem.\n",
    "\n",
    "The code below defines these control parameters and optimizes the problem before printing the solution and objective values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8b00ade4",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Control parameters.\n",
    "p.controls.feastol = 1e-9               # Set the feasibility tolerance to a very small value.\n",
    "p.controls.timelimit = 3                # Set a time limit in seconds.\n",
    "#p.controls.worklimit = 0.6             # Alternatively, set a work limit instead of time as a deterministic stopping criterion.\n",
    "p.controls.nlpsolver = xp.constants.NLPSOLVER_GLOBAL        # Set the NLP solver to global.\n",
    "\n",
    "# Solve the problem.\n",
    "p.optimize()\n",
    "\n",
    "# Print a solution summary.\n",
    "xsol = p.getSolution(x)\n",
    "ysol = p.getSolution(y)\n",
    "rsol = p.getSolution(r)\n",
    "print(f\"Sum of radii for N = {N} is {p.attributes.objval}\")\n",
    "for i in CIRCLES:\n",
    "    print(f\"{i}: x = {xsol[i]}, y = {ysol[i]}, r = {rsol[i]}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bfa4ee35",
   "metadata": {},
   "source": [
    "## Visualization"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7d5441fd",
   "metadata": {},
   "source": [
    "The code cell below generates a plot of the circles defined by the solution inside the unit square, with a blue line and label indicating the radius of each circle."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "90bb8e2b",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "N = 6, Sum of radii = 1.2028389108720672\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAbAAAAGiCAYAAACGUJO6AAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAAdARJREFUeJztnQm4TXX3x5dr1iANMieVNFLEWyqppJTmEoUmTTRQKRXSRJH6V6SU0iCaSBFJKGlEs5ei4k1ECpmH/X8+Z7Vz3e5wzrlnnz2tz/Ps57rXvffsu8/ev/Vba33XWiUcx3HEMAzDMEJGjt8nYBiGYRjpYAbMMAzDCCVmwAzDMIxQYgbMMAzDCCVmwAzDMIxQYgbMMAzDCCVmwAzDMIxQYgbMMAzDCCVmwAzDMIxQYgbMMAzDiIcBe//996V169ZSrVo1KVGihIwZM6bIn5k6daocfvjhUrZsWdl3333l2WefTfd8DcMwDCM9A7ZmzRqpX7++DBo0KKnv//HHH+XUU0+V5s2byxdffCE33HCDXH755TJx4sRUX9owDMMw/qFEcZr54oGNHj1azjzzzAK/55ZbbpFx48bJN99888/XLrjgAvnzzz9lwoQJ6b60YRiGEXNKef0CH330kZx44onbfa1ly5YJT6wgNmzYkDhctm7dKitWrJDddtstYTQNwzCM8ICftHr16kTqKScnJzwGbMmSJbLnnntu9zU+X7Vqlaxbt07Kly//r5/p27ev9OnTx+tTMwzDMLLIokWLpEaNGuExYOnQo0cP6dat2z+fr1y5UmrVqiXz5s2TXXfd1Z+TWrsWBYvIO++IvPuuyIoVIpxLw4Yihx6qxyGHiFSunP5rbN4s8sMPIl9+KfLVV9uOLVv0d7doIXLSSSIHHUT89l8/vmnTJpkyZUoi31i6dOni/b0RZtP69TLl/fel+cyZUpow9vz5+h877rjtfXQ/8rCVKZO5e+jnn7d/b+fMEVm/Xt9PXo/3l+OAA/J9j7NFrO+lrVtFvvhCZNIkfd7nztWv87y79wbHwQfLpsqVZcq0acW/Tjz77r3x9dd68O81a/T/Dzts271Rt66v90Y6EEGrW7eu7LTTTpJJPDdgVapUkaVLl273NT7feeed8/W+ALUiR14wXoQRswq5u8cfF3nuOZG//lLj0aGDyOmnizRuLFKyZGZfD2+1adNtn//xh8jbb4u8+abIY4+J3HefPkTXXCNy4YUiuW4IFp0KFSokrlHsFp2iWL1aBOHQm2/KpsmTpcKgQbLbq69KaTYFt98u0qiRyD77iGQwvPEvuHdr1hQ5+ujtF67vvhP5/HNdLB95ROSee/T7uMc4mjXjoZBsErt7iU3E3/dH4vjtN5Hddxc57TS9P5o00fckj+HI6HXi2WdNyW1I2dB+9JHIW2+JPPigCJGpvfcWad1a743jjsv8GuQhGU8BOcWAHx89enSh39O9e3fn4IMP3u5rbdu2dVq2bJn066xcuTLxWsuXL3eywqZNjjNypOMceyx/pOPsuafj9OzpON9/7/jKhg2O8/bbjnPWWY6Tk+M4O+3kONdc4zjffZf4740bNzpjxoxJfDQcx9m61XGmT3ecdu0cp0wZfS8POcTZ2LOnXqf1653AwXs8aZLjXHut4+y1l55zxYr6+d/vczaIzb00f77j3Hyz4+y2m17revVYtPS+2bw5WNeJ+3XCBH3ma9TQ861Vy3Huu89xli51ggxrN2s4a3kmSdmArV692pk9e3bi4IQGDhyY+PfPP/+c+P9bb73Vad++/T/fv2DBAqdChQrOzTff7MyZM8cZNGiQU7JkSWcCb0TQDBgL3muv6U3MzdGsmeOMGqWLStBYuFCNKsa1RAnHad/e2fjDD/FYdIpi9WrHGTLEcQ49VN/Hffd1nAce4GYM1+LM/fjll47To4fjVK6sf0vz5o7zyiv8EZ6+dGiuUTpgmN5803FOOUWfnV12cZxu3dLaIPh2nbg3PvnEcS691HHKlXOc0qXxDBzngw/0/wJGYAzYlClTEieS9+jYsWPi//nYjIU/z880aNDAKVOmjFOnTh3nmWeeSek1s2LA3nvPcRo31kXipJMc5/PPnVCAcR08OGHINlasqA/T4sVOLMFD7tJFPVM81DPOcJyJEx1ny5bwL87svkeMcJyjj9Z7tGpVx+nd23F++82TlwvlNSoK/pZBg7Z5tg0bOs6wYY6zZk24r9PvvzvOwIG6UePvYuPGxjtAhiwwBswPPDVgLPZnnqlvPAZs8mQnlPz1l7Oxb199mPbYQ41anoU7svAeXnWV45QqpZ7K7bc7zt8RgcAuOsUBr4y/d4cd1Fj36aNeZwYJ/TXKDc8Bxn+ffdTjIqT86afRu05btjjOO+84zsknbzPQhKMjbMDi2wuRt/iFF1SUQZJ01CiRjz8WOf54CSU77CBy44367/PPV5EHAoUff5TIsnKlJtj33VffPwQuP/2kIohatSSyoIBDWMTf2qmTyL33qgAFkc/GjX6fXbCecVSmKIXbtRM58EBV+L74osgRR0jkyMnRZx7R19SpIohK+JwDkVAEiacB+/VXEbqHtG8vcsopIt9+q4t+yKSpBfLQQyoBRsGEYnHwYH2Yo6QYQ5FVp47+rddfL7JggcjNN4sUoGyNJKjkuA7ffy9y6ql6HerV0wUaBVucQfp+wgn6fLO5++ADkbFj9XmIA82aicyYITJ6tMgvv6jBvuACJOASJeJnwPC2qKn45BN9c3nYsy3NzwZ0P6EE4KKLRDp3FmnTZltNSZjhoaxfnx5luunASON57bKLxBa8zWHDtG6Ia8N73rz5tvq2OEHNJEa9QQOqZtVoYbxyly7EhRIldKPOffH00yLvvade6EsvRWZDGy8D9swzWjdByIlQQiE9HCMBNWJDhoi8/rrI+PH6EFMsGUbWrRO56Sb9GypV0oeSMFq1an6fWXAgHM6mbMoUXbwJNT76aHy8MbyuY49VT/yqq/QZp14qKpGVdClVSuTSSzXSxMaWcOo550TCG4uHAaNYlM4evIkUIbMTydPeKtKcdZZ6Ln/+qaGE6dMlVHDu7KjJ8dx/v8iHH+pO0sgfNmkYeO73666LvjeW2+tatkw75hBarlDB7zMLFnvsobniV17RNSAC3lj0DRhNgc89VzscsBt98snMtQYKE+zGP/tMb1qEKuzUw/De5fa6Zs/W3XWIOg/4Bm2xuN/zemMhXqzyhU41rVrpfXH11doCKo7hwlRgPcztjV1+uT5rISTaBoywE2FClEhvvCHSpUu8wwkk/WlXhEd23nkiI0dKYFmyRD0HFl3X66I/oJGeN3bJJeqNsWDRkzEK0EeSFk8o7BAtDRxoXleq3hjDhVFjs6kNYUgxugaMHQVx3mnTtI8YKi1DvU+EKyxk9FJ89VUJHDNnaqgTmTgJePO6iu+NEX59+WUVNZAn+t//JNSQ0/3Pf/R+JrKA4tBInY4ddY1ExUs/0FmzJEzkRDYm3rat5rp4YPPMI4s9JHURtCCr5TqNGyeBAa+QEFDVqrqzzt3c1CgeeN14suSJWKxQ5IYNQqAPPKBNdvEu+RsopzDSh40Az1qVKvrs4ZmFhGgasFtv1ZDha6+Z8SoIPJrhw3UhwJAxvsFPUMpRlIxBJUbPrtAUhpkHoQMey377qQEghBSmjemVV2oJxW23aR43w+M5Ykv16ip+Ib3AesAmIQQEch5YsWDsyYABIg8/bGHDZDwx4t+MbznjDJFPP9U8WbbZtElrl1BH9e+vHUXinKv0GmbWTZ6s9YHkxujWcuedwb7mKIlRVRL+xugS+jIyC00AWA/waNkk0DCgZ89A3xfRMmAUJ19xxTb5sFE0dCnAWyXnRIgJkUc25z+Rq6TImpwG9WpRr80LCuSOUORSE0nEAsETYpkgLlYYLzY45GtHjND7xfAG3v+771ZjRkSE1mR8HsT7IlIGjAF0uL/0PaN1UkAveCDZay8Nt6JEwvuh5CAbsMPjPaNvG0aUtj9G9uAZYaddrpzIDTfoZoLIRZCeHcKGF1+s9yceOveL4T233aZDVClj4SOeWACJhgEjsUsNCKEobvQsT6+NBMcco8Wg9NMjnOi1qoudHbkujBcqUVOR+QfvOc8MzxDeN2HcoDzX5LwotuUw45VdbrxRNzV4Ytwf3btL0IiGAUM1Q/gJmTBKGiM9qJMjMX7ZZSrq8CpBTkgIsQa1O4xvN+PlP7ReYgNI6J3wUa9efp+Rdtmnhx9iI/peGv54YuvWqadOpCZg4dvwGzCK71h4yd9wGMUbx8CCQdcGdlv0GvRiV414gPIGNh0nnZT51zDS49prdcdN3R0bwRo1/DsXNlKErfr00fZvhn/cdZfWiSH4Qb16+OESFMIvo2cxJGY/aJDfZxINUCCRzKcJMG2IMg35ScQDTzyhjVaNYEHOo2tXFXb4Bd4/o45oRHDHHf6dh6Gwvj71lDaLJr0QoI4d4TZgNKQk50XjTlqjGJmBXMhRR2kD5Ex2MqewnHwLB0pRI5hQA0S3Dsj2QNTly3WRRB1J6JCogOE/5curV0z4/+yzA9M7Mbx3B6Eo4rLM9qItkpE5WDRYxGiMmqmqfEIQhHjpb0idnhHs+kC3wJlc5erV2XldcnDcI3/9papUSjyM4FCjhhoxunYw8T0AjaHDa8BQrjFmo29f26V5AcXNhPgI4RR3TD0L4Omni+y6qxpEFkgj2LgDQumZSA1WNmaKMZiU3pdEVRAMGMFsO/XkkzpANQAtp8K58lMb0qOH7uZNBOCtCowQEvHvdGGXRvJ34UIVbmDEjPDAQoVSlGJWL8Hbv+ceVb1R0mEEl44dVRWKeM7nfFg4DRgPFPNs8L6CVHQZNQ45RDvW9+unm4Z0oHMCO2qaB9s4lPDBBhFJPQaMKQFegIdPsTL3h4k2wsFjj2nki3y5j6HEcBowlGxHHqmzgAxvQZHGQMR0Otb/+qtKs8mjoCgzwgmFrGxmMDLFDSfnBxvRb77RvFsch82GkT320HWYnJiPocTwGbB587QAliSi4T3UfLBR4GZNBXZlFMfS2SFbrakMb+A9xLj897+ZDyXmDh0GqL7ISAI66fgcSgyfAaM+iY7pXDwjO7BZmDhR5IcfUgsdkvNy3y8j3NSvr4XFeEuZCiUiDGGc/YEHWugw7KHEG27w5eXDZcCoPSCXQqsjGpAa2YFdFuILio9TDR1a/7rogHCKLi2ZCiUSesIY4t1b6DC8ocR+/XQQrVc50sgYMKbJ/vmnVukb2YPNAj3QEGMkk7ClCShhp0cfzcbZGdkOJc6dq42fiwMGEK+LUg1KNozw0qGDSL16vnRvCZcBI/e1zz4acjCyC3VcSOq/+67w75s9WzuHI8HfbbdsnZ2RLfDAUJ7RbmzFivR/D6UZ3E/UfhnhplQpfR/ffVePLBI+A8ZCatL57EPN3Y47al6rqDDT/vtrmMmIriqRsgpCR+lApw0axLJzP/jgTJ+d4QcMoqXImec/i7L6cBmwJUusAaxfMA+ImqDCDBjNfxF74H1Zt43oUrmyhokJEdOpI1UYmvnHH9pp3ogGJUrohoY2U0zOzhLhMmDMpzr6aL/PIr6wefjkk/xDR+y6iIEfcYQ2+zSiDY2e8cjvvDO1n1u7VvNnhCGtXVS0aNZMpGVLLYvIkhcWLgPWoIEmkg1/oEM9N2Z+aiMKGj/9VHdhFuKNPjvvrCIMVMFz5iT/c6jVVq7UiQRGNDc2X32lfWqzQE7oalEM/2DEBV5wXgOGUWMn3qKFyPHH+3V2RrahUL1mzeSLm7lPmNvXqpXI3nt7fXaGH5x4oq4TqTY+iIUBo52N4R8ULNItIa8BYy4bQwiZ4mzEKy+KJ0XOI5lODJ99JjJrlnXRiTI5f/dHfOUVkd9+8/7lJExE1IDVrq3CPSKkfExV3MV9cvLJVBeUkuuuay4ffFCiQPEXIWoaY7jTMvIDASFRQEru/kXDhpqozQ27rbp1zfuKI9wsCHaefrro7+U+wfPiJjSifU+ULJncPRErA0a4IqLQlIC2cAwtplsP6aRkQTuBgvW77zbLtdfOlg4dSiZmA+aF9CEzQAsr1Xj99SLSjFjZn35Sa+gqQylwZtdlc9niR6VK2nGFlmGFTSz4/XfNfxF2ZHEzosuuuyZ3T2SAcK04MRAHVK+uRe0//5z8z7z8sq4LsN9+f0rVqo5Mm5Z/xAcnqSDviygQ9YgDBxYxldVtF+UWpLIDt7qv+EJIsKiJBdykLGbMhjOiz5VX6gzA/BaiDGLFOgGDht9sVo87TgcZFzTbb889teSK78XbqlJFPwLqZO6dVOnUSeSBB1SnUSBVq24zYISD6I/IzLDCYpJGtCGs3LixhghpNJAf1A8is6Z3nhF9GjfW3Tjvu4epBTNgAYFWg0TgaDP30EPbnnPCitkAR6pWrSTutdwGbPx4LWQlfGjEG7wwvPAFC0Tq1Nn+/9iJERvv39+vszP8iJaxmcGAsaB5RLhCiBGGHBjlNO+8ozktRH0896Sc8jvcPDjtBongkYpyIfyIMUoFmmi88YYKSjjctne0NvxX/U+FCiKLF2vu66CDbI6TIXLeedr0mXrAvHBT07zXuujEi9OT7J9aDMwDC2AZBQ4NNaIYlGQ8MNYO8qW0qPv++11k8eISiWhNKrz44r83UNQj/isyyH8gY1y2THMexLoNg00NNy87btpM5Yav0fPQar/i2z/1iis8eQkzYFmCGk666BT0f+vWiaxZs62YnZptyqsOO6zo392rl+av6tUrJRs3Hi5DhmyRjRv5t9aYEvVjbiAwXHn5cpFVqzREfeyxGj7MD84nX0ViqYoi85dIhd9/lxK2qzZy77hRE5GYdScRINywjU48KVtWQ0VmwMIPxovNSGFdmvJSkICjYFBp7pRvHj2/galEAVE2cxQmOPw3X4ksEPlrjymyA8law4DTTtMpy2+/LXLRRfo1agYxaPyfET9OO03k0ksLKCotPpYDM9LnlFOs9svYBq4+G5rcEwvovoEbb3nSeHLkkRpiIqnvAeaBZTFF4Nb+esWmTZtk4sSJ0rJlSyntZdNjum4s/kUqnP2Sd69hhBNCytRibNigISTajqEG4t9G/NhvP63LIaHuAWbAsgTahx128PY1qAMrV25L4nU8bdq/eolIznqRFid6+CJGKDn1VJGePXXsDglWDBi7cCOe5ORoIt8jD8ziP0ZqbN6s+n4SZLiVhpG3X2n58mq4SPwioabQ2YgvDRuKfPmlJ7/aDJiRGm7X8QMO8PtMjCBCUSISWgwYixYqRDNg8aZhQ2015gFmwIzUoNMCNGrk95kYQV6wMGCMTiGWTQ2YEV8aereBMQNmpMZHH6Wr8TfitGDRE23ePG3MaQIOX8cs7bef7iHef1/SGrN0//2MatLzYOpFKpMyEuyzj2eN2MNlwAqqBDayhxkwoyjwzpFOk/+qVs3vs5G4j1n6/nuRZ54RadduW8PvZMcs8dr0aOZ1+XeXLnqkBC/gFrbH2oB9+63fZ2CQ10DmaAIOoyDIjyLkYG6c2/zZ8HXM0hFH6F4i1TFLOE4YPbdLEPXIBTc4KATGZ0jcZfTUEuATG/5BF3rCQoZRmJCDWkFWWzNgvo9Zyh2STHXMEnqcrl21jSVzKjF2BYUiC6VyZfGCcBkwj2oJjCRZuVJk/XqNaRsZg4WFhQGnhZ6YzHwk/JNKrqNDB5H58/X3EPKhBCu/XMc556i+gmqIvN19mHYyfLh2g9pvv5Jy/vnFWB4wXERMzICFeszSjz/qlPYfflAP7rHH9Jzo05oS5oH97YEZ/sHK59b6GBnPdZAk/+UXTZgT0km2zaSb65gwQTs3nXWWLjx5i9ndXAc7aXbyuZk0SfMk1B/TOKFPH0defPFAOf/8NP8gVlqspOXAMnJfkJ+iyQn3BRueojyw3XKNWXK9MCK6qY5ZYmISj7v7NrK5uvZanY5Tpoz/BixcOTBUTQxQNPzhrbf0I0+UEalcB6nNo4/eNo37lFO2ytSp6SQ7/qZiRf2ItM3I6Jgl3iM8sC/yOTBeeccsARsbNkepjlliNumHH25rg8cSQHQ4JeMF7Jok7h5YyZJ6Bd2n1cgub76pH9PK4hpBznWgfCf0yI6d3/3SSzmybl1JWbFiU3qbZ9dw4QYYGYEOXfvuq4GQZEqr7r9fpH17ldFjcF54YZtXzggmNjruUkq7Skb8MWaJx5tRXs8/r948xg9hKRsg9FsjRqRx8h71tgvX3UWchE7XZsD8KWAmEA6W14hcroMF66abdPoF+0R3JE/a9seVTZPUiymFzQBMRXTtKgAxQu7GxP1aYTC+Ke+AbPfnCCXn/tytjsnve/H6OPL7v2RZs7W8eEG4DFiLFiK9e6s/W9hwLcMb74vVjLyGGbDI5Trgmmv0gOnTHdltt3Wy885pLhFuyOiPPySuFDUDMF60E5GrY54Dw4AxpoEJr0Z2IZtLOwC251YDFrlch1sh4S68ffrkyFlnfZ/+H7Lzzult1Q0jsh4Y9Udkmp94QmMuRnYglvHBByIXX7wtjGhkJEzk/g4ibe5a362b1t8gVWYSRVGQz+jUSfMjpBqGDlWVGMfdd6vDfPnl+r1NmogsX665DgQjyO2femqb8URCT07t/PMdOfTQhbJmzQHppS8cNjkVRFZtkQqOZ52EJO4zAMPC74Oek73+DlvG14ABMQ56osyZYx3Rs8Xjj2vM6qCDLCnvUZjoqKP+/bV0u3UxkisvN9zw768tXiwycqQeeenbtyTz4CV9jsf9ErlV5K8u3s/Ci+sMwLCwvoQ3udBwhRDh7LO1qptF1fAepHDPPadbfJ5IlAaGYRipkLdqPkOktRoNGjRIateuLeXKlZMmTZrIp0V0mHz44Ydl//33l/Lly0vNmjWla9eusp6ODumAlpN4CC0DzD/3HrS3xLauuEK9L+Y7Gf8KE0Xt+OOPTTJy5FuJj2n9ji/ny1+yg/x1+Q2WMjUk0S7GA1KOB40aNUq6desmQ4YMSRgvjFPLli1l7ty5UjmfflcjRoyQW2+9VYYNGyZHHXWUzJs3Ty6++GIpUaKEDBw4ML2zRkZP35v/+z+R229P73cYRYNghmISvN6aNTXBggrRiHyYiDxYuXJbEn9bWjmwstwna0X+WioSw/yXUcAgXL89MIxOp06d5JJLLpEDDzwwYcgqVKiQMFD5MWPGDGnatKm0a9cu4bWddNJJ0rZt2yK9tkJhMUWq9cADWslpeAPSNiap3nWXfo4HxsqG6sAwCsOd2xFjGb2RC+o8/PbANm7cKDNnzpQePXr887WcnBw58cQT5aP8KuESyemj5IUXXkgYrMaNG8uCBQtk/Pjx0p4S8QLYsGFD4nBZhWQq8UxsShz/VOK99JIasXvuSeXPiCzutfnnGhU39zVgwDZ5G7+T2p5y5bRkv1IlCSsZvU4RpdjXiHuE7sSoXCJ8ne1eSo5NNAL3gBKOk/x2evHixVK9evWEV3XkkUf+8/Xu3bvLtGnT5BM6gebDI488IjfddJPwUps3b5arrrpKHi9EhHHnnXdKnz598g1H4u0ZhmEY4WHt2rWJKNzKlStlZ7dGMAN4romeOnWq3HfffTJ48OBEzuyHH36Q66+/Xu6++27pSXOvfMDDI8+W2wND/NG8eXPZLfdkT7LFtC9Ag4xSLuawC5w0aZK0aNFCShen9xitHNigIJahkMiFDrM0TXvjjX+3M4/jdYowxb5GzN1gE0rhO5XVfIwgdi8lwcKF8jtdpj0gJQO2++67S8mSJWVpnoQcn1fJ3U00FxgpwoWX/11Jecghh8iaNWvkiiuukNtvvz0RgsxL2bJlE0deuEG2u0kIYxFCvOAC7TqZ9uyHaPGv65QKVLISNqQNBO0gcv8eunxScUvLhgg8rMW6TjEh7WuE0aJAjo8MKqOGMMLYvVQIX3whpT3qiZmSiKNMmTLSsGFDmTx58j9f27p1a+Lz3CHFvK5jXiOFEYQUopcFg9FiSl/nzp5JNWMFLcmZxfH00/+u0CX/xSwOt+eQYSQzududI2fEk5kzgzMPjNDe0KFDZfjw4TJnzhy5+uqrEx4VqkTo0KHDdiKP1q1bJ/JdI0eOlB9//DHhbuOV8XXXkBVbx8yiizHEiBnpw04ZcQwKT7rJ5gczGGjhYBiFwT1CR2FmeXz+ud9nY/htwEg9eEDKObA2bdrIsmXLpFevXrJkyRJp0KCBTJgwQfb828IuXLhwO4/rjjvuSNR88fGXX36RPfbYI2G87r333sz9FdSfDRqkoUQau7mN34zkobCcFl1cS8KyBYEBs6GiRlFwjxx+uG4wzQOLL46j7/+VV26bJ5hB0hJxdOnSJXEUJNrY7gVKlZLevXsnDk+huS+vTa9ERtrS9NdI/iaj08ZXX4m8/37hzf0OPtiTG9GIEJRg4M1zr7DhQfRDB5eICjmMQmCAGbW6Hnlg0WpsR2cOFIl0jkh1JG2coSMK41fJexWlFmIULIuTFagaBcG8FzZF3CuM8qUW7Jtv/D4rww/c0iozYEnAyNJXXtHePmecYbOIkmHCBAr5RG69VUOIReHOMp81y/NTM0IKOS8EPwceqFPU8ejHj/f7rAw/4H3HE8+nzWAmiJYBA2axE7L4/nuRM8/U3I6RPzNmiJx7rkirVsl3M6lbVzcIltcwCoJ7g4FmtB6jHKZlS5GxY/0+KyPbEDZ+6y2R00/37CWiZ8Bcd5ULx0RAFuhcbakM2Ta695RT1KNinn2y+Qm+j+S8KcuMwgwYoUOX1q01lORRPzwjoNBekPyXGbA0oFPEmDEi776rnphHhXSh5MMPdfwuxaUY+lTbc2H0zAMzChJwzJ27LdQMePioEceN8/PMjGyD14063aMuHNE2YEDoggWawly8DetcrzFprgstuCZO1I4bqdK4sciCBVYPZuS/60bAwT2SO6yPuMrCiPFi7Fj1vj0cghttAwZ4Gu+8oyooHqpvv5VYwqJCfddpp2mR8ttvp2e8AANIKJHNgWHkXbTowIGAIzcsZJMm2RDauDBnjnrivO8eEn0DBtSEkfNBfIAqKm51TIRPGV9Dlw26pBBaLU5Xf8aqcE1tR23k3STxbJHzIGSYG5oMkIseMcKvszOyydChIjReP+kkT18mHgYM9t5bVXctWqjEnsLqOMzwQY3ZrJnIa6/p/DQ6oGTCpWeRIr9opQqGC4Xw1F/ml7SnrRS7cTrm2EDUaLN2rcgzz4hcdpmWU3hIfAwYUI/y6qsMHNOFvEkTfeiiCF3lH35Y5czk/lBksgvOFCxG7KgJCxkG4JEz6+nYY/P/f7rk8LwVMPzWiAgjR4owwJL2UR4TLwMGeB+9eol8/LF6YMh9mXkVJW+MThl4XV276mgUFo3cqrBMQJPWAw6IXzjWKNyAIZaioUBB+Wime9N824gmjqNeNvdBnTqev1z8DJgLhotaJrpQMHjvsMN0MQ5zeMMd2920qSoEUV/SXovcnxfghXHNomT8jfRYtEifp8KS9mwemXRAtxwbfRRNPvtMu/TgbWeB+BowoEsAHSg+/VSlvsTujzlG66TCJtIYMGBbvzE3VFNQKCdT0Hpq2TLzwgzto8lGqSjV2cUXqyEbMiRbZ2ZkEzbMtWuLnHxyVl4u3gbMhc4S772n0nJECSjskJvztSB7ZH/+qTcM7Z3oZUgTYyBE6pXXlRvya9T3WEgo3uCBP/mkKl3JgRWlYCU38uCDVpcZNb76SoViRLWyNHnADJgLsl92DXSYQOr7448iJ5yg9SyPPKLGIijMnq25rerVRW66Sb1G6i4eeij754K3x4Tu//43+69tBAN6jzKBmfBgMtx+u4qM+vb1+syMbHLbbZr3yuI8RjNgeSG80batFj4zX4yw3I03qrG48EKRl1/elmvKFniBGKj779dibDxGushT04VsGYOLqMIP6DVJ+NVCQvEFD5yoRbIjM7hf2Hg99piNPYoKH3ygrcJIyZQunbWXTWugZWw8MpR8HOwuifFTS4WxoMs2vRaJ9yOYYFwA+bRMQpKb3BzKLowVSXLUXRgqQjWoudjFUm/B+XDTVKqkB91GatTQQsK8BaWZhr+bHReLGKUJ2QhdGsGBjdWUKakXKHfrpmo1BFQ8W0Z4cRxNYSCEO//8rL60GbBkqFpV5I479GDHiGgBw8Iukvg/xuOQQ1SqTo/BmjX1ZzhoZomByQ/a6mAc581TWT8KHjw/jFd+yj6+hiHjZ/idxJk3b9aD/2N8ATFo8lIIOzgvzoGpuDTu5fw42ClnssCQnEa/fiIvvqiTnY34wMaFWU9u/jVZaGPWs6fIDTdohCNv6ykjPLz5pjaJoLeqh30P88MMWKrQUaBzZz0wEl9+qXkzDkZGDBumhsQFD4jENcYEo8NuhZ9DLJLfmJeKFfVhxtNiR4NhxABxsFAUlhzld9Idg4PRFRg6Dowe5zd8uBo7zgODRikBIUlqNjC66ULvu7POUiOGyqygOiAjWvzyi8hTT2nSPp0IBJsd8rbXXqsF8Vle/IwMdd1gA0J/VbocZRkzYMWhfHntrcjhgvFCWk4dlmtAMCY//CDy3XfqbZFDw6ChHiT8SEcQQpLktgry1pLBNRyML8gvDs1wz6+/1nod1+g++6yeM8aSMgIO/p1q6PGuu9SzQ43WpUv6f4MRHnjPCRkTDkwHjB73CwsfH6+6KtNnaGRDuPG//2n+y+t0RT6YAcs0eEhVqqjXxcgRvCHk+atWqfd20UVqJMitZTpvVhSEDTFuuefzoK7E9SckiiSfnASCFfJ7iFlQOCZzY+LRdeigXU3wwmjbZUQXOo2Tu2LCAVGDdKE7B57YzTerCpgaIiM8wo1HHtEaVDbjPmA+e6ZhN0IdFmG1885Tz4tc2RdfiPz0kyqv6NCcbeNVELvsItKmjeavyL1R+8Z5Y9QwsoQwyXNggIsC44dBHDgwG2du+An5YMLamei40L+/bvho/oowyQhH6PDSSzXffv31vp2GGbBMwEOHp0UeCMNFE91zztkWriNZTdGvDy52ShB2bN5c8xL0U+RvqldP5Lrr1CtjseJvKio/yI6MMKoRTRAb0RSbEGImxEAUP+PNsXkilGgEn9tv1806Of8sFS3nhxmw4houOi8juiCOj7eFNJjkNp4W+a2wgrGlkJuF6uefNVHLHDHyXCRsC+ooTkwckNQb0ZVMc89TzpEpcocSrSg+2EyerOkGnnGfQocuZsDSfYiZ8oyKjzwRNVnvv6+tVEhEpzvpOKjgfTGCBkNGIfeKFRo6OPNMFabkZvfd1Yg9+qju1I1o8cIL6imR+8r0zhvPHS+eeX1//JHZ321kBiIzpBjYsPsYOnQxA5YqFBfjmbRsqSpEEplvvZW82CHMEGLk5qXbNAsZBpsc2SWXbN9RAW+NejjEHPmVChjhBGUt4WQ60px6auZ/Pxs/xESEn9kYUvJhBIdVq1SAxiaVyJOPoUMXM2DJQliQxRvJO2IHHjSGRNJCJ25Qr8MiRqiHUML48Vq3Rj0QyV0MHfJ8pkEj7DCiEXWgYJ2cF8ozr9hnHx23Qv71llu8ex0j9XQJCmryXqx9dPwJAGbAknlwKQBGJo7BYmGmeBmZedQ9rmTqzqj5IqxAUpewIZ4X42jwzFBj0r/RQonhB4+bSMMTT6hi0EuIcCAkQs3Ks2f4T8+e+v7jeSHsCghmwIryujBUhMJwnekx2LFjIFznQEHNF8aKUgHCC4RTKW6lw4KFEqMVOuQ5yAZsjOixibCDptqGfzz7rMh992nek649AcIMWFFeF90qcJmfe877nWfY2X9/zQlS1/P441owTaKXUKKrTjTCBV1aqM/yOnSYF6IbKHqpRWQ2H732jOwzcqS+/2wkyG0HDDNg+SUqaUya2+sqasqssQ28U2501xvjOtImi3AQYSgjXNBZBcUtG7psb+AIUVO6QQNqdv4IqIzs8frrmvfC82ZDGsCUiRmw3FDHdeSRKhMePdq8rkx4Y9SK0KiVZsGEhGwRChfkopC30z3GDypU0NwL0RCk2+RXDe8ZNUpHozDvj2LlgDZaDuZZ+VWcR2d2xpLQVZ4aJ6P43hhDN2n0SYspHgK8WXIqRrChVALatdORJ36CvJ7WZuRTKV+xnJi3DB+u7zsHUZPiNBj3GDNg5LtQz/FgYMAwXgFS2USCVq30ujIbbflyVZnRGd8IJkxQYPEC2qIFIXSEEaMpNgX0PKumTvRGKk/DAsL+9DlEvBFg4wXxNmAUSlLbgsIKsQGhioDUN0SOAw7QXT3CDurHjj3WGrcGEer46OnpGq2gNJ12w4kMT6SFFYssuVYrds4MDNelzpW6zXvu0Z6UAQ0b5ibY5tVLNm7Uan8Uhs88ow+E4S1sDqilQxzDbrppU1WXBWGHb6hXTOicqeDkgfHEggYGdehQbY7dtauKrJhCbhvP9GFKBu27GP+EaIZ/h4Tgm1ivHlSUhnhcKG3MeGUPQhLkxMiFffyxhoTMEwvGho6EvdsajaGmQYUNDzWG5MUQBTFQlvlkRurQw5WoyOrV2qA7RMYrngaMglqMFztMwhEmkfdnAcLz5WHBiNGOi1yk4d8zQfgItegbb2jZQxggl0qXFzZFSO2RettmKPn3nO4aTJagaw7XMYTTM+JlwFAYMrzRNV5+SYMNxQ1XsPMjnGhGLPusW6dhQ7wZ3o+wPRP0TmQTRL0S8+qQ2hMSMwqGXDSTNPr10w46vPe77SZhJCd2zShpPEuNF7s3w39YNPGCMWKEgthkGNmBkSUoRAkjEdYNWJuglBSKQ4aoB0ktJ56EeWP5h4nxulBb47UybBcDRvPtkBIfEUfv3trlmgGNYX1QowphK4wYi+jhh6uHvMcefp9VtJkzR71fyhomTND+lV7mnNesUcVg7oPNCuFkFlP3YDHl4w47pK6AZCgmAhSGYuKN8axTIsPwzbgzfbpeE953jBb1mSE2XPEyYFSVIw3FZSb/ZQQzJ0bfOzotICDAU2b6s5F5uLYocGvU0NwHYbh0JfcoFSlML+wjRezpQFiralWRatW2/5j33/RpzOuNIUihfx/5HRpwIw+nG0zc+Ppr7UOKMIc8IV4XCs6IUCoW8V4GLlKYybwqI5hQc8LCSmyeHAbqRFp52YYjc5BjpNHyrbdqg1y6LOy8c3KGit6WLH40t549WweYrly5/fdhSHIbG9o/uUYGw5Lbw3IPcL2x3N4Zvxvj5xrCefO0AwefEwrLDRL62rXVe+f+YaFGGES9IeNf6Oc4YoR2uMfzCGm+JyV+/lk9reefF6lTR0sNaA0VgtquVIi2AVu6VMMkhBCeesrqjYIOixzhLBYhFrdzztHOAHfcYSNsigtGCI/kxRd1dttdd+W/mBHqw1ix8atVS/OSfE4+iea6eMVsLsgn5/WMKlb0/hnDCK9Ysb1x4yMTDzCwdOjAAHL/kAvDmGG0XGNGDRkF0FddJVK5skSOBQt0yCxeKH1c6ehPD9IIhAvjZcCQibIAcjMjFChf3u8zMpJh7701d0E+g8UTA0ZynmJzpj4b6eU/aA3ENF124hdcsO3/+Bq5RwrK8a7IkWCsMEYYA5pbX321GgI8KoyYn2Ag8aA48pN9k2/76iv9WzjyGjW8NYw3B5tb5tZhkMO8ud2yRTd+gwdrg4BddlHvix6W5BIjTKlIizYockRhRazfCA/Nm+vsKUI+dLN/+mnd+TNUj7Zf5o0l73XhbbEjZzNA6UjduhoCJOfIgafF9ST8hpCDRQ9jxfe9+652ow/T7p0wJio7joKMGn056eBBEwOOKlW0pRyeGVGAsLBsmXaKx9si7M77xrNCqRBtt+KAEwJWrlxJgZCzfPny5H7g448dJyfHce6914kTGzdudMaMGZP4GHq2bnWc005znCpVHGfhQse57jqCR47TtKnjzJtXrF8dqetUEB984Dj77ec45co5Tr9+jvPWW45z9dWOU6OGXseKFR3nggscZ8QIx1mxIn7XaOlSx3nqKcc56ijHKVlSr0mJEo5Tp47j3HST4yxalNSvyfp1mj/fcR56yHGOP95xSpXS9/fiix3n00+dIMPazRrOWp5JopXRc3dbtIZCyWaijfBCSIecBe8nOQy8iGnTNN+BN4YYwTra51/bhReFN0UokKbJKHARbRBmIqzO6CB274QTUSPGsY8g+S8mDaN6pZEtIWq8NkKqzD9DschAVsKttNfyq2kwqRDCuygJCZmiGL3lFi0xeOQRPV/OnXZQMSR6IURyJiQyCRUEfBSAUQSIA3hIO3TQVkfkLAgF8TCjpHvsMZVH05087mFFwoWMPiHMimFnA8BzQD0dGwAaKJPDCnOux8uwI5teDow+dYjcW6geKcHh4P5izBLhbVfpmOmxSxgrZO9uqJODz1FoYkxPPVVzd3QbCVOo00OitcIT22Znjmw2hH29jHxA7fbyy6oaw6tAWYU3hrAAdSIlEuyYWbgpho7bAs3ihpHHkNOQFTBU5A9pEUV+x0gelJkIiDgAdSMbA+5B8mZ8jnEDhGEIizAqKP3YLLiqTO5Tt2QA44fQgveKA+/XVU/mVlNSmsD/8/28hxhJxDd4V/w77pu0SBswdi8sZhY6jGYokQeamW3UtQC7X9SKCHXwxvDOUJOxeWGXHHVDRkiLwvwHHlDDxcKLweLeR7AR9b8/W2CgkKLT0QPPDKUfnWMI4SF8QewCGCDaofGRUoSiQOXplh+45QpuLRuFxqaajpkBY5FjpAIKKwsdRgsedLwsdrmoEHPH+8lbkNNBMYcho8cldX+0zSG0mEyhbpiYP1/kpptU+s5unYWOa4I3au23suOZLVqk9WQMfUTFyeaJe5NBoJQZsKGgVs31uviYu3ibcGBcVIIeEw0RBzcMiWpaxlj7oWhCfgIvDCOVt2s93gZ5AWp+MGZMf8Zbq15dDRl5hDDD30sIa//9RfbdV+saCQ0SOuTeJ6Rqxit7IPAgbEjI79ln9WtEf/CkaJZLh/+99tIOGLxnbKjw1Pic7zHjlTGiYcAGDhRZtUrzAEY0If5PnoswDoXN+YEhY74RoUVa6VDXw2LPpobWQuQy8GDCYrTogIHijB07tT10Wufvc1s5MdTR8iL+gbeF1wWMdEHhyT2GUhAxGWuS4SnhN2C//abhJZLWcWzWGScQaZDnQlVX1KgMvC8WEQwZUwgIJeK94cE0aaLfQ/6M8E6Q8rjMZurcWe9l8rnkuFgI8UBpiouH2aCB32dq5AWvn1wZBcUIjshPYsjwjnlfDU8IvwFjV058mkXNiDZ4WCwM5B0wSsmAEozO5DQKZnQInRdQdAFhR5LoeDcYCry7dDunp8OSJZrLInKAccbTOvlkrc8ih8K5k99CpUatj0mngw8trlBCo1Z0W1URRqQxdZA2SxEh3GoHdtcMrqNtVBw6TBsqpaceBtECnepTaXO0444a8qGoF4OGt8OYCQpFMSSueoydM0aOA2VY7oa1qeQvCAOSo8otmSYM6Pbo4+uA5Br1GcaLUBRfp/YND5JcihE+8KBpIE4Ym3uV/DyGrW9fvX9NJZoRwm3AkLayKJGwN+IDUnkWfHr7FWfcCtJld5Aju2NGdrgjQzgwcnRpyCt/xpDtuadKqV1lGcbKHQVCEh/vCkNEgXFu2Ghx7hgo10jSQJfibJSU7NrpmIFgxYhGaPG117RGlRA2njb5WCIJTZv6fXahJ7wGjEWCRpa0g4l4x2UjD+SGyIVRk5OpeWGIIVhsOJDfA0bJnUuV24vi3+Re3cJUDsLYdHTAmKEQpIt73sGLHGy4XBjxQaEqHR8wpLQ14u8yogd5V0LU77yjhgwjxmaFomhrNh5DA0Z7F/q+kTA14gfyeLp0YAQy3dLHhTAPoyk4MGyZAm+PLu+ElpBVE7485RQLK0Ud3t+WLTX3yvp1883qaXMvIMO39z9GIg523+QM0h2HboQbhBmIHhglESYwuOy+6ZiB2vDLL0VatbLFK07grdNE+ZtvVHpPFIl7gMa8RgwM2Gef6cEu3Ign5J/ofkAhaTKte/wGr4tyDyTwKAwZMvngg9YyKM7g2ZMGwQOn2B5vjM/zFuobETNgKA8JvbBrMeILQwipkUJ2HiaviwJly3UZLqxj5o1lz4ANGjRIateuLeXKlZMmTZrIpxSEFsKff/4pnTt3lqpVq0rZsmWlbt26Mh6FV7riDeLHnTpZF4K4g8Sd3JHbzifIXtfvv+tcKfO6jKK8MUYGmTfmjQEbNWqUdOvWTXr37i2zZs2S+vXrS8uWLeU3VFn5sHHjRmnRooX89NNP8uqrr8rcuXNl6NChUp1OCelAJwKkyeRADINdK3VcjKgIutdlsmmjKPC+GNuCuta8scwbsIEDB0qnTp3kkksukQMPPFCGDBkiFSpUkGHsFvKBr69YsULGjBkjTZs2TXhuzZo1Sxi+tBg7VkccUN1uGBSFAjvXoMDIF6T+ub0ua+BqpOKN0XnF9caYbfj2236fVfhl9HhTM2fOlB652jbl5OTIiSeeKB8xCycfxo4dK0ceeWQihPjGG2/IHnvsIe3atZNbbrlFShYQAtywYUPicFn1d1PMTRs3yiaKPfG+/BrxHWA2UY+U62MsoIvFscdq8e+FF/p7nQgZ0haK/ndI/N1wYQjfj1jeS0G7TsjtUamS62XNo4Cfnq8hVKxu8ug+KuE4yQdZFy9enAj9zZgxI2GUXLp37y7Tpk2TT6g2z0O9evUS4cMLL7xQrrnmGvnhhx8SH6+77rpEGDI/7rzzTumTT2f5ESNGJLw9wzAMIzysXbs24bisXLlSds7gjD7PC5m3bt0qlStXlieffDLhcTVs2FB++eUX6d+/f4EGDA+PPFtuD6xmzZrS/KuvZLfhw3Ukhg2tzHeXM2nSpETOsXQqPQLDDvkmOh0wRoVda7avE/0NL7hAu3MgKGHkSciJ7b0U5OtESypKhxB4jBihHV9Cwu+E0z0gJSuw++67J4zQ0qVLt/s6n1cp4GKiPOSNzR0uPOCAA2TJkiWJkGQZZurkAaUiR15KT5okpZm4ayquQuF6x2rRIUdAyyZ6I6ZQWpGR60RDYLrZ0ybq/fc1PxshYncvBfk6sUni/jrzTN2wjRmjE8lDgFfXJiURB8YGD2oySsBcHhaf5w4p5gbhBmFDvs9l3rx5CcOWn/EqlO++E2nePLWfMaIPOQHuiwLysJ5A5J1BqhhM1IV0kY+Y8TICCM2faeJACcmxx6pgKMakrEIktIcMfvjw4TJnzhy5+uqrZc2aNQlVInTo0GE7kQf/jwrx+uuvTxiucePGyX333ZcQdaS1aLiznAwjN9wXyI+pE/Sa9eu1dx2jMuhnhzKWLvWGkQ2Idk2ZItKunU41oFQjprPGUk4ktWnTRpYtWya9evVKhAEbNGggEyZMkD0ZLyFMOl+YUCa6kLuaOHGidO3aVQ499NCECARjhgoxZej2feCBqf+cEQ8DhjIV2bE7cdkL6EZPjc7s2SIvvJC08tEwMgoplqefFqEcCb0AnTzoSBOzjVRaSoguXbokjvyYymiIPBBe/JgQS3GhI7jF4438OOQQFfYwx8srA8YAVQQaeGDUdx1xhDevYxjJhs6vv17XRfKw3JuMa4nRcN9w9UI89FC/z8AIKnjnGDEMmBegfCXnAMztMuNlBIWTThKZNk1k0SKR445DVSdxwQyYEa0wohcGDJk+xgsjidKQBLphBG1tnDZNu780aybyyy8SB3JCJ5c2jII4/HDNBWSy6p/fx4JQqZIuEOn28DQMrzngAN1g0SuWe5aQd8QJlwGrWdPvMzCCzF57qRorUyEUjBchGWrMyO2GqHDUiCn77qtGjLIljNjChRJlwmXAdtrJ7zOIDETB6IfMpA8+9uuX2s/TdIKB2JQ+4RjzzOTHX3/pFHWGJ9OjNC/9++vPIy496yxG70j6UEzsKgWLy9y5IhTNs2mi7pE/wDDC8nBPm6YiD4QdixdLVAmXAQthE8sgw1g1pny8955I374iRYx1245bbxX5z39Evv9eG2dTkpJf5A7RKBUT9GDOy6RJ+rPUH1OjTgrr9tuL8QfhKUFxH1gEGzz4e+yhqi4aBhtGmKhZUx9smqKzEStg3FXYCZcBMzyBtE69eqmFzF9+WeSqq/TfCPKwHWz68itXwRbk533RaJuRWa5jTVOLYjUWwODQsqw4Hpgrld9xR7W6/E7DCCN7761GjLDGiSeqwCNiWEdcIyGy494m3bN6tcgxx+T/fdSq0/qP78Xbyp0SImqRargdj2vwYJElS/R3v/iivv6KFWk6PRTQc1LpGjByZ+xWMYKEDS3nZYSd/fbTe5l8GHJ7crkRSsWYAYsx1D6y5pPueeihbc4GYcVsQPvCm24SOe00tRnkwKBYgwbIg6UTQiTUQocNFFxMeK5RoxgnYRgB4sADNZrAzrR9e5HXX9cHPwKYAYt5DgwRB/d269YaOcOTKsoDo9AfI4Pn5DopP/0kUqtW6ufAdAgOoFkLdqNY44IqVxZZtiz1HptXX601ZMRBrc7LiBr16+sIltNPF2GMFcMxI4AZMCMRHmf9vuMOkTfeSM4DO+88kSFDGD6qzbGpmyRKkSpE+3CacHx69dK+pMWCpFuuad5J8cgjqiZ57jlv+ygahp+cdpqqtVBgIf0lBBNyzICFEBwGFvy8kJdav76krFlTdMtIfgeN2/leoB8om7Tp00UOO6zoc8DYdOqkZSe81tChIhs36sHmDqN0+eX6vdiE5csZTKqCEZpaPPXUNuNJyQrnzrgjmry755QeO4is3ywVnCRFq6gM+eOJZRJeMYwo0727NrzmQSM/RvF/iCnhOCxlwYaJzBUrVpTly5fLbjFqVFkQLPCI5IyCof5shx0KnqI7fvx4abX//lKaOXbUA7z1libijO2vUatWNtAyatdp3TrdRZID+PxzzQ1kYSIzA5FXrlwpOxcrR7A90cjkGUY6tG2rDy9jKMx4GXGhfHmR0aN1/BDCpVRD7gHCQoghpEIF9TDy2w0ye61ly5bh2Q1mmosuSsQqK1QYW/D3uMP/SMDRQiS/IjXDiDI1aqgRI3GNioqYfggbRZgBCyHcZ/mFx8gjlSu3JfF/cbVfImtEym8VKexZJEnXuLHIs8+K1K2bxXMzjADxn/+IPPmkyMUXawL8uuskbFgI0YgWqFsYe1IQVEtT9AYULRtGnOnYUeTGG0W6dtXebiHDDJgRLQgLFtRBg8aNSCfJfRmGodx/v0iLFvpchKxnohkwI1q4hWX55b2QDtO0ceBAP87MMIJJyZJaA0lugnxY8IXp/2AGzIgOFKFRcOZ2pc/No4+KfPihyLBhqoIxDGP7DjaDBom89prIK69IWDADZkQH6logrwdG6PC22zRJTf2LYRj/5vzzRc49V6Rz59CEEsNlwGjZYBgF4Xahz23AcocO77vPt1MzjFAwaJB+DEkoMVwGLILzbIwM4nahzx1CzB06LKg1h2EYoQwlhsuA/fij32dgBBlChRgpt92YhQ4NI9KhxHAZsK+/9vsMjCDDOBQ6ETPryEKHhhH5UGK4DBhdlA2jMAPGmGew0KFhRD6UGC4DZh6YURB//ikyf75Io0YWOjSMmIQSw2XAWJjy62JrGLNm6Uc8MB44lIgWOjSMzIQSb75Zgki4DBix2GTGBRvxDB8SKly4UHu6DRhgoUPDyEQokbHrzz8fyAhYuAwY4+KZX28YeWEwX4MGIrffrl22zzzT7zMyjGjQqZPI3nvrsxUwwmXAjjhCZOJEv8/CCBoM5sPrYjglnli/fqGcbWQYgaRMGZF77hF5802R6dMlSITLgNExecoUkdWr/T4TI0igNvzjD/XOTz5Zh/QZhpE52rTRCMettwZKVh8+A0bD1nfe8ftMjCAxdqxIxYoiixaJ9O3r99kYRvTIydFni83iuHESFMJlwGrVEjn4YF2wDAPYDY4Zo+Oo27XTXaJhGJmnZUuR444T6dFDGwUEgHAZMDj9dN0BBOQCGj4zd67IggUiGzaI3H2332djGNGlRAnNL3/zjciIERIEwmnAaOr70Ud+n4kRBEaO1I+XXy5Sp47fZ2MY0aZJE5GzzhLp1Us3jT4TPgOGEpEi1Zde8vtMjCCEDx9/XCfK9unj99kYRjy4916tt3ziCb/PJIQGjGQiu20K60yNGG+Q9dLihpY3SOgNw/CeAw7QRtlI631eg8NnwOCKK0TWrhV58UW/z8TwE/odsqFx290YhpEdevfW0pXhw8VPwmnAatTQXNjgwYGqSTCyyC+/iHz7rchRR4lUquT32RhGvKhZU3NhPq/B4TRg7pwaenNRl2DEj5tu0o/9+/t9JoYRT665RmTOHJFp03w7hfAasOOPF6lbV3cARryghGL0aJEqVbTvoWEY2YeON+TDfFyDw2vAyH1cfbXIq6+qIsaID//3fyrhvf56v8/EMOJdF3bNNbqZXLzYl1MIrwED1IjkP2j3b8TH+6JgmckEbhjRMAx/aN9en8WnnvLl5cNtwHbcUaRnT1XCfPed32djZANqT5i+3LGjSKlSfp+NYcSbihVFLrpIn0vauWWZcBswV1K/114qqTaizfr1OpOI0AUyXsMw/IdUDiFE6jKzTE4kZtUQUnrjDZEZM/w+G8NLqPfC+6KpaLVqfp+NYRhQv75I06a+iDnCb8CgbVuRQw8N3KwaI4OsXCly113671tu8ftsDMPIDWKOyZNF/vtfySY5kZpV88EHIq+/7vfZGF5A25q//hLZbz8bWGkYQeOcc0T22ENkyJCsvmw0DBiccop252AnsHy532djZJKPPxZ58EH9N+8vOTDDMIIDSsSLLxZ54YWsjrqKjgFjUcP6o4S59lq/z8bIFOvW6YPBqJStW7V9jWEYwePss7M+6io6BgwYs/LIIzojykKJ0QC14Y8/itSrp8liFKeGYQSPxo1FKlcWGTs2ay8ZLQMGF16ooUSknRZKjEboECNGz0veV8MwgqtFaN3aDFixsFBitEKHjRrpFFjk82bADCPYtG4tMneuyLx5WXm56BkwN5T46KMaSsyyKsbIAJRCXHmlyM8/izzzjMi4cVr3dfjhfp+ZYRiFceKJIuXKZa2oOZoGzA0ldumiXpiP7f6NNBg4UCduP/20drsmJMHOjhCFYRjBZYcd1IhlKYyYE/mF8NhjtUYBIYARfCZMEOneXYvS27XTeUPz51v40DDCwumni0yfropEj4m2AStdWuTll7Xh5BlnaCGsEVyInV9wgUirVlq4DIQiKlTQ+W+GYQSf007Tkpfx4z1/qWgbMNhtN3Vn8cA6dNALawQPV6RBruvFF0VKltSv896ddJLG1Q3DCIcGoXHjrIQRo2/A4KCDdFEcM0ZnSFm/xOApDimCXLZMb/qdd9av//abFkVa+NAwwkXr1poOYPCsh8TDgAGLIMrEhx7SGWJGMOAGp7vGJ5/oRIF99932f2+/rR9PPdW30zMMI831lpSNxwK6eE0E7NxZZ0rhhdG7ywyZv2zcKHL++XqTI5U/5pjt/x/vC++Z6n7DMMLDIYeIVKqkzQhIAXhEvAwY3Hij7voZjEixc58+1hzWD9hInHuuyKRJGtrNT6Qxc6ZIw4Z+nJ1hGMWBNZVnl2c4aCHEQYMGSe3ataVcuXLSpEkT+fTTT5P6uZEjR0qJEiXkzDPPFF9hevP99+sgTGZLWU4su6xdqzHy995TlSGTBPLzzr76ygyYYYSVRo1EPv88WAZs1KhR0q1bN+ndu7fMmjVL6tevLy1btpTfSLgXwk8//SQ33XSTHJM3TOQX1Br93/+J9O+v6kQ8AsN7fvlF5LjjNDyIzLag8MI336gRMwNmGOGkYUORxYtFliwJTghx4MCB0qlTJ7nkkksSnw8ZMkTGjRsnw4YNk1spPs2HLVu2yIUXXih9+vSRDz74QP5EMl0IGzZsSBwuq1atSnzctGlT4sgYNPwlv8KMqRYtREaMEKlSRcKKe20yeo0yCbsxOqQgkZ8yRaRBAw3j5sesWVrVf/DBBX9PVK9TALBrlBx2nQqB57t8+cRzvwlvzANKOE7y8bONGzdKhQoV5NVXX90uDNixY8eEUXoDFVk+4K199dVXMnr0aLn44osT3zuGvEcB3HnnnQljl5cRI0YkXt8wDMMID2vXrpV27drJypUrZWe3TCbbHtjy5csT3tSee+653df5/L///W++PzN9+nR5+umn5Ysvvkj6dXr06JEIU+b2wGrWrCnNmzeX3ShM9oKlS9U7IO/CTDE6QoQMdoGTJk2SFi1aSGm6kAQBprMyDoUShosu0jKGMmWK/jnCjPRBfPzxeFyngGHXKDnsOhUBjk758vL7Y49J6FSIq1evlvbt28vQoUNl9913T/rnypYtmzjywg3i2U1So4bIO++IXHWV5sSoS0LoQRgrZHh6nVLhf/8TIdSMWIO+lNddl5zik9wX4UY2ER7+HYG5TgHGrlFy2HUqgEMPTTTm9urapGTAMEIlS5aUpXgrueDzKvnkjubPn58Qb7RGcfY3W/9u5VSqVCmZO3eu7LPPPhIYMJrDhmny8eabtZKcz2kIbCQPUWnGoHTtKrLjjiITJ2qH6mQxAYdhRIOGDUX69dOuOn6rEMuUKSMNGzaUyZMnb2eQ+PzII4/81/fXq1dPvv7660T40D1OP/30RCiQfxMWDBx4CIxh+fJLYqMizZqJXH+9yJo1fp9ZOMDrohnvZZdphw2MUSrGC6gdYXQKSWDDMMJLw783oV9/HYwQIrkpRBuNGjWSxo0by8MPPyxr1qz5R5XYoUMHqV69uvTt2zdRJ3YwKrJc7LLLLomPeb8eOOrWFXn/fc2HUTdGpwjzxpL3ut56K/0WUBiwAw/ULvSGYYSX2rW1IwfagiDUgbVp00YGDBggvXr1kgYNGiQ8qQkTJvwj7Fi4cKH8+uuvEgmQe7Mg5/bGECLYbLHtmTFDr01ur6s4/Qtnz7bpy4YRpY4c334bnE4cXbp0kZ9//jlRq/XJJ58kunG4TJ06VZ599tkCf5b/K0xCH2hvbMgQEcKn+++vggSP4rqhgZuSOWtNmyIV1VwX7z07ruKwaJHI3ntn6iwNw/ATnmUaGHhAfLrRZ8Ibu/JKkR9+0P6Jzz0nggDlzjuRW0qsWLhQ1YUojIhtM6qGwuNMNO1Edo9IiLlghmGEn2rVgiHiMERl9T166Jh7OnmgsMGQUevk0S4jMBBKxYjjkdIGilZc1P+1a6eii0zAjY5SlaF4hmGEn6pVddafB5gBSxcKqh94QD2y887TOqe99tIO69Q9RaVBMC29aLF19NGqCkSc0auXGnDUmskUJaeCmz81D8wwokG1arop9QAzYJkogB40SL0vFIt4JCecoF0k8FB+/11CCQYK9SWlDnQoKVdO5LXXRH7+Wb+O0tALaP4J5oEZRjSo6t2zbAYsU9Dfi6bA5IQY0Ii3wuBMmgUjvR8wQGTePAks7JAYPodxYhgdk5EHD9bw4Jw5Iu++K3L22VSge3seeGCEI22IpWFEg2reRVPiN9AyG7JRDBYHYgTmXY0dq9Of6e5B/ohx28jMjzjC31ZVy5eLfPihniOhQc6X0Cjnhjjl5JOzf354YBgvrw2lYRjZgefZo6HBtkp4CbVjl1+uB0Mc8WIwZigY8cjwNOrV0zoJxg3wEc/NC6OBsaJAOPdBOBAoC6D/I4aVjiooLv0CD8zCh4YRHUqVog+hlsdk+ldn/Dca+UNXCQwEB+E6KtNpWOsak1GjtP8fRo3qdRZx3nRCk7vuqp4R9VV8Ts9GjooV9QDmEbn5NoZ00tIJY8DBv12FJD9PkTBiE4wmXmCQ+lFyvibgMIzoeWGLzIBFA7fPH4aqTh01IBgVcmQ//aTeEh9TVe4wPO6ll1QRudNO+vsPOkh7EdK6Cw+P18qU5N2rEKL1QDQCBo8Se0YesXXrtAyygPm9BVaHEORAG8XvIb2cX1c6mvywt6QccvNm1YI9+aTuXQv7v8CTZwRXpjADli3wrmi59Nln27wuJPiud4ZoAs+DOVh8xAPj4x57qFSdO/aPP0T++kul7RwrV4qsWKEH3hdhSkAxiAHkQBWJ4SIsyOuwE3K9tiCCB0YzYMMIGARJ2FsRzKBV5/HHixx2WHI/i7H7z390wAVLAB3XMEh5p4zwyE+froYS6CNOOhpBc2H/F3g8mnRvBsxLMCxvv615Lz7SsQMjwl2PUALDwkEOKhN5J8KIFBgvWCCyZMk2Q0moki0fnh3xaIwkI24IZ7K1DAp4nJy3Rze7YWSC6tU1dU0KOVkD9vLL2/arRO0xRoiV8w5qyD0GEU+LIRhuxUph/xd4PFIVmwHLNCj5KPzFaH3wgd5p5JpQIJ52mrZf8lokgeKnVi092OoBhdV4ZBg4zg2JP1s4JPMYsrZtNdzoJxhgrldonkojjhDUIODBPpA96Q03HCe9epX6l9COqBntQflebu3c+zL2jXRkKyhY07ixGkiWCx7XZP4v0Hg0WSLAyZAQgXHAWGEEKPyl1RRKQgqcEVAQM0BGz3bNL4UfTxdNNTt31qcKb+yVVzQm8vjjmiPjiWSryFPiB4RJwSbbGgGkTRvNOxE+vPZaje6Tan744any+eeb5YsvZLuDxywdyBjw8+yF8fSeeCK5/ws0HpXFmAErDmy/WPzZCpGRJVxHeynyONRV0TeQeEMQQY1IRhhJP+c7cqQaYp5SWmLRLgrjm03YpoLVgBkBzYFR0//OO5rTomeB64E1alQqsRfMfbRsqT+HgJhbmui4C8EQAiSFgbFCLPL886n9XyDxaONuBiwdGB1C816ME/0A99tPZNIkjS3ccENIZEF5ngYMF0F5nko6bjz0kBoypFM8bdn0wMyAGQGGvBV9vO+4I3kPjHapTGMCAjIIQRihlxdCg64Wi5QwQRL2x0X9X+Dx6Jm2lSIVUP7hcd17r6oBMV7MBSNsGBUIJRL6pMs+E5bvu0+9M57Y22/3tsWTeWCxAWffXYzDcr7I5xFOQLduIvXrayOb9etLJr5eWOSbgEanTio25vuGDtVIPcfdd6vomH4Hn36q05pcI4Unx6PI7y/s/4LOmi3lvPnFTghYuXIlrd2d5cuX+3MCmzc7zrPPOk6tWo6Tk+M4l1/uOIsWOUFj48aNzpgxYxIfM8ZffznOPfc4zs47O86OOzpO796Os2qV4wkLF7JOOM7bbzuhu04Rw+trxG2lZsGOeBwrE2s4a3kmsRBiUZDXQup+8cWqf2UKMdsnutDHAcQoeF5I88npufPPSAhwX2YSj/qlGYYRTSxWU1i4EN+eBZtAM53amzSR2EImmv6NhExvvFHkggs0CE99WabCim7o0M2FGZEFVTVR+LCzadMmmThxorRs2VJKm3q2QH4f/ILs1V0yjhmwgrwuPK65c1WsgeTIbk4F6RSGi4PxMWiKyZmdf37xPSj3Gru5MCOycKv4OYghU3Crliu3JfG32BJRMOtz1okXWAgxN2RUkRbhaXE30sGC+i27M/8NsirCqc2bqzfG58UdG24emGFEky1bPPm1ZsBcWHxbtBC5/371uj75JEQaVZ8gdIgnRvEzEnxyheiH08UMmGFEk83ePNNmwIDRJgg0qOOaOtW8rlTB+5o9Ww1a06Zq1NLBQoiGEU02mwHzhtdfFznqKJ25RYUhC7CROqgyaad1xhmaD6PwJdVxMBgwOpbSZd8wjOiwerUnvza+BgwJ+F13iZxzjnaGZ05BUb1djMJhzsOLL4r07Styzz16bVORmpHZp6KT1laGYUSHpUs9+bU5sU0oUhZProuFlk4THnVLjh0YIVSbtMmePFmHJjHHLFkwYAy1NAwjOizJ1Qgyg+TEMhbbsaO2SaKRLUW6VkCbeRgdg7CDAmiMWLIKRfPADCN6LCumQrkAcmJnvC66SLtI4HW1b+/3GUUbxscgisEgIbdP5iZm0p8ZMMOIFksthFg8EBRceqnIq6+q7BvlnJGd5sB4YswfO+mkosOJFkI0jOh1NfojhTRCCuTERrBBN3UEBhzulGIjO+y/v8i774osWqRDkgoTdmDAGGHLTW8YRvhZ4k3+Kz4GbOBAkSefFHn6aZ17Zfjjibkz05gxVpDEnhCixze9YRhZxMOISvQN2Ntvi3Tvrso4+hsa/ubE8IDHjNk22Cg/DwwsD2YY0eBX757laBswdvv06WvVSuXyhv+0bq0DQanBy69jh+uBmQEzjGjw66+edTaKrgEjaUhXiOrVdddfsqTfZ2S44A2zsaCcgRZUece2cLObkMMwosHixZ5Nco+mASO/0ratyrYpqN15Z7/PyMgNdXfkIxnFwiYjt7ye/6tSxTwww4gKv/5qBiwlHnlEZOJErfXad1+/z8bIDzqfkAtbt04VormnO++9t8i8eX6enWEYmYJnuWZN8YLoGbDvvxe57TaRLl207sgIdgNghmG+9tr2+TDGsjBU1DCM8Lftmz1bpH59T359TiSLlVGy9evn99kYyUDn+nPPFencWeS337YZMFpQeVT8aBhGlmCq/dq1Iocc4smvz4lc6JCu8sOGRWNeeVzAC4NrrtFQIgYMZs3y9bQMwygmTLUHM2BJhg6vvVakWTO/z8ZIBRK8uUOJdeuK7LijhRENI+zMnKk6hJ128uTXR8eAkfMidMgsKiO8oUTexzVrtOjZ3b0ZhhFeA9bw74iKB0TDgDF36p13RPr3t9BhmHnwQZFVq7T1lwk5DCMaAo6GZsAKhpwJhbFNmliT3rDDRGzEHAMGaBjRhByGEe5OSGvXijRq5NlLhN+AkTch1ITq0AZThh/ymDk5Ih99pJ+bkMMwwsnMvyMohx/u2UvkhH5AJROVGdFx3HF+n42RCWgldfPNWoROONjCiIYRbgFHxYqevUS4Ddgzz2iVtwk3okXXriK77moGzDDCzExvBRzhNmAULT/wgE5WRrFmRAcMF6FECps//tjvszEMI4ACjnAbMCb8/vCDyHXX+X0mhhcwu61cOZGFC1XMYRhGePjoIxVwHHOMpy8TXgP2+ONa3d20qd9nYngBEwQuukj/TdNfwzDCw9ix2qCgcWNPXyacBmzRIr1AtB4y5WF0cb1rcp2GYYSHsWN1eC2KYg8JpwF78knNk1x4od9nYngJHvY++4h8843VgxlGmBr4zp0rcvrpnr9U+AzYxo0iQ4eKdOjgWX8tI2CKRHjqKb/PxDCMZHjzTc1fn3iieE34DNiUKSJLl4pcfrnfZ2Jkg06dREqWFHn2Wb/PxDCMZA0YxouhtR6TE8rYau3ang1IMwJGmTIiBx2kbWnwvg3DCC6//64jrbIQPgyfAaPvIdadi2PijfjQrp3W/eWe2mwYRvAYP16f1dNOy8rLhcuAffedKhCzZN2NgHDVVfqR3KdhGMGOkDVurKOtskC4DBgjU6gP8rg4zggY9FKrUUPkk0/UCzcMI3hs2CAyYUJWHYzwGbBTTtG8iBEvqClZv95aSxlGUJk6VeSvv8yAFcjXX4u0auX3WRh+cO21+vH//s/vMzEMIz/IUSOwO/hgyRbhMmDgcWsSI6DUqydStqzIuHEimzb5fTaGYeTmzz9FXnpJ5LLLsiqwC5cBo/sGk3qN+MFDQWdrQhRvvOH32RiGkZvnntMylyzX5+aErrWQx721jADD0NLSpUUGD/b7TAzDcEFYxTN5zjkiVapINgmfATPiCx4Y4UO6scyZ4/fZGIYBPI/0PqS5epZJy4ANGjRIateuLeXKlZMmTZrIp59+WuD3Dh06VI455hipVKlS4jjxxBML/f5CMQMWb9zheMjqGadjGIb/4H3RLceH8qaUDdioUaOkW7du0rt3b5k1a5bUr19fWrZsKb8xPTcfpk6dKm3btpUpU6bIRx99JDVr1pSTTjpJfvnll9TP1gxYvKlVS2S33bSN2PDhmg8zDMM/WMeZ1+fTaKtSqf7AwIEDpVOnTnLJJZckPh8yZIiMGzdOhg0bJrfeeuu/vv/FF1/c7vOnnnpKXnvtNZk8ebJ0oKN8PmzYsCFxuKxatSrxcVO1arLJFGgF4l6bvNdov/1KJUrnypcXWbdOpGPHrdK9+9akfy97k0svLSkLFpRI/J5HH90ixxzz74Ji7EmbNiVl1qwSsnmzyLJlm7f7//79c+T553MSv6NcOUceemirHHFEioXJiHgwZDNnquqJyc0Zuk7GNuwaJUfsr9PTT4vsuqvIBRcUqg726vqUcJzkWxts3LhRKlSoIK+++qqceeaZ/3y9Y8eO8ueff8obSajDVq9eLZUrV5ZXXnlFTiugX9add94pffr0+dfXR4wYkXh9IzU6dWohPXp8InXqrJLffy8n1157vNx55wypW/fPpH7+0UcbyO67r5O2befK99/vIv36NZYnnpgkpUptf+ts2pQjc+bsKjvuuFHuuONoGTFi/D//t2DBztK3bxN55JH3pHz5LTJ1ag0ZN66O9O//fsb/XsMwgsXatWulXbt2snLlStmZbkp+GLDFixdL9erVZcaMGXLkkUf+8/Xu3bvLtGnT5BNa/RTBNddcIxMnTpRvv/02kUNL1gMj9Pjrr7/KboSQjAJ3OZMmTZIWLVpIadR6uTywV17ZLA0a6OdNm5aUrl23yrnnJvfWV6pUSubM2fyPwOioo0rK3XdvlRNOyP/nf/pJ5IgjSm3ngX35JQX6peTTTzfLnnsSNs+RKVNKyCuvbEntj7zhBpFZs0R69xY5+2ztztKkSUauk7ENu0bJEevr9MYbOpfxo49EDjyw0G/9/fffpWrVqhk3YCmHEItDv379ZOTIkYm8WEHGC8qWLZs48sINErubJA3yu076NZ1KsmIF43pyEp2ZCsq7YmQmTtTpCHj/NWtu+317781mJifx+/J//W2v6dKokc6mrFu3dCLiwNv7/vt8T4pp2MqVRX7+WeSkk0SqVSOGLXL00an9jn/O0+6norBrlByxvE6DB+uDncRoK6+uTUoGbPfdd5eSJUvKUgZK5oLPqxSh/x8wYEDCgL377rty6KGHpne2Rtq0aaMldKhdH3pIZI899OtffJGd1//xR5HXXxf54Qe1O489pufE6KCUoMs19x+Bg6uvFrntNpG+fbGwHp25YRj/goUD+fyIEeInKW1/y5QpIw0bNkwIMFy2bt2a+Dx3SDEvDzzwgNx9990yYcIEaYTFNrLOqFFaOkXEDa0NbSVXr5ZEWDG/o2VL/TkitqVKiSxZsn2IEB1FKrz2mopIMV6ABujDD9OYUckvYN4QyhKmNROOyCdfahiGh7BxRFB13nniJynL6JHQU9s1fPhwmTNnjlx99dWyZs2af1SJKAt79Ojxz/fff//90rNnz4RKkdqxJUuWJI6/TALtC0z6xnG54w6RnXbSjVR+B+FDF+5RInXw2WeqnG3WLLXXrVNHDZb7tr/1lt7/KQ8WqFRJP/7xhxov/pBnnrHCZsPIFtOmibz9tsi99+ru1kdSfvU2bdrIsmXLpFevXglD1KBBg4RntSdJExFZuHCh5ORq9/T4448n1Ivnnnvudr+HOjLUhkZyEDFbu7bw7yFXtX59SVmzZlseyv1Z5PN8Hbp107A14bvDDiv6tXv1Umdn33319zJXEs+J4+67NarntkBDT7F8OcIbkerVRY49ltIJTVnNmCFy+OGa/6KtJQpc95ySZjO50Qoiq7ZIBUekBMMuiYliyHDzDMPwDhYTQjhE0mgd5TMpqRD9AhVixYoVZfny5bFVIbLQ77ij32cRLPDmMISJRqIdO+qssCQUiSjHxo8fL61atYpf4j1J7BolR+yu05gxImedJfLuuyInnJD0j6FCREORaRViuHohGkZ+XHihziBiZxj8/ZhhhJPNmzX3RR4iBePlJf4GMI2koX67qLQhu0Fq7GjtFdnd4OzZIsccLfLhDKlQ4W/5bsmSIvfdp5NgUam4ChTDMDLH889rrpmIR0AwAxYSaDOWCJcVkQMrV25L4vuiar+kNLLFtSI75Yjkbr1GV5emTdULa9HCxu4YRiahaJTmASi6AqQkt6fcCBeu7j6vhcbC9+unEkpqBgzDyByDBtG9QOSeeyRImAEzwsWyZfrRrcTODR058MRQJKZcYGYYRr6sXKkheqTG1L4ECDNgRrhgF0jxGP2o8oMHjbYfaP0Nwyg+/ftrHQ71NAHDDJgRLn79VQvPCpo9RLsPJPU9e+r3GoaRPog2BgzQRqZuG50AYQbMCJ8HhgErDB44vLQrrzRZvWEURzbPvL3atTUsH0DCZcBsMTLwqoraCVLsTu+rN99komq2zswwosXAgdo7jlZtTMMNIOEyYNY/0UjGAwMGrrZtK3LddRZKNIx0QofkvOg7V0ijdr8JlwFbsMDvMzD8ZMsWkfnzdSBZMjz6qIUSDaM4oUOanQaYcBkwZoAY8YVpnHQ0piNwMlgo0TAiGTp0MQNmhIeZM/VjsgYMLJRoGJELHYbTgH35pd9nYPhtwPbbT6RixdR+zkKJhhGp0GE4DdjcuSIbNvh9FoZffP65SMOGqf+chRINI1Khw3AaMLrVWhgxvgIO+hymY8DyhhKXLs302RlGuJkTrtBhOA0YDVwZ62vED8LHCDiOOCL930EokXHQdOowDENhfPq556q6NyShw3AaMHYGb73l91kYfsD7ziTX4uwOCSW+/vo2MYjlw4y4s3WryEUXifzvfyKjR4cmdBhOA3bSSSJTp2p3ZCNejB0rcsopKsYoDhjAhx7Sfz/1VEZOzTBCS8+eujl86SWRevUkbITLgDHKmjzYxIl+n4mRTdgd4jUxcTkTsOOEW24Ree+9zPxOwwgbL72k0xvuv1+kVSsJI+EyYDVqiNSvr7txIz6wQyxZUj2wTHLssTph1jq8GHFj5kyRSy/VzdxNN0lYCZcBA3bh48drzYIRD9iwYGwqVcrs7332WZ0rxj21enVmf7dhBJUlS0TOOEPk0EN1bl5Bo4lCQDgN2B9/iEyZ4veZGNmA95owX+vWmf/du+yixnHRIt2JktA2jCizYYPIWWfpvY5oo1w5CTPhM2DUAR10kBamGtEHL4mHrV07b37/AQdoLoAi5wBOnDWMjOE4IlddJTJ7tsiYMYEcUBl9A4a7e801Im+8ocl9I7pguB5/XGtU9tzTu9chgd2vn8i994qMGuXd6xiGn/zf/+mGkLBh48YSBcJnwIBwD/UKvBFGdJk8WeT773XD4jU33yxy4YUil1wiMn26969nGNlkzBiRG2/U+7x9e4kK4TRgFLTyJjz5pMrqjWgyeLDIIYeING2aHc+eDdF//qMe2aefev+ahpENxo8XOf98kXPOEenbV6JEOA0YXH21qmnYWRjRA2EFAgu8r2yppPDqeU3UWS1baq7AMMIexTj7bC1BoZE15SgRIrwGjJ35MceIPPigtQSKarx+hx00rJdNdtxRZNw4HdvSooXIN99k9/UNI1N88IGqtps3F3n5Ze0lGzHCa8AA1dgnn6igw4iW9/XYY9oZe6edsv/6zBuj20vNmtr9hUnQhhEmPv5YQ+FNmmj/T5pYR5BwGzAWlxNOELntNh23YUSDPn3UcGHA/IKi6XfeEdljD5FmzcwTM8LlebVoIdKggZaHhKxBb3wMGCB/ZpbNc8/5fSZGJuC9ZKDeHXeoWMdPMF4UzFMvc9xxIrNm+Xs+hpFMzuvkk1UmP2GChuEjTPgNWKNG2s+ud2+R9ev9PhujuGC4CN1RcBkEdt9dO4HUqSNy/PEasjaMIPL22yKnnqpt1+gfGnHjFQ0DBvfcI7J4scigQX6fiVEcMA7E6++6K1gxe8KJ774rcvDBGrYmRGMYQWLMGO1viPfFvyMcNoyeAatbV+SKK0TuvFPk55/9PhsjHajnw+ti2kC2lYfJQDiTkAwTockvWMjaCAKOIzJwoNZ4nXmmyCuvBGvz5zHRMGBAgR7NWTt1Mll9GGEu0ddfiwwbFtxaFST2hGnathXp2FG7Gph4yPCL9eu1c4zbYYOenhGUysfDgCF9ZsLupEnWYipsfPmlhoF79BA5/HAJNOxuMbJMdWbne9ppIn/+6fdZGXHj119VWETvzhdeUDFbUDd+HhIdAwZ0T7jsMt2RWCgxPKHDiy/WrvCMNw8DdAa54Qb1xqi3odZm7ly/z8qIC599puK1RYtE3n8/mCH3LBEtAwZ05iCUePnlFkoMU+iQLtllykioOOkk7ZmYk6NGjOJnw/CSESNUZYhS9/PPNScbY6JnwNxQIqoxQjxGcKHre1hChwVByym8sKOP1s4H3HO2cTIyDbnWW29Vb6tNG5GpU0WqVpW4Ez0D5oYSu3fXw3bFwWThQm0ySqf5sIQOC9s00c6MRDrhaxLrVpNoZIpVq1Qi37+/Rpgo9A/5JOVMEU0D5oamqIlgtzJvnt9nY+RmzRp9ICm0RPYbttBhfpBAJ5FOQp3EOgn2n37y+6yMsEMLM0b8EK1gLArt1bI1nSEERNeAsaAQL8bNpiOzKcWCAeE1PBQGVeK10K4pShDiIbGOSozCZyZKM1naMFJh82bdhDdsqAaLXCuRJSMmBswN7TDfaelSrd2xmh3/ufde9bqef17nbkUREusIUzBmzDOj8Nm8MSNVr4vQOh7XzJnarMGImQFzk+zMwqGzeOfOlmD3E4wWDyWtos46SyINnTueeELrEn/4wbwxIzWva906FQfRoMHyXTE2YMAOGGUiC0rXrmbE/IC8EPVelDfQsDcu0DvRvDGjKL79VuTII7f3umIukU+GeBgwIO8yeLBO+r3lFjNi2YQGvSzgHEOGxC8JndcbY5o418HuQQOvCy+LMpK1a83rSpH4GDC4+mqRhx9WOap5Ytlh5EiR888XOffcYPc5zKY31q6d3ot8bt5YfHG9LiIS5nWlRbwMGFx//TZPjEXEhB3eQb2K63m9+KJIqVJ+n1EwvTGELZQWGPFgxQqNApnXVWziZ8AAw4U3QNNfa8bqXdeASy/V6QAYsjh7XoV5Y+QEEbXss49urOgNaUQTjBW1grzXzC6kA415XcUingbMzYkx34ndD5JVK3bODCtXbt81AOUdvQKN/L0xutrTCJgany5dtKkxYzFMrRgd2JTgde+7r0ivXiLt24vMn6/zC83rKhbxXllQhFEgiKigcWM1aEb6UJzsdg0YN866BiRL7doiw4eLfPWVyIEHao4MKTX3o+VpwwubEEp4DjpIoz7HHy/y3/+KPPKIyJ57+n12kSDeBix3M1Z68p16qsiAAbZopAN1dmwCeGg/+UTbeBmpQa0YhfcffKBttk45RRc9rqcRLmgmzvNAKzvWmNmztc1YnTp+n1mkMAOWu2PHTTdpQ1bGZNg8seT46y8NfRECY6QIi+3++/t9VuGGzvYYsTffFFm+XL1aGh9/953fZ2YUBREd8ptEd5iOPG2aRiPq1/f7zCKJGTAXRAb3369hG9x8dsPErc0bKxhGOtAOCkEM5Qk0G2UWm1F8CL0iMPriC5HnnhOZNUtDUSyMY8Zo/ZARrEL9o47STdzixSKjR4vMmKGzuwzPMAOWFzwJepFdcIHIVVeZN1aQ13XttSLNm4vUqKG5G8oTTKzhzcaKpD9CD1pxIbenDdfee+sstSVL/D7D+EINX58++u8rrhCpUEHk1VdVXXrmmZb/zQK24hQUUkRi73pj1OqgqIv7jCe80ddeU6/r6afV68ILQ11leEvZsiIXXaS7erwx8mPUDjGZl80WHfAtWuA95HjfflukdWvNZ9GiDj77TPNe55xjJSNZxAxYMt4YO2AKD0nGEi6LY/jmvfc0PEJHDa6DeV3+cdhhIk8+KfLLLyo6QiDQrJlutKglYwCikVl+/11LQ7j3mbz9v//pe4BnDNYt3hds9UnGG6PokAQ6Me7LLlMPhDxEHHa87PYx5CecoCGRyZN1yrV5Xf5DvpFNBFECdv+IZ667TqR6dZErr1Rl6IYNfp9luAuPEXd16KDXlJZPqJU/+kifC4rQCRsavmEGLFnYYZGoJVTAzUwegj5m1Hls3CiRAsOMCu6887QeiRwgsX3KDZB1G8GCjQUbDMK75GWov2OTwcaDgaG8j+TP8CKMwmEQKekDhuDutpsW5aMspFsKXheCGlShlt8KBGbAUqVRI+1jx0FegjqPvfbSCntu8DCzerV2zsDDRD315ZcaJiGMSmzfHtrgg6gGYcGPP2qYl9D3woXqRVSurKFG8rkUnRu6WUN0QT9KQuTVqql4i44yiGQIEeLhdu8evenhEcAMWLpQ60GNBzc/NTq0BKKjAv8mdBOWPBkPMFJthn3y8FLTRZwfA82DSy9Da8IbPthskBO7/XatzUPazQgXQuKEwogo1KunCzOdU8Jyv2YCwqqEXAm3IsRgw0aPwlq11MP67Td9tm+80XJbAcdWpuJCvRg5Mh4AKu35N6GbSpVUKUYogq4ULBxBgZAnDyiFssT4CRHS2oYRMxgslG1GtKhaVd9bDnI7LOC897SwQpxATz6KbQkZc9ApPSrGik0mTXPdg8/pT8h9zvPJgWdKRMWIvgEbNGiQ9O/fX5YsWSL169eXRx99VBrTNqUAXnnlFenZs6f89NNPst9++8n9998vrVDyRImddtJ+Z4QfeEhc4zBihHowxx2nhakIQdjxZfNhoTs8zYqJ5VNsTHkASjV2nMiBeYA5vzJlsndOhn8gPHAXbmTh5HURJnDfojYljIxBo6kwuTU8Odew0auRDhNhM1ZI2ykE529gSgJiDAy2hcXjZcBGjRol3bp1kyFDhkiTJk3k4YcflpYtW8rcuXOlMjH2PMyYMUPatm0rffv2ldNOO01GjBghZ555psyaNUsOxnuJGjwQ5Mk4yEWQf8CYcdCmioeJBYC/3V0UOAhV0J08Ew8UcXvCgu5DjMya4mPgvGiZxeKFIbUHON5QBkHuhyN3LpR75o8/dPTHlClq1Ag3Y9i4bwi9EXLGs3M/uv9mM5fp+4rXZuwRIgvCoXzM/W9yeuRq8zNWfMRYlS+f2XMyfKeE46SmBcdoHXHEEfLYY48lPt+6davUrFlTrr32WrmVGVB5aNOmjaxZs0beeuutf772n//8Rxo0aJAwgvmxYcOGxOGycuVKqVWrlsybN0923XVXCS0UQs+Zo8l1doYc5JncoZo8YGwCqlTRj4T1UELhGfFQstjwdpGv4CDRTCcGYvZLl8qmlStlSr9+0rxzZym9bp2KS1hsONhFcwQplOkTmzZtkilTpkjz5s2ldFC9iaBdIzqAMEGYexdDgeHg3lu69N/DOLmPuXc53PuY5sTcw/wuIhL82y345f5372kOPidCsGyZ/n5eh3/nbSSAoeRZ4fcjXmFTiKFiJE2WxpTYvZQcK1askLp168qff/4pFTO5BjkpsGHDBqdkyZLO6NGjt/t6hw4dnNNPPz3fn6lZs6bz0EMPbfe1Xr16OYceemiBr9O7d2+Mqh122GGHHRE65s+f72SSlEKIy5cvly1btsieeWbZ8Pl/8STygTxZft/P1wuiR48eiTClC1Z7r732koULF2bWekeMVatWJbzhRYsWyc6EI418setUNHaNksOuU3K4UbRMR9ACqUIsW7Zs4sgLxstukqLhGtl1Khq7TkVj1yg57DolR06GW8+l9Nt23313KVmypCwlLp0LPq9CLDof+Hoq328YhmEYGTdgZcqUkYYNG8pk+uH9DSIOPj+Stkr5wNdzfz9MmjSpwO83DMMwDE9CiOSmOnbsKI0aNUrUfiGjR2V4ySWXJP6/Q4cOUr169YRsHq6//npp1qyZPPjgg3LqqafKyJEj5fPPP5cnaVGUJIQTe/funW9Y0diGXafksOtUNHaNksOuk7/XKWUZPSChdwuZkcM/8sgjCXk9HHfccVK7dm159tlntytkvuOOO/4pZH7ggQeiV8hsGIZhZJW0DJhhGIZh+I018zUMwzBCiRkwwzAMI5SYATMMwzBCiRkwwzAMI5QExoAxogX1Yrly5RKKxk8Z/VEIKBvr1auX+P5DDjlExjMmJAakcp2GDh0qxxxzjFSqVClxnHjiiUVe1yiQ6r3kQolHiRIlEtMS4kCq14mWbp07d5aqVasm5NA0Z43Dc5fqdaK0aP/995fy5csn2kx17dpV1udtRBwx3n//fWndurVUq1Yt8QyNGTOmyJ+ZOnWqHH744Yl7ad99991OuZ40TgAYOXKkU6ZMGWfYsGHOt99+63Tq1MnZZZddnKVLl+b7/R9++GGiqfADDzzgfPfdd84dd9zhlC5d2vn666+dKJPqdWrXrp0zaNAgZ/bs2c6cOXOciy++2KlYsaLzv//9z4kqqV4jlx9//NGpXr26c8wxxzhnnHGGE3VSvU408m7UqJHTqlUrZ/r06YnrNXXqVOeLL75wokyq1+nFF190ypYtm/jINZo4caJTtWpVp2vXrlk/92wyfvx45/bbb3def/31RNPevA3f87JgwQKnQoUKTrdu3RJr+KOPPppY0ydMmJDS6wbCgDVu3Njp3LnzP59v2bLFqVatmtO3b998v//88893Tj311O2+1qRJE+fKK690okyq1ykvmzdvdnbaaSdn+PDhTlRJ5xpxXY466ijnqaeecjp27BgLA5bqdXr88cedOnXqOBs3bnTiRKrXie89/vjjt/sai3TTpk2duCBJGLDu3bs7Bx100HZfa9OmjdOyZcuUXsv3EOLGjRtl5syZifBW7oaPfP4RU2Lzga/n/n5gqGZB3x8F0rlOeVm7dm1iflGoZ6p5cI3uuuuuxDDWyy67TOJAOtdp7NixifZvhBCZJsEw2vvuuy8xnSKqpHOdjjrqqMTPuGHGBQsWJMKs1rjBmzXc92702RrREnbSuU55ueWWWxIx6rw3Tpyv0fTp0+Xpp5+WL5hgHRPSuU4sxO+9955ceOGFiQX5hx9+kGuuuSaxIaJFUBRJ5zq1a9cu8XNHH3000S3ZvHmzXHXVVXLbbbdl6azDQUFrOONp1q1bl8gfJoPvHpiRHfr165cQKYwePTqRjDZEVq9eLe3bt0+IXZi0YBQMTbvxUulhSkNvJq3ffvvtBU5VjysIE/BMBw8eLLNmzZLXX39dxo0bJ3fffbffpxZJfPfAbESLd9fJZcCAAQkD9u6778qhhx4qUSXVazR//vxEf07UU7kXaihVqpTMnTtX9tlnH4ka6dxLKA9Lly6d+DmXAw44ILGTJtTGpIqokc516tmzZ2JTdPnllyc+RyFNs/MrrrgiYfAzPQ8rrBS0hjNTLVnvC3y/mjaixbvrBDROZvc3YcKExASBKJPqNaIM4+uvv06ED93j9NNPl+bNmyf+jQQ6iqRzLzVt2jQRNnQNPMybNy9h2KJovNK9TuSZ8xop1+hb21kP1nAnIFJVpKfPPvtsQlJ5xRVXJKSqS5YsSfx/+/btnVtvvXU7GX2pUqWcAQMGJOThvXv3jo2MPpXr1K9fv4QE+NVXX3V+/fXXf47Vq1c7USXVa5SXuKgQU71OCxcuTChYu3Tp4sydO9d56623nMqVKzv33HOPE2VSvU6sRVynl156KSEVf+edd5x99tknoZyOMqtXr06U63BgVgYOHJj4988//5z4f64R1yqvjP7mm29OrOGU+4RWRg/UAdSqVSux4CJd/fjjj//5v2bNmiUWlty8/PLLTt26dRPfjxxz3LhxThxI5TrttddeiZsp78FDFmVSvZfiaMDSuU4zZsxIlKuwoCOpv/feexMlCFEnleu0adMm584770wYrXLlyjk1a9Z0rrnmGuePP/5wosyUKVPyXWvca8NHrlXen2nQoEHiunI/PfPMMym/ro1TMQzDMEKJ7zkwwzAMw0gHM2CGYRhGKDEDZhiGYYQSM2CGYRhGKDEDZhiGYYQSM2CGYRhGKDEDZhiGYYQSM2CGYRhGKDEDZhiGYYQSM2CGYRhGKDEDZhiGYUgY+X8XOVtkBnY4AgAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print(f\"N = {N}, Sum of radii = {p.attributes.objval}\")\n",
    "\n",
    "# Create a plot\n",
    "fig, ax = plt.subplots()\n",
    "\n",
    "# Plot each circle\n",
    "for (x, y, r) in zip(xsol, ysol, rsol):\n",
    "    # Draw the circle with red edge\n",
    "    circle = plt.Circle((x, y), r, edgecolor='red', facecolor='none')\n",
    "    ax.add_patch(circle)\n",
    "\n",
    "    # Draw the radius line in blue\n",
    "    ax.plot([x, x + r], [y, y], color='blue')\n",
    "\n",
    "    # Add the radius length as label\n",
    "    ax.text(x + r / 2, y, f'R={r:.2f}', color='blue', fontsize=8, ha='center', va='bottom')\n",
    "\n",
    "# Set the limits and aspect ratio\n",
    "ax.set_xlim(0, 1)\n",
    "ax.set_ylim(0, 1)\n",
    "ax.set_aspect('equal', 'box')\n",
    "ax.grid(True)\n",
    "\n",
    "# Show the plot\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.12.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
