blackoilmodel.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_MODEL_HH
29 #define EWOMS_BLACK_OIL_MODEL_HH
30 
31 #include <opm/material/densead/Math.hpp>
32 
33 #include "blackoilproblem.hh"
34 #include "blackoilindices.hh"
39 #include "blackoilratevector.hh"
41 #include "blackoillocalresidual.hh"
42 #include "blackoilnewtonmethod.hh"
43 #include "blackoilproperties.hh"
47 
51 
52 #include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
53 #include <opm/common/Unused.hpp>
54 #include <opm/common/ErrorMacros.hpp>
55 #include <opm/common/Exceptions.hpp>
56 
57 #include <sstream>
58 #include <string>
59 
60 namespace Ewoms {
61 template <class TypeTag>
63 
64 template <class TypeTag>
66 }
67 
68 namespace Ewoms {
69 namespace Properties {
72  VtkBlackOil,
73  VtkBlackOilSolvent,
74  VtkBlackOilPolymer,
75  VtkComposition));
76 
78 SET_TYPE_PROP(BlackOilModel, LocalResidual,
80 
83 
86 
89 
92 
95 
98 
101 
104 
108 
111  Ewoms::BlackOilIndices<GET_PROP_VALUE(TypeTag, EnableSolvent)?1:0, GET_PROP_VALUE(TypeTag, EnablePolymer)?1:0, /*PVOffset=*/0>);
112 
114 SET_PROP(BlackOilModel, FluidSystem)
115 {
116 private:
117  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
118  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
119 
120 public:
121  typedef Opm::FluidSystems::BlackOil<Scalar> type;
122 };
123 
124 // by default, the ECL solvent module is disabled
125 SET_BOOL_PROP(BlackOilModel, EnableSolvent, false);
126 SET_BOOL_PROP(BlackOilModel, EnablePolymer, false);
127 // by default, ebos formulates the conservation equations in terms of mass not surface volumes
128 SET_BOOL_PROP(BlackOilModel, BlackoilConserveSurfaceVolume, false);
129 } // namespace Properties
130 
193 template<class TypeTag >
194 class BlackOilModel
195  : public MultiPhaseBaseModel<TypeTag>
196 {
197  typedef typename GET_PROP_TYPE(TypeTag, Model) Implementation;
198  typedef MultiPhaseBaseModel<TypeTag> ParentType;
199 
200  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
201  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
202  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
203  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
204  typedef typename GET_PROP_TYPE(TypeTag, Discretization) Discretization;
205  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
206  typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
207 
208  enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
209  enum { numComponents = FluidSystem::numComponents };
210  enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
211 
212  static const bool compositionSwitchEnabled = Indices::gasEnabled;
213  static const bool waterEnabled = Indices::waterEnabled;
214 
215  typedef BlackOilSolventModule<TypeTag> SolventModule;
216  typedef BlackOilPolymerModule<TypeTag> PolymerModule;
217 
218 
219 public:
220  BlackOilModel(Simulator& simulator)
221  : ParentType(simulator)
222  {}
223 
227  static void registerParameters()
228  {
230 
233 
234  // register runtime parameters of the VTK output modules
237  }
238 
242  void finishInit()
243  {
244  maxOilSaturation_.resize(this->numGridDof(), 0.0);
246 
247  Dune::FMatrixPrecision<Scalar>::set_singular_limit(1e-35);
248  }
249 
253  static std::string name()
254  { return "blackoil"; }
255 
259  std::string primaryVarName(int pvIdx) const
260  {
261  std::ostringstream oss;
262 
263  if (pvIdx == Indices::waterSaturationIdx)
264  oss << "saturation_" << FluidSystem::phaseName(FluidSystem::waterPhaseIdx);
265  else if (pvIdx == Indices::pressureSwitchIdx)
266  oss << "pressure_switching";
267  else if (static_cast<int>(pvIdx) == Indices::compositionSwitchIdx)
268  oss << "composition_switching";
269  else if (SolventModule::primaryVarApplies(pvIdx))
270  return SolventModule::primaryVarName(pvIdx);
271  else if (PolymerModule::primaryVarApplies(pvIdx))
272  return PolymerModule::primaryVarName(pvIdx);
273  else
274  assert(false);
275 
276  return oss.str();
277  }
278 
282  std::string eqName(int eqIdx) const
283  {
284  std::ostringstream oss;
285 
286  if (Indices::conti0EqIdx <= eqIdx && eqIdx < Indices::conti0EqIdx + numComponents)
287  oss << "conti_" << FluidSystem::phaseName(eqIdx - Indices::conti0EqIdx);
288  else if (SolventModule::eqApplies(eqIdx))
289  return SolventModule::eqName(eqIdx);
290  else if (PolymerModule::eqApplies(eqIdx))
291  return PolymerModule::eqName(eqIdx);
292  else
293  assert(false);
294 
295  return oss.str();
296  }
297 
301  Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
302  {
303  // do not care about the auxiliary equations as they are supposed to scale
304  // themselves
305  if (globalDofIdx >= this->numGridDof())
306  return 1.0;
307 
308  // saturations are always in the range [0, 1]!
309  if (Indices::waterSaturationIdx == pvIdx)
310  return 1.0;
311 
312  // oil pressures usually are in the range of 100 to 500 bars for typical oil
313  // reservoirs (which is the only relevant application for the black-oil model).
314  else if (Indices::pressureSwitchIdx == pvIdx)
315  return 1.0/300e5;
316 
317  // deal with primary variables stemming from the solvent module
318  else if (SolventModule::primaryVarApplies(pvIdx))
319  return SolventModule::primaryVarWeight(pvIdx);
320 
321  // deal with primary variables stemming from the polymer module
322  else if (PolymerModule::primaryVarApplies(pvIdx))
323  return PolymerModule::primaryVarWeight(pvIdx);
324 
325  // if the primary variable is either the gas saturation, Rs or Rv
326  assert(Indices::compositionSwitchIdx == pvIdx);
327 
328  auto pvMeaning = this->solution(0)[globalDofIdx].primaryVarsMeaning();
329  if (pvMeaning == PrimaryVariables::Sw_po_Sg)
330  return 1.0; // gas saturation
331  else if (pvMeaning == PrimaryVariables::Sw_po_Rs)
332  return 1.0/250.; // gas dissolution factor
333  else {
334  assert(pvMeaning == PrimaryVariables::Sw_pg_Rv);
335  return 1.0/0.025; // oil vaporization factor
336  }
337 
338  }
339 
343  Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx OPM_UNUSED) const
344  {
345  // do not care about the auxiliary equations as they are supposed to scale
346  // themselves
347  if (globalDofIdx >= this->numGridDof())
348  return 1.0;
349 
350  else if (SolventModule::eqApplies(eqIdx))
351  return SolventModule::eqWeight(eqIdx);
352 
353  else if (PolymerModule::eqApplies(eqIdx))
354  return PolymerModule::eqWeight(eqIdx);
355 
356  // it is said that all kilograms are equal!
357  return 1.0;
358  }
359 
368  template <class DofEntity>
369  void serializeEntity(std::ostream& outstream, const DofEntity& dof)
370  {
371  unsigned dofIdx = static_cast<unsigned>(asImp_().dofMapper().index(dof));
372 
373  // write phase state
374  if (!outstream.good()) {
375  OPM_THROW(std::runtime_error, "Could not serialize degree of freedom " << dofIdx);
376  }
377 
378  // write the primary variables
379  const auto& priVars = this->solution(/*timeIdx=*/0)[dofIdx];
380  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx)
381  outstream << priVars[eqIdx] << " ";
382 
383  // write the pseudo primary variables
384  outstream << priVars.primaryVarsMeaning() << " ";
385  outstream << priVars.pvtRegionIndex() << " ";
386 
387  if (maxOilSaturation_.size() > 0)
388  outstream << maxOilSaturation_[dofIdx] << " ";
389 
390  SolventModule::serializeEntity(*this, outstream, dof);
391  PolymerModule::serializeEntity(*this, outstream, dof);
392 
393  }
394 
403  template <class DofEntity>
404  void deserializeEntity(std::istream& instream,
405  const DofEntity& dof)
406  {
407  unsigned dofIdx = static_cast<unsigned>(asImp_().dofMapper().index(dof));
408 
409  // read in the "real" primary variables of the DOF
410  auto& priVars = this->solution(/*timeIdx=*/0)[dofIdx];
411  for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
412  if (!instream.good())
413  OPM_THROW(std::runtime_error,
414  "Could not deserialize degree of freedom " << dofIdx);
415  instream >> priVars[eqIdx];
416  }
417 
418  // read the pseudo primary variables
419  unsigned primaryVarsMeaning;
420  instream >> primaryVarsMeaning;
421 
422  unsigned pvtRegionIdx;
423  instream >> pvtRegionIdx;
424 
425  if (maxOilSaturation_.size() > 0)
426  instream >> maxOilSaturation_[dofIdx];
427 
428  if (!instream.good())
429  OPM_THROW(std::runtime_error,
430  "Could not deserialize degree of freedom " << dofIdx);
431 
432  SolventModule::deserializeEntity(*this, instream, dof);
433  PolymerModule::deserializeEntity(*this, instream, dof);
434 
435  typedef typename PrimaryVariables::PrimaryVarsMeaning PVM;
436  priVars.setPrimaryVarsMeaning(static_cast<PVM>(primaryVarsMeaning));
437  priVars.setPvtRegionIndex(pvtRegionIdx);
438  }
439 
447  template <class Restarter>
448  void deserialize(Restarter& res)
449  {
450  ParentType::deserialize(res);
451 
452  // set the PVT indices of the primary variables. This is also done by writing
453  // them into the restart file and re-reading them, but it is better to calculate
454  // them from scratch because the input could have been changed in this regard...
455  ElementContext elemCtx(this->simulator_);
456  auto elemIt = this->gridView().template begin</*codim=*/0>();
457  auto elemEndIt = this->gridView().template end</*codim=*/0>();
458  for (; elemIt != elemEndIt; ++ elemIt) {
459  elemCtx.updateStencil(*elemIt);
460  for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timIdx=*/0); ++dofIdx) {
461  unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timIdx=*/0);
462  updatePvtRegionIndex_(this->solution(/*timeIdx=*/0)[globalDofIdx],
463  elemCtx,
464  dofIdx,
465  /*timeIdx=*/0);
466  }
467  }
468 
469  this->solution(/*timeIdx=*/1) = this->solution(/*timeIdx=*/0);
470  }
471 
479  Scalar maxOilSaturation(unsigned globalDofIdx) const
480  { return maxOilSaturation_[globalDofIdx]; }
481 
486  void setMaxOilSaturation(const Scalar value, unsigned globalDofIdx)
487  { maxOilSaturation_[globalDofIdx] = value; }
488 
498  {
499  if (maxOilSaturation_.size() > 0) {
500  unsigned nGridDofs = this->numGridDof();
501  assert(maxOilSaturation_.size() == nGridDofs);
502  for (unsigned dofIdx = 0; dofIdx < nGridDofs; ++dofIdx) {
503  const PrimaryVariables& priVars = this->solution(/*timeIdx=*/0)[dofIdx];
504  Scalar So = 0.0;
505  switch (priVars.primaryVarsMeaning()) {
506  case PrimaryVariables::Sw_po_Sg:
507  So = 1.0;
508  if( waterEnabled)
509  So -= priVars[Indices::waterSaturationIdx];
510  if( compositionSwitchEnabled )
511  So -= priVars[Indices::compositionSwitchIdx];
512  break;
513  case PrimaryVariables::Sw_pg_Rv:
514  So = 0.0;
515  break;
516  case PrimaryVariables::Sw_po_Rs:
517  So = 1.0;
518  if (waterEnabled)
519  So -= priVars[Indices::waterSaturationIdx];
520  break;
521  }
522 
523  maxOilSaturation_[dofIdx] = std::max(maxOilSaturation_[dofIdx], So);
524  }
525  }
526  }
527 
528 /*
529  // hack: this interferes with the static polymorphism trick
530 protected:
531  friend ParentType;
532  friend Discretization;
533 */
534 
535  template <class Context>
536  void supplementInitialSolution_(PrimaryVariables& priVars,
537  const Context& context,
538  unsigned dofIdx,
539  unsigned timeIdx)
540  { updatePvtRegionIndex_(priVars, context, dofIdx, timeIdx); }
541 
542  void registerOutputModules_()
543  {
544  ParentType::registerOutputModules_();
545 
546  // add the VTK output modules which make sense for the blackoil model
547  SolventModule::registerOutputModules(*this, this->simulator_);
548  PolymerModule::registerOutputModules(*this, this->simulator_);
549 
550  this->addOutputModule(new Ewoms::VtkBlackOilModule<TypeTag>(this->simulator_));
551  this->addOutputModule(new Ewoms::VtkCompositionModule<TypeTag>(this->simulator_));
552  }
553 
554 private:
555  Implementation& asImp_()
556  { return *static_cast<Implementation*>(this); }
557  const Implementation& asImp_() const
558  { return *static_cast<const Implementation*>(this); }
559 
560  template <class Context>
561  void updatePvtRegionIndex_(PrimaryVariables& priVars,
562  const Context& context,
563  unsigned dofIdx,
564  unsigned timeIdx)
565  {
566  unsigned regionIdx = context.problem().pvtRegionIndex(context, dofIdx, timeIdx);
567  priVars.setPvtRegionIndex(regionIdx);
568  }
569 
570  std::vector<Scalar> maxOilSaturation_;
571 };
572 } // namespace Ewoms
573 
574 #endif
#define SET_BOOL_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant boolean value.
Definition: propertysystem.hh:361
void deserializeEntity(std::istream &instream, const DofEntity &dof)
Reads the current solution variables for a degree of freedom from a restart file. ...
Definition: blackoilmodel.hh:404
Base class for all problems which use the black-oil model.
Definition: blackoilmodel.hh:65
Definition: baseauxiliarymodule.hh:37
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilmodule.hh:123
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkcompositionmodule.hh:105
Contains the quantities which are are constant within a finite volume in the black-oil model...
std::string eqName(int eqIdx) const
Given an equation index, return a human readable name.
Definition: blackoilmodel.hh:282
void setMaxOilSaturation(const Scalar value, unsigned globalDofIdx)
Sets an elements maximum oil phase saturation observed during the simulation (used for restarting fro...
Definition: blackoilmodel.hh:486
void serializeEntity(std::ostream &outstream, const DofEntity &dof)
Write the current solution for a degree of freedom to a restart file.
Definition: blackoilmodel.hh:369
Implements a boundary vector for the fully implicit black-oil model.
Implements a vector representing mass, molar or volumetric rates for the black oil model...
Definition: blackoilratevector.hh:50
VTK output module for the black oil model&#39;s parameters.
Definition: vtkblackoilmodule.hh:89
#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 std::string name()
Definition: blackoilmodel.hh:253
Scalar eqWeight(unsigned globalDofIdx, unsigned eqIdx OPM_UNUSED) const
Returns the relative weight of an equation.
Definition: blackoilmodel.hh:343
static void registerOutputModules(Model &model, Simulator &simulator)
Register all solvent specific VTK and ECL output modules.
Definition: blackoilsolventmodules.hh:500
A newton solver which is specific to the black oil model.
Definition: blackoilnewtonmethod.hh:45
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
A fully-implicit black-oil flow model.
Definition: blackoilmodel.hh:62
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition: blackoilmodel.hh:227
Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
Returns the relative weight of a primary variable for calculating relative errors.
Definition: blackoilmodel.hh:301
Scalar maxOilSaturation(unsigned globalDofIdx) const
Returns an elements maximum oil phase saturation observed during the simulation.
Definition: blackoilmodel.hh:479
#define INHERITS_FROM(...)
Syntactic sugar for NEW_TYPE_TAG.
Definition: propertysystem.hh:230
static void registerParameters()
Register all run-time parameters for the black-oil solvent module.
Definition: blackoilsolventmodules.hh:488
static void registerOutputModules(Model &model, Simulator &simulator)
Register all polymer specific VTK and ECL output modules.
Definition: blackoilpolymermodules.hh:377
Provides a Darcy flux module for the blackoil model.
Definition: blackoildarcyfluxmodule.hh:48
Declares the properties required by the black oil model.
This file contains the default flux module of the blackoil model.
void finishInit()
Apply the initial conditions to the model.
Definition: multiphasebasemodel.hh:159
Calculates the local residual of the black oil model.
Base class for all problems which use the black-oil model.
Definition: blackoilproblem.hh:44
VTK output module for the fluid composition.
void deserialize(Restarter &res)
Deserializes the state of the model.
Definition: blackoilmodel.hh:448
Implements a boundary vector for the fully implicit black-oil model.
Definition: blackoilboundaryratevector.hh:44
Implements a vector representing mass, molar or volumetric rates for the black oil model...
static void registerParameters()
Register all run-time parameters for the black-oil polymer module.
Definition: blackoilpolymermodules.hh:365
The multi-dimensional Newton method.
Definition: newtonmethod.hh:56
Represents the primary variables used by the black-oil model.
Calculates the local residual of the black oil model.
Definition: blackoillocalresidual.hh:43
This template class contains the data which is required to calculate the fluxes of the fluid phases o...
#define SET_PROP(EffTypeTagName, PropTagName)
Set a property for a specific type tag.
Definition: propertysystem.hh:297
The primary variable and equation indices for the black-oil model.
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition: multiphasebasemodel.hh:147
VTK output module for the black oil model&#39;s parameters.
The primary variable and equation indices for the black-oil model.
Definition: blackoilindices.hh:39
void finishInit()
Apply the initial conditions to the model.
Definition: blackoilmodel.hh:242
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
void updateMaxOilSaturations()
Update the maximum oil saturation observed during the simulation for all elements.
Definition: blackoilmodel.hh:497
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
Definition: multiphasebasemodel.hh:47
The primary variable and equation indices for the black-oil model.
Represents the primary variables used by the black-oil model.
Definition: blackoilprimaryvariables.hh:58
This template class contains the data which is required to calculate the fluxes of the fluid phases o...
Definition: blackoilextensivequantities.hh:52
VTK output module for the fluid composition.
Definition: vtkcompositionmodule.hh:77
#define SET_TYPE_PROP(EffTypeTagName, PropTagName,...)
Set a property which defines a type.
Definition: propertysystem.hh:377
std::string primaryVarName(int pvIdx) const
Given an primary variable index, return a human readable name.
Definition: blackoilmodel.hh:259
A newton solver which is specific to the black oil model.
Contains the classes required to extend the black-oil model by polymer.
Contains the quantities which are are constant within a finite volume in the black-oil model...
Definition: blackoilintensivequantities.hh:54