Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <cmath>
       7             : #include <cstddef>
       8             : #include <limits>
       9             : #include <optional>
      10             : #include <pup.h>
      11             : #include <string>
      12             : #include <type_traits>
      13             : #include <unordered_set>
      14             : #include <vector>
      15             : 
      16             : #include "DataStructures/DataBox/ObservationBox.hpp"
      17             : #include "DataStructures/DataBox/TagName.hpp"
      18             : #include "DataStructures/DataBox/ValidateSelection.hpp"
      19             : #include "DataStructures/DataVector.hpp"
      20             : #include "DataStructures/ExtractPoint.hpp"
      21             : #include "DataStructures/Tensor/Tensor.hpp"
      22             : #include "DataStructures/Tensor/TypeAliases.hpp"
      23             : #include "Domain/Tags.hpp"
      24             : #include "Options/String.hpp"
      25             : #include "ParallelAlgorithms/EventsAndTriggers/Event.hpp"
      26             : #include "Utilities/ErrorHandling/Error.hpp"
      27             : #include "Utilities/Serialization/CharmPupable.hpp"
      28             : #include "Utilities/TMPL.hpp"
      29             : #include "Utilities/TypeTraits/IsA.hpp"
      30             : 
      31             : /// \cond
      32             : namespace Frame {
      33             : struct Inertial;
      34             : }  // namespace Frame
      35             : namespace Parallel {
      36             : template <typename Metavariables>
      37             : class GlobalCache;
      38             : }  // namespace Parallel
      39             : namespace domain::Tags {
      40             : template <size_t Dim, typename Frame>
      41             : struct Coordinates;
      42             : }  // namespace domain::Tags
      43             : /// \endcond
      44             : 
      45             : namespace Events {
      46             : /*!
      47             :  * \brief ERROR if tensors get too big.
      48             :  *
      49             :  * The magnitudes of the components of the specified tensors are
      50             :  * checked, and if any exceed the specified threshold an ERROR is
      51             :  * thrown, terminating the evolution.
      52             :  *
      53             :  * Any `Tensor` in the `db::DataBox` can be checked but must be listed
      54             :  * in the `Tensors` template parameter. Any additional compute tags
      55             :  * that hold a `Tensor` can also be added to the `Tensors` template
      56             :  * parameter. Finally, `Variables` and other non-tensor compute tags
      57             :  * used to calculate tensors can be listed in the
      58             :  * `NonTensorComputeTags`.
      59             :  *
      60             :  * \note The `NonTensorComputeTags` are intended to be used for `Variables`
      61             :  * compute tags like `Tags::DerivCompute`
      62             :  */
      63             : template <size_t Dim, typename Tensors, typename NonTensorComputeTags>
      64           1 : class ErrorIfDataTooBig : public Event {
      65             :  public:
      66             :   /// \cond
      67             :   explicit ErrorIfDataTooBig(CkMigrateMessage* /*unused*/) {}
      68             :   using PUP::able::register_constructor;
      69             :   WRAPPED_PUPable_decl_template(ErrorIfDataTooBig);  // NOLINT
      70             :   /// \endcond
      71             : 
      72           0 :   struct VariablesToCheck {
      73           0 :     static constexpr Options::String help = "Subset of variables to check";
      74           0 :     using type = std::vector<std::string>;
      75           0 :     static size_t lower_bound_on_size() { return 1; }
      76             :   };
      77             : 
      78           0 :   struct Threshold {
      79           0 :     static constexpr Options::String help = "Threshold at which to ERROR";
      80           0 :     using type = double;
      81           0 :     static type lower_bound() { return 0.0; }
      82             :   };
      83             : 
      84           0 :   using options = tmpl::list<VariablesToCheck, Threshold>;
      85             : 
      86           0 :   static constexpr Options::String help = "ERROR if tensors get too big";
      87             : 
      88           0 :   ErrorIfDataTooBig() = default;
      89             : 
      90           0 :   ErrorIfDataTooBig(const std::vector<std::string>& variables_to_check,
      91             :                     const double threshold, const Options::Context& context)
      92             :       : variables_to_check_(variables_to_check.begin(),
      93             :                             variables_to_check.end()),
      94             :         threshold_(threshold) {
      95             :     db::validate_selection<Tensors>(variables_to_check, context);
      96             :   }
      97             : 
      98           0 :   using compute_tags_for_observation_box =
      99             :       tmpl::append<Tensors, NonTensorComputeTags>;
     100             : 
     101           0 :   using return_tags = tmpl::list<>;
     102           0 :   using argument_tags =
     103             :       tmpl::list<::Tags::ObservationBox,
     104             :                  ::domain::Tags::Coordinates<Dim, Frame::Inertial>>;
     105             : 
     106             :   template <typename ComputeTagsList, typename DataBoxType,
     107             :             typename Metavariables, typename ArrayIndex,
     108             :             typename ParallelComponent>
     109           0 :   void operator()(const ObservationBox<ComputeTagsList, DataBoxType>& box,
     110             :                   const tnsr::I<DataVector, Dim>& coordinates,
     111             :                   Parallel::GlobalCache<Metavariables>& /*cache*/,
     112             :                   const ArrayIndex& /*array_index*/,
     113             :                   const ParallelComponent* const /*component*/,
     114             :                   const ObservationValue& /*observation_value*/) const {
     115             :     tmpl::for_each<Tensors>([&](const auto tensor_tag_v) {
     116             :       using tensor_tag = tmpl::type_from<decltype(tensor_tag_v)>;
     117             :       const std::string tag_name = db::tag_name<tensor_tag>();
     118             :       if (variables_to_check_.count(tag_name) != 0) {
     119             :         const auto check_components = [&](const auto& tensor) {
     120             :           for (const auto& component : tensor) {
     121             :             for (size_t point = 0; point < component.size(); ++point) {
     122             :               if (std::abs(component[point]) > threshold_) {
     123             :                 ERROR_NO_TRACE(
     124             :                     tag_name
     125             :                     << " too big with value " << extract_point(tensor, point)
     126             :                     << " at position\n"
     127             :                     << extract_point(coordinates, point) << "\nwith ElementId: "
     128             :                     << get<::domain::Tags::Element<Dim>>(box).id());
     129             :               }
     130             :             }
     131             :           }
     132             :         };
     133             :         const auto& tensor = get<tensor_tag>(box);
     134             :         if constexpr (tt::is_a_v<std::optional,
     135             :                                  std::decay_t<decltype(tensor)>>) {
     136             :           if (tensor.has_value()) {
     137             :             check_components(*tensor);
     138             :           }
     139             :         } else {
     140             :           check_components(tensor);
     141             :         }
     142             :       }
     143             :     });
     144             :   }
     145             : 
     146           0 :   using is_ready_argument_tags = tmpl::list<>;
     147             : 
     148             :   template <typename Metavariables, typename ArrayIndex, typename Component>
     149           0 :   bool is_ready(Parallel::GlobalCache<Metavariables>& /*cache*/,
     150             :                 const ArrayIndex& /*array_index*/,
     151             :                 const Component* const /*meta*/) const {
     152             :     return true;
     153             :   }
     154             : 
     155           1 :   bool needs_evolved_variables() const override { return true; }
     156             : 
     157             :   // NOLINTNEXTLINE(google-runtime-references)
     158           0 :   void pup(PUP::er& p) override {
     159             :     Event::pup(p);
     160             :     p | variables_to_check_;
     161             :     p | threshold_;
     162             :   }
     163             : 
     164             :  private:
     165           0 :   std::unordered_set<std::string> variables_to_check_{};
     166           0 :   double threshold_ = std::numeric_limits<double>::signaling_NaN();
     167             : };
     168             : 
     169             : /// \cond
     170             : template <size_t Dim, typename Tensors, typename NonTensorComputeTags>
     171             : PUP::able::PUP_ID
     172             :     ErrorIfDataTooBig<Dim, Tensors,
     173             :                       NonTensorComputeTags>::my_PUP_ID = 0;  // NOLINT
     174             : /// \endcond
     175             : }  // namespace Events