dune-fem  2.4.1-rc
implicit.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_IMPLICIT_HH
2 #define DUNE_FEM_SOLVER_RUNGEKUTTA_IMPLICIT_HH
3 
4 //- system includes
5 #include <sstream>
6 #include <vector>
7 
8 //- dune-common includes
9 #include <dune/common/exceptions.hh>
10 
11 //- dune-fem includes
15 
16 namespace DuneODE
17 {
18 
19  // ImplicitRungeKuttaSolver
20  // ------------------------
21 
23  template< class HelmholtzOperator, class NonlinearSolver, class TimeStepControl = ImplicitRungeKuttaTimeStepControl >
25  : public BasicImplicitRungeKuttaSolver< HelmholtzOperator, NonlinearSolver, TimeStepControl >
26  {
29 
30  public:
31  typedef HelmholtzOperator HelmholtzOperatorType;
33 
34  typedef typename TimeStepControlType::TimeProviderType TimeProviderType;
37 
46  ImplicitRungeKuttaSolver ( HelmholtzOperatorType &helmholtzOp,
47  TimeProviderType &timeProvider,
48  int order,
49  const ParametersType& tscParam,
50  const NonlinearSolverParametersType& nlsParam )
51  : BaseType( helmholtzOp,
52  butcherTable( tscParam.selectedSolver( order ) ),
53  TimeStepControlType( timeProvider, tscParam ),
54  nlsParam )
55  {}
56 
57  ImplicitRungeKuttaSolver ( HelmholtzOperatorType &helmholtzOp,
58  TimeProviderType &timeProvider,
59  int order,
61  : BaseType( helmholtzOp,
62  butcherTable( ParametersType( parameter ).selectedSolver( order ) ),
63  TimeStepControlType( timeProvider, ParametersType( parameter ) ),
64  NonlinearSolverParametersType( parameter ) )
65  {}
66 
74  ImplicitRungeKuttaSolver ( HelmholtzOperatorType &helmholtzOp,
75  TimeProviderType &timeProvider,
76  const ParametersType& tscParam = ParametersType(),
77  const NonlinearSolverParametersType& nlsParam = NonlinearSolverParametersType() )
78  : BaseType( helmholtzOp,
79  butcherTable( tscParam.selectedSolver( 1 ) ),
80  TimeStepControlType( timeProvider, tscParam ),
81  nlsParam )
82  {}
83 
84  ImplicitRungeKuttaSolver ( HelmholtzOperatorType &helmholtzOp,
85  TimeProviderType &timeProvider,
87  : BaseType( helmholtzOp,
88  butcherTable( ParametersType( parameter ).selectedSolver( 1 ) ),
89  TimeStepControlType( timeProvider, ParametersType( parameter ) ),
90  NonlinearSolverParametersType( parameter ) )
91  {}
92 
93  protected:
94  static SimpleButcherTable< double > butcherTable ( const int solverId )
95  {
96  switch( solverId )
97  {
98  case 1:
100  case 2:
101  return gauss2ButcherTable();
102  case 3:
103  return implicit3ButcherTable();
104  case 4:
105  return implicit34ButcherTable();
106  default:
107  DUNE_THROW( NotImplemented, "Implicit Runge-Kutta method with id " << solverId << " not implemented." );
108  }
109  }
110  };
111 
112 } // namespace DuneODE
113 
114 #endif // #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_IMPLICIT_HH
static SimpleButcherTable< double > butcherTable(const int solverId)
Definition: implicit.hh:94
ImplicitRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const ParametersType &tscParam=ParametersType(), const NonlinearSolverParametersType &nlsParam=NonlinearSolverParametersType())
constructor
Definition: implicit.hh:74
BaseType::ParametersType ParametersType
Definition: implicit.hh:35
BaseType::TimeStepControlType TimeStepControlType
Definition: implicit.hh:32
Implicit RungeKutta ODE solver.
Definition: implicit.hh:24
Definition: butchertable.hh:15
NonlinearSolver::ParametersType NonlinearSolverParametersType
Definition: basicimplicit.hh:70
ImplicitRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
Definition: implicit.hh:84
SimpleButcherTable< double > implicitEulerButcherTable()
Definition: butchertable.cc:58
Implicit RungeKutta ODE solver.
Definition: basicimplicit.hh:52
Definition: multistep.hh:16
ImplicitRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, int order, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
Definition: implicit.hh:57
SimpleButcherTable< double > gauss2ButcherTable()
Definition: butchertable.cc:73
SimpleButcherTable< double > implicit3ButcherTable()
Definition: butchertable.cc:44
HelmholtzOperator HelmholtzOperatorType
Definition: implicit.hh:31
static ParameterContainer & container()
Definition: io/parameter.hh:190
TimeStepControlType::TimeProviderType TimeProviderType
Definition: implicit.hh:34
SimpleButcherTable< double > implicit34ButcherTable()
Definition: butchertable.cc:24
BaseType::NonlinearSolverParametersType NonlinearSolverParametersType
Definition: implicit.hh:36
TimeStepControlType::ParametersType ParametersType
Definition: basicimplicit.hh:69
ImplicitRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, int order, const ParametersType &tscParam, const NonlinearSolverParametersType &nlsParam)
constructor
Definition: implicit.hh:46