Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <array>
       7             : #include <cstddef>
       8             : #include <optional>
       9             : #include <unordered_map>
      10             : #include <utility>
      11             : #include <vector>
      12             : 
      13             : #include "Domain/Block.hpp"
      14             : #include "Domain/Creators/Tags/Domain.hpp"
      15             : #include "Domain/Creators/Tags/InitialExtents.hpp"
      16             : #include "Domain/Creators/Tags/InitialRefinementLevels.hpp"
      17             : #include "Domain/Domain.hpp"
      18             : #include "Domain/ElementMap.hpp"
      19             : #include "Domain/ExcisionSphere.hpp"
      20             : #include "Domain/InterfaceLogicalCoordinates.hpp"
      21             : #include "Domain/Structure/CreateInitialMesh.hpp"
      22             : #include "Domain/Structure/ElementId.hpp"
      23             : #include "Domain/Structure/InitialElementIds.hpp"
      24             : #include "Evolution/DiscontinuousGalerkin/Initialization/QuadratureTag.hpp"
      25             : #include "Evolution/Systems/CurvedScalarWave/Worldtube/Tags.hpp"
      26             : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
      27             : #include "NumericalAlgorithms/Spectral/Quadrature.hpp"
      28             : #include "Utilities/Gsl.hpp"
      29             : #include "Utilities/TMPL.hpp"
      30             : 
      31             : namespace CurvedScalarWave::Worldtube::Initialization {
      32             : 
      33             : /*!
      34             :  * \brief Initializes a map of the grid coordinates centered on the worldtube of
      35             :  * all faces that abut the worldtube with corresponding ElementIds.
      36             :  *
      37             :  * \details The worldtube singleton computes an internal solution in the grid
      38             :  * frame and uses this map to evaluate compute it at the grid coordinates of
      39             :  * each element face abutting the worldtube each time step. This data is sent to
      40             :  * the corresponding element where it is used to apply pointwise boundary
      41             :  * conditions. The `ElementFacesGridCoordinates` holds a map of all the element
      42             :  * ids abutting the worldtube with the corresonding grid coordinates.
      43             :  *
      44             :  * \warning This currently assumes that initial domain remains the same and
      45             :  * there is no AMR. To support this, the worldtube could send the
      46             :  * coefficients of its internal solution to each element which can evaluate
      47             :  * it on their current grid in the boundary conditions.
      48             :  *
      49             :  *  DataBox changes:
      50             :  * - Adds: nothing
      51             :  * - Removes: nothing
      52             :  * - Modifies: Tags::ElementFacesCoordinatesMap<Dim>
      53             :  */
      54             : template <size_t Dim>
      55           1 : struct InitializeElementFacesGridCoordinates {
      56           0 :   using return_tags = tmpl::list<Tags::ElementFacesGridCoordinates<Dim>>;
      57           0 :   using simple_tags = tmpl::list<Tags::ElementFacesGridCoordinates<Dim>>;
      58           0 :   using compute_tags = tmpl::list<>;
      59           0 :   using simple_tags_from_options =
      60             :       tmpl::list<::domain::Tags::InitialExtents<Dim>,
      61             :                  ::domain::Tags::InitialRefinementLevels<Dim>,
      62             :                  evolution::dg::Tags::Quadrature>;
      63           0 :   using const_global_cache_tags = tmpl::list<>;
      64           0 :   using mutable_global_cache_tags = tmpl::list<>;
      65           0 :   using argument_tags = tmpl::flatten<
      66             :       tmpl::list<simple_tags_from_options, ::domain::Tags::Domain<Dim>,
      67             :                  Tags::ExcisionSphere<Dim>>>;
      68           0 :   static void apply(
      69             :       const gsl::not_null<std::unordered_map<
      70             :           ElementId<Dim>, tnsr::I<DataVector, Dim, Frame::Grid>>*>
      71             :           element_faces_grid_coords,
      72             :       const std::vector<std::array<size_t, Dim>>& initial_extents,
      73             :       const std::vector<std::array<size_t, Dim>>& initial_refinement,
      74             :       const Spectral::Quadrature& quadrature, const Domain<Dim>& domain,
      75             :       const ::ExcisionSphere<Dim>& excision_sphere) {
      76             :     const auto& blocks = domain.blocks();
      77             :     const auto& worldtube_grid_coords = excision_sphere.center();
      78             :     const auto& neighboring_blocks = excision_sphere.abutting_directions();
      79             :     for (const auto& [block_id, _] : neighboring_blocks) {
      80             :       const auto element_ids =
      81             :           initial_element_ids(block_id, initial_refinement.at(block_id));
      82             :       for (const auto element_id : element_ids) {
      83             :         const auto direction = excision_sphere.abutting_direction(element_id);
      84             :         if (direction.has_value()) {
      85             :           const auto mesh = ::domain::Initialization::create_initial_mesh(
      86             :               initial_extents, element_id, quadrature);
      87             :           const auto face_mesh = mesh.slice_away(direction.value().dimension());
      88             :           const auto& current_block = blocks.at(block_id);
      89             :           const ElementMap<Dim, Frame::Grid> element_map{
      90             :               element_id,
      91             :               current_block.is_time_dependent()
      92             :                   ? current_block.moving_mesh_logical_to_grid_map().get_clone()
      93             :                   : current_block.stationary_map().get_to_grid_frame()};
      94             :           const auto face_logical_coords =
      95             :               interface_logical_coordinates(face_mesh, direction.value());
      96             :           auto faces_grid_coords = element_map(face_logical_coords);
      97             :           for (size_t i = 0; i < 3; ++i) {
      98             :             faces_grid_coords.get(i) -= worldtube_grid_coords.get(i);
      99             :           }
     100             :           element_faces_grid_coords->operator[](element_id) =
     101             :               std::move(faces_grid_coords);
     102             :         }
     103             :       }
     104             :     }
     105             :   }
     106             : };
     107             : }  // namespace CurvedScalarWave::Worldtube::Initialization