3.4
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
multidomain/facet/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 *****************************************************************************/
24#ifndef DUMUX_DISCRETIZATION_BOX_FACET_COUPLING_FOURIERS_LAW_HH
25#define DUMUX_DISCRETIZATION_BOX_FACET_COUPLING_FOURIERS_LAW_HH
26
27#include <cmath>
28#include <vector>
29
30#include <dune/common/exceptions.hh>
31#include <dune/common/fvector.hh>
32#include <dune/common/float_cmp.hh>
33
34#include <dumux/common/math.hh>
37
42
43namespace Dumux {
44
51template<class TypeTag>
53: public FouriersLawImplementation<TypeTag, DiscretizationMethod::box>
54{
56
59 using FVElementGeometry = typename GridGeometry::LocalView;
60 using SubControlVolume = typename GridGeometry::SubControlVolume;
61 using SubControlVolumeFace = typename GridGeometry::SubControlVolumeFace;
62 using Extrusion = Extrusion_t<GridGeometry>;
63 using GridView = typename GridGeometry::GridView;
64 using Element = typename GridView::template Codim<0>::Entity;
65 using CoordScalar = typename GridView::ctype;
66 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
67
68 static constexpr int dim = GridView::dimension;
69 static constexpr int dimWorld = GridView::dimensionworld;
70
71public:
72
76 template<class Problem, class ElementVolumeVariables, class ElementFluxVarsCache>
77 static Scalar flux(const Problem& problem,
78 const Element& element,
79 const FVElementGeometry& fvGeometry,
80 const ElementVolumeVariables& elemVolVars,
81 const SubControlVolumeFace& scvf,
82 const ElementFluxVarsCache& elemFluxVarCache)
83 {
84 // if this scvf is not on an interior boundary, use the standard law
85 if (!scvf.interiorBoundary())
86 return ParentType::flux(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarCache);
87
88 static const Scalar xi = getParamFromGroup<Scalar>(problem.paramGroup(), "FacetCoupling.Xi", 1.0);
89 if ( !Dune::FloatCmp::eq(xi, 1.0, 1e-6) )
90 DUNE_THROW(Dune::NotImplemented, "Xi != 1.0 cannot be used with the Box-Facet-Coupling scheme");
91
92 // get some references for convenience
93 const auto& fluxVarCache = elemFluxVarCache[scvf];
94 const auto& shapeValues = fluxVarCache.shapeValues();
95 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
96 const auto& insideVolVars = elemVolVars[insideScv];
97
98 // on interior Neumann boundaries, evaluate the flux using the facet thermal conductivity
99 const auto bcTypes = problem.interiorBoundaryTypes(element, scvf);
100 if (bcTypes.hasOnlyNeumann())
101 {
102 // compute tpfa flux from integration point to facet centerline
103 const auto& facetVolVars = problem.couplingManager().getLowDimVolVars(element, scvf);
104
105 // interpolate temperature to scvf integration point
106 Scalar T = 0.0;
107 for (const auto& scv : scvs(fvGeometry))
108 T += elemVolVars[scv].temperature()*shapeValues[scv.indexInElement()][0];
109
110 using std::sqrt;
111 // If this is a surface grid, use the square root of the facet extrusion factor
112 // as an approximate average distance from scvf ip to facet center
113 using std::sqrt;
114 const auto a = facetVolVars.extrusionFactor();
115 auto gradT = scvf.unitOuterNormal();
116 gradT *= dim == dimWorld ? 0.5*a : 0.5*sqrt(a);
117 gradT /= gradT.two_norm2();
118 gradT *= (facetVolVars.temperature() - T);
119
120 return -1.0*Extrusion::area(scvf)
121 *insideVolVars.extrusionFactor()
122 *vtmv(scvf.unitOuterNormal(),
123 facetVolVars.effectiveThermalConductivity(),
124 gradT);
125 }
126
127 // on interior Dirichlet boundaries use the facet temperature and evaluate flux
128 else if (bcTypes.hasOnlyDirichlet())
129 {
130 // create vector with nodal temperatures
131 std::vector<Scalar> temperatures(element.subEntities(dim));
132 for (const auto& scv : scvs(fvGeometry))
133 temperatures[scv.localDofIndex()] = elemVolVars[scv].temperature();
134
135 // substitute with facet temperatures for those scvs touching this facet
136 for (const auto& scvfJ : scvfs(fvGeometry))
137 if (scvfJ.interiorBoundary() && scvfJ.facetIndexInElement() == scvf.facetIndexInElement())
138 temperatures[ fvGeometry.scv(scvfJ.insideScvIdx()).localDofIndex() ]
139 = problem.couplingManager().getLowDimVolVars(element, scvfJ).temperature();
140
141 // evaluate gradT at integration point
142 Dune::FieldVector<Scalar, dimWorld> gradT(0.0);
143 for (const auto& scv : scvs(fvGeometry))
144 gradT.axpy(temperatures[scv.localDofIndex()], fluxVarCache.gradN(scv.indexInElement()));
145
146 // apply matrix thermal conductivity and return the flux
147 return -1.0*Extrusion::area(scvf)
148 *insideVolVars.extrusionFactor()
149 *vtmv(scvf.unitOuterNormal(),
150 insideVolVars.effectiveThermalConductivity(),
151 gradT);
152 }
153
154 // mixed boundary types are not supported
155 else
156 DUNE_THROW(Dune::NotImplemented, "Mixed boundary types are not supported");
157 }
158
159 // compute transmissibilities ti for analytical jacobians
160 template<class Problem, class ElementVolumeVariables, class FluxVarCache>
161 static std::vector<Scalar> calculateTransmissibilities(const Problem& problem,
162 const Element& element,
163 const FVElementGeometry& fvGeometry,
164 const ElementVolumeVariables& elemVolVars,
165 const SubControlVolumeFace& scvf,
166 const FluxVarCache& fluxVarCache,
167 unsigned int phaseIdx)
168 {
169 DUNE_THROW(Dune::NotImplemented, "transmissibilty computation for BoxFacetCouplingFouriersLaw");
170 }
171};
172
173} // end namespace Dumux
174
175#endif
Define some often used mathematical functions.
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
The available discretization methods in Dumux.
Helper classes to compute the integration elements.
The reference frameworks and formulations available for splitting total fluxes into a advective and d...
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:863
Definition: adapt.hh:29
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
Specialization of Fourier's Law for the box method.
Definition: flux/box/fourierslaw.hh:45
static Scalar flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache)
Returns the heat flux within the porous medium (in J/s) across the given sub-control volume face.
Definition: flux/box/fourierslaw.hh:65
Fourier's law for the box scheme in the context of coupled models where coupling occurs across the fa...
Definition: multidomain/facet/box/fourierslaw.hh:54
static Scalar flux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVarsCache &elemFluxVarCache)
Compute the conductive heat flux across a sub-control volume face.
Definition: multidomain/facet/box/fourierslaw.hh:77
static std::vector< Scalar > calculateTransmissibilities(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const FluxVarCache &fluxVarCache, unsigned int phaseIdx)
Definition: multidomain/facet/box/fourierslaw.hh:161
Declares all properties used in Dumux.
This file contains the data which is required to calculate energy fluxes due to molecular diffusion w...