update dask

day2_hpcenv/4_parallel-programming/2_dask/1-Introduction-to-Dask.ipynb

-{
- "cells": [
-  {
-   "cell_type": "markdown",
-   "id": "linear-bangkok",
-   "metadata": {},
-   "source": [
-    "### HIGH THROUGHPUT COMPUTING WITH DASK\n",
-    "\n",
-    "**Organisers:** Alan O’Cais, David Swenson  \n",
-    "**Website:** https://www.cecam.org/workshop-details/1022\n",
-    "\n",
-    "**Synopsis:**\n",
-    "High-throughput (task-based) computing is a flexible approach to parallelisation. It involves splitting a problem into loosely-coupled tasks. A scheduler then orchestrates the parallel execution of those tasks, allowing programs to adaptively scale their resource usage. E-CAM has extended the data-analytics framework Dask with a capable and efficient library to handle such workloads. This workshop will be held as a series of virtual seminars/tutorials on tools in the Dask HPC ecosystem.\n",
-    "\n",
-    "**Programme:**\n",
-    "- 21 January 2021, 3pm CET (2pm UTC): Dask - a flexible library for parallel computing in Python\n",
-    "  - YouTube link: https://youtu.be/Tl8rO-baKuY\n",
-    "  - GitHub Repo: https://github.com/jacobtomlinson/dask-video-tutorial-2020\n",
-    "\n",
-    "- 4 February 2021, 3pm CET (2pm UTC): Dask-Jobqueue - a library that integrates Dask with standard HPC queuing systems, such as SLURM or PBS\n",
-    "  - YouTube link: https://youtu.be/iNxhHXzmJ1w\n",
-    "  - GitHub Repo: https://github.com/ExaESM-WP4/workshop-Dask-Jobqueue-cecam-2021-02\n",
-    "\n",
-    "- 11 February 2021, 3pm CET (2pm UTC) : Jobqueue-Features - a library that enables functionality aimed at enhancing scalability\n",
-    "  - YouTube link: https://youtu.be/FpMua8iJeTk\n",
-    "  - GitHub Repo: https://github.com/E-CAM/jobqueue_features_workshop_materials"
-   ]
-  }
- ],
- "metadata": {
-  "kernelspec": {
-   "display_name": "Python 3 (ipykernel)",
-   "language": "python",
-   "name": "python3"
-  },
-  "language_info": {
-   "codemirror_mode": {
-    "name": "ipython",
-    "version": 3
-   },
-   "file_extension": ".py",
-   "mimetype": "text/x-python",
-   "name": "python",
-   "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython3",
-   "version": "3.8.5"
-  }
- },
- "nbformat": 4,
- "nbformat_minor": 5
-}

day2_hpcenv/4_parallel-programming/2_dask/1_dask_MonteCarloPi.ipynb

