gradient_descent.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "9f4f1972",
   "metadata": {},
   "source": [
    "## ECE6960 Homework2 Yufeng Yang"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a13ab0a9",
   "metadata": {},
   "source": [
    "### Package Preparation\n",
    "Import numpy for matrix computing.\n",
    "Import matplotlib for plotting"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "id": "fe6c05e4",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5af9b4f9",
   "metadata": {},
   "source": [
    "### Data Generation\n",
    "Specify m and n, sigma(variance for noise) for the data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "id": "74a07122",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Specify m, n\n",
    "m = 500\n",
    "n = 50\n",
    "# Specify the mean and variance of noise \\delta, which follows normal distribution.\n",
    "mu = 0\n",
    "# This can be changed.\n",
    "sigma = 0.8\n",
    "np.random.seed(1234)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "id": "512c273c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# define a function to generate experimental data.\n",
    "# Input args: m,n are matrix shape; mu and sigma are for noise data.\n",
    "def data_generator(m,n,mu, sigma):\n",
    "    # A is a matrix with size R^{m*n}\n",
    "    \n",
    "    A = np.empty([m,n])\n",
    "    # ground truth \n",
    "    z = 10*np.random.rand(n,1)\n",
    "    # create noise vector with the length equals to m\n",
    "    for i in range(m):\n",
    "        for j in range(n):\n",
    "            A[i,j] = np.random.normal(mu, 1)\n",
    "        norm_A_i = np.linalg.norm(A[i,:],2)\n",
    "        # scale each row of A such that the norm of each row is equal to 1.\n",
    "        A[i,:] = A[i,:]/norm_A_i\n",
    "    delta = np.random.normal(mu, sigma, size=(m,1))\n",
    "    return A, z, delta\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "id": "c6ffc318",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(500, 50)\n",
      "(50, 1)\n",
      "(500, 1)\n",
      "(500, 1)\n",
      "1.0\n"
     ]
    }
   ],
   "source": [
    "A, z, delta = data_generator(m,n,mu,sigma)\n",
    "# generate observations\n",
    "b = np.matmul(A,z)+delta\n",
    "# check about A,z,b's shape\n",
    "print(A.shape)\n",
    "print(z.shape)\n",
    "print(delta.shape)\n",
    "print(b.shape)\n",
    "print(np.linalg.norm(A[1,:],2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 74,
   "id": "f70f5081",
   "metadata": {},
   "outputs": [],
   "source": [
    "#print(A)\n",
    "#print(\"---***---\")\n",
    "#print(z)\n",
    "#print(\"---***---\")\n",
    "#print(delta)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 75,
   "id": "613169e2",
   "metadata": {},
   "outputs": [],
   "source": [
    "# define a function compute least square loss function value.\n",
    "def least_square(A,x,b):\n",
    "    residual = np.matmul(A,x)-b\n",
    "    loss = 0.5*np.linalg.norm(residual,2)**2\n",
    "    return loss"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 76,
   "id": "d5b6e986",
   "metadata": {},
   "outputs": [],
   "source": [
    "# define a function compute gradient\n",
    "def least_square_grad(A,x,b):\n",
    "    residual = np.matmul(A,x)-b\n",
    "    # gradient is A^T(Ax-b)\n",
    "    grad = np.matmul(np.transpose(A),residual)\n",
    "    return grad"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 77,
   "id": "56d4e639",
   "metadata": {},
   "outputs": [],
   "source": [
    "# generate initial guess\n",
    "x0 = z+5*np.random.normal(1,1,size = (z.shape[0],z.shape[1]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 78,
   "id": "075c22ef",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-1.83199015]\n",
      " [12.31107652]\n",
      " [ 9.38727108]\n",
      " [14.80615973]\n",
      " [15.73165837]\n",
      " [ 5.7759575 ]\n",
      " [-2.12012345]\n",
      " [16.60972677]\n",
      " [19.2831095 ]\n",
      " [17.68120053]\n",
      " [ 1.69597812]\n",
      " [ 8.35374792]\n",
      " [ 8.30220288]\n",
      " [19.42893025]\n",
      " [ 8.47250201]\n",
      " [17.88888907]\n",
      " [ 3.40437308]\n",
      " [ 3.85503629]\n",
      " [17.79712773]\n",
      " [14.90010191]\n",
      " [14.46398133]\n",
      " [15.71520615]\n",
      " [11.98141918]\n",
      " [16.57216393]\n",
      " [21.50076316]\n",
      " [10.4507935 ]\n",
      " [13.94751482]\n",
      " [ 8.16827851]\n",
      " [ 2.29515193]\n",
      " [ 9.26120032]\n",
      " [ 5.08596502]\n",
      " [ 8.05424058]\n",
      " [13.38487195]\n",
      " [ 9.17213981]\n",
      " [19.9165741 ]\n",
      " [13.3769188 ]\n",
      " [13.00540409]\n",
      " [ 7.64899932]\n",
      " [12.69516724]\n",
      " [17.20771868]\n",
      " [-0.71844652]\n",
      " [-0.86499653]\n",
      " [ 9.23951829]\n",
      " [ 2.24459299]\n",
      " [16.95257829]\n",
      " [ 6.92885652]\n",
      " [-0.77029597]\n",
      " [12.2173633 ]\n",
      " [ 4.24657174]\n",
      " [14.26345398]]\n",
      "[[1.9151945 ]\n",
      " [6.22108771]\n",
      " [4.37727739]\n",
      " [7.85358584]\n",
      " [7.79975808]\n",
      " [2.72592605]\n",
      " [2.76464255]\n",
      " [8.01872178]\n",
      " [9.58139354]\n",
      " [8.75932635]\n",
      " [3.5781727 ]\n",
      " [5.00995126]\n",
      " [6.83462935]\n",
      " [7.12702027]\n",
      " [3.70250755]\n",
      " [5.61196186]\n",
      " [5.03083165]\n",
      " [0.1376845 ]\n",
      " [7.72826622]\n",
      " [8.82641191]\n",
      " [3.64885984]\n",
      " [6.15396178]\n",
      " [0.75381242]\n",
      " [3.68824006]\n",
      " [9.33140102]\n",
      " [6.51378143]\n",
      " [3.97202578]\n",
      " [7.88730143]\n",
      " [3.16836122]\n",
      " [5.68098653]\n",
      " [8.6912739 ]\n",
      " [4.36173424]\n",
      " [8.02147642]\n",
      " [1.43766825]\n",
      " [7.04260971]\n",
      " [7.04581308]\n",
      " [2.18792106]\n",
      " [9.24867629]\n",
      " [4.42140755]\n",
      " [9.09315959]\n",
      " [0.59809223]\n",
      " [1.84287084]\n",
      " [0.47355279]\n",
      " [6.74880944]\n",
      " [5.9462478 ]\n",
      " [5.33310163]\n",
      " [0.43324063]\n",
      " [5.6143308 ]\n",
      " [3.29668446]\n",
      " [5.02966833]]\n"
     ]
    }
   ],
   "source": [
    "print(x0)\n",
    "print(z)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "id": "d12f9a6a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# specify learning rate\n",
    "lr = 0.1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "id": "5f2f8192",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "12068.560241773142\n",
      "(50, 1)\n"
     ]
    }
   ],
   "source": [
    "print(least_square(A,x0,b))\n",
    "print(least_square_grad(A,x0,b).shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 81,
   "id": "8e9ebd21",
   "metadata": {},
   "outputs": [],
   "source": [
    "# define a stopping creteria\n",
    "stop_criteria = 1e-8\n",
    "norm_grad = np.linalg.norm(least_square_grad(A,x0,b),2)\n",
    "err_ratio = np.linalg.norm(x0-z,2)/np.linalg.norm(z,2)\n",
    "obj_val = np.empty(0)\n",
    "estimation_err = np.empty(0)\n",
    "#loss=east_square(A,x0,b)\n",
    "#raining_loss = np.append(training_loss, loss)\n",
    "loop = 0\n",
    "#print(obj_val)\n",
    "#print(estimation_err)\n",
    "while norm_grad > stop_criteria:\n",
    "    \n",
    "    grad = least_square_grad(A,x0,b)\n",
    "    # gradient update\n",
    "    x0 = x0 - lr*grad\n",
    "    norm_grad = np.linalg.norm(grad,2)\n",
    "    loss_val = least_square(A,x0,b)\n",
    "    # record training loss and estimation error.\n",
    "    #print(obj_val)\n",
    "    obj_val = np.append(obj_val, loss_val)\n",
    "    \n",
    "    err_ratio = np.linalg.norm(x0-z,2)/np.linalg.norm(z,2)\n",
    "    estimation_err = np.append(estimation_err, err_ratio)\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 82,
   "id": "7eae098f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkQAAAHHCAYAAABeLEexAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABCZElEQVR4nO3deXQUZdr+8as7SXdCQhLWhLCLyL4JihEVkQwBUUH5qYwMy4gyYjKIuIzMKAouQXQAWZRxHAEdGJfxFRVfUTaDYlhkeVkNqAEiEKJCEtYs3fX7A1PQBCSRTleF/n7O6UO6qrrrroJjLu96niqHYRiGAAAAgpjT6gIAAACsRiACAABBj0AEAACCHoEIAAAEPQIRAAAIegQiAAAQ9AhEAAAg6BGIAABA0CMQAQCAoEcgAoLInDlz5HA4tGvXLqtL8VFSUqJHH31UDRs2lNPpVP/+/a0u6azsev5OVxVqBOwo1OoCAOD111/XCy+8oNGjR+vyyy9Xo0aNLK3nueeeU+vWrW0bzCrL/PnzlZubq9GjR1tdChBwdIgAWG7ZsmWqX7++pkyZosGDB6t79+6W1vPcc89pwYIFZZYPHjxYx48fV+PGjQNfVADMnz9fU6dOtboMwBJ0iABYLjc3V7GxsVaXcV4hISEKCQmxugwAlYAOEQC9/PLLatOmjdxutxISEpSSkqK8vDyfbXbu3KkBAwYoPj5e4eHhatCggQYOHKj8/Hxzm8WLF+uaa65RbGysoqKi1KJFC/31r38953537dolh8Oh5cuXa+vWrXI4HHI4HPr888/1+eefmz+f7TNz5swxlw0bNkxRUVHau3ev+vfvr6ioKNWpU0cPP/ywPB6Pz+e9Xq9eeukltWvXTuHh4apTp4569+6tr7/+WpLkcDh09OhRzZ0716xn2LBhks49Pqc85+/6669X27ZttW3bNvXo0UPVqlVT/fr1NWnSpHP/xZzG4XAoNTVV8+bNU4sWLRQeHq7OnTtrxYoV5fr8+Wq8/vrr9fHHH2v37t3mcTdp0qRc3w1cDOgQAUHuqaee0vjx45WUlKSRI0cqMzNTr7zyitauXauVK1cqLCxMRUVFSk5OVmFhof785z8rPj5ee/fu1cKFC5WXl6eYmBht3bpVN910k9q3b68JEybI7Xbr22+/1cqVK8+57zp16ujNN9/Us88+qyNHjigtLU2S1KpVK23fvr1Cx+HxeJScnKyuXbvqxRdf1JIlS/T3v/9dzZo108iRI83thg8frjlz5qhPnz665557VFJSoi+++EKrVq1Sly5d9Oabb+qee+7RlVdeqREjRkiSmjVrdkHnr9ShQ4fUu3dv3Xbbbbrjjjv03//+V3/5y1/Url079enT57zHmJ6errffflujRo2S2+3Wyy+/rN69e2vNmjVq27btBdX4t7/9Tfn5+frhhx80ZcoUSVJUVNR5awIuGgaAoDF79mxDkpGVlWUYhmHk5uYaLpfL6NWrl+HxeMztZsyYYUgyXn/9dcMwDGPDhg2GJOPdd98953dPmTLFkGT8+OOPFa6re/fuRps2bXyWLV++3JBkLF++3Gd5VlaWIcmYPXu2uWzo0KGGJGPChAk+23bq1Mno3Lmz+X7ZsmWGJGPUqFFlavB6vebPkZGRxtChQ8ts81vPX+kxSjLeeOMNc1lhYaERHx9vDBgwoOxJOYMkQ5Lx9ddfm8t2795thIeHG7feeqtfauzbt6/RuHHj89YCXIy4ZAYEsSVLlqioqEijR4+W03nqPwf33nuvoqOj9fHHH0uSYmJiJEmffvqpjh07dtbvKh0D9MEHH8jr9VZu4edw3333+by/9tpr9f3335vv33vvPTkcDj355JNlPutwOCq8v/Kev1JRUVH6wx/+YL53uVy68sorfWr8NYmJiercubP5vlGjRurXr58+/fTTMpcGf2uNQLAiEAFBbPfu3ZKkFi1a+Cx3uVy65JJLzPVNmzbVmDFj9Nprr6l27dpKTk7WzJkzfcYP3XnnnerWrZvuuecexcXFaeDAgXrnnXcCFo5KxwOdrkaNGjp06JD5/rvvvlNCQoJq1qzpl32W9/yVatCgQZngdWaNv6Z58+Zlll122WU6duyYfvzxR7/UCAQrAhGAcvn73/+uTZs26a9//auOHz+uUaNGqU2bNvrhhx8kSREREVqxYoWWLFmiwYMHa9OmTbrzzjv1u9/97pzdi19zro7Nub6rKsz+OleNhmEEuBIAZyIQAUGs9H46mZmZPsuLioqUlZVV5n477dq10+OPP64VK1boiy++0N69ezVr1ixzvdPpVM+ePTV58mRt27ZNzz77rJYtW6bly5dXuLYaNWpIUpnZWhfS0WjWrJn27dungwcP/up25b18VtHzd6F27txZZtmOHTtUrVq1Mt2x31Ljb7lsCFwsCERAEEtKSpLL5dK0adN8uhT/+te/lJ+fr759+0qSCgoKVFJS4vPZdu3ayel0qrCwUJLOGjI6duwoSeY2FdG4cWOFhISUmVb+8ssvV/i7Sg0YMECGYWj8+PFl1p1+/JGRkWWC2NmU9/z5S0ZGhtavX2++z87O1gcffKBevXqds/tUkRojIyN9LoMCwYRp90AQq1OnjsaOHavx48erd+/euuWWW5SZmamXX35ZV1xxhTkAeNmyZUpNTdXtt9+uyy67TCUlJXrzzTcVEhKiAQMGSJImTJigFStWqG/fvmrcuLFyc3P18ssvq0GDBrrmmmsqXFtMTIxuv/12TZ8+XQ6HQ82aNdPChQuVm5v7m4+3R48eGjx4sKZNm6adO3eqd+/e8nq9+uKLL9SjRw+lpqZKkjp37qwlS5Zo8uTJSkhIUNOmTdW1a9fffP78pW3btkpOTvaZdi/prAHvt9TYuXNnvf322xozZoyuuOIKRUVF6eabb/brMQC2ZeUUNwCBdeaU7FIzZswwWrZsaYSFhRlxcXHGyJEjjUOHDpnrv//+e+Puu+82mjVrZoSHhxs1a9Y0evToYSxZssTcZunSpUa/fv2MhIQEw+VyGQkJCcbvf/97Y8eOHeet62zT7g3DMH788UdjwIABRrVq1YwaNWoYf/rTn4wtW7acddp9ZGRkmc8/+eSTxpn/mSspKTFeeOEFo2XLlobL5TLq1Klj9OnTx1i3bp25zTfffGNcd911RkREhCHJnIL/W8/frx3j0KFDyzXVXZKRkpJi/Pvf/zaaN29uuN1uo1OnTmVuS3AhNR45csS46667jNjYWEMSU/ARVByGwWg+ALA7h8OhlJQUzZgxw+pSgIsSY4gAAEDQIxABAICgRyACAABBj1lmAFAFMNwTqFx0iAAAQNAjEAEAgKDHJbNy8Hq92rdvn6pXr86t7QEAqCIMw9Dhw4eVkJAgp/PXe0AEonLYt2+fGjZsaHUZAADgN8jOzlaDBg1+dRsCUTlUr15d0skTGh0dbXE1AACgPAoKCtSwYUPz9/ivIRCVQ+llsujoaAIRAABVTHmGuzCoGgAABD0CEQAACHoEIgAAEPQIRAAAIOgRiAAAQNAjEAEAgKBHIAIAAEGPQAQAAIIegQgAAAQ9AhEAAAh6BCIAABD0CEQAACDoEYgs5PEa2p9/XHt+PmZ1KQAABDWedm+h3MMnlJi2TKFOh7597karywEAIGjRIbKQK+Tk6S/xGvJ4DYurAQAgeBGILOQOCzF/LirxWlgJAADBjUBkodIOkUQgAgDASgQiC4WFOORwnPy50OOxthgAAIIYgchCDofD7BIVFtMhAgDAKgQii7lCT/4VFHkIRAAAWIVAZDF36MmB1XSIAACwDoHIYm46RAAAWI5AZLHSQFRYzKBqAACsQiCyGGOIAACwHoHIYqc6RAQiAACsQiCyGB0iAACsRyCymDnLrIQxRAAAWIVAZDGzQ8SjOwAAsAyByGLmGCICEQAAliEQWYwOEQAA1iMQWYwOEQAA1iMQWcxFIAIAwHIEIouVzjLjkhkAANYhEFnsVIeIafcAAFiFQGQxVwiDqgEAsJqlgWjFihW6+eablZCQIIfDoQULFvisNwxD48aNU7169RQREaGkpCTt3LnTZ5uDBw9q0KBBio6OVmxsrIYPH64jR474bLNp0yZde+21Cg8PV8OGDTVp0qTKPrRyc4cxhggAAKtZGoiOHj2qDh06aObMmWddP2nSJE2bNk2zZs3S6tWrFRkZqeTkZJ04ccLcZtCgQdq6dasWL16shQsXasWKFRoxYoS5vqCgQL169VLjxo21bt06vfDCC3rqqaf06quvVvrxlQcdIgAArBdq5c779OmjPn36nHWdYRiaOnWqHn/8cfXr10+S9MYbbyguLk4LFizQwIEDtX37di1atEhr165Vly5dJEnTp0/XjTfeqBdffFEJCQmaN2+eioqK9Prrr8vlcqlNmzbauHGjJk+e7BOcrOIO49EdAABYzbZjiLKyspSTk6OkpCRzWUxMjLp27aqMjAxJUkZGhmJjY80wJElJSUlyOp1avXq1uc11110nl8tlbpOcnKzMzEwdOnTorPsuLCxUQUGBz6uyuOkQAQBgOdsGopycHElSXFycz/K4uDhzXU5OjurWreuzPjQ0VDVr1vTZ5mzfcfo+zpSWlqaYmBjz1bBhwws/oHNgDBEAANazbSCy0tixY5Wfn2++srOzK21fjCECAMB6tg1E8fHxkqQDBw74LD9w4IC5Lj4+Xrm5uT7rS0pKdPDgQZ9tzvYdp+/jTG63W9HR0T6vykKHCAAA69k2EDVt2lTx8fFaunSpuaygoECrV69WYmKiJCkxMVF5eXlat26duc2yZcvk9XrVtWtXc5sVK1aouLjY3Gbx4sVq0aKFatSoEaCjOTdXCHeqBgDAapYGoiNHjmjjxo3auHGjpJMDqTdu3Kg9e/bI4XBo9OjReuaZZ/Thhx9q8+bNGjJkiBISEtS/f39JUqtWrdS7d2/de++9WrNmjVauXKnU1FQNHDhQCQkJkqS77rpLLpdLw4cP19atW/X222/rpZde0pgxYyw6al+lHaIiD4EIAACrWDrt/uuvv1aPHj3M96UhZejQoZozZ44effRRHT16VCNGjFBeXp6uueYaLVq0SOHh4eZn5s2bp9TUVPXs2VNOp1MDBgzQtGnTzPUxMTH67LPPlJKSos6dO6t27doaN26cLabcS6fGEBUWM+0eAACrOAzDMKwuwu4KCgoUExOj/Px8v48n+ianQL2nfqHaUS59/fjv/PrdAAAEs4r8/rbtGKJgcapDxCUzAACsQiCymHmnasYQAQBgGQKRxU6/DxFXLwEAsAaByGKu0FN/Bcw0AwDAGgQii7lPC0TcnBEAAGsQiCxWeslM4uaMAABYhUBkMafTcWqmGYEIAABLEIhsoHQcER0iAACsQSCygdJxRIUl3K0aAAArEIhsgA4RAADWIhDZgJtABACApQhENuAKZVA1AABWIhDZgDv05OM76BABAGANApENuBhUDQCApQhENuDmkhkAAJYiENkAY4gAALAWgcgGmGUGAIC1CEQ24PplUDUdIgAArEEgsoHSZ5nRIQIAwBoEIhtwhzHLDAAAKxGIbIAOEQAA1iIQ2cCpDhGBCAAAKxCIbMBNhwgAAEsRiGzAHVY6y4wxRAAAWIFAZAOMIQIAwFoEIhsoHUNU5CEQAQBgBQKRDZR2iAqLCUQAAFiBQGQDdIgAALAWgcgGXCG/DKqmQwQAgCUIRDZQ+nDXQjpEAABYgkBkA67SQFTMtHsAAKxAILKB0g4RY4gAALAGgcgGTnWICEQAAFiBQGQDLjpEAABYikBkA+7Q0llmjCECAMAKBCIbYAwRAADWIhDZgDntvsQrwzAsrgYAgOBDILKB0jFEhiGVeAlEAAAEGoHIBkrHEEk88R4AACsQiGygtEMknbxsBgAAAotAZAMhTodCnQ5JdIgAALACgcgmzJszljD1HgCAQCMQ2YQ59Z4OEQAAAUcgsgnXaVPvAQBAYBGIbMK8WzWBCACAgCMQ2QRjiAAAsA6ByCYYQwQAgHUIRDbBGCIAAKxDILIJVwgdIgAArEIgsgl3GIOqAQCwCoHIJugQAQBgHQKRTbjDSgMRs8wAAAg0ApFNuEMYVA0AgFUIRDZxqkNEIAIAINAIRDbhokMEAIBlCEQ2UTrLrMhDIAIAINAIRDZhdoiKGVQNAECgEYhswnx0Bx0iAAACjkBkE+ajO4oJRAAABBqByCZKO0SFdIgAAAg4ApFNuEJ/eXQHHSIAAALO1oHI4/HoiSeeUNOmTRUREaFmzZrp6aeflmEY5jaGYWjcuHGqV6+eIiIilJSUpJ07d/p8z8GDBzVo0CBFR0crNjZWw4cP15EjRwJ9OL+KMUQAAFjH1oHo+eef1yuvvKIZM2Zo+/btev755zVp0iRNnz7d3GbSpEmaNm2aZs2apdWrVysyMlLJyck6ceKEuc2gQYO0detWLV68WAsXLtSKFSs0YsQIKw7pnE6NIWKWGQAAgRZqdQG/5quvvlK/fv3Ut29fSVKTJk30n//8R2vWrJF0sjs0depUPf744+rXr58k6Y033lBcXJwWLFiggQMHavv27Vq0aJHWrl2rLl26SJKmT5+uG2+8US+++KISEhKsObgzuOgQAQBgGVt3iK6++motXbpUO3bskCT93//9n7788kv16dNHkpSVlaWcnBwlJSWZn4mJiVHXrl2VkZEhScrIyFBsbKwZhiQpKSlJTqdTq1evPut+CwsLVVBQ4POqbOYlM+5UDQBAwNm6Q/TYY4+poKBALVu2VEhIiDwej5599lkNGjRIkpSTkyNJiouL8/lcXFycuS4nJ0d169b1WR8aGqqaNWua25wpLS1N48eP9/fh/CrzkhmBCACAgLN1h+idd97RvHnzNH/+fK1fv15z587Viy++qLlz51bqfseOHav8/HzzlZ2dXan7kyT3L7PM6BABABB4tu4QPfLII3rsscc0cOBASVK7du20e/dupaWlaejQoYqPj5ckHThwQPXq1TM/d+DAAXXs2FGSFB8fr9zcXJ/vLSkp0cGDB83Pn8ntdsvtdlfCEZ2beR+iEgZVAwAQaLbuEB07dkxOp2+JISEh8npPdlGaNm2q+Ph4LV261FxfUFCg1atXKzExUZKUmJiovLw8rVu3ztxm2bJl8nq96tq1awCOonwYQwQAgHVs3SG6+eab9eyzz6pRo0Zq06aNNmzYoMmTJ+vuu++WJDkcDo0ePVrPPPOMmjdvrqZNm+qJJ55QQkKC+vfvL0lq1aqVevfurXvvvVezZs1ScXGxUlNTNXDgQNvMMJMYQwQAgJVsHYimT5+uJ554Qvfff79yc3OVkJCgP/3pTxo3bpy5zaOPPqqjR49qxIgRysvL0zXXXKNFixYpPDzc3GbevHlKTU1Vz5495XQ6NWDAAE2bNs2KQzonxhABAGAdh3H6bZ9xVgUFBYqJiVF+fr6io6MrZR8Hjxbp8qcXS5K+e+5GhTgdlbIfAACCRUV+f9t6DFEwKR1DJNElAgAg0AhENuE6LRAx0wwAgMAiENlEqNOh0qtkdIgAAAgsApFNOBwOZpoBAGARApGNlM40IxABABBYBCIbcXFzRgAALEEgshFXCI/vAADACgQiG3GH0SECAMAKBCIbOdUhIhABABBIBCIbcYfx+A4AAKxAILIRNx0iAAAsQSCyEXMMkYdB1QAABBKByEbMMUTFdIgAAAgkApGNnOoQEYgAAAgkApGN0CECAMAaBCIbKX10Bx0iAAACi0BkI+bDXYsZVA0AQCARiGzEXRqI6BABABBQBCIb4eGuAABYg0BkI6VjiLgxIwAAgUUgshE6RAAAWINAZCPmGCICEQAAARX6Wz60dOlSLV26VLm5ufJ6fX95v/76634pLBid6hAxywwAgECqcCAaP368JkyYoC5duqhevXpyOByVUVdQctEhAgDAEhUORLNmzdKcOXM0ePDgyqgnqLkZQwQAgCUqPIaoqKhIV199dWXUEvQYQwQAgDUqHIjuuecezZ8/vzJqCXrmozsIRAAABFSFL5mdOHFCr776qpYsWaL27dsrLCzMZ/3kyZP9VlywOTWGiEHVAAAEUoUD0aZNm9SxY0dJ0pYtW3zWMcD6wjCGCAAAa1Q4EC1fvrwy6oCYZQYAgFUu6MaMP/zwg3744Qd/1RL0GEMEAIA1KhyIvF6vJkyYoJiYGDVu3FiNGzdWbGysnn766TI3aUTF8OgOAACsUeFLZn/729/0r3/9SxMnTlS3bt0kSV9++aWeeuopnThxQs8++6zfiwwWTLsHAMAaFQ5Ec+fO1WuvvaZbbrnFXNa+fXvVr19f999/P4HoApgdIo9XhmEwSB0AgACp8CWzgwcPqmXLlmWWt2zZUgcPHvRLUcGqtEMk0SUCACCQKhyIOnTooBkzZpRZPmPGDHXo0MEvRQUr12mBqMhDIAIAIFAqfMls0qRJ6tu3r5YsWaLExERJUkZGhrKzs/W///u/fi8wmLhCTusQFXulcAuLAQAgiFS4Q9S9e3ft2LFDt956q/Ly8pSXl6fbbrtNmZmZuvbaayujxqDhcDh8xhEBAIDAqHCHSJISEhIYPF1J3CFOFZV4VVjM4zsAAAiUcgWiTZs2qW3btnI6ndq0adOvbtu+fXu/FBas3GFOHS6kQwQAQCCVKxB17NhROTk5qlu3rjp27CiHwyHDMMps53A45PHQ2bgQpeOICosJRAAABEq5AlFWVpbq1Klj/ozK4w775fEddIgAAAiYcgWixo0bmz/v3r1bV199tUJDfT9aUlKir776ymdbVBwdIgAAAq/Cs8x69Ohx1hsw5ufnq0ePHn4pKpi5w0pnmXHpEQCAQKlwIDrXIyV+/vlnRUZG+qWoYEaHCACAwCv3tPvbbrtN0smB08OGDZPb7TbXeTwebdq0SVdffbX/KwwypzpEBCIAAAKl3IEoJiZG0skOUfXq1RUREWGuc7lcuuqqq3Tvvff6v8IgY3aIeJYZAAABU+5ANHv2bElSkyZN9Mgjj6hatWqVVlQwc4eenGVGIAIAIHAqPIZoyJAh2rt3b5nlO3fu1K5du/xRU1AzH91BIAIAIGAqHIiGDRumr776qszy1atXa9iwYf6oKai5Q0svmTHLDACAQKlwINqwYYO6detWZvlVV12ljRs3+qOmoEaHCACAwKtwIHI4HDp8+HCZ5fn5+Ty2ww8YQwQAQOBVOBBdd911SktL8wk/Ho9HaWlpuuaaa/xaXDCiQwQAQOCVe5ZZqeeff17XXXedWrRooWuvvVaS9MUXX6igoEDLli3ze4HBxsUYIgAAAq7CHaLWrVtr06ZNuuOOO5Sbm6vDhw9ryJAh+uabb9S2bdvKqDGouOkQAQAQcBXuEElSQkKCnnvuOX/XAp0+y4xABABAoPymQJSXl6c1a9YoNzdXXq/vL+4hQ4b4pbBgRYcIAIDAq3Ag+uijjzRo0CAdOXJE0dHRPg96dTgcBKIL5KJDBABAwFV4DNFDDz2ku+++W0eOHFFeXp4OHTpkvg4ePFgZNQaV0mn3dIgAAAicCgeivXv3atSoUTzLrJIw7R4AgMCrcCBKTk7W119/XRm1nNXevXv1hz/8QbVq1VJERITatWvns3/DMDRu3DjVq1dPERERSkpK0s6dO32+4+DBgxo0aJCio6MVGxur4cOH68iRIwE7horg0R0AAARehccQ9e3bV4888oi2bdumdu3aKSwszGf9Lbfc4rfiDh06pG7duqlHjx765JNPVKdOHe3cuVM1atQwt5k0aZKmTZumuXPnqmnTpnriiSeUnJysbdu2KTw8XJI0aNAg7d+/X4sXL1ZxcbH++Mc/asSIEZo/f77favUXxhABABB4DsMwjIp8wOk8d1PJ4XD49fEdjz32mFauXKkvvvjirOsNw1BCQoIeeughPfzww5JOPkIkLi5Oc+bM0cCBA7V9+3a1bt1aa9euVZcuXSRJixYt0o033qgffvhBCQkJ562joKBAMTExys/PV3R0tN+O72zWZB3UHf/I0CW1I7Xs4esrdV8AAFzMKvL7u8KXzLxe7zlf/n6W2YcffqguXbro9ttvV926ddWpUyf985//NNdnZWUpJydHSUlJ5rKYmBh17dpVGRkZkqSMjAzFxsaaYUiSkpKS5HQ6tXr16rPut7CwUAUFBT6vQKFDBABA4FU4EAXS999/r1deeUXNmzfXp59+qpEjR2rUqFGaO3euJCknJ0eSFBcX5/O5uLg4c11OTo7q1q3rsz40NFQ1a9Y0tzlTWlqaYmJizFfDhg39fWjnxI0ZAQAIvAqPIZowYcKvrh83btxvLuZMXq9XXbp0Me+K3alTJ23ZskWzZs3S0KFD/bafM40dO1Zjxowx3xcUFAQsFJ2aZcagagAAAqXCgej999/3eV9cXKysrCyFhoaqWbNmfg1E9erVU+vWrX2WtWrVSu+9954kKT4+XpJ04MAB1atXz9zmwIED6tixo7lNbm6uz3eUlJTo4MGD5ufP5Ha75Xa7/XUYFUKHCACAwKtwINqwYUOZZQUFBRo2bJhuvfVWvxRVqlu3bsrMzPRZtmPHDjVu3FiS1LRpU8XHx2vp0qVmACooKNDq1as1cuRISVJiYqLy8vK0bt06de7cWZK0bNkyeb1ede3a1a/1+oPZIfJ4ZRiGz53AAQBA5fDLGKLo6GiNHz9eTzzxhD++zvTggw9q1apVeu655/Ttt99q/vz5evXVV5WSkiLp5Ky20aNH65lnntGHH36ozZs3a8iQIUpISFD//v0lnewo9e7dW/fee6/WrFmjlStXKjU1VQMHDizXDLNAc4ecvFO1YUjFngpNAAQAAL/Rb3q469nk5+crPz/fX18nSbriiiv0/vvva+zYsZowYYKaNm2qqVOnatCgQeY2jz76qI4ePaoRI0YoLy9P11xzjRYtWmTeg0iS5s2bp9TUVPXs2VNOp1MDBgzQtGnT/Fqrv7jDTmXUIo/X7BgBAIDKU+H7EJ0ZJAzD0P79+/Xmm2+qe/futrzZ4YUK5H2IvF5Dl/z1fyVJ6x5PUq0oa8YyAQBQ1VXk93eFO0RTpkzxee90OlWnTh0NHTpUY8eOrejX4QxOp0NhIQ4VewwVeRhYDQBAIJQrEG3atElt27aV0+lUVlZWZdcU9FwhThV7PDzgFQCAACnXAJVOnTrpp59+kiRdcskl+vnnnyu1qGDnDjs5sJqp9wAABEa5AlFsbKzZGdq1a5e8Xn5RVyZXSOnNGTnPAAAEQrkumQ0YMEDdu3dXvXr15HA41KVLF4X8Mj38TN9//71fCwxGpTPNCrlbNQAAAVGuQPTqq6/qtttu07fffqtRo0bp3nvvVfXq1Su7tqBV2iHikhkAAIFR7llmvXv3liStW7dODzzwAIGoEp3qEBGIAAAIhApPu589e3Zl1IHTMIYIAIDA4jbINuQOZZYZAACBRCCyIfMBrwQiAAACgkBkQ+5QZpkBABBIBCIbokMEAEBgVTgQzZ07Vx9//LH5/tFHH1VsbKyuvvpq7d6926/FBStXKLPMAAAIpAoHoueee04RERGSpIyMDM2cOVOTJk1S7dq19eCDD/q9wGBUOqiaDhEAAIFR4Wn32dnZuvTSSyVJCxYs0IABAzRixAh169ZN119/vb/rC0puLpkBABBQFe4QRUVFmQ93/eyzz/S73/1OkhQeHq7jx4/7t7ogxaBqAAACq8Idot/97ne655571KlTJ+3YsUM33nijJGnr1q1q0qSJv+sLSgyqBgAgsCrcIZo5c6YSExP1448/6r333lOtWrUknXykx+9//3u/FxiM3AyqBgAgoCrcIYqNjdWMGTPKLB8/frxfCgIdIgAAAq3CHaJFixbpyy+/NN/PnDlTHTt21F133aVDhw75tbhgxaM7AAAIrAoHokceeUQFBQWSpM2bN+uhhx7SjTfeqKysLI0ZM8bvBQYj7kMEAEBgVfiSWVZWllq3bi1Jeu+993TTTTfpueee0/r1680B1rgwzDIDACCwKtwhcrlcOnbsmCRpyZIl6tWrlySpZs2aZucIF4YxRAAABFaFO0TXXHONxowZo27dumnNmjV6++23JUk7duxQgwYN/F5gMGIMEQAAgVXhDtGMGTMUGhqq//73v3rllVdUv359SdInn3yi3r17+73AYESHCACAwKpwh6hRo0ZauHBhmeVTpkzxS0FgDBEAAIFW4UAkSR6PRwsWLND27dslSW3atNEtt9yikJAQvxYXrMwOkYcOEQAAgVDhQPTtt9/qxhtv1N69e9WiRQtJUlpamho2bKiPP/5YzZo183uRwcYVwiUzAAACqcJjiEaNGqVmzZopOztb69ev1/r167Vnzx41bdpUo0aNqowag054GPchAgAgkCrcIUpPT9eqVatUs2ZNc1mtWrU0ceJEdevWza/FBSvXL5ce6RABABAYFe4Qud1uHT58uMzyI0eOyOVy+aWoYOemQwQAQEBVOBDddNNNGjFihFavXi3DMGQYhlatWqX77rtPt9xyS2XUGHRKxxB5vIZKGFgNAEClq3AgmjZtmpo1a6bExESFh4crPDxc3bp106WXXqqXXnqpMmoMOqUdIomZZgAABEKFxxDFxsbqgw8+0M6dO/XNN99Iklq1aqVLL73U78UFq9IOkXRyHFE1rkQCAFCpftN9iCSpefPmat68uT9rwS9CQ5wKcTrk8RqMIwIAIADKFYjGjBlT7i+cPHnyby4Gp7hCnDru9TDTDACAAChXINqwYUO5vszhcFxQMTjFHebU8WIPj+8AACAAyhWIli9fXtl14Ayl44i4ZAYAQOWr8CwzBAb3IgIAIHAIRDbF88wAAAgcApFNuUN5fAcAAIFCILIpVyiXzAAACBQCkU2VBiI6RAAAVD4CkU25zQ4R0+4BAKhsBCKbctMhAgAgYAhENlU6qJoxRAAAVD4CkU0xhggAgMAhENkUY4gAAAgcApFN0SECACBwCEQ25eY+RAAABAyByKa4MSMAAIFDILIpZpkBABA4BCKbYgwRAACBQyCyKfPGjB4CEQAAlY1AZFPmGKJipt0DAFDZCEQ2VTqGiA4RAACVj0BkU6c6RAQiAAAqG4HIplwhjCECACBQCEQ25Q7j0R0AAAQKgcim3CFMuwcAIFAIRDZ1qkNEIAIAoLJVqUA0ceJEORwOjR492lx24sQJpaSkqFatWoqKitKAAQN04MABn8/t2bNHffv2VbVq1VS3bl098sgjKikpCXD1FeMK+WWWGYEIAIBKV2UC0dq1a/WPf/xD7du391n+4IMP6qOPPtK7776r9PR07du3T7fddpu53uPxqG/fvioqKtJXX32luXPnas6cORo3blygD6FC6BABABA4VSIQHTlyRIMGDdI///lP1ahRw1yen5+vf/3rX5o8ebJuuOEGde7cWbNnz9ZXX32lVatWSZI+++wzbdu2Tf/+97/VsWNH9enTR08//bRmzpypoqIiqw7pvFyMIQIAIGCqRCBKSUlR3759lZSU5LN83bp1Ki4u9lnesmVLNWrUSBkZGZKkjIwMtWvXTnFxceY2ycnJKigo0NatW8+6v8LCQhUUFPi8Ao1ZZgAABE6o1QWcz1tvvaX169dr7dq1Zdbl5OTI5XIpNjbWZ3lcXJxycnLMbU4PQ6XrS9edTVpamsaPH++H6n+70g5RsceQ12vI6XRYWg8AABczW3eIsrOz9cADD2jevHkKDw8P2H7Hjh2r/Px885WdnR2wfZdyh4WYP3NzRgAAKpetA9G6deuUm5uryy+/XKGhoQoNDVV6erqmTZum0NBQxcXFqaioSHl5eT6fO3DggOLj4yVJ8fHxZWadlb4v3eZMbrdb0dHRPq9AK+0QSQysBgCgstk6EPXs2VObN2/Wxo0bzVeXLl00aNAg8+ewsDAtXbrU/ExmZqb27NmjxMRESVJiYqI2b96s3Nxcc5vFixcrOjparVu3DvgxlVdYiEOOX66SMbAaAIDKZesxRNWrV1fbtm19lkVGRqpWrVrm8uHDh2vMmDGqWbOmoqOj9ec//1mJiYm66qqrJEm9evVS69atNXjwYE2aNEk5OTl6/PHHlZKSIrfbHfBjKi+HwyFXiFOFJV4GVgMAUMlsHYjKY8qUKXI6nRowYIAKCwuVnJysl19+2VwfEhKihQsXauTIkUpMTFRkZKSGDh2qCRMmWFh1+bhDTwYiOkQAAFQuh2EYhtVF2F1BQYFiYmKUn58f0PFEXZ5Zop+OFOqTB65Vq3qBH8cEAEBVVpHf37YeQxTs3KHcnBEAgEAgENlYaSBilhkAAJWLQGRjLjpEAAAEBIHIxk51iJhlBgBAZSIQ2RgdIgAAAoNAZGPu0JOP72AMEQAAlYtAZGN0iAAACAwCkY0xhggAgMAgENmYi2n3AAAEBIHIxswbM3oIRAAAVCYCkY2ZHaJiAhEAAJWJQGRjpbPM6BABAFC5CEQ2RocIAIDAIBDZ2KkxRMwyAwCgMhGIbIwOEQAAgUEgsjHGEAEAEBgEIhujQwQAQGAQiGzMHcJ9iAAACAQCkY25w3h0BwAAgUAgsjFXCA93BQAgEAhENnaqQ0QgAgCgMhGIbMwV8sssMwIRAACVikBkY6UdIgIRAACVi0BkY6VjiLhkBgBA5SIQ2RhjiAAACAwCkY2d6hAx7R4AgMpEILIxdxiDqgEACAQCkY2dPobIMAyLqwEA4OJFILKx0jFEklTsIRABAFBZCEQ2VtohkhhHBABAZSIQ2Zg79NRfD+OIAACoPAQiG3M4HNyLCACAACAQ2ZwrlLtVAwBQ2QhENld62YwOEQAAlYdAZHO1o9ySpL15xyyuBACAixeByOZaJ0RLkrbuLbC4EgAALl4EIptr80sg2rIv3+JKAAC4eBGIbK5NQowkaes+OkQAAFQWApHNlV4y++HQceUfK7a4GgAALk4EIpuLiQhTo5rVJElbuWwGAEClIBBVAaXjiLhsBgBA5SAQVQFt658cR8TAagAAKgeBqApoTYcIAIBKRSCqAkovmX334xEdKyqxuBoAAC4+BKIqoG71cNWt7pZhSNv3H7a6HAAALjoEoiqitEu0jXFEAAD4HYGoijAHVvMIDwAA/I5AVEWYU+/30yECAMDfCERVROkjPDJzDquoxGtxNQAAXFwIRFVEgxoRig4PVbHH0M5cBlYDAOBPBKIqwuFw8KBXAAAqCYGoCmlb/5dxRHsZRwQAgD8RiKoQOkQAAFQOAlEVYt6LaH+BPF7D4moAALh4EIiqkEvqRCk8zKljRR7t+vmo1eUAAHDRIBBVISFOh1rVO9kl2sI4IgAA/IZAVMW0/WUc0TbGEQEA4DcEoirGvGM1gQgAAL8hEFUxpTPNtuzLl2EwsBoAAH8gEFUxl8VHKdTpUN6xYu3LP2F1OQAAXBRsHYjS0tJ0xRVXqHr16qpbt6769++vzMxMn21OnDihlJQU1apVS1FRURowYIAOHDjgs82ePXvUt29fVatWTXXr1tUjjzyikpKSQB6K37hDQ9Q8rrokBlYDAOAvtg5E6enpSklJ0apVq7R48WIVFxerV69eOnr01JTzBx98UB999JHeffddpaena9++fbrtttvM9R6PR3379lVRUZG++uorzZ07V3PmzNG4ceOsOCS/aMs4IgAA/MphVKGBKD/++KPq1q2r9PR0XXfddcrPz1edOnU0f/58/b//9/8kSd98841atWqljIwMXXXVVfrkk0900003ad++fYqLi5MkzZo1S3/5y1/0448/yuVynXe/BQUFiomJUX5+vqKjoyv1GMtjzsosPfXRNvVsWVf/GnaF1eUAAGBLFfn9besO0Zny809eIqpZs6Ykad26dSouLlZSUpK5TcuWLdWoUSNlZGRIkjIyMtSuXTszDElScnKyCgoKtHXr1rPup7CwUAUFBT4vO2lTn0d4AADgT1UmEHm9Xo0ePVrdunVT27ZtJUk5OTlyuVyKjY312TYuLk45OTnmNqeHodL1pevOJi0tTTExMearYcOGfj6aC9OqXrQcDimn4IR+OlJodTkAAFR5VSYQpaSkaMuWLXrrrbcqfV9jx45Vfn6++crOzq70fVZElDtUTWtFSqJLBACAP1SJQJSamqqFCxdq+fLlatCggbk8Pj5eRUVFysvL89n+wIEDio+PN7c5c9ZZ6fvSbc7kdrsVHR3t87Kb1ubAamaaAQBwoWwdiAzDUGpqqt5//30tW7ZMTZs29VnfuXNnhYWFaenSpeayzMxM7dmzR4mJiZKkxMREbd68Wbm5ueY2ixcvVnR0tFq3bh2YA6kEbUvHEe2lQwQAwIUKtbqAX5OSkqL58+frgw8+UPXq1c0xPzExMYqIiFBMTIyGDx+uMWPGqGbNmoqOjtaf//xnJSYm6qqrrpIk9erVS61bt9bgwYM1adIk5eTk6PHHH1dKSorcbreVh3dB2tAhAgDAb2wdiF555RVJ0vXXX++zfPbs2Ro2bJgkacqUKXI6nRowYIAKCwuVnJysl19+2dw2JCRECxcu1MiRI5WYmKjIyEgNHTpUEyZMCNRhVIrSR3js+vmYDp8oVvXwMIsrAgCg6qpS9yGyit3uQ1Tq6rSl2pd/Qm+PuEpdL6lldTkAANjKRXsfIvhqncD9iAAA8AcCURXWtv7JtLuFcUQAAFwQAlEVVjqOaBsdIgAALgiBqAor7RDtzD2iE8Uei6sBAKDqIhBVYfHR4aoZ6ZLHaygz57DV5QAAUGURiKowh8Oh9g1OXjabumSHij1eiysCAKBqIhBVcaOTLlN4mFPLM3/U2P/ZLO6iAABAxRGIqriODWM1867LFeJ06L/rftDfP9thdUkAAFQ5BKKLQM9WcXru1raSpBnLv9UbGbusLQgAgCqGQHSRuPOKRnrod5dJkp78cKv+d/N+iysCAKDqIBBdRFJvuFR/uKqRDEMa/dZGrfr+Z6tLAgCgSiAQXUQcDofG39JWvdvEq8jj1b1vfK1vcrhpIwAA50MgusiEOB2aOrCjrmxSU4dPlGjo62u0N++41WUBAGBrBKKLUHhYiP45pIsui4vSgYJCDfnXau08cJgp+QAAnIPD4LfkeRUUFCgmJkb5+fmKjo62upxy259/XLe9/JX255+QJNWPjdD1Leqo+2V11O3S2op0h1pcIQAAlaciv78JROVQVQORJH334xFN+GibMr77WUWn3ck6LMShK5rUVPfL6uj6FnV1ad0ohTgdFlYKAIB/EYj8rCoHolLHikq06vuf9Xnmj/o880ftOXjMZ73DIcVEhKlmpEs1q7lUI9KlWpEn/6xZzaWo8FCFhTgVFuKQK8R58ufQU+9DQ5xyOiSnwyHHL3+efJ0c7O1wSA6d/Fk6+XPpfh2/vHNcYB670M8DAKwT4nSoXkyEX7+TQORnF0MgOlPWT0f1eWau0nf8qIzvflZhCc9BAwBYp251t9b8Lcmv31mR398MIglSTWtHqmntpvpjt6Yq8XiVd7xYh44W6eejRT5/HjxWpINHi3SsyKNij/fkq8RQocer4hLvqWWek7naaxi/vCTjlz+9hiGv15CZvI1Tf5Tm8TNT+Zkx3SizRcUR/QHAvtxh1s7zIhBBoSFO1Y5yq3aUW82tLgYAAAsw7R4AAAQ9AhEAAAh6BCIAABD0CEQAACDoEYgAAEDQIxABAICgRyACAABBj0AEAACCHoEIAAAEPQIRAAAIegQiAAAQ9AhEAAAg6BGIAABA0CMQAQCAoBdqdQFVgWEYkqSCggKLKwEAAOVV+nu79Pf4ryEQlcPhw4clSQ0bNrS4EgAAUFGHDx9WTEzMr27jMMoTm4Kc1+vVvn37VL16dTkcDr9+d0FBgRo2bKjs7GxFR0f79buDCefRPziP/sF59A/Oo38E83k0DEOHDx9WQkKCnM5fHyVEh6gcnE6nGjRoUKn7iI6ODrp/qJWB8+gfnEf/4Dz6B+fRP4L1PJ6vM1SKQdUAACDoEYgAAEDQIxBZzO1268knn5Tb7ba6lCqN8+gfnEf/4Dz6B+fRPziP5cOgagAAEPToEAEAgKBHIAIAAEGPQAQAAIIegQgAAAQ9ApGFZs6cqSZNmig8PFxdu3bVmjVrrC7J9lasWKGbb75ZCQkJcjgcWrBggc96wzA0btw41atXTxEREUpKStLOnTutKdam0tLSdMUVV6h69eqqW7eu+vfvr8zMTJ9tTpw4oZSUFNWqVUtRUVEaMGCADhw4YFHF9vTKK6+offv25s3uEhMT9cknn5jrOYe/zcSJE+VwODR69GhzGefy/J566ik5HA6fV8uWLc31nMPzIxBZ5O2339aYMWP05JNPav369erQoYOSk5OVm5trdWm2dvToUXXo0EEzZ8486/pJkyZp2rRpmjVrllavXq3IyEglJyfrxIkTAa7UvtLT05WSkqJVq1Zp8eLFKi4uVq9evXT06FFzmwcffFAfffSR3n33XaWnp2vfvn267bbbLKzafho0aKCJEydq3bp1+vrrr3XDDTeoX79+2rp1qyTO4W+xdu1a/eMf/1D79u19lnMuy6dNmzbav3+/+fryyy/NdZzDcjBgiSuvvNJISUkx33s8HiMhIcFIS0uzsKqqRZLx/vvvm++9Xq8RHx9vvPDCC+ayvLw8w+12G//5z38sqLBqyM3NNSQZ6enphmGcPGdhYWHGu+++a26zfft2Q5KRkZFhVZlVQo0aNYzXXnuNc/gbHD582GjevLmxePFio3v37sYDDzxgGAb/HsvrySefNDp06HDWdZzD8qFDZIGioiKtW7dOSUlJ5jKn06mkpCRlZGRYWFnVlpWVpZycHJ/zGhMTo65du3Jef0V+fr4kqWbNmpKkdevWqbi42Oc8tmzZUo0aNeI8noPH49Fbb72lo0ePKjExkXP4G6SkpKhv374+50zi32NF7Ny5UwkJCbrkkks0aNAg7dmzRxLnsLx4uKsFfvrpJ3k8HsXFxfksj4uL0zfffGNRVVVfTk6OJJ31vJaugy+v16vRo0erW7duatu2raST59Hlcik2NtZnW85jWZs3b1ZiYqJOnDihqKgovf/++2rdurU2btzIOayAt956S+vXr9fatWvLrOPfY/l07dpVc+bMUYsWLbR//36NHz9e1157rbZs2cI5LCcCERDEUlJStGXLFp+xBii/Fi1aaOPGjcrPz9d///tfDR06VOnp6VaXVaVkZ2frgQce0OLFixUeHm51OVVWnz59zJ/bt2+vrl27qnHjxnrnnXcUERFhYWVVB5fMLFC7dm2FhISUGeF/4MABxcfHW1RV1Vd67jiv5ZOamqqFCxdq+fLlatCggbk8Pj5eRUVFysvL89me81iWy+XSpZdeqs6dOystLU0dOnTQSy+9xDmsgHXr1ik3N1eXX365QkNDFRoaqvT0dE2bNk2hoaGKi4vjXP4GsbGxuuyyy/Ttt9/y77GcCEQWcLlc6ty5s5YuXWou83q9Wrp0qRITEy2srGpr2rSp4uPjfc5rQUGBVq9ezXk9jWEYSk1N1fvvv69ly5apadOmPus7d+6ssLAwn/OYmZmpPXv2cB7Pw+v1qrCwkHNYAT179tTmzZu1ceNG89WlSxcNGjTI/JlzWXFHjhzRd999p3r16vHvsbysHtUdrN566y3D7XYbc+bMMbZt22aMGDHCiI2NNXJycqwuzdYOHz5sbNiwwdiwYYMhyZg8ebKxYcMGY/fu3YZhGMbEiRON2NhY44MPPjA2bdpk9OvXz2jatKlx/Phxiyu3j5EjRxoxMTHG559/buzfv998HTt2zNzmvvvuMxo1amQsW7bM+Prrr43ExEQjMTHRwqrt57HHHjPS09ONrKwsY9OmTcZjjz1mOBwO47PPPjMMg3N4IU6fZWYYnMvyeOihh4zPP//cyMrKMlauXGkkJSUZtWvXNnJzcw3D4ByWB4HIQtOnTzcaNWpkuFwu48orrzRWrVpldUm2t3z5ckNSmdfQoUMNwzg59f6JJ54w4uLiDLfbbfTs2dPIzMy0tmibOdv5k2TMnj3b3Ob48ePG/fffb9SoUcOoVq2aceuttxr79++3rmgbuvvuu43GjRsbLpfLqFOnjtGzZ08zDBkG5/BCnBmIOJfnd+eddxr16tUzXC6XUb9+fePOO+80vv32W3M95/D8HIZhGNb0pgAAAOyBMUQAACDoEYgAAEDQIxABAICgRyACAABBj0AEAACCHoEIAAAEPQIRAAAIegQiAFVOkyZNNHXqVKvLqHTXX3+9Ro8ebXUZQFAgEAE4p2HDhql///7m+0D/gp4zZ45iY2PLLF+7dq1GjBhRqfv+/PPP5XA41KZNG3k8Hp91sbGxmjNnTqXuH0BgEYgABFxRUdEFfb5OnTqqVq2an6r5dd9//73eeOONgOwrEDwej7xer9VlALZDIAJQLsOGDVN6erpeeuklORwOORwO7dq1S5K0ZcsW9enTR1FRUYqLi9PgwYP1008/mZ+9/vrrlZqaqtGjR6t27dpKTk6WJE2ePFnt2rVTZGSkGjZsqPvvv19HjhyRdLJD88c//lH5+fnm/p566ilJZS+Z7dmzR/369VNUVJSio6N1xx136MCBA+b6p556Sh07dtSbb76pJk2aKCYmRgMHDtThw4fPe9x//vOf9eSTT6qwsPCs63ft2iWHw6GNGzeay/Ly8uRwOPT555+bx+JwOPTpp5+qU6dOioiI0A033KDc3Fx98sknatWqlaKjo3XXXXfp2LFjPt9fUlKi1NRUxcTEqHbt2nriiSd0+hOXCgsL9fDDD6t+/fqKjIxU165dzf1Kp7psH374oVq3bi232609e/ac97iBYEMgAlAuL730khITE3Xvvfdq//792r9/vxo2bKi8vDzdcMMN6tSpk77++mstWrRIBw4c0B133OHz+blz58rlcmnlypWaNWuWJMnpdGratGnaunWr5s6dq2XLlunRRx+VJF199dWaOnWqoqOjzf09/PDDZeryer3q16+fDh48qPT0dC1evFjff/+97rzzTp/tvvvuOy1YsEALFy7UwoULlZ6erokTJ573uEePHq2SkhJNnz79t54601NPPaUZM2boq6++UnZ2tu644w5NnTpV8+fP18cff6zPPvuszH7mzp2r0NBQrVmzRi+99JImT56s1157zVyfmpqqjIwMvfXWW9q0aZNuv/129e7dWzt37jS3OXbsmJ5//nm99tpr2rp1q+rWrXvBxwJcdCx+uCwAGxs6dKjRr18/8/2ZTyE3DMN4+umnjV69evksy87ONiQZmZmZ5uc6dep03v29++67Rq1atcz3s2fPNmJiYsps17hxY2PKlCmGYRjGZ599ZoSEhBh79uwx12/dutWQZKxZs8YwDMN48sknjWrVqhkFBQXmNo888ojRtWvXc9ayfPlyQ5Jx6NAhY9asWUbNmjWNvLw8wzAMIyYmxpg9e7ZhGIaRlZVlSDI2bNhgfvbQoUOGJGP58uU+37VkyRJzm7S0NEOS8d1335nL/vSnPxnJycnm++7duxutWrUyvF6vuewvf/mL0apVK8MwDGP37t1GSEiIsXfvXp/ae/bsaYwdO9Y8h5KMjRs3nvNYARgGHSIAF+T//u//tHz5ckVFRZmvli1bSjrZlSnVuXPnMp9dsmSJevbsqfr166t69eoaPHiwfv755zKXjX7N9u3b1bBhQzVs2NBc1rp1a8XGxmr79u3msiZNmqh69erm+3r16ik3N7dc+xg+fLhq1aql559/vtx1nU379u3Nn+Pi4lStWjVdcsklPsvOrOmqq66Sw+Ew3ycmJmrnzp3yeDzavHmzPB6PLrvsMp/zn56e7nPuXS6Xz74BlBVqdQEAqrYjR47o5ptvPmtYqFevnvlzZGSkz7pdu3bppptu0siRI/Xss8+qZs2a+vLLLzV8+HAVFRX5fdB0WFiYz3uHw1HuwcWhoaF69tlnNWzYMKWmpvqsczpP/n+lcdq4nuLi4vPW4HA4Lqgm6eS5DwkJ0bp16xQSEuKzLioqyvw5IiLCJ1QBKItABKDcXC5XmSnol19+ud577z01adJEoaHl/0/KunXr5PV69fe//90MFe+8885593emVq1aKTs7W9nZ2WaXaNu2bcrLy1Pr1q3LXc/53H777XrhhRc0fvx4n+V16tSRJO3fv1+dOnWSJJ8B1hdq9erVPu9XrVql5s2bKyQkRJ06dZLH41Fubq6uvfZav+0TCEZcMgNQbk2aNNHq1au1a9cu/fTTT/J6vUpJSdHBgwf1+9//XmvXrtV3332nTz/9VH/84x9/NcxceumlKi4u1vTp0/X999/rzTffNAdbn76/I0eOaOnSpfrpp5/OeiktKSlJ7dq106BBg7R+/XqtWbNGQ4YMUffu3dWlSxe/Hv/EiRP1+uuv6+jRo+ayiIgIXXXVVZo4caK2b9+u9PR0Pf74437b5549ezRmzBhlZmbqP//5j6ZPn64HHnhAknTZZZdp0KBBGjJkiP7nf/5HWVlZWrNmjdLS0vTxxx/7rQYgGBCIAJTbww8/rJCQELVu3Vp16tTRnj17lJCQoJUrV8rj8ahXr15q166dRo8erdjYWLPzczYdOnTQ5MmT9fzzz6tt27aaN2+e0tLSfLa5+uqrdd999+nOO+9UnTp1NGnSpDLf43A49MEHH6hGjRq67rrrlJSUpEsuuURvv/2234//hhtu0A033KCSkhKf5a+//rpKSkrUuXNnjR49Ws8884zf9jlkyBAdP35cV155pVJSUvTAAw/43JRy9uzZGjJkiB566CG1aNFC/fv319q1a9WoUSO/1QAEA4dx+oVvAACAIESHCAAABD0CEQAACHoEIgAAEPQIRAAAIOgRiAAAQNAjEAEAgKBHIAIAAEGPQAQAAIIegQgAAAQ9AhEAAAh6BCIAABD0CEQAACDo/X/VQ6AX7aNskwAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# make a plot.\n",
    "plt.plot(obj_val)\n",
    "plt.xlabel('Iteration Number')\n",
    "plt.ylabel('loss function')\n",
    "plt.title('loss function plot')\n",
    "plt.savefig(\"Visualization of loss function\", dpi = 150)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 83,
   "id": "73e9067d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAHHCAYAAABXx+fLAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAABNRElEQVR4nO3deXhU9d3+8XtmksxkIQEM2RAJi2WTEEskv6gISiSgjxaVCpRHMFVcwDVaK1IJuAWUIioIFitQHkXUurRWUUwJFY2iLAIKCIjsCYslgUASkjm/P2AOjAmYgcmchHm/rmsuMmfOnPnMgcvcflebYRiGAAAAgojd6gIAAAACjQAEAACCDgEIAAAEHQIQAAAIOgQgAAAQdAhAAAAg6BCAAABA0CEAAQCAoEMAAgAAQYcABMDUu3dv9e7d2+oyvBQUFMhms6mgoMDqUs46N998s5KTk60uA7AEAQhoBGbPni2bzXbSxxdffFHna3333XcaN26cfvzxx/or+DS8+OKLmj17ttVloI6eeuopvfvuu1aXAZy2EKsLAFB3jz32mNq0aVPjePv27et8je+++07jx49X7969a/zf/8cff3ymJZ62F198UbGxsbr55pu9jl922WU6fPiwwsLCrCkMtXrqqac0cOBADRgwwOpSgNNCAAIakf79+ystLa3ert8QQ4bdbpfL5bK6DJ+VlZUpMjKyxnG3263Kysoz+k4nuzaAuqMLDDjLvP766+revbuaNGmi6Ohode3aVc8995yko11pv/3tbyVJl19+udmF5hlf8/MxQJ7xN2+88YbGjx+vli1bqkmTJho4cKBKSkpUUVGh++67T3FxcYqKilJ2drYqKiq86pk1a5auuOIKxcXFyel0qnPnzpo+fbrXOcnJyfr222+1ePFisyZPHScbA/Tmm2+qe/fuCg8PV2xsrP73f/9XO3bs8Drn5ptvVlRUlHbs2KEBAwYoKipKLVq00IMPPqjq6uo63c8PP/xQPXv2VGRkpJo0aaKrr75a3377ba2fs2nTJl111VVq0qSJhg4dKkmy2Wy666679Oqrr6pLly5yOp1asGCBJGnFihXq37+/oqOjFRUVpT59+tTozvR0fy5evFgjR45UXFyczj333JPW67lf8+fP1yOPPKKEhARFRkbq2muv1bZt237x+5aVlemBBx5Qq1at5HQ61aFDB02aNEmGYZjn2Gw2lZWVac6cOebf189b7oCGjhYgoBEpKSnR3r17vY7ZbDadc845kqSFCxdqyJAh6tOnjyZOnChJWrt2rT777DPde++9uuyyy3TPPffo+eef1yOPPKJOnTpJkvnnyeTl5Sk8PFwPP/ywNm7cqBdeeEGhoaGy2+3673//q3HjxumLL77Q7Nmz1aZNG40dO9Z87/Tp09WlSxdde+21CgkJ0T//+U+NHDlSbrdbo0aNkiRNmTJFd999t6KiojRmzBhJUnx8/EnrmT17trKzs3XRRRcpLy9PxcXFeu655/TZZ59pxYoVatq0qXludXW1srKylJ6erkmTJumTTz7Rn//8Z7Vr10533nnnKb/33LlzNXz4cGVlZWnixIk6dOiQpk+frksvvVQrVqzw6kKsqqpSVlaWLr30Uk2aNEkRERHma//+97/1xhtv6K677lJsbKwZ+Hr27Kno6Gg99NBDCg0N1UsvvaTevXtr8eLFSk9P96pl5MiRatGihcaOHauysrJT1i1JTz75pGw2m/74xz9q9+7dmjJlijIzM7Vy5UqFh4fX+h7DMHTttddq0aJFuuWWW5SamqqPPvpIf/jDH7Rjxw49++yz5n259dZb1aNHD912222SpHbt2v1iTUCDYgBo8GbNmmVIqvXhdDrN8+69914jOjraqKqqOum13nzzTUOSsWjRohqv9erVy+jVq5f5fNGiRYYk44ILLjAqKyvN40OGDDFsNpvRv39/r/dnZGQYrVu39jp26NChGp+TlZVltG3b1utYly5dvD775zV46q2srDTi4uKMCy64wDh8+LB53vvvv29IMsaOHWseGz58uCHJeOyxx7yueeGFFxrdu3ev8VknOnDggNG0aVNjxIgRXseLioqMmJgYr+Oez3n44YdrXEeSYbfbjW+//dbr+IABA4ywsDBj06ZN5rGdO3caTZo0MS677DLzmOfv/tJLLz3l36uH5361bNnSKC0tNY+/8cYbhiTjueee86r7xL+vd99915BkPPHEE17XHDhwoGGz2YyNGzeaxyIjI43hw4f/Yj1AQ0UXGNCITJs2TQsXLvR6fPjhh+brTZs2VVlZmRYuXOjXzx02bJhCQ0PN5+np6TIMQ7///e+9zktPT9e2bdtUVVVlHjuxtcHTgtWrVy/98MMPKikp8bmWr7/+Wrt379bIkSO9xtFcffXV6tixo/71r3/VeM8dd9zh9bxnz5764YcfTvk5Cxcu1P79+zVkyBDt3bvXfDgcDqWnp2vRokU13nOyFqVevXqpc+fO5vPq6mp9/PHHGjBggNq2bWseT0xM1O9+9zstWbJEpaWlXtcYMWKEHA7HKWs+0bBhw9SkSRPz+cCBA5WYmKgPPvjgpO/54IMP5HA4dM8993gdf+CBB2QYhte/NaCxowsMaER69OhxykHQI0eO1BtvvKH+/furZcuW6tu3r2688Ub169fvjD73vPPO83oeExMjSWrVqlWN4263WyUlJWa33Geffabc3FwVFhbq0KFDXueXlJSY16qrLVu2SJI6dOhQ47WOHTtqyZIlXsdcLpdatGjhdaxZs2b673//e8rP2bBhgyTpiiuuqPX16Ohor+chISEnHZvz85l7e/bs0aFDh2r9Dp06dZLb7da2bdvUpUuXk17jl5x//vlez202m9q3b3/K5Q+2bNmipKQkr+DkqcnzOnC2IAABZ5G4uDitXLlSH330kT788EN9+OGHmjVrloYNG6Y5c+ac9nVP1vJwsuPGsQGzmzZtUp8+fdSxY0dNnjxZrVq1UlhYmD744AM9++yzcrvdp11TXfnSanIiT21z585VQkJCjddDQrz/8+l0OmW3196ofrIxN77wxzUAHEcAAs4yYWFhuuaaa3TNNdfI7XZr5MiReumll/Too4+qffv2stlsAavln//8pyoqKvSPf/zDqxWptu6jutbVunVrSdL69etrtM6sX7/efP1MeQb1xsXFKTMz0y/X9GjRooUiIiK0fv36Gq+tW7dOdru9RuuarzwtWB6GYWjjxo1KSUk56Xtat26tTz75RAcOHPBqBVq3bp35ukcg/x0B9YExQMBZZN++fV7P7Xa7+QvPMz3ds37M/v37670eT+uLccIU6pKSEs2aNavGuZGRkXWqKS0tTXFxcZoxY4bXlPsPP/xQa9eu1dVXX33mhUvKyspSdHS0nnrqKR05cqTG63v27DntazscDvXt21fvvfeeV5dUcXGxXnvtNV166aU1uth89be//U0HDhwwn7/11lvatWuX+vfvf9L3XHXVVaqurtbUqVO9jj/77LOy2Wxe763r3xfQUNECBDQiH374ofl/4ye6+OKL1bZtW91666366aefdMUVV+jcc8/Vli1b9MILLyg1NdUcx5GamiqHw6GJEyeqpKRETqfTXKfH3/r27Wu2SN1+++06ePCgZs6cqbi4OO3atcvr3O7du2v69Ol64okn1L59e8XFxdU6/iY0NFQTJ05Udna2evXqpSFDhpjT4JOTk3X//ff7pfbo6GhNnz5dN910k379619r8ODBatGihbZu3ap//etfuuSSS2oEBV888cQTWrhwoS699FKNHDlSISEheumll1RRUaGnn376jOtv3ry5Lr30UmVnZ6u4uFhTpkxR+/btNWLEiJO+55prrtHll1+uMWPG6Mcff1S3bt308ccf67333tN9993nNdW9e/fu+uSTTzR58mQlJSWpTZs2NabuAw2apXPQANTJqabBSzJmzZplGIZhvPXWW0bfvn2NuLg4IywszDjvvPOM22+/3di1a5fX9WbOnGm0bdvWcDgcXlPMTzYN/s0336y1nq+++srreG5uriHJ2LNnj3nsH//4h5GSkmK4XC4jOTnZmDhxovHKK68YkozNmzeb5xUVFRlXX3210aRJE0OSWcfPp8F7zJ8/37jwwgsNp9NpNG/e3Bg6dKixfft2r3OGDx9uREZG1rifnjrrYtGiRUZWVpYRExNjuFwuo127dsbNN99sfP3117/4OYZxdBr8qFGjan1t+fLlRlZWlhEVFWVEREQYl19+ufH55597nXOye32qeiUZ8+bNM0aPHm3ExcUZ4eHhxtVXX21s2bLF69yfT4M3jKPT/++//34jKSnJCA0NNc4//3zjmWeeMdxut9d569atMy677DIjPDzckMSUeDQ6NsM4oW0aANCoFRQU6PLLL9ebb76pgQMHWl0O0GAxBggAAAQdAhAAAAg6BCAAABB0GAMEAACCDi1AAAAg6BCAAABA0GEhxFq43W7t3LlTTZo0Ybl3AAAaCcMwdODAASUlJZ10bz4PAlAtdu7cecb78AAAAGts27ZN55577inPIQDVwrMJ4LZt2854Px4AABAYpaWlatWqlddmvidDAKqFp9srOjqaAAQAQCNTl+ErDIIGAABBhwAEAACCDgEIAAAEHQIQAAAIOgQgAAAQdAhAAAAg6BCAAABA0CEAAQCAoEMAAgAAQYcABAAAgg4BCAAABB0CEAAACDpshhpAZRVV+u+hSrlCHYqNclpdDgAAQYsWoAD665LNunTiIv354/VWlwIAQFAjAAVQeKhDklR+xG1xJQAABDcCUAC5Qo/e7vIj1RZXAgBAcGsQAWjatGlKTk6Wy+VSenq6li5detJz3377baWlpalp06aKjIxUamqq5s6d63XOzTffLJvN5vXo169ffX+NX+Q81gJ0mAAEAIClLB8EPX/+fOXk5GjGjBlKT0/XlClTlJWVpfXr1ysuLq7G+c2bN9eYMWPUsWNHhYWF6f3331d2drbi4uKUlZVlntevXz/NmjXLfO50Wj/o2GV2gRGAAACwkuUtQJMnT9aIESOUnZ2tzp07a8aMGYqIiNArr7xS6/m9e/fWddddp06dOqldu3a69957lZKSoiVLlnid53Q6lZCQYD6aNWsWiK9zSowBAgCgYbA0AFVWVmrZsmXKzMw0j9ntdmVmZqqwsPAX328YhvLz87V+/XpddtllXq8VFBQoLi5OHTp00J133ql9+/ad9DoVFRUqLS31etQHxgABANAwWNoFtnfvXlVXVys+Pt7reHx8vNatW3fS95WUlKhly5aqqKiQw+HQiy++qCuvvNJ8vV+/frr++uvVpk0bbdq0SY888oj69++vwsJCORyOGtfLy8vT+PHj/ffFTsLTBVZRRQsQAABWsnwM0Olo0qSJVq5cqYMHDyo/P185OTlq27atevfuLUkaPHiweW7Xrl2VkpKidu3aqaCgQH369KlxvdGjRysnJ8d8XlpaqlatWvm9blfIsUHQlbQAAQBgJUsDUGxsrBwOh4qLi72OFxcXKyEh4aTvs9vtat++vSQpNTVVa9euVV5enhmAfq5t27aKjY3Vxo0baw1ATqczIIOkzS6wKgIQAABWsnQMUFhYmLp37678/HzzmNvtVn5+vjIyMup8HbfbrYqKipO+vn37du3bt0+JiYlnVO+ZYhYYAAANg+VdYDk5ORo+fLjS0tLUo0cPTZkyRWVlZcrOzpYkDRs2TC1btlReXp6ko+N10tLS1K5dO1VUVOiDDz7Q3LlzNX36dEnSwYMHNX78eN1www1KSEjQpk2b9NBDD6l9+/Ze0+St4DphFphhGLLZbJbWAwBAsLI8AA0aNEh79uzR2LFjVVRUpNTUVC1YsMAcGL1161bZ7ccbqsrKyjRy5Eht375d4eHh6tixo/7v//5PgwYNkiQ5HA6tWrVKc+bM0f79+5WUlKS+ffvq8ccft3wtIE8XmHR0ILQnEAEAgMCyGYZhWF1EQ1NaWqqYmBiVlJQoOjrab9c9Uu3W+WM+lCStHHulmkaE+e3aAAAEO19+f1u+EGIwCXXYFWI/2u3FYogAAFiHABRgDIQGAMB6BKAAYyo8AADWIwAFmJPFEAEAsBwBKMCO7wfGGCAAAKxCAAqw8LBjY4DoAgMAwDIEoADz7AdWwSBoAAAsQwAKMM8ssMMEIAAALEMACjDGAAEAYD0CUICxDhAAANYjAAXYiRuiAgAAaxCAAszTBcYYIAAArEMACjBmgQEAYD0CUIAxBggAAOsRgALMXAiRMUAAAFiGABRgzhA2QwUAwGoEoAAzF0JkM1QAACxDAAowcwxQFV1gAABYhQAUYMdXgqYFCAAAqxCAAiw8lGnwAABYjQAUYGyGCgCA9QhAAcZmqAAAWI8AFGDOEBZCBADAagSgADu+ECIBCAAAqxCAAoxp8AAAWI8AFGCuYytBV1a5Ve02LK4GAIDgRAAKME8LkCRVsB0GAACWIAAF2IkBiJlgAABYgwAUYA67TWEOVoMGAMBKBCALOI+tBcRiiAAAWIMAZAFzJhgBCAAASxCALMBq0AAAWIsAZAE2RAUAwFoEIAscXwyRAAQAgBUIQBZwHdsP7HAlXWAAAFiBAGQBZyjT4AEAsBIByAJ0gQEAYC0CkAXCzWnwdIEBAGAFApAFXHSBAQBgKQKQBVgIEQAAaxGALEAAAgDAWgQgC7gYAwQAgKUIQBZgDBAAANZqEAFo2rRpSk5OlsvlUnp6upYuXXrSc99++22lpaWpadOmioyMVGpqqubOnet1jmEYGjt2rBITExUeHq7MzExt2LChvr9GnZkLIRKAAACwhOUBaP78+crJyVFubq6WL1+ubt26KSsrS7t37671/ObNm2vMmDEqLCzUqlWrlJ2drezsbH300UfmOU8//bSef/55zZgxQ19++aUiIyOVlZWl8vLyQH2tU6ILDAAAa1kegCZPnqwRI0YoOztbnTt31owZMxQREaFXXnml1vN79+6t6667Tp06dVK7du107733KiUlRUuWLJF0tPVnypQp+tOf/qTf/OY3SklJ0d/+9jft3LlT7777bgC/2cl5usAqWAgRAABLWBqAKisrtWzZMmVmZprH7Ha7MjMzVVhY+IvvNwxD+fn5Wr9+vS677DJJ0ubNm1VUVOR1zZiYGKWnp9fpmoEQziwwAAAsFWLlh+/du1fV1dWKj4/3Oh4fH69169ad9H0lJSVq2bKlKioq5HA49OKLL+rKK6+UJBUVFZnX+Pk1Pa/9XEVFhSoqKsznpaWlp/V96srTBcYYIAAArGFpADpdTZo00cqVK3Xw4EHl5+crJydHbdu2Ve/evU/renl5eRo/frx/izyF45uhMgYIAAArWNoFFhsbK4fDoeLiYq/jxcXFSkhIOOn77Ha72rdvr9TUVD3wwAMaOHCg8vLyJMl8ny/XHD16tEpKSszHtm3bzuRr/SIWQgQAwFqWBqCwsDB1795d+fn55jG32638/HxlZGTU+Tput9vswmrTpo0SEhK8rllaWqovv/zypNd0Op2Kjo72etQnzzR4WoAAALCG5V1gOTk5Gj58uNLS0tSjRw9NmTJFZWVlys7OliQNGzZMLVu2NFt48vLylJaWpnbt2qmiokIffPCB5s6dq+nTp0uSbDab7rvvPj3xxBM6//zz1aZNGz366KNKSkrSgAEDrPqaXsLDjgagClqAAACwhOUBaNCgQdqzZ4/Gjh2roqIipaamasGCBeYg5q1bt8puP95QVVZWppEjR2r79u0KDw9Xx44d9X//938aNGiQec5DDz2ksrIy3Xbbbdq/f78uvfRSLViwQC6XK+DfrzaeafAMggYAwBo2wzAMq4toaEpLSxUTE6OSkpJ66Q77b1mlLnx8oSRp45P9FeKwfDkmAAAaPV9+f/Ob1wKeQdCSVF7FOCAAAAKNAGQBZ8jx285MMAAAAo8AZAG73WaGoMOVBCAAAAKNAGQRTzcY+4EBABB4BCCLuFgNGgAAyxCALMJq0AAAWIcAZJHjO8LTAgQAQKARgCziZEd4AAAsQwCyiCvEMwaIAAQAQKARgCzCGCAAAKxDALKIOQaIlaABAAg4ApBFzGnwLIQIAEDAEYAsQhcYAADWIQBZxAxArAQNAEDAEYAs4mQlaAAALEMAskg4XWAAAFiGAGQRFwshAgBgGQKQRTwLIVbQBQYAQMARgCzCLDAAAKxDALJIeBizwAAAsAoByCLOkGNjgFgIEQCAgCMAWcTFNHgAACxDALIICyECAGAdApBFPAGIWWAAAAQeAcgiLIQIAIB1CEAW8YwBYiFEAAACjwBkkRPXATIMw+JqAAAILgQgi7iOTYN3G9KRagIQAACBRACyiCvs+K1nJhgAAIFFALJImMMum+3oz+UshggAQEARgCxis9nMbjAWQwQAILAIQBYyV4OmCwwAgIAiAFmIHeEBALAGAchCxxdDpAsMAIBAIgBZyHksALEYIgAAgUUAstDxHeEJQAAABBIByELHZ4ERgAAACCQCkIU8LUDsCA8AQGARgCwUHsYYIAAArEAAshBdYAAAWIMAZCEn0+ABALAEAchCrAQNAIA1CEAWCmclaAAALEEAshBbYQAAYI0GEYCmTZum5ORkuVwupaena+nSpSc9d+bMmerZs6eaNWumZs2aKTMzs8b5N998s2w2m9ejX79+9f01fHZ8IUTGAAEAEEiWB6D58+crJydHubm5Wr58ubp166asrCzt3r271vMLCgo0ZMgQLVq0SIWFhWrVqpX69u2rHTt2eJ3Xr18/7dq1y3zMmzcvEF/HJ7QAAQBgDZ8C0JEjR9SuXTutXbvWbwVMnjxZI0aMUHZ2tjp37qwZM2YoIiJCr7zySq3nv/rqqxo5cqRSU1PVsWNHvfzyy3K73crPz/c6z+l0KiEhwXw0a9bMbzX7C9PgAQCwhk8BKDQ0VOXl5X778MrKSi1btkyZmZnHC7LblZmZqcLCwjpd49ChQzpy5IiaN2/udbygoEBxcXHq0KGD7rzzTu3bt++k16ioqFBpaanXIxBcLIQIAIAlfO4CGzVqlCZOnKiqqqoz/vC9e/equrpa8fHxXsfj4+NVVFRUp2v88Y9/VFJSkleI6tevn/72t78pPz9fEydO1OLFi9W/f39VV9ceNPLy8hQTE2M+WrVqdfpfygeuEMYAAQBghRBf3/DVV18pPz9fH3/8sbp27arIyEiv199++22/FfdLJkyYoNdff10FBQVyuVzm8cGDB5s/d+3aVSkpKWrXrp0KCgrUp0+fGtcZPXq0cnJyzOelpaUBCUGMAQIAwBo+B6CmTZvqhhtu8MuHx8bGyuFwqLi42Ot4cXGxEhISTvneSZMmacKECfrkk0+UkpJyynPbtm2r2NhYbdy4sdYA5HQ65XQ6ff8CZ8gTgCqqaAECACCQfA5As2bN8tuHh4WFqXv37srPz9eAAQMkyRzQfNddd530fU8//bSefPJJffTRR0pLS/vFz9m+fbv27dunxMREf5XuFyyECACANXwOQB579uzR+vXrJUkdOnRQixYtTus6OTk5Gj58uNLS0tSjRw9NmTJFZWVlys7OliQNGzZMLVu2VF5eniRp4sSJGjt2rF577TUlJyebY4WioqIUFRWlgwcPavz48brhhhuUkJCgTZs26aGHHlL79u2VlZV1ul+3XnjWAWIQNAAAgeVzACorK9Pdd9+tv/3tb3K7j3bdOBwODRs2TC+88IIiIiJ8ut6gQYO0Z88ejR07VkVFRUpNTdWCBQvMgdFbt26V3X58rPb06dNVWVmpgQMHel0nNzdX48aNk8Ph0KpVqzRnzhzt379fSUlJ6tu3rx5//HFLurlOhTFAAABYw2YYhuHLG26//XZ98sknmjp1qi655BJJ0pIlS3TPPffoyiuv1PTp0+ul0EAqLS1VTEyMSkpKFB0dXW+fs/tAuXo8eXT9os15V8lms9XbZwEAcLbz5fe3zy1Af//73/XWW2+pd+/e5rGrrrpK4eHhuvHGG8+KABQonhYg6ehA6BOfAwCA+uPzOkCHDh2qsW6PJMXFxenQoUN+KSpYhJ8QeOgGAwAgcHwOQBkZGcrNzfVaEfrw4cMaP368MjIy/Frc2S7UYZfDfrTbi8UQAQAIHJ+7wKZMmaJ+/frp3HPPVbdu3SRJ33zzjVwulz766CO/F3i2c4XYVVZZTQsQAAAB5HMA6tq1qzZs2KBXX31V69atkyQNGTJEQ4cOVXh4uN8LPNu5Qh1HA1AVAQgAgEDxKQAdOXJEHTt21Pvvv68RI0bUV01B5fhUeLrAAAAIFEt3g8cJiyFW0gIEAECgWLobPE5oAaILDACAgGnUu8GfDcwNURkEDQBAwFi6GzyOd4ExBggAgMDxKQBVVVXp8ssvV9++fZWQkFBfNQUVz2KIbIgKAEDg+DQGKCQkRHfccYcqKirqq56g42RDVAAAAs7nQdA9evTQihUr6qOWoOQKYRo8AACB5vMYoJEjR+qBBx7Q9u3b1b179xqDoFNSUvxWXDA4PgaIFiAAAALF5wA0ePBgSdI999xjHrPZbDIMQzabTdXV/CL3BdPgAQAIPJ8D0ObNm+ujjqDlGQRdzkKIAAAEjM8BqHXr1vVRR9BiGjwAAIHn8yBoSZo7d64uueQSJSUlacuWLZKO7hL/3nvv+bW4YEAXGAAAgedzAJo+fbpycnJ01VVXaf/+/eaYn6ZNm2rKlCn+ru+sxzR4AAACz+cA9MILL2jmzJkaM2aMHA6HeTwtLU2rV6/2a3HB4PhCiHSBAQAQKD4HoM2bN+vCCy+scdzpdKqsrMwvRQUTpsEDABB4PgegNm3aaOXKlTWOL1iwQJ06dfJHTUHFsxAim6ECABA4Ps8Cy8nJ0ahRo1ReXi7DMLR06VLNmzdPeXl5evnll+ujxrOaOQiaLjAAAALG5wB06623Kjw8XH/605906NAh/e53v1NSUpKee+45c5FE1J3ZBcYsMAAAAsbnACRJQ4cO1dChQ3Xo0CEdPHhQcXFx/q4raHhagA6zECIAAAFzWgHIIyIiQhEREf6qJSi5mAYPAEDAndZCiPCf411gjAECACBQCEAW87QAVVa55XYbFlcDAEBwIABZzLMQosRAaAAAAoUAZDHXiQGIqfAAAATEaQ2Czs/PV35+vnbv3i232/uX9iuvvOKXwoKFw25TqMOmI9UGA6EBAAgQnwPQ+PHj9dhjjyktLU2JiYmy2Wz1UVdQcYU4dKS6igAEAECA+ByAZsyYodmzZ+umm26qj3qCkjPUoQMVVXSBAQAQID6PAaqsrNTFF19cH7UErfCwo38Nh2kBAgAgIHwOQLfeeqtee+21+qglaLEhKgAAgeVzF1h5ebn+8pe/6JNPPlFKSopCQ0O9Xp88ebLfigsW5mrQTIMHACAgfA5Aq1atUmpqqiRpzZo1Xq8xIPr0mKtBMwYIAICA8DkALVq0qD7qCGpsiAoAQGCd0UKI27dv1/bt2/1VS9CiCwwAgMDyOQC53W499thjiomJUevWrdW6dWs1bdpUjz/+eI1FEVE3x3eE5/4BABAIPneBjRkzRn/96181YcIEXXLJJZKkJUuWaNy4cSovL9eTTz7p9yLPdq4QzxggWoAAAAgEnwPQnDlz9PLLL+vaa681j6WkpKhly5YaOXIkAeg0eFqAmAYPAEBg+NwF9tNPP6ljx441jnfs2FE//fSTX4oKNuFhxwZBE4AAAAgInwNQt27dNHXq1BrHp06dqm7duvmlqGBzvAuMMUAAAASCzwHo6aef1iuvvKLOnTvrlltu0S233KLOnTtr9uzZeuaZZ06riGnTpik5OVkul0vp6elaunTpSc+dOXOmevbsqWbNmqlZs2bKzMyscb5hGBo7dqwSExMVHh6uzMxMbdiw4bRqCwSnOQiaFiAAAALB5wDUq1cvff/997ruuuu0f/9+7d+/X9dff73Wr1+vnj17+lzA/PnzlZOTo9zcXC1fvlzdunVTVlaWdu/eXev5BQUFGjJkiBYtWqTCwkK1atVKffv21Y4dO8xznn76aT3//POaMWOGvvzyS0VGRiorK0vl5eU+1xcIx6fB0wIEAEAg2AzDMKwsID09XRdddJHZreZ2u9WqVSvdfffdevjhh3/x/dXV1WrWrJmmTp2qYcOGyTAMJSUl6YEHHtCDDz4oSSopKVF8fLxmz56twYMH/+I1S0tLFRMTo5KSEkVHR5/ZF6yDV7/cojHvrFFmp3i9PDyt3j8PAICzkS+/v+s0C2zVqlW64IILZLfbtWrVqlOem5KSUudCKysrtWzZMo0ePdo8ZrfblZmZqcLCwjpd49ChQzpy5IiaN28uSdq8ebOKioqUmZlpnhMTE6P09HQVFhbWGoAqKipUUVFhPi8tLa3zd/CHcM8sMBZCBAAgIOoUgFJTU1VUVKS4uDilpqbKZrOptoYjm82m6uq6/xLfu3evqqurFR8f73U8Pj5e69atq9M1/vjHPyopKckMPEVFReY1fn5Nz2s/l5eXp/Hjx9e5bn9zMQYIAICAqlMA2rx5s1q0aGH+3FBMmDBBr7/+ugoKCuRyuU77OqNHj1ZOTo75vLS0VK1atfJHiXXCZqgAAARWnQJQ69atzZ+3bNmiiy++WCEh3m+tqqrS559/7nXuL4mNjZXD4VBxcbHX8eLiYiUkJJzyvZMmTdKECRP0ySefeHW7ed5XXFysxMREr2t6drH/OafTKafTWee6/c0VQgsQAACB5PMssMsvv7zWBQ9LSkp0+eWX+3StsLAwde/eXfn5+eYxt9ut/Px8ZWRknPR9Tz/9tB5//HEtWLBAaWneg4bbtGmjhIQEr2uWlpbqyy+/POU1reRiIUQAAALK560wDMOQzWarcXzfvn2KjIz0uYCcnBwNHz5caWlp6tGjh6ZMmaKysjJlZ2dLkoYNG6aWLVsqLy9PkjRx4kSNHTtWr732mpKTk81xPVFRUYqKipLNZtN9992nJ554Queff77atGmjRx99VElJSRowYIDP9QXC8RYgusAAAAiEOgeg66+/XtLRgc4333yzV5dRdXW1Vq1apYsvvtjnAgYNGqQ9e/Zo7NixKioqUmpqqhYsWGAOYt66davs9uMNVdOnT1dlZaUGDhzodZ3c3FyNGzdOkvTQQw+prKxMt912m/bv369LL71UCxYsOKNxQvXJMwaIvcAAAAiMOq8D5GmRmTNnjm688UaFh4ebr4WFhSk5OVkjRoxQbGxs/VQaQIFeB2jn/sO6eMK/FeqwacOTV9X75wEAcDby+zpAkjRr1ixJUnJysh588MHT6u5C7TzT4I9UG6qqdivE4fPQLAAA4AOfxwDl5ubWRx1BzbMQonR0O4woAhAAAPXK5wAkSW+99ZbeeOMNbd26VZWVlV6vLV++3C+FBRNnyPHAU36kWlHO0/prAQAAdeRzU8Pzzz+v7OxsxcfHa8WKFerRo4fOOecc/fDDD+rfv3991HjWs9ttCgvxLIbIQGgAAOqbzwHoxRdf1F/+8he98MILCgsL00MPPaSFCxfqnnvuUUlJSX3UGBRcIawGDQBAoPgcgLZu3WpOdw8PD9eBAwckSTfddJPmzZvn3+qCSHgYq0EDABAoPgeghIQEcyXo8847T1988YWko3uE1XFGPWrBhqgAAASOzwHoiiuu0D/+8Q9JR9cGuv/++3XllVdq0KBBuu666/xeYLBgNWgAAALH5+lGf/nLX+R2H/0lPWrUKJ1zzjn6/PPPde211+r222/3e4HB4viO8LQAAQBQ33wOQHa73WtrisGDB2vw4MF+LSoYOUPZEBUAgEA5rQVnysvLtWrVKu3evdtsDfK49tpr/VJYsAlnDBAAAAHjcwBasGCBhg0bpr1799Z4zWazqbqaX+Cnw+wCq2IMEAAA9c3nQdB33323fvvb32rXrl1yu91eD8LP6fPMAmNHeAAA6p/PAai4uFg5OTmKj4+vj3qC1vFZYAQgAADqm88BaODAgSooKKiHUoKbZyFEBkEDAFD/fB4DNHXqVP32t7/Vp59+qq5duyo0NNTr9XvuucdvxQUTZyhbYQAAECg+B6B58+bp448/lsvlUkFBgWw2m/mazWYjAJ0musAAAAgcnwPQmDFjNH78eD388MNe6wHhzBzfCoMWIAAA6pvPCaayslKDBg0i/PgZK0EDABA4PqeY4cOHa/78+fVRS1BjIUQAAALH5y6w6upqPf300/roo4+UkpJSYxD05MmT/VZcMDG7wKoIQAAA1DefA9Dq1at14YUXSpLWrFnj9dqJA6LhGxezwAAACBifA9CiRYvqo46g56QLDACAgGEkcwMRzm7wAAAETJ1agK6//nrNnj1b0dHRuv7660957ttvv+2XwoLN8b3A6AIDAKC+1SkAxcTEmON7oqOjGetTD5gGDwBA4NQpAM2aNcv8efbs2fVVS1BjJWgAAALH5zFAV1xxhfbv31/jeGlpqa644gp/1BSUXCeMATIMw+JqAAA4u/kcgAoKClRZWVnjeHl5uT799FO/FBWMPIOg3YZ0pJoABABAfarzNPhVq1aZP3/33XcqKioyn1dXV2vBggVq2bKlf6sLIp7d4KWjiyGGhTBBDwCA+lLnAJSamiqbzSabzVZrV1d4eLheeOEFvxYXTJwhdtlskmEcHQcU7Qr95TcBAIDTUucAtHnzZhmGobZt22rp0qVq0aKF+VpYWJji4uLkcDjqpchgYLPZ5Ayxq/yIm6nwAADUszoHoNatW0uS3G5+OdcXV6hD5UfcLIYIAEA983mgyZw5c/Svf/3LfP7QQw+padOmuvjii7Vlyxa/Fhds2BEeAIDA8DkAPfXUUwoPD5ckFRYWaurUqXr66acVGxur+++/3+8FBhNzR3i6wAAAqFc+b4a6bds2tW/fXpL07rvvauDAgbrtttt0ySWXqHfv3v6uL6g4Q1gNGgCAQPC5BSgqKkr79u2TJH388ce68sorJUkul0uHDx/2b3VBxsWGqAAABITPLUBXXnmlbr31Vl144YX6/vvvddVVV0mSvv32WyUnJ/u7vqDCGCAAAALD5xagadOmKSMjQ3v27NHf//53nXPOOZKkZcuWaciQIX4vMJh4NkRlGjwAAPXL5xagpk2baurUqTWOjx8/3i8FBTNzEHQVLUAAANSn09pv4dNPP9X//u//6uKLL9aOHTskSXPnztWSJUv8WlywcdEFBgBAQPgcgP7+978rKytL4eHhWr58uSoqKiRJJSUleuqpp/xeYDDxdIEdrqQLDACA+uRzAHriiSc0Y8YMzZw5U6Ghx/eruuSSS7R8+XK/Fhds6AIDACAwfA5A69ev12WXXVbjeExMjPbv3+9zAdOmTVNycrJcLpfS09O1dOnSk5777bff6oYbblBycrJsNpumTJlS45xx48aZm7Z6Hh07dvS5LivQBQYAQGD4HIASEhK0cePGGseXLFmitm3b+nSt+fPnKycnR7m5uVq+fLm6deumrKws7d69u9bzDx06pLZt22rChAlKSEg46XW7dOmiXbt2mY/GMjbJFcJK0AAABILPAWjEiBG699579eWXX8pms2nnzp169dVX9eCDD+rOO+/06VqTJ0/WiBEjlJ2drc6dO2vGjBmKiIjQK6+8Uuv5F110kZ555hkNHjxYTqfzpNcNCQlRQkKC+YiNjfWpLqt4xgDRAgQAQP3yeRr8ww8/LLfbrT59+ujQoUO67LLL5HQ69eCDD+ruu++u83UqKyu1bNkyjR492jxmt9uVmZmpwsJCX8vysmHDBiUlJcnlcikjI0N5eXk677zzzuiagRAeRhcYAACB4HMAstlsGjNmjP7whz9o48aNOnjwoDp37qyoqCifrrN3715VV1crPj7e63h8fLzWrVvna1mm9PR0zZ49Wx06dNCuXbs0fvx49ezZU2vWrFGTJk1qfU9FRYU5m02SSktLT/vzz8TxLjACEAAA9cnnAOQRFhamzp07+7MWv+jfv7/5c0pKitLT09W6dWu98cYbuuWWW2p9T15eXoNYyNFpdoExBggAgPp0Wgsh+kNsbKwcDoeKi4u9jhcXF59ygLOvmjZtql/96le1Dtz2GD16tEpKSszHtm3b/Pb5vmAaPAAAgWFZAAoLC1P37t2Vn59vHnO73crPz1dGRobfPufgwYPatGmTEhMTT3qO0+lUdHS018MKEcfGAJVVVFny+QAABIvT7gLzh5ycHA0fPlxpaWnq0aOHpkyZorKyMmVnZ0uShg0bppYtWyovL0/S0YHT3333nfnzjh07tHLlSkVFRal9+/aSpAcffFDXXHONWrdurZ07dyo3N1cOh6NRbNQaH+2SJBWVlFtcCQAAZzdLA9CgQYO0Z88ejR07VkVFRUpNTdWCBQvMgdFbt26V3X68kWrnzp268MILzeeTJk3SpEmT1KtXLxUUFEiStm/friFDhmjfvn1q0aKFLr30Un3xxRdq0aJFQL/b6UiMORqASsurdLCiSlFOS/96AAA4a9kMwzCsLqKhKS0tVUxMjEpKSgLeHdY19yMdqKjSwvsv0/nxtc9aAwAANfny+9uyMUCoXWLTo61AO+kGAwCg3hCAGpjEmHBJ0q79hy2uBACAsxcBqIFJogUIAIB6RwBqYJJoAQIAoN4RgBqYxKbHAhAtQAAA1BsCUAOTFOPpAqMFCACA+kIAamDMFqD95WKFAgAA6gcBqIHxLIZ4+Ei1Sg4fsbgaAADOTgSgBsYV6lDzyDBJ0s79jAMCAKA+EIAaIE8r0C7GAQEAUC8IQA2QZzHEnUyFBwCgXhCAGiAWQwQAoH4RgBogtsMAAKB+EYAaIFqAAACoXwSgBijJXA2aFiAAAOoDAagB8swCKyopl9vNYogAAPgbAagBio92yWaTjlQb2ltWYXU5AACcdQhADVCow664Jk5JR7fEAAAA/kUAaqDMmWCMAwIAwO8IQA2UOROMFiAAAPyOANRA0QIEAED9IQA1UJ6ZYKwFBACA/xGAGijPWkDsBwYAgP8RgBooc0d4xgABAOB3BKAGquWxFqDdB8pVVe22uBoAAM4uBKAGKjbKqVCHTW5DKj7AYogAAPgTAaiBstttio/2dIMxDggAAH8iADVgScemwjMTDAAA/yIANWCJTWkBAgCgPhCAGrDjiyHSAgQAgD8RgBqw49th0AIEAIA/EYAaMFqAAACoHwSgBsxcDJH9wAAA8CsCUAPm2Q5j78FKlR+ptrgaAADOHgSgBqxZRKicIUf/ioroBgMAwG8IQA2YzWYzt8TYSTcYAAB+QwBq4I6vBUQLEAAA/kIAauCOzwSjBQgAAH8hADVwScdmgrEdBgAA/kMAauASj40BYjsMAAD8hwDUwB1fC4gWIAAA/IUA1MB51gJiOwwAAPyHANTAeVqASsurVFZRZXE1AACcHQhADVwTV6iaOEMkMRMMAAB/sTwATZs2TcnJyXK5XEpPT9fSpUtPeu63336rG264QcnJybLZbJoyZcoZX7MxSDR3hWccEAAA/mBpAJo/f75ycnKUm5ur5cuXq1u3bsrKytLu3btrPf/QoUNq27atJkyYoISEBL9cszFgLSAAAPzL0gA0efJkjRgxQtnZ2ercubNmzJihiIgIvfLKK7Wef9FFF+mZZ57R4MGD5XQ6/XLNxsAzEHoHLUAAAPiFZQGosrJSy5YtU2Zm5vFi7HZlZmaqsLAwoNesqKhQaWmp16Mh8SyGyFpAAAD4h2UBaO/evaqurlZ8fLzX8fj4eBUVFQX0mnl5eYqJiTEfrVq1Oq3Pry/mYoisBQQAgF9YPgi6IRg9erRKSkrMx7Zt26wuycvx7TBoAQIAwB9CrPrg2NhYORwOFRcXex0vLi4+6QDn+rqm0+k86ZiihuD4dhjlMgxDNpvN4ooAAGjcLGsBCgsLU/fu3ZWfn28ec7vdys/PV0ZGRoO5ZkPgWQzx8JFqlRw+YnE1AAA0fpa1AElSTk6Ohg8frrS0NPXo0UNTpkxRWVmZsrOzJUnDhg1Ty5YtlZeXJ+noIOfvvvvO/HnHjh1auXKloqKi1L59+zpdszFyhTrUPDJMP5VVauf+cjWNCLO6JAAAGjVLA9CgQYO0Z88ejR07VkVFRUpNTdWCBQvMQcxbt26V3X68kWrnzp268MILzeeTJk3SpEmT1KtXLxUUFNTpmo1VYoxLP5VValfJYXVOira6HAAAGjWbYRiG1UU0NKWlpYqJiVFJSYmioxtG2Lh1ztf6ZG2xHh9wgW76f62tLgcAgAbHl9/fzAJrJJKashYQAAD+QgBqJJJYCwgAAL8hADUSnplgO2kBAgDgjBGAGglPCxCLIQIAcOYIQI2EpwWoqKRcbjfj1gEAOBMEoEYiPtolm006Um1ob1mF1eUAANCoEYAaiVCHXXFNjm7XsWs/A6EBADgTBKBGJDHGMxOMcUAAAJwJAlAj4lkLaCctQAAAnBECUCNCCxAAAP5BAGpEzLWAWAwRAIAzQgBqRMzVoFkMEQCAM0IAakTYDgMAAP8gADUinkHQxaXlOlhRZXE1AAA0XgSgRqRFlFNtYiPlNqSC9butLgcAgEaLANSI2Gw2ZXVJkCQtWFNkcTUAADReBKBGpt8FRwPQonW7VX6k2uJqAABonAhAjUxKyxglRLtUVlmtzzfttbocAAAaJQJQI2O325TVJV4S3WAAAJwuAlAjlHWsG2zhd8WqqnZbXA0AAI0PAagR6pHcXM0iQvXfQ0e09MefrC4HAIBGhwDUCIU47MrsdLQb7ONviy2uBgCAxocA1Eh5ZoMtWFMkt9uwuBoAABoXAlAjdUn7WEWGOVRUWq5VO0qsLgcAgEaFANRIuUIdurxjnCRmgwEA4CsCUCPmWRX6o2+LZBh0gwEAUFcEoEbs8o5xCnPYtXlvmTbsPmh1OQAANBoEoEYsyhminufHSqIbDAAAXxCAGrmsC9gcFQAAXxGAGrnMTvGy26TvdpVq20+HrC4HAIBGgQDUyDWPDFN6m3MkHR0MDQAAfhkB6CzQj24wAAB8QgA6C/Q9tjv8sq3/1e7ScourAQCg4SMAnQUSY8LVrVVTGYb08XfsDQYAwC8hAJ0l+p2wKCIAADg1AtBZIutYN1jhpn0qOXTE4moAAGjYCEBnibYtotQhvomq3Iby19ENBgDAqRCAziKeViBmgwEAcGoEoLOIZ1Xo/2zYo/+WVVpcDQAADRcB6CzSOTFabVtEqvyIW799qVDb/8vK0AAA1IYAdBax2Wx6ceivlRDt0sbdB3X9i5/r250lVpcFAECDQwA6y3RMiNY7oy5Wh/gm2n2gQjfOKNR/vt9jdVkAADQoBKCzUGJMuN64I0P/r21zlVVW6/ezv9Jby7ZbXRYAAA0GAegsFRMeqjm/76FruyWpym3owTe/0dR/b5BhGFaXBgCA5RpEAJo2bZqSk5PlcrmUnp6upUuXnvL8N998Ux07dpTL5VLXrl31wQcfeL1+8803y2azeT369etXn1+hQXKGODRlUKru6NVOkjTp4+/1yDtrVFXttrgyAACsZXkAmj9/vnJycpSbm6vly5erW7duysrK0u7du2s9//PPP9eQIUN0yy23aMWKFRowYIAGDBigNWvWeJ3Xr18/7dq1y3zMmzcvEF+nwbHbbXq4f0c99psustmkeUu36ra5y7R2VymtQQCAoGUzLP4tmJ6erosuukhTp06VJLndbrVq1Up33323Hn744RrnDxo0SGVlZXr//ffNY//v//0/paamasaMGZKOtgDt379f77777mnVVFpaqpiYGJWUlCg6Ovq0rtEQLVhTpHtfX6GKqqMtQK3PiVDfzvHK6pKgX5/XTHa7zeIKAQA4fb78/ra0BaiyslLLli1TZmamecxutyszM1OFhYW1vqewsNDrfEnKysqqcX5BQYHi4uLUoUMH3Xnnndq3b99J66ioqFBpaanX42zU74IEvXF7hjI7xSksxK4t+w5p5qebNXBGoXo8la/Rb69WwfrdqqyiiwwAcHYLsfLD9+7dq+rqasXHx3sdj4+P17p162p9T1FRUa3nFxUd3/6hX79+uv7669WmTRtt2rRJjzzyiPr376/CwkI5HI4a18zLy9P48eP98I0avm6tmurl4ReprKJKi7/fo4++LdK/1+7W3oMVmrd0q+Yt3arIMIdaNHEqPCxEEWEORYQ55Ap1mD+Hh4Yo0ulQRFiIopwORTpDjv189HikM0QOu01ut6Eqt6Fq80+3qt1SldstGZLDblOIw64Qu00Ou02hDvuxP22y22xyG0ffe/RPeT13G5JNkt1mk91+7E+bTXbb0fWQ7DbJkGQYkmEYMnT0/YZx/E9Jstkkm2xH/zzx52OvGcYJ19HRNx295vF7arP97M9j1zjRied7rlMbm07dCneq957uNQHACs0iQxXXxGXZ51sagOrL4MGDzZ+7du2qlJQUtWvXTgUFBerTp0+N80ePHq2cnBzzeWlpqVq1ahWQWq0S6QzRVV0TdVXXRFVWuVX4wz599G2RFn5XrD0HKlS2j1WkAQD1Z2TvdnqoX0fLPt/SABQbGyuHw6HiYu/dy4uLi5WQkFDrexISEnw6X5Latm2r2NhYbdy4sdYA5HQ65XQ6T+MbnB3CQuzq9asW6vWrFnriNxdow+6DKi0/okOV1TpcWXX0zyPVOlxZrUPmo0plFdUqq6hSWWXV0T8rqnXw2PNqt3GsZed4C4/DbjN/lnRCy5ChI9Vurz+rDUN2m00Om032Y++x22xy2CXHsZl9ko61Bh1tETLMFiLJfayJ6OetQpLn+dHvfmILj2Sc8PxoW8vRlqDjbSi2Y9c4sYXneOuO4fXc834P71ah2lplvFt3DOPn7/GdP64BAPUh0mltG4ylnx4WFqbu3bsrPz9fAwYMkHR0EHR+fr7uuuuuWt+TkZGh/Px83XfffeaxhQsXKiMj46Sfs337du3bt0+JiYn+LP+sZLfb1CGhidVlAABQryyfBp+Tk6OZM2dqzpw5Wrt2re68806VlZUpOztbkjRs2DCNHj3aPP/ee+/VggUL9Oc//1nr1q3TuHHj9PXXX5uB6eDBg/rDH/6gL774Qj/++KPy8/P1m9/8Ru3bt1dWVpYl3xEAADQslo8BGjRokPbs2aOxY8eqqKhIqampWrBggTnQeevWrbLbj+e0iy++WK+99pr+9Kc/6ZFHHtH555+vd999VxdccIEkyeFwaNWqVZozZ47279+vpKQk9e3bV48//nhQd3MBAIDjLF8HqCE6W9cBAgDgbNZo1gECAACwAgEIAAAEHQIQAAAIOgQgAAAQdAhAAAAg6BCAAABA0CEAAQCAoEMAAgAAQYcABAAAgg4BCAAABB0CEAAACDqWb4baEHm2RystLbW4EgAAUFee39t12eaUAFSLAwcOSJJatWplcSUAAMBXBw4cUExMzCnPYTf4Wrjdbu3cuVNNmjSRzWbz67VLS0vVqlUrbdu2jZ3mzwD30T+4j/7BffQP7uOZC/Z7aBiGDhw4oKSkJNntpx7lQwtQLex2u84999x6/Yzo6Oig/Mfpb9xH/+A++gf30T+4j2cumO/hL7X8eDAIGgAABB0CEAAACDoEoABzOp3Kzc2V0+m0upRGjfvoH9xH/+A++gf38cxxD+uOQdAAACDo0AIEAACCDgEIAAAEHQIQAAAIOgQgAAAQdAhAATRt2jQlJyfL5XIpPT1dS5cutbqkBu0///mPrrnmGiUlJclms+ndd9/1et0wDI0dO1aJiYkKDw9XZmamNmzYYE2xDVheXp4uuugiNWnSRHFxcRowYIDWr1/vdU55eblGjRqlc845R1FRUbrhhhtUXFxsUcUN0/Tp05WSkmIuMJeRkaEPP/zQfJ176LsJEybIZrPpvvvuM49xH+tm3LhxstlsXo+OHTuar3MffxkBKEDmz5+vnJwc5ebmavny5erWrZuysrK0e/duq0trsMrKytStWzdNmzat1teffvppPf/885oxY4a+/PJLRUZGKisrS+Xl5QGutGFbvHixRo0apS+++EILFy7UkSNH1LdvX5WVlZnn3H///frnP/+pN998U4sXL9bOnTt1/fXXW1h1w3PuuedqwoQJWrZsmb7++mtdccUV+s1vfqNvv/1WEvfQV1999ZVeeuklpaSkeB3nPtZdly5dtGvXLvOxZMkS8zXuYx0YCIgePXoYo0aNMp9XV1cbSUlJRl5enoVVNR6SjHfeecd87na7jYSEBOOZZ54xj+3fv99wOp3GvHnzLKiw8di9e7chyVi8eLFhGEfvW2hoqPHmm2+a56xdu9aQZBQWFlpVZqPQrFkz4+WXX+Ye+ujAgQPG+eefbyxcuNDo1auXce+99xqGwb9FX+Tm5hrdunWr9TXuY93QAhQAlZWVWrZsmTIzM81jdrtdmZmZKiwstLCyxmvz5s0qKiryuqcxMTFKT0/nnv6CkpISSVLz5s0lScuWLdORI0e87mXHjh113nnncS9Porq6Wq+//rrKysqUkZHBPfTRqFGjdPXVV3vdL4l/i77asGGDkpKS1LZtWw0dOlRbt26VxH2sKzZDDYC9e/equrpa8fHxXsfj4+O1bt06i6pq3IqKiiSp1nvqeQ01ud1u3Xfffbrkkkt0wQUXSDp6L8PCwtS0aVOvc7mXNa1evVoZGRkqLy9XVFSU3nnnHXXu3FkrV67kHtbR66+/ruXLl+urr76q8Rr/FusuPT1ds2fPVocOHbRr1y6NHz9ePXv21Jo1a7iPdUQAAoLIqFGjtGbNGq+xAqi7Dh06aOXKlSopKdFbb72l4cOHa/HixVaX1Whs27ZN9957rxYuXCiXy2V1OY1a//79zZ9TUlKUnp6u1q1b64033lB4eLiFlTUedIEFQGxsrBwOR40R+MXFxUpISLCoqsbNc9+4p3V311136f3339eiRYt07rnnmscTEhJUWVmp/fv3e53PvawpLCxM7du3V/fu3ZWXl6du3brpueee4x7W0bJly7R79279+te/VkhIiEJCQrR48WI9//zzCgkJUXx8PPfxNDVt2lS/+tWvtHHjRv491hEBKADCwsLUvXt35efnm8fcbrfy8/OVkZFhYWWNV5s2bZSQkOB1T0tLS/Xll19yT3/GMAzdddddeuedd/Tvf/9bbdq08Xq9e/fuCg0N9bqX69ev19atW7mXv8DtdquiooJ7WEd9+vTR6tWrtXLlSvORlpamoUOHmj9zH0/PwYMHtWnTJiUmJvLvsa6sHoUdLF5//XXD6XQas2fPNr777jvjtttuM5o2bWoUFRVZXVqDdeDAAWPFihXGihUrDEnG5MmTjRUrVhhbtmwxDMMwJkyYYDRt2tR47733jFWrVhm/+c1vjDZt2hiHDx+2uPKG5c477zRiYmKMgoICY9euXebj0KFD5jl33HGHcd555xn//ve/ja+//trIyMgwMjIyLKy64Xn44YeNxYsXG5s3bzZWrVplPPzww4bNZjM+/vhjwzC4h6frxFlghsF9rKsHHnjAKCgoMDZv3mx89tlnRmZmphEbG2vs3r3bMAzuY10QgALohRdeMM477zwjLCzM6NGjh/HFF19YXVKDtmjRIkNSjcfw4cMNwzg6Ff7RRx814uPjDafTafTp08dYv369tUU3QLXdQ0nGrFmzzHMOHz5sjBw50mjWrJkRERFhXHfddcauXbusK7oB+v3vf2+0bt3aCAsLM1q0aGH06dPHDD+GwT08XT8PQNzHuhk0aJCRmJhohIWFGS1btjQGDRpkbNy40Xyd+/jLbIZhGNa0PQEAAFiDMUAAACDoEIAAAEDQIQABAICgQwACAABBhwAEAACCDgEIAAAEHQIQAAAIOgQgIMj17t1b9913n9Vl1GCz2fTuu+9aXUa9GTdunFJTU60uAwhaBCAgyL399tt6/PHHzefJycmaMmVKwD7/ZEFg165dXjteA4A/hVhdAABrNW/evF6uW1lZqbCwsNN+P7tWn54jR44oNDTU6jKABo8WICDIndgF1rt3b23ZskX333+/bDabbDabed6SJUvUs2dPhYeHq1WrVrrnnntUVlZmvp6cnKzHH39cw4YNU3R0tG677TZJ0h//+Ef96le/UkREhNq2batHH31UR44ckSTNnj1b48eP1zfffGN+3uzZsyXV7AJbvXq1rrjiCoWHh+ucc87RbbfdpoMHD5qv33zzzRowYIAmTZqkxMREnXPOORo1apT5WbXxtD7NnTtXycnJiomJ0eDBg3XgwAGv7/XzFrHU1FSNGzfOfG6z2fTSSy/pf/7nfxQREaFOnTqpsLBQGzduVO/evRUZGamLL75YmzZtqlHDSy+9pFatWikiIkI33nijSkpKvF5/+eWX1alTJ7lcLnXs2FEvvvii+dqPP/4om82m+fPnq1evXnK5XHr11Ve1ZcsWXXPNNWrWrJkiIyPVpUsXffDBBye9D0AwIgABML399ts699xz9dhjj2nXrl3atWuXJGnTpk3q16+fbrjhBq1atUrz58/XkiVLdNddd3m9f9KkSerWrZtWrFihRx99VJLUpEkTzZ49W999952ee+45zZw5U88++6wkadCgQXrggQfUpUsX8/MGDRpUo66ysjJlZWWpWbNm+uqrr/Tmm2/qk08+qfH5ixYt0qZNm7Ro0SLNmTNHs2fPNgPVyWzatEnvvvuu3n//fb3//vtavHixJkyY4PO984S/lStXqmPHjvrd736n22+/XaNHj9bXX38twzBq1Ltx40a98cYb+uc//6kFCxZoxYoVGjlypPn6q6++qrFjx+rJJ5/U2rVr9dRTT+nRRx/VnDlzvK7z8MMP695779XatWuVlZWlUaNGqaKiQv/5z3+0evVqTZw4UVFRUT5/J+CsZvFmrAAs9vPduFu3bm08++yzXufccsstxm233eZ17NNPPzXsdrtx+PBh830DBgz4xc975plnjO7du5vPc3NzjW7dutU4T5LxzjvvGIZhGH/5y1+MZs2aGQcPHjRf/9e//mXY7XajqKjIMAzDGD58uNG6dWujqqrKPOe3v/2tMWjQoJPWkpuba0RERBilpaXmsT/84Q9Genq6+by2+9GtWzcjNzfXq9Y//elP5vPCwkJDkvHXv/7VPDZv3jzD5XJ5fbbD4TC2b99uHvvwww8Nu91u7trdrl0747XXXvP67Mcff9zIyMgwDMMwNm/ebEgypkyZ4nVO165djXHjxp30ewMwDMYAAfhF33zzjVatWqVXX33VPGYYhtxutzZv3qxOnTpJktLS0mq8d/78+Xr++ee1adMmHTx4UFVVVYqOjvbp89euXatu3bopMjLSPHbJJZfI7XZr/fr1io+PlyR16dJFDofDPCcxMVGrV68+5bWTk5PVpEkTr/fs3r3bp/okKSUlxfzZU0/Xrl29jpWXl6u0tNT8/uedd55atmxpnpORkWF+pyZNmmjTpk265ZZbNGLECPOcqqoqxcTEeH32z+/7PffcozvvvFMff/yxMjMzdcMNN3jVB4AuMAB1cPDgQd1+++1auXKl+fjmm2+0YcMGtWvXzjzvxIAiSYWFhRo6dKiuuuoqvf/++1qxYoXGjBmjysrKeqnz54N/bTab3G73Gb3HbrfLMAyvc2obV3TidTxjp2o79kv1eHjGN82cOdPrvq9Zs0ZffPGF17k/v++33nqrfvjhB910001avXq10tLS9MILL9Tpc4FgQQsQAC9hYWGqrq72OvbrX/9a3333ndq3b+/TtT7//HO1bt1aY8aMMY9t2bLlFz/v5zp16qTZs2errKzM/GX/2WefyW63q0OHDj7V5KsWLVqYY6EkqbS0VJs3b/bLtbdu3aqdO3cqKSlJkvTFF1+Y3yk+Pl5JSUn64YcfNHToUJ+v3apVK91xxx264447NHr0aM2cOVN33323X+oGzga0AAHwkpycrP/85z/asWOH9u7dK+noTK7PP/9cd911l1auXKkNGzbovffeqzGo9+fOP/98bd26Va+//ro2bdqk559/Xu+8806Nz9u8ebNWrlypvXv3qqKiosZ1hg4dKpfLpeHDh2vNmjVatGiR7r77bt10001md1N9ueKKKzR37lx9+umnWr16tYYPH+7VzXYmPN/pm2++0aeffqp77rlHN954o7kEwPjx45WXl6fnn39e33//vVavXq1Zs2Zp8uTJp7zufffdp48++kibN2/W8uXLtWjRIrObEsBRBCAAXh577DH9+OOPateunVq0aCHp6PiWxYsX6/vvv1fPnj114YUXauzYsWbLxclce+21uv/++3XXXXcpNTVVn3/+uTk7zOOGG25Qv379dPnll6tFixaaN29ejetEREToo48+0k8//aSLLrpIAwcOVJ8+fTR16lT/ffGTGD16tHr16qX/+Z//0dVXX60BAwZ4dfudifbt2+v666/XVVddpb59+yolJcVrmvutt96ql19+WbNmzVLXrl3Vq1cvzZ49W23atDnldaurqzVq1Ch16tRJ/fr1069+9Suv6wKQbMbPO7cBAADOcrQAAQCAoEMAAgAAQYcABAAAgg4BCAAABB0CEAAACDoEIAAAEHQIQAAAIOgQgAAAQNAhAAEAgKBDAAIAAEGHAAQAAIIOAQgAAASd/w+j5eOFOVCqCQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# make a plot\n",
    "plt.plot(estimation_err)\n",
    "plt.xlabel('iteration numbers')\n",
    "plt.ylabel('estimation error')\n",
    "plt.title('Estimation error plot')\n",
    "plt.savefig('Visualization of estimation error', dpi = 150)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 84,
   "id": "fdca2de8",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0.35560305 0.14735522 0.08134191 0.05784741 0.05121708 0.04884181\n",
      " 0.04852684 0.04815635 0.04826833 0.04813881 0.0482223  0.04816008\n",
      " 0.04820469 0.04817373 0.04819594 0.04818076 0.04819151 0.04818419\n",
      " 0.04818932 0.04818583 0.04818825 0.0481866  0.04818774 0.04818697\n",
      " 0.0481875  0.04818714 0.04818738 0.04818722 0.04818733 0.04818725\n",
      " 0.04818731 0.04818727 0.04818729 0.04818728 0.04818729 0.04818728\n",
      " 0.04818729 0.04818728 0.04818729 0.04818728 0.04818728 0.04818728\n",
      " 0.04818728 0.04818728 0.04818728 0.04818728 0.04818728 0.04818728\n",
      " 0.04818728 0.04818728 0.04818728 0.04818728 0.04818728 0.04818728\n",
      " 0.04818728 0.04818728]\n"
     ]
    }
   ],
   "source": [
    "print(estimation_err)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.18"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}