Selecting Solvers and Solver Parameters

Solvers are used when the field is defined implicitly, i.e., when the there is a relation between the field values at various locations and what one is given. For example, in an electrostatic simulation, the potential is found by solving Poisson’s equation,

\[- \nabla^2 \phi = \rho\]

This partial differential quation is then discretized to obtain a large linear equation that in the problem interior relates a linear combination of the values of the potentials at a node and nearby nodes to the value of the potential at that node. At the boundaries, the value of of the potential is directly related to the desired boundary value. Consequently, Vorpal must solve a large linear system, where the number of independent values, i.e., the length of the vector, is the number of field values in the problem.

As an example, in a Poisson solve for a 10 cell by 10 cell problem, because the potential must also be known at the node above the last cell, there are \(11 \times 11 = 121\) values of potential to solve for. The matrix for this solve therefore has \(121 \times 121 = 14641\) elements. One can see that these matrices become very large as the problem size increases.

Vorpal gives the user a fair number of choices for solving these problems. The first choice is whether to use a direct solver, which uses methods like LU decomposition, or an iterative solver, which finds the solution through successive matrix operations that converge to the solution. This first choice generally depends on the size of the problem. For problems that are too large, one cannot hold the matrix in memory, and so one generally goes with an iterative method. This is not definitive, however, as for smaller problems, an iterative method may still get one to solution more rapidly.

When using a direct solver, it is important to ensure that the matrix is not singular. In almost all cases this is true, but with fully periodic boundary conditions, the Laplace matrix is singular, as can be seen by the fact that the constant function is in its null space. One can make the matrix nonsingular by applying a boundary condition at a single cell. However, this problem has no solution when the total charge in the system does not vanish. Thus, one must ensure that the system is overall neutral. This can be enforced in Vorpal by adding a neutralizing background charge density.

Iterative solvers have no problem with periodic boundary conditions, and with them one need not impose any single-cell bounary condition. Iterative solvers works very well on the Poisson problem. However, for iterative solvers one much choose a tolerance, which is a measure of the residual reduction compared to the initial residual, is the stopping criterion for the iterative solver. Values of \(10^{-8}\) often are sufficient for giving meaningful results. If the solver does not converge, increase the tolerance. If the resulting potentials miss a lot of small scale structure, reduce the tolerance.

For interative solvers, one must also choose a preconditioner. Preconditioners transform the linear system used in the Poisson solver into systems with more favorable convergence behavior. For a simple fully periodic system, a pre-conditioner may not be necessary. For most other cases, using a pre-conditioner significantly improves the convergence behavior. Multi-grid preconditioning tends to yield the best convergence behavior and is especially good for Poisson’s equation.

For a comparison of different pre-conditioners and solvers for electrostatic simulations, see P. Messmer et.al. [messmer2004parallel].

The Visual Setup user interface by default sets up what are believed to be the best solvers, preconditioners, and parameters. However, it also allows the freedom to try others. Experimentation may be needed to arrive at the fastest running simulation.