12#ifndef DUMUX_MUTLIDOMAIN_COUPLING_JACOBIAN_PATTERN_HH 
   13#define DUMUX_MUTLIDOMAIN_COUPLING_JACOBIAN_PATTERN_HH 
   16#include <dune/common/indices.hh> 
   17#include <dune/istl/matrixindexset.hh> 
   27template<
bool isImplicit, 
class CouplingManager, 
class GridGeometryI, 
class GridGeometryJ, std::size_t i, std::size_t j,
 
   31                                                Dune::index_constant<i> domainI,
 
   32                                                const GridGeometryI& gridGeometryI,
 
   33                                                Dune::index_constant<j> domainJ,
 
   34                                                const GridGeometryJ& gridGeometryJ)
 
   36    const auto numDofsI = gridGeometryI.numDofs();
 
   37    const auto numDofsJ = gridGeometryJ.numDofs();
 
   38    Dune::MatrixIndexSet pattern;
 
   39    pattern.resize(numDofsI, numDofsJ);
 
   44        for (
const auto& elementI : elements(gridGeometryI.gridView()))
 
   46            const auto& stencil = couplingManager.
couplingStencil(domainI, elementI, domainJ);
 
   47            const auto globalI = gridGeometryI.elementMapper().index(elementI);
 
   48            for (
const auto globalJ : stencil)
 
   49                pattern.add(globalI, globalJ);
 
 
   65template<
bool isImplicit, 
class CouplingManager, 
class GridGeometryI, 
class GridGeometryJ, std::size_t i, std::size_t j,
 
   67                                    GridGeometryI::isCellCenter()), 
int> = 0>
 
   69                                                Dune::index_constant<i> domainI,
 
   70                                                const GridGeometryI& gridGeometryI,
 
   71                                                Dune::index_constant<j> domainJ,
 
   72                                                const GridGeometryJ& gridGeometryJ)
 
   74    Dune::MatrixIndexSet pattern(gridGeometryI.numDofs(), gridGeometryJ.numDofs());
 
   76    for (
const auto& elementI : elements(gridGeometryI.gridView()))
 
   78        const auto ccGlobalI = gridGeometryI.elementMapper().index(elementI);
 
   79        for (
auto&& faceGlobalJ : couplingManager.couplingStencil(domainI, elementI, domainJ))
 
   80                pattern.add(ccGlobalI, faceGlobalJ);
 
   91template<
bool isImplicit, 
class CouplingManager, 
class GridGeometryI, 
class GridGeometryJ, std::size_t i, std::size_t j,
 
   93                                    GridGeometryI::isFace()), 
int> = 0>
 
   95                                                Dune::index_constant<i> domainI,
 
   96                                                const GridGeometryI& gridGeometryI,
 
   97                                                Dune::index_constant<j> domainJ,
 
   98                                                const GridGeometryJ& gridGeometryJ)
 
  100    Dune::MatrixIndexSet pattern(gridGeometryI.numDofs(), gridGeometryJ.numDofs());
 
  102    auto fvGeometry = 
localView(gridGeometryI);
 
  103    for (
const auto& elementI : elements(gridGeometryI.gridView()))
 
  105        fvGeometry.bindElement(elementI);
 
  108        for (
auto&& scvf : scvfs(fvGeometry))
 
  110            const auto faceGlobalI = scvf.dofIndex();
 
  111            for (
auto&& globalJ : couplingManager.couplingStencil(domainI, scvf, domainJ))
 
  112                pattern.add(faceGlobalI, globalJ);
 
  124template<
bool isImplicit, 
class CouplingManager, 
class GridGeometryI, 
class GridGeometryJ, std::size_t i, std::size_t j,
 
  127                                                Dune::index_constant<i> domainI,
 
  128                                                const GridGeometryI& gridGeometryI,
 
  129                                                Dune::index_constant<j> domainJ,
 
  130                                                const GridGeometryJ& gridGeometryJ)
 
  132    Dune::MatrixIndexSet pattern(gridGeometryI.numDofs(), gridGeometryJ.numDofs());
 
  134    auto fvGeometry = 
localView(gridGeometryI);
 
  135    for (
const auto& elementI : elements(gridGeometryI.gridView()))
 
  137        fvGeometry.bindElement(elementI);
 
  138        for (
const auto& scv : scvs(fvGeometry))
 
  140            const auto globalI = scv.dofIndex();
 
  141            const auto& stencil = couplingManager.couplingStencil(domainI, elementI, scv, domainJ);
 
  142            for (
const auto globalJ : stencil)
 
  144                assert(globalJ < gridGeometryJ.numDofs());
 
  145                pattern.add(globalI, globalJ);
 
  147                if (gridGeometryI.isPeriodic())
 
  149                    if (gridGeometryI.dofOnPeriodicBoundary(globalI))
 
  151                        const auto globalIP = gridGeometryI.periodicallyMappedDof(globalI);
 
  153                        if (globalI > globalIP)
 
  154                            pattern.add(globalIP, globalJ);
 
  169template<
bool isImplicit, 
class CouplingManager, 
class GridGeometryI, 
class GridGeometryJ, std::size_t i, std::size_t j,
 
  170         typename std::enable_if_t<DiscretizationMethods::isCVFE<typename GridGeometryI::DiscretizationMethod>, 
int> = 0>
 
  172                                                Dune::index_constant<i> domainI,
 
  173                                                const GridGeometryI& gridGeometryI,
 
  174                                                Dune::index_constant<j> domainJ,
 
  175                                                const GridGeometryJ& gridGeometryJ)
 
  177    Dune::MatrixIndexSet pattern;
 
  180    if constexpr (isImplicit)
 
  182        pattern.resize(gridGeometryI.numDofs(),  gridGeometryJ.numDofs());
 
  183        auto fvGeometry = 
localView(gridGeometryI);
 
  184        for (
const auto& elementI : elements(gridGeometryI.gridView()))
 
  186            fvGeometry.bindElement(elementI);
 
  187            const auto& stencil = couplingManager.couplingStencil(domainI, elementI, domainJ);
 
  188            for (
const auto& scv : scvs(fvGeometry))
 
  190                for (
const auto globalJ : stencil)
 
  191                    pattern.add(scv.dofIndex(), globalJ);
 
The interface of the coupling manager for multi domain problems.
Definition multidomain/couplingmanager.hh:37
const CouplingStencilType< i, j > & couplingStencil(Dune::index_constant< i > domainI, const Element< i > &elementI, Dune::index_constant< j > domainJ) const
returns an iterable container of all indices of degrees of freedom of domain j that couple with / inf...
Definition multidomain/couplingmanager.hh:102
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
Dune::MatrixIndexSet getCouplingJacobianPattern(const CouplingManager &couplingManager, Dune::index_constant< i > domainI, const GridGeometryI &gridGeometryI, Dune::index_constant< j > domainJ, const GridGeometryJ &gridGeometryJ)
Helper function to generate coupling Jacobian pattern (off-diagonal blocks) for cell-centered schemes...
Definition couplingjacobianpattern.hh:30
The available discretization methods in Dumux.
constexpr CCMpfa ccmpfa
Definition method.hh:146
constexpr CCTpfa cctpfa
Definition method.hh:145
constexpr Staggered staggered
Definition method.hh:149
constexpr FCStaggered fcstaggered
Definition method.hh:151