25#ifndef DUMUX_DISCRETIZATION_BOX_FICKS_LAW_HH
26#define DUMUX_DISCRETIZATION_BOX_FICKS_LAW_HH
28#include <dune/common/fvector.hh>
29#include <dune/common/fmatrix.hh>
43template<
class TypeTag,
class DiscretizationMethod, ReferenceSystemFormulation referenceSystem>
50template <
class TypeTag, ReferenceSystemFormulation referenceSystem>
58 using FVElementGeometry =
typename GridGeometry::LocalView;
59 using SubControlVolume =
typename GridGeometry::SubControlVolume;
60 using SubControlVolumeFace =
typename GridGeometry::SubControlVolumeFace;
67 using Element =
typename GridView::template Codim<0>::Entity;
70 enum { dim = GridView::dimension} ;
71 enum { dimWorld = GridView::dimensionworld} ;
74 numPhases = ModelTraits::numFluidPhases(),
75 numComponents = ModelTraits::numFluidComponents()
77 using DimWorldMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
78 using ComponentFluxVector = Dune::FieldVector<Scalar, numComponents>;
85 {
return referenceSystem; }
93 static ComponentFluxVector
flux(
const Problem& problem,
94 const Element& element,
95 const FVElementGeometry& fvGeometry,
96 const ElementVolumeVariables& elemVolVars,
97 const SubControlVolumeFace& scvf,
99 const ElementFluxVariablesCache& elemFluxVarsCache)
102 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
103 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
106 const auto& fluxVarsCache = elemFluxVarsCache[scvf];
107 const auto& shapeValues = fluxVarsCache.shapeValues();
111 for (
auto&& scv : scvs(fvGeometry))
112 rho +=
massOrMolarDensity(elemVolVars[scv], referenceSystem, phaseIdx)*shapeValues[scv.indexInElement()][0];
114 ComponentFluxVector componentFlux(0.0);
115 for (
int compIdx = 0; compIdx < numComponents; compIdx++)
117 if constexpr (!FluidSystem::isTracerFluidSystem())
118 if (compIdx == FluidSystem::getMainComponent(phaseIdx))
121 const auto D = averageDiffusionCoefficient_(phaseIdx, compIdx, insideVolVars, outsideVolVars, problem, scvf);
124 const auto massOrMoleFrac = [&](
const SubControlVolume& scv){
return massOrMoleFraction(elemVolVars[scv], referenceSystem, phaseIdx, compIdx); };
125 componentFlux[compIdx] = discreteFlux_(fvGeometry, scvf, fluxVarsCache, massOrMoleFrac, D, rho);
128 if constexpr (!FluidSystem::isTracerFluidSystem())
129 if (BalanceEqOpts::mainComponentIsBalanced(phaseIdx))
130 componentFlux[FluidSystem::getMainComponent(phaseIdx)] -= componentFlux[compIdx];
132 return componentFlux;
136 static std::array<std::vector<Scalar>, numComponents>
138 const Element& element,
139 const FVElementGeometry& fvGeometry,
140 const ElementVolumeVariables& elemVolVars,
141 const SubControlVolumeFace& scvf,
142 const FluxVarCache& fluxVarCache,
146 const auto& shapeValues = fluxVarCache.shapeValues();
147 for (
auto&& scv : scvs(fvGeometry))
148 rho +=
massOrMolarDensity(elemVolVars[scv], referenceSystem, phaseIdx)*shapeValues[scv.indexInElement()][0];
150 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
151 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
153 std::array<std::vector<Scalar>, numComponents> ti;
154 for (
int compIdx = 0; compIdx < numComponents; compIdx++)
156 if constexpr (!FluidSystem::isTracerFluidSystem())
157 if (compIdx == FluidSystem::getMainComponent(phaseIdx))
160 const auto D = averageDiffusionCoefficient_(phaseIdx, compIdx, insideVolVars, outsideVolVars, problem, scvf);
162 ti[compIdx].resize(fvGeometry.numScv());
163 for (
auto&& scv : scvs(fvGeometry))
164 ti[compIdx][scv.indexInElement()] = -rho*
vtmv(scvf.unitOuterNormal(), D, fluxVarCache.gradN(scv.indexInElement()))*Extrusion::area(scvf);
171 static Scalar averageDiffusionCoefficient_(
const int phaseIdx,
const int compIdx,
172 const VolumeVariables& insideVV,
const VolumeVariables& outsideVV,
173 const Problem& problem,
174 const SubControlVolumeFace& scvf)
177 auto [insideD, outsideD] = diffusionCoefficientsAtInterface_(phaseIdx, compIdx, insideVV, outsideVV);
180 insideD *= insideVV.extrusionFactor();
181 outsideD *= outsideVV.extrusionFactor();
187 static std::pair<Scalar, Scalar>
188 diffusionCoefficientsAtInterface_([[maybe_unused]]
const int phaseIdx,
const int compIdx,
189 const VolumeVariables& insideVV,
const VolumeVariables& outsideVV)
191 if constexpr (!FluidSystem::isTracerFluidSystem())
193 const auto mainCompIdx = FluidSystem::getMainComponent(phaseIdx);
194 const auto insideD = insideVV.effectiveDiffusionCoefficient(phaseIdx, mainCompIdx, compIdx);
195 const auto outsideD = outsideVV.effectiveDiffusionCoefficient(phaseIdx, mainCompIdx, compIdx);
196 return { std::move(insideD), std::move(outsideD) };
200 const auto insideD = insideVV.effectiveDiffusionCoefficient(0, 0, compIdx);
201 const auto outsideD = outsideVV.effectiveDiffusionCoefficient(0, 0, compIdx);
202 return { std::move(insideD), std::move(outsideD) };
206 template<
class EvaluateVariable,
class Tensor>
207 static Scalar discreteFlux_(
const FVElementGeometry& fvGeometry,
208 const SubControlVolumeFace& scvf,
209 const FluxVarCache& fluxVarsCache,
211 const Tensor& D,
const Scalar preFactor)
213 Dune::FieldVector<Scalar, dimWorld> gradX(0.0);
214 for (
auto&& scv : scvs(fvGeometry))
216 return -1.0*preFactor*
vtmv(scvf.unitOuterNormal(), D, gradX)*Extrusion::area(scvf);
Define some often used mathematical functions.
Helper classes to compute the integration elements.
The available discretization methods in Dumux.
The reference frameworks and formulations available for splitting total fluxes into a advective and d...
Container storing the diffusion coefficients required by Fick's law. Uses the minimal possible contai...
A free function to average a Tensor at an interface.
VolumeVariables::PrimaryVariables::value_type massOrMoleFraction(const VolumeVariables &volVars, ReferenceSystemFormulation referenceSys, const int phaseIdx, const int compIdx)
returns the mass or mole fraction to be used in Fick's law based on the reference system
Definition: referencesystemformulation.hh:66
VolumeVariables::PrimaryVariables::value_type massOrMolarDensity(const VolumeVariables &volVars, ReferenceSystemFormulation referenceSys, const int phaseIdx)
evaluates the density to be used in Fick's law based on the reference system
Definition: referencesystemformulation.hh:55
ReferenceSystemFormulation
The formulations available for Fick's law related to the reference system.
Definition: referencesystemformulation.hh:45
Dune::DenseMatrix< MAT >::value_type vtmv(const Dune::DenseVector< V1 > &v1, const Dune::DenseMatrix< MAT > &M, const Dune::DenseVector< V2 > &v2)
Evaluates the scalar product of a vector v2, projected by a matrix M, with a vector v1.
Definition: math.hh:863
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition: extrusion.hh:177
Scalar faceTensorAverage(const Scalar T1, const Scalar T2, const Dune::FieldVector< Scalar, dim > &normal)
Average of a discontinuous scalar field at discontinuity interface (for compatibility reasons with th...
Definition: facetensoraverage.hh:41
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:150
forward declaration of the method-specific implemetation
Definition: flux/box/fickslaw.hh:44
static std::array< std::vector< Scalar >, numComponents > calculateTransmissibilities(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const FluxVarCache &fluxVarCache, const int phaseIdx)
Definition: flux/box/fickslaw.hh:137
static constexpr ReferenceSystemFormulation referenceSystemFormulation()
Definition: flux/box/fickslaw.hh:84
static ComponentFluxVector flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const int phaseIdx, const ElementFluxVariablesCache &elemFluxVarsCache)
Returns the diffusive fluxes of all components within a fluid phase across the given sub-control volu...
Definition: flux/box/fickslaw.hh:93
Container storing the diffusion coefficients required by Fick's law. Uses the minimal possible contai...
Definition: fickiandiffusioncoefficients.hh:44
Declares all properties used in Dumux.