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/Spectral/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