12#ifndef DUMUX_DISCRETIZATION_STAGGERED_FV_ELEMENT_GEOMETRY_HH 
   13#define DUMUX_DISCRETIZATION_STAGGERED_FV_ELEMENT_GEOMETRY_HH 
   33template<
class GG, 
bool enableGr
idGeometryCache>
 
   48    using GridView = 
typename GG::GridView;
 
   53    using Element = 
typename GridView::template Codim<0>::Entity;
 
   57    using ParentType::ParentType;
 
   61    template<
class CellCenterOrFaceFVGr
idGeometry>
 
   63    : ParentType(gridGeometry.actualGridGeometry()) {}
 
 
   66    using ParentType::scvf;
 
   69        return this->gridGeometry().scvf(eIdx, localScvfIdx);
 
 
 
   82class StaggeredFVElementGeometry<GG, false>
 
   84    using ThisType = StaggeredFVElementGeometry<GG, false>;
 
   85    using GridView = 
typename GG::GridView;
 
   90    using Element = 
typename GridView::template Codim<0>::Entity;
 
   92    using SubControlVolume = 
typename GG::SubControlVolume;
 
   94    using SubControlVolumeFace = 
typename GG::SubControlVolumeFace;
 
   96    using GridGeometry = GG;
 
  100    template<
class CellCenterOrFaceFVGr
idGeometry>
 
  101    StaggeredFVElementGeometry(
const CellCenterOrFaceFVGridGeometry& gridGeometry)
 
  102    : gridGeometryPtr_(&gridGeometry.actualGridGeometry()) {}
 
  105    StaggeredFVElementGeometry(
const GridGeometry& gridGeometry)
 
  106    : gridGeometryPtr_(&gridGeometry) {}
 
  109    const SubControlVolumeFace& scvf(GridIndexType eIdx, LocalIndexType localScvfIdx)
 const 
  111        return scvf(this->gridGeometry().localToGlobalScvfIndex(eIdx, localScvfIdx));
 
  116    const SubControlVolume& scv(GridIndexType scvIdx)
 const 
  118        if (scvIdx == scvIndices_[0])
 
  121            return neighborScvs_[findLocalIndex_(scvIdx, neighborScvIndices_)];
 
  126    const SubControlVolumeFace& scvf(GridIndexType scvfIdx)
 const 
  128        auto it = std::find(scvfIndices_.begin(), scvfIndices_.end(), scvfIdx);
 
  129        if (it != scvfIndices_.end())
 
  130            return scvfs_[std::distance(scvfIndices_.begin(), it)];
 
  132            return neighborScvfs_[findLocalIndex_(scvfIdx, neighborScvfIndices_)];
 
  140    friend inline Dune::IteratorRange<typename std::array<SubControlVolume, 1>::const_iterator>
 
  141    scvs(
const ThisType& g)
 
  143        using IteratorType = 
typename std::array<SubControlVolume, 1>::const_iterator;
 
  144        return Dune::IteratorRange<IteratorType>(g.scvs_.begin(), g.scvs_.end());
 
  152    friend inline Dune::IteratorRange<typename std::vector<SubControlVolumeFace>::const_iterator>
 
  153    scvfs(
const ThisType& g)
 
  155        using IteratorType = 
typename std::vector<SubControlVolumeFace>::const_iterator;
 
  156        return Dune::IteratorRange<IteratorType>(g.scvfs_.begin(), g.scvfs_.end());
 
  160    std::size_t numScv()
 const 
  161    { 
return scvs_.size(); }
 
  164    std::size_t numScvf()
 const 
  165    { 
return scvfs_.size(); }
 
  172    StaggeredFVElementGeometry bind(
const Element& element) &&
 
  174        this->bind_(element);
 
  175        return std::move(*
this);
 
  179    void bind(
const Element& element) &
 
  180    { this->bind_(element); }
 
  187    StaggeredFVElementGeometry bindElement(
const Element& element) &&
 
  189        this->bindElement_(element);
 
  190        return std::move(*
this);
 
  194    void bindElement(
const Element& element) &
 
  195    { this->bindElement_(element); }
 
  199    { 
return static_cast<bool>(element_); }
 
  203    { 
return *element_; }
 
  206    const GridGeometry& gridGeometry()
 const 
  207    { 
return *gridGeometryPtr_; }
 
  210    bool hasBoundaryScvf()
 const 
  211    { 
return hasBoundaryScvf_; }
 
  214    typename SubControlVolume::Traits::Geometry geometry(
const SubControlVolume& scv)
 const 
  215    { 
return gridGeometryPtr_->element(scv.dofIndex()).geometry(); }
 
  218    typename SubControlVolumeFace::Traits::Geometry geometry(
const SubControlVolumeFace& scvf)
 const 
  220        const auto insideElementIndex = scvf.insideScvIdx();
 
  221        const auto elementGeometry = (insideElementIndex != scvIndices_[0]) ?
 
  222            element_->geometry() :
 
  223            gridGeometryPtr_->element(insideElementIndex).geometry();
 
  224        const auto center = elementGeometry.center();
 
  227        auto lowerLeft = elementGeometry.corner(0);
 
  228        auto upperRight = elementGeometry.corner(elementGeometry.corners()-1);
 
  234        auto inPlaneAxes = std::move(std::bitset<SubControlVolumeFace::Traits::dimWorld>{}.set());
 
  237        return {lowerLeft, upperRight, inPlaneAxes};
 
  242    void bindElement_(
const Element& element)
 
  246        scvfs_.reserve(
element.subEntities(1));
 
  247        scvfIndices_.reserve(
element.subEntities(1));
 
  248        makeElementGeometries_();
 
  253    void bind_(
const Element& element)
 
  255        bindElement_(element);
 
  257        neighborScvs_.reserve(
element.subEntities(1));
 
  258        neighborScvfIndices_.reserve(
element.subEntities(1));
 
  259        neighborScvfs_.reserve(
element.subEntities(1));
 
  261        std::vector<GridIndexType> handledNeighbors;
 
  262        handledNeighbors.reserve(
element.subEntities(1));
 
  263        for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
 
  265            if (intersection.neighbor())
 
  267                const auto outside = intersection.outside();
 
  268                const auto outsideIdx = gridGeometry().elementMapper().index(outside);
 
  271                if ( std::find(handledNeighbors.begin(), handledNeighbors.end(), outsideIdx) == handledNeighbors.end() )
 
  273                    makeNeighborGeometries_(outside, outsideIdx);
 
  274                    handledNeighbors.push_back(outsideIdx);
 
  281    void makeElementGeometries_()
 
  283        const auto& 
