12#ifndef DUMUX_CVFE_ELEMENT_BOUNDARY_TYPES_HH 
   13#define DUMUX_CVFE_ELEMENT_BOUNDARY_TYPES_HH 
   40    template<
class Problem, 
class FVElementGeometry>
 
   42                const typename FVElementGeometry::Element& element,
 
   43                const FVElementGeometry& fvGeometry)
 
   45        bcTypes_.resize(fvGeometry.numScv());
 
   47        hasDirichlet_ = 
false;
 
   50        for (
const auto& scv : scvs(fvGeometry))
 
   52            const auto scvIdxLocal = scv.localDofIndex();
 
   53            bcTypes_[scvIdxLocal].reset();
 
   55            if (fvGeometry.gridGeometry().dofOnBoundary(scv.dofIndex()))
 
   57                bcTypes_[scvIdxLocal] = problem.boundaryTypes(element, scv);
 
   58                hasDirichlet_ = hasDirichlet_ || bcTypes_[scvIdxLocal].hasDirichlet();
 
   59                hasNeumann_ = hasNeumann_ || bcTypes_[scvIdxLocal].hasNeumann();
 
 
   69    { 
return hasDirichlet_; }
 
 
   76    { 
return hasNeumann_; }
 
 
   83    template<
class FVElementGeometry>
 
   84    const BoundaryTypes& 
get(
const FVElementGeometry& fvGeometry, 
const typename FVElementGeometry::SubControlVolumeFace& scvf)
 const 
   86        assert(scvf.boundary());
 
   87        const auto localDofIdx = fvGeometry.scv(scvf.insideScvIdx()).localDofIndex();
 
   88        assert(localDofIdx < bcTypes_.size());
 
   89        return bcTypes_[localDofIdx];
 
 
   97    template<
class FVElementGeometry>
 
   98    const BoundaryTypes& 
get(
const FVElementGeometry&, 
const typename FVElementGeometry::SubControlVolume& scv)
 const 
  100        const auto localDofIdx = scv.localDofIndex();
 
  101        assert(localDofIdx < bcTypes_.size());
 
  102        return bcTypes_[localDofIdx];
 
 
  106    std::vector<BoundaryTypes> bcTypes_;
 
  107    bool hasDirichlet_ = 
false;
 
  108    bool hasNeumann_ = 
false;
 
 
This class stores an array of BoundaryTypes objects.
Definition cvfe/elementboundarytypes.hh:28
BoundaryTypes BoundaryTypes
Definition cvfe/elementboundarytypes.hh:30
bool hasDirichlet() const
Returns whether the element has a vertex which contains a Dirichlet value.
Definition cvfe/elementboundarytypes.hh:68
const BoundaryTypes & get(const FVElementGeometry &, const typename FVElementGeometry::SubControlVolume &scv) const
Definition cvfe/elementboundarytypes.hh:98
bool hasNeumann() const
Returns whether the element potentially features a Neumann boundary segment.
Definition cvfe/elementboundarytypes.hh:75
const BoundaryTypes & get(const FVElementGeometry &fvGeometry, const typename FVElementGeometry::SubControlVolumeFace &scvf) const
Definition cvfe/elementboundarytypes.hh:84
void update(const Problem &problem, const typename FVElementGeometry::Element &element, const FVElementGeometry &fvGeometry)
Update the boundary types for all vertices of an element.
Definition cvfe/elementboundarytypes.hh:41
The available discretization methods in Dumux.