diff --git a/applications/test/nearWallDist-wave/testWallDist2.C b/applications/test/nearWallDist-wave/testWallDist2.C
index 364a73279f..45690a33f4 100644
--- a/applications/test/nearWallDist-wave/testWallDist2.C
+++ b/applications/test/nearWallDist-wave/testWallDist2.C
@@ -29,7 +29,7 @@ Description
 
 #include "fvCFD.H"
 #include "wallFvPatch.H"
-#include "MeshWave.H"
+#include "FaceCellWave.H"
 #include "wallPoint.H"
 
 
@@ -158,7 +158,7 @@ int main(int argc, char *argv[])
         forAll(allCellInfo, cellI)
         {
             scalar dist = allCellInfo[cellI].distSqr();
-            if (allCellInfo[cellI].valid())
+            if (allCellInfo[cellI].valid(wallDistCalc.data()))
             {
                 wallDistUncorrected[cellI] = Foam::sqrt(dist);
             }
@@ -180,7 +180,7 @@ int main(int argc, char *argv[])
                 const label meshFaceI = patchField.patch().start() + patchFaceI;
 
                 scalar dist = allFaceInfo[meshFaceI].distSqr();
-                if (allFaceInfo[meshFaceI].valid())
+                if (allFaceInfo[meshFaceI].valid(wallDistCalc.data()))
                 {
                     patchField[patchFaceI] = Foam::sqrt(dist);
                 }
diff --git a/applications/utilities/mesh/advanced/selectCells/selectCells.C b/applications/utilities/mesh/advanced/selectCells/selectCells.C
index 827449c02c..8f86a145d5 100644
--- a/applications/utilities/mesh/advanced/selectCells/selectCells.C
+++ b/applications/utilities/mesh/advanced/selectCells/selectCells.C
@@ -299,7 +299,7 @@ label selectOutsideCells
         mesh,
         outsideFaces.shrink(),
         outsideFacesInfo.shrink(),
-        mesh.globalData().nTotalCells()+1  // max iterations
+        mesh.globalData().nTotalCells()+1   // max iterations
     );
 
     // Now regionCalc should hold info on cells that are reachable from
diff --git a/src/OpenFOAM/algorithms/MeshWave/FaceCellWave.C b/src/OpenFOAM/algorithms/MeshWave/FaceCellWave.C
index 3a9155614f..a653afd688 100644
--- a/src/OpenFOAM/algorithms/MeshWave/FaceCellWave.C
+++ b/src/OpenFOAM/algorithms/MeshWave/FaceCellWave.C
@@ -35,75 +35,16 @@ License
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
-template <class Type>
-const Foam::scalar Foam::FaceCellWave<Type>::geomTol_ = 1e-6;
+template <class Type, class TrackingData>
+const Foam::scalar
+Foam::FaceCellWave<Type, TrackingData>::geomTol_ = 1e-6;
 
-template <class Type>
-Foam::scalar Foam::FaceCellWave<Type>::propagationTol_ = 0.01;
+template <class Type, class TrackingData>
+const Foam::scalar
+Foam::FaceCellWave<Type, TrackingData>::propagationTol_ = 0.01;
 
-// Write to ostream
-template <class Type>
-Foam::Ostream& Foam::FaceCellWave<Type>::writeFaces
-(
-    const label nFaces,
-    const labelList& faceLabels,
-    const List<Type>& faceInfo,
-    Ostream& os
-)
-{
-    // Write list contents depending on data format
-    if (os.format() == IOstream::ASCII)
-    {
-        os << nFaces;
-
-        for (label i = 0; i < nFaces; i++)
-        {
-            os << ' ' << faceLabels[i];
-        }
-        for (label i = 0; i < nFaces; i++)
-        {
-            os << ' ' << faceInfo[i];
-        }
-    }
-    else
-    {
-        os << nFaces;
-
-        for (label i = 0; i < nFaces; i++)
-        {
-            os << faceLabels[i];
-        }
-        for (label i = 0; i < nFaces; i++)
-        {
-            os << faceInfo[i];
-        }
-    }
-    return os;
-}
-
-
-// Read from istream
-template <class Type>
-Foam::Istream& Foam::FaceCellWave<Type>::readFaces
-(
-    label& nFaces,
-    labelList& faceLabels,
-    List<Type>& faceInfo,
-    Istream& is
-)
-{
-    is >> nFaces;
-
-    for (label i = 0; i < nFaces; i++)
-    {
-        is >> faceLabels[i];
-    }
-    for (label i = 0; i < nFaces; i++)
-    {
-        is >> faceInfo[i];
-    }
-    return is;
-}
+template <class Type, class TrackingData>
+Foam::label Foam::FaceCellWave<Type, TrackingData>::dummyTrackData_ = 12345;
 
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
@@ -113,8 +54,8 @@ Foam::Istream& Foam::FaceCellWave<Type>::readFaces
 // Updates:
 //      - changedCell_, changedCells_, nChangedCells_,
 //      - statistics: nEvals_, nUnvisitedCells_
-template <class Type>
-bool Foam::FaceCellWave<Type>::updateCell
+template <class Type, class TrackingData>
+bool Foam::FaceCellWave<Type, TrackingData>::updateCell
 (
     const label cellI,
     const label neighbourFaceI,
@@ -125,7 +66,7 @@ bool Foam::FaceCellWave<Type>::updateCell
 {
     nEvals_++;
 
-    bool wasValid = cellInfo.valid();
+    bool wasValid = cellInfo.valid(td_);
 
     bool propagate =
         cellInfo.updateCell
@@ -134,7 +75,8 @@ bool Foam::FaceCellWave<Type>::updateCell
             cellI,
             neighbourFaceI,
             neighbourInfo,
-            tol
+            tol,
+            td_
         );
 
     if (propagate)
@@ -146,7 +88,7 @@ bool Foam::FaceCellWave<Type>::updateCell
         }
     }
 
-    if (!wasValid && cellInfo.valid())
+    if (!wasValid && cellInfo.valid(td_))
     {
         --nUnvisitedCells_;
     }
@@ -160,8 +102,8 @@ bool Foam::FaceCellWave<Type>::updateCell
 // Updates:
 //      - changedFace_, changedFaces_, nChangedFaces_,
 //      - statistics: nEvals_, nUnvisitedFaces_
-template <class Type>
-bool Foam::FaceCellWave<Type>::updateFace
+template <class Type, class TrackingData>
+bool Foam::FaceCellWave<Type, TrackingData>::updateFace
 (
     const label faceI,
     const label neighbourCellI,
@@ -172,7 +114,7 @@ bool Foam::FaceCellWave<Type>::updateFace
 {
     nEvals_++;
 
-    bool wasValid = faceInfo.valid();
+    bool wasValid = faceInfo.valid(td_);
 
     bool propagate =
         faceInfo.updateFace
@@ -181,7 +123,8 @@ bool Foam::FaceCellWave<Type>::updateFace
             faceI,
             neighbourCellI,
             neighbourInfo,
-            tol
+            tol,
+            td_
         );
 
     if (propagate)
@@ -193,7 +136,7 @@ bool Foam::FaceCellWave<Type>::updateFace
         }
     }
 
-    if (!wasValid && faceInfo.valid())
+    if (!wasValid && faceInfo.valid(td_))
     {
         --nUnvisitedFaces_;
     }
@@ -207,8 +150,8 @@ bool Foam::FaceCellWave<Type>::updateFace
 // Updates:
 //      - changedFace_, changedFaces_, nChangedFaces_,
 //      - statistics: nEvals_, nUnvisitedFaces_
-template <class Type>
-bool Foam::FaceCellWave<Type>::updateFace
+template <class Type, class TrackingData>
+bool Foam::FaceCellWave<Type, TrackingData>::updateFace
 (
     const label faceI,
     const Type& neighbourInfo,
@@ -218,7 +161,7 @@ bool Foam::FaceCellWave<Type>::updateFace
 {
     nEvals_++;
 
-    bool wasValid = faceInfo.valid();
+    bool wasValid = faceInfo.valid(td_);
 
     bool propagate =
         faceInfo.updateFace
@@ -226,7 +169,8 @@ bool Foam::FaceCellWave<Type>::updateFace
             mesh_,
             faceI,
             neighbourInfo,
-            tol
+            tol,
+            td_
         );
 
     if (propagate)
@@ -238,7 +182,7 @@ bool Foam::FaceCellWave<Type>::updateFace
         }
     }
 
-    if (!wasValid && faceInfo.valid())
+    if (!wasValid && faceInfo.valid(td_))
     {
         --nUnvisitedFaces_;
     }
@@ -248,8 +192,11 @@ bool Foam::FaceCellWave<Type>::updateFace
 
 
 // For debugging: check status on both sides of cyclic
