Parallel Computation With yt¶

How do I run my yt job on a subset of available processes¶

You can set the communicator keyword in the enable_parallelism() call to a specific MPI communicator to specify a subset of availble MPI processes. If none is specified, it defaults to COMM_WORLD.

Creating Parallel and Serial Sections in a Script¶

Many yt operations will automatically run in parallel (see the next section for a full enumeration), however some operations, particularly ones that print output or save data to the filesystem, will be run by all processors in a parallel script. For example, in the script above the lines print v,c and p.save() will be run on all 16 processors. This means that your terminal output will contain 16 repetitions of the output of the print statement and the plot will be saved to disk 16 times (overwritten each time).

yt provides two convenience functions that make it easier to run most of a script in parallel but run some subset of the script on only one processor. The first, is_root(), returns True if run on the ‘root’ processor (the processor with MPI rank 0) and False otherwise. One could rewrite the above script to take advantage of is_root() like so:

import yt
yt.enable_parallelism()

ds = yt.load("RD0035/RedshiftOutput0035")
v, c = ds.find_max("density")
p = yt.ProjectionPlot(ds, "x", "density")
if yt.is_root():
    print v, c
    p.save()

The second function, only_on_root() accepts the name of a function as well as a set of parameters and keyword arguments to pass to the function. This is useful when the serial component of your parallel script would clutter the script or if you like writing your scripts as a series of isolated function calls. I can rewrite the example from the beginning of this section once more using only_on_root() to give you the flavor of how to use it:

import yt
yt.enable_parallelism()

def print_and_save_plot(v, c, plot, print=True):
    if print:
       print v, c
    plot.save()

ds = yt.load("RD0035/RedshiftOutput0035")
v, c = ds.find_max("density")
p = yt.ProjectionPlot(ds, "x", "density")
yt.only_on_root(print_and_save_plot, v, c, plot, print=True)

Spatial Decomposition¶

During this process, the index will be decomposed along either all three axes or along an image plane, if the process is that of projection. This type of parallelism is overall less efficient than grid-based parallelism, but it has been shown to obtain good results overall.

The following operations use spatial decomposition:

• Halo Finding
• Volume Rendering: Making 3D Photorealistic Isocontoured Images

Grid Decomposition¶

The alternative to spatial decomposition is a simple round-robin of data chunks, which could be grids, octs, or whatever chunking mechanism is used by the code frontend begin used. This process allows yt to pool data access to a given data file, which ultimately results in faster read times and better parallelism.

The following operations use chunk decomposition:

• Projections (see Available Objects)
• Slices (see Available Objects)
• Cutting planes (see Available Objects)
• Derived Quantities (see Processing Objects: Derived Quantities)
• 1-, 2-, and 3-D profiles (see Profiles and Histograms)
• Isocontours & flux calculations (see 3D Surfaces and Sketchfab)

Parallelization over Multiple Objects and Datasets¶

If you have a set of computational steps that need to apply identically and independently to several different objects or datasets, a so-called embarrassingly parallel task, yt can do that easily. See the sections below on Parallelizing over Multiple Objects and Parallelization over Multiple Datasets (including Time Series).

for dataset in dataset_series:
    <process>

yt.enable_parallelism()
for dataset in dataset_series.piter():
    <process>

yt.enable_parallelism()
my_dictionary = {}
for dataset in dataset_series.piter(storage=my_dictionary):
    <process>
    sto.result = <some information processed for this dataset>
    sto.result_id = <some identfier for this dataset>

print my_dictionary

Chunk Decomposition¶

In general, these types of parallel calculations scale very well with number of processors. They are also fairly memory-conservative. The two limiting factors is therefore the number of chunks in the dataset, and the speed of the disk the data is stored on. There is no point in running a parallel job of this kind with more processors than chunks, because the extra processors will do absolutely nothing, and will in fact probably just serve to slow down the whole calculation due to the extra overhead. The speed of the disk is also a consideration - if it is not a high-end parallel file system, adding more tasks will not speed up the calculation if the disk is already swamped with activity.

