dune-fem  2.8-git
semiimplicit.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_SEMIIMPLICIT_HH
2 #define DUNE_FEM_SOLVER_RUNGEKUTTA_SEMIIMPLICIT_HH
3 
4 //- system includes
5 #include <vector>
6 
7 //- dune-common includes
8 #include <dune/common/dynmatrix.hh>
9 #include <dune/common/dynvector.hh>
10 #include <dune/common/exceptions.hh>
11 
12 //- dune-fem includes
16 
17 namespace DuneODE
18 {
19 
20  // SemiImplicitRungeKuttaSourceTerm
21  // --------------------------------
22 
23  template< class ExplicitOperator >
25  {
27 
28  public:
29  typedef ExplicitOperator ExplicitOperatorType;
30 
31  typedef typename ExplicitOperatorType::DestinationType DestinationType;
32 
33  template< class ButcherTable >
35  const ButcherTable &butcherTable,
36  const Dune::DynamicMatrix< double > &implicitA )
37  : explicitOp_( explicitOp ),
38  alpha_( butcherTable.A() ),
39  gamma_( butcherTable.stages() ),
40  c_( butcherTable.c() ),
41  uex_( "SIRK u-explicit", explicitOp_.space() ),
42  limiter_( explicitOp_.hasLimiter() )
43  {
44  Dune::DynamicMatrix< double > Ainv( implicitA );
45  Ainv.invert();
46  alpha_.rightmultiply( Ainv );
47 
48  for( int i = 0; i < butcherTable.stages(); ++i )
49  {
50  gamma_[ i ] = 1.0;
51  for( int j = 0; j < i; ++j )
52  gamma_[ i ] -= alpha_[ i ][ j ];
53  }
54  }
55 
56  bool operator() ( double time, double timeStepSize, int stage,
57  const DestinationType &u, const std::vector< DestinationType * > &update,
58  DestinationType &source )
59  {
60  uex_.assign( u );
61  uex_ *= gamma_[ stage ];
62  for( int k = 0; k < stage; ++k )
63  uex_.axpy( alpha_[ stage ][ k ], *update[ k ] );
64  explicitOp_.setTime( time + c_[ stage ]*timeStepSize );
65  explicitOp_( uex_, source );
66 
67  return true;
68  }
69 
70  void limit( DestinationType& update, const double time )
71  {
72  if( limiter_ )
73  {
74  // set correct time
75  explicitOp_.setTime( time );
76  // copy given function
77  uex_.assign( update );
78  // apply limiter
79  explicitOp_.limit( uex_, update );
80  }
81  }
82 
83  double initialTimeStepEstimate ( double time, const DestinationType &u ) const
84  {
85  explicitOp_.setTime( time );
86  explicitOp_.initializeTimeStepSize( u );
87  return explicitOp_.timeStepEstimate();
88  }
89 
90  double timeStepEstimate () const
91  {
92  return explicitOp_.timeStepEstimate();
93  }
94 
95  private:
96  ExplicitOperatorType &explicitOp_;
97  Dune::DynamicMatrix< double > alpha_;
98  Dune::DynamicVector< double > gamma_, c_;
99  DestinationType uex_;
100  const bool limiter_ ;
101  };
102 
103 
104 
105  // SemiImplicitRungeKuttaSolver
106  // ----------------------------
107 
109  template< class ExplicitOperator, class HelmholtzOperator, class NonlinearSolver >
111  : public BasicImplicitRungeKuttaSolver< HelmholtzOperator, NonlinearSolver, ImplicitRungeKuttaTimeStepControl, SemiImplicitRungeKuttaSourceTerm< ExplicitOperator > >
112  {
115 
116  public:
117  typedef ExplicitOperator ExplicitOperatorType;
118  typedef HelmholtzOperator HelmholtzOperatorType;
121 
125 
136  HelmholtzOperatorType &helmholtzOp,
137  TimeProviderType &timeProvider, int order,
138  const ParameterType& tscParams,
139  const NonlinearSolverParameterType& nlsParams )
140  : BaseType( helmholtzOp,
141  defaultButcherTable( order, false ),
142  TimeStepControlType( timeProvider, tscParams ),
143  SourceTermType( explicitOp, defaultButcherTable( order, true ), defaultButcherTable( order, false ).A() ),
144  nlsParams )
145  {}
146 
148  HelmholtzOperatorType &helmholtzOp,
149  TimeProviderType &timeProvider,
150  int order,
152  : BaseType( helmholtzOp,
153  defaultButcherTable( order, false ),
154  TimeStepControlType( timeProvider, parameter ),
155  SourceTermType( explicitOp, defaultButcherTable( order, true ), defaultButcherTable( order, false ).A() ),
156  NonlinearSolverParameterType( parameter ) )
157  {}
158 
168  HelmholtzOperatorType &helmholtzOp,
169  TimeProviderType &timeProvider,
170  const ParameterType& tscParams,
171  const NonlinearSolverParameterType& nlsParams )
172  : BaseType( helmholtzOp,
173  defaultButcherTable( 1, false ),
174  TimeStepControlType( timeProvider, tscParams ),
175  SourceTermType( explicitOp, defaultButcherTable( 1, true ), defaultButcherTable( 1, false ).A() ),
176  nlsParams )
177  {}
178 
180  HelmholtzOperatorType &helmholtzOp,
181  TimeProviderType &timeProvider,
183  : BaseType( helmholtzOp,
184  defaultButcherTable( 1, false ),
185  TimeStepControlType( timeProvider, parameter ),
186  SourceTermType( explicitOp, defaultButcherTable( 1, true ), defaultButcherTable( 1, false ).A() ),
187  NonlinearSolverParameterType( parameter ) )
188  {}
189 
190 
192  HelmholtzOperatorType &helmholtzOp,
193  TimeProviderType &timeProvider,
194  const SimpleButcherTable< double >& explButcherTable,
195  const SimpleButcherTable< double >& implButcherTable,
197  : BaseType( helmholtzOp,
198  implButcherTable,
199  TimeStepControlType( timeProvider, parameter ),
200  SourceTermType( explicitOp, explButcherTable, implButcherTable.A() ),
201  NonlinearSolverParameterType( parameter ) )
202  {}
203  protected:
204  static SimpleButcherTable< double > defaultButcherTable ( int order, bool expl )
205  {
206  switch( order )
207  {
208  case 1:
209  return semiImplicitEulerButcherTable( expl );
210  case 2:
211  return semiImplicitSSP222ButcherTable( expl );
212  case 3:
213  return semiImplicit33ButcherTable( expl );
214  case 4:
215  return semiImplicitIERK45ButcherTable( expl );
216  default:
217  DUNE_THROW( NotImplemented, "Semi-implicit Runge-Kutta method of order " << order << " not implemented." );
218  }
219  }
220  };
221 
222 } // namespace DuneODE
223 
224 #endif // #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_IMPLICIT_HH
Definition: multistep.hh:17
SimpleButcherTable< double > semiImplicit33ButcherTable(bool expl)
Definition: butchertable.cc:256
SimpleButcherTable< double > semiImplicitEulerButcherTable(bool expl)
Definition: butchertable.cc:190
SimpleButcherTable< double > semiImplicitIERK45ButcherTable(bool expl)
Definition: butchertable.cc:335
SimpleButcherTable< double > semiImplicitSSP222ButcherTable(bool expl)
Definition: butchertable.cc:286
static ParameterContainer & container()
Definition: io/parameter.hh:193
Implicit RungeKutta ODE solver.
Definition: basicimplicit.hh:54
Definition: butchertable.hh:17
Definition: semiimplicit.hh:25
ExplicitOperatorType::DestinationType DestinationType
Definition: semiimplicit.hh:31
ExplicitOperator ExplicitOperatorType
Definition: semiimplicit.hh:29
double timeStepEstimate() const
Definition: semiimplicit.hh:90
bool operator()(double time, double timeStepSize, int stage, const DestinationType &u, const std::vector< DestinationType * > &update, DestinationType &source)
Definition: semiimplicit.hh:56
double initialTimeStepEstimate(double time, const DestinationType &u) const
Definition: semiimplicit.hh:83
SemiImplicitRungeKuttaSourceTerm(ExplicitOperatorType &explicitOp, const ButcherTable &butcherTable, const Dune::DynamicMatrix< double > &implicitA)
Definition: semiimplicit.hh:34
void limit(DestinationType &update, const double time)
Definition: semiimplicit.hh:70
Implicit RungeKutta ODE solver.
Definition: semiimplicit.hh:112
SemiImplicitRungeKuttaSolver(ExplicitOperatorType &explicitOp, HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
Definition: semiimplicit.hh:179
ExplicitOperator ExplicitOperatorType
Definition: semiimplicit.hh:117
BaseType::ParameterType ParameterType
Definition: semiimplicit.hh:123
SemiImplicitRungeKuttaSolver(ExplicitOperatorType &explicitOp, HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, int order, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
Definition: semiimplicit.hh:147
BaseType::TimeStepControlType TimeStepControlType
Definition: semiimplicit.hh:119
SemiImplicitRungeKuttaSolver(ExplicitOperatorType &explicitOp, HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const ParameterType &tscParams, const NonlinearSolverParameterType &nlsParams)
constructor
Definition: semiimplicit.hh:167
SemiImplicitRungeKuttaSolver(ExplicitOperatorType &explicitOp, HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const SimpleButcherTable< double > &explButcherTable, const SimpleButcherTable< double > &implButcherTable, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
Definition: semiimplicit.hh:191
HelmholtzOperator HelmholtzOperatorType
Definition: semiimplicit.hh:118
BaseType::NonlinearSolverParameterType NonlinearSolverParameterType
Definition: semiimplicit.hh:124
static SimpleButcherTable< double > defaultButcherTable(int order, bool expl)
Definition: semiimplicit.hh:204
SemiImplicitRungeKuttaSolver(ExplicitOperatorType &explicitOp, HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, int order, const ParameterType &tscParams, const NonlinearSolverParameterType &nlsParams)
constructor
Definition: semiimplicit.hh:135
BaseType::SourceTermType SourceTermType
Definition: semiimplicit.hh:120
TimeStepControlType::TimeProviderType TimeProviderType
Definition: semiimplicit.hh:122
Definition: timestepcontrol.hh:23
Definition: timestepcontrol.hh:166
general base for time providers
Definition: timeprovider.hh:36