PrecomputeCceDependencies.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/ComplexDataVector.hpp"
10 #include "DataStructures/SpinWeighted.hpp"
11 #include "Evolution/Systems/Cce/IntegrandInputSteps.hpp"
12 #include "Evolution/Systems/Cce/OptionTags.hpp"
13 #include "Evolution/Systems/Cce/Tags.hpp"
15 #include "NumericalAlgorithms/Spectral/SwshCollocation.hpp"
16 #include "NumericalAlgorithms/Spectral/SwshTags.hpp"
17 #include "Utilities/Gsl.hpp"
18 #include "Utilities/TMPL.hpp"
19 #include "Utilities/VectorAlgebra.hpp"
20 
21 namespace Cce {
22 
23 namespace detail {
24 // A convenience function for computing the spin-weighted derivatives of \f$R\f$
25 // divided by \f$R\f$, which appears often in Jacobians to transform between
26 // Bondi coordinates and the numerical coordinates used in CCE.
27 template <typename DerivKind>
28 void angular_derivative_of_r_divided_by_r_impl(
31  Spectral::Swsh::Tags::derivative_spin_weight<DerivKind>>*>
32  d_r_divided_by_r,
33  const SpinWeighted<ComplexDataVector, 0>& boundary_r, size_t l_max,
34  size_t number_of_radial_points) noexcept;
35 
36 } // namespace detail
37 
38 /*!
39  * \brief A set of procedures for computing the set of inputs to the CCE
40  * integrand computations that can be computed before any of the intermediate
41  * integrands are evaluated.
42  *
43  * \details The template specializations of this template are
44  * compatible with acting as a the mutator in a \ref DataBoxGroup
45  * `db::mutate_apply` operation. For flexibility in defining the \ref
46  * DataBoxGroup structure, the tags for `Tensor`s used in these functions are
47  * also organized into type lists:
48  * - type alias `integration_independent_tags`: with a subset of
49  * `Cce::pre_computation_tags`, used for both input and output.
50  * - type alias `boundary_values`: with a subset of
51  * `Cce::pre_computation_boundary_tags`, used only for input.
52  * - type alias `pre_swsh_derivatives` containing hypersurface quantities. For
53  * this struct, it will only ever contain `Cce::Tags::BondiJ`, and is used as
54  * input.
55  *
56  * The `BoundaryPrefix` tag allows easy switching between the
57  * regularity-preserving version and standard CCE
58  *
59  */
60 template <template <typename> class BoundaryPrefix, typename Tag>
62 
63 /// Computes \f$1 - y\f$ for the CCE system.
64 template <template <typename> class BoundaryPrefix>
65 struct PrecomputeCceDependencies<BoundaryPrefix, Tags::OneMinusY> {
66  using boundary_tags = tmpl::list<>;
67  using pre_swsh_derivative_tags = tmpl::list<>;
68  using integration_independent_tags = tmpl::list<>;
69 
70  using return_tags = tmpl::list<Tags::OneMinusY>;
71  using argument_tags = tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>;
72 
73  static void apply(
75  one_minus_y,
76  const size_t l_max, const size_t number_of_radial_points) noexcept {
77  const size_t number_of_angular_points =
79  const DataVector one_minus_y_collocation =
80  1.0 - Spectral::collocation_points<Spectral::Basis::Legendre,
81  Spectral::Quadrature::GaussLobatto>(
82  number_of_radial_points);
83  // iterate through the angular 'chunks' and set them to their 1-y value
84  for (size_t i = 0; i < number_of_radial_points; ++i) {
85  ComplexDataVector angular_view{
86  get(*one_minus_y).data().data() + number_of_angular_points * i,
87  number_of_angular_points};
88  angular_view = one_minus_y_collocation[i];
89  }
90  }
91 };
92 
93 /// Computes a volume version of Bondi radius of the worldtube \f$R\f$ from its
94 /// boundary value (by repeating it over the radial dimension)
95 template <template <typename> class BoundaryPrefix>
96 struct PrecomputeCceDependencies<BoundaryPrefix, Tags::BondiR> {
97  using boundary_tags = tmpl::list<BoundaryPrefix<Tags::BondiR>>;
98  using pre_swsh_derivative_tags = tmpl::list<>;
99  using integration_independent_tags = tmpl::list<>;
100 
101  using return_tags = tmpl::list<Tags::BondiR>;
102  using argument_tags =
103  tmpl::append<boundary_tags, tmpl::list<Tags::NumberOfRadialPoints>>;
104 
105  static void apply(
107  const Scalar<SpinWeighted<ComplexDataVector, 0>>& boundary_r,
108  const size_t number_of_radial_points) noexcept {
109  fill_with_n_copies(make_not_null(&get(*r).data()), get(boundary_r).data(),
110  number_of_radial_points);
111  }
112 };
113 
114 /// Computes \f$\partial_u R / R\f$ from its boundary value (by repeating it
115 /// over the radial dimension).
116 template <template <typename> class BoundaryPrefix>
117 struct PrecomputeCceDependencies<BoundaryPrefix, Tags::DuRDividedByR> {
118  using boundary_tags = tmpl::list<BoundaryPrefix<Tags::DuRDividedByR>>;
119  using pre_swsh_derivative_tags = tmpl::list<>;
120  using integration_independent_tags = tmpl::list<>;
121 
122  using return_tags = tmpl::list<Tags::DuRDividedByR>;
123  using argument_tags =
124  tmpl::append<boundary_tags, tmpl::list<Tags::NumberOfRadialPoints>>;
125 
126  static void apply(
128  du_r_divided_by_r,
130  boundary_du_r_divided_by_r,
131  const size_t number_of_radial_points) noexcept {
132  fill_with_n_copies(make_not_null(&get(*du_r_divided_by_r).data()),
133  get(boundary_du_r_divided_by_r).data(),
134  number_of_radial_points);
135  }
136 };
137 
138 /// Computes \f$\eth R / R\f$ by differentiating and repeating the boundary
139 /// value of \f$R\f$.
140 template <template <typename> class BoundaryPrefix>
141 struct PrecomputeCceDependencies<BoundaryPrefix, Tags::EthRDividedByR> {
142  using boundary_tags = tmpl::list<BoundaryPrefix<Tags::BondiR>>;
143  using pre_swsh_derivative_tags = tmpl::list<>;
144  using integration_independent_tags = tmpl::list<>;
145 
146  using return_tags = tmpl::list<Tags::EthRDividedByR>;
147  using argument_tags =
148  tmpl::append<boundary_tags,
149  tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>>;
150 
151  static void apply(
153  eth_r_divided_by_r,
154  const Scalar<SpinWeighted<ComplexDataVector, 0>>& boundary_r,
155  const size_t l_max, const size_t number_of_radial_points) noexcept {
156  detail::angular_derivative_of_r_divided_by_r_impl<
157  Spectral::Swsh::Tags::Eth>(make_not_null(&get(*eth_r_divided_by_r)),
158  get(boundary_r), l_max,
159  number_of_radial_points);
160  }
161 };
162 
163 /// Computes \f$\eth \eth R / R\f$ by differentiating and repeating the boundary
164 /// value of \f$R\f$.
165 template <template <typename> class BoundaryPrefix>
166 struct PrecomputeCceDependencies<BoundaryPrefix, Tags::EthEthRDividedByR> {
167  using boundary_tags = tmpl::list<BoundaryPrefix<Tags::BondiR>>;
168  using pre_swsh_derivative_tags = tmpl::list<>;
169  using integration_independent_tags = tmpl::list<>;
170 
171  using return_tags = tmpl::list<Tags::EthEthRDividedByR>;
172  using argument_tags =
173  tmpl::append<boundary_tags,
174  tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>>;
175 
176  static void apply(
178  eth_eth_r_divided_by_r,
179  const Scalar<SpinWeighted<ComplexDataVector, 0>>& boundary_r,
180  const size_t l_max, const size_t number_of_radial_points) noexcept {
181  detail::angular_derivative_of_r_divided_by_r_impl<
183  make_not_null(&get(*eth_eth_r_divided_by_r)), get(boundary_r), l_max,
184  number_of_radial_points);
185  }
186 };
187 
188 /// Computes \f$\eth \bar{\eth} R / R\f$ by differentiating and repeating the
189 /// boundary value of \f$R\f$.
190 template <template <typename> class BoundaryPrefix>
191 struct PrecomputeCceDependencies<BoundaryPrefix, Tags::EthEthbarRDividedByR> {
192  using boundary_tags = tmpl::list<BoundaryPrefix<Tags::BondiR>>;
193  using pre_swsh_derivative_tags = tmpl::list<>;
194  using integration_independent_tags = tmpl::list<>;
195 
196  using return_tags = tmpl::list<Tags::EthEthbarRDividedByR>;
197  using argument_tags =
198  tmpl::append<boundary_tags,
199  tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>>;
200 
201  static void apply(
203  eth_ethbar_r_divided_by_r,
204  const Scalar<SpinWeighted<ComplexDataVector, 0>>& boundary_r,
205  const size_t l_max, const size_t number_of_radial_points) noexcept {
206  detail::angular_derivative_of_r_divided_by_r_impl<
208  make_not_null(&get(*eth_ethbar_r_divided_by_r)), get(boundary_r), l_max,
209  number_of_radial_points);
210  }
211 };
212 
213 /// Computes \f$K = \sqrt{1 + J \bar{J}}\f$.
214 template <template <typename> class BoundaryPrefix>
215 struct PrecomputeCceDependencies<BoundaryPrefix, Tags::BondiK> {
216  using boundary_tags = tmpl::list<>;
217  using pre_swsh_derivative_tags = tmpl::list<Tags::BondiJ>;
218  using integration_independent_tags = tmpl::list<>;
219 
220  using return_tags = tmpl::list<Tags::BondiK>;
221  using argument_tags = tmpl::push_front<pre_swsh_derivative_tags, Tags::LMax>;
222 
223  static void apply(
225  const size_t /*l_max*/,
226  const Scalar<SpinWeighted<ComplexDataVector, 2>>& j) noexcept {
227  get(*k).data() = sqrt(1.0 + get(j).data() * conj(get(j)).data());
228  }
229 };
230 
231 /*!
232  * \brief Convenience routine for computing all of the CCE inputs to integrand
233  * computation that do not depend on intermediate integrand results. It should
234  * be executed before moving through the hierarchy of integrands.
235  *
236  * \details Provided a \ref DataBoxGroup with the appropriate tags (including
237  * `Cce::pre_computation_boundary_tags`, `Cce::pre_computation_tags`,
238  * `Cce::Tags::BondiJ` and `Tags::LMax`), this function will
239  * apply all of the necessary mutations to update the
240  * `Cce::pre_computation_tags` to their correct values for the current values
241  * for the remaining (input) tags.
242  *
243  * The `BoundaryPrefix` template template parameter is to be passed a prefix
244  * tag associated with the boundary value prefix used in the computation (e.g.
245  * `Cce::Tags::BoundaryValue`), and allows easy switching between the
246  * regularity-preserving version and standard CCE.
247  */
248 template <template <typename> class BoundaryPrefix, typename DataBoxType>
250  const gsl::not_null<DataBoxType*> box) noexcept {
251  tmpl::for_each<pre_computation_tags>([&box](auto x) {
252  using integration_independent_tag = typename decltype(x)::type;
253  using mutation =
255  db::mutate_apply<mutation>(box);
256  });
257 }
258 } // namespace Cce
A set of procedures for computing the set of inputs to the CCE integrand computations that can be com...
Definition: PrecomputeCceDependencies.hpp:61
The set of utilities for performing Cauchy characteristic evolution and Cauchy characteristic matchin...
Definition: BoundaryComputeAndSendToEvolution.hpp:28
constexpr size_t number_of_swsh_collocation_points(const size_t l_max) noexcept
Convenience function for determining the number of spin-weighted spherical harmonic collocation value...
Definition: SwshCollocation.hpp:25
const DataVector & collocation_points(const size_t num_points) noexcept
Collocation points.
Definition: Spectral.cpp:274
void mutate_all_precompute_cce_dependencies(const gsl::not_null< DataBoxType *> box) noexcept
Convenience routine for computing all of the CCE inputs to integrand computation that do not depend o...
Definition: PrecomputeCceDependencies.hpp:249
Definition: Determinant.hpp:11
void fill_with_n_copies(const gsl::not_null< VectorType *> result, const VectorType &to_copy, const size_t times_to_copy) noexcept
Creates or fills a vector with data from to_repeat copied times_to_repeat times in sequence...
Definition: VectorAlgebra.hpp:54
Make a spin-weighted type T with spin-weight Spin. Mathematical operators are restricted to addition...
Definition: SpinWeighted.hpp:24
Defines classes and functions used for manipulating DataBox&#39;s.
Stores a collection of complex function values.
Definition: ComplexDataVector.hpp:47
Definition: DataBoxTag.hpp:27
constexpr auto apply(F &&f, const DataBox< BoxTags > &box, Args &&... args) noexcept
Apply the invokable f with argument Tags TagsList from DataBox box
Definition: DataBox.hpp:1444
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:30
Stores a collection of function values.
Definition: DataVector.hpp:42
Wraps the template metaprogramming library used (brigand)
Defines functions and classes from the GSL.
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:37
gsl::not_null< T * > make_not_null(T *ptr) noexcept
Construct a not_null from a pointer. Often this will be done as an implicit conversion, but it may be necessary to perform the conversion explicitly when type deduction is desired.
Definition: Gsl.hpp:879
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:45
Tensor< T, Symmetry<>, index_list<> > Scalar
Scalar type.
Definition: TypeAliases.hpp:21
Require a pointer to not be a nullptr
Definition: Gsl.hpp:182