Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <cstddef>
       7             : #include <limits>
       8             : 
       9             : #include "DataStructures/Tensor/TypeAliases.hpp"
      10             : #include "DataStructures/VariablesTag.hpp"
      11             : #include "Domain/Tags.hpp"
      12             : #include "Evolution/Systems/GrMhd/ValenciaDivClean/PrimitiveFromConservativeOptions.hpp"
      13             : #include "Evolution/Systems/GrMhd/ValenciaDivClean/Tags.hpp"
      14             : #include "Options/String.hpp"
      15             : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
      16             : #include "PointwiseFunctions/Hydro/Tags.hpp"
      17             : #include "Utilities/TMPL.hpp"
      18             : 
      19             : /// \cond
      20             : class DataVector;
      21             : namespace EquationsOfState {
      22             : template <bool IsRelativistic, size_t ThermodynamicDim>
      23             : class EquationOfState;
      24             : }  // namespace EquationsOfState
      25             : namespace PUP {
      26             : class er;
      27             : }  // namespace PUP
      28             : template <size_t Dim>
      29             : class Mesh;
      30             : namespace gsl {
      31             : template <typename T>
      32             : class not_null;
      33             : }  // namespace gsl
      34             : template <typename TagsList>
      35             : class Variables;
      36             : /// \endcond
      37             : 
      38             : namespace grmhd::ValenciaDivClean {
      39             : /*!
      40             :  * \brief Reduces oscillations inside an element in an attempt to guarantee a
      41             :  * physical solution of the conserved variables for which the primitive
      42             :  * variables can be recovered.
      43             :  *
      44             :  * The algorithm uses the conditions of FixConservatives on \f$\tilde{D}\f$ and
      45             :  * \f$\tilde{\tau}\f$ to reduce oscillations inside an element. Oscillations are
      46             :  * reduced by rescaling the conserved variables about the mean to bring them
      47             :  * into the required range. When rescaling \f$\tilde{D}\f$ because it is
      48             :  * negative, it is important to also rescale \f$\tilde{\tau}\f$ and
      49             :  * \f$\tilde{S}_i\f$ by the same amount. At least, this is what is observed in
      50             :  * the cylindrical blast wave test problem.
      51             :  *
      52             :  * This currently doesn't use the check on \f$\tilde{S}^2\f$, but instead checks
      53             :  * that the primitive variables can be recovered. If the primitives cannot be
      54             :  * recovered then we flatten to the mean values in the element.
      55             :  */
      56             : template <typename RecoverySchemesList>
      57           1 : class Flattener {
      58             :  public:
      59             :   /// \brief Require that the mean of TildeD is positive, otherwise terminate
      60             :   /// the simulation.
      61           1 :   struct RequirePositiveMeanTildeD {
      62           0 :     using type = bool;
      63           0 :     static constexpr Options::String help = {
      64             :         "Require that the mean of TildeD is positive, otherwise terminate the "
      65             :         "simulation."};
      66             :   };
      67             : 
      68             :   /// \brief Require that the mean of TildeYe is positive, otherwise terminate
      69             :   /// the simulation.
      70           1 :   struct RequirePositiveMeanTildeYe {
      71           0 :     using type = bool;
      72           0 :     static constexpr Options::String help = {
      73             :         "Require that the mean of TildeYe is positive, otherwise terminate the "
      74             :         "simulation."};
      75             :   };
      76             : 
      77             :   /// \brief Require that the mean of TildeTau is physical, otherwise terminate
      78             :   /// the simulation.
      79           1 :   struct RequirePhysicalMeanTildeTau {
      80           0 :     using type = bool;
      81           0 :     static constexpr Options::String help = {
      82             :         "Require that the mean of TildeTau is physical, otherwise terminate "
      83             :         "the simulation."};
      84             :   };
      85             : 
      86             :   /// \brief If true, then the primitive variables are updated at the end of the
      87             :   /// function.
      88             :   ///
      89             :   /// This is useful when not using any pointwise conserved variable fixing,
      90             :   /// treating the case that the means do not satisfy the bounds as an error.
      91           1 :   struct RecoverPrimitives {
      92           0 :     using type = bool;
      93           0 :     static constexpr Options::String help = {
      94             :         "If true, then the primitive variables are updated at the end of the "
      95             :         "function."};
      96             :   };
      97             : 
      98           0 :   using options =
      99             :       tmpl::list<RequirePositiveMeanTildeD, RequirePositiveMeanTildeYe,
     100             :                  RequirePhysicalMeanTildeTau, RecoverPrimitives>;
     101           0 :   static constexpr Options::String help = {
     102             :       "Reduces oscillations (flattens) the conserved variables according to "
     103             :       "the variable fixing procedure described in Foucart's thesis.\n"};
     104             : 
     105           0 :   Flattener(bool require_positive_mean_tilde_d,
     106             :             bool require_positive_mean_tilde_ye,
     107             :             bool require_physical_mean_tilde_tau, bool recover_primitives);
     108             : 
     109           0 :   Flattener() = default;
     110           0 :   Flattener(const Flattener& /*rhs*/) = default;
     111           0 :   Flattener& operator=(const Flattener& /*rhs*/) = default;
     112           0 :   Flattener(Flattener&& /*rhs*/) = default;
     113           0 :   Flattener& operator=(Flattener&& /*rhs*/) = default;
     114           0 :   ~Flattener() = default;
     115             : 
     116             :   // NOLINTNEXTLINE(google-runtime-references)
     117           0 :   void pup(PUP::er& p);
     118             : 
     119           0 :   using return_tags =
     120             :       tmpl::list<Tags::TildeD, Tags::TildeYe, Tags::TildeTau,
     121             :                  Tags::TildeS<Frame::Inertial>,
     122             :                  ::Tags::Variables<hydro::grmhd_tags<DataVector>>>;
     123           0 :   using argument_tags = tmpl::list<
     124             :       Tags::TildeB<>, Tags::TildePhi,
     125             :       gr::Tags::SqrtDetSpatialMetric<DataVector>,
     126             :       gr::Tags::SpatialMetric<DataVector, 3>,
     127             :       gr::Tags::InverseSpatialMetric<DataVector, 3>, domain::Tags::Mesh<3>,
     128             :       domain::Tags::DetInvJacobian<Frame::ElementLogical, Frame::Inertial>,
     129             :       hydro::Tags::GrmhdEquationOfState,
     130             :       grmhd::ValenciaDivClean::Tags::PrimitiveFromConservativeOptions>;
     131             : 
     132           0 :   void operator()(
     133             :       gsl::not_null<Scalar<DataVector>*> tilde_d,
     134             :       gsl::not_null<Scalar<DataVector>*> tilde_ye,
     135             :       gsl::not_null<Scalar<DataVector>*> tilde_tau,
     136             :       gsl::not_null<tnsr::i<DataVector, 3>*> tilde_s,
     137             :       gsl::not_null<Variables<hydro::grmhd_tags<DataVector>>*> primitives,
     138             :       const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b,
     139             :       const Scalar<DataVector>& tilde_phi,
     140             :       const Scalar<DataVector>& sqrt_det_spatial_metric,
     141             :       const tnsr::ii<DataVector, 3, Frame::Inertial>& spatial_metric,
     142             :       const tnsr::II<DataVector, 3, Frame::Inertial>& inv_spatial_metric,
     143             :       const Mesh<3>& mesh,
     144             :       const Scalar<DataVector>& det_logical_to_inertial_inv_jacobian,
     145             :       const EquationsOfState::EquationOfState<true, 3>& eos,
     146             :       const grmhd::ValenciaDivClean::PrimitiveFromConservativeOptions&
     147             :           primitive_from_conservative_options) const;
     148             : 
     149             :  private:
     150             :   template <typename LocalRecoverySchemesList>
     151             :   // NOLINTNEXTLINE(readability-redundant-declaration)
     152           0 :   friend bool operator==(const Flattener<LocalRecoverySchemesList>& lhs,
     153             :                          const Flattener<LocalRecoverySchemesList>& rhs);
     154             : 
     155           0 :   bool require_positive_mean_tilde_d_ = false;
     156           0 :   bool require_positive_mean_tilde_ye_ = false;
     157           0 :   bool require_physical_mean_tilde_tau_ = false;
     158           0 :   bool recover_primitives_ = false;
     159             : };
     160             : 
     161             : template <typename RecoverySchemesList>
     162           0 : bool operator!=(const Flattener<RecoverySchemesList>& lhs,
     163             :                 const Flattener<RecoverySchemesList>& rhs);
     164             : }  // namespace grmhd::ValenciaDivClean