| 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::FacetCouplingManager | ||
| 11 | */ | ||
| 12 | #ifndef DUMUX_BOX_FACETCOUPLING_MANAGER_HH | ||
| 13 | #define DUMUX_BOX_FACETCOUPLING_MANAGER_HH | ||
| 14 | |||
| 15 | #include <algorithm> | ||
| 16 | #include <cassert> | ||
| 17 | |||
| 18 | #include <dumux/common/properties.hh> | ||
| 19 | #include <dumux/common/numeqvector.hh> | ||
| 20 | #include <dumux/discretization/method.hh> | ||
| 21 | #include <dumux/discretization/elementsolution.hh> | ||
| 22 | #include <dumux/multidomain/couplingmanager.hh> | ||
| 23 | #include <dumux/multidomain/facet/couplingmanager.hh> | ||
| 24 | |||
| 25 | namespace Dumux { | ||
| 26 | |||
| 27 | /*! | ||
| 28 | * \ingroup FacetCoupling | ||
| 29 | * \brief Manages the coupling between bulk elements and lower dimensional elements | ||
| 30 | * where the coupling occurs across the facets of the bulk grid. This implementation | ||
| 31 | * is to be used in conjunction with models using the box scheme in the bulk domain. | ||
| 32 | * | ||
| 33 | * \tparam MDTraits The multidomain traits containing the types on all sub-domains | ||
| 34 | * \tparam CouplingMapper Class containing maps on the coupling between dofs of different grids | ||
| 35 | * \tparam bulkDomainId The domain id of the bulk problem | ||
| 36 | * \tparam lowDimDomainId The domain id of the lower-dimensional problem | ||
| 37 | */ | ||
| 38 | template<class MDTraits, class CouplingMapper, std::size_t bulkDomainId, std::size_t lowDimDomainId> | ||
| 39 | class FacetCouplingManager<MDTraits, CouplingMapper, bulkDomainId, lowDimDomainId, DiscretizationMethods::Box> | ||
| 40 | : public virtual CouplingManager< MDTraits > | ||
| 41 | { | ||
| 42 | using ParentType = CouplingManager< MDTraits >; | ||
| 43 | |||
| 44 | // convenience aliases and instances of the two domain ids | ||
| 45 | using BulkIdType = typename MDTraits::template SubDomain<bulkDomainId>::Index; | ||
| 46 | using LowDimIdType = typename MDTraits::template SubDomain<lowDimDomainId>::Index; | ||
| 47 | static constexpr auto bulkId = BulkIdType(); | ||
| 48 | static constexpr auto lowDimId = LowDimIdType(); | ||
| 49 | |||
| 50 | // the sub-domain type tags | ||
| 51 | template<std::size_t id> using SubDomainTypeTag = typename MDTraits::template SubDomain<id>::TypeTag; | ||
| 52 | |||
| 53 | // further types specific to the sub-problems | ||
| 54 | template<std::size_t id> using PrimaryVariables = GetPropType<SubDomainTypeTag<id>, Properties::PrimaryVariables>; | ||
| 55 | template<std::size_t id> using Problem = GetPropType<SubDomainTypeTag<id>, Properties::Problem>; | ||
| 56 | template<std::size_t id> using NumEqVector = Dumux::NumEqVector<PrimaryVariables<id>>; | ||
| 57 | template<std::size_t id> using ElementBoundaryTypes = GetPropType<SubDomainTypeTag<id>, Properties::ElementBoundaryTypes>; | ||
| 58 | template<std::size_t id> using LocalResidual = typename MDTraits::template SubDomain<id>::LocalResidual; | ||
| 59 | |||
| 60 | template<std::size_t id> using GridGeometry = GetPropType<SubDomainTypeTag<id>, Properties::GridGeometry>; | ||
| 61 | template<std::size_t id> using FVElementGeometry = typename GridGeometry<id>::LocalView; | ||
| 62 | template<std::size_t id> using SubControlVolume = typename GridGeometry<id>::SubControlVolume; | ||
| 63 | template<std::size_t id> using SubControlVolumeFace = typename GridGeometry<id>::SubControlVolumeFace; | ||
| 64 | template<std::size_t id> using GridView = typename GridGeometry<id>::GridView; | ||
| 65 | template<std::size_t id> using Element = typename GridView<id>::template Codim<0>::Entity; | ||
| 66 | template<std::size_t id> using GridIndexType = typename GridView<id>::IndexSet::IndexType; | ||
| 67 | |||
| 68 | template<std::size_t id> using GridVariables = GetPropType<SubDomainTypeTag<id>, Properties::GridVariables>; | ||
| 69 | template<std::size_t id> using GridVolumeVariables = typename GridVariables<id>::GridVolumeVariables; | ||
| 70 | template<std::size_t id> using ElementVolumeVariables = typename GridVolumeVariables<id>::LocalView; | ||
| 71 | template<std::size_t id> using VolumeVariables = typename ElementVolumeVariables<id>::VolumeVariables; | ||
| 72 | template<std::size_t id> using GridFluxVariablesCache = typename GridVariables<id>::GridFluxVariablesCache; | ||
| 73 | template<std::size_t id> using ElementFluxVariablesCache = typename GridFluxVariablesCache<id>::LocalView; | ||
| 74 | |||
| 75 | // extract corresponding grid ids from the mapper | ||
| 76 | static constexpr int bulkDim = GridView<bulkDomainId>::dimension; | ||
| 77 | static constexpr int lowDimDim = GridView<lowDimDomainId>::dimension; | ||
| 78 | static constexpr auto bulkGridId = CouplingMapper::template gridId<bulkDim>(); | ||
| 79 | static constexpr auto lowDimGridId = CouplingMapper::template gridId<lowDimDim>(); | ||
| 80 | |||
| 81 | static constexpr bool lowDimUsesBox = GridGeometry<lowDimId>::discMethod == DiscretizationMethods::box; | ||
| 82 | |||
| 83 | /*! | ||
| 84 | * \brief The coupling context of the bulk domain. Contains all data of the lower- | ||
| 85 | * dimensional domain which is required for the computation of a bulk element | ||
| 86 | * residual. This boils down to the geometries of all lower-dimensional elements | ||
| 87 | * connected to a given bulk element, as well as their element volume variables. | ||
| 88 | */ | ||
| 89 | struct BulkCouplingContext | ||
| 90 | { | ||
| 91 | bool isSet; | ||
| 92 | GridIndexType< bulkId > elementIdx; | ||
| 93 | std::vector< FVElementGeometry<lowDimId> > lowDimFvGeometries; | ||
| 94 | std::vector< std::vector<VolumeVariables<lowDimId>> > lowDimElemVolVars; | ||
| 95 | |||
| 96 | 2435900 | void reset() | |
| 97 | { | ||
| 98 |
2/2✓ Branch 0 taken 164288 times.
✓ Branch 1 taken 1188176 times.
|
2435900 | lowDimFvGeometries.clear(); |
| 99 | 2435900 | lowDimElemVolVars.clear(); | |
| 100 | 2435900 | isSet = false; | |
| 101 | 2435900 | } | |
| 102 | }; | ||
| 103 | |||
| 104 | /*! | ||
| 105 | * \brief The coupling context of the lower-dimensional (codim 1) domain. Contains | ||
| 106 | * all data of the bulk domain which is required for computation of element | ||
| 107 | * residuals in the lower-dimensional domain. This is essentially everything | ||
| 108 | * that is necessary to compute the fluxes in bulk domain entering a given | ||
| 109 | * lower-dimensional element. Thus, we store and bind the local views of the | ||
| 110 | * neighboring elements. | ||
| 111 | * | ||
| 112 | * \note We need unique ptrs here because the local views have no default constructor. | ||
| 113 | */ | ||
| 114 | struct LowDimCouplingContext | ||
| 115 | { | ||
| 116 | private: | ||
| 117 | // For dim == dimworld in bulk domain, we know there is always going | ||
| 118 | // to be two bulk neighbors for a lower-dimensional element. Choose | ||
| 119 | // the container type depending on this | ||
| 120 | static constexpr int bulkDimWorld = GridView<bulkId>::dimensionworld; | ||
| 121 | static constexpr bool bulkIsSurfaceGrid = bulkDim != bulkDimWorld; | ||
| 122 | |||
| 123 | template< class T > | ||
| 124 | using ContainerType = typename std::conditional< bulkIsSurfaceGrid, | ||
| 125 | std::vector< T >, | ||
| 126 | std::array< T, 2 > >::type; | ||
| 127 | |||
| 128 | public: | ||
| 129 | bool isSet; | ||
| 130 | GridIndexType< lowDimId > elementIdx; | ||
| 131 | ContainerType< std::unique_ptr<FVElementGeometry<bulkId>> > bulkFvGeometries; | ||
| 132 | ContainerType< std::unique_ptr<ElementVolumeVariables<bulkId>> > bulkElemVolVars; | ||
| 133 | ContainerType< std::unique_ptr<ElementFluxVariablesCache<bulkId>> > bulkElemFluxVarsCache; | ||
| 134 | std::unique_ptr< LocalResidual<bulkId> > bulkLocalResidual; | ||
| 135 | ContainerType< ElementBoundaryTypes<bulkId> > bulkElemBcTypes; | ||
| 136 | |||
| 137 | 106202 | void reset() | |
| 138 | { | ||
| 139 | 106202 | isSet = false; | |
| 140 |
12/12✓ Branch 0 taken 101758 times.
✓ Branch 1 taken 6 times.
✓ Branch 2 taken 101758 times.
✓ Branch 3 taken 6 times.
✓ Branch 4 taken 101758 times.
✓ Branch 5 taken 6 times.
✓ Branch 6 taken 8753 times.
✓ Branch 7 taken 175 times.
✓ Branch 8 taken 8753 times.
✓ Branch 9 taken 175 times.
✓ Branch 10 taken 8753 times.
✓ Branch 11 taken 175 times.
|
1274175 | auto resetPtr = [&] (auto&& uniquePtr) { uniquePtr.reset(); }; |
| 141 | 318606 | std::for_each(bulkFvGeometries.begin(), bulkFvGeometries.end(), resetPtr); | |
| 142 | 318606 | std::for_each(bulkElemVolVars.begin(), bulkElemVolVars.end(), resetPtr); | |
| 143 | 318606 | std::for_each(bulkElemFluxVarsCache.begin(), bulkElemFluxVarsCache.end(), resetPtr); | |
| 144 |
2/2✓ Branch 0 taken 55311 times.
✓ Branch 1 taken 23 times.
|
106202 | bulkLocalResidual.reset(); |
| 145 | 106202 | } | |
| 146 | |||
| 147 | // resize containers for surface grids | ||
| 148 | template<bool s = bulkIsSurfaceGrid, std::enable_if_t<s, int> = 0> | ||
| 149 | 1248 | void resize(std::size_t numEmbedments) | |
| 150 | { | ||
| 151 | 1248 | bulkFvGeometries.resize(numEmbedments); | |
| 152 | 1248 | bulkElemVolVars.resize(numEmbedments); | |
| 153 | 1248 | bulkElemFluxVarsCache.resize(numEmbedments); | |
| 154 | 1248 | bulkElemBcTypes.resize(numEmbedments); | |
| 155 | 1248 | } | |
| 156 | |||
| 157 | // no memory reservation necessary for non-surface grids | ||
| 158 | template<bool s = bulkIsSurfaceGrid, std::enable_if_t<!s, int> = 0> | ||
| 159 | ✗ | void resize(std::size_t numEmbedments) | |
| 160 | ✗ | {} | |
| 161 | }; | ||
| 162 | |||
| 163 | public: | ||
| 164 | |||
| 165 | //! types used for coupling stencils | ||
| 166 | template<std::size_t i, std::size_t j = (i == bulkId) ? lowDimId : bulkId> | ||
| 167 | using CouplingStencilType = typename CouplingMapper::template Stencil< CouplingMapper::template gridId<GridView<j>::dimension>() >; | ||
| 168 | |||
| 169 | //! the type of the solution vector | ||
| 170 | using SolutionVector = typename MDTraits::SolutionVector; | ||
| 171 | |||
| 172 | /*! | ||
| 173 | * \brief Initialize the coupling manager. | ||
| 174 | * | ||
| 175 | * \param bulkProblem The problem to be solved on the bulk domain | ||
| 176 | * \param lowDimProblem The problem to be solved on the lower-dimensional domain | ||
| 177 | * \param couplingMapper The mapper object containing the connectivity between the domains | ||
| 178 | * \param curSol The current solution | ||
| 179 | */ | ||
| 180 | 27 | void init(std::shared_ptr< Problem<bulkId> > bulkProblem, | |
| 181 | std::shared_ptr< Problem<lowDimId> > lowDimProblem, | ||
| 182 | std::shared_ptr< CouplingMapper > couplingMapper, | ||
| 183 | const SolutionVector& curSol) | ||
| 184 | { | ||
| 185 | 27 | couplingMapperPtr_ = couplingMapper; | |
| 186 | |||
| 187 | // set the sub problems | ||
| 188 |
2/4✓ Branch 1 taken 23 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 23 times.
✗ Branch 4 not taken.
|
27 | this->setSubProblem(bulkProblem, bulkId); |
| 189 |
2/4✓ Branch 1 taken 23 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 23 times.
✗ Branch 4 not taken.
|
27 | this->setSubProblem(lowDimProblem, lowDimId); |
| 190 | |||
| 191 | // copy the solution vector | ||
| 192 | 27 | ParentType::updateSolution(curSol); | |
| 193 | |||
| 194 | // determine all bulk elements that couple to low dim elements | ||
| 195 | 54 | const auto& bulkMap = couplingMapperPtr_->couplingMap(bulkGridId, lowDimGridId); | |
| 196 | 108 | bulkElemIsCoupled_.assign(bulkProblem->gridGeometry().gridView().size(0), false); | |
| 197 | 81 | std::for_each( bulkMap.begin(), | |
| 198 | bulkMap.end(), | ||
| 199 | 7428 | [&] (const auto& entry) { bulkElemIsCoupled_[entry.first] = true; }); | |
| 200 | 27 | } | |
| 201 | |||
| 202 | /*! | ||
| 203 | * \brief The coupling stencil of a given bulk domain element. | ||
| 204 | */ | ||
| 205 | 338458 | const CouplingStencilType<bulkId>& couplingStencil(BulkIdType domainI, | |
| 206 | const Element<bulkId>& element, | ||
| 207 | LowDimIdType domainJ) const | ||
| 208 | { | ||
| 209 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 328383 times.
|
338458 | const auto eIdx = this->problem(domainI).gridGeometry().elementMapper().index(element); |
| 210 | |||
| 211 |
4/4✓ Branch 0 taken 11989 times.
✓ Branch 1 taken 316394 times.
✓ Branch 2 taken 11989 times.
✓ Branch 3 taken 316394 times.
|
676916 | if (bulkElemIsCoupled_[eIdx]) |
| 212 | { | ||
| 213 | 25738 | const auto& map = couplingMapperPtr_->couplingMap(bulkGridId, lowDimGridId); | |
| 214 | 12869 | auto it = map.find(eIdx); | |
| 215 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 11989 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 11989 times.
|
25738 | assert(it != map.end()); |
| 216 | 25738 | return it->second.couplingStencil; | |
| 217 | } | ||
| 218 | |||
| 219 | 651178 | return getEmptyStencil(lowDimId); | |
| 220 | } | ||
| 221 | |||
| 222 | /*! | ||
| 223 | * \brief The coupling stencil of the lower-dimensional domain with the bulk domain. | ||
| 224 | */ | ||
| 225 | 6578 | const CouplingStencilType<lowDimId>& couplingStencil(LowDimIdType domainI, | |
| 226 | const Element<lowDimId>& element, | ||
| 227 | BulkIdType domainJ) const | ||
| 228 | { | ||
| 229 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 6118 times.
|
6578 | const auto eIdx = this->problem(lowDimId).gridGeometry().elementMapper().index(element); |
| 230 | |||
| 231 | 13156 | const auto& map = couplingMapperPtr_->couplingMap(lowDimGridId, bulkGridId); | |
| 232 | 6578 | auto it = map.find(eIdx); | |
| 233 |
2/4✓ Branch 0 taken 6118 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 6118 times.
✗ Branch 3 not taken.
|
13156 | if (it != map.end()) return it->second.couplingStencil; |
| 234 | ✗ | else return getEmptyStencil(bulkId); | |
| 235 | } | ||
| 236 | |||
| 237 | /*! | ||
| 238 | * \brief returns true if a bulk scvf flux depends on data in the facet domain. | ||
| 239 | * \note for box this is the case if an scvf lies on an interior boundary | ||
| 240 | */ | ||
| 241 | bool isCoupled(const Element<bulkId>& element, | ||
| 242 | const SubControlVolumeFace<bulkId>& scvf) const | ||
| 243 | { return scvf.interiorBoundary(); } | ||
| 244 | |||
| 245 | /*! | ||
| 246 | * \brief returns true if a bulk scvf coincides with a facet element. | ||
| 247 | * \note for box, this is always true for coupled scvfs | ||
| 248 | */ | ||
| 249 | bool isOnInteriorBoundary(const Element<bulkId>& element, | ||
| 250 | const SubControlVolumeFace<bulkId>& scvf) const | ||
| 251 | { return isCoupled(element, scvf); } | ||
| 252 | |||
| 253 | /*! | ||
| 254 | * \brief returns the vol vars of a lower-dimensional element coinciding with a bulk scvf. | ||
| 255 | */ | ||
| 256 | 4823120 | const VolumeVariables<lowDimId>& getLowDimVolVars(const Element<bulkId>& element, | |
| 257 | const SubControlVolumeFace<bulkId>& scvf) const | ||
| 258 | { | ||
| 259 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 4217968 times.
|
4823120 | const auto eIdx = this->problem(bulkId).gridGeometry().elementMapper().index(element); |
| 260 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 4217968 times.
|
4823120 | assert(bulkContext_.isSet); |
| 261 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 4217968 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4217968 times.
|
9646240 | assert(bulkElemIsCoupled_[eIdx]); |
| 262 | |||
| 263 | 9646240 | const auto& map = couplingMapperPtr_->couplingMap(bulkGridId, lowDimGridId); | |
| 264 | 9646240 | const auto& couplingData = map.find(eIdx)->second; | |
| 265 | |||
| 266 | // search the low dim element idx this scvf is embedded in | ||
| 267 | // and find out the local index of the scv it couples to | ||
| 268 | unsigned int coupledScvIdx; | ||
| 269 | 15093148 | auto it = std::find_if( couplingData.elementToScvfMap.begin(), | |
| 270 | couplingData.elementToScvfMap.end(), | ||
| 271 | 7849420 | [&] (auto& dataPair) | |
| 272 | { | ||
| 273 | 4236604 | const auto& scvfList = dataPair.second; | |
| 274 | 12709812 | auto it = std::find(scvfList.begin(), scvfList.end(), scvf.index()); | |
| 275 | 16946416 | coupledScvIdx = std::distance(scvfList.begin(), it); | |
| 276 | 12709812 | return it != scvfList.end(); | |
| 277 | } ); | ||
| 278 | |||
| 279 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 4217968 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4217968 times.
|
9646240 | assert(it != couplingData.elementToScvfMap.end()); |
| 280 | 4823120 | const auto lowDimElemIdx = it->first; | |
| 281 | 9646240 | const auto& s = map.find(bulkContext_.elementIdx)->second.couplingElementStencil; | |
| 282 | 14469360 | const auto& idxInContext = std::distance( s.begin(), std::find(s.begin(), s.end(), lowDimElemIdx) ); | |
| 283 |
2/4✗ Branch 4 not taken.
✓ Branch 5 taken 4217968 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 4217968 times.
|
19292480 | assert(std::find(s.begin(), s.end(), lowDimElemIdx) != s.end()); |
| 284 | |||
| 285 | if (lowDimUsesBox) | ||
| 286 | 14469360 | return bulkContext_.lowDimElemVolVars[idxInContext][coupledScvIdx]; | |
| 287 | else | ||
| 288 | return bulkContext_.lowDimElemVolVars[idxInContext][0]; | ||
| 289 | } | ||
| 290 | |||
| 291 | /*! | ||
| 292 | * \brief returns the lower-dimensional element coinciding with a bulk scvf. | ||
| 293 | */ | ||
| 294 | const Element<lowDimId> getLowDimElement(const Element<bulkId>& element, | ||
| 295 | const SubControlVolumeFace<bulkId>& scvf) const | ||
| 296 | { | ||
| 297 | const auto lowDimElemIdx = getLowDimElementIndex(element, scvf); | ||
| 298 | return this->problem(lowDimId).gridGeometry().element(lowDimElemIdx); | ||
| 299 | } | ||
| 300 | |||
| 301 | /*! | ||
| 302 | * \brief returns the index of the lower-dimensional element coinciding with a bulk scvf. | ||
| 303 | */ | ||
| 304 | const GridIndexType<lowDimId> getLowDimElementIndex(const Element<bulkId>& element, | ||
| 305 | const SubControlVolumeFace<bulkId>& scvf) const | ||
| 306 | { | ||
| 307 | const auto eIdx = this->problem(bulkId).gridGeometry().elementMapper().index(element); | ||
| 308 | |||
| 309 | assert(bulkElemIsCoupled_[eIdx]); | ||
| 310 | const auto& map = couplingMapperPtr_->couplingMap(bulkGridId, lowDimGridId); | ||
| 311 | const auto& couplingData = map.at(eIdx); | ||
| 312 | |||
| 313 | // search the low dim element idx this scvf is embedded in | ||
| 314 | auto it = std::find_if( couplingData.elementToScvfMap.begin(), | ||
| 315 | couplingData.elementToScvfMap.end(), | ||
| 316 | [&] (auto& dataPair) | ||
| 317 | { | ||
| 318 | const auto& scvfList = dataPair.second; | ||
| 319 | auto it = std::find(scvfList.begin(), scvfList.end(), scvf.index()); | ||
| 320 | return it != scvfList.end(); | ||
| 321 | } ); | ||
| 322 | |||
| 323 | assert(it != couplingData.elementToScvfMap.end()); | ||
| 324 | return it->first; | ||
| 325 | } | ||
| 326 | |||
| 327 | /*! | ||
| 328 | * \brief Evaluates the coupling element residual of a bulk domain element with respect | ||
| 329 | * to a dof in the lower-dimensional domain (dofIdxGlobalJ). This is essentially | ||
| 330 | * the fluxes across the bulk element facets that coincide with the lower-dimensional | ||
| 331 | * element whose dof idx is dofIdxGlobalJ. | ||
| 332 | */ | ||
| 333 | template< class BulkLocalAssembler > | ||
| 334 | typename LocalResidual<bulkId>::ElementResidualVector | ||
| 335 | 66160 | evalCouplingResidual(BulkIdType, | |
| 336 | const BulkLocalAssembler& bulkLocalAssembler, | ||
| 337 | LowDimIdType, | ||
| 338 | GridIndexType<lowDimId> dofIdxGlobalJ) | ||
| 339 | { | ||
| 340 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 61680 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 61680 times.
|
132320 | const auto& map = couplingMapperPtr_->couplingMap(bulkGridId, lowDimGridId); |
| 341 | |||
| 342 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 61680 times.
|
66160 | assert(bulkContext_.isSet); |
| 343 |
3/6✗ Branch 0 not taken.
✓ Branch 1 taken 61680 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 61680 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 61680 times.
|
198480 | assert(bulkElemIsCoupled_[bulkContext_.elementIdx]); |
| 344 |
1/2✗ Branch 2 not taken.
✓ Branch 3 taken 61680 times.
|
198480 | assert(map.find(bulkContext_.elementIdx) != map.end()); |
| 345 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 61680 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 61680 times.
|
66160 | assert(bulkContext_.elementIdx == this->problem(bulkId).gridGeometry().elementMapper().index(bulkLocalAssembler.element())); |
| 346 | |||
| 347 | 66160 | const auto& fvGeometry = bulkLocalAssembler.fvGeometry(); | |
| 348 | 132320 | typename LocalResidual<bulkId>::ElementResidualVector res(fvGeometry.numScv()); | |
| 349 | 66160 | res = 0.0; | |
| 350 | |||
| 351 | // compute fluxes across the coupling scvfs | ||
| 352 | 198480 | const auto& couplingScvfs = map.find(bulkContext_.elementIdx)->second.dofToCouplingScvfMap.at(dofIdxGlobalJ); | |
| 353 |
4/4✓ Branch 0 taken 127632 times.
✓ Branch 1 taken 61680 times.
✓ Branch 2 taken 127632 times.
✓ Branch 3 taken 61680 times.
|
480112 | for (auto scvfIdx : couplingScvfs) |
| 354 | { | ||
| 355 | 140816 | const auto& scvf = fvGeometry.scvf(scvfIdx); | |
| 356 | 281632 | const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx()); | |
| 357 | 933040 | res[insideScv.localDofIndex()] += evalBulkFluxes_( bulkLocalAssembler.element(), | |
| 358 | bulkLocalAssembler.fvGeometry(), | ||
| 359 | bulkLocalAssembler.curElemVolVars(), | ||
| 360 | bulkLocalAssembler.elemBcTypes(), | ||
| 361 | bulkLocalAssembler.elemFluxVarsCache(), | ||
| 362 | bulkLocalAssembler.localResidual(), | ||
| 363 | std::array<GridIndexType<bulkId>, 1>({scvfIdx}) ); | ||
| 364 | } | ||
| 365 | |||
| 366 | 66160 | return res; | |
| 367 | } | ||
| 368 | |||
| 369 | /*! | ||
| 370 | * \brief Evaluates the coupling element residual of a lower-dimensional domain element | ||
| 371 | * with respect to a dof in the bulk domain (dofIdxGlobalJ). This is essentially | ||
| 372 | * the fluxes across the facets of the neighboring bulk element that coincide with | ||
| 373 | * the given element. | ||
| 374 | */ | ||
| 375 | template< class LowDimLocalAssembler > | ||
| 376 | typename LocalResidual<lowDimId>::ElementResidualVector | ||
| 377 | 121218 | evalCouplingResidual(LowDimIdType, | |
| 378 | const LowDimLocalAssembler& lowDimLocalAssembler, | ||
| 379 | BulkIdType, | ||
| 380 | GridIndexType<bulkId> dofIdxGlobalJ) | ||
| 381 | { | ||
| 382 | // make sure this is called for the element for which the context was set | ||
| 383 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 115938 times.
|
121218 | assert(lowDimContext_.isSet); |
| 384 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 115938 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 115938 times.
|
121218 | assert(this->problem(lowDimId).gridGeometry().elementMapper().index(lowDimLocalAssembler.element()) == lowDimContext_.elementIdx); |
| 385 | |||
| 386 | // evaluate sources for all scvs in lower-dimensional element | ||
| 387 | 363654 | typename LocalResidual<lowDimId>::ElementResidualVector res(lowDimLocalAssembler.fvGeometry().numScv()); | |
| 388 | 121218 | res = 0.0; | |
| 389 |
4/4✓ Branch 0 taken 236772 times.
✓ Branch 1 taken 115938 times.
✓ Branch 2 taken 236772 times.
✓ Branch 3 taken 115938 times.
|
746892 | for (const auto& scv : scvs(lowDimLocalAssembler.fvGeometry())) |
| 390 | 929824 | res[scv.localDofIndex()] -= evalSourcesFromBulk(lowDimLocalAssembler.element(), | |
| 391 | lowDimLocalAssembler.fvGeometry(), | ||
| 392 | lowDimLocalAssembler.curElemVolVars(), | ||
| 393 | scv); | ||
| 394 | 121218 | return res; | |
| 395 | } | ||
| 396 | |||
| 397 | /*! | ||
| 398 | * \brief Computes the sources in a lower-dimensional element stemming from the bulk domain. | ||
| 399 | */ | ||
| 400 | 520856 | NumEqVector<lowDimId> evalSourcesFromBulk(const Element<lowDimId>& element, | |
| 401 | const FVElementGeometry<lowDimId>& fvGeometry, | ||
| 402 | const ElementVolumeVariables<lowDimId>& elemVolVars, | ||
| 403 | const SubControlVolume<lowDimId>& scv) const | ||
| 404 | { | ||
| 405 | // make sure the this is called for the element of the context | ||
| 406 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 389160 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 389160 times.
|
520856 | assert(this->problem(lowDimId).gridGeometry().elementMapper().index(element) == lowDimContext_.elementIdx); |
| 407 | |||
| 408 | 520856 | NumEqVector<lowDimId> sources(0.0); | |
| 409 | 1041712 | const auto& map = couplingMapperPtr_->couplingMap(lowDimGridId, bulkGridId); | |
| 410 | 520856 | auto it = map.find(lowDimContext_.elementIdx); | |
| 411 |
4/4✓ Branch 0 taken 120512 times.
✓ Branch 1 taken 268648 times.
✓ Branch 2 taken 120512 times.
✓ Branch 3 taken 268648 times.
|
1041712 | if (it == map.end()) |
| 412 | ✗ | return sources; | |
| 413 | |||
| 414 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 389160 times.
|
520856 | assert(lowDimContext_.isSet); |
| 415 |
6/6✓ Branch 0 taken 777456 times.
✓ Branch 1 taken 389160 times.
✓ Branch 2 taken 777456 times.
✓ Branch 3 taken 389160 times.
✓ Branch 4 taken 777456 times.
✓ Branch 5 taken 389160 times.
|
4690296 | for (unsigned int i = 0; i < it->second.embedments.size(); ++i) |
| 416 | { | ||
| 417 | 2085152 | const auto& embedment = it->second.embedments[i]; | |
| 418 | |||
| 419 | // list of scvfs in the bulk domain whose fluxes enter this scv | ||
| 420 | // if low dim domain uses tpfa, this is all scvfs lying on this element | ||
| 421 | // if it uses box, it is the one scvf coinciding with the given scv | ||
| 422 | 1042576 | const auto& coincidingScvfs = embedment.second; | |
| 423 |
1/4✗ Branch 3 not taken.
✓ Branch 4 taken 777456 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
|
4170304 | const auto& scvfList = lowDimUsesBox ? std::vector<GridIndexType<lowDimId>>{ coincidingScvfs[scv.localDofIndex()] } |
| 424 | : coincidingScvfs; | ||
| 425 | |||
| 426 |
5/12✓ Branch 1 taken 777456 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 777456 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 777456 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 3456 times.
✓ Branch 10 taken 774000 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
|
2513888 | sources += evalBulkFluxes_(this->problem(bulkId).gridGeometry().element(embedment.first), |
| 427 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 777456 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 777456 times.
|
2085152 | *(lowDimContext_.bulkFvGeometries[i]), |
| 428 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 777456 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 777456 times.
|
2085152 | *(lowDimContext_.bulkElemVolVars[i]), |
| 429 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 777456 times.
|
1042576 | lowDimContext_.bulkElemBcTypes[i], |
| 430 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 777456 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 777456 times.
|
2085152 | *(lowDimContext_.bulkElemFluxVarsCache[i]), |
| 431 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 777456 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 777456 times.
|
2085152 | *lowDimContext_.bulkLocalResidual, |
| 432 | scvfList); | ||
| 433 | } | ||
| 434 | |||
| 435 | 400344 | return sources; | |
| 436 | } | ||
| 437 | |||
| 438 | /*! | ||
| 439 | * \brief For the assembly of the element residual of a bulk domain element | ||
| 440 | * we need to prepare all variables of lower-dimensional domain elements | ||
| 441 | * that are coupled to the given bulk element | ||
| 442 | */ | ||
| 443 | template< class Assembler > | ||
| 444 | 2594260 | void bindCouplingContext(BulkIdType, const Element<bulkId>& element, const Assembler& assembler) | |
| 445 | { | ||
| 446 | // clear context | ||
| 447 | 2594260 | bulkContext_.reset(); | |
| 448 | |||
| 449 | // set index in context in any case | ||
| 450 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1297130 times.
|
2594260 | const auto bulkElemIdx = this->problem(bulkId).gridGeometry().elementMapper().index(element); |
| 451 | 2594260 | bulkContext_.elementIdx = bulkElemIdx; | |
| 452 | |||
| 453 | // if element is coupled, actually set the context | ||
| 454 |
6/6✓ Branch 0 taken 164311 times.
✓ Branch 1 taken 1132819 times.
✓ Branch 2 taken 164311 times.
✓ Branch 3 taken 1132819 times.
✓ Branch 4 taken 164311 times.
✓ Branch 5 taken 1132819 times.
|
7782780 | if (bulkElemIsCoupled_[bulkElemIdx]) |
| 455 | { | ||
| 456 | 657244 | const auto& map = couplingMapperPtr_->couplingMap(bulkGridId, lowDimGridId); | |
| 457 | |||
| 458 |
2/4✗ Branch 1 not taken.
✓ Branch 2 taken 164311 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 164311 times.
|
985866 | auto it = map.find(bulkElemIdx); assert(it != map.end()); |
| 459 | 328622 | const auto& elementStencil = it->second.couplingElementStencil; | |
| 460 | 657244 | bulkContext_.lowDimFvGeometries.reserve(elementStencil.size()); | |
| 461 | 657244 | bulkContext_.lowDimElemVolVars.reserve(elementStencil.size()); | |
| 462 | |||
| 463 | // local view on the bulk fv geometry | ||
| 464 |
4/8✗ Branch 0 not taken.
✓ Branch 1 taken 164311 times.
✓ Branch 3 taken 164311 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 164311 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 164311 times.
✗ Branch 10 not taken.
|
985866 | const auto bulkFvGeometry = localView(this->problem(bulkId).gridGeometry()).bindElement(element); |
| 465 | |||
| 466 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 164311 times.
|
328622 | const auto& ldGridGeometry = this->problem(lowDimId).gridGeometry(); |
| 467 | 657244 | auto fvGeom = localView(ldGridGeometry); | |
| 468 |
4/4✓ Branch 0 taken 167974 times.
✓ Branch 1 taken 164311 times.
✓ Branch 2 taken 167974 times.
✓ Branch 3 taken 164311 times.
|
1657762 | for (const auto lowDimElemIdx : elementStencil) |
| 469 | { | ||
| 470 |
2/4✓ Branch 1 taken 167974 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 152000 times.
✗ Branch 5 not taken.
|
639948 | const auto& ldSol = assembler.isImplicit() ? this->curSol(lowDimId) : assembler.prevSol()[lowDimId]; |
| 471 |
1/2✓ Branch 1 taken 167974 times.
✗ Branch 2 not taken.
|
335948 | const auto elemJ = ldGridGeometry.element(lowDimElemIdx); |
| 472 |
1/2✓ Branch 1 taken 167974 times.
✗ Branch 2 not taken.
|
335948 | fvGeom.bindElement(elemJ); |
| 473 | |||
| 474 |
4/10✓ Branch 1 taken 167974 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 167974 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 167974 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 167974 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
|
1343792 | std::vector<VolumeVariables<lowDimId>> elemVolVars(fvGeom.numScv()); |
| 475 |
2/4✓ Branch 1 taken 167974 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 167974 times.
✗ Branch 5 not taken.
|
671896 | const auto& coupledScvfIndices = it->second.elementToScvfMap.at(lowDimElemIdx); |
| 476 |
1/2✓ Branch 1 taken 167974 times.
✗ Branch 2 not taken.
|
335948 | makeCoupledLowDimElemVolVars_(element, bulkFvGeometry, elemJ, fvGeom, ldSol, coupledScvfIndices, elemVolVars); |
| 477 | |||
| 478 | 335948 | bulkContext_.isSet = true; | |
| 479 |
1/2✓ Branch 1 taken 167974 times.
✗ Branch 2 not taken.
|
335948 | bulkContext_.lowDimFvGeometries.emplace_back( std::move(fvGeom) ); |
| 480 |
1/2✓ Branch 1 taken 167974 times.
✗ Branch 2 not taken.
|
335948 | bulkContext_.lowDimElemVolVars.emplace_back( std::move(elemVolVars) ); |
| 481 | } | ||
| 482 | } | ||
| 483 | 2594260 | } | |
| 484 | |||
| 485 | /*! | ||
| 486 | * \brief For the assembly of the element residual of a bulk domain element | ||
| 487 | * we need to prepare the local views of one of the neighboring bulk | ||
| 488 | * domain elements. These are used later to compute the fluxes across | ||
| 489 | * the faces over which the coupling occurs | ||
| 490 | * \note Since the low-dim coupling residua are fluxes stemming from | ||
| 491 | * the bulk domain, we have to prepare the bulk coupling context | ||
| 492 | * for the neighboring element (where fluxes are calculated) as well. | ||
| 493 | */ | ||
| 494 | template< class Assembler > | ||
| 495 | 106202 | void bindCouplingContext(LowDimIdType, const Element<lowDimId>& element, const Assembler& assembler) | |
| 496 | { | ||
| 497 | // reset contexts | ||
| 498 | 106202 | bulkContext_.reset(); | |
| 499 | 106202 | lowDimContext_.reset(); | |
| 500 | |||
| 501 | // set index in context in any case | ||
| 502 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 55334 times.
|
106202 | const auto lowDimElemIdx = this->problem(lowDimId).gridGeometry().elementMapper().index(element); |
| 503 | 106202 | lowDimContext_.elementIdx = lowDimElemIdx; | |
| 504 | |||
| 505 | 212404 | const auto& map = couplingMapperPtr_->couplingMap(lowDimGridId, bulkGridId); | |
| 506 | 106202 | auto it = map.find(lowDimElemIdx); | |
| 507 | |||
| 508 | // if element is coupled, actually set the context | ||
| 509 |
2/4✓ Branch 0 taken 55334 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 55334 times.
✗ Branch 3 not taken.
|
212404 | if (it != map.end()) |
| 510 | { | ||
| 511 | // first bind the low dim context for the first neighboring bulk element | ||
| 512 | // this will set up the lower-dimensional data on this current lowdim element | ||
| 513 | // which will be enough to compute the fluxes in all neighboring elements | ||
| 514 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 55334 times.
|
106202 | const auto& bulkGridGeom = this->problem(bulkId).gridGeometry(); |
| 515 | 318574 | const auto bulkElem = bulkGridGeom.element(it->second.embedments[0].first); | |
| 516 |
1/2✓ Branch 1 taken 55318 times.
✗ Branch 2 not taken.
|
106202 | bindCouplingContext(bulkId, bulkElem, assembler); |
| 517 | |||
| 518 | // reserve space in the context and bind local views of neighbors | ||
| 519 |
1/2✓ Branch 1 taken 1232 times.
✗ Branch 2 not taken.
|
106202 | const auto& embedments = it->second.embedments; |
| 520 |
1/2✓ Branch 1 taken 1232 times.
✗ Branch 2 not taken.
|
106202 | const auto numEmbedments = embedments.size(); |
| 521 |
1/2✓ Branch 1 taken 1232 times.
✗ Branch 2 not taken.
|
106202 | lowDimContext_.resize(numEmbedments); |
| 522 | |||
| 523 | 212372 | auto bulkFvGeom = localView(bulkGridGeom); | |
| 524 |
0/2✗ Branch 1 not taken.
✗ Branch 2 not taken.
|
424808 | auto bulkElemVolVars = localView(assembler.gridVariables(bulkId).curGridVolVars()); |
| 525 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 55334 times.
|
424808 | auto bulkElemFluxVarsCache = localView(assembler.gridVariables(bulkId).gridFluxVarsCache()); |
| 526 | |||
| 527 |
2/2✓ Branch 0 taken 110556 times.
✓ Branch 1 taken 55334 times.
|
318526 | for (unsigned int i = 0; i < numEmbedments; ++i) |
| 528 | { | ||
| 529 |
2/4✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 100000 times.
✗ Branch 5 not taken.
|
412324 | const auto& bulkSol = assembler.isImplicit() ? this->curSol(bulkId) : assembler.prevSol()[bulkId]; |
| 530 |
2/4✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 110556 times.
✗ Branch 5 not taken.
|
636844 | const auto curBulkElem = bulkGridGeom.element(embedments[i].first); |
| 531 | |||
| 532 |
1/2✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
|
212324 | bulkFvGeom.bind(curBulkElem); |
| 533 |
1/2✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
|
212324 | bulkElemVolVars.bind(curBulkElem, bulkFvGeom, bulkSol); |
| 534 |
1/2✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
|
212324 | bulkElemFluxVarsCache.bind(curBulkElem, bulkFvGeom, bulkElemVolVars); |
| 535 | |||
| 536 | 212324 | lowDimContext_.isSet = true; | |
| 537 |
3/6✗ Branch 0 not taken.
✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 110556 times.
✓ Branch 5 taken 110556 times.
✗ Branch 6 not taken.
|
424648 | lowDimContext_.bulkElemBcTypes[i].update(this->problem(bulkId), curBulkElem, bulkFvGeom); |
| 538 |
4/8✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 110556 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 110556 times.
✓ Branch 8 taken 110556 times.
✗ Branch 9 not taken.
|
636972 | lowDimContext_.bulkFvGeometries[i] = std::make_unique< FVElementGeometry<bulkId> >( std::move(bulkFvGeom) ); |
| 539 |
4/8✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 110556 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 110556 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 110556 times.
|
212324 | lowDimContext_.bulkElemVolVars[i] = std::make_unique< ElementVolumeVariables<bulkId> >( std::move(bulkElemVolVars) ); |
| 540 |
4/8✓ Branch 1 taken 110556 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 110556 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 110556 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 110556 times.
|
212324 | lowDimContext_.bulkElemFluxVarsCache[i] = std::make_unique< ElementFluxVariablesCache<bulkId> >( std::move(bulkElemFluxVarsCache) ); |
| 541 | } | ||
| 542 | |||
| 543 | // finally, set the local residual | ||
| 544 |
3/6✓ Branch 1 taken 55334 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 55334 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 55334 times.
|
212404 | lowDimContext_.bulkLocalResidual = std::make_unique< LocalResidual<bulkId> >(assembler.localResidual(bulkId)); |
| 545 | } | ||
| 546 | 106202 | } | |
| 547 | |||
| 548 | /*! | ||
| 549 | * \brief After deflecting the solution of the lower-dimensional domain, | ||
| 550 | * we have to update the element volume variables object if the context. | ||
| 551 | */ | ||
| 552 | template< class BulkLocalAssembler > | ||
| 553 | 85032 | void updateCouplingContext(BulkIdType domainI, | |
| 554 | const BulkLocalAssembler& bulkLocalAssembler, | ||
| 555 | LowDimIdType domainJ, | ||
| 556 | GridIndexType<lowDimId> dofIdxGlobalJ, | ||
| 557 | const PrimaryVariables<lowDimId>& priVarsJ, | ||
| 558 | unsigned int pvIdxJ) | ||
| 559 | { | ||
| 560 | // communicate deflected solution | ||
| 561 |
1/2✓ Branch 0 taken 80552 times.
✗ Branch 1 not taken.
|
85032 | ParentType::updateCouplingContext(domainI, bulkLocalAssembler, domainJ, dofIdxGlobalJ, priVarsJ, pvIdxJ); |
| 562 | |||
| 563 | // Since coupling only occurs via the fluxes, the context does not | ||
| 564 | // have to be updated in explicit time discretization schemes, where | ||
| 565 | // they are strictly evaluated on the old time level | ||
| 566 | if (!bulkLocalAssembler.isImplicit()) | ||
| 567 | return; | ||
| 568 | |||
| 569 | // skip the rest if context is empty | ||
| 570 |
1/2✓ Branch 0 taken 80552 times.
✗ Branch 1 not taken.
|
85032 | if (bulkContext_.isSet) |
| 571 | { | ||
| 572 | 170064 | const auto& map = couplingMapperPtr_->couplingMap(bulkGridId, lowDimGridId); | |
| 573 | 85032 | const auto& couplingEntry = map.at(bulkContext_.elementIdx); | |
| 574 | 85032 | const auto& couplingElemStencil = couplingEntry.couplingElementStencil; | |
| 575 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 80552 times.
|
85032 | const auto& ldGridGeometry = this->problem(lowDimId).gridGeometry(); |
| 576 | |||
| 577 | // find the low-dim elements in coupling stencil, where this dof is contained in | ||
| 578 |
1/4✓ Branch 1 taken 80552 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
|
259576 | const auto couplingElements = [&] () |
| 579 | { | ||
| 580 | if (lowDimUsesBox) | ||
| 581 | { | ||
| 582 |
0/2✗ Branch 0 not taken.
✗ Branch 1 not taken.
|
161104 | std::vector< Element<lowDimId> > lowDimElems; |
| 583 | 246904 | std::for_each( couplingElemStencil.begin(), couplingElemStencil.end(), | |
| 584 | 81392 | [&] (auto lowDimElemIdx) | |
| 585 | { | ||
| 586 | 81392 | auto element = ldGridGeometry.element(lowDimElemIdx); | |
| 587 |
8/8✓ Branch 1 taken 114624 times.
✓ Branch 2 taken 48 times.
✓ Branch 3 taken 114624 times.
✓ Branch 4 taken 48 times.
✓ Branch 6 taken 10752 times.
✓ Branch 7 taken 768 times.
✓ Branch 8 taken 10752 times.
✓ Branch 9 taken 768 times.
|
126192 | for (unsigned int i = 0; i < element.geometry().corners(); ++i) |
| 588 | { | ||
| 589 | 250752 | const auto dofIdx = ldGridGeometry.vertexMapper().subIndex(element, i, lowDimDim); | |
| 590 |
4/4✓ Branch 0 taken 76352 times.
✓ Branch 1 taken 38272 times.
✓ Branch 3 taken 4224 times.
✓ Branch 4 taken 6528 times.
|
125376 | if (dofIdxGlobalJ == dofIdx) { lowDimElems.emplace_back( std::move(element) ); break; } |
| 591 | } | ||
| 592 | } ); | ||
| 593 | 161104 | return lowDimElems; | |
| 594 | } | ||
| 595 | // dof index = element index for cc schemes | ||
| 596 | else | ||
| 597 | return std::vector<Element<lowDimId>>( {ldGridGeometry.element(dofIdxGlobalJ)} ); | ||
| 598 | } (); | ||
| 599 | |||
| 600 | // update all necessary vol vars in context | ||
| 601 |
4/4✓ Branch 0 taken 80576 times.
✓ Branch 1 taken 80552 times.
✓ Branch 2 taken 80576 times.
✓ Branch 3 taken 80552 times.
|
340152 | for (const auto& element : couplingElements) |
| 602 | { | ||
| 603 | // find index in coupling context | ||
| 604 | 170112 | const auto eIdxGlobal = ldGridGeometry.elementMapper().index(element); | |
| 605 | 255168 | auto it = std::find(couplingElemStencil.begin(), couplingElemStencil.end(), eIdxGlobal); | |
| 606 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 80576 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 80576 times.
|
170112 | const auto idxInContext = std::distance(couplingElemStencil.begin(), it); |
| 607 |
3/6✗ Branch 0 not taken.
✓ Branch 1 taken 80576 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 80576 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 80576 times.
|
255168 | assert(it != couplingElemStencil.end()); |
| 608 | |||
| 609 |
1/2✓ Branch 1 taken 80576 times.
✗ Branch 2 not taken.
|
85056 | auto& elemVolVars = bulkContext_.lowDimElemVolVars[idxInContext]; |
| 610 |
1/2✓ Branch 1 taken 80576 times.
✗ Branch 2 not taken.
|
85056 | const auto& fvGeom = bulkContext_.lowDimFvGeometries[idxInContext]; |
| 611 |
1/2✓ Branch 1 taken 80576 times.
✗ Branch 2 not taken.
|
85056 | const auto& coupledScvfIndices = couplingEntry.elementToScvfMap.at(eIdxGlobal); |
| 612 |
3/6✓ Branch 1 taken 80576 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 80576 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 80576 times.
✗ Branch 8 not taken.
|
255168 | makeCoupledLowDimElemVolVars_(bulkLocalAssembler.element(), bulkLocalAssembler.fvGeometry(), |
| 613 | element, fvGeom, this->curSol(lowDimId), coupledScvfIndices, elemVolVars); | ||
| 614 | } | ||
| 615 | } | ||
| 616 | } | ||
| 617 | |||
| 618 | /*! | ||
| 619 | * \brief Update the coupling context for a derivative bulk -> bulk. | ||
| 620 | * Here, we simply have to update the solution. | ||
| 621 | */ | ||
| 622 | template< class BulkLocalAssembler > | ||
| 623 | ✗ | void updateCouplingContext(BulkIdType domainI, | |
| 624 | const BulkLocalAssembler& bulkLocalAssembler, | ||
| 625 | BulkIdType domainJ, | ||
| 626 | GridIndexType<bulkId> dofIdxGlobalJ, | ||
| 627 | const PrimaryVariables<bulkId>& priVarsJ, | ||
| 628 | unsigned int pvIdxJ) | ||
| 629 | { | ||
| 630 | // communicate deflected solution | ||
| 631 | 8225600 | ParentType::updateCouplingContext(domainI, bulkLocalAssembler, domainJ, dofIdxGlobalJ, priVarsJ, pvIdxJ); | |
| 632 | ✗ | } | |
| 633 | |||
| 634 | /*! | ||
| 635 | * \brief After deflecting the solution of the bulk domain, we have to update | ||
| 636 | * the element volume variables of the neighboring bulk elements stored | ||
| 637 | * in the context. | ||
| 638 | */ | ||
| 639 | template< class LowDimLocalAssembler > | ||
| 640 | 158288 | void updateCouplingContext(LowDimIdType domainI, | |
| 641 | const LowDimLocalAssembler& lowDimLocalAssembler, | ||
| 642 | BulkIdType domainJ, | ||
| 643 | GridIndexType<bulkId> dofIdxGlobalJ, | ||
| 644 | const PrimaryVariables<bulkId>& priVarsJ, | ||
| 645 | unsigned int pvIdxJ) | ||
| 646 | { | ||
| 647 | // communicate deflected solution | ||
| 648 |
1/2✓ Branch 0 taken 153008 times.
✗ Branch 1 not taken.
|
158288 | ParentType::updateCouplingContext(domainI, lowDimLocalAssembler, domainJ, dofIdxGlobalJ, priVarsJ, pvIdxJ); |
| 649 | |||
| 650 | // Since coupling only occurs via the fluxes, the context does not | ||
| 651 | // have to be updated in explicit time discretization schemes, where | ||
| 652 | // they are strictly evaluated on the old time level | ||
| 653 | if (!lowDimLocalAssembler.isImplicit()) | ||
| 654 | return; | ||
| 655 | |||
| 656 | // skip the rest if context is empty | ||
| 657 |
1/2✓ Branch 0 taken 153008 times.
✗ Branch 1 not taken.
|
158288 | if (lowDimContext_.isSet) |
| 658 | { | ||
| 659 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 153008 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 153008 times.
|
158288 | assert(lowDimContext_.elementIdx == this->problem(lowDimId).gridGeometry().elementMapper().index(lowDimLocalAssembler.element())); |
| 660 | |||
| 661 | 316576 | const auto& map = couplingMapperPtr_->couplingMap(lowDimGridId, bulkGridId); | |
| 662 | 158288 | auto it = map.find(lowDimContext_.elementIdx); | |
| 663 | |||
| 664 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 153008 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 153008 times.
|
316576 | assert(it != map.end()); |
| 665 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 153008 times.
|
158288 | const auto& embedments = it->second.embedments; |
| 666 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 153008 times.
|
158288 | const auto& bulkGridGeometry = this->problem(bulkId).gridGeometry(); |
| 667 | |||
| 668 | // TODO more efficient (only one dof update per embedment) | ||
| 669 | // update the elem volvars in context of those elements that contain globalJ | ||
| 670 | 158288 | unsigned int embedmentIdx = 0; | |
| 671 |
4/4✓ Branch 0 taken 305920 times.
✓ Branch 1 taken 153008 times.
✓ Branch 2 taken 305920 times.
✓ Branch 3 taken 153008 times.
|
1109360 | for (const auto& embedment : embedments) |
| 672 | { | ||
| 673 |
1/2✓ Branch 1 taken 304384 times.
✗ Branch 2 not taken.
|
631424 | const auto elementJ = bulkGridGeometry.element(embedment.first); |
| 674 | |||
| 675 |
5/5✓ Branch 0 taken 4608 times.
✓ Branch 1 taken 1511344 times.
✓ Branch 2 taken 4608 times.
✓ Branch 3 taken 1206960 times.
✓ Branch 4 taken 304384 times.
|
1571056 | for (unsigned int i = 0; i < elementJ.subEntities(bulkDim); ++i) |
| 676 | { | ||
| 677 |
2/4✓ Branch 1 taken 1205424 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1205424 times.
✗ Branch 5 not taken.
|
2507616 | const auto dofIdx = bulkGridGeometry.vertexMapper().subIndex(elementJ, i, bulkDim); |
| 678 |
2/2✓ Branch 0 taken 156032 times.
✓ Branch 1 taken 1054000 times.
|
1253808 | if (dofIdx == dofIdxGlobalJ) |
| 679 | { | ||
| 680 | // element contains the deflected dof | ||
| 681 |
2/4✓ Branch 1 taken 155648 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 155648 times.
✗ Branch 5 not taken.
|
324096 | const auto& fvGeom = *lowDimContext_.bulkFvGeometries[embedmentIdx]; |
| 682 |
4/8✓ Branch 1 taken 155648 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 155648 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 155648 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 155648 times.
✗ Branch 11 not taken.
|
648192 | (*lowDimContext_.bulkElemVolVars[embedmentIdx]).bindElement(elementJ, fvGeom, this->curSol(bulkId)); |
| 683 | } | ||
| 684 | } | ||
| 685 | |||
| 686 | // keep track of embedments | ||
| 687 | 317248 | embedmentIdx++; | |
| 688 | } | ||
| 689 | } | ||
| 690 | } | ||
| 691 | |||
| 692 | /*! | ||
| 693 | * \brief After deflecting the solution of the lower-dimensional domain has been deflected | ||
| 694 | * during the assembly of the element residual of a lower-dimensional element, we | ||
| 695 | * have to communicate this to the volume variables stored in the context as well | ||
| 696 | * as the transmissibilities. | ||
| 697 | */ | ||
| 698 | template< class LowDimLocalAssembler > | ||
| 699 | 42688 | void updateCouplingContext(LowDimIdType domainI, | |
| 700 | const LowDimLocalAssembler& lowDimLocalAssembler, | ||
| 701 | LowDimIdType domainJ, | ||
| 702 | GridIndexType<lowDimId> dofIdxGlobalJ, | ||
| 703 | const PrimaryVariables<lowDimId>& priVarsJ, | ||
| 704 | unsigned int pvIdxJ) | ||
| 705 | { | ||
| 706 | // communicate deflected solution | ||
| 707 |
1/2✓ Branch 0 taken 40512 times.
✗ Branch 1 not taken.
|
42688 | ParentType::updateCouplingContext(domainI, lowDimLocalAssembler, domainJ, dofIdxGlobalJ, priVarsJ, pvIdxJ); |
| 708 | |||
| 709 | // Since coupling only occurs via the fluxes, the context does not | ||
| 710 | // have to be updated in explicit time discretization schemes, where | ||
| 711 | // they are strictly evaluated on the old time level | ||
| 712 | if (!lowDimLocalAssembler.isImplicit()) | ||
| 713 | return; | ||
| 714 | |||
| 715 | // skip the rest if context is empty | ||
| 716 |
1/2✓ Branch 0 taken 40512 times.
✗ Branch 1 not taken.
|
42688 | if (lowDimContext_.isSet) |
| 717 | { | ||
| 718 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 40512 times.
|
42688 | assert(bulkContext_.isSet); |
| 719 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 40512 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 40512 times.
|
42688 | assert(lowDimContext_.elementIdx == this->problem(lowDimId).gridGeometry().elementMapper().index(lowDimLocalAssembler.element())); |
| 720 | |||
| 721 | // update the corresponding vol vars in the bulk context | ||
| 722 | 85376 | const auto& bulkMap = couplingMapperPtr_->couplingMap(bulkGridId, lowDimGridId); | |
| 723 | 42688 | const auto& bulkCouplingEntry = bulkMap.at(bulkContext_.elementIdx); | |
| 724 | 42688 | const auto& couplingElementStencil = bulkCouplingEntry.couplingElementStencil; | |
| 725 | 128064 | auto it = std::find(couplingElementStencil.begin(), couplingElementStencil.end(), lowDimContext_.elementIdx); | |
| 726 |
4/8✗ Branch 0 not taken.
✓ Branch 1 taken 40512 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 40512 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 40512 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 40512 times.
|
170752 | assert(it != couplingElementStencil.end()); |
| 727 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 40512 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 40512 times.
|
85376 | const auto idxInContext = std::distance(couplingElementStencil.begin(), it); |
| 728 | |||
| 729 | // get neighboring bulk element from the bulk context (is the same element as first entry in low dim context) | ||
| 730 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 40512 times.
|
85248 | const auto& bulkElement = this->problem(bulkId).gridGeometry().element(bulkContext_.elementIdx); |
| 731 |
2/4✓ Branch 1 taken 40448 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 27680 times.
✗ Branch 5 not taken.
|
72480 | const auto& bulkFvGeometry = *lowDimContext_.bulkFvGeometries[0]; |
| 732 | |||
| 733 |
1/2✓ Branch 1 taken 40448 times.
✗ Branch 2 not taken.
|
42688 | auto& elemVolVars = bulkContext_.lowDimElemVolVars[idxInContext]; |
| 734 |
1/2✓ Branch 1 taken 40448 times.
✗ Branch 2 not taken.
|
42688 | const auto& fvGeom = bulkContext_.lowDimFvGeometries[idxInContext]; |
| 735 | 42688 | const auto& element = lowDimLocalAssembler.element(); | |
| 736 |
1/2✓ Branch 1 taken 40448 times.
✗ Branch 2 not taken.
|
42688 | const auto& scvfIndices = bulkCouplingEntry.elementToScvfMap.at(lowDimContext_.elementIdx); |
| 737 | |||
| 738 |
2/4✓ Branch 1 taken 40448 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 40448 times.
✗ Branch 5 not taken.
|
85376 | makeCoupledLowDimElemVolVars_(bulkElement, bulkFvGeometry, element, fvGeom, |
| 739 | this->curSol(lowDimId), scvfIndices, elemVolVars); | ||
| 740 | } | ||
| 741 | } | ||
| 742 | |||
| 743 | //! Empty stencil to be returned for elements that aren't coupled | ||
| 744 | template<std::size_t id, std::enable_if_t<(id == bulkId || id == lowDimId), int> = 0> | ||
| 745 | const typename CouplingMapper::template Stencil<id>& | ||
| 746 | ✗ | getEmptyStencil(Dune::index_constant<id>) const | |
| 747 | 652098 | { return std::get<(id == bulkId ? 0 : 1)>(emptyStencilTuple_); } | |
| 748 | |||
| 749 | private: | ||
| 750 | |||
| 751 | /*! | ||
| 752 | * \brief Prepares the volume variables in a lower-dimensional | ||
| 753 | * element coupled to the given bulk element. | ||
| 754 | */ | ||
| 755 | template<class LowDimSolution, class ScvfIdxStorage> | ||
| 756 | 455238 | void makeCoupledLowDimElemVolVars_(const Element<bulkId>& bulkElement, | |
| 757 | const FVElementGeometry<bulkId>& bulkFvGeometry, | ||
| 758 | const Element<lowDimId>& lowDimElement, | ||
| 759 | const FVElementGeometry<lowDimId>& lowDimFvGeometry, | ||
| 760 | const LowDimSolution& lowDimSol, | ||
| 761 | const ScvfIdxStorage& coupledBulkScvfIndices, | ||
| 762 | std::vector<VolumeVariables<lowDimId>>& elemVolVars) | ||
| 763 | { | ||
| 764 | 895564 | const auto lowDimElemSol = elementSolution(lowDimElement, lowDimSol, lowDimFvGeometry.gridGeometry()); | |
| 765 | |||
| 766 | // for tpfa, make only one volume variables object for the element | ||
| 767 | // for the update, use the first (and only - for tpfa) scv of the element | ||
| 768 | if (!lowDimUsesBox) | ||
| 769 | elemVolVars[0].update(lowDimElemSol, this->problem(lowDimId), lowDimElement, *scvs(lowDimFvGeometry).begin()); | ||
| 770 | // for box, make vol vars either: | ||
| 771 | // - in each scv center, which geometrically coincides with the scvf integration point (Neumann coupling) | ||
| 772 | // - at the vertex position (Dirichlet coupling) | ||
| 773 | else | ||
| 774 | { | ||
| 775 | // recall that the scvfs in the coupling map were ordered such | ||
| 776 | // that the i-th scvf in the map couples to the i-th lowdim scv | ||
| 777 |
4/4✓ Branch 0 taken 585580 times.
✓ Branch 1 taken 289062 times.
✓ Branch 2 taken 585580 times.
✓ Branch 3 taken 289062 times.
|
1835864 | for (unsigned int i = 0; i < coupledBulkScvfIndices.size(); ++i) |
| 778 | { | ||
| 779 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 585580 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 585580 times.
|
1850776 | const auto& scvf = bulkFvGeometry.scvf(coupledBulkScvfIndices[i]); |
| 780 |
3/4✗ Branch 0 not taken.
✓ Branch 1 taken 585580 times.
✓ Branch 2 taken 375448 times.
✓ Branch 3 taken 2400 times.
|
925388 | const auto bcTypes = this->problem(bulkId).interiorBoundaryTypes(bulkElement, scvf); |
| 781 |
4/4✓ Branch 0 taken 583180 times.
✓ Branch 1 taken 2400 times.
✓ Branch 2 taken 583180 times.
✓ Branch 3 taken 2400 times.
|
1850776 | if (bcTypes.hasOnlyNeumann()) |
| 782 |
1/2✗ Branch 2 not taken.
✓ Branch 3 taken 583180 times.
|
1845976 | FacetCoupling::makeInterpolatedVolVars(elemVolVars[i], this->problem(lowDimId), lowDimSol, |
| 783 | lowDimFvGeometry, lowDimElement, lowDimElement.geometry(), | ||
| 784 | scvf.ipGlobal()); | ||
| 785 | else | ||
| 786 |
3/6✗ Branch 0 not taken.
✓ Branch 1 taken 2400 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2400 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2400 times.
|
7200 | elemVolVars[i].update( lowDimElemSol, this->problem(lowDimId), lowDimElement, lowDimFvGeometry.scv(i) ); |
| 787 | } | ||
| 788 | } | ||
| 789 | 455238 | } | |
| 790 | |||
| 791 | //! evaluates the bulk-facet exchange fluxes for a given facet element | ||
| 792 | template<class BulkScvfIndices> | ||
| 793 | 905088 | NumEqVector<bulkId> evalBulkFluxes_(const Element<bulkId>& elementI, | |
| 794 | const FVElementGeometry<bulkId>& fvGeometry, | ||
| 795 | const ElementVolumeVariables<bulkId>& elemVolVars, | ||
| 796 | const ElementBoundaryTypes<bulkId>& elemBcTypes, | ||
| 797 | const ElementFluxVariablesCache<bulkId>& elemFluxVarsCache, | ||
| 798 | const LocalResidual<bulkId>& localResidual, | ||
| 799 | const BulkScvfIndices& scvfIndices) const | ||
| 800 | { | ||
| 801 | |||
| 802 | 905088 | NumEqVector<bulkId> coupledFluxes(0.0); | |
| 803 | 4397808 | for (const auto& scvfIdx : scvfIndices) | |
| 804 | 2320144 | coupledFluxes += localResidual.evalFlux(this->problem(bulkId), | |
| 805 | elementI, | ||
| 806 | fvGeometry, | ||
| 807 | elemVolVars, | ||
| 808 | elemBcTypes, | ||
| 809 | elemFluxVarsCache, | ||
| 810 | fvGeometry.scvf(scvfIdx)); | ||
| 811 | 905088 | return coupledFluxes; | |
| 812 | } | ||
| 813 | |||
| 814 | std::shared_ptr<CouplingMapper> couplingMapperPtr_; | ||
| 815 | |||
| 816 | //! store bools for all bulk elements that indicate if they | ||
| 817 | //! are coupled, so that we don't have to search in the map every time | ||
| 818 | std::vector<bool> bulkElemIsCoupled_; | ||
| 819 | |||
| 820 | //! empty stencil to return for non-coupled elements | ||
| 821 | using BulkStencil = typename CouplingMapper::template Stencil<bulkId>; | ||
| 822 | using LowDimStencil = typename CouplingMapper::template Stencil<lowDimId>; | ||
| 823 | std::tuple<BulkStencil, LowDimStencil> emptyStencilTuple_; | ||
| 824 | |||
| 825 | //! The coupling contexts of the two domains | ||
| 826 | BulkCouplingContext bulkContext_; | ||
| 827 | LowDimCouplingContext lowDimContext_; | ||
| 828 | }; | ||
| 829 | |||
| 830 | } // end namespace Dumux | ||
| 831 | |||
| 832 | #endif | ||
| 833 |