From c6873cc44cd427f20f7317a79d2e3f132d7c92cc Mon Sep 17 00:00:00 2001
From: Kilian Weishaupt <kilian.weishaupt@iws.uni-stuttgart.de>
Date: Thu, 15 Nov 2018 10:51:11 +0100
Subject: [PATCH] [staggeredFF] Add hasBoundaryScvf() to fvGeometry

---
 .../freeflow/elementvolumevariables.hh        | 40 ++++++++++---------
 .../staggered/fvelementgeometry.hh            | 11 +++++
 .../staggered/fvgridgeometry.hh               |  8 ++++
 3 files changed, 41 insertions(+), 18 deletions(-)

diff --git a/dumux/discretization/staggered/freeflow/elementvolumevariables.hh b/dumux/discretization/staggered/freeflow/elementvolumevariables.hh
index 420b66fed7..70106b9fca 100644
--- a/dumux/discretization/staggered/freeflow/elementvolumevariables.hh
+++ b/dumux/discretization/staggered/freeflow/elementvolumevariables.hh
@@ -91,7 +91,8 @@ public:
               const SolutionVector& sol)
     {
         // the last parameter {} is needed for the PassKey pattern which restricts access to the ElementVolVars class
-        gridVolVars().updateBoundary_(element, fvGeometry, sol, {});
+        if (fvGeometry.hasBoundaryScvf())
+            gridVolVars().updateBoundary_(element, fvGeometry, sol, {});
     }
 
     //! function to prepare the vol vars within the element
@@ -181,24 +182,27 @@ 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())
-                continue;
-
-            volumeVariables_.resize(localIdx+1);
-            volVarIndices_.resize(localIdx+1);
-
-            auto boundaryPriVars = gridVolVars().getBoundaryPriVars(problem, sol, element, scvf);
-            auto elemSol = elementSolution<FVElementGeometry>(std::move(boundaryPriVars));
-            volumeVariables_[localIdx].update(elemSol,
-                                              problem,
-                                              element,
-                                              scvI);
-            volVarIndices_[localIdx] = scvf.outsideScvIdx();
-            ++localIdx;
+            // 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;
+
+                volumeVariables_.resize(localIdx+1);
+                volVarIndices_.resize(localIdx+1);
+
+                auto boundaryPriVars = gridVolVars().getBoundaryPriVars(problem, sol, element, scvf);
+                auto elemSol = elementSolution<FVElementGeometry>(std::move(boundaryPriVars));
+                volumeVariables_[localIdx].update(elemSol,
+                                                  problem,
+                                                  element,
+                                                  scvI);
+                volVarIndices_[localIdx] = scvf.outsideScvIdx();
+                ++localIdx;
+            }
         }
     }
 
diff --git a/dumux/discretization/staggered/fvelementgeometry.hh b/dumux/discretization/staggered/fvelementgeometry.hh
index 75d83d8ce3..e34dea10ff 100644
--- a/dumux/discretization/staggered/fvelementgeometry.hh
+++ b/dumux/discretization/staggered/fvelementgeometry.hh
@@ -208,6 +208,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:
 
     //! create scvs and scvfs of the bound element
@@ -238,6 +242,9 @@ private:
                                     geometryHelper);
                 scvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
                 scvfCounter++;
+
+                if (intersection.boundary())
+                    hasBoundaryScvf_ = true;
             }
         }
     }
@@ -301,6 +308,8 @@ private:
         neighborScvfIndices_.clear();
         neighborScvs_.clear();
         neighborScvfs_.clear();
+
+        hasBoundaryScvf_ = false;
     }
 
     const Element* elementPtr_; //!< the element to which this fvgeometry is bound
@@ -318,6 +327,8 @@ private:
 
     std::vector<IndexType> neighborScvfIndices_;
     std::vector<SubControlVolumeFace> neighborScvfs_;
+
+    bool hasBoundaryScvf_ = false;
 };
 
 
diff --git a/dumux/discretization/staggered/fvgridgeometry.hh b/dumux/discretization/staggered/fvgridgeometry.hh
index c22674ed6f..8535d30fd5 100644
--- a/dumux/discretization/staggered/fvgridgeometry.hh
+++ b/dumux/discretization/staggered/fvgridgeometry.hh
@@ -281,6 +281,7 @@ public:
         scvfs_.reserve(numScvf);
         scvfIndicesOfScv_.resize(numScvs);
         localToGlobalScvfIndices_.resize(numScvs);
+        hasBoundaryScvf_.resize(numScvs, false);
 
         // Build the scvs and scv faces
         IndexType scvfIdx = 0;
@@ -328,6 +329,8 @@ public:
                                         geometryHelper);
                     localToGlobalScvfIndices_[eIdx][localFaceIndex] = scvfIdx;
                     scvfsIndexSet.push_back(scvfIdx++);
+
+                    hasBoundaryScvf_[eIdx] = true;
                 }
             }
 
@@ -398,6 +401,10 @@ public:
         return FaceFVGridGeometry<ThisType>(this);
     }
 
+    //! Returns whether one of the geometry's scvfs lies on a boundary
+    bool hasBoundaryScvf(IndexType eIdx) const
+    { return hasBoundaryScvf_[eIdx]; }
+
 private:
 
     // mappers
@@ -409,6 +416,7 @@ private:
     std::vector<std::vector<IndexType>> scvfIndicesOfScv_;
     std::vector<std::vector<IndexType>> localToGlobalScvfIndices_;
     IndexType numBoundaryScvf_;
+    std::vector<bool> hasBoundaryScvf_;
 };
 
 /*!
-- 
GitLab