24#ifndef DUMUX_DISCRETIZATION_CC_TPFA_FOURIERS_LAW_HH
25#define DUMUX_DISCRETIZATION_CC_TPFA_FOURIERS_LAW_HH
38template<
class TypeTag, DiscretizationMethod discMethod>
39class FouriersLawImplementation;
45template <
class TypeTag>
52 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
55 using Element =
typename GridView::template Codim<0>::Entity;
59 static const int dim = GridView::dimension;
60 static const int dimWorld = GridView::dimensionworld;
63 class TpfaFouriersLawCacheFiller
68 template<
class FluxVariablesCacheFiller>
69 static void fill(FluxVariablesCache& scvfFluxVarsCache,
70 const Problem& problem,
71 const Element& element,
72 const FVElementGeometry& fvGeometry,
73 const ElementVolumeVariables& elemVolVars,
74 const SubControlVolumeFace& scvf,
75 const FluxVariablesCacheFiller& fluxVarsCacheFiller)
77 scvfFluxVarsCache.updateHeatConduction(problem, element, fvGeometry, elemVolVars, scvf);
82 class TpfaFouriersLawCache
85 using Filler = TpfaFouriersLawCacheFiller;
87 void updateHeatConduction(
const Problem& problem,
88 const Element& element,
89 const FVElementGeometry& fvGeometry,
90 const ElementVolumeVariables& elemVolVars,
91 const SubControlVolumeFace &scvf)
93 tij_ = Implementation::calculateTransmissibility(problem, element, fvGeometry, elemVolVars, scvf);
96 const Scalar& heatConductionTij()
const
111 static Scalar
flux(
const Problem& problem,
112 const Element& element,
113 const FVElementGeometry& fvGeometry,
114 const ElementVolumeVariables& elemVolVars,
115 const SubControlVolumeFace& scvf,
116 const ElementFluxVarsCache& elemFluxVarsCache)
119 Scalar tij = elemFluxVarsCache[scvf].heatConductionTij();
122 const auto tInside = elemVolVars[scvf.insideScvIdx()].temperature();
123 const auto tOutside = scvf.numOutsideScvs() == 1 ? elemVolVars[scvf.outsideScvIdx()].temperature()
124 : branchingFacetTemperature_(problem, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf, tInside, tij);
126 return tij*(tInside - tOutside);
131 const Element& element,
132 const FVElementGeometry& fvGeometry,
133 const ElementVolumeVariables& elemVolVars,
134 const SubControlVolumeFace& scvf)
138 const auto insideScvIdx = scvf.insideScvIdx();
139 const auto& insideScv = fvGeometry.scv(insideScvIdx);
140 const auto& insideVolVars = elemVolVars[insideScvIdx];
142 const auto insideLambda = insideVolVars.effectiveThermalConductivity();
146 if (scvf.boundary() || scvf.numOutsideScvs() > 1)
148 tij = scvf.area()*ti;
153 const auto outsideScvIdx = scvf.outsideScvIdx();
154 const auto& outsideScv = fvGeometry.scv(outsideScvIdx);
155 const auto& outsideVolVars = elemVolVars[outsideScvIdx];
157 const auto outsideLambda = outsideVolVars.effectiveThermalConductivity();
169 tij = scvf.area()*(ti * tj)/(ti + tj);
178 static Scalar branchingFacetTemperature_(
const Problem& problem,
179 const Element& element,
180 const FVElementGeometry& fvGeometry,
181 const ElementVolumeVariables& elemVolVars,
182 const ElementFluxVarsCache& elemFluxVarsCache,
183 const SubControlVolumeFace& scvf,
184 Scalar insideTemperature,
187 Scalar sumTi(insideTi);
188 Scalar sumTempTi(insideTi*insideTemperature);
190 for (
unsigned int i = 0; i < scvf.numOutsideScvs(); ++i)
192 const auto outsideScvIdx = scvf.outsideScvIdx(i);
193 const auto& outsideVolVars = elemVolVars[outsideScvIdx];
194 const auto& flippedScvf = fvGeometry.flipScvf(scvf.index(), i);
195 const auto& outsideFluxVarsCache = elemFluxVarsCache[flippedScvf];
197 auto outsideTi = outsideFluxVarsCache.heatConductionTij();
199 sumTempTi += outsideTi*outsideVolVars.temperature();
201 return sumTempTi/sumTi;
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
The available discretization methods in Dumux.
DiscretizationMethod
The available discretization methods in Dumux.
Definition: method.hh:37
Tensor::field_type computeTpfaTransmissibility(const SubControlVolumeFace &scvf, const SubControlVolume &scv, const Tensor &T, typename 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:47
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:149
forward declaration of the method-specific implementation
Definition: flux/fourierslaw.hh:37
Fourier's law for cell-centered finite volume schemes with two-point flux approximation.
Definition: flux/cctpfa/fourierslaw.hh:47
static Scalar calculateTransmissibility(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf)
Compute transmissibilities.
Definition: flux/cctpfa/fourierslaw.hh:130
static Scalar flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVarsCache &elemFluxVarsCache)
Compute the heat condution flux assuming thermal equilibrium.
Definition: flux/cctpfa/fourierslaw.hh:111
TpfaFouriersLawCache Cache
export the type for the corresponding cache
Definition: flux/cctpfa/fourierslaw.hh:108
Declares all properties used in Dumux.
Free functions to evaluate the transmissibilities associated with flux evaluations across sub-control...