13#ifndef DUMUX_DISCRETIZATION_BOX_DISPERSION_FLUX_HH 
   14#define DUMUX_DISCRETIZATION_BOX_DISPERSION_FLUX_HH 
   16#include <dune/common/fvector.hh> 
   17#include <dune/common/fmatrix.hh> 
   18#include <dune/common/float_cmp.hh> 
   31template<
class TypeTag, 
class DiscretizationMethod, ReferenceSystemFormulation referenceSystem>
 
   38template <
class TypeTag, ReferenceSystemFormulation referenceSystem>
 
   46    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   47    using SubControlVolume = 
typename GridGeometry::SubControlVolume;
 
   48    using SubControlVolumeFace = 
typename GridGeometry::SubControlVolumeFace;
 
   52    using ElementFluxVariablesCache = 
typename GridFluxVariablesCache::LocalView;
 
   53    using FluxVarCache = 
typename GridFluxVariablesCache::FluxVariablesCache;
 
   58    using Element = 
typename GridView::template Codim<0>::Entity;
 
   60    using Indices = 
typename ModelTraits::Indices;
 
   62    static constexpr int dim = GridView::dimension;
 
   63    static constexpr int dimWorld = GridView::dimensionworld;
 
   65    static constexpr int numPhases = ModelTraits::numFluidPhases();
 
   66    static constexpr int numComponents = ModelTraits::numFluidComponents();
 
   68    using DimWorldMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
 
   69    using ComponentFluxVector = Dune::FieldVector<Scalar, numComponents>;
 
   70    using HeatFluxScalar = Scalar;
 
   78    { 
return referenceSystem; }
 
 
   87                                                           const Element& element,
 
   88                                                           const FVElementGeometry& fvGeometry,
 
   89                                                           const ElementVolumeVariables& elemVolVars,
 
   90                                                           const SubControlVolumeFace& scvf,
 
   92                                                           const ElementFluxVariablesCache& elemFluxVarsCache)
 
   94        ComponentFluxVector componentFlux(0.0);
 
   96        const auto& fluxVarsCache = elemFluxVarsCache[scvf];
 
   97        const auto& shapeValues = fluxVarsCache.shapeValues();
 
  100        Scalar rhoMassOrMole(0.0);
 
  101        for (
auto&& scv : scvs(fvGeometry))
 
  104            rhoMassOrMole += rho * shapeValues[scv.indexInElement()][0];
 
  107        const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
 
  108        const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
 
  110        for (
int compIdx = 0; compIdx < numComponents; compIdx++)
 
  112            if constexpr (!FluidSystem::isTracerFluidSystem())
 
  113                if (compIdx == FluidSystem::getMainComponent(phaseIdx))
 
  117            const auto& dispersionTensor = [&]()
 
  120                    ModelTraits::CompositionalDispersionModel::compositionalDispersionTensor(problem, scvf, fvGeometry,
 
  121                                                                                             elemVolVars, elemFluxVarsCache,
 
  124                if (Dune::FloatCmp::eq(insideVolVars.extrusionFactor(), outsideVolVars.extrusionFactor(), 1e-6))
 
  125                    return insideVolVars.extrusionFactor()*tensor;
 
  128                    const auto insideTensor = insideVolVars.extrusionFactor() * tensor;
 
  129                    const auto outsideTensor = outsideVolVars.extrusionFactor() * tensor;
 
  137            Dune::FieldVector<Scalar, dimWorld> gradX(0.0);
 
  138            for (
auto&& scv : scvs(fvGeometry))
 
  140                const auto x = 
massOrMoleFraction(elemVolVars[scv], referenceSystem, phaseIdx, compIdx);
 
  141                gradX.axpy(x, fluxVarsCache.gradN(scv.indexInElement()));
 
  145            componentFlux[compIdx] = -1.0 * rhoMassOrMole * 
vtmv(scvf.unitOuterNormal(), dispersionTensor, gradX)*Extrusion::area(fvGeometry, scvf);
 
  146            if constexpr (!FluidSystem::isTracerFluidSystem())
 
  147                if (BalanceEqOpts::mainComponentIsBalanced(phaseIdx))
 
  148                    componentFlux[FluidSystem::getMainComponent(phaseIdx)] -= componentFlux[compIdx];
 
  150        return componentFlux;
 
 
  158                                                const Element& element,
 
  159                                                const FVElementGeometry& fvGeometry,
 
  160                                                const ElementVolumeVariables& elemVolVars,
 
  161                                                const SubControlVolumeFace& scvf,
 
  163                                                const ElementFluxVariablesCache& elemFluxVarsCache)
 
  165        const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
 
  166        const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
 
  169        const auto& dispersionTensor = [&]()
 
  172                ModelTraits::ThermalDispersionModel::thermalDispersionTensor(problem, scvf, fvGeometry,
 
  173                                                                             elemVolVars, elemFluxVarsCache,
 
  176            if (Dune::FloatCmp::eq(insideVolVars.extrusionFactor(), outsideVolVars.extrusionFactor(), 1e-6))
 
  177                return insideVolVars.extrusionFactor()*tensor;
 
  180                const auto insideTensor = insideVolVars.extrusionFactor() * tensor;
 
  181                const auto outsideTensor = outsideVolVars.extrusionFactor() * tensor;
 
  189        const auto& fluxVarsCache = elemFluxVarsCache[scvf];
 
  190        Dune::FieldVector<Scalar, GridView::dimensionworld> gradTemp(0.0);
 
  191        for (
auto&& scv : scvs(fvGeometry))
 
  192            gradTemp.axpy(elemVolVars[scv].temperature(), fluxVarsCache.gradN(scv.indexInElement()));
 
  195        return -1.0*
vtmv(scvf.unitOuterNormal(), dispersionTensor, gradTemp)*Extrusion::area(fvGeometry, scvf);
 
 
 
static constexpr ReferenceSystemFormulation referenceSystemFormulation()
Definition box/dispersionflux.hh:77
static ComponentFluxVector compositionalDispersionFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const int phaseIdx, const ElementFluxVariablesCache &elemFluxVarsCache)
Returns the dispersive fluxes of all components within a fluid phase across the given sub-control vol...
Definition box/dispersionflux.hh:86
static HeatFluxScalar thermalDispersionFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const int phaseIdx, const ElementFluxVariablesCache &elemFluxVarsCache)
Returns the thermal dispersive flux across the given sub-control volume face.
Definition box/dispersionflux.hh:157
Definition box/dispersionflux.hh:32
Defines all properties used in Dumux.
Helper classes to compute the integration elements.
A free function to average a Tensor at an interface.
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:880
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.
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:29
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition extrusion.hh:166
Traits of a flux variables type.
Definition flux/traits.hh:32
static constexpr bool hasStationaryVelocityField()
Definition flux/traits.hh:33