{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Quickly Simulate an RBC Model\n", "\n", "Here I demonstrate how how relatively straightforward it is to appoximate, solve, and simulate a DSGE model using `linearsolve`. In the example that follows, I describe the procedure more carefully.\n", "\n", "\\begin{align}\n", "C_t^{-\\sigma} & = \\beta E_t \\left[C_{t+1}^{-\\sigma}(\\alpha A_{t+1} K_{t+1}^{\\alpha-1} + 1 - \\delta)\\right]\\\\\n", "C_t + K_{t+1} & = A_t K_t^{\\alpha} + (1-\\delta)K_t\\\\\n", "\\log A_{t+1} & = \\rho_a \\log A_{t} + \\epsilon_{t+1}\n", "\\end{align}\n", "\n", "In the block of code that immediately follows, I input the model, solve for the steady state, compute the log-linear approximation of the equilibirum conditions, and compute some impulse responses following a shock to technology $A_t$." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtMAAAEACAYAAABxtvotAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAB1oUlEQVR4nO29d3wc1bn//z676l2yLVuSi1yxjRu2MTZVQCgGAyaU4FSSkIQEbkLaDfnldUlyv2k3JDeEJJRAyKVDQgkQehPGBuMq926vm3rvbXd+f4xW2t2ZlXakLSPpeb9e85LmmbIfjWbPPnvmc56jNE1DEARBEARBEATrOGItQBAEQRAEQRCGK5JMC4IgCIIgCMIgkWRaEARBEARBEAaJJNOCIAiCIAiCMEgkmRYEQRAEQRCEQSLJtCAIgiAIgiAMkpCSaaXU5Uqp/UqpQ0qpO022K6XUvT3bdyilFvtse0QpVamU2hVwTI5S6m2l1MGen9lD/3MEQRAEabMFQRCix4DJtFLKCfwFWAnMBdYopeYG7LYSmNmzfB2432fb/wGXm5z6TuBdTdNmAu/2rAuCIAhDQNpsQRCE6BJKz/Qy4JCmaUc0TesEngGuCdjnGuAxTWcDkKWUygPQNG0tUGty3muAR3t+fxRYPQj9giAIgj/SZguCIESRUJLpAuCEz/rJnpjVfQIZr2laGUDPz9wQtAiCIAj9I222IAhCFAklmVYmscA5yEPZRxAEQYg80mYLgiBEkbgQ9jkJTPJZnwiUDmKfQCqUUnmappX1PF6sNNtJKSUNvCAIwxZN08wS10gibbYgCMIQsNpuh9IzvQmYqZSaqpRKAG4CXg7Y52Xgiz0jxJcDDd7Hgf3wMvClnt+/BLwUbEdN02yxvPqqxk9/qi8XXPDT3t8PHIi9Nu/y05/+NOYaRJfoEl36EiOkze5ZXli3kwt++lMu+OlPmXLBBVzw05+y6nc/5Z3D78Rcm93vXdElukartsEwYDKtaVo3cDvwJrAX+IemabuVUrcqpW7t2e014AhwCHgI+Jb3eKXU08DHwGlKqZNKqa/2bPoNcIlS6iBwSc+6renu7vu9vt7V+/vu3dHXEgyXyxVrCaaILmuILmvYVVcskDa7j06fRru9vh4Atxt2V+0e9IdmuLHrvSu6rCG6rGNnbVYJxeaBpmmvoTe+vrEHfH7XgNuCHLsmSLwGuDhkpTbAN5n2Ze9eWLUK4kK6moIgCJFF2mydji5jo+3xQG1bLWXNZeSn58dAlSAIIw2ZAdECXV19vy9adHPv7x0dcOhQ9PWYcfPNN8dagimiyxqiyxp21SXElnafRnvCokWAnkwD7KrcZXJE9LHrvSu6rCG6rGNnbVZRdnnUFQyllGYXjU8+CQcPmm+bPx+uuy66egRBsDdKKbToD0CMKXZqs//8yoc8t+Vdv1hcHJx7LmQmZnLH8jtQalT9ewRBGIDBtNvSM20B355pl6vYb9v+/f7bY0VxcXGsJZhipquwsBCllCyyDPulsLAw6u8pYWA6fBrl+h5/prdnuqGjgZONJ2Ogyp/h1GbbgXDoks8eWZQKb7stLl8LBPNMA3R26r3WcwMn7RWCcuzYMdsMAhKEoaCU9G7akU6TRtvjAU0DpXSrx6TMSSZHCiMZ+ewRILzttvRMW8C3XS4sLDJs32UDC15RUVGsJZhiV12CIIxcfJPpLJ9eKG/v9J6qPXg0T5RV+WPXtlF0CULoSDJtgYFsHAcO6IMRBUEQhNjT0W3eaHuT6abOJo43HI+iIkEQRiKSTFvAt2c60DPt3X7gQPT0mDGSfW6CIAhW6PJptOt9atp6fDqjY13Vw65to+gShNCRZNoC/XmmvdjB6iEI4aKwsJB33nkn1jIEYVB0us0bbd9kem/V3phbPQRBGN5IMm0BX5uHmWca9HrT7e3R0WOGXf1kdtUlCMLIpdPH5mHmmQZo6WrhaN3RKKryx65to+gShNCRZNoCZj3Tqan+62437NsXHT2CIAhCcLpMeqbjScXt9o/F2uohCMLwRkrjhYjH49+b4XIVM3VqEfPmwSef+O+7axf0TLYVdYqLi235zd2Krp/9LDIarJ73N7/5DQ899BCVlZVMmjSJX/7yl1x77bUR0RbIz4p/FpnzFg3+vPv27WPlypX8+te/5qabbgqfqGDY5EY4ceIE3/nOd/jwww/xeDysWbOGP//5z5HRJoQV32S63uUiq7CQ8SzA4/nYb7+91XtZ5VmF0+GMtsQR0WZHk0jrskmzM+ywy2dWrD63pWc6RMwqecTFwbx5xviRI9DaGnlNQmSZPn06H374IQ0NDfz0pz/l85//PGVlZbGWFRO2bt3KpZdeyp/+9KfoJNI2we12s2rVKqZMmYLL5eLUqVOj6u8f7nSaVPMYxxycKsEv1t7dzuG6w9GSJQj9UlpaynXXXce4ceOYOnUq9957b7/7b9y4kRUrVpCVlUVeXh633347nZ2dUVJrL2L1uS3JdIgEWjwKC4uIj4eJEyEz03+bxwN790ZPmy927EkA++rqjxtuuIH8/HwcDgef+cxnmDlzJhs3boy1rKjz4YcfcvXVV/Poo4+yatWqWMuJKhs3bqS0tJS7776b1NRUkpKSOPfcc2MtSwiRLo+xznQ8KUxOmWXYd3fl7mjJ8sOubaPoig0ej4errrqKhQsXcurUKd59913uuece3nzzzaDHOJ1O/vCHP1BdXc3HH3/Mu+++y3333RdF1fYhVp/bkkyHSLCeaaXg9NON26Sqx/DnscceY9GiRWRlZZGVlcWuXbuorq6Otayo88ADD3D22Wdz4YUXxlpK1Dlx4gRTpkwhLk4cccORLrex4XYQx9RU4yPFfdX76PaEULJJECLIpk2bqKqq4q677iIhIYFp06bxta99jWeeeSboMUuWLGH58uXExcVRWFjIN77xDT744IMoqrYPsfrclmQ6RAJ7pl2uYryfr2bJtMsFzc0Rl2XArjU47aorGMeOHeNrX/saf/7zn6mpqaG+vp558+aNyiloH3jgAY4fP853v/vdWEuJOpMmTeL48eN0h1IXU7AdgZ5p0JPp/KQZJDoT/fbtcHdwsOZgNOUB9m0bRVdsOHbsGKWlpb3JYFZWFr/61a+oqKgIesyBAwdYtWoVEyZMICMjg//v//v/RmXHTyw/t6W7JUTMPkvj4/Wf+fmQnQ11dX3bNA327IFly6KjbyRhh4EaLS0tKKUYN24cAH//+9/ZFcXHDUMZKBhu0tPTeeONN7j44ou58847+c1vfhOdF7bBjbBs2TLy8vK48847+fnPf47T6WTLli2cc845sZYmhIBZT7ODeDR3HLPHzmZ7xXa/bburdjNn3JxoyRNsSKybnUmTJjF16lQOHgz9i903v/lNzjjjDJ5++mnS09O55557eO655yKo0ogdPrNi+bktPdMhEmjzKCws6u2ZVsp8IGIsrB529ZPZVVcw5s6dy/e//31WrFjB+PHj2blz56hOoLKysnj77bd5/fXX+a//+q9Yy4kaTqeTV155hUOHDjF58mQmTpzIs88+G2tZQgh4PBpdHmOdaQdxdHfDvFxjo72/ej+d7ugO3LJr2yi6YsOyZcvIyMjgf/7nf2hra8PtdrNr1y42bdoU9JimpiYyMjJIS0tj37593H///VFUbB9i+bmt7P7YWiml2UHj0aPw6KP+scJCuPlm/ffycnjgAeNxd9wBWVmR1TZcUUqNStuEMPIIdi/3xFUMJMUMu7TZnV1uLv3l//OLKRxcwF1cdBGcc66b3330O9q62/z2+fScT7Ng/IJoShWijN0/e0pLS/n+97/P+++/T0dHB6eddhq/+MUv+NSnPmW6/9q1a/n617/OyZMnOeOMM7jwwgt57733WLduXZSVDy/C2W5Lz3SI9OeZBhg/HsaONR5XUhJRWQbs6iezqy5BEEYm7Z3+jXa9y4Wjx9nY3Q1Oh9PU0lFSXhINeb3YtW0UXbEjPz+fp59+mvLycurq6tiwYUPQRBrg/PPPZ9++fTQ3N/Phhx/y3//935JIRxlJpkMkWDUPL0rBwoXGfUpKdP+0IAiCED1aO8wreUBf58jC8cZG+2jdUerb6yMpTRCEEYYk0yESrM60LwsX6km1L/X1emWPaGFXP5lddQmCMDLp6PJvtLMKC3GgN9rezpHJmZPJSc7x209DY3u5/8DESGLXtlF02YuVK1eSlpZmWH71q1/FWpqAVPMIGbNqHoGlZzMyYPp0OHTIP15SAlOnRkyaIAiCEECgzQOMPdNKKRZNWMR7R9/z26+kvITzp5yPCuwdEYQY8frrr8dagtAP0jMdIoE2j0DPtJdFi4yxPXugoyMisgzY1U9mV12CIIxM2jr9G+1Az7SXheMXovBPmuva6zjWcCziGsG+baPoEoTQkWQ6RPqrM+3L7NmQlOQf6+qC3bGZqVYQBGFU0mHSM+0MsHkAZCZlMi17mmHfbWXbIqZNEISRhdg8QqS/OtO+xMXB/PkQWBJy2zZYvDhy+rzY1U9mpmvKlCnyGFUYEUyZMiXWEoQAAnumdc+0sWca4Iy8Mzhcd9gvtqdqD1fMvILEOP+ZEsPNcGqz7UA4dMlnjwDhbbclmQ6RUDzTXs44w5hMnzgB1dXm5fNGK65ojswUBGFUETgAEYyeaS+zx84mKS6J9u723liXp4vdVbtZnBeFXhAhqshnjxBuxOYRImZ1ps1sHgB5eZCba4xHo+a0Xf1kossaossadtUlxI6ObrM600abB0CcI475ufMN54hGzWm73ruiyxqiyzp21mYVSaZDZKA6074opfdOB7J9O3g84dUlCIIgGGnvHLjOtC+LJiwyxI43HKemtSbc0gRBGGHIdOIh8uKLejLsy+rV5tU7AFpa4Pe/NybPn/sczJwZCYWCINgNmU48djzx7hYe/vAVv1geizmNq8nJgW9/239/TdO4f/P9VLZU+sXPm3weF0+7ONJyBUGwCTKdeASx4pkGSE2FWbOM8WhPLy4IgjAaCbR5QP89096a04Fsr9iOR5NHioIgBEeS6RAJtc60L2a91vv2QVtb2GQZsKsHSXRZQ3RZw666hNjR0RVanWlfFoxfgEP5fyw2djRypO5IRDSCfe9d0WUN0WUdO2uziiTTIRJqnWlfZs7Ue6h9cbth587w6RIEQRCMdJr2TJsPQPSSlpDGzByjD09qTguC0B/imQ6RRx6B48f9Y1/+MgxUpvDNN+Hjj/1jeXnwjW+EV58gCPZDPNOx47fPvcNru9b5xaZyMVM4D6Xgrrv0weKB7KvexzO7nvGLOZWTH5z9A5LjkyMpWRAEGyCe6QhipZqHL2ZVPcrKoKJi6JoEQRAEc9pNGm2vzUPTgldWmpkzk9R4/0eKbs3Nzkp5pCgIgjmSTIeIlTrTvuTmQn6+Mb4tQk8N7epBEl3WEF3WsKsuIXZ0mdSZ9k4nDsGtHk6HkwXjFxjikao5bdd7V3RZQ3RZx87arCLJdIgMtmcazHund+zQ/dOCIAhC+Ol0B6/mAcEHIYJ5zenSplIqmuWRoiAIRkJKppVSlyul9iulDiml7jTZrpRS9/Zs36GUWjzQsUqpRUqpDUqpEqXUZqXUsvD8SZEhsOEtLCwKOZmeN8+YeLe2woED4dHmS1FRUfhPGgZElzVElzXsqitWSJsNnd3+PSBZhYUhJ9Pj08aTn258pBiJ3mm73ruiyxqiyzp21maVAZNppZQT+AuwEpgLrFFKzQ3YbSUws2f5OnB/CMf+Fvi5pmmLgLt61m3LYKp5eElOhtmzjfEtW4amSRAEIRBps3X6q+YBwW0eXoLVnO729JOFC4IwKgmlZ3oZcEjTtCOapnUCzwDXBOxzDfCYprMByFJK5Q1wrAZk9PyeCZQO8W+JKIOpM+2LWc3pQ4egtnZIsgzY1YMkuqwhuqxhV10xQtpsoNMdvM409N8zDTA/dz5O5fSLtXa1sqtyV9g0gn3vXdFlDdFlHTtrs0ooyXQBcMJn/WRPLJR9+jv2DuBupdQJ4HfAj0NWHWU8HnN/s5Vketo0yM42xjdtGrwuQRAEE0Z9mw3QZeKZDmUAopfk+GROzz3dEN90ShptQRD8CSUdNKu1F1hENNg+/R37TeC7mqY9r5S6Efgb8CkzATfffDOFhYUAZGVlsWjRol6vjfebTSTX9R4Mfd3l0rfPmFGEUtbOd+aZ8Ne/6uuFhfr2558vxuGASy4Jj15vLJrXZzive2N20WP3dW/MLnrstn7PPfdQUlLS217FiFHfZkNfMl3vcuk6ejzT3ja8u3vg8y0rWMbLb74MQOEi/e9Z/+F6ssqzuOHKG8Ki1xuzyz1s93VvzC567L7ujdlFj+96UVGRLfSUlJRQX18PgKunvbDKgJO2KKVWAD/TNO2ynvUfA2ia9muffR4EijVNe7pnfT9QBBQGO1Yp1QBkaZqmKaUU0KBpWgYB2GECgNZW+G2AOzApCe40DOvpn7Y2+P3vjY8Xr7oKliwZmkZBEOxHLCZtkTZb54a7/0xVS7Vf7Ey+RSq5AHz2szBrVv/n0DSNh7Y+RGmTv6Nl4fiFXDvn2rDqFQTBHkRq0pZNwEyl1FSlVAJwE/BywD4vA1/sGSG+HL2RLRvg2FLggp7fLwIOWhEeTcy8dceOFVs+T3IyzJ9vjG/cqE8iEA6837rshuiyhuiyhl11xYhR32YDhoGCumc6dJsH6B+qywqMRUt2Ve6ipbNlyBrBvveu6LKG6LKOnbVZZUCbh6Zp3Uqp24E3ASfwiKZpu5VSt/ZsfwB4DbgCOAS0Al/u79ieU38N+KNSKg5oRx9RbkvMkmmnc3DnWrbMOGFLRYU+VflAU5MLgiAMhLTZOl0mVTesDED0cvq403nr8Fu0drX2xtyam61lWzlvynlD1ikIwvBnQJtHrLHDI8OKCrj/fv/YuHFw222DO9/f/gYnTvjH5s2D668f3PkEQbAnsbB5xBo7tNkAl/3i13R0d/jFzuFHxJMMWLPXvXPkHdYdX+cXy0zM5DvLv4NDydxngjCSiJTNY9QzlBrTZiwzmepgzx5oahr8OQVBEIQ+zOpBW6nm4cvS/KWogLGZDR0N7K/eP2h9giCMHCSZDgGzRvfQoeJBn2/uXEhL8495POGZxMWuHiTRZQ3RZQ276hJiQ7fbg9vjX8+03uVC0efPC9XmAZCVlMWsMcbRiptKh14mz673ruiyhuiyjp21WUWS6RAIp2fae6zZ48XNm83rWQuCIAih095pbLQVTr/eZSvJNGA6EPFI3RGqWqos6xMEYWQhnukQ2LcPnnnGP3baabBmzeDP2dgI99yj90j7cv31un9aEIThj3imY0N1QyvX/8G/nmkcyZzLj3rXzzkHLrkk9HNqmsZfNv2F6lb/cnvLCpZxxcwrhqRXEAT7IJ7pCGFm87Ay+6EZGRkwZ44xvnHj0M4rCIIw2mntMDbajoDiVVZ7ppVSnJl/piFeUl5iGOgoCMLoQpLpEDBrdPfvLx7yec0GIh4/DuXlgz+nXT1IossaossadtUlxIYOE5tHfcDMZlYGIHpZOGEhCc4Ev1inu5PtFdutn6wHu967ossaoss6dtZmFUmmQyDcnmkvkyfD+PHGuPROC4IgDJ6OLrMa0/6NttWeaYCkuCQWjl9oiG86tYlYW1sEQYgd4pkOgY8+grfe8o8tXw6XXz70c2/eDP/+t38sPh6+9z19xkRBEIYv4pmODZv2n+SHTz/sF0snnyU+88zMnQs33mj93JUtldy36T5D/EsLv8TU7KnWTygIgq0Qz3SECHedaV8WLICkJP9YVxeUlITn/IIgCKMN855p/0Z7MDYPgNzUXAqzCg3xjafkkaIgjFYkmQ4Bs0Z3587isJw7IQEWLTLGN26EwXTu2NWDJLqsIbqsYVddQmwwS6ZrXAf91gdj8/BiViZvX/U+GtobLJ/Lrveu6LKG6LKOnbVZRZLpEIiUZ9rLmcYB4tTVwX6ZXEsQBMEybZ1m1TyG7pn2MnvsbDISM/xiGhqfnPpk8CcVBGHYIp7pEHj1VdgUMNHVFVeYV+MYLE88AYcO+ccmToSvfhXUqHJcCsLIQTzTseH5dTv50zvP+8Vymcdcru9dnzABbr118K+x9tha3jv6nl8swZnAd5d/l+R4GfAiCMMV8UxHiEjUmQ5k+XJj7ORJvVSeIAiCEDqdJt3OQ60zHciSvCXEO/x92J3uTjaXbh7aiQVBGHZIMh0CZo3utm3FYX2N6dP1npJA1q2zdh67epBElzVElzXsqkuIDe0mNo8q116/9aEm06kJqZyRd4YhvuHkBrrcoY9utOu9K7qsIbqsY2dtVpFkOgTMGt1w90wrBeeea4wfPAgVFeF9LUEQhJGMec+0v2d6sNU8fDl70tk4lP/HaEtXCyXlJUM/uSAIwwbxTIfA44/D4cP+sc99DmbODO/reDzwpz/pgw99mT8frrsuvK8lCELkEc90bPjjSx/w4rb3/WKTOY9pXNy7npgIP/7x0F/rhb0vsKNih18sOymb/zjrPwyJtiAI9kc80xEiknWmfXE44JxzjPFdu4wJtiAIgmCOWc+0M0x1pgM5Z5Kx0a5rr2N35e7wvIAgCLZHkukQMEumN20qjshrLVwIqan+MU3TZ2EMBbt6kESXNUSXNeyqS4gNZsl0hcu/99jj0ZehMj5tPLPGzDLE159YH9IU43a9d0WXNUSXdeyszSqSTIeAWQ9GOOtM+xIfb17ZY9s2aGmJzGsKgiCMJDpMGm2FsdEe6iBEL2a90+XN5RyuO2yytyAIIw3xTIfAvfdCba1/7D/+A8aMiczrtbfDH/4AHR3+8fPPh4suisxrCoIQfsQzHRt+/Oi/+PhoiV/sNK4hD//qGz/8ofFJ4GDQNI1Htj3CicYTfvHCrEJuXnTz0F9AEISoIZ7pCBGNOtO+JCXB0qXG+MaNxgRbEARB8CeUOtMQvp5ppRTnTjaWY3LVuzjRcMLkCEEQRhKSTIeAWYP70UfFEX3N5cuNVpL2dtiypf/j7OpBEl3WEF3WsKsuITZ0dht7QMpcWw2xcCXTALPGzCI3NdcQX39ifb/H2fXeFV3WEF3WsbM2q0gyHQLRqDMdSHo6LFpkjH/8cXg/AARBEEYaXe7QeqbDVdED9N5pM+/0vup9VLVUhe+FBEGwHeKZHgBNg5//3Bi/6y69lF0kqamBP/9Z1+DL1VfD4sWRfW1BEIaOeKZjw1f+9HeO1Bzziy3iZrIo9IvdcgtMnBi+13V73Nz7yb00dDT4v/aERayevTp8LyQIQsQQz3QEcLuNMYcj8ok06AMc58wxxtevD09JJ0EQhJFIp8l03vHOyHmmvTgdTs6edLYhvqNiBw3tDSZHCIIwEpBkegCCTdgSLa+P2RTjNTWwf7/5/nb1IIkua4gua9hVlxAbuk1sHlWnNhpi4bR5eDkj7wxS4lP8Yh7Nw4aTG0z3t+u9K7qsIbqsY2dtVpFkegCiXckjkPx8mDbNGF+71mj/EARBEMx7phMTI1dn2pcEZwLLCpYZ4ptLN9PSKZMFCMJIRDzTA1BXB3/8o38sKwvuuCN6Go4cgcceM8ZvvBHmzo2eDkEQrCGe6dhw1a9/T1NHk1/sxknfo/JEhl/s05+GBQvC//qtXa384eM/0OXxT+qXT1zO5TMuD/8LCoIQNsQzHQFiUckjkKlToaDAGH//ffFOC4IgBNLtMTbcKUmR90z3vlZ8CkvzjZMFbC7dTGNHY2ReVBCEmCHJ9AAEs3lE0+ujFFx4oTFeVQU7d/rH7OpBEl3WEF3WsKsuITZ0mdg8XIeN9Z4jWWb03MnnkuBM8H89TzcfuD7wi9n13hVd1hBd1rGzNqtIMj0AwQYgRpvp02HKFGO8uNi84oggCMJoxOPRTHumk0w805EYgOglNSGVswrOMsS3lW+jtq02ci8sCELUEc/0ABw+DI8/7h+bOhW+9KXoazl2DP7+d2N81Srz6ccFQYgt4pmOPu2d3Vz+q1/4xZwOJz85979Yu9Z/3wsvhAsuiJyWtq42/vjJH2nvbveLLxi/gE/P+XTkXlgQhEEjnukIYAfPtJcpU2DGDGN87VqZFVEQBAGgtd2kxrQjzrTdjnS7mRyfbFp3emfFTipbKiP74oIgRA1Jpgcg1nWmA7noImOssRE2bdJ/t6sHSXRZQ3RZw666hOjT3mVstOOc8ezYUWyIR9Lm4WX5xOWkxqf6xTQ03j/6PmDfe1d0WUN0WcfO2qwiyfQAxLrOdCD5+eazIq5bB52d0dcjCIJgJzo6jcl0rHqmQa87fe5k4+xbe6v3UtpUGnkBgiBEnJCSaaXU5Uqp/UqpQ0qpO022K6XUvT3bdyilFodyrFLqP3q27VZK/Xbof074CWbzKCoqiroWLxdeqFf48KWlBTZsiK2u/hBd1hBd1rCrrlgxmtvs1g7zqcTPOqvIEI+WPe7MgjPJSMwwxN87+p5t713RZQ3RZR07a7PKgMm0UsoJ/AVYCcwF1iilAqcKWQnM7Fm+Dtw/0LFKqQuBa4AFmqadDvwuHH9QuLFLNQ9fcnNh/nxj/KOPoK0t+noEQbAPo73N7jCxecQ74017pqNh8wCIc8Rx/pTzDfFDtYc43nA8OiIEQYgYofRMLwMOaZp2RNO0TuAZ9AbVl2uAxzSdDUCWUipvgGO/CfxG07QOAE3TbDkaww51ps0oKgJHwH+vvR3uu684FnIGJNbXKxiiyxqia1gwqtvstk7znult24oN8WgO3D5jwhlkJ2Ub4n/5x1+wY1Utu76nRJc17KoL7K3NKqEk0wXACZ/1kz2xUPbp79hZwHlKqU+UUh8opc60Ijxa2LFnGiAnB844wxjfu1e3fAiCMGoZ1W22Wc90gjMep7HMdNR6pkEvz1dUWGSIV7RUcKTuSPSECIIQdkJJps1q7QV+jQ62T3/HxgHZwHLgh8A/lAp0AsceO3qmvZx/PoYPiIkTi/jww9jo6Q87XC8zRJc1RNewYFS32eY2jzjOPbfIEI92SdH54+czLmWcX6xwUSHvHn3Xdr3Tdn1PiS5r2FUX2FubVUKpS3ESmOSzPhEIHIIcbJ+Efo49CbzQU91/o1LKA4wFqgIF3HzzzRQWFgKQlZXFokWLev8J3scEkVrfsqUYlwsKe3oUXK5iduygt2GO9Ov3t56ZCYmJxezZ46/vxAlYsULfHkt9si7ro239nnvuoaSkpLe9ihGjus3evOEj6l0usnpev97loqJF6/VMu1z6/oWFRXR3R/cecSgHqaWpbHJtonCRrs9V4sKFi/Mmn8eccXNifg/LuqyPtvWSkhLq6+sBcLlcDIYBZ0BUSsUBB4CLgVPAJuCzmqbt9tnnSuB24ArgLOBeTdOW9XesUupWIF/TtLuUUrOAd4HJgVNnxXo2rZdfhq1b/WNXXQVNTcW9/4xY0twMf/xj3+NKl6uYwsIiFiyAT9togq3iYntcr0BElzVElzViMQPiaG+zn/2ghPvf/5df7KzChZyVkc2OHUV+8XHj4LbboqcNQNM0/rrlr5Q1lwF6Ml24qJAxyWP41pnfwulwRldQEOz6nhJd1rCrLrCvtojMgKhpWjd6o/smsBf4h7dh7WlcAV4DjgCHgIeAb/V3bM8xjwDTlFK70Ae5fCmmLXAQ7FZnOpC0NFi+3BjfsQNOnoy+HkEQYstob7PNPdOxqzMdiFKKi6YaZ9+qaavhk1OfRF+QIAhDZsCe6VgT816OZ/VBfb7ccAOcfnps9JjR3g733gutrf7xggK45RZjTWpBEKJDLHqmY02s2+y/vv4xT33ypl/sotPO4o4rV/K//+u/b3o6fP/7URTXg6ZpPLr9UVz1Lr94ojORb5/1bVITUs0PFAQh4kSkZ3q0Y9dqHr4kJZlPM37qlN5DLQiCMFroNGm0k+JjW2c6EKUUl8+4HIX/53WHu4P3jr4XG1GCIAwaSaYHwK51pgNZvBjGj+8bXOPlnXegoyM2mnyx2/XyIrqsIboEu2OWTCfExfHRR8WGeCxsHl4mpE1gSf4SXCUuv/jWsq2UNZXFRpQPdn1PiS5r2FUX2FubVSSZHoBgpfHshsMBl19ujDc1wbp10dcjCIIQCzpMekAS4uNM60x3d0MsnY4XFl5IgjPBL6ah8cahN2xXKk8QhOCIZ3oAHngAysv9Y9/4BuTlxUbPQJh5vOPi9BHr2dmx0SQIoxXxTEefnz/1Gu8f2OgX+8KKlXz1srP4xS+MHSQ/+UlsrXsbTm7gjUNvGOI3zL2B03NtNDhHEEYJ4pmOAHav5hHIpZca9XV3w1tvxUaPIAhCNOl0Gx8nJvY0inap6OHLmflnMjZlrCH+1uG36HLHyNQtCIIlJJkegGADEO3q9dm+vZgVK4zxvXvh6NHo6/Fi1+sluqwhugS702nSA5KUEE9xcbGtBiF6+XDth1w+w+jRa+ho4KMTH8VAkY5d31Oiyxp21QX21mYVSaYHYLh4pn057zy95FMgb7wBHk/09QiCIEQLs57phHi90Tazc8S6ZxpgRs4MZubMNMTXHV9HQ3tDDBQJgmAF8UwPwK9+BZ2d/rE779TL0dmZ7dvhxReN8SuvhDPPjL4eQRiNiGc6+nzz/ifYW3HIL/bjVZ/lsqWz+MtfoCpg8vNvfQtyc6MoMAjVrdXct+k+PJp/j8f83PlcN/e6GKkShNGHeKYjwHCoM23GggUwcaIx/v770NYWfT2CIAjRoMukZzo5QW+07Wjz8DI2ZSxnFZxliO+s3MnxhuMxUCQIQqhIMt0PbrfRFqGUXobOrl4fry6lzEvltbZCLKTb/XrZDdFlDbvqEqKPWTKdGB9HcXGxLW0evvfuBYUXkBKfYtgnFqXy7PqeEl3WsKsusLc2q0gy3Q/B/NLDZXruiRNh4UJjfONGKIv9nACCIAhhp7Pb2NWcmGDfah6+JMUlcfHUiw3x0qZSNpVuioEiQRBCQTzT/dDSAnff7R9LSYH//M+YyBkUjY3w5z8bfd8TJsDXvobpRAaCIIQH8UxHn+vv/hPVLTV+sQe/ejunTRrLU0/BgQP++990E8yeHUWBA+DRPPx1y18pb/af4CDBmcBtZ95GZlJmjJQJwuhgMO22zetSxJbhVmPajIwMOP98fVpxX8rL4eOP4dxzY6Mr1rR3t1PXVkdzZzNt3W20dbX1/mztaqWtu4327na6Pd24PW48mge35vb7XdM0nA4nTuUkzhFn+D3eEU9yfDJJcUkkx/X89FlPTUglPSGdpLgk1HB53CEINqfbxOaRlGDvah6+OJSDlTNW8veSv/vFO92dvHrwVdbMWyPthSDYjGGWGkaX/sriFRcXU1RUFFU9oWCma8UK2LkTKioC94U5c2DMmNjoijRNHU1UtlRS115HXVsd9e31vb+3deujMF0lLgoXFQ7+RdxD1xnniCM9IZ20hDTSE9NJT0jn0NZDXHThRWQlZZGdnE1yXLItPkCH030vjE46TSY6SUqI66kzXWTYFusBiGb37pSsKSzJW8KWsi1+8QM1B9hdtZt5ufNiossOiC5r2FUX2FubVSSZ7ofhWskjEKcTrr4aHn4YfJ++dnfDK6/Al740fHzgwWjqaKK0qZTSplLKmssobSqlubM51rJCotvTrSf57XW9MVeZi5o9fY+qE52JZCdn68l1UjY5yTmMTRnL2JSxpCWk2SLRFgQ70O0JXs1jOPRMe7lk+iUcqDlAU2eTX/z1g68zLXua6UBFQRBig3im++HECfjb3/xjBQW613g48uaburUjkKuvhsWLo69nsGiaRnVrNYdqD3G0/uiwSpwjQaIzsTex9i65qblkJ2fjUDLGOJaIZzr6XPTz/zbUan7rJ/9FQryTN96ADRv897/sMkxnjbUDe6v28uzuZw3xRRMWsXr26ugLEoRRgHimw8xI6Zn2cuGF+rTi9fX+8bfegpkzzWdNtAutXa0cqTvC4drDHK47TGNHY6wl2YYOdwenmk5xqumUXzzeEU9uai7j08YzPnU849PGk5uaKz1awoils8ttSKQdykFCvD7S2s51ps2YM24Oc8bOYW/1Xr94SXkJ83PnMz1neoyUCYLgiyTT/TBSPNNeEhLgqqvg8cf94+3t8PrrcOONsdEVjJrWGnZW7uRgzUFKm0rRCF9vl1M5yUrKonRnKUvOXkJyXDLJ8cmGn/GOeJwOJw7lwKmcfr8rpXB73Lg1d+9ARd/fO9wdtHe3097dTluXPqDRO7CxrauNps4mmjqa6PIYP82H7OUGujxdpkl2RmIGeWl55KXn9f5MT0gPySoyHO97YfTQ3mlstOMceqOt15kuMmyPtc1joHv3iplXcLT+KO3d7X7xfx/4N98885skOBNioitWiC5r2FUX2FubVSSZ7oeRUM0jkOnT9drT27f7x/fs0Xut58yJjS4vLZ0t7KrcxY6KHYYk0CpO5SQ3NZdxqeN6vcbZydlkJ2WTnpiOQzkobiumaG5ReMQPAk3T6HB30NzZTFNHU2+C/eHJD5mUM4m6dn3gpJkPdLA0djTS2NHI/pr9vbHU+FTy0/PJS88jPz2fgvQC0hNt/KhCEExo7TA22t5kGuxfZ9qM9MR0Lp1+KS/vf9kvXtdeR7GrmEunXxojZYIgeBHPdD/s2AEvvOAfmz8frrsuJnLCRmsr/OUveh1tX9LT4bbbICkpunq63F3sr9nPjoodHKo9ZHhMGwrexNk3IcxNzfX7IB2uaJpGc2dzbyWSuvY6alprqG6tprq12rRnOxxkJGYwMWMiBekFTMyYSF56XsR6wUYq4pmOLscq6vnS/ff4xbKSMvnXnd8F4JNP9Kdwvpx5Jlx5ZZQEDhJN03h0+6O46l1+cYXilsW3UJBREBthgjACEc90mBmJPdOgTzyzciU895x/vKkJ3n5bt4JEg9q2Wj4+8TE7KnbQ4e6wdKxTOZmUOYnp2dOZlj2N8WnjR0TibIZSSi+Zl5jO5MzJfts0TaOxo7E3sa5qraKqpYqKlgrDY2GrNHY0sqdqD3uq9ug6UOSm5jIpcxKTMiYxKXMS2UnZUklEsA2mNg9nX7swnKp5+KKU4qpZV3H/5vv9nlJpaLy8/2W+vuTrOB0yA5cgxIqRmX2EiZHmmfbl9NP1nvfA2cC2bNF73wsLI6frVOMp1p9Yz96qvZZ80GOSxzA9ZzozcmZQmFUYll7S4f5/VEqRmZRJZlKm32Akb5Jd0VJBZUslFc0VVLRUUN1aPaief9A/uD9Z/wkViyrYXLoZ0O0hvsl1fnp+TL7U2PX/KESXtk5jD0i8j2c6J6fIsD3WyXSo9+6YlDEUFRbxzhH/GbgqWipYf2I95085Pya6oo3osoZddYG9tVlFkul+GGnVPHxRSn+06XIZpxp/8UX45jfDa/fQNI2DNQdZf2K94VFlf0xIm8CC8QuYM3YO2cnZ4RM0wvFNsmeNmdUb7/Z0U9lSqdfjbiqjrLmMiuYK3NrgZp9p6WphX/U+9lXvA/QnBvnp+UzJmsLkzMlMzpxMUlyUfUPCqKXDpGc6Ia6v0R5u1TwCWTFxBbsqdxmmGi92FTMjZwb56fkxUiYIoxvxTPdDcbG++HL++XDRRbFQExk2boTXXjPG583TveFDfYLv9rjZVbmL9SfWU9lSGdIxGYkZzM+dz4LxCxifNn5oAoQBcXvcVLZU9k52c6rxFBUtFYPuwfbFaw3xJtdTMqeMqoGN4pmOLmt3HOWuFx71i80YW8jDt98MwMGD8OST/sfMmAGf/3yUBIaB0qZSHtrykOGp3pjkMXxj6TdkXIMgDBHxTIeZ/mweI4WlS2HXLjh+3D++a5dee3rhwsGf+2DNQd48/CbVrdUD7pvgTGDuuLksGL+AwqxCmWwkijgdTr1MXnoei/P02Xu63F2UNZdxqvEUJxtPcqrpFPXt9ZbPraFR0aJbTDae2gjoH/pTsqZQmFXIlMwpZCZlhvPPEUYxpjaPuP6reQynnmmA/PR8zpl8DuuOr/OL17TV8OahN7nqtCgNehEEoRfJWPqhP5tHcWCXtU2wqsvhgE9/GhITjdtefRVqa61rqGqp4okdT/Dkzid7E2lXict03/SEdC6ZdgnfW/E9Vs9ezbTsaVFNpEfK/zHcxDvjmZw5mRWTVnDD6Tdwx/I7+MHZP2BawzTOmXQOkzMnD9obXdNWw9ayrbyw9wX+sOEP/HHDH/nXvn9RUl4yqIQdYn+9BHvQYdJoJzj1RluvM208xg6eaatcWHgheWl5hviWsi3srdprcoR17PqeEl3WsKsusLc2q4ywftbwMlKreQSSlQWrVsHzz/vHOzv10oBf+YqedA9Ea1crH7g+YFPppgEtAuNSxnH2pLOZP37+iK3CMdJIS0hjcuZkiqYXAbo9pKy5jBMNJzjReILjDccHNa17XXsddeV1lJSXAJCdlE1hViFTs6dSmFVIRmJGGP8KYSTT2WWWTA/vOtNmOB1Orpt7HQ9uftBQGvPl/S9TkFEg7xtBiCLime6HF180Tm6yejUsWhQLNZHnhRf0Ch+BXHCBPhV5MNweN5tLN1PsKqatu63f15iSOYVzJp/DzJyZUlJthKFpGnXtdRxvOM6x+mMcbzhOTVvNkM87JnlMb3I9NWsqqQmpYVAbHcQzHV2eeHcLD3/4il/snOmL+eUXrgaguhr+/Gf/Y3Jy4NvfjpbC8LKldAuvHHjFEJ+WPY0vLPiCtLGCMAjEMx1mRoNn2pcrrtC90/X1/vG1a/WZEydPNh5zouEEL+1/aUBf9JTMKXxq2qeYlDkpfIIFW6GUIic5h5zkHBZNWARAc2dzb3J9rOEYFc0VlqeFr2mroaathi1lWwDITc1lWvY0pmZNZUrWFKkWIvRibvMY/nWmg7E4bzEHaw/2VtPxcqTuCB+f/JizJ50dI2WCMLoQz3Q/9GfzsKvXZyi6kpJ0/3RgZ4am6b3W7T5zgLg9bt47+h6PbHuk30Q6KymLG0+/kcL6Qlsm0iPx/xhJrOpKS0hj7ri5rJy5kluX3sp/nvOfrJm3hrMnnU1+ej4K6z1nlS2VbDi5gad3Pc3/rPsfHtryEL9/6vccqTtCl3uYjSYTwkqHSaOdENdXZ9qONo+hvNeVUlx92tWkJxgr5Lx75F3KmspioiuSiC5r2FUX2FubVUZwP+vQGcl1poMxebJu6wi8x+vr9RJ6n/60PsDwhb0vUNYcvKFOcCZw/pTzWT5xOXGOOCpVaGXxhJFNcnwyp409jdPGngZAe3c7JxpO4Kp34ap3UdpUaqnnWkPjVNMpXJUumrY3EeeIY1LGJKZlT2Na9jTy0vOkMswootOk0U6MHzl1ps1IiU9h9ezVPL7jcb+4W3Pz/N7n+caSbxDvHOEfXIIQY8Qz3Q9/+xucOOEf+/KXYcqUmMiJGh4P/P3vxr9dQ2Nm0UZcjrf9prT1RaFYNGERF029aFTVExbCQ3t3O8cbjuOqd3G07ijlzeWWbSG+JMUlUZhV2Jtcj0keE1UfqXimo8tvn3uH13b5l4y7fsnF3H7VeYDetv33f/sfoxTcddfQa+rHmjcPvcnHJz82xJfmL2XVrFUxUCQIwxPxTIeZ0eaZ9uItl/fAA9DRocc6aGQf/2Ldh0dYsgRSUozHFaQXsGrWKvLSjSWbBCEUkuKSmDVmVu+sjW1dbRxrOMbRuqMcrT8a8sQ/Xtq72/1maMxIzOhNrKdlTyMtIS3sf4MQO9r7sXmA3rY5HHpS7UXT9HWnMxoKI8fF0y7maP1Rw+yIm0s3MzVrKqfnnh4jZYIw8pHnn/0wGupMByM7Wx+QCFDJbjZxH3Ucwe3WJ3TxvTYO5eDCwgv56uKvBk2kR/r1CjeiSyc5PpnZY2ezcuZKvnXmt/jB2T/g+rnXsyRvCdlJfdPLB6tjHkhjRyMl5SW8sPcFfvfR77hv0328cegNDtQcoNPdGaG/QogWXSaNdlJ8X51psJ/VI1zvqThHHNfNuc601OhL+1+y/EVU2iBriC7r2FmbVUZBP+vgGS11poMxf77G89veZY/L/7Fpayvs3w9z58K41LFcO/taCjIKYqRSGE2kJaQxL3ce83LnAVDfXs/RuqO87HqZ1IRUy3WuK1sqewc0OpWTiRkTmZY9jek508lPzxe/9TCj023mmfZvtOPj9Rr6vsR6EGK4GJc6jsumX8arB1/1i3e6O3l217N8bcnXpPqNIEQA8Uz3w+9+B80Bn83f+x5kjIJa+J3uTl7Y+wK7yvexZQu0mZSPvnLRMu646hIZ3CLYAk3TqGqt4kjdEY7WHcVV76LD3THo83n91tOzpzMtexo5yTmW/dbimY4u3334H2w7uccvdsel17P67Hm963/4AzQ0+B/3ne/oT+NGApqm8c89/2RP1R7DtlljZrFm3hqpPy0I/RAxz7RS6nLgj4ATeFjTtN8EbFc9268AWoGbNU3bGuKxPwDuBsZpmtZ/seIoMxqreQA0tDfw9K6nKW8uJy4O5s2DrVvB7da3J5DObK6hdfsMjs/Xa1ALQqxRSpGbmktuai7LJy7H7XFT2lTKkbojHKk7wsnGk7g1d8jnC/RbZyZmMj1nem+NaztPHjNa2+yBqnnAyKs1HYhSimtOu4aqliqqWqv8th2oOcAHxz6gqLAoNuIEYYQy4DNMpZQT+AuwEpgLrFFKzQ3YbSUws2f5OnB/KMcqpSYBlwDHh/yXRIDRVmca4GTjSR7a+pDfIJbUVDhNr2TGWOZwJt8khxloGjz3HNTVRV5XpBBd1hhOupwOJ5MyJ3FB4QV8+Ywv86Nzf8Tn5n+OFRNXMD51vOXXaOhoYGvZVp7b8xx3f3Q3D25+kLcPv83h2sO2qm89qtvsfmwe/XmmY5lMR+I9lRiXyGfmfYZEZ6Lx9VzFHKg5EBNd4UB0WcOuusDe2qwSSs/0MuCQpmlHAJRSzwDXAL7PkK4BHut5trdBKZWllMoDCgc49g/AfwIvheFvCSseT19PrC8j2TO9s2InL+1/ybTsXW4uTHdcROOu8/wm2mhrg2efha9+dXT02gvDlwRnAjPHzGTmmJmAPjvj0bqjHKk7wuG6wzR2NFo6X1lzGWXNZaw/sZ44RxyTMyf3WkImpE2I5aP0UdlmA3SafKlJTvBvmOw2ADFSjE0Zy7VzruWZXc8Ytr2w9wW+tvhrjEkZEwNlgjDyCCU1LAB8Kw6fBM4KYZ+C/o5VSl0NnNI0bbsd/VvBEmmv1KKioqjqCZXB6NI0jfdd77P22FrT7fGOeK6dcy2zx8zliVY4csR/e3k5/PvfsHp18FqtI+l6RQPRZY3B6EpLSGP++PnMHz8fTdOoaavRE+vaw5b91t2e7l47CegTaUzLnmZZU5gYlW029N8z7b1H7GbziOR7avbY2Zw/5XxD297e3c6zu5/llsW3kOBMiLquoSC6rGFXXWBvbVYJJZk2azUDR5cE28c0rpRKAX4CXBrC68eE0VLJw6N5eHHvi+ys3Gm6PT0hnTXz15Cfng/A9dfDgw8aB/Bs3w75+XBW4Ee2IAwDlFKMTRnL2JSxLCtYhkfzcKrxFIfrDvf6rT2aZ+AT9dDa1cquyl0RVNwvo7LNBujsNjbcgdU87GbziDRFhUWUNZVxsPagX7yypZKX9r3E9XOvlwGJgjBEQkkPTwKTfNYnAqUh7pMQJD4dmAp4ezgmAluVUss0TfOvOA/cfPPNFBYWApCVlcWiRYt6v9F4PTfhXl+8WF93ufT1wsIi4uP9PT5FRUURe/3Brt9zzz0hXx+P5uH/Pfr/OFp/lMJFhfrf21Ovt3BRIfnp+UysmciBLQfIL9KT6Y0bi5k0CVpaiuju9r8+b74JR48WM2GC8fVGwvWK5ro3Zhc9o/F6TcqcxOFth5nGND537udw1bv41xv/orSplOw5eukH3/cLwIbnNlB+qJysCVnEkFHZZhcVFdHt6abe5dJft+f1d2z5hCpXeu8++/cXc/y43maB3oatXw+zZkVen9l6pN9Taz9Yy9jusdSk11DbVut3z+6u2k35P8qZP36+tEE2bIPCsW7X6+V7rWKtp6SkhPr6egBcPe2HVQYsjaeUigMOABcDp4BNwGc1Tdvts8+VwO3oI8PPAu7VNG1ZKMf2HO8ClpqNDI9VmaWaGvjTn/xj2dl6CSXQ/wHef4adCFWX2+Pmhb0vsLtqt+n208edzurZq4OWvSspgX/9yxhPTtb902PHDk5XtBFd1hBdOo0djb2WkCN1R2jpajHd7+cX/jzqpfFGa5sNsOrXv6O5w7+e6VP/8X3yx6T33iMvvAA7dvgfd+21sHBhFIX6EK17t6K5goe3PkyXx7/3XqG44fQbmDvOf4yqvNetIbqsY1dtgymNF1KdaaXUFcA96KWSHtE07ZdKqVsBNE17oKfM0p+By9HLLH1Z07TNwY41Ob8LmzXMFRVw//3+sXHj4Lbboi4l7Lg9bp7f+7xpHVKAosIiLphywYCP/l57DTZuNMazsuCWWyBNZmoWRgGaplHZUtlrCTlWf6w3YYlFMg2js80GuOwXv6aj29/r/sL3fkRORnLv+ssv66U+fbnqKliyJBoKY8uuyl08t+c5QzzOEccXF36RyZmTY6BKEOxFxOpMa5r2GvBaQOwBn981wDTNNDvWZJ/CUHREk5FaY3qgRPrKmVdyZsGZIZ3rssv0wYfHA4pk1dfDU0/BzTdDgvnYFkEYMSilGJ82nvFp4zl70tl0e7o50XCidzBiLBiNbTZgWokoJWngOtMjsZqHGfNy53Gq8RQfn/zYL97t6ebpnU/z1cVfZWzK2CBHC4IQDJkrNwgDDUD09frYif50uT1untvzXNBEetWsVSEn0gBOJ3zmM5CTY9xWWgr//KdeYnAgXbFEdFlDdA1MnCOOqdlTuXjaxbGWMqrodntwewLKMClIiHMCffeI3QYgRvvevWT6JZw25jRDvK27jSd3PElzZ3NMdIWK6LKGXXWBvbVZRZLpIJg1rsO5mofb4+afe/7J3uq9ptuvmnUVS/OXWj5vaip87nOQkmLcdvAgvPoq2HzGekEQRgDtncZGO84Rh8Ph/7R2tNSZDoZDObhu7nUUpBcYttW11/HUzqfodHfGQJkgDF9C8kzHklj57/btg2cCat2fdhqsWRN1KUPGm0h7p0UO5OrTrmZx3uIhvcbJk/B//2f+JeTii+G884Z0ekEYlgzGezfciVWbXd3QyvV/+K1fLDk+mdd/8iO/2Lp18M47/seecw5cckmkFdqLls4WHt76MHXtxilsZ42ZxU3zbsKhpL9NGH0Mpt2Wd0oQRkqdaU3TeGn/S6aJtEKFJZEGmDgRrrvOfNKWd9/V61ALgiBEimA904aYzWwesSI1IZXPL/g8KfHGx4oHag7w6oFXsXtnmyDYBUmmgzDQAES7en0Cda09tpYdFTsM+4UzkfYyZw6sXGm+7Z57ig0zJ9qB4fJ/tAuiS7ArbR3GHpAEn9Ke/XmmY2nziOW9OyZlDGvmrTH90vH868+z7vi6GKjqH7u+10WXdeyszSqSTAdhJHimd1bs5H3X+4a4QnHN7Gs4I++MsL/msmVw9tnGuKbBs8/qAxMFQRDCTUeXsdGOdxobbbtNJx5rJmVO4ro516EwPlZ89+i7bC3banKUIAi+iGc6CB99BG+95R9bvhwuvzzqUgbF8YbjPFryKG7NbdgW7h7pQDQNnn8edpnMppycDF/8IuTlRezlBcE2iGc6emzaf5IfPv2wX2xiVj5P3PF1v9ju3XqlIV/mzoUbb4y0QnvzyclPeP3Q64a4QrF69moWTojRrDaCEGXEMx1GhnOd6dq2Wp7Z9YxpIn3e5PMimkiD7ptevRqmTDFua2uDxx/XJ8URBEEIF2Y90wkmM7jazeZhF86aeBYrJq4wxDU0/rXvX+ys2BkDVYIwPJBkOgjDtc70m++8yVM7n6K1q9Ww7fRxp3PR1IuioiMuDm66SZ81EsDlKu7d1toKjz0GVVVRkdIvdv0/ii5r2FWXED3Mk+m+Rtt7j9jN5mGne/fS6ZcyL3ceAK4SV29cQ+PFfS8GnaMgmtjpevkiuqxjZ21WkWQ6CMPRM+32uHnf9T7VrYYZfpmYMZHVs1cPOEV4OPFaOsaMMW5raYFHH4Vqo1RBEATLtHUae0DMPNNSzSM4SimunX0ts8fONmzzaB6e2/Nc0BKrgjCaEc90EF59FTZt8o9dcYU+wM6OaJrGy/tfZlv5NsO2rKQsbll8C2kJaTFQBo2Neg3q2lrjtvR0+PKXzWdRFIThjnimo8fz63byp3ee94stmTSP33/1er9YWRk8+KD/sRMmwK23Rlrh8KHb082zu57lYO1BwzancvKZeZ9h1phZMVAmCJFnMO22zftaY8dwqzO97vg600Q60ZnIZ+d/duiJtKZBR4duem5vN/5sb9e7d7xLV1fv7xnd3dzi8LBuP7S1YpgScfNGOOc8B6mZcfoc5d4lLq7vZ0JC35KY6L+elKR3gycn2/ufJAhCxOg06V5OMGkPpGd6YOIccXxm3md4eufTHK477LfNrbl5dtezrJm/hhk5M2KkUBDshWQeQRjI5lFcXExRUVHU9PTH4drDvHv0XUD3uRUuKgT0aWNvPP1GclNz+z+BpkFDA9TXQ1OT3pXc2Oj/e3MzeDyD1rjR5eK8yYWUlOg5uR+NsPPfsGiRnhcPibi4vsTau6SmBl2KP/mEooui4yO3gp3uL19El2BX2k1sHr7JtPcesdsARLveu+vWruOm827iqZ1PcbT+qN82t+bmmV3P8Nn5n2Va9rSo6rLr9RJd1rGzNqtIMh2E4VLNo7mzmRf2vmC67cqZVzI9Z3pfoL0damp0o3JNTd/vtbVR+TRJTtYT5m3boLPTf1t7O5SUwIIFkGKckCt0urv1LwFNTaHtf+yY/sLp6eZLRoa+DDnLFwQhkpj3TBsbbbsNQLQz8c541sxfw5M7nuRYwzG/bd2ebp7a+RQ3nn6jWD6EUY94poPw+ONw2P/pFp/7HMycGXUpQdE0jSd2PGF4DAdw7rglfCr5dH2WFO9SXx99kSa0tur5a2BCDfoH3YIFeh5rKxITITNTXzIy+n7PytKX9HRwyHhewR/xTEePP728lue3vucXu3L+efzwuov9Yu3t8Jvf+B+bmAg//nGkFQ5fOro7eGLHE5xoPGHY5lAOrjntGqlDLYwYxDMdRoZDz/RHJz7SE2lNI622mazyejKqmyhsTeDMJA3UllhLNCUlBRYu1BPqwA7xri49Pn++nqPaho4OqKzUFzMcDv/k2rtkZ+ujK1NT9QLcgiBEhHaTp2tJJo223Wwew4HEuEQ+v+DzPLb9MU41nfLb5tE8vLjvRdq621g+cXmMFApCbJFkOgi29kxrGqWunex96zFOP1VDVkU98R264J0VTaxeellkSuAFDvbz/p6UpC8JCfpFMlmKP/mEIp95xlOBWVX6TIltbaA0Dw5PN8rjxqG52dPt5tIF3Uyd7Nb/GZ2dekLb2dm3eNe9AyHb2iz7uotdLooKC4d+bTweqKvTFzPi4/sS6+zsvt9zcvSkO6BX265eMtEl2JWBBiB67xGnU/9e69t57vHoSyweLtn13g3UlRiXyBcWfoHHtz9uSKgB3jj0Bq1drVxYeGFES7AOl+tlF+yqC+ytzSqSTAfBdtU82tvhwAE4fJiuQwc4vv99ZnS3G3abkjWFBGeC9fMnJuoFoTMz/b3C3t/T0/VkebAcOwaTJ/uFxk2GGwp1S42ZA+WxnXD1VDhjaYivoWl6cu1baaS1VS9qbbY0Nw/+77FKV1fwnm2Hwz+5zsmBU6f0xDwzU+wjghACHSaNtlk1D6X0tjxw9+7uoTVxo4GkuCS+uPCLPLv7WY7UHTFsX3tsLa1drVwx8wocStotYfQgnukg3HuvsS7yf/yH+QQkEaO1Ffbvhz174MgRcLvRNI291XupbDEmZZMyJvkPOAxEKT1RGzMGxo71/xlDG0JjIzzxRHAHxSWXwDnnROjFu7r0pNpbvSRwaWjQt7mNU7NHBaez738WuIh1xPaIZzp6/PjRf/Hx0RK/2DeKrmFN0RmGff/nf/Tv2r788If6W0oYmG5PNy/sfSHojIhzx83l03M+TZxD+uuE4Yd4psNIzHqmm5th3z49gXa5DLaF8uZy00Q6PSGdqdlT+wJK6QlXfn7fMmGCLbteMjL0iVuefBJOnjRuf/tt/XvFpz4VgdzRa7/Izg6+j6bpPdnexLqhoW+pr9eXVuP07WHB7dbnXTebez0pyf9LkXfJzpZ628KoI9Q606C/7QOTaanoETpxjjiun3s9rx54lS1lxrE5e6r20N7dzmdO/wyJcYkxUCgI0UU+cYMQVc+0261bOLZs0UuIBOnVaelsCToj1Zz8BThmzKK4ro6iq66CvDzdumETBrpe3qnH//EPOHTIuH39ej1nXb06vANBQ/o/KgVpafpSUGC+T0eHnlzX1fUl2F4P9SBKD4bk5W5v1799BH4DUUpPqMeN80+yx40bcok/u3rc7KpLiB6d3cb3WGK80TMN9pq4xa737kC6HMrBqlmrSIlP4cPjHxq2H6k7wiPbHmHN/DVkJWVFTVesEF3WsbM2q0gyHYSoVPOoq4OtW/XCywP4dz2ahz1Ve/Boek+1x6FoHJdBfV42Z523hpRFn9ItAcXFEI4BdTEgIQHWrIEXX4Rdu4zbd+/W89I1a/TebFuRmAi5ufoSiLdn2ze59v6sqQl/r7am6eeurdVtQr6kpfUl1r4/09PFMiIMa7rcxkY7OcG80TZry6Wih3WUUlw87WJS4lN48/Cbhu0VLRU8tOUhPjPvM0zOnGxyBkEYGYhn2gRNg5//3Bi/664wjAVzu/UEZ8sW3Qcd4t92qPYQh7srqZ48jppJY2jIzcQd72TRhEWsnr16iKLshabB66/Dxo3m29PS4KabYOLE6OqKGO3tfclvTU3fz5oa47PoSJGY6J9ge5esLEmyh4B4pqPHV/70d47U+E8s8ovrb+bceYWGfR96SB/j68stt4ygNiUGbC/fzkv7X+rt8PHFqZysmrWKM/KM/nVBsBvimQ4TZmPNHI4hJtIdHbBpE2zYYK2KRGYmVVPG8Up3Iw3jpvslNmOSx3DFzCuGIMqeKAUrV+r1qIuLjdubm+H//g+uvlqf4GXYk5TU52sPpLW1L7EOXML5XLqjw9wyEh9v7MkeN04fFClVRgQb0ek2sXkkmH/ESa3p8LNwwkKS4pJ4bs9zdHn8L6Zbc/PS/peobKnkkumXSKUPYcQhybQJoVg8Qvb6tLbCJ5/oS7uxlJ0pWVkwbx7MmYN7wnj+ufWvNLT4+xqcysn1c683lMGzqwfJqi6loKhIH1v30kvG/0l3N7zwgl4B5OKLB995avvrlZKiL5Mm+e+gabpHu7q6b/FODx/qVOqh0NUFZWX6go+X2+k09mJ7k2ynM3yvHyJ2/T8K0aPbzObh03D73iN2mlLcrvfuYHSdNvY0vnLGV3h619M0djQatn988mOqW6u5bu51JMUNbvzGSLpe0cCuusDe2qwiybQJYank0dwMH3+s90abzZsdiMMBs2fDkiUwbVpvdvjx8XWm1TsunnYxeel5FkUNP+bP1/OzZ54xzxHXrdMLXXz607Yabxl5lOqbZXHGDP9tHR19CXZVVd/vtbWWJ7UJitsNFRX64ovDoX8D8ibXubl9SbZUGBEiiJln2ncAoi92GoA40shLz+PrS77OM7ue4WSjsTzTwdqD/G3r31gzfw05yTkxUCgI4Uc80ybU1cEf/+gfy8qCO+4I4eCGBr30xNatobXO2dl6Ar1okW4G9tXRVsd9m+4zPDLLS8vja0u+NqoelTU26gl1aan59nHj4IYbzMf/CT243XpC7Ztke39G+hm3w6En1IFJ9pgxIzrJFs909Ljq17+nqcP/G/cTt32PieOMo5Wfe844yPnTnx4htjGb0O3p5pX9r7C9Yrvp9qS4JFbPXs3ssbOjrEwQ+kc802FiUD3THR16N+lHHw08wYdSMGcOLF0KU6eaehQ0TePfB/5tSKQViqtOu2pUJdLQV4v6pZfMK31UVemDiq64Qv9eImPmTHA6+5LZOXP64l7LiG+C7V1CtSYNhMfT10O+d29f3DuRkDe59i5jx47oJFsIP90eY+dFUhDPtJ1sHiOVOEccq2evJjc1l3eOvIOG/xes9u52ntn1DMsnLueSaZfgdETfHiYI4WJ0ZWQhMlCNadC9PoCeiJSUwJ/+BB9+2H8i7XTC4sX6VIo33uhn5whkV+UuDtcdNsTPmngW+ekmA9UCddmMcOiKj4frroOLLjLf3tWlJ9svvhiasyZcuiJBVHV5LSMzZ8KKFfrIzq9+FX70I/j+9/UC4CtXwtKlFLe3h3eaOE3Tvd5798LatfD88/DAA/DLX+rvqWeegXffhR07dN92kB50u/4fhegxUDLte4/YyeZh13s3HLqUUpwz+RzWzF9jGN/jZcPJDTyy7RHq2uqipisSiC7r2FmbVaTrx4SQa0wfPw5vvBHce+AlLk63cpx9NmRmDvj6bV1tvHHoDUM8MzGTi6YGySRHCUrB+efrnZcvvGCeW+3Yof9LbrgBxo+PvsYRg1J6/en0dP2LH+hWpKIivW62txe7srKvJ9tKpZr+8CbZNTX6jKC+mrKy/K0i48ZJKYZRjsej0WVSzSM5MfSeabmFIsesMbO4ZfEtPL3zaerajUnzqaZTPLjlQa457RrmjJtjcgZBsDfimTbh8GF4/HH/2NSp8KUv9aw0NOhzXJv5DXxJSIBly2D5coMfuj9e3v8yW8u2GuJr5q3htLGnhXyekU5lJfzzn+YzbYP+HWblSv1hgNg+okRbW19i7Ztkh7PCSDC8SXZgGb8hzvo4FMQzHR3aO7u5/Fe/8Is5HU7eveu/TPd/7z39QYgvF14IF1wQKYUCQGtXKy/te4n9NfuD7rOsYBmXTr+UOIf09QmxQTzTYaJfm8fmzfDmm/13Y8TF6Y/Lzz5bnyfbAsfqj5km0nPGzpFEOoDcXPja1+C113SnTSDd3fDKK3D0KKxaFdOcavSQnAyTJ+uLL75Jtm+yHc4k2zuN+8GD/vH0dKMfe9y48NpVhJjS2m5sj/tLxuxk8xhNpMSncNO8m9hwcgNvH3nbdIKXjac2cqLhBNfPvZ4xKWNioFIQrCOeaROC2jxOnoR//xu6uih2ucwPnjsXbrtNL35sMZHu9nTzyoFXDPFEZyIrZ64M6Rx29SBFSldCAqxerS/BpnvftQvuuw8OHYqerqEy4nR5k+wlS+Dyy+ELX9D92Hfeqfuzr75a/wI6Y8ag5ooP+n4EPWE/ckSv9f7vf+sz/tx9N/z2t/D3v+vfuDZs0B9JNTSEPCupYB/au4yNdrzTv0HwvXftZPMYce/1AVBKsWLSCr5yxlfISsoy3aesuYwHNj/AhpMbCHzKMdqu11Cxqy6wtzarSM+0CUGrefh6NwOZMEFPEgoLB/2664+vp7q12hC/eNrFZCRaTzBGE4sWQUGBbvuoNJblprERnnhCt3xcdtkoq0ltZ5KS9AlpAiel6egw9mJXV+s9z+GitRWOHdMXXxISjLM+jh0rsz7amI5Ok2RaeqZtzcSMiXxjyTd4af9L7Ks2frZ2ebp449Ab7K3ayzWzr5Ga1IKtEc+0CZs36x1YvixeDFe3Petf1gv0R8UXXQRnnDGkD9qa1hru33y/YUT6xIyJfOWMr4y6UniDpasLXn9dL/MdjMxMuOaavjF1wjCis9NoF6mq0pPsSLcTTqeeUPsm2N4lwbxSgXimo8NuVyW3/d99frHctLH84we3m+6/bZte+ceXhQvh2msjpVAIhqZpbDy1kbcOv4VbM6+GFe+I55Lpl3Bm/pkoGQAjRBjxTIeJ4DaPGuOGz3zG6A8dBG8dfsuQSDuUg6tmjb6a0kMhPl53DEydqn8h6ugw7tPQAI89ppf5vuQS6aUeViQk6I8gCgr8411dxhrZVVX6JDXhSuzc7r7zBn6pzsz0T669ixAVOkKwefhtkzrTtkEpxVkTz2JS5iSe2/MctW21hn26PF28dvA19lTt4ZrTriE7OTsGSgUhOCFlaUqpy5VS+5VSh5RSd5psV0qpe3u271BKLR7oWKXU3UqpfT37v6iUygrLXxQGTG0eTk3/YO6h16MZhg/M4w3HTUc3nz3pbManWavtZlcPUrR1zZ8P3/xm/73PmzfD979fzJEjUZMVMvJ/tEbx+vWQl6dPYXfxxXDTTXo995/8RL8RbrhBL+l3+un6yFVnmCeIaGjQPdeffAKvvgqPPgq//314X8MCo63NNk+m/fuKpM60NaKtKz89n1uX3sqygmVB93HVu/jhQz9k06lNBi91rJH/o3XsrM0qA/ZMK6WcwF+AS4CTwCal1Muapu3x2W0lMLNnOQu4HzhrgGPfBn6saVq3Uup/gB8DPwrfnzZ4zBrV5M4G44bkZMuDDAPRNI23D79tiGcmZnLBFKnTNBSysvRxblu2wFtvmU/k0tys91LPnw+XXqoXfhBGEHFxerHxwILjHg/U1fX1NPtOrR7qjD82ZTS22W2dxh6QhH56ps2SaakzHXsSnAlcMfMK5oydw0v7X6K+vd6wT7enm1cPvkpJeQlXzrqy30nMBCFahGLzWAYc0jTtCIBS6hngGsC3Yb4GeKzHKLdBKZWllMoDCoMdq2naWz7HbwCuH+ofEy7MGtWkFn+LR1FhIYwZM+QCxvuq93Gi8YQhfuHUC/t9TBmMoqKiIemJFLHSpZRu55g+HV5+WS+T50thoa5r5044cEDvvFy2LPwdl1aR/6M1LOtyOPT375gxMHt2X1zT9NGqvpYR7zToLS1h1RxBRl2bHUrPtO89Yiebx4h5T4WRqdlT+ebSb/L2kbfZXLrZb1vhokJAn+jloS0PsTR/KRdNvYjk+KF1bA0V+T9ax87arBJKMl0A+GZ7J9F7MgbapyDEYwG+AjwbgpaoYNaoJjab+KXHDK0Gpkfz8O7Rdw3x8anjWTB+wZDOLfiTna3Pir1pkz7fjtkXpo4OvYR4SQlccQVMmRJ1mUKsUUr3P2dm6t/AfGlt7UusfXuyozH40Rqjrs3uNEmmE8y6n3uwk81DMCcxLpFVs1Yxd9xcXtr3Eg0dDYZ9NDQ2lW5iT9UeLpl+CQvHL5QBikJMCCWZNrszAz85gu0z4LFKqZ8A3cCTwQTcfPPNFPaUnMvKymLRokW932i8nptwrm/bBqCvu1z69qScNn27Tz3boosuGtLrbS3byuaP9G/d3m/brhIXM6bO6B10aPX899xzT8Svz2DWvbFY6lEKWluLmT8famqKOHYMNmy4hwkTFvX2ULtcxbhcUFFRxMKFkJRUTHLy6LxeZuuj+v5KSaG4x2BfdOmlfdvdbormzYPqau65/35Kdu2iMC1NT75jw6hrs7eU9BWRr+9poxMKF/rt7z2muLiYujoIbOPHjYucvv7WR/V7KoT149uPc7r7dDomdrClbAsbntvAhBkT/D4zAVq6WthatpWs8ixyknNG7fUaLveX77WKtZ6SkhLqe8quuvqbs6AfBiyNp5RaAfxM07TLetZ/DKBp2q999nkQKNY07eme9f3oLVVhf8cqpb4E3ApcrGma6SdPLMosPf+8/tjfly8nPMGUzr4Gu9jlouiHP9QHNA2CTncn935yL82dzX7xqVlT+eLCLw7623VxcXHvTWIn7KZL0/Qe6AceKCYvryjofklJcP75uvWjn46usGO36+VFdFlA01AOR9RL443GNvuxdzbzyDr/eqbnTV/C//vCVb3rvvdIbS3ce6//ObKz4TvfibRSI7a8d7GnLle9i3ueuYes2VlB91EoluYv5YLCC0hLSIuaNjteL7CvLrCvtsGUxgslmY4DDgAXA6eATcBnNU3b7bPPlcDtwBXojwTv1TRtWX/HKqUuB/4XuEDTtKp+Xj/qDfOzJuWkv9H+R/KS6vyDt96qT9YyCD5wfcD7rvcN8a8v+boMqIgibW3w7rv6IMX+brPMTL2c+IIFQ7bJC6OIWNSZHo1t9l9f/5inPnnTL3bRaWdx1xrzmWMbG+F//9c/lp6uT8op2Bu3x80npz6h2FVMpzv4YOEEZwLnTDqHFZNWkOA0rwMvCGZEpM50z8jt24E3ASfwSE/DemvP9geA19Ab5UNAK/Dl/o7tOfWfgUTg7Z5e2A2apt1qRXykCPTOOTzdxLfWQ1LAjjmDm5GppbOF9SfWG+LzcudJIh1lkpNh1Sp9zp1XX4XSUvP9GhrgxRfh44/hU5/S7bSSVAt2ZDS22Z0mhucks1GGPdhpOnHBGk6Hk7Mnnc283Hm8dfgtdlXuMt2v093J+6732VS6iaLCIs6YcAZOhzPKaoXRQkh1pjVNe03TtFmapk3XNO2XPbEHehplNJ3berbP1zRtc3/H9sRnaJo2SdO0RT2LLRplMDaqSW11OJV/T0txZWXQWc8G4oNjHxi+UTuUg4umXjSo8/np8vEg2Qm76yoogFtu0RPr/qodlpfr05I//jiUlUVel90QXcOD0dZmmyXTgQMQfe8ROw1AtOu9a3ddGYkZXD/3er648IuMTQk+30NzZzP/PvBv7tt0H3ur9kasPrXdr5cdsbM2q8jUeiYENqopbTXGmcIzMgZ17tq2WkOpH4Cl+UvJSR5cT7cQHhwOvYze7bfr08f31/N85Ag8+CA895xe1EEQhNjRYdKtPJhqHvYqyiKEwrTsaXxz6Tf51LRPkegMPp1tTVsNz+5+lke2PcLh2sO2m/RFGN4M6JmONbHw3z3wgN4D6WXS8fXcmP22/4QeS5fq3ZgW+efuf7K7ardfLMGZwHfO+g6pCamDVCxEgspKeOcdvf50fygFc+boAxUHaaEXRiix8EzHmli02T9/6jXeP7DRL/aFFSv56mVmVf10fvELY8fJT35ibgERhgctnS2sPbaWzaWbcWvufvedmDGR86ecz8ycmVJOT/AjIp7p0UhgJ0dKWw2OMQE7DaLG9KnGU4ZEGuCcSedIIm1DcnPhs5+FY8f0GRRPnTLfT9Ngzx59mT1bT6rzxfouCFGj0230aCQOUH4nLs6YTHd3SzI9nElNSGXlzJWcNfEs3jv6XlA/NcDJxpM8tfMp8tLyOH/K+cweO1uSamHQiM3DBMOs4a01hhnxig8dwgqapvH2EeO04WkJaayYtMKqxKDY1YM0nHVNmaL7qW+8ceAxp/v2wV//Ck8+CSdORFZXLBBdgh3pNLF5JCX4Z8WB94hdphS36707nHXlJOdw/dzr+drirzE1a2q/+5Y1l/Hs7me5f/P97KrchUfzRExXLLCrLrC3NqtIz7QJpj3TQ/RMu+pduOpdhnhRYZGU7RkGKAVz58Jpp8HWrfDBB9DcHHz/gwf1ZcoUWL5cP85wDwmCEBbMeqYT4vv/eLPTlOJCZCjIKOCLC7/IodpDvHPkHSpaKoLuW9lSyXN7niMnOYflE5ezaMIi+WwWQkY80yb86lfQ2VNsw9ndwXnrfs255/r0ZDgcurkusLu6Hx7f/jiH6w77xcamjOVbZ36rd7ZDYfjQ1aUn1evWQVPTwPtnZ8NZZ+kl+BKDj5ERRhjimY4O33rgSfaUH/SL/XjVZ7ls6aygx/zlL8bBw9/6lm7vEkYemqaxr3ofa4+tpax54FJMSXFJLMlbwrKCZWQmZUZBoWAXxDMdJnx7J1LaaoCAXsXsbEuJdGlTqSGRBriw8EJJpIcp8fF6crxkCWzbpifVDQ3B96+rgzfegPff1yuFnHUWZGVFTa4gjGg6u43+jOSE/s3PUmt6dKGUYs64OcweO5uDtQf5wPUBp5qCDIQB2rvbWX9iPR+f/Jg5Y+ewYtIKJmZMjKJiYTghmVwAbjd4fCxTya16Mu03LmHMGEten3XH1xliY5LHMGfcnEGqDI5dPUgjVVdcHJx5Jnz723D11fr3rP7o6NAnfvnjH/WZNg8dMi/HNVKvV6Swqy4hOnSZDUCMD15nGuxTa9qu9+5I1aWUYtaYWdyy+Ba+sOALTM6c3O/+Hs3D7qrdPLz1YR7e+jAl5SV0uY3fukbq9YokdtZmFemZDiBYjenAZDpUqlur2Vu11xA/Z/I50is9gnA69R7nRYtg505Yv14vrRcMTdOnrN+7V++hXrxYt4D4lV8UBCEkzHqmExMGruYRiHimRw9KKabnTGda9jRc9S7WHV9n+gTZl5ONJznZeJI3Dr3BgvELWJy3mAlpUg9VEM+0gZYWuPvuvvU5e55nYt1OzjnHZ6dVq/Q60yHw8v6X2Vq21S+WnpDOd5Z/hziHfJcZqWiaPrHLhg36QMRQcDhg5kzdOjJjhgxYHAmIZzo6XH/3n6huqfGLPfjV2zltUvCZ8Z5+Gvbv94/ddJNe3lIYnVS2VLLh5AZ2VOyg2xPaN6uC9AKW5C9hXu48GbA4QhDPdBgIqZJHiD3TjR2NbC/fboivmLRCEukRjlIwfbq+VFfrSfX27f17Mj0e/cN9/369WMzChbBgAYwbFz3dgjAc6TaxeSRJz7RgkdzUXK4+7Wounnoxm0s3s6l0E82d/ZRtAk41neLU/lO8cegN5uXOY8H4BUzJnCI1q0cZ0vcVgF9jqmkkt5on06F4fT4+8bFhFqbkuGSW5C0Zss5g2NWDNJp1jR2rP8z47nfhU58Krarijh3FfPihXnHggQfgo4+gsTHiUgdkNP8fBfvSaeJhHYxnWupM9zGadaUmpHJB4QXcsfwOrp19bUhWjgNbDrC1bCv/V/J/3LPhHt4+/DYVzcFL8UULu/4fwd7arCLdowH4JtMJXS3EuTv8k+n4+JCMra1drWwp22KILytYRmKc1EYbjaSkwLnnwooVuvVjy5bgAxB9KS/Xl7ffhsJCvbd6zhxISoqKbEGwPWaP5FMSrVfzkJ5pwZc4RxwLJyxkwfgFnGo6xZbSLeyq3EWXp/9vXQ0dDaw/sZ71J9aTm5rLgvELmJ87X0rsjWDEMx3AiRPwt7/pv2fWH+OMkr+Tnq77WAGYMAFuvXXA8xS7iil2FfvF4h3xfHfFd0mJTwmvaGHY0tCgl9bbutVaz7PTqVtI5szRJ4RJkVvKlohnOjpc9PP/Nsxc99ZP/ouE+OAlTN94Q7df+XLZZfqXXUEIRkd3Bzsrd7KldEtI9ap9mZgxkTlj5zBn3BxykgeYTleIGeKZDgNmNab9SkqH4JfudHfyyclPDPEl+UskkRb8yMyEoiI4/3y9l3rrVjhwwL88oxlut77fgQP6QMUpU/QZGmfPloogwuiis8ttSKQdytFvIg32sXkIw4vEuESW5i9laf5SSptK2Vq2lR0VO+h0dw54rLcayNtH3mZC2oTexHpcyjjxWA9zxDMdgG9j6q0x7Wfz6Emm+/P6bCndQlt3m1/MoRysmBj5Lg+7epBEV/84HDBrll5N4Hvfg9zcYgoKQjvW44GjR+HVV+H3v9efrKxbBxUVA1tIrGKX6xWIXXUJkae90+jNMBvgHXiP2MXmYdd7V3QNTH56PqtmreIHZ/+AybWTmZkzM+SSt+XN5bzvep/7Nt3Hnzf+mbcPv82x+mO4Pe6BD7aAna5XIHbWZhXpmQ5gwNkPB+iZdnvcfHzyY0N8wfgF4pcSQiItTe9lLiqCmhrYsUOvXV1bG9rxJ07oyzvv6IMdZ87Ul2nTIEEqNwkjjFCTacM+Us1DCBMJzgSm5UyjaEERLZ0t7K7azY6KHZxsPBnS8TVtNb0e66S4JKZlT2Nmzkxm5MwgPVEeNQ4HxDMdwI4d8MIL+u9nbvwLqa1V5ObqyQ0At9wCE4NPKbqtbBsv7X/JL6ZQ3LbsNsamBK95Kgj9oWlw6pR+f+7ZA839V2syxenU7SAzZuh+69zcgMmIhLAjnunIc6yini/df49fLCspk3/d+d1+j9u4EV57zT925plw5ZVhFiiMWmrbatlZsZNdlbuoaq0a1Dny0vKYOWYm07OnMzFjIk5H//YlYeiIZzoMeG0eSvOQ3KZ3BYbaM+3RPKZTh88eO1sSaWFIKKV/h5s4ES6/HE6e7JtBsb4+tHO43fpEMkeO6OupqXp1kGnTYOpUfSp0Sa6F4YZpz7RTeqaF2JOTnMMFhRdwQeEFVLVUsa96H3ur91LaVBryOcqayyhrLmPtsbXEO+KZkjWFqVlTmZo9lQlpE2QmZZsg/4UAvI1pYnsDjp4a0b3JdEoKJCcD5l6ffdX7qGmrMcTPnXxuJKSaYlcPkuiyRn+6HA6YPFmvPPCd78A3vqEPYLQ6uUtLC+zeDa+8AvfeC/fcAy+9pE8uU1dn7rcejtdLGNm0dRpHDcaH4Jm2SzJt13tXdFljIF3jUsdx3pTz+PqSr3PH8ju4fMblTM6cjCL0HowuTxeHag/x9pG3+euWv3L3+rt5dtezbDy1kfLmcsNA3FB0xRI7a7OK9EwH4G1MU3yS4t5qHgP4pT868ZEhNi17GgUZIY4kEwSLKAV5efpy0UX6bIsHDuh1rI8dG7gqiC/eMn3btunr6el60u5dxo+PzN8gCEOhw6RnOiGu/xrTYD4AUap5CNEgKymL5ROXs3zicpo7mzlYc5CDtQc5XHuYDndHyOdp625jb/Ve9lbvBSDRmcikzElMzpzM5MzJFKRL7hEtxDMdQHGxvhSc/ISZh14HdJ/p1KnAokWwerXpcWVNZTy45UFD/IsLv8i07GmRkisIQeno0C0dBw/qS1PT0M6XkNBnNSko0Je0tPBoHamIZzryrN1xlLteeNQvNmNsIQ/ffnO/xx08CE8+6R+bMQM+//kwCxSEEHF73JxoPNGbXFe2VA7pfA7lID89n4kZEylIL6Ago4DspGwpwzcA4pkOA2Y90702j356pjeVbjLEJqRNYGrW1HDKE4SQSUzUJ3WZM0e3bFRU6LWsjxyB48etP9Lu7PT3XINeJ9ubWBcU6D3kiTLBpxBFTG0eg/RMS8+0EEucDieFWYUUZhVyyfRLqG+v51DtIY7UHcFV76K1q9XS+Tyap7e2tZfkuGQKMgooSC8gPz2fgowC0hKkV2SoSDIdgDfB8NaYBvNkuri4mKKiIgDau9vZWbHTcK4z88+M+jdAX112QnRZI9y6lNIn75wwQZ/SvLtbH8R49KieHJ86FZolxOUqprCwT1dDg77s2dP3Ojk5fa81YYKeYEe6B9uu/0ch8nSYfCs0s3kE3iN2qjNtx3tXdFkjErqykrJ6J4jRNI2KlgqO1h3lSN0RjjUcC2miGFeJi8JFhb3rbd1tHKo9xKHaQ72x9IR0JqRNYELaBPLS85iQNiEqPdh2/V8OBkmmA/D2TFjpmd5evp0uj3+XRqIzkfnj50dCoiAMmbg4vZJHYSFceKFuCTl+HFwu/WdpqV79wyqaptfGrqnRBzd6SUvTE+vx4/WSfLm5MHaseUIjCFbo7DJJpqWahzDCUEr1JrwrJq3A7XFT2lTK0fqjHG84zomGE5b81r40dTbRVNvEwdqDvbFEZyLj08YzIW0Cuam55KbmMi5lHMnxyeH6k0YU4pkO4MUXYee2bs5b+0sU+uvOnq0nAvzkJ4ZPf03TuG/TfYYakssKlnHFzCuiJVsQwkpXl55QHz+uLydOQHt7eF9DKb0cnze5zs3Vv6+OGTNyJpcRz3TkeeLdLTz84St+sXOmL+aXX7i63+NqauBPf/KP5eTAt78dboWCEHk8mofKlkqONxzvXRo7GsP+OukJ6b3JdW5qLmNTxjImZQwp8Slhf61YIZ7pMNDVBclttb2JNPT0TGdmmnajHWs4ZlqMfWn+0kjKFISIEh+vD7ydMkVf1zSorNTtIN6lstJatZBANE2f1bG2Fvbt89+Wmakn1WPH6os3yc7ICKj7Lox6TG0e0jMtjDIcytHbc72sYBkA9e31nGw8yanGU5xqOkVZU5nhKbpVmjqbaOps4nDdYb94SnyKnlgnj+lNsMemjCUrKSukGUmHOyP/L7RId7e/Xxp6PrwDLB5er8+mU8aBh4VZheSm5kZSZlDs6kESXdawmy6ldIvG3r3FXHVVEaAPSCwv90+w6+rC83peL7bvYEfQy1RmZek9iDk5es92Tg7s2VPMlVcWiW1kFNJp6pk2rzPt+56ySzJtt/e6F9FlDTvqykrKomRDCZcVXQb09V6fajxFaVMpp5pOUdlSaVqf2iqtXa29PeK+KBQZiRnkJOeQnZxNTnJO77Ljkx1cevGlQ35tOyDJdADd3f5+aTBPpgGaO5t76zv6Ir3SwmggIaGvBrWX9na9akh5ub6UlUFV1eD812a43X2ebF9cLigp0b3ZWVn6kpnZ93tWlt6rLZVGRh4dJiU4EkP4ViV1poXRhm/v9RKWANDt6aaqpYry5nLKmssoby6nvLk8pMGNoaCh0dDRQENHA0frj/ptc+10sTVxK1lJWYYlMzGTzKRMkuOSh0UpP/FMB/C3v0HK2y+RV76tN7ZoEWTddDksX+6379pja3nv6Ht+sbSENL67/Ls4HU4EQdAT4KoqPcmurOxbGhqiryUxUU+yMzL6lsxMfYIa75KcHL5p1cUzHXl++9w7vLZrnV/s+iUXc/tV5/V7nMcD//3f/jGl4K67wvf/F4ThiKZp1LXXUd5cTmVLJZUtlVS1VFHTVhOWXmwrxDniyEjMICMxg8zEzN7fMxIzSE9MJz0hndSE1LBOqy6e6TAQas+0R/OwpXSL4fgzJpwhibQg+OB09pXJ86WjQ0+yvcl1dbW+NDSYT2UeDjo6+l6vP72+yXVamr6kphp/N7MKCNGl3aQ72czmEYjDoS++vn9N07/8yf9VGM0opXqtGHPHze2Nd3u6qWmt6U2wq1urqW6tpratFrcWpsePAXR7uqltq6W2rTa4XhRpCWm9yXVaQlrvkpqQqv+M138mOBMi0tMtTUYAZp5ppxNDMv3Ey0/QkOnftaZQLMlfEmmJ/WJH3xaILquMBl2JiX0zKvrS1aUPSqyu1u0c3p+1tdDWZn6uwPrXQ8Hthvp6fRmIxEQ9sU5NhZQU408h8nSZGJ2TTDwcZvdufLz+BcuX7u7oJtOj4b0eTkSXNcKpK84Rx/i08YxPG+8X92ge6tvre5Prmtaa3iS7qTP41LuBNbAHi4bWOzByIOId8aQmpJISn0JqfM/PgPXBIMl0AJ7WdhK6WvxiKs6hmy592F+9n/hM/wZ75piZZCX57ycIgjXi4/XBjuPHG7e1temDHL1VQLy/V1Xpj+aj7Vrr6NCX2uCdJkKE6XQbk+nE+NA+2uLizJNpQRBCx6EcvT3Zs8bM8tvW5e6irr2O2rZa6trqenuZ69rrOM7xIGeMHF2eLurb66lvrw/reSWZDiC+scYQU2Ny/Opx1bXVkTA9AQ3/T+4z88+MuL6BsOM3YxBdVhFd5iQn60t+vn/8y18uwu2Gxsa+XmXfpbFRX8I1EFKwD53dodk8zO5dO0wpHuv3VDBElzVElznxzvjemtSBaMs0mjube5Nb36Who4HGjsawDYSMNJJMB2CWTDvH+Vs8tpRtMSTS2UnZzMiZEVFtgiAEx+nUS+VlZ5tv1zRoadGT6oaGvgS7sRGamvqWzuHRdgs9mJXGS0oIrUaiXaYUF4TRiFJK9zknpjMpc5Jhu6ZptHe309jRSGNHY2+C3djRSFOHbuto6miirTuI/y+KhJRMK6UuB/4IOIGHNU37TcB21bP9CqAVuFnTtK39HauUygGeBQoBF3CjpmlhqlI7eMySaYdPMt3t6WZr2VaD12dJ/hJblG8ZDb6tcCK6rDGcdSnVN4AwsGfbl44O/+S6uVlPwn1/NjdDa+vQJq2JJKOpze4K0eZhdo/Yodb0cH5PxQLRZQ276oKBtSmlSI5PJjk+2eDT9qXb092bXDd3NtPUof9s6WrRf3a29K53eyLzBh8wmVZKOYG/AJcAJ4FNSqmXNU3b47PbSmBmz3IWcD9w1gDH3gm8q2nab5RSd/as/yh8f5p1PB5IajHpmc7tS6b3VO2htauV8kPlvcm0Uzk5Y8IZ0ZLZLyUlJbZ844gua4gua4RTV2Kivowd2/9+mqYn1K2tepJt9jMWjKY2G6DTbfRlJJv0TJvdI3aoNT0a3lPhRHRZw666IHza4hxxZCdnk50c5LFkD5qm0eHuoKWzhdauVlq6en4GrA9KQwj7LAMOaZp2BEAp9QxwDeDbMF8DPNZTXHSDUipLKZWH3oMR7NhrgKKe4x8FignSMP/9K7dZ+qOGQlxtLRX0PedVCl6uzqZl534AypvLAWhvbu/dZ+64uaQmpEZNY3/Uh1KCIAaILmuILmvEQpdSfZU8xo2L+sv3R8zb7NseeCq8f1E/lDdWGWJmPdNm94hZz/Rbb0W3EsuHH9aTG5sJc/tFdFlDdFkn+toUkNSzGCfiA9BHx33B8plDSaYLgBM+6yfRezIG2qdggGPHa5pWBqBpWplSKuglVUd2hyAzPLgBX/eNQ8EhTzmdNcGfZp5ZEPuBh4IgCD3EvM3eXX5gcMrDhJVqHoGcPBlmMQNQUwMHYnu5TBFd1hBd1rGzNquEMmWMmRE4sABVsH1COdbWdMc56UxOMMTry+sByE3NZVKG0TgfK1wuV6wlmCK6rCG6rGFXXTFiVLfZACmJRv+G2T0SwqzjEae+3hVrCaaILmuILuvYWZtlNE3rdwFWAG/6rP8Y+HHAPg8Ca3zW9wN5/R3r3afn9zxgf5DX12SRRRZZhusyUBsb7gVps2WRRRZZhrRYbXdDeRa2CZiplJoKnAJuAj4bsM/LwO09/rqzgIaex4BV/Rz7MvAl4Dc9P18ye3Gr86MLgiCMcqTNFgRBiCIDJtOapnUrpW4H3kQvlfSIpmm7lVK39mx/AHgNvcTSIfQyS1/u79ieU/8G+IdS6qvAceCGsP5lgiAIoxBpswVBEKKL0qI9/64gCIIgCIIgjBBCGYAYE5RSlyul9iulDvXUNLUNSimXUmqnUqpEKbU5hjoeUUpVKqV2+cRylFJvK6UO9vzsv/Bi9HT9TCl1quealSilroiBrklKqfeVUnuVUruVUt/picf0mvWjK6bXTCmVpJTaqJTa3qPr5z3xWF+vYLrscI85lVLblFL/7lmP+fsxmti13ZY2e1C67PB+kjbbmi5pswenb8jtti17ppU+ccABfCYOQB8ss6ffA6OEUsoFLNU0rTrGOs4HmtHrxc7rif0WqPWZWCFb07SoTqwQRNfPgGZN034XTS0BuvLQB1BtVUqlA1uA1cDNxPCa9aPrRmJ4zZRSCkjVNK1ZKRUPrAO+A3ya2F6vYLouJ/b32PeApUCGpmmr7PB+jBZ2brelzR6Urp8R+/eTtNnWdEmbPTh9Q2637doz3TvpgKZpnYB34gDBB03T1gK1AeFr0CdUoOfn6mhqgqC6Yo6maWVaz5TJmqY1AXvR6+rG9Jr1oyumaDrNPavxPYtG7K9XMF0xRSk1EbgSeNgnHPP3YxSRdnsApM22hrTZ1pA22zrharftmkwHm1DALmjAW0qpLUqpr8daTAB+EysAdpr76Hal1I6eR4oxfdytlCoEzgA+wUbXLEAXxPia9Tz+KgEqgbc1TbPF9QqiC2J7ve4B/hPw+MRifq2iiJ3bbWmzB4e02dZ0gbTZVnRB7O+xewhDu23XZNruEweco2naYmAlcFvPIzKhf+4HpgOLgDLg97ESopRKA54H7tA0rTFWOgIx0RXza6ZpmlvTtEXARGCZUmpetDWYEURXzK6XUmoVUKlp2pZovaYNsXO7LW22dWLe/niRNjt0pM0OnXC223ZNpk8CvtMKTgRKY6TFgKZppT0/K4EX0R9v2oWKHj+X19dVGWM9AGiaVtHzZvIADxGja9bj13oeeFLTtBd6wjG/Zma67HLNerTUA8XoHreYXy8zXTG+XucAV/d4c58BLlJKPYGNrlUUsG27LW22dezS/kibPTikzQ6JsLXbdk2meycdUEoloE8c8HKMNQGglErtGXCAUioVuBTY1f9RUcU7sQL0M7FCtPHemD1cSwyuWc8giL8BezVN+1+fTTG9ZsF0xfqaKaXGKaWyen5PBj4F7CP218tUVyyvl6ZpP9Y0baKmaYXo7dV7mqZ9Hpu+HyOELdttabMHR6zbnx4N0mZb0yVttgXC2m5rUZ7qNtQFfUKBA8Bh4Cex1uOjaxqwvWfZHUttwNPoj0a60HuFvgqMAd4FDvb8zLGJrseBncCOnhs1Lwa6zkV/7LwDKOlZroj1NetHV0yvGbAA2Nbz+ruAu3risb5ewXTF/B7r0VEE/NsO1yoGf7vt2m1pswetK+bvJ2mzLeuSNnvwGofUbtuyNJ4gCIIgCIIgDAfsavMQBEEQBEEQBNsjybQgCIIgCIIgDBJJpgVBEARBEARhkEgyLQiCIAiCIAiDRJJpQRAEQRAEQRgkkkwLgiAIgiAIwiCRZFoQBEEQBEEQBokk04IgCIIgCIIwSP5/dcgSU1eOo0IAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Import numpy, pandas, linearsolve, matplotlib.pyplot\n", "import numpy as np\n", "import pandas as pd\n", "import linearsolve as ls\n", "import matplotlib.pyplot as plt\n", "plt.style.use('classic')\n", "%matplotlib inline\n", "\n", "# Input model parameters\n", "parameters = pd.Series(dtype=float)\n", "parameters['alpha'] = .35\n", "parameters['beta'] = 0.99\n", "parameters['delta'] = 0.025\n", "parameters['rhoa'] = .9\n", "parameters['sigma'] = 1.5\n", "parameters['A'] = 1\n", "\n", "# Funtion that evaluates the equilibrium conditions\n", "def equations(variables_forward,variables_current,parameters):\n", " \n", " # Parameters \n", " p = parameters\n", " \n", " # Variables\n", " fwd = variables_forward\n", " cur = variables_current\n", " \n", " # Household Euler equation\n", " euler_eqn = p.beta*fwd.c**-p.sigma*(p.alpha*fwd.a*fwd.k**(p.alpha-1)+1-p.delta) - cur.c**-p.sigma\n", " \n", " # Goods market clearing\n", " market_clearing = cur.c + fwd.k - (1-p.delta)*cur.k - cur.a*cur.k**p.alpha\n", " \n", " # Exogenous technology\n", " technology_proc = p.rhoa*np.log(cur.a) - np.log(fwd.a)\n", " \n", " # Stack equilibrium conditions into a numpy array\n", " return np.array([\n", " euler_eqn,\n", " market_clearing,\n", " technology_proc\n", " ])\n", "\n", "# Initialize the model\n", "model = ls.model(equations = equations,\n", " n_states=2,\n", " n_exo_states = 1,\n", " var_names=['a','k','c'],\n", " parameters = parameters)\n", "\n", "# Compute the steady state numerically\n", "guess = [1,1,1]\n", "model.compute_ss(guess)\n", "\n", "# Find the log-linear approximation around the non-stochastic steady state and solve\n", "model.approximate_and_solve()\n", "\n", "# Compute impulse responses and plot\n", "model.impulse(T=41,t0=5,shocks=None)\n", "\n", "fig = plt.figure(figsize=(12,4))\n", "ax1 =fig.add_subplot(1,2,1)\n", "model.irs['e_a'][['a','k','c']].plot(lw='5',alpha=0.5,grid=True,ax=ax1).legend(loc='upper right',ncol=3)\n", "ax2 =fig.add_subplot(1,2,2)\n", "model.irs['e_a'][['e_a','a']].plot(lw='5',alpha=0.5,grid=True,ax=ax2).legend(loc='upper right',ncol=2)" ] } ], "metadata": { "anaconda-cloud": {}, "celltoolbar": "Raw Cell Format", "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.9.7" } }, "nbformat": 4, "nbformat_minor": 1 }