1#ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_SEMIIMPLICIT_HH
2#define DUNE_FEM_SOLVER_RUNGEKUTTA_SEMIIMPLICIT_HH
8#include <dune/common/dynmatrix.hh>
9#include <dune/common/dynvector.hh>
10#include <dune/common/exceptions.hh>
23 template<
class ExplicitOperator >
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() )
44 Dune::DynamicMatrix< double > Ainv( implicitA );
46 alpha_.rightmultiply( Ainv );
48 for(
int i = 0; i < butcherTable.stages(); ++i )
51 for(
int j = 0; j < i; ++j )
52 gamma_[ i ] -= alpha_[ i ][ j ];
56 bool operator() (
double time,
double timeStepSize,
int stage,
57 const DestinationType &u,
const std::vector< DestinationType * > &update,
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 );
75 explicitOp_.setTime( time );
77 uex_.assign( update );
79 explicitOp_.limit( uex_, update );
85 explicitOp_.setTime( time );
86 explicitOp_.initializeTimeStepSize( u );
87 return explicitOp_.timeStepEstimate();
92 return explicitOp_.timeStepEstimate();
97 Dune::DynamicMatrix< double > alpha_;
98 Dune::DynamicVector< double > gamma_, c_;
100 const bool limiter_ ;
109 template<
class ExplicitOperator,
class HelmholtzOperator,
class NonlinearSolver >
111 :
public BasicImplicitRungeKuttaSolver< HelmholtzOperator, NonlinearSolver, ImplicitRungeKuttaTimeStepControl, SemiImplicitRungeKuttaSourceTerm< ExplicitOperator > >
200 SourceTermType( explicitOp, explButcherTable, implButcherTable.A() ),
217 DUNE_THROW( NotImplemented,
"Semi-implicit Runge-Kutta method of order " << order <<
" not implemented." );
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
NonlinearSolver::ParameterType NonlinearSolverParameterType
Definition: basicimplicit.hh:70
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 ¶meter=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 ¶meter=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 ¶meter=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