4 #ifndef DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH 5 #define DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH 10 #include <dune/common/ios_state.hh> 33 linear_solver_time(0.0),
34 linear_solver_iterations(0),
35 nonlinear_solver_iterations(0)
55 template<
class T,
class IGOS,
class PDESOLVER,
class TrlV,
class TstV = TrlV>
58 typedef typename PDESOLVER::Result PDESolverResult;
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;
121 setStepNumber(res.successful.timesteps+1);
144 T
apply (T time, T dt, TrlV& xold, TrlV& xnew)
147 ios_base_all_saver format_attribute_saver(std::cout);
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();
225 res.total.timesteps += 1;
228 PDESolverResult pderes = pdesolver.result();
239 for (
unsigned i=1; i<method->s(); ++i)
delete x[i];
249 res.total.timesteps += 1;
254 res.successful.timesteps += 1;
255 if (verbosityLevel>=1){
256 std::ios_base::fmtflags oldflags = std::cout.flags();
257 std::cout <<
"::: timesteps " << std::setw(6) << res.successful.timesteps
258 <<
" (" << res.total.timesteps <<
")" << std::endl;
259 std::cout <<
"::: nl iterations " << std::setw(6) << res.successful.nonlinear_solver_iterations
260 <<
" (" << res.total.nonlinear_solver_iterations <<
")" << std::endl;
261 std::cout <<
"::: lin iterations " << std::setw(6) << res.successful.linear_solver_iterations
262 <<
" (" << res.total.linear_solver_iterations <<
")" << std::endl;
263 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
264 << res.successful.assembler_time <<
" (" << res.total.assembler_time <<
")" << std::endl;
265 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
266 << res.successful.linear_solver_time <<
" (" << res.total.linear_solver_time <<
")" << std::endl;
267 std::cout.flags(oldflags);
285 T
apply (T time, T dt, TrlV& xold, F& f, TrlV& xnew)
291 ios_base_all_saver format_attribute_saver(std::cout);
293 std::vector<TrlV*> x(1);
296 if (verbosityLevel>=1){
297 std::ios_base::fmtflags oldflags = std::cout.flags();
298 std::cout <<
"TIME STEP [" << method->name() <<
"] " 299 << std::setw(6) << step
301 << std::setw(12) << std::setprecision(4) << std::scientific
304 << std::setw(12) << std::setprecision(4) << std::scientific
307 << std::setw(12) << std::setprecision(4) << std::scientific
310 std::cout.flags(oldflags);
314 igos.preStep(*method,time,dt);
317 for (
unsigned r=1; r<=method->s(); ++r)
319 if (verbosityLevel>=2){
320 std::ios_base::fmtflags oldflags = std::cout.flags();
321 std::cout <<
"STAGE " 324 << std::setw(12) << std::setprecision(4) << std::scientific
325 << time+method->d(r)*dt
327 std::cout.flags(oldflags);
342 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
346 igos.interpolate(r,*x[r-1],f,*x[r]);
350 pdesolver.apply(*x[r]);
355 PDESolverResult pderes = pdesolver.result();
364 res.total.timesteps += 1;
367 PDESolverResult pderes = pdesolver.result();
378 for (
unsigned i=1; i<method->s(); ++i)
delete x[i];
388 res.total.timesteps += 1;
393 res.successful.timesteps += 1;
394 if (verbosityLevel>=1){
395 std::ios_base::fmtflags oldflags = std::cout.flags();
396 std::cout <<
"::: timesteps " << std::setw(6) << res.successful.timesteps
397 <<
" (" << res.total.timesteps <<
")" << std::endl;
398 std::cout <<
"::: nl iterations " << std::setw(6) << res.successful.nonlinear_solver_iterations
399 <<
" (" << res.total.nonlinear_solver_iterations <<
")" << std::endl;
400 std::cout <<
"::: lin iterations " << std::setw(6) << res.successful.linear_solver_iterations
401 <<
" (" << res.total.linear_solver_iterations <<
")" << std::endl;
402 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
403 << res.successful.assembler_time <<
" (" << res.total.assembler_time <<
")" << std::endl;
404 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
405 << res.successful.linear_solver_time <<
" (" << res.total.linear_solver_time <<
")" << std::endl;
406 std::cout.flags(oldflags);
416 PDESOLVER& pdesolver;
425 #endif // DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH OneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, PDESOLVER &pdesolver_)
construct a new one step scheme
Definition: implicitonestep.hh:75
void setMethod(const TimeSteppingParameterInterface< T > &method_)
redefine the method to be used; can be done before every step
Definition: implicitonestep.hh:132
T apply(T time, T dt, TrlV &xold, TrlV &xnew)
do one step;
Definition: implicitonestep.hh:144
OneStepMethodPartialResult total
Definition: implicitonestep.hh:41
OneStepMethodPartialResult successful
Definition: implicitonestep.hh:42
OneStepMethodResult()
Definition: implicitonestep.hh:43
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:285
int linear_solver_iterations
Definition: implicitonestep.hh:27
double assembler_time
Definition: implicitonestep.hh:25
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
double linear_solver_time
Definition: implicitonestep.hh:26
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: implicitonestep.hh:84
int nonlinear_solver_iterations
Definition: implicitonestep.hh:28
void setResult(const OneStepMethodResult &result_)
Set a new result.
Definition: implicitonestep.hh:118
Definition: implicitonestep.hh:22
const PDESOLVER & getPDESolver() const
Access to the (non) linear solver.
Definition: implicitonestep.hh:96
OneStepMethodPartialResult()
Definition: implicitonestep.hh:30
PDESOLVER & getPDESolver()
Access to the (non) linear solver.
Definition: implicitonestep.hh:102
Do one step of a time-stepping scheme.
Definition: implicitonestep.hh:56
void setStepNumber(int newstep)
change number of current step
Definition: implicitonestep.hh:93
OneStepMethodResult Result
Definition: implicitonestep.hh:61
Definition: implicitonestep.hh:39
const Result & result() const
Definition: implicitonestep.hh:107
unsigned int timesteps
Definition: implicitonestep.hh:24