All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Pages
ExplicitArraysSatDerivativesFluidState.hpp
1 /*
2  Copyright 2015 Andreas Lauser
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 3 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 
20 #ifndef OPM_EXPLICIT_ARRAYS_SAT_DERIVATIVES_FLUID_STATE_HEADER_INCLUDED
21 #define OPM_EXPLICIT_ARRAYS_SAT_DERIVATIVES_FLUID_STATE_HEADER_INCLUDED
22 
23 #include <opm/material/densead/Evaluation.hpp>
24 #include <opm/material/densead/Math.hpp>
25 
26 #include <opm/core/props/BlackoilPhases.hpp>
27 
28 #include <array>
29 
30 namespace Opm
31 {
32 
42 {
43 public:
44  enum { numPhases = BlackoilPhases::MaxNumPhases };
45  enum { numComponents = 3 };
46 
47  typedef Opm::DenseAd::Evaluation<double, numPhases> Evaluation;
48  typedef Evaluation Scalar;
49 
51  : phaseUsage_(phaseUsage)
52  {
53  globalSaturationArray_ = 0;
54 
55  // initialize the evaluation objects for the saturations
56  for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
57  saturation_[phaseIdx] = Evaluation::createVariable(0.0, phaseIdx);
58  }
59  }
60 
67  void setIndex(unsigned arrayIdx)
68  {
69  int np = phaseUsage_.num_phases;
70 
71  // copy the saturations values from the global value. the derivatives do not need
72  // to be modified for these...
73  for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
74  if (!phaseUsage_.phase_used[phaseIdx]) {
75  saturation_[phaseIdx].setValue(0.0);
76  }
77  else {
78  saturation_[phaseIdx].setValue(globalSaturationArray_[np*arrayIdx + phaseUsage_.phase_pos[phaseIdx]]);
79  }
80  }
81  }
82 
91  void setSaturationArray(const double* saturations)
92  { globalSaturationArray_ = saturations; }
93 
97  const Evaluation& saturation(int phaseIdx) const
98  { return saturation_[phaseIdx]; }
99 
100  // TODO (?) temperature, pressure, composition, etc
101 
102 private:
103  const PhaseUsage phaseUsage_;
104  const double* globalSaturationArray_;
105 
106  std::array<Evaluation, numPhases> saturation_;
107 };
108 
109 } // namespace Opm
110 
111 #endif
This is a fluid state which translates global arrays and translates them to a subset of the fluid sta...
Definition: ExplicitArraysSatDerivativesFluidState.hpp:41
void setIndex(unsigned arrayIdx)
Sets the currently used array index.
Definition: ExplicitArraysSatDerivativesFluidState.hpp:67
void setSaturationArray(const double *saturations)
Set the array containing the phase saturations.
Definition: ExplicitArraysSatDerivativesFluidState.hpp:91
Definition: BlackoilPhases.hpp:36
const Evaluation & saturation(int phaseIdx) const
Returns the saturation of a phase for the current cell index.
Definition: ExplicitArraysSatDerivativesFluidState.hpp:97