dune-pdelab 2.7-git
Loading...
Searching...
No Matches
implicitonestep.hh
Go to the documentation of this file.
1// -*- tab-width: 2; indent-tabs-mode: nil -*-
2// vi: set et ts=2 sw=2 sts=2:
3
4#ifndef DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH
5#define DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH
6
7#include <iostream>
8#include <iomanip>
9
10#include <dune/common/ios_state.hh>
12
13namespace Dune {
14 namespace PDELab {
15
21 // Status information of Newton's method
38
46
48
55 template<class T, class IGOS, class PDESOLVER, class TrlV, class TstV = TrlV>
57 {
58 typedef typename PDESOLVER::Result PDESolverResult;
59
60 public:
62
64
76 IGOS& igos_, PDESOLVER& pdesolver_)
77 : method(&method_), igos(igos_), pdesolver(pdesolver_), verbosityLevel(1), step(1), res()
78 {
79 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
80 verbosityLevel = 0;
81 }
82
84 void setVerbosityLevel (int level)
85 {
86 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
87 verbosityLevel = 0;
88 else
89 verbosityLevel = level;
90 }
91
93 void setStepNumber(int newstep) { step = newstep; }
94
96 const PDESOLVER & getPDESolver() const
97 {
98 return pdesolver;
99 }
100
102 PDESOLVER & getPDESolver()
103 {
104 return pdesolver;
105 }
106
107 const Result& result() const
108 {
109 return res;
110 }
111
113
118 void setResult (const OneStepMethodResult& result_)
119 {
120 res = result_;
122 }
123
125
133 {
134 method = &method_;
135 }
136
144 T apply (T time, T dt, TrlV& xold, TrlV& xnew)
145 {
146 // save formatting attributes
147 ios_base_all_saver format_attribute_saver(std::cout);
148
149 // do statistics
150 OneStepMethodPartialResult step_result;
151
152 std::vector<TrlV*> x(1); // vector of pointers to all steps
153 x[0] = &xold; // initially we have only one
154
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
159 << " time (from): "
160 << std::setw(12) << std::setprecision(4) << std::scientific
161 << time
162 << " dt: "
163 << std::setw(12) << std::setprecision(4) << std::scientific
164 << dt
165 << " time (to): "
166 << std::setw(12) << std::setprecision(4) << std::scientific
167 << time+dt
168 << std::endl;
169 std::cout.flags(oldflags);
170 }
171
172 // prepare assembler
173 igos.preStep(*method,time,dt);
174
175 // loop over all stages
176 for (unsigned r=1; r<=method->s(); ++r)
177 {
178 if (verbosityLevel>=2){
179 std::ios_base::fmtflags oldflags = std::cout.flags();
180 std::cout << "STAGE "
181 << r
182 << " time (to): "
183 << std::setw(12) << std::setprecision(4) << std::scientific
184 << time+method->d(r)*dt
185 << "." << std::endl;
186 std::cout.flags(oldflags);
187 }
188
189 // prepare stage
190 igos.preStage(r,x);
191
192 // get vector for current stage
193 if (r==method->s())
194 {
195 // last stage
196 x.push_back(&xnew);
197 if (r>1) xnew = *(x[r-1]); // if r=1 then xnew has already initial guess
198 }
199 else
200 {
201 // intermediate step
202 x.push_back(new TrlV(igos.trialGridFunctionSpace()));
203 if (r>1)
204 *(x[r]) = *(x[r-1]); // use result of last stage as initial guess
205 else
206 *(x[r]) = xnew;
207 }
208
209 // solve stage
210 try {
211 pdesolver.apply(*x[r]);
212 }
213 catch (...)
214 {
215 // time step failed -> accumulate to total only
216 PDESolverResult pderes = pdesolver.result();
217 step_result.assembler_time += pderes.assembler_time;
218 step_result.linear_solver_time += pderes.linear_solver_time;
219 step_result.linear_solver_iterations += pderes.linear_solver_iterations;
220 step_result.nonlinear_solver_iterations += pderes.iterations;
221 res.total.assembler_time += step_result.assembler_time;
225 res.total.timesteps += 1;
226 throw;
227 }
228 PDESolverResult pderes = pdesolver.result();
229 step_result.assembler_time += pderes.assembler_time;
230 step_result.linear_solver_time += pderes.linear_solver_time;
231 step_result.linear_solver_iterations += pderes.linear_solver_iterations;
232 step_result.nonlinear_solver_iterations += pderes.iterations;
233
234 // stage cleanup
235 igos.postStage();
236 }
237
238 // delete intermediate steps
239 for (unsigned i=1; i<method->s(); ++i) delete x[i];
240
241 // step cleanup
242 igos.postStep();
243
244 // update statistics
245 res.total.assembler_time += step_result.assembler_time;
249 res.total.timesteps += 1;
250 res.successful.assembler_time += step_result.assembler_time;
254 res.successful.timesteps += 1;
255 if (verbosityLevel>=1){
256 std::ios_base::fmtflags oldflags = std::cout.flags();
257 std::cout << "::: timesteps " << std::setw(6) << res.successful.timesteps
258 << " (" << res.total.timesteps << ")" << std::endl;
259 std::cout << "::: nl iterations " << std::setw(6) << res.successful.nonlinear_solver_iterations
260 << " (" << res.total.nonlinear_solver_iterations << ")" << std::endl;
261 std::cout << "::: lin iterations " << std::setw(6) << res.successful.linear_solver_iterations
262 << " (" << res.total.linear_solver_iterations << ")" << std::endl;
263 std::cout << "::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
264 << res.successful.assembler_time << " (" << res.total.assembler_time << ")" << std::endl;
265 std::cout << "::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
266 << res.successful.linear_solver_time << " (" << res.total.linear_solver_time << ")" << std::endl;
267 std::cout.flags(oldflags);
268 }
269
270 step++;
271 return dt;
272 }
273
284 template<typename F>
285 T apply (T time, T dt, TrlV& xold, F& f, TrlV& xnew)
286 {
287 // do statistics
288 OneStepMethodPartialResult step_result;
289
290 // save formatting attributes
291 ios_base_all_saver format_attribute_saver(std::cout);
292
293 std::vector<TrlV*> x(1); // vector of pointers to all steps
294 x[0] = &xold; // initially we have only one
295
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
300 << " time (from): "
301 << std::setw(12) << std::setprecision(4) << std::scientific
302 << time
303 << " dt: "
304 << std::setw(12) << std::setprecision(4) << std::scientific
305 << dt
306 << " time (to): "
307 << std::setw(12) << std::setprecision(4) << std::scientific
308 << time+dt
309 << std::endl;
310 std::cout.flags(oldflags);
311 }
312
313 // prepare assembler
314 igos.preStep(*method,time,dt);
315
316 // loop over all stages
317 for (unsigned r=1; r<=method->s(); ++r)
318 {
319 if (verbosityLevel>=2){
320 std::ios_base::fmtflags oldflags = std::cout.flags();
321 std::cout << "STAGE "
322 << r
323 << " time (to): "
324 << std::setw(12) << std::setprecision(4) << std::scientific
325 << time+method->d(r)*dt
326 << "." << std::endl;
327 std::cout.flags(oldflags);
328 }
329
330 // prepare stage
331 igos.preStage(r,x);
332
333 // get vector for current stage
334 if (r==method->s())
335 {
336 // last stage
337 x.push_back(&xnew);
338 }
339 else
340 {
341 // intermediate step
342 x.push_back(new TrlV(igos.trialGridFunctionSpace()));
343 }
344
345 // set boundary conditions and initial value
346 igos.interpolate(r,*x[r-1],f,*x[r]);
347
348 // solve stage
349 try {
350 pdesolver.apply(*x[r]);
351 }
352 catch (...)
353 {
354 // time step failed -> accumulate to total only
355 PDESolverResult pderes = pdesolver.result();
356 step_result.assembler_time += pderes.assembler_time;
357 step_result.linear_solver_time += pderes.linear_solver_time;
358 step_result.linear_solver_iterations += pderes.linear_solver_iterations;
359 step_result.nonlinear_solver_iterations += pderes.iterations;
360 res.total.assembler_time += step_result.assembler_time;
364 res.total.timesteps += 1;
365 throw;
366 }
367 PDESolverResult pderes = pdesolver.result();
368 step_result.assembler_time += pderes.assembler_time;
369 step_result.linear_solver_time += pderes.linear_solver_time;
370 step_result.linear_solver_iterations += pderes.linear_solver_iterations;
371 step_result.nonlinear_solver_iterations += pderes.iterations;
372
373 // stage cleanup
374 igos.postStage();
375 }
376
377 // delete intermediate steps
378 for (unsigned i=1; i<method->s(); ++i) delete x[i];
379
380 // step cleanup
381 igos.postStep();
382
383 // update statistics
384 res.total.assembler_time += step_result.assembler_time;
388 res.total.timesteps += 1;
389 res.successful.assembler_time += step_result.assembler_time;
393 res.successful.timesteps += 1;
394 if (verbosityLevel>=1){
395 std::ios_base::fmtflags oldflags = std::cout.flags();
396 std::cout << "::: timesteps " << std::setw(6) << res.successful.timesteps
397 << " (" << res.total.timesteps << ")" << std::endl;
398 std::cout << "::: nl iterations " << std::setw(6) << res.successful.nonlinear_solver_iterations
399 << " (" << res.total.nonlinear_solver_iterations << ")" << std::endl;
400 std::cout << "::: lin iterations " << std::setw(6) << res.successful.linear_solver_iterations
401 << " (" << res.total.linear_solver_iterations << ")" << std::endl;
402 std::cout << "::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
403 << res.successful.assembler_time << " (" << res.total.assembler_time << ")" << std::endl;
404 std::cout << "::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
405 << res.successful.linear_solver_time << " (" << res.total.linear_solver_time << ")" << std::endl;
406 std::cout.flags(oldflags);
407 }
408
409 step++;
410 return dt;
411 }
412
413 private:
415 IGOS& igos;
416 PDESOLVER& pdesolver;
417 int verbosityLevel;
418 int step;
419 Result res;
420 };
421
423 } // end namespace PDELab
424} // end namespace Dune
425#endif // DUNE_PDELAB_INSTATIONARY_IMPLICITONESTEP_HH
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:285
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
PDESOLVER & getPDESolver()
Access to the (non) linear solver.
Definition implicitonestep.hh:102
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
const PDESOLVER & getPDESolver() const
Access to the (non) linear solver.
Definition implicitonestep.hh:96
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition implicitonestep.hh:84
Base parameter class for time stepping scheme parameters.
Definition onestepparameter.hh:44