1 #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_EXPLICIT_HH 2 #define DUNE_FEM_SOLVER_RUNGEKUTTA_EXPLICIT_HH 61 template<
class DestinationImp>
68 typedef typename DestinationType :: DiscreteFunctionSpaceType
SpaceType;
74 std::vector< std::vector<double> >
a;
75 std::vector<double>
b;
76 std::vector<double>
c;
77 std::vector<DestinationType*>
Upd;
91 a(0),b(0),c(0), Upd(0),
99 for (
int i=0; i<ord_; ++i)
109 a[0][0]=0.; a[0][1]=0.; a[0][2]=0.; a[0][3]=0.;
110 a[1][0]=1.0; a[1][1]=0.; a[1][2]=0.; a[1][3]=0.;
111 a[2][0]=0.25; a[2][1]=0.25; a[2][2]=0.; a[2][3]=0.;
112 a[3][0]=1./6.; a[3][1]=1./6.; a[3][2]=2./3.; a[3][3]=0.;
113 b[0]=1./6.; b[1]=1./6.; b[2]=2./3.; b[3]=0.;
114 c[0]=0.; c[1]=1.0; c[2]=0.5; c[3]=1.0;
117 a[0][0]=0.; a[0][1]=0.; a[0][2]=0.;
118 a[1][0]=1.0; a[1][1]=0.; a[1][2]=0.;
119 a[2][0]=0.25; a[2][1]=0.25; a[2][2]=0.;
120 b[0]=1./6.; b[1]=1./6.; b[2]=2./3.;
121 c[0]=0.; c[1]=1; c[2]=0.5;
124 a[0][0]=0.; a[0][1]=0.;
125 a[1][0]=1.0; a[1][1]=0.;
134 default : std::cerr <<
"Runge-Kutta method of this order not implemented" 140 for (
int i=0; i<ord_; ++i)
142 Upd.push_back(
new DestinationType(
"URK",op_.space()) );
144 Upd.push_back(
new DestinationType(
"Ustep",op_.space()) );
158 for(
size_t i=0; i<Upd.size(); ++i)
172 tp_.provideTimeStepEstimate( op_.timeStepEstimate() );
177 void solve(DestinationType& U0, MonitorType& monitor )
185 DUNE_THROW(InvalidStateException,
"ExplicitRungeKuttaSolver wasn't initialized before first call!");
189 const double dt = tp_.deltaT();
191 const double t = tp_.time();
200 tp_.provideTimeStepEstimate( op_.timeStepEstimate() );
202 for (
int i=1; i<ord_; ++i)
204 (Upd[ord_])->assign(U0);
205 for (
int j=0; j<i ; ++j)
207 (Upd[ord_])->
axpy((a[i][j]*dt), *(Upd[j]));
211 op_.setTime( t + c[i]*dt );
214 op_( *(Upd[ord_]), *(Upd[i]) );
217 tp_.provideTimeStepEstimate( op_.timeStepEstimate() );
221 for (
int j=0; j<ord_; ++j)
223 U0.axpy((b[j]*dt), *(Upd[j]));
229 out <<
"ExplRungeKutta, steps: " << ord_
231 <<
"\\\\" <<std::endl;
247 #endif // #ifndef DUNE_FEM_SOLVER_RUNGEKUTTA_EXPLICIT_HH
std::vector< std::vector< double > > a
Definition: explicit.hh:74
TimeProviderBase & tp_
Definition: explicit.hh:238
void description(std::ostream &out) const
print description of ODE solver to out stream
Definition: explicit.hh:227
Interface class for ODE Solver.
Definition: odesolverinterface.hh:11
DestinationImp DestinationType
Definition: explicit.hh:66
void initialize(const DestinationType &U0)
apply operator once to get dt estimate
Definition: explicit.hh:163
~ExplicitRungeKuttaSolver()
destructor
Definition: explicit.hh:156
void solve(DestinationType &U0, MonitorType &monitor)
solve the system
Definition: explicit.hh:177
Definition: multistep.hh:16
void reset()
Definition: odesolverinterface.hh:34
void axpy(const T &a, const T &x, T &y)
Definition: space/basisfunctionset/functor.hh:37
DestinationType::DiscreteFunctionSpaceType SpaceType
Definition: explicit.hh:68
general base for time providers
Definition: timeprovider.hh:35
const int ord_
Definition: explicit.hh:78
Definition: coordinate.hh:4
ExplicitRungeKuttaSolver(OperatorType &op, TimeProviderBase &tp, const int pord, const Dune::Fem::ParameterReader ¶meter=Dune::Fem::Parameter::container())
Definition: explicit.hh:148
OdeSolverInterface< DestinationImp >::MonitorType MonitorType
Definition: explicit.hh:70
static ParameterContainer & container()
Definition: io/parameter.hh:190
SpaceOperatorInterface< DestinationImp > OperatorType
Definition: explicit.hh:67
bool initialized_
Definition: explicit.hh:240
std::vector< double > b
Definition: explicit.hh:75
std::vector< DestinationType * > Upd
Definition: explicit.hh:77
ExplicitRungeKuttaSolver(OperatorType &op, TimeProviderBase &tp, const int pord, bool verbose)
constructor
Definition: explicit.hh:87
Definition: odesolverinterface.hh:17
Exlicit RungeKutta ODE solver.
Definition: explicit.hh:62
std::vector< double > c
Definition: explicit.hh:76
OperatorType & op_
Definition: explicit.hh:236