diff --git a/dumux/discretization/cellcentered/gridvolumevariables.hh b/dumux/discretization/cellcentered/gridvolumevariables.hh
index b9b250ffd78686829fa1a33cfd07ec6a5724831d..21ad3d9edb5b61404b83ea8dbebdb64240af2861 100644
--- a/dumux/discretization/cellcentered/gridvolumevariables.hh
+++ b/dumux/discretization/cellcentered/gridvolumevariables.hh
@@ -80,25 +80,28 @@ public:
volumeVariables_[scv.dofIndex()].update(elemSol, problem(), element, scv);
}
- // handle the boundary volume variables
- for (auto&& scvf : scvfs(fvGeometry))
+ if (fvGeometry.hasBoundaryScvf())
{
- // if we are not on a boundary, skip the rest
- if (!scvf.boundary())
- continue;
-
- // check if boundary is a pure dirichlet boundary
- const auto bcTypes = problem().boundaryTypes(element, scvf);
- if (bcTypes.hasOnlyDirichlet())
+ // handle the boundary volume variables
+ for (auto&& scvf : scvfs(fvGeometry))
{
- const auto insideScvIdx = scvf.insideScvIdx();
- const auto& insideScv = fvGeometry.scv(insideScvIdx);
- const auto dirichletPriVars = elementSolution<typename FVGridGeometry::LocalView>(problem().dirichlet(element, scvf));
-
- volumeVariables_[scvf.outsideScvIdx()].update(dirichletPriVars,
- problem(),
- element,
- insideScv);
+ // if we are not on a boundary, skip the rest
+ if (!scvf.boundary())
+ continue;
+
+ // check if boundary is a pure dirichlet boundary
+ const auto bcTypes = problem().boundaryTypes(element, scvf);
+ if (bcTypes.hasOnlyDirichlet())
+ {
+ const auto insideScvIdx = scvf.insideScvIdx();
+ const auto& insideScv = fvGeometry.scv(insideScvIdx);
+ const auto dirichletPriVars = elementSolution<typename FVGridGeometry::LocalView>(problem().dirichlet(element, scvf));
+
+ volumeVariables_[scvf.outsideScvIdx()].update(dirichletPriVars,
+ problem(),
+ element,
+ insideScv);
+ }
}
}
}
diff --git a/dumux/discretization/cellcentered/tpfa/elementvolumevariables.hh b/dumux/discretization/cellcentered/tpfa/elementvolumevariables.hh
index 9612e6c1ea50e0020e9f6b669b233a7e9d0d72ea..15f70efdbe7c0446035c3cf71e6a6456e0925f6e 100644
--- a/dumux/discretization/cellcentered/tpfa/elementvolumevariables.hh
+++ b/dumux/discretization/cellcentered/tpfa/elementvolumevariables.hh
@@ -152,28 +152,31 @@ public:
++localIdx;
}
- // Update boundary volume variables
- for (auto&& scvf : scvfs(fvGeometry))
+ if (fvGeometry.hasBoundaryScvf())
{
- // if we are not on a boundary, skip to the next scvf
- if (!scvf.boundary())
+ // Update boundary volume variables
+ for (auto&& scvf : scvfs(fvGeometry))
+ {
+ // if we are not on a boundary, skip to the next scvf
+ if (!scvf.boundary())
continue;
- // check if boundary is a pure dirichlet boundary
- const auto bcTypes = problem.boundaryTypes(element, scvf);
- if (bcTypes.hasOnlyDirichlet())
- {
- const auto dirichletPriVars = elementSolution<FVElementGeometry>(problem.dirichlet(element, scvf));
-
- volumeVariables_.resize(localIdx+1);
- volVarIndices_.resize(localIdx+1);
- volumeVariables_[localIdx].update(dirichletPriVars,
- problem,
- element,
- scvI);
- volVarIndices_[localIdx] = scvf.outsideScvIdx();
- ++localIdx;
- }
+ // check if boundary is a pure dirichlet boundary
+ const auto bcTypes = problem.boundaryTypes(element, scvf);
+ if (bcTypes.hasOnlyDirichlet())
+ {
+ const auto dirichletPriVars = elementSolution<FVElementGeometry>(problem.dirichlet(element, scvf));
+
+ volumeVariables_.resize(localIdx+1);
+ volVarIndices_.resize(localIdx+1);
+ volumeVariables_[localIdx].update(dirichletPriVars,
+ problem,
+ element,
+ scvI);
+ volVarIndices_[localIdx] = scvf.outsideScvIdx();
+ ++localIdx;
+ }
+ }
}
//! Check if user added additional DOF dependencies, i.e. the residual of DOF globalI depends
diff --git a/dumux/discretization/cellcentered/tpfa/fvelementgeometry.hh b/dumux/discretization/cellcentered/tpfa/fvelementgeometry.hh
index c1e4b4388b1de6558445f1bf0c03ac832179b029..3ccaf0ef026651cf5dd796658dd1d49838ab4bad 100644
--- a/dumux/discretization/cellcentered/tpfa/fvelementgeometry.hh
+++ b/dumux/discretization/cellcentered/tpfa/fvelementgeometry.hh
@@ -156,6 +156,10 @@ public:
const FVGridGeometry& fvGridGeometry() const
{ return *fvGridGeometryPtr_; }
+ //! Returns whether one of the geometry's scvfs lies on a boundary
+ bool hasBoundaryScvf() const
+ { return fvGridGeometry().hasBoundaryScvf(scvIndices_[0]); }
+
private:
const Element* elementPtr_;
@@ -351,6 +355,10 @@ public:
const FVGridGeometry& fvGridGeometry() const
{ return *fvGridGeometryPtr_; }
+ //! Returns whether one of the geometry's scvfs lies on a boundary
+ bool hasBoundaryScvf() const
+ { return hasBoundaryScvf_; }
+
private:
IndexType findFlippedScvfIndex_(IndexType insideScvIdx, IndexType globalOutsideScvIdx)
@@ -412,11 +420,15 @@ private:
scvIndices.resize(scvfNeighborVolVarIndices.size() + 1);
scvIndices[0] = eIdx;
std::copy(scvfNeighborVolVarIndices.begin(), scvfNeighborVolVarIndices.end(), scvIndices.begin()+1);
+
+ const bool onBoundary = intersection.boundary() && !intersection.neighbor();
+ hasBoundaryScvf_ = (hasBoundaryScvf_ || onBoundary);
+
scvfs_.emplace_back(intersection,
intersection.geometry(),
scvFaceIndices[scvfCounter],
scvIndices,
- intersection.boundary() && !intersection.neighbor());
+ onBoundary);
scvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
scvfCounter++;
@@ -536,6 +548,8 @@ private:
neighborScvs_.clear();
neighborScvfs_.clear();
flippedNeighborScvfIndices_.clear();
+
+ hasBoundaryScvf_ = false;
}
const Element* elementPtr_; //!< the element to which this fvgeometry is bound
@@ -555,6 +569,8 @@ private:
std::vector<IndexType> neighborScvfIndices_;
std::vector<SubControlVolumeFace> neighborScvfs_;
std::vector<std::vector<IndexType>> flippedNeighborScvfIndices_;
+
+ bool hasBoundaryScvf_ = false;
};
} // end namespace Dumux
diff --git a/dumux/discretization/cellcentered/tpfa/fvgridgeometry.hh b/dumux/discretization/cellcentered/tpfa/fvgridgeometry.hh
index be19deafec40f9c5f412c2b2c0021a12dfec0e4f..9d59a04246533e3b67f09611aa6e304795290d69 100644
--- a/dumux/discretization/cellcentered/tpfa/fvgridgeometry.hh
+++ b/dumux/discretization/cellcentered/tpfa/fvgridgeometry.hh
@@ -172,6 +172,7 @@ public:
scvs_.resize(numScvs);
scvfs_.reserve(numScvf);
scvfIndicesOfScv_.resize(numScvs);
+ hasBoundaryScvf_.resize(numScvs, false);
// Build the scvs and scv faces
IndexType scvfIdx = 0;
@@ -254,6 +255,8 @@ public:
ScvfGridIndexStorage({eIdx, this->gridView().size(0) + numBoundaryScvf_++}),
true);
scvfsIndexSet.push_back(scvfIdx++);
+
+ hasBoundaryScvf_[eIdx] = true;
}
}
@@ -314,6 +317,10 @@ public:
const ConnectivityMap &connectivityMap() const
{ return connectivityMap_; }
+ //! Returns whether one of the geometry's scvfs lies on a boundary
+ bool hasBoundaryScvf(IndexType eIdx) const
+ { return hasBoundaryScvf_[eIdx]; }
+
private:
// find the scvf that has insideScvIdx in its outsideScvIdx list and outsideScvIdx as its insideScvIdx
IndexType findFlippedScvfIndex_(IndexType insideScvIdx, IndexType outsideScvIdx)
@@ -338,6 +345,7 @@ private:
std::vector<SubControlVolumeFace> scvfs_;
std::vector<std::vector<IndexType>> scvfIndicesOfScv_;
IndexType numBoundaryScvf_;
+ std::vector<bool> hasBoundaryScvf_;
//! needed for embedded surface and network grids (dim < dimWorld)
std::vector<std::vector<IndexType>> flipScvfIndices_;