16#ifndef DUMUX_DISCRETIZATION_CVFE_DARCYS_LAW_HH 
   17#define DUMUX_DISCRETIZATION_CVFE_DARCYS_LAW_HH 
   34template<
class Scalar, 
class Gr
idGeometry>
 
   37    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   38    using SubControlVolumeFace = 
typename FVElementGeometry::SubControlVolumeFace;
 
   40    using GridView = 
typename GridGeometry::GridView;
 
   41    using Element = 
typename GridView::template Codim<0>::Entity;
 
   43    static constexpr int dim = GridView::dimension;
 
   44    static constexpr int dimWorld = GridView::dimensionworld;
 
   58    template<
class Problem, 
class ElementVolumeVariables, 
class ElementFluxVarsCache>
 
   59    static Scalar 
flux(
const Problem& problem,
 
   60                       const Element& element,
 
   61                       const FVElementGeometry& fvGeometry,
 
   62                       const ElementVolumeVariables& elemVolVars,
 
   63                       const SubControlVolumeFace& scvf,
 
   65                       const ElementFluxVarsCache& elemFluxVarCache)
 
   67        const auto& fluxVarCache = elemFluxVarCache[scvf];
 
   68        const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
 
   69        const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
 
   70        const auto& insideVolVars = elemVolVars[insideScv];
 
   71        const auto& outsideVolVars = elemVolVars[outsideScv];
 
   73        auto insideK = insideVolVars.permeability();
 
   74        auto outsideK = outsideVolVars.permeability();
 
   77        insideK *= insideVolVars.extrusionFactor();
 
   78        outsideK *= outsideVolVars.extrusionFactor();
 
   83        const auto& shapeValues = fluxVarCache.shapeValues();
 
   86        Dune::FieldVector<Scalar, dimWorld> gradP(0.0);
 
   88        for (
auto&& scv : scvs(fvGeometry))
 
   90            const auto& volVars = elemVolVars[scv];
 
   93                rho += volVars.density(phaseIdx)*shapeValues[scv.indexInElement()][0];
 
   96            gradP.axpy(volVars.pressure(phaseIdx), fluxVarCache.gradN(scv.indexInElement()));
 
  100            gradP.axpy(-rho, problem.spatialParams().gravity(scvf.center()));
 
  103        return -1.0*
vtmv(scvf.unitOuterNormal(), K, gradP)*Extrusion::area(fvGeometry, scvf);
 
 
  107    template<
class Problem, 
class ElementVolumeVariables, 
class FluxVarCache>
 
  109                                                           const Element& element,
 
  110                                                           const FVElementGeometry& fvGeometry,
 
  111                                                           const ElementVolumeVariables& elemVolVars,
 
  112                                                           const SubControlVolumeFace& scvf,
 
  113                                                           const FluxVarCache& fluxVarCache)
 
  115        const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
 
  116        const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
 
  117        const auto& insideVolVars = elemVolVars[insideScv];
 
  118        const auto& outsideVolVars = elemVolVars[outsideScv];
 
  120        auto insideK = insideVolVars.permeability();
 
  121        auto outsideK = outsideVolVars.permeability();
 
  124        insideK *= insideVolVars.extrusionFactor();
 
  125        outsideK *= outsideVolVars.extrusionFactor();
 
  129        std::vector<Scalar> ti(fvGeometry.numScv());
 
  130        for (
const auto& scv : scvs(fvGeometry))
 
  131            ti[scv.indexInElement()] =
 
  132                -1.0*Extrusion::area(fvGeometry, scvf)*
vtmv(scvf.unitOuterNormal(), K, fluxVarCache.gradN(scv.indexInElement()));
 
 
 
  139template<
class TypeTag, 
class DiscretizationMethod>
 
  140class DarcysLawImplementation;
 
  146template<
class TypeTag, 
class DM>
 
  148: 
public CVFEDarcysLaw<GetPropType<TypeTag, Properties::Scalar>, GetPropType<TypeTag, Properties::GridGeometry>>
 
 
Darcy's law for control-volume finite element schemes.
Definition flux/cvfe/darcyslaw.hh:36
static std::vector< Scalar > calculateTransmissibilities(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const FluxVarCache &fluxVarCache)
Definition flux/cvfe/darcyslaw.hh:108
static Scalar flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const int phaseIdx, const ElementFluxVarsCache &elemFluxVarCache)
Returns the advective flux of a fluid phase across the given sub-control volume face.
Definition flux/cvfe/darcyslaw.hh:59
forward declaration of the method-specific implementation
Definition flux/ccmpfa/darcyslaw.hh:27
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
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition parameters.hh:149
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
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.