                    Finite Element Discretization Library
                                   __
                       _ __ ___   / _|  ___  _ __ ___
                      | '_ ` _ \ | |_  / _ \| '_ ` _ \
                      | | | | | ||  _||  __/| | | | | |
                      |_| |_| |_||_|   \___||_| |_| |_|

                               https://mfem.org

This directory contains some sophisticated solvers for discrete problems
generated by MFEM.

# Block Solvers

The miniapp block-solvers compares the relative performance of some solvers for
the discrete saddle point problems arising from the mixed finite element
discretization of the second order scalar elliptic problem
                   -div(K grad p) = r.h.s.
The user can customize the problem settings by providing the mesh, the
coefficient K, and the essential/natural boundary assignment.

The discrete saddle point problem has a block structure of the form

                   [ M  B^T ] [u] = [g]
                   [ B   0  ] [p] = [f]

The solvers for the above block system include:

1. The divergence free solver (both coupled and decoupled modes).

   The idea of the solver is to exploit a multilevel decomposition of the
   Raviart-Thomas finite element space to find a particular flux solution
   satisfying the divergence constraint B u_p = f, and then solve the remaining
   divergence-free component of the flux u in the kernel space of the discrete
   divergence operator B. In our implementation, these two flux components and
   the pressure p can be solved one by one, or altogether simultaneously by
   switching the flag DFSParameters::coupled_solve (see div_free_solver.hpp).

   The multilevel decompositions of the Raviart-Thomas and discrete L2 space are
   given by DFSData, see the documentation in div_free_solver.hpp for more
   details. In block-solvers.cpp, the data are generated geometrically through
   mesh refinement and collection of "TransferOperators" of the spaces, This
   process is done in DFSSpaces. In general, the multilevel decompositions can
   also be obtained via element agglomeration (e.g., use ParElag).

   For more details see

   [1] Vassilevski, Multilevel Block Factorization Preconditioners (Appendix F.3),
       Springer, 2008.
   [2] Voronin, Lee, Neumuller, Sepulveda, Vassilevski, Space-time discretizations
       using constrained first-order system least squares (CFOSLS).
       J. Comput. Phys. 373: 863-876, 2018.

2. MINRES preconditioned by a block diagonal preconditioner.

   This solver is the one used in examples/ex5p.cpp.

   First, define the approximate Schur complement S = B diag(M)^{-1} B^T. Then,
   BoomerAMG is used as a preconditioner for S, denoted by AMG(S). Lastly, the
   block diagonal preconditioner is defined as

             P^{-1} = [ diag(M)^{-1}     0    ]
                      [  0             AMG(S) ]

# Low-Order Refined Solvers

The miniapp `lor` (and its parallel counterpart `plor`) demonstrate the use of
the `LOR` and `LORSolver` classes for the low-order refined preconditioning of
high-order diffusion problems. Preconditioners for H1, H(curl), H(div), and L2
finite element spaces are constructed.

In serial, a direct solver is used for the LOR system if MFEM is compiled with
SuiteSparse support (otherwise a simple Gauss-Seidel smoother is constructed
using the low-order system).

In parallel, hypre's scalable AMG solvers are used for the low-order systems:
`HypreBoomerAMG` is used for H1 and L2 spaces, `HypreAMS` is used for H(curl)
(and H(div) in 2D), and `HypreADS` is used for H(div) in 3D.
