3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
box/fourierslaw.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
25#ifndef DUMUX_DISCRETIZATION_BOX_FOURIERS_LAW_HH
26#define DUMUX_DISCRETIZATION_BOX_FOURIERS_LAW_HH
27
28#include <dumux/common/math.hh>
32
33namespace Dumux {
34
35// forward declaration
36template<class TypeTag, DiscretizationMethod discMethod>
37class FouriersLawImplementation;
38
43template <class TypeTag>
45{
48 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
49 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
51 using ElementVolumeVariables = typename GetPropType<TypeTag, Properties::GridVolumeVariables>::LocalView;
52 using ElementFluxVariablesCache = typename GetPropType<TypeTag, Properties::GridFluxVariablesCache>::LocalView;
54 using Element = typename GridView::template Codim<0>::Entity;
55
56public:
57 static Scalar flux(const Problem& problem,
58 const Element& element,
59 const FVElementGeometry& fvGeometry,
60 const ElementVolumeVariables& elemVolVars,
61 const SubControlVolumeFace& scvf,
62 const ElementFluxVariablesCache& elemFluxVarsCache)
63 {
64 // get inside and outside diffusion tensors and calculate the harmonic mean
65 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
66 const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
67 const auto& insideVolVars = elemVolVars[insideScv];
68 const auto& outsideVolVars = elemVolVars[outsideScv];
69
70 // effective diffusion tensors
71 auto insideLambda = Deprecated::template effectiveThermalConductivity<ThermalConductivityModel>(
72 insideVolVars, problem.spatialParams(), element, fvGeometry, insideScv);
73 auto outsideLambda = Deprecated::template effectiveThermalConductivity<ThermalConductivityModel>(
74 outsideVolVars, problem.spatialParams(), element, fvGeometry, outsideScv);
75
76 // scale by extrusion factor
77 insideLambda *= insideVolVars.extrusionFactor();
78 outsideLambda *= outsideVolVars.extrusionFactor();
79
80 // the resulting averaged diffusion tensor
81 const auto lambda = problem.spatialParams().harmonicMean(insideLambda, outsideLambda, scvf.unitOuterNormal());
82
83 // evaluate gradTemp at integration point
84 const auto& fluxVarsCache = elemFluxVarsCache[scvf];
85
86 // compute the temperature gradient with the shape functions
87 Dune::FieldVector<Scalar, GridView::dimensionworld> gradTemp(0.0);
88 for (auto&& scv : scvs(fvGeometry))
89 gradTemp.axpy(elemVolVars[scv].temperature(), fluxVarsCache.gradN(scv.indexInElement()));
90
91 // comute the heat conduction flux
92 return -1.0*vtmv(scvf.unitOuterNormal(), lambda, gradTemp)*scvf.area();
93 }
94};
95
96} // end namespace Dumux
97
98#endif
Helpers for deprecation.
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
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
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: 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: box/fourierslaw.hh:57
Declares all properties used in Dumux.