NormalVectors.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
6 #include <cstddef>
7 
8 #include "DataStructures/DataVector.hpp"
11 #include "Utilities/Gsl.hpp"
12 #include "Utilities/TMPL.hpp"
13 #include "Utilities/TypeTraits/CreateHasTypeAlias.hpp"
14 
15 namespace TestHelpers::evolution::dg::detail {
16 CREATE_HAS_TYPE_ALIAS(inverse_spatial_metric_tag)
17 CREATE_HAS_TYPE_ALIAS_V(inverse_spatial_metric_tag)
18 
19 template <bool HasInverseSpatialMetricTag = false>
20 struct inverse_spatial_metric_tag {
21  template <typename System>
22  using f = tmpl::list<>;
23 };
24 
25 template <>
26 struct inverse_spatial_metric_tag<true> {
27  template <typename System>
28  using f = typename System::inverse_spatial_metric_tag;
29 };
30 
31 // On input `inv_spatial_metric` is expected to have components on the interval
32 // [0, 1]. The components are rescaled by 0.01, and 1 is added to the diagonal.
33 // This is to give an inverse spatial metric of the form:
34 // \delta^{ij} + small^{ij}
35 // This is done to give a physically reasonable inverse spatial metric
36 template <size_t Dim>
37 void adjust_inverse_spatial_metric(
38  const gsl::not_null<tnsr::II<DataVector, Dim>*> inv_spatial_metric) {
39  for (size_t i = 0; i < Dim; ++i) {
40  for (size_t j = i; j < Dim; ++j) {
41  inv_spatial_metric->get(i, j) *= 0.01;
42  }
43  }
44  for (size_t i = 0; i < Dim; ++i) {
45  inv_spatial_metric->get(i, i) += 1.0;
46  }
47 }
48 
49 // On input the `unit_normal_covector` is the unnormalized normal covector. On
50 // output `unit_normal_covector` is the normalized (and hence actually unit)
51 // normal covector, and `unit_normal_vector` is the unit normal vector. The
52 // inverse spatial metric is used for computing the magnitude of the
53 // unnormalized normal vector.
54 template <size_t Dim>
55 void normalize_vector_and_covector(
56  const gsl::not_null<tnsr::i<DataVector, Dim>*> unit_normal_covector,
57  const gsl::not_null<tnsr::I<DataVector, Dim>*> unit_normal_vector,
58  const tnsr::II<DataVector, Dim>& inv_spatial_metric) {
59  for (size_t i = 0; i < Dim; ++i) {
60  unit_normal_vector->get(i) =
61  inv_spatial_metric.get(i, 0) * get<0>(*unit_normal_covector);
62  for (size_t j = 1; j < Dim; ++j) {
63  unit_normal_vector->get(i) +=
64  inv_spatial_metric.get(i, j) * unit_normal_covector->get(j);
65  }
66  }
67 
68  const DataVector normal_magnitude =
69  sqrt(get(dot_product(*unit_normal_covector, *unit_normal_vector)));
70  for (auto& t : *unit_normal_covector) {
71  t /= normal_magnitude;
72  }
73  for (auto& t : *unit_normal_vector) {
74  t /= normal_magnitude;
75  }
76 }
77 } // namespace TestHelpers::evolution::dg::detail
get
constexpr Tag::type & get(Variables< TagList > &v) noexcept
Return Tag::type pointing into the contiguous array.
Definition: Variables.hpp:660
sqrt
auto sqrt(const TensorExpression< T, typename T::type, tmpl::list<>, tmpl::list<>, tmpl::list<>> &t)
Returns the tensor expression representing the square root of a tensor expression that evaluates to a...
Definition: SquareRoot.hpp:88
dot_product
void dot_product(const gsl::not_null< Scalar< DataType > * > dot_product, const Tensor< DataType, Symmetry< 1 >, index_list< Index >> &vector_a, const Tensor< DataType, Symmetry< 1 >, index_list< Index >> &vector_b) noexcept
Compute the Euclidean dot product of two vectors or one forms.
Definition: DotProduct.hpp:25
CREATE_HAS_TYPE_ALIAS
#define CREATE_HAS_TYPE_ALIAS(ALIAS_NAME)
Generate a type trait to check if a class has a type alias with a particular name,...
Definition: CreateHasTypeAlias.hpp:27
DotProduct.hpp
cstddef
DataVector
Stores a collection of function values.
Definition: DataVector.hpp:46
Gsl.hpp
Tensor.hpp
TMPL.hpp
gsl::not_null
Require a pointer to not be a nullptr
Definition: ReadSpecThirdOrderPiecewisePolynomial.hpp:13