Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <array>
       7             : #include <cstddef>
       8             : #include <pup.h>
       9             : #include <string>
      10             : #include <vector>
      11             : 
      12             : #include "DataStructures/DataBox/DataBox.hpp"
      13             : #include "DataStructures/DataBox/DataBoxTag.hpp"
      14             : #include "DataStructures/DataBox/ValidateSelection.hpp"
      15             : #include "DataStructures/DataVector.hpp"
      16             : #include "DataStructures/Tensor/Tensor.hpp"
      17             : #include "Domain/Amr/Flag.hpp"
      18             : #include "Domain/Tags.hpp"
      19             : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
      20             : #include "Options/Context.hpp"
      21             : #include "Options/ParseError.hpp"
      22             : #include "Options/String.hpp"
      23             : #include "ParallelAlgorithms/Amr/Criteria/Criterion.hpp"
      24             : #include "Utilities/Algorithm.hpp"
      25             : #include "Utilities/TMPL.hpp"
      26             : 
      27             : /// \cond
      28             : template <size_t>
      29             : class ElementId;
      30             : /// \endcond
      31             : 
      32             : namespace amr::Criteria {
      33             : 
      34             : namespace TruncationError_detail {
      35             : /*!
      36             :  * \brief Apply the truncation error criterion to a single tensor component
      37             :  *
      38             :  * The `result` is the current decision in each dimension based on the previous
      39             :  * tensor components. This function will update the flags if necessary. It takes
      40             :  * the "max" of the current and new flags, where the "highest" flag is
      41             :  * `Flag::IncreaseResolution`, followed by `Flag::DoNothing`, and then
      42             :  * `Flag::DecreaseResolution`.
      43             :  */
      44             : template <typename ValueType, size_t Dim>
      45             : void max_over_components(
      46             :     gsl::not_null<std::array<Flag, Dim>*> result,
      47             :     const gsl::not_null<std::array<DataVector, Dim>*> power_monitors_buffer,
      48             :     const ValueType& tensor_component, const Mesh<Dim>& mesh,
      49             :     double target_abs_truncation_error, double target_rel_truncation_error);
      50             : }  // namespace TruncationError_detail
      51             : 
      52             : /*!
      53             :  * \brief Refine the grid towards the target truncation error
      54             :  *
      55             :  * - If any tensor component has a truncation error above the target value, the
      56             :  *   element will be p-refined.
      57             :  * - If all tensor components still satisfy the target even with one mode
      58             :  *   removed, the element will be p-coarsened.
      59             :  *
      60             :  * For details on how the truncation error is computed see
      61             :  * `PowerMonitors::truncation_error`.
      62             :  *
      63             :  * \tparam Dim Spatial dimension of the grid
      64             :  * \tparam TensorTags List of tags of the tensors to be monitored
      65             :  */
      66             : template <size_t Dim, typename TensorTags>
      67           1 : class TruncationError : public Criterion {
      68             :  public:
      69           0 :   struct VariablesToMonitor {
      70           0 :     using type = std::vector<std::string>;
      71           0 :     static constexpr Options::String help = {
      72             :         "The tensors to monitor the truncation error of."};
      73           0 :     static size_t lower_bound_on_size() { return 1; }
      74             :   };
      75           0 :   struct AbsoluteTargetTruncationError {
      76           0 :     static std::string name() { return "AbsoluteTarget"; }
      77           0 :     using type = double;
      78           0 :     static type lower_bound() { return 0.0; }
      79           0 :     static constexpr Options::String help = {
      80             :         "The absolute target truncation error."};
      81             :   };
      82           0 :   struct RelativeTargetTruncationError {
      83           0 :     static std::string name() { return "RelativeTarget"; }
      84           0 :     using type = double;
      85           0 :     static type lower_bound() { return 0.0; }
      86           0 :     static constexpr Options::String help = {
      87             :         "The relative target truncation error."};
      88             :   };
      89             : 
      90           0 :   using options = tmpl::list<VariablesToMonitor, AbsoluteTargetTruncationError,
      91             :                              RelativeTargetTruncationError>;
      92             : 
      93           0 :   static constexpr Options::String help = {
      94             :       "Refine the grid towards the target truncation error"};
      95             : 
      96           0 :   TruncationError() = default;
      97             : 
      98           0 :   TruncationError(std::vector<std::string> vars_to_monitor,
      99             :                   double target_abs_truncation_error,
     100             :                   double target_rel_truncation_error,
     101             :                   const Options::Context& context = {});
     102             : 
     103             :   /// \cond
     104             :   explicit TruncationError(CkMigrateMessage* msg);
     105             :   using PUP::able::register_constructor;
     106             :   WRAPPED_PUPable_decl_template(TruncationError);  // NOLINT
     107             :   /// \endcond
     108             : 
     109           0 :   std::string observation_name() override { return "TruncationError"; }
     110             : 
     111           0 :   using compute_tags_for_observation_box = tmpl::list<>;
     112             : 
     113           0 :   using argument_tags = tmpl::list<::Tags::DataBox>;
     114             : 
     115             :   template <typename DbTagsList, typename Metavariables>
     116           0 :   std::array<Flag, Dim> operator()(const db::DataBox<DbTagsList>& box,
     117             :                                    Parallel::GlobalCache<Metavariables>& cache,
     118             :                                    const ElementId<Dim>& element_id) const;
     119             : 
     120           0 :   void pup(PUP::er& p) override;
     121             : 
     122             :  private:
     123           0 :   std::vector<std::string> vars_to_monitor_{};
     124           0 :   double target_abs_truncation_error_{};
     125           0 :   double target_rel_truncation_error_{};
     126             : };
     127             : 
     128             : // Out-of-line definitions
     129             : /// \cond
     130             : 
     131             : template <size_t Dim, typename TensorTags>
     132             : TruncationError<Dim, TensorTags>::TruncationError(
     133             :     std::vector<std::string> vars_to_monitor,
     134             :     const double target_abs_truncation_error,
     135             :     const double target_rel_truncation_error, const Options::Context& context)
     136             :     : vars_to_monitor_(std::move(vars_to_monitor)),
     137             :       target_abs_truncation_error_(target_abs_truncation_error),
     138             :       target_rel_truncation_error_(target_rel_truncation_error) {
     139             :   db::validate_selection<TensorTags>(vars_to_monitor_, context);
     140             : }
     141             : 
     142             : template <size_t Dim, typename TensorTags>
     143             : TruncationError<Dim, TensorTags>::TruncationError(CkMigrateMessage* msg)
     144             :     : Criterion(msg) {}
     145             : 
     146             : template <size_t Dim, typename TensorTags>
     147             : template <typename DbTagsList, typename Metavariables>
     148             : std::array<Flag, Dim> TruncationError<Dim, TensorTags>::operator()(
     149             :     const db::DataBox<DbTagsList>& box,
     150             :     Parallel::GlobalCache<Metavariables>& /*cache*/,
     151             :     const ElementId<Dim>& /*element_id*/) const {
     152             :   auto result = make_array<Dim>(Flag::Undefined);
     153             :   const auto& mesh = db::get<domain::Tags::Mesh<Dim>>(box);
     154             :   std::array<DataVector, Dim> power_monitors_buffer{};
     155             :   // Check all tensors and all tensor components in turn
     156             :   tmpl::for_each<TensorTags>(
     157             :       [&result, &box, &mesh, &power_monitors_buffer, this](const auto tag_v) {
     158             :         // Stop if we have already decided to refine every dimension
     159             :         if (result == make_array<Dim>(Flag::IncreaseResolution)) {
     160             :           return;
     161             :         }
     162             :         using tag = tmpl::type_from<std::decay_t<decltype(tag_v)>>;
     163             :         const std::string tag_name = db::tag_name<tag>();
     164             :         // Skip if this tensor is not being monitored
     165             :         if (not alg::found(vars_to_monitor_, tag_name)) {
     166             :           return;
     167             :         }
     168             :         const auto& tensor = db::get<tag>(box);
     169             :         for (const auto& tensor_component : tensor) {
     170             :           TruncationError_detail::max_over_components(
     171             :               make_not_null(&result), make_not_null(&power_monitors_buffer),
     172             :               tensor_component, mesh, target_abs_truncation_error_,
     173             :               target_rel_truncation_error_);
     174             :         }
     175             :       });
     176             :   return result;
     177             : }
     178             : 
     179             : template <size_t Dim, typename TensorTags>
     180             : void TruncationError<Dim, TensorTags>::pup(PUP::er& p) {
     181             :   p | vars_to_monitor_;
     182             :   p | target_abs_truncation_error_;
     183             :   p | target_rel_truncation_error_;
     184             : }
     185             : 
     186             : template <size_t Dim, typename TensorTags>
     187             : PUP::able::PUP_ID TruncationError<Dim, TensorTags>::my_PUP_ID = 0;  // NOLINT
     188             : /// \endcond
     189             : 
     190             : }  // namespace amr::Criteria