Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 : #include <memory>
8 : #include <string>
9 : #include <unordered_map>
10 :
11 : #include "DataStructures/DataBox/Protocols/Mutator.hpp"
12 : #include "DataStructures/DataVector.hpp"
13 : #include "DataStructures/Tensor/Tensor.hpp"
14 : #include "DataStructures/Variables.hpp"
15 : #include "Domain/Block.hpp"
16 : #include "Domain/CoordinateMaps/CoordinateMap.hpp"
17 : #include "Domain/CoordinateMaps/Tags.hpp"
18 : #include "Domain/Creators/Tags/Domain.hpp"
19 : #include "Domain/Domain.hpp"
20 : #include "Domain/ElementMap.hpp"
21 : #include "Domain/FunctionsOfTime/FunctionOfTime.hpp"
22 : #include "Domain/FunctionsOfTime/Tags.hpp"
23 : #include "Domain/Structure/Element.hpp"
24 : #include "Domain/Structure/ElementId.hpp"
25 : #include "Domain/Tags.hpp"
26 : #include "Evolution/DgSubcell/Mesh.hpp"
27 : #include "Evolution/DgSubcell/Tags/Coordinates.hpp"
28 : #include "Evolution/DgSubcell/Tags/DidRollback.hpp"
29 : #include "Evolution/DgSubcell/Tags/Inactive.hpp"
30 : #include "Evolution/DgSubcell/Tags/Mesh.hpp"
31 : #include "Evolution/DgSubcell/Tags/OnSubcellFaces.hpp"
32 : #include "Evolution/Initialization/InitialData.hpp"
33 : #include "NumericalAlgorithms/Spectral/Basis.hpp"
34 : #include "NumericalAlgorithms/Spectral/LogicalCoordinates.hpp"
35 : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
36 : #include "NumericalAlgorithms/Spectral/Quadrature.hpp"
37 : #include "PointwiseFunctions/AnalyticData/Tags.hpp"
38 : #include "PointwiseFunctions/InitialDataUtilities/InitialData.hpp"
39 : #include "PointwiseFunctions/InitialDataUtilities/Tags/InitialData.hpp"
40 : #include "Time/Tags/Time.hpp"
41 : #include "Utilities/CallWithDynamicType.hpp"
42 : #include "Utilities/ErrorHandling/Assert.hpp"
43 : #include "Utilities/Gsl.hpp"
44 : #include "Utilities/ProtocolHelpers.hpp"
45 : #include "Utilities/TMPL.hpp"
46 :
47 : namespace evolution::dg::subcell {
48 :
49 : /*!
50 : * \brief Allocate or assign background general relativity quantities on
51 : * cell-centered and face-centered FD grid points, for evolution systems run on
52 : * a curved spacetime without solving Einstein equations (e.g. ValenciaDivclean,
53 : * ForceFree),
54 : *
55 : * \warning This mutator assumes that the GR analytic data or solution
56 : * specifying background spacetime metric is time-independent.
57 : *
58 : */
59 : template <typename System, typename Metavariables, bool ComputeOnlyOnRollback>
60 1 : struct BackgroundGrVars : tt::ConformsTo<db::protocols::Mutator> {
61 0 : static constexpr size_t volume_dim = System::volume_dim;
62 :
63 0 : using gr_vars_tag = typename System::spacetime_variables_tag;
64 0 : using inactive_gr_vars_tag =
65 : evolution::dg::subcell::Tags::Inactive<gr_vars_tag>;
66 0 : using subcell_faces_gr_tag = evolution::dg::subcell::Tags::OnSubcellFaces<
67 : typename System::flux_spacetime_variables_tag, volume_dim>;
68 :
69 0 : using GrVars = typename gr_vars_tag::type;
70 0 : using InactiveGrVars = typename inactive_gr_vars_tag::type;
71 0 : using SubcellFaceGrVars = typename subcell_faces_gr_tag::type;
72 :
73 0 : using argument_tags = tmpl::list<
74 : ::Tags::Time, domain::Tags::FunctionsOfTime,
75 : domain::Tags::Domain<volume_dim>, domain::Tags::Element<volume_dim>,
76 : domain::Tags::ElementMap<volume_dim, Frame::Grid>,
77 : domain::CoordinateMaps::Tags::CoordinateMap<volume_dim, Frame::Grid,
78 : Frame::Inertial>,
79 : evolution::dg::subcell::Tags::Mesh<volume_dim>,
80 : evolution::dg::subcell::Tags::Coordinates<3, Frame::Inertial>,
81 : subcell::Tags::DidRollback, evolution::initial_data::Tags::InitialData>;
82 :
83 0 : using return_tags =
84 : tmpl::list<gr_vars_tag, inactive_gr_vars_tag, subcell_faces_gr_tag>;
85 :
86 : template <typename T>
87 0 : static void apply(
88 : const gsl::not_null<GrVars*> active_gr_vars,
89 : const gsl::not_null<InactiveGrVars*> inactive_gr_vars,
90 : const gsl::not_null<SubcellFaceGrVars*> subcell_face_gr_vars,
91 : const double time,
92 : const std::unordered_map<
93 : std::string,
94 : std::unique_ptr<::domain::FunctionsOfTime::FunctionOfTime>>&
95 : functions_of_time,
96 : const Domain<volume_dim>& domain, const Element<volume_dim>& element,
97 : const ElementMap<volume_dim, Frame::Grid>& logical_to_grid_map,
98 : const domain::CoordinateMapBase<Frame::Grid, Frame::Inertial, volume_dim>&
99 : grid_to_inertial_map,
100 : const Mesh<volume_dim>& subcell_mesh,
101 : const tnsr::I<DataVector, volume_dim, Frame::Inertial>&
102 : subcell_inertial_coords,
103 : const bool did_rollback, const T& solution_or_data) {
104 : const size_t num_subcell_pts = subcell_mesh.number_of_grid_points();
105 :
106 : if (gsl::at(*subcell_face_gr_vars, 0).number_of_grid_points() != 0) {
107 : // Evolution phase
108 :
109 : // Check if the mesh is actually moving i.e. block coordinate map is
110 : // time-dependent. If not, we can skip the evaluation of GR variables
111 : // since they may stay with their values assigned at the initialization
112 : // phase.
113 : const auto& element_id = element.id();
114 : const size_t block_id = element_id.block_id();
115 : const Block<volume_dim>& block = domain.blocks()[block_id];
116 :
117 : if (block.is_time_dependent()) {
118 : if (did_rollback or not ComputeOnlyOnRollback) {
119 : if (did_rollback) {
120 : // Right after rollback, subcell GR vars are stored in the
121 : // `inactive` one.
122 : ASSERT(inactive_gr_vars->number_of_grid_points() == num_subcell_pts,
123 : "The size of subcell GR variables ("
124 : << inactive_gr_vars->number_of_grid_points()
125 : << ") is not equal to the number of FD grid points ("
126 : << subcell_mesh.number_of_grid_points() << ").");
127 :
128 : cell_centered_impl(inactive_gr_vars, time, subcell_inertial_coords,
129 : solution_or_data);
130 :
131 : } else {
132 : // In this case the element didn't rollback but started from FD.
133 : // Therefore subcell GR vars are in the `active` one.
134 : ASSERT(active_gr_vars->number_of_grid_points() == num_subcell_pts,
135 : "The size of subcell GR variables ("
136 : << active_gr_vars->number_of_grid_points()
137 : << ") is not equal to the number of FD grid points ("
138 : << subcell_mesh.number_of_grid_points() << ").");
139 :
140 : cell_centered_impl(active_gr_vars, time, subcell_inertial_coords,
141 : solution_or_data);
142 : }
143 : if constexpr (not std::is_same_v<
144 : typename SubcellFaceGrVars::value_type::tags_list,
145 : tmpl::list<>>) {
146 : face_centered_impl(subcell_face_gr_vars, time, functions_of_time,
147 : logical_to_grid_map, grid_to_inertial_map,
148 : subcell_mesh, solution_or_data);
149 : }
150 : }
151 : }
152 : } else {
153 : // Initialization phase
154 : (*inactive_gr_vars).initialize(num_subcell_pts);
155 :
156 : fd::verify_subcell_mesh(subcell_mesh);
157 : if constexpr (not std::is_same_v<
158 : typename SubcellFaceGrVars::value_type::tags_list,
159 : tmpl::list<>>) {
160 : const size_t num_face_centered_mesh_grid_pts =
161 : (subcell_mesh.extents(0) + 1) * subcell_mesh.extents(1) *
162 : subcell_mesh.extents(2);
163 : for (size_t d = 0; d < volume_dim; ++d) {
164 : gsl::at(*subcell_face_gr_vars, d)
165 : .initialize(num_face_centered_mesh_grid_pts);
166 : }
167 : face_centered_impl(subcell_face_gr_vars, time, functions_of_time,
168 : logical_to_grid_map, grid_to_inertial_map,
169 : subcell_mesh, solution_or_data);
170 : }
171 :
172 : cell_centered_impl(inactive_gr_vars, time, subcell_inertial_coords,
173 : solution_or_data);
174 : }
175 : }
176 :
177 : private:
178 : template <typename Vars, typename T>
179 0 : static void cell_centered_impl(
180 : const gsl::not_null<Vars*> background_gr_vars, const double time,
181 : const tnsr::I<DataVector, volume_dim, Frame::Inertial>& inertial_coords,
182 : const T& solution_or_data) {
183 : GrVars temp{background_gr_vars->data(), background_gr_vars->size()};
184 :
185 : using derived_classes =
186 : tmpl::at<typename Metavariables::factory_creation::factory_classes,
187 : evolution::initial_data::InitialData>;
188 : call_with_dynamic_type<void, derived_classes>(
189 : &solution_or_data, [&temp, &inertial_coords,
190 : &time](const auto* const solution_or_data_ptr) {
191 : temp.assign_subset(evolution::Initialization::initial_data(
192 : *solution_or_data_ptr, inertial_coords, time,
193 : typename GrVars::tags_list{}));
194 : });
195 : }
196 :
197 : template <typename T>
198 0 : static void face_centered_impl(
199 : const gsl::not_null<SubcellFaceGrVars*> face_centered_gr_vars,
200 : const double time,
201 : const std::unordered_map<
202 : std::string,
203 : std::unique_ptr<::domain::FunctionsOfTime::FunctionOfTime>>&
204 : functions_of_time,
205 : const ElementMap<volume_dim, Frame::Grid>& logical_to_grid_map,
206 : const domain::CoordinateMapBase<Frame::Grid, Frame::Inertial, volume_dim>&
207 : grid_to_inertial_map,
208 : const Mesh<volume_dim>& subcell_mesh, const T& solution_or_data) {
209 : const size_t comp_dim = fd::get_computational_dim(subcell_mesh);
210 : fd::verify_subcell_mesh(subcell_mesh);
211 :
212 : for (size_t dim = 0; dim < comp_dim; ++dim) {
213 : const auto basis = subcell_mesh.basis();
214 : auto quadrature = subcell_mesh.quadrature();
215 : auto extents = subcell_mesh.extents().indices();
216 :
217 : gsl::at(extents, dim) = subcell_mesh.extents(0) + 1;
218 : gsl::at(quadrature, dim) = Spectral::Quadrature::FaceCentered;
219 :
220 : const Mesh<volume_dim> face_centered_mesh{extents, basis, quadrature};
221 : const auto face_centered_logical_coords =
222 : logical_coordinates(face_centered_mesh);
223 : const auto face_centered_inertial_coords = grid_to_inertial_map(
224 : logical_to_grid_map(face_centered_logical_coords), time,
225 : functions_of_time);
226 :
227 : using derived_classes =
228 : tmpl::at<typename Metavariables::factory_creation::factory_classes,
229 : evolution::initial_data::InitialData>;
230 : call_with_dynamic_type<void, derived_classes>(
231 : &solution_or_data,
232 : [&face_centered_gr_vars, &face_centered_inertial_coords, &dim,
233 : &time](const auto* const solution_or_data_ptr) {
234 : gsl::at(*face_centered_gr_vars, dim)
235 : .assign_subset(evolution::Initialization::initial_data(
236 : *solution_or_data_ptr, face_centered_inertial_coords, time,
237 : typename SubcellFaceGrVars::value_type::tags_list{}));
238 : });
239 : }
240 : }
241 : };
242 :
243 : } // namespace evolution::dg::subcell
|
