12#ifndef DUMUX_PNM_2P_FLUXVARIABLESCACHE_HH 
   13#define DUMUX_PNM_2P_FLUXVARIABLESCACHE_HH 
   17#include <dune/common/reservedvector.hh> 
   27template<
class AdvectionType, 
int maxNumCorners = 4>
 
   30    using Scalar = 
typename AdvectionType::Scalar;
 
   31    static constexpr auto numPhases = 2;
 
   32    using NumCornerVector = Dune::ReservedVector<Scalar, maxNumCorners>;
 
   38    template<
class Problem, 
class Element, 
class FVElementGeometry,
 
   39             class ElementVolumeVariables, 
class SubControlVolumeFace>
 
   41                const Element& element,
 
   42                const FVElementGeometry& fvGeometry,
 
   43                const ElementVolumeVariables& elemVolVars,
 
   44                const SubControlVolumeFace& scvf,
 
   47        const auto eIdx = fvGeometry.gridGeometry().elementMapper().index(element);
 
   48        throatCrossSectionShape_ = fvGeometry.gridGeometry().throatCrossSectionShape(eIdx);
 
   49        throatShapeFactor_ = fvGeometry.gridGeometry().throatShapeFactor(eIdx);
 
   50        pc_ = std::max(elemVolVars[0].capillaryPressure(), elemVolVars[1].capillaryPressure());
 
   51        pcEntry_ = problem.spatialParams().pcEntry(element, elemVolVars);
 
   52        pcSnapoff_ = problem.spatialParams().pcSnapoff(element, elemVolVars);
 
   53        throatInscribedRadius_ = problem.spatialParams().throatInscribedRadius(element, elemVolVars);
 
   54        throatLength_ = problem.spatialParams().throatLength(element, elemVolVars);
 
   56        poreToPoreDistance_ = element.geometry().volume();
 
   59        using FluidSystem = 
