deal.II version GIT relicensing-2250-g88cb8ba3cb 2024-12-13 12:20:00+00:00
\(\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\}}\)
Loading...
Searching...
No Matches
time_stepping.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2014 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_time_stepping_h
16#define dealii_time_stepping_h
17
18
19#include <deal.II/base/config.h>
20
22
23#include <functional>
24#include <vector>
25
27
32namespace TimeStepping
33{
144
145
146
166
167
168
173 template <typename VectorType>
175 {
176 public:
180 virtual ~TimeStepping() = default;
181
192 virtual double
194 std::vector<std::function<VectorType(const double, const VectorType &)>>
195 &F,
196 std::vector<std::function<
197 VectorType(const double, const double, const VectorType &)>> &J_inverse,
198 double t,
199 double delta_t,
200 VectorType &y) = 0;
201
205 struct Status
206 {};
207
211 virtual const Status &
212 get_status() const = 0;
213 };
214
215
216
220 template <typename VectorType>
221 class RungeKutta : public TimeStepping<VectorType>
222 {
223 public:
227 virtual ~RungeKutta() override = default;
228
232 virtual void
234
246 double
248 std::vector<std::function<VectorType(const double, const VectorType &)>>
249 &F,
250 std::vector<std::function<
251 VectorType(const double, const double, const VectorType &)>> &J_inverse,
252 double t,
253 double delta_t,
254 VectorType &y) override;
255
267 virtual double
269 const std::function<VectorType(const double, const VectorType &)> &f,
270 const std::function<
271 VectorType(const double, const double, const VectorType &)>
273 double t,
274 double delta_t,
275 VectorType &y) = 0;
276
277 protected:
281 unsigned int n_stages;
282
286 std::vector<double> b;
287
291 std::vector<double> c;
292
296 std::vector<std::vector<double>> a;
297 };
298
299
300
305 template <typename VectorType>
306 class ExplicitRungeKutta : public RungeKutta<VectorType>
307 {
308 public:
309 using RungeKutta<VectorType>::evolve_one_time_step;
310
317
322
326 void
327 initialize(const runge_kutta_method method) override;
328
342 double
344 const std::function<VectorType(const double, const VectorType &)> &f,
345 const std::function<
346 VectorType(const double, const double, const VectorType &)>
348 double t,
349 double delta_t,
350 VectorType &y) override;
351
359 double
361 const std::function<VectorType(const double, const VectorType &)> &f,
362 double t,
363 double delta_t,
364 VectorType &y);
365
369 struct Status : public TimeStepping<VectorType>::Status
370 {
372 : method(invalid)
373 {}
374
376 };
377
381 const Status &
382 get_status() const override;
383
384 private:
388 void
390 const std::function<VectorType(const double, const VectorType &)> &f,
391 const double t,
392 const double delta_t,
393 const VectorType &y,
394 std::vector<VectorType> &f_stages) const;
395
400 };
401
402
403
409 template <typename VectorType>
410 class LowStorageRungeKutta : public RungeKutta<VectorType>
411 {
412 public:
413 using RungeKutta<VectorType>::evolve_one_time_step;
414
421
426
430 void
431 initialize(const runge_kutta_method method) override;
432
444 double
446 const std::function<VectorType(const double, const VectorType &)> &f,
447 const std::function<
448 VectorType(const double, const double, const VectorType &)>
450 double t,
451 double delta_t,
452 VectorType &y) override;
453
463 double
465 const std::function<VectorType(const double, const VectorType &)> &f,
466 double t,
467 double delta_t,
468 VectorType &solution,
469 VectorType &vec_ri,
470 VectorType &vec_ki);
471
478 void
479 get_coefficients(std::vector<double> &a,
480 std::vector<double> &b,
481 std::vector<double> &c) const;
482
486 struct Status : public TimeStepping<VectorType>::Status
487 {
489 : method(invalid)
490 {}
491
493 };
494
498 const Status &
499 get_status() const override;
500
501 private:
505 void
507 const std::function<VectorType(const double, const VectorType &)> &f,
508 const double t,
509 const double factor_solution,
510 const double factor_ai,
511 const VectorType &current_ri,
512 VectorType &vec_ki,
513 VectorType &solution,
514 VectorType &next_ri) const;
515
520 };
521
522
523
528 template <typename VectorType>
529 class ImplicitRungeKutta : public RungeKutta<VectorType>
530 {
531 public:
532 using RungeKutta<VectorType>::evolve_one_time_step;
533
540
547 const unsigned int max_it = 100,
548 const double tolerance = 1e-6);
549
553 void
554 initialize(const runge_kutta_method method) override;
555
567 double
569 const std::function<VectorType(const double, const VectorType &)> &f,
570 const std::function<
571 VectorType(const double, const double, const VectorType &)>
573 double t,
574 double delta_t,
575 VectorType &y) override;
576
581 void
583 const double tolerance);
584
589 struct Status : public TimeStepping<VectorType>::Status
590 {
592 : method(invalid)
593 , n_iterations(numbers::invalid_unsigned_int)
594 , norm_residual(numbers::signaling_nan<double>())
595 {}
596
598 unsigned int n_iterations;
600 };
601
605 const Status &
606 get_status() const override;
607
608 private:
612 void
614 const std::function<VectorType(const double, const VectorType &)> &f,
615 const std::function<
616 VectorType(const double, const double, const VectorType &)>
618 double t,
619 double delta_t,
620 VectorType &y,
621 std::vector<VectorType> &f_stages);
622
626 void
628 const std::function<void(const VectorType &, VectorType &)> &get_residual,
629 const std::function<VectorType(const VectorType &)>
631 VectorType &y);
632
636 void
638 const std::function<VectorType(const double, const VectorType &)> &f,
639 double t,
640 double delta_t,
641 const VectorType &new_y,
642 const VectorType &y,
643 VectorType &tendency,
644 VectorType &residual) const;
645
649 unsigned int max_it;
650
654 double tolerance;
655
660 };
661
662
663
668 template <typename VectorType>
669 class EmbeddedExplicitRungeKutta : public RungeKutta<VectorType>
670 {
671 public:
672 using RungeKutta<VectorType>::evolve_one_time_step;
673
680
686 const double coarsen_param = 1.2,
687 const double refine_param = 0.8,
688 const double min_delta = 1e-14,
689 const double max_delta = 1e100,
690 const double refine_tol = 1e-8,
691 const double coarsen_tol = 1e-12);
692
697 {
698 free_memory();
699 }
700
704 void
706
710 void
711 initialize(const runge_kutta_method method) override;
712
726 double
728 const std::function<VectorType(const double, const VectorType &)> &f,
729 const std::function<
730 VectorType(const double, const double, const VectorType &)>
732 double t,
733 double delta_t,
734 VectorType &y) override;
735
743 double
745 const std::function<VectorType(const double, const VectorType &)> &f,
746 double t,
747 double delta_t,
748 VectorType &y);
749
753 void
755 const double refine_param,
756 const double min_delta,
757 const double max_delta,
758 const double refine_tol,
759 const double coarsen_tol);
760
767 struct Status : public TimeStepping<VectorType>::Status
768 {
770 : method(invalid)
771 {}
772
775 unsigned int n_iterations;
778 };
779
783 const Status &
784 get_status() const override;
785
786 private:
790 void
792 const std::function<VectorType(const double, const VectorType &)> &f,
793 const double t,
794 const double delta_t,
795 const VectorType &y,
796 std::vector<VectorType> &f_stages);
797
803
809
814
819
825
831
836 bool last_same_as_first = false;
837
841 std::vector<double> b1;
842
846 std::vector<double> b2;
847
852 VectorType *last_stage = nullptr;
853
858 };
859} // namespace TimeStepping
860
862
863#endif
double evolve_one_time_step(const std::function< VectorType(const double, const VectorType &)> &f, double t, double delta_t, VectorType &y)
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 Status & get_status() const override
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)
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
EmbeddedExplicitRungeKutta(const runge_kutta_method method, const double coarsen_param=1.2, const double refine_param=0.8, const double min_delta=1e-14, const double max_delta=1e100, const double refine_tol=1e-8, const double coarsen_tol=1e-12)
const Status & get_status() const override
void initialize(const runge_kutta_method method) override
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
double evolve_one_time_step(const std::function< VectorType(const double, const VectorType &)> &f, double t, double delta_t, VectorType &y)
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
ExplicitRungeKutta(const runge_kutta_method method)
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)
const Status & get_status() const override
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 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
void set_newton_solver_parameters(const unsigned int max_it, const double tolerance)
ImplicitRungeKutta(const runge_kutta_method method, const unsigned int max_it=100, const double tolerance=1e-6)
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 &current_ri, VectorType &vec_ki, VectorType &solution, VectorType &next_ri) const
void initialize(const runge_kutta_method method) override
LowStorageRungeKutta(const runge_kutta_method method)
void get_coefficients(std::vector< double > &a, std::vector< double > &b, std::vector< double > &c) const
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)
const Status & get_status() const 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
std::vector< double > b
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
virtual ~RungeKutta() override=default
virtual void initialize(const runge_kutta_method method)=0
std::vector< std::vector< double > > a
std::vector< double > c
virtual 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)=0
virtual 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)=0
virtual ~TimeStepping()=default
virtual const Status & get_status() const =0
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
@ LOW_STORAGE_RK_STAGE9_ORDER5
@ LOW_STORAGE_RK_STAGE3_ORDER3
@ LOW_STORAGE_RK_STAGE7_ORDER4
@ LOW_STORAGE_RK_STAGE5_ORDER4
embedded_runge_kutta_time_step exit_delta_t