1 #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_BASICROW_HH
2 #define DUNE_FEM_SOLVER_RUNGEKUTTA_BASICROW_HH
14 #include <dune/common/dynmatrix.hh>
15 #include <dune/common/dynvector.hh>
30 bool operator() (
double time,
double timeStepSize,
int stage,
const T &u,
const std::vector< T * > &update, T &source )
36 void limit( T& update,
const double time ) {}
98 template<
class HelmholtzOperator,
class NonlinearSolver,
class TimeStepControl,
class SourceTerm = NoROWRungeKuttaSourceTerm >
128 template<
class ButcherTable >
131 const ButcherTable &butcherTable,
141 alpha_( butcherTable.A() ),
145 c_( butcherTable.c() ),
149 linAbsTol_( parameter.linAbsTolParameter( ) ),
155 setup( butcherTable );
158 template<
class ButcherTable >
161 const ButcherTable &butcherTable,
175 template<
class ButcherTable >
178 const ButcherTable &butcherTable,
185 template<
class ButcherTable >
188 const ButcherTable &butcherTable,
200 template<
class ButcherTable >
203 const ButcherTable &butcherTable,
209 template<
class ButcherTable >
212 const ButcherTable &butcherTable,
218 template<
class ButcherTable >
219 void setup(
const ButcherTable& butcherTable )
223 std::cout <<
"ROW method of order=" << butcherTable.order() <<
" with " <<
stages_ <<
" stages" << std::endl;
227 for(
int i = 0; i <
stages(); ++i )
229 std::ostringstream name;
230 name <<
"RK stage " << i;
235 for(
int i = 0; i <
stages(); ++i )
246 for(
int i = 0; i <
stages(); ++i )
257 const double helmholtzEstimate =
helmholtzOp_.timeStepEstimate();
259 double sourceTermEstimate =
sourceTerm_.initialTimeStepEstimate( time, U0 );
261 if( sourceTermEstimate < 0.0 ) sourceTermEstimate = helmholtzEstimate ;
263 timeStepControl_.initialTimeStepSize( helmholtzEstimate, sourceTermEstimate );
273 typename HelmholtzOperatorType::JacobianOperatorType jOp(
"jacobianOperator", U.space(), U.space() );
277 assert( timeStepSize > 0.0 );
279 for(
int s = 0; s <
stages(); ++s )
285 for(
int k = 0; k < s; ++k )
288 updateStage *= (
gamma_[s]*timeStepSize);
289 for(
int k = 0; k < s; ++k )
292 rhs_.assign( updateStage );
295 const double stageTime = time +
c_[ s ]*timeStepSize;
310 jInv(
rhs_, updateStage );
311 monitor.linearSolverIterations_ += jInv.iterations();
316 jInv(
rhs_, updateStage );
317 monitor.linearSolverIterations_ += jInv.iterations();
331 for(
int s = 0; s <
stages(); ++s )
335 Uerr.axpy( -1.0, U );
336 const double errorU = Uerr.scalarProductDofs( Uerr );
337 const double normU = U.scalarProductDofs( U );
339 if( normU > 0 && errorU > 0 )
343 std::cout << std::scientific <<
"Error in RK = " << error <<
" norm " << errorU <<
" " << normU << std::endl;
349 for(
int s = 0; s <
stages(); ++s )
353 monitor.error_ = error;
363 out <<
"Generic " <<
stages() <<
"-stage implicit Runge-Kutta solver.\\\\" << std::endl;
370 typedef typename DestinationType :: ConstDofIteratorType ConstDofIteratorType ;
371 const ConstDofIteratorType uend = U.dend();
373 for( ConstDofIteratorType u = U.dbegin(), uerr = Uerr.dbegin(); u != uend; ++u, ++uerr )
376 double uerrval = *uerr ;
379 double norm =
std::abs( uval - uerrval );
Double abs(const Double &a)
Definition: double.hh:872
static double sqrt(const Double &v)
Definition: double.hh:887
static constexpr T max(T a)
Definition: utility.hh:77
Definition: multistep.hh:17
T getValue(const std::string &key) const
get mandatory parameter
Definition: reader.hh:149
static bool verbose()
obtain the cached value for fem.verbose
Definition: io/parameter.hh:437
static DUNE_EXPORT ParameterContainer & container()
Definition: io/parameter.hh:192
Definition: io/parameter.hh:544
Interface class for ODE Solver.
Definition: odesolverinterface.hh:21
Monitor MonitorType
monitor type
Definition: odesolverinterface.hh:59
virtual void solve(DestinationType &u)
solve where is the internal operator.
Definition: odesolverinterface.hh:75
DestinationImp DestinationType
type of destination
Definition: odesolverinterface.hh:62
Definition: odesolverinterface.hh:27
Definition: basicrow.hh:28
double timeStepEstimate() const
Definition: basicrow.hh:45
void limit(T &update, const double time)
Definition: basicrow.hh:36
double initialTimeStepEstimate(double time, const T &u) const
Definition: basicrow.hh:39
bool operator()(double time, double timeStepSize, int stage, const T &u, const std::vector< T * > &update, T &source)
Definition: basicrow.hh:30
Definition: basicrow.hh:54
virtual int maxLinearIterationsParameter() const
Definition: basicrow.hh:89
const Dune::Fem::ParameterReader & parameter() const
Definition: basicrow.hh:72
const std::string keyPrefix_
Definition: basicrow.hh:57
ROWSolverParameter(const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: basicrow.hh:67
virtual double linAbsTolParameter() const
Definition: basicrow.hh:74
Dune::Fem::ParameterReader parameter_
Definition: basicrow.hh:59
virtual bool linearSolverVerbose() const
Definition: basicrow.hh:84
ROWSolverParameter(const std::string keyPrefix, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: basicrow.hh:62
virtual double linReductionParameter() const
Definition: basicrow.hh:79
ROW RungeKutta ODE solver.
Definition: basicrow.hh:101
void solve(DestinationType &U, MonitorType &monitor)
solve the system
Definition: basicrow.hh:269
Dune::DynamicVector< double > beta_
Definition: basicrow.hh:395
BasicROWRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const ButcherTable &butcherTable, const TimeStepControlType &timeStepControl, const ParametersType ¶meter, const NonlinearSolverParametersType &nlsParam)
constructor
Definition: basicrow.hh:176
void description(std::ostream &out) const
print description of ODE solver to out stream
Definition: basicrow.hh:361
SourceTerm SourceTermType
Definition: basicrow.hh:112
int stages_
Definition: basicrow.hh:392
TimeStepControl TimeStepControlType
Definition: basicrow.hh:111
SourceTerm sourceTerm_
Definition: basicrow.hh:390
std::vector< DestinationType * > update_
Definition: basicrow.hh:398
NonlinearSolver::ParametersType NonlinearSolverParametersType
Definition: basicrow.hh:115
BasicROWRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const ButcherTable &butcherTable, const TimeStepControlType &timeStepControl, const SourceTermType &sourceTerm, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: basicrow.hh:159
Dune::DynamicMatrix< double > alpha2_
Definition: basicrow.hh:394
HelmholtzOperatorType & helmholtzOp_
Definition: basicrow.hh:387
Dune::DynamicMatrix< double > alpha_
Definition: basicrow.hh:394
const PreconditionOperatorType * preconditioner_
Definition: basicrow.hh:404
BaseType::MonitorType MonitorType
Definition: basicrow.hh:106
BaseType::DestinationType DestinationType
Definition: basicrow.hh:107
NonlinearSolverType nonlinearSolver_
Definition: basicrow.hh:388
int stages() const
Definition: basicrow.hh:359
ROWSolverParameter ParametersType
Definition: basicrow.hh:114
const double linAbsTol_
Definition: basicrow.hh:400
const double linReduction_
Definition: basicrow.hh:400
TimeStepControl timeStepControl_
Definition: basicrow.hh:389
DestinationType temp_
Definition: basicrow.hh:397
NonlinearSolver NonlinearSolverType
Definition: basicrow.hh:110
BasicROWRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const ButcherTable &butcherTable, const TimeStepControlType &timeStepControl, const SourceTermType &sourceTerm, const ParametersType ¶meter, const NonlinearSolverParametersType &nlsParam)
constructor
Definition: basicrow.hh:129
double delta_
Definition: basicrow.hh:393
const bool linVerbose_
Definition: basicrow.hh:401
const int maxLinearIterations_
Definition: basicrow.hh:402
Dune::Fem::TimeProviderBase TimeProviderType
Definition: basicrow.hh:119
void initialize(const DestinationType &U0)
apply operator once to get dt estimate
Definition: basicrow.hh:251
BasicROWRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const ButcherTable &butcherTable, const ParametersType ¶meter, const NonlinearSolverParametersType &nlsParam)
constructor
Definition: basicrow.hh:201
HelmholtzOperator::SpaceOperatorType::PreconditionOperatorType PreconditionOperatorType
Definition: basicrow.hh:117
HelmholtzOperator HelmholtzOperatorType
Definition: basicrow.hh:109
double infNorm(const DestinationType &U, const DestinationType &Uerr) const
Definition: basicrow.hh:368
BasicROWRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const ButcherTable &butcherTable, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: basicrow.hh:210
BasicROWRungeKuttaSolver(HelmholtzOperatorType &helmholtzOp, TimeProviderType &timeProvider, const ButcherTable &butcherTable, const TimeStepControlType &timeStepControl, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: basicrow.hh:186
Dune::DynamicVector< double > gamma_
Definition: basicrow.hh:395
DestinationType rhs_
Definition: basicrow.hh:397
Dune::DynamicVector< double > c_
Definition: basicrow.hh:395
~BasicROWRungeKuttaSolver()
destructor
Definition: basicrow.hh:244
void setup(const ButcherTable &butcherTable)
Definition: basicrow.hh:219
Definition: timestepcontrol.hh:166
general base for time providers
Definition: timeprovider.hh:36