| Line | Branch | Exec | Source |
|---|---|---|---|
| 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 | // SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder | ||
| 5 | // SPDX-License-Identifier: GPL-3.0-or-later | ||
| 6 | // | ||
| 7 | /*! | ||
| 8 | * \file | ||
| 9 | * \ingroup FacetCoupling | ||
| 10 | * \copydoc Dumux::CCTpfaFacetCouplingFouriersLaw | ||
| 11 | */ | ||
| 12 | #ifndef DUMUX_DISCRETIZATION_CC_TPFA_FACET_COUPLING_FOURIERS_LAW_HH | ||
| 13 | #define DUMUX_DISCRETIZATION_CC_TPFA_FACET_COUPLING_FOURIERS_LAW_HH | ||
| 14 | |||
| 15 | #include <array> | ||
| 16 | #include <cmath> | ||
| 17 | |||
| 18 | #include <dune/common/float_cmp.hh> | ||
| 19 | #include <dune/common/fvector.hh> | ||
| 20 | |||
| 21 | #include <dumux/common/math.hh> | ||
| 22 | #include <dumux/common/parameters.hh> | ||
| 23 | #include <dumux/common/properties.hh> | ||
| 24 | |||
| 25 | #include <dumux/discretization/method.hh> | ||
| 26 | #include <dumux/discretization/extrusion.hh> | ||
| 27 | #include <dumux/discretization/cellcentered/tpfa/computetransmissibility.hh> | ||
| 28 | |||
| 29 | #include <dumux/flux/cctpfa/fourierslaw.hh> | ||
| 30 | |||
| 31 | namespace Dumux { | ||
| 32 | |||
| 33 | //! Forward declaration of the implementation | ||
| 34 | template<class TypeTag, | ||
| 35 | bool isNetwork> | ||
| 36 | class CCTpfaFacetCouplingFouriersLawImpl; | ||
| 37 | |||
| 38 | /*! | ||
| 39 | * \ingroup FacetCoupling | ||
| 40 | * \brief Fouriers's law for cell-centered finite volume schemes with two-point flux approximation | ||
| 41 | * in the context of coupled models where the coupling occurs across the facets of the bulk | ||
| 42 | * domain elements with a lower-dimensional domain defined on these facets. | ||
| 43 | */ | ||
| 44 | template<class TypeTag> | ||
| 45 | using CCTpfaFacetCouplingFouriersLaw = | ||
| 46 | CCTpfaFacetCouplingFouriersLawImpl< TypeTag, ( int(GetPropType<TypeTag, Properties::GridGeometry>::GridView::dimension) | ||
| 47 | < int(GetPropType<TypeTag, Properties::GridGeometry>::GridView::dimensionworld) ) >; | ||
| 48 | |||
| 49 | /*! | ||
| 50 | * \ingroup FacetCoupling | ||
| 51 | * \brief Specialization of CCTpfaFacetCouplingFouriersLawImpl for dim=dimWorld | ||
| 52 | */ | ||
| 53 | template<class TypeTag> | ||
| 54 | class CCTpfaFacetCouplingFouriersLawImpl<TypeTag, /*isNetwork*/false> | ||
| 55 | : public FouriersLawImplementation<TypeTag, DiscretizationMethods::CCTpfa> | ||
| 56 | { | ||
| 57 | using Implementation = CCTpfaFacetCouplingFouriersLawImpl<TypeTag, false>; | ||
| 58 | using ParentType = FouriersLawImplementation<TypeTag, DiscretizationMethods::CCTpfa>; | ||
| 59 | |||
| 60 | using Scalar = GetPropType<TypeTag, Properties::Scalar>; | ||
| 61 | using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; | ||
| 62 | using FVElementGeometry = typename GridGeometry::LocalView; | ||
| 63 | using SubControlVolume = typename GridGeometry::SubControlVolume; | ||
| 64 | using SubControlVolumeFace = typename GridGeometry::SubControlVolumeFace; | ||
| 65 | using Extrusion = Extrusion_t<GridGeometry>; | ||
| 66 | |||
| 67 | using GridView = typename GridGeometry::GridView; | ||
| 68 | using Element = typename GridView::template Codim<0>::Entity; | ||
| 69 | |||
| 70 | /*! | ||
| 71 | * \brief The cache class used with this specialization of Fourier's law. | ||
| 72 | */ | ||
| 73 | class FacetCouplingFouriersLawCache | ||
| 74 | { | ||
| 75 | public: | ||
| 76 | //! export the corresponding filler class | ||
| 77 | using Filler = typename ParentType::Cache::Filler; | ||
| 78 | |||
| 79 | //! we store the transmissibilities associated with the interior | ||
| 80 | //! cell, outside cell, and the fracture facet in an array. Access | ||
| 81 | //! to this array should be done using the following indices: | ||
| 82 | static constexpr int insideTijIdx = 0; | ||
| 83 | static constexpr int outsideTijIdx = 1; | ||
| 84 | static constexpr int facetTijIdx = 2; | ||
| 85 | |||
| 86 | //! Export transmissibility storage type | ||
| 87 | using HeatConductionTransmissibilityContainer = std::array<Scalar, 3>; | ||
| 88 | |||
| 89 | //! update subject to a given problem | ||
| 90 | template< class Problem, class ElementVolumeVariables > | ||
| 91 | void updateHeatConduction(const Problem& problem, | ||
| 92 | const Element& element, | ||
| 93 | const FVElementGeometry& fvGeometry, | ||
| 94 | const ElementVolumeVariables& elemVolVars, | ||
| 95 | const SubControlVolumeFace &scvf) | ||
| 96 | { | ||
| 97 | 1706616 | tij_ = Implementation::calculateTransmissibility(problem, element, fvGeometry, elemVolVars, scvf); | |
| 98 | } | ||
| 99 | |||
| 100 | //! We use the same name as in the TpfaFouriersLawCache so | ||
| 101 | //! that this cache and the law implementation for non-coupled | ||
| 102 | //! models can be reused here on facets that do not lie on an | ||
| 103 | //! interior boundary, i.e. do not coincide with a facet element | ||
| 104 | Scalar heatConductionTij() const | ||
| 105 | 2422656 | { return tij_[insideTijIdx]; } | |
| 106 | |||
| 107 | //! returns the transmissibility associated with the inside cell | ||
| 108 | Scalar heatConductionTijInside() const | ||
| 109 | 72576 | { return tij_[insideTijIdx]; } | |
| 110 | |||
| 111 | //! returns the transmissibility associated with the outside cell | ||
| 112 | Scalar heatConductionTijOutside() const | ||
| 113 | 72576 | { return tij_[outsideTijIdx]; } | |
| 114 | |||
| 115 | //! returns the transmissibility associated with the outside cell | ||
| 116 | Scalar heatConductionTijFacet() const | ||
| 117 | 72576 | { return tij_[facetTijIdx]; } | |
| 118 | |||
| 119 | private: | ||
| 120 | HeatConductionTransmissibilityContainer tij_; | ||
| 121 | }; | ||
| 122 | |||
| 123 | public: | ||
| 124 | //! export the type for the corresponding cache | ||
| 125 | using Cache = FacetCouplingFouriersLawCache; | ||
| 126 | |||
| 127 | //! state the discretization method this implementation belongs to | ||
| 128 | using DiscretizationMethod = DiscretizationMethods::CCTpfa; | ||
| 129 | static constexpr DiscretizationMethod discMethod{}; | ||
| 130 | |||
| 131 | //! Compute the conductive heat flux | ||
| 132 | template< class Problem, class ElementVolumeVariables, class ElementFluxVarsCache > | ||
| 133 | 1247616 | static Scalar flux(const Problem& problem, | |
| 134 | const Element& element, | ||
| 135 | const FVElementGeometry& fvGeometry, | ||
| 136 | const ElementVolumeVariables& elemVolVars, | ||
| 137 | const SubControlVolumeFace& scvf, | ||
| 138 | const ElementFluxVarsCache& elemFluxVarsCache) | ||
| 139 | { | ||
| 140 |
6/6✓ Branch 0 taken 1211328 times.
✓ Branch 1 taken 36288 times.
✓ Branch 2 taken 1211328 times.
✓ Branch 3 taken 36288 times.
✓ Branch 4 taken 1211328 times.
✓ Branch 5 taken 36288 times.
|
3742848 | if (!problem.couplingManager().isOnInteriorBoundary(element, scvf)) |
| 141 | 1211328 | return ParentType::flux(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache); | |
| 142 | |||
| 143 | // get inside/outside volume variables | ||
| 144 | 36288 | const auto& fluxVarsCache = elemFluxVarsCache[scvf]; | |
| 145 | 72576 | const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()]; | |
| 146 | 72576 | const auto& facetVolVars = problem.couplingManager().getLowDimVolVars(element, scvf); | |
| 147 | 72576 | const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()]; | |
| 148 | |||
| 149 | // the inside and outside temperatures | ||
| 150 |
1/2✓ Branch 0 taken 36288 times.
✗ Branch 1 not taken.
|
36288 | const Scalar tInside = insideVolVars.temperature(); |
| 151 |
1/2✓ Branch 0 taken 36288 times.
✗ Branch 1 not taken.
|
36288 | const Scalar tFacet = facetVolVars.temperature(); |
| 152 |
1/2✓ Branch 0 taken 36288 times.
✗ Branch 1 not taken.
|
36288 | const Scalar tOutside = outsideVolVars.temperature(); |
| 153 | |||
| 154 |
1/2✓ Branch 0 taken 36288 times.
✗ Branch 1 not taken.
|
36288 | Scalar flux = fluxVarsCache.heatConductionTijInside()*tInside |
| 155 |
1/2✓ Branch 0 taken 36288 times.
✗ Branch 1 not taken.
|
36288 | + fluxVarsCache.heatConductionTijFacet()*tFacet; |
| 156 | |||
| 157 |
1/2✓ Branch 0 taken 36288 times.
✗ Branch 1 not taken.
|
36288 | if (!scvf.boundary()) |
| 158 | 72576 | flux += fluxVarsCache.heatConductionTijOutside()*tOutside; | |
| 159 | return flux; | ||
| 160 | } | ||
| 161 | |||
| 162 | // The flux variables cache has to be bound to an element prior to flux calculations | ||
| 163 | // During the binding, the transmissibility will be computed and stored using the method below. | ||
| 164 | template< class Problem, class ElementVolumeVariables > | ||
| 165 | static typename Cache::HeatConductionTransmissibilityContainer | ||
| 166 | 1706616 | calculateTransmissibility(const Problem& problem, | |
| 167 | const Element& element, | ||
| 168 | const FVElementGeometry& fvGeometry, | ||
| 169 | const ElementVolumeVariables& elemVolVars, | ||
| 170 | const SubControlVolumeFace& scvf) | ||
| 171 | { | ||
| 172 | typename Cache::HeatConductionTransmissibilityContainer tij; | ||
| 173 |
6/6✓ Branch 0 taken 1650672 times.
✓ Branch 1 taken 55944 times.
✓ Branch 2 taken 1650672 times.
✓ Branch 3 taken 55944 times.
✓ Branch 4 taken 1650672 times.
✓ Branch 5 taken 55944 times.
|
5119848 | if (!problem.couplingManager().isCoupled(element, scvf)) |
| 174 | { | ||
| 175 | //! use the standard Fourier's law and only compute one transmissibility | ||
| 176 | 1650672 | tij[Cache::insideTijIdx] = ParentType::calculateTransmissibility(problem, element, fvGeometry, elemVolVars, scvf); | |
| 177 | 1650672 | return tij; | |
| 178 | } | ||
| 179 | |||
| 180 | //! xi factor for the coupling conditions | ||
| 181 |
5/8✓ Branch 0 taken 1 times.
✓ Branch 1 taken 55943 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 1 times.
✗ Branch 10 not taken.
|
55944 | static const Scalar xi = getParamFromGroup<Scalar>(problem.paramGroup(), "FacetCoupling.Xi", 1.0); |
| 182 | |||
| 183 |
2/2✓ Branch 0 taken 27972 times.
✓ Branch 1 taken 27972 times.
|
55944 | const auto insideScvIdx = scvf.insideScvIdx(); |
| 184 |
2/2✓ Branch 0 taken 27972 times.
✓ Branch 1 taken 27972 times.
|
55944 | const auto& insideScv = fvGeometry.scv(insideScvIdx); |
| 185 | 55944 | const auto& insideVolVars = elemVolVars[insideScvIdx]; | |
| 186 | 55944 | const auto wIn = Extrusion::area(fvGeometry, scvf) | |
| 187 | 111888 | *computeTpfaTransmissibility(fvGeometry, scvf, insideScv, | |
| 188 | insideVolVars.effectiveThermalConductivity(), | ||
| 189 | insideVolVars.extrusionFactor()); | ||
| 190 | |||
| 191 | // proceed depending on the interior BC types used | ||
| 192 | 55944 | const auto iBcTypes = problem.interiorBoundaryTypes(element, scvf); | |
| 193 | |||
| 194 | // neumann-coupling | ||
| 195 |
2/4✓ Branch 0 taken 55944 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 55944 times.
✗ Branch 3 not taken.
|
111888 | if (iBcTypes.hasOnlyNeumann()) |
| 196 | { | ||
| 197 | 111888 | const auto& facetVolVars = problem.couplingManager().getLowDimVolVars(element, scvf); | |
| 198 | 55944 | const auto wFacet = 2.0*Extrusion::area(fvGeometry, scvf)*insideVolVars.extrusionFactor() | |
| 199 | 55944 | /facetVolVars.extrusionFactor() | |
| 200 | 223776 | *vtmv(scvf.unitOuterNormal(), | |
| 201 | facetVolVars.effectiveThermalConductivity(), | ||
| 202 | scvf.unitOuterNormal()); | ||
| 203 | |||
| 204 | // The fluxes across this face and the outside face can be expressed in matrix form: | ||
| 205 | // \f$\mathbf{C} \bar{\mathbf{u}} + \mathbf{D} \mathbf{u} + \mathbf{E} \mathbf{u}_\gamma\f$, | ||
| 206 | // where \f$\gamma$\f denotes the domain living on the facets and \f$\bar{\mathbf{u}}$\f are | ||
| 207 | // intermediate face unknowns in the matrix domain. Equivalently, flux continuity reads: | ||
| 208 | // \f$\mathbf{A} \bar{\mathbf{u}} = \mathbf{B} \mathbf{u} + \mathbf{M} \mathbf{u}_\gamma\f$. | ||
| 209 | // Combining the two, we can eliminate the intermediate unknowns and compute the transmissibilities | ||
| 210 | // that allow the description of the fluxes as functions of the cell and Dirichlet temperatures only. | ||
| 211 |
1/2✓ Branch 0 taken 55944 times.
✗ Branch 1 not taken.
|
55944 | if (!scvf.boundary()) |
| 212 | { | ||
| 213 | 55944 | const auto outsideScvIdx = scvf.outsideScvIdx(); | |
| 214 | 55944 | const auto& outsideVolVars = elemVolVars[outsideScvIdx]; | |
| 215 |
2/2✓ Branch 0 taken 27972 times.
✓ Branch 1 taken 27972 times.
|
55944 | const auto wOut = -1.0*Extrusion::area(fvGeometry, scvf) |
| 216 |
4/4✓ Branch 0 taken 27972 times.
✓ Branch 1 taken 27972 times.
✓ Branch 2 taken 27972 times.
✓ Branch 3 taken 27972 times.
|
167832 | *computeTpfaTransmissibility(fvGeometry, scvf, fvGeometry.scv(outsideScvIdx), |
| 217 | outsideVolVars.effectiveThermalConductivity(), | ||
| 218 | outsideVolVars.extrusionFactor()); | ||
| 219 | |||
| 220 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 55944 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 55944 times.
|
55944 | if ( !Dune::FloatCmp::eq(xi, 1.0, 1e-6) ) |
| 221 | { | ||
| 222 | // optimized implementation: factorization obtained using sympy | ||
| 223 | // see CCTpfaFacetCouplingDarcysLaw for more details | ||
| 224 | ✗ | const Scalar factor = wIn * wFacet / ( wIn * wOut * ( 2.0 * xi - 1.0 ) + wFacet * ( xi * ( wIn + wOut ) + wFacet ) ); | |
| 225 | ✗ | tij[Cache::insideTijIdx] = factor * ( wOut * xi + wFacet ); | |
| 226 | ✗ | tij[Cache::outsideTijIdx] = factor * ( wOut * ( 1.0 - xi ) ); | |
| 227 | ✗ | tij[Cache::facetTijIdx] = factor * ( - wOut - wFacet ); | |
| 228 | } | ||
| 229 | else | ||
| 230 | { | ||
| 231 | 55944 | tij[Cache::insideTijIdx] = wFacet*wIn/(wIn+wFacet); | |
| 232 | 111888 | tij[Cache::facetTijIdx] = -tij[Cache::insideTijIdx]; | |
| 233 | 111888 | tij[Cache::outsideTijIdx] = 0.0; | |
| 234 | } | ||
| 235 | } | ||
| 236 | else | ||
| 237 | { | ||
| 238 | ✗ | tij[Cache::insideTijIdx] = wFacet*wIn/(wIn+wFacet); | |
| 239 | ✗ | tij[Cache::facetTijIdx] = -tij[Cache::insideTijIdx]; | |
| 240 | ✗ | tij[Cache::outsideTijIdx] = 0.0; | |
| 241 | } | ||
| 242 | } | ||
| 243 | ✗ | else if (iBcTypes.hasOnlyDirichlet()) | |
| 244 | { | ||
| 245 | ✗ | tij[Cache::insideTijIdx] = wIn; | |
| 246 | ✗ | tij[Cache::outsideTijIdx] = 0.0; | |
| 247 | ✗ | tij[Cache::facetTijIdx] = -wIn; | |
| 248 | } | ||
| 249 | else | ||
| 250 | ✗ | DUNE_THROW(Dune::NotImplemented, "Interior boundary types other than pure Dirichlet or Neumann"); | |
| 251 | |||
| 252 | return tij; | ||
| 253 | } | ||
| 254 | }; | ||
| 255 | |||
| 256 | /*! | ||
| 257 | * \ingroup FacetCoupling | ||
| 258 | * \brief Specialization of CCTpfaFacetCouplingFouriersLawImpl for dim<dimWorld | ||
| 259 | */ | ||
| 260 | template<class TypeTag> | ||
| 261 | class CCTpfaFacetCouplingFouriersLawImpl<TypeTag, /*isNetwork*/true> | ||
| 262 | : public FouriersLawImplementation<TypeTag, DiscretizationMethods::CCTpfa> | ||
| 263 | { | ||
| 264 | using Implementation = CCTpfaFacetCouplingFouriersLawImpl<TypeTag, true>; | ||
| 265 | using ParentType = FouriersLawImplementation<TypeTag, DiscretizationMethods::CCTpfa>; | ||
| 266 | |||
| 267 | using Scalar = GetPropType<TypeTag, Properties::Scalar>; | ||
| 268 | using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; | ||
| 269 | using FVElementGeometry = typename GridGeometry::LocalView; | ||
| 270 | using SubControlVolume = typename GridGeometry::SubControlVolume; | ||
| 271 | using SubControlVolumeFace = typename GridGeometry::SubControlVolumeFace; | ||
| 272 | using Extrusion = Extrusion_t<GridGeometry>; | ||
| 273 | |||
| 274 | using GridView = typename GridGeometry::GridView; | ||
| 275 | using Element = typename GridView::template Codim<0>::Entity; | ||
| 276 | |||
| 277 | /*! | ||
| 278 | * \brief The cache class used with this specialization of Fourier's law. | ||
| 279 | */ | ||
| 280 | class FacetCouplingFouriersLawCache | ||
| 281 | { | ||
| 282 | public: | ||
| 283 | //! export the corresponding filler class | ||
| 284 | using Filler = typename ParentType::Cache::Filler; | ||
| 285 | |||
| 286 | //! we store the transmissibilities associated with the interior | ||
| 287 | //! cell and the fracture facet in an array. Access to this array | ||
| 288 | //! should be done using the following indices: | ||
| 289 | static constexpr int insideTijIdx = 0; | ||
| 290 | static constexpr int facetTijIdx = 1; | ||
| 291 | |||
| 292 | //! Export transmissibility storage type | ||
| 293 | using HeatConductionTransmissibilityContainer = std::array<Scalar, 2>; | ||
| 294 | |||
| 295 | //! update subject to a given problem | ||
| 296 | template< class Problem, class ElementVolumeVariables > | ||
| 297 | void updateHeatConduction(const Problem& problem, | ||
| 298 | const Element& element, | ||
| 299 | const FVElementGeometry& fvGeometry, | ||
| 300 | const ElementVolumeVariables& elemVolVars, | ||
| 301 | const SubControlVolumeFace &scvf) | ||
| 302 | { | ||
| 303 | 1706616 | tij_ = Implementation::calculateTransmissibility(problem, element, fvGeometry, elemVolVars, scvf); | |
| 304 | } | ||
| 305 | |||
| 306 | //! We use the same name as in the TpfaFouriersLawCache so | ||
| 307 | //! that this cache and the law implementation for non-coupled | ||
| 308 | //! models can be reused here on facets that do not lie on an | ||
| 309 | //! interior boundary, i.e. do not coincide with a facet element | ||
| 310 | Scalar heatConductionTij() const | ||
| 311 | 2422656 | { return tij_[insideTijIdx]; } | |
| 312 | |||
| 313 | //! returns the transmissibility associated with the inside cell | ||
| 314 | Scalar heatConductionTijInside() const | ||
| 315 | 72576 | { return tij_[insideTijIdx]; } | |
| 316 | |||
| 317 | //! returns the transmissibility associated with the outside cell | ||
| 318 | Scalar heatConductionTijFacet() const | ||
| 319 | 72576 | { return tij_[facetTijIdx]; } | |
| 320 | |||
| 321 | private: | ||
| 322 | HeatConductionTransmissibilityContainer tij_; | ||
| 323 | }; | ||
| 324 | |||
| 325 | public: | ||
| 326 | //! export the type for the corresponding cache | ||
| 327 | using Cache = FacetCouplingFouriersLawCache; | ||
| 328 | |||
| 329 | //! state the discretization method this implementation belongs to | ||
| 330 | using DiscretizationMethod = DiscretizationMethods::CCTpfa; | ||
| 331 | static constexpr DiscretizationMethod discMethod{}; | ||
| 332 | |||
| 333 | //! Compute the conductive heat flux | ||
| 334 | template< class Problem, class ElementVolumeVariables, class ElementFluxVarsCache > | ||
| 335 | 1247616 | static Scalar flux(const Problem& problem, | |
| 336 | const Element& element, | ||
| 337 | const FVElementGeometry& fvGeometry, | ||
| 338 | const ElementVolumeVariables& elemVolVars, | ||
| 339 | const SubControlVolumeFace& scvf, | ||
| 340 | const ElementFluxVarsCache& elemFluxVarsCache) | ||
| 341 | { | ||
| 342 |
6/6✓ Branch 0 taken 1211328 times.
✓ Branch 1 taken 36288 times.
✓ Branch 2 taken 1211328 times.
✓ Branch 3 taken 36288 times.
✓ Branch 4 taken 1211328 times.
✓ Branch 5 taken 36288 times.
|
3742848 | if (!problem.couplingManager().isOnInteriorBoundary(element, scvf)) |
| 343 | 1211328 | return ParentType::flux(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache); | |
| 344 | |||
| 345 | // get inside/facet volume variables | ||
| 346 | 36288 | const auto& fluxVarsCache = elemFluxVarsCache[scvf]; | |
| 347 | 72576 | const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()]; | |
| 348 | 72576 | const auto& facetVolVars = problem.couplingManager().getLowDimVolVars(element, scvf); | |
| 349 | |||
| 350 | 72576 | return fluxVarsCache.heatConductionTijInside()*insideVolVars.temperature() | |
| 351 | 108864 | + fluxVarsCache.heatConductionTijFacet()*facetVolVars.temperature(); | |
| 352 | } | ||
| 353 | |||
| 354 | // The flux variables cache has to be bound to an element prior to flux calculations | ||
| 355 | // During the binding, the transmissibility will be computed and stored using the method below. | ||
| 356 | template< class Problem, class ElementVolumeVariables > | ||
| 357 | static typename Cache::HeatConductionTransmissibilityContainer | ||
| 358 | 1706616 | calculateTransmissibility(const Problem& problem, | |
| 359 | const Element& element, | ||
| 360 | const FVElementGeometry& fvGeometry, | ||
| 361 | const ElementVolumeVariables& elemVolVars, | ||
| 362 | const SubControlVolumeFace& scvf) | ||
| 363 | { | ||
| 364 | typename Cache::HeatConductionTransmissibilityContainer tij; | ||
| 365 |
6/6✓ Branch 0 taken 1650672 times.
✓ Branch 1 taken 55944 times.
✓ Branch 2 taken 1650672 times.
✓ Branch 3 taken 55944 times.
✓ Branch 4 taken 1650672 times.
✓ Branch 5 taken 55944 times.
|
5119848 | if (!problem.couplingManager().isCoupled(element, scvf)) |
| 366 | { | ||
| 367 | //! use the standard Fourier's law and only compute one transmissibility | ||
| 368 | 1650672 | tij[Cache::insideTijIdx] = ParentType::calculateTransmissibility(problem, element, fvGeometry, elemVolVars, scvf); | |
| 369 | 1650672 | return tij; | |
| 370 | } | ||
| 371 | |||
| 372 | //! xi factor for the coupling conditions | ||
| 373 |
5/8✓ Branch 0 taken 1 times.
✓ Branch 1 taken 55943 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 1 times.
✗ Branch 10 not taken.
|
55944 | static const Scalar xi = getParamFromGroup<Scalar>(problem.paramGroup(), "FacetCoupling.Xi", 1.0); |
| 374 | |||
| 375 | // On surface grids only xi = 1.0 can be used, as the coupling condition | ||
| 376 | // for xi != 1.0 does not generalize for surface grids where the normal | ||
| 377 | // vectors of the inside/outside elements have different orientations. | ||
| 378 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 55944 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 55944 times.
|
55944 | if (Dune::FloatCmp::ne(xi, 1.0, 1e-6)) |
| 379 | ✗ | DUNE_THROW(Dune::InvalidStateException, "Xi != 1.0 cannot be used on surface grids"); | |
| 380 | |||
| 381 |
2/2✓ Branch 0 taken 27972 times.
✓ Branch 1 taken 27972 times.
|
55944 | const auto insideScvIdx = scvf.insideScvIdx(); |
| 382 |
2/2✓ Branch 0 taken 27972 times.
✓ Branch 1 taken 27972 times.
|
55944 | const auto& insideScv = fvGeometry.scv(insideScvIdx); |
| 383 | 55944 | const auto& insideVolVars = elemVolVars[insideScvIdx]; | |
| 384 | 55944 | const auto wIn = Extrusion::area(fvGeometry, scvf) | |
| 385 | 111888 | *computeTpfaTransmissibility(fvGeometry, scvf, insideScv, | |
| 386 | insideVolVars.effectiveThermalConductivity(), | ||
| 387 | insideVolVars.extrusionFactor()); | ||
| 388 | |||
| 389 | // proceed depending on the interior BC types used | ||
| 390 | 55944 | const auto iBcTypes = problem.interiorBoundaryTypes(element, scvf); | |
| 391 | |||
| 392 | // neumann-coupling | ||
| 393 |
2/4✓ Branch 0 taken 55944 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 55944 times.
✗ Branch 3 not taken.
|
111888 | if (iBcTypes.hasOnlyNeumann()) |
| 394 | { | ||
| 395 | // Here we use the square root of the facet extrusion factor | ||
| 396 | // as an approximate average distance from scvf ip to facet center | ||
| 397 | using std::sqrt; | ||
| 398 | 111888 | const auto& facetVolVars = problem.couplingManager().getLowDimVolVars(element, scvf); | |
| 399 | 55944 | const auto wFacet = 2.0*Extrusion::area(fvGeometry, scvf)*insideVolVars.extrusionFactor() | |
| 400 | 55944 | /sqrt(facetVolVars.extrusionFactor()) | |
| 401 | 223776 | *vtmv(scvf.unitOuterNormal(), | |
| 402 | facetVolVars.effectiveThermalConductivity(), | ||
| 403 | scvf.unitOuterNormal()); | ||
| 404 | |||
| 405 | 55944 | tij[Cache::insideTijIdx] = wFacet*wIn/(wIn+wFacet); | |
| 406 | 167832 | tij[Cache::facetTijIdx] = -tij[Cache::insideTijIdx]; | |
| 407 | } | ||
| 408 | ✗ | else if (iBcTypes.hasOnlyDirichlet()) | |
| 409 | { | ||
| 410 | ✗ | tij[Cache::insideTijIdx] = wIn; | |
| 411 | ✗ | tij[Cache::facetTijIdx] = -wIn; | |
| 412 | } | ||
| 413 | else | ||
| 414 | ✗ | DUNE_THROW(Dune::NotImplemented, "Interior boundary types other than pure Dirichlet or Neumann"); | |
| 415 | |||
| 416 | 55944 | return tij; | |
| 417 | } | ||
| 418 | }; | ||
| 419 | |||
| 420 | } // end namespace Dumux | ||
| 421 | |||
| 422 | #endif | ||
| 423 |