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 <memory>
       9             : #include <string>
      10             : 
      11             : #include "DataStructures/Tensor/Tensor.hpp"
      12             : #include "IO/Exporter/PointwiseInterpolator.hpp"
      13             : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
      14             : #include "Options/String.hpp"
      15             : #include "PointwiseFunctions/InitialDataUtilities/Background.hpp"
      16             : #include "PointwiseFunctions/InitialDataUtilities/InitialData.hpp"
      17             : #include "PointwiseFunctions/InitialDataUtilities/InitialGuess.hpp"
      18             : #include "Utilities/Serialization/CharmPupable.hpp"
      19             : #include "Utilities/TMPL.hpp"
      20             : #include "Utilities/TaggedTuple.hpp"
      21             : 
      22             : /// \cond
      23             : namespace PUP {
      24             : class er;
      25             : }  // namespace PUP
      26             : /// \endcond
      27             : 
      28             : /*!
      29             :  * \brief Load numeric data from volume data files
      30             :  *
      31             :  * This class loads all requested tags from volume data files using
      32             :  * `spectre::Exporter::interpolate_to_points`. This is an easy and useful
      33             :  * alternative to `importers::ElementDataReader` to load numeric data with the
      34             :  * following advantages:
      35             :  *
      36             :  * - No need to work with parallel components, actions, or phases to load
      37             :  *   numeric data. Also, no need to parse numeric data options separately from
      38             :  *   analytic data options. Just include this class in the list of analytic data
      39             :  *   classes for the executable.
      40             :  * - The data is interpolated and returned in serial, just like analytic data.
      41             :  * - The data can be re-interpolated to any set of points on request, which is
      42             :  *   an easy way to handle AMR or other domain changes.
      43             :  *
      44             :  * However, it also comes with the following caveats:
      45             :  *
      46             :  * - The volume data files must have datasets with the same names as the tags
      47             :  *   being loaded. For example, if the tag being loaded is `gr::Tags::Shift`,
      48             :  *   then the volume data files must have datasets named `Shift_x`, `Shift_y`,
      49             :  *   and `Shift_z`. If you need processing of the data before it can be used,
      50             :  *   e.g. to load some datasets and compute other derived quantities from them,
      51             :  *   then consider writing a custom class or use `importers::ElementDataReader`.
      52             :  * - Reading in data on the same grid that it was written on may not give you
      53             :  *   the same data back, at least not on Gauss-Lobatto grids. This is because
      54             :  *   grid points on element boundaries are disambiguated by
      55             :  *   `block_logical_coordinates` and `element_logical_coordinates`, so a
      56             :  *   boundary point may be written by one element but read back in from its
      57             :  *   neighbor. To avoid this, consider using `importers::ElementDataReader`
      58             :  *   which has functionality to read in data on the same grid that it was
      59             :  *   written on. This is not possible with this class because it is not aware of
      60             :  *   the element structure (it operates in a pointwise manner).
      61             :  * - Large datasets may not fit in memory. Each element will open the H5 files
      62             :  *   and interpolate data from them, so this may not fit in memory if the
      63             :  *   datasets are large and/or many elements are doing this at the same time. To
      64             :  *   avoid this, consider using `importers::ElementDataReader` which reads one
      65             :  *   H5 file at a time on the node level.
      66             :  */
      67           1 : class NumericData {
      68             :  public:
      69           0 :   struct FileGlob {
      70           0 :     static constexpr Options::String help =
      71             :         "Path or glob pattern to the data file";
      72           0 :     using type = std::string;
      73             :   };
      74           0 :   struct Subgroup {
      75           0 :     static constexpr Options::String help = {
      76             :         "The subgroup within the file, excluding extensions"};
      77           0 :     using type = std::string;
      78             :   };
      79           0 :   struct ObservationStep {
      80           0 :     static constexpr Options::String help =
      81             :         "The observation step at which to read data";
      82           0 :     using type = int;
      83             :   };
      84           0 :   struct ExtrapolateIntoExcisions {
      85           0 :     static constexpr Options::String help = {
      86             :         "Whether to extrapolate data into excised regions"};
      87           0 :     using type = bool;
      88             :   };
      89           0 :   using options =
      90             :       tmpl::list<FileGlob, Subgroup, ObservationStep, ExtrapolateIntoExcisions>;
      91           0 :   static constexpr Options::String help =
      92             :       "Numeric data loaded from volume data files";
      93             : 
      94           0 :   NumericData() = default;
      95           0 :   NumericData(const NumericData&) = default;
      96           0 :   NumericData& operator=(const NumericData&) = default;
      97           0 :   NumericData(NumericData&&) = default;
      98           0 :   NumericData& operator=(NumericData&&) = default;
      99           0 :   ~NumericData() = default;
     100             : 
     101           0 :   NumericData(std::string file_glob, std::string subgroup, int observation_step,
     102             :               bool extrapolate_into_excisions);
     103             : 
     104           0 :   const std::string& file_glob() const { return file_glob_; }
     105             : 
     106           0 :   const std::string& subgroup() const { return subgroup_; }
     107             : 
     108           0 :   int observation_step() const { return observation_step_; }
     109             : 
     110           0 :   bool extrapolate_into_excisions() const {
     111             :     return extrapolate_into_excisions_;
     112             :   }
     113             : 
     114             :   template <typename DataType, size_t Dim, typename... RequestedTags>
     115           0 :   tuples::TaggedTuple<RequestedTags...> variables(
     116             :       const tnsr::I<DataType, Dim>& x,
     117             :       tmpl::list<RequestedTags...> /*meta*/) const {
     118             :     return spectre::Exporter::interpolate_to_points<
     119             :         tmpl::list<RequestedTags...>>(
     120             :         file_glob_, subgroup_,
     121             :         spectre::Exporter::ObservationStep{observation_step_}, x,
     122             :         extrapolate_into_excisions_);
     123             :   }
     124             : 
     125             :   template <size_t Dim, typename... RequestedTags>
     126           0 :   tuples::TaggedTuple<RequestedTags...> variables(
     127             :       const tnsr::I<DataVector, Dim>& x, const Mesh<Dim>& /*mesh*/,
     128             :       const InverseJacobian<DataVector, Dim, Frame::ElementLogical,
     129             :                             Frame::Inertial>& /*inv_jacobian*/,
     130             :       tmpl::list<RequestedTags...> /*meta*/) const {
     131             :     return variables(x, tmpl::list<RequestedTags...>{});
     132             :   }
     133             : 
     134             :   // NOLINTNEXTLINE
     135           0 :   void pup(PUP::er& p);
     136             : 
     137             :  protected:
     138           0 :   std::string file_glob_{};
     139           0 :   std::string subgroup_{};
     140           0 :   int observation_step_{};
     141           0 :   bool extrapolate_into_excisions_{};
     142             : };
     143             : 
     144           0 : bool operator==(const NumericData& lhs, const NumericData& rhs);
     145           0 : bool operator!=(const NumericData& lhs, const NumericData& rhs);
     146             : 
     147             : namespace elliptic::analytic_data {
     148             : 
     149             : /*!
     150             :  * \brief Load numeric data from volume data files
     151             :  *
     152             :  * \see ::NumericData
     153             :  */
     154           1 : class NumericData : public elliptic::analytic_data::Background,
     155             :                     public elliptic::analytic_data::InitialGuess,
     156             :                     public ::NumericData {
     157             :  public:
     158             :   using ::NumericData::NumericData;
     159             : 
     160           0 :   explicit NumericData(CkMigrateMessage* m);
     161             :   using PUP::able::register_constructor;
     162           0 :   WRAPPED_PUPable_decl_template(NumericData);
     163             : 
     164             :   // NOLINTNEXTLINE
     165           0 :   void pup(PUP::er& p) override;
     166             : };
     167             : 
     168           0 : bool operator==(const NumericData& lhs, const NumericData& rhs);
     169           0 : bool operator!=(const NumericData& lhs, const NumericData& rhs);
     170             : 
     171             : }  // namespace elliptic::analytic_data
     172             : 
     173             : namespace evolution::initial_data {
     174             : 
     175             : /*!
     176             :  * \brief Load numeric data from volume data files
     177             :  *
     178             :  * \see ::NumericData
     179             :  */
     180           1 : class NumericData : public evolution::initial_data::InitialData,
     181             :                     public ::NumericData {
     182             :  public:
     183             :   using ::NumericData::NumericData;
     184             : 
     185           0 :   explicit NumericData(CkMigrateMessage* m);
     186             :   using PUP::able::register_constructor;
     187           0 :   WRAPPED_PUPable_decl_template(NumericData);
     188             : 
     189           0 :   std::unique_ptr<evolution::initial_data::InitialData> get_clone()
     190             :       const override;
     191             : 
     192             :   // NOLINTNEXTLINE
     193           0 :   void pup(PUP::er& p) override;
     194             : };
     195             : 
     196           0 : bool operator==(const NumericData& lhs, const NumericData& rhs);
     197           0 : bool operator!=(const NumericData& lhs, const NumericData& rhs);
     198             : 
     199             : }  // namespace evolution::initial_data