Skip to content

Instantly share code, notes, and snippets.

@fehiepsi
Last active April 5, 2019 16:45
Show Gist options
  • Save fehiepsi/13770a6d72c9e0edc37bc71f9b5bd381 to your computer and use it in GitHub Desktop.
Save fehiepsi/13770a6d72c9e0edc37bc71f9b5bd381 to your computer and use it in GitHub Desktop.
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"\n",
"import jax.numpy as np\n",
"from jax import lax, random, jit\n",
"\n",
"import numpyro.distributions as dist\n",
"from numpyro.mcmc import hmc_kernel"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"algo = \"HMC\"\n",
"true_mean, true_std = 1., 2.\n",
"warmup_steps, num_samples = 1000, 8000\n",
"\n",
"def potential_fn(z):\n",
" return 0.5 * np.sum(((z - true_mean) / true_std) ** 2)\n",
"\n",
"def kinetic_fn(r, m_inv):\n",
" return 0.5 * np.sum(m_inv * r ** 2)\n",
"\n",
"init_kernel, sample_kernel = hmc_kernel(potential_fn, kinetic_fn, algo)\n",
"init_samples = np.array(0.)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"HMCState(z=array(1.025145, dtype=float32), z_grad=array(0.00628626, dtype=float32), potential_energy=array(7.903422e-05, dtype=float32), num_steps=array(4, dtype=int32), accept_prob=array(0.86167413, dtype=float32), step_size=array(1.3247993, dtype=float32), inverse_mass_matrix=array([2.7177565], dtype=float32), rng=array([ 80214108, 1281723661], dtype=uint32))"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"hmc_state = init_kernel(init_samples,\n",
" num_steps=10,\n",
" num_warmup_steps=warmup_steps)\n",
"hmc_state"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJztnXd4VGX2x78nCSEk1JCASAtNkCZgkCqCgNIU3bVgl9Uf6tpXlwVR0VUUdde26rosiq4VewEFAUFApYQOAtJCL6H3kuT9/TF3kpnJ3P7eMjfn8zx5MnPnznvP3Pvec8973nPOS0IIMAzDMMEhyWsBGIZhGLmwYmcYhgkYrNgZhmECBit2hmGYgMGKnWEYJmCwYmcYhgkYrNgZhmECBit2hmGYgMGKnWEYJmCkeHHQrKwskZOT48WhGYZhEpZFixbtFUJk6+3niWLPyclBXl6eF4dmGIZJWIhos5H9DLtiiOhtItpDRCsjtmUS0TQiWqf8r2FFWIZhGEYeZnzs7wDoF7NtBIAZQohmAGYo7xmGYRgPMazYhRCzAeyP2TwYwLvK63cBXCFJLoZhGMYidqNiagshdiqvdwGobbM9hmEYxibSwh1FqLC7anF3IhpGRHlElFdQUCDrsAzDMEwMdhX7biKqAwDK/z1qOwohxgkhcoUQudnZutE6DMMwjEXsKvZvANyivL4FwNc222MYhmFsYibc8SMAvwJoTkTbiOg2AGMB9CWidQD6KO8ZhmECybYDxzFrrapjwjcYTlASQlyn8lFvSbIwDMP4mktemo3jp4uQP3ag16JowrViGIZhDHL8dJHXIhiCFTvDMEzAYMXOMAwTMFixMwzDBAxW7AzDMAGDFTvDMEzAYMXOMAwTMFixMwzDBAxW7AzDMAGDFTvDMEzAYMXOMAwTMFixMwzDBAxW7AzDMAGDFTvDMEzAYMXOMAwTMFixMwzDBAxW7AzDMAFDimInogeJaBURrSSij4goTUa7DMMwjHlsK3YiqgvgPgC5QojWAJIBDLHbLsMwDGMNWa6YFACViCgFQDqAHZLaZRiGYUxiW7ELIbYD+AeALQB2AjgkhPjBbrsMwzCMNWS4YmoAGAygEYCzAWQQ0Y1x9htGRHlElFdQUGD3sAzDMIwKMlwxfQBsEkIUCCHOAPgCQNfYnYQQ44QQuUKI3OzsbAmHZRiGYeIhQ7FvAdCZiNKJiAD0BrBaQrsMwzCMBWT42OcD+AzAYgArlDbH2W2XYRiGsUaKjEaEEKMBjJbRFsMwDGMPzjxlGIYJGOVasTcaORl3vJfntRiMyxQXC9z70RIs2XLAa1EYH7Fl33EcPnnGazGkUK4VuxDA1FW7vRaDcZk9R07h22U7cOf7i7wWJbAUFwt8vXQ7ioqFpe8fOHYa039z997s8cJMDHp1rqvHdIpyrdgZhnGGjxduxf0fL8X78zZb+v6w9/Jw+//ysP/YacmSabNl/3FXj+cU5UKxCyHQ6Znp+DRvq9eiMEy5YO/RU1H/zZK/L6RgC4uKpclUnigXir2wWGD34VMY+cUKr0VhfICANfcAwyQK5UKxu8XSrQcxaTnXP0sUCOS1CAzjCFLi2JkQV7z+MwBgUNuzPZaEYYKBkbHVtgPHkVW5ItIqJDsuT6JQriz2wmKBHQdPeC0G4zEigT0xExduwT9/WOu1GI5jZizV/bmZGPYeRzhFUq4UOwB0Hfuj1yIwNpixejdemLpGSluUgJ6Yv32+Av/6cb3XYviO2b9zxdhIyoViT8D7l1Hhtnfz8PrMDV6LwcRw/HQh/jt7I4qVuPVEHhUFgXKh2BmGcZbnp6zFmO9W47uVO6O2s1HlDazYGd8hhIBw0OQLijFZXCxcT+BR48jJQgDA8dNFUttly98arNgdYNPeY4b33XHwBE5IvhkSndHfrEKjkd95LYbveW3menR4ahp2HTrptSil8xWSFHEizn/4iYRX7JOX70TOiMlYuf2Q6j7kci/p9Y9ZhvftOvZH3Pz2fOeEcQEhBHJGTMZL03439b19R09hX5zMxP/9ai0N3ShB0RnTV4dqqew+7L1iT1JOanlN/tpYcBS/bNjrtRglJL5iXxFKCBr0r8Qt3rMw3z9VBtftPoKeL8zEAQtD/FdmrDO1//lPT8f5T083fRy7eKF6Tp4J9qgsnOxlseZXGRLNBXPxP3/C9f8NGWg7D53AVo9rziS8YjeCk/5avyCEwBkJdTXemLUB+fuOY+baPRKk8jduWe4bCo6ixWNT8Pmibar7HDx+Gp/F+by4WOD46cK435HVrb9dtgM5IybbcumEB8Xl4FbTpcuzP+LC52d6KoMUxU5E1YnoMyJaQ0SriaiLjHbd4ojBGszFxQKvz1zves3mPYdP4vZ38zTlfH3mejQb9X1g6kkHid93HQEATNMoQ3vfx0vx8KfLsH7P0ajtL03/HS0fn+rodf1EKY63dvcRy22E3Z2yXDFueE9PnC7CniPeu7GcQJbF/gqAKUKIFgDOg4uLWRup97H1gHa26WClFIAeM9bswQtT1+Kpb38ztL8sXpmxDtNX78ZXS9Xr0HySF7L2rLhQ7OKUlXa6sPxU9tuj+Mljf/OXS7YDAA4dL6vY/TTBGJZFlismzML8/XIbjOC6/87DBWNmONa+l9hW7ERUDUAPAG8BgBDitBDioN12ZXH8dKHuZObGAmNRLOGbbuPeY/h+xU6dvb1l79FTtqwRK8raqqI5U1SMnYfKPnxve3ehtQZ1SETX3N8n/Ya/fbY8apuffkbJpZcs1L0fLcF2h8qALN0aX02Fw23X2RjBeI0Mi70RgAIAE4hoCRGNJ6IMCe1K4dQZ+Vbfos0HcNcHi6X4tJ0i9+npmtaIEAIFR8pGpBjVzQVHTiHfRFinFk9+uwpdnv2xjFU6Z11plMG033YjZ8RkKccL41a0lB1VF9aT037bjYkq6wn4wXJPKnHFyGdvnH7qJBsKjuLDBVvQ96XZ+GW9fyJdzCBDsacA6ADg30KI9gCOARgRuxMRDSOiPCLKKyiQWNfBB53aLLPW7sF/Z2+M2mZ1CTGrjJu9ER3HTLesnDuOmY6eMSOhSGNt1to9yBkx2ZC1NXNNqD8cOaXuRx43214ZgaJigZvemo8Fm5wb2uuhpoCFEFizK3GtQyDCFeNAPy4slmdAfTB/M2au0Q8MWLn9MABg0z45xovbyFDs2wBsE0KEg7E/Q0jRRyGEGCeEyBVC5GZnZ0s4bHxyRkzGG7NKiyTJtGZi27I66rx1wkKM+S56GmLs965NSwAAflKKJqkpXrWfNmn5DkNWzISf8wEA3XxSdG334ZOYs24v7v94ia9cGADwdcTcyZVv/IzJy/3t5isqFmXCN2XkJ63ZdRhzlVFa5NzZrkOnsEzFbWKWUV+uxNB3nHHx+Qnbil0IsQvAViJqrmzqDcDd2cUYnp+yFs8oilPtJj52Kn4ImRlkJmNMWbXLwAHVj1csSosv9Xt5tnVBdB6E93y4BNePj06o8pme9C1qly8/wio8VViMuz9c7JJEwNtzN0W5vNTYuv847ngvDyfPFOHO9xehxWNToj4viYpRfqOVe6Pfy3Nw41tlk/Xu/nCx4QAHGRh98Ash8Mv6vb6cs5EVFXMvgA+IaDmAdgCekdSuLmp6aFyMqyOW2LAyr9m6X8VyFsJQpME2JfJn8oqdpob1svvk+Dna510Lv90favHjbqKmNA6dOIPNEtwEf59kzAZ7atJvmLpqN2atLYgbtlkaFSPif5CghDPb48X4f7V0O64fPx+f5pXmH2jNu81dtxf/9788Vx4EUhS7EGKp4mZpK4S4QgjhaCrlim2H8LPBSQ0n+1Xk9XHqYr0xawM+WrBF9fNNe49FlVMwmmSid16s/p6nJzvjUrK7jJ2RX5MzYjJGfB6KPFm98zBaPj4V3yyTt9Shuo/dfFuDX5uLwye1HzxTVu4sCaMEQkpn0vIdqsrHzhkO6jKDv2zYByB+jH/YGNt6oDTLtNmo71XbGvrOAkz7bTfOFCWIYneby16bixvGz8fB49ox2zeMn4fN+9xJ7d1QED0CWLn9EL5cUjaTcMLPm0peGyn+paXUgVBdGivlFNSUiRM3qBACXy3ZbimKaNuB+Ncv3rk1gt6v+3hhKPJk1Y7Q5Nksn2bg5uv069OFxbjz/cUYMm5eybZmo77HPR8uwSvTzZV+iOTtiP5rld2HT+KJb1ahUKM/uG3ol+1n0TeIDMMtrNDdqKeTkIo9TLu/T9MMWft5/T68MNX6MmKrdhzC77uPYPXOw7r7Fsb4Swb9ay4enLiszH5PRiQ3uVE/ZPfhk3EnSAuVTjZ7XWmE0pdLtpV5QBlBr9N/v3IXHpi4FK9ZWPmn38tz4m6Pd24dwWP3kNXDh5XHtjjJeVoT5jeMn4dnYyb27/5wMX5Q3C9qUUVTlTmiUwaSyh75YgXe+SUfc30SSlhYVIzuz5WWADB7zs3q/D2HnQ/fTGjFboR4ev/g8dO6daPz9x7DwFfn4pKXZqP/KyHl4sVg02yn+XXjvqj3nZ6ZUSYy5XRhMRYoGX2RcxEPTlymmrRhhwPKyGqPSjyylnV21OYk94FjpzHsf3lxMze1CIt0SkKugtMuVa1Rllnr8Of1+/CfmPkpI1E64YdFvESzMIu3HMCOgydKjCAtyYycs7nr9lqub7N8W2k/f8rAPMOtExaWKewVPutl5hV8QIrXAjhNPIu+3d+n6X5vX0xqfryywEIAe46c9JV/0cikcO7T+r/fdxg4xacKi3DkZCGyKlcs2TZuzkb88NtuZFWpqPFNdSYv34nXr7f01bh8snArBratg4yKpbeebLVwpqgYa3Z6Exe/Yrv66PYPb/wCAOhxjpxw53gRNEaJfCDMMBDXDiAqjDqSZdt8k2hfQvAtdkntxPNjCwhcMGYGOo5xpvTs0VOFjqRT6026AU55IEKtfrxgC26Kc1POXLtHfVFiAwLd/cES5OqUAY580O+NUwu+dD/tY+Xl7zedjDN/034M/3w5Ln9trm6i1ISfN+FPFuKtj54qRL+XZxsKD5QVGx6vTS0jdrcvFgbR0wxlP/9oQfzMX4n5U9IIvmJ3MEHJLnpqoaMHtcrDHDGg/MPo/Y7wiKaoWODRr1ZgxBcr4sZOP/71Ktz89gIzYpZwpqi4ZOGJeHw4PzQJHTkfoOYPnvDzJvzlE3Uf/i/r9+KqN3/Fm0o2bGFRMc5/ahq+Xro97v7hfhPOndhQcAzX/OdX9R+D0FzMjwYtyUj++ukybIiofaSlYOOVlJDNxoKjuObNX6PyRozcR36Pknzim1WWZXTjtwVfsXstgA1OJMjiDKO+XBH1/uSZoridd/m2Q3h/XtkoHxkdPVwF0Qh6LtHnp2hPuO9QLM6w2+vIyULsO3Yao79ZZfp4xcVCcwez7lujBe2cJnxNJy3fiQX5+zWToE6eKcIOAyPTU4Xm74f//ZqPpo98V2aCP7LLxT/H2if+nV/yTcviJgmv2PV0QpKDj8dbLFqXehQVC7w9Vz+szC+1pD/Jiw49fG7KmqibRT9m3sBBdNpwIjZbNvHOQ+NHvsOrFqKFzPJORJjivmOncfx0IXZ6tKRe+FoNnbAQB4+fxh3vLUJXA6UnhoybByEEWj0+xXBBuNHfrEJhsSjTx45rGE1GH6Z67pyjpwoxfs5GFBeLmBGL872yHEyeOte27CXtzhQVo6hY4NtlO+JmBL7500bc1CUHQCjCIDwZZYSiYoHkpPgnY+ehE6hTrZLqd8M3Uf7YgYaOtfdo9MSzU0EDP6/fi25NszT3iT20lZvq2KlCpCQTKqYkazeuw0kLlUaNRLUY+UlPRITZzv69AJf9a26UywaQ9xD8Zf1ezZIdkceds25vSd0iPZZsOYhfNuzDMYOLv783b7Nq37vvoyUlr63OY4VHAWrn/5nvVuPD+VvQsGYG0iqU2tBuGBsJb7E7hdn730pkTKzFfdm/5qLFY1NUb4rIDmh2Ym2LxhqMXZ4tay05kUkru4LhDePVoyKEEPj3rA1SFnpuNXoqBr9WOhkZe6XDZ+rg8TPSSwsbwcqlilXqMrl+/Hz8d469RCa1u2mCiQSpN2Y6PxLS4tCJUIjtuj1HsG63uyVMEt5i11fA1p6PboSmxibfhBWfE0O168bNw7xHehve30hWrFH0fo6dn1tcLJAUMxJp+fgUdGhQQzUBRs8Kjje3EflQ2q08kP0YvxyL3yW0krTmF4zep3pzNk6Q8Bb7QZ3EExXvQ1wWbd6vEwcuLxFEC0NRAybbPHqqEDsOnjBsUf5z2u8mj1DKtxLrq4RR+72xyTQAcPx0keGsxoPHT5uu9Bm+UdVccTPX7pG2RKGRZ4dWfykqFtinEdYpg3tsVKO0s85qLA/FRDKZmeeJ+p7vH4f6JLxi18uUNHNB//jvX9HnxZ9Mf08LKzVNnPDBnSoswiUv2SjnC1h2M4SHpFb5bNE21VvtuSlr0PVZc+tWRrrN2v19GlqNnmr4u5EVH9XcVUMnLMStHtX8jqeUjLrArKqzSQ7Uj9+hGuuufnd8vjj6XrMzopIyGPPw+ZDwil13mO9xXMSDE5cpK+To15spwcBTxay75kyRsJ2eb5Wx36/R/FzvGj386TIs2qw+Ua2uBOJj1SLbf+w0Wj5e+hDQamWThZo78ZChG7TmIvRYvydxV3ZSK2Ghh9M6g+PYJeCHBKXPF29Hv5fn4P15m40dx9phbLFW8sSm0UgHo7i9dGA8Yie7SxaV8JGvXZZSOnaqEJv2HkOfF+2N8qzgxPJ6ZvqPgDB0r7/5k73lGp2EFbvCwvzSFO95G/fFjSO3ev9OWRmqfPfoVysN7W9E5v2SfLhhLrWz6pJNnPYBG0FvVBEPI5Z/Xv5+x1dEenXGOmyUNEIIc9Nb89ErZk3bWDqOmW7YILjz/UWGj52rUaLDqnH13BTz11ePcHa22sjAS1994BX7dysMLDkH4Oo3S1O8h4ybF9dv+MF8YxZ3LFqp7gwMxyXLQC1m2Yr1FX7Qa7nFRn6xQvUzWfzw225c/M+fbLvaIiOhFm/RryNTcOQUnlDJtrWDHaNluUpBrh/irPokC7WAC6O6xwmkKXYiSiaiJUQ0SVabfkMrLVqmX84HXgcAwPyN+5C/VzveudFI9+O27SLLcxJu5qIXZqrus87FJRjthqiasarDxJaJ9prfJcWL+8i7ZgmZcez3A1gNoKrENsslzzswbLTCtRGr76iR6DeAGWJ/a/h9vIJpsnIR3CjUlUho1XuXhZf182UhxWInonoABgIYL6M9NV6Zvs6TzD63MVNZkXGPxVtiI3P89VTbcfBEIGKwtVipUe8d8FdtIC+RZbG/DGA4gCpqOxDRMADDAKBBgwaWDjJ+TtlkFIYxi1XVF5slGVsTJ+oYHgxljNRgDzqFxcV4YeoavD7TvxErbjx8bSt2IhoEYI8QYhER9VTbTwgxDsA4AMjNzQ22WcEEkp2SF4hYJzHrkgnxt8/lTFYnuotRhsXeDcDlRDQAQBqAqkT0vhDiRglt66JXUsAt3BoCj/7aWMgko47ZEgJWMOJj72szE5gxTqGHCXqxuOFjt63YhRAjAYwEAMVif9gtpe4nDhxz5wHz7q/WQi6ZUu54z3z0h1nOFBWXi/mgRGH7wRNobbB0xPo9RzExL/4yeIlC4OPY3cLK6i6MN8hcR/YulRDB4y7G5jNyuTeiVrsTuFFSQGrZXiHELACzZLZpBjtV5uzixqoojDpuuFfi8f1K75JQGEaNQFnsTlSZMwqrdW/ZvE99IRGGKW8ESrF7yYMTl3otQrnGD0XCGMYvsGKXxD7JRbkYcyR6eBrDyIQVO8MwTMBgxc4EgqCn0jOMGVixM4GAXTEMUwordiYQrJe80ATDJDIJpdiP+CQlmPEfSwwsDMEw5YWEUuwMwzCJjhs5L6zYGYZh3MQFzc6KnWEYJmCwYmcYhgkYrNgZhmECBit2hmEYN3Eh54IVO8MwTMBgxc4wDOMmHBXDMAzDmMW2Yiei+kQ0k4h+I6JVRHS/DMEYhmEYa8hYGq8QwENCiMVEVAXAIiKaJoT4TULbDMMwjElsW+xCiJ1CiMXK6yMAVgOoa7ddhmEYxhpSfexElAOgPYD5cT4bRkR5RJRXUFAg87AMwzAJA7kweypNsRNRZQCfA3hACHE49nMhxDghRK4QIjc7O1vWYRmGYZgYpCh2IqqAkFL/QAjxhYw2GYZhgogbq33JiIohAG8BWC2EeNG+SAzDMIwdZFjs3QDcBOBiIlqq/A2Q0C7DMAxjAdvhjkKIuXCndjzDMAxjAM48ZRiGcZHDJ844fgxW7AzDMC7yzbKdjh+DFTvDMIybiASIimEYhmGM40I5dlbsDMMwQYMVO8MwjIu44Ilhxc4wDOMmCZF5yjAMwxiHLXaGYZiAwZOnDMMwAaOYwx0ZhmEYs7BiZxiGCRis2BmGYdyEJ08ZhmGCBU+eMgzDBAzBk6cMwzDBImHi2ImoHxGtJaL1RDRCRpsMwzBBJCnJ+XWJZKx5mgzgdQD9AbQEcB0RtbTbLsMwTBC5sFmW48eQYbFfAGC9EGKjEOI0gI8BDJbQLsMwDGMBGYq9LoCtEe+3KdsYhmEYD3Bt8pSIhhFRHhHlFRQUuHVYhmEYX0FIAB87gO0A6ke8r6dsi0IIMU4IkSuEyM3OzpZwWIZhGCYeMhT7QgDNiKgREaUCGALgGwntMgzDBA5y3mBHit0GhBCFRHQPgKkAkgG8LYRYZVsyhmGYAOKCXrev2AFACPEdgO9ktMUwDMPYgzNPGYZhAgYrdoZhmIDBip1hGCZgsGJnGIYJGKzYGYZh3MSFsBhW7AzDMAGDFTvDMEzAYMXOMAzjIhWSnVe7rNgZhmFchBU7wzBMwKiQnBjVHRmGYRiDnFO7iuPHYMXOMIzjZFWu6LUIvoFdMQzDBIJxN5/vtQjlClbsDMMwAYMVO8M4TNt61eJuf/6qti5L4h1u1CBnSmHFzjAecU1uff2dXOSfV5+H3IY1HGn73DpVHWmXiQ8rdoZhAADJSYT6memOtc24hy3FTkQvENEaIlpORF8SUXVZgjEMw/iFV69r77UIprBrsU8D0FoI0RbA7wBG2hcp8Rk14FyvRWAYX5Eo9voFjTLjbm9UM8NlSexhS7ELIX4QQhQqb+cBqGdfJHWG9WjsZPPSaFIrsTpBELgm19GuF8Vzf2yD+3s3c+14TtPn3FqOH4MoMVT7e7ddEHe7gHBZEnvI9LH/CcD3Etsrw70XN8WtXXN097u9eyMnxWB8yKODWrp2rGs7NvBEsWdmpDp+jMHtznak3cRQ68FBV7ET0XQiWhnnb3DEPqMAFAL4QKOdYUSUR0R5BQUFloStklYBT1zeCh1ztGfuRw301hUiEuvhXsLV59dDlYopXosRCIb3a17yWpZSS3FhArJnc+etd8Z5dBW7EKKPEKJ1nL+vAYCIbgUwCMANQqirNCHEOCFErhAiNzs725bQ113QQPNzIkK9GpVsHaM88sLV56Fdg/jz306FwSUqes/u27vLdxsmiDcjLokseyJiyzwjon4AhgO4SAhxXI5IcvCyIyWqxa5Fks/vTD+fc1mipaY4E50cDkV08hInio9dDT/3r3jY7SmvAagCYBoRLSWiNyXIpIuRk6y1T05NZ2J1E4Fv7+mu+bnqDZjY96V0kgi4vlMD1Eiv4MrxujSuibF/cCZT9dGBLTG0Ww4GtKkjpb02deNn2jLuYTcqpqkQor4Qop3yd6cswazStFZl3X06NPCHW6F2Vfcr3mVX8a7K3sfDOjvXuMsWFRHhmSvblCjDRlnqkVAynonD+zVHzcrOTJ7WyEjF6MtaSas6mKNxLvyOEaMx1YXqjHbxv4QmiXcTzXjoIkPfTU1OwjNXtpErkM/wKmyrTrU0dG5c07H2iz0eK2c5pHQjqZIWPTp4ZUg7Ke1qTI0xcbi4hf8nmBNSsZvthk2yo614te93bFQD13fSnpg1guzbZECbs6S15dU9/NwfnS145dbP+vLPXaPeN1NGiPVrmHPvzRney/Sx61aPDggY3K4unrjMfpin7HMX9AfF6zd0wPInLimz3Y2oJaMkpGK3QoOIGhixN0gYSiBH8rt/ip9IEcufezYx1a7aTWn3zMhcXGDDMwPQv3X0w07PYm8goQZKlYopaB/jxrulaw4+v6srLmoeHemlN1dotnZK7NzHwLYhF5AMFSpbD7uh1q9sX1daW5FL1aldt8gaOslJhKppZedWLm1d1gDrc25t+wJaICEVe7v66pMzahdm9vBeWD+mP9Y81Q/392mGxxxMaNGzWHo1txfu2bVJTVx0Ttk2bujUoEwC1/B+LQy1+eHtnXT3+fcNHfAvizUzZBozyUmEC5tF/3495XRFhCL4S99zNPcd1Nb4JCIR4Xy9UFClU151vvHs2GqVtCdljUQpVUkzGPRmUxO7NYEcicxw5ll/1R89aSWHad0T42/JtSSTXRJSsTetZW3NwJTkJKRVSEaF5CRc2qrsk/Q8jQeGHtdqlGCd+kCPqPdpFZJ123vzxtIVZ2KTRtSUWIPMdNwXkREZtmoih4hq93DXplkAoi3DyAxeIqB/mzq47DxrmYmxlq5sZM4d1KmWpnIMe0SqYj29bDQ6UGs3o25FtXN3hcEs1NpV458vJ5GZhas2gjeNjzxQCanYtQi7U14ZYtyyvLhFLXx//4X4S9/mcT/X84d2bVITg9uX3gR61zd802pZmc3PKn14XZNbH8tGl/XpxW074nXYqpv1156lsnnk/5Qdgx2r+OL9rLQKpccM7+5JuQlFuEgr26rbb7TiUzdyHY0eQ60prRK+Mx/uiepqlrrDXWzCrR3RXTFEAGDSvdohvE7jxxD9wCn2ZrVDk1m6w+MYzq1TVdXvGdnB7724acnrW7o0tCChNfSG5lrU05nY69bUuWgVtzD6vEr3sGRCpAIIvz5bZXRQ5rvK/5oxi0LH+9lm8zSs6OFGWRkl8xax/n+r/u+3DLoterWoFXUuvY6ICuOnQmGBUuz39W5meLkxs31hzJWtcV/vZrivd8g/v25Mf9ytKPnYIW9k22aVZqSVaYV4ykOP/97snh/QjJ/ZDLV04vPNWFVqSVpaTcR+J6ppeFesAAAUX0lEQVR/aRz8l5G9jQsGdUs90p3QSkkQUlM0sasZyRzFTRjaETlZ1iaqrYYRFhV7q1D1Rkb9WsmLajNKwiv29NRSf3XnxplITzVvkRm552/o1BB/6XsOKiQn4bbujVAhOQm1qqQhf+xADGpb6oapHGMRNqyZUeYG0+oINdJTdWXSUlJ6/vt497CRcyYrYugfV58npZ1YkuKMtiJ/azjktUm2+eSZRwYYm4DWw86QPTbpJ/ZB0rdl6ZyR3mFifcpqalGvHdmGMhHhpWvN9w8zer3FWfrzcx1U6iXpoXY+Kto01qyQ8Io9akJKkvJZ8lhfwz7tMGGFfEX7s6PCp6yi1VfVOpBASLGHF/qIlOKnv/bE3L+Zj512EyddQoPa1sE393TD5RYmf4foFJ0zDsV5pc2EoR2RP3agqivOjHK9o0djPNT3HLwYozzV2mhk4SFotudfd0F00MGV7eshI1U/uCAS2QX/PrjdXIY0+9gdQEaKfGy/rpGRatqnfW6dqvjw/zrhsUEt0bN5Lc1MxKqVQhby1TqLQ1itaRNuP5KGNTN0fe2ANyVhwhOrGQZGDvljBwIwby0SEdrWq26oGJVT54CirRBDnFfPvPWodmqSkgj39m6GqmkV8NmdXUqs1/QYRXqWEuXS4iztBajDv2dIx1LlbNaIf9Zm/Zsm2RmoXTUNNxuY7/rAQEgvAFRKTUbXJiEjY2i3HMMjNiN98lGXSoonvGKvU039aT1xWGdMdLI+SQxdm2ShYkoykpMI9/dRj5WuVCEFa5/uh4dUonCA0H0/5YEeWPXkpVJllDl8HnNla919kpNIN6po7B9CZRy0Ime6N83ChKEdNdvZ+MyAqPd+saTiiaE3ulT7tHFWyKWktyaBHrk5mZh0b3cserQPMmLch0O75QAoVfB6tKtfPW5ehSyevLxViXtkvMp8kJEQ4m5Ns3QT1cL3x3u3dcLap/th9GWtMKyHdpJf5LU6p3bZWlWRn99+oTurwCWsYv/iz13x5o0dNG/eTo1ropOD9UmMQCjrdweAiinJuoonrUJymZtO71h6WJ25jydrL41FGf6gREacU7tKmbC5sOIAgKzKFZFiICv1mo71NY8HlPWzm32ItQ/7Vi08ELS+EhYjUjyrD5029aphzvBeuKmzuoXaonbIEm9ZR9viTklOKhNlAwB3XNQE+WMHokZGKjY+MwD3RUSCGSP6x9mN0qmUmlwSsVbVRnQYYHyEn5xEqJhiziWkdm/Fm/9xmoRV7B0a1EC/1tEZglZvFqdPe70a6Zg4rDMevkQ74xEwWJJYpQPFbpVpscabYNVyV/1RiX6Jl5U4+rJWJSnxoyNqncT7VWayQO1wX+9m+M+N56NJdgZu7GQ+jNXIMyTSSk+OuDgj+pubnK2fmV7iUhrYtg6yq1TETV0aorEywdqrRS1M/8tFGNwuOuwwMvbbKElJhP/r0TjKF65XDqFJdkbJSlzXXdAgbmbnc39sgy9i6u6Y4WxlAlhrZKxHnWppJUsBhis2Wrln9L6Tqcy/RRo0TpOwij3R6NS4puZwsWNODfQ5txYqpDi76IEVV0z19Apxw0i1RhN6x2lfP2QdN8hMt/RglX1+MlKTUatqGmY81BP1M9N169ZbIayM6lavhBoRmZOxPu5I9EIRa1dNw8JRfdAku3JU/4otX71uTH90s6DYgVBVyUhfePjUV1KOF1vegIjwN+VhpXadru3YwFb57PTUFOSPHWh6MjxSno+HdTaVyGiVcG5NFxe9B+V2gUvDdTRc4pM7uoCI0OP5mbbbsuNHjzfUvbtn06gU7qzKqdh79LT1gwC4rXsj9GpRC02yK2PLfv3Ft2TGWvdvfRZenPZ71LZYBdSmXtnyEloSqD1nKiSX2umdGmfi87u6ljzUtJC14lDTWpXRIDNdahG2MK9e1x4TF25F67plXT6yrtb5DWvg0lZn4dO8rVLas9ONnrisJd6bt1lzn3hzJ9fk1keLs6riPAPXXRaBstjN3ArV01Px5OWtNPe586ImuLCZNSuniTIsbnW2sfoz4Rs5fPG1YstjO4/aGrDxOplev376Cv0J0ZkP98TCUX0099GvhUIRseXqFo0TS6o1q12lJOrhxWvOwwWNMnHV+eq1fsyiVXvk/IY1HPO5xjtV397TParukB1irePaVdNwX+9mmtfI7i/9/K6utrKuZXJrt0aY8VBPQ/s+1PccfH13NwChPuymUgckKXYieoiIBBFZ04IeoVf8Z0T/FnjvNmMhUrF0bZqFaQ/2KBOnq8cLV7XF13d3MxXGWaZuisa+epavkZuoSloFqSsxtTy7KhY80hs3SKiFH4kRGWtXTcMnd3SJW1TqnaEdo/zfRpXU01e0NjVJrRV9YtfyrZSaLK1Ozz+vMZA8pAjcVHlYR4ZrLnvcXG5I5IpbNyqTxVYSzIBSt5FMLsjJjFLYkbdW31a1XVfmkdj2RxBRfQCXANhiXxzbspjav50SBTG0mzOFoZrVNl+FMq1Csm6HMKo03Aj3m/HQRUhPTUaXZ3+01U4tFeVm5ycYCYHTomfzWujZvBYOnzxj+DsD29ZBrxa1cKqwCICxpLm+LWvjvdsuwE1vLbAsq9fE9rUuTWpi5sM9oyJiqqVXwMyHe6KwqFizrbCCjFw7dXC7umUmg6O/o31PLHw0NMKUmSz7yZ1dAABTVu6S2KocZDzKXwIwHL4qWmmMrMoVkT92ILpbdLe4jdq8QGyfDsfqNtdIn9Zao9MMTbIrx80lCGcDWo1vfmpwKzS0uei4J5UsLRySqGx9+ZLPbIojk9JCZMZK5jbKyihjbDXKyjBs8FgxTC5uUQu/jry4zPZ4IceyF9bxkwK0pdiJaDCA7UKIZQb2HUZEeUSUV1BQYOewUUTO9PslIcUoZkcY/7lJ21ca1mPdmmbh23u6l1l0A3Cv8zWsmYGFo/pgWA9rCRk3dcnBT3/thSGKK8tstc5I2lus/WEGmX3P7hDeqedZSnISnr+qLT6703qYYjxmPHQRXrW4gEssnRtnRhka9/QyG4Nvnshr75cKj7quGCKaDiBeebJRAB5ByA2jixBiHIBxAJCbmyvt1991URN8OH8Lth884dl6nkYJKyerE7JqxFMq8aI6gNKb3o1noAw/fNcmWSVlBKyw9ul+UTHjbhG2BuOVd1DjqcGt8NjXqwwtjOIV12gsKGOVJtmVo9YlbpCZjrW7jxhaJSpM7L3fuXEmBrSpg5u75ODhS+NneMtWwpEyeL3Mpm6vE0LEDX8gojYAGgFYplie9QAsJqILhBCuOZ2Skgh1qqVh+8ETbh3SMu0b1MCap/rZ9v2qYew+cFGzO8ilrc7CW3M3Yf2eo5r7mc0elEVqShL+PrgVejTLxgMTlxr6TmQ5XTuXJ9FGrrG8f3snLNp8wNJ9ElaoHw/rEvdzJ4w/P55uy64YIcQKIUQtIUSOECIHwDYAHdxU6mHChfYdCNWVjh2l3qZuNWRVTsWDKtl2Rjpt2P2TJlnh3d69UZmkGCfJzEjF9L9c5NrxrHBzl5yocrt+H1H6hewqFdEvzsLQ/kdEvPL2YieAKtQnXI/ZiZhnP1ElrQLyHu1bpv6NmV/dOCsDD/Y5B+Nu1o9t/nOvphjQ5ixcayBk89FBLX2naNvVd3ad1XjEu6GNdsvw2qGR8fyJ8jBIFDnVuERZA9mK2y5S73jtggkjLf1Ssdo9IRz9YMYnV14hItzfp5n+jghZxW/cICe5xW0m39c9ym/rNDJu6PqZ6ZgzvBfOrl4JHy+0Hj2c6ErWCle0r4vxczfh4nP1VmGKf3JevKYdHhvU0lBBOtWWfXTeA2Wxe1BEzRfw86wsrc6upun2ulapId7MhPtIM/FL0tC7fma6bpEtNbzsB173wdZ1qyF/7EDLD/PUlKSSEZNZ/Hj7BUSxl2+L/ZYuOTi3TlX88XxriwiXRwa3q4v8sQNVE6NkcodSz7upiyMIxj3OUeLyB7RxpxKpEfxVCcsipT52b+XwijrVK+H7+y80/b0v/twVuw6ddECiYKK5Dq3Gp/1an2UrZJPxNw1qpmP9mP5ISU7CuNkbvRYHQEAUO/vYraFWNrVKWgqOnCx0WZoQowac60pCkVP4yc/KuIcd37wTBEKxl3dXjIzFsyNZ8EifknPqNv9nMVPVa8pp10tojKwBnKj46zFjkfI+eSo7CadSqrkl+RjnLHUzk7IP9A7lNzSyWAGxPBB5naxOUicCgVDsb954Pm7o1MDV8DaGcRbzSqdPy9rIHzswbsErt/A6McdrbuseqhR7tk5JcKcJhFnWtFZljLmyjddiWOavKrUsEpGsyhWx9+gpr8VwnfLuiinnP7+EazrWxzUd5dfTMUsgFHsiE7RoiTnDe6GwWLveNsMwzsKKnZFKpdRkAN4U3vID0n3t5duzIZ3yErUUCB87U37x26Lksijvrh3GHsG8K5hyw8yHe2Lf0dOOHiNcH71hVnDD45hgwYo9gZkzvBe2HfB/HXonyapcEVmV5S2sHY/01BS8dUsu2mmsbMQGdmIQjtp5dOC5HkviLKzYE5j6memon8lWpBv0Pre2q8d7Z2hHfDB/i+MPLVmEcykSJUkwqC68MMH+dQzjEuEsRrvrlYZpdXY1PJNAIbxjrmyNRtkZqotyM+7Cip1hJNCmXjVMe7BHuU2Sq1m5Iv7Wr4XXYugSXuUrcsHrIMKKnQkUL1zVFg1repNS30wp38r4l9u7N0b7BjXQMSfTa1EcxbZiJ6J7AdwNoAjAZCHEcNtSMYxFrs71PuuP8S9JSRR4pQ7YVOxE1AvAYADnCSFOEZHeulQMwzCMw9hNULoLwFghxCkAEELssS8SwzAMYwe7iv0cABcS0Xwi+omIOqrtSETDiCiPiPIKCgpsHpZhGIZRQ9cVQ0TTAZwV56NRyvczAXQG0BHAJ0TUWIiyFRmEEOMAjAOA3NzcclKxgWEYxn10FbsQoo/aZ0R0F4AvFEW+gIiKAWQBYJOcYRjGI+y6Yr4C0AsAiOgcAKkA9toVimEYhrGO3XDHtwG8TUQrAZwGcEs8NwzDMAzjHrYUuxDiNIAbJcnCMAzDSIC8MLCJqADAZotfz4I/3T0slzlYLnOwXObxq2x25GoohNAtyOOJYrcDEeUJIXK9liMWlsscLJc5WC7z+FU2N+TiFZQYhmECBit2hmGYgJGIin2c1wKowHKZg+UyB8tlHr/K5rhcCedjZxiGYbRJRIudYRiG0SChFDsR9SOitUS0nohGuHC8t4loj5KAFd6WSUTTiGid8r+Gsp2I6FVFtuVE1CHiO7co+68joltsylSfiGYS0W9EtIqI7veDXEp7aUS0gIiWKbI9qWxvpBSKW09EE4koVdleUXm/Xvk8J6Ktkcr2tUR0qQTZkoloCRFN8otMSpv5RLSCiJYSUZ6yzQ/XsjoRfUZEa4hoNRF18VouImqunKfw32EiesBruZT2HlT6/Eoi+ki5F7zrY0KIhPgDkAxgA4DGCJUuWAagpcPH7AGgA4CVEdueBzBCeT0CwHPK6wEAvkdowfrOAOYr2zMBbFT+11Be17AhUx0AHZTXVQD8DqCl13IpbRKAysrrCgDmK8f8BMAQZfubAO5SXv8ZwJvK6yEAJiqvWyrXtyKARsp1T7Yp218AfAhgkvLec5mUdvMBZMVs88O1fBfA7crrVADV/SBXhHzJAHYBaOi1XADqAtgEoFJE37rVyz5m+wS79QegC4CpEe9HAhjpwnFzEK3Y1wKoo7yuA2Ct8vo/AK6L3Q/AdQD+E7E9aj8J8n0NoK8P5UoHsBhAJ4SSMVJiryOAqQC6KK9TlP0o9tpG7mdRlnoAZgC4GMAk5RieyhTRTj7KKnZPryWAaggpKvKTXDGyXALgZz/IhZBi34rQgyJF6WOXetnHEskVEz55YbYp29ymthBip/J6F4Dayms1+RyTWxnCtUfIMvaFXIrLYymAPQCmIWR1HBRCFMY5TokMyueHANR0QLaXAQwHUKy8r+kDmcIIAD8Q0SIiGqZs8/paNkKoQusExX01nogyfCBXJEMAfKS89lQuIcR2AP8AsAXAToT6zCJ42McSSbH7DhF6rHoSVkRElQF8DuABIcRhv8glhCgSQrRDyEq+AICnS9cT0SAAe4QQi7yUQ4PuQogOAPoDuJuIekR+6NG1TEHIBflvIUR7AMcQcnF4LRcAQPFVXw7g09jPvJBL8ekPRuiBeDaADAD93JQhlkRS7NsBRK5UXE/Z5ja7iagOACj/w8sBqsknXW4iqoCQUv9ACPGFX+SKRAhxEMBMhIag1YkoXHAu8jglMiifVwOwT7Js3QBcTkT5AD5GyB3ziscylaBYexChZSW/ROhh6PW13AZgmxBivvL+M4QUvddyhekPYLEQYrfy3mu5+gDYJIQoEEKcAfAFQv3Osz6WSIp9IYBmykxzKkJDsW88kOMbAOFZ9FsQ8nGHt9+szMR3BnBIGR5OBXAJEdVQnuyXKNssQUQE4C0Aq4UQL/pFLkW2bCKqrryuhJDvfzVCCv4qFdnCMl8F4EfF4voGwBAleqARgGYAFliRSQgxUghRTwiRg1Cf+VEIcYOXMoUhogwiqhJ+jdA1WAmPr6UQYheArUTUXNnUG8BvXssVwXUodcOEj++lXFsAdCaidOX+DJ8v7/qYjIkMt/4QmuX+HSG/7SgXjvcRQj6zMwhZMbch5AubAWAdgOkAMpV9CcDrimwrAORGtPMnAOuVv6E2ZeqO0FBzOYClyt8Ar+VS2msLYIki20oAjyvbGysddD1Cw+eKyvY05f165fPGEW2NUmReC6C/pOvZE6VRMZ7LpMiwTPlbFe7TPrmW7QDkKdfyK4SiR/wgVwZC1m21iG1+kOtJAGuUfv8eQpEtnvUxzjxlGIYJGInkimEYhmEMwIqdYRgmYLBiZxiGCRis2BmGYQIGK3aGYZiAwYqdYRgmYLBiZxiGCRis2BmGYQLG/wPgbTDyx3wefgAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"hmc_states = lax.scan(lambda state, i: sample_kernel(state), hmc_state, np.arange(num_samples))\n",
"plt.plot(hmc_states.z);"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"HMCState(z=array(0.5759466, dtype=float32), z_grad=array(-0.10601334, dtype=float32), potential_energy=array(0.02247766, dtype=float32), num_steps=array(3, dtype=int32), accept_prob=array(0.316412, dtype=float32), step_size=array(1.2542776, dtype=float32), inverse_mass_matrix=array([3.831791], dtype=float32), rng=array([ 80214108, 1281723661], dtype=uint32))"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"hmc_state = init_kernel(init_samples,\n",
" num_warmup_steps=warmup_steps)\n",
"hmc_state"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJztnXecFdXZx3/PdnYpUpbeQUAEEUSkiSJKBA0mMTFqYomJqDGWxPc1EI0lMWrUxJLkVbGGaFAjGgsqggqC0pEO0ntbet1ly3n/uDN3586dcs7MmXIv5+vHD3vnzpx57syZZ57znOc8DzHGoFAoFIrsISdqARQKhUIhF6XYFQqFIstQil2hUCiyDKXYFQqFIstQil2hUCiyDKXYFQqFIstQil2hUCiyDCmKnYh+TUTLiWgZEU0goiIZ7SoUCoVCHN+KnYhaAbgdQF/GWA8AuQCu9NuuQqFQKLyRJ7GdOkRUCaAYwHannZs0acLat28v6dQKhUJxcrBgwYI9jLFSt/18K3bG2DYiegLAZgDHAXzKGPvUvB8RjQYwGgDatm2L+fPn+z21QqFQnFQQ0Sae/WS4YhoCuAxABwAtAZQQ0U/N+zHGxjHG+jLG+paWur5wFAqFQuERGZOnFwLYwBgrY4xVAngHwEAJ7SoUCoXCAzIU+2YA/YmomIgIwDAAKyW0q1AoFAoP+FbsjLE5AN4GsBDAUq3NcX7bVSgUCoU3pETFMMbuB3C/jLYUCoVC4Q+18lShUCiyDKXYFQqFIstQil2hyEDW7DqMT5btjFoMRUxRil2hyEB+9uo83PzagqjFUMQUpdgVigxk6/7jUYugiDFKsSsUCkWWoRS7QqFQZBlKsSsUCkWWoRR7ljN3wz7MWFMWtRgKhSJElGLPcq54fhaueWlu1GIoXNh1qBzPTlsHxljUoiiyAKXYFYoY8Kt/L8SfP1mFb3cdjloURRagFLtCEQOOVFQDAKprlMWu8I9S7AqFQpFlKMWuUMQI5WJXyEApdoUiBlDUAiiyCqXYFRnD4i0H0H7MJBw8Vhm1KCclf/hgBaas2BW1GAoOlGJXZAyX/eMrAMC1r6jwzSh4+asNuHH8/KjFUHCgFHuGU1ldg1e+2oCq6pqoRQmNrfuORS1CbFBx7worpCh2IjqFiN4molVEtJKIBshoV+HOq19txIMfrMD4WZuiFiU0KrPwJUYcTvbVuw5jXdmRlG1Kr8ujoqoac9bvjVoMKciy2J8G8AljrBuAXgBWSmpX4cKh8oS/+UhFVcSShEeUsd5rdx/BP75YG8m5hz/5JYb9ZXrKNqXX5fGHD1bgx+NmY3UWLBLzrdiJqAGAIQBeAgDG2AnG2AG/7SrEOJksNy8/ddXOQ3hxxnrf575y3Cw8PvlbqS/SHQePY/n2QwDs7+Prc6xHZFuUW0oaukLff/RExJL4R4bF3gFAGYBXiOgbInqRiEoktOuJKSt24fzHv8jK4boV2Rwmd6SiCq/N3iTFj3zpMzPx0CT/A8mKykS/kjlqeHPeluTfzOa1dc+7yyy3n//ENDz04QppsiiyAxmKPQ9AHwDPMsZ6AzgKYIx5JyIaTUTziWh+WVlw2Qbv/e9SbNx7DHuPuL91527YhzETl5wUE1BlhyuiFkGYB99fjnv/uwxfr0v1ex47US3cVpVJEZdXVnsKm8zJSbxKayQqdr/d78WZG+QIogAAEM+ER8yRodi3AtjKGJujfX4bCUWfAmNsHGOsL2Osb2lpqYTTWpOr3ZSqGneL/YrnZ+GNeVuQDek57Cw9nf3HMm94uVcbEh/3oMjt0F/iI5+ZgV5/+FT4+KCf+ZPAxog92WDo+VbsjLGdALYQUVdt0zAAkY0N9bdtJt6br9ftwXpT1IMrWWBduCHzJ+ov8fVlR+U1GjAnqk4Ot2LUUBY5NmVFxdwG4HUiWgLgTAAPS2pXmLxc3WLn1+xxeUNf/cIcXGCKepBFEF32cHkldh8uD6DlBEHclxpJbQbVY9aVHcEB0+jq42U7AjqbIgjKK6vxzGdrIn0h58lohDG2CEBfGW35JXveufFn+JNfYsfBcmx89JJAz6Nb7PUK83DYZzSKX72u96+gjIHfvLUYLRsU4euxw5LbYmJ3KDh5bvo6PDV1DeoX5eH6QR0ikSHrVp7qhrqIgs+G58bt4Q/CY7PjYHDWuhVXndPWdxtxt9gBYLvpup4E3raM4GhFFY6dcDcs9Dmhcs1iv3H8fJzpYT7HD1mn2PVJRJGHgedZP1xeieXbD3qUKjiSFiT3nplLjgQNJ8v6lWlFuzWVDVEaYXDg2Als3R9cXP/p90/GGQ/wK2i9j0xZsQsHQk5cl32KPWmxy30Ybnh1Hi55Zqb0IfimvUexeItaz2WF+UrL0G8iFvvXa/dg1c5DJhn0yfnwxnkVlfKigrKZIY99gcF//iLQc/DM3c3S0hJkvI89Tnh53txCBQFg3sb9yfZlGlDnPT4NADz7qanW6StHoBiiv6RlXHYRxX71i4kIXqt7E+bV/sMHagESD4fK45FWY8nWxMh+S4CjBzeyzmLX4VHWyX0zWCeu2cUXHpmJo/kg7ousNQth9hm/E8YKPk5oq9Vl3VqZi9hEyUKLnWn/BtR+MM16ZtLSRCjc8SwcrievtfZSkvJykhQVI2sSVhEfFtm4RCuqqvHXKauF24uyj2SdxX6yPm4iqXFenrkB7cdMwlXjZgcnEAcXP/Ulznl4qu33QaRBKK+S8wI8WfvZyYBZH7+9YCueny6eQC7KFe3Zp9i1iylyTUVerHFZzGSm2iWFgtHY/eesjQBqJ3miYtXOw9h1yF551y3MBQAU5Mrrpv/79hIp7UjtBzHtU0ZkpnXINKJ0qXgl+xS7B1tqzgZ+BRfULTZ2nvVlR5LpWJduPYj2YyZh4x7nJfDVAsohNwCHe00Nk/4A6BZ7vqbYjZFOXhXrzDX+EtDply5MXVyUH/1jWiFppGPmzXmb0WHspFhnYxVxc64x5HKPcl4r+h4jmaTFLvDkVcQgF8d0g8K54C/Tce5jibCtiQu3AgA+X7Xb8fiqauffa4yFDqLDnfvYF+j9xylS29y41z6qwKti9f/uCT8XURxymOhZLWXz0KSVYAw4VhGfEYHZOMzN4VeTZUdqR6BRDsSyT7F7OKZBnXz+9gO6WZU+Xy4iEzVOC30mL9+JlTsO2X5vx7YDx3HwuPsijBNVNcLWmf6SNood1cSU7vKS5avnIRMjmnjRb6OA7gydXIHrb3wJK4tdIroSOCwQ01pckBuUONz4tSRF9Fyug/V1078WYMTTM/wJ40CXez/G0CemCR1jdW2iMob2aysIn5oqHiXhlShLAeqwgAa1+m+TsapYGqbL7fS8OBHlL8o6xb5HK7Bx+xvfcB8jlAnScNc/WbYDT0z+ll84B152KZbg1u9FHv2ol6hv3X9caP+CvPQFSvJSA3hryG2JOGMMz3y2BnuPuEf2uEkQB1dhUCMk3X8dxW+89fWFuO8968pURry6oaJ8zrJOsetscvDPmnGziIwTR8b+ffNrC/F3SYWN1+9xXmjk9ly5PXhGBSYxyCRQWjQoAgDUyU9fbuFlktwKr/rKzRiYu2Ef/jplNX47UU4UjhFzmoMwkK3YGWN4f/H25OcoKnxNWroD42dZ15I14nU0oSx2Ab5eu0d6sVm3Pvvnj+VY5WZ4fNKijH1nKf77zba07cafKDIZFCUFeQk5rZS4NIvd63EuAujGwtSVu7Fg036PZ7Hm4qdSXWUNi/nniLzCM6g9WlGFKSt2cbU3bXUZbp9QO6q2MjbKK6tTlL+RDXuOBhZ6bG61fhH/9TW+A6J0oGXGE65RVV2Dq1+cg5+8OMd9ZxeWbq3N1OjWQTa4WNNeEQkPdDUatKYmzN2MO99clP614VQik0FRYr4tzOK74ye81S6tbcfb4ydy2N8+XyPc/sFjldjMOeqUOeSvqKrG7sPlOFReia/X7Ulu57lOd09cghvHz8c6iypgnX/3Ee56a3Hys7mYiJV9+8cPV+D2Cd9gtmm9xYJN+zH0iWn412x3a9sL5p8q8uKMy6OVUYpdv94rPERtlFdWp1j6N46fn9auHZWGUEKZRkLq8+ivS7gPlWu/j9rHzsOXq8uS/mmrn6Zb8SOe/tJT7dLadoI/Ls/FR2v17UVPTseQx/kyFe7zMIKtrmF4a/4W3PfeMjz68ark9t+8uRj9/vQZbn19Ia5+YU5SAfPYIPM37gOQsNzNVNWwZOiuFVZdcqeWl/6IKRBCf3Hc995yd6E8IMvNp6JiOPGjVK94flZKnHUcdBtPfDKvRen24BmbEbVSV2w/hBe+FF9S7ZXjJ6px7ctzcdS02jH1NyT+tYp1P1xeyVUQwdymG8a5GDH3ist9tuiMuwP2Ob82exPufnsJxs/ahOemr0tun7x8JwAkU0nrqWd5fOz6KmIvq1RFHseg0xi7PUu3vLbA9jfGxWjKLMXu4026ZGtqkQyvERaXPDMDN7w6L3bLjPntdXFGPjMDf/popY8WxKiySY9gvP9OiqbnA5/inIc/S36+ctws231F+tSjHwtcA0MHm7pyF/74oX3q3bBUwbJtB5M+azsrvyg/EfpbXpl6D0QKRZzgWKfgx0gLarGUjtnwMYv68bKd+HTFzkBl8Is0xU5EuUT0DRF9KKvNIDG+WUUe7vV7juLzVbvxrMHKCZJ92gO1dKtz9SY3KzwD0pHYYu2Kcca4jmH2+n1Cbdux34cv/yWHcFY3I++Cbk1d2+cxNC7928yUCUsr9AnrpHKmhBvzmpdq57Xmbdxn6W7Ruealua6y8LB6d2J5/rTVqauug0iJ8fVaw3yC6TvPz06WrDy9A0CgZl1w/m1xjEPxz1fxRQKkC+Euj97h3rGIdDHibrGH08se+2QVV+jaPxzCRO0ktXLFhImIoSiSBsBt385N67q28ca8LdznA+z7m9VvvP6VudhrsPB/9NwsnH7/ZKHzmeG5f1v2JdY7fLIs1TqWYbGbqxsZJ2LTLXaPnS3TfexE1BrAJQBelNFeGKRYOD6VxA2vznffCcDMNXvw3qJaBS3T8HANmmHAPe8uRc/7J3v2A/L45v9v2jqMfWep636PCyzssnywpC1Q4t83qNWRZUeci4KP45jf2MOxEMoI74tn39ETjiMeI/WLwinvIOMuOLmL0qKxBPpITFzs0iz2pwDcDSDQ5WMyrTRjJXgvzXoJk/vpS3NwxxuLLIfNdv2Bd1Usj9J5fc5mHK6okh7iN3PNnpTP5hTC2w4cxxyBFMFWD9aRiir837Ra95eVj/1fszai/ZhJ3OcBxKyxoCbGXpu92XcbopLx/hRzzLyRSUt2YMfB2lXELRrUEZSCD/OtDuIFazyH/sjtOHgcFVXVGZl737diJ6JLAexmjC1w2W80Ec0novllZd5SpwblTnDTc1b96Itvvad/1ZWS8byXn9Xact9OpSUAgGb1Cx3b7NaiXspnPe2vjoyXol0T5rTH5gdv6OPT8GMfRT0YgElLtqdtM/NcwMUQ4mKNycD8U+Zt3CccEHDrvxfih8/WTkyLXB+eCVY7ZNwHpyZ0w2fAI5/j6hfmeDaE+ndo7Ok4Gciw2AcBGEVEGwG8AeACInrNvBNjbBxjrC9jrG9paamnE4le34kLtgpZinYEloTJ0KzdIogRPVoAAL57RkvL7/t3bAQAaFRcgNcMfsJtB47jiGGCy/hS9Gp52kWimFszt+/nIdapNKUlZoylxF97ReShNbt2dx+2d6FE8RLYfbgcnyzb4brfF6t2p4Q4Agm/+ROfiq+w3nag1mIXsaR53HV2NCwp8HysHcbnw9gjFmzab2lEPP3ZGox9Jz1dhPEK1A3JNWWFb8XOGBvLGGvNGGsP4EoAnzPGfupbMgnc9Z/FXJai20jg63VyKw15eU28OHMDelpMWOm+Ugbg3v8uM2xHSgSEnzh2Hd58IbKVGmMs5bcBiZetWTl5yZSg/6L3Fm1zXQ5v9kv3+9NnNnuGP2/GAFz70lzc/NpC1xj+n706L22NAAC86zJBHxXmXtewOKHYT29Z31N7e49UYOpK+3udz7E0e33ZUUyYmz5hnZJSIEuiYgLHfJ12HypH+zGTUiYkPbUb8g2oLd8ndmKnavXm30BEXHnVL/jLNCzb5hxKqWObWM2kyWtqGP766bc4XO4xPNAhlUDyHBYbvRSk0K/bHW8swo3j5ztavAGHT/uCMWCbljVTJFtpXPAi8vLth7iyZwJAZXUNRj49A9NXl+G6V+bijjdq025s3X8Mk5fXKvqSglRLe8bq1DkkXsKKRLNCqmJnjE1jjF0qs01T+ymf1+xOLC2+441FgSTUChqeFwqv9euWKMv4t9FVsr7sKJ7krMA+/Mkvuaz9z1btxjOfr8XDH6W6SmQmbbIaPXgaKZiaufm1hQBgWeg7LqsK7dDFi2rxnJ+zivQN475O0VX//Doxmb7v6AmUHa7Aih2HMGbiEmzak2qgmNswXz67RGTWGMs3ChwmmYy22I2P2bC/TJPWbtDIfJPb1eAkMvsN7c/pVxo7dTdhbmq0x7oy57qttfKY/enp+1jNe3jT69a/3qrQd8z1Og4JFJexorqGJesZeIFHOS/detDSCHM60qtB8KYW27/9gCn/v+k+mucGzH3CbiW0G198u1s4SksW0Xn3PZB2fw33I+gOKZOFmw5gQKfGKd3n0HF/D6XrAiWHHcSspXQFx6vwzItC+EmXj8ffX13DXKvfyI5jL6+sRmFejqV1X1FVjcI8sWpdvPcmZXLco4c/6Pw0APDdv89Er9YN0raL5Hs37rnaUDw6+T1jICLufpm2m484diPvLIxuziKjLHbzBZdV5Ddsi10v3GF8aL0W7Eh2XlPvI9Mm47duuTCcEL1W6y1SuLqeg+PB4vE2bD9wHF+tdfaPivwet5fE1v3H0O33n+DfczdbKpWxE8UjQby88KPy7e485LzQSmexRXoMr8pz4eYDrm29t2ib80vLdK+8iFJVXaPNKcVjfiOzFLuJovxU8b1bhOEi4qt121d/uZkVnfkw0VHJoi3pD4xdO04v2Av+Mr32WM5Oz7MXj4U3eflOy9z9e1MqyQssULLY9tLMDcnMlxs13+1HS60nYPXcJyLwpiQ2J+0y8+dP/IeGuiGSKMyMfj/nbtiHK56fZVnwfMnWA9h39ASqqp3vmfnbF2ZswLUv2ddwMPdfLy+ZuycuQc8HPo1NTqbMcsWYbplZ6c3dwLf02aLhJDsOHse8jfsxqpd13Dgvuw+VY86GffiuQzsy+0C6gjL7DQ3fmH2KFoJ87x9fYeOjl6Rtv/6VeXjtF+d4lNHTYdZRMRwm+/TV1ovIznpoqmPbIuhZG7s2r5e06O1+p6wRphXG0E+r8z87LZykdV5hDJi9fi+u1CatjX5x/eeM+vtXaN+42DJVc2pbDOb+b5x/cCt072XEo7td4hKQlFEWu5ti8FrxzXgjBzzyOW6f8E1KnVMv9Hv4M9w24RvLvM023pM0Pli8Hf+atZHrfCL9yU+ptpkWro0g4tZTP6fvU22x0TwpN2ONe5iayMumd9tTbL+79uW5yQlduzajDH8Lm7fmb0GVwMI0xpBU6kDqS9DYvdyUOuD8LFh1VfNCs3RXIP99C2wxoyAZZrH756u1e7hK63mcCBfC7UG/zSXFKuAWFWM4l9PkqetZ5KCnhHXDLM+0b3en7WO+P4u2HPCUVtfqHti59NwmPpP3AsxSgYQ1TI+Darn77SXYuu8YCvNzcct5nVz3t7oPnUpLsK7sKHq0Sp9sNbLtwPGUlduihbfNBsABH6HTYQdi2JFZit1w0V6csR4rdxw2fe/ehpVSdyq9JsKRiiq8NnsTRp/bMbnNypoNImzOStrU9LnBdTjen+NWIk7HfD9W70qfgDU/vMu38y2ySj9Z+iY7l55bn0iuArbZLSxjLi7K5ZnPEwEBp3KkHX7eIs/P93u3whOfrrb0txsZ9OjnaN+4OPnZ78+/fcI3wq5YosR5rUaSUZBZit3w90OT0lO/i76pk+16jLow8/BHK/HvOZtTOhmQqFxz6d9mJj8nq55L7APmh9msQt3CHVftPJQWgfHijPWYZDMRmHIuTs3Ou4jMrEALLSx9r/c6/Vzp2C3J5z0lg/U14VkJLJMTVTX4v2lrMXpIR/edA4THArbqZ/r15kkbzOOiAbwtMuM5hpC473FxxWSUj92Nzfv4bi4PXqye9xclVqgZIxRqGEtLeTB7/V4wxnxZcDU1DM9PX5eS6MsJt1Nd/NQMXPF8agm5hyatxDcW4WRmeB+WO99c5L4TkCaslQtH1gNklcfcd8sRP9v6hO6EuZvx1NQ1ltZwmLQ6JZh0voD1aEB/dGWNjHl0gb7GYe1u+/DeKSt2od+fpnpPtSFARil2t+v7n/n2VdAd27XY5kVv6ErWaHF+vXZvmuJ75ONVuOutxej/iH0SKTe+WrcHj3y8Kql43SzYOIwQbXPNmDCL2tgim58sw+iSZ2ZyZUQE3K9h8r5HvEL1v5qBcVQbeVREHAbcoI515lIZrLFQpPp98JO3XdSw0/PzWHkSdG4cPx+7D1e4lrmUQWYpdhdTyGuSpqkWWf38+CmNE3t2FrVbqTs3Kkxxy2kz+ab9nX6PX6V/SLIFYpanvoViECos7YJ5qG8f1eJMahhd9PkHwnqZH62o4iqHKIJ58l+EQY9+jsrqGl+5fYK8dmHcloxS7G5XhDfqwsyb89PTb/qxCI2H3vnmIm53iQiPmBRbuiJP/eyUcsEqhFEE2UN98wvcajQyb6P3kE0z5vbt7tf+o85pK1ZsD9eHzkvQOW6++/eZOPtPU913DIn9xyqx/+gJS0PPWPHJiRgMcH2RUYrd7WL/8Kw28s7lx2I3HTvdR7UlO8wJtdxib2/990LpMoTFGouoGB56uoTJ6ZjT3P7PfxZb7vf795ZZbtfRX0hB6dFuzeu57xQB612Su3mdD/GTJjl9iVICL6LIVvJhjKQySrG74dVit8KXD9d07DZzdrkASMuIGPgZa7l2QDup7Zk7PldkjsVjvJQzz3y1yxJ1HTFftfw7MLx7M6H9DxxLjDCidgp5jWDSKyWd06GR8LF+lWeKYST5VsqK6HIioxS76+SVhAtWu+DHe1tfrpFvobvx2CepOaXDmHnXSYZvesS8IMjLlQ/jYbHKSmjEXL5PNqPObCW0/wszNgQkiRh+J7rbNioOPTZfuWJCxG3yVGYMqZ+WPlzCF2XhhcU2ybnMPP2Zt2yRXvC7VH79nlRXi/kh5lksEsaDeM2A9gDsozxe12rOztmwD2WHvaeRtoOnZJsVbrlRgqai0lt6DrtCMVzH+ugR3VvUD9RdEsYipsxS7C7X45CEKkr6oxOGBegF3pzZmZLp0oqYXvpk37CbjKxjKKl282sLAjg/v2I3ljv8l6HIeRTc/gbn+gUbahhQLpi7yU8furhHc+8Hc+CWnVIGGbvy1Aq3fNk8kLY2OCYLyNKI4wvHvz/T+fu4/GZdCrv4aKsVsrK48uw2QtEtUbgD7bBaBMaDfr0nLtyKiQu9rVHxdF4mbvGfUpzPnbZYJDmaV3z3RCJqQ0RfENEKIlpORHfIEMwKNz+bXdzqFxZJpNxYt/sItuw7Fqqvmoeo6lk64VeklTsO4cMltXUleQptmPHzSueuK8uco16CnCQXjcnmzcsTV8orq3GfSxSSEwzWRTi4jzf0Oaci8joi17syhGdYhsVeBeAuxthCIqoHYAERTWGMrZDQdgrlLr46u77/wSL+YrR6E9e+PJf7mDDh9c+FOdnk18f+m7cS4YWXntHSsr24Wewyi1rz3icisXj0IHO/h8Hrczb7Ggl+5GOey1t/5r/eYTybvi12xtgOxthC7e/DAFYCEJu+5+TRj+0rkjshchnjXrCYd4J41U7xaj2ekdRP9eIKXix2P2zZz5fqQM8DIrOP8BpvOST3hRJ3/L7MnXK2uOHl1CK3JuPi2ImoPYDeANwTnntg/zHnSAO7C/Yu5/L9F2esDzxkzS/7XFY/RoEsi9rOlcHTvh8JvuV8CY77Us4K2ykrdiUzSPJeO8p4Gzxc/IwiGcSVr8i9CaPgijTFTkR1AUwEcCdjLG1tNRGNJqL5RDS/rMzbxI5roitPrSbYduC4YwKfuPDgB9I9XL6RZYG8+tXGRHum7UG7JEUn3f1Mfq3ccQg3jp+Pe/+b8B/zjsBEXTGZzvhZ/iJ5fM1PMiasfEUSjoVRxEeKYieifCSU+uuMsXes9mGMjWOM9WWM9S0tLfV0niAf8Gx7ZgZ3bhLauWTdFj2ZmNkHOXVlepI2M37u34geLYT29xOupueh0TNd8r4UCWJ+85PpJWDF4q0+Jk7hwWIXuN5WSe1kIyMqhgC8BGAlY+yv/kWyx3XSwYfpmJ+bUSH9rvhN7CWCLFdMp9JEbu2wnWGl9QqF9q+QGK7GOxlORCe9shYhdB+7wL5OtXNlIUObDQJwDYALiGiR9v9ICe2m4Z4P2zsyYuBPVl7RXCh+0ZN2eXmwGhan52znRfTen96yvudz6egFxY8KZP5UPZSfvu0aej5W3BEjNrEdxn30He7IGJuJkPpc3ItJ6MV3gybqFKknqmqkJlzTYRZ/8VLlw3EpGvPtp4CD8cgvV5dxh9XSSRYV4xc/9wgQD0kMI9GfCBnlf3DzsfuJD/1KgusiDKUOmItUh0+Xez8OpF39/nm5jXe/vcTzeXmUgLFveV1JCaT6Ymet38t/HOLpigmy7J0fvF4rvSh1kDZiGC/ojFLsbsWA/dyM2yZ84+NohUy83Ec/5d94DHajUeFvZEgWf4keqXDDq2tVL0o9Yc5mqfKYzxE0GaXY3YjaFaOQQ9j30c2COlpRJW2C2KhvhBa1gMXSFRM3F4SOV8WeQwTGEnWJgyKM25hVij3IpedP/KhXYG1nIvrkn0yiei+7PWjlldUSFbvRYhdZhq4sdhHcRvd25BCFsoAoaLJKsQdJZQgZ2TKJy5/9Wn6jTP8n3AfLTWF6iWu24/Y3al1+IpZbdQ2LpY89rjjV+HWCdF9MgISxhlgpdk7eCTFt6MlO6C41F41ZU8OkWez6wqTEafkf8GrGhC38uNOptCRqEdIgCmHkGMINkQ1gAAAgAElEQVQLOqPysbsRZGd2e6jimE4302BgaD9mUujndXvODpVX4qoXZtfuL+nBFMovIuiLyQR3gt+QxCAgUOhl+IIgqyx2UavqlGL+pb1uffCjZcGVw1NEy5QVuwMJZRXVa3HLIOiXWCp2Cv7aqclTQUTvh9CqP5ebMf7raMuPKbzj1m8OHA8mo6aIayU/Nydt7+KCXNv9M2EAmRPD1d4huNhVuKMoQb5p3ayLuRv3BXdyC9o2Kg71fADQr32jQNsP8v49feWZnk/8/HQ56XrNiFhu42dtSvHJ5+c6vxYywRUTP7WemPeQsVgxarJLsQt2ZhGrJm7DRt50r1IJ+BKMeWdpyue/xDTEVNZlcKsvYCYlBt5FioxwxcRQ+xypqAq8SI1aeSqIaGcWUY4x0+uRLAyJ2SWQRlRui+mrxeoSGJU5EdC2sX1USSZM5ufG7aHygJeJVuWKUcSKLHgOLRGddJdlcfnJ+d2wuADjb+iH5356lnXbPuQKiziupBVl0Rbved+D5KRW7B2b8MfRHi7nn2h14sVr+0ppJwqytTibqHErKxzOTztv3TQApfUKcXGP5jZte27aliZ1vadGtsLv3GnHGMTB24366+TbT2yrqBhBCrVUsi9dx6c869XJx68v7MK1b+uGcrLYxdGvyIuMDukUyREVUcUtb9zLV0QbSEz+Gq9/28bOk+ey02tMvGUA3v/VYKlt+p23ioOZYfcTOjW1f+molaeC6H25nYPv0Qxv35JViCNuk7BhI1KGLqxLpSvB7i38F9Dwyg/Pau34/WVnthJceSpPsZ/Woj7OatcILSWn6PWt2GPwLNllGnG6/MpiF0SPihG5cLz6Os9gas+4e6iIWKbzRd8ZvSJD9BiGLiddMbwv7yAUynd7teQ4L397Mu31t28eILG1WvyOXuPQlV79ekPUIliSVYq9UisyTABe4PBlE/gfUmOVHT/1UTNFsXduWjdtm4whpMjvD9ti533pBOG6MZ76zdH9XfdxQ6aIJYXBZB7x25/i8CjJmnuTTVYp9scnfwsgoawv6t5Matt5ubJcMVKa4eZHLkN8EWQ8SHFcFl9+ohqAv2IdfjH+1NY2i8/0EcVNQzq6thdkCmtZ+LHYr+rXNhaT+V4MtYxxxRDRxUT0LRGtJaIxMtqUgVs1cCL+G2O00nkry1sR9jLqq85p6+m4oCYU4+AXNfPizMRwOuiFKbzYxaATETY8MhJjRnRzbWPt7iOyxZKOn9HrL8/vFAuL3e4pWb7dWz54WfhW7ESUC+AfAEYA6A7gKiLq7rddXzJp/57fpan7vh4mT3OJ8PtLu3vqWGG7YmRaFDKUssh7LTSLvbI6nBM5YHyROk+8Edd9WLkzWsXCg58XaRS1Vm8fdqqUdjIlKqYfgLWMsfWMsRMA3gBwmYR2PcOrfwhAFWcBjUGdGyf/bt6gCD8f3AFntGogLFvYrhin05XYhB7mEFlmM5QhelQTgDLPE4RcxjZlvPszwBPjqyh7Tg7fC04ml/dpJaWdTHHFtAKwxfB5q7YtBSIaTUTziWh+WZnYUmpRRN6IvJMf555amrbNy7PjtTPWDWgCywq7l88qCVagyDLysOLLY6EEfchQlJ/+GG/dH89apDIJ2xNjNfqdtS6eCcNCmzxljI1jjPVljPUtLU1Xkn74fu/U94h+/d184UTE7fOWldfCq8U+955hHs/nxRVjfcyuQ94tLLe2rYiDvrUiiBeBn2yMU359nkRJFHZYRQfpkXhGjKN7KzIlV8w2AG0Mn1tr26RjV0rLvLBEV9an1HEupEFIDWN0wuoF8IPe1kOzJ39sn5XQ60InpyXKTjjp0aMnrH3LVoecJmnxTnuXFZMpSFag3ZrXk9sgB7df0Bmv/OxsoWNK6xUK7d+mUTEK8rIqwC2W8K6adnvxZ0p2x3kATiWiDkRUAOBKAO9LaDeNN28agGHd0idEzddJ/3jtgHaubfpZUdqztXXUzfd7t8bqh0ZYfud18jRMf2KQE7xXn+N+T6xoJ/JCsOHjO87F367q7bsdES7t1RJDu7pP4huVQZGHl7joHWtev0j4HHHDbzdtVl/sBWqnKw4er0z57CZXRljsjLEqAL8CMBnASgBvMcaW+23XiiZ1CzGwc5O07WalpyumvNwcXHiaczy7HyXmdKidBRWEzjTGNf/24tRQOC9ugyAneEVepEb3RFGe/xwzRMS1wtMNEbcJ768136fBFv3cCbt+bHZT6sz+nTfXXpzw+yxdeoZYX7BbmNjrwU9TPscivl5GI4yxjxhjXRhjnRhjf5LRpgjmy8irO4j8vT29HBuENdzMYH3VsZhIEyWIWPs2jergo9vPTblmZ7QWjyrKNHhvt/lV8YtzO0g5/5AuYi+Ik4XF9w3H70ae5mr4ecHtxZ8pUTGhYhUpkaaHuMMdyddFjstim5TKOhJkCuJXNaiTj+4t66dc79NbOvvtYxGtYsGWffGLOLFbaRpFCot/33hOKOdpWs+7O6lBcT5ycwh/vzpc1xyQOT72yLFzxbgf6O8ieznSTVn1auO8WtZSDoffUFIo7sIIwmLX74lR1r7tnGuospS//Wl5Wdk5xeE7r9lgkfXwG58Fu+ADJ7yIMbBTOKOE/xneNZTzAMD8ey8EkAiYiGISXpTwgqMlYdXhzc8sr2L3++h4sYaOnXCOm/eif1ItdvHj09vz30iTugXYc6S2pqdVi5ef1RpntWuI85+YZtmGTIs9T1Ii/DaN6ghZ7CKumBl3D0VFVSJSSfQONKlbaFku0RjJ0cAlSszMD/q0wv2Xni4oSXgUSnA78t6fJnUTE61//XGiKHr7MZPs28wWH3vkmO4O72Ul8qcIvRxbXsmfaEpv3y0k0866a9kgusiHNJlsZGzuIKPRSl+9y1/uk2GnuUem8PD0lcEM3RlLhC12burNGrRbzFXHoNhFXRdN6haiQbHYyyBMZKhPr0YMT8RdlGSFYk+fPK3dsnqXcz6KsN+uXZrVxcZHL8E1/Z07xmM/PAOL7hvO1abx9xp/Te92DT05MPYdPZG2bVAn50UXZkRernbYWexXe0hs5ifVspGSArFBLv+ow+yKETpN8uhLeqYWMskhQkNNOf/qgs5ijUbAzwa1597XOCL0iten/6x2DS23D+vWFA+Min6UkxWK3fzWJcOv2rzPvvyY/8lT8WOKtdVr+tAuXaYEJQV5qKft++uL+Mr3eRaKA56Mgk5i2GU99vJibVKSqL3Ju7hMJuKl/bz5k0Svi/4CuWFwajSN8dloYRodWdUw7dikJG2/oLjQYhR1/3f5lWIlZ56nILB6YRfk5uCl68+2rGUQNhmn2I9WuOd24fax+w539LBcX/v3e72dY2iJEpOYGx+9BLcOdba0jBERKRJJLM+VJ2jxmq+N3T1xtNhtG4/Oh1nHQrFfP7C97f68hbL9zifofcA8SUxUex3N7rGvxwxLW0j396v74NoB7RP7+xPJFb/zHl5zzhtz67hNqtcrysOQLlZ5otLPHZMgOQAZqNh5rl2cfew6bi8fkaaNCcKM3a2GMdsVhoMcIhcKJSxPN/+8fcesh82O1yGG8Y71i9J9zkEMvb26YsyjmByqvYxmHVaQl5O2kC5FYXno32LZO1Pv71PaxCQv1bxvTRNGo8Mt+uixy8/A+Bv6pW236ppxqo6WcYqdp8iF3QUed81ZKZ9J+88rXu6jfozdsfpkYqOS9GGymbaNijH3nmGO0Q4lhXloapF7ZNy1Z9nmtO5Y6n8oaf556y3SAFvtx0WECl80JwuvqObdRK+LPnma3vcpuRqZt8Sd39BSXjzq5SReu8FxU/79l693L6PJc+441fPNOMXOs9TaSml2bFKSVmXdt8XuyRWTHs9t5HcjT8MzV/XGOR3dJyvzcikt0oFXouKCPHx422BMvnNI2ncy5hl547CddqtjmKjs0qz2ZSNafDpKeJWk33eVfry5hGNhXg6uPqctNj56CfJzczDxlgH4x9V9bNupqja69eyvb7/21msQRO6ILrN+b1s1FCue0dhijsALF3QTX31qdbvismARyEDF3rd9I2x4ZKTjPlbX168S5z0P7zFGnfSMITFVnfxcjOLMZ5Kv+yht4tj1B8fcCfUFFg1LCtDVYrHFsm3+867zunOcHgZjnpPGJbWjDjt/chzht9hNOwr+tBobi918jc5q1wiXnJEaOfPoD3qKnQyyVpcmZB49pBM++NVgnG3zsrDjjNan4H+/09U2QiVIrMJLeW6ZXe4e2WScYgfc34xWrhgist3uWQ7PR6bKOKpXS9Qrykvb7oZunRmP4BlFiCaY4qHAZOY/IyGLop3i/vngDrjwtKbJQt2jOYo783BOBzHFwgO3Yk/T62K9q6rG+8vuyn614aO8+XvsJtNFnqfqpMxAT495g24d2hl3fye8Fag61ha79b4bH70k+feVZ7ex3kkyGanY3bBU4Ei/8ET+oti9vBSSPnaX73nQJ8rs5AhzZGguHG52e8mkcd1CvHjd2fjnrE0AgA8Wb7fd969X9MLSB8TXAzjBszRfdNm5WVGIRnzolcDMk6eifYCIQpvCaNc4cR07+ZzT6dIsgiX+FteoRYPw67DakZWK3aovE6U/uDnEN7N+Rd/W1ufx4oox+diTES0s9XsebtQs1RSL3coVY/qJYTy3YTpJjtkUDAESE4b1LCJZrJCUdQBA7f31OhFp7Jc/EViQZc7zU8G50tmquEfQhsFvL+6Gd385EGfY1DWIA3ZBDA0tthtX6epFeMwLxsLyw2enYrfysYPSE+IDePmrDa7t2UUTeLlFVj52wBBr7HBHJtzYP+Wznk86DnM2nBkEYo2scLUF916Y7Bv8rpjUHY0W+0Pf6+F6vK6YzSUceWufTrp9MN755cCUbUHfwjoFuejd1r9/PKi+9tD3etgGMVgtrjKKoY9gG5akGhVhPRdZqtitfOnAhj3pIXc7Dpa7tldjY9VbneePlznHNOtHmI/VH2yn+24X1mi08s8zLKYIK2xN55vfX4R/3tAP4645S2jk8cvzO3Htd9+l3fG/gv7UXYfc768OfySP836N6xb6foBTFp1xNJYMd/T4RDetV4Q+mpL1U0Tc7tipvwmuLqvMtCCL7691253f1b42s52OMRNVZG5WKnYriCgtzwfvg6zHzk++cwjeGF1rNe87ml7c+Rpt1Z6THICDxe7Dbw/AMmbd7Hxx62wXdRcP/yIQGpYU4LwupRh+enPHN9RN53VMsXjuvribbU1XY73QGwZ3sFyF66SIvt3pnCvISBBBNl4nT0VXyzMPrrwgsPNsltqk0BDhw9sGW38h8ScbU0aIPotRX3sjGZe2V6deUV5ywmh9mXvmPwLQv2OjtG086A+ZOTTw0HH39AZAwoVy1QuzU86Z9LVrn2sfTHvs+pndMXY+djf8uJjsPhsZO+I07nYLOYLqC/JyAVjfC5FnU+bKQb0pr6Mm0cnTWsPAWo6wsA1zlSCHbf+X+BuNTYm+6OPkfvRlsRPR40S0ioiWENG7RBTaLMik285N/r3ewsViJifHe9yzXZkyc4Hl+kXW78kBFpkRdV962opDn52DiPDYD88A4CVhVQIZCk60BS8K8I+a71m/viUWv1fEijL+7vaNi/Guyedc26Y7+nm9xrHbuf9sj+dw5YniZ52GmRwCPr/rPPzp+4l75rVv8tCnrXc1ZLTSRS3wlN9uF6nmRSgP+HXFTAHQgzF2BoDVAMb6F4mPtgalyqOvc4jShlZEZLus3khjm5lxPaTP2v2RIG2UYBPuqD/YTp3J9kGj1D9/2Kc17rqoCx4c1UNrWwwZlofoMNaLL9L4Ip36m/Mw/e6h4o0YMA4OiAhneqhmpXP7sFMBAB1NoZEFeTmYPTa9kLT594sut7dz5Xm5lX78wnbHEhE6ltZFxyZ1Pctl92wYty55YDgmjO5vuR8PforWGOXr2SoRlz8ogDUjPPhS7IyxTxlj+hh4NgDruMCA4bHECekvACK+BRl2SqooPxfP/fQsvP4L61V4Gx4ZmRbJYtdW0hXj8FOMHeeL/znfcrueFfK2YaemFUlwi7/u0KTEVQZb2czXVrwJX3RuWtcyFbLIbzH2I4K/0LSLujfDxkcvSQu1JAD166SP7NJ87B61q1liLz8h+ZJwuYsiE9m1QQP6v/7mklK3135RvygfhXnOowF9YZtbW+JrAGr/PrPNKVh8/3CMNIU7hoXMydMbAHwssT1uujZ3LooMAIu3HrR8AfDcPKd9Lu7RPBlqZjUicOvA+hBa5DEuLshNKmGzfE7nc6unql8e8wPdppH4wgvRhyKo4AERMawyN8omx2YFdNoCJb8ZsiLCTupkzVvts5dXpp37RrStB10i12rb9WeeiJYilImrYieiqUS0zOL/ywz73IPE7NXrDu2MJqL5RDS/rKxMjvR625z7pRXkAN/Nc9vDi3FlVsD66jsn/7bo5KmO/vLQh4fm3BrrHx6JVX+8GH+4rAdObVo3xc0FwDZiJVUG87WVo9l5Lq1TVIzTi84Y2gbUrgsAgIoqvrAUPXPilF+nJ1Ozlsf6PpoXCFXFVLHPuHsoPr/LPnTR7l6k/WYPOlNfqWpGZKJ546OXoJizCpaXVbuy2/SK6y9kjF3o9D0RXQ/gUgDDmMMTxhgbB2AcAPTt21dqr+VpbGCnxumK3eYhM+M2maiv9hPJNqcv/R6o+eBe+3k/LNt+SDgtLMA/rO3fsTHm/m4YmppytOfkEIpycjGocxNM+c15eHLK6pTveZ4bczUb0Q58/6juuOfdZbbfB/FAmC0q44jOqRarkb7tEy9JkUIk5v70lx/1Sll/AADVNWLxjvauPO8Xzuqat2lUnL7RgK40C3Jz0KA4H2WHK1La4umrc+8ZBjCg38OfccnpNS+7G+WV9iuarYhRUIzvqJiLAdwNYBRjzL4GXdBwaJ5ESgHTNs5b4dYXG9TJxyM/6InXfm6f8c7s3y7Kz8Vnd52Hv2nJshrXLUx7uNPk8CifcT+zUpfFnA37fB3/k3PaeT5W1jLtYydqQyadpm2MpzO7GNwgi32tKvScqBZTVrVRMamtD+kiPnnnZbWsHomiR2Qt/8N3MGvMBcnva1Np6J/taVqvKKWfNq1XiAe+2912f31lOM8K3ajQJ+LDSingN4797wAKAUzRBJ7NGLvZt1SC8PTDHCLLsmE88LwArurnnM/jPzcPxDpTvL3X5Eeio9r9xyq1/eS8yOKGkyvGarT11k0DUqpO6ejXCai9Vs3qF6K8siYtHUVyP85r9cxVvXH7hG9sZTJTyekK0jm/a1O8v3h72ojPS4hvbYSWy36Gy67/pobFiVGruXh4euACv1xz73F0GiA/Nyclg6IVI3o0T3MxuiG8QMlhd7u0ykHhS7EzxuJf9tyA+Ub95Jx2WLv7MD5cssPlOP/nblRSgEYl/tLC+l2gcaTCWjmltWd6pPUC3NcPbI9Xv94IIPGgOo2AZfdfmUuz+9mk5zWKrCuBWWOGgQHo9LuPHI9x+70DtVj7wvxcrmtTJeiKefxHZ+B/v9M1rSarrwlAgZvotiuZRjZhGw//uLpPWoI0I/07NsLs9amjTpkiXtWvLZZsXYq2Lq4sWWTsylMjq3e5Lxs3K3X9Db9pn/vipvhYsP4E8TIM/O3F3ZKFtx8YdToqq2vw+pzN+P2l3XFR92YY/OcvAKTXq5S1vDrMS2/Ms6LnzXZSBgC/BaZHuVitp7CiUtAVU5iXa+n/DrLvGiWsXZBlLXcy4orDFRMEbtfBHJYMyKvQBCQUu9uoXiZZkSvm67V7XfcpyLW+s5f3ae0aBhinHBCAlYLmk8+L1XvL+Z0s80zn5RBaN6xVJN8zVYYx5wV3WsTlFdlRCEZrinc9AW/7dkv+bc8Rry7njiavXRezCgWOE1ahyW7x8HEmKxQ7TwUds89Ppyg/F49rEz52xKUP+nXF8Gbtc1Rqpu8aFudjUOf0lAk5OYSbhnTETdq9CTLnttOvEgmF61zKV7ChR6vaRW28qQN4cgEZaSmpaIOXruspfJf7WLHJZlmE8SKJiZoAkOGumPO6lGLj3qNcoWbGDtfSFMrm1hnDmvDgxWtIsGhUmHk5vBXf3GdfnWjsyNPw6fKd2qdo4rJFQuHIlFLAjj99vwcmLtyq7ZfYxnsWfgtfzvUKq/Rj7XVwlltEnLEjuqE/R1F3RToZrdj/eUM/AMDLM92LZegd7p1fDkybwDDHYJuJi1q3D3d0lrAwLwcVVTUCUUAJRvRozi2bbVvEZ9Ha4XRYPW1S1yklrMhCH2ORCnNCt15tTsHiLQcAJEZ5p7Woj5U7DiWvqdvIgCcXUMr+kt6Dnix27V+zO82pbbLxxdxyfids3nss7RievnjTeXx5+hXpZLRi1/lJ/7bYfbgCz01fZ7uP/qD0sajY4vYQxcxgx+GK1BS1buJNvnMInpu+Dr19JLXyih6/f3639IozIljdg/O7luKxy8/AqDNbpn+pIZJb3jgy69w01S3zzi0DHZW3m5tLzxnzgz7hVKnX8dJ3h3Ytxd8+X4Nhpzlfu9RUFol/ze9RfWWuxdHigim4yQrFXpiXizEjujkqdifcEi7FZaLHa9Hq9k1K8OjlzvMIQdGxtC4W3XdRIHkziAhX2FR9b1xSgL1HT6BHS/ckbzzk5hByDcrIrMjdjIO6hXlY+YeL7fOVm5BmsXvou73bNsT6h0e6HrvTUJ2Kf1I4daGSIhiyQrHzIGuCLUps83BItn6S/lJJl+WUYnlhY6KIKBAvyka/9jweH3OMuY6Vx0Nv7gKfIx2v8LwQ3py3pXZ/fRLZzcdu+lcRDFkRFeOXhhEqHhHsJgJlWz9xGaHIIO7v7Ae+2x2NLeYI9KRgvQKMJvLLYEOucVFjIIu6WCzJOsXetZl1yJpTh+vQpAQf/MqmnmKMEF204peY60RH/CgOnuIrsrh+UAfL7ed1KcUrPzsbtw6N7wTi+V3TRxNufUYp9HDIKlfMh7cNRuuGdg+lc5fryVFwI2qCymLnhfhIIg8i4O2bB9imhzXCUxjFL0MtFGecsKoD4DaJXFvrV2n4IMkqi71HqwZp/lyRNLijetlHV8SBSi1/SDdTUe0wraCoHsgwXCoEQt/2jdJyozsecxLrJ7L4W3Y8v8IbWaXYrfi7lhaXRzF4yYUeJrrFbs5MGJSyjYN/OsznX0TZyFpAlC3U+tj5FigFfV/fu3UQ7hl5WsBnSSVOc1Px1mQSGNKlFCN7NscDo9zLYXVp5i2NblgUabkrWpncTfInTxP/WimvsBWaUp/xxTiyuVpLcNW9BZ9LM2gl2KvNKbiRI9VItpJVPnYrivJz8X8/OYtr318M7ogz2zTEFc/PClgqb/Rs3QBPX3kmLjQtHJGu2CWUC5SN6G+s8eAD9/KbTmZf8XdOr12ZPPz05q450RXhkfWKXYScHLLN1R0XLjszfeViFMol7pb093u3wkszN1gW1LAjTkPpTICI8P6vBgmFC6tLHA5KsWcBgT0sFto77hbqwt9fhKqaGjQpKcRdw7ugiKMQtxfiMP8QB0Szdur9JyckJ3B+LoUeJhwHst7HfjIgW9VmslXVqKQATesVISeHuKrR//nynsm/PbliMvhaRUlYBsLHdwxJ1mE9mVAWu8KWk8HOMSb+EouKUXgh7Bdh56Z10blpvIMigkCKxU5EdxERIyLxkugK38ifPE3AW5gjSPTkYe2buC8a8kKeTWUtN/RrUxu/Hf21ygSS4Y5qpBMovi12ImoDYDiAzf7FUXghzEIgDUsSE2Uik5K8NCopwL6jJ1K29WjVAK/87GwMCKjgQqrFLn4dvV7631zUhauQSbai9HqwyHg6nwRwN4D3JLSl8EBQE4RW3Dq0E0rrFeJ7FtE5fpl0+2Cs3nUE1708N2V7kEvrc12KSQTF7cNOjeS8UaP71ssrnYvbKPzhS7ET0WUAtjHGFqtQsegIbIGShXehMC8X1/RvJ/eEGi0a1EGLBnXw7E/6BOZ6MZPr8eIpx4s3dh9O5HA35nLPVK4d0A7jZ22KWgxLXBU7EU0FYFUj7R4Av0PCDeMKEY0GMBoA2rZtKyBi+Ey4sT827zuK305cGrUokRB1SOOIni1CO5d/gyS+Bk0Th5KBUXGkvMp9J4VvXBU7Y+xCq+1E1BNABwC6td4awEIi6scY22nenzE2DsA4AOjbt2+sDZ4BnRpjQKfGGaPYA8sVE0ir8SKbB5rz77V8dCPF7/W+95LTMqYwTpR4dsUwxpYCSDo/iWgjgL6MsT0S5FIIkM3KKa4kVlsedS34rDDj73r94tyTN/+LCCqOPQtQC5TC59mf9sEny3Ym5wL0YtUKZ1TfCgdpip0x1l5WW4p4oEa89jStV4RrB7RPfg6z6pIbM387FGWHK6IWw5K2jYoBAMMiquUqkzi/o5TF7sDpLetj+fZDUYvhTkA97GRYdBPnh9MrrRsWo3XD4qjFsOS0FvUx7X/OR7vG8ZQvW1CK3YGJtwxEZXX8422jjmJRKEQIK5Q1bIoLwltP4oZS7A4U5eeGuvjHK8pvqVCEjzlU9qHv9YhIknRUdkdFGrWFiSMWRKHIEHq0qp9WbzlKlGLPAqRHxUhuT6HIdro2qx+1CCkoV0wWoNI5eEfmtXtw1OkpdUAVmYlerEWEOGRCNaIUuyIN9Z7wxnUD20ctgkICjUrEXSrxUuvKFZMVKD3snxYNiqIWQZFh3DHsVOgLj+OW5kAp9ixAdvGCltpiG30xycnA6S3j5SNVxJ+GJQV48sdnAgBq4qXXlSsmG5Adxz68ezO8/otzAituEUdiZnApMoTaCLJ4dSCl2LMIWeqdiDCo88lR5VC5sRR+qC0jGakYaShXTDaQdMUoNaVQhEmyKE3Mpk+VYs8ClD5XKKJBt9TjZlQpxZ4F6F1KpQZXKMKlcd1EaGSn0roRS5KK8rFnEXGzGjKJeA2kFZnCwHQy/FAAAAfUSURBVE5NMP6GfhjYKV6BBkqxZwG6QlcWu0IRPkO6lEYtQhrKFZNF5CiLXRjdUlcvRUU2oRR7FqCvestVil0Y/dopN5Yim1CKPQtger4ipZuE0ReWqEunyCZ8K3Yiuo2IVhHRciJ6TIZQCjFIu4vN66t8J6LoS8GVG0uRTfiaPCWioQAuA9CLMVZBRJlfoTYDqV+Uj7/8qBcGn3pyrBaVie6KyVFjV0UW4Tcq5hYAjzLGKgCAMbbbv0gKL1x+VuuoRchI4rrARKHwg187pQuAc4loDhFNJ6Kz7XYkotFENJ+I5peVlfk8rUIhhxrlY1dkIa4WOxFNBdDc4qt7tOMbAegP4GwAbxFRR2aR6owxNg7AOADo27evWg+iiAVM+dgVWYirYmeMXWj3HRHdAuAdTZHPJaIaAE0AKJNckREkfexKryuyCL+umP8CGAoARNQFQAGAPX6FUijCQlnsimzE7+TpywBeJqJlAE4AuM7KDaNQxJVkSTOl1xVZhC/Fzhg7AeCnkmRRKEJHWeyKbERF7yoUUAa7IrtQil1xUhO3yjcKhQyUYlcooKpQKbILpdgVJzVqql+RjSjFrlAAIOVlV2QRSrErFApFlqEUu+KkJldbclqQpx4FRfagap4qTmq+17sV1pYdwa1DO0ctikIhDaXYFSc1+bk5GDvitKjFUCikosafCoVCkWUoxa5QKBRZhlLsCoVCkWUoxa5QKBRZhlLsCoVCkWUoxa5QKBRZhlLsCoVCkWUoxa5QKBRZBkVRyY6IygBs8nh4E8SzrqqSSwwllxhKLnHiKpsfudoxxkrddopEsfuBiOYzxvpGLYcZJZcYSi4xlFzixFW2MORSrhiFQqHIMpRiVygUiiwjExX7uKgFsEHJJYaSSwwllzhxlS1wuTLOx65QKBQKZzLRYlcoFAqFAxml2InoYiL6lojWEtGYEM73MhHtJqJlhm2NiGgKEa3R/m2obSciekaTbQkR9TEcc522/xoius6nTG2I6AsiWkFEy4nojjjIpbVXRERziWixJtuD2vYORDRHk+FNIirQthdqn9dq37c3tDVW2/4tEX1Hgmy5RPQNEX0YF5m0NjcS0VIiWkRE87VtcbiXpxDR20S0iohWEtGAqOUioq7addL/P0REd0Ytl9ber7U+v4yIJmjPQnR9jDGWEf8DyAWwDkBHAAUAFgPoHvA5hwDoA2CZYdtjAMZof48B8Gft75EAPgZAAPoDmKNtbwRgvfZvQ+3vhj5kagGgj/Z3PQCrAXSPWi6tTQJQV/s7H8Ac7ZxvAbhS2/4cgFu0v38J4Dnt7ysBvKn93V27v4UAOmj3PdenbL8B8G8AH2qfI5dJa3cjgCambXG4l/8E8Avt7wIAp8RBLoN8uQB2AmgXtVwAWgHYAKCOoW9dH2Uf832Bw/ofwAAAkw2fxwIYG8J52yNVsX8LoIX2dwsA32p/Pw/gKvN+AK4C8Lxhe8p+EuR7D8BFMZSrGMBCAOcgsRgjz3wfAUwGMED7O0/bj8z31rifR1laA/gMwAUAPtTOEalMhnY2Il2xR3ovATRAQlFRnOQyyTIcwFdxkAsJxb4FiRdFntbHvhNlH8skV4x+8XS2atvCphljbIf2904AzbS/7eQLTG5tCNcbCcs4FnJpLo9FAHYDmIKE1XGAMVZlcZ6kDNr3BwE0DkC2pwDcDaBG+9w4BjLpMACfEtECIhqtbYv6XnYAUAbgFc199SIRlcRALiNXApig/R2pXIyxbQCeALAZwA4k+swCRNjHMkmxxw6WeK1GElZERHUBTARwJ2PsUFzkYoxVM8bORMJK7gegWxRy6BDRpQB2M8YWRCmHA4MZY30AjABwKxENMX4Z0b3MQ8IF+SxjrDeAo0i4OKKWCwCg+apHAfiP+bso5NJ8+pch8UJsCaAEwMVhymAmkxT7NgBtDJ9ba9vCZhcRtQAA7d/d2nY7+aTLTUT5SCj11xlj78RFLiOMsQMAvkBiCHoKEemF043nScqgfd8AwF7Jsg0CMIqINgJ4Awl3zNMRy5REs/bAGNsN4F0kXoZR38utALYyxuZon99GQtFHLZfOCAALGWO7tM9Ry3UhgA2MsTLGWCWAd5Dod5H1sUxS7PMAnKrNNBcgMRR7PwI53gegz6Jfh4SPW99+rTYT3x/AQW14OBnAcCJqqL3Zh2vbPEFEBOAlACsZY3+Ni1yabKVEdIr2dx0kfP8rkVDwP7SRTZf5hwA+1yyu9wFcqUUPdABwKoC5XmRijI1ljLVmjLVHos98zhj7SZQy6RBRCRHV0/9G4h4sQ8T3kjG2E8AWIuqqbRoGYEXUchm4CrVuGP38Ucq1GUB/IirWnk/9ekXXx2RMZIT1PxKz3KuR8NveE8L5JiDhM6tEwor5ORK+sM8ArAEwFUAjbV8C8A9NtqUA+hrauQHAWu3/n/mUaTASQ80lABZp/4+MWi6tvTMAfKPJtgzAfdr2jloHXYvE8LlQ216kfV6rfd/R0NY9mszfAhgh6X6ej9qomMhl0mRYrP2/XO/TMbmXZwKYr93L/yIRPRIHuUqQsG4bGLbFQa4HAazS+v2/kIhsiayPqZWnCoVCkWVkkitGoVAoFBwoxa5QKBRZhlLsCoVCkWUoxa5QKBRZhlLsCoVCkWUoxa5QKBRZhlLsCoVCkWUoxa5QKBRZxv8DgvCKTN+NMlUAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"hmc_states = lax.scan(lambda state, i: sample_kernel(state), hmc_state, np.arange(num_samples))\n",
"plt.plot(hmc_states.z);"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"HMCState(z=array(0.48313957, dtype=float32), z_grad=array(-0.1292151, dtype=float32), potential_energy=array(0.03339309, dtype=float32), num_steps=array(5, dtype=int32), accept_prob=array(0.99635506, dtype=float32), step_size=array(0.8920079, dtype=float32), inverse_mass_matrix=array([4.1756053], dtype=float32), rng=array([3589372389, 3207553271], dtype=uint32))"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"hmc_state = init_kernel(init_samples,\n",
" num_warmup_steps=warmup_steps,\n",
" rng=random.PRNGKey(4))\n",
"hmc_state"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJztnXecHMWxx391UemUTzmcBEjiCELSIQRCAoQIIvphsLENNibIYEx62DwRTDJgPWxj/GzAFsHGBJGzAIECGBAITlkoJ5R1p3wKJ53u+v2xM7uzsxO6Z3rC7vaXjz7szc501870VHdXV1cRYwwKhUKhyB0KohZAoVAoFHJRil2hUChyDKXYFQqFIsdQil2hUChyDKXYFQqFIsdQil2hUChyDKXYFQqFIsdQil2hUChyDKXYFQqFIscoiqLSjh07soqKiiiqVigUiqxl1qxZWxlj5W7nSVHsRHQLgKsBMAALAPycMVZvd35FRQWqq6tlVK1QKBR5AxF9x3Oeb1MMEXUHcCOAKsbY0QAKAVzqt1yFQqFQeEOWjb0IQHMiKgLQAsBGSeUqFAqFQhDfip0xtgHAHwGsBbAJwC7G2Ed+y1UoFAqFN2SYYtoBuBBAHwDdALQkossszhtLRNVEVF1bW+u3WoVCoVDYIMMUMxrAasZYLWOsAcAbAE4yn8QYm8AYq2KMVZWXuy7qKhQKhcIjMhT7WgDDiKgFERGA0wEsllCuQqFQKDwgw8Y+E8BrAGYj4epYAGCC33IVCoVC4Q0pXjGMsXsYYwMYY0czxi5njB2QUW6+0NjE8Mo369DYpNIUKhQK/6iQAjHgxZnf4bbX5+PZGWuiFkWhUOQASrHHgJ37GgAA2/cejFgShUKRCyjFHgMKCggA0MiUKUahUPhHKfYYUEAJxd6kbOwKhUICSrHHgELtKajFU4VCIQOl2GNAcsSeJ3qdMYbXZ61HfUNj1KIoFDmJUuwxgJKKPT80+2fLt+LWV+fh9++rfWwKRRAoxa4InT0HDgEAtuxW2x0UiiBQij0GUNQChIzmBAQG/zOULbvrlUlHoTChFLsiAuStKZzw0FT84rlZ/gtSKHIIpdgVoUOSpyifLlNhoBUKI0qxKyIjT9aKFYrQUYpdoVAocgyl2BUKhSLHUIo9Bsi2OcedPPu5OQFjLOmmqog/SrHHCKaMzoqY8swXa3D0PZOxcef+qEVRcKAUewxQI1hvqI4wPCYv3AwAWLd9X8SSKHhQij0G5Kt62nfQ39Re6fUQSW4qU2QDSrErImPGym2+rldKJnxUZ5odSFHsRNSWiF4joiVEtJiITpRRbr6gTDHeyJegaXFAb6MywkAogqdIUjl/AfAhY+xiIioB0EJSuYochCS5AT02fYWUchTu5JvnVrbje8RORG0AjATwNAAwxg4yxnb6LTcfmb9hF3btb4hUhprd9Zi3Ljse34fagp4ieEgZ2bMKGaaYPgBqAfyTiOYQ0VNE1FJCuXmDPoKds3YnLn96ZqSynP6nT3HhY19EKgMvyhITHvqIPV+SwWQ7MhR7EYDBAJ5gjA0CsBfAOPNJRDSWiKqJqLq2VgVtsmP++l2R1l+XRZtQlL03fJRJJjuQodjXA1jPGNOHmq8hoejTYIxNYIxVMcaqysvLJVQrxpcrt6Fi3CS1wSKHUDliFQprfCt2xthmAOuIqL926HQAi/yWK5sXv14LAPhmzfaIJVHIQqn18AjC7DV77Q7MWbtDfsEKaV4xNwB4QfOIWQXg55LKlYaaQeYOew8cwlH3TI5aDIVPLnp8BgBgzfhzI5Yk95Ci2BljcwFUySgrKFKLP/Eb5+Wb3dJvKIAtu+slSaJQ5CZ5s/M0z3RnBje9NAd3vbUgajEA+DehKNO6QuFM3ih2naAH7B8u3Iw356wPthIPvD13I57/am3UYgAACn1OUeIW/Gvyt5sx9MEpeGvOBk/Xr966F3uzyBspl5i7bieqc3DdLW8Uu6zdjm5c+/ws3PLyvFDqioptew5g1z7vG6n8Poq4jdh/8dws1NQdwM0vz/V0/Wl//CTy/QtubNt7IGoRAuF7j32Bi//+ZdRiSCdvFLtOzAZ7WcmQB6Zg4P0feb7e7hls3Lkf89e773rdse+g57rjyuy18d7tu2zLHgC5ee9zkbxR7KKDRL8hZUXIN/u/Xd960vhpuOBvzrteX61eh0snfCVfqICZuWob/m/qchw41Bi1KL4IYu/AwUNN0svMd/JGsYuEunhzznpU3j0Zy7fUBSpSLrFh536c/9fPsW2P+5Tdj2fSZ8u3er42Sn444Ss88vEy3PP2t56un7lqG5ZF2B67t20OAOjdQX60ENmeag2NTXhn3sbYrcWESd4odj2IEU8jmrq4BgCweHO8FHtNXT1210cbJMyOpz9bjQUbduFNjgVE/RGUFok3v2x3DfUaMuKHE77CmX/+j2Rp+OnQqgRAdswu/zp1OW6cOAeTv90StSiRkT+KXaBFFmgnh9bjcwo39MGpGPnw9ICFCQPv9zVsxfLdtr24880F0kwQcRxD3vTSHFz3/CzHc/RXIYiOVbZ5Z7O2z2HX/vxdD8gbxa5z22vz0dDobNMLezPTb99ayH3uTh/eKEEiEpCryYeSKAh5yH7DxDl4YeZaLNgQbXC2IHl77kZ84BICWX++FEDXKvs9y2MLTJKsVexvz92A+9/lD0ljbI4ra/c4npsasSf+3n+wEVf96xvU1Mnd8Tj+gyW45O8zpJYZNTxupX5evLDcVjPqDajcbLEDBylmk+S1U11Up05o+Php+Nu05XIrdmDrngMY8NsPMDekXAdZq9hvemkunvlitadr3UYdyTRgWgs58u4PMXVJDYY+ONVTfXb8/dOV+GaNexCkpiaG+971tugWR/yE2y2QoGEZY9ixNx7TdC8Kc8qiLXhn3kYp9fMsdhsJol8NbGbsIOuGnfvxx4+WBVOvBV+u3Ib6hiY8+dmqUOrLWsUuinEK79o4Y5YsZtGm3fjnF2ukl/v16u2oGDcJa7bu9V2WyLvpx6QqQ7E8O2MNBv3uY6z2+Lu37TmA9xds8nSteYTu5VZc/e9q3Dhxjqf6M+ThPS/Al2FviK7F+ULWK3beqayXxdO4BAwLyq78xuxE6IMvV23L+M6riYBHUj/mBxk23ulLE4le1mzzptiv+Xc1fvnCbGyXMOoP2hSzZuteDLr/I6zfsc/ye9676UXK3fUNrmZPQL4pRgYLN+xCxbhJ0lxMw9YkOaDYxa/hHLDHZsheENBTcrp3MenTMgjqXthhdR/W7UgkaznksgjPVb7vEpx5pXodduxrwNtzrU03Vh4pW/ccwFaTicZLB3TLS3Nx+p8+dT1vZwy9VyZpM7KPF8l1mQxrhSj7FTvneSKDXkqaYuKh3cL2BAGCVThWMyH+QEwRLZ5aVOvlHi0x7Y0IugNNtmWbiqxs9VUPTEHVA1MszxeRd+qSGq7z3HYbiyLjnuprOU1xC0zESdYrdn5zCb+NXZ/ux2XUGpQqc7oPQZoIrIrmDcRU38C3JX/zrnopI2orzIvrXmGMYeLX8iNuPvP5ajz+yQoA7m153XZrE40Zkd968FAT5oXk/eGEn/cmZY7lv2bq4i2oGDcJtXXRB0zLesXO2+C8jNjjQlDyON27IAcqfhQiz87WrXsOYNjvp+L3HyzxXpEDsmZ0X67ahnveke/tdP97i/Dwh0sBGPdkWJ/77JffCZXN85vHf7AEFz4mdxQugpOMu/Y3cOU9TnbeAs/4XzPWAEg4O2TIFPIoUVZqvMjwtsDprCnjMlLXicJ326vS4hE16EXpnVoEwk+W1uC351VmfL+iRlvQ8yhGchemz7nUvgPBBwTTozGan+f6HfuwaRf/vgyR9rBoUzCbud6euwElhQUYc0xXAMCeA4fQvLgQhTY+sFbvzYj/nYbd9RxeOCQ+a0+1C6diw3mXs37E/i6nP6+I//MSbSXcapRjXmzaE0KChGhMMfzlzFi5NTla4SHqfnMDx4iNB7/vKO/1dT7iA+nJVWZ9l75fYsTD03GJQBxyvT1EOei56aW5uO6F2QAStu+j75mMO94QywrGpdSR0hciPze5OzcGM35pip2IColoDhG9J6tMO4whATbu5Bt1kICNXbcPbrdILnDRE6mdom/P3YCj75mMxRZTL5m49fKNTQxj/vIZpgiu4Mt6SX/8ZCpJRB3PixOScpBRjdNIlbFwwjsfc6/32Pc6xtC4izftFn72Iqfv2Bts2IuauvqkPK/MWhdIHbq+WMXhrqnzxYptaddGicwR+00AFksszxZjKAHeKaJILzq4V1sAQL/OZRnfGReFPtX8ob/dGLBid/l+9/4GLN60G7e+6i1zk16+ccOO14QKj3zsvpvPbIrhXcDjx/qOrdu+D1c/+43baTYlpk5myf8zjA/Iju+FF2euxe02I1jjLR/zl8881xH1bAtA2j237KAkCKnri/fmi29E06996rNV+NNHibWOhsZw75wUxU5EPQCcC+ApGeW58ZXFhhoR3N7ndi20EKUuvUFRYeL7xoB3WLg1Cb9TP7184/byE38/zZcJgKc+nREhRawc/8ESTFlscMHz7dUC/Ftw8VGEGSu2omLcJO7z73hzga2Xjd91jTjFtGGMTx4/r4WMax+YtBh/nZbwTrpXWyTftT+cIH6yRuyPArgNQCh7yIyPlNsrxkP5btcUartlgu6Nw3qpzLV4WT8Y0CVzlmMmsh29kmbIuvhBJ6CW2Wn4vePJWUpEz+7LlanBHGPB7zAp8BOUyOJS/V36z7Ja7+UK4FuxE9F5AGoYY44BnYloLBFVE1F1ba28H2d8wPPX77T1czaOvnlXpt1OK9ZG7EH5S+vwx/NgaGpi+PPHy6RsdzeaH9bv2IeKcZPQ764PHK85o7Izh5y+RZOCV/Wg78r8n9fne7q+Ytwkrq32H37rHEpXBN8xzzkvP9TYhKUBZHr60ZOpdIgMzm1oeQ2/XTwIcsXGPhzABUS0BsBLAEYR0fPmkxhjExhjVYyxqvLycl8VWt22DTv344K/fYG7OGKbu4069O/dNzJp57vW6A83RZjchALg8xVb8ZepyzH4dx+7llt3IDEtfGz6CutyDb9/lLY13C0/5SEOBWI84/MQU91NMtlLD/mcaa3Z5n1tYLJEpc2DrJG2WymPTV8ppR53OVKSmMN367Hzva45AfzedlbkhFcMY+x2xlgPxlgFgEsBTGOMXeZbMgGq12zHj7UefaFNQgTjzXZrnPr3V/6r2nYX2bvzNuJAaEl4XTV7kkMc9v7h46fh0SmpUf36Hdbuf8b2yZtw+ECD+3nJjhOEy56e6XJ2cBz0OdPiNSlZbUtvbGSeFIBX99qQBuz4bntmYDXZ5htzcV7DdzthDP3wzy9WC0UCNbuWRkFW+rGnPVfGMPa5WfjOZfSU5tXg0s6ML8EDkzKTedwwcQ5umDgHL32zjqs8v3CXz3nehp378eiU5bYbO5K4fG01va8/5L7pxs/9OqWf+2yPOKdSbpm03OCNI2IOqAV4n+XV7PaW7MXvugbv9VbrTbLfDzdTjGzue3cRLn6CPyGOn9G+LKQqdsbYJ4yx82SW6QbvSCR9dOR8kTEglZWt2uuDa2hswi9fmCUcCtQs7XLT9akt7mLofrdAIordFyvSzSJutkIrxb7/II9i9/5W8qxp8ZrI/NpCuUfBVgHEmDflpHtZiMJbl1tn5VZOg8XMLiZLKr7g2p8RI7JyxG6EKP29sVsYNR51a5z7DQuwn0m0AS/csAvvL9iM3wja/szynmHKVm/320QWzK75dzUenZKeKsyLqYAnlouoWWD/wcbkorjMLdn6iP3Cv32Om18ST1zhZxQ8fWmNJ/OI8f7OXruDOyomr6yNhvMqxk3CgvUJ0ybvT7WaBcnPacpCX4AX+Q0FRPjOY6x/WWSlYhdR0lbIbhPcAzdKLXKKlc93xZ4Dh7DXMGLu7+LB4oabCvXqVSJ61ZF3f4gTHpKblhBIdXzz1u/CWzbxynmud8NqZjB33U7fCu+ix2e4RsXs27ElAGBQr3ZcZZpFemtuoiPRn7WbzFbrFmEo4aA90+x+98xV2zJmOUTA0s3yPYNEyErFnubH7mHnqVtDCzoFo/B2boHzjSnTeDxUnHDbTCEi176Dh3BAs797McUEsbHD94Ki7+uD13gjtTWJwzu14jrf/D7pnde67YkF9sdtPKiS11v8JAaGXfsapClfZiFn4HtJDJ+/27YXv3tvEaYu3oIfTvgqY/GWKLV5MSqyUrEb4bWTGqfwUSXQ8GpFiMrv+x+fOifeFZGr8u7JGHjfR3hw0iLXUZ+TwvNyC2t211uGLeAZMSeDX4Fh2Za6tJ2gvCNuu+cuMzSy3T1zS7SRWY5zueZEIdOWbMGE/zi7ODY1AQPv/wj/83oq3MGuffYd9Q0T5+D0P33iIGT474Sxvuuen42nP1+NT7SQIitr080uBEpuXoyKrFfsNXX12MexYOfXfCNKUxPD+ws2WS5GiXYsccnkZEbUlJDI0r4aSzY5T1N5nicvq7fuxdCHplqGLbBSdgcONVp6njCWmSbNaUa0a38DRj/yKZZs3m3bGfneNGTAKm8tIH+BuFlxYdrfV/6rGg+9vyQ5G7NCd8F9c8765DF9D4XOsi11SXv+u/M2YmXtXscRvvnOhek7ntnu0/8uIKAwYmf2rFfsr1SvT1vstBuZGF9C6e5XWoF6cuiEXOvwyxdm4/mZiW3h1784O5kCzK1+80viajrilE8Ut0GHVxtxg4tCO+qeybamF6f3pb6hEW/OWZ981owxnPbHT2zPt5Li+hfmYKiFPZ+51G3m8+VbsaJmD85+9DPstYm7LnNRcff+TK+Nz5dvNYzYna9v3zIRH+k5UxgDXcbubZsDAC4a1N3y+n99sca2bP1xG0Uwp3s888//wfl/+zztmN0+AwaW0aajmNXa3VsiCj03r5msV+y8PP15yg4W1Aj4v19Jebv89u3EDtiV2vZm865HJ/ZIdq16bdZ695MscPNAsbuLbqFseYKm7fQQXfLV6nW45eV5eI/THdVqNjVlcWJUnqk4mNDo13jrgthAw8NlT89MRiN1a/NNWjiKB99fbDqe+H+v9i0AAF3aNLO8vt5hY5p+L2Up36hDUpjfC3MHTQSUFilTTOiE0TD0xZztFrZE1xABpobjV97NAplyjLj5jDObd9ko7659DRmmjfcXJLbTy+5g9R20NZw5J50mDuZ7/tj0lZixkt/1lacLsBuxvzlnvbQELlvqEveex2HASp6MDs6DDGaT06Zd+/GHyUtdr7OT+YOFmyM1TupRT/WZ2FrT+g0BKC1KN1m1KEn/O2hyTrF7md42NjGp9k4jVmYQt5rMSsFNAb43z3k24PWnuY1QrXI7JupLVTjw/o8yTBt6v1VSaN/8rOpevGl32qYqu3p5TSZlzewzQ5pv2ZTFW4T2NHClCLToGOet24lbXp6HO98Uywxkh55+b8J/3BbCWZoPe1JG7Zi+Td7Yno0LoEu32OckMJd700tzufY7ODVbs6g8IS8WOeRNEDFXlpeVAgB2awr+q1Xpewms3rcrh/fhLl8GOajYxa854aEpGPKAe9AsO5zahJfRNlFiVKN7criVcYeLEvAaL77EZTppjLhnxO0ZnHx4R0/yjPnLZ2nrKTrNTQt6vPe8rRZ33wo9kJTo80spGHfNbjUI0UMBbxEMHWDXkdRp5W1zifbJYN3R6M9St3f//oMlyU1IT36W6iz0WdjymsyFcfPPlBJjyVTmDye4p/m75eW5tt+9Ptu9ozFTZJjS3vduKin53HU7M+6D0f1x/vqdCJqcU+xWmEP5Ghva9KU12LrnIHY6uF/xYNfjW422eUYHJ/5+WtKTw+9cYqvHEL5Vvfk2tWTAKbAsk1grbeQtvj/A/oLvPfYF3psvvmlJz6eavm/Cup4nPrV3E5QV+rVHO23Rc7D1omcSG1OMeR1i576G5H2xilC5ZXemGUyfDYs6ioiMos1umKLw7OB9/JN0H37jAvA/TYvHt7xsv7t8h09dw0POKXa9MWzfezAZt+Tk/013dTMq238LJGG2rRP2W5ztpohPOkyNzSNePxtZVtbuwYszrbPquOE12UDtnnr8Zcpy+87OwkvCjIgS2FN/CHe9tQD7tA78lepEcDanBT23+gFghYe43mT6PwBssMnLu6rWIhKicI188uheLXYwZJpM7LDbDGRnzkyayEwyueEkzVOfO5uWkmVIXFB7+MOlWpmJv716NPrJ4cGLvZExS9Efox6PfMa4URnR9YzPWlbsETvbfloqNo2DjU0Z3gdGjI3x3Xkb8fpsb14tAPCYKWgUT5CulBze6rzl5XlYsGEXRvSzNrnoHausl25/QyOe/yrVeen6ZbOLOSPIRXRju9I9bUSw80sPiibGbBfDzey2cUX9r8e/sC7bVC7vbV+6uQ6dykrRu0PLjO/scgi8bvIAY4xPAb/0zTrcf+HRruZHICW/2WWTF6/XCdUReA1hY2o11z6fmdjJeIqsWyykI1xONg58bpg4J7nDzcjcdXx2ujdMi1TmjSHOeNN8urujXaRAXaEGtF4tgHwBXq5eh5rd9aHn0LGrz7wr0g7G+EfsD0yyHpTMX2+dC0EvV3QQdcnfv8Qpf/jE8ju7ssyLxCJP+C2OBd1EoYlSRUbeRs+ZMPYuZaViF5kiWzW2IAL/i3jVrHIJ2s/jBuh1AWbog3KDaVmZlPRbscAm6Yl++60WQsPk2udnJ23idth5/tjxxCcrcdWz1VLc204f0Ml3GTrum9wyN/0kjvvHOJudvXaHlE7PTqk2mKYHIu/6bYKpDkU6qvJWpanrQuj2c9YUw3uOrN5zwG8/lFMQ+EwEsvyc/cphZVLSs83c925mkhIgXiESvlnNF/ZWBD/5Zo333Lx93466+gYpAc1Ei3BN1KKhz9wamxguejwzYUXrZkXYrW3K4w0UllCqmRKb1y283JaNLp19UgbO8j5bXmvadSsskjBZOWJ3gqeHTj/F/12WraZ4fPGttpDLJij1y9Nh2CUll43fCJiy8dLpHXPvR/jjR+4bfpzrFV9zsPK0Oa1/ZoYrJxNPs+KCtE7pua++sz3XCG8sFi+T848c8tEu21KXfEK8I/bLn/46TS95dUoQIfcUu+BZcUg8a4anMRoVwKkWL1Oc4fl9Z/z5P7j9DbGpsTdZ7IUxJx7hZcPO/aHPSdzMk6u37cXJ/zsNNbvrLUP4MsaEOpWvVm2z3GDWuXVmyAEnM2UhUdpAhjdTEa9uNJbN+6473YVVtXuS7XfBBn5zqLGZqRG7B7iUouS3zqm8yq6thcvj2j3LgG/WbMcZj3wa2Og2KK8RXgUy8et1wQhgIMpZiRWbbFwj/TJp/ias37Ef7y/YZKlYGIPlzbDr+GastPba0U831uF0LwoKyNNu8TA8S6wwjtKXbeFf60s3/6oRuyU92zv75LoRhFcMT13c13Be9Lv3FmF5zR4sF2hgQnIEpPbC2KCRraS5tgbUOK0UEoOglxKzbh36UaPyclLchQXkrRPkHX17KNup6IRlX7xQ4z3ICndHIupJRNOJaBERfUtEN8kQzIm+He2zwTAwPDplmeP12/akFrc+WiTuY2xVp0xEG6PsnJJe5eDFy8affOD3HywO3IRjV35jE0vGPuE53/a4vnnHVLYdhURyQgzY8MwXq5MLsjzrKW/MXu+ycc6/Ug5jriFjxH4IwK2MsUoAwwBcT0SVEsq1xc1GZbSNXnFSRcb3um/7jBXyElXLRLSjiNn6X3YRo3vnlrEqaOas3cF97kYbk5HV5h2ngYd5IdHqVD8p9f4weSle1nYi8wwo/vuVebZmJj+kb4qUXnwGvhU7Y2wTY2y29rkOwGIALoEp/OHUa5obhtPi2Ocmxb56614M+Z14MDDZI1seRW08JbAReyClKqJEVsA6193QhlfUqT2bvVv+bDHb/uULszOOibj77hV0DXY6n+DtfQ8jv60RqTZ2IqoAMAjATIvvxhJRNRFV19Zm7qSURYZidzjX3D9M/HqtaxQ8cXn82eN4ym0IKEN72I0xCoLrFL2VG6WTlqjEVrdOT+rCbYrhcBGxMpcG2TSdyibyVnfYs2ppip2IWgF4HcDNjLGM7XqMsQmMsSrGWFV5eXjueXYPYfbaHdIWMayCObnV7wTPNcZZy0mHeQuDqwBe9Zhdyg2RmDx2hK3kZSpL46tlF1oCCK5jtYI3Wqlbp+xF4jB/JyBJsRNRMRJK/QXG2BsyyvRKZsYX6xv6zert0l4cp2mpl75DNJWdU8IKhTPGuOfvcqbU42Hsc5kxingw7loNQhUw2LdJ3h2XPBi3zTttUNrkMbuXHUP7tLeXSYInza9enIPFgmEmgPDNmjK8YgjA0wAWM8Ye8S8SR50O3200NRTHKZDpSQcxQvLSUfNE9jPKGpRbYh5YYtK4YeKcqEXA8oA9hhizDzH9N4uIiY42eYd6iNPGLhurHal6J8Mbo8VJXK/xjbLRxj4cwOUARhHRXO3fORLKlYLd/XxvfuZGjbjoMa5GwLkBxJccsbkjwZFvnZco78zb6GmEany1wjZD2CIjELxHjLHJwujofAcBY4x9jmjXfFywvosLNuzC6CM7B147z+p93/KWaXZ60ReBdxu2KHF5H4MkNkrHgknznXPZesHLz7XLT+o0AEnboBTikD2ugxFjO/OaqlKEnDfOOt1D84g9iN7JLSzsOcd0wbRbT007xvPch/XtgIVaWNywkzLkEjHW64HgpSPjSRRthtcrJgziEA/KeAcO2WSgkknOK3anHtz8wKNoflY++TwvX3FBQeBTunxQeuZ73cUiiFUusffAIZQ1E5uoe7EPn3FUZ3x/cA8A8dlAx6vfDwUwoja2sw6G2OxBkZWKXaQHdvZJjb4rt1rsiXqEoxMHKeostrnLxHyrRZVetvF/01agZYnYb/TSHB/6r2MwSksUstAm4UoQyBiMzF7rLYmNE2/MTpmzrKJryiYrFbsIIs85CjVvtT9js08XsN4dWvi6XicOG5Qmf+s/lo8T5ny4MejrA0c0bKzdrNepeTQrLkzey16S2iMPbpuL8oWcV+xxXhwDrCO91XEsuDoGgiqvAAAgAElEQVSZmGT95J37GvDIx8tCXfwyM2NlPOP5ZDO1ps7MDbvHv3GX8/qR3oFE2X7yldyedwLxsCc44DWbShj9lZ72bgjnjr0gCCM/pBGRGNvZSoPg4p3dzM19ETDx7HiTZMvAcU0tzs57ksnKEXuXNvwLXMbHPG7MgLTvFm1M98/9h0Vi5qDxmk3F6VWpqZO7my8OJhlFdNitJa7fsc/xOopgxO5UlTLFxJzSokK05MwCbzTFtG9RkvbdpAXy/YRF4U0IbMZJ2dY3yF3VN8rYubX1iv5R3cQzRencd8FROH9gN8vvXKMIKgLHb0CzMJ0BrN6LODhJhE1WKnbG+B9W2HGQxfGo2CVL4QRP5+Pn3rYsLcKoAdmVtzWfsAseOn2pc5RWff0oBLftJCLRXMNElkMDL9mp2MG41WHYKalCI6L4G2XNiqWXH0Zy36D5x+VDohYhMLya4qIwxcTVahjGpiQjWaXYv1mzHa/PWp94eB4yjns1e/glCBt1mFuniwpT962iQws8/P1jpZbvtb+9aFCg+VwUGl5b2q79iT0IE79ZK08YF5zetS9WRLdD+2BAORPsyCrF/vbcDUlPDe54PoYHHdWA/enPV9t+N/Frb40+zJGJeaYzrG+HjHP8eBwUEOGcY7oKX/fID4/zXKdscmDSYUttnZh7pM5bcxNhkJ3yFYQBwVtYBJn4Se/nhaxS7AQSHv2GEG/HlQcmLZZe5l1vLZReph1FBalmcmyPtpbn+HWJLC0qxGe3nearjCjJ5QU6cwpJXooCmCG3beFsCrTTDlHvZxF1MfVLdil2Sjw4xhj/4imMI/bsffkqu6Z7nchOUODE/oZG9O7QAv07l+H60w63nPncee6RnsvXZwQlRXKbY3FheM87e1tWcARh+nSy16+oqYutjV2ZYhzQE8k6ZYExY2wHZuWYTQzu3RYvXn1CJHX/4B9fgjHgyK5lKCwgy5en2EcWJ/1ZylYDFw/pmfzcVWDvgxfCfnGj4C7BztsqDpJfnPT2Pz5dhQUWcWniMJ4L2xSUXYqdCA2NTfj3l99h5z6+4FBG0002p5BjDDjp8Ohym67dvi8wt7UCB83eqcx7JLyRR6TuV4dWJQ5n+kdfKBThqZ9WBSBJcIiaZIIYsTuNyL/waDIKk44Bt0OdrNJ0RMA+wSTBuRKmIg4/Q88JKtsjh5L/z1QELTg3olmWaxiqxWE7uTmq3+jK4BO9yGTfAbF3bwVHmr9fjOwrVKaTrdypVcbFRPP+TSNCqSe7FLuHl9M4Yo9rdpVsQ/ZLoitgqwGenzgjxil4mNPxcptZxq1n9AtPiAAQDVVxXE/rhXYjI/vxb0z7QVUPx7Z3yGEUF5d3v1NZOPH+s0uxe3g5c2bEHqPfIVsUXaFbLW43SrL/8JruvFBSVMA15PCzDhEHRDvZYYe1l1r/XedVOo7Y7dwy9bW5fEJKSyOis4loKRGtIKJxMsq0rMfDNcaGkG8PNyhkb7jSFbrV83UahbmWa/i8drtzwCo//OuK49P+tlMwUW2Qk4XojHk4x5qQSIlFNgv3PETt7ggAV5/cJ7S6fCt2IioE8BiAMQAqAfyIiCr9lmtdl/g1Xp9nt4C9KHR+ckKvUOoJg0d+MNDTdakRe+Z3fl7IsNxbS4sLsHjTbtfzvIZojguit1P2ugaBPLeH6NU60KNd89DqkjFiHwpgBWNsFWPsIICXAFwoodwMvLyoRtuayMPtHtJDqOjQkvPMODTNBHaD6D4deX9LOuRkipE0Yg+aEyx245oJwv0vDjQrtlYjPD9XxLyj72MRhcj7SF8mYXbsMhR7dwDrDH+v146lQURjiaiaiKpra52jwtkhelsKyPvO09YBBLuygvddj0PD1JGd7DdpirG4F7vr3bNJ2Zfr+VLRmtCuhbsbW5wG7BcNlhdnx0+YaN3Tigcib2bAPQcO4cAhMY+eIAjz8Ye2msMYm8AYq2KMVZWXewzRKjoVpPSpm0gGeqfmM+KI8P3J4zTYkx2pjkz/NyI6Yj+2Rxvf8ohCBAzt475QGCcb+5ijU7F52rls09cR1ak8UR3rBDpuAqFrG/GZ9B8mL8VFj88Qvk42Ye58l6HYNwDoafi7h3ZMOqI2O0K6gm4u4BPtVNP5x1onhfDD+zeOwGM/Hiy93CCws3N6bbgFyRG7/4af5rse0ntE4FPashT77aZMYFbcf+FRjt/PWbsj+Zn3vvMuQLdulsi4yWNmEeksCggY7DEm0fodzvlZwyDMwZkMxf4NgCOIqA8RlQC4FMA7EsrNQPS98Dp1AxJZ1sOksltrnHusU4TD+Iz2/GwaskJv8F70XqkpvkwUd4lXMcoasV1S1dP1nB8NdV6U32nYKSt7IvH4TxKx6XkmWyKLodkc6ylsfCt2xtghAL8CMBnAYgCvMMa+9VuuFV5W5b2svZ3avxz3uYx48pnDO5VhaIU8H+XkiN2DWjaHEE7blGQqL6gsNrxSW43YJ988MpD63BZqjd/yKkxek43+M3nSV4q8n9mu1sPc/SzFxs4Ye58x1o8xdhhj7EEZZVohcmOG9W0PkDd3uXvOPwodWznEKJH4fPhHe/LqlMGTP8uMc+JVxJRXjPi1T1yWbr4qSNfsaXz6m2DCAvPKbaVs+3cpyzjmPHPjq8/NA8N4n3hvO++bpLfpthwLyiJ+LnF7B+JMVm2FE3mweifgx5vkmSuqcM2IzE0FHy7c7L3QmPDXHw2y/W72b89wvd5Kb3h98fyMZFqUFKX9bRwV9++cqTSjpIDzbdNt1Ha43S+eNH3G58cdKZVzeB3UGjER4baz+uOwcm9utflEdil2kXNJ30osrtn1ekYN6Izbzh6A7w/ukfb9tr0Hhct0q0unr0dfcFHOH9gNr117ouV3PP7WMvPH6opApEw75dezXcrc0q1tc9x5jvc48WYm2ChMXrn5F09dznP5Wo9D7xSrxThT5JWf91US8dcWfT17tm+BF68ZJnZRiMjOKeCVeEjBi8CLT5T45zdWTHFhAf5k2lHpNX/Db87qn7HwaP5Jdi+jaJV6XlKnjDNVFe0td74SR6uwUgZeR95OfuxWfHjzCEz79ak2ZTn/7Yczj+pieZz3Zea9P7JkdtKvQZpiRDpotzKt1nLibJI5rb9HV27JZJViFxqxa//pNvZnrxzKf21ADad9y5IML47MyuXUNZzT195K8fOM2K1O8bpbNxlSgPP8AV1a266B8AwWhx/uvkvUDquNPbzBvXgHsjNc4oq7PR6eDiTdFMMnGO96lYgpxm1GbSVaXHfw/vOK49EvJua/7FLsIjZ27Vx9xC7iouf2YvC+COZFMLNfvX7MiB8Th9WW/qBeASsx27f0lkTAKaSAKG73b+VD5+D5q7xnompVmjAB3XHOANww6nCUNSvizs7Ea6JYs83ZX1zGMzXKInv3s0gbLnRZeLBU7AbZlz84hrsuIwO1mfFxPdt6brdmTo5g46Id2aXYRTcoGfzYZSo43rK8hCWwK5vnXRlo2HXJb821MKnweF1IHDU5xWM3c4ZLcgo3sQoLyFcHol9ZXFiAW8/sjwX3nsW950HWPZOykSutPL5reGOai/zOO85x32yVUb6hoXgNhVxACZPec1cNxdd3nO56Pk9aQLdfHWZM+KxS7CJTPKKEytKnj15G+16/FznP/JL6eWfvu+BoAImUa7zlvDknc5Mw31Tem6BWESD1kngU1mHlrRy/D3oTi15+lLF7XH+hdoKTiF4WT3nXq0QeQV+X52mFrA5yQJfWKGtWjCKXzmFAlzLXDV9xI6sUu5i7Y4LUC+ivMfzue0cb5OBcBLM4ZlYIom30uaus1wo6lZWiTYtirBl/blrKNTdZN+zM3GrNM7Lg6WR1s4URq4VGscU2Z9nMcu054D2ImBXEoTRtrzV89hPTxt3GLlbGbWfxjZp/flIF13lBm8Bl2NhFO2beQdqkBZuk1emHLFPsol4x5GnEbsXlw3qnyvZXVBojjkhfRbcbLevHj+hkvThjnpKm/PjFWxPPJTzP4p7zM8Pyp+LCGMviFs0V8/17dMpyeYXD3301/s5XfmHtamrHTacfkSGDG04iGjtAtw1ROredzdcByNxhaVUW734AWSzZXMc9+FhVuzdgafjIKsUuAkH3Y0/9zX2tzcn9OiemjdeeepincqzKNS94yjYDRRlR0Er56+IwBhzT3cOoNeLwxWSQn4f3bjg5da3WCstKi9Ls8pNuPBljXZI69zVsypHRERqjO/LQrU0z7ra072BqlnTtKXzvygCLHbgAUGHhECBzfSdMwmy6WaXYxezkBGNIAbHRvvW5+i7HNs35FkXdRi5f3Z65aPMzm+lu0nPEri7TFw2NiRjZJR4Wl6Qt8rmULaokAeuXw8lXXzbfbUuMyN5faD/lNnJ09zbJ9qInpPixae/AUd3aYPSRzovCXnC6rV0EM4QZg4aZMceQqTOYv3i80d6/cQReHms9g7Hq/KNwd+Sp0rXjC9EWk12KXWDcvUtriPqtFGkKds/nf84egK5tmtmOLtwwy2/1crllIbLrdH5yQu+0vw8e0hS7h51wzUsK8Y/Lh+DIrq2Fr9WptLnWamOM3+Y+9+4zk5/NNvhxHCFuRVi4IZECb87ancLXlhQVYPmDYyxlkhlHnmcQI9p57zton6iitykLWFFB5jN2orJba7RpUYw/XmKxsG5RQNgpBk/s2yHrZgnZpdgF7m1BYsDuycberMh6lHHiYR3w5e2no0VJET68eQT+culxjuU41Wketdldw7uT8job81AngeQiRs46qgs6tvLu39u1TTOblzLx/1alRckfI2Kvtusw7OBPPZjO4Z3EvTV4KC4ssFS8bm57xlvE3ZYd7muQaspoZhJRwhcP6ZFxjOfq60/LbPtP/ERebgM9PIlflCnGBtGbS0RobNT92EU8L9wZ0KU1LjzOOb1YRo2GA2NHWNtUzSODq4b3SbuU91f0LW+F+y88Co//ZDBevOYE/P0y98BQMilyiLtw17lH4q3rT/I0Yv/eIOd7ruuy4yu8JWTQ4TV5yUJoDYjz7IOSM13ZYa7F2En5vV92139w04hkrKNyi13IbSSa5wrI396HKHAOIxczRO7tN2t2oG2LYmzcVQ8gkf3lGM7pbmALjsywYcqmCrcpn0gD++mJFQDgGIK4vKwUtXUHuMp78Rr+HZvGYFxmrtY6NVEbu5X7pBVrxp/LV6ANbVsUo1+XMiyv2ZPxnZcpOc+MxK1YLztFF2/abS8TXxF8OPw+0ft12bBeeP6rtcm/CYRjurfBgg270s4zmgmtaj+yi3czoplmxYXce2iKCwkNNh2qcne0wY8bVU1dPfe5vIujfrD7LXYNKBkoS7IcIusFZaX89+XnJ/exVEBpJoXUUa4yeRJlZIRs8HjDjO6FsnBqv8YO27iD+JIhPVDZtTXGHJ0KQOb2k8IeWzrNDERlOd4c9IuAV689EXPvdg8lbaRdyxLH8NM8LU5319Vdp3kw52Q92yZwXNDk1Ig9EULA8DesP0cGh++2205U2TPCoGYnhUQ4qpvzDElkF+eTP63CoF72YWh17MpyC0WQJhfk3pfKbq3x1artjuYpI2WGUBTfH9IDf7jEnCWKr5zCAnJMBt6suABXn+zsZsmDeUZilE50xJ7R/pEYMTuFbbB75sVO99ul0U2+eSTqGxILxkd14x/9m39u/y5l+PDbzVqVyivGkk27nEfdTo0o22xkOpkNXe7vcLpnZsUsGpahX+cyLPnd2a7n8jT3Myo7O2e1cqtDdKeh55oymfDTKrzyixPRktOU5CaH1bHzB3bDEFOiZ6f0hYwxLPndGPz6rP7CMrnBu/nMyoHA3J/yvLcinfAJfdq7C6Z9PbBnW7zzq+G4YZT32VtU3jRZpdjnrXN2MTM/3x37Ur63Udxfp0iOnuWR/DvKtIQVZRZK58LjunkuVxfTaaSVXDwNcCDj+TYbHtD4i47Bi1efoB0XL6t1s2IM7cOfI9bosmml2Kxk+OuPBqF5yAnYdVnMMrYzpMSzU8x3nnMkHvqvYzKOmxUhz+0+j3PnLJA5WLFDr/fYHm2FOg7zmUZnp6zxiiGiPxDREiKaT0RvEpH7XNlXffbfHV/RzrF3j6LfdFJYXmcQsjsoXRn85uzMkZsvM5DNuVZue3KnqHZl8ddhDjR26dBeOOnwREhWp3vw1e2n48vbRyX/lhE4yqo+t7bj9PXr152IW0b34w5V+86vhjt+f3xF+7RHXVJYgJ7tU2shohatDMXOcX2RQIwBpwVl0XqtyJTfsJKURYunHwM4mjF2LIBlAG73L5I9Tg364YsHOjaioX28J1fwinmzDBEJb5gy20iD6qBkZkQCwKVHk3FXvNfiXoeHN/Spn1Ulf7m5TTlNrbu0aZa2eFbiMdXWfRccZaiP/zrdFdQpAuaQ3u1x0+gjuO+LmyfSwUNNaUpwQNf0xXhRU4R5dMxzuZeE9Wa+ve8sU70e237Aa2K8+FLsjLGPGGP6/uGvAGTuMJCI8R5V3zU67buiAnJURO1ahrftXMeqvdXVJ24X7wPXFXsypICPluI0ZeUp1quiN+Zx3bY35Vp58xlHoHWzIlQKLE65IWNU1LZFie39FlFUXmOFH54W6I2/vouH9MCa8eeiW9tE5yJDqbj93uJCcqzH7ju7KJ0ZNnaO39+iNN0EpYdusHpX9HrNs0Tz+ocsS6nx/mWNKcbElQA+kFheBsbnZLWQ5tjAIjDGOD3IQ5ybR8yjET+/wipnp+iCqBeM+UnnGtZJTjqsI+bfe5anhCRhYZ4xuaY2BPD1HafjsPKWtnF/RPCjnGV0cm6K3WwGyQhL7bM+nt9fWlSYFqZBf9etntUBLdSGU4gEKzl4MXcmUaXxc22lRDSFiBZa/LvQcM6dAA4BeMGhnLFEVE1E1bW1tZ6EdVLOjLl5xXiq0hdOL5aTG5qRQ2ZTjMXvaMmZ9s+vLVvkHtp5sHgdxfJivl9e+c4mPd2wvgmT3o2jDre9tlPrZph666lptmavhNVsZ4wbham3npJZv4sARYXpM+WMkbipgGF9nReRvdu2M49ZtTW9ebilM2zvEE5DLOFP6nOs3B0ZY6MZY0db/HsbAIjoCgDnAfgJc5CcMTaBMVbFGKsqL/eWydvpoRcUOCdViIuzo96gnB7x+zeOSC5upVL7Udr/dW4Z3Q/Tf3MqV93WTyeYO2MXI0RkoUsE/X69Nmu95fei75TdIEG3AZeFNMvwlcZP4NJubZtb2uZ3OUR1TNRBJuWV/r25GbiFas60sfP9CKuooVbonnVuni5Os0jHzFSmv4MeyNjh1yvmbAC3AbiAMeacgVcCM1Zus/3ObZXfqYGcfLj3JLTODSTTjMIzla/s1hq/1IJ6mUe+5p9xfJ926FTGF+jLaZHJ6qsMrxgfnYB+j3u2b+5ypjfswsMOP7wDhvZpj9sFc2vaNRf9ectYsOOSw+a4aDA0rxw45GyyYIylyZhpihFcPM1Y0+C8znAij3tiUCYSc7FRKXa/O0//BqAUwMea4vyKMXatb6k8oMdKt8PpMfrZZeio1i2VpeYJ4qIYrjq5D64e0RdPfbZKu876PJERsPXv5I9+yXPO13ecjvqGpozjh5W3xOcrtroXIJkWJUXC2YoA+9+qH5dk8UnjyuF9cNoAU0YtGzne+dVwNDKG/nd9KF+QdAk4TrFfINy6hy8Okc7g3u1wSr9yfLqsVqud793U8w8A7jkIzqzsjDs5klN7wSyv4+7XAPHrFXM4Y6wnY+w47V8kSl2ne1v70aDe9p690jpnqFec3m+n/KZueoF3CsrrjwwA5x6T8Ioxxmhv0zxxvXFTia1MHHV0at0MvThiukSNa8hl21g+wY3Y7z6/MiNVop1JqKiwAKU24aVl0q2tWNhn84BFNOdss+JCPHvl0GTYZN6B9dtzNyY/O42SW5QUYsJPqzJiyANIbkJzw/gT+3Rsicd+PDgZS94s79Y9B5Ofy8u875wWJat2no44wrvJRH9RT+mXad8P0w1JtP826w9zoxWJG15UWIDJN4/El+NSm2huHn0Exl90TFqQKZ2+Lkk/Vj50DnfdOnFZ63ALuWw3idOn8E1BDNljiHng8L3juqFHu/QBlNMz9WrxSK0tifOHS461/a7IYXZ+kgeTLCGRM9YqljyQnibwgoHed3KLklWKvbUp6uKEy9NjjDs2oii8Yiy6jCcuG4IfVPXwnACisIDw+nUneZapf5cydDDY7ZsVF+LSob0sZwhFhQX4xOCqaD4lynyqQdC3Y0tcobko2sV1KQjQFGOFH1Ow1SBGFPOM4dFLB2WYtpwWT72uyyQ38nFeblTY5llP2nkB27zNCW+MM7sw41VllWI3Y/bLduqNg7qnTrZyq6/6dynDwxcP5FaKUccuS08mnFuK3Mjcu8/AtF+finu1XZ92Ef0uHdoLfctb4pKqQPfiJfGzYD12ZF98c+do9xMd68+km8HkyVj6OeZNcLYblFw6xlW1idyyizbVcUjJn6npnGOcw+jee34lTjrMeZd6Lwc31guP645Zd43GNSP6AAhvAGAmuxS7y01yWkCNQiXJig1uJsh48Ved3Mf2Ozv5/XgV8cCTENnv7KGtaY3BbnTVrW1zTLv11DTlFiRubebRHx6XMXNNXUsoLyvFned4Xyi0s/FfNDhhymJgyXv17JVD8SuTf7/5aj0XL6+3yFotebgbvF4ud51b6fj9FcP74MVrhjme8/DFx+Lm0VrER4tqO7QqxZ3nVmLN+HNDjQ9jJLsUuwWn9k9Nu8xhS43wToP8jnCMyHyoxrKCyse5Zvy5+O15zg3fzIxxo/DUz6q4zw+qnes2c6cF9GzErdl+b1B3yx3FRq4ZyR9z/V8/Pz6tI7Cr/4wjU/Ht9VNalBS65hO46uS+OLp7a+7IobzthXfELmNw1ay4EGdp99ytuDA3JRnJesX+5E+rsFAL4PPb8yptN0DwPk/RleugH1vSPTLUJV5r7KblTqF5vfK3Hw8SOn9Qz0RgUS+d3oc3jxC+JizCDoVxav9OaR2B3YAoFZkz5WVldaZ5xN+rQwu8d8OItHUeGfCO2P3cT2NmK97YTWHtdzCT9Yq9uLAgGYGupKgA/W1SvTnd/6B61SCV8TNXVGG6YWEzDPws/vz4hN4oKSpwHV3qnHdsakTHU6t+r71YZAZIzI8pm/iuT6cEe+QHx+HaUw7D4F6ZM2a/C4a8r+ba7Xz7I2Xdz9ROcGeisrFnVWo8HuwaQhRBwGTodV1q8+8aNYA/1Zss/NzB/l3KsOyBMdJkMaPfH9meB98fHM4iqR1RL57rmF2NCwwj9i5tmqUF4TISlvi8sZdktQ/GubHPuHEqTLJKsfMsotkRnFeM9fFWpUV4Y84G3+XH5cWOO0FMulY+dE4MRszRCPDLUw9LmlG+vH1Uxga2sEyEvKXzLp7LeJ7GGb7bgFHfBcuTiF0mWWWKufmMfigrLUrLUiOD84+Vv3Hg6O7xnd57JapOhmeUJZrAhIfCAgrV99iKqKq/7ewByXyoXdtkrqPwiuVXfl4zKb/7sA+BPGRDOqiN2H8uIYSzCFk1Yu/etjkWmDKd+KWya2v84PieOKJzK0f/VDuIMh/ya9eeiH5dynDsvR9JkjIeRGLOEiTXZjhx/Tm899lvm+EdsTvtYYmSsSP7YvGm3a47nWWTVSN2P9g1RH2KNKhXO2kr9VUV7S3DfnoxF6Rs7DHwiono3eEZjcXh/gRB1DMGO3Qzjdtt9ys+r7dLGLugR2rrDOVlpYY1HedrerRrgVevPQntBGI6ySCrRux+sBs5fLBws89yw8nZmZtqK5NT+pWjviE9VKxb3s104qkIvSLr17x1/XC0KpXolqovnvKd5pkbHBKaGAlDsd88uh8uHdoL3ds2x7f7dwVenx/yR7EH9NzJyhaTo4Tx8lhF33zxGr6oewly61l4TdFm5jjNz18WYdnYeXdZF4UQHrewgJIb4ILywpJF/phiQqpnwb1nSi1PT8zbPIBNQKLIUjK8nNi3A47t0cYyxKoZo59CLhFTvZFSaK6mGJ8/gPNyO5PNaf39B0JzIqaPJ59G7KlH8J/fnIaRf5geSD2yU6Z9f3APbN51ANeMtI/hEhZhr09NHOscs8MIr81TIQcvt9kupo0TvIOJqCKNxrW95Y9iN3yWmQgi6OdaVFiAm/SAQxHDG48jGrzH744zRsUx/denYs1WvqBYQZNcPHUZshvl5911nHY953nHV7TH7LU7M44HZZiLu/U150wxdt4Rce1Zs4mg8kTKQI/O6MVlNc4YZ5p9OrbEaQM6RShNCmOsGMfzfHa1vCP24QFHGDXTunliTHxkSLlnRcm5EfuBQ9ZbeOO6yJFNxHnEPqxvB0y4fAhO7R8PxSeLuN7xsBZPea8PMpS1Fb07tMSL15yAqt7tQ62XFykjdiK6lYgYEYXbbVoQtq3N2PD0nKJ2xH365kbcMyadeVSXtHyuuUBsxyOCcp3vMS0c74i/0iYxSpCcdFjH2LY331IRUU8AZwJY618c/3TUNhndMrpf+JXH9SWURMz1ek4SticSL7rC5R2rdGsjlhQ7WQ+nhtLvkp94UrmEDFPMnwHcBuBtCWVJoyXnZoyrHTIGKdKJq5LJZeJ6x1PRHYOdhvL+/qLCAvzmrP4YfWR61NNsnyV7xZdiJ6ILAWxgjM3LVhu27ID/uUxpTKedOU1MXyvR992rfhUZTFx/Gt8u1XzAVbET0RQAVn5KdwK4AwkzjCtENBbAWADo1auXgIhihJ1p6AdVPfHCzFhYoQInWztvO6b/+tS46s0kcQ28xrk/ybf0OdbkQsNVsTPGLJOAEtExAPoA0EfrPQDMJqKhjLGMACyMsQkAJgBAVVVV4NqXVwn57Qjuv/BoHNO9Dca9scA9/6GvmhSy6dPRfUdr1ESVWs2NsLxi/Jr/4nn3gsezKYYxtgBA0pA1nScAAAndSURBVLeMiNYAqGKMbZUgV9ZQWEBorhZsFAHRIeSogLyokXS8yTmjqegAJ8gB0REeEisrFEaKCuP5iiajjgY8JFYL9t6Q1moYYxVxGq3zNoemALPNfnTLSLx+3Ym4aFC4QfYViqBJ5gkIuh6l1z0Rz+GAD3SbJO9mGhl63a5OIsKQmO5MUyiyAaXXvZFzIQX0bOW8299leNGMGtAZg3q1td0UVaqF3I1r+i5F/Jh3z5mBzialEbAtJu67neNKzin27w3qjhdmrsXJnEGBZLw7bZoX481fDrf9ftyYAWjTvBjnHuscckCh0Ak79okown7sHjsAvy62uZoy0Y2cU+zHV7THmvHnph37jZZt3YowHnyb5sUYN2ZA4PUoFGGTn2oz/uScjd2K7zksXsbVT1ihiDOpJOtu5ylTShTk3IjdCqumNaBLGZZsrpNiisl1bjz9CBzXs03UYihihPJWiTd5MWK3aoQXHJcII6pG7O789xn9MGpAZ/cTFQpJnB6ThCLZSl6M2K02ORTwBrtQKBQZ6CaW3fUNgZT/xGVDsOfAId/lNObplDw/RuwOx9SIXaEQZ/2OfQCA77btC6T8kqICtJcQTqFvefzjAQVBXih2K81eENKWaIUiF6mr5xtND+2T2KB3Sr9oTCt3nVuJF64+IZK6oyQvTDFWK/O6JaZRaXaFQpiGJuvcwmYG9myLlQ+dE9lGo2bFhaEnuo4DeTFit2pTx1ckRhIjjygPWRqFIvs5qhu/l5TaPRo++TFit1g8HdizLZY+cDZKi1TIXYVCFD2HaXmZykAWR/JDsdscV0pdofBHtlgyx40ZgIoOLaIWIzTyQ7GrmaBCIZfkO5Udmv3aUw6LWoRQyQsbu9rWrFDIpURLANKrff6MgrOJ/Bix50X3pVCER9sWJXjyp1UY0rtd1KIoLMgPxR61AApFDnJGpQozEVdyeixbXJhQ6X5jOisUCkU24VuxE9ENRLSEiL4loodlCCULZVtXKBT5iC9TDBGdBuBCAAMZYweIKJYh2ZR6VygU+YTfEft1AMYzxg4AAGOsxr9I8lGWGIVCkU/4Vez9AIwgoplE9CkRHS9DKNkok4xCocgnXE0xRDQFQBeLr+7Urm8PYBiA4wG8QkR9mUUiUSIaC2AsAPTq1cuPzMKoEbtCocgnXBU7Y2y03XdEdB2ANzRF/jURNQHoCKDWopwJACYAQFVVVXZsV1MoFIosxK8p5i0ApwEAEfUDUAJgq1+hpKFG6gqFIg/xu0HpGQDPENFCAAcB/MzKDKNQKBSK8PCl2BljBwFcJkkW6agBu0KhyEdyeudpixIVllehUOQfOR0r5tVrT8LUxVvQrFgpeIVCkT/ktGI/vFMrHN6pVdRiKBQKRajktClGoVAo8hGl2BUKhSLHUIpdoVAocgyl2BUKhSLHUIpdoVAocgyl2BUKhSLHUIpdoVAocgyl2BUKhSLHoChidhFRLYDvPF7eEXGKIJlCySWGkksMJZc4cZXNj1y9GWPlbidFotj9QETVjLGqqOUwo+QSQ8klhpJLnLjKFoZcyhSjUCgUOYZS7AqFQpFjZKNinxC1ADYoucRQcomh5BInrrIFLlfW2dgVCoVC4Uw2jtgVCoVC4UBWKXYiOpuIlhLRCiIaF0J9zxBRjZbTVT/Wnog+JqLl2v/baceJiP5Pk20+EQ02XPMz7fzlRPQznzL1JKLpRLSIiL4lopviIJdWXjMi+pqI5mmy3acd70NEMzUZXiaiEu14qfb3Cu37CkNZt2vHlxLRWRJkKySiOUT0Xlxk0spcQ0QLiGguEVVrx+LwLNsS0WtEtISIFhPRiVHLRUT9tfuk/9tNRDdHLZdW3i1am19IRBO1dyG6NsYYy4p/AAoBrATQF0AJgHkAKgOucySAwQAWGo49DGCc9nkcgP/VPp8D4AMkUq0OAzBTO94ewCrt/+20z+18yNQVwGDtcxmAZQAqo5ZLK5MAtNI+FwOYqdX5CoBLteN/B3Cd9vmXAP6ufb4UwMva50rt+ZYC6KM990Kfsv03gBcBvKf9HblMWrlrAHQ0HYvDs3wWwNXa5xIAbeMgl0G+QgCbAfSOWi4A3QGsBtDc0LauiLKN+b7BYf0DcCKAyYa/bwdwewj1ViBdsS8F0FX73BXAUu3zPwD8yHwegB8B+IfheNp5EuR7G8AZMZSrBYDZAE5AYjNGkfk5ApgM4ETtc5F2HpmfrfE8j7L0ADAVwCgA72l1RCqToZw1yFTskT5LAG2QUFQUJ7lMspwJ4Is4yIWEYl+HREdRpLWxs6JsY9lkitFvns567VjYdGaMbdI+bwbQWftsJ19gcmtTuEFIjIxjIZdm8pgLoAbAx0iMOnYyxg5Z1JOUQft+F4AOAcj2KIDbADRpf3eIgUw6DMBHRDSLiMZqx6J+ln0A1AL4p2a+eoqIWsZALiOXApiofY5ULsbYBgB/BLAWwCYk2swsRNjGskmxxw6W6FYjcSsiolYAXgdwM2Nsd1zkYow1MsaOQ2KUPBTAgCjk0CGi8wDUMMZmRSmHAyczxgYDGAPgeiIaafwyomdZhIQJ8gnG2CAAe5EwcUQtFwBAs1VfAOBV83dRyKXZ9C9EokPsBqAlgLPDlMFMNin2DQB6Gv7uoR0Lmy1E1BUAtP/XaMft5JMuNxEVI6HUX2CMvREXuYwwxnYCmI7EFLQtEemJ0431JGXQvm8DYJtk2YYDuICI1gB4CQlzzF8ilimJNtoDY6wGwJtIdIZRP8v1ANYzxmZqf7+GhKKPWi6dMQBmM8a2aH9HLddoAKsZY7WMsQYAbyDR7iJrY9mk2L8BcIS20lyCxFTsnQjkeAeAvor+MyRs3Prxn2or8cMA7NKmh5MBnElE7bSe/UztmCeIiAA8DWAxY+yRuMilyVZORG21z82RsP0vRkLBX2wjmy7zxQCmaSOudwBcqnkP9AFwBICvvcjEGLudMdaDMVaBRJuZxhj7SZQy6RBRSyIq0z8j8QwWIuJnyRjbDGAdEfXXDp0OYFHUchn4EVJmGL3+KOVaC2AYEbXQ3k/9fkXXxmQsZIT1D4lV7mVI2G3vDKG+iUjYzBqQGMVchYQtbCqA5QCmAGivnUsAHtNkWwCgylDOlQBWaP9+7lOmk5GYas4HMFf7d07UcmnlHQtgjibbQgB3a8f7ag10BRLT51LteDPt7xXa930NZd2pybwUwBhJz/NUpLxiIpdJk2Ge9u9bvU3H5FkeB6Bae5ZvIeE9Ege5WiIxum1jOBYHue4DsERr988h4dkSWRtTO08VCoUix8gmU4xCoVAoOFCKXaFQKHIMpdgVCoUix1CKXaFQKHIMpdgVCoUix1CKXaFQKHIMpdgVCoUix1CKXaFQKHKM/weAwoCD9+EqwAAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"hmc_states = lax.scan(lambda state, i: sample_kernel(state), hmc_state, np.arange(num_samples))\n",
"plt.plot(hmc_states.z);"
]
}
],
"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.6.8"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment