XM_40017/notebooks/mpi_tutorial.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "8d800917",
   "metadata": {},
   "source": [
    "# Tutorial: Using MPI in Julia\n",
    "Message Passing Interface (MPI) is a standardized and portable library specification for communication between parallel processes in distributed memory systems. Julia offers a convenient way to work with MPI for creating efficient parallel and distributed applications. In this tutorial, you will learn how to use MPI from Julia to perform parallel computing tasks.\n",
    "\n",
    "## MPI launches separate Julia instances\n",
    "When you run an MPI-enabled Julia script, MPI takes care of spawning multiple instances of the Julia executable, each acting as a separate process. These workers can communicate with each other using MPI communication functions. This enables parallel processing and distributed computation. Here's a summary of how it works:\n",
    "\n",
    "-- TODO: insert picture here --\n",
    "\n",
    "-  **MPI Spawns Processes**: The `mpiexec` command launches multiple instances of the Julia executable, creating separate worker processes. In this example, 4 Julia workers are spawned.\n",
    "\n",
    "-  **Worker Communication**: These workers can communicate with each other using MPI communication functions, allowing them to exchange data and coordinate actions.\n",
    "\n",
    "- **Parallel Tasks**: The workers execute parallel tasks simultaneously, working on different parts of the computation to potentially speed up the process.\n",
    "\n",
    "\n",
    " \n",
    "\n",
    "\n",
    "## Installing MPI.jl and MPIClusterManagers Packages\n",
    "To use MPI in Julia, you'll need the MPI.jl package, and if you intend to run MPI programs in a Jupyter Notebook, you'll also need the MPIClusterManagers package. These packages provide the necessary bindings to the MPI library and cluster management capabilities. To install the packages, open a terminal and run the following commands:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3cb5f151",
   "metadata": {},
   "outputs": [],
   "source": [
    "using Pkg\n",
    "Pkg.add(\"MPI\")\n",
    "Pkg.add(\"MPIClusterManagers\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ed45a4b2",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-block alert-info\">\n",
    "    <b>Tip:</b>\n",
    "The package <code>MPI.jl</code> is the Julia interface to MPI. Note that it is not a MPI library by itself. It is just a thin wrapper between MPI and Julia. To use this interface, you need an actual MPI library installed in your system such as OpenMPI or MPICH. Julia downloads and installs a MPI library for you, but it is also possible to use a MPI library already available in your system. This is useful, e.g., when running on HPC clusters. See the documentation of <code>MPI.jl</code> for further details.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7a36916e",
   "metadata": {},
   "source": [
    "## Writing a HelloWorld MPI Program in Julia\n",
    "Let's start by creating a simple MPI program that prints a message along with the rank of each worker. \n",
    "\n",
    "Create a new Julia script, for example, `mpi_hello_world.jl`:\n",
    "\n",
    "```julia\n",
    "using MPI\n",
    "\n",
    "# Initialize MPI\n",
    "MPI.Init()\n",
    "\n",
    "# Get the default communicator (MPI_COMM_WORLD) for all processes\n",
    "comm = MPI.COMM_WORLD\n",
    "\n",
    "# Get the number of processes in this communicator\n",
    "nranks = MPI.Comm_size(comm)\n",
    "\n",
    "# Get the rank of the current process within the communicator\n",
    "rank = MPI.Comm_rank(comm)\n",
    "\n",
    "# Print a message with the rank of the current process\n",
    "println(\"Hello, I am process $rank of $nranks processes!\")\n",
    "\n",
    "# Finalize MPI\n",
    "MPI.Finalize()\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6caa8d74",
   "metadata": {},
   "source": [
    "### MPI Communicators\n",
    "In MPI, a **communicator** is a context in which a group of processes can communicate with each other. `MPI_COMM_WORLD` is one of the MPI standard communicators, it represents all processes in the MPI program. Custom communicators can also be created to group processes based on specific requirements or logical divisions. \n",
    "\n",
    "The **rank** of a processor is a unique identifier assigned to each process within a communicator. It allows processes to distinguish and address each other in communication operations. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "19f41e38",
   "metadata": {},
   "source": [
    "## Running the HelloWorld MPI Program\n",
    "\n",
    "To run MPI applications in parallel, you need a launcher like `mpiexec`. MPI codes written in Julia are not an exception to this rule. From the system terminal, you can run\n",
    "```\n",
    "$ mpiexec -np 4 mpi_hello_world.jl\n",
    "```\n",
    "In this command, `-np 4` specifies the desired number of processes. \n",
    "But it will probably not work since the version of `mpiexec` needs to match with the MPI version we are using from Julia. You can find the path to the `mpiexec` binary you need to use with these commands\n",
    "\n",
    "```julia\n",
    "julia> using MPI\n",
    "julia> MPI.mpiexec_path\n",
    "```\n",
    "\n",
    "and then try again\n",
    "```\n",
    "$ /path/to/my/mpiexec -np 4 julia mpi_hello_world.jl\n",
    "```\n",
    "with your particular path.\n",
    "\n",
    "However, this is not very convenient. Don't worry if you could not make it work! A more elegant way to run MPI code is from the Julia REPL directly, by using these commands:\n",
    "```julia\n",
    "julia> using MPI\n",
    "julia> mpiexec(cmd->run(`$cmd -np 4 julia mpi_hello_world.jl`))\n",
    "```\n",
    "\n",
    "Now, you should see output from 4 ranks.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0592e58c",
   "metadata": {},
   "source": [
    "## Running MPI Programs in Jupyter Notebook with MPIClusterManagers\n",
    "If you want to run your MPI code from a Jupyter Notebook, you can do so using the `MPIClusterManagers` package.\n",
    "\n",
    "1. Load the packages and start an MPI cluster with the desired number of workers:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cf66dd39",
   "metadata": {},
   "outputs": [],
   "source": [
    "using MPIClusterManagers\n",
    "# Distributed package is needed for addprocs()\n",
    "using Distributed\n",
    "\n",
    "manager = MPIWorkerManager(4)\n",
    "addprocs(manager)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d40fe3ee",
   "metadata": {},
   "source": [
    "2. Run your MPI code inside a `@mpi_do` block to execute it on the cluster workers:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0a51d1f2",
   "metadata": {},
   "outputs": [],
   "source": [
    "@mpi_do manager begin\n",
    "    using MPI\n",
    "    comm = MPI.COMM_WORLD\n",
    "    rank = MPI.Comm_rank(comm)\n",
    "    println(\"Hello from process $rank\")\n",
    "end\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "38ed88c1",
   "metadata": {},
   "source": [
    "MPI is automatically initialized and finalized within the `@mpi_do` block.\n",
    "\n",
    "3. Remove processes when done:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e0b53cc1",
   "metadata": {},
   "outputs": [],
   "source": [
    "rmprocs(manager)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5466a650",
   "metadata": {},
   "source": [
    "## Point-to-Point Communication with MPI\n",
    "MPI provides point-to-point communication using blocking send and receiving functions `MPI.send`, `MPI.recv`; or their non-blocking versions `MPI.Isend`, and `MPI.Irecv!`. These functions allow individual processes to send and receive data between each other.\n",
    "\n",
    "### Blocking communication\n",
    "\n",
    "Let's demonstrate how to send and receive with an example:\n",
    "\n",
    "```julia\n",
    "using MPI\n",
    "\n",
    "MPI.Init()\n",
    "\n",
    "comm = MPI.COMM_WORLD\n",
    "rank = MPI.Comm_rank(comm)\n",
    "\n",
    "# Send and receive messages using blocking MPI.send and MPI.recv\n",
    "if rank == 0\n",
    "    data = \"Hello from process $rank !\"\n",
    "    MPI.send(data, comm, dest=1)\n",
    "elseif rank == 1\n",
    "    received_data = MPI.recv(comm, source=0)\n",
    "    println(\"Process $rank received: $received_data\")\n",
    "end\n",
    "\n",
    "MPI.Finalize()\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d4dfe654",
   "metadata": {},
   "source": [
    "In this example, process 0 sends a message using `MPI.send`, and process 1 receives it using `MPI.recv`.\n",
    "\n",
    "### Non-blocking communication\n",
    "\n",
    "To demonstrate asynchronous communication, let's modify the example using `MPI.Isend` and `MPI.Irecv!`:\n",
    "\n",
    "```julia\n",
    "using MPI\n",
    "\n",
    "MPI.Init()\n",
    "\n",
    "comm = MPI.COMM_WORLD\n",
    "rank = MPI.Comm_rank(comm)\n",
    "\n",
    "# Asynchronous communication using MPI.Isend and MPI.Irecv!\n",
    "if rank == 0\n",
    "    data = \"Hello from process $rank !\"\n",
    "    request = MPI.Isend(data, comm, dest=1)\n",
    "    # Other computation can happen here\n",
    "    MPI.Wait(request)\n",
    "elseif rank == 1\n",
    "    received_data = Array{UInt8}(undef, 50)  # Preallocate buffer\n",
    "    request = MPI.Irecv!(received_data, comm, source=0)\n",
    "    # Other computation can happen here\n",
    "    MPI.Wait(request)\n",
    "    println(\"Process $rank received: $(String(received_data))\")\n",
    "end\n",
    "\n",
    "MPI.Finalize()\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "024db538",
   "metadata": {},
   "source": [
    "In this example, process 0 uses `MPI.Isend` to send the message asynchronously. This function returns immediately, allowing the sender process to continue its execution. However, the actual sending of data is done asynchronously in the background. Similar to `MPI.Isend`, `MPI.Irecv!` returns immediately, allowing the receiver process to continue executing. \n",
    "\n",
    "<div class=\"alert alert-block alert-warning\">\n",
    "<b>Important:</b> In asynchronous communication, always use <code>MPI.Wait()</code> to ensure the communication is finished before accessing the send or receive buffer.\n",
    "</div>\n",
    "\n",
    "\n",
    "## Collective Communication with MPI\n",
    "MPI provides collective communication functions for communication involving multiple processes. Let's explore some of these functions:\n",
    "\n",
    "- MPI.Gather: Gathers data from all processes to a single process.\n",
    "- MPI.Scatter: Distributes data from one process to all processes.\n",
    "- MPI.Bcast: Broadcasts data from one process to all processes.\n",
    "- MPI.Barrier: Synchronizes all processes.\n",
    "\n",
    "\n",
    "Let's illustrate the use of `MPI.Gather` and `MPI.Scatter` with an example:\n",
    "\n",
    "```julia\n",
    "# TODO: check if this runs correctly\n",
    "using MPI\n",
    "using Random\n",
    "\n",
    "MPI.Init()\n",
    "\n",
    "comm = MPI.COMM_WORLD\n",
    "rank = MPI.Comm_rank(comm)\n",
    "size = MPI.Comm_size(comm)\n",
    "\n",
    "# Root processor generates random data\n",
    "data = rand(rank == 0 ? size * 2 : 0)\n",
    "\n",
    "# Scatter data to all processes\n",
    "local_data = Vector{Float64}(undef, 2)\n",
    "MPI.Scatter!(data, local_data, comm, root=0)\n",
    "\n",
    "# Compute local average\n",
    "local_average = sum(local_data) / length(local_data)\n",
    "\n",
    "# Gather local averages at the root processor\n",
    "gathered_averages = Vector{Float64}(undef, size)\n",
    "MPI.Gather!(local_average, gathered_averages, comm, root=0)\n",
    "\n",
    "if rank == 0\n",
    "    # Compute global average of sub-averages\n",
    "    global_average = sum(gathered_averages) / size\n",
    "    println(\"Global average: $global_average\")\n",
    "end\n",
    "\n",
    "MPI.Finalize()\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e65cb53f",
   "metadata": {},
   "outputs": [],
   "source": [
    "using MPI\n",
    "using Random\n",
    "\n",
    "# TODO: check if this runs correctly\n",
    "\n",
    "MPI.Init()\n",
    "\n",
    "comm = MPI.COMM_WORLD\n",
    "rank = MPI.Comm_rank(comm)\n",
    "size = MPI.Comm_size(comm)\n",
    "\n",
    "# Root processor generates random data\n",
    "data = rand(rank == 0 ? size * 2 : 0)\n",
    "\n",
    "# Scatter data to all processes\n",
    "local_data = Vector{Float64}(undef, 2)\n",
    "MPI.Scatter!(data, local_data, comm, root=0)\n",
    "\n",
    "# Compute local average\n",
    "local_average = sum(local_data) / length(local_data)\n",
    "\n",
    "# Gather local averages at the root processor\n",
    "gathered_averages = Vector{Float64}(undef, size)\n",
    "MPI.Gather!(local_average, gathered_averages, comm, root=0)\n",
    "\n",
    "if rank == 0\n",
    "    # Compute global average of sub-averages\n",
    "    global_average = sum(gathered_averages) / size\n",
    "    println(\"Global average: $global_average\")\n",
    "end\n",
    "\n",
    "MPI.Finalize()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dfd5da9e",
   "metadata": {},
   "source": [
    "In this example, the root processor generates random data and then scatters it to all processes using MPI.Scatter. Each process calculates the average of its local data, and then the local averages are gathered using MPI.Gather. The root processor computes the global average of all sub-averages and prints it."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5e8f6e6a",
   "metadata": {},
   "source": [
    "# License\n",
    "\n",
    "TODO: replace link to website\n",
    "\n",
    "This notebook is part of the course [Programming Large Scale Parallel Systems](http://localhost:8000/) at Vrije Universiteit Amsterdam and may be used under a [CC BY 4.0](https://creativecommons.org/licenses/by/4.0/) license."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c9364808",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.9.1",
   "language": "julia",
   "name": "julia-1.9"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.9.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}