diff --git a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/pointEdgeStructuredWalk.H b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/pointEdgeStructuredWalk.H
index 55d003a16a..e6852dffa2 100644
--- a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/pointEdgeStructuredWalk.H
+++ b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/pointEdgeStructuredWalk.H
@@ -37,11 +37,7 @@ SourceFiles
 #define pointEdgeStructuredWalk_H
 
 #include "point.H"
-#include "label.H"
-#include "scalar.H"
 #include "tensor.H"
-#include "pTraits.H"
-
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -75,10 +71,12 @@ class pointEdgeStructuredWalk
     // Private Member Functions
 
         //- Evaluate distance to point.
+        template<class TrackingData>
         inline bool update
         (
             const pointEdgeStructuredWalk& w2,
-            const scalar tol
+            const scalar tol,
+            TrackingData& td
         );
 
 public:
@@ -111,80 +109,105 @@ public:
 
         // Needed by meshWave
 
-            //- Check whether still contains original (invalid) value.
-            inline bool valid() const;
+            //- Check whether origin has been changed at all or
+            //  still contains original (invalid) value.
+            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 pointEdgeStructuredWalk&,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             ) const;
 
             //- Convert origin to relative vector to leaving point
             //  (= point coordinate)
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyPatch& patch,
                 const label patchPointI,
-                const point& pos
+                const point& pos,
+                TrackingData& td
             );
 
             //- Convert relative origin to absolute by adding entering point
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyPatch& patch,
                 const label patchPointI,
-                const point& pos
+                const point& pos,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to origin
-            inline void transform(const tensor& rotTensor);
+            template<class TrackingData>
+            inline void transform
+            (
+                const tensor& rotTensor,
+                TrackingData& td
+            );
 
             //- Influence of edge on point
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const polyMesh& mesh,
                 const label pointI,
                 const label edgeI,
                 const pointEdgeStructuredWalk& edgeInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same point.
             //  Merge new and old info.
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const polyMesh& mesh,
                 const label pointI,
                 const pointEdgeStructuredWalk& newPointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same point.
             //  No information about current position whatsoever.
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const pointEdgeStructuredWalk& newPointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of point on edge.
+            template<class TrackingData>
             inline bool updateEdge
             (
                 const polyMesh& mesh,
                 const label edgeI,
                 const label pointI,
                 const pointEdgeStructuredWalk& pointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const pointEdgeStructuredWalk&, TrackingData&)
+            const;
+
 
     // Member Operators
 
         //Note: Used to determine whether to call update.
         inline bool operator==(const pointEdgeStructuredWalk&) const;
-
         inline bool operator!=(const pointEdgeStructuredWalk&) const;
 
 
