![]() |
Reference documentation for deal.II version 9.6.0
|
#include <deal.II/base/time_stepping.h>
Classes | |
struct | Status |
Public Member Functions | |
ImplicitRungeKutta ()=default | |
ImplicitRungeKutta (const runge_kutta_method method, const unsigned int max_it=100, const double tolerance=1e-6) | |
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 |
void | set_newton_solver_parameters (const unsigned int max_it, const double tolerance) |
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_stages (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, std::vector< VectorType > &f_stages) |
void | newton_solve (const std::function< void(const VectorType &, VectorType &)> &get_residual, const std::function< VectorType(const VectorType &)> &id_minus_tau_J_inverse, VectorType &y) |
void | compute_residual (const std::function< VectorType(const double, const VectorType &)> &f, double t, double delta_t, const VectorType &new_y, const VectorType &y, VectorType &tendency, VectorType &residual) const |
Private Attributes | |
unsigned int | max_it |
double | tolerance |
Status | status |
This class is derived from RungeKutta and implement the implicit methods. This class works only for Diagonal Implicit Runge-Kutta (DIRK) methods.
Definition at line 529 of file time_stepping.h.
|
default |
Default constructor. initialize(runge_kutta_method) and set_newton_solver_parameters(unsigned int,double) need to be called before the object can be used.
TimeStepping::ImplicitRungeKutta< VectorType >::ImplicitRungeKutta | ( | const runge_kutta_method | method, |
const unsigned int | max_it = 100, | ||
const double | tolerance = 1e-6 ) |
Constructor. This function calls initialize(runge_kutta_method) and initialize the maximum number of iterations and the tolerance of the Newton solver.
|
overridevirtual |
Initialize the implicit 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
(I-\tau J)^{-1} where I is the identity matrix, \tau is given, and J is the Jacobian \frac{\partial f}{\partial y} . The input parameters this function receives 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 >.
void TimeStepping::ImplicitRungeKutta< VectorType >::set_newton_solver_parameters | ( | const unsigned int | max_it, |
const double | tolerance ) |
Set the maximum number of iterations and the tolerance used by the Newton solver.
|
overridevirtual |
Return the status of the current object.
Implements TimeStepping::TimeStepping< VectorType >.
|
private |
Compute the different stages needed.
|
private |
Newton solver used for the implicit stages.
|
private |
Compute the residual needed by the Newton solver.
|
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 |
Maximum number of iterations of the Newton solver.
Definition at line 649 of file time_stepping.h.
|
private |
Tolerance of the Newton solver.
Definition at line 654 of file time_stepping.h.
|
private |
Status structure of the object.
Definition at line 659 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.