![]() |
ryujin 2.1.1 revision 9dcb748690310d6a540ebb8b066d1a0834fc7604
|
#include <source/hyperbolic_module.h>
Public Types | |
using | HyperbolicSystem = typename Description::HyperbolicSystem |
using | HyperbolicSystemView = typename HyperbolicSystem::template View< dim, Number > |
using | state_type = typename HyperbolicSystemView::state_type |
using | flux_type = typename HyperbolicSystemView::flux_type |
using | scalar_type = typename OfflineData< dim, Number >::scalar_type |
using | vector_type = typename HyperbolicSystemView::vector_type |
using | precomputed_vector_type = typename HyperbolicSystemView::precomputed_vector_type |
using | precomputed_initial_vector_type = typename HyperbolicSystemView::precomputed_initial_vector_type |
Public Member Functions | |
HyperbolicModule (const MPI_Comm &mpi_communicator, std::map< std::string, dealii::Timer > &computing_timer, const OfflineData< dim, Number > &offline_data, const HyperbolicSystem &hyperbolic_system, const InitialValues< Description, dim, Number > &initial_values, const std::string &subsection="/HyperbolicModule") | |
void | prepare () |
Static Public Attributes | |
static constexpr unsigned int | problem_dimension |
static constexpr unsigned int | n_precomputed_values |
static constexpr unsigned int | n_precomputation_cycles |
static constexpr unsigned int | n_precomputed_initial_values |
Functons for performing explicit time steps | |
bool | precompute_only_ |
IDViolationStrategy | id_violation_strategy_ |
template<int stages> | |
Number | step (const vector_type &old_U, std::array< std::reference_wrapper< const vector_type >, stages > stage_U, std::array< std::reference_wrapper< const precomputed_vector_type >, stages > stage_precomputed, const std::array< Number, stages > stage_weights, vector_type &new_U, precomputed_vector_type &new_precomputed, Number tau=Number(0.)) const |
void | apply_boundary_conditions (vector_type &U, Number t) const |
void | cfl (Number new_cfl) const |
auto & | cfl () const |
auto & | offline_data () const |
auto & | hyperbolic_system () const |
auto & | precomputed_initial () const |
auto & | alpha () const |
auto & | n_restarts () const |
auto & | n_warnings () const |
Explicit forward Euler time-stepping for hyperbolic systems with convex limiting.
This module is described in detail in [12], Alg. 1.
Definition at line 75 of file hyperbolic_module.h.
using ryujin::HyperbolicModule< Description, dim, Number >::HyperbolicSystem = typename Description::HyperbolicSystem |
Definition at line 81 of file hyperbolic_module.h.
using ryujin::HyperbolicModule< Description, dim, Number >::HyperbolicSystemView = typename HyperbolicSystem::template View<dim, Number> |
Definition at line 86 of file hyperbolic_module.h.
using ryujin::HyperbolicModule< Description, dim, Number >::state_type = typename HyperbolicSystemView::state_type |
Definition at line 98 of file hyperbolic_module.h.
using ryujin::HyperbolicModule< Description, dim, Number >::flux_type = typename HyperbolicSystemView::flux_type |
Definition at line 103 of file hyperbolic_module.h.
using ryujin::HyperbolicModule< Description, dim, Number >::scalar_type = typename OfflineData<dim, Number>::scalar_type |
Shorthand typedef for dealii::LinearAlgebra::distributed::Vector<Number>.
Definition at line 108 of file hyperbolic_module.h.
using ryujin::HyperbolicModule< Description, dim, Number >::vector_type = typename HyperbolicSystemView::vector_type |
Definition at line 113 of file hyperbolic_module.h.
using ryujin::HyperbolicModule< Description, dim, Number >::precomputed_vector_type = typename HyperbolicSystemView::precomputed_vector_type |
Typedef for a MultiComponentVector storing precomputed values.
Definition at line 130 of file hyperbolic_module.h.
using ryujin::HyperbolicModule< Description, dim, Number >::precomputed_initial_vector_type = typename HyperbolicSystemView::precomputed_initial_vector_type |
Typedef for a MultiComponentVector storing precomputed initial_values.
Definition at line 142 of file hyperbolic_module.h.
ryujin::HyperbolicModule< Description, dim, Number >::HyperbolicModule | ( | const MPI_Comm & | mpi_communicator, |
std::map< std::string, dealii::Timer > & | computing_timer, | ||
const OfflineData< dim, Number > & | offline_data, | ||
const HyperbolicSystem & | hyperbolic_system, | ||
const InitialValues< Description, dim, Number > & | initial_values, | ||
const std::string & | subsection = "/HyperbolicModule< Description, dim, Number >" |
||
) |
void ryujin::HyperbolicModule< Description, dim, Number >::prepare |
Prepare time stepping. A call to prepare()
allocates temporary storage and is necessary before any of the following time-stepping functions can be called.
Definition at line 75 of file hyperbolic_module.template.h.
Number ryujin::HyperbolicModule< Description, dim, Number >::step | ( | const vector_type & | old_U, |
std::array< std::reference_wrapper< const vector_type >, stages > | stage_U, | ||
std::array< std::reference_wrapper< const precomputed_vector_type >, stages > | stage_precomputed, | ||
const std::array< Number, stages > | stage_weights, | ||
vector_type & | new_U, | ||
precomputed_vector_type & | new_precomputed, | ||
Number | tau = Number(0.) |
||
) | const |
Given a reference to a previous state vector old_U
perform an explicit euler step (and store the result in new_U
). The function returns the computed maximal time step size tau_max according to the CFL condition.
The time step is performed with either tau_max (if tau
is set to 0), or tau (if tau
is nonzero). Here, tau_max is the computed maximal time step size and tau
is the last parameter of the function.
The function takes an optional array of states stage_U
together with a an array of weights stage_weights
to construct a modified high-order flux. The standard high-order flux reads (cf [12], Eq. 12):
\begin{align} \newcommand{\bF}{{\boldsymbol F}} \newcommand{\bU}{{\boldsymbol U}} \newcommand\bUni{\bU^n_i} \newcommand\bUnj{\bU^n_j} \newcommand{\polf}{{\mathbb f}} \newcommand\Ii{\mathcal{I}(i)} \newcommand{\bc}{{\boldsymbol c}} \sum_{j\in\Ii} \frac{m_{ij}}{m_{j}} \; \frac{m_{j}}{\tau_n}\big( \tilde\bU_j^{H,n+1} - \bU_j^{n}\big) \;=\; \bF^n_i + \sum_{j\in\Ii}d_{ij}^{H,n}\big(\bUnj-\bUni\big), \qquad\text{with}\quad \bF^n_i\;:=\; \sum_{j\in\Ii}\Big(-(\polf(\bUni)+\polf(\bUnj)) \cdot\bc_{ij}\Big). \end{align}
Instead, the function assembles the modified high-order flux:
\begin{align} \newcommand{\bF}{{\boldsymbol F}} \newcommand{\bU}{{\boldsymbol U}} \newcommand\bUnis{\bU^{s,n}_i} \newcommand\bUnjs{\bU^{s,n}_j} \newcommand{\polf}{{\mathbb f}} \newcommand\Ii{\mathcal{I}(i)} \newcommand{\bc}{{\boldsymbol c}} \tilde{\bF}^n_i\;:=\; \big(1-\sum_{s=\{1:\text{stages}\}}\omega_s\big)\bF^n_i \;+\; \sum_{s=\{1:stages\}}\omega_s \bF^{s,n}_i \qquad\text{with}\quad \bF^{s,n}_i\;:=\; \sum_{j\in\Ii}\Big(-(\polf(\bUnis)+\polf(\bUnjs)) \cdot\bc_{ij}\Big). \end{align}
where \(\omega_s\) denotes the weigths for the given stages \(\bU^{s,n}\).
new_U
. It is best to simply call HyperbolicModule::apply_boundary_conditions() on the appropriate vector immediately after performing a time step. Definition at line 150 of file hyperbolic_module.template.h.
void ryujin::HyperbolicModule< Description, dim, Number >::apply_boundary_conditions | ( | vector_type & | U, |
Number | t | ||
) | const |
This function postprocesses a given state U
to conform with all prescribed boundary conditions at time t
. This implies that on slip (and no-slip) boundaries the normal momentum is set to zero; on Dirichlet boundaries the appropriate state at time t
is substituted; and on "flexible" boundaries depending on the fact whether we have supersonic or subsonic inflow/outflow the appropriate Riemann invariant is prescribed. See [11] for details.
U
Definition at line 1012 of file hyperbolic_module.template.h.
References ryujin::do_nothing, LIKWID_MARKER_START, and LIKWID_MARKER_STOP.
|
inline |
Sets the relative CFL number used for computing an appropriate time-step size to the given value. The CFL number must be a positive value. If chosen to be within the interval \((0,1)\) then the low-order update and limiting stages guarantee invariant domain preservation.
Definition at line 260 of file hyperbolic_module.h.
|
inline |
Returns the relative CFL number used for computing an appropriate time-step size.
Definition at line 270 of file hyperbolic_module.h.
|
inline |
Return a reference to the OfflineData object
Definition at line 275 of file hyperbolic_module.h.
|
inline |
Return a reference to the HyperbolicSystem object
Definition at line 280 of file hyperbolic_module.h.
|
inline |
Return a reference to the precomputed initial data vector
Definition at line 285 of file hyperbolic_module.h.
|
inline |
Return a reference to alpha vector storing indicator values. Note that the values stored in alpha correspond to the last step executed by this class. This value can be recomputed for a given state vector by setting precompute_only_ to true and calling the step() function.
Definition at line 293 of file hyperbolic_module.h.
|
inline |
The number of restarts issued by the step() function.
Definition at line 298 of file hyperbolic_module.h.
|
inline |
The number of ID violation warnings encounterd in the step() function.
Definition at line 304 of file hyperbolic_module.h.
|
staticconstexpr |
Definition at line 92 of file hyperbolic_module.h.
|
staticconstexpr |
Definition at line 118 of file hyperbolic_module.h.
|
staticconstexpr |
Definition at line 124 of file hyperbolic_module.h.
|
staticconstexpr |
Definition at line 136 of file hyperbolic_module.h.
|
mutable |
Definition at line 307 of file hyperbolic_module.h.
|
mutable |
Definition at line 310 of file hyperbolic_module.h.