4 #ifndef DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH
5 #define DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH
10 #include <dune/common/ios_state.hh>
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;
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();
228 for (
unsigned i=1; i<r; ++i)
delete x[i];
234 PDESolverResult pderes = pdesolver.result();
245 for (
unsigned i=1; i<method->
s(); ++i)
delete x[i];
261 if (verbosityLevel>=1){
262 std::ios_base::fmtflags oldflags = std::cout.flags();
269 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
271 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
273 std::cout.flags(oldflags);
291 T
apply (T time, T dt, TrlV& xold, F& f, TrlV& xnew)
297 ios_base_all_saver format_attribute_saver(std::cout);
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 igos.interpolate(r,*x[r-1],f,*x[r]);
356 pdesolver.apply(*x[r]);
361 PDESolverResult pderes = pdesolver.result();
373 for (
unsigned i=1; i<r; ++i)
delete x[i];
379 PDESolverResult pderes = pdesolver.result();
390 for (
unsigned i=1; i<method->
s(); ++i)
delete x[i];
406 if (verbosityLevel>=1){
407 std::ios_base::fmtflags oldflags = std::cout.flags();
414 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
416 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
418 std::cout.flags(oldflags);
428 PDESOLVER& pdesolver;
const std::string s
Definition: function.hh:843
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.
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Definition: implicitonestep.hh:23
int linear_solver_iterations
Definition: implicitonestep.hh:27
double linear_solver_time
Definition: implicitonestep.hh:26
int nonlinear_solver_iterations
Definition: implicitonestep.hh:28
unsigned int timesteps
Definition: implicitonestep.hh:24
OneStepMethodPartialResult()
Definition: implicitonestep.hh:30
double assembler_time
Definition: implicitonestep.hh:25
Definition: implicitonestep.hh:40
OneStepMethodResult()
Definition: implicitonestep.hh:43
OneStepMethodPartialResult total
Definition: implicitonestep.hh:41
OneStepMethodPartialResult successful
Definition: implicitonestep.hh:42
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
const PDESOLVER & getPDESolver() const
Access to the (non) linear solver.
Definition: implicitonestep.hh:96
void setStepNumber(int newstep)
change number of current step
Definition: implicitonestep.hh:93
OneStepMethodResult Result
Definition: implicitonestep.hh:61
void setMethod(const TimeSteppingParameterInterface< T > &method_)
redefine the method to be used; can be done before every step
Definition: implicitonestep.hh:132
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 Result & result() const
Definition: implicitonestep.hh:107
PDESOLVER & getPDESolver()
Access to the (non) linear solver.
Definition: implicitonestep.hh:102
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: implicitonestep.hh:84