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 <optional>
9 :
10 : #include "DataStructures/DataBox/Tag.hpp"
11 : #include "DataStructures/DataVector.hpp"
12 : #include "DataStructures/Tensor/EagerMath/DotProduct.hpp"
13 : #include "DataStructures/Tensor/EagerMath/Magnitude.hpp"
14 : #include "DataStructures/Tensor/Tensor.hpp"
15 : #include "DataStructures/Variables.hpp"
16 : #include "Domain/CoordinateMaps/CoordinateMap.hpp"
17 : #include "Domain/Structure/Direction.hpp"
18 : #include "Domain/Structure/DirectionMap.hpp"
19 : #include "Evolution/DiscontinuousGalerkin/Actions/ComputeTimeDerivativeHelpers.hpp"
20 : #include "Evolution/DiscontinuousGalerkin/NormalVectorTags.hpp"
21 : #include "Utilities/Gsl.hpp"
22 : #include "Utilities/TMPL.hpp"
23 :
24 : namespace evolution::dg::Actions::detail {
25 : struct OneOverNormalVectorMagnitude : db::SimpleTag {
26 : using type = Scalar<DataVector>;
27 : };
28 :
29 : template <size_t Dim>
30 : struct NormalVector : db::SimpleTag {
31 : using type = tnsr::I<DataVector, Dim, Frame::Inertial>;
32 : };
33 :
34 : /*!
35 : * \brief Computes the normal covector, magnitude of the unnormalized normal
36 : * covector, and in curved spacetimes the normal vector.
37 : *
38 : * The `fields_on_face` argument is used as a temporary buffer via the
39 : * `OneOverNormalVectorMagnitude` tag and also used to return the normalized
40 : * normal vector via the `NormalVector<Dim>` tag. The `fields_on_face` argument
41 : * also serves as the input for the inverse spatial metric when the spacetime is
42 : * curved. The `fields_on_face` argument must have tags `NormalVector<Dim>` and
43 : * `OneOverNormalVectorMagnitude`. If the system specifies
44 : * `System::inverse_spatial_metric_tag` then this tag must also be in the
45 : * Variables.
46 : */
47 : template <typename System, size_t Dim, typename FieldsOnFaceTags>
48 : void unit_normal_vector_and_covector_and_magnitude_impl(
49 : const gsl::not_null<Scalar<DataVector>*> face_normal_magnitude,
50 : const gsl::not_null<tnsr::i<DataVector, Dim, Frame::Inertial>*>
51 : unit_normal_covector,
52 : const gsl::not_null<Variables<FieldsOnFaceTags>*> fields_on_face,
53 : const tnsr::i<DataVector, Dim, Frame::Inertial>&
54 : unnormalized_normal_covector) {
55 : if constexpr (has_inverse_spatial_metric_tag_v<System>) {
56 : using inverse_spatial_metric_tag =
57 : typename System::inverse_spatial_metric_tag;
58 : auto& normal_vector = get<NormalVector<Dim>>(*fields_on_face);
59 : const auto& inverse_spatial_metric =
60 : get<inverse_spatial_metric_tag>(*fields_on_face);
61 : // Compute unnormalized normal vector
62 : for (size_t i = 0; i < Dim; ++i) {
63 : normal_vector.get(i) = inverse_spatial_metric.get(i, 0) *
64 : get<0>(unnormalized_normal_covector);
65 : for (size_t j = 1; j < Dim; ++j) {
66 : normal_vector.get(i) += inverse_spatial_metric.get(i, j) *
67 : unnormalized_normal_covector.get(j);
68 : }
69 : }
70 : // Perform normalization
71 : dot_product(face_normal_magnitude, normal_vector,
72 : unnormalized_normal_covector);
73 : get(*face_normal_magnitude) = sqrt(get(*face_normal_magnitude));
74 : auto& one_over_normal_vector_magnitude =
75 : get<OneOverNormalVectorMagnitude>(*fields_on_face);
76 : get(one_over_normal_vector_magnitude) = 1.0 / get(*face_normal_magnitude);
77 : for (size_t i = 0; i < Dim; ++i) {
78 : unit_normal_covector->get(i) = unnormalized_normal_covector.get(i) *
79 : get(one_over_normal_vector_magnitude);
80 : normal_vector.get(i) *= get(one_over_normal_vector_magnitude);
81 : }
82 : } else {
83 : magnitude(face_normal_magnitude, unnormalized_normal_covector);
84 : auto& one_over_normal_vector_magnitude =
85 : get<OneOverNormalVectorMagnitude>(*fields_on_face);
86 : get(one_over_normal_vector_magnitude) = 1.0 / get(*face_normal_magnitude);
87 : for (size_t i = 0; i < Dim; ++i) {
88 : unit_normal_covector->get(i) = unnormalized_normal_covector.get(i) *
89 : get(one_over_normal_vector_magnitude);
90 : }
91 : }
92 : }
93 :
94 : /*!
95 : * \brief Computes the normal vector, covector, and their magnitude on the face
96 : * in the given direction.
97 : *
98 : * The normal covector and its magnitude are returned using the
99 : * `normal_covector_quantities` function argument. The `fields_on_face` argument
100 : * is used as a temporary buffer via the `OneOverNormalVectorMagnitude` tag and
101 : * also used to return the normalized normal vector via the `NormalVector<Dim>`
102 : * tag. The `fields_on_face` argument also serves as the input for the inverse
103 : * spatial metric when the spacetime is curved. The `fields_on_face` argument
104 : * must have tags `NormalVector<Dim>` and `OneOverNormalVectorMagnitude`. If the
105 : * system specifies `System::inverse_spatial_metric_tag` then this tag must also
106 : * be in the Variables.
107 : */
108 : template <typename System, size_t Dim, typename FieldsOnFaceTags>
109 : void unit_normal_vector_and_covector_and_magnitude(
110 : const gsl::not_null<
111 : DirectionMap<Dim, std::optional<Variables<tmpl::list<
112 : evolution::dg::Tags::MagnitudeOfNormal,
113 : evolution::dg::Tags::NormalCovector<Dim>>>>>*>
114 : normal_covector_quantities,
115 : const gsl::not_null<Variables<FieldsOnFaceTags>*> fields_on_face,
116 : const Direction<Dim>& direction,
117 : const std::unordered_map<Direction<Dim>, tnsr::i<DataVector, Dim>>&
118 : unnormalized_normal_covectors,
119 : const domain::CoordinateMapBase<Frame::Grid, Frame::Inertial, Dim>&
120 : moving_mesh_map) {
121 : const bool mesh_is_moving = not moving_mesh_map.is_identity();
122 : const auto& unnormalized_normal_covector =
123 : unnormalized_normal_covectors.at(direction);
124 :
125 : if (auto& normal_covector_quantity =
126 : normal_covector_quantities->at(direction);
127 : detail::has_inverse_spatial_metric_tag_v<System> or mesh_is_moving or
128 : not normal_covector_quantity.has_value()) {
129 : if (not normal_covector_quantity.has_value()) {
130 : normal_covector_quantity =
131 : Variables<tmpl::list<evolution::dg::Tags::MagnitudeOfNormal,
132 : evolution::dg::Tags::NormalCovector<Dim>>>{
133 : fields_on_face->number_of_grid_points()};
134 : }
135 : detail::unit_normal_vector_and_covector_and_magnitude_impl<System>(
136 : make_not_null(&get<evolution::dg::Tags::MagnitudeOfNormal>(
137 : *normal_covector_quantity)),
138 : make_not_null(&get<evolution::dg::Tags::NormalCovector<Dim>>(
139 : *normal_covector_quantity)),
140 : fields_on_face, unnormalized_normal_covector);
141 : }
142 : }
143 : } // namespace evolution::dg::Actions::detail
|
