SpECTRE
v2024.05.11
|
A Rusanov/local Lax-Friedrichs Riemann solver. More...
#include <Rusanov.hpp>
Public Types | |
using | options = tmpl::list<> |
using | dg_package_field_tags = tmpl::list< Tags::TildeE, Tags::TildeB, Tags::TildePsi, Tags::TildePhi, Tags::TildeQ, ::Tags::NormalDotFlux< Tags::TildeE >, ::Tags::NormalDotFlux< Tags::TildeB >, ::Tags::NormalDotFlux< Tags::TildePsi >, ::Tags::NormalDotFlux< Tags::TildePhi >, ::Tags::NormalDotFlux< Tags::TildeQ >, AbsCharSpeed > |
using | dg_package_data_temporary_tags = tmpl::list< gr::Tags::Lapse< DataVector >, gr::Tags::Shift< DataVector, 3 > > |
using | dg_package_data_primitive_tags = tmpl::list<> |
using | dg_package_data_volume_tags = tmpl::list<> |
using | dg_boundary_terms_volume_tags = tmpl::list<> |
Public Types inherited from ForceFree::BoundaryCorrections::BoundaryCorrection | |
using | creatable_classes = tmpl::list< Rusanov > |
Public Member Functions | |
Rusanov (const Rusanov &)=default | |
Rusanov & | operator= (const Rusanov &)=default |
Rusanov (Rusanov &&)=default | |
Rusanov & | operator= (Rusanov &&)=default |
void | pup (PUP::er &p) override |
std::unique_ptr< BoundaryCorrection > | get_clone () const override |
Public Member Functions inherited from ForceFree::BoundaryCorrections::BoundaryCorrection | |
BoundaryCorrection (const BoundaryCorrection &)=default | |
BoundaryCorrection & | operator= (const BoundaryCorrection &)=default |
BoundaryCorrection (BoundaryCorrection &&)=default | |
BoundaryCorrection & | operator= (BoundaryCorrection &&)=default |
virtual std::unique_ptr< BoundaryCorrection > | get_clone () const =0 |
Static Public Member Functions | |
static double | dg_package_data (gsl::not_null< tnsr::I< DataVector, 3, Frame::Inertial > * > packaged_tilde_e, gsl::not_null< tnsr::I< DataVector, 3, Frame::Inertial > * > packaged_tilde_b, gsl::not_null< Scalar< DataVector > * > packaged_tilde_psi, gsl::not_null< Scalar< DataVector > * > packaged_tilde_phi, gsl::not_null< Scalar< DataVector > * > packaged_tilde_q, gsl::not_null< tnsr::I< DataVector, 3, Frame::Inertial > * > packaged_normal_dot_flux_tilde_e, gsl::not_null< tnsr::I< DataVector, 3, Frame::Inertial > * > packaged_normal_dot_flux_tilde_b, gsl::not_null< Scalar< DataVector > * > packaged_normal_dot_flux_tilde_psi, gsl::not_null< Scalar< DataVector > * > packaged_normal_dot_flux_tilde_phi, gsl::not_null< Scalar< DataVector > * > packaged_normal_dot_flux_tilde_q, gsl::not_null< Scalar< DataVector > * > packaged_abs_char_speed, const tnsr::I< DataVector, 3, Frame::Inertial > &tilde_e, const tnsr::I< DataVector, 3, Frame::Inertial > &tilde_b, const Scalar< DataVector > &tilde_psi, const Scalar< DataVector > &tilde_phi, const Scalar< DataVector > &tilde_q, const tnsr::IJ< DataVector, 3, Frame::Inertial > &flux_tilde_e, const tnsr::IJ< DataVector, 3, Frame::Inertial > &flux_tilde_b, const tnsr::I< DataVector, 3, Frame::Inertial > &flux_tilde_psi, const tnsr::I< DataVector, 3, Frame::Inertial > &flux_tilde_phi, const tnsr::I< DataVector, 3, Frame::Inertial > &flux_tilde_q, const Scalar< DataVector > &lapse, const tnsr::I< DataVector, 3, Frame::Inertial > &shift, const tnsr::i< DataVector, 3, Frame::Inertial > &normal_covector, const tnsr::I< DataVector, 3, Frame::Inertial > &normal_vector, const std::optional< tnsr::I< DataVector, 3, Frame::Inertial > > &, const std::optional< Scalar< DataVector > > &normal_dot_mesh_velocity) |
static void | dg_boundary_terms (gsl::not_null< tnsr::I< DataVector, 3, Frame::Inertial > * > boundary_correction_tilde_e, gsl::not_null< tnsr::I< DataVector, 3, Frame::Inertial > * > boundary_correction_tilde_b, gsl::not_null< Scalar< DataVector > * > boundary_correction_tilde_psi, gsl::not_null< Scalar< DataVector > * > boundary_correction_tilde_phi, gsl::not_null< Scalar< DataVector > * > boundary_correction_tilde_q, const tnsr::I< DataVector, 3, Frame::Inertial > &tilde_e_int, const tnsr::I< DataVector, 3, Frame::Inertial > &tilde_b_int, const Scalar< DataVector > &tilde_psi_int, const Scalar< DataVector > &tilde_phi_int, const Scalar< DataVector > &tilde_q_int, const tnsr::I< DataVector, 3, Frame::Inertial > &normal_dot_flux_tilde_e_int, const tnsr::I< DataVector, 3, Frame::Inertial > &normal_dot_flux_tilde_b_int, const Scalar< DataVector > &normal_dot_flux_tilde_psi_int, const Scalar< DataVector > &normal_dot_flux_tilde_phi_int, const Scalar< DataVector > &normal_dot_flux_tilde_q_int, const Scalar< DataVector > &abs_char_speed_int, const tnsr::I< DataVector, 3, Frame::Inertial > &tilde_e_ext, const tnsr::I< DataVector, 3, Frame::Inertial > &tilde_b_ext, const Scalar< DataVector > &tilde_psi_ext, const Scalar< DataVector > &tilde_phi_ext, const Scalar< DataVector > &tilde_q_ext, const tnsr::I< DataVector, 3, Frame::Inertial > &normal_dot_flux_tilde_e_ext, const tnsr::I< DataVector, 3, Frame::Inertial > &normal_dot_flux_tilde_b_ext, const Scalar< DataVector > &normal_dot_flux_tilde_psi_ext, const Scalar< DataVector > &normal_dot_flux_tilde_phi_ext, const Scalar< DataVector > &normal_dot_flux_tilde_q_ext, const Scalar< DataVector > &abs_char_speed_ext, dg::Formulation dg_formulation) |
Static Public Attributes | |
static constexpr Options::String | help |
A Rusanov/local Lax-Friedrichs Riemann solver.
Let \(U\) be the evolved variables, \(F^i\) the corresponding fluxes, and \(n_i\) be the outward directed unit normal to the interface. Denoting \(F := n_i F^i\), the Rusanov boundary correction is
\begin{align*} G_\text{Rusanov} = \frac{F_\text{int} - F_\text{ext}}{2} - \frac{\text{max}\left(|\lambda_\text{int}|, |\lambda_\text{ext}|\right)}{2} \left(U_\text{ext} - U_\text{int}\right), \end{align*}
where "int" and "ext" stand for interior and exterior, and \(\lambda\) is the characteristic/signal speed. The minus sign in front of the \(F_{\text{ext}}\) is necessary because the outward directed normal of the neighboring element has the opposite sign, i.e. \(n_i^{\text{ext}}=-n_i^{\text{int}}\).
For the GRFFE system the largest characteristic speeds \(\lambda\) of our interest are given as
\begin{align*} \lambda_{\pm} = -\beta^i n_i \pm \alpha. \end{align*}
which correspond to fast mode waves.
dg_boundary_terms
function returns \(G - F_\text{int}\)
|
overridevirtual |
|
staticconstexpr |