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 PDESolverResult pderes = pdesolver.result();
239 for (
unsigned i=1; i<method->
s(); ++i)
delete x[i];
255 if (verbosityLevel>=1){
256 std::ios_base::fmtflags oldflags = std::cout.flags();
263 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
265 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
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();
367 PDESolverResult pderes = pdesolver.result();
378 for (
unsigned i=1; i<method->
s(); ++i)
delete x[i];
394 if (verbosityLevel>=1){
395 std::ios_base::fmtflags oldflags = std::cout.flags();
402 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
404 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
406 std::cout.flags(oldflags);
416 PDESOLVER& pdesolver;
425 #endif // DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH