{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "8835f784",
   "metadata": {},
   "source": [
    "# `evans` — Combined trip distribution + assignment (Evans 1976)\n",
    "\n",
    "**What.** Evans fuses trip *distribution* and *assignment* into one convex program:\n",
    "the OD matrix is endogenous, distributed by a doubly-constrained entropy/gravity model\n",
    "at the equilibrium travel costs, while the flows satisfy Wardrop UE. `evans` alternates\n",
    "an all-or-nothing assignment, an entropy distribution update, and a line search to the\n",
    "joint optimum (`[evans1976derivation]`,\n",
    "[docs/REFERENCES.md](../../docs/REFERENCES.md)).\n",
    "\n",
    "**Why it is in the benchmark.** It replaces the ad-hoc feedback loop between separate\n",
    "distribution and assignment steps with one convergent program (ADR-007). See the\n",
    "[model compendium](../../docs/MODELS.md) and\n",
    "[docs/ARCHITECTURE.md](../../docs/ARCHITECTURE.md) (P1).\n",
    "\n",
    "**Scope.** Runs on the built-in symmetric Evans anchor (2 origins × 2 destinations,\n",
    "β = 0.5) whose equilibrium collapses to a binary logit split, and certifies the gap."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "62e0cf60",
   "metadata": {},
   "source": [
    "## How this notebook is graded\n",
    "\n",
    "**A notebook never claims a number it does not compute in that cell.** Every scored\n",
    "quantity below is recomputed live by the P1 `Evaluator` from the flows the model\n",
    "emitted, in the cell where it is claimed. Model self-reports are shown only as\n",
    "provenance and diffed against the certificate as an honesty check, exactly as the\n",
    "harness treats them ([README](../../README.md), *Certified, not self-reported*)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "b2e434eb",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-21T13:45:43.167281Z",
     "iopub.status.busy": "2026-07-21T13:45:43.166896Z",
     "iopub.status.idle": "2026-07-21T13:45:45.128781Z",
     "shell.execute_reply": "2026-07-21T13:45:45.127643Z"
    }
   },
   "outputs": [],
   "source": [
    "# Setup. `evans` is a core model: a plain `pip install -e .` suffices — no\n",
    "# optional extra, so no guard cell. The inline backend is Agg-based (headless CI\n",
    "# renders into the notebook); NEVER matplotlib.use(\"Agg\") in-kernel — it silently\n",
    "# suppresses inline figure capture.\n",
    "%matplotlib inline\n",
    "import numpy as np\n",
    "\n",
    "from tabench import (\n",
    "    Budget,\n",
    "    Evaluator,\n",
    "    EvansCombinedModel,\n",
    "    RngBundle,\n",
    "    Trace,\n",
    "    evans_symmetric_scenario,\n",
    "    viz,\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "efbfb52c",
   "metadata": {},
   "source": [
    "## The scenario\n",
    "\n",
    "Two origins (1, 2) and two destinations (3, 4), each producing/attracting 10 trips,\n",
    "one congestible link per OD pair. By symmetry the doubly-constrained gravity collapses\n",
    "to a binary logit split at dispersion β = 0.5 — a scalar fixed point. Content-hashed (P2)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "466140ec",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-21T13:45:45.132183Z",
     "iopub.status.busy": "2026-07-21T13:45:45.131773Z",
     "iopub.status.idle": "2026-07-21T13:45:45.137186Z",
     "shell.execute_reply": "2026-07-21T13:45:45.136664Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "scenario      : evans\n",
      "content hash  : 69a9ca0453081aee…\n",
      "zones         : 4  (2 origins, 2 destinations)\n",
      "total trips   : 20.0  (gravity dispersion beta = 0.5)\n"
     ]
    }
   ],
   "source": [
    "scenario = evans_symmetric_scenario()\n",
    "net = scenario.network\n",
    "\n",
    "print(f\"scenario      : {scenario.name}\")\n",
    "print(f\"content hash  : {scenario.content_hash()[:16]}…\")\n",
    "print(f\"zones         : {scenario.demand.n_zones}  (2 origins, 2 destinations)\")\n",
    "print(f\"total trips   : {scenario.combined_demand.total}  (gravity dispersion beta = {scenario.combined_demand.beta})\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6c92508d",
   "metadata": {},
   "source": [
    "## Solve\n",
    "\n",
    "The model contract ([CONTRIBUTING.md](../../CONTRIBUTING.md)): a model receives\n",
    "`(scenario, budget, rng, trace)`, records checkpoints, and respects the budget.\n",
    "Budgets are hardware-free (iterations / shortest-path calls; wall-clock is recorded\n",
    "but never the ranking axis, P7). Whatever the model writes into `self_report` is\n",
    "provenance, not a score."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "a5d6a444",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-21T13:45:45.139512Z",
     "iopub.status.busy": "2026-07-21T13:45:45.139134Z",
     "iopub.status.idle": "2026-07-21T13:45:45.145766Z",
     "shell.execute_reply": "2026-07-21T13:45:45.145279Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "model              : evans\n",
      "emitted flows      : [6.917388 3.082612 3.082612 6.917388]\n",
      "self-reported gap  : 0.000e+00  (provenance only)\n",
      "self distribution gap: 4.441e-16  (provenance only)\n"
     ]
    }
   ],
   "source": [
    "model = EvansCombinedModel()\n",
    "bundle = model.solve(scenario, Budget(iterations=200), RngBundle(0), Trace())\n",
    "\n",
    "final = bundle.final\n",
    "print(f\"model              : {model.name}\")\n",
    "print(f\"emitted flows      : {np.round(final.link_flows, 6)}\")\n",
    "print(f\"self-reported gap  : {final.self_report['relative_gap']:.3e}  (provenance only)\")\n",
    "print(f\"self distribution gap: {final.self_report['distribution_gap']:.3e}  (provenance only)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "22d30e7c",
   "metadata": {},
   "source": [
    "## Certify (P1)\n",
    "\n",
    "The harness certifies the relative gap of the joint program (assignment + distribution)\n",
    "from the emitted flows. We recompute the symmetric equilibrium split in-cell with\n",
    "`brentq` — the gravity collapses to `p = trips / (1 + exp(β(c_near(p) − c_far(q))))`,\n",
    "`c_near = 1 + 0.1 p`, `c_far = 3 + 0.1 q`, `q = trips − p` — rather than quoting it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "607e0817",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-21T13:45:45.147995Z",
     "iopub.status.busy": "2026-07-21T13:45:45.147699Z",
     "iopub.status.idle": "2026-07-21T13:45:45.153806Z",
     "shell.execute_reply": "2026-07-21T13:45:45.153303Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "certified relative gap : 0.000e+00\n",
      "feasible               : 1\n",
      "logit split p, q       : 6.9174, 3.0826  (recomputed via brentq)\n"
     ]
    }
   ],
   "source": [
    "from scipy.optimize import brentq\n",
    "\n",
    "evaluator = Evaluator(scenario)\n",
    "metrics = evaluator.evaluate(final.link_flows)\n",
    "gap = metrics[\"relative_gap\"]\n",
    "print(f\"certified relative gap : {gap:.3e}\")\n",
    "print(f\"feasible               : {metrics['feasible']:.0f}\")\n",
    "\n",
    "assert metrics[\"feasible\"] == 1.0\n",
    "assert abs(gap) < 1e-10\n",
    "\n",
    "# Honesty diff (P1).\n",
    "assert np.isclose(final.self_report[\"relative_gap\"], gap, rtol=1e-9, atol=1e-12)\n",
    "\n",
    "# Analytic anchor RECOMPUTED: the symmetric logit split via brentq.\n",
    "trips = 10.0\n",
    "beta = scenario.combined_demand.beta\n",
    "\n",
    "def _split(p):\n",
    "    c_near = 1.0 + 0.1 * p\n",
    "    c_far = 3.0 + 0.1 * (trips - p)\n",
    "    return p - trips / (1.0 + np.exp(beta * (c_near - c_far)))\n",
    "\n",
    "p = brentq(_split, 0.0, trips)\n",
    "q = trips - p\n",
    "ref_flows = np.array([p, q, q, p])  # link order 1->3, 1->4, 2->3, 2->4\n",
    "print(f\"logit split p, q       : {p:.4f}, {q:.4f}  (recomputed via brentq)\")\n",
    "assert np.allclose(final.link_flows, ref_flows, atol=1e-3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "037888f4",
   "metadata": {},
   "source": [
    "## Visualize\n",
    "\n",
    "Both figures come from `tabench.viz`. Left/top: the Evans equilibrium link flows.\n",
    "Right/bottom: the emitted flows against the `brentq`-recomputed symmetric split —\n",
    "on-diagonal means the combined program reached the joint optimum."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "d91c725a",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-21T13:45:45.155990Z",
     "iopub.status.busy": "2026-07-21T13:45:45.155709Z",
     "iopub.status.idle": "2026-07-21T13:45:45.407278Z",
     "shell.execute_reply": "2026-07-21T13:45:45.406523Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAANUAAADWCAYAAABG1MctAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAgmklEQVR4nO3deXxMZ///8ddkspAJWQhiiyWxJgi3NUn17i1FLKXEVkWtlVJKY7u19G4ttX7vXy213L2rRC2hllZLLNEERS0tIqJB7BFrSUIymfP7I7dhKsuEk5nB5/l45FGZc53rfM6k71xnTq5zjkavz1IQQqjGztoFCPGikVAJoTIJlRAqk1AJoTIJlRAqk1AJoTIJlRAqk1AJoTIJlRAqk1C9wKZNm0ZQUJCqfbq6uvH9998DkJycjKurG7///nuR1pSSksIbb3TCy6s8lStXfqIOWyOhes60a9eOcePGmdV2+PDhbNq0qchqqVixIomJp6hTp06RbQNg/vwFpKSkEBcXy6FDh4p0W2qwt3YBQn2KopCdnY2Li0uRbker1VK2bNki3QbA2bNnadCgPtWrVy/ybamhUCOVwWBg9uw5+PvXo2zZcgQGBrJhw0bjstq167B06X9M1vntt99wc3Pn/PnzAMybN4/mzVvg5VWeOnXqMmrUaO7du2dsHxkZSeXKldm+fQeNGzehfPkKvPlmF65evWpsExsby9///prxcOD111sb+89NXu2Tk5Nxc3Pn8OEjJu0XLFiAn58fBoOB2NhYXF3d2L59B0FBwZQtW4727TuQmppKdHQ0jRs3oWLFSgwYMJD09HRjH+3atSMiIoJx48ZRubI3Pj6+fP31MtLS0ggPD6dChYo0aBBAdHS0ybbj4+Pp0qUr5ctXwMfHl8GDB3Pjxg0Ahg4dSlzcHhYu/BJXVzdcXd1ITk421hgdHc0rr7TE07MM+/bty/VQa/ny5TRt2gxPzzLUqFGTDz+MKPDnnpe/Hv49rCMmZjctW75KuXJehIS8zunTp/Ps48yZs9SrV58PP4xAUZ6c2+3v78+mTZv49ttVuLq6MXTo0Fz7OXHiBO3bd6Bs2XJUqVKV998fYfz/Kj4+Hjc3d65fvw7AzZu3cHNz5513+hvXnzlzJq1btwHg1q3bDBw4iGrVqlO2bDkCAhqyYsUKs9+XQoVq9uw5rFq1irlz5/LLL78QHh7O4MGDiYuLw87Oji5duhAVFWWyzpo1a2nWrJnxWNjOzo7PP5/OL7/sY+HCBfz88898/PEkk3XS0zP44osvWLRoEVu2/MDFixeZOPEjAPR6PW+99RaBgYHs2bOH6Oho+vXri0ajybXm/Np7e3vz6quvEhkZabJOZGQkvXr1ws7u0dszffp0Zs2aybZtW7l06RL9+vVjwYKFLF26hDVrVrNz504WLVps0s+3366iVKlS7Nq1kyFDBjNq1Cj69u1HkyZN+fnn3bz22msMHjzEGMbbt2/ToUNH6tWrR0zMLtati+LatVT69u1nrKFJkyb07duXxMRTJCaeomLFisbtTZ48mcmTJ3PgwAHq1vV74r1YuvQ/fPhhBP369WXfvr18++23VKtW1bh86NChtGvXLtf3sTA+/fRTpkz5jJiYXdjba3nvvWG5tjt+/Dht2rQhLCyMWbNm5voz3LVrF61ataJz584kJp5i+vTpT7RJS0vjzTe74Obmxq5dO1m27GtiYmKIiMj5hVG7dm08PDyIi9sDwL59e/Hw8GDPnj3GPuLi9hh/AU2ZMoWEhASioqI4cGA/c+bMoVSpUua/AXp9lmLOV1raPcXT01PZu3ePyevvvvuu0rdvX0Wvz1IOHz6kuLi4KGfPnlH0+iwlM/OB4uvrqyxatCjPfqOiopRKlSoZv1+27GtFp9MpiYmJxtcWLlyoVK1aVdHrs5Rr164pOp1OiYnZZVbdBbVfs2aNUqFCBSUt7Z6i12cpv/56UHFxcVGSkv5Q9PosZdeunYpOp1O2b99uXGfGjBmKTqdTTp8+bXxt2LBhSseOHY3ft279utKqVSvj9w8e3FfKlCmjDBjQ3/japUsXFZ1Op+zbt1fR67OUadOmKh06dDCpLzn5nKLT6ZSEhJPGfj/8cLRJm4c1bty4weT1Tz/9l9KsWVPj9z4+PsqkSZPyfK8mTpxoUl9uXzqdTtmw4TtFr89SkpL+UHQ6nXL48KE836stW35QdDqdcu/eXZOa9uyJUypUqKDMnTunwJ9ht25hyqBBA/OsY+nSpUqFChWUO3dum2y3RIkSyuXLlxS9Pkvp0aOHMnLkSEWvz1IiIj5Uxo8fr1SoUEGJjz+hZGSkK56ensq2bdsUvT5L6dq1izJkyGCz/v/K7cvskerMmTOkp6fTqVNnypevYPxatWoVZ8+eBaBevXrUrFmTtWuj/pf+OFJTU+ncudNjv3li6NChI7Vq1aZChYoMGTKEmzdvmhw6OTs7m/wGLVeuLKmpqQB4eLjz1lu9ePPNLnTv3p2FCxcaDw0vXLhgUtusWbPzbQ/Qvn07tFotmzfnnElauXIlwcHBeHt7m+y/n19d47/LlPHE2dmZqlWrmLz2sMaH6tZ9tI5Wq8XDw506dR7vpwyAcb1jx44TGxtrsg+NGzcBML7H+QkICMhzWWpqKleuXKFly1fybDN58iQWLVpU4HYK8vh79fAz1+PvzYULF+nUqTNjx45h+PDhz7y9xMRT+Pv7odPpjK81bdoUg8FgPPQMCgokLi4OyBmVWrZ8hcDAFsTGxnH48GGysrJo1qwpAAMGDGDduvUEBQXx0Ucfs3///kLVY3ao0tLSAFizZjWxsbHGr/379/PNN8uM7cLCwoiKWgvA2rVRtGr1Dzw8PICcY/Du3bvj51eX5cu/YffuGGbNmglAZmamsQ8HB9PzJxqNxuR4e8GCBURHb6NJk6asX/8djRr9jYMHD+Ll5WVSW//+/fNtD+Do6EjPnj2IjIwkMzOTtWujePvt3k/sv739o5o0Gk0eNRpMXsutzeOvPTzcMRgMxve4TZs2JvsQGxvL4cOHadGixRM1/ZWzs3Oey4oVK1bg+mr563sFj/YRoHTp0jRq1JCoqHX8+eefFqkpKCiIhIQEkpKSOHXqFM2aNSMoKIi4uDji4vYQEBBgfP9CQkI4fvwY4eHhXL16lY4d3+Cf/5xo9rbMDlXNmjVxcnLi4sWLVK9ezeTr8eP6sLCuxMef5MiRo2zcuJGwsG7GZUePHsVgMDBlyhQaN26Mj4+PyahRGPXr12f06FFER2+jdu3arF27Fnt7e5O6PDzc823/UJ8+fYiJiWHp0qVkZ+vp0KHDU9X0rOrXr09CQgLe3pWfeI8f/hZ2cHAkOzu70H2XKFGCypUrs3v3z2qXXWjFixdj9erVFCvmxJtvduHu3bvP1F+NGjU5duy48Rc/wP79+7Gzs8PX1xfIOWpwc3Nj5syZ+Pv74+LiQlBQEHv27CEuLu6JEzqlS5emV69eLFmymGnTprFs2TLMZXaoSpQowfDhwxk/fgIrV67kzJmzHD16lEWLFrFy5UpjO29vb5o2bcrw4cMwGAyEhrY1LqtWrRpZWVksWrSIs2fPsWrVKr766r9mFwtw7tw5Jk/+hAMHDnD+/Hl27NhJUlISNWrUfOr2NWvWpHHjxkyaNJkuXbpSvHjxQtWklkGDBnLr1i369x/AoUOHOXPmLNu37yA8PNwYpMqVK/Prr4dITk7mxo0bJiNAQcaPH8e8efP48ssvSUpKMv78Hpo8+ROGDBmi+n7lRqfTsWbNGuzt7enaNczkDHBhdesWRrFixXj33aHEx8fz888/ExExhh49uhsPsTUaDYGBLVizZq0xQH5+fjx48IDdu3cTGBho7G/KlCn88MMPJCWd4eTJk2zd+hM1atQwu55Cnf2bOPGfjBkTwZw5c2nSpAldunRl69ZtT3z+CAsL49ix47Rv387kf1B/f3+mTp3C//3fv2nevDlr1qxl0qSPC1MCzs7OnD6dyNtv96FRo78xcuQIBg0aSP/+7zxT+7ff7k1mZia9ez956GcpXl5ebNu2lezsbDp37kyLFi0YP348rq6uxjOR778/HK1WS9OmzahWrToXLlwwu/9evXoxbdo0li79D02bNqN79x4kJSUZl6ekXOXixYuq71deXFxciIpai6IodOvWzWSkKQxnZ2fWr1/HrVu3+PvfX6NPn760bNmSmTNnmrQLDAwkOzub4OCcUNnZ2REY2AKNRmP8PAU5Hwk++eRfBAYG0rZtKFqtlq++Mv1TUX40cuOXHDNmzGDDhg3s3bvX2qWI59xLP03p3r17xMfHs3jxEgYPtsyhj3ixvfShioiIoGXLVwkKCsr1rJ8QhSWHf0Ko7KUfqYRQm4RKCJVJqIRQmYRKCJXJRYoWVm/iaRy0uV+mUlSyshV+/8zXott8mUmoLMxBqyF2ovlTXtQQ/FmiRbf3spPDPyFUJqESQmUSKiFUJqESQmUSKiFUJqESQmUSKiFUJqESQmUSqufIqpUr8Cpdkh+32OaN+UUOCdVz4sL5ZCKXf02jvzW2dimiABKq54DBYGD0yOF8Nm0mjo5O1i5HFEBC9RxYtGAejZs0pX6DvO9AK2yHTKi1cQkn4/nh+418t/kna5cizCShsnH79+3lwvnztGiSM0qlXkshYtT7XLt6lb79B1q5OpEbCZWN69t/oEl43uwYyqB3w2kb2t6KVYn8yGcqIVQmI9VzZv2mLdYuQRRARiohVCahEkJlEiohVCahEkJlEiohVCahEkJlEiohVCahEkJlEiohVCahEkJlMk3JwrKyFYvf2zwrWx6WaUkyUlmYpZ/4Ya1tvszkmb9CqExGKiFUJqESQmUSKiFUJmf/bNzZyzd5f9YGEs9fR2unIaxVPT7q/w/s7eVHZ6vkRIUN23HgNF3GfsNff0Bl3F04uWa0BMtGyeGfjXqQmUWPiZFPBArg2q17dJsQafGahHkkVDZq7rexZGZl57l8169J3LiTZsGKhLkkVDYoM0vPovX7821jUBS+WLPHQhWJwpBQ2aAVPx7h5p/pBbZbvH6/jFY2SEJlYzKz9MyO3G1W27T7mTJa2SAJlY1Z8eMRLl67Y3Z7Ga1sj4TKhhRmlHpIRivbI6GyIYUdpR6S0cq2SKhsxNOMUg/JaGVbJFQ24mlHqYdktLIdMk3JwhpN+uOJ1xRDNskHv0P/4NlC4VbJj9JVG+W67NAnPs/UtzCfjFQWltul7em3Lj1zoAD+vHIaRXmyf7mc3rJkRqaFOWg1xE6sYfLapC//YPGJZ+/boH/Akt6u1K5azuR1S98T42UnI5UNcC7moGJfjqr1JZ6OhMoGNPOvqko/5T1dqVjGTZW+xNOTUNmAlo18aB/s90x9aDQapg3rgFYrP1Jrk5+AjVg4oRsRff5BuVIlCrWe1s4Of5/yrJs5gNeb1y6i6kRhyIkKG2Gv1TKq92uM6v0a12+nsWHXb3y08Ic822vtNBxYHoFbieLyOcrGSKhsUGk3He4lnfNto9FoKO/paqGKRGHI4Z8QKpNQCaEyCZUQKpPPVBZw7tw5hg0bxrVrqZy/rif9gzicdTrj8nVrV/PFv+dgMBjo0bM34cNH5Cw4uxXSr0HWPbBzAAddzn8rtYQSZa20N6IgMlJZQHh4OBMmTODAgf1UCFuMo5OTcdmNGzeYMX0KG77/iV2xvxAXt5s/Tp/OWehRC+r0hrp9oZg7VAyCsgGQHG2lPRHmkFAVsZMnT2Jv70CLFi0A0BZzNbkJ5vnks/j61sDNzR2tVkuz5oH8uGVzzkJXb9BowMEZintCZhpkZwLyaBxbJod/RSwpKQkXFx3du3fn8uUr3HQJAqYbl1epWo2EhJNcuXIZd3cPdm6Ppq6fH+Ufn2CRnQk34uHWadDYgU8Hi++HMJ+Eqojp9dns3buPuLhYPD098W4Yyu6YnbR89TUA3N09+HTq57zzdk8cnZyoW9cPO632UQeKAsnbwasJlG0IN07Cpb1Qs5N1dqgAmzdvJjg4GDc3N2uXYjVy+FfEypf3IiAggIoVK+Lk5IRzlUBOHDtm0qZtaHt+2r6bTT9so0zZclSrVv3Rwkt7QeuUEyiAUrXh7kUUfYYF98J8s2bNxte3Bq+80pKJEz9i27Zt3Lt3z9plWZSEqog1bNiQ1NRUbt26jcFgIOPSEXxrmF5PdT01FYBrKSls2rCeTm92zVmQegzSU8Cr8aPGt5PAvhhoi1lqFwpl9+4YEhMTGTMmgszMB0yaNJmqVavx+uutrV2axcjhXxGzt7dn0qSPCQ0NRVEUHNwCCGndlrd6dGH23HmU8/Ji/NjRnEo4iVar5eNPPsPd3QM4Dxd2g6MLnIgEFLAvDsXcoHoHNBrbPVnh7u5GjRo1uHLlKlevpnD58mUMBoO1y7IYCZUFhISEEBISAjy6R0XkqnXG5Uu++ib3FRsOK/La1DZgwED27NlDqVIetGzZkp49ezBv3heULFnS2qVZjIRKqGrXrl2ULFmSVq1CCA4OpkWL5jg75z85+EUjn6mEqs6cSWLFiuWULVuGpUuX4O9fj9at2zB16lRrl2YxMlIJ1fn5+VGlShV8fHyoWrUqkZGRHDx4kAkTJli7NIuQUAlVTZ78CXFxcfz222/4+voQHBzM/PnzCQoKtnZpFiOhEqr6888/ee+99wgODqJ06dLWLscqJFQviczMTAAcHYv20vs5c2Yb/33jxg1KlSpVpNuzRXKi4iWgMWQR1vNthr031BiuopKRkcHIkR9QrpwXPj6+lCvnxciRH5CW9vLc511C9TJQsnGy03P79m30en2RbmrChH/yxx+n2bRpI6dOJbB58yaSkpKYOPGjIt2uLZHDP6GqH3/8kb179+Lh4Q5AmTJlWLbsa5o3b8HcuXOsXJ1lSKgsLCtbMeve5ndTrua7PNtgXj8A2c/wgIJfDx5k9epVXL16BWdnZ7p0Dct3Hp+iKNjZmU6h0mjscn1wwotKQmVhDlrz5uw5FPCTKczcP3szt/lXR48cYenSxQx/fwS1a9UmPSODO3du57tO69at6dOnL5MnT6Zy5UokJ5/n008/pU2bNk9Vw/NIQmVh5j4nasVP93gvnyeB2GnM7ys9PZ1+fc1qamL16m/p2jWMunVzrph0cXHBxcUl33WmTp3C2LFjadOmDVlZWTg6OtK1axemTPms8AU8pyRUNsrxKUeXgmRkmF6H5eDggIPDk08duX//PmfOnCEg4CYj3h9GenoGtWvX5p3+A3B3d8+zfxcXF+bPn8+8efO4fv06pUuXtukZ9UVBQvWSGfruYJPvu4Z1o1u37k+0S0tLQ1EUDh48wMSPPqaESwkWL1nMF//v33w8aXKB29FoNHh6eqpV9nNFQvWSWfjlYooXL278PrdRCqBYsZyLINu2DcXTswwA3bp1Z8T7w7h//75xOUDlyt5mjUbJyeeeofLnh4SqiCUnJ9O0aTMURUFRFDp16sTixYtM2qxYEckHH3yAwWCgefPmfP/9/+6mdPUwXPkl59929lCzW85Fis+gePHiZl2KodPp8p5m9JczeStXRj5TTS8aCVUR8/Ly4vjxY5QuXZrU1FRq1qzF2LFjqF790X0oIiIiWLJkMaGhoVSpUpVNmzZBce+cQFULBdcqkLAm50pg3zcsVnurViH8+OOPNGgQgIuLC1FRa/Hz86fYYyMdwKRJk9mxYzsA06dPZ9y4cRar0RbJjIoi5ujoaPyNf/duzg1QHr+0/OjR3zAYDHTq1AlHR0deffVV/vvfrx918PAGL4asnEvr/8cSf/bp1Kkz/v7+RESMZujQIWQ+eMDw4e8/0e706dPGmRrz5s0v+sJsnIxUFpCcnEyTJk24f/8BnTt3wtfX17js5MmTlCjx6EFv3t7exMTEcPXm3Zw70iZvh+QdOYd/tR6dUMi2wD0f7LRa+vTtR5++/fJtFxwcRFBQMNWrVycjI4O33uqda7vIyBVFUKXtkVBZgLe3NykpKcTHxxMS8jrx8fHUqVMn33XW7fgdrv4KVVqDRw1IWA1/bIYanY1tfo2/wN/qVCrq8gv01VdfsXHjRpKTk9m2bRv+/s/2qNXnnYTKgurUqUPFihWJjFxp/GNo7dq1uXv3rrFNcnIyru6l2HfsBOjv5wQKoHQ9uBxn0t/0b3YRNb2PxerPi5OTE926dQPg9u078pnK2gW86OLj47ly5QoAFy5c4OzZszRt2sS4vEGD+mg0GjZs2EBmZiYxMTEUrxAAjiVByYa7F3Ma3koEB9PnAUfvP82v8Rcsti/meJlmTuRFQlXEDh8+QoMGAZQtW5ZGjf5Gx44d6dixI9WqVefo0d8A+Pzz6QwaNBgvr/LUrutH7Dm7nM9QZRrAH5vgyALIuAZVXn+i/+nf7LLwHomCaPT6rJdn+vBzYNTczfxn04FCrbNj/uB8P1vlzP17G4Cvly0v9C3DMh88YPToUdy9+ydfL1teqHVfRjJS2ZCL1+7wzZZDhV6vqEer1atXvbRTjp6GhMqGzIn8mSx9dqHXK8xnq4yMDNLT041fWVlZ+bY/k5TE0aNHeaNTJ7P6X7Jkaa6vjxgx0qz1XwRy9s9GPO0o9ZC5ZwLNnVALkJ2dzaJFCxkwcKDZFxnOnz+fMmU8eeONRzM/Ro0azYkT+VzH8oKRUNmIpx2lHno4WhX0dytzJ9QCbNq0kSpVqlKnTl1OnDhuVh1RUWvp2PENSpUqRVBQEGPGjOHw4cNs3LjBrPVfBBIqC3v4gILH6R+kce7Ar8/cd4ePtlDe7x9PvK7Jvk+1//3b3Am1V69cIXrbVmbMmFWoGnx8fFi+/Bveeqs3r7wSTELCKTZt2oirq2uh+nmeSagsLCtbeeKS+rQbF0F59mlH6bcuYcjOwk5rOvron+IeFQkJJ7lz5w4jRgzP6UOfzf37GQzo349x4yfg6/voGVvHj5uOYk5OTgwZMoQvv/ySpUuXcPHiRS5evIif38sx00JCZWEOWg2xE00f+jb+iwS+fnIAKzxFYX5PVxrUqGjy8iuf/F7orpq3CMS/Xn3j94mJp/hy4UJmzJyN618eixMUFIxGo8n1c1eHDh2BnIsWb926Weg6nkcSKhtQqWzel6dbqy8nJyecnJyM35csWRKNhlzvOHv79i1VtvmikFPqNqBDSz+KO+V9wsDBQcucMb34eGinfPsJaVoTj5JF8yyounX95A+/ZpKRygZUKuvO5yPeIGLuBh5kPXkHWY1GQxP/6lxJvZ1nH76VPZk6vGPul7VrtManxdvbF+2P/PLly0yZMoWjR48arx976PfffyvSbdsKCZWNCGsVQD3fCkRtP8KJpCtcunabxPOpebbXAC3qV6O0m45GdSrTO7RxnqOdYufAvPkLgaJ/QMHgwYMpXtyZkSNHvnRPUHxIQmVDanqX4Z8Dcu7+um7HUYZ9vjbPtlqtHVEzB5jdd1GH6aGjR3/jzJkki23PFslnqudAdnbO6fZrN/+0ciUFq1WrFikpKdYuw6pkpHqOpNy4Y+0SCtShQwd69uzJwIGDKFPGdBJuaGiolaqyLAmVUNXSpTkTamfPnm3yukajkVAJ9Zw7d45hw4Zx7Voq56/rSf8gDmedzrh83drVfPHvORgMBnr07E348BE5C85uhbQUsvXptG2+DOeSHnAnCyq9CjoP6+xMAY4dK/wfml808pnKAsLDw5kwYQIHDuynQthiHB/7o+qNGzeYMX0KG77/iV2xvxAXt5s/Tp/OWehRC2r3RFOtNdV8a/F6t3BwrQrnd1ppT4Q5ZKQqYidPnsTe3oEWLVoAoC3mavK3ovPJZ/H1rYGbW85MiGbNA/lxy2bK+70Grt4A2LlXx7emF3dvXwddObh2xPI7ko/WrduwdetPwKMpS7mJjf3ZkmVZjYSqiCUlJeHioqN79+5cvnyFmy5BwHTj8ipVq5GQcJIrVy7j7u7Bzu3R1PXzo/xjc0+V7Ez27Y6m67BpsOE7cK325IasaODAR6f2w8OHWrES2yChKmJ6fTZ79+4jLi4WT09PvBuGsjtmJy1ffQ0Ad3cPPp36Oe+83RNHJyfq1vXDTqt91IGiYDi3jU4D+/H7vmh4cAd8/26lvcldWFiY8d+9evWyYiW2QUJVxMqX9yIgIICKFXNmjjtXCeTEsWPGUAG0DW1P29D2AMydPcM4pQiAS3vRaIuhKAoJR+Kgekewy3ueoDVs2bLFrHZy9k+oomHDhqSmpnLr1m1cXUuScekIvjWGmbS5nppKaU9PrqWksGnDetZv2sLOw+ch9RhkpIJrZbZv+Y7eH8xg5jc7rLQneRs7tuCbZ8opdaEae3t7Jk36mNDQUBRFwcEtgJDWbXmrRxdmz51HOS8vxo8dzamEk2i1Wj7+5DPc3T2A8zlP+XB0Qbl4gTMODiz6ZBDcV0Cjhbo9rL1rRnIa3ZSEygJCQkIICQkBHl1OH7lqnXH5kq++yX3FhjkjmlZrR8zX/2Tb3mP8a+GGIq1VPDv5O5UQKpNQCaEyCZUQKpNQCaEyCZUQKpNQCaEyCZUQKpNQCaEyCZUQKpNQCaEymaZkYVnZCsGfJRbY7m7K1XyX6w3m9fNwm8JyZKSysL8+8eNF3ebLTEYqCzv0iY9Z7dZEpzFoat7L7e00ZvclLEtGKiFUJqESQmUSKiFUJqESQmUSKiFUJqESQmUSKiFUJqESQmUSqudAdraBvUdPc+jEWWuXIswgMyqeE2Nmr7J2CcJMMlIJoTIJlRAqk1AJoTIJlRAqk1AJoTI5+2dhDx9QUJC7KSn5LtcbFLP7AvOv4xLPTkYqCzP30naH4iXyXe5YvKTq2xTqkJHKwhy0GmIn1iiwXZa+OvW77+LW3Yxcl78dUofP3iu4H8Dse1kIdchIZaMc7LXMG9ct1ye9+1TyZGy/VlaoSphDQmXDXmtcg/WzBtKykQ+l3XSUK1WCoV2DWD9rICV0xaxdnsiDHP7ZuGb+VVg17R1rlyEKQUYqIVQmoRJCZRIqIVQmoRJCZRIqIVQmoRJCZRIqGzdxfASNA/zwKl2S48d+t3Y5wgwSKhvXrkMnNv6wlYqVKlu7FGEm+eOvjWveItDaJYhCkpFKCJVJqIRQmYRKCJVJqIRQmYTKxkWMGkFD/1pcuXyJnt0607xxfWuXJAogZ/9s3Mw5/7Z2CaKQZKQSQmUSKiFUJqESQmUSKiFUJqESQmUSKiFUJqESQmUSKiFUJqESQmUyo8LCsrIVi9/bXB5QYFkyUlmYg/bJe6O/iNt8mWn0+iz5NSaEimSkEkJlEiohVCahEkJlEiohVCahEkJlEiohVCahEkJlEiohVPb/AXmF5K4R3fBSAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 80x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcIAAAHXCAYAAAAvA1ptAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAABJlUlEQVR4nO3dd3wUZf4H8M/2bMumBwgJHU4UKRaQIhYUToqA0uQsFBFBUSwnKkrxEH8qnsfRROygCAgHdk5EQRBQPMBCBwkRSLLp2Z2tM78/YpaEBNhNdjNbPu/XK69Xdmd25vvsJvnkmXnmGYXH45ZAREQUo5RyF0BERCQnBiEREcU0BiEREcU0BiEREcU0BiEREcU0BiEREcU0BiEREcU0BiEREcU0BiEREcU0BiFRHU2bNg3333+/X+uuWLECPXv2DHFFdfP22++gbdt2aNIkA3v37pW7HKIGp5a7ACKSj9vtxhNPPIF169aie/fucpdDJAv2CImilMfjgSRdeCrh3NxcOBwOtG/fvk77cLvddXodUThhEFLU6tChA+bNewXXXXc9Gjdugttuux2FhUV45JFHkZWVhc6du2Dnzp2+9cvKyjBlykNo27Yd2rZth4cfngqbzeZbvm3bNlxzTXc0aZKB0aP/hrKy8mr7O3bsOEaMGIGWLVvhsssuw0svvQRRFC9a58KFCzFgwMBqz3300VpceeVVAIA9e/bgxhv7oGnTTLRo0RIjRow477YslgQsXboU3bpdg8aNm6C8vPy8de3duxdXXXU1AKB9+0vRsWMnAEB5eTkee+xxXHrpZWjVqjXuu+8+lJSUAABOnDgBiyUBy5cvR6dOnXHJJe19NQ4YMADNmjVHp06d8fbb7/hqmjt3LkaMGIHHHnscWVlZuPTSy/DRR2t9y0VRxJIlS3DllVchI6MpOnfugq+++goAIEmSb1lWVhb69++PgwcPXvQ9JQoEg5Ci2rp1a7F8+Xs4cGA//vjjD/Tp0wfXXdcbx48fx7Bht2Pq1Km+dadNm4Zjx45hx47v8f3323H48CE8+eRTAICiomKMGjUK9957L7KzT+BvfxuNVatW+V5rt9tx662D0Lt3bxw4sB+ff/45PvpoLZYvX37RGocNG4YdO3YgJyfH99yHH37oC7zHH/87+vXrh+zsEzhwYD+mTJlywe2tXr0G69atRU7OSahUqvPW1bFjR+zY8T0A4LfffsXevXsAAJMnP4CioiJs2/Yd9u3bC7fbg8cff7zaPj7//HN8881m7Nu3F7m5uRg8eAjGjh2HY8eO4v33V2Du3Ln45ptvfetv2vQ1unfvjuPHj2P69KcxZcoUlJWVAQCWLl2KRYsW4/XXX0dOzkls2LAemZmZAIBly97Ae++9hw8/XIljx45h4MCBGDFiJFwu10XfVyJ/MQgpqo0dOw5NmzaFxWLBTTfdhKSkJAwaNAgqlQpDhw7Fb7/th8vlgiiKWLVqNWbOnIGkpCQkJyfj2WefxcqVKyGKIr788gs0atQYY8eOgVqtxl//+ldce+21vv18+eVGWCwJmDRpErRaLTIzMzFx4kSsXr3mojWmpaXhuuuuw6pVqwEA+fn52Lx5M0aOrAhCjUaNkydP4vTp09DpdOjRo8cFt/fQQ1PQuHFj6HS6gOuyWq3YsGEDXn75ZSQkJMBoNOLpp5/C2rXr4PV6fes98cQTSEhIgMFgwMqVH6JHj+4YOnQIVCoV2rdvj9GjR2P16tW+9Tt27OhbPnJkRZAdOXIUAPDGG2/iySenoXPnTlAoFMjMzES7du0AAMuWLcNTTz2FVq1aQa1WY+LEiXA4HPjxxx8v+r4S+YuDZSiqpaWl+r43GPTVHuv1ekiSBLvdDpfLBZfLhaysLN/y5s2bw+l0oqCgAKdPn/H1UiplZmbC6XQAALKzs7F///5qrxdFCRkZGX7VOXLkSLz00kt45JGpWLNmDbp2vdq3vwULFuL//u8F9O59HRISEjBhwr2YMGHCebfVtGlT3/eB1nXiRDZEUUTHjpdXe16pVCI3N7fKPs6+F9nZ2di48b/V9uH1irjmmmt8j9PT03zfKxQK6PVxKC+v6BGePHkSrVq1qrWe7OxsTJhwH1Sqs/+zu1xunDp1qvbGE9UBg5AIQEpKCrRaLbKzs5GWVvFHOzs7GzqdDsnJyWjcuBFOnjxZ7TU5OTlITU0BAGRkZKBTp07YtOmrOu2/f/9bMHXqVPzvf3uwcuWHGD9+nG9Zy5Yt8Nprr0GSJOzYsQO33joYV111NTp37lTrtpTKs6ERaF1Nm2ZAqVTiwIEDMBgMNZafOHHiz30oqu1jwIABeOutN/3ax7kyMzNx7NgxXH311TWWZWRk4IUX5qJPnz512jaRP3holAgV4TFs2O2YPfs5FBYWobCwELNmzcaIESOgVCpx8819cfr0abz99jvweDz48ssvsWXLFt/r+/Xri7y8PLz++jI4HA54vV4cPnwYW7du9Wv/er0egwYNwnPPPYeDBw9i8ODBvmUffPAB8vLyoFAoYLFYoFQqq/WQLiTQutLT09G/f388/vjjKCgoAFAxsvTjjz8+7z5GjhyBLVu2YP369XC73XC73di3bx927/7JrxrHjLkHL7zwf9i3bx8kScLJkyd9A2LuvXc85sx5HocPHwYAlJaW4tNPP/WdXyQKBgYh0Z9eeOEFZGVloWvXrujatRtatmyJ55+fAwBISkrE+++vwJIlS5CV1Qzvvvsuhg0b5nutyWTC+vXr8e2336JDh8vRokULjBs3Hrm5eX7vf9Sokdi0aRP69+8Ps9nse/6bb75Bjx490aRJBkaNugPPPTcbl19++QW2dFZd6lq8eBEsFguuu+56NG2aiX79/oo9e85/oX2TJk2wdu1HeOutt9G2bTu0bt0Gjz32uN9hNXHiRIwbNxb33DMGGRlNceutg3HyZMXAoQkTJuCOO+7A3/52J5o2zcTVV3f167wrUSAUHo/7whcaERERRTH2CImIKKYxCImIKKYxCImIKKYxCImIKKYxCImIKKYxCImIKKZFXBBWTol1sdvLEBER+SPiglAQBNxz950QBKHe2yosyA9CReGD7QlvbE94Y3vCWyDtEQQ7XE6n3+tHXBAGkz/3ioskbE94Y3vCG9sT3vxtjyDYUFpSBIfT/85STAchERFFj4oQLIZeb4DZbPH7dQxCIiKKeIJgPxuC8QlQKBQXf9GfeBsmIiKKeBqNBgajCSZTfEAhCLBHSEREEczpdEAURajVGpjNloBDEGAQEhFRhBIEO4qLCiAItnptJyYPjeYUuvHlz2X4PVeAyZiPVmla3HyZCfF6ldylERGRHyrOCRZBrzfAYDDVa1sxFYQF5R68/LkVu47ZAVQMx1UqK24euuzbQgzsFI8J1ydBowq8a01ERA2jaggGOjCmNjEThAXlHkxZfhpnStwwxSnRr4MZzcwO6E0WfHfYjq0HbVi7uwSnS9yYOTgdaoYhEVFYEr3eoIUgIFMQlpWVYfasmb7HLpcTubm5WLbsTZjM5pDsc94XVpwpcaNjlh6zh6TBFKeCNT8XKakmXH+JCSesLjz+4Rl8f8SOj34swYiuCSGpg4iI6sbj8UCtVsNoMkOSpKCEICBTEJrNZrz08jzf4w0b1uO3334NWQjmFLqx86gdRp3SF4LnapaixdODUvHI+6ex/qdS3H6VBSplYG/yypUrsWTJa3C7XejduzfUag3atGmNO++8EwAwefJk9OvXD5deehkmTZoEm80GtVqNf/3rVVx++eVYsWIFNm78L4qLi3HixAmMGzcWDz74IMrLy3HXXXfj9OnTAIDnnnsOffrcWP83hogoTPzvhIBFmwrw6x9OuL0SNCoFLs3QYdKNyejcTA9RFFFgzUViUgq0Wl3QQhAIk0Ojm7/ehFF3jK51mdvthtvt9j2uyxyjG3+tOA/Yt4O51hCsdHnTOLRK0+Jongt7sx3o0lzv9z4OHjyITz/9DF999V+o1Wrcd999aNGiBf7zn/W488474fF48O2332LevHkQRRHr1/8HOp0Ov/zyC55+ejrWr/8PAODXX3/FN99shsfjwZVXXoX77rsPmzZtQlJSEtau/QiSJKGsrCzg94CIKBx5PCImvXcKu45V/9vu9kr44biAMctycEUzLWbfooTZaIJGow16DbIH4cGDB1BebsMVV1xZ6/J169ZizepVNZ4vsObBrvcvqE6cqbjOpJlZgDU/1/e81+up9hgAWiV5cfiMiCM5Bcgyavxux6effIJdO3eiZ8+eAABBcKB161Y4evQIjh45hL1796FTp04oLytBSUkJZs6ajQMHDkGlUqKwsBDW/FyUlZWiW7eucPw5FDg1JQUH9v+GJo3TsXXrVjz22KO4+aab0KVLZ1jza/5DUFt7IhnbE97YnvAWKe2ZusaGfae8UAK4pLESo67UoXWKCkesXnzwoxP7T4v48YQL09Yr8cptTriseX5vOyU13a/1ZA/CrzdtQu/evaFS1d5TGzJkKAYMGOh7LAgC7p84AckpaTAYDH7tw2jMh1JZBoMpASmpZ4fZVpwjrP5GqXUV6yYkWJCS6v+hWoPRiDFjx2DatGnVni8rs+H7Hbuwc+cujBgxAimp6Vj6+hto06Yd3n33PdhsNnTocDlSUtNhNsfDYknw1aTV6WBJSESzZs3w/ffb8cUXX+DFl+Zh+PBhmDBhQo0aamtPJGN7whvbE94ioT3//aUMP58qhVIBzBqSjkGd433LLmkJDLhKwsqtOXjxKwf2nRKxJ8+Imy4L/ik0WS+odwgCvv9+O66/4YbzrqPRaGAwGHxfej97gVW1TqvoSn936MIXXbq9EnYcsf/5Gl1A++jduzfWrl2HwsJCAEB+fj7OnDmDwYMHY82aj/DNN9/gpptuAlAxWKhx40ZQKBR4//33L7rt06dPw2Aw4I477sD999+Pfft+Dqg2IqJw9MaWIgBAzzaGaiEIwDcYZkTPDPRoY6i2frDJGoTbt29Ds2bNkZHRNKT7uekyE7RqBb47bMfxfNd51/t8XxkKbV60z4hDy7TAjkNfcskleOyxRzFo0CB0794dw4YNR1FREbp06Yxjx47h6quvRlxcHABg/PjxeOutt9CzZ09fcF7Ir7/+huuvvwE9e/bE4sWL8cADkwOqjYgoHB3Krbhn4JSbU6o9Lwh2FFhzIYpeKJVKPPTn8sr1g03WQ6Nff/01buzTJ+T7McepMLhLPFbtKsHfPzyNJwekoXOzON9yt1fCZ3vLsHBTAQDgjm7+376jquHDh2P48OE1nt+3b2+1x61atcL27dt9jysPp44eXX3A0LfffgMAaNasGUeJElFUcXlEiCKgVAJt0s8egat6sbxCUdFXa5Oug0oBiCLgcImI0wa3DydrEP5jzvMNtq9x1ybhVLEH3x2y4fEPT6NFqhatk7xQ6/Kx46gdRTYvAGDi9cm4prWxweoiIopFaiWgQEW4Fdk8SDSqzztjTJHNA1GqWF8bgtSSfbBMQ1GrFHj21jSs212KdbtLcTzfhaO5Z6dYuzQjDqO6WRiCREQNQKlUIt2ixpkSDxZ9XYBpt6Scd9q0hZsKIAFoZFFDqQz+Gb2YCUIAUCkVuP0qC4ZeEY+fcxw4fLIAiYkJaJmqRYvU4F+bQkRE5zfs6nj8+7+F2PBTGe7oloimyalQqzXVQvB4vgsf/6/Mt34oxORtmJRKBTpm6XFdWw1ubG9iCBIRyeDO7olIMSnh8Ei4Y3E2lm0tg9MtAag4F7j46wKMWpwNp0dCskGBO7snhqSOmOoREhFR+PC6Bbx6mxYPrnaiyC7htc2FeGNLEeI0CjjcEjzeilBMNqmwYJgeWnVo+m4x2SMkIiJ5CYINpSXFyEo14dNHWmBkNwssBiU8XgnlDhEer4QEgxJ3dLPg46nNkWoOXVyxR0hERA3K6XSgtKS42sCYaf3TMK1/Gk4WuJBf5kGqWY3M5LOnrewhrIdBSEREDUqr1cEcn/DntYLV7yKRmaytFoANgYdGiYioQQiCHS6XCwqFAgaDMai3UqoPBiEREYVc5cXyTkcoD3LWDYOQiIhCquqMMSZz3aawDCUGIRERhcz5pk0LJxwsQ0REIaNWa2AwmGAyx4dlCAIMQiIiCgGn0wGtVgeNRgONJvwOh1bFQ6NERBRUgmBDcVEBBPuFb4YeLhiEREQUNJUzxuj1BugNkXE3HwYhEREFRdUQDNeBMbVhEBIRUVB4Pd6IC0GAg2WIiKievF4PVCo1jCYzAERUCALsERIRUT0Igh3W/Fzf1GmRFoIAg5CIiOqo6sXyGo1G7nLqjEFIREQBi4QZY/zFICQiooBIkgSbrSwqQhDgYBkiIgqAJElQKBRISkqN2HOC52KPkIiI/CIINhQU5EEURSiVyqgIQYBBSEREfqi8WF6r0UZNAFZiEBIR0QVF6owx/mIQEhHReXm9nqgOQYCDZYiI6AJUKjUSk1Kh0WiiMgQB9giJiKgWgmBHeVkpJEmCVht95wWrYhASEVE1lRfLi6JX7lIaBIOQiIh8omnGGH8xCImICADgdDpiLgQBDpYhIqI/abU6mM0W6A3GmAlBgD1CIqKYJwh2uN0Vt1EyGE0xFYIAg5CIKKZVXCxfBIcgyF2KbBiEREQxquqMMSZzvNzlyIZBSEQUgypGh0b3jDH+4mAZIqIYpFKpYTAYYTJbYjoEAQYhEVFMcTod0Gp10Gq10Gq1cpcTFnholIgoRgiCHcVFBRAEu9ylhBUGIRFRDKg6Y4xeb5C7nLDCICQiinKxOG1aIBiERERRzuN2MQQvQLbBMm63G+++8zb27t0DjUaLZs2bY8qUh+Qqh4go6ni9XqhUKpjMFgBgCJ6HbEG4YsVyKBQK/Gv+AigUChQXFclVChFR1BFFEdb8M0hKToVGw9GhFyJLEDocDmz+ehMWL1nq+w8lITGx1nXdbjfcbrfvsRDD0wAREflDEGyQJBF6vQFqtUbucsKeLEGYm3sGJpMJ69Z+hJ9//hlarRbDhg9Hhw6X11h33bq1WLN6VY3nC6x5sOv19arD6/XAmp9br22EE7YnvLE94S1a2iOKIiRJBAA4nU64rHkyVxQcdfl8UlLT/VpPliD0ekXk5+ejaWYmRv/tThw/fgz/eG425r3yKhISEqqtO2TIUAwYMND3WBAE3D9xApJT0mAw1G8IsDU/1+83KhKwPeGN7Qlv0dAeSZJQYM2FVmuA0+lEalojuUsKmlB+PrIEYUpKChQKJXr17AUAaNGiJdLS0pCdfaJGEGo0Gmg07NoTEV2IJElQKBRISk6FQqGMmp5gQ5Dl8on4+Hh06NABe/buAQDk5eYiLy8PTTOaylEOEVFEEwQ7CgvyIYoilEoVR4cGSLZRo/dOmIAlixdhxfLlUCoVuHfCfUhKTparHCKiiFT1YnkGYN3IFoTp6Y0wY+ZsuXZPRBTxOGNMcHBmGSKiCOTxeBiCQcLbMBERRSC1Wo3EpBRoNFqGYD2xR0hEFEEEwQZbeRkAQKvVMQSDgEFIRBQhBMGG0pJieL0eSJIkdzlRg0FIRBQBKkOQ5wSDj0FIRBTmnA4HQzCEOFiGiCjMabRamMzxMBhMDMEQYI+QiChMCYIdHo8bSqUSRqOZIRgiDEIiojBUebG8INjlLiXqMQiJiMJM1RljTKZ4ucuJegxCIqIwwmnTGh4HyxARhRGVUgW9wQiz2cIQbCAMQiKiMOByOaHRaKHV6aDV6eQuJ6bw0CgRkcwEwYaiQiscHBgjCwYhEZGMqs4YE6c3yF1OTGIQEhHJhNOmhQcGIRGRTNwuF0MwDHCwDBFRA/N6vVCpVDDHJwAAQ1Bm7BESETUgQbDDas2F2+2GQqFgCIYBBiERUQPxXSwfp4dazQNy4YJBSETUADhjTPhiEBIRhZgkiSgvK2UIhin2zYmIQkiSJCgUSiQlp0KpVDIEwxB7hEREIVI5Y4wkilCpVAzBMMUgJCIKgcqL5dVqNcAADGsMQiKiIOOMMZGFQUhEFEQej5shGGE4WIaIKIjUag0SEpOh1eoYghGCPUIioiAQBDvstnIAgE4XxxCMIAxCIqJ6qrxY3uNxQ5IkucuhADEIiYjqgTPGRD4GIRFRHTkcAkMwCnCwDBFRHWm1OphM8TAYTQzBCMYeIRFRgByCHR6PB0qlEkaTmSEY4RiEREQBEAQbSkqK4BBscpdCQcIgJCLyU9UZY4ymeLnLoSBhEBIR+aFidChnjIlGHCxDROQHpVIJvcEIs9nCEIwyDEIiogtwuZzQaLTQ6eKg08XJXQ6FAA+NEhGdhyDYUVRohcMhyF0KhRCDkIioFlVnjImL08tdDoUQg5CI6BycNi22MAiJiKqQJAkup4MhGENkGywzedJEqNUaaLVaAMCQIUPRvUcPucohIoIoeqFUqhBvSQQAhmCMkHXU6NSpj6B5ixZylkBEBAAQRRHW/FwkJadCrdbIXQ41oLC/fMLtdsPtdvseCwJHbxFRcAmCDZIkQq83QKUK+z+LFGQKj8cty10kJ0+aCL1eD0kCWrdujdGj/4Z4i6XGeqtWfYg1q1fVeH7evHnQ6+s3ksvr9UTVDz3bE97YnvAkiiIkSQQAKJWqqDkcGi2fT6W6tCclNd2v9WQLQmt+PlJSU+HxeLBy5Qc4mX0CTz41vcZ6tfUI7584AW+/8x4MBkM9a8j1+42KBGxPeGN7wo8oiiiw5kKni4PT6URqWiO5SwqaaPh8qgple2QbNZqSmgoAUKvV6N9/APbv31/rehqNBgaDwfdV314gERFQMTpUqVQiKTmVo0NjnCxB6HA4YLOdvYXJtu+2ogUHzRBRAxEEO4qLCiBJElQqNUMwxslyALmkpBjzXn4ZouiFJAHp6emY/MAUOUohohhT9WJ5IkCmIExPb4QXX3pZjl0TUQzjjDFUG84sQ0Qxwe12MwSpVtEztpaI6AI0Gg0SEpKg1cUxBKkaBiERRTVBsAESoDcYoeNdJKgWPDRKRFFLEGwoLSmG2+2SuxQKYwxCIopKlSFYeU6Q6HwYhEQUdRwOoVoI8pwgXQjPERJR1NFqtTCazDAazQxBuij2CIkoajgcArzeinsKmkzxDEHyC4OQiKKCINhRUlwIwW67+MpEVTAIiSjiVZ0xxmgyy10ORRgGIRFFNE6bRvXFwTJEFNEUCgX0eiPM8RaGINUJg5CIIpLb7YJarUFcnB5xnDGG6oGHRoko4giCDYUF+XA6HXKXQlGAQUhEEaXqjDE6XZzc5VAUYBASUcQ4d9o0nhOkYGAQElFEkCQJToeDIUhBx8EyRBT2RFGEUqmEJSEJABiCFFTsERJRWBMEO6z5Z+DxeKBQKBiCFHQMQiIKW5UXy8fF6aFSqeQuh6IUg5CIwhJnjKGGwiAkorAjiiLKSksYgtQgOFiGiMKKJElQKpVISk6FSqViCFLIsUdIRGFDEGwoKS6EJElQq9UMQWoQDEIiCguVF8srlfyzRA2LP3FEJDvOGENyYhASkazcbhdDkGTFwTJEJCu1WgNLQhJ0ujiGIMki4B7hxx9/jOLi4hCUQkSxRBDscAh2KBQKxMXpGYIkm4CD8OWX56FNm7a49tremD79GWzcuBHl5eWhqI2IolTlxfIul0vuUogCD8Jvv/0Ghw4dwt///jhcLidmzJiJFi1a4uab+4aiPiKKMtVnjLHIXQ5R3c4RJiYmoG3btjh9+gzOnMnFqVOnIIpisGsjoijjcAicNo3CTsBBOG7ceGzbtg3JyUno3bs3Ro0aiQUL/o34+PhQ1EdEUUSj0cJoNMNoMjMEKWwEHISbN29GfHw8+vS5Cb169UL37tfAYDCEojYiihIOhwCtRguVSgWTmf80U3gJ+BzhsWNHsXz5e0hPT8OyZa+jQ4fL0bdvPzz//POhqI+IIlzltGl2wSZ3KUS1qtM5wssuuwzNmzdH69at0aJFC6xYsQI//PADnnrqqWDXR0QRrOqMMUajWe5yiGoVcBDOnDkL3333Hfbu3Ys2bVqjV69eWLhwIXr27BWK+ogoQnHaNIoUAQdhaWkpJk+ejF69eiIlJSUUNRFRNJDAEKSIEHAQvvLKPN/3BQUFSE5ODmpBRBTZ3G43NBoN9AYj9DDKXQ7RRQU8WEYQBDz88FQ0atQYrVu3QaNGjfHww1Nhs/FEOFGsEwQ7Cgvy4HQIcpdC5LeAg/Cpp57GkSOHsWHDehw8eAAff7wBR48exfTpz4SiPiKKEFVnjNHq4uQuh8hvAR8a/fzzz7F9+3YkJSUCANLS0vDOO2/jmmu645//fCXoBRJR+Ks+bRrPCVJkCbhHKEkSlMrqP+QKhRKSJAWtKCKKHJIkwSHYGYIUsQIOwr59++Kuu+7GTz/9D1arFbt3/4QxY8agX79+dSpg8+avMXzYbdi1a2edXk9E8hFFEQqFAgmJyQxBilgBB+Hzz89BZmZT9OvXD23atMUtt9yCjIwmmDPnHwHvPC8vD5u++gpt2rQN+LVEJC9RFFFgzYXX64FCoWAIUsQKOAhNJhMWLlyI3NwzOHToIM6cOY2FCxfCbA5s1ghRFPHakkUYO3YcNBrNeddzu92w2+2+L0HgaDQiuQmCDZIkQqeLg1Kpkrsconqp0xRrAKBQKJCamlrnHX/yycdo1+4vaNmq1QXXW7duLdasXlXj+QJrHux6fZ33DwBerwfW/Nx6bSOcsD3hLVraI4oiJKnitmtOpxMua57MFQVHtHw+ldgeICU13a/1/ArCrKxmfh32OHHid792mp2djZ07dmDW7Ocuuu6QIUMxYMBA32NBEHD/xAlITkmr910vrPm5fr9RkYDtCW/R0B5R9MKanwu93gCn04nUtEZylxQ00fD5VMX2+M+vIHz//RVB3emB/b8hPz8PD015AABQXFyMpa+dRHFREW7uW33QjUajueChUyJqGBUjxlVISk6FSqWOmp4gkV9BOGPGTGza9BUA4IUXXsC0adPqtdOb+/arFngzZzyLW/r3x9VXd63XdokoNATBDqfTAYslEWo1/zGl6OLXYJnDhw/D4/EAABYsWBjSgogovFReLK/kqFCKUn71CHv16omePXuhVatWEAQBo0f/rdb1VqxYXqciZs6aXafXEVFoccYYigV+BeGbb76J9evX48SJE9i4cSM6dLgs1HURkcxcLidDkGKCX0Go0+kwfPhwAEBxcUm9zxESUfjTaLSItyQiLk7PEKSoFvB1hHWZQYaIIocg2KFUKKCL00Ovr98lSkSRIOCZZYgoegmCDaUlRXC6nHKXQtRgGIREBKAyBIsrzgmaLXKXQ9RgGIREBIdgPxuCHBhDMSbgIHz99WW1Pv/QQw/XtxYikolGo4XRaGYIUkwKOAgXLlyI9evXV3vukUcexa+//hq0ooioYTgcAkRRhEqthskczxCkmBTwqNE1a1Zj0KBbkZycjJ49e+Lvf/87fvrpJ6xf/58QlEdEoVJ5sbzJFA+jKbDbqBFFk4CDsHXr1njvvXcxevTfcO21vXDgwEFs2LAeFgtPrhNFiqozxhiMJrnLIZKVX0H4yy+/VHus0+lw3333YcmSJVi27HXk5OQgJycHl13GGWeIwh2nTSOqzq8g7NmzFxQKBSRJqrFs4MBBACpu1FtUVBjc6ogo6CRJZAgSVeFXEBYXF4W6DiIKMY/HDbVaA4PBBEmSGIJEf+J1hEQxQBBsKLDmwel0AABDkKiKgAfLnDp1CnPmzMGePXtQVlZebdm+fXuDVhgRBUfVGWO0Wp3c5RCFnYCDcMKECdDrDXj44YdhMHBCXqJwVm3aNJ4TJKpVwEG4Z89eHDt2FFqtNhT1EFGQSJIEwW5nCBJdRMDnCP/yl78gNzc3FLUQUZBIkgiFQoHExGSGINFFBNwjHDhwIEaNGoXx4+9FWlpqtWW33HJL0AojoroRBDvKy0qRlJwKlUoldzlEYS/gIFy2rGLS7Xnz5lV7XqFQMAiJZFb1YnmlkoPCifwRcBD+/PO+UNRBRPXEGWOI6ob/MhJFAdHrRWkpR4cS1YVfPcK+ffvhyy+/AHB2urXabN26JXiVEZHflCoVkpNSoVKrGYJEAfIrCMePH+f7ftKk+0NWDBEFRhBscLlciI9PgFqjkbscoojkVxAOGzbM9/0dd9wRsmKIyH9VL5YnorrzKwg/++wzvzbGUaNEDYMzxhAFj19B+MQT0y66Di+fIGoYLpeTIUgURH4FIS+ZIAofGo0W8fEJiNMbGIJEQcDLJ4gihCDY4XI6oVAooDcYGYJEQcIgJIoAlRfLO5yC3KUQRR0GIVGYqzZjjNkidzlEUYdBSBTGHJw2jSjkAp5rlIgajlqjgcFogskUzxAkChG/gvBC06pVxSnWiILD6XBAo9VCrdbwcChRiPkVhJxWjajhVF4sbzLHw2g0y10OUdTzKwg5rRpRcBzJdeLfXxXg1z8ccHkk6LVK9GhtwOQ+yUg2qavNGGMwmOQulygm1Okc4fLly7Fq1WpYrfnYvn07vvvuO+Tl5WPo0CHBro8oKng8Ip5cb8MPJ0qrPV8qiFi7uxT/+V8pRl5lxN1XShwYQ9TAAh41+tJLL2PRokW47bbbkJOTAwBo1KgR5s+fH/TiiKLF+Lf/wK4TXigAtEnX4umBqVh8VxM8cGMSGieoIYrA+zttWLpdZAgSNbCAg/Ddd9/F6tWrcffddwGo+GVt2bIljh8/HuzaiKLCmh+KseeEA0oFMGtoOlY/0AzDrk7ANW2MGH9dMj5+KBMju1UMiFn7PweOW10yV0wUWwIOQrvdjkaNGgGA779Wt9sNnU4X3MqIosS724oBANe1UWNQ5/hqywTBjgJrLh65yYJWaVpIABb8t7DhiySKYQEH4VVXXYlly5ZVe+6995aja9euQSuKKFqUCh5kF7ihUADju8dVW1Z1xhiNRosxvRIBALuO2eUolShmBTxYZu7cFzBo0CCsWPE+bDYbbrrpZuTl5WH9+v+EoDyiyHaq2AMAiNMokB5/9v9OoZYZY3q2rbjBrsMjyVIrUawKOAhbtGiOXbt24ssvNyI7OxsZGRno168vjEZjKOojimgmXUX4ebyAKIoAAEmSYLeV1xgdml/mBQAoOU6GqEHV6fIJvV6PwYNvrdeO//HcbBQXF0GhUEKv12PM2LFo0aJlvbZJFG6aJKih1yoguCRs+NmNcWkSFAoFEpNSoFAoqo0Off2bQt9riKjh+PUbN3nyZL82tnDhQr93PPWRR329yF07d2LRwgV46eVX/H49USRQKpW44RITPt1bhuW7nOh3yWk0aZwOpVJVbb1j+U5s3m8DANzZI0GGSolil1+DZeLj431fCoUSq1evQW5uHrRaHfLy8rFmzUc1frEvpuqhVLvdhspLMc7ldrtht9t9X4LA+7FRZHno5hTo1AoUCcDdywVs3n92MIzHI+KDHcUYveQk3F4JafFqDOkSf4GtEVGwKTwed0Bn5u+88y6MHn0H+vXr53vuyy+/xHvvLcfy5e8FtPMF/56PX3/9BQDw5JNPI6tZsxrrrFr1IdasXlXj+Xnz5kGv1we0v3N5vR6oVNFzGIrtCU+iKGLfHy48tcEFZ8VpQJi0gEGrQKlDgqNiPA0S9AosGWlEqjky7o4WLZ9PJbYnvNWlPSmp6X6tF3AQNm2aiezsE1Aqz/6yer1eNGvWHDk5JwMqstI332zG99u34cmnptdY5na74Xa7fY8FQcD9Eyfg7Xfeg8FgqNP+Klnzc/1+oyIB2xN+vF4vrPlnoNcbcCxXwLKdEnYeFeD2nv21M2gVuLG9CU/0T4UpLrAjK3KKhs+nKrYnvIWyPQH/u5CVlYn33lv+58wyFVasWIHMzMw6F3Hdddfj9aVLUVZWBrO5+mz7Go0GGo2mztsmkpNKpUJScirUag2aulxYcGc67C4RO47YUCqISI9Xo2srfbV/LImoYQUchC+99BJGjboDixcvRmZmJk6ePIlTp07hgw/e93sbNpsNTqcTSUlJAIBdu3bCbDbBZOJs+xQdBMEOj9sFk9kCjUZbbZlBq8QN7Xl7JaJwEXAQ9ujRA/v27cUXX3yBM2dy0bhxI9x8c18kJib4vQ273YZX5s2Dy+WCUqlAfHw8pk17ihMNU1SoerE8EYW/Op1JTUhIwMiRI1FQUIDk5OSAX5+amoa5L/xfXXZNFNZqmzGGiMJbnSbdfuihh9GoUWO0bt0GjRo1xsMPT4XNZgtFfUQRw+V0MgSJIlDAQfj009Nx9OgRbNiwHgcPHsDHH2/A0aNHMX36M6GojyhiaLRamOMTGIJEESbgQ6Off/45tm/fjqSkipny09LS8M47b+Oaa7rjn//kzDAUewTBDpVKBa1WB4OBc+4SRZqAe4SSJEF5zqzACoUSksQZ8yn2CIINpSVFcDo44xFRpAo4CPv27Yu77robP/30P1itVuze/RPGjBlTbaYZolhQEYLF0OsNMJktcpdDRHUUcBA+//wcZGY2Rb9+/dCmTVvccsstyMhogjlz/hGK+ojCUsXo0GIOjCGKAgGfIzSZTFi4cCEWLFgAq9WKlJQU/hGgmKNWa2AwmGAyx/PnnyjC1XlGVq/XC51Oh7KyMt9z8fGcNZ+im9PpgFar+3PqPx4OJYoGAQfhDz/8gIcffhj79x/wDZCRpIqbjRYVFQa9QKJwUXmxvNlsgcHI6QCJokXAQThx4v24/fbb8Oabb9b7NkhEkaLqjDF6XiJBFFUCDsL8/HxMmzaN50UoZnDaNKLoFvCo0WHDhuGzzz4LRS1EYcnr8TAEiaJYwD3C6dOno0+fPvjXv+YjNTW12rIVK5YHrTAiuXk9HqjUahhNFbdMYggSRaeAg3DChAnQarXo1q0bDAaeI6ToVHmxfGJSCrRandzlEFEIBRyE27Ztw8GDB2rcSZ4oWlSdMebcm+oSUfQJ+Bxhu3btUF5eHopaiGRXNQR5TpAoNgTcIxw4cCCGDx+BcePGIS2t+jnCW265JWiFETU0SZJgs5UzBIliTMBB+NZbbwEA5s2bV+15hULBIKSIVTkpRFJSKhQKBUOQKIYEHIQ//7wvFHUQyUYQbLDZypGUlAqlMuCzBUQU4fhbTzGt8pygVqNlL5AoRjEIKWZxYAwRAQxCilFer4chSEQA6nEbJqJIplKpkZiUCo1GwxAkinHsEVJMEQQbystKIEkStFqeFyQiBiHFkMpzgqIoyl0KEYURBiHFBA6MIaLzYRBS1HM6HQxBIjovDpahqKfV6mA2W6A3GBmCRFQDe4QUtQTBDrfbBYVCAYPRxBAkoloxCCkqVZwTLIJDEOQuhYjCHIOQok7VgTEmc7zc5RBRmGMQUlQRBDsHxhBRQDhYhqKKSqWGwWCEyWxhCBKRXxiEFBWcTge0Wh20Wi20Wq3c5RBRBOGhUYp4gmBDcVEBBMEudylEFIEYhBTRqg6M0esNcpdDRBGIQUgRi9OmEVEwMAgpYnncboYgEdUbB8tQxPF6vVCpVDCZLQDAECSiemGPkCKKINhgzT/jmzqNIUhE9cUgpIhR9ZygWq2RuxwiihIMQooIHBhDRKHCIKSwJ0kSbOVlDEEiCglZBsu4XC68+uor+CMnB1qtFvHxFtx77wQ0atxYjnIojEmSBIVCgaTkVCgUSoYgEQWdbKNG+/S5CZ07d4FCocAXn3+GJUsWY+as2XKVQ2FIFEUUFuQjMSkFSqVK7nKIKErJcmhUq9WiS5crfP/dt2nbFvn5ebWu63a7YbfbfV8C7y8XEwTBBkkSodFo2AskopAKi+sIP/v0U1x55VW1Llu3bi3WrF5V4/kCax7sen299uv1emDNz63XNsJJtLRHFEVIkggAcDqdcFlr/ycp0kTL51OJ7QlvbA+Qkpru13oKj8ct1aWoYFm79iPs/vFHPDtjJnQ6XY3lbrcbbrfb91gQBNw/cQLefuc9GAz1m1vSmp/r9xsVCaKhPR6PBwXWXOj1BjidTqSmNZK7pKCJhs+nKrYnvLE9/pO1R7hhw3rs2rkTzzw7o9YQBACNRgONhteMxQq1Wo3EpBRoNNqo6QkSUXiT7fKJTz7egG3ffYfpzzwLo9EoVxkUJgTBBlt5GQBAq9XxvCARNRhZeoQFBQV49913kJ6ejlkzZwCo6Pk9P/cFOcohmVW9WL7ycgkiooYiSxAmJydj1eqP5Ng1hRnOGENEcuPMMiQbp8PBECQi2YXF5RMUmzRaLUzmeBgMJoYgEcmGPUJqcIJgh8fjhlKphNFoZggSkawYhNSgKs4JFkEQ7HKXQkQEgEFIDajqwBiTKV7ucoiIADAIqYEIgp0DY4goLHGwDDUIlVIFvcEIs9nCECSisMIgpJByuZzQaLTQ6nTQnmcaPSIiOfHQKIWMINhQVGiFgwNjiCiMMQgpJKoOjInT1+8uIUREocQgpKDjtGlEFEkYhBR0bpeLIUhEEYODZShovF4vVCoVzPEJAMAQJKKIwB4hBYUg2GC15sLtdkOhUDAEiShiMAip3nznBOP0UKt5kIGIIguDkOqFA2OIKNIxCKnOJElEeVkZQ5CIIhqPY1GdSJIEhUKJpORUKJVKhiARRSz2CClglTPGSKIIlUrFECSiiMYgpIBUnhNUq9UAA5CIogCDkPzGgTFEFI0YhOQXj8fNECSiqMTBMuQXtVqDhMRkaLU6hiARRRX2COmCBMEGu60cAKDTxTEEiSjqMAjpvCrPCXo8bkiSJHc5REQhwSCkWnFgDBHFCgYh1eBwCAxBIooZHCxDNWi1OphM8TAYTQxBIop67BGSj0Oww+PxQKlUwmgyMwSJKCYwCAlAxTnBkpIiOASb3KUQETUoBiFVGxhjNMXLXQ4RUYNiEMY4QbBzYAwRxTQOlolxSqUSeoMRZrOFIUhEMYlBGKNcLic0Gi10ujjodHFyl0NEJBseGo1BlfcTdDgEuUshIpIdgzDGVB0YExenl7scIiLZMQhjCKdNIyKqiUEYIyRJgsvpZAgSEZ2Dg2VigCh6oVSqEG9JBACGIBFRFewRRjlBsMGanwuPxw2FQsEQJCI6B4MwilWeE4yL00OlYuefiKg2DMIoxYExRET+YRBGIVEUUV5WyhAkIvKDbMfL3nzzDez+8Qfk5+fjxRdfRvMWLeQqJapIkgSlUomk5FQolSqGIBHRRcjWI+zWrRtmPzcHqampcpUQdURRRHFRASRJgkqlZggSEflBth5h+/aX+rWe2+2G2+32PRYETgtWG0GwQZJEqFQquUshIoooYT+UcN26tVizelWN5wusebDr6zdFmNfrgTU/t17bCAeiKEKSRACA0+mEy5onc0XBES2fTyW2J7yxPeGtLu1JSU33a72wD8IhQ4ZiwICBvseCIOD+iROQnJIGg8FQr21b83P9fqPCldvtRmFBHvR6A5xOJ1LTGsldUtBEw+dTFdsT3tie8BbK9oR9EGo0Gmg0GrnLCFsajQYJCUnQ6uKipidIRNSQwj4IqXaCYAMkQG8wQse7SBAR1ZlsQbj0tSX46afdKC4uxpw5zyEuTo9/L1goVzkRperF8noY5S6HiCiiyRaEE+6bKNeuI9q5M8YQEVH9cGaZCOJwCJw2jYgoyHiOMIJotVoYTWYYjWaGIBFRkLBHGAEcDgFeb8U9BU2meIYgEVEQMQjDnCDYUFJcCMFuk7sUIqKoxCAMY1UHxhhNZrnLISKKSgzCMCUIdg6MISJqABwsE6YUCgX0eiPM8RaGIBFRCDEIw4zb7YJarUFcnB5xnDGGiCjkeGg0jAiCDYUF+XA6HXKXQkQUMxiEYaLqwBidLk7ucoiIYgaDMAycO20azwkSETUcBqHMJEmC0+FgCBIRyYSDZWQkiiKUSiUsCUkAwBAkIpIBe4QyEQQbrPln4PF4oFAoGIJERDJhEMqg8pxgXJweKpVK7nKIiGIag7CBcWAMEVF4YRA2IFEUUVZayhAkIgojHCzTQCRJglKpRFJyKlQqFUOQiChMsEfYAATBhuLiQkiSBLVazRAkIgojDMIQqzwnqFLyrSYiCkf86xxCHBhDRBT+GIQh4na7GIJERBGAg2VCRK3WwJKQBJ0ujiFIRBTGGIRBJgg2KKBAnN7A+wkSEUUAHhoNospzgi6XS+5SiIjITwzCIKk+MMYidzlEROQnBmEQOBwCB8YQEUUoniMMAo1GC6PRDKPJzBAkIoow7BHWg8MhQPR6oVKpYDLHMwSJiCIQg7COBMGGkuJC2AWb3KUQEVE9MAjroOrAGKPRLHc5RERUDwzCAHHaNCKi6MIgrAOGIBFR9OCoUT+53S5oNFro9Ubo9Ua5yyEioiBhj9APgmBDYUE+HA5B7lKIiCjIGIQXUfWcoE4XJ3c5REQUZAzCC+DAGCKi6McgPA9JkuAQBIYgEVGU42CZWoiiCKVSiYTEZABgCBIRRTH2CM8hCDYUWHPh9XqgUCgYgkREUY5BWEXlOUGdLg5KpUrucoiIqAEwCP/EgTFERLFJtnOEp0+fwsIFC1BWVgqDwYBJkx9AZmaWLLWIohdlpSUMQSKiGCRbj3Dpa6+hT58++Nf8Bbj11iFYtHCBLHVIkgSlUoWk5FSGIBFRDJIlCEtKSnDs2FH0urY3AKBrt26wWgtw5vTpBq1DFEWUlBRBkiSo1RqGIBFRDJLl0GiB1YqEhESoVBUDUhQKBVJSUmC1WtGoceNq67rdbrjdbt9jQQjONGeCYIMkiVAy/IiIYlrYX0e4bt1arFm9qsbzBdY82PX6Om1TFEVIkggAcDqdcFnz6lVjuPB6PbDm58pdRtCwPeGN7QlvbA+Qkpru13qyBGFySgqKi4vg9XqhUqkgSRKsVitSUlJqrDtkyFAMGDDQ91gQBNw/cQKSU9JgMBgC3rfL5URRoRV6vQFOpxOpaY3q1ZZwYs3P9fuDjwRsT3hje8Ib2+M/Wc4RWiwWtGjRElu3fAsA2LljB5KTk2scFgUAjUYDg8Hg+9LXsRd4dntaxFsSOTCGiIgAyHhodMKE+7Bw4QKsW7cWer0BkyZNDun+BMEOpUIBXZween3gPUkiIopOsgVhk4wMzHl+boPsy3exvMEIXVz9epRERBRdon5mmWozxpgtcpdDRERhJqqD0CHYOW0aERFdUNhfPlEfGo0WRqMZRpOZIUhERLWKyh6hwyFAFEWo1GqYzPEMQSIiOq+oC0JBsKGkuBCC3SZ3KUREFAGiKgirDowxGE1yl0NERBEgaoKQ9xMkIqK6iJoglCSJIUhERAGL+FGjHo8barUGBoMJkiQxBImIKCAR3SMUBDsKrHlwOh0AwBAkIqKARWyP0OEQ4HY5oNcboNXq5C6HiIgiVMT2CEtLinhOkIiI6i3ieoSSJP35nQIqtbZed6wXBAF2uz04hYUBtie8sT3hje0Jb3Vtj16vv2hnKeKC0OGoOB/497//XeZKiIgo3L39znsXvYm7wuNxSxdcI8yIooiioiLExcXV65Bo5Z3uFy9ZWu+b/YYDtie8sT3hje0Jb/VpT1T2CJVKJZKTk4O2Pb1ef9H/FiIJ2xPe2J7wxvaEt1C1J2IHyxAREQUDg5CIiGJazAahRqPB7cOGQ6PRyF1KULA94Y3tCW9sT3gLdXsibrAMERFRMMVsj5CIiAhgEBIRUYxjEBIRUUyLuOsIg+X06VNYuGAByspKYTAYMGnyA8jMzJK7rDp58803sPvHH5Cfn48XX3wZzVu0kLukenG5XHj11VfwR04OtFot4uMtuPfeCWjUuLHcpdXJP56bjeLiIigUSuj1eowZOxYtWrSUu6x627z5ayxetBCPPf53XH11V7nLqbPJkyZCrdZAq9UCAIYMGYruPXrIXFXdud1uvPvO29i7dw80Gi2aNW+OKVMekrusOikrK8PsWTN9j10uJ3Jzc7Fs2Zswmc1B20/MBuHS115Dnz59cN31N2DH999j0cIFmPvCi3KXVSfdunXDrbcOxrPPPC13KUHTp89N6Ny5CxQKBb74/DMsWbIYM2fNlrusOpn6yKMwGo0AgF07d2LRwgV46eVXZK6qfvLy8rDpq6/Qpk1buUsJiqlTH4n4fyArrVixHAqFAv+avwAKhQLFRUVyl1RnZrMZL708z/d4w4b1+O23X4MagkCMHhotKSnBsWNH0eva3gCArt26wWotwJnTp2WurG7at780qLPtyE2r1aJLlyt80yK1adsW+fl5MldVd5UhCAB2uw1AZN8tRRRFvLZkEcaOHRc1w/OjhcPhwOavN2HkqDt8vz8JiYkyVxU8m7/ehBtuuDHo243JHmGB1YqEhESoVCoAFTf0TUlJgdVqjdjDb9Hss08/xZVXXiV3GfWy4N/z8euvvwAAnnwysnvun3zyMdq1+wtatmoldylBs2DBfEgS0Lp1a4we/TfEWyxyl1QnublnYDKZsG7tR/j555+h1WoxbPhwdOhwudyl1dvBgwdQXm7DFVdcGfRtx2SPkCLH2rUf4cyZM7hj9N/kLqVeHnhwChYvWYoRI0dhxYr35C6nzrKzs7Fzxw4Mve12uUsJmlmznsPL8/6J/3vxJZjj47Fw4b/lLqnOvF4R+fn5aJqZiRf+70WMGTsWr/7zFRQXF8tdWr19vWkTevfu7evABFNMBmFySgqKi4vg9XoBVNzj0Gq1IiUlRebKqKoNG9Zj186deOrp6dDpdHKXExTXXXc9fvnlV5SVlcldSp0c2P8b8vPz8NCUBzB50kQcPnwIS19bgo1ffiF3aXWWkpoKAFCr1ejffwD2798vc0V1l5KSAoVCiV49ewEAWrRoibS0NGRnn5C5svpxCAK+/347rr/hhpBsPyYPjVosFrRo0RJbt3yL666/ATt37EBycjIPi4aRTz7egG3ffYdnnp1R7RxbpLHZbHA6nUhKSgIA7Nq1E2azCSaTSebK6ubmvv1wc99+vsczZzyLW/r3j9hRow6HA16v1/cztu27rWgRwYNm4uPj0aFDB+zZuwddulyBvNxc5OXloWlGU7lLq5ft27ehWbPmyAhRO2J2irVTf/yBhQsXoLy8DHq9AZMmTUZWs2Zyl1UnS19bgp9+2o3i4mKYzWbExenx7wUL5S6rzgoKCnD/xAlIT09HXFzFvcc0Gg2en/uCzJUFLj8/D6/MmweXywWlUoH4+HjceefdUTNCMdKDMDf3DOa9/DJE0QtJAtLT03HPmLFIS0uTu7Q6y809gyWLF6G0tAxKpQK33T4M3bpdI3dZ9TL96adwY58+uP760PQIYzYIiYiIgBg9R0hERFSJQUhERDGNQUhERDGNQUhERDGNQUhERDGNQUhERDGNQUhERDGNQUh0jq5du+GLLyqmDFuxYgV69uzp1zJ/rFq1CjfddHPwivXDmjUf4e677/E9tlgSsG/fvgatIVAPPjgFb7/9zgXX6dChAz755BMAgb+vgbwHDocDo0ePRlZWFq6//gZs3boVWVmRee9Sql1MTrFG59e/f3/s2vVDtdvr6HQ6HD9+TMaqGtbOnTvqtMwfw4cPx/Dhw+u1jUCIoojZs2fjgw/eb7B9nuv++++HxWLBCy/4PzPQY489iptv7otRo0b6Nc9sKN/X9evX4/DhIzh8+DB0Oh22bt0akv2QfNgjpBpmzZqJU6f+8H3FUgiGktvtbvB9bty4EYmJibj00kvr9Ho5agaAZs2aoXXrVli/fr0s+6/qxIkTaN26VdRM/E41MQjJb0888QQmT55c7bl//vNVDBkyFACwadPX6N37OmRmZqFt23Z45JFHIQiCb90OHTrg1Vf/hRtv7IOMjKa45ZZbkJOTA6DiDiDPPjsDbdq0RdOmmejS5QrfIcjarFnzEbp3746srCxcd9312Llzp29ZcXEx7r77HmRlZeHKK6/CkiVLYLEkVKuj8pAaAHzyySfo0KHDeZdXVduyWbNmo3nzFrj00svw+uvLfM/PnTsXw4ePwNSpj6BZs+aYOXNmjcOp5x6iW7RoEfr3719t+dKlS3H11V3RuHETTJgwAUVFxbjnnjFo2jQTPXv2wqFDh877Pn322ee49tprazy/bds2XHHFlcjKysI994xBSUkJgIo/+hZLApYvX45OnTrjkkvaAwD27NmDAQMGoFmz5ujUqXO1w5Zz587FiBEj8NhjjyMrKwuXXnoZPvpoLQBgyZIlWLVqNZYtewNNmmSga9duAC7+GQFA79698dlnn5+3bVWd+75e6GftXEeOHEHHjp2wdOnSGsuefvppvPjiS/jiiy/RpEkGnn/++RrrlJWVYcqUh9C2bTu0bdsODz88FTabDUDFne9nzJgJoOJnvGXLVhgzZqzvtdde2zsswj7WMQjJbyNGjMSGDR9XC7cPP/wQo0aNBADo9XGYP/9f+P334/jyyy+xdetWLDhn8u9Vqz7EG28sw9GjR2AwGDFnzhwAwNdfb8aaNWuwZcu3yMk5ifXr/4PWrVvXWsfGjRvxzDPPYNGiRfj999/xyCNTMWLESBQWFgKoCOySkhLs27cPH3+8AR98sDIUbwcA4Lff9kOhUODQoYN46603MWvWLGzbts23/KuvvsKVV16Bo0eP4Omn63ZD3k8//QxffPEFfvppN77+ejP6978F9903Ab//fhyXX94Bzz777Hlf+/PPP6Nt2zY1nl+58kN88snH2LdvH4qLi/Hkk09WW/7555/jm282Y9++vcjNzcXgwUMwduw4HDt2FO+/vwJz587FN99861t/06av0b17dxw/fhzTpz+NKVOmoKysDBMnTsTw4cMwfvw4nDr1h+/Qsj+fUbt27fDzzz/X6T0Dzv+zVtXu3bsxcOAgzJgxAxMmTKixfM6cOXj00UfQr19fnDr1B5566qka60ybNg3Hjh3Djh3f4/vvt+Pw4UN48smK9Xr16uU7lPrLL7/AbDb7fj6Kiorx66+/olevmv+oUMNiEFINs2bNRlZWlu/r1lsHAwC6dOmMJk2a4NNPPwMA7N27Fzk5ORgwYAAAoHv37ujYsSNUKhVatGiOMWPuwXfffVdt2+PGjUfz5s0RFxeH4cOHYc+evQAAjUYNh8OB/fsPwO12IzMz87xB+PrryzBlyoPo1KkTlEolBg0ahLZt22Ljxo3wer1Yu3YdnnlmOhISEtC4cWNMmTIlRO8UYDQa8eST06DVanH11Vdj2LBh1f6ot29/CUaPHg21Wg2DwVCnfTz44INISkpE48aN0aNHD/zlL5fgmmuugVqtxuDBg7F37/kHfVTckSS+xvMPPfQQGjdujISEBEyf/jRWr14DURR9y5944gkkJCTAYDBg5coP0aNHdwwdOgQqlQrt27fH6NGjsXr1at/6HTt29C0fOXIkXC4Xjhw5WmtN/n5GZnN8vW4oe76ftUpffbUJd9wxGkuWLMbQoUPqtA9RFLFq1WrMnDkDSUlJSE5OxrPPPouVK1dCFEX06tUTe/fuRWlpKbZs2YJBgwYhOTkZBw4cwHffbcUll1yCpKTEOreRgoODZaiGGTOexaRJk2pdNmLECKxcuRK3334bVq5ciUGDBvr+wO/e/RNmzZqF3377DQ6HAx6Pp0aYpaefvb2NwWBAeXnFDWqvvfZaPPnkk5gzZw4OHTqE3r174x//eA7NmzevUUN2djZmz34Oc+fO9T3ndntw6tRpFBQUwOVyITMz07csKyuzxjaCpVGjRtUGFmVmZlbrETZtWv/7p6Wlpfq+Nxj0sFgsvsd6vR42W/l5X5uQkICystIaz1d9TzIzM+FyuWC1Wn3PNW16dnl2djY2bvxvtZGSXq+Ia645e2ufqp+rQqGAXh/n+2zP5e9nVFZWioSEhPO27WLO97NWafHixejduzd69+5d531YrVa4XK5q703z5s3hdDpRUFCA1NRUtGnTBt9//z22bNmC8ePHw+12YevWrX/+nLM3GA7YI6SADB8+DN9++y1OnTqFNWs+wogRI3zLxo0bh169emHv3j3IyTmJZ599BoD/d/m6997x2LTpK/zyy8/Q6XR44oknal0vIyMD//jHP5Cdne37On36FB55ZCqSk5Oh0Whw8uRJ3/onT1Y/N2Q0mqod3j1zJtfvGs915syZagNKcnJy0LjKDZ6Vygv/ihmNxqDVUpsOHTrg0KHDNZ7Pzj77/uTk5ECr1SIlJcX3nFKp8H2fkZGBAQMGVHu///gjB2vWrIY/zn0P/PmMAODgwYPVzt0G27Jlr+PgwYN4/PHH67yNlJQUaLVaZGdn+57Lzs6GTqdDcnIyAKBXr57YvHkzdu36Addccw2uvfZabNmyFVu2bEWvXr3q3Q6qPwYhBaRp06bo1q0bHnjgQWi12moDMcrKypCQYIHRaMTBgwfxxhtv+r3d3bt/ws6dO+FyuaDX62EwGKBW137A4t57x2P+/Pn43//2QJIk2O12bN78Df744w+oVCoMGTIEc+Y8j+LiYpw+fRrz58+v9vqOHS/HmjVr4HA4cPz471i2bFmt+/GHzWbDiy++CJfLhR9//BGrV6/G8OHD/H59x46XY+XKD+HxeLBv3z58+OGHda6lNn/9a79ah/vPnz8fp0+fRnFxMebMeR633Tb0vKE9cuQIbNmyBevXr4fb7Ybb7ca+ffuwe/dPftWQmpqG33//HZJU8U+RP58RAGzZsgX9+vUNoLWBSUxMxIYNG7Br1w+YOvURX32BUCqVGDbsdsye/RwKC4tQWFiIWbNmY8SIEb73s1evXlix4n20bt0aJpMJPXr0wLfffosjR46ge/fuwW4W1QGDkGqYMWMmmjTJqPZVORAFAEaNGolNmzZh+PDhUCjO9hxeffWfmD//32jSJANTp07FbbcN9XufZWVlePTRR9GiRUu0adMWZ86cOe91Z3/9618xc+ZMTJkyBc2aNcPll3fEkiWLfee4XnzxRRiNRnTocDkGDhyEkSNHVHv99OnTUVJSilatWmP8+PEYOXJkIG9PNe3bXwKPx4u2bdvhrrvuxjPPTK91lOb5vPjii/jhh13IymqGGTNmYtSoUXWupTY333wzCgoK8Ntvv1V7fsSI4RgwYCA6dLgcJpPpgtf4NWnSBGvXfoS33nobbdu2Q+vWbfDYY4+jrKz2Q5/nuuuuu3Dq1Gk0a9bc94f/Yp9RdnY2Dh06jMGDBwfW4AAlJiZg/fr/YO/evZgy5aE6heELL7yArKwsdO3aFV27dkPLli3x/PNnB+b07NkLZWVlvp8Li8WC1q1bo1OnToiPr3n+lhoe71BPUW/fvn3o1etalJQUy12KLNasWYNPPvkUb7/9ltylnNe5n9GUKQ+hS5cuuOeeu+UtjGICg5CiXqwHYSTgZ0Ry4qFRIiKKaewREhFRTGOPkIiIYhqDkIiIYhqDkIiIYhqDkIiIYhqDkIiIYhqDkIiIYhqDkIiIYhqDkIiIYtr/AwbliL0BrdXhAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 480x460 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "display(viz.plot_network_flows(net, final.link_flows))\n",
    "display(viz.plot_flow_scatter((\"Evans equilibrium (brentq)\", ref_flows), {\"evans\": final.link_flows}))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e3596374",
   "metadata": {},
   "source": [
    "## Takeaways & pointers\n",
    "\n",
    "- **Certified, not self-reported.** The joint gap came from `Evaluator`; the\n",
    "  self-reports (assignment + distribution) were only diffed against it.\n",
    "- **Endogenous demand.** The OD matrix is an output, distributed by gravity at the\n",
    "  equilibrium costs — the split recomputed here via `brentq`.\n",
    "- **Where next.** Fixed-demand UE: [`bfw`](05-bfw.ipynb); elastic (total) demand:\n",
    "  [`fw-elastic`](13-fw-elastic.ipynb); ADR-007 in the\n",
    "  [model compendium](../../docs/MODELS.md)."
   ]
  }
 ],
 "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.10.12"
  },
  "tabench": {
   "covers": [],
   "requires_extra": null,
   "track": "static",
   "unit": "evans"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}