All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
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
This is a stand-alone version of boost::alignment::aligned_allocator from Boost 1.58.
unsigned globalSpaceIndex(unsigned dofIdx, unsigned timeIdx) const
Return the global spatial index for a sub-control volume.
Definition: fvbaseelementcontext.hh:358
Scalar dofTotalVolume(unsigned dofIdx, unsigned timeIdx) const
Return the total volume associated with a degree of freedom.
Definition: fvbaseelementcontext.hh:384
void updateStencilTopology(const Element &elem)
Update the topological part of the stencil, but nothing else.
Definition: fvbaseelementcontext.hh:172
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
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
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:189
const GridView & gridView() const
Return a reference to the grid view.
Definition: fvbaseelementcontext.hh:295
#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
bool enableStorageCache() const
Returns true iff the cache for the storage term ought to be used for this context.
Definition: fvbaseelementcontext.hh:530
Problem & problem()
Return the object which specifies the pysical setup of the simulation.
Definition: simulator.hh:202
int stashedDofIdx() const
Return the (local) index of the DOF for which the primary variables were stashed. ...
Definition: fvbaseelementcontext.hh:477
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
size_t numPrimaryDof(unsigned timeIdx) const
Return the number of primary degrees of freedom of the current element.
Definition: fvbaseelementcontext.hh:313
bool onBoundary() const
Returns whether the current element is on the domain&#39;s boundary.
Definition: fvbaseelementcontext.hh:391
unsigned focusDofIndex() const
Returns the degree of freedom on which the simulator is currently &quot;focused&quot; on.
Definition: fvbaseelementcontext.hh:271
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.
const Simulator & simulator() const
Return a reference to the simulator.
Definition: fvbaseelementcontext.hh:277
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 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
#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 Model & model() const
Return a reference to the model.
Definition: fvbaseelementcontext.hh:289
FvBaseElementContext(const Simulator &simulator)
The constructor.
Definition: fvbaseelementcontext.hh:95
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
const PrimaryVariables & primaryVars(unsigned dofIdx, unsigned timeIdx) const
Return the primary variables for a given local index.
Definition: fvbaseelementcontext.hh:456
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
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
const GradientCalculator & gradientCalculator() const
Return a reference to the gradient calculation class of the chosen spatial discretization.
Definition: fvbaseelementcontext.hh:510
PrimaryVariables & primaryVars(unsigned dofIdx, unsigned timeIdx)
Return the primary variables for a given local index.
Definition: fvbaseelementcontext.hh:448
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 Element & element() const
Return the current element.
Definition: fvbaseelementcontext.hh:301
const IntensiveQuantities * thermodynamicHint(unsigned dofIdx, unsigned timeIdx) const
Return the thermodynamic hint for a given local index.
Definition: fvbaseelementcontext.hh:426
const Stencil & stencil(unsigned timeIdx OPM_UNUSED) const
Return the current finite element geometry.
Definition: fvbaseelementcontext.hh:336
size_t numDof(unsigned timeIdx) const
Return the number of sub-control volumes of the current element.
Definition: fvbaseelementcontext.hh:307
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
void setFocusDofIndex(unsigned dofIdx)
Sets the degree of freedom on which the simulator is currently &quot;focused&quot; on.
Definition: fvbaseelementcontext.hh:263
const Problem & problem() const
Return a reference to the problem.
Definition: fvbaseelementcontext.hh:283
void updateAllExtensiveQuantities()
Compute the extensive quantities of all sub-control volume faces of the current element for all time ...
Definition: fvbaseelementcontext.hh:235
Scalar dofVolume(unsigned dofIdx, unsigned timeIdx) const
Return the element-local volume associated with a degree of freedom.
Definition: fvbaseelementcontext.hh:371
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