13#ifndef DUMUX_DISCRETIZATION_STAGGERED_FICKS_LAW_HH 
   14#define DUMUX_DISCRETIZATION_STAGGERED_FICKS_LAW_HH 
   17#include <dune/common/fvector.hh> 
   33template<
class TypeTag, 
class DiscretizationMethod, ReferenceSystemFormulation referenceSystem>
 
   40template <
class TypeTag, ReferenceSystemFormulation referenceSystem>
 
   45    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   46    using SubControlVolumeFace = 
typename GridGeometry::SubControlVolumeFace;
 
   49    using GridView = 
typename GridGeometry::GridView;
 
   50    using Element = 
typename GridView::template Codim<0>::Entity;
 
   52    using Indices = 
typename ModelTraits::Indices;
 
   55    static constexpr int numComponents = ModelTraits::numFluidComponents();
 
   56    static constexpr int numPhases = ModelTraits::numFluidPhases();
 
   58    using NumEqVector = Dune::FieldVector<Scalar, numComponents>;
 
   60    static_assert(ModelTraits::numFluidPhases() == 1, 
"Only one phase supported!");
 
   69    { 
return referenceSystem; }
 
 
   83    template<
class Problem, 
class ElementVolumeVariables>
 
   84    static NumEqVector 
flux(
const Problem& problem,
 
   85                            const Element& element,
 
   86                            const FVElementGeometry& fvGeometry,
 
   87                            const ElementVolumeVariables& elemVolVars,
 
   88                            const SubControlVolumeFace &scvf)
 
   90        NumEqVector 
flux(0.0);
 
   95        if (scvf.boundary() && problem.boundaryTypes(element, scvf).isOutflow(Indices::conti0EqIdx + 1))
 
   98        const int phaseIdx = 0;
 
  100        const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
 
  101        const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
 
  102        const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
 
  104        const Scalar insideDistance = (insideScv.dofPosition() - scvf.ipGlobal()).two_norm();
 
  105        const Scalar insideDensity = 
massOrMolarDensity(insideVolVars, referenceSystem, phaseIdx);
 
  107        for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
 
  109            if (compIdx == FluidSystem::getMainComponent(phaseIdx))
 
  112            const Scalar massOrMoleFractionInside = 
massOrMoleFraction(insideVolVars, referenceSystem, phaseIdx, compIdx);
 
  113            const Scalar massOrMoleFractionOutside =  
massOrMoleFraction(outsideVolVars, referenceSystem, phaseIdx, compIdx);
 
  114            const Scalar insideDiffCoeff = getEffectiveDiffusionCoefficient_(insideVolVars, phaseIdx, compIdx) * insideVolVars.extrusionFactor();
 
  119                const Scalar outsideDensity = 
massOrMolarDensity(outsideVolVars, referenceSystem, phaseIdx);
 
  120                const Scalar avgDensity = 0.5*(insideDensity + outsideDensity);
 
  121                flux[compIdx] = avgDensity * insideDiffCoeff
 
  122                                * (massOrMoleFractionInside - massOrMoleFractionOutside) / insideDistance;
 
  126                const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
 
  127                const Scalar outsideDiffCoeff = getEffectiveDiffusionCoefficient_(outsideVolVars, phaseIdx, compIdx)
 
  128                                      * outsideVolVars.extrusionFactor();
 
  129                const Scalar outsideDistance = (outsideScv.dofPosition() - scvf.ipGlobal()).two_norm();
 
  130                const Scalar outsideDensity = 
massOrMolarDensity(outsideVolVars, referenceSystem, phaseIdx);
 
  132                const Scalar avgDensity = 0.5*(insideDensity + outsideDensity);
 
  133                const Scalar avgDiffCoeff = 
harmonicMean(insideDiffCoeff, outsideDiffCoeff, insideDistance, outsideDistance);
 
  135                flux[compIdx] = avgDensity * avgDiffCoeff
 
  136                                * (massOrMoleFractionInside - massOrMoleFractionOutside) / (insideDistance + outsideDistance);
 
  141        const Scalar cumulativeFlux = std::accumulate(
flux.begin(), 
flux.end(), 0.0);
 
  142        flux[FluidSystem::getMainComponent(0)] = -cumulativeFlux;
 
  144        flux *= Extrusion::area(fvGeometry, scvf);
 
 
  150    static Scalar getEffectiveDiffusionCoefficient_(
const VolumeVariables& volVars, 
const int phaseIdx, 
const int compIdx)
 
  152        return volVars.effectiveDiffusionCoefficient(phaseIdx, FluidSystem::getMainComponent(phaseIdx), compIdx);
 
 
Container storing the diffusion coefficients required by Fick's law. Uses the minimal possible contai...
Definition fickiandiffusioncoefficients.hh:32
DiscretizationMethods::Staggered DiscretizationMethod
Definition flux/staggered/freeflow/fickslaw.hh:63
static NumEqVector flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf)
Returns the diffusive fluxes of all components within a fluid phase across the given sub-control volu...
Definition flux/staggered/freeflow/fickslaw.hh:84
FickianDiffusionCoefficients< Scalar, numPhases, numComponents > DiffusionCoefficientsContainer
Definition flux/staggered/freeflow/fickslaw.hh:75
static constexpr ReferenceSystemFormulation referenceSystemFormulation()
Definition flux/staggered/freeflow/fickslaw.hh:68
static constexpr DiscretizationMethod discMethod
Definition flux/staggered/freeflow/fickslaw.hh:65
FluxVariablesCaching::EmptyDiffusionCache Cache
Definition flux/staggered/freeflow/fickslaw.hh:73
forward declaration of the method-specific implementation
Definition flux/box/fickslaw.hh:32
Defines all properties used in Dumux.
Helper classes to compute the integration elements.
Container storing the diffusion coefficients required by Fick's law. Uses the minimal possible contai...
Classes related to flux variables caching.
constexpr Scalar harmonicMean(Scalar x, Scalar y, Scalar wx=1.0, Scalar wy=1.0) noexcept
Calculate the (weighted) harmonic mean of two scalar values.
Definition math.hh:57
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:54
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:43
ReferenceSystemFormulation
The formulations available for Fick's law related to the reference system.
Definition referencesystemformulation.hh:33
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
Define some often used mathematical functions.
The available discretization methods in Dumux.
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition extrusion.hh:166
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Definition fluxvariablescaching.hh:55