notebooks/matching_gaussian_prior.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "GJcMJWZDBy2h"
      },
      "source": [
        "# Matching q=Normal(loc, softplus(uscale)) to p=Normal(0,1) by minimizing KL(q|p)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 116,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "uJd1yRdoJgA7",
        "outputId": "ac8f8121-db81-47a0-f47e-eb001f0b38c8"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Requirement already satisfied: torch in /Users/tomasz.kusmierczyk/anaconda3/lib/python3.8/site-packages (1.9.0)\n",
            "Requirement already satisfied: typing-extensions in /Users/tomasz.kusmierczyk/anaconda3/lib/python3.8/site-packages (from torch) (4.4.0)\n",
            "Note: you may need to restart the kernel to use updated packages.\n"
          ]
        }
      ],
      "source": [
        "%load_ext autoreload\n",
        "%autoreload 2\n",
        "%pip install -q torch"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 117,
      "metadata": {
        "id": "PkWZ_bKBJYeT"
      },
      "outputs": [],
      "source": [
        "import torch\n",
        "import numpy as np"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 118,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "uiEnOgwGJaEn",
        "outputId": "c8936bc9-58b4-4ecf-fcea-88d08592c3c4"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Using CPU\n"
          ]
        }
      ],
      "source": [
        "if torch.cuda.is_available():\n",
        "    torch.set_default_tensor_type('torch.cuda.FloatTensor')\n",
        "    ENV = torch.cuda\n",
        "    print(\"Using GPU\")\n",
        "else:\n",
        "    torch.set_default_tensor_type('torch.FloatTensor')\n",
        "    ENV = torch\n",
        "    print(\"Using CPU\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 119,
      "metadata": {
        "id": "JVdgdm7Cd5UI"
      },
      "outputs": [],
      "source": [
        "from matplotlib import pyplot as plt\n",
        "plt.rcParams['figure.figsize'] = (5, 5) # Width and height\n",
        "\n",
        "import seaborn as sns; \n",
        "sns.set(); \n",
        "sns.set(font_scale=2); \n",
        "sns.set_style(\"white\"); "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 120,
      "metadata": {},
      "outputs": [],
      "source": [
        "import sys\n",
        "sys.path.append(\"../\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 121,
      "metadata": {},
      "outputs": [],
      "source": [
        "import reparameterized as r\n",
        "from reparameterized import sampling\n",
        "from reparameterized import densities"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Create q"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 122,
      "metadata": {
        "id": "18IvegmjJYea"
      },
      "outputs": [],
      "source": [
        "parameter_shape = (3,)\n",
        "parameter_init_value = torch.randn(torch.Size(parameter_shape))\n",
        "\n",
        "sampler, variational_params, aux_objs = sampling.create_factorized_gaussian_sampler(parameter_init_value,\n",
        "    loc_initalization=lambda parameter_init_value: parameter_init_value.clone().detach()+20,\n",
        "    uscale_initialization=lambda parameter_init_value: torch.randn_like(parameter_init_value)*20,\n",
        ")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 123,
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "{'loc': tensor([18.6674, 18.8212, 21.1151], requires_grad=True),\n",
              " 'unnormalized_scale': tensor([  0.7450,   1.3804, -22.8532], requires_grad=True)}"
            ]
          },
          "execution_count": 123,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "variational_params"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 124,
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "(torch.Size([10240, 3]), torch.Size([10240]))"
            ]
          },
          "execution_count": 124,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "samples, nlls = sampler(n_samples=10240)\n",
        "samples.shape, nlls.shape"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 125,
      "metadata": {},
      "outputs": [
        {
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWUAAAFTCAYAAADhph0LAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3X1cFWX+//EXCIKAiaiomBYiYGqYkqj9MsxMdy1Xt0w2b1KTr9/W1NYyLG9Wkc1NuyO0stayvMkbyHRJE9M1MzI1zHS9QQE1M+8QBREBgfn94fecRA5yUI4O+X4+Hj3Ynetz5swM4/tc55prBifDMAxERMQUnG/2BoiIyG8UyiIiJqJQFhExEYWyiIiJKJRFRExEoSwiYiIK5Wrql19+ITg42OZ/bdq0ISwsjMcee4y33nqLY8eOlbsey2v2799/A7e+crKzs3nllVd48MEHadOmDQ888AATJkzg6NGjVfo+M2bMIDg4mJdeeqnU8qo8Runp6ZV+ja33f+mllwgODmbGjBnXvU32KCkp4eDBg6WWbdmyheDgYDp27HhDtuFWoVD+HWjTpg3t27e3/teyZUu8vLzYs2cPc+bM4ZFHHmHNmjU3ezOvSXZ2Nk8++STz588nOzuboKAg8vPz+eyzz+jbty/79u272Ztol1OnTvG3v/2NqKiom70plbZr1y769evHsmXLbvam3BJcbvYGyPV7++23uf3228ssP378ODNnzmTVqlWMGzeOOnXq0Llz51I1q1evBqBp06Y3ZFsra/LkyaSnpxMeHs6bb76Jl5cXBQUFTJ06leXLl/P888+TmJhIjRo1HLYNVXGMNm3axJdffknr1q1vyvtfj0WLFrF79+4yPeKQkBBWr16Ni4tipCqpp/w71qhRI9544w169uzJxYsXmThxIhcvXixVExAQQEBAADVr1rxJW1m+9PR01q5di4eHBzNnzsTLywsANzc3/vGPfxAQEEB6ejpfffWVQ7fjZh+jm/3+5alVqxYBAQHccccdN3tTflcUyr9zTk5OTJkyBTc3N44ePcoXX3xxszfJbv/+978xDINu3brh7e1dqq1GjRo89thjwG89SZHfA4XyLaBevXo8+OCDAGzYsKFU29UuIq1bt45du3bxzDPPEBYWRrt27XjyySf59ttvAcjLy+O1116jW7dutGnThm7duvHWW2+V6Y0PHjy43IuSV/53+UW2nTt3AtCuXTub+3XPPfcAkJKSUqnjsWXLFiIjI+nUqRPt27dn+PDh1veypbwLfVu2bOHZZ5+le/futGnThs6dOzN8+PAyHxLdunXj5ZdfBmD37t0EBwfTrVs34LcLto8++ihpaWlERERw9913c//997Nw4cKrvr/FTz/9xNNPP027du249957efrpp9m4cWOZuuXLlxMcHGz9MLvSwoULCQ4OZvDgwaW27fPPPwfgo48+KvU7utqFvoKCAubNm0e/fv1o164dbdu25dFHHyU2Npbs7Owy9ZZzbs2aNezbt4/Ro0fTqVMn7r77bh555BHmzJlDYWGhze3+vdFg0C2iXbt2rFmzhh9++MHu13z99desWLECV1dX7rzzTo4cOcL27dsZMWIE7733Hm+88QYHDhygWbNmNG7cmJ9//pk5c+aQmZnJK6+8Yl1PUFAQRUVFdr3nnXfeaf3fhw8fBrA5Xg7g5+cHQGZmJufPn8fT07PC9S9atIiYmBgMw8DX1xc/Pz+2bt3KgAEDaNWqlV3bCJCYmEhUVBQlJSU0atSI4OBgMjMz+fbbb/n222/ZtWsX48ePBy5diHV1deXQoUN4eHjQsmVLGjRoUGp9586dY/jw4eTk5NCiRQsyMjIICAiocDt++OEHFixYAFw6zpmZmSQnJ5OcnMyYMWN49tln7d6nK7m5udG+fXsOHz7M6dOnadiwIU2aNCn1O7IlKyuLoUOHkpqaipOTE82bN6dmzZocOHCA9957j5UrVzJ37lyb+/f9998zbtw4APz9/alVqxZpaWm89dZb/PTTT7z33nvXvD/VhiHV0pEjR4ygoCAjKCjIOHLkSIX1a9eutdYXFhZal1uWpaamWpeNHz/eunz06NHGuXPnDMMwjHPnzhl9+/Y1goKCjJYtWxoPPvigsXv3buvrPvzwQyMoKMi46667jDNnzlz3Pt5zzz1GUFCQ8eOPP9psz83NrdQx2L9/v9GqVSsjODjYmD9/vlFSUmIYhmFkZmYagwcPtq5r/PjxpV535TEqLi427rvvPiMoKMhYtWpVqdrPP//cCA4ONlq2bFlqmz777DMjKCjI+POf/1yq/vLfY48ePYzMzEzDMAzjzJkz1u2r6HfUr18/49ixY4ZhGEZJSYmxYMECIzg42AgODjZSUlIq3AaLBQsWGEFBQcagQYNKLbe816uvvlpq+ffff28EBQUZYWFhpZY/9dRTRlBQkPHoo48aaWlp1uXHjx+3tvXo0cPIz8+3uT8jRowwTp06ZW375JNPrG0//fSTzW3/PdHwxS3i8l6kra+PttSpU4fp06dbL7B5eXnx5JNPApfmrU6dOrVU7/Kpp56iZs2aFBcXc+DAgeve5vz8fADc3d1ttl++vKCgoML1ffTRRxQVFdG3b18GDx6Mk5MTcGl4Jy4ujttuu82u7Tp9+jSZmZnUqVOHP/7xj6Xa+vbtS//+/XnkkUfIzc21a30WTz/9NPXq1QPA29vbun1X4+npybvvvkujRo2AS9cQBg0axOOPP45hGMybN69S23C9fvjhB77//nvc3NyYM2dOqd5ww4YNeeedd2jUqBGHDh3is88+K/N6b29v3n77berXr29d9tRTT9GsWTMAduzY4fiduMkUyreIy8d57fnHDtC+fXtrIFtYhgxcXFzo1KlTqTYXFxfrBbnz589fz+YCVDjNraSkpFLr27RpEwB/+tOfyrR5e3vTvXt3u9ZTt25dateuTXZ2NhMmTCjzATRt2jRef/11WrZsWants4yRV0b37t3LDIUA1nHj5ORkiouLK73ea/X1118DEB4eTpMmTcq0e3l5WbfNUnu5sLAwmx/C/v7+AJX+oKuONKZ8i7j8ZK5du7Zdr2nYsGGZZa6ursClf1y2pmhZ2o3L/nZCTEwMe/bsses9w8PDeeaZZ4BLU64uXrxYbi/48gs/5fWmLS5cuMCpU6cAyh2rDQ4OtmsbXVxcGDNmDK+88grLly9n+fLlNG7cmP/3//4f4eHhdOnShVq1atm1rsvZCteK3HXXXTaXBwYGApc+HE+dOmXtSTvaoUOHrrpdgHWutqX2crbOOfjt91vZD+LqSKF8i8jIyAAuXTSzd76rh4dHlbz3/v372b59u121l8959fb2Jicnh7Nnz9qsvXy5j4/PVdebk5Nj/d/l7Ze9wxdw6Sv1HXfcwccff8zWrVs5duwYCQkJJCQk4OnpSWRkJCNHjrR7fXDpwlpllbcvly+/cOFCpdd7rSzfkK520dWybba+TVk+1Mtj3AJ/KEmhfIv48ccfgWv7iny9LLMDKqt58+b8/PPP5T7j4tdffwUu9TAr6plePs85Ly/P5rcFyxi2vcLDwwkPD+fcuXNs2bKF7777jg0bNvDrr7/y9ttv4+npyZAhQyq1zsoqL3AvDzx7P2yqIrwtgXu1YYZz586VqpXSNKZ8Czh27BibN28GKHNhyszatGkDXJqHa4vlok/btm0rXJebmxuNGzcGYO/evTZrLN8mKlJYWMj+/fut66lduzbdu3fn73//O+vXr+fPf/4zcOnmF0e78iFBFpZtq1u3rvXioWWMvrz5vpbhnethGfst7xjDpbnagPXinZSmUL4FxMTEUFJSQosWLaw3kVQHDz/8MABfffVVmSGM4uJi600Nti7c2fLQQw8BEB8fX6YtLy/P7oc2ffXVV/Tu3ZsXXnihzNdpZ2dn6wXQy8c/nZ0d809t3bp1NnulixcvBij1+65Tpw4AR48eLTNOX1RUZPOGE/jtwrA9Qwddu3YFYOPGjTa/4eTm5rJixQoA7r///grXdytSKP+OHTp0iDFjxrB+/XpcXV2JiYlx6IN7qlrLli0JDw/n/PnzjBkzhjNnzgCXpr9NmjSJ9PR0/P39reFdkeHDh+Ph4cG6deuIi4uz3tCSk5PD888/b3dPsWvXrnh6epKens706dNLfe3/9ddf+fDDDwF44IEHrMstX9VPnjxZpXemZWZmMnbsWOuYeXFxMe+88w5r1qzBzc2NyMhIa21ISAg1atQgLy+Pt956yzor49y5c0ycONHmhTf4bXzYMlx0Nffeey8dO3aksLCQZ555ptSjSk+cOMGzzz7LiRMnaNq0KU888cS17vbvmsaUfweee+65UhfvCgoKOHnypDVkvLy8eO2112jfvv3N2sRrNm3aNAYMGMCWLVt48MEHad68Ob/88gvZ2dnUrl2b2bNn290L9fPzY8aMGbzwwgu88847LFmyhMaNG5ORkcGFCxfo2rWrzWlaV/L09GTmzJmMGjWK+fPn89lnn9GsWTMKCws5fPgwRUVFtG7dmv/5n/+xviYwMBAnJydOnTpFz549adSokbU3ez0eeughNmzYQHh4OP7+/pw4cYLMzExcXFz45z//WWqmiY+PD0899RTz5s1j3rx5fPHFF/j6+nLw4EEKCgr461//avOOOcuslKSkJHr37k1YWBiTJ08ud5veeOMNnn76afbv388jjzxCQEAArq6uHDhwgKKiIpo0acLs2bPLTLeUSxTKvwP//e9/S/1/FxcXateuTWhoKF26dCEiIqLC2Qlm1ahRIz777DPeeecd/vOf/7B//35q167No48+yujRoyu85fdKPXr0YOnSpcyZM4eUlBTS09Np2bIlzz77LL/88otdoQyX5gcvXLiQTz75hB07dnDgwAHc3d1p1aoVf/zjHxk0aFCpD0p/f3/+8Y9/8P7773Ps2DEuXrxIZmZmpbbdlgcffJDBgwcTFxfHnj17cHNz46GHHmLkyJHWMfnLjR8/noCAAJYsWUJaWhoFBQWEhYUxcuRIatasaTOU//znP5OWlkZiYiKHDh2q8MJhgwYNWLZsGYsWLWL16tXWce+AgAD+8Ic/MHDgQOtQipTlZNwKc0xERKoJjSmLiJiIQllExEQUyiIiJqJQFhExEc2+4NLttf/9739p0KBBtZrHKyLmVlxczKlTp2jTpk2FD82yUChzaUrZwIEDb/ZmiMjv1KJFi7j33nvtqlUo89sjExctWnTDHnEoIr9/x48fZ+DAgZV6LKtCmd8e1NKoUaNy/x6ciMi1qsywqC70iYiYiEJZRMREFMoiIiaiUBYRMRGFsoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWW4phReLq6RGxFF0m7XcUmq61qD3CyuvWpP4Rp8btDUiZamnLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFrlCRVPiNGVOHElT4kSuUNG0OU2ZE0dST1lExESuOZSXL19OcHAwP/zwg832gwcP8vzzzxMeHk7btm3p3bs3CxcupKSkxGZ9Tk4Or732Gj179iQkJIRu3brx6quvkpuba7O+uLiYJUuW0LdvX9q1a0fnzp0ZO3YsBw8evNZdEhG56a4plH/88UdiYmLKbd+3bx/9+vVj1apV+Pn50aVLF44fP05MTAxRUVFl6nNzcxk0aBBz587FycmJrl274uTkxLx584iIiODcuXNlXjNp0iSmTJnC8ePHuf/++2nSpAmrV6/mscceY8+ePdeyWyIiN12lQ3nt2rUMHz6cvLw8m+2GYRAVFUVubi4zZ85k8eLFzJ49m6SkJIKDg0lMTCQpKanUa2JjY0lNTaV///6sXr2auLg4kpKS6NOnD2lpacTGxpbZhuXLl9O6dWu++uorZs2aRUJCAtHR0eTl5fHSSy9hGEZld01E5KazO5SPHz9OVFQUo0ePpqSkhPr169usS05OJjU1lbCwMPr0+e2CiI+PD1OmTAFgwYIF1uU5OTnEx8fj5eXF+PHjcXa+tEkuLi5MmTKFOnXqkJCQUOpD4KOPPgLgpZdeonbt2tblf/nLX7jvvvtITU1ly5Yt9u6aiIhp2B3KsbGxrFy5kjZt2rB06VKaN29us27Tpk0AdO/evUxbaGgo9erVIyUlxTpWvG3bNvLz8+nUqRNeXl6l6j09PencuTP5+fls27YNuBTiO3bswNvbm3vvvbfMe1je95tvvrF310RETMPuUG7evDkzZswgPj6e4ODgcuvS0tIACAoKstnu7+9PSUkJ6enppeoDAwPLfV+A1NRUANLT0zEMg4CAAGuv2lb9/v377dkt+Z3RHGKp7uyepzxixAi76k6ePAlAgwYNbLZblmdmZgJw6tQpu+pPnz5tV72vr2+perm1aI6xVHdVPk/5woULALi7u9tstyy3jBFbftaqVatK6t3c3ErViYhUJ1UeypYhBScnJ5vtllkRlp+OrhcRqU6qPJQ9PDwAyM/Pt9leUFBQqs7eekvPuKL6wsLCUvUiItVJlYeyZUzXMmZ8pSvHhKu6vqIxbRERM6vyULbMorDMqricYRhkZGRQo0YNAgICKqwHrLM0LDM+WrRogbOzs3X5lTIyMoDyZ3+IiJhZlYdyly5dAFi/fn2Ztu3bt5OVlUVoaKh1TnKHDh1wd3dn8+bNZS7OnT9/ns2bN+Ph4UFoaCiA9X+fPn2a7du3l3mPdevWARAeHl6l+yUiciNUeSiHhYURGBhIcnIyy5Ytsy7PysoiOjoagGHDhlmXe3h40LdvX7Kzs4mOjqaoqAiAoqIipk2bRk5ODhEREaVuLBkwYAAA0dHRZGVlWZcvXbqU7777jtatW9OxY8eq3jUREYer8ucpOzs7M336dIYMGcLkyZNJSEjA19eXrVu3kp2dTf/+/enWrVup14wdO5YtW7awYsUKUlJSaNWqFXv27OHIkSO0atWK0aNHl6rv1asXa9eu5csvv+QPf/gDYWFhnDhxgp07d3LbbbcxY8aMqt4tEZEbwiHPUw4JCSE+Pp6ePXty+PBhkpOT8fPzIzo6mqlTp5ap9/b2ZsmSJQwePJiioiI2bNiAs7MzkZGRzJ8/H09PzzKvef3113n55Zfx9fVl48aNnDhxgkceeYT4+Phy7w4UETG7a+4pX/5QIVtatGhBXFyc3evz9vZm0qRJTJo0ya56FxcXhg4dytChQ+1+DxERs9NfHhERMRGFsoiIiSiURURMRKEsUkn2PB5UjxCVa1XlU+JEfu8qejwo6BGhcu3UUxYRMRGFsoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWUTERBTKIiImolCWakN3ycmtQHf0SbWhO+nkVqCesoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWUTERBTKIiImolAWETERhbKIiIkolEVETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYc+5H7lypUsWrSI/fv3U1JSgr+/P4899hiDBg2iRo0apWoPHjzIrFmzSElJ4ezZszRr1oyIiAgGDBiAs3PZz46cnBzef/991q1bx7Fjx6hfvz49evRg1KhReHl5OXK3REQcxmE95ZkzZxIVFcXevXtp3749HTt25Oeff2b69OmMGTMGwzCstfv27aNfv36sWrUKPz8/unTpwvHjx4mJiSEqKqrMunNzcxk0aBBz587FycmJrl274uTkxLx584iIiODcuXOO2i0REYdySE85NTWVjz76CB8fHz799FP8/f0BOHHiBE8++STr1q1j7dq19OzZE8MwiIqKIjc3l5kzZ9Knz6U/55OVlcXQoUNJTEzk4YcfpmfPntb1x8bGkpqaSv/+/YmOjsbZ2ZmioiImTJjAypUriY2NZfLkyY7YNRERh3JIT/m7777DMAz+9Kc/WQMZoGHDhgwYMACAbdu2AZCcnExqaiphYWHWQAbw8fFhypQpACxYsMC6PCcnh/j4eLy8vBg/frx1aMPFxYUpU6ZQp04dEhISyMvLc8SuiYg4lENC2cnJCbjUM77SmTNnAPD29gZg06ZNAHTv3r1MbWhoKPXq1SMlJYXc3FzgUpjn5+fTqVOnMmPHnp6edO7cmfz8fGvoi4hUJw4J5S5duuDk5MSaNWv44IMPyMrKIicnh4SEBObPn0+dOnV4/PHHAUhLSwMgKCjI5rr8/f0pKSkhPT29VH1gYKDN+ubNmwOXhlBERKobh4RyQEAAMTExuLu788Ybb9C5c2c6dOjAxIkTad26NUuXLqVx48YAnDx5EoAGDRrYXJdleWZmJgCnTp2yq/706dNVt0MiIjeIw2ZftG/fns6dO+Ph4UGnTp2477778PT0ZNeuXXz66afW2RcXLlwAwN3d3eZ6LMstY8SWn7Vq1bKrXkSkOnHI7IsdO3bw9NNP06RJExITE7n99tuBS2PMo0aNYv78+Xh5efHcc89ZL9RZxqGvZAlvy8/K1ouIVCcO6SlPnz6d8+fP88orr1gDGS7NvnjzzTdxcXHh448/5sKFC3h4eACQn59vc10FBQUA1jp768vrSYuImFmVh3J+fj47d+6kdu3ahISElGlv2rQp/v7+5OXlcfjwYXx9fYHfxoyvdOUYcmXrRUSqkyoP5XPnzmEYRpnbqC9nabt48aJ1FoVlVsXlDMMgIyODGjVqEBAQAHDVesA6SyM4OPjad0JE5Cap8lCuV68e3t7enD17lp07d5ZpP3HiBOnp6bi6utK8eXO6dOkCwPr168vUbt++naysLEJDQ61zkjt06IC7uzubN28uczHv/PnzbN68GQ8PD0JDQ6t618TBCi8W3+xNELnpqjyUnZ2d6devHwATJ04sdQNJVlYW48aN4+LFizz++ON4enoSFhZGYGAgycnJLFu2rFRtdHQ0AMOGDbMu9/DwoG/fvmRnZxMdHU1RUREARUVFTJs2jZycHCIiIvRQomqopmsNer+wstz/RG4FDpl9MWbMGHbu3MnWrVt5+OGH6dChA05OTvz000/k5ORwzz33MH78eOBSiE+fPp0hQ4YwefJkEhIS8PX1ZevWrWRnZ9O/f3+6detWav1jx45ly5YtrFixgpSUFFq1asWePXs4cuQIrVq1YvTo0Y7YLRERh3PI7As3Nzc++ugjJkyYQIsWLUhJSWHr1q00btyYF154gQULFlhnUQCEhIQQHx9Pz549OXz4MMnJyfj5+REdHc3UqVPLrN/b25slS5YwePBgioqK2LBhA87OzkRGRjJ//nw8PT0dsVsiIg7nsOcpu7q6MmTIEIYMGWJXfYsWLYiLi7N7/d7e3kyaNIlJkyZd6yaKOEzhxWJqupZ/sbuidrl1OfQh9yK3Ksv4eHkS3+hTbpvc2vTnoERETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWyiIiJKJRFRExEoSwiYiIKZRERE1Eoi4iYiEJZRMREFMoiIiaiUBYRMRGFsoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWUTERBTKIiImolAWETERhbLITVB4sbhKauT3x+Vmb4DIraimaw16v7DyqjWJb/S5QVsjZqKesoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWUTERBw6Je7o0aO88847fPvtt2RlZVG3bl26du3KmDFjaNCgQanagwcPMmvWLFJSUjh79izNmjUjIiKCAQMG4Oxc9rMjJyeH999/n3Xr1nHs2DHq169Pjx49GDVqFF5eXo7cLRERh3FYT3nXrl306dOHzz77jDp16hAeHo6zszPLli1jwIABZGdnW2v37dtHv379WLVqFX5+fnTp0oXjx48TExNDVFRUmXXn5uYyaNAg5s6di5OTE127dsXJyYl58+YRERHBuXPnHLVbIiIO5ZBQLiwsZNy4cZw7d45JkyaRmJjIO++8w9q1a+nZsyc///wzs2bNAsAwDKKiosjNzWXmzJksXryY2bNnk5SURHBwMImJiSQlJZVaf2xsLKmpqfTv35/Vq1cTFxdHUlISffr0IS0tjdjYWEfsloiIwzkklFevXs2hQ4fo3bs3gwcPti53c3Pj5Zdfpn79+hw8eBCA5ORkUlNTCQsLo0+f3+5g8vHxYcqUKQAsWLDAujwnJ4f4+Hi8vLwYP368dWjDxcWFKVOmUKdOHRISEsjLy3PEromIOJRDQnnt2rUADBs2rExb48aNSU5O5sMPPwRg06ZNAHTv3r1MbWhoKPXq1SMlJYXc3FwAtm3bRn5+Pp06dSozduzp6Unnzp3Jz89n27ZtVbpPIiI3gkMu9O3ZswdXV1datmzJsWPHSExM5Oeff8bb25sePXoQEhJirU1LSwMgKCjI5rr8/f05ffo06enptG3b1lofGBhos7558+YApKamEh4eXpW7JSLicFUeyoWFhRw7doxGjRqxZs0aJk6cyIULF6zt//rXvxg+fLj1At7JkycByszGsLAsz8zMBODUqVN21Z8+fboK9kZE5Maq8uELyzBDdnY248ePp3v37qxZs4Zt27bx1ltv4e3tzYcffsjSpUsBrIHt7u5uc32W5ZYxYsvPWrVq2VUvIlKdVHkoFxQUAJfCtmPHjrz++uv4+/tz22230atXL2bMmAHAO++8g2EY1gt1Tk5ONtdnGEapn5WtFxGpTqo8lC/vwT755JNl2rt27UrDhg05ceIEhw8fxsPDA4D8/Hyb67OEvKXO3vryetJyc+iB7SL2qfIx5dq1a+Pq6srFixe5/fbbbdb4+flx4sQJzpw5g6+vL3v37iUzM5OAgIAytVeOIfv6+gK/jTFXVC/moIe6i9inynvKNWrUsIbriRMnbNZYAtXHx8c6i8Iyq+JyhmGQkZFRap1XqwdIT08HIDg4+Dr2QkTk5nDIPOUHHngAgDVr1pRpy8jI4OjRo/j6+tK0aVO6dOkCwPr168vUbt++naysLEJDQ61zkjt06IC7uzubN28uczHv/PnzbN68GQ8PD0JDQ6t6t0REHM4hofyXv/wFDw8PVqxYQWJionV5dnY2kyZNoqSkhIEDB+Ls7ExYWBiBgYEkJyezbNkya21WVhbR0dFA6ZtQPDw86Nu3L9nZ2URHR1NUVARAUVER06ZNIyddsNgSAAAbJklEQVQnh4iICD2USESqJYfcPNKkSRNeeeUVXnzxRcaNG8e8efPw9fVlx44dnDlzhk6dOjF8+HDg0myK6dOnM2TIECZPnkxCQgK+vr5s3bqV7Oxs+vfvT7du3Uqtf+zYsWzZsoUVK1aQkpJCq1at2LNnD0eOHKFVq1aMHj3aEbslIuJwDntKXK9evUhISKBnz578+uuvfPfdd9SrV48XXniBuXPn4urqaq0NCQkhPj6enj17cvjwYZKTk/Hz8yM6OpqpU6eWWbe3tzdLlixh8ODBFBUVsWHDBpydnYmMjGT+/Pl4eno6ardERBzKoc9Tvuuuu4iLi7OrtkWLFnbXwqVgnjRpEpMmTbrWzRMRMR395RERERNRKIuImIhCWUTERBTKIiImolAWETERhbKIiIkolEVETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWyiIiJKJRFRExEoSxiUoUXi6+rXaonhz7kXkSuXU3XGvR+YWW57Ylv9LmBWyM3inrKIiImolAWETERhbKIiIkolEVETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWyiIiJKJSlSujhOCJVQw8kkiqhh+eIVA31lEVETOSGhfLZs2e5//77CQ4Ottl+8OBBnn/+ecLDw2nbti29e/dm4cKFlJSU2KzPycnhtddeo2fPnoSEhNCtWzdeffVVcnNzHbkbIiIOdcNCOTo6mlOnTtls27dvH/369WPVqlX4+fnRpUsXjh8/TkxMDFFRUWXqc3NzGTRoEHPnzsXJyYmuXbvi5OTEvHnziIiI4Ny5c47eHRERh7ghofzFF1+wevVqm22GYRAVFUVubi4zZ85k8eLFzJ49m6SkJIKDg0lMTCQpKanUa2JjY0lNTaV///6sXr2auLg4kpKS6NOnD2lpacTGxt6I3RIRqXIOD+UTJ04QExNDu3btqFGjRpn25ORkUlNTCQsLo0+f3y4G+fj4MGXKFAAWLFhgXZ6Tk0N8fDxeXl6MHz8eZ+dLu+Di4sKUKVOoU6cOCQkJ5OXlOXjPRESqnsNDeeLEiRQUFDBjxgyb7Zs2bQKge/fuZdpCQ0OpV68eKSkp1rHibdu2kZ+fT6dOnfDy8ipV7+npSefOncnPz2fbtm1VvCciIo7n0FD+9NNP2bRpE+PGjeOOO+6wWZOWlgZAUFCQzXZ/f39KSkpIT08vVR8YGGizvnnz5gCkpqZe17aLiNwMDgvln3/+mddee41OnToxcODAcutOnjwJQIMGDWy2W5ZnZmYCWC8WVlR/+vTpa9twEZGbyCGhXFxcTFRUFE5OTvzzn//Eycmp3NoLFy4A4O7ubrPdstwyRmz5WatWLbvqRUSqE4eE8ty5c/nxxx95+eWX8fPzu/oG/N+FuvKC2zCMUj8rWy8iUp1UeSjv27ePWbNmER4ezhNPPFFhvYeHBwD5+fk22wsKCkrV2VtfXk9aRMTMqvzZF2+99RYXL16kqKiIcePGlWqz3J1nWT5hwgR8fX3Zu3cvmZmZBAQElFnflWPIvr6+wG9jzBXVi4hUJ1Ueypax3OTk5HJrEhMTAfjb3/5GYGAgGzduJC0tjY4dO5aqMwyDjIwMatSoYQ1sy6wLyyyMK1lmaZR3O7eIiJlVeShffqPHlVq1akVxcXGp6WpdunRh7ty5rF+/vswsje3bt5OVlUVYWJh1TnKHDh1wd3dn8+bN5OXlWYczAM6fP8/mzZvx8PAgNDS0ivdMRMTxbvpT4sLCwggMDCQ5OZlly5ZZl2dlZREdHQ3AsGHDrMs9PDzo27cv2dnZREdHU1RUBEBRURHTpk0jJyeHiIiIMjeWiIhUBzf9ecrOzs5Mnz6dIUOGMHnyZBISEvD19WXr1q1kZ2fTv39/unXrVuo1Y8eOZcuWLaxYsYKUlBRatWrFnj17OHLkCK1atWL06NE3aW9ERK7PTe8pA4SEhBAfH0/Pnj05fPgwycnJ+Pn5ER0dzdSpU8vUe3t7s2TJEgYPHkxRUREbNmzA2dmZyMhI5s+fj6en543fCRGRKnBDe8p79uwpt61FixbExcXZvS5vb28mTZrEpEmTqmLTRERMwRQ9ZRERuUShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWySDVVeLH4utrFnG76bdYicm1qutag9wsry21PfKNPuW1iXuopi4iYiEJZRMREFMoiIiaiUBYRMRGFsoiIiSiURURMRKEsdtGcV5EbQ/OUxS6aEytyY6inLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWyiIiJKJRFRExEoSwiYiIKZRERE1Eoi4iYiEJZRMREFMoiIiaiUBYRMRGFsoiIiTjsecrFxcUsXryYzz//nIyMDIqLi2natCm9evUiMjISNze3UvUHDx5k1qxZpKSkcPbsWZo1a0ZERAQDBgzA2bnsZ0dOTg7vv/8+69at49ixY9SvX58ePXowatQovLy8HLVbIiIO5ZCecnFxMSNHjiQmJoaMjAzatm1LWFgYJ0+eJC4ujsGDB3PhwgVr/b59++jXrx+rVq3Cz8+PLl26cPz4cWJiYoiKiiqz/tzcXAYNGsTcuXNxcnKia9euODk5MW/ePCIiIjh37pwjdktExOEcEsrx8fF8/fXXBAcHs2bNGj7++GPmzp1LUlIS7dq146effuLdd98FwDAMoqKiyM3NZebMmSxevJjZs2eTlJREcHAwiYmJJCUllVp/bGwsqamp9O/fn9WrVxMXF0dSUhJ9+vQhLS2N2NhYR+yWiIjDOSSUP//8cwAmTJhAw4YNrct9fHyYOnUqAKtWrQIgOTmZ1NRUwsLC6NOnT6naKVOmALBgwQLr8pycHOLj4/Hy8mL8+PHWoQ0XFxemTJlCnTp1SEhIIC8vzxG7JiLiUA4J5bp169K8eXNCQkLKtN15550AnDx5EoBNmzYB0L179zK1oaGh1KtXj5SUFHJzcwHYtm0b+fn5dOrUqczYsaenJ507dyY/P59t27ZV5S6JiNwQDgnlOXPm8OWXX+Lh4VGmbdeuXQA0atQIgLS0NACCgoJsrsvf35+SkhLS09NL1QcGBtqsb968OQCpqanXsQciIjfHDZ0SZxgGcXFxAPTo0QP4rcfcoEEDm6+xLM/MzATg1KlTdtWfPn26irZaROTGuaGh/Oabb7J161bq169PZGQkgHUWhru7u83XWJZbxogtP2vVqmVXvYhIdXLDQvntt9/mgw8+oGbNmsTGxuLj43NpA/7vQp2Tk5PN1xmGUepnZetFblWFF4urpEZuLIfdPGJRVFTEtGnTWLp0KW5ubsyaNYsOHTpY2y3jzvn5+TZfX1BQUKrO3vryetIit4qarjXo/cLKq9YkvtHnqu1y4zm0p3z+/HmeeeYZli5dym233caHH35IeHh4qRpfX1/gtzHjK105hlzZehGR6sRhoZydnc3gwYPZtGkTjRs3ZtGiRaV6yBaWWRSWWRWXMwyDjIwMatSoQUBAQIX1gHWWRnBwcJXsh4jIjeSQUC4sLGTEiBHs3r2bFi1asGTJknKnvHXp0gWA9evXl2nbvn07WVlZhIaGWuckd+jQAXd3dzZv3lzmYt758+fZvHkzHh4ehIaGVvFeiYg4nkNCOS4ujh07dtC4cWMWLFhgnZNsS1hYGIGBgSQnJ7Ns2TLr8qysLKKjowEYNmyYdbmHhwd9+/YlOzub6OhoioqKgN/GrnNycoiIiNBDiSpBF3tEzKPKL/SdPXvWelu0j48P06dPL7f29ddfx9nZmenTpzNkyBAmT55MQkICvr6+bN26lezsbPr370+3bt1KvW7s2LFs2bKFFStWkJKSQqtWrdizZw9HjhyhVatWjB49uqp363dNF4REzKPKQ3nnzp3WmRG7d+9m9+7d5da+/vrrAISEhBAfH09cXBxbtmzhwIED3HHHHTz//PM88cQTZV7n7e3NkiVLmD17NuvWrWPDhg00btyYyMhInnnmGTw9Pat6t0REbogqD+UHHnjgmm5xbtGihfVuP3t4e3szadIkJk2aVOn3EhExK/3lERERE1Eoi4iYiEJZRMREFMoiIiaiUBYRMRGFsoiIiSiURURMRKEscgur6BZ73YJ/4zn8ecoiYl4V3WKv2+tvPPWUbwHq7YhUH+op3wLUGxKpPtRTFhExEYWyiIiJKJRFRExEoSwiYiIKZRERE1Eoi4iYiEJZRMREFMoiIiaiUBaRctlzN6juGK1auqNPRMpV0d2goDtCq5p6yiIiJqJQFhExEYWyiIiJKJRFRExEoVzN6cq3yO+LZl9Uc7o6LvL7op6yiIiJKJRF5Lroj69WLQ1fiMh10Z8bq1rqKZucehkit5Zq31P+7rvvmDNnDqmpqVy8eJHWrVszYsQIunTpcrM3rUqoFyJya6nWPeXly5czbNgwfvzxR0JCQmjXrh0//vgjkZGRLF269GZvnoighxpVVrXtKZ88eZIpU6ZQu3ZtPv30U4KCggDYuXMnw4YN45VXXqFr1640bNjwJm+pyK1N0zYrp9r2lBcuXEhhYSFDhw61BjJASEgIkZGRFBQUqLcsItVOtQ3lTZs2AdC9e/cybQ8//DAA33zzzQ3dJhGR61UtQ9kwDNLS0nB2dqZ58+Zl2u+8806cnZ1JS0vDMIybsIX20TiayCWa6/ybajmmnJ2dTWFhIT4+PtSsWbNMu4uLC3Xr1uX06dOcP38eLy+vq66vuPjSL/z48eN2b8PFohJcXa7+mVZYVELNCmoiX/nqqu1zJz7Mxbysctt/+eWXq7bbU3Mj1lFdtvP3tI7qsp0AJ08cu+q/hbkTH77q683KkimWjLGHk2HmrmQ5jh07RteuXWnSpAn/+c9/bNZ069aNo0eP8s0331R4se+HH35g4MCBjthUEREWLVrEvffea1dttewpOztXPOpSmc+aNm3asGjRIho0aECNGjWuZ9NERKyKi4s5deoUbdq0sfs11TKUPTw8ACgoKCi3xtJWq1atCtfn7u5u96eYiEhl3HHHHZWqr5YX+ry8vPDw8ODMmTMUFRWVaS8qKuLMmTO4ublx22233YQtFBG5NtUylJ2cnGjRogXFxcUcOnSoTPvBgwcpKSkpNX9ZRKQ6qJahDFifbbFu3boybZZl4eHhN3SbRESuV7UN5cceeww3Nzf+9a9/8d///te6fNeuXcydOxd3d3cGDBhwE7dQRKTyquWUOItFixYxbdo0XF1d6dSpE4ZhsGXLFoqKipgxYwZ9+uh+ehGpXqp1KANs2LCBuXPnsmfPHmrWrElwcDB//etf6dy5883eNBGRSqv2oSwi8ntSLecp32jLly/n5ZdftnlXjuXuwvK0b9+exYsX2/U++fn5fPLJJ/z73//ml19+oXbt2nTt2pUxY8bg6+t7PbvgUOUdn8GDB7N169YKXz9q1ChGjx5dYd2QIUP4/vvvy21fu3ZtpeeEOkJxcTGLFy/m888/JyMjg+LiYpo2bUqvXr2IjIzEzc2tVP3BgweZNWsWKSkpnD17lmbNmhEREcGAAQPsulHKojqcP5U9Nhs3bmT+/Pns2rWLvLw8GjRoQJcuXRg5ciSNGjWy+32ry7kDCuUK/fjjj8TExJTbvmfPHgCCg4NtTsHz9/e3630uXrzIyJEjSU5OpnHjxoSHh5ORkUF8fDxff/01y5Ytw8/P79p2woGudnzuu+++cm9xz8vLY/369QDcdddddr3Xvn378PDw4KGHHrLZ7unpadd6HKm4uJiRI0fy9ddf4+HhQdu2bXFxceGnn34iLi6OjRs38sknn1hvatq3bx8DBw4kNzeX9u3bc/fdd7NlyxZiYmLYsWMHr7/+ul3vWx3On8oemw8++IA33ngDZ2dnQkJCqFevHnv37mXp0qV89dVXLFy4kICAALveuzqcO1aGlCspKclo166dERQUZAQFBRnbtm0rUzNr1iwjKCjIWLly5XW910cffWQEBQUZI0aMMAoKCqzL33zzTSMoKMj43//93+tavyPYc3zK8+KLLxpBQUHGP//5T7vqf/nlFyMoKMgYOnTotW7uDbF48WIjKCjI6N27t3H8+HHr8tOnTxsRERFGUFCQ8frrrxuGYRglJSVG7969jaCgIGPFihWlai3L16xZY9f7VofzpzLH5sCBA8Zdd91l3HPPPcb27duttYWFhcbUqVONoKAgo3///na9b3U5dyyq7ZQ4Rzp+/DhRUVGMHj2akpIS6tevX26tpafcunXra34/wzCYN28eTk5OTJ48udST75577jn8/f3ZsGEDR44cueb3qEqVOT62JCYmsnLlSoKCgnj++eftes3evXuB6zvON8Lnn38OwIQJE0p9S/Dx8WHq1KkArFq1CoDk5GRSU1MJCwsrNVPIx8eHKVOmALBgwYIK37O6nD+VOTYrV66kuLiYYcOG0a5dO2utq6srEyZMwMfHhx07dnD06NEK37e6nDsWCmUbYmNjWblyJW3atGHp0qU2n9lssXfvXjw8POweprAlNTWVEydO0LJlS26//fZSbc7OznTr1g0wz0P7K3N8rnT+/HlmzJgBwNSpU20+etWWqvjwuxHq1q1L8+bNCQkJKdN25513Apf+lBlc/Q81hIaGUq9ePVJSUsjNzb3qe1aX86cyx8bV1ZXg4GA6dOhQptbV1dW6n5b6q6ku546FxpRtaN68OTNmzOBPf/rTVS+0nD17ll9//ZXWrVszb948Vq5cyeHDh6lduzYPPvggo0aNsutvBKalpQEQGBhY7vYA7N+//xr2purZe3xsmTNnDqdOnaJXr16Ehoba/TrLP6zs7GyGDx/O7t27KSgooE2bNqb66+Vz5swpt23Xrl0A1gtUlt97eY8D8Pf35/Tp06Snp9O2bdty11tdzp/KHJsxY8YwZswYm7V5eXnWfbbnYl91OXcs1FO2YcSIEfTt27fCwLF8Ldq9ezdvvfUW9erVo2PHjhQXF7Ns2TIef/xxMjIyKny/U6dOAdCgQQOb7Zblp0+frsxuOIy9x+dKZ8+eZcGCBTg5OfHss89W6rWWYz1lyhROnjxJhw4duP3229m6dSuRkZF88sknlVrfjWYYBnFxcQD06NED+K2XV9HvPTMz86rrrm7nz5VsHZur+de//kVeXh533303jRs3rrC+up076ilfB8sncGBgIO+99x5NmzYFLn2ST548mS+++IJx48axfPnyq64nLy8PKP8xo+7u7qXqqqvFixdz4cIFunXrRosWLex+XVZWFsePH8fFxYUZM2bw6KOPWttWr17Niy++yIwZMwgLC7N7JseN9uabb7J161bq169PZGQkABcuXAB++/1eyd7fe3U/f2wdm/Js3LiR999/H2dnZ1588cUK110dzx31lK/D0KFDWbduHfPnz7cGMlx63vM//vEPGjZsyO7du9mxY8dV12PpcTo5OdlsN/7v/h6jGt/nU1xczKJFiwAq/Id3JR8fHzZv3syqVatK/aMC6NWrFwMHDrTOfzWjt99+mw8++ICaNWsSGxuLj48PUHW/9+p8/pR3bGz5+uuvGT16NMXFxYwdO5aOHTtWuP7qeO4olK9DjRo1aNq0qc0TqVatWnTq1Am4NLxxNZaH9ufn59tsr8wD+81q27ZtnDp1ittvv71SY8kWPj4+1otBV3rwwQeBio/zjVZUVMTf//533n33Xdzc3Jg9e3apC1f2/t4tdeWpjudPRcfmSgkJCTz77LMUFBTw7LPPMmLECLvfq7qdOxq+cCDLVDHL19TyWO62Km/ssKIxw+rgq68u/VHMXr16Vfm6LcelvFC6Gc6fP89zzz3Hpk2buO2223j33XfLhI6vry979+4lMzPT5k0Q9v7eq9v5Y8+xuVxsbCzvvfceTk5OvPzyywwdOrTKtsWM5456ytdh9uzZjBkzhtTUVJvtv/zyC1DxFWLL1XfLFeUrpaenl6qrjjZu3AjAww9X/q8Sf/fdd7z44ot8/PHHNtvtPc43SnZ2NoMHD2bTpk00btyYRYsW2Qwdy2wJW793wzDIyMigRo0aFd61Vp3OH3uPDVw6BhMnTuS9996jZs2avPnmm5UO5Op27oBC+bqkpqaSlJTEl19+Wabt9OnTJCcn4+rqWuHYV0BAAE2aNGHPnj0cO3asVFtJSQn/+c9/cHJyMt3UHXudOXOGI0eOUKtWLVq1alXp1+fn5/Pvf/+b+fPn2/zzXytWrADg/vvvv+5tvV6FhYWMGDGC3bt306JFC5YsWVJuGFp+n5bbzS+3fft2srKyCA0NxcvL66rvWV3On8ocG4BXX32VhIQEvLy8+PDDD6/pW1Z1OncsFMrXISIiAoB58+aRkpJiXX7+/HkmTJhAbm4u/fr1K/W1MSsri/T0dH799ddS6/rLX/5CcXExEydOLHWV/O233+bQoUM8/PDDNGvWzMF75BiWOah33XUXLi5XHzGzdXzuv/9+mjRpwtGjR3nttdcoLi62tn322Wd8+eWXNGjQgH79+jlmByohLi6OHTt20LhxYxYsWHDVHlhYWBiBgYEkJyezbNky6/KsrCyio6MBGDZsWKnXVOfzpzLH5ptvvuHjjz/GxcWF999/n7CwsArXX93PHQs9utMOlqed2XpK3Kuvvsq8efNwdnamffv21K1blx9++IEzZ85w7733Mnfu3FIXWGbNmsXs2bMJCwsrdQttYWEhQ4cOJSUlhQYNGtC+fXsOHjzI/v378fPzY+nSpaZ50teVrnZ8AD799FOio6Pp06cPM2fOvOq6yjs+27dvZ/jw4eTl5dGsWTNatmzJkSNHrHdUfvTRR6Vux70Zzp49S3h4OPn5+bRu3fqqdzpaHjS0c+dOhgwZQl5eHm3btsXX15etW7eSnZ1N//79yzzsqbqeP5U9Nv379+enn36iYcOGVw3kv/71r9bhnep87lxOF/qu00svvUTbtm1ZuHAhe/bsoaSkhGbNmhEZGcmQIUNwdXW1az01a9bkww8/5IMPPuCLL75gw4YNNGjQgIiICEaNGmXaQLZHVlYWcH3jdu3bt+fzzz/nvffeIzk5mQ0bNlC3bl0ee+wxRo4cWWpK4s2yc+dO6wWj3bt3X/WKviWUQ0JCiI+PJy4uji1btnDgwAHuuOMOnn/+eZ544gm739vs509ljk1MTIz129WJEydITEwst/aJJ56ocMy9Opw7l1NPWUTERDSmLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWyiIiJKJRFRExEoSwiYiIKZRERE/n/zPcFquFDeXgAAAAASUVORK5CYII=",
            "text/plain": [
              "<Figure size 360x360 with 1 Axes>"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "plt.hist(samples[:,0].cpu().detach().numpy(), bins=30);\n",
        "plt.title(\"Dim=0 distribution\");"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Create p"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 126,
      "metadata": {},
      "outputs": [],
      "source": [
        "calc_p_nll = densities.create_gaussian_nll(parameter_shape)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 127,
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "torch.Size([10240])"
            ]
          },
          "execution_count": 127,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "p_nlls = calc_p_nll(samples)\n",
        "p_nlls.shape"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "HdEt2uU1g4qM"
      },
      "source": [
        "## Optimize"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 128,
      "metadata": {
        "id": "Gc6q2sg2JYeg"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "epoch=0: loss= 591.3127\n",
            "epoch=1: loss= 581.8459\n",
            "epoch=2: loss= 582.8884\n",
            "epoch=3: loss= 577.7267\n",
            "epoch=4: loss= 569.5424\n",
            "epoch=5: loss= 563.8992\n",
            "epoch=6: loss= 558.2569\n",
            "epoch=7: loss= 554.4103\n",
            "epoch=8: loss= 544.3384\n",
            "epoch=9: loss= 548.1510\n",
            "epoch=100: loss= 177.9917\n",
            "epoch=200: loss= 34.2460\n",
            "epoch=300: loss= 4.5056\n",
            "epoch=400: loss= 0.4106\n",
            "epoch=500: loss= 0.0527\n",
            "epoch=600: loss= 0.0097\n",
            "epoch=700: loss=-0.0003\n",
            "epoch=800: loss=-0.0129\n",
            "epoch=900: loss= 0.0169\n",
            "epoch=1000: loss= 0.0097\n",
            "epoch=1100: loss= 0.0074\n",
            "epoch=1200: loss=-0.0042\n",
            "epoch=1300: loss= 0.0014\n",
            "epoch=1400: loss= 0.0052\n",
            "epoch=1500: loss= 0.0152\n",
            "epoch=1600: loss=-0.0080\n",
            "epoch=1700: loss= 0.0076\n",
            "epoch=1800: loss=-0.0056\n",
            "epoch=1900: loss= 0.0077\n",
            "epoch=2000: loss=-0.0014\n",
            "epoch=2100: loss=-0.0047\n",
            "epoch=2200: loss=-0.0041\n",
            "epoch=2300: loss= 0.0223\n",
            "epoch=2400: loss= 0.0060\n",
            "epoch=2500: loss= 0.0097\n",
            "epoch=2600: loss= 0.0025\n",
            "epoch=2700: loss= 0.0108\n",
            "epoch=2800: loss= 0.0070\n",
            "epoch=2900: loss= 0.0125\n",
            "epoch=3000: loss=-0.0089\n",
            "epoch=3100: loss= 0.0030\n",
            "epoch=3200: loss=-0.0016\n",
            "epoch=3300: loss=-0.0031\n",
            "epoch=3400: loss= 0.0016\n",
            "epoch=3500: loss=-0.0085\n",
            "epoch=3600: loss= 0.0076\n",
            "epoch=3700: loss= 0.0110\n",
            "epoch=3800: loss= 0.0008\n",
            "epoch=3900: loss= 0.0980\n",
            "epoch=4000: loss= 0.0012\n",
            "epoch=4100: loss= 0.0030\n",
            "epoch=4200: loss= 0.0003\n",
            "epoch=4300: loss= 0.0152\n",
            "epoch=4400: loss=-0.0032\n",
            "epoch=4500: loss=-0.0120\n",
            "epoch=4600: loss= 0.0147\n",
            "epoch=4700: loss= 0.0061\n",
            "epoch=4800: loss=-0.0088\n",
            "epoch=4900: loss= 0.0133\n",
            "epoch=5000: loss=-0.0141\n",
            "epoch=5100: loss= 0.0153\n",
            "epoch=5200: loss= 0.0004\n",
            "epoch=5300: loss=-0.0047\n",
            "epoch=5400: loss=-0.0030\n",
            "epoch=5500: loss= 0.0290\n",
            "epoch=5600: loss= 0.0290\n",
            "epoch=5700: loss= 0.0172\n",
            "epoch=5800: loss= 0.0044\n",
            "epoch=5900: loss= 0.0127\n",
            "epoch=6000: loss= 0.0203\n",
            "epoch=6100: loss= 0.0082\n",
            "epoch=6200: loss= 0.0014\n",
            "epoch=6300: loss= 0.0231\n",
            "epoch=6400: loss= 0.0074\n",
            "epoch=6500: loss= 0.0132\n",
            "epoch=6600: loss= 0.0337\n",
            "epoch=6700: loss=-0.0120\n",
            "epoch=6800: loss= 0.0180\n",
            "epoch=6900: loss= 0.0165\n",
            "epoch=7000: loss= 0.0129\n",
            "epoch=7100: loss=-0.0085\n",
            "epoch=7200: loss= 0.0532\n",
            "epoch=7300: loss= 0.0048\n",
            "epoch=7400: loss= 0.0066\n",
            "epoch=7500: loss= 0.0001\n",
            "epoch=7600: loss= 0.0198\n",
            "epoch=7700: loss= 0.0111\n",
            "epoch=7800: loss= 0.0279\n",
            "epoch=7900: loss=-0.0007\n",
            "epoch=8000: loss=-0.0047\n",
            "epoch=8100: loss= 0.0216\n",
            "epoch=8200: loss= 0.0051\n",
            "epoch=8300: loss=-0.0011\n",
            "epoch=8400: loss= 0.0348\n",
            "epoch=8500: loss= 0.0238\n",
            "epoch=8600: loss= 0.0135\n",
            "epoch=8700: loss= 0.0103\n",
            "epoch=8800: loss=-0.0141\n",
            "epoch=8900: loss= 0.0191\n",
            "epoch=9000: loss=-0.0027\n",
            "epoch=9100: loss= 0.0106\n",
            "epoch=9200: loss= 0.0041\n",
            "epoch=9300: loss= 0.0039\n",
            "epoch=9400: loss= 0.0367\n",
            "epoch=9500: loss= 0.0009\n",
            "epoch=9600: loss= 0.0065\n",
            "epoch=9700: loss= 0.0318\n",
            "epoch=9800: loss= 0.0255\n",
            "epoch=9900: loss= 0.0060\n"
          ]
        }
      ],
      "source": [
        "optimized_parameters = variational_params.values()\n",
        "optimizer = torch.optim.Adam(optimized_parameters, lr=0.1) \n",
        "n_epochs = 10000\n",
        "n_posterior_samples = 100\n",
        "\n",
        "for e in range(n_epochs):\n",
        "    optimizer.zero_grad()\n",
        "\n",
        "    samples, q_nlls = sampler(n_samples=n_posterior_samples)\n",
        "    p_nlls = calc_p_nll(samples)\n",
        "    KLD = p_nlls-q_nlls\n",
        "\n",
        "    loss_vi = KLD.mean()\n",
        "    loss_vi.backward()\n",
        "    optimizer.step()\n",
        "\n",
        "    if e<10 or e%100==0:\n",
        "        print(f\"epoch={e}: loss={loss_vi: .4f}\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 129,
      "metadata": {},
      "outputs": [
        {
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWUAAAFTCAYAAADhph0LAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XtYVVXi//E3F0UBk7ygYlqAHAwNUxK1yTDHdKZytDKZVDSLaZpSZyzTMhxFykmzMsyyGdIZL3nN9EealmVlZGqYl1Ax8JpXEAUBAYH9+8PvOYkc5KAc3ejn9Tw92F7rLNbC4+cs1l57bxfDMAxERMQUXK91B0RE5DcKZRERE1Eoi4iYiEJZRMREFMoiIiaiUBYRMRGFcg3166+/EhwcbPe/tm3bEh4eziOPPMLbb7/N0aNHK2zH+po9e/Zcxd5XTXZ2Nq+99hr33Xcfbdu25d5772Xs2LEcPny4Wr/P5MmTCQ4O5qWXXipzvDp/Runp6VV+jb3v/9JLLxEcHMzkyZOvuE+OKC0tZd++fWWObdy4keDgYDp16nRV+nCjUChfB9q2bUuHDh1s/7Vu3Rpvb2927tzJzJkzefDBB1m9evW17uZlyc7O5vHHH2fOnDlkZ2djsVgoKCjg448/pm/fvuzevftad9EhGRkZ/OMf/2D06NHXuitVtmPHDvr168fixYuvdVduCO7XugNy5d555x1uueWWcsePHTvGlClTWLlyJaNGjaJ+/fp06dKlTJ1Vq1YB0KJFi6vS16oaN24c6enpRERE8NZbb+Ht7U1hYSETJkxg2bJlPP/88yQmJuLm5ua0PlTHz2j9+vV89tlntGnT5pp8/ysxf/58UlJSys2IQ0NDWbVqFe7uipHqpJnydaxp06a8+eab9OrVi3PnzvHKK69w7ty5MnUCAwMJDAykdu3a16iXFUtPT+fzzz/H09OTKVOm4O3tDYCHhwevvvoqgYGBpKen88UXXzi1H9f6Z3Stv39F6tatS2BgILfeeuu17sp1RaF8nXNxcWH8+PF4eHhw+PBhPv3002vdJYf9v//3/zAMg+7du+Pj41OmzM3NjUceeQT4bSYpcj1QKN8AGjZsyH333QfAunXrypRd6iTS2rVr2bFjB8888wzh4eG0b9+exx9/nO+++w6A/Px83njjDbp3707btm3p3r07b7/9drnZeFRUVIUnJS/+78KTbNu3bwegffv2dsd15513ApCcnFyln8fGjRuJjo6mc+fOdOjQgaeeesr2veyp6ETfxo0bee655+jRowdt27alS5cuPPXUU+U+JLp3787LL78MQEpKCsHBwXTv3h347YTtQw89RFpaGpGRkdxxxx3cc889zJs375Lf32rbtm08+eSTtG/fnrvuuosnn3ySb775ply9ZcuWERwcbPswu9i8efMIDg4mKiqqTN8++eQTAGbNmlXm7+hSJ/oKCwuZPXs2/fr1o3379rRr146HHnqIadOmkZ2dXa6+9T23evVqdu/ezfDhw+ncuTN33HEHDz74IDNnzqSoqMhuv683Wgy6QbRv357Vq1fz448/Ovyar7/+muXLl1OrVi1uu+02Dh06xJYtW3j66ad5//33efPNN/nll19o2bIlzZo14+DBg8ycOZPMzExee+01WzsWi4Xi4mKHvudtt91m+/OBAwcA7K6XA/j5+QGQmZlJXl4eXl5elbY/f/584uLiMAwDX19f/Pz82LRpEwMGDCAkJMShPgIkJiYyevRoSktLadq0KcHBwWRmZvLdd9/x3XffsWPHDsaMGQOcPxFbq1Yt9u/fj6enJ61bt6Zx48Zl2jtz5gxPPfUUOTk5tGrVir179xIYGFhpP3788Ufmzp0LnP85Z2ZmkpSURFJSEiNGjOC5555zeEwX8/DwoEOHDhw4cICTJ0/SpEkTmjdvXubvyJ6srCyeeOIJUlNTcXFxISAggNq1a/PLL7/w/vvvs2LFChISEuyO74cffmDUqFEA+Pv7U7duXdLS0nj77bfZtm0b77///mWPp8YwpEY6dOiQYbFYDIvFYhw6dKjS+p9//rmtflFRke249Vhqaqrt2JgxY2zHhw8fbpw5c8YwDMM4c+aM0bdvX8NisRitW7c27rvvPiMlJcX2ug8//NCwWCzG7bffbpw6deqKx3jnnXcaFovF+Omnn+yW5+bmVulnsGfPHiMkJMQIDg425syZY5SWlhqGYRiZmZlGVFSUra0xY8aUed3FP6OSkhLj7rvvNiwWi7Fy5coydT/55BMjODjYaN26dZk+ffzxx4bFYjEefvjhMvUv/Hvs2bOnkZmZaRiGYZw6dcrWv8r+jvr162ccPXrUMAzDKC0tNebOnWsEBwcbwcHBRnJycqV9sJo7d65hsViMQYMGlTlu/V6vv/56meM//PCDYbFYjPDw8DLHBw8ebFgsFuOhhx4y0tLSbMePHTtmK+vZs6dRUFBgdzxPP/20kZGRYSv73//+Zyvbtm2b3b5fT7R8cYO4cBZp79dHe+rXr8+kSZNsJ9i8vb15/PHHgfP7VidMmFBmdjl48GBq165NSUkJv/zyyxX3uaCgAIA6derYLb/weGFhYaXtzZo1i+LiYvr27UtUVBQuLi7A+eWd+Ph4brrpJof6dfLkSTIzM6lfvz5//OMfy5T17duX/v378+CDD5Kbm+tQe1ZPPvkkDRs2BMDHx8fWv0vx8vLivffeo2nTpsD5cwiDBg3i0UcfxTAMZs+eXaU+XKkff/yRH374AQ8PD2bOnFlmNtykSRNmzJhB06ZN2b9/Px9//HG51/v4+PDOO+/QqFEj27HBgwfTsmVLALZu3er8QVxjCuUbxIXrvI78Ywfo0KGDLZCtrEsG7u7udO7cuUyZu7u77YRcXl7elXQXoNJtbqWlpVVqb/369QD86U9/Klfm4+NDjx49HGrn5ptvpl69emRnZzN27NhyH0ATJ05k6tSptG7dukr9s66RV0WPHj3KLYUAtnXjpKQkSkpKqtzu5fr6668BiIiIoHnz5uXKvb29bX2z1r1QeHi43Q9hf39/gCp/0NVEWlO+QVz4Zq5Xr55Dr2nSpEm5Y7Vq1QLO/+Oyt0XLWm5c8OyEuLg4du7c6dD3jIiI4JlnngHOb7k6d+5chbPgC0/8VDSbtjp79iwZGRkAFa7VBgcHO9RHd3d3RowYwWuvvcayZctYtmwZzZo143e/+x0RERF07dqVunXrOtTWheyFa2Vuv/12u8eDgoKA8x+OGRkZtpm0s+3fv/+S/QJse7WtdS9k7z0Hv/39VvWDuCZSKN8g9u7dC5w/aebofldPT89q+d579uxhy5YtDtW9cM+rj48POTk5nD592m7dC483aNDgku3m5OTY/lzRuBxdvoDzv1Lfeuut/Pe//2XTpk0cPXqUpUuXsnTpUry8vIiOjubZZ591uD04f2Ktqioay4XHz549W+V2L5f1N6RLnXS19s3eb1PWD/WKGDfAg5IUyjeIn376Cbi8X5GvlHV3QFUFBARw8ODBCu9xceTIEeD8DLOymemF+5zz8/Pt/rZgXcN2VEREBBEREZw5c4aNGzfy/fffs27dOo4cOcI777yDl5cXQ4YMqVKbVVVR4F4YeI5+2FRHeFsD91LLDGfOnClTV8rSmvIN4OjRo2zYsAGg3IkpM2vbti1wfh+uPdaTPu3atau0LQ8PD5o1awbArl277Nax/jZRmaKiIvbs2WNrp169evTo0YN//vOffPnllzz88MPA+YtfnO3imwRZWft28803204eWtfoK9rva13euRLWtd+KfsZwfq82YDt5J2UplG8AcXFxlJaW0qpVK9tFJDXB/fffD8AXX3xRbgmjpKTEdlGDvRN39vz+978HYMmSJeXK8vPzHb5p0xdffEHv3r154YUXyv067erqajsBeuH6p6urc/6prV271u6sdMGCBQBl/r7r168PwOHDh8ut0xcXF9u94AR+OzHsyNJBt27dAPjmm2/s/oaTm5vL8uXLAbjnnnsqbe9GpFC+ju3fv58RI0bw5ZdfUqtWLeLi4px6457q1rp1ayIiIsjLy2PEiBGcOnUKOL/9LSYmhvT0dPz9/W3hXZmnnnoKT09P1q5dS3x8vO2ClpycHJ5//nmHZ4rdunXDy8uL9PR0Jk2aVObX/iNHjvDhhx8CcO+999qOW39VP3HiRLVemZaZmcnIkSNta+YlJSXMmDGD1atX4+HhQXR0tK1uaGgobm5u5Ofn8/bbb9t2ZZw5c4ZXXnnF7ok3+G192LpcdCl33XUXnTp1oqioiGeeeabMrUqPHz/Oc889x/Hjx2nRogWPPfbY5Q77uqY15evA3//+9zIn7woLCzlx4oQtZLy9vXnjjTfo0KHDteriZZs4cSIDBgxg48aN3HfffQQEBPDrr7+SnZ1NvXr1ePfddx2ehfr5+TF58mReeOEFZsyYwcKFC2nWrBl79+7l7NmzdOvWze42rYt5eXkxZcoUhg0bxpw5c/j4449p2bIlRUVFHDhwgOLiYtq0acNf/vIX22uCgoJwcXEhIyODXr160bRpU9ts9kr8/ve/Z926dURERODv78/x48fJzMzE3d2df/3rX2V2mjRo0IDBgwcze/ZsZs+ezaeffoqvry/79u2jsLCQv/3tb3avmLPuSlmzZg29e/cmPDyccePGVdinN998kyeffJI9e/bw4IMPEhgYSK1atfjll18oLi6mefPmvPvuu+W2W8p5CuXrwM8//1zm/93d3alXrx5hYWF07dqVyMjISncnmFXTpk35+OOPmTFjBl999RV79uyhXr16PPTQQwwfPrzSS34v1rNnTxYtWsTMmTNJTk4mPT2d1q1b89xzz/Hrr786FMpwfn/wvHnz+N///sfWrVv55ZdfqFOnDiEhIfzxj39k0KBBZT4o/f39efXVV/nggw84evQo586dIzMzs0p9t+e+++4jKiqK+Ph4du7ciYeHB7///e959tlnbWvyFxozZgyBgYEsXLiQtLQ0CgsLCQ8P59lnn6V27dp2Q/nhhx8mLS2NxMRE9u/fX+mJw8aNG7N48WLmz5/PqlWrbOvegYGB/OEPf2DgwIG2pRQpz8W4EfaYiIjUEFpTFhExEYWyiIiJKJRFRExEoSwiYiLafcH5y2t//vlnGjduXKP28YqIuZWUlJCRkUHbtm0rvWmWlUKZ81vKBg4ceK27ISLXqfnz53PXXXc5VFehzG+3TJw/f/5Vu8WhiFz/jh07xsCBA6t0W1aFMr/dqKVp06YVPg9ORORyVWVZVCf6RERMRKEsImIiCmURERNRKIuImMhlh/KyZcsIDg7mxx9/tFu+b98+nn/+eSIiImjXrh29e/dm3rx5FT74MCcnhzfeeINevXoRGhpK9+7def311yt8rExJSQkLFy6kb9++tG/fni5dujBy5MgKn8QgIlITXFYo//TTT8TFxVVYvnv3bvr168fKlSvx8/Oja9euHDt2jLi4OEaPHl2ufm5uLoMGDSIhIQEXFxe6deuGi4sLs2fPJjIy0vZMrwvFxMQwfvx4jh07xj333EPz5s1ZtWoVjzzyiMNPThYRMR2jitasWWO0b9/esFgshsViMTZv3lymvLS01Ojdu7dhsViM5cuX246fPHnSdnz16tVlXhMXF2dYLBYjJibGKCkpMQzDMM6dO2e8+OKLhsViMSZOnFiuDxaLxXj44YeNnJwc2/EFCxYYFovF6N27t1FaWurwmA4dOmRYLBbj0KFDDr9GRKQyl5MtDs+Ujx07xujRoxk+fDilpaU0atTIbr2kpCRSU1MJDw+nT58+tuMNGjRg/PjxQNmnG+fk5LBkyRK8vb0ZM2aM7SkS7u7ujB8/nvr167N06VLy8/Ntr5k1axYAL730UpmnEv/5z3/m7rvvJjU1lY0bNzo6NBER03A4lKdNm8aKFSto27YtixYtIiAgwG699evXA+efzHCxsLAwGjZsSHJysm2tePPmzRQUFNC5c+dyj4fx8vKiS5cuFBQUsHnzZuB8iG/duhUfHx+7ly1av++3337r6NBEREzD4VAOCAhg8uTJLFmyxPbMLnvS0tIAsFgsdsv9/f0pLS21PVDRWj8oKKjC7wuQmpoKQHp6OoZhEBgYaPfZbNb6e/bscWRYIlVWdK6kWuqI2OPwZdZPP/20Q/VOnDgBUOG13tbj1ueTWR/uWVn9kydPOlTf19e3TH2R6la7lhu9X1hxyTqJb/a5ZLlIRap9n7L1cesV3abOety6Rmz9Wrdu3Wqp7+HhUaaeiEhNUu2hbF1ScHFxsVtu/N9zWq1fnV1fRKQmqfZQ9vT0BM7fON6ewsLCMvUcrW+dGVdWv6ioqEx9karSerBcS9V+605fX1927dpFZmYmgYGB5covXhO2rgFb15ivtH5la9oilalszVjrxeJM1T5Ttu6isO6quJBhGOzduxc3NzdbYF+qPmDbpWHd8dGqVStcXV1txy+2d+9eoOLdHyIiZlbtody1a1cAvvzyy3JlW7ZsISsri7CwMNue5I4dO1KnTh02bNhQ7uRcXl4eGzZswNPTk7CwMADbn0+ePMmWLVvKfY+1a9cCEBERUa3jEhG5Gqo9lMPDwwkKCiIpKYnFixfbjmdlZREbGwvA0KFDbcc9PT3p27cv2dnZxMbGUlxcDEBxcTETJ04kJyeHyMjIMheWDBgwAIDY2FiysrJsxxctWsT3339PmzZt6NSpU3UPTUTE6ap9TdnV1ZVJkyYxZMgQxo0bx9KlS/H19WXTpk1kZ2fTv39/unfvXuY1I0eOZOPGjSxfvpzk5GRCQkLYuXMnhw4dIiQkhOHDh5ep/8ADD/D555/z2Wef8Yc//IHw8HCOHz/O9u3buemmm5g8eXJ1D0tE5Kpwyv2UQ0NDWbJkCb169eLAgQMkJSXh5+dHbGwsEyZMKFffx8eHhQsXEhUVRXFxMevWrcPV1ZXo6GjmzJmDl5dXuddMnTqVl19+GV9fX7755huOHz/Ogw8+yJIlSyq8OlBExOwue6Z84U2F7GnVqhXx8fEOt+fj40NMTAwxMTEO1Xd3d+eJJ57giSeecPh7iBSdK6F2LccfYilytelp1nJD0SXSYnZ6HJSIiIkolEVETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQluuK7oUsNZ0uHpHriu6FLDWdZsoiIiaiUBZxgsqWUbTMIhXR8oWIE2gZRS6XZsoiIiaiUBYRMRGFsoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWUTERBTKIiImolAWETERhbKIiIkolEVETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWyiIiJKJRFRExEoSwiYiIKZRERE3F3ZuMrVqxg/vz57Nmzh9LSUvz9/XnkkUcYNGgQbm5uZeru27eP6dOnk5yczOnTp2nZsiWRkZEMGDAAV9fynx05OTl88MEHrF27lqNHj9KoUSN69uzJsGHD8Pb2duawREScxmkz5SlTpjB69Gh27dpFhw4d6NSpEwcPHmTSpEmMGDECwzBsdXfv3k2/fv1YuXIlfn5+dO3alWPHjhEXF8fo0aPLtZ2bm8ugQYNISEjAxcWFbt264eLiwuzZs4mMjOTMmTPOGpaIiFM5ZaacmprKrFmzaNCgAR999BH+/v4AHD9+nMcff5y1a9fy+eef06tXLwzDYPTo0eTm5jJlyhT69OkDQFZWFk888QSJiYncf//99OrVy9b+tGnTSE1NpX///sTGxuLq6kpxcTFjx45lxYoVTJs2jXHjxjljaCIiTuWUmfL333+PYRj86U9/sgUyQJMmTRgwYAAAmzdvBiApKYnU1FTCw8NtgQzQoEEDxo8fD8DcuXNtx3NycliyZAne3t6MGTPGtrTh7u7O+PHjqV+/PkuXLiU/P98ZQxMRcSqnhLKLiwtwfmZ8sVOnTgHg4+MDwPr16wHo0aNHubphYWE0bNiQ5ORkcnNzgfNhXlBQQOfOncutHXt5edGlSxcKCgpsoS8iUpM4JZS7du2Ki4sLq1ev5t///jdZWVnk5OSwdOlS5syZQ/369Xn00UcBSEtLA8Bisdhty9/fn9LSUtLT08vUDwoKsls/ICAAOL+EIiJS0zgllAMDA4mLi6NOnTq8+eabdOnShY4dO/LKK6/Qpk0bFi1aRLNmzQA4ceIEAI0bN7bblvV4ZmYmABkZGQ7VP3nyZPUNSETkKnHa7osOHTrQpUsXPD096dy5M3fffTdeXl7s2LGDjz76yLb74uzZswDUqVPHbjvW49Y1YuvXunXrOlRfRKQmccrui61bt/Lkk0/SvHlzEhMTueWWW4Dza8zDhg1jzpw5eHt78/e//912os66Dn0xa3hbv1a1vohITeKUmfKkSZPIy8vjtddeswUynN998dZbb+Hu7s5///tfzp49i6enJwAFBQV22yosLASw1XO0fkUzaRERM6v2UC4oKGD79u3Uq1eP0NDQcuUtWrTA39+f/Px8Dhw4gK+vL/DbmvHFLl5Drmp9uX4UnSu51l0QcbpqX744c+YMhmGUu4z6Qtayc+fOERQUxDfffENaWhqdOnUqU88wDPbu3YubmxuBgYHAb7surLswLmbdpREcHHzFYxFzqV3Ljd4vrLhkncQ3+1yyXMTsqn2m3LBhQ3x8fDh9+jTbt28vV378+HHS09OpVasWAQEBdO3aFYAvv/yyXN0tW7aQlZVFWFiYbU9yx44dqVOnDhs2bCh3Mi8vL48NGzbg6elJWFhYdQ9NRMTpqj2UXV1d6devHwCvvPJKmQtIsrKyGDVqFOfOnePRRx/Fy8uL8PBwgoKCSEpKYvHixWXqxsbGAjB06FDbcU9PT/r27Ut2djaxsbEUFxcDUFxczMSJE8nJySEyMlI3JRKRGskpuy9GjBjB9u3b2bRpE/fffz8dO3bExcWFbdu2kZOTw5133smYMWOA8yE+adIkhgwZwrhx41i6dCm+vr5s2rSJ7Oxs+vfvT/fu3cu0P3LkSDZu3Mjy5ctJTk4mJCSEnTt3cujQIUJCQhg+fLgzhiUi4nRO2X3h4eHBrFmzGDt2LK1atSI5OZlNmzbRrFkzXnjhBebOnWvbRQEQGhrKkiVL6NWrFwcOHCApKQk/Pz9iY2OZMGFCufZ9fHxYuHAhUVFRFBcXs27dOlxdXYmOjmbOnDl4eXk5Y1giIk7ntPsp16pViyFDhjBkyBCH6rdq1Yr4+HiH2/fx8SEmJoaYmJjL7aKIiOnoySMiIiaiUBYRMRGFsoiIiSiURURMRKEscg04csm4Liu/MTn1adYiYp8uGZeKaKYsImIiCmURERNRKIuImIhCWUTERBTKIiImolAWETERhbKIiIkolEVETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWyiIiJKJRFRExEoSwiYiIKZRERE1Eoi4iYiEJZRMREFMoiIiaiUBYRMRGFsoiIiSiURURMRKEsImIiCmURkyo6V3JF5VIzuV/rDoiIfbVrudH7hRUVlie+2ecq9kauFs2URURMxKkz5cOHDzNjxgy+++47srKyuPnmm+nWrRsjRoygcePGZeru27eP6dOnk5yczOnTp2nZsiWRkZEMGDAAV9fynx05OTl88MEHrF27lqNHj9KoUSN69uzJsGHD8Pb2duawREScxmkz5R07dtCnTx8+/vhj6tevT0REBK6urixevJgBAwaQnZ1tq7t792769evHypUr8fPzo2vXrhw7doy4uDhGjx5dru3c3FwGDRpEQkICLi4udOvWDRcXF2bPnk1kZCRnzpxx1rBERJzKKaFcVFTEqFGjOHPmDDExMSQmJjJjxgw+//xzevXqxcGDB5k+fToAhmEwevRocnNzmTJlCgsWLODdd99lzZo1BAcHk5iYyJo1a8q0P23aNFJTU+nfvz+rVq0iPj6eNWvW0KdPH9LS0pg2bZozhiUi4nROCeVVq1axf/9+evfuTVRUlO24h4cHL7/8Mo0aNWLfvn0AJCUlkZqaSnh4OH36/HbiokGDBowfPx6AuXPn2o7n5OSwZMkSvL29GTNmjG1pw93dnfHjx1O/fn2WLl1Kfn6+M4YmIuJUTgnlzz//HIChQ4eWK2vWrBlJSUl8+OGHAKxfvx6AHj16lKsbFhZGw4YNSU5OJjc3F4DNmzdTUFBA586dy60de3l50aVLFwoKCti8eXO1jklE5Gpwyom+nTt3UqtWLVq3bs3Ro0dJTEzk4MGD+Pj40LNnT0JDQ21109LSALBYLHbb8vf35+TJk6Snp9OuXTtb/aCgILv1AwICAEhNTSUiIqI6hyUi4nTVHspFRUUcPXqUpk2bsnr1al555RXOnj1rK//Pf/7DU089ZTuBd+LECYByuzGsrMczMzMByMjIcKj+yZMnq2E0IiJXV7UvX1iXGbKzsxkzZgw9evRg9erVbN68mbfffhsfHx8+/PBDFi1aBGAL7Dp16thtz3rcukZs/Vq3bl2H6ouI1CTVHsqFhYXA+bDt1KkTU6dOxd/fn5tuuokHHniAyZMnAzBjxgwMw7CdqHNxcbHbnmEYZb5Wtb6ISE1S7aF84Qz28ccfL1ferVs3mjRpwvHjxzlw4ACenp4AFBQU2G3PGvLWeo7Wr2gmLSJiZtUeyvXq1aNWrVoA3HLLLXbr+Pn5AXDq1Cl8fX2B39aML3bxGnJV64uI1CTVHspubm4EBgYCcPz4cbt1rIHaoEED2y4K666KCxmGwd69e8u0ean6AOnp6QAEBwdfwShERK4Np+xTvvfeewFYvXp1ubK9e/dy+PBhfH19adGiBV27dgXgyy+/LFd3y5YtZGVlERYWZtuT3LFjR+rUqcOGDRvKnczLy8tjw4YNeHp6EhYWVt3DEhFxOqeE8p///Gc8PT1Zvnw5iYmJtuPZ2dnExMRQWlrKwIEDcXV1JTw8nKCgIJKSkli8eLGtblZWFrGxsUDZi1A8PT3p27cv2dnZxMbGUlxcDEBxcTETJ04kJyeHyMhI3ZRIRGokp1w80rx5c1577TVefPFFRo0axezZs/H19WXr1q2cOnWKzp0789RTTwHnd1NMmjSJIUOGMG7cOJYuXYqvry+bNm3jiMC2AAAY9klEQVQiOzub/v3707179zLtjxw5ko0bN7J8+XKSk5MJCQlh586dHDp0iJCQEIYPH+6MYYmIOJ3T7hL3wAMPsHTpUnr16sWRI0f4/vvvadiwIS+88AIJCQm2k4EAoaGhLFmyhF69enHgwAGSkpLw8/MjNjaWCRMmlGvbx8eHhQsXEhUVRXFxMevWrcPV1ZXo6GjmzJmDl5eXs4YlIuJUTr2f8u233058fLxDdVu1auVwXTgfzDExMcTExFxu98Rkis6VULuW27Xuhsg1pcdBiWno8UciehyUiIipKJRFRExEoSwiYiIKZRERE1Eoi4iYiEJZRMREFMoiIiaiUBYRMRGFsoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWUTERBTKIiImolAWETERhbJIDVV0ruSKysWc9OQRkRpKT2q5PmmmLFeFZm0ijtFMWa6KymZ1oJmdCGimLCJiKgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYWyiIiJKJRFRExEoSwiYiIKZRERE1Eoi4iYiEJZRMREFMoiIiZy1UL59OnT3HPPPQQHB9st37dvH88//zwRERG0a9eO3r17M2/ePEpLS+3Wz8nJ4Y033qBXr16EhobSvXt3Xn/9dXJzc505DBERp7pqoRwbG0tGRobdst27d9OvXz9WrlyJn58fXbt25dixY8TFxTF69Ohy9XNzcxk0aBAJCQm4uLjQrVs3XFxcmD17NpGRkZw5c8bZwxERcYqrEsqffvopq1atsltmGAajR48mNzeXKVOmsGDBAt59913WrFlDcHAwiYmJrFmzpsxrpk2bRmpqKv3792fVqlXEx8ezZs0a+vTpQ1paGtOmTbsawxIRqXZOD+Xjx48TFxdH+/btcXNzK1eelJREamoq4eHh9Onz203OGzRowPjx4wGYO3eu7XhOTg5LlizB29ubMWPG4Op6fgju7u6MHz+e+vXrs3TpUvLz8508MhGR6uf0UH7llVcoLCxk8uTJdsvXr18PQI8ePcqVhYWF0bBhQ5KTk21rxZs3b6agoIDOnTvj7e1dpr6XlxddunShoKCAzZs3V/NIREScz6mh/NFHH7F+/XpGjRrFrbfeardOWloaABaLxW65v78/paWlpKenl6kfFBRkt35AQAAAqampV9R3EZFrwWmhfPDgQd544w06d+7MwIEDK6x34sQJABo3bmy33Ho8MzMTwHaysLL6J0+evLyOi4hcQ04J5ZKSEkaPHo2Liwv/+te/cHFxqbDu2bNnAahTp47dcutx6xqx9WvdunUdqi8iUpM4JZQTEhL46aefePnll/Hz87t0B/7vRF1FwW0YRpmvVa0vIlKTVHso7969m+nTpxMREcFjjz1WaX1PT08ACgoK7JYXFhaWqedo/Ypm0iIiZuZe3Q2+/fbbnDt3juLiYkaNGlWmzHp1nvX42LFj8fX1ZdeuXWRmZhIYGFiuvYvXkH19fYHf1pgrqy8iUpNUeyhb13KTkpIqrJOYmAjAP/7xD4KCgvjmm29IS0ujU6dOZeoZhsHevXtxc3OzBbZ114V1F8bFrLs0KrqcW0TEzKo9lC+80ONiISEhlJSUlNmu1rVrVxISEvjyyy/L7dLYsmULWVlZhIeH2/Ykd+zYkTp16rBhwwby8/NtyxkAeXl5bNiwAU9PT8LCwqp5ZCIiznfN7xIXHh5OUFAQSUlJLF682HY8KyuL2NhYAIYOHWo77unpSd++fcnOziY2Npbi4mIAiouLmThxIjk5OURGRpa7sEREpCao9plyVbm6ujJp0iSGDBnCuHHjWLp0Kb6+vmzatIns7Gz69+9P9+7dy7xm5MiRbNy4keXLl5OcnExISAg7d+7k0KFDhISEMHz48Gs0GhGRK3PNZ8oAoaGhLFmyhF69enHgwAGSkpLw8/MjNjaWCRMmlKvv4+PDwoULiYqKori4mHXr1uHq6kp0dDRz5szBy8vr6g9CRKQaXNWZ8s6dOyssa9WqFfHx8Q635ePjQ0xMDDExMdXRNRERUzDFTFlERM5TKIuImIhCWUTERBTKIiImolAWETERhbKIiIkolEVETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFrlOFZ0rqZY6cnVd8/spi4hz1K7lRu8XVlyyTuKbfa5Sb8RRmilLtdCMS6R6aKYs1aKyWZlmZCKO0UxZRMREFMoiIiaiUBYRMRGFsoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWUTERBTKIiImolAWETERhbKIiIkolEVETEShLCJiIgplcYjulyxydeh+yuIQ3S9Z5OrQTFlExEQUyiIiJqJQFrmBVXauQOcSrj6nrSmXlJSwYMECPvnkE/bu3UtJSQktWrTggQceIDo6Gg8PjzL19+3bx/Tp00lOTub06dO0bNmSyMhIBgwYgKtr+c+OnJwcPvjgA9auXcvRo0dp1KgRPXv2ZNiwYXh7eztrWCLXFZ0rMB+nzJRLSkp49tlniYuLY+/evbRr147w8HBOnDhBfHw8UVFRnD171lZ/9+7d9OvXj5UrV+Ln50fXrl05duwYcXFxjB49ulz7ubm5DBo0iISEBFxcXOjWrRsuLi7Mnj2byMhIzpw544xhiYg4nVNCecmSJXz99dcEBwezevVq/vvf/5KQkMCaNWto374927Zt47333gPAMAxGjx5Nbm4uU6ZMYcGCBbz77rusWbOG4OBgEhMTWbNmTZn2p02bRmpqKv3792fVqlXEx8ezZs0a+vTpQ1paGtOmTXPGsEREnM4pofzJJ58AMHbsWJo0aWI73qBBAyZMmADAypUrAUhKSiI1NZXw8HD69OlTpu748eMBmDt3ru14Tk4OS5YswdvbmzFjxtiWNtzd3Rk/fjz169dn6dKl5OfnO2NoIiJO5ZRQvvnmmwkICCA0NLRc2W233QbAiRMnAFi/fj0APXr0KFc3LCyMhg0bkpycTG5uLgCbN2+moKCAzp07l1s79vLyokuXLhQUFLB58+bqHJKIyFXhlFCeOXMmn332GZ6enuXKduzYAUDTpk0BSEtLA8Bisdhty9/fn9LSUtLT08vUDwoKsls/ICAAgNTU1CsYgYjItXFVt8QZhkF8fDwAPXv2BH6bMTdu3Njua6zHMzMzAcjIyHCo/smTJ6up1yIiV89VDeW33nqLTZs20ahRI6KjowFsuzDq1Klj9zXW49Y1YuvXunXrOlRfRKQmuWqh/M477/Dvf/+b2rVrM23aNBo0aHC+A/93os7FxcXu6wzDKPO1qvVFRGoSp9+QqLi4mIkTJ7Jo0SI8PDyYPn06HTt2tJVb150LCgrsvr6wsLBMPUfrVzSTFhExM6fOlPPy8njmmWdYtGgRN910Ex9++CERERFl6vj6+gK/rRlf7OI15KrWFxGpSZwWytnZ2URFRbF+/XqaNWvG/Pnzy8yQray7KKy7Ki5kGAZ79+7Fzc2NwMDASusDtl0awcHB1TIOEZGrySmhXFRUxNNPP01KSgqtWrVi4cKFFW5569q1KwBffvllubItW7aQlZVFWFiYbU9yx44dqVOnDhs2bCh3Mi8vL48NGzbg6elJWFhYNY9KRMT5nBLK8fHxbN26lWbNmjF37lzbnmR7wsPDCQoKIikpicWLF9uOZ2VlERsbC8DQoUNtxz09Penbty/Z2dnExsZSXFwM/LZ2nZOTQ2RkpG5KJCI1UrWf6Dt9+rTtsugGDRowadKkCutOnToVV1dXJk2axJAhQxg3bhxLly7F19eXTZs2kZ2dTf/+/enevXuZ140cOZKNGzeyfPlykpOTCQkJYefOnRw6dIiQkBCGDx9e3cMSEbkqqj2Ut2/fbtsZkZKSQkpKSoV1p06dCkBoaChLliwhPj6ejRs38ssvv3Drrbfy/PPP89hjj5V7nY+PDwsXLuTdd99l7dq1rFu3jmbNmhEdHc0zzzyDl5dXdQ9LROSqqPZQvvfeey/rEudWrVrZrvZzhI+PDzExMcTExFT5e4mImJWePCIiYiIKZdEjf0RMxOlX9In5VfZIINBjgUSuFs2URURMRKEsImIiCmURERNRKIuImIhCWUTERBTKIlIhR7ZLaktl9dKWOBGpkLZLXn2aKYuImIhCWUTERBTKIiImolAWETERhbKIiIkolEVETEShLCJiIgplERETUSiLiJiIQllExEQUyiIiJqJQFhExEYXyDUB38RKpOXSXuBtAZXf60l2+RMxDM2URuSKV/Sam39SqRjNlEbki+k2semmmLCJiIgplERETUSiLiJiIQllExEQUyjWczmyLXF+0+6KG09OGxeyKzpVQu5bbFde5USiURcSpNHGoGi1fiIiYSI0P5e+//57BgwfTqVMnOnToQFRUFOvXr7/W3RIRuSw1OpSXLVvG0KFD+emnnwgNDaV9+/b89NNPREdHs2jRomvdPRGRKquxa8onTpxg/Pjx1KtXj48++giLxQLA9u3bGTp0KK+99hrdunWjSZMm17inIlKZyk703UgnAmtsKM+bN4+ioiL++te/2gIZIDQ0lOjoaKZNm8aiRYsYMWLENezlpTnyRis8V4LHDfJmlBuX7p/xmxobytZ14x49epQru//++5k2bRrffvutqUPZ0bPSerPKje5G2lZXI0PZMAzS0tJwdXUlICCgXPltt92Gq6sraWlpGIaBi4vLNeiliFSXG2lbXY0M5ezsbIqKimjQoAG1a9cuV+7u7s7NN9/MyZMnycvLw9vb+5LtlZScvyru2LFj1drPc8Wl1HK/9LnUc/lZlyz/9ddfL1mnsnKztFFT+nk9tVFT+lldbezbf/CS/94c+fdY3ayZYs0YR7gYhmE4q0POcvToUbp160bz5s356quv7Nbp3r07hw8f5ttvv630ZN+PP/7IwIEDndFVERHmz5/PXXfd5VDdGjlTdnWt/NOuKp81bdu2Zf78+TRu3Bg3t5q/JiUi5lBSUkJGRgZt27Z1+DU1MpQ9PT0BKCwsrLCOtaxu3bqVtlenTh2HP8VERKri1ltvrVL9GnnxiLe3N56enpw6dYri4uJy5cXFxZw6dQoPDw9uuumma9BDEZHLUyND2cXFhVatWlFSUsL+/fvLle/bt4/S0tIy+5dFRGqCGhnKAF27dgVg7dq15cqsxyIiIq5qn0RErlSNDeVHHnkEDw8P/vOf//Dzzz/bju/YsYOEhATq1KnDgAEDrmEPRUSqrkZuibOaP38+EydOpFatWnTu3BnDMNi4cSPFxcVMnjyZPn2uj83kInLjqNGhDLBu3ToSEhLYuXMntWvXJjg4mL/97W906dLlWndNRKTKanwoi4hcT2rkPuUbWUFBAbNmzeKzzz7j4MGDuLi4EBAQwMMPP8zAgQMdurAGfrsqsiIdOnRgwYIF1dRr8/j++++ZOXMmqampnDt3jjZt2vD000/bThw74vjx48yYMYOkpCQyMjJo1qwZf/rTn/jLX/5i97L/61VJSQkLFizgk08+Ye/evZSUlNCiRQseeOABoqOj8fDwqLSN4uJi2rdvT1FRkd3yJk2a8O2331Z3101NoVyD5OXlMXjwYH7++Wd8fHzo2LEj586dY9u2bbz66qts2LCB6dOnO3RV4s6dOwEIDg62u3XQ39+/2vt/rS1btoyXX36Z2rVr07lzZ0pLS9m4cSPR0dFMnDiRyMjISts4duwYkZGRHDt2jJCQENq0acOWLVuIj4/nhx9+YNasWdSqVesqjObaKikp4dlnn+Xrr7/G09OTdu3a4e7uzrZt24iPj+ebb77hf//7X6UXb6WlpVFUVETLli1p165duXIfHx9nDcG8DKkxpk6dalgsFmPIkCFGdna27fihQ4eMnj17GhaLxViwYIFDbU2fPt2wWCzGihUrnNVdUzl+/LjRtm1bIywszEhNTbUd37Ztm9GhQwfjjjvuMI4dO1ZpO3/9618Ni8VizJgxw3YsLy/PeOKJJwyLxWJ8+OGHTum/2SxYsMCwWCxG7969y/zcTp48aURGRhoWi8WYOnVqpe0sW7bMsFgsxnvvvefM7tYoNXZL3I3ok08+ASAuLq7MlYq33HILL774IgArV650qC3rTLlNmzbV3Etzsj4U4YknnrD7UITCwsJKHyG2d+9evv76a1q2bMkzzzxjO+7p6clrr72Gm5sb8+bNc9oYzMT6Xhw7dmyZG341aNCACRMmAI69F2+096EjFMo1RF5eHrfddhuhoaG0aNGiXLl1ueHEiRMOtbdr1y48PT2vy2UKeyp7KAJQ6drld999h2EY3HfffeXW7v38/AgJCeHw4cOkpaVVU6/N6+abbyYgIIDQ0NByZbfddhvg2Htx165dgEL5QlpTriG8vLwuOQvbsWMHAE2bNq20rdOnT3PkyBHatGnD7NmzWbFiBQcOHKBevXrcd999DBs27Lp6tqFRTQ9FsIZtUFCQ3fKAgAB27NjBnj17aNWqVfUNwIRmzpxZYZmj70XDMNi1axeNGzfmq6++YtGiRaSnp+Ph4cHdd9/NsGHD7P59Xe80U74OFBUV8f777wPQs2fPSutbZycpKSm8/fbbNGzYkE6dOlFSUsLixYt59NFH2bt3r1P7fDVZH4rg4+NzyYcinD17lry8vArbsc78fH197ZY3btwYgMzMzGrodc1kGAbx8fFA5e/FQ4cOkZubS0ZGBv/85z/x8PCgU6dOeHh4sHLlSvr160dycvLV6LapaKZ8Db3wwgukpKRUWu/+++/nhRdesFtmGAZjx45l//79BAYG8thjj1XannUdLygoiPfff9+2HJKfn8+4ceP49NNPGTVqFMuWLavCaMzr7NmzwKVv41qnTh2ASz6pxtqOtW5FbeTn5192X2u6t956i02bNtGoUSOio6MvWdf6PmzSpAkffPABt99+O3B+m9ybb77JrFmzGDlyJF988YVD2+uuFwrla+jIkSPs27ev0noZGRl2j5eUlBATE0NiYiL169dn+vTpDu2TfeKJJ+jZsydeXl40aNDAdtzT05NXX32VzZs3k5KSwtatW7nzzjsdH5BJVddDEaztVLS8YW3DkbauR++88w7//ve/qV27NtOmTSvz3rKnV69efP3117i6upZZLnN3d+fFF19k48aNpKSksHbtWh588EFnd980FMrX0JVcnJGfn8/zzz/PunXr8PHx4cMPPyQwMNCh17q5udk9WQjnZ5OdO3dmxYoVpKSkXBehXF0PRbC2U1BQcNltXI+Ki4uZOHEiixYtwsPDg+nTp9OxY8dKX+fi4kKzZs3slrm6uhIREUFKSgo///yzQlnMLTMzk6effpqUlBSaNm1KQkJChSefLkejRo2A335dr+kufiiCu3vZt72jD0WwriVXtGZs/Y2mojXn61FeXh5///vfWb9+PTfddBPvvfeeQ4HsCOv7sKIPweuVTvTVMIcPHyYyMpKUlBQsFguLFi2qciC/++67jBgxgtTUVLvlv/76K+DYTo6aoLoeimD9OVe05S09PR3ghnm4QnZ2NlFRUaxfv55mzZoxf/78KgXy/Pnz+cc//sH3339vt/x6ex86SqFcg5w+fZqhQ4fy66+/0rFjRz766KPLesOmpqayZs0aPvvss3JlJ0+eJCkpiVq1atGpU6fq6LYpVMdDEaxtfPXVV5SWlpYpO3LkCLt27aJ58+bX/XY4OL/jx/rbWqtWrVi4cGGVP4wOHTrEZ599ZrsQ5UKFhYWsWbMGgN/97nfV0ueaQqFcg8TGxnLgwAFCQkJISEigXr16lb7m4MGDpKenc+bMGdsx6z0eZs+eXWbLUV5eHmPHjiU3N5d+/frZtnhdD6r6UIQjR46Qnp5OVlaW7ViLFi3o2rUr+/bt45133rEdz8/PJyYmhpKSEoYOHXp1BnSNxcfHs3XrVpo1a8bcuXMrnRzY+3n269cPNzc3EhMTbQEMcO7cOeLi4jh8+DD33ntvlZ4EfT3QrTtriPT0dB588EEMw6Bjx44V/iNo0KABY8eOtf1/9+7dOXz4MP/617945JFHbMdff/11Zs+ejaurKx06dODmm2/mxx9/5NSpU9x1110kJCRcdyesqvJQhKioKDZt2sSwYcMYPny47fihQ4d4/PHHycjIwGKx4O/vz5YtW8jIyODee+/l/fffL7dmfb05ffo0ERERFBQU0KZNm0te4DF16lSg4p/nnDlzmDRpEoZhcMcdd+Dn58e2bds4duwYAQEBzJs3j4YNGzp9TGZyfb97riObNm2ybbXavHlzhfWaN29eJpQr8tJLL9GuXTvmzZvHzp07KS0tpWXLlkRHRzNkyJDr8k5nAwcOxM/Pj4SEBJKTk6lduzYdOnSo0kMRWrRowZIlS4iPj+fbb7/lwIEDtGjRgsGDBzNkyJDrPpABtm/fbjv5lpKScsm99tZQrsjgwYMJCgoiISGB7du3k5qaip+fH8888wxPP/00Xl5e1dr3mkAzZRERE9GasoiIiSiURURMRKEsImIiCmURERNRKIuImIhCWUTERBTKIiImolAWETERhbKIiIkolEVETOT/A7cVU9f6Vj+EAAAAAElFTkSuQmCC",
            "text/plain": [
              "<Figure size 360x360 with 1 Axes>"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "samples, nlls = sampler(n_samples=10240)\n",
        "plt.hist(samples[:,0].cpu().detach().numpy(), bins=30);\n",
        "plt.title(\"Dim=0 distribution\");"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "um_HOFQK-lEz"
      },
      "outputs": [],
      "source": []
    }
  ],
  "metadata": {
    "accelerator": "GPU",
    "colab": {
      "collapsed_sections": [
        "-C2DLYqynn54"
      ],
      "provenance": []
    },
    "gpuClass": "standard",
    "kernelspec": {
      "display_name": "Python 3.8.5 ('base')",
      "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.8.5"
    },
    "vscode": {
      "interpreter": {
        "hash": "72dff520818572b45c401a256658ee25fa99b55260abeb0336ef1b1e43b37267"
      }
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}