to an output file. The storage array, used for programming convenience in collecting the 2-D average from all the relevant cores, became excessively large at high resolution and processor numbers, and we abandoned the previous method in favour of saving chunks of the column density array into separate files, each maintained by the root processors in the $z$-direction. A similar method was implemented for $y$-averages and $x$-averages.

![Code speed as a function of simulation time (top panel) and maximum particle number on any core (bottom panel) for $256^3$ resolution on 2048 cores. Standard domain decomposition (SDD) quickly becomes unbalanced with particles and achieves only the speed of the particle-laden mid-plane cores. With the Particle Block Domain Decomposition (PBDD) scheme the speed stays close to its optimal value, and the particle number per core (bottom panel) does not rise much beyond $10^4$.[]{data-label="f:speed_t"}](speed_t.eps){width="8.7cm"}

The above-mentioned global arrays had been used in the code for programming convenience and did not imply excessive memory usage at moderate resolution and processor numbers. Purging those arrays in favour of loops or smaller 1-D arrays was relatively straight-forward.

### Particle migration {#s:parmig}

At the end of a sub-time-step each processor checks if any particles have left its 