4 #ifndef DUNE_PDELAB_MULTISTEP_METHOD_HH
5 #define DUNE_PDELAB_MULTISTEP_METHOD_HH
10 #include <dune/common/ios_state.hh>
11 #include <dune/common/shared_ptr.hh>
12 #include <dune/common/timer.hh>
30 template<
class T,
class MGOS,
class PDESolver,
class TrialV,
53 MGOS& mgos_, PDESolver& pdeSolver_)
54 : parameters(¶meters_), mgos(mgos_), pdeSolver(pdeSolver_),
57 if(mgos.trialGridFunctionSpace().gridView().comm().rank() > 0)
64 if(mgos.trialGridFunctionSpace().gridView().comm().rank() == 0)
70 parameters = ¶meters_;
84 template<
class OldValues>
85 T
apply(T time, T dt,
const OldValues& oldValues, TrialV& xnew)
91 ios_base_all_saver format_attribute_saver(std::cout);
94 std::cout <<
"TIME STEP [" << parameters->
name() <<
"] "
95 << std::setw(6) << step
97 << std::setw(12) << std::setprecision(4) << std::scientific
100 << std::setw(12) << std::setprecision(4) << std::scientific
103 << std::setw(12) << std::setprecision(4) << std::scientific
109 std::cout <<
"== prepare assembler" << std::endl;
112 mgos.preStep(time,dt, oldValues);
114 std::cout <<
"== prepare assembler (" << subTimer.elapsed() <<
"s)"
119 std::cout <<
"== apply solver" << std::endl;
122 pdeSolver.apply(xnew);
124 std::cout <<
"== apply solver (" << subTimer.elapsed() <<
"s)"
129 std::cout <<
"== cleanup assembler" << std::endl;
134 std::cout <<
"== cleanup assembler (" << subTimer.elapsed() <<
"s)"
140 std::cout <<
"== time step done (" << allTimer.elapsed() <<
"s)"
161 template<
typename OldValues,
typename F>
162 T
apply(T time, T dt,
const OldValues& oldValues, F&
f, TrialV& xnew)
168 ios_base_all_saver format_attribute_saver(std::cout);
171 std::cout <<
"TIME STEP [" << parameters->
name() <<
"] "
172 << std::setw(6) << step
174 << std::setw(12) << std::setprecision(4) << std::scientific
177 << std::setw(12) << std::setprecision(4) << std::scientific
180 << std::setw(12) << std::setprecision(4) << std::scientific
186 std::cout <<
"== prepare assembler" << std::endl;
189 mgos.preStep(time, dt, oldValues);
191 std::cout <<
"== prepare assembler (" << subTimer.elapsed() <<
"s)"
196 std::cout <<
"== setup result vector" << std::endl;
200 mgos.interpolate(*oldValues[0],f,xnew);
202 std::cout <<
"== setup result vector (" << subTimer.elapsed() <<
"s)"
207 std::cout <<
"== apply solver" << std::endl;
210 pdeSolver.apply(xnew);
212 std::cout <<
"== apply solver (" << subTimer.elapsed() <<
"s)"
217 std::cout <<
"== cleanup assembler" << std::endl;
222 std::cout <<
"== cleanup assembler (" << subTimer.elapsed() <<
"s)"
228 std::cout <<
"== time step done (" << allTimer.elapsed() <<
"s)"
244 shared_ptr<const TrialV>
apply(T time, T dt)
250 ios_base_all_saver format_attribute_saver(std::cout);
253 std::cout <<
"TIME STEP [" << parameters->
name() <<
"] "
254 << std::setw(6) << step
256 << std::setw(12) << std::setprecision(4) << std::scientific
259 << std::setw(12) << std::setprecision(4) << std::scientific
262 << std::setw(12) << std::setprecision(4) << std::scientific
268 std::cout <<
"== prepare assembler" << std::endl;
271 mgos.preStep(step, time, dt);
273 std::cout <<
"== prepare assembler (" << subTimer.elapsed() <<
"s)"
278 std::cout <<
"== setup result vector" << std::endl;
282 xnew(
new TrialV(*mgos.getCache()->getUnknowns(step-1)));
284 std::cout <<
"== setup result vector (" << subTimer.elapsed() <<
"s)"
289 std::cout <<
"== apply solver" << std::endl;
292 pdeSolver.apply(*xnew);
294 std::cout <<
"== apply solver (" << subTimer.elapsed() <<
"s)"
299 std::cout <<
"== cleanup assembler" << std::endl;
304 std::cout <<
"== cleanup assembler (" << subTimer.elapsed() <<
"s)"
308 mgos.getCache()->setUnknowns(step, xnew);
313 std::cout <<
"== time step done (" << allTimer.elapsed() <<
"s)"
335 shared_ptr<const TrialV>
apply(T time, T dt, F&
f)
341 ios_base_all_saver format_attribute_saver(std::cout);
344 std::cout <<
"TIME STEP [" << parameters->
name() <<
"] "
345 << std::setw(6) << step
347 << std::setw(12) << std::setprecision(4) << std::scientific
350 << std::setw(12) << std::setprecision(4) << std::scientific
353 << std::setw(12) << std::setprecision(4) << std::scientific
359 std::cout <<
"== prepare assembler" << std::endl;
362 mgos.preStep(step, time, dt);
364 std::cout <<
"== prepare assembler (" << subTimer.elapsed() <<
"s)"
369 std::cout <<
"== setup result vector" << std::endl;
372 shared_ptr<TrialV> xnew(
new TrialV(mgos.trialGridFunctionSpace()));
375 mgos.interpolate(*mgos.getCache()->getUnknowns(step-1),
f,*xnew);
377 std::cout <<
"== setup result vector (" << subTimer.elapsed() <<
"s)"
382 std::cout <<
"== apply solver" << std::endl;
385 pdeSolver.apply(*xnew);
387 std::cout <<
"== apply solver (" << subTimer.elapsed() <<
"s)"
392 std::cout <<
"== cleanup assembler" << std::endl;
397 std::cout <<
"== cleanup assembler (" << subTimer.elapsed() <<
"s)"
401 mgos.getCache()->setUnknowns(step, xnew);
406 std::cout <<
"== time step done (" << allTimer.elapsed() <<
"s)"
417 #endif // DUNE_PDELAB_MULTISTEP_METHOD_HH
Do one step of a multi-step time-stepping scheme.
Definition: method.hh:32
shared_ptr< const TrialV > apply(T time, T dt, F &f)
do one step (with caching)
Definition: method.hh:335
T apply(T time, T dt, const OldValues &oldValues, F &f, TrialV &xnew)
do one step;
Definition: method.hh:162
T apply(T time, T dt, const OldValues &oldValues, TrialV &xnew)
do one step;
Definition: method.hh:85
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: method.hh:62
virtual std::string name() const =0
Return name of the scheme.
shared_ptr< const TrialV > apply(T time, T dt)
do one step (with caching)
Definition: method.hh:244
const F & f
Definition: common/constraints.hh:145
void setMethod(const Parameters ¶meters_)
redefine the method to be used; can be done before every step
Definition: method.hh:69
MultiStepMethod(const Parameters ¶meters_, MGOS &mgos_, PDESolver &pdeSolver_)
construct a new multi-step scheme
Definition: method.hh:52
Base parameter class for multi step time schemes.
Definition: parameter.hh:26