SpECTRE Documentation Coverage Report
Current view: top level - ParallelAlgorithms/LinearSolver/Multigrid - ElementsAllocator.hpp Hit Total Coverage
Commit: 9a905b0737f373631c1b8e8389b8f26e67fa5313 Lines: 1 4 25.0 %
Date: 2024-03-28 09:03:18
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 <algorithm>
       7             : #include <charm++.h>
       8             : #include <cstddef>
       9             : #include <optional>
      10             : #include <unordered_map>
      11             : #include <unordered_set>
      12             : #include <vector>
      13             : 
      14             : #include "Domain/Creators/Tags/Domain.hpp"
      15             : #include "Domain/Creators/Tags/InitialExtents.hpp"
      16             : #include "Domain/Creators/Tags/InitialRefinementLevels.hpp"
      17             : #include "Domain/ElementDistribution.hpp"
      18             : #include "Domain/Structure/ElementId.hpp"
      19             : #include "Domain/Tags/ElementDistribution.hpp"
      20             : #include "Elliptic/DiscontinuousGalerkin/Tags.hpp"
      21             : #include "NumericalAlgorithms/Convergence/Tags.hpp"
      22             : #include "Parallel/GlobalCache.hpp"
      23             : #include "Parallel/Local.hpp"
      24             : #include "Parallel/Printf.hpp"
      25             : #include "Parallel/Protocols/ArrayElementsAllocator.hpp"
      26             : #include "Parallel/Section.hpp"
      27             : #include "Parallel/Tags/Section.hpp"
      28             : #include "ParallelAlgorithms/LinearSolver/Multigrid/Hierarchy.hpp"
      29             : #include "ParallelAlgorithms/LinearSolver/Multigrid/Tags.hpp"
      30             : #include "Utilities/ErrorHandling/Error.hpp"
      31             : #include "Utilities/PrettyType.hpp"
      32             : #include "Utilities/ProtocolHelpers.hpp"
      33             : #include "Utilities/Serialization/Serialize.hpp"
      34             : #include "Utilities/System/ParallelInfo.hpp"
      35             : #include "Utilities/TaggedTuple.hpp"
      36             : 
      37             : namespace LinearSolver::multigrid {
      38             : 
      39             : /*!
      40             :  * \brief A `Parallel::protocols::ArrayElementsAllocator` that creates array
      41             :  * elements to cover the initial computational domain multiple times at
      42             :  * different refinement levels, suitable for the
      43             :  * `LinearSolver::multigrid::Multigrid` algorithm
      44             :  *
      45             :  * The `initial_element_ids` function and the option-created
      46             :  * `domain::Tags::Domain` and `domain::Tags::InitialRefinementLevels` determine
      47             :  * the initial set of element IDs in the domain. This is taken as the finest
      48             :  * grid in the multigrid hierarchy. Coarser grids are determined by successively
      49             :  * applying `LinearSolver::multigrid::coarsen`, up to
      50             :  * `LinearSolver::multigrid::Tags::MaxLevels` grids.
      51             :  *
      52             :  * Array elements are created for all element IDs on all grids, meaning they all
      53             :  * share the same parallel component, action list etc. Elements are connected to
      54             :  * their neighbors _on the same grid_ by the
      55             :  * `domain::Initialization::create_initial_element` function (this is
      56             :  * independent of the multigrid code, the function is typically called in an
      57             :  * initialization action). Elements are connected to their parent and children
      58             :  * _across grids_ by the multigrid-tags set here and in the multigrid
      59             :  * initialization actions (see `LinearSolver::multigrid::parent_id` and
      60             :  * `LinearSolver::multigrid::child_ids`).
      61             :  *
      62             :  * This allocator also creates two sets of sections (see `Parallel::Section`):
      63             :  * - `Parallel::Tags::Section<ElementArray,
      64             :  *   LinearSolver::multigrid::Tags::MultigridLevel>`: One section per grid. All
      65             :  *   elements are part of a section with this tag.
      66             :  * - `Parallel::Tags::Section<ElementArray,
      67             :  *   LinearSolver::multigrid::Tags::IsFinestGrid>`: A single section that holds
      68             :  *   only the elements on the finest grid. Holds `std::nullopt` on all other
      69             :  *   elements.
      70             :  *
      71             :  * The elements are distributed on processors using the
      72             :  * `domain::BlockZCurveProcDistribution` for every grid independently. An
      73             :  * unordered set of `size_t`s can be passed to the `apply` function which
      74             :  * represents physical processors to avoid placing elements on.
      75             :  */
      76             : template <size_t Dim, typename OptionsGroup>
      77           1 : struct ElementsAllocator
      78             :     : tt::ConformsTo<Parallel::protocols::ArrayElementsAllocator> {
      79             :   template <typename ElementArray>
      80           0 :   using array_allocation_tags =
      81             :       tmpl::list<Parallel::Tags::Section<ElementArray, Tags::MultigridLevel>,
      82             :                  Parallel::Tags::Section<ElementArray, Tags::IsFinestGrid>>;
      83             : 
      84             :   template <typename ElementArray, typename Metavariables,
      85             :             typename... InitializationTags>
      86           0 :   static void apply(Parallel::CProxy_GlobalCache<Metavariables>& global_cache,
      87             :                     const tuples::TaggedTuple<InitializationTags...>&
      88             :                         original_initialization_items,
      89             :                     const tuples::tagged_tuple_from_typelist<
      90             :                         typename ElementArray::
      91             :                             array_allocation_tags>& /*array_allocation_items*/
      92             :                     = {},
      93             :                     const std::unordered_set<size_t>& procs_to_ignore = {}) {
      94             :     // Copy the initialization items so we can adjust them on each refinement
      95             :     // level
      96             :     auto initialization_items =
      97             :         deserialize<tuples::TaggedTuple<InitializationTags...>>(
      98             :             serialize<tuples::TaggedTuple<InitializationTags...>>(
      99             :                 original_initialization_items)
     100             :                 .data());
     101             :     auto& local_cache = *Parallel::local_branch(global_cache);
     102             :     auto& element_array =
     103             :         Parallel::get_parallel_component<ElementArray>(local_cache);
     104             :     const auto& domain = get<domain::Tags::Domain<Dim>>(local_cache);
     105             :     const auto& initial_extents =
     106             :         get<domain::Tags::InitialExtents<Dim>>(initialization_items);
     107             :     auto& initial_refinement_levels =
     108             :         get<domain::Tags::InitialRefinementLevels<Dim>>(initialization_items);
     109             :     auto& children_refinement_levels =
     110             :         get<Tags::ChildrenRefinementLevels<Dim>>(initialization_items);
     111             :     auto& parent_refinement_levels =
     112             :         get<Tags::ParentRefinementLevels<Dim>>(initialization_items);
     113             :     const auto& quadrature =
     114             :         Parallel::get<elliptic::dg::Tags::Quadrature>(local_cache);
     115             :     const std::optional<domain::ElementWeight>& element_weight =
     116             :         get<domain::Tags::ElementDistribution>(local_cache);
     117             :     std::optional<size_t> max_levels =
     118             :         get<Tags::MaxLevels<OptionsGroup>>(local_cache);
     119             :     const size_t number_of_procs =
     120             :         Parallel::number_of_procs<size_t>(local_cache);
     121             :     if (max_levels == 0) {
     122             :       ERROR_NO_TRACE(
     123             :           "The 'MaxLevels' option includes the finest grid, so '0' is not a "
     124             :           "valid value. Set the option to '1' to effectively disable "
     125             :           "multigrid.");
     126             :     }
     127             :     const size_t num_iterations =
     128             :         get<Convergence::Tags::Iterations<OptionsGroup>>(local_cache);
     129             :     if (UNLIKELY(num_iterations == 0)) {
     130             :       Parallel::printf(
     131             :           "%s is disabled (zero iterations). Only creating a single grid.\n",
     132             :           pretty_type::name<OptionsGroup>());
     133             :       max_levels = 1;
     134             :     }
     135             :     const auto& blocks = domain.blocks();
     136             :     size_t multigrid_level = 0;
     137             :     do {
     138             :       // Store the current grid as child grid before coarsening it
     139             :       children_refinement_levels = initial_refinement_levels;
     140             :       initial_refinement_levels = parent_refinement_levels;
     141             :       // Construct coarsened (parent) grid
     142             :       if (not max_levels.has_value() or multigrid_level < *max_levels - 1) {
     143             :         parent_refinement_levels =
     144             :             LinearSolver::multigrid::coarsen(initial_refinement_levels);
     145             :       }
     146             :       // Create element IDs for all elements on this level
     147             :       std::vector<ElementId<Dim>> element_ids{};
     148             :       for (const auto& block : blocks) {
     149             :         const std::vector<ElementId<Dim>> block_element_ids =
     150             :             initial_element_ids(block.id(),
     151             :                                 initial_refinement_levels[block.id()],
     152             :                                 multigrid_level);
     153             :         element_ids.insert(element_ids.begin(), block_element_ids.begin(),
     154             :                            block_element_ids.end());
     155             :       }
     156             :       // Create an array section for this refinement level
     157             :       std::vector<CkArrayIndex> array_indices(element_ids.size());
     158             :       std::transform(
     159             :           element_ids.begin(), element_ids.end(), array_indices.begin(),
     160             :           [](const ElementId<Dim>& local_element_id) {
     161             :             return Parallel::ArrayIndex<ElementId<Dim>>(local_element_id);
     162             :           });
     163             :       using MultigridLevelSection =
     164             :           Parallel::Section<ElementArray, Tags::MultigridLevel>;
     165             :       get<Parallel::Tags::Section<ElementArray, Tags::MultigridLevel>>(
     166             :           initialization_items) = MultigridLevelSection{
     167             :           multigrid_level, MultigridLevelSection::cproxy_section::ckNew(
     168             :                                element_array.ckGetArrayID(),
     169             :                                array_indices.data(), array_indices.size())};
     170             :       using FinestGridSection =
     171             :           Parallel::Section<ElementArray, Tags::IsFinestGrid>;
     172             :       get<Parallel::Tags::Section<ElementArray, Tags::IsFinestGrid>>(
     173             :           initialization_items) =
     174             :           multigrid_level == 0
     175             :               ? std::make_optional(FinestGridSection{
     176             :                     true, FinestGridSection::cproxy_section::ckNew(
     177             :                               element_array.ckGetArrayID(),
     178             :                               array_indices.data(), array_indices.size())})
     179             :               : std::nullopt;
     180             :       // Create the elements for this refinement level and distribute them among
     181             :       // processors
     182             :       const size_t num_of_procs_to_use =
     183             :           static_cast<size_t>(sys::number_of_procs()) - procs_to_ignore.size();
     184             :       // Distributed with weighted space filling curve
     185             :       if (element_weight.has_value()) {
     186             :         const std::unordered_map<ElementId<Dim>, double> element_costs =
     187             :             domain::get_element_costs(blocks, initial_refinement_levels,
     188             :                                       initial_extents, element_weight.value(),
     189             :                                       quadrature);
     190             :         const domain::BlockZCurveProcDistribution<Dim> element_distribution{
     191             :             element_costs,   num_of_procs_to_use,
     192             :             blocks,          initial_refinement_levels,
     193             :             initial_extents, procs_to_ignore};
     194             : 
     195             :         for (const auto& element_id : element_ids) {
     196             :           const size_t target_proc =
     197             :               element_distribution.get_proc_for_element(element_id);
     198             :           element_array(element_id)
     199             :               .insert(global_cache, initialization_items, target_proc);
     200             :         }
     201             :       } else {
     202             :         // Distributed with round-robin
     203             :         size_t which_proc = 0;
     204             :         for (const auto& element_id : element_ids) {
     205             :           while (procs_to_ignore.find(which_proc) != procs_to_ignore.end()) {
     206             :             which_proc = which_proc + 1 == number_of_procs ? 0 : which_proc + 1;
     207             :           }
     208             : 
     209             :           element_array(element_id)
     210             :               .insert(global_cache, initialization_items, which_proc);
     211             : 
     212             :           which_proc = which_proc + 1 == number_of_procs ? 0 : which_proc + 1;
     213             :         }
     214             :       }
     215             :       Parallel::printf(
     216             :           "%s level %zu has %zu elements in %zu blocks distributed on %d "
     217             :           "procs.\n",
     218             :           pretty_type::name<OptionsGroup>(), multigrid_level,
     219             :           element_ids.size(), blocks.size(), num_of_procs_to_use);
     220             :       ++multigrid_level;
     221             :     } while (initial_refinement_levels != parent_refinement_levels);
     222             :     element_array.doneInserting();
     223             :   }
     224             : };
     225             : 
     226             : }  // namespace LinearSolver::multigrid

Generated by: LCOV version 1.14