dune-fem  2.4.1-rc
odesolver.hh
Go to the documentation of this file.
1 #ifndef RUNGEKUTTA_ODE_SOLVER_HH
2 #define RUNGEKUTTA_ODE_SOLVER_HH
3 
4 // inlcude all used headers before, that they don not appear in DuneODE
5 #ifdef HUNDERT
6 #endif
7 //- system includes
8 #include <iostream>
9 #include <cmath>
10 #include <vector>
11 #include <pthread.h>
12 #include <cassert>
13 #include <sys/times.h>
14 
15 //- dune-common includes
16 #include <dune/common/nullptr.hh>
17 
18 //- dune-fem includes
19 #include <dune/fem/io/parameter.hh>
25 
26 // include headers of PARDG
27 #include "pardg.hh"
28 
29 namespace DuneODE
30 {
31 
32  using namespace Dune;
33  using namespace Fem;
34  using namespace std;
35 
36 #ifdef USE_PARDG_ODE_SOLVER
37  // see rungekutta/timestepcontrol.hh for the adaptive control of the clf for implicit solvers
39  {
42 
44  : ImplicitRungeKuttaSolverParameters( "fem.ode.", parameter )
45  {}
46 
47  ODEParameters( const std::string keyPrefix, const ParameterReader &parameter = Parameter::container() )
48  : ImplicitRungeKuttaSolverParameters( keyPrefix, parameter )
49  {}
50 
51  // choice of linear solver for the implicit ODE solver
52  virtual PARDG::IterativeLinearSolver *linearSolver(PARDG::Communicator & comm) const
53  {
54  static const std::string methodTypeTable[] = { "gmres", "cg", "bicgstab" };
55  const int method = parameter_.getEnum( keyPrefix_ + "linearsolver", methodTypeTable, 0 );
56 
57  PARDG::IterativeLinearSolver *solver = nullptr;
58  switch( method )
59  {
60  case 0:
61  int cycles;
62  if( parameter_.exists( "fem.ode.gmrescycles" ) )
63  {
64  std::cerr << "WARNING: deprecated key, use '" << keyPrefix_ << "gmres.cycles' instead!" << std::endl;
65  cycles = parameter_.getValue< int >( "fem.ode.gmrescycles", 15 );
66  }
67  else
68  cycles = parameter_.getValue< int >( keyPrefix_ + "gmres.cycles", 15 );
69  solver = new PARDG::GMRES( comm, cycles );
70  break;
71 
72  case 1:
73  solver = new PARDG::CG( comm );
74  break;
75 
76  case 2:
77  solver = new PARDG::BICGSTAB( comm );
78  break;
79  }
80  if( !solver )
81  DUNE_THROW( InvalidStateException, "Unable to create linear solver." );
82 
83  // tolerance for the linear solver
84  double tol = parameter_.getValue< double >( keyPrefix_ + "solver.tolerance" , 1e-8 );
85  static const std::string errorTypeTable[] = { "absolute", "relative" };
86  // errormeassure used in the linear solver
87  int errorType = parameter_.getEnum( keyPrefix_ + "solver.errormeasure", errorTypeTable, 0 );
88  solver->set_tolerance(tol,(errorType==1));
89  // max iterations that the linear solver should do
90  int maxIter = parameter_.getValue< int >( keyPrefix_ + "solver.iterations" , 1000 );
91  solver->set_max_number_of_iterations(maxIter);
92  return solver;
93  }
94 
95  // overload clone method
96  virtual ODEParameters* clone () const { return new ODEParameters( *this ); }
97  };
98 
99 
100  template <class Operator>
101  class OperatorWrapper : public PARDG::Function
102  {
103  // type of discrete function
104  typedef typename Operator::DestinationType DestinationType;
105  public:
108  : op_(op)
109  {}
110 
113  void operator()(const double *u, double *f, int i = 0)
114  {
115  // set actual time of iteration step
116  op_.setTime( this->time() );
117 
118  // call operator apply
119  op_( u, f );
120  }
121 
123  int dim_of_argument(int i = 0) const
124  {
125  assert( i == 0 );
126  return op_.size();
127  }
128 
130  int dim_of_value(int i = 0) const
131  {
132  assert( i == 0 );
133  return op_.size();
134  }
135 
137  const Operator& op() const { return op_; }
138 
139  protected:
140  // operator to call
142  };
143 
144 
145  template <class Operator>
146  class LimiterWrapper : public PARDG::Limiter
147  {
148  // type of discrete function
149  typedef typename Operator::DestinationType DestinationType;
150  public:
153  : op_(op)
154  {}
155 
156  void operator()( double* ) { abort(); }
157 
160  void operator()(const double *u, double *f)
161  {
162  // call operator apply
163  op_.limit( u, f );
164  }
165 
167  const Operator& op() const { return op_; }
168 
169  protected:
170  // operator to call
172  };
173 
174 
180  /* \brief Base class for ParDG ode solvers */
181  template<class DestinationImp>
182  class ParDGOdeSolverBase : public OdeSolverInterface<DestinationImp>
183  {
185  public:
187  typedef DestinationImp DestinationType;
188 
191 
194  protected:
197  const int order ) :
198  timeProvider_( tp ),
199  comm_(PARDG::Communicator::instance()),
200  order_( order ),
201  initialized_( false ),
202  odeSolver_( 0 )
203  {
204  }
205 
207  {
208  if( !initialized_ )
209  {
210  if( !odeSolver_ )
211  odeSolver_ = createOdeSolver();
212  initialized_ = true;
213  }
214  }
215 
216  public:
217  // initialize time step size
218  void initialize ( const DestinationType& U0 )
219  {
220  // initialized dt on first call
221  if ( ! initialized_ )
222  {
223  // create instance of ode solver
224  initializeOdeSolver();
225 
226  // apply operator once
227  spaceOperator().initializeTimeStepSize( U0 );
228 
229  // set initial time step
230  timeProvider_.provideTimeStepEstimate( spaceOperator().timeStepEstimate() );
231  }
232  }
233 
234  protected:
236  ~ParDGOdeSolverBase() { delete odeSolver_; odeSolver_ = 0; }
237 
239  PARDG::ODESolver& odeSolver()
240  {
241  assert( odeSolver_ );
242  return *odeSolver_;
243  }
244 
245  public:
246  void description(std::ostream& out) const
247  {
248  out << name() << ", steps: " << order_;
249  // would be nice to add info on cfl number, but TimeProviderBase doesn't have this
250  //<< ", cfl: " << this->timeProvider_.cfl() << "\\\\" <<std::endl;
251  }
252 
253  protected:
255  virtual std::string name() const = 0;
256 
258  virtual PARDG::ODESolver* createOdeSolver() = 0;
259 
261  virtual const OperatorType& spaceOperator() const = 0;
262 
263  protected:
265  PARDG::Communicator & comm_;
266  const int order_;
268 
269  private:
270  // use method odeSolver for access
271  PARDG::ODESolver* odeSolver_;
272  };
273 
275  //
276  // --ExplicitOdeSolver
277  //
279  template<class DestinationImp>
281  public ParDGOdeSolverBase< DestinationImp >
282  {
284 
285  protected:
286  using BaseType :: order_;
287  using BaseType :: comm_;
288  using BaseType :: initialized_;
289  using BaseType :: timeProvider_;
290  using BaseType :: odeSolver ;
291 
292  public:
293  using BaseType :: initialize;
294  using BaseType :: solve ;
295 
299 
301  ExplicitOdeSolver(OperatorType& op,
302  TimeProviderBase &tp,
303  const int order, bool verbose = false) :
304  BaseType(tp, order ),
305  expl_( op ),
306  verbose_( verbose )
307  {}
308 
310  virtual ~ExplicitOdeSolver() {}
311 
313  void solve( DestinationType& U0,
314  MonitorType& monitor )
315  {
316  // no nonlinear system to solve for time step update
317  monitor.reset() ;
318 
319  // initialize
320  if( ! initialized_ )
321  {
322  DUNE_THROW(InvalidStateException,"ExplicitOdeSolver wasn't initialized before first call!");
323  }
324 
325  // get dt
326  const double dt = timeProvider_.deltaT();
327 
328  // should be larger then zero
329  assert( dt > 0.0 );
330 
331  // get time
332  const double time = timeProvider_.time();
333 
334  // get leakPointer
335  double* u = U0.leakPointer();
336 
337  // call ode solver
338  const bool convergence = odeSolver().step(time, dt, u);
339 
340  // set time step estimate of operator
341  timeProvider_.provideTimeStepEstimate( expl_.op().timeStepEstimate() );
342 
343  assert( convergence );
344  if(! convergence )
345  {
346  timeProvider_.invalidateTimeStep();
347  if( comm_.id() == 0 )
348  std::cerr << "No Convergence of ExplicitOdeSolver! \n";
349  }
350  }
351 
352  protected:
354  virtual std::string name() const { return "ExplicitOdeSolver"; }
355 
357  const OperatorType& spaceOperator() const { return expl_.op(); }
358 
360  PARDG::ODESolver* createOdeSolver()
361  {
362  PARDG::ODESolver* odeSolver = 0;
363  switch ( order_ )
364  {
365  case 1 : odeSolver = new PARDG::ExplicitEuler(comm_, expl_); break;
366  case 2 : odeSolver = new PARDG::ExplicitTVD2 (comm_, expl_); break;
367  case 3 : odeSolver = new PARDG::ExplicitTVD3 (comm_, expl_); break;
368  case 4 : odeSolver = new PARDG::ExplicitRK4 (comm_, expl_); break;
369  default : odeSolver = new PARDG::ExplicitBulirschStoer(comm_, expl_, 7);
370  if( comm_.id() == 0 )
371  {
372  std::cerr << "Runge-Kutta method of this order not implemented.\n"
373  << "Using 7-stage Bulirsch-Stoer scheme.\n"
374  << std::endl;
375  }
376  }
377 
378  // only set output when general verbose mode is enabled
379  // (basically to avoid output on every rank)
380  if( verbose_ && Parameter :: verbose() )
381  {
382  odeSolver->DynamicalObject::set_output(cout);
383  }
384  assert( odeSolver );
385  return odeSolver;
386  }
387 
388  protected:
390  const bool verbose_;
391  };
392 
394  //
395  // --ImplicitOdeSolver
396  //
398  template<class DestinationImp>
400  public ParDGOdeSolverBase< DestinationImp >
401  {
403  protected:
404  using BaseType :: order_;
405  using BaseType :: comm_;
406  using BaseType :: initialized_;
407  using BaseType :: timeProvider_;
408  using BaseType :: odeSolver ;
409 
410  public:
411  using BaseType :: initialize;
412  using BaseType :: solve ;
413 
417 
418  ImplicitOdeSolver(OperatorType& op,
419  TimeProviderBase& tp,
420  const int order,
421  const ODEParameters& parameter ) :
422  BaseType(tp, order),
423  impl_( op ),
424  linsolver_( 0 ),
425  param_( parameter.clone() ),
426  verbose_( param_->verbose() ),
427  cfl_( param_->cflStart() ),
428  cflMax_( param_->cflMax() )
429  {
430  }
431 
432  ImplicitOdeSolver(OperatorType& op,
433  TimeProviderBase& tp,
434  const int order,
435  const ParameterReader &parameter = Parameter::container() ) :
436  BaseType(tp, order),
437  impl_( op ),
438  linsolver_( 0 ),
439  param_( new ODEParameters( parameter ) ),
440  verbose_( param_->verbose() ),
441  cfl_( param_->cflStart() ),
442  cflMax_( param_->cflMax() )
443  {
444  }
445 
447  {
448  delete linsolver_; linsolver_ = 0;
449  delete param_; param_ = 0;
450  }
451 
452  protected:
454  virtual std::string name() const { return "ImplicitOdeSolver"; }
455 
457  const ODEParameters& parameter() const { assert( param_ ); return *param_; }
458 
460  virtual const OperatorType& spaceOperator() const { return impl_.op(); }
461 
463  PARDG::ODESolver* createOdeSolver()
464  {
465  linsolver_ = parameter().linearSolver( comm_ );
466  assert( linsolver_ );
467 
468  PARDG::DIRK* odeSolver = 0;
469  switch ( order_ )
470  {
471  case 1: odeSolver = new PARDG::ImplicitEuler(comm_, impl_); break;
472  case 2: odeSolver = new PARDG::Gauss2(comm_, impl_); break;
473  case 3: odeSolver = new PARDG::DIRK3 (comm_, impl_); break;
474  case 4: odeSolver = new PARDG::DIRK34 (comm_, impl_); break;
475  default : if( comm_.id() == 0 )
476  {
477  std::cerr << "Runge-Kutta method of this order not implemented"
478  << std::endl;
479  }
480  abort();
481  }
482 
483  assert( odeSolver );
484 
485  odeSolver->set_linear_solver( *linsolver_ );
486  odeSolver->set_tolerance( parameter().tolerance() );
487  odeSolver->set_max_number_of_iterations( parameter().iterations() );
488 
489  // only set output when general verbose mode is enabled
490  // (basically to avoid output on every rank)
491  if( verbose_ == ODEParameters :: fullVerbosity &&
493  {
494  odeSolver->IterativeSolver::set_output(cout);
495  odeSolver->DynamicalObject::set_output(cout);
496  linsolver_->IterativeSolver::set_output(cout);
497  }
498  return odeSolver;
499  }
500 
501  virtual int numberOfIterations()
502  {
503  return static_cast<PARDG::DIRK&> (odeSolver()).number_of_iterations();
504  }
505 
506  public:
508  void solve( DestinationType& U0, MonitorType& monitor )
509  {
510  // initialize
511  if( ! initialized_ )
512  {
513  DUNE_THROW(InvalidStateException,"ImplicitOdeSolver wasn't initialized before first call!");
514  }
515 
516  const double dt = timeProvider_.deltaT();
517  assert( dt > 0.0 );
518  const double time = timeProvider_.time();
519 
520  // get pointer to solution
521  double* u = U0.leakPointer();
522 
523  const bool convergence =
524  odeSolver().step( time, dt, u,
525  monitor.newtonIterations_,
526  monitor.linearSolverIterations_,
527  monitor.maxNewtonIterations_,
528  monitor.maxLinearSolverIterations_ );
529 
530  assert( linsolver_->number_of_iterations() == monitor.linearSolverIterations_ );
531 
532  double factor( 1 );
533  bool changed =
534  parameter().cflFactor( impl_.op().timeStepEstimate(), spaceOperator().timeStepEstimate(),
535  monitor.linearSolverIterations_, convergence, factor );
536 
537  if( (factor >= std::numeric_limits< double >::min()) &&
539  {
540  // only apply when factor is small or max cfl was not reached yet
541  if( factor < 1.0 || cfl_ <= cflMax_ )
542  {
543  cfl_ *= factor;
544  }
545  else
546  changed = false ;
547  }
548  else
549  DUNE_THROW( InvalidStateException, "invalid cfl factor: " << factor );
550 
551  if (convergence)
552  {
553  // timeProvider_.provideTimeStepEstimate( cfl_ * spaceOperator().timeStepEstimate() );
554  timeProvider_.provideTimeStepEstimate( timeStepEstimate(cfl_) );
555 
556  if( changed && (verbose_ >= ODEParameters :: cflVerbosity ) && (MPIManager::rank() == 0) )
557  {
558  derr << " New cfl number is: "<< cfl_ << ", iterations per time step("
559  << "ILS: " << monitor.linearSolverIterations_
560  << ", Newton: " << monitor.newtonIterations_
561  << ")"
562  << std::endl;
563  }
564  }
565  else
566  {
567  if( factor >= 1.0 )
568  DUNE_THROW( InvalidStateException, "Invalid cfl factor (no convergence): " << factor );
569  timeProvider_.provideTimeStepEstimate( factor * dt );
570  timeProvider_.invalidateTimeStep();
571 
572  // output only in verbose mode
573  if( (verbose_ >= ODEParameters :: noConvergenceVerbosity) && (MPIManager::rank() == 0) )
574  derr << "No convergence: New cfl number is " << cfl_ << std::endl;
575  }
576 
577  linsolver_->reset_number_of_iterations();
578  }
579 
580  protected:
581  virtual double timeStepEstimate(double cfl)
582  {
583  return cfl*spaceOperator().timeStepEstimate();
584  }
585 
587  PARDG::IterativeLinearSolver* linsolver_;
589  const int verbose_;
590  double cfl_;
591  const double cflMax_;
592  }; // end ImplicitOdeSolver
593 
594 
596  //
597  // --SemiImplicitOdeSolver
598  //
600  template<class DestinationImp>
602  public ImplicitOdeSolver< DestinationImp >
603  {
604  typedef ImplicitOdeSolver< DestinationImp > BaseType;
605  protected:
606  using BaseType :: order_;
607  using BaseType :: comm_;
608  using BaseType :: initialized_;
609  using BaseType :: timeProvider_;
610  using BaseType :: odeSolver ;
611  using BaseType :: linsolver_ ;
612  using BaseType :: parameter ;
613  using BaseType :: impl_ ;
614  using BaseType :: verbose_ ;
615  using BaseType :: cfl_;
616 
617  public:
620 
621  SemiImplicitOdeSolver(OperatorType& explOp,
622  OperatorType& implOp,
623  TimeProviderBase& tp,
624  const int order,
625  const ODEParameters& parameter = ODEParameters() ) :
626  BaseType( implOp, tp, order, parameter ),
627  expl_( explOp ),
628  limiter_( explOp )
629  {
630  }
631 
632  // initialize time step size
633  void initialize ( const DestinationType& U0 )
634  {
635  // initialized dt on first call
636  if( !initialized_ )
637  {
638  BaseType::initializeOdeSolver();
639 
640  // apply operators once
641  expl_.op().initializeTimeStepSize( U0 );
642  impl_.op().initializeTimeStepSize( U0 );
643 
644  // set initial time step
645  double dtEst = std::numeric_limits< double >::max();
646  parameter().initTimeStepEstimate( expl_.op().timeStepEstimate(), impl_.op().timeStepEstimate(), dtEst, cfl_ );
647  timeProvider_.provideTimeStepEstimate( dtEst );
648  }
649  }
650 
651 
652  protected:
654  virtual std::string name() const { return "SemiImplicitOdeSolver"; }
655 
657  const OperatorType& spaceOperator() const { return expl_.op(); }
658 
660  PARDG::ODESolver* createOdeSolver()
661  {
662  linsolver_ = parameter().linearSolver( comm_ );
663  assert( linsolver_ );
664 
665  PARDG::SIRK* odeSolver = 0;
666  switch ( order_ )
667  {
668  case 1: odeSolver = new PARDG::SemiImplicitEuler(comm_, impl_, expl_); break;
669  case 2: odeSolver = new PARDG::IMEX_SSP222 (comm_, impl_, expl_); break;
670  case 3: odeSolver = new PARDG::SIRK33 (comm_, impl_, expl_); break;
671  //odeSolver = new PARDG::IMEX_ARK34 (comm_, impl_, expl_); break;
672  case 4: odeSolver = new PARDG::IERK45 (comm_, impl_, expl_); break;
673  //odeSolver = new PARDG::IMEX_ARK46 (comm_, impl_, expl_); break;
674  default : if( comm_.id() == 0 )
675  {
676  std::cerr << "Runge-Kutta method of this order not implemented"
677  << std::endl;
678  }
679  abort();
680  }
681 
682  assert( odeSolver );
683 
684  odeSolver->set_linear_solver(*linsolver_);
685  odeSolver->set_tolerance( parameter().tolerance() );
686  odeSolver->set_max_number_of_iterations( parameter().iterations() );
687 
688  if( expl_.op().hasLimiter() )
689  odeSolver->set_expl_limiter( limiter_ );
690 
691  // only set output when general verbose mode is enabled
692  // (basically to avoid output on every rank)
693  if( verbose_ == ODEParameters :: fullVerbosity &&
695  {
696  odeSolver->IterativeSolver::set_output(cout);
697  odeSolver->DynamicalObject::set_output(cout);
698  // linsolver_->set_output(cout);
699  }
700  return odeSolver;
701  }
702 
703  virtual int numberOfIterations()
704  {
705  return static_cast<PARDG::SIRK&> (odeSolver()).number_of_iterations();
706  }
707 
708  protected:
709  virtual double timeStepEstimate(double cfl)
710  {
711  // take the minimum of the explicit part and the implicit part scaled by the cfl number
712  // return spaceOperator().timeStepEstimate();
713  return std::min( spaceOperator().timeStepEstimate(),
714  cfl * impl_.op().timeStepEstimate() );
715  }
718 
719  }; // end SemiImplicitOdeSolver
720 
721 #endif // USE_PARDG_ODE_SOLVER
722 
727 } // namespace DuneODE
728 
729 #undef USE_EXTERNAL_BLAS
730 #endif // #ifndef RUNGEKUTTA_ODE_SOLVER_HH
PARDG::ODESolver & odeSolver()
return reference to ode solver
Definition: odesolver.hh:239
ODESpaceOperatorInterface for Operators that work with PARDG ODE solvers of the type where is a dis...
Definition: spaceoperatorif.hh:28
const ODEParameters & parameter() const
return reference to parameter class
Definition: odesolver.hh:457
Definition: odesolver.hh:146
const Operator & op() const
return reference to real operator
Definition: odesolver.hh:137
virtual ODEParameters * clone() const
Definition: odesolver.hh:96
PARDG::ODESolver * createOdeSolver()
create implicit ode solver
Definition: odesolver.hh:660
virtual std::string name() const
return name of ode solver
Definition: odesolver.hh:354
void initialize(const DestinationType &U0)
initialize solver
Definition: odesolver.hh:218
Interface class for ODE Solver.
Definition: odesolverinterface.hh:11
const int verbose_
Definition: odesolver.hh:589
static int rank()
Definition: mpimanager.hh:116
OperatorWrapper(Operator &op)
constructor
Definition: odesolver.hh:107
static constexpr T min(T a)
Definition: utility.hh:81
Definition: odesolver.hh:38
~ParDGOdeSolverBase()
destructor
Definition: odesolver.hh:236
BaseType::OperatorType OperatorType
Definition: odesolver.hh:296
ImplicitOdeSolver(OperatorType &op, TimeProviderBase &tp, const int order, const ParameterReader &parameter=Parameter::container())
Definition: odesolver.hh:432
PARDGSpaceOperatorInterface< DestinationType > OperatorType
type of discretization operator
Definition: odesolver.hh:190
virtual const OperatorType & spaceOperator() const
for initialization
Definition: odesolver.hh:460
void initialize(const DestinationType &U0)
initialize solver
Definition: odesolver.hh:633
DestinationImp DestinationType
type of destination function
Definition: odesolver.hh:187
static constexpr T max(T a)
Definition: utility.hh:65
int linearSolverIterations_
Definition: odesolverinterface.hh:26
virtual int numberOfIterations()
Definition: odesolver.hh:501
OperatorWrapper< OperatorType > expl_
Definition: odesolver.hh:716
const OperatorType & spaceOperator() const
for initialization
Definition: odesolver.hh:357
bool initialized_
Definition: odesolver.hh:267
PARDG::IterativeLinearSolver * linsolver_
Definition: odesolver.hh:587
Definition: odesolver.hh:182
LimiterWrapper(Operator &op)
constructor
Definition: odesolver.hh:152
BaseType::OperatorType OperatorType
Definition: odesolver.hh:414
virtual int numberOfIterations()
Definition: odesolver.hh:703
Definition: multistep.hh:16
Definition: odesolver.hh:280
BaseType::MonitorType MonitorType
Definition: odesolver.hh:416
virtual PARDG::IterativeLinearSolver * linearSolver(PARDG::Communicator &comm) const
Definition: odesolver.hh:52
Operator & op_
Definition: odesolver.hh:171
ODEParameters(const std::string keyPrefix, const ParameterReader &parameter=Parameter::container())
Definition: odesolver.hh:47
void reset()
Definition: odesolverinterface.hh:34
LimiterWrapper< OperatorType > limiter_
Definition: odesolver.hh:717
OperatorWrapper< OperatorType > impl_
Definition: odesolver.hh:586
BaseType::DestinationType DestinationType
Definition: odesolver.hh:297
int maxLinearSolverIterations_
Definition: odesolverinterface.hh:28
virtual double timeStepEstimate(double cfl)
Definition: odesolver.hh:581
int dim_of_value(int i=0) const
return size of destination
Definition: odesolver.hh:130
BaseType::OperatorType OperatorType
Definition: odesolver.hh:618
Definition: odesolver.hh:101
general base for time providers
Definition: timeprovider.hh:35
const double cflMax_
Definition: odesolver.hh:591
abstract operator
Definition: operator.hh:25
const Operator & op() const
return reference to real operator
Definition: odesolver.hh:167
Definition: odesolver.hh:399
virtual double timeStepEstimate(double cfl)
Definition: odesolver.hh:709
Dune::Fem::ParameterReader parameter_
Definition: timestepcontrol.hh:40
virtual std::string name() const
return name of ode solver
Definition: odesolver.hh:454
Definition: coordinate.hh:4
PARDG::ODESolver * createOdeSolver()
create implicit ode solver
Definition: odesolver.hh:463
ExplicitOdeSolver(OperatorType &op, TimeProviderBase &tp, const int order, bool verbose=false)
constructor
Definition: odesolver.hh:301
Operator & op_
Definition: odesolver.hh:141
void initializeOdeSolver()
Definition: odesolver.hh:206
ODEParameters(const ParameterReader &parameter=Parameter::container())
Definition: odesolver.hh:43
const OperatorType & spaceOperator() const
for initialization
Definition: odesolver.hh:657
BaseType::DestinationType DestinationType
Definition: odesolver.hh:619
int dim_of_argument(int i=0) const
return size of argument
Definition: odesolver.hh:123
static ParameterContainer & container()
Definition: io/parameter.hh:190
BaseType::DestinationType DestinationType
Definition: odesolver.hh:415
STL namespace.
BaseType::MonitorType MonitorType
Definition: odesolver.hh:298
const int order_
Definition: odesolver.hh:266
PARDG::ODESolver * createOdeSolver()
create explicit ode solver
Definition: odesolver.hh:360
void solve(DestinationType &U0, MonitorType &monitor)
solve system
Definition: odesolver.hh:313
ImplicitOdeSolver(OperatorType &op, TimeProviderBase &tp, const int order, const ODEParameters &parameter)
Definition: odesolver.hh:418
static bool verbose()
obtain the cached value for fem.verbose
Definition: io/parameter.hh:444
int newtonIterations_
Definition: odesolverinterface.hh:25
double cfl_
Definition: odesolver.hh:590
OdeSolverInterface< DestinationImp >::MonitorType MonitorType
type of monitor class
Definition: odesolver.hh:193
int maxNewtonIterations_
Definition: odesolverinterface.hh:27
Definition: odesolverinterface.hh:17
const std::string keyPrefix_
Definition: timestepcontrol.hh:29
TimeProviderBase & timeProvider_
Definition: odesolver.hh:264
Definition: odesolver.hh:601
virtual ~ImplicitOdeSolver()
Definition: odesolver.hh:446
ParDGOdeSolverBase(TimeProviderBase &tp, const int order)
constructor
Definition: odesolver.hh:196
void operator()(double *)
Definition: odesolver.hh:156
PARDG::Communicator & comm_
Definition: odesolver.hh:265
void solve(DestinationType &U0, MonitorType &monitor)
solve
Definition: odesolver.hh:508
void operator()(const double *u, double *f, int i=0)
Definition: odesolver.hh:113
Definition: timestepcontrol.hh:21
virtual std::string name() const
return name of ode solver
Definition: odesolver.hh:654
const ODEParameters * param_
Definition: odesolver.hh:588
const bool verbose_
Definition: odesolver.hh:390
void operator()(const double *u, double *f)
Definition: odesolver.hh:160
#define PARDG
Definition: pardg.hh:26
SemiImplicitOdeSolver(OperatorType &explOp, OperatorType &implOp, TimeProviderBase &tp, const int order, const ODEParameters &parameter=ODEParameters())
Definition: odesolver.hh:621
virtual ~ExplicitOdeSolver()
destructor
Definition: odesolver.hh:310
OperatorWrapper< OperatorType > expl_
Definition: odesolver.hh:389
void description(std::ostream &out) const
print description of ODE solver to out stream
Definition: odesolver.hh:246