diff --git a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/pointEdgeStructuredWalkI.H b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/pointEdgeStructuredWalkI.H
index d60e9b52f7..47adcdb775 100644
--- a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/pointEdgeStructuredWalkI.H
+++ b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/pointEdgeStructuredWalkI.H
@@ -29,13 +29,15 @@ License
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Update this with w2.
+template<class TrackingData>
 inline bool Foam::pointEdgeStructuredWalk::update
 (
     const pointEdgeStructuredWalk& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    if (!valid())
+    if (!valid(td))
     {
         // current not yet set. Walked from w2 to here (=point0)
         dist_ = w2.dist_ + mag(point0_-w2.previousPoint_);
@@ -105,17 +107,20 @@ inline const Foam::vector& Foam::pointEdgeStructuredWalk::data() const
 }
 
 
-inline bool Foam::pointEdgeStructuredWalk::valid() const
+template<class TrackingData>
+inline bool Foam::pointEdgeStructuredWalk::valid(TrackingData& td) const
 {
     return previousPoint_ != vector::max;
 }
 
 
 // Checks for cyclic points
+template<class TrackingData>
 inline bool Foam::pointEdgeStructuredWalk::sameGeometry
 (
     const pointEdgeStructuredWalk& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 ) const
 {
     scalar diff = Foam::mag(dist() - w2.dist());
@@ -138,18 +143,25 @@ inline bool Foam::pointEdgeStructuredWalk::sameGeometry
 }
 
 
+template<class TrackingData>
 inline void Foam::pointEdgeStructuredWalk::leaveDomain
 (
     const polyPatch& patch,
     const label patchPointI,
-    const point& coord
+    const point& coord,
+    TrackingData& td
 )
 {
     previousPoint_ -= coord;
 }
 
 
-inline void Foam::pointEdgeStructuredWalk::transform(const tensor& rotTensor)
+template<class TrackingData>
+inline void Foam::pointEdgeStructuredWalk::transform
+(
+    const tensor& rotTensor,
+    TrackingData& td
+)
 {
     previousPoint_ = Foam::transform(rotTensor, previousPoint_);
 }
@@ -157,11 +169,13 @@ inline void Foam::pointEdgeStructuredWalk::transform(const tensor& rotTensor)
 
 // Update absolute geometric quantities. Note that distance (dist_)
 // is not affected by leaving/entering domain.
+template<class TrackingData>
 inline void Foam::pointEdgeStructuredWalk::enterDomain
 (
     const polyPatch& patch,
     const label patchPointI,
-    const point& coord
+    const point& coord,
+    TrackingData& td
 )
 {
     // back to absolute form
@@ -170,18 +184,20 @@ inline void Foam::pointEdgeStructuredWalk::enterDomain
 
 
 // Update this with information from connected edge
+template<class TrackingData>
 inline bool Foam::pointEdgeStructuredWalk::updatePoint
 (
     const polyMesh& mesh,
     const label pointI,
     const label edgeI,
     const pointEdgeStructuredWalk& edgeInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     if (inZone())
     {
-        return update(edgeInfo, tol);
+        return update(edgeInfo, tol, td);
     }
     else
     {
@@ -191,17 +207,19 @@ inline bool Foam::pointEdgeStructuredWalk::updatePoint
 
 
 // Update this with new information on same point
+template<class TrackingData>
 inline bool Foam::pointEdgeStructuredWalk::updatePoint
 (
     const polyMesh& mesh,
     const label pointI,
     const pointEdgeStructuredWalk& newPointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     if (inZone())
     {
-        return update(newPointInfo, tol);
+        return update(newPointInfo, tol, td);
     }
     else
     {
@@ -211,29 +229,33 @@ inline bool Foam::pointEdgeStructuredWalk::updatePoint
 
 
 // Update this with new information on same point. No extra information.
+template<class TrackingData>
 inline bool Foam::pointEdgeStructuredWalk::updatePoint
 (
     const pointEdgeStructuredWalk& newPointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(newPointInfo, tol);
+    return update(newPointInfo, tol, td);
 }
 
 
 // Update this with information from connected point
+template<class TrackingData>
 inline bool Foam::pointEdgeStructuredWalk::updateEdge
 (
     const polyMesh& mesh,
     const label edgeI,
     const label pointI,
     const pointEdgeStructuredWalk& pointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
     if (inZone())
     {
-        return update(pointInfo, tol);
+        return update(pointInfo, tol, td);
     }
     else
     {
@@ -242,6 +264,17 @@ inline bool Foam::pointEdgeStructuredWalk::updateEdge
 }
 
 
+template <class TrackingData>
+inline bool Foam::pointEdgeStructuredWalk::equal
+(
+    const pointEdgeStructuredWalk& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 inline bool Foam::pointEdgeStructuredWalk::operator==
diff --git a/src/fvMotionSolver/motionDiffusivity/inversePointDistance/inversePointDistanceDiffusivity.C b/src/fvMotionSolver/motionDiffusivity/inversePointDistance/inversePointDistanceDiffusivity.C
index 7ed2a34bba..486052b9db 100644
--- a/src/fvMotionSolver/motionDiffusivity/inversePointDistance/inversePointDistanceDiffusivity.C
+++ b/src/fvMotionSolver/motionDiffusivity/inversePointDistance/inversePointDistanceDiffusivity.C
@@ -91,6 +91,9 @@ void Foam::inversePointDistanceDiffusivity::correct()
     List<pointEdgePoint> pointWallDist(mesh.nPoints());
     List<pointEdgePoint> edgeWallDist(mesh.nEdges());
 
+    int dummyTrackData = 0;
+
+
     {
         // Seeds
         List<pointEdgePoint> seedInfo(nPatchEdges);
@@ -108,7 +111,7 @@ void Foam::inversePointDistanceDiffusivity::correct()
             {
                 label pointI = meshPoints[i];
 
-                if (!pointWallDist[pointI].valid())
+                if (!pointWallDist[pointI].valid(dummyTrackData))
                 {
                     // Not yet seeded
                     seedInfo[nPatchEdges] = pointEdgePoint
@@ -135,7 +138,8 @@ void Foam::inversePointDistanceDiffusivity::correct()
 
             pointWallDist,
             edgeWallDist,
-            mesh.globalData().nTotalPoints() // max iterations
+            mesh.globalData().nTotalPoints(),// max iterations
+            dummyTrackData
         );
     }
 
diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C
index 292ec5d142..d7abe804f7 100644
--- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C
+++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C
@@ -756,6 +756,9 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo
     // Distance to wall
     List<pointData> pointWallDist(mesh.nPoints());
 
+    // Dummy additional info for PointEdgeWave
+    int dummyTrackData = 0;
+
 
     // 1. Calculate distance to points where displacement is specified.
     {
@@ -783,7 +786,8 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo
             wallInfo,
             pointWallDist,
             edgeWallDist,
-            mesh.globalData().nTotalPoints()    // max iterations
+            mesh.globalData().nTotalPoints(),   // max iterations
+            dummyTrackData
         );
     }
 
@@ -813,7 +817,7 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo
                 {
                     label pointI = e[ep];
 
-                    if (!pointMedialDist[pointI].valid())
+                    if (!pointMedialDist[pointI].valid(dummyTrackData))
                     {
                         maxPoints.append(pointI);
                         maxInfo.append
@@ -857,7 +861,7 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo
                 {
                     label pointI = meshPoints[i];
 
-                    if (!pointMedialDist[pointI].valid())
+                    if (!pointMedialDist[pointI].valid(dummyTrackData))
                     {
                         maxPoints.append(pointI);
                         maxInfo.append
@@ -888,7 +892,8 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo
 
             pointMedialDist,
             edgeMedialDist,
-            mesh.globalData().nTotalPoints()    // max iterations
+            mesh.globalData().nTotalPoints(),   // max iterations
+            dummyTrackData
         );
 
         // Extract medial axis distance as pointScalarField
diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointData.C b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointData.C
index 6d97d1cd00..1d7f5090ae 100644
--- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointData.C
+++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointData.C
@@ -32,19 +32,20 @@ Foam::Ostream& Foam::operator<<(Ostream& os, const pointData& wDist)
     if (os.format() == IOstream::ASCII)
     {
         return os
-            << wDist.origin() << token::SPACE << wDist.distSqr()
+            << static_cast<const pointEdgePoint&>(wDist)
             << token::SPACE << wDist.s() << token::SPACE << wDist.v();
     }
     else
     {
         return os
-            << wDist.origin() << wDist.distSqr() << wDist.s() << wDist.v();
+            << static_cast<const pointEdgePoint&>(wDist)
+            << wDist.s() << wDist.v();
     }
 }
 
 Foam::Istream& Foam::operator>>(Istream& is, pointData& wDist)
 {
-    return is >> wDist.origin_ >> wDist.distSqr_ >> wDist.s_ >> wDist.v_;
+    return is >> static_cast<pointEdgePoint&>(wDist) >> wDist.s_ >> wDist.v_;
 }
 
 // ************************************************************************* //
diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointData.H b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointData.H
index 050b9d92c7..7719079821 100644
--- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointData.H
+++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointData.H
@@ -25,10 +25,10 @@ Class
     Foam::pointData
 
 Description
-    Holds information regarding nearest wall point. Used in pointEdgeWave.
-    (so not standard meshWave)
-
-    To be used in wall distance calculation.
+    Variant of pointEdgePoint with some transported additional data. 
+    WIP - should be templated on data like wallDistData.
+    Passive vector v_ is not a coordinate (so no enterDomain/leaveDomain
+    transformation needed)
 
 SourceFiles
     pointDataI.H
@@ -39,9 +39,10 @@ SourceFiles
 #ifndef pointData_H
 #define pointData_H
 
-#include "point.H"
-#include "label.H"
-#include "tensor.H"
+#include "pointEdgePoint.H"
+//#include "point.H"
+//#include "label.H"
+//#include "tensor.H"
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -49,51 +50,22 @@ SourceFiles
 namespace Foam
 {
 
-// Class forward declarations
-class polyPatch;
-class polyMesh;
-
 /*---------------------------------------------------------------------------*\
                            Class pointData Declaration
 \*---------------------------------------------------------------------------*/
 
 class pointData
+:
+    public pointEdgePoint
 {
     // Private data
 
-        //- position of nearest wall center
-        point origin_;
-
-        //- normal distance (squared) from point to origin
-        scalar distSqr_;
-
         //- additional information.
         scalar s_;
 
         //- additional information.
         vector v_;
 
-
-    // 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.
-        inline bool update
-        (
-            const point&,
-            const pointData& w2,
-            const scalar tol
-        );
-
-        //- Combine current with w2. Update distSqr, origin if w2 has smaller
-        //  quantities and returns true.
-        inline bool update
-        (
-            const pointData& w2,
-            const scalar tol
-        );
-
 public:
 
     // Constructors
@@ -118,10 +90,6 @@ public:
 
         // Access
 
-            inline const point& origin() const;
-
-            inline scalar distSqr() const;
-
             inline scalar s() const;
 
             inline const vector& v() const;
@@ -129,81 +97,60 @@ public:
 
         // Needed by meshWave
 
-            //- Check whether origin has been changed at all or
-            //  still contains original (invalid) value.
-            inline bool valid() const;
-
-            //- Check for identical geometrical data. Used for cyclics checking.
-            inline bool sameGeometry(const pointData&, const scalar tol)
-             const;
-
-            //- Convert origin to relative vector to leaving point
-            //  (= point coordinate)
-            inline void leaveDomain
-            (
-                const polyPatch& patch,
-                const label patchPointI,
-                const point& pos
-            );
-
-            //- Convert relative origin to absolute by adding entering point
-            inline void enterDomain
-            (
-                const polyPatch& patch,
-                const label patchPointI,
-                const point& pos
-            );
-
             //- Apply rotation matrix to origin
-            inline void transform(const tensor& rotTensor);
+            template<class TrackingData>
+            inline void transform
+            (
+                const tensor& rotTensor,
+                TrackingData& td
+            );
 
             //- Influence of edge on point
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const polyMesh& mesh,
                 const label pointI,
                 const label edgeI,
                 const pointData& edgeInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same point.
             //  Merge new and old info.
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const polyMesh& mesh,
                 const label pointI,
                 const pointData& newPointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same point.
             //  No information about current position whatsoever.
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const pointData& newPointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of point on edge.
+            template<class TrackingData>
             inline bool updateEdge
             (
                 const polyMesh& mesh,
                 const label edgeI,
                 const label pointI,
                 const pointData& pointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
-
-    // Member Operators
-
-        //Note: Used to determine whether to call update.
-        inline bool operator==(const pointData&) const;
-
-        inline bool operator!=(const pointData&) const;
-
-
     // IOstream Operators
 
         friend Ostream& operator<<(Ostream&, const pointData&);
diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointDataI.H b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointDataI.H
index a04038c247..c2228eca82 100644
--- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointDataI.H
+++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/pointData/pointDataI.H
@@ -26,118 +26,12 @@ License
 #include "polyMesh.H"
 #include "transform.H"
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-// Update this with w2 if w2 nearer to pt.
-inline bool Foam::pointData::update
-(
-    const point& pt,
-    const pointData& w2,
-    const scalar tol
-)
-{
-    scalar dist2 = magSqr(pt - w2.origin());
-
-    if (!valid())
-    {
-        distSqr_ = dist2;
-        origin_ = w2.origin();
-        s_ = w2.s();
-        v_ = w2.v();
-
-        return true;
-    }
-
-
-//    if (v_ != w2.v())
-//    {
-//        return false;
-//    }
-
-
-    scalar diff = distSqr_ - dist2;
-
-    if (diff < 0)
-    {
-        // already nearer to pt
-        return false;
-    }
-
-    if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
-    {
-        // don't propagate small changes
-        return false;
-    }
-    else
-    {
-        // update with new values
-        distSqr_ = dist2;
-        origin_ = w2.origin();
-        s_ = w2.s();
-        v_ = w2.v();
-
-        return true;
-    }
-}
-
-
-// Update this with w2 (information on same point)
-inline bool Foam::pointData::update
-(
-    const pointData& w2,
-    const scalar tol
-)
-{
-    if (!valid())
-    {
-        // current not yet set so use any value
-        distSqr_ = w2.distSqr();
-        origin_ = w2.origin();
-        s_ = w2.s();
-        v_ = w2.v();
-        return true;
-    }
-
-
-//    if (v_ != w2.v())
-//    {
-//        return false;
-//    }
-
-
-    scalar diff = distSqr_ - w2.distSqr();
-
-    if (diff < 0)
-    {
-        // already nearer to pt
-        return false;
-    }
-
-    if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
-    {
-        // don't propagate small changes
-        return false;
-    }
-    else
-    {
-        // update with new values
-        distSqr_ =  w2.distSqr();
-        origin_ = w2.origin();
-        s_ = w2.s();
-        v_ = w2.v();
-
-        return true;
-    }
-}
-
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Null constructor
 inline Foam::pointData::pointData()
 :
-    origin_(point::max),
-    distSqr_(GREAT),
+    pointEdgePoint(),
     s_(GREAT),
     v_(point::max)
 {}
@@ -152,8 +46,7 @@ inline Foam::pointData::pointData
     const vector& v
 )
 :
-    origin_(origin),
-    distSqr_(distSqr),
+    pointEdgePoint(origin, distSqr),
     s_(s),
     v_(v)
 {}
@@ -162,8 +55,7 @@ inline Foam::pointData::pointData
 // Construct as copy
 inline Foam::pointData::pointData(const pointData& wpt)
 :
-    origin_(wpt.origin()),
-    distSqr_(wpt.distSqr()),
+    pointEdgePoint(wpt),
     s_(wpt.s()),
     v_(wpt.v())
 {}
@@ -171,18 +63,6 @@ inline Foam::pointData::pointData(const pointData& wpt)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline const Foam::point& Foam::pointData::origin() const
-{
-    return origin_;
-}
-
-
-inline Foam::scalar Foam::pointData::distSqr() const
-{
-    return distSqr_;
-}
-
-
 inline Foam::scalar Foam::pointData::s() const
 {
     return s_;
@@ -195,157 +75,143 @@ inline const Foam::vector& Foam::pointData::v() const
 }
 
 
-inline bool Foam::pointData::valid() const
-{
-    return origin_ != point::max;
-}
-
-
-// Checks for cyclic points
-inline bool Foam::pointData::sameGeometry
+template <class TrackingData>
+inline void Foam::pointData::transform
 (
-    const pointData& w2,
-    const scalar tol
-) const
-{
-    scalar diff = Foam::mag(distSqr() - w2.distSqr());
-
-    if (diff < SMALL)
-    {
-        return true;
-    }
-    else
-    {
-        if ((distSqr() > SMALL) && ((diff/distSqr()) < tol))
-        {
-            return true;
-        }
-        else
-        {
-            return false;
-        }
-    }
-}
-
-
-inline void Foam::pointData::leaveDomain
-(
-    const polyPatch& patch,
-    const label patchPointI,
-    const point& coord
+    const tensor& rotTensor,
+    TrackingData& td
 )
 {
-    origin_ -= coord;
-}
-
-
-inline void Foam::pointData::transform(const tensor& rotTensor)
-{
-    origin_ = Foam::transform(rotTensor, origin_);
-}
-
-
-// Update absolute geometric quantities. Note that distance (distSqr_)
-// is not affected by leaving/entering domain.
-inline void Foam::pointData::enterDomain
-(
-    const polyPatch& patch,
-    const label patchPointI,
-    const point& coord
-)
-{
-    // back to absolute form
-    origin_ += coord;
+    pointEdgePoint::transform(rotTensor, td);
+    v_ = Foam::transform(rotTensor, v_);
 }
 
 
 // Update this with information from connected edge
+template <class TrackingData>
 inline bool Foam::pointData::updatePoint
 (
     const polyMesh& mesh,
     const label pointI,
     const label edgeI,
     const pointData& edgeInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return
-        update
+    if
+    (
+        pointEdgePoint::updatePoint
         (
-            mesh.points()[pointI],
+            mesh,
+            pointI,
+            edgeI,
             edgeInfo,
-            tol
-        );
+            tol,
+            td
+        )
+    )
+    {
+        s_ = edgeInfo.s_;
+        v_ = edgeInfo.v_;
+        return true;
     }
-
+    else
+    {
+        return false;
+    }
+}
 
 // Update this with new information on same point
+template <class TrackingData>
 inline bool Foam::pointData::updatePoint
 (
     const polyMesh& mesh,
     const label pointI,
     const pointData& newPointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return
-        update
+    if
+    (
+        pointEdgePoint::updatePoint
         (
-            mesh.points()[pointI],
+            mesh,
+            pointI,
             newPointInfo,
-            tol
-        );
+            tol,
+            td
+        )
+    )
+    {
+        s_ = newPointInfo.s_;
+        v_ = newPointInfo.v_;
+        return true;
+    }
+    else
+    {
+        return false;
+    }
 }
 
 
 // Update this with new information on same point. No extra information.
+template <class TrackingData>
 inline bool Foam::pointData::updatePoint
 (
     const pointData& newPointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(newPointInfo, tol);
+    if (pointEdgePoint::updatePoint(newPointInfo, tol, td))
+    {
+        s_ = newPointInfo.s_;
+        v_ = newPointInfo.v_;
+        return true;
+    }
+    else
+    {
+        return false;
+    }
 }
 
 
 // Update this with information from connected point
+template <class TrackingData>
 inline bool Foam::pointData::updateEdge
 (
     const polyMesh& mesh,
     const label edgeI,
     const label pointI,
     const pointData& pointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
+
 )
 {
-    const pointField& points = mesh.points();
-
-    const edge& e = mesh.edges()[edgeI];
-
-    const point edgeMid(0.5*(points[e[0]] + points[e[1]]));
-
-    return
-        update
+    if
+    (
+        pointEdgePoint::updateEdge
         (
-            edgeMid,
+            mesh,
+            edgeI,
+            pointI,
             pointInfo,
-            tol
-        );
-}
-
-
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
-
-inline bool Foam::pointData::operator==(const pointData& rhs) const
-{
-    return origin() == rhs.origin();
-}
-
-
-inline bool Foam::pointData::operator!=(const pointData& rhs) const
-{
-    return !(*this == rhs);
+            tol,
+            td
+        )
+    )
+    {
+        s_ = pointInfo.s_;
+        v_ = pointInfo.v_;
+        return true;
+    }
+    else
+    {
+        return false;
+    }
 }
 
 
diff --git a/src/meshTools/PointEdgeWave/PointEdgeWave.C b/src/meshTools/PointEdgeWave/PointEdgeWave.C
index 304675d71f..a11d6601be 100644
--- a/src/meshTools/PointEdgeWave/PointEdgeWave.C
+++ b/src/meshTools/PointEdgeWave/PointEdgeWave.C
@@ -32,36 +32,25 @@ License
 #include "PstreamCombineReduceOps.H"
 #include "debug.H"
 #include "typeInfo.H"
+#include "globalMeshData.H"
+#include "pointFields.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
-template <class Type>
-Foam::scalar Foam::PointEdgeWave<Type>::propagationTol_ = 0.01;
+template <class Type, class TrackingData>
+Foam::scalar Foam::PointEdgeWave<Type, TrackingData>::propagationTol_ = 0.01;
 
-
-// Offset labelList. Used for transferring from one cyclic half to the other.
-template <class Type>
-void Foam::PointEdgeWave<Type>::offset(const label val, labelList& elems)
-{
-    forAll(elems, i)
-    {
-        elems[i] += val;
-    }
-}
+template <class Type, class TrackingData>
+Foam::label Foam::PointEdgeWave<Type, TrackingData>::dummyTrackData_ = 12345;
 
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-// Gets point-point correspondence. Is
-// - list of halfA points (in cyclic patch points)
-// - list of halfB points (can overlap with A!)
-// - for every patchPoint its corresponding point
 // Handle leaving domain. Implementation referred to Type
-template <class Type>
-void Foam::PointEdgeWave<Type>::leaveDomain
+template <class Type, class TrackingData>
+void Foam::PointEdgeWave<Type, TrackingData>::leaveDomain
 (
-    const polyPatch& meshPatch,
-    const primitivePatch& patch,
+    const polyPatch& patch,
     const labelList& patchPointLabels,
     List<Type>& pointInfo
 ) const
@@ -74,17 +63,16 @@ void Foam::PointEdgeWave<Type>::leaveDomain
 
         const point& pt = patch.points()[meshPoints[patchPointI]];
 
-        pointInfo[i].leaveDomain(meshPatch, patchPointI, pt);
+        pointInfo[i].leaveDomain(patch, patchPointI, pt, td_);
     }
 }
 
 
 // Handle entering domain. Implementation referred to Type
-template <class Type>
-void Foam::PointEdgeWave<Type>::enterDomain
+template <class Type, class TrackingData>
+void Foam::PointEdgeWave<Type, TrackingData>::enterDomain
 (
-    const polyPatch& meshPatch,
-    const primitivePatch& patch,
+    const polyPatch& patch,
     const labelList& patchPointLabels,
     List<Type>& pointInfo
 ) const
@@ -97,15 +85,16 @@ void Foam::PointEdgeWave<Type>::enterDomain
 
         const point& pt = patch.points()[meshPoints[patchPointI]];
 
-        pointInfo[i].enterDomain(meshPatch, patchPointI, pt);
+        pointInfo[i].enterDomain(patch, patchPointI, pt, td_);
     }
 }
 
 
 // Transform. Implementation referred to Type
-template <class Type>
-void Foam::PointEdgeWave<Type>::transform
+template <class Type, class TrackingData>
+void Foam::PointEdgeWave<Type, TrackingData>::transform
 (
+    const polyPatch& patch,
     const tensorField& rotTensor,
     List<Type>& pointInfo
 ) const
@@ -116,19 +105,22 @@ void Foam::PointEdgeWave<Type>::transform
 
         forAll(pointInfo, i)
         {
-            pointInfo[i].transform(T);
+            pointInfo[i].transform(T, td_);
         }
     }
     else
     {
         FatalErrorIn
         (
-            "PointEdgeWave<Type>::transform(const tensorField&, List<Type>&)"
-        )   << "Parallel cyclics not supported" << abort(FatalError);
+            "PointEdgeWave<Type, TrackingData>::transform"
+            "(const tensorField&, List<Type>&)"
+        )   << "Non-uniform transformation on patch " << patch.name()
+            << " not supported for point fields"
+            << abort(FatalError);
 
         forAll(pointInfo, i)
         {
-            pointInfo[i].transform(rotTensor[i]);
+            pointInfo[i].transform(rotTensor[i], td_);
         }
     }
 }
@@ -139,19 +131,18 @@ void Foam::PointEdgeWave<Type>::transform
 // Updates:
 //      - changedPoint_, changedPoints_, nChangedPoints_,
 //      - statistics: nEvals_, nUnvisitedPoints_
-template <class Type>
-bool Foam::PointEdgeWave<Type>::updatePoint
+template <class Type, class TrackingData>
+bool Foam::PointEdgeWave<Type, TrackingData>::updatePoint
 (
     const label pointI,
     const label neighbourEdgeI,
     const Type& neighbourInfo,
-    const scalar tol,
     Type& pointInfo
 )
 {
     nEvals_++;
 
-    bool wasValid = pointInfo.valid();
+    bool wasValid = pointInfo.valid(td_);
 
     bool propagate =
         pointInfo.updatePoint
@@ -160,7 +151,8 @@ bool Foam::PointEdgeWave<Type>::updatePoint
             pointI,
             neighbourEdgeI,
             neighbourInfo,
-            tol
+            propagationTol_,
+            td_
         );
 
     if (propagate)
@@ -172,7 +164,7 @@ bool Foam::PointEdgeWave<Type>::updatePoint
         }
     }
 
