![]() |
deal.II version GIT relicensing-2838-gd85d4b70e9 2025-03-13 22:40:00+00:00
|
#include <deal.II/base/time_stepping.h>
Classes | |
struct | Status |
Public Member Functions | |
LowStorageRungeKutta ()=default | |
LowStorageRungeKutta (const runge_kutta_method method) | |
void | initialize (const runge_kutta_method method) override |
double | evolve_one_time_step (const std::function< VectorType(const double, const VectorType &)> &f, const std::function< VectorType(const double, const double, const VectorType &)> &id_minus_tau_J_inverse, double t, double delta_t, VectorType &y) override |
double | evolve_one_time_step (const std::function< VectorType(const double, const VectorType &)> &f, double t, double delta_t, VectorType &solution, VectorType &vec_ri, VectorType &vec_ki) |
void | get_coefficients (std::vector< double > &a, std::vector< double > &b, std::vector< double > &c) const |
const Status & | get_status () const override |
double | evolve_one_time_step (std::vector< std::function< VectorType(const double, const VectorType &)> > &F, std::vector< std::function< VectorType(const double, const double, const VectorType &)> > &J_inverse, double t, double delta_t, VectorType &y) override |
Protected Attributes | |
unsigned int | n_stages |
std::vector< double > | b |
std::vector< double > | c |
std::vector< std::vector< double > > | a |
Private Member Functions | |
void | compute_one_stage (const std::function< VectorType(const double, const VectorType &)> &f, const double t, const double factor_solution, const double factor_ai, const VectorType ¤t_ri, VectorType &vec_ki, VectorType &solution, VectorType &next_ri) const |
Private Attributes | |
Status | status |
The LowStorageRungeKutta class is derived from RungeKutta and implements a specific class of explicit methods. The main advantages of low-storage methods are the reduced memory consumption and the reduced memory access.
Definition at line 410 of file time_stepping.h.
|
default |
Default constructor. This constructor creates an object for which you will want to call initialize(runge_kutta_method)
before it can be used.
TimeStepping::LowStorageRungeKutta< VectorType >::LowStorageRungeKutta | ( | const runge_kutta_method | method | ) |
Constructor. This function calls initialize(runge_kutta_method).
|
overridevirtual |
Initialize the explicit Runge-Kutta method.
Implements TimeStepping::RungeKutta< VectorType >.
|
overridevirtual |
This function is used to advance from time t
to t+ delta_t
. f
is the function \( f(t,y) \) that should be integrated, the input parameters are the time t and the vector y and the output is value of f at this point. id_minus_tau_J_inverse
is a function that computes \(
inv(I-\tau J)\) where \( I \) is the identity matrix, \( \tau \) is given, and \( J \) is the Jacobian \( \frac{\partial f}{\partial y} \). The input parameters are the time, \( \tau \), and a vector. The output is the value of function at this point. evolve_one_time_step returns the time at the end of the time step.
Implements TimeStepping::RungeKutta< VectorType >.
double TimeStepping::LowStorageRungeKutta< VectorType >::evolve_one_time_step | ( | const std::function< VectorType(const double, const VectorType &)> & | f, |
double | t, | ||
double | delta_t, | ||
VectorType & | solution, | ||
VectorType & | vec_ri, | ||
VectorType & | vec_ki | ||
) |
This function is used to advance from time t
to t+ delta_t
. This function is similar to the one derived from RungeKutta, but does not required id_minus_tau_J_inverse because it is not used for explicit methods. evolve_one_time_step returns the time at the end of the time step. Note that vec_ki holds the evaluation of the differential operator, and vec_ri holds the right-hand side for the differential operator application.
void TimeStepping::LowStorageRungeKutta< VectorType >::get_coefficients | ( | std::vector< double > & | a, |
std::vector< double > & | b, | ||
std::vector< double > & | c | ||
) | const |
Get the coefficients of the scheme. Note that here vector a
is not the conventional definition in terms of a Butcher tableau but merely one of the sub-diagonals. More details can be found in step-67 and the references therein.
|
overridevirtual |
Return the status of the current object.
Implements TimeStepping::TimeStepping< VectorType >.
|
private |
Compute one stage of low storage rk.
|
overridevirtualinherited |
This function is used to advance from time t
to t+ delta_t
. F
is a vector of functions \( f(t,y) \) that should be integrated, the input parameters are the time t and the vector y and the output is value of f at this point. J_inverse
is a vector functions that compute the inverse of the Jacobians associated to the implicit problems. The input parameters are the time, \( \tau \), and a vector. The output is the value of function at this point. This function returns the time at the end of the time step. When using Runge-Kutta methods, F
and J_inverse
can only contain one element.
Implements TimeStepping::TimeStepping< VectorType >.
|
private |
Status structure of the object.
Definition at line 519 of file time_stepping.h.
|
protectedinherited |
Number of stages of the Runge-Kutta method.
Definition at line 281 of file time_stepping.h.
|
protectedinherited |
Butcher tableau coefficients.
Definition at line 286 of file time_stepping.h.
|
protectedinherited |
Butcher tableau coefficients.
Definition at line 291 of file time_stepping.h.
|
protectedinherited |
Butcher tableau coefficients.
Definition at line 296 of file time_stepping.h.