28 #ifndef EWOMS_MULTI_PHASE_BASE_MODEL_HH 29 #define EWOMS_MULTI_PHASE_BASE_MODEL_HH 31 #include <opm/material/densead/Math.hpp> 40 #include <opm/material/fluidmatrixinteractions/NullMaterial.hpp> 41 #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp> 42 #include <opm/material/heatconduction/DummyHeatConductionLaw.hpp> 43 #include <opm/common/Unused.hpp> 46 template <
class TypeTag>
51 namespace Properties {
83 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
84 typedef Opm::NullMaterialTraits<Scalar, FluidSystem::numPhases> Traits;
87 typedef Opm::NullMaterial<Traits> type;
101 Opm::DummyHeatConductionLaw<
typename GET_PROP_TYPE(TypeTag, Scalar)>);
106 HeatConductionLawParams,
119 template <
class TypeTag>
120 class MultiPhaseBaseModel :
public GET_PROP_TYPE(TypeTag, Discretization)
122 typedef typename GET_PROP_TYPE(TypeTag, Discretization) ParentType;
123 typedef typename GET_PROP_TYPE(TypeTag, Model) Implementation;
124 typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
125 typedef typename GET_PROP_TYPE(TypeTag, ThreadManager) ThreadManager;
126 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
127 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
128 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
129 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
130 typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
131 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
133 typedef typename GridView::template Codim<0>::Iterator ElementIterator;
134 typedef typename GridView::template Codim<0>::Entity Element;
137 enum { numComponents = FluidSystem::numComponents };
140 MultiPhaseBaseModel(Simulator& simulator)
141 : ParentType(simulator)
149 ParentType::registerParameters();
161 ParentType::finishInit();
166 Scalar minDofVolume = 1e100;
167 for (
unsigned globalDofIdx = 0; globalDofIdx < this->numGridDof(); ++ globalDofIdx)
168 if (this->isLocalDof(globalDofIdx))
169 minDofVolume = std::min(minDofVolume, this->dofTotalVolume(globalDofIdx));
170 minDofVolume = this->gridView().comm().min(minDofVolume);
172 Scalar newtonTolerance =
EWOMS_GET_PARAM(TypeTag, Scalar, NewtonRawTolerance);
173 newtonTolerance /= std::sqrt(minDofVolume);
174 this->newtonMethod().setTolerance(newtonTolerance);
194 assert(0 <= phaseIdx && phaseIdx < numPhases);
207 ElementContext elemCtx(this->simulator_);
208 ElementIterator elemIt = threadedElemIt.beginParallel();
211 for (; !threadedElemIt.isFinished(elemIt); elemIt = threadedElemIt.increment()) {
212 const Element& elem = *elemIt;
213 if (elem.partitionType() != Dune::InteriorEntity)
216 elemCtx.updateStencil(elem);
217 elemCtx.updateIntensiveQuantities(0);
219 const auto& stencil = elemCtx.stencil(0);
221 for (
unsigned dofIdx = 0; dofIdx < elemCtx.numDof(0); ++dofIdx) {
222 const auto& scv = stencil.subControlVolume(dofIdx);
223 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, 0);
226 this->localResidual(threadId).addPhaseStorage(tmp,
231 tmp *= scv.volume()*intQuants.extrusionFactor();
240 storage = this->gridView_.comm().sum(storage);
243 void registerOutputModules_()
245 ParentType::registerOutputModules_();
253 const Implementation& asImp_()
const 254 {
return *
static_cast<const Implementation *
>(
this); }
This file contains the necessary classes to calculate the velocity out of a pressure potential gradie...
bool phaseIsConsidered(unsigned phaseIdx OPM_UNUSED) const
Returns true iff a fluid phase is used by the model.
Definition: multiphasebasemodel.hh:182
#define SET_BOOL_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant boolean value.
Definition: propertysystem.hh:361
Specifies a flux module which uses the Darcy relation.
Definition: darcyfluxmodule.hh:65
Definition: baseauxiliarymodule.hh:37
The base class for the vertex centered finite volume discretization scheme.
Definition: vcfvdiscretization.hh:49
VTK output module for quantities which make sense for all models which deal with multiple fluid phase...
Definition: vtkmultiphasemodule.hh:99
Provides an STL-iterator like interface to iterate over the enties of a GridView in OpenMP threaded a...
Definition: threadedentityiterator.hh:41
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
Defines the common properties required by the porous medium multi-phase models.
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
#define SET_TAG_PROP(EffTypeTagName, PropTagName, ValueTypeTagName)
Define a property containing a type tag.
Definition: propertysystem.hh:436
#define INHERITS_FROM(...)
Syntactic sugar for NEW_TYPE_TAG.
Definition: propertysystem.hh:230
#define SET_INT_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant integer value.
Definition: propertysystem.hh:345
static unsigned threadId()
Return the index of the current OpenMP thread.
Definition: threadmanager.hh:119
The base class for the vertex centered finite volume discretization scheme.
void finishInit()
Apply the initial conditions to the model.
Definition: multiphasebasemodel.hh:159
#define SET_SPLICES(TypeTagName,...)
Define splices for a given type tag.
Definition: propertysystem.hh:213
VTK output module for the temperature in which assume thermal equilibrium.
Definition: vtktemperaturemodule.hh:59
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtktemperaturemodule.hh:82
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
#define SET_PROP(EffTypeTagName, PropTagName)
Set a property for a specific type tag.
Definition: propertysystem.hh:297
void globalPhaseStorage(EqVector &storage, unsigned phaseIdx)
Compute the total storage inside one phase of all conservation quantities.
Definition: multiphasebasemodel.hh:192
This class calculates the pressure potential gradients and the filter velocities for multi-phase flow...
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkmultiphasemodule.hh:137
static void registerParameters()
Register all run-time parameters for the immiscible model.
Definition: multiphasebasemodel.hh:147
A base class for fully-implicit multi-phase porous-media flow models which assume multiple fluid phas...
Definition: multiphasebasemodel.hh:47
The base class for the problems of ECFV discretizations which deal with a multi-phase flow through a ...
Definition: multiphasebaseproblem.hh:59
#define SET_TYPE_PROP(EffTypeTagName, PropTagName,...)
Set a property which defines a type.
Definition: propertysystem.hh:377
The base class for the problems of ECFV discretizations which deal with a multi-phase flow through a ...