The best advice for these sort of calculations is to run with just a few processors and go from there, seeing if it the runtime improves noticeably.

Projections, Slices, and Cutting Planes

Projections, slices and cutting planes are the most common methods of creating two-dimensional representations of data. All three have been parallelized in a chunk-based fashion.

• Projections: projections are parallelized utilizing a quad-tree approach. Data is loaded for each processor, typically by a process that consolidates open/close/read operations, and each grid is then iterated over and cells are deposited into a data structure that stores values corresponding to positions in the two-dimensional plane. This provides excellent load balancing, and in serial is quite fast. However, the operation by which quadtrees are joined across processors scales poorly; while memory consumption scales well, the time to completion does not. As such, projections can often be done very fast when operating only on a single processor! The quadtree algorithm can be used inline (and, indeed, it is for this reason that it is slow.) It is recommended that you attempt to project in serial before projecting in parallel; even for the very largest datasets (Enzo 1024^3 root grid with 7 levels of refinement) in the absence of IO the quadtree algorithm takes only three minutes or so on a decent processor.
• Slices: to generate a slice, chunks that intersect a given slice are iterated over and their finest-resolution cells are deposited. The chunks are decomposed via standard load balancing. While this operation is parallel, it is almost never necessary to slice a dataset in parallel, as all data is loaded on demand anyway. The slice operation has been parallelized so as to enable slicing when running in situ.
• Cutting planes: cutting planes are parallelized exactly as slices are. However, in contrast to slices, because the data-selection operation can be much more time consuming, cutting planes often benefit from parallelism.

Object-Based¶

Like chunk decomposition, it does not help to run with more processors than the number of objects to be iterated over. There is also the matter of the kind of work being done on each object, and whether it is disk-intensive, cpu-intensive, or memory-intensive. It is up to the user to figure out what limits the performance of their script, and use the correct amount of resources, accordingly.

Disk-intensive jobs are limited by the speed of the file system, as above, and extra processors beyond its capability are likely counter-productive. It may require some testing or research (e.g. supercomputer documentation) to find out what the file system is capable of.

If it is cpu-intensive, it’s best to use as many processors as possible and practical.

For a memory-intensive job, each processor needs to be able to allocate enough memory, which may mean using fewer than the maximum number of tasks per compute node, and increasing the number of nodes. The memory used per processor should be calculated, compared to the memory on each compute node, which dictates how many tasks per node. After that, the number of processors used overall is dictated by the disk system or CPU-intensity of the job.

Domain Decomposition¶

The various types of analysis that utilize domain decomposition use them in different enough ways that they are be discussed separately.

Halo-Finding

Halo finding, along with the merger tree that uses halo finding, operates on the particles in the volume, and is therefore mostly chunk-agnostic. Generally, the biggest concern for halo finding is the amount of memory needed. There is subtle art in estimating the amount of memory needed for halo finding, but a rule of thumb is that the HOP halo finder is the most memory intensive (HaloFinder()), and Friends of Friends (FOFHaloFinder()) being the most memory-conservative. It has been found that parallelHF() needs roughly 1 MB of memory per 5,000 particles, although recent work has improved this and the memory requirement is now smaller than this. But this is a good starting point for beginning to calculate the memory required for halo-finding. For more information, see Halo Finding.

Volume Rendering

The simplest way to think about volume rendering, is that it load-balances over the i/o chunks in the dataset. Each processor is given roughly the same sized volume to operate on. In practice, there are just a few things to keep in mind when doing volume rendering. First, it only uses a power of two number of processors. If the job is run with 100 processors, only 64 of them will actually do anything. Second, the absolute maximum number of processors is the number of chunks. In order to keep work distributed evenly, typically the number of processors should be no greater than one-eighth or one-quarter the number of processors that were used to produce the dataset. For more information, see Volume Rendering: Making 3D Photorealistic Isocontoured Images.