-    if (!wasValid && pointInfo.valid())
+    if (!wasValid && pointInfo.valid(td_))
     {
         --nUnvisitedPoints_;
     }
@@ -186,18 +178,17 @@ bool Foam::PointEdgeWave<Type>::updatePoint
 // Updates:
 //      - changedPoint_, changedPoints_, nChangedPoints_,
 //      - statistics: nEvals_, nUnvisitedPoints_
-template <class Type>
-bool Foam::PointEdgeWave<Type>::updatePoint
+template <class Type, class TrackingData>
+bool Foam::PointEdgeWave<Type, TrackingData>::updatePoint
 (
     const label pointI,
     const Type& neighbourInfo,
-    const scalar tol,
     Type& pointInfo
 )
 {
     nEvals_++;
 
-    bool wasValid = pointInfo.valid();
+    bool wasValid = pointInfo.valid(td_);
 
     bool propagate =
         pointInfo.updatePoint
@@ -205,7 +196,8 @@ bool Foam::PointEdgeWave<Type>::updatePoint
             mesh_,
             pointI,
             neighbourInfo,
-            tol
+            propagationTol_,
+            td_
         );
 
     if (propagate)
@@ -217,7 +209,7 @@ bool Foam::PointEdgeWave<Type>::updatePoint
         }
     }
 
