Parallelization of the Hydro Routines

As this is a parallel code, we divide up this task among processors in a rather straightforward way. In everything which follows, ``NPROCS'' represents the number of processors used in the calculation. First, we break up the $N$ particles into NPROCS subsets, each with approximately the same number of particles. As an example, for a simulation with 10,000 particles and four processors, we define the lower and upper limits `` ${\tt n\_lower}$'' and `` ${\tt n\_upper}$'' for each MPI process such that the rank 0 process has values of $1$ and $2500$, respectively, rank 1 has $2501$ and $5000$, etc. Each processor is responsible for the piece of the total particle summation defined by ${\tt n\_lower} \le i \le {\tt n\_upper}$. Note that since there is no restriction on particle $j$, our ``gather-scatter'' method can contribute density contributions to particles within the processor's particle range as well as to particles outside of the range, defeating efforts to use HPF-style ``DO INDEPENDENT'' formulations. After the partial sums are calculated, we use a ``MPI_ALLREDUCE'' command to total up the partial sums and redistribute all of the values to all of the processors.

The routine we use to calculate the ranges of particles fow which each processor is responsible is ``GRAVQUANT'', which also calculates the sectioning of the 3-d gravity for input into our FFT algorithm, discussed below in Sec. 3.3.3. It is ABSOLUTELY necessary to call GRAVQUANT when starting the SPH code, without exceptions. See Sec. 4 on initial conditions for more on this topic.