37    using Scalar = 
typename MDTraits::Scalar;
 
   41    static constexpr auto freeFlowMassIndex = 
typename MDTraits::template SubDomain<0>::Index();
 
   42    static constexpr auto porousMediumIndex = 
typename MDTraits::template SubDomain<1>::Index();
 
   48    using FreeFlowMassTypeTag = 
typename MDTraits::template SubDomain<freeFlowMassIndex>::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<freeFlowMassIndex>::Geometry::GlobalCoordinate;
 
   73    struct FreeFlowMassCouplingContext
 
   75        Element<porousMediumIndex> element;
 
   76        VolumeVariables<porousMediumIndex> volVars;
 
   77        FVElementGeometry<porousMediumIndex> fvGeometry;
 
   79        std::size_t freeFlowMassScvfIdx;
 
   80        std::size_t porousMediumScvfIdx;
 
   81        mutable VelocityVector velocity; 
 
   84    struct PorousMediumCouplingContext
 
   86        Element<freeFlowMassIndex> element;
 
   87        VolumeVariables<freeFlowMassIndex> volVars;
 
   88        FVElementGeometry<freeFlowMassIndex> fvGeometry;
 
   90        std::size_t porousMediumScvfIdx;
 
   91        std::size_t freeFlowMassScvfIdx;
 
   92        mutable VelocityVector velocity; 
 
   95    using CouplingMapper = DarcyDarcyBoundaryCouplingMapper<MDTraits>; 
 
  105    void init(std::shared_ptr<Problem<freeFlowMassIndex>> freeFlowMassProblem,
 
  106              std::shared_ptr<Problem<porousMediumIndex>> darcyProblem,
 
  109        this->
setSubProblems(std::make_tuple(freeFlowMassProblem, darcyProblem));
 
  112        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];
 
  149        auto& context = std::get<1-j>(couplingContext_);
 
  150        for (
auto& c : context)
 
  152            if (c.dofIdx == dofIdxGlobalJ)
 
  154                const auto elemSol = 
elementSolution(c.element, this->curSol(domainJ), this->problem(domainJ).gridGeometry());
 
  155                const auto& scv = *scvs(c.fvGeometry).begin();
 
  156                c.volVars.update(elemSol, this->
problem(domainJ), c.element, scv);
 
 
  166    template<std::
size_t i>
 
  168                                const FVElementGeometry<i>& fvGeometry,
 
  169                                const SubControlVolumeFace<i> scvf)
 const 
  171        auto& contexts = std::get<i>(couplingContext_);
 
  173        if (contexts.empty() || couplingContextBoundForElement_[i] != scvf.insideScvIdx())
 
  174            bindCouplingContext_(domainI, fvGeometry);
 
  176        for (
const auto& context : contexts)
 
  178            const auto expectedScvfIdx = domainI == 
freeFlowMassIndex ? context.freeFlowMassScvfIdx : context.porousMediumScvfIdx;
 
  179            if (scvf.index() == expectedScvfIdx)
 
  183        DUNE_THROW(Dune::InvalidStateException, 
"No coupling context found at scvf " << scvf.center());
 
 
  194    template<std::
size_t i, std::
size_t j>
 
  196                                           const Element<i>& element,
 
  197                                           Dune::index_constant<j> domainJ)
 const 
  199        const auto eIdx = this->
problem(domainI).gridGeometry().elementMapper().index(element);
 
  200        return couplingMapper_.couplingStencil(domainI, eIdx, domainJ);
 
 
  208    template<std::
size_t i>
 
  209    bool isCoupled(Dune::index_constant<i> domainI, 
const SubControlVolumeFace<i>& scvf)
 const 
  210    { 
return couplingMapper_.isCoupled(domainI, scvf); }
 
 
  216    template<std::
size_t i>
 
  217    bool isCoupledElement_(Dune::index_constant<i> domainI, std::size_t eIdx)
 const 
  221    template<std::
size_t i>
 
  222    VolumeVariables<i> volVars_(Dune::index_constant<i> domainI,
 
  223                                const Element<i>& element,
 
  224                                const SubControlVolume<i>& scv)
 const 
  226        VolumeVariables<i> volVars;
 
  228            element, this->
curSol(domainI), this->
problem(domainI).gridGeometry()
 
  230        volVars.update(elemSol, this->
problem(domainI), element, scv);
 
  237    template<std::
size_t i>
 
  238    void bindCouplingContext_(Dune::index_constant<i> domainI, 
const Element<i>& element)
 const 
  240        const auto fvGeometry = 
localView(this->
problem(domainI).gridGeometry()).bindElement(element);
 
  241        bindCouplingContext_(domainI, fvGeometry);
 
  247    template<std::
size_t i>
 
  248    void bindCouplingContext_(Dune::index_constant<i> domainI, 
const FVElementGeometry<i>& fvGeometry)
 const 
  250        auto& context = std::get<domainI>(couplingContext_);
 
  253        const auto eIdx = this->
problem(domainI).gridGeometry().elementMapper().index(fvGeometry.element());
 
  256        if (!isCoupledElement_(domainI, eIdx))
 
  259        couplingContextBoundForElement_[domainI] = eIdx;
 
  261        for (
const auto& scvf : scvfs(fvGeometry))
 
  265                const auto otherElementIdx = couplingMapper_.outsideElementIndex(domainI, scvf);
 
  266                constexpr auto domainJ = Dune::index_constant<1-domainI>();
 
  267                const auto& otherGridGeometry = this->
problem(domainJ).gridGeometry();
 
  268                const auto& otherElement = otherGridGeometry.element(otherElementIdx);
 
  269                auto otherFvGeometry = 
localView(otherGridGeometry).bindElement(otherElement);
 
  274                    volVars_(domainJ, otherElement, *std::begin(scvs(otherFvGeometry))),
 
  275                    std::move(otherFvGeometry),
 
  278                    couplingMapper_.flipScvfIndex(domainI, scvf),
 
  285    mutable std::tuple<std::vector<FreeFlowMassCouplingContext>, std::vector<PorousMediumCouplingContext>> couplingContext_;
 
  286    mutable std::array<std::size_t, 2> couplingContextBoundForElement_;
 
  288    CouplingMapper couplingMapper_;
 
 
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
The interface of the coupling manager for multi domain problems.