12#ifndef DUMUX_MULTIDOMAIN_FREEFLOW_COUPLING_MANAGER_CVFE_HH 
   13#define DUMUX_MULTIDOMAIN_FREEFLOW_COUPLING_MANAGER_CVFE_HH 
   20#include <dune/common/exceptions.hh> 
   21#include <dune/common/indices.hh> 
   22#include <dune/common/float_cmp.hh> 
   23#include <dune/geometry/referenceelements.hh> 
   60    template<std::
size_t id> 
using SubDomainTypeTag = 
typename Traits::template SubDomain<id>::TypeTag;
 
   63    template<std::
size_t id> 
using GridView = 
typename GridGeometry<id>::GridView;
 
   64    template<std::
size_t id> 
using Element = 
typename GridView<id>::template Codim<0>::Entity;
 
   65    template<std::
size_t id> 
using ElementSeed = 
typename GridView<id>::Grid::template Codim<0>::EntitySeed;
 
   66    template<std::
size_t id> 
using FVElementGeometry = 
typename GridGeometry<id>::LocalView;
 
   67    template<std::
size_t id> 
using SubControlVolume = 
typename FVElementGeometry<id>::SubControlVolume;
 
   68    template<std::
size_t id> 
using SubControlVolumeFace = 
typename FVElementGeometry<id>::SubControlVolumeFace;
 
   69    template<std::
size_t id> 
using GridVariables = 
typename Traits::template SubDomain<id>::GridVariables;
 
   70    template<std::
size_t id> 
using ElementVolumeVariables = 
typename GridVariables<id>::GridVolumeVariables::LocalView;
 
   71    template<std::
size_t id> 
using GridFluxVariablesCache = 
typename GridVariables<id>::GridFluxVariablesCache;
 
   75    using Scalar = 
typename Traits::Scalar;
 
   76    using SolutionVector = 
typename Traits::SolutionVector;
 
   78    using CouplingStencilType = std::vector<std::size_t>;
 
   80    using GridVariablesTuple = 
typename Traits::template TupleOfSharedPtr<GridVariables>;
 
   82    using FluidSystem = 
typename VolumeVariables<freeFlowMassIndex>::FluidSystem;
 
   84    using VelocityVector = 
typename SubControlVolumeFace<freeFlowMassIndex>::GlobalPosition;
 
   85    using ShapeValue = 
typename Dune::FieldVector<Scalar, 1>;
 
   87    static_assert(std::is_same_v<VelocityVector, typename SubControlVolumeFace<freeFlowMomentumIndex>::GlobalPosition>);
 
   89    struct MomentumCouplingContext
 
   91        FVElementGeometry<freeFlowMassIndex> fvGeometry;
 
   92        ElementVolumeVariables<freeFlowMassIndex> curElemVolVars;
 
   93        ElementVolumeVariables<freeFlowMassIndex> prevElemVolVars;
 
   97    struct MassAndEnergyCouplingContext
 
   99        MassAndEnergyCouplingContext(FVElementGeometry<freeFlowMomentumIndex>&& f, 
const std::size_t i)
 
  100        : fvGeometry(std::move(f))
 
  104        FVElementGeometry<freeFlowMomentumIndex> fvGeometry;
 
  108    using MomentumDiscretizationMethod = 
typename GridGeometry<freeFlowMomentumIndex>::DiscretizationMethod;
 
  109    using MassDiscretizationMethod = 
typename GridGeometry<freeFlowMassIndex>::DiscretizationMethod;
 
  113    static constexpr auto pressureIdx = VolumeVariables<freeFlowMassIndex>::Indices::pressureIdx;
 
  121    void init(std::shared_ptr<Problem<freeFlowMomentumIndex>> momentumProblem,
 
  122              std::shared_ptr<Problem<freeFlowMassIndex>> massProblem,
 
  123              GridVariablesTuple&& gridVariables,
 
  124              const SolutionVector& 
curSol)
 
  126        this->momentumCouplingContext_().clear();
 
  127        this->massAndEnergyCouplingContext_().clear();
 
  129        this->