+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Dask local cluster example"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## What is Dask? (https://docs.dask.org/en/latest/)\n",
+    "\n",
+    "* combine a blocked algorithm approach\n",
+    "* with dynamic and memory aware task scheduling\n",
+    "* to realise a parallel out-of-core NumPy clone\n",
+    "* optimized for interactive computational workloads\n",
+    "\n",
+    "-----------------------------------"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### WORKSHOP on DASK - HIGH THROUGHPUT COMPUTING WITH DASK\n",
+    "\n",
+    "**Organisers:** Alan O’Cais, David Swenson  \n",
+    "**Website:** https://www.cecam.org/workshop-details/1022\n",
+    "\n",
+    "**Synopsis:**\n",
+    "High-throughput (task-based) computing is a flexible approach to parallelisation. It involves splitting a problem into loosely-coupled tasks. A scheduler then orchestrates the parallel execution of those tasks, allowing programs to adaptively scale their resource usage. E-CAM has extended the data-analytics framework Dask with a capable and efficient library to handle such workloads. This workshop will be held as a series of virtual seminars/tutorials on tools in the Dask HPC ecosystem.\n",
+    "\n",
+    "**Programme:**\n",
+    "- 21 January 2021, 3pm CET (2pm UTC): Dask - a flexible library for parallel computing in Python\n",
+    "  - YouTube link: https://youtu.be/Tl8rO-baKuY\n",
+    "  - GitHub Repo: https://github.com/jacobtomlinson/dask-video-tutorial-2020\n",
+    "\n",
+    "- 4 February 2021, 3pm CET (2pm UTC): Dask-Jobqueue - a library that integrates Dask with standard HPC queuing systems, such as SLURM or PBS\n",
+    "  - YouTube link: https://youtu.be/iNxhHXzmJ1w\n",
+    "  - GitHub Repo: https://github.com/ExaESM-WP4/workshop-Dask-Jobqueue-cecam-2021-02\n",
+    "\n",
+    "- 11 February 2021, 3pm CET (2pm UTC) : Jobqueue-Features - a library that enables functionality aimed at enhancing scalability\n",
+    "  - YouTube link: https://youtu.be/FpMua8iJeTk\n",
+    "  - GitHub Repo: https://github.com/E-CAM/jobqueue_features_workshop_materials\n",
+    "  \n",
+    "------------------------------------"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Example problem: Monte-Carlo estimate of $\\pi$\n",
+    "\n",
+    "<img src=\"https://upload.wikimedia.org/wikipedia/commons/8/84/Pi_30K.gif\" width=\"25%\" align=left alt=\"PI monte-carlo estimate\"/>\n",
+    "\n",
+    "## Problem description\n",
+    "\n",
+    "Suppose we want to estimate the number $\\pi$ using a [Monte-Carlo method](https://en.wikipedia.org/wiki/Pi#Monte_Carlo_methods), i.e. obtain a numerical estimate based on a random sampling approach, and that we want at least single precision floating point accuracy.\n",
+    "\n",
+    "We take advantage of the fact that the area of a quarter circle with unit radius is $\\pi/4$ and that hence the probability of a randomly chosen point inside a unit square to lie within that circle is $\\pi/4$ as well.\n",
+    "\n",
+    "So for N randomly chosen pairs $(x, y)$ with $x\\in[0, 1)$ and $y\\in[0, 1)$ we count the number $N_{circ}$ of pairs that also satisfy $(x^2 + y^2) < 1$ and estimage $\\pi \\approx 4 \\cdot N_{circ} / N$."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Monte-Carlo estimate with NumPy on a single CPU\n",
+    "\n",
+    "* NumPy is the fundamental package for scientific computing with Python (https://numpy.org/).\n",
+    "* It contains a powerful n-dimensional array object and useful random number capabilities."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "import numpy"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "def calculate_pi_single(size_in_bytes):\n",
+    "    \n",
+    "    \"\"\"Calculate pi using a Monte Carlo method.\"\"\"\n",
+    "    \n",
+    "    rand_array_shape = (int(size_in_bytes / 8 / 2), 2)\n",
+    "    \n",
+    "    # 2D random array with positions (x, y)\n",
+    "    xy = numpy.random.uniform(low=0.0, high=1.0, size=rand_array_shape)\n",
+    "    \n",
+    "    # check if position (x, y) is in unit circle\n",
+    "    xy_inside_circle = (xy ** 2).sum(axis=1) < 1\n",
+    "\n",
+    "    # pi is the fraction of points in circle x 4\n",
+    "    pi = 4 * xy_inside_circle.sum() / xy_inside_circle.size\n",
+    "\n",
+    "    print(f\"\\nfrom {xy.nbytes / 1e9} GB randomly chosen positions\")\n",
+    "    print(f\"   pi estimate: {pi}\")\n",
+    "    print(f\"   pi error: {abs(pi - numpy.pi)}\\n\")\n",
+    "    \n",
+    "    return pi"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Let's calculate...\n",
+    "\n",
+    "Observe how the error decreases with an increasing number of randomly chosen positions!"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "%time pi = calculate_pi_single(size_in_bytes=10_000_000) # 10 MB\n",
+    "%time pi = calculate_pi_single(size_in_bytes=100_000_000) # 100 MB\n",
+    "%time pi = calculate_pi_single(size_in_bytes=1_000_000_000) # 1 GB"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Are we already better than single precision floating point resolution?"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "numpy.finfo(numpy.float32)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## We won't be able to scale the problem to several Gigabytes or Terabytes!"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Problems\n",
+    "\n",
+    "* slowness of the numpy-only single CPU approach! (we could scale the problem using the [multiprocessing](https://docs.python.org/3.8/library/multiprocessing.html) and/or [threading](https://docs.python.org/3.8/library/threading.html) libraries)\n",
+    "* frontend/login node compute resources are shared and CPU, memory (and IO bandwidth) user demands will collide"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Monte-Carlo estimate with Dask on multiple CPUs\n",
+    "\n",
+    "We define a Dask cluster with 8 CPUs and 24 GB of memory."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "import dask.distributed"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "cluster = dask.distributed.LocalCluster(\n",
+    "    n_workers=1, threads_per_worker=8, memory_limit=24e9,\n",
+    "    ip=\"0.0.0.0\"\n",
+    ")\n",
+    "\n",
+    "client = dask.distributed.Client(cluster)\n",
+    "client"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Use dask.array for randomly chosen positions"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "import numpy, dask.array"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "def calculate_pi_dask(size_in_bytes, number_of_chunks):\n",
+    "    \n",
+    "    \"\"\"Calculate pi using a Monte Carlo method.\"\"\"\n",
+    "    \n",
+    "    array_shape = (int(size_in_bytes / 8 / 2), 2)\n",
+    "    chunk_size = (int(array_shape[0] / number_of_chunks), 2)\n",
+    "    \n",
+    "    # 2D random positions array using dask.array\n",
+    "    xy = dask.array.random.uniform(\n",
+    "        low=0.0, high=1.0, size=array_shape,\n",
+    "        # specify chunk size, i.e. task number\n",
+    "        chunks=chunk_size )\n",
+    "  \n",
+    "    xy_inside_circle = (xy ** 2).sum(axis=1) < 1\n",
+    "\n",
+    "    pi = 4 * xy_inside_circle.sum() / xy_inside_circle.size\n",
+    "    \n",
+    "    # start Dask calculation\n",
+    "    pi = pi.compute()\n",
+    "\n",
+    "    print(f\"\\nfrom {xy.nbytes / 1e9} GB randomly chosen positions\")\n",
+    "    print(f\"   pi estimate: {pi}\")\n",
+    "    print(f\"   pi error: {abs(pi - numpy.pi)}\\n\")\n",
+    "    display(xy)\n",
+    "    \n",
+    "    return pi"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Let's calculate again...\n",
+    "Observe the wall time decreases of the 1 Gigabyte and 10 Gigabyte random sample $\\pi$ estimates!"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "%time pi = calculate_pi_dask(size_in_bytes=1_000_000_000, number_of_chunks=10) # 1 GB"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "%time pi = calculate_pi_dask(size_in_bytes=10_000_000_000, number_of_chunks=100) # 10 GB"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Let's go larger than memory...\n",
+    "Because Dask splits the computation into single managable tasks, we can scale up easily!"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "tags": []
+   },
+   "outputs": [],
+   "source": [
+    "%time pi = calculate_pi_dask(size_in_bytes=100_000_000_000, number_of_chunks=250) # 100 GB"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Are we now better than single precision floating point resolution?\n",
+    "Not at all, if we require an order of magnitude better..."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "numpy.finfo(numpy.float32)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## We could increase the local cluster CPU resources...\n",
+    "However, the above Dask cluster size is always limited by the memory/CPU resources of a single compute node."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# %time pi = calculate_pi(size_in_bytes=1_000_000_000_000, number_of_chunks=2_500) # 1 TB"
+   ]
+  }
+ ],
+ "metadata": {
+  "anaconda-cloud": {},
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.10.4"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}