All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
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 "fvbaseproperties.hh"
32 
34 #include <ewoms/io/restart.hh>
35 #include <ewoms/disc/common/restrictprolong.hh>
36 
37 #include <opm/common/Unused.hpp>
38 #include <opm/common/ErrorMacros.hpp>
39 #include <opm/common/Exceptions.hpp>
40 
41 #include <dune/common/fvector.hh>
42 
43 #include <iostream>
44 #include <limits>
45 #include <string>
46 
47 namespace Ewoms {
48 
58 template<class TypeTag>
60 {
61 private:
62  typedef typename GET_PROP_TYPE(TypeTag, Problem) Implementation;
63  typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
64 
65  static const int vtkOutputFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
67 
68  typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
69  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
70  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
71  typedef typename GET_PROP_TYPE(TypeTag, ThreadManager) ThreadManager;
72  typedef typename GET_PROP_TYPE(TypeTag, NewtonMethod) NewtonMethod;
73 
74  typedef typename GET_PROP_TYPE(TypeTag, VertexMapper) VertexMapper;
75  typedef typename GET_PROP_TYPE(TypeTag, ElementMapper) ElementMapper;
76 
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;
81 
82  enum {
83  dim = GridView::dimension,
84  dimWorld = GridView::dimensionworld
85  };
86 
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;
90 
91  typedef typename GridView::Grid::ctype CoordScalar;
92  typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
93 
94 public:
95  // the default restriction and prolongation for adaptation is simply an empty one
97 
98 private:
99  // copying a problem is not a good idea
100  FvBaseProblem(const FvBaseProblem& ) = delete;
101 
102 public:
111  : gridView_(simulator.gridView())
112  , elementMapper_(gridView_)
113  , vertexMapper_(gridView_)
114  , boundingBoxMin_(std::numeric_limits<double>::max())
115  , boundingBoxMax_(-std::numeric_limits<double>::max())
116  , simulator_(simulator)
117  , defaultVtkWriter_(0)
118  {
119  // calculate the bounding box of the local partition of the grid view
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]);
126  }
127  }
128 
129  // communicate to get the bounding box of the whole domain
130  for (unsigned i = 0; i < dim; ++i) {
131  boundingBoxMin_[i] = gridView_.comm().min(boundingBoxMin_[i]);
132  boundingBoxMax_[i] = gridView_.comm().max(boundingBoxMax_[i]);
133  }
134 
135  if (enableVtkOutput_())
136  defaultVtkWriter_ = new VtkMultiWriter(gridView_, asImp_().name());
137  }
138 
139  ~FvBaseProblem()
140  { delete defaultVtkWriter_; }
141 
146  static void registerParameters()
147  {
148  Model::registerParameters();
149  EWOMS_REGISTER_PARAM(TypeTag, Scalar, MaxTimeStepSize,
150  "The maximum size to which all time steps are limited to [s]");
151  EWOMS_REGISTER_PARAM(TypeTag, Scalar, MinTimeStepSize,
152  "The minimum size to which all time steps are limited to [s]");
153  EWOMS_REGISTER_PARAM(TypeTag, unsigned, MaxTimeStepDivisions,
154  "The maximum number of divisions by two of the timestep size "
155  "before the simulation bails out");
156  }
157 
163  void finishInit()
164  { }
165 
170  void prefetch(const Element& elem OPM_UNUSED) const
171  {
172  // do nothing by default
173  }
174 
178  void gridChanged()
179  {
180  elementMapper_.update();
181  vertexMapper_.update();
182 
183  if (enableVtkOutput_())
184  defaultVtkWriter_->gridChanged();
185  }
186 
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"); }
202 
213  template <class Context>
214  void constraints(Constraints& constraints OPM_UNUSED,
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"); }
219 
232  template <class Context>
233  void source(RateVector& rate OPM_UNUSED,
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"); }
238 
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"); }
255 
271  template <class Context>
272  Scalar extrusionFactor(const Context& context OPM_UNUSED,
273  unsigned spaceIdx OPM_UNUSED,
274  unsigned timeIdx OPM_UNUSED) const
275  { return asImp_().extrusionFactor(); }
276 
277  Scalar extrusionFactor() const
278  { return 1.0; }
279 
285  {}
286 
291  { }
292 
297  { }
298 
303  { }
304 
309  { }
310 
317  void endTimeStep()
318  { }
319 
325  void endEpisode()
326  {
327  std::cerr << "The end of episode " << simulator().episodeIndex() + 1 << " has been "
328  << "reached, but the problem does not override the endEpisode() method. "
329  << "Doing nothing!\n";
330  }
331 
335  void finalize()
336  {
337  const auto& executionTimer = simulator().executionTimer();
338 
339  Scalar executionTime = executionTimer.realTimeElapsed();
340  Scalar setupTime = simulator().setupTimer().realTimeElapsed();
341  Scalar prePostProcessTime = simulator().prePostProcessTimer().realTimeElapsed();
342  Scalar localCpuTime = executionTimer.cpuTimeElapsed();
343  Scalar globalCpuTime = executionTimer.globalCpuTimeElapsed();
344  Scalar writeTime = simulator().writeTimer().realTimeElapsed();
345  Scalar linearizeTime = simulator().linearizeTimer().realTimeElapsed();
346  Scalar solveTime = simulator().solveTimer().realTimeElapsed();
347  Scalar updateTime = simulator().updateTimer().realTimeElapsed();
348  unsigned numProcesses = static_cast<unsigned>(this->gridView().comm().size());
349  unsigned threadsPerProcess = ThreadManager::maxThreads();
350  if (gridView().comm().rank() == 0) {
351  std::cout << std::setprecision(3)
352  << "Simulation of problem '" << asImp_().name() << "' finished.\n"
353  << "\n"
354  << "------------------------ Timing receipt ------------------------\n"
355  << "Setup time: " << setupTime << " seconds" << Simulator::humanReadableTime(setupTime)
356  << ", " << setupTime/(executionTime + setupTime)*100 << "%\n"
357  << "Simulation time: " << executionTime << " seconds" << Simulator::humanReadableTime(executionTime)
358  << ", " << executionTime/(executionTime + setupTime)*100 << "%\n"
359  << " Linearization time: " << linearizeTime << " seconds" << Simulator::humanReadableTime(linearizeTime)
360  << ", " << linearizeTime/executionTime*100 << "%\n"
361  << " Linear solve time: " << solveTime << " seconds" << Simulator::humanReadableTime(solveTime)
362  << ", " << solveTime/executionTime*100 << "%\n"
363  << " Newton update time: " << updateTime << " seconds" << Simulator::humanReadableTime(updateTime)
364  << ", " << updateTime/executionTime*100 << "%\n"
365  << " Pre/postprocess time: " << prePostProcessTime << " seconds" << Simulator::humanReadableTime(prePostProcessTime)
366  << ", " << prePostProcessTime/executionTime*100 << "%\n"
367  << " Output write time: " << writeTime << " seconds" << Simulator::humanReadableTime(writeTime)
368  << ", " << writeTime/executionTime*100 << "%\n"
369  << "First process' simulation CPU time: " << localCpuTime << " seconds" << Simulator::humanReadableTime(localCpuTime) << "\n"
370  << "Number of processes: " << numProcesses << "\n"
371  << "Threads per processes: " << threadsPerProcess << "\n"
372  << "Total CPU time: " << globalCpuTime << " seconds" << Simulator::humanReadableTime(globalCpuTime) << "\n"
373  << "\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
377  << "%\n"
378  << "\n"
379  << "Our simulation hours are 24/7. Thank you for choosing us.\n"
380  << "----------------------------------------------------------------\n"
381  << "\n"
382  << std::flush;
383  }
384  }
385 
391  {
392  unsigned maxFails = EWOMS_GET_PARAM(TypeTag, unsigned, MaxTimeStepDivisions);
393  Scalar minTimeStepSize = EWOMS_GET_PARAM(TypeTag, Scalar, MinTimeStepSize);
394 
395  // if the time step size of the simulator is smaller than
396  // the specified minimum size and we're not going to finish
397  // the simulation or an episode, try with the minimum size.
398  if (simulator().timeStepSize() < minTimeStepSize &&
399  !simulator().episodeWillBeOver() &&
400  !simulator().willBeFinished())
401  {
402  simulator().setTimeStepSize(minTimeStepSize);
403  }
404 
405  for (unsigned i = 0; i < maxFails; ++i) {
406  bool converged = model().update();
407  if (converged)
408  return;
409 
410  Scalar dt = simulator().timeStepSize();
411  Scalar nextDt = dt / 2;
412  if (nextDt < minTimeStepSize)
413  break; // give up: we can't make the time step smaller anymore!
414  simulator().setTimeStepSize(nextDt);
415 
416  // update failed
417  if (gridView().comm().rank() == 0)
418  std::cout << "Newton solver did not converge with "
419  << "dt=" << dt << " seconds. Retrying with time step of "
420  << nextDt << " seconds\n" << std::flush;
421  }
422 
423  OPM_THROW(std::runtime_error,
424  "Newton solver didn't converge after "
425  << maxFails << " time-step divisions. dt="
426  << simulator().timeStepSize());
427  }
428 
435  {
436  Scalar dtNext = std::min(EWOMS_GET_PARAM(TypeTag, Scalar, MaxTimeStepSize),
437  newtonMethod().suggestTimeStepSize(simulator().timeStepSize()));
438 
439  if (dtNext < simulator().maxTimeStepSize()
440  && simulator().maxTimeStepSize() < dtNext*2)
441  {
442  dtNext = simulator().maxTimeStepSize()/2 * 1.01;
443  }
444 
445  return dtNext;
446  }
447 
457  {
458  return simulator().timeStepIndex() > 0 &&
459  (simulator().timeStepIndex() % 10 == 0);
460  }
461 
470  bool shouldWriteOutput() const
471  { return true; }
472 
479  { model().advanceTimeLevel(); }
480 
488  std::string name() const
489  { return "sim"; }
490 
494  const GridView& gridView() const
495  { return gridView_; }
496 
501  const GlobalPosition& boundingBoxMin() const
502  { return boundingBoxMin_; }
503 
508  const GlobalPosition& boundingBoxMax() const
509  { return boundingBoxMax_; }
510 
514  const VertexMapper& vertexMapper() const
515  { return vertexMapper_; }
516 
520  const ElementMapper& elementMapper() const
521  { return elementMapper_; }
522 
527  { return simulator_; }
528 
532  const Simulator& simulator() const
533  { return simulator_; }
534 
538  Model& model()
539  { return simulator_.model(); }
540 
544  const Model& model() const
545  { return simulator_.model(); }
546 
551  { return model().newtonMethod(); }
552 
556  const NewtonMethod& newtonMethod() const
557  { return model().newtonMethod(); }
558  // \}
559 
565  {
566  return RestrictProlongOperator();
567  }
568 
577  {
578  return 0;
579  }
580 
594  template <class Restarter>
595  void serialize(Restarter& res)
596  {
597  if (enableVtkOutput_())
598  defaultVtkWriter_->serialize(res);
599  }
600 
611  template <class Restarter>
612  void deserialize(Restarter& res)
613  {
614  if (enableVtkOutput_())
615  defaultVtkWriter_->deserialize(res);
616  }
617 
624  void writeOutput(bool verbose = true)
625  {
626  if (verbose && gridView().comm().rank() == 0)
627  std::cout << "Writing visualization results for the current time step.\n"
628  << std::flush;
629 
630  // calculate the time _after_ the time was updated
631  Scalar t = simulator().time() + simulator().timeStepSize();
632 
633  if (enableVtkOutput_())
634  defaultVtkWriter_->beginWrite(t);
635 
636  model().prepareOutputFields();
637 
638  if (enableVtkOutput_()) {
639  model().appendOutputFields(*defaultVtkWriter_);
640  defaultVtkWriter_->endWrite();
641  }
642  }
643 
649  { return defaultVtkWriter_; }
650 
651 private:
652  bool enableVtkOutput_() const
653  { return EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput); }
654 
656  Implementation& asImp_()
657  { return *static_cast<Implementation *>(this); }
658 
660  const Implementation& asImp_() const
661  { return *static_cast<const Implementation *>(this); }
662 
663  // Grid management stuff
664  const GridView gridView_;
665  ElementMapper elementMapper_;
666  VertexMapper vertexMapper_;
667  GlobalPosition boundingBoxMin_;
668  GlobalPosition boundingBoxMax_;
669 
670  // Attributes required for the actual simulation
671  Simulator& simulator_;
672  mutable VtkMultiWriter *defaultVtkWriter_;
673 };
674 
675 } // namespace Ewoms
676 
677 #endif
const Model & model() const
Returns numerical model used for the problem.
Definition: fvbaseproblem.hh:544
void finishInit()
Called by the Ewoms::Simulator in order to initialize the problem.
Definition: fvbaseproblem.hh:163
int episodeIndex() const
Returns the index of the current episode.
Definition: simulator.hh:452
void endEpisode()
Called when the end of an simulation episode is reached.
Definition: fvbaseproblem.hh:325
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 serialize(Restarter &res)
Write the multi-writer&#39;s state to a restart file.
Definition: vtkmultiwriter.hh:378
std::string name() const
The problem name.
Definition: fvbaseproblem.hh:488
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
Simplifies multi-threaded capabilities.
Definition: threadmanager.hh:52
const GlobalPosition & boundingBoxMin() const
The coordinate of the corner of the GridView&#39;s bounding box with the smallest values.
Definition: fvbaseproblem.hh:501
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
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:189
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
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:63
const ElementMapper & elementMapper() const
Returns the mapper for elements to indices.
Definition: fvbaseproblem.hh:520
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
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
const NewtonMethod & newtonMethod() const
Returns object which implements the Newton method.
Definition: fvbaseproblem.hh:556
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 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
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
void deserialize(Restarter &res)
Read the multi-writer&#39;s state from a restart file.
Definition: vtkmultiwriter.hh:412
const Simulator & simulator() const
Returns Simulator object used by the simulation.
Definition: fvbaseproblem.hh:532
Declare the properties used by the infrastructure code of the finite volume discretizations.
FvBaseProblem(Simulator &simulator)
Definition: fvbaseproblem.hh:110
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 endTimeStep()
Called by the simulator after each time integration.
Definition: fvbaseproblem.hh:317
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
void deserialize(Restarter &res)
This method restores the complete state of the problem from disk.
Definition: fvbaseproblem.hh:612
const GridView & gridView() const
The GridView which used by the problem.
Definition: fvbaseproblem.hh:494
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
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
const Ewoms::Timer & solveTimer() const
Returns a reference to the timer object which measures the time needed by the solver.
Definition: simulator.hh:304
Definition: restrictprolong.hh:140
const VertexMapper & vertexMapper() const
Returns the mapper for vertices to indices.
Definition: fvbaseproblem.hh:514
void timeIntegration()
Called by Ewoms::Simulator in order to do a time integration on the model.
Definition: fvbaseproblem.hh:390
bool shouldWriteOutput() const
Returns true if the current solution should be written to disk (i.e.
Definition: fvbaseproblem.hh:470
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
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
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
Scalar timeStepSize() const
Returns the time step length so that we don&#39;t miss the beginning of the next episode or cross the en...
Definition: simulator.hh:347
double realTimeElapsed() const
Return the real time [s] elapsed during the periods the timer was active since the last reset...
Definition: timer.hh:121
void finalize()
Called after the simulation has been run sucessfully.
Definition: fvbaseproblem.hh:335
bool shouldWriteRestartFile() const
Returns true if a restart file should be written to disk.
Definition: fvbaseproblem.hh:456
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
Simplifies writing multi-file VTK datasets.
void gridChanged()
Updates the internal data structures after mesh refinement.
Definition: vtkmultiwriter.hh:126
void writeOutput(bool verbose=true)
Write the relevant secondary variables of the current solution into an VTK output file...
Definition: fvbaseproblem.hh:624
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
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
Scalar nextTimeStepSize()
Called by Ewoms::Simulator whenever a solution for a time step has been computed and the simulation t...
Definition: fvbaseproblem.hh:434
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
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
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:75
const GlobalPosition & boundingBoxMax() const
The coordinate of the corner of the GridView&#39;s bounding box with the largest values.
Definition: fvbaseproblem.hh:508
Model & model()
Returns numerical model used for the problem.
Definition: fvbaseproblem.hh:538
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
int timeStepIndex() const
Returns number of time steps which have been executed since the beginning of the simulation.
Definition: simulator.hh:360
Scalar time() const
Return the number of seconds of simulated time which have elapsed since the start time...
Definition: simulator.hh:254
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
NewtonMethod & newtonMethod()
Returns object which implements the Newton method.
Definition: fvbaseproblem.hh:550