typename ElementVolumeVariables::VolumeVariables::FluidSystem;
 
   60        const auto& spatialParams = problem.spatialParams();
 
   61        nPhaseIdx_ = 1 - spatialParams.template wettingPhase<FluidSystem>(element, elemVolVars);
 
   64        surfaceTension_ = 0.5*(elemVolVars[0].surfaceTension() + elemVolVars[1].
surfaceTension());
 
   66        const auto& cornerHalfAngles = spatialParams.cornerHalfAngles(element);
 
   67        wettingLayerArea_.clear(); wettingLayerArea_.resize(cornerHalfAngles.size());
 
   68        const Scalar totalThroatCrossSectionalArea = spatialParams.throatCrossSectionalArea(element, elemVolVars);
 
   72            const Scalar theta = spatialParams.contactAngle(element, elemVolVars);
 
   73            for (
int i = 0; i< cornerHalfAngles.size(); ++i)
 
   77            throatCrossSectionalArea_[
wPhaseIdx()] = std::min(
 
   78                std::accumulate(wettingLayerArea_.begin(), wettingLayerArea_.end(), 0.0),
 
   79                totalThroatCrossSectionalArea
 
   81            throatCrossSectionalArea_[
nPhaseIdx()] = totalThroatCrossSectionalArea - throatCrossSectionalArea_[
wPhaseIdx()];
 
   85            for (
int i = 0; i< cornerHalfAngles.size(); ++i)
 
   86                wettingLayerArea_[i] = 0.0;
 
   88            throatCrossSectionalArea_[
wPhaseIdx()] = totalThroatCrossSectionalArea;
 
   89            throatCrossSectionalArea_[
nPhaseIdx()] = 0.0;
 
   92        for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
 
   94            singlePhaseCache_.fill(problem, element, fvGeometry, scvf, elemVolVars, *
this, phaseIdx);
 
   95            nonWettingPhaseCache_.fill(problem, element, fvGeometry, scvf, elemVolVars, *
this, phaseIdx);
 
   96            wettingLayerCache_.fill(problem, element, fvGeometry, scvf, elemVolVars, *
this, phaseIdx);
 
   99        for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
 
  101            transmissibility_[phaseIdx] = AdvectionType::calculateTransmissibility(
 
  102                problem, element, fvGeometry, scvf, elemVolVars, *
this, phaseIdx
 
 
  111    { 
return throatCrossSectionShape_; }
 
 
  117    { 
return throatShapeFactor_; }
 
 
  123    { 
return transmissibility_[phaseIdx]; }
 
 
  129    { 
return throatCrossSectionalArea_[phaseIdx]; }
 
 
  135    { 
return throatCrossSectionalArea_[0] + throatCrossSectionalArea_[1]; }
 
 
  141    { 
return throatLength_; }
 
 
  147    { 
return throatInscribedRadius_; }
 
 
  159    { 
return pcSnapoff_; }
 
 
  171    { 
return surfaceTension_; }
 
 
  183    { 
return surfaceTension_ / pc_;}
 
 
  190    { 
return wettingLayerArea_[cornerIdx]; }
 
 
  196    { 
return 1 - nPhaseIdx_; }
 
 
  202    { 
return nPhaseIdx_; }
 
 
  208    { 
return singlePhaseCache_; }
 
 
  214    { 
return nonWettingPhaseCache_; }
 
 
  220    { 
return wettingLayerCache_; }
 
 
  226    { 
return poreToPoreDistance_; }
 
 
  230    Scalar throatShapeFactor_;
 
  231    std::array<Scalar, numPhases> transmissibility_;
 
  232    std::array<Scalar, numPhases> throatCrossSectionalArea_;
 
  233    Scalar throatLength_;
 
  234    Scalar throatInscribedRadius_;
 
  238    Scalar surfaceTension_;
 
  240    NumCornerVector wettingLayerArea_;
 
  241    std::size_t nPhaseIdx_;
 
  242    Scalar poreToPoreDistance_;
 
  244    typename AdvectionType::Transmissibility::SinglePhaseCache singlePhaseCache_;
 
  245    typename AdvectionType::Transmissibility::NonWettingPhaseCache nonWettingPhaseCache_;
 
  246    typename AdvectionType::Transmissibility::WettingLayerCache wettingLayerCache_;
 
 
Flux variables cache for the two-phase-flow PNM Store data required for flux calculation.
Definition porenetwork/2p/fluxvariablescache.hh:29
Scalar throatCrossSectionalArea() const
Returns the throats's total cross-sectional area.
Definition porenetwork/2p/fluxvariablescache.hh:134
Scalar pc() const
Returns the capillary pressure within the throat.
Definition porenetwork/2p/fluxvariablescache.hh:164
bool invaded() const
Definition porenetwork/2p/fluxvariablescache.hh:176
Scalar poreToPoreDistance() const
Returns the throats's pore-to-pore-center distance.
Definition porenetwork/2p/fluxvariablescache.hh:225
Scalar wettingLayerCrossSectionalArea(const int cornerIdx) const
Returns the cross-sectional area of a wetting layer within one of the throat's corners.
Definition porenetwork/2p/fluxvariablescache.hh:189
const auto & nonWettingPhaseFlowVariables() const
Returns the throats's cached flow variables for the nonwetting phase.
Definition porenetwork/2p/fluxvariablescache.hh:213
Scalar curvatureRadius() const
Returns the curvature radius within the throat.
Definition porenetwork/2p/fluxvariablescache.hh:182
Scalar throatCrossSectionalArea(const int phaseIdx) const
Returns the throats's cross-sectional area for a given phaseIdx.
Definition porenetwork/2p/fluxvariablescache.hh:128
Scalar transmissibility(const int phaseIdx) const
Returns the throats's transmissibility.
Definition porenetwork/2p/fluxvariablescache.hh:122
Scalar pcSnapoff() const
Returns the throats's snap-off capillary pressure.
Definition porenetwork/2p/fluxvariablescache.hh:158
static bool constexpr isSolDependent
Definition porenetwork/2p/fluxvariablescache.hh:36
Scalar throatLength() const
Returns the throats's length.
Definition porenetwork/2p/fluxvariablescache.hh:140
Throat::Shape throatCrossSectionShape() const
Returns the throats's cross-sectional shape.
Definition porenetwork/2p/fluxvariablescache.hh:110
Scalar throatShapeFactor() const
Returns the throats's shape factor.
Definition porenetwork/2p/fluxvariablescache.hh:116
const auto & singlePhaseFlowVariables() const
Returns the throats's cached flow variables for single-phase flow.
Definition porenetwork/2p/fluxvariablescache.hh:207
Scalar pcEntry() const
Returns the throats's entry capillary pressure.
Definition porenetwork/2p/fluxvariablescache.hh:152
Scalar throatInscribedRadius() const
Returns the throats's inscribed radius.
Definition porenetwork/2p/fluxvariablescache.hh:146
std::size_t nPhaseIdx() const
Returns the index of the nonwetting phase.
Definition porenetwork/2p/fluxvariablescache.hh:201
Scalar surfaceTension() const
Returns the surface tension within the throat.
Definition porenetwork/2p/fluxvariablescache.hh:170
std::size_t wPhaseIdx() const
Returns the index of the wetting phase.
Definition porenetwork/2p/fluxvariablescache.hh:195
void update(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, bool invaded)
Definition porenetwork/2p/fluxvariablescache.hh:40
const auto & wettingLayerFlowVariables() const
Returns the throats's cached flow variables for the wetting phase.
Definition porenetwork/2p/fluxvariablescache.hh:219
constexpr Scalar wettingLayerCrossSectionalArea(const Scalar curvatureRadius, const Scalar contactAngle, const Scalar cornerHalfAngle) noexcept
Return the cross-sectional area of a wetting layer residing in a corner of a throat.
Definition throatproperties.hh:243
Shape
Collection of different pore-throat shapes.
Definition throatproperties.hh:26
Definition discretization/porenetwork/fvelementgeometry.hh:24
This file contains functions related to calculate pore-throat properties.