12#ifndef DUMUX_DISCRETIZATION_CC_TPFA_FOURIERS_LAW_NONEQUILIBRIUM_HH 
   13#define DUMUX_DISCRETIZATION_CC_TPFA_FOURIERS_LAW_NONEQUILIBRIUM_HH 
   24template<
class TypeTag, 
class DiscretizationMethod>
 
   31template <
class TypeTag>
 
   34    using Implementation = FouriersLawNonEquilibriumImplementation<TypeTag, DiscretizationMethods::CCTpfa>;
 
   38    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   39    using SubControlVolumeFace = 
typename GridGeometry::SubControlVolumeFace;
 
   41    using GridView = 
typename GridGeometry::GridView;
 
   42    using ElementVolumeVariables = 
typename GetPropType<TypeTag, Properties::GridVolumeVariables>::LocalView;
 
   43    using Element = 
typename GridView::template Codim<0>::Entity;
 
   44    using ElementFluxVarsCache = 
typename GetPropType<TypeTag, Properties::GridFluxVariablesCache>::LocalView;
 
   46    static constexpr int dim = GridView::dimension;
 
   47    static constexpr int dimWorld = GridView::dimensionworld;
 
   53    static constexpr auto numEnergyEq = numEnergyEqSolid + numEnergyEqFluid;
 
   54    static constexpr auto sPhaseIdx = ModelTraits::numFluidPhases();
 
   57    using DiscretizationMethod = DiscretizationMethods::CCTpfa;
 
   59    static constexpr DiscretizationMethod discMethod{};
 
   61    using Cache = FluxVariablesCaching::EmptyHeatConductionCache;
 
   67    static Scalar flux(
const Problem& problem,
 
   68                       const Element& element,
 
   69                       const FVElementGeometry& fvGeometry,
 
   70                       const ElementVolumeVariables& elemVolVars,
 
   71                       const SubControlVolumeFace& scvf,
 
   73                       const ElementFluxVarsCache& elemFluxVarsCache)
 
   76        Scalar tOutside = 0.0;
 
   78        if (phaseIdx < numEnergyEqFluid)
 
   80            tInside += elemVolVars[scvf.insideScvIdx()].temperatureFluid(phaseIdx);
 
   81            tOutside += elemVolVars[scvf.outsideScvIdx()].temperatureFluid(phaseIdx);
 
   85            tInside += elemVolVars[scvf.insideScvIdx()].temperatureSolid();
 
   86            tOutside += elemVolVars[scvf.outsideScvIdx()].temperatureSolid();
 
   89        Scalar tij = calculateTransmissibility(problem, element, fvGeometry, elemVolVars, scvf, phaseIdx);
 
   90        return tij*(tInside - tOutside);
 
   94    static Scalar calculateTransmissibility(
const Problem& problem,
 
   95                                            const Element& element,
 
   96                                            const FVElementGeometry& fvGeometry,
 
   97                                            const ElementVolumeVariables& elemVolVars,
 
   98                                            const SubControlVolumeFace& scvf,
 
  101        const auto insideScvIdx = scvf.insideScvIdx();
 
  102        const auto& insideScv = fvGeometry.scv(insideScvIdx);
 
  103        const auto& insideVolVars = elemVolVars[insideScvIdx];
 
  104        const auto computeLambda = [&](
const auto& v){
 
  105            if constexpr (numEnergyEq == 1)
 
  106                return v.effectiveThermalConductivity();
 
  107            else if constexpr (numEnergyEqFluid == 1)
 
  108                return (phaseIdx != sPhaseIdx)
 
  109                        ? v.effectiveFluidThermalConductivity()
 
  110                        : v.effectiveSolidThermalConductivity();
 
  112                return v.effectivePhaseThermalConductivity(phaseIdx);
 
  115        const auto insideLambda = computeLambda(insideVolVars);
 
  119        if (scvf.boundary() || scvf.numOutsideScvs() > 1)
 
  120            return Extrusion::area(fvGeometry, scvf)*ti;
 
  123            const auto outsideScvIdx = scvf.outsideScvIdx();
 
  124            const auto& outsideScv = fvGeometry.scv(outsideScvIdx);
 
  125            const auto& outsideVolVars = elemVolVars[outsideScvIdx];
 
  126            const auto outsideLambda = computeLambda(outsideVolVars);
 
  133                tj = 
computeTpfaTransmissibility(fvGeometry, fvGeometry.flipScvf(scvf.index()), outsideScv, outsideLambda, outsideVolVars.extrusionFactor());
 
  139                return Extrusion::area(fvGeometry, scvf)*(ti * tj)/(ti + tj);
 
Definition box/fourierslawnonequilibrium.hh:30
Defines all properties used in Dumux.
Helper classes to compute the integration elements.
Classes related to flux variables caching.
Tensor::field_type computeTpfaTransmissibility(const FVElementGeometry &fvGeometry, const typename FVElementGeometry::SubControlVolumeFace &scvf, const typename FVElementGeometry::SubControlVolume &scv, const Tensor &T, typename FVElementGeometry::SubControlVolume::Traits::Scalar extrusionFactor)
Free function to evaluate the Tpfa transmissibility associated with the flux (in the form of flux = T...
Definition tpfa/computetransmissibility.hh:36
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:310
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
The available discretization methods in Dumux.
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition extrusion.hh:166
Free functions to evaluate the transmissibilities associated with flux evaluations across sub-control...