fvbaseelementcontext.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_ELEMENT_CONTEXT_HH
29 #define EWOMS_FV_BASE_ELEMENT_CONTEXT_HH
30 
31 #include "fvbaseproperties.hh"
32 
34 
35 #include <opm/common/Unused.hpp>
36 #include <opm/common/ErrorMacros.hpp>
37 #include <opm/common/Exceptions.hpp>
38 
39 #include <dune/common/fvector.hh>
40 
41 #include <vector>
42 
43 namespace Ewoms {
44 
51 template<class TypeTag>
53 {
54  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) Implementation;
55 
56  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
57  typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
58  typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
59  typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
60 
61  // the history size of the time discretization in number of steps
62  enum { timeDiscHistorySize = GET_PROP_VALUE(TypeTag, TimeDiscHistorySize) };
63 
64  struct DofStore_ {
65  IntensiveQuantities intensiveQuantities[timeDiscHistorySize];
66  PrimaryVariables priVars[timeDiscHistorySize];
67  const IntensiveQuantities *thermodynamicHint[timeDiscHistorySize];
68  };
69  typedef std::vector<DofStore_> DofVarsVector;
70  typedef std::vector<ExtensiveQuantities> ExtensiveQuantitiesVector;
71 
72  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
73  typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
74  typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
75  typedef typename GET_PROP_TYPE(TypeTag, Stencil) Stencil;
76  typedef typename GET_PROP_TYPE(TypeTag, GradientCalculator) GradientCalculator;
77  typedef typename GET_PROP_TYPE(TypeTag, SolutionVector) SolutionVector;
78 
79  typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
80  typedef typename GridView::template Codim<0>::Entity Element;
81 
82  static const unsigned dim = GridView::dimension;
83  static const unsigned numEq = GET_PROP_VALUE(TypeTag, NumEq);
84 
85  typedef typename GridView::ctype CoordScalar;
86  typedef Dune::FieldVector<CoordScalar, dim> GlobalPosition;
87 
88  // we don't allow copies of element contexts!
90 
91 public:
96  : gridView_(simulator.gridView())
97  , stencil_(gridView_, simulator.model().dofMapper() )
98  {
99  // remember the simulator object
100  simulatorPtr_ = &simulator;
101  enableStorageCache_ = EWOMS_GET_PARAM(TypeTag, bool, EnableStorageCache);
102  stashedDofIdx_ = -1;
103  focusDofIdx_ = -1;
104  }
105 
106  static void *operator new(size_t size) {
107  return Ewoms::aligned_alloc(alignof(FvBaseElementContext), size);
108  }
109 
110  static void operator delete(void *ptr) {
111  Ewoms::aligned_free(ptr);
112  }
113 
121  void updateAll(const Element& elem)
122  {
123  asImp_().updateStencil(elem);
124  asImp_().updateAllIntensiveQuantities();
125  asImp_().updateAllExtensiveQuantities();
126  }
127 
134  void updateStencil(const Element& elem)
135  {
136  // remember the current element
137  elemPtr_ = &elem;
138 
139  // update the stencil. the center gradients are quite expensive to calculate and
140  // most models don't need them, so that we only do this if the model explicitly
141  // enables them
142  stencil_.update(elem);
143 
144  // resize the arrays containing the flux and the volume variables
145  dofVars_.resize(stencil_.numDof());
146  extensiveQuantities_.resize(stencil_.numInteriorFaces());
147  }
148 
155  void updatePrimaryStencil(const Element& elem)
156  {
157  // remember the current element
158  elemPtr_ = &elem;
159 
160  // update the finite element geometry
161  stencil_.updatePrimaryTopology(elem);
162 
163  dofVars_.resize(stencil_.numPrimaryDof());
164  }
165 
172  void updateStencilTopology(const Element& elem)
173  {
174  // remember the current element
175  elemPtr_ = &elem;
176 
177  // update the finite element geometry
178  stencil_.updateTopology(elem);
179  }
180 
186  {
187  if (!enableStorageCache_) {
188  // if the storage cache is disabled, we need to calculate the storage term
189  // from scratch, i.e. we need the intensive quantities of all of the history.
190  for (unsigned timeIdx = 0; timeIdx < timeDiscHistorySize; ++ timeIdx)
191  asImp_().updateIntensiveQuantities(timeIdx);
192  }
193  else
194  // if the storage cache is enabled, we only need to recalculate the storage
195  // term for the most recent point of history (i.e., for the current iterative
196  // solution)
197  asImp_().updateIntensiveQuantities(/*timeIdx=*/0);
198  }
199 
206  void updateIntensiveQuantities(unsigned timeIdx)
207  { updateIntensiveQuantities_(timeIdx, numDof(timeIdx)); }
208 
215  void updatePrimaryIntensiveQuantities(unsigned timeIdx)
216  { updateIntensiveQuantities_(timeIdx, numPrimaryDof(timeIdx)); }
217 
228  void updateIntensiveQuantities(const PrimaryVariables& priVars, unsigned dofIdx, unsigned timeIdx)
229  { asImp_().updateSingleIntQuants_(priVars, dofIdx, timeIdx); }
230 
236  { asImp_().updateExtensiveQuantities(/*timeIdx=*/0); }
237 
245  void updateExtensiveQuantities(unsigned timeIdx)
246  {
247  gradientCalculator_.prepare(/*context=*/asImp_(), timeIdx);
248 
249  for (unsigned fluxIdx = 0; fluxIdx < numInteriorFaces(timeIdx); fluxIdx++) {
250  extensiveQuantities_[fluxIdx].update(/*context=*/asImp_(),
251  /*localIndex=*/fluxIdx,
252  timeIdx);
253  }
254  }
255 
263  void setFocusDofIndex(unsigned dofIdx)
264  { focusDofIdx_ = dofIdx; }
265 
271  unsigned focusDofIndex() const
272  { return focusDofIdx_; }
273 
277  const Simulator& simulator() const
278  { return *simulatorPtr_; }
279 
283  const Problem& problem() const
284  { return simulatorPtr_->problem(); }
285 
289  const Model& model() const
290  { return simulatorPtr_->model(); }
291 
295  const GridView& gridView() const
296  { return gridView_; }
297 
301  const Element& element() const
302  { return *elemPtr_; }
303 
307  size_t numDof(unsigned timeIdx) const
308  { return stencil(timeIdx).numDof(); }
309 
313  size_t numPrimaryDof(unsigned timeIdx) const
314  { return stencil(timeIdx).numPrimaryDof(); }
315 
320  size_t numInteriorFaces(unsigned timeIdx) const
321  { return stencil(timeIdx).numInteriorFaces(); }
322 
327  size_t numBoundaryFaces(unsigned timeIdx) const
328  { return stencil(timeIdx).numBoundaryFaces(); }
329 
336  const Stencil& stencil(unsigned timeIdx OPM_UNUSED) const
337  { return stencil_; }
338 
347  const GlobalPosition& pos(unsigned dofIdx, unsigned timeIdx OPM_UNUSED) const
348  { return stencil_.subControlVolume(dofIdx).globalPos(); }
349 
358  unsigned globalSpaceIndex(unsigned dofIdx, unsigned timeIdx) const
359  { return stencil(timeIdx).globalSpaceIndex(dofIdx); }
360 
361 
371  Scalar dofVolume(unsigned dofIdx, unsigned timeIdx) const
372  { return stencil(timeIdx).subControlVolume(dofIdx).volume(); }
373 
384  Scalar dofTotalVolume(unsigned dofIdx, unsigned timeIdx) const
385  { return model().dofTotalVolume(globalSpaceIndex(dofIdx, timeIdx)); }
386 
391  bool onBoundary() const
392  { return element().hasBoundaryIntersections(); }
393 
404  const IntensiveQuantities& intensiveQuantities(unsigned dofIdx, unsigned timeIdx) const
405  {
406 #ifndef NDEBUG
407  assert(0 <= dofIdx && dofIdx < numDof(timeIdx));
408 
409  if (enableStorageCache_ && timeIdx != 0)
410  OPM_THROW(std::logic_error,
411  "If caching of the storage term is enabled, only the intensive quantities "
412  "for the most-recent substep (i.e. time index 0) are available!");
413 #endif
414 
415  return dofVars_[dofIdx].intensiveQuantities[timeIdx];
416  }
417 
426  const IntensiveQuantities *thermodynamicHint(unsigned dofIdx, unsigned timeIdx) const
427  {
428  assert(0 <= dofIdx && dofIdx < numDof(timeIdx));
429  return dofVars_[dofIdx].thermodynamicHint[timeIdx];
430  }
434  IntensiveQuantities& intensiveQuantities(unsigned dofIdx, unsigned timeIdx)
435  {
436  assert(0 <= dofIdx && dofIdx < numDof(timeIdx));
437  return dofVars_[dofIdx].intensiveQuantities[timeIdx];
438  }
439 
448  PrimaryVariables& primaryVars(unsigned dofIdx, unsigned timeIdx)
449  {
450  assert(0 <= dofIdx && dofIdx < numDof(timeIdx));
451  return dofVars_[dofIdx].priVars[timeIdx];
452  }
456  const PrimaryVariables& primaryVars(unsigned dofIdx, unsigned timeIdx) const
457  {
458  assert(0 <= dofIdx && dofIdx < numDof(timeIdx));
459  return dofVars_[dofIdx].priVars[timeIdx];
460  }
461 
469  { return stashedDofIdx_ != -1; }
470 
477  int stashedDofIdx() const
478  { return stashedDofIdx_; }
479 
485  void stashIntensiveQuantities(unsigned dofIdx)
486  {
487  assert(0 <= dofIdx && dofIdx < numDof(/*timeIdx=*/0));
488 
489  intensiveQuantitiesStashed_ = dofVars_[dofIdx].intensiveQuantities[/*timeIdx=*/0];
490  priVarsStashed_ = dofVars_[dofIdx].priVars[/*timeIdx=*/0];
491  stashedDofIdx_ = static_cast<int>(dofIdx);
492  }
493 
499  void restoreIntensiveQuantities(unsigned dofIdx)
500  {
501  dofVars_[dofIdx].priVars[/*timeIdx=*/0] = priVarsStashed_;
502  dofVars_[dofIdx].intensiveQuantities[/*timeIdx=*/0] = intensiveQuantitiesStashed_;
503  stashedDofIdx_ = -1;
504  }
505 
510  const GradientCalculator& gradientCalculator() const
511  { return gradientCalculator_; }
512 
521  const ExtensiveQuantities& extensiveQuantities(unsigned fluxIdx, unsigned timeIdx OPM_UNUSED) const
522  { return extensiveQuantities_[fluxIdx]; }
523 
530  bool enableStorageCache() const
531  { return enableStorageCache_; }
532 
536  void setEnableStorageCache(bool yesno)
537  { enableStorageCache_ = yesno; }
538 
539 private:
540  Implementation& asImp_()
541  { return *static_cast<Implementation*>(this); }
542 
543  const Implementation& asImp_() const
544  { return *static_cast<const Implementation*>(this); }
545 
546 protected:
552  void updateIntensiveQuantities_(unsigned timeIdx, size_t numDof)
553  {
554  // update the intensive quantities for the whole history
555  const SolutionVector& globalSol = model().solution(timeIdx);
556 
557  // update the non-gradient quantities
558  for (unsigned dofIdx = 0; dofIdx < numDof; dofIdx++) {
559  unsigned globalIdx = globalSpaceIndex(dofIdx, timeIdx);
560  const PrimaryVariables& dofSol = globalSol[globalIdx];
561  dofVars_[dofIdx].priVars[timeIdx] = dofSol;
562 
563  dofVars_[dofIdx].thermodynamicHint[timeIdx] =
564  model().thermodynamicHint(globalIdx, timeIdx);
565 
566  const auto *cachedIntQuants = model().cachedIntensiveQuantities(globalIdx, timeIdx);
567  if (cachedIntQuants) {
568  dofVars_[dofIdx].intensiveQuantities[timeIdx] = *cachedIntQuants;
569  }
570  else {
571  updateSingleIntQuants_(dofSol, dofIdx, timeIdx);
572  model().updateCachedIntensiveQuantities(dofVars_[dofIdx].intensiveQuantities[timeIdx],
573  globalIdx,
574  timeIdx);
575  }
576  }
577  }
578 
579  void updateSingleIntQuants_(const PrimaryVariables& priVars, unsigned dofIdx, unsigned timeIdx)
580  {
581 #ifndef NDEBUG
582  if (enableStorageCache_ && timeIdx != 0)
583  OPM_THROW(std::logic_error,
584  "If caching of the storage term is enabled, only the intensive quantities "
585  "for the most-recent substep (i.e. time index 0) are available!");
586 #endif
587 
588  dofVars_[dofIdx].priVars[timeIdx] = priVars;
589  dofVars_[dofIdx].intensiveQuantities[timeIdx].update(/*context=*/asImp_(), dofIdx, timeIdx);
590  }
591 
592  IntensiveQuantities intensiveQuantitiesStashed_;
593  PrimaryVariables priVarsStashed_;
594 
595  GradientCalculator gradientCalculator_;
596 
597  std::vector<DofStore_, Ewoms::aligned_allocator<DofStore_, alignof(DofStore_)> > dofVars_;
598  std::vector<ExtensiveQuantities, Ewoms::aligned_allocator<ExtensiveQuantities, alignof(ExtensiveQuantities)> > extensiveQuantities_;
599 
600  const Simulator *simulatorPtr_;
601  const Element *elemPtr_;
602  const GridView gridView_;
603  Stencil stencil_;
604 
605  int stashedDofIdx_;
606  int focusDofIdx_;
607  bool enableStorageCache_;
608 };
609 
610 } // namespace Ewoms
611 
612 #endif
bool enableStorageCache() const
Returns true iff the cache for the storage term ought to be used for this context.
Definition: fvbaseelementcontext.hh:530
This is a stand-alone version of boost::alignment::aligned_allocator from Boost 1...
void updateStencilTopology(const Element &elem)
Update the topological part of the stencil, but nothing else.
Definition: fvbaseelementcontext.hh:172
size_t numDof(unsigned timeIdx) const
Return the number of sub-control volumes of the current element.
Definition: fvbaseelementcontext.hh:307
void updateAllIntensiveQuantities()
Compute the intensive quantities of all sub-control volumes of the current element for all time indic...
Definition: fvbaseelementcontext.hh:185
void setEnableStorageCache(bool yesno)
Specifies if the cache for the storage term ought to be used for this context.
Definition: fvbaseelementcontext.hh:536
const ExtensiveQuantities & extensiveQuantities(unsigned fluxIdx, unsigned timeIdx OPM_UNUSED) const
Return a reference to the extensive quantities of a sub-control volume face.
Definition: fvbaseelementcontext.hh:521
void updatePrimaryIntensiveQuantities(unsigned timeIdx)
Compute the intensive quantities of all sub-control volumes of the current element for a single time ...
Definition: fvbaseelementcontext.hh:215
Definition: baseauxiliarymodule.hh:37
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:189
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
void updateAll(const Element &elem)
Construct all volume and extensive quantities of an element from scratch.
Definition: fvbaseelementcontext.hh:121
void stashIntensiveQuantities(unsigned dofIdx)
Stash the intensive quantities for a degree of freedom on internal memory.
Definition: fvbaseelementcontext.hh:485
const Stencil & stencil(unsigned timeIdx OPM_UNUSED) const
Return the current finite element geometry.
Definition: fvbaseelementcontext.hh:336
Problem & problem()
Return the object which specifies the pysical setup of the simulation.
Definition: simulator.hh:202
void updateIntensiveQuantities(const PrimaryVariables &priVars, unsigned dofIdx, unsigned timeIdx)
Compute the intensive quantities of a single sub-control volume of the current element for a single t...
Definition: fvbaseelementcontext.hh:228
unsigned focusDofIndex() const
Returns the degree of freedom on which the simulator is currently "focused" on.
Definition: fvbaseelementcontext.hh:271
const Model & model() const
Return a reference to the model.
Definition: fvbaseelementcontext.hh:289
void updateIntensiveQuantities_(unsigned timeIdx, size_t numDof)
Update the first &#39;n&#39; intensive quantities objects from the primary variables.
Definition: fvbaseelementcontext.hh:552
Declare the properties used by the infrastructure code of the finite volume discretizations.
int stashedDofIdx() const
Return the (local) index of the DOF for which the primary variables were stashed. ...
Definition: fvbaseelementcontext.hh:477
const IntensiveQuantities * thermodynamicHint(unsigned dofIdx, unsigned timeIdx) const
Return the thermodynamic hint for a given local index.
Definition: fvbaseelementcontext.hh:426
const Problem & problem() const
Return a reference to the problem.
Definition: fvbaseelementcontext.hh:283
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
bool haveStashedIntensiveQuantities() const
Returns true if no intensive quanties are stashed.
Definition: fvbaseelementcontext.hh:468
const GlobalPosition & pos(unsigned dofIdx, unsigned timeIdx OPM_UNUSED) const
Return the position of a local entities in global coordinates.
Definition: fvbaseelementcontext.hh:347
const IntensiveQuantities & intensiveQuantities(unsigned dofIdx, unsigned timeIdx) const
Return a reference to the intensive quantities of a sub-control volume at a given time...
Definition: fvbaseelementcontext.hh:404
FvBaseElementContext(const Simulator &simulator)
The constructor.
Definition: fvbaseelementcontext.hh:95
const Element & element() const
Return the current element.
Definition: fvbaseelementcontext.hh:301
This class stores an array of IntensiveQuantities objects, one intensive quantities object for each o...
Definition: fvbaseelementcontext.hh:52
IntensiveQuantities & intensiveQuantities(unsigned dofIdx, unsigned timeIdx)
Return a reference to the intensive quantities of a sub-control volume at a given time...
Definition: fvbaseelementcontext.hh:434
void updateIntensiveQuantities(unsigned timeIdx)
Compute the intensive quantities of all sub-control volumes of the current element for a single time ...
Definition: fvbaseelementcontext.hh:206
size_t numPrimaryDof(unsigned timeIdx) const
Return the number of primary degrees of freedom of the current element.
Definition: fvbaseelementcontext.hh:313
PrimaryVariables & primaryVars(unsigned dofIdx, unsigned timeIdx)
Return the primary variables for a given local index.
Definition: fvbaseelementcontext.hh:448
bool onBoundary() const
Returns whether the current element is on the domain&#39;s boundary.
Definition: fvbaseelementcontext.hh:391
void updateExtensiveQuantities(unsigned timeIdx)
Compute the extensive quantities of all sub-control volume faces of the current element for a single ...
Definition: fvbaseelementcontext.hh:245
void restoreIntensiveQuantities(unsigned dofIdx)
Restores the intensive quantities for a degree of freedom from internal memory.
Definition: fvbaseelementcontext.hh:499
const GridView & gridView() const
Return a reference to the grid view.
Definition: fvbaseelementcontext.hh:295
size_t numInteriorFaces(unsigned timeIdx) const
Return the number of non-boundary faces which need to be considered for the flux apporixmation.
Definition: fvbaseelementcontext.hh:320
const GradientCalculator & gradientCalculator() const
Return a reference to the gradient calculation class of the chosen spatial discretization.
Definition: fvbaseelementcontext.hh:510
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:75
void updateStencil(const Element &elem)
Compute the finite volume geometry for an element.
Definition: fvbaseelementcontext.hh:134
void updatePrimaryStencil(const Element &elem)
Update the primary topological part of the stencil, but nothing else.
Definition: fvbaseelementcontext.hh:155
size_t numBoundaryFaces(unsigned timeIdx) const
Return the number of boundary faces which need to be considered for the flux apporixmation.
Definition: fvbaseelementcontext.hh:327
unsigned globalSpaceIndex(unsigned dofIdx, unsigned timeIdx) const
Return the global spatial index for a sub-control volume.
Definition: fvbaseelementcontext.hh:358
void setFocusDofIndex(unsigned dofIdx)
Sets the degree of freedom on which the simulator is currently "focused" on.
Definition: fvbaseelementcontext.hh:263
Scalar dofTotalVolume(unsigned dofIdx, unsigned timeIdx) const
Return the total volume associated with a degree of freedom.
Definition: fvbaseelementcontext.hh:384
const Simulator & simulator() const
Return a reference to the simulator.
Definition: fvbaseelementcontext.hh:277
void updateAllExtensiveQuantities()
Compute the extensive quantities of all sub-control volume faces of the current element for all time ...
Definition: fvbaseelementcontext.hh:235
const PrimaryVariables & primaryVars(unsigned dofIdx, unsigned timeIdx) const
Return the primary variables for a given local index.
Definition: fvbaseelementcontext.hh:456
Scalar dofVolume(unsigned dofIdx, unsigned timeIdx) const
Return the element-local volume associated with a degree of freedom.
Definition: fvbaseelementcontext.hh:371