-template <class Type>
-void Foam::FaceCellWave<Type>::checkCyclic(const polyPatch& patch) const
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::checkCyclic
+(
+    const polyPatch& patch
+) const
 {
     const cyclicPolyPatch& nbrPatch =
         refCast<const cyclicPolyPatch>(patch).neighbPatch();
@@ -259,10 +206,22 @@ void Foam::FaceCellWave<Type>::checkCyclic(const polyPatch& patch) const
         label i1 = patch.start() + patchFaceI;
         label i2 = nbrPatch.start() + patchFaceI;
 
-        if (!allFaceInfo_[i1].sameGeometry(mesh_, allFaceInfo_[i2], geomTol_))
+        if
+        (
+           !allFaceInfo_[i1].sameGeometry
+            (
+                mesh_,
+                allFaceInfo_[i2],
+                geomTol_,
+                td_
+            )
+        )
         {
-            FatalErrorIn("FaceCellWave<Type>::checkCyclic(const polyPatch&)")
-                << "problem: i:" << i1 << "  otheri:" << i2
+            FatalErrorIn
+            (
+                "FaceCellWave<Type, TrackingData>"
+                "::checkCyclic(const polyPatch&)"
+            )   << "problem: i:" << i1 << "  otheri:" << i2
                 << "   faceInfo:" << allFaceInfo_[i1]
                 << "   otherfaceInfo:" << allFaceInfo_[i2]
                 << abort(FatalError);
@@ -270,8 +229,11 @@ void Foam::FaceCellWave<Type>::checkCyclic(const polyPatch& patch) const
 
         if (changedFace_[i1] != changedFace_[i2])
         {
-            FatalErrorIn("FaceCellWave<Type>::checkCyclic(const polyPatch&)")
-                << " problem: i:" << i1 << "  otheri:" << i2
+            FatalErrorIn
+            (
+                "FaceCellWave<Type, TrackingData>"
+                "::checkCyclic(const polyPatch&)"
+            )   << " problem: i:" << i1 << "  otheri:" << i2
                 << "   faceInfo:" << allFaceInfo_[i1]
                 << "   otherfaceInfo:" << allFaceInfo_[i2]
                 << "   changedFace:" << changedFace_[i1]
@@ -283,8 +245,8 @@ void Foam::FaceCellWave<Type>::checkCyclic(const polyPatch& patch) const
 
 
 // Check if has cyclic patches
-template <class Type>
-bool Foam::FaceCellWave<Type>::hasCyclicPatch() const
+template <class Type, class TrackingData>
+bool Foam::FaceCellWave<Type, TrackingData>::hasCyclicPatch() const
 {
     forAll(mesh_.boundaryMesh(), patchI)
     {
@@ -298,8 +260,8 @@ bool Foam::FaceCellWave<Type>::hasCyclicPatch() const
 
 
 // Copy face information into member data
-template <class Type>
-void Foam::FaceCellWave<Type>::setFaceInfo
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::setFaceInfo
 (
     const labelList& changedFaces,
     const List<Type>& changedFacesInfo
@@ -309,13 +271,13 @@ void Foam::FaceCellWave<Type>::setFaceInfo
     {
         label faceI = changedFaces[changedFaceI];
 
-        bool wasValid = allFaceInfo_[faceI].valid();
+        bool wasValid = allFaceInfo_[faceI].valid(td_);
 
         // Copy info for faceI
         allFaceInfo_[faceI] = changedFacesInfo[changedFaceI];
 
         // Maintain count of unset faces
-        if (!wasValid && allFaceInfo_[faceI].valid())
+        if (!wasValid && allFaceInfo_[faceI].valid(td_))
         {
             --nUnvisitedFaces_;
         }
@@ -329,8 +291,8 @@ void Foam::FaceCellWave<Type>::setFaceInfo
 
 
 // Merge face information into member data
-template <class Type>
-void Foam::FaceCellWave<Type>::mergeFaceInfo
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::mergeFaceInfo
 (
     const polyPatch& patch,
     const label nFaces,
@@ -347,7 +309,7 @@ void Foam::FaceCellWave<Type>::mergeFaceInfo
 
         Type& currentWallInfo = allFaceInfo_[meshFaceI];
 
-        if (currentWallInfo != neighbourWallInfo)
+        if (!currentWallInfo.equal(neighbourWallInfo, td_))
         {
             updateFace
             (
@@ -364,8 +326,8 @@ void Foam::FaceCellWave<Type>::mergeFaceInfo
 // Construct compact patchFace change arrays for a (slice of a) single patch.
 // changedPatchFaces in local patch numbering.
 // Return length of arrays.
-template <class Type>
-Foam::label Foam::FaceCellWave<Type>::getChangedPatchFaces
+template <class Type, class TrackingData>
+Foam::label Foam::FaceCellWave<Type, TrackingData>::getChangedPatchFaces
 (
     const polyPatch& patch,
     const label startFaceI,
@@ -394,8 +356,8 @@ Foam::label Foam::FaceCellWave<Type>::getChangedPatchFaces
 
 
 // Handle leaving domain. Implementation referred to Type
-template <class Type>
-void Foam::FaceCellWave<Type>::leaveDomain
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::leaveDomain
 (
     const polyPatch& patch,
     const label nFaces,
@@ -410,14 +372,14 @@ void Foam::FaceCellWave<Type>::leaveDomain
         label patchFaceI = faceLabels[i];
 
         label meshFaceI = patch.start() + patchFaceI;
-        faceInfo[i].leaveDomain(mesh_, patch, patchFaceI, fc[meshFaceI]);
+        faceInfo[i].leaveDomain(mesh_, patch, patchFaceI, fc[meshFaceI], td_);
     }
 }
 
 
 // Handle entering domain. Implementation referred to Type
-template <class Type>
-void Foam::FaceCellWave<Type>::enterDomain
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::enterDomain
 (
     const polyPatch& patch,
     const label nFaces,
@@ -432,14 +394,14 @@ void Foam::FaceCellWave<Type>::enterDomain
         label patchFaceI = faceLabels[i];
 
         label meshFaceI = patch.start() + patchFaceI;
-        faceInfo[i].enterDomain(mesh_, patch, patchFaceI, fc[meshFaceI]);
+        faceInfo[i].enterDomain(mesh_, patch, patchFaceI, fc[meshFaceI], td_);
     }
 }
 
 
 // Transform. Implementation referred to Type
-template <class Type>
-void Foam::FaceCellWave<Type>::transform
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::transform
 (
     const tensorField& rotTensor,
     const label nFaces,
@@ -452,22 +414,22 @@ void Foam::FaceCellWave<Type>::transform
 
         for (label faceI = 0; faceI < nFaces; faceI++)
         {
-            faceInfo[faceI].transform(mesh_, T);
+            faceInfo[faceI].transform(mesh_, T, td_);
         }
     }
     else
     {
         for (label faceI = 0; faceI < nFaces; faceI++)
         {
-            faceInfo[faceI].transform(mesh_, rotTensor[faceI]);
+            faceInfo[faceI].transform(mesh_, rotTensor[faceI], td_);
         }
     }
 }
 
 
 // Offset mesh face. Used for transferring from one cyclic half to the other.
-template <class Type>
-void Foam::FaceCellWave<Type>::offset
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::offset
 (
     const polyPatch&,
     const label cycOffset,
@@ -483,8 +445,8 @@ void Foam::FaceCellWave<Type>::offset
 
 
 // Tranfer all the information to/from neighbouring processors
-template <class Type>
-void Foam::FaceCellWave<Type>::handleProcPatches()
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::handleProcPatches()
 {
     // Send all
 
@@ -532,7 +494,11 @@ void Foam::FaceCellWave<Type>::handleProcPatches()
             }
 
             UOPstream toNeighbour(procPatch.neighbProcNo(), pBufs);
-            writeFaces(nSendFaces, sendFaces, sendFacesInfo, toNeighbour);
+            //writeFaces(nSendFaces, sendFaces, sendFacesInfo, toNeighbour);
+            toNeighbour
+                << SubList<label>(sendFaces, nSendFaces)
+                << SubList<Type>(sendFacesInfo, nSendFaces);
+
         }
     }
 
@@ -550,26 +516,19 @@ void Foam::FaceCellWave<Type>::handleProcPatches()
                 refCast<const processorPolyPatch>(patch);
 
             // Allocate buffers
-            label nReceiveFaces = 0;
             labelList receiveFaces(patch.size());
             List<Type> receiveFacesInfo(patch.size());
 
             {
                 UIPstream fromNeighbour(procPatch.neighbProcNo(), pBufs);
-                readFaces
-                (
-                    nReceiveFaces,
-                    receiveFaces,
-                    receiveFacesInfo,
-                    fromNeighbour
-                );
+                fromNeighbour >> receiveFaces >> receiveFacesInfo;
             }
 
             if (debug)
             {
                 Pout<< " Processor patch " << patchI << ' ' << patch.name()
                     << " communicating with " << procPatch.neighbProcNo()
-                    << "  Receiving:" << nReceiveFaces
+                    << "  Receiving:" << receiveFaces.size()
                     << endl;
             }
 
@@ -579,7 +538,7 @@ void Foam::FaceCellWave<Type>::handleProcPatches()
                 transform
                 (
                     procPatch.reverseT(),
-                    nReceiveFaces,
+                    receiveFaces.size(),
                     receiveFacesInfo
                 );
             }
@@ -588,7 +547,7 @@ void Foam::FaceCellWave<Type>::handleProcPatches()
             enterDomain
             (
                 patch,
-                nReceiveFaces,
+                receiveFaces.size(),
                 receiveFaces,
                 receiveFacesInfo
             );
@@ -597,7 +556,7 @@ void Foam::FaceCellWave<Type>::handleProcPatches()
             mergeFaceInfo
             (
                 patch,
-                nReceiveFaces,
+                receiveFaces.size(),
                 receiveFaces,
                 receiveFacesInfo
             );
@@ -607,10 +566,8 @@ void Foam::FaceCellWave<Type>::handleProcPatches()
 
 
 // Transfer information across cyclic halves.
-// Note: Cyclic is two patches in one: one side from 0..size/2 and other
-// side from size/2 .. size.
-template <class Type>
-void Foam::FaceCellWave<Type>::handleCyclicPatches()
+template <class Type, class TrackingData>
+void Foam::FaceCellWave<Type, TrackingData>::handleCyclicPatches()
 {
     forAll(mesh_.boundaryMesh(), patchI)
     {
@@ -696,17 +653,19 @@ void Foam::FaceCellWave<Type>::handleCyclicPatches()
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Set up only. Use setFaceInfo and iterate() to do actual calculation.
-template <class Type>
-Foam::FaceCellWave<Type>::FaceCellWave
+template <class Type, class TrackingData>
+Foam::FaceCellWave<Type, TrackingData>::FaceCellWave
 (
     const polyMesh& mesh,
     UList<Type>& allFaceInfo,
-    UList<Type>& allCellInfo
+    UList<Type>& allCellInfo,
+    TrackingData& td
 )
 :
     mesh_(mesh),
     allFaceInfo_(allFaceInfo),
     allCellInfo_(allCellInfo),
+    td_(td),
     changedFace_(mesh_.nFaces(), false),
     changedFaces_(mesh_.nFaces()),
     nChangedFaces_(0),
@@ -723,20 +682,22 @@ Foam::FaceCellWave<Type>::FaceCellWave
 
 // Iterate, propagating changedFacesInfo across mesh, until no change (or
 // maxIter reached). Initial cell values specified.
-template <class Type>
-Foam::FaceCellWave<Type>::FaceCellWave
+template <class Type, class TrackingData>
+Foam::FaceCellWave<Type, TrackingData>::FaceCellWave
 (
     const polyMesh& mesh,
     const labelList& changedFaces,
     const List<Type>& changedFacesInfo,
     UList<Type>& allFaceInfo,
     UList<Type>& allCellInfo,
-    const label maxIter
+    const label maxIter,
+    TrackingData& td
 )
 :
     mesh_(mesh),
     allFaceInfo_(allFaceInfo),
     allCellInfo_(allCellInfo),
+    td_(td),
     changedFace_(mesh_.nFaces(), false),
     changedFaces_(mesh_.nFaces()),
     nChangedFaces_(0),
@@ -759,7 +720,7 @@ Foam::FaceCellWave<Type>::FaceCellWave
     {
         FatalErrorIn
         (
-            "FaceCellWave<Type>::FaceCellWave"
+            "FaceCellWave<Type, TrackingData>::FaceCellWave"
             "(const polyMesh&, const labelList&, const List<Type>,"
             " UList<Type>&, UList<Type>&, const label maxIter)"
         )
@@ -775,15 +736,15 @@ Foam::FaceCellWave<Type>::FaceCellWave
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 
-template <class Type>
-Foam::label Foam::FaceCellWave<Type>::getUnsetCells() const
+template <class Type, class TrackingData>
+Foam::label Foam::FaceCellWave<Type, TrackingData>::getUnsetCells() const
 {
     return nUnvisitedCells_;
 }
 
 
-template <class Type>
-Foam::label Foam::FaceCellWave<Type>::getUnsetFaces() const
+template <class Type, class TrackingData>
+Foam::label Foam::FaceCellWave<Type, TrackingData>::getUnsetFaces() const
 {
     return nUnvisitedFaces_;
 }
@@ -791,8 +752,8 @@ Foam::label Foam::FaceCellWave<Type>::getUnsetFaces() const
 
 
 // Propagate cell to face
-template <class Type>
-Foam::label Foam::FaceCellWave<Type>::faceToCell()
+template <class Type, class TrackingData>
+Foam::label Foam::FaceCellWave<Type, TrackingData>::faceToCell()
 {
     const labelList& owner = mesh_.faceOwner();
     const labelList& neighbour = mesh_.faceNeighbour();
@@ -808,7 +769,7 @@ Foam::label Foam::FaceCellWave<Type>::faceToCell()
         label faceI = changedFaces_[changedFaceI];
         if (!changedFace_[faceI])
         {
-            FatalErrorIn("FaceCellWave<Type>::faceToCell()")
+            FatalErrorIn("FaceCellWave<Type, TrackingData>::faceToCell()")
                 << "Face " << faceI
                 << " not marked as having been changed"
                 << abort(FatalError);
@@ -823,7 +784,7 @@ Foam::label Foam::FaceCellWave<Type>::faceToCell()
         label cellI = owner[faceI];
         Type& currentWallInfo = allCellInfo_[cellI];
 
-        if (currentWallInfo != neighbourWallInfo)
+        if (!currentWallInfo.equal(neighbourWallInfo, td_))
         {
             updateCell
             (
@@ -841,7 +802,7 @@ Foam::label Foam::FaceCellWave<Type>::faceToCell()
             cellI = neighbour[faceI];
             Type& currentWallInfo2 = allCellInfo_[cellI];
 
-            if (currentWallInfo2 != neighbourWallInfo)
+            if (!currentWallInfo2.equal(neighbourWallInfo, td_))
             {
                 updateCell
                 (
@@ -876,8 +837,8 @@ Foam::label Foam::FaceCellWave<Type>::faceToCell()
 
 
 // Propagate cell to face
-template <class Type>
-Foam::label Foam::FaceCellWave<Type>::cellToFace()
+template <class Type, class TrackingData>
+Foam::label Foam::FaceCellWave<Type, TrackingData>::cellToFace()
 {
     const cellList& cells = mesh_.cells();
 
@@ -891,7 +852,7 @@ Foam::label Foam::FaceCellWave<Type>::cellToFace()
         label cellI = changedCells_[changedCellI];
         if (!changedCell_[cellI])
         {
-            FatalErrorIn("FaceCellWave<Type>::cellToFace()") << "Cell " << cellI
+            FatalErrorIn("FaceCellWave<Type, TrackingData>::cellToFace()") << "Cell " << cellI
                 << " not marked as having been changed"
                 << abort(FatalError);
         }
@@ -906,7 +867,7 @@ Foam::label Foam::FaceCellWave<Type>::cellToFace()
             label faceI = faceLabels[faceLabelI];
             Type& currentWallInfo = allFaceInfo_[faceI];
 
-            if (currentWallInfo != neighbourWallInfo)
+            if (!currentWallInfo.equal(neighbourWallInfo, td_))
             {
                 updateFace
                 (
@@ -952,8 +913,8 @@ Foam::label Foam::FaceCellWave<Type>::cellToFace()
 
 
 // Iterate
-template <class Type>
-Foam::label Foam::FaceCellWave<Type>::iterate(const label maxIter)
+template <class Type, class TrackingData>
+Foam::label Foam::FaceCellWave<Type, TrackingData>::iterate(const label maxIter)
 {
     if (hasCyclicPatches_)
     {
diff --git a/src/OpenFOAM/algorithms/MeshWave/FaceCellWave.H b/src/OpenFOAM/algorithms/MeshWave/FaceCellWave.H
index 457fb7a413..6ad8a1ae32 100644
--- a/src/OpenFOAM/algorithms/MeshWave/FaceCellWave.H
+++ b/src/OpenFOAM/algorithms/MeshWave/FaceCellWave.H
@@ -71,7 +71,7 @@ TemplateName(FaceCellWave);
                            Class FaceCellWave Declaration
 \*---------------------------------------------------------------------------*/
 
-template <class Type>
+template<class Type, class TrackingData = int>
 class FaceCellWave
 :
     public FaceCellWaveName
@@ -87,6 +87,9 @@ class FaceCellWave
         //- Information for all cells
         UList<Type>& allCellInfo_;
 
+        //- Additional data to be passed into container
+        TrackingData& td_;
+
         //- Has face changed
         boolList changedFace_;
 
@@ -115,27 +118,6 @@ class FaceCellWave
         label iter_;
 
 
-    // Static Functions
-
-        //- Write faces info
-        static Ostream& writeFaces
-        (
-            const label nFaces,
-            const labelList& faceLabels,
-            const List<Type>& faceInfo,
-            Ostream& os
-        );
-
-        //- Read faces info
-        static Istream& readFaces
-        (
-            label& nFaces,
-            labelList& faceLabels,
-            List<Type>& faceInfo,
-            Istream& is
-        );
-
-
     // Private Member Functions
 
         //- Disallow default bitwise copy construct
@@ -250,7 +232,12 @@ class FaceCellWave
       // Private static data
 
             static const scalar geomTol_;
-            static scalar propagationTol_;
+            static const scalar propagationTol_;
+
+            //- Used as default trackdata value to satisfy default template
+            //  argument.
+            static label dummyTrackData_;
+
 
 public:
 
@@ -277,7 +264,8 @@ public:
         (
             const polyMesh&,
             UList<Type>& allFaceInfo,
-            UList<Type>& allCellInfo
+            UList<Type>& allCellInfo,
+            TrackingData& td = dummyTrackData_
         );
 
         //- Construct from mesh and list of changed faces with the Type
@@ -290,7 +278,8 @@ public:
             const List<Type>& changedFacesInfo,
             UList<Type>& allFaceInfo,
             UList<Type>& allCellInfo,
-            const label maxIter
+            const label maxIter,
+            TrackingData& td = dummyTrackData_
         );
 
 
@@ -310,6 +299,12 @@ public:
                 return allCellInfo_;
             }
 
+            //- Additional data to be passed into container
+            const TrackingData& data() const
+            {
+                return td_;
+            }
+
             //- Access mesh
             const polyMesh& mesh() const
             {
diff --git a/src/OpenFOAM/algorithms/MeshWave/MeshWave.C b/src/OpenFOAM/algorithms/MeshWave/MeshWave.C
index 4e887dfc53..461cf22495 100644
--- a/src/OpenFOAM/algorithms/MeshWave/MeshWave.C
+++ b/src/OpenFOAM/algorithms/MeshWave/MeshWave.C
@@ -26,17 +26,25 @@ License
 #include "MeshWave.H"
 #include "polyMesh.H"
 
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+template <class Type, class TrackingData>
+Foam::label Foam::MeshWave<Type, TrackingData>::dummyTrackData_ = 12345;
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Iterate, propagating changedFacesInfo across mesh, until no change (or
 // maxIter reached).
-template <class Type>
-Foam::MeshWave<Type>::MeshWave
+template <class Type, class TrackingData>
+Foam::MeshWave<Type, TrackingData>::MeshWave
 (
     const polyMesh& mesh,
     const labelList& changedFaces,
     const List<Type>& changedFacesInfo,
-    const label maxIter
+    const label maxIter,
+    TrackingData& td
 )
 :
     allFaceInfo_(mesh.nFaces()),
@@ -48,21 +56,23 @@ Foam::MeshWave<Type>::MeshWave
         changedFacesInfo,
         allFaceInfo_,
         allCellInfo_,
-        maxIter
+        maxIter,
+        td
     )
 {}
 
 
 // Iterate, propagating changedFacesInfo across mesh, until no change (or
 // maxIter reached). Initial cell values specified.
-template <class Type>
-Foam::MeshWave<Type>::MeshWave
+template <class Type, class TrackingData>
+Foam::MeshWave<Type, TrackingData>::MeshWave
 (
     const polyMesh& mesh,
     const labelList& changedFaces,
     const List<Type>& changedFacesInfo,
     const List<Type>& allCellInfo,
-    const label maxIter
+    const label maxIter,
+    TrackingData& td
 )
 :
     allFaceInfo_(mesh.nFaces()),
@@ -74,7 +84,8 @@ Foam::MeshWave<Type>::MeshWave
         changedFacesInfo,
         allFaceInfo_,
         allCellInfo_,
-        maxIter
+        maxIter,
+        td
     )
 {}
 
diff --git a/src/OpenFOAM/algorithms/MeshWave/MeshWave.H b/src/OpenFOAM/algorithms/MeshWave/MeshWave.H
index 13d6561224..a93d01cfab 100644
--- a/src/OpenFOAM/algorithms/MeshWave/MeshWave.H
+++ b/src/OpenFOAM/algorithms/MeshWave/MeshWave.H
@@ -53,7 +53,7 @@ TemplateName(MeshWave);
                            Class MeshWave Declaration
 \*---------------------------------------------------------------------------*/
 
-template <class Type>
+template <class Type, class TrackingData = int>
 class MeshWave
 :
     public MeshWaveName
@@ -67,7 +67,7 @@ class MeshWave
         List<Type> allCellInfo_;
 
         //- Wave calculation engine.
-        FaceCellWave<Type> calc_;
+        FaceCellWave<Type, TrackingData> calc_;
 
     // Private Member Functions
 
@@ -77,6 +77,12 @@ class MeshWave
         //- Disallow default bitwise assignment
         void operator=(const MeshWave&);
 
+      // Private static data
+
+            //- Used as default trackdata value to satisfy default template
+            //  argument.
+            static label dummyTrackData_;
+
 public:
 
     // Constructors
@@ -89,7 +95,8 @@ public:
             const polyMesh& mesh,
             const labelList& initialChangedFaces,
             const List<Type>& changedFacesInfo,
-            const label maxIter
+            const label maxIter,
+            TrackingData& td = dummyTrackData_
         );
 
         //- Construct from mesh, list of changed faces with the Type
@@ -102,7 +109,8 @@ public:
             const labelList& initialChangedFaces,
             const List<Type>& changedFacesInfo,
             const List<Type>& allCellInfo,
-            const label maxIter
+            const label maxIter,
+            TrackingData& td = dummyTrackData_
         );
 
 
@@ -120,6 +128,11 @@ public:
             return allCellInfo_;
         }
 
+        //- Additional data to be passed into container
+        const TrackingData& data() const
+        {
+            return calc_.data();
+        }
 
         //- Iterate until no changes or maxIter reached. Returns number of
         //  unset cells (see getUnsetCells)
diff --git a/src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.C b/src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.C
index 1d77ff45ee..9c904bbc6b 100644
--- a/src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.C
+++ b/src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.C
@@ -25,35 +25,13 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "directionInfo.H"
-#include "hexMatcher.H"
-#include "meshTools.H"
+//#include "hexMatcher.H"
+//#include "meshTools.H"
 #include "polyMesh.H"
 
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-bool Foam::directionInfo::equal(const edge& e, const label v0, const label v1)
-{
-    return
-        (e.start() == v0 && e.end() == v1)
-     || (e.start() == v1 && e.end() == v0);
-}
-
-
-Foam::point Foam::directionInfo::eMid
-(
-    const primitiveMesh& mesh,
-    const label edgeI
-)
-{
-    const edge& e = mesh.edges()[edgeI];
-
-    return
-        0.5
-      * (mesh.points()[e.start()] + mesh.points()[e.end()]);
-}
-
-
 // Find edge among edgeLabels that uses v0 and v1
 Foam::label Foam::directionInfo::findEdge
 (
@@ -67,9 +45,7 @@ Foam::label Foam::directionInfo::findEdge
     {
         label edgeI = edgeLabels[edgeLabelI];
 
-        const edge& e = mesh.edges()[edgeI];
-
-        if (equal(e, v0, v1))
+        if (mesh.edges()[edgeI] == edge(v0, v1))
         {
             return edgeI;
         }
@@ -225,159 +201,6 @@ Foam::label Foam::directionInfo::edgeToFaceIndex
 }
 
 
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-// Update this cell with neighbouring face information
-bool Foam::directionInfo::updateCell
-(
-    const polyMesh& mesh,
-    const label thisCellI,
-    const label neighbourFaceI,
-    const directionInfo& neighbourInfo,
-    const scalar        // tol
-)
-{
-    if (index_ >= -2)
-    {
-        // Already determined.
-        return false;
-    }
-
-    if (hexMatcher().isA(mesh, thisCellI))
-    {
-        const face& f = mesh.faces()[neighbourFaceI];
-
-        if (neighbourInfo.index() == -2)
-        {
-            // Geometric information from neighbour
-            index_ = -2;
-        }
-        else if (neighbourInfo.index() == -1)
-        {
-            // Cut tangential to face. Take any edge connected to face
-            // but not used in face.
-
-            // Get first edge on face.
-            label edgeI = mesh.faceEdges()[neighbourFaceI][0];
-
-            const edge& e = mesh.edges()[edgeI];
-
-            // Find face connected to face through edgeI and on same cell.
-            label faceI =
-                meshTools::otherFace
-                (
-                    mesh,
-                    thisCellI,
-                    neighbourFaceI,
-                    edgeI
-                );
-
-            // Find edge on faceI which is connected to e.start() but not edgeI.
-            index_ =
-                meshTools::otherEdge
-                (
-                    mesh,
-                    mesh.faceEdges()[faceI],
-                    edgeI,
-                    e.start()
-                );
-        }
-        else
-        {
-            // Index is a vertex on the face. Convert to mesh edge.
-
-            // Get mesh edge between f[index_] and f[index_+1]
-            label v0 = f[neighbourInfo.index()];
-            label v1 = f[(neighbourInfo.index() + 1) % f.size()];
-
-            index_ = findEdge(mesh, mesh.faceEdges()[neighbourFaceI], v0, v1);
-        }
-    }
-    else
-    {
-        // Not a hex so mark this as geometric.
-        index_ = -2;
-    }
-
-
-    n_ = neighbourInfo.n();
-
-    return true;
-}
-
-
-// Update this face with neighbouring cell information
-bool Foam::directionInfo::updateFace
-(
-    const polyMesh& mesh,
-    const label thisFaceI,
-    const label neighbourCellI,
-    const directionInfo& neighbourInfo,
-    const scalar    // tol
-)
-{
-    // Handle special cases first
-
-    if (index_ >= -2)
-    {
-        // Already determined
-        return false;
-    }
-
-    // Handle normal cases where topological or geometrical info comes from
-    // neighbouring cell
-
-    if (neighbourInfo.index() >= 0)
-    {
-        // Neighbour has topological direction (and hence is hex). Find cut
-        // edge on face.
-        index_ =
-            edgeToFaceIndex
-            (
-                mesh,
-                neighbourCellI,
-                thisFaceI,
-                neighbourInfo.index()
-            );
-    }
-    else
-    {
-        // Neighbour has geometric information. Use.
-        index_ = -2;
-    }
-
-
-    n_ = neighbourInfo.n();
-
-    return true;
-}
-
-
-// Merge this with information on same face
-bool Foam::directionInfo::updateFace
-(
-    const polyMesh& mesh,
-    const label,    // thisFaceI
-    const directionInfo& neighbourInfo,
-    const scalar    // tol
-)
-{
-    if (index_ >= -2)
-    {
-        // Already visited.
-        return false;
-    }
-    else
-    {
-        index_ = neighbourInfo.index();
-
-        n_ = neighbourInfo.n();
-
-        return true;
-    }
-}
-
-
 // * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
 
 Foam::Ostream& Foam::operator<<
diff --git a/src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.H b/src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.H
index 814e89a16e..a42c44d41f 100644
--- a/src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.H
+++ b/src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.H
@@ -86,12 +86,6 @@ class directionInfo
 
     // Private Member Functions
 
-        //- edge uses two labels
-        static bool equal(const edge& e, const label, const label);
-
-        //- Calculate mid point of edge.
-        static point eMid(const primitiveMesh& mesh, const label edgeI);
-
         //- Find edge among edgeLabels that uses v0 and v1
         static label findEdge
         (
@@ -158,71 +152,90 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value.
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData& td) const;
 
             //- Check for identical geometrical data. Used for cyclics checking.
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const directionInfo&,
-                const scalar
+                const scalar,
+                TrackingData& td
             ) const;
 
             //- Convert any absolute coordinates into relative to (patch)face
             //  centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor&
+                const tensor&,
+                TrackingData& td
             );
 
             //- Influence of neighbouring face.
-            bool updateCell
+            template<class TrackingData>
+            inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const directionInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell.
-            bool updateFace
+            template<class TrackingData>
+            inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const directionInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face.
-            bool updateFace
+            template<class TrackingData>
+            inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const directionInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const directionInfo&, TrackingData& td) const;
+
     // Member Operators
 
         // Needed for List IO
diff --git a/src/dynamicMesh/meshCut/directions/directionInfo/directionInfoI.H b/src/dynamicMesh/meshCut/directions/directionInfo/directionInfoI.H
index 1f174e4bf3..93ebcc92ce 100644
--- a/src/dynamicMesh/meshCut/directions/directionInfo/directionInfoI.H
+++ b/src/dynamicMesh/meshCut/directions/directionInfo/directionInfoI.H
@@ -24,9 +24,8 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "polyMesh.H"
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
+#include "meshTools.H"
+#include "hexMatcher.H"
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
@@ -60,18 +59,21 @@ inline Foam::directionInfo::directionInfo(const directionInfo& w2)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline bool Foam::directionInfo::valid() const
+template<class TrackingData>
+inline bool Foam::directionInfo::valid(TrackingData& td) const
 {
     return index_ != -3;
 }
 
 
 // No geometric data so never any problem on cyclics
+template<class TrackingData>
 inline bool Foam::directionInfo::sameGeometry
 (
     const polyMesh&,
     const directionInfo& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
  const
 {
@@ -80,12 +82,14 @@ inline bool Foam::directionInfo::sameGeometry
 
 
 // index_ is already offset in face
+template<class TrackingData>
 inline void Foam::directionInfo::leaveDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {}
 
@@ -94,12 +98,14 @@ inline void Foam::directionInfo::leaveDomain
 // (Note: f[0] on other domain is connected to f[0] in this domain,
 //        f[1]      ,,                         f[size-1]   ,,
 // etc.)
+template<class TrackingData>
 inline void Foam::directionInfo::enterDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {
     if (index_ >= 0)
@@ -112,14 +118,184 @@ inline void Foam::directionInfo::enterDomain
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::directionInfo::transform
 (
     const polyMesh&,
-    const tensor& rotTensor
+    const tensor& rotTensor,
+    TrackingData& td
 )
 {}
 
 
+// Update this cell with neighbouring face information
+template<class TrackingData>
+inline bool Foam::directionInfo::updateCell
+(
+    const polyMesh& mesh,
+    const label thisCellI,
+    const label neighbourFaceI,
+    const directionInfo& neighbourInfo,
+    const scalar,       // tol
+    TrackingData& td
+)
+{
+    if (index_ >= -2)
+    {
+        // Already determined.
+        return false;
+    }
+
+    if (hexMatcher().isA(mesh, thisCellI))
+    {
+        const face& f = mesh.faces()[neighbourFaceI];
+
+        if (neighbourInfo.index() == -2)
+        {
+            // Geometric information from neighbour
+            index_ = -2;
+        }
+        else if (neighbourInfo.index() == -1)
+        {
+            // Cut tangential to face. Take any edge connected to face
+            // but not used in face.
+
+            // Get first edge on face.
+            label edgeI = mesh.faceEdges()[neighbourFaceI][0];
+
+            const edge& e = mesh.edges()[edgeI];
+
+            // Find face connected to face through edgeI and on same cell.
+            label faceI =
+                meshTools::otherFace
+                (
+                    mesh,
+                    thisCellI,
+                    neighbourFaceI,
+                    edgeI
+                );
+
+            // Find edge on faceI which is connected to e.start() but not edgeI.
+            index_ =
+                meshTools::otherEdge
+                (
+                    mesh,
+                    mesh.faceEdges()[faceI],
+                    edgeI,
+                    e.start()
+                );
+        }
+        else
+        {
+            // Index is a vertex on the face. Convert to mesh edge.
+
+            // Get mesh edge between f[index_] and f[index_+1]
+            label v0 = f[neighbourInfo.index()];
+            label v1 = f[(neighbourInfo.index() + 1) % f.size()];
+
+            index_ = findEdge(mesh, mesh.faceEdges()[neighbourFaceI], v0, v1);
+        }
+    }
+    else
+    {
+        // Not a hex so mark this as geometric.
+        index_ = -2;
+    }
+
+
+    n_ = neighbourInfo.n();
+
+    return true;
+}
+
+
+// Update this face with neighbouring cell information
+template<class TrackingData>
+inline bool Foam::directionInfo::updateFace
+(
+    const polyMesh& mesh,
+    const label thisFaceI,
+    const label neighbourCellI,
+    const directionInfo& neighbourInfo,
+    const scalar,   // tol
+    TrackingData& td
+)
+{
+    // Handle special cases first
+
+    if (index_ >= -2)
+    {
+        // Already determined
+        return false;
+    }
+
+    // Handle normal cases where topological or geometrical info comes from
+    // neighbouring cell
+
+    if (neighbourInfo.index() >= 0)
+    {
+        // Neighbour has topological direction (and hence is hex). Find cut
+        // edge on face.
+        index_ =
+            edgeToFaceIndex
+            (
+                mesh,
+                neighbourCellI,
+                thisFaceI,
+                neighbourInfo.index()
+            );
+    }
+    else
+    {
+        // Neighbour has geometric information. Use.
+        index_ = -2;
+    }
+
+
+    n_ = neighbourInfo.n();
+
+    return true;
+}
+
+
+// Merge this with information on same face
+template<class TrackingData>
+inline bool Foam::directionInfo::updateFace
+(
+    const polyMesh& mesh,
+    const label,    // thisFaceI
+    const directionInfo& neighbourInfo,
+    const scalar,   // tol
+    TrackingData& td
+)
+{
+    if (index_ >= -2)
+    {
+        // Already visited.
+        return false;
+    }
+    else
+    {
+        index_ = neighbourInfo.index();
+
+        n_ = neighbourInfo.n();
+
+        return true;
+    }
+}
+
+
+template <class TrackingData>
+inline bool Foam::directionInfo::equal
+(
+    const directionInfo& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 inline bool Foam::directionInfo::operator==(const Foam::directionInfo& rhs)
diff --git a/src/dynamicMesh/meshCut/wallLayerCells/wallLayerCells.C b/src/dynamicMesh/meshCut/wallLayerCells/wallLayerCells.C
index 0f82f6777f..c4b991f6fb 100644
--- a/src/dynamicMesh/meshCut/wallLayerCells/wallLayerCells.C
+++ b/src/dynamicMesh/meshCut/wallLayerCells/wallLayerCells.C
@@ -171,7 +171,7 @@ Foam::wallLayerCells::wallLayerCells
         {
             const wallNormalInfo& info = faceInfo[faceI];
 
-            if (info.valid())
+            if (info.valid(regionCalc.data()))
             {
                 const face& f = mesh.faces()[faceI];
 
@@ -215,7 +215,7 @@ Foam::wallLayerCells::wallLayerCells
     {
         const wallNormalInfo& info = cellInfo[cellI];
 
-        if (info.valid() && !usesCoupledPatch(cellI))
+        if (info.valid(regionCalc.data()) && !usesCoupledPatch(cellI))
         {
             refineCells.append(refineCell(cellI, info.normal()));
         }
diff --git a/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfo.C b/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfo.C
index abb1dfc951..ce0add8f25 100644
--- a/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfo.C
+++ b/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfo.C
@@ -25,11 +25,6 @@ License
 
 #include "wallNormalInfo.H"
 
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-Foam::point Foam::wallNormalInfo::greatVector(GREAT, GREAT, GREAT);
-
-
 // * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
 
 Foam::Ostream& Foam::operator<<
diff --git a/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfo.H b/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfo.H
index dfdf9bd8c9..f4553ec47c 100644
--- a/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfo.H
+++ b/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfo.H
@@ -65,15 +65,11 @@ class wallNormalInfo
     // Private Member Functions
 
         //- Evaluate distance to point. Update normal_
-        inline bool update(const wallNormalInfo& w2);
+        template<class TrackingData>
+        inline bool update(const wallNormalInfo& w2, TrackingData& td);
 
 public:
 
-    // Static data members
-
-        //- initial point far away.
-        static point greatVector;
-
     // Constructors
 
         //- Construct null
@@ -96,71 +92,90 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value.
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData& td) const;
 
             //- Check for identical geometrical data. Used for cyclics checking.
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const wallNormalInfo&,
-                const scalar
+                const scalar,
+                TrackingData& td
             ) const;
 
             //- Convert any absolute coordinates into relative to (patch)face
             //  centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor&
+                const tensor&,
+                TrackingData& td
             );
 
             //- Influence of neighbouring face.
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const wallNormalInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const wallNormalInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const wallNormalInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const wallNormalInfo&, TrackingData& td) const;
+
     // Member Operators
 
         // Needed for List IO
diff --git a/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfoI.H b/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfoI.H
index 91569e3239..762a399353 100644
--- a/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfoI.H
+++ b/src/dynamicMesh/meshCut/wallLayerCells/wallNormalInfo/wallNormalInfoI.H
@@ -28,9 +28,14 @@ License
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Update this with w2 if not yet set.
-inline bool Foam::wallNormalInfo::update(const wallNormalInfo& w2)
+template<class TrackingData>
+inline bool Foam::wallNormalInfo::update
+(
+    const wallNormalInfo& w2,
+    TrackingData& td
+)
 {
-    if (!w2.valid())
+    if (!w2.valid(td))
     {
         FatalErrorIn
         (
@@ -39,7 +44,7 @@ inline bool Foam::wallNormalInfo::update(const wallNormalInfo& w2)
 
         return false;
     }
-    else if (valid())
+    else if (valid(td))
     {
         // our already set. Stop any transfer
         return false;
@@ -58,7 +63,7 @@ inline bool Foam::wallNormalInfo::update(const wallNormalInfo& w2)
 // Null constructor
 inline Foam::wallNormalInfo::wallNormalInfo()
 :
-    normal_(greatVector)
+    normal_(point::max)
 {}
 
 
@@ -84,18 +89,21 @@ inline const Foam::vector& Foam::wallNormalInfo::normal() const
 }
 
 
-inline bool Foam::wallNormalInfo::valid() const
+template<class TrackingData>
+inline bool Foam::wallNormalInfo::valid(TrackingData& td) const
 {
-    return normal_ != greatVector;
+    return normal_ != point::max;
 }
 
 
 // No geometric data so never any problem on cyclics
+template<class TrackingData>
 inline bool Foam::wallNormalInfo::sameGeometry
 (
     const polyMesh&,
     const wallNormalInfo& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 ) const
 {
     return true;
@@ -103,74 +111,97 @@ inline bool Foam::wallNormalInfo::sameGeometry
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::wallNormalInfo::leaveDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {}
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::wallNormalInfo::transform
 (
     const polyMesh&,
-    const tensor& rotTensor
+    const tensor& rotTensor,
+    TrackingData& td
 )
 {}
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::wallNormalInfo::enterDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {}
 
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool Foam::wallNormalInfo::updateCell
 (
     const polyMesh&,
     const label thisCellI,
     const label neighbourFaceI,
     const wallNormalInfo& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(neighbourWallInfo);
+    return update(neighbourWallInfo, td);
 }
 
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool Foam::wallNormalInfo::updateFace
 (
     const polyMesh&,
     const label thisFaceI,
     const label neighbourCellI,
     const wallNormalInfo& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(neighbourWallInfo);
+    return update(neighbourWallInfo, td);
 }
 
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool Foam::wallNormalInfo::updateFace
 (
     const polyMesh&,
     const label thisFaceI,
     const wallNormalInfo& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(neighbourWallInfo);
+    return update(neighbourWallInfo, td);
+}
+
+
+template <class TrackingData>
+inline bool Foam::wallNormalInfo::equal
+(
+    const wallNormalInfo& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
 }
 
 
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.C
index c9919351fc..5370d79f73 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.C
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.C
@@ -2165,6 +2165,9 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
     // refinementLevel data on seed faces
     DynamicList<refinementData> seedFacesInfo(mesh_.nFaces()/100);
 
+    // Dummy additional info for FaceCellWave
+    int dummyTrackData = 0;
+
 
     // Additional buffer layer thickness by changing initial count. Usually
     // this happens on boundary faces. Bit tricky. Use allFaceInfo to mark
@@ -2173,7 +2176,7 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
     {
         label faceI = facesToCheck[i];
 
-        if (allFaceInfo[faceI].valid())
+        if (allFaceInfo[faceI].valid(dummyTrackData))
         {
             // Can only occur if face has already gone through loop below.
             FatalErrorIn
@@ -2247,7 +2250,7 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
     forAll(faceNeighbour, faceI)
     {
         // Check if face already handled in loop above
-        if (!allFaceInfo[faceI].valid())
+        if (!allFaceInfo[faceI].valid(dummyTrackData))
         {
             label own = faceOwner[faceI];
             label nei = faceNeighbour[faceI];
@@ -2262,7 +2265,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
                     faceI,
                     own,
                     allCellInfo[own],
-                    FaceCellWave<refinementData>::propagationTol()
+                    FaceCellWave<refinementData, int>::propagationTol(),
+                    dummyTrackData
                 );
                 seedFaces.append(faceI);
                 seedFacesInfo.append(allFaceInfo[faceI]);
@@ -2275,7 +2279,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
                     faceI,
                     nei,
                     allCellInfo[nei],
-                    FaceCellWave<refinementData>::propagationTol()
+                    FaceCellWave<refinementData, int>::propagationTol(),
+                    dummyTrackData
                 );
                 seedFaces.append(faceI);
                 seedFacesInfo.append(allFaceInfo[faceI]);
@@ -2289,7 +2294,7 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
     for (label faceI = mesh_.nInternalFaces(); faceI < mesh_.nFaces(); faceI++)
     {
         // Check if face already handled in loop above
-        if (!allFaceInfo[faceI].valid())
+        if (!allFaceInfo[faceI].valid(dummyTrackData))
         {
             label own = faceOwner[faceI];
 
@@ -2301,7 +2306,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
                 faceI,
                 own,
                 allCellInfo[own],
-                FaceCellWave<refinementData>::propagationTol()
+                FaceCellWave<refinementData, int>::propagationTol(),
+                dummyTrackData
             );
             seedFaces.append(faceI);
             seedFacesInfo.append(faceData);
@@ -2310,11 +2316,12 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
 
 
     // face-cell-face transport engine
-    FaceCellWave<refinementData> levelCalc
+    FaceCellWave<refinementData, int> levelCalc
     (
         mesh_,
         allFaceInfo,
-        allCellInfo
+        allCellInfo,
+        dummyTrackData
     );
 
     while (true)
@@ -2418,7 +2425,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
                             faceI,
                             cellI,
                             cellInfo,
-                            FaceCellWave<refinementData>::propagationTol()
+                            FaceCellWave<refinementData, int>::propagationTol(),
+                            dummyTrackData
                         );
 
                         if (faceData.count() > allFaceInfo[faceI].count())
@@ -2627,6 +2635,9 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
     // Initial information about (distance to) cellLevel on all faces
     List<refinementDistanceData> allFaceInfo(mesh_.nFaces());
 
+    // Dummy additional info for FaceCellWave
+    int dummyTrackData = 0;
+
 
     // Mark cells with wanted refinement level
     forAll(cellsToRefine, i)
@@ -2643,7 +2654,7 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
     // Mark all others with existing refinement level
     forAll(allCellInfo, cellI)
     {
-        if (!allCellInfo[cellI].valid())
+        if (!allCellInfo[cellI].valid(dummyTrackData))
         {
             allCellInfo[cellI] = refinementDistanceData
             (
@@ -2660,14 +2671,13 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
     // refinementLevel data on seed faces
     DynamicList<refinementDistanceData> seedFacesInfo(mesh_.nFaces()/100);
 
-
     const pointField& cc = mesh_.cellCentres();
 
     forAll(facesToCheck, i)
     {
         label faceI = facesToCheck[i];
 
-        if (allFaceInfo[faceI].valid())
+        if (allFaceInfo[faceI].valid(dummyTrackData))
         {
             // Can only occur if face has already gone through loop below.
             FatalErrorIn
@@ -2685,7 +2695,7 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
 
         label ownLevel =
         (
-            allCellInfo[own].valid()
+            allCellInfo[own].valid(dummyTrackData)
           ? allCellInfo[own].originLevel()
           : cellLevel_[own]
         );
@@ -2709,7 +2719,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
                 faceI,
                 own,        // not used (should be nei)
                 neiData,
-                FaceCellWave<refinementDistanceData>::propagationTol()
+                FaceCellWave<refinementDistanceData, int>::propagationTol(),
+                dummyTrackData
             );
         }
         else
@@ -2718,7 +2729,7 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
 
             label neiLevel =
             (
-                allCellInfo[nei].valid()
+                allCellInfo[nei].valid(dummyTrackData)
               ? allCellInfo[nei].originLevel()
               : cellLevel_[nei]
             );
@@ -2732,7 +2743,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
                     faceI,
                     nei,
                     refinementDistanceData(level0Size, cc[nei], neiLevel+1),
-                    FaceCellWave<refinementDistanceData>::propagationTol()
+                    FaceCellWave<refinementDistanceData, int>::propagationTol(),
+                    dummyTrackData
                 );
                 allFaceInfo[faceI].updateFace
                 (
@@ -2740,7 +2752,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
                     faceI,
                     own,
                     refinementDistanceData(level0Size, cc[own], ownLevel+1),
-                    FaceCellWave<refinementDistanceData>::propagationTol()
+                    FaceCellWave<refinementDistanceData, int>::propagationTol(),
+                    dummyTrackData
                 );
             }
             else
@@ -2752,7 +2765,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
                     faceI,
                     nei,
                     refinementDistanceData(level0Size, cc[nei], neiLevel),
-                    FaceCellWave<refinementDistanceData>::propagationTol()
+                    FaceCellWave<refinementDistanceData, int>::propagationTol(),
+                    dummyTrackData
                 );
                 allFaceInfo[faceI].updateFace
                 (
@@ -2760,7 +2774,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
                     faceI,
                     own,
                     refinementDistanceData(level0Size, cc[own], ownLevel),
-                    FaceCellWave<refinementDistanceData>::propagationTol()
+                    FaceCellWave<refinementDistanceData, int>::propagationTol(),
+                    dummyTrackData
                 );
             }
         }
@@ -2775,13 +2790,13 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
     forAll(faceNeighbour, faceI)
     {
         // Check if face already handled in loop above
-        if (!allFaceInfo[faceI].valid())
+        if (!allFaceInfo[faceI].valid(dummyTrackData))
         {
             label own = faceOwner[faceI];
 
             label ownLevel =
             (
-                allCellInfo[own].valid()
+                allCellInfo[own].valid(dummyTrackData)
               ? allCellInfo[own].originLevel()
               : cellLevel_[own]
             );
@@ -2790,7 +2805,7 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
 
             label neiLevel =
             (
-                allCellInfo[nei].valid()
+                allCellInfo[nei].valid(dummyTrackData)
               ? allCellInfo[nei].originLevel()
               : cellLevel_[nei]
             );
@@ -2805,7 +2820,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
                     faceI,
                     own,
                     refinementDistanceData(level0Size, cc[own], ownLevel),
-                    FaceCellWave<refinementDistanceData>::propagationTol()
+                    FaceCellWave<refinementDistanceData, int>::propagationTol(),
+                    dummyTrackData
                 );
                 seedFacesInfo.append(allFaceInfo[faceI]);
             }
@@ -2818,7 +2834,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
                     faceI,
                     nei,
                     refinementDistanceData(level0Size, cc[nei], neiLevel),
-                    FaceCellWave<refinementDistanceData>::propagationTol()
+                    FaceCellWave<refinementDistanceData, int>::propagationTol(),
+                    dummyTrackData
                 );
                 seedFacesInfo.append(allFaceInfo[faceI]);
             }
@@ -2829,14 +2846,15 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement2
     seedFacesInfo.shrink();
 
     // face-cell-face transport engine
-    FaceCellWave<refinementDistanceData> levelCalc
+    FaceCellWave<refinementDistanceData, int> levelCalc
     (
         mesh_,
         seedFaces,
         seedFacesInfo,
         allFaceInfo,
         allCellInfo,
-        mesh_.globalData().nTotalCells()+1
+        mesh_.globalData().nTotalCells()+1,
+        dummyTrackData
     );
 
 
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementData.H b/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementData.H
index 5473ce04d2..171de17109 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementData.H
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementData.H
@@ -26,7 +26,7 @@ Class
 
 Description
     Transfers refinement levels such that slow transition between levels is
-    maintained. Used in meshWave.
+    maintained. Used in FaceCellWave.
 
 SourceFiles
     refinementDataI.H
@@ -108,71 +108,91 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value.
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData& td) const;
 
             //- Check for identical geometrical data. Used for cyclics checking.
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const refinementData&,
-                const scalar
+                const scalar,
+                TrackingData& td
             ) const;
 
             //- Convert any absolute coordinates into relative to (patch)face
             //  centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor&
+                const tensor&,
+                TrackingData& td
             );
 
             //- Influence of neighbouring face.
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const refinementData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const refinementData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const refinementData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const refinementData&, TrackingData& td) const;
+
+
     // Member Operators
 
         // Needed for List IO
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDataI.H b/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDataI.H
index c8ddde6223..0325942ee0 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDataI.H
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDataI.H
@@ -23,10 +23,6 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Null constructor
@@ -51,18 +47,21 @@ inline Foam::refinementData::refinementData
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline bool Foam::refinementData::valid() const
+template<class TrackingData>
+inline bool Foam::refinementData::valid(TrackingData& td) const
 {
     return count_ != -1;
 }
 
 
 // No geometric data so never any problem on cyclics
+template<class TrackingData>
 inline bool Foam::refinementData::sameGeometry
 (
     const polyMesh&,
     const refinementData&,
-    const scalar
+    const scalar,
+    TrackingData& td
 ) const
 {
     return true;
@@ -70,47 +69,55 @@ inline bool Foam::refinementData::sameGeometry
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::refinementData::leaveDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {}
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::refinementData::transform
 (
     const polyMesh&,
-    const tensor& rotTensor
+    const tensor& rotTensor,
+    TrackingData& td
 )
 {}
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::refinementData::enterDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {}
 
 
 // Update cell with neighbouring face information
+template<class TrackingData>
 inline bool Foam::refinementData::updateCell
 (
     const polyMesh&,
     const label thisCellI,
     const label neighbourFaceI,
     const refinementData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    if (!valid())
+    if (!valid(td))
     {
         FatalErrorIn("refinementData::updateCell") << "problem"
             << abort(FatalError);
@@ -161,17 +168,19 @@ inline bool Foam::refinementData::updateCell
 
 
 // Update face with neighbouring cell information
+template<class TrackingData>
 inline bool Foam::refinementData::updateFace
 (
     const polyMesh&,
     const label thisFaceI,
     const label neighbourCellI,
     const refinementData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     // From cell to its faces.
-    if (!valid())
+    if (!valid(td))
     {
         refinementCount_ = neighbourInfo.refinementCount();
         count_ = neighbourInfo.count();
@@ -194,16 +203,18 @@ inline bool Foam::refinementData::updateFace
 
 
 // Update face with coupled face information
+template<class TrackingData>
 inline bool Foam::refinementData::updateFace
 (
     const polyMesh&,
     const label thisFaceI,
     const refinementData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     // From face to face (e.g. coupled faces)
-    if (!valid())
+    if (!valid(td))
     {
         refinementCount_ = neighbourInfo.refinementCount();
         count_ = neighbourInfo.count();
@@ -225,6 +236,17 @@ inline bool Foam::refinementData::updateFace
 }
 
 
+template <class TrackingData>
+inline bool Foam::refinementData::equal
+(
+    const refinementData& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 inline bool Foam::refinementData::operator==(const Foam::refinementData& rhs)
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDistanceData.H b/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDistanceData.H
index 39a6d3721b..d838ef7178 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDistanceData.H
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDistanceData.H
@@ -68,11 +68,13 @@ class refinementDistanceData
     // Private Member Functions
 
         //- Updates with neighbouring data. Returns true if something changed.
+        template<class TrackingData>
         inline bool update
         (
             const point&,
             const refinementDistanceData& neighbourInfo,
-            const scalar tol
+            const scalar tol,
+            TrackingData&
         );
 
 public:
@@ -137,71 +139,94 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value.
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData&) const;
 
             //- Check for identical geometrical data. Used for cyclics checking.
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const refinementDistanceData&,
-                const scalar
+                const scalar,
+                TrackingData&
             ) const;
 
             //- Convert any absolute coordinates into relative to (patch)face
             //  centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData&
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData&
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor&
+                const tensor&,
+                TrackingData&
             );
 
             //- Influence of neighbouring face.
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const refinementDistanceData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData&
             );
 
             //- Influence of neighbouring cell.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const refinementDistanceData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData&
             );
 
             //- Influence of different value on same face.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const refinementDistanceData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData&
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal
+            (
+                const refinementDistanceData&,
+                TrackingData&
+            ) const;
+
     // Member Operators
 
         // Needed for List IO
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDistanceDataI.H b/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDistanceDataI.H
index d5d9afb45e..a501d293b9 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDistanceDataI.H
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/refinementDistanceDataI.H
@@ -57,16 +57,18 @@ inline Foam::label Foam::refinementDistanceData::wantedLevel(const point& pt)
 }
 
 
+template<class TrackingData>
 inline bool Foam::refinementDistanceData::update
 (
     const point& pos,
     const refinementDistanceData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    if (!valid())
+    if (!valid(td))
     {
-        if (!neighbourInfo.valid())
+        if (!neighbourInfo.valid(td))
         {
             FatalErrorIn("refinementDistanceData::update(..)")
                 << "problem" << abort(FatalError);
@@ -142,40 +144,47 @@ inline Foam::refinementDistanceData::refinementDistanceData
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline bool Foam::refinementDistanceData::valid() const
+template<class TrackingData>
+inline bool Foam::refinementDistanceData::valid(TrackingData& td) const
 {
     return level0Size_ != -1;
 }
 
 
 // No geometric data so never any problem on cyclics
+template<class TrackingData>
 inline bool Foam::refinementDistanceData::sameGeometry
 (
     const polyMesh&,
     const refinementDistanceData&,
-    const scalar
+    const scalar,
+    TrackingData& td
 ) const
 {
     return true;
 }
 
 
+template<class TrackingData>
 inline void Foam::refinementDistanceData::leaveDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {
     origin_ -= faceCentre;
 }
 
 
+template<class TrackingData>
 inline void Foam::refinementDistanceData::transform
 (
     const polyMesh&,
-    const tensor& rotTensor
+    const tensor& rotTensor,
+    TrackingData& td
 )
 {
     origin_ = Foam::transform(rotTensor, origin_);
@@ -183,12 +192,14 @@ inline void Foam::refinementDistanceData::transform
 
 
 // Update absolute geometric quantities.
+template<class TrackingData>
 inline void Foam::refinementDistanceData::enterDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {
     // back to absolute form
@@ -197,49 +208,80 @@ inline void Foam::refinementDistanceData::enterDomain
 
 
 // Update cell with neighbouring face information
+template<class TrackingData>
 inline bool Foam::refinementDistanceData::updateCell
 (
     const polyMesh& mesh,
     const label thisCellI,
     const label neighbourFaceI,
     const refinementDistanceData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const point& pos = mesh.cellCentres()[thisCellI];
 
-    return update(pos, neighbourInfo, tol);
+    return update(pos, neighbourInfo, tol, td);
 }
 
 
 // Update face with neighbouring cell information
+template<class TrackingData>
 inline bool Foam::refinementDistanceData::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const label neighbourCellI,
     const refinementDistanceData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const point& pos = mesh.faceCentres()[thisFaceI];
 
-    return update(pos, neighbourInfo, tol);
+    return update(pos, neighbourInfo, tol, td);
 }
 
 
 // Update face with coupled face information
+template<class TrackingData>
 inline bool Foam::refinementDistanceData::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const refinementDistanceData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const point& pos = mesh.faceCentres()[thisFaceI];
 
-    return update(pos, neighbourInfo, tol);
+    return update(pos, neighbourInfo, tol, td);
+}
+
+
+template <class TrackingData>
+inline bool Foam::refinementDistanceData::equal
+(
+    const refinementDistanceData& rhs,
+    TrackingData& td
+) const
+{
+    if (!valid(td))
+    {
+        if (!rhs.valid(td))
+        {
+            return true;
+        }
+        else
+        {
+            return false;
+        }
+    }
+    else
+    {
+        return operator==(rhs);
+    }
 }
 
 
@@ -251,24 +293,10 @@ inline bool Foam::refinementDistanceData::operator==
 )
  const
 {
-    if (!valid())
-    {
-        if (!rhs.valid())
-        {
-            return true;
-        }
-        else
-        {
-            return false;
-        }
-    }
-    else
-    {
-        return
-            level0Size_ == rhs.level0Size_
-         && origin_ == rhs.origin_
-         && originLevel_ == rhs.originLevel_;
-    }
+    return
+        level0Size_ == rhs.level0Size_
+     && origin_ == rhs.origin_
+     && originLevel_ == rhs.originLevel_;
 }
 
 
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.C b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.C
index dd803cec72..57f48b9935 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.C
+++ b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/fvcSmooth.C
@@ -31,8 +31,6 @@ License
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-Foam::scalar Foam::smoothData::maxRatio = 1.0;
-
 void Foam::fvc::smooth
 (
     volScalarField& field,
@@ -40,7 +38,7 @@ void Foam::fvc::smooth
 )
 {
     const fvMesh& mesh = field.mesh();
-    smoothData::maxRatio = 1 + coeff;
+    scalar maxRatio = 1 + coeff;
 
     DynamicList<label> changedFaces(mesh.nFaces()/100 + 100);
     DynamicList<smoothData> changedFacesInfo(changedFaces.size());
@@ -54,12 +52,12 @@ void Foam::fvc::smooth
         const label nbr = neighbour[facei];
 
         // Check if owner value much larger than neighbour value or vice versa
-        if (field[own] > smoothData::maxRatio*field[nbr])
+        if (field[own] > maxRatio*field[nbr])
         {
             changedFaces.append(facei);
             changedFacesInfo.append(smoothData(field[own]));
         }
-        else if (field[nbr] > smoothData::maxRatio*field[own])
+        else if (field[nbr] > maxRatio*field[own])
         {
             changedFaces.append(facei);
             changedFacesInfo.append(smoothData(field[nbr]));
@@ -98,16 +96,19 @@ void Foam::fvc::smooth
     // Set initial field on faces
     List<smoothData> faceData(mesh.nFaces());
 
+    smoothData::trackData td;
+    td.maxRatio = maxRatio;
 
     // Propagate information over whole domain
-    FaceCellWave<smoothData > smoothData
+    FaceCellWave<smoothData, smoothData::trackData> smoothData
     (
         mesh,
         changedFaces,
         changedFacesInfo,
         faceData,
         cellData,
-        mesh.globalData().nTotalCells()    // max iterations
+        mesh.globalData().nTotalCells(),   // max iterations
+        td
     );
 
     forAll(field, celli)
@@ -130,7 +131,6 @@ void Foam::fvc::spread
 )
 {
     const fvMesh& mesh = field.mesh();
-    smoothData::maxRatio = 1;
 
     DynamicList<label> changedFaces(mesh.nFaces()/100 + 100);
     DynamicList<smoothData> changedFacesInfo(changedFaces.size());
@@ -208,12 +208,16 @@ void Foam::fvc::spread
     changedFaces.shrink();
     changedFacesInfo.shrink();
 
+    smoothData::trackData td;
+    td.maxRatio = 1.0;
+
     // Propagate information over whole domain
-    FaceCellWave<smoothData> smoothData
+    FaceCellWave<smoothData, smoothData::trackData> smoothData
     (
         mesh,
         faceData,
-        cellData
+        cellData,
+        td
     );
 
     smoothData.setFaceInfo(changedFaces, changedFacesInfo);
@@ -311,7 +315,7 @@ void Foam::fvc::sweep
 
     forAll(field, celli)
     {
-        if (cellData[celli].valid())
+        if (cellData[celli].valid(sweepData.data()))
         {
             field[celli] = max(field[celli], cellData[celli].value());
         }
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/smoothData.H b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/smoothData.H
index d027f9f2a2..f51124e920 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/smoothData.H
+++ b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/smoothData.H
@@ -47,27 +47,39 @@ namespace Foam
 
 class smoothData
 {
+
+public:
+
+    //- Class used to pass additional data in
+    class trackData
+    {
+    public:
+
+        //- Cut off distance
+        scalar maxRatio;
+    };
+
+
+private:
+
     scalar value_;
 
     // Private Member Functions
 
         //- Update - gets information from neighbouring face/cell and
         //  uses this to update itself (if necessary) and return true
+        template<class TrackingData>
         inline bool update
         (
             const smoothData& svf,
             const scalar scale,
-            const scalar tol
+            const scalar tol,
+            TrackingData& td
         );
 
 
 public:
 
-    // Static data members
-
-        //- Field fraction
-        static scalar maxRatio;
-
 
     // Constructors
 
@@ -93,72 +105,91 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData& td) const;
 
             //- Check for identical geometrical data
             //  Used for cyclics checking
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const smoothData&,
-                const scalar
+                const scalar,
+                TrackingData& td
             ) const;
 
             //- Convert any absolute coordinates into relative to
             //  (patch)face centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor&
+                const tensor&,
+                TrackingData& td
             );
 
             //- Influence of neighbouring face
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const smoothData& svf,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const smoothData& svf,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const smoothData& svf,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const smoothData&, TrackingData& td) const;
+
 
         // Member Operators
 
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/smoothDataI.H b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/smoothDataI.H
index 4799ede9ec..4694f593aa 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/smoothDataI.H
+++ b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/smoothDataI.H
@@ -25,14 +25,16 @@ License
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+template<class TrackingData>
 inline bool Foam::smoothData::update
 (
     const smoothData& svf,
     const scalar scale,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    if (!valid() || (value_ < VSMALL))
+    if (!valid(td) || (value_ < VSMALL))
     {
         // My value not set - take over neighbour
         value_ = svf.value()/scale;
@@ -73,90 +75,116 @@ inline Foam::smoothData::smoothData(const scalar value)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline bool Foam::smoothData::valid() const
+template<class TrackingData>
+inline bool Foam::smoothData::valid(TrackingData& td) const
 {
     return value_ > -SMALL;
 }
 
 
+template<class TrackingData>
 inline bool Foam::smoothData::sameGeometry
 (
     const polyMesh&,
     const smoothData&,
-    const scalar
+    const scalar,
+    TrackingData& td
 ) const
 {
     return true;
 }
 
 
+template<class TrackingData>
 inline void Foam::smoothData::leaveDomain
 (
     const polyMesh&,
     const polyPatch&,
     const label,
-    const point&
+    const point&,
+    TrackingData& td
 )
 {}
 
 
+template<class TrackingData>
 inline void Foam::smoothData::transform
 (
     const polyMesh&,
-    const tensor&
+    const tensor&,
+    TrackingData& td
 )
 {}
 
 
+template<class TrackingData>
 inline void Foam::smoothData::enterDomain
 (
     const polyMesh&,
     const polyPatch&,
     const label,
-    const point&
+    const point&,
+    TrackingData& td
 )
 {}
 
 
+template<class TrackingData>
 inline bool Foam::smoothData::updateCell
 (
     const polyMesh&,
     const label,
     const label,
     const smoothData& svf,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     // Take over info from face if more than deltaRatio larger
-    return update(svf, maxRatio, tol);
+    return update(svf, td.maxRatio, tol, td);
 }
 
 
+template<class TrackingData>
 inline bool Foam::smoothData::updateFace
 (
     const polyMesh&,
     const label,
     const label,
     const smoothData& svf,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     // Take over information from cell without any scaling (scale = 1.0)
-    return update(svf, 1.0, tol);
+    return update(svf, 1.0, tol, td);
 }
 
 
 // Update this (face) with coupled face information.
+template<class TrackingData>
 inline bool Foam::smoothData::updateFace
 (
     const polyMesh&,
     const label,
     const smoothData& svf,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     // Take over information from coupled face without any scaling (scale = 1.0)
-    return update(svf, 1.0, tol);
+    return update(svf, 1.0, tol, td);
+}
+
+
+template <class TrackingData>
+inline bool Foam::smoothData::equal
+(
+    const smoothData& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
 }
 
 
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/sweepData.H b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/sweepData.H
index 1bfec88fdf..ff64713e0b 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/sweepData.H
+++ b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/sweepData.H
@@ -54,11 +54,13 @@ class sweepData
 
         //- Update - gets information from neighbouring face/cell and
         //  uses this to update itself (if necessary) and return true
+        template<class TrackingData>
         inline bool update
         (
             const sweepData& svf,
             const point& position,
-            const scalar tol
+            const scalar tol,
+            TrackingData& td
         );
 
 
@@ -94,72 +96,91 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData& td) const;
 
             //- Check for identical geometrical data
             //  Used for cyclics checking
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const sweepData&,
-                const scalar
+                const scalar,
+                TrackingData& td
             ) const;
 
             //- Convert any absolute coordinates into relative to
             //  (patch)face centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor&
+                const tensor&,
+                TrackingData& td
             );
 
             //- Influence of neighbouring face
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const sweepData& svf,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const sweepData& svf,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const sweepData& svf,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const sweepData&, TrackingData& td) const;
+
 
         // Member Operators
 
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/sweepDataI.H b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/sweepDataI.H
index a042886085..0c9f8e3e60 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcSmooth/sweepDataI.H
+++ b/src/finiteVolume/finiteVolume/fvc/fvcSmooth/sweepDataI.H
@@ -27,14 +27,16 @@ License
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+template<class TrackingData>
 inline bool Foam::sweepData::update
 (
     const sweepData& svf,
     const point& position,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    if (!valid())
+    if (!valid(td))
     {
         operator=(svf);
         return true;
@@ -85,51 +87,60 @@ inline Foam::sweepData::sweepData(const scalar value, const point& origin)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline bool Foam::sweepData::valid() const
+template<class TrackingData>
+inline bool Foam::sweepData::valid(TrackingData& td) const
 {
     return value_ > -SMALL;
 }
 
 
+template<class TrackingData>
 inline bool Foam::sweepData::sameGeometry
 (
     const polyMesh&,
     const sweepData&,
-    const scalar
+    const scalar,
+    TrackingData& td
 ) const
 {
     return true;
 }
 
 
+template<class TrackingData>
 inline void Foam::sweepData::leaveDomain
 (
     const polyMesh&,
     const polyPatch&,
     const label,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {
     origin_ -= faceCentre;
 }
 
 
+template<class TrackingData>
 inline void Foam::sweepData::transform
 (
     const polyMesh&,
-    const tensor& rotTensor
+    const tensor& rotTensor,
+    TrackingData& td
 )
 {
     origin_ = Foam::transform(rotTensor, origin_);
 }
 
 
+template<class TrackingData>
 inline void Foam::sweepData::enterDomain
 (
     const polyMesh&,
     const polyPatch&,
     const label,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {
     // back to absolute form
@@ -137,42 +148,59 @@ inline void Foam::sweepData::enterDomain
 }
 
 
+template<class TrackingData>
 inline bool Foam::sweepData::updateCell
 (
     const polyMesh& mesh,
     const label thisCellI,
     const label,
     const sweepData& svf,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(svf, mesh.cellCentres()[thisCellI], tol);
+    return update(svf, mesh.cellCentres()[thisCellI], tol, td);
 }
 
 
+template<class TrackingData>
 inline bool Foam::sweepData::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const label,
     const sweepData& svf,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(svf, mesh.faceCentres()[thisFaceI], tol);
+    return update(svf, mesh.faceCentres()[thisFaceI], tol, td);
 }
 
 
 // Update this (face) with coupled face information.
+template<class TrackingData>
 inline bool Foam::sweepData::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const sweepData& svf,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(svf, mesh.faceCentres()[thisFaceI], tol);
+    return update(svf, mesh.faceCentres()[thisFaceI], tol, td);
+}
+
+
+template <class TrackingData>
+inline bool Foam::sweepData::equal
+(
+    const sweepData& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
 }
 
 
diff --git a/src/finiteVolume/fvMesh/wallDist/wallPointYPlus/wallPointYPlus.H b/src/finiteVolume/fvMesh/wallDist/wallPointYPlus/wallPointYPlus.H
index 3c2ccce7cd..2ca7e7197c 100644
--- a/src/finiteVolume/fvMesh/wallDist/wallPointYPlus/wallPointYPlus.H
+++ b/src/finiteVolume/fvMesh/wallDist/wallPointYPlus/wallPointYPlus.H
@@ -32,6 +32,9 @@ Description
     the damping function becomes 1, since y gets initialized to GREAT and
     yStar to 1.
 
+    Note: should feed the additional argument (yPlusCutoff) through as a
+    template argument into FaceCellWave
+
 SourceFiles
     wallPointYPlusI.H
     wallPointYPlus.C
@@ -58,16 +61,20 @@ class wallPointYPlus
     public wallPointData<scalar>
 {
 
+private:
+
     // Private Member Functions
 
         //- Evaluate distance to point. Update distSqr, origin from whomever
         //  is nearer pt. Return true if w2 is closer to point,
         //  false otherwise.
+        template<class TrackingData>
         inline bool update
         (
             const point&,
             const wallPointYPlus& w2,
-            const scalar tol
+            const scalar tol,
+            TrackingData& td
         );
 
 public:
@@ -94,39 +101,45 @@ public:
 
     // Member Functions
 
-        // Needed by meshWave
+        // Needed by FaceCellWave
 
             //- Influence of neighbouring face.
             //  Calls update(...) with cellCentre of cellI
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh& mesh,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const wallPointYPlus& neighbourWallInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell.
             //  Calls update(...) with faceCentre of faceI
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh& mesh,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const wallPointYPlus& neighbourWallInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face.
             //  Merge new and old info.
             //  Calls update(...) with faceCentre of faceI
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh& mesh,
                 const label thisFaceI,
                 const wallPointYPlus& neighbourWallInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 };
 
diff --git a/src/finiteVolume/fvMesh/wallDist/wallPointYPlus/wallPointYPlusI.H b/src/finiteVolume/fvMesh/wallDist/wallPointYPlus/wallPointYPlusI.H
index a9bab14083..6e9784264a 100644
--- a/src/finiteVolume/fvMesh/wallDist/wallPointYPlus/wallPointYPlusI.H
+++ b/src/finiteVolume/fvMesh/wallDist/wallPointYPlus/wallPointYPlusI.H
@@ -31,11 +31,13 @@ namespace Foam
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool wallPointYPlus::update
 (
     const point& pt,
     const wallPointYPlus& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     scalar dist2 = magSqr(pt - w2.origin());
@@ -103,13 +105,15 @@ inline wallPointYPlus::wallPointYPlus
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool wallPointYPlus::updateCell
 (
     const polyMesh& mesh,
     const label thisCellI,
     const label neighbourFaceI,
     const wallPointYPlus& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const vectorField& cellCentres = mesh.primitiveMesh::cellCentres();
@@ -118,19 +122,22 @@ inline bool wallPointYPlus::updateCell
     (
         cellCentres[thisCellI],
         neighbourWallInfo,
-        tol
+        tol,
+        td
     );
 }
 
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool wallPointYPlus::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const label neighbourCellI,
     const wallPointYPlus& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const vectorField& faceCentres = mesh.faceCentres();
@@ -139,18 +146,21 @@ inline bool wallPointYPlus::updateFace
     (
         faceCentres[thisFaceI],
         neighbourWallInfo,
-        tol
+        tol,
+        td
     );
 }
 
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool wallPointYPlus::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const wallPointYPlus& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const vectorField& faceCentres = mesh.faceCentres();
@@ -159,7 +169,8 @@ inline bool wallPointYPlus::updateFace
     (
         faceCentres[thisFaceI],
         neighbourWallInfo,
-        tol
+        tol,
+        td
     );
 }
 
diff --git a/src/fvMotionSolver/motionDiffusivity/inverseFaceDistance/inverseFaceDistanceDiffusivity.C b/src/fvMotionSolver/motionDiffusivity/inverseFaceDistance/inverseFaceDistanceDiffusivity.C
index 62587d26ff..6464f87341 100644
--- a/src/fvMotionSolver/motionDiffusivity/inverseFaceDistance/inverseFaceDistanceDiffusivity.C
+++ b/src/fvMotionSolver/motionDiffusivity/inverseFaceDistance/inverseFaceDistanceDiffusivity.C
@@ -115,7 +115,7 @@ void Foam::inverseFaceDistanceDiffusivity::correct()
         mesh,
         changedFaces,
         faceDist,
-        mesh.globalData().nTotalCells()+1 // max iterations
+        mesh.globalData().nTotalCells()+1   // max iterations
     );
 
     const List<wallPoint>& faceInfo = waveInfo.allFaceInfo();
diff --git a/src/meshTools/cellClassification/cellInfo.H b/src/meshTools/cellClassification/cellInfo.H
index 4fbf65176c..a51e255f4c 100644
--- a/src/meshTools/cellClassification/cellInfo.H
+++ b/src/meshTools/cellClassification/cellInfo.H
@@ -63,13 +63,15 @@ class cellInfo
         //- Update current cell/face type with neighbouring
         //  type. Return true if information needs to propagate,
         //  false otherwise.
+        template<class TrackingData>
         inline bool update
         (
             const cellInfo& w2,
             const label thisFaceI,
             const label thisCellI,
             const label neighbourFaceI,
-            const label neighbourCellI
+            const label neighbourCellI,
+            TrackingData& td
         );
 
 public:
@@ -100,71 +102,90 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value.
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData& td) const;
 
             //- Check for identical geometrical data. Used for cyclics checking.
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const cellInfo&,
-                const scalar
+                const scalar,
+                TrackingData& td
             ) const;
 
             //- Convert any absolute coordinates into relative to (patch)face
             //  centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch& patch,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch& patch,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor& rotTensor
+                const tensor& rotTensor,
+                TrackingData& td
             );
 
             //- Influence of neighbouring face.
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const cellInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const cellInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const cellInfo& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const cellInfo&, TrackingData& td) const;
+
     // Member Operators
 
         //Note: Used to determine whether to call update.
diff --git a/src/meshTools/cellClassification/cellInfoI.H b/src/meshTools/cellClassification/cellInfoI.H
index 8bdc6c8aac..15971fdb69 100644
--- a/src/meshTools/cellClassification/cellInfoI.H
+++ b/src/meshTools/cellClassification/cellInfoI.H
@@ -29,13 +29,15 @@ License
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Update this with w2 information
+template<class TrackingData>
 inline bool Foam::cellInfo::update
 (
     const cellInfo& w2,
     const label thisFaceI,
     const label thisCellI,
     const label neighbourFaceI,
-    const label neighbourCellI
+    const label neighbourCellI,
+    TrackingData& td
 )
 {
     if
@@ -113,18 +115,21 @@ inline Foam::cellInfo::cellInfo(const cellInfo& w2)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline bool Foam::cellInfo::valid() const
+template<class TrackingData>
+inline bool Foam::cellInfo::valid(TrackingData& td) const
 {
     return type_ != cellClassification::NOTSET;
 }
 
 
 // No geometric data so never any problem on cyclics
+template<class TrackingData>
 inline bool Foam::cellInfo::sameGeometry
 (
     const polyMesh&,
     const cellInfo& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
  const
 {
@@ -133,44 +138,52 @@ inline bool Foam::cellInfo::sameGeometry
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::cellInfo::leaveDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {}
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::cellInfo::transform
 (
     const polyMesh&,
-    const tensor& rotTensor
+    const tensor& rotTensor,
+    TrackingData& td
 )
 {}
 
 
 // No geometric data.
+template<class TrackingData>
 inline void Foam::cellInfo::enterDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {}
 
 
 // Update this with neighbour information
+template<class TrackingData>
 inline bool Foam::cellInfo::updateCell
 (
     const polyMesh&,
     const label thisCellI,
     const label neighbourFaceI,
     const cellInfo& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     return update
@@ -179,19 +192,22 @@ inline bool Foam::cellInfo::updateCell
         -1,
         thisCellI,
         neighbourFaceI,
-        -1
+        -1,
+        td
     );
 }
 
 
 // Update this with neighbour information
+template<class TrackingData>
 inline bool Foam::cellInfo::updateFace
 (
     const polyMesh&,
     const label thisFaceI,
     const label neighbourCellI,
     const cellInfo& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     return update
@@ -200,17 +216,20 @@ inline bool Foam::cellInfo::updateFace
         thisFaceI,
         -1,
         -1,
-        neighbourCellI
+        neighbourCellI,
+        td
     );
 }
 
 // Update this with neighbour information
+template<class TrackingData>
 inline bool Foam::cellInfo::updateFace
 (
     const polyMesh&,
     const label thisFaceI,
     const cellInfo& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     return update
@@ -219,11 +238,23 @@ inline bool Foam::cellInfo::updateFace
         thisFaceI,
         -1,
         -1,
-        -1
+        -1,
+        td
     );
 }
 
 
+template <class TrackingData>
+inline bool Foam::cellInfo::equal
+(
+    const cellInfo& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 inline bool Foam::cellInfo::operator==(const Foam::cellInfo& rhs) const
diff --git a/src/meshTools/cellDist/patchWave/patchDataWave.C b/src/meshTools/cellDist/patchWave/patchDataWave.C
index 988b8a0be5..c65f59a0f3 100644
--- a/src/meshTools/cellDist/patchWave/patchDataWave.C
+++ b/src/meshTools/cellDist/patchWave/patchDataWave.C
@@ -24,7 +24,6 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "patchDataWave.H"
-#include "MeshWave.H"
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
@@ -93,7 +92,7 @@ Foam::label Foam::patchDataWave<TransferType>::getValues
 
         scalar dist = wpn.distSqr();
 
-        if (cellInfo[cellI].valid())
+        if (cellInfo[cellI].valid(waveInfo.data()))
         {
             distance_[cellI] = Foam::sqrt(dist);
 
@@ -138,7 +137,7 @@ Foam::label Foam::patchDataWave<TransferType>::getValues
 
             scalar dist = faceInfo[meshFaceI].distSqr();
 
-            if (faceInfo[meshFaceI].valid())
+            if (faceInfo[meshFaceI].valid(waveInfo.data()))
             {
                 // Adding SMALL to avoid problems with /0 in the turbulence
                 // models
diff --git a/src/meshTools/cellDist/patchWave/patchDataWave.H b/src/meshTools/cellDist/patchWave/patchDataWave.H
index ec562f7229..714e3b745b 100644
--- a/src/meshTools/cellDist/patchWave/patchDataWave.H
+++ b/src/meshTools/cellDist/patchWave/patchDataWave.H
@@ -45,6 +45,7 @@ SourceFiles
 #include "cellDistFuncs.H"
 #include "FieldField.H"
 #include "UPtrList.H"
+#include "MeshWave.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -54,7 +55,7 @@ namespace Foam
 // Forward declaration of classes
 class polyMesh;
 class wallPoint;
-template<class Type> class MeshWave;
+//template<class Type, class TrackingData> class MeshWave;
 
 /*---------------------------------------------------------------------------*\
                            Class patchDataWave Declaration
diff --git a/src/meshTools/cellDist/patchWave/patchWave.C b/src/meshTools/cellDist/patchWave/patchWave.C
index 826430b003..7d104a8382 100644
--- a/src/meshTools/cellDist/patchWave/patchWave.C
+++ b/src/meshTools/cellDist/patchWave/patchWave.C
@@ -26,7 +26,6 @@ License
 #include "patchWave.H"
 #include "polyMesh.H"
 #include "wallPoint.H"
-#include "MeshWave.H"
 #include "globalMeshData.H"
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
@@ -82,7 +81,7 @@ Foam::label Foam::patchWave::getValues(const MeshWave<wallPoint>& waveInfo)
     {
         scalar dist = cellInfo[cellI].distSqr();
 
-        if (cellInfo[cellI].valid())
+        if (cellInfo[cellI].valid(waveInfo.data()))
         {
             distance_[cellI] = Foam::sqrt(dist);
         }
@@ -112,7 +111,7 @@ Foam::label Foam::patchWave::getValues(const MeshWave<wallPoint>& waveInfo)
 
             scalar dist = faceInfo[meshFaceI].distSqr();
 
-            if (faceInfo[meshFaceI].valid())
+            if (faceInfo[meshFaceI].valid(waveInfo.data()))
             {
                 // Adding SMALL to avoid problems with /0 in the turbulence
                 // models
diff --git a/src/meshTools/cellDist/patchWave/patchWave.H b/src/meshTools/cellDist/patchWave/patchWave.H
index b9e5dc50cd..b3f718cc52 100644
--- a/src/meshTools/cellDist/patchWave/patchWave.H
+++ b/src/meshTools/cellDist/patchWave/patchWave.H
@@ -39,7 +39,7 @@ SourceFiles
 
 #include "cellDistFuncs.H"
 #include "FieldField.H"
-
+#include "MeshWave.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -49,7 +49,6 @@ namespace Foam
 // Forward declaration of classes
 class polyMesh;
 class wallPoint;
-template<class Type> class MeshWave;
 
 /*---------------------------------------------------------------------------*\
                            Class patchWave Declaration
diff --git a/src/meshTools/cellDist/wallPoint/wallPoint.H b/src/meshTools/cellDist/wallPoint/wallPoint.H
index 5f2386828d..8215795896 100644
--- a/src/meshTools/cellDist/wallPoint/wallPoint.H
+++ b/src/meshTools/cellDist/wallPoint/wallPoint.H
@@ -76,11 +76,13 @@ class wallPoint
         //- Evaluate distance to point. Update distSqr, origin from whomever
         //  is nearer pt. Return true if w2 is closer to point,
         //  false otherwise.
+        template<class TrackingData>
         inline bool update
         (
             const point&,
             const wallPoint& w2,
-            const scalar tol
+            const scalar tol,
+            TrackingData& td
         );
 
 
@@ -122,71 +124,90 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value.
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData& td) const;
 
             //- Check for identical geometrical data. Used for cyclics checking.
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const wallPoint&,
-                const scalar
+                const scalar,
+                TrackingData& td
             ) const;
 
             //- Convert any absolute coordinates into relative to (patch)face
             //  centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor&
+                const tensor&,
+                TrackingData& td
             );
 
             //- Influence of neighbouring face.
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const wallPoint& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const wallPoint& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const wallPoint& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const wallPoint&, TrackingData& td) const;
+
 
     // Member Operators
 
diff --git a/src/meshTools/cellDist/wallPoint/wallPointData.C b/src/meshTools/cellDist/wallPoint/wallPointData.C
index 246cba03ab..0de93818aa 100644
--- a/src/meshTools/cellDist/wallPoint/wallPointData.C
+++ b/src/meshTools/cellDist/wallPoint/wallPointData.C
@@ -24,9 +24,6 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "wallPointData.H"
-#include "point.H"
-#include "scalar.H"
-
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -47,6 +44,7 @@ Ostream& operator<<
     return os << wDist.data();
 }
 
+
 template <class Type>
 Istream& operator>>
 (
@@ -59,6 +57,7 @@ Istream& operator>>
     return is >> wDist.data_;
 }
 
+
 // ************************************************************************* //
 
 } // End namespace Foam
diff --git a/src/meshTools/cellDist/wallPoint/wallPointData.H b/src/meshTools/cellDist/wallPoint/wallPointData.H
index 3f8d75f2b8..a826f0003a 100644
--- a/src/meshTools/cellDist/wallPoint/wallPointData.H
+++ b/src/meshTools/cellDist/wallPoint/wallPointData.H
@@ -74,11 +74,13 @@ class wallPointData
         //- Evaluate distance to point. Update distSqr, origin from whomever
         //  is nearer pt. Return true if w2 is closer to point,
         //  false otherwise.
+        template<class TrackingData>
         inline bool update
         (
             const point&,
             const wallPointData<Type>& w2,
-            const scalar tol
+            const scalar tol,
+            TrackingData& td
         );
 
 
@@ -112,35 +114,41 @@ public:
 
             //- Influence of neighbouring face.
             //  Calls update(...) with cellCentre of cellI
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh& mesh,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const wallPointData<Type>& neighbourWallInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell.
             //  Calls update(...) with faceCentre of faceI
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh& mesh,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const wallPointData<Type>& neighbourWallInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face.
             //  Merge new and old info.
             //  Calls update(...) with faceCentre of faceI
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh& mesh,
                 const label thisFaceI,
                 const wallPointData<Type>& neighbourWallInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
     // Member Operators
diff --git a/src/meshTools/cellDist/wallPoint/wallPointDataI.H b/src/meshTools/cellDist/wallPoint/wallPointDataI.H
index f059f6b5e9..9bd47575dd 100644
--- a/src/meshTools/cellDist/wallPoint/wallPointDataI.H
+++ b/src/meshTools/cellDist/wallPoint/wallPointDataI.H
@@ -21,8 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Description
-
 \*---------------------------------------------------------------------------*/
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -32,19 +30,20 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-
 // Update this with w2 if w2 nearer to pt.
 template <class Type>
+template <class TrackingData>
 inline bool wallPointData<Type>::update
 (
     const point& pt,
     const wallPointData<Type>& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     scalar dist2 = magSqr(pt - w2.origin());
 
-    if (!valid())
+    if (!valid(td))
     {
         // current not yet set so use any value
         distSqr() = dist2;
@@ -122,13 +121,15 @@ inline Type& wallPointData<Type>::data()
 
 // Update this with w2 if w2 nearer to pt.
 template <class Type>
+template <class TrackingData>
 inline bool wallPointData<Type>::updateCell
 (
     const polyMesh& mesh,
     const label thisCellI,
     const label,
     const wallPointData<Type>& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const vectorField& cellCentres = mesh.primitiveMesh::cellCentres();
@@ -137,20 +138,23 @@ inline bool wallPointData<Type>::updateCell
     (
         cellCentres[thisCellI],
         neighbourWallInfo,
-        tol
+        tol,
+        td
     );
 }
 
 
 // Update this with w2 if w2 nearer to pt.
 template <class Type>
+template <class TrackingData>
 inline bool wallPointData<Type>::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const label,
     const wallPointData<Type>& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const vectorField& faceCentres = mesh.faceCentres();
@@ -159,19 +163,22 @@ inline bool wallPointData<Type>::updateFace
     (
         faceCentres[thisFaceI],
         neighbourWallInfo,
-        tol
+        tol,
+        td
     );
 }
 
 
 // Update this with w2 if w2 nearer to pt.
 template <class Type>
+template <class TrackingData>
 inline bool wallPointData<Type>::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const wallPointData<Type>& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     const vectorField& faceCentres = mesh.faceCentres();
@@ -180,7 +187,8 @@ inline bool wallPointData<Type>::updateFace
     (
         faceCentres[thisFaceI],
         neighbourWallInfo,
-        tol
+        tol,
+        td
     );
 }
 
diff --git a/src/meshTools/cellDist/wallPoint/wallPointI.H b/src/meshTools/cellDist/wallPoint/wallPointI.H
index b999790e76..e9cdc63d1e 100644
--- a/src/meshTools/cellDist/wallPoint/wallPointI.H
+++ b/src/meshTools/cellDist/wallPoint/wallPointI.H
@@ -29,11 +29,13 @@ License
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool Foam::wallPoint::update
 (
     const point& pt,
     const wallPoint& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     //Already done in calling algorithm
@@ -45,7 +47,7 @@ inline bool Foam::wallPoint::update
 
     scalar dist2 = magSqr(pt - w2.origin());
 
-    if (!valid())
+    if (!valid(td))
     {
         // current not yet set so use any value
         distSqr_ = dist2;
@@ -128,18 +130,21 @@ inline Foam::scalar& Foam::wallPoint::distSqr()
 }
 
 
-inline bool Foam::wallPoint::valid() const
+template<class TrackingData>
+inline bool Foam::wallPoint::valid(TrackingData& td) const
 {
     return origin_ != point::max;
 }
 
 
 // Checks for cyclic faces
+template<class TrackingData>
 inline bool Foam::wallPoint::sameGeometry
 (
     const polyMesh&,
     const wallPoint& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
  const
 {
@@ -163,22 +168,26 @@ inline bool Foam::wallPoint::sameGeometry
 }
 
 
+template<class TrackingData>
 inline void Foam::wallPoint::leaveDomain
 (
     const polyMesh&,
     const polyPatch&,
     const label,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {
     origin_ -= faceCentre;
 }
 
 
+template<class TrackingData>
 inline void Foam::wallPoint::transform
 (
     const polyMesh&,
-    const tensor& rotTensor
+    const tensor& rotTensor,
+    TrackingData& td
 )
 {
     origin_ = Foam::transform(rotTensor, origin_);
@@ -187,12 +196,14 @@ inline void Foam::wallPoint::transform
 
 // Update absolute geometric quantities. Note that distance (distSqr_)
 // is not affected by leaving/entering domain.
+template<class TrackingData>
 inline void Foam::wallPoint::enterDomain
 (
     const polyMesh&,
     const polyPatch&,
     const label,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData& td
 )
 {
     // back to absolute form
@@ -201,13 +212,15 @@ inline void Foam::wallPoint::enterDomain
 
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool Foam::wallPoint::updateCell
 (
     const polyMesh& mesh,
     const label thisCellI,
     const label neighbourFaceI,
     const wallPoint& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     return
@@ -215,19 +228,22 @@ inline bool Foam::wallPoint::updateCell
         (
             mesh.cellCentres()[thisCellI],
             neighbourWallInfo,
-            tol
+            tol,
+            td
         );
     }
 
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool Foam::wallPoint::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const label neighbourCellI,
     const wallPoint& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     return
@@ -235,17 +251,20 @@ inline bool Foam::wallPoint::updateFace
         (
             mesh.faceCentres()[thisFaceI],
             neighbourWallInfo,
-            tol
+            tol,
+            td
         );
 }
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool Foam::wallPoint::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const wallPoint& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     return
@@ -253,11 +272,23 @@ inline bool Foam::wallPoint::updateFace
         (
             mesh.faceCentres()[thisFaceI],
             neighbourWallInfo,
-            tol
+            tol,
+            td
         );
 }
 
 
+template <class TrackingData>
+inline bool Foam::wallPoint::equal
+(
+    const wallPoint& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 inline bool Foam::wallPoint::operator==(const Foam::wallPoint& rhs) const
diff --git a/src/parallel/decompose/decompositionMethods/structuredDecomp/structuredDecomp.C b/src/parallel/decompose/decompositionMethods/structuredDecomp/structuredDecomp.C
index 0e837a88f8..1c1e41ca74 100644
--- a/src/parallel/decompose/decompositionMethods/structuredDecomp/structuredDecomp.C
+++ b/src/parallel/decompose/decompositionMethods/structuredDecomp/structuredDecomp.C
@@ -155,7 +155,7 @@ Foam::labelList Foam::structuredDecomp::decompose
     bool haveWarned = false;
     forAll(finalDecomp, cellI)
     {
-        if (!cellData[cellI].valid())
+        if (!cellData[cellI].valid(deltaCalc.data()))
         {
             if (!haveWarned)
             {
diff --git a/src/parallel/decompose/decompositionMethods/structuredDecomp/topoDistanceData.H b/src/parallel/decompose/decompositionMethods/structuredDecomp/topoDistanceData.H
index d5c6e840d8..37d5f7f331 100644
--- a/src/parallel/decompose/decompositionMethods/structuredDecomp/topoDistanceData.H
+++ b/src/parallel/decompose/decompositionMethods/structuredDecomp/topoDistanceData.H
@@ -97,71 +97,90 @@ public:
 
             //- Check whether origin has been changed at all or
             //  still contains original (invalid) value.
-            inline bool valid() const;
+            template<class TrackingData>
+            inline bool valid(TrackingData& td) const;
 
             //- Check for identical geometrical data. Used for cyclics checking.
+            template<class TrackingData>
             inline bool sameGeometry
             (
                 const polyMesh&,
                 const topoDistanceData&,
-                const scalar
+                const scalar,
+                TrackingData& td
             ) const;
 
             //- Convert any absolute coordinates into relative to (patch)face
             //  centre
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Reverse of leaveDomain
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyMesh&,
                 const polyPatch&,
                 const label patchFaceI,
-                const point& faceCentre
+                const point& faceCentre,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to any coordinates
+            template<class TrackingData>
             inline void transform
             (
                 const polyMesh&,
-                const tensor&
+                const tensor&,
+                TrackingData& td
             );
 
             //- Influence of neighbouring face.
+            template<class TrackingData>
             inline bool updateCell
             (
                 const polyMesh&,
                 const label thisCellI,
                 const label neighbourFaceI,
                 const topoDistanceData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of neighbouring cell.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const label neighbourCellI,
                 const topoDistanceData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same face.
+            template<class TrackingData>
             inline bool updateFace
             (
                 const polyMesh&,
                 const label thisFaceI,
                 const topoDistanceData& neighbourInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const topoDistanceData&, TrackingData& td) const;
+
     // Member Operators
 
         // Needed for List IO
diff --git a/src/parallel/decompose/decompositionMethods/structuredDecomp/topoDistanceDataI.H b/src/parallel/decompose/decompositionMethods/structuredDecomp/topoDistanceDataI.H
index 5fcf99511c..7ae06e331b 100644
--- a/src/parallel/decompose/decompositionMethods/structuredDecomp/topoDistanceDataI.H
+++ b/src/parallel/decompose/decompositionMethods/structuredDecomp/topoDistanceDataI.H
@@ -53,18 +53,21 @@ inline Foam::topoDistanceData::topoDistanceData
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline bool Foam::topoDistanceData::valid() const
+template <class TrackingData>
+inline bool Foam::topoDistanceData::valid(TrackingData& td) const
 {
     return distance_ != -1;
 }
 
 
 // No geometric data so never any problem on cyclics
+template <class TrackingData>
 inline bool Foam::topoDistanceData::sameGeometry
 (
     const polyMesh&,
     const topoDistanceData&,
-    const scalar
+    const scalar,
+    TrackingData&
 ) const
 {
     return true;
@@ -72,44 +75,52 @@ inline bool Foam::topoDistanceData::sameGeometry
 
 
 // No geometric data.
+template <class TrackingData>
 inline void Foam::topoDistanceData::leaveDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData&
 )
 {}
 
 
 // No geometric data.
+template <class TrackingData>
 inline void Foam::topoDistanceData::transform
 (
     const polyMesh&,
-    const tensor& rotTensor
+    const tensor& rotTensor,
+    TrackingData&
 )
 {}
 
 
 // No geometric data.
+template <class TrackingData>
 inline void Foam::topoDistanceData::enterDomain
 (
     const polyMesh&,
     const polyPatch& patch,
     const label patchFaceI,
-    const point& faceCentre
+    const point& faceCentre,
+    TrackingData&
 )
 {}
 
 
 // Update cell with neighbouring face information
+template <class TrackingData>
 inline bool Foam::topoDistanceData::updateCell
 (
     const polyMesh&,
     const label thisCellI,
     const label neighbourFaceI,
     const topoDistanceData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData&
 )
 {
     if (distance_ == -1)
@@ -126,13 +137,15 @@ inline bool Foam::topoDistanceData::updateCell
 
 
 // Update face with neighbouring cell information
+template <class TrackingData>
 inline bool Foam::topoDistanceData::updateFace
 (
     const polyMesh& mesh,
     const label thisFaceI,
     const label neighbourCellI,
     const topoDistanceData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData&
 )
 {
     // From cell to its faces.
@@ -150,12 +163,14 @@ inline bool Foam::topoDistanceData::updateFace
 
 
 // Update face with coupled face information
+template <class TrackingData>
 inline bool Foam::topoDistanceData::updateFace
 (
     const polyMesh&,
     const label thisFaceI,
     const topoDistanceData& neighbourInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData&
 )
 {
     // From face to face (e.g. coupled faces)
@@ -171,6 +186,17 @@ inline bool Foam::topoDistanceData::updateFace
 }
 
 
+template <class TrackingData>
+inline bool Foam::topoDistanceData::equal
+(
+    const topoDistanceData& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 inline bool Foam::topoDistanceData::operator==
diff --git a/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDelta.C b/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDelta.C
index d97395d864..7d10bfc8d6 100644
--- a/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDelta.C
+++ b/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDelta.C
@@ -37,9 +37,6 @@ namespace Foam
 defineTypeNameAndDebug(smoothDelta, 0);
 addToRunTimeSelectionTable(LESdelta, smoothDelta, dictionary);
 
-scalar smoothDelta::deltaData::maxDeltaRatio = 1.2;
-
-
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Fill changedFaces (with face labels) and changedFacesInfo (with delta)
@@ -69,12 +66,12 @@ void smoothDelta::setChangedFaces
 
         // Check if owner delta much larger than neighbour delta or vice versa
 
-        if (ownDelta > deltaData::maxDeltaRatio * neiDelta)
+        if (ownDelta > maxDeltaRatio_ * neiDelta)
         {
             changedFaces.append(faceI);
             changedFacesInfo.append(deltaData(ownDelta));
         }
-        else if (neiDelta > deltaData::maxDeltaRatio * ownDelta)
+        else if (neiDelta > maxDeltaRatio_ * ownDelta)
         {
             changedFaces.append(faceI);
             changedFacesInfo.append(deltaData(neiDelta));
@@ -108,7 +105,6 @@ void smoothDelta::setChangedFaces
 
 void smoothDelta::calcDelta()
 {
-    deltaData::maxDeltaRatio = maxDeltaRatio_;
     const volScalarField& geometricDelta = geometricDelta_();
 
     // Fill changed faces with info
@@ -130,14 +126,15 @@ void smoothDelta::calcDelta()
 
 
     // Propagate information over whole domain.
-    FaceCellWave<deltaData> deltaCalc
+    FaceCellWave<deltaData, scalar> deltaCalc
     (
         mesh_,
         changedFaces,
         changedFacesInfo,
         faceDeltaData,
         cellDeltaData,
-        mesh_.globalData().nTotalCells()+1   // max iterations
+        mesh_.globalData().nTotalCells()+1,  // max iterations
+        maxDeltaRatio_
     );
 
     forAll(delta_, cellI)
diff --git a/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDelta.H b/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDelta.H
index 627c23801e..7a302b9379 100644
--- a/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDelta.H
+++ b/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDelta.H
@@ -76,23 +76,19 @@ class smoothDelta
 
                 //- Update. Gets information from neighbouring face/cell and
                 //  uses this to update itself (if nessecary) and return true.
+                template<class TrackingData>
                 inline bool update
                 (
                     const deltaData& w2,
                     const scalar scale,
-                    const scalar tol
+                    const scalar tol,
+                    TrackingData& td
                 );
 
 
         public:
 
-            // Static data members
-
-                //- delta fraction
-                static scalar maxDeltaRatio;
-
-
-                // Constructors
+            // Constructors
 
                 //- Construct null
                 inline deltaData();
@@ -115,72 +111,91 @@ class smoothDelta
 
                     //- Check whether origin has been changed at all or
                     //  still contains original (invalid) value.
-                    inline bool valid() const;
+                    template<class TrackingData>
+                    inline bool valid(TrackingData& td) const;
 
                     //- Check for identical geometrical data.
                     //  Used for cyclics checking.
+                    template<class TrackingData>
                     inline bool sameGeometry
                     (
                         const polyMesh&,
                         const deltaData&,
-                        const scalar
+                        const scalar,
+                        TrackingData& td
                     ) const;
 
                     //- Convert any absolute coordinates into relative to
                     //  (patch)face centre
+                    template<class TrackingData>
                     inline void leaveDomain
                     (
                         const polyMesh&,
                         const polyPatch&,
                         const label patchFaceI,
-                        const point& faceCentre
+                        const point& faceCentre,
+                        TrackingData& td
                     );
 
                     //- Reverse of leaveDomain
+                    template<class TrackingData>
                     inline void enterDomain
                     (
                         const polyMesh&,
                         const polyPatch&,
                         const label patchFaceI,
-                        const point& faceCentre
+                        const point& faceCentre,
+                        TrackingData& td
                     );
 
                     //- Apply rotation matrix to any coordinates
+                    template<class TrackingData>
                     inline void transform
                     (
                         const polyMesh&,
-                        const tensor&
+                        const tensor&,
+                        TrackingData& td
                     );
 
                     //- Influence of neighbouring face.
+                    template<class TrackingData>
                     inline bool updateCell
                     (
                         const polyMesh&,
                         const label thisCellI,
                         const label neighbourFaceI,
                         const deltaData& neighbourInfo,
-                        const scalar tol
+                        const scalar tol,
+                        TrackingData& td
                     );
 
                     //- Influence of neighbouring cell.
+                    template<class TrackingData>
                     inline bool updateFace
                     (
                         const polyMesh&,
                         const label thisFaceI,
                         const label neighbourCellI,
                         const deltaData& neighbourInfo,
-                        const scalar tol
+                        const scalar tol,
+                        TrackingData& td
                     );
 
                     //- Influence of different value on same face.
+                    template<class TrackingData>
                     inline bool updateFace
                     (
                         const polyMesh&,
                         const label thisFaceI,
                         const deltaData& neighbourInfo,
-                        const scalar tol
+                        const scalar tol,
+                        TrackingData& td
                     );
 
+                    //- Same (like operator==)
+                    template<class TrackingData>
+                    inline bool equal(const deltaData&, TrackingData& td) const;
+
                 // Member Operators
 
                     // Needed for List IO
diff --git a/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDeltaDeltaDataI.H b/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDeltaDeltaDataI.H
index 0508610d7a..cf7ebbc526 100644
--- a/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDeltaDeltaDataI.H
+++ b/src/turbulenceModels/LES/LESdeltas/smoothDelta/smoothDeltaDeltaDataI.H
@@ -30,16 +30,17 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-
 // Update this with w2 if applicable
+template<class TrackingData>
 inline bool smoothDelta::deltaData::update
 (
     const smoothDelta::deltaData& w2,
     const scalar scale,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    if (!valid() || (delta_ < VSMALL))
+    if (!valid(td) || (delta_ < VSMALL))
     {
         // My delta not set. Take over neighbour.
         delta_ = w2.delta()/scale;
@@ -82,94 +83,120 @@ inline smoothDelta::deltaData::deltaData(const scalar delta)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline bool smoothDelta::deltaData::valid() const
+template<class TrackingData>
+inline bool smoothDelta::deltaData::valid(TrackingData& td) const
 {
     return delta_ > -SMALL;
 }
 
 
 // Checks for cyclic faces
+template<class TrackingData>
 inline bool smoothDelta::deltaData::sameGeometry
 (
     const polyMesh&,
     const deltaData&,
-    const scalar
+    const scalar,
+    TrackingData& td
 ) const
 {
     return true;
 }
 
 
+template<class TrackingData>
 inline void smoothDelta::deltaData::leaveDomain
 (
     const polyMesh&,
     const polyPatch&,
     const label,
-    const point&
+    const point&,
+    TrackingData& td
 )
 {}
 
 
+template<class TrackingData>
 inline void smoothDelta::deltaData::transform
 (
     const polyMesh&,
-    const tensor&
+    const tensor&,
+    TrackingData& td
 )
 {}
 
 
 // Update absolute geometric quantities.
+template<class TrackingData>
 inline void smoothDelta::deltaData::enterDomain
 (
     const polyMesh&,
     const polyPatch&,
     const label,
-    const point&
+    const point&,
+    TrackingData& td
 )
 {}
 
 
 // Update this (cellI) with face information.
+template<class TrackingData>
 inline bool smoothDelta::deltaData::updateCell
 (
     const polyMesh&,
     const label,
     const label,
     const deltaData& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     // Take over info from face if more than deltaRatio larger.
-    return update(neighbourWallInfo, maxDeltaRatio, tol);
+    return update(neighbourWallInfo, td, tol, td);
 }
 
 
 // Update this (face) with cell information.
+template<class TrackingData>
 inline bool smoothDelta::deltaData::updateFace
 (
     const polyMesh&,
     const label,
     const label,
     const deltaData& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     // Take over information from cell without any scaling (scale = 1.0)
-    return update(neighbourWallInfo, 1.0, tol);
+    return update(neighbourWallInfo, 1.0, tol, td);
 }
 
 
 // Update this (face) with coupled face information.
+template<class TrackingData>
 inline bool smoothDelta::deltaData::updateFace
 (
     const polyMesh&,
     const label,
     const deltaData& neighbourWallInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     // Take over information from coupled face without any scaling (scale = 1.0)
-    return update(neighbourWallInfo, 1.0, tol);
+    return update(neighbourWallInfo, 1.0, tol, td);
+}
+
+
+template <class TrackingData>
+inline bool smoothDelta::deltaData::equal
+(
+    const deltaData& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
 }