28 #ifndef EWOMS_FV_BASE_PROBLEM_HH 29 #define EWOMS_FV_BASE_PROBLEM_HH 35 #include <ewoms/disc/common/restrictprolong.hh> 37 #include <opm/common/Unused.hpp> 38 #include <opm/common/ErrorMacros.hpp> 39 #include <opm/common/Exceptions.hpp> 41 #include <dune/common/fvector.hh> 58 template<
class TypeTag>
62 typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
63 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
65 static const int vtkOutputFormat =
GET_PROP_VALUE(TypeTag, VtkOutputFormat);
68 typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
69 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
74 typedef typename GET_PROP_TYPE(TypeTag, VertexMapper) VertexMapper;
75 typedef typename GET_PROP_TYPE(TypeTag, ElementMapper) ElementMapper;
77 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
78 typedef typename GET_PROP_TYPE(TypeTag, BoundaryRateVector) BoundaryRateVector;
79 typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
80 typedef typename GET_PROP_TYPE(TypeTag, Constraints) Constraints;
83 dim = GridView::dimension,
84 dimWorld = GridView::dimensionworld
87 typedef typename GridView::template Codim<0>::Entity Element;
88 typedef typename GridView::template Codim<dim>::Entity Vertex;
89 typedef typename GridView::template Codim<dim>::Iterator VertexIterator;
91 typedef typename GridView::Grid::ctype CoordScalar;
92 typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
112 , elementMapper_(gridView_)
113 , vertexMapper_(gridView_)
114 , boundingBoxMin_(std::numeric_limits<double>::max())
115 , boundingBoxMax_(-std::numeric_limits<double>::max())
117 , defaultVtkWriter_(0)
120 VertexIterator vIt = gridView_.template begin<dim>();
121 const VertexIterator vEndIt = gridView_.template end<dim>();
122 for (; vIt!=vEndIt; ++vIt) {
123 for (
unsigned i=0; i<dim; i++) {
124 boundingBoxMin_[i] = std::min(boundingBoxMin_[i], vIt->geometry().corner(0)[i]);
125 boundingBoxMax_[i] = std::max(boundingBoxMax_[i], vIt->geometry().corner(0)[i]);
130 for (
unsigned i = 0; i < dim; ++i) {
131 boundingBoxMin_[i] = gridView_.comm().min(boundingBoxMin_[i]);
132 boundingBoxMax_[i] = gridView_.comm().max(boundingBoxMax_[i]);
135 if (enableVtkOutput_())
140 {
delete defaultVtkWriter_; }
148 Model::registerParameters();
150 "The maximum size to which all time steps are limited to [s]");
152 "The minimum size to which all time steps are limited to [s]");
154 "The maximum number of divisions by two of the timestep size " 155 "before the simulation bails out");
170 void prefetch(
const Element& elem OPM_UNUSED)
const 180 elementMapper_.update();
181 vertexMapper_.update();
183 if (enableVtkOutput_())
196 template <
class Context>
197 void boundary(BoundaryRateVector& values OPM_UNUSED,
198 const Context& context OPM_UNUSED,
199 unsigned spaceIdx OPM_UNUSED,
200 unsigned timeIdx OPM_UNUSED)
const 201 { OPM_THROW(std::logic_error,
"Problem does not provide a boundary() method"); }
213 template <
class Context>
215 const Context& context OPM_UNUSED,
216 unsigned spaceIdx OPM_UNUSED,
217 unsigned timeIdx OPM_UNUSED)
const 218 { OPM_THROW(std::logic_error,
"Problem does not provide a constraints() method"); }
232 template <
class Context>
234 const Context& context OPM_UNUSED,
235 unsigned spaceIdx OPM_UNUSED,
236 unsigned timeIdx OPM_UNUSED)
const 237 { OPM_THROW(std::logic_error,
"Problem does not provide a source() method"); }
249 template <
class Context>
250 void initial(PrimaryVariables& values OPM_UNUSED,
251 const Context& context OPM_UNUSED,
252 unsigned spaceIdx OPM_UNUSED,
253 unsigned timeIdx OPM_UNUSED)
const 254 { OPM_THROW(std::logic_error,
"Problem does not provide a initial() method"); }
271 template <
class Context>
273 unsigned spaceIdx OPM_UNUSED,
274 unsigned timeIdx OPM_UNUSED)
const 275 {
return asImp_().extrusionFactor(); }
277 Scalar extrusionFactor()
const 328 <<
"reached, but the problem does not override the endEpisode() method. " 329 <<
"Doing nothing!\n";
342 Scalar localCpuTime = executionTimer.cpuTimeElapsed();
343 Scalar globalCpuTime = executionTimer.globalCpuTimeElapsed();
348 unsigned numProcesses =
static_cast<unsigned>(this->
gridView().comm().size());
350 if (
gridView().comm().rank() == 0) {
351 std::cout << std::setprecision(3)
352 <<
"Simulation of problem '" << asImp_().name() <<
"' finished.\n" 354 <<
"------------------------ Timing receipt ------------------------\n" 356 <<
", " << setupTime/(executionTime + setupTime)*100 <<
"%\n" 358 <<
", " << executionTime/(executionTime + setupTime)*100 <<
"%\n" 360 <<
", " << linearizeTime/executionTime*100 <<
"%\n" 362 <<
", " << solveTime/executionTime*100 <<
"%\n" 364 <<
", " << updateTime/executionTime*100 <<
"%\n" 366 <<
", " << prePostProcessTime/executionTime*100 <<
"%\n" 368 <<
", " << writeTime/executionTime*100 <<
"%\n" 370 <<
"Number of processes: " << numProcesses <<
"\n" 371 <<
"Threads per processes: " << threadsPerProcess <<
"\n" 374 <<
"Note 1: If not stated otherwise, all times are wall clock times\n" 375 <<
"Note 2: Taxes and administrative overhead are " 376 << (executionTime - (linearizeTime+solveTime+updateTime+prePostProcessTime+writeTime))/executionTime*100
379 <<
"Our simulation hours are 24/7. Thank you for choosing us.\n" 380 <<
"----------------------------------------------------------------\n" 392 unsigned maxFails =
EWOMS_GET_PARAM(TypeTag,
unsigned, MaxTimeStepDivisions);
393 Scalar minTimeStepSize =
EWOMS_GET_PARAM(TypeTag, Scalar, MinTimeStepSize);
398 if (
simulator().timeStepSize() < minTimeStepSize &&
405 for (
unsigned i = 0; i < maxFails; ++i) {
406 bool converged =
model().update();
411 Scalar nextDt = dt / 2;
412 if (nextDt < minTimeStepSize)
418 std::cout <<
"Newton solver did not converge with " 419 <<
"dt=" << dt <<
" seconds. Retrying with time step of " 420 << nextDt <<
" seconds\n" << std::flush;
423 OPM_THROW(std::runtime_error,
424 "Newton solver didn't converge after " 425 << maxFails <<
" time-step divisions. dt=" 436 Scalar dtNext = std::min(
EWOMS_GET_PARAM(TypeTag, Scalar, MaxTimeStepSize),
439 if (dtNext <
simulator().maxTimeStepSize()
440 &&
simulator().maxTimeStepSize() < dtNext*2)
479 {
model().advanceTimeLevel(); }
495 {
return gridView_; }
502 {
return boundingBoxMin_; }
509 {
return boundingBoxMax_; }
515 {
return vertexMapper_; }
521 {
return elementMapper_; }
527 {
return simulator_; }
533 {
return simulator_; }
539 {
return simulator_.
model(); }
545 {
return simulator_.
model(); }
551 {
return model().newtonMethod(); }
557 {
return model().newtonMethod(); }
594 template <
class Restarter>
597 if (enableVtkOutput_())
611 template <
class Restarter>
614 if (enableVtkOutput_())
626 if (verbose &&
gridView().comm().rank() == 0)
627 std::cout <<
"Writing visualization results for the current time step.\n" 633 if (enableVtkOutput_())
636 model().prepareOutputFields();
638 if (enableVtkOutput_()) {
639 model().appendOutputFields(*defaultVtkWriter_);
649 {
return defaultVtkWriter_; }
652 bool enableVtkOutput_()
const 656 Implementation& asImp_()
657 {
return *
static_cast<Implementation *
>(
this); }
660 const Implementation& asImp_()
const 661 {
return *
static_cast<const Implementation *
>(
this); }
664 const GridView gridView_;
665 ElementMapper elementMapper_;
666 VertexMapper vertexMapper_;
667 GlobalPosition boundingBoxMin_;
668 GlobalPosition boundingBoxMax_;
671 Simulator& simulator_;
672 mutable VtkMultiWriter *defaultVtkWriter_;
void finishInit()
Called by the Ewoms::Simulator in order to initialize the problem.
Definition: fvbaseproblem.hh:163
Scalar maxTimeStepSize() const
Aligns the time step size to the episode boundary and to the end time of the simulation.
Definition: simulator.hh:403
void endEpisode()
Called when the end of an simulation episode is reached.
Definition: fvbaseproblem.hh:325
void serialize(Restarter &res)
Write the multi-writer's state to a restart file.
Definition: vtkmultiwriter.hh:378
Simplifies multi-threaded capabilities.
Definition: threadmanager.hh:52
void beginWrite(double t)
Called whenever a new time step must be written.
Definition: vtkmultiwriter.hh:135
unsigned markForGridAdaptation()
Mark grid cells for refinement or coarsening.
Definition: fvbaseproblem.hh:576
void serialize(Restarter &res)
This method writes the complete state of the problem to the harddisk.
Definition: fvbaseproblem.hh:595
Definition: baseauxiliarymodule.hh:37
void constraints(Constraints &constraints OPM_UNUSED, const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED) const
Evaluate the constraints for a control volume.
Definition: fvbaseproblem.hh:214
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:189
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:63
int timeStepIndex() const
Returns number of time steps which have been executed since the beginning of the simulation.
Definition: simulator.hh:360
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
void boundary(BoundaryRateVector &values OPM_UNUSED, const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED) const
Evaluate the boundary conditions for a boundary segment.
Definition: fvbaseproblem.hh:197
void endWrite(bool onlyDiscard=false)
Finalizes the current writer.
Definition: vtkmultiwriter.hh:340
void endIteration()
Called by the simulator after each Newton-Raphson update.
Definition: fvbaseproblem.hh:308
void beginIteration()
Called by the simulator before each Newton-Raphson iteration.
Definition: fvbaseproblem.hh:302
void beginEpisode()
Called at the beginning of an simulation episode.
Definition: fvbaseproblem.hh:290
RestrictProlongOperator restrictProlongOperator()
return restriction and prolongation operator
Definition: fvbaseproblem.hh:564
Load or save a state of a problem to/from the harddisk.
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:68
const Ewoms::Timer & solveTimer() const
Returns a reference to the timer object which measures the time needed by the solver.
Definition: simulator.hh:304
const Model & model() const
Returns numerical model used for the problem.
Definition: fvbaseproblem.hh:544
void deserialize(Restarter &res)
Read the multi-writer's state from a restart file.
Definition: vtkmultiwriter.hh:412
Declare the properties used by the infrastructure code of the finite volume discretizations.
FvBaseProblem(Simulator &simulator)
Definition: fvbaseproblem.hh:110
void source(RateVector &rate OPM_UNUSED, const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED) const
Evaluate the source term for all phases within a given sub-control-volume.
Definition: fvbaseproblem.hh:233
void endTimeStep()
Called by the simulator after each time integration.
Definition: fvbaseproblem.hh:317
Scalar extrusionFactor(const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED) const
Return how much the domain is extruded at a given sub-control volume.
Definition: fvbaseproblem.hh:272
void deserialize(Restarter &res)
This method restores the complete state of the problem from disk.
Definition: fvbaseproblem.hh:612
double realTimeElapsed() const
Return the real time [s] elapsed during the periods the timer was active since the last reset...
Definition: timer.hh:121
Definition: restrictprolong.hh:140
void timeIntegration()
Called by Ewoms::Simulator in order to do a time integration on the model.
Definition: fvbaseproblem.hh:390
const Ewoms::Timer & prePostProcessTimer() const
Returns a reference to the timer object which measures the time needed for pre- and postprocessing of...
Definition: simulator.hh:290
void beginTimeStep()
Called by the simulator before each time integration.
Definition: fvbaseproblem.hh:296
static unsigned maxThreads()
Return the maximum number of threads of the current process.
Definition: threadmanager.hh:113
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
The multi-dimensional Newton method.
Definition: newtonmethod.hh:56
void advanceTimeLevel()
Called by the simulator after everything which can be done about the current time step is finished an...
Definition: fvbaseproblem.hh:478
void setTimeStepSize(Scalar value)
Set the current time step size to a given value.
Definition: simulator.hh:331
bool shouldWriteRestartFile() const
Returns true if a restart file should be written to disk.
Definition: fvbaseproblem.hh:456
const VertexMapper & vertexMapper() const
Returns the mapper for vertices to indices.
Definition: fvbaseproblem.hh:514
const Ewoms::Timer & writeTimer() const
Returns a reference to the timer object which measures the time needed to write the visualization out...
Definition: simulator.hh:318
void finalize()
Called after the simulation has been run sucessfully.
Definition: fvbaseproblem.hh:335
Simplifies writing multi-file VTK datasets.
const Ewoms::Timer & updateTimer() const
Returns a reference to the timer object which measures the time needed to the solutions of the non-li...
Definition: simulator.hh:311
const GridView & gridView() const
The GridView which used by the problem.
Definition: fvbaseproblem.hh:494
void prefetch(const Element &elem OPM_UNUSED) const
Allows to improve the performance by prefetching all data which is associated with a given element...
Definition: fvbaseproblem.hh:170
void gridChanged()
Updates the internal data structures after mesh refinement.
Definition: vtkmultiwriter.hh:126
int episodeIndex() const
Returns the index of the current episode.
Definition: simulator.hh:452
void writeOutput(bool verbose=true)
Write the relevant secondary variables of the current solution into an VTK output file...
Definition: fvbaseproblem.hh:624
static void registerParameters()
Registers all available parameters for the problem and the model.
Definition: fvbaseproblem.hh:146
Simulator & simulator()
Returns Simulator object used by the simulation.
Definition: fvbaseproblem.hh:526
const Ewoms::Timer & linearizeTimer() const
Returns a reference to the timer object which measures the time needed for linarizing the solutions...
Definition: simulator.hh:297
Scalar nextTimeStepSize()
Called by Ewoms::Simulator whenever a solution for a time step has been computed and the simulation t...
Definition: fvbaseproblem.hh:434
Scalar timeStepSize() const
Returns the time step length so that we don't miss the beginning of the next episode or cross the en...
Definition: simulator.hh:347
std::string name() const
The problem name.
Definition: fvbaseproblem.hh:488
bool shouldWriteOutput() const
Returns true if the current solution should be written to disk (i.e.
Definition: fvbaseproblem.hh:470
static std::string humanReadableTime(Scalar timeInSeconds, bool isAmendment=true)
Given a time step size in seconds, return it in a format which is more easily parsable by humans...
Definition: simulator.hh:707
const Ewoms::Timer & executionTimer() const
Returns a reference to the timer object which measures the time needed to run the simulation...
Definition: simulator.hh:283
const NewtonMethod & newtonMethod() const
Returns object which implements the Newton method.
Definition: fvbaseproblem.hh:556
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:75
Model & model()
Returns numerical model used for the problem.
Definition: fvbaseproblem.hh:538
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:648
const Ewoms::Timer & setupTimer() const
Returns a reference to the timer object which measures the time needed to set up and initialize the s...
Definition: simulator.hh:276
Base class for all problems which use a finite volume spatial discretization.
Definition: fvbaseproblem.hh:59
void initialSolutionApplied()
Callback used by the model to indicate that the initial solution has been determined for all degrees ...
Definition: fvbaseproblem.hh:284
void gridChanged()
Handle changes of the grid.
Definition: fvbaseproblem.hh:178
const GlobalPosition & boundingBoxMin() const
The coordinate of the corner of the GridView's bounding box with the smallest values.
Definition: fvbaseproblem.hh:501
const Simulator & simulator() const
Returns Simulator object used by the simulation.
Definition: fvbaseproblem.hh:532
const ElementMapper & elementMapper() const
Returns the mapper for elements to indices.
Definition: fvbaseproblem.hh:520
const GlobalPosition & boundingBoxMax() const
The coordinate of the corner of the GridView's bounding box with the largest values.
Definition: fvbaseproblem.hh:508
Scalar time() const
Return the number of seconds of simulated time which have elapsed since the start time...
Definition: simulator.hh:254
void initial(PrimaryVariables &values OPM_UNUSED, const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED) const
Evaluate the initial value for a control volume.
Definition: fvbaseproblem.hh:250
NewtonMethod & newtonMethod()
Returns object which implements the Newton method.
Definition: fvbaseproblem.hh:550