SpECTRE Documentation Coverage Report
Current view: top level - ParallelAlgorithms/LinearSolver/Multigrid - Multigrid.hpp Hit Total Coverage
Commit: b5f497991094937944b0a3f519166bb54739d08a Lines: 3 17 17.6 %
Date: 2024-03-28 18:20:13
Legend: Lines: hit not hit

          Line data    Source code
       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <cstddef>
       7             : 
       8             : #include "DataStructures/DataBox/Prefixes.hpp"
       9             : #include "IO/Observer/Actions/RegisterWithObservers.hpp"
      10             : #include "IO/Observer/Helpers.hpp"
      11             : #include "ParallelAlgorithms/Actions/Goto.hpp"
      12             : #include "ParallelAlgorithms/LinearSolver/AsynchronousSolvers/ElementActions.hpp"
      13             : #include "ParallelAlgorithms/LinearSolver/Multigrid/ElementActions.hpp"
      14             : #include "ParallelAlgorithms/LinearSolver/Multigrid/ElementsRegistration.hpp"
      15             : #include "ParallelAlgorithms/LinearSolver/Multigrid/ObserveVolumeData.hpp"
      16             : #include "Utilities/TMPL.hpp"
      17             : 
      18             : /// Items related to the multigrid linear solver
      19             : ///
      20             : /// \see `LinearSolver::multigrid::Multigrid`
      21             : namespace LinearSolver::multigrid {
      22             : 
      23             : /// A label indicating the pre-smoothing step in a V-cycle multigrid algorithm,
      24             : /// i.e. the smoothing step before sending the residual to the coarser (parent)
      25             : /// grid
      26           1 : struct VcycleDownLabel {};
      27             : 
      28             : /// A label indicating the post-smoothing step in a V-cycle multigrid algorithm,
      29             : /// i.e. the smoothing step before sending the correction to the finer (child)
      30             : /// grid
      31           1 : struct VcycleUpLabel {};
      32             : 
      33             : /*!
      34             :  * \brief A V-cycle geometric multgrid solver for linear equations \f$Ax = b\f$
      35             :  *
      36             :  * This linear solver iteratively corrects an initial guess \f$x_0\f$ by
      37             :  * restricting the residual \f$b - Ax\f$ to a series of coarser grids, solving
      38             :  * for a correction on the coarser grids, and then prolongating (interpolating)
      39             :  * the correction back to the finer grids. The solves on grids with different
      40             :  * scales can very effectively solve large-scale modes in the solution, which
      41             :  * Krylov-type linear solvers such as GMRES or Conjugate Gradients typically
      42             :  * struggle with. Therefore, a multigrid solver can be an effective
      43             :  * preconditioner for Krylov-type linear solvers (see
      44             :  * `LinearSolver::gmres::Gmres` and `LinearSolver::cg::ConjugateGradient`). See
      45             :  * \cite Briggs2000jp for an introduction to multigrid methods.
      46             :  *
      47             :  * \par Grid hierarchy
      48             :  * This geometric multigrid solver relies on a strategy to coarsen the
      49             :  * computational grid in a way that removes small-scale modes. We currently
      50             :  * h-coarsen the domain, meaning that we create multigrid levels by successively
      51             :  * combining two elements into one along every dimension of the grid. We only
      52             :  * p-coarsen the grid in the sense that we choose the smaller of the two
      53             :  * polynomial degrees when combining elements. This strategy follows
      54             :  * \cite Vincent2019qpd. See `LinearSolver::multigrid::ElementsAllocator` and
      55             :  * `LinearSolver::multigrid::coarsen` for the code that creates the multigrid
      56             :  * hierarchy.
      57             :  *
      58             :  * \par Inter-mesh operators
      59             :  * The algorithm relies on operations that project data between grids. Residuals
      60             :  * are projected from finer to coarser grids ("restriction") and solutions are
      61             :  * projected from coarser to finer grids ("prolongation"). We use the standard
      62             :  * \f$L_2\f$-projections implemented in
      63             :  * `Spectral::projection_matrix_child_to_parent` and
      64             :  * `Spectral::projection_matrix_parent_to_child`, where "child" means the finer
      65             :  * grid and "parent" means the coarser grid. Note that the residuals \f$b -
      66             :  * Ax\f$ may or may not already include a mass matrix and Jacobian factors,
      67             :  * depending on the implementation of the linear operator \f$A\f$. To account
      68             :  * for this, provide a tag as the `ResidualIsMassiveTag` that holds a `bool`.
      69             :  * See section 3 in \cite Fortunato2019jl for details on the inter-mesh
      70             :  * operators.
      71             :  *
      72             :  * \par Smoother
      73             :  * On every level of the grid hierarchy the multigrid solver relies on a
      74             :  * "smoother" that performs an approximate linear solve on that level. The
      75             :  * smoother can be any linear solver that solves the `smooth_fields_tag` for the
      76             :  * source in the `smooth_source_tag`. Useful smoothers are asynchronous linear
      77             :  * solvers that parallelize well, such as `LinearSolver::Schwarz::Schwarz`. The
      78             :  * multigrid algorithm doesn't assume anything about the smoother, but requires
      79             :  * only that the `PreSmootherActions` and the `PostSmootherActions` passed to
      80             :  * `solve` leave the `smooth_fields_tag` in a state that represents an
      81             :  * approximate solution to the `smooth_source_tag`, and that
      82             :  * `db::add_tag_prefix<LinearSolver::Tags::OperatorAppliedTo,
      83             :  * smooth_fields_tag>` is left up-to-date as well.
      84             :  * The smoother can assume that, on entry, these two tags represent the initial
      85             :  * fields and the linear operator applied to the initial fields, respectively.
      86             :  * Here's an example of setting up a smoother for the multigrid solver:
      87             :  *
      88             :  * \snippet Test_MultigridAlgorithm.cpp setup_smoother
      89             :  *
      90             :  * The smoother can be used to construct an action list like this:
      91             :  *
      92             :  * \snippet Test_MultigridAlgorithm.cpp action_list
      93             :  *
      94             :  * \par Algorithm overview
      95             :  * Every iteration of the multigrid algorithm performs a V-cycle over the grid
      96             :  * hierarchy. One V-cycle consists of first "going down" the grid hierarchy,
      97             :  * from the finest to successively coarser grids, smoothing on every level
      98             :  * ("pre-smoothing", can be controlled by the option
      99             :  * `LinearSolver::multigrid::Tags::EnablePreSmoothing`),
     100             :  * and then "going up" the grid hierarchy again, smoothing on
     101             :  * every level again ("post-smoothing"). When going down, the algorithm projects
     102             :  * the remaining residual of the smoother to the next-coarser grid, setting it
     103             :  * as the source for the smoother on the coarser grid. When going up again, the
     104             :  * algorithm projects the solution of the smoother to the next-finer grid,
     105             :  * adding it to the solution on the finer grid as a correction. The bottom-most
     106             :  * coarsest grid (the "tip" of the V-cycle) may skip the post-smoothing, so the
     107             :  * result of the pre-smoother is immediately projected up to the finer grid
     108             :  * (controlled by the
     109             :  * `LinearSolver::multigrid::Tags::EnablePostSmoothingAtBottom` option). On the
     110             :  * top-most finest grid (the "original" grid that represents the overall
     111             :  * solution) the algorithm applies the smoothing and the corrections from the
     112             :  * coarser grids directly to the solution fields.
     113             :  *
     114             :  * \par AMR
     115             :  * AMR is not yet fully supported by the multigrid solver. When AMR is enabled,
     116             :  * only a single multigrid level can be used (the finest grid). To support AMR
     117             :  * fully, we have to send AMR decisions to coarser grids to refine those as
     118             :  * well.
     119             :  */
     120             : template <typename Metavariables, typename FieldsTag, typename OptionsGroup,
     121             :           typename ResidualIsMassiveTag,
     122             :           typename SourceTag =
     123             :               db::add_tag_prefix<::Tags::FixedSource, FieldsTag>>
     124           1 : struct Multigrid {
     125           0 :   static constexpr size_t Dim = Metavariables::volume_dim;
     126           0 :   using fields_tag = FieldsTag;
     127           0 :   using options_group = OptionsGroup;
     128           0 :   using source_tag = SourceTag;
     129             : 
     130           0 :   using operand_tag = FieldsTag;
     131             : 
     132           0 :   using smooth_source_tag = source_tag;
     133           0 :   using smooth_fields_tag = fields_tag;
     134             : 
     135           0 :   using component_list = tmpl::list<
     136             :       detail::ElementsRegistrationComponent<Metavariables, OptionsGroup>>;
     137             : 
     138           0 :   using observed_reduction_data_tags = observers::make_reduction_data_tags<
     139             :       tmpl::list<async_solvers::reduction_data>>;
     140             : 
     141           0 :   using initialize_element = tmpl::list<
     142             :       async_solvers::InitializeElement<FieldsTag, OptionsGroup, SourceTag>,
     143             :       detail::InitializeElement<Dim, FieldsTag, OptionsGroup, SourceTag>>;
     144             : 
     145           0 :   using amr_projectors =
     146             :       tmpl::list<initialize_element,
     147             :                  detail::ProjectMultigridSections<Metavariables>>;
     148             : 
     149           0 :   using register_element =
     150             :       tmpl::list<detail::RegisterElement<Dim, OptionsGroup>,
     151             :                  async_solvers::RegisterElement<FieldsTag, OptionsGroup,
     152             :                                                 SourceTag, Tags::IsFinestGrid>,
     153             :                  observers::Actions::RegisterWithObservers<
     154             :                      detail::RegisterWithVolumeObserver<OptionsGroup>>>;
     155             : 
     156             :   template <typename ApplyOperatorActions, typename PreSmootherActions,
     157             :             typename PostSmootherActions, typename Label = OptionsGroup>
     158           0 :   using solve = tmpl::list<
     159             :       async_solvers::PrepareSolve<FieldsTag, OptionsGroup, SourceTag, Label,
     160             :                                   Tags::IsFinestGrid, false>,
     161             :       detail::ReceiveResidualFromFinerGrid<Dim, FieldsTag, OptionsGroup,
     162             :                                            SourceTag>,
     163             :       detail::PreparePreSmoothing<FieldsTag, OptionsGroup, SourceTag>,
     164             :       // No need to apply the linear operator here:
     165             :       // - On the finest grid, the operator applied to the fields should have
     166             :       //   already been computed at this point, either applied to the initial
     167             :       //   fields before the multigrid solver is invoked, or by the smoother at
     168             :       //   the end of the previous V-cycle.
     169             :       // - On coarser grids, the initial fields are zero, so the operator
     170             :       //   applied to them is also zero.
     171             :       PreSmootherActions,
     172             :       detail::SkipPostSmoothingAtBottom<FieldsTag, OptionsGroup, SourceTag>,
     173             :       detail::SendResidualToCoarserGrid<FieldsTag, OptionsGroup,
     174             :                                         ResidualIsMassiveTag, SourceTag>,
     175             :       detail::ReceiveCorrectionFromCoarserGrid<Dim, FieldsTag, OptionsGroup,
     176             :                                                SourceTag>,
     177             :       ApplyOperatorActions, ::Actions::Label<detail::PostSmoothingBeginLabel>,
     178             :       PostSmootherActions,
     179             :       detail::SendCorrectionToFinerGrid<FieldsTag, OptionsGroup, SourceTag>,
     180             :       detail::ObserveVolumeData<FieldsTag, OptionsGroup, SourceTag>,
     181             :       async_solvers::CompleteStep<FieldsTag, OptionsGroup, SourceTag, Label,
     182             :                                   Tags::IsFinestGrid, false>>;
     183             : };
     184             : 
     185             : }  // namespace LinearSolver::multigrid

Generated by: LCOV version 1.14