24#ifndef DUMUX_DISCRETIZATION_CC_TPFA_FOURIERS_LAW_HH
25#define DUMUX_DISCRETIZATION_CC_TPFA_FOURIERS_LAW_HH
37template<
class TypeTag, DiscretizationMethod discMethod>
38class FouriersLawImplementation;
44template <
class TypeTag>
51 using FVElementGeometry =
typename GridGeometry::LocalView;
52 using SubControlVolumeFace =
typename GridGeometry::SubControlVolumeFace;
54 using GridView =
typename GridGeometry::GridView;
56 using VolumeVariables =
typename ElementVolumeVariables::VolumeVariables;
57 using Element =
typename GridView::template Codim<0>::Entity;
61 static const int dim = GridView::dimension;
62 static const int dimWorld = GridView::dimensionworld;
65 class TpfaFouriersLawCacheFiller
70 template<
class FluxVariablesCacheFiller>
71 static void fill(FluxVariablesCache& scvfFluxVarsCache,
72 const Problem& problem,
73 const Element& element,
74 const FVElementGeometry& fvGeometry,
75 const ElementVolumeVariables& elemVolVars,
76 const SubControlVolumeFace& scvf,
77 const FluxVariablesCacheFiller& fluxVarsCacheFiller)
79 scvfFluxVarsCache.updateHeatConduction(problem, element, fvGeometry, elemVolVars, scvf);
84 class TpfaFouriersLawCache
87 using Filler = TpfaFouriersLawCacheFiller;
89 void updateHeatConduction(
const Problem& problem,
90 const Element& element,
91 const FVElementGeometry& fvGeometry,
92 const ElementVolumeVariables& elemVolVars,
93 const SubControlVolumeFace &scvf)
95 tij_ = Implementation::calculateTransmissibility(problem, element, fvGeometry, elemVolVars, scvf);
98 const Scalar& heatConductionTij()
const
123 static Scalar
flux(
const Problem& problem,
124 const Element& element,
125 const FVElementGeometry& fvGeometry,
126 const VolumeVariables& insideVolVars,
127 const VolumeVariables& outsideVolVars,
128 const SubControlVolumeFace& scvf,
129 const ElementFluxVarsCache& elemFluxVarsCache)
131 if constexpr (isMixedDimensional_)
132 if (scvf.numOutsideScv() != 1)
133 DUNE_THROW(Dune::Exception,
"This flux overload requires scvf.numOutsideScv() == 1");
136 Scalar tij = elemFluxVarsCache[scvf].heatConductionTij();
139 const auto tInside = insideVolVars.temperature();
140 const auto tOutside = outsideVolVars.temperature();
142 return tij*(tInside - tOutside);
152 static Scalar
flux(
const Problem& problem,
153 const Element& element,
154 const FVElementGeometry& fvGeometry,
155 const ElementVolumeVariables& elemVolVars,
156 const SubControlVolumeFace& scvf,
157 const ElementFluxVarsCache& elemFluxVarsCache)
160 Scalar tij = elemFluxVarsCache[scvf].heatConductionTij();
163 const auto tInside = elemVolVars[scvf.insideScvIdx()].temperature();
164 const auto tOutside = scvf.numOutsideScvs() == 1 ? elemVolVars[scvf.outsideScvIdx()].temperature()
165 : branchingFacetTemperature_(problem, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf, tInside, tij);
167 return tij*(tInside - tOutside);
172 const Element& element,
173 const FVElementGeometry& fvGeometry,
174 const ElementVolumeVariables& elemVolVars,
175 const SubControlVolumeFace& scvf)
179 const auto insideScvIdx = scvf.insideScvIdx();
180 const auto& insideScv = fvGeometry.scv(insideScvIdx);
181 const auto& insideVolVars = elemVolVars[insideScvIdx];
183 const auto insideLambda = insideVolVars.effectiveThermalConductivity();
187 if (scvf.boundary() || scvf.numOutsideScvs() > 1)
189 tij = Extrusion::area(scvf)*ti;
194 const auto outsideScvIdx = scvf.outsideScvIdx();
195 const auto& outsideScv = fvGeometry.scv(outsideScvIdx);
196 const auto& outsideVolVars = elemVolVars[outsideScvIdx];
198 const auto outsideLambda = outsideVolVars.effectiveThermalConductivity();
200 if constexpr (dim == dimWorld)
210 tij = Extrusion::area(scvf)*(ti * tj)/(ti + tj);
219 static Scalar branchingFacetTemperature_(
const Problem& problem,
220 const Element& element,
221 const FVElementGeometry& fvGeometry,
222 const ElementVolumeVariables& elemVolVars,
223 const ElementFluxVarsCache& elemFluxVarsCache,
224 const SubControlVolumeFace& scvf,
225 Scalar insideTemperature,
228 Scalar sumTi(insideTi);
229 Scalar sumTempTi(insideTi*insideTemperature);
231 for (
unsigned int i = 0; i < scvf.numOutsideScvs(); ++i)
233 const auto outsideScvIdx = scvf.outsideScvIdx(i);
234 const auto& outsideVolVars = elemVolVars[outsideScvIdx];
235 const auto& flippedScvf = fvGeometry.flipScvf(scvf.index(), i);
236 const auto& outsideFluxVarsCache = elemFluxVarsCache[flippedScvf];
238 auto outsideTi = outsideFluxVarsCache.heatConductionTij();
240 sumTempTi += outsideTi*outsideVolVars.temperature();
242 return sumTempTi/sumTi;
245 static constexpr bool isMixedDimensional_ =
static_cast<int>(GridView::dimension) <
static_cast<int>(GridView::dimensionworld);
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
The available discretization methods in Dumux.
Helper classes to compute the integration elements.
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 Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition: extrusion.hh:177
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:150
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:46
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:171
static Scalar flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVarsCache &elemFluxVarsCache)
Returns the heat flux within the porous medium (in J/s) across the given sub-control volume face.
Definition: flux/cctpfa/fourierslaw.hh:152
TpfaFouriersLawCache Cache
export the type for the corresponding cache
Definition: flux/cctpfa/fourierslaw.hh:110
static Scalar flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const VolumeVariables &insideVolVars, const VolumeVariables &outsideVolVars, const SubControlVolumeFace &scvf, const ElementFluxVarsCache &elemFluxVarsCache)
Returns the heat flux within the porous medium (in J/s) across the given sub-control volume face.
Definition: flux/cctpfa/fourierslaw.hh:123
Declares all properties used in Dumux.
Free functions to evaluate the transmissibilities associated with flux evaluations across sub-control...