25#ifndef DUMUX_DISCRETIZATION_BOX_FOURIERS_LAW_HH
26#define DUMUX_DISCRETIZATION_BOX_FOURIERS_LAW_HH
36template<
class TypeTag, DiscretizationMethod discMethod>
37class FouriersLawImplementation;
43template <
class TypeTag>
49 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
53 using Element =
typename GridView::template Codim<0>::Entity;
56 static Scalar
flux(
const Problem& problem,
57 const Element& element,
58 const FVElementGeometry& fvGeometry,
59 const ElementVolumeVariables& elemVolVars,
60 const SubControlVolumeFace& scvf,
61 const ElementFluxVariablesCache& elemFluxVarsCache)
64 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
65 const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
66 const auto& insideVolVars = elemVolVars[insideScv];
67 const auto& outsideVolVars = elemVolVars[outsideScv];
70 auto insideLambda = insideVolVars.effectiveThermalConductivity();
71 auto outsideLambda = outsideVolVars.effectiveThermalConductivity();
74 insideLambda *= insideVolVars.extrusionFactor();
75 outsideLambda *= outsideVolVars.extrusionFactor();
78 const auto lambda = problem.spatialParams().harmonicMean(insideLambda, outsideLambda, scvf.unitOuterNormal());
81 const auto& fluxVarsCache = elemFluxVarsCache[scvf];
84 Dune::FieldVector<Scalar, GridView::dimensionworld> gradTemp(0.0);
85 for (
auto&& scv : scvs(fvGeometry))
86 gradTemp.axpy(elemVolVars[scv].
temperature(), fluxVarsCache.gradN(scv.indexInElement()));
89 return -1.0*
vtmv(scvf.unitOuterNormal(), lambda, gradTemp)*scvf.area();
Define some often used mathematical functions.
The available discretization methods in Dumux.
DiscretizationMethod
The available discretization methods in Dumux.
Definition: method.hh:37
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:840
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
std::string temperature() noexcept
I/O name of temperature for equilibrium models.
Definition: name.hh:51
forward declaration of the method-specific implementation
Definition: flux/fourierslaw.hh:37
static Scalar flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache)
Definition: flux/box/fourierslaw.hh:56
Declares all properties used in Dumux.