Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <memory>
       7             : #include <optional>
       8             : 
       9             : #include "DataStructures/DataBox/Prefixes.hpp"
      10             : #include "DataStructures/Tensor/TypeAliases.hpp"
      11             : #include "Evolution/BoundaryCorrection.hpp"
      12             : #include "Evolution/Systems/GrMhd/ValenciaDivClean/Tags.hpp"
      13             : #include "NumericalAlgorithms/DiscontinuousGalerkin/Formulation.hpp"
      14             : #include "Options/String.hpp"
      15             : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
      16             : #include "PointwiseFunctions/Hydro/Tags.hpp"
      17             : #include "Utilities/Gsl.hpp"
      18             : #include "Utilities/Serialization/CharmPupable.hpp"
      19             : #include "Utilities/TMPL.hpp"
      20             : 
      21             : /// \cond
      22             : class DataVector;
      23             : namespace gsl {
      24             : template <typename T>
      25             : class not_null;
      26             : }  // namespace gsl
      27             : namespace PUP {
      28             : class er;
      29             : }  // namespace PUP
      30             : /// \endcond
      31             : 
      32             : namespace grmhd::ValenciaDivClean::BoundaryCorrections {
      33             : /*!
      34             :  * \brief A Rusanov/local Lax-Friedrichs Riemann solver
      35             :  *
      36             :  * Let \f$U\f$ be the state vector of evolved variables, \f$F^i\f$ the
      37             :  * corresponding fluxes, and \f$n_i\f$ be the outward directed unit normal to
      38             :  * the interface. Denoting \f$F := n_i F^i\f$, the %Rusanov boundary correction
      39             :  * is
      40             :  *
      41             :  * \f{align*}
      42             :  * G_\text{Rusanov} = \frac{F_\text{int} - F_\text{ext}}{2} -
      43             :  * \frac{\text{max}\left(\{|\lambda_\text{int}|\},
      44             :  * \{|\lambda_\text{ext}|\}\right)}{2} \left(U_\text{ext} - U_\text{int}\right),
      45             :  * \f}
      46             :  *
      47             :  * where "int" and "ext" stand for interior and exterior, and
      48             :  * \f$\{|\lambda|\}\f$ is the set of characteristic/signal speeds. The minus
      49             :  * sign in front of the \f$F_{\text{ext}}\f$ is necessary because the outward
      50             :  * directed normal of the neighboring element has the opposite sign, i.e.
      51             :  * \f$n_i^{\text{ext}}=-n_i^{\text{int}}\f$. The characteristic/signal speeds
      52             :  * are given in the documentation for
      53             :  * grmhd::ValenciaDivClean::characteristic_speeds(). Since the fluid is
      54             :  * travelling slower than the speed of light, the speeds we are interested in
      55             :  * are
      56             :  *
      57             :  * \f{align*}{
      58             :  *   \lambda_{\pm}&=\pm\alpha-\beta^i n_i,
      59             :  * \f}
      60             :  *
      61             :  * which correspond to the divergence cleaning field.
      62             :  *
      63             :  * \note
      64             :  * - In the strong form the `dg_boundary_terms` function returns
      65             :  * \f$G - F_\text{int}\f$
      66             :  * - It may be possible to use the slower speeds for the magnetic field and
      67             :  * fluid part of the system in order to make the flux less dissipative for
      68             :  * those variables.
      69             :  */
      70           1 : class Rusanov final : public evolution::BoundaryCorrection {
      71             :  private:
      72           0 :   struct AbsCharSpeed : db::SimpleTag {
      73           0 :     using type = Scalar<DataVector>;
      74             :   };
      75             : 
      76             :  public:
      77           0 :   using options = tmpl::list<>;
      78           0 :   static constexpr Options::String help = {
      79             :       "Computes the Rusanov or local Lax-Friedrichs boundary correction term "
      80             :       "for the GRMHD system."};
      81             : 
      82           0 :   Rusanov() = default;
      83           0 :   Rusanov(const Rusanov&) = default;
      84           0 :   Rusanov& operator=(const Rusanov&) = default;
      85           0 :   Rusanov(Rusanov&&) = default;
      86           0 :   Rusanov& operator=(Rusanov&&) = default;
      87           0 :   ~Rusanov() override = default;
      88             : 
      89             :   /// \cond
      90             :   explicit Rusanov(CkMigrateMessage* /*unused*/);
      91             :   using PUP::able::register_constructor;
      92             :   WRAPPED_PUPable_decl_template(Rusanov);  // NOLINT
      93             :   /// \endcond
      94           0 :   void pup(PUP::er& p) override;  // NOLINT
      95             : 
      96           0 :   std::unique_ptr<BoundaryCorrection> get_clone() const override;
      97             : 
      98           0 :   using dg_package_field_tags =
      99             :       tmpl::list<Tags::TildeD, Tags::TildeYe, Tags::TildeTau,
     100             :                  Tags::TildeS<Frame::Inertial>, Tags::TildeB<Frame::Inertial>,
     101             :                  Tags::TildePhi, ::Tags::NormalDotFlux<Tags::TildeD>,
     102             :                  ::Tags::NormalDotFlux<Tags::TildeYe>,
     103             :                  ::Tags::NormalDotFlux<Tags::TildeTau>,
     104             :                  ::Tags::NormalDotFlux<Tags::TildeS<Frame::Inertial>>,
     105             :                  ::Tags::NormalDotFlux<Tags::TildeB<Frame::Inertial>>,
     106             :                  ::Tags::NormalDotFlux<Tags::TildePhi>, AbsCharSpeed>;
     107           0 :   using dg_package_data_temporary_tags = tmpl::list<
     108             :       gr::Tags::Lapse<DataVector>, gr::Tags::Shift<DataVector, 3>,
     109             :       hydro::Tags::SpatialVelocityOneForm<DataVector, 3, Frame::Inertial>>;
     110           0 :   using dg_package_data_primitive_tags =
     111             :       tmpl::list<hydro::Tags::RestMassDensity<DataVector>,
     112             :                  hydro::Tags::ElectronFraction<DataVector>,
     113             :                  hydro::Tags::Temperature<DataVector>,
     114             :                  hydro::Tags::SpatialVelocity<DataVector, 3>>;
     115           0 :   using dg_package_data_volume_tags =
     116             :       tmpl::list<hydro::Tags::GrmhdEquationOfState>;
     117           0 :   using dg_boundary_terms_volume_tags = tmpl::list<>;
     118             : 
     119           0 :   static double dg_package_data(
     120             :       gsl::not_null<Scalar<DataVector>*> packaged_tilde_d,
     121             :       gsl::not_null<Scalar<DataVector>*> packaged_tilde_ye,
     122             :       gsl::not_null<Scalar<DataVector>*> packaged_tilde_tau,
     123             :       gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*> packaged_tilde_s,
     124             :       gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*> packaged_tilde_b,
     125             :       gsl::not_null<Scalar<DataVector>*> packaged_tilde_phi,
     126             :       gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_d,
     127             :       gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_ye,
     128             :       gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_tau,
     129             :       gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*>
     130             :           packaged_normal_dot_flux_tilde_s,
     131             :       gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
     132             :           packaged_normal_dot_flux_tilde_b,
     133             :       gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_phi,
     134             :       gsl::not_null<Scalar<DataVector>*> packaged_abs_char_speed,
     135             : 
     136             :       const Scalar<DataVector>& tilde_d, const Scalar<DataVector>& tilde_ye,
     137             :       const Scalar<DataVector>& tilde_tau,
     138             :       const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s,
     139             :       const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b,
     140             :       const Scalar<DataVector>& tilde_phi,
     141             : 
     142             :       const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_d,
     143             :       const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_ye,
     144             :       const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_tau,
     145             :       const tnsr::Ij<DataVector, 3, Frame::Inertial>& flux_tilde_s,
     146             :       const tnsr::IJ<DataVector, 3, Frame::Inertial>& flux_tilde_b,
     147             :       const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_phi,
     148             : 
     149             :       const Scalar<DataVector>& lapse,
     150             :       const tnsr::I<DataVector, 3, Frame::Inertial>& shift,
     151             :       const tnsr::i<DataVector, 3,
     152             :                     Frame::Inertial>& /*spatial_velocity_one_form*/,
     153             : 
     154             :       const Scalar<DataVector>& /*rest_mass_density*/,
     155             :       const Scalar<DataVector>& /*electron_fraction*/,
     156             :       const Scalar<DataVector>& /*temperature*/,
     157             :       const tnsr::I<DataVector, 3, Frame::Inertial>& /*spatial_velocity*/,
     158             : 
     159             :       const tnsr::i<DataVector, 3, Frame::Inertial>& normal_covector,
     160             :       const tnsr::I<DataVector, 3, Frame::Inertial>& normal_vector,
     161             :       const std::optional<tnsr::I<DataVector, 3, Frame::Inertial>>&
     162             :       /*mesh_velocity*/,
     163             :       const std::optional<Scalar<DataVector>>& normal_dot_mesh_velocity,
     164             :       const EquationsOfState::EquationOfState<true, 3>&
     165             :       /*equation_of_state*/);
     166             : 
     167           0 :   static void dg_boundary_terms(
     168             :       gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_d,
     169             :       gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_ye,
     170             :       gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_tau,
     171             :       gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*>
     172             :           boundary_correction_tilde_s,
     173             :       gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
     174             :           boundary_correction_tilde_b,
     175             :       gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_phi,
     176             :       const Scalar<DataVector>& tilde_d_int,
     177             :       const Scalar<DataVector>& tilde_ye_int,
     178             :       const Scalar<DataVector>& tilde_tau_int,
     179             :       const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s_int,
     180             :       const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b_int,
     181             :       const Scalar<DataVector>& tilde_phi_int,
     182             :       const Scalar<DataVector>& normal_dot_flux_tilde_d_int,
     183             :       const Scalar<DataVector>& normal_dot_flux_tilde_ye_int,
     184             :       const Scalar<DataVector>& normal_dot_flux_tilde_tau_int,
     185             :       const tnsr::i<DataVector, 3, Frame::Inertial>&
     186             :           normal_dot_flux_tilde_s_int,
     187             :       const tnsr::I<DataVector, 3, Frame::Inertial>&
     188             :           normal_dot_flux_tilde_b_int,
     189             :       const Scalar<DataVector>& normal_dot_flux_tilde_phi_int,
     190             :       const Scalar<DataVector>& abs_char_speed_int,
     191             :       const Scalar<DataVector>& tilde_d_ext,
     192             :       const Scalar<DataVector>& tilde_ye_ext,
     193             :       const Scalar<DataVector>& tilde_tau_ext,
     194             :       const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s_ext,
     195             :       const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b_ext,
     196             :       const Scalar<DataVector>& tilde_phi_ext,
     197             :       const Scalar<DataVector>& normal_dot_flux_tilde_d_ext,
     198             :       const Scalar<DataVector>& normal_dot_flux_tilde_ye_ext,
     199             :       const Scalar<DataVector>& normal_dot_flux_tilde_tau_ext,
     200             :       const tnsr::i<DataVector, 3, Frame::Inertial>&
     201             :           normal_dot_flux_tilde_s_ext,
     202             :       const tnsr::I<DataVector, 3, Frame::Inertial>&
     203             :           normal_dot_flux_tilde_b_ext,
     204             :       const Scalar<DataVector>& normal_dot_flux_tilde_phi_ext,
     205             :       const Scalar<DataVector>& abs_char_speed_ext,
     206             :       dg::Formulation dg_formulation);
     207             : };
     208             : 
     209           0 : bool operator==(const Rusanov& lhs, const Rusanov& rhs);
     210           0 : bool operator!=(const Rusanov& lhs, const Rusanov& rhs);
     211             : }  // namespace grmhd::ValenciaDivClean::BoundaryCorrections