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