{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "879b6016",
   "metadata": {},
   "source": [
    "# Factor Model as Constraint\n",
    "\n",
    "In this notebook, we explore how to formulate a portfolio optimization problem where the risk estimator is not given explicitly but is instead represented by a factor model, as is common in US equity models [1]."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "92d02774",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-01-31T10:04:06.254231Z",
     "iopub.status.busy": "2025-01-31T10:04:06.253981Z",
     "iopub.status.idle": "2025-01-31T10:04:07.022972Z",
     "shell.execute_reply": "2025-01-31T10:04:07.022230Z"
    },
    "nbsphinx": "hidden"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Requirement already satisfied: numpy in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (2.2.2)\r\n",
      "Requirement already satisfied: scipy in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (1.15.1)\r\n",
      "Requirement already satisfied: gurobipy in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (11.0.3)\r\n",
      "Requirement already satisfied: pandas in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (2.2.3)\r\n",
      "Requirement already satisfied: matplotlib in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (3.10.0)\r\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Requirement already satisfied: python-dateutil>=2.8.2 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from pandas) (2.9.0.post0)\r\n",
      "Requirement already satisfied: pytz>=2020.1 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from pandas) (2025.1)\r\n",
      "Requirement already satisfied: tzdata>=2022.7 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from pandas) (2025.1)\r\n",
      "Requirement already satisfied: contourpy>=1.0.1 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from matplotlib) (1.3.1)\r\n",
      "Requirement already satisfied: cycler>=0.10 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from matplotlib) (0.12.1)\r\n",
      "Requirement already satisfied: fonttools>=4.22.0 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from matplotlib) (4.55.8)\r\n",
      "Requirement already satisfied: kiwisolver>=1.3.1 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from matplotlib) (1.4.8)\r\n",
      "Requirement already satisfied: packaging>=20.0 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from matplotlib) (24.2)\r\n",
      "Requirement already satisfied: pillow>=8 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from matplotlib) (11.1.0)\r\n",
      "Requirement already satisfied: pyparsing>=2.3.1 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from matplotlib) (3.2.1)\r\n",
      "Requirement already satisfied: six>=1.5 in /opt/hostedtoolcache/Python/3.11.11/x64/lib/python3.11/site-packages (from python-dateutil>=2.8.2->pandas) (1.17.0)\r\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Note: you may need to restart the kernel to use updated packages.\n"
     ]
    }
   ],
   "source": [
    "# Install dependencies\n",
    "%pip install numpy scipy gurobipy pandas matplotlib"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "7bea3d05",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-01-31T10:04:07.025206Z",
     "iopub.status.busy": "2025-01-31T10:04:07.024977Z",
     "iopub.status.idle": "2025-01-31T10:04:07.652479Z",
     "shell.execute_reply": "2025-01-31T10:04:07.651787Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import gurobipy as gp\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "fb3cd77a",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-01-31T10:04:07.654719Z",
     "iopub.status.busy": "2025-01-31T10:04:07.654468Z",
     "iopub.status.idle": "2025-01-31T10:04:07.662706Z",
     "shell.execute_reply": "2025-01-31T10:04:07.662149Z"
    },
    "nbsphinx": "hidden"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Set parameter WLSAccessID\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Set parameter WLSSecret\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Set parameter LicenseID to value 2443533\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WLS license 2443533 - registered to Gurobi GmbH\n"
     ]
    }
   ],
   "source": [
    "# Hidden cell to avoid licensing messages\n",
    "# when docs are generated.\n",
    "with gp.Model():\n",
    "    pass"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "96d8f302",
   "metadata": {},
   "source": [
    "## Background\n",
    "\n",
    "The *standard mean-variance (Markowitz) portfolio selection model* determines optimal investments that balance risk and expected return. In this notebook, we maximize the expected return of $n$ assets subject to a limit on the variance, whose standard formulation reads:\n",
    "\n",
    "\\begin{alignat*}{2}\n",
    "\\tag{1}\n",
    "\\max_{x \\in \\mathbf{R}^{n}} \\quad        & \\mu^\\top x \\\\\n",
    "\\mbox{s.t.} \\quad & 1_n^\\top x = 1 \\\\\n",
    "& x^\\top\\Sigma x \\leq \\bar{\\sigma}^2 \\\\\n",
    "& 0 \\le x \\le 1\\\\\n",
    "\\end{alignat*}\n",
    "\n",
    "Here $\\mu \\in \\mathbf{R}^n$ is the return estimator, $\\Sigma \\in \\mathbf{R}^{n,n}$ is used as an estimator for the variance, and $1_n$ is the vector of $n$ ones. While it is certainly possible that the matrix $\\Sigma$ is given explicitly (e.g., as the covariance matrix of a time series), it is often expressed implicitly through a *factor model*.  In that case, the matrix $\\Sigma$ takes the form\n",
    "\n",
    "$$\n",
    "\\Sigma = X \\Sigma_0 X^T + D\n",
    "$$\n",
    "\n",
    "where\n",
    "\n",
    "- $X \\in \\mathbf{R}^{n,k}$ is the factor exposure matrix,\n",
    "- $\\Sigma_0 \\in \\mathbf{R}^{k,k}$ is the symmetric positive definite (SPD) factor covariance matrix, and\n",
    "- $D \\in \\mathbf{R}^{n,n}$ is the *diagonal* matrix of specific risks with $d_{ii} > 0$ for all $i$.\n",
    "\n",
    "Effectively, this splits the risk into two sources:  One that arises from common risk factors (macroeconomic conditions, market trends, etc.) and a specific risk that is uncorrelated among the assets.  The number of factors, $k$, is typically much smaller than the number $n$ of assets.\n",
    "\n",
    "The important observation from a computational point of view is the following: The factor model data is of size $nk + n$ (with $k \\ll n$), while the covariance matrix $\\Sigma$ is on the order of $n^2$. As we will see, optimization models using the factor model admit a smaller representation and typically offer improved computational performance."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "204fa209",
   "metadata": {},
   "source": [
    "## A synthetic factor model\n",
    "\n",
    "Leading industry factor models are commercial products that are not necessary for our demonstration purposes. We are more interested in a qualitative comparison of solver performance using *some* factor model than in the actual economic meaning.  For this reason, our first step will be to create a simple, synthetic multivariate factor model that we can use in the study that follows.\n",
    "\n",
    "The following function uses an uncorrelated factor covariance matrix $\\Sigma_0$ that is used for sampling a multivariate normal distribution on `num_factors` factors along `time_steps` sampling points for `num_assets` assets:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "a9d07248",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-01-31T10:04:07.664571Z",
     "iopub.status.busy": "2025-01-31T10:04:07.664395Z",
     "iopub.status.idle": "2025-01-31T10:04:07.668690Z",
     "shell.execute_reply": "2025-01-31T10:04:07.668290Z"
    }
   },
   "outputs": [],
   "source": [
    "def factor_model(num_assets, num_factors, timesteps):\n",
    "    # Generate random factor model, risk is X * sigma0 * X.T + cov(u)\n",
    "    sigma0 = np.diag(1 + np.arange(num_factors) + np.random.rand(num_factors))\n",
    "    X = np.random.normal(size=(num_assets, num_factors))\n",
    "    alpha = np.random.normal(loc=1, size=(num_assets, 1))\n",
    "    u = np.random.multivariate_normal(\n",
    "        np.zeros(num_assets), np.eye(num_assets), timesteps\n",
    "    ).T\n",
    "    d = np.diag(np.cov(u))  # NOTE: This is the _diagonal_ of D!\n",
    "\n",
    "    # Time series in factor space\n",
    "    TS_factor = np.random.multivariate_normal(\n",
    "        np.zeros(num_factors), sigma0, timesteps\n",
    "    ).T\n",
    "\n",
    "    # Estimate mu from time series in full space\n",
    "    mu = np.mean(alpha + X @ TS_factor + u, axis=1)\n",
    "\n",
    "    return X, sigma0, d, mu"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "70fdf914",
   "metadata": {},
   "source": [
    "We skip the details of this statistical procedure for constructing $X$, $\\Sigma_0$, $D$, and $\\mu$; more background and details can be found in [2, Sect. 6.6].  The essential point is that we now have a synthetic factor model that acts similarly to a commercial factor model regarding solver performance."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9c7749c0",
   "metadata": {},
   "source": [
    "## Taking advantage of the factor model structure\n",
    "\n",
    "Ignoring the constant factor $\\frac{1}{2}\\gamma$ for a moment, the objective function we want to maximize includes\n",
    "\n",
    "$$\n",
    "x^T \\Sigma x = x^T (X \\Sigma_0 X + D) x = x^T X \\Sigma_0 X^T x + x^T D x\n",
    "$$\n",
    "\n",
    "The first term, $x^T X \\Sigma_0 X^T x$, would result in a quadratic function having $(n+1)\\frac n2$ terms.  Since $\\Sigma_0$ is SPD, it admits a *Cholesky factorization* $\\Sigma_0 = LL^T$ where $L \\in \\mathbf{R}^{k,k}$ is a triangular matrix.  This allows us to rewrite the first term as\n",
    "\n",
    "$$\n",
    "x^T X \\Sigma_0 X^T x = x^T \\underbrace{X L}_{B} \\underbrace{L^T X^T}_{B^T} x = \\underbrace{x^T B}_{y^T} \\underbrace{B^T x}_{y} = y^T y\n",
    "$$\n",
    "\n",
    "where we have substituted $y = B^T x$ in the last step. Also, note that $B = X L \\in \\mathbf{R}^{n, k}$ comprises only $nk$ elements.\n",
    "\n",
    "The second term, $x^T D x$, comprises only $n$ terms (i.e., $\\sum_i d_{ii} x_{i}^2$) and can be used as is.  Putting this together, the standard optimization model (1) can be rewritten as\n",
    "\n",
    "\\begin{alignat}{2}\n",
    "\\tag{2}\n",
    "\\max_{x \\in \\mathbf{R}^{n}, y \\in \\mathbf{R}^k} \\quad        & \\mu^\\top x \\\\\n",
    "\\mbox{s.t.} \\quad & 1_n^\\top x = 1 \\\\\n",
    "& y = B^T x\\\\\n",
    "& y^T y + \\sum_{i=1}^n d_{ii} x_i^2 \\leq \\bar{\\sigma}^2\\\\\n",
    "& 0 \\le x \\le 1.\n",
    "\\end{alignat}\n",
    "\n",
    "Note that the $k$ variables $y$ do not have bound constraints.  Model (2) contains the much smaller matrix $B$ instead of $\\Sigma$ in model (1), at the expense of $k$ additional optimization variables.  Generally, form (2) is advantageous for solver performance; we shall compare both in the next step."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "de62456d",
   "metadata": {},
   "source": [
    "## Comparing the two optimization models\n",
    "\n",
    "In order to compare the solution times for models (1) and (2) above, we will define two auxiliary functions to build an optimization model in either form.  Let's start with a function for the traditional form (1):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "6cef3785",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-01-31T10:04:07.670493Z",
     "iopub.status.busy": "2025-01-31T10:04:07.670317Z",
     "iopub.status.idle": "2025-01-31T10:04:07.673687Z",
     "shell.execute_reply": "2025-01-31T10:04:07.673287Z"
    }
   },
   "outputs": [],
   "source": [
    "def build_sigma_model(model, v_max, sigma, mu):\n",
    "    x = model.addMVar(mu.shape, lb=0.0, ub=1.0)\n",
    "    model.addConstr(x.sum() == 1)\n",
    "    model.addConstr(x @ sigma @ x <= v_max)\n",
    "    model.setObjective(mu @ x, gp.GRB.MAXIMIZE)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ff6c5f4a",
   "metadata": {},
   "source": [
    "The second function build the equivalent model (2):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "ea68c2f5",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-01-31T10:04:07.675466Z",
     "iopub.status.busy": "2025-01-31T10:04:07.675296Z",
     "iopub.status.idle": "2025-01-31T10:04:07.679101Z",
     "shell.execute_reply": "2025-01-31T10:04:07.678721Z"
    }
   },
   "outputs": [],
   "source": [
    "def build_factor_model(model, v_max, B, d, mu):\n",
    "    x = model.addMVar(B.shape[0], lb=0.0, ub=1.0)  # n variables\n",
    "    y = model.addMVar(B.shape[1], lb=-float(\"inf\"), ub=float(\"inf\"))  # k variables\n",
    "    model.addConstr(x.sum() == 1)\n",
    "    model.addConstr(B.T @ x - y == 0)\n",
    "    model.addConstr(y @ y + (d * x**2).sum() <= v_max)\n",
    "    model.setObjective(x @ mu, gp.GRB.MAXIMIZE)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7376cd78",
   "metadata": {},
   "source": [
    "We are now ready to run a small benchmark.  We will keep a fixed risk aversion coefficient $\\gamma$, and solve both models over a range of data with an increasing number of assets:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "a619c048",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-01-31T10:04:07.680944Z",
     "iopub.status.busy": "2025-01-31T10:04:07.680775Z",
     "iopub.status.idle": "2025-01-31T10:04:21.946536Z",
     "shell.execute_reply": "2025-01-31T10:04:21.945869Z"
    }
   },
   "outputs": [],
   "source": [
    "np.random.seed(0xACAC)  # Fix seed for reproducibility\n",
    "\n",
    "num_factors = 72  # USE4 has 72 factors, too (see [1])\n",
    "timesteps = 700\n",
    "v_max = 4\n",
    "\n",
    "problem_dims = np.linspace(100, 1000, 16, endpoint=True, dtype=int)\n",
    "time_sigma = []\n",
    "time_factor = []\n",
    "\n",
    "for num_assets in problem_dims:\n",
    "    X, sigma0, d, mu = factor_model(num_assets, num_factors, timesteps)\n",
    "\n",
    "    with gp.Model() as m:\n",
    "        m.Params.OutputFlag = 0\n",
    "        sigma = X @ sigma0 @ X.T + np.diag(d)\n",
    "        build_sigma_model(m, v_max, sigma, mu)\n",
    "        m.optimize()\n",
    "        time_sigma.append(m.RunTime)\n",
    "\n",
    "    with gp.Model() as m:\n",
    "        m.Params.OutputFlag = 0\n",
    "        L = np.linalg.cholesky(sigma0)\n",
    "        B = X @ L\n",
    "        build_factor_model(m, v_max, B, d, mu)\n",
    "        m.optimize()\n",
    "        time_factor.append(m.RunTime)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "f2e53d36",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2025-01-31T10:04:21.948697Z",
     "iopub.status.busy": "2025-01-31T10:04:21.948427Z",
     "iopub.status.idle": "2025-01-31T10:04:22.075447Z",
     "shell.execute_reply": "2025-01-31T10:04:22.074786Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAYAAABB4NqyAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlHJYcgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAbKRJREFUeJzt3XtcU/X/B/DXmFxFrsodBS95xXsa3k2U1Lxnaqaoqd1MjdKkC95KTfOeZWnevr9M01ArTSUSr3iX8n7FUAS8IgIKYzu/P06bTAaMsXE29no+HnuwnfPh7L19gL35XGWCIAggIiIisiI2UgdAREREVN6YABEREZHVYQJEREREVocJEBEREVkdJkBERERkdZgAERERkdVhAkRERERWp5LUAZgjlUqFW7duoUqVKpDJZFKHQ0RERHoQBAGPHj2Cn58fbGyKb+NhAqTDrVu3EBgYKHUYREREZIAbN24gICCg2DJMgHSoUqUKAPENdHFxkTga86RQKLB7925069YNtra2Uodj9Vgf5oX1YV5YH+bFlPWRmZmJwMBAzed4cZgA6aDu9nJxcWECVASFQgEnJye4uLjwD4oZYH2YF9aHeWF9mJfyqA99hq9wEDQRERFZHSZAREREZHWYABEREZHV4RigMlAqlVAoFFKHIQmFQoFKlSrhyZMnUCqVJnkOW1tbyOVyk1ybiIisGxMgAwiCgLS0NGRkZEgdimQEQYCPjw9u3Lhh0rWS3Nzc4OPjw/WYiIjIqJgAGUCd/Hh5ecHJyckqP5xVKhWysrLg7Oxc4mJThhAEATk5Obh9+zYAwNfX1+jPQURE1osJUCkplUpN8uPp6Sl1OJJRqVTIy8uDg4ODSRIgAHB0dAQA3L59G15eXuwOIyIio+Eg6FJSj/lxcnKSOBLroH6frXWsFRERmQYTIANZY7eXFPg+ExGRKTABIiIionKhVAJ798qwb58/9u6VwUSTiPXCBIiIiIhMLiYGCAoCunathAULWqJr10oIChKPS4EJkESUSiA+HvjpJ/GrlFnwiBEj0LdvX+kCICKiCi0mBnjlFeDmTe3jKSnicSmSIM4Ck0BMDDBhgvYPQkAAsHgx0L9/+cezePFiCIJQ/k9MREQVnlIpfubp+pgRBEAmAyZOBPr0Acpzsi9bgMqZOWbBrq6ucHNzK/8nJiKiCm///sKfeQUJAnDjhliuPDEBMgJBALKzS75lZgLjxxedBQNilpyZqd/1Sttos3nzZoSEhMDR0RGenp4ICwtDdnZ2oS6wR48eYejQoahcuTJ8fX2xcOFCdOrUCRMnTtSUqVmzJr766itERETA2dkZNWrUwK+//oo7d+6gT58+cHZ2RuPGjXH8+HHN99y7dw9DhgyBv78/nJycEBISgp9++ql0L4KIiCxKaqpxyxkLEyAjyMkBnJ1Lvrm6ii09RREEMUt2ddXvejk5+seYmpqKIUOGYNSoUTh//jzi4+PRv39/nV1fkZGROHjwIH799VfExsZi//79OHnyZKFy33zzDdq0aYNTp06hZ8+eGDZsGIYPH47XX38dJ0+eRK1atTB8+HDNczx58gQtWrTA9u3bcebMGYwdOxbDhg3D0aNH9X8hRERkUfRdyL+8F/znGCArkZqaivz8fPTv3x81atQAAISEhBQq9+jRI6xduxbr169Hly5dAACrV6+Gn59fobJdu3bFm2++CRsbG0RHR+Pbb7/F888/j4EDBwIAPvroI4SGhiI9PR0+Pj7w9/fHhx9+qPn+9957D7t27cLPP/+MVq1ameJlExGRxNq3F8e5pqTo7rmQycTz7duXb1xsATICJycgK6vk244d+l1vxw79rleaxaibNGmCLl26ICQkBAMHDsSKFSvw4MGDQuWuXbsGhUKhlZC4urqibt26hco2bNhQc9/b2xuAdlKlPqbez0upVGLmzJkICQmBh4cHnJ2dsWvXLiQnJ+v/QoiIyKLI5eIkH13Ua90uWlS+A6ABJkBGIZMBlSuXfOvWTcxyi1rcWCYDAgPFcvpcrzSLJMvlcsTGxuKPP/5AgwYNsHTpUtStWxdJSUkGv25bW9sCscuKPKZSqQAA8+bNw+LFi/HRRx9hz549SExMRHh4OPLy8gyOgYiIzF///sDYsYWPBwQAmzdLMwOaCVA5KpgFP5u8lEcWLJPJ0LZtW0yfPh2nTp2CnZ0dtmzZolWmZs2asLW1xbFjxzTHHj58iEuXLpX5+Q8ePIg+ffrg9ddfR5MmTVCzZk2jXJeIiMzfqVPi1zFjlIiMPI7Y2HwkJUmT/ABMgMpd//5ituvvr33c1FnwkSNHMGvWLBw/fhzJycmIiYnBnTt3UL9+fa1yVapUQUREBCZNmoQ9e/bg7NmzeOONN2BjY1Pmfbnq1KmD2NhYHDp0COfPn8ebb76J9PT0Ml2TiIjM35kzwNGjQKVKQHS0Ch06pKBjR6Hcu70K4iBoCfTvLy74tH+/OO3P11cc/GXKHwQXFxfs27cPixYtQmZmJmrUqIH58+eje/fu2Lhxo1bZBQsW4K233sLLL78MFxcXTJ48GTdu3ICDg0OZYvj0009x7do1hIeHw8nJCWPHjkXfvn3x8OHDMl2XiIjM2+rV4teXXwb+Gx4qOUkToH379mHevHk4ceIEUlNTsWXLlmK3ZBgxYgTWrl1b6HiDBg1w9uxZAMC0adMwffp0rfN169bFhQsXjBp7WcnlQKdO5fd89evXx86dO3WeW7NmjdbjKlWq4Mcff9Q8zs7OxvTp0zG2QAfutWvXkJmZqfV9z06pDwoK0jrm4eGBrVu3GvgKiIjIEuXlAevWifffeEPaWAqStAssOzsbTZo0wbJly/Qqv3jxYqSmpmpuN27cgIeHh2batVrDhg21yh04cMAU4VdYp06dwk8//YSrV6/i5MmTGDp0KACgT58+EkdGRESW5vffgbt3xd6Ol16SOpqnJG0B6t69O7p37653eVdXV7i6umoeb926FQ8ePMDIkSO1ylWqVAk+Pj56Xzc3Nxe5ubmax+qWDYVCAYVCoVVWoVBAEASoVCrN7KaKRqVS4auvvsLFixdhZ2eH5s2bY+/evfDw8NC8ZnXLjvq9MGUsgiBAoVBALmVnsZlT/5w++/NK0mB9mBfWh7RWrpQDsMHrryshCCqT1kdprmnRY4B++OEHhIWFaRb2U7t8+TL8/Pzg4OCA0NBQzJ49G9WrVy/yOrNnzy7UbQYAu3fvhtMzi+2ok6usrKwKO327Vq1aiIuLK3T82S4vQFw40ZTy8vLw+PFj7Nu3D/n5+SZ9roogNjZW6hCoANaHeWF9lL979xywa1c3AEBQ0B7s2JGtOWeK+sgpxRYJMsFMtgGXyWQljgEq6NatW6hevTrWr1+PV199VXP8jz/+QFZWFurWrYvU1FRMnz4dKSkpOHPmDKpUqaLzWrpagAIDA3H37l24uLholX3y5Alu3LiBoKCgMg8KtmSCIODRo0eoUqVKmWeHFefJkye4fv06AgMDrfr9LolCoUBsbCy6du2qtRYTSYP1YV5YH9KZM8cG0dFytG2rwp49SgCmrY/MzExUrVoVDx8+LPT5/SyLbQFau3Yt3NzcCiVMBbvUGjdujNatW6NGjRr4+eef8UYRo6/s7e1hb29f6LitrW2hylEqlZDJZLCxsYGNjfWuIqDu9lK/F6ainn6vqy6oML5P5oX1YV5YH+VLEAD1vKXRo21ga6v9WWGK+ijN9SzyE1wQBKxatQrDhg2DnZ1dsWXd3Nzw3HPP4cqVK+UUHREREe3bB1y9Km7e/corUkdTmEUmQHv37sWVK1eKbNEpKCsrC1evXoVveW8zS0REZMVWrRK/Dh4sJkHmRtIEKCsrC4mJiUhMTAQAJCUlITExUbM5ZlRUFIYPH17o+3744Qe0bt0ajRo1KnTuww8/xN69e3H9+nUcOnQI/fr1g1wux5AhQ0z6WoiIiEiUmQls2iTeHzVK2liKImkCdPz4cTRr1gzNmjUDAERGRqJZs2aIjo4GAKSmphbaKfzhw4f45Zdfimz9uXnzJoYMGYK6devi1VdfhaenJw4fPoxq1aqZ9sVYAEEQMHbsWHh4eEAmk2kSTyIiImPasAF4/BioXx944QWpo9FN0kHQnTp1KrR6cEHPrlAMiGsBFTfNbcOGDcYIzfSUyvLdCwPAzp07sWbNGsTHx6NmzZqoWrVqma43Z84c7Ny5k4kUERFpUXd/jRpVePNvc2Gxs8AsWkwMMGECcPPm02MBAeJW8SbcFlc9FqpNmzYmew5D5OXllTiYnYiILMPZs8CRI+LGp8OGSR1N0SxyELRFi4kRh8MXTH4AICVFPB4TY5KnHTFiBN577z0kJydDJpMhKCgIO3fuRLt27eDm5gZPT0+8/PLLuHr1qtb3qbsUPTw8ULlyZbRs2RJHjhzBmjVr8OWXX+Lvv/+GTCaDTCbTtNglJyejT58+cHZ2houLC1599VWtXd+nTZuGpk2bYuXKlQgODub6PkREFYi69cecNj7VhS1AxiAIgD6rTyqVwPjxYnld15DJxJahsDD9usOcnPRuW1y8eDFq1aqF77//HseOHYNcLse+ffsQGRmJxo0bIysrC9HR0ejXrx8SExNhY2ODrKwsdOzYEf7+/vj111/h4+ODkydPQqVSYdCgQTh16hT27NmDP//8E4DYPalSqTTJz969e5Gfn493330XgwYNQnx8vCaeK1eu4JdffkFMTAy3uCAiqiAKbnxqroOf1ZgAGUNOjnHm+AmC2DJUYL+zYmVlAZUr61XU1dUVVapUgVwu1+yTNmDAAK0yq1atQrVq1XDu3Dk0atQI69evx507d3Ds2DF4eHgAAGrXrg1AXAixcuXKhfZdi42NxenTp5GUlITAwEAAwLp169CwYUMcO3YMzz//PACx22vdunUcnE5EVIGoNz718QFKsdWnJNgFZsUuX76MIUOGoGbNmnBxcUFQUBAAaGbeJSYmolmzZprkRx/nz59HYGCgJvkBgAYNGsDNzQ3nz5/XHKtRowaTHyKiCuaHH8SvERHiGCBzZubhWQgnJ7E1piT79gE9epRcbscOoEMH/Z63DHr16oUaNWpgxYoV8PPzg0qlQqNGjTSbvDo6Opbp+sWprGfLFRERWYaUFGDnTvG+uXd/AUyAjEMm068rqls3cbZXSorucUAymXi+WzeTT4m/d+8eLl68iBUrVqB9+/YAgAMHDmiVady4MVauXIn79+/rbAWytbWFUqnUOla/fn3cuHEDN27c0LQCnTt3DhkZGWjQoIGJXg0REUlt3TpApQLatQOee07qaErGLrDyJJeLU92BwoOX1Y8XLTJ58gMA7u7u8PT0xPfff48rV67gr7/+QmRkpFaZIUOGwMfHB3379sXBgwdx7do1/PLLL0hISAAAVK9eXbN69927d5Gbm4uwsDCEhIRg6NChOHnyJI4ePYrhw4ejY8eOaNmypclfFxERlT9BeDr7S49dqswCE6Dy1r8/sHkz4O+vfTwgQDxuwnWACrKxscGGDRtw4sQJNGrUCO+//z7mzZunVcbOzg67d++Gl5cXevTogZCQEMyZM0cza6t3794IDw9H586dUa1aNfz000+QyWTYtm0b3N3d0aFDB4SFhaFmzZrYuHFjubwuIiIqf/v3A1eumO/Gp7qwC0wK/fsDffqU+0rQEydOxMSJEzWPw8LCcO7cOa0yz67MXaNGDWzevLnQtVQqFezt7bFp0ybY2Gjn0dWrV8e2bduKjGPatGmYNm1a6V8AERGZJfXgZ3Pd+FQXJkBSkcuBTp2kjoKIiKhMLGHjU13YBUZEREQGs4SNT3VhAkREREQGs4SNT3VhAkREREQGsZSNT3VhAmSgZwcLk2nwfSYiMl+WsvGpLkyASsnW1hYAkKPP5qdUZur3Wf2+ExGRecjLA/73P/G+JQ1+VuMssFKSy+Vwc3PD7du3AQBOTk6QWVKnp5GoVCrk5eXhyZMnhabBG4MgCMjJycHt27fh5ubGHeOJiMzM778Dd+5YxsanujABMoB693N1EmSNBEHA48eP4ejoaNIE0M3NTWu3eSIiMg/q7i9L2PhUFwsMWXoymQy+vr7w8vKCQqGQOhxJKBQK7Nu3Dx06dDBZ95StrS1bfoiIzFBKCvDHH+J9S+z+ApgAlYlcLrfaD2i5XI78/Hw4ODhwfA4RkZWxtI1PdeEgaCIiItKbJW58qgsTICIiItKbJW58qgsTICIiItKbeuPTQYMsZ+NTXZgAERERkV4Kbnxqyd1fABMgIiIi0tPGjeLGp/XqWdbGp7owASIiIiK9qLu/3njDsjY+1YUJEBEREZXIkjc+1YUJEBEREZXIkjc+1YUJEBERERXL0jc+1YUJEBERERXL0jc+1YUJEBERERXL0jc+1YUJEBERERWpImx8qgsTICIiIipSRdj4VBcmQERERKRTwY1PK1LrD8AEiIiIiIpQcOPTgQOljsa4JE2A9u3bh169esHPzw8ymQxbt24ttnx8fDxkMlmhW1pamla5ZcuWISgoCA4ODmjdujWOHj1qwldBRERUMalbfyx941NdJE2AsrOz0aRJEyxbtqxU33fx4kWkpqZqbl5eXppzGzduRGRkJKZOnYqTJ0+iSZMmCA8Px+3bt40dPhERUYVVkTY+1UXSyWzdu3dHdwMWFPDy8oKbm5vOcwsWLMCYMWMwcuRIAMDy5cuxfft2rFq1ClOmTClLuERERFZj40YgJ6dibHyqi0XO5m/atClyc3PRqFEjTJs2DW3btgUA5OXl4cSJE4iKitKUtbGxQVhYGBISEoq8Xm5uLnJzczWPMzMzAQAKhQIKhcJEr8Kyqd8Xvj/mgfVhXlgf5oX1YZiVK+UAbDBihBL5+SqjXdeU9VGaa1pUAuTr64vly5ejZcuWyM3NxcqVK9GpUyccOXIEzZs3x927d6FUKuH9zCYl3t7euHDhQpHXnT17NqZPn17o+O7du+Hk5GT011GRxMbGSh0CFcD6MC+sD/PC+tBfcnIVHD36IuRyFXx8YrFjR27J31RKpqiPnJwcvctaVAJUt25d1K1bV/O4TZs2uHr1KhYuXIj/qTcpMUBUVBQiIyM1jzMzMxEYGIhu3brBxcWlTDFXVAqFArGxsejatStsbW2lDsfqsT7MC+vDvLA+Sm/yZHGIcM+ewGuvdTHqtU1ZH+oeHH1YVAKkS6tWrXDgwAEAQNWqVSGXy5Genq5VJj09HT4+PkVew97eHvb29oWO29ra8pelBHyPzAvrw7ywPswL60M/eXnAjz+K90ePtoGtrWnmS5miPkpzPYtfBygxMRG+vr4AADs7O7Ro0QJxcXGa8yqVCnFxcQgNDZUqRCIiIouxfXvF2/hUF0lbgLKysnDlyhXN46SkJCQmJsLDwwPVq1dHVFQUUlJSsG7dOgDAokWLEBwcjIYNG+LJkydYuXIl/vrrL+zevVtzjcjISERERKBly5Zo1aoVFi1ahOzsbM2sMCIiIiraDz+IXyvSxqe6SPrSjh8/js6dO2seq8fhREREYM2aNUhNTUVycrLmfF5eHj744AOkpKTAyckJjRs3xp9//ql1jUGDBuHOnTuIjo5GWloamjZtip07dxYaGE1ERETabt16uvFpRW83kDQB6tSpEwRBKPL8mjVrtB5PnjwZkydPLvG648aNw7hx48oaHhERkVVZu/bpxqcF5hxVSBY/BoiIiIjKriJvfKoLEyAiIiKq0Buf6sIEiIiIiCr0xqe6MAEiIiKychV941NdmAARERFZuYq+8akuTICIiIisnLr76403AJlM2ljKCxMgIiIiK3buHHD4MCCXA8OGSR1N+WECREREZMXUrT8vvwxY05rBFXiRayIiItJFqRSnvd+4AaxcKR6zlsHPakyAiIiIrEhMDDBhAnDz5tNjNjbAkyfSxSQFJkBERERWIiYGeOUVcdXnglQqcf0fuRzo31+a2MobxwARERFZAaVSbPkpZgtOTJwolrMGTICIiIiswP792t1ezxIEcUzQ/v3lF5OUmAARERFZgdRU45azdEyAiIiIrICvr3HLWTomQERERFagfXsgIKDo8zIZEBgolrMGTICIiIisgFwOREToPqfe/mLRIrGcNWACREREZAVu33666KGzs/a5gABg82brmQIPcB0gIiKiCk8QgNGjgfR0oFEjce+vY8fEAc++vmK3l7W0/KgxASIiIqrgvvsO+O03wM4O+PFHoHJloFMnqaOSFrvAiIiIKrALF4DISPH+l18CjRtLG4+5YAJERERUQeXlAa+9Bjx+DHTtCowfL3VE5oMJEBERUQX12WfAqVOApyewZo246SmJ+FYQERFVQHv2APPmifdXrgT8/KSNx9wwASIiIqpg7t8Hhg0TZ3+NGQP07St1ROaHCRAREVEFIgjAW28BKSnAc88BCxdKHZF5YgJERERUgaxbB2zaBFSq9HTKOxXGBIiIiKiCuHoVGDdOvD9jBtCypbTxmDMmQERERBVAfj7w+utAVhbQoQMwebLUEZk3JkBEREQVwOefi1tcuLoC//uf9W1tUVpMgIiIiCzcoUPAzJni/eXLgerVpY3HEjABIiIismCZmcDQoYBKJU59HzxY6ogsAxMgIiIiC/bee8D160BQEPD111JHYzmYABEREVmojRvFae82NsD//R/g4iJ1RJaDCRAREZEFSk4WFzwEgE8+Adq2lTYeS8MEiIiIyMIolcDw4UBGBtC6tbjpKZWOpAnQvn370KtXL/j5+UEmk2Hr1q3Flo+JiUHXrl1RrVo1uLi4IDQ0FLt27dIqM23aNMhkMq1bvXr1TPgqiIiIyte8ecDevYCzs7jas62t1BFZHkkToOzsbDRp0gTLli3Tq/y+ffvQtWtX7NixAydOnEDnzp3Rq1cvnDp1Sqtcw4YNkZqaqrkdOHDAFOETERGVuxMnnrb4LFkC1KolbTyWqlJpv+H8+fPYsGED9u/fj3///Rc5OTmoVq0amjVrhvDwcAwYMAD29vZ6Xat79+7o3r273s+9aNEircezZs3Ctm3b8Ntvv6FZs2aa45UqVYKPj4/e1yUiIrIE2dnAa6+Jqz6/8gowYoTUEVkuvROgkydPYvLkyThw4ADatm2L1q1bo1+/fnB0dMT9+/dx5swZfPLJJ3jvvfcwefJkTJw4Ue9EyFAqlQqPHj2Ch4eH1vHLly/Dz88PDg4OCA0NxezZs1G9mFWhcnNzkZubq3mcmZkJAFAoFFAoFKYJ3sKp3xe+P+aB9WFeWB/mpSLVx8SJNrh0SQ5/fwFff52P/HypIyo9U9ZHaa4pEwRB0KdgcHAwJk2ahNdeew1ubm5FlktISMDixYvRuHFjfPzxx/oHIpNhy5Yt6Nu3r97fM3fuXMyZMwcXLlyAl5cXAOCPP/5AVlYW6tati9TUVEyfPh0pKSk4c+YMqlSpovM606ZNw/Tp0wsdX79+PZycnPSOh4iIyFSOHPHB7NmtIZMJmD79EBo3vit1SGYnJycHr732Gh4+fAiXEtYE0DsBUigUsC3FKKvSli9tArR+/XqMGTMG27ZtQ1hYWJHlMjIyUKNGDSxYsABvvPGGzjK6WoACAwNx9+7dEt9Aa6VQKBAbG4uuXbuWqp7JNFgf5oX1YV4qQn2kpgLNm1fCvXsyREYqMWeOSuqQDGbK+sjMzETVqlX1SoD07gIrbZCm/CHbsGEDRo8ejU2bNhWb/ACAm5sbnnvuOVy5cqXIMvb29jq762xtbS32l6W88D0yL6wP88L6MC+WWh8qFTBmDHDvHtC0KTBrlhy2tpa/06kp6qM01zPqLLDjx49j3759xrxkIT/99BNGjhyJn376CT179iyxfFZWFq5evQpfX1+TxkVERGQKS5cCu3cDDg7A+vWAiYfXWo1SzwIrzrBhw3Dp0iUolUq9ymdlZWm1zCQlJSExMREeHh6oXr06oqKikJKSgnXr1gEQu70iIiKwePFitG7dGmlpaQAAR0dHuLq6AgA+/PBD9OrVCzVq1MCtW7cwdepUyOVyDBkyxJgvlYiIyOROnwY++ki8P38+UL++tPFUJEZtAYqLi8O1a9f0Ln/8+HE0a9ZMM4U9MjISzZo1Q3R0NAAgNTUVycnJmvLff/898vPz8e6778LX11dzmzBhgqbMzZs3MWTIENStWxevvvoqPD09cfjwYVSrVs1Ir5KIiMj0njwRp7zn5gI9ewJvvy11RBWLUVuA/Pz8SlW+U6dOKG4M9po1a7Qex8fHl3jNDRs2lCoGIiIiczRlCnDmDODlBaxaBchkUkdUsRjUArRjx45CW1AAwK5du/DHH3+UOSgiIiJrtmsXsHixeH/1ajEJIuMyKAGaMmWKznE+giBgypQpZQ6KiIjIWt2583SF53HjgB49JA2nwjIoAbp8+TIaNGhQ6Hi9evWKnW5ORERERRMEYPRoIC0NaNAAmDtX6ogqLoPGALm6uuLatWsICgrSOn7lyhVUrlzZGHERERFVeEolsH+/uNChry9w4QLw66+AnZ24y7ujo9QRVlwGJUB9+vTBxIkTsWXLFtT6bxvaK1eu4IMPPkDv3r2NGiAREVFFFBMDTJgA3LxZ+Nzs2eKih2Q6BnWBzZ07F5UrV0a9evUQHByM4OBg1K9fH56envjqq6+MHSMREVGFEhMj7uauK/kBgGL27yYjMbgL7NChQ4iNjcXff/8NR0dHNG7cGB06dDB2fERERBWKUim2/BS1CoxMBkRGAv36AXLL3/HCbBm8DpBMJkO3bt3QoUMH2NvbQ8YFCoiIiEq0f3/RLT+AmBjduCGW69Sp3MKyOgZ1galUKsycORP+/v5wdnZGUlISAOCzzz7DDz/8YNQAiYiIKpLUVOOWI8MYlAB9/vnnWLNmDebOnQs7OzvN8UaNGmHlypVGC46IiKiicXLSrxz38DYtgxKgdevW4fvvv8fQoUMhL9BB2aRJE1y4cMFowREREVUUubnAV18Bw4cXX04mAwIDgfbtyycua2VQApSSkoLatWsXOq5SqaBQKMocFBERUUUhCMCmTeJO7pMmAZmZgHoZvWeHz6ofL1rEAdCmZlAC1KBBA+zfv7/Q8c2bN2t2diciIrJ2R44A7doBr74KJCWJ3VqrVwNXrgC//AL4+2uXDwgANm8G+veXJl5rYtAssOjoaERERCAlJQUqlQoxMTG4ePEi1q1bh99//93YMRIREVmU69eBqChgwwbxsZMTMHky8OGHgHrDhP79gT59tFeCbt+eLT/lxeCVoH/77TfMmDEDlStXRnR0NJo3b47ffvsNXbt2NXaMREREFiEzU1zFeeFCccyPTCZubDpzZuHWHkBMdjjVXRoGrwPUvn17xMbGGjMWIiIii5SfD6xcCURHi7u5A8CLLwLz53NLC3Nl0BigGzdu4GaBVZyOHj2KiRMn4vvvvzdaYEREROZOEIA//gCaNAHefltMfurWFTc0/fNPJj/mzKAE6LXXXsOePXsAAGlpaQgLC8PRo0fxySefYMaMGUYNkIiIyBydPg2EhwM9egDnzgGensDSpeLxXr0Kz/Ai82JQAnTmzBm0atUKAPDzzz8jJCQEhw4dwo8//og1a9YYMz4iIiKzkpYGjBkjtu7ExgJ2duLg5itXgHHjAFtbqSMkfRg0BkihUMDe3h4A8Oeff6J3794AgHr16iGVa3cTEVEFlJMDLFgAzJkDZGeLxwYOFB/XrCltbFR6BrUANWzYEMuXL8f+/fsRGxuLl156CQBw69YteHp6GjVAIiIiU1Mqgfh44KefxK9K5dNzKhXwv/+JY3s++0xMflq3Bg4cAH7+mcmPpTKoBejLL79Ev379MG/ePERERKBJkyYAgF9//VXTNUZERGQJYmKACRO0d2gPCAAWLwaqVgUiI4ETJ8TjNWqILT6DBnGMj6UzKAHq1KkT7t69i8zMTLi7u2uOjx07Fk767vJGREQksZgY4JVXxNlcBaWkAAMGPH1cpQrwySdiouTgUL4xkmkYvA6QXC7XSn4AIEi9uQkREZGZUyrFhObZ5AfQPvbmm8CMGYCXV/nFRqan9xigl156CYcPHy6x3KNHj/Dll19i2bJlZQqMiIjIlPbv1+72KsrgwUx+KiK9W4AGDhyIAQMGwNXVFb169ULLli3h5+cHBwcHPHjwAOfOncOBAwewY8cO9OzZE/PmzTNl3ERERGWi76RlTm6umPROgN544w28/vrr2LRpEzZu3Ijvv/8eDx8+BADIZDI0aNAA4eHhOHbsGOrXr2+ygImIiIzB19e45ciylGoMkL29PV5//XW8/vrrAICHDx/i8ePH8PT0hC1XfiIiIgvSvr042yslRfc4IJlMPN++ffnHRqZn0DpAaq6urvDx8WHyQ0REFkcuF6e6F5X8AMCiRWI5qnjKlAARERFZMjc33ccDAoDNm4H+/cs1HCpHBk+DJyIismRKJfDBB+L9d94Rt7VITRXH/LRvz5afio4JEBERWaV164DERMDVFZg+XVz1mawHu8CIiMjqZGWJKzsD4v5eTH6sj8EJUEZGBlauXImoqCjcv38fAHDy5EmkpKQYLTgiIiJTmDdP7O6qVQsYN07qaEgKBnWB/fPPPwgLC4OrqyuuX7+OMWPGwMPDAzExMUhOTsa6deuMHScREZFR3LwpJkAAMHcuYG8vbTwkDYNagCIjIzFixAhcvnwZDgV2hevRowf27dtntOCIiIiM7ZNPgMePxYHO/fpJHQ1JxaAE6NixY3jzzTcLHff390daWpre19m3bx969eoFPz8/yGQybN26tcTviY+PR/PmzWFvb4/atWtjzZo1hcosW7YMQUFBcHBwQOvWrXH06FG9YyIioorr+HFx8DMALFjwdL0fsj4GJUD29vbIzMwsdPzSpUuoVq2a3tfJzs5GkyZN9N44NSkpCT179kTnzp2RmJiIiRMnYvTo0di1a5emzMaNGxEZGYmpU6fi5MmTaNKkCcLDw3H79m294yIioopHEIDISPH+sGFAy5bSxkPSMmgMUO/evTFjxgz8/PPPAMS9wJKTk/HRRx9hwIABel+ne/fu6N69u97lly9fjuDgYMyfPx8AUL9+fRw4cAALFy5EeHg4AGDBggUYM2YMRo4cqfme7du3Y9WqVZgyZYrO6+bm5iI3N1fzWJ3cKRQKKBQKveOzJur3he+PeWB9mBfWh3lR18Mvv6iwfz/g6Chg2rR8sHqkYcrfj9Jc06AEaP78+XjllVfg5eWFx48fo2PHjkhLS0NoaCi++OILQy6pl4SEBISFhWkdCw8Px8SJEwEAeXl5OHHiBKKiojTnbWxsEBYWhoSEhCKvO3v2bEyfPr3Q8d27d8PJyck4wVdQsbGxUodABbA+zAvrw3woFDaIjFQAsEevXpdw+vQFnD4tdVTWzRS/Hzk5OXqXNSgBcnV1RWxsLA4cOIB//vkHWVlZaN68eaHkxNjS0tLg7e2tdczb2xuZmZl4/PgxHjx4AKVSqbPMhQsXirxuVFQUItXtohBbgAIDA9GtWze4uLgY90VUEAqFArGxsejatSv3gjMDrA/zwvowLwqFAuPGXUNamjN8fAQsX14Tzs41pQ7Lapny90PX8JyilGkl6Hbt2qFdu3ZluYRZsLe3h72OeZC2trb841UCvkfmhfVhXlgf5uHuXWDjxroAgC++kMHdnXViDkzx+1Ga6xmcAB07dgx79uzB7du3oVKptM4tWLDA0MsWy8fHB+np6VrH0tPT4eLiAkdHR8jlcsjlcp1lfHx8TBITERGZty++sEFOjhyNGwuIiOC0LxIZlADNmjULn376KerWrQtvb2/ICswjlJlwTmFoaCh27NihdSw2NhahoaEAADs7O7Ro0QJxcXHo27cvAEClUiEuLg7juNQnEZHVuXABWL5cnPA8b54Scjm3wCSRQT8JixcvxqpVqzBixIgyPXlWVhauXLmieZyUlITExER4eHigevXqiIqKQkpKimZl6bfeegtff/01Jk+ejFGjRuGvv/7Czz//jO3bt2uuERkZiYiICLRs2RKtWrXCokWLkJ2drZkVRkRE1mPyZECplKFVq1R07swNv+gpgxIgGxsbtG3btsxPfvz4cXTu3FnzWD0QOSIiAmvWrEFqaiqSk5M154ODg7F9+3a8//77WLx4MQICArBy5UrNFHgAGDRoEO7cuYPo6GikpaWhadOm2LlzZ6GB0UREVLHFxQG//QZUqiQgIuIcgA5Sh0RmxKAE6P3338eyZcuwaNGiMj15p06dIAhCked1rfLcqVMnnDp1qtjrjhs3jl1eRERWTKl8uujhW2+p4O+fJW1AZHYMSoA+/PBD9OzZE7Vq1UKDBg0KjbqOiYkxSnBERESGWLMG+OcfwM0N+OQTFY4ckToiMjcGJUDjx4/Hnj170LlzZ3h6epp04DMREVFpPHoEfPqpeD86GvD0lDYeMk8GJUBr167FL7/8gp49exo7HiIiojKZOxdISwNq1wbefVfqaMhcGbQZqoeHB2rVqmXsWIiIiMrkxg3gq6/E+3PnAnZ20sZD5sugBGjatGmYOnVqqfbcICIiMrWPPwaePAE6dAD+Ww6OSCeDusCWLFmCq1evwtvbG0FBQYUGQZ88edIowREREenr2DHg//5PvL9gAcDhqVQcgxKgvkyriYjIjAjC02nvw4cDLVpIGw+ZP4MSoKlTpxo7DiIiIoPFxAAHDgCOjsAXX0gdDVkCg8YAERERmYvcXHHLC0D8GhAgbTxkGfRuAfLw8MClS5dQtWpVuLu7F7v2z/37940SHBERUUmWLgWuXQN8fYFJk6SOhiyF3gnQwoULUaVKFc19Ln5IRERSu3MHmDlTvD9rFlC5srTxkOXQOwGKiIjQ3C/rLvBERETGMH06kJkJNGsmDn4m0pdBY4Dkcjlu375d6Pi9e/cgl8vLHBQREVFJzp8Hli8X78+fD9hwVCuVgkE/LkXt4J6bmws7LrtJRETlYNIkcdf3Pn2Azp2ljoYsTammwS9ZsgQAIJPJsHLlSjg7O2vOKZVK7Nu3D/Xq1TNuhERERM+IjQW2bwcqVRK3vCAqrVIlQAsXLgQgtgAtX75cq7vLzs4OQUFBWK5ujyQiIjIBpRL44APx/rvvAs89J208ZJlKlQAlJSUBADp37oyYmBi4u7ubJCgiIqKirF4NnD4NuLsD0dFSR0OWyqCVoPfs2WPsOIiIiEr06BHw6afi/ehowMND2njIcnHMPBERWYwvvwTS04E6dYB33pE6GrJkTICIiMgiJCeL090BYN48gJOOqSyYABERkUWIigKePAE6dQJ695Y6GrJ0TICIiMjsHTkCrF8PyGRiKxB3Y6KyMmgQNABkZGTg6NGjuH37NlQqlda54VyPnIiIjEQQgMhI8X5EBNC8ubTxUMVgUAL022+/YejQocjKyoKLi4vWxqgymYwJEBERGc3mzcChQ4CTE/D551JHQxWFQV1gH3zwAUaNGoWsrCxkZGTgwYMHmtv9+/eNHSMREVmpJ0+Ajz4S70+eDPj7SxsPVRwGJUApKSkYP348nJycjB0PERGRxtKlQFIS4OcHfPih1NFQRWJQAhQeHo7jx48bOxYiIiKNO3eednnNmgVUrixtPFSxGDQGqGfPnpg0aRLOnTuHkJAQ2Nraap3vzfmJRERkAKUS2L8fSE0FNm4EMjPFQc/DhkkdGVU0BiVAY8aMAQDMmDGj0DmZTAalUlm2qIiIyOrExAATJgA3b2of79sXsOGiLWRkBiVAz057JyIiKouYGOCVV8Qp78+aOhVo2BDo37/846KKizk1ERFJSqkUW350JT9qEyeK5YiMxeAEaO/evejVqxdq166N2rVro3fv3ti/f78xYyMiIiuwf3/hbq+CBAG4cUMsR2QsBiVA//d//4ewsDA4OTlh/PjxGD9+PBwdHdGlSxesX7/e2DESEVEFlppq3HJE+jBoDNAXX3yBuXPn4v3339ccGz9+PBYsWICZM2fitddeM1qARERUcalUwOHD+pX19TVtLGRdDGoBunbtGnr16lXoeO/evZGUlFTmoIiIqOK7fFnc2X3JkuLLyWRAYCDQvn25hEVWwqAEKDAwEHFxcYWO//nnnwgMDCz19ZYtW4agoCA4ODigdevWOHr0aJFlO3XqBJlMVujWs2dPTZkRI0YUOv/SSy+VOi4iIjI+pVLc0b1xY3FcT+XKwOjRYqLz7C7v6seLFgFyebmHShWYQV1gH3zwAcaPH4/ExES0adMGAHDw4EGsWbMGixcvLtW1Nm7ciMjISCxfvhytW7fGokWLEB4ejosXL8LLy6tQ+ZiYGOTl5Wke37t3D02aNMHAgQO1yr300ktYvXq15rG9vX2p4iIiIuM7dw4YNQo4ckR8HBYGrFgBBAUB3bsXXgcoIEBMfjgFnozNoATo7bffho+PD+bPn4+ff/4ZAFC/fn1s3LgRffr0KdW1FixYgDFjxmDkyJEAgOXLl2P79u1YtWoVpkyZUqi8h4eH1uMNGzbAycmpUAJkb28PHx+fUsVCRESmoVAA8+YB06cDeXmAiwuwYIGYDKlbefr3B/r0eboStK+v2O3Flh8yBYMSIADo168f+vXrV6Ynz8vLw4kTJxAVFaU5ZmNjg7CwMCQkJOh1jR9++AGDBw9G5Wc2iYmPj4eXlxfc3d3x4osv4vPPP4enp6fOa+Tm5iI3N1fzODMzEwCgUCigUChK+7Ksgvp94ftjHlgf5oX1oS0xERg7thISE8VMp0cPFb7+WomAACA/v3D5tm2f3lepxFtZsD7MiynrozTXlAlCcUtPmdatW7fg7++PQ4cOITQ0VHN88uTJ2Lt3L46o20iLcPToUbRu3RpHjhxBq1atNMfVrULBwcG4evUqPv74Yzg7OyMhIQFyHf9KTJs2DdOnTy90fP369dzxnojIQAqFDJs21cUvv9SBUmkDZ+c8jB59Gh073iw01ofIGHJycvDaa6/h4cOHcHFxKbas3i1AHh4euHTpEqpWrQp3d3fIivnpvX//vv7RlsEPP/yAkJAQreQHAAYPHqy5HxISgsaNG6NWrVqIj49Hly5dCl0nKioKkZGRmseZmZkIDAxEt27dSnwDrZVCoUBsbCy6du1aaDNcKn+sD/PC+gCOH5dhzBg5zp4VPyv69lVhyRIZfHwaA2hcrrGwPsyLKetD3YOjD70ToIULF6JKlSqa+8UlQPqqWrUq5HI50tPTtY6np6eXOH4nOzsbGzZs0Lkh67Nq1qyJqlWr4sqVKzoTIHt7e52DpG1tbfnLUgK+R+aF9WFerLE+Hj8Gpk0DvvpK7LqqVg1YtgwYONAGUu++ZI31Yc5MUR+luZ7eCVBERITm/ogRI0oVUFHs7OzQokULxMXFoW/fvgDEjVbj4uIwbty4Yr9306ZNyM3Nxeuvv17i89y8eRP37t2DL1fRIiIymYMHxUHNly6Jj197DVi8GKhaVdq4iHQxKB2Xy+W4fft2oeP37t3TOcamOJGRkVixYgXWrl2L8+fP4+2330Z2drZmVtjw4cO1Bkmr/fDDD+jbt2+hgc1ZWVmYNGkSDh8+jOvXryMuLg59+vRB7dq1ER4eXqrYiIioZNnZ4mal7duLyY+vL7BtG/Djj0x+yHwZNAusqHHTubm5sLOzK9W1Bg0ahDt37iA6OhppaWlo2rQpdu7cCW9vbwBAcnIybGy087SLFy/iwIED2L17d6HryeVy/PPPP1i7di0yMjLg5+eHbt26YebMmVwLiIjIyPbsAd54A1BvAjBypDi93c1N0rCISlSqBGjJf+uVy2QyrFy5Es7OzppzSqUS+/btQ7169UodxLhx44rs8oqPjy90rG7dukUmYY6Ojti1a1epYyAiIv1lZgIffQQsXy4+DgwUFzRkQztZilIlQAsXLgQgtgAtX75cq7vLzs4OQUFBWK7+bSAiIoulVBa9IOHOncDYscCNG+Ljt98G5swRFzckshSlSoDUG5127twZMTExcHd3N0lQREQknZgY3VtSfPGF2OW1Zo14rGZN4IcfxA1NiSyNQWOA9uzZY+w4iIjIDMTEAK+8Ajw7yuDmTUA9GVgmExOkzz8XNzIlskQGJUCjRo0q9vyqVasMCoaIiKSjVIqJTXH7A1SqBPz1l9glRmTJDEqAHjx4oPVYoVDgzJkzyMjIwIsvvmiUwIiIqHzt36/d7aVLfr6YKBFZOoMSoC1bthQ6plKp8Pbbb6NWrVplDoqIiMpfaqpxyxGZM6OtS25jY4PIyEjNTDEiIrIs+i6Wz0X1qSIw6sYsV69eRX5+vjEvSURE5aR9e3G2V1FkMnG9H47/oYrAoC6wgjunA+K6QKmpqdi+fbvWnmFERGQ55HKgRw/g++8Ln1Pvf71o0dP1gIgsmUEJ0KlTp7Qe29jYoFq1apg/f36JM8SIiMg8nToFrF0r3nd1BR4+fHouIEBMfvr3lyQ0IqPjOkBERITMTODVV4HcXKBXL3E9oAMHdK8ETVQRGJQAERFRxSEIwJgxwJUrQPXq4krPlSpxhWeq2PROgJo1awaZuhO4BCdPnjQ4ICIiKl/ffgv8/LOY9GzcCHh4SB0RkenpnQD17dvXhGEQEZEUTp4E3n9fvD93LvDCC9LGQ1Re9E6Apk6daso4iIionD18CAwcCOTlAb17AxMnSh0RUfkp0xigEydO4Pz58wCAhg0bolmzZkYJioiITEsQgNGjgWvXgBo1xHE/eo5yIKoQDEqAbt++jcGDByM+Ph5ubm4AgIyMDHTu3BkbNmxAtWrVjBkjEREZ2TffAJs3A7a24vgfd3epIyIqXwatBP3ee+/h0aNHOHv2LO7fv4/79+/jzJkzyMzMxPjx440dIxERGdGJE4B6Pdu5c4FWraSNh0gKBrUA7dy5E3/++Sfq16+vOdagQQMsW7YM3bp1M1pwRERkXBkZT8f99O0LTJggdURE0jCoBUilUsHW1rbQcVtbW6hUqjIHRURExicIwBtvAElJQFAQsGoVx/2Q9TIoAXrxxRcxYcIE3Lp1S3MsJSUF77//Prp06WK04IiIyHi+/lpc4ZnjfogMTIC+/vprZGZmIigoCLVq1UKtWrUQHByMzMxMLF261NgxEhFRGR07BnzwgXj/q6+A55+XNh4iqRk0BigwMBAnT57En3/+iQsXLgAA6tevj7CwMKMGR0REZZeRIe7zpVCIm5m+957UERFJz+B1gGQyGbp27YquXbsaMx4iIjIiQQBGjQKuXweCg4EffuC4HyKglF1gCQkJ+P3337WOrVu3DsHBwfDy8sLYsWORm5tr1ACJiMhwS5YAW7Y8Hffz39JtRFavVAnQjBkzcPbsWc3j06dP44033kBYWBimTJmC3377DbNnzzZ6kEREVHpHjwKTJon3588HWraUNh4ic1KqBCgxMVFrlteGDRvQunVrrFixApGRkViyZAl+/vlnowdJRESl8+DB03E/AwYA48ZJHRGReSlVAvTgwQN4e3trHu/duxfdu3fXPH7++edx48YN40VHRESlJgjAyJHAv/8CNWty3A+RLqVKgLy9vZGUlAQAyMvLw8mTJ/HCCy9ozj969EjnAolERFR+Fi0Ctm0D7OzEcT+urlJHRGR+SpUA9ejRA1OmTMH+/fsRFRUFJycntG/fXnP+n3/+Qa1atYweJBER6efIEWDyZPH+ggVAixbSxkNkrko1DX7mzJno378/OnbsCGdnZ6xduxZ2dnaa86tWreJeYEREErl/Hxg0CMjPF/f7eucdqSMiMl+lSoCqVq2Kffv24eHDh3B2doZcLtc6v2nTJjg7Oxs1QCIiKlnBcT+1agErVnDcD1FxDFoI0bWIDmUPD48yBUNERIZZuBD49VeO+yHSl0F7gRERkfk4fBj46CPx/qJFQPPmkoZDZBGYABERWbD798X1fvLzxa9vvSV1RESWgQkQEZGFUqmAiAjgxg2gdm2O+yEqDbNIgJYtW4agoCA4ODigdevWOHr0aJFl16xZA5lMpnVzcHDQKiMIAqKjo+Hr6wtHR0eEhYXh8uXLpn4ZRETlasEC4PffAXt7cdyPi4vUERFZDskToI0bNyIyMhJTp07FyZMn0aRJE4SHh+P27dtFfo+LiwtSU1M1t3///Vfr/Ny5c7FkyRIsX74cR44cQeXKlREeHo4nT56Y+uUQEZWLQ4eAKVPE+4sWAc2aSRoOkcWRPAFasGABxowZg5EjR6JBgwZYvnw5nJycsGrVqiK/RyaTwcfHR3MruD2HIAhYtGgRPv30U/Tp0weNGzfGunXrcOvWLWzdurUcXhERkWnduwcMHgwoleLXN9+UOiIiy2PQNHhjycvLw4kTJxAVFaU5ZmNjg7CwMCQkJBT5fVlZWahRowZUKhWaN2+OWbNmoWHDhgCApKQkpKWlISwsTFPe1dUVrVu3RkJCAgYPHlzoerm5ucjNzdU8zszMBAAoFAooFIoyv86KSP2+8P0xD6wP82LK+lCpgOHD5bhxwwa1awv4+ut85Ocb/WkqFP5+mBdT1kdprilpAnT37l0olUqtFhxA3HPswoULOr+nbt26WLVqFRo3boyHDx/iq6++Qps2bXD27FkEBAQgLS1Nc41nr6k+96zZs2dj+vTphY7v3r0bTk5Ohrw0qxEbGyt1CFQA68O8GKM+lErg3DlPPHjgAHf3J7h40R07djSEra0S7767DwcOZBohUuvA3w/zYor6yMnJ0buspAmQIUJDQxEaGqp53KZNG9SvXx/fffcdZs6cadA1o6KiEBkZqXmcmZmJwMBAdOvWDS4cVaiTQqFAbGwsunbtyg1wzQDrw7wYqz62bJEhMlKOlJSCU7sEAMDixQJGj25XxkitA38/zIsp60Pdg6MPSROgqlWrQi6XIz09Xet4eno6fHx89LqGra0tmjVrhitXrgCA5vvS09Ph6+urdc2mTZvqvIa9vT3s7e11Xpu/LMXje2ReWB/mpSz1ERMjju8RhGfPiMmQl1clsKpLh78f5sUU9VGa60k6CNrOzg4tWrRAXFyc5phKpUJcXJxWK09xlEolTp8+rUl2goOD4ePjo3XNzMxMHDlyRO9rEhFJSakEJkzQlfyIZDLg/ffFckRkGMm7wCIjIxEREYGWLVuiVatWWLRoEbKzszFy5EgAwPDhw+Hv74/Zs2cDAGbMmIEXXngBtWvXRkZGBubNm4d///0Xo0ePBiDOEJs4cSI+//xz1KlTB8HBwfjss8/g5+eHvn37SvUyiYj0tn8/cPNm0ecFQVz8cP9+oFOncguLqEKRPAEaNGgQ7ty5g+joaKSlpaFp06bYuXOnZhBzcnIybGyeNlQ9ePAAY8aMQVpaGtzd3dGiRQscOnQIDRo00JSZPHkysrOzMXbsWGRkZKBdu3bYuXNnoQUTiYjMUWqqccsRUWGSJ0AAMG7cOIwbN07nufj4eK3HCxcuxMKFC4u9nkwmw4wZMzBjxgxjhUhEVC6SksRVnfVRYJgjEZWSWSRARETW7sQJYN48YNMmca2f4shkQEAA0L59+cRGVBFJvhI0EZG1EgRg1y6gSxegZUtg40Yx+enWDYiOFhOdZzc3VT9etAiQy8s9ZKIKgy1ARETlTKEANmwAvvoK+Ocf8ZhcLk57//BDQL1iR5Mm4mywggOiAwLE5Kd///KOmqhiYQJERFROHj0CVqwQE5gbN8RjlSsDY8YAEycCNWpol+/fH+jTR5ztlZoqjvlp354tP0TGwASIiMjEUlOBJUuAb78FHj4Uj3l7A+PHA2+/Dbi7F/29cjmnuhOZAhMgIiITuXABWLwY+N//gLw88dhzz4ndXMOGAVyZg0g6TICIiIzs0CEZZs1qhaNHny7L36YNMGkS0Ls3YMPpJ0SSYwJERKQHpbL4sTgqFbBtmziVPSGhEgBfyGQC+vSRYdIkMQEiIvPBBIiIqAQxMbpnYy1eDPToAaxbB8yfD1y6JJ6zsxPQseO/WLDAH40acfNNInPEBIiIqBgxMcArrxTemDQlBRgwAHBxATIzxWNubsA77wBvvZWPkyf/Rt26/uUeLxHphwkQEVERituVXX0sM1NsDfrgA+CNN4AqVcR1fojIvDEBIiIqQkm7squtXg2EhZk+HiIyHs5FICIqgr67rd+5Y9o4iMj4mAARERWhShX9ynFXdiLLwy4wIiIdYmPFVZqLw13ZiSwXW4CIiArIyhITn27dxPE/3t7ice7KTlSxMAEiIvpPfDwQEgIsXy4+fvdd4OpV4JdfAP9nZrQHBACbN3NXdiJLxS4wIrJ6OTlAVJS4YSkAVK8OrFoFdOkiPuau7EQVDxMgIrJqhw4BERHAlSvi4zFjgK++Ehc4LIi7shNVLOwCIyKr9OSJuDlpu3Zi8uPvD/zxB/D994WTHyKqeNgCRERW5+hRsdXnwgXxcUSEOJjZzU3KqIioPLEFiIisRm4u8Mkn4s7sFy6IM7y2bQPWrGHyQ2Rt2AJERFbh1Cmxpef0afHxkCHA0qWAp6e0cRGRNNgCREQVmkIBTJ8OtGolJj9Vq4rT19evZ/JDZM3YAkREFdaZM2Krz8mT4uMBA4BvvgG8vKSNi4ikxxYgIqpw8vOB2bOBFi3E5MfdXWzx2bSJyQ8RidgCREQVyoULwIgRwJEj4uNevYDvvuOGpUSkjS1ARFQhKJXAggVAs2Zi8uPqKs7u2raNyQ8RFcYWICKyKEpl4S0pkpKAkSOBAwfEMuHhwMqV4n5dRES6MAEiIosREwNMmCDu0q7m5ibu5ZWXBzg7i61Ao0cX3r2diKggJkBEZBFiYoBXXgEEQft4Rob4tVEj4LffgKCg8o6MiCwRxwARkdlTKsWWn2eTn4IePgQCA8svJiKybEyAiMjs7d+v3e2ly40bYjkiIn0wASIis5WfL3Z9jRunX/nUVNPGQ0QVB8cAEZHZSU0VZ3F99x2QkqL/93G6OxHpiwkQEZkFQRC7sL75BvjlF7H1BwCqVQNGjQLWrQPS0nSPA5LJxCnv7duXb8xEZLnMogts2bJlCAoKgoODA1q3bo2jR48WWXbFihVo37493N3d4e7ujrCwsELlR4wYAZlMpnV76aWXTP0yiMgAjx4B334LNG4MdOwIbNwoJj9t2gD/93/i2J45c4CvvxbLPzu9Xf140SJALi/X0InIgkmeAG3cuBGRkZGYOnUqTp48iSZNmiA8PBy3b9/WWT4+Ph5DhgzBnj17kJCQgMDAQHTr1g0pz7STv/TSS0hNTdXcfvrpp/J4OUSkp7NnxbE9/v7AO++IG5c6OQFjxgCnTgEHDwJDhwL29mL5/v3FXdz9/bWvExAgHu/fv/xfAxFZLsm7wBYsWIAxY8Zg5MiRAIDly5dj+/btWLVqFaZMmVKo/I8//qj1eOXKlfjll18QFxeH4cOHa47b29vDx8dHrxhyc3ORm5ureZyZmQkAUCgUUCgUpX5N1kD9vvD9MQ+WUh8KBbBtmwzffWeDvXuf/v9Vp46At95SYdgwFdzcnpZ9Vq9eQI8ewIEDMs1K0O3aCZDLdZeXiqXUh7Uw6/pQKiE7cECztLnQrl2Fb8o0ZX2U5pqSJkB5eXk4ceIEoqKiNMdsbGwQFhaGhIQEva6Rk5MDhUIBDw8PrePx8fHw8vKCu7s7XnzxRXz++efw9PTUeY3Zs2dj+vTphY7v3r0bTk5OpXhF1ic2NlbqEKgAc62Pe/ccsHt3DezeHYQHDxwAADY2Alq1SkX37klo3PguZDLg0CH9r+niAmRnA7t2mShoIzDX+rBW5lYfvgkJCFm5Eo737mmOPfb0xOnRo5EaGiphZOXDFPWRk5Ojd1mZIBS3tJhp3bp1C/7+/jh06BBCC1T25MmTsXfvXhxRb+dcjHfeeQe7du3C2bNn4eAg/mHdsGEDnJycEBwcjKtXr+Ljjz+Gs7MzEhISINeRWetqAQoMDMTdu3fh4uJihFda8SgUCsTGxqJr166wtbWVOhyrJ0V9KJW6W2LUBAHYu1eG5cttsG2bDEqlOFjH21vAqFEqjB6tqrALF/L3w7yYY33ItmyBfPBgQBBQcFib8N+gNuWGDRD69ZMmOBMzZX1kZmaiatWqePjwYYmf35J3gZXFnDlzsGHDBsTHx2uSHwAYPHiw5n5ISAgaN26MWrVqIT4+Hl26dCl0HXt7e9irBxoUYGtraza/LOaK75F5Ka/60LUnV0AAsHgxEBYmztj65hvg/Pmn59u3B959F+jXTwY7OzmAit3MD/D3w9yYTX0olcAHH+ic0igTBEAmQ6UPPwQGDKjQ3WGmqI/SXE/SBKhq1aqQy+VIT0/XOp6enl7i+J2vvvoKc+bMwZ9//onGjRsXW7ZmzZqoWrUqrly5ojMBIiL9FbUnV0qK+Pfa3h5QN6hWrgwMGyYOcg4JKf9YicxSSUubC8LTpc07dSq3sKyNpLPA7Ozs0KJFC8TFxWmOqVQqxMXFaXWJPWvu3LmYOXMmdu7ciZYtW5b4PDdv3sS9e/fgy1XSiMqkuD251Mdyc4F69cRp67duiVPcmfyQVRME4MoV4H//E/8bGDFCv+/j0uYmJXkXWGRkJCIiItCyZUu0atUKixYtQnZ2tmZW2PDhw+Hv74/Zs2cDAL788ktER0dj/fr1CAoKQlpaGgDA2dkZzs7OyMrKwvTp0zFgwAD4+Pjg6tWrmDx5MmrXro3w8HDJXidRRaDPnlyA2P3VubPp4yEyS1lZwNGjQEICcPiweLt7t/TX4T/tJiV5AjRo0CDcuXMH0dHRSEtLQ9OmTbFz5054e3sDAJKTk2Fj87Sh6ttvv0VeXh5eeeUVretMnToV06ZNg1wuxz///IO1a9ciIyMDfn5+6NatG2bOnKlznA8R6SczE3hmFYoi/fd/CZHlUCrFDF89qr99e/3G3wgCcOmSmOQkJIi3M2cAlUq7nJ0d0KIFEBoKtGoFvP9+0UubA4CzM9C2bdlfFxVJ8gQIAMaNG4dxRex2GB8fr/X4+vXrxV7L0dERu8x5XiyRBREE4MABYNUq4OefAX1nmPIfV7IoxY3qf3aFzczMwq079+8Xvmb16sALL4gJT2go0LTp01U9AcDWVhxMJ5PpToKyssSVQNetAwpM8iHjMYsEiIjMy61b4t/dVauAy5efHn/uOSA9XfwM4J5cVCEUN6r/lVfEPVaqVHnaunP2bOGy9vZAy5ZioqNOevz8in9e9dLmzyZegYHibIJly4BNm8QWqW3bgGfWuqOyYwJERACAvDxg+3Yx6dmx42kLvrMzMGiQuCFpaCiwZYvuf1y5JxdZHH1G9U+YUPhcUJB2606TJmIXV2n17w/06aO76613b6BvX7EJtm1b4I8/xOclo2ECRGTlzp0Tk55164A7d54eb9dOTHoGDhSTILWi/nENCBCTH+7JRWZPEIDkZHFWlj6j+kNCgO7dxaTnhReM28crl+ue6t65s7ghXvfuwIULYqK1fTvQvLnxntvKMQEiskKZmeKu66tWiUMY1Hx8gIgIYORIoG7dor+/uH9cicxKfj5w5gwC4uNhEx8P/P03kJgIPHig/zWiooAhQ0wVYdEaNRJ/QXv0AP75B+jQQfzv46WXyj+WCogJEJGVEAQxYVm1ShxaoB7QLJcDL78MvPGG+M9mJT3/KhT1jytRmRg6GwsQN4f75x8xwTl1Svx6+jRsnzxBi2fLVqoE1KgBXL1a8nWlHNXv7w/s2yeOC4qLE39Zv/tO/IWlMmECRGThlEpxz619+/xRubIMnTtrf16kpDwd0HzlytPj9eqJXVzDhoktP0SSK81srDt3niY66mTn0qXC088BCM7OuB8YCLfOnSFv0UKckdWwoZgEBQWJvyTmPKrf1VUcmDd6tNhtN3q0uFL01KlPB99RqTEBIrJgTz8vKgFoiQULxL/X8+eLf9tXrRLHTj47oPmNN8ShDPzbSUUqS0uMIUqajfXhh+JsK3XSk5Ki+zo+PkCzZmKS89/X/OrVcWDnTvTo0QPyZ/eKWrzYMkb129kBa9eK0+u/+AKYPl0cx/Tdd+KUeio1JkBEFqqoz4ubN8Ukp6CiBjQT6VSalpiyEARxvZv0dODtt4ufjTVvXuFzdeoUSnZ0NmcqFEXHYEmj+mUy4PPPxSTo7beB1avFNSs2bRKn6lOpMAEiskDFzd5Vs7ERN5wePVpcv4dILyW1xGzerDspEARxYNndu9q3e/eKP5aXp39sPXqIA9WaNQMaNzbeh76ljeofO1ZcZ2jQIGDXLqBjR3GGGFcgLRUmQEQW5MEDcVLI+vUlz95VqcTPCyY/pDd91sUZORKIjRV/GJ9NbHJzDXteW9viW2nUXn/ddLOxLG1U/8svA/HxQM+eYpfgCy8AO3cC9etLHZnFYAJEZKZUKuDiReDQIXEB2kOHgPPnS3cNbiZNesvIEAeNlZRZZ2YCy5cXfd7eHqhatfDN07PoY0eP6rd7Lls4tD3/vPjHoXt3ccn2Nm2AX3+VftC2hWACRGQkZR0z+uiR+DmgTngSEsTPpGfVqQMEBwO7d5d8TX5ekE4PHgAnTwInTjy96TMdXK1fP7G1RFdi4+RU+tH17duLY27MfTaWOapVS/yj0bu3+EcjLEycKfbqq1JHZvaYABEZQWnHjAqC+HmjbtlJSABOny48g9fRUdw4OjRU/OfuhReAatXEZMsSZu+SEZQ1s37wQDvROXECuHZNd1lvb3FAcknGjzdud5FcbjmzscxR1ariGkFDh4p71QwaJP4xev99TvUsBhMgojLSZ8zoSy8Bx49rJzwFt51Qq1FDTHTUCU/jxrpnuPLzQgLlPS0cKH1mff9+4WQnKUn3tYODgRYtnt6aNwfc3KTLrC1pNpY5cnQUZ4O9/z6wdKk4A+Lff4EFC/iHoAhMgIjKQJ8xo4MHiy07SqX2eTs78XNHnfDos4F0Qfy8KEflNS382ecsLrNevVpMxAomO9ev675WzZqFk52idheXMrO2tNlY5kb9n1GNGuK6SUuWiD+z//d/YoJEWpgAEZXB/v0ljxlVT27x9dVu3WneXBwvWhbqz4s9e/Lxxx+J6N69KTp3rsTPC2MydFp4WeTni91MxWXWI0bo/t5atQonO+7u+j+31Jm1pc3GMjcymdj6ExAADB8u/vyGhYmDoz09pY7OrDABIjJQcjKwZo1+ZRcvBt57zzTd8XI50LGjgOzsFHTs2KTiJz/l2RVVUhOfTAaMGyeuNfDkibgOTnY2ZA8fIvDwYdjcuCFODc/OFm//ndd8LepYdnbxizyp+fmJr//ZbqyyYkuM5Rs0SKy3Pn3Efvc2bcRl4WvWlDoys8EEiEhPgiDus7htG7B1q7j0hr4aN+ZYRKMwRVeUIIhTu+/cEW+3bz+9f+JE8U18giAmCCEhWocrAWhuWDSl89VXXBeHitahA3DwoDhN/tIlsfl5+3agZUupIzMLTICIipGfL/4TvG2beCs4xEImE/+pOnNG/Py0utlY5rZXlLorShDENQUKJjIl3S/NasS6VK4sjqlxcgIqV4bKyQl3srNRrUYN2Dg7i+f/O6d1v6hjf/8tvqaScJ0DKkmDBuKsi549xX3UOnYUB0uHh1t9Cx8TIKJnZGWJq8tv2yb+s3T//tNzDg5At25iq/LLLwNeXk8/l61qNlZ5DwrOzRW7mkoabe7lZXhC4+wsrjFQrZp4nWrVgMePgQ0bSv7e33/Xai1RKhQ4vGMHevToARtDNqoMDua6OGQ8fn7Avn3iH6rdu4FevcQd5h88eFrG1IP6zRATICKIS5/89pvYtfXnn9or+nt6in8v+vQBunYV/0kvSOoxo+XOWIOC1RthpqYCaWniraj7t2+XPCZGodDeIbxyZe1kRv21qPu6ZskolcCBA+WfiHCdAzK2KlXERD08HNizRzv5AUw7qL8gpRKyvXvhv28fZJUriyuAS/RzzASIKhx9e2YuXnw6nufwYe3PmJo1gb59xaSnTRugUgm/KZKOGS3PPyj6DAqeMEFcov/27ZITm5wc48b3+efiflHVqoldSWUlZSJidZk1mZyNjbhlhi4Ff3/79DHNz/R/LceVbt5ES0Bco0jClicmQFShFNcz07evuNXE1q1i4nPhgvb3tmwp/t737Qs0bFj6QcuSjBktrz8oSqW42eXvv5c8KPjmTaB6df2vXaUK4OMjZo0+PoXv+/iIKxcPGFDytdq2FddAMSYpExHOxiJjKmndDvXvb5Uq4s+aupW0atWn9wve1MednUv+gynFchIlYAJEFUZRv183b4qfnW5u2ntrVaoEvPii+PnSu7f4eVYm5joouCgKhdhKk54u3gref/Z2927hfTqKI5MVndA8m9w826eoS0iItGNipExEOBuLjEXf3ZEfPxb/6Shqy5Rn2dvrTozUNw8P4J13im85njjRdC1PRWACRBVCcT0zahkZ4j82PXqIrTzdu4vjAI2ivAcF67ME9ZtviomLerbTs0lNwdHd+pDJxDcwM7PksrGxQJcupbt+ccxhTAwTEbJ0+s4aXLtWXFBTPVuy4E39N0V9e/JEHDR582bJq8IWRRCAGzfEfzDK8XeMCRBZpLw8cQzP2bPiNPT4eP1+9375RRzIbFTGbtrNzxeztYwMcaDigweF7589W/ILvntXTIKKI5eL/6F5exd/Uw8Ulsn02yvKFH/EOCaGqGzat9evJXXoUP3+mRAEcdFOXYlRwYTp/HngypWSr6dvC5WRMAEikzFGj1B+vjhm7+zZp8nO2bPisfx83d9jAyXaYz98kYpU+GI/2kMF8Ynv3i3ji3qWPoOC335b7G/LzNSdzDx7/9Ej48XXtCnQrFnRSY2npzgwsjS4VxSRZTJ2S6pMJo7/cXYWl24oSny8ODmjJOW8rhUTIDKJ0vYIKZXiptXqBEf99eLFopd0cXERBys3aiRuLLpsGdAPMViMCQjE0ye+gQBMwGJsQf+y/X6pVMC9e8CtW09vBw+WPKjw9m3xQ7u0nJ3FPZzc3MSvBe9nZOi3D8fChcZvjZG6JYZdUUSGk+L3V9+Wp3Je14oJEBldST1C33wj/qwXbNU5f17sStalcmVxMdNGjcSER530+Ps//adFqQRyN8Tgu3uvANB+Yn+kYDNewVuem9G+vY5fbkEAHj7UTmyevaWkiC0O6p1NSys4WJxb/2wio+uxm5t4K24BPaVSXLDIGgcFE1HZlPfvrzmM4dOBCRAZVcEeoUJdUYLYFfX227q/18EBqF+/cKJTvXrJvTRyKLEYEwAIeLaoDQSoACzNHQP5/Cvi+jPPJjiPH+v/IqtVE1dW9fMTf3l37Cj5e1atMm6rhTn8QWFLDJHlKu/fX6lbjnVgAmQFlHlKnP5mP3KupsKpli9C3mkPuV3ZPxjz82W4fBn4919xfNvly8CRI+LPdkldUcHBQOvW2slOzZp6fF7n5OiejXD8OJzuFd0VZQPAPus+8NFHRV/b3f1pYlPUzcdH7G9TUyr1GxRsipYYM/yDQkRUpP9anvL37EHiH3+gaffuqMSVoK1Eea8TA+Dw5BhUXzABTZVPPyBvfRiA5MjFeGFuyR+QeXniBqCXL4tJjjrRuXKlEpKSXoZKVbhpph9isBlFd0W9gs0Y+EV/DBksiGNZ7twBbt8Bzt0tfrrl3btlXzm4bVtxaednExtfX91bIZRE6pYYM/uDQkRULLkcQseOSMnORpOOHSX9W8UEqLzExECYMAGyAv+pCwEBkJlwCfDDk2PQal7hRMRHmQKfea/gMDbjhbn9kZcnrnelneCIt+vXda9/J4MARzxGVYcsNKiRg7oB2ajllwPbx5kYsHksZBDw7LqgNhAgANiAwcA77sDw+0VP5SqOnV3hxbZyc8XBRyX5/PMKOSjYXP6gEBFZCiZA5SEmBsKAVyA8kxQIN1OAAa9A9ovxlwBXPlGgxvzxKG5MTJ15Y/DZylQ8yXgMRyEHlZGNyshGS2Sj43/3nZCDKjbZcLfNRhV5NpyEbNgrc2Cb919LzBMAF/+76UEGwA4KIOP204NVquhePbSox1WqFF52XcquKICDgomILAwTIFNTKpEzdgIcikxEZMgd9Q4cPD3FVoycnKe3x48191VZOcjLyIEi8zEUmTlQPsqBkJUDIScHssePYfMkB/K8HFTKy4Fdfg7sVU9Q3IxvGwCeuI+ZD8aV/BpUAHKLOe/oKE7VcnISW3Ru3Sr5mrNmAcOGiUmNg0PJ5UsidVeUOgYOCiYisghmkQAtW7YM8+bNQ1paGpo0aYKlS5eiVatWRZbftGkTPvvsM1y/fh116tTBl19+iR49emjOC4KAqVOnYsWKFcjIyEDbtm3x7bffok6dOuXxcrQo4/eXMDhXgMPD9BI/OG0AOPx3M6abfs/Ds009OHg4QeZcWUxk1MlM5cratwLHFHZ22LV/P8L79YOtvf3TC+q74FVoqBE233qG1F1RRERkMSRPgDZu3IjIyEgsX74crVu3xqJFixAeHo6LFy/Cy8urUPlDhw5hyJAhmD17Nl5++WWsX78effv2xcmTJ9GoUSMAwNy5c7FkyRKsXbsWwcHB+OyzzxAeHo5z587BwRitDaVwMT4VDfQodws+uA1v5MCp0O0xHDX3822dACcnyCo7Qe7sCHkVJ9i6OsHOzQn27k5w9BRvslMnMWhDvxKf9+6kuQiY2Kn0L0yhgNLRsfD8dKkXvGJXFBER6UHyBGjBggUYM2YMRo4cCQBYvnw5tm/fjlWrVmHKlCmFyi9evBgvvfQSJk2aBACYOXMmYmNj8fXXX2P58uUQBAGLFi3Cp59+ij7/rb67bt06eHt7Y+vWrRg8eHD5vTgAqfDVKwGa3/wn2Id30qyB5+YGeLgBtQqsi+fqqj0DuzjKPH/c2hQAH2UKbFA4EVFBhlR5AELeMXIiwq4oIiKyAJImQHl5eThx4gSioqI0x2xsbBAWFoaEhASd35OQkIDIyEitY+Hh4di6dSsAICkpCWlpaQgLC9Ocd3V1RevWrZGQkKAzAcrNzUVu7tNBLpn/7XatUCigMHTl3/8I7UJxAwHwR9GJyE0EoPusUHR8seTn0jscGZA0YQF8FgyCCjKt51b9NxT7+oT58JKpoFLomOZVYhwKra9aevWCbMMGyCMjIUtJ0RwW/P2hnD8fQq9ehq+oTDoVWx9U7lgf5oX1YV5MWR+luaakCdDdu3ehVCrh7e2tddzb2xsXLlzQ+T1paWk6y6elpWnOq48VVeZZs2fPxvTp0wsd3717N5ycnPR7MUVQKoFPnediddbQIhOR6Cpfonf2Tr0WFC6VDnb48dpsvPjrUvirniYit2z8saf3OLh0sMOOMj5pbGys7hP29sCSJfA8dw4ODx7gibs77jVoILbOGP2FklqR9UGSYH2YF9aHeTFFfeSUYq04ybvAzEFUVJRWq1JmZiYCAwPRrVs3uLi4lPn6+T/IMHCQHRZhotbKyDcRgPexEINW9kGvXjrGyxhDjx5Q5k3AieUHNStBN3yrLQaXcSVohUKB2NhYdO3aFbbF7VnVq1eZnof0o3d9ULlgfZgX1od5MWV9qHtw9CFpAlS1alXI5XKkp6drHU9PT4ePj4/O7/Hx8Sm2vPpreno6fAts/Z2eno6mTZvqvKa9vT3sC85k+o+tra1RKufVV4FKlQag3fi+CE55ujfW9YD2WLBYbvLJSba2tmjxQReTXZt/UMwH68O8sD7MC+vDvJiiPkpzvRK2mDQtOzs7tGjRAnFxcZpjKpUKcXFxCA0N1fk9oaGhWuUBsRlNXT44OBg+Pj5aZTIzM3HkyJEir1ke+vcHrv0rx7Q9ndB7/RBM29MJV6+bPvkhIiKiwiTvAouMjERERARatmyJVq1aYdGiRcjOztbMChs+fDj8/f0xe/ZsAMCECRPQsWNHzJ8/Hz179sSGDRtw/PhxfP/99wAAmUyGiRMn4vPPP0edOnU00+D9/PzQt29fqV4mAE5OIiIiMheSJ0CDBg3CnTt3EB0djbS0NDRt2hQ7d+7UDGJOTk6GTYG1Ztq0aYP169fj008/xccff4w6depg69atmjWAAGDy5MnIzs7G2LFjkZGRgXbt2mHnzp3lvgYQERERmSfJEyAAGDduHMaN070lQ3x8fKFjAwcOxMCBA4u8nkwmw4wZMzBjxgxjhUhEREQViKRjgIiIiIikwASIiIiIrA4TICIiIrI6TICIiIjI6jABIiIiIqvDBIiIiIisDhMgIiIisjpmsQ6QuREEcWPS0myqZm0UCgVycnKQmZnJvXXMAOvDvLA+zAvrw7yYsj7Un9vqz/HiMAHS4dGjRwCAwMBAiSMhIiKi0nr06BFcXV2LLSMT9EmTrIxKpcKtW7dQpUoVyGQyqcMxS5mZmQgMDMSNGzfg4uIidThWj/VhXlgf5oX1YV5MWR+CIODRo0fw8/PT2kZLF7YA6WBjY4OAgACpw7AILi4u/INiRlgf5oX1YV5YH+bFVPVRUsuPGgdBExERkdVhAkRERERWhwkQGcTe3h5Tp06Fvb291KEQWB/mhvVhXlgf5sVc6oODoImIiMjqsAWIiIiIrA4TICIiIrI6TICIiIjI6jABIiIiIqvDBIg0Zs+ejeeffx5VqlSBl5cX+vbti4sXL2qVefLkCd599114enrC2dkZAwYMQHp6ulaZ5ORk9OzZE05OTvDy8sKkSZOQn59fni+lwpkzZw5kMhkmTpyoOca6KH8pKSl4/fXX4enpCUdHR4SEhOD48eOa84IgIDo6Gr6+vnB0dERYWBguX76sdY379+9j6NChcHFxgZubG9544w1kZWWV90uxeEqlEp999hmCg4Ph6OiIWrVqYebMmVp7QLE+TGffvn3o1asX/Pz8IJPJsHXrVq3zxnrv//nnH7Rv3x4ODg4IDAzE3LlzjfciBKL/hIeHC6tXrxbOnDkjJCYmCj169BCqV68uZGVlacq89dZbQmBgoBAXFyccP35ceOGFF4Q2bdpozufn5wuNGjUSwsLChFOnTgk7duwQqlatKkRFRUnxkiqEo0ePCkFBQULjxo2FCRMmaI6zLsrX/fv3hRo1aggjRowQjhw5Ily7dk3YtWuXcOXKFU2ZOXPmCK6ursLWrVuFv//+W+jdu7cQHBwsPH78WFPmpZdeEpo0aSIcPnxY2L9/v1C7dm1hyJAhUrwki/bFF18Inp6ewu+//y4kJSUJmzZtEpydnYXFixdryrA+TGfHjh3CJ598IsTExAgAhC1btmidN8Z7//DhQ8Hb21sYOnSocObMGeGnn34SHB0dhe+++84or4EJEBXp9u3bAgBh7969giAIQkZGhmBrayts2rRJU+b8+fMCACEhIUEQBPGXwsbGRkhLS9OU+fbbbwUXFxchNze3fF9ABfDo0SOhTp06QmxsrNCxY0dNAsS6KH8fffSR0K5duyLPq1QqwcfHR5g3b57mWEZGhmBvby/89NNPgiAIwrlz5wQAwrFjxzRl/vjjD0EmkwkpKSmmC74C6tmzpzBq1CitY/379xeGDh0qCALrozw9mwAZ673/5ptvBHd3d62/Vx999JFQt25do8TNLjAq0sOHDwEAHh4eAIATJ05AoVAgLCxMU6ZevXqoXr06EhISAAAJCQkICQmBt7e3pkx4eDgyMzNx9uzZcoy+Ynj33XfRs2dPrfccYF1I4ddff0XLli0xcOBAeHl5oVmzZlixYoXmfFJSEtLS0rTqxNXVFa1bt9aqEzc3N7Rs2VJTJiwsDDY2Njhy5Ej5vZgKoE2bNoiLi8OlS5cAAH///TcOHDiA7t27A2B9SMlY731CQgI6dOgAOzs7TZnw8HBcvHgRDx48KHOc3AyVdFKpVJg4cSLatm2LRo0aAQDS0tJgZ2cHNzc3rbLe3t5IS0vTlCn4gas+rz5H+tuwYQNOnjyJY8eOFTrHuih/165dw7fffovIyEh8/PHHOHbsGMaPHw87OztERERo3lNd73nBOvHy8tI6X6lSJXh4eLBOSmnKlCnIzMxEvXr1IJfLoVQq8cUXX2Do0KEAwPqQkLHe+7S0NAQHBxe6hvqcu7t7meJkAkQ6vfvuuzhz5gwOHDggdShW6caNG5gwYQJiY2Ph4OAgdTgE8Z+Cli1bYtasWQCAZs2a4cyZM1i+fDkiIiIkjs76/Pzzz/jxxx+xfv16NGzYEImJiZg4cSL8/PxYH6QXdoFRIePGjcPvv/+OPXv2ICAgQHPcx8cHeXl5yMjI0Cqfnp4OHx8fTZlnZyKpH6vLUMlOnDiB27dvo3nz5qhUqRIqVaqEvXv3YsmSJahUqRK8vb1ZF+XM19cXDRo00DpWv359JCcnA3j6nup6zwvWye3bt7XO5+fn4/79+6yTUpo0aRKmTJmCwYMHIyQkBMOGDcP777+P2bNnA2B9SMlY772p/4YxASINQRAwbtw4bNmyBX/99VehpscWLVrA1tYWcXFxmmMXL15EcnIyQkNDAQChoaE4ffq01g92bGwsXFxcCn14UNG6dOmC06dPIzExUXNr2bIlhg4dqrnPuihfbdu2LbQsxKVLl1CjRg0AQHBwMHx8fLTqJDMzE0eOHNGqk4yMDJw4cUJT5q+//oJKpULr1q3L4VVUHDk5ObCx0f4Ik8vlUKlUAFgfUjLWex8aGop9+/ZBoVBoysTGxqJu3bpl7v4CwGnw9NTbb78tuLq6CvHx8UJqaqrmlpOToynz1ltvCdWrVxf++usv4fjx40JoaKgQGhqqOa+eet2tWzchMTFR2Llzp1CtWjVOvTaCgrPABIF1Ud6OHj0qVKpUSfjiiy+Ey5cvCz/++KPg5OQk/N///Z+mzJw5cwQ3Nzdh27Ztwj///CP06dNH59TfZs2aCUeOHBEOHDgg1KlTh9OuDRARESH4+/trpsHHxMQIVatWFSZPnqwpw/ownUePHgmnTp0STp06JQAQFixYIJw6dUr4999/BUEwznufkZEheHt7C8OGDRPOnDkjbNiwQXBycuI0eDI+ADpvq1ev1pR5/Pix8M477wju7u6Ck5OT0K9fPyE1NVXrOtevXxe6d+8uODo6ClWrVhU++OADQaFQlPOrqXieTYBYF+Xvt99+Exo1aiTY29sL9erVE77//nut8yqVSvjss88Eb29vwd7eXujSpYtw8eJFrTL37t0ThgwZIjg7OwsuLi7CyJEjhUePHpXny6gQMjMzhQkTJgjVq1cXHBwchJo1awqffPKJ1pRp1ofp7NmzR+fnRUREhCAIxnvv//77b6Fdu3aCvb294O/vL8yZM8dor0EmCAWWzSQiIiKyAhwDRERERFaHCRARERFZHSZAREREZHWYABEREZHVYQJEREREVocJEBEREVkdJkBERERkdZgAERERkdVhAkRERERWhwkQEZWLESNGQCaTYc6cOVrHt27dCplMJlFUpiGTybB161apwyCiYjABIqJy4+DggC+//BIPHjyQOhQisnJMgIio3ISFhcHHxwezZ88uttwvv/yChg0bwt7eHkFBQZg/f36pnufevXsYMmQI/P394eTkhJCQEPz0009aZTZv3oyQkBA4OjrC09MTYWFhyM7OBgDEx8ejVatWqFy5Mtzc3NC2bVv8+++/mu/dtm0bmjdvDgcHB9SsWRPTp09Hfn4+ACAoKAgA0K9fP8hkMs3jv//+G507d0aVKlXg4uKCFi1a4Pjx46V6XURkPEyAiKjcyOVyzJo1C0uXLsXNmzd1ljlx4gReffVVDB48GKdPn8a0adPw2WefYc2aNXo/z5MnT9CiRQts374dZ86cwdixYzFs2DAcPXoUAJCamoohQ4Zg1KhROH/+POLj49G/f38IgoD8/Hz07dsXHTt2xD///IOEhASMHTtW0023f/9+DB8+HBMmTMC5c+fw3XffYc2aNfjiiy8AAMeOHQMArF69GqmpqZrHQ4cORUBAAI4dO4YTJ05gypQpsLW1NfStJKKyMtq+8kRExYiIiBD69OkjCIIgvPDCC8KoUaMEQRCELVu2CAX/FL322mtC165dtb530qRJQoMGDcr0/D179hQ++OADQRAE4cSJEwIA4fr164XK3bt3TwAgxMfH67xOly5dhFmzZmkd+9///if4+vpqHgMQtmzZolWmSpUqwpo1a8r0GojIeNgCRETl7ssvv8TatWtx/vz5QufOnz+Ptm3bah1r27YtLl++DKVSqdf1lUolZs6ciZCQEHh4eMDZ2Rm7du1CcnIyAKBJkybo0qULQkJCMHDgQKxYsUIzLsnDwwMjRoxAeHg4evXqhcWLFyM1NVVz7b///hszZsyAs7Oz5jZmzBikpqYiJyenyJgiIyMxevRohIWFYc6cObh69aper4WITIMJEBGVuw4dOiA8PBxRUVEmuf68efOwePFifPTRR9izZw8SExMRHh6OvLw8AGJXXGxsLP744w80aNAAS5cuRd26dZGUlARA7L5KSEhAmzZtsHHjRjz33HM4fPgwACArKwvTp09HYmKi5nb69GlcvnwZDg4ORcY0bdo0nD17Fj179sRff/2FBg0aYMuWLSZ5/URUMiZARCSJOXPm4LfffkNCQoLW8fr16+PgwYNaxw4ePIjnnnsOcrlcr2sfPHgQffr0weuvv44mTZqgZs2auHTpklYZmUyGtm3bYvr06Th16hTs7Oy0EpJmzZohKioKhw4dQqNGjbB+/XoAQPPmzXHx4kXUrl270M3GRvyTamtrq7O16rnnnsP777+P3bt3o3///li9erVer4eIjI8JEBFJIiQkBEOHDsWSJUu0jn/wwQeIi4vDzJkzcenSJaxduxZff/01PvzwQ02ZLl264Ouvvy7y2nXq1EFsbCwOHTqE8+fP480330R6errm/JEjRzBr1iwcP34cycnJiImJwZ07d1C/fn0kJSUhKioKCQkJ+Pfff7F7925cvnwZ9evXBwBER0dj3bp1mD59Os6ePYvz589jw4YN+PTTTzXXDwoKQlxcHNLS0vDgwQM8fvwY48aNQ3x8PP79918cPHgQx44d01yTiCQg9SAkIrIOBQdBqyUlJQl2dnbCs3+KNm/eLDRo0ECwtbUVqlevLsybN0/rfI0aNYSpU6cW+Vz37t0T+vTpIzg7OwteXl7Cp59+KgwfPlzz/OfOnRPCw8OFatWqCfb29sJzzz0nLF26VBAEQUhLSxP69u0r+Pr6CnZ2dkKNGjWE6OhoQalUaq6/c+dOoU2bNoKjo6Pg4uIitGrVSvj+++8153/99Vehdu3aQqVKlYQaNWoIubm5wuDBg4XAwEDBzs5O8PPzE8aNGyc8fvzYgHeSiIxBJgiCIHUSRkRERFSe2AVGREREVocJEBEREVkdJkBERERkdZgAERERkdVhAkRERERWhwkQERERWR0mQERERGR1mAARERGR1WECRERERFaHCRARERFZHSZAREREZHX+H0Atopqici2yAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fix, ax = plt.subplots()\n",
    "ax.plot(problem_dims, time_sigma, c=\"blue\", marker=\"o\", label=\"sigma\")\n",
    "ax.plot(problem_dims, time_factor, c=\"red\", marker=\"o\", label=\"factor\")\n",
    "ax.legend()\n",
    "plt.xlabel(\"No. assets\")\n",
    "plt.ylabel(\"Solution time (sec.)\")\n",
    "plt.grid(True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ed3cf872",
   "metadata": {},
   "source": [
    "From the data, it is clear that the solve times of (2) are much smaller both in absolute and asymptotic terms.  The data may appear somewhat nonsmooth over the problem dimensions.  These are artifacts from multi-threading and different code paths in the computational kernels used for dense matrix multiplication."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "afddf4c0",
   "metadata": {},
   "source": [
    "## Takeaways\n",
    "\n",
    "- Incorporating explicit factors to model risk in a gurobipy model is straightforward.\n",
    "- Using factor risk models instead of the derived full covariance matrix can greatly improve solution time."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b32f8363",
   "metadata": {},
   "source": [
    "## References\n",
    "\n",
    "[1] Menchero, J., Orr, D., Wang, J.: The Barra US equity model (USE4), methodology notes. English,\n",
    "MSCI (May (2011)\n",
    "\n",
    "[2] Gérard Cornuéjols, Javier Peña, and Reha Tütüncü. Optimization Methods in Finance. Cambridge University Press, 2 edition, 2018. doi:10.1017/9781107297340."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}