tutorials/composite_mtbo.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Composite Bayesian Optimization with Multi-Task Gaussian Processes\n",
    "\n",
    "In this tutorial, we'll be describing how to perform multi-task Bayesian optimization over composite functions. In these types of problems, there are several related outputs, and an overall easy to evaluate objective function that we wish to maximize.\n",
    "\n",
    "**Multi-task Bayesian Optimization** was first proposed by [Swersky et al, NeurIPS, '13](https://papers.neurips.cc/paper/2013/hash/f33ba15effa5c10e873bf3842afb46a6-Abstract.html) in the context of fast hyper-parameter tuning for neural network models; however, we demonstrate a more advanced use-case of **[composite Bayesian optimization](https://proceedings.mlr.press/v97/astudillo19a.html)** where the overall function that we wish to optimize is a cheap-to-evaluate (and known) function of the outputs. In general, we expect that using more information about the function should yield improved performance when attempting to optimize it, particularly if the metric function itself is quickly varying.\n",
    "\n",
    "See [the composite BO tutorial w/ HOGP](https://github.com/pytorch/botorch/blob/main/tutorials/composite_bo_with_hogp.ipynb) for a more technical introduction. In general, we suggest using MTGPs for unstructured task outputs and the HOGP for matrix / tensor structured outputs.\n",
    "\n",
    "\n",
    "We will use a Multi-Task Gaussian process ([MTGP](https://papers.nips.cc/paper/2007/hash/66368270ffd51418ec58bd793f2d9b1b-Abstract.html)) with an ICM kernel to model all of the outputs in this problem. MTGPs can be easily accessed in Botorch via the `botorch.models.KroneckerMultiTaskGP` model class (for the \"block design\" case of fully observed outputs at all inputs). Given $T$ tasks (outputs) and $n$ data points, they assume that the responses, $Y \\sim \\mathbb{R}^{n \\times T},$ are distributed as $\\text{vec}(Y) \\sim \\mathcal{N}(f, D)$ and $f \\sim \\mathcal{GP}(\\mu_{\\theta}, K_{XX} \\otimes K_{T}),$ where $D$ is a (diagonal) noise term."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import torch\n",
    "import os\n",
    "import time\n",
    "\n",
    "from botorch.test_functions import Hartmann\n",
    "from botorch.models import SingleTaskGP, KroneckerMultiTaskGP\n",
    "from gpytorch.mlls import ExactMarginalLogLikelihood\n",
    "from botorch.sampling import IIDNormalSampler\n",
    "from botorch.optim import optimize_acqf\n",
    "from botorch.optim.fit import fit_gpytorch_torch\n",
    "from botorch.acquisition.objective import GenericMCObjective\n",
    "from botorch.acquisition.monte_carlo import qExpectedImprovement\n",
    "from botorch.models.transforms.outcome import Standardize\n",
    "from botorch.models.transforms.input import Normalize\n",
    "\n",
    "SMOKE_TEST = os.environ.get(\"SMOKE_TEST\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Set device, dtype and random seed"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "torch.random.manual_seed(2010)\n",
    "\n",
    "if torch.cuda.is_available():\n",
    "    torch.cuda.set_device(\"cuda:0\")\n",
    "tkwargs = {\n",
    "    \"dtype\": torch.double,\n",
    "    \"device\": torch.device(\"cuda:0\" if torch.cuda.is_available() else \"cpu\"),\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Problem Definition\n",
    "\n",
    "The function that we wish to optimize is based off of a contextual version of the Hartmann-6 test function, where following [Feng et al, NeurIPS, '20](https://proceedings.neurips.cc/paper/2020/hash/faff959d885ec0ecf70741a846c34d1d-Abstract.html) we convert the sixth task dimension into a task indicator. We assume that we evaluate all $20$ contexts at once."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from botorch.test_functions import Hartmann\n",
    "from torch import Tensor\n",
    "\n",
    "class ContextualHartmann6(Hartmann):\n",
    "    def __init__(self, num_tasks: int = 20, noise_std = None, negate = False):\n",
    "        super().__init__(dim=6, noise_std = noise_std, negate = negate)\n",
    "        self.task_range = torch.linspace(0, 1, num_tasks).unsqueeze(-1)\n",
    "        self._bounds = [(0.0, 1.0) for _ in range(self.dim - 1)]\n",
    "        self.bounds = torch.tensor(self._bounds).t()\n",
    "        \n",
    "    def evaluate_true(self, X: Tensor) -> Tensor:\n",
    "        batch_X = X.unsqueeze(-2)\n",
    "        batch_dims = X.ndim - 1\n",
    "        \n",
    "        expanded_task_range = self.task_range\n",
    "        for _ in range(batch_dims):\n",
    "            expanded_task_range = expanded_task_range.unsqueeze(0)\n",
    "        task_range = expanded_task_range.repeat(*X.shape[:-1], 1, 1).to(X)\n",
    "        concatenated_X = torch.cat(\n",
    "            (batch_X.repeat(*[1]*batch_dims, self.task_range.shape[0], 1), task_range), dim=-1\n",
    "        )\n",
    "        return super().evaluate_true(concatenated_X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use `GenericMCObjective` to define the differentiable function that we are optimizing. Here, it is defined as \n",
    "$$g(f) = \\sum_{i=1}^T \\cos(f_i^2 + f_i w_i)$$\n",
    "where $w$ is a weight vector (drawn randomly once at the start of the optimization). As this function is a non-linear function of the outputs $f,$ we cannot compute acquisition functions via computation of the posterior mean and variance, but rather have to compute posterior samples and evaluate acquisitions with Monte Carlo sampling. \n",
    "\n",
    "For greater than $10$ or so tasks, it is computationally challenging to sample the posterior over all tasks jointly using conventional approaches, except that [Maddox et al, '21](https://arxiv.org/abs/2106.12997) have devised an efficient method for exploiting the structure in the posterior distribution of the MTGP, enabling efficient MC based optimization of objectives using MTGPs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "problem = ContextualHartmann6(noise_std = 0.001, negate=True).to(**tkwargs)\n",
    "\n",
    "# we choose 20 random weights\n",
    "weights = torch.randn(20, **tkwargs)\n",
    "\n",
    "def callable_func(samples, X=None):\n",
    "    res = -torch.cos((samples**2) + samples * weights)\n",
    "    return res.sum(dim=-1)\n",
    "\n",
    "objective = GenericMCObjective(callable_func)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "bounds = problem.bounds"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Define helper functions\n",
    "\n",
    "Define helper functions used for optimizing the acquisition function and for constructing the batch expected improvement acquisition, which we optimize for both the batch GP and MTGP."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def optimize_acqf_and_get_candidate(acq_func, bounds, batch_size):\n",
    "    \"\"\"Optimizes the acquisition function, and returns a new candidate and a noisy observation.\"\"\"\n",
    "    # optimize\n",
    "    candidates, _ = optimize_acqf(\n",
    "        acq_function=acq_func,\n",
    "        bounds=bounds,\n",
    "        q=batch_size,\n",
    "        num_restarts=10,\n",
    "        raw_samples=512,  # used for intialization heuristic\n",
    "        options={\"batch_limit\": 5, \"maxiter\": 200, \"init_batch_limit\": 5},\n",
    "    )\n",
    "    # observe new values \n",
    "    new_x = candidates.detach()\n",
    "    return new_x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def construct_acqf(model, objective, num_samples, best_f):\n",
    "    sampler = IIDNormalSampler(num_samples=num_samples)\n",
    "    qEI = qExpectedImprovement(\n",
    "        model=model, \n",
    "        best_f=best_f,\n",
    "        sampler=sampler, \n",
    "        objective=objective,\n",
    "    )\n",
    "    return qEI"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Set environmental parameters, we use 20 initial data points and optimize for 20 steps with a batch size of 3 candidate points at each evaluation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "if SMOKE_TEST:\n",
    "    n_init = 5\n",
    "    n_steps = 1\n",
    "    batch_size = 2\n",
    "    num_samples = 4\n",
    "    n_trials = 4\n",
    "    verbose = False\n",
    "else:\n",
    "    n_init = 20\n",
    "    n_steps = 20\n",
    "    batch_size = 3\n",
    "    num_samples = 128\n",
    "    n_trials = 3\n",
    "    verbose = True"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## BO Loop\n",
    "\n",
    "Warning... this optimization loop can take a while, especially on the CPU. We compare to both random sampling and a batch GP fit in a composite manner on every output. The batch GP does not take into account any correlations between the different tasks."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch  0: best_value (random, mtgp, batch) = (-9.52, -9.52, -9.52), batch time = 11.39, mtgp time = 50.09\n",
      "Batch  1: best_value (random, mtgp, batch) = (-9.52, -9.52, -9.52), batch time = 11.67, mtgp time = 22.61\n",
      "Batch  2: best_value (random, mtgp, batch) = (-9.52, -9.52, -9.52), batch time = 12.64, mtgp time = 28.46"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch  3: best_value (random, mtgp, batch) = (-9.52, -9.52, -9.52), batch time = 9.27, mtgp time = 43.40"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch  4: best_value (random, mtgp, batch) = (-9.52, -5.68, -9.52), batch time = 4.93, mtgp time = 30.97\n",
      "Batch  5: best_value (random, mtgp, batch) = (-9.52, -5.68, -9.52), batch time = 6.38, mtgp time = 13.53"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch  6: best_value (random, mtgp, batch) = (-9.52, -5.68, -9.52), batch time = 6.23, mtgp time = 42.83\n",
      "Batch  7: best_value (random, mtgp, batch) = (-9.24, -4.61, -9.52), batch time = 6.24, mtgp time = 14.94\n",
      "Batch  8: best_value (random, mtgp, batch) = (-9.24, -4.47, -9.52), batch time = 9.92, mtgp time = 21.80\n",
      "Batch  9: best_value (random, mtgp, batch) = (-9.24, -4.47, -6.24), batch time = 14.26, mtgp time = 30.95"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch 10: best_value (random, mtgp, batch) = (-9.24, -3.26, -4.54), batch time = 12.88, mtgp time = 26.41\n",
      "Batch 11: best_value (random, mtgp, batch) = (-9.24, -3.26, -3.57), batch time = 15.02, mtgp time = 27.82\n",
      "Batch 12: best_value (random, mtgp, batch) = (-9.24, -3.26, -2.50), batch time = 13.59, mtgp time = 27.51\n",
      "Batch 13: best_value (random, mtgp, batch) = (-9.24, -2.75, -2.50), batch time = 24.31, mtgp time = 33.26\n",
      "Batch 14: best_value (random, mtgp, batch) = (-9.24, -2.75, -2.50), batch time = 11.34, mtgp time = 28.75\n",
      "Batch 15: best_value (random, mtgp, batch) = (-9.24, -2.75, -2.50), batch time = 11.85, mtgp time = 39.74\n",
      "Batch 16: best_value (random, mtgp, batch) = (-9.24, -2.35, -2.31), batch time = 13.90, mtgp time = 29.78\n",
      "Batch 17: best_value (random, mtgp, batch) = (-9.24, -2.35, -2.31), batch time = 11.93, mtgp time = 29.11\n",
      "Batch 18: best_value (random, mtgp, batch) = (-9.24, -2.35, -2.31), batch time = 11.93, mtgp time = 30.58\n",
      "Batch 19: best_value (random, mtgp, batch) = (-8.54, -2.35, -2.31), batch time = 6.24, mtgp time = 27.86\n",
      "Batch  0: best_value (random, mtgp, batch) = (-6.34, -6.34, -6.34), batch time = 12.88, mtgp time = 36.25\n",
      "Batch  1: best_value (random, mtgp, batch) = (-4.86, -6.34, -6.34), batch time = 12.11, mtgp time = 33.74\n",
      "Batch  2: best_value (random, mtgp, batch) = (-4.86, -6.34, -6.34), batch time = 12.00, mtgp time = 46.60\n",
      "Batch  3: best_value (random, mtgp, batch) = (-4.86, -6.34, -6.34), batch time = 12.92, mtgp time = 36.88\n",
      "Batch  4: best_value (random, mtgp, batch) = (-4.86, -6.34, -6.34), batch time = 6.23, mtgp time = 30.78\n",
      "Batch  5: best_value (random, mtgp, batch) = (-4.86, -3.68, -6.34), batch time = 8.35, mtgp time = 27.63\n",
      "Batch  6: best_value (random, mtgp, batch) = (-4.86, -2.66, -4.93), batch time = 7.50, mtgp time = 49.98\n",
      "Batch  7: best_value (random, mtgp, batch) = (-4.86, -2.66, -4.47), batch time = 11.87, mtgp time = 37.50\n",
      "Batch  8: best_value (random, mtgp, batch) = (-4.86, -2.66, -2.81), batch time = 11.41, mtgp time = 41.93\n",
      "Batch  9: best_value (random, mtgp, batch) = (-4.86, -2.66, -2.81), batch time = 10.96, mtgp time = 44.75\n",
      "Batch 10: best_value (random, mtgp, batch) = (-4.86, -2.66, -2.54), batch time = 12.96, mtgp time = 34.05\n",
      "Batch 11: best_value (random, mtgp, batch) = (-4.86, -2.66, -2.49), batch time = 13.49, mtgp time = 40.42\n",
      "Batch 12: best_value (random, mtgp, batch) = (-4.86, -2.66, -2.26), batch time = 13.78, mtgp time = 37.96\n",
      "Batch 13: best_value (random, mtgp, batch) = (-4.86, -2.35, -2.26), batch time = 11.03, mtgp time = 49.89\n",
      "Batch 14: best_value (random, mtgp, batch) = (-4.86, -2.35, -2.26), batch time = 10.55, mtgp time = 33.78\n",
      "Batch 15: best_value (random, mtgp, batch) = (-4.86, -2.32, -2.26), batch time = 11.93, mtgp time = 31.45\n",
      "Batch 16: best_value (random, mtgp, batch) = (-4.86, -2.32, -2.26), batch time = 11.73, mtgp time = 25.93"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch 17: best_value (random, mtgp, batch) = (-4.86, -2.32, -2.26), batch time = 11.65, mtgp time = 41.50\n",
      "Batch 18: best_value (random, mtgp, batch) = (-4.86, -2.32, -2.26), batch time = 12.54, mtgp time = 36.79"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch 19: best_value (random, mtgp, batch) = (-4.86, -2.32, -2.26), batch time = 10.62, mtgp time = 52.57\n",
      "Batch  0: best_value (random, mtgp, batch) = (-5.52, -5.52, -5.52), batch time = 11.71, mtgp time = 39.39"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch  1: best_value (random, mtgp, batch) = (-5.52, -5.52, -5.52), batch time = 14.40, mtgp time = 34.61"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch  2: best_value (random, mtgp, batch) = (-5.52, -5.52, -5.52), batch time = 14.71, mtgp time = 26.74"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch  3: best_value (random, mtgp, batch) = (-5.52, -5.52, -5.52), batch time = 12.25, mtgp time = 31.19\n",
      "Batch  4: best_value (random, mtgp, batch) = (-5.52, -2.74, -3.40), batch time = 11.82, mtgp time = 31.47"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/xarfuse/uid-228567/ef01a148-seed-nspid4026531836_cgpid45489687-ns-4026531840/gpytorch/utils/cholesky.py:44: NumericalWarning:\n",
      "\n",
      "A not p.d., added jitter of 1.0e-06 to the diagonal\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Batch  5: best_value (random, mtgp, batch) = (-5.52, -2.74, -3.29), batch time = 10.88, mtgp time = 43.98\n",
      "Batch  6: best_value (random, mtgp, batch) = (-5.52, -2.74, -3.29), batch time = 11.60, mtgp time = 32.99\n",
      "Batch  7: best_value (random, mtgp, batch) = (-5.52, -2.74, -2.97), batch time = 11.25, mtgp time = 37.71\n",
      "Batch  8: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.48), batch time = 11.48, mtgp time = 31.68\n",
      "Batch  9: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.48), batch time = 11.50, mtgp time = 29.20\n",
      "Batch 10: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.48), batch time = 12.28, mtgp time = 35.44\n",
      "Batch 11: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.48), batch time = 11.04, mtgp time = 31.78\n",
      "Batch 12: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.48), batch time = 9.10, mtgp time = 46.31\n",
      "Batch 13: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.48), batch time = 10.81, mtgp time = 55.77\n",
      "Batch 14: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.48), batch time = 12.39, mtgp time = 50.37\n",
      "Batch 15: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.40), batch time = 10.58, mtgp time = 26.95\n",
      "Batch 16: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.35), batch time = 11.29, mtgp time = 39.82\n",
      "Batch 17: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.31), batch time = 10.38, mtgp time = 42.04\n",
      "Batch 18: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.31), batch time = 10.90, mtgp time = 18.43\n",
      "Batch 19: best_value (random, mtgp, batch) = (-5.52, -2.33, -2.31), batch time = 10.75, mtgp time = 47.58"
     ]
    }
   ],
   "source": [
    "\n",
    "mtgp_trial_objectives = []\n",
    "batch_trial_objectives = []\n",
    "rand_trial_objectives = []\n",
    "\n",
    "for trial in range(n_trials):\n",
    "    init_x = (bounds[1] - bounds[0]) * torch.rand(n_init, bounds.shape[1], **tkwargs) + bounds[0]\n",
    "\n",
    "    init_y = problem(init_x)\n",
    "\n",
    "    mtgp_train_x, mtgp_train_y = init_x, init_y\n",
    "    batch_train_x, batch_train_y = init_x, init_y\n",
    "    rand_x, rand_y = init_x, init_y\n",
    "\n",
    "    best_value_mtgp = objective(init_y).max()\n",
    "    best_value_batch = best_value_mtgp\n",
    "    best_random = best_value_mtgp\n",
    "    \n",
    "    for iteration in range(n_steps):\n",
    "        # we empty the cache to clear memory out\n",
    "        torch.cuda.empty_cache()\n",
    "\n",
    "        mtgp_t0 = time.time()\n",
    "        mtgp = KroneckerMultiTaskGP(\n",
    "            mtgp_train_x, \n",
    "            mtgp_train_y, \n",
    "        )\n",
    "        mtgp_mll = ExactMarginalLogLikelihood(mtgp.likelihood, mtgp)\n",
    "        fit_gpytorch_torch(mtgp_mll, options={\"maxiter\": 3000, \"lr\": 0.01, \"disp\": False})\n",
    "        mtgp_acqf = construct_acqf(mtgp, objective, num_samples, best_value_mtgp)\n",
    "        new_mtgp_x = optimize_acqf_and_get_candidate(mtgp_acqf, bounds, batch_size)\n",
    "        mtgp_t1 = time.time()\n",
    "\n",
    "        batch_t0 = time.time()\n",
    "        batchgp = SingleTaskGP(\n",
    "            batch_train_x, \n",
    "            batch_train_y, \n",
    "        )\n",
    "        batch_mll = ExactMarginalLogLikelihood(batchgp.likelihood, batchgp)\n",
    "        fit_gpytorch_torch(batch_mll, options={\"maxiter\": 3000, \"lr\": 0.01, \"disp\": False})\n",
    "        batch_acqf = construct_acqf(batchgp, objective, num_samples, best_value_batch)\n",
    "        new_batch_x = optimize_acqf_and_get_candidate(batch_acqf, bounds, batch_size)\n",
    "        batch_t1 = time.time()\n",
    "\n",
    "        mtgp_train_x = torch.cat((mtgp_train_x, new_mtgp_x), dim=0)\n",
    "        batch_train_x = torch.cat((batch_train_x, new_batch_x), dim=0)\n",
    "\n",
    "        mtgp_train_y = torch.cat((mtgp_train_y, problem(new_mtgp_x)), dim=0)\n",
    "        batch_train_y = torch.cat((batch_train_y, problem(new_batch_x)), dim=0)\n",
    "\n",
    "        best_value_mtgp = objective(mtgp_train_y).max()\n",
    "        best_value_batch = objective(batch_train_y).max()\n",
    "        \n",
    "        new_rand_x = (bounds[1] - bounds[0]) * torch.rand(batch_size, bounds.shape[1], **tkwargs) + bounds[0]\n",
    "        rand_x = torch.cat((rand_x, new_rand_x))\n",
    "        rand_y = torch.cat((rand_y, problem(new_rand_x)))\n",
    "        best_random = objective(rand_y).max()\n",
    "\n",
    "        if verbose:\n",
    "            print(\n",
    "                f\"\\nBatch {iteration:>2}: best_value (random, mtgp, batch) = \"\n",
    "                f\"({best_random:>4.2f}, {best_value_mtgp:>4.2f}, {best_value_batch:>4.2f}), \"\n",
    "                f\"batch time = {batch_t1-batch_t0:>4.2f}, mtgp time = {mtgp_t1-mtgp_t0:>4.2f}\", end=\"\"\n",
    "            )\n",
    "        else:\n",
    "            print(\".\", end=\"\")\n",
    "            \n",
    "    mtgp_trial_objectives.append(objective(mtgp_train_y).detach().cpu())\n",
    "    batch_trial_objectives.append(objective(batch_train_y).detach().cpu())\n",
    "    rand_trial_objectives.append(objective(rand_y).detach().cpu())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plot Results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we plot the results, where we see that the MTGP tends to outperform both the batch GP and the random baseline. The optimization procedure seems to have a good deal less noise than the batch GP.\n",
    "\n",
    "However, as demonstrated above, optimizing the acquisition function and fitting the MTGP tend to take a bit longer."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "mtgp_results = torch.stack(mtgp_trial_objectives)[:, n_init:].cummax(1).values\n",
    "batch_results = torch.stack(batch_trial_objectives)[:, n_init:].cummax(1).values\n",
    "random_results = torch.stack(rand_trial_objectives)[:, n_init:].cummax(1).values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'Best Objective Achieved')"
      ]
     },
     "execution_count": 12,
     "metadata": {
      "bento_obj_id": "140024665143952"
     },
     "output_type": "execute_result"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "W0720 165607.890 font_manager.py:1280] findfont: Font family ['STIXGeneral'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.892 font_manager.py:1280] findfont: Font family ['STIXGeneral'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.895 font_manager.py:1280] findfont: Font family ['STIXGeneral'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.897 font_manager.py:1280] findfont: Font family ['STIXNonUnicode'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.898 font_manager.py:1280] findfont: Font family ['STIXNonUnicode'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.900 font_manager.py:1280] findfont: Font family ['STIXNonUnicode'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.901 font_manager.py:1280] findfont: Font family ['STIXSizeOneSym'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.902 font_manager.py:1280] findfont: Font family ['STIXSizeTwoSym'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.904 font_manager.py:1280] findfont: Font family ['STIXSizeThreeSym'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.905 font_manager.py:1280] findfont: Font family ['STIXSizeFourSym'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.906 font_manager.py:1280] findfont: Font family ['STIXSizeFiveSym'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.908 font_manager.py:1280] findfont: Font family ['cmsy10'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.909 font_manager.py:1280] findfont: Font family ['cmr10'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.910 font_manager.py:1280] findfont: Font family ['cmtt10'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.911 font_manager.py:1280] findfont: Font family ['cmmi10'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.913 font_manager.py:1280] findfont: Font family ['cmb10'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.914 font_manager.py:1280] findfont: Font family ['cmss10'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.915 font_manager.py:1280] findfont: Font family ['cmex10'] not found. Falling back to DejaVu Sans.\n",
      "W0720 165607.920 font_manager.py:1280] findfont: Font family ['DejaVu Sans Display'] not found. Falling back to DejaVu Sans.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZUAAAEPCAYAAACKplkeAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOydd5grV32w3zMzatvb7b6+1/W6YmNw/YwNNjIgAkQEY1MT4AsQkg8IJZSYYpKQAAmmhNDBAUyxDcIhCINcsDHFgAvuvdx+7/amOjPn+0OzdyWtdjXSzqy02vM+j55djWZGZ6TR+Z1fF1JKFAqFQqHwAq3ZA1AoFApF+6CEikKhUCg8QwkVhUKhUHiGEioKhUKh8AwlVBQKhULhGUqogGz0MTa8v+FjW/HRTtfTTtfSbtfTTteyxq+nKkqoLAPbsps9BE9pp+tpp2uhza6nna4FdT0LUEJFoVAoFJ6hhIpCoVAoPEMJFYVCoVB4hhIqCoVCofAMJVQUCoVC4RlKqCgUCoXCM5RQUSgUCoVnKKGiUCgUCs8wmj0Ar4jG4icAnwfOBiaArwGXp5IJq9ljU6wdCrZN1rbpNtrmp6VQ1EVb3PnRWLwfuAF4AHgZcBTwH44mdlmzx6dYOzySyRDSNCVUFGuWdjF/vRWIAC9PJROpVDLxJeBy4F3RWLyn2YNTrA325fOMFgqMFQrNHopC0TTaRai8CPh5KpmYKtn2fUfQnN/EcSnWCBnL4olMBoCsbTNrKaurYm3SLjr6ccBNpRtSycTOaCyedl77yWIHjg3vb7iAmmkWGNm/t6FjW5F2up6VvBYpJffnC8zI+fvosZkptnhoAlPfTeuyVq9naOPmqtvbRaj0O875Ssad1xZlYN3Ght90ZP/eRT/Y1Ug7Xc9KXsvT2SxaNkupnVUaBkNdXZ69h/puWhd1PeW0i/mLRer7i6Xq/iv8Rcr2/+inTJOnc7mF2y2Lgt1eJdEVCje0i1AZB/qqbO9dRINRrAAHpttbqEyZJg9lMlBNeErJmGk2Y1gKRVNpF/PXQ47v5BDRWHwr0Om8plhhZnKS/ZOSjW0We5exLA4UChwsFMjWcMaPmSYbgsEVG5tC0Qq0i6byM+AF0Vi8u2TbJUAGuKWJ41qzjM5IJjLto6mMFgrcNTPDH6an2ZnN1hQoAOOFAvYaMAEqFKW0i6byJeDtwI+isfgngCOBjwKfrggzVqwQo7OSvAnpvKQjKJo9nIYpSMnjmQwH8/m6jzWlZNKy6FeJkIo1RFtoKqlkYhy4ENCd8OHLgSuAjzR7bGsR25aMp4sr9Lm/q5HRQoE7pqcbEihzqERIxVqjbZZQqWTiAeCCZo9DAeMZyVzqz2RGsqVaCEUTKEjpyhxlAzuzWQ4sQ5jMMWqaHLXssygUq4e2ESqK1mF0dn7inmgRTWWkUOChdHrFfRxZyyJtWXTo+oq+r0LRLNrC/KVoLUZn5ifumRwUrOYKln25HA80QaDMMapCixVrCCVUFJ6SNyXT2fJtk02MAnsym+XRxXJJVgjlV1GsJZRQUXhKqelrjmY466WUPJxOsyubdbG3v0xaFgUPhVrWttmTyzGjilYqWhDlU1F4SjWhMpHx7/2klOzM5RbU4tlZKIAHjnZPkJKnMhkiLvwqAjCEICgEAU0r/hWCWdtmNptltFA4JEwMITixs5NeFbKsaCHU3ajwlGpCZSojkVIihPf5KuOmydNVtJEp26aVkvn3LUfACcFUPk9PxXWaUnLv7CzHd3QwGAgsf5AKhQcsKlSisfjhLs8RSCUTj3s3JMVqZSYnyVVxH1g2TGWhN+L9ex5YC/6KJUxntpTcn06zIxJRJWEULcFSmspTdVT4VfGSirKor0om0pLeyPI0lUcPWjw1Ov8elpQ8amepVgs4M6MRmSl/JajD1n5BT7j+cWRNyBUkuiYI6KBrYLSKR9LxH5lSsiUUavZoFGucpYRKtOT/o4B3At8AHnaEyAnAa4GPrcA4FauAaqavOcYzkm3LOPeTIzZPjpSff1IWsBd5S8nCJVHehMeHJf2dki19GgEXQiFnwr4pm/H0wpMKAdoKVaDpkILuvqLPZTEez2Q4WCjQqWl06joduk6nphHUWkX6KdYCiwqVVDJx49z/0Vj8MuDSVDJxT8kuP47G4tcDn3Fa9yrWMKWlWaqxnLDiXeM2jx5cqI9MysZMX+OzMJW1OaxPMNBRfZrOW7B/ymZsdnHrk5SwUik4B2YF1kGb7QMaoSWWgtOmyXTFNq2KL0sXgiHDYF0wSK+u++LvUqxN3Drqnw08WGX7vcBpHo9JsQopLc1SjVyhseKS+yZtHtq/8MQ5aZGRjYfUWhY8PSoZnpEEKo23sugDarUCw+kcPLTfZnOfYF2X+8+xWtKnLSX78nn25fMENY11gQBDgQBGhXA5FI1Wh7aTKUiGpyUHpyWTGbnAhj47oXOEZbGuSzDUJQjorSHQbCkZN012ZfLcM5wjW7FikICFxJSSAjYFKTGRpNOzhMf3LrKvjYnEkrKqmbYVuGiwjxdt7PXsfG6Fyi7g76Kx+GdSyUTpJ/1WYI9no1GsWpYyfc0xmalPqByctrlvr111cm9US6kkvbBpY0tjS9g9LpnKSPoW0bLqx2IfFrB4To8QEBQaASEICY2g0AhrGiGhYSB4IDvDY9ksM7lideqlyIo0vx+1YLT4PBIUdARhObLFlsWJ3LaLf6Us1nrLSousbZORFhlpk5cLp/aCbTNsFpiSJrZocCWRqdQPVw/5JwNNESofBr4HfCAai+90BPHhwCDwes9Go2h5pJRVBciIiy6PExnJJpf37sFpm3v2VBcoUkrPhMpqZSoLU9mVVqUWaoYTMs+vrRFGZR0h0wHALOlIYQJpj4a4HFpDYVr1uBIqqWTi2mgs/kfgYmArEAKuA36SSib+5P8wFa1CpgB37mxMkXeTWW/bkkcO2uwcW3zfWSxM14GJCj+wpeQBe4q77fGWNesomoPr5MdUMvEU8KloLG6kkglVIW+Nks43PpnP5sC0JMYido50XnL3bouZGpVVJutZFbchUsL+McHoZHPeP2vkeWpghHRoldkOXWDnNay0QaBgEK7mR5ICIQXCFmAX/7dNiaaX7yvk/L5z+80d32o84+guT8/nSqhEY3ED+AfgLcA6oCMai3cCnwLemUom1vavfA2RWYbVScqiCWyoipN5z0TRIb+Usx8nN2Vart01zXQafvIbnZ0HvA0T1rvyhDakEXrFokGAHjLROkz0iPPoMBFV3j67t4PCVAvlyUiQBQ27oCFNDVnQkJZAVpnYB3SDszaFuODoEGft0Nna7/7zHdm/l6GNmz0e/OrFrabyfuBNjhD5lLOtw4kK+zfgXT6OUdFCLEdTAXj4gM3TFaYt05ZMurSpT1NYs4avx3YLfvo7nUzOq9WuJLg+Q9exE4Q3zzZ8FjunMXnnejK7uup2TGzqFTzzMJ2zj9Tqmsgr0QTomsBwklJ1zckhcjGcsAHHbdRbJgptteNWqLwGeFkqmbgnGot/kqI5bDgai18K3KqEytohs0yddDYHs7nGxcLEGnTQmxbcck8Hf3qi/lJ9WtBCiyzU7AIDWbqOmSDQt7wvNLOnk8k71mNny8cmgPXdgk29YkHlAdvKc9ymCKceprGhR2NTr+DETRraSmWSKnzF7V16uJOTUslOJwJMsUZYrqYyR0HaZLHISps8litnr4Rl5aY0ipQwPg3TaUEmB+mcIJ0thiOblv8T4b4RwfDkwoKRmwZtBrrnn1uhAmZPBiuSLz468siAP5+Xlg0Q3tPP0Fg3x28uji1swGH9GtsGBIf1aYQD1T+bqYkZevqKA982KNixQVV5aifcCpW9wOnA7yu2v9TJYVGsEWbzkl22+1r2lSJIIslJG6uFjVhSwsgk7DygsfOgYNdBQTrbOqtoISTnnWJz5vE2mlYMsb7HnuQee2JZn+qgCLJehBdYjIJodAidDgw6nb+Bbg3tOMGRWoRANQeLC3Zs0Ng2qErItBtuhcpngP+JxuJfBvRoLP4OJ5P+lcB7fR6jokXIFiR5WzLTBo7ybA5uu1fj0d0ahYrLMS3Im60jRErp7ZS89FyLLUNF8SGl5Pf2GA/bjSXfCeBw0cHxWg/rRKiuci0btXBDAkUTcPIWjU29SqC0I27zVL4QjcVHgL930pQ+DDwC/FUqmfiB/8NUtALpPNgtrGG45aGdgtQfdGZbSPtww/HbbF5whkU4WMxw745IfpEe4RF7cSe7DvRogQVlxHUh2GKEOSXYTY9efy+WXi3AEYHyUvuagIAOxpzDXK9eyTndaSmB0sa4DSk+0hEeSoCsYdIFuaqFynQaUn/UeWRX/RNa0JCs75d0hCESgo6QJBKCYECuSCJ2f8cMWzd1MqdIdHfY3JAb5pFCuSmyQ9N4Tm8vmwIBNgWDDAYC6B4WixRC0GcYHBeJEGiw+vFIK2TPK3zDrfnrsWgsfhtwJXB1KpmY8XlcihYkk5dVCnWsDFOzcN+TGpOz7iZIs9CJUVIpUkp4eKcgV3B3fCggOWyd5PANxceGfkkzK8inZy1sJNOywAQFnkpPsTNfnnzYq+u8ZdMmNnncrEsIQb9hMBQIMBgIEFAVjRVL4FaoXAhcAnwC+Hw0Fv+xI2BuqCgwqWhj0nlW3MFuWfCHhzR+fa9Goa5Iq3DNPTRNcvaJNs88xqZ0nhRAOIhnQsSWkjGZZ2G93ir7AjlsstIig/OQJhOBHDPm2KJnGDIM/mbTJgaW0Vb48HCYXr3CUCYE3bq+oHqxQrEYbn0qNwM3R2PxvwVe4AiYa4CZaCz+nVQy8X7/h6poNum8rFpG3S+e3i/4xR90Rqe8n9A2D9q86CyLdX2en7qMGVkgae6jen/KOljiI9gSDPKWjRvpNurPY5ljQzDI9nBtQaxQ1KKuuzCVTFhAEkhGY/FTnez69zoZ94o2J1NYGU1lNgM33alz/1Pe25sCRjEc91nH2itizrrbmli+QFmCo8Jh3rRxI5FlXEy3YXBMJOLpuBRrl7qESjQWP94JI34lcCxwM/BX/g1P0SrkTYlp+S9U9o0Krv1l9ciscFBy5vE2kVDtMeTzWYLB8pV3MADbNxad7SuBLSW7pfucHjf0aDqbQkE2BAIcFYlwUkdH1c6ObglqGicu8xwKRSluo78uc0xeJziZ9d8AvptKJvb5P0RFK5B2qnn4Gf314NOCn/5Wr5qlfspRNuefarkWCOnZHB2dzS1uOCxz5Eu0FB1Bv6jtRA+hEUEnLPRDf41sjm29fZyysXETVyVCCE7s6FA97BWe4vYOfbPTpOtVqWTiPp/HpGhBMoWiMPFDU5ESfnOfxq/uWViuY32/5AWnW2xZt/riQXbL8tjZw0UHzzHWNXSutLQ4vNfbciY7IpG6/TDSluSmbcxl1G/LjEpmtNWfQDtHs65HCAh2agQ6WmtR4PaO2pZKJmQ0FjeisfgRqWTiSZ/HpWgx5jSV5QgVKWEmQ3l5ewm/uqe6/+TM4y3OP3VlfB9lSIkXbpDdFeVsDiMCVmOfXxDoDYKs83hDCIYCgQWZ8h2axpAewK51PgmFtE12yiY7aZGfrt6Nsx6mxiVyqn26ZTT7eoJdGp3rdDqHDDSj+WZMt0IlHI3FPwu8wYl6DEVj8X7gKuDVqWRiwudxKprMXCHJRs1f6Sz84CaDA+O1b3pNk7zwDItnHNUE7SQj0XdbiGWsxAGmAiaTR85XVBYSDn80gGE3tqLt0ySZGffHasCgk1dSLfnRxmY33vp7FM0hP2OTn7GZeKpAZEAn1O1uFRbq0Ql2er9icytU/sWp9fVy4Gpnm+08rnCEjaKNSTvmL7vBufZPj2uuBEo4KHn5eRaHb1h5gSJGbPQD1sIqmIsQMCCwyG/y0Z7y9pWbs0H6q9UscTUw6NTcDUoDBgIBhjzOpFe0PlJCetQiPeouRbl/G00VKi8FLkolE09EY3FJMbx4MhqL/1+g6T3qo7H4XwHfrPLS36SSiS81YUhtR2aZ5q/9o7UnuMEeySuea9LfXXNXz9B1MCyJ2G0hZiSHimSVWMCkLMqZgF4s7x4JCEKBhX1CSvlZd7lQOVlG2LiMeldpFzUsBHD0MsqnKBRe4FaobAaq+VHGgRWcAmpyAZTp9E80cSxtg2lJ8o7lxUI2pK5U9lPvChfLnghZdDhuHbB5/okmYRvwoPe6lgFh1naMbOvXCIzayICAfm9W9jlsnjLKS6jsyPufB9Kl60qgKJqO69pfwEXAzyu2v8F5rVX4g6pLVhtZ4mmtVerczNqMjVoU9pmYGQs9bSLqdBZbEsamysuHvPXoDOHKYKb9dZ12SQI50ENLmwE6QmBk3RRPqY/HAllKo6IHLYMh27tQ4MXoW0ZGvULhFW7vwn8GfhSNxa8FjGgs/hngmcD/AV7t8xjXNGbW5uADORd7Lo6UxaghKYtmncWid6bGbdL9C0vITqRtzBGJhUQ0oKWM5wR2SZ2RLsNeKFCaQL9PoZgPB8tNXzsKYYTPtYx1oLuybpdC0QTc1v66OhqL7wTeCdwPPMfpp3JuKpn4nf/DdM3j0Vh8EHgc+HQqmfhyrQPGhvdjW43Fj5pmgZH9exs61i2FtGRs/8o4rS3LZGp8bMH20VnIzGoUsMk1EI+/b6Y84W8gaJLL+Rt5ZNv2ku/REQAra+N1FXaJ5KHe8rNum5Skc4010ZrDtm3SM4ufo1fXmM4vb/GxUix2n61WVuv1yE5BTlu42HE7rw1t3Fx1u2t92REel7rdf4XZB3zIaXesA68CvhSNxTtSycQVSx04sG5jw286sn/voh+sV2QnLQr9KzNZTI2P0dM/sGD7hLSJ2BINi5Bd/1gmJ8pvs/URQSjkr48hl8ss/h4CNvVrBH1Y2O/Wc6T1+UVASAp2BPoxFunX7pb0zDQdXYu7L7eEw3SuEk1lsftstbJar6d/KED3xoVVrZc7ry0qVKKx+IdSycQ/Of9/bKmTpJKJDzc8gurv3QtsqrVfKpl4yPn78wp/z8+isXgIuCwai382lUz4V9HPZ2omp60AOUc5afRDHMmWm5nWhZt7TV0hfBEoAA9VmL6OLoQxfDZ9BYRYNQJF0f4spam8Bvgn5//XLbGfdNoLe8nFwFdd7LfUr/Vap/Dl9tUcBSZbQBzmTacfeoMu7ZFcuVAZDDfxooR/vhSAhwPlJreViPrqVQ56RQux6N2YSiaOK/n/iBUbUfH9vgZ8zaPTNX+pvwxss7nDl/Z8OHEjvVSkhJGKisNDoeYJle5wMd/ED6aEyT6jPIv+2IL/JZFV1Jeilai39P0Q0FG5PZVM7PR0VN7wF8AI8HSzB7IcZLP69zrkSt6/EVEwWRCYcl6oRHRJZ7PmQJ+1lErT12FmkC7pr1kqrGmEVG6KooVwW/r+z4CvA0MVLwlHE2iqQTcai//QcdLf44zlEufx9tXsT6EFfCr5Ek2pkbpflVrKYKi8da9rRFHDqBKsUh0TQhU+yI7g0lnwy+HBQIZfRMpL4O0o+G/6UlqKotVwe0d+GrjRKSCZdbH/SvMw8EZgqyPoHgBen0omvt3sgS2XpmsqJRHEjUjnSif9kEsnfTAAXUFBQIeAUfxbjyxKGzYdXf6v4G0kN0QmuTWyMNz3+BUwfS3oKa9QNJl6yrS8PpVMtGQThFQy8UHgg80ehx80W1MpFyoNaCq5ysiv2qIpHISN3drKl7yvk1lh8YOuUZ4ILAyzPifbxQardkOu5dCl6xit/iEp1hxuhcqfnBDfXT6PR1FBK2kqjUR/1euk7wjChh6tMROZz0gkY5rFbiPHLiPP/cEM01r5F6RJeGG6j7NzXb6PR5m+FK3IUnkqpdkvHwL+KxqLf94Jzy2bGVLJxKoN2W11mq+plPpU6kPKhZrKUuavjhBs6G4tgTKmmTwSyPBoIMsuI09aW/xT6LI1Lp0ZYrvpXRtjXQi6dJ2QrtEdKDqJhPNQZVkUrchSS53dFeG4AnhxlW1Nd9S3M03VVOR8ODENhBTPmJAtqawY0CQ9gern6ArDum7NtzTBHDYH9UJZocfFyAubx40cjwQzDOvuLL7bCiEunRmke5nRXgKIaBpduk6XrhPRNBCCqbRBT9Bfc5pC4QVLCZXnreA4FIuwXE0lk3dvtErnwcjP723OFaGcG0ud5q8FTvqQTWcYQhUtTzUNesP+iJMJzeQ34Wn+GJolL7zX+nQJZ2e7iWZ60ZcQiVtDIXqUuUrhI1ZeMvJwjsyYu5Xo/j6dUI/G5tMirDvOO+16qeTHW0qfR2PxLiCQSibGnefbgIlUMuFB9wvFYjSiqdi2ZGxWMjwjyRZcHOCQmdWI5BYz79TvURmtYvrqCAp6fBIgpRwMWPyhc5R7g2lsD98uKAVbzCBbzSBbzRCHm0E6XWgngVay6SnaivSYyZ7fZ9h7ZxYzU//CyQiKlREqpURj8WcD1wN/C/zA2fwXwAejsfgLU8nEHz0bkaKMappKOi+ZzUkMHQK6E3arQcGGkRnJ6KzE8ths1kg48XAVJ73rPBMXjGomP+0YZ7TCRGUjGde9+QA0CdvMEMcWwhxdCLPBCtCIkc5oI6EyO2ySHmn8803P2OS6VkdFZTc063qsvGT/PVlGH823VN0Qt/r4fwBfAf6nZNvnnez6K5xS+AofKK39lc5L9k9JJtOld9DK3E0N5ahU0VS8FCo/7hzjySrhvIvRZWsMWLVveQEM2gGOzYc52gwTlssP220XofJ4aoanbvWiYUC7GTja7Xoax61QeSZwQSqZOLQ8SSUThWgs/m/AP/g3vNXHngkbu3IGFsXe5qGAIGxA0HA3wUi76NPIOMJkIt285Uhj4cSVQsVG82hytZDsNNwJlCHL4NxMN6fmO32vGFwNXYiaHTZXA4W0zdO3ed2BRuElgQ7B+pPCGC5MzOFejVCXztAO70xf1CFUZp1qv49XbD8OPO9ztGrZN2lz/97aa3pNg5BRO0PcLkjS+yxydfhF/KJeJ33GhFlz/gp1IekPSjSPVJUJzawZybW1EOS8bA87CuGGTFZe0S7+lImdhZaomq1YSPdmg8POjLDh5DC6y949/dsCdG9e2E9lubgVKt8BktFY/D+BJx0LwXGOj+U7no9qFVKwJA8fcPeLs23I5F3sl5MtIVBowPxVafoaCBVNX16Zv0Yq/CiHmUFeMTPfKEnOpFnX0evNmy2TdjF9TTxdfjN2DOl0DNUfQm0WChgB7yezZtHM6wn36mw8JUzPYUbLaMNuhcoHgGknCXKuqOQI8AXg4z6Ob9XwyAG7LKfDE1poVVh/OHH1THrdo/t+WC+f4DaaAYbs+R922m6d8iXtI1TKV0JHXtDJhpPrr29W7JTY5+HImku7Xc9ycduj3gQ+BnwsGov3AGIulNh53iLr6eYwNivZM+G9v0ParRPSUW/iY9VMeoFn2fIjWrkEX2e37sq3HYSKlZdM7yn/zPu2t+5nrmgedWdjpZKJKYrC5BTH/PUqYPHm2W2ObUse3O9T2nuT636V4kXio5eRX5XmryEXUV3Noh18KpO7y/0pkQGdULcqpKFYSL1NugJOq9+/Bc5y+pe8x7/htT5PjEpm/QpRbyHzV7060wLzV9j2zPRFFfNXKwuVdtBUKv0pfduUlqKojtvkx8OBtwJvAjqBIPDSVDLxU/+H2Lqk8/BU2r+Zv6XMX3WIlbwFk4V5TUUgGQxJz0rZZ4TNbElhR11Cn62Eip9MVvhTlFBRLMaSv8RoLH4R8HdAzCl//1Hgu8BTwL0rN8zW5JERDdvPPkwtZf5yz2iufBLtC0oMzUt/SvmqedAylqy71WxWu/nLtiSTu5Q/ReGOWsu7652Q4dNSycQ9cxujsbj/I1sFzOYFER+FymrVVBY66YsiySufSmXl4KEWdtLTBprKzD4Tq6TQaLBLIzKg/CmK6tQSKjc4jvj10Vj860CiVbs/tiWt5FOpQ75VOunXhYoHe+VTGVlF/pR2yKav5k9Z7dek8I8lrdypZOIi4ASnB/xXgL3RWPzTjjBqnWV0u7JKNZX9mfLbatDRVDwzfy2I/GpdTWW1m76okp+i/CmKpajpOk0lE4+mkol3AFuAjwDPd5z1343G4pdEY/HWXSaucprdSrgUy6VQyVvw1Ez5bbW5wzF/eWT/qvSprGthTWW1m76klCryS1EXruNxUslEOpVMfDGVTDwDuADYB3xb9a33kRbRVOopJvnEjI4l5yfSvqDNkIfmLxu5oNR9K/tUVrtQSQ9bFEoKmeohQdfG1hXiiubT0N2RSiZ+CfwyGotvAd7s/bAUVJS9byb1iLZHJssduMf2WIfMXl4oKuOaVVZIstPWiHhQmt4vVrtQWaClHB5AeJnFqmg7lrXkSCUTexyTmMIPWsT85Va22RIem6oQKr3zF+HFXLTQSd+6WgptIVSUP0VRH627xFO0TEixWyf9nrRGukSNCOuSrZ3zIkl4cLdVOulb2Z9CGzjqKzWVXiVUFDVQQqWVaRHzl1uhUmn6OqrbKvOj6B5MsMMVTvqhFs6kZ5VrKtlJi+xEyaJAh54tSqgolqa1f5FrnRbRVNyO4pFK01dPuf3OG/PX6gknZpULlUotpWdLwHUDKEVz0HTo3hxwFb4f7PFHp3AtVKKxuAacB2xPJRNXOts6UsmE6vzoA9KWLeOodxNOPJoTjJZk0mtIjqoQKl7Mr5U+FWX+8g8VSry6EAIGjwkR6W9utQO3BSW3AT8DjnXcx1dGY/HtwO3RWPzCVDJxn/9DXWO0iEDBZS+VRytMX9u6bMIV97a+zIVRRtjM+FxIsjBsk7m3gD27fC1RAA8Y3hSgKOQtAsFJT87lltHHKpz0qt5XS9N7eKDpAoU6NJVPA7cD5zj5KaSSiaeisfhXnNcu8neYa49WcdLjMk+llulLCJZd8rEy6XHA9raQpJ2VjH8/64lAmSPjaQifXz0WXCCgd6sSKq1K55BOjw/95hvB7drxbOCdqWRiosLE/q/AGT6NbW3TIuHEuFCa0ibsmi2/lY7prfCn+BL55e2PKPuI6alAaSe6NhgEIiqupxUJdWkMHBVs9jAO4fYu6QYyVbYHgebrW+1IC5m/amkqj03pyBKNYX3Ypi9Yfow/OSremh9wl3QAACAASURBVL7yT7eQJG8xjrygs9lDUFRBDwoGdwRbKiHV7a/yduDdjmYCRT9LV4lZTOExrWT+quWof7SG6QuPhMqCkvceaipSSvJPl0vyrucE0HsbH3iXbtBreLPmyszOEulc+YldCEHvtgCRPr1kW7H8faMECxDqbh+tp2nXI4rBE0awtT5Lt0LlPcD10Vj8LUAoGovfDRwFzAIv9HmMa5MW0lSWGoppw+PTi2fRz+GJpuJjIUlrTGLPlAhPAzrPCiCWUbBsXTDIYMAbwTc1nqGnP+LJuZaLHhJsOKnxRkL6fo2hjX52t1tZ2u16lourX2Uqmbg7GosfDbzOiQCzgS8C300lE9N+Dc55z/cCZwEnAb9KJRPPXWTfE4DPO/6fCeBrwOWpZGJV2jSk1TqaSmnw13hOMFLS2XE4q5G35593GTabIgvF0HKFio1kzMdCkpWmr+BWfVkChVUeTrwURqg9r0vhDW5Dii8H/juVTHzR/yGVcaLTyvh3jv+mKtFYvN9pKPYA8DJHi/oPx2d02coO2SNaSlMpSpVfHzC4ef/SDsFjSgpIlrJcR/2EZmH6WEgyVylUti3/3Ks58XEpjFBrmVsUrYXbu+P1wKPRWPxX0Vj8TdFYvMfncc3xk1QysTWVTFwM3L/Efm8FIsDLU8lEKpVMfAm4HHjXCo7VU1rJp2IjyVtw64HamkE10xcsP5x42MdCklJK8jvLxx3atnxfSLsKFV1pKoolcCVUUsnEEcC5wB3Ax4B90Vj8u9FY/AXRWNy3OyyVTLhdr78I+HkqmZgq2fZ9R9Cc79Pw/KVFjHZzkV97M1pZn5Rq9Adtjuiq/pXpy7R/jWiVTnrv/CnmQRtZEtsoQmBsUJrKYijzl2Ip6mnS9dtUMvFO4DDgJcAk8I0WadJ1HPBQ6YZUMrETSDuvrT5aRFOZExG7K/JQ+oI2R3dbhx7PGixwyRE5jEXuqOX6VBaUZ/Ewk74y6iu4VV92iKYOaG0qVPRge16Xwhvq/mWmkgkZjcWngSlgDNjmz9Dqot9xzlcy7ry2KGPD+7GtxhwYlgVTE6MNHVvz3JM2coUT8SzbJjNbHndRwCanmeyc7i7bfsZAmlP6FmZ45xZJ+s7pkrQTfGAjyWqSnFb8mxWypn1sl5Yte949a5LOLB4jYts26Rl3MSSZx0W5gW5TgfRMYalDahIUgqm8dxnwlmUyNT7m2fmWQ2hSMJNpXLCYZoGR/Xs9HVMzWavXM7Rxc9Xt9RSUfBbwSuBiR1u5wclbSdRxjl5gU639UsnEQ7X2qUK1GVjUKrI7sG5jA29VRH/qAF09gw0fvxT5yTyWtbI2sMzsNJHOcuGhYRG0cuzLBggOZeh95jB6R4GHdXhULC30Sl/VBAgEtpDkahznhi2Bbjr0xf0q6ZlpOrq6F3390BgtyfS+8pqoXcd2ElhGHgZAp6bRE/EuBHhqfIye/gHPztcoQsD6rRHEMrSwkf17F52QViPqespxG/31OLAduAv4nBNKfLCB97sY+KqL/eq9Y8eBvirbexfRYFqeVgkptp0KxBkL1p+5H6Oz6NtobB3vzTXpEvo9Mn8V9tvIkrqJWgcYQ8s37xhe1KVpQfSgWJZAUbQ/bn+Z3we+3aAGcYhUMvE1J3/Eax6q9J1EY/GtQGelr2XV0CIhxdLxp2gR65BAaTZHFsKeFZJckJ+yTfdk0mxXJ72K/FLUYlGhEo3Fz08lE7c4T28ENkdj8ao6USqZuMm3EbrjZ8B7o7F4d0ky5iVOvbJbahzbmrSKo15Kdqd1jK68i73dE7IFEVnMNQlLzXXEyKBlcGGm17NxVBMqXtC2iY/KSa+owVKayvVOSC6O/0QuYpaSfhWVjMbiHU7yI8AWoCcai7/CeZ4saRD2JeDtwI+isfgngCOBjwKfrggzXjXIFgkptpHsntUxtpQbvLZPhzlnfzWLYznCuWUOGxDoWvF5SAo0D0vWN4o0Jfk9FZFfHgmVdtVUjHB7mvUU3rGUUNlR8v8RKzCWaqwHrqnYNvf8COApiprSeDQWvxD4T+Anjh/lCkewrE5aRFOZMSUjOY2eCk2lLxcgbLufgLuEhgf+eU/J77GhxKKn9Qj0Pm+EQbsKFRVOrKjFokLFyfOY412pZOIdlfs40VyfdzLuPSeVTDzl1mmfSiYeAC7wYxwrjZSt00r46XTx49e7yzWV3rx7R7nQlp9R7weVpq+QR/4U2tn8FW7P61J4x5IzQzQW7wMGgDdHY/HPVJngTwD+wi+hsmZpUKBIVz0ai1+icDnNPz1b3M/oKhcqPQX3QqWFWj2UUVmaxYt6X3PobSpUlKaiqEWtmeE1wGeczPvHqrwuHCe+wksa8KfMyAIT0l2gr4HGBi3kSrA8PStAyIVCpQ5NpTS6VspiiXlrXGKO21jjEmvK9jzazTQFeSO75D6FvRX+lMO98adobSpUhFAlWhS1WXJmSCUTX4jG4lcBBxbpQz8L3O3f8NYm9RWTlEzIAjPSfbivic2kLNAnlq44bEnYk9bQIyZCl2y/P8Czfhlm4IBB0K4vW/yAMxdJc6XqmgnMOt5IHxDoPd5oKu2ao6IFREt1GFS0JjWXm6lkYiIai5+aSiYejMbi4VQykWXeNDaZSiZazP3aBrhctUskYzJPpoFQsRlpEhY64SUC9/ZkIG8Lwl05zro+wrNvajxDvNVvEq+ivmik9tEqQWkpCje4XVLlo7H4w8CflWx7E/CI00hL4SFusultJMMy15BAmWPczh/qlVKNJ2YhbNq8+q7csgRKqyNC0HWWd6X0A22qqajER4Ub3C6qPuc0yipNJLwSONqJ/nqRT+NbU+SkxZQsIE0bWcM/kpY25jIdERaSCZlnQISqvn7wgMnbHpxmyLu6iOBM4nq/htEv0Ps19D6BCHg7YeWzGYLh2oJQGEVfiubhhNmu4cRKU1G4wa1QOQfYNGf6omgWG43G4u8C2qc8Z5NJYzEi8wjLRvcp+1F/0kLfaUFBIqxipeVxK4NmFSv75rX5Umnnj1sEKuTWbJfNwRcINtRRLLG3A/o7nNW7ABHE9/pRcgYiXc0xRLWvUGlPDUzhLW5/dZZTnLEynGZ961SpWv0U5pJTfPpEjbsKBO9c6NCXJb5zq8S5XWkQ2n+4yc9eN83zx4cg537i1MMCbQ3lN7RtjorSVBQucCtUrnNKoPwb8KQTSnwc8AHnNYUHmI5/Q3idTS8lgbtMAnc1XhDy/tOz3BJPYxvQc6A+DUCssQVuO4YTo3wqCpe4nR3e4fhVfuQ494WzuP0O8E6fx7hmKOCDpiIlgTtMAn9qTKCkDfj9S2a576wcCIiYGgFZn5TwOwo1XMUxbguNDk1DE8VsnHq7MNpSYpY85sS8cISGIQSa45QPCEHQ+RsSom01FZX4qHCDK6GSSiZmgDdGY/F3OMUaAZ4oqQis8IBD5q9luFOkhOmCYCwvGMvClvsKHPF4uUAxDXjwhDC5kMDWwNYEtg6WbaFpxdDax6Z0pkyNye1Z+s6e7zhYT3mWOXSf56Ijw+EFPpqpbMbTJlm2LIqVdm0RXAs9IND8/iIVbUE9nR814FnA9lQycaWzraOkUrBiGUgpD5m/DmkqGYnxkInILm0OM+1ikuJMQSNrgXQy5QfyFkdMlEeRZXXBN4/uZlcwMO9Msap03QoWH529jWfSz+HnRKzjv9OfNSxM5lCmL4Vb3HZ+3Ob0LDnWmYKujMbi24Hbo7H4halk4j7/h9remKV1u+yiyhH+WQ5tvLZ/JQAc7UK9yeiCbxzbw55O94LB6G685tccfqZttGv2equhnPQKt7j9RX4auB0YmltHOxWEv+K8plgmpTknwpaISelKoLglrQu+XqdAARY052rE/OWnT6Vd/RethhIqCre4nSHOBo5PJROT0Vi8dKb7V5Wn4g2F0sx2C8SUdwIlFxTccUYHR/ZKjlyiu7xlFtCN+UBiAezuz1Oa+9iI+ctPU3y7Rlq1Gsr8pXCL2xmi22nNW0nQr66Paw2ztIGKDdp0uVCxNmhYR8x/1PsygnvGyr++TkPy7CGTkD7fo0AaYB2u86yIXTOsLJfLECpJrjeF5FvBErOahO56zV+iWN3WL9o10bDVUImPCre4nSFuB97taCZQ9LN0lZjFFMukTFOxQUyXCwDrMA3zxOLXtXtW49uPh7A2zE+oEV3yhmOyGCHDsyLA0wGzrINOl6ljyPomcb/DiZX5a2VQ5i+FW9wKlfcA10dj8bcAoWgsfjdwlFP6/oU+j3FNUKj0qVRoKrK7+KOeyAuueSqEVTK5a0Lyiu05BkLeJk1OBctDkRuL/PJwQFVQmsrKoMxfCre40mlTycTdTvHITzhJkDc6guYY5zXFMqn0qWgVPhXZo7EvLfjvx0LMmuU/8Bcflmdbl/e1XSa9ECo+W02UUPEfTUflqChc43qWcBIgv+jvcNYu84mPEuRCTeUh2+BHj4UwK8xP56wvcMqAP8UnpwLlQqW3kXBipamseoyw8qco3LPoLBGNxVOpZCLq/H9rjfNIYAT4bCqZqLWvogJbSqzSxMcMiBI5UTDg6n2hBR7v43tNnrfRXQvhRlDmLwWqPIuiTpaaJZ4s+f9xF+c6ErgK2OrBuNqa/VM2T43aWI5yYiI5aDs/XFPQs0dwUsn+wwF9gUA5Y6jA8zcXfI2sWmj+qj/Qz0+hIlRI8YpgrKEK04rls6hQSSUTby75/w21ThSNxXVgysvBtQt5U/LgAZt791j8aY/NwelqDvX5r+LUUatMqIyF5idzDckLD8tz2qC/jd4LwiZjlAQPNBJOrBIf2wJDaSqKOqin9tcLgJcDhzt9VXYCV6WSid9TFDwW0OnraFuYx4YtfnqfydNjElkhM6ZzkkIdMmAgV77zmJMjENGLUV5+OOUrqTR9dRd0NOqfXDQfpYoyfa0MutJUFHXgtvbX3zlRX38CnnAsD88F/i4ai781lUx81f+htib7Jm2uuavAHTu9m+gHcuXnGg/pDIVsLjkiR7/HYcOL4UXkFz5n0yuhsjKoxEdFPdTTT+XVqWTi+6Ubo7H4a4GPAGtOqByYsrnq7gC/fjpHvT21dA12rNfo6yhOihlpkZ+L/ipIDnusfELfstHmvGOyhFawdsFUoFxbalSo+JpNr4pJrggq8VFRD25nik3AtVW2fx/4ssdjammms5Iv3JLnK7flyRQqG+4uzmCn4BlbNJ6xReeEjRrhwPwPdac9y6zTk16M2XSmyif0E7fZSH1lJ9AFkV8N+FPwW1Px79QKB00HzVjel1goFJiYmDj0PJs3GR4e9mB0rUG7XU8mm6NQKBAIuJ/fSnH7u/wdcApwR8X2k4A/NvTOq5Qf3V3gszfnq7525JDgFacG2NhT/iPUdUFvePG+H2apEyYr0Uo61EgBsnPlV4qVQqWR6sSsUk1F08GItJYWFMhBsGvlx7RcLaVQKDA6Osr69evRnO+rUMgTCAQ9GmHzabfryfV0Mzo6yuDgYEOCZak8lfNKnn4b+O9oLP5N4CGnp8rxwBuAjzc6+NXIq54d4Eu/yvPU6Lwg2NAtuPi0AM8+XGuoYVRpiRZtvNyfIjuF760TJwMFnujJkJF5dL1Yk3g8uPzmXPjc3KqxdVRtIgM6g0eHXOy5chj7NYY2hps9jLqZmJgoEyiK1kfTNNavX8/IyAjr16+v+/ilZopfVtn2qSrbrnLMYGuCoCF4/0Uh3vq9LD0hyZ+fGuT8Y3SMBqOcLCnLagdrYxVCpdtfgTJjmPzv4SPkDoUP5xbso9nQaTbm0PHTaueXpqIyyL1FCZTVh6Y1tkCmRp6KqzshGouvOdP2S042mMiECJvj9A90LOtchYpy9NpE9UKSfnHfwEyJQKlOT8FoKJwYv81fPp03oEJoFYqGqSdPpRc4winJ8rhTC4xUMmHWPrq90DTBX54V5LrfLf9clUJFVJi/bB+FSka3eLg3XXO/HRONpx+txpBipakoFI1TU6hEY/FNwBeAP3MacgkgH43FrwXemUomRlZmqO2JWZEpuZKayoN9s1ja/PtFCoKTx3vK9hnIBdicbtC/4GODLkP4d3IjojQVhaJRlhQqjnbyO2AS+L/AA45gOQF4F/CbaCx+2pzW4gfRWPxo4L3AWU602a9SycRzq+z3V8A3q5zib1LJxJf8Gt9yWWD+mqzwqfT4s2ouCJsH+8u/tuOHw5w03eXZe6xGLUUPCFXmfQX4xQPNMXBcdMKas9avOLU+4fcB9wAvSyUTpbPd7dFY/FvAj4F/AD7s4xhPBGKOcHMTt3dBRevjJ3wc27IpEypSLuhN75f565G+NDl9/r2CluDosbCnIVWrse6XKp6oqMa733cZ99x7Px/7yAc5+8zTF7z+rvf+I/fe/wDvfdf/4zOfn1/DWpaFlBLDmJ9q//7tbyN6YXFd/Mijj3P1tQnuvf8BZmZmCYdDbN92OLEXRrnweecfOubnqZv49ys+XxbiGw6F2LZtK6++9GJOf9Yzfbz6+qglVF4GXFohUMCp9RWNxS9zor/8FCo/SSUT11HURq4Fhmrs/wc/NSevKW3OJaYlomQBJwOAD5GtNpL7KrWU8S4CtrcTqp9BP/75U5RQUVRnoL+f5PWpBUJl77797Nm7D4AN69eTvO7qQ6998tOfY+/efXzm3/91wfl+9evf8q+fvIJLLo7z1je/kaHBAaamprjpltv4zOe/xK7de/ir17267Jif/Oh76HoxEjOdzvCTn17Phz76L1zxqY9z/HHH+nTl9VHrZ39YKpm4d4nX7/G71H01gdZOmLIkR2WkSjixD5Pn4z0ZZkvKsOg2nLAMZ/xi+KmpKCe9YqU5+6zT+eMddzEyOla2/Rc33MyZZzyrrnOl0xmu+Ox/8dIXv5C/fO2rGBocAKCnp4c/f0mMD33gPWzZvGnJc3R0RLjk4jibNm7gtt94EDXkEbV+QSIaiy9lEAk4iZCtxOPRWNyMxuIPR2PxtzR7MEshpcQs0VS00YrIrx7vJ06J5N6B6bJtx0x2ErG8LyzmazixTydX4cSKxRjo7+fkk07gFzfcdGibZVncePOtXPDc85Y8tpI77ryb6ZkZXvHyl1V9/YzTn0X0wue5Opdt2xj6ChYGrEEt89dDwPnADYu8/nzgER/G1Qj7gA8Bv3eCCV4FfCkai3ekkokrFjtobHg/ttWYMmRZMDUx2vCATSSzWvbQ8859Fa9HbHK5zMIDl8Hu7jwToXkbm5Bw7AGDXD6DbXv7fgEJac0fRTNrGEwtYV+zLJOp8bFFX1+MwJQgbbaeYDHNAiP79zZ7GHWTzZsUChVljaTEtPzrWLoUhUJj96OUEsu2iF74XL511Q94RfwlCCH4wx/vImAYnHB80fRkmoWy65W2jZRywWewa/duIuEwvb3dCz+fKliW6Yw/j20XBcjMzCz/m/w5wyOj/J9zznR1HpcXS6GQJ5vJLHnPDW3cXHV7LaHyPeAz0Vj8BalkYk/pC05U1meBRSfsajgRZUvrdUWz10P1nDeVTPwc+HnJpp9FY/EQcFk0Fv/sYma0gXUb63mbMvSnDtDVM9jw8Rlp0mHP54mEptLA/ISv9QcIhWpHq0gkO7uyC0qrVOPJnmzZ8yOmIwyJLghBLpchFIrUfR2LEYlAh0/1qvrCYTqWWJ1NjY/R0z9Q93k3bIu0ZPTXyP69i/6IW5nh4eEFdbEKhTyG7leRnaUJBBqtYSfQNZ3zzzuXL37lG9z3wMOcduozuOGmW3j+hc8lGCw6Pw0jUHa9wslMr/wMgsEgtrQxjMChzPXx8Qle81eHeiNiWRaf/NfLOeXkk9D14rhffslfHno9Eglz5PZtfOJfPsqxxxzd0HVVY66WWTgSYWjdurqPr/UJfx54KfBwNBb/HvCgc8xJwCuBG4Ev1vmeF7ssle/FL/taZ5zbWzEKrNRJD6CN1h/5ZQrJrRvHeaqnMQ3j5DHvQogr8bNBlx/RX7qhwokVSxMMBLjweefzs5+nOPrI7dz+hzv4v294bd3n2XrYFnK5PHv27uOwLcXFQn9/X5mTPxqLUzFFlDnqW5UlhUoqmTCjsfhFTk7Kq4HXAAXgYeA9wH+lkom6uomkkomvAV9b9sjrY2U6Wy1CYU8eEdTQB3REyaRVkBU5KpV1v2r4VLK6xQ1bxjgYaUzt3TIbYjDnX3XV1ZanopIeFW6IvSDK29/1PnYcczQnn3QC69bVCkhdyGnPPIWBgX6u+t41vO8971jwum2v3vikmrpgKpkoAJ9wHquNvwBGgKebNQBpS6wJC2lbmAcL6P0G+qCOFtLKNZW8RCvJUalV8n4yYJI6bISpYGNxEroNzxrucbFn4/gVUqyzeBuB5aCaUSnccMQR29i2bSvf/cG1vP1tjcUCBQIB3v+ed/Lhyz9OJBLm4pe/jE2bNpLJZrnv/gf53g+uZWhwkA0b6q8S3GxaPr00Got3OMmPAFuAnmgs/grneTKVTKSd/X7oOOnvceadS5zH25sZlmxnbOYUEmmDOWpijpponRoFLY8215xr1H3J+wPhHDdsGVtQCLIvZ7BtprZPRLcFW33WUnARWtgovlUnbrEeKu3Mas9sf9ELonz9ym9zzjlnNmwIeeapz+A/P/sprvnhj3nP+z/M5NQkhmFw+NbDOOesM3jJi19IZ6f3of5+I6RsqmWoJtFYfDvw5CIvH5FKJp5y9vu4o5lsdfwxDwCfSSUT367xFg1/ANf97gCRGo76woEC5sHqJSkOyOwhE9iusRl2JOZtpbuPKvDjt0xXPa4am2aDXLB3kJDd+MToqaNewLYBzZfS952axvbI0uNsxFE/eEyQzqHWnOxWs6N+XYWzt92aWrXr9VT77iqouuptzV9QCY7QqGmXSCUTHwQ+uDKjco89s7iSZDkCfTxY4ICVZQfzq5LJAfdmraMmI5y7vx/dk9gGb+gI+tdLxb8+Kq3z+SkUq5WWFyqrGWlJ7Ex1oSKR2I6StKczR89YeUTH1KA7i92pI908c7Qb0UICBaDTR/+Ef4mPyvylUCwXt4247lpk+0A0Fm+5UN1WwZ61FzWuleohw+E8vaPlX8XkwNJCJWJqnLevn9NGe1pOoAgBXcFVGE5stNbnqFCsRmqVvj8VOA04IRqLv6GKGepYYPWFJ6wQ9uwSpq+S6sQj4Tw9Y+UOuefODHLew6tz5dwRXH0lWpTpS6Hwhlrmry3A25waX1+v8nrayapXVMGeWdwvMlfzK6tbTAesBeYvPzs++k2Xz6G5SqgoFK1LreTHnwI/jcbiu1LJhK/ViNuFJ+1Z8tIGU6Kny6O+BMW6OsadBbTHLSIFSVDAG7U+AoX5Sc2vkvcrgaZBh4+mL3xLfFydWqFC0Wq4ctSnkomt0Vi8L5VMTFA0i3UCFwKPpZKJB3wf5SrBlpLsXN7JjI2scKhIQH/UxLh7XtgYgFHh2vKr5P1K4LfpC6WpKBQtjVtH/cvmstKjsXjQ6cL4beCOkkTENU++xE8iZqt46G1J4K7abVTtwZVfNQeERrfU6RGBkkf9wYF+m74EoCuholC0LG5njQ87vhWAVwDdwGbgDODfncKNa55cqVCZWShU9EcttCrbSzF7oXDKykd6R9AJSJ2IKK8eW0CSke5yZjQNIj6bvvxqI6zCiRUKb3A7ex3jlMEHeDFwdSqZmI3G4r8EjvJxfKuK3Nzkm5OIQoXwsCWBP5VrKenjBN+7eL7nhy4Fr9i7AU00R1OpRo8IuBYqXSFvSksvhR+mLxVOvPLcOjHRlPc9r6+vKe+7lnA7e+WBQDQW14ALgF8428PNrgDcSsxpKtVMX/pjFtp0ScFIDXadY5PploceXYbRFIECEFjkVgigERHuSm13hvwfu/KnKNYK0VicG276ZbOHUTduNZXfOH1TLEcQ3eJsfyuwVA/7NUXeqeO1wPRlSwJ3l2sp5rE6e9eX90BZl2le/SADwWItvtxoK7oOkRXou+RHiRYlVBS1ePf7LuO++x881MskYBgMDQ1y4QXn84r4S1ydY3hkhDvuvJsXXvR8n0fbXNwKlbc7QmUj8MpUMlGIxuJDwEccc9iax5bykKNepMuFSjUtxXyGwUi4vA/KumxzhMpipq9Dr6MRFvqhyLZqdK7Q0P3wNqlwYoUbnnf+ubz/vX8PTlfGO+/6Ex/950/QEQnz5y/9s5rH//o3t/PLW29TQoX5oo4vqtg2Eo3FN8+Vnl/rHIr8ykiEWSJUqmkpx+hY3YKRcLluMJRtTovVoAsraK8IFIWKKIYNV9KzQo5upakoWgFd1zn92adx2JbN7NpV7LS+a/cevvjlr/PQw49i2TZHHbGdv3nLGznm6KP4+je/zdU//DFSSmIveyWf/uQ/c9yOY7n+Fzfy/at/yPDIKBs2rOM1l17Mhc87/9D7ZDJZPv6JT3P77/9IKBzipS9+Ea999SubeOW1cf0Ljcbi26Ox+Eejsfg3Sjaf7M+wVh85p5pXpT9Ff7xCSxFgnmIwHjKxtPntEVOj02xOm9CAC/f6nLbS1yHY2KMteARXKGDND7GrhIqiXnK5HL+85Tb27d/Pec85B4CP/csnCYfDXPXfX+Hqq77Bxo3rufyfi70N3/SG13HhBedzwvE7SF53NcftOJZ77r2fz33hy/zd297M//zwu/z1G/+ST/7H57jn3vsPvc+P/+en/FnsBSSu+Q6vvvRi/vs73+PJp5rWc9AVrqaCaCx+nuOcv8/pT/9Gp8/JzdFY/LWpZOJH/g+1tck5/hSttNR9lbwU61gd2a0x3FHuTxnKBSGgQWXU2ApQy/w1x7pAAKPB1sVe4YemosKJFW64+ZbbuPW23wJgmiaGYfC6V1/C8ccdC8DnrvgEmhCEQsVyGOefdy6pG3/J2Ng4AwP9C87345/8lNOf9UyefdqpAJx95ul85B//gb6+3kP7nHXm6Tzj5BMBeMHzn8cXvvhVdu7czRHbt63INTeC2/Xly9XZUwAAHtZJREFUvwDvTCUTX4rG4hkck1g0Fn+Nk8OihAo22LJMU9GfWKilzOWgDPeVm74Gu8NYER3j0drJkV6zWORXJUf260wInWmrsRbGXuC1QqTCiRVuKfWpmKbJ0zt38bkvfJmHH3mUj37o/dxzz31c9f1r2L17L9lcjrkGiPlC9RCYPXv38cxnlBt7zjn7zLLnmzdtPPR/MFi0O2dzOc+vzUvcLtFOLCkoWbqUvs7JYVnz5KRVjPoq+XSMByq0lGOKWgrAcKB8xT8kQhASyBWe4DQEmgvzV1+HoCeisS4YJKhpjT+EQBeNF+v3OqRYmb4UjWAYBkcdeQR//ca/5Ne/vZ1du/dw+T9/ghOOP45vX/llktddzT999B9rnqeWXUKswnJNbhd+EugBRiu2HwG0tthcAWwpKSDLtBIxZqMPl98yhWcUP+6CYTMpKjQVUVyFyA6BmFo5E5ibDHVNgy19xf0imsYxNVr5usWWsvioJ9nJ4x+ZroSKwgN++7vfUzBNXv+aS+noKP4+Hnvs8SWP2bJ5E7t27S7blrrxl2xYv+6QyWs14laoJIGvRGPx91D0sfQ6fVY+BfzU3yG2PoeSHkvyU4yHK7SUzRqyt6iljPSYZZNoDwYhJ8FQdgqYWpFhg0vT18YeQdBDDUrTIdSjE+zSCHVrBLs0NN3byX1k/wRDGzs8PadCMYeUkr379vPNb13Ftm1bOfmkohC46+57OOvMZ/Pr397O7/94JwAHDw6zccN6wqEQo6NjTE1NEQwWI7nef9nl3Pab33H2madz95/u5YrPfoGP/9OHm3x1y8OtUHkHcBUwJ3rHnNp+Pwf+3sfxrQpyWMVQ4jknuykxHiv3O5g75iO7hjurmL4cZMfKrpxrOenDAVjfNT+mcK9G//bGk1KEBoZyjK95VmO5lFJHvRCCgf4+zjj9WVx6cZx169bx6ktewac/9wVs2+bcc87inz7yQS776L/w4Y99nMs//AGef8H5/Po3t/P6N/0N//Dud3DOWWfw9r99C1/9xrf4109ewfp1Q/z9O97Gqaes7qBaMedMckM0Fj8O2AHYwCOpZOJhX0e3MjRsa7rudweI9AxywM4ycSCDdrCoseiPmoRunTdvyTBkLg2Dsxq/8ehxntbn03tO1wY4Xu9xdpYYD5nl/YZd0hMBy4aCDQULpIs29+u1EEGKAi8zO02ks7vs9aPXC7odIaAbgg2nhDCCrS8URvbvZWjj5mYPwzNW6/UMDw+zbt26sm2FQp5AoHnVI7ymXa+n2ndXQdUVcF3BNKlk4iHgoWgsfjTQGY3FRSqZWPO1v3LYiOlS01eFlnKMcUigAIzo5W6oUk0FIZARUbXK8VKEAjDUVT7ZmzaYNYTTsWEDzfFTTE/YdPfNn0MT5VWH+48KrAqBolAomkdNoRKNxd8CHA98L5VM3B6Nxa8B/sJZ4d8bjcWfn0omRlZmuK1JLm8iMkUhIMZt9APlKsKc6ctGcrCnwGyJGqIBAxXl5mVH/UKlJ7Jw0WBoxcdiBMW8FgJgBaFzkX4oXRsMOgZWviS/QqFYXSw5S0Rj8Q8AHwIeAd4QjcU/7vRROddRff7Zef0dKzfk1sKSEmvGRkOyuzNH190mW5j3n0xstbnjmClGwnnGQgWsikm+XwTRK/wasrM+v4qmQVcDfUxCLhMJAxFB37bmlJBRKBSri1qzymuBeCqZOBV4s1NA8q2pZOI3qWTi18BfA7UrqbUxeSzEtOTBvllu2jDK4APlk/vt56Z5sH+W4chCgUKl6ctBRkRdjUm6Qo1F2roRKkLA4DEhz6OzFApFe1JrVtkCpJz/rwNCwKGe9Klk4jGncvGaJWtbWLMmdw5NcdS9QcKZ+Y8002HzxEmLlzXRgB1ad5UXin4Vt/Q2WGU37EKo9B0eINip/CgKhcIdtYzkwVQyYVMUINloLJ5LJROVrt81vYTNzVg80p0mr0tOvL1c63joWTmsCqtRBJ0hEWRQhDhS66RLVDcryQ6xoIR+NcJBCNRZh1ILg9YhGOgJ0BnU0QyBpgvEiKB/3XwUi9Cgc0j5URQKhXvUjLFMstMF7u+foe+gxpYnKxzuxxo8c6QbTQr6cwaDRojQ1oXmrmq4zVfpqTcjXEDoGAM9rLGpJ3wo8gsgKwVd69UtoVAoGqemphKNxb+1xHN8qka+anjEmmEmYHHRDZ1l262NGsfRU1bYxl6v4SJ1BFwKFV2Dzjod9MaQhhYSRDStTKAoFAqFF9QSKrcBW0ue/6ri+dw+a5JC1uKe3inOSEU49u5yDaQ0g36OurLlDYEMC0R2cRNYd7hOB70GgU1F/0in3pzeLQqFor1ZUqikkonnrtxQVh+P5rNsfEDjjBvKCyzaPQJre8WkLajL+c6cX2UxoSKgu04HfWCdhuZoNp0+9CVRKNyy87fNaRh7+Nmrrx7cN//7Km68+Ra+c+VXmj0UVygD+jJ4eleW8xPlN6kMQe6iIFQUYJQhUZZV7wbZIYpV1qrQ4fT0co0GgY3zByhNRaFwz7vfdxn33f8gesnvJmAYbNmymT9/aYyLnn9BU8fXSiih0iB33jXJaddG0OS8oLB1yEeDh6oRl9JIociljuleJPN9MQIbNERg/pgOJVQUiroobdKF01L41tt+w79f8Z9EwhGec+7ZTR1fq6CESgNM7i5w5+XjBAolAkVICs8LYW+oPlk3VH04KJABMV/92CEShI46hIrQi0JlDl0UHfUKhaJxQqEQ0Qufx4033cItv/o1zzn3bB548CG++o1v8eRTTyOE4PjjdvD//uav2eR0cLz0dW/ilX/x5+zctZtbf/UbEPDc887l/73tzQghsG2bb33n+/w8dRPpTIb/c/YZ9PT0lL3vgQMH+eJXv8GDDz5CwSxw5BHb+es3vp4dxxb7Jb7j3R/ghON3kM/nueHGXxIMBXnjX76WrYdt5rP/+SX27z/I0UcfyQff9y7WDQ15/rm09MwSjcUvjsbi/xONxfdEY/GZaCx+RzQWf9Ui+54QjcVvjMbi6Wgsvjcai38sGot7vhzPjFtc/4H96DPl2w88V2BtW/zt6vWnHDquQhgFArCh202vxnmMDRrCUFqKQuEHhUIBIQSFQoEPXf5xjj7qSK7+7pV855tfxjItPnXF5w/ta+g61/zwOs4849lc870ruez97+EnP72e3//hDgBuvPlWrvnRdbzvPe/g2u9dydlnncH1v7jh0PGWZfEP//hRDN3gK//1Ga668iscv+NY3v2+yzh4cLj4HobODTfezMknnsAPf/At/vwlL+a/vvx1fvTj/+VT//ZPXPWtr5JOZ7jmh9f58nm0tFAB3gXMOD1bXgrcDHw3Gvv/7Z15mBTVtcB/PTM9C8OwzABhEQUFUUyeihpQ44pX5QYl18SYfDFxSYzkKSoYYxITt6iJikseGjce6jMuCYn1iFCo5RpFXzY1ibK5IGhgMkBQZBkYpvv9MaeHmqZ7prvpganJ+X3ffN19l6p7qmruqXvuvee4KeFCxrq+wDPi5HIScC1wKXBNsRu0dP4nbFrVdv/n347ZQq+9K7PWSZbFIE9zVWvdkB+w0lIY2KuEfAYZsbK2oxR0kl5RisLmzZv53dz5vLlwESeecBzxeJyHZt3Dt795FuXxONXV1Rx5xFgWL3m7Tb3R+4/i8LGHUVpaypiDD6SmpifvL/8AgBd+/zKHHXIwBx34GeLxOEcdeTgHjN6/te6f/vwaK1eu4jvnn0vv3r2oqqriG2d+hbKyMl74/faFuIMHD+LYYz5HWVkZnztyHI2NjZw6cQJ9evemV00Nh445iA8//EenXJeubv46Jc0D8nPGusGibGaE0icDVcBpge+tBwJjXS/gamPdTZJWFIafXs2sFfUc9EyLEll8yBaax8SJrc+uNAodpSBKJVkTIxaDwbUlVGYZZPSJx9mjvHyHEUxVXSnVfdveZt2foij5Ew7ShYxQRo7Yhyu+fymHHToGgJcWvIo3Zy6r6v9JU1MTiUSC5ua2L6GDBw9q87uivIItW1rCYTQ0rOaQgw9skz90yGCWL18BwMpV9fTsWU1dbW1rfjweZ+CnPsXKVfWtaZ8aMGD78SsqJK1/m7QtW7O7kNoZurRSyeJS/3VxvR9mAvBUmvJ4DLgROAZ4olhtemH9el4+cRPrejUz/K04r07azBdX9G63zk5Fc6yI0TysjL3rYvTOopyGVlYyvDL7SElRlJ0nPFGfSCSYetkV1PSs5vCxhwHw5luLmH7bDCafdw4T7UlUVFQwb/7T3D7jrjbHae+lrqmpaYfNZ82JtlumMwVWTCaTxEL1Skp2PEdsF1kourRSycIR4oo/zH7Ac+GEwPdWGOs2SV5WpfKv1fUkmnPd5w6HJ5NMr6tl9jEbmfvZDRzYUEVzY2PWQI2JCmiqBDbmfIodGFyTJLYlyfq2sb2IAcPL4tQ0bqIYAW22bWtiTf3KIhxp99OdZCHC8jRu3UZTU9obcTJJ87Ztu6U9O7QlR5LJJIlEok39S6acz5RLLsd/MsCebHjzrYXU1PTk1Iknt55rydK325w3CTQnmtscJ0myNa1fvzpW1de3yV++4gOScoyBAwewceMmVtXX06+uZbSydetW6v/5T8YfdzRNTVt3aOu20GcqrTnRTDKZzHw9JL1x8+Z2n7lskUgjpVSMdeNlzuTctKy+wEcZqqyTvKzU9s/fyfJx8nfjgmWUNVZRkWW+JFkWo3nvUuJprlQqy2Fon5LcHEHGoCJDubJYjP2rq+lbVrxbGNWQtZnoTrIQYXlWr169Q6jdpqatlBbxuc2HQsP+xmIxSkpK2tQfPmwY55z1Ne6b9RCHjDmYQYMGsnHjJt5btoJ99h7GvPlP876Yrdat+5gBA/oTA0pLStscJ0asNe2IcZ9l1oO/ZPGSd9hv1EheefWPvP3Ou1T36EE8Xs7Yww5l6B5DmDnrIS6eMpnSkhL+5+FfAzHGH38M8Xj5Dm0tC32m0kpLSonFYhmvRyqccGVVFf3aDyeckV16Z411vYFBHZWTsMXpdYcBjwBzAt97IEO1TFvPYzsTg74jTh1YwR96lbNh0TZizWmniUFijxIIK5QYDKiBQb1KyDA6pW88Tk2OK7MGxOO6ikuJLFHc2Z4JN2kiLy14lRtvvp1bbr6eE084nu/98CrK42WceMLx/OSqK7j08h8xeco0Ztx2U4fHm3SKZc3atVz30+k0NjYybuxhnDbpFOY9+TQApaWlXHvlD/jFvbP42lnfprw8zr4jR3D79BvazLPsTmKZ7HOdhbHuW8B9HZULfK9Nl2usqwUWyEqwYwPf25iW3wDcGfjeNWnpG4BrAt+7uZ3TFXwB1tSvpKSuPy+v+ISGvzaxJWSeah5cSrJ2uw2zogz2rC2hZxYnxXtVVrLXbp4XierbcCa6kyxEWJ7Vq1fTP+1tN/Um3F3orvJkundpZDTR7NKRSuB7M4GZ+dQx1vUA5gLlwOfTFYqwWOZOwvWGAtWS12nUxuPsO6iSnqWwbtE2NjQmobaE2NDtCqU8tRQ4wy0oicUY1aMH/eP/1s6eFUXpJnTpORVjXRkwGxgJHBn4XkOWovOBy4x1NYHvfSJpZwCbgRc7u537VFayrnYbyRHQc22Cyn1LiWXSIGlUlJRwQHU1PdWMpShKN6FLKxXgF4AFLgZqjXXjQnmvB76XMjjdDVwEPG6suxHYG7gauLWYe1SyES8pYURlJYv6b6KsLpaTQukbj7NfVRVx3YioKEo3oqsrlRPl8+cZ8oYD79NiVlsnK8PukOXDHwG3iWLZJfQvL6ehqYm1TU1Zy5TEYgyIxxlcUaGjE0VRuiVdWqkEvjcsj7ILgd3qf3pEVRUbmpt32LMSj8UYWF7OwPJy4rqbXVGUbkyXVipRo6KkhLFpHkUV5d+dRCJBiZp5I0Uikci4cz8X9E4ritJp9OnTh4aGBhKJ3L1WKLuXRCJBQ0MDffu2u288KzpSURSl04jH49TV1bFmzZpW31SNmzdTWVXVYd2o0N3k2bxxA4OG7EG8wG0OqlQURelU4vE4A0Jec9fUryzI/UdXpdvJ09xUsEJBzV+KoihKMVGloiiKohQNVSqKoihK0VCloiiKohQNVSqKoihK0VCloiiKohSNXRpPRVEURene6EhFURRFKRqqVBRFUZSioUpFURRFKRqqVBRFUZSioUpFURRFKRqqVBRFUZSioUpFURRFKRrq+r4AjHWjgRnA4cBHwEzgmsD30iMJdymMdSOAy4BxwKeBlwLfOzZDuS4vn7HudODrwCFAb2AJMD3wvUfTynV5WWhp55eAacAooBpYDjwE3BT43tZQuUjIk8JYN0TuTTVQE/jehlBeJGQx1p0N3J8h6zuB790dKhcJeWhpaxnwXeCbwJ7AamB24HtTQ2UKkkeVSp4Y6/oCzwALgUnAPsAtMur70e5uXwccAFjg/4DyTAUiJN80YBkwFVgjcj1irOsX+N4MoiULQB3wPHCz/AN/FrgaGAhcSPTkSXEzsEGUSisRleV4YHPo93upLxGU535gPHANsBgYCoxOZe6MPKpU8mcyUAWcFvjeeiAw1vUCrjbW3SRpXZUnAt+bQ8tD8xugX4YyUZHvlMD31oR+P2esGyzKZoakRUUWAt+7Jy3peWnrBca6KYHvJaMkDy3P2FHAycANolzCREoW4U/hkVYakZHHWHcy8BXgwMD3FmYpVrA8OqeSPxOAp9Iu6mNyA47Zje3qkMD3cgkUHgn50hRKiteBwaHfkZClHdamjSgjI4+xrlSU+7UykkwnMrLkSJTkORd4rh2Fws7Io0olf/aT4WIrge+tADZJXtSJsnxHAEtDvyMni7Gu1FjXw1j3OeAi4C4ZpRAxeSYDlcCdWfKjJEuKd41124x1S4x156flRUmescBSY90dxrr1xrpNxrrHZaSfomB5VKnkT1+xeaezTvKiTiTlM9aNF9vvLaHkKMqyUf5eAl6UhRUpIiGPsa4O+AkwLfC9pizFIiGLsAr4sSwMOQX4A3C3sW5qqEyU5BkInA0cJGawc2TBi2esi0mZguXROZXCyOTaOZYlPYpESj5j3TDgEWBO4HsPpGVHShYZbfWQiforgTuA/wzlR0Ge64E/BL7nd1AuCrIQ+N5TwFOhpPnGugrgR8a6n4fMypGQR9oUAyYFvreWlv+hVfISczzwrJQrSB5VKvmzDuiTIb13Fs0eNSIln7GuFpgPrADOTMuOlCy0dGCvydeXjXVrgAeNdbcEvvduFOQx1h0gNvujjXWptvaQz97GuubA9zZHQZYO+A3wZWCYrAKLkjzrgPdSCkV4GdgqK8Ce3Rl5VKnkz+J0m6KxbqgsmVycvVpkiIx8xroewFyZzP584Hsb04pERpYspBTMcODdiMgzEogDr2bI+xD4b+BbEZElF1Jv7VGSZxFQkSE9BqRGXQXLo3Mq+TMfOMlYVxNKO0PWr7+4G9tVLCIhn2zemi2d2ITA9xoyFIuELO1wpHwuk88oyPMycFza342SZ0NLi6MgS3t8UVa1LZffUZJnLvAfxrrwloKj5WXgr/K7YHk08mOeyKaghcCb8s+yN3ArcHvge11xk1Mr8mZv5eelQC/gKvntB763KSryGevuBc4DLgb+mJb9euB7W6IiCy3yPCmbzd4CmkWhXArMDXzvK0T42QvtSG/dUR8lWYx1v5Vn7G9AqXSuZwIXpW20jYo8vaSd/5A9RDXS5sWB7xl2Uh4dqeRJ4HvrZCdqKfCE7Ei9LdQ5d2UGyNv9bHHVMjr0ewDRku9E+fy5mFrCf4OIliwAf5IVObOBX8sqox/IiiOInjztEjFZlsg80W/l/owGvpFSKERMHtl7crzMmzwmy76flTmiVJmC5dGRiqIoilI0dKSiKIqiFA1VKoqiKErRUKWiKIqiFA1VKoqiKErRUKWiKIqiFA1VKoqiKErRUDctSpfHWHcs8Hzge7Ecind2W0bJfoV9gC/m4DSxS2CsawQmZ3C42W0w1h0NPA2MDnzvvRyqKJ2AKhUlJ4x1L4grh6MC31uQlvcALRumzt5tDdx1nCdO9fpJbIk2hHaPb8lQ96PA9wbuikYa684EXkl1roHvVe6Cc5YBFwBfE79RMeAD4HHg5sD3Pu7M8we+93uJ4aLsRlSpKPmwBrjXWHdw4Htbd3djCkViRpQEvtdcQPU6YHkG55Xp9Ax8b1uBTdwpRL7bgK+G46h38jnjgA8MAaYCLwBNwGfEFf4bxrojAt9b1Vnnbyd2i7ILUaWi5MNMCYT1A3HbsAMS22QZYALfe0bSRgBvA8cFvveCsW6BuFOpBL4hTup+KO4w7havvK8BXw1878PQsU8QtyzDJMLjJYHvvSh5VcBNwOniSfVt4PrA934r+VcDE4FAIipOEl9b6e0/GvgpcIAEy3pRgk3VG+t84CQgJuak0wPfeyLfiyht+Vbge3uE0q4Dzgx8b5ixbqTId5IE6Toc+BdwVeB797O9E79B3Lj0lFDK08Rf07/Ec7NvrHsy8L1TjXVJ4LzA92ZK/fOBC+VarxLXMNeKzzQjZqSxEhL40+Inamrge/OyiHWx+CvbTyIEpvirsW4S8AbwC8DJ+dPbUyZK6JyUic5Yd4EoqCES2uCXwM8C32tKmUSBsyQw293GumclbWTge+/k8EzUiXzj5Rp+ANwS+N59+d5TZTs6Ua/kw1Yx/3zfWLf/ThynSZTJS/LmP0OUxSXi1XaoOLv8blq9KdLRDpC680KeVm+V6HXjJDLddcBjEpY3xV7yzPcFnktvlCi/Z8VcM1iO1x8IjHWxwPcs8JCYlSoLUSh5XB8keuJUuRYPAHdJ/BhEsX9BTJJ9xTvwU+KraZSUsYHvnZpBzrOB6eKwsq90ul+Xaxg+/9XiD6oPsAC4PxQZMJ0zgUfTFAq0mKWa5XwTxVFhhxjrzhA/UymleYa4zL88regp4qk6k0+qjp6J68WMOUqcKk4BbjPWjc6ljUpmdKSi5EXge68Y6+4H7jPWHRWKn54v7wS+9ytaOpDHpQO9M/C91ZL2lJhOwlyfGrkY666SN+2TjHW/E2V3UuB770vZxyX929LhAtTKMbKZ7iYD7wa+lwpJvMlYd4WEjz1UnD7uSh4MfO9NWuSdLSFtR4jH3Ckiy1K2j3T+Ju7LO+JC4OHA956W328Y6+6QSIZTQuXuCHxvOdvv0dmiZDOFGRglI4lsLJT+Zrg4MuyIacA9ge+l4rK8YaybLiO360LlHgx87yNpY2uiuGzv6JkYKPFDGiV6Y2Csq9mJZ1pRpaIUyPcl0M9k4K4Cj/F+6Htqwnt5WlpVWp1FqS+B760z1n0E7CkrsUrF3BPuEErSgkWtFQ+t2RgRPoewRD73yVOpbAh3csKCwPfG53GMd0LfN8tnD2NdaqFA63yJzPE8yvYOtT1GyIgrzBJ5Wx/Q0fmzHHNbB/1JqXzm2ueMAsYY6y4LpcVoka88lJZtziiXZ+LHwByg3lj3HPCkeO1t7xlROkCVipI3ge+tN9ZdCDwgb34dkcnMmsgxLUz6xHoMaAzVOzwUjjcTuSwuSDfvpNqe79trvhP1uV4jQm0sdIl1Mkc5O7ofYd6RuZds7C/HW5olP13+BHB54Hu3ZiocUtjZ7mmHz0Tge3+X+asjJJTCVOBKY9248Fyekh86p6IUROB7nsxL3JGWlXqjDb9N7lmk07bO48gkax+ZwH1PFM7B4cLGuj1lAjhXlkisjEznzNYZFsLmtOtDPtdIzD1rwuFejXXlxrpLjXXDczjEUlmIEGZ/iT2eybSVC48CXzbW7ZueYawrDQUcS8U3b+zgGVma4X4OMNb1zLE9HT4Txro+sgrwpcD3fixKcRPwpRzPoWRARyrKznCBmIs+Ca2kagDWAhPE9FAt5YrBFTJC+hi4VmzzTwa+t9FYN1PmBF6TuYUjAU/O/ViOx78bmGqs+64ELqqTFVZ/DHzv9SLJgFyz/sa6QwPf+7OxbgxwgnS0uXIncIGsSFsMfE/mIWaF/q/3M9b9RQIuhfkv4B6Zp3keOEiu032B7yUzmO1y4XbgNGC+se4ieR7CS4prgc+nXYOTjXV3iZnqcjGhhY/3oLHuf4HfyeKNXwF/EbNruwS+t6G9Z8JY9ysxZz5hrLtGTF6fEbPi24VcAKUFHakoBRP43kqZXxkSSksC3wQmGOvel9VUqQh5hb7ExOXt/j4ZHa2WFT0TQ/tFpgHzxC6+UcpeGfhergqFwPeWyVLjr8qy3FeBd0MhmIvFXOBeYK6xbrGYXW7N8/rcAPxGOu91osQnBL63ThY7PCzx4OekVwx872FZLXWPdKa/lBHnDwsVSBY/HAc8IntkPhalMg/4EDgkNekvTJFJ+w9kIcQcGSmVyfEek0n5nwEbZGL9FVkhmCtZnwl5Tr8gI9Nl8mL0qORnWzat5IBGflQUpeiIO5vFwKmduPRa6YKoUlEUpVMw1s2TjaoTgQ91x/u/B2r+UhSlszhXJsz/nrZcXOnG6EhFURRFKRo6UlEURVGKhioVRVEUpWioUlEURVGKhioVRVEUpWioUlEURVGKxv8D0gnTLbtXA7YAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "bento_obj_id": "140024665028592",
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(mtgp_results.mean(0))\n",
    "plt.fill_between(\n",
    "    torch.arange(n_steps * batch_size), \n",
    "    mtgp_results.mean(0) - 2. * mtgp_results.std(0) / (n_trials ** 0.5),\n",
    "    mtgp_results.mean(0) + 2. * mtgp_results.std(0) / (n_trials ** 0.5),\n",
    "    alpha = 0.3, label = \"MTGP\",\n",
    ")\n",
    "\n",
    "plt.plot(batch_results.mean(0))\n",
    "plt.fill_between(\n",
    "    torch.arange(n_steps * batch_size), \n",
    "    batch_results.mean(0) - 2. * batch_results.std(0) / (n_trials ** 0.5),\n",
    "    batch_results.mean(0) + 2. * batch_results.std(0) / (n_trials ** 0.5),\n",
    "    alpha = 0.3, label = \"Batch\"\n",
    ")\n",
    "\n",
    "plt.plot(random_results.mean(0))\n",
    "plt.fill_between(\n",
    "    torch.arange(n_steps * batch_size), \n",
    "    random_results.mean(0) - 2. * random_results.std(0) / (n_trials ** 0.5),\n",
    "    random_results.mean(0) + 2. * random_results.std(0) / (n_trials ** 0.5),\n",
    "    alpha = 0.3, label = \"Random\"\n",
    ")\n",
    "\n",
    "plt.legend(loc = \"lower right\", fontsize = 15)\n",
    "plt.xlabel(\"Number of Function Queries\")\n",
    "plt.ylabel(\"Best Objective Achieved\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anp_metadata": {
   "path": "notebooks/fbsource/fbcode/pytorch/botorch/tutorials/composite_mtbo.ipynb"
  },
  "bento_stylesheets": {
   "bento/extensions/flow/main.css": true,
   "bento/extensions/kernel_selector/main.css": true,
   "bento/extensions/kernel_ui/main.css": true,
   "bento/extensions/new_kernel/main.css": true,
   "bento/extensions/system_usage/main.css": true,
   "bento/extensions/theme/main.css": true
  },
  "disseminate_notebook_info": {
   "backup_notebook_id": "2967299390260046"
  },
  "kernelspec": {
   "display_name": "ae",
   "language": "python",
   "name": "bento_kernel_ae"
  },
  "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.8.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}