Scheme of integration

Equations (8 - 9) are solved using finite differences with a constant step in space $\Delta$x = $\Delta$y = $\Delta$z = 1 and a constant step in the expansion parameter $\Delta$a . We use the ``leap-frog'' scheme to advance coordinates and velocities from one moment to another. (In the following we drop tildas for all dimensionless variables.) At any moment a_n = a₀ + n$\Delta$a , we have the coordinates x _n and the potential $\phi_{n}^{}$ . Velocities p _{n - 1/2} are defined at a_{n - 1/2} = a_n - $\Delta$a/2 . The coordinates and the velocities for the next moment are found using:

 

$$\nabla \phi_{n}^{} \Delta$$ p _{n + 1/2} =

$${F(a_{n+1/2})\over a_{n+1/2}^2} \Delta$$ x _{n + 1} = (14)

In order to solve eq.(9) we approximate the Laplacian operator using the 7-point ``crest'' template:

 

$$\nabla^{2}_{} \phi \approx \phi_{i\pm 1,j,k}^{} \phi_{i,j\pm 1,k}^{} \phi_{i,j,k\pm 1}^{} \phi_{i,j,k}^{}$$ (15)

where (i,j,k) = 1,...,N_grid . This leads to a large system of linear equations relating unknown variables $\phi_{i,j,k}^{}$ with known right-hand side of the discrete form of the Poisson equation 3$\Omega_{0}^{}$($\rho_{i,j,k}^{}$ - 1)/2a . The system of equations is solved exactly by the FFT technique.

The density on the mesh $\rho_{i,j,k}^{}$ is obtained from particle positions using the Cloud-In-Cell method. In order to find the ``acceleration'' g = - $\nabla$$\phi$ for each particle, the gravitational potential is differentiated on the mesh:

 

$$\phi_{i+1,j,k}^{} \phi_{i-1,j,k}^{}$$ (16)

Then the acceleration is interpolated to the position of the particle using a three-linear interpolation. This scheme for the force interpolation (the same interpolation as in the density assignment) is very important because it does not produce a force acting on the particle itself. (Thus, an isolated point does not produce a force at the position of the particle). While this might sound like a natural condition for any realistic method, only two methods - PM and TREE - do not have this self-force. In P ³M the effect is minimized. In the case of multigrid methods the self-force cannot be excluded - only minimized. Typically this is achieved by placing extended buffer zones around regions with high resolution (e.g. Kravtsov 1997). No precautions were made in AP ³M method, which might result in spurious effects in regions were multi-level grids are introduced.

Thus, the main scheme of the PM method consists of the following four blocks repeated every time step:

• Find density on the mesh using the Cloud-In-Cell technique.
• Solve the Poisson equation using two three-dimensional FFTs.
• Advance velocities and coordinates of the particles.
• Advance time and print results.