diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.H b/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.H
index 28057107c48f58e968def9a62b7f445e90f265ca..23dd4193495c1a94400ae2eb4121210d8ff0347b 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.H
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.H
@@ -822,6 +822,7 @@ public:
             const labelList& cellPoints
             (
                 const label celli,
+                labelHashSet&,
                 DynamicList<label>&
             ) const;
 
@@ -876,6 +877,7 @@ public:
             const labelList& cellEdges
             (
                 const label celli,
+                labelHashSet&,
                 DynamicList<label>&
             ) const;
 
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCellPoints.C b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCellPoints.C
index f7c0e4699a0b341442972dffb793ab54b8db7e0c..0f96e519ccae60ebd73040d7d6f9e4fccdc7fdb3 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCellPoints.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCellPoints.C
@@ -58,6 +58,7 @@ const Foam::labelListList& Foam::primitiveMesh::cellPoints() const
 const Foam::labelList& Foam::primitiveMesh::cellPoints
 (
     const label celli,
+    labelHashSet& set,
     DynamicList<label>& storage
 ) const
 {
@@ -70,27 +71,22 @@ const Foam::labelList& Foam::primitiveMesh::cellPoints
         const faceList& fcs = faces();
         const labelList& cFaces = cells()[celli];
 
-        labelSet_.clear();
+        set.clear();
 
-        forAll(cFaces, i)
+        for (const label facei : cFaces)
         {
-            const labelList& f = fcs[cFaces[i]];
-
-            forAll(f, fp)
-            {
-                labelSet_.insert(f[fp]);
-            }
+            set.insert(fcs[facei]);
         }
 
         storage.clear();
-        if (labelSet_.size() > storage.capacity())
+        if (set.size() > storage.capacity())
         {
-            storage.setCapacity(labelSet_.size());
+            storage.setCapacity(set.size());
         }
 
-        forAllConstIter(labelHashSet, labelSet_, iter)
+        for (const label pointi : set)
         {
-            storage.append(iter.key());
+            storage.append(pointi);
         }
 
         return storage;
@@ -100,10 +96,8 @@ const Foam::labelList& Foam::primitiveMesh::cellPoints
 
 const Foam::labelList& Foam::primitiveMesh::cellPoints(const label celli) const
 {
-    return cellPoints(celli, labels_);
+    return cellPoints(celli, labelSet_, labels_);
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
 // ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C
index 3cf689595ab870c71bd8f006206347c99894497b..cc91132dc6783a0a14b7826f1e286da5f37760ed 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C
@@ -630,6 +630,7 @@ const Foam::labelList& Foam::primitiveMesh::faceEdges(const label facei) const
 const Foam::labelList& Foam::primitiveMesh::cellEdges
 (
     const label celli,
+    labelHashSet& set,
     DynamicList<label>& storage
 ) const
 {
@@ -641,28 +642,22 @@ const Foam::labelList& Foam::primitiveMesh::cellEdges
     {
         const labelList& cFaces = cells()[celli];
 
-        labelSet_.clear();
+        set.clear();
 
-        forAll(cFaces, i)
+        for (const label facei : cFaces)
         {
-            const labelList& fe = faceEdges(cFaces[i]);
-
-            forAll(fe, feI)
-            {
-                labelSet_.insert(fe[feI]);
-            }
+            set.insert(faceEdges(facei));
         }
 
         storage.clear();
-
-        if (labelSet_.size() > storage.capacity())
+        if (set.size() > storage.capacity())
         {
-            storage.setCapacity(labelSet_.size());
+            storage.setCapacity(set.size());
         }
 
-        forAllConstIter(labelHashSet, labelSet_, iter)
+        for (const label edgei : set)
         {
-            storage.append(iter.key());
+            storage.append(edgei);
         }
 
         return storage;
@@ -672,7 +667,7 @@ const Foam::labelList& Foam::primitiveMesh::cellEdges
 
 const Foam::labelList& Foam::primitiveMesh::cellEdges(const label celli) const
 {
-    return cellEdges(celli, labels_);
+    return cellEdges(celli, labelSet_, labels_);
 }
 
 
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/combineFaces.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/combineFaces.C
index 10e5ee81c631ba9196338feb9cc2af0e2365a7a6..6499498e3e2438289ac3e88781b302d5a6c99144 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/combineFaces.C
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/combineFaces.C
@@ -305,6 +305,7 @@ Foam::labelListList Foam::combineFaces::getMergeSets
     // Lists of faces that can be merged.
     DynamicList<labelList> allFaceSets(boundaryCells.size() / 10);
     // Storage for on-the-fly cell-edge addressing.
+    labelHashSet set;
     DynamicList<label> storage;
 
     // On all cells regionise the faces
@@ -314,7 +315,7 @@ Foam::labelListList Foam::combineFaces::getMergeSets
 
         const cell& cFaces = mesh_.cells()[celli];
 
-        const labelList& cEdges = mesh_.cellEdges(celli, storage);
+        const labelList& cEdges = mesh_.cellEdges(celli, set, storage);
 
         // Region per face
         Map<label> faceRegion(cFaces.size());
diff --git a/src/mesh/snappyHexMesh/meshRefinement/meshRefinementProblemCells.C b/src/mesh/snappyHexMesh/meshRefinement/meshRefinementProblemCells.C
index d85faeb02c9bb3f120f5ffe953833279c82e1968..ffd42013dc3244d9f84e3c2f8afa6bfab5e1bb06 100644
--- a/src/mesh/snappyHexMesh/meshRefinement/meshRefinementProblemCells.C
+++ b/src/mesh/snappyHexMesh/meshRefinement/meshRefinementProblemCells.C
@@ -718,10 +718,11 @@ Foam::labelList Foam::meshRefinement::markFacesOnProblemCells
     // On-the-fly addressing storage.
     DynamicList<label> dynFEdges;
     DynamicList<label> dynCPoints;
+    labelHashSet pSet;
 
     forAll(cellLevel, celli)
     {
-        const labelList& cPoints = mesh_.cellPoints(celli, dynCPoints);
+        const labelList& cPoints = mesh_.cellPoints(celli, pSet, dynCPoints);
 
         // Get number of anchor points (pointLevel <= cellLevel)
 
@@ -729,10 +730,8 @@ Foam::labelList Foam::meshRefinement::markFacesOnProblemCells
         label nNonAnchorBoundary = 0;
         label nonBoundaryAnchor = -1;
 
-        forAll(cPoints, i)
+        for (const label pointi : cPoints)
         {
-            label pointi = cPoints[i];
-
             if (pointLevel[pointi] <= cellLevel[celli])
             {
                 // Anchor point