| 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 DualNetworkCoupling | ||
| 10 | * \ingroup DarcyDarcyCoupling | ||
| 11 | * \brief Coupling manager for equal-dimension boundary coupling | ||
| 12 | */ | ||
| 13 | |||
| 14 | #ifndef DUMUX_MULTIDOMAIN_DUAL_NETWORK_COUPLINGMANAGER_HH | ||
| 15 | #define DUMUX_MULTIDOMAIN_DUAL_NETWORK_COUPLINGMANAGER_HH | ||
| 16 | |||
| 17 | #include <iostream> | ||
| 18 | #include <vector> | ||
| 19 | #include <tuple> | ||
| 20 | |||
| 21 | #include <dune/common/indices.hh> | ||
| 22 | #include <dune/common/exceptions.hh> | ||
| 23 | #include <dune/common/reservedvector.hh> | ||
| 24 | |||
| 25 | #include <dumux/common/properties.hh> | ||
| 26 | #include <dumux/common/math.hh> | ||
| 27 | #include <dumux/common/enumerate.hh> | ||
| 28 | #include <dumux/common/numeqvector.hh> | ||
| 29 | #include <dumux/common/dimensionlessnumbers.hh> | ||
| 30 | #include <dumux/discretization/elementsolution.hh> | ||
| 31 | #include <dumux/discretization/method.hh> | ||
| 32 | #include <dumux/multidomain/couplingmanager.hh> | ||
| 33 | #include "extendedsourcestencil.hh" | ||
| 34 | |||
| 35 | #include "couplingmapper.hh" | ||
| 36 | |||
| 37 | namespace Dumux::PoreNetwork { | ||
| 38 | |||
| 39 | /*! | ||
| 40 | * \ingroup DualNetworkCoupling | ||
| 41 | * \ingroup DarcyDarcyCoupling | ||
| 42 | * \brief Coupling manager for dual network approach for pore network models | ||
| 43 | * \note Concept and algorithms described in Koch et al (2021) https://doi.org/10.1007/s11242-021-01602-5 | ||
| 44 | */ | ||
| 45 | template<class MDTraits> | ||
| 46 | class PNMHeatTransferCouplingManager | ||
| 47 | : public CouplingManager<MDTraits> | ||
| 48 | { | ||
| 49 | using ParentType = CouplingManager<MDTraits>; | ||
| 50 | using ThisType = PNMHeatTransferCouplingManager<MDTraits>; | ||
| 51 | |||
| 52 | using Scalar = typename MDTraits::Scalar; | ||
| 53 | using SolutionVector = typename MDTraits::SolutionVector; | ||
| 54 | |||
| 55 | template<std::size_t i> using SubDomainTypeTag = typename MDTraits::template SubDomain<i>::TypeTag; | ||
| 56 | template<std::size_t i> using Problem = GetPropType<SubDomainTypeTag<i>, Properties::Problem>; | ||
| 57 | template<std::size_t i> using PrimaryVariables = GetPropType<SubDomainTypeTag<i>, Properties::PrimaryVariables>; | ||
| 58 | template<std::size_t i> using NumEqVector = Dumux::NumEqVector<PrimaryVariables<i>>; | ||
| 59 | template<std::size_t i> using GridVolumeVariables = GetPropType<SubDomainTypeTag<i>, Properties::GridVolumeVariables>; | ||
| 60 | template<std::size_t i> using ElementVolumeVariables = typename GridVolumeVariables<i>::LocalView; | ||
| 61 | template<std::size_t i> using ElementFluxVariablesCache = typename GetPropType<SubDomainTypeTag<i>, Properties::GridFluxVariablesCache>::LocalView; | ||
| 62 | template<std::size_t i> using VolumeVariables = typename GridVolumeVariables<i>::VolumeVariables; | ||
| 63 | template<std::size_t i> using GridGeometry = typename MDTraits::template SubDomain<i>::GridGeometry; | ||
| 64 | template<std::size_t i> using FVElementGeometry = typename GridGeometry<i>::LocalView; | ||
| 65 | template<std::size_t i> using SubControlVolumeFace = typename GridGeometry<i>::SubControlVolumeFace; | ||
| 66 | template<std::size_t i> using SubControlVolume = typename GridGeometry<i>::SubControlVolume; | ||
| 67 | template<std::size_t i> using GridView = typename GridGeometry<i>::GridView; | ||
| 68 | template<std::size_t i> using Element = typename GridView<i>::template Codim<0>::Entity; | ||
| 69 | template<std::size_t i> using Indices = typename GetPropType<SubDomainTypeTag<i>, Properties::ModelTraits>::Indices; | ||
| 70 | |||
| 71 | template<std::size_t id> using GridVariables = typename MDTraits::template SubDomain<id>::GridVariables; | ||
| 72 | using GridVariablesTuple = typename MDTraits::template TupleOfSharedPtr<GridVariables>; | ||
| 73 | |||
| 74 | using CouplingMapper = DualNetworkCouplingMapper<Scalar>; | ||
| 75 | |||
| 76 | template<std::size_t i> | ||
| 77 | static constexpr auto domainIdx() | ||
| 78 | { return typename MDTraits::template SubDomain<i>::Index{}; } | ||
| 79 | |||
| 80 | using CouplingStencil = std::vector<std::size_t>; | ||
| 81 | using DimLessNum = DimensionlessNumbers<Scalar>; | ||
| 82 | |||
| 83 | |||
| 84 | public: | ||
| 85 | static constexpr auto solidDomainIdx = typename MDTraits::template SubDomain<0>::Index(); | ||
| 86 | static constexpr auto voidDomainIdx = typename MDTraits::template SubDomain<1>::Index(); | ||
| 87 | |||
| 88 | private: | ||
| 89 | |||
| 90 | using VoidElementSolution = std::decay_t<decltype(elementSolution(std::declval<Element<voidDomainIdx>>(), std::declval<SolutionVector>()[voidDomainIdx], std::declval<GridGeometry<voidDomainIdx>>()))>; | ||
| 91 | using SolidElementSolution = std::decay_t<decltype(elementSolution(std::declval<Element<solidDomainIdx>>(), std::declval<SolutionVector>()[solidDomainIdx], std::declval<GridGeometry<solidDomainIdx>>()))>; | ||
| 92 | |||
| 93 | 2086288 | struct CouplingContextForOneConnection | |
| 94 | { | ||
| 95 | std::size_t connectionGlobalId; | ||
| 96 | std::vector<FVElementGeometry<solidDomainIdx>> solidFVGeometry; // only one needed, but FVElementGeometry can't be default constructed | ||
| 97 | VolumeVariables<solidDomainIdx> solidVolVars; | ||
| 98 | std::vector<FVElementGeometry<voidDomainIdx>> voidFVGeometry; | ||
| 99 | std::vector<ElementVolumeVariables<voidDomainIdx>> voidElemVolVars; | ||
| 100 | std::vector<ElementFluxVariablesCache<voidDomainIdx>> voidElemFluxVarsCache; | ||
| 101 | }; | ||
| 102 | |||
| 103 | using CouplingContextForOnePore = std::pair<std::size_t, std::vector<CouplingContextForOneConnection>>; | ||
| 104 | using CouplingContextForOneElement = Dune::ReservedVector<CouplingContextForOnePore, 2>; | ||
| 105 | |||
| 106 | struct ElementCouplingContext | ||
| 107 | { | ||
| 108 | std::size_t boundElementIndex() const | ||
| 109 | { return boundElementIndex_; } | ||
| 110 | |||
| 111 | ✗ | std::size_t boundDomainId() const | |
| 112 | ✗ | { return domaindId_; } | |
| 113 | |||
| 114 | auto& data() | ||
| 115 | 818168 | { return data_; } | |
| 116 | |||
| 117 | template<std::size_t i> | ||
| 118 | 2534262 | void bind(Dune::index_constant<i> domainI, | |
| 119 | const Element<i>& element, | ||
| 120 | const ThisType& couplingManager) | ||
| 121 | { | ||
| 122 | // do nothing if context is already bound correctly | ||
| 123 | 5068524 | const auto eIdx = couplingManager.gridGeometry(domainI).elementMapper().index(element); | |
| 124 |
5/6✓ Branch 0 taken 1267118 times.
✓ Branch 1 taken 13 times.
✓ Branch 2 taken 1237405 times.
✓ Branch 3 taken 29713 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 1237405 times.
|
2534262 | if (domaindId_ == i && boundElementIndex_ == eIdx && initialized_) |
| 125 | return; | ||
| 126 | |||
| 127 | // do nothing if the element does not couple at all | ||
| 128 | // (this check is probably more expensive, so we do it after the above one) | ||
| 129 | constexpr auto domainJ = Dune::index_constant<1-i>{}; | ||
| 130 |
2/4✓ Branch 1 taken 29726 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 29726 times.
✗ Branch 4 not taken.
|
59452 | if (couplingManager.couplingStencil(domainI, element, domainJ).empty()) |
| 131 | return; | ||
| 132 | |||
| 133 | 59452 | data_.clear(); | |
| 134 | |||
| 135 | // each element has two pores for which we need to bind some connection coupling context | ||
| 136 |
2/2✓ Branch 0 taken 59452 times.
✓ Branch 1 taken 29726 times.
|
178356 | for (int localVertexIdx = 0; localVertexIdx < 2; ++localVertexIdx) |
| 137 | { | ||
| 138 | 237808 | const auto vIdx = couplingManager.gridGeometry(domainI).vertexMapper().subIndex(element, localVertexIdx, 1); | |
| 139 |
2/2✓ Branch 1 taken 54 times.
✓ Branch 2 taken 59398 times.
|
118904 | if (!couplingManager.isCoupledPore(domainI, vIdx)) |
| 140 | 108 | continue; | |
| 141 | |||
| 142 | 118796 | initialized_ = true; | |
| 143 | 118796 | domaindId_ = i; | |
| 144 | 118796 | boundElementIndex_ = eIdx; | |
| 145 | |||
| 146 | 118796 | const auto& connections = [&] | |
| 147 | { | ||
| 148 | if constexpr(domainI == solidDomainIdx) | ||
| 149 | 235548 | return couplingManager.couplingMapper().solidToVoidConnections(vIdx); | |
| 150 | else | ||
| 151 | 120840 | return couplingManager.couplingMapper().voidToSolidConnections(vIdx); | |
| 152 | 237592 | }(); | |
| 153 | 118796 | const auto numConnections = std::distance(connections.begin(), connections.end()); | |
| 154 |
2/4✗ Branch 1 not taken.
✓ Branch 2 taken 59398 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 59398 times.
|
118796 | assert(numConnections == (domainI == solidDomainIdx ? couplingManager.couplingMapper().solidToVoidConnectionIds().at(vIdx).size() : couplingManager.couplingMapper().voidToSolidConnectionIds().at(vIdx).size())); |
| 155 | |||
| 156 | 237592 | data_.push_back(CouplingContextForOnePore{}); | |
| 157 | |||
| 158 | // iterate over all solid/void connections for the given pore | ||
| 159 | 118796 | data_.back().first = vIdx; | |
| 160 | 237592 | data_.back().second.reserve(numConnections); | |
| 161 |
4/4✓ Branch 0 taken 521572 times.
✓ Branch 1 taken 59398 times.
✓ Branch 2 taken 521572 times.
✓ Branch 3 taken 59398 times.
|
1161940 | for (const auto& connection : connections) |
| 162 | { | ||
| 163 | 1043144 | CouplingContextForOneConnection context; | |
| 164 | 1043144 | context.connectionGlobalId = connection.id; | |
| 165 | |||
| 166 | 1043144 | const auto& solidElement = couplingManager.solidGridGeometry_->element(connection.someSolidElementIdx); | |
| 167 | 1043144 | auto solidFVGeometry = localView(*couplingManager.solidGridGeometry_); | |
| 168 |
1/2✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
|
1043144 | solidFVGeometry.bindElement(solidElement); |
| 169 |
1/2✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
|
1043144 | context.solidFVGeometry.push_back(solidFVGeometry); |
| 170 | |||
| 171 |
2/2✓ Branch 0 taken 1043144 times.
✓ Branch 1 taken 521572 times.
|
3129432 | for (const auto& solidScv : scvs(solidFVGeometry)) |
| 172 | { | ||
| 173 |
2/2✓ Branch 0 taken 521572 times.
✓ Branch 1 taken 521572 times.
|
2086288 | if (solidScv.dofIndex() == connection.solidVertexIdx) |
| 174 |
1/2✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
|
1043144 | context.solidVolVars = couplingManager.volVars(solidDomainIdx, solidElement, solidScv); |
| 175 | } | ||
| 176 | |||
| 177 | 1043144 | auto voidFVGeometry = localView(couplingManager.gridGeometry(voidDomainIdx)); | |
| 178 |
3/8✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 521572 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 521572 times.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
|
2086288 | auto voidElemVolVars = localView(couplingManager.gridVariables(voidDomainIdx).curGridVolVars()); |
| 179 |
5/10✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 521572 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 521572 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 502484 times.
✓ Branch 10 taken 19088 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
|
2086288 | auto voidElemFluxVarsCache = localView(couplingManager.gridVariables(voidDomainIdx).gridFluxVarsCache()); |
| 180 | |||
| 181 |
1/2✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
|
1043144 | const auto numConvectionVoidElems = connection.convectionVoidElementIdx.size(); |
| 182 |
1/2✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
|
1043144 | context.voidFVGeometry.reserve(numConvectionVoidElems); |
| 183 |
1/2✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
|
1043144 | context.voidElemVolVars.reserve(numConvectionVoidElems); |
| 184 |
1/2✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
|
1043144 | context.voidElemFluxVarsCache.reserve(numConvectionVoidElems); |
| 185 |
4/4✓ Branch 0 taken 1233344 times.
✓ Branch 1 taken 521572 times.
✓ Branch 2 taken 1233344 times.
✓ Branch 3 taken 521572 times.
|
5596120 | for (const auto voidElemIdx : connection.convectionVoidElementIdx) |
| 186 | { | ||
| 187 |
1/2✓ Branch 1 taken 1233344 times.
✗ Branch 2 not taken.
|
2466688 | const auto voidElem = couplingManager.gridGeometry(voidDomainIdx).element(voidElemIdx); |
| 188 |
1/2✓ Branch 1 taken 1233344 times.
✗ Branch 2 not taken.
|
2466688 | voidFVGeometry.bindElement(voidElem); |
| 189 |
2/4✓ Branch 1 taken 1233344 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1233344 times.
✗ Branch 5 not taken.
|
4933376 | voidElemVolVars.bindElement(voidElem, voidFVGeometry, couplingManager.curSol(voidDomainIdx)); |
| 190 |
1/2✓ Branch 1 taken 1233344 times.
✗ Branch 2 not taken.
|
2466688 | voidElemFluxVarsCache.bindElement(voidElem, voidFVGeometry, voidElemVolVars); |
| 191 | |||
| 192 |
1/2✓ Branch 1 taken 1233344 times.
✗ Branch 2 not taken.
|
2466688 | context.voidFVGeometry.push_back(voidFVGeometry); |
| 193 |
1/2✓ Branch 1 taken 1233344 times.
✗ Branch 2 not taken.
|
2466688 | context.voidElemVolVars.push_back(voidElemVolVars); |
| 194 |
1/2✓ Branch 1 taken 1233344 times.
✗ Branch 2 not taken.
|
2466688 | context.voidElemFluxVarsCache.push_back(voidElemFluxVarsCache); |
| 195 | } | ||
| 196 |
4/6✓ Branch 1 taken 521572 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 521572 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 502484 times.
✓ Branch 7 taken 19088 times.
|
2086288 | data_.back().second.push_back(std::move(context)); |
| 197 | } | ||
| 198 | } | ||
| 199 | } | ||
| 200 | |||
| 201 | 1253354 | const auto& operator[](const std::size_t dofIdx) const | |
| 202 | { | ||
| 203 |
1/2✓ Branch 0 taken 1879782 times.
✗ Branch 1 not taken.
|
4386490 | for (const auto& d : data_) |
| 204 | { | ||
| 205 |
2/2✓ Branch 0 taken 1253354 times.
✓ Branch 1 taken 626428 times.
|
1879782 | if (d.first == dofIdx) |
| 206 | 1253354 | return d.second; | |
| 207 | } | ||
| 208 | ✗ | DUNE_THROW(Dune::InvalidStateException, "No connection found"); | |
| 209 | } | ||
| 210 | |||
| 211 | private: | ||
| 212 | |||
| 213 | bool initialized_ = false; | ||
| 214 | std::size_t domaindId_; | ||
| 215 | std::size_t boundElementIndex_; | ||
| 216 | mutable Dune::ReservedVector<CouplingContextForOnePore, 2> data_; | ||
| 217 | }; | ||
| 218 | |||
| 219 | public: | ||
| 220 | |||
| 221 | //! export traits | ||
| 222 | using MultiDomainTraits = MDTraits; | ||
| 223 | //! export stencil types | ||
| 224 | using CouplingStencils = std::unordered_map<std::size_t, CouplingStencil>; | ||
| 225 | |||
| 226 | /*! | ||
| 227 | * \brief Methods to be accessed by main | ||
| 228 | */ | ||
| 229 | // \{ | ||
| 230 | |||
| 231 | template<class HostGridView, class HostGridData, class VoidGridView, class SolidGridView> | ||
| 232 | 1 | void init(std::shared_ptr<Problem<solidDomainIdx>> solidProblem, | |
| 233 | std::shared_ptr<Problem<voidDomainIdx>> voidProblem, | ||
| 234 | const HostGridView& hostGridView, | ||
| 235 | const HostGridData& hostGridData, | ||
| 236 | const VoidGridView& voidGridView, | ||
| 237 | const SolidGridView& solidGridView, | ||
| 238 | const SolutionVector& curSol) | ||
| 239 | { | ||
| 240 | 2 | voidGridGeometry_ = &voidProblem->gridGeometry(); | |
| 241 | 2 | solidGridGeometry_ = &solidProblem->gridGeometry(); | |
| 242 |
2/4✗ Branch 5 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 1 times.
|
5 | couplingMapper_ = std::make_unique<CouplingMapper>(hostGridView, hostGridData, voidProblem->gridGeometry(), solidProblem->gridGeometry()); |
| 243 | 2 | this->setSubProblems(std::make_tuple(solidProblem, voidProblem)); | |
| 244 | 1 | this->updateSolution(curSol); | |
| 245 | |||
| 246 | 1 | const auto& someVoidElement = voidGridGeometry_->element(0); | |
| 247 | 1 | const auto& someSolidElement = solidGridGeometry_->element(0); | |
| 248 | 2 | voidElementSolution_.update(someVoidElement, curSol[voidDomainIdx], *voidGridGeometry_); | |
| 249 | 2 | solidElementSolution_.update(someSolidElement, curSol[solidDomainIdx], *solidGridGeometry_); | |
| 250 | 1 | setupExtendedStencil(); | |
| 251 | 1 | } | |
| 252 | |||
| 253 | // \} | ||
| 254 | |||
| 255 | /*! | ||
| 256 | * \brief Methods to be accessed by the assembly | ||
| 257 | */ | ||
| 258 | // \{ | ||
| 259 | |||
| 260 | /*! | ||
| 261 | * \brief returns an iterable container of all indices of degrees of freedom of domain j | ||
| 262 | * that couple with / influence the element residual of the given element of domain i | ||
| 263 | * | ||
| 264 | * \param domainI the domain index of domain i | ||
| 265 | * \param element the coupled element of domain í | ||
| 266 | * \param domainJ the domain index of domain j | ||
| 267 | */ | ||
| 268 | template<std::size_t i, std::size_t j> | ||
| 269 | ✗ | const CouplingStencil& couplingStencil(Dune::index_constant<i> domainI, | |
| 270 | const Element<i>& element, | ||
| 271 | Dune::index_constant<j> domainJ) const | ||
| 272 | { | ||
| 273 | ✗ | const auto eIdx = gridGeometry(domainI).elementMapper().index(element); | |
| 274 | |||
| 275 | ✗ | auto getStencils = [this, domainI]() -> const auto& | |
| 276 | { | ||
| 277 | ✗ | return (domainI == solidDomainIdx) ? couplingMapper().solidToVoidStencils() : couplingMapper().voidToSolidStencils(); | |
| 278 | }; | ||
| 279 | |||
| 280 | ✗ | const auto& stencils = getStencils(); | |
| 281 | |||
| 282 | ✗ | return (stencils.count(eIdx)) ? stencils.at(eIdx) : emptyStencil_; | |
| 283 | } | ||
| 284 | |||
| 285 | /*! | ||
| 286 | * \brief Returns whether a given solid grain/void pore body is coupled to the other domain | ||
| 287 | */ | ||
| 288 | template<std::size_t i> | ||
| 289 | 2613832 | bool isCoupledPore(Dune::index_constant<i> domainI, const std::size_t dofIdx) const | |
| 290 | { | ||
| 291 | 10455328 | const auto& isCoupledDof = [&] | |
| 292 | { | ||
| 293 | if constexpr(domainI == solidDomainIdx) | ||
| 294 | 1711518 | return couplingMapper().isCoupledSolidDof(); | |
| 295 | else //voidDomainIdx | ||
| 296 | 902314 | return couplingMapper().isCoupledVoidDof(); | |
| 297 | 2613832 | }(); | |
| 298 | 5227664 | return isCoupledDof[dofIdx]; | |
| 299 | } | ||
| 300 | |||
| 301 | /*! | ||
| 302 | * \brief Returns summed conductive flux between void and solid for one void pore body or one solid grain | ||
| 303 | */ | ||
| 304 | template<std::size_t i> | ||
| 305 | 25504 | Scalar conductionSource(Dune::index_constant<i> domainI, | |
| 306 | const Element<i>& element, | ||
| 307 | const FVElementGeometry<i>& fvGeometry, | ||
| 308 | const ElementVolumeVariables<i>& elemVolVars, | ||
| 309 | const SubControlVolume<i>& scv) const | ||
| 310 | { | ||
| 311 | 25504 | const auto& solidSol = this->curSol(solidDomainIdx); | |
| 312 | 25504 | const auto& voidSol = this->curSol(voidDomainIdx); | |
| 313 | |||
| 314 | 25504 | Scalar source = 0.0; | |
| 315 | 25504 | const auto& connections = [&] | |
| 316 | { | ||
| 317 | if constexpr(domainI == solidDomainIdx) | ||
| 318 | 52344 | return couplingMapper().solidToVoidConnections(scv.dofIndex()); | |
| 319 | else //voidDomainIdx | ||
| 320 | 24168 | return couplingMapper().voidToSolidConnections(scv.dofIndex()); | |
| 321 | 51008 | }(); | |
| 322 | |||
| 323 | // iterate over all connection throats | ||
| 324 |
4/4✓ Branch 0 taken 111842 times.
✓ Branch 1 taken 12752 times.
✓ Branch 2 taken 111842 times.
✓ Branch 3 taken 12752 times.
|
472872 | for (const auto& connection : connections) |
| 325 | { | ||
| 326 | 223684 | const Scalar t = getConnectionTransmissiblity(domainI, connection, elemVolVars, scv); | |
| 327 | const Scalar deltaT = [&] | ||
| 328 | { | ||
| 329 | if constexpr(domainI == solidDomainIdx) | ||
| 330 | 301104 | return solidSol[scv.dofIndex()][Indices<solidDomainIdx>::temperatureIdx] - voidSol[connection.voidVertexIdx][Indices<voidDomainIdx>::temperatureIdx]; | |
| 331 | else //voidDomainIdx | ||
| 332 | 146264 | return voidSol[scv.dofIndex()][Indices<voidDomainIdx>::temperatureIdx] - solidSol[connection.solidVertexIdx][Indices<solidDomainIdx>::temperatureIdx]; | |
| 333 | 223684 | }(); | |
| 334 | |||
| 335 | 223684 | source -= t * deltaT; | |
| 336 | } | ||
| 337 | |||
| 338 | 76512 | source /= (scv.volume() * fvGeometry.gridGeometry().coordinationNumber()[scv.dofIndex()]); | |
| 339 | 25504 | return source; | |
| 340 | } | ||
| 341 | |||
| 342 | template<class Connection, class Context> | ||
| 343 | 12939279 | Scalar convectiveHeatFluxForOneConnection(const Connection& connection, const Context& context) const | |
| 344 | { | ||
| 345 | 12939279 | const auto& solidSol = this->curSol(solidDomainIdx); | |
| 346 | 12939279 | const auto& voidSol = this->curSol(voidDomainIdx); | |
| 347 | 12939279 | const Scalar tSolid = solidSol[connection.solidVertexIdx]; | |
| 348 | 12939279 | Scalar resultConvection = 0.0; | |
| 349 | |||
| 350 | // iterate over all convection void elements | ||
| 351 |
4/4✓ Branch 0 taken 31948674 times.
✓ Branch 1 taken 12939279 times.
✓ Branch 2 taken 31948674 times.
✓ Branch 3 taken 12939279 times.
|
115654464 | for (const auto& [localConvectionVoidElementIdx, convectionVoidElementIdx] : enumerate(connection.convectionVoidElementIdx)) |
| 352 | { | ||
| 353 | 31948674 | const auto& voidElement = voidGridGeometry_->element(convectionVoidElementIdx); | |
| 354 | 31948674 | const auto& voidFVGeometry = context.voidFVGeometry[localConvectionVoidElementIdx]; | |
| 355 | 31948674 | const auto& voidElemVolVars = context.voidElemVolVars[localConvectionVoidElementIdx]; | |
| 356 | 31948674 | const auto& voidElemFluxVarsCache = context.voidElemFluxVarsCache[localConvectionVoidElementIdx]; | |
| 357 | |||
| 358 | 63897348 | const Scalar distance = [&, convectionVoidElementIdx = convectionVoidElementIdx, solidVertexIdx = connection.solidVertexIdx] | |
| 359 | { | ||
| 360 | 31948674 | const auto& throatCenter = this->problem(voidDomainIdx).spatialParams().throatCenter(convectionVoidElementIdx); | |
| 361 |
1/2✓ Branch 0 taken 61174665 times.
✗ Branch 1 not taken.
|
93123339 | for (const auto& solidScv : scvs(context.solidFVGeometry[0])) |
| 362 | { | ||
| 363 | 122349330 | if (solidScv.dofIndex() == solidVertexIdx) | |
| 364 | 95846022 | return (solidScv.dofPosition() - throatCenter).two_norm(); | |
| 365 | } | ||
| 366 | ✗ | DUNE_THROW(Dune::InvalidStateException, "No solid scv found"); | |
| 367 |
3/4✗ Branch 1 not taken.
✓ Branch 2 taken 31948674 times.
✓ Branch 3 taken 31948674 times.
✓ Branch 4 taken 29225991 times.
|
125072013 | }(); |
| 368 | |||
| 369 | using VoidFluxVariables = GetPropType<SubDomainTypeTag<voidDomainIdx>, Properties::FluxVariables>; | ||
| 370 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 31948674 times.
|
31948674 | VoidFluxVariables fluxVars; |
| 371 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 31948674 times.
|
31948674 | const auto& scvf = voidFVGeometry.scvf(0/*localScvfIdx*/); //only one scvf per element -> localScvfIdx = 0 |
| 372 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 31948674 times.
|
31948674 | fluxVars.init(this->problem(voidDomainIdx), voidElement, voidFVGeometry, voidElemVolVars, scvf, voidElemFluxVarsCache); |
| 373 | |||
| 374 | static constexpr auto phaseIdx = 0; | ||
| 375 | 159743370 | const Scalar flux = fluxVars.advectiveFlux(phaseIdx, [phaseIdx = phaseIdx](const auto& volVars){ return volVars.mobility(phaseIdx);}); | |
| 376 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
|
31948674 | const Scalar velocity = flux / voidElemFluxVarsCache[scvf].throatCrossSectionalArea(phaseIdx); |
| 377 | |||
| 378 | 31948674 | const auto tFluidMean = [&] | |
| 379 | { | ||
| 380 | ✗ | Scalar result = 0.0; | |
| 381 | ✗ | const auto numScv = voidFVGeometry.numScv(); | |
| 382 | assert(numScv == 2); | ||
| 383 | ✗ | for (const auto& voidScv : scvs(voidFVGeometry)) | |
| 384 | ✗ | result += 1.0/numScv * voidSol[voidScv.dofIndex()][Indices<voidDomainIdx>::temperatureIdx]; | |
| 385 | ✗ | return result; | |
| 386 | }; | ||
| 387 | |||
| 388 | 31948674 | const auto tFluidUpstream = [&] | |
| 389 | { | ||
| 390 | ✗ | const auto upstreamDofIdx = flux > 0.0 ? voidFVGeometry.scv(scvf.insideScvIdx()).dofIndex() : voidFVGeometry.scv(scvf.outsideScvIdx()).dofIndex(); | |
| 391 | ✗ | return voidSol[upstreamDofIdx][Indices<voidDomainIdx>::temperatureIdx]; | |
| 392 | }; | ||
| 393 | |||
| 394 | enum class FluidTemperatureMode {mean, self, upstream}; | ||
| 395 |
3/4✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
|
31948675 | static const auto fluidTemperatureMode = [&] |
| 396 | { | ||
| 397 |
3/6✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
|
1 | static const auto mode = getParam<std::string>("DualNetwork.FluidTemperatureMode", "mean"); |
| 398 | 2 | std::cout << "Using FluidTemperatureMode " << mode << std::endl; | |
| 399 |
1/2✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
|
1 | if (mode == "mean") |
| 400 | return FluidTemperatureMode::mean; | ||
| 401 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 1 times.
|
1 | else if (mode == "self") |
| 402 | return FluidTemperatureMode::self; | ||
| 403 | ✗ | else if (mode == "upstream") | |
| 404 | return FluidTemperatureMode::upstream; | ||
| 405 | else | ||
| 406 | ✗ | DUNE_THROW(Dune::IOError, mode << " is not a valid FluidTemperatureMode"); | |
| 407 |
1/2✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
|
2 | }(); |
| 408 | |||
| 409 | 63897348 | const Scalar tFluid = [&, voidVertexIdx = connection.voidVertexIdx] | |
| 410 | { | ||
| 411 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 31948674 times.
|
31948674 | if (fluidTemperatureMode == FluidTemperatureMode::mean) |
| 412 | ✗ | return tFluidMean(); | |
| 413 |
1/2✓ Branch 0 taken 31948674 times.
✗ Branch 1 not taken.
|
31948674 | else if (fluidTemperatureMode == FluidTemperatureMode::self) |
| 414 | 63897348 | return voidSol[voidVertexIdx][Indices<voidDomainIdx>::temperatureIdx]; | |
| 415 | else | ||
| 416 | ✗ | return tFluidUpstream(); | |
| 417 | 31948674 | }(); | |
| 418 | |||
| 419 |
8/8✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 31948673 times.
✓ Branch 4 taken 1 times.
✓ Branch 5 taken 31948673 times.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 31948673 times.
|
127794696 | const Scalar meanKinematicViscosity = 0.5*(voidElemVolVars[0].viscosity(phaseIdx)/voidElemVolVars[0].density(phaseIdx) |
| 420 |
8/8✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 31948673 times.
✓ Branch 4 taken 1 times.
✓ Branch 5 taken 31948673 times.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 31948673 times.
|
127794696 | + voidElemVolVars[1].viscosity(phaseIdx)/voidElemVolVars[1].density(phaseIdx)); |
| 421 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
|
31948674 | const Scalar characteristicLength = 2.0*voidElemFluxVarsCache[scvf].throatInscribedRadius(); |
| 422 | using std::abs; | ||
| 423 |
4/4✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 31948673 times.
|
63897348 | const Scalar Re = DimLessNum::reynoldsNumber(abs(velocity), characteristicLength, meanKinematicViscosity); |
| 424 | |||
| 425 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
|
31948674 | static const Scalar fixedLambda = getParam<Scalar>("DualNetwork.FixedConvectionLambda", -1.0); |
| 426 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
|
31948674 | static const Scalar lambaFactor = getParam<Scalar>("DualNetwork.ConvectionLambaFactor", 0.9); |
| 427 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
|
31948674 | static const Scalar lambaExponent = getParam<Scalar>("DualNetwork.ConvectionLambaExponent", 0.4); |
| 428 | using std::pow; | ||
| 429 |
1/2✓ Branch 0 taken 31948674 times.
✗ Branch 1 not taken.
|
31948674 | const Scalar lambda = fixedLambda > 0.0 ? fixedLambda : 1.0 + lambaFactor*pow(Re, lambaExponent); //approximation: see eq.30 in Koch et al (2021) https://doi.org/10.1007/s11242-021-01602-5 |
| 430 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
|
31948674 | const Scalar selfA = connection.connectionArea / connection.convectionVoidElementIdx.size(); |
| 431 | |||
| 432 | 95846022 | const auto neighborA = [&] | |
| 433 | { | ||
| 434 | // get the area associated to the other void dof | ||
| 435 |
1/2✓ Branch 0 taken 47696623 times.
✗ Branch 1 not taken.
|
79645297 | for (const auto& voidScv : scvs(voidFVGeometry)) |
| 436 | { | ||
| 437 |
2/2✓ Branch 0 taken 31948674 times.
✓ Branch 1 taken 15747949 times.
|
47696623 | if (voidScv.dofIndex() != connection.voidVertexIdx) |
| 438 | { | ||
| 439 | 95846022 | const auto& otherConnections = couplingMapper().voidToSolidConnections(voidScv.dofIndex()); | |
| 440 |
2/4✓ Branch 0 taken 177395633 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 177395633 times.
✗ Branch 3 not taken.
|
177395633 | for (const auto& otherConn : otherConnections) |
| 441 | { | ||
| 442 |
2/2✓ Branch 0 taken 31948674 times.
✓ Branch 1 taken 145446959 times.
|
177395633 | if (otherConn.solidVertexIdx == connection.solidVertexIdx) |
| 443 | { | ||
| 444 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 31948674 times.
|
31948674 | if (otherConn.voidVertexIdx != voidScv.dofIndex()) |
| 445 | ✗ | DUNE_THROW(Dune::InvalidStateException, "Void dofs don't match"); | |
| 446 | |||
| 447 | 63897348 | return otherConn.connectionArea/otherConn.convectionVoidElementIdx.size(); | |
| 448 | } | ||
| 449 | } | ||
| 450 | } | ||
| 451 | } | ||
| 452 | ✗ | DUNE_THROW(Dune::InvalidStateException, "No neighbor area found"); | |
| 453 | }; | ||
| 454 | |||
| 455 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
|
31948674 | static const bool useAvgA = getParam<bool>("DualNetwork.UseAverageConvectionArea", false); |
| 456 |
1/2✓ Branch 0 taken 31948674 times.
✗ Branch 1 not taken.
|
31948674 | const Scalar A = useAvgA ? 0.5*(selfA + neighborA()) : selfA; |
| 457 | |||
| 458 |
2/2✓ Branch 0 taken 12458322 times.
✓ Branch 1 taken 19490352 times.
|
31948674 | const Scalar deltaT = (elementCouplingContext_.boundDomainId() == voidDomainIdx) ? (tSolid - tFluid) : (tFluid - tSolid); |
| 459 | |||
| 460 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 31948673 times.
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
|
31948674 | static const int verbose = getParam<int>("DualNetwork.SourceVerboseForDof", -1); |
| 461 |
1/6✗ Branch 0 not taken.
✓ Branch 1 taken 31948674 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
31948674 | if (verbose >= 0 && (connection.voidVertexIdx == verbose || connection.solidVertexIdx == verbose)) |
| 462 | { | ||
| 463 | ✗ | std::cout << "\n" << std::endl; | |
| 464 | ✗ | const auto domain = elementCouplingContext_.boundDomainId() == solidDomainIdx ? "solid" : "void"; | |
| 465 | ✗ | std::cout << "At " << domain << ", dof " << verbose << ", flow elem " << convectionVoidElementIdx << ", globalId " << connection.id << std::endl; | |
| 466 | ✗ | std::cout << std::scientific << std::setprecision(15) << "velocity " << velocity << ", Re " << Re << ", delTaT " << deltaT << ", result " << lambda*A/distance*deltaT << std::endl; | |
| 467 | ✗ | std::cout << std::scientific << std::setprecision(15) << "A " << A << ", distance " << distance << std::endl; | |
| 468 | ✗ | std::cout << std::endl; | |
| 469 | } | ||
| 470 | |||
| 471 | 31948674 | resultConvection += lambda*A/distance * deltaT; | |
| 472 | } | ||
| 473 | |||
| 474 | 12939279 | return resultConvection; | |
| 475 | } | ||
| 476 | |||
| 477 | /*! | ||
| 478 | * \brief Returns summed conductive heat fluxes for one void pore body coupled to solid grains (or the other way around) | ||
| 479 | * that occur due to convection in void throats | ||
| 480 | */ | ||
| 481 | template<std::size_t i> | ||
| 482 | 2493956 | Scalar convectionSource(Dune::index_constant<i> domainI, | |
| 483 | const Element<i>& element, | ||
| 484 | const FVElementGeometry<i>& fvGeometry, | ||
| 485 | const ElementVolumeVariables<i>& elemVolVars, | ||
| 486 | const SubControlVolume<i>& scv) const | ||
| 487 | { | ||
| 488 | 4987912 | bindCouplingContext(domainI, element); | |
| 489 | |||
| 490 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 1246976 times.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
|
2493956 | static const int verbose = getParam<int>("DualNetwork.SourceVerboseForDof", -1); |
| 491 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1246978 times.
|
2493956 | if (scv.dofIndex() == verbose) |
| 492 | ✗ | std::cout << "Start Source at elemn " << fvGeometry.gridGeometry().elementMapper().index(element) << " *******************************" << std::endl; | |
| 493 | |||
| 494 | // iterate over all connection throats | ||
| 495 | 2493956 | const auto& connections = [&] | |
| 496 | { | ||
| 497 | if constexpr (domainI == solidDomainIdx) | ||
| 498 | 4895766 | return couplingMapper().solidToVoidConnections(scv.dofIndex()); | |
| 499 | else | ||
| 500 | 2586102 | return couplingMapper().voidToSolidConnections(scv.dofIndex()); | |
| 501 | 4987912 | }(); | |
| 502 | |||
| 503 | 2493956 | Scalar source = 0.0; | |
| 504 | 2493956 | const auto& contextForPore = elementCouplingContext_[scv.dofIndex()]; | |
| 505 | |||
| 506 |
4/4✓ Branch 0 taken 12932091 times.
✓ Branch 1 taken 1246978 times.
✓ Branch 2 taken 12932091 times.
✓ Branch 3 taken 1246978 times.
|
59210232 | for (const auto& [localConnectionIdx, connection] : enumerate(connections)) |
| 507 | 51728364 | source += convectiveHeatFluxForOneConnection(connection, contextForPore[localConnectionIdx]); | |
| 508 | |||
| 509 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1246978 times.
|
2493956 | if (scv.dofIndex() == verbose) |
| 510 | ✗ | std::cout << std::scientific << std::setprecision(15) << "total conv source " << source << "\n\n ********************" << std::endl; | |
| 511 | |||
| 512 | 4987912 | source /= (scv.volume() * fvGeometry.gridGeometry().coordinationNumber()[scv.dofIndex()]); | |
| 513 | |||
| 514 | 2493956 | return source; | |
| 515 | } | ||
| 516 | |||
| 517 | //TODO: this should not be in the general coupling manager | ||
| 518 | template<std::size_t i, std::size_t j> | ||
| 519 | Scalar sourceWithFixedTransmissibility(Dune::index_constant<i> domainI, | ||
| 520 | const Element<i>& element, | ||
| 521 | const FVElementGeometry<i>& fvGeometry, | ||
| 522 | const ElementVolumeVariables<i>& elemVolVars, | ||
| 523 | const SubControlVolume<i>& scv, | ||
| 524 | Dune::index_constant<j> domainJ) const | ||
| 525 | { | ||
| 526 | const auto& voidSol = this->curSol(voidDomainIdx); | ||
| 527 | const auto& solidSol = this->curSol(solidDomainIdx); | ||
| 528 | |||
| 529 | Scalar source = 0.0; | ||
| 530 | |||
| 531 | // iterate over all connection throats | ||
| 532 | for (const auto& connection : couplingMapper().voidToSolidConnections(scv.dofIndex())) | ||
| 533 | { | ||
| 534 | const Scalar t = this->problem(voidDomainIdx).getInternalReferenceHeatTransmissibilityCoupling(); | ||
| 535 | const Scalar deltaT = [&] | ||
| 536 | { | ||
| 537 | if constexpr(domainI == solidDomainIdx) | ||
| 538 | return solidSol[scv.dofIndex()][Indices<solidDomainIdx>::temperatureIdx] - voidSol[connection.voidVertexIdx][Indices<voidDomainIdx>::temperatureIdx]; | ||
| 539 | else //voidDomainIdx | ||
| 540 | return voidSol[scv.dofIndex()][Indices<voidDomainIdx>::temperatureIdx] - solidSol[connection.solidVertexIdx][Indices<solidDomainIdx>::temperatureIdx]; | ||
| 541 | }(); | ||
| 542 | |||
| 543 | source -= t * deltaT; | ||
| 544 | } | ||
| 545 | |||
| 546 | source /= (scv.volume() * fvGeometry.gridGeometry().coordinationNumber()[scv.dofIndex()]); | ||
| 547 | return source; | ||
| 548 | } | ||
| 549 | |||
| 550 | template<std::size_t i, class Connection> | ||
| 551 | 223684 | Scalar getConnectionTransmissiblity(Dune::index_constant<i> domainI, | |
| 552 | const Connection& connection, | ||
| 553 | const ElementVolumeVariables<i>& elemVolVars, | ||
| 554 | const SubControlVolume<i>& scv) const | ||
| 555 | { | ||
| 556 | static constexpr bool isVoid = (domainI == voidDomainIdx); | ||
| 557 | |||
| 558 | 111842 | const auto poreRadiusVoid = [&] | |
| 559 | { | ||
| 560 | 111842 | static const bool useExactPoreRadiusVoid = getParam<bool>("DualNetwork.UseExactPoreRadiusVoid", false); | |
| 561 | 111842 | if (useExactPoreRadiusVoid) | |
| 562 | { | ||
| 563 | using std::sqrt; | ||
| 564 | ✗ | static const Scalar R = getParam<Scalar>("DualNetwork.SphereRadius", 50e-6); | |
| 565 | ✗ | static const Scalar overlapFactor = getParam<Scalar>("DualNetwork.OverlapFactor"); | |
| 566 | ✗ | static const Scalar dx = overlapFactor*R; | |
| 567 | ✗ | static const Scalar r = sqrt(3.0) * dx - R; | |
| 568 | ✗ | return r; | |
| 569 | } | ||
| 570 | else | ||
| 571 | 111842 | return gridGeometry(voidDomainIdx).poreInscribedRadius(connection.voidVertexIdx); | |
| 572 | 335526 | }(); | |
| 573 | |||
| 574 | 111842 | const auto poreRadiusSolid = [&] | |
| 575 | { | ||
| 576 | 111842 | static const bool useExactPoreRadiusSolid = getParam<bool>("DualNetwork.UseExactPoreRadiusSolid", false); | |
| 577 | 111842 | if (useExactPoreRadiusSolid) | |
| 578 | { | ||
| 579 | ✗ | static const Scalar R = getParam<Scalar>("DualNetwork.SphereRadius", 50e-6); | |
| 580 | ✗ | return R; | |
| 581 | } | ||
| 582 | else | ||
| 583 | 335526 | return this->problem(solidDomainIdx).spatialParams().poreExtendedRadius(connection.solidVertexIdx); | |
| 584 | 335526 | }(); | |
| 585 | |||
| 586 | 75276 | const Scalar fluidThermalConductivity = [&] | |
| 587 | { | ||
| 588 | if constexpr (isVoid) | ||
| 589 | 73132 | return elemVolVars[scv].effectiveThermalConductivity(); | |
| 590 | else | ||
| 591 | { | ||
| 592 | 283100 | const auto& voidElement = voidGridGeometry_->element(connection.someVoidElementIdx); | |
| 593 | 75276 | auto voidFVGeometry = localView(*voidGridGeometry_); | |
| 594 | 75276 | voidFVGeometry.bindElement(voidElement); | |
| 595 | |||
| 596 | 132548 | for (const auto& voidScv : scvs(voidFVGeometry)) | |
| 597 | { | ||
| 598 | 132548 | if (voidScv.dofIndex() == connection.voidVertexIdx) | |
| 599 | 150552 | return volVars(voidDomainIdx, voidElement, voidScv).effectiveThermalConductivity(); | |
| 600 | } | ||
| 601 | |||
| 602 | ✗ | DUNE_THROW(Dune::InvalidStateException, "No scv found"); | |
| 603 | } | ||
| 604 | 374236 | }(); | |
| 605 | |||
| 606 | 36566 | const Scalar solidThermalConductivity = [&] | |
| 607 | { | ||
| 608 | if constexpr (!isVoid) //solid | ||
| 609 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 75275 times.
|
150552 | return elemVolVars[scv].effectiveThermalConductivity(); |
| 610 | else | ||
| 611 | { | ||
| 612 | 141938 | const auto& solidElement = solidGridGeometry_->element(connection.someSolidElementIdx); | |
| 613 | 36566 | auto solidFVGeometry = localView(*solidGridGeometry_); | |
| 614 | 36566 | solidFVGeometry.bindElement(solidElement); | |
| 615 | |||
| 616 | 68806 | for (const auto& solidScv : scvs(solidFVGeometry)) | |
| 617 | { | ||
| 618 | 68806 | if (solidScv.dofIndex() == connection.solidVertexIdx) | |
| 619 | 73132 | return volVars(solidDomainIdx, solidElement, solidScv).effectiveThermalConductivity(); | |
| 620 | } | ||
| 621 | |||
| 622 | ✗ | DUNE_THROW(Dune::InvalidStateException, "No scv found"); | |
| 623 | } | ||
| 624 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 75275 times.
|
296816 | }(); |
| 625 | |||
| 626 | 223684 | const Scalar kappa = fluidThermalConductivity / solidThermalConductivity; | |
| 627 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 111840 times.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
|
223684 | static const Scalar Nu = getParam<Scalar>("DualNetwork.Nu", 1.0); |
| 628 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 111840 times.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
|
223684 | static const Scalar Bi = getParam<Scalar>("DualNetwork.Bi", 1.0); |
| 629 | |||
| 630 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 111840 times.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
|
223684 | static const bool useExactConnectionLength = getParam<bool>("DualNetwork.UseExactConnectionLength", false); |
| 631 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 111842 times.
|
223684 | const Scalar length = useExactConnectionLength ? poreRadiusSolid + poreRadiusVoid : connection.connectionLength; |
| 632 | |||
| 633 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 111840 times.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
|
223684 | static const bool useExactConnectionAreaSphere = getParam<bool>("DualNetwork.UseExactConnectionAreaSphere", false); |
| 634 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 111840 times.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
|
223684 | static const Scalar connectionAreaShapeFactor = getParam<Scalar>("DualNetwork.ConnectionAreaShapeFactor", 0.9); |
| 635 | 111842 | const Scalar area = [&]() | |
| 636 | { | ||
| 637 | 111842 | if (useExactConnectionAreaSphere) | |
| 638 | { | ||
| 639 | ✗ | static const Scalar R = getParam<Scalar>("DualNetwork.SphereRadius", 50e-6); | |
| 640 | ✗ | static const Scalar overlapFactor = getParam<Scalar>("DualNetwork.OverlapFactor"); | |
| 641 | ✗ | static const auto dx = overlapFactor*R; | |
| 642 | ✗ | static const auto h = R - dx; | |
| 643 | ✗ | static const auto interfacialArea = 4*M_PI*R*R - 6*(2*M_PI*R*h); | |
| 644 | ✗ | assert(std::abs(length - std::sqrt(3.0) * dx) < 1e-14); // TODO remove | |
| 645 | ✗ | return interfacialArea/8.0*connectionAreaShapeFactor; | |
| 646 | } | ||
| 647 | else | ||
| 648 | 111842 | return connection.connectionArea*connectionAreaShapeFactor; | |
| 649 | 223684 | }(); | |
| 650 | //distance weighted harmonic mean of thermal conductivities | ||
| 651 | 223684 | const Scalar meanThermalConductivity = (fluidThermalConductivity*Bi*(poreRadiusSolid + poreRadiusVoid))/(poreRadiusSolid*kappa + poreRadiusVoid*Bi/Nu); | |
| 652 | 223684 | return area / length * meanThermalConductivity; | |
| 653 | } | ||
| 654 | |||
| 655 | // \} | ||
| 656 | |||
| 657 | /*! | ||
| 658 | * \brief Return the volume variables of domain i for a given element and scv | ||
| 659 | */ | ||
| 660 | template<std::size_t i> | ||
| 661 | 150552 | VolumeVariables<i> volVars(Dune::index_constant<i> domainI, | |
| 662 | const Element<i>& element, | ||
| 663 | const SubControlVolume<i>& scv) const | ||
| 664 | { | ||
| 665 |
0/2✗ Branch 1 not taken.
✗ Branch 2 not taken.
|
150552 | VolumeVariables<i> volVars; |
| 666 |
0/4✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
301104 | const auto elemSol = elementSolution(element, this->curSol(domainI), gridGeometry(domainI)); |
| 667 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 75276 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
|
150552 | volVars.update(elemSol, this->problem(domainI), element, scv); |
| 668 | 150552 | return volVars; | |
| 669 | } | ||
| 670 | |||
| 671 | template<std::size_t i, std::size_t j, class LocalAssemblerI> | ||
| 672 | 973616 | decltype(auto) evalCouplingResidual(Dune::index_constant<i> domainI, | |
| 673 | const LocalAssemblerI& localAssemblerI, | ||
| 674 | Dune::index_constant<j> domainJ, | ||
| 675 | std::size_t dofIdxGlobalJ) | ||
| 676 | { | ||
| 677 | static_assert(i != j, "A domain cannot be coupled to itself!"); | ||
| 678 | |||
| 679 | 973616 | typename LocalAssemblerI::LocalResidual::ElementResidualVector residual; | |
| 680 | |||
| 681 | 973616 | const auto& element = localAssemblerI.element(); | |
| 682 | 973616 | const auto& fvGeometry = localAssemblerI.fvGeometry(); | |
| 683 | 973616 | const auto& curElemVolVars = localAssemblerI.curElemVolVars(); | |
| 684 | |||
| 685 | 1947232 | residual.resize(fvGeometry.numScv()); | |
| 686 |
2/2✓ Branch 0 taken 973616 times.
✓ Branch 1 taken 486808 times.
|
3894464 | for (const auto& scv : scvs(fvGeometry)) |
| 687 | { | ||
| 688 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 973616 times.
|
1947232 | auto couplingSource = this->problem(domainI).source(element, fvGeometry, curElemVolVars, scv); |
| 689 | 3894464 | couplingSource *= -scv.volume()*curElemVolVars[scv].extrusionFactor(); | |
| 690 | 3894464 | residual[scv.indexInElement()] = couplingSource; | |
| 691 | } | ||
| 692 | |||
| 693 | 973616 | return residual; | |
| 694 | } | ||
| 695 | |||
| 696 | //! Bind the coupling context for a low dim element TODO remove Assembler | ||
| 697 | template<std::size_t i, class Assembler = int> | ||
| 698 | ✗ | void bindCouplingContext(Dune::index_constant<i> domainI, const Element<i>& element, const Assembler& assembler = 0) const | |
| 699 |
5/14✓ Branch 5 taken 1007 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 2184 times.
✗ Branch 9 not taken.
✓ Branch 11 taken 2184 times.
✗ Branch 12 not taken.
✓ Branch 14 taken 1007 times.
✗ Branch 15 not taken.
✓ Branch 17 taken 1007 times.
✗ Branch 18 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
|
1267131 | { elementCouplingContext_.bind(domainI, element, *this); } |
| 700 | |||
| 701 | /*! | ||
| 702 | * \brief Update the coupling context | ||
| 703 | */ | ||
| 704 | template<std::size_t i, class LocalAssemblerI, std::size_t j, class PriVars> | ||
| 705 | 1636336 | void updateCouplingContext(Dune::index_constant<i> domainI, | |
| 706 | const LocalAssemblerI& localAssemblerI, | ||
| 707 | Dune::index_constant<j> domainJ, | ||
| 708 | const std::size_t dofIdxGlobalJ, | ||
| 709 | const PriVars& priVars, | ||
| 710 | int pvIdxJ) | ||
| 711 | { | ||
| 712 | 5941728 | this->curSol(domainJ)[dofIdxGlobalJ][pvIdxJ] = priVars[pvIdxJ]; | |
| 713 | |||
| 714 | // each element has two pores for which we need to bind some connection coupling context | ||
| 715 |
2/2✓ Branch 0 taken 1635728 times.
✓ Branch 1 taken 818168 times.
|
8180464 | for (auto& context : elementCouplingContext_.data()) |
| 716 | { | ||
| 717 |
4/4✓ Branch 0 taken 17045136 times.
✓ Branch 1 taken 1635728 times.
✓ Branch 2 taken 17045136 times.
✓ Branch 3 taken 1635728 times.
|
43904640 | for (auto& connInfo : context.second) |
| 718 | { | ||
| 719 |
1/2✓ Branch 0 taken 13877856 times.
✗ Branch 1 not taken.
|
34090272 | const auto& staticConnectionInfo = couplingMapper().connectionInfo()[connInfo.connectionGlobalId]; |
| 720 | |||
| 721 | // when the solid dof is deflected, just updated the solid volVars | ||
| 722 | if constexpr (domainJ == solidDomainIdx) | ||
| 723 | { | ||
| 724 | 6334560 | const auto& solidElement = gridGeometry(domainJ).element(staticConnectionInfo.someSolidElementIdx); | |
| 725 |
2/2✓ Branch 0 taken 6334560 times.
✓ Branch 1 taken 3167280 times.
|
31672800 | for (const auto& scv : scvs(connInfo.solidFVGeometry.front())) |
| 726 | { | ||
| 727 |
6/6✓ Branch 0 taken 556664 times.
✓ Branch 1 taken 5777896 times.
✓ Branch 2 taken 556664 times.
✓ Branch 3 taken 5777896 times.
✓ Branch 4 taken 556664 times.
✓ Branch 5 taken 5777896 times.
|
38007360 | solidElementSolution_[scv.localDofIndex()] = this->curSol(domainJ)[scv.dofIndex()]; |
| 728 |
2/2✓ Branch 0 taken 556664 times.
✓ Branch 1 taken 5777896 times.
|
12669120 | if (scv.dofIndex() == dofIdxGlobalJ) |
| 729 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 556664 times.
|
1113328 | connInfo.solidVolVars.update(solidElementSolution_, this->problem(domainJ), solidElement, scv); |
| 730 | } | ||
| 731 | } | ||
| 732 | else // deflection of void dofs | ||
| 733 | { | ||
| 734 |
3/6✓ Branch 0 taken 13877856 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 13877856 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 13877856 times.
✗ Branch 5 not taken.
|
83267136 | assert(staticConnectionInfo.convectionVoidElementIdx.size() == connInfo.voidFVGeometry.size()); |
| 735 |
4/4✓ Branch 0 taken 34084704 times.
✓ Branch 1 taken 13877856 times.
✓ Branch 2 taken 34084704 times.
✓ Branch 3 taken 13877856 times.
|
123680832 | for (int voidElemLocalIdx = 0; voidElemLocalIdx < staticConnectionInfo.convectionVoidElementIdx.size(); ++voidElemLocalIdx) |
| 736 | { | ||
| 737 | 68169408 | const auto eIdx = staticConnectionInfo.convectionVoidElementIdx[voidElemLocalIdx]; | |
| 738 | 68169408 | const auto& element = gridGeometry(voidDomainIdx).element(eIdx); | |
| 739 | 68169408 | const auto& fvGeometry = connInfo.voidFVGeometry[voidElemLocalIdx]; | |
| 740 | 68169408 | auto& elemVolVars = connInfo.voidElemVolVars[voidElemLocalIdx]; | |
| 741 | |||
| 742 |
2/2✓ Branch 0 taken 68169408 times.
✓ Branch 1 taken 34084704 times.
|
272677632 | for (const auto& scv : scvs(fvGeometry)) |
| 743 | { | ||
| 744 |
2/2✓ Branch 0 taken 4219776 times.
✓ Branch 1 taken 63949632 times.
|
136338816 | if (scv.dofIndex() == dofIdxGlobalJ) |
| 745 | { | ||
| 746 | 25318656 | voidElementSolution_[scv.localDofIndex()] = this->curSol(voidDomainIdx)[scv.dofIndex()]; | |
| 747 |
1/2✗ Branch 1 not taken.
✓ Branch 2 taken 4219776 times.
|
8439552 | getVolVarAccess_(domainJ, gridVars_(voidDomainIdx).curGridVolVars(), elemVolVars, scv).update(voidElementSolution_, this->problem(voidDomainIdx), element, scv); |
| 748 | } | ||
| 749 | } | ||
| 750 | |||
| 751 |
2/2✓ Branch 0 taken 34084704 times.
✓ Branch 1 taken 34084704 times.
|
204508224 | for (const auto& scvf : scvfs(fvGeometry)) |
| 752 | { | ||
| 753 | if constexpr (getPropValue<SubDomainTypeTag<j>, Properties::EnableGridFluxVariablesCache>()) | ||
| 754 | gridVars_(voidDomainIdx).gridFluxVarsCache().cache(eIdx, scvf.index()).update(problem(voidDomainIdx), element, fvGeometry, elemVolVars, scvf); | ||
| 755 | else | ||
| 756 |
3/6✗ Branch 0 not taken.
✓ Branch 1 taken 34084704 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 34084704 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 34084704 times.
|
204508224 | connInfo.voidElemFluxVarsCache[voidElemLocalIdx][scvf].update(this->problem(voidDomainIdx), element, fvGeometry, elemVolVars, scvf); |
| 757 | } | ||
| 758 | } | ||
| 759 | } | ||
| 760 | } | ||
| 761 | } | ||
| 762 | 1636336 | } | |
| 763 | |||
| 764 | const auto& couplingContext() const | ||
| 765 |
2/4✓ Branch 1 taken 4362 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2014 times.
✗ Branch 5 not taken.
|
6376 | { return elementCouplingContext_; } |
| 766 | |||
| 767 | |||
| 768 | template<std::size_t i> | ||
| 769 | ✗ | const auto& gridGeometry(Dune::index_constant<i> domainI) const | |
| 770 | { | ||
| 771 | if constexpr (i == voidDomainIdx) | ||
| 772 | ✗ | return *voidGridGeometry_; | |
| 773 | else | ||
| 774 | ✗ | return *solidGridGeometry_; | |
| 775 | } | ||
| 776 | |||
| 777 | const CouplingMapper& couplingMapper() const | ||
| 778 |
16/32✓ Branch 0 taken 3695232 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 3695232 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 10182624 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 10182624 times.
✗ Branch 7 not taken.
✓ Branch 51 taken 3594 times.
✗ Branch 52 not taken.
✓ Branch 54 taken 3594 times.
✗ Branch 55 not taken.
✓ Branch 57 taken 3594 times.
✗ Branch 58 not taken.
✓ Branch 60 taken 3594 times.
✗ Branch 61 not taken.
✓ Branch 63 taken 3594 times.
✗ Branch 64 not taken.
✓ Branch 66 taken 3594 times.
✗ Branch 67 not taken.
✗ Branch 68 not taken.
✓ Branch 69 taken 1 times.
✗ Branch 70 not taken.
✓ Branch 71 taken 1 times.
✓ Branch 73 taken 4362 times.
✗ Branch 74 not taken.
✓ Branch 76 taken 4362 times.
✗ Branch 77 not taken.
✓ Branch 79 taken 2014 times.
✗ Branch 80 not taken.
✓ Branch 82 taken 2014 times.
✗ Branch 83 not taken.
|
103392822 | { return *couplingMapper_; } |
| 779 | |||
| 780 | /*! | ||
| 781 | * \brief set the pointers to the grid variables | ||
| 782 | * \param gridVariables A tuple of shared pointers to the grid variables | ||
| 783 | */ | ||
| 784 | void setGridVariables(GridVariablesTuple&& gridVariables) | ||
| 785 | 1 | { gridVariables_ = gridVariables; } | |
| 786 | |||
| 787 | /*! | ||
| 788 | * \brief set a pointer to one of the grid variables | ||
| 789 | * \param gridVariables a pointer to the grid variables | ||
| 790 | * \param domainIdx the domain index of the grid variables | ||
| 791 | */ | ||
| 792 | template<class GridVariables, std::size_t i> | ||
| 793 | void setGridVariables(std::shared_ptr<GridVariables> gridVariables, Dune::index_constant<i> domainIdx) | ||
| 794 | { std::get<i>(gridVariables_) = gridVariables; } | ||
| 795 | |||
| 796 | /*! | ||
| 797 | * \brief Return a reference to the grid variables of a sub problem | ||
| 798 | * \param domainIdx The domain index | ||
| 799 | */ | ||
| 800 | template<std::size_t i> | ||
| 801 | 1043144 | const GridVariables<i>& gridVariables(Dune::index_constant<i> domainIdx) const | |
| 802 | { | ||
| 803 |
2/4✓ Branch 0 taken 1043144 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1043144 times.
✗ Branch 3 not taken.
|
2086288 | if (std::get<i>(gridVariables_)) |
| 804 | 3129432 | return *std::get<i>(gridVariables_); | |
| 805 | else | ||
| 806 | ✗ | DUNE_THROW(Dune::InvalidStateException, "The gridVariables pointer was not set. Use setGridVariables() before calling this function"); | |
| 807 | } | ||
| 808 | |||
| 809 | /*! | ||
| 810 | * \brief extend the jacobian pattern of the diagonal block of domain i | ||
| 811 | * by those entries that are not already in the uncoupled pattern | ||
| 812 | * \note Such additional dependencies can arise from the coupling, e.g. if a coupling source | ||
| 813 | * term depends on a non-local average of a quantity of the same domain | ||
| 814 | */ | ||
| 815 | template<std::size_t id, class JacobianPattern> | ||
| 816 | ✗ | void extendJacobianPattern(Dune::index_constant<id> domainI, JacobianPattern& pattern) const | |
| 817 | { | ||
| 818 | ✗ | extendedSourceStencil_.extendJacobianPattern(*this, domainI, pattern); | |
| 819 | ✗ | } | |
| 820 | |||
| 821 | /*! | ||
| 822 | * \brief evaluate additional derivatives of the element residual of a domain with respect | ||
| 823 | * to dofs in the same domain that are not in the regular stencil (per default this is not the case) | ||
| 824 | * \note Such additional dependencies can arise from the coupling, e.g. if a coupling source | ||
| 825 | * term depends on a non-local average of a quantity of the same domain | ||
| 826 | * \note This is the same for box and cc | ||
| 827 | */ | ||
| 828 | template<std::size_t i, class LocalAssemblerI, class JacobianMatrixDiagBlock, class GridVariables> | ||
| 829 | ✗ | void evalAdditionalDomainDerivatives(Dune::index_constant<i> domainI, | |
| 830 | const LocalAssemblerI& localAssemblerI, | ||
| 831 | const typename LocalAssemblerI::LocalResidual::ElementResidualVector&, | ||
| 832 | JacobianMatrixDiagBlock& A, | ||
| 833 | GridVariables& gridVariables) | ||
| 834 | { | ||
| 835 | // std::cout << "calling evalAdditionalDomainDerivatives " << std::endl; | ||
| 836 | 25528 | extendedSourceStencil_.evalAdditionalDomainDerivatives(*this, domainI, localAssemblerI, this->curSol(domainI), A, gridVariables); //TODO: call without curSol, but currently protected .../dumux/dumux/multidomain/couplingmanager.hh | |
| 837 | ✗ | } | |
| 838 | |||
| 839 | //! Return a reference to an empty stencil | ||
| 840 | std::vector<std::size_t>& emptyStencil() | ||
| 841 | { return emptyStencil_; } | ||
| 842 | |||
| 843 | |||
| 844 | //! Return a reference to an empty stencil | ||
| 845 | template<std::size_t i> | ||
| 846 | ✗ | const std::vector<std::size_t>& emptyStencil(Dune::index_constant<i> domainI) const | |
| 847 |
1/2✓ Branch 0 taken 8736 times.
✗ Branch 1 not taken.
|
10920 | { return emptyStencil_; } |
| 848 | |||
| 849 | template<std::size_t i> | ||
| 850 | ✗ | const auto& gridView(Dune::index_constant<i> domainI) const | |
| 851 | { | ||
| 852 | 4 | return gridGeometry(domainI).gridView(); | |
| 853 | } | ||
| 854 | |||
| 855 | |||
| 856 | 1 | void setupExtendedStencil() | |
| 857 | { | ||
| 858 | 1 | extendedSourceStencil_.clear(); | |
| 859 | |||
| 860 |
4/4✓ Branch 2 taken 1007 times.
✓ Branch 3 taken 1 times.
✓ Branch 4 taken 1007 times.
✓ Branch 5 taken 1 times.
|
2016 | for (const auto& element : elements(gridView(voidDomainIdx))) |
| 861 | { | ||
| 862 | 2014 | const auto eIdx = gridGeometry(voidDomainIdx).elementMapper().index(element); | |
| 863 | 2014 | const auto vIdx0 = gridGeometry(voidDomainIdx).vertexMapper().subIndex(element, 0, 1); | |
| 864 | 2014 | const auto vIdx1 = gridGeometry(voidDomainIdx).vertexMapper().subIndex(element, 1, 1); | |
| 865 | |||
| 866 |
9/12✓ Branch 3 taken 1007 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1007 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 27915 times.
✓ Branch 9 taken 1007 times.
✓ Branch 10 taken 27915 times.
✓ Branch 11 taken 1007 times.
✓ Branch 12 taken 9632 times.
✓ Branch 13 taken 18283 times.
✓ Branch 15 taken 27915 times.
✗ Branch 16 not taken.
|
29929 | for (const auto& is : intersections(gridView(voidDomainIdx), element)) |
| 867 | { | ||
| 868 |
4/4✓ Branch 0 taken 9632 times.
✓ Branch 1 taken 18283 times.
✓ Branch 2 taken 9632 times.
✓ Branch 3 taken 18283 times.
|
55830 | if (is.neighbor()) |
| 869 | { | ||
| 870 |
1/2✓ Branch 1 taken 9632 times.
✗ Branch 2 not taken.
|
9632 | const auto& outsideElement = is.outside(); |
| 871 |
2/2✓ Branch 0 taken 19264 times.
✓ Branch 1 taken 9632 times.
|
28896 | for (int i = 0; i < 2; ++i) |
| 872 | { | ||
| 873 |
2/4✓ Branch 1 taken 19264 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 19264 times.
✗ Branch 5 not taken.
|
38528 | const auto outsideVertexIdx = gridGeometry(voidDomainIdx).vertexMapper().subIndex(outsideElement, i, 1); |
| 874 |
2/2✓ Branch 0 taken 9632 times.
✓ Branch 1 taken 9632 times.
|
19264 | if (outsideVertexIdx != vIdx0 && outsideVertexIdx != vIdx1) |
| 875 |
3/6✓ Branch 1 taken 9632 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 9632 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 9632 times.
✗ Branch 8 not taken.
|
19264 | extendedSourceStencil_.stencil()[eIdx].push_back(outsideVertexIdx); |
| 876 | } | ||
| 877 | } | ||
| 878 | } | ||
| 879 | |||
| 880 | 2014 | removeDuplicates_(extendedSourceStencil_.stencil()[eIdx]); | |
| 881 | } | ||
| 882 | 1 | } | |
| 883 | |||
| 884 | |||
| 885 | protected: | ||
| 886 | |||
| 887 | template<std::size_t i> | ||
| 888 | ✗ | VolumeVariables<i>& getVolVarAccess_(Dune::index_constant<i> domainIdx, GridVolumeVariables<i>& gridVolVars, ElementVolumeVariables<i>& elemVolVars, const SubControlVolume<i>& scv) | |
| 889 | { | ||
| 890 | if constexpr (getPropValue<SubDomainTypeTag<i>, Properties::EnableGridVolumeVariablesCache>()) | ||
| 891 | return gridVolVars.volVars(scv); | ||
| 892 | else | ||
| 893 |
4/8✗ Branch 0 not taken.
✓ Branch 1 taken 1489024 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1489024 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2730752 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 2730752 times.
|
8439552 | return elemVolVars[scv]; |
| 894 | } | ||
| 895 | |||
| 896 | /*! | ||
| 897 | * \brief Return a reference to the grid variables of a sub problem | ||
| 898 | * \param domainIdx The domain index | ||
| 899 | */ | ||
| 900 | template<std::size_t i> | ||
| 901 | 4219776 | GridVariables<i>& gridVars_(Dune::index_constant<i> domainIdx) | |
| 902 | { | ||
| 903 |
2/4✓ Branch 0 taken 4219776 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4219776 times.
✗ Branch 3 not taken.
|
8439552 | if (std::get<i>(gridVariables_)) |
| 904 | 12659328 | return *std::get<i>(gridVariables_); | |
| 905 | else | ||
| 906 | ✗ | DUNE_THROW(Dune::InvalidStateException, "The gridVariables pointer was not set. Use setGridVariables() before calling this function"); | |
| 907 | } | ||
| 908 | |||
| 909 | //! Removes duplicate entries from the coupling stencils | ||
| 910 | ✗ | void removeDuplicates_(std::vector<std::size_t>& stencil) | |
| 911 | { | ||
| 912 | ✗ | std::sort(stencil.begin(), stencil.end()); | |
| 913 | ✗ | stencil.erase(std::unique(stencil.begin(), stencil.end()), stencil.end()); | |
| 914 | ✗ | } | |
| 915 | |||
| 916 | |||
| 917 | bool convectiveHeatTransfer_; | ||
| 918 | |||
| 919 | std::vector<std::size_t> emptyStencil_; | ||
| 920 | std::unique_ptr<const CouplingMapper> couplingMapper_; | ||
| 921 | const GridGeometry<voidDomainIdx>* voidGridGeometry_; | ||
| 922 | const GridGeometry<solidDomainIdx>* solidGridGeometry_; | ||
| 923 | |||
| 924 | VoidElementSolution voidElementSolution_; | ||
| 925 | SolidElementSolution solidElementSolution_; | ||
| 926 | mutable ElementCouplingContext elementCouplingContext_; | ||
| 927 | |||
| 928 | /*! | ||
| 929 | * \brief A tuple of std::shared_ptrs to the grid variables of the sub problems | ||
| 930 | */ | ||
| 931 | GridVariablesTuple gridVariables_; | ||
| 932 | |||
| 933 | //! the extended source stencil object | ||
| 934 | PoreNetwork::PNMHeatExtendedSourceStencil<ThisType> extendedSourceStencil_; | ||
| 935 | }; | ||
| 936 | |||
| 937 | }; // end namespace Dumux::PoreNetwork | ||
| 938 | |||
| 939 | #endif | ||
| 940 |