finance/simulation/iron_condor.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# _*Pricing Iron Condor Option*_ \n",
    "\n",
    "The latest version of this notebook is available on https://github.com/Qiskit/qiskit-iqx-tutorials.\n",
    "\n",
    "***\n",
    "### Contributors\n",
    "Stefan Woerner<sup>[1]</sup>, Daniel Egger<sup>[1]</sup>\n",
    "### Affliation\n",
    "- <sup>[1]</sup>IBMQ"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Introduction\n",
    "<br>\n",
    "Suppose a <a href=\"http://www.theoptionsguide.com/iron-condor.aspx\">iron condor option</a> with strike prices $K_1 < K_2 < K_3 < K_4$ and an underlying asset whose spot price at maturity $S_T$ follows a given random distribution.\n",
    "The corresponding payoff function is defined as:\n",
    "<br>\n",
    "<br>\n",
    "$$ F(S_T) = \n",
    "\\begin{cases}\n",
    "0         ,& S_T < K_1 \\\\\n",
    "S_T - K_1 ,& K_1 \\leq S_T < K_2 \\\\\n",
    "K_2 - K_1 ,& K_2 \\leq S_T < K_3 \\\\\n",
    "K_3 - S_T ,& K_3 \\leq S_T < K_4 \\\\\n",
    "0         ,& S_T \\geq K_4. \n",
    "\\end{cases}$$\n",
    "<br>\n",
    "In the following, a quantum algorithm based on amplitude estimation is used to estimate the expected payoff, i.e., the fair price before discounting, for the option:\n",
    "<br>\n",
    "<br>\n",
    "$$\\mathbb{E}\\left[ F(S_T) \\right].$$\n",
    "<br>\n",
    "The approximation of the objective function is explained in detail in the following paper:<br>\n",
    "<a href=\"https://arxiv.org/abs/1806.06893\">Quantum Risk Analysis. Woerner, Egger. 2018.</a>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "import numpy as np\n",
    "\n",
    "from qiskit import BasicAer\n",
    "from qiskit.aqua.algorithms import AmplitudeEstimation\n",
    "from qiskit.aqua.components.uncertainty_models import LogNormalDistribution\n",
    "from qiskit.aqua.components.uncertainty_problems import UnivariateProblem\n",
    "from qiskit.aqua.components.uncertainty_problems import UnivariatePiecewiseLinearObjective as PwlObjective"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Uncertainty Model\n",
    "\n",
    "We construct a circuit factory to load a log-normal random distribution into a quantum state.\n",
    "The distribution is truncated to a given interval $[low, high]$ and discretized using $2^n$ grid points, where $n$ denotes the number of qubits used.\n",
    "The unitary operator corresponding to the circuit factory implements the following: \n",
    "$$\\big|0\\rangle_{n} \\mapsto \\big|\\psi\\rangle_{n} = \\sum_{i=0}^{2^n-1} \\sqrt{p_i}\\big|i\\rangle_{n},$$\n",
    "where $p_i$ denote the probabilities corresponding to the truncated and discretized distribution and where $i$ is mapped to the right interval using the affine map:\n",
    "$$ \\{0, \\ldots, 2^n-1\\} \\ni i \\mapsto \\frac{high - low}{2^n - 1} * i + low \\in [low, high].$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# number of qubits to represent the uncertainty\n",
    "num_uncertainty_qubits = 3\n",
    "\n",
    "# parameters for considered random distribution\n",
    "S = 2.0 # initial spot price\n",
    "vol = 0.4 # volatility of 40%\n",
    "r = 0.05 # annual interest rate of 4%\n",
    "T = 40 / 365 # 40 days to maturity\n",
    "\n",
    "# resulting parameters for log-normal distribution\n",
    "mu = ((r - 0.5 * vol**2) * T + np.log(S))\n",
    "sigma = vol * np.sqrt(T)\n",
    "mean = np.exp(mu + sigma**2/2)\n",
    "variance = (np.exp(sigma**2) - 1) * np.exp(2*mu + sigma**2)\n",
    "stddev = np.sqrt(variance)\n",
    "\n",
    "# lowest and highest value considered for the spot price; in between, an equidistant discretization is considered.\n",
    "low  = np.maximum(0, mean - 3*stddev)\n",
    "high = mean + 3*stddev\n",
    "\n",
    "# construct circuit factory for uncertainty model\n",
    "uncertainty_model = LogNormalDistribution(num_uncertainty_qubits, mu=mu, sigma=sigma, low=low, high=high)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaAAAAE0CAYAAABqwecMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xu8ZXP9x/HXO8JkGJPLkGSQSJSaqUyUGUT49SNiKvWLlPRL+pVLUjGRcgn5kTQpk27ThZ9yi3EZcme6TcbIYMgowgzmai6f3x/fdVizZu+z9z57n72WOe/n47Ef56zv+q61PnufffZnr/X9ru9XEYGZmVm3vaLsAMzMbGByAjIzs1I4AZmZWSmcgMzMrBROQGZmVgonIDMzK4UTkJmZlcIJyMzMSlF6ApK0jaTrJc2X9LikkySt0mCbN0n6fVZ/kaRHJV0oaaNCvQmSosZj6/59VmZm1siqZR5c0lDgOmAasA+wBXAmKTF+tZdNhwAPAxcDjwObAScCIyS9PSKW5OpOBw4pbD+zmfjWW2+9GD58eDNV+8W8efNYc801Szt+PVWNCxxbX1Q1LnBsfVF2XFOmTHkqItZvqnJElPYAvgzMBtbOlR0LzM+XNbmv9wIBvC1XNgG4p6/xjRgxIsp04403lnr8eqoaV4Rj64uqxhXh2Pqi7Lha+cwt+xLcnsA1EfFcrmwiMAjYucV9PZ39XK0TgZmZWf8qOwFtTbpE9qKIeJR0BtSwnUbSKyStJmkr4FTgbuCuQrVtJD2XtRXdIqnVxGZmZv1AUeJo2JIWA8dExHcK5Y8BF0fE8Q22/z2wR7Y4BdgrIp7Mrf888AKpjWl94ChgBLBTRBQTVc82hwGHAQwbNmzExIkT+/LUOmLu3LkMHjy4tOPXU9W4wLH1RVXjAsfWF2XHNWbMmCkRMbKpys1eq+uPB7AY+HyN8lnAKU1svyXwTuCjpDOpKcAavdQfROq8cFkz8bkNqLaqxhXh2PqiqnFFOLa+KDsuXkZtQLOBdWqUDwHmNNo4Ih6IiDsj4qekM6G3Ah/ppf4C4CrgbX0L18zMOqXsBDSdQluPpE2ANSm0DTUSEY8AzwCbN1O9lX2bmVnnlZ2Argb2kLRWrmwssAC4qZUdZR0R1iVdYqtXZxCp592U1kM1M7NOKvVGVOAC4EjgUkmnkc5exgFnRa5rtqQZwE0RcWi2/G1gCXAn6VLdG0n3Dz1I6saNpCHAFcBPgRnAesAXgI2BA7vw3MzMrBelJqCImC1pV+A84HJSMjmblITyVgXyw/PcA3yO1FttDeBR4BLgWxExL6uzCPg3aUSFDYCFwO3AzhFxT388HzMza17ZZ0BExDRglwZ1hheWJ5Kd6fSyzUJgv3bjMwMYftyVbe/jqO2WcHCb+5l56t5tx2FWFWW3AZmZ2QDlBGRmZqVwAjIzs1I4AZmZWSmcgMzMrBROQGZmVgonIDMzK4UTkJmZlcIJyMzMSuEEZGZmpXACMjOzUjgBmZlZKZyAzMysFE5AZmZWCicgMzMrhROQmZmVwgnIzMxKUXoCkrSNpOslzZf0uKSTJK3SYJs3Sfp9Vn+RpEclXShpoxp195E0VdJCSdMkje2/Z2NmZs0qdUpuSUOB64BpwD7AFsCZpMT41V42HQI8DFwMPA5sBpwIjJD09ohYku1/J+AS4HzgSGAv4BeSZkfEtf3ypMzMrCmlJiDgcGAQsF9EPAdMkrQ2ME7S6VnZCiLiNuC2XNFkSY8B1wJvBv6YlX8NuDkijsyWb5T0JuCErK6ZmZWk7EtwewLXFBLNRFJS2rnFfT2d/VwNQNLqwBjgV4V6E4FRkoa0Hq6ZmXVK2Qloa2B6viAiHgXmZ+t6JekVklaTtBVwKnA3cFe2egvglcX9A/eRnvcb2gvdzMzaUXYCGgrMqVE+O1vXyFXAIlKSeTXwHxGxLLdvaux/dmG9mZmVQBFR3sGlxcDREXFOoXwWMCEivtJg+y1JiWdLUqeFecCOEbFQ0o7ALcD2EfGXwjZ/B3aPiEk19nkYcBjAsGHDRkycOLGdp9iWuXPnMnjw4NKOX09V44L+i23qrGfb3sewQfDEgvb2sd3Gnb9yPBD/np1Q1djKjmvMmDFTImJkM3XL7oQwG1inRvkQap8ZLSciHsh+vVPSH0g94z4C/IiXznSK++9Zrrn/iBgPjAcYOXJkjB49ulEY/Wby5MmUefx6qhoX9F9sBx93Zdv7OGq7JZw5tb1/uZkHjW47jqKB+PfshKrGVtW4ain7Etx0Cm09kjYB1mTFtpteRcQjwDPA5lnRg8Di4v6z5WWksyAzMytJ2QnoamAPSWvlysYCC4CbWtlR1hFhXdJZEBGxCLgROKBQdSxwe0S0f03FzMz6rOxLcBeQbhC9VNJppLOXccBZ+a7ZkmYAN0XEodnyt4ElwJ2kS2lvBI4lnfXkG21OJt0j9B3gMtKNqHsB7+vfp2VmZo2UegYUEbOBXYFVgMuBrwNnk0Y1yFs1q9PjHuDdwA+BK0lJ7BJgh4iYl9v/LcAHgd2Aa4D/BD7iURDMzMpX9hkQETEN2KVBneGF5Yksf6bT27aXkc5+zFY6wzvUOaLdThYzT9277Ths4Cm7DcjMzAYoJyAzMyuFE5CZmZXCCcjMzErhBGRmZqVwAjIzs1I4AZmZWSmcgMzMrBROQGZmVgonIDMzK4UTkJmZlcIJyMzMSuEEZGZmpXACMjOzUrQ8HYOk7YB3ABsCa5Cmwf47cFs2v4+ZmVlDTSUgSZsDnwEOAoYBy0gzkS4C1gFeBSyTdBNwIfDLiFjWLxGbmdlKoeElOEkXAvcC2wMnAW8F1oiI9SPitRExGNgAeD8wFTgduE/STv0XtpmZvdw1cwa0ENg6Ih6pVyEingKuBq6W9EXgAGDjzoRoZmYro4ZnQBFxRG/Jp0b9ZRHxy4j4ZTP1JW0j6XpJ8yU9LukkSas02Obtki6SNCPb7n5JJ0pao1BvnKSo8Xhfs8/HzMz6R8udEPIkbQvsDAi4KSKmtrj9UOA6YBqwD7AFcCYpMX61l03HZnVPAx4A3gycnP3cv1D3WaCYcO5rJU4zM+u8PicgSZ8BTgGuB9YEzpB0VESc38JuDgcGAftFxHPAJElrA+MknZ6V1XJaRPw7tzxZ0kLg+5I2LZyxLYmIO1qIyczMuqCZTgivqrPqS8CoiDggIvYCPgt8pcXj7wlcU0g0E0lJaed6GxWST48/ZT83aDEGMzMrQTM3ov5d0kE1ykXqjt0j+nD8rYHp+YKIeBSYn61rxbuyeO4vlK8j6SlJiyX9SdJ+fYjTzMw6TBG95w1J7wHOBl4AjoyIu7Pyz5K6ZV9Pug9oV+DYiDi36YNLi4FjIuI7hfLHgIsj4vgm97Mh8Ffgqog4OFf+UdIZ0Z+BwcCngb2A/SPi0jr7Ogw4DGDYsGEjJk6c2OzT6bi5c+cyePDg0o5fT1Xjgv6LbeqsZ9vex7BB8MSC9vax3cZDlluualydMhDfa+0qO64xY8ZMiYiRzdRtmIAAJAn4JCnhTAK+FBH/lPQWXrpUdnNE/LmVQLMEdHREnFMonwVMiIiGl/QkrUbqyPBaYERvozFkz+M2YFBEbN9o3yNHjox77rmnUbV+M3nyZEaPHl3a8eupalzQf7ENP+7Ktvdx1HZLOHNqW/1+mHnq3sstVzWuThmI77V2lR2XpKYTUFNjwUXyA2Ar4AlgqqTjgekR8b/Zo6Xkk5lNGkmhaAhppIVeZQnlYuBNwF6NhgKKlG0vBd7cqKu3mZn1r5YGI42I5yLiGGAH4J3AdEkfbOP40ym09UjahNSrbnrNLZZ3Nqn79j4R0Uz9Hn1przIzsw5qqhecpG9IujNrxB8PLIyIfYBPASdKuim7HNeqq4E9JK2VKxsLLABuahDXl4HPAR+NiFuaOVh2xvQB4C8RsbQP8ZqZWYc0c+H3h8A2pHt+5pMa6CdJ2iYirpO0PWmg0kmSLouIw1o4/gXAkcClkk4DNgfGAWflu2ZLmkG60fXQbPkjwDeBCcAsSTvk9vlgTzftbHDUS0hnU2uSEuYOwL4txGhmZv2gmQS0J3BAREwCkHQr8DRpJIIZ2ZnEeZJ+RkoeTYuI2ZJ2Bc4DLie1+5xdYz+rAvk2m92znwdnj7xDSIkJYAbwP8BGpC7afwT2joirW4nTzMw6r5kENB34mKQppIFJPw3MAx7LV8o6AHy+1QAiYhqwS4M6wwvLB7Ni4qm13aGtxmNmZt3RTAL6OOmM4ilS4/3DpDOihf0Yl5mZreQaJqCIuB8YJWlNYDXPempmZp3Q9N1nETGPdOnNzMysbc10w/5YqzdtSnq9pHf3PSwzM1vZNXMj6lHAg5JO7u1eH0nrSjpI0uWkkak36lSQZma28mmmDWh7SWNJN31+RdJc0oRuTwGLSEPpbAa8jjS0zk+BwyNiVr9FbWZmL3tNtQFl02v/UtIWwG7A24ANSTd3PgHcDNwKTI6Ixf0Uq5mZrURaGgI3Ih4EHuynWMzMbABpaTBSMzOzTnECMjOzUjgBmZlZKZyAzMysFC0lIEn/IclJy8zM2tZqMvktaf6d0yS9sT8CMjOzgaHVBLQFMB44EPibpNslfUrS2p0PzczMVmYtJaCImBkRJ0bEZsB7SRO+nQ38U9JPJI3pjyDNzGzl0+f2nIi4ISI+BrwBmAIcBFwn6WFJX5DU0k2uZmY2sPQ5AUnaWdIE4H5gW+C7pKmyfw18Hbi4yf1sI+l6SfMlPS7ppEajb0t6u6SLJM3Itrtf0omS1qhRd0dJd0pakCXHI1t9rmZm1nktnaVI2pQ0Q+rHgeHAZOAw4NKIWJRVu17S7aRBSRvtbyhwHTAN2IfUxnQmKTF+tZdNx2Z1TwMeAN4MnJz93D+3/9cD1wBXAF8G3gGcJWl+RFzYzHM2M7P+0eplsoeAx0lTdP8oIh6uU+9e4K4m9nc4MAjYLyKeAyZlHRrGSTo9K6vltIj4d255sqSFwPclbRoRj2Tlx2TxfjQilgA3SHodcKKkH0ZENBGjmZn1g1Yvwb0f2DQivtZL8iEi/h4RzXRI2BO4ppBoJpKS0s697P/fNYr/lP3coLD/S7Pkk9//a0mXDc3MrCStJqCRpGkYViBpI0kntLi/rYHp+YKIeBSYn61rxbuAZaQ2KSStCWxS3D9pLqOeY5uZWUlaTUAnks4eanlNtr4VQ4E5NcpnZ+uaImlD4CvAT3JnU+tkP4v7n507tpmZlaTVNiAB9dpNXstLH+6tqLW/3o6zfEVpNeBXwFzgC03uv265pMNIHSsYNmwYkydPbiaMfjF37txSj19PVeOC/ovtqO2WNK7UwLBB7e+n+NyqGlenDMT3WruqGlctDROQpJ5eb5A+tL8nqdg5YA1gO+DaFo8/m5fOVPKGUPvMqBibSN293wTsGBH5BNizfXH/QwvrlxMR40mjPTBy5MgYPXp0ozD6zeTJkynz+PVUNS7ov9gOPu7Ktvdx1HZLOHNqe7fHzTxo9HLLVY2rUwbie61dVY2rlmbedfOBp7PfBTwLPFOo8wJwNXB+i8efTqEtRtImpKm+i203tZxN6r793ogotiXNk/SP4v5zy83s38zM+knDBBQRvybdXIqki4CTeusB16KrgWMkrRURz2dlY4EFwE29bSjpy8DngAMj4pZe9v8BSV+NiKW5/f8D+Fvb0ZuZWZ+1OhbcIR1MPgAXAIuASyXtlrW/jAPOynfNzkY8+GFu+SPAN0mX32ZJ2iH3WD+3/zNIbVM/kTRG0rHAp0lJ1PcAmZmVqNTx2iJitqRdgfOAy0ntMmeTklDeqkB+eJ7ds58HZ4+8Q0g3yhIRMyS9DziLdDb0L+Aoj4JgZla+Zjoh3AUcHBHTJN1Ng95pEfGOVgKIiGnALg3qDC8sH8yKiafetreQhuAxM7MKaeYM6F5Sm0zP7750ZWZmbWumE8Ihud8P7tdozMxswOjzdAxmZmbtaKYNqGG7T16rbUBmZjYwNdsG5HYfMzPrqGbagA7uQhxmZjbAuA3IzMxKUfp9QGZmNjD5PiAzMyuF7wMyM7NStDwWXDYB3MGk4W02Av4J3An8OCJe6Gh0Zma20mqpE4KkNwIPAN8FtgWWZj+/C8yQtE3HIzQzs5VSq2dA40kT0r07Ih7tKZT0OuBK0vQK7+lceGZmtrJqNQGNBD6cTz4AEfGopBOAn3csMhtwhndoeul2p6meeerebcdhZo21eh/QTGCNOuvWAB6ts87MzGw5rSag44BvSHpnvlDSDsBJwJc6FZiZma3c+jIY6drAbZKeBJ4ENsgeTwPHA5f1Q5xmZraS6ctgpPf2UyxmZjaAlD4YadZ1+1xgFDAHuBD4ekQs7WWb1YBTgB1IHSPWiAjVqDcB+HiNXbwxIqa3H72ZmfVVyzeidpKkocB1wDRgH2AL4ExS29RXe9n0VcAngbuA24Bdeqk7HTikUDazbxGbmVmnlJqAgMOBQcB+EfEcMEnS2sA4SadnZSuIiDmSXh0RIekIek9A8yLijs6HbmZm7Wh5OgZJYyVdJ+lRSU8WHy3ubk/gmkKimUhKSjv3tmFEeFBUM7OXsVaH4vkI8GNgBvBa4HfAFdl+ngPOa/H4W5Mukb0ou8l1frauE7aR9JykRZJukdRrYjMzs+5QKycSkv4E/AY4FVgMjIyIP0paC5gE/CYivt3C/hYDx0TEdwrljwEXR8TxTezjCODcOp0QPg+8QGpjWh84ChgB7BQRd9XZ32HAYQDDhg0bMXHixGafTsfNnTuXwYMHl3b8evorrqmznm17H8MGwRMLGtfrzXYbD1mhrKqxVTWuTqnq/wBUN7ay4xozZsyUiBjZTN1W24C2BG6NiKWSlpLuCSIinpd0GnA20HQCytTKgKpT3tqOI85ZbqfSlaRkdDywb51txpPGvGPkyJExevTodsPos8mTJ1Pm8evpr7jaHUIH0lA8Z05tr2lz5kGjVyiramxVjatTqvo/ANWNrapx1dJqG9CzwOrZ77OAN+bWCVi3xf3NBtapUT6E1CW7oyJiAXAV8LZO79vMzFrT6teee4A3A9eQ2n9OkLSEdJnrBNK8QK2YTqGtR9ImwJoU2oY6zB0YzMxK1moC+hawafb7Cdnv5wOrAHeTtZ204GrgGElrRcTzWdlY0hTgN7W4r4YkDSL1vJvS6X2bmVlrWkpA2f00d2S/zwH2kbQ6sHq9e3YauAA4Erg0a0PaHBgHnJXfn6QZwE0RcWiubE/SmdL22fIHs1V3R8QjkoaQeuj9lNRrbz3gC8DGwIF9iNXMzDqoY1NyS2p5Su6ImC1pV1L37ctJ7T5nk5JQMc5VCmXf46WzMYBfZz8PASYAi4B/k0ZU2ABYCNwO7BwR97QSp5mZdV5LCSibkvv3wGtIl7GeJE3J/V/A1yS9LyKmtbLPrH5vIxkQEcObKSusXwjs10osZmbWPZ6S28zMStFqN+yRwAm1puQmdUp4e6cCMzOzlZun5DYzs1K0egnuOOBMSQ9HxIv3/OSm5D6mk8GZ2cvX8A6N0tDuaA8zT9277Tisf3hKbjMzK4Wn5DYzs1KUPiW3mZkNTH0aAlfSa4BRwKtJl97uiIjHOxmYmZmt3Fq9EXUV4FzgUyw/MsFSSeOBz0XEsg7GZ2ZmK6lWu2F/HfgEqbPBcNLU2cOz5U+w4hA6ZmZmNbV6Ce6/gK8WZj19FDhDUpAGFj2hU8GZmdnKq9UzoA2Av9ZZ99dsvZmZWUOtJqC/Ax+qs+5DwP3thWNmZgNFq5fgvgFMzAYf/Q3wBOms5wBgDPWTk5mZ2XJanZDuV5LmkDojnAO8ElhMmprhfRExqfMhmpnZyqjpBCTplaRJ6P4WEaMkvYI0y+hT7nptZmataqUNaClwA/BGgIhYFhFPOvmYmVlfNJ2AskTzADCskwFI2kbS9ZLmS3pc0knZDa+9bbOapDMk/UHSgqwLeL26+0iaKmmhpGmSxnYyfjMz65tWe8F9BThB0nadOLikocB1pMFO9yFN6XAUqY2pN68CPgnMB27rZf87AZcANwJ7kmZt/YWk3dsO3szM2tJqL7ivAusCf5Y0i9QLbrmzj4h4Rwv7O5w0msJ+EfEcMEnS2sA4SadnZSuIiDmSXh0RIekIYJc6+/8acHNEHJkt3yjpTaSbZa9tIU4zM+uwVhPQvcDfOnj8PYFrColmInAasDNweb0NI6LuZTcASauTuoYfWVg1EbhI0pCIeLZPUZuZWdta7YZ9cIePvzWpY0P+GI9Kmp+tq5uAmrAFqZv49EL5faRLj28A7m5j/2Zm1oamEpCkQcBepIFH/wlcHxFPdOD4Q4E5NcpnZ+va3Tc19j+7sN7MzEqgBleykLQ5qaPA8Fzxc8CBEdFWO4qkxcDREXFOoXwWMCEivtLEPo4Azo0IFcp3BG4Bto+Iv+TKtyQNKbR7rRtnJR0GHAYwbNiwERMnTmz9iXXI3LlzGTx4cGnHr6e/4po6q/0rosMGwRML2tvHdhsPWaGsqrFVNS6odmydMND+P5s1ZsyYKRExspm6zZwBnQ4sA95NGvFgM+B84PvZ7+2YDaxTo3wItc+MWt03Nfbfs1xz/xExHhgPMHLkyBg9enSbYfTd5MmTKfP49fRXXAcfd2Xb+zhquyWcObVP8yy+aOZBo1coq2psVY0Lqh1bJwy0/8/+0Ew37FGkKRhujYiFEXEf8GngdZI2avP400ltPS+StAmwJiu23bTqQdIwQVsXyrcmJdS/t7l/MzNrQzMJaCPgoULZg4CADds8/tXAHpLWypWNBRYAN7Wz44hYRLr/54DCqrHA7e4BZ2ZWrmbPbXtvKOq7C0jdpC+VdBqwOWlW1bPyXbMlzQBuiohDc2V7ks6Uts+WP5itujsiHsl+PxmYLOk7wGWkjhR7Ae/rp+djZmZNajYBXSNpSY3y64vlEdH0pHQRMVvSrsB5pC7Xc4CzWXFq71WB4vA83wM2zS3/Ovt5CDAh2/8tWWL6BvAZ4GHgI+12njAzs/Y1k4AaDYvTloiYRv2RDHrqDG+mrM62l5HOfszMrEIaJqCI6NcEZGZmA1Org5GamZl1hBOQmZmVwgnIzMxK4QRkZmalcAIyM7NSOAGZmVkpnIDMzKwUTkBmZlYKJyAzMyuFE5CZmZXCCcjMzErhBGRmZqVwAjIzs1I4AZmZWSmcgMzMrBROQGZmVgonIDMzK0XpCUjSNpKulzRf0uOSTpK0ShPbDZF0kaTZkp6V9DNJ6xbqTJAUNR5b998zMjOzZjSckrs/SRoKXAdMA/YBtgDOJCXGrzbY/JfAVsAngWXAacBlwLsL9aYDhxTKZrYTt5mZta/UBAQcDgwC9ouI54BJktYGxkk6PStbgaRRwB7AzhFxc1Y2C7hT0m4RcV2u+ryIuKN/n4aZmbWq7EtwewLXFBLNRFJS2rnBdk/0JB+AiLgLeDhbZ2ZmFVd2AtqadInsRRHxKDA/W9f0dpn7amy3jaTnJC2SdIuk3hKbmZl1iSKivINLi4FjIuI7hfLHgIsj4vg6200iXVrbt1D+U2DziHhXtvx54AVSG9P6wFHACGCn7Iyp1r4PAw4DGDZs2IiJEye28QzbM3fuXAYPHlza8evpr7imznq27X0MGwRPLGhvH9ttPGSFsqrGVtW4oNqxdcJA+/9s1pgxY6ZExMhm6pbdBgRQKwOqTnlL20XEOcutlK4kJaPjgX2pISLGA+MBRo4cGaNHj24QRv+ZPHkyZR6/nv6K6+Djrmx7H0dtt4Qzp7b3tp550OgVyqoaW1XjgmrH1gkD7f+zP5R9CW42sE6N8iHAnD5st05v20XEAuAq4G0txGhmZv2g7AQ0nUKbjaRNgDWp3cZTd7tMvbahovKuO5qZGVB+Aroa2EPSWrmyscAC4KYG220oaaeeAkkjgc2zdTVJGkTqJTelnaDNzKx9ZSegC4BFwKWSdss6AIwDzsp3zZY0Q9IPe5Yj4nbgGuBiSftJ2hf4GXBLzz1A2UgJf5D0aUm7ShoL3AhsDHyzW0/QzMxqK7UTQkTMlrQrcB5wOan95mxSEspbFSgOz/OhrO6PSIn0CuDI3PpFwL9JIypsACwEbifdvHpPR5+ImZm1rPRecBExDdilQZ3hNcrmkIbYKQ6z07N+IbBfB0I0s5XM8A710Gunp9/MU/duO4aXu7IvwZmZ2QDlBGRmZqVwAjIzs1I4AZmZWSmcgMzMrBROQGZmVgonIDMzK4UTkJmZlcIJyMzMSlH6SAjWXVW4Axx8F7iZ+QzIzMxK4gRkZmalcAIyM7NSOAGZmVkpnIDMzKwUTkBmZlYKJyAzMytF6QlI0jaSrpc0X9Ljkk6SVJx+u9Z2QyRdJGm2pGcl/UzSujXq7SNpqqSFkqZJGts/z8TMzFpR6o2okoYC1wHTgH2ALYAzSYnxqw02/yWwFfBJYBlwGnAZ8O7c/ncCLgHOB44E9gJ+IWl2RFzb0SdjZtamgXajeNkjIRwODAL2i4jngEmS1gbGSTo9K1uBpFHAHsDOEXFzVjYLuFPSbhFxXVb1a8DNEXFktnyjpDcBJwBOQGZmJSr7EtyewDWFRDORlJR2brDdEz3JByAi7gIeztYhaXVgDPCrwrYTgVGShrQfvpmZ9VXZCWhrYHq+ICIeBeZn65reLnNfbrstgFfWqHcf6Xm/oQ/xmplZh5SdgIYCc2qUz87WtbNdz89ivdmF9WZmVgJFRHkHlxYDR0fEOYXyWcCEiPhKne0mAXMj4gOF8p8BwyNiR0k7ArcA20fEX3J1tgT+DuweEZNq7Psw4LBscSvg/j4/wfatBzxV4vHrqWpc4Nj6oqpxgWPri7Lj2jQi1m+mYtmdEGYD69QoH0LtM5z8drWe4Dq57Wbnyop1qLf/iBgPjO/l2F0j6Z6IGFl2HEVVjQscW19UNS5wbH1R1bhqKfsS3HQKbT2SNgHWpHYbT93tMvm2oQeBxTXqbU3qtv33PsRrZmYdUnYCuhrYQ9JaubKxwALgpgbbbZjd5wOApJHA5tk6ImIRcCNwQGHbscDtEfFs++GbmVlflZ2ALgAWAZdK2i1rfxkHnJXvmi1phqQf9ixHxO3ANcDFkvaTtC/wM+CW3D1AACcDoyV9R9JoSaeTbkY9qd+fWWdU4lI65KO6AAAYxUlEQVRgDVWNCxxbX1Q1LnBsfVHVuFZQaicESEPxAOcBo0jtMhcC4yJiaa7OTGByRBycK1sHOBv4ACmRXgEcGRHLNb5lyekbwJak+4TGRcTEfnxKZmbWhNITkJmZDUxlX4IzM7MBygnIzMxK4QRkZmalcAIyM7NSlD0SggGSROrNtzfwRuDVwFLgCeAO0rBEpdw4m90YvBcg4NcR8bSk1wJHkwZ8nQmMj4ipXYzpS8BV3TxmsyQNAlaNiOdzZesDRwDbkG6C/jNwfjfvRZP0VmBQRNyWK3sf8OVcXH8h9RK9rfZeuif7n3g/8DYggHtIf/NK9JrKpo15CtglIm4pMYZdgNWAKyNiXvZe+yzpnsiHSP+bj5cRXzPcC65k2RvmKmAEKeEsAjYm/dNdTXojbQWcHBEndzm2d5DmTRoMLAGeIc3DdBUpQd4LbAtsCOwWEX/oUlzLSK/PdODnwC8jYkY3jt2IpKuAByLi89nyKNLfcRkwhZTIRwAvkD687u1SXHcAl0fEKdnyJ0i3PNwI3JDFtStpQsf9I+K33Ygri+U24NCIuC9bHkp6340A5mbVBpO+jO2RT+79HNd/97J6EHAGcA7wAEBEnN+NuAAkvR64HtgkK3oY2B2YRBpu7EHS58YCYEREPNat2FoSEX6U+AB+QXoDb5crew3we+CSbHln0j/iJ7oc2yTSB9Q6pKktzgMeA34LvDKrszrpA/bGLsa1DPgW6d6vRaRkeDfwBWDjkv+eTwH75JbvIH1QrJUrG0Ia6eOaLsb1HGkA3p7lGcC5NepdAPyly6/ZMuAdueUfkr7svC9X9j7S+I5ndzmupdnPWo/8uqVdfs1+RTpjfT3pislPss+R23rea6RBSf8CfL+bsbX0PMoOYKA/SDff7l+jfHj2Bt8oWz6+hA+Gp4E9c8sbZP9suxfq7Q081cW4XvzAIk2rcVj2Ib8ke0zOytYt4e85H3hPbvmF4uuVe83mdfl9lk9Ai0kzChfr7QYs7PJrVkxA/wb+p0a9o4FHuhjXZcA/gUPIrhbl1q2Txf2ebsVTOP7jwIG55U2zePYr1DsE+HsZMTbzcCeE8omUaIqWZut6Zm69k+5PohfZI79MoazWctdExOyIGB8RuwKvBY4iXRO/AHhc0pVdDulvpJl4ezxB+oZatC4pWXXLH4CDcsv3ArVGTH47MKsrEdW3DqnNp2gK6XJvV0TEvsDHgWOAu7MpXl5c3a046hgK/Cu33PM3e6RQ7yHS/0UluRNC+a4DviHprxHxELx4Dfx/SW+wns4Hg4FuD6A6BTha0q3APNJZ2CzgM5JuiIilklYF/pv0wVuqiPgX6Zr8OZI2BT5MGny2m04FfibpH8DFwCnAGZKeJl12EylBfYsVp4vvT8cDt0p6BXAuqfPBjyW9mnTGCKkN6H+A47oYV4/9swGFof50K+uRLiV2TURcK+nNpNflSkm/J52JdaUdqhdPks56eiwFvk/6wpO3AeXHWl/Zp2AD/UH6dvI30iWRGcA0UsPhHJa//HU6qbG9m7GNJH0YLM5iehp4C+la80PA5aTGz0XAmC7Gtdwlm6o9gE+SPiifBe7Kfl+aPZZkP/8PeFWX49oeuJ0a7RfZz6eBz5fwetVqX/lRjXrfB/5Q4t91Q9KXirnAmdnrVtYluMtqvUY16p0LXFfWa9bo4V5wFSBpFeBA0of7GqRE9POIeKbUwICsy/V/kM6WL4mIf0raEDiW1MvmEeDCiPhjF2M6EfhBVLl7qbQu6ezrHaQPrleQGtbvA66IiCklxvZG4J014rotIhaXFVcjkj4FPBgRN5QcxyjSQMhbAXtHCd3WJQ0jfYF5uEG9L5Lajq/vTmStcQIyM7NSuA2oQiS9iTRj61BSI+ccYHp06V6RVklaJXLTZpRN0hqkm2OXATPK/jaftY9tTu7G4oh4tMyYXm6yG1KJEr8pZzcXKyLm58q2J7sRu8yz2Ze9sq8B+hEAnyBdyqp1z8FS0mgDh5QU236k681XAe/PysZmMS3N4v5Ul2P6KLl7okhfpE4l9SrraWt5HjiupNdsBPA7UrvZ0sJjFmlCxK62/2Rx/Qepu/pU4JfUaL8gXZrr9j0tu5O7Tyor2xf4I6nNbDGpV9zeXY5rCKmtbnEWxw+AVYAfF/4/bwXWK+O91sRz2L/bf89WHu6GXTJJnyM1rl4BjCb1Wnll9tiAdBPqFcAFkj7b5dgOBH5D6n20GPhldh3+J6QPsiNJN75dIGmPLoZ2POkG2B6nZbF8C3gP6TU7EzhR0vFdjAtJu5Nek9cA3yHNynst6YNqHHAW6UPhtqy3Y7fiei/pBuI1SL3xXg/cKOnMnrOMEl1NGoIKAEkfAC4FFpJ65H2ZdD/Vb7PXt1tOJo0M8UXSl8R3kXou7kK6MXYYKakPz+pai9wGVDJJDwEXRMTpDeodCxweEZt3JzKQdDcwJSIOz5YPAn4EnBcRR+XqXQRsEhG7dSmu+aQegjdly08Cp0TEOYV6RwOfi4hNa+ymv2KbAvwtIj5eKP8c6R6lzUn3Kd0G3BERvQ330sm4biENEXRIruwTpO7+k4APR8RCSe8kdUZYpRtxZXEsA3aIiLuy5T8CsyLi/YV6VwFrRsTOXYrrYeCbEfGDbPmtpFsTDomIH+fqfQo4PiI260Zc2TF/1GTVTYHR3fx7tsJnQOXbkNRVt5G76OJNeJmtSGdAPa4gnZkVb+68lDSgZbc8Szor6zGENORI0V9IZ5HdtA3w0xrlPwVeB2wVEQtJZ0cf6GJc2xbjiogfkc4WdwBuyO4JqoJtSVcFisaTBiftlg146T48yMZ8I42zljeD2vct9aePky5dbtfg0bUvX33hBFS+vwKfym4QrCm7RPKprG43Bemad4+egSHnFOrNJd293i2/I90gu1q2fB3pptOiD5Pu+O+mJ0nd6YveQno9e24mfoSXRrnohoXAmsXCSA3oO5I+QG8DuvYtvhhK7vdneem9ljeP7n5mPUxK0D3eTWr3eVeh3o5AtzuXPADcEBFv7+1BOTcVN8294Mp3FGng0WmSLiWN8DyH9A+5DqlX3AdIN6y+r8uxPUL6Rn8NQKSRD0aR7hnJ25zlhwXpb18mDS3zN0kXkm6IPU3Strx0V/8uwFtJQ/p303jgZElrktrJXiANb/MV0oCtPfcubU53P7T+CuxJSt7LiYiHsmFmrgImdDGmvGskLcl+H0K6aXZyoc7WpLHZuuUC0qga25GS4oGk994JkgaTzrDfRhoEt9ttQHewYiKsJUijb1SSE1DJIuLWrEvnsaSxujYpVPkHqZH2jIgonvr3t0spjCMVEXfWqPcRoGtzokTEM5J2IH2of5GXLrONyh4vkNo13h0Rd3crriy2U7I2jeOAE3uKSaOe/0+u6mLgm10M7RLgeEmvjho3OEfEk5J2JvX66kpbXs7Xa5Q9WaNsf9Lo7F0REedlVyY+TDozPDYiLpD0GKntrGc8vwuAb3crrsy5pF6CjdzE8mMTVoo7IVSMpFfx0uWsOZG796CqJL2OFGtXx+nKHX84y9/V/2CUfw/QK0n3iawBPFTWa2P9I7ssvl5E/LvsWF7OnIDMzKwUvgRXEdnU1xsA90fECg2wktYD9oqIi7seXA3ZNfA/Agd1+zJXVae9zsVSuWnMm5WNE3dARJzU5eO+LKaXrjFV+BRSvF3/Jp+NHr4/6X02ISKmS3oL6ZJmz/vsuxHx+27H1iyfAZVM0uqk7rH7ZUXLSCPufjH/4VnS/Rl79bJ6TdLd9MeRTcUQEVd1Ka5KTnudxVLJacybJWl/4Fddfp9Vcnrpqk4VnsWyB6nzzTOk3oHrA/uQ2m2nkb6AjSB1gNk/Ii7rVmwtKXsohoH+AE4g9Xr7FGn6g8+T5vR4ANgyV6+MIVIqOSUxFZ32OjtuVacxf12Tj8NLeJ9VcnppKjpVeHbcW4FfA6tky8dncfywUO8npBueuxZbS8+j7AAG+oPU7fqIQtmGwM2kqYlHZWVlJKApvDQl8aaFx5uzf9ADe8q6GFclp73OjlnlacyL49LVenT1y0QWWyWnl66RgCoxVXh2zGdJZ9A9y0OzeHcp1Nud1EGoa7G18nAbUPk2oXCDaUT8S9KupG8v12VD4HTz/oceI0lnZqeR7rs5OrL5RyT13ET5r4goTgPc33qmvb45W67KtNdQ3WnMnwduAC5sUG8n0i0B3fRymV66ElOFZxaw/I3FPb8PKtR7Fekm5EpyAirf48CWvPRhCkCkbsQfkvQd0ql21zsfRPoKNV7Sr4BvAH+VdF72e5mqOu01VHca87uAIRFRHEZpOdmUFt1W5emlKzlVOOkS3AmSHsiO/W1S28+XJN0cEc9nXxKPJb0nq6nsU7CB/iAN7jm5QZ0vU8KlkRpxvJl0d/rjpLaqMqckruq011WdxvxrwD+aqPceutg2lR2zktNLU+GpwkntZTNz7/kHSW14Pf8LU0nJejawfTdja+XhXnAly75djQVOjYine6n3EeC9kRvNuCySPgScTrocMjoibm6wSX/FUclpr6s4jXmVvdynly5rqvDsVogdSZ1dro+IBdmN7J/kpffZz6NLvQb7wgnI+iS7jLQmMDcqNCuqmb18OAGZmVkpPB3Dy4SkH0j6Ydlx1FLV2KoaF1Q3NknXSarUJa4eVY2tqnFBtWMD94J7ORlDdb8wVDW2qsYF1Y1NVDMuqG5sVY0Lqh2bL8GZmVk5KpsZbXmS1simPaicqsZW1bigurFJemUV44LqxlbVuKDasYET0MvJ3qT7R6qoqrFVNS4oITZJn5X0oKTnJd0p6WM1qr2t23FVObaqxlX12JrlBGQ2AGT3bp1LGrj166SbiSdI+k02vYVje5nEVfXYWuE2oJJJavbmtfWBbaK7w+RXMraqxgXVjU3SPcANEXFsrmxX4GekO+r3jjRvURnTflQytqrGVfXYWuEzoPK9BxhGGrKlt0e3x8CqcmxVjavKsW1FmpPoRdmIAjuQpq64XdIWXY6pR1Vjq2pcUO3YmuZu2OX7G2kW1LG9VZL0QdIEcN1U1diqGhdUN7ZnSYNmLiciZkp6F3Alaf6dk7sYU4+qxlbVuKDasTXNZ0Dlu5P0raWRIPXp76aqxlbVuKC6sU0B9q0ZSMRsYFfSVAP/28WYelQ1tqrGBdWOrWlOQOU7HfhcE/WuAjbr51iKqhpbVeOC6sb2U2BzSbXmTSIiFgD/SZov6NEuxgXVja2qcUG1Y2uaOyGYmVkpfAZkZmalcAIyM7NSOAGZmVkpnIDMzKwUTkDWK0kHS5qSjTc1W9KfJJ3VT8c6UNLBTdQbJylyj8clXdLsjXeSJmR3kpeu2eec1e153g/UWT8jWz+uv2Jocb/Lvc6dPo6kV0g6IntPLpD0nKR7Jf2vpD51cVfyF0kfr7N+gqThddZ9VxWc46nKnICsLklfJnXjvAbYD/gv4Lek7p394UDg4CbrPguMyh5HA9sD10tas4ltT27hOP2tlecMsBDYTNLIfKGktwObZuv7O4ZmFV/nTh9nIvAN4FLSe/LjpO7t74q+d+89EBgK/LwP254BHCTp9X089oDjkRCsN0cA34+I43Nll0v6elkB5SyJiDuy3++Q9CjwB2Av4NfFypJWAVaJiBci4sEuxtlp84A/Ah8i3WjY40PADcCIMoLq0a3XWdKewAHAXhFxdW7V//X17CdzJPCTiFicO9aqpGT6MeA1wIclPQh8PSJeHM0iG4XgFuAzwFFtxDBg+AzIerMO8K9iYf7bZc9lFkn7SpouaaGkWyRtU9wuuwQzVdIiSf+QdEr2z42kCcD+wM65S2vjWoh1SvZzeI247iWdGbwzv64Q23sk3ShprqRnJU2W9Nbc+p0k3SRpvqSnlabUXqu3gCSNkvS77BLhPEl/lnRQ/rXr43OeCBzY80Gb/TwwK+9YDNlr8JvC/kZndbbNv5aNXud6x5G0t6RlkjYrHGezrLze2fbO2c8VBn/t69lPdubyLuA3hVWfB44ljSpwFfAJ4EfAujV2cwnpLMifrU3wGZD15o/A57Kziysi4uk69TYFzgK+BiwgDQ9/jaQtI2IhgKTdSWOfXQwcA7yZ9K1yXeDw7PfXkZLef2f7fayFWIdnP/9VKDsdOAl4gjrzokgaDUwCbiRdxpkH7AhsDPxJ0o7A9cBlwAezmE8lXar5YC8xbQrcClxA+mDeEbhI0rKI+AV9f86XAt8DdiKd9b2bNLr2/5EuA3UjhrzhNH6d6x3nn6SpBD4OjMvVPxj4N4UBN3PmZT/PkHRmRDzSYsy17Jrt9y+F8p1JI0+fnn2xujUiZtbZx22kwWi3q7EfK4oIP/yo+SAliYdI45YtA+4lfcisnaszIVv/rlzZpsAS4PBc2R3AjYX9HwssBV6bLf8GmNxEXOOAp0hfoFYF3kBKHs8BGxXi2r7G9hOAe3LLt5MuZ6nO8f5QI/Zdsv1v2+RrqSzW75M+zHrKm3rO+eed/f5b4LvZ7+cDl2W/PwWM60QMwGTgN4Wy0fnn3eLrXO843yAlLeXinAl8u5fXYkPgr9mxgzQI7PHA4Dbe7+OBu2uUfx/4R3bMCcDwXvaxavbe/1Rf4xhID58mWl0R8VfgjaQG3vNJHwxfA+6RNDhX9cmIuC233SOkS2LvgBfbBd7Gim0zvyRdBh7Vh/DWBRZnj/uBzYGxEfHPXJ1ZEfHn3naSdVp4J/DjyD5BCutflcX3K0mr9jyAW7Jj121zkTRUqUfWI7lYDyMlzHZNBD4oaXXSWdgKl9+6EEOPhq9zAz8ifWkZnS2PyZYvqrdBRPwLeCuwB+lscB3gFOA2SavBiz04/5w9FmWXiP+s1KvzlTV2uyEpgRedQjozepj0v3B0dlZcK64lwJxsX9aAE5D1KiIWRcTlEXFERGwDfBLYEjg0V+3JGps+CWyU/b4e8ErS5Zm8nuWaAyo28CzwdmAk8FrSt9KrC3WKx6tlKCmx/rOX9auQEvDi3GMR6Tlt0su+JwBjSZfFds/i/RGwRhNxNfI7YDDpw3FN4PISYujRzOtcV0Q8RDrbOiQrOgS4KyLubbDd0oi4NiL+m3R57yLSpa9R2foJEbE96cvPEmDHiNg+IkZErpNBzhqkv2vxOI9m+/0A6YrATsAtqn87wiI6+/qutNwGZC2JiB9KOh3YOle8QY2qG5Au2UH6Vrm4Rr1h2c9n+hDKkohodC9PM43Rs0mXFzeqs35Otp9x1G6PeLzWRpLWAPYGjoiIC3LlHfnSFxHzJF0BfAH4dUTMK9bpQAwLgdUKZbW+LHRiROMLgR8odf3fjxZ7kUXEMknXkpJX8cN/S2B21G/D7PEMdc5csoT1e6WpsMeRpkI4W9J3sgSVtw59e08POD4DsrokrZBYJK1PmnEx/613A6VJsHrqvI70rfMuSN9USZfkDijs7kDSh//t2fILdPmbY/bBfSfwXz29ymqsvwPYKiLuqfGomYCA1UlnTi9+o856zRV7dbXznL9HOvO5oM76dmN4jOW/aAC8t0+R9n4cSB0rXiBdSnwFdS4pAkgaVmfVfwLzSX/PvLfQXIeA+6kxRUat9wVwd/bz1YW66wOvAv7exPEGPJ8BWW+mSvotcC3pktqmpJs+5wM/ztV7CviJpJ5ecCdl9Sfk6pxI6hl3EenDZTtSz6gfRERPr6vpwD6S9iV9+D3eywd8Jx0HXAdcLWk86Xr/KFID+hWkzhLXS1pGakh/nnTJZ2/gKxGxwodNRDwr6W7gBEnPkRLtcaRLh2vnqvb5OUfEZNKlq3rr243h/4BDJZ1NmmFzDKnNpa/qPteIWCjpZ8BngV9ExJxe9vMrSc8DvyJ1VtgAOAjYh9T4X9z2LaQOC43cSnqt1o+If+fKfy7pT8DNpMudI0hnnrOA+wr7GEk6I7wNa6zsXhB+VPdB+jC4lnSZaSHpn/3nwNa5OhNIPcj2I33rW0T6R16hdxipLWIq6ZvuY6T2i1Vz69cjfeg9Q3bZq05c48h6g/US+wRyPbAarSN1tb2ZlFznkHrVbZ9b/07g96SedvOAaaSu50N6ieH1pPtU5pEmBTu2GHuzz7mF571cL7h2YwC+TOoB9jxpErT/ZMVecE29zo2eK7BbVr5bg+f4iexv8Vj2XnqGlCBH16l/OfChJt7vqwFPAx8rlH8gO96/SEn8OVLif2uNfZxDocekH/UfnpDO2pLdYLhtRIxsVNesN1nb4lhgs4hY1sH9PgrsERHFs5Vadc8BXh8Re9dZP4GUOGfWWLcK8AhwXET8tK2gBwhfgjOzUknaCtiGNITN1zucfIaSbtJttk3mDOB+SW+IGpdWGziAdAm6bvuVLc+dEMysbN8nXdq9ijTcTcdExOyIGBSpI0wz9R8j3WJQr1fkZaRLtLUIODTSvUDWBF+CMzOzUvgMyMzMSuEEZGZmpXACMjOzUjgBmZlZKZyAzMysFE5AZmZWCicgMzMrhROQmZmV4v8Bzzme+z7FdgYAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot probability distribution\n",
    "x = uncertainty_model.values\n",
    "y = uncertainty_model.probabilities\n",
    "plt.bar(x, y, width=0.2)\n",
    "plt.xticks(x, size=15, rotation=90)\n",
    "plt.yticks(size=15)\n",
    "plt.grid()\n",
    "plt.xlabel('Spot Price at Maturity $S_T$ (\\$)', size=15)\n",
    "plt.ylabel('Probability ($\\%$)', size=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Payoff Function\n",
    "\n",
    "The payoff function equals zero as long as the spot price at maturity $S_T$ is less than the strike price $K_1$, then increases linearly, and is bounded by $K_2$.\n",
    "The implementation uses two comparators, that flip an ancilla qubit each from $\\big|0\\rangle$ to $\\big|1\\rangle$ if $S_T \\geq K_1$ and $S_T \\leq K_2, and these ancillas are used to control the linear part of the payoff function.\n",
    "\n",
    "The linear part itself is then approximated as follows.\n",
    "We exploit the fact that $\\sin^2(y + \\pi/4) \\approx y + 1/2$ for small $|y|$.\n",
    "Thus, for a given approximation scaling factor $c_{approx} \\in [0, 1]$ and $x \\in [0, 1]$ we consider\n",
    "$$ \\sin^2( \\pi/2 * c_{approx} * ( x - 1/2 ) + \\pi/4) \\approx \\pi/2 * c_{approx} * ( x - 1/2 ) + 1/2 $$ for small $c_{approx}$.\n",
    "\n",
    "We can easily construct an operator that acts as \n",
    "$$\\big|x\\rangle \\big|0\\rangle \\mapsto \\big|x\\rangle \\left( \\cos(a*x+b) \\big|0\\rangle + \\sin(a*x+b) \\big|1\\rangle \\right),$$\n",
    "using controlled Y-rotations.\n",
    "\n",
    "Eventually, we are interested in the probability of measuring $\\big|1\\rangle$ in the last qubit, which corresponds to\n",
    "$\\sin^2(a*x+b)$.\n",
    "Together with the approximation above, this allows to approximate the values of interest.\n",
    "The smaller we choose $c_{approx}$, the better the approximation.\n",
    "However, since we are then estimating a property scaled by $c_{approx}$, the number of evaluation qubits $m$ needs to be adjusted accordingly.\n",
    "\n",
    "For more details on the approximation, we refer to:\n",
    "<a href=\"https://arxiv.org/abs/1806.06893\">Quantum Risk Analysis. Woerner, Egger. 2018.</a>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set the strike price (should be within the low and the high value of the uncertainty)\n",
    "strike_price_1 = 1.438\n",
    "strike_price_2 = 1.896\n",
    "strike_price_3 = 2.126\n",
    "strike_price_4 = strike_price_3 + strike_price_2 - strike_price_1\n",
    "\n",
    "# set the approximation scaling for the payoff function\n",
    "c_approx = 0.25\n",
    "\n",
    "# setup piecewise linear objective fcuntion\n",
    "breakpoints = [uncertainty_model.low, strike_price_1, strike_price_2, strike_price_3, strike_price_4]\n",
    "slopes = [0, 1, 0, -1, 0]\n",
    "offsets = [0, 0, strike_price_2 - strike_price_1, strike_price_2 - strike_price_1, 0]\n",
    "f_min = 0\n",
    "f_max = strike_price_2 - strike_price_1\n",
    "iron_condor_objective = PwlObjective(\n",
    "    uncertainty_model.num_target_qubits, \n",
    "    uncertainty_model.low, \n",
    "    uncertainty_model.high,\n",
    "    breakpoints,\n",
    "    slopes,\n",
    "    offsets,\n",
    "    f_min,\n",
    "    f_max,\n",
    "    c_approx\n",
    ")\n",
    "\n",
    "# construct circuit factory for payoff function\n",
    "iron_condor = UnivariateProblem(\n",
    "    uncertainty_model,\n",
    "    iron_condor_objective\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAE+CAYAAACnXJZvAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XeYVPX1x/H3AVQ6YsPOWmJBUxRLVNSAKBI7FlAUFRRRQNjENM0vtpBiCktRAUUQRcHeoqIoYC8QYwyKjaZUaSIu0vb8/jh3ZRh2d2a2zHfuzHk9zzy7c+fe3c8Ow5y599tEVXHOOeeqq17oAM455+LNC4lzzrka8ULinHOuRryQOOecqxEvJM4552rEC4lzzrka8ULicpaI3CQimnBbKCKPish+ATMdJyL/FpHvRESjbU1EZIKILI9yXlbJsWOT/p7y291Z/SM25+ktImdXsH2uiPw9RCYXTw1CB3Auha+BU6Pv9wVuBV4SkUNU9dsAeUYCS4FOwLpo29XAGUAPYAHweRXHzwIuT9q2tJYzpqs38D/giaTt5wDLsx/HxZUXEpfrNqrqW9H3b4nIfOBV4OfAwwHyHASMUtVpSds+VtVH0zj+24S/Jyep6nuhM7h48UtbLm5mRF+LAETkGBF5Krrs9a2I/EdEupfvLCI7RJehLk38IWLmiMg/E7Z1EJG3o/2XiMgdItI0euxn0aWs+sCQ6JLUWBGZC/QCDiu/VFXdP6z8d4jIoUnbp4rIIwn3x4rIdBE5WUT+G/3dr4nIIUnH1ReR34nIJyKyTkS+FJGx5T8TaAtcmnCJ7bLosa0ubYnIBSLyQfRzvhCRQSLSIOHxy6Kf8UMReTHKNEtEulT3+XDx4YXExU1R9HVx9LU18DpwBXZ56VFgjIhcCKCqK4DH2fpy0s+inzUGQETaAM8Dy4BzgRuBi4DyN/B/A8dE3/8j+v5W7DLQs9glq2MS9qmUiDRIvKX1V29tb+BvwCDgQmAX4CERkYR9RgI3Aw8BpwO/BJpEj10TZX42Ife/Ksl7CjARew7OAoYB1wHDK9j9AeAp7Hn5FJggIntW8290MeGXtlzOS3iz3Re4A/gGmAygqhMS9hPgFWBP4Ergweih0cALIrKvqs6Otl0OzFDVD6L7fwDmAWeq6qbo560AJorIMar6JnZpDWBu4uUpEfkKaJXmJau2wIakv+8HqvpZGscm2gE4TlU/jX5GPaxgHgjMEpGDsDOlAao6NOG4iQCq+qGIfAt8lUbuW4Cpqlp+Vvd89Dz8WUT+qKpfJuw7WFXviTLNAJZgRWxEhn+fixE/I3G5bkfsjXcD8DFWTLqq6iIAEWkpIkNFZF7Cfr2BAxJ+xktYkbg0OqYZ0IXobCRyFPB4eRGJPApsBNrV4t/zEXBk0u2LavycueVFJPJh9LX803/76OvYavzs74lIfeBwtm6Pmoi9fySfgb1Q/o2qLsc6EvgZSZ7zMxKX674GOgKKXc5aqFtOWT0W+Cl2melDYDXWi+qs8h1UVUVkDNBTRG4CLsBe+w8k/JzdsE/PJBy3SUSWY5/+a0upqk6vhZ+zKun++uhrw+jrjljD/uoa/p6dgG1Iem4S7ic/NxXlaojLa15IXK7bWNkbr4g0BE4D+qnqiITtFZ1pj8HaPdoDlwFPqOrKhMcXYe0MiT+/PvaGvKImf0AGvou+bpu0fQes7SYTy4EmItK8hsVkGXaWt0vS9lbR12w9Ny6H+aUtF2fbYb2oysdzlF+2OjN5R1X9ArvscjN2qWpM0i5vA+dExaNcF+zD1mu1G7tS5W0NB5dvEJG9sHaPTL0cfe1RxT4pzxaiS30zgPOTHroAKAPerEY2l2f8jMTFlqp+LSLvAn8QkdXYG9tvscthzSs4ZDR2rf9L4MWkx/4IvAc8ISJ3Ytf1/wpMihra65yqfhn9PbeKSCn2Qe96qvGpX1U/FpFRwD9EZBesE8L2wHmq2i3abRbQSUQ6YWcwc6J2jWQ3ApOiy4MTgB9ilxLvSmpodwXKz0hc3F0EzAHGAUOwBvJxlez7DNZ4fq+qliU+oKozgc7YJZzHsMLyIHBe3cSu1EXAfOB+4E9Yj6mPq/mzrsHOwC7GuvmWAGsTHv8j1vj/EPAu1n16K6r6AtANOAJ4GhiIdYHuV81cLs+IL7XrCoWI/BwrJgdUo7utc64SXkhc3hOR3YEfYAPp5qvq6YEjOZdX/NKWKwS9sbEk3wH9A2dxLu/4GYlzzrka8TMS55xzNVIQ3X932mknLSoqqtax3377LU2aNEm9Y46IU944ZYV45Y1TVohX3jhlhZrlnTFjxjJV3Tnljqqa97e2bdtqdU2ZMqXax4YQp7xxyqoar7xxyqoar7xxyqpas7zAdE3jPdYvbTnnnKsRLyTOOedqxAuJc865GvFC4pxzrka8kDjnnKsRLyTOpTJ+PBQVcWKHDlBUZPdzVZyyurxREONInKu28eOhd28oLUUA5s2z+wDdu4dMtrU4ZXV5xQuJc1W54QYoLd1yW2kp9O0LH1d3dvc6MnRoxVlvuMELiatTXkicq8r8+RVv//pr+OMfs5sllcrmzavsb3CulngbiXNV2Xvvire3bg1lZbl1a906s7/BuVrihcS5qlx33dbbGjeGQYOynyWVQYMsW6JGjXIzq8srXkicq8rChfZ1991REfvUP2pUbrY5dO9u2Vq35vuLXGeckZtZXV7xQuJcZUpLYeRIOOccWLCAaS+/DHPn5vYbc/fuMHcu06ZMgZNOgtdfhw0bQqdyec4LiXOVue8+WLECiotDJ6me4mJYsAAeeSR0EpfnvJA4V5GyMigpgbZtoV270Gmqp3NnOOAAGDy48h5dztUCLyTOVeSFF2DWLBg4EERCp6meevVgwAB49114443QaVwe80LiXEUGD4bddoMLLgidpGYuvRRatrSzK+fqiBcS55LNnGlnJP36wbbbhk5TM02a2DQpjz1mHQWcqwNeSJxLVlICDRtunqcq7vr2tctzw4aFTuLylBcS5xJ99ZX11urRA3baKXSa2rHXXnD++XD33fDNN6HTuDzkhcS5RCNHwrp11kidT4qLYfVqGDMmdBKXh7yQOFdu3Tq4/Xbo1AnatAmdpnYddRQceywMGQKbNoVO4/KMFxLnyj30ECxeHN8BiKkMHAizZ8Mzz4RO4vKMFxLnwAbsDR4MBx8Mp5wSOk3dOOccmyts8ODQSVye8ULiHMArr8B778V7AGIqDRpA//4wbZr9rc7VEi8kzoF1+d1xR7jkktBJ6lavXja2xAcoulrkhcS5zz+HJ5+EPn1s/Y58tv320LMnPPggLFoUOo3LE15InBs61C77XHNN6CTZce21sHEj3HFH6CQuT2S9kIhIGxF5SURKRWShiNwiIvUzOL6eiMwQERWR0+syqysAX38N99wDXbvC7ruHTpMd++9vC16NGAFr14ZO4/JAVguJiLQEJgMKnAXcAvwSuDmDH3MFsEftp3MFafRoWLPGGtkLSXExLFsG48eHTuLyQLbPSPoAjYAuqvqiqo7AisgvRKR5qoOjQjQIuKFuY7qCsHGjXdY6/nhbd6SQnHgi/OQn1ujua5W4Gsp2IekMTFLV1QnbJmDF5cQ0jr8VeB14qQ6yuULz5JMwb17+DkCsioidhc2cCS++GDqNi7lsF5KDgFmJG1R1PlAaPVYpEfkRcDlwXZ2lc4Vl8GDYZx8488zQScLo1g1atfKuwK7GRLN4WisiG4BfqWpJ0vYvgXGqen0Vx04D3lbVX4tIETAHOENVK5zvQUR6A70BWrVq1XbChAnVyrxmzRqaNm1arWNDiFPekFmbzZpF26uv5rO+ffnyvPPSOiYfn9vW48axz5gxvDN2LKWtW2chWcXy8bnNFTXJ2759+xmqekTKHVU1azdgAzCggu0LgEFVHNcNWAw0j+4XYQ32p6fze9u2bavVNWXKlGofG0Kc8gbNetFFqs2aqX79ddqH5OVzu2SJ6nbbqV51VZ3mSSUvn9scUZO8wHRN4z0225e2VgLbV7C9BbCqogNEZBvgb8BfgXoisj1Q3jDfRESa1UVQl8cWLLAJGq+4Apqn7OOR33bZBS6+GMaNg+XLQ6dxMZXtQjKLpLYQEdkLaEJS20mCJsCewD+xQrQSeD96bALgkwa5zAwfDmVlNu+Us0b3tWth1KjQSVxMZbuQPAd0SjqL6AqsBaZVcswaoH3S7cLoseuB7nUT1eWl0lJbvOrss62h3cGhh0LHjlZg168PncbFULYLyQhgHfCYiHSMGsRvAv6pCV2CReQzERkNoKobVXVq4g14K9r1A1V9O7t/gou1ceNg5crC7PJbleJiWLgQHnkkdBIXQ1ktJKq6EjgJqA88jQ1GHAzcmLRrg2gf52pPWZl1dW3bFo47LnSa3HLqqXDggdYl2gcougw1yPYvVNUPgQ4p9ilK8fhcIE8XjXB1ZtIk+PhjuP/+/F1zpLrq1bN16q+5Bl5/Hdq1C53IxYjP/usKx+DBNjHj+eeHTpKbevSAli19gKLLmBcSVxj+9z+bCqRvX9h229BpclOTJnDVVfD44zBnTug0Lka8kLjCMGSILVp11VWhk+S2vn3tMtewYaGTuBjxQuLy31dfwX332aWbHXcMnSa37bmnXfq7+25YvTr1/s7hhcQVghEjYN06a0x2qRUXwzffwJgxoZO4mPBC4vLbunW2pOypp8LBB4dOEw9HHmndo4cMgU2bQqdxMeCFxOW3iRNh8WIfgJipgQOtwf2pp0IncTHghcTlL1Xr8tumDZx8cug08XL22dC6tXcFdmnxQuLy1yuvwH/+Y5+ufQBiZho0gGuvtefw3/8OncblOC8kLn8NHmy9tC6+OHSSeOrVC5o2tefRuSp4IXH56bPP7Pp+nz42fsRlrkUL6NnT2pkWLgydxuUwLyQuPw0bZpdn+vYNnSTerr0WNm60nm/OVcILics/X38N99wD3brBbruFThNv++0HZ55pY3HWrg2dxuUoLyQu/9x9N6xZY43sruaKi20Z3vvvD53E5SgvJC6/bNxol7VOOAEOPzx0mvxwwglw2GHWFdjXKnEV8ELi8ssTT8C8eT4AsTaJ2Nndhx/CCy+ETuNykBcSl18GD4Z994UzzgidJL906wa77uoDFF2FvJC4/PHOO/DGG9bTqL6v1Fyrtt3WesA9/zx89FHoNC7HeCFx+aOkBJo3t7EPrvZddRVst52flbiteCFx+eHLL+Hhh200drNmodPkp513hksugXHjYNmy0GlcDvFC4vLD8OFQVmaXtVzdGTgQvvsORo0KncTlEC8kLv6+/dbe2M45B4qKQqfJb4ccYjMpDx8O69eHTuNyhBcSF3/jxsHKlT4AMVuKi2HRInjoodBJXI7wQuLirazMGn+POMJW9XN1r1MnOOgg62rtAxQdXkhc3D3/PHzyiX1K9jVHsqNePRgwwNYpee210GlcDvBC4uJt8GDYfXc477zQSQpLjx6www6+VokDvJC4OPvgA5g8Gfr1swFzLnsaN7ZxJU88AbNnh07jAvNC4uJryBBbtKp379BJClPfvjaDwLBhoZO4wLyQuHhautSmNe/Rw5bTddm3xx5wwQUwejSsXh06jQvIC4mLpxEjYN067/IbWnExfPONLSTmClbWC4mItBGRl0SkVEQWisgtIlLlDHsicoiIPB/tv05E5ovI3SLiy98VonXrbOnXzp2tG6oL54gjoF07GDoUNm0KncYFktVCIiItgcmAAmcBtwC/BG5OcWgLYA5wHdAJuBHoCDwrIg3qLLDLTRMmwJIlfjaSKwYOhDlz4MknQydxgWT7TbgP0AjooqqrgRdFpDlwk4jcFm3biqq+AbyRsGmqiHwJvAD8CPh3Hed2uULVupyWT9Xhwjv7bJuapqQEunQJncYFkO1LW52BSUkFYwJWXE7M8Gctj756v89CMm0avP++fQr2AYi5oX59myzz1VdhxozQaVwA2S4kBwGzEjeo6nygNHqsSiJST0S2FZEDgb8A7wLv1EVQl6MGD4addoLu3UMncYl69oSmTX2AYoESzeJcOSKyAfiVqpYkbf8SGKeq16c4/nmsjQRgBvBzVV1ayb69gd4ArVq1ajthwoRqZV6zZg1Nmzat1rEhxClvplkbLVjAUZdcwryLL2ZugMWr8vm5rQ37Dx/O7k88wVsPPsj6nXfO6Fh/butOTfK2b99+hqoekXJHVc3aDdgADKhg+wJgUBrH/wA4GrgYO7OZATRMdVzbtm21uqZMmVLtY0OIU96Ms/bvr7rNNqoLF9ZJnlTy+rmtDZ9/riqiev31GR/qz23dqUleYLqm8d6e7UtbK4HtK9jeAliV6mBV/VRV31bV+7Ezk8OAi2o3ostJq1bZWIVu3WA37/Wdk/bdF846y8b4lJaGTuOyKNuFZBZJbSEishfQhKS2k1RUdR6wAti31tK53HX33baAVXFx6CSuKsXFsGIF3Hdf6CQui7JdSJ4DOolI4qLaXYG1wLRMflDU4L4jNr7E5bONG20+pxNPhMMOC53GVeX44+Hww60rcFlZ6DQuS7JdSEYA64DHRKRj1CB+E/BPTegSLCKficjohPt/F5G/iMg5ItJeRK4BJgGfY92HXT57/HGYP98HIMaBiP07zZoFL7wQOo3LkqwWElVdCZwE1Aeexka0D8ZGqidqEO1TbjpwPDAa+BdwLfAo8FNV/baOY7vQBg+26+9nnBE6iUtH167WjuVdgQtG1qcXUdUPgQ4p9ilKuj8BP/MoTG+/DW++aVPG169ySjaXK7bd1qaY//3vYeZMm4XA5TWf/dfltpISaN4cLr88dBKXiauugoYN7QOAy3teSFzu+uILePhhuOIKaNYs9f4ud+y0E1xyifXeWrYsdBpXx1IWEhHpISK+cpDLvttvt0ka+/cPncRVx8CB8N13MHJk6CSujqVzRjIG2A9ARDaJyFF1G8k5bMzIqFFwzjk2s6yLnzZt4JRT7APB+vWh07g6lE4hWQnsHn0v2FoiztWte++FlSt9AGLcFRfDokUwcWLoJK4OpdNrazJwn4h8jBWRsSJSaZdbVfUzFlczZWXWSHvkkXDssaHTuJro1AkOPti6Al98sU/9n6fSKSQ9gWuAA4HDsZHkX9VlKFfgnnsOPvkExo/3N564E4EBA6BPH1uv5IQTQidydSBlIVHVUuDvACLSEbhBVd+v62CugA0eDHvsAeefHzqJqw2XXALXX2//rl5I8lI6vbY2iciR0d2pQIXL4TpXKz74AF56Cfr1g222CZ3G1YbGje2M5Mkn4fPPQ6dxdSCdxvb1wHbR9z2AzFascS4TJSXQqBH07h06iatN11xjMxMMGxY6iasD6bSRfAjcJCJPYL22zhORylbMUlW9s9bSucKydKm1i1x+OeywQ+g0rjbtsYfNwTV6NNx8M7RoETqRq0XpFJL+wEhsckUFrqtiXwW8kLjqGTEC1q2zxlmXf4qL7YPC6NHwi1+ETuNqUcpLW6r6hqr+UFW3wc5Ifqqq9Sq5+ax6rnrWrYM77oDOneGgg1Lv7+KnbVto1w6GDrU1ZlzeyHSurfbYpS7nateDD8KSJT4AMd8VF8O8edbw7vJGRtPIq+o0ABE5GmgH7IAtd/uaqr5d+/FcQVC1RvZDDoGOHUOncXXprLNgn32sK/C554ZO42pJRoVERJoADwOdgE3Acmy52/oi8jxwfjTuxLn0TZ0K778Pd93lAxDzXf36cO21dmby7rs2e4GLvUwvbd0GHAN0Axqq6m5Aw+j+McBfazeeKwiDB9u04927h07isqFnT1sWoKQkdBJXSzItJOcCv1HVh1W1DEBVy1T1YeC3gA9Fdpn59FN45hm4+mobP+LyX/Pm0KsXPPQQLFgQOo2rBZkWkhbAF5U89gXQvGZxXMEZOhQaNLBC4grHtdfa5Jy33x46iasFmRaS94GrRba8kB3dvzp63Lm0NFizBsaMgQsvhN12Cx3HZdM++1jD+8iRUOrNqnGXaSG5HmtonyUifxGRYhH5M/ARcEr0uHNVGz8eioo47owzbAGrAw4InciFUFwMK1bAXntxYocOtoDZ+PGhU7lqyLT778sichjwB6w9ZDdgEfA20EVVfYyJq9r48TaPVmkp35/W/ulP9ibije2FZf5866W3YoW9FubN2zzHmr8WYiXTMxJU9UNV7aaq+6lq4+jrRV5EXFpuuGHrSxmlpbbdFZYbbrAxRIn8tRBLGRUSEfm7iLSpqzCuAMyfn9l2l7/8tZA3qtP99wMReUdE+oiIT+HpMrP33pltd/nLXwt5I6NCoqr7AB2BWcDfgEUi8kC0cqJzqQ0aZKObEzVubNtdYRk0yP7tE/lrIZaq00YyRVV7YA3t/YE9gUkiMk9EbhaRfWs7pMsjJ5xg4weaNUNFoHVrGDXKG1cLUffu9m/fujXft5T84x/+WoihjAtJOVVdo6qjgRuB14G9gN8Bn4jIkyLSupYyunwyfLj11Pnvf5n28sswd66/cRSy7t1h7lzeHTPG7i9bFjaPq5ZqFRIRKRKRG0VkNvACsAbrDtwMOBMoAibUVkiXJ7791j6Bduli3X2di5QWFUGnTjbSff360HFchjLttXWJiLwMfAZcCowB9lHVn6vqo6q6TlWfBa4FKluO1xWqe++FVatg4MDQSVwuGjgQFi+GiRNDJ3EZyvSMZBSwGOikqvuq6q2q+mUF+30C/LGiHyAibUTkJREpFZGFInKLiFS5sqKIHCkiY0Tks+i4j6MzooYZ5nehlJXBkCE2bfixx4ZO43JRp05w8ME2G3Ty+BKX0zIa2Q7srqorU+2kqouAm5O3i0hLYDK2yuJZwH7AP7CC9vsqfmTXaN+/Ap8CPwJujb766jhx8Nxz8MknNrLd1xxxFRGBAQOgTx949VXrmOFiIdMpUlIWkRT6AI2w6VRWAy+KSHPgJhG5LdpWkb+q6lcJ96eKyHfASBFprarzapjL1bWSEth9dzjfVxpwVbjkErj+ejsr8UISGxk3totIVxGZLCLzRWRp8i3F4Z2BSUkFYwJWXE6s7KCkIlLuvejrLhn9AS77PvgAJk+Gfv1gm21Cp3G5rHFjOyN58kmYPTt0GpemTBvbLwLuxRrb9wSeAp6Jfs5qYHiKH3EQNpjxe6o6HyiNHsvEsUAZ8HGGx7lsKymxRavKJ+RzrirXXGODVocODZ3EpSnTM5JfYW0TfaP7d6hqT2AfYBlWEKrSElhVwfaV0WNpEZFdgRuA+6q4HOZywdKl1i7SowfsuGPoNC4O9tgDunaFe+6B1f7fOw5EM+gdISJrgNNVdaqIbABOVtWp0WPnAINVtaiK4zcA16nqkKTtC4Cxqppy2k8R2RZrsN8TaFtZu42I9AZ6A7Rq1arthAnVG9ayZs0amjZtWq1jQ8i1vK3vvZd9xo7lnXvvpTRpDqVcy5pKnPLGKStsnbfZxx/Ttk8fPrvmGr7MsXa1uD+3mWjfvv0MVU09lENV074BC7CuvwBzgasTHusCfJPi+KXAjRVsXwP8Ko3fL1ibynLgoHRzt23bVqtrypQp1T42hJzK+913qq1aqXbuXOHDOZU1DXHKG6esqpXkbddOtahIdePGrOepSl48t2kCpmsa77GZXtqajnW5BWsf+YOIXCkil2KTOL6d4vhZJLWFiMheQBOS2k4qMRjrNnyWqqazvwtpwgRYssQHILrqGTjQptB58snQSVwKmRaSPwPliwX8AXgHuAMb4b6M6FJSFZ4DOolIs4RtXYG1wLSqDhSR32GTRF6sqq9lmNtlm6o1srdpAyefHDqNi6Ozz7apdEpKQidxKaRVSESkkYicCxwHNBCRVqq6SlXPApoC26vq0aqaqr/eCGAd8JiIdIzaMW4C/qkJjebRCPbRCfcvAv4EjAMWiMhPE247Z/D3umyZNg3+8x/7VOkDEF111K8P/fvb4MQZM0KncVVIWUiiaeFnAg9jl6/uAz4WkVMA1ObXSqtrhVrD+ElAfeBpbPT7YGwG4UQNon3KnRJ9vQx4M+l2Wjq/22VZSYn10rr44tBJXJz16gVNm9oARZez0jkjuQ0br3E80Bg4BBsMOLI6v1BtzfcOqtpIVXdT1f9T1U1J+xSp6mUJ9y9TVankNrY6OVwd+uwzeOopG1jWqFHoNC7OWrSwYjJxIixcGDqNq0Q6heQY4Peq+rqqfqeqHwFXAXuLyG51G8/F0tCh0KAB9O2bel/nUunfHzZtsinmXU5Kp5DsBiS3fXyOdcXdtdYTuXhbtcoGknXrBrv55wxXC/bbD846C0aOhNJUY55dCOn22vI5nV16Ro+2Bay8y6+rTQMHwvLlcP/9oZO4CqRbSCYlTcy4KNr+UoaTNrp8tnEjDBtms7YefnjoNC6fnHACHHaYdeLwtUpyTjrTyG+1rohzFXriCZg3z3vYuNonYmcll14KL7xgi2C5nJGykKiqFxKXnpIS2GcfOPPM0ElcPuraFX7zG3udeSHJKRmvR+Jchd59F15/Ha691gaSOVfbttvOpph//nn46KPQaVwCLySudpSUQLNm0LNn6CQun/XpYwVlyJDU+7qs8ULiam7BAnjoIbjiCmjePHQal8923tmW4x03znpxuZzghcTV3PDhUFZmA8ecq2sDBsDatTauxOUELySuZkpL7T/02WdbQ7tzde3QQ21G6dtvh/XrQ6dxeCFxNTVuHKxc6QMQXXYNHGhzbz38cOgkDi8kribKyqzRs21baNcudBpXSE49FQ480Aco5ggvJK76Jk2CWbN8zRGXffXqWVvJ9OnW7dwF5YXEVV9JiU3MeMEFoZO4QtSjB7Rs6Sso5gAvJK56Zs60qSr69YNttw2dxhWiJk2gd294/HGYMyd0moLmhcRVT0kJNGxo/5GdC6VfP7vMNXx46CQFzQuJy9xXX8F999mlhZ12Cp3GFbI994Tzz4e774ZvvgmdpmB5IXGZGzkS1q2zxk7nQhs4EFavhjFjQicpWF5IXGbWr7eBYJ06QZs2odM4B0cdBccea13RN20KnaYgeSFxmZk4ERYv9gGILrcMHAizZ8PTT4dOUpC8kLj0qdqiVQcf7OtBuNxyzjmw997eFTgQLyQufa++Cu+95wMQXe5p0MCo1nPKAAAfqUlEQVQmDZ02zV6jLqu8kLj0lZTADjvAxReHTuLc1q64wsaW+FlJ1nkhcemZPdvWZO/TBxo3Dp3Gua1tvz1cfjk8+CAsWhQ6TUHxQuLSM3SoLaF7zTWhkzhXuQEDYONGuPPO0EkKihcSl9rq1XDPPdC1K+yxR+g0zlVu//3hjDOskKxdGzpNwfBC4lIbPdpGDXuXXxcHAwfCsmXwwAOhkxQMLySuaps22WWtdu3giCNCp3EutZ/9DH78Y+uq7muVZIUXEle1J5+EuXP9bMTFh4i9XmfOhJdeCp2mIHghcVUrKYGiIluT3bm4uPBC2GUXOytxdS7rhURE2ojISyJSKiILReQWEamf4phtReRvIvKqiKwVET9fzYYZM2wQ4rXXWo8t5+Jiu+2sh+Gzz9oqnq5OZbWQiEhLYDKgwFnALcAvgZtTHNoYuAIoBd6oy4wuQUkJNG0KPXuGTuJc5q6+2hZdGzo0dJK8l+0zkj5AI6CLqr6oqiOwIvILEWle2UGqugrYQVU7AY9nJ2qBW7gQJkyAXr2gRYvQaZzL3C67QPfucO+9sGJF6DR5LduFpDMwSVVXJ2ybgBWXE6s6UNW7X2TV7bdbj63+/UMnca76Bg6E0lK4667QSfJatgvJQcAWFyxVdT52yeqgLGdxlSkttcWrzjoL9tsvdBrnqu9HP4IOHWDYMNiwIXSavNUgy7+vJbCqgu0ro8dqjYj0BnoDtGrViqlTp1br56xZs6bax4ZQG3l3e/ppDly+nPdOPJGv6/BvL8TnNlvilBXqNu+OHTvyw5df5sNbb2Vphw41/nn+3FZAVbN2AzYAAyrYvgAYlObP6Ed0pSvdW9u2bbW6pkyZUu1jQ6hx3rIy1YMPVj3sMPu+DhXcc5tFccqqWsd5N21S/cEPVI86qlZe04X03ALTNY332Gxf2loJbF/B9hZUfKbisu2FF+Cjj6C42NcccfmhXj2bzPGdd+Ctt0KnyUvZLiSzSGoLEZG9gCYktZ24QEpKYNddbYJG5/LFpZfaNPM+QLFOZLuQPAd0EpFmCdu6AmuBaVnO4pJ99BE8/zz07Wv9753LF02bwpVXwqOPwrx5odPknWwXkhHAOuAxEekYNYjfBPxTE7oEi8hnIjI68UAR6Swi5wE/ie6fF91aZy9+nhsyxEYEX3VV6CTO1b5+/exy7fDhoZPknawWElVdCZwE1AeexgYjDgZuTNq1QbRPojuBh4Fe0f2Ho1v7uspbUJYvh3Hj4JJLYOedQ6dxrvbtvTece66NKVmzJnSavJLt7r+o6odAlX3wVLUonW2uFo0caQsBDRgQOolzdae4GB56CMaOtTMUVyt89l8H69fbSPaTT4ZDDw2dxrm689OfwtFH22XcsrLQafKGFxIHDz9sc2sVF4dO4lzdKy6Gzz6Df/0rdJK84YWk0Klal98DD4ROnUKnca7unXsu7LWXdwWuRV5ICt3rr8P06Ta5XT1/ObgC0KCBtY9MmQLvvx86TV7wd45CV1ICLVtaby3nCsWVV0Ljxvb6dzXmhaSQzZkDjz9u40aaNAmdxrnsadkSLrsMHngAliwJnSb2vJAUsmHD7HJW376hkziXfQMGWI/FO+8MnST2vJAUqtWr4e674fzzYc89Q6dxLvsOOABOO80KyXffhU4Ta15ICtWYMfDNN9bI7lyhKi6GpUvhwQdDJ4k1LySFaNMmGDoUjj0WjjoqdBrnwunQAX74Q+sK7Kt5V5sXkkL09NMwe7YPQHROxM7KP/jAugO7avFCUohKSqB1azj77NBJnAvvootsolIfoFhtXkgKzXvvwbRp0L+/DcxyrtA1bAhXXw3PPAOffho6TSx5ISk0JSU2ZqRXr9T7Olcorr7aFnMbMiR0kljyQlJIFi2y3ik9e9qyo845s+uucOGF1ptx5crQaWLHC0khueMO2LgRrr02dBLncs/AgVBaauOrXEa8kBSKtWthxAg44wzYf//QaZzLPT/5CfzsZzbjw8aNodPEiheSQjF+PCxb5l1+natKcTF88QU89ljoJLHihaQQlK858uMfw4knhk7jXO467TTYbz+fFThDXkgKweTJMHOmfdoSCZ3GudxVv75N5vjmm/D226HTxIYXkkJQUgKtWkG3bqGTOJf7Lr8cWrTwAYoZ8EKS72bNgmefhWuuge22C53GudzXtClccQU88oi1l7iUvJDku6FDrYD06RM6iXPx0b+/tS0OHx46SSx4IclnK1bAvfdC9+6wyy6h0zgXH61bQ5cuMGoUfPtt6DQ5zwtJPhs1ygZY+ZojzmWuuBhWrbIPY65KXkjy1YYNdlp+0km23oJzLjPHHANHHmnzb5WVhU6T07yQ5KtHHoEFC3wAonPVJWL/fz75BJ57LnSanOaFJB+pWtfFAw6Azp1Dp3Euvs47D/bYwwcopuCFJB+9+Sa8+64NrKrn/8TOVds220C/fjao94MPQqfJWf4uk49KSmya+B49QidxLv5694ZGjXytkipkvZCISBsReUlESkVkoYjcIiL10ziuhYiMEZGVIvK1iIwXkR2zkTlOtlu8GB591F78TZuGjuNc/O2wA1x6Kdx/PyxdGjpNTspqIRGRlsBkQIGzgFuAXwI3p3H4ROBnwBXAZcCRwBN1kTPO9njiCWsk7NcvdBTn8seAAbBunS3F4LaS7TOSPkAjoIuqvqiqI7Ai8gsRaV7ZQSJyDNAJuFRVH1XVx4GLgXYi0rFOko4fD0VFnNihAxQV2f1cNn487L03e02caCPZX3kldCLn8sdBB8GPfgQ33xyv94QsvYdlu5B0Biap6uqEbROw4lLV/OadgSWq+v27o6q+A8yJHqtd48fbpaF58xBVmDfP7ufqC6c87xdfIGCDEHM5r3NxM348fPwxlJXF6z0hS+9hDerkp1buIODlxA2qOl9ESqPHnq7iuFkVbP8oeqx23XCDvRknKi2Fyy6DP/2p1n9djX3yydYrupWW2t/RvXuYTM7lkxtusEtbifw94XvZLiQtgVUVbF8ZPVad4/at6AAR6Q30BmjVqhVTp05NO+SJ8+dT0aodunEjX+28c9o/J1t2/vDDivPOn8+0DP7ubFuzZk1G/y6hxSlvnLJC7uf194QUVDVrN2ADMKCC7QuAQVUc9yLweAXbxwOvp/q9bdu21Yy0bq1qw/q2vLVundnPyZa45Y1MmTIldISMxClvnLKqxiBv3P6P1VJeYLqm8d6e7TaSlcD2FWxvQcVnHKmO2z7FcdUzaBA0brzltsaNbXsuilte5+Imbv/Hspw324VkFkltGiKyF9CEittAKj0uUlnbSc10724z57ZujYrYlNKjRuVue0Pc8joXN3H7P5blvNkuJM8BnUSkWcK2rsBaYFqK43YVkXblG0TkCKx9pG5mU+veHebOZdrLL8Pcubn7gikXt7zOxU3c/o9lMW+2C8kIYB3wmIh0jBrEbwL+qQldgkXkMxEZXX5fVd8EJgHjRKSLiJyNtY+8pqqTs/oXOOec20JWC4mqrgROAupjXX1vBgYDNybt2iDaJ1E37KzlHmAcMAM4py7zOuecSy3b3X9R1Q+BDin2Kapg2yrg8ujmnHMuR/jsv84552rEC4lzzrkaERtzkt9E5CtgXjUP3wlYVotx6lqc8sYpK8Qrb5yyQrzyxikr1Cxva1VNOXS/IApJTYjIdFU9InSOdMUpb5yyQrzyxikrxCtvnLJCdvL6pS3nnHM14oXEOedcjXghSW1U6AAZilPeOGWFeOWNU1aIV944ZYUs5PU2EuecczXiZyTOOedqxAuJc865GvFC4pxzrka8kDjnnKsRLyTOOedqJOuz/7raEa0s+XNAgIdVdbmI7AlcB+wHzAVGqeoH4VKCiPwGeDZ0jnSJSCOggap+k7BtZ6Af0AYoA/4D3KGqX4dJ6Vxu8e6/ERERbH2T04CDgR2ATcAS4C1grKp+Ei7hZiJyFPAC0BTYCKwAOgHPYplnAocCuwIdVfXVQFERkTJAsSWRHwAmqupnofKkIiLPAp+q6oDo/jHYKpxl2Bo4ArQF1gMdVHVmwKyHAY1U9Y2EbacCv2Nz0XsfuClxn1wR/Z87Azgce41Mxz505PSbkog0x+au6qCqr4XOA99n6gBsC/xLVb+NPgD1xVaSnY19sFxYJ78/x//NsiJ6wp/F3iCWYKs47oG9uJ/D/iEOBG5V1VtD5SwnIi9iZ5PnAN9ii4Odjb3RnaeqG0RkO+AJoKGqtg+YtQz4K/BD4GQs97+xovKQqi4Ila0iIrIM6KWqT0b338Ke47PLz1JEpAXwFPCdqnYKmPUt4GlVHRTd7wncDUwBXsaK3knA8cC55X9ToKxvYM/rR9H9ltiHobbAmmi3ptiHtk6JZ4QhiMg1VTzcCPgbMAT4FEBV78hGroqIyP7AS8Be0aY5wCnAi8D2wOfY+9daoK2qflnrIVS14G/Ag9gL4ocJ23YHngceje6fiL3ge+ZA3uVA54T7u2CfPk9J2u80YFngrGXAUdH3LYHe0Yt+Y3SbGm3bMfTzGmUsBU5IuL8++XlNeG6/DZx1dWI24DNgWAX7jQDez5XXQXR/NHYmfWrCtlOBlcDgHHgdlGFn92WV3BIf2xQ460PYmef+2JWU+6L3szeAZtE+O0X7jKyLDN7YbjoDv9WE6/hqp4B9gLNFZDdVnQb8CRgQKGMijW6J90naVtH9oFR1paqOUtWTgD2BX2Kn4iOAhSLyr6ABzf+AxDO4Jdh/zmQ7YkUnpLKk+62BRyrY7xHsE2kuORO4RVWfL98QfT8I6BIs1WZPAUuBXkB9Va1XfsNeDwL8LNqWvCx4trUDBqnqZ6q6Avg91k76d43O7FR1GVDClq/tWuOFxAj2CSPZpuixFtH9t4EDshWqCjOA60SkmYjUA64HFgBXi0h9ABFpAFyDvTHmHFVdrKpDVPVYYB/gRuwsMLS/AL8VkZ7RczgI+JuInCwi24rIdlE7xJ+xT4IhvQp0T7g/E6houvAjsddHLtkeaxNJNgNr2wtKVc8GLgV+BbwrIsclPhwmVaVaAosT7pf/WyevwTQb+wBX67zXlpkM/FFE/quqs+H7a7hDsX+g8kb2pkAu9NS5Abv+uQK7PFSKNbQ9AnwqIuWN7btjlwtymqrOw97A/5IDWR4Tkf7Yp7fBwMfYB4nyT86Kfbh4CnuTCel64PXow8QwrJH9XhHZAbtkCNZGMhD4bZCEWzpXRMoL3UqgogWTdsIu2QWnqi+IyI+w5+9fIvI81isyaPtNBZZiZ6PlNgEjsbPpRLtQR9m9sR2Ius0+j53+z8Oui++DNbpfqKrPRfvdhq0Y1jVU1nJR5tOxDwOPquoiEdkV+DWb/467VfXfAWMiIjcCd2kd9RapKyKyI9AVOAr7hFwPK9wfAc+o6oyA8b4nIj8B7gSOZnORI+H7ldglpCFhEpqo00WysaraM2m/kUAbVT0+O8nSE/3fug277DYSKy7tVfWVoMEAEXkCWJH8XFaw3zDgYFXtWOsZvJCY6JLQBcCPgYZYw+UD0TVH53KaiByMFZPkoveGqm4ImS0TInIl8Lmqvhw6S0Wi7uCDsQ9rp2kOdKsWkVZAY1Wdk2K/X2CdLl6q9QxeSPKPiNRX1YrafHKGiDTEGgTLgM9y8c0uaiPZl4QxRao6P2wq53KPN7YnEZFDRORcEblCRHpF3x8SOlcyEekiIk+IyLMicka0rauIzAXWi8i86NNdUCJycTS+ofx+AxH5C/aJ+b9YZ4AVIpIL1/ABEJG2IvIUdj35I+B1bHzDHBFZICK3iEjjoCHziERC56iIiDRK/rcWkZ9E7wttQ+XKOSH7P+fSDeiJtStU1Hd8EzblyOWhc0ZZL4hyvQY8iTW2X4m17YzGRrM+GOXuFDjrh8DVCff/EeX9P+A4rOviTdhgqetz4Lk9BWsbm471zLoJG5S6Psr8S6x31H+AljmQ93RsXM4HwEQSxsAk7HM04cc6nEI0piFh29nY4NSNwIboOT8t9HMaZWsBPB7l2gjcBdQH7k16X3gd2Cl03jT/pnPr6nUQ/I/LhRvQP3rB3I6NAt4petHUj75vBwyP3mD65kDed4ERCfe7R9n+kbTfGGBy4KylwIkJ95cCAyrY7zpgXg48tzOAeyt5jczFzuIbRm+AdwTOenLCm9nwKPumqFhLwn65UEg2seWAxHOiN+M3on/766LvN1LBANAAeYdi06D0B3pEHx4eBb6IiuLO2PizBcCdofOm+TfVWSHxNhJARGZjb8y3pdjv10AfVd03O8kqzbEa6KKqk6P7LbDeOR01oZEyuuQ1UlWDjc8QkUVAP1V9NLq/DjtLmpq038nAU6raKPspt8ixFjhTVV9M2t4Sm1HgEFX9SER6AH9V1d1C5IwyvYbNC3Z5wrae2Jvgi1iPw+9E5Gis0T3YwLmo19ZPVfWd6P6/gQWqekbSfs8CTVT1xAAxE3PMAf6kqndF9w/DCvXlqnpvwn5XYmfS+4RJCiJyT5q7tsYGUdb668DbSMyuwDtp7PcOOTBYCuvamfhiKJ+raFXSfmuwgV8hPYUNntw2uj8ZuLCC/S7EPvWFthTruZfsx9jzXj6OaB6bB6qGcihwf+IGVb0Hm87np8DL0ZiSXHQo1o022ShsEsfQdmHz+DGI5tTC5q1K9BkVj4fJpkuxs6Qfpri1ruwH1JQPSDT/Ba4UkVdUtaL+7uUzlV4Z7RvaPGx210kAqrop6pb4UdJ++7LliNcQfoeNwP6fiNwNPA38VUQOZfOguQ7AYdhMsKGNAm4VkSZY28N6bGT4DcAU3TweZl8gdA+u74AmyRtVdUY0EnsSdrnopiznqkzi5Y+v2fwBKNG35MYH3DlYQZ4W3T8euxR3LNY2We44wr8OPgXeUdUeVe0kIudh7Wi1zguJ+SU2IPFDEXkMm/J8FfbC3x44CLumuye5MVL8MZKmOlDVtyvY7yK2fNFnnaquEJGfYm/Ev8A+6QEcE93WY5dhjlfVd8Ok3ExVB0WXYX6LTdsC9jp4EBuEVm4DNvdaSP/FrtM/lfyAqs6OismzwNgs56rMJBHZGH3fAvgJmz9MlDsIWJTNUJUYAQwRkR9iRe8C7EPRH0SkKTYB4uFAMRB6RvC3sAKXSuKA1VrlbSQREdkPGxV+KpunYy73BdZz52+qmnxqm7NEZG9glarmxJQTACJSxJaD5j7X3BxDsg02zqUhMDuXnsNyInIVNk3KYVrJwNnozOpxrP0s2Cf9aIaDZJ+q6gNJ+02NtudC1/VrsUuu22CzRIwQkQuxNqjySTtHAb8J+RqOuiEfp6pDU+y3E9bGN62q/aqVwQvJ1qJ+4+VtC6tUNfQsr865HBFd5t5JVb8KnSVXeCHJM9Fp97+B7rlwqUhiuHStxGQZY+dyhReSBNEbyC7Ax6q6VUNgdGr4c1Udl/VwW+b4eRUPN8Ea1H5LNIW8qj6bjVwVkRgtXQvxWsY4XdE8XOer6i2BcwRdDramojORxKWBZ2B/R/A3UbFZlc/F/j+NVdVZIvJj4GY2f/i5XRPWf6lVoQfJ5MIN2A54GHuj2IQ1pI4GWiTtF3xgV5QjTqu3LQPOSrj/FtYbqlnCthZY75hJOfDcvogtVbs9dm18OPAlNoPANgmvl+ewXlzBX79p/E11NhAtgwz7Y70Ny1+Xn2NvcLOxYv0uNn38EmDPHHjO3sBmyi2/3zLKWBblXM3mAZXNQuWMsnXCPogtjp7X1dgCViuxwaq3R//vNmFLRtd+htD/YLlwA/6A9dK6ElsYaED0gv4U+EHCfrlSSGZgPVsux/qGJ95+FL3ALyjfFjhrbJaujXLEaRnjvdO89Qn9uiUHloPNMG9slgaOisXD2EqOYB0wVgKjk/a7D3irTjKE/gfLhRvW3bdf0rZdgVeAr4Bjom25UkgEW+d8KTZtwz4Jj7WI/hNsNedSoKzvADcm3P8C6FbBfj2Ar3Ig77KkN4udo+fz5KT9fp4DhaT87DPVLRfOTBcCFyTcbx3l6pK03+XAJznwOkguJF8BAyvYL/jUPlj35I4J91tG+Tsk7XcK1nmo1jP4OBKzF0kDDVV1sYichFXxySLSndzo347aq2KUiDwE/BH4r4gMj77PNX8BxovIF8A4Ni9duxy7nCXYaXguLF0Lm5cxfh0bHJe4jPHLaoM/c2UZ42+Al4G7U+zXDuvaHlLw5WBrKJeXBl7LlgNTy79Pnm6oMTaItdZ5ITELgR9gZyDfU+sb3k1ESrBTx6CN7MlUdRXQT0RGYX3bPwX+Sg6tKa3xWroW4rWM8TtYO96/qtopWvsltODLwVZDXJYGfh0bKPlplOXv2Kzbv4lm6/gmmo/v11jhq3Xea4vvJz3bV1V/VsU+v8M+TasGnPyuKiLSDVsOdE9scrbgy4CWk5gsXQuxWsb4/4Deqpo8gDZ5vxOAm1W1fXaSVZgh+HKwmZAYLQ0sIvtjc9iVvw7mYmf5j2AzBcwDirAPRu1V9T+1nsELyfdd57oCf1HV5VXsdxF2rfzyyvYJLbrs0gRYozm+SqIrHJIDy8HWBcmRpYGj8WPHYT0NX1LVtdHA6ivY/OHnAVX9sk5+vxcS55xzNZELs2y6OiIid4nI6NA50hGnrBC/vM7VJW9sz4CI3AXUU9VeobOkqT3x+bAQp6wQo7wiMhm7+nBS6CypxCkrxCtvXWb1QpKZ2Lx5AKjq/qEzpCtOWSF2eYX4vG7jlBXilbfOsnobSR6Lun3uoqqhF95JKU5ZIX55natLcamkOUFEGkZrfMTFadhKb3EQp6wQo7wisk1cXrdxygrxyluXWb2QZCY2bx6uMIhIXxH5XES+EZG3ReSSCnY7nBx43cYpK8Qrb+is3kYSQyKSbp/1ikbiZlWcskK88kYDUIdhywC/h40jGCsiZwGXqOrakPkSxSkrxCtvLmT1NhIyfvNoE3pku9i61x9j0yBUZQ/g6JB545QV4pVXRKYDL6vqrxO2nQSMx0Y3n6a2KNfRwBueNX1xypsLWb2QEK83DwAR+Q+2+FbXFPudB0wM/CKPTdYoR2zyisg3wBmqOjVpexG2Xkp9oDM2H1ToN7vYZIV45c2FrN5GYv4H/E9Vz6/qBvwzdNDI28BP09ivfELEkOKUFeKV92vszWELqjoXOBabEv8N4MjsxqpQnLJCvPIGz+pnJHw/8dqpqto6xX7nYmt4By3AIrIfcIiqPpViv0ZYF9XkqbqzJk5ZoxyxySsiTwLfqOrFlTzeCJu4rzOBJxuNU9YoT2zy5kJWLyTE683DuXIicj5QDJyuqisq2ac+cCc22eg+2cyXlCM2WaMsscmbC1m9kDjnnKsRbyNxzjlXI15InHPO1YgXEldQROQyEZkRjQBeKSLviUid9MYTkQNE5CYR2T6NfW8SEU24LRSRR6P2u1THXhYd07R2kjuXGS8krmCILZd8NzAJ6AL0AJ4EzqyjX3kAcCOQspBEvgaOiW7XAT8BXhKRJimO+1d0TGk1czpXIz5Fiisk/YCRqnp9wranReTmUIGSbFTVt6Lv3xKR+cCrwM+Bh5N3jnri1FfVr4CvshfTuS35GYkrJNsDi5M3akLXRREpii4TXSQi90WXwJaKyI3Jx4lIh2iCvO9EZImI3FF+eUlEfgY8He06J/qZczPMOyP6WhT9zLEiMl1EzhaRmcB3wNEVXdoSkUYicpuIzBORdSIyR0T+nJT/ChGZGT0+T0R+jXPV4GckrpD8G+gffdJ/RlWXV7Hv34BngPOAE4AbRWSZqt4OICJtgOeBF4Fzgb2AvwD7AqdGv+s64O/YZbRFwLoM8xZFXxcnbbsNuAVYgs3mukU7iogIdsnuGOBWrCDtARyfsM+vgD9FP2sq0Ba4VURKVXV4hjldoVNVv/mtIG7Aj4DZ2PQmZcBM7A25ecI+RdHjLyQdexewAFtqGWAC8Cl2aal8nwuiY4+J7p8e3S9KI9tN2FQWDaLbAcAUYDWwW7TP2Ojn/STp2Mui7U2j+52i+2dW8ruaA2uAG5O234IVrfqp8vrNb4k3v7TlCoaq/hc4GGtcvwObK+v/gOkV9Hh6POn+Y8DuwJ7R/aOAx1V1U8I+jwIbgXbVjLgjsCG6fYyd3XRV1UUJ+yxQ1f+k+DkdgBVa+UwNxwBNgIdFpEH5DXgZaMXmv9G5tPilLVdQVHUd1nbxNICI9MJ6cvUChiTsujTp0PL7uwHzo69Lkn72JhFZDuxQzXhfAx2xs4nFwEJVTZ56YslWR21tR+xSWmXKJ/ibWcnjewE+DZBLmxcSV9BUdbSI3AYclPTQLpXcX5TwdYt9ol5UOwIVzneUho2qOj3FPunMabQcK3SVKc93OhUXpo/T+B3Ofc8vbbmCISLJxQER2RlowdZvqOck3S9vMP8yuv82cE5UPBL3aQC8Ft1fH31tWIPY1fESsIOInF7J428Ca4HdVXV6BbdvshfV5QM/I3GF5INoyu0XsEtVrbGeVaXAvUn7HhItL/Ao1murFzBAVcuix/+ILWv6hIjcibUr/BWYpKpvRvuUf7K/SkQmAKWq+kHd/GlbeBEbdPmAiNyC9SDbDThBVa9S1VUichMwRERaA69gHyoPANqranIRda5KXkhcIbkFOAsYirVjLMYW/OmqqnOS9v01dunnUWy8xq3A991iVXWmiHTGutA+hvWuejA6rnyfeSJyHXAt0B87mymqiz8skaqqiJwTZR6ILRG9EHggYZ/bRGQhNv34L7G/8RNgYl3nc/nHp5F3LkG0POkcbOnSZ8KmcS4evI3EOedcjXghcc45VyN+acs551yN+BmJc865GvFC4pxzrka8kDjnnKsRLyTOOedqxAuJc865Gvl/T2D1lNDCAdoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot exact payoff function (evaluated on the grid of the uncertainty model)\n",
    "x = uncertainty_model.values\n",
    "def payoff(x):\n",
    "    if x <= strike_price_1:\n",
    "        return 0\n",
    "    elif x < strike_price_2:\n",
    "        return x - strike_price_1\n",
    "    elif x < strike_price_3:\n",
    "        return strike_price_2 - strike_price_1\n",
    "    elif x < strike_price_4:\n",
    "        return strike_price_2 - strike_price_1 + strike_price_3 - x\n",
    "    else:\n",
    "        return 0\n",
    "y = [payoff(x_) for x_ in x]\n",
    "plt.plot(x, y, 'ro-')\n",
    "plt.grid()\n",
    "plt.title('Payoff Function', size=15)\n",
    "plt.xlabel('Spot Price', size=15)\n",
    "plt.ylabel('Payoff', size=15)\n",
    "plt.xticks(x, size=15, rotation=90)\n",
    "plt.yticks(size=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "exact expected value:\t0.3569\n"
     ]
    }
   ],
   "source": [
    "# evaluate exact expected value (normalized to the [0, 1] interval)\n",
    "exact_value = np.dot(uncertainty_model.probabilities, y)\n",
    "print('exact expected value:\\t%.4f' % exact_value)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evaluate Expected Payoff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set number of evaluation qubits (=log(samples))\n",
    "m = 5\n",
    "\n",
    "# construct amplitude estimation \n",
    "ae = AmplitudeEstimation(m, iron_condor)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# result = ae.run(quantum_instance=BasicAer.get_backend('qasm_simulator'), shots=100)\n",
    "result = ae.run(quantum_instance=BasicAer.get_backend('statevector_simulator'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Exact value:    \t0.3569\n",
      "Estimated value:\t0.3428\n",
      "Probability:    \t0.9697\n"
     ]
    }
   ],
   "source": [
    "print('Exact value:    \\t%.4f' % exact_value)\n",
    "print('Estimated value:\\t%.4f' % result['estimation'])\n",
    "print('Probability:    \\t%.4f' % result['max_probability'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZwAAAEPCAYAAAB2s3LUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAHjVJREFUeJzt3XvUHFWd7vHvw0UIt3APDDJE8BJhZHnkFWEOI4ncwzlyEU0WembFQaMeFWYWOiAiBHRYBkdAF+MClk44nNGEGWA4wyWEcHnDXQgSxEkCBg13EZyXxBiIQH7nj13BSr39dlffqvN2ns9avd7uXbt27dqp9K+rau9digjMzMy6bZNeV8DMzDYODjhmZlYJBxwzM6uEA46ZmVXCAcfMzCrhgGNmZpVwwDGrQ9IMSYPZ+0FJM5pcf6KkKJY1Qt4bJT1WZ/mlkoYkbVFy2++UFJKObqbOZt3igGO24ZgN/IWk/YoLJG0KnARcFxFrKq+ZWQc44JhtOP4fsBqYWmPZJGAcKSiZjUoOOGYtknSwpP+Q9LykP0haJOmTrZYXEauAG4EpNRZPBV4E7sy2vYekWZJ+LelVSU9IOk/S5nXqu1l2ie3zhfRvSfpNIW0vSVdnl/BWS5or6V2t7psZwGa9roDZhiwiZuTeTyws3gu4F7gMeA3478AsSWsjYna2ziCgYll1zAY+IemAiHgYIAsiJwA/jog3s3y7AC8Dfwu8AkwAzgV2Br7Y5G6uR9LO2X69CEzP9u0sYL6k9/iSnrXKAcesRRExZ917SQLuAt4OfJbWL33NJQWQqcDDWdpRwI75MiNiEbAot/17gVeByySdFhFvtLh9gNOBLYDDIuKVrPz7gOXANODyNsq2jZgvqZm1SNIOkr4v6Sng9ew1HXh3q2VmZw//TjrLUZY8BXgKeCC37U0knS5piaRXs23/H2AMKei143BgHrAquwy3GbAC+Bkw0GbZthFzwDFr3ZWkYPAd4Ejgg8A/A1u2We5s4M+BgyVtCRwHzI71p3Y/HZgJ/BvwUeBA4NRsWbvb3xn4JH8KouteHwb2bLNs24j5kppZC7JAcCzwpYi4LJfeiR9xd5Dun0wFdge2Zfgluo8DcyLinNy2929Q7pvAG8DbCuk7Fj7/F/AIcEGNMlY22IbZiBxwzFqzBbAp8NYNdEnbks422nrIVES8KenfSEFlD2BJRPy8kG1MftuZuj3kIiIkPQe8N1fnTYGPFLLeTjqreswdBKyTHHDMWhARKyQ9BJwjaSWwFjiTdK9juw5sYjbwJVLvtHNqLJ8PfEHSQuBXwF8D40uU++/AdEmPku4LfRbYqpDnH4GTgTskXQo8D+wGHAoMRsS/Nr03ZjjgmLXjZOAK4Crgd8ClpC/vL3Wg7PtJvcLGA3NqLD8X2Il02SuAa4C/A65vUO45pHs0FwB/BL4PLAY+sy5DRPxW0kHAPwCXANsDLwB3AyNOvWPWiKp+xLSkdwJfBQ4C/gK4u8b4hlrrjSUd/MeTOjvcCJwaEb8r5DsO+BbwLtIvv/Mi4upO7oOZmTWvF73U9gMmA09kr7KuBiaSfolNI/UIWu/XnKRDgGtJo7GPAW4CZks6st1Km5lZe3pxhrNJRKzN3l8D7NzoDEfSwcB9wKERcVeWdiDwU+CIiLgtS5sHbB4RH8mtezOwXUQc0o39MTOzcio/w1kXbJp0DPDiumCTlfMg8OtsGdmU7ZOA4g3NOaTxDGNbq7GZmXXCaBn4OQFYWiN9SbYMYB9g8xr5lpD2s+XR32Zm1r7R0kttB9L8UkVDwN65PNTIN1RYvh5J00nTkTBmzJgD9tyzMwOp165dyyabjJZ43jtup3LcTuW4ncrpZDs98cQTL0fELmXyjpaAA7UH06lGevGzRkhPiRFXkLq2MjAwEAsXLmynjm8ZHBxk4sSJHSmrn7mdynE7leN2KqeT7ZTNJVjKaPkpMEQaC1C0PX86oxnKpRXzQO0zJDMzq8hoCThL+dO9mrz8vZ0nSRMMFvNNII0Cb6YLtpmZddhoCThzgd2ycTYASBog3b+ZC29N634naf6pvCnA/RGxoqK6mplZDZXfw5G0FWngJ6SJCbeTdFL2+eaIWC1pGbAgIk4BiIj7szE2V0n6CumMZSZwz7oxOJlvAoOSLiENCp2cvY7u+o6ZmVldveg0sCvpGR556z6/gzR/1GakmXjzpgIXk5438tbUNvkMEXFPFry+BXyBNE7n5Ii4tYP1NzOzFlQecCJiOX/qOTZSnvE10l4BPp296q17PY0nMDSznPFn3tQwz/JvH1tBTayfjZZ7OGZmNso54JiZWSUccMzMrBIOOGZmVgkHHDMzq4QDjpmZVcIBx8zMKuGAY2ZmlXDAMTOzSjjgmJlZJRxwzMysEg44ZmZWCQccMzOrhAOOmZlVwgHHzMwq4YBjZmaVcMAxM7NKOOCYmVklHHDMzKwSDjhmZlYJBxwzM6uEA46ZmVXCAcfMzCrhgGNmZpVwwDEzs0o44JiZWSUccMzMrBIOOGZmVgkHHDMzq4QDjpmZVcIBx8zMKuGAY2ZmlXDAMTOzSjjgmJlZJRxwzMysEg44ZmZWCQccMzOrhAOOmZlVwgHHzMwq4YBjZmaVcMAxM7NKOOCYmVklKg84kvaVdLuk1ZKel3S+pE0brDNDUozw+lou35Uj5JnQ/T0zM7N6NqtyY5J2AG4DFgPHAfsA3yUFvrPrrPpD4JZC2vHAGcDcQvpS4NOFtOWt1djMzDql0oADfB4YA5wYESuB+ZK2A2ZIujBLGyYingWezadJ+gawNCIWFbL/ISIe6ELdzcysDVVfUjsGmFcILHNIQejQsoVI2hE4Apjd2eqZmVm3VB1wJpAueb0lIp4GVmfLyjoJ2JwUrIr2lbRS0hpJ90gqHcjMzKx7FBHVbUx6HfhqRFxSSH8WuCoizipZzh3A2Ig4oJB+GvBH0j2iXYDTgQOAQyLiwRHKmg5MBxg3btwBc+bUimHNW7VqFdtss01Hyupnbqdyut1Ojz23omGe9+0xtmvb7xQfT+V0sp0mTZr0cEQMlMlb9T0cgFoRTiOkD88o7U66/HbGsIIjvlfIexMp+JxF6mQwvDIRVwBXAAwMDMTEiRPLVKOhwcFBOlVWP3M7ldPtdpp25k0N8yz/ZPe23yk+nsrpVTtVfUltCNi+RvpY4JWSZXyCFKCubpQxIl4FbgY+ULaCZmbWHVUHnKUU7tVI2hPYmsK9nTqmAvdExDNNbLe664ZmZlZT1QFnLnCUpG1zaVOAV4EFjVaWNB44iJK90ySNIfWMe7jZipqZWWdVHXAuA9YA10k6PLthPwO4KN9VWtIyST+qsf5U4A3gmuICSWMl3S3pc5IOkzQFuBPYA7igC/tiZmZNqLTTQEQMSToMuBS4gXTf5mJS0CnWq9Z0N1OB2yPipRrL1gAvkWYs2BV4DbgfODQiFnZkB8zMrGWV91KLiMXARxrkGT9C+vvrrPMacGJblTMzs67xbNFmZlYJBxwzM6uEA46ZmVXCAcfMzCrhgGNmZpVwwDEzs0o44JiZWSUccMzMrBIOOGZmVommAo6kWtPNmJmZNdTsGc5zki6U9N6u1MbMzPpWswHncuAk4BeSfippuqTtulAvMzPrM00FnIg4NyL2Bo4AHgcuAl6Q9GNJh3ejgmZm1h9a6jQQEXdExF8DuwFfBt4DzJO0XNIMSX/WyUqamdno124vtQHgw6THRg8BdwOfAZZJ+lSbZZuZWR9pOuBI2kvSuZKeBG4Hdgf+BviziPhfwF6kez3f6WhNzcxsVGvqAWyS7iCd0TwLXAnMioin8nki4k1JPwFO61Qlzcxs9Gv2iZ8vA5OB+RERdfItAt7Rcq3MzKzvNHtJ7VLgvlrBRtI2kj4MEBGvF898zMxs49ZswLkT2HeEZe/JlpuZmQ3TbMBRnWXbAKvbqIuZmfWxhvdwsstkE3NJn5F0dCHblsCxwGOdq5qZmfWTMp0GPkQa3AkQwMeBNwp5/ggsBb7auaqZmVk/aRhwIuI7ZGNqJP0aOCEiFnW7YmZm1l+a6hYdEe7qbGZmLSlzD2cycE9ErMze1xURN3ekZmZm1lfKnOHcCBwEPJi9D0burRaAH9JmZmbDlAk47wBeyL03MzNrWplOA0/Vem9mZtaMMvdwtmqmwIjw4E8zMxumzCW1VaR7M2X5Ho6ZmQ1TJuD8Dc0FHDMzs2HK3MO5soJ6mJlZn2v3EdNmZmallOk08CAwLSIWS3qIBpfXIuLATlXOzMz6R5l7OP8JvJp77/s5ZmbWtDL3cD6dez+tq7UxM7O+1fI9HCW7SKr3UDYzMzOghYAjabKk+4DXgN8Ar0m6T9KxHa+dmZn1jaYCjqTPATeQBoOeRnoY22nZ5//IlpuZmQ3T1PNwgLOAKyLiC4X0yyRdBnwduLwjNTMzs77S7CW1nYDrRlh2LbBjowIk7SvpdkmrJT0v6XxJdafDkTReUtR4zamR9zhJj0l6TdJiSVNK7ZmZmXVVs2c4dwKHAvNrLDsUuKveypJ2AG4DFgPHAfsA3yUFvrNLbP8rwL25zy8Xyj+EFPh+AJwKTAZmSxqKiFtLlG9mZl1SZuDnvrmP3wd+KGkn4Hrgt8CuwAnAMcBnGhT3eWAMcGJErATmS9oOmCHpwiytnscj4oE6y78B3BURp2af75S0H3AO4IBjZtZDZc5wfsH6gz0FfC57FZ/+eQv1Z4s+BphXCCxzgJmkM6QbStSnJklbAJNIZzZ5c4BZksZGxIpWyzczs/aUCTiTOri9CcAd+YSIeFrS6mxZo4AzS9KOpDOr2cDXI2LdLAj7AJsDSwvrLCFdsns38FB71Tczs1aVmWlgQQe3twPwSo30oWzZSNYA/0S6LLYSmAicQQoyx+XKpkb5Q4Xl65E0HZgOMG7cOAYHB+vVv7RVq1Z1rKx+5nYqp9vtdPr73miYZzT8O/l4KqdX7dRsp4G3SNoE2LKYXuKJn7XmYtMI6evKfAH4Ui5pUNKLwA8kvT8iFtUpXyOkryv7CuAKgIGBgZg4cWL92pc0ODhIp8rqZ26ncrrdTtPOvKlhnuWf7N72O8XHUzm9aqdmB35K0hmSlgGvA7+v8apnCNi+RvpYap/51HNN9vcDubKpUf66z82Wb2ZmHdTsOJxTgTOBH5HOHP4BOB94AlhOdmmqjqWkezVvkbQnsDXD7700EoW/T5KC4IRCvgnA2qyOZmbWI80GnM8C5wIXZp+vj4jzgP1IAeNdDdafCxwladtc2hTS4w+avVd0Uvb3YYCIWEMaJ/TxQr4pwP3uoWZm1lvN3sN5B7AoIt6U9DrZ5aqIWCvpB8APSWdAI7mMdJZ0naSZwN7ADOCifFfp7JLdgog4Jfs8A9iWNOhzJfBh4KvAdRHx81z53yTd37mENE5ocvY6usn9NDOzDmv2DOd3wDbZ+6eB/5ZbtgNpUOeIImIIOIw0VucG4DzgYtJZU95mrD+eZylpnM4s4GbgZOA72d98+feQznwOB+YBHwVO9iwDZma91+wZzr3AB0lf+j8hzRCwI/BH4IvA7Y0KiIjFwEca5Blf+DyHNICzoYi4nnR2Y2ZmG5BmA84MYI/s/QWkS2rTSGc284Evd6piZmbWX5oKOBHxOPB49n4N6Vk4p3WhXmZm1mfaGfj5dmB34PmIeK5zVTIzs37UyiOmvyDpGeAp4KfA05KelfS/O147MzPrG83ONHAOcClpPM2xwED2dy7w/Wy5mZnZMM1eUvsicEFEfKOQfks2t9kXSTMPmJmZrafZS2pjGPmpnguoMZmnmZkZNB9wrgdOHGHZx4Ab26uOmZn1qzKPmJ6c+zgXuFDSeIY/Yno/4O87X0UzM+sHZe7h3MjwR0nvARxVI++/kJ7EaWZmtp4yAecdXa+FmZn1vTKPmH6qioqYmVl/a3qmAUmbkToIHALsCPwXcDfpUQGNH4xuZmYbpaYCjqRdgVuB/UlP+HwROJg0/uZRSUdGxEudrqSZmY1+zXaLvgjYCfhQROwdEQdHxN7Ah7L0izpdQTMz6w/NBpzJwBkR8VA+Mfv8NdI0N2ZmZsM0G3C2AH4/wrLfA29rrzpmZtavmg04DwBnSNo6n5h9PiNbbmZmNkyzvdROB+4EnpF0K6nTwK6kQaACJna0dmZm1jeaOsOJiEXAu4ArgF2AI0gB5zLgXRHxaMdraGZmfaH0GY6kzYEDgV9HxJndq5KZmfWjZs5w3gTuAN7bpbqYmVkfKx1wImIt8EtgXPeqY2Zm/arZXmpfB86R9L5uVMbMzPpXs73UzibNKLBI0nOkXmqRzxARB3aobmZm1keaDTi/yF5mZmZNKRVwJI0hTWvzC+A3wG0R8WI3K2ZmZv2lzCOm9wZuA8bnkldK+kRE3NqtipmZWX8p02ngQmAt8FfAVsB+wCPA5V2sl5mZ9ZkyAedg4OyIuDciXouIJcDngD+XtHt3q2dmZv2iTMDZHfhVIe1J0txpu3W8RmZm1pfKjsOJxlnMzMxGVrZb9DxJb9RIv72YHhG7tl8tMzPrN2UCznldr4WZmfW9hgEnIhxwzMysbc3OpWZmZtYSBxwzM6uEA46ZmVXCAcfMzCrhgGNmZpVwwDEzs0pUHnAk7SvpdkmrJT0v6XxJmzZY54OSZklalq33uKRzJW1ZyDdDUtR4Hd3dvTIzs0aafQBbWyTtQHrUwWLgOGAf4LukwHd2nVWnZHlnAr8E9ge+mf39WCHvCqAYYJa0W3czM2tPpQEH+DwwBjgxIlYC8yVtB8yQdGGWVsvMiHgp93lQ0mvA5ZL2ioincsveiIgHulN9MzNrVdWX1I4B5hUCyxxSEDp0pJUKwWadR7K/nrvNzGwUqDrgTACW5hMi4mlgdbasGX9JejDc44X07SW9LOl1SY9IOrHl2pqZWccooronD0h6HfhqRFxSSH8WuCoizipZzm7Az4GbI2JaLv1TpDOeRcA2pAfFTQY+FhHXjVDWdGA6wLhx4w6YM2dOs7tV06pVq9hmm206UlY/czuV0+12euy5FQ3zvG+PsV3bfqf4eCqnk+00adKkhyNioEzeXgScr0TE9wrpzwFXRsTXS5TxNlLHg7cDB0TEUJ28Au4DxkTE+xuVPTAwEAsXLmyUrZTBwUEmTpzYkbL6mdupnG630/gzb2qYZ/m3j+3a9jvFx1M5nWwnSaUDTtWX1IaA7WukjwVeabRyFkCuAvYDJtcLNgCRoul1wP6Nul6bmVl3Vd1LbSmFezWS9gS2pnBvZwQXk7pTHxERZfKv4yeWmpn1WNVnOHOBoyRtm0ubArwKLKi3oqSvAV8GPhUR95TZWHZGdALwaES82VqVzcysE6o+w7kMOBW4TtJMYG9gBnBRvqu0pGXAgog4Jft8MnABcCXwnKSDcmU+ua7btKQFwLWks6Wtgc8CBwHHd3e3zMyskUoDTkQMSToMuBS4gXTf5mJS0CnWK3/P5cjs77TslfdpUiACWAb8LbA7qcv0z4BjI2JuJ+pvZmatq/oMh4hYDHykQZ7xhc/TGB5oaq13ShtVMzOzLvJs0WZmVgkHHDMzq4QDjpmZVcIBx8zMKuGAY2ZmlXDAMTOzSjjgmJlZJRxwzMysEg44ZmZWCQccMzOrhAOOmZlVwgHHzMwq4YBjZmaVcMAxM7NKOOCYmVklHHDMzKwSDjhmZlYJBxwzM6uEA46ZmVXCAcfMzCrhgGNmZpVwwDEzs0o44JiZWSUccMzMrBIOOGZmVgkHHDMzq4QDjpmZVcIBx8zMKuGAY2ZmlXDAMTOzSjjgmJlZJRxwzMysEg44ZmZWCQccMzOrhAOOmZlVwgHHzMwq4YBjZmaVcMAxM7NKOOCYmVklHHDMzKwSDjhmZlYJBxwzM6tE5QFH0r6Sbpe0WtLzks6XtGmJ9cZKmiVpSNIKST+WtFONfMdJekzSa5IWS5rSnT0xM7NmVBpwJO0A3AYEcBxwPnA6cF6J1a8GJgKfAaYBHwSuL5R/CHAtcCdwDHATMFvSkR3ZATMza9lmFW/v88AY4MSIWAnMl7QdMEPShVnaMJIOBo4CDo2Iu7K054CfSjo8Im7Lsn4DuCsiTs0+3ylpP+Ac4Nbu7ZaZNWP8mTc1zLP828dWUBOrUtUB5xhgXiGwzAFmAocCN9RZ78V1wQYgIh6U9Ots2W2StgAmAacW1p0DzJI0NiJWdGg/zLqq+IV8+vveYFouzV/GNhpVHXAmAHfkEyLiaUmrs2UjBZwJwNIa6UuyZQD7AJvXyLeEdOnw3cBDrVXbNmaNfo37y3/D5X+7DUvVAWcH4JUa6UPZslbW2zuXhxr5hgrL1yNpOjA9+7hK0uN16tGMnYGXO1RWPxv17aSZ3d/GqYV2qmKbRVVvs8XtNXU89aIdNxCd/H+3V9mMVQccSB0GijRCeivrFT+rzvpExBXAFQ223TRJCyNioNPl9hu3Uzlup3LcTuX0qp2q7hY9BGxfI30stc9gGq23fW69oVxaMQ8Nyjczsy6rOuAs5U/3XACQtCewNbXv0Yy4XiZ/b+dJ4PUa+SYAa4EnWqivmZl1SNUBZy5wlKRtc2lTgFeBBQ3W2y0bZwOApAHS/Zu5ABGxhjT+5uOFdacA9/egh1rHL9P1KbdTOW6nctxO5fSknRTR6NZJBzeWBn4uBn5B6gq9N3ARcElEnJ3LtwxYEBGn5NJuIfU0+wrpjGUm8NuI+KtcnkOAQeBS0qDQyVn+oyPC43DMzHqo0jOciBgCDgM2JXWBPg+4GDi3kHWzLE/eVNJZ0D8DVwEPAycUyr8HOAk4HJgHfBQ42cHGzKz3Kj3DMTOzjZdni27Ak4021kobSfpg1j7LsvUel3SupC0L+WZIihqvo7u7V53XYjuNH2H/59TIO+qPJWi5nUY6TkLS13L5rhwhT61OSRs0Se+UdLmkRyW9KWmw5Ho9+27qxTicUSM32ehi0mSj+wDfJQXqs+usCmmy0feQJhtdd8/peqB4z+la4AekKXkmkyYbHRotlwHbaKMpWd6ZwC+B/YFvZn8/Vsi7AigGmCXt1r1KbR5LkO5F3pv7vN6gvX44lqCtdvohcEsh7XjgDLKORTlLgU8X0pa3VuOe2o/07/wA8LYm1uvdd1NE+DXCC/gaaXzPdrm0vwdW59NqrHcwaaDph3NpB2Zph+fS5gF3FNa9Gbin1/teQRvtUiNtetZGe+XSZgAv93o/e9hO47M2+R8Nyh/1x1I77TRCWTcBSwppVwILe72fHWqrTXLvrwEGS6zT0+8mX1Krb6TJRseQJhutt96wyUaBdZONkpts9F8L684BDpY0tv3qV6KlNoqIl2okP5L93bVz1dtgtHosNdRHxxJ0qJ0k7QgcAczubPU2HBGxtoXVevrd5IBT37BJQyPiadKvrXrXfDs12eho0Gob1fKXpFP84nx220t6WdLrkh6RdGLLte2ddttpVnad/gVJF0kak1vWL8cSdO54OonUJsPudQH7SlopaY2keyS1FfBHmZ5+Nzng1NeNyUZ3yOWhRr66k41ugFpto/VI2g34OvB/C79ul5EuqXyCdG/neeDaURh0Wm2nNcA/AaeQhhRcDnyB9b9I++VYgg4dT6RhFD+LiOIMI4+QHvr4P4FPkoZfzJd0YAt1HY16+t3kTgONbVCTjW6gWm2jlFF6G+n0fRXwd+sVHPEvhbw3APeRHqp3XSuV7aGm2ykiXgC+lEsalPQi8ANJ74+IRXXKH43HErR/PO1Ouvx2xrCCI75XyHsTqYPCWaROBhuDnn03+QynPk822lirbQSAJJEG8u4HTI40OHhEke5eXgfsX6Z7+gakrXYquCb7+4Fc2dQof7QdS9CZdvoE6cvx6kYZI+JV0s3wDzTK2yd6+t3kgFOfJxttrNU2WudiUvfX4yKiTP51Rtuv9nbbKS8Kf/vlWILOtNNUUm+qZ5rY7mg7nlrV0+8mB5z6NqbJRlvVahuRDcj7MvCpSNMSNZSdEZ0APBoRb7ZW5Z5ouZ1qOCn7+zD01bEEbbaTpPHAQZTsnZZ1vjiGrC03Ar39bup1X/IN+UW6OfYCMJ80P9t00n2GbxXyLQN+VEi7BfgVcCLp2vDjwN2FPIcAbwCXABOBC0m/II7s9b53u42Ak0m/KmeRviDyr11y+RaQBp4dSQo0N2dt9NFe73tF7TSDNPDxxGy980lfvtf227HUTjvl0s8k/TqvNc5rLHA38DlSB4wppEGTa4CBXu97C221FenHx0nA/cB/5j5vNVI79fK7qeeNtqG/gH2BO7L/5C+QRsNvWsizHLiykLZ99mX6CrAS+Amwc43yjyfNnr2GdEo7tdf7XEUbkQbgxQivabl8P8r+c7wK/CH7wjim1/tcYTtNBRaSZlv4Y/YFcj6wRT8eS622Uy59EXDLCOVuSbr/90zWRiuyL9+Der3PLbbT+Dr/h8aP1E69/G7y5J1mZlYJ38MxM7NKOOCYmVklHHDMzKwSDjhmZlYJBxwzM6uEA46ZmVXCAcfMzCrhgGNmZpX4/3wikDlpOrGeAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZwAAAEPCAYAAAB2s3LUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XmYXFWd//H3h0WIBMKagAwSQTHDMqKEzWGkIyBbMIgg62gUDCqy+CATQISwTJyAbBPGgfwcCQ5CcIBhhgAJa4LgBkiYIBAMGlYTWQIhJuzf3x/nNhTVVdW3qqtudVc+r+epp7vOPefW91RV17fvueeeUkRgZmbWaiu1OwAzM1sxOOGYmVkhnHDMzKwQTjhmZlYIJxwzMyuEE46ZmRXCCceaStIESVHldkTOfWyR7WftsvKx2X4Gtyb6fHH0cZ/XSpqVo94qkk6Q9JCk5ZIWS7pZ0i4NPm5/eU7Hlr0n/iJppqRP5WjblbXZuohYrfmccKwVXgF2rnCbkbP9FsAZQPkH/U3ZfpY1J8yG42gpSSsDNwATgf8F9gHGAm8DsyQd1sBu+8tz2u2z2eMeDWwA3CXpQ720+V3W5okWx2Ytskq7A7CO9FZE/LrZO42I54Hnm73ffuhYYF9g74goTdL/I2kaMEXS7Ih4tq8P1Mbn9L6IWAog6X7gSeBw4LzyipIErBYRS4Cmv6+sOD7CsbaQdIqk+ZJek7RI0gxJG0rqAm7Mqv0pG0JZkLV53/CPpOHZ/UMkXS5piaRnuofuJP2TpOckPS9pkqSVSh5/hKRpkp6WtEzS77MhrJWy7VXjyLZ/OGv/UtZ+pqSPl/Vxk2wYbLmkBZKOyvn0HA/cVZZsun0PWB04suRxFkj6oaTvS1ooaamkn0ka0ltfKg2pSVpf0hWSXsz6NkvSyLK+dT/md7LnfHH2fNR9NBgRT5OS3vBs3xMkvSBpF0n3Aa8BB1UaUpO0cvZeelzS61ksU8tiHSPp/uy9tlDSuZJWrTdO6zsf4VhLSOrx3oqIt7JtXwZOBcYDvwfWIw2xrEEaNvku8EPgAODPwOu9PNwk4GfAF4GvAVdI+iSwaXZ/O+Ac4EFgWtZmY2Be1u5VYFvgTGAQ8INacUhaF7gHeBH4Bmk46mTgdklbRMTy7L/y/wHWJyWH17L9rwv8ocbztgnpg/fCStsj4glJc4HPlG06FJgPfB3YCDgX+DFwUK2+VHED8NGszQvASaQhr09GxPySel8C/g8YB/wNcAFpGPBbNfbdg6Q1Sc/LwpLiDwJXZP14HHgu61e5y4AvZ/VmZ/s5sGTfXwKuzuqdCmxOen1XyvpnRYoI33xr2g2YAESV2/CsziXAdTX2Mbq0fkn52Kx8cHZ/eHb/8pI6awFvkj7UVy4p/y1wTZXHE+mfr1OBP+aI42xSslm3pGwd0rmrY7L7+2RtdyypsynwFjCrRt93ytqNqVHnBuDRkvsLgJe6n5es7HDgHeBv63xO98ru71pSZw3SEchlZY/5BLBKSdlFwMJe3h/djzcke843Aa7Jnpdty95DY8radmXlW2f3R2T3j6vxuj5Z+v7Iyr8GLAfWa/ffy4p28xGOtcIrwO4Vyp/Lfs4BjpR0Jumk9QMR8XYfHu+O7l8iYomk54HZZfucD3y4+46k1YFTSB/MHwZWLdm2SmRHY1XsDtwGLCk5knsVeADoHnraAVgUEb8pie1JSQ800L88bovsnEjmeuBKYHvg0Tr2swPwfETM7i6IiL9Kmg6Uz5C7q+x5egQYKukDEfFGL4/zcsnvLwBfi4g5JWUB3NLLPkZlP6dW2b4F6bX9edkR952kYcmtSUdFVhAnHGuFtyLi/hrbfwKsSRqKOR14UdK/AxMaTDwvl91/o0rZ6iX3JwFHkYa5fpfVHwOcltVbSnXrk45EDq6wrTv5bQj8pcL2v5D6Xk33RIBNa9TZtKRe6X7fFWlYbymVh6Fq2QhYVKF8EWm4qlSl51jAB7Lfa/kMaSjyBeDpiHinbPviHElrPeCvkSYTVLJ+9vPmKts36WX/1mROOFa47MPlQuDC7JzF4cA/kz5ELy0ojIOAyRFxbneBpH1ztn2JNF357ArbXs1+LgSGVtg+lDScU1FEPJ2d0P888K/l2yV9hPSfefljDy2rNwgYTDpfU48/l+8rM4zU72Z5sOyIrFye7015EVhD0lpVkk53vONI5+/K/SnHY1gTeZaatVVEPB0R/0Ia8toyK+7+z3b1yq2aYhAlJ86Vrn05pKxOtTjuALYCfh8R95fd5mV17gOGSdqx5DE+DPR6gSNwMbCbpM9V2HZOFvd/lJXvofdfvHkA6UO7+0gz73P6G9Kw2LuTEiR9kDRN+54csRfpzuznl6tsn0f6J2Z4hdfp/oh4sZgwrZuPcKwVVpG0U4XypyPiWUmXkf77/DXpfM8o4GOkWWuQPigAjla67mRZRMxtcoy3AcdImp/FcgywWlmdanFcABwB3ClpMulDbRiwK3BPRFxNGsZ5CPgvSeNJs9TOovIwW7nJpPNE/y3ph8As0jDckaST//8YPa/BWQ7cJOk80rDYecB/R8QjvfTlfSJipqR7gWsknUw6ivguKUH3uEamnSJinqQpwPmShgJ3ky5sPTAiDomIdySdCPynpLVI54TeADYD9s/qFX3B64qt3bMWfOusG7VnqZ2W1RkL3Ev6oF9Gmlp7ZNl+TiTNMHoLWFDSrtIstdFlbRcAPywrmwrcX3J/GPDfwBLS+YlzSVOK391/tTiy8g8Bl2dtX88e80pgq5I6HyatrrA828fRwLXUmKVW0nYV4DvZc7McWEz6wNylQt0FwPnZc78I+CtpKvDa9T6nWdkGwE+zx1xOOrG+fY7nuMe+KsSap84E4IUK5V2UzFLLylYmm11ISibP0HNW2t7AL7LnZQlp0so5lMyw862Ym7IXpDCSPkqa178TaSz6FxHRlaPdENK0y/1JQ4HTSdMhXyyrN4b0ZvoY6U14ZkRc08w+mPUn2TmfayPC15VYv9aOczhbka5ReDy75XUN6T+co0j/JW1Puh7hXUoLG14H3EX6r+Ym4OoqY+FmZlagdhzhrBTZFEhJ1wLr93aEI2ln4Jeki9Huzsp2IJ3g3CMibs/KZgKrRsRnS9reDKwVEQ2tsmvW3/kIxwaKwo9woud8+zz2Jl1Ed3fJfn5Lmta4N4Ck1Ugnn39e1nYasHP3ulJmnSYihjvZ2EAwUKZFjwAeq1D+aLYN0hpJq1ao9yipn1u0LDozM+vVQJkWvQ49r2qGNItms5I6VKi3uGz7+0gaR7owjEGDBm23ySb5Lj5+5513WGmlgZKvK+uEPoD70WprPp5Otb66Re//s/XXPtTL/cjv8ccffyEiNshTd6AkHKh85bEqlJffV432RMQUYArAyJEj4/77a63I8p5Zs2bR1dWVq25/1Ql9APej5ZT9Cc2bV7se/bgPdXI/8pP0ZN66AyWFL6byty6uzXtHNItLysrrQOUjJDMzK8hASTiP8d65mlKl53aeIC1LX15vBGmZ9nqmYJuZWZMNlIRzC7Bhdp0NANk3EG6WbSMiXiddf3NQWduDgV9FxCsFxWpmZhUUfg4nWwhwn+zuxsBakrq/oe/miFiWrW81OyKOBIiIX2XX2PxU0ndJRyyTSOtW3V6y+7OBWZIuIl0Uuk9226vlHTMzs5raMWlgKPBfZWXd9z9CWqNpFdIaSaUOIS1p/xNKlrYprRAR92TJ6xzgm6TrdA6LiFubGL/ZiqXgi8OtcxWecCJiAe/NHKtWZ3iFspeBr2a3Wm1voGzJGzPrafjJNzXcdsG/5P3qILP3DJRzOGZmNsA54ZhZTTdOPZ4bpx7f7jCsAwykCz/NrA22WfREu0OwDuEjHDMzK4QTjpmZFcIJx8zMCuGEY2ZmhXDCMTOzQniWmpnVdNUn9mx3CNYhnHDMrKZT9zq23SFYh/CQmpmZFcIJx8xq2nrhfLZeOL/dYVgH8JCamdU0/YoTABg+fnqbI7GBzkc4ZmZWCCccMzMrhBOOmZkVwgnHzMwK4YRjZmaFcMIxM7NCeFq0mdU0+isXtTsE6xBOOGZW08MbfrTdIViH8JCamZkVwgnHzGqaOGMyE2dMbncY1gGccMyspsMemslhD81sdxjWAZxwzMysEE44ZmZWCCccMzMrhBOOmZkVwgnHzMwK4Qs/zaymucM2b3cI1iGccMyspv3GXtzuEKxDeEjNzMwK4YRjZmaFcMIxs5oWTBrNgkmj2x2GdQAnHDMzK4QTjpmZFcIJx8zMCuGEY2ZmhXDCMTOzQjjhmJlZIbzSgJnVdMqe3253CNYhnHDMrKart92r3SFYhyh8SE3SlpLukLRM0nOSzpK0ci9tJkiKKrdTSupNrVJnROt7ZmZmtRR6hCNpHeB24BFgDLA5cD4p8Z1Wo+mPgRllZfsD44FbysofA75aVragsYjN7NA56U/PRzrWV0UPqX0DGAQcEBFLgNskrQVMkHRuVtZDRDwDPFNaJun7wGMRMaes+l8j4tctiN1shfSDmZcATjjWd0UPqe0NzCxLLNNISWjXvDuRtC6wB3B1c8MzM7NWKTrhjCANeb0rIp4ClmXb8joQWJWUrMptKWmJpNcl3SMpdyIzM7PWUUQU92DSm8BJEXFRWfkzwE8j4tSc+7kTGBIR25WVHw+8QTpHtAFwIrAdsEtE/LbKvsYB4wCGDRu23bRplXJYT0uXLmXw4MG56vZXndAHcD8aNffZV3LVO/aI/QGYfOUN75Zts/GQinX9WvQvRfRj1KhRD0TEyDx12zEtulKGU5XynhWljUjDb+N77Dji4rK6N5GSz6mkSQY9g4mYAkwBGDlyZHR1deUJg1mzZpG3bn/VCX0A96NRY0++KVe9Y7Of58997+NiweFdFev6tehf+ls/ih5SWwysXaF8CPByzn18iZSgrumtYkQsB24GPpU3QDMza42iE85jlJ2rkbQJsAZl53ZqOAS4JyKeruNxixs3NDOziooeUrsFOEnSmhHxalZ2MLAcmN1bY0nDgZ2Ab+V5MEmDSDPjHmgkWDOD4eOntzsE6xBFH+FcCrwOXC9p9+yE/QTggtKp0pLmS/qPCu0PAd4Cri3fIGmIpF9IOlrSbpIOBu4CNgYmtqAvZmZWh0KPcCJisaTdgEuAG0nnbS4kJZ3yuCotd3MIcEdEPF9h2+vA86QVC4YCrwG/AnaNiPub0gEzM2tY4bPUIuIR4LO91BlepXzbGm1eAw7oU3Bm1sONU48HYL+xF/dS06w2rxZtZjVts+iJdodgHcJfwGZmZoVwwjEzs0I44ZiZWSGccMzMrBBOOGZmVgjPUjOzmq76xJ7tDsE6hBOOmdV06l7H9l7JLAcPqZmZWSHqSjiSKi03Y2YdbOuF89l64fx2h2EdoN4jnGclnSvpb1sSjZn1O9OvOIHpV5zQ7jCsA9SbcC4DDgQelvQbSeMkrdWCuMzMrMPUlXAi4oyI2AzYA5gHXAD8WdLPJO3eigDNzKwzNDRpICLujIgvAxuSvvL848BMSQskTZD0oWYGaWZmA19fZ6mNBD5D+troxcAvgKOA+ZKO6OO+zcysg9SdcCRtKukMSU8AdwAbAV8DPhQR/whsSjrXc15TIzUzswGtrgs/Jd1JOqJ5BpgKXB4RT5bWiYi3JV0FHN+sIM3MbOCrd6WBF4B9gNsiImrUmwN8pOGozKzfGP2Vi9odgnWIehPOJcDvKiUbSYOBT0XE3RHxJvBkj9ZmNuA8vOFH2x2CdYh6z+HcBWxZZdvHs+1mZmY91JtwVGPbYGBZH2Ixs35o4ozJTJwxud1hWAfodUhN0meArpKioyTtVVZtdWBfYG7zQjOz/uCwh2YCXjXa+i7POZwdSRd3AgRwEPBWWZ03gMeAk5oXmpmZdZJeE05EnEd2TY2kPwFfiIg5rQ7MzMw6S12z1CLCU53NzKwhec7h7APcExFLst9rioibmxKZmZl1lDxHONOBnYDfZr8H1WerBeAvaTMzsx7yJJyPAH8u+d3MViBzh23e7hCsQ+SZNPBkpd/NbMWw39iL2x2CdYg853A+WM8OI8IXf5qZWQ95htSWks7N5OVzOGZm1kOehPM16ks4ZtZBFkwaDcDw8dPbHIkNdHnO4UwtIA4zM+twff2KaTMzs1zyTBr4LTA2Ih6RdB+9DK9FxA7NCs7MzDpHnnM4vweWl/zu8zlmZla3POdwvlry+9iWRmNmZh2r4XM4SjaQVOtL2czMzIA6V4uGdxfzPA3YLmv/lqQHgH+OiJuaHJ+Ztdkpe3673SFYh6gr4Ug6GvgRcAdwPPAXYChwAPC/kr4VEZc1PUoza5urty3/gl+zxtR7hHMqMCUivllWfqmkS4HvAU44ZmbWQ73ncNYDrq+y7Tpg3d52IGlLSXdIWibpOUlnSaq5HI6k4ZKiwm1ahbpjJM2V9JqkRyQdnKtnZlbRoXNmcOicGe0OwzpAvUc4dwG7ArdV2LYrcHetxpLWAW4HHgHGAJsD55MS32k5Hv+7wL0l918o2/8upMT3I+A4YB/gakmLI+LWHPs3szI/mHkJ4KE167s8F35uWXL3X4EfS1oPuIH3zuF8AdgbOKqX3X0DGAQcEBFLgNskrQVMkHRuVlbLvIj4dY3t3wfujojjsvt3SdoKOB1wwjEza6M8RzgP8/6LPQUcnd3Kv/1zBrVXi94bmFmWWKYBk0hHSDfmiKciSasBo0hHNqWmAZdLGhIRrzS6fzMz65s8CWdUEx9vBHBnaUFEPCVpWbatt4RzuaR1SUdWVwPfi4juVRA2B1YFHitr8yhpyG4L4L6+hW9mZo3Ks9LA7CY+3jrAyxXKF2fbqnkd+DfSsNgSoAsYT0oyY0r2TYX9Ly7b/j6SxgHjAIYNG8asWbNqxf+upUuX5q7bX3VCH8D9aNSJ27zVcP1qcfq16F/6Wz/qvvCzm6SVgNXLy3N842eltdhUpbx7n38GSq8+myVpEfAjSdtGxJwa+1eV8u59TwGmAIwcOTK6urpqR98dwKxZ5K3bX3VCH8D9aNTYk/Ndp31s9vP8ue99XCw4vKtiXb8W/Ut/60dd06Kz5WzGS5oPvAm8WuFWy2Jg7QrlQ6h85FPLtdnPT5Xsmwr7775f7/7NzKyJ6j3COQ44GTgX+GfgHOBt4BDgA8DEXto/RjpX8y5JmwBr0PPcS2+i7OcTpCQ4AigdBhwBvAM8Xuf+zQx/06c1T70Xfn4dOIOUcABuiIgzga1ICeNjvbS/BdhT0polZQeTvv6g3nNFB2Y/HwCIiNdJ1wkdVFbvYOBXnqFmZtZe9R7hfASYExFvS3qTbLgqIt6R9CPgx6QjoGouJR0lXS9pErAZMAG4oHSqdDZkNzsijszuTwDWJF30uQT4DHAScH1E/F/J/s8mnd+5iHSd0D7ZzVesmZm1Wb1HOC8Cg7PfnwI+WbJtHdJFnVVFxGJgN9K1OjcCZwIXko6aSq3C+6/neYx0nc7lwM3AYcB52c/S/d9DOvLZHZgJfB44zKsMmDXuxqnHc+PU49sdhnWAeo9w7gW2J33oX0VaIWBd4A3gGNIq0jVFxCPAZ3upM7zs/jTSBZy9iogbSEc3ZtYE2yx6ot0hWIeoN+FMADbOfp9IGlIbSzqyuY33ZlCamZm9T10JJyLmAfOy318nfSeOj7XNzKxXfbnw82+AjYDnIuLZ5oVkZmadqN5JA0j6pqSngSeB3wBPSXpG0reaHp2ZmXWMelcaOB24hHQ9zb7AyOznLcC/ZtvNzMx6qHdI7RhgYkR8v6x8Rra22THAWU2JzMz6has+sWe7Q7AOUW/CGUT1b/WcjWepmXWcU/fyn7U1R73ncG4ADqiy7YuAF10yM7OK8nzF9D4ld28BzpU0nJ5fMb0V8E/ND9HM2mnrhfMBeHjDj7Y5Ehvo8gypTafnV0lvDFQa2L2S9E2cZtYhpl9xAuBVo63v8iScj7Q8CjMz63h5vmL6ySICMTOzzlb3SgOSViFNENgFWBd4CfgF6asC6vuSdDMzW2HUlXAkDQVuBf4OWAAsAnYmXX/zkKTPRcTzzQ7SzMwGvnqnRV8ArAfsGBGbRcTOEbEZsGNWfkGzAzQzs85Qb8LZBxgfEfeVFmb3TyEtc2NmZtZDvedwVgNerbLtVeADfQvHzPqb0V+5qN0hWIeoN+H8Ghgv6c6I+Gt3oaQ1gPHZdjPrIL7g05ql3oRzInAX8LSkW0mTBoaSLgIV0NXU6MzMrGPUdQ4nIuYAHwOmABsAe5ASzqXAxyLioaZHaGZtNXHGZCbOmNzuMKwD5D7CkbQqsAPwp4g4uXUhmVl/cthDMwGvGm19V88RztvAncDftigWMzPrYLkTTkS8A/wBGNa6cMzMrFPVex3O94DTJW3TimDMzKxz1TtL7TTSigJzJD1LmqUWpRUiYocmxWZmZh2k3oTzcHYzMzOrS66EI2kQaVmbh4GFwO0RsaiVgZlZ/zB32ObtDsE6RJ6vmN4MuB0YXlK8RNKXIuLWVgVmZv3DfmMvbncI1iHyTBo4F3gH+Afgg8BWwIPAZS2My8zMOkyehLMzcFpE3BsRr0XEo8DRwIclbdTa8MzMrFPkSTgbAX8sK3uCtHbahk2PyMz6lQWTRrNg0uh2h2EdIO91ONF7FTMzs+ryToueKemtCuV3lJdHxNC+h2VmZp0mT8I5s+VRmJlZx+s14USEE46ZmfVZvWupmZmZNcQJx8zMClHvWmpmtoI5Zc9vtzsE6xBOOGZW09Xb7tXuEKxDeEjNzMwK4YRjZjUdOmcGh86Z0e4wrAMUnnAkbSnpDknLJD0n6SxJK/fSZntJl0uan7WbJ+kMSauX1ZsgKSrcPCZg1qAfzLyEH8y8pN1hWAco9ByOpHVIX3XwCDAG2Bw4n5T4TqvR9OCs7iTgD8DfAWdnP79YVvcVoDzBPNrX2M3MrG+KnjTwDWAQcEBELAFuk7QWMEHSuVlZJZMi4vmS+7MkvQZcJmnTiHiyZNtbEfHr1oRvZmaNKnpIbW9gZllimUZKQrtWa1SWbLo9mP302m1mZgNA0QlnBPBYaUFEPAUsy7bV49OkL4abV1a+tqQXJL0p6UFJBzQcrZmZNY0iivvmAUlvAidFxEVl5c8AP42IU3PuZ0Pg/4CbI2JsSfkRpCOeOcBg0hfF7QN8MSKur7KvccA4gGHDhm03bdq0XH1ZunQpgwcPzlW3v+qEPoD70ai5z76Sq96xR+wPwOQrb3i3bJuNh1Ss69eifymiH6NGjXogIkbmqduOCz8rZThVKe9ZUfoA8HNgKfCd9+044sqyujcCvwROByomnIiYAkwBGDlyZHR1deUJg1mzZpG3bn/VCX0A96NRY0++KVe9Y7Of58997+NiweFdFev6tehf+ls/ik44i4G1K5QPAV7urbEkAT8FtgL+PiIW16ofESHpemCSpJUj4u0GYjZboQ0fP73dIViHKDrhPEbZuRpJmwBrUHZup4oLSdOp94iIPPW7+RtLzczarOhJA7cAe0pas6TsYGA5MLtWQ0mnkI7uj4iIe/I8WHZE9AXgIR/dmJm1V9FHOJcCxwHXS5oEbAZMAC4onSotaT4wOyKOzO4fBkwEpgLPStqpZJ9PdE+bljQbuI50tLQG8HVgJ2D/1nbLrHPdOPV4APYbe3GbI7GBrtCEExGLJe0GXALcSDpvcyEp6ZTHVbrczeeyn2OzW6mvkhIRwHzgBGAj0pTp3wH7RsQtzYjfbEW0zaIn2h2CdYjCZ6lFxCPAZ3upM7zs/lh6JppK7Y7sQ2hmZtZCXi3azMwK4YRjZmaFcMIxM7NCOOGYmVkh2rG0jZkNIFd9Ys92h2AdwgnHzGo6da9je69kloOH1MzMrBBOOGZW09YL57P1wvntDsM6gIfUzKym6VecAHjVaOs7H+GYmVkhnHDMzKwQTjhmZlYIJxwzMyuEE46ZmRXCCcfMzArhadFmVtPor1zU7hCsQzjhmFlND2/40XaHYB3CQ2pmZlYIJxwzq2nijMlMnDG53WFYB3DCMbOaDntoJoc9NLPdYVgHcMIxM7NCOOGYmVkhnHDMzKwQTjhmZlYIJxwzMyuEL/w0s5rmDtu83SFYh3DCMbOa9ht7cbtDsA7hITUzMyuEE46ZmRXCCcfMalowaTQLJo1udxjWAZxwzMysEE44ZmZWCCccMzMrhBOOmZkVwgnHzMwK4YRjZmaF8EoDZlbTKXt+u90hWIdwwjGzmq7edq92h2AdwkNqZmZWCCccM6vp0DkzOHTOjHaHYR3AQ2pmVtMPZl4CeGjN+s5HOGZmVojCE46kLSXdIWmZpOcknSVp5Rzthki6XNJiSa9I+pmk9SrUGyNprqTXJD0i6eDW9MTMzOpRaMKRtA5wOxDAGOAs4ETgzBzNrwG6gKOAscD2wA1l+98FuA64C9gbuAm4WtLnmtIBMzNrWNHncL4BDAIOiIglwG2S1gImSDo3K+tB0s7AnsCuEXF3VvYs8BtJu0fE7VnV7wN3R8Rx2f27JG0FnA7c2rpumVktw0++qeG2C/5l3yZGYu1UdMLZG5hZllimAZOAXYEba7Rb1J1sACLit5L+lG27XdJqwCjguLK204DLJQ2JiFea1A+zppv77CuMbeCD2R/INlAUnXBGAHeWFkTEU5KWZduqJZwRwGMVyh/NtgFsDqxaod6jpKHDLYD7GgvbVjT+j3zg8WvW/xWdcNYBXq5Qvjjb1ki7zUrqUKHe4rLt7yNpHDAuu7tU0rwacZRaH3ghZ93+qhP6AP2sH5rUcNOG+tGHx8u3/+5fSr71s8ZjNv21aHX/qjxGv3pP9UER/dg0b8V2XIcTFcpUpbyRduX3VaU8FUZMAab08tg9H1i6PyJG1tuuP+mEPoD70Z90Qh/A/WiVoqdFLwbWrlA+hMpHML21W7uk3eKSsvI69LJ/MzNrsaITzmO8d84FAEmbAGtQ+RxN1XaZ0nM7TwBvVqg3AngHeLyBeM3MrEmKTji3AHtKWrOk7GBgOTC7l3YbZtfZACBpJOn8zS0AEfE66fqbg8raHgz8qgUz1OoehuuHOqEP4H70J53QB3A/WkIRvZ06aeKDpQs/HwEeJk2F3gy4ALgoIk4rqTcfmB0RR5aUzSDNNPusb93LAAAGwUlEQVQu6YhlEvCXiPiHkjq7ALOAS0gXhe6T1d8rInwdjplZGxV6hBMRi4HdgJVJU6DPBC4EziirukpWp9QhpKOgnwA/BR4AvlC2/3uAA4HdgZnA54HDnGzMzNqv0CMcMzNbcXm16JwkfV3SH7JFQR+QtFud7T8p6W1JbZ3b30g/JB0t6TZJi7KFU+8tYn26Vi/0WpRG+iFp+6wP87N28ySdIWn1ouKuEFNDr0dJ+5Wy91xIGt17i+brSx8kHSDpPknLJb0oaYakNVodc5VYGv3bGCnp1iz+lyTdLmnHImIGICJ86+VGGs57m7RW2yjSkN5yYOuc7QXcCywEXhho/QCeIp183B/YA7iCdB7t8y2MdR3gOdJir3uQ1uH7K3BOjrYzgD8BXyQNuz4O/KJNz3lD/QB+CNwNfJ20aO1xwCvAdQOpH2X7GJf9DQQweiD1gbRo8GukBYe7svfVZGDIQOkHsAnp8pA7gX2z26zsfbVpIbEX/WQNxBswD/hJyf2VgLnAlTnb/yMwH5jY5oTTUD+A9SuU/RK4q4WxnkK6tmqtkrJ/ApaVllVot3P2gfaZkrIdsrLd2/CcN9qPDSqUjcv6selA6UdJ3XWA54Ej25hwGn0t1gdeBb5edMxN7sc3SP9wrl32urwNfLOI2D2k1gtJm5Fmx/28uywi3gH+i7RwaG/t1yTNqPsu8EaLwuxVX/oREZWGAR8EhjYzxjLVFnodRFrotVa7Hgu9ko54en29WqChfkTE8xWKH8x+tvJ5r6bR16Pb2aSj/DtaEFtejfbhS9nPK1oVWJ0a7ceqwFvA0pKypVmZKrZoMiec3nVfSFppUdB1JW3QS/vTgUcj4oZe6rVaX/tRbmfSFPdW6bFga0Q8RfovrtJFwFXbZUoXei1So/2o5NOkocy86/01U8P9kPR3wFdJ/3S1U6N92JH0nB8p6RlJb0r6jaRPty7Umhrtx3VZnfMlDZU0lDRLeDHpH8+Wc8LpXUOLggJI+jhwDHBCC+KqV8P9KCfpa8AngR81Ia5qWrHQa+4+NlFT4pG0IfA94D+jyvdGtVhf+jEZ+LeImN/0qOrTaB82BD4OnAaMB/YjnTOZIWlYs4PMoaF+RMRzpHO3XwQWZbcDgD2rHFE3XTsW72w7SUOAjXqrFxGl/0XUtSho5mJgakTMrS/CfArsR+ljbkf6ALk4Iu7K06YPWr3Qa1H6FI+kD5CGQpcC32liXPWqux+SDiF9WO/XqqDq1MhrsRIwGDgoImYASPol8CTwbdIknKI18lpsBFxLuobxqKz4GOAmSZ/OjpJaaoVMOKTlb/5fjnri/YuCli6PU3NRUEl7A38PfFtSd93V0yatDSyPtBxPX7S8H+/bSToPdBNpHP7E/GE2pC8LvVYaHixd6LVIjfYDSG8W0mzCrYC/j3TxdDvU3Q9JqwLnkc5hrpS979fKNq8hac2IeLUVwVbR6GvxUvZzVndBRCyR9ACwZdOiy6/RfpxE+sw/MCLeBJB0J/AH0nBn+ZdXNt0KOaQWET+OCPV2y6p3Hx1UWhT0pRqHoh8n/Vf0B9IbZDHpcHzd7PeTBkg/AMjGe2eS/qs7JCLe7mv8vWjlQq9FarQf3S4ExgBjyo5Ui9ZIP9YA/oa0fFX338BD2bZpvDcJoiiNvhaPko4cyk+si3ROrWiN9mME8PvuZAMQEW8Avyd9gWXLrZAJpx4R8UfSdRzvLgoqaaXs/i01ml5LGi8tvV0BLMl+/88WhVxRH/qBpMHAzdnd0RGxrFVxlmjZQq8Fa7QfSDoFOBY4ItKyTe3USD+W0vNv4NBs26nA4a0JtapGX4vppOQyqrsgG87ejvcSaJEa7ceTwNbZEC0AklYDtgYWtCDOnto5n3yg3Eh/JG+TThqOAqZSdsEkaTriW8CuNfYzgfZeh9NQP4BbSVO6DwN2Kr21MNZ1gD8Dt5HWxhtH+gA7p6zefOA/yspmAH8knRDdnzTDqJ0Xftbdj+y5DuDy8uecCtfo9Nd+VNjPcNp74Wej76kbsrZfIV0wOZt0XdE6A6UfpAT5JmlYfF9gNCl5vQl8opDYi36yBuqNdMX3fOB14HfAbmXbu7I/pK4a+5hAGxNOo/3I7le8tTjWLUlXRS/P/sDOBlYuq7OANDGjtGzt7IP6ZdIR5VVUuHi1wOe87n6Q/hmo9ryPHSj9qLCPtiWcPr6nBgP/DryYtb0d2GYgvaeyst1IK1i8lN1m1/rMavbNi3eamVkhfA7HzMwK4YRjZmaFcMIxM7NCOOGYmVkhnHDMzKwQTjhmZlYIJxwzMyuEE46ZmRXi/wPqkbRjapwSQAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot estimated values for \"a\"\n",
    "plt.bar(result['values'], result['probabilities'], width=0.5/len(result['probabilities']))\n",
    "plt.xticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.yticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.title('\"a\" Value', size=15)\n",
    "plt.ylabel('Probability', size=15)\n",
    "plt.ylim((0,1))\n",
    "plt.grid()\n",
    "plt.show()\n",
    "\n",
    "# plot estimated values for option price (after re-scaling and reversing the c_approx-transformation)\n",
    "plt.bar(result['mapped_values'], result['probabilities'], width=1/len(result['probabilities']))\n",
    "plt.plot([exact_value, exact_value], [0,1], 'r--', linewidth=2)\n",
    "plt.xticks(size=15)\n",
    "plt.yticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.title('Estimated Option Price', size=15)\n",
    "plt.ylabel('Probability', size=15)\n",
    "plt.ylim((0,1))\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}