{ "cells": [ { "cell_type": "markdown", "id": "03d685b8", "metadata": {}, "source": [ "# Example 5, gradients\n", "\n", "For this example we will show how this code can allow gradients involving the contours of a function can be to be taken with respect to the parameters of a function." ] }, { "cell_type": "code", "execution_count": 1, "id": "28b92ba7", "metadata": {}, "outputs": [], "source": [ "%load_ext autoreload\n", "%autoreload 2" ] }, { "cell_type": "code", "execution_count": 2, "id": "1dac839a", "metadata": {}, "outputs": [], "source": [ "import jax\n", "import jax.numpy as jnp\n", "\n", "jax.config.update(\"jax_enable_x64\", True)\n", "\n", "import matplotlib.pyplot as plt\n", "import matplotlib.colors as colors\n", "\n", "from jax_zero_contour import ZeroSolver" ] }, { "cell_type": "markdown", "id": "da577e9d", "metadata": {}, "source": [ "To demonstrate how gradients work let's write a function that takes in a position along with a couple of free parameters. The shape of the zero contours will change depending on the values of these parameters." ] }, { "cell_type": "code", "execution_count": 3, "id": "760fa33d", "metadata": {}, "outputs": [], "source": [ "def f(p, A, B):\n", " return A * ((p + B)**4 - 30 * p**2 - 20 * p).sum(axis=0) + 2.0" ] }, { "cell_type": "markdown", "id": "2b8b7c69", "metadata": {}, "source": [ "Next we need to write a function that takes in these parameters and calculates some quantity of interest using the contours. For this example we will calculate the centroid of each contour. Inside this function we just need to `jax.tree_util.Partial` our input parameters into the function before proceeding as normal.\n", "\n", "```{note}\n", "For plotting reasons we will also return the calculated contours, when taking the jacobian later we will use teh `has_aux=True` to indicate the function returns this extra data.\n", "```" ] }, { "cell_type": "code", "execution_count": 4, "id": "0c99629c", "metadata": {}, "outputs": [], "source": [ "zs = ZeroSolver()" ] }, { "cell_type": "code", "execution_count": 5, "id": "866cd0f5", "metadata": {}, "outputs": [], "source": [ "def centroids(A, B):\n", " _f = jax.tree_util.Partial(f, A=A, B=B)\n", " init_guess = jnp.array([[-2.0, 0.0], [4.0, 4.0]])\n", " paths, _ = zs.zero_contour_finder(\n", " _f,\n", " init_guess\n", " )\n", " return jnp.nanmean(paths['path'], axis=1), paths" ] }, { "cell_type": "markdown", "id": "da49b7c9", "metadata": {}, "source": [ "As a test we can see what the output looks like for a given set of values for `A` and `B`." ] }, { "cell_type": "code", "execution_count": 10, "id": "c3057540", "metadata": {}, "outputs": [], "source": [ "center, paths = centroids(0.01, 0.5)" ] }, { "cell_type": "code", "execution_count": 24, "id": "25e30610", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaMAAAGdCAYAAAC/5RwpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABnsUlEQVR4nO3dd3xUVfrH8c9MJpn03klIo/fepQuoqNiwYAEVG1hWdy27ru76W2VXXXUtK5YVu9hFBBWkV+klQIAAISEhnfRkMuX+/rhJAKUkJDN3yvN+veaVkATuEyDznXPuOc/RKYqiIIQQQmhIr3UBQgghhISREEIIzUkYCSGE0JyEkRBCCM1JGAkhhNCchJEQQgjNSRgJIYTQnISREEIIzRm0LuBcbDYbeXl5BAUFodPptC5HCCFECymKQmVlJfHx8ej1Zx//OHUY5eXlkZiYqHUZQgghWiknJ4eEhISzft6pwygoKAhQv4ng4GCNqxFCCNFSFRUVJCYmNj2fn41Th1Hj1FxwcLCEkRBCuLDz3WqRBQxCCCE0J2EkhBBCcxJGQgghNCdhJIQQQnMSRkIIITQnYSSEEEJzEkZCCCE0J2EkhBBCcxJGQgghNCdhJIQQQnMSRkIIITQnYSSEEEJzTt0o1ZNZrDaOl9dRVGWiuNJEaXU9VSYLVSYLdWYbeh3odKDX6dDpdOh14O2lJ8zfhzB/b8ICfE5739tLXncIIZyXhJHGFEXhSHE16XkV7MkrZ39+JVnF1Rw7UYvFprTJNfQ6iA/1IynCn/bh/rQPDyApwp/UqADSogIlqIQQmpMw0kBOaQ3LMwrZeLiETUdKKamuP+PX+XjpiQoyEhVkJCLAh0BfAwFGA37eXigK2BQ1rGyKgk1RMJltnKgxU1ZTT2lNPSeq6ymrNWNT4NiJWo6dqGUdJb+7RqfYQLrHhdAtPpju8cF0iw/G30f+awghHEeecRwku6SGb7Yf4+c9Bew7XnHa54wGPV3jgunRLpiuccGkRgaSHOlPTJAven3rjlu32hSKq0xkl9ZwtKSG7NIaskuqOVpaQ2ZBFZUmC+m5FaTnnqzJS6+jR3wwg1LCGZisPsICfFpVhxBCnItOUZS2mQuyg4qKCkJCQigvL3fJw/XMVhuLdh1n/uZsNh4ubfq4XgcDk8MZ2SmKwSnh9EoIxcfg+KkyRVHIKa1lT145e49XsKdhqrCgwvS7r+0YHcjwDpGM6hTFkNQI/Hy8HF6vEML1NPd5XMLIDqpMFuZvyuZ/a49wvLwOUBcbjOgQyZV92jGuS7RTjzSOnahhc1Ypm46cYHNWKZmFVad93segZ3BKOKM7RzO6cxRpUYEaVSqEcHYSRhowW23M35TNK78cbLoPFBlo5JYhSVw3IIH4UD+NK7wwpdX1bDpSwqoDxazaX0heQ8A2SosK4JIecVzSM5ZuccHnPV5YCOE5JIwcbH1mMU9+l87h4moAkiP8uWdUGlP6tsPX232mtBRFIbOwipX7i1h1oIhfj5Rgtp78L9Q+3J9LesZyWc84erYLkWASnkdRQLEBOtDLSlUJIwcprzUzZ/E+5m/OASAiwIeHxnfkhkHtPWLJdEWdmeX7Cvkx/Tgr9xdhstiaPpcaFcA1/RKY0rcd7Vx0VCgENiucyFIflflQmdfwNh8qj0NNCZiqoL4KLCbglKdULyMYfMEYBP7h6iMoDkKTICwJIjpCTDfwCdDom7M/CSMH2J59gtmfbie3rBaAW4cm8aeJnQny9da4Mm1Umyys3F/E4vTjLNtXQJ35ZDANSQ3n6n4JXNozjkCjLOIUTqoiD/J2QNE+KMyAogwoPgCWuvP+1gung/BUSBgAScMheYT6azeZVZAwsiNFUZi3Los5P+7DbFVIjvDnhet6MzA5XOvSnEZlnZkf0/P5dlsuGw6f3NsU4OPFFX3aMW1we3q0C9GwQuHxLPWQvxuObYKcXyFnM1QcO/PXGnwhLAWC49WRTVCs+giOB/8I8AkEY6D6dTov0HupIypLnfowVUBNKVQXqyOrE0fVkVZRBlQV/P56YcnQZTJ0vRwSBrn0dJ+EkZ1YbQpPf5/OxxuzAbisZxz/vKanx46GmiO3rJbvtufy9dZjTffUAHonhHDT4PZc3jteNtkK+1MUKNwLh5ZD5jLI3vD7EY/OC6K7QnQ3iOqsvh/VRQ0HvZ3u/VYVQv4uyN4IWesgdwtYT9kIH5YM/W6FPtPUAHQxEkZ2UGe28sBn21mytwCdDv5yaVfuGJEiN+mbSVEUNh4u5dNN2fyUfrxp4UOQr4EbBiZy27BkEsL8Na5SuJW6CshcqobPoeXqPZ5T+YWpI4/Ehkd8P3WEoyVTFRxaBvt+gP0/Qn2l+nG9AXrfACMehog0bWtsAQmjNmayWLn7o62s3F+Ej0HPK9f34dKecZrW5MqKq0x8tfUYn23K5mhJDaB2fpjUPZbbR6TQPylM4wqFy6opVZ/E932vBtCpowyDHyQPh7RxkDZWHf0484vJ+mrY8x1sfV+dTgTQ6aHX9TDuKXWa0MlJGLUhi9XG7E+389OefHy99cybPoihaRGa1eNObDaFlQcK+d/aI6zLPHlvqU9iKHePTGVi99hWt0QSHqCuHPYugD3fwpHVYLOc/FxEB+g0CTqMg/bDwNtXuzpbI2cTrH4RDv6s/to7AEb+EYbOAoNR29rOQcKoDT353W4+3piNj5ee/00fwEUdozSrxZ3tO17Be2uPsGBHHvVWdSVeh+hA7hudxuW94z1iqbxoAZsVDq+EnZ/BvoWn3/+J7g7droBuV6r3fJx59NNSudvgp8fVRRcAsT3h6nchuou2dZ2FhFEb+XjjUZ78Lh2dDt66uT8TurveDURXU1Rp4sMNWXywPouKOvUVbkKYH/eMSuPa/glutYlYXIDSw7DtQ9j5uboyrVFkZ+h9PXS9EiI7aFefIygK7Pocfv6zus/J4AsXPwOD7nK64JUwagNbj5Zy/VsbsdgUHp3UmftGu/l/cCdTWWfmo41H+d+aI03tlWKCjcwe25HrByRq0lxWaMRmhYNLYfM7kPnLyY/7hkLPa6HPTeriAyd7Ira7ygJYcN/Jv5O7VkJ8X01L+i0Jo1aqMlm45D+rySmtZXKvOF67sa+smtNIbb2Vzzdn89bqw02NZxPC/HhwXEeu6tsOg0zfua+aUnUUtOU9KDt68uMdxqvLnTtNcur7JQ6hKLDpHXW/0ri/qh+zmJzm70XCqJUe+2oXn2/JoV2oHz8+dBHBso9IcyaLlc835/Da8kyKKtVjLlKjAnj44k5c2iNOFjq4kxNZsOEN2PYRWNQOJ/iGQt+bYeAdaocC8XuKAuv+A+lfw/QfwFf7jeUSRq2w/lAxN73zKzodzJ85hMGpsnLOmdTWW/loYxZvrjzEiRozoG6gfXJyN+mC4eryd8PaV9RVcYpV/VhsL/VeSI9rwEf2oZ1TdQn8dzBUF0HSCLj5a81XD0oYXSCz1cZlr67hQEEVtwxJ4v+m9HDIdUXLVdaZeW9tFm+vPkR1vfrENal7LI9f0oXkSPdtPOmWsn+FVf9SN3s2ShsLwx+ClJGedy+oNY7vhHmXqZtlu0yGqR/ar3tEM0gYXaB5647w94V7CQ/wYfkjowj1d95D8ISqqNLEy78cYP6mbGwKeHvpuGVIMg+O60iIv0yvOrXcbbDi2ZM34HV66H4VDH8Q4nprW5srO7IGPr5a3fA78lEY+xfNSpEwugC19VYuen45xVX1PHtVD6YNTrL7NUXb2Z9fyXOL97HqQBGgHufxxKVduaZfO1l84mzyd8OK52D/YvXXOi/oO01tdROeom1t7mLXF/DNTEAH076CjuM1KaO5z+OyDOkUH288SnFVPYnhfkwdkKh1OaKFOscG8cHtg/jg9kF0iA6kpLqeP365k6lvbSAjv0Lr8gRAWTZ8dQfMHaEGkU4PvW+E+7fAFa9JELWlXlNhwB2AooZSVZHWFZ2ThFEDk8XKW6sPA3D/mI6y29+FjeoUxeIHLuKxSV3w8/Zic9YJLnt1Lf/4YS9VJsv5/wDR9kyVsOwZeH0gpH+lfqz71XDfr3DVXFkdZy+T5kBMT6gthZ+f0Lqac5Jn3AY/7s6nuMpEbLAvV/Vrp3U5opV8DHruHZ3GL4+MYlL3WKw2hXfXHmHiy6tZc9C5XyG6FZtV3Sf0Wn9Y82+1ZU/yRXD3arhuHkR10rpC92YwwhWvqiPQ3V+qG4edlMPC6J///Cc6nY6HHnrIUZdskQ82ZAEwbbBnHBfuKdqF+jH3lv7MmzGQxHA/cstqueV/m3jim11U1Jm1Ls+95W6Dd8bA9/erGzLDU+H6T+C2hbI4wZHa9YPB96rv//wX9QWCE3LIs+7mzZt566236NWrlyMu12L78yvZnl2Gt5eOGwa117ocYQdjOkfz04MjuW2ouijls005THx5NSv3F2pcmRsyVcKPj8O749RlxsYQmPCsOiXXdbIs09bC6MfUs5uK98PO+VpXc0Z2D6OqqiqmTZvGO++8Q1iYc55Rs2BHLgCjO0cTFeQcLTRE2wswGvj7lT2Yf9cQkiL8OV5ex/R5m3nim93U1jvnq0WXk7EI3hgMv74Jig16Xgf3b4Vhs8Eg2yQ04xuirlQEWPlPsDrfrIDdw2jWrFlcdtlljB9//mWFJpOJioqK0x72pigKC3aonX+v7OP8B1WJ1huSGsGPD17EjOHJ6HTw2aZsrnxjLQcLKrUuzXVVl8Dnt8D8m6AiF0KT1N3/17wLgXLkilMYNBMCoqE8Wz140MnYNYzmz5/Ptm3bmDNnTrO+fs6cOYSEhDQ9EhPtv7x6T14FuWW1+Hl7Ma5LjN2vJ5yDv4+Bpy/vzoe3DyIy0MiBgiouf30tn2/Oxom33jmnQ8vhzWHqE5zeACP+APdtVJuZCufh7QcDblff3zhX21rOwG5hlJOTw4MPPsgnn3yCr2/zeiM98cQTlJeXNz1ycnLsVV6TxnsGIzpG4ucj5+R4mos6RvHjgxdxUcdI6sw2Hvt6Nw/M3yFLwJvDYlJviH90FVTlQ2QnuHMZjP+b9JBzVgNuB723eoT58Z1aV3Mau4XR1q1bKSwspF+/fhgMBgwGA6tWreLVV1/FYDBgtf5+jt5oNBIcHHzaw95W7FeX+Y7pHG33awnnFBVk5IMZg3hsUhe89DoW7szjqjfWkVVcrXVpzqswA94ZCxteV3894A64axXE99G0LHEeQTHQ5VL1/d1falvLb9gtjMaNG8fu3bvZsWNH02PAgAFMmzaNHTt24OWl/Sikpt7CjpwyAEZ2itS2GKEpvV7HvaPT+OLuocQEGzlYWMUVr6+VPUlnkv6NumS7IB38I+DG+TD5JRkNuYqe16lv078Bm03bWk5htzAKCgqiR48epz0CAgKIiIigRw/n6IS9I7sMq00hPsSXhDD5QRLQPymMhbNH0Ld9KBV1FmbM28yXW+w/XewSbFZY+hR8NQPMNZAyCu7dAJ0v0boy0RIdLgZjsLrQJHeL1tU08ejdnVuOngCgv5yBI04RHezL/LuGMKVPPBabwp++2sVryw569sIGcy18cat6cBvAsAfg5m/UaR/hWrx9IW2M+n7msnN/rQMZHHmxlStXOvJy57U7txyAPomh2hYinI7R4MVLU/sQF+rHmysP8e+lB6gxW3l0YmfP6wBeUwqfToVjm8HLB6a8CT2v1boq0Rpp42DvAnUl5Bjn6Fnn0SOjxk7O3eIce6S5cA16vY7HJnXhqcndAHhz5SGeXbTPs0ZI1cXwweVqEPmGwq0LJIjcQdpY9W3uFqh3joU6HhtGlXVmckprAegaF6RxNcKZ3T4ipenE33fXHuH15ZkaV+QgVUXw/mR1oUJgDNz+MyQN07oq0RZCEyEoXu2S4SRLvD02jLKKawCIDDTKaa7ivG4ZksQzV3YH4N9LD/DFZjdf1FBfA59eB0X7IDAWpi+C6C5aVyXaUrt+6ttjzrGIwWPD6GipOjRNipBVdKJ5bh2azH2j0wB44tvdbD1aqnFFdmKzqoex5W0Hv3A1iCI7al2VaGvxfdW3+bu1raOBx4ZRdqk6MkoKlzASzfeniZ2Z3CsOq03hgc92UFZTr3VJbW/1C5DxA3gZ4cbPILKD1hUJe4jqrL4tOahtHQ08NowKyusAiAttXqsiIQB0Oh1zru5JUoQ/uWW1/Plb53hV2WZyt8Kq59X3r3gV2g/Rth5hPxENo93iTHCCRTkeG0ZFVSYAooMkjETLBPl68/qN/fDS61i8O58Nh0q0LqltmOvg23tAsUKPa6D3DVpXJOwpPEV9W1+pLt/XmOeGUaUaRpGBcn6RaLmeCSHcOEjtKv/c4n3YbNq/smy1HR9D8QF15dylL2pdjbA3g1Ft5wRqo1uNeWwYVdapXZlD/Lw1rkS4qofGdyLQaGB3bjnLM1z8xFirGdY2dFe46BHwl64kHiEwVn1beVzbOvDgMKozq13D/Xw89q9AtFJkoJGpA9TR0cJdeRpX00p7F6iHrgVEQb9bta5GOEpAQ4Poau2nmj32mbi2IYx8vbXvHi5c1+TecQD8sreg6QWOSzq4VH3b92b1EDbhGYwNG/7rtT/l2HPDqF7CSLRe38RQ4kN8qa63srWh8a5LOrpOfZsyUts6hGMZG1qhmaq0rQMPDqM6s3qOh5+EkWgFnU5HWnQgAMcbtgu4nIo8KM8BnRckDta6GuFIjaNgc622deChYWS1KdRb1TCSkZForZhgdXtAQYWLhlHjq2LfYPAJ0LYW4Vi6xgjQfjWoR4aRhx0AIOws2FddkVlea9a4kgul/ROR0EhjGCnan/jqkWGk1+sw6NVIqrdo/48gXFt+hTrFERvsqhuoG16eWc1OsRNfOJAThFAjjwwjAB+D+q1LGInWOnZCDaOEMBddhRaWpB6aV18FpYe1rkY4kqVhatmg/QspCSOrCy/HFZqrqbewP19dFpsa5aL3WwxGiOujvn9ss6alCAdrDCMnWM7vsWHk7aV+6yYZGYlWWLm/CJPFRvtwf9KiArUu58IlDlLfNu43Ep7B1LC/yAkWrnhsGPl4yTSdaL1Fu9U2Kpf0iEWnc+GlMY1Hie9doJ7wKjxDdcO/dUCUtnXgwWEU7OfqK6CE1nJKa/g5XW0weXnveI2raaX4vtBuANjMsO19rasRjtL4wiMgWts68OAwighQjxovrXbDw9GEQ7y1+hAWm8JFHSPp0S5E63Jab9BM9e36153iSAFhZ4pyysgoUtta8OAwCpcwEq2QXVLDF5uPATBrjJuchNrjWojuBnVl6mmvwr3VV4OlofOCTNNppzGMSiSMRAspisLfFu6h3mpjRIdIBqe4yXELXgaY8A/1/U1vQ8FebesR9lXdcOyJtz8YtV984/FhVFolYSRaZsneApZnFOLtpePvV3Z37YULv9VhHHS+DGwW+PZudSOsJ1kx5+Sx67+16nn18+6iLEd9GxSnbR0NPD6MSqpNGlciXElpdT1/+TYdgJkXpbr2cu6zmfwS+IVB/i5Y7WEnvuq9YMWzvw+kVc+rH9e7US/L4gPq28iO2tbRwGPDKD5U3XGcW+aizS2FwymKwpPf7aa4ykSH6EAeGOccP8RtLij25LHjq5+HI6u1rceRRj0KY/5yeiA1BtGYv6ifdxclmepbJwkjg9YFaCUxzB+AY6U1KIriXlMtwi4+35zD4t35GPQ6Xrm+j3t3fO95LRxaDjs+ga/ugHvWqCHlCRoDZ8Wz6kIOa737BRFA8UH1bYRzhJHHjowSGsKo0mSRvUbivHbmlPHUgj0APDKhs3ss5T6fS19UV9dVF8IXt4HFg6a0Rz2q9uuz1qtv3S2I4GQYOcnIyGPDyM/Hi8hAIwA5pdofLCWcV2l1Pfd9so16q40J3WK4Z1Sq1iU5ho8/TP0QjCGQsxEWPeI5Xb1XPX8yiKz1Z1/U4KrMteqBigCRnbStpYHHhhFA+3C1OWDOiRqNKxHOympTeHD+dnLLakmJDODFqb09a0o3siNc+z/13JvtH8GGN7SuyP5OvUf016Lf30NyB/m7AUXdX+QfoXU1gIeHUWK4OlWXXSphJM7spaX7WXOwGD9vL+be3L/pID2P0vFiuPj/1PeX/AXSv9G2Hns602KFMy1qcHW5W9W37fqDk7y48tgFDEDTstyDBVUaVyKc0bfbj/HGikMA/POannSODdK4Ig0NnQVlR9XNsN/eDYHRkDxC66rans165sUKjb+2ucmRM01hNEDbOk7h0WHUKUZ9cjlQUKlxJcLZ/Hq4hMe+2g3A3aNSubJPO40r0phOB5P+CRV5kPEDzL8Jbv8ZortqXVnbGvPE2T/nTosYmsKon7Z1nMKjp+k6xTSMjAorsdk85MasOK9DRVXc9dFW6q02Lu0Zy2MTu2hdknPQe8E170LiEKgrh4+vgfJcrasSLVVTevJE3/i+2tZyCo8Oo6SIAHwMeurMNlnEIAB15dzt72+mvNZMn8RQXpraB73eOebUnYK3H9z4mboCqyIXPpoi5x+5msZRUXga+DtPX0WPDiMvvY4ODfeNDsh9I49XZ7Zy14dbOFpSQ0KYH+/eNsC9N7ZeKP9wuPkbCE5QW8p8fBXUlmldlWiuI6vUt+2HalvHb3h0GAFNN6X3Ha/QuBKhJZtN4U9f7WLL0RME+Rp4f8bApn1o4gxCE+HWBeqhbPm74ZNrwSQv6FzC4ZXq27QxmpbxWx4fRo076XfnlmtcidDSS0sPsHBnHga9jrdu7k+HaA9eOddckR3g1u/ANxSObYbPblA3UwrnVV3csMcISBmpbS2/4fFh1CtBDaNdx8q0LURo5tNfs3l9hdo08rmrezKsg/anXrqMmO7qlJ1PIGStaWgbJMeyOK3GKbro7uryfCfi8WHUPT4YvQ4KKkwUVEgHb0+zbF8BT36nvlJ8YFxHpg5I1LgiF5TQH276HAy+cPBn+PYu99mP424ap+hSR2tZxRl5fBj5+xjo2DAlszOnTNtihEPtzClj9qfbsSlwbf8E/jDeORpGuqTkEXD9J6D3hj3fwoJZEkjORlEgc7n6fuoobWs5A48PIzg5VSf3jTzH0ZJqbn9/M7VmKxd1jGTO1T09q+ecPXQcD9e+Bzov2PkZfH+/BJIzOb4TKo6px4w72f0ikDACoHdiKADbs8s0rUM4Rml1PdPnbaakup5uccG8eXN/vL3kR6FNdLtC3Rir81LPQvr+AbDZtK5KgNo5A9Sj5b39tK3lDOQnEBiQHAbAtuwTWKzyg+POauut3PHBZo4UV9Mu1I/3Zwwk0OjRXbHaXo+r4Zp3GgLp44YRkvxcaS5jkfq2y2Rt6zgLCSOgU3QQQb4Gauqt7JX9Rm6r8TiI7dllBPsa+OD2gUQH+2pdlnvqcQ1c/bZ69MSOj2GhjJA0VXIICveC3gCdJmpdzRlJGAF6vY4BSeroaHPWCY2rEfagKArPLNzDkr0F+Hjpefe2gbKXyN56XgtXv3PyLKQfHpRA0krjFF3yCPAL07aWs5AwajAwRe3RtPlIqcaVCHt4Z81hPthwFICXr+/DoBTn6cnl1npeC1c1jJC2fQg/PCSBpIXGM6icdIoO7BxGc+bMYeDAgQQFBREdHc2UKVPYv3+/PS95wQYmN4RRVimKpxyt7CG+35nHc4szAHjysq5c1itO44o8TK/r4Kq3GgLpA1j0BwkkRyraD8d3qFN03a/WupqzsmsYrVq1ilmzZrFx40aWLl2K2WxmwoQJVFdX2/OyF6RXQghGg56S6noOFkqPLXex4VAJf/xiJwC3D0/hzotSNa7IQ/WaClPmAjrY+j4sfkQCyVF2zlffdrgYApzjiPEzsesyop9++um0X7///vtER0ezdetWRo50rnXuRoMXg1LCWXOwmLUHi5sO3hOu60BBJXd9tKXpXKInL3Ozg+BcTe/rAQW+vQe2vAfo4LJ/O82x127JZoPdX6rv975e21rOw6H3jMrL1U2l4eHOOV8/oqEn2drMYo0rEa2VX17H9Pc2UVlnYWBymJxL5Cx63wBT3gR0sOV/sOgRtTOAsI+j66A8B4zB0OkSras5J4dtsLDZbDz00EMMHz6cHj16nPFrTCYTJpOp6dcVFY5dZj28IYw2Hi7BbLXJRkgXVVlnZvq8TeSV15EWFcA7t8q5RE6lz42AAt/dpwaSTg+XviAjJHtonKLrdiV4O/c2Boc9286aNYv09HTmz59/1q+ZM2cOISEhTY/ERMc2rewWF0xEgA819VbpxuCi6i027v14Gxn5lUQFGXl/xiBC/X20Lkv8Vp+b4Mo3AB1sfgd+fFRGSG2ttgzSv1bf7zNN01KawyFhNHv2bH744QdWrFhBQkLCWb/uiSeeoLy8vOmRk5PjiPKa6PW6puMDZKrO9SiKwuNf72JtZjH+Pl7Mmz6QxHB/rcsSZ9N3GlzxGqCDTW/Dj49JILWlnfPBUgvR3aD9EK2rOS+7hpGiKMyePZtvv/2W5cuXk5KScs6vNxqNBAcHn/ZwtBEd1NUmaw8WOfzaonX+veQA32zPxUuv47/T+jUdnCicWL9bGgIJ2PQW/PSEBFJbUJSGRSLAgNtdYgrUrveMZs2axaeffsqCBQsICgoiPz8fgJCQEPz8nK9RH8CIjlEA7DxWTkWdmWBfb40rEs3xya9Hmw7Im3NVT0Z3dq6Dw8Q59LsFUNQedr++qT5xTnzOJZ5AndbR9VC8H7wDoJdzr6JrZNeR0Ztvvkl5eTmjR48mLi6u6fH555/b87Kt0i7Uj9TIAKw2hfWZJVqXI5ph2b4C/vpdOgAPje/I1IFyQJ7L6XcrXP4f9f2N/4Wf/yIjpNbY8j/1bc9rwdfxM0wXwq4jI1ftZDCqcxSHi6tZtq+AST1itS5HnMOuYycPyJs6IIEHx8kBeS6r/3Q1gH54CDa+oY6MJvxDRkgtVZ4Lexeo7w+4XdtaWkDWLp/Bxd1iAFiWUYjV5pqB6gnyy+u484Mt1JqtjOoUxbNXyQF5Lm/ADJj8svr+htdh6V9lhNRSv84FmwWSRkB8H62raTYJozMYlBxOiJ83pdX1bD0qXbydUW29lZkfbqGw0kTnmCDemNZP9oW5iwG3w2Uvqe+vfw2WPiWB1Fx1FWq7JYBh92taSkvJT+8ZGLz0jO2i3gBfujdf42rEbymKwh+/2snu3HLCA3x497YBckCeuxl4h9oqCGD9q/DL0xJIzbHtQzBVQGQn6DhB62paRMLoLBqn6pbuLXDZe1/u6rXlmSzadRxvLx1vTusne4nc1cA74dIX1ffX/Qd++ZsE0rlYzbDxTfX9obNB71pP765VrQON7BSFj5eerJIaMqWLt9P4cfdxXlp6AIB/TOnB4FTn7UIs2sCgmacE0iuw7BkJpLNJ/xoqjkFAlMss5z6VhNFZBBoNDGvYALtkb4HG1QiA9Nxy/vDFDgDuGJHC9QPba1uQcIxBM+GSF9T3174kgXQmVgusbvg7GnKv0/ehOxMJo3M4dapOaKuo0sTMD7dQZ7YxqlMUT1zSReuShCMNvgsm/Ut9f+1LsOzvEkinSv8KSjLBLxwG3aV1NRdEwugcxndVw2hHThmFFXUaV+O5LFYbsz/dxvGGLtyv3dQXg6yc8zxD7jklkF6We0iNrBZY9bz6/rD7weiaZ7HJT/Q5xAT70jsxFIBf9hVqW4wH+9dPGfx6pJRAo4G3bhkgLZo82ZB7Tk7ZrXtFln2DOioqPeTSoyKQMDqvCU1TdbLEWwuLdh3nnTVHAHjxul50iA7UuCKhucF3nVzUsP5VWPKk5waS1XxyVDT8ATC67s+HhNF5NIbRuswSqkwWjavxLJmFlfzpq50A3DMqjUk94jSuSDiNQTNP7kPa8Lrn9rLb9oE6KvKPhIEzta6mVSSMzqNDdCDJEf7UW22sPiDHSjhKZZ2Zuz7aSk29lWFpEfxxQietSxLOZuCdJzs1bHwDfv6zZwVSXQWsmKO+P/pxlx4VgYTReel0OllV52CKovDoV7s4XFRNXIgvr94oCxbEWQy8Aya/or6/8b/w0+OeE0jrXoGaYojooDaZdXHyE94MF3dTO3cvzyjEbLVpXI37e2fNYX5Mz8fbSz0kLzLQqHVJwpkNmHHy+Ilf53rGibHlubDhDfX9i58BL9df1CNh1Az9k8IID/ChvNbM5qxSrctxa9uyT/D8T/sBeOry7vRtH6ZxRcIl9J9+yhHmb8HiP7l3IC3/P7DUQfth0PlSratpExJGzeCl153SOFWm6uylvNbMA59tx2JTuLx3PDcPlg4LogX63XoykDa/A4seAZsbzmRkb4Sdn6nvu9F5TxJGzSSNU+1LURT+/M1ujp2opX24P89e1UPOJhIt1+8WuPINQKeedrrYzQLJaoYfHlbf73szJPTXtp42JGHUTBd1jMRo0HPsRC0Z+ZVal+N2Pt2UzaLdxzHodbx2Y1/Z2CouXN9pMOW/qIH0Hiz6g/sE0q9zoXAP+IXB+Ge0rqZNSRg1k7+PgYs6RgIyVdfWDhRU8szCvQA8NqlLU9cLIS5Yn5vgqrmATj1s7vv7wWbVuqrWKT92cin3xf8HAe7VsV7CqAUap+qWSDeGNlNvsfGHz3dgsqgNUO8YkaJ1ScJd9L4Brn4bdHrY8TF8M1Od5nJVPz0O5mpIHAJ9pmldTZuTMGqBsV1i0OkgPbeCvLJarctxC68tP8ievApC/b154dpe6PVyn0i0oV5T4dp5oDeo5/18cRtYTFpX1XK7v4J9C0HnBZNfcrmD85rD/b4jO4oKMtKvYanxL/tkqq61tmWf4I0VmQA8O6Un0cGudwaLcAHdp8ANn4KXEfYvgs9uhPoaratqvoo8dWUgwMg/QUx3beuxEwmjFpog3RjaRG29lUe+2IlNgSl94rmsl/SdE3bUaSJM+wK8/eHQMvjkWjC5wEIkRYEFs6GuDOL7wsg/al2R3UgYtVDjfaONh0uoqHPh+WeN/eunDI4UVxMb7Mvfr+ihdTnCE6SOhpu/AWMwHF0HH06B2hNaV3Vum99Vw9PgC1e97RadFs5GwqiFUqMCSYsKwGxVWLlfGqdeiG3ZJ/hgQxYAz1/bixB/9/0BE04maSjcukBdGp27Bd6/HCqddEFScSYs+av6/vi/Q5R7NwuWMLoAjb3qZKqu5eotNp74ejeKAtf0S2BkpyitSxKepl0/mL4IAqKhYDf8bwKUHNK6qtPV18AXt4KlFlJGufShec0lYXQBGqfqVmYUUm9xk810DvL26kPsL6gkIsCHJy/rqnU5wlPFdIc7foawFCg7qgZS7jatq1IpirpgoXCPGphXv+2Wq+d+y/2/QzvomxhKZKCRSpNFGqe2wKGiKl5dpq6ee+ryboQF+GhckfBo4alwxxKI660exfD+ZDi0XOuqYNPbsPNTdX/Ute9BUKzWFTmEhNEF0Ot1jO6sTi8tzyjUuBrXoCgKf/l2N/VWGyM7RXFF73itSxICAqPVKbuUUeqG0k+mwtYPtKvn4C/q5lZQ7xOlXKRdLQ4mYXSBGrt4r5AwapYf0/PZeLgUo0HPs1OkCapwIsYgmPYl9LgWbGZY+IB6JpLV4tg6ju+Cr2aAYoM+N8Ow+x17fY1JGF2gizpGYtDrOFxczZHiaq3LcWp1ZivPLtoHwD2j0kgM99e4IiF+w2CEa96FMX9Rf/3rXHUvkqOWfhdnwkdXgakCkoarXRY87AWbhNEFCvL1ZlBKOCBTdefz7prD5JbVEhfiyz2j0rQuR4gz0+lg1KMw9SN1c+zhFfDOWCjab9/rlh6Gj6ao961ie8GNn6nh6GEkjFpBpurOL7+8jjdWqMtmH7+kC34+XhpXJMR5dLtCXdgQkqgGxdtjYPsn9jk5Nj8d3psE5TkQ0VHdlOsb0vbXcQESRq0wpiGMfj1SQpXJwfPLLuL5nzOoNVvpnxQmixaE64jtCTNXQPJF6sKGBffB13dCTRuunj20HN6/FKoKIKaHupAi0HP33UkYtUJqZABJEf6YrQprDxZrXY7TySys5NvtuQA8NbmbLFoQriUwSu3WMPavarfs9K/g9YGw68vWjZJsNlj5L/joaqgrh8TBahAFxbRd7S5IwqgVdDodYzqro6OV+2Wq7rde+eUgigITu8fIgXnCNem91OakdyyBqK7qfZ1v7oT3L4OstS3/83K3wbxJsPI5QIF+t8Gt34NfaFtX7nIkjFppVMN+o3WHZGR0qgMFlSzafRyAh8a7d08t4QESBsDdq2Hsk+pRFEfXqYH0/mT1rCFT1dl/r80GRzfAN3fDO2Mg51d1gcSUN+GKV8Fbjk4BMGhdgKsblByOQa8jp7SWnNIaWbbc4D8No6JLesTSNS5Y63KEaD2Dj3qeUO8bYc1LsO1DyFqjPgy+6sbZiA4QlqwGTMVxqDgGh1aoCxQa9boBxj8NwXIP9VQSRq0UYDTQJzGULUdPsOFQiYQRsO94BYt2H0enk1GRcEMhCeo+oBEPqd0a9nyjrro7+LP6OBOfIOh6OQy8ExL6O7RcVyFh1AaGpUWw5egJ1h0qZurARK3L0dzcVepS7kt7xtE5NkjjaoSwk9D2MO6v6tRd/i7I2QQnstSHpU4d+QTFq01ZO14M3n5aV+zUJIzawNC0SF5dnsn6QyUoiuLRq8byympZtEu9V3SvbHAVnkCnU5utxvXWuhKXJgsY2kDf9qEYDXqKKk0c9vDWQB9syMJiUxiSGk6Pdp65eU8I0XISRm3A19uLng1PvNuzy7QtRkM2m8LihhV0d45I1bgaIYQrkTBqI33bhwKwI8dBjRWdkF6v46cHR/L8Nb2aWiUJIURzSBi1kT6JYQDsyCnTthANzV11iF3HyrmmfwJ6vefeNxNCtJyEURvp0zAyyjheSZ3Zqm0xGsgtq+WfP2Zw07sbKa4yaV2OEMLFSBi1kfgQX6KCjFhsCum55VqX43A/NtwrGpgcTkyw7CgXQrSMhFEb0el09Gnov7bbA8Poh4bl3JN7xWlciRDCFUkYtaFOMYEAZBaeo0+VGyqoqGNHThk6HUzqEat1OUIIF+SQMHrjjTdITk7G19eXwYMHs2nTJkdc1uHSotQwOlTkWWG06kARAL0SQokOkik6IUTL2T2MPv/8cx5++GGefvpptm3bRu/evZk4cSKFhe535EKH6MaRkWdtfF21Xw2j0Z0892AwIUTr2D2MXnrpJWbOnMmMGTPo1q0bc+fOxd/fn/fee8/el3a41IaRUXGVifIas8bVOIaiKKxvOD5jZKdIjasRQrgqu4ZRfX09W7duZfz48ScvqNczfvx4NmzY8LuvN5lMVFRUnPZwJYFGA3Eh6jRVpodM1R0pruZEjRkfg56e7UK1LkcI4aLsGkbFxcVYrVZiYk4/TjcmJob8/Pzfff2cOXMICQlpeiQmul4H7JTIAACySz1jqm5bQ/ujXu1C8DHIehghxIVxqmePJ554gvLy8qZHTk7O+X+Tk4kKMgJQXFmvcSWO0binSo4VF0K0hl2PkIiMjMTLy4uCgoLTPl5QUEBs7O+XABuNRoxGoz1LsrvIwIYw8pAuBI3L2DvHyLlFQogLZ9eRkY+PD/3792fZsmVNH7PZbCxbtoyhQ4fa89KaiQj0AaDIQ8LoYGElAB0a9lgJIcSFsPvheg8//DC33XYbAwYMYNCgQbzyyitUV1czY8YMe19aEydHRu4/TVdntlJQoYZuSkSAxtUIIVyZ3cPo+uuvp6ioiKeeeor8/Hz69OnDTz/99LtFDe4iqjGMKt1/ZFTU8D36GPSE+ntrXI0QwpU55Njx2bNnM3v2bEdcSnONT8rlte6/z6iwsg6A6CCjRx+1LoRoPadaTecOfL29ADBZbBpXYn9lDRt7wwN8NK5ECOHqJIza2Mkwcv8zjRoD19fgpXElQghXJ2HUxgwNJ5xarIrGldhfY+AaveW/kRCideRZpI0pDRnkCadu29x/JlII4SASRm3M1pBGeg+4od84IvKE+2NCCPuSMGpjZqv6xOzl5f5h1HivyGR2//tjQgj7kjBqY9X16hNzgI9DVs1rKshX/R7LPGAZuxDCviSM2li1yQJAgNH9V5jFBKvHZRRWmFAU91+wIYSwHwmjNnaiRm0DFOLn/h0JooPVbhO1ZisVdRaNqxFCuDIJozaWX652JWgcNbgzfx8DEQ0bXrNLajSuRgjhyiSM2lhBhRpGsR4QRgAdotVu3QcKKjWuRAjhyiSM2tjxhpFRbIhnhFGnhnOMDhRKGAkhLpyEURtrGhl5SBh1jw8GYGdOmbaFCCFcmoRRG8tvCKM4Dwmj/klhAOzIKWvaYyWEEC0lYdSG6sxW8srUMEoI89e4GsdIiwokxM+bOrON9NxyrcsRQrgoCaM2lFlYhdWmEOrvTXSQUetyHEKv1zEkNRyAVQeKNK5GCOGqJIza0P589SZ+55ggjzpsbmyXaABWZBRqXIkQwlVJGLWh/Q3Lm7vEBmlciWON6ayG0c5j5RQ23DMTQoiWkDBqQ/uOVwDQJS5Y40ocKzrYl77tQwH4YddxbYsRQrgkCaM21DRN52EjI4Are8cD8P3OPI0rEUK4IgmjNpJbVkthpQmDXudx03QAl/WKR69Tl3hnFlZpXY4QwsVIGLWRLVmlAHRvF4K/Bxwf8VtRQcamhQyfbcrWuBohhKuRMGojmxvCaGDDJlBPdNPg9gB8tfUYdXLgnhCiBSSM2sjmIycAGJAcrnEl2hnVKZp2oX6U15r5bnuu1uUIIVyIhFEbKK8xNy3rHpDsuSMjL72OGcOTAXh7zWFsNjlwTwjRPBJGbaBxii41MoDIQM/ovHA2NwxqT5CvgcNF1fyyr0DrcoTQhNWmyIuxFvK8O+12sOag2gZnaFqExpVoL9Bo4OYhSby58hAvLT3A+K4x6PWe041CeBarTWHXsTLWZRazNrOYw0XVVJks1NRb8fXW0ykmiC6xQQxNi+CSHnH4entpXbLTkjBqA4092UZ1itK4Eudw98hUPt5wlIz8ShbuyuPKPu20LkmINlVZZ+aLLcd4f/0Rckprz/g1dWYbu46Vs+tYOV9sOcbfF+7l2n4J3HlRqsccMdMSEkatdLSkmqySGgx6HcM6RGpdjlMI9ffhrpGp/HvpAV5eeoBLe8bh7SUzwsL11ZmtvLvmMHNXHabKZAEgyGhgeIdIhneMpE9CKCF+3gQYvaios5BxvILdueUs2JFHblkt7649wmebsvnDxZ2YPiwZg/xcNJEwaqXVDaOiAclhBBrlr7PR7SNSeH99FlklNXy55VjTsm8hXJGiKCzenc9zi/eRW6aOhDpEB3L78BSu6tsOP5/fT79FBBpJiQzgkp5xPDKhMyv3F/LGiky2ZZfxj0X7+GZbLm9M60dKZICjvx2nJLHcSo1TdCNliu40AUYDs8Z0AOA/yw7IviPhsooqTdzz8VZmfbqN3LJa4kN8+c8NfVjy0EhuGtz+jEH0W156HeO6xvDVPcP459U9CfX3Zu/xCq54bS2/7JWFPiBh1Cp1ZivrMksAuV90JjcNbk+7UD8KKkz8d0Wm1uUI0WI/pedz8cur+HlPAQa9jgfGdWTZI6O5sk+7C1qYo9fruGFQe5Y8NJIBSWFUmizc+eEW3l59yA7VuxYJo1ZYe7CYWrOV+BBfunlYp+7m8PX24snLugIwd9VhDhdJzzrhGurMVp5akM49H2+lrMZMt7hgvp89gocv7tSskdD5RAf78unMIdw2NAmA5xZn8MLPGSiK5y4HlzBqhSV78wGY0D3Wow7Ta4lJPWIZ1SmKequNpxbs8egfNuEacstquXbuej7ccBRQV4d+N2s43eLb9gWnj0HP36/swaOTOgPwxopD/N8P+zz2Z0TC6AJZbQq/7FNPNp3QLUbjapyXTqfj71d0x8egZ21mMYt2y3lHwnn9eriEK15bS3puBeEBPsybMZAnLu2Kj8F+T5X3je7AP6b0AOC9dUd4eekBu13LmUkYXaAtWaWUVtcT4ufNwBTP7UfXHMmRAdw3Og2AZxbupbLOrHFFQpxOURQ+2niUae/+Skl1fcO03PCmU4zt7eYhSTxzZXcAXl2eybtrDjvkus5EwugCLWlYATOua7TsoWmGe0alkRzhT2GliWcX7dO6HCGa1Fts/PnbdP76XToWm8LkXnF8fe8wEsL8HVrHrUOT+dNEdcru2cX7+Cnds2YR5Fn0AthsCosbppsmdIvVuBrX4OvtxT+v6YVOB/M357A8Q5azCu2V15q57b1NfLYpG50OHpvUhddu7NsmixQuxH2j07h1aBKKAg99voMdOWWa1KEFCaMLsOXoCY6X1xFkNDC6syzpbq4hqRHcMTwFgMe+3s2J6nqNKxKeLK+slqlzN7DhcAkBPl7877YB3Ds6TdPFSDqdjqcmd2N05yjqzDbu/XgrJVUmzepxJAmjC7BwZx6grqKTxoct88eJnekQHUhRpYknF6RrXY7wUBn5FVz93/XsL6gkKsjI53cPZWwX51iIZPDS89qNfUmNCuB4eR0Pfb4Dqwd0AJcwaiGL1dY0RXd57ziNq3E9vt5evDy1Dwa9jkW7jrNghzaH8L289ACvLjt4xs+9uuygx65o8gTrM4u57s0N5FfU0SE6kG/vG0aPdiFal3WaIF9v3pzWH19vPWsOFp/1/6o7kTBqofWHSiipric8wIfh0hj1gvRMCGH2WLVV0F++TedIcbXDa/DS63jpDIH06rKDvLT0AF5y7IVbWrAjl9vmbaLSZGFQcjhf3TPU4QsVmqtzbBDPXdUTgFeXH2xqPeauJIxa6PuGKbpLe8bKKrpWmD2mA4NSwqkyWZj1yTaH9657YFxHHr6402mB1BhED1/ciQfGdXRoPcK+FEXhzZWHeHD+DsxWhUt7xvLhHYMI9ffRurRzurpfAjcNbq8uaJi/nePlZz6uwh3Is2kL1Jmt/Jyudl24vFe8xtW4tsZ58YgAH/Yer+CZH/Y6vIZTA6nTX36UIHJTVpvC09/v4V8/ZQBwx4gUXr+xn8vc731qcjd6tgvhRI2ZR7/a5bYdGiSMWmDVgSIqTRbiQnwZmCwbXVsrJtiXV27og04Hn/6arcn9owfGdcTHS0+91YaPl16CyM3Uma3c+/FWPtxwFJ0OnrysK3+d3M2lTh/29fbi5ev7YDSo948+3nhU65LsQsKoBRqn6Cb3inOp/8zO7KKOUdzfcNTEn7/ZzcGCSode/9VlB5uCqN5q84gbxZ6itLqem97ZyJK9BfgY9Lx+Yz/uvChV67IuSIfoQB6/pAugNlXN0uA+q71JGDVTtcnCsn3qRs3Le8sUXVt6cHwnhqVFUF1vZeaHWyirccz+o1PvER149pLf3UMSriu7pIZr3lzPtuwygn0NfHzHYC7r5dqrX28bmszQ1AhqzVYe/sL9lntLGDXTL/sKqDPbSI7wp6eTLQN1dV56Ha/d2JeEMD+ySmqY/el2LFabXa95psUKZ1rUIFzPrmNlXP3mOo4UV9Mu1I+v7x3GIDfoH6nX63hxam+CjAa2ZZfx9mr36l8nYdRM3+9Qp+gu7x0vx0XYQUSgkXduHYC/jxdrM4t5drF9+9dZbcoZFys0BpK7ver0FCsyCrn+rY0UV6nNTr+5bxgdY4K0LqvNtAv146+XdwPglV8OkF1So3FFbcduYZSVlcUdd9xBSkoKfn5+pKWl8fTTT1Nf73otYMpq6ll9UF3jf4VM0dlN17hgXpraG4B567L4YnOO3a71h3OsmntgXEf+cHEnu11b2Mf8Tdnc+eEWas1WLuoYyed3DyEm2Ffrstrcdf0TGJYWgcli48kF6W6zus5uYZSRkYHNZuOtt95iz549vPzyy8ydO5c///nP9rqk3fyUno/ZqtAlNsitXmU5o0k94nhovBoSf/luN+szizWuSDg7RVF4aekBHv9mN1abwjX9Enhv+kCCfL21Ls0udDod/5jSAx+DntUHili4yz26e9stjCZNmsS8efOYMGECqampXHHFFfzxj3/km2++sdcl7WbhrpNTdML+Hhjbkct6xWG2Ktz90VYy8iu0Lkk4qXqLjT9+uavpHt8DYzvw4nW93H5DempUILMbVqE+s3Av5TWuf0aYQ//FysvLCQ8/+41Ek8lERUXFaQ+tFVbWseFQCSBTdI6i1+v493W9GZQcTqXJwox5m91657m4MBV1Zma8v4mvtx3DS6/juat68vCEzh5zT/fuUamkRQVQXGXiXz9naF1OqzksjDIzM3nttde4++67z/o1c+bMISQkpOmRmJjoqPLOatGu49gU6JMYSmK4c/awcke+3l68fWt/0ho6F8+Yt5kKOSFWNMgrq+W6NzewLrMEfx8v3r1tADcNbq91WQ5lNHjxbEPvus82ZbMnr1zjilqnxWH0+OOPo9PpzvnIyDg9pXNzc5k0aRLXXXcdM2fOPOuf/cQTT1BeXt70yMmx3w3s5mo8LkJGRY4X6u/D+zMGERVkJCO/kns+2urwHnbC+ezJK+eq/65rOv7hi7uHOux4cGczJDWCyb3iUBR1us6VFzPolBZWX1RURElJyTm/JjU1FR8ftQFhXl4eo0ePZsiQIbz//vvo9c3Pv4qKCkJCQigvLyc4OLglZbaJnNIaLnp+BTod/PrEOKLdcGWOK0jPLef6tzZQXW/l4m4x/HdaP7e/JyDObNWBIu77eCvV9VY6Rgcyb8ZAp+267Si5ZbWMfXElJouN/07rx6U9nWtzb3Ofxw0t/YOjoqKIimre6aa5ubmMGTOG/v37M2/evBYFkTP4oWGVypCUCAkiDfVoF8I7tw1g+rzNLN1bwJ++3MlLU/tISyYP8/nmbP78bTpWm8LQ1Ajm3tKfED/3XDHXEu1C/bh7ZCqvLs/kucX7GNsl2mWawJ7KbumQm5vL6NGjad++PS+++CJFRUXk5+eTn59vr0u2ucZedFf0kSk6rQ1Li+TNaf0w6HV8tyOPv7rR/gpxboqi8O8l+3nsa3Xp9lV92/HB7YMkiE5xz+g0YoN9OXailv+tPaJ1ORfEbmG0dOlSMjMzWbZsGQkJCcTFxTU9XEFmYSX7jldg0Ou4pEes1uUIYFzXGF66Xu3y/cmv2fzzxwwJJDdnslh55IudvLY8E4D7x3bgpam98TG41iyLvfn7GHjsks4AvLnyECVVJo0rajm7/YtOnz4dRVHO+HAF3+9Up+hGdopy+gO4PMkVveObTr98a/Vh/vXTfpf5PyVaprjKxLR3fuWb7bl46XX865qePOJBS7db6sre7egeH0yVycJ/Vx7SupwWk5cXZ6AoCj/sOnlchHAuNw5qz98a+nPNXXVIRkhuKCO/gitfX8eWoycI8jUwb/pArh/oWUu3W0qv1/HoJPWYiY82HCW3zLX25kkYnUFGfiWHi6rxMei5uFuM1uWIM5g+PIVnruwOqCOk5xbvk0ByE7/sLeCa/64nt6yW5Ah/vr1vOCM7NW/RlKcb2TGSoakR1FttvLz0gNbltIiE0Rk0jopGd4py2/5W7uDWocn835QeALyz5gj/94MEkitTFIW3Vh1i5kdbqK63MjQ1gu9mDadDdKDWpbkMnU7Ho5PUe0ffbDvGAQcfVtkaEka/oU7RqfeLJstGV6d3y5Aknr1KDaT31h3hz9/uluMfXJDJYuVPX+1izo8ZKArcNLg9H94xSO7XXoC+7cOY1D0Wm4JLjY4kjH5jT14FR0tq8PXWM66LZ+7qdjXTBifxr2t6otfBZ5tyuP+zbZgs0qnBVRRW1DHtnV/5ausx9Dr42+XdeHZKD9nY3AoPT1CPQPkxPZ/9+a4xOpJ/7d9o7NA9tks0AcYW7wkWGrl+YHtev6kfPl56Fu/O584PtlBtsmhdljiPrUdLmfzaWnWhgtHAvBmDmD48RVbMtVKnmCAu7aluSXl9RabG1TSPhNEpFEVhUeMUXS+ZonM1l/aM473pA/H38WLNwWKmvfsrJ6pd7zBHT6AoCh9uyOL6tzZSWGmiU0wgC2YPZ5QsVGgzs8eo54L9sCuPzMIqjas5PwmjU+w8Vs6xE7X4+3h5bONFVzeiYySf3DmYUH9vduSUcc3c9W51NLM7qDNbeeTLnTy1YA8Wm8JlveL49r7hpEbJQoW21C0+mAndYlAUeMMFRkcSRqf4oaH9z7iuMfj5uF5vJ6Hq2z6ML+4eSnyIL4eLqrnqv+vYln1C67IEavPhq/+7nm+2qRtZ/3JpV16/sa9MidvJA+PU0dGCHblkFVdrXM25SRg1sNkUFu9unKKTja6urlNMEN/OGk73+GBKquu58e2NTf++QhurDhQx+bW17D1eQUSADx/dMYiZI1Pl/pAd9WgXwpjOUdgUeHftYa3LOScJowbbc8rIK68j0GiQeWs3ERPsyxd3D2Vcl2hMFhv3fbKNt1Ydkr1IDmax2njh5wymz9tEea2Z3omhLLx/BMPSIrUuzSPcNTINgC+3HHPqnnUSRg1+3qN2Ex/X1TXbr4szCzAaePvWAdw2NAmAOT9m8PAXO+WQPgfJK6vlhrc38saKQygKTBvcni/uHkJ8qJ/WpXmMIanh9GwXgsli4+ON2VqXc1YSRqgrexrDaGJ36dDtbrz0Ov5+ZQ/+dnk3vPQ6vt2ey7Vz17tc7y5Xs3RvAZe+uqZp2fbrN/Xl2at6YjTIiz1H0ul03DUyFYAPN2Q57QsxCSPgYGEVR0tq8DHopQeWG5s+PIWP7hhEeIAP6bkVXP7aWjYcOvepxaLl6sxW/vb9HmZ+uIWyGjM924XwwwMjZLuEhi7pEUu7UD9Kquv5ausxrcs5IwkjYEnDqGhEh0gCZVWPWxuWFsn3s4fTo10wpdX13Py/X3l3zWG5j9RG0nPLmfzaWt5fnwXA7cNT+OreoSRFBGhbmIczeOm5Y0QKAO+vz3LK/+8SRsCSvQUATJAO3R4hIcyfr+4ZxlV922G1Kfxj0T7u/GALpbJB9oJZbQpvrMhkyhvryCysIjLQyLzpA3nq8m4yLeckrhuQgL+PF5mFVWw47HwzAh4fRnlltew6Vo5Op+4vEp7B19uLl6b25v+m9MDHoGdZRiGX/Gc1G53wh9TZHS2p5vq3NvDCz/ux2BQmdo9hyR9GMkZ6OzqVIF9vrurbDlDPO3I2Hh9Gv+xTR0X924cRFWTUuBrhSDqdjluGJPHdfcNJjQqgoMLETe9s5OWlB6TzdzNYrDbeWX2Yia+sZsvREwQaDbx4XW/m3tyf8ADptu2Mbh2aDKizQfnlddoW8xseH0ZL9qhhJIfoea5u8cH8cP8IruufgE2B/yw7yLVz17tEPy+t7DtewdVvrufZxfuoM9sYmhrBjw9exLX9E2QTqxPrHBvEoJRwrDaFTzc51zJvjw6jyjpz07SMhJFn8/cx8MJ1vXnl+j4EGQ1szy7j0lfXMHfVIRklncJksfLvJfu5/LW17DpWTpCvgX9d05NPZw4mMdxf6/JEM9zasOfus03ZmK02jas5yaPDaF1mCRabQkpkgDRpFABM6duOn/8wklGdoqi32Pjnjxlc/eZ6DrrQiZn2sjyjgIkvr+a15ZlN94Z+eXgU1w9sL6MhFzKxeyyRgT4UVZpYtb9I63KaeHQYrT6o/kNI+x9xqvhQP96fMZDnr+1FkK+BnTllXPbqWl5eesBpNwza0+GiKmbM28Tt728hq6SGqCAjb07rx1u3DCAm2Ffr8kQLeXvpmdJHXcjw5dYcjas5yWPDSFGUplcFEkbit3Q6HVMHJKqrwjpHUW+18Z9lBxn371X8uPu4U+7TaGtVJgtzftzHxFdWs2J/Ed5eOu4elcqKP47mkp7STNiVXTcgEYBl+wqdpl+dx4bRoaJqcstq8THoGZwarnU5wknFhfjx3vSB/HdaP9qF+pFbVsu9n2xj2ru/csBNp+7qzFbeXXOYUc+v4K1VhzFbFUZ3juLnh0byxCVdZWO4G+gcG0SvhBAsNoXvduRpXQ4AHvu/atUBdVQ0OCUcfx+P/WsQzaDT6bi0ZxxjOkfz5qpDzF11iPWHSrjkP2uYOiCRWWPSSAhz/Zv39RYbn2/J4fXlBymoUF8tp0QG8NfJXRnbRRb4uJvr+iew61g5X27J4fbhyZrf9/PYZ+HGMJIpOtFcfj5ePHxxJ67rn8Czi/bx0558PtuUzVdbc7huQCL3jXbNUDJZrCzYnseryw9y7ITaPLZdqB8PjOvANf0SMHh57ASKW7uidzv+b9E+MvIrycivpGtcsKb1eGQY1Zmt/NqwpFvCSLRUYrg/c2/pz6Yjpfxn2QHWZZbw6a/ZfLklh2v7q6HkCsucCyvr+GRjNp/8epTiKrUVUlSQkdljOnDDoERp4+PmQvy9Gd0piiV7C1i067iEkRY2Hi7BZLERH+JLh2hZ0i0uzKCUcD65cwibs0r5zy8HWZtZzGebsvl8czYTusVy27BkhqSGaz798Vu7j5Uzb90RFu7Kw2xVF2LEBvty+4hkbhmSjJ+PhJCnuKxXHEv2FvDDrjwemdBJ0/+rHhlG6xuODbioY5TTPVEI1zMwOZyP7xzM5qxSXl12kDUHi/lpTz4/7cknLSqAa/oncHXfBGJDtFsGXVxlYuHOPL7Zlsvu3PKmj/drH8qM4SlM6hGLt0zHeZzxXWPw9daTVVLDnrwKerQL0awWjwyjxq4LQ9MiNK5EuJOByeF8dMdgDhRU8sH6LL7Zlsuhomqe/2k/L/y8n4FJ4UzoHsPE7rEOmcY7dqKGX/YW8POeAjZllTZ1kjDodVzWK44Zw1Pokxhq9zqE8wowGhjbJZrFu/NZuCtP0zDSKU68YaKiooKQkBDKy8sJDm6b+czKOjO9/74EmwIbnhhLXIgcfyzso7LOzI+78/lq6zE2ZZWe9rmUyACGpkUwJDWCXu1CaB/uj15/4aN0q03hSHE1u3PL2HSklPWHSjhaUnPa1/ROCOHqfglM7hVHRKA0BRaqxbuPc98n20gI82PNo2PafLaouc/jHjcy2pJ1ApsCSRH+EkTCroJ8vZk6MJGpAxPJLatlyZ58fkrPZ3NWKUeKqzlSXM2nv6rNKoOMBjrHBtE+wp/EMH/iQnwJ8fMmxM8bXx8v1JeMCnVmG2U1Zspq6ykor+NoaQ1HS2o4WFBJdf3p3SG89Dr6tw9jQvcYJnSLpX2E8y+qEI43pnM0vt56jp2oZd/xSrrFa7OQwePCqHGKbkiKTNEJx2kX6seM4SnMGJ5CRZ2ZTYdL2XC4hC1ZpezLr6TSZGHL0RNsOXrigq/h5+1Ft/hg+iaGMqxDBAOTwwny9W7D70K4Iz8fL0Z0iOSXfYUszyiQMHKUpjBKk64LQhvBvt6M7xbD+IZO8WarjczCKg4WVpFTWsOxEzUUVJioqDVTUWem1mxFhw6dDowGPaF+PoT4exMZaCQpwp+kcH/SogNJiwrEqxVTfcJzjesawy/7CvllXyGzx3bUpAaPCqM6s5U9eRUADEiSMBLOwdtLT9e4YM33eQjPNaazeirvzmNlFFWaNDlo1KPWcqbnlmOxKUQFGUkIk/tFQggBEBviS492wSgKrNhfqEkNHhVG27PLAOibGCr7i4QQ4hTjGvoPLt8nYWR323PUm8N924dpXIkQQjiXkQ2t0TYeKcGmwenGnhVGjSOj9qGa1iGEEM6mV0II/j5elNWY2Zdf4fDre0wYFVTUcby8Dr0Oemq4y1gIIZyRt5eegcnqwq4NDS3THMljwmhPntqPq0N0IAFyOJgQQvzOsIYWaY1bYBzJY8Jo33H1VE5ZPiuEEGfW2K/z18OlWKw2h17bY8Job8P+om4SRkIIcUbd40MI8jVQabKQkV/p0Gt7TBjtO66GkYyMhBDizLz0uqZO7jtyyhx6bY8Io5p6C0dKqgEJIyGEOJdeCeoCr13Hyhx6XY8IowMFVSgKRAYaNWlzIYQQrqJXQigAu46Vn/sL25hHhNGR4ioAOkQHaFyJEEI4t94NYXSgoJKaeovDrushYaQeMpYSKWEkhBDnEhviS3SQEZtCU2NpR3BIGJlMJvr06YNOp2PHjh2OuORpjhSr94uSIySMhBDifBrvG+3JddxUnUPC6NFHHyU+Pt4RlzqjrIYwkpGREEKcX4foIAAONzx3OoLdw+jHH39kyZIlvPjii/a+1BkpitI0MpIwEkKI80uLUp8rDxVVOeyadu2LU1BQwMyZM/nuu+/w9/c/79ebTCZMJlPTrysqWj9feaLGTJVJvQmXGH7+GoQQwtOlRgUCcKjQDUZGiqIwffp07rnnHgYMGNCs3zNnzhxCQkKaHomJia2uI7+8DoDIQB98vb1a/ecJIYS7axwZ5VfUNb2Yt7cWh9Hjjz+OTqc75yMjI4PXXnuNyspKnnjiiWb/2U888QTl5eVNj5ycnJaW9zsFFWoYRQf5tvrPEkIITxDq70NEgA8AR4ocMzpq8TTdI488wvTp08/5NampqSxfvpwNGzZgNJ6+yXTAgAFMmzaNDz744He/z2g0/u7rW6sxjGJDJIyEEKK5kiMDKKmu52hpNT0T7H/sTovDKCoqiqioqPN+3auvvso//vGPpl/n5eUxceJEPv/8cwYPHtzSy16wggr1HlRMsHReEEKI5mp8ziysMJ3nK9uG3RYwtG/f/rRfBwaqN8TS0tJISEiw12V/p6BSpumEEKKlGp8zCysdE0Zu34GhqOEvUnrSCSFE80U3jIyKHBRGDjvyNDk5GUVRHHW5JlV16kqQIF853VUIIZrr5MioziHXc/uRUXVDo79AOWpcCCGaLTrIsSMjtw+jxjXyEkZCCNF84Q1Lu0/U1Dvkeu4fRg3TdAESRkII0Wy+3mo8mCw2h1zP7cOo2iT3jIQQoqWMBrVjjcksYdQmas1WAPykFZAQQjSbsWlkZHXI4jO3DyOdTqd1CUII4XIaR0Y2BSw2CaNWa4wiB/xdCiGE2zAaTsaDI+4buX8YNaSRgqSREEI0l5f+5KyS1Sojo1ZrnKbTYL+tEEK4LMspAeRtsP/tDvcPo4a3NkkjIYRotnrryak5by/7R4Xbh1HjX6LZAcNMIYRwF3UNK5G99DoMehkZtVqInzcAZQ7aRSyEEO6gvNYMQLCvwSGrkt0+jEL9G8Ko4S9WCCHE+VU0PGc2vqC3N48Jo/IaCSMhhGiu4ip1NimsoUedvXlAGDm22Z8QQriD/PJaAGKDHXMwqfuHUdM9IxkZCSFEc+U3HDceI2HUNiIahphFVY45k0MIIdxBTmkNAAlhfg65ntuHUUK4P3DyL1YIIcT5HSqqAiAtOtAh13P7MEoMU8Po2IlajSsRQgjXYLUpHC6uBqBDlIRRm0iOPDkyqnfQIVFCCOHKDhVVUW+xEeDjRXyoTNO1idhgX4KMBiw2hcPFVVqXI4QQTm9HThkAPdqFnNYw1Z7cPox0Oh2dYoMA2J9fqXE1Qgjh/HY2hFHvxFCHXdPtwwiga5waRum55RpXIoQQzm/D4RIA+rUPc9g1PSKM+iSqf6Hbs8u0LUQIIZxcXlkth4uq0etgaFqEw67rEWHUr30oALtzy2URgxBCnMPag8UA9EoIdVhfOvCQMEqJDCAiwAeTxdZ0Y04IIcTv/Zh+HICxXaIdel2PCCOdTsewDpEArD1YpHE1QgjhnMprzKzNVEdGl/aMc+i1PSKMAC5qCKNVDUNQIYQQp/sx/Thmq0KnmEA6OKjzQiOPCaNRnaMAdcni8XLpxiCEEL/16aZsAK7pl+Dwa3tMGMUE+zIgSV1V93N6vsbVCCGEc0nPLWfXsXJ8vPRc21/CyK4m9YgFYOGu4xpXIoQQzuXdNYcBmNgjlohAo8Ovb3D4FTV0ee94nlu8j61HT5BZWOXwOVEhzsZitVFYaSK/oo7yWjOVdRaq6iyYrTbMVhs2RcFLr8fbS4ePl54Ao4EgXwMhft5EBRmJCjJiNHhp/W0IF5VVXM33O/MAuHtkqiY1eFQYxQT7MrZLNL/sK+SLLTn8+dKuWpckPIjFaiOrpJo9eRXsz68kt6yWvLJa8srqyK+ow2pTWvXnhwf4kBzhT0pkIKlRAaREnnz4ektQibN7fUUmNgXGdI6iR7sQTWrwqDACuH5g+6Ywemh8R/x9PO6vQDiAzaZwoLCSzVkn2JtXzt68CjLyKzGdY9O1t5eO6CBfQv29Cfb1JsBowGjQY/DS4aXTYbEpWGw2TGYblSYLlXUWymvqKa6qp95qo7S6ntLqerb9ptOIl15Hx+hAeieE0isxhN4JoXSODcLby6Nm6cVZpOeW8/W2YwDcP66jZnV43DPx2C7RtA/3J7u0hi+3HOO2YclalyTcgMVqY+/xCjYdKeXXI6Vszio941H3/j5edIkNoktcMEnh/rQL8yM+1I92oX5EBhovqEOyoiiU15rJLaslq7iGI8VVHC6u5khxNYeLqimvNZORX0lGfiWfb8kBwMegp1e7EIakRjA0LYJ+7cPw85HRk6dRFIVnFu5FUeCK3vEO7UX3WzpFUVo3N2BHFRUVhISEUF5eTnBwcJv9uR9tPMpfv0snIcyPFX8cLa8QxQUpqKhj1f4iVuwvZO3BYipNltM+7+ftRf+kMHonhtAtLoRu8WoA6R3Ukr9RfnkdO4+VsetYGTtzytl1rIyKutNr9fHS0ycxlCFpEYzoEEm/9qEY5OfC7X2xJYdHv9qFr7ee5Y+MtsvZRc19HvfIMKozWxnxr+UUV9Xzjyk9uHlIUpv92cJ92WwK23PKWLavgJX7i9h7vOK0zwf5GhiUHM6gFPXRo12IU77QURSFrJIaNh8pZcPhEjYcKiG/ou60rwnyNXBRx0hGd4pmVOcoYoJ9NapW2MuxEzVMemUNVSYLj03qwr2j0+xyHQmj8/hgfRZPf7+HyEAjq/40mgCjx81YimawWG1sOlLKj+n5/Lwnn8JKU9PndDq1meToTlGM7hxFr4RQhx1E1pYUReFoSQ0bD5ew/lAJaw4WceI3U4xd44IZ3zWai7vF0LNdCDqd632f4iSL1cYt/9vEhsMl9Gsfypf3DLPb/10Jo/Oot9gY/9IqsktruH9sBx6Z0LlN/3zhukwWK+syi/kpPZ+lewtOe2IOMhoY3SWasV2iGNkxSpP9GPZmtSnsOlbGyv1FrDxQxK5jZZz6LBEX4sv4rjFM6B7D4JQIfAzON/oT5/bc4n28vfow/j5eLHrgIlIiA+x2LQmjZvhx93Hu/WQb3l46Fj9wER1jgtr8GsI11FtsrDpQxA+78li2r5CqU+7/hPl7M6FbLJN6xDKsQ4TH7ecpqTKx6kARS/cWsOpAETX11qbPNYbzhG4xjO4cRZCv444cEBdmwY5cHpy/A4A3burHZb3s2xBVwqgZFEXhzg+2sCyjkIHJYXx+11CH31wW2rHaFH49XMKCHXn8mH78tJv60UFGJvVQA2hQcrjczG9QZ7ay/lAxS/cWsHRvIcVVJ6ctvb10DEmNYEL3WC7uGkNsiNxncjarDhRx5webMVsV7h2dxmOTutj9mhJGzXTsRA0TXl5NTb2VRyd15r7RHexyHeEcFEVhR04Z3+/MY9Gu46fdA4oOMjK5VzyX9Yqjb2KovDA5j8YFHUv3FrBkbz6Hi6pP+3zvxFAmdIthYvcYOkTLrIPWNmeVcsv/fqXObOOyXnG8ekNfh9zjlDBqgc83Z/PY17vx0uv48p6hmq61F/ZxoKCS73fk8f3OPLJLa5o+HuLnzSU9YrmiTzyDUyJccgGCszhUVKUG0558tuecfp8pNSqACd1imdA9hj4JEvSOtnJ/Ifd+vI1as5XRnaN4+5YBDrvXJ2HUAoqicP9n2/lh13HiQnxZMGs40bKU1eXllNawcFce3+/IIyO/sunjft5eXNwthit6xzOyU5TcgLeDwoo6ftlXyJK9+azPLKHeerLzRFSQkYu7xTChWwzD0iLl79/Ovtueyx+/3InFpnBRx0jevmWAQzc4Sxi19Fp1Zqa8vo7DxdX0Tgjh87uHSj8vF1RYWcfiXcf5fmfeaW1xvL10jOoUxRV92jG+a7S0gXKgyjozK/cXsWRvASsyTl8cIgsg7MdstfGvHzN4d+0RQO2w8OJ1vR0e/hJGFyCruJop/11HWY2ZCd1ieGNaP6fctChOd6K6np/25LNwZx4bD5fQ2G9Up4OhqRFc0TueST1iCfX30bZQgcliZePhUn7eoy6bL6o8fQHE0LRIxnaOYkyXaJIi7Lfc2N3lltXyh/k72JRVCsA9o9J4dGJnTaZHJYwu0K+HS7jlf5uot6o3+f5zfR9ZSeWEKuvMLN1bwMKdeaw5WIzllI7XfduHMrlXPJN7xUnnACdmsynsOFbGkj3qfabDxacvgEiNCmBM52jGdolmYHK4TOc1g6IofLYph+cW76PKZCHQaODF63o3neWmBQmjVli2r4B7Pt6K2aowuVcc/57a2+P2ljij2noryzMKWbgzj+X7C6k/pQN2t7hgLu+tBlBiuL+GVYoLlVlYybJ9hazYX8iWrBOnvcAI8PFiRMdIxnSOZmSnKLv0UHN16bnlPPPDXjYdUUdD/ZPCePG63nbd0NocEkat9POefGZ9sg2LTWFoagRv3dqfYJnPdjiTxcqaA8Us3JXH0r0Fp224TI0K4Ire8UzuFS8HJbqZijozaw8WsyKjkBX7i07bzwSQHOHPsA6RDEuLYGhqhFt2wmiu7JIa/rPsIN9sP4aigK+3nj9O6MyM4SlOsTrUKcJo0aJFPPPMM+zatQtfX19GjRrFd9991+zfr2UYAaw+UMS9H2+lut5Kp5hA5t7cn9QoedKzt9p6K6sPFvHznnx+2Vtw2mbUhDA/Lu8dz+W94ukaFyQ90jyAzaawJ6+CFfvVUdOuY+W/O4iwS2wQw9IiGd4hgkEp4R6xECIjv4K5Kw+xcNfxpr+PKX3i+dOkLrRzopGj5mH09ddfM3PmTJ577jnGjh2LxWIhPT2dqVOnNvvP0DqMQB36znh/M0WVJgKNBp6/theX9rRv+wxPVF5rZkVGIT+l57PqQBG15pMjoMbNqJf3jqNPYqgEkIerqDOz6XAp6w+VsP5Q8WnL9gH0OugSG0z/pLCmR0KYn1v8v6mtt7J493E+25TNlqMnmj4+qlMUD43vSF8n3COpaRhZLBaSk5P5+9//zh133HHBf44zhBGo59bc/+n2ppUpUwck8JfLuhHi5/6vvuypsLKOX/YW8tOefDYcKsZsPflfsV2oHxO7xzKxewwDksOdYrpBOKfiKhMbGoJp/aESjpbU/O5rooKM9E0MpWe7EHq0C6F7fLDL7CWsMllYfaCIRbuPsyKjsGmq2kuvY2L3GO4b3UGzo8KbQ9Mw2rRpE4MHD+a9997j1VdfJT8/nz59+vDCCy/Qo0ePs/4+k8mEyXRybriiooLExETNwwjUlusvLNnP26sPoygQG+zLM1d25+JuMW7xissRrDaFncfKWNlwH2B3bvlpn+8YHcikHrFM7B5L9/hg+XsVF+R4eS3bjpax9egJtmarx76f+kKnUVSQkR7xwXSJC6ZDVCBp0YGkRgVofm+4ymRhV04ZW46eYO3BYrZln76YIzHcjxsGtue6/gkuEaiahtH8+fO58cYbad++PS+99BLJycn8+9//ZsmSJRw4cIDw8PAz/r6//e1v/P3vf//dx50hjBptzirl0a92caRhGerwDhE8eVk3usY5R33O5kR1PasPFrEio5BVB35/Tk7vhBAmNgRQmtyPE3ZQZ7ay65h6wu2evArSc8s5VFSF7SzPfNFBRtKiAmkf7k9cqC/xoX7Eh/ip74f4tUn3AkVRKKtRj4rPKa0hs7CKA4VV7M+v4GBhFb99Vk6K8GdS91gu7RlHrwTXOk/KLmH0+OOP869//eucX7Nv3z62bdvGtGnTeOutt7jrrrsAddSTkJDAP/7xD+6+++4z/l5nHhmdqrbeyqvLD/K/NUeot9rQ6+Dy3vHMHtPB44+hMFms7MguY92hEtYeLGJHTtlpP/RBvgZGdlQPoxvVOYroIOd/ZSfcT229lX35FezJLedAQRWHitRHQYXpvL/X38eLYF9vgv0MDW+9CfY14O2lR6cDvU6HTkdTYNSYLFTUWaioNVNRZ6ai1kJZbT11ZttZr9Eu1I8+7UMZlhbBRR2iaB/hutsV7BJGRUVFlJSUnPNrUlNTWbduHWPHjmXNmjWMGDGi6XODBw9m/PjxPPvss826nrPcMzqb7JIa/vnTPhbvzm/62MTuMdw2LJmhqREu9erlQpmtNvbkVbD+UDEbDpWwOav0dz9kXWKDGN05mjGdo+iXFCZdLYTTqqgzc7iomkOFVRw7Ucvx8lryyus4XlbL8fK601oZtYXIQCPtwvxIiwqgU0wQHaMD6ZkQ4lYv0pr7PN6iBl1RUVFERUWd9+v69++P0Whk//79TWFkNpvJysoiKSmpJZd0au0j/PnvtP6k55bz+vJMftqTz897Cvh5TwFpUQHcMLA9k3vHERfiPMssW6u8xsy2nBNszTrBlqOl7MwpP23lG0BkoA9D09Q9IKNkg6JwIcG+3vRJDKVPYugZP19RZ6asunGEo74tr1VHO2abDUVRp+AUBRRAURpGUqeNotRRVUywr/S/PIXdlnY/9NBDfPXVV7z33nskJSXxwgsvsHDhQjIyMggLa97yQ2cfGf3WwYJKPtiQxbfbcqluWPGi08HA5HAmdY9lTJdozXdDt0SVycLehjn29NxydueWk1n0+/nsED9vBqeEMywtgmEdIukYHegRo0IhxPlpvs/IbDbzxBNP8NFHH1FbW8vgwYN55ZVX6N69e7P/DFcLo0ZVJgvfbc9lwY5cNmedOO1zSRH+DE2NYEByOAOSwkiK8Nf8idtksXKkuJrMwioOFlSRWVTFvuMVHCmu/l3wgLr7vX9SOAOSwxiQFEZaVKCcTyOEOCPNw6gtuGoYnSqvrJbFu4+zYn8hm46U/m6JaZCvga6xwXSNCyIlMoDkyACSIgKIDfZtszNH6sxWiipN5JXVkldeS15ZHbllteSeqOVoSTXZpTVnXVkUF+JL9/iQhv0ZwfRKCCUqyHNbrwghWkbCyAlVmSxsPFTC5qOlbMk6we5j5acdOvZbQb4GogKNBPoa8PfxIsDHQIDRQIDRC6PBC6tNwWKzYbEqDe8r1JmtlNeq89hlNebzrtppFOxroEN0YNOjU0wQPdqFEOnBPb+EEK1nlwUMonUCjQbGd4thfLcYAOotNg41TIntL6gkq7iaoyU1ZJfWUFNvpbLOQmVd26ze8fHSE9+4Z6LhkRDqR0KYHx2iA4kKMmo+XSiE8FwSRhryMejpGhf8uw2ziqJQZbJQUGGiuMpETb2FapO16W21yYLJYsNLr8Og1+Hl1fBWr8do0BPq702onw+h/t6E+HkT6u9NoNEgYSOEcFoSRk5Ip9MR5OtNkK+3HI0ghPAIsvtQCCGE5iSMhBBCaE7CSAghhOYkjIQQQmhOwkgIIYTmJIyEEEJoTsJICCGE5iSMhBBCaE7CSAghhOYkjIQQQmhOwkgIIYTmJIyEEEJoTsJICCGE5py6a3fjuX8VFRUaVyKEEOJCND5/n+8cV6cOo8rKSgASExM1rkQIIURrVFZWEhISctbPO/Wx4zabjby8PIKCgux6MFxFRQWJiYnk5OS4xfHm5+JJ3yt41vfrSd8ryPfrKhRFobKykvj4ePT6s98ZcuqRkV6vJyEhwWHXCw4Odql/5NbwpO8VPOv79aTvFeT7dQXnGhE1kgUMQgghNCdhJIQQQnMSRoDRaOTpp5/GaDRqXYrdedL3Cp71/XrS9wry/bobp17AIIQQwjPIyEgIIYTmJIyEEEJoTsJICCGE5iSMhBBCaE7C6AwWLVrE4MGD8fPzIywsjClTpmhdkl2ZTCb69OmDTqdjx44dWpdjF1lZWdxxxx2kpKTg5+dHWloaTz/9NPX19VqX1mbeeOMNkpOT8fX1ZfDgwWzatEnrkuxizpw5DBw4kKCgIKKjo5kyZQr79+/XuiyH+Oc//4lOp+Ohhx7SupQ2J2H0G19//TW33HILM2bMYOfOnaxbt46bbrpJ67Ls6tFHHyU+Pl7rMuwqIyMDm83GW2+9xZ49e3j55ZeZO3cuf/7zn7UurU18/vnnPPzwwzz99NNs27aN3r17M3HiRAoLC7Uurc2tWrWKWbNmsXHjRpYuXYrZbGbChAlUV1drXZpdbd68mbfeeotevXppXYp9KKKJ2WxW2rVrp7z77rtal+IwixcvVrp06aLs2bNHAZTt27drXZLDPP/880pKSorWZbSJQYMGKbNmzWr6tdVqVeLj45U5c+ZoWJVjFBYWKoCyatUqrUuxm8rKSqVjx47K0qVLlVGjRikPPvig1iW1ORkZnWLbtm3k5uai1+vp27cvcXFxXHLJJaSnp2tdml0UFBQwc+ZMPvroI/z9/bUux+HKy8sJDw/XuoxWq6+vZ+vWrYwfP77pY3q9nvHjx7NhwwYNK3OM8vJyALf4tzybWbNmcdlll532b+xuJIxOcfjwYQD+9re/8eSTT/LDDz8QFhbG6NGjKS0t1bi6tqUoCtOnT+eee+5hwIABWpfjcJmZmbz22mvcfffdWpfSasXFxVitVmJiYk77eExMDPn5+RpV5Rg2m42HHnqI4cOH06NHD63LsYv58+ezbds25syZo3UpduURYfT444+j0+nO+Wi8pwDwl7/8hWuuuYb+/fszb948dDodX375pcbfRfM093t97bXXqKys5IknntC65FZp7vd7qtzcXCZNmsR1113HzJkzNapctIVZs2aRnp7O/PnztS7FLnJycnjwwQf55JNP8PX11bocu/KIdkBFRUWUlJSc82tSU1NZt24dY8eOZc2aNYwYMaLpc4MHD2b8+PE8++yz9i611Zr7vU6dOpWFCxeedk6U1WrFy8uLadOm8cEHH9i71DbR3O/Xx8cHgLy8PEaPHs2QIUN4//33z3m+iquor6/H39+fr7766rSVn7fddhtlZWUsWLBAu+LsaPbs2SxYsIDVq1eTkpKidTl28d1333HVVVfh5eXV9DGr1YpOp0Ov12MymU77nCtz6vOM2kpUVBRRUVHn/br+/ftjNBrZv39/UxiZzWaysrJISkqyd5ltornf66uvvso//vGPpl/n5eUxceJEPv/8cwYPHmzPEttUc79fUEdEY8aMaRrxukMQAfj4+NC/f3+WLVvWFEY2m41ly5Yxe/ZsbYuzA0VRuP/++/n2229ZuXKl2wYRwLhx49i9e/dpH5sxYwZdunThsccec5sgAg8Jo+YKDg7mnnvu4emnnyYxMZGkpCReeOEFAK677jqNq2tb7du3P+3XgYGBAKSlpTn0QENHyc3NZfTo0SQlJfHiiy9SVFTU9LnY2FgNK2sbDz/8MLfddhsDBgxg0KBBvPLKK1RXVzNjxgytS2tzs2bN4tNPP2XBggUEBQU13RcLCQnBz89P4+raVlBQ0O/uhQUEBBAREeF298gkjH7jhRdewGAwcMstt1BbW8vgwYNZvnw5YWFhWpcmWmHp0qVkZmaSmZn5u7B1h5nq66+/nqKiIp566iny8/Pp06cPP/300+8WNbiDN998E4DRo0ef9vF58+Yxffp0xxck2oRH3DMSQgjh3Nxj0lwIIYRLkzASQgihOQkjIYQQmpMwEkIIoTkJIyGEEJqTMBJCCKE5CSMhhBCakzASQgihOQkjIYQQmpMwEkIIoTkJIyGEEJqTMBJCCKG5/wdkl+K+yEovkAAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.plot(*paths['path'].T)\n", "plt.plot(*center[0].T, 'x', color='C0')\n", "plt.plot(*center[1].T, 'x', color='C1')\n", "plt.gca().set_aspect(1);" ] }, { "cell_type": "markdown", "id": "1e7c2bfa", "metadata": {}, "source": [ "As with any jax function we can evaluate the jacobian with respect to the inputs:" ] }, { "cell_type": "code", "execution_count": 22, "id": "80b15cdf", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(Array([[ 17.77701254, 17.91892068],\n", " [-39.5424976 , -39.5436664 ]], dtype=float64, weak_type=True),\n", " Array([[-0.75541393, -0.75484016],\n", " [-0.2485473 , -0.24775855]], dtype=float64, weak_type=True))" ] }, "execution_count": 22, "metadata": {}, "output_type": "execute_result" } ], "source": [ "grad, paths = jax.jacobian(centroids, argnums=(0, 1), has_aux=True)(0.01, 0.5)\n", "grad" ] }, { "cell_type": "code", "execution_count": null, "id": "de4f8bbf", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "lensing", "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.10.9" } }, "nbformat": 4, "nbformat_minor": 5 }