4#ifndef DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH
5#define DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH
11#include <dune/pdelab/instationary/onestepparameter.hh>
22 struct OneStepMethodPartialResult
24 unsigned int timesteps;
25 double assembler_time;
26 double linear_solver_time;
27 int linear_solver_iterations;
28 int nonlinear_solver_iterations;
30 OneStepMethodPartialResult() :
33 linear_solver_time(0.0),
34 linear_solver_iterations(0),
35 nonlinear_solver_iterations(0)
39 struct OneStepMethodResult
41 OneStepMethodPartialResult total;
42 OneStepMethodPartialResult successful;
43 OneStepMethodResult() : total(), successful()
55 template<
class T,
class IGOS,
class PDESOLVER,
class TrlV,
class TstV = TrlV>
58 typedef typename PDESOLVER::Result PDESolverResult;
61 typedef OneStepMethodResult Result;
76 IGOS& igos_, PDESOLVER& pdesolver_)
77 : method(&method_), igos(igos_), pdesolver(pdesolver_), verbosityLevel(1), step(1), res()
79 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
86 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
89 verbosityLevel = level;
107 const Result& result()
const
144 T
apply (T time, T dt, TrlV& xold, TrlV& xnew)
150 OneStepMethodPartialResult step_result;
152 std::vector<TrlV*> x(1);
155 if (verbosityLevel>=1){
156 std::ios_base::fmtflags oldflags = std::cout.flags();
157 std::cout <<
"TIME STEP [" << method->
name() <<
"] "
158 << std::setw(6) << step
160 << std::setw(12) << std::setprecision(4) << std::scientific
163 << std::setw(12) << std::setprecision(4) << std::scientific
166 << std::setw(12) << std::setprecision(4) << std::scientific
169 std::cout.flags(oldflags);
173 igos.preStep(*method,time,dt);
176 for (
unsigned r=1; r<=method->
s(); ++r)
178 if (verbosityLevel>=2){
179 std::ios_base::fmtflags oldflags = std::cout.flags();
180 std::cout <<
"STAGE "
183 << std::setw(12) << std::setprecision(4) << std::scientific
184 << time+method->
d(r)*dt
186 std::cout.flags(oldflags);
197 if (r>1) xnew = *(x[r-1]);
202 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
211 pdesolver.apply(*x[r]);
216 PDESolverResult pderes = pdesolver.result();
217 step_result.assembler_time += pderes.assembler_time;
218 step_result.linear_solver_time += pderes.linear_solver_time;
219 step_result.linear_solver_iterations += pderes.linear_solver_iterations;
220 step_result.nonlinear_solver_iterations += pderes.iterations;
221 res.total.assembler_time += step_result.assembler_time;
222 res.total.linear_solver_time += step_result.linear_solver_time;
223 res.total.linear_solver_iterations += step_result.linear_solver_iterations;
224 res.total.nonlinear_solver_iterations += step_result.nonlinear_solver_iterations;
225 res.total.timesteps += 1;
228 for (
unsigned i=1; i<r; ++i)
delete x[i];
234 PDESolverResult pderes = pdesolver.result();
235 step_result.assembler_time += pderes.assembler_time;
236 step_result.linear_solver_time += pderes.linear_solver_time;
237 step_result.linear_solver_iterations += pderes.linear_solver_iterations;
238 step_result.nonlinear_solver_iterations += pderes.iterations;
245 for (
unsigned i=1; i<method->
s(); ++i)
delete x[i];
251 res.total.assembler_time += step_result.assembler_time;
252 res.total.linear_solver_time += step_result.linear_solver_time;
253 res.total.linear_solver_iterations += step_result.linear_solver_iterations;
254 res.total.nonlinear_solver_iterations += step_result.nonlinear_solver_iterations;
255 res.total.timesteps += 1;
256 res.successful.assembler_time += step_result.assembler_time;
257 res.successful.linear_solver_time += step_result.linear_solver_time;
258 res.successful.linear_solver_iterations += step_result.linear_solver_iterations;
259 res.successful.nonlinear_solver_iterations += step_result.nonlinear_solver_iterations;
260 res.successful.timesteps += 1;
261 if (verbosityLevel>=1){
262 std::ios_base::fmtflags oldflags = std::cout.flags();
263 std::cout <<
"::: timesteps " << std::setw(6) << res.successful.timesteps
264 <<
" (" << res.total.timesteps <<
")" << std::endl;
265 std::cout <<
"::: nl iterations " << std::setw(6) << res.successful.nonlinear_solver_iterations
266 <<
" (" << res.total.nonlinear_solver_iterations <<
")" << std::endl;
267 std::cout <<
"::: lin iterations " << std::setw(6) << res.successful.linear_solver_iterations
268 <<
" (" << res.total.linear_solver_iterations <<
")" << std::endl;
269 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
270 << res.successful.assembler_time <<
" (" << res.total.assembler_time <<
")" << std::endl;
271 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
272 << res.successful.linear_solver_time <<
" (" << res.total.linear_solver_time <<
")" << std::endl;
273 std::cout.flags(oldflags);
291 T
apply (T time, T dt, TrlV& xold, F& f, TrlV& xnew)
294 OneStepMethodPartialResult step_result;
299 std::vector<TrlV*> x(1);
302 if (verbosityLevel>=1){
303 std::ios_base::fmtflags oldflags = std::cout.flags();
304 std::cout <<
"TIME STEP [" << method->
name() <<
"] "
305 << std::setw(6) << step
307 << std::setw(12) << std::setprecision(4) << std::scientific
310 << std::setw(12) << std::setprecision(4) << std::scientific
313 << std::setw(12) << std::setprecision(4) << std::scientific
316 std::cout.flags(oldflags);
320 igos.preStep(*method,time,dt);
323 for (
unsigned r=1; r<=method->
s(); ++r)
325 if (verbosityLevel>=2){
326 std::ios_base::fmtflags oldflags = std::cout.flags();
327 std::cout <<
"STAGE "
330 << std::setw(12) << std::setprecision(4) << std::scientific
331 << time+method->
d(r)*dt
333 std::cout.flags(oldflags);
348 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
352 TrlV
const *
const init_guess = (r==1) ? &xnew : x[r-1];
354 igos.interpolate(r,*init_guess,f,*x[r]);
358 pdesolver.apply(*x[r]);
363 PDESolverResult pderes = pdesolver.result();
364 step_result.assembler_time += pderes.assembler_time;
365 step_result.linear_solver_time += pderes.linear_solver_time;
366 step_result.linear_solver_iterations += pderes.linear_solver_iterations;
367 step_result.nonlinear_solver_iterations += pderes.iterations;
368 res.total.assembler_time += step_result.assembler_time;
369 res.total.linear_solver_time += step_result.linear_solver_time;
370 res.total.linear_solver_iterations += step_result.linear_solver_iterations;
371 res.total.nonlinear_solver_iterations += step_result.nonlinear_solver_iterations;
372 res.total.timesteps += 1;
375 for (
unsigned i=1; i<r; ++i)
delete x[i];
381 PDESolverResult pderes = pdesolver.result();
382 step_result.assembler_time += pderes.assembler_time;
383 step_result.linear_solver_time += pderes.linear_solver_time;
384 step_result.linear_solver_iterations += pderes.linear_solver_iterations;
385 step_result.nonlinear_solver_iterations += pderes.iterations;
392 for (
unsigned i=1; i<method->
s(); ++i)
delete x[i];
398 res.total.assembler_time += step_result.assembler_time;
399 res.total.linear_solver_time += step_result.linear_solver_time;
400 res.total.linear_solver_iterations += step_result.linear_solver_iterations;
401 res.total.nonlinear_solver_iterations += step_result.nonlinear_solver_iterations;
402 res.total.timesteps += 1;
403 res.successful.assembler_time += step_result.assembler_time;
404 res.successful.linear_solver_time += step_result.linear_solver_time;
405 res.successful.linear_solver_iterations += step_result.linear_solver_iterations;
406 res.successful.nonlinear_solver_iterations += step_result.nonlinear_solver_iterations;
407 res.successful.timesteps += 1;
408 if (verbosityLevel>=1){
409 std::ios_base::fmtflags oldflags = std::cout.flags();
410 std::cout <<
"::: timesteps " << std::setw(6) << res.successful.timesteps
411 <<
" (" << res.total.timesteps <<
")" << std::endl;
412 std::cout <<
"::: nl iterations " << std::setw(6) << res.successful.nonlinear_solver_iterations
413 <<
" (" << res.total.nonlinear_solver_iterations <<
")" << std::endl;
414 std::cout <<
"::: lin iterations " << std::setw(6) << res.successful.linear_solver_iterations
415 <<
" (" << res.total.linear_solver_iterations <<
")" << std::endl;
416 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
417 << res.successful.assembler_time <<
" (" << res.total.assembler_time <<
")" << std::endl;
418 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
419 << res.successful.linear_solver_time <<
" (" << res.total.linear_solver_time <<
")" << std::endl;
420 std::cout.flags(oldflags);
430 PDESOLVER& pdesolver;
Do one step of a time-stepping scheme.
Definition: implicitonestep.hh:57
T apply(T time, T dt, TrlV &xold, F &f, TrlV &xnew)
do one step; This is a version which interpolates constraints at the start of each stage
Definition: implicitonestep.hh:291
void setStepNumber(int newstep)
change number of current step
Definition: implicitonestep.hh:93
void setMethod(const TimeSteppingParameterInterface< T > &method_)
redefine the method to be used; can be done before every step
Definition: implicitonestep.hh:132
PDESOLVER & getPDESolver()
Access to the (non) linear solver.
Definition: implicitonestep.hh:102
void setResult(const OneStepMethodResult &result_)
Set a new result.
Definition: implicitonestep.hh:118
T apply(T time, T dt, TrlV &xold, TrlV &xnew)
do one step;
Definition: implicitonestep.hh:144
OneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, PDESOLVER &pdesolver_)
construct a new one step scheme
Definition: implicitonestep.hh:75
const PDESOLVER & getPDESolver() const
Access to the (non) linear solver.
Definition: implicitonestep.hh:96
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: implicitonestep.hh:84
Utility class for storing and resetting stream attributes.
Definition: ios_state.hh:34
virtual std::string name() const =0
Return name of the scheme.
virtual R d(int r) const =0
Return entries of the d Vector.
virtual unsigned s() const =0
Return number of stages of the method.
Utility class for storing and resetting stream attributes.
Dune namespace.
Definition: alignedallocator.hh:13