Line data    Source code

       1            0 : // Distributed under the MIT License.
       2              : // See LICENSE.txt for details.
       3              : 
       4              : #pragma once
       5              : 
       6              : #include <algorithm>
       7              : #include <cmath>
       8              : #include <cstddef>
       9              : #include <pup.h>
      10              : #include <string>
      11              : 
      12              : #include "DataStructures/DataVector.hpp"
      13              : #include "Evolution/Systems/GeneralizedHarmonic/Bbh/CompletionCriteria.hpp"
      14              : #include "Evolution/Systems/GeneralizedHarmonic/Bbh/CompletionSingleton.hpp"
      15              : #include "Evolution/Systems/GeneralizedHarmonic/Constraints.hpp"
      16              : #include "Evolution/Systems/GeneralizedHarmonic/Tags.hpp"
      17              : #include "Options/String.hpp"
      18              : #include "Parallel/ArrayCollection/IsDgElementCollection.hpp"
      19              : #include "Parallel/GlobalCache.hpp"
      20              : #include "Parallel/Reduction.hpp"
      21              : #include "ParallelAlgorithms/EventsAndTriggers/Event.hpp"
      22              : #include "Time/Tags/Time.hpp"
      23              : #include "Utilities/ErrorHandling/Error.hpp"
      24              : #include "Utilities/Serialization/CharmPupable.hpp"
      25              : #include "Utilities/TMPL.hpp"
      26              : 
      27            0 : namespace gh::bbh::Events {
      28              : /*!
      29              :  * \brief Computes per-element Linf norms of the generalized-harmonic
      30              :  * constraints, latching BBH completion criteria when thresholds are exceeded.
      31              :  *
      32              :  * \details This event is intended to run once common-horizon finding is
      33              :  * active (typically behind a separation-based trigger). It monitors the gauge
      34              :  * and three-index constraints, using a reduction to determine if their Linf
      35              :  * norms exceed completion thresholds.
      36              :  */
      37            1 : class CheckConstraintThresholds : public Event {
      38            0 :   using ReductionData = Parallel::ReductionData<
      39              :       Parallel::ReductionDatum<double, funcl::AssertEqual<>>,
      40              :       Parallel::ReductionDatum<double, funcl::Max<>>,
      41              :       Parallel::ReductionDatum<double, funcl::Max<>>>;
      42              : 
      43              :  public:
      44              :   /// \cond
      45              :   explicit CheckConstraintThresholds(CkMigrateMessage* /*unused*/) {}
      46              :   using PUP::able::register_constructor;
      47              :   WRAPPED_PUPable_decl_template(CheckConstraintThresholds);  // NOLINT
      48              :   /// \endcond
      49              : 
      50            0 :   using compute_tags_for_observation_box =
      51              :       tmpl::list<gh::Tags::GaugeConstraintCompute<3, Frame::Inertial>,
      52              :                  gh::Tags::ThreeIndexConstraintCompute<3, Frame::Inertial>>;
      53            0 :   using options = tmpl::list<>;
      54            0 :   static constexpr Options::String help =
      55              :       "Checks local Linf norms of constraints against BBH completion "
      56              :       "thresholds and forwards reduced maxima to the BBH completion singleton "
      57              :       "reduction callback.";
      58            0 :   static std::string name() { return "BbhCheckConstraintThresholds"; }
      59              : 
      60            0 :   CheckConstraintThresholds() = default;
      61              : 
      62            0 :   using return_tags = tmpl::list<>;
      63            0 :   using argument_tags = tmpl::list<
      64              :       ::Tags::Time, gh::Tags::GaugeConstraint<DataVector, 3, Frame::Inertial>,
      65              :       gh::Tags::ThreeIndexConstraint<DataVector, 3, Frame::Inertial>>;
      66              : 
      67              :   template <typename Metavariables, typename ArrayIndex, typename Component>
      68            0 :   void operator()(
      69              :       const double time,
      70              :       const tnsr::a<DataVector, 3, Frame::Inertial>& gauge_constraint,
      71              :       const tnsr::iaa<DataVector, 3, Frame::Inertial>& three_index_constraint,
      72              :       Parallel::GlobalCache<Metavariables>& cache,
      73              :       const ArrayIndex& array_index, const Component* const /*component*/,
      74              :       const ObservationValue& /*observation_value*/) const {
      75              :     const double local_gauge_linf = local_linf_norm(gauge_constraint);
      76              :     const double local_three_index_linf =
      77              :         local_linf_norm(three_index_constraint);
      78              :     if constexpr (Parallel::is_dg_element_collection_v<Component>) {
      79              :       ERROR(
      80              :           "BbhCheckConstraintThresholds currently requires array components "
      81              :           "(not DgElementCollection).");
      82              :     } else {
      83              :       const auto& self_proxy =
      84              :           Parallel::get_parallel_component<Component>(cache)[array_index];
      85              :       auto& reduction_target_proxy = Parallel::get_parallel_component<
      86              :           gh::bbh::CompletionSingleton<Metavariables>>(cache);
      87              :       Parallel::contribute_to_reduction<
      88              :           gh::bbh::Actions::ProcessConstraintMaxima>(
      89              :           ReductionData{time, local_gauge_linf, local_three_index_linf},
      90              :           self_proxy, reduction_target_proxy);
      91              :     }
      92              :   }
      93              : 
      94            0 :   using is_ready_argument_tags = tmpl::list<>;
      95              :   template <typename Metavariables, typename ArrayIndex, typename Component>
      96            0 :   bool is_ready(Parallel::GlobalCache<Metavariables>& /*cache*/,
      97              :                 const ArrayIndex& /*array_index*/,
      98              :                 const Component* const /*meta*/) const {
      99              :     return true;
     100              :   }
     101              : 
     102            1 :   bool needs_evolved_variables() const override { return true; }
     103              : 
     104              :  private:
     105              :   template <typename TensorType>
     106            0 :   static double local_linf_norm(const TensorType& tensor) {
     107              :     double result = 0.0;
     108              :     for (size_t storage_index = 0; storage_index < tensor.size();
     109              :          ++storage_index) {
     110              :       const auto& component = tensor[storage_index];
     111              :       result = std::max(result, max(abs(component)));
     112              :     }
     113              :     return result;
     114              :   }
     115              : };
     116              : }  // namespace gh::bbh::Events