-    if (!wasValid && pointInfo.valid())
+    if (!wasValid && pointInfo.valid(td_))
     {
         --nUnvisitedPoints_;
     }
@@ -231,19 +223,18 @@ bool Foam::PointEdgeWave<Type>::updatePoint
 // Updates:
 //      - changedEdge_, changedEdges_, nChangedEdges_,
 //      - statistics: nEvals_, nUnvisitedEdge_
-template <class Type>
-bool Foam::PointEdgeWave<Type>::updateEdge
+template <class Type, class TrackingData>
+bool Foam::PointEdgeWave<Type, TrackingData>::updateEdge
 (
     const label edgeI,
     const label neighbourPointI,
     const Type& neighbourInfo,
-    const scalar tol,
     Type& edgeInfo
 )
 {
     nEvals_++;
 
-    bool wasValid = edgeInfo.valid();
+    bool wasValid = edgeInfo.valid(td_);
 
     bool propagate =
         edgeInfo.updateEdge
@@ -252,7 +243,8 @@ bool Foam::PointEdgeWave<Type>::updateEdge
             edgeI,
             neighbourPointI,
             neighbourInfo,
-            tol
+            propagationTol_,
+            td_
         );
 
     if (propagate)
@@ -264,7 +256,7 @@ bool Foam::PointEdgeWave<Type>::updateEdge
         }
     }
 
-    if (!wasValid && edgeInfo.valid())
+    if (!wasValid && edgeInfo.valid(td_))
     {
         --nUnvisitedEdges_;
     }
@@ -274,9 +266,9 @@ bool Foam::PointEdgeWave<Type>::updateEdge
 
 
 // Check if patches of given type name are present
-template <class Type>
+template <class Type, class TrackingData>
 template <class PatchType>
-Foam::label Foam::PointEdgeWave<Type>::countPatchType() const
+Foam::label Foam::PointEdgeWave<Type, TrackingData>::countPatchType() const
 {
     label nPatches = 0;
 
@@ -291,199 +283,110 @@ Foam::label Foam::PointEdgeWave<Type>::countPatchType() const
 }
 
 
-// Collect changed patch points
-template <class Type>
-void Foam::PointEdgeWave<Type>::getChangedPatchPoints
-(
-    const primitivePatch& patch,
-
-    DynamicList<Type>& patchInfo,
-    DynamicList<label>& patchPoints,
-    DynamicList<label>& owner,
-    DynamicList<label>& ownerIndex
-) const
-{
-    const labelList& meshPoints = patch.meshPoints();
-    const faceList& localFaces = patch.localFaces();
-    const labelListList& pointFaces = patch.pointFaces();
-
-    forAll(meshPoints, patchPointI)
-    {
-        label meshPointI = meshPoints[patchPointI];
-
-        if (changedPoint_[meshPointI])
-        {
-            patchInfo.append(allPointInfo_[meshPointI]);
-            patchPoints.append(patchPointI);
-
-            label patchFaceI = pointFaces[patchPointI][0];
-
-            const face& f = localFaces[patchFaceI];
-
-            label index = findIndex(f, patchPointI);
-
-            owner.append(patchFaceI);
-            ownerIndex.append(index);
-        }
-    }
-
-    patchInfo.shrink();
-    patchPoints.shrink();
-    owner.shrink();
-    ownerIndex.shrink();
-}
-
-
-// Update overall for changed patch points
-template <class Type>
-void Foam::PointEdgeWave<Type>::updateFromPatchInfo
-(
-    const polyPatch& meshPatch,
-    const primitivePatch& patch,
-    const labelList& owner,
-    const labelList& ownerIndex,
-    List<Type>& patchInfo
-)
-{
-    const faceList& localFaces = patch.localFaces();
-    const labelList& meshPoints = patch.meshPoints();
-
-    // Get patch and mesh points.
-    labelList changedPatchPoints(patchInfo.size());
-    labelList changedMeshPoints(patchInfo.size());
-
-    forAll(owner, i)
-    {
-        label faceI = owner[i];
-
-        const face& f = localFaces[faceI];
-
-        label index = (f.size() - ownerIndex[i]) % f.size();
-
-        changedPatchPoints[i] = f[index];
-        changedMeshPoints[i] = meshPoints[f[index]];
-    }
-
-    // Adapt for entering domain
-    enterDomain(meshPatch, patch, changedPatchPoints, patchInfo);
-
-    // Merge received info
-    forAll(patchInfo, i)
-    {
-        updatePoint
-        (
-            changedMeshPoints[i],
-            patchInfo[i],
-            propagationTol_,
-            allPointInfo_[changedMeshPoints[i]]
-        );
-    }
-}
-
-
-
 // Transfer all the information to/from neighbouring processors
