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;
132 T
apply (T time, T dt, TrlV& xold, TrlV& xnew)
135 ios_base_all_saver format_attribute_saver(std::cout);
140 std::vector<TrlV*> x(1);
143 if (verbosityLevel>=1){
144 std::ios_base::fmtflags oldflags = std::cout.flags();
145 std::cout <<
"TIME STEP [" << method->name() <<
"] " 146 << std::setw(6) << step
148 << std::setw(12) << std::setprecision(4) << std::scientific
151 << std::setw(12) << std::setprecision(4) << std::scientific
154 << std::setw(12) << std::setprecision(4) << std::scientific
157 std::cout.flags(oldflags);
161 igos.preStep(*method,time,dt);
164 for (
unsigned r=1; r<=method->s(); ++r)
166 if (verbosityLevel>=2){
167 std::ios_base::fmtflags oldflags = std::cout.flags();
168 std::cout <<
"STAGE " 171 << std::setw(12) << std::setprecision(4) << std::scientific
172 << time+method->d(r)*dt
174 std::cout.flags(oldflags);
185 if (r>1) xnew = *(x[r-1]);
190 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
199 pdesolver.apply(*x[r]);
204 PDESolverResult pderes = pdesolver.result();
213 res.total.timesteps += 1;
216 PDESolverResult pderes = pdesolver.result();
227 for (
unsigned i=1; i<method->s(); ++i)
delete x[i];
237 res.total.timesteps += 1;
242 res.successful.timesteps += 1;
243 if (verbosityLevel>=1){
244 std::ios_base::fmtflags oldflags = std::cout.flags();
245 std::cout <<
"::: timesteps " << std::setw(6) << res.successful.timesteps
246 <<
" (" << res.total.timesteps <<
")" << std::endl;
247 std::cout <<
"::: nl iterations " << std::setw(6) << res.successful.nonlinear_solver_iterations
248 <<
" (" << res.total.nonlinear_solver_iterations <<
")" << std::endl;
249 std::cout <<
"::: lin iterations " << std::setw(6) << res.successful.linear_solver_iterations
250 <<
" (" << res.total.linear_solver_iterations <<
")" << std::endl;
251 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
252 << res.successful.assembler_time <<
" (" << res.total.assembler_time <<
")" << std::endl;
253 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
254 << res.successful.linear_solver_time <<
" (" << res.total.linear_solver_time <<
")" << std::endl;
255 std::cout.flags(oldflags);
273 T
apply (T time, T dt, TrlV& xold, F& f, TrlV& xnew)
279 ios_base_all_saver format_attribute_saver(std::cout);
281 std::vector<TrlV*> x(1);
284 if (verbosityLevel>=1){
285 std::ios_base::fmtflags oldflags = std::cout.flags();
286 std::cout <<
"TIME STEP [" << method->name() <<
"] " 287 << std::setw(6) << step
289 << std::setw(12) << std::setprecision(4) << std::scientific
292 << std::setw(12) << std::setprecision(4) << std::scientific
295 << std::setw(12) << std::setprecision(4) << std::scientific
298 std::cout.flags(oldflags);
302 igos.preStep(*method,time,dt);
305 for (
unsigned r=1; r<=method->s(); ++r)
307 if (verbosityLevel>=2){
308 std::ios_base::fmtflags oldflags = std::cout.flags();
309 std::cout <<
"STAGE " 312 << std::setw(12) << std::setprecision(4) << std::scientific
313 << time+method->d(r)*dt
315 std::cout.flags(oldflags);
330 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
334 igos.interpolate(r,*x[r-1],f,*x[r]);
338 pdesolver.apply(*x[r]);
343 PDESolverResult pderes = pdesolver.result();
352 res.total.timesteps += 1;
355 PDESolverResult pderes = pdesolver.result();
366 for (
unsigned i=1; i<method->s(); ++i)
delete x[i];
376 res.total.timesteps += 1;
381 res.successful.timesteps += 1;
382 if (verbosityLevel>=1){
383 std::ios_base::fmtflags oldflags = std::cout.flags();
384 std::cout <<
"::: timesteps " << std::setw(6) << res.successful.timesteps
385 <<
" (" << res.total.timesteps <<
")" << std::endl;
386 std::cout <<
"::: nl iterations " << std::setw(6) << res.successful.nonlinear_solver_iterations
387 <<
" (" << res.total.nonlinear_solver_iterations <<
")" << std::endl;
388 std::cout <<
"::: lin iterations " << std::setw(6) << res.successful.linear_solver_iterations
389 <<
" (" << res.total.linear_solver_iterations <<
")" << std::endl;
390 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
391 << res.successful.assembler_time <<
" (" << res.total.assembler_time <<
")" << std::endl;
392 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
393 << res.successful.linear_solver_time <<
" (" << res.total.linear_solver_time <<
")" << std::endl;
394 std::cout.flags(oldflags);
404 PDESOLVER& pdesolver;
void setMethod(const TimeSteppingParameterInterface< T > &method_)
redefine the method to be used; can be done before every step
Definition: implicitonestep.hh:120
OneStepMethodResult Result
Definition: implicitonestep.hh:61
int linear_solver_iterations
Definition: implicitonestep.hh:27
OneStepMethodResult()
Definition: implicitonestep.hh:43
OneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, PDESOLVER &pdesolver_)
construct a new one step scheme
Definition: implicitonestep.hh:75
Do one step of a time-stepping scheme.
Definition: implicitonestep.hh:56
Definition: adaptivity.hh:27
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: implicitonestep.hh:84
T apply(T time, T dt, TrlV &xold, TrlV &xnew)
do one step;
Definition: implicitonestep.hh:132
double assembler_time
Definition: implicitonestep.hh:25
const Result & result() const
Definition: implicitonestep.hh:107
PDESOLVER & getPDESolver()
Access to the (non) linear solver.
Definition: implicitonestep.hh:102
double linear_solver_time
Definition: implicitonestep.hh:26
OneStepMethodPartialResult successful
Definition: implicitonestep.hh:42
OneStepMethodPartialResult()
Definition: implicitonestep.hh:30
Definition: implicitonestep.hh:39
OneStepMethodPartialResult total
Definition: implicitonestep.hh:41
void setStepNumber(int newstep)
change number of current step
Definition: implicitonestep.hh:93
unsigned int timesteps
Definition: implicitonestep.hh:24
const PDESOLVER & getPDESolver() const
Access to the (non) linear solver.
Definition: implicitonestep.hh:96
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:273
int nonlinear_solver_iterations
Definition: implicitonestep.hh:28
Definition: implicitonestep.hh:22