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