Line data Source code
1 1 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : /// \file
5 : /// Defines tags related to domain quantities
6 :
7 : #pragma once
8 :
9 : #include <array>
10 : #include <cstddef>
11 : #include <memory>
12 : #include <string>
13 : #include <unordered_map>
14 : #include <unordered_set>
15 : #include <vector>
16 :
17 : #include "DataStructures/DataBox/Subitems.hpp"
18 : #include "DataStructures/DataBox/Tag.hpp"
19 : #include "DataStructures/DataBox/TagName.hpp"
20 : #include "DataStructures/Tensor/EagerMath/Determinant.hpp"
21 : #include "DataStructures/Tensor/EagerMath/DeterminantAndInverse.hpp"
22 : #include "DataStructures/Tensor/TypeAliases.hpp"
23 : #include "DataStructures/Variables.hpp"
24 : #include "Domain/Structure/Direction.hpp"
25 : #include "Domain/Structure/Element.hpp"
26 : #include "Utilities/GetOutput.hpp"
27 : #include "Utilities/Gsl.hpp"
28 : #include "Utilities/NoSuchType.hpp"
29 : #include "Utilities/Requires.hpp"
30 : #include "Utilities/TMPL.hpp"
31 : #include "Utilities/TypeTraits.hpp"
32 : #include "Utilities/TypeTraits/IsA.hpp"
33 :
34 : /// \cond
35 : class DataVector;
36 : template <size_t VolumeDim>
37 : class Domain;
38 : template <size_t VolumeDim, typename Frame>
39 : class ElementMap;
40 : template <size_t VolumeDim>
41 : class Mesh;
42 : /// \endcond
43 :
44 : namespace domain {
45 : /// \ingroup ComputationalDomainGroup
46 : /// \brief %Tags for the domain.
47 : namespace Tags {
48 : /// \ingroup DataBoxTagsGroup
49 : /// \ingroup ComputationalDomainGroup
50 : /// The ::Element associated with the DataBox
51 : template <size_t VolumeDim>
52 1 : struct Element : db::SimpleTag {
53 0 : using type = ::Element<VolumeDim>;
54 : };
55 :
56 : /// \ingroup DataBoxTagsGroup
57 : /// \ingroup ComputationalDomainGroup
58 : /// \brief The computational grid of the Element in the DataBox
59 : /// \details The corresponding interface tag uses Mesh::slice_through to compute
60 : /// the mesh on the face of the element.
61 : template <size_t VolumeDim>
62 1 : struct Mesh : db::SimpleTag {
63 0 : using type = ::Mesh<VolumeDim>;
64 : };
65 :
66 : /// \ingroup DataBoxTagsGroup
67 : /// \ingroup ComputationalDomainGroup
68 : /// The coordinate map from the ElementLogical frame to the TargetFrame
69 : template <size_t VolumeDim, typename TargetFrame = Frame::Inertial>
70 1 : struct ElementMap : db::SimpleTag {
71 0 : static constexpr size_t dim = VolumeDim;
72 0 : using target_frame = TargetFrame;
73 0 : using source_frame = Frame::ElementLogical;
74 :
75 0 : static std::string name() {
76 : return "ElementMap(" + get_output(TargetFrame{}) + ")";
77 : }
78 0 : using type = ::ElementMap<VolumeDim, TargetFrame>;
79 : };
80 :
81 : /// \ingroup DataBoxTagsGroup
82 : /// \ingroup ComputationalDomainGroup
83 : /// The coordinates in a given frame.
84 : template <size_t Dim, typename Frame>
85 1 : struct Coordinates : db::SimpleTag {
86 0 : static std::string name() { return get_output(Frame{}) + "Coordinates"; }
87 0 : using type = tnsr::I<DataVector, Dim, Frame>;
88 : };
89 :
90 : /// \ingroup DataBoxTagsGroup
91 : /// \ingroup ComputationalDomainGroup
92 : /// The coordinates in the target frame of `MapTag`. The `SourceCoordsTag`'s
93 : /// frame must be the source frame of `MapTag`
94 : template <class MapTag, class SourceCoordsTag,
95 : template <size_t, class> class CoordinatesTag = Coordinates>
96 1 : struct MappedCoordinates
97 : : CoordinatesTag<MapTag::dim, typename MapTag::target_frame>,
98 : db::ComputeTag {
99 0 : using base = CoordinatesTag<MapTag::dim, typename MapTag::target_frame>;
100 0 : using return_type = typename base::type;
101 0 : using argument_tags = tmpl::list<MapTag, SourceCoordsTag>;
102 0 : static constexpr auto function(
103 : const gsl::not_null<return_type*> target_coords,
104 : const typename MapTag::type& element_map,
105 : const tnsr::I<DataVector, MapTag::dim, typename MapTag::source_frame>&
106 : source_coords) {
107 : if (source_coords.get(0).size() == 0) {
108 : // When subcell is not supported on an element due to topology, the
109 : // subcell mesh is
110 : // {0, Basis::Uninitialized, Quadrature::Uninitialized}. It is never
111 : // used, but the compute items still need to run
112 : return;
113 : }
114 : *target_coords = element_map(source_coords);
115 : }
116 : };
117 :
118 : /// \ingroup DataBoxTagsGroup
119 : /// \ingroup ComputationalDomainGroup
120 : /// \brief The inverse Jacobian from the source frame to the target frame.
121 : ///
122 : /// Specifically, \f$\partial x^{\bar{i}} / \partial x^i\f$, where \f$\bar{i}\f$
123 : /// denotes the source frame and \f$i\f$ denotes the target frame.
124 : template <size_t Dim, typename SourceFrame, typename TargetFrame>
125 1 : struct InverseJacobian : db::SimpleTag {
126 0 : static std::string name() {
127 : return "InverseJacobian(" + get_output(SourceFrame{}) + "," +
128 : get_output(TargetFrame{}) + ")";
129 : }
130 0 : using type = ::InverseJacobian<DataVector, Dim, SourceFrame, TargetFrame>;
131 : };
132 :
133 : /// \ingroup DataBoxTagsGroup
134 : /// \ingroup ComputationalDomainGroup
135 : /// Computes the inverse Jacobian of the map held by `MapTag` at the coordinates
136 : /// held by `SourceCoordsTag`. The coordinates must be in the source frame of
137 : /// the map.
138 : template <typename MapTag, typename SourceCoordsTag>
139 1 : struct InverseJacobianCompute
140 : : InverseJacobian<MapTag::dim, typename MapTag::source_frame,
141 : typename MapTag::target_frame>,
142 : db::ComputeTag {
143 0 : using base = InverseJacobian<MapTag::dim, typename MapTag::source_frame,
144 : typename MapTag::target_frame>;
145 0 : using return_type = typename base::type;
146 0 : using argument_tags = tmpl::list<MapTag, SourceCoordsTag>;
147 0 : static constexpr auto function(
148 : const gsl::not_null<return_type*> inv_jacobian,
149 : const typename MapTag::type& element_map,
150 : const tnsr::I<DataVector, MapTag::dim, typename MapTag::source_frame>&
151 : source_coords) {
152 : *inv_jacobian = element_map.inv_jacobian(source_coords);
153 : }
154 : };
155 :
156 : /// \ingroup DataBoxTagsGroup
157 : /// \ingroup ComputationalDomainGroup
158 : /// \brief The Jacobian from the source frame to the target frame.
159 : ///
160 : /// Specifically, \f$\partial x^{i} / \partial \xi^{\bar{i}}\f$, where
161 : /// \f$\xi^\bar{i}\f$ denotes the source frame and \f$x^i\f$ denotes the target
162 : /// frame.
163 : template <size_t Dim, typename SourceFrame, typename TargetFrame>
164 1 : struct Jacobian : db::SimpleTag {
165 0 : static std::string name() {
166 : return "Jacobian(" + get_output(SourceFrame{}) + "," +
167 : get_output(TargetFrame{}) + ")";
168 : }
169 0 : using type = ::Jacobian<DataVector, Dim, SourceFrame, TargetFrame>;
170 : };
171 :
172 : /// \ingroup DataBoxTagsGroup
173 : /// \ingroup ComputationalDomainGroup
174 : /// Computes the Jacobian of the map from the `InverseJacobian<Dim, SourceFrame,
175 : /// TargetFrame>` tag.
176 : template <size_t Dim, typename SourceFrame, typename TargetFrame>
177 1 : struct JacobianCompute : Jacobian<Dim, SourceFrame, TargetFrame>,
178 : db::ComputeTag {
179 0 : using base = Jacobian<Dim, SourceFrame, TargetFrame>;
180 0 : using return_type = typename base::type;
181 0 : using argument_tags =
182 : tmpl::list<InverseJacobian<Dim, SourceFrame, TargetFrame>>;
183 0 : static constexpr auto function(
184 : const gsl::not_null<return_type*> jacobian,
185 : const ::InverseJacobian<DataVector, Dim, SourceFrame, TargetFrame>&
186 : inv_jac) {
187 : *jacobian = determinant_and_inverse(inv_jac).second;
188 : }
189 : };
190 :
191 : /// \ingroup DataBoxTagsGroup
192 : /// \ingroup ComputationalDomainGroup
193 : /// \brief The determinant of the inverse Jacobian from the source frame to the
194 : /// target frame.
195 : template <typename SourceFrame, typename TargetFrame>
196 1 : struct DetInvJacobian : db::SimpleTag {
197 0 : using type = Scalar<DataVector>;
198 0 : static std::string name() {
199 : return "DetInvJacobian(" + get_output(SourceFrame{}) + "," +
200 : get_output(TargetFrame{}) + ")";
201 : }
202 : };
203 :
204 : /// \ingroup DataBoxTagsGroup
205 : /// \ingroup ComputationalDomainGroup
206 : /// Computes the determinant of the inverse Jacobian.
207 : template <size_t Dim, typename SourceFrame, typename TargetFrame>
208 1 : struct DetInvJacobianCompute : db::ComputeTag,
209 : DetInvJacobian<SourceFrame, TargetFrame> {
210 0 : using base = DetInvJacobian<SourceFrame, TargetFrame>;
211 0 : using return_type = typename base::type;
212 0 : using argument_tags =
213 : tmpl::list<InverseJacobian<Dim, SourceFrame, TargetFrame>>;
214 0 : static void function(const gsl::not_null<return_type*> det_inv_jac,
215 : const ::InverseJacobian<DataVector, Dim, SourceFrame,
216 : TargetFrame>& inv_jac) {
217 : determinant(det_inv_jac, inv_jac);
218 : }
219 : };
220 :
221 : /// \ingroup DataBoxTagsGroup
222 : /// \ingroup ComputationalDomainGroup
223 : /// \brief The determinant of the Jacobian from the source frame to the target
224 : /// frame.
225 : template <typename SourceFrame, typename TargetFrame>
226 1 : struct DetJacobian : db::SimpleTag {
227 0 : using type = Scalar<DataVector>;
228 0 : static std::string name() {
229 : return "DetJacobian(" + get_output(SourceFrame{}) + "," +
230 : get_output(TargetFrame{}) + ")";
231 : }
232 : };
233 :
234 : /// \ingroup DataBoxTagsGroup
235 : /// \ingroup ComputationalDomainGroup
236 : /// \brief The inverse Jacobian times the determinant of the Jacobian.
237 : ///
238 : /// This quantity is divergence-free analytically. See
239 : /// `::dg::metric_identity_det_jac_times_inv_jac` for more information.
240 : template <size_t Dim, typename SourceFrame, typename TargetFrame>
241 1 : struct DetTimesInvJacobian : db::SimpleTag {
242 0 : using type = ::InverseJacobian<DataVector, Dim, SourceFrame, TargetFrame>;
243 0 : static std::string name() {
244 : return "DetTimesInvJacobian(" + get_output(SourceFrame{}) + "," +
245 : get_output(TargetFrame{}) + ")";
246 : }
247 : };
248 :
249 : /// @{
250 : /// \ingroup DataBoxTagsGroup
251 : /// \ingroup ComputationalDomainGroup
252 : /// The set of directions to neighboring Elements
253 : template <size_t VolumeDim>
254 1 : struct InternalDirections : db::SimpleTag {
255 0 : static constexpr size_t volume_dim = VolumeDim;
256 0 : using type = std::unordered_set<::Direction<VolumeDim>>;
257 : };
258 :
259 : template <size_t VolumeDim>
260 0 : struct InternalDirectionsCompute : InternalDirections<VolumeDim>,
261 : db::ComputeTag {
262 0 : static constexpr size_t volume_dim = VolumeDim;
263 0 : using base = InternalDirections<VolumeDim>;
264 0 : using return_type = std::unordered_set<::Direction<VolumeDim>>;
265 0 : using argument_tags = tmpl::list<Element<VolumeDim>>;
266 0 : static void function(const gsl::not_null<return_type*> directions,
267 : const ::Element<VolumeDim>& element) {
268 : for (const auto& direction_neighbors : element.neighbors()) {
269 : directions->insert(direction_neighbors.first);
270 : }
271 : }
272 : };
273 : /// @}
274 :
275 : /// @{
276 : /// \ingroup DataBoxTagsGroup
277 : /// \ingroup ComputationalDomainGroup
278 : /// The set of directions which correspond to external boundaries.
279 : /// Used for representing data on the interior side of the external boundary
280 : /// faces.
281 : template <size_t VolumeDim>
282 1 : struct BoundaryDirectionsInterior : db::SimpleTag {
283 0 : static constexpr size_t volume_dim = VolumeDim;
284 0 : using type = std::unordered_set<::Direction<VolumeDim>>;
285 : };
286 :
287 : template <size_t VolumeDim>
288 0 : struct BoundaryDirectionsInteriorCompute
289 : : BoundaryDirectionsInterior<VolumeDim>,
290 : db::ComputeTag {
291 0 : static constexpr size_t volume_dim = VolumeDim;
292 0 : using base = BoundaryDirectionsInterior<VolumeDim>;
293 0 : using return_type = std::unordered_set<::Direction<VolumeDim>>;
294 0 : using argument_tags = tmpl::list<Element<VolumeDim>>;
295 0 : static void function(const gsl::not_null<return_type*> directions,
296 : const ::Element<VolumeDim>& element) {
297 : *directions = element.external_boundaries();
298 : }
299 : };
300 : /// @}
301 :
302 : /// @{
303 : /// \ingroup DataBoxTagsGroup
304 : /// \ingroup ComputationalDomainGroup
305 : /// The set of directions which correspond to external boundaries. To be used
306 : /// to represent data which exists on the exterior side of the external boundary
307 : /// faces.
308 : template <size_t VolumeDim>
309 1 : struct BoundaryDirectionsExterior : db::SimpleTag {
310 0 : static constexpr size_t volume_dim = VolumeDim;
311 0 : using type = std::unordered_set<::Direction<VolumeDim>>;
312 : };
313 :
314 : template <size_t VolumeDim>
315 0 : struct BoundaryDirectionsExteriorCompute
316 : : BoundaryDirectionsExterior<VolumeDim>,
317 : db::ComputeTag {
318 0 : static constexpr size_t volume_dim = VolumeDim;
319 0 : using base = BoundaryDirectionsExterior<VolumeDim>;
320 0 : using return_type = std::unordered_set<::Direction<VolumeDim>>;
321 0 : using argument_tags = tmpl::list<Element<VolumeDim>>;
322 0 : static constexpr auto function(const gsl::not_null<return_type*> directions,
323 : const ::Element<VolumeDim>& element) {
324 : *directions = element.external_boundaries();
325 : }
326 : };
327 : /// @}
328 :
329 : /// \ingroup DataBoxTagsGroup
330 : /// \ingroup ComputationalDomainGroup
331 : /// \brief Tag which is either a SimpleTag for quantities on an
332 : /// interface, base tag to a compute item which acts on tags on an interface, or
333 : /// base tag to a compute item which slices a tag from the volume to an
334 : /// interface.
335 : ///
336 : /// The contained object will be a map from ::Direction to the item type of
337 : /// `Tag`, with the set of directions being those produced by `DirectionsTag`.
338 : ///
339 : /// If a SimpleTag is desired on the interface, then this tag can be added to
340 : /// the DataBox directly. If a ComputeTag which acts on Tags on the interface is
341 : /// desired, then the tag should be added using `InterfaceCompute`. If a
342 : /// ComputeTag which slices a TensorTag or a VariablesTag in the volume to an
343 : /// Interface is desired, then it should be added using `Slice`. In all cases,
344 : /// the tag can then be retrieved using `Tags::Interface<DirectionsTag, Tag>`.
345 : ///
346 : /// \tparam DirectionsTag the item of directions
347 : /// \tparam Tag the tag labeling the item
348 : ///
349 : /// \see InterfaceCompute, Slice
350 : template <typename DirectionsTag, typename Tag>
351 1 : struct Interface : db::SimpleTag {
352 : static_assert(db::is_simple_tag_v<DirectionsTag>);
353 : static_assert(db::is_simple_tag_v<Tag>);
354 0 : static std::string name() {
355 : return "Interface<" + db::tag_name<DirectionsTag>() + ", " +
356 : db::tag_name<Tag>() + ">";
357 : };
358 0 : using tag = Tag;
359 0 : using type = std::unordered_map<::Direction<DirectionsTag::volume_dim>,
360 : typename Tag::type>;
361 : };
362 :
363 : /// \ingroup DataBoxTagsGroup
364 : /// \ingroup ComputationalDomainGroup
365 : /// ::Direction to an interface
366 : template <size_t VolumeDim>
367 1 : struct Direction : db::SimpleTag {
368 0 : using type = ::Direction<VolumeDim>;
369 : };
370 :
371 : } // namespace Tags
372 : } // namespace domain
373 :
374 : namespace db {
375 : namespace detail {
376 : template <typename DirectionsTag, typename VariablesTag>
377 : struct InterfaceSubitemsImpl {
378 : using type = tmpl::transform<
379 : typename VariablesTag::type::tags_list,
380 : tmpl::bind<domain::Tags::Interface, tmpl::pin<DirectionsTag>, tmpl::_1>>;
381 :
382 : using tag = domain::Tags::Interface<DirectionsTag, VariablesTag>;
383 :
384 : template <typename Subtag>
385 : static void create_item(
386 : const gsl::not_null<typename tag::type*> parent_value,
387 : const gsl::not_null<typename Subtag::type*> sub_value) {
388 : sub_value->clear();
389 : for (auto& direction_vars : *parent_value) {
390 : const auto& direction = direction_vars.first;
391 : auto& parent_vars = get<typename Subtag::tag>(direction_vars.second);
392 : auto& sub_var = (*sub_value)[direction];
393 : for (auto vars_it = parent_vars.begin(), sub_var_it = sub_var.begin();
394 : vars_it != parent_vars.end(); ++vars_it, ++sub_var_it) {
395 : sub_var_it->set_data_ref(&*vars_it);
396 : }
397 : }
398 : }
399 :
400 : // The `return_type` can be anything for Subitems because the DataBox figures
401 : // out the correct return type, we just use the `return_type` type alias to
402 : // signal to the DataBox we want mutating behavior.
403 : using return_type = NoSuchType;
404 :
405 : template <typename Subtag>
406 : static void create_compute_item(
407 : const gsl::not_null<typename Subtag::type*> sub_value,
408 : const typename tag::type& parent_value) {
409 : for (const auto& direction_vars : parent_value) {
410 : const auto& direction = direction_vars.first;
411 : const auto& parent_vars =
412 : get<typename Subtag::tag>(direction_vars.second);
413 : auto& sub_var = (*sub_value)[direction];
414 : auto sub_var_it = sub_var.begin();
415 : for (auto vars_it = parent_vars.begin();
416 : vars_it != parent_vars.end(); ++vars_it, ++sub_var_it) {
417 : // clang-tidy: do not use const_cast
418 : // The DataBox will only give out a const reference to the
419 : // result of a compute item. Here, that is a reference to a
420 : // const map to Tensors of DataVectors. There is no (publicly
421 : // visible) indirection there, so having the map const will
422 : // allow only allow const access to the contained DataVectors,
423 : // so no modification through the pointer cast here is
424 : // possible.
425 : sub_var_it->set_data_ref(const_cast<DataVector*>(&*vars_it)); // NOLINT
426 : }
427 : }
428 : }
429 : };
430 : } // namespace detail
431 :
432 : template <typename DirectionsTag, typename VariablesTag>
433 0 : struct Subitems<domain::Tags::Interface<DirectionsTag, VariablesTag>,
434 : Requires<tt::is_a_v<Variables, typename VariablesTag::type>>>
435 : : detail::InterfaceSubitemsImpl<DirectionsTag, VariablesTag> {};
436 :
437 : } // namespace db
|
