Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include "DataStructures/SpinWeighted.hpp"
7 : #include "DataStructures/Tensor/TypeAliases.hpp"
8 : #include "Evolution/Systems/Cce/OptionTags.hpp"
9 : #include "Evolution/Systems/Cce/Tags.hpp"
10 : #include "NumericalAlgorithms/SpinWeightedSphericalHarmonics/SwshTags.hpp"
11 :
12 : /// \cond
13 : class ComplexDataVector;
14 : class Matrix;
15 : /// \endcond
16 :
17 : namespace Cce {
18 :
19 : /*!
20 : * \brief Provides access to a lazily cached integration matrix for the \f$Q\f$
21 : * and \f$W\f$ equations in CCE hypersurface evaluation.
22 : *
23 : * \details The provided matrix acts on the integrand collocation points and
24 : * solves the equation,
25 : *
26 : * \f[
27 : * (1 - y) \partial_y f + 2 f = g,
28 : * \f]
29 : *
30 : * for \f$f\f$ given integrand \f$g\f$.
31 : */
32 1 : const Matrix& precomputed_cce_q_integrator(size_t number_of_radial_grid_points);
33 :
34 : /*!
35 : * \brief A utility function for evaluating the \f$Q\f$ and \f$W\f$ hypersurface
36 : * integrals during CCE evolution.
37 : *
38 : * \details Computes and returns by `not_null` pointer the solution to the
39 : * equation
40 : *
41 : * \f[
42 : * (1 - y) \partial_y f + 2 f = A + (1 - y) B,
43 : * \f]
44 : *
45 : * where \f$A\f$ is provided as `pole_of_integrand` and \f$B\f$ is provided as
46 : * `regular_integrand`. The value `one_minus_y` is required for determining the
47 : * integrand and `l_max` is required to determine the shape of the spin-weighted
48 : * spherical harmonic mesh.
49 : */
50 1 : void radial_integrate_cce_pole_equations(
51 : gsl::not_null<ComplexDataVector*> integral_result,
52 : const ComplexDataVector& pole_of_integrand,
53 : const ComplexDataVector& regular_integrand,
54 : const ComplexDataVector& boundary, const ComplexDataVector& one_minus_y,
55 : size_t l_max, size_t number_of_radial_points);
56 :
57 : namespace detail {
58 : // needed because the standard transpose utility cannot create an arbitrary
59 : // ordering of blocks of data. This returns by pointer the configuration useful
60 : // for the linear solve step for H integration
61 : void transpose_to_reals_then_imags_radial_stripes(
62 : gsl::not_null<DataVector*> result, const ComplexDataVector& input,
63 : size_t number_of_radial_points, size_t number_of_angular_points);
64 : } // namespace detail
65 :
66 : /// @{
67 : /*!
68 : * \brief Computational structs for evaluating the hypersurface integrals during
69 : * CCE evolution. These are compatible with use in `db::mutate_apply`.
70 : *
71 : * \details
72 : * The integral evaluated and the corresponding inputs required depend on the
73 : * CCE quantity being computed. In any of these, the only mutated tag is `Tag`,
74 : * where the result of the integration is placed. The supported `Tag`s act in
75 : * the following ways:
76 : * - If the `Tag` is `Tags::BondiBeta` or `Tags::BondiU`, the integral to be
77 : * evaluated is simply \f[ \partial_y f = A, \f] where \f$A\f$ is retrieved with
78 : * `Tags::Integrand<Tag>`.
79 : * - If the `Tag` is `Tags::BondiQ` or `Tags::BondiW`, the integral to be
80 : * evaluated is \f[ (1 - y) \partial_y f + 2 f = A + (1 - y) B, \f] where
81 : * \f$A\f$ is retrieved with `Tags::PoleOfIntegrand<Tag>` and \f$B\f$ is
82 : * retrieved with `Tags::RegularIntegrand<Tag>`.
83 : * - If `Tag` is `Tags::BondiH`, the integral to be evaluated is:
84 : *
85 : * \f[
86 : * (1 - y) \partial_y f + L f + L^\prime \bar{f} = A + (1 - y) B,
87 : * \f]
88 : *
89 : * for \f$f\f$, where \f$A\f$ is retrieved with `Tags::PoleOfIntegrand<Tag>`,
90 : * \f$B\f$ is retrieved with `Tags::RegularIntegrand<Tag>`, \f$L\f$ is retrieved
91 : * with `Tags::LinearFactor<Tag>`, and \f$L^\prime\f$ is retrieved with
92 : * `Tags::LinearFactorForConjugate<Tag>`. The presence of \f$L\f$ and
93 : * \f$L^\prime\f$ ensure that the only current method we have for evaluating the
94 : * \f$H\f$ hypersurface equation is a direct linear solve, rather than the
95 : * spectral matrix multiplications which are available for the other integrals.
96 : *
97 : * In each case, the boundary value at the world tube for the integration is
98 : * retrieved from `BoundaryPrefix<Tag>`.
99 : *
100 : * Additional type aliases `boundary_tags` and `integrand_tags` are provided for
101 : * template processing of the required input tags necessary for these functions.
102 : * These type aliases are `tmpl::list`s with the subsets of `argument_tags` from
103 : * specific other parts of the CCE computation. Because they play different
104 : * roles, and have different extents, it is better for tag management to give
105 : * separated lists for the dependencies.
106 : */
107 : template <template <typename> class BoundaryPrefix, typename Tag>
108 1 : struct RadialIntegrateBondi {
109 0 : using boundary_tags = tmpl::list<BoundaryPrefix<Tag>>;
110 0 : using integrand_tags = tmpl::list<Tags::Integrand<Tag>>;
111 :
112 0 : using return_tags = tmpl::list<Tag>;
113 0 : using argument_tags =
114 : tmpl::append<integrand_tags, boundary_tags,
115 : tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>>;
116 0 : static void apply(
117 : gsl::not_null<
118 : Scalar<SpinWeighted<ComplexDataVector, Tag::type::type::spin>>*>
119 : integral_result,
120 : const Scalar<SpinWeighted<ComplexDataVector, Tag::type::type::spin>>&
121 : integrand,
122 : const Scalar<SpinWeighted<ComplexDataVector, Tag::type::type::spin>>&
123 : boundary,
124 : size_t l_max, size_t number_of_radial_points);
125 : };
126 :
127 : template <template <typename> class BoundaryPrefix>
128 0 : struct RadialIntegrateBondi<BoundaryPrefix, Tags::BondiQ> {
129 0 : using boundary_tags = tmpl::list<BoundaryPrefix<Tags::BondiQ>>;
130 0 : using integrand_tags = tmpl::list<Tags::PoleOfIntegrand<Tags::BondiQ>,
131 : Tags::RegularIntegrand<Tags::BondiQ>>;
132 0 : using integration_independent_tags = tmpl::list<Tags::OneMinusY>;
133 :
134 0 : using return_tags = tmpl::list<Tags::BondiQ>;
135 0 : using argument_tags =
136 : tmpl::append<integrand_tags, boundary_tags, integration_independent_tags,
137 : tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>>;
138 0 : static void apply(
139 : gsl::not_null<Scalar<SpinWeighted<ComplexDataVector, 1>>*>
140 : integral_result,
141 : const Scalar<SpinWeighted<ComplexDataVector, 1>>& pole_of_integrand,
142 : const Scalar<SpinWeighted<ComplexDataVector, 1>>& regular_integrand,
143 : const Scalar<SpinWeighted<ComplexDataVector, 1>>& boundary,
144 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
145 : size_t l_max, size_t number_of_radial_points);
146 : };
147 :
148 : template <template <typename> class BoundaryPrefix>
149 0 : struct RadialIntegrateBondi<BoundaryPrefix, Tags::BondiW> {
150 0 : using boundary_tags = tmpl::list<BoundaryPrefix<Tags::BondiW>>;
151 0 : using integrand_tags = tmpl::list<Tags::PoleOfIntegrand<Tags::BondiW>,
152 : Tags::RegularIntegrand<Tags::BondiW>>;
153 0 : using integration_independent_tags = tmpl::list<Tags::OneMinusY>;
154 :
155 0 : using return_tags = tmpl::list<Tags::BondiW>;
156 0 : using argument_tags =
157 : tmpl::append<integrand_tags, boundary_tags, integration_independent_tags,
158 : tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>>;
159 0 : static void apply(
160 : gsl::not_null<Scalar<SpinWeighted<ComplexDataVector, 0>>*>
161 : integral_result,
162 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& pole_of_integrand,
163 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& regular_integrand,
164 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& boundary,
165 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
166 : size_t l_max, size_t number_of_radial_points);
167 : };
168 :
169 : template <template <typename> class BoundaryPrefix>
170 0 : struct RadialIntegrateBondi<BoundaryPrefix, Tags::BondiH> {
171 0 : using boundary_tags = tmpl::list<BoundaryPrefix<Tags::BondiH>>;
172 0 : using integrand_tags =
173 : tmpl::list<Tags::PoleOfIntegrand<Tags::BondiH>,
174 : Tags::RegularIntegrand<Tags::BondiH>,
175 : Tags::LinearFactor<Tags::BondiH>,
176 : Tags::LinearFactorForConjugate<Tags::BondiH>>;
177 0 : using integration_independent_tags = tmpl::list<Tags::OneMinusY>;
178 :
179 0 : using return_tags = tmpl::list<Tags::BondiH>;
180 0 : using argument_tags =
181 : tmpl::append<integrand_tags, boundary_tags, integration_independent_tags,
182 : tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>>;
183 0 : static void apply(
184 : gsl::not_null<Scalar<SpinWeighted<ComplexDataVector, 2>>*>
185 : integral_result,
186 : const Scalar<SpinWeighted<ComplexDataVector, 2>>& pole_of_integrand,
187 : const Scalar<SpinWeighted<ComplexDataVector, 2>>& regular_integrand,
188 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& linear_factor,
189 : const Scalar<SpinWeighted<ComplexDataVector, 4>>&
190 : linear_factor_of_conjugate,
191 : const Scalar<SpinWeighted<ComplexDataVector, 2>>& boundary,
192 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
193 : size_t l_max, size_t number_of_radial_points);
194 : };
195 : /// @}
196 : } // namespace Cce
|