4 #ifndef DUNE_PDELAB_INSTATIONARY_EXPLICITONESTEP_HH
5 #define DUNE_PDELAB_INSTATIONARY_EXPLICITONESTEP_HH
10 #include <dune/common/ios_state.hh>
63 template<
class R,
class IGOS>
72 CFLTimeController (R cfl_, R target_,
const IGOS& igos_) : cfl(cfl_), target(target_), igos(igos_)
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 : method(&method_), igos(igos_), ls(ls_), verbosityLevel(1), step(1), D(igos),
132 DUNE_THROW(
Exception,
"explicit one step method called with implicit scheme");
133 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
150 : method(&method_), igos(igos_), ls(ls_), verbosityLevel(1), step(1), D(igos),
151 tc(&tc_), allocated(false)
154 DUNE_THROW(
Exception,
"explicit one step method called with implicit scheme");
159 if (allocated)
delete tc;
165 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
168 verbosityLevel = level;
186 DUNE_THROW(
Exception,
"explicit one step method called with implicit scheme");
197 T
apply (T time, T dt, TrlV& xold, TrlV& xnew)
199 DefaultLimiter limiter;
200 return apply(time,dt,xold,xnew,limiter);
203 template<
typename Limiter>
204 T
apply (T time, T dt, TrlV& xold, TrlV& xnew, Limiter& limiter)
207 ios_base_all_saver format_attribute_saver(std::cout);
209 mytag <<
"ExplicitOneStepMethod::apply(): ";
211 std::vector<TrlV*> x(1);
213 if(verbosityLevel>=4)
214 std::cout << mytag <<
"Creating residual vectors alpha and beta..."
216 TstV alpha(igos.testGridFunctionSpace()), beta(igos.testGridFunctionSpace());
217 if(verbosityLevel>=4)
219 <<
"Creating residual vectors alpha and beta... done."
222 if (verbosityLevel>=1){
223 std::ios_base::fmtflags oldflags = std::cout.flags();
224 std::cout <<
"TIME STEP [" << method->
name() <<
"] "
225 << std::setw(6) << step
227 << std::setw(12) << std::setprecision(4) << std::scientific
230 << std::setw(12) << std::setprecision(4) << std::scientific
233 << std::setw(12) << std::setprecision(4) << std::scientific
236 std::cout.flags(oldflags);
240 if(verbosityLevel>=4)
241 std::cout << mytag <<
"Preparing assembler..." << std::endl;
242 igos.preStep(*method,time,dt);
243 if(verbosityLevel>=4)
244 std::cout << mytag <<
"Preparing assembler... done." << std::endl;
247 for(
unsigned r=1; r<=method->
s(); ++r)
250 stagetag <<
"stage " << r <<
": ";
251 if (verbosityLevel>=4)
252 std::cout << stagetag <<
"Start." << std::endl;
254 if (verbosityLevel>=2){
255 std::ios_base::fmtflags oldflags = std::cout.flags();
256 std::cout <<
"STAGE "
259 << std::setw(12) << std::setprecision(4) << std::scientific
260 << time+method->
d(r)*dt
262 std::cout.flags(oldflags);
270 if (r>1) xnew = *(x[r-1]);
275 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
283 if (verbosityLevel>=4) std::cout <<
"assembling D, alpha, beta ..." << std::endl;
289 limiter.prestage(*x[r-1]);
291 if(verbosityLevel>=4)
292 std::cout << stagetag <<
"Assembling residual..." << std::endl;
293 igos.explicit_jacobian_residual(r,x,D,alpha,beta);
294 if(verbosityLevel>=4)
295 std::cout << stagetag <<
"Assembling residual... done."
302 newdt = std::min(newdt, dt);
304 if (verbosityLevel>=4){
305 std::ios_base::fmtflags oldflags = std::cout.flags();
306 std::cout << stagetag
308 << std::setw(12) << std::setprecision(4) << std::scientific
311 << std::setw(12) << std::setprecision(4) << std::scientific
314 std::cout.flags(oldflags);
317 if (verbosityLevel>=2 && newdt!=dt)
319 std::ios_base::fmtflags oldflags = std::cout.flags();
320 std::cout <<
"changed dt to "
321 << std::setw(12) << std::setprecision(4) << std::scientific
324 std::cout.flags(oldflags);
330 if (verbosityLevel>=4)
331 std::cout << stagetag
332 <<
"Combining residuals with selected dt..."
335 if (verbosityLevel>=4)
336 std::cout << stagetag
337 <<
"Combining residuals with selected dt... done."
341 if (verbosityLevel>=4)
342 std::cout << stagetag <<
"Solving diagonal system..."
344 ls.apply(D,*x[r],alpha,0.99);
345 if (verbosityLevel>=4)
346 std::cout << stagetag <<
"Solving diagonal system... done."
350 limiter.poststage(*x[r]);
353 if (verbosityLevel>=4)
354 std::cout << stagetag <<
"Cleanup..." << std::endl;
356 if (verbosityLevel>=4)
357 std::cout << stagetag <<
"Cleanup... done" << std::endl;
359 if (verbosityLevel>=4)
360 std::cout << stagetag <<
"Finished." << std::endl;
364 for(
unsigned i=1; i<method->
s(); ++i)
delete x[i];
367 if (verbosityLevel>=4)
368 std::cout << mytag <<
"Cleanup..." << std::endl;
370 if (verbosityLevel>=4)
371 std::cout << mytag <<
"Cleanup... done." << std::endl;
392 const TimeSteppingParameterInterface<T> *method;
398 TimeControllerInterface<T> *tc;
CFLTimeController(R cfl_, const IGOS &igos_)
Definition: explicitonestep.hh:69
Default time controller; just returns given dt.
Definition: explicitonestep.hh:44
ExplicitOneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, LS &ls_, TC &tc_)
construct a new one step scheme
Definition: explicitonestep.hh:149
virtual std::string name() const =0
Return name of the scheme.
virtual RealType suggestTimestep(RealType time, RealType givendt)
Return name of the scheme.
Definition: explicitonestep.hh:82
void setTarget(R target_)
Definition: explicitonestep.hh:75
virtual unsigned s() const =0
Return number of stages of the method.
R RealType
Definition: explicitonestep.hh:67
limit time step to maximum dt * CFL number
Definition: explicitonestep.hh:64
R RealType
Definition: explicitonestep.hh:29
R RealType
Definition: explicitonestep.hh:47
void setMethod(const TimeSteppingParameterInterface< T > &method_)
redefine the method to be used; can be done before every step
Definition: explicitonestep.hh:182
virtual ~TimeControllerInterface()
every abstract base class has a virtual destructor
Definition: explicitonestep.hh:36
Do one step of an explicit time-stepping scheme.
Definition: explicitonestep.hh:109
ExplicitOneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, LS &ls_)
construct a new one step scheme
Definition: explicitonestep.hh:127
Base class for all PDELab exceptions.
Definition: exceptions.hh:17
virtual RealType suggestTimestep(RealType time, RealType givendt)=0
Return name of the scheme.
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: explicitonestep.hh:163
Definition: adaptivity.hh:27
void setStepNumber(int newstep)
change number of current step
Definition: explicitonestep.hh:172
virtual bool implicit() const =0
Return true if method is implicit.
Controller interface for adaptive time stepping.
Definition: explicitonestep.hh:26
virtual RealType suggestTimestep(RealType time, RealType givendt)
Return name of the scheme.
Definition: explicitonestep.hh:51
CFLTimeController(R cfl_, R target_, const IGOS &igos_)
Definition: explicitonestep.hh:72
virtual R d(int r) const =0
Return entries of the d Vector.
T apply(T time, T dt, TrlV &xold, TrlV &xnew, Limiter &limiter)
Definition: explicitonestep.hh:204
~ExplicitOneStepMethod()
Definition: explicitonestep.hh:157
T apply(T time, T dt, TrlV &xold, TrlV &xnew)
do one step;
Definition: explicitonestep.hh:197
Insert standard boilerplate into log messages.
Definition: logtag.hh:180