{ "cells": [ { "cell_type": "markdown", "id": "optional-abraham", "metadata": {}, "source": [ "# Introduction\n", "\n", "The dask library provides parallel versions of many operations available in numpy and pandas. It does this by breaking up an array into chunks.\n", "\n", "The full documentation for the Dask library is available at https://docs.dask.org/en/latest/" ] }, { "cell_type": "code", "execution_count": 1, "id": "applicable-bunny", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\n", "\n", "\n", "\n", "\n", "
\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
Array Chunk
Bytes 800 B 200 B
Shape (10, 10) (5, 5)
Count 4 Tasks 4 Chunks
Type float64 numpy.ndarray
\n", "		\n", "\n", "\n", " \n", " \n", " \n", " \n", "\n", " \n", " \n", " \n", " \n", "\n", " \n", " \n", "\n", " \n", " 10\n", " 10\n", "\n", "
" ], "text/plain": [ "dask.array" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import dask.array\n", "\n", "dask.array.ones((10,10), chunks=(5,5))" ] }, { "cell_type": "markdown", "id": "aging-lightweight", "metadata": {}, "source": [ "Here I've made a 10x10 array, that's made up of 4 5x5 chunks.\n", "\n", "Dask will run operations on different chunks in parallel. It evaluates arrays lazily - data is only computed for the chunks needed, if you don't explicitly load the data (e.g. by using `.compute()`, or plotting or saving the data) it builds up a graph of operations that can be run later when the actual values are needed." ] }, { "cell_type": "code", "execution_count": 2, "id": "hybrid-pakistan", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAa0AAAGcCAYAAABqcP8/AAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdd3hUdfY/8PedkkIKIQESCKFKQg8ovQgEBFmUBTEqCMIqgqsr9hX0a3sU0d31i/501xVBIY0EQbpIS5AiLWBCTYHQCS2NZNImM+f3B9/Mmk2buVM+907O63nyB1Pufc8c7jl3Zu7ckYiIwBhjjKmARnQAxhhjzFo8tBhjjKkGDy3GGGOqobP1DqtXr3ZGDrcydOhQtGvXTnSMGi5fvowDBw6IjqFokiQhOjpadIxaDhw4gMuXL4uOoWhhYWEYMmSI6Bi1cL9snM39kmxgNpsJAP818rd69WpbnlaXSExMFP68qOFPiaZOnSr8eVH6X3R0tOgy1cL90ro/W/ulrLcHV69eDSLiv//6M5vNcp5OlxL9HCn1LzExUXRpGhQdHS38OVLq39SpU0WXp0HcL+v+k9sv+TMtxhhjqsFDizHGmGrw0GKMMaYaPLQYY4ypBg8txhhjqsFDizHGmGrw0GKMMaYaPLQYY4ypBg8txhhjqsFDizHGmGrw0GKMMaYaPLQYY4ypBg8txhhjquHSoXXx4kVMmjQJd+7cqXF5WloaJk6ciICAAPj5+WHs2LHYv3+/w9Zr7fIXLFiApKQkh63XnYiqHQD89NNPCA8Ph05X/8+/ce3qV1/tAOueW7l4u7OPqLpZu3xRdXPZ0EpLS0P//v0xbtw4+Pv7Wy4/dOgQhg4dCj8/P5w5cwbnz59H586dMWrUKGzfvt3u9dqy/GeffRYLFy7EO++8Y/d63Ymo2p07dw6TJk3CwoULcePGjQZvy7WrW321s+W5lYO3O/uIqpsqtjmyQfWPmtn6o11FRUXUrl07mjdvXo3LTSYT9ezZk9q0aUOlpaWWy6uqqigiIoLCwsKovLzcpnXZu/y0tDSSJImSkpJsXp/c58cVqn8E0laiakdENG3aNFq8eDEZjUYKDQ0lrVbb4O3tqZ3c58cVpk6dKutHDuurHZHtz60tXL3dyX1+nM3R/ZLIuXWTs3wR/dIlQ+vtt98mnU5HV69erXF5SkoKAaAXX3yx1n3ef/99AkBr1qyxaV2OWH50dDS1a9eOjEajTetzx6ElqnZEVKPhWbuByq2dOw6t+mpHJO+5tZartzt3G1qi6iZ3+a7ul05/e5CIsGzZMgwaNAht27atcV1ycjIAoH///rXuV33Zrl27ZK9b7vKnTJmCK1euYMuWLbLX7Q5E1g4AvL29bb4P1+6uhmoHyHturcXbnXwi6yZ3+a6um9OHVnp6Om7cuIHIyMha12VkZAAA2rVrV+u60NBQAEBWVpbsdctdft++fQEA27Ztk71udyCydnJx7e5qqHbOxtudfCLrJper6+b0oXXy5EkAdf8HLiwsBAD4+PjUus7X1xcAUFBQIHvdcpdfvWFVZ2+qRNZOLq7dXQ3Vztl4u5NPZN3kcnXdnD60cnNzAQDNmze36X5EBACQJMnhmRpbvr+/PyRJsmRvqpRau4Zw7e6SWztn4+2uYUqtW0NcXTenD63y8nIAgF6vr3VdQEAAAMBgMNS6rvqy6tvIYc/ydTodysrKZK/bHYisnT24dg3Xztl4u5NPZN3s4cq6OX1oeXl5AQCMRmOt67p16wYAuHLlSq3rrl69CgAIDw+XvW57ll9VVeX0Dz2VTmTt7MG1a7h2zsbbnXwi62YPV9bN6UOrTZs2AICioqJa140ePRoAcPTo0VrXVV82ZswY2euWu/w7d+6AiCzZmyqRtZOLa3dXQ7VzNt7u5BNZN7lcXTenD61evXoBqHuva+TIkejRowfWrFljeVkMACaTCYmJiQgLC8PEiRNlr1vu8qv3BquzN1UiaycX1+6uhmrnbLzdySeybnK5um5OH1qRkZFo3bo10tPTa69co8Hy5cuRn5+PP/3pT7h+/Try8vLwwgsvIDs7G99++63l5XK1GTNmQJIknD9/vtF1y1k+cPcUKgAwbtw4mY/aPYisnVxcu7saqp0cvN25hsi6yeXqujl9aEmShDlz5uDQoUO4du1aresHDx6MX3/9FUVFRYiIiEDHjh2RnZ2N3bt3Y/z48bVun5ubC19fX7Rv396q9du6fABYt24dQkNDhbxSUBLRtdu8eTMkSYIkSbh69SpMJpPl38uWLavzPly7uxqrna3PLW93riG6bqrY5lxx2o3CwkIKDQ2t81xatigoKCBvb2+aM2eOXctpSPW5tFatWmXzfd3xNE5NpXbueBqnplI7dzuNU1Opm6LPPUhEdOzYMQoKCqKvvvrK5vtWr3vmzJkUHBxMubm5spbRmHPnzlHnzp3prbfeknV/dxxaRE2jdu44tIiaRu3cbWgRNY26Kfbcg9X69euH1NRUbN26tc7fh2nMjRs3kJOTg127diEkJMQJCYFvvvkGixYtwqJFi5yyfLXi2qkX106duG71c84viNWjY8eO2Lx5s6z7hoSEYN++fQ5OVNOnn37q1OWrGddOvbh26sR1q5tLf7mYMcYYswcPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqYbqh1ZpaanoCEwmrp06GY1GGI1G0TGYDO6wzcn6aZIDBw6AiBydRZaVK1di+vTp0Ov1oqOowurVq0VHAABUVVVh1apVmDlzpugoAICDBw+KjtCgK1euKKZ2x44dg0ajQd++fUVHAXD3ubH25+RF4H7pWLKG1pIlSxydwy4//fST6Aiq8fjjj4uOUIPc3wtqag4cOIADBw6IjqFYSh5a3C8dSyKl7ALI8MMPP+Cxxx7DY489hqSkJNFxmA0effRRrF27FmvXrsUjjzwiOg6zUmlpKYKCggAAt2/fho+Pj+BEzFru0i9V/ZlWQkICAGDjxo0oKSkRnIZZq7i42PIKKz4+XnAaZouNGzeioqICFRUV2LJli+g4zAbu0i9VO7Tu3LljeZlbWVmJjRs3Ck7ErLV+/XrLB/mbNm1CUVGR4ETMWrGxsdBoNNBqtYiLixMdh1nJnfqlaofW2rVrUVVVBQCQJAmxsbGCEzFrxcXFQZIkAIDJZMKGDRsEJ2LWKCgowPbt22EymVBVVYWtW7ciPz9fdCxmBXfql6odWrGxsTUa344dO5CXlyc4FWvM7du3sWvXLphMJstlat6AmpIffvgBZrPZ8m8iwrp16wQmYtZyp36pyqF18+ZN/PLLLzUaHwCsWbNGUCJmrf8+bNtsNiMlJQU3btwQlIhZ6793LoiIdzhUwN36pSqHVmJiIjSamtGJCDExMYISMWvFxMTU+Z0VtW5ATcW1a9fw66+/1nilZTabsWfPHly9elVgMtYYd+uXqhxaMTExtfYazGYzDhw4gCtXrghKxRpz+fJlHD58uEbjA+7WTq0bUFOxatWqWo0PALRaLX744QcBiZi13K1fqm5o5eTk4NixY3Xuret0OiQmJgpIxayRkJAArVZb63IiwpEjR3DhwgXXh2JWqavxAXc/H1m5cqWARMwa7tgvVTe0Vq1aBZ2u7hN5VFVV8R67gtXX+AD1bkBNwdmzZ3H8+PE6Gx8RIS0tDdnZ2QKSsca4Y79U3dCKjY2t92SdRIQTJ07g9OnTLk7FGpORkYHTp0/Xew42o9HIe+wKFRcX1+C56vR6Pe9wKJQ79ktVDa3jx48jMzOzwdt4eHgo5sSi7D/i4+MbPUlnRkYGTp486aJEzFpxcXENntXdaDRixYoVrgvErOKu/VJVQ8uavbnKykreY1egmJgYq37OgvfYleXYsWM4d+5co7fLyclBWlqaCxIxa7lrv5R1lndRjEYjoqOja1x24MABhIeHW07iWS03Nxdt2rRxZTxWj2vXrmHQoEEYNGiQ5bK8vDxkZ2dj8ODBNW5bWVnp6nisAefPn6+1zWVlZQEAwsPDa1yek5OjmJ8rYe7bL1V9lncigkajwerVq2sVhylbUlISpk2bVuvwd6Z80dHR0Gg0qj5TeFPkLv1SVW8PMsYYa9p4aDHGGFMNHlqMMcZUg4cWY4wx1eChxRhjTDV4aDHGGFMNHlqMMcZUg4cWY4wx1eChxRhjTDV4aDHGGFMNHlqMMcZUg4cWY4wx1eChxRhjTDV4aDHGGFMNHlqMMcZUg4cWY4wx1eChxRhjTDV4aDHGGFMNHlqMMcZUg4cWY4wx1eChxRhjTDV4aDHGGFMNHlqMMcZUg4cWY4wx1eChxRhjTDV4aDHGGFMNHlqMMcZUg4cWY4wx1eChxRhjTDV4aDHGGFMNHlqMMcZUQyc6gK3Ky8tx69YtGAwGlJSUAACys7Nx4sQJ+Pj4oGXLlvD39xecktXlzp07uH37NgwGA7KysgAAR48ehY+PD3x9fdGyZUt4eXkJTsn+m9lstmxzBQUFuH37NiRJwtGjR9GiRQv4+PigVatW0Gh4H1hp3LFfSkREokPU5+rVq0hOTsaRI0eQmZmJrKwsXLx4EY1FDgoKQnh4OLp3745evXph9OjR6NOnD29ULmI2m5Geno6UlBScPHkSZ86cQVZWFvLz8xu8n0ajQfv27REREYGIiAgMGDAAUVFRaNu2rYuSM4PBgH379mHfvn2WumVlZaGioqLB+3l6eiIiIsKy3Q0fPhzDhw9Hs2bNXJScNZV+qaihRUTYt28fVq1ahV27diErKwuenp7o168funfvjvDwcISHh6NNmzaWvfOAgAAYDAbL382bN5GTk4OMjAxkZWXh2LFjuH37NoKCgjBy5Eg88sgjmDJlCm9MDmYwGPDjjz9i3bp1+OWXX5Cfn49WrVqhX79+iIiIQLdu3dCpUye0bt0aPj4+lr/CwkKUlJTAYDAgNzcXWVlZyMzMxJkzZ/Dbb7+hsrISERERGDt2LKZNm4ahQ4dCkiTRD9etXLx4EXFxcdi6dSsOHz4Mo9GIbt26oVevXggPD0dERAQ6duwIHx8fBAQEwMfHB8DdmlfX78KFC5ZGefLkSWRmZsLDwwODBg3ChAkT8OSTT6J9+/aCH6l7abL9khTg0qVL9O6771KnTp0IAEVGRtJbb71FO3bsoNLSUruWbTKZKC0tjZYsWUITJkwgvV5Pfn5+NHv2bNq3b5+DHkHT9csvv9CsWbPI19eXPDw8aOLEifT5559Teno6mc1mu5ZtMBho+/bttGDBAurduzcBoC5dutD7779PV65ccdAjaJrKysrou+++o5EjR5IkSdS6dWuaM2cOxcfHU25urt3Lv3btGsXFxdHTTz9NLVu2JI1GQ6NHj6YVK1ZQeXm5Ax5B09XU+6XQoZWdnU3PPPMMeXh4UJs2bej111+n48ePO3WdN2/epC+++IL69+9PAGjEiBG0detWp67T3ZjNZtqyZQsNHTqUANCAAQPoyy+/pFu3bjl1vWlpafTqq69SSEgIeXp60ty5c+ncuXNOXae7KS4upn/84x/Upk0b8vT0pOjoaNq4cSNVVlY6bZ2VlZW0YcMGmjp1Knl6elJoaCgtWbKESkpKnLZOd8T98i4hQ+v69es0a9Ys0mq11LVrV1q+fDlVVFS4PMe+fftowoQJlsb766+/ujyD2uzdu5fuvfdekiSJJk6cSPv373d5hoqKClq6dCl16dKFdDodPf3003Tz5k2X51ATo9FIn3/+OQUFBZGvry+9/vrrdO3aNZfnuHr1Kr3yyivk4+NDrVq1oq+++oqqqqpcnkNNuF/W5NKhVVVVRV9++SUFBARQhw4dKCEhQRH/YVNTU2ns2LEkSRI988wzTn/FoEY3b96k2bNnkyRJ9OCDD9Jvv/0mOhIZjUaKiYmhdu3aUYsWLejrr78mk8kkOpbi7N+/n/r27Uuenp60cOFCun37tuhIdOvWLXrjjTdIr9fTvffeS4cOHRIdSXG4X9bNZUMrJyeHBg4cSB4eHrRw4UIyGAyuWrXVEhMTKTQ0lIKCgmjDhg2i4yjGjz/+SIGBgRQWFkZr1qwRHaeW4uJiev3110mv19PQoUPpwoULoiMpQnl5Ob300kskSRI98MADlJmZKTpSLadPn6bRo0eTRqOh1157zalvU6oJ98v6uWRorVmzhgICAigyMpJOnz7tilXKdufOHXrmmWdIkiR69dVXm/RGVF5eTi+++CJJkkRz586l4uJi0ZEadOLECerVqxcFBgY2+Z2Oc+fOUf/+/cnPz49iY2NFx2nU999/T76+vjR48OAmv9PB/bJhTh1aJpOJXnvtNQJA8+bNo7KyMmeuzqFiY2MtG1FT/LwkNzeXBgwYQH5+fpSQkCA6jtVKS0stG9HChQvtPoJRjbZu3UrNmzenvn37KvLVVX3OnDlDvXv3phYtWtDOnTtFx3E57pfWcdrQqqyspCeffJI8PDxUsadXl4yMDOrSpQuFh4dTTk6O6Dguc+7cOerSpQt17dpVVU3v91asWEF6vZ5mzZpFRqNRdByXiY+PJ71eTzNnzlRV06tWWlpK06ZNI09PT0pKShIdx2W4X1rPKUPLYDDQgw8+SH5+frRjxw5nrMJlrl+/Tv369aM2bdo4/fBSJUhLS6OQkBC677776MaNG6Lj2GXr1q3k4+NDEydOtPv7K2rwxRdfkEajoddff13VrzBNJhO9/PLLpNFo6KuvvhIdx+m4X9rG4UOrsrKS/vCHP1DLli3pyJEjjl68EEVFRTRy5EgKDg6m7Oxs0XGcJisri1q3bk1RUVF0584d0XEc4uDBgxQYGEgPPfSQW7/i+uabb0iSJPrb3/4mOorDLFq0iCRJomXLlomO4jTcL23n0KFlNptp9uzZ1KxZM7f7zpPBYKAhQ4ZQ586dhXy/xdmuXr1KHTt2pIEDByr+gAtbHTp0iHx9fWnGjBmqfgVSnw0bNpBOp6MPP/xQdBSHe+edd0ir1SryqFV7cb+Ux6FD68033yQPDw/avn27IxerGDdv3qTw8HDq16+fIg9Blau4uJj69OlD3bt3V8R3eJxhy5YtpNfr6X/+539ER3Goffv2kZeXF/3lL38RHcVp5s6dS97e3nTgwAHRURyK+6U8DhtaGzZsIEmSaOXKlY5apCLl5ORQYGAgPf3006KjOMzMmTOpVatWdPHiRdFRnOrbb78lSZJoy5YtoqM4xK1bt6hdu3b0xz/+0a2/VF1VVUUTJkyg9u3bU15enug4DsH9Uj6HDK1Lly5RUFAQzZkzxxGLU7zNmze7zX+45cuXkyRJtHHjRtFRXOKpp56iwMBA1X8XyGw206RJkygsLMxtXx3/Xl5eHnXo0IEmTpyo+rd4uV/ax+6hZTKZaNiwYdSrV68mcYRWtVdffZV8fX1VfcLWzMxMatasGS1YsEB0FJcpLi6mbt260ahRo1Td/P73f/+X9Ho9HTx4UHQUl9m7dy/pdDpVH1HI/dL+fmn30Fq6dCnpdDpKT0+3O4yaVFRUUM+ePWnixImio8j2wAMPUGRkpFsfVVeX1NRU0mq1tGLFCtFRZLl8+TL5+vrSBx98IDqKy7399tvk7++v2oOhuF/a3y/tGlq3b9+mli1b0quvvmp3EDXas2cPSZJE69atEx3FZomJiSRJEv3yyy+iowjxwgsvUFBQkCrfWps6dSrdc889qvzysL1KS0upU6dONH36dNFRbMb90jH90q6hNXfuXAoNDXWb7/TIMWPGDOrQoYOqftjOYDBQaGioWx1MYquCggIKDg5W3VF3O3bsIAD0888/i44izMaNGwkA7d69W3QUm3C/dEy/lD20Ll68SHq9npYvXy575e7g6tWr5OXlRV9//bXoKFb74osvqFmzZnT9+nXRUYT697//TZ6enqr6FeThw4fThAkTRMcQbsyYMTR69GjRMazG/fIuR/RL2UPrhRdeoLCwMCE/RqY0f/7zn6l9+/aqOCN8ZWUldejQgV566SXRUYQrLy+n0NBQeuWVV0RHscru3bsJAO3du1d0FOGSk5MJgGJ+Ar4x3C//w95+KWtoXb9+nby9vVV9FI8jnT9/nvR6vSoOgV+2bBl5eHjQ5cuXRUdRhM8//5x8fHxU8cOf48ePp5EjR4qOoRjDhg1TxYFQ3C9rsrdfyhpan376KQUGBjbJD4LrM336dBoyZIjoGI3q378/zZo1S3QMxTAYDBQQEECfffaZ6CgNOnfuHEmSRJs3bxYdRTHWr19PGo1G8d+5435Zmz39UgMZ4uLi8MQTT8DLy0vO3d3SrFmzcODAAWRmZoqOUq8zZ84gNTUVs2bNEh1FMZo1a4ZHH30UcXFxoqM0KCYmBq1bt8b48eNFR1GMiRMnolWrVoiPjxcdpUHcL2uzp1/aPLSOHj2KEydOYObMmTavzJ2NHTsWoaGhSEhIEB2lXjExMWjfvj1GjhwpOoqizJw5E7/99htOnDghOkqdiAhxcXGYMWMGdDqd6DiKodPpEB0djZUrV4qOUi/ul3Wzp1/aPLSSkpLQtWtXDB482OaVuTONRoMnnngCiYmJoqPUKzExEU8++SQ0GlkvsN3WiBEj0LFjRyQlJYmOUqfU1FScO3cOM2bMEB1FcWbOnImsrCykp6eLjlIn7pd1s6df2ty9du3ahQkTJti8oqZgwoQJyMrKwsWLF0VHqeXcuXO4cOEC164OkiRh/Pjx2LVrl+goddq5cyfatGmDyMhI0VEUZ8CAAWjVqhV27twpOkqduF/WT26/tGloFRYWIj09HaNHj7ZpJU3FsGHD4OXlhd27d4uOUktKSgqaNWuGgQMHio6iSFFRUUhNTcWdO3dER6klJSUFUVFRkCRJdBTFkSQJo0aNQkpKiugotXC/bJjcfmnT0Nq9ezeISFGfiSQmJkKSJEiSZNUHnbbe3hZeXl4YMmQIkpOTHbpcR0hJScHw4cPh6ekpOgoAZdUNAEaPHg2TyYQ9e/Y4fNn2qKysxP79+zFq1CjRUSyUVruoqCjs2bMHVVVVDl+2PbhfNkxuv7RpaB09ehTdu3dHixYtbFqJMz3xxBMgIowZM8Ypt7fV0KFDcfToUacs2x5Hjx7FkCFDRMewUFrdWrVqhfDwcKSmpjpl+XKdOXMGpaWlGDp0qOgoFkqr3dChQ1FcXKy4I3e5XzZOTr+0aWhlZGQgIiLCphU0NRERETh79ixMJpPoKBZGoxE5OTno1q2b6CiKFh4ejqysLNExasjMzIRWq0WXLl1ER1Gsrl27QqPRKK523C8bJ6df2jS0MjMzuQiNiIiIQEVFhaIOxsjJyYHRaOTaNSIiIkJxe+uZmZno1KmTYt7WVSJvb2+EhYUpsna8zTVMTr+0emgREc6ePYvw8HBZ4ZqK6v+kStqAsrKyIEkSunbtKjqKokVERChubz0rK4u3OSt069ZNUdsc90vryOmXVg8tg8GAsrIyBAcH256sAVVVVUhKSsIDDzyAkJAQeHt7o3fv3vjiiy9gNptr3T4jIwOTJ09G8+bN4ePjgxEjRmDfvn31Lt/W29urefPm8PDwQF5entPWYavbt2+jWbNm8PX1ddgy3a1uANC6dWuUlJSgoqLCqeuxxe3bt3mbs0KrVq0Utc1xv7SOrH5p7fmerl27RgBoz549ss4XVZ9NmzYRAPr4448pPz+fbt26Rf/v//0/0mg09Prrr9e4bXZ2NgUEBFBoaCht376diouL6fjx4zRu3Djq2LEjeXp62nV7R2nZsiX961//csqy5fjiiy8oJCTEoct0x7pVnzlcSSfPHT58OM2fP9+hy3TH2j3//PM0atQopyxbDu6X1rO1X1o9tLKysggA/fbbb7KC1WfTpk11/mebMWMG6fV6KioqslwWHR1NAGjNmjU1bnv16lXy9PSs9aTaentH6dSpE33yySdOWbYcixYtonvuucehy3THuqWmphIAOnfunFOWL0dkZCS99dZbDl2mO9buzTffpPvuu88py5aD+6X1bO2XNr09CNw9wagjPfTQQ3V+MTAyMhJGoxGnTp2yXPbzzz8DQK2ThrZt27bO945tvb2j+Pj4WJ4vJSgtLYWPj49Dl+mudQOguNrxNtc4Hx8flJSUOG35tuJ+aT1b+6XVQ6v6i2Xl5eW2p2pAUVER3n33XfTu3RstWrSwfJHtjTfeAHB3owWAiooKFBcXw8vLq87PZlq3bl3j37be3pHKysrg7e3ttOXbytPTk+tmhbKyMgBQVO28vLwc/hmbu9ZOaXUDuF9aw9baWT20/Pz8AMDhezMPP/wwPvzwQzz77LPIysqC2WwGEWHJkiUA7h6FA9xtvH5+figvL68zQ35+fo1/23p7RyouLrY8X0rg5+eH4uJihy7TXesGQFG18/X15dpZQYnbHMD90hq21s7moeXIDchkMmH//v0ICQnB/Pnz0apVK8v51ar3en+v+sST1S9jq92+fbvOQyZtvb2jKHEDcuTG4851A5Q1tLh21lHiNgdwv7SGzbWz9sMvk8lEer2e4uPjbfqQrTFRUVEEgP72t7/RrVu3qLS0lJKTk6l9+/YEgHbs2GG57dmzZykwMLDG0S2nTp2i8ePHU+vWrWt9UGjr7R3BYDAQANqwYYPDly3X2rVrSaPRUHl5ucOW6W51IyJauXIleXp6ktlsdsry5Xj88cdpypQpDl2mO9buoYceoieffNIpy5aD+6V15PRLq4cWEVF4eDi9//77NgdryK1bt2jevHkUFhZGer2egoODafbs2bRgwQICQABqHBWUmZlJkydPJn9/f/L29qYBAwbQ5s2bacyYMZbbP/PMM7Jvb6+0tDQCQKdOnXLYMu2Vnp7u8EzuVjciorfffpt69uzp0GXayxmZ3LF24eHh9N577zl0mfbiftk4Of1SIvq/N0GtMGnSJPj5+Sn+561FWr16NaZPnw6DwaCYU++UlZXB19cXa9aswZQpU0THUazo6GiYzWasXbtWdBSL2NhYPPvsszAYDNBqtaLjKJLRaISPjw9iYmLwxBNPiI5jwf2ycXL6pU3nHoyIiEBGRoascE1FVlaW4s4V51G6rhEAACAASURBVO3tjfbt2yvqNDdKpMRTJinxXJZKc/78eRiNRkXWjvtlw+T0S5uGVmRkJE6dOmU5rJLVdvjwYUX+wmxkZCSOHDkiOoZilZSU4MyZM4qrXY8ePaDX63H48GHRURTr8OHD8PT0RPfu3UVHqYH7ZePk9EubhlZUVBQqKirw66+/2rSSpsJkMmHv3r2K/KXSUaNGYffu3XWen4zB8iOCSvqxReDuIe/9+/dX5K9hK0VycjIGDx6sqO9pAdwvGyO3X9o0tKq/Ga3En7ZWgtTUVBQWFjrtB9PsMXr0aOTn5yM9PV10FEVKSUlBz549ERISIjpKLVFRUYr8NWylSElJQVRUlOgYtXC/bJjcfmnT0ALubkDbtm2z9W5Nwvbt29G2bVtF/thinz590KpVK2zfvl10FEXavn27IhsfcHeby87ORk5OjugoipOZmYkLFy4ounbcL+smt1/aPLSmTp2Ko0eP8geMdUhISMDUqVNFx6iTJEl45JFHEBsbKzqK4pw6dQrHjx/Ho48+KjpKne6//34EBwcjISFBdBTFiY+PR9u2bTFkyBDRUerE/bJ+cvulrFdaYWFhiIuLs3ll7uzIkSPIyMjAzJkzRUep18yZM3Hq1CmkpaWJjqIoK1euRIcOHTB8+HDRUeqk0+nw+OOPIyYmBjZ8Q8XtERHi4+Mxffp0xX4dgPtl3ezplzYPLY1Gg+nTpyM2NpY/1P+d2NhYdOvWDQMGDBAdpV5Dhw7FPffcg5iYGNFRFMNkMiEhIQEzZ860nBJHiZ566ilkZ2fj0KFDoqMoxr59+5CTk6PoHUXul3Wzq1/K+RZzZmYmaTQaWr16tZy7u528vDzy8/Ojv//976KjNGrx4sXUvHlzKiwsFB1FERISEkir1dLZs2dFR2lUv379KDo6WnQMxZgyZQoNGDBAdIxGcb+syd5+KWtoEd39wbDIyEhFnadNlPfee48CAwPpzp07oqM0qqioiAICAmjRokWiowhnNpupd+/eijpnXUNWr15NGo1GUacIE+XUqVOk0Who3bp1oqNYhfvlf9jbL2UPrbS0NJIkiTZv3ix3EW6hqKiIWrRoQR988IHoKFZ7++23KSgoiIqLi0VHEWrt2rUkSRKdOHFCdBSrmEwm6tmzJ82ePVt0FOGmTZtGPXr0IJPJJDqKVbhf3uWIfil7aBERTZo0iXr37k1Go9Gexajam2++SQEBAZSfny86itVu3bpFfn5+9O6774qOIkxFRQV1796dHn30UdFRbBIbG0s6nY7S0tJERxHmyJEjpNFoKCkpSXQUm3C/dEy/tGtoZWdnk5eXF3322Wf2LEa1MjMzydPTk7766ivRUWz22WefkaenJ2VkZIiOIsTixYvJ29ubcnJyREexidlspqFDh9LAgQNV8yrDkUwmEw0aNIhGjBihurfauF86pl/aNbSIiN59913y8/OjK1eu2Lso1YmKiqJ7772XqqqqREexmdFopMjISHrggQdER3G5S5cuka+vr2o/1zt+/DjpdDpavny56Cgu969//Yt0Oh2lp6eLjiIL90v7+6XdQ6u0tJS6dOlC48ePb1J7fl999RVptVo6dOiQ6Ciy7d27lzQaDS1dulR0FJepqqqiqKgo6tatm0N/FNPVXn75ZWrRogWdP39edBSXyc7OJn9/f/rrX/8qOops3C/t75d2Dy0iosOHD5OHh4dq91xtlZ6eTt7e3qo6+KI+b731Fnl5edGxY8dER3GJ9957jzw9Peno0aOio9ilrKyM+vbtSwMGDKCKigrRcZyuvLyc7r33XurTpw+VlpaKjmMX7pf2ccjQIiJasmQJ6XQ6+uWXXxy1SEUqLCyke+65h6KiotxiT8loNNL9999PERERqjhk3x47d+4krVZL//znP0VHcYgzZ86Qr68vvfzyy6KjON3zzz9P/v7+lJ2dLTqKQ3C/lM9hQ8tsNtOUKVMoKCiITp8+7ajFKkp5eTmNHj2a2rZtS7m5uaLjOMyVK1coODiYHnjgAbfdaz9x4gS1aNGCHn/8cdFRHCouLo4kSXKbQVyXzz//nCRJcqsv53K/lM9hQ4vo7vu1w4cPp9DQULpw4YIjFy2cyWSi6Oho8vf3d8u30tLT0ykgIICeeOIJt3gF+XuXL1+m9u3b05AhQ8hgMIiO43AfffSRKg8Bt8aqVatIo9Go4mwztuJ+KY9DhxbR3VN09OzZk3r06EHXrl1z9OKFMJlM9Oyzz5KXlxft3r1bdByn2blzJ3l6etK8efPcZnBduXKFunXrRn369KGCggLRcZxm/vz55Onp6VZfXt24cSN5eHjQa6+9JjqK03C/tJ3DhxbR3T3biIgI6tSpE2VmZjpjFS5TXl5Ojz76KHl5edHGjRtFx3G6devWkaenJz322GOqPrqO6O5nPh06dHCrhlAfk8lETz/9NOn1elq5cqXoOHZbvnw56XQ6mjdvnuq+j2Ur7pe2ccrQIrp71oWBAwdSq1atVHtYeEFBAY0aNYoCAgLc/gPT30tOTiZ/f38aM2aMak+s++uvv1JQUBANHTqU8vLyRMdxCbPZTAsWLCBJkujTTz8VHUcWs9lMH330EUmSRO+8847oOC7D/dJ6ThtaREQlJSU0YcIEVZ414vDhw9SpUycKDQ1V7RcZ7XHs2DFq06YNdenShVJTU0XHsZrZbKbPP/+cPDw86OGHH3bLz7Aas2TJEtJqtRQdHa2qnY6CggKaMmUKabVa1fULR+B+aR2nDi2iu1/mfP/990mr1dKjjz6q+I3IbDbTkiVLyMPDg8aNG0c3btwQHUmY3NxcGjNmDHl6etKXX36p+Ldp8vLyaPLkyaTT6eijjz5ym8/l5Ni1axeFhISoZqfj4MGD1LFjR2rbtq1bf27cGO6XjXP60KqWkpJCbdu2pZCQEFq5cqUiG2B6ejoNGzaMtFotvffee0266VUzm830ySefkFarpREjRij2jOirV6+m4OBgCg4Oph07doiOowg3b96k8ePHk06no/nz51NRUZHoSLUUFhbS/PnzSafTUVRUlFt9lcQe3C/r57KhRXR3T3jevHmk0Who9OjRinnbLS8vj1555RXS6XQ0ePBgtzyk3V6pqak0YMAA0uv19PrrryvmrPbHjh2jESNGkFarpeeff96tjxCUw2Qy0dKlSykoKIhCQ0MpPj5eETtjJpOJVq5cSSEhIdSqVStavny5IhuzSNwv6+bSoVXt4MGDdN9995EkSTR58mQ6cuSIiBh048YN+utf/0p+fn7UsmVLWrp0qSI2aKUymUz09ddfU2BgIPn7+9PChQvp5s2bQrIcPHiQHn74YZIkiQYOHKiKt8BEunXrFj3zzDOk0Wioe/fuFBMTI+QnMiorK2nFihUUERFBWq2W5s2b12QOlJGL+2VNQoYW0d23ndavX08DBgwgABQVFUUxMTFUUlLi1PWaTCbavXs3Pf3009SsWTMKDg6mv/3tb03+BxFtUVRURIsXL6bWrVtTs2bNaM6cOfTLL784fU/5zp07tHLlSho9ejQBoMGDB9OmTZt4D90GZ86coaeeeop0Oh116tSJPvzwQ5ecdPfcuXP0wQcfUMeOHUmv19Of/vQn1R/e7UrcL/9D2ND6vW3bttEf//hH8vDwIF9fX5o5cyatWrWKrl+/7pDll5aW0s6dO+nNN9+k9u3bEwC699576Z///KfqT74pksFgoC+//JL69u1LAKhjx460cOFC2rVrF5WVlTlkHbm5uRQfH09PPvkkNWvWjDw9PWnKlCm0c+dOhyy/qcrJyaGXXnqJWrduTZIk0f33309ffPEFnThxwiE7AWazmY4fP05Lliyh4cOHkyRJFBwcTK+88orbnf3B1Zp6v5SIiKAQeXl5SExMRFJSEg4cOACTyYSePXti2LBh6NatGyIiIhAREYHg4GD4+PgAAOLi4jBjxgwAQGVlJfLz85GdnY3MzExkZWUhNTUVBw4cQHl5OcLDwzFlyhTMnDkTPXv2FPlQ3c7JkycRExODdevW4ezZs/Dy8sLQoUPRv39/hIeHIyIiAvfccw+CgoKg1+sBAAkJCZg+fToAwGAw4Pr168jKykJGRgYyMzOxb98+nDp1Cnq9HkOGDMETTzyBxx9/HIGBgSIfqlupqqrCzz//jPj4eGzbtg0FBQUIDg7G/fffjz59+iA8PBzh4eFo37695XlPT0+HVqtFr169AAD5+fm4dOmSZZs7fvw49uzZg5s3byIwMBDjx4/HjBkzMG7cOOh0OpEP16001X6pqKH1eyUlJdizZw+Sk5ORmpqKjIwM3Lhxw3K9JEkICAhAcXExfHx8UFZWhsrKSsv1vr6+CA8PR+/evREVFYWoqCi0a9dOxENpci5duoTk5GSkpKTgxIkTyMrKgsFgsFzv6ekJT09PlJaWws/PD4WFhfj9f8OQkBB069YNAwYMQFRUFEaMGGHZ6JjzmEwm/Pbbb0hOTsa+fftw5swZXLhwAVVVVZbb+Pr6WmolSRJKSkos1+l0OnTq1Andu3fHiBEjEBUVhb59+0Kj0bj8sTQ1TalfKnZo1aWwsBBZWVm4desWSkpKkJ6ejsWLF2PcuHGYOXMmfH190aJFC3Tu3BlhYWGi47LfuXz5MnJyclBQUICSkhLExMRgx44deOuttxAZGQlfX1+0bNkSERERaN68uei47P9UVlbi3LlzuHLlCgoLC3Hnzh289tprkCQJ//jHP+Dv74+AgACEhYWhc+fO8PDwEB2Z/R937Zeqeq0eEBCAgQMHWv598OBBAMDx48cxffp03qNTsLCwMMuGYTab8frrrwMAysvL8dhjj4mMxhrg4eGB7t27o3v37gCA/fv3o6ioCADQo0cPDBkyRGQ81gB37ZfqTI27jS8uLg4AcP36dezfv19wImatPXv2WN66iImJgclkEpyIWSshIQF6vR4eHh5YtWqV6DjMSu7UL1U7tJKTk3H79m0AgF6v5w1IRVatWmV5G+n27dvYs2eP4ETMGlVVVUhMTITRaERlZSXi4uJqfN7FlMud+qVqh1ZCQoKl8RmNRiQkJMBoNApOxRpjNBqRmJho+RBYr9cjISFBcCpmjR07diA/P9/y74KCAiQnJwtMxKzlTv1SlUOrsrISa9asqXH0S1FREXbs2CEwFbPGzz//jDt37lj+bTQakZSUhIqKCoGpmDWq3xqsptfrER8fLzARs4a79UtVDq0tW7bUONQW4A1ILeLj42s0PuDu4brbtm0TlIhZo7y8HOvWrauxd240GrFmzRqUlZUJTMYa4279UpVDKz4+vtaXFI1GI3788cca3wdiylJaWooNGzbUeltCq9WqdgNqKjZt2oTS0tJal5eVlWHr1q0CEjFruVu/VN3QKi4uxubNm+t8P7aiogJbtmwRkIpZY8OGDXW+DVhVVYUNGzbU2htkyhEXFwetVlvrcq1WazkqjSmPO/ZL1Q2t/36L4vd4A1K2uLi4er8bYjQasWHDBhcnYta4c+cOtm7dWueRglVVVdi8ebPlu1tMWdyxX6puaMXFxUGSpDqvq6qqwtatW2sc4cSUoaCgANu3b6/3O1mSJCE2NtbFqZg11qxZA7PZXO/1JpMJ69evd2EiZi137JeqGlq3bt1CcnJyg19GJSKsW7fOhamYNX744YdGG9/OnTst3yVhytHYzgTvcCiTu/ZLVQ2tNWvWNHr2BLPZzBuQAsXFxaGx01yaTCasXbvWRYmYNa5fv449e/Y0uN2ZTCbs3r0bN2/edGEy1hh37ZeqGlrWfIubiLBv3z7k5ua6IBGzxrVr17B///5GhxZgXY2Z66xevbrBV8jVTCYTVq9e7YJEzFru2i9VdZb3/0ZE0Gg0WL16NaKjo0XHYTZISkrCtGnTrGqITFmio6Oh0WiQlJQkOgqzgbv0S1W90mKMMda08dBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6rBQ4sxxphq8NBijDGmGjy0GGOMqQYPLcYYY6ohERHZcod+/fo5K4sseXl58PPzg4eHh+goFn//+98xduxY0TFq2L59O958803RMSwqKipQUlKCoKAg0VEsJEnCsWPHRMeo5Y033sDOnTtFx7AoLi4GAPj5+QlO8h/jxo3Dp59+KjpGLdwvG2drv9TZsnAiQlpaGiZPnoxu3brZHK4p+OSTT1BQUCA6Ri0FBQVIS0vDggULREdRpDNnzmDDhg2iY9Tp/PnzKC8vx+TJk0VHUaR169bh/PnzomPUwv2ycXL6pU1Dq9r06dMRHR0t565ujYjwySefiI7RoMWLF4uOoEhJSUmKHVoA0Lt3b65dPbKzs0VHaBD3y7rJ7Zf8mRZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHVcOnQunjxIiZNmoQ7d+7UuDwtLQ0TJ05EQEAA/Pz8MHbsWOzfv99h67V2+QsWLEBSUpLD1utORNUOAH766SeEh4dDp6v/Rwm4dvWrr3aAdc+tXLzd2UdU3axdvqi6uWxopaWloX///hg3bhz8/f0tlx86dAhDhw6Fn58fzpw5g/Pnz6Nz584YNWoUtm/fbvd6bVn+s88+i4ULF+Kdd96xe73uRFTtzp07h0mTJmHhwoW4ceNGg7fl2tWtvtrZ8tzKwdudfUTVTRXbHNnAbDYTAFq9erUtd6OioiJq164dzZs3r8blJpOJevbsSW3atKHS0lLL5VVVVRQREUFhYWFUXl5u07rsXX5aWhpJkkRJSUk2r0/u8+MKiYmJZGO5iUhc7YiIpk2bRosXLyaj0UihoaGk1WobvL09tZP7/LjC1KlTKTo62ub71Vc7ItufW1u4eruT+/w4m6P7JZFz6yZn+SL6pUuG1ttvv006nY6uXr1a4/KUlBQCQC+++GKt+7z//vsEgNasWWPTuhyx/OjoaGrXrh0ZjUab1ueOQ0tU7YioRsOzdgOVWzt3HFr11Y5I3nNrLVdvd+42tETVTe7yXd0vnf72IBFh2bJlGDRoENq2bVvjuuTkZABA//79a92v+rJdu3bJXrfc5U+ZMgVXrlzBli1bZK/bHYisHQB4e3vbfB+u3V0N1Q6Q99xai7c7+UTWTe7yXV03pw+t9PR03LhxA5GRkbWuy8jIAAC0a9eu1nWhoaEAgKysLNnrlrv8vn37AgC2bdsme93uQGTt5OLa3dVQ7ZyNtzv5RNZNLlfXzelD6+TJkwDq/g9cWFgIAPDx8al1na+vLwCgoKBA9rrlLr96w6rO3lSJrJ1cXLu7Gqqds/F2J5/Iusnl6ro5fWjl5uYCAJo3b27T/YgIACBJksMzNbZ8f39/SJJkyd5UKbV2DeHa3SW3ds7G213DlFq3hri6bk4fWuXl5QAAvV5f67qAgAAAgMFgqHVd9WXVt5HDnuXrdDqUlZXJXrc7EFk7e3DtGq6ds/F2J5/IutnDlXVz+tDy8vICABiNxlrXdevWDQBw5cqVWtddvXoVABAeHi573fYsv6qqyukfeiqdyNrZg2vXcO2cjbc7+UTWzR6urJvTh1abNm0AAEVFRbWuGz16NADg6NGjta6rvmzMmDGy1y13+Xfu3AERWbI3VSJrJxfX7q6GaudsvN3JJ7Jucrm6bk4fWr169QJQ917XyJEj0aNHD6xZs8byshgATCYTEhMTERYWhokTJ8pet9zlV+8NVmdvqkTWTi6u3V0N1c7ZeLuTT2Td5HJ13Zw+tCIjI9G6dWukp6fXXrlGg+XLlyM/Px9/+tOfcP36deTl5eGFF15AdnY2vv32W8vL5WozZsyAJEk4f/58o+uWs3zg7ilUAGDcuHEyH7V7EFk7ubh2dzVUOzl4u3MNkXWTy9V1c/rQkiQJc+bMwaFDh3Dt2rVa1w8ePBi//vorioqKEBERgY4dOyI7Oxu7d+/G+PHja90+NzcXvr6+aN++vVXrt3X5ALBu3TqEhoYKeaWgJKJrt3nzZkiSBEmScPXqVZhMJsu/ly1bVud9uHZ3NVY7W59b3u5cQ3TdVLHNueK0G4WFhRQaGlrnubRsUVBQQN7e3jRnzhy7ltOQ6nNprVq1yub7uuNpnJpK7dzxNE5NpXbudhqnplI3RZ97kIjo2LFjFBQURF999ZXN961e98yZMyk4OJhyc3NlLaMx586do86dO9Nbb70l6/7uOLSImkbt3HFoETWN2rnb0CJqGnVT7LkHq/Xr1w+pqanYunVrnb8P05gbN24gJycHu3btQkhIiBMSAt988w0WLVqERYsWOWX5asW1Uy+unTpx3ernnF8Qq0fHjh2xefNmWfcNCQnBvn37HJyopk8//dSpy1czrp16ce3UietWN5f+cjFjjDFmDx5ajDHGVIOHFmOMMdXgocUYY0w1eGgxxhhTDR5ajDHGVIOHFmOMMdXgocUYY0w1eGgxxhhTDR5ajDHGVIOHFmOMMdXgocUYY0w1eGgxxhhTDdUPrWPHjomOwGTi2qnTtWvXcP36ddExmAzusM3J+mmSv/71r/j4448dnUWW7OxsdOzYEXq9XnQUVejXr5/oCAAAo9GICxcuoGvXrqKjAAAKCwtFR2jQjh07FFO7W7duAQBatWolOMldFy5cwAMPPCA6Rr24XzqWTUNLkiQsWLDAWVlsVlhYiLS0NAQHB2PgwIGi4wAAHnzwQXTv3l10jFp69OihqNodOnQIp06dwpAhQ9C8eXPRcQDc/f+tRI888ohihjsAfPfdd9BoNHjwwQdFR7Ho06eP6Ai1cL9snKx+Ket3khVi0aJFBID69OkjOgqzUc+ePQkAffLJJ6KjMBucPXuWJEkiAJSVlSU6DrOBu/RLVX+mFRsbCwA4fvw4MjMzBadh1jpz5gxOnToFAFi5cqXgNMwW8fHx0Ol00Ov1SEpKEh2H2cBd+qVqh9aJEyeQkZEBAPDw8MCqVasEJ2LWWrVqleU99d8PMKZ8sbGxMBqNMBqN+P7770XHYVZyp36p2qH1+8ZXWVmJFStWiA3ErBYTEwOj0Qjg7gaUmJgoOBGzxrFjx3D27FnLv3NycpCeni4wEbOWO/VLVQ4tIrLs8VW7ePEijh49KjAVs8bhw4dx8eJFy7+rNyAiEpiKWWPVqlXw8PCw/Fuv16t6j72pcLd+qcqhdeDAAVy5cqXGZWp/ydtU/H6Pr9qVK1dw+PBhQYmYNYgI8fHxqKystFxmNBoRExPDOxwK5279UpVD67/3+IC7e+wxMTEwmUyCUrHGmM1mxMfH19jjA9S9ATUVe/bsQW5ubq3Lc3NzsX//fgGJmLXcrV+qbmiZTCYkJCTU2OOrduvWLezdu1dAKmaNlJQUyxdTf6+yshJxcXGq3ICairoaH8BvESqdO/ZL1Q2tnTt3Ij8/v87r9Ho9EhISXJyIWSshIaHeb+Ln5eUhJSXFxYmYNYxGIxITE+tsfEajEQkJCbVePTNlcMd+qbqh1VDja2jjYmJVVlbihx9+qLe5qXUDagq2b9+OoqKieq8vLCzErl27XJiIWcsd+6WqhlZ5eTnWrl3b4F5dSUkJtm/f7sJUzBpbt25FSUlJvdcbjUYkJSWhvLzchamYNRpqfADvcCiVu/ZLVQ2tLVu2oLS0tMHbaLVaxMfHuygRs1ZCQgK0Wm2DtyktLcXWrVtdlIhZo7S0FOvWrWuw8RmNRqxZswZlZWUuTMYa4679UlVDa9WqVY0eXltVVYWNGzfCYDC4KBVrTElJCTZt2oSqqqpGb8sf6ivLpk2brBpGZWVl2Lx5swsSMWu5a7+USEVfsti7dy8qKios/yYijBs3Du+88w7uv//+GrcdMGCAYs4e3tQVFRXhyJEjNS7bvXs3Pv7441pvTXh5eWH48OGujMcakJWVhUuXLtW47MMPP4QkSfif//mfGpe3b98e4eHhrozHGuCu/VJVQ+u/ERE0Gg1Wr16N6Oho0XGYDZKSkjBt2jSYzWbRUZiNoqOjodFo+IS5KuMu/VJVbw8yxhhr2nhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlSDhxZjjDHV4KHFGGNMNXhoMcYYUw0eWowxxlRDIiISHaIhFRUVyMrKQmZmJrKysnDp0iUUFRXBYDDAYDDg8uXLCAwMRIsWLeDn54fmzZujZcuWCA8PR/fu3REeHo7AwEDRD6NJysvLQ1ZWFjIyMpCZmYm8vDwUFRWhuLgYBQUFyM/PR1hYGHx9feHj4wN/f3906NAB4eHhiIiIQNeuXeHp6Sn6YTQ5ZrMZly5dsmx32dnZlm2uqKgIN27cAAAEBwejefPm8PHxQfPmzREeHm6pXfv27SFJkuBH0vQ0hX6puKGVl5eH3bt3IyUlBcnJycjMzITZbIZWq0WHDh3QoUMHBAQEwMfHBz4+PggICEBpaSlKSkpQUlKCwsJC3Lx5E5mZmSgvLwcAhISE4P7778fo0aMRFRWF8PBwwY/SPWVkZFjqtnfvXktz8/b2Rnh4uKXJ+fn5wcfHB82aNUNhYaFlgyosLMSFCxdwOW25ngAAIABJREFU8eJFS80jIiIQFRWF0aNHY9SoUYrfoNSosrIShw4dQnJyMlJSUnD48GGUlZUBAFq2bImIiAjLNhcQEABfX18AsGxvBoMBBQUFyMrKwu3btwEAzZo1w8CBAy3b3KBBg6DX64U9RnfVFPulIobWzZs3sWrVKsTFxeHYsWOQJAn33nsvoqKiMGDAAFl73b/fWzx+/DhSUlKwZ88elJSUoF27dnjsscfw1FNPITIy0omPzP399ttviImJwQ8//ICrV6/Cz88PI0eOxKhRo9CnTx+Eh4fbvNddvbeYlZWFQ4cOISUlBb/99huICP3798eTTz6JadOmoVWrVk58ZO6tsrISP/30E2JiYrBt2zaUlpaiQ4cOiIqKwvDhw9GjRw9Ze935+fnIzMzE6dOnsXfvXqSkpODSpUvw8fHBgw8+iKeeegoTJkzgAWaHJt8vSRCTyUQ//vgjPfTQQ6TT6cjf35+efvpp2rBhAxUWFjplnUajkfbv30/vvvsudenShQBQnz596LPPPnPaOt1Rfn4+/f3vf6devXoRAOratSu9//77dODAATIajU5ZZ0FBAa1fv55mz55Nfn5+pNfr6eGHH6b169eTyWRyyjrd0fHjx+mFF16goKAg0mg0NGbMGFq6dCmdO3fOaevMzs6mb775hkaPHk0ajYZatmxJL774Ip08edJp63Q33C//w+VDq7KyklasWEHdunUjjUZDDz74IMXHx1NpaalLc5jNZtq3bx/NnTuXmjdvTs2bN6e33nqLbt686dIcanL9+nV68803yc/PjwICAui5556j/fv3uzyHwWCguLg4Gj9+PGk0GurRowfFxMQ4bWC6g4MHD9KkSZNIkiTq1q0bLV68mC5duuTyHBcvXqRFixZRREQEaTQamjJlCh05csTlOdSC+2VtLhtaZrOZli9fTh07diS9Xk+zZ8+mM2fOuGr1DSoqKqKPP/6YWrduTc2aNaNXX32VCgoKRMdSjLy8PJo/fz55e3tTcHAwffrpp3Tnzh3RsYiI6PTp0/TUU0+RTqejzp0708qVK8lsNouOpRipqak0duxYAkCDBg2iDRs2KOL5MZvNtG7dOhowYAABoPHjx1NaWproWIrB/bJ+LhlaaWlpNGTIENLpdPTnP/+ZLly44IrV2sxgMNAXX3xBrVu3ppCQEIqNjVXEBi6K2Wym77//nlq1akUhISH05ZdfunwPz1o5OTk0d+5c0mq1NGLECDp+/LjoSEIVFBTQX/7yF9JqtTRs2DDauXOn6Ej12rZtm6U/vPzyy1RUVCQ6klDcLxvm1KFlMBjo5ZdfJp1OR0OHDlXNnlR+fj79+c9/Jo1GQ6NGjaLMzEzRkVzu9OnTNGLECNJqtfTiiy+q5jO/Y8eO0aBBg0iv19Nrr72m2CHrTAkJCRQSEkKtW7emFStWqGLHq/qVRcuWLalt27b0ww8/iI7kctwvreO0oXXy5Enq0aMHBQYG0rJly1Sx4fy3I0eO0H333Ue+vr4UGxsrOo7LfPfdd9SsWTMaMGAAHT16VHQcm5lMJlq6dCm1aNGCevfurZi3VZytpKSEZs2aRZIk0XPPPUf5+fmiI9ksLy+Pnn32WZIkiebMmdNkdjq4X1rPKUNr5cqV5OPjQwMGDKCcnBxnrMJljEYjvffee6TRaGjmzJlUUlIiOpLTlJaWWhrG/PnzqbKyUnQku1y6dImGDRtG3t7etHTpUtFxnOr06dPUu3dvCgoKos2bN4uOY7eNGzdSYGAgde/e3e3f6uV+aRuHDq2qqiqaO3cuSZJEf/3rX93qaK7qjahPnz505coV0XEc7tKlS9SzZ09q2bIlbdmyRXQch6msrKTXXnuNJEmi559/nqqqqkRHcrgff/yRvL29adiwYUKOCHSWCxcu0ODBg8nHx4c2btwoOo7Dcb+Ux2FDq6ysjCZPnkze3t60fv16Ry1WUS5cuEDdu3enDh06UEZGhug4DnP69GkKCwujXr16uVXT+701a9aQl5cXPfroo1ReXi46jsN8++23pNVq6bnnnnOrpletsrKS5syZQzqdjr777jvRcRyG+6V8DhlaxcXFNHbsWAoICKA9e/Y4YpGKlZ+fT8OGDaPAwED69ddfRcex2+HDh6lVq1Y0aNAgun37tug4TnXgwAEKCgqiUaNGqebAkoZ88sknJEkSvfnmm6KjOF31Y33vvfdER7Eb90v72D20iouLaeDAgdS2bVs6ceKEIzIpnsFgoIkTJ5Kfnx+lpqaKjiPb4cOHydfXlx5++OEm84F3eno6tWnThoYMGUIGg0F0HNk++OAD0mg09K9//Ut0FJf58ssvSaPR0EcffSQ6imzcL+3vl3YNrcrKSnrwwQepZcuWbvV2mTUqKytpwoQJqn3s2dnZFBwcTBMnTlT9ARe2qn7sY8eOpYqKCtFxbPbvf/+bJEly+4NL6qLmx8790jH9UvbQMpvNNGvWLNW/2rCHwWCgoUOHUlhYmKo+C7p69Sp17NiRBg0a5NZHQzbkyJEj5OvrS9OnT1fVuQvXr19PWq2WFi1aJDqKMO+99x5ptVpVfZeL+6Xj+qXsofXOO++Qp6enor9p7wq3bt2iiIgIuu+++1TxAX9ZWRn17duXevbsSXl5eaLjCPXzzz+TXq+nDz74QHQUqxw5coQ8PT3pxRdfFB1FuOeee468vLzo2LFjoqNYhfvlXY7ol7KG1vbt20mj0dA333wja6XuJisri/z9/ekvf/mL6CiNeu655yggIED13wdxlC+//JK0Wi0lJyeLjtKgwsJC6ty5M40bN05VrwydpaqqikaPHk1du3ZVzHkw68P9siZ7+6XNQ+v69evUpk0beuyxx2St0F0lJSWRJEm0du1a0VHqtXr1asVnFGHatGkUHBxMubm5oqPUSw0ZXU0NvUgNGUWwp1/aNLTMZjONGTOGwsPDFb93I0L1q5iLFy+KjlJLTk4O+fv781tLdSgsLKQuXbrQ+PHjFXn6nK+//pq0Wi2lpKSIjqI4O3bsII1GQ8uWLRMdpRbulw2T2y9tGlrff/89abXaJvtBYmPKysqoW7duNHnyZNFRavnDH/5AvXr1UsXnbiIcPHiQNBoNxcfHi45Sw7Vr1yy/XcTq9sYbb1CLFi3oxo0boqPUwP2yYXL7pdVDKz8/n1q3bk3z58+3OVxTsnv3bpIkiTZt2iQ6isXatWtJkiTeU2/Ec889R8HBwYr6LbVp06ZR+/btm+xRntYwGAzUsWNHmj17tugoFtwvrSOnX1o9tObOnau4DVqpnnjiCerQoYMivrxavUHPmjVLdBTFq240L730kugoRPSfDdodz7vnaErbMeN+aT1b+6VVQ+vEiROk0WgoLi7OrnBNxdWrV8nPz48WL14sOgq9//77FBAQoLi3TpTqu+++I61WK/zLn2bz/2/vzKOiOu83/r3DLCwzMiCbqMiAAopLlFJBo8igASLGAwnHgmLUNiatqU2TpuSkrUnVNlpbq1abejwnKijpMXWpqUkUYXABFHfrNoOIEmVRUcYZUGCY5/dHDvRHZOcuw8z9nMMfITPv81y/vM/3ztz3vteKcePGYe7cuYL6GEgkJCQgIiJC8OuSYl72jt7mZY+aVlpaGsaMGSMute0FH374Iby8vGAymQTzYDab4eXlhY8//lgwDwONlpYWhIeHY/HixYL62L9/PxiGcZitftjgwoULYBgGX331laA+xLzsPb3Jy26bVllZGaRSKXJyclgx5yg8fPgQSqUSGzZsEMzD2rVr4ebmhgcPHgjmYSCyY8cOyGQylJeXC+YhOjraJhf02Dovv/wyoqKiBNMX87Jv9CYvu21ab7zxBoKDg+3ysQdc88tf/hLDhg0TZMXes2fP4O/vj1//+te8aw90mpqaoNFo8LOf/UwQ/SNHjoCIUFJSIoj+QKa4uBhEhIKCAkH0xbzsOz3Nyy6b1qNHj+Ds7Czeyd1H7t27B5lMJshZV1ZWFuRyuXgzah/ZsmULXF1dYTQaedeeM2cO4uLieNe1F2JiYpCcnMy7rpiX/aOneSmhLtizZw8xDEPz5s3r6mUineDv70/x8fGUnZ3Nu/bOnTspKSmJ/Pz8eNe2B+bPn08A6F//+hevuvfv36dvvvmGfvzjH/Oqa08sWbKEDh06RA8fPuRVV8zL/tHTvOyyaWVlZVFKSgq5u7uzas6RyMjIoCNHjlBVVRVvmvfu3aOCggLKyMjgTdPecHd3p6SkJN5POHJycsjFxYXmzp3Lq6498eqrr5JcLqc9e/bwqivmZf/pSV522rRu3rxJxcXFYvD1k1deeYVUKhXl5OTwppmVlUXu7u6UmJjIm6Y9kpGRQceOHaPy8nLeNLOzsyk1NZVcXV1507Q33NzcKDk5mdcTDjEv2aEnedlp09q3bx95e3vTzJkzOTHnKDg7O1NKSgrt27ePN839+/dTamoqKRQK3jTtkYSEBPLw8KADBw7wonf79m06f/48paWl8aJnz6Snp9OpU6fo3r17vOiJeckOPcnLTpuWTqcjrVZLTk5OnJhzJGbNmkUlJSX05MkTzrWMRiOdP3+eZs2axbmWvSOTySg2NpZ0Oh0vevn5+eTs7ExTp07lRc+emT59OikUCsrPz+dFT8xL9uguLztsWs3NzXTy5EmKjY3l1JyjoNVqqaWlhQoLCznXKigoIAAUExPDuZYjEBsbSwUFBWSxWDjX0ul09OKLL5KzszPnWvaOq6srRUVF8XLCIeYlu3SXlx02rZKSEjKbzaTVajk15yj4+PhQeHg4LxNIp9PRhAkTyMvLi3MtR0Cr1ZLJZKLz589zrqXT6cTgY5HY2FjKy8vjXEfMS3bpLi87bFqFhYU0dOhQGjlyJKfmHIkZM2bQiRMnONcpLCwUP2WxyOjRo8nX15dOnjzJqU55eTndu3dPrB2LzJgxgyoqKuju3buc6oh5yT5d5WWHTev69es0duxYTk05GuHh4XTjxg1ONQCQXq8Xa8cyfNSudfzw8HBOdRyJcePGERFxXjsxL9mnqznXYdPS6/UUGhrKqSlHIyQkhOrq6uj+/fucaVRWVpLJZKKQkBDONByR0NBQ0uv1nGro9Xry8/MjtVrNqY4j4enpSYMHD+aldmJesktXedlh0zIYDGLwsUzrHzWXE6h1bHECsUtISAgvwSfOOfbho3ZiXrJPV3n5XNOqra2l2tpasQgs4+/vTyqVitMJZDAYSK1Wk4+PD2cajkhoaCjV1NSQ0WjkTMNgMIgnGxwQFhbG6ZwT85IbusrLDpsWEXEWfLW1tfTuu+9ScHAwyeVy8vDwoMTExHYrRQ4cOEAMw7T93L59m+bNm0dqtZoGDx5MSUlJVFZW9tzYDx48oOXLl1NgYCDJ5XLy9vamlJQUunjxYrvXNTY20ooVKygsLIxcXV3J09OT5syZQwcPHqSWlhZOjpthGPLy8mr79+WC2tpaThuWo9au9d90oNbOUetGROTt7c153YjEvGSbLvPy+zvonj17FkSEsrIy1nfxraqqgkajga+vL7788ksYjUbo9XqkpKSAYRhs27at3evnzp0LIsLcuXNRVFQEs9mM3NxcuLi4IDIyst1rKysrMWLECPj6+uLQoUMwmUy4cuUKYmJi4OzsjKKiorbX/uQnP4G7uzuOHDmChoYGVFdX41e/+hWIiNPHdY8fPx6/+c1vOBs/MzMTkyZN4mRsR66dXq8HEeHixYucjA8AQUFBnDzp2pHrBgCrVq1CSEgIZ+OLealj/bhb6Swvn2taOp0ORMTJ49kXLVoEIsLnn3/e7vetz35ycXFBdXV12+9bi/Dll1+2e/1rr70GImr3cMPXX38dRITdu3e3e21VVRUUCgUiIiLafqfRaDBlypTn/IWEhHBahKlTp2L58uWcjb9s2TLExMRwMrYj166yshJEhBMnTnAyPgD4+Pjgb3/7G+vjOnLdAGDDhg0YMmQIZ+OLeanrzyF2SWd5+dzXg2azmYiIVCoVCx/y2rN//34iIpo9e3a73ysUCoqLi6OnT5/S4cOHn3tfZGRku/8ePnw4EX23Wq6VAwcOkEQioaSkpHav9fPzo/DwcDp37lzb/RoJCQlUVFRES5cupVOnTrV9xNXr9TRjxoz+HWQXqFQqMplMnI1vMplIqVRyMrYj1651LnBdO3HOsQ/Xc07Myxn9O8gu6Kx2zzWt1u1qpFIpqwYaGxvJaDSSs7NzhwX29fUlIqLq6urn/t/3t/qXy+VERGS1WtuNbbVayd3dvd33uwzDtO1mUFpaSkREW7ZsoaysLLp16xbFxcXRoEGDKCEhoe2PhCtkMhk1NzdzNr7FYmG9bkRi7WQyGRHRgKudo9eNiJ85RyTmJRd0VrvnmpabmxsR/e8Mgi0UCgW5u7vTs2fPOuyeNTU1RER9emihQqEgtVpNUqmUmpubCd997fncT+sWOQzDUEZGBh09epTq6urowIEDBIBSUlJo/fr1/TvQLuDqbLoVpVJJ9fX1rI/r6LVrPWauayfOOfbhes6Jecl/7Z5rWq0vYrsIRETJyclERHTo0KF2v29sbKS8vDxycXGh+Pj4Po2dkpJCFoulw00W165dSwEBAW1nRWq1uu1ua5lMRrNmzWpbgfN9b2zCR9Pi6qsQR67dQG1aRI5dNyLu55yYl/zXrtOmxUX4ffLJJ6TRaOidd96h//znP2QymchgMFB6ejpVVVXRxo0b2z729mXs4OBgWrJkCX399ddkNBrp0aNHtHXrVlq5ciX9+c9/bvcR/q233qLLly9TY2Mj3b9/n/70pz8RAE43veRjAnHVtBy5dnw0La5q58h1I/qudoMGDeJsfDEvBcjL76/MqKioABGhsLCQlRUg3+fhw4d45513oNFoIJPJ4O7ujvj4eOTl5bW9pri4GETU7qd16eP3fz979uy299XW1uLdd99FUFAQZDIZvL298dJLLyE3N7edh4sXL+LNN9/E6NGj4erqCk9PT0RFRWHbtm2wWq2cHDcA+Pr6YuPGjZyNv379egwdOpSz8R21dseOHQMRobKykpPxASA6Ohq/+MUvOBnbUesGAG+//TamTZvG2fhiXvKfl881LYvFAmdnZ2RnZ3NmxhF58uQJGIZ5bjkqm+zbtw8SiQT19fWcaTgin332GVxcXNDS0sKZRlpaGubMmcPZ+I5KQkICFi5cyNn4Yl5yQ1d5+dzXg05OThQUFMT5fl2OhsFgIACcbtUTGhpKVquVbt68yZmGI9K6L6BE0umDvvsNH5vyOiJcb2Yr5iU3dJWXHc5CcQKxj16vJ7lcThqNhjONkSNHklQqFWvHMnzs4h0aGkq3bt3idHm2o9HY2EgVFRW81E6cc+zSVV522LTCwsI4fwaNo6HX6ykoKIiT+6hakcvlFBgYKNaOZfhoWmFhYWSxWDrcI06kb5SWllJLSwsvtRPnHLt0lZcdNq2IiAi6du0ap7taOxpFRUUUERHBuU5ERAQVFxdzruMoPHr0iPR6Pee1Gz16NLm4uFBRURGnOo5EUVERKZVKzpuWmJfs01Vedti0ZsyYQQDo+PHjnBpzFJqamqioqIjT5aGtxMbG0okTJ8SvmVhCp9MRwzA0ffp0TnUUCgVNnTq13e7dIv0jPz+fpk2b1rajCVeIecku3eVlh01r8ODBNG7cOHECsURRURE1NDTw0rS0Wi2ZzWY6c+YM51qOgE6no0mTJpGHhwfnWrGxsZSfn8+5jiMAgAoKCtp2deASMS/Zpbu87HQ5lFarpby8PM6MORI6nY40Gg0FBgZyrjVq1CgaPny4WDuWyM/P5yX4iL6bc5WVleL1ERa4cuUK1dTU8HKiSCTmJZt0l5edNq3ExES6fPkyGQwGrrw5DHv37qXExETe9BITE2nv3r286dkr165do+vXr/NWux/84Afk7e0t1o4F9u7dS35+fvTCCy/woifmJXt0l5edNq24uDgaNmwY7dq1ixNjjsK5c+fo6tWrlJGRwZvm/Pnz6dKlS3T58mXeNO2RrKwsCggI4Px6VitSqZR+9KMf0c6dOwkAL5r2CADatWsXzZ8/n5ycnHjRFPOSHXqSl502LYlEQmlpaZSVlSVOoH6QnZ1No0aNosmTJ/OmOW3aNAoKCqLs7GzeNO0Nq9VKOTk5lJGRwelNxd8nIyODSktLxWuS/aCwsJDKysp4PVEU85IdepKXXc7GhQsX0p07d8RVMX2kubmZPv/8c8rIyCCGYXjTZRiGFixYQLt3727bqVmkd+h0Ovr2229pwYIFvOpGRkbSmDFjaOfOnbzq2hPZ2dk0fvx4mjBhAq+6Yl72jx7nZXd7QEVFReGVV15hc1sph2H79u2QyWSoqKjgXbu8vBxSqVTcE62PvPzyy3jxxRcF0V6/fj2USiUePnwoiP5ApqamBq6urti0aZMg+mJe9p2e5mW3Tevf//43GIbB5cuXWTPnCFgsFoSGhmLJkiWCeViwYAFGjx7N6Uav9siFCxfAMAy+/vprQfTNZjO8vb2xYsUKQfQHMh988AF8fHwE2zRazMu+0Zu87LZpWa1WTJw4EWlpaayYcxRycnLg5OQEvV4vmIdr165BIpFg3759gnkYiKSkpGDixImcPnahO1auXAl3d3fU1dUJ5mGgUVdXB3d3d6xZs0YwD2Je9o3e5GW3TQsA/vnPf8LJyQnXrl3rtzlHoLm5GWPHjrWJP9xXX30VEydOFD9t9ZBLly7ZRKN//PgxBg0ahNWrVwvqYyDxu9/9Dp6ennjy5ImgPsS87B29zcseNS2LxYJJkyYhLi6uX+YchY0bN0KhUAj6KauVK1euQCaT4dNPPxXais1jtVoxbdo0TJ48WdBPWa384Q9/gJubG+7cuSO0FZunrKwMLi4uWLdundBWxLzsJb3Nyx41LQAoKSmBRCJBTk5On805AtXV1VCr1fjtb38rtJU23n//fXh4eKCmpkZoKzbNjh07IJFIcPr0aaGtAAAaGxsRFhaG5ORkoa3YPElJSQgPD0dTU5PQVgCIedlT+pKXPW5aALB06VL4+fnh8ePHvTbnKKSnpyMgIABms1loK23U19cjMDAQixYtEtqKzWI0GjFkyBD8/Oc/F9pKO3Jzc0FEnD7xeqCzf/9+MAyD/Px8oa20Q8zL7ulLXvaqadXW1sLb25vTx1cPZPbt29fpI6KFZu/evWAYBgcPHhTaik2Snp5uswGTlpaGYcOG4cGDB0JbsTlqamrg7+9vk5kk5mXX9DUve9W0AOCrr74CwzDYvn17b99q19y5cweenp546623hLbSKYsXL4aHhwfKy8uFtmJTbN26FRKJBIcPHxbaSoc8evQIgYGBSExMtIlrbbZCS0sL4uPjMWLECNTW1gptp0PEvOyY/uRlr5sW8N01EmdnZ1y8eLEvb7c7mpqaEB0djXHjxqGhoUFoO53y9OlTTJgwAT/84Q/R2NgotB2b4L///S9cXV1t6hpkR5w+fRpyuRxr164V2orNsHr1ashkMhQXFwttpUvEvGxPf/OyT02rqakJU6ZMQVhYmM2e4fDJm2++CZVKZROrBbvj6tWrcHNzw7Jly4S2IjgPHjxASEgIpk+fDovFIrSdblm3bh2kUilyc3OFtiI433zzDaRSKTZu3Ci0lW4R87I9/c3LPjUtALh79y4CAgIwZcoUwe4+twV+//vfw8nJCQcOHBDaSo/54osvIJFIHPoeILPZjMmTJ0Oj0aCyslJoOz3CarUiPT0dKpUKZ8+eFdqOYJSUlECpVOL1118fMF+Xinn5HWzkZZ+bFgCUlpbCx8cHs2fPRnNzc3+GGpBs3boVDMNg27ZtQlvpNf/4xz8GrPf+0tTUhMTERHh5eeHGjRtC2+kVTU1NSEhIGJDe2aC0tBS+vr4DMnPEvGQnL/vVtADg1KlTcHNzQ1pams3cI8EH2dnZcHJywqpVq4S20mdWrFgBJycn7N69W2grvNHY2Ih58+ZBqVTizJkzQtvpEyaTCZGRkdBoNCgrKxPaDm+UlpZixIgRiIqKGrCfVsS87H9e9rtpAcDRo0ehVCqRkJBgU/cnccVf/vIXMAyDzMxMoa30m/feew8SiQQbNmwQ2grnmEwmvPTSS1CpVDZ3T09vefDgASZNmgQ/Pz9cuHBBaDucc+7cOfj6+iIyMnLA734v5mX/YKVpAcCZM2fg4+ODyMhI3L9/n61hbQqr1YqPPvoIDMPY1SquDRs2QCKRYPny5QPmGkFvqa2tRXR0NHx9fXHu3Dmh7bCCyWRCfHw8lEoljhw5IrQdzsjPz8egQYOg1WphNBqFtsMKYl72HdaaFgDo9XoEBgYiODjYboKhlcePHyM5ORkymQy7du0S2g7rbN++HVKpFKmpqXYTDK2UlJRAo9EgKCgIpaWlQtthlWfPniE1NRUKhQJbtmwR2g6rWK1WbNiwAXK5HOnp6XZ3m4aYl32D1aYFAFVVVdBqtVAoFNi8eTPbwwtCa+j5+/ujoKBAaDuccfToUfj5+WHkyJF2MYmsVis2btwIuVyOWbNm2e3eiy0tLfjoo4/g5ORkNycdraHXeh3EXp9SIOZl72G9aQHf7XLcurTxtddeG7Bh0dzcjHXr1kEulyM+Pn7AHkdvqKqqQlxcHBQKBf76178OiPuXOqKqqgrJycmQSqVYvXq13Ybe/ycvLw9DhgxBcHAwCgsLhbbTZ44fPw6NRoOhQ4fi2LFjQtvhHDEvewcnTasVnU6HgIAAeHh4YMuWLQMqAE+ePIkJEyZAoVDgk08+cYjQa6WlpQWrVq2CXC7HpEmTcOrUKaEt9RiLxYJNmzbB3d0dgYGBOH78uNCWeKW6uhqJiYlgGAZLliwZUPsVVldXY+HChWAYBnPmzLHbaz2dIeZlz+C0aQHf3cSZmZkJmUyGiIgInDx5kmvJfnHv3j0sXrwYDMMgPj4eBoNBaEuCcf36dcTFxUEikeCNN95AVVWV0Ja65NhnxjeCAAAE9klEQVSxY5g4cSLkcjk+/PDDAbssmg327NmDYcOGwdPTE3//+99tenl1Y2MjNm3aBLVajYCAAMEfwCkkYl52D+dNq5WrV68iNjYWRAStVou8vDy+pHtEeXk5fvrTn0KhUGD48OH44osvhLZkM+Tk5MDf3x/Ozs5YtmyZzT2UMDc3FzExMSAizJw50yFvuu0Ik8mE9957D3K5HAEBAdi8eTOePn0qtK02GhoasGnTJgwbNgxyuRyZmZkOsQS8J4h52Tm8Na1W8vPzERcXByJCdHQ0cnJyBNtk1mq1orCwEIsWLYJMJoNGo8Gnn36KZ8+eCeLHlnn69Ck2b96MESNGQCaTYfHixSgqKhLMT0NDA3bv3o2oqCgQEWbNmmXXi2T6w507d/D222/D2dkZfn5++OMf/4hvv/1WMD8VFRVYvXo1fH194erqiuXLlwvqx5YR8/J5eG9arRQXF2Pu3LmQSqUYNGgQlixZgoKCAl6uHd28eRMff/wxRo4cCSLC+PHjsWPHDofcWqW3NDU14bPPPsPYsWNBRBg1ahRWrlyJW7duca5tsViQn5+PxYsXY9CgQZDJZEhOTraZJw3bOlVVVXj//ffh6ekJiUSCuLg47Ny5E0+ePOFc22g0Yvv27dBqtZBIJPDy8sIHH3wwYBcd8I2Yl/9DsKbVSm1tLbZu3YqpU6eCiKBSqTBz5kysWbMGZ8+eZeVm1/v372PPnj1Yvnw5IiIiQETw9PTE0qVLceLECRaOwjG5cuUKMjMz4efnByJCUFAQli5dip07d+Lu3busaJSVlWHr1q1ITU2Fp6cniAhjxozBmjVrUF1dzYqGo9HY2IiDBw8iNTUVcrkcTk5OiIiIQGZmJnJzc1n5CrG5uRlnz57FmjVrMHPmTCgUCigUCiQlJWHPnj12d88VX4h5CTAAQDbC9evX6ciRI5SXl0fHjx8no9FIKpWKQkJCKCQkhMLCwkij0ZBSqSSVSkVqtZrc3NyooaGBTCYTmc1mqq+vp7t375LBYCCDwUA3btyg6upqkkqlFBkZSVqtlrRaLcXExJCTk5PQh2wXWCwWOnbsGOXl5ZFOp6OzZ8+SxWKhIUOGUGhoKIWGhlJISAgNHTqUlEolubm5kUqlIldXV6qvr6e6urq2+t26dYtu3LjRVj+z2UxqtZpiYmJIq9VSfHw8hYaGCn3IdkNtbS0dPnyY8vPzSafT0a1bt0gqlVJQUFBb7UaOHEmenp7k7u5OSqWSlEolERGZzWYym81kNBrp0aNHVFpaSnq9nvR6PZWXl5PFYqHg4OC2ORcfH08eHh4CH7H94Kh5aVNN6//T0tJC58+fp0uXLrVNhBs3blBFRQU1NjZ2+j6JREJDhgxpK1xoaCiFh4dTdHQ0qVQqHo/AcTGZTFRUVERXr14lg8FAer2eDAYDVVdXk9Vq7fR9zs7OFBAQQGFhYW2N7oUXXqCJEyfazISxd27fvk2nT59um28Gg4Fu3rxJRqOxy/ep1WoKDg6m0NDQtvpFRUVRQEAAT84dG0fKS5ttWl3R3NxMZrOZ6urqyGw2k5ubW9sZoKurq9D2RLqgoaGh7Qy9vr6elEolqdVqUqlUJJVKhbYn0gWt881sNhMRtc05tVotsDORrrC3vByQTUtERERExDGRCG1ARERERESkp4hNS0RERERkwCAloltCmxAREREREekJ/wcvGgR1ROBKywAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data = dask.array.ones((10,10), chunks=(5,5))\n", "\n", "(data + 5).visualize()" ] }, { "cell_type": "code", "execution_count": 3, "id": "altered-liberal", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "6.0" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "(data + 5)[4,6].compute()" ] }, { "cell_type": "markdown", "id": "driving-wonder", "metadata": {}, "source": [ "Since I've only asked for a single element of the array, only that one chunk needs to be calculated. This is less important for this tiny example, but it does get important when you've got a huge dataset that you only need a subset of." ] }, { "cell_type": "markdown", "id": "quick-clerk", "metadata": {}, "source": [ "## Creating Dask arrays manually\n", "\n", "Something to keep in mind with Dask arrays is that you can't directly write to an array element" ] }, { "cell_type": "code", "execution_count": 4, "id": "advisory-privacy", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "ERROR Item assignment with not supported\n" ] } ], "source": [ "try:\n", " data[4, 6] = 3\n", "except Exception as e:\n", " print('ERROR', e)" ] }, { "cell_type": "markdown", "id": "foreign-advocate", "metadata": {}, "source": [ "This keeps us from trying to do calculations with loops, which Dask can't parallelise and are really inefficient in Python. Instead you should use whole-array operations.\n", "\n", "Say we want to calculate the area of grid cells, knowing the latitude and longitude of grid point centres.\n", "\n", "For each grid point we need to calculate\n", "\n", "$$A = R^2\\int_{\\phi_0}^{\\phi_1}\\int_{\\theta_0}^{\\theta_1}\\cos\\theta d\\theta d\\phi = R^2(\\phi_1 - \\phi_0)(\\cos \\theta_1 - \\cos \\theta_0)$$\n", "\n", "Rather than doing this in a loop over each grid point, we can use numpy array operations, computing the width and height of each cell and then performing an outer vector product to create a 2d array" ] }, { "cell_type": "code", "execution_count": 5, "id": "fewer-hello", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\n", "\n", "\n", "\n", "\n", "
\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
Array Chunk
Bytes 1.20 MB 80.00 kB
Shape (300, 500) (100, 100)
Count 90 Tasks 15 Chunks
Type float64 numpy.ndarray
\n", "		\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", " \n", "\n", " \n", " 500\n", " 300\n", "\n", "
" ], "text/plain": [ "dask.array" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import numpy\n", "from numpy import deg2rad\n", "\n", "# Grid centres\n", "lon = dask.array.linspace(0, 360, num=500, endpoint=False, chunks=100)\n", "lat = dask.array.linspace(-90, 90, num=300, endpoint=True, chunks=100)\n", "\n", "# Grid spacing\n", "dlon = (lon[1] - lon[0]).compute()\n", "dlat = (lat[1] - lat[0]).compute()\n", "\n", "# Grid edges in radians\n", "lon0 = deg2rad(lon - dlon/2)\n", "lon1 = deg2rad(lon + dlon/2)\n", "lat0 = deg2rad((lat - dlat/2).clip(-90, 90))\n", "lat1 = deg2rad((lat + dlat/2).clip(-90, 90))\n", "\n", "# Compute cell dimensions\n", "cell_width = lon1 - lon0\n", "cell_height = numpy.sin(lat1) - numpy.sin(lat0)\n", "\n", "# Area\n", "R_earth = 6_371_000\n", "A = R_earth**2 * numpy.outer(cell_height, cell_width)\n", "A" ] }, { "cell_type": "markdown", "id": "affecting-arbitration", "metadata": {}, "source": [ "You can see that even though only numpy operations were used we still ended up with a Dask array at the end. An array much larger than would fit in memory can be created this way, with the values only needing to be computed if they're actually used.\n", "\n", "Plotting Dask arrays works just like numpy arrays:" ] }, { "cell_type": "code", "execution_count": 6, "id": "olympic-chase", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAV0AAAEDCAYAAACWDNcwAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAYQUlEQVR4nO3df6yc1X3n8fdnrq8NBaJADMRgd82mpCtCG2gtulpWK6dJEzZFpdmKiKhNqRat+wdpEjVRsVNp06ay1lm1tJGSVr0tKEYNEG8ShIWyIYYtYpEAx06BYBwaN3hTxxau+bGYVTD3zv3sH89jMoH7Y67v3DPP3OfzkkYzc+bMM+eg5HuPv+fHI9tEREQZnWE3ICKiTRJ0IyIKStCNiCgoQTcioqAE3YiIghJ0IyIKStCNiGVJ0q2Sjkp6so+6/0rS/ZKekPSApLVL1a4E3YhYrr4IXNVn3T8FbrP988BngP+2VI1K0I2IZcn2g8DzvWWS3ibpG5L2Svrfkv5N/dElwP31678HrlmqdiXoRkSbTAC/Z/sXgU8Cf1mXPw78Rv36A8BZkt6yFA1YsRQXjYhoGklnAv8O+B+SThavqp8/CXxe0u8ADwI/BKaWoh0JuhHRFh3gRduXvf4D24eB/wSvBeffsP1/l6oRERHLnu2XgGckXQugyjvr16slnYyHW4Bbl6od8wZdSadJ2i3pcUn7JP1xXX6OpF2Svlc/n93znS2SDkh6WtL7lqrxERGzkXQH8DDws5IOSboB+E3gBkmPA/v48YTZRuBpSf8InA9sXbJ2zXe0o6rkxxm2X5Y0DjwEfIxqKP687W2SNgNn275J0iXAHcAVwAXAfcDbbXeXqhMREaNi3pGuKy/Xb8frh6n+Qmyvy7cDv16/vga40/YJ288AB6gCcERE6/U1kSZpDNgL/AzwBduPSjrf9hEA20cknVdXvxB4pOfrh+qy119zE7AJYIyxX/wp3nTqvYiI1jjOC8dsn7uYa7zvXWf4uef7+8f33idO3Gu7300W8+or6NapgcskvRm4S9Klc1TXDGVvyGHYnqBaM8ebdI5/Se/upykR0XL3+Sv/Z7HXeO75Lrvv/em+6o6t+d7qxf5erwUtGbP9oqQHqLbWPStpTT3KXQMcrasdAtb1fG0tcHiu60qis+q0hTQlItrqlcVfwsA004u/0CmYN+hKOheYrAPu6cB7gM8CO4HrgW318931V3YCt0u6mWoi7WJg9zw/AmNZvRYRZRgzOaS5/X5GumuA7XVetwPssH2PpIeBHfUyjB8A1wLY3idpB/AU1Y6OG+ddudARndNPX0Q3IqI1/t9gLtPYka7tJ4DLZyh/DpgxEWt7Kwta55aRbkSUY0x3SHdCb8Y2YAlWNKMpEdEO02+c3y+iGZFOgpUrh92KiGgJA912B11gxdiwWxERLZKRboJuRBRiYLLtOV2vbEZTImL5M056wSuyeiEiCjF0hxNzmxF0LTGdkW5EFFLtSBscSQeB40AXmLK9Yba6zYh0Aq+Y6ciGiIilILozHhOzKO+yfWy+Ss0IuoA7CboRUUY1kTacmNOIoOuOmF6V1QsRUUa1TrfvoLta0p6e9xP1KYmvv+Q3JRn46xk+f00jgi5kpBsRZU33P9I9NleOtnal7cP1ueK7JH3X9oMzVWxG0E1ONyIKWuBId/7rVXcTxvZRSXdR3S2nuUHXEt3xBN2IKMOI7oBuhi7pDKBj+3j9+r3AZ2ar34igi8BjCboRUc4C0gvzOZ/qjjpQxdTbbX9jtsqNCbrdnHcTEYUY8aoHM3lv+/vAO/ut34igazKRFhHlVJsjhrMLthFBt0ovDLsREdEmS7A5oi+NCbqZSIuIUmzRdYtHuiYj3Ygoa7rtI90h/dGJiBaqJtKGE/4aE3Snx4fdiIhoi0ykkfRCRJTVbfWBN4LpRrQkItpgkDvSFqoZoS453YgobLrNqxcg6YWIKKc68KbFQdeZSIuIgoyYHNJIrxFBN+mFiCjJprmbIyStA24D3kp1L7cJ25+T9EfAfwH+pa76Kdtfr7+zBbiB6iZtH7V973y/M6QlcxHRSmr05ogp4BO2vy3pLGCvpF31Z39u+097K0u6BLgOeAdwAXCfpLfb7s76CxnpRkRBpsEjXdtHgCP16+OS9gMXzvGVa4A7bZ8AnpF0gOoU9Yfn/J2xId2EPiJaaSQm0iStBy4HHgWuBD4i6beBPVSj4ReoAvIjPV87xNxBOut0I6Ioo0EeYr4gfYc6SWcCXwU+bvslSX8F/AnVSP1PgD8D/jPMmCh5wzBW0iZgE8DY2WdDRroRUUh1C/YGn70gaZwq4H7J9tcAbD/b8/nfAPfUbw8B63q+vhY4/Ppr1rcongBY9dPrPKQ/OhHRSmruebqqbvxzC7Df9s095WvqfC/AB4An69c7gdsl3Uw1kXYxsHvuHwGPZ6QbEWWYZu9IuxL4MPAdSY/VZZ8CPiTpMqr2HwR+F8D2Pkk7gKeoVj7cOOfKBagu0UnQjYhyGjvStf0QM+dpvz7Hd7YCW/tuhcArEnQjogxbjR7pLj2BMpEWEYVUE2lt3gYMoATdiCil5fdIQ6azYnrYrYiIlqgm0hqa0y1F2QYcEQWNxI60pSJBZywj3YgoYyR2pC0pwViCbkQU1OobUwrT6SToRkQZNkxOtznoCsZXzLN/IiJiQKr0QouDLpix7EiLiIIauyOtBAFjSnohIspo/ZIxyaxMeiEiihlsekHSGNW54j+0ffVcdRsSdGEsE2kRUdCA75H2MWA/8Kb5KjYj6GJWjU0NuxkR0RLV6oXBnL0gaS3wq1SHfP3+fPWbEXQFKzLSjYhCBrw54i+APwDO6qdyM4IuZkUm0iKioAWkF1ZL2tPzfqK+8w2SrgaO2t4raWM/F2tE0O0ApyW9EBGFLHD1wjHbG2b57Erg1yS9HzgNeJOkv7P9W7NdrBFBF5lOjnaMiIIGsXrB9hZgC0A90v3kXAEXGhJ0BazoZMlYRJRhi6k270jryJw+NjnsZkREiwx6c4TtB4AH5qvXiKAL0CHphYgoo/U70jpMc/rYq8NuRkS0SKuDrgTjWacbEYXkEHOq0W5ERCkD3gbct0YE3Q5mVSfrdCOiDBumWn2IOWZcWTIWEeW0Or0gmfGMdCOikNbndIU5TQm6EVGO2x10oZMDbyKioMZOpElaB9wGvBWYpjph53OSzgG+DKwHDgIftP1C/Z0twA1AF/io7Xvn/A3MacqOtIgow252TncK+ITtb0s6C9graRfwO8D9trdJ2gxsBm6SdAlwHfAO4ALgPklvtz3rTFnukRYRZYluU1cv2D4CHKlfH5e0H7gQuAbYWFfbTrXn+Ka6/E7bJ4BnJB0ArgAenu03pKxeiIiyRiKnK2k9cDnwKHB+HZCxfUTSeXW1C4FHer52qC57/bU2AZsAVl+wMumFiChmJM5ekHQm8FXg47ZfkmZt8EwfvOE0m/rk9QmAt/3cGc5EWkQU4yqvOwx9BV1J41QB90u2v1YXPytpTT3KXQMcrcsPAet6vr4WODzn9TFj2QYcEQU1efWCgFuA/bZv7vloJ3A9sK1+vrun/HZJN1NNpF0M7J7zNzCndZJeiIgy3OSJNKp7AH0Y+I6kx+qyT1EF2x2SbgB+AFwLYHufpB3AU1QrH26ca+UC1KsXMtKNiIIam16w/RAz52kB3j3Ld7ZS3QO+L1m9EBGljcTqhaUiYJwE3Ygow2550AVnc0REFNX4JWNLScB4gm5EFNTYnG4JWTIWESUZMd3g1QtLLiPdiChtWPcfb0zQHcst2COilEykJehGRGFtz+mOK0E3Ispp/Ug3IqIUA9PTLQ66AoYzjxgRrWSgzSNdSayc/ajIiIiBa/U6XchINyIKa3PQre6RlpFuRJSiTKR1hnSgcES0VJtHuhERRRk8oNULkk4DHgRWUcXUr9j+9Gz1GxF0q6Mdk9WNiJIG9q/rE8Av2365vrXZQ5L+p+1HZqrciKAbEVHcgNILtg28XL8drx+zXr0RQVcoE2kRUVb/QXe1pD097yfqu5m/RtIYsBf4GeALth+d7WKNCLoAnaQXIqKUhW2OOGZ7w5yXq+4DeZmkNwN3SbrU9pMz1U2ki4hWqm7ZM/9jYdf0i8ADwFWz1WnESLfaBpz0QkQUNLjVC+cCk7ZflHQ68B7gs7PVb0TQjYgobYAHG64Bttd53Q6ww/Y9s1VuSNAVY0qmIyIKMYNcvfAEcHm/9RsSdCMiSlK7TxmLiCiuzduAq4m0pBcioqAh3Qu3EUE3IqKoIR5iPu/wUtKtko5KerKn7I8k/VDSY/Xj/T2fbZF0QNLTkt63VA2PiFgMub/HoPXzb/ovMvNC3z+3fVn9+DqApEuA64B31N/5y3oZRUREs7jPx4DNG3RtPwg83+f1rgHutH3C9jPAAeCKRbQvImJZWczs1UckPVGnH86uyy4E/rmnzqG67A0kbZK0R9Kef3muu4hmREQsXJPTCzP5K+BtwGXAEeDP6vKZMtMzNtv2hO0Ntjec+5ZkICKiIFNtA+7nMWCnFHRtP2u7a3sa+Bt+nEI4BKzrqboWOLy4JkZELIGm5nRnImlNz9sPACdXNuwErpO0StJFwMXA7sU1MSJi8IaVXph3na6kO4CNVAf5HgI+DWyUdBnV34GDwO8C2N4naQfwFDAF3FifMxkR0SxN3ZFm+0MzFN8yR/2twNbFNCoiYsk1NeiWYGB6WHvyIqJ1lip10I9GBN2IiOKWYGVCPxJ0I6KVWj7SNV0nvRARBbU76EZEFNT2nG41kTak/wIR0U5tDroREaWp7YeYZ8lYRLRBI4KuMV0nvRARBSW9EBFRSNsn0iATaRFRWJuDroFugm5ElNTmoAsZ6UZEOaLlqxcMmUiLiHKS0yULxiKirDYHXePkdCOirFYHXcNkYm5EFNT69EJERFFtDrpGTHo4BwpHRAu55asXALok6EZEQQMa6UpaB9wGvJVqTcCE7c/NVr8RQbfaHJGgGxHlDDCnOwV8wva3JZ0F7JW0y/ZTM1VuTNCddGfYzYiINhlQ0LV9BDhSvz4uaT9wIdDkoCu6JOhGRCFmIUF3taQ9Pe8nbE/MVFHSeuBy4NHZLtaIoAswnYm0iChELCi9cMz2hnmvKZ0JfBX4uO2XZqvXiKBrxKseG3YzIqJFBrlOV9I4VcD9ku2vzVW3IUEXJknQjYiCBrd6QcAtwH7bN89XvxlB12IyI92IKGlwI90rgQ8D35H0WF32Kdtfn6nyvEFX0q3A1cBR25fWZecAXwbWAweBD9p+of5sC3AD0AU+avve+X6jWjKWibSIKGSAp4zZfgj6X/Paz0j3i8DnqRb/nrQZuN/2Nkmb6/c3SboEuA54B3ABcJ+kt9vuzv0TopslYxFRUlO3Adt+sF4G0esaYGP9ejvwAHBTXX6n7RPAM5IOAFcAD8/1G9OIVzy+oIZHRCzGqG0DPr9eEIztI5LOq8svBB7pqXeoLnsDSZuATQCrL1jJdEa6EVHQcjllbKa8xoxdqxcXTwBc9HNnOkvGIqKYhW2OGKhTDbrPSlpTj3LXAEfr8kPAup56a4HD813MFieSXoiIkkYs6O4Erge21c9395TfLulmqom0i4Hd812sukda0gsRUcYCd6QNVD9Lxu6gmjRbLekQ8GmqYLtD0g3AD4BrAWzvk7SD6qCHKeDG+VcuVDvSMpEWESVpejhRt5/VCx+a5aN3z1J/K7B1IY0wZCItIsoZwZzuQFV3jshEWkSU09j0Qgm2eGU66YWIKKjVQTcj3YgorN0jXXKebkQU1u6gm3W6EVFQ2+8GbMPkdFYvREQZjV6nW8I0HX7UXTnsZkREm7ih63RLmc4t2COioFaPdKv0QlYvREQhbd8cMY14pZuJtIgop9UTaVhZMhYRRbU66E4Dr3Qb0ZSIaAPT7ok0I6Zy4E1EFNT6ibSprNONiJJaHXQRJ5JeiIhCWr85woZuRroRUYrd3EPMy1CCbkSU1eaR7rTh1W42R0REOa1OL1Qj3azTjYhCqvNkh/LTjQi6NkxOZaQbEQW1eaRrxHRyuhFRULvTC86dIyKirFavXrChO5n0QkQU0vZTxgA8pMMnIqJ9qs0Rg4m6km4FrgaO2r50vvoNCbrCWb0QESUNbqD3ReDzwG39VG5G0DW4m4m0iChnUCNd2w9KWt9v/QYF3Yx0I6KQheV0V0va0/N+wvbEqf70ooKupIPAcaALTNneIOkc4MvAeuAg8EHbL8x5IYOmEnQjopQFnb1wzPaGQf3yIEa677J9rOf9ZuB+29skba7f3zT3JQTJ6UZEScvoEPNrgI316+3AA8wXdA0kvRARpXh0b9dj4JuSDPx1nec43/YRANtHJJ3Xz4U63UW2JCJiIQa3ZOwOqoHmakmHgE/bvmW2+osNulfaPlwH1l2SvruAhm4CNgGMnX12RroRUdaAsgu2P7SQ+osKurYP189HJd0FXAE8K2lNPcpdAxyd5bsTwATAqnXrPKyhfkS0k6aHE3ROOehKOgPo2D5ev34v8BlgJ3A9sK1+vnveaxk6kxnpRkQhZpCbIxZkMSPd84G7JJ28zu22vyHpW8AOSTcAPwCunfdKQ0xqR0T7CA9sc8RCnXLQtf194J0zlD8HvHuh19PUqbYkIuIUjFrQHaiMdCOitDYHXQHKkrGIKGVEc7qDY+hMDrsREdEmI7d6YaCSXoiIotzy9IKhk4m0iCjFtDvoQnK6EVFY23O6SS9EREkjt053oDKRFhGltTnoZslYRBRlQ7f1qxeGdD/kiGinNo90MYy9OuxGRESrtDnoiox0I6IgA0OKOY0IuhnpRkRZBrc9p9vNSDciCjHtnkiTF3Q75IiIxWtzThfD2IkE3YgoqNVBl0ykRURJOfAGTSXoRkQhBtp9tKPpTGZLWkQU1OaRLhnpRkRRLd8GLJvOqzlQNyIKMbj163SncrZjRBTU7h1phm5yuhFRULtzukYnkl6IiELstq9eAKYy0o2Igto+0uXVnHgTEaUYDyml2ZygO5X0QkQU0vqjHYe4Zi4iWmq5LRmTdBXwOWAM+Fvb22atPG2mf/SjpWpKRMRPMODlNNKVNAZ8AfgV4BDwLUk7bT814xeGeJO4iGghL79DzK8ADtj+PoCkO4FrgBmDrm08lXuwR0Q5y20i7ULgn3veHwJ+qbeCpE3ApvrtiV1TX35yidoyTKuBY8NuxIAtxz5B+jVKfnaxFzjOC/fe56+s7rP6QP/7LVXQ1QxlP5FAsT0BTABI2mN7wxK1ZWiWY7+WY58g/RolkvYs9hq2rxpEW05FZ4muewhY1/N+LXB4iX4rImJkLFXQ/RZwsaSLJK0ErgN2LtFvRUSMjCVJL9iekvQR4F6qJWO32t43x1cmlqIdDbAc+7Uc+wTp1ygZ6T7JQ9p/HBHRRkuVXoiIiBkk6EZEFDT0oCvpKklPSzogafOw29MvSbdKOirpyZ6ycyTtkvS9+vnsns+21H18WtL7htPq+UlaJ+nvJe2XtE/Sx+ryke2bpNMk7Zb0eN2nP67LR7ZPJ0kak/QPku6p3y+HPh2U9B1Jj51cHrYc+vUa20N7UE2y/RPwr4GVwOPAJcNs0wLa/h+AXwCe7Cn778Dm+vVm4LP160vqvq0CLqr7PDbsPszSrzXAL9SvzwL+sW7/yPaNat34mfXrceBR4N+Ocp96+vb7wO3APcvof4MHgdWvKxv5fp18DHuk+9p2YduvAie3Czee7QeB519XfA2wvX69Hfj1nvI7bZ+w/QxwgKrvjWP7iO1v16+PA/updhiObN9cebl+O14/zAj3CUDSWuBXgb/tKR7pPs1h2fRr2EF3pu3CFw6pLYNwvu0jUAUv4Ly6fCT7KWk9cDnVyHCk+1b/M/wx4Ciwy/bI9wn4C+APgN6TW0a9T1D9QfympL31cQGwPPoFDP883Xm3Cy8TI9dPSWcCXwU+bvslaaYuVFVnKGtc32x3gcskvRm4S9Klc1RvfJ8kXQ0ctb1X0sZ+vjJDWaP61ONK24clnQfskvTdOeqOUr+A4Y90l9t24WclrQGon4/W5SPVT0njVAH3S7a/Vhcvi77ZfhF4ALiK0e7TlcCvSTpIlZb7ZUl/x2j3CQDbh+vno8BdVOmCke/XScMOusttu/BO4Pr69fXA3T3l10laJeki4GJg9xDaNy9VQ9pbgP22b+75aGT7JunceoSLpNOB9wDfZYT7ZHuL7bW211P9/+Z/2f4tRrhPAJLOkHTWydfAe4EnGfF+/YRhz+QB76eaIf8n4A+H3Z4FtPsO4AgwSfXX9gbgLcD9wPfq53N66v9h3cengf847PbP0a9/T/XPsyeAx+rH+0e5b8DPA/9Q9+lJ4L/W5SPbp9f1byM/Xr0w0n2iWsn0eP3YdzImjHq/eh/ZBhwRUdCw0wsREa2SoBsRUVCCbkREQQm6EREFJehGRBSUoBsRUVCCbkREQf8f7gnKGepsFm8AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "plt.pcolormesh(A)\n", "plt.colorbar()" ] }, { "cell_type": "markdown", "id": "activated-state", "metadata": {}, "source": [ "And we can check the total area is about $5.1 \\times 10^{14}$" ] }, { "cell_type": "code", "execution_count": 7, "id": "cleared-auditor", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "'5.100645e+14'" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "'%e'%(A.sum().compute())" ] }, { "cell_type": "markdown", "id": "sustained-snake", "metadata": {}, "source": [ "## Creating Dask arrays from NetCDF files\n", "\n", "The most common way of creating a Dask array is to read them from a netcdf file with Xarray. You can give `open_dataset()` and `open_mfdataset()` a `chunks` parameter, which is how large chunks should be in each dimension of the file.\n", "\n", "If you use `open_mfdataset()`, by default each input file will be its own chunk." ] }, { "cell_type": "code", "execution_count": 8, "id": "outdoor-collectible", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'tas' (time: 1980, lat: 144, lon: 192)>\n",
       "dask.array<open_dataset-0e410aadf9156f5cf0b6cd8c2849d202tas, shape=(1980, 144, 192), dtype=float32, chunksize=(1, 144, 192), chunktype=numpy.ndarray>\n",
       "Coordinates:\n",
       "  * time     (time) datetime64[ns] 1850-01-16T12:00:00 ... 2014-12-16T12:00:00\n",
       "  * lat      (lat) float64 -89.38 -88.12 -86.88 -85.62 ... 86.88 88.12 89.38\n",
       "  * lon      (lon) float64 0.9375 2.812 4.688 6.562 ... 353.4 355.3 357.2 359.1\n",
       "    height   float64 ...\n",
       "Attributes:\n",
       "    standard_name:  air_temperature\n",
       "    long_name:      Near-Surface Air Temperature\n",
       "    comment:        near-surface (usually, 2 meter) air temperature\n",
       "    units:          K\n",
       "    cell_methods:   area: time: mean\n",
       "    cell_measures:  area: areacella\n",
       "    history:        2019-11-08T06:41:45Z altered by CMOR: Treated scalar dime...\n",
       "    _ChunkSizes:    [  1 144 192]
" ], "text/plain": [ "\n", "dask.array\n", "Coordinates:\n", " * time (time) datetime64[ns] 1850-01-16T12:00:00 ... 2014-12-16T12:00:00\n", " * lat (lat) float64 -89.38 -88.12 -86.88 -85.62 ... 86.88 88.12 89.38\n", " * lon (lon) float64 0.9375 2.812 4.688 6.562 ... 353.4 355.3 357.2 359.1\n", " height float64 ...\n", "Attributes:\n", " standard_name: air_temperature\n", " long_name: Near-Surface Air Temperature\n", " comment: near-surface (usually, 2 meter) air temperature\n", " units: K\n", " cell_methods: area: time: mean\n", " cell_measures: area: areacella\n", " history: 2019-11-08T06:41:45Z altered by CMOR: Treated scalar dime...\n", " _ChunkSizes: [ 1 144 192]" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import xarray\n", "\n", "path = 'https://dapds00.nci.org.au/thredds/dodsC/fs38/publications/CMIP6/CMIP/CSIRO-ARCCSS/ACCESS-CM2/historical/r1i1p1f1/Amon/tas/gn/v20191108/tas_Amon_ACCESS-CM2_historical_r1i1p1f1_gn_185001-201412.nc'\n", "ds = xarray.open_dataset(path, chunks={'time': 1})\n", "ds.tas" ] }, { "cell_type": "markdown", "id": "promotional-transsexual", "metadata": {}, "source": [ "There are a few ways to turn a Dask xarray.DataArray back into a numpy array. `.load()` will compute the dask data and retain metadata, `.values` will compute the dask data and return a numpy array, and `.data` will return the dask array itself." ] }, { "cell_type": "code", "execution_count": 9, "id": "settled-yellow", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.DataArray 'tas' ()>\n",
       "array(295.87695312)\n",
       "Coordinates:\n",
       "    time     datetime64[ns] 1850-01-16T12:00:00\n",
       "    lat      float64 -26.88\n",
       "    lon      float64 186.6\n",
       "    height   float64 2.0\n",
       "Attributes:\n",
       "    standard_name:  air_temperature\n",
       "    long_name:      Near-Surface Air Temperature\n",
       "    comment:        near-surface (usually, 2 meter) air temperature\n",
       "    units:          K\n",
       "    cell_methods:   area: time: mean\n",
       "    cell_measures:  area: areacella\n",
       "    history:        2019-11-08T06:41:45Z altered by CMOR: Treated scalar dime...\n",
       "    _ChunkSizes:    [  1 144 192]
" ], "text/plain": [ "\n", "array(295.87695312)\n", "Coordinates:\n", " time datetime64[ns] 1850-01-16T12:00:00\n", " lat float64 -26.88\n", " lon float64 186.6\n", " height float64 2.0\n", "Attributes:\n", " standard_name: air_temperature\n", " long_name: Near-Surface Air Temperature\n", " comment: near-surface (usually, 2 meter) air temperature\n", " units: K\n", " cell_methods: area: time: mean\n", " cell_measures: area: areacella\n", " history: 2019-11-08T06:41:45Z altered by CMOR: Treated scalar dime...\n", " _ChunkSizes: [ 1 144 192]" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ds.tas[0,50,99].load()" ] }, { "cell_type": "code", "execution_count": 10, "id": "accessible-reputation", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array(295.87695312)" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ds.tas[0,50,99].values" ] }, { "cell_type": "code", "execution_count": 11, "id": "blond-malta", "metadata": {}, "outputs": [ { "data": { "text/html": [ "\n", "\n", "\n", "\n", "\n", "
\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
Array Chunk
Bytes 218.97 MB 110.59 kB
Shape (1980, 144, 192) (1, 144, 192)
Count 1981 Tasks 1980 Chunks
Type float32 numpy.ndarray
\n", "		\n", "\n", "\n", " \n", " \n", " \n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", " \n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", " \n", " \n", "\n", " \n", " \n", "\n", " \n", " \n", " \n", "\n", " \n", " \n", " \n", "\n", " \n", " \n", "\n", " \n", " 192\n", " 144\n", " 1980\n", "\n", "
" ], "text/plain": [ "dask.array" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ds.tas.data" ] }, { "cell_type": "markdown", "id": "naval-impossible", "metadata": {}, "source": [ "## Distributed Dask\n", "\n", "Without any special setup, Dask will run operations in threaded mode. You can configure it to run in distributed mode instead with" ] }, { "cell_type": "code", "execution_count": 9, "id": "revolutionary-antenna", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n" ] }, { "data": { "text/html": [ "\n", "\n", "\n", "\n", "\n", "
\n", "

Client

\n", "\n", "		\n", "

Cluster

\n", "
    \n", "
  • Workers: 4
    • \n", "
  • Cores: 12
    • \n", "
  • Memory: 17.13 GB
    • \n", "
\n", "
" ], "text/plain": [ "" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" }, { "name": "stderr", "output_type": "stream", "text": [ "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n", "distributed.comm.tcp - WARNING - Could not set timeout on TCP stream: [Errno 92] Protocol not available\n" ] } ], "source": [ "import dask.distributed\n", "import tempfile\n", "\n", "try:\n", " client\n", "except NameError:\n", " dask_worker_dir = tempfile.TemporaryDirectory()\n", " \n", " client = dask.distributed.Client(\n", " local_directory = dask_worker_dir.name,\n", " )\n", "client" ] }, { "cell_type": "markdown", "id": "mental-algebra", "metadata": {}, "source": [ "This will by default ask for all resources available on your computer.\n", "\n", "```{hint}\n", "The `try: ... except NameError:` structure is to make sure only on Dask client is created, in case you execute the notebook cell more than once. If you're writing a python script rather than using a notebook it's not needed.\n", "```\n", "\n", "```{warning}\n", "It's important to set the `local_directory` parameter, otherwise Dask will store temporary files in the current working directory which can be a problem if filesystem quotas are enabled.\n", "```\n", "\n", "Other useful options are:\n", " * `n_workers`: Number of distributed processes\n", " * `threads_per_worker`: Number of shared memory threads within each process\n", " * `memory_limit`: Memory available to each process (e.g. `'4gb'`)\n", " \n", "```{warning}\n", "If you're using a shared system, be polite and don't take over the whole system with your Dask cluster, set reasonable limits. If running on NCI's Gadi supercomputer, `climtas.nci.GadiClient()` will inspect the PBS resources requested by `qsub` and set up the cluster using those limits\n", "```\n", "\n", "You can follow the dashboard link displayed by Jupyter to get an interactive view of what the Dask cluster is doing.\n", "\n", "To stop the Dask cluster run" ] }, { "cell_type": "code", "execution_count": 7, "id": "pleased-representative", "metadata": {}, "outputs": [], "source": [ "client.close()" ] }, { "cell_type": "markdown", "id": "assumed-yahoo", "metadata": {}, "source": [ "This isn't normally needed, Dask will clean itself up at the end of your script automatically, but it can be helpful if you're experimenting with different cluster sizes." ] }, { "cell_type": "code", "execution_count": null, "id": "north-december", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python [conda env:analysis3] *", "language": "python", "name": "conda-env-analysis3-py" }, "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.6" } }, "nbformat": 4, "nbformat_minor": 5 }