setSubProblems(std::make_tuple(momentumProblem, massProblem));
 
  130        gridVariables_ = gridVariables;
 
  133        computeCouplingStencils_();
 
 
  137    void init(std::shared_ptr<Problem<freeFlowMomentumIndex>> momentumProblem,
 
  138              std::shared_ptr<Problem<freeFlowMassIndex>> massProblem,
 
  139              GridVariablesTuple&& gridVariables,
 
  140              const SolutionVector& 
curSol,
 
  141              const SolutionVector& prevSol)
 
  143        init(momentumProblem, massProblem, std::forward<GridVariablesTuple>(gridVariables), 
curSol);
 
 
  149    void init(std::shared_ptr<Problem<freeFlowMomentumIndex>> momentumProblem,
 
  150              std::shared_ptr<Problem<freeFlowMassIndex>> massProblem,
 
  151              GridVariablesTuple&& gridVariables,
 
  154        this->momentumCouplingContext_().clear();
 
  155        this->massAndEnergyCouplingContext_().clear();
 
  157        this->
setSubProblems(std::make_tuple(momentumProblem, massProblem));
 
  158        gridVariables_ = gridVariables;
 
  161        computeCouplingStencils_();
 
 
  174    Scalar 
pressure(
const Element<freeFlowMomentumIndex>& element,
 
  175                    const FVElementGeometry<freeFlowMomentumIndex>& fvGeometry,
 
  176                    const SubControlVolumeFace<freeFlowMomentumIndex>& scvf,
 
  177                    const bool considerPreviousTimeStep = 
false)
 const 
  179        assert(!(considerPreviousTimeStep && !this->isTransient_));
 
 
  190    Scalar 
pressure(
const Element<freeFlowMomentumIndex>& element,
 
  191                    const FVElementGeometry<freeFlowMomentumIndex>& fvGeometry,
 
  192                    const SubControlVolume<freeFlowMomentumIndex>& scv,
 
  193                    const bool considerPreviousTimeStep = 
false)
 const 
  195        assert(!(considerPreviousTimeStep && !this->isTransient_));
 
 
  206    Scalar 
density(
const Element<freeFlowMomentumIndex>& element,
 
  207                   const FVElementGeometry<freeFlowMomentumIndex>& fvGeometry,
 
  208                   const SubControlVolumeFace<freeFlowMomentumIndex>& scvf,
 
  209                   const bool considerPreviousTimeStep = 
false)
 const 
  211        assert(!(considerPreviousTimeStep && !this->isTransient_));
 
  212        bindCouplingContext_(Dune::index_constant<freeFlowMomentumIndex>(), element, fvGeometry.elementIndex());
 
  216            const auto eIdx = fvGeometry.elementIndex();
 
  217            const auto& scv = this->momentumCouplingContext_()[0].fvGeometry.scv(eIdx);
 
  219            const auto& volVars = considerPreviousTimeStep ?
 
  220                this->momentumCouplingContext_()[0].prevElemVolVars[scv]
 
  221                : this->momentumCouplingContext_()[0].curElemVolVars[scv];
 
  223            return volVars.density();
 
  229            using ShapeValue = 
typename Dune::FieldVector<Scalar, 1>;
 
  230            const auto& localBasis = this->momentumCouplingContext_()[0].fvGeometry.feLocalBasis();
 
  231            std::vector<ShapeValue> shapeValues;
 
  232            const auto ipLocal = element.geometry().local(scvf.ipGlobal());
 
  233            localBasis.evaluateFunction(ipLocal, shapeValues);
 
  236            for (
const auto& scv : scvs(this->momentumCouplingContext_()[0].fvGeometry))
 
  238                const auto& volVars = considerPreviousTimeStep ?
 
  239                    this->momentumCouplingContext_()[0].prevElemVolVars[scv]
 
  240                    : this->momentumCouplingContext_()[0].curElemVolVars[scv];
 
  241                rho += volVars.density()*shapeValues[scv.indexInElement()][0];
 
  247            DUNE_THROW(Dune::NotImplemented,
 
  248                "Density interpolation for discretization scheme " << MassDiscretizationMethod{}
 
 
  255    Scalar 
density(
const Element<freeFlowMomentumIndex>& element,
 
  256                   const FVElementGeometry<freeFlowMomentumIndex>& fvGeometry,
 
  257                   const SubControlVolume<freeFlowMomentumIndex>& scv,
 
  258                   const bool considerPreviousTimeStep = 
false)
 const 
  260        assert(!(considerPreviousTimeStep && !this->isTransient_));
 
  261        bindCouplingContext_(Dune::index_constant<freeFlowMomentumIndex>(), element, scv.elementIndex());
 
  265            const auto eIdx = scv.elementIndex();
 
  266            const auto& scvI = this->momentumCouplingContext_()[0].fvGeometry.scv(eIdx);
 
  268            const auto& volVars = considerPreviousTimeStep ?
 
  269                this->momentumCouplingContext_()[0].prevElemVolVars[scvI]
 
  270                : this->momentumCouplingContext_()[0].curElemVolVars[scvI];
 
  272            return volVars.density();
 
  278            using ShapeValue = 
typename Dune::FieldVector<Scalar, 1>;
 
  279            const auto& localBasis = this->momentumCouplingContext_()[0].fvGeometry.feLocalBasis();
 
  280            std::vector<ShapeValue> shapeValues;
 
  281            const auto ipLocal = element.geometry().local(scv.dofPosition());
 
  282            localBasis.evaluateFunction(ipLocal, shapeValues);
 
  285            for (
const auto& scvI : scvs(this->momentumCouplingContext_()[0].fvGeometry))
 
  287                const auto& volVars = considerPreviousTimeStep ?
 
  288                    this->momentumCouplingContext_()[0].prevElemVolVars[scvI]
 
  289                    : this->momentumCouplingContext_()[0].curElemVolVars[scvI];
 
  290                rho += volVars.density()*shapeValues[scvI.indexInElement()][0];
 
  295            DUNE_THROW(Dune::NotImplemented,
 
  296                "Density interpolation for discretization scheme " << MassDiscretizationMethod{}
 
 
  304                              const FVElementGeometry<freeFlowMomentumIndex>& fvGeometry,
 
  305                              const SubControlVolumeFace<freeFlowMomentumIndex>& scvf,
 
  306                              const bool considerPreviousTimeStep = 
false)
 const 
  308        assert(!(considerPreviousTimeStep && !this->isTransient_));
 
  309        bindCouplingContext_(Dune::index_constant<freeFlowMomentumIndex>(), element, fvGeometry.elementIndex());
 
  313            const auto eIdx = fvGeometry.elementIndex();
 
  314            const auto& scv = this->momentumCouplingContext_()[0].fvGeometry.scv(eIdx);
 
  315            const auto& volVars = considerPreviousTimeStep ?
 
  316                this->momentumCouplingContext_()[0].prevElemVolVars[scv]
 
  317                : this->momentumCouplingContext_()[0].curElemVolVars[scv];
 
  318            return volVars.viscosity();
 
  324            using ShapeValue = 
typename Dune::FieldVector<Scalar, 1>;
 
  325            const auto& localBasis = this->momentumCouplingContext_()[0].fvGeometry.feLocalBasis();
 
  326            std::vector<ShapeValue> shapeValues;
 
  327            const auto ipLocal = element.geometry().local(scvf.ipGlobal());
 
  328            localBasis.evaluateFunction(ipLocal, shapeValues);
 
  331            for (
const auto& scv : scvs(this->momentumCouplingContext_()[0].fvGeometry))
 
  333                const auto& volVars = considerPreviousTimeStep ?
 
  334                    this->momentumCouplingContext_()[0].prevElemVolVars[scv]
 
  335                    : this->momentumCouplingContext_()[0].curElemVolVars[scv];
 
  336                mu += volVars.viscosity()*shapeValues[scv.indexInElement()][0];
 
  342            DUNE_THROW(Dune::NotImplemented,
 
  343                "Viscosity interpolation for discretization scheme " << MassDiscretizationMethod{}
 
 
  351                              const FVElementGeometry<freeFlowMomentumIndex>& fvGeometry,
 
  352                              const SubControlVolume<freeFlowMomentumIndex>& scv,
 
  353                              const bool considerPreviousTimeStep = 
false)
 const 
  355        assert(!(considerPreviousTimeStep && !this->isTransient_));
 
  356        bindCouplingContext_(Dune::index_constant<freeFlowMomentumIndex>(), element, fvGeometry.elementIndex());
 
  360            const auto eIdx = fvGeometry.elementIndex();
 
  361            const auto& scvI = this->momentumCouplingContext_()[0].fvGeometry.scv(eIdx);
 
  362            const auto& volVars = considerPreviousTimeStep ?
 
  363                this->momentumCouplingContext_()[0].prevElemVolVars[scvI]
 
  364                : this->momentumCouplingContext_()[0].curElemVolVars[scvI];
 
  365            return volVars.viscosity();
 
  371            using ShapeValue = 
typename Dune::FieldVector<Scalar, 1>;
 
  372            const auto& localBasis = this->momentumCouplingContext_()[0].fvGeometry.feLocalBasis();
 
  373            std::vector<ShapeValue> shapeValues;
 
  374            const auto ipLocal = element.geometry().local(scv.dofPosition());
 
  375            localBasis.evaluateFunction(ipLocal, shapeValues);
 
  378            for (
const auto& scvI : scvs(this->momentumCouplingContext_()[0].fvGeometry))
 
  380                const auto& volVars = considerPreviousTimeStep ?
 
  381                    this->momentumCouplingContext_()[0].prevElemVolVars[scvI]
 
  382                    : this->momentumCouplingContext_()[0].curElemVolVars[scvI];
 
  383                mu += volVars.viscosity()*shapeValues[scvI.indexInElement()][0];
 
  389            DUNE_THROW(Dune::NotImplemented,
 
  390                "Viscosity interpolation for discretization scheme " << MassDiscretizationMethod{}
 
 
  398                                const SubControlVolumeFace<freeFlowMassIndex>& scvf)
 const 
  403        bindCouplingContext_(Dune::index_constant<freeFlowMassIndex>(), element, eIdx);
 
  405        const auto& fvGeometry = this->massAndEnergyCouplingContext_()[0].fvGeometry;
 
  406        const auto& localBasis = fvGeometry.feLocalBasis();
 
  408        std::vector<ShapeValue> shapeValues;
 
  409        const auto ipLocal = element.geometry().local(scvf.ipGlobal());
 
  410        localBasis.evaluateFunction(ipLocal, shapeValues);
 
  413        VelocityVector velocity(0.0);
 
  414        for (
const auto& scv : scvs(fvGeometry))
 
  415            velocity.axpy(shapeValues[scv.localDofIndex()][0], this->curSol(
freeFlowMomentumIndex)[scv.dofIndex()]);
 
 
  423    VelocityVector 
elementVelocity(
const FVElementGeometry<freeFlowMassIndex>& fvGeometry)
 const 
  425        bindCouplingContext_(Dune::index_constant<freeFlowMassIndex>(), fvGeometry.element());
 
  427        const auto& momentumFvGeometry = this->massAndEnergyCouplingContext_()[0].fvGeometry;
 
  428        const auto& localBasis = momentumFvGeometry.feLocalBasis();
 
  431        VelocityVector velocity(0.0);
 
  432        std::vector<ShapeValue> shapeValues;
 
  433        localBasis.evaluateFunction(referenceElement(fvGeometry.element()).position(0,0), shapeValues);
 
  435        for (
const auto& scv : scvs(momentumFvGeometry))
 
  436            velocity.axpy(shapeValues[scv.localDofIndex()][0], this->curSol(
freeFlowMomentumIndex)[scv.dofIndex()]);
 
 
  445    template<std::
size_t j>
 
  446    const CouplingStencilType& 
couplingStencil(Dune::index_constant<freeFlowMomentumIndex> domainI,
 
  447                                               const Element<freeFlowMomentumIndex>& elementI,
 
  448                                               const SubControlVolume<freeFlowMomentumIndex>& scvI,
 
  449                                               Dune::index_constant<j> domainJ)
 const 
  450    { 
return emptyStencil_; }
 
 
  466    const CouplingStencilType& 
couplingStencil(Dune::index_constant<freeFlowMassIndex> domainI,
 
  467                                              const Element<freeFlowMassIndex>& elementI,
 
  468                                              Dune::index_constant<freeFlowMomentumIndex> domainJ)
 const 
  471        return massAndEnergyToMomentumStencils_[eIdx];
 
 
  482    const CouplingStencilType& 
couplingStencil(Dune::index_constant<freeFlowMomentumIndex> domainI,
 
  483                                               const Element<freeFlowMomentumIndex>& elementI,
 
  484                                               Dune::index_constant<freeFlowMassIndex> domainJ)
 const 
  487        return momentumToMassAndEnergyStencils_[eIdx];
 
 
  516    template<std::
size_t i, std::
size_t j, 
class LocalAssemblerI>
 
  518                               const LocalAssemblerI& localAssemblerI,
 
  519                               Dune::index_constant<j> domainJ,
 
  520                               std::size_t dofIdxGlobalJ,
 
  521                               const PrimaryVariables<j>& priVarsJ,
 
  524        this->
curSol(domainJ)[dofIdxGlobalJ][pvIdxJ] = priVarsJ[pvIdxJ];
 
  530                bindCouplingContext_(domainI, localAssemblerI.element());
 
  533                const auto& deflectedElement = 
problem.gridGeometry().element(dofIdxGlobalJ);
 
  535                const auto& fvGeometry = momentumCouplingContext_()[0].fvGeometry;
 
  536                const auto& scv = fvGeometry.scv(dofIdxGlobalJ);
 
  538                if constexpr (ElementVolumeVariables<freeFlowMassIndex>::GridVolumeVariables::cachingEnabled)
 
  539                    gridVars_(
freeFlowMassIndex).curGridVolVars().volVars(scv).update(std::move(elemSol), 
problem, deflectedElement, scv);
 
  541                    momentumCouplingContext_()[0].curElemVolVars[scv].update(std::move(elemSol), 
problem, deflectedElement, scv);
 
  549                bindCouplingContext_(domainI, localAssemblerI.element());
 
  552                const auto& deflectedElement = 
problem.gridGeometry().element(this->momentumCouplingContext_()[0].eIdx);
 
  554                const auto& fvGeometry = this->momentumCouplingContext_()[0].fvGeometry;
 
  556                for (
const auto& scv : scvs(fvGeometry))
 
  558                    if(scv.dofIndex() == dofIdxGlobalJ)
 
  560                        if constexpr (ElementVolumeVariables<freeFlowMassIndex>::GridVolumeVariables::cachingEnabled)
 
  561                            this->gridVars_(
freeFlowMassIndex).curGridVolVars().volVars(scv).update(std::move(elemSol), 
problem, deflectedElement, scv);
 
  563                           this->momentumCouplingContext_()[0].curElemVolVars[scv].update(std::move(elemSol), 
problem, deflectedElement, scv);
 
  569            DUNE_THROW(Dune::NotImplemented,
 
  570                "Context update for discretization scheme " << MassDiscretizationMethod{}
 
 
  600    template<std::
size_t i, 
class AssembleElementFunc>
 
  603        if (elementSets_.empty())
 
  604            DUNE_THROW(Dune::InvalidStateException, 
"Call computeColorsForAssembly before assembling in parallel!");
 
  611        for (
const auto& elements : elementSets_)
 
  615                const auto element = grid.entity(elements[eIdx]);
 
  616                assembleElement(element);
 
 
  622    void bindCouplingContext_(Dune::index_constant<freeFlowMomentumIndex> domainI,
 
  623                              const Element<freeFlowMomentumIndex>& elementI)
 const 
  626        if (momentumCouplingContext_().empty())
 
  629            bindCouplingContext_(domainI, elementI, momentumCouplingContext_()[0].fvGeometry.gridGeometry().elementMapper().index(elementI));
 
  632    void bindCouplingContext_(Dune::index_constant<freeFlowMomentumIndex> domainI,
 
  633                              const Element<freeFlowMomentumIndex>& elementI,
 
  634                              const std::size_t eIdx)
 const 
  636        if (momentumCouplingContext_().empty())
 
  639            fvGeometry.bind(elementI);
 
  648                prevElemVolVars.bindElement(elementI, fvGeometry, (*prevSol_)[
freeFlowMassIndex]);
 
  650            momentumCouplingContext_().emplace_back(MomentumCouplingContext{std::move(fvGeometry), std::move(curElemVolVars), std::move(prevElemVolVars), eIdx});
 
  652        else if (eIdx != momentumCouplingContext_()[0].eIdx)
 
  654            momentumCouplingContext_()[0].eIdx = eIdx;
 
  655            momentumCouplingContext_()[0].fvGeometry.bind(elementI);
 
  656            momentumCouplingContext_()[0].curElemVolVars.bind(elementI, momentumCouplingContext_()[0].fvGeometry, this->
curSol(
freeFlowMassIndex));
 
  659                momentumCouplingContext_()[0].prevElemVolVars.bindElement(elementI, momentumCouplingContext_()[0].fvGeometry, (*prevSol_)[
freeFlowMassIndex]);
 
  663    void bindCouplingContext_(Dune::index_constant<freeFlowMassIndex> domainI,
 
  664                              const Element<freeFlowMassIndex>& elementI)
 const 
  667        if (massAndEnergyCouplingContext_().empty())
 
  668            bindCouplingContext_(domainI, elementI, this->
problem(
freeFlowMassIndex).gridGeometry().elementMapper().index(elementI));
 
  670            bindCouplingContext_(domainI, elementI, massAndEnergyCouplingContext_()[0].fvGeometry.gridGeometry().elementMapper().index(elementI));
 
  673    void bindCouplingContext_(Dune::index_constant<freeFlowMassIndex> domainI,
 
  674                              const Element<freeFlowMassIndex>& elementI,
 
  675                              const std::size_t eIdx)
 const 
  677        if (massAndEnergyCouplingContext_().empty())
 
  680            auto fvGeometry = 
localView(gridGeometry);
 
  681            fvGeometry.bindElement(elementI);
 
  682            massAndEnergyCouplingContext_().emplace_back(std::move(fvGeometry), eIdx);
 
  684        else if (eIdx != massAndEnergyCouplingContext_()[0].eIdx)
 
  686            massAndEnergyCouplingContext_()[0].eIdx = eIdx;
 
  687            massAndEnergyCouplingContext_()[0].fvGeometry.bindElement(elementI);
 
  695    template<std::
size_t i>
 
  696    const GridVariables<i>& gridVars_(Dune::index_constant<i> domainIdx)
 const 
  698        if (std::get<i>(gridVariables_))
 
  699            return *std::get<i>(gridVariables_);
 
  701            DUNE_THROW(Dune::InvalidStateException, 
"The gridVariables pointer was not set. Use setGridVariables() before calling this function");
 
  708    template<std::
size_t i>
 
  709    GridVariables<i>& gridVars_(Dune::index_constant<i> domainIdx)
 
  711        if (std::get<i>(gridVariables_))
 
  712            return *std::get<i>(gridVariables_);
 
  714            DUNE_THROW(Dune::InvalidStateException, 
"The gridVariables pointer was not set. Use setGridVariables() before calling this function");
 
  718    void computeCouplingStencils_()
 
  722        auto momentumFvGeometry = 
localView(momentumGridGeometry);
 
  723        auto massFvGeometry = 
localView(massGridGeometry);
 
  725        massAndEnergyToMomentumStencils_.clear();
 
  726        massAndEnergyToMomentumStencils_.resize(massGridGeometry.gridView().size(0));
 
  728        momentumToMassAndEnergyStencils_.clear();
 
  729        momentumToMassAndEnergyStencils_.resize(momentumGridGeometry.gridView().size(0));
 
  731        assert(massAndEnergyToMomentumStencils_.size() == momentumToMassAndEnergyStencils_.size());
 
  733        for (
const auto& element : elements(momentumGridGeometry.gridView()))
 
  735            momentumFvGeometry.bindElement(element);
 
  736            massFvGeometry.bindElement(element);
 
  737            const auto eIdx = momentumFvGeometry.elementIndex();
 
  739            for (
const auto& scv : scvs(momentumFvGeometry))
 
  740                massAndEnergyToMomentumStencils_[eIdx].push_back(scv.dofIndex());
 
  742            for (
const auto& scv : scvs(massFvGeometry))
 
  743                momentumToMassAndEnergyStencils_[eIdx].push_back(scv.dofIndex());
 
  747    CouplingStencilType emptyStencil_;
 
  748    std::vector<CouplingStencilType> momentumToMassAndEnergyStencils_;
 
  749    std::vector<CouplingStencilType> massAndEnergyToMomentumStencils_;
 
  751    std::vector<MomentumCouplingContext>& momentumCouplingContext_()
 const 
  752    { 
return momentumCouplingContextImpl_; }
 
  754    std::vector<MassAndEnergyCouplingContext>& massAndEnergyCouplingContext_()
 const 
  755    { 
return massAndEnergyCouplingContextImpl_; }
 
  757    mutable std::vector<MassAndEnergyCouplingContext> massAndEnergyCouplingContextImpl_;
 
  758    mutable std::vector<MomentumCouplingContext> momentumCouplingContextImpl_;
 
  761    GridVariablesTuple gridVariables_;
 
  763    const SolutionVector* prevSol_;
 
  766    std::deque<std::vector<ElementSeed<freeFlowMomentumIndex>>> elementSets_;
 
 
  772template<class Traits, class DiscretizationMethod = typename Detail::MomentumDiscretizationMethod<Traits>::type>
 
  776template<
class Traits, 
class D>
 
  778{ 
using type = std::false_type; };
 
 
The interface of the coupling manager for free flow systems.
Definition couplingmanager_cvfe.hh:50
VelocityVector faceVelocity(const Element< freeFlowMassIndex > &element, const SubControlVolumeFace< freeFlowMassIndex > &scvf) const
Returns the velocity at a given sub control volume face.
Definition couplingmanager_cvfe.hh:397
void init(std::shared_ptr< Problem< freeFlowMomentumIndex > > momentumProblem, std::shared_ptr< Problem< freeFlowMassIndex > > massProblem, GridVariablesTuple &&gridVariables, const SolutionVector &curSol)
Methods to be accessed by main.
Definition couplingmanager_cvfe.hh:121
void init(std::shared_ptr< Problem< freeFlowMomentumIndex > > momentumProblem, std::shared_ptr< Problem< freeFlowMassIndex > > massProblem, GridVariablesTuple &&gridVariables, const typename ParentType::SolutionVectorStorage &curSol)
use as binary coupling manager in multi model context
Definition couplingmanager_cvfe.hh:149
Scalar effectiveViscosity(const Element< freeFlowMomentumIndex > &element, const FVElementGeometry< freeFlowMomentumIndex > &fvGeometry, const SubControlVolume< freeFlowMomentumIndex > &scv, const bool considerPreviousTimeStep=false) const
Returns the effective viscosity at a given sub control volume.
Definition couplingmanager_cvfe.hh:350
static constexpr auto freeFlowMomentumIndex
Definition couplingmanager_cvfe.hh:53
static constexpr auto pressureIdx
Definition couplingmanager_cvfe.hh:113
const CouplingStencilType & couplingStencil(Dune::index_constant< freeFlowMassIndex > domainI, const Element< freeFlowMassIndex > &elementI, Dune::index_constant< freeFlowMomentumIndex > domainJ) const
returns an iterable container of all indices of degrees of freedom of domain j that couple with / inf...
Definition couplingmanager_cvfe.hh:466
const CouplingStencilType & couplingStencil(Dune::index_constant< freeFlowMomentumIndex > domainI, const Element< freeFlowMomentumIndex > &elementI, Dune::index_constant< freeFlowMassIndex > domainJ) const
returns an iterable container of all indices of degrees of freedom of domain j that couple with / inf...
Definition couplingmanager_cvfe.hh:482
Scalar density(const Element< freeFlowMomentumIndex > &element, const FVElementGeometry< freeFlowMomentumIndex > &fvGeometry, const SubControlVolumeFace< freeFlowMomentumIndex > &scvf, const bool considerPreviousTimeStep=false) const
Returns the density at a given sub control volume face.
Definition couplingmanager_cvfe.hh:206
Scalar effectiveViscosity(const Element< freeFlowMomentumIndex > &element, const FVElementGeometry< freeFlowMomentumIndex > &fvGeometry, const SubControlVolumeFace< freeFlowMomentumIndex > &scvf, const bool considerPreviousTimeStep=false) const
Returns the effective viscosity at a given sub control volume face.
Definition couplingmanager_cvfe.hh:303
VelocityVector elementVelocity(const FVElementGeometry< freeFlowMassIndex > &fvGeometry) const
Returns the velocity at the element center.
Definition couplingmanager_cvfe.hh:423
Scalar density(const Element< freeFlowMomentumIndex > &element, const FVElementGeometry< freeFlowMomentumIndex > &fvGeometry, const SubControlVolume< freeFlowMomentumIndex > &scv, const bool considerPreviousTimeStep=false) const
Returns the density at a given sub control volume.
Definition couplingmanager_cvfe.hh:255
Scalar pressure(const Element< freeFlowMomentumIndex > &element, const FVElementGeometry< freeFlowMomentumIndex > &fvGeometry, const SubControlVolume< freeFlowMomentumIndex > &scv, const bool considerPreviousTimeStep=false) const
Returns the pressure at a given sub control volume.
Definition couplingmanager_cvfe.hh:190
void assembleMultithreaded(Dune::index_constant< i > domainId, AssembleElementFunc &&assembleElement) const
Execute assembly kernel in parallel.
Definition couplingmanager_cvfe.hh:601
Scalar pressure(const Element< freeFlowMomentumIndex > &element, const FVElementGeometry< freeFlowMomentumIndex > &fvGeometry, const SubControlVolumeFace< freeFlowMomentumIndex > &scvf, const bool considerPreviousTimeStep=false) const
Returns the pressure at a given sub control volume face.
Definition couplingmanager_cvfe.hh:174
void init(std::shared_ptr< Problem< freeFlowMomentumIndex > > momentumProblem, std::shared_ptr< Problem< freeFlowMassIndex > > massProblem, GridVariablesTuple &&gridVariables, const SolutionVector &curSol, const SolutionVector &prevSol)
use as regular coupling manager in a transient setting
Definition couplingmanager_cvfe.hh:137
void computeColorsForAssembly()
Compute colors for multithreaded assembly.
Definition couplingmanager_cvfe.hh:579
static constexpr auto freeFlowMassIndex
Definition couplingmanager_cvfe.hh:54
const CouplingStencilType & couplingStencil(Dune::index_constant< freeFlowMomentumIndex > domainI, const Element< freeFlowMomentumIndex > &elementI, const SubControlVolume< freeFlowMomentumIndex > &scvI, Dune::index_constant< j > domainJ) const
The coupling stencil of domain I, i.e. which domain J DOFs the given domain I element's residual depe...
Definition couplingmanager_cvfe.hh:446
typename ParentType::SolutionVectorStorage SolutionVectorStorage
Definition couplingmanager_cvfe.hh:58
void attachSolution(const SolutionVectorStorage &curSol)
Attach a solution vector stored outside of this class.
Definition multidomain/couplingmanager.hh:310
void setSubProblems(const std::tuple< std::shared_ptr< SubProblems >... > &problems)
set the pointers to the sub problems
Definition multidomain/couplingmanager.hh:275
const Problem< i > & problem(Dune::index_constant< i > domainIdx) const
Return a reference to the sub problem.
Definition multidomain/couplingmanager.hh:297
SubSolutionVector< i > & curSol(Dune::index_constant< i > domainIdx)
the solution vector of the subproblem
Definition multidomain/couplingmanager.hh:326
void updateSolution(const SolutionVector &curSol)
Updates the entire solution vector, e.g. before assembly or after grid adaption Overload might want t...
Definition multidomain/couplingmanager.hh:207
CouplingManager()
Default constructor.
Definition multidomain/couplingmanager.hh:70
typename Traits::template TupleOfSharedPtr< SubSolutionVector > SolutionVectorStorage
the type in which the solution vector is stored in the manager
Definition multidomain/couplingmanager.hh:59
Coloring schemes for shared-memory-parallel assembly.
Defines all properties used in Dumux.
Element solution classes and factory functions.
free functions for the evaluation of primary variables inside elements.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
PrimaryVariables evalSolution(const Element &element, const typename Element::Geometry &geometry, const typename FVElementGeometry::GridGeometry &gridGeometry, const CVFEElementSolution< FVElementGeometry, PrimaryVariables > &elemSol, const typename Element::Geometry::GlobalCoordinate &globalPos, bool ignoreState=false)
Interpolates a given box element solution at a given global position. Uses the finite element cache o...
Definition evalsolution.hh:152
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
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)
updates all data and variables that are necessary to evaluate the residual of the element of domain i...
Definition couplingmanager_cvfe.hh:517
void parallelFor(const std::size_t count, const FunctorType &functor)
A parallel for loop (multithreading)
Definition parallel_for.hh:160
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
The available discretization methods in Dumux.
The interface of the coupling manager for multi domain problems.
A linear system assembler (residual and Jacobian) for finite volume schemes with multiple domains.
Distance implementation details.
Definition cvfelocalresidual.hh:25
constexpr FCDiamond fcdiamond
Definition method.hh:152
constexpr CCTpfa cctpfa
Definition method.hh:145
constexpr Box box
Definition method.hh:147
auto computeColoring(const GridGeometry &gg, int verbosity=1)
Compute iterable lists of element seeds partitioned by color.
Definition coloring.hh:239
Parallel for loop (multithreading)
Type trait that is specialized for coupling manager supporting multithreaded assembly.
Definition multistagemultidomainfvassembler.hh:78
std::false_type type
Definition couplingmanager_cvfe.hh:778
Definition couplingmanager_cvfe.hh:773