blackoilprimaryvariables.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_BLACK_OIL_PRIMARY_VARIABLES_HH
29 #define EWOMS_BLACK_OIL_PRIMARY_VARIABLES_HH
30 
31 #include "blackoilproperties.hh"
34 
36 
37 #include <dune/common/fvector.hh>
38 
39 #include <opm/material/constraintsolvers/NcpFlash.hpp>
40 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
41 #include <opm/material/fluidstates/SimpleModularFluidState.hpp>
42 #include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
43 #include <opm/common/Valgrind.hpp>
44 
45 namespace Ewoms {
46 template <class TypeTag, bool enableSolvent>
48 
49 template <class TypeTag, bool enablePolymer>
51 
57 template <class TypeTag>
59 {
60  typedef FvBasePrimaryVariables<TypeTag> ParentType;
61  typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) Implementation;
62 
63  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
64  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
65  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
66  typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
67  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
68  typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
69  typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
70 
71  // number of equations
72  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
73 
74  // primary variable indices
75  enum { waterSaturationIdx = Indices::waterSaturationIdx };
76  enum { pressureSwitchIdx = Indices::pressureSwitchIdx };
77  enum { compositionSwitchIdx = Indices::compositionSwitchIdx };
78 
79  static const bool compositionSwitchEnabled = Indices::gasEnabled;
80  static const bool waterEnabled = Indices::waterEnabled;
81 
82  // phase indices from the fluid system
83  enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
84  enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
85  enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
86  enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
87 
88  // component indices from the fluid system
89  enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
90  enum { enableSolvent = GET_PROP_VALUE(TypeTag, EnableSolvent) };
91  enum { enablePolymer = GET_PROP_VALUE(TypeTag, EnablePolymer) };
92  enum { gasCompIdx = FluidSystem::gasCompIdx };
93  enum { waterCompIdx = FluidSystem::waterCompIdx };
94  enum { oilCompIdx = FluidSystem::oilCompIdx };
95 
96  typedef typename Opm::MathToolbox<Evaluation> Toolbox;
97  typedef Dune::FieldVector<Scalar, numComponents> ComponentVector;
100 
101 
102  static_assert(numPhases == 3, "The black-oil model assumes three phases!");
103  static_assert(numComponents == 3, "The black-oil model assumes three components!");
104 
105 public:
106  enum PrimaryVarsMeaning {
107  Sw_po_Sg, // threephase case
108  Sw_po_Rs, // water + oil case
109  Sw_pg_Rv, // water + gas case
110  };
111 
113  : ParentType()
114  {
115  Opm::Valgrind::SetUndefined(*this);
116  pvtRegionIdx_ = 0;
117  }
118 
123  : ParentType(value)
124  {
125  Opm::Valgrind::SetUndefined(primaryVarsMeaning_);
126  pvtRegionIdx_ = 0;
127  }
128 
132  BlackOilPrimaryVariables(const BlackOilPrimaryVariables& value) = default;
133 
142  void setPvtRegionIndex(unsigned value)
143  { pvtRegionIdx_ = static_cast<unsigned short>(value); }
144 
148  unsigned pvtRegionIndex() const
149  { return pvtRegionIdx_; }
150 
155  PrimaryVarsMeaning primaryVarsMeaning() const
156  { return primaryVarsMeaning_; }
157 
162  void setPrimaryVarsMeaning(PrimaryVarsMeaning newMeaning)
163  { primaryVarsMeaning_ = newMeaning; }
164 
168  template <class FluidState>
169  void assignMassConservative(const FluidState& fluidState,
170  const MaterialLawParams& matParams,
171  bool isInEquilibrium = false)
172  {
173  typedef typename std::remove_reference<typename FluidState::Scalar>::type ConstEvaluation;
174  typedef typename std::remove_const<ConstEvaluation>::type FsEvaluation;
175  typedef typename Opm::MathToolbox<FsEvaluation> FsToolbox;
176 
177 #ifndef NDEBUG
178  // make sure the temperature is the same in all fluid phases
179  for (unsigned phaseIdx = 1; phaseIdx < numPhases; ++phaseIdx) {
180  Opm::Valgrind::CheckDefined(fluidState.temperature(0));
181  Opm::Valgrind::CheckDefined(fluidState.temperature(phaseIdx));
182 
183  assert(fluidState.temperature(0) == fluidState.temperature(phaseIdx));
184  }
185 #endif // NDEBUG
186 
187  // for the equilibrium case, we don't need complicated
188  // computations.
189  if (isInEquilibrium) {
190  assignNaive(fluidState);
191  return;
192  }
193 
194  // If your compiler bails out here, you're probably not using a suitable black
195  // oil fluid system.
196  typename FluidSystem::template ParameterCache<Scalar> paramCache;
197  paramCache.setRegionIndex(pvtRegionIdx_);
198  paramCache.setMaxOilSat(FsToolbox::value(fluidState.saturation(oilPhaseIdx)));
199 
200  // create a mutable fluid state with well defined densities based on the input
201  typedef Opm::NcpFlash<Scalar, FluidSystem> NcpFlash;
202  typedef Opm::CompositionalFluidState<Scalar, FluidSystem> FlashFluidState;
203  FlashFluidState fsFlash;
204  fsFlash.setTemperature(FsToolbox::value(fluidState.temperature(/*phaseIdx=*/0)));
205  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
206  fsFlash.setPressure(phaseIdx, FsToolbox::value(fluidState.pressure(phaseIdx)));
207  fsFlash.setSaturation(phaseIdx, FsToolbox::value(fluidState.saturation(phaseIdx)));
208  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
209  fsFlash.setMoleFraction(phaseIdx, compIdx, FsToolbox::value(fluidState.moleFraction(phaseIdx, compIdx)));
210  }
211 
212  paramCache.updateAll(fsFlash);
213  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
214  if (!FluidSystem::phaseIsActive(phaseIdx))
215  continue;
216 
217  Scalar rho = FluidSystem::template density<FlashFluidState, Scalar>(fsFlash, paramCache, phaseIdx);
218  fsFlash.setDensity(phaseIdx, rho);
219  }
220 
221  // calculate the "global molarities"
222  ComponentVector globalMolarities(0.0);
223  for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
224  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
225  if (!FluidSystem::phaseIsActive(phaseIdx))
226  continue;
227 
228  globalMolarities[compIdx] +=
229  fsFlash.saturation(phaseIdx) * fsFlash.molarity(phaseIdx, compIdx);
230  }
231  }
232 
233  // use a flash calculation to calculate a fluid state in
234  // thermodynamic equilibrium
235 
236  // run the flash calculation
237  NcpFlash::template solve<MaterialLaw>(fsFlash, matParams, paramCache, globalMolarities);
238 
239  // use the result to assign the primary variables
240  assignNaive(fsFlash);
241  }
242 
243  template <class FluidState, class SolventContainer>
244  void assignMassConservative(const FluidState& fluidState,
245  const MaterialLawParams& matParams,
246  Scalar solSat,
247  bool isInEquilibrium = false)
248  {
249  assignMassConservative(fluidState, matParams, isInEquilibrium);
250 
251  // set the primary variables of the solvent module
252  SolventModule::assignPrimaryVars(*this, solSat);
253 
254  // set the primary variables of the polymer module
255  PolymerModule::assignPrimaryVars(*this, solSat);
256  }
257 
261  template <class FluidState>
262  void assignNaive(const FluidState& fluidState)
263  {
264  typedef typename std::remove_reference<typename FluidState::Scalar>::type ConstEvaluation;
265  typedef typename std::remove_const<ConstEvaluation>::type FsEvaluation;
266  typedef typename Opm::MathToolbox<FsEvaluation> FsToolbox;
267 
268  bool gasPresent = (fluidState.saturation(gasPhaseIdx) > 0.0);
269  bool oilPresent = (fluidState.saturation(oilPhaseIdx) > 0.0);
270 
271  // determine the meaning of the primary variables
272  if ((gasPresent && oilPresent) || (!gasPresent && !oilPresent))
273  // gas and oil: both hydrocarbon phases are in equilibrium (i.e., saturated
274  // with the "protagonist" component of the other phase.)
275  primaryVarsMeaning_ = Sw_po_Sg;
276  else if (oilPresent) {
277  // only oil: if dissolved gas is enabled, we need to consider the oil phase
278  // composition, if it is disabled, the gas component must stick to its phase
279  if (FluidSystem::enableDissolvedGas())
280  primaryVarsMeaning_ = Sw_po_Rs;
281  else
282  primaryVarsMeaning_ = Sw_po_Sg;
283  }
284  else {
285  assert(gasPresent);
286  // only gas: if vaporized oil is enabled, we need to consider the gas phase
287  // composition, if it is disabled, the oil component must stick to its phase
288  if (FluidSystem::enableVaporizedOil())
289  primaryVarsMeaning_ = Sw_pg_Rv;
290  else
291  primaryVarsMeaning_ = Sw_po_Sg;
292  }
293 
294  // assign the actual primary variables
295  if (primaryVarsMeaning() == Sw_po_Sg) {
296  if (waterEnabled)
297  (*this)[waterSaturationIdx] = FsToolbox::value(fluidState.saturation(waterPhaseIdx));
298  (*this)[pressureSwitchIdx] = FsToolbox::value(fluidState.pressure(oilPhaseIdx));
299  if( compositionSwitchEnabled )
300  (*this)[compositionSwitchIdx] = FsToolbox::value(fluidState.saturation(gasPhaseIdx));
301  }
302  else if (primaryVarsMeaning() == Sw_po_Rs) {
303  const auto& Rs = Opm::BlackOil::getRs_<FluidSystem, Scalar, FluidState>(fluidState, pvtRegionIdx_);
304 
305  if (waterEnabled)
306  (*this)[waterSaturationIdx] = FsToolbox::value(fluidState.saturation(waterPhaseIdx));
307  (*this)[pressureSwitchIdx] = FsToolbox::value(fluidState.pressure(oilPhaseIdx));
308  if( compositionSwitchEnabled )
309  (*this)[compositionSwitchIdx] = Rs;
310  }
311  else {
312  assert(primaryVarsMeaning() == Sw_pg_Rv);
313 
314  const auto& Rv = Opm::BlackOil::getRv_<FluidSystem, Scalar, FluidState>(fluidState, pvtRegionIdx_);
315  if (waterEnabled)
316  (*this)[waterSaturationIdx] = FsToolbox::value(fluidState.saturation(waterPhaseIdx));
317 
318  (*this)[pressureSwitchIdx] = FsToolbox::value(fluidState.pressure(gasPhaseIdx));
319  if( compositionSwitchEnabled )
320  (*this)[compositionSwitchIdx] = Rv;
321  }
322  }
323 
337  bool adaptPrimaryVariables(const Problem& problem, unsigned globalDofIdx, Scalar eps = 0.0)
338  {
339  static const Scalar thresholdWaterFilledCell = 1.0 - eps;
340 
341  // this function accesses quite a few black-oil specific low-level functions
342  // directly for better performance (instead of going the canonical way through
343  // the IntensiveQuantities). The reason is that most intensive quantities are not
344  // required to be able to decide if the primary variables needs to be switched or
345  // not, so it would be a waste to compute them.
346  Scalar Sw = 0.0;
347  if (waterEnabled)
348  Sw = (*this)[Indices::waterSaturationIdx];
349 
350  if (primaryVarsMeaning() == Sw_po_Sg) {
351 
352  // special case for cells with almost only water
353  if (Sw >= thresholdWaterFilledCell) {
354 
355  // make sure water saturations does not exceed 1.0
356  if (waterEnabled)
357  (*this)[Indices::waterSaturationIdx] = 1.0;
358  // the hydrocarbon gas saturation is set to 0.0
359  if (compositionSwitchEnabled)
360  (*this)[Indices::compositionSwitchIdx] = 0.0;
361 
362  return false;
363  }
364 
365  // phase equilibrium, i.e., both hydrocarbon phases are present.
366  Scalar Sg = 0.0;
367  if (compositionSwitchEnabled)
368  Sg = (*this)[Indices::compositionSwitchIdx];
369 
370  Scalar So = 1.0 - Sw - Sg - solventSaturation();
371 
372  Scalar So2 = 1.0 - Sw - solventSaturation();
373  if (Sg < -eps && So2 > 0.0 && FluidSystem::enableDissolvedGas()) {
374  // the hydrocarbon gas phase disappeared and some oil phase is left,
375  // i.e., switch the primary variables to { Sw, po, Rs }.
376  //
377  // by a "lucky" coincidence the pressure switching variable already
378  // represents the oil phase pressure, so we do not need to change
379  // this. For the gas dissolution factor, we use the low-level blackoil
380  // PVT objects to calculate the mole fraction of gas saturated oil.
381  Scalar po = (*this)[Indices::pressureSwitchIdx];
382  Scalar T = asImp_().temperature_();
383  Scalar SoMax = problem.model().maxOilSaturation(globalDofIdx);
384  Scalar RsSat = FluidSystem::oilPvt().saturatedGasDissolutionFactor(pvtRegionIdx_,
385  T,
386  po,
387  So2,
388  SoMax);
389  setPrimaryVarsMeaning(Sw_po_Rs);
390  if (compositionSwitchEnabled)
391  (*this)[Indices::compositionSwitchIdx] = RsSat;
392 
393  // because more than one primary variable switch can occur at a time,
394  // call this method recursively
395  asImp_().adaptPrimaryVariables(problem, globalDofIdx, eps);
396 
397  return true;
398  }
399 
400  Scalar Sg2 = 1.0 - Sw - solventSaturation();
401  if (So < -eps && Sg2 > 0.0 && FluidSystem::enableVaporizedOil()) {
402  // the oil phase disappeared and some hydrocarbon gas phase is still
403  // present, i.e., switch the primary variables to { Sw, pg, Rv }.
404  Scalar po = (*this)[Indices::pressureSwitchIdx];
405 
406  // we only have the oil pressure readily available, but we need the gas
407  // pressure, i.e. we must determine capillary pressure
408  Scalar pC[numPhases] = { 0.0 };
409  const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
410  computeCapillaryPressures_(pC, /*So=*/0.0, Sg2 + solventSaturation(), Sw, matParams);
411  Scalar pg = po + (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
412 
413  // we start at the Rv value that corresponds to that of oil-saturated
414  // hydrocarbon gas
415  Scalar T = asImp_().temperature_();
416  Scalar SoMax = problem.model().maxOilSaturation(globalDofIdx);
417  Scalar RvSat =
418  FluidSystem::gasPvt().saturatedOilVaporizationFactor(pvtRegionIdx_,
419  T,
420  pg,
421  0.0,
422  SoMax);
423 
424  setPrimaryVarsMeaning(Sw_pg_Rv);
425  (*this)[Indices::pressureSwitchIdx] = pg;
426  if (compositionSwitchEnabled)
427  (*this)[Indices::compositionSwitchIdx] = RvSat;
428 
429  // because more than one primary variable switch can occur at a time,
430  // call this method recursively
431  asImp_().adaptPrimaryVariables(problem, globalDofIdx, eps);
432 
433  return true;
434  }
435 
436  return false;
437  }
438  else if (primaryVarsMeaning() == Sw_po_Rs) {
439 
440  assert(compositionSwitchEnabled);
441 
442  // special case for cells with almost only water
443  if (Sw >= thresholdWaterFilledCell) {
444  // switch back to phase equilibrium mode if the oil phase vanishes (i.e.,
445  // the water-only case)
446  setPrimaryVarsMeaning(Sw_po_Sg);
447  if (waterEnabled)
448  (*this)[Indices::waterSaturationIdx] = 1.0; // water saturation
449 
450  (*this)[Indices::compositionSwitchIdx] = 0.0; // hydrocarbon gas saturation
451 
452  // because more than one primary variable switch can occur at a time,
453  // call this method recursively
454  asImp_().adaptPrimaryVariables(problem, globalDofIdx, eps);
455 
456  return true;
457  }
458 
459  // Only the oil and the water phases are present. The hydrocarbon gas phase
460  // appears as soon as more of the gas component is present in the oil phase
461  // than what saturated oil can hold.
462  Scalar T = asImp_().temperature_();
463  Scalar po = (*this)[Indices::pressureSwitchIdx];
464  Scalar So = 1.0 - Sw - solventSaturation();
465  Scalar SoMax = std::max(So, problem.model().maxOilSaturation(globalDofIdx));
466  Scalar RsSat =
467  FluidSystem::oilPvt().saturatedGasDissolutionFactor(pvtRegionIdx_, T, po, So, SoMax);
468  Scalar Rs = (*this)[Indices::compositionSwitchIdx];
469  if (Rs > RsSat*(1.0 + eps)) {
470  // the gas phase appears, i.e., switch the primary variables to { Sw, po,
471  // Sg }.
472  setPrimaryVarsMeaning(Sw_po_Sg);
473  (*this)[Indices::compositionSwitchIdx] = 0.0; // hydrocarbon gas saturation
474 
475  // because more than one primary variable switch can occur at a time,
476  // call this method recursively
477  asImp_().adaptPrimaryVariables(problem, globalDofIdx, eps);
478 
479  return true;
480  }
481 
482  return false;
483  }
484  else {
485  assert(primaryVarsMeaning() == Sw_pg_Rv);
486  assert(compositionSwitchEnabled);
487 
488  Scalar pg = (*this)[Indices::pressureSwitchIdx];
489  Scalar Sg = 1.0 - Sw - solventSaturation();
490 
491  // special case for cells with almost only water
492  if (Sw >= thresholdWaterFilledCell) {
493  // switch to phase equilibrium mode because the hydrocarbon gas phase
494  // disappears. here we need the capillary pressures to calculate the oil
495  // phase pressure using the gas phase pressure
496  Scalar pC[numPhases] = { 0.0 };
497  const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
498  computeCapillaryPressures_(pC,
499  /*So=*/0.0,
500  /*Sg=*/Sg + solventSaturation(),
501  Sw,
502  matParams);
503  Scalar po = pg + (pC[oilPhaseIdx] - pC[gasPhaseIdx]);
504 
505  setPrimaryVarsMeaning(Sw_po_Sg);
506  if (waterEnabled)
507  (*this)[Indices::waterSaturationIdx] = 1.0;
508 
509  (*this)[Indices::pressureSwitchIdx] = po;
510  (*this)[Indices::compositionSwitchIdx] = 0.0; // hydrocarbon gas saturation
511 
512  // because more than one primary variable switch can occur at a time,
513  // call this method recursively
514  asImp_().adaptPrimaryVariables(problem, globalDofIdx, eps);
515 
516  return true;
517  }
518 
519  // Only the gas and the water phases are present. The oil phase appears as
520  // soon as more of the oil component is present in the hydrocarbon gas phase
521  // than what saturated gas contains. Note that we use the blackoil specific
522  // low-level PVT objects here for performance reasons.
523  Scalar T = asImp_().temperature_();
524  Scalar SoMax = problem.model().maxOilSaturation(globalDofIdx);
525  Scalar RvSat =
526  FluidSystem::gasPvt().saturatedOilVaporizationFactor(pvtRegionIdx_,
527  T,
528  pg,
529  /*So=*/Scalar(0.0),
530  SoMax);
531 
532  Scalar Rv = (*this)[Indices::compositionSwitchIdx];
533  if (Rv > RvSat*(1.0 + eps)) {
534  // switch to phase equilibrium mode because the oil phase appears. here
535  // we also need the capillary pressures to calculate the oil phase
536  // pressure using the gas phase pressure
537  Scalar pC[numPhases] = { 0.0 };
538  const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
539  computeCapillaryPressures_(pC,
540  /*So=*/0.0,
541  /*Sg=*/Sg + solventSaturation(),
542  Sw,
543  matParams);
544  Scalar po = pg + (pC[oilPhaseIdx] - pC[gasPhaseIdx]);
545 
546  setPrimaryVarsMeaning(Sw_po_Sg);
547  (*this)[Indices::pressureSwitchIdx] = po;
548  (*this)[Indices::compositionSwitchIdx] = Sg; // hydrocarbon gas saturation
549 
550  // because more than one primary variable switch can occur at a time,
551  // call this method recursively
552  asImp_().adaptPrimaryVariables(problem, globalDofIdx, eps);
553 
554  return true;
555  }
556 
557  return false;
558  }
559 
560  assert(false);
561  return false;
562  }
563 
564  BlackOilPrimaryVariables& operator=(const BlackOilPrimaryVariables& other) = default;
565  BlackOilPrimaryVariables& operator=(Scalar value)
566  {
567  for (unsigned i = 0; i < numEq; ++i)
568  (*this)[i] = value;
569 
570  return *this;
571  }
572 
573 private:
574  Implementation& asImp_()
575  { return *static_cast<Implementation*>(this); }
576 
577  const Implementation& asImp_() const
578  { return *static_cast<const Implementation*>(this); }
579 
580  Scalar solventSaturation() const
581  {
582  if (!enableSolvent)
583  return 0.0;
584 
585  return (*this)[Indices::solventSaturationIdx];
586  }
587 
588  Scalar polymerConcentration() const
589  {
590  if (!enablePolymer)
591  return 0.0;
592 
593  return (*this)[Indices::polymerConcentrationIdx];
594  }
595 
596  template <class Container>
597  void computeCapillaryPressures_(Container& result,
598  Scalar So,
599  Scalar Sg,
600  Scalar Sw,
601  const MaterialLawParams& matParams) const
602  {
603  typedef Opm::SimpleModularFluidState<Scalar,
604  numPhases,
605  numComponents,
606  FluidSystem,
607  /*storePressure=*/false,
608  /*storeTemperature=*/false,
609  /*storeComposition=*/false,
610  /*storeFugacity=*/false,
611  /*storeSaturation=*/true,
612  /*storeDensity=*/false,
613  /*storeViscosity=*/false,
614  /*storeEnthalpy=*/false> SatOnlyFluidState;
615  SatOnlyFluidState fluidState;
616  fluidState.setSaturation(waterPhaseIdx, Sw);
617  fluidState.setSaturation(oilPhaseIdx, So);
618  fluidState.setSaturation(gasPhaseIdx, Sg);
619 
620  MaterialLaw::capillaryPressures(result, matParams, fluidState);
621  }
622 
623  // the standard blackoil model is isothermal
624  Scalar temperature_() const
625  { return FluidSystem::surfaceTemperature; }
626 
627  PrimaryVarsMeaning primaryVarsMeaning_;
628  unsigned short pvtRegionIdx_;
629 };
630 
631 
632 } // namespace Ewoms
633 
634 #endif
Definition: baseauxiliarymodule.hh:37
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
Contains the classes required to extend the black-oil model by solvents.
static void assignPrimaryVars(PrimaryVariables &priVars, Scalar polymerConcentration)
Assign the polymer specific primary variables to a PrimaryVariables object.
Definition: blackoilpolymermodules.hh:511
unsigned pvtRegionIndex() const
Return the index of the region which should be used for PVT properties.
Definition: blackoilprimaryvariables.hh:148
Declares the properties required by the black oil model.
bool adaptPrimaryVariables(const Problem &problem, unsigned globalDofIdx, Scalar eps=0.0)
Adapt the interpretation of the switching variables to be physically meaningful.
Definition: blackoilprimaryvariables.hh:337
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilprimaryvariables.hh:47
Represents the primary variables used by the a model.
Definition: fvbaseprimaryvariables.hh:48
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
PrimaryVarsMeaning primaryVarsMeaning() const
Return the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:155
void setPvtRegionIndex(unsigned value)
Set the index of the region which should be used for PVT properties.
Definition: blackoilprimaryvariables.hh:142
void setPrimaryVarsMeaning(PrimaryVarsMeaning newMeaning)
Set the interpretation which should be applied to the switching primary variables.
Definition: blackoilprimaryvariables.hh:162
Represents the primary variables used by the a model.
void assignMassConservative(const FluidState &fluidState, const MaterialLawParams &matParams, bool isInEquilibrium=false)
Set the primary variables from an arbitrary fluid state in a mass conservative way.
Definition: blackoilprimaryvariables.hh:169
Represents the primary variables used by the black-oil model.
Definition: blackoilprimaryvariables.hh:58
void assignNaive(const FluidState &fluidState)
Directly retrieve the primary variables from an arbitrary fluid state.
Definition: blackoilprimaryvariables.hh:262
BlackOilPrimaryVariables(Scalar value)
Constructor with assignment from scalar.
Definition: blackoilprimaryvariables.hh:122
static void assignPrimaryVars(PrimaryVariables &priVars, Scalar solventSaturation)
Assign the solvent specific primary variables to a PrimaryVariables object.
Definition: blackoilsolventmodules.hh:615
Contains the classes required to extend the black-oil model by polymer.