fvbaseproblem.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_FV_BASE_PROBLEM_HH
29#define EWOMS_FV_BASE_PROBLEM_HH
30
31#include <dune/common/fvector.hh>
32
36
39
41
42#include <functional>
43#include <iomanip>
44#include <iostream>
45#include <limits>
46#include <memory>
47#include <string>
48
49namespace Opm::Properties {
50
51template <class TypeTag, class MyTypeTag>
52struct NewtonMethod;
53
54} // namespace Opm::Properties
55
56namespace Opm {
57
67template<class TypeTag>
69{
70private:
71 using Implementation = GetPropType<TypeTag, Properties::Problem>;
73
74 static constexpr auto vtkOutputFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
76
82
85
90
91 enum {
92 dim = GridView::dimension,
93 dimWorld = GridView::dimensionworld
94 };
95
96 using Element = typename GridView::template Codim<0>::Entity;
97 using Vertex = typename GridView::template Codim<dim>::Entity;
98 using VertexIterator = typename GridView::template Codim<dim>::Iterator;
99
100 using CoordScalar = typename GridView::Grid::ctype;
101 using GlobalPosition = Dune::FieldVector<CoordScalar, dimWorld>;
102
103public:
104 // the default restriction and prolongation for adaptation is simply an empty one
106
107private:
108 // copying a problem is not a good idea
109 FvBaseProblem(const FvBaseProblem&) = delete;
110
111public:
119 explicit FvBaseProblem(Simulator& simulator)
120 : nextTimeStepSize_(0.0)
121 , gridView_(simulator.gridView())
122 , elementMapper_(gridView_, Dune::mcmgElementLayout())
123 , vertexMapper_(gridView_, Dune::mcmgVertexLayout())
124 , boundingBoxMin_(std::numeric_limits<double>::max())
125 , boundingBoxMax_(-std::numeric_limits<double>::max())
126 , simulator_(simulator)
127 {
128 // calculate the bounding box of the local partition of the grid view
129 for (const auto& vertex : vertices(gridView_)) {
130 for (unsigned i = 0; i < dim; ++i) {
131 boundingBoxMin_[i] = std::min(boundingBoxMin_[i], vertex.geometry().corner(0)[i]);
132 boundingBoxMax_[i] = std::max(boundingBoxMax_[i], vertex.geometry().corner(0)[i]);
133 }
134 }
135
136 // communicate to get the bounding box of the whole domain
137 for (unsigned i = 0; i < dim; ++i) {
138 boundingBoxMin_[i] = gridView_.comm().min(boundingBoxMin_[i]);
139 boundingBoxMax_[i] = gridView_.comm().max(boundingBoxMax_[i]);
140 }
141
142 if (enableVtkOutput_()) {
143 const bool asyncVtkOutput =
144 simulator_.gridView().comm().size() == 1 &&
145 Parameters::Get<Parameters::EnableAsyncVtkOutput>();
146
147 // asynchonous VTK output currently does not work in conjunction with grid
148 // adaptivity because the async-IO code assumes that the grid stays
149 // constant. complain about that case.
150 const bool enableGridAdaptation = Parameters::Get<Parameters::EnableGridAdaptation>();
151 if (asyncVtkOutput && enableGridAdaptation) {
152 throw std::runtime_error("Asynchronous VTK output currently cannot be used "
153 "at the same time as grid adaptivity");
154 }
155
156 defaultVtkWriter_ = std::make_unique<VtkMultiWriter>(asyncVtkOutput,
157 gridView_,
158 asImp_().outputDir(),
159 asImp_().name());
160 }
161 }
162
167 static void registerParameters()
168 {
169 Model::registerParameters();
170 Parameters::Register<Parameters::MaxTimeStepSize<Scalar>>
171 ("The maximum size to which all time steps are limited to [s]");
172 Parameters::Register<Parameters::MinTimeStepSize<Scalar>>
173 ("The minimum size to which all time steps are limited to [s]");
174 Parameters::Register<Parameters::MaxTimeStepDivisions>
175 ("The maximum number of divisions by two of the timestep size "
176 "before the simulation bails out");
177 Parameters::Register<Parameters::EnableAsyncVtkOutput>
178 ("Dispatch a separate thread to write the VTK output");
179 Parameters::Register<Parameters::ContinueOnConvergenceError>
180 ("Continue with a non-converged solution instead of giving up "
181 "if we encounter a time step size smaller than the minimum time "
182 "step size.");
183 }
184
193 { return true; }
194
203 {
204 if (Parameters::Get<Parameters::EnableStorageCache>() && asImp_().recycleFirstIterationStorage()) {
205 return 1;
206 }
207 return 2;
208 }
209
218 std::string outputDir() const
219 {
220 return simulatorOutputDir();
221 }
222
229 static std::string helpPreamble(int, const char **argv)
230 {
231 std::string desc = Implementation::briefDescription();
232 if (!desc.empty()) {
233 desc = desc + "\n";
234 }
235
236 return "Usage: " + std::string(argv[0]) + " [OPTIONS]\n" + desc;
237 }
238
245 static std::string briefDescription()
246 { return ""; }
247
248 // TODO (?): detailedDescription()
249
266 static int handlePositionalParameter(std::function<void(const std::string&,
267 const std::string&)>,
268 std::set<std::string>&,
269 std::string& errorMsg,
270 int,
271 const char** argv,
272 int paramIdx,
273 int)
274 {
275 errorMsg = std::string("Illegal parameter \"") + argv[paramIdx] + "\".";
276 return 0;
277 }
278
285 {}
286
291 void prefetch(const Element&) const
292 {
293 // do nothing by default
294 }
295
300 {
301 elementMapper_.update(gridView_);
302 vertexMapper_.update(gridView_);
303
304 if (enableVtkOutput_()) {
305 defaultVtkWriter_->gridChanged();
306 }
307 }
308
318 template <class Context>
319 void boundary(BoundaryRateVector&,
320 const Context&,
321 unsigned,
322 unsigned) const
323 { throw std::logic_error("Problem does not provide a boundary() method"); }
324
335 template <class Context>
336 void constraints(Constraints&,
337 const Context&,
338 unsigned,
339 unsigned) const
340 { throw std::logic_error("Problem does not provide a constraints() method"); }
341
354 template <class Context>
355 void source(RateVector&,
356 const Context&,
357 unsigned,
358 unsigned) const
359 { throw std::logic_error("Problem does not provide a source() method"); }
360
371 template <class Context>
372 void initial(PrimaryVariables&,
373 const Context&,
374 unsigned,
375 unsigned) const
376 { throw std::logic_error("Problem does not provide a initial() method"); }
377
393 template <class Context>
394 Scalar extrusionFactor(const Context&,
395 unsigned,
396 unsigned) const
397 { return asImp_().extrusionFactor(); }
398
399 Scalar extrusionFactor() const
400 { return 1.0; }
401
407 {}
408
413 {}
414
419 {}
420
425 {}
426
431 {}
432
440 {}
441
448 {
449 std::cerr << "The end of episode " << simulator().episodeIndex() + 1 << " has been "
450 << "reached, but the problem does not override the endEpisode() method. "
451 << "Doing nothing!\n";
452 }
453
457 void finalize()
458 {
459 const auto& executionTimer = simulator().executionTimer();
460
461 const Scalar executionTime = executionTimer.realTimeElapsed();
462 const Scalar setupTime = simulator().setupTimer().realTimeElapsed();
463 const Scalar prePostProcessTime = simulator().prePostProcessTimer().realTimeElapsed();
464 const Scalar localCpuTime = executionTimer.cpuTimeElapsed();
465 const Scalar globalCpuTime = executionTimer.globalCpuTimeElapsed();
466 const Scalar writeTime = simulator().writeTimer().realTimeElapsed();
467 const Scalar linearizeTime = simulator().linearizeTimer().realTimeElapsed();
468 const Scalar solveTime = simulator().solveTimer().realTimeElapsed();
469 const Scalar updateTime = simulator().updateTimer().realTimeElapsed();
470 const unsigned numProcesses = static_cast<unsigned>(this->gridView().comm().size());
471 const unsigned threadsPerProcess = ThreadManager::maxThreads();
472 if (gridView().comm().rank() == 0) {
473 std::cout << std::setprecision(3)
474 << "Simulation of problem '" << asImp_().name() << "' finished.\n"
475 << "\n"
476 << "------------------------ Timing ------------------------\n"
477 << "Setup time: " << setupTime << " seconds"
478 << humanReadableTime(setupTime)
479 << ", " << setupTime / (executionTime + setupTime) * 100 << "%\n"
480 << "Simulation time: " << executionTime << " seconds"
481 << humanReadableTime(executionTime)
482 << ", " << executionTime / (executionTime + setupTime) * 100 << "%\n"
483 << " Linearization time: " << linearizeTime << " seconds"
484 << humanReadableTime(linearizeTime)
485 << ", " << linearizeTime / executionTime * 100 << "%\n"
486 << " Linear solve time: " << solveTime << " seconds"
487 << humanReadableTime(solveTime)
488 << ", " << solveTime / executionTime * 100 << "%\n"
489 << " Newton update time: " << updateTime << " seconds"
490 << humanReadableTime(updateTime)
491 << ", " << updateTime / executionTime * 100 << "%\n"
492 << " Pre/postprocess time: " << prePostProcessTime << " seconds"
493 << humanReadableTime(prePostProcessTime)
494 << ", " << prePostProcessTime / executionTime * 100 << "%\n"
495 << " Output write time: " << writeTime << " seconds"
496 << humanReadableTime(writeTime)
497 << ", " << writeTime / executionTime * 100 << "%\n"
498 << "First process' simulation CPU time: " << localCpuTime << " seconds"
499 << humanReadableTime(localCpuTime) << "\n"
500 << "Number of processes: " << numProcesses << "\n"
501 << "Threads per processes: " << threadsPerProcess << "\n"
502 << "Total CPU time: " << globalCpuTime << " seconds"
503 << humanReadableTime(globalCpuTime) << "\n"
504 << "\n"
505 << "----------------------------------------------------------------\n"
506 << std::endl;
507 }
508 }
509
515 {
516 const unsigned maxFails = asImp_().maxTimeIntegrationFailures();
517 Scalar minTimeStep = asImp_().minTimeStepSize();
518
519 std::string errorMessage;
520 for (unsigned i = 0; i < maxFails; ++i) {
521 bool converged = model().update();
522 if (converged) {
523 return;
524 }
525
526 const Scalar dt = simulator().timeStepSize();
527 Scalar nextDt = dt / 2.0;
528 if (dt < minTimeStep * (1.0 + 1e-9)) {
529 if (asImp_().continueOnConvergenceError()) {
530 if (gridView().comm().rank() == 0) {
531 std::cout << "Newton solver did not converge with minimum time step of "
532 << dt << " seconds. Continuing with unconverged solution!\n"
533 << std::flush;
534 }
535 return;
536 }
537 else {
538 errorMessage = "Time integration did not succeed with the minumum time step size of " +
539 std::to_string(double(minTimeStep)) + " seconds. Giving up!";
540 break; // give up: we can't make the time step smaller anymore!
541 }
542 }
543 else if (nextDt < minTimeStep) {
544 nextDt = minTimeStep;
545 }
546 simulator().setTimeStepSize(nextDt);
547
548 // update failed
549 if (gridView().comm().rank() == 0) {
550 std::cout << "Newton solver did not converge with "
551 << "dt=" << dt << " seconds. Retrying with time step of "
552 << nextDt << " seconds\n" << std::flush;
553 }
554 }
555
556 if (errorMessage.empty()) {
557 errorMessage = "Newton solver didn't converge after " +
558 std::to_string(maxFails) + " time-step divisions. dt=" +
559 std::to_string(double(simulator().timeStepSize()));
560 }
561 throw std::runtime_error(errorMessage);
562 }
563
567 Scalar minTimeStepSize() const
568 { return Parameters::Get<Parameters::MinTimeStepSize<Scalar>>(); }
569
575 { return Parameters::Get<Parameters::MaxTimeStepDivisions>(); }
576
583 { return Parameters::Get<Parameters::ContinueOnConvergenceError>(); }
584
588 void setNextTimeStepSize(Scalar dt)
589 { nextTimeStepSize_ = dt; }
590
596 Scalar nextTimeStepSize() const
597 {
598 if (nextTimeStepSize_ > 0.0) {
599 return nextTimeStepSize_;
600 }
601
602 Scalar dtNext = std::min(Parameters::Get<Parameters::MaxTimeStepSize<Scalar>>(),
603 newtonMethod().suggestTimeStepSize(simulator().timeStepSize()));
604
605 if (dtNext < simulator().maxTimeStepSize() &&
606 simulator().maxTimeStepSize() < dtNext * 2)
607 {
608 dtNext = simulator().maxTimeStepSize() / 2 * 1.01;
609 }
610
611 return dtNext;
612 }
613
623 {
624 return simulator().timeStepIndex() > 0 &&
625 (simulator().timeStepIndex() % 10 == 0);
626 }
627
636 bool shouldWriteOutput() const
637 { return true; }
638
645 { model().advanceTimeLevel(); }
646
654 std::string name() const
655 { return "sim"; }
656
660 const GridView& gridView() const
661 { return gridView_; }
662
667 const GlobalPosition& boundingBoxMin() const
668 { return boundingBoxMin_; }
669
674 const GlobalPosition& boundingBoxMax() const
675 { return boundingBoxMax_; }
676
680 const VertexMapper& vertexMapper() const
681 { return vertexMapper_; }
682
686 const ElementMapper& elementMapper() const
687 { return elementMapper_; }
688
692 Simulator& simulator()
693 { return simulator_; }
694
698 const Simulator& simulator() const
699 { return simulator_; }
700
704 Model& model()
705 { return simulator_.model(); }
706
710 const Model& model() const
711 { return simulator_.model(); }
712
716 NewtonMethod& newtonMethod()
717 { return model().newtonMethod(); }
718
722 const NewtonMethod& newtonMethod() const
723 { return model().newtonMethod(); }
724 // \}
725
731 {
733 }
734
743 {
744 return 0;
745 }
746
760 template <class Restarter>
761 void serialize(Restarter& res)
762 {
763 if (enableVtkOutput_()) {
764 defaultVtkWriter_->serialize(res);
765 }
766 }
767
778 template <class Restarter>
779 void deserialize(Restarter& res)
780 {
781 if (enableVtkOutput_()) {
782 defaultVtkWriter_->deserialize(res);
783 }
784 }
785
792 void writeOutput(bool verbose = true)
793 {
794 if (!enableVtkOutput_()) {
795 return;
796 }
797
798 if (verbose && gridView().comm().rank() == 0) {
799 std::cout << "Writing visualization results for the current time step.\n"
800 << std::flush;
801 }
802
803 // calculate the time _after_ the time was updated
804 const Scalar t = simulator().time() + simulator().timeStepSize();
805
806 defaultVtkWriter_->beginWrite(t);
807 model().prepareOutputFields();
808 model().appendOutputFields(*defaultVtkWriter_);
809 defaultVtkWriter_->endWrite(false);
810 }
811
817 { return *defaultVtkWriter_; }
818
819protected:
821
822 bool enableVtkOutput_() const
823 { return Parameters::Get<Parameters::EnableVtkOutput>(); }
824
825private:
827 Implementation& asImp_()
828 { return *static_cast<Implementation*>(this); }
829
831 const Implementation& asImp_() const
832 { return *static_cast<const Implementation*>(this); }
833
834 // Grid management stuff
835 const GridView gridView_;
836 ElementMapper elementMapper_;
837 VertexMapper vertexMapper_;
838 GlobalPosition boundingBoxMin_;
839 GlobalPosition boundingBoxMax_;
840
841 // Attributes required for the actual simulation
842 Simulator& simulator_;
843 std::unique_ptr<VtkMultiWriter> defaultVtkWriter_{};
844};
845
846} // namespace Opm
847
848#endif
Definition: restrictprolong.hh:142
Base class for all problems which use a finite volume spatial discretization.
Definition: fvbaseproblem.hh:69
static int handlePositionalParameter(std::function< void(const std::string &, const std::string &)>, std::set< std::string > &, std::string &errorMsg, int, const char **argv, int paramIdx, int)
Handles positional command line parameters.
Definition: fvbaseproblem.hh:266
std::string name() const
The problem name.
Definition: fvbaseproblem.hh:654
const Simulator & simulator() const
Returns Simulator object used by the simulation.
Definition: fvbaseproblem.hh:698
void writeOutput(bool verbose=true)
Write the relevant secondary variables of the current solution into an VTK output file.
Definition: fvbaseproblem.hh:792
void boundary(BoundaryRateVector &, const Context &, unsigned, unsigned) const
Evaluate the boundary conditions for a boundary segment.
Definition: fvbaseproblem.hh:319
void finalize()
Called after the simulation has been run sucessfully.
Definition: fvbaseproblem.hh:457
Simulator & simulator()
Returns Simulator object used by the simulation.
Definition: fvbaseproblem.hh:692
const Model & model() const
Returns numerical model used for the problem.
Definition: fvbaseproblem.hh:710
static void registerParameters()
Registers all available parameters for the problem and the model.
Definition: fvbaseproblem.hh:167
void source(RateVector &, const Context &, unsigned, unsigned) const
Evaluate the source term for all phases within a given sub-control-volume.
Definition: fvbaseproblem.hh:355
bool recycleFirstIterationStorage() const
Return if the storage term of the first iteration is identical to the storage term for the solution o...
Definition: fvbaseproblem.hh:192
unsigned maxTimeIntegrationFailures() const
Returns the maximum number of subsequent failures for the time integration before giving up.
Definition: fvbaseproblem.hh:574
unsigned markForGridAdaptation()
Mark grid cells for refinement or coarsening.
Definition: fvbaseproblem.hh:742
void serialize(Restarter &res)
This method writes the complete state of the problem to the harddisk.
Definition: fvbaseproblem.hh:761
void beginIteration()
Called by the simulator before each Newton-Raphson iteration.
Definition: fvbaseproblem.hh:424
unsigned intensiveQuantityHistorySize() const
Returns the required history size for intensive quantities cache.
Definition: fvbaseproblem.hh:202
void endTimeStep()
Called by the simulator after each time integration.
Definition: fvbaseproblem.hh:439
FvBaseProblem(Simulator &simulator)
Definition: fvbaseproblem.hh:119
const GlobalPosition & boundingBoxMax() const
The coordinate of the corner of the GridView's bounding box with the largest values.
Definition: fvbaseproblem.hh:674
RestrictProlongOperator restrictProlongOperator()
return restriction and prolongation operator
Definition: fvbaseproblem.hh:730
void timeIntegration()
Called by Opm::Simulator in order to do a time integration on the model.
Definition: fvbaseproblem.hh:514
Model & model()
Returns numerical model used for the problem.
Definition: fvbaseproblem.hh:704
Scalar extrusionFactor() const
Definition: fvbaseproblem.hh:399
void endEpisode()
Called when the end of an simulation episode is reached.
Definition: fvbaseproblem.hh:447
void beginTimeStep()
Called by the simulator before each time integration.
Definition: fvbaseproblem.hh:418
static std::string helpPreamble(int, const char **argv)
Returns the string that is printed before the list of command line parameters in the help message.
Definition: fvbaseproblem.hh:229
void deserialize(Restarter &res)
This method restores the complete state of the problem from disk.
Definition: fvbaseproblem.hh:779
void prefetch(const Element &) const
Allows to improve the performance by prefetching all data which is associated with a given element.
Definition: fvbaseproblem.hh:291
Scalar extrusionFactor(const Context &, unsigned, unsigned) const
Return how much the domain is extruded at a given sub-control volume.
Definition: fvbaseproblem.hh:394
bool shouldWriteRestartFile() const
Returns true if a restart file should be written to disk.
Definition: fvbaseproblem.hh:622
Scalar nextTimeStepSize() const
Called by Opm::Simulator whenever a solution for a time step has been computed and the simulation tim...
Definition: fvbaseproblem.hh:596
Scalar minTimeStepSize() const
Returns the minimum allowable size of a time step.
Definition: fvbaseproblem.hh:567
void beginEpisode()
Called at the beginning of an simulation episode.
Definition: fvbaseproblem.hh:412
void initial(PrimaryVariables &, const Context &, unsigned, unsigned) const
Evaluate the initial value for a control volume.
Definition: fvbaseproblem.hh:372
std::string outputDir() const
Determine the directory for simulation output.
Definition: fvbaseproblem.hh:218
Scalar nextTimeStepSize_
Definition: fvbaseproblem.hh:820
const VertexMapper & vertexMapper() const
Returns the mapper for vertices to indices.
Definition: fvbaseproblem.hh:680
void gridChanged()
Handle changes of the grid.
Definition: fvbaseproblem.hh:299
const GlobalPosition & boundingBoxMin() const
The coordinate of the corner of the GridView's bounding box with the smallest values.
Definition: fvbaseproblem.hh:667
static std::string briefDescription()
Returns a human readable description of the problem for the help message.
Definition: fvbaseproblem.hh:245
void endIteration()
Called by the simulator after each Newton-Raphson update.
Definition: fvbaseproblem.hh:430
void setNextTimeStepSize(Scalar dt)
Impose the next time step size to be used externally.
Definition: fvbaseproblem.hh:588
bool continueOnConvergenceError() const
Returns if we should continue with a non-converged solution instead of giving up if we encounter a ti...
Definition: fvbaseproblem.hh:582
void constraints(Constraints &, const Context &, unsigned, unsigned) const
Evaluate the constraints for a control volume.
Definition: fvbaseproblem.hh:336
void finishInit()
Called by the Opm::Simulator in order to initialize the problem.
Definition: fvbaseproblem.hh:284
const NewtonMethod & newtonMethod() const
Returns object which implements the Newton method.
Definition: fvbaseproblem.hh:722
const GridView & gridView() const
The GridView which used by the problem.
Definition: fvbaseproblem.hh:660
bool enableVtkOutput_() const
Definition: fvbaseproblem.hh:822
EmptyRestrictProlong RestrictProlongOperator
Definition: fvbaseproblem.hh:105
NewtonMethod & newtonMethod()
Returns object which implements the Newton method.
Definition: fvbaseproblem.hh:716
bool shouldWriteOutput() const
Returns true if the current solution should be written to disk (i.e. as a VTK file)
Definition: fvbaseproblem.hh:636
void initialSolutionApplied()
Callback used by the model to indicate that the initial solution has been determined for all degrees ...
Definition: fvbaseproblem.hh:406
void advanceTimeLevel()
Called by the simulator after everything which can be done about the current time step is finished an...
Definition: fvbaseproblem.hh:644
VtkMultiWriter & defaultVtkWriter() const
Method to retrieve the VTK writer which should be used to write the default ouput after each time ste...
Definition: fvbaseproblem.hh:816
const ElementMapper & elementMapper() const
Returns the mapper for elements to indices.
Definition: fvbaseproblem.hh:686
static unsigned maxThreads()
Return the maximum number of threads of the current process.
Definition: threadmanager.hpp:66
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:65
Declare the properties used by the infrastructure code of the finite volume discretizations.
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: fvbaseprimaryvariables.hh:141
auto Get(bool errorIfNotRegistered=true)
Retrieve a runtime parameter.
Definition: parametersystem.hpp:187
Definition: blackoilmodel.hh:79
Definition: blackoilboundaryratevector.hh:39
std::string humanReadableTime(double timeInSeconds, bool isAmendment=true)
Given a time step size in seconds, return it in a format which is more easily parsable by humans.
std::string simulatorOutputDir()
Determine and check the configured directory for simulation output.
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
Specify the maximum size of a time integration [s].
Definition: fvbaseparameters.hh:106