-template <class Type>
-void Foam::PointEdgeWave<Type>::handleProcPatches()
+template <class Type, class TrackingData>
+void Foam::PointEdgeWave<Type, TrackingData>::handleProcPatches()
 {
-    // 1. Send all point info on processor patches. Send as
-    // face label + offset in face.
+    // 1. Send all point info on processor patches.
 
-    forAll(mesh_.boundaryMesh(), patchI)
+    PstreamBuffers pBufs(Pstream::nonBlocking);
+
+    DynamicList<Type> patchInfo;
+    DynamicList<label> thisPoints;
+    DynamicList<label> nbrPoints;
+
+    forAll(mesh_.globalData().processorPatches(), i)
     {
-        const polyPatch& patch = mesh_.boundaryMesh()[patchI];
+        label patchI = mesh_.globalData().processorPatches()[i];
+        const processorPolyPatch& procPatch =
+            refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchI]);
 
-        if (isA<processorPolyPatch>(patch))
+        patchInfo.clear();
+        patchInfo.reserve(procPatch.nPoints());
+        thisPoints.clear();
+        thisPoints.reserve(procPatch.nPoints());
+        nbrPoints.clear();
+        nbrPoints.reserve(procPatch.nPoints());
+
+        // Get all changed points in reverse order
+        const labelList& neighbPoints = procPatch.neighbPoints();
+        forAll(neighbPoints, thisPointI)
         {
-            // Get all changed points in relative addressing
-
-            DynamicList<Type> patchInfo(patch.nPoints());
-            DynamicList<label> patchPoints(patch.nPoints());
-            DynamicList<label> owner(patch.nPoints());
-            DynamicList<label> ownerIndex(patch.nPoints());
-
-            getChangedPatchPoints
-            (
-                patch,
-                patchInfo,
-                patchPoints,
-                owner,
-                ownerIndex
-            );
-
-            // Adapt for leaving domain
-            leaveDomain(patch, patch, patchPoints, patchInfo);
-
-            const processorPolyPatch& procPatch =
-                refCast<const processorPolyPatch>(patch);
-
-            if (debug)
+            label meshPointI = procPatch.meshPoints()[thisPointI];
+            if (changedPoint_[meshPointI])
             {
-                Pout<< "Processor patch " << patchI << ' ' << patch.name()
-                    << " communicating with " << procPatch.neighbProcNo()
-                    << "  Sending:" << patchInfo.size() << endl;
-            }
-
-            {
-                OPstream toNeighbour
-                (
-                    Pstream::blocking,
-                    procPatch.neighbProcNo()
-                );
-
-                toNeighbour << owner << ownerIndex << patchInfo;
+                patchInfo.append(allPointInfo_[meshPointI]);
+                thisPoints.append(thisPointI);
+                nbrPoints.append(neighbPoints[thisPointI]);
             }
         }
+
+        // Adapt for leaving domain
+        leaveDomain(procPatch, thisPoints, patchInfo);
+
+        if (debug)
+        {
+            Pout<< "Processor patch " << patchI << ' ' << procPatch.name()
+                << " communicating with " << procPatch.neighbProcNo()
+                << "  Sending:" << patchInfo.size() << endl;
+        }
+
+        UOPstream toNeighbour(procPatch.neighbProcNo(), pBufs);
+        toNeighbour << nbrPoints << patchInfo;
     }
 
 
+    pBufs.finishedSends();
+
     //
     // 2. Receive all point info on processor patches.
     //
 
-    forAll(mesh_.boundaryMesh(), patchI)
+    forAll(mesh_.globalData().processorPatches(), i)
     {
-        const polyPatch& patch = mesh_.boundaryMesh()[patchI];
+        label patchI = mesh_.globalData().processorPatches()[i];
+        const processorPolyPatch& procPatch =
+            refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchI]);
+
+        List<Type> patchInfo;
+        labelList patchPoints;
 
-        if (isA<processorPolyPatch>(patch))
         {
-            const processorPolyPatch& procPatch =
-                refCast<const processorPolyPatch>(patch);
+            UIPstream fromNeighbour(procPatch.neighbProcNo(), pBufs);
+            fromNeighbour >> patchPoints >> patchInfo;
+        }
 
-            List<Type> patchInfo;
-            labelList owner;
-            labelList ownerIndex;
+        if (debug)
+        {
+            Pout<< "Processor patch " << patchI << ' ' << procPatch.name()
+                << " communicating with " << procPatch.neighbProcNo()
+                << "  Received:" << patchInfo.size() << endl;
+        }
+
+        // Apply transform to received data for non-parallel planes
+        if (!procPatch.parallel())
+        {
+            transform(procPatch, procPatch.forwardT(), patchInfo);
+        }
+
+        // Adapt for entering domain
+        enterDomain(procPatch, patchPoints, patchInfo);
+
+        // Merge received info
+        const labelList& meshPoints = procPatch.meshPoints();
+        forAll(patchInfo, i)
+        {
+            label meshPointI = meshPoints[patchPoints[i]];
+
+            if (!allPointInfo_[meshPointI].equal(patchInfo[i], td_))
             {
-                IPstream fromNeighbour
+                updatePoint
                 (
-                    Pstream::blocking,
-                    procPatch.neighbProcNo()
+                    meshPointI,
+                    patchInfo[i],
+                    allPointInfo_[meshPointI]
                 );
-
-                fromNeighbour >> owner >> ownerIndex >> patchInfo;
             }
-
-            if (debug)
-            {
-                Pout<< "Processor patch " << patchI << ' ' << patch.name()
-                    << " communicating with " << procPatch.neighbProcNo()
-                    << "  Received:" << patchInfo.size() << endl;
-            }
-
-            // Apply transform to received data for non-parallel planes
-            if (!procPatch.parallel())
-            {
-                transform(procPatch.forwardT(), patchInfo);
-            }
-
-            updateFromPatchInfo
-            (
-                patch,
-                patch,
-                owner,
-                ownerIndex,
-                patchInfo
-            );
         }
     }
 
@@ -508,20 +411,21 @@ void Foam::PointEdgeWave<Type>::handleProcPatches()
 
     // Combine on master. Reduce operator has to handle a list and call
     // Type.updatePoint for all elements
-    combineReduce(sharedData, listUpdateOp<Type>());
+    combineReduce(sharedData, listUpdateOp<Type>(propagationTol_, td_));
 
     forAll(pd.sharedPointLabels(), i)
     {
         label meshPointI = pd.sharedPointLabels()[i];
 
-        // Retrieve my entries from the shared points
-        if (sharedData[pd.sharedPointAddr()[i]].valid())
+        // Retrieve my entries from the shared points.
+        const Type& nbrInfo = sharedData[pd.sharedPointAddr()[i]];
+
+        if (!allPointInfo_[meshPointI].equal(nbrInfo, td_))
         {
             updatePoint
             (
                 meshPointI,
-                sharedData[pd.sharedPointAddr()[i]],
-                propagationTol_,
+                nbrInfo,
                 allPointInfo_[meshPointI]
             );
         }
@@ -529,11 +433,14 @@ void Foam::PointEdgeWave<Type>::handleProcPatches()
 }
 
 
