Parallelisation

Minim can be parallelised at two levels:

  • MPI distributes a state across multiple processors, which may be on different nodes of a cluster. This is the primary form of parallelisation and is used to study large systems.

  • OpenMP shares the work of a single processor between multiple threads on the same node. This is used to speed up a few key operations without any communication cost.

Both are enabled by default when the library is built with make.

Running with MPI

Parallelisation is enabled at compile time with the -DPARALLEL flag, which the Makefile includes by default, so no extra setup is needed. Since MPI requires initialising, user programs must call mpi.init before creating any states:

int main(int argc, char **argv) {
  mpi.init(&argc, &argv);
  ...
}

The number of mpi ranks and the current rank can be reference by mpi.size and mpi.rank, respectively.

A program is then run in the usual way with mpirun:

mpirun -np 4 ./run.exe

Everything proceeds identically whether the program is run on one processor or many: the energy and gradient functions always return the total values over the whole system, and the coordinates can always be read and written in full. The parallelisation is abstracted away by the state’s Communicator, which handles all of the communication between processors in the background.

The degrees of freedom of each State are distributed over the processors automatically when the state is created; there is no need to manually split the data. All of the built-in minimisers and potentials work in parallel without modification. By default a state is spread over all of the MPI ranks, but a subset can also be used by passing the ranks argument to its constructor, which allows different states to be minimised simultaneously on different processors.

For grid potentials, the number of processors along each grid dimension can be controlled with setCommArray, e.g.:

PhaseField pot;
pot.setGridSize({nx, ny, 1});
pot.setCommArray({2, 2, 1});

which splits the grid over a 2×2 array of processors (and so requires at least four ranks). By default it will attempt to automatically set the number along each dimension to make the sub-grids as cubic as possible.

Using OpenMP

OpenMP is enabled at compile time with the -fopenmp flag, which the Makefile includes by default. It is used to parallelise some of the more intensive computations that are local to each processor, such as element loops in the phase field potentials, the vector operations used by L-BFGS, and the dot products of grid data.

The number of threads is controlled in the usual way with the OMP_NUM_THREADS environment variable:

OMP_NUM_THREADS=4 mpirun -np 4 ./run.exe

which runs four MPI ranks with four OpenMP threads each.

Running under Slurm

On a cluster managed by Slurm, the three modes of running are:

MPI only — one MPI rank per core, no OpenMP threads. Use --ntasks to request the number of ranks and launch with srun:

#!/bin/bash
#SBATCH --nodes=1
#SBATCH --ntasks=32
#SBATCH --time=01:00:00

srun ./run.exe

OpenMP only — a single MPI rank with several threads. Use --ntasks=1 with --cpus-per-task for the thread count, and set OMP_NUM_THREADS to match:

#!/bin/bash
#SBATCH --nodes=1
#SBATCH --ntasks=1
#SBATCH --cpus-per-task=16
#SBATCH --time=01:00:00

export OMP_NUM_THREADS=$SLURM_CPUS_PER_TASK
srun ./run.exe

Hybrid MPI + OpenMP — several ranks, each with several threads. The total core count is --ntasks times --cpus-per-task. Set OMP_NUM_THREADS to the number of threads per rank, and bind each rank (and its threads) to its own cores so that threads belonging to different ranks do not interfere:

#!/bin/bash
#SBATCH --nodes=2
#SBATCH --ntasks=8
#SBATCH --cpus-per-task=8
#SBATCH --time=01:00:00

export OMP_NUM_THREADS=$SLURM_CPUS_PER_TASK
export OMP_PROC_BIND=true
export OMP_PLACES=cores
srun --cpubind=cores ./run.exe

The same layouts apply when launching with mpirun instead of srun (e.g. mpirun -np 32 ./run.exe); inside a Slurm job, mpirun will usually detect the allocated ranks automatically.

Beware of oversubscription: requesting more ranks or threads than the allocated cores will slow the program down. In particular, leaving OMP_NUM_THREADS unset may cause each rank to spawn as many threads as there are cores on the node, i.e. ntasks-fold more than intended.

How the parallelisation is implemented

Data distribution

Each state holds a Communicator, created automatically from its potential, which defines how the degrees of freedom are divided between the processors. Each processor is assigned one block: a set of degrees of freedom specific to that processor. To compute the energy and gradient of its block, a processor may also require degrees of freedom owned by its neighbours, so it additionally stores a halo of degrees of freedom from the neighbouring blocks. Together, the block plus halo form the proc data held by that processor.

There are two types of communicator.

CommGrid is used by potentials defined on a structured grid, where the grid is split into a commArray of sub-grids, one per processor, each surrounded by a halo of the given haloWidth.

CommUnstructured is used by element-based potentials such as BarAndHinge, where the degrees of freedom are split into equal-sized blocks in the order they are given. For these, the energy and gradient are broken into elements, each depending on only a few degrees of freedom, and each element is assigned to the processor that holds the most of its degrees of freedom. The remaining degrees of freedom of that element then lie in the halo. Defining a potential in terms of elements therefore allows the halo regions to be determined automatically. For quantities that cannot be split into elements, such as a system-wide volume constraint, a potential may instead define a blockEnergyGradient function that receives the communicator and performs its own communication.

When creating a state, the communicator sets up the block and halo layout, remaps the elements to local indices, and creates MPI derived datatypes describing the halo and edge regions of the local data.

Communication

All communication is handled by the communicator. The key operations are:

  • communicate sends the edges of each block into the halo regions of the neighbouring processors, updating stale halo data (e.g. after randomly perturbing the coordinates).

  • communicateAccumulate instead adds the values of a processor’s halo onto the corresponding entries of the neighbours’ blocks. This is how gradient contributions from elements that span two blocks are combined.

  • gather and scatter convert between the global data and the local proc data.

  • Reductions such as sum, norm and dotProduct compute a result from the local blocks and sum it over all processors, so every processor receives the total.

Since halo entries are duplicated between processors, operations that combine data across processors (like dotProduct) only count each block’s entries once, excluding the halo.

The state’s energy and gradient functions build on these operations. The block functions compute a single processor’s contribution without any communication, so summing the block energies over the processors gives the total. The plain energy / gradient functions do this summing (and the equivalent gather for the gradient) internally, which is why they can be used as if the program were serial. See the State documentation for the full set of functions and their parallel behaviour.

Where MPI communication can be avoided, it is. For example, the gradient returned by procGradient is correct over the whole proc region (block plus halo), so the minimisers can take their steps using only local data and only need to communicate when the next gradient is required.

Hybrid MPI + OpenMP

OpenMP is used within a single processor to parallelise operations that require no communication, such as loops over the local energy elements and the vector arithmetic in L-BFGS. It is applied with #pragma omp parallel for directives, so it requires no changes to user code. The dot product over grid data is a special case: rather than excluding the halo entries one by one, each processor sums over its interior and the result is reduced over the processors with MPI.