28#ifndef EWOMS_FV_BASE_PROBLEM_HH
29#define EWOMS_FV_BASE_PROBLEM_HH
31#include <dune/common/fvector.hh>
35#include <opm/models/discretization/common/restrictprolong.hh>
40#include <opm/models/utils/simulatorutils.hpp>
47namespace Opm::Properties {
49template <
class TypeTag,
class MyTypeTag>
65template<
class TypeTag>
90 dim = GridView::dimension,
91 dimWorld = GridView::dimensionworld
94 using Element =
typename GridView::template Codim<0>::Entity;
95 using Vertex =
typename GridView::template Codim<dim>::Entity;
96 using VertexIterator =
typename GridView::template Codim<dim>::Iterator;
98 using CoordScalar =
typename GridView::Grid::ctype;
99 using GlobalPosition = Dune::FieldVector<CoordScalar, dimWorld>;
118 : nextTimeStepSize_(0.0)
125 , defaultVtkWriter_(0)
131 for (
unsigned i=0; i<dim; i++) {
132 boundingBoxMin_[i] = std::min(boundingBoxMin_[i],
vIt->geometry().corner(0)[i]);
133 boundingBoxMax_[i] = std::max(boundingBoxMax_[i],
vIt->geometry().corner(0)[i]);
138 for (
unsigned i = 0; i < dim; ++i) {
139 boundingBoxMin_[i] = gridView_.comm().min(boundingBoxMin_[i]);
140 boundingBoxMax_[i] = gridView_.comm().max(boundingBoxMax_[i]);
143 if (enableVtkOutput_()) {
145 simulator_.gridView().comm().size() == 1 &&
146 Parameters::Get<Parameters::EnableAsyncVtkOutput>();
151 bool enableGridAdaptation = Parameters::Get<Parameters::EnableGridAdaptation>();
153 throw std::runtime_error(
"Asynchronous VTK output currently cannot be used "
154 "at the same time as grid adaptivity");
156 std::string
outputDir = asImp_().outputDir();
164 {
delete defaultVtkWriter_; }
172 Model::registerParameters();
173 Parameters::Register<Parameters::MaxTimeStepSize<Scalar>>
174 (
"The maximum size to which all time steps are limited to [s]");
175 Parameters::Register<Parameters::MinTimeStepSize<Scalar>>
176 (
"The minimum size to which all time steps are limited to [s]");
177 Parameters::Register<Parameters::MaxTimeStepDivisions>
178 (
"The maximum number of divisions by two of the timestep size "
179 "before the simulation bails out");
180 Parameters::Register<Parameters::EnableAsyncVtkOutput>
181 (
"Dispatch a separate thread to write the VTK output");
182 Parameters::Register<Parameters::ContinueOnConvergenceError>
183 (
"Continue with a non-converged solution instead of giving up "
184 "if we encounter a time step size smaller than the minimum time "
220 std::string
desc = Implementation::briefDescription();
225 "Usage: "+std::string(
argv[0]) +
" [OPTIONS]\n"
257 const std::string&)>,
258 std::set<std::string>&,
291 elementMapper_.update(gridView_);
292 vertexMapper_.update(gridView_);
294 if (enableVtkOutput_())
307 template <
class Context>
312 {
throw std::logic_error(
"Problem does not provide a boundary() method"); }
324 template <
class Context>
329 {
throw std::logic_error(
"Problem does not provide a constraints() method"); }
343 template <
class Context>
348 {
throw std::logic_error(
"Problem does not provide a source() method"); }
360 template <
class Context>
365 {
throw std::logic_error(
"Problem does not provide a initial() method"); }
382 template <
class Context>
386 {
return asImp_().extrusionFactor(); }
388 Scalar extrusionFactor()
const
438 std::cerr <<
"The end of episode " <<
simulator().episodeIndex() + 1 <<
" has been "
439 <<
"reached, but the problem does not override the endEpisode() method. "
440 <<
"Doing nothing!\n";
448 const auto& executionTimer =
simulator().executionTimer();
451 Scalar setupTime =
simulator().setupTimer().realTimeElapsed();
454 Scalar
globalCpuTime = executionTimer.globalCpuTimeElapsed();
461 if (
gridView().comm().rank() == 0) {
462 std::cout << std::setprecision(3)
463 <<
"Simulation of problem '" << asImp_().name() <<
"' finished.\n"
465 <<
"------------------------ Timing ------------------------\n"
466 <<
"Setup time: " << setupTime <<
" seconds"
468 <<
", " << setupTime/(
executionTime + setupTime)*100 <<
"%\n"
475 <<
" Linear solve time: " <<
solveTime <<
" seconds"
478 <<
" Newton update time: " <<
updateTime <<
" seconds"
484 <<
" Output write time: " <<
writeTime <<
" seconds"
487 <<
"First process' simulation CPU time: " <<
localCpuTime <<
" seconds"
494 <<
"----------------------------------------------------------------\n"
505 unsigned maxFails = asImp_().maxTimeIntegrationFailures();
509 for (
unsigned i = 0; i <
maxFails; ++i) {
510 bool converged =
model().update();
519 std::cout <<
"Newton solver did not converge with minimum time step of "
520 <<
dt <<
" seconds. Continuing with unconverged solution!\n"
526 "Time integration did not succeed with the minumum time step size of "
537 std::cout <<
"Newton solver did not converge with "
538 <<
"dt=" <<
dt <<
" seconds. Retrying with time step of "
539 <<
nextDt <<
" seconds\n" << std::flush;
544 "Newton solver didn't converge after "
545 +std::to_string(
maxFails)+
" time-step divisions. dt="
546 +std::to_string(
double(
simulator().timeStepSize()));
554 {
return Parameters::Get<Parameters::MinTimeStepSize<Scalar>>(); }
561 {
return Parameters::Get<Parameters::MaxTimeStepDivisions>(); }
569 {
return Parameters::Get<Parameters::ContinueOnConvergenceError>(); }
575 { nextTimeStepSize_ =
dt; }
584 if (nextTimeStepSize_ > 0.0)
585 return nextTimeStepSize_;
609 return simulator().timeStepIndex() > 0 &&
630 {
model().advanceTimeLevel(); }
646 {
return gridView_; }
653 {
return boundingBoxMin_; }
660 {
return boundingBoxMax_; }
666 {
return vertexMapper_; }
672 {
return elementMapper_; }
678 {
return simulator_; }
684 {
return simulator_; }
690 {
return simulator_.model(); }
696 {
return simulator_.model(); }
702 {
return model().newtonMethod(); }
708 {
return model().newtonMethod(); }
745 template <
class Restarter>
748 if (enableVtkOutput_())
762 template <
class Restarter>
765 if (enableVtkOutput_())
777 if (!enableVtkOutput_())
781 std::cout <<
"Writing visualization results for the current time step.\n"
788 model().prepareOutputFields();
789 model().appendOutputFields(*defaultVtkWriter_);
799 {
return defaultVtkWriter_; }
802 Scalar nextTimeStepSize_;
804 bool enableVtkOutput_()
const
805 {
return Parameters::Get<Parameters::EnableVtkOutput>(); }
809 Implementation& asImp_()
810 {
return *
static_cast<Implementation *
>(
this); }
813 const Implementation& asImp_()
const
814 {
return *
static_cast<const Implementation *
>(
this); }
817 const GridView gridView_;
818 ElementMapper elementMapper_;
819 VertexMapper vertexMapper_;
820 GlobalPosition boundingBoxMin_;
821 GlobalPosition boundingBoxMax_;
824 Simulator& simulator_;
825 mutable VtkMultiWriter *defaultVtkWriter_;
Definition restrictprolong.hh:142
Base class for all problems which use a finite volume spatial discretization.
Definition fvbaseproblem.hh:67
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:256
std::string name() const
The problem name.
Definition fvbaseproblem.hh:639
const Simulator & simulator() const
Returns Simulator object used by the simulation.
Definition fvbaseproblem.hh:683
void writeOutput(bool verbose=true)
Write the relevant secondary variables of the current solution into an VTK output file.
Definition fvbaseproblem.hh:775
void boundary(BoundaryRateVector &, const Context &, unsigned, unsigned) const
Evaluate the boundary conditions for a boundary segment.
Definition fvbaseproblem.hh:308
void finalize()
Called after the simulation has been run sucessfully.
Definition fvbaseproblem.hh:446
Simulator & simulator()
Returns Simulator object used by the simulation.
Definition fvbaseproblem.hh:677
const Model & model() const
Returns numerical model used for the problem.
Definition fvbaseproblem.hh:695
static void registerParameters()
Registers all available parameters for the problem and the model.
Definition fvbaseproblem.hh:170
void source(RateVector &, const Context &, unsigned, unsigned) const
Evaluate the source term for all phases within a given sub-control-volume.
Definition fvbaseproblem.hh:344
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:195
unsigned maxTimeIntegrationFailures() const
Returns the maximum number of subsequent failures for the time integration before giving up.
Definition fvbaseproblem.hh:560
unsigned markForGridAdaptation()
Mark grid cells for refinement or coarsening.
Definition fvbaseproblem.hh:727
void serialize(Restarter &res)
This method writes the complete state of the problem to the harddisk.
Definition fvbaseproblem.hh:746
void beginIteration()
Called by the simulator before each Newton-Raphson iteration.
Definition fvbaseproblem.hh:413
void endTimeStep()
Called by the simulator after each time integration.
Definition fvbaseproblem.hh:428
FvBaseProblem(Simulator &simulator)
Definition fvbaseproblem.hh:117
const GlobalPosition & boundingBoxMax() const
The coordinate of the corner of the GridView's bounding box with the largest values.
Definition fvbaseproblem.hh:659
RestrictProlongOperator restrictProlongOperator()
return restriction and prolongation operator
Definition fvbaseproblem.hh:715
void timeIntegration()
Called by Opm::Simulator in order to do a time integration on the model.
Definition fvbaseproblem.hh:503
Model & model()
Returns numerical model used for the problem.
Definition fvbaseproblem.hh:689
void endEpisode()
Called when the end of an simulation episode is reached.
Definition fvbaseproblem.hh:436
void beginTimeStep()
Called by the simulator before each time integration.
Definition fvbaseproblem.hh:407
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:217
void deserialize(Restarter &res)
This method restores the complete state of the problem from disk.
Definition fvbaseproblem.hh:763
void prefetch(const Element &) const
Allows to improve the performance by prefetching all data which is associated with a given element.
Definition fvbaseproblem.hh:281
Scalar extrusionFactor(const Context &, unsigned, unsigned) const
Return how much the domain is extruded at a given sub-control volume.
Definition fvbaseproblem.hh:383
bool shouldWriteRestartFile() const
Returns true if a restart file should be written to disk.
Definition fvbaseproblem.hh:607
Scalar nextTimeStepSize() const
Called by Opm::Simulator whenever a solution for a time step has been computed and the simulation tim...
Definition fvbaseproblem.hh:582
Scalar minTimeStepSize() const
Returns the minimum allowable size of a time step.
Definition fvbaseproblem.hh:553
void beginEpisode()
Called at the beginning of an simulation episode.
Definition fvbaseproblem.hh:401
void initial(PrimaryVariables &, const Context &, unsigned, unsigned) const
Evaluate the initial value for a control volume.
Definition fvbaseproblem.hh:361
std::string outputDir() const
Determine the directory for simulation output.
Definition fvbaseproblem.hh:206
const VertexMapper & vertexMapper() const
Returns the mapper for vertices to indices.
Definition fvbaseproblem.hh:665
void gridChanged()
Handle changes of the grid.
Definition fvbaseproblem.hh:289
const GlobalPosition & boundingBoxMin() const
The coordinate of the corner of the GridView's bounding box with the smallest values.
Definition fvbaseproblem.hh:652
static std::string briefDescription()
Returns a human readable description of the problem for the help message.
Definition fvbaseproblem.hh:235
void endIteration()
Called by the simulator after each Newton-Raphson update.
Definition fvbaseproblem.hh:419
void setNextTimeStepSize(Scalar dt)
Impose the next time step size to be used externally.
Definition fvbaseproblem.hh:574
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:568
void constraints(Constraints &, const Context &, unsigned, unsigned) const
Evaluate the constraints for a control volume.
Definition fvbaseproblem.hh:325
void finishInit()
Called by the Opm::Simulator in order to initialize the problem.
Definition fvbaseproblem.hh:274
const NewtonMethod & newtonMethod() const
Returns object which implements the Newton method.
Definition fvbaseproblem.hh:707
const GridView & gridView() const
The GridView which used by the problem.
Definition fvbaseproblem.hh:645
NewtonMethod & newtonMethod()
Returns object which implements the Newton method.
Definition fvbaseproblem.hh:701
bool shouldWriteOutput() const
Returns true if the current solution should be written to disk (i.e.
Definition fvbaseproblem.hh:621
void initialSolutionApplied()
Callback used by the model to indicate that the initial solution has been determined for all degrees ...
Definition fvbaseproblem.hh:395
void advanceTimeLevel()
Called by the simulator after everything which can be done about the current time step is finished an...
Definition fvbaseproblem.hh:629
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:798
const ElementMapper & elementMapper() const
Returns the mapper for elements to indices.
Definition fvbaseproblem.hh:671
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:66
void endWrite(bool onlyDiscard=false) override
Finalizes the current writer.
Definition vtkmultiwriter.hh:384
void deserialize(Restarter &res)
Read the multi-writer's state from a restart file.
Definition vtkmultiwriter.hh:437
void gridChanged()
Updates the internal data structures after mesh refinement.
Definition vtkmultiwriter.hh:165
void beginWrite(double t) override
Called whenever a new time step must be written.
Definition vtkmultiwriter.hh:174
void serialize(Restarter &res)
Write the multi-writer's state to a restart file.
Definition vtkmultiwriter.hh:403
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.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
std::string humanReadableTime(double timeInSeconds, bool isAmendment)
Given a time step size in seconds, return it in a format which is more easily parsable by humans.
Definition simulatorutils.cpp:45
std::string simulatorOutputDir()
Determine and check the configured directory for simulation output.
Definition simulatorutils.cpp:119
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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:235
Load or save a state of a problem to/from the harddisk.
Specify the maximum size of a time integration [s].
Definition fvbaseparameters.hh:106
Simplifies writing multi-file VTK datasets.