-template <class Type>
-void Foam::PointEdgeWave<Type>::handleCyclicPatches()
+template <class Type, class TrackingData>
+void Foam::PointEdgeWave<Type, TrackingData>::handleCyclicPatches()
 {
-    // 1. Send all point info on cyclic patches. Send as
-    // face label + offset in face.
+    // 1. Send all point info on cyclic patches.
+
+    DynamicList<Type> nbrInfo;
+    DynamicList<label> nbrPoints;
+    DynamicList<label> thisPoints;
 
     forAll(mesh_.boundaryMesh(), patchI)
     {
@@ -544,31 +451,44 @@ void Foam::PointEdgeWave<Type>::handleCyclicPatches()
             const cyclicPolyPatch& cycPatch =
                 refCast<const cyclicPolyPatch>(patch);
 
-            const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
+            nbrInfo.clear();
+            nbrInfo.reserve(cycPatch.nPoints());
+            nbrPoints.clear();
+            nbrPoints.reserve(cycPatch.nPoints());
+            thisPoints.clear();
+            thisPoints.reserve(cycPatch.nPoints());
 
-            DynamicList<Type> nbrInfo(nbrPatch.nPoints());
-            DynamicList<label> nbrPoints(nbrPatch.nPoints());
-            DynamicList<label> nbrOwner(nbrPatch.nPoints());
-            DynamicList<label> nbrIndex(nbrPatch.nPoints());
+            // Collect nbrPatch points that have changed
+            {
+                const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
+                const edgeList& pairs = cycPatch.coupledPoints();
+                const labelList& meshPoints = nbrPatch.meshPoints();
 
-            getChangedPatchPoints
-            (
-                nbrPatch,
-                nbrInfo,
-                nbrPoints,
-                nbrOwner,
-                nbrIndex
-            );
+                forAll(pairs, pairI)
+                {
+                    label thisPointI = pairs[pairI][0];
+                    label nbrPointI = pairs[pairI][1];
+                    label meshPointI = meshPoints[nbrPointI];
+
+                    if (changedPoint_[meshPointI])
+                    {
+                        nbrInfo.append(allPointInfo_[meshPointI]);
+                        nbrPoints.append(nbrPointI);
+                        thisPoints.append(thisPointI);
+                    }
+                }
+
+                // nbr : adapt for leaving domain
+                leaveDomain(nbrPatch, nbrPoints, nbrInfo);
+            }
 
-            // nbr : adapt for leaving domain
-            leaveDomain(nbrPatch, nbrPatch, nbrPoints, nbrInfo);
 
             // Apply rotation for non-parallel planes
 
             if (!cycPatch.parallel())
             {
                 // received data from half1
-                transform(cycPatch.forwardT(), nbrInfo);
+                transform(cycPatch, cycPatch.forwardT(), nbrInfo);
             }
 
             if (debug)
@@ -578,19 +498,24 @@ void Foam::PointEdgeWave<Type>::handleCyclicPatches()
                     << endl;
             }
 
-            // Half1: update with data from halfB
-            updateFromPatchInfo
-            (
-                cycPatch,
-                cycPatch,
-                nbrOwner,
-                nbrIndex,
-                nbrInfo
-            );
+            // Adapt for entering domain
+            enterDomain(cycPatch, thisPoints, nbrInfo);
 
-            if (debug)
+            // Merge received info
+            const labelList& meshPoints = cycPatch.meshPoints();
+            forAll(nbrInfo, i)
             {
-                //checkCyclic(patch);
+                label meshPointI = meshPoints[thisPoints[i]];
+
+                if (!allPointInfo_[meshPointI].equal(nbrInfo[i], td_))
+                {
+                    updatePoint
+                    (
+                        meshPointI,
+                        nbrInfo[i],
+                        allPointInfo_[meshPointI]
+                    );
+                }
             }
         }
     }
@@ -601,22 +526,24 @@ void Foam::PointEdgeWave<Type>::handleCyclicPatches()
 
 // Iterate, propagating changedPointsInfo across mesh, until no change (or
 // maxIter reached). Initial point values specified.
-template <class Type>
-Foam::PointEdgeWave<Type>::PointEdgeWave
+template <class Type, class TrackingData>
+Foam::PointEdgeWave<Type, TrackingData>::PointEdgeWave
 (
     const polyMesh& mesh,
     const labelList& changedPoints,
     const List<Type>& changedPointsInfo,
 
-    List<Type>& allPointInfo,
-    List<Type>& allEdgeInfo,
+    UList<Type>& allPointInfo,
+    UList<Type>& allEdgeInfo,
 
-    const label maxIter
+    const label maxIter,
+    TrackingData& td
 )
 :
     mesh_(mesh),
     allPointInfo_(allPointInfo),
     allEdgeInfo_(allEdgeInfo),
+    td_(td),
     changedPoint_(mesh_.nPoints(), false),
     changedPoints_(mesh_.nPoints()),
     nChangedPoints_(0),
