CSD_2D_model_v2.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load Some Packages\n",
    "import numpy as np\n",
    "import imageio\n",
    "import scipy as sp\n",
    "from scipy import signal\n",
    "import matplotlib.pyplot as plt\n",
    "import os, glob, ntpath"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "maxIter = 5000\n",
    "\n",
    "# Given the range of the rectangle\n",
    "lenX = 250\n",
    "lenY = 250\n",
    "\n",
    "deltaX = 1\n",
    "deltaY = 1\n",
    "deltaT = 0.1\n",
    "# Need to use CFL condition to verify the upper boundary of deltaT"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Given Some Coefficients in Partial(K)/Partial(t)\n",
    "D = 0.75\n",
    "\n",
    "# Given Some Coefficients in F(K)\n",
    "A = -0.3\n",
    "Kr = 0.03\n",
    "Ktheta = 0.2\n",
    "Km = 1.0\n",
    "\n",
    "# Given Some Coefficients in Partial(r)/Partial(t)\n",
    "B = 0.0001\n",
    "C = 10.0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Given Initial values\n",
    "k_init = np.full([lenX, lenY], Kr, dtype=np.float64)\n",
    "kt_init = np.zeros([lenX, lenY], dtype=np.float64)\n",
    "r_init = np.zeros([lenX, lenY], dtype=np.float64)\n",
    "rt_init = np.zeros([lenX, lenY], dtype=np.float64)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # The number of the random points you want to pick\n",
    "# N = 10\n",
    "\n",
    "# # Pick N random point(s) and raise its or them value(s)\n",
    "# for n in range(N):\n",
    "#     a = np.random.randint(0, lenX)\n",
    "#     b = np.random.randint(0, lenY)\n",
    "#     k_init[a, b] = k_init[a, b] + np.random.uniform(0,1-Kr)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 1.,  1.,  1.],\n",
       "       [ 1., -8.,  1.],\n",
       "       [ 1.,  1.,  1.]])"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Definie the filter of the laplacian operator\n",
    "filter_laplacian = np.array([[1/deltaX/deltaY, 1/deltaY**2, 1/deltaX/deltaY],\n",
    "                           [1/deltaX**2, -4/deltaX/deltaY - 2/deltaX**2 - 2/deltaY**2, 1/deltaX**2],\n",
    "                           [1/deltaX/deltaY, 1/deltaY**2, 1/deltaX/deltaY]])\n",
    "filter_laplacian"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def Dirichlet_bd( maxIter, N, clamping_ticks, K, Kt, R, Rt, top, bottom, left, right):\n",
    "    \n",
    "    random_points = np.zeros([N,2])\n",
    "    random_points[:,0] = np.random.randint(0, lenX, N)\n",
    "    random_points[:,1] = np.random.randint(0, lenY, N)\n",
    "    \n",
    "    foucs_point_loc = np.random.randint(0, N)\n",
    "    foucs_point = random_points[foucs_point_loc,:]\n",
    "    foucs_point_value = np.zeros([maxIter, 2])\n",
    "    \n",
    "    # Fixed the boundary values\n",
    "    K[:,0] = left\n",
    "    K[:,-1] = right\n",
    "    K[0,:] = top\n",
    "    K[-1,:] = bottom\n",
    "    \n",
    "    for i in range(maxIter):\n",
    "\n",
    "        if i == clamping_ticks:\n",
    "            for p in range(N):\n",
    "                K[int(random_points[p,0]), int(random_points[p,1])] = K[int(random_points[p,0]), int(random_points[p,1])] + np.random.uniform(0,1-Kr)\n",
    "#                K[int(random_points[p,0]), int(random_points[p,1])] = 0.9\n",
    "\n",
    "        # Fick's second law\n",
    "        R = R + deltaT*Rt\n",
    "        Kt = D*sp.signal.convolve2d(K, filter_laplacian, mode = 'same') + A*(K - Kr)*(K - Ktheta)*(K - Km)*(K + 0.1) - R*K\n",
    "        K = K + deltaT*Kt\n",
    "\n",
    "        # Fixed the boundary values\n",
    "        K[:,0] = left\n",
    "        K[:,-1] = right\n",
    "        K[0,:] = top\n",
    "        K[-1,:] = bottom\n",
    "\n",
    "        # plot the images\n",
    "        plt.imsave('./images/{}.png'.format(i), K, cmap='gray')\n",
    "\n",
    "        # Fick's second law\n",
    "        Rt = B*((K-Kr) - C*R)\n",
    "\n",
    "        \n",
    "        foucs_point_value[i,0] = i\n",
    "        foucs_point_value[i,1] = K[int(foucs_point[0]),int(foucs_point[1])]\n",
    "        \n",
    "    return foucs_point_value"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def Neumann_bd( maxIter, N, clamping_ticks, K, Kt, R, Rt, top, bottom, left, right):\n",
    "    \n",
    "    random_points = np.zeros([N,2])\n",
    "    random_points[:,0] = np.random.randint(0, lenX, N)\n",
    "    random_points[:,1] = np.random.randint(0, lenY, N)\n",
    "    \n",
    "    foucs_point_loc = np.random.randint(0, N)\n",
    "    foucs_point = random_points[foucs_point_loc,:]\n",
    "    foucs_point_value = np.zeros([maxIter, 2])\n",
    "    \n",
    "    \n",
    "    for i in range(maxIter):\n",
    "\n",
    "        if i == clamping_ticks:\n",
    "            for p in range(N):\n",
    "#                 K[int(random_points[p,0]), int(random_points[p,1])] = K[int(random_points[p,0]), int(random_points[p,1])] + np.random.uniform(0,1-Kr)\n",
    "                K[int(random_points[p,0]), int(random_points[p,1])] = 0.9\n",
    "        \n",
    "        # Fick's second law\n",
    "        R = R + deltaT*Rt\n",
    "\n",
    "        Kt = sp.signal.convolve2d(K, filter_laplacian, mode = 'same')\n",
    "        \n",
    "        # Fixed the differential boundary values\n",
    "        Kt[:,0] = left\n",
    "        Kt[:,-1] = right\n",
    "        Kt[0,:] = top\n",
    "        Kt[-1,:] = bottom\n",
    "        \n",
    "        Kt = D*Kt + A*(K - Kr)*(K - Ktheta)*(K - Km)*(K + 0.1) - R*K\n",
    "        K = K + deltaT*Kt\n",
    "        \n",
    "        # plot the images\n",
    "        plt.imsave('./images/{}.png'.format(i), K, cmap='gray')\n",
    "\n",
    "        # Fick's second law\n",
    "        Rt = B*((K-Kr) - C*R)\n",
    "\n",
    "        foucs_point_value[i,0] = i\n",
    "        foucs_point_value[i,1] = K[int(foucs_point[0])+5,int(foucs_point[1])+5]\n",
    "        \n",
    "    return foucs_point_value"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def Periodic_bd( maxIter, N, clamping_ticks, K, Kt, R, Rt):\n",
    "\n",
    "    random_points = np.zeros([N,2])\n",
    "    random_points[:,0] = np.random.randint(0, lenX, N)\n",
    "    random_points[:,1] = np.random.randint(0, lenY, N)\n",
    "    \n",
    "    foucs_point_loc = np.random.randint(0, N)\n",
    "    foucs_point = random_points[foucs_point_loc,:]\n",
    "    foucs_point_value = np.zeros([maxIter, 2])\n",
    "    \n",
    "#     raise_value_ticks = (1 - Kr) / (clamping_ticks*0.01)\n",
    "    \n",
    "    \n",
    "    for i in range(maxIter):\n",
    "\n",
    "        if i == clamping_ticks:\n",
    "            for p in range(N):\n",
    "#                 K[int(random_points[p,0]), int(random_points[p,1])] = K[int(random_points[p,0]), int(random_points[p,1])] + np.random.uniform(0,1-Kr)\n",
    "                K[int(random_points[p,0]), int(random_points[p,1])] = 0.9\n",
    "            \n",
    "            \n",
    "            \n",
    "#         if i >= clamping_ticks*0.995 and i <= clamping_ticks*1.005:\n",
    "#             for p in range(N):\n",
    "#                 K[int(random_points[p,0]), int(random_points[p,1])] = K[int(random_points[p,0]), int(random_points[p,1])] + raise_value_ticks\n",
    "\n",
    "#         l = 100\n",
    "\n",
    "#         if i >= clamping_ticks-l and i <= clamping_ticks+l:\n",
    "# #             for p in range(N):\n",
    "# #                 K[int(random_points[p,0]), int(random_points[p,1])] = np.exp((i-clamping_ticks)*4/l)/(np.exp((i-clamping_ticks)*4/l) + 1)\n",
    "# #                 if K[int(random_points[p,0]), int(random_points[p,1])] <= Kr:\n",
    "# #                     K[int(random_points[p,0]), int(random_points[p,1])] = Kr\n",
    "# #                 elif K[int(random_points[p,0]), int(random_points[p,1])] >= 1:\n",
    "# #                     K[int(random_points[p,0]), int(random_points[p,1])] = 1\n",
    "#             for p in range(N):\n",
    "#                 K[int(random_points[p,0]), int(random_points[p,1])] = K[int(random_points[p,0]), int(random_points[p,1])] + \\\n",
    "#                                 (Km - K[int(random_points[p,0]), int(random_points[p,1])])/2/l\n",
    "        \n",
    "        \n",
    "        \n",
    "        # Fick's second law\n",
    "        R = R + deltaT*Rt\n",
    "        Kt = D*sp.signal.convolve2d(K, filter_laplacian, mode = 'same', boundary = 'wrap') + A*(K - Kr)*(K - Ktheta)*(K - Km)*(K + 0.1) - R*K\n",
    "        K = K + deltaT*Kt\n",
    "        \n",
    "        # plot the images\n",
    "        plt.imsave('./images/{}.png'.format(i), K, cmap='gray')\n",
    "\n",
    "        # Fick's second law\n",
    "        Rt = B*((K-Kr) - C*R)\n",
    "\n",
    "        foucs_point_value[i,0] = i\n",
    "        foucs_point_value[i,1] = K[int(foucs_point[0])+5,int(foucs_point[1])+5]\n",
    "        \n",
    "    return foucs_point_value"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def compute_fick_law(boundary, N = 1, clamping_ticks = 100, top = 0.0, bottom = 0.0, left = 0.0, right = 0.0):\n",
    "    '''N: The number of the random points you want to pick'''\n",
    "    '''clamping_ticks: '''\n",
    "    \n",
    "    # remove other images in the dir\n",
    "    for image in sorted(glob.glob('./images/*.png'), key=lambda k: int(ntpath.basename(k).replace('.png', ''))):\n",
    "        os.remove(image)\n",
    "    \n",
    "    # Fixed the initial values\n",
    "    K = k_init\n",
    "    Kt = kt_init\n",
    "    R = r_init\n",
    "    Rt = rt_init\n",
    "    \n",
    "    if boundary == 'Diri':\n",
    "        point_value = Dirichlet_bd( maxIter, N, clamping_ticks, K, Kt, R, Rt, top, bottom, left, right)\n",
    "        \n",
    "    elif boundary == 'Neum':\n",
    "        point_value = Neumann_bd( maxIter, N, clamping_ticks, K, Kt, R, Rt, top, bottom, left, right)\n",
    "        \n",
    "    elif boundary == 'Perid':\n",
    "        point_value = Periodic_bd( maxIter, N, clamping_ticks, K, Kt, R, Rt)\n",
    "        \n",
    "    return point_value"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Diri, Neum, Perid\n",
    "# point_value = compute_fick_law('Diri', N = 1, clamping_ticks = 200, top = Kr, bottom = Kr, left = Kr, right = Kr)\n",
    "# point_value = compute_fick_law('Neum', N = 1, clamping_ticks = 200)\n",
    "point_value = compute_fick_law('Perid', N = 1, clamping_ticks = 200)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load all timesteps images\n",
    "images = []\n",
    "for image in sorted(glob.glob('./images/*.png'), key=lambda k: int(ntpath.basename(k).replace('.png', ''))):\n",
    "    images.append(imageio.imread(image))\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x1c5451f860>"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQsAAAD8CAYAAABgtYFHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJztfV2sNclV3ap77/fNg20EjmN7GE9kEw1RTKSAZRlLRsgRCsFWpIEHIhMpHoHF8GALkIjEYB6wxAuJ+FEQkZVBWIwjwLEEyKPIiTEWEXnAYBsZ/+AYD+DgYUYefizGEZLnO/dWHu7Z91tnnb3rp7vPOX3P7SW1uruqurq6umrV2ruqz0k5ZyxYsGBBDSeHLsCCBQuuBxayWLBgQRMWsliwYEETFrJYsGBBExayWLBgQRMWsliwYEETdkYWKaXvSCl9NqX0RErpkV3dZ8GCBftB2sU6i5TSKYA/AfAvATwJ4CMAvifn/MeT32zBggV7wa6UxWsAPJFz/rOc83MA3gPgwR3da8GCBXvA2Y7yvQ/AF+j8SQDfHCU+OTnJJycnqKmcHhW0rEzdH1JKk6Yvxe87btcY0+Z7487Pz/865/wP20u3iV2RhVf7G6VPKT0M4OH1MZ7//OdfJsr5arNzDucwjY/i9Lg3rCd+7hjSUTUsOudwDfP20TEAnJychOl4awmv5avhpbpg2DUAcHFxUUwLbLdRvcbOOVzTev1DNy8cAL70pS/932ohC9gVWTwJ4H46fxmApzhBzvlRAI8CwNnZWR7C7jlnpJQ2OjCfR8de2lJYT/yccQiiGEoWHhHYuXZ4JgIvTwA4PT3dKmdEHBpXqyNOb4gGqouLi61n9NqghRtJmPK29mdbSukqDcfZeQuJtWJXZPERAA+klF4B4C8BvAnAvy1dYJVtlcmVMbRz6vXRy1EcG2G0EHEvUejxULKISMM6fkk9RIpCrwUuyaKFKJQkpjJReMTXDmxtXRUFh3nntleSUDLhuDHYCVnknFcppbcB+ACAUwDvyjl/Okpvo4U9XCXvqwrgvcVZfl46vl9Jcdw07IoovLAWsigRgbWVFpJQpaHnBk916PVRPUXwlIUNhOyfMwLgcnth5+fnV+dMGnyNksJ1URbIOb8fwPtb00dM7hFIRATe8VDCOBZ1UWvgY4nCS1ciBz7uJQvbPOVwcnJypQoi0ohIJAqz8Fp9eVDTw/LySMLCS8Rxenq6FW/nfL1uUykjYIdk0YsW2VciDj6OSOCmEcY+iGLIvoUsWswQYJMISiQREUwU5tWFEkcJPKLbdUoIHnlExOEpC43zwPmNxezIwptCNRLgCmkxPzjvFrNj7p1/LmhRFLpvURheWKu6qCmLVtVhYVrWCCcnJ67UZ2LxCMHO1USJ1ENJWWgcK4trYYb0gBtB6eGYDIBtpRH5M2qKo1a266guhqiKUhpvtK0pjRIReGlaVQY7KzkuCi9dw+f6TCWVwefm8wDK0/qsHCycTZSLiwtXXZTiSqpjKkVhmAVZAJtTQx4ikmglCKDs0/DSlMJ64ueGqI5biUDDakqCw1oVROTEtE4fkUFreM356SmNUt0pmCAMp6enG7MZSh52LyUKbbNqgniqwwjk6MwQGxWYcVs2uzbn7JovikhhMPFwuJbxuhBGa4P20pdGUD6PiKGFDCzeIwQL486rqsA2vcbSedd5pMFxWmZ9HkOr3yKaHtX2y+ri4uLi6lzNCG6jnOb09BTn5+db+Vj61Wp1dd1YzIIsgG2fRau9pUqiNWxI5V0HwqgRRQ+RKIl4RBHFl0iDwz3C8OI8c0PjIxLheFUWnurQ+OiZS1BCsL23EtNTF2p28DXq3wA2fR5WTiOTqXwXsyELqyB9YS3OmpJiiBSClyY69q6/jvAaetT4S2pD09VIo4UMNFydkuqg5DAjj4hElCy8e3nltfspauoiUhUWp+e2t7Zl7VJNiMhPwfe1czvWPjUGsyALVRW1dDW0+C8sv1LYdSOG3kah6VvIJJLokXLo2XTkB3zzgcM4HZshTBIAcHZ2thFn9/BIi++vKzpb6tlTFB5ZsKpgNaHnBiYU25+fn2+U0VMWnGYMZkEWAK7kUklFlF4SV5KdR4ShJDAlKRyKYGoNuIdIItLwpLh2Mg3rIQnt5N7o75GCKhAmhBbVEREG58MYqiyYHDyTxDq1Ekbk4Dw/P8fZ2RlWq5UbZ8qC/YFjMBuy0IZolcoVX2rwShJ83EsY11ldtKKkGFrDNH4MSXidF/CVheebODk5wenp6cZ1aoJ4zlO7rlVZePXmtQ81HyJlYW2czQb1WSj5cBvWdRZ2DwszZdEzWESYBVnYC/Qqh8/5gT3zQeNaCMJL75WvR4nsm2B6VUXtnMN6yUNVhhfXoi4AhOaG58y0cCUJu8ZTI9G5Rxp6HNW5+iIszCMLVhCmLKztqKqw+0eDp4arsjg/P3fL24NZkAVwaVNeXFxgtVpt2W/c+TzG9hQIpzOw6rBzTs/HNYKYC2H0EkUp3lMKXholiMhBGHVEb8QHtsnBMyE0nPPQcE916P1LZOGRhldfBn7f2iaVIPTYyEKnT60fWJyZHhy/Wq2upk8tTNd0HJ2yALDlkImYtDQ6ep2aVUQpTK8ZQwi7Jox9EoV2HlUQLerBUxCtRMHx3PGZMJQ8mARYXWgZ1AxRM6VEnh6iAQ3ARqeOzBAjDmuDSsbq6ORrvUHTFm6NxSzIAtj8PYuIMPgllSQZIyICD4fu/PvAFCNMD0mURuyIKLiDR3uPFM7OzkL/xK1bt7ZUBO+VMCLVVKvHyE9hYezgNNOACUJ9FUYalk43a5O8etPKZ6bIUZkh/FKskjzCKBFEqROoirD0Qwikdl67fl/w6iOqr0ilReqCj3u2Fr8Ed3yO5z2TRUldcJxHMJZ/tPJTNwb7LxiRr8KOeUWlmQrsr2AisHLyYOnVv7VvLiMvAzfTZSxmQRbAtvfZIwy1vfQ4ck4qMWheJSfoTYI2Rm/P6SJCaiENHs09lRBtrAQ8gvHUh6csPFWhKgjAFpGUSFZ9Fp5PzZSF+iZOTk42lIXFcfu3YyMVJQt9D+zzO7oVnLwHNgmDG5pnlnjngL/+ws5bFcFY9bALAiopqVr6VkWmaTyC8EbfkqpQH0APUahKUN+FmiGaNjJ3bMEWsK0qPJ+FQdWFqgrbK2nox1/aqU0NsLLw/BNmWrA/wsrN5DQVZkEWLEu9F8AKg5evlohD849Ig+NvsrrwCKREMEM2Hrkjv0PU2fncMy0sPMr31q1bG8SiCkdNJCUPr14iZeGZH9ZOlSjUwWk+CvZZWHvkmRF+D6X3NmU7ng1Z2HQQEDtybKUaEHucdcrK8jd4ki0iCr6XHutL2Je66FEC3nmJCDx1wMetROApCW9U95SAdm7PH3F2drZxrZGEXeORjJKFR0j6HKw+tE5q0NWZ3NmVJLwpUJ0ytThvKtWIZbVa4ezs7IqQ+DP1ozZDWEGokvAUhdexDdEIUFISQzr/HBRJ7dlr6TWsVVGU4rkj2vvl63rNDvZNeH4KVRh2HpHTULKI2hWwSRbeIOiRhT0DrzcyX4bdj01yGzx5hab2iylIwjALsrBKArZnQ3hvL4MJRU2T0gvU+GhaVTu9l2boc44hk9ZRrXZ9jSD0vDSi1tSGEgSw6Rfg914yNZQINA2TgcallHD79u2QlOy6iNwiteXBmy5Vf4WnNJRQjCT0Q7GU7vorjChu3bq1QSiMoyMLwJ8F0b01LmZSb3QzeFNKNRLoIYY5KAnGWFXhjZw1ZRGFa5poCpQVQWR2eMTgEUZEFkwInk9Dw1qUxVCyMOXgKQwmAp065bbM/g2+l5aL00yBWZCFvRggVhZccZrGKtHzLVj+HnrUwlTEMDSfsaqiJ++IRFpIgUdmPfekfjR9WiIKjxQszMySyHQpmTlROfV5Su8jcmyqsvDIAsCGqmACYbIw80NNDY8ceKZkLGZBFsDdh7KXp0RhcQw1Q5RVWzC1utgF8bQ8S0uH7z1vUSalLZL93FlVWbQoiogUbt265c6KnJ2duQ5TNUM8VcFhXr1471Id7kwcJycnGwThmSHW5m2/Wq02yEJJweJswRfj6MwQezk8ZaREETk0vZcYKQwPESnMzcTYFcaoCk5fMk80HXD3mwtWFiVzhDcmilKcqpEof1tnwWVRM6T0LEB5QZZHFqowPGXBJodXn2yWGJQweLn4WMyGLGyFXWnxCbD9c2E8rQr4/9FQUwaWno8tzptWrRFLq7rge9XSjUmjI2GJcDlNi3LQTlXzB3A4d2BVDaoSOO6ee+5xFYaFe9eZg9MzPbzVnUD8WxdDyMJTEmyaXFxcXC0f4GlRO7bB0zabSmW/Bavx1Wp1pdY53RjMgiyATfPDYBXNL1LXQHimB3dqDo86ZkQSUdpounUIaveaIv+W/CKZ7aXrVRt2rNI+UhER4ahqsPUWLapDiUnz9ZQOn5fUEtBGFqYoVFWYmgawYYJYOdh/4RG/kYPCpl+tDGMxC7KwSsk5X7GhepAtnaXlSrZGxaaKKobIHIlUQU1NtObT+vylKd7atS1hUZqITGpkEKkQdQp6xx4peOclH4USBZ97hGHqwTNHSl+rGlkYLFzrS8lCHZvWnu2eqjJskDQfhZXFM0f0HdhHYt7P/h/dbAhw18Hp/YgHf0BjEhHYVAJe47Vz3s91KrREaFPlPyRO05VMFc7L9p6jMHIkcue2jhqpBI8YouOSyWNmSOSnUDPEIwoGT516ysI6vpog1q5ZTfA9dXqU78+EwOYHcLdfHc0KTlYWrBZUZbCzU1WGHnPnYxLQl+yNDpGa6FUPveTTSxS19DXV4R1HiiJKy3sdlXkU9gjDOqySBpsGJZPC21vnb5k98ZQFE4Tn4BxKFkYO1ibUBOFvQvRn/lVZ6LsoqYfT08vFjEf1ibo9lEcS7KxsnSmJZLPlG8Hr4K2d/pDKhMswNE0L+ZSIJCIPe5+RI1RJQzu65xT1iKNVWegCLf6mREnC81kYWOXqyO2pCh6A2PxgZcFrKGyzNAZdrWl9x+7LZeD+NBazIIuULmdDgE3VoGmsIlRqnZ2dXTmPjEg8pcF5ecdA20yIl670bIcgkJIS0DAN79mskavpwB08pbTVeb3Oz8rA1IGaHGdnZ7h9+7Ybfs8994Tmye3bt7fMECuX982IPdOUysKWZnszI9ameSZktVpdxZkvw8LtLwCsH5yenl6l40WMOeer8LGYDVkw07PTUiWcpbMK4U6rjk8L55fsdXivPFPOeOwLY8yY0rUthMF7JSLtdMC2L8NzdKr/wjNLPP+F7iMzRONU8egzjSUL69A8M8LmiNULmyF37ty5+pI0upct2uK2agrcVMfRzIYA2BiV7MHUf6E/YsrHlgfvNUwJwFMwQNtvAGg+cyWVHgIpmRiaX8n8iAjEIwTtuF6njUjAUygtzs+ILDwzRMnCEI3UNkjxQMcbsPk3hJ7K0PvZdGv0LrnPWP7ch1gtjcEsyMLkoNlr7Oxk/wU3OiUNVRW2bzElojLVHJ49zzc1mQxVESWTrHStpxj4mDuXhUeKQyW/dtLSCk0lB3ZmtqiMVmURlb9Ud6YegM1/RLc4XaVZIwu+T6Qs+J58HyOQ2SiLlNLnAXwZwDmAVc751SmlFwL4bwBeDuDzAP5NzvlLtbx4lOHZD5NrTAzc4SN56MV7yiJ4rsGq4VAqY6wJ4l0fKYuagvOku/owgO3P0j214ZkdLc7N0pqLVp8FlzciWw/cdtWMZoLgY/7hGq++dbEiw2Y6PAUzG7JY41/knP+azh8B8KGc80+llB5Zn/9oKQOVndzAeFGWpbU4M004D09d8LWqLIaaIZzfXDGEQKIRNAqLJLuliVSF+jE8teFtUcePSMM7ZnVhyiLyVXiqwiNLzwzgjQc+b0GWN7hx/mpq8H2YEHT2cO5myIMAXr8+fgzA/0KFLIDtUcarZDZJuKOb48h7iXbeoyz4Os8UsTJEaWt5jUVLub3zKLw3r0iO897eVY+qKKkLTddLGt6iLPWTqN9iSrLg0d0jBv6YzFMvvE7Cc5wqUZyenm6YI3NYlJUB/FZKKQP4LznnRwG8JOf8NADknJ9OKb3YuzCl9DCAhwHgBS94wdVy3GiBiT2wMqepC90zk3qLtqgc7v2ukykyhQnSsmnH4ZHY63DqB1AfAXdk3W7durWx53DbNFynTvnY0kf+jxJJ2LFXdwyPMNhnwcrC4lRd8Idi3KYBbBAd1yf7QHSbgiQMY8nidTnnp9aE8MGU0v9pvXBNLI8CwEtf+tJckrQazl+T6gvWBStccUM6f+26MaQyBK2qokSCLfloGq1PLw+vg3n7ktznPFhpsJ+hRVlEfg2vw3nOzR5VwWB/gSpRVg7erIcuLuSPw/iYVYqqB1bhPElw8K9Oc85PrffPpJR+E8BrAHwxpXTvWlXcC+CZWj7WALjzayXoiMbqgjsrX+/dxzq07qNyDYmbE0qNu0QYEVF48HwVHjF4ZMBh3miv5gfPXHgqwdtYzShZ6OYNVFEdeMqCySEijmjUV7IA7pKEmuBsbrACubi4uDJZ2G8xxXLvwV6PlNLzUkovsGMA3w7gUwAeB/DQOtlDAN7XVJBg5IlGnRap6I0MHO+NlpomOq+hN/1QjL2P95wtdaXvxOK4c+kWvVPvfUf+DHVKthKFKhTvGlUzJTXkyf4ojT5PTR1x2ugbGX52rY+IFMdiDN28BMBvrhvKGYBfzTn/z5TSRwC8N6X0FgB/AeC7axlZhUZszzYZqwZPHnqmCN+npCg4/S5Uw1g1sisC8giV79eiTpgkvNG41omUHEpbSTV4aUodsqQqSgNNrT65jXntzWsHph4U7MTk51MHJ68G9baxGEwWOec/A/DPnfC/AfBtvfl57OyNSpaWzQ+Ns3i20/hl8cts6cDetVH8nFHq9LW0JeJQ0mbCiN6NNyJH6iIikpoJ4imPUn4eSQwhC68uuf1w2+V6MZ8Dw8iAySGa+bBn0rCDk8XUsIdTE8JrQPqJOpOGEoQni69rZ69hKuUR5aOkwEQAYKuuvXfH8fp+SqTgEYTGeyaFJ+15Xxqg9JnH1qGF6TQqf/hlYNVg6Zgw7Bp9togojoYsuCJsD2wvOuEHZ2Zl0uAGGP3OhTYGdj5xmQyeIqk5Rzl9Lay1jlrCvHhvZKyNoqU4NTOUFPh6bbBsW7f8gE3r1KpOpepai+jvEK3teQpoCkRqzmtH2p65Xm/fvr3103l27jkvLT/7InuKwXAWZAFsOziB7RHKwrTzc1o7tjyjD3BUGh6DspgKJeLU88jM0PPIFLE81G/QY5bU/Bw1HwXfx3veCGpK9NattWPec77RQOmZYpGvQvvUGMyCLDzZyhXAFeH90I02PDVFvPv1jvhj4/eJMaOip0K8PHtMEAv3TEsb3fW6WsMvdRwmh5Y8W4nC8yfwcam+OJzTavu2vIwgWC2zGVIjB+95x2IWZAHEzq/aZtdqPhyuakTXZUzZyQ9FGjWC6CWQqG6VJPRdeNe3bEo2nqJgYolUiIbbNZ5ibSUKbxbDe8c1c5bDdN0F740kLI2Rn84KRs+tg633foZgNmThjUQtJMFqg1WFKgxFj7oYokRK6L12KqUwNL/WPLgz6nnLu9SO7Q0ealYoqZTUB3e8HqJQgogIQ0mgVn88aHF98Q87scpQFeIRoxIEpxmLWZCFPpg+NKexY92zZON8o/upvakvLkpby6s3flfoVRql+i3lp+9HCcPLjzsHk4SSQ6QWIuKITIyIeJQwFDpNyWF2XKqLGmlE70DNj9LaiYg0vecfi1mQBbC92lJHE1YP7MuIXrY3utm+pfNa+mOaZo0aZ09D4rSl0Soa2Vo6rhcXqYiITOxeLWatByWHaNN6UdLk+Jrq8J7ZVIXFqzmlz+7Vfe87jjArsiixIpOGRyJcIRxu59GPivBLvK4kMCU8FaANTUcqr5HyeYkEdDTkrbYCM1r6zMukOSwybxRMBvyVqLUhdjJyvfHGJlXtflzXSjI8BWrnJycnG99/8PoKmyq11ZwWd/APyaaCx/SeCeKl072xMedteyabVt/Egkt4ctZr5PpualsLsdTiW+6t5Y2gqyN5szDel+qI1ZR+KR3VsWdKR89leZeIOFI6QzALsgCGv1wglsOeohhbxqnI5FDE1CNHW9KW6r6Wdyup1IjDI4/S4KJhBs/M8H4jU39fk/NVXxvXB5suGs9lYn+FHetg553X6mgsZkMWBo85vXUTUzw859k7OzFXZeLVS29dja3bSH14ebc0Zo8UNLxnqyEiC/1JPG6XPMKzn8HK2bqKcswzReQ6FWZDFvpwWgnRSMXXMNNGaVrKcegOfyi01lNUv55Tz0PkeKt17lobiTpPC7xpUSYM/vMfjyzYzFCnpK0i5rUTJXNEUVMLUR1652MwG7IwtHZ0lWG919fixmDfhLOr56jl3VLvPe9Tz0uk4RFOiSRaSESJwvthXftDH/v5O4O1RXY22n30I8mIKJhEWoiRjz3Tp1dR1TAbsoheMO912nTf5bupimMsdPoUiDtsLTxq/KWOovctQRWFkgYTBf/pFX/bweRgaqPFwWnl43ZWuqZHlU2B2ZGFsmQ0QvF1pTy9c7tP1Pl3OVLPFbV6LHXiKKzlnfH7Lnn0S+Xi42gBUq3zRGso9L9HeX/nzp0NJXB2dnZFKKYszATR36Oo1U9EiN76oqHvrhfj51OuOVoqvHbtXDHmmWrXTjEV14uIQEr+rJY68AYNXanpTZ2WplV5atVbxBXdN3qOCNHU8y7a5tGQxSEa74JttEyZDsGciNkjj2jtxSExdZ0dTQ/btw9jgY/aexjamebWCT3TaapvMKbC1HV2NGQxFNFHQT3XzhVjnql27SHIOfouIypLz9oGL4xNMiUH7/OE6BuVVqds6bkj2LOXvluZCrNxcI6163rvtaAPu64zlfS99+P00bWtS62jWQZdth39QlXJf9DiE4oIe0zfmOL9zYYsGNGvHHsoVUIpj57rbgJKHam0kChyDkZxGq7LpjVdiTg0zvLyflK/ZREUz5BZWo8gvHyNTPjjtt6ZPc3XAxNqLa3FTdWeZ0MWXsPQMM/rbOn4Gs23dM+e8u0i7RRo6Qi7yLul3j1noIV713ty2nMk8l/0aTuxNLzGwWtfnllgHRzA1RecFsdfsOoKTjVN+MeCbVq1tgy79NzR7Ipt/HWsV2dHpSwiorDj0uhj8MjDO6+V46ai9dmVQJjMeTakZXSMBodSvHcekUbrM7FyYrXA5MFptcOrv4LVhf7OqOXRihJJcLym8c7HYDZkYfDUhMbzfqp7TpF+DkTjKYFe5TFWqUSKQePsvNaYvdFSw3u2FnPE9hFRGCloeOT87FnzMWTjaz1TZQrMhixKjFl76Eh1TO2xPyRBTXnfHi98LW2p7qO/47O8Wzu9Nv5aGlUuJbXhOR/53CME9mN46WuE4dWpV86aWdHynEenLLyHahlBLI736hyNKj8qxwIfVn9qfphs168uLb6FCLy0NYJQcmD/hZdey1uaFTHwDymllDb8JPp9kpJBzy9lcX3yXp/DO/cWh5XqeQxmQRbA3cYYOcFKRBHlZ/AqqiWPmwzuUNq5vFG1lleN6FsdeV7aKA2TF2886wHUv2/hj7msjUZTosD2b3+2OjZbnrlWHyW1MRazIQv7gk8f2rzOWjEW5jnLuMGp01NZ2zu2c6+SrzO5aKf3OoyXpqYmtGNyuP2svXZwfn+c3uL4/z9V2nOcfajFx7bZbAUf8/PydxUaZ4gUi9YfXxets6iZH/o5vP6GBn/typuXnj9408/ph2IWZFEjAk5jx7rvmX/2iEHLUstjqvhdQTt5LV6fWeV6rS69Ed1Lw3t+39YpPd9DNHIqGemmU6emCLzylNSFHasqqZGFHnt1V1MR0YdrunmEo3mMxSzIAthebBJVpKZhEyM69tCjGqZWGL3X1jp+77W9+bXmoUQR+SOid2mkocpSOxKn006iHcY6qx3zn1B5BKkoEUbNfNHrtf74mJ+tRA4eIbAaiYhiCp/FbL4N8Vi1ZbNrNR8OLymPqUf+QyqJMfFReq0vz4Hs1WMr6fN799KrNC+NuF64dr6oTKU68pSCN8PRMutRqxsAWx2eF1y1PDdf47X7oZiFsvCIQk0RfnC7xvaej4Lz9u7XEjZl/D4xlRLRY4aN1Hzu+S2896TvWf/wNxpZSyOuEQN3VCujmiEMnfEo1V2pTnumo21fel5PLamfolRHXp5jMQuyALD1gMB2hVqYqg9Oa8eWJ58zdHRccBfaYUrn9i6soylR8Lld6xF8SWZHJGH+ipLfwkwPawtqVlia1h8b7k3j1a3ta51bTYsovoU89mKGpJTelVJ6JqX0KQp7YUrpgymlz633X7MOTymln08pPZFS+kRK6VWtBYkaR0kyRmTgKRA7Lp1zWBQ3FkPzmIrQvOf3njeqPy8vryF6ClCJvkQINRXRstU6UTQDF73/MeA8PRVQI4Wxz7wvM+SXAfwCgHdT2CMAPpRz/qmU0iPr8x8F8AYAD6y3bwbwzvW+iJwzVqvV1lRPND2kfynHMs1rnCXS4bQ1BaJhtRewa8XSqgA8aR3JbU+Sa73ZaGx7u47NAO4APIVp19l79Gx9AFitVlfXGOzY/rpPy6wflnGH5A/AeG9pS7//6ZWh5d3ovkSWFxcXV+2f64d/79P6CG8WfufOnY1jL+1YVJVFzvl3AfytBD8I4LH18WMAvpPC350v8WEAX51SurelILXOrWpBZaz3UjR/PW4dPUqjaxQ2R7SSnJc22ttxVP/Ru/FG2hazI5LgpQEmkutemDeojN1Kz+uVozRIlp5Rn1fvMxZDfRYvyTk/vW4AT6eUXrwOvw/AFyjdk+uwp0uZsTKIJJTu7TptvFwp2im89F78rjr/2HwjNTAWUblUmUQKBdh0cPKxpbEwUyMlQoj8D5rGVMlqtbpSLJyGZzB0URI/s3akkrpoeQdR++R2bvf1SKNEDJpW1Xip3sZiagenV4tuS0wpPQzgYQD4qq/6qg0y8EaZiKEj5QGUP4XWvUcspfMa9qU2xhKIXm/l1s7hNXqeaWACiJZZK6G0qAp1YlqHWK1WG6s22azRac4WuwfFAAAgAElEQVTVaoWzs7PQ7LKy8314OpTriOsiIk6uJz7OeXvlak1hrFYrl0h0cNV0Xn5jMZQsvphSunetKu4F8Mw6/EkA91O6lwF4yssg5/wogEcB4KUvfWn25pI5zJNxkdLgMLmnSxQRhsbNCTVVUFoL4BGGl57Vgs6KWJynOlhF8Ll1El0ybiMyn3OnjsgCuOsD8QYNbmNMDnzMedo+amO21/YWKWbP9DA/Q0lxcBz3F/Z9MJmMxVCyeBzAQwB+ar1/H4W/LaX0Hlw6Nv8ur82VGiISiMIj5cH52T5SFq3lKl03NN+hqCmJWseuxUf30vr08lATw5QHfz/idRaN98iEnaE6+itZGPjYOgt/Ns8KgdWQEoX303h6zPlxPVmbZGXhqSolAiYK/XMjz/zgze7DhDJF26ySRUrp1wC8HsCLUkpPAvgJXJLEe1NKbwHwFwC+e538/QDeCOAJAH8P4HtbCmEVF0ksrVwdESwP23ukwnHRsUcsnPexQDu5R3hKEpFy4HPPFLE87N2xwojMDFUU1mG5XahPIqWEO3fuXD2HdurT09ON0TXnfDUbos8WqQtD7YtbJQjelCw8x2ZJUdhsByuIkunB52NRJYuc8/cEUd/mpM0A3jqkIDbdc35+OU3kMShXIle4hSuZaGP1pGFEBJEaia7pUSBj4XX2lvOIJCIFYlCC4A6kI6ilNz8Bh5sj0qBOR76vdWQL43fI34/Y72SqIjk/P8fZ2dlVu7Hrzs/Pr37yzvZ2rX7hGk3reirDUxUlslBVcXFxsdHutf1/5Stf2VAYvD333HPFKdWppk5nsYJTJZk2DmZGjmN55akJ9V0MURW1cs8ZSg4t6ZUkLFzz4/fAJoeO0hwGbC4TV6WhMxrepo5M3rhD6HMbYdhmhGOOT/1lK1M1/PuZlm+L+aZb7blq64yYPFrSefcYi1mQBYCNxuA9JFd41NntXPdevB4zWtJEOBSBDCEGz9TQNBquJofnBFXzhMOAuzKeicYIgNuBRwgAtkwQg8548DtXX4WFsS+DFYaaVFqGUn3Zs7HiUrPDUxlqfqiq9s498ojIZCxmQRY556uKUhbmClXFUVITkTli+1qnjhTJ0OebGkPJQQmgJR9VFkoYqgxqvyOhKoJ9EtZprON7m9fwVVlo+XW1Jh9bGjVJ7NjK2OK3KPkrrJ0rcURk4Tk3I2XBYZG/YyxmQRYArphWZZTHyB5hWB6898I0bqgZMkZ97BM9pBI9f+T/UIXhqQsmClYXTBZ27ikM60R8HS8FZ/LgZeD8vi8utpd2e2Sh5ojns2ATKqpDi4/IolVZmP/O4mqk4c2WHJ0Z4hFCyRQxX4WqB1UVERHUSGGskjgUxpoi0bU1hVbyWVinsI7GnUn9HfaedbqUSUR9FgxPWeSccXZ2tjGiG3GYE9TKw0QBYOO/PyInZ1RXXhuMfAnqs7BnN/LwiCIyObzvS/Y2G7IvMAF4ZgeTBLDtIa8pBz6P4j2UlEkL5kQ4ETkoUWiZPb8Ek4HGMQFwuHZyTquEYGpC/RMaxuUzv4S+bzYxVGEY1MHJz6k/uBsRhnffyGehW2RKeEqiRB6Rs/OoZkNs2qi0co1ZuMTUTDAeqXgE4RELx2na2vGhoB2awzSNd+wpDFVv3iIqb+Ti5dl2ztOYPNprBzNTwaZfbc+EYlOj5+fnV0Rx584d3Lp1C6vVauN/R40odNr09PT06q8FjRT4Z/ynVBbcdj2yMJNDlcLFxQWee+451yw5Pz+/mjr1zBNekjAWsyELr7N7vgt1IEW+ixY1EZFHpCJq595z9dbDUDOiNb6mLlrK7CkMNkM8ZcG+CX6/1gFZYZhfgv0R/ByRuWFkYOe2KTF5//2hi7+sHPrXg0PJQv0UXhuPTIiaj6Lm/DwqZQFgSz1w5XJlazzgk4OnFKbu3FNCCaqHNFrzLzXwlvtx2bz6jgjEoGaKmiEMUw18rUEbvpXFCIHD2eQwNcLqxcglcmxq3URLv/W+an5YOVrMEFbOETF44aWVn1O07VmQBZNARBieSRGZHJFJUVMUmtbbl56hdN5zrYW1duCSeqhdV1IapWuAzZWdPKKqs1OPDbrgyfJkRWFhTBDRSMk+CH3HphKUKKzdeCaIKh4unx3z/aJ7q7Lw/HMXF3dXcHKH9/wZHmGUpl1tG4tZkAWADaJQ80MJAahPhep5qfNqfC1tLU0PavcaqzAiMojS1ZSNV09cf54D1CMK7ogcpnFMHEYQem738UwQJiszQXhvBMMmiIWpSVJycLaQReRbs/bOPgsljNo0qadO9H5jMQuyyDlvSC8lDGPiUqWXNrtHq6rQ815V0fPcLWl2ZZJ46qJUPiUSjyy0c0YEwU5JjjcVEf2pssUxOXAcj+J2bOnMBGGHJxODOjU9nwY/f40sPP+a+ix4007PSqRmZnhOUQ0fi1mQBYCriuE9NzaVbyrtaiRRuu/UqqKVBKZCjVB2qS5YQXAeRgZMHhZmvgMOZ/+GjfoRvIbPZTdyYGemPYuaIPxRnPdnxlN+SMYDn2eKmLJoMUMi0ijNJI7FLMgi57wxJcSsymrDUxRa6fxiaoojIhU9rpW9dN5yTUv6Vh9E73lp31Ju9VEYmAR49LY9mxDcka0T88hujZ9nO6zDmzowh+jZ2Zk7RWqKhB2Zmib64rSmKrw68tQFk4WaIh5ZqInChGBt374sjWZROG4sZkMWEQlElWtxdn2p05dIwtK1lHGqZx163dTmSJS3d85773pgc7ZDP8ZiQmGHp5GANWa7b0k2q5li9+d1HGZu2J59Ft63HzobosrCytZKFmqCeGaHkgWrAiWNlm8/Iueo7cdiFmQBbDs4PcJQgrB9i2LwzBGPMHpIZEpTYgrUOn1LesD/4Rs+Z7AKiaS4miC2N+lvxwxvJGRS4DAjjtVqtfXdh91fHZi8MaFEZog961iyKPnf1CnJprdHFpG5wcecdixmQRZqagAICQPY/EOimrnh3SsyM6IRtFeFlJ5zDMaqCzU1SnnzuZKGXhfVJZsm3HHZn8H+C0b0VSfHXVxcbHzbYR2cycLUxHPPPbexUjMii1Zl4dUHP7+aHxbmKQo7ZrLQr1OVLDjNnTt3NvpQpDDGYhZkAWz/LDqw3XlbCMLS95gfEamUMOSafWCsutAwjzSiPKK6ZfMDwEaH5nCGhnF+TEC8fFwd3Wxu8LXqu2CzqYcsSqiZIRFZqNmtZkjk02Dy8EjjqJSFZ8959p36KVoIw7tfFLYrVTEVdqEutLNHRBF5/0vnqia40ZYaMI+EOmKzw1Lvk1LamCLlVZpGEkwaXBe9zs1IWZQGLp46VcJgZcGkwbMc3u9hRIrjKH0WXFFqt1kF895TGCW14bF8dO7tS8elsJ74HvSqhSj9EAXh5WfKQX0Ddm4d9/z8/KqzssrwNsvXzo1seEaFPyHnqVDzXeiMh34gZtfyL2zptyDR8m4rA4PJLyIL7rgRWWhbv7jYnClRn4aSBcdzfxqLWZAFsG2GREoiIgVO46kK77yEoWmHxA9Br8IoEUMtn0h5WBw/X2SOWEdmyc+koSiF6/15Fab3P6Y8Net9A2JplOSA7SXprT6LaMDSgY733OE5zPpBtJTAMzuYaI6eLJgoOKy0eSiRRklV3EQwGTBBeGaLodRh9JyJwxAtwLJ37v2ZMStFJgYLU1OD1YiqF1605RGEkochcsCqutC2pkvUua17fYBNC89EUZLxZhSnMEGAmZCFsmBEFC2EUVMdfE/v2Ctbz3nt+n3BUwyt6iIiCiUTTl8qh8FGdj72fBnA3bUU3pSrznh4z+ORHCsLVhKWt0cQrFpKqkKfl9ubkpvFcRtXwvCUhae+lSxKamQsZkEWAKpE4cm6HqIw9BBDb/x1gEciY/IoqTpgkyAsvcZ5hMFOUSUKHlnZ3AC21QOXlf0aTAa8alT9FD2rN/kZlTTsXGdsdKD0fBamSCLlEf2gDpPSFG13FmShkgpoNz0iguC8ba/kUSOUseSxa3Kpdfye+JKaUCURzZxE9yipCwAhYXhEAeBqJNWvTC1P/RRd1YcqCyUQ3ayM/Lx6rHVg99d64Dat7ZtnStQM8fwRGu+Rz9EqC37hWpE9PosWZVFq3F58b8fflwoZqxQ8wtC4KNyO9Vzz1zDrnNyhlDBUWSjB6LoKzgfAxsIvi7P7qopQZWH5MGHwKtMhysKOo3btmdscZn/NWFIhSgwaPxazIAuTWfxgLUxcIo1ISejeIxKvfKXzWvpdo0QYqgC89N65XeMpC+8els4bbbkj8nQq34MJxCsbd14AWx2eZ1qAu2Ri8Tad6pkZOl1q4byPiCKqm6gtem24RBR2DTsqI1URmSBH5eAE/KmmFmIYQxR6/9rxMaFFRdTCgDqhlOrPOrw6PTmOVQVw1+/AJGLp9blYqbC/gu/LROUpCW/Wo7QUncvCdeBtXGYeKL0B0jNRjARUYUVxYzEbslDWHEISwHiiGItDkUtJXbTE19J65MDhFqdhNUUCbJKGdSB1SlqYtRMlCiUR7xlspoNVA9+3RVnY/UooKdhS+1ZVwOrCW/nJJFLzhUyBWZBFS0UOURfe3rtvFHbdVEUPIXjpI5KIFAjgEwTHt5CFQUmDV30atBMbUSiJKGmwGWL3YdPHrrVNfSCMocrC4vTc9rUO7ykHJhXNC9gchMdiFmQBxLZcLznwcYkMStdF15bCrwNayMALL5FQRBrWUaN3oCM4gxs2d1zt6BbGhKIrMpk09L5MGOpv4Xh+Pg4rIRq41PTwFICqCiYL7ReenyLqS2MxG7LoJYjoZbSGDS3jmPh9oKYuetSHEgbgr+704vU4yl+dnzzNyTMfOtJ7pKBlNscl58FLzTnOysvKwnPY1lSFocVvoc+pPgmOi1Z+quqwfFRZTNE2Z0EWQ0hCyaBGCiWVEcVHaYfE7xM9hKDp9dronPdAH0lwXgrrpMDd0d2Ug6cstPPybId1FI9gOH9P6ejsiIWp0oieyWtbSgAeaUTqIToGNglB03DeYzELsgDqDk6OA/ypqYgsai+QUSMVD3MiihZ4ZBJ1ek6vBNGbP3dAvodull5nLTxCYxMiCgc2VYbGRV+Wlr44rZEFH3P71XDuxGp6RGpBScEjCzs/utmQHhXB6fVaTePt9bilbGPiD4WaumhRH14HjTptjVC8MO9cyQTY/gRelQawqTA8BVGK07rwvjeJwDM4ESJy8AZCT3kMjfP60BhUySKl9C4A/xrAMznnf7YOeweA7wfwV+tkb885v38d92MA3gLgHMAP5pw/0FKQ0gPzuR1rmIZr3t6xV4YFdWgn8oiipkQ03MtT49ivwSThOUFrcRbvEYUqK45TREThhWsb5+eMiMDOa+qB4w5GFgB+GcAvAHi3hP9czvmnOSCl9EoAbwLwDQC+FsBvp5S+PudcneiNyKGmKqIwPa6de5V5XVWFYYi6iIggOu4hilZfBrCtJrSMll+kLDRO1YQplFKY5sv+lBpUGWl4pDIiU6WkLPhcN/VfjEGVLHLOv5tSenljfg8CeE/O+SsA/jyl9ASA1wD4vYb7FJVEiTw4Dz0eqjauO1EY9kEYmk6JYojzs1ZW9XOUlAWvxvQIIfqFcVYZlndPWfW4RVFEYZGysHwjJbFvZRHhbSmlNwP4KIAfyTl/CcB9AD5MaZ5ch20hpfQwgIeBTbuvRhaaRtNG8YwaUdw0jCUMIHaGalhrefSzdlYYqjqUMFRZANs/ssOzHrzCUcN7CKL2TLaPplVrSqNFeXBajzzGYChZvBPATwLI6/3PAPg+AF7NuqXMOT8K4FEAODs7yyUCGEIEUeW0VNqxqApDTV20pimlL5FDqx9D84/uExGJ9+m73Z+nUT0ThNO1+CpaZkO4LN5zaXwraXB4iRi868ZgEFnknL9oxymlXwTw39enTwK4n5K+DMBTLXmWHJsczmEaH8XpcW9YT/xcoQogSuM5GUsqovXeHmFwnIYbWFFoOg3XhVTsY/CIwfNNcFiNPNWf0VIPvNdreJpUr6k5/mvhU2AQWaSU7s05P70+/S4An1ofPw7gV1NKP4tLB+cDAP6gNd8SUXAaPV6Ioh29CsK7JsrDIwGNi1QHo9Vk0fuxI1PJR1WHIfJjcDyXi1H6QGuswrC91/k1LlrktXeySCn9GoDXA3hRSulJAD8B4PUppW/EpYnxeQA/sC7gp1NK7wXwxwBWAN6aO2dCIgJoVRZ6HN2vJewmooUI9Lx0bKgpCyWRUvm8D8VYWXD+ngPUCwd8Yphi5aPXNr18PVPD9hFR9MSNRZpDJzk9Pc3Pe97zmlRFq8rw0njnUdiQNNcFLcrCS1ML845LYaV9dMxhJROjJX20Z0TfgrT6LIDy/7bycdT+S+Etaez82Wef/VjO+dVuwRswixWcQJt5sRDFNGgxRbw0tbCSYtAwQ8n5WVMaWh5VFp76sHxLez7elbLoIYtaXE2R6P2GYnZkEamHVsJoPe8t17FhSsKI4kudv9dJ2lrOiChKcS3l9ODVTVS+luMhRNEbNwazIQtgOr9EK1HUKvBYicIwFWGUzkvHhlLHnQKsLDxiKKmZqH6iMg4ljDH7G0cWkZSKwqI0ip4X2pvmGNAygnppSgqhx5TQ+7SO/JHpUAv34qKw6LwHrW22pd337u14qrY8C7JoNTW8a2p59l7XmubY0KIyWq+L1ERP+khhtJJPr9qpqR4PHnmWyhOd19q9dvhWkojihmIWZAHUH6700FOZHa1pjhU1wojiWwgA8JeEW7ge10yDnvcUEULkcygpKC996b61sF514YXtmiQMsyOL1uMobAxRLKhjKGHoeemY0WrSDPFzeCTCZagphx5l4aXpafO9qiO6dgxmQxaGFhMkihtLFAuhtGFXhAG0qQ/DGILwwvn+pbR6Tcs9a9f1KOpeMukpawmzIYsSSfRU9th733SUZLemG9rpIvPDiyuFtyiOGkpkxvfVZ2rJtzW8pe2PSXN0ymIoISxqYjdoIY0Wez5yGmqa0nv0OnOJSLw8e2c7xsyAlDCULHridkEUwIzIwjCWJGpxQ9LdZNSkeC1NZJoA8desXh6avuT07InbN1ra8hSEEd1rDGZFFgtRzBO7IIwoHRCvoKyZL1Ou+dgVhrbxMSSxmCEdcWPSLrhEK2EAbWZJKX1kYpTiFa2mzlSO01YM9WUMOS/dbwhmQRa7srFq91qwG4w1S7x0kVrg63oURMtU7a4wJWEMTTMEsyALwxBn5aIo9odWc8LSAu0qwwtvUR6l60p5eSThhXtoVVg1jCENL+xGKAtg2IMu/on9o+aM7M2rZb1GS+dvDWs1RTRey9eKlrRTk0brfXsxG7IwTEUQQ9MvaEMraXgdvSW+ZeTvMRlqplFUPq8sU2BoO58qfAhmRRYLUVw/TGWa1OJLqmFIeGu8l64VUyqQQ5KEYTZkMSVRLCSxX/QQRkv6XgephQNlk6GmbKI0XroxGGOajIkbi9mQRYRFTVwP7IIwgD7TpSeuJ02EVlXSgjGD5dSKPMJsyWIhieuHXudnS/opTJeW6xm95e/BvpyjQ+5Xw+zIYun01x9DVAbQRhpRul4nZeu99ompZvd2Vf5ZkcWumXrB/jBkirVnhmUsKfSQx66wb5UxFrMhi8WBeZwYQxql61rNh9YZjjm0qV31gaP3WUSYw0td0I9e04SvA9rIpoVkovSt10yFfanoo/ZZlLAQxfXGEJUx9Npe4tBrIuzS+TllXke9zqKEhSSOC/skDb5m6D2jvHaFofe4MT6LCAtRLPAwlHDm4Nj0MKad76uPzJYsFpI4foxRGJrH0HzGrLMYgkObJ2MwO7JYSOLmYQrS2EU+c8QhyzYbspjzC1qwH0zlW9i3WtgV5tYnZkMWCxYwhk61lvIzzJk4dkUQU+S7kMWC2WIqsyLK13Bo8tilgrix6ywW3ExMrTK8/A37II65mRetOKklSCndn1L6nZTSZ1JKn04p/dA6/IUppQ+mlD633n/NOjyllH4+pfRESukTKaVX7fohFhw/cs576WR2H97mlF/vfadElSwArAD8SM75nwJ4LYC3ppReCeARAB/KOT8A4EPrcwB4A4AH1tvDAN45aYkX3GgcYlT2Onzrdgjs6r5Vssg5P51z/sP18ZcBfAbAfQAeBPDYOtljAL5zffwggHfnS3wYwFenlO6dvOQLFizYwi4JqkVZXCGl9HIA3wTg9wG8JOf8NHBJKABevE52H4Av0GVPrsMWLJgEhxy154p91EmzgzOl9HwAvw7gh3POzxYcQV7E1lOklB7GpZmyYMGCa4AmZZFSuoVLoviVnPNvrIO/aObFev/MOvxJAPfT5S8D8JTmmXN+NOf86pzzq4cWfsHNxqIuLrGvemiZDUkAfgnAZ3LOP0tRjwN4aH38EID3Ufib17MirwXwd2auLFgwJQ69PmIu2Fc9pBorpZS+BcD/BvBJABfr4Lfj0m/xXgD/CMBfAPjunPPfrsnlFwB8B4C/B/C9OeePVu6xDBELurGQxV00qouPjVHyVbLYBxayWNCLhSi20dCXR5FF12zIggVzwEIUPnZdL8ty7wXXCocgiqm/gN0lSn8APRYLWSy4FtibE2/i+/T+4vjU910+JFtwo7Broti3Whnzl4mHxEIWC2aLXXXiufk8dvn/JVMqjIUsFswS+zAH5ohdEUdKaXR+syGLKR5mwXFgyo59XUjCw9z+MW02ZAHMr3IW7AdTdei55eNhaLs+lKOUMSuyYOxyCmjBPDBFpxybx6Gdm4Yh7XzffWS2ZGFYzJMFHoZ28rmaJWNU9b5IY/ZkASwq49gwpsMOufYQPpAxbXXo1OquzfhrQRaGRWVcb+yTJOZOKj3teMhguSzKwkIY1xX7MBv2QShTYYjDcohymLK/XDuyABaz5DphLiRxCEfoUN9Dy/U9xDFVf5kNWYyRWr3XLdg9dmkGlNJNkcdU6L1H658f1dLtqi/MhiwMi2q4/tiFOVBLMzZ+aNpetKqFUtraILmrQXR2ZGHoJY2FZA6PuZHErs2ToT4ED0M7fc0ncaN8FkNIYyGM/WNqohhKAocij57rasTgpSkRR4/SGIPZkMWUDLmojP1iDkTRG94a35uO0WJGRGlL5NFCHEe/zmJqSbWojN1jSl9AL0kcijhaUcqv10kZkUBv+BjMiiyA+kMuZsl8MNWoPFWnn4JQWuOHwtpiTV1E5BGpin2QxmzIQjt1qx3WUgmLWTIdpjQ5WjtyKwm0hA0lj560Q9dHjG2nOecmkhmK2ZAFEKuAUiX2ksZCGMMxVYcaQxIt6fZFGhFayaRHPXjtXMNa0ozBLMiCK2fXKmAhjN1jCt/DWEIYQw67MkEAf/SPwjXMOwdi0ojChmIWZAHEZkgrSy4KY7cYO+041J/QQwq9BDKkbLvAUKJoVRZHuc7Ce7DWsJa4IekW7IYoaunGKIeeuJ6wfaGXGGqkMdWzzIosgHZFUVMZC2FMg30oit6ZCTufijDmRBY1omiJY5I4yqnTFkkVXTeGMBbEGNthhpgFLSZG73FPXFSuUjkZ+2pzJaKIfCJjMTuyYEQPbnF67ZCpKM8ZdJMxxGSohU+hALzj1rDW45bzIdB25bUzC4vS1shhV2qCMRuyMLQ+eI8zZzFLpsUuiKLF1GgJ6yGYUngUFsW3DEgtqnlIG2whkilIbzZkUWJKL5yvm5IwFpSxa6Jo6fRD41rCWp+lJb6nrY4liShuKqIAZkYW3gNqZbYoijGEcZNJZUiHiMKHqoZafE1JTKk6ojAvrjSdz/FRu9brhpBIq9k+FLMgC35ZrTZYrRIjqbcQho+ho89QoqilbyWF3vDa/VrLXHoGRgtJtLT7MaRxlGZIpCz0QaMK5zjNtxZ2kzG0M/R0qqGmRnRdK1m0kk/r87R2Om+NQ6tvoUQiPZjSBAFmRhZAm3e3RBKtFduiTG4CoeyKKGojuZe2Z98axvfoVR1ReUvh0UKoHqKoqekeEhlDNooqWaSU7gfwbgAvBXAB4NGc839KKb0DwPcD+Kt10rfnnN+/vubHALwFwDmAH8w5f6DhPlcP1fuAJaeod873HGLKHAOmIgkNG6sgWjq/xg25prSPwoaC2483RcphGt8Sd3Fx0XzNGLQoixWAH8k5/2FK6QUAPpZS+uA67udyzj/NiVNKrwTwJgDfAOBrAfx2Sunrc87npZuU5JdKuVp6L0/vPAorle0YcGii6FUOUbyGn5ycFPPrIQuD5dlSF4DfKa0z83U9ZFFr5ycnJ8X2uTefRc75aQBPr4+/nFL6DID7Cpc8COA9OeevAPjzlNITAF4D4PdK96k9TIkkNJ+FMMahtWHViMJLO1Y5AD4plAihFM5gYhiqLGqEqmqA4y8uLtyB0ghjSifoEHT5LFJKLwfwTQB+H8DrALwtpfRmAB/Fpfr4Ei6J5MN02ZNwyCWl9DCAh9en/+/ZZ5/9GwB/3Vn+Q+FFuD5lBaS8LXbuAXGt63bm+CdjLm4mi5TS8wH8OoAfzjk/m1J6J4CfBJDX+58B8H0APCrean0550cBPEr5fzTn/Oq+4h8G16mswPUq73UqK3C9yptS+uiY632DbPsmt3BJFL+Sc/4NAMg5fzHnfJ5zvgDwi7g0NYBLJXE/Xf4yAE+NKeSCBQsOjypZpEuD6pcAfCbn/LMUfi8l+y4An1ofPw7gTSmle1JKrwDwAIA/mK7ICxYsOARazJDXAfh3AD6ZUvr4OuztAL4npfSNuDQxPg/gBwAg5/zplNJ7AfwxLmdS3lqbCVnj0XqS2eA6lRW4XuW9TmUFrld5R5U13RQP/4IFC8ahyWexYMGCBQcni5TSd6SUPptSeiKl9Mihy+MhpfT5lNInU0ofN49ySumFKaUPppQ+t95/zYHK9q6U0jMppU9RmFu2dImfX9f1J1JKr5pJed+RUvrLdf1+PKX0Ror7sXV5P5tS+ld7Luv9KaXfSSl9JqX06ZTSD2EF0qkAAAKCSURBVK3DZ1e/hbJOV7e84GPfG4BTAH8K4OsA3AbwRwBeecgyBeX8PIAXSdh/BPDI+vgRAP/hQGX7VgCvAvCpWtkAvBHA/8Dl9PZrAfz+TMr7DgD/3kn7ynWbuAfAK9Zt5XSPZb0XwKvWxy8A8CfrMs2ufgtlnaxuD60sXgPgiZzzn+WcnwPwHlyuAL0OeBDAY+vjxwB85yEKkXP+XQB/K8FR2R4E8O58iQ8D+GqZ1do5gvJGuFoNnHP+cwC2GngvyDk/nXP+w/XxlwHY6uXZ1W+hrBG66/bQZHEfgC/QubvacwbIAH4rpfSx9cpTAHhJvlwKj/X+xQcr3Taiss25vt+2lu7vIpNuNuWV1cuzrl8pKzBR3R6aLJpWe84Ar8s5vwrAGwC8NaX0rYcu0EDMtb7fCeAfA/hGXH6H9DPr8FmUV1cvl5I6YXstr1PWyer20GRxLVZ75pyfWu+fAfCbuJRrXzSJud4/c7gSbiEq2yzrO894NbC3ehkzrd9dr7Q+NFl8BMADKaVXpJRu4/LT9scPXKYNpJSely4/zUdK6XkAvh2Xq1UfB/DQOtlDAN53mBK6iMr2OIA3r732rwXwdyanD4m5rgaOVi9jhvW7l5XW+/LWFry4b8Sl5/ZPAfz4ocvjlO/rcOk1/iMAn7YyAvgHAD4E4HPr/QsPVL5fw6W8vIPL0eItUdlwKT3/87quPwng1TMp739dl+cT60Z8L6X/8XV5PwvgDXsu67fgUpp/AsDH19sb51i/hbJOVrfLCs4FCxY04dBmyIIFC64JFrJYsGBBExayWLBgQRMWsliwYEETFrJYsGBBExayWLBgQRMWsliwYEETFrJYsGBBE/4/F8HUSTnK+UMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10a700828>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(images[3000], cmap='gray')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "# compute the gif\n",
    "imageio.mimsave('images.gif', images)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x1c605b9898>]"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAD8CAYAAABZ/vJZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3XucXFWZ7//Pt+/pTncn6TS5dZJOyIVcgABNiArqgGDAkcgIQ5hRUfGHCPjSUc8xzMzhOJwz85Mz4+BxYBQUBKNCGFCJM0FEcIwyENNAQm4kNEkgTW6de1/Sl6p6zh+1q1N0+lJ9rV3dz/v1qlftWrX23k9Vd++n11p77yUzwznnnOutrHQH4JxzLjN5AnHOOdcnnkCcc871iScQ55xzfeIJxDnnXJ94AnHOOdcnnkCcc871iScQ55xzfeIJxDnnXJ/kpDuAwTR+/HirrKxMdxjOOZdRXn755UNmVt5TvWGdQCorK6murk53GM45l1EkvZVKPe/Ccs451yeeQJxzzvWJJxDnnHN9klICkbRU0nZJNZJWdPJ+vqRVwfvrJFUG5YslbQgeGyVdE5RPlfRbSdskbZH0paRtjZP0rKQ3guexQbkkfSfYx2uSzh+IL8A551zf9JhAJGUD9wFXAvOBGyTN71DtJuComc0C7gHuDso3A1VmtghYCtwvKQeIAF81s3nAEuC2pG2uAJ4zs9nAc8Frgv3PDh43A9/tw+d1zjk3QFJpgSwGasxsp5m1Ao8ByzrUWQY8Eiw/AVwmSWbWZGaRoLwAMAAz22dmrwTL9cA2YEon23oE+FhS+Y8s7iVgjKRJvfiszjnnBlAqCWQKsCfpdS2nDvan1QkSxnGgDEDSRZK2AJuAW5ISCsH7lcB5wLqgaIKZ7Qu2tQ84oxdxIOlmSdWSquvq6lL4eM455/oilQSiTso6zoPbZR0zW2dmC4ALgTskFbSvJI0GngS+bGYnBiAOzOwBM6sys6ry8h6vg8k4z2zZz29fP5juMJxzLqUEUgtMTXpdAeztqk4wxlEKHEmuYGbbgEZgYVAvl3jy+ImZ/Syp6oFE11TwnDhaphLHsPbim4f5/MqX+czD61m383C6w3HOjXCpJJD1wGxJMyTlAcuB1R3qrAZuDJavBZ43MwvWyQGQNB2YC+yWJOBBYJuZ/XM327oReCqp/FPB2VhLgOOJrq6RYtX6txmVm01JQQ73r92Z7nCccyNcj7cyMbOIpNuBZ4Bs4CEz2yLpLqDazFYTTwYrJdUQb3ksD1a/GFghqQ2IAbea2SFJFwOfBDZJ2hDU/WszWwN8E3hc0k3A28B1wftrgKuAGqAJ+Ex/P3ymWb/7KJeedQYVY0fx0Au7OH6yjdJRuekOyzk3QqV0L6zgwL6mQ9mdScvNnDrQJ9dZCazspPwPdD6mgZkdBi7rpNyA21KJdzg61NDCO8dO8pn3VXLetLHcv3Ynv9tRx9XnTk53aM65EcqvRM8Qbx1uAmDG+CLOrShldH6Oj4M459LKE0iG2HMknkCmlxWSk51FVeVY1u060sNazjk3eDyBZIhEC6RibCEAi2eMo+ZgA4caWtIZlnNuBPMEkiHePtLExJICCnKzAVhcOQ6AV946ms6wnHMjmCeQDHGwvpkJpe3XYLJgcilZgs3vHE9jVM65kcwTSIY43NDK+KK89tej8rKZfUYxmzyBOOfSxBNIhjjS2Mq4pAQCsHBKKZveOU78DGfnnBtankAygJnFE8jodyeQs6eUcKihlf0nmtMUmXNuJPMEkgHqWyK0RmOML8p/V/nZFaUAbKr1bizn3NDzBJIBjjS0ApzWhTV/UjCQvrenGxk759zA8wSSAQ43xhNIWYcurFF52VSOL+L1fZ5AnHNDzxNIBjjWFE8gYwrzTntv3sQSXt9fP9QhOeecJ5BM0NASn8SxuOD0e1+eNbGYt480tddxzrmh4gkkA9Q3Bwkkv5MEMqkEgO3eCnHODTFPIBngVAvk9Lk/zppYDMDr+30cxDk3tDyBZICG5gjZWaIg9/QfV8XYUYzOz+H1fd4Ccc4NrZQSiKSlkrZLqpG0opP38yWtCt5fJ6kyKF8saUPw2CjpmqR1HpJ0UNLmDttalbTO7sSMhZIqJZ1Meu97/fngmaS+uY3R+TnEZwJ+N0mcNbHYu7Ccc0OuxxkJJWUD9wGXA7XAekmrzWxrUrWbgKNmNkvScuBu4HpgM1AVTIs7Cdgo6ZdmFgEeBu4FfpS8PzO7Pmnf3wKSr5J708wW9eFzZrT6lgijOxn/SDhrUjFPbdiLmXWaZJxzbjCk0gJZDNSY2U4zawUeA5Z1qLMMeCRYfgK4TJLMrClIFgAFQPtNm8xsLfH50zul+JHwz4FHU/okw1hDc6TTM7ASzppYQn1zhL3H/ZYmzrmhk0oCmQLsSXpdG5R1WidIGMeBMgBJF0naAmwCbklKKD25BDhgZm8klc2Q9Kqk30m6pLOVJN0sqVpSdV1dXYq7CreGHlog84IzsfyCQufcUEolgXTWJ9Lx9q9d1jGzdWa2ALgQuENSQSd1O3MD72597AOmmdl5wFeAn0oqOW2nZg+YWZWZVZWXl6e4q3BraIkwupsWyNz2M7F8HMQ5N3RSSSC1wNSk1xXA3q7qSMoBSunQPWVm24BGYGFPOwy28WfAqqT1W8zscLD8MvAmMCeF+DNeQ3P3LZDR+TlMG1fINm+BOOeGUCoJZD0wW9IMSXnAcmB1hzqrgRuD5WuB583MgnVyACRNB+YCu1PY54eA182sNlEgqTwY0EfSTGA2sDOFbWW8hpYIRXndn+9w1sRib4E454ZUjwkkGLO4HXgG2AY8bmZbJN0l6eqg2oNAmaQa4t1LiVN9LyZ+5tUG4OfArWZ2CEDSo8CLwFxJtZJuStrtck4fPH8/8JqkjcQH6m8xsy4H4YeTk21RRuVld1vnrEkl7KxroLktOkRROedGuh5P4wUwszXAmg5ldyYtNwPXdbLeSmBlF9u8oZv9fbqTsieBJ1OJd7hpaYtRkNt9Apk3sZiYwY4D9ZxTMWaIInPOjWR+JXrIRWNGazTW6VXoyRJnYvk4iHNuqHgCCblEl9SoHlog08YVUpSXzTa/pYlzboh4Agm5k0EC6akLKytLzJ1YzFZvgTjnhognkJBLtQUC8W6sbftOYNbxMh3nnBt4nkBCLpFA8nsYA4F4AqlvjvDOsZODHZZzznkCCbvmthiQWgtk/uTEQLqPgzjnBp8nkJBLjIH0dB0IxC8mlGDrXh8Hcc4NPk8gIdec4iA6QGFeDpVlRX4qr3NuSHgCCbmTrakPogPMm1TMNp/e1jk3BDyBhFxzJD4G0tOFhAnzJpbw1uGm9nnUnXNusHgCCbnm1tS7sODUFenbvRXinBtknkBCrjnSywQSnIm11c/Ecs4NMk8gIdfbMZDJpQWUFOT4QLpzbtB5Agm5xHUgqbZAJDFvUomfyuucG3SeQELuZFuUvOwssrM6mzW4c/MmlbB9fz3RmN/SxDk3eDyBhFxzWzSl25gkmz+5hJNtUd463DhIUTnnXIoJRNJSSdsl1Uha0cn7+ZJWBe+vk1QZlC+WtCF4bJR0TdI6D0k6KGlzh219Q9I7SetdlfTeHcE+tkv6cF8/dCZpboumPP6RMD84E2uzd2M55wZRjwkkmIf8PuBKYD5wg6T5HardBBw1s1nAPcDdQflmoMrMFgFLgfsTc6QDDwdlnbnHzBYFjzVBHPOJT3W7IFjvXxNzpA9nzW3RlMc/EuZMKCYvJ4vN7xwfpKiccy61FshioMbMdppZK/AYsKxDnWXAI8HyE8BlkmRmTcGc6gAFQHunvJmtBXozp/ky4DEzazGzXUBNENuw1hqNkZfTuy6svJws5k0q4bXaY4MUlXPOpZZApgB7kl7XBmWd1gkSxnGgDEDSRZK2AJuAW5ISSndul/Ra0M01thdxIOlmSdWSquvq6lLYVbi1RmLkZfd+qOqcKaVsfucEMR9Id84NklSOTJ2d/tPxqNRlHTNbZ2YLgAuBOyQV9LC/7wJnAouAfcC3ehEHZvaAmVWZWVV5eXkPuwq/1qiR28sWCMDZFaU0tETY5QPpzrlBksqRqRaYmvS6AtjbVZ1gjKOUDt1TZrYNaAQWdrczMztgZlEziwHf51Q3VSpxDDutkSj5fWmBVJQCsKnWx0Gcc4MjlSPTemC2pBmS8ogPZK/uUGc1cGOwfC3wvJlZsE4OgKTpwFxgd3c7kzQp6eU1xAfiE/tYHpzxNQOYDfwxhfgzWlvUyM1J/RqQhFnloynIzeI1TyDOuUGS01MFM4tIuh14BsgGHjKzLZLuAqrNbDXwILBSUg3xlsfyYPWLgRWS2oAYcKuZHQKQ9CjwQWC8pFrgf5rZg8D/kbSIePfUbuDzQRxbJD0ObAUiwG1mFh2ILyHMWiMxSgp6/DGdJic7iwWTS9n0jg+kO+cGR0pHpuBU2jUdyu5MWm4GrutkvZXAyi62eUMX5Z/sJo6/B/4+lZiHi7Y+nIWVcPaUUlat30M0Zr26kt0551LhV6KHXGskRm4fxkAgPg5ysi3Km3UNAxyVc855Agm9vlwHkpAYSPdxEOfcYPAEEnJ9vQ4EYMb40RTlZfsFhc65QeEJJOT6MwaSnSUWTillwx5PIM65gecJJOT6MwYCcMH0sWzde6J9YirnnBsonkBCrj9jIBBPIJGYeTeWc27AeQIJMTOLX0jYjxbIedPitxJ7+e2jAxWWc84BnkBCrTUan842vx8tkHFFecwsL+KVtzyBOOcGlieQEGuLxu8VmZvdv4sAL5g2llfePoaZ35nXOTdwPIGEWGsk3gLp62m8CRdMH8uRxlZ2H24aiLCccw7wBBJqbUEXVl9u557s/OnBOIh3YznnBpAnkBAbqBbIrPLRFBfkeAJxzg0oTyAhlhhE789pvABZWeL8aWN51c/Ecs4NIE8gITZQLRCAqulj2X6gnmNNrf3elnPOgSeQUGsboBYIwJIzyzCDdbuO9FzZOedSkNKRSdJSSdsl1Uha0cn7+ZJWBe+vk1QZlC+WtCF4bJR0TdI6D0k6KGlzh239o6TXJb0m6eeSxgTllZJOJm3ve/354Jkg0QLpz4WECedUlFKQm8VLOw/3e1vOOQcpJBBJ2cB9wJXAfOAGSfM7VLsJOGpms4B7gLuD8s1AlZktApYC9yemuAUeDso6ehZYaGbnADuAO5Lee9PMFgWPW1L5gJlsoMZAAPJzsqmaPo4X3/QE4pwbGKkcmRYDNWa208xagceAZR3qLAMeCZafAC6TJDNrMrNIUF5AfJpaAMxsLfHpb9/FzH6dtM5LQEXKn2aYGcgWCMCSmeN4fX89Rxt9HMQ513+pHJmmAHuSXtcGZZ3WCQ7+x4EyAEkXSdoCbAJuSUoOqfgs8HTS6xmSXpX0O0mX9GI7GSlxJXp/bmWS7D1nlgGwbpe3Qpxz/ZfKkamz+2h0vCdGl3XMbJ2ZLQAuBO6QVJBKYJL+BogAPwmK9gHTzOw84CvATyWVdLLezZKqJVXX1dWlsqvQGugWyNlTxjAqN5uXdvpAunOu/1I5MtUCU5NeVwB7u6oTjHGU0qF7ysy2AY3Awp52KOlG4E+Bv7TgBk5m1mJmh4Pll4E3gTkd1zWzB8ysysyqysvLU/h44TWQZ2EltlNVOdbHQZxzAyKVI9N6YLakGZLygOXA6g51VgM3BsvXAs+bmQXr5ABImg7MBXZ3tzNJS4GvA1ebWVNSeXkwoI+kmcBsYGcK8WesUy2Q/t1MMdmSmWVsP1DP4YaWAdumc25k6jGBBGMWtwPPANuAx81si6S7JF0dVHsQKJNUQ7x7KXGq78XARkkbgJ8Dt5rZIQBJjwIvAnMl1Uq6KVjnXqAYeLbD6brvB16TtJH4QP0tZjas+2Laz8IaoC4sODUO8l/eCnHO9VNOz1XAzNYAazqU3Zm03Axc18l6K4GVXWzzhi7KZ3VR/iTwZCrxDheR6MCOgQCcM6WU0lG5rN1Rx0fPnTxg23XOjTx+JXqIRWLxcxWyB7ALKyc7i4tnj+d3O+p8fhDnXL94Agmx9gmlsgb2x/SBOeUcrG/h9f31A7pd59zI4gkkxKKxeBdWdtbAtUAA3j87fnba2h2ZfZqzcy69PIGE2EBNadvRxNICzppYzO88gTjn+sETSIhFYjGys4Q0sAkE4t1Y63cfobGlNzcGcM65UzyBhFgkZgPefZXwgTnltEXNLyp0zvWZJ5AQi0SN3EFKIBdUjqUwL5vfbj84KNt3zg1/nkBCLBozcgbwGpBk+TnZfGBOOb/ZdoBYzE/ndc71nieQEGuLxsgZpBYIwBULJnDgRAuvvXN80PbhnBu+PIGEWCRq5AzwGVjJLp07gews8est+wdtH8654csTSIhFYkbOAF9EmKy0MJclM8fx660HBm0fzrnhyxNIiEVisUFtgQBcPm8CNQcbeLOuYVD345wbfjyBhFgkaoM6BgJw+YKJADzrrRDnXC95AgmxSCw2oHfi7cyUMaNYOKWEX232cRDnXO94AgmxSHTwLiRMduXCSWzYc4w9R5p6ruyccwFPICHWNojXgSS7OpgX5N9f2zfo+3LODR8pHZ0kLZW0XVKNpBWdvJ8vaVXw/jpJlUH54mBWwQ2SNkq6JmmdhyQdlLS5w7bGSXpW0hvB89igXJK+E+zjNUnn9+eDZ4JobHCvA0mYOq6Q86aNYfXGjlPdO+dc13pMIME85PcBVwLzgRskze9Q7SbgaDCb4D3A3UH5ZqDKzBYBS4H7E3OkAw8HZR2tAJ4zs9nAc5yaHvdK4vOgzwZuBr6bygfMZG1DMIiecPW5k9m27wQ1B32OEOdcalJpgSwGasxsp5m1Ao8ByzrUWQY8Eiw/AVwmSWbWFMypDlAAtN8zw8zWAp3NaZ68rUeAjyWV/8jiXgLGSJqUQvwZKxId/EH0hI+cPYksweqN3o3lnEtNKkenKcCepNe1QVmndYKEcRwoA5B0kaQtwCbglqSE0pUJZrYv2NY+4IxexDGsRAfxbrwdnVFSwJKZZfxy416f6tY5l5JUEkhnR7COR5gu65jZOjNbAFwI3CGpoHch9ioOJN0sqVpSdV1dZk+Y1Ba1AZ9MqjvLFk1m16FGNuw5NmT7dM5lrlQSSC0wNel1BdBxtLW9TjDGUUqH7ikz2wY0Agt72N+BRNdU8Jy433gqcWBmD5hZlZlVlZeX97CrcIvEYoN6K5OOrjp7EqNys3m8unbI9umcy1ypHJ3WA7MlzZCUBywHVneosxq4MVi+FnjezCxYJwdA0nRgLrC7h/0lb+tG4Kmk8k8FZ2MtAY4nurqGq0jMyB7CFkhxQS5XnT2JX27cS1Orz1TonOtejwkkGLO4HXgG2AY8bmZbJN0l6eqg2oNAmaQa4CucOnPqYmCjpA3Az4FbzewQgKRHgReBuZJqJd0UrPNN4HJJbwCXB68B1gA7gRrg+8Ct/fjcGWEwJ5TqyvUXTqWhJcLTm/zKdOdc93J6rgJmtob4ATy57M6k5Wbguk7WWwms7GKbN3RRfhi4rJNyA25LJd7hYjAnlOrKhZVjmTG+iFXVe/j4BRVDum/nXGbxK9FDbLAnlOqMJK6rquCPu46w61DjkO7bOZdZPIGEWCQ2uBNKdeXj51eQnSUeW//2kO/bOZc5PIGEWCQ6tGdhJUwoKeBD885g1fo9NLdFh3z/zrnM4AkkxOIzEg59CwTg0++dwbGmNlZv8PtjOec65wkkxOJzoqfnR7Rk5jjmTijm4f/a7VemO+c65QkkxOITSqWnBSKJG99bydZ9J1i/+2haYnDOhZsnkJCKxYyYMWT3wurMx86bTOmoXB75r91pi8E5F16eQEKqLRYDGLK78XamMC+H5RdO5Vdb9vtshc6503gCCaloLD7ukM4WCMCn31dJluCBtTvTGodzLnw8gYRUWzSeQNJ1FlbCpNJR/Nl5FTxevYe6+pa0xuKcCxdPICEViaa/Cyvh8x+YSWs0xkMv7Ep3KM65EEn/0cl1KixdWAAzy0dz1cJJ/PjFtzjR3JbucJxzIeEJJKTaggSSrtN4O/rCB8+kviXCwy/sTncozrmQ8AQSUtH2MZBw/IgWTinl8vkT+P7anRxrak13OM65EAjH0cmdJnEabzpuptiVr14xh4bWCPf7GVnOOTyBhFYkZC0QgLMmlnD1uZP54Qu7OFjfnO5wnHNpltLRSdJSSdsl1Uha0cn7+ZJWBe+vk1QZlC+WtCF4bJR0TU/blPT7pHX2SvpFUP5BSceT3ruzYxzDSSSELRCAv/rQHNqixn3P16Q7FOdcmvWYQCRlA/cBVwLzgRskze9Q7SbgqJnNAu4B7g7KNwNVZrYIWArcLymnu22a2SVmtihY50XgZ0n7+X3iPTO7q4+fOSNEQnIdSEeV44v486oKfvrHt33CKedGuFRaIIuBGjPbaWatwGPAsg51lgGPBMtPAJdJkpk1BXOqAxQAidu69rhNScXApcAvevuhhoNTLZDwdGEl/NWH5pCXncXf/8e2dIfinEujVI5OU4A9Sa9rg7JO6wQJ4zhQBiDpIklbgE3ALcH7qWzzGuA5MzuRVPaeoCvsaUkLUog9YyVaILkha4EAnFFSwG2XzuI32w7w+zfq0h2Ocy5NUkkgnR3BOk4Q0WUdM1tnZguAC4E7JBWkuM0bgEeTXr8CTDezc4F/oYuWiaSbJVVLqq6ry9yDWyREFxJ25rPvm8G0cYX8r3/f2n7VvHNuZEklgdQCU5NeVwAdp6lrryMpBygFjiRXMLNtQCOwsKdtSioj3s31H0nrnzCzhmB5DZAraXzHYM3sATOrMrOq8vLyFD5eOLVFw9uFBVCQm81fXzWPHQca+PFLb6U7HOdcGqRydFoPzJY0Q1IesBxY3aHOauDGYPla4Hkzs2CdHABJ04G5wO4Utnkd8O9m1n6uqKSJkhQsLw5iP9yrT5tBErcyCdsgerIPL5jAJbPH80+/3sH+435ar3MjTY8JJBizuB14BtgGPG5mWyTdJenqoNqDQJmkGuArQOK03IuBjZI2AD8HbjWzQ11tM2m3y3l39xXEE9NmSRuB7wDLbRjPtdp+N96QncabTBL/+2MLaYvGuPOpzekOxzk3xHJSqRR0Ga3pUHZn0nIz8VZDx/VWAitT3WbSex/spOxe4N5U4h0OIiGYUCoV08uK+PKH5nD3r17nV5v3s3ThxHSH5JwbIuE+Oo1gYbobb08+d8kM5k0q4c6nNvvdep0bQTyBhFRb+2m84f8R5WZncffHz+ZQQwt/t3prusNxzg2R8B+dRqhoSG9l0pVzKsZw+5/M4slXalmzaV+6w3HODQFPICEVlilte+OLl83m3IpS/vrnmzhwws/Kcm648wQSUpGQXwfSmdzsLO65fhEtbTG+9m8bicWG7Ulyzjk8gYRW4kr0TOnCSphZPpq//dN5/P6NQz5viHPDnCeQkIpkwIWEXfmLxdP4yDmT+MdnXufFN4fttZ7OjXieQEKqvQsrA87C6kgSd3/8HGaML+KLj77i4yHODVOZd3QaITK5BQIwOj+H733iAppao9z2k1dojfgNF50bbjyBhFQkamQJsjI0gQDMnlDMNz9+DtVvHeVvf7GJYXznGedGpJRuZeKGXlssllFnYHXl6nMnU3Ognu88X8OM8aP5wgfPTHdIzrkB4gkkpKJRC+VkUn3xV5fPYdfhJu7+1evMGF/I0oWT0h2Sc24AZP6/uMNUJGYZcR+sVEjiH689h/OmjeHLqzZQvftIzys550LPE0hItUVjob8Tb28U5Gbz/U9VMal0FJ95eD1b9h5Pd0jOuX4aPkeoYSY6jFogCeNH5/Pjz13E6Pwcbnzoj+ysa0h3SM65fvAEElJtURtWLZCEKWNGsfKmizCDT/xgHbVHm9IdknOuj1I6QklaKmm7pBpJKzp5P1/SquD9dZIqg/LFkjYEj42Srulpm5IelrQrab1FQbkkfSeo/5qk8/v74cMsEotl3G1MUjXrjNE88tnF1LdEuP7+l3jrcGO6Q3LO9UGPCURSNnAfcCUwH7hB0vwO1W4CjprZLOAe4O6gfDNQZWaLgKXA/ZJyUtjmfzOzRcFjQ1B2JTA7eNwMfLf3HzdzDKdB9M4snFLKo//fEhpb40nkTe/Oci7jpNICWQzUmNlOM2sFHgOWdaizDHgkWH4CuEySzKwpmP8coABIXEmWyjY7Wgb8yOJeAsZIGrbng0aisYyYTKo/Fk4p5bGbl9AWjXH9/S+xfX99ukNyzvVCKkeoKcCepNe1QVmndYKEcRwoA5B0kaQtwCbgluD9nrb590E31T2S8nsRx7ARjdmw7cJKdtbEElZ9/j1kCa773n+xbqfffNG5TJFKAunsKNbxnhRd1jGzdWa2ALgQuENSQQ/bvAM4K6g/Dvh6L+JA0s2SqiVV19XVdbJKZmiLWsbeB6u3Zp0xmie/8F7GF+fzyQf/yC837k13SM65FKSSQGqBqUmvK4COf+HtdSTlAKXAu64WM7NtQCOwsLttmtm+oJuqBfgh8e6uVOPAzB4wsyozqyovL0/h44VTZJjcyiRVU8cV8rMvvJdzp5byxUdf5YG1b/q9s5wLuVSOUOuB2ZJmSMoDlgOrO9RZDdwYLF8LPG9mFqyTAyBpOjAX2N3dNhPjGpIEfIz4QHxiH58KzsZaAhw3s2E7+XZkBLVAEsYU5rHypov4yNmT+Ic1r/Pfn3iN5rZousNyznWhx3thmVlE0u3AM0A28JCZbZF0F1BtZquBB4GVkmqItzyWB6tfDKyQ1AbEgFvN7BBAZ9sM1vmJpHLiXVYbgFuC8jXAVUAN0AR8pn8fPdwiMaMgd+S0QBIKcrP5lxvO48zyIr7zfA3bD9TzvU9cwOQxo9IdmnOuAw3nboKqqiqrrq5Odxh9suzePzCmMI9HPru458rD1DNb9vPVxzdSkJvFvX9xPktmlqU7JOdGBEkvm1lVT/VG3r+4GSJ+JfrI6sLq6MMLJvKL295LSUEuf/H9l/j2b3a0z9TonEs/TyAhFYnFMnI624E264xinrr9fSxbNIVv/+YNbvj+S7xz7GS6w3LO4QkktCLRkXEdSCqKC3K55/pF3HP9uWzde4Irv72WX27c62dpOZdfxoPCAAATqUlEQVRmnkBCqi02vG7nPhCuOa+CNV+6hBnlo/nio69yy49f5mB9c7rDcm7E8iNUSI3E03hTMb2siCdveQ8rrjyL326v4/J/XsuTL9d6a8S5NPAEElJtURtRFxL2Rk52Frd84Eye/tIlzD5jNF/9t418yucXcW7I+REqpCKx2Ig/C6snZ5aP5vHPv4dvfHQ+G94+xoe/vZa7f/U6jS2Rnld2zvWbJ5CQindh+Y+nJ1lZ4tPvm8HzX/sgyxZN4bv/+SaXfet3PLXhHe/Wcm6Q+REqpNqiMXJzvAWSqvLifP7punODmzLm8aXHNnD1vS/whzcOpTs054YtTyAhFYnZsJ8PZDBcMH0sT912Md+67lyONLbyiQfX8YkfrGNT7fF0h+bcsONHqBAysxEzH8hgyM4SH7+ggue++gH+9iPz2LL3OB+99w/csvJlNr/jicS5gdLjzRTd0GuLxvvu/TqQ/inIzeZzl8zkzy+cyg/W7uSHL+zmV1v28ydzy7n90llcMH1cukN0LqP5ESqE2oL7Pfl1IAOjpCCXr1wxlz+suJSvXTGHDXuO8fHvvsgND7zEf24/SCzmg+3O9YW3QEIoErRA/DqQgVU6KpfbL53NZ943g0f/+DYPrN3Jp3+4njPLi/j0+2bw8fOnUJjnfxLOpcqPUCHUFou3QPw6kMFRlJ/D5y6ZyR++fin3XH8uRfk5/I9fbGbJPzzH/79mG3uONKU7ROcygv+7FULtLRA/C2tQ5eVkcc15FXxs0RRefusoP3xhNz/4wy4e+P1OLp41nusvnMrl8yeQn5Od7lCdC6WUjlCSlkraLqlG0opO3s+XtCp4f52kyqB8saQNwWOjpGt62qaknwTlmyU9JCk3KP+gpONJ27uzvx8+rNrHQLwFMiQkUVU5jvv+8nx+/9//hC9fNoeddY3c/tNXWfIPz/G//n0rOw7UpztM50KnxxaIpGzgPuByoBZYL2m1mW1NqnYTcNTMZklaDtwNXE98PvOqYFrcScBGSb8ErJtt/gT4RLDdnwKfA74bvP69mf1p/z5y+EViibOwPIEMtcljRvGlD83m9ktn8ULNIVat38OPXtzNg3/YxfxJJVy9aDIfPXcyU3yKXedS6sJaDNSY2U4ASY8By4DkBLIM+Eaw/ARwrySZWXJncgHxxNHtNs1sTWIFSX8EKnr7oTJdpP0sLO/CSpfsLPH+OeW8f045hxtaeGrDXlZv3Ms3n36dbz79Oosrx/HRRZP5yNmTGFeUl+5wnUuLVBLIFGBP0uta4KKu6gStjeNAGXBI0kXAQ8B04JPB+z1uM+i6+iTwpaTi90jaCOwFvmZmW1KIP+Ocug7EWyBhUDY6n89ePIPPXjyDtw438suNe/nFhr38j19s5hurt7C4chxXLJjA5fMnUDG2MN3hOjdkUkkgnR3FOp4432UdM1sHLJA0D3hE0tMpbvNfgbVm9vvg9SvAdDNrkHQV8Atg9mnBSjcDNwNMmzat808UcpGYt0DCanpZEbdfOpvb/mQW2/bV8x+b9vLs1gP83S+38ne/3MqCySVcMX8iVyyYwFkTi5H8nwA3fKWSQGqBqUmvK4i3ADqrUyspBygFjiRXMLNtkhqBhT1tU9L/BMqBzyetfyJpeY2kf5U03szedbc8M3sAeACgqqoqI68Qa2+B5HgCCStJzJ9cwvzJJfy3D5/FrkONPLt1P7/ecoBvP7eDe36zg0mlBbx/djkfmFvO+84cT2lhbrrDdm5ApZJA1gOzJc0A3gGWA3/Roc5q4EbgReBa4Hkzs2CdPUG31XRgLrAbONbVNiV9DvgwcJmZxRI7kDQROBBsdzHxM8gO9+1jh1tiDCTXr0TPGDPGF3Hz+8/k5vefSV19C89tO8DvdtSxZvM+VlXvIUuwaOqY9nGVc6aU+oWiLuP1mECCg//twDNANvCQmW2RdBdQbWargQeBlZJqiLc8lgerXwyskNQGxIBbEy2GzrYZrPM94C3gxaD5/zMzu4t4YvqCpAhwElhuw3TCh8RZWH6AyUzlxfksXzyN5YunEYnG2LDnGGt31PG7Nw7xf597g2//5g1G5+dwwfSxXDRzHEtmlnH2lFK/95nLOBqmx2Ag3oVVXV2d7jB67T+3H+TTP1zPz259L+dPG5vucNwAOtLYygs1h1i36zDrdh7hjYPxaXgL87LjCWXGOC6sHMfZFaV+WxWXNpJeNrOqnur5b2gIJa5E9/lAhp9xRXl89Nz4tSQAhxpa+OOuI6zbeZh1u47wT7/eAcRPIz5rYjHnTRvDeVPHct60McwYX+SD8i5UPIGEUPtZWH4a77A3fnQ+V509iavOngTA0cZWXt1zlFffPsarbx/jF6/u5ccvvQ3AmMJcFk0dw7kVY1gwuYQFU0qZXFrgScWljSeQEPLrQEausUV5XHrWBC49awIA0ZjxZl0Dr74dTyqvvH2UtTvqSNyBfkxhbjyZTC5tf54xvohsPwHDDQFPICHk14G4hOwsMWdCMXMmFHP9hfHrmppaI2zbV8/WvcfZsvcEW/ae4OEXdtManL03KjebOROLmXPGaOZMKGb2hPjzJG+tuAHmCSSE2trnA/E/dne6wrz4GVwXTD91gkVbNEbNwQa27D3B5neOs31/Pb/dfpB/e7m2vU5xfg6zJoxmzhmnksrM8iIml44iy1ssrg88gYRQayT+n2SeX0joUpSbncW8SSXMm1TCtRecun3ckcZWdhyo540D9ew40MCOA/U8u+0Aq6pP3UkoPyeL6WWFVJYVMaO8iBllRcwYH3+UF+d7q8V1yRNICLUECaQg1+ehcP0zriiPJTPLWDKz7F3lhxpa2HGgnt2Hmth1qIFdh5p4s66B324/2N4CBijKy6ZyfBGVZUVUjB1FxbhCKsaOYurY+LP/jo5snkBCqCUSBeL/GTo3GMaPzmf86Hzee+a7yyPRGHuPNbPrcCO7DzWy61AjOw81snXfCZ7deqB9nCV5O1PHjaJibCFTxwbP40YxqXQUE0sLGJ3vh5jhzH+6IdTSFnRh+ZXJbojlZGcxrayQaWWFfGBO+bvei8WMg/Ut1B5tovboSfYciT/XHmti455jPL1pX/tdFBKK83OYWFrAxNICJpUWMLF0VPBcwMSSeFnpqFzvJstQnkBCqDkSJT8ny/+oXKhkZak9GVRVnv5+NGYcONHMniNN7D/RzL7jzewPHvtONLPjQB0H61voePOLgtwsJpYUcEZxAeOL8ygfnU95cfwxPmm5rCjfxwVDxhNICLW0xbz7ymWc7CwxecwoJnczW2NbNEZdfUt7ctl3/CQHgmRzqKGF7fvr+UP9IU40Rzpdf0xh7mkJZlxRHuOK8hhbmBcs5zK2MI8xhXl+Pcwg8wQSQi2RGPk+OOmGodzsrB6TDEBzW5TDja3U1bdQV9/CoYaW9uXE6w17jlFX30JTa7TTbUhQOiqXcYV5jG1PMLmMLcprLxtXmEdpYS4lBbmUjoo/CnK99Z8qTyAh1BJ0YTk3UhXkZjNlzKiU5p5vbotytKmVI42tHG1s40hTK0cbg9dNp57fOXaSze8c50hj62knAyTLy86iZFQuJaNy2pNK6ah3J5nS4P2SpPeKC3Ioys8ZUXdV9gQSQi2RmJ8e6VyKCnKzmVQaP/MrFWZGY2u0PcmcaG7j+MlTjxMnI8Fz/PWRxlZ2HWpsL4v1cAPz/JwsigtyGJ2fw+jEc34uo/Ozg9e5p95/V534cnF+PBGNys0O/QWenkBCqLnVWyDODRZJ7QfsqeN6N4d9LGY0tkZOTzbNbTS2RGhojtDQEqE+abmhJcLeYyfbl+ub2951rU13RuVmU5iXTWF+NkV5OYzKiz8X5iXKcyjMjT8XJcqC96eOK2ThlNK+fEUp8wQSQsdPtlE6yqc/dS5ssrJEcUEuxQW5VPRjqp6WSPRUsmmOxJNP8DjRHKGpJUJTa5Sm1giNrVFOtkZpbIlwsi3+fKihJXg/XqezcaCPnjuZf7nhvH582p6llEAkLQX+L/HZA39gZt/s8H4+8CPgAuLTzF5vZruDqWcfSFQDvmFmP+9um8E0t48B44BXgE+aWWtX++jrBw+z4yfbmHXG6HSH4ZwbJPk52eSPzqZsdP6AbC8WM5ojURpbgmTTGqEwb/C7wXvsJ5GUDdwHXAnMB26QNL9DtZuAo2Y2C7gHuDso3wxUmdkiYClwv6ScHrZ5N3CPmc0Gjgbb7m4fw463QJxzvZGVJQrzcigvzmdaWSHzJpUwvaxo8PebQp3FQI2Z7TSzVuKtg2Ud6iwDHgmWnwAukyQzazKzxAndBUCi46/TbSp+7tylwTYItvmx7vaR6gfNFGbGMU8gzrkMkEoX1hRgT9LrWuCiruqYWUTScaAMOCTpIuAhYDrx7qiIpK62WQYcS0o6tcG2u91HKh+0N17ff4Iv/vTVgd5sSpojUVojMSp6ObjnnHNDLZUE0tl/+R1PIeiyjpmtAxZImgc8Iunpbup3t69U4kDSzcDNANOmTetklZ4V5GQze0L6xiAumV3O1edMTtv+nXMuFakkkFpgatLrCmBvF3VqJeUApcCR5Apmtk1SI7Cwm20eAsZIyglaIcn76nEfwX4eIBi4r6qqSu1cuQ4qxxfxr395QV9Wdc65ESOVMZD1wGxJMyTlAcuB1R3qrAZuDJavBZ43MwvWyQGQNB2YC+zuaptmZsBvg20QbPOp7vbRq0/rnHNuwPTYAgnGG24HniF+yu1DZrZF0l1AtZmtBh4EVkqqId4qWB6sfjGwQlIbEANuNbNDAJ1tM1jn68Bjkv438GqwbbrZh3POuTTQcP4nvqqqyqqrq9MdhnPOZRRJL5tZVU/1/H4Zzjnn+sQTiHPOuT7xBOKcc65PPIE455zrE08gzjnn+mRYn4UlqQ54qx+bGM8g3CplAHhcveNx9Y7H1TvDMa7pZlbeU6VhnUD6S1J1KqeyDTWPq3c8rt7xuHpnJMflXVjOOef6xBOIc865PvEE0r0Heq6SFh5X73hcveNx9c6IjcvHQJxzzvWJt0Ccc871iSeQTkhaKmm7pBpJK4Z431Ml/VbSNklbJH0pKP+GpHckbQgeVyWtc0cQ63ZJHx7E2HZL2hTsvzooGyfpWUlvBM9jg3JJ+k4Q12uSzh+kmOYmfScbJJ2Q9OV0fF+SHpJ0UNLmpLJefz+SbgzqvyHpxs72NQBx/aOk14N9/1zSmKC8UtLJpO/te0nrXBD8/GuC2Ps1pXQXcfX65zbQf69dxLUqKabdkjYE5UP5fXV1bEjf75iZ+SPpQfz28m8CM4E8YCMwfwj3Pwk4P1guBnYA84FvAF/rpP78IMZ8YEYQe/YgxbYbGN+h7P8AK4LlFcDdwfJVQGL2ySXAuiH62e0nPn3ykH9fwPuB84HNff1+gHHAzuB5bLA8dhDiugLICZbvToqrMrleh+38EXhPEPPTwJWDEFevfm6D8ffaWVwd3v8WcGcavq+ujg1p+x3zFsjpFgM1ZrbTzFqBx4BlQ7VzM9tnZq8Ey/XANk7NC9+ZZcBjZtZiZruAGuKfYagsAx4Jlh8BPpZU/iOLe4n4TJOTBjmWy4A3zay7i0cH7fsys7WcPktmb7+fDwPPmtkRMzsKPAssHei4zOzXFp/1E+Al4rN/dimIrcTMXrT4UehHSZ9lwOLqRlc/twH/e+0urqAV8efAo91tY5C+r66ODWn7HfMEcropwJ6k17V0fwAfNJIqgfOAdUHR7UFT9KFEM5WhjdeAX0t6WfG55wEmmNk+iP+CA2ekIa6E5bz7Dzvd3xf0/vtJx/f2WeL/qSbMkPSqpN9JuiQomxLEMhRx9ebnNtTf1yXAATN7I6lsyL+vDseGtP2OeQI5XWf9lEN+qpqk0cCTwJfN7ATwXeBMYBGwj3gzGoY23veZ2fnAlcBtkt7fTd0h/R4Vnxr5auDfgqIwfF/d6SqOof7e/gaIAD8JivYB08zsPOArwE8llQxhXL39uQ31z/MG3v1PypB/X50cG7qs2kUMAxabJ5DT1QJTk15XAHuHMgBJucR/QX5iZj8DMLMDZhY1sxjwfU51uwxZvGa2N3g+CPw8iOFAomsqeD441HEFrgReMbMDQYxp/74Cvf1+hiy+YPD0T4G/DLpZCLqIDgfLLxMfX5gTxJXczTUocfXh5zaU31cO8GfAqqR4h/T76uzYQBp/xzyBnG49MFvSjOC/2uXA6qHaedDH+iCwzcz+Oak8efzgGiBxhshqYLmkfEkzgNnEB+8GOq4iScWJZeKDsJuD/SfO4rgReCoprk8FZ4IsAY4nmtmD5F3/Gab7+0rS2+/nGeAKSWOD7psrgrIBJWkp8HXgajNrSiovl5QdLM8k/v3sDGKrl7Qk+B39VNJnGci4evtzG8q/1w8Br5tZe9fUUH5fXR0bSOfvWH/OChiuD+JnL+wg/t/E3wzxvi8m3px8DdgQPK4CVgKbgvLVwKSkdf4miHU7/TzTo5u4ZhI/w2UjsCXxvQBlwHPAG8HzuKBcwH1BXJuAqkH8zgqBw0BpUtmQf1/EE9g+oI34f3k39eX7IT4mURM8PjNIcdUQ7wdP/I59L6j78eDnuxF4Bfho0naqiB/Q3wTuJbgQeYDj6vXPbaD/XjuLKyh/GLilQ92h/L66Ojak7XfMr0R3zjnXJ96F5Zxzrk88gTjnnOsTTyDOOef6xBOIc865PvEE4pxzrk88gTjnnOsTTyDOOef6xBOIc865Pvl/SO5DlTLUUvcAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1c545280f0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(point_value[:2000,1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}