24#ifndef DUMUX_DISCRETIZATION_BOX_FACET_COUPLING_FOURIERS_LAW_HH
25#define DUMUX_DISCRETIZATION_BOX_FACET_COUPLING_FOURIERS_LAW_HH
30#include <dune/common/exceptions.hh>
31#include <dune/common/fvector.hh>
32#include <dune/common/float_cmp.hh>
51template<
class TypeTag>
59 using FVElementGeometry =
typename GridGeometry::LocalView;
60 using SubControlVolume =
typename GridGeometry::SubControlVolume;
61 using SubControlVolumeFace =
typename GridGeometry::SubControlVolumeFace;
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;
68 static constexpr int dim = GridView::dimension;
69 static constexpr int dimWorld = GridView::dimensionworld;
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)
85 if (!scvf.interiorBoundary())
86 return ParentType::flux(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarCache);
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");
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];
99 const auto bcTypes = problem.interiorBoundaryTypes(element, scvf);
100 if (bcTypes.hasOnlyNeumann())
103 const auto& facetVolVars = problem.couplingManager().getLowDimVolVars(element, scvf);
107 for (
const auto& scv : scvs(fvGeometry))
108 T += elemVolVars[scv].temperature()*shapeValues[scv.indexInElement()][0];
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);
120 return -1.0*Extrusion::area(scvf)
121 *insideVolVars.extrusionFactor()
122 *
vtmv(scvf.unitOuterNormal(),
123 facetVolVars.effectiveThermalConductivity(),
128 else if (bcTypes.hasOnlyDirichlet())
131 std::vector<Scalar> temperatures(element.subEntities(dim));
132 for (
const auto& scv : scvs(fvGeometry))
133 temperatures[scv.localDofIndex()] = elemVolVars[scv].temperature();
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();
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()));
147 return -1.0*Extrusion::area(scvf)
148 *insideVolVars.extrusionFactor()
149 *
vtmv(scvf.unitOuterNormal(),
150 insideVolVars.effectiveThermalConductivity(),
156 DUNE_THROW(Dune::NotImplemented,
"Mixed boundary types are not supported");
160 template<
class Problem,
class ElementVolumeVariables,
class FluxVarCache>
162 const Element& element,
163 const FVElementGeometry& fvGeometry,
164 const ElementVolumeVariables& elemVolVars,
165 const SubControlVolumeFace& scvf,
166 const FluxVarCache& fluxVarCache,
167 unsigned int phaseIdx)
169 DUNE_THROW(Dune::NotImplemented,
"transmissibilty computation for BoxFacetCouplingFouriersLaw");
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
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...