Reference documentation for deal.II version Git e8f99f7e31 2020-07-08 12:50:43 +0200
\(\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 
33 namespace TimeStepping
34 {
55  {
69  };
70 
71 
72 
80  {
84  };
85 
86 
87 
92  template <typename VectorType>
94  {
95  public:
99  virtual ~TimeStepping() = default;
100 
111  virtual double
112  evolve_one_time_step(
113  std::vector<std::function<VectorType(const double, const VectorType &)>>
114  & F,
115  std::vector<std::function<
116  VectorType(const double, const double, const VectorType &)>> &J_inverse,
117  double t,
118  double delta_t,
119  VectorType & y) = 0;
120 
124  struct Status
125  {};
126 
130  virtual const Status &
131  get_status() const = 0;
132  };
133 
134 
135 
139  template <typename VectorType>
140  class RungeKutta : public TimeStepping<VectorType>
141  {
142  public:
146  virtual ~RungeKutta() override = default;
147 
151  virtual void
152  initialize(const runge_kutta_method method) = 0;
153 
165  double
166  evolve_one_time_step(
167  std::vector<std::function<VectorType(const double, const VectorType &)>>
168  & F,
169  std::vector<std::function<
170  VectorType(const double, const double, const VectorType &)>> &J_inverse,
171  double t,
172  double delta_t,
173  VectorType &y) override;
174 
186  virtual double
187  evolve_one_time_step(
188  const std::function<VectorType(const double, const VectorType &)> &f,
189  const std::function<
190  VectorType(const double, const double, const VectorType &)>
191  & id_minus_tau_J_inverse,
192  double t,
193  double delta_t,
194  VectorType &y) = 0;
195 
196  protected:
200  unsigned int n_stages;
201 
205  std::vector<double> b;
206 
210  std::vector<double> c;
211 
215  std::vector<std::vector<double>> a;
216  };
217 
218 
219 
224  template <typename VectorType>
225  class ExplicitRungeKutta : public RungeKutta<VectorType>
226  {
227  public:
229 
235  ExplicitRungeKutta() = default;
236 
241 
245  void
246  initialize(const runge_kutta_method method) override;
247 
259  double
260  evolve_one_time_step(
261  const std::function<VectorType(const double, const VectorType &)> &f,
262  const std::function<
263  VectorType(const double, const double, const VectorType &)>
264  & id_minus_tau_J_inverse,
265  double t,
266  double delta_t,
267  VectorType &y) override;
268 
276  double
277  evolve_one_time_step(
278  const std::function<VectorType(const double, const VectorType &)> &f,
279  double t,
280  double delta_t,
281  VectorType &y);
282 
286  struct Status : public TimeStepping<VectorType>::Status
287  {
289  : method(invalid)
290  {}
291 
293  };
294 
298  const Status &
299  get_status() const override;
300 
301  private:
305  void
306  compute_stages(
307  const std::function<VectorType(const double, const VectorType &)> &f,
308  const double t,
309  const double delta_t,
310  const VectorType & y,
311  std::vector<VectorType> &f_stages) const;
312 
317  };
318 
319 
320 
325  template <typename VectorType>
326  class ImplicitRungeKutta : public RungeKutta<VectorType>
327  {
328  public:
330 
336  ImplicitRungeKutta() = default;
337 
344  const unsigned int max_it = 100,
345  const double tolerance = 1e-6);
346 
350  void
351  initialize(const runge_kutta_method method) override;
352 
364  double
365  evolve_one_time_step(
366  const std::function<VectorType(const double, const VectorType &)> &f,
367  const std::function<
368  VectorType(const double, const double, const VectorType &)>
369  & id_minus_tau_J_inverse,
370  double t,
371  double delta_t,
372  VectorType &y) override;
373 
378  void
379  set_newton_solver_parameters(const unsigned int max_it,
380  const double tolerance);
381 
386  struct Status : public TimeStepping<VectorType>::Status
387  {
389  : method(invalid)
390  , n_iterations(numbers::invalid_unsigned_int)
391  , norm_residual(numbers::signaling_nan<double>())
392  {}
393 
395  unsigned int n_iterations;
397  };
398 
402  const Status &
403  get_status() const override;
404 
405  private:
409  void
410  compute_stages(
411  const std::function<VectorType(const double, const VectorType &)> &f,
412  const std::function<
413  VectorType(const double, const double, const VectorType &)>
414  & id_minus_tau_J_inverse,
415  double t,
416  double delta_t,
417  VectorType & y,
418  std::vector<VectorType> &f_stages);
419 
423  void
424  newton_solve(
425  const std::function<void(const VectorType &, VectorType &)> &get_residual,
426  const std::function<VectorType(const VectorType &)>
427  & id_minus_tau_J_inverse,
428  VectorType &y);
429 
433  void
434  compute_residual(
435  const std::function<VectorType(const double, const VectorType &)> &f,
436  double t,
437  double delta_t,
438  const VectorType &new_y,
439  const VectorType &y,
440  VectorType & tendency,
441  VectorType & residual) const;
442 
449 
453  unsigned int max_it;
454 
458  double tolerance;
459 
464  };
465 
466 
467 
472  template <typename VectorType>
473  class EmbeddedExplicitRungeKutta : public RungeKutta<VectorType>
474  {
475  public:
477 
483  EmbeddedExplicitRungeKutta() = default;
484 
490  const double coarsen_param = 1.2,
491  const double refine_param = 0.8,
492  const double min_delta = 1e-14,
493  const double max_delta = 1e100,
494  const double refine_tol = 1e-8,
495  const double coarsen_tol = 1e-12);
496 
501  {
502  free_memory();
503  }
504 
508  void
509  free_memory();
510 
514  void
515  initialize(const runge_kutta_method method) override;
516 
528  double
529  evolve_one_time_step(
530  const std::function<VectorType(const double, const VectorType &)> &f,
531  const std::function<
532  VectorType(const double, const double, const VectorType &)>
533  & id_minus_tau_J_inverse,
534  double t,
535  double delta_t,
536  VectorType &y) override;
537 
545  double
546  evolve_one_time_step(
547  const std::function<VectorType(const double, const VectorType &)> &f,
548  double t,
549  double delta_t,
550  VectorType &y);
551 
555  void
556  set_time_adaptation_parameters(const double coarsen_param,
557  const double refine_param,
558  const double min_delta,
559  const double max_delta,
560  const double refine_tol,
561  const double coarsen_tol);
562 
569  struct Status : public TimeStepping<VectorType>::Status
570  {
573  unsigned int n_iterations;
575  double error_norm;
576  };
577 
581  const Status &
582  get_status() const override;
583 
584  private:
588  void
589  compute_stages(
590  const std::function<VectorType(const double, const VectorType &)> &f,
591  const double t,
592  const double delta_t,
593  const VectorType & y,
594  std::vector<VectorType> &f_stages);
595 
601 
606  double refine_param;
607 
611  double min_delta_t;
612 
616  double max_delta_t;
617 
622  double refine_tol;
623 
628  double coarsen_tol;
629 
635 
639  std::vector<double> b1;
640 
644  std::vector<double> b2;
645 
651 
656  };
657 } // namespace TimeStepping
658 
660 
661 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:191
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
Tensor< 2, dim, Number > F(const Tensor< 2, dim, Number > &Grad_u)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
std::vector< double > b
std::vector< double > c
std::vector< std::vector< double > > a
embedded_runge_kutta_time_step exit_delta_t
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
T signaling_nan()
embedded_runge_kutta_time_step
Definition: time_stepping.h:79