Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
time_stepping.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2014 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_time_stepping_h
17 #define dealii_time_stepping_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 
24 #include <functional>
25 #include <vector>
26 
28 
36 namespace TimeStepping
37 {
58  {
72  };
73 
74 
75 
83  {
87  };
88 
89 
90 
95  template <typename VectorType>
97  {
98  public:
102  virtual ~TimeStepping() = default;
103 
114  virtual double
115  evolve_one_time_step(
116  std::vector<std::function<VectorType(const double, const VectorType &)>>
117  & F,
118  std::vector<std::function<
119  VectorType(const double, const double, const VectorType &)>> &J_inverse,
120  double t,
121  double delta_t,
122  VectorType & y) = 0;
123 
127  struct Status
128  {};
129 
133  virtual const Status &
134  get_status() const = 0;
135  };
136 
137 
138 
145  template <typename VectorType>
146  class RungeKutta : public TimeStepping<VectorType>
147  {
148  public:
152  virtual ~RungeKutta() override = default;
153 
157  virtual void
158  initialize(const runge_kutta_method method) = 0;
159 
171  double
173  std::vector<std::function<VectorType(const double, const VectorType &)>>
174  & F,
175  std::vector<std::function<
176  VectorType(const double, const double, const VectorType &)>> &J_inverse,
177  double t,
178  double delta_t,
179  VectorType &y) override;
180 
192  virtual double
194  const std::function<VectorType(const double, const VectorType &)> &f,
195  const std::function<
196  VectorType(const double, const double, const VectorType &)>
197  & id_minus_tau_J_inverse,
198  double t,
199  double delta_t,
200  VectorType &y) = 0;
201 
202  protected:
206  unsigned int n_stages;
207 
211  std::vector<double> b;
212 
216  std::vector<double> c;
217 
221  std::vector<std::vector<double>> a;
222  };
223 
224 
225 
230  template <typename VectorType>
231  class ExplicitRungeKutta : public RungeKutta<VectorType>
232  {
233  public:
235 
241  ExplicitRungeKutta() = default;
242 
247 
251  void
252  initialize(const runge_kutta_method method) override;
253 
265  double
267  const std::function<VectorType(const double, const VectorType &)> &f,
268  const std::function<
269  VectorType(const double, const double, const VectorType &)>
270  & id_minus_tau_J_inverse,
271  double t,
272  double delta_t,
273  VectorType &y) override;
274 
282  double
284  const std::function<VectorType(const double, const VectorType &)> &f,
285  double t,
286  double delta_t,
287  VectorType &y);
288 
292  struct Status : public TimeStepping<VectorType>::Status
293  {
295  : method(invalid)
296  {}
297 
299  };
300 
304  const Status &
305  get_status() const override;
306 
307  private:
311  void
313  const std::function<VectorType(const double, const VectorType &)> &f,
314  const double t,
315  const double delta_t,
316  const VectorType & y,
317  std::vector<VectorType> &f_stages) const;
318 
323  };
324 
325 
326 
331  template <typename VectorType>
332  class ImplicitRungeKutta : public RungeKutta<VectorType>
333  {
334  public:
336 
342  ImplicitRungeKutta() = default;
343 
350  const unsigned int max_it = 100,
351  const double tolerance = 1e-6);
352 
356  void
357  initialize(const runge_kutta_method method) override;
358 
370  double
372  const std::function<VectorType(const double, const VectorType &)> &f,
373  const std::function<
374  VectorType(const double, const double, const VectorType &)>
375  & id_minus_tau_J_inverse,
376  double t,
377  double delta_t,
378  VectorType &y) override;
379 
384  void
385  set_newton_solver_parameters(const unsigned int max_it,
386  const double tolerance);
387 
392  struct Status : public TimeStepping<VectorType>::Status
393  {
395  : method(invalid)
398  {}
399 
401  unsigned int n_iterations;
403  };
404 
408  const Status &
409  get_status() const override;
410 
411  private:
415  void
417  const std::function<VectorType(const double, const VectorType &)> &f,
418  const std::function<
419  VectorType(const double, const double, const VectorType &)>
420  & id_minus_tau_J_inverse,
421  double t,
422  double delta_t,
423  VectorType & y,
424  std::vector<VectorType> &f_stages);
425 
429  void
430  newton_solve(
431  const std::function<void(const VectorType &, VectorType &)> &get_residual,
432  const std::function<VectorType(const VectorType &)>
433  & id_minus_tau_J_inverse,
434  VectorType &y);
435 
439  void
441  const std::function<VectorType(const double, const VectorType &)> &f,
442  double t,
443  double delta_t,
444  const VectorType &new_y,
445  const VectorType &y,
446  VectorType & tendency,
447  VectorType & residual) const;
448 
455 
459  unsigned int max_it;
460 
464  double tolerance;
465 
470  };
471 
472 
473 
478  template <typename VectorType>
479  class EmbeddedExplicitRungeKutta : public RungeKutta<VectorType>
480  {
481  public:
483 
489  EmbeddedExplicitRungeKutta() = default;
490 
496  const double coarsen_param = 1.2,
497  const double refine_param = 0.8,
498  const double min_delta = 1e-14,
499  const double max_delta = 1e100,
500  const double refine_tol = 1e-8,
501  const double coarsen_tol = 1e-12);
502 
507  {
508  free_memory();
509  }
510 
514  void
515  free_memory();
516 
520  void
521  initialize(const runge_kutta_method method) override;
522 
534  double
536  const std::function<VectorType(const double, const VectorType &)> &f,
537  const std::function<
538  VectorType(const double, const double, const VectorType &)>
539  & id_minus_tau_J_inverse,
540  double t,
541  double delta_t,
542  VectorType &y) override;
543 
551  double
553  const std::function<VectorType(const double, const VectorType &)> &f,
554  double t,
555  double delta_t,
556  VectorType &y);
557 
561  void
563  const double refine_param,
564  const double min_delta,
565  const double max_delta,
566  const double refine_tol,
567  const double coarsen_tol);
568 
575  struct Status : public TimeStepping<VectorType>::Status
576  {
579  unsigned int n_iterations;
581  double error_norm;
582  };
583 
587  const Status &
588  get_status() const override;
589 
590  private:
594  void
596  const std::function<VectorType(const double, const VectorType &)> &f,
597  const double t,
598  const double delta_t,
599  const VectorType & y,
600  std::vector<VectorType> &f_stages);
601 
607 
612  double refine_param;
613 
617  double min_delta_t;
618 
622  double max_delta_t;
623 
628  double refine_tol;
629 
634  double coarsen_tol;
635 
641 
645  std::vector<double> b1;
646 
650  std::vector<double> b2;
651 
657 
662  };
663 } // namespace TimeStepping
664 
666 
667 #endif
TimeStepping::ImplicitRungeKutta::evolve_one_time_step
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
TimeStepping::RungeKutta::evolve_one_time_step
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
Physics::Elasticity::Kinematics::F
Tensor< 2, dim, Number > F(const Tensor< 2, dim, Number > &Grad_u)
TimeStepping::RungeKutta
Definition: time_stepping.h:146
TimeStepping::RK_THIRD_ORDER
@ RK_THIRD_ORDER
Definition: time_stepping.h:60
TimeStepping::EmbeddedExplicitRungeKutta::b2
std::vector< double > b2
Definition: time_stepping.h:650
TimeStepping::CASH_KARP
@ CASH_KARP
Definition: time_stepping.h:70
TimeStepping::EmbeddedExplicitRungeKutta::~EmbeddedExplicitRungeKutta
~EmbeddedExplicitRungeKutta() override
Definition: time_stepping.h:506
TimeStepping::EmbeddedExplicitRungeKutta::get_status
const Status & get_status() const override
TimeStepping::FORWARD_EULER
@ FORWARD_EULER
Definition: time_stepping.h:59
TimeStepping::runge_kutta_method
runge_kutta_method
Definition: time_stepping.h:57
TimeStepping::EmbeddedExplicitRungeKutta::refine_param
double refine_param
Definition: time_stepping.h:612
TimeStepping::EmbeddedExplicitRungeKutta::Status::method
runge_kutta_method method
Definition: time_stepping.h:577
TimeStepping::EmbeddedExplicitRungeKutta::min_delta_t
double min_delta_t
Definition: time_stepping.h:617
TimeStepping::TimeStepping::Status
Definition: time_stepping.h:127
TimeStepping::BACKWARD_EULER
@ BACKWARD_EULER
Definition: time_stepping.h:62
TimeStepping::ExplicitRungeKutta::status
Status status
Definition: time_stepping.h:322
TimeStepping::ExplicitRungeKutta
Definition: time_stepping.h:231
TimeStepping::EmbeddedExplicitRungeKutta::Status::delta_t_guess
double delta_t_guess
Definition: time_stepping.h:580
TimeStepping
Definition: time_stepping.h:36
TimeStepping::DELTA_T
@ DELTA_T
Definition: time_stepping.h:84
VectorType
TimeStepping::RungeKutta::~RungeKutta
virtual ~RungeKutta() override=default
TimeStepping::ImplicitRungeKutta::compute_stages
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)
TimeStepping::ExplicitRungeKutta::ExplicitRungeKutta
ExplicitRungeKutta()=default
TimeStepping::EmbeddedExplicitRungeKutta::free_memory
void free_memory()
Physics::Elasticity::Kinematics::e
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
TimeStepping::ExplicitRungeKutta::evolve_one_time_step
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
TimeStepping::EmbeddedExplicitRungeKutta::max_delta_t
double max_delta_t
Definition: time_stepping.h:622
TimeStepping::HEUN_EULER
@ HEUN_EULER
Definition: time_stepping.h:66
TimeStepping::EmbeddedExplicitRungeKutta::initialize
void initialize(const runge_kutta_method method) override
TimeStepping::EmbeddedExplicitRungeKutta::Status::error_norm
double error_norm
Definition: time_stepping.h:581
TimeStepping::ImplicitRungeKutta::initialize
void initialize(const runge_kutta_method method) override
TimeStepping::RungeKutta::c
std::vector< double > c
Definition: time_stepping.h:216
TimeStepping::SDIRK_TWO_STAGES
@ SDIRK_TWO_STAGES
Definition: time_stepping.h:65
TimeStepping::RungeKutta::a
std::vector< std::vector< double > > a
Definition: time_stepping.h:221
TimeStepping::ImplicitRungeKutta::Status::method
runge_kutta_method method
Definition: time_stepping.h:400
TimeStepping::EmbeddedExplicitRungeKutta::coarsen_param
double coarsen_param
Definition: time_stepping.h:606
TimeStepping::RungeKutta::n_stages
unsigned int n_stages
Definition: time_stepping.h:206
TimeStepping::ImplicitRungeKutta::Status::Status
Status()
Definition: time_stepping.h:394
TimeStepping::invalid
@ invalid
Definition: time_stepping.h:71
TimeStepping::ExplicitRungeKutta::Status
Definition: time_stepping.h:292
TimeStepping::ExplicitRungeKutta::Status::Status
Status()
Definition: time_stepping.h:294
double
TimeStepping::ImplicitRungeKutta::Status
Definition: time_stepping.h:392
TimeStepping::ExplicitRungeKutta::get_status
const Status & get_status() const override
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
numbers::signaling_nan
T signaling_nan()
Definition: signaling_nan.h:229
TimeStepping::RungeKutta::b
std::vector< double > b
Definition: time_stepping.h:211
TimeStepping::ExplicitRungeKutta::compute_stages
void compute_stages(const std::function< VectorType(const double, const VectorType &)> &f, const double t, const double delta_t, const VectorType &y, std::vector< VectorType > &f_stages) const
TimeStepping::ImplicitRungeKutta::compute_residual
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
TimeStepping::ExplicitRungeKutta::Status::method
runge_kutta_method method
Definition: time_stepping.h:298
TimeStepping::EmbeddedExplicitRungeKutta::status
Status status
Definition: time_stepping.h:661
TimeStepping::MAX_DELTA_T
@ MAX_DELTA_T
Definition: time_stepping.h:86
TimeStepping::DOPRI
@ DOPRI
Definition: time_stepping.h:68
TimeStepping::ImplicitRungeKutta::set_newton_solver_parameters
void set_newton_solver_parameters(const unsigned int max_it, const double tolerance)
TimeStepping::EmbeddedExplicitRungeKutta::b1
std::vector< double > b1
Definition: time_stepping.h:645
TimeStepping::ImplicitRungeKutta::status
Status status
Definition: time_stepping.h:469
numbers
Definition: numbers.h:207
TimeStepping::EmbeddedExplicitRungeKutta
Definition: time_stepping.h:479
TimeStepping::EmbeddedExplicitRungeKutta::Status::exit_delta_t
embedded_runge_kutta_time_step exit_delta_t
Definition: time_stepping.h:578
TimeStepping::EmbeddedExplicitRungeKutta::last_stage
VectorType * last_stage
Definition: time_stepping.h:656
TimeStepping::ExplicitRungeKutta::initialize
void initialize(const runge_kutta_method method) override
TimeStepping::ImplicitRungeKutta::Status::n_iterations
unsigned int n_iterations
Definition: time_stepping.h:401
TimeStepping::ImplicitRungeKutta::Status::norm_residual
double norm_residual
Definition: time_stepping.h:402
TimeStepping::RK_CLASSIC_FOURTH_ORDER
@ RK_CLASSIC_FOURTH_ORDER
Definition: time_stepping.h:61
TimeStepping::RungeKutta::initialize
virtual void initialize(const runge_kutta_method method)=0
TimeStepping::ImplicitRungeKutta::skip_linear_combi
bool skip_linear_combi
Definition: time_stepping.h:454
TimeStepping::EmbeddedExplicitRungeKutta::Status::n_iterations
unsigned int n_iterations
Definition: time_stepping.h:579
TimeStepping::ImplicitRungeKutta::ImplicitRungeKutta
ImplicitRungeKutta()=default
TimeStepping::EmbeddedExplicitRungeKutta::EmbeddedExplicitRungeKutta
EmbeddedExplicitRungeKutta()=default
TimeStepping::BOGACKI_SHAMPINE
@ BOGACKI_SHAMPINE
Definition: time_stepping.h:67
TimeStepping::FEHLBERG
@ FEHLBERG
Definition: time_stepping.h:69
TimeStepping::CRANK_NICOLSON
@ CRANK_NICOLSON
Definition: time_stepping.h:64
TimeStepping::EmbeddedExplicitRungeKutta::evolve_one_time_step
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
numbers::invalid_unsigned_int
static const unsigned int invalid_unsigned_int
Definition: types.h:191
TimeStepping::EmbeddedExplicitRungeKutta::refine_tol
double refine_tol
Definition: time_stepping.h:628
signaling_nan.h
TimeStepping::ImplicitRungeKutta::max_it
unsigned int max_it
Definition: time_stepping.h:459
config.h
TimeStepping::embedded_runge_kutta_time_step
embedded_runge_kutta_time_step
Definition: time_stepping.h:82
TimeStepping::EmbeddedExplicitRungeKutta::last_same_as_first
bool last_same_as_first
Definition: time_stepping.h:640
TimeStepping::ImplicitRungeKutta
Definition: time_stepping.h:332
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
TimeStepping::MIN_DELTA_T
@ MIN_DELTA_T
Definition: time_stepping.h:85
TimeStepping::IMPLICIT_MIDPOINT
@ IMPLICIT_MIDPOINT
Definition: time_stepping.h:63
TimeStepping::EmbeddedExplicitRungeKutta::Status
Definition: time_stepping.h:575
TimeStepping::EmbeddedExplicitRungeKutta::set_time_adaptation_parameters
void set_time_adaptation_parameters(const double coarsen_param, const double refine_param, const double min_delta, const double max_delta, const double refine_tol, const double coarsen_tol)
TimeStepping::ImplicitRungeKutta::newton_solve
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)
TimeStepping::EmbeddedExplicitRungeKutta::coarsen_tol
double coarsen_tol
Definition: time_stepping.h:634
TimeStepping::ImplicitRungeKutta::tolerance
double tolerance
Definition: time_stepping.h:464
TimeStepping::ImplicitRungeKutta::get_status
const Status & get_status() const override
TimeStepping::EmbeddedExplicitRungeKutta::compute_stages
void compute_stages(const std::function< VectorType(const double, const VectorType &)> &f, const double t, const double delta_t, const VectorType &y, std::vector< VectorType > &f_stages)