26 #ifndef EWOMS_DIFFUSION_MODULE_HH
27 #define EWOMS_DIFFUSION_MODULE_HH
32 #include <dune/common/fvector.hh>
35 namespace Properties {
45 template <
class TypeTag,
bool enableDiffusion>
51 template <
class TypeTag>
55 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
56 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
57 typedef typename FluidSystem::ParameterCache ParameterCache;
70 template <
class Context>
72 int spaceIdx,
int timeIdx)
79 template <
class TypeTag>
83 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
84 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
85 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
88 enum { numPhases = FluidSystem::numPhases };
89 enum { numComponents = FluidSystem::numComponents };
90 enum { conti0EqIdx = Indices::conti0EqIdx };
92 typedef Opm::MathToolbox<Evaluation> Toolbox;
105 template <
class Context>
107 int spaceIdx,
int timeIdx)
109 const auto &extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
111 const auto &fluidStateI = context.intensiveQuantities(extQuants.interiorIndex(), timeIdx).fluidState();
112 const auto &fluidStateJ = context.intensiveQuantities(extQuants.exteriorIndex(), timeIdx).fluidState();
114 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
116 Evaluation rhoMolar = fluidStateI.molarDensity(phaseIdx);
117 rhoMolar += Toolbox::value(fluidStateJ.molarDensity(phaseIdx));
120 for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
122 flux[conti0EqIdx + compIdx] +=
124 * extQuants.moleFractionGradientNormal(phaseIdx, compIdx)
125 * extQuants.effectiveDiffusionCoefficient(phaseIdx, compIdx);
137 template <
class TypeTag,
bool enableDiffusion>
143 template <
class TypeTag>
147 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
148 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
150 typedef typename FluidSystem::ParameterCache ParameterCache;
159 OPM_THROW(std::logic_error,
"Method tortuosity() does not make sense "
160 "if diffusion is disabled");
169 OPM_THROW(std::logic_error,
"Method diffusionCoefficient() does not "
170 "make sense if diffusion is disabled");
179 OPM_THROW(std::logic_error,
"Method effectiveDiffusionCoefficient() "
180 "does not make sense if diffusion is "
189 template <
class Flu
idState>
191 ParameterCache ¶mCache,
192 const ElementContext &elemCtx,
201 template <
class TypeTag>
205 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
206 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
208 typedef typename FluidSystem::ParameterCache ParameterCache;
210 enum { numPhases = FluidSystem::numPhases };
211 enum { numComponents = FluidSystem::numComponents };
219 {
return diffusionCoefficient_[phaseIdx][compIdx]; }
226 {
return tortuosity_[phaseIdx]; }
233 {
return tortuosity_[phaseIdx] * diffusionCoefficient_[phaseIdx][compIdx]; }
240 template <
class Flu
idState>
242 ParameterCache ¶mCache,
243 const ElementContext &elemCtx,
247 const auto &intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
249 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
250 if (!elemCtx.model().phaseIsConsidered(phaseIdx))
256 tortuosity_[phaseIdx] =
257 1.0 / (intQuants.porosity() * intQuants.porosity())
258 * std::pow(std::max(0.0001,
260 * intQuants.fluidState().saturation(phaseIdx)),
263 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
264 diffusionCoefficient_[phaseIdx][compIdx] =
265 FluidSystem::diffusionCoefficient(fluidState,
274 Scalar tortuosity_[numPhases];
275 Scalar diffusionCoefficient_[numPhases][numComponents];
284 template <
class TypeTag,
bool enableDiffusion>
290 template <
class TypeTag>
294 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
301 void update_(
const ElementContext &elemCtx,
int faceIdx,
int timeIdx)
304 template <
class Context,
class Flu
idState>
306 const FluidState &fluidState)
318 OPM_THROW(std::logic_error,
"Method moleFractionGradient() does not "
319 "make sense if diffusion is disabled.");
331 OPM_THROW(std::logic_error,
"Method effectiveDiffusionCoefficient() "
332 "does not make sense if diffusion is "
340 template <
class TypeTag>
344 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
347 enum { dimWorld = GridView::dimensionworld };
351 typedef Dune::FieldVector<Scalar, dimWorld> DimVector;
358 void update_(
const ElementContext &elemCtx,
int faceIdx,
int timeIdx)
360 const auto& gradCalc = elemCtx.gradientCalculator();
363 DimVector temperatureGrad;
365 const auto &face = elemCtx.stencil(timeIdx).interiorFace(faceIdx);
366 const auto &extQuants = elemCtx.extensiveQuantities(faceIdx, timeIdx);
368 const auto &intQuantsInside = elemCtx.intensiveQuantities(extQuants.interiorIndex(), timeIdx);
369 const auto &intQuantsOutside = elemCtx.intensiveQuantities(extQuants.exteriorIndex(), timeIdx);
371 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
372 if (!elemCtx.model().phaseIsConsidered(phaseIdx))
376 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
379 DimVector moleFractionGradient(0.0);
380 gradCalc.calculateGradient(moleFractionGradient,
383 moleFractionCallback);
385 moleFractionGradientNormal_[phaseIdx][compIdx] =
386 face.normal()*moleFractionGradient;
387 Valgrind::CheckDefined(moleFractionGradientNormal_[phaseIdx][compIdx]);
391 effectiveDiffusionCoefficient_[phaseIdx][compIdx] =
392 (intQuantsInside.effectiveDiffusionCoefficient(phaseIdx, compIdx)
394 intQuantsOutside.effectiveDiffusionCoefficient(phaseIdx, compIdx))
397 Valgrind::CheckDefined(effectiveDiffusionCoefficient_[phaseIdx][compIdx]);
402 template <
class Context,
class Flu
idState>
404 const FluidState &fluidState)
406 const auto &stencil = context.stencil(timeIdx);
407 const auto &face = stencil.boundaryFace(bfIdx);
409 const auto &elemCtx = context.elementContext();
410 int insideScvIdx = face.interiorIndex();
411 const auto &insideScv = stencil.subControlVolume(insideScvIdx);
413 const auto &intQuantsInside = elemCtx.intensiveQuantities(insideScvIdx, timeIdx);
414 const auto &fluidStateInside = intQuantsInside.fluidState();
417 DimVector distVec = face.integrationPos();
418 distVec -= context.element().geometry().global(insideScv.localGeometry().center());
420 Scalar dist = distVec * face.normal();
426 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
427 if (!elemCtx.model().phaseIsConsidered(phaseIdx))
430 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
431 DimVector moleFractionGradient(0.0);
435 moleFractionGradientNormal_[phaseIdx][compIdx] =
436 (fluidState.moleFraction(phaseIdx, compIdx)
438 fluidStateInside.moleFraction(phaseIdx, compIdx))
440 Valgrind::CheckDefined(moleFractionGradientNormal_[phaseIdx][compIdx]);
444 effectiveDiffusionCoefficient_[phaseIdx][compIdx] =
445 intQuantsInside.effectiveDiffusionCoefficient(phaseIdx, compIdx);
447 Valgrind::CheckDefined(effectiveDiffusionCoefficient_[phaseIdx][compIdx]);
460 {
return moleFractionGradientNormal_[phaseIdx][compIdx]; }
470 {
return effectiveDiffusionCoefficient_[phaseIdx][compIdx]; }
473 Scalar moleFractionGradientNormal_[numPhases][numComponents];
474 Scalar effectiveDiffusionCoefficient_[numPhases][numComponents];
static void registerParameters()
Register all run-time parameters for the diffusion module.
Definition: diffusionmodule.hh:98
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: diffusionmodule.hh:46
Scalar diffusionCoefficient(int phaseIdx, int compIdx) const
Returns the molecular diffusion coefficient for a component in a phase.
Definition: diffusionmodule.hh:167
Declare the properties used by the infrastructure code of the finite volume discretizations.
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
Scalar effectiveDiffusionCoefficient(int phaseIdx, int compIdx) const
The effective diffusion coeffcient of a component in a fluid phase at the face's integration point...
Definition: diffusionmodule.hh:469
void update_(const ElementContext &elemCtx, int faceIdx, int timeIdx)
Update the quantities required to calculate the diffusive mass fluxes.
Definition: diffusionmodule.hh:358
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:485
void update_(const ElementContext &elemCtx, int faceIdx, int timeIdx)
Update the quantities required to calculate the diffusive mass fluxes.
Definition: diffusionmodule.hh:301
void setPhaseIndex(short phaseIdx)
Set the index of the fluid phase for which a mole fraction should be returned.
Definition: quantitycallbacks.hh:429
static void addDiffusiveFlux(RateVector &flux, const Context &context, int spaceIdx, int timeIdx)
Adds the diffusive mass flux flux to the flux vector over a flux integration point.
Definition: diffusionmodule.hh:71
void updateBoundary_(const Context &context, int bfIdx, int timeIdx, const FluidState &fluidState)
Definition: diffusionmodule.hh:305
Scalar moleFractionGradientNormal(int phaseIdx, int compIdx) const
The the gradient of the mole fraction times the face normal.
Definition: diffusionmodule.hh:316
Callback class for a mole fraction of a component in a phase.
Definition: quantitycallbacks.hh:402
Scalar moleFractionGradientNormal(int phaseIdx, int compIdx) const
The the gradient of the mole fraction times the face normal.
Definition: diffusionmodule.hh:459
void updateBoundary_(const Context &context, int bfIdx, int timeIdx, const FluidState &fluidState)
Definition: diffusionmodule.hh:403
NEW_PROP_TAG(Grid)
The type of the DUNE grid.
Scalar effectiveDiffusionCoefficient(int phaseIdx, int compIdx) const
Returns the effective molecular diffusion coefficient of the porous medium for a component in a phase...
Definition: diffusionmodule.hh:177
void update_(FluidState &fs, ParameterCache ¶mCache, const ElementContext &elemCtx, int dofIdx, int timeIdx)
Update the quantities required to calculate diffusive mass fluxes.
Definition: diffusionmodule.hh:190
Scalar tortuosity(int phaseIdx) const
Returns the tortuousity of the sub-domain of a fluid phase in the porous medium.
Definition: diffusionmodule.hh:157
Definition: baseauxiliarymodule.hh:35
void update_(FluidState &fluidState, ParameterCache ¶mCache, const ElementContext &elemCtx, int dofIdx, int timeIdx)
Update the quantities required to calculate diffusive mass fluxes.
Definition: diffusionmodule.hh:241
Provides the quantities required to calculate diffusive mass fluxes.
Definition: diffusionmodule.hh:285
This method contains all callback classes for quantities that are required by some extensive quantiti...
static void registerParameters()
Register all run-time parameters for the diffusion module.
Definition: diffusionmodule.hh:63
Provides the volumetric quantities required for the calculation of molecular diffusive fluxes...
Definition: diffusionmodule.hh:138
Scalar effectiveDiffusionCoefficient(int phaseIdx, int compIdx) const
Returns the effective molecular diffusion coefficient of the porous medium for a component in a phase...
Definition: diffusionmodule.hh:232
Scalar tortuosity(int phaseIdx) const
Returns the tortuousity of the sub-domain of a fluid phase in the porous medium.
Definition: diffusionmodule.hh:225
Scalar effectiveDiffusionCoefficient(int phaseIdx, int compIdx) const
The effective diffusion coeffcient of a component in a fluid phase at the face's integration point...
Definition: diffusionmodule.hh:329
Scalar diffusionCoefficient(int phaseIdx, int compIdx) const
Returns the molecular diffusion coefficient for a component in a phase.
Definition: diffusionmodule.hh:218
static void addDiffusiveFlux(RateVector &flux, const Context &context, int spaceIdx, int timeIdx)
Adds the mass flux due to molecular diffusion to the flux vector over the flux integration point...
Definition: diffusionmodule.hh:106
void setComponentIndex(short compIdx)
Set the index of the component for which the mole fraction should be returned.
Definition: quantitycallbacks.hh:436