28#ifndef OPM_DIFFUSION_MODULE_HH
29#define OPM_DIFFUSION_MODULE_HH
34#include <opm/material/common/Valgrind.hpp>
36#include <dune/common/fvector.hh>
46template <
class TypeTag,
bool enableDiffusion>
52template <
class TypeTag>
70 template <
class Context>
81template <
class TypeTag>
90 enum { numPhases = FluidSystem::numPhases };
91 enum { numComponents = FluidSystem::numComponents };
92 enum { conti0EqIdx = Indices::conti0EqIdx };
94 using Toolbox = Opm::MathToolbox<Evaluation>;
107 template <
class Context>
109 unsigned spaceIdx,
unsigned timeIdx)
111 const auto& extQuants = context.extensiveQuantities(spaceIdx, timeIdx);
113 const auto& fluidStateI = context.intensiveQuantities(extQuants.interiorIndex(), timeIdx).fluidState();
114 const auto& fluidStateJ = context.intensiveQuantities(extQuants.exteriorIndex(), timeIdx).fluidState();
116 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
118 Evaluation rhoMolar = fluidStateI.molarDensity(phaseIdx);
119 rhoMolar += Toolbox::value(fluidStateJ.molarDensity(phaseIdx));
122 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
124 flux[conti0EqIdx + compIdx] +=
126 * extQuants.moleFractionGradientNormal(phaseIdx, compIdx)
127 * extQuants.effectiveDiffusionCoefficient(phaseIdx, compIdx);
139template <
class TypeTag,
bool enableDiffusion>
145template <
class TypeTag>
159 throw std::logic_error(
"Method tortuosity() does not make sense "
160 "if diffusion is disabled");
169 throw std::logic_error(
"Method diffusionCoefficient() does not "
170 "make sense if diffusion is disabled");
179 throw std::logic_error(
"Method effectiveDiffusionCoefficient() "
180 "does not make sense if diffusion is disabled");
188 template <
class Flu
idState>
190 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>&,
191 const ElementContext&,
200template <
class TypeTag>
208 enum { numPhases = FluidSystem::numPhases };
209 enum { numComponents = FluidSystem::numComponents };
217 {
return diffusionCoefficient_[phaseIdx][compIdx]; }
224 {
return tortuosity_[phaseIdx]; }
231 {
return tortuosity_[phaseIdx] * diffusionCoefficient_[phaseIdx][compIdx]; }
238 template <
class Flu
idState>
240 typename FluidSystem::template ParameterCache<typename FluidState::Scalar>& paramCache,
241 const ElementContext& elemCtx,
245 using Toolbox = Opm::MathToolbox<Evaluation>;
247 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
248 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
249 if (!elemCtx.model().phaseIsConsidered(phaseIdx))
256 const Evaluation& base =
259 * intQuants.fluidState().saturation(phaseIdx));
260 tortuosity_[phaseIdx] =
261 1.0 / (intQuants.porosity() * intQuants.porosity())
262 * Toolbox::pow(base, 10.0/3.0);
264 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
265 diffusionCoefficient_[phaseIdx][compIdx] =
266 FluidSystem::diffusionCoefficient(fluidState,
275 Evaluation tortuosity_[numPhases];
276 Evaluation diffusionCoefficient_[numPhases][numComponents];
285template <
class TypeTag,
bool enableDiffusion>
291template <
class TypeTag>
308 template <
class Context,
class Flu
idState>
325 throw std::logic_error(
"The method moleFractionGradient() does not "
326 "make sense if diffusion is disabled.");
339 throw std::logic_error(
"The method effectiveDiffusionCoefficient() "
340 "does not make sense if diffusion is disabled.");
347template <
class TypeTag>
355 enum { dimWorld = GridView::dimensionworld };
356 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
357 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
359 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
360 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
367 void update_(
const ElementContext& elemCtx,
unsigned faceIdx,
unsigned timeIdx)
369 const auto& gradCalc = elemCtx.gradientCalculator();
372 const auto& face = elemCtx.stencil(timeIdx).interiorFace(faceIdx);
373 const auto& normal = face.normal();
374 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, timeIdx);
376 const auto& intQuantsInside = elemCtx.intensiveQuantities(extQuants.interiorIndex(), timeIdx);
377 const auto& intQuantsOutside = elemCtx.intensiveQuantities(extQuants.exteriorIndex(), timeIdx);
379 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
380 if (!elemCtx.model().phaseIsConsidered(phaseIdx))
384 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
387 DimEvalVector moleFractionGradient(0.0);
388 gradCalc.calculateGradient(moleFractionGradient,
391 moleFractionCallback);
393 moleFractionGradientNormal_[phaseIdx][compIdx] = 0.0;
394 for (
unsigned i = 0; i < normal.size(); ++i)
395 moleFractionGradientNormal_[phaseIdx][compIdx] +=
396 normal[i]*moleFractionGradient[i];
397 Opm::Valgrind::CheckDefined(moleFractionGradientNormal_[phaseIdx][compIdx]);
401 effectiveDiffusionCoefficient_[phaseIdx][compIdx] =
402 (intQuantsInside.effectiveDiffusionCoefficient(phaseIdx, compIdx)
404 intQuantsOutside.effectiveDiffusionCoefficient(phaseIdx, compIdx))
407 Opm::Valgrind::CheckDefined(effectiveDiffusionCoefficient_[phaseIdx][compIdx]);
412 template <
class Context,
class Flu
idState>
416 const FluidState& fluidState)
418 const auto& stencil = context.stencil(timeIdx);
419 const auto& face = stencil.boundaryFace(bfIdx);
421 const auto& elemCtx = context.elementContext();
422 unsigned insideScvIdx = face.interiorIndex();
423 const auto& insideScv = stencil.subControlVolume(insideScvIdx);
425 const auto& intQuantsInside = elemCtx.intensiveQuantities(insideScvIdx, timeIdx);
426 const auto& fluidStateInside = intQuantsInside.fluidState();
429 DimVector distVec = face.integrationPos();
430 distVec -= context.element().geometry().global(insideScv.localGeometry().center());
432 Scalar dist = distVec * face.normal();
438 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
439 if (!elemCtx.model().phaseIsConsidered(phaseIdx))
442 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
445 moleFractionGradientNormal_[phaseIdx][compIdx] =
446 (fluidState.moleFraction(phaseIdx, compIdx)
448 fluidStateInside.moleFraction(phaseIdx, compIdx))
450 Opm::Valgrind::CheckDefined(moleFractionGradientNormal_[phaseIdx][compIdx]);
454 effectiveDiffusionCoefficient_[phaseIdx][compIdx] =
455 intQuantsInside.effectiveDiffusionCoefficient(phaseIdx, compIdx);
457 Opm::Valgrind::CheckDefined(effectiveDiffusionCoefficient_[phaseIdx][compIdx]);
470 {
return moleFractionGradientNormal_[phaseIdx][compIdx]; }
480 {
return effectiveDiffusionCoefficient_[phaseIdx][compIdx]; }
483 Evaluation moleFractionGradientNormal_[numPhases][numComponents];
484 Evaluation effectiveDiffusionCoefficient_[numPhases][numComponents];
void update_(const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate the diffusive mass fluxes.
Definition: diffusionmodule.hh:303
const Evaluation & effectiveDiffusionCoefficient(unsigned, unsigned) const
The effective diffusion coeffcient of a component in a fluid phase at the face's integration point.
Definition: diffusionmodule.hh:336
const Evaluation & moleFractionGradientNormal(unsigned, unsigned) const
The the gradient of the mole fraction times the face normal.
Definition: diffusionmodule.hh:322
void updateBoundary_(const Context &, unsigned, unsigned, const FluidState &)
Definition: diffusionmodule.hh:309
void updateBoundary_(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Definition: diffusionmodule.hh:413
const Evaluation & moleFractionGradientNormal(unsigned phaseIdx, unsigned compIdx) const
The the gradient of the mole fraction times the face normal.
Definition: diffusionmodule.hh:469
const Evaluation & effectiveDiffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
The effective diffusion coeffcient of a component in a fluid phase at the face's integration point.
Definition: diffusionmodule.hh:479
void update_(const ElementContext &elemCtx, unsigned faceIdx, unsigned timeIdx)
Update the quantities required to calculate the diffusive mass fluxes.
Definition: diffusionmodule.hh:367
Provides the quantities required to calculate diffusive mass fluxes.
Definition: diffusionmodule.hh:286
Scalar diffusionCoefficient(unsigned, unsigned) const
Returns the molecular diffusion coefficient for a component in a phase.
Definition: diffusionmodule.hh:167
Scalar effectiveDiffusionCoefficient(unsigned, unsigned) const
Returns the effective molecular diffusion coefficient of the porous medium for a component in a phase...
Definition: diffusionmodule.hh:177
void update_(FluidState &, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > &, const ElementContext &, unsigned, unsigned)
Update the quantities required to calculate diffusive mass fluxes.
Definition: diffusionmodule.hh:189
Scalar tortuosity(unsigned) const
Returns the tortuousity of the sub-domain of a fluid phase in the porous medium.
Definition: diffusionmodule.hh:157
Evaluation tortuosity(unsigned phaseIdx) const
Returns the tortuousity of the sub-domain of a fluid phase in the porous medium.
Definition: diffusionmodule.hh:223
void update_(FluidState &fluidState, typename FluidSystem::template ParameterCache< typename FluidState::Scalar > ¶mCache, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the quantities required to calculate diffusive mass fluxes.
Definition: diffusionmodule.hh:239
Evaluation effectiveDiffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
Returns the effective molecular diffusion coefficient of the porous medium for a component in a phase...
Definition: diffusionmodule.hh:230
Evaluation diffusionCoefficient(unsigned phaseIdx, unsigned compIdx) const
Returns the molecular diffusion coefficient for a component in a phase.
Definition: diffusionmodule.hh:216
Provides the volumetric quantities required for the calculation of molecular diffusive fluxes.
Definition: diffusionmodule.hh:140
static void addDiffusiveFlux(RateVector &, const Context &, unsigned, unsigned)
Adds the diffusive mass flux flux to the flux vector over a flux integration point.
Definition: diffusionmodule.hh:71
static void registerParameters()
Register all run-time parameters for the diffusion module.
Definition: diffusionmodule.hh:63
static void addDiffusiveFlux(RateVector &flux, const Context &context, unsigned spaceIdx, unsigned timeIdx)
Adds the mass flux due to molecular diffusion to the flux vector over the flux integration point.
Definition: diffusionmodule.hh:108
static void registerParameters()
Register all run-time parameters for the diffusion module.
Definition: diffusionmodule.hh:100
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: diffusionmodule.hh:47
Callback class for a mole fraction of a component in a phase.
Definition: quantitycallbacks.hh:425
void setComponentIndex(unsigned compIdx)
Set the index of the component for which the mole fraction should be returned.
Definition: quantitycallbacks.hh:460
void setPhaseIndex(unsigned phaseIdx)
Set the index of the fluid phase for which a mole fraction should be returned.
Definition: quantitycallbacks.hh:453
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:242
This method contains all callback classes for quantities that are required by some extensive quantiti...