A parallel multigrid solver develped by James Bordner

---

Grids: Computational grids used are 3D Cartesian meshes, and PDEs are discretized using finite differences or finite volumes. Parallelism is obtained by partitioning the grid into subregions, and assigning grid points within each subregion to a processor in a 3D Cartesian processor grid. MPI is used for transferring data between processors.

Solvers: The main solver in MGMPI is multigrid. Multigrid methods are capable of solving three-dimensional elliptic boundary value problems in time linearly proportional to the number of grid points. In addition to a multigrid V-cycle method, MGMPI will include numerous Krylov subspace solvers, including CG, BiCG, BiCG-STAB, CGNR, CGNE, and CGS. These will be able to be used alone, or with multigrid as a preconditioner.

Robustness:To increase the robustness of the multigrid algorithms, MGMPI will include support for semi-coarsening, line and plane smoothers, operator-based prolongation and restriction, and Galerkin coarsening. The user will have detailed control over the construction of the multigrid cycling strategy, including multiple semi-coarsened grids, alternating line or plane smoothing directions, etc.

Efficiency: In addition to robustness, efficiency is important as well. Multiple matrix storage schemes are implemented, allowing the matrix to be tailored to specific PDE characteristics (such as symmetry or constant coefficients), as well as to characteristics of the computer memory hierarchy. Furthermore, blocking and fusion techniques will be used to improve performance of the smoothing operators on cache-based machines. And, while the bulk of the code is written in C++, Fortran versions of many of the computational kernels are provided, which typically run at about twice the speed of the corresponding C++ versions.

Adaptivity: A long-term goal is to add adaptivity to MGMPI at several levels. One, to add adaptive mesh refinement (AMR) capabilities so that the computational grid can adapt to varying resolution scales within the PDE domain. Two, to make the algorithms adaptable to the PDE being solved; for example, by the automatic modification of the number of smoothings, smoothing algorithm, the grid hierarchy, etc., during a solve, or between sequences of solves. And three, to make the implementation parameters (blocking size, etc.) adaptable to the machine on which the solution is being computed, allowing the code efficiency to increase as more iterations are computed. These levels of adaptability will be especially useful in implicit schemes for time-dependent simulations where a sequence of related elliptic PDEs is solved.