Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include "Evolution/Systems/Cce/OptionTags.hpp"
       7             : #include "Evolution/Systems/Cce/Tags.hpp"
       8             : #include "NumericalAlgorithms/Spectral/Spectral.hpp"
       9             : #include "Parallel/AlgorithmExecution.hpp"
      10             : #include "Parallel/GlobalCache.hpp"
      11             : #include "Time/Tags/TimeStepId.hpp"
      12             : #include "Utilities/Gsl.hpp"
      13             : #include "Utilities/TMPL.hpp"
      14             : 
      15             : namespace Cce::Actions {
      16             : 
      17             : /*!
      18             :  * \ingroup ActionsGroup
      19             :  * \brief Given initial boundary data for the Klein-Gordon variable \f$\Psi\f$,
      20             :  * computes its initial hypersurface data.
      21             :  *
      22             :  * \details This action is to be called after boundary data has been received,
      23             :  * but before the time-stepping evolution loop. So, it should be either late in
      24             :  * an initialization phase or early (before a `Actions::Goto` loop or similar)
      25             :  * in the `Evolve` phase.
      26             :  */
      27           1 : struct InitializeKleinGordonFirstHypersurface {
      28           0 :   using const_global_cache_tags =
      29             :       tmpl::list<Tags::LMax, Tags::NumberOfRadialPoints>;
      30             : 
      31             :   template <typename DbTags, typename... InboxTags, typename Metavariables,
      32             :             typename ArrayIndex, typename ActionList,
      33             :             typename ParallelComponent>
      34           0 :   static Parallel::iterable_action_return_t apply(
      35             :       db::DataBox<DbTags>& box,
      36             :       const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
      37             :       Parallel::GlobalCache<Metavariables>& /*cache*/,
      38             :       const ArrayIndex& /*array_index*/, const ActionList /*meta*/,
      39             :       const ParallelComponent* const /*meta*/) {
      40             :     // In some contexts, this action may get re-run (e.g. self-start procedure)
      41             :     // In those cases, we do not want to alter the existing hypersurface data,
      42             :     // so we just exit. However, we do want to re-run the action each time
      43             :     // the self start 'reset's from the beginning
      44             :     if (db::get<::Tags::TimeStepId>(box).slab_number() > 0 or
      45             :         not db::get<::Tags::TimeStepId>(box).is_at_slab_boundary()) {
      46             :       return {Parallel::AlgorithmExecution::Continue, std::nullopt};
      47             :     }
      48             : 
      49             :     db::mutate<Tags::KleinGordonPsi>(
      50             :         [](const gsl::not_null<Scalar<SpinWeighted<ComplexDataVector, 0>>*>
      51             :                kg_psi_volume,
      52             :            const Scalar<SpinWeighted<ComplexDataVector, 0>>& boundary_kg_psi,
      53             :            const size_t l_max, const size_t number_of_radial_points) {
      54             :           const DataVector one_minus_y_collocation =
      55             :               1.0 -
      56             :               Spectral::collocation_points<Spectral::Basis::Legendre,
      57             :                                            Spectral::Quadrature::GaussLobatto>(
      58             :                   number_of_radial_points);
      59             : 
      60             :           const size_t boundary_size =
      61             :               Spectral::Swsh::number_of_swsh_collocation_points(l_max);
      62             : 
      63             :           for (size_t i = 0; i < number_of_radial_points; i++) {
      64             :             ComplexDataVector angular_view_kg_psi{
      65             :                 get(*kg_psi_volume).data().data() + boundary_size * i,
      66             :                 boundary_size};
      67             : 
      68             :             // this gives psi = (boundary value) / r
      69             :             angular_view_kg_psi =
      70             :                 get(boundary_kg_psi).data() * one_minus_y_collocation[i] / 2.;
      71             :           }
      72             :         },
      73             :         make_not_null(&box),
      74             :         db::get<Tags::BoundaryValue<Tags::KleinGordonPsi>>(box),
      75             :         db::get<Tags::LMax>(box), db::get<Tags::NumberOfRadialPoints>(box));
      76             : 
      77             :     return {Parallel::AlgorithmExecution::Continue, std::nullopt};
      78             :   }
      79             : };
      80             : }  // namespace Cce::Actions