element = *element_;
 
  284        const auto eIdx = gridGeometry().elementMapper().index(element);
 
  285        scvs_[0] = SubControlVolume(
element.geometry(), eIdx);
 
  286        scvIndices_[0] = eIdx;
 
  288        const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
 
  289        const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
 
  291        typename GridGeometry::GeometryHelper geometryHelper(element, gridGeometry().gridView());
 
  294        for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
 
  296            const auto& scvfNeighborVolVarIndex = neighborVolVarIndices[scvfCounter];
 
  298            if (intersection.neighbor() || intersection.boundary())
 
  300                geometryHelper.updateLocalFace(gridGeometry().intersectionMapper(), intersection);
 
  301                std::vector<GridIndexType> scvIndices{eIdx, scvfNeighborVolVarIndex};
 
  302                scvfs_.emplace_back(intersection,
 
  303                                    intersection.geometry(),
 
  304                                    scvFaceIndices[scvfCounter],
 
  307                scvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
 
  310                if (intersection.boundary())
 
  311                    hasBoundaryScvf_ = 
true;
 
  317    void makeNeighborGeometries_(
const Element& element, 
const GridIndexType eIdx)
 
  322        neighborScvs_.emplace_back(
element.geometry(), eIdx);
 
  323        neighborScvIndices_.push_back(eIdx);
 
  325        const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
 
  326        const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
 
  328        typename GridGeometry::GeometryHelper geometryHelper(element, gridGeometry().gridView());
 
  331        for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
 
  333            const auto& scvfNeighborVolVarIndex = neighborVolVarIndices[scvfCounter];
 
  334            geometryHelper.updateLocalFace(gridGeometry().intersectionMapper(), intersection);
 
  336            if (intersection.neighbor())
 
  339                if (intersection.outside() == *element_)
 
  341                    std::vector<GridIndexType> scvIndices{eIdx, scvfNeighborVolVarIndex};
 
  342                    neighborScvfs_.emplace_back(intersection,
 
  343                                                intersection.geometry(),
 
  344                                                scvFaceIndices[scvfCounter],
 
  348                    neighborScvfIndices_.push_back(scvFaceIndices[scvfCounter]);
 
  352            else if (intersection.boundary())
 
  357    const LocalIndexType findLocalIndex_(
const GridIndexType idx,
 
  358                                         const std::vector<GridIndexType>& indices)
 const 
  360        auto it = std::find(indices.begin(), indices.end(), idx);
 
  361        assert(it != indices.end() && 
"Could not find the scv/scvf! Make sure to properly bind this class!");
 
  362        return std::distance(indices.begin(), it);
 
  368        scvfIndices_.clear();
 
  371        neighborScvIndices_.clear();
 
  372        neighborScvfIndices_.clear();
 
  373        neighborScvs_.clear();
 
  374        neighborScvfs_.clear();
 
  376        hasBoundaryScvf_ = 
false;
 
  379    std::optional<Element> element_; 
 
  380    const GridGeometry* gridGeometryPtr_;  
 
  383    std::array<GridIndexType, 1> scvIndices_;
 
  384    std::array<SubControlVolume, 1> scvs_;
 
  386    std::vector<GridIndexType> scvfIndices_;
 
  387    std::vector<SubControlVolumeFace> scvfs_;
 
  389    std::vector<GridIndexType> neighborScvIndices_;
 
  390    std::vector<SubControlVolume> neighborScvs_;
 
  392    std::vector<GridIndexType> neighborScvfIndices_;
 
  393    std::vector<SubControlVolumeFace> neighborScvfs_;
 
  395    bool hasBoundaryScvf_ = 
false;
 
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models This builds up th...
Definition discretization/cellcentered/tpfa/fvelementgeometry.hh:53
const SubControlVolumeFace & scvf(GridIndexType eIdx, LocalIndexType localScvfIdx) const
Definition discretization/staggered/fvelementgeometry.hh:67
StaggeredFVElementGeometry(const CellCenterOrFaceFVGridGeometry &gridGeometry)
Definition discretization/staggered/fvelementgeometry.hh:62
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition discretization/staggered/fvelementgeometry.hh:55
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition discretization/staggered/fvelementgeometry.hh:53
Stencil-local finite volume geometry (scvs and scvfs) for staggered models This builds up the sub con...
Definition discretization/staggered/fvelementgeometry.hh:34
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models This builds up th...
static std::size_t normalAxis(const Vector &v)
Returns the normal axis index of a unit vector (0 = x, 1 = y, 2 = z)
Definition normalaxis.hh:26
Corners::value_type center(const Corners &corners)
The center of a given list of corners.
Definition center.hh:24
Base class for the finite volume geometry vector for face-centered staggered models This builds up th...
typename GridView::IndexSet::IndexType GridIndex
Definition indextraits.hh:27
unsigned int LocalIndex
Definition indextraits.hh:28