Complex-valued coarse solver

[sebastiangrimberg]

The multigrid smoothers have recently been upgraded to smooth using the true complex-valued system matrix for frequency domain applications which, while more expensive per iteration, improves convergence by roughly halving the number of iterations in most observed cases. At high frequencies, the solver is expected to perform worse due to the indefinite nature of the system matrix but even limiting the number of multigrid levels to 1 (with "MGMaxLevels": 1) and using a sparse direct solver does not converge (or takes a large number of iterations). It can be tested that if the system matrix is real (remove lumped ports, absorbing boundary conditions), the solver converges in a constant number of iterations with 1 multigrid level (1 or 2 iterations based on stability of the direct solve and the desired precision), and performs very well even at large-ish frequencies with > 1 multigrid level. These experiments can be validated with the cpw example.

The conclusion is that we probably stand to gain in terms of convergence if the coarse solve (sparse direct solve) is performed on the complex-valued matrix instead of the real approximation A_r + A_i as is done currently. There are two relatively straightforward ways to go about this:

1. Link to complex-valued STRUMPACK, SuperLU_DIST, and MUMPS and upgrade the MFEM interfaces to support input for a pair of HypreParMatrix representing the real and imaginary parts.
2. Continue to use the real-valued libraries, but pass the 2x2 block matrix [A_r -A_i; A_i A_r] for the equivalent real formulation to the solver. This will be a more expensive factorization and solve by ~2x than the complex-valued one, but is perhaps easier to implement for testing.

[sebastiangrimberg]

This also allows us to more thoroughly investigate complex-shifted Laplacian approaches for the high-frequency regime (it is widely stated that the real shift -omega^2 M -> +omega^2 M hurts convergence at high frequency).

[sebastiangrimberg]

We can also keep an eye here on the eigenvalue bounds computed by the Chebyshev smoother at each level. For very high frequencies, the coarser levels may become negative definite after which and sort of shifted-Laplacian approach is counter productive.