{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Optimal control problem 2\n",
    "\n",
    "B. A. Conway and K. M. Larson : Collocation Versus Differential Inclusion in Direct Optimization"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use single shooting for transcription"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import simulatetraj as a\n",
    "import casadi as cs\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* define the integrator\n",
    "* integrator grid and transcription grid will match for single shooting"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "sym_dict=a.get_symbols(n_x=2,n_u=1,n_p=0,N=50,cs_type=cs.SX)\n",
    "x=sym_dict['x']\n",
    "u=sym_dict['u']\n",
    "f=cs.vertcat(x[1],-x[1]+u)\n",
    "_=a.scale_ode(sym_dict=sym_dict,f=f)\n",
    "_=a.integrator(sym_dict=sym_dict,options=['rk',{}])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "******************************************************************************\n",
      "This program contains Ipopt, a library for large-scale nonlinear optimization.\n",
      " Ipopt is released as open source code under the Eclipse Public License (EPL).\n",
      "         For more information visit https://github.com/coin-or/Ipopt\n",
      "******************************************************************************\n",
      "\n",
      "This is Ipopt version 3.14.11, running with linear solver MUMPS 5.4.1.\n",
      "\n",
      "Number of nonzeros in equality constraint Jacobian...:       52\n",
      "Number of nonzeros in inequality constraint Jacobian.:        0\n",
      "Number of nonzeros in Lagrangian Hessian.............:       50\n",
      "\n",
      "Total number of variables............................:       52\n",
      "                     variables with only lower bounds:        0\n",
      "                variables with lower and upper bounds:        0\n",
      "                     variables with only upper bounds:        0\n",
      "Total number of equality constraints.................:        3\n",
      "Total number of inequality constraints...............:        0\n",
      "        inequality constraints with only lower bounds:        0\n",
      "   inequality constraints with lower and upper bounds:        0\n",
      "        inequality constraints with only upper bounds:        0\n",
      "\n",
      "iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls\n",
      "   0  6.7346939e-01 4.81e-02 2.75e-02  -1.0 0.00e+00    -  0.00e+00 0.00e+00   0\n",
      "   1  5.7790126e-01 4.44e-16 1.39e-17  -2.5 3.43e-01    -  1.00e+00 1.00e+00f  1\n",
      "\n",
      "Number of Iterations....: 1\n",
      "\n",
      "                                   (scaled)                 (unscaled)\n",
      "Objective...............:   5.7790126446789747e-01    5.7790126446789747e-01\n",
      "Dual infeasibility......:   1.3877787807814457e-17    1.3877787807814457e-17\n",
      "Constraint violation....:   4.4408920985006262e-16    4.4408920985006262e-16\n",
      "Variable bound violation:   0.0000000000000000e+00    0.0000000000000000e+00\n",
      "Complementarity.........:   0.0000000000000000e+00    0.0000000000000000e+00\n",
      "Overall NLP error.......:   4.4408920985006262e-16    4.4408920985006262e-16\n",
      "\n",
      "\n",
      "Number of objective function evaluations             = 2\n",
      "Number of objective gradient evaluations             = 2\n",
      "Number of equality constraint evaluations            = 2\n",
      "Number of inequality constraint evaluations          = 0\n",
      "Number of equality constraint Jacobian evaluations   = 2\n",
      "Number of inequality constraint Jacobian evaluations = 0\n",
      "Number of Lagrangian Hessian evaluations             = 1\n",
      "Total seconds in IPOPT                               = 0.005\n",
      "\n",
      "EXIT: Optimal Solution Found.\n",
      "      solver  :   t_proc      (avg)   t_wall      (avg)    n_eval\n",
      "       nlp_f  |        0 (       0)   3.00us (  1.50us)         2\n",
      "       nlp_g  |        0 (       0) 184.00us ( 92.00us)         2\n",
      "  nlp_grad_f  |        0 (       0)  29.00us (  9.67us)         3\n",
      "  nlp_hess_l  |   1.00ms (  1.00ms) 878.00us (878.00us)         1\n",
      "   nlp_jac_g  |   1.00ms (333.33us)   1.19ms (396.33us)         3\n",
      "       total  |   9.00ms (  9.00ms)   5.84ms (  5.84ms)         1\n",
      "Objective value: Paper - 0.577678 | code - 0.5779012644678975\n"
     ]
    }
   ],
   "source": [
    "nlp=cs.Opti()\n",
    "X0=nlp.variable(sym_dict['n_x'],1)\n",
    "U =nlp.variable(sym_dict['n_u'],sym_dict['N'])\n",
    "obj=cs.sumsqr(U)*2/sym_dict['N']\n",
    "\n",
    "x0=cs.DM([0,0])\n",
    "a.simulate(sym_dict=sym_dict,t0=0,tf=2,x0=x0,u=U)\n",
    "X=sym_dict['x_res']\n",
    "\n",
    "nlp.minimize(obj)\n",
    "nlp.subject_to(X0-x0==0)\n",
    "nlp.subject_to(cs.horzcat(1,-2.694528)@X[:,-1]+1.155356==0)\n",
    "nlp.solver('ipopt')\n",
    "nlp.set_initial(U,cs.np.linspace(0,1,sym_dict['N']))\n",
    "sol=nlp.solve()\n",
    "obj_val=sol.value(nlp.f)\n",
    "print(f'Objective value: Paper - 0.577678 | code - {obj_val}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* plot results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGhCAYAAACzurT/AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy80BEi2AAAACXBIWXMAAA9hAAAPYQGoP6dpAABcBUlEQVR4nO3deXwU9f3H8dfuZnMfHIEEQiDc9yEgEU+QS7EqtiqKB6Jiq9JqaWultiDqT9RSxYOKF2JRBK1XVUQxGkVB0ABy3zchF0fubDa78/tjIBC5kpBldjfv5+Oxj52dzMx+vpkk+87Md75jMwzDQERERMQidqsLEBERkfpNYUREREQspTAiIiIillIYEREREUspjIiIiIilFEZERETEUgojIiIiYimFEREREbGUwoiIiIhYSmFERERELFWrMDJ9+nRSUlIIDw8nNTWVZcuWnXTZWbNmYbPZqjzCw8NrXbCIiIgElxqHkXnz5jF+/HgmTZrE8uXL6dmzJ8OGDSMnJ+ek68TGxrJv377Kx86dO8+oaBEREQketpreKC81NZVzzz2XF154AQCv10tycjK///3vefDBB49bftasWdx///0cOnSo1kV6vV4yMzOJiYnBZrPVejsiIiJy9hiGQWFhIc2bN8duP/nxj5CabLS8vJyMjAwmTJhQOc9utzN48GCWLFly0vWKiopo1aoVXq+X3r178/jjj9O1a9eTLu9yuXC5XJWv9+7dS5cuXWpSqoiIiPiJ3bt306JFi5N+vUZhJC8vD4/HQ0JCQpX5CQkJbNiw4YTrdOzYkZkzZ9KjRw/y8/OZOnUq559/PmvXrj1pYVOmTGHy5MnHzX/11VeJjIysSckiIiJikZKSEu68805iYmJOuVyNTtNkZmaSlJTE4sWL6d+/f+X8Bx54gG+++YalS5eedhtut5vOnTtz44038uijj55wmV8eGSkoKCA5OZm8vDxiY2OrW261alm4cCFDhgzB6XTW2Xb9SbC3Ue0LfMHeRrUv8AV7G33ZvoKCAuLj48nPzz/l53eNjozEx8fjcDjIzs6uMj87O5vExMRqbcPpdHLOOeewZcuWky4TFhZGWFjYCdf1xQ+Cr7brT4K9jWpf4Av2Nqp9gS/Y2+iL9lV3ezW6miY0NJQ+ffqQlpZWOc/r9ZKWllblSMmpeDweVq9eTbNmzWry1iIiIhKkanRkBGD8+PGMHj2avn370q9fP6ZNm0ZxcTFjxowB4NZbbyUpKYkpU6YA8Mgjj3DeeefRrl07Dh06xD//+U927tzJnXfeWbctERERkYBU4zAycuRIcnNzmThxIllZWfTq1YsFCxZUdmrdtWtXlct3Dh48yNixY8nKyqJhw4b06dOHxYsX1/nVMR6PB7fbXaN13G43ISEhlJWV4fF46rQef1FRUYHdbqeGV3CLiIicNTUOIwDjxo1j3LhxJ/xaenp6ldfPPPMMzzzzTG3eptqKiorYs2dPjT9wDcMgMTGR3bt3B+34JYZh0KxZM/bu3UtSUhKhoaFWlyQiIlJFrcKIP/F4POzZs4fIyEiaNGlSo1Dh9XopKioiOjr6lIOxBDKPx0N+fj7FxcVs376d9u3bB21bRUQkMAV8GHG73RiGQZMmTYiIiKjRul6vl/LycsLDw4P2A9rr9eJ2u4mNjWX37t2V7RUREfEXQfMJHKynWepKsIYtEREJfPqEEhEREUspjIiIiIilFEZERETEUgojAWjfvn2MGjWKDh06YLfbuf/++60uSUREpNYURgKQy+WiSZMm/P3vf6dnz55WlyMiIgHM/v00Ome+C4X7rKvBsnf2EcMwKCmvqPajtNxTo+VP9ajJoGu5ubkkJiby+OOPV85bvHgxoaGhVe79cyIpKSk8++yz3HrrrcTFxdX6eyUiIvVceQn2H16gQ/bH2LJWWVZGwI8z8kulbg9dJn5uyXuve2QYkaHV+5Y2adKEmTNnMmLECIYOHUrHjh255ZZbGDduHIMGDfJxpSIiIsDqd7CVHaI4tCmhbQdbVkbQhZFAMnz4cMaOHctNN91E3759iYqKqrzBoIiIiE8ZBix9GYDt8YPoaHdYVkrQhZEIp4N1jwyr1rJer5fCgkJiYmPqZFCwCGfNd+TUqVPp1q0b7777LhkZGYSFhZ1xHSIiIqe183vIWYvhjGRX44vpaGEpQRdGbDZbtU+VeL1eKkIdRIaGWDZC6datW8nMzMTr9bJjxw66d+9uSR0iIlLPLH0JAG+3a3ETZWkpQRdGAkl5eTk333wzI0eOpGPHjtx5552sXr2apk2bWl2aiIgEs/w9sOFTALx974SfdlhajsKIhR566CHy8/N57rnniI6OZv78+dx+++188sknp1135cqVABQVFZGbm8vKlSsJDQ2lS5cuPq5aREQC3o+vgeGBlIugaRdgh6XlKIxYJD09nWnTpvH1118TGxsLwOzZs+nZsycvvvgid9999ynXP+eccyqnMzIymDNnDq1atWLHjh2+LFtERAKduwyWv2FO97vL2loOUxixyIABA3C73VXmpaSkkJ+fX631azKmiYiISKW170PJfohtAR2Hg9f6z5OgG/RMRERETsIwKjuucu4d4PCPYxIKI36oa9euREdHn/Dx1ltvWV2eiIgEqj0/wr6V4AiD3qOtrqaSf0QiqWL+/PnHncI5IiEh4SxXIyIiQePIUZHu10FUY2trOYbCiB9q1aqV1SWIiEiwKcyCdR+a06n+0XH1CJ2mERERqQ9+eh28FZB8HjTzrzu+K4yIiIgEu4pyyHjdnPazoyKgMCIiIhL81n0ERdkQ0ww6X2V1NcdRGBEREQl2yw53XO17Ozic1tZyAgojIiIiwWzvcvOSXrvTry7nPZbCiIiISDBb9rL53PUaiPHP4SEURgLQ+++/z5AhQ2jSpAmxsbH079+fzz//3OqyRETE3xTnwZr3zOnU31pbyykojASgb7/9liFDhjB//nwyMjIYOHAgV155JStWrLC6NBER8ScZs8BTDs17Q4u+VldzUsE36JlhgLukest6veay5Q6w10Euc0aCzVatRXNzc+nevTt/+MMf+Nvf/gbA4sWLGTBgAJ999hmDBg066brTpk2r8vrxxx/no48+4uOPP65yN18REanHPG748TVz2o+PikAwhhF3CTzevFqL2oEGdfnef8uE0KhqLdqkSRNmzpzJiBEjGDp0KB07duSWW25h3LhxpwwiJ+L1eiksLKRRo0a1qVpERILR2g+hMBOiE8z+In4s+MJIABk+fDhjx47lpptuom/fvkRFRTFlypQab2fq1KkUFRVx/fXX+6BKEREJOIYBS543p/uNhZAwa+s5jeALI85I8whFNXi9XgoKC4mNicFeV6dpamjq1Kl069aNd999l4yMDMLCavYDM2fOHCZPnsxHH31E06ZNa/z+IiIShHYuhn0/Q0g49Lnd6mpOK/jCiM1W7VMleL3g9JjL10UYqYWtW7eSmZmJ1+tlx44ddO/evdrrzp07lzvvvJN3332XwYMH+7BKEREJKEumm889b/Sru/OeTPCFkQBSXl7OzTffzMiRI+nYsSN33nknq1evrtYRjrfffpvbb7+duXPncsUVV5yFakVEJCDs3wob55vT591jbS3VpDBioYceeoj8/Hyee+45oqOjmT9/PrfffjuffPLJKdebM2cOo0eP5tlnnyU1NZWsrCwAIiIiiIuLOxuli4iIv/rhRcCA9sOgSQerq6kWjTNikfT0dKZNm8bs2bOJjY3Fbrcze/ZsFi1axIsvvnjKdV9++WUqKiq49957adasWeXjvvvuO0vVi4iIXyo5ACvfMqf732ttLTWgIyMWGTBgAG63u8q8lJQU8vPzT7tuenq6j6oSEZGAljHLHOIioTu0vtjqaqpNR0ZERESCQUX50fvQ9L+32oNw+gOFET/UtWtXoqOjT/h46623rC5PRET80doPoHAfRCdCt99YXU2N6DSNH5o/f/5xp3COSEjwzzsuioiIhQwDlrxgTvcbCyGh1tZTQwojfqhVq1ZWlyAiIoFkx3eQtQpCIqCv/w9y9ktBc5rGMAyrS/Br+v6IiASxI4Oc9RoFkYF3n7KADyMOhwMwBxCTkyspMe9k7HQ6La5ERETqVN4W2PSZOR0gg5z9UsCfpgkJCSEyMpLc3FycTmeN7jHj9XopLy+nrKysbu5N44c8Hg+FhYUUFhbSsGHDyvAmIiJB4od/m88dLof4dtbWUksBH0ZsNhvNmjVj+/bt7Ny5s0brGoZBaWkpERER2ALoEqiaMAyD4uJimjVrRmJiotXliIhIXSo5ACvnmNP9A/OoCARBGAEIDQ2lffv2NT5V43a7+fbbb7n44ouD9vRFRUUFX331Fb169QrawCUiUm/9NBMqSiGxO6RcZHU1tRYUYQTAbrcTHh5eo3UcDgcVFRWEh4cHbRhxu93qvCoiEowqXMcMcjYuoAY5+6Xg7CghIiIS7Na8D0XZ5iBnXX9tdTVnRGFEREQk0BjG0ct5U+8KuEHOfklhREREJNBsTYPs1eCMhD5jrK7mjNUqjEyfPp2UlBTCw8NJTU1l2bJl1Vpv7ty52Gw2RowYUZu3FREREYDvppnPfW4LyEHOfqnGYWTevHmMHz+eSZMmsXz5cnr27MmwYcPIyck55Xo7duzgz3/+MxddFLi9fUVERCy3JwN2LAJ7iHl33iBQ4zDy9NNPM3bsWMaMGUOXLl2YMWMGkZGRzJw586TreDwebrrpJiZPnkybNm3OqGAREZF67ftnzOfu10FcC2trqSM1urS3vLycjIwMJkyYUDnPbrczePBglixZctL1HnnkEZo2bcodd9zBokWLTvs+LpcLl8tV+bqgoAAwL1M92d1sa+PItupym/4m2Nuo9gW+YG+j2hf4/KqN+zcTsv4TbIA79V6og5p82b7qbrNGYSQvLw+Px3PcbewTEhLYsGHDCdf57rvveO2111i5cmW132fKlClMnjz5uPlffPEFkZGRNSm5WhYuXFjn2/Q3wd5GtS/wBXsb1b7A5w9t7LnrNVIw2Bd7Dst+3AZsq7Nt+6J9R+6Ldjo+HfSssLCQW265hVdeeYX4+PhqrzdhwgTGjx9f+bqgoIDk5GSGDh1KbGxsndXndrtZuHAhQ4YMCepBz4K5jWpf4Av2Nqp9gc9v2li4j5Dp5lmIJiMeY3hyap1s1pftO3Jm43RqFEbi4+NxOBxkZ2dXmZ+dnX3C+55s3bqVHTt2cOWVV1bO83q95huHhLBx40batm173HphYWGEhYUdN9/pdPrkB8FX2/Unwd5GtS/wBXsb1b7AZ3kbf3oFPOWQfB4hbS6s8837on3V3V6NOrCGhobSp08f0tLSKud5vV7S0tLo37//cct36tSJ1atXs3LlysrHVVddxcCBA1m5ciXJyck1eXsREZH6qfQQ/PS6OX3hHy0txRdqfJpm/PjxjB49mr59+9KvXz+mTZtGcXExY8aYg67ceuutJCUlMWXKFMLDw+nWrVuV9Rs0aABw3HwRERE5iZ9mQnkhNOkM7YdaXU2dq3EYGTlyJLm5uUycOJGsrCx69erFggULKju17tq1C7tdA7uKiIjUCXcZ/PCiOX3BfRCEn7G16sA6btw4xo0bd8Kvpaenn3LdWbNm1eYtRURE6qef50BxDsS2gO7XWl2NTwRfvBIREQkWXg8sft6cPn8cOIKzk7DCiIiIiL9a/z84sA0iGkLvW62uxmcURkRERPyRYRy9IV6/uyA0ytJyfElhRERExB9tS4d9KyEkAvr91upqfEphRERExB99P8187n0rRDW2tBRfUxgRERHxN5krzCMjNgf0v9fqanxOYURERMTffP+s+dztN9CwlbW1nAUKIyIiIv4kbzOs/dCcvuA+S0s5WxRGRERE/MmifwEGdBwOifXj1ikKIyIiIv7iwHZY9Y45ffFfrK3lLFIYERER8RffPQ2GB9oNhqTeVldz1iiMiIiI+INDu2Hl2+b0xQ9YW8tZpjAiIiLiD76fBl43tL4YWqZaXc1ZpTAiIiJitYJ9sHy2OV3PjoqAwoiIiIj1Fj8HHhe07A8pF1pdzVmnMCIiImKlolz46XVz+uI/g81mbT0WUBgRERGx0pIXoKIUmveGtoOsrsYSCiMiIiJWKTkAP75qTl/yQL08KgIKIyIiItb54UUoL4LE7tDhMqursYzCiIiIiBXK8mHpS+b0xX+pt0dFQGFERETEGktfBlc+NOkMna60uhpLKYyIiIicba5C+GG6OX3xn8Fevz+O63frRURErPDja1B6EBq3g67XWF2N5RRGREREzqbyEvNyXoCL/gR2h7X1+AGFERERkbNp+RtQnAsNWkL366yuxi8ojIiIiJwt7lL4/llz+sLx4HBaW4+fUBgRERE5W356HQr3QVwy9BpldTV+Q2FERETkbCgvhu+eNqcv/guEhFlbjx9RGBERETkblr1i9hVpmKKjIr+gMCIiIuJrrsKjfUUu+av6ivyCwoiIiIivLZ0BpQfMcUW6X291NX5HYURERMSXSg/B4ufN6QETwBFiaTn+SGFERETEl374t3lTvCadoeuvra7GLymMiIiI+ErJAVjyb3N64IR6fw+ak9F3RURExFcWPwflhZDQvd7fmfdUFEZERER8oSgXlr5kTg/8m46KnIK+MyIiIr7w/TRwl0Dz3tDxcqur8WsKIyIiInWtMAt+fNWcHvgQ2GzW1uPnFEZERETq2qKnoaIMklOh3SCrq/F7CiMiIiJ1KX8PZLxuTuuoSLUojIiIiNSlb6eCpxxaXQitL7a6moCgMCIiIlJXDu6EFbPN6Ut1VKS6FEZERETqyrdPgbcC2gyEVudbXU3AUBgRERGpC7mbYOUcc3rgQ9bWEmAURkREROrCV4+A4YWOwyH5XKurCSgKIyIiImdqz0+w/mOw2WHQRKurCTgKIyIiImfCMODLh83pnqOgaWdLywlECiMiIiJnYksa7FgEjjDzzrxSYwojIiIiteX1Hj0q0m8sxLWwtJxApTAiIiJSW2veg+zVEBYLF/3J6moClsKIiIhIbVSUw1ePmtMX3AeRjaytJ4DVKoxMnz6dlJQUwsPDSU1NZdmyZSdd9v3336dv3740aNCAqKgoevXqxezZs2tdsIiIiF/IeB0O7YToBDjvbqurCWg1DiPz5s1j/PjxTJo0ieXLl9OzZ0+GDRtGTk7OCZdv1KgRDz30EEuWLGHVqlWMGTOGMWPG8Pnnn59x8SIiIpZwFcI3T5nTl/wVQqOsrSfA1TiMPP3004wdO5YxY8bQpUsXZsyYQWRkJDNnzjzh8gMGDOCaa66hc+fOtG3blvvuu48ePXrw3XffnXHxIiIillgyHUryoFFb6H2r1dUEvJCaLFxeXk5GRgYTJhy9dMlutzN48GCWLFly2vUNw+Crr75i48aNPPnkkyddzuVy4XK5Kl8XFBQA4Ha7cbvdNSn5lI5sqy636W+CvY1qX+AL9jaqfYHvuDYW5xKy+DlsQMWAv2F4AW/gtt+X+7C627QZhmFUd6OZmZkkJSWxePFi+vfvXzn/gQce4JtvvmHp0qUnXC8/P5+kpCRcLhcOh4N///vf3H777Sd9n4cffpjJkycfN3/OnDlERkZWt1wREZE6123Pm7TN/YKDka35tsPDujPvKZSUlDBq1Cjy8/OJjY096XI1OjJSWzExMaxcuZKioiLS0tIYP348bdq0YcCAASdcfsKECYwfP77ydUFBAcnJyQwdOvSUjakpt9vNwoULGTJkCE6ns86260+CvY1qX+AL9jaqfYGvShuLMwl58WsAYkb8i+GtL7a4ujPny3145MzG6dQojMTHx+NwOMjOzq4yPzs7m8TExJOuZ7fbadeuHQC9evVi/fr1TJky5aRhJCwsjLCwsOPmO51On/yw+2q7/iTY26j2Bb5gb6PaF/icTifORU+Zp2TaDCSkwyCrS6pTvtiH1d1ejTqwhoaG0qdPH9LS0irneb1e0tLSqpy2OR2v11ulT4iIiIjfy14Dq94xpwc/bGkpwabGp2nGjx/P6NGj6du3L/369WPatGkUFxczZswYAG699VaSkpKYMmUKAFOmTKFv3760bdsWl8vF/PnzmT17Ni+++GLdtkRERMSHHF8/BhjQ9dfQvJfV5QSVGoeRkSNHkpuby8SJE8nKyqJXr14sWLCAhIQEAHbt2oXdfvSAS3FxMffccw979uwhIiKCTp068eabbzJy5Mi6a4WIiIgPNSlYg33rl2APgUv/bnU5QadWHVjHjRvHuHHjTvi19PT0Kq8fe+wxHnvssdq8jYiIiPW8HrrufducPncsNG5rbT1BSPemEREROQXbqreJK9uNER4HlzxgdTlBSWFERETkZFyFONIfB8B74Z91MzwfURgRERE5me+fxVacQ1FYAt6+d1hdTdBSGBERETmR/D2w+HkA1jUfCY5QiwsKXgojIiIiJ5L2KFSU4W3Zn31xfayuJqgpjIiIiPzS3uWwai4AnsGP6v4zPqYwIiIicizDgM8fMqd73ADNellaTn2gMCIiInKs9R/DrsUQEgGD/mF1NfWCwoiIiMgRFS5YONGcPv/3ENfC2nrqCYURERGRI5a9Age3Q3QCXHCf1dXUGwojIiIiACUH4NunzOlL/wFh0dbWU48ojIiIiAB88ySU5UNCd+g1yupq6hWFERERkbzN8OOr5vSwx8DusLaeekZhREREZOFE8FZAh8ugzQCrq6l3FEZERKR+25IGG+eDzQFDHrW6mnpJYUREROqvinL47K/mdOpvoUkHa+uppxRGRESk/lo6A/ZvhqgmMOBBq6uptxRGRESkfirMMq+gARg8GcLjrK2nHlMYERGR+mnhJCgvghbnQs8bra6mXlMYERGR+mfX0sN35bXB5U+BXR+HVtJ3X0RE6hevB+b/2ZzufQsk9ba2HlEYERGRemb5G5C1yuwjMmiS1dUICiMiIlKflByAtEfM6YEPQVS8tfUIoDAiIiL1yVePQelBaNoV+t5hdTVymMKIiIjUD/tWQcbr5vTwp8ARYm09UklhREREgp9hwGcPgOGFrr+GlAutrkiOoTAiIiLBb/W7sGsJOCNh6GNWVyO/oDAiIiLBzVUIX/zDnL7oTxCXZG09chyFERERCW7f/hOKsqBhazj/91ZXIyegMCIiIsErdyMs+bc5fdkTEBJmbT1yQgojIiISnAwDPvkjeN3Q4TLoeJnVFclJKIyIiEhwWjkHdn5vdlod/k+rq5FTUBgREZHgU7wfvvi7OT3gQWjQ0tp65JQURkREJPgsnAilB8yRVs+7x+pq5DQURkREJLjs+B5WvmlOXzkNHE5Ly5HTUxgREZHgUVFudloF6HMbJPeztBypHoUREREJHoufg7yNENUEBj9sdTVSTQojIiISHA5sNwc4Axj2OEQ0tLYeqTaFERERCXyGAfP/DBVl0PoS6H6d1RVJDSiMiIhI4Fv7AWz5EhyhcMXTYLNZXZHUgMKIiIgEtrJ8WPCgOX3RnyC+nbX1SI0pjIiISGD76jEoyoZGbeGC+62uRmpBYURERALX3gxY9oo5/atnwBlubT1SKwojIiISmDwV8PH9gAE9RkKbS6yuSGpJYURERALTkuchaxWEx8HQx6yuRs6AwoiIiASevC3w9RRzetgUiG5qbT1yRhRGREQksHi98L9x4HFB20HQa5TVFckZUhgREZHA8tNrsGsJhEabN8LTmCIBT2FEREQCx6Fd8OXD5vTgh6FBSyurkTqiMCIiIoHBMODj+6C8CFr2h753WF2R1BGFERERCQw/vw1bvwJHGFz1PNj1ERYsarUnp0+fTkpKCuHh4aSmprJs2bKTLvvKK69w0UUX0bBhQxo2bMjgwYNPubyIiMhxCrNhwQRzeuAEiG9vbT1Sp2ocRubNm8f48eOZNGkSy5cvp2fPngwbNoycnJwTLp+ens6NN97I119/zZIlS0hOTmbo0KHs3bv3jIsXEZF6Yv6foOwQNOsF/X9vdTVSx2ocRp5++mnGjh3LmDFj6NKlCzNmzCAyMpKZM2eecPm33nqLe+65h169etGpUydeffVVvF4vaWlpZ1y8iIjUA2s/hPUfgz0Ern4BHCFWVyR1rEZ7tLy8nIyMDCZMmFA5z263M3jwYJYsWVKtbZSUlOB2u2nUqNFJl3G5XLhcrsrXBQUFALjdbtxud01KPqUj26rLbfqbYG+j2hf4gr2Nat8ZKjlAyPw/YwM8/e/D27gTnOXvpfbhmW/7dGyGYRjV3WhmZiZJSUksXryY/v37V85/4IEH+Oabb1i6dOlpt3HPPffw+eefs3btWsLDT3xDo4cffpjJkycfN3/OnDlERkZWt1wREQlw5+x8iZYHvqcgPIlvOj6C1+60uiSpgZKSEkaNGkV+fj6xsbEnXe6sHut64oknmDt3Lunp6ScNIgATJkxg/Pjxla8LCgoq+5qcqjE15Xa7WbhwIUOGDMHpDM4f8GBvo9oX+IK9jWpf7dm2fEnIiu8xsBF5w2tcltS3TrdfXdqHtXfkzMbp1CiMxMfH43A4yM7OrjI/OzubxMTEU647depUnnjiCb788kt69OhxymXDwsIICws7br7T6fTJD4KvtutPgr2Nal/gC/Y2qn01VJYPn/0ZANt59xCS0v80K/ie9mHttlkdNerAGhoaSp8+fap0Pj3SGfXY0za/9NRTT/Hoo4+yYMEC+va1JtmKiEgAWTABCvZAwxS49CGrqxEfq/FpmvHjxzN69Gj69u1Lv379mDZtGsXFxYwZMwaAW2+9laSkJKZMMe+m+OSTTzJx4kTmzJlDSkoKWVlZAERHRxMdHV2HTRERkaCw4VNY+RZgg2tegtAoqysSH6txGBk5ciS5ublMnDiRrKwsevXqxYIFC0hISABg165d2I8ZFe/FF1+kvLyca6+9tsp2Jk2axMMPP3xm1YuISHApyoX//cGcvuAP0PI8a+uRs6JWHVjHjRvHuHHjTvi19PT0Kq937NhRm7cQEZH6xjDgk/uhJA+adoGBOj1TX2hgfxER8Q+r5sGGT8DuNE/PhBx/IYMEJ4URERGxXv4emP8Xc3rAg9Ds1FddSnBRGBEREWt5vfDhPeAqgBbnwgX3W12RnGUKIyIiYq0fX4Xt30BIBIyYoXvP1EMKIyIiYp28LbBwojk99FGIb2dtPWIJhREREbGGpwI++C1UlEKbAdD3DqsrEosojIiIiDW+fwb2/gRhcXD1dLDrI6m+0p4XEZGzb9/PkP6EOT38KYhrYW09YimFERERObvcZfDB78BbAZ2vhB4jra5ILKYwIiIiZ9fCf0DOOohqAr+aBjab1RWJxRRGRETk7Nn4GSx72ZweMQOi4q2tR/yCwoiIiJwdBfvMwc0AzrsX2g+2th7xGwojIiLie16veRlv6QFI7AGDJ1ldkfgRhREREfG9xc+ao6w6I+HamboJnlShMCIiIr61JwO+esycvvxJiG9vbT3idxRGRETEd1yF8N4d5mW8Xa6Gc26xuiLxQwojIiLiO5/+GQ5uh7hkuPJZXcYrJ6QwIiIivrHqHVg1F2x2+PUrENHQ6orETymMiIhI3TuwHT4Zb05f/AC06m9tPeLXFEZERKRuedxmP5HyQkg+Dy7+i9UViZ9TGBERkbqVPgX2Zph34/3NK+AIsboi8XMKIyIiUne2pMGip83pq56FBi2trUcCgsKIiIjUjfy98P5YwIDeo6HrNVZXJAFCYURERM6cxw3/HQMl+yGxO1z+lNUVSQBRGBERkTNm//oR2L3U7Cdy/X/AGW51SRJA1KtIRETOSLNDP+HY/qL5YsR0aNTG2oIk4OjIiIiI1N7B7Zyz8xVzuv846HyltfVIQNKRERERqR13GSHv3Y7NW4q3RT/sgx+2uiIJUDoyIiIitbPgr9iyV+MKicFzzavgcFpdkQQohREREam5n+dBxiwMbGS0+h3ENre6IglgOk0jIiI1k7MePrkfAO+FfyK3uLu19UjA05ERERGpPlcRvDMa3CXQZgDei3TfGTlzCiMiIlI9hmEeEcnbCDHN4Nevgt1hdVUSBBRGRESkepbOgNXvgs0B174O0U2srkiChMKIiIic3rZv4POHzOmhj0Gr/tbWI0FFYURERE7t4E549zYwPNDzRjjvbqsrkiCjMCIiIidXXgJzb4LSA9D8HPjVM2CzWV2VBBmFEREROTHDgI/uhezVENUERr4Jzgirq5IgpDAiIiIn9v2zsPZ9sIeYd+KNa2F1RRKkFEZEROR4m7+ELx82py9/Elqdb2k54jvzV2fx0U47FR6vZTUojIiISFX7t8J7twMG9B4Nfe+wuiLxkQ1ZBTz4wRq+yrTz3+WZltWhMCIiIke5Cs0Oq2X50KIfDP+nOqwGqfwSN7+dnUGp20uHOC/X9rbu/kIKIyIiYvJ64YPfQe56iE6EkbMhJMzqqsQHvF6D++etYOf+Elo0COe29l5CHNZFAoURERExLfoXbPgEHKHmlTMxiVZXJD4y7ctNfL0xl7AQOy/c2Isop7X1KIyIiAis+x98/Zg5fcXTkHyutfWIz3y+NovnvtoCwBO/6U7X5rEWV6QwIiIiezPg/bvM6dTfQe9brK1HfGZLThF/eudnAMZckMI15/jH5doKIyIi9dmh3fD2jVBRCu2HwbDHra5IfKSwzM1vZ/9EkauC1NaN+NvwzlaXVElhRESkviorgDkjoSgbErrBta+B3WF1VeIDXq/Bn975ma25xSTGhvPCqN44Leyw+kv+U4mIiJw9ngr47+2QsxaiE2DUPAiLsboq8ZF/p2/hi3XZhDrszLilD01i/OsqKYUREZH66PMJsGUhhETAjXM11HsQ+3pjDv9auAmAR0d0pVdyA2sLOgGFERGR+mbpS7DsZcAGv3kFknpbXZH4yM79xdz39goMA0altmTkuS2tLumEahVGpk+fTkpKCuHh4aSmprJs2bKTLrt27Vp+85vfkJKSgs1mY9q0abWtVUREztSmz2HBg+b0kMnQ+Upr6xGfKSxzc9d/Migoq6B3ywZMurKL1SWdVI3DyLx58xg/fjyTJk1i+fLl9OzZk2HDhpGTk3PC5UtKSmjTpg1PPPEEiYkaQEdExDJZq81+IoYXet8K5//B6orER9weL/e8tZyN2YU0jQnjxZv7EBbiv52TaxxGnn76acaOHcuYMWPo0qULM2bMIDIykpkzZ55w+XPPPZd//vOf3HDDDYSF+VeHGRGReqMwy7xyprwIWl9sDmyme84EJcMwmPS/tSzanEeE08Fro88lITbc6rJOKaQmC5eXl5ORkcGECRMq59ntdgYPHsySJUvqrCiXy4XL5ap8XVBQAIDb7cbtdtfZ+xzZVl1u098EexvVvsAX7G30i/aVF+OYMxJ7wV6Mxu2puGYmeAHvmdfkF+3zsUBr42vf72DO0l3YbPDMdd3plBB5ytp92b7qbrNGYSQvLw+Px0NCQkKV+QkJCWzYsKEmmzqlKVOmMHny5OPmf/HFF0RGRtbZ+xyxcOHCOt+mvwn2Nqp9gS/Y22hV+2xGBalbp5FQuApXSAzfJtxFydeL6/x9gn3/QWC08ef9Nl7fZAdsjGjlwbX9J+Zvr966vmhfSUlJtZarURg5WyZMmMD48eMrXxcUFJCcnMzQoUOJja27MfTdbjcLFy5kyJAhOJ0W3yXIR4K9jWpf4Av2NlraPsOL4+Nx2AtXYTgjcdz0XwYk9anTtwj2/QeB08ZVe/L568wfMfByU79kJv2qE7ZqnIrzZfuOnNk4nRqFkfj4eBwOB9nZ2VXmZ2dn12nn1LCwsBP2L3E6nT75QfDVdv1JsLdR7Qt8wd5GS9r3xT9g9Ttgc2C7/j+EpJzns7cK9v0H/t3GvYdK+d2clZS5vQzo2ITJV3cjpIYjrPqifdXdXo0qDQ0NpU+fPqSlpVXO83q9pKWl0b9//5pVKCIivrNkOix+zpy++gVoP8TaesRnCsvc3P76j+QWuuiUGMPzN55T4yBitRqfphk/fjyjR4+mb9++9OvXj2nTplFcXMyYMWMAuPXWW0lKSmLKlCmA2el13bp1ldN79+5l5cqVREdH065duzpsioiIALDqXfj8b+b04Ieh1yhLyxHfqfB4uXfOispLeGfedi4x4f559OZUahxGRo4cSW5uLhMnTiQrK4tevXqxYMGCyk6tu3btwm4/msgyMzM555xzKl9PnTqVqVOncskll5Cenn7mLRARkaO2fgUf3m1Op94NF9xvaTniO0cu4f12U27lJbzNG0RYXVat1KoD67hx4xg3btwJv/bLgJGSkoJhGLV5GxERqYnMFTDvFvOS3a6/hmGPayyRIPbqou28dfgS3mdv6EX3FnFWl1RrgXVSSURETmz/Vnjz2sODml0C18wAu/7EB6v3l+/h/+avB+Ch4Z0Z2jWwRzjXT6qISKArzIY3fw0leZDYA0a+CSEa8TpYLVyXzV/+uwqA2y9ozR0Xtra4ojOnMCIiEsjK8uGta+HgDmjQCm76L4TX3XhM4l9+2Lafe+csx+M1+E3vFvz9is7VGkvE3ymMiIgEKlcRvHUdZK2CyHi45QOISTj9ehKQ1uzN5843fqK8wsvgzgk8+Zvu2O2BH0RAYUREJDC5S2HujbB7KYTHmUGkcVurqxIf2ZpbxOiZyyhyVXBem0a8MCrwxhI5leBpiYhIfVFRDu/cCtu/hdBouPl9aNbD6qrERzIPlXLLq0vZX1xO96Q4Xrm1L+FOh9Vl1SmFERGRQOKpgPfugM1fQEgEjJoHLfpaXZX4yIHicm55bSmZ+WW0aRLFrDGBOajZ6SiMiIgECq8XProH1v8PHKFww1uQcqHVVYmPFLkquO31ZWzNLaZZXDiz70ilcXRwXiWlMCIiEggMAz79I6yaBzYHXDcL2g2yuirxkTK3h7Fv/MSqPfk0igpl9h2pJAXo6KrVoTAiIuLvDMO810zGLMAGv34ZOl1hdVXiI+UVXsbNWcGSbfuJDgvhjTH9aNc02uqyfEphRETE3331GPzwb3P66heg+7XW1iM+U17h5d45y/lyfTahIXZeubVvQA/zXl21ujeNiIicJYv+BYummtPDp8I5N1tbj/hMeYWXe96qGkT6t21sdVlnhcKIiIi/WvQvSHvEnB7yCPQba2094jOuCg/3vrWcL9fnEHY4iFzcoYnVZZ01CiMiIv7om6fg6/8zpwc+BBfcZ2094jOuCg/3vLmctA1mEHl1dF8ual9/gggojIiI+BfDgPQp8M2T5utBE+GiP1lbk/iMq8LD3W8u56vDQeS10edyYft4q8s66xRGRET8hWGYp2W+e9p8PeQRHREJYmVuD3e/mcHXG3PrdRABhREREf9gGLDwH7D4efP1sMeh/73W1iQ+U+b28Ls3M0jfmEu40wwiF7Srn0EEFEZERKx3ZByRI5fvXv5PSL3L2prEZ8rcHn47O4NvNplBZOboczm/HgcRUBgREbGWYcBnD8Cyl83Xv3oG+t5ubU3iMyXlFfx2dgaLNueZQeS2czm/bf0OIqAwIiJiHa8XPh0PGa8DNrjqOeh9q9VViY8cLC5nzKwfWbn7EBFOBzNvO7fejCNyOgojIiJW8Hrg4/tgxWzABiP+Db1GWV2V+Mi+/FJueW0ZW3KKiItwMvO2c+nTqqHVZfkNhRERkbOtwgXv3wXrPgSbHa55GXpcZ3VV4iNbcoq49bWlZOaXkRgbzn/u6EeHhBiry/IrCiMiImeTqwjm3QTb0sHuhN+8Cl1HWF2V+MjK3YcY8/oyDpa4aRMfxX/u6EeLhpFWl+V3FEZERM6W4v0w5zrYmwHOKLjhLWg70OqqxEcWbc7lt7MzKCn30KNFHK/fdi6No8OsLssvKYyIiJwNBXvh7esgbxNENIKb/gst+lhdlfjIJ6sy+eO8lbg9Bhe2i2fGLX2IDtNH7snYrS5ARCTYRZdlEvLGcDOIxCbB7QsURILY7CU7+P3bK3B7DK7o3ozXbuurIHIa+u6IiPiQLXMFF27+P2wVhdC4PdzyATRItros8QHDMJj25WaeTdsMwM3ntWTyVd1w2G0WV+b/FEZERHxlWzqOuaMIqSjG26wX9pvfgygNcBWMytwe/vreKj5amQnAfYPac//g9thsCiLVoTAiIuIL6z6C9+7E5iknN7oLDW76AHtUI6urEh/IKSzjrv9ksHL3IULsNh65uhujUltaXVZAURgREalrP8yAzyeA4cXb6Up+CBvBZWEaVyIYrcss4M43fiQzv4y4CCcv3tS73t9npjbUgVVEpK54PTD/AVjwVzC80GcMnmtexWt3Wl2Z+MDCddlcO2MxmflltImP4sN7L1AQqSUdGRERqQuuInjvTtj0mfl68GS44D6oqLC2LqlzhmHw8rfbeGLBBgwDLmjXmH+P6kNcpEJnbSmMiIicqYJ98PZI2PczOMLg1y9B12usrkp8oLzCy0MfrObdjD0A3JTakoev6orToRMNZ0JhRETkTGStgTnXm4OaRcbDjW9Dcj+rqxIfOFBczu/nrmLZjgPYbTDxV10YfX6KrpipAwojIiK1teVLeOc2KC+E+A4w6h1o1NrqqsQHdhfBP19ayp6DpcSEhfD8qHMY0LGp1WUFDYUREZHa+GkmfPpnMDyQchGMnA0RuiV8sDEMg7d/3M0zaxx4jFJaNorktdF9aa+77tYphRERkZrweuDLSbD4efN1zxvhyucgJNTauqTOlZRX8PcP1vD+ir2AjUGdmvD09eeoo6oPKIyIiFRX6SHzipktC83XA/4GlzwA6jMQdLbkFHHPWxlsyi7CYbdxRYsKpo7qRWiogogvKIyIiFRHzgaYOwoObIWQcLjqBehxndVViQ98/HMmD763iuJyD01iwph2fXfy1v2gjqo+pDAiInI66z+GD34H5UUQlwwj34TmvayuSupYeYWXx+evZ9biHQCc16YRz914Dg3DHcxfZ21twU5hRETkZLxeSJ8C3z5lvk65CK6bpZvdBaG9h0q5963lrNx9CIB7B7blj4M7EOKw43a7rS2uHlAYERE5kbJ8eP8u2LTAfJ16Nwx9FBzqMxBsPlmVyUMfrCG/1E1chJNnRvbk0k4JVpdVryiMiIj8Uu4ms3/I/s3miKpXPgu9brS6Kqlj+aVuJn20hg9XZgLQs0UcL4zqTXKjSIsrq38URkREjrXxM3hvrDmQWWyS2T8kqbfVVUkdW7wljz+9+zP78stw2G3cO7Adv7+0nYZ1t4jCiIgIgMcNXz0G308zX7c8H65/A6I1ymYwKXN7+OfnG3ntu+0ApDSO5OmRvejdUgPWWUlhREQkfw/893bYvdR83e8uGPa4+ocEmbWZ+fxx3ko2ZRcBMCq1JX+/ojORofootJr2gIjUbxsXwIe/g9KDEBYLVz0PXUdYXZXUIY/X4OVvt/H0wo24PQbx0WE8dW13dVL1IwojIlI/VZRD2mRY8oL5uvk5cO3rutFdkNmSU8iD763mp50HARjaJYEpv+5O4+gwiyuTYymMiEj9c3CneVpm70/m69S7YchkCNEHVLAoc3v4d/pWXkzfgttjEBXqYNJVXbmuTwuNpOqHFEZEpH5Z/wl8dI85jkh4HFz9b+j8K6urkjq0ZOt+HvpgNdvyigEY1Kkpj4zoRlKDCIsrk5NRGBGR+sFdBl8+DEtfNF8n9YVrZ0LDVpaWJXXnYHE5j89fz7sZewBoGhPGw1d15fJuiToa4udqdUH19OnTSUlJITw8nNTUVJYtW3bK5d999106depEeHg43bt3Z/78+bUqVkSkVvatglcGHg0i/cfBmM8URIKEYRh8sGIPg5/+pjKI3HxeS7780yUM795MQSQA1DiMzJs3j/HjxzNp0iSWL19Oz549GTZsGDk5OSdcfvHixdx4443ccccdrFixghEjRjBixAjWrFlzxsWLiJySpwK+nQqvXAo56yCqCYx6B4b9H4SEWl2d1IGd+4u5deYy/jjvZ/YXl9MhIZr37u7PYyO6ExuuS7MDRY1P0zz99NOMHTuWMWPGADBjxgw+/fRTZs6cyYMPPnjc8s8++yyXXXYZf/nLXwB49NFHWbhwIS+88AIzZsw4w/LPTNbuLRgHtrLrQAlhoU5C7HbsdnDYbEen7ea002FTuhYJJPu3mnfa3XP4yG3nK+FX03STuyBRWObm3+lbee277ZRXeAkNsfOHS9tx18VtCQ3RKKqBpkZhpLy8nIyMDCZMmFA5z263M3jwYJYsWXLCdZYsWcL48eOrzBs2bBgffvjhSd/H5XLhcrkqXxcUFADgdrvr7u6JhkHe23dzZflKXnk+g2cqrsXFyf9TcththDvtRDgdhDsdRFSZdhAR6iA2PISY8BBiwp3EhoeYryOcxISZ07ERThpFOglzOuqmDdVw5PsVrHedVPsCX5230TCwr3gD+5cTsblLMMJi8Ax9AqP79WCzwVn+Xgb7Pjzb7fN4Df67fC/PfLmF/cXlAJzfthGTr+xMSuMoMDy43Z46fU/twzPf9unUKIzk5eXh8XhISKg6UExCQgIbNmw44TpZWVknXD4rK+uk7zNlyhQmT5583PwvvviCyMi6uYGR3VtOhCcSh83gdyGfMNSewQMVd/GTt+MJl/d4DYpdHopdZ/5DHuYwiA6BaCfEOA2iDk9HOw1inBAbCnFOg7hQCHeYfz/P1MKFC898I35M7Qt8ddHGMPchztn1KgkFqwDIje7MilZjKd0TA3s+O+Ptn4lg34dno32b8m18sMNOZon5R7FJuMHVrbx0a5jDuqU5rPPx+2sf1lxJSUm1lvPLq2kmTJhQ5WhKQUEBycnJDB06lNjY2Dp7H/ewK/jhv0/RL3sObYr38W7oI3j7/ZaKiyfgCYnA4zXweA3cHoOyCg9lbg+l5V7K3B7KKjyUur2UlXsodXsoKfdQWFZBQZn78HMFhWUVFJa5KTg8v6C0ggqvgctjw+WB/S6AUyeNyFAHTWPCKh8JsWEkxIbTLC6c5nHhJDUIp1FU6ElPIbndbhYuXMiQIUNwOoPv/KnaF/jqqo229R/h+GwSttKDGI4wvJf+gwbn3sVAm7WH7IN9H56N9u3YX8wTCzaRtiEXgNjwEMYNbMtN/ZLPyikZ7cPaO3Jm43RqFEbi4+NxOBxkZ2dXmZ+dnU1iYuIJ10lMTKzR8gBhYWGEhR0/+JDT6azzb1R23Dl4RtyL/atJ2Fa+hWPZDBxbvoCrXoCUC+r0vQzDoKC0gv3FLg4Ul7O/uNx8LnJVTucVucgpcJFdUEZBWQUl5R527C9hx/6Tp8uwEDtJDSJo3iDi6HPDCFo2iqR5rBOv4ZvvnT9R+wJfrdtYsA/m/xk2fGK+btYT2zUv42jaibN3QvT0gn0f+qJ9+SVunvtqM/9ZsgO3x8Bht3FzakvuH9yBhlFnvwOy9mHttlkdNQojoaGh9OnTh7S0NEaMGAGA1+slLS2NcePGnXCd/v37k5aWxv333185b+HChfTv378mb+1bEQ1gxL+h6zXwvz/AgW0wa7h5s6xBkyAsuk7exmazERfpJC7SSZsmp1++pLyiMphkF7rIKSgju6CMffllZB4qZe+hUnIKXbgqvGzLK64c4OeXnDYHL2z9npTGUSQ3iqRlo0haNTafWzaOJCzEn/5ki1ST1wsZr5tjh7gKwB4CF46Hi/+iK2UCXJGrglnfb+eVRdvJLzX7HAzs2ISHruhMu6YxFlcnvlDj0zTjx49n9OjR9O3bl379+jFt2jSKi4srr6659dZbSUpKYsqUKQDcd999XHLJJfzrX//iiiuuYO7cufz000+8/PLLdduSutB+CNz7A3zxd1j+H1j2Mmz63LxxVptLzno5kaEhpMSHkBIfddJlyiu8ZOWXsfdwOMk8VMreg6XsOVTCrgMlZB4qw+2FrbnFbM09PqzYbJDUIILW8VGVj5T4KNrER5HUIIIQh3qlix/K2QAf3we7fzBfJ/WFq56DhK7W1iVnpNhVwX+W7OTlb7dysMQMIR0Sonnoii5c0qEa/8FJwKpxGBk5ciS5ublMnDiRrKwsevXqxYIFCyo7qe7atQu7/egH2Pnnn8+cOXP4+9//zt/+9jfat2/Phx9+SLdu3equFXUpPM4MH11GmH/sDu2E/1wFPUfBkEcg2r9+IUJD7LRsbB7hOJGSMhdvf7SANj1SySxwsetACbv2m0Fl5/4SilwV7DlYyp6DpSzanFdlXafDRnKjSNo2iaZd0+hjnqOI0fX7YoUKFyx6Ghb9C7xuCI2GQRPh3DvBriN8gaq03MPsH3bw0jfbKq+QaRMfxX2D2/OrHs1x2DWsQrCrVQfWcePGnfS0THp6+nHzrrvuOq677rravJV12g2Ce5bAwonw00z4eQ5s/BQu/Qf0vT1g/vA5HXbiw+HCdo2PO3dnGAZ5ReXs2F/M9lzzNM+OvGK25xWzfX8x5RVetuUWsy23mIXrqvb7SYgNqxJS2idE0yEhhsan6EwrckZ2LoGP/wB5m8zXHS6DK/4FcS2srUtqrczt4a2lu3gxfSt5ReZwDq0aR3LfoPZc1bO5jszWI355NY3fCIuBXz1jHhX5dDxkrTI7yq14E371NCT1sbrCM2Kz2WgSE0aTmDDOTWlU5Wter8G+gjK25RaxNaeIrbnFbMkpYmtuETmFLrILzMfirfurrNcw0kn7hBg6HA4n7Zqaz/G6XbfUVskBSHvE7B8CENUUhj9lHr1U8A1Ixa4K3vlpNy+mbyWn0AwhLRpG8IdB7fn1OUkKIfWQwkh1JJ8Ld6XDj6/BV4/CvpXwyiDoO8Y8UhLZ6HRbCDh2u42kw1fnXNS+6qmp/FK3GVIOB5QtOUVszilk14ESDpa4Wbb9AMu2H6iyTqOoUDokRNMpMZaOiTF0SIihY2IM0WH6EZST8FSYAeTr/4PSg+a83qNhyGSIaGhtbVIrOQVlzFq8g7eW7qrsmJrUIIJxl7bjN71baOTUekyfBNVld0DqXdDlalj4D1g1zzx9s+4jsy9Jz1Fgrx+/SHERTs5p2ZBzWlb9QCgt97A11wwmm7KL2JxtPu8+WMKB4nJ+2HaAH7ZVDSktGkbQ8XAw6ZgYQ+dmsbSOj8Kp/4zqt+2L4LO/Qs5a83XTrjD8n3V+ub2cHRuyCnjl2+387+e9uD0GACmNI7nzojZc17eFrugThZEai0mAX78M59xinrLJ3QAf3QvLZ8Nljwf8qZszERHqoFtSHN2S4qrMLy33sCWniA1ZBWzKLmRDViGbsgvJLnBVdp5N23D0RouhDjttm0bTOTGGTs1i6JQYS6dmMTSJDlN/lGB3aJd5Ndu6j8zXEQ1h4EPQZww49OcqkBiGwaLNebyyaFuVzvF9WzVk7MVtGNw5QR1TpZJ+u2ur9UXwu+/gh39D+pPmJYavXArdrzNP3ejW5JUiQh10bxFH9xZVQ8rB4nI2ZhdWBpSNhx9FrgrW7ytg/b4CWHF0+cZRoZXhpHOzWDo3i6FVw/Cz3BrxBYfXhf2bKfDDdKgoA5sd+t4BA/8WlKdBg1lpuYePV2Uy87vtbMgqBMBug8u6JXLnRW3o3VKn2OR4CiNnwuGEC+6Dbr+BtEdh1VxY/S6s+x+c9ztzAKaIBlZX6bcaRoVyXpvGnNemceU8wzDYc7CUDVmFbNhXwIasQtZnFbAjr5j9xeV8v2U/32852mk2xG6jSZiDtOLVdE2Ko9PhkNI0RiElIBgGtrXvM2jdgzjch0/hpVwElz+pMUMCzPp9hby7PJMPV+6lsKwCMG9ncX3fZG6/oPVJhx8QAYWRuhHXAn79Epx3t3mIecci+P5Z89TNJX81LwXWiJDVYrOZY5skN4pkSJejN1gsLfewOafw8BGTwsojJwVlFewrtfG/Vfv436p9lcvHR4fSuVksXZodOYoSS5sm6oviV7Z+DWmTCclcQQhgxCVjG/Z/0PkqXSUTIIpdFXy4fA8vrXKw65g7tyc3imBUv1aM6teSuEiNSSSnpzBSl5r3gtEfw+Yv4It/QN5GWPBXWPYSDJ4Mna/UH9laigh10KNFA3q0aFA5zzAMduUV8uan6UQndWRTTjHr9xWwfX8xeUXlLNqcV+VcdWiInQ4J0XROPBpQujSL1R/Ls21PBqQ9DNu/BcAIjWJD42G0u+VZnJF1dyNM8Z3Ve/KZs2wX/1u5l+JyD2DD6bAxtEsiN/RL5oK28djVH0RqQGGkrtls0GEYtB0EK/4DXz9u3uvmnVsg+Ty49O9mfxM5YzabjeYNIujW0GD4gDaVg7qVlFewMauwyhGUDYf7oqzZW8CavVXvItk8LrwynBzpi5LSOEp/TOtazgbz0vgjN7RzhELfO6jo/wc2ffMj7ZwR1tYnp7TnYAmfrNrHRyszzf5ch6U0jqRHVCF/u/FSEhvWzX28pP5RGPEVR4h5eqb7dfD9c7D4ebOT6xu/Ms+JD5igyxR9JDI05LhLj71eg90HS1i/r4B1x4SUPQdLycwvIzO/rMoVPRFOR+Wlxl2axdCpmTk+SqyGwa+5Q7sg/Qn4+W0wvGbn1J43woAHoUFLcLutrlBOIqewjPmr9vG/nzNZvutQ5fzQEDuXd0vkhnNb0ic5hs8++4zGGthQzoDCiK+FxcClD5kDpC36l3kDvh2LYNYiM5QM/Bu0Ot/qKoOe3W6jVeMoWjWO4rJuzSrnF5S52XBMOFm/r4CN2YWUuj2s3H2IlbsPVdlOi4YRh6/mMYNKp8QYWjWO0iWKJ1Kwz+w79dNr4DHvN0KnX5lXmzXtZG1tclL5JW4+W7OPj1dlsmTrfrzmsCDYbHBe68Zc2bM5w7sn0iDS7AfnVpiUOqAwcrbENjfvo3HhH80bfR0JJa9fDq0vhgF/g1b9ra6y3okNd9KvdSP6tT56+ajHa7A9r/jw6R2zw+yGfQVk5pdVjovy5fqj9+oJd9rNEWUTjg7c1jGxHg+Bf2C7GUJWvnU0hKRcBIMfhhZ9LS1NTizzUClfbcghbX02323JqxyYDKBXcgOu6tmcK3o0IyFWV6mJbyiMnG1xLcz72lz4R/juafOKm+3fmo/Wl5inbxRKLOWw22jX1LwJ4JU9m1fOzy9xsz6rgA2Hr+jZkGUeRSlze1m1J59Ve/KrbCc+OpSOiea4KB0TYuiQGEP7ptFEBesQ+Dnr4btnYPV/wfCY85LPgwF/hTYD1Xnbj3i8Bj/vOcRX63NI25BTpQ8IQKfEGK7s2ZyrejYnuZEuyRXfC9K/igGgQbJ5E74jR0pWvAnbvzEfyefB+eOg4/CAuTtwfRAX6TxuXBSP12Dn/mKzw2xWIRuzCtiYVcjOAyXkFZWT94txUeDoEPgdEo/eULBtk2jCnQG6r/dmmD/DRzqmgtmB+6I/qV+UHyksc7Nocx5p63NI35jD/uLyyq/ZbNC7ZUMu7dSUoV0SaJ8QY2GlUh8pjFitQUu4chpcNB6+nQor55gdXef9AI3awHn3QK+bIFT/nfgjh91GmybRtGkSzeXdj/ZFKSmvYFN2ERsPn+Y5cr+e3MITD4Fvt0GrxlGVR2TaH35u2cBPT/UYBuz4zuwHte3rwzNt5uXrF42H5udYWp5AkauCn3YcYMm2/fyw7QBr9ubj8R49/RITFsLFHZswqFNTBnRsSqMojYUk1lEY8RcNWsJVz5kdWpe9bN4h+MA28/43X/+fOTR2v7vMe+OI34sMDaFXcgN6JTeoMv9AcTmbsgvZnF1oDoWfVcTG7ELyS91szytme14xC9dlV1mnUZiD9/Iy6JAQezj4RNGmSZQ19+pxl5qnYZa9BFmrzXk2B/S43jzK16Tj2a1HKhW7Kvhp50F+2LafJVv3s/oX4QOgdXwUgzo15dLOTTk3pZEGARS/oTDib2ISYdBE8xD3irfMe3Uc3AGLpsLi58w/+v3HQdPOVlcqtdDoJEPg5xa62JxTxJYc867Hm7OL2JpbRF5ROQdcNr7dvJ9vN1c93RMTHkKb+CgzoBx5bhJFSuMoIkLr+JTPod3mVTEZb0Dp4WHbQyKg1yjzlgi6F9NZ5fUabMsrZtWeQ6zak8/K3YdYszefil+Ej+RGEZzX2vx5O69tY5IaaCwX8U8KI/4qNApS74Jz7zDPxS9+Hvb8aPYtWfGmeXVCn9vMw+IhfnooX6rFZrPRNDacprHhXNAuvsrXsg8V89bHX9KkbQ+27S9hW24x2/KK2HOwlMKyCn7ek8/Pv+g4C5AQG0arxlGkNI48/BxFq8aRtGocSUx1x0o5cipm2Uuw4VNzjBCAuJbQ707zztW6iZ3PGYbBvvwyVu05xMrd+azac4jVe/IpdFUct2xSgwj6tzXDR2rrRup8KgFDYcTf2R3Q5WrzsWspLHne/GDYsch8RDQy/zvtcxvEt7e6WqljjaJCaRsLw89tUTnCLECZ28PO/SVsyy1iW14xW3OLzKCSW0RBWQXZBS6yC1ws237guG3GR4fSoqF5/58WDSNo0TCC5IbmdPMGEYR7S2DNe7D0ZchZe3TF1hdDv99Cx8vVsdpHDh4+jbcpp4gt2WY/o03ZhVU6mx4R7rTTtXkcPVs0oGdyHL1bNlT4kIClMBJIWqaaj/w95tGR5f+Bgr2w5AXzoaMl9Ub44RFiOyYef9XDoZJyduwvYef+YnbkHX7eX8yuI1f4HH4cO6CbDS/n2ddzreMbhjt+JAIXAOX2cLY2u5K8LrcSndydZnERNMGOokjtuSo87D1Yyu6Dpew6UMLWHDNwbMouIq/IdcJ1HHYbHRNi6Jlsho8eLRrQISGaEPX5kCChMBKI4lqYQ2lf9GfY8iVkvG7enO+XR0t6jITGGumyvmkQGUqvyNDjOs+CeXnnzv0l7DlYwp6DpRRnbaH1nv/Rr+BzEo2jV/ds9TZjjudS3vVcQsHWaNhaCCwGzA/GpjFhJMaF0zQmjCYxYTSJDic+JpQm0ebr+MPPAXu5ci0ZhkGhq4LcQhc5BS72HiplZ14hP2yxM/vVZew5WEZ2YRmGcfJtJDWIqLzku32Cefl3+6Yxdd8PSMSPKIwEMkcIdLzMfJzkaElIfAc6OLvBgU6QoCsd6ruYcCfd4h10y/kKtsyBnd9Vfs0Ii8XVcQR7Uq5hU0gnmuWXcX1+GfsKysjKP/woKMPjNfsw7MsvO/37hYXQMCqUhpFO4iJDaRDhrDod5SQ61M7WAli/r5AGUeFEhTmICgshLMR+9q8WOoZhGJS6PRSWVVBQ6qagrIKCMnfl6/1F5eQWlZnBo9BF7uGHq8J7gq3ZgUOVryKcDpIbmafHWsdHHR5zJoZ2TaOJDtZB8UROQT/1weLI0ZKL/2IeLVnxJmz6HFveJjqzCV58H5r3hu7XQtdfQ2yz029TgoeryDx6tv5/sOlzcJcc/oIN2gyAc27G1ukKwp0RtAPanWQzHq9BXpGLffllZOWXVn4A5xa5yC0sJ7fIRd7h1+UVXgpdFRS6Kth1fNeVXwjhubVLqsxx2G1EhTqIDgsxw4nTjtNhJ9RhJzTEfHYenjafbZXhpeqRh6MvDAPcHgNXhQdXhZfyCu8vps3XxS4PBaXu465Oqa6YsBDiY8Jo3iCcpLhwinN2Mei8XrRuEkNyo0gaR4VaGrRE/I3CSLCxO6DDMPNRlk/F2v+xP/0lmhatw5a5HDKXw+cPQcqFZjDpOByim1pdtfhCWb4ZPNZ9ZAbUimOOZDRqYw6m1/MGM8hWk8NuIyE23LxHyQlOAx1x7OmKQyXlHCx2c6jUzaGScg6VuDlUWs7BEjf5JW4OFLvIOVAAzjCKXR5K3eZQ8h6vcfhoxPFXjZxNdhvERjiJCQ8hNtxJbLg53Tj66GmpJsecrmoSE1bllIrb7Wb+/J0M79GsSidkETlKYSSYhcdh9LiBH/bEMvzivjg3fQpr/gu7lx7tX/Lx/ZDUxzzV0+FySOiqe4gEstKDsPEzM4Bs/erojeoAGraGLldB56shqbdP97PNZqv84D4d88N6PsOHD8DpdOLxGpSUV1Ds8lDkqqDYVUFxeQWuCi/uCi/lHi9ujxd3hYHLc3Re+S9OjxzbumObeuRoSliIg7CQI9NH54WG2IkJD6kMH5GhDh3FEPExhZH6IrqpOW5J6l1wcKd56ea6D2Hfz7D3J/Px1WPmGBIdL4MOl5lHT3RVjn/zemHfCtjyFWxNg93Ljt6kDiC+g3lZeOerILF7QARNh91GTLiz+uOhiEjAUxipjxq2Mu8fctF4KMiETQtg4wLzJn35u8zh6Je9DKHR0HageTfh1hebH2wB8GEW9Ar2mUc9tqbB1q+Pjoh6REI3M3x0uRqa6moqEfF/CiP1XWxz6Hu7+SgvMQPJxs/MvgZFWbD+Y/MBENXUPFqScqE5pkl8e4WTsyDUXYBt43zYu8wMH8cORAYQFmuGxXaDzLvlamh2EQkwCiNyVGikObpmx8sPH/5fCVvSYMe35uH/4hxY+775AIhOMINJqwugRV9o2gUcOrR+RrxeyNtk3rl59zJCdi3h8gPbYM2xC9nMu+IeCR8t+ur7LiIBTWFETsxuNzs5JvWGS/4CFS7YmwHbD3d83b0MirLNvidr3jPXcYRBYjfzg/LII76jOR6KHM8wzDFhstdC1irze7p7GZQdqlzkyHEnI74jtpap5imzNgMhqvEJNykiEoj0KSHVExIGrc43H/wV3GVmONmxCHYuNo+ilOWb8/ZmHLNeBDTrAc16mbeXj29v9j2JTqhfp3jKiyFnA2SvMcNH9lpz+pjgUSkkwjzakdyPiuZ9+WLdIYZcdZ0uCxWRoKUwIrXjDIeUC8wHmP/lH9wOmSsOP1aaj/JC81Li3Uurrh8WawaTxu2PBpTG7cwxL8Jjz3Zr6kZZARzcAYd2ms9HHge2w4FtHDv4ViWbw2x7QtfDASTVvOrl8GkXw+3GvWn+2WuDiIgFFEakbths5kBajdpAt9+Y87xeOLDVDCf7foa8zbB/s/kB7So4/ijKEaHRENPM7Fwb2/wX04kQ3gDC48xAg4+PrhiGecSnOA+Kc6Hk8HNxnvkoyoZDu8w2/fKqll+KamqGjoSu5hUvCV3No0W6fFpE6jmFEfEdu/3wUY/20OP6o/MrXOaRgrxNZkDJ22xOH9hqfvCXF5mhZf/m075FSEg4wwgjZOdkM6CEx0JYDNid5mi0NsfhZ9sx0w6w2c0RSSvKzKHR3WXgLoWK0sPTJeaj5AB43dVvc2RjaJhS9dGgFTTtrJFuRUROQmFEzr6QMPPDuWnn479WXmyOo1GYaT4X7IXCfeZ4KIX7oDDraGABbBVlhFMGB/J9W3NYrBk0oppAVLz5iIw3XzdIPho6AvUUk4iIhRRGxL+ERkF8O/NxKl4PuApwFx3guy/nc1G/HoRUlJj9NlwF5tcNz+Fn7+Fpb9V5IeHgjDD7v4REHJ6OODw/0pyObGw+nOFnp/0iIvWQwogEJrsDIhpCSDQFkS0xWp4PutpERCQg2a0uQEREROo3hRERERGxlMKIiIiIWEphRERERCylMCIiIiKWUhgRERERSymMiIiIiKUURkRERMRSCiMiIiJiKYURERERsZTCiIiIiFhKYUREREQspTAiIiIilgqIu/YahgFAQUFBnW7X7XZTUlJCQUEBziC942uwt1HtC3zB3ka1L/AFext92b4jn9tHPsdPJiDCSGFhIQDJyckWVyIiIiI1VVhYSFxc3Em/bjNOF1f8gNfrJTMzk5iYGGw2W51tt6CggOTkZHbv3k1sbGydbdefBHsb1b7AF+xtVPsCX7C30ZftMwyDwsJCmjdvjt1+8p4hAXFkxG6306JFC59tPzY2Nih/wI4V7G1U+wJfsLdR7Qt8wd5GX7XvVEdEjlAHVhEREbGUwoiIiIhYql6HkbCwMCZNmkRYWJjVpfhMsLdR7Qt8wd5GtS/wBXsb/aF9AdGBVURERIJXvT4yIiIiItZTGBERERFLKYyIiIiIpRRGRERExFIKIyIiImKpoAsj06dPJyUlhfDwcFJTU1m2bNkpl3/33Xfp1KkT4eHhdO/enfnz51f5umEYTJw4kWbNmhEREcHgwYPZvHmzL5twSjVp3yuvvMJFF11Ew4YNadiwIYMHDz5u+dtuuw2bzVblcdlll/m6GadUkzbOmjXruPrDw8OrLBPI+3DAgAHHtc9ms3HFFVdULuNP+/Dbb7/lyiuvpHnz5thsNj788MPTrpOenk7v3r0JCwujXbt2zJo167hlavp77Ss1bd/777/PkCFDaNKkCbGxsfTv35/PP/+8yjIPP/zwcfuvU6dOPmzFqdW0jenp6Sf8Gc3KyqqyXKDuwxP9ftlsNrp27Vq5jD/twylTpnDuuecSExND06ZNGTFiBBs3bjztelZ/FgZVGJk3bx7jx49n0qRJLF++nJ49ezJs2DBycnJOuPzixYu58cYbueOOO1ixYgUjRoxgxIgRrFmzpnKZp556iueee44ZM2awdOlSoqKiGDZsGGVlZWerWZVq2r709HRuvPFGvv76a5YsWUJycjJDhw5l7969VZa77LLL2LdvX+Xj7bffPhvNOaGathHMIYyPrX/nzp1Vvh7I+/D999+v0rY1a9bgcDi47rrrqiznL/uwuLiYnj17Mn369Gotv337dq644goGDhzIypUruf/++7nzzjurfGDX5mfCV2ravm+//ZYhQ4Ywf/58MjIyGDhwIFdeeSUrVqyoslzXrl2r7L/vvvvOF+VXS03beMTGjRurtKFp06aVXwvkffjss89Wadfu3btp1KjRcb+D/rIPv/nmG+69915++OEHFi5ciNvtZujQoRQXF590Hb/4LDSCSL9+/Yx777238rXH4zGaN29uTJky5YTLX3/99cYVV1xRZV5qaqrx29/+1jAMw/B6vUZiYqLxz3/+s/Lrhw4dMsLCwoy3337bBy04tZq275cqKiqMmJgY44033qicN3r0aOPqq6+u61JrraZtfP311424uLiTbi/Y9uEzzzxjxMTEGEVFRZXz/G0fHgEYH3zwwSmXeeCBB4yuXbtWmTdy5Ehj2LBhla/P9HvmK9Vp34l06dLFmDx5cuXrSZMmGT179qy7wupQddr49ddfG4Bx8ODBky4TTPvwgw8+MGw2m7Fjx47Kef68D3NycgzA+Oabb066jD98FgbNkZHy8nIyMjIYPHhw5Ty73c7gwYNZsmTJCddZsmRJleUBhg0bVrn89u3bycrKqrJMXFwcqampJ92mr9Smfb9UUlKC2+2mUaNGVeanp6fTtGlTOnbsyN13383+/fvrtPbqqm0bi4qKaNWqFcnJyVx99dWsXbu28mvBtg9fe+01brjhBqKioqrM95d9WFOn+x2si++ZP/F6vRQWFh73O7h582aaN29OmzZtuOmmm9i1a5dFFdZer169aNasGUOGDOH777+vnB9s+/C1115j8ODBtGrVqsp8f92H+fn5AMf9zB3LHz4LgyaM5OXl4fF4SEhIqDI/ISHhuHOXR2RlZZ1y+SPPNdmmr9Smfb/017/+lebNm1f5gbrsssv4z3/+Q1paGk8++STffPMNl19+OR6Pp07rr47atLFjx47MnDmTjz76iDfffBOv18v555/Pnj17gODah8uWLWPNmjXceeedVeb70z6sqZP9DhYUFFBaWlonP/f+ZOrUqRQVFXH99ddXzktNTWXWrFksWLCAF198ke3bt3PRRRdRWFhoYaXV16xZM2bMmMF7773He++9R3JyMgMGDGD58uVA3fzt8heZmZl89tlnx/0O+us+9Hq93H///VxwwQV069btpMv5w2dhSJ1sRfzeE088wdy5c0lPT6/SwfOGG26onO7evTs9evSgbdu2pKenM2jQICtKrZH+/fvTv3//ytfnn38+nTt35qWXXuLRRx+1sLK699prr9G9e3f69etXZX6g78P6Ys6cOUyePJmPPvqoSn+Kyy+/vHK6R48epKam0qpVK9555x3uuOMOK0qtkY4dO9KxY8fK1+effz5bt27lmWeeYfbs2RZWVvfeeOMNGjRowIgRI6rM99d9eO+997JmzRpL+yBVV9AcGYmPj8fhcJCdnV1lfnZ2NomJiSdcJzEx8ZTLH3muyTZ9pTbtO2Lq1Kk88cQTfPHFF/To0eOUy7Zp04b4+Hi2bNlyxjXX1Jm08Qin08k555xTWX+w7MPi4mLmzp1brT9sVu7DmjrZ72BsbCwRERF18jPhD+bOncudd97JO++8c9zh8F9q0KABHTp0CIj9dzL9+vWrrD9Y9qFhGMycOZNbbrmF0NDQUy7rD/tw3LhxfPLJJ3z99de0aNHilMv6w2dh0ISR0NBQ+vTpQ1paWuU8r9dLWlpalf+cj9W/f/8qywMsXLiwcvnWrVuTmJhYZZmCggKWLl160m36Sm3aB2YP6EcffZQFCxbQt2/f077Pnj172L9/P82aNauTumuitm08lsfjYfXq1ZX1B8M+BPOyO5fLxc0333za97FyH9bU6X4H6+Jnwmpvv/02Y8aM4e23365ySfbJFBUVsXXr1oDYfyezcuXKyvqDYR+CeZXKli1bqvUPgZX70DAMxo0bxwcffMBXX31F69atT7uOX3wW1kk3WD8xd+5cIywszJg1a5axbt0646677jIaNGhgZGVlGYZhGLfccovx4IMPVi7//fffGyEhIcbUqVON9evXG5MmTTKcTqexevXqymWeeOIJo0GDBsZHH31krFq1yrj66quN1q1bG6WlpX7fvieeeMIIDQ01/vvf/xr79u2rfBQWFhqGYRiFhYXGn//8Z2PJkiXG9u3bjS+//NLo3bu30b59e6OsrOyst682bZw8ebLx+eefG1u3bjUyMjKMG264wQgPDzfWrl1buUwg78MjLrzwQmPkyJHHzfe3fVhYWGisWLHCWLFihQEYTz/9tLFixQpj586dhmEYxoMPPmjccsstlctv27bNiIyMNP7yl78Y69evN6ZPn244HA5jwYIFlcuc7nvmz+176623jJCQEGP69OlVfgcPHTpUucyf/vQnIz093di+fbvx/fffG4MHDzbi4+ONnJycs94+w6h5G5955hnjww8/NDZv3mysXr3auO+++wy73W58+eWXlcsE8j484uabbzZSU1NPuE1/2od33323ERcXZ6Snp1f5mSspKalcxh8/C4MqjBiGYTz//PNGy5YtjdDQUKNfv37GDz/8UPm1Sy65xBg9enSV5d955x2jQ4cORmhoqNG1a1fj008/rfJ1r9dr/OMf/zASEhKMsLAwY9CgQcbGjRvPRlNOqCbta9WqlQEc95g0aZJhGIZRUlJiDB061GjSpInhdDqNVq1aGWPHjrXkD8SxatLG+++/v3LZhIQEY/jw4cby5curbC+Q96FhGMaGDRsMwPjiiy+O25a/7cMjl3n+8nGkTaNHjzYuueSS49bp1auXERoaarRp08Z4/fXXj9vuqb5nZ1NN23fJJZeccnnDMC9lbtasmREaGmokJSUZI0eONLZs2XJ2G3aMmrbxySefNNq2bWuEh4cbjRo1MgYMGGB89dVXx203UPehYZiXsUZERBgvv/zyCbfpT/vwRG0Dqvxe+eNnoe1w8SIiIiKWCJo+IyIiIhKYFEZERETEUgojIiIiYimFEREREbGUwoiIiIhYSmFERERELKUwIiIiIpZSGBERERFLKYyIiIiIpRRGRERExFIKIyIiImKp/weOLXcn0QkK0wAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAi8AAAGdCAYAAADaPpOnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy80BEi2AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA0m0lEQVR4nO3df3xU1Z3/8fckmUyIEpAvkABG8FdV/AEKmzyi64o1P6Qsj7KPtUXxB6VKt2pabGrRdJGQ8lX8QRFrWbNaKVqLoG3VdrWENBhoS4QK+F21yoqNYoUErasTGExuMuf7B83ImF8zYW7m3juv5+ORB8ydc0/OJ3cm88k5557jM8YYAQAAuERashsAAAAQD5IXAADgKiQvAADAVUheAACAq5C8AAAAVyF5AQAArkLyAgAAXIXkBQAAuEpGshuQaOFwWPv27dPQoUPl8/mS3RwAABADY4xaW1s1duxYpaX13bfiueRl3759ys/PT3YzAADAALz33ns68cQT+yzjueRl6NChko4En5OTk9C6LcvSxo0bVVpaKr/fn9C6nYD43M/rMXo9Psn7MRKf+9kVYzAYVH5+fuRzvC+eS166hopycnJsSV6ys7OVk5PjyRcl8bmf12P0enyS92MkPvezO8ZYpnwwYRcAALgKyQsAAHAVkhcAAOAqnpvzEgtjjDo6OtTZ2RnXeZZlKSMjQ59++mnc5w6m9PR0ZWRkcKs4AMCTUi55aW9v1/79+xUKheI+1xijvLw8vffee45PDLKzszVmzBhlZmYmuykAACRUSiUv4XBYTU1NSk9P19ixY5WZmRlXEhIOh3Xw4EEdf/zx/S6gkyzGGLW3t+uDDz5QU1OTTj/9dMe2FQCAgUip5KW9vV3hcFj5+fnKzs6O+/xwOKz29nZlZWU5OiEYMmSI/H6/3n333Uh7AQDwCud+AtvIyYlHoqRCjACA1MQnHAAAcBWSFwAA4CokLwAAwFVSasIuAADonTFGh62+1zGzrA61dR4pmywkLwAAQMYYXVHTqB3v/m8MpTNUVtapZC0llvLDRsYYhdo7Yv463N4ZV/m+vuLJWidMmKCVK1dGHZs8ebKWLFmS2B8IACAlHbY6Y0xcki/le14OW52auLg2Kd/7zz8oU3Zmyl8CAIDDvLyoWNmZ6T0+Z1mWams3aoi/5+cHA5+cAAAgSnZmeq9/XFs+o0C6krpNTsonL0P86frzD8piKhsOh9UabNXQnKEJWQQumVkrAABulfLJi8/ni3noJhwOq+Pv2ehgr2CblpbWbY6MZVmD2gYAAJwg5SfsusWoUaO0f//+yONgMKimpqYktggAgOQgeXGJL37xi/rZz36m3//+93r11Vc1d+5cpacz7AQASD0pP2zkFpWVlWpqatI///M/a9iwYVq6dCk9LwCAlETy4hI5OTlat25d1LG5c+cmqTUAACQPw0YAAMBVSF4AAICrkLwAAABXIXkBAACukpLJSzK38R4sqRAjACA1pdTdRn6/X5IUCoU0ZMiQJLfGXqFQSNJnMQMAUpsxRoetzl6fD7X3/pzTpFTykp6eruHDh+vAgQOSpOzs7Lg2lgqHw2pvb9enn3466NsDxMoYo1AopAMHDmj48OEsZAcAkDFGV9Q0ase7/5vspiRESiUvkpSXlydJkQQmHsYYHT58WEOGDEnqbpqxGD58eCRWAEBqO2x1xpy4TB1/guM3Dk655MXn82nMmDEaPXp03BsbWpalLVu26J/+6Z8cPRzj9/vpcQEA9OjlRcXKzuz9M2KIP93xf6CnXPLSJT09Pe4P+PT0dHV0dCgrK8vRyQsAAL3JzkxXdqa7P/6dOXEDAACgF7YmL1u2bNHMmTM1duxY+Xw+Pfvss32W/9WvfqWSkhKNGjVKOTk5KioqUm1trZ1NBAAALmNr8nLo0CFNmjRJq1atiqn8li1bVFJSohdeeEE7duzQpZdeqpkzZ2rXrl12NhMAALiIrYNe06dP1/Tp02Muv3LlyqjHd911l5577jn95je/0fnnn5/g1gEAADdy9JyXcDis1tZWjRgxItlNAQAADuHo6cbLly/XwYMH9dWvfrXXMm1tbWpra4s8DgaDko7c1hzvrdD96aov0fU6BfG5n9dj9Hp8kvdjJL7ksKyOo/5vyfINfAsZu2KMpz6fGaRNcHw+n5555hnNmjUrpvJr167V/Pnz9dxzz6m4uLjXckuWLFF1dXWP52dnZw+0uQAAeEZbp7Rw+5H+insLOhRw4FJgoVBIc+bM0SeffKKcnJw+yzoyeVm3bp2+/vWv6+mnn9aMGTP6LNtTz0t+fr4+/PDDfoOPl2VZqqurU0lJiSfXeSE+9/N6jF6PT/J+jMSXHKH2Dk1aukmS9P/u+OIxrfNiV4zBYFAjR46MKXlx3LDRk08+qa9//etat25dv4mLJAUCAQUCgW7H/X6/bS8cO+t2AuJzP6/H6PX4JO/HSHyDy28+WzH3SNuO/eM/0THGU5etycvBgwe1Z8+eyOOmpia98sorGjFihE466SRVVlbq/fff1+OPPy7pyFDP3Llz9cADD6iwsFDNzc2SpCFDhmjYsGF2NhUAALiErXcbvfzyyzr//PMjtzlXVFTo/PPP1+LFiyVJ+/fv1969eyPlH374YXV0dOjmm2/WmDFjIl8LFiyws5kAALiaMUah9o4+vjqT3cSEsrXnZdq0aeprSs2aNWuiHjc0NNjZHAAAPMcYoytqGmPeNdoLHL3OCwAA6NthqzPmxGXq+BM0xO/AW43i5LgJuwAAYGBeXlSs7Mzek5Mh/nT5fL5en3cLkhcAADwiOzP9mG6DdguGjQAAgKuQvAAAAFcheQEAAK5C8gIAAFyF5AUAALgKyQsAAHAVkhcAAOAqJC8AAMBVSF4AAICreH8ZPgAAXMwYo8NW77tCe23H6FiQvAAA4FCpuGN0LBg2AgDAoVJxx+hY0PMCAIALpMqO0bEgeQEAwAVSZcfoWDBsBAAAXIXkBQAAuArJCwAAcBWSFwAA4CokLwAAwFVIXgAAgKuQvAAAAFcheQEAAK7CajcAACQJmy4ODMkLAABJwKaLA8ewEQAAScCmiwNHzwsAAEnGpovxIXkBACDJ2HQxPgwbAQAAVyF5AQAArkLyAgAAXMXW5GXLli2aOXOmxo4dK5/Pp2effbbfcxoaGnTBBRcoEAjotNNO05o1a+xsIgAAcBlbk5dDhw5p0qRJWrVqVUzlm5qaNGPGDF166aV65ZVXdMstt+iGG25QbW2tnc0EAAAuYuvU5unTp2v69Okxl6+pqdHJJ5+sH/7wh5Kks846S3/4wx90//33q6yszK5mAgAAF3HUfVmNjY0qLi6OOlZWVqZbbrml13Pa2trU1tYWeRwMBiVJlmXJsqyEtq+rvkTX6xTE535ej9Hr8Unej5H4ji7bEXWe5TO2tSuR7LqG8dTnqOSlublZubm5Ucdyc3MVDAZ1+PBhDRkypNs5y5YtU3V1dbfjGzduVHZ2ti3trKurs6VepyA+9/N6jF6PT/J+jMQntXVKXR/DtbUbFXDZArqJvoahUCjmso5KXgaisrJSFRUVkcfBYFD5+fkqLS1VTk5OQr+XZVmqq6tTSUmJ/H5/Qut2AuJzP6/H6PX4JO/HmErxZWRk9Lnp4uH2Tmn7ZklSWVmpaxaps+sado2cxMJRP6m8vDy1tLREHWtpaVFOTk6PvS6SFAgEFAgEuh33+/22vTHsrNsJiM/9vB6j1+OTvB+j1+PLyMjQVY++HPPeRUd+Ho76SO5Xoq9hPHU5ap2XoqIi1dfXRx2rq6tTUVFRkloEAED82HTRXrameQcPHtSePXsij5uamvTKK69oxIgROumkk1RZWan3339fjz/+uCTpm9/8pn784x9r4cKF+vrXv65Nmzbpqaee0vPPP29nMwEAsA2bLiaercnLyy+/rEsvvTTyuGtuyty5c7VmzRrt379fe/fujTx/8skn6/nnn9d3vvMdPfDAAzrxxBP1k5/8hNukAQCuxaaLiWfrT3PatGkypvdbv3paPXfatGnatWuXja0CAABu5qg5LwAAAP0heQEAAK5C8gIAAFyF5AUAALgKyQsAAHAVkhcAAOAqJC8AAMBVWDUHAIA4GWN63HTRsjrU1vn3TRdhG5IXAADiYIzRFTWNfexdlBHZLRr2YNgIAIA4sOli8tHzAgDAAH1+00XLslRbu1FlZaXy+/1sumgTkhcAAAbo85suWj6jQLqUnZkhv5+PWLswbAQAAFyF5AUAALgKyQsAAHAVkhcAAOAqJC8AAMBVSF4AAICrkLwAAABX4SZ0AACO0tu+RV1C7FuUdCQvAAD8Xf/7FsEJGDYCAODv2LfIHeh5AQCgB5/ft+jz2LcoeUheAADowef3LYJzMGwEAABcheQFAAC4CskLAABwFZIXAADgKiQvAADAVUheAACAq5C8AAAAV+EGdgBAymDfIm+wPXlZtWqV7rvvPjU3N2vSpEl68MEHVVBQ0Gv5lStX6qGHHtLevXs1cuRIXXHFFVq2bJmysrLsbioAwMPYt8g7bB02Wr9+vSoqKlRVVaWdO3dq0qRJKisr04EDB3osv3btWt1+++2qqqrSG2+8oUcffVTr16/X97//fTubCQBIAexb5B229rysWLFC8+fP17x58yRJNTU1ev7557V69Wrdfvvt3cpv3bpVF110kebMmSNJmjBhgq666ipt27bNzmYCAFIM+xa5m23JS3t7u3bs2KHKysrIsbS0NBUXF6uxsbHHcy688EI98cQT2r59uwoKCvSXv/xFL7zwgq699tpev09bW5va2toij4PBoCTJsixZlpWgaBSp8+h/vYb43M/rMXo9Psn7MSYzPsvqiPzf7wvL7+t98KGjo6PX5/r+Ht6+fpJ9McZTn88YYxL63f9u3759GjdunLZu3aqioqLI8YULF2rz5s299qb86Ec/0q233ipjjDo6OvTNb35TDz30UK/fZ8mSJaquru52fO3atcrOzj72QAAAntDWKS3cfuRv9nsLOhRgVMhRQqGQ5syZo08++UQ5OTl9lnXU3UYNDQ2666679B//8R8qLCzUnj17tGDBAi1dulR33HFHj+dUVlaqoqIi8jgYDCo/P1+lpaX9Bh8vy7JUV1enkpIS+f3+hNbtBMTnfl6P0evxSd6PMZnxhdo7tHD7JklSWVmpLTtGe/36SfbF2DVyEgvbkpeRI0cqPT1dLS0tUcdbWlqUl5fX4zl33HGHrr32Wt1www2SpHPPPVeHDh3SN77xDf37v/+70tK6d/EFAgEFAoFux/1+v20vHDvrdgLicz+vx+j1+CTvx5iM+PzmszksR76/fX+/e/36SYmPMZ66bLvbKDMzU1OmTFF9fX3kWDgcVn19fdQw0tFCoVC3BCU9/Ui/nk2jWwAAwGVsHTaqqKjQ3LlzNXXqVBUUFGjlypU6dOhQ5O6j6667TuPGjdOyZcskSTNnztSKFSt0/vnnR4aN7rjjDs2cOTOSxAAAgNRma/Iye/ZsffDBB1q8eLGam5s1efJkbdiwQbm5uZKkvXv3RvW0LFq0SD6fT4sWLdL777+vUaNGaebMmbrzzjvtbCYAAHAR2yfslpeXq7y8vMfnGhoaohuTkaGqqipVVVXZ3SwAgIf0t+y/xNL/XuKou40AAIgXy/6nHnaVBgC4WjzL/kss/e8F9LwAADyjv2X/JZb+9wKSFwCAZ2Rnptuy+BychWEjAADgKiQvAADAVUheAACAq5C8AAAAVyF5AQAArsKUbACAo/W3ei4r56YekhcAgGOxei56wrARAMCx4lk9l5VzUwc9LwAAV+hv9VxWzk0dJC8AAFdg9Vx0YdgIAAC4CskLAABwFZIXAADgKiQvAADAVUheAACAq5C8AAAAV+GeMwBA0rD0PwaC5AUAkBQs/Y+BYtgIAJAULP2PgaLnBQCQdCz9j3iQvAAAko6l/xEPho0AAICrkLwAAABXIXkBAACuQvICAABcheQFAAC4ClO7AQC26Gn1XMvqUFunFGrvkGX4+xkDQ/ICAEi4vlfPzdDC7ZsGvU3wDtvT3lWrVmnChAnKyspSYWGhtm/f3mf5jz/+WDfffLPGjBmjQCCgL3zhC3rhhRfsbiYAIIFYPRd2srXnZf369aqoqFBNTY0KCwu1cuVKlZWVaffu3Ro9enS38u3t7SopKdHo0aP1i1/8QuPGjdO7776r4cOH29lMAICNjl4917Is1dZuVFlZqfx+vyRWz0X8bE1eVqxYofnz52vevHmSpJqaGj3//PNavXq1br/99m7lV69erY8++khbt26NvKgnTJhgZxMBADY7evVcy2cUSJeyMzPk9zNzAQNj2yunvb1dO3bsUGVlZeRYWlqaiouL1djY2OM5v/71r1VUVKSbb75Zzz33nEaNGqU5c+botttuU3p6z12KbW1tamtrizwOBoOSjmT3lmUlMCJF6kt0vU5BfO7n9Ri9Hp/knRgtq+Oo/1uyfCby/6P/9RqvxyfZF2M89fmMMSah3/3v9u3bp3Hjxmnr1q0qKiqKHF+4cKE2b96sbdu2dTvnzDPP1DvvvKOrr75aN910k/bs2aObbrpJ3/72t1VVVdXj91myZImqq6u7HV+7dq2ys7MTFxAAIGZtndLC7Uf+Pr63oEMBprSgH6FQSHPmzNEnn3yinJycPss6qs8uHA5r9OjRevjhh5Wenq4pU6bo/fff13333ddr8lJZWamKiorI42AwqPz8fJWWlvYbfLwsy1JdXZ1KSkoiw1peQnzu5/UYvR6f5J0YQ+0dkTuKyspKPxs28kh8vfF6fJJ9MXaNnMTCtuRl5MiRSk9PV0tLS9TxlpYW5eXl9XjOmDFj5Pf7o4aIzjrrLDU3N6u9vV2ZmZndzgkEAgoEAt2O+/1+2144dtbtBMTnfl6P0evxSe6P0W8+m4B7JJbojxu3x9cfr8cnJT7GeOqy7VbpzMxMTZkyRfX19ZFj4XBY9fX1UcNIR7vooou0Z88ehcPhyLH/+Z//0ZgxY3pMXAAAyWGMUai9o4+vzv4rAQbI1mGjiooKzZ07V1OnTlVBQYFWrlypQ4cORe4+uu666zRu3DgtW7ZMknTjjTfqxz/+sRYsWKBvfetbeuutt3TXXXfp29/+tp3NBADEoe8F6AD72Zq8zJ49Wx988IEWL16s5uZmTZ48WRs2bFBubq4kae/evUpL+6zzJz8/X7W1tfrOd76j8847T+PGjdOCBQt022232dlMAEAcWIAOyWb7hN3y8nKVl5f3+FxDQ0O3Y0VFRXrppZdsbhUAIBGOXoCuJyxABzs46m4jAIC7HL0AHTBY2NITAAC4CskLAABwFZIXAADgKiQvAADAVZhlBQCIYozRYav3ReZYgA7JRvICAIhgATq4AcNGAIAIFqCDG9DzAgDoEQvQwalIXgAAPWIBOjgVw0YAAMBVSF4AAICrkLwAAABXIXkBAACuwkwsAEghLEAHLyB5AYAUwQJ08AqGjQAgRbAAHbyCnhcASEEsQAc3I3kBgBTEAnRwM4aNAACAq5C8AAAAVyF5AQAArsKAJwB4BGu4IFWQvACAB7CGC1IJw0YA4AGs4YJUQs8LAHgMa7jA60heAMBjWMMFXsewEQAAcBWSFwAA4CokLwAAwFUYFAUAF2ANF+AzJC8A4HCs4QJEG5Rho1WrVmnChAnKyspSYWGhtm/fHtN569atk8/n06xZs+xtIAA4GGu4ANFs73lZv369KioqVFNTo8LCQq1cuVJlZWXavXu3Ro8e3et577zzjm699VZdfPHFdjcRAFyDNVyAQeh5WbFihebPn6958+Zp4sSJqqmpUXZ2tlavXt3rOZ2dnbr66qtVXV2tU045xe4mAoBrdK3h0tsXiQtSga09L+3t7dqxY4cqKysjx9LS0lRcXKzGxsZez/vBD36g0aNH6/rrr9fvf//7Pr9HW1ub2traIo+DwaAkybIsWZZ1jBFE66ov0fU6BfG5n9dj9Hp8Us8xWlZH1POWzwx6uxLF69fQ6/FJ9sUYT322Ji8ffvihOjs7lZubG3U8NzdXb775Zo/n/OEPf9Cjjz6qV155JabvsWzZMlVXV3c7vnHjRmVnZ8fd5ljU1dXZUq9TEJ/7eT1Gr8cnRcfY1il1/bqurd2ogAemtHj9Gno9PinxMYZCoZjLOupuo9bWVl177bV65JFHNHLkyJjOqaysVEVFReRxMBhUfn6+SktLlZOTk9D2WZaluro6lZSUyO/3J7RuJyA+9/N6jF6PT+o5xlB7hxZu3yRJKisrdfXS/16/hl6PT7Ivxq6Rk1jY+g4YOXKk0tPT1dLSEnW8paVFeXl53cq//fbbeueddzRz5szIsXA4fKShGRnavXu3Tj311KhzAoGAAoFAt7r8fr9tLxw763YC4nM/r8fotfiOXsPFMj61dR75V8b392OfTU88Ert7k5cuXruGn+f1+KTExxhPXba+AzIzMzVlyhTV19dHbncOh8Oqr69XeXl5t/JnnnmmXn311ahjixYtUmtrqx544AHl5+fb2VwAGHQ9r+GSEelpAdCd7el7RUWF5s6dq6lTp6qgoEArV67UoUOHNG/ePEnSddddp3HjxmnZsmXKysrSOeecE3X+8OHDJanbcQDwAtZwAeJne/Iye/ZsffDBB1q8eLGam5s1efJkbdiwITKJd+/evUpLY4slAHh5UbH8vrBqazeqrKy0Wzc6a7gARwzKwGl5eXmPw0SS1NDQ0Oe5a9asSXyDAMCBsjPT5felKZAuZWdmeGJuC2AHujwAAICrkLwAAABXIXkBAACuwoAqANjo6DVcehJq7/05AD0jeQEAm/S8hguAY8WwEQDYhDVcAHvQ8wIAg+DlRcXKzuw9OWENFyB2JC8AMAiyM9NdvaEi4CQMGwEAAFcheQEAAK5CHyYADBC3QQPJQfICAAPAbdBA8jBsBAADwG3QQPLQ8wIAx4jboIHBRfICAMeI26CBwcWwEQAAcBWSFwAA4Cr0cwJAD7gNGnAukhcA+BxugwacjWEjAPgcboMGnI2eFwDoA7dBA85D8gIAfeA2aMB5GDYCAACuwp8TAFIOdxIB7kbyAiClcCcR4H4MGwFIKdxJBLgfPS8AUhZ3EgHuRPICIGVxJxHgTgwbAQAAVyF5AQAArkJ/KQBP4TZowPtIXgB4BrdBA6lhUIaNVq1apQkTJigrK0uFhYXavn17r2UfeeQRXXzxxTrhhBN0wgknqLi4uM/yANCF26CB1GB7z8v69etVUVGhmpoaFRYWauXKlSorK9Pu3bs1evTobuUbGhp01VVX6cILL1RWVpbuuecelZaW6vXXX9e4cePsbi4Aj+A2aMC7bO95WbFihebPn6958+Zp4sSJqqmpUXZ2tlavXt1j+Z///Oe66aabNHnyZJ155pn6yU9+onA4rPr6erubCsBDum6D7u2LxAVwL1t7Xtrb27Vjxw5VVlZGjqWlpam4uFiNjY0x1REKhWRZlkaMGNHj821tbWpra4s8DgaDkiTLsmRZ1jG0vruu+hJdr1MQn/t5Pcb+4rOsjqiyls8MSrsSKdWvodt5PT7Jvhjjqc9njLHt3b1v3z6NGzdOW7duVVFRUeT4woULtXnzZm3btq3fOm666SbV1tbq9ddfV1ZWVrfnlyxZourq6m7H165dq+zs7GMLAICjGCO1h3t/vj0sLXr5yN9k9xZ0KMCUFsA1QqGQ5syZo08++UQ5OTl9lnX03UZ333231q1bp4aGhh4TF0mqrKxURUVF5HEwGFR+fr5KS0v7DT5elmWprq5OJSUl8vv9Ca3bCYjP/bwcozFGsx/epl1/DcZUvqys1JWr53r5GkrE5wV2xdg1chILW9/ZI0eOVHp6ulpaWqKOt7S0KC8vr89zly9frrvvvlu/+93vdN555/VaLhAIKBAIdDvu9/tte+HYWbcTEJ/7eTHGUHtHzInL1PEnKCc7y9XzWrx4DY9GfO6X6BjjqcvW5CUzM1NTpkxRfX29Zs2aJUmRybfl5eW9nnfvvffqzjvvVG1traZOnWpnEwG40Eu3XaKc43rujZW4kwjwOtv7VCsqKjR37lxNnTpVBQUFWrlypQ4dOqR58+ZJkq677jqNGzdOy5YtkyTdc889Wrx4sdauXasJEyaoublZknT88cfr+OOPt7u5AFxgCBsqAinN9nf/7Nmz9cEHH2jx4sVqbm7W5MmTtWHDBuXm5kqS9u7dq7S0z+7Yfuihh9Te3q4rrrgiqp6qqiotWbLE7uYCAACHG5Q/XcrLy3sdJmpoaIh6/M4779jfIACOw55EAGJFvyuApGNPIgDxGJS9jQCgL/HsSXTyUMOeRECKo+cFgKP0tSeRZVl6sW4jdxIBKY7kBYCjZPdxJ5HlMyJvAUDyAsB2TMYFkEgkLwBsxWRcAInGhF0AtopnMu7U8ScwGRdAv+h5ATBo+pqMK7GsP4DYkLwAGDR9TcYFgFgxbAQAAFyFP4EAHBPuJAIw2EheAAwYdxIBSAaGjQAMGHcSAUgGel4AJAR3EgEYLCQvABKCO4kADBZ+0wDoFZNxATgRyQuAHjEZF4BTMWEXQI+YjAvAqeh5AdAvJuMCcBKSFwD9YjIuACfhtxGQopiMC8CtSF6AFMRkXABuxoRdIAUxGReAm9HzAqQ4JuMCcBuSFyDFMRkXgNvwGwvwIGOMQu0dvT7PZFwAbkbyAniMMdKVP/mTdu79ONlNAQBbMGEX8Jj2sGJOXJiMC8CN6HkBPIzJuAC8iOQFcJm+FpezrA61hz97zGRcAF7EbzXARWJbXI63NQBvY84L4CIsLgcAg/Qn2qpVq3TfffepublZkyZN0oMPPqiCgoJeyz/99NO644479M477+j000/XPffcoy996UuD0VTANXqaz2JZlmprN6qsrFQ52VnMZwHgSbb3vKxfv14VFRWqqqrSzp07NWnSJJWVlenAgQM9lt+6dauuuuoqXX/99dq1a5dmzZqlWbNm6bXXXrO7qYCrdM1n+fxXIF3KzswgcQHgWbYnLytWrND8+fM1b948TZw4UTU1NcrOztbq1at7LP/AAw/o8ssv1/e+9z2dddZZWrp0qS644AL9+Mc/trupQNJ1LS7X+xeLywGArcNG7e3t2rFjhyorKyPH0tLSVFxcrMbGxh7PaWxsVEVFRdSxsrIyPfvss3Y2FUg6dnoGgNjYmrx8+OGH6uzsVG5ubtTx3Nxcvfnmmz2e09zc3GP55ubmHsu3tbWpra0t8jgYDEo6MvZvWdaxNL+brvoSXa9TEF9yhdo7Yk5cppw0XBkKd4vF6TEeK6/HJ3k/RuJzP7tijKc+199TuWzZMlVXV3c7vnHjRmVnZ9vyPevq6myp1ymILznaOqWut+T/ndqhzD4GdTPTPtRvf/vbXp93aoyJ4vX4JO/HSHzul+gYQ6FQzGVtTV5Gjhyp9PR0tbS0RB1vaWlRXl5ej+fk5eXFVb6ysjJqmCkYDCo/P1+lpaXKyck5xgiiWZaluro6lZSUyO/3J7RuJyA+e/W1uJwkHW7vlLZvliTNnF46oMXlkh2j3bwen+T9GInP/eyKsWvkJBa2Ji+ZmZmaMmWK6uvrNWvWLElSOBxWfX29ysvLezynqKhI9fX1uuWWWyLH6urqVFRU1GP5QCCgQCDQ7bjf77fthWNn3U5AfIkX73yWI20c+NuTa+h+Xo+R+Nwv0THGU5ftw0YVFRWaO3eupk6dqoKCAq1cuVKHDh3SvHnzJEnXXXedxo0bp2XLlkmSFixYoEsuuUQ//OEPNWPGDK1bt04vv/yyHn74YbubCtiGxeUAIHFsT15mz56tDz74QIsXL1Zzc7MmT56sDRs2RCbl7t27V2lpnw3uX3jhhVq7dq0WLVqk73//+zr99NP17LPP6pxzzrG7qcCgYLNEADg2gzJht7y8vNdhooaGhm7HvvKVr+grX/mKza0CEqe/+SxHr8/CZokAcGz4DQocI9ZnAYDBxcaMwDFiPgsADC56XoAEYj4LANiP5AXoB/NZAMBZ+C0L9IH5LADgPMx5AfrAfBYAcB56XoAYMZ8FAJyB5AUpjfksAOA+/CZGymI+CwC4E3NekLKYzwIA7kTPCyDmswCAm5C8wLM+P5/FsjrU1imF2jvkNz7mswCAS/HbGp7U+3yWDC3cvikpbQIAJAZzXuBJzGcBAO+i5wWuFM8tzl3zWSzLUm3tRpWVlcrv90eeZz4LALgLyQtcJ95bnLvms1g+o0C6lJ2ZIb+flz4AuBXDRnAdhoQAILXx5yccZyBDQr1hSAgAvIfkBY4y0CEhAEDqYNgIjsKQEACgP/zJCsdiSAgA0BOSFwwqdnEGABwrPhkwaNjFGQCQCMx5waBhPgsAIBHoeUHCcIszAGAwkLwgIbjFGQAwWBg2QkIwJAQAGCz86YuYMCQEAHAKkhf0iyEhAICTMGyEfjEkBABwEv48RmRIyLI61NYphdo75DefDeswJAQAcBKSlxTXfUgoQwu3b+q1PENCAIBk41PI42KZaMuQEADATWxLXj766CN961vf0m9+8xulpaXpX//1X/XAAw/o+OOP77V8VVWVNm7cqL1792rUqFGaNWuWli5dqmHDhtnVTE+Ld6LtS7ddoi0v1qusrFR+v7/b8wwJAQCcwLbk5eqrr9b+/ftVV1cny7I0b948feMb39DatWt7LL9v3z7t27dPy5cv18SJE/Xuu+/qm9/8pvbt26df/OIXdjXT0+KdaDviuEwF0qXszAz5/XTKAQCcyZZPqDfeeEMbNmzQn/70J02dOlWS9OCDD+pLX/qSli9frrFjx3Y755xzztEvf/nLyONTTz1Vd955p6655hp1dHQoI4MP089L9NorHR0dCW0fAAB2sCUjaGxs1PDhwyOJiyQVFxcrLS1N27Zt07/8y7/EVM8nn3yinJycPhOXtrY2tbW1RR4Hg0FJkmVZsixrgBH0rKu+RNc7EMYYXfmTP2nn3o9jKu/3heX39X5nfEdHh6Pis4PX45O8H6PX45O8HyPxuZ9dMcZTn88YYxL63SXdddddeuyxx7R79+6o46NHj1Z1dbVuvPHGfuv48MMPNWXKFF1zzTW68847ey23ZMkSVVdXdzu+du1aZWdnx994hzBGag/3/nx7WFr0cmy558lDjRac3SmmqwAAnCoUCmnOnDmRjou+xNXzcvvtt+uee+7ps8wbb7wRT5U9CgaDmjFjhiZOnKglS5b0WbayslIVFRVR5+bn56u0tLTf4ONlWZbq6upUUlLS44TWRIm3V+Wl2y7RkASsvTJY8SWL1+OTvB+j1+OTvB8j8bmfXTF2jZzEIq7k5bvf/a6+9rWv9VnmlFNOUV5eng4cOBB1vKOjQx999JHy8vL6PL+1tVWXX365hg4dqmeeeabfH0wgEFAgEOh23O/32/bCsbNu6cgicbEmLlPHn6Dc4ccl9C4gu+NLNq/HJ3k/Rq/HJ3k/RuJzv0THGE9dcSUvo0aN0qhRo/otV1RUpI8//lg7duzQlClTJEmbNm1SOBxWYWFhr+cFg0GVlZUpEAjo17/+tbKysuJpnmuwySEAAANny4Tds846S5dffrnmz5+vmpoaWZal8vJyXXnllZE7jd5//31ddtllevzxx1VQUKBgMKjS0lKFQiE98cQTCgaDkS6kUaNGKT3dG4ujsckhAADHxrZPxZ///OcqLy/XZZddFlmk7kc/+lHkecuytHv3boVCIUnSzp07tW3bNknSaaedFlVXU1OTJkyYYFdTE4oVbQEAsJdtycuIESN6XZBOkiZMmKCjb3SaNm2abLjxaVDF26vCkBAAAPFjPCJGxhiF2nvedblLvL0q/+e4TJITAADiRPISo8NWpyYt3aT+dl3uQq8KAAD2IHmxAb0qAADYh+QlRkP86fp/d3xRtbUbe911+eiyJC4AANiD5CVGPp9P2ZkZ7LoMAECS9b5THwAAgAORvAAAAFcheQEAAK5C8gIAAFyF5AUAALgKyQsAAHAVkhcAAOAqJC8AAMBVSF4AAICrkLwAAABXIXkBAACuQvICAABcheQFAAC4iue2RjbGSJKCwWDC67YsS6FQSMFgUH6/P+H1JxvxuZ/XY/R6fJL3YyQ+97Mrxq7P7a7P8b54LnlpbW2VJOXn5ye5JQAAIF6tra0aNmxYn2V8JpYUx0XC4bD27dunoUOHyufzJbTuYDCo/Px8vffee8rJyUlo3U5AfO7n9Ri9Hp/k/RiJz/3sitEYo9bWVo0dO1ZpaX3PavFcz0taWppOPPFEW79HTk6OZ1+UEvF5gddj9Hp8kvdjJD73syPG/npcujBhFwAAuArJCwAAcBWSlzgEAgFVVVUpEAgkuym2ID7383qMXo9P8n6MxOd+TojRcxN2AQCAt9HzAgAAXIXkBQAAuArJCwAAcBWSFwAA4CopnbysWrVKEyZMUFZWlgoLC7V9+/Y+yz/99NM688wzlZWVpXPPPVcvvPBC1PPGGC1evFhjxozRkCFDVFxcrLfeesvOEPoVT4yPPPKILr74Yp1wwgk64YQTVFxc3K381772Nfl8vqivyy+/3O4wehVPfGvWrOnW9qysrKgyTruG8cQ3bdq0bvH5fD7NmDEjUsZJ12/Lli2aOXOmxo4dK5/Pp2effbbfcxoaGnTBBRcoEAjotNNO05o1a7qVifd9bad4Y/zVr36lkpISjRo1Sjk5OSoqKlJtbW1UmSVLlnS7hmeeeaaNUfQu3vgaGhp6fI02NzdHlXPzNezpPebz+XT22WdHyjjlGi5btkz/8A//oKFDh2r06NGaNWuWdu/e3e95TvgsTNnkZf369aqoqFBVVZV27typSZMmqaysTAcOHOix/NatW3XVVVfp+uuv165duzRr1izNmjVLr732WqTMvffeqx/96EeqqanRtm3bdNxxx6msrEyffvrpYIUVJd4YGxoadNVVV+nFF19UY2Oj8vPzVVpaqvfffz+q3OWXX679+/dHvp588snBCKebeOOTjqwIeXTb33333ajnnXQN443vV7/6VVRsr732mtLT0/WVr3wlqpxTrt+hQ4c0adIkrVq1KqbyTU1NmjFjhi699FK98soruuWWW3TDDTdEfbgP5DVhp3hj3LJli0pKSvTCCy9ox44duvTSSzVz5kzt2rUrqtzZZ58ddQ3/8Ic/2NH8fsUbX5fdu3dHtX/06NGR59x+DR944IGo2N577z2NGDGi2/vQCddw8+bNuvnmm/XSSy+prq5OlmWptLRUhw4d6vUcx3wWmhRVUFBgbr755sjjzs5OM3bsWLNs2bIey3/1q181M2bMiDpWWFho/u3f/s0YY0w4HDZ5eXnmvvvuizz/8ccfm0AgYJ588kkbIuhfvDF+XkdHhxk6dKh57LHHIsfmzp1rvvzlLye6qQMSb3w//elPzbBhw3qtz2nX8Fiv3/3332+GDh1qDh48GDnmpOt3NEnmmWee6bPMwoULzdlnnx11bPbs2aasrCzy+Fh/ZnaKJcaeTJw40VRXV0ceV1VVmUmTJiWuYQkSS3wvvviikWT+93//t9cyXruGzzzzjPH5fOadd96JHHPqNTxw4ICRZDZv3txrGad8FqZkz0t7e7t27Nih4uLiyLG0tDQVFxersbGxx3MaGxujyktSWVlZpHxTU5Oam5ujygwbNkyFhYW91mmngcT4eaFQSJZlacSIEVHHGxoaNHr0aJ1xxhm68cYb9be//S2hbY/FQOM7ePCgxo8fr/z8fH35y1/W66+/HnnOSdcwEdfv0Ucf1ZVXXqnjjjsu6rgTrt9A9PceTMTPzGnC4bBaW1u7vQffeustjR07Vqeccoquvvpq7d27N0ktHJjJkydrzJgxKikp0R//+MfIcS9ew0cffVTFxcUaP3581HEnXsNPPvlEkrq93o7mlM/ClExePvzwQ3V2dio3NzfqeG5ubrex1y7Nzc19lu/6N5467TSQGD/vtttu09ixY6NehJdffrkef/xx1dfX65577tHmzZs1ffp0dXZ2JrT9/RlIfGeccYZWr16t5557Tk888YTC4bAuvPBC/fWvf5XkrGt4rNdv+/bteu2113TDDTdEHXfK9RuI3t6DwWBQhw8fTshr3mmWL1+ugwcP6qtf/WrkWGFhodasWaMNGzbooYceUlNTky6++GK1trYmsaWxGTNmjGpqavTLX/5Sv/zlL5Wfn69p06Zp586dkhLze8tJ9u3bp9/+9rfd3odOvIbhcFi33HKLLrroIp1zzjm9lnPKZ6HndpVGYtx9991at26dGhoaoia1XnnllZH/n3vuuTrvvPN06qmnqqGhQZdddlkymhqzoqIiFRUVRR5feOGFOuuss/Sf//mfWrp0aRJblniPPvqozj33XBUUFEQdd/P1SzVr165VdXW1nnvuuag5IdOnT4/8/7zzzlNhYaHGjx+vp556Stdff30ymhqzM844Q2eccUbk8YUXXqi3335b999/v372s58lsWX2eOyxxzR8+HDNmjUr6rgTr+HNN9+s1157LWnzp+KVkj0vI0eOVHp6ulpaWqKOt7S0KC8vr8dz8vLy+izf9W88ddppIDF2Wb58ue6++25t3LhR5513Xp9lTznlFI0cOVJ79uw55jbH41ji6+L3+3X++edH2u6ka3gs8R06dEjr1q2L6Zdgsq7fQPT2HszJydGQIUMS8ppwinXr1umGG27QU0891a2L/vOGDx+uL3zhC664hj0pKCiItN1L19AYo9WrV+vaa69VZmZmn2WTfQ3Ly8v1X//1X3rxxRd14okn9lnWKZ+FKZm8ZGZmasqUKaqvr48cC4fDqq+vj/rL/GhFRUVR5SWprq4uUv7kk09WXl5eVJlgMKht27b1WqedBhKjdGSW+NKlS7VhwwZNnTq13+/z17/+VX/72980ZsyYhLQ7VgON72idnZ169dVXI2130jU8lviefvpptbW16Zprrun3+yTr+g1Ef+/BRLwmnODJJ5/UvHnz9OSTT0bd5t6bgwcP6u2333bFNezJK6+8Emm7V66hdOROnj179sT0R0SyrqExRuXl5XrmmWe0adMmnXzyyf2e45jPwoRN/XWZdevWmUAgYNasWWP+/Oc/m2984xtm+PDhprm52RhjzLXXXmtuv/32SPk//vGPJiMjwyxfvty88cYbpqqqyvj9fvPqq69Gytx9991m+PDh5rnnnjP//d//bb785S+bk08+2Rw+fHjQ4zMm/hjvvvtuk5mZaX7xi1+Y/fv3R75aW1uNMca0traaW2+91TQ2Npqmpibzu9/9zlxwwQXm9NNPN59++qnj46uurja1tbXm7bffNjt27DBXXnmlycrKMq+//nqkjJOuYbzxdfnHf/xHM3v27G7HnXb9Wltbza5du8yuXbuMJLNixQqza9cu8+677xpjjLn99tvNtddeGyn/l7/8xWRnZ5vvfe975o033jCrVq0y6enpZsOGDZEy/f3MBlu8Mf785z83GRkZZtWqVVHvwY8//jhS5rvf/a5paGgwTU1N5o9//KMpLi42I0eONAcOHHB8fPfff7959tlnzVtvvWVeffVVs2DBApOWlmZ+97vfRcq4/Rp2ueaaa0xhYWGPdTrlGt54441m2LBhpqGhIer1FgqFImWc+lmYssmLMcY8+OCD5qSTTjKZmZmmoKDAvPTSS5HnLrnkEjN37tyo8k899ZT5whe+YDIzM83ZZ59tnn/++ajnw+GwueOOO0xubq4JBALmsssuM7t37x6MUHoVT4zjx483krp9VVVVGWOMCYVCprS01IwaNcr4/X4zfvx4M3/+/KT9UjEmvvhuueWWSNnc3FzzpS99yezcuTOqPqddw3hfo2+++aaRZDZu3NitLqddv67bZj//1RXT3LlzzSWXXNLtnMmTJ5vMzExzyimnmJ/+9Kfd6u3rZzbY4o3xkksu6bO8MUduDx8zZozJzMw048aNM7NnzzZ79uwZ3MD+Lt747rnnHnPqqaearKwsM2LECDNt2jSzadOmbvW6+Roac+TW4CFDhpiHH364xzqdcg17iktS1PvKqZ+Fvr8HAAAA4AopOecFAAC4F8kLAABwFZIXAADgKiQvAADAVUheAACAq5C8AAAAVyF5AQAArkLyAgAAXIXkBQAAuArJCwAAcBWSFwAA4CokLwAAwFX+P0QDgAX6yjwDAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.plot(cs.np.linspace(0,2,50+1),sol.value(X[0,:]).flatten(),label='x_1')\n",
    "plt.plot(cs.np.linspace(0,2,50+1),sol.value(X[1,:]).flatten(),label='x_2')\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "plt.figure()\n",
    "plt.step(cs.np.linspace(0,2,50+1),cs.np.hstack((cs.np.nan,sol.value(U))),label='u')\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the single shooting method is not robust and may be prone to failure. Multiple shooting methods are recommended."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "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.11.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}