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 <set>
       8             : #include <sstream>
       9             : 
      10             : #include "DataStructures/DataBox/DataBox.hpp"
      11             : #include "DataStructures/LinkedMessageId.hpp"
      12             : #include "IO/Logging/Verbosity.hpp"
      13             : #include "NumericalAlgorithms/SphericalHarmonics/StrahlkorperFunctions.hpp"
      14             : #include "Parallel/GlobalCache.hpp"
      15             : #include "Parallel/Printf/Printf.hpp"
      16             : #include "ParallelAlgorithms/ApparentHorizonFinder/FastFlow.hpp"
      17             : #include "ParallelAlgorithms/ApparentHorizonFinder/Protocols/Callback.hpp"
      18             : #include "ParallelAlgorithms/ApparentHorizonFinder/Tags.hpp"
      19             : #include "Utilities/ErrorHandling/Error.hpp"
      20             : #include "Utilities/Numeric.hpp"
      21             : #include "Utilities/PrettyType.hpp"
      22             : #include "Utilities/ProtocolHelpers.hpp"
      23             : 
      24             : namespace ah::callbacks {
      25             : /*!
      26             :  * \brief Callback when an apparent horizon find fails. The template \p Ignore
      27             :  * says whether to ignore the failure or raise an ERROR.
      28             :  */
      29             : template <typename HorizonMetavars, bool Ignore>
      30           1 : struct FailedHorizonFind : tt::ConformsTo<ah::protocols::Callback> {
      31             :  private:
      32           0 :   using Fr = typename HorizonMetavars::frame;
      33             : 
      34             :  public:
      35             :   template <typename DbTags, typename Metavariables>
      36           0 :   static void apply(db::DataBox<DbTags>& box,
      37             :                     const Parallel::GlobalCache<Metavariables>& /*cache*/,
      38             :                     const FastFlow::Status failure_reason) {
      39             :     std::ostringstream os{};
      40             : 
      41             :     const auto& time = db::get<ah::Tags::CurrentTime>(box).value();
      42             :     const auto& all_storage = db::get<ah::Tags::Storage<Fr>>(box);
      43             :     const auto& current_time_storage = all_storage.at(time);
      44             :     const auto& current_iteration = current_time_storage.current_iteration;
      45             : 
      46             :     os << pretty_type::name<HorizonMetavars>()
      47             :        << ": Horizon find failed at time " << time
      48             :        << ". Reason = " << failure_reason
      49             :        << ". Number of compute coords retries = "
      50             :        << current_iteration.compute_coords_retries << ".";
      51             : 
      52             :     if (failure_reason == FastFlow::Status::InterpolationFailure) {
      53             :       const auto& block_coord_holders =
      54             :           current_iteration.block_coord_holders.value();
      55             : 
      56             :       std::vector<size_t> missing_indices{};
      57             :       missing_indices.reserve(block_coord_holders.size());
      58             :       for (size_t i = 0; i < block_coord_holders.size(); i++) {
      59             :         if (not block_coord_holders[i].has_value()) {
      60             :           missing_indices.push_back(i);
      61             :         }
      62             :       }
      63             : 
      64             :       // Get the actual points
      65             :       const auto& strahlkorper = current_iteration.strahlkorper;
      66             :       const auto& fast_flow = db::get<ah::Tags::FastFlow>(box);
      67             :       const size_t l_mesh = fast_flow.current_l_mesh(strahlkorper);
      68             :       const auto prolonged_strahlkorper =
      69             :           ylm::Strahlkorper<Fr>(l_mesh, l_mesh, strahlkorper);
      70             :       const auto coords = ylm::cartesian_coords(prolonged_strahlkorper);
      71             : 
      72             :       // Now output some information about them
      73             :       os << "\n Invalid points (in " << pretty_type::name<Fr>()
      74             :          << " frame) are:\n";
      75             :       for (const size_t index : missing_indices) {
      76             :         os << " (" << get<0>(coords)[index] << "," << get<1>(coords)[index]
      77             :            << "," << get<2>(coords)[index] << ")\n";
      78             :       }
      79             :     }
      80             : 
      81             :     if constexpr (Ignore) {
      82             :       const ::Verbosity verbosity = db::get<ah::Tags::Verbosity>(box);
      83             :       if (verbosity >= ::Verbosity::Quiet) {
      84             :         Parallel::printf("%s\n", os.str());
      85             :       }
      86             :     } else {
      87             :       ERROR(os.str());
      88             :     }
      89             :   }
      90             : };
      91             : }  // namespace ah::callbacks