@@ -632,7 +559,7 @@ Foam::PointEdgeWave<Type>::PointEdgeWave
     {
         FatalErrorIn
         (
-            "PointEdgeWave<Type>::PointEdgeWave"
+            "PointEdgeWave<Type, TrackingData>::PointEdgeWave"
             "(const polyMesh&, const labelList&, const List<Type>,"
             " List<Type>&, List<Type>&, const label maxIter)"
         )   << "size of pointInfo work array is not equal to the number"
@@ -645,7 +572,7 @@ Foam::PointEdgeWave<Type>::PointEdgeWave
     {
         FatalErrorIn
         (
-            "PointEdgeWave<Type>::PointEdgeWave"
+            "PointEdgeWave<Type, TrackingData>::PointEdgeWave"
             "(const polyMesh&, const labelList&, const List<Type>,"
             " List<Type>&, List<Type>&, const label maxIter)"
         )   << "size of edgeInfo work array is not equal to the number"
@@ -671,7 +598,7 @@ Foam::PointEdgeWave<Type>::PointEdgeWave
     {
         FatalErrorIn
         (
-            "PointEdgeWave<Type>::PointEdgeWave"
+            "PointEdgeWave<Type, TrackingData>::PointEdgeWave"
             "(const polyMesh&, const labelList&, const List<Type>,"
             " List<Type>&, List<Type>&, const label maxIter)"
         )   << "Maximum number of iterations reached. Increase maxIter." << endl
@@ -683,17 +610,19 @@ Foam::PointEdgeWave<Type>::PointEdgeWave
 }
 
 
-template <class Type>
-Foam::PointEdgeWave<Type>::PointEdgeWave
+template <class Type, class TrackingData>
+Foam::PointEdgeWave<Type, TrackingData>::PointEdgeWave
 (
     const polyMesh& mesh,
-    List<Type>& allPointInfo,
-    List<Type>& allEdgeInfo
+    UList<Type>& allPointInfo,
+    UList<Type>& allEdgeInfo,
+    TrackingData& td
 )
 :
     mesh_(mesh),
     allPointInfo_(allPointInfo),
     allEdgeInfo_(allEdgeInfo),
+    td_(td),
     changedPoint_(mesh_.nPoints(), false),
     changedPoints_(mesh_.nPoints()),
     nChangedPoints_(0),
@@ -709,31 +638,31 @@ Foam::PointEdgeWave<Type>::PointEdgeWave
 
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
-template <class Type>
-Foam::PointEdgeWave<Type>::~PointEdgeWave()
+template <class Type, class TrackingData>
+Foam::PointEdgeWave<Type, TrackingData>::~PointEdgeWave()
 {}
 
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 
-template <class Type>
-Foam::label Foam::PointEdgeWave<Type>::getUnsetPoints() const
+template <class Type, class TrackingData>
+Foam::label Foam::PointEdgeWave<Type, TrackingData>::getUnsetPoints() const
 {
     return nUnvisitedPoints_;
 }
 
 
-template <class Type>
-Foam::label Foam::PointEdgeWave<Type>::getUnsetEdges() const
+template <class Type, class TrackingData>
+Foam::label Foam::PointEdgeWave<Type, TrackingData>::getUnsetEdges() const
 {
     return nUnvisitedEdges_;
 }
 
 
 // Copy point information into member data
-template <class Type>
-void Foam::PointEdgeWave<Type>::setPointInfo
+template <class Type, class TrackingData>
+void Foam::PointEdgeWave<Type, TrackingData>::setPointInfo
 (
     const labelList& changedPoints,
     const List<Type>& changedPointsInfo
@@ -743,13 +672,13 @@ void Foam::PointEdgeWave<Type>::setPointInfo
     {
         label pointI = changedPoints[changedPointI];
 
-        bool wasValid = allPointInfo_[pointI].valid();
+        bool wasValid = allPointInfo_[pointI].valid(td_);
 
         // Copy info for pointI
         allPointInfo_[pointI] = changedPointsInfo[changedPointI];
 
         // Maintain count of unset points
-        if (!wasValid && allPointInfo_[pointI].valid())
+        if (!wasValid && allPointInfo_[pointI].valid(td_))
         {
             --nUnvisitedPoints_;
         }
@@ -766,8 +695,8 @@ void Foam::PointEdgeWave<Type>::setPointInfo
 
 
 // Propagate information from edge to point. Return number of points changed.
-template <class Type>
-Foam::label Foam::PointEdgeWave<Type>::edgeToPoint()
+template <class Type, class TrackingData>
+Foam::label Foam::PointEdgeWave<Type, TrackingData>::edgeToPoint()
 {
     for
     (
@@ -780,7 +709,7 @@ Foam::label Foam::PointEdgeWave<Type>::edgeToPoint()
 
         if (!changedEdge_[edgeI])
         {
-            FatalErrorIn("PointEdgeWave<Type>::edgeToPoint()")
+            FatalErrorIn("PointEdgeWave<Type, TrackingData>::edgeToPoint()")
                 << "edge " << edgeI
                 << " not marked as having been changed" << nl
                 << "This might be caused by multiple occurences of the same"
@@ -797,14 +726,13 @@ Foam::label Foam::PointEdgeWave<Type>::edgeToPoint()
         {
             Type& currentWallInfo = allPointInfo_[e[eI]];
 
-            if (currentWallInfo != neighbourWallInfo)
+            if (!currentWallInfo.equal(neighbourWallInfo, td_))
             {
                 updatePoint
                 (
                     e[eI],
                     edgeI,
                     neighbourWallInfo,
-                    propagationTol_,
                     currentWallInfo
                 );
             }
@@ -843,8 +771,8 @@ Foam::label Foam::PointEdgeWave<Type>::edgeToPoint()
 
 
 // Propagate information from point to edge. Return number of edges changed.
-template <class Type>
-Foam::label Foam::PointEdgeWave<Type>::pointToEdge()
+template <class Type, class TrackingData>
+Foam::label Foam::PointEdgeWave<Type, TrackingData>::pointToEdge()
 {
     const labelListList& pointEdges = mesh_.pointEdges();
 
@@ -859,7 +787,7 @@ Foam::label Foam::PointEdgeWave<Type>::pointToEdge()
 
         if (!changedPoint_[pointI])
         {
-            FatalErrorIn("PointEdgeWave<Type>::pointToEdge()")
+            FatalErrorIn("PointEdgeWave<Type, TrackingData>::pointToEdge()")
                 << "Point " << pointI
                 << " not marked as having been changed" << nl
                 << "This might be caused by multiple occurences of the same"
@@ -877,14 +805,13 @@ Foam::label Foam::PointEdgeWave<Type>::pointToEdge()
 
             Type& currentWallInfo = allEdgeInfo_[edgeI];
 
-            if (currentWallInfo != neighbourWallInfo)
+            if (!currentWallInfo.equal(neighbourWallInfo, td_))
             {
                 updateEdge
                 (
                     edgeI,
                     pointI,
                     neighbourWallInfo,
-                    propagationTol_,
                     currentWallInfo
                 );
             }
@@ -912,8 +839,11 @@ Foam::label Foam::PointEdgeWave<Type>::pointToEdge()
 
 
 // Iterate
-template <class Type>
-Foam::label Foam::PointEdgeWave<Type>::iterate(const label maxIter)
+template <class Type, class TrackingData>
+Foam::label Foam::PointEdgeWave<Type, TrackingData>::iterate
+(
+    const label maxIter
+)
 {
     if (nCyclicPatches_ > 0)
     {
@@ -975,4 +905,5 @@ Foam::label Foam::PointEdgeWave<Type>::iterate(const label maxIter)
     return iter;
 }
 
+
 // ************************************************************************* //
diff --git a/src/meshTools/PointEdgeWave/PointEdgeWave.H b/src/meshTools/PointEdgeWave/PointEdgeWave.H
index 9423fbbd71..c3e70db195 100644
--- a/src/meshTools/PointEdgeWave/PointEdgeWave.H
+++ b/src/meshTools/PointEdgeWave/PointEdgeWave.H
@@ -58,13 +58,9 @@ SourceFiles
 #ifndef PointEdgeWave_H
 #define PointEdgeWave_H
 
-#include "label.H"
 #include "boolList.H"
 #include "scalarField.H"
-#include "pointFields.H"
-#include "tensor.H"
-#include "primitivePatch.H"
-#include "PtrList.H"
+#include "tensorField.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -73,7 +69,7 @@ namespace Foam
 
 // Forward declaration of classes
 class polyMesh;
-
+class polyPatch;
 
 /*---------------------------------------------------------------------------*\
                         Class PointEdgeWaveName Declaration
@@ -86,7 +82,7 @@ TemplateName(PointEdgeWave);
                            Class PointEdgeWave Declaration
 \*---------------------------------------------------------------------------*/
 
-template <class Type>
+template <class Type, class TrackingData = int>
 class PointEdgeWave
 :
     public PointEdgeWaveName
@@ -98,6 +94,10 @@ class PointEdgeWave
         //  up to Type implementation)
         static scalar propagationTol_;
 
+        //- Used as default trackdata value to satisfy default template
+        //  argument.
+        static label dummyTrackData_;
+
 
     // Private data
 
@@ -105,10 +105,13 @@ class PointEdgeWave
         const polyMesh& mesh_;
 
         //- Wall information for all points
-        List<Type>& allPointInfo_;
+        UList<Type>& allPointInfo_;
 
         //- Information on all mesh edges
-        List<Type>& allEdgeInfo_;
+        UList<Type>& allEdgeInfo_;
+
+        //- Additional data to be passed into container
+        TrackingData& td_;
 
         //- Has point changed
         boolList changedPoint_;
@@ -127,9 +130,6 @@ class PointEdgeWave
         //- Number of cyclic patches
         label nCyclicPatches_;
 
-        //- For every cyclic patch two primitivePatches
-        PtrList<primitivePatch> cycHalves_;
-
         //- Number of evaluations
         label nEvals_;
 
@@ -140,15 +140,10 @@ class PointEdgeWave
 
     // Private Member Functions
 
-        //- Add value to all elements of labelList
-        static void offset(const label val, labelList& elems);
-
-
         //- Adapt pointInfo for leaving domain
         void leaveDomain
         (
-            const polyPatch& meshPatch,
-            const primitivePatch& patch,
+            const polyPatch&,
             const List<label>& patchPointLabels,
             List<Type>& pointInfo
         ) const;
@@ -156,8 +151,7 @@ class PointEdgeWave
         //- Adapt pointInfo for entering domain
         void enterDomain
         (
-            const polyPatch& meshPatch,
-            const primitivePatch& patch,
+            const polyPatch&,
             const List<label>& patchPointLabels,
             List<Type>& pointInfo
         ) const;
@@ -165,6 +159,7 @@ class PointEdgeWave
         //- Transform. Implementation referred to Type
         void transform
         (
+            const polyPatch& patch,
             const tensorField& rotTensor,
             List<Type>& pointInfo
         ) const;
@@ -176,7 +171,6 @@ class PointEdgeWave
             const label pointI,
             const label neighbourEdgeI,
             const Type& neighbourInfo,
-            const scalar tol,
             Type& pointInfo
         );
 
@@ -186,7 +180,6 @@ class PointEdgeWave
         (
             const label pointI,
             const Type& neighbourInfo,
-            const scalar tol,
             Type& pointInfo
         );
 
@@ -197,7 +190,6 @@ class PointEdgeWave
             const label edgeI,
             const label neighbourPointI,
             const Type& neighbourInfo,
-            const scalar tol,
             Type& edgeInfo
         );
 
@@ -207,26 +199,6 @@ class PointEdgeWave
             template <class PatchType>
             label countPatchType() const;
 
-            //- Get info on patch points
-            void getChangedPatchPoints
-            (
-                const primitivePatch& patch,
-                DynamicList<Type>& patchInfo,
-                DynamicList<label>& patchPoints,
-                DynamicList<label>& owner,
-                DynamicList<label>& ownerIndex
-            ) const;
-
-            //- Merge data from patch into overall data
-            void updateFromPatchInfo
-            (
-                const polyPatch& meshPatch,
-                const primitivePatch& patch,
-                const labelList& owner,
-                const labelList& ownerIndex,
-                List<Type>& patchInfo
-            );
-
             //- Merge data from across processor boundaries
             void handleProcPatches();
 
@@ -270,9 +242,10 @@ public:
             const polyMesh& mesh,
             const labelList& initialPoints,
             const List<Type>& initialPointsInfo,
-            List<Type>& allPointInfo,
-            List<Type>& allEdgeInfo,
-            const label maxIter
+            UList<Type>& allPointInfo,
+            UList<Type>& allEdgeInfo,
+            const label maxIter,
+            TrackingData& td = dummyTrackData_
         );
 
         //- Construct from mesh. Use setPointInfo and iterate() to do
@@ -280,8 +253,9 @@ public:
         PointEdgeWave
         (
             const polyMesh& mesh,
-            List<Type>& allPointInfo,
-            List<Type>& allEdgeInfo
+            UList<Type>& allPointInfo,
+            UList<Type>& allEdgeInfo,
+            TrackingData& td = dummyTrackData_
         );
 
 
@@ -291,18 +265,24 @@ public:
 
     // Member Functions
 
-        //- Get allPointInfo
-        const List<Type>& allPointInfo() const
+        //- Access allPointInfo
+        UList<Type>& allPointInfo() const
         {
             return allPointInfo_;
         }
 
-        //- Get allEdgeInfo
-        const List<Type>& allEdgeInfo() const
+        //- Access allEdgeInfo
+        UList<Type>& allEdgeInfo() const
         {
             return allEdgeInfo_;
         }
 
+        //- Additional data to be passed into container
+        const TrackingData& data() const
+        {
+            return td_;
+        }
+
         //- Get number of unvisited edges, i.e. edges that were not (yet)
         //  reached from walking across mesh. This can happen from
         //  - not enough iterations done
@@ -340,19 +320,29 @@ public:
 \*---------------------------------------------------------------------------*/
 
 //- List update operation
-template <class Type>
+template <class Type, class TrackingData = int>
 class listUpdateOp
 {
+    //- Additional data to be passed into container
+
+    const scalar tol_;
+
+    TrackingData& td_;
 
 public:
+    listUpdateOp(const scalar tol, TrackingData& td)
+    :
+        tol_(tol),
+        td_(td)
+    {}
 
     void operator()(List<Type>& x, const List<Type>& y) const
     {
         forAll(x, i)
         {
-            if (y[i].valid())
+            if (y[i].valid(td_))
             {
-                x[i].updatePoint(y[i], PointEdgeWave<Type>::propagationTol());
+                x[i].updatePoint(y[i], tol_, td_);
             }
         }
     }
diff --git a/src/meshTools/PointEdgeWave/pointEdgePoint.H b/src/meshTools/PointEdgeWave/pointEdgePoint.H
index a876ab2e87..81f4019bd4 100644
--- a/src/meshTools/PointEdgeWave/pointEdgePoint.H
+++ b/src/meshTools/PointEdgeWave/pointEdgePoint.H
@@ -73,19 +73,23 @@ class pointEdgePoint
         //- 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 pointEdgePoint& w2,
-            const scalar tol
+            const scalar tol,
+            TrackingData& td
         );
 
         //- Combine current with w2. Update distSqr, origin if w2 has smaller
         //  quantities and returns true.
+        template<class TrackingData>
         inline bool update
         (
             const pointEdgePoint& w2,
-            const scalar tol
+            const scalar tol,
+            TrackingData& td
         );
 
 
@@ -97,7 +101,7 @@ public:
         inline pointEdgePoint();
 
         //- Construct from origin, distance
-        inline pointEdgePoint(const point& origin, const scalar distSqr);
+        inline pointEdgePoint(const point&, const scalar);
 
         //- Construct as copy
         inline pointEdgePoint(const pointEdgePoint&);
@@ -116,76 +120,102 @@ 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.
-            inline bool sameGeometry(const pointEdgePoint&, const scalar tol)
-             const;
+            template<class TrackingData>
+            inline bool sameGeometry
+            (
+                const pointEdgePoint&,
+                const scalar tol,
+                TrackingData& td
+            ) const;
 
             //- Convert origin to relative vector to leaving point
             //  (= point coordinate)
+            template<class TrackingData>
             inline void leaveDomain
             (
                 const polyPatch& patch,
                 const label patchPointI,
-                const point& pos
+                const point& pos,
+                TrackingData& td
             );
 
             //- Convert relative origin to absolute by adding entering point
+            template<class TrackingData>
             inline void enterDomain
             (
                 const polyPatch& patch,
                 const label patchPointI,
-                const point& pos
+                const point& pos,
+                TrackingData& td
             );
 
             //- Apply rotation matrix to origin
-            inline void transform(const tensor& rotTensor);
+            template<class TrackingData>
+            inline void transform
+            (
+                const tensor& rotTensor,
+                TrackingData& td
+            );
 
             //- Influence of edge on point
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const polyMesh& mesh,
                 const label pointI,
                 const label edgeI,
                 const pointEdgePoint& edgeInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same point.
             //  Merge new and old info.
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const polyMesh& mesh,
                 const label pointI,
                 const pointEdgePoint& newPointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of different value on same point.
             //  No information about current position whatsoever.
+            template<class TrackingData>
             inline bool updatePoint
             (
                 const pointEdgePoint& newPointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
             //- Influence of point on edge.
+            template<class TrackingData>
             inline bool updateEdge
             (
                 const polyMesh& mesh,
                 const label edgeI,
                 const label pointI,
                 const pointEdgePoint& pointInfo,
-                const scalar tol
+                const scalar tol,
+                TrackingData& td
             );
 
+            //- Same (like operator==)
+            template<class TrackingData>
+            inline bool equal(const pointEdgePoint&, TrackingData& td) const;
+
 
     // Member Operators
 
-        //Note: Used to determine whether to call update.
+        // Needed for List IO
         inline bool operator==(const pointEdgePoint&) const;
-
         inline bool operator!=(const pointEdgePoint&) const;
 
 
diff --git a/src/meshTools/PointEdgeWave/pointEdgePointI.H b/src/meshTools/PointEdgeWave/pointEdgePointI.H
index 38247e1583..43bdd088a2 100644
--- a/src/meshTools/PointEdgeWave/pointEdgePointI.H
+++ b/src/meshTools/PointEdgeWave/pointEdgePointI.H
@@ -29,16 +29,18 @@ License
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Update this with w2 if w2 nearer to pt.
+template<class TrackingData>
 inline bool Foam::pointEdgePoint::update
 (
     const point& pt,
     const pointEdgePoint& 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;
@@ -72,13 +74,15 @@ inline bool Foam::pointEdgePoint::update
 
 
 // Update this with w2 (information on same point)
+template<class TrackingData>
 inline bool Foam::pointEdgePoint::update
 (
     const pointEdgePoint& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    if (!valid())
+    if (!valid(td))
     {
         // current not yet set so use any value
         distSqr_ = w2.distSqr();
@@ -155,17 +159,20 @@ inline Foam::scalar Foam::pointEdgePoint::distSqr() const
 }
 
 
-inline bool Foam::pointEdgePoint::valid() const
+template<class TrackingData>
+inline bool Foam::pointEdgePoint::valid(TrackingData& td) const
 {
     return origin_ != point::max;
 }
 
 
 // Checks for cyclic points
+template<class TrackingData>
 inline bool Foam::pointEdgePoint::sameGeometry
 (
     const pointEdgePoint& w2,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 ) const
 {
     scalar diff = Foam::mag(distSqr() - w2.distSqr());
@@ -188,18 +195,25 @@ inline bool Foam::pointEdgePoint::sameGeometry
 }
 
 
+template<class TrackingData>
 inline void Foam::pointEdgePoint::leaveDomain
 (
     const polyPatch& patch,
     const label patchPointI,
-    const point& coord
+    const point& coord,
+    TrackingData& td
 )
 {
     origin_ -= coord;
 }
 
 
-inline void Foam::pointEdgePoint::transform(const tensor& rotTensor)
+template<class TrackingData>
+inline void Foam::pointEdgePoint::transform
+(
+    const tensor& rotTensor,
+    TrackingData& td
+)
 {
     origin_ = Foam::transform(rotTensor, origin_);
 }
@@ -207,11 +221,13 @@ inline void Foam::pointEdgePoint::transform(const tensor& rotTensor)
 
 // Update absolute geometric quantities. Note that distance (distSqr_)
 // is not affected by leaving/entering domain.
+template<class TrackingData>
 inline void Foam::pointEdgePoint::enterDomain
 (
     const polyPatch& patch,
     const label patchPointI,
-    const point& coord
+    const point& coord,
+    TrackingData& td
 )
 {
     // back to absolute form
@@ -220,78 +236,74 @@ inline void Foam::pointEdgePoint::enterDomain
 
 
 // Update this with information from connected edge
+template<class TrackingData>
 inline bool Foam::pointEdgePoint::updatePoint
 (
     const polyMesh& mesh,
     const label pointI,
     const label edgeI,
     const pointEdgePoint& edgeInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return
-        update
-        (
-            mesh.points()[pointI],
-            edgeInfo,
-            tol
-        );
+    return update(mesh.points()[pointI], edgeInfo, tol, td);
 }
 
 
 // Update this with new information on same point
+template<class TrackingData>
 inline bool Foam::pointEdgePoint::updatePoint
 (
     const polyMesh& mesh,
     const label pointI,
     const pointEdgePoint& newPointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return
-        update
-        (
-            mesh.points()[pointI],
-            newPointInfo,
-            tol
-        );
+    return update(mesh.points()[pointI], newPointInfo, tol, td);
 }
 
 
 // Update this with new information on same point. No extra information.
+template<class TrackingData>
 inline bool Foam::pointEdgePoint::updatePoint
 (
     const pointEdgePoint& newPointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    return update(newPointInfo, tol);
+    return update(newPointInfo, tol, td);
 }
 
 
 // Update this with information from connected point
+template<class TrackingData>
 inline bool Foam::pointEdgePoint::updateEdge
 (
     const polyMesh& mesh,
     const label edgeI,
     const label pointI,
     const pointEdgePoint& pointInfo,
-    const scalar tol
+    const scalar tol,
+    TrackingData& td
 )
 {
-    const pointField& points = mesh.points();
-
     const edge& e = mesh.edges()[edgeI];
+    return update(e.centre(mesh.points()), pointInfo, tol, td);
+}
 
-    const point edgeMid(0.5*(points[e[0]] + points[e[1]]));
 
-    return
-        update
-        (
-            edgeMid,
-            pointInfo,
-            tol
-        );
+template <class TrackingData>
+inline bool Foam::pointEdgePoint::equal
+(
+    const pointEdgePoint& rhs,
+    TrackingData& td
+) const
+{
+    return operator==(rhs);
 }