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 <memory> 9 : #include <string> 10 : #include <unordered_map> 11 : #include <utility> 12 : 13 : #include "DataStructures/DataBox/Protocols/Mutator.hpp" 14 : #include "Domain/CoordinateMaps/CoordinateMap.hpp" 15 : #include "Domain/CoordinateMaps/Tags.hpp" 16 : #include "Domain/ElementMap.hpp" 17 : #include "Domain/FunctionsOfTime/FunctionOfTime.hpp" 18 : #include "Domain/FunctionsOfTime/Tags.hpp" 19 : #include "Evolution/DgSubcell/Tags/Coordinates.hpp" 20 : #include "Evolution/DgSubcell/Tags/Inactive.hpp" 21 : #include "Evolution/DgSubcell/Tags/Mesh.hpp" 22 : #include "Evolution/DgSubcell/Tags/OnSubcellFaces.hpp" 23 : #include "Evolution/Initialization/Tags.hpp" 24 : #include "Evolution/Systems/ScalarAdvection/Tags.hpp" 25 : #include "Evolution/Systems/ScalarAdvection/VelocityField.hpp" 26 : #include "NumericalAlgorithms/Spectral/Mesh.hpp" 27 : #include "NumericalAlgorithms/Spectral/Quadrature.hpp" 28 : #include "Time/Tags/Time.hpp" 29 : #include "Utilities/ErrorHandling/Assert.hpp" 30 : #include "Utilities/Gsl.hpp" 31 : #include "Utilities/ProtocolHelpers.hpp" 32 : #include "Utilities/TMPL.hpp" 33 : 34 : /// \cond 35 : namespace domain { 36 : namespace Tags { 37 : template <size_t Dim, typename Frame> 38 : struct Coordinates; 39 : template <size_t VolumeDim> 40 : struct Mesh; 41 : } // namespace Tags 42 : } // namespace domain 43 : /// \endcond 44 : 45 : namespace ScalarAdvection::subcell { 46 : /*! 47 : * \brief Allocate and set the velocity field needed for evolving 48 : * ScalarAdvection system when using a DG-subcell hybrid scheme. 49 : * 50 : * Uses: 51 : * - DataBox: 52 : * * `Tags::Time` 53 : * * `evolution::dg::subcell::Tags::Mesh<Dim>` 54 : * * `domain::Tags::ElementMap<Dim, Frame::Grid>` 55 : * * `Tags::CoordinateMap<Dim, Frame::Grid, Frame::Inertial>` 56 : * * `domain::Tags::FunctionsOfTime` 57 : * * `evolution::dg::subcell::Tags::Coordinates<Dim, Frame::ElementLogical>` 58 : * 59 : * DataBox changes: 60 : * - Adds: 61 : * * `evolution::dg::subcell::Tags::Inactive<velocity_field>` 62 : * * `evolution::dg::subcell::Tags::OnSubcellFaces<velocity_field, Dim>` 63 : * where `velocity_field` is `ScalarAdvection::Tags::VelocityField<Dim>`. 64 : * 65 : * - Removes: nothing 66 : * - Modifies: nothing 67 : * 68 : * \note This mutator is meant to be used with 69 : * `Initialization::Actions::AddSimpleTags`. 70 : */ 71 : template <size_t Dim> 72 1 : struct VelocityAtFace : tt::ConformsTo<db::protocols::Mutator> { 73 0 : using simple_tags_from_options = tmpl::list<::Tags::Time>; 74 : 75 0 : using velocity_field = ::ScalarAdvection::Tags::VelocityField<Dim>; 76 0 : using subcell_velocity_field = 77 : ::evolution::dg::subcell::Tags::Inactive<velocity_field>; 78 0 : using subcell_faces_velocity_field = 79 : ::evolution::dg::subcell::Tags::OnSubcellFaces<velocity_field, Dim>; 80 : 81 0 : using vars = typename velocity_field::type; 82 0 : using subcell_vars = typename subcell_velocity_field::type; 83 0 : using face_vars = typename subcell_faces_velocity_field::type::value_type; 84 : 85 0 : using return_tags = 86 : tmpl::list<subcell_velocity_field, subcell_faces_velocity_field>; 87 : 88 0 : using argument_tags = tmpl::list< 89 : ::Tags::Time, 90 : evolution::dg::subcell::Tags::Mesh<Dim>, 91 : domain::Tags::ElementMap<Dim, Frame::Grid>, 92 : domain::CoordinateMaps::Tags::CoordinateMap<Dim, Frame::Grid, 93 : Frame::Inertial>, 94 : domain::Tags::FunctionsOfTime, 95 : evolution::dg::subcell::Tags::Coordinates<Dim, Frame::ElementLogical>>; 96 : 97 0 : static void apply( 98 : const gsl::not_null<subcell_vars*> cell_centered_vars, 99 : const gsl::not_null<std::array<face_vars, Dim>*> face_centered_vars, 100 : const double initial_time, const Mesh<Dim>& subcell_mesh, 101 : const ElementMap<Dim, Frame::Grid>& logical_to_grid_map, 102 : const domain::CoordinateMapBase<Frame::Grid, Frame::Inertial, Dim>& 103 : grid_to_inertial_map, 104 : const std::unordered_map< 105 : std::string, 106 : std::unique_ptr<domain::FunctionsOfTime::FunctionOfTime>>& 107 : functions_of_time, 108 : const tnsr::I<DataVector, Dim, Frame::ElementLogical>& 109 : subcell_logical_coordinates) { 110 : // check if the subcell mesh is isotropic 111 : ASSERT(Mesh<Dim>(subcell_mesh.extents(0), subcell_mesh.basis(0), 112 : subcell_mesh.quadrature(0)) == subcell_mesh, 113 : "The subcell mesh must have isotropic basis, quadrature. and " 114 : "extents but got " 115 : << subcell_mesh); 116 : 117 : const size_t num_grid_points = subcell_mesh.number_of_grid_points(); 118 : const auto cell_centered_inertial_coords = 119 : grid_to_inertial_map(logical_to_grid_map(subcell_logical_coordinates), 120 : initial_time, functions_of_time); 121 : 122 : // Set cell-centered vars. Need to first do without prefix then move into 123 : // prefixed Variables. 124 : vars no_prefix_cell_centered_vars{num_grid_points}; 125 : // Note: This assumes the velocity field is the hard-coded one in the 126 : // compute tag 127 : ::ScalarAdvection::Tags::VelocityFieldCompute<Dim>::function( 128 : make_not_null(&no_prefix_cell_centered_vars), 129 : cell_centered_inertial_coords); 130 : *cell_centered_vars = std::move(no_prefix_cell_centered_vars); 131 : 132 : // Set face-centered vars. 133 : for (size_t dim = 0; dim < Dim; ++dim) { 134 : const auto basis = make_array<Dim>(subcell_mesh.basis(0)); 135 : auto quadrature = make_array<Dim>(subcell_mesh.quadrature(0)); 136 : auto extents = make_array<Dim>(subcell_mesh.extents(0)); 137 : gsl::at(extents, dim) = subcell_mesh.extents(0) + 1; 138 : gsl::at(quadrature, dim) = Spectral::Quadrature::FaceCentered; 139 : const Mesh<Dim> face_centered_mesh{extents, basis, quadrature}; 140 : const auto face_centered_logical_coords = 141 : logical_coordinates(face_centered_mesh); 142 : const auto face_centered_inertial_coords = grid_to_inertial_map( 143 : logical_to_grid_map(face_centered_logical_coords), initial_time, 144 : functions_of_time); 145 : 146 : ::ScalarAdvection::Tags::VelocityFieldCompute<Dim>::function( 147 : make_not_null(&(gsl::at(*face_centered_vars, dim))), 148 : face_centered_inertial_coords); 149 : } 150 : } 151 : }; 152 : } // namespace ScalarAdvection::subcell