4#ifndef DUNE_PDELAB_INSTATIONARY_EXPLICITONESTEP_HH
5#define DUNE_PDELAB_INSTATIONARY_EXPLICITONESTEP_HH
11#include <dune/pdelab/common/logtag.hh>
12#include <dune/pdelab/instationary/onestepparameter.hh>
63 template<
class R,
class IGOS>
69 CFLTimeController (R cfl_,
const IGOS& igos_) : cfl(cfl_), target(1e100), igos(igos_)
72 CFLTimeController (R cfl_, R target_,
const IGOS& igos_) : cfl(cfl_), target(target_), igos(igos_)
75 void setTarget (R target_)
84 RealType suggested = cfl*igos.suggestTimestep(givendt);
85 if (time+2.0*suggested<target)
87 if (time+suggested<target)
88 return 0.5*(target-time);
108 template<
class T,
class IGOS,
class LS,
class TrlV,
class TstV = TrlV,
class TC = SimpleTimeController<T> >
111 typedef typename TrlV::ElementType Real;
112 typedef typename IGOS::template MatrixContainer<Real>::Type M;
128 "Constructor with default assigned reduction is no longer supported. "
129 "Use the new constructor by passing a suitable reduction to the linear "
130 "system. In particular, use a dummy reduction (i.e. 0.99) on systems "
131 "known to have diagonal mass matrices (e.g. FV or DG) and where the "
132 "linear solver performs a block inversion on the first iteration "
133 "(essentially any). Otherwise, specify a proper reduction suitable "
137 : method(&method_), igos(igos_), ls(ls_), verbosityLevel(1), step(1), D(igos),
142 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
170 , ls_reduct{ ls_reduct_ }
174 "explicit one step method called with implicit scheme");
175 if (igos.trialGridFunctionSpace().gridView().comm().rank() > 0)
192 "Constructor with default assigned reduction is no longer supported. "
193 "Use the new constructor by passing a suitable reduction to the linear "
194 "system. In particular, use a dummy reduction (i.e. 0.99) on systems \n"
195 "known to have diagonal mass matrices (e.g. FV or DG) and where the \n"
196 "linear solver performs a block inversion on the first iteration \n"
197 "(essentially any). Otherwise, specify a proper reduction suitable \n"
201 : method(&method_), igos(igos_), ls(ls_), verbosityLevel(1), step(1), D(igos),
202 tc(&tc_), allocated(false), ls_reduct{0.99}
233 , ls_reduct{ ls_reduct_ }
237 "explicit one step method called with implicit scheme");
242 if (allocated)
delete tc;
248 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
251 verbosityLevel = level;
268 void setReduction(
const double& ls_reduction_) { ls_reduct = ls_reduction_; }
292 T
apply (T time, T dt, TrlV& xold, TrlV& xnew)
294 DefaultLimiter limiter;
295 return apply(time,dt,xold,xnew,limiter);
306 template<
typename Limiter>
307 T
apply (T time, T dt, TrlV& xold, TrlV& xnew, Limiter& limiter)
312 mytag <<
"ExplicitOneStepMethod::apply(): ";
314 std::vector<TrlV*> x(1);
316 if(verbosityLevel>=4)
317 std::cout << mytag <<
"Creating residual vectors alpha and beta..."
319 TstV alpha(igos.testGridFunctionSpace()), beta(igos.testGridFunctionSpace());
320 if(verbosityLevel>=4)
322 <<
"Creating residual vectors alpha and beta... done."
325 if (verbosityLevel>=1){
326 std::ios_base::fmtflags oldflags = std::cout.flags();
327 std::cout <<
"TIME STEP [" << method->
name() <<
"] "
328 << std::setw(6) << step
330 << std::setw(12) << std::setprecision(4) << std::scientific
333 << std::setw(12) << std::setprecision(4) << std::scientific
336 << std::setw(12) << std::setprecision(4) << std::scientific
339 std::cout.flags(oldflags);
343 if(verbosityLevel>=4)
344 std::cout << mytag <<
"Preparing assembler..." << std::endl;
345 igos.preStep(*method,time,dt);
346 if(verbosityLevel>=4)
347 std::cout << mytag <<
"Preparing assembler... done." << std::endl;
350 for(
unsigned r=1; r<=method->
s(); ++r)
353 stagetag <<
"stage " << r <<
": ";
354 if (verbosityLevel>=4)
355 std::cout << stagetag <<
"Start." << std::endl;
357 if (verbosityLevel>=2){
358 std::ios_base::fmtflags oldflags = std::cout.flags();
359 std::cout <<
"STAGE "
362 << std::setw(12) << std::setprecision(4) << std::scientific
363 << time+method->
d(r)*dt
365 std::cout.flags(oldflags);
373 if (r>1) xnew = *(x[r-1]);
378 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
386 if (verbosityLevel>=4) std::cout <<
"assembling D, alpha, beta ..." << std::endl;
392 limiter.prestage(*x[r-1]);
394 if(verbosityLevel>=4)
395 std::cout << stagetag <<
"Assembling residual..." << std::endl;
396 igos.explicit_jacobian_residual(r,x,D,alpha,beta);
397 if(verbosityLevel>=4)
398 std::cout << stagetag <<
"Assembling residual... done."
407 if (verbosityLevel>=4){
408 std::ios_base::fmtflags oldflags = std::cout.flags();
409 std::cout << stagetag
411 << std::setw(12) << std::setprecision(4) << std::scientific
414 << std::setw(12) << std::setprecision(4) << std::scientific
417 std::cout.flags(oldflags);
420 if (verbosityLevel>=2 && newdt!=dt)
422 std::ios_base::fmtflags oldflags = std::cout.flags();
423 std::cout <<
"changed dt to "
424 << std::setw(12) << std::setprecision(4) << std::scientific
427 std::cout.flags(oldflags);
433 if (verbosityLevel>=4)
434 std::cout << stagetag
435 <<
"Combining residuals with selected dt..."
438 if (verbosityLevel>=4)
439 std::cout << stagetag
440 <<
"Combining residuals with selected dt... done."
444 if (verbosityLevel>=4)
445 std::cout << stagetag <<
"Solving linear system..."
447 ls.apply(D,*x[r],alpha,ls_reduct);
448 if (verbosityLevel>=4)
449 std::cout << stagetag <<
"Solving linear system... done."
453 limiter.poststage(*x[r]);
456 if (verbosityLevel>=4)
457 std::cout << stagetag <<
"Cleanup..." << std::endl;
459 if (verbosityLevel>=4)
460 std::cout << stagetag <<
"Cleanup... done" << std::endl;
462 if (verbosityLevel>=4)
463 std::cout << stagetag <<
"Finished." << std::endl;
467 for(
unsigned i=1; i<method->
s(); ++i)
delete x[i];
470 if (verbosityLevel>=4)
471 std::cout << mytag <<
"Cleanup..." << std::endl;
473 if (verbosityLevel>=4)
474 std::cout << mytag <<
"Cleanup... done." << std::endl;
488 template <
typename F>
489 T
apply (T time, T dt, TrlV& xold, F& f, TrlV& xnew)
491 DefaultLimiter limiter;
492 return apply(time,dt,xold,f,xnew,limiter);
504 template<
typename F,
typename Limiter>
505 T
apply (T time, T dt, TrlV& xold, F&f, TrlV& xnew, Limiter& limiter)
510 mytag <<
"ExplicitOneStepMethod::apply(): ";
512 std::vector<TrlV*> x(1);
514 if(verbosityLevel>=4)
515 std::cout << mytag <<
"Creating residual vectors alpha and beta..."
517 TstV alpha(igos.testGridFunctionSpace()), beta(igos.testGridFunctionSpace());
518 if(verbosityLevel>=4)
520 <<
"Creating residual vectors alpha and beta... done."
527 if(verbosityLevel>=4)
528 std::cout << mytag <<
"Creating residual vector and update for residual formulation of linear problem per stage"
530 TrlV residual(igos.testGridFunctionSpace());
531 TrlV update(igos.testGridFunctionSpace());
532 if(verbosityLevel>=4)
533 std::cout << mytag <<
"Creating residual vector and update for residual... done."
537 if (verbosityLevel>=1){
538 std::ios_base::fmtflags oldflags = std::cout.flags();
539 std::cout <<
"TIME STEP [" << method->
name() <<
"] "
540 << std::setw(6) << step
542 << std::setw(12) << std::setprecision(4) << std::scientific
545 << std::setw(12) << std::setprecision(4) << std::scientific
548 << std::setw(12) << std::setprecision(4) << std::scientific
551 std::cout.flags(oldflags);
555 if(verbosityLevel>=4)
556 std::cout << mytag <<
"Preparing assembler..." << std::endl;
557 igos.preStep(*method,time,dt);
558 if(verbosityLevel>=4)
559 std::cout << mytag <<
"Preparing assembler... done." << std::endl;
562 for(
unsigned r=1; r<=method->
s(); ++r)
565 stagetag <<
"stage " << r <<
": ";
566 if (verbosityLevel>=4)
567 std::cout << stagetag <<
"Start." << std::endl;
569 if (verbosityLevel>=2){
570 std::ios_base::fmtflags oldflags = std::cout.flags();
571 std::cout <<
"STAGE "
574 << std::setw(12) << std::setprecision(4) << std::scientific
575 << time+method->
d(r)*dt
577 std::cout.flags(oldflags);
589 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
593 if (verbosityLevel>=4) std::cout <<
"assembling D, alpha, beta ..." << std::endl;
599 limiter.prestage(*x[r-1]);
602 TrlV
const *
const init_guess = (r==1) ? &xnew : x[r-1];
604 igos.interpolate(r,*init_guess,f,*x[r]);
606 if(verbosityLevel>=4)
607 std::cout << stagetag <<
"Assembling residual..." << std::endl;
608 igos.explicit_jacobian_residual(r,x,D,alpha,beta);
609 if(verbosityLevel>=4)
610 std::cout << stagetag <<
"Assembling residual... done."
619 if (verbosityLevel>=4){
620 std::ios_base::fmtflags oldflags = std::cout.flags();
621 std::cout << stagetag
623 << std::setw(12) << std::setprecision(4) << std::scientific
626 << std::setw(12) << std::setprecision(4) << std::scientific
629 std::cout.flags(oldflags);
632 if (verbosityLevel>=2 && newdt!=dt)
634 std::ios_base::fmtflags oldflags = std::cout.flags();
635 std::cout <<
"changed dt to "
636 << std::setw(12) << std::setprecision(4) << std::scientific
639 std::cout.flags(oldflags);
645 if (verbosityLevel>=4)
646 std::cout << stagetag
647 <<
"Combining residuals with selected dt..."
650 if (verbosityLevel>=4)
651 std::cout << stagetag
652 <<
"Combining residuals with selected dt... done."
659 using Backend::native;
660 native(D).mv(native(*x[r]), native(residual));
662 auto cc = igos.trialConstraints();
664 if (verbosityLevel>=4)
665 std::cout << stagetag <<
"Solving linear system..."
667 ls.apply(D, update, residual, ls_reduct);
668 if (verbosityLevel>=4)
669 std::cout << stagetag <<
"Solving linear system... done."
674 limiter.poststage(*x[r]);
677 if (verbosityLevel>=4)
678 std::cout << stagetag <<
"Cleanup..." << std::endl;
680 if (verbosityLevel>=4)
681 std::cout << stagetag <<
"Cleanup... done" << std::endl;
683 if (verbosityLevel>=4)
684 std::cout << stagetag <<
"Finished." << std::endl;
688 for(
unsigned i=1; i<method->
s(); ++i)
delete x[i];
691 if (verbosityLevel>=4)
692 std::cout << mytag <<
"Cleanup..." << std::endl;
694 if (verbosityLevel>=4)
695 std::cout << mytag <<
"Cleanup... done." << std::endl;
716 const TimeSteppingParameterInterface<T> *method;
722 TimeControllerInterface<T> *tc;
limit time step to maximum dt * CFL number
Definition: explicitonestep.hh:65
virtual RealType suggestTimestep(RealType time, RealType givendt) override
Return name of the scheme.
Definition: explicitonestep.hh:82
Base class for all PDELab exceptions.
Definition: exceptions.hh:19
Do one step of an explicit time-stepping scheme.
Definition: explicitonestep.hh:110
T apply(T time, T dt, TrlV &xold, F &f, TrlV &xnew, Limiter &limiter)
do one step;
Definition: explicitonestep.hh:505
ExplicitOneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, LS &ls_, TC &tc_)
construct a new one step scheme
Definition: explicitonestep.hh:200
T apply(T time, T dt, TrlV &xold, F &f, TrlV &xnew)
do one step;
Definition: explicitonestep.hh:489
T apply(T time, T dt, TrlV &xold, TrlV &xnew)
do one step;
Definition: explicitonestep.hh:292
ExplicitOneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, LS &ls_, double ls_reduct_)
Construct a new Explicit One Step Method object.
Definition: explicitonestep.hh:158
T apply(T time, T dt, TrlV &xold, TrlV &xnew, Limiter &limiter)
do one step;
Definition: explicitonestep.hh:307
void setMethod(const TimeSteppingParameterInterface< T > &method_)
redefine the method to be used; can be done before every step
Definition: explicitonestep.hh:278
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: explicitonestep.hh:246
void setStepNumber(int newstep)
change number of current step
Definition: explicitonestep.hh:255
ExplicitOneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, LS &ls_, TC &tc_, double ls_reduct_)
Construct a new Explicit One Step Method object.
Definition: explicitonestep.hh:220
ExplicitOneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, LS &ls_)
construct a new one step scheme
Definition: explicitonestep.hh:136
void setReduction(const double &ls_reduction_)
Set the Reduction for linear system soltion.
Definition: explicitonestep.hh:268
Insert standard boilerplate into log messages.
Definition: logtag.hh:172
Default time controller; just returns given dt.
Definition: explicitonestep.hh:45
virtual RealType suggestTimestep(RealType time, RealType givendt) override
Return name of the scheme.
Definition: explicitonestep.hh:51
Controller interface for adaptive time stepping.
Definition: explicitonestep.hh:27
virtual ~TimeControllerInterface()
every abstract base class has a virtual destructor
Definition: explicitonestep.hh:36
virtual RealType suggestTimestep(RealType time, RealType givendt)=0
Return name of the scheme.
Utility class for storing and resetting stream attributes.
Definition: ios_state.hh:34
#define DUNE_THROW(E, m)
Definition: exceptions.hh:218
void set_constrained_dofs(const CG &cg, typename XG::ElementType x, XG &xg)
construct constraints from given boundary condition function
Definition: constraints.hh:796
constexpr auto min
Function object that returns the smaller of the given values.
Definition: hybridutilities.hh:506
virtual bool implicit() const =0
Return true if method is implicit.
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