13#ifndef DUMUX_MULTIDOMAIN_BOUNDARY_FFPM_FFMMOMENTUMPM_COUPLINGMANAGER_STAGGERED_TPFA_HH
14#define DUMUX_MULTIDOMAIN_BOUNDARY_FFPM_FFMMOMENTUMPM_COUPLINGMANAGER_STAGGERED_TPFA_HH
19#include <dune/common/float_cmp.hh>
20#include <dune/common/exceptions.hh>
34template<
class MDTraits>
38 using Scalar =
typename MDTraits::Scalar;
43 static constexpr auto porousMediumIndex =
typename MDTraits::template SubDomain<1>::Index();
48 using FreeFlowMomentumTypeTag =
typename MDTraits::template SubDomain<freeFlowMomentumIndex>::TypeTag;
49 using PorousMediumTypeTag =
typename MDTraits::template SubDomain<porousMediumIndex>::TypeTag;
51 using CouplingStencils = std::unordered_map<std::size_t, std::vector<std::size_t> >;
52 using CouplingStencil = CouplingStencils::mapped_type;
55 template<std::
size_t id>
56 using SubDomainTypeTag =
typename MDTraits::template SubDomain<id>::TypeTag;
64 template<std::
size_t id>
using FVElementGeometry =
typename GridGeometry<id>::LocalView;
66 template<std::
size_t id>
using Element =
typename GridView<id>::template Codim<0>::Entity;
68 template<std::
size_t id>
using SubControlVolumeFace =
typename FVElementGeometry<id>::SubControlVolumeFace;
69 template<std::
size_t id>
using SubControlVolume =
typename FVElementGeometry<id>::SubControlVolume;
71 using VelocityVector =
typename Element<freeFlowMomentumIndex>::Geometry::GlobalCoordinate;
73 struct FreeFlowMomentumCouplingContext
75 Element<porousMediumIndex> element;
76 VolumeVariables<porousMediumIndex> volVars;
77 FVElementGeometry<porousMediumIndex> fvGeometry;
78 std::size_t freeFlowMomentumScvfIdx;
79 std::size_t porousMediumScvfIdx;
80 std::size_t porousMediumDofIdx;
81 VelocityVector gravity;
84 struct PorousMediumCouplingContext
86 Element<freeFlowMomentumIndex> element;
87 FVElementGeometry<freeFlowMomentumIndex> fvGeometry;
88 std::size_t porousMediumScvfIdx;
89 std::size_t freeFlowMomentumScvfIdx;
90 std::size_t freeFlowMomentumDofIdx;
94 using CouplingMapper = FFMomentumPMCouplingMapperStaggeredCCTpfa<MDTraits, FFMomentumPMCouplingManagerStaggeredCCTpfa<MDTraits>>;
104 void init(std::shared_ptr<Problem<freeFlowMomentumIndex>> freeFlowMomentumProblem,
105 std::shared_ptr<Problem<porousMediumIndex>> porousMediumProblem,
108 this->
setSubProblems(std::make_tuple(freeFlowMomentumProblem, porousMediumProblem));
111 couplingMapper_.
update(*
this);
125 template<std::
size_t i,
class Assembler>
126 void bindCouplingContext(Dune::index_constant<i> domainI,
const Element<i>& element,
const Assembler& assembler)
const
128 bindCouplingContext_(domainI, element);
134 template<std::
size_t i, std::
size_t j,
class LocalAssemblerI>
136 const LocalAssemblerI& localAssemblerI,
137 Dune::index_constant<j> domainJ,
138 std::size_t dofIdxGlobalJ,
139 const PrimaryVariables<j>& priVarsJ,
142 this->
curSol(domainJ)[dofIdxGlobalJ][pvIdxJ] = priVarsJ[pvIdxJ];
154 auto& context = std::get<porousMediumIndex>(couplingContext_);
155 for (
auto& c : context)
157 if (c.freeFlowMomentumDofIdx == dofIdxGlobalJ)
159 const auto& scvf = localAssemblerI.fvGeometry().scvf(c.porousMediumScvfIdx);
160 c.faceVelocity =
faceVelocity(localAssemblerI.element(), scvf);
169 auto& context = std::get<freeFlowMomentumIndex>(couplingContext_);
170 for (
auto& c : context)
172 if (c.porousMediumDofIdx == dofIdxGlobalJ)
174 const auto& ggJ = c.fvGeometry.gridGeometry();
175 const auto& scv = *scvs(c.fvGeometry).begin();
176 const auto elemSol =
elementSolution(c.element, this->curSol(domainJ), ggJ);
177 c.volVars.update(elemSol, this->
problem(domainJ), c.element, scv);
188 template<std::
size_t i>
190 const FVElementGeometry<i>& fvGeometry,
191 const SubControlVolumeFace<i> scvf)
const
193 auto& contexts = std::get<i>(couplingContext_);
195 if (contexts.empty() || couplingContextBoundForElement_[i] != scvf.insideScvIdx())
196 bindCouplingContext_(domainI, fvGeometry);
198 for (
const auto& context : contexts)
200 const auto expectedScvfIdx = domainI ==
freeFlowMomentumIndex ? context.freeFlowMomentumScvfIdx : context.porousMediumScvfIdx;
201 if (scvf.index() == expectedScvfIdx)
205 DUNE_THROW(Dune::InvalidStateException,
"No coupling context found at scvf " << scvf.center());
216 const CouplingStencil&
couplingStencil(Dune::index_constant<porousMediumIndex> domainI,
217 const Element<porousMediumIndex>& element,
218 Dune::index_constant<freeFlowMomentumIndex> domainJ)
const
220 const auto eIdx = this->
problem(domainI).gridGeometry().elementMapper().index(element);
233 const CouplingStencil&
couplingStencil(Dune::index_constant<freeFlowMomentumIndex> domainI,
234 const Element<freeFlowMomentumIndex>& elementI,
235 const SubControlVolume<freeFlowMomentumIndex>& scvI,
236 Dune::index_constant<porousMediumIndex> domainJ)
const
238 return couplingMapper_.
couplingStencil(domainI, elementI, scvI, domainJ);
247 template<
class LocalAssemblerI, std::
size_t j>
249 const LocalAssemblerI& localAssemblerI,
250 const SubControlVolume<freeFlowMomentumIndex>& scvI,
251 Dune::index_constant<j> domainJ,
252 std::size_t dofIdxGlobalJ)
const
254 return localAssemblerI.evalLocalResidual();
262 bool isCoupledLateralScvf(Dune::index_constant<freeFlowMomentumIndex> domainI,
const SubControlVolumeFace<freeFlowMomentumIndex>& scvf)
const
268 template<std::
size_t i>
269 bool isCoupled(Dune::index_constant<i> domainI,
const SubControlVolumeFace<i>& scvf)
const
274 return couplingMapper_.
isCoupled(domainI, scvf);
282 bool isCoupled(Dune::index_constant<freeFlowMomentumIndex> domainI,
283 const SubControlVolume<freeFlowMomentumIndex>& scv)
const
284 {
return couplingMapper_.
isCoupled(domainI, scv); }
290 const SubControlVolumeFace<porousMediumIndex>& scvf)
const
293 auto velocity = scvf.unitOuterNormal();
295 std::for_each(velocity.begin(), velocity.end(), [](
auto& v){ v = abs(v); });
299 couplingMapper_.
outsideDofIndex(Dune::index_constant<porousMediumIndex>(), scvf)
307 template<std::
size_t i>
308 VolumeVariables<i> volVars_(Dune::index_constant<i> domainI,
309 const Element<i>& element,
310 const SubControlVolume<i>& scv)
const
312 VolumeVariables<i> volVars;
314 element, this->
curSol(domainI), this->
problem(domainI).gridGeometry()
316 volVars.update(elemSol, this->
problem(domainI), element, scv);
323 template<std::
size_t i>
324 bool isCoupledElement_(Dune::index_constant<i> domainI, std::size_t eIdx)
const
330 template<std::
size_t i>
331 void bindCouplingContext_(Dune::index_constant<i> domainI,
const Element<i>& element)
const
333 const auto fvGeometry =
localView(this->
problem(domainI).gridGeometry()).bindElement(element);
334 bindCouplingContext_(domainI, fvGeometry);
340 template<std::
size_t i>
341 void bindCouplingContext_(Dune::index_constant<i> domainI,
const FVElementGeometry<i>& fvGeometry)
const
343 auto& context = std::get<domainI>(couplingContext_);
346 const auto eIdx = this->
problem(domainI).gridGeometry().elementMapper().index(fvGeometry.element());
349 if (!isCoupledElement_(domainI, eIdx))
352 couplingContextBoundForElement_[domainI] = eIdx;
354 for (
const auto& scvf : scvfs(fvGeometry))
356 if (couplingMapper_.
isCoupled(domainI, scvf))
361 constexpr auto domainJ = Dune::index_constant<1-i>();
362 const auto& otherGridGeometry = this->
problem(domainJ).gridGeometry();
363 const auto& otherElement = otherGridGeometry.element(otherElementIdx);
364 auto otherFvGeometry =
localView(otherGridGeometry).bindElement(otherElement);
369 volVars_(domainJ, otherElement, *std::begin(scvs(otherFvGeometry))),
370 std::move(otherFvGeometry),
374 this->
problem(domainJ).spatialParams().gravity(scvf.center())
381 constexpr auto domainJ = Dune::index_constant<1-i>();
382 const auto& otherGridGeometry = this->
problem(domainJ).gridGeometry();
383 const auto& otherElement = otherGridGeometry.element(otherElementIdx);
384 auto otherFvGeometry =
localView(otherGridGeometry).bindElement(otherElement);
386 const auto otherScvfIdx = couplingMapper_.
flipScvfIndex(domainI, scvf);
387 const auto& otherScvf = otherFvGeometry.scvf(otherScvfIdx);
388 const auto& otherScv = otherFvGeometry.scv(otherScvf.insideScvIdx());
392 std::move(otherFvGeometry),
403 mutable std::tuple<std::vector<FreeFlowMomentumCouplingContext>, std::vector<PorousMediumCouplingContext>> couplingContext_;
404 mutable std::array<std::size_t, 2> couplingContextBoundForElement_;
406 CouplingMapper couplingMapper_;
The interface of the coupling manager for multi domain problems.
Definition: multidomain/couplingmanager.hh:48
void attachSolution(SolutionVectorStorage &curSol)
Attach a solution vector stored outside of this class.
Definition: multidomain/couplingmanager.hh:322
void setSubProblems(const std::tuple< std::shared_ptr< SubProblems >... > &problems)
set the pointers to the sub problems
Definition: multidomain/couplingmanager.hh:287
const Problem< i > & problem(Dune::index_constant< i > domainIdx) const
Return a reference to the sub problem.
Definition: multidomain/couplingmanager.hh:309
SubSolutionVector< i > & curSol(Dune::index_constant< i > domainIdx)
the solution vector of the subproblem
Definition: multidomain/couplingmanager.hh:338
typename Traits::template TupleOfSharedPtr< SubSolutionVector > SolutionVectorStorage
the type in which the solution vector is stored in the manager
Definition: multidomain/couplingmanager.hh:71
Coupling manager for Stokes and Darcy domains with equal dimension specialization for staggered-cctpf...
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:37
bool isCoupled(Dune::index_constant< i > domainI, const SubControlVolumeFace< i > &scvf) const
Returns whether a given scvf is coupled to the other domain.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:269
bool isCoupled(Dune::index_constant< freeFlowMomentumIndex > domainI, const SubControlVolume< freeFlowMomentumIndex > &scv) const
If the boundary entity is on a coupling boundary.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:282
bool isCoupledLateralScvf(Dune::index_constant< freeFlowMomentumIndex > domainI, const SubControlVolumeFace< freeFlowMomentumIndex > &scvf) const
Returns whether a given scvf is coupled to the other domain.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:262
void bindCouplingContext(Dune::index_constant< i > domainI, const Element< i > &element, const Assembler &assembler) const
Methods to be accessed by the assembly.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:126
const auto & couplingContext(Dune::index_constant< i > domainI, const FVElementGeometry< i > &fvGeometry, const SubControlVolumeFace< i > scvf) const
Access the coupling context needed for the Stokes domain.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:189
const CouplingStencil & couplingStencil(Dune::index_constant< freeFlowMomentumIndex > domainI, const Element< freeFlowMomentumIndex > &elementI, const SubControlVolume< freeFlowMomentumIndex > &scvI, Dune::index_constant< porousMediumIndex > domainJ) const
returns an iterable container of all indices of degrees of freedom of domain j that couple with / inf...
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:233
void init(std::shared_ptr< Problem< freeFlowMomentumIndex > > freeFlowMomentumProblem, std::shared_ptr< Problem< porousMediumIndex > > porousMediumProblem, SolutionVectorStorage &curSol)
Methods to be accessed by main.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:104
typename ParentType::SolutionVectorStorage SolutionVectorStorage
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:45
static constexpr auto porousMediumIndex
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:43
void updateCouplingContext(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ, const PrimaryVariables< j > &priVarsJ, int pvIdxJ)
Update the coupling context.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:135
auto faceVelocity(const Element< porousMediumIndex > &element, const SubControlVolumeFace< porousMediumIndex > &scvf) const
Returns the velocity at a given sub control volume face.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:289
const CouplingStencil & couplingStencil(Dune::index_constant< porousMediumIndex > domainI, const Element< porousMediumIndex > &element, Dune::index_constant< freeFlowMomentumIndex > domainJ) const
The coupling stencils.
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:216
static constexpr auto freeFlowMomentumIndex
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:42
decltype(auto) evalCouplingResidual(Dune::index_constant< freeFlowMomentumIndex > domainI, const LocalAssemblerI &localAssemblerI, const SubControlVolume< freeFlowMomentumIndex > &scvI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ) const
evaluate the coupling residual special interface for fcstaggered methods
Definition: ffmomentumpm/couplingmanager_staggered_cctpfa.hh:248
void update(const CouplingManager &couplingManager)
Main update routine.
Definition: couplingmapper_staggered_cctpfa.hh:81
bool isCoupled(Dune::index_constant< i > domainI, const SubControlVolumeFace< i > &scvf) const
If the boundary entity is on a coupling boundary.
Definition: couplingmapper_staggered_cctpfa.hh:259
const std::vector< std::size_t > & couplingStencil(Dune::index_constant< CouplingManager::porousMediumIndex > domainI, const std::size_t eIdxI, Dune::index_constant< CouplingManager::freeFlowMomentumIndex > domainJ) const
returns an iterable container of all indices of degrees of freedom of domain j that couple with / inf...
Definition: couplingmapper_staggered_cctpfa.hh:205
std::size_t flipScvfIndex(Dune::index_constant< i > domainI, const SubControlVolumeFace< i > &scvf) const
Return the scvf index of the flipped scvf in the other domain.
Definition: couplingmapper_staggered_cctpfa.hh:289
bool isCoupledElement(Dune::index_constant< i >, std::size_t eIdx) const
Return if an element residual with index eIdx of domain i is coupled to domain j.
Definition: couplingmapper_staggered_cctpfa.hh:245
std::size_t outsideDofIndex(Dune::index_constant< i > domainI, const SubControlVolumeFace< i > &scvf) const
Return the outside element index (the element index of the other domain)
Definition: couplingmapper_staggered_cctpfa.hh:313
std::size_t outsideElementIndex(Dune::index_constant< i > domainI, const SubControlVolumeFace< i > &scvf) const
Return the outside element index (the element index of the other domain)
Definition: couplingmapper_staggered_cctpfa.hh:301
bool isCoupledLateralScvf(Dune::index_constant< CouplingManager::freeFlowMomentumIndex > domainI, const SubControlVolumeFace< CouplingManager::freeFlowMomentumIndex > &scvf) const
If the boundary entity is on a coupling boundary.
Definition: couplingmapper_staggered_cctpfa.hh:270
Defines all properties used in Dumux.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:26
auto elementSolution(const Element &element, const SolutionVector &sol, const GridGeometry &gg) -> std::enable_if_t< GridGeometry::discMethod==DiscretizationMethods::cctpfa||GridGeometry::discMethod==DiscretizationMethods::ccmpfa, CCElementSolution< typename GridGeometry::LocalView, std::decay_t< decltype(std::declval< SolutionVector >()[0])> > >
Make an element solution for cell-centered schemes.
Definition: cellcentered/elementsolution.hh:101
decltype(auto) evalCouplingResidual(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ) const
evaluates the element residual of a coupled element of domain i which depends on the variables at the...
Definition: multidomain/couplingmanager.hh:249
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
The interface of the coupling manager for multi domain problems.
The local element solution class for staggered methods.