diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtra.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtra.C
new file mode 100644
index 0000000000..7ab4c43084
--- /dev/null
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtra.C
@@ -0,0 +1,173 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "PrimitivePatchExtra.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+// Construct from components
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+PrimitivePatchExtra
+(
+    const ListType<Face>& faces,
+    const pointField& points
+)
+:
+    PrimitivePatch<Face, ListType, PointField, PointType>(faces, points),
+    sortedEdgeFacesPtr_(NULL),
+    edgeOwnerPtr_(NULL)
+{}
+
+
+// Construct as copy
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+PrimitivePatchExtra
+(
+    const PrimitivePatchExtra<Face, ListType, PointField, PointType>& pp
+)
+:
+    PrimitivePatch<Face, ListType, PointField, PointType>(pp),
+    sortedEdgeFacesPtr_(NULL),
+    edgeOwnerPtr_(NULL)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+~PrimitivePatchExtra()
+{
+    clearOut();
+}
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+void Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+clearOut()
+{
+    PrimitivePatch<Face, ListType, PointField, PointType>::clearOut();
+    clearTopology();
+}
+
+
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+void Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+clearTopology()
+{
+    PrimitivePatch<Face, ListType, PointField, PointType>::clearTopology();
+    deleteDemandDrivenData(sortedEdgeFacesPtr_);
+    deleteDemandDrivenData(edgeOwnerPtr_);
+}
+
+
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+const Foam::labelListList&
+Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+sortedEdgeFaces() const
+{
+    if (!sortedEdgeFacesPtr_)
+    {
+        calcSortedEdgeFaces();
+    }
+
+    return *sortedEdgeFacesPtr_;
+}
+
+
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+const Foam::labelList&
+Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+edgeOwner() const
+{
+    if (!edgeOwnerPtr_)
+    {
+        calcEdgeOwner();
+    }
+
+    return *edgeOwnerPtr_;
+}
+
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "PrimitivePatchExtraAddressing.C"
+#include "PrimitivePatchExtraCleanup.C"
+#include "PrimitivePatchExtraSearch.C"
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtra.H b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtra.H
new file mode 100644
index 0000000000..094dfb7340
--- /dev/null
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtra.H
@@ -0,0 +1,211 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    Foam::PrimitivePatchExtra
+
+Description
+    PrimitivePatch with some extra functionality.
+
+SourceFiles
+    PrimitivePatchExtra.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef PrimitivePatchExtra_H
+#define PrimitivePatchExtra_H
+
+#include "PrimitivePatch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class PrimitivePatchExtra Declaration
+\*---------------------------------------------------------------------------*/
+
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType=point
+>
+class PrimitivePatchExtra
+:
+    public PrimitivePatch<Face, ListType, PointField, PointType>
+{
+
+public:
+
+    // Public typedefs
+
+        typedef Face FaceType;
+        typedef ListType<Face> FaceListType;
+        typedef PointField PointFieldType;
+
+private:
+
+    // Private typedefs
+    typedef PrimitivePatch<Face, ListType, PointField, PointType> TemplateType;
+
+    // Private data
+
+    // Demand driven private data
+
+        //- Edge-face addressing (sorted)
+        mutable labelListList* sortedEdgeFacesPtr_;
+
+        //- Label of face that 'owns' edge (i.e. e.vec() is righthanded walk
+        //  along face)
+        mutable labelList* edgeOwnerPtr_;
+
+
+    // Private Member Functions
+
+        //- Calculate sorted edgeFaces
+        void calcSortedEdgeFaces() const;
+
+        //- Calculate owner
+        void calcEdgeOwner() const;
+
+protected:
+    // Protected Member Functions
+
+    // Edit
+
+        //- Fill faceZone with currentZone for every face reachable
+        //  from faceI without crossing edge marked in borderEdge.
+        //  Note: faceZone has to be sized nFaces before calling this fun.
+        void markZone
+        (
+            const boolList& borderEdge,
+            const label faceI,
+            const label currentZone,
+            labelList& faceZone
+        ) const;
+
+        //- (size and) fills faceZone with zone of face.
+        //  Zone is area reachable by edge crossing without crossing borderEdge
+        //  (bool for every edge in surface). Returns number of zones.
+        label markZones
+        (
+            const boolList& borderEdge,
+            labelList& faceZone
+        ) const;
+
+        //- Determine the mapping for a sub mesh.
+        //  Only include faces for which boolList entry is true
+        //  Sets: pointMap: from new to old localPoints
+        //        faceMap:  new to old faces
+        void subsetMap
+        (
+            const boolList& include,
+            labelList& pointMap,
+            labelList& faceMap
+        ) const;
+
+public:
+
+    // Constructors
+
+        //- Construct from components
+        PrimitivePatchExtra
+        (
+            const ListType<Face>& faces,
+            const pointField& points
+        );
+
+        //- Construct as copy
+        PrimitivePatchExtra
+        (
+            const PrimitivePatchExtra<Face, ListType, PointField, PointType>&
+        );
+
+
+    // Destructor
+
+        virtual ~PrimitivePatchExtra();
+
+        void clearOut();
+
+        void clearTopology();
+
+
+    // Member Functions
+
+    // Access functions for demand driven data
+
+    // Topological data; no mesh required.
+
+        //- Return edge-face addressing sorted
+        //  (for edges with more than 2 faces) according to the
+        //  angle around the edge.
+        //  Orientation is anticlockwise looking from
+        //  edge.vec(localPoints())
+        const labelListList& sortedEdgeFaces() const;
+
+        //- If 2 face neighbours: label of face where ordering of edge
+        //  is consistent with righthand walk.
+        //  If 1 neighbour: label of only face.
+        //  If >2 neighbours: undetermined.
+        const labelList& edgeOwner() const;
+
+
+    // Addressing into mesh
+
+
+    // Other patch operations
+
+
+    // Check
+
+        //- Check triply (or more) connected edges. 
+        //  Return list of faces sharing these edges.
+        void checkEdges(const bool verbose) const;
+
+        //- Check orientation (normals) and normals of neighbouring faces
+        boolList checkOrientation(const bool verbose) const;
+
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "PrimitivePatchExtra.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraAddressing.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraAddressing.C
new file mode 100644
index 0000000000..1773014a4d
--- /dev/null
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraAddressing.C
@@ -0,0 +1,210 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Description
+    Contains fix for PrimitivePatch addressing (which doesn't work if surface
+    is non-manifold). Should be moved into PrimitivePatch.
+
+\*---------------------------------------------------------------------------*/
+
+#include "PrimitivePatchExtra.H"
+#include "HashTable.H"
+#include "SortableList.H"
+#include "transform.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+void Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+calcSortedEdgeFaces() const
+{
+    if (sortedEdgeFacesPtr_)
+    {
+        FatalErrorIn("PrimitivePatchExtra<Face, ListType, PointField>::calcSortedEdgeFaces()")
+            << "sortedEdgeFacesPtr_ already set"
+            << abort(FatalError);
+    }
+
+    const labelListList& eFaces = TemplateType::edgeFaces();
+    const edgeList& edgeLst = TemplateType::edges();
+    const pointField& locPointLst = TemplateType::localPoints();
+    const List<FaceType>& locFaceLst = TemplateType::localFaces();
+
+    // create the lists for the various results. (resized on completion)
+    sortedEdgeFacesPtr_ = new labelListList(eFaces.size());
+    labelListList& sortedEdgeFaces = *sortedEdgeFacesPtr_;
+
+    forAll(eFaces, edgeI)
+    {
+        const labelList& myFaceNbs = eFaces[edgeI];
+
+        if (myFaceNbs.size() > 2)
+        {
+            // Get point on edge and normalized direction of edge (= e2 base
+            // of our coordinate system)
+            const edge& e = edgeLst[edgeI];
+
+            const point& edgePt = locPointLst[e.start()];
+
+            vector e2 = e.vec(locPointLst);
+            e2 /= mag(e2) + VSMALL;
+
+            // Get opposite vertex for 0th face
+            const Face& f = locFaceLst[myFaceNbs[0]];
+
+            label fp0 = findIndex(f, e[0]);
+            label fp1 = f.fcIndex(fp0);
+            label vertI = (f[fp1] != e[1] ? f[fp1] : f.fcIndex(fp1));
+
+            // Get vector normal both to e2 and to edge from opposite vertex
+            // to edge (will be x-axis of our coordinate system)
+            vector e0 = e2 ^ (locPointLst[vertI] - edgePt);
+            e0 /= mag(e0) + VSMALL;
+
+            // Get y-axis of coordinate system
+            vector e1 = e2 ^ e0;
+
+            SortableList<scalar> faceAngles(myFaceNbs.size());
+
+            // e0 is reference so angle is 0
+            faceAngles[0] = 0;
+
+            for (label nbI = 1; nbI < myFaceNbs.size(); nbI++)
+            {
+                // Get opposite vertex
+                const FaceType& f = locFaceLst[myFaceNbs[nbI]];
+                label fp0 = findIndex(f, e[0]);
+                label fp1 = f.fcIndex(fp0);
+                label vertI = (f[fp1] != e[1] ? f[fp1] : f.fcIndex(fp1));
+
+                vector vec = e2 ^ (locPointLst[vertI] - edgePt);
+                vec /= mag(vec) + VSMALL;
+
+                faceAngles[nbI] = pseudoAngle
+                (
+                    e0,
+                    e1,
+                    vec
+                );
+            }
+
+            faceAngles.sort();
+
+            sortedEdgeFaces[edgeI] = IndirectList<label>
+            (
+                myFaceNbs,
+                faceAngles.indices()
+            );
+        }
+        else
+        {
+            // No need to sort. Just copy.
+            sortedEdgeFaces[edgeI] = myFaceNbs;
+        }
+    }
+}
+
+
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+void Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+calcEdgeOwner() const
+{
+    if (edgeOwnerPtr_)
+    {
+        FatalErrorIn("PrimitivePatchExtra<Face, ListType, PointField>::calcEdgeOwner()")
+            << "edgeOwnerPtr_ already set"
+            << abort(FatalError);
+    }
+
+    // create the owner list
+    edgeOwnerPtr_ = new labelList
+    (
+        TemplateType::nEdges()
+    );
+    labelList& edgeOwner = *edgeOwnerPtr_;
+
+    const edgeList& edgeLst = TemplateType::edges();
+    const labelListList& eFaces = TemplateType::edgeFaces();
+    const List<FaceType>& locFaceLst = TemplateType::localFaces();
+
+
+    forAll(edgeLst, edgeI)
+    {
+        const edge& e = edgeLst[edgeI];
+
+        const labelList& myFaces = eFaces[edgeI];
+
+        if (myFaces.size() == 1)
+        {
+            edgeOwner[edgeI] = myFaces[0];
+        }
+        else
+        {
+            // Find the first face whose vertices are aligned with the edge.
+            // (in case of multiply connected edge the best we can do)
+            edgeOwner[edgeI] = -1;
+
+            forAll(myFaces, i)
+            {
+                const FaceType& f = locFaceLst[myFaces[i]];
+
+                if (f.findEdge(e) > 0)
+                {
+                    edgeOwner[edgeI] = myFaces[i];
+                    break;
+                }
+            }
+
+            if (edgeOwner[edgeI] == -1)
+            {
+                FatalErrorIn("PrimitivePatchExtra<Face, ListType, PointField>::calcEdgeOwner()")
+                    << "Edge " << edgeI << " vertices:" << e
+                    << " is used by faces " << myFaces
+                    << " vertices:"
+                    << IndirectList<FaceType>(locFaceLst, myFaces)()
+                    << " none of which use the edge vertices in the same order"
+                    << nl << "I give up" << abort(FatalError);
+            }
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraCleanup.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraCleanup.C
new file mode 100644
index 0000000000..433802ae72
--- /dev/null
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraCleanup.C
@@ -0,0 +1,248 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "PrimitivePatchExtra.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+// Check/fix edges with more than two faces
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+void Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+checkEdges
+(
+    const bool verbose
+) const
+{
+    const labelListList& eFaces =
+       PrimitivePatch<Face, ListType, PointField>::edgeFaces();
+
+    const edgeList& edgeLst =
+        PrimitivePatch<Face, ListType, PointField>::edges();
+
+    forAll (eFaces, edgeI)
+    {
+        const labelList& myFaces = eFaces[edgeI];
+
+        if (myFaces.size() == 0)
+        {
+            FatalErrorIn("PrimitivePatchExtra::checkEdges(bool verbose)")
+                << "Edge " << edgeI << " with vertices " << edgeLst[edgeI]
+                << " has no edgeFaces"
+                << exit(FatalError);
+        }
+        else if (myFaces.size() > 2)
+        {
+            WarningIn
+            (
+                "PrimitivePatchExtra::checkEdges(bool verbose)"
+            )   << "Edge " << edgeI << " with vertices " << edgeLst[edgeI]
+                << " has more than 2 faces connected to it : " << myFaces
+                << endl;
+        }
+    }
+}
+
+
+// Check normals and orientation
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+Foam::boolList
+Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+checkOrientation
+(
+    const bool verbose
+) const
+{
+    const edgeList& edgeLst = TemplateType::edges();
+    const labelListList& faceEs = TemplateType::faceEdges();
+    const List<FaceType>& faceLst = TemplateType::faces();
+    const label numEdges = TemplateType::nEdges();
+    const pointField& pointLst = TemplateType::points();
+    const vectorField& normLst = TemplateType::faceNormals();
+
+    // Check edge normals, face normals, point normals.
+    forAll (faceEs, faceI)
+    {
+        const labelList& edgeLabels = faceEs[faceI];
+
+        if (edgeLabels.size() < 3)
+        {
+            FatalErrorIn("PrimitivePatchExtra::checkOrientation(bool)")
+                << "face " << faceLst[faceI]
+                << " has fewer than 3 edges. Edges:" << edgeLabels
+                << exit(FatalError);
+        }
+
+        bool valid = true;
+        forAll (edgeLabels, i)
+        {
+            if (edgeLabels[i] < 0 || edgeLabels[i] >= numEdges)
+            {
+                WarningIn
+                (
+                    "PrimitivePatchExtra::checkOrientation(bool)"
+                )   << "edge number " << edgeLabels[i] << " on face " << faceI
+                    << " out of range"
+                    << "\nThis usually means that the input surface has "
+                    << "edges with more than 2 faces connected.\n"
+                    << endl;
+                valid = false;
+            }
+        }
+        if (!valid)
+        {
+            continue;
+        }
+
+
+        //
+        //- Compute normal from 3 points, use the first as the origin
+        //
+        const FaceType& f = faceLst[faceI];
+        const point p0(pointLst[f[0]]);
+        const point p1(pointLst[f[1]]);
+        const point p2(pointLst[f[f.size()-1]]);
+
+        const vector pointNormal((p1 - p0) ^ (p2 - p0));
+        if ((pointNormal & normLst[faceI]) < 0)
+        {
+            FatalErrorIn("PrimitivePatchExtra::checkOrientation(bool)")
+                << "Normal calculated from points not consistent with"
+                " faceNormal" << endl
+                << "face: " << f << endl
+                << "points: " << p0 << ' ' << p1 << ' ' << p2 << endl
+                << "pointNormal:" << pointNormal << endl
+                << "faceNormal:" << normLst[faceI]
+                << exit(FatalError);
+        }
+    }
+
+
+    const labelListList& eFaces = TemplateType::edgeFaces();
+    const pointField& locPointsLst = TemplateType::localPoints();
+
+    // Storage for holding status of edge. True if normal flips across this
+    // edge
+    boolList borderEdge(numEdges, false);
+
+    forAll (edgeLst, edgeI)
+    {
+        const labelList& neighbours = eFaces[edgeI];
+
+        if (neighbours.size() == 2)
+        {
+            // Two faces, A and B. Check if edge orientation is
+            // anticlockwise on both.
+            const labelList& fEdgesA = faceEs[neighbours[0]];
+            const labelList& fEdgesB = faceEs[neighbours[1]];
+
+            // Get next edge after edgeI
+            label nextEdgeA = fEdgesA.fcIndex(findIndex(fEdgesA, edgeI));
+            label nextEdgeB = fEdgesB.fcIndex(findIndex(fEdgesB, edgeI));
+
+            // Now check if nextEdgeA and nextEdgeB have any common points
+            if
+            (
+                edgeLst[nextEdgeA].start() == edgeLst[nextEdgeB].start()
+             || edgeLst[nextEdgeA].start() == edgeLst[nextEdgeB].end()
+             || edgeLst[nextEdgeA].end() == edgeLst[nextEdgeB].start()
+             || edgeLst[nextEdgeA].end() == edgeLst[nextEdgeB].end()
+            )
+            {
+                borderEdge[edgeI] = true;
+                if (verbose)
+                {
+                    // just list first three points
+                    // to simplify generating the message
+                    WarningIn("PrimitivePatchExtra::checkOrientation(bool)")
+                        << "face orientation incorrect." << nl
+                        << "edge neighbours:" << neighbours << nl
+                        << "face " << neighbours[0] << " has edges "
+                        << fEdgesA << nl
+                        << "    with points " << nl
+                        << "    " << edgeLst[fEdgesA[0]].start() << ' '
+                        << edgeLst[fEdgesA[0]].end() << nl
+                        << "    " << edgeLst[fEdgesA[1]].start() << ' '
+                        << edgeLst[fEdgesA[1]].end() << nl
+                        << "    " << edgeLst[fEdgesA[2]].start() << ' '
+                        << edgeLst[fEdgesA[2]].end()
+                        << endl
+                        << "face " << neighbours[1] << " has edges "
+                        << fEdgesB << nl
+                        << "    with points " << nl
+                        << "    " << edgeLst[fEdgesB[0]].start() << ' '
+                        << edgeLst[fEdgesB[0]].end() << nl
+                        << "    " << edgeLst[fEdgesB[1]].start() << ' '
+                        << edgeLst[fEdgesB[1]].end() << nl
+                        << "    " << edgeLst[fEdgesB[2]].start() << ' '
+                        << edgeLst[fEdgesB[2]].end() << nl
+                        << endl;
+                }
+            }
+        }
+        else if (neighbours.size() != 1)
+        {
+            if (verbose)
+            {
+                const edge& e = edgeLst[edgeI];
+                WarningIn("PrimitivePatchExtra::checkOrientation(bool)")
+                    << "Wrong number of edge neighbours." << endl
+                    << "Edge:" << e
+                    << "with points:" << locPointsLst[e.start()]
+                    << ' ' << locPointsLst[e.end()]
+                    << " has neighbours:" << neighbours << endl;
+            }
+            borderEdge[edgeI] = true;
+        }
+    }
+
+    return borderEdge;
+}
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraSearch.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraSearch.C
new file mode 100644
index 0000000000..c9d4d75fa6
--- /dev/null
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatchExtra/PrimitivePatchExtraSearch.C
@@ -0,0 +1,235 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "PrimitivePatchExtra.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+// Finds area, starting at faceI, delimited by borderEdge. Marks all visited
+// faces (from face-edge-face walk) with currentZone.
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+void Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::markZone
+(
+    const boolList& borderEdge,
+    const label faceI,
+    const label currentZone,
+    labelList& faceZone
+) const
+{
+    // List of faces whose faceZone has been set.
+    labelList changedFaces(1, faceI);
+
+    const labelListList& faceEs = TemplateType::faceEdges();
+    const labelListList& eFaces = TemplateType::edgeFaces();
+
+    while (true)
+    {
+        // Pick up neighbours of changedFaces
+        DynamicList<label> newChangedFaces(2*changedFaces.size());
+
+        forAll (changedFaces, i)
+        {
+            label faceI = changedFaces[i];
+
+            const labelList& fEdges = faceEs[faceI];
+
+            forAll(fEdges, i)
+            {
+                label edgeI = fEdges[i];
+
+                if (!borderEdge[edgeI])
+                {
+                    const labelList& eFaceLst = eFaces[edgeI];
+
+                    forAll(eFaceLst, j)
+                    {
+                        label nbrFaceI = eFaceLst[j];
+
+                        if (faceZone[nbrFaceI] == -1)
+                        {
+                            faceZone[nbrFaceI] = currentZone;
+                            newChangedFaces.append(nbrFaceI);
+                        }
+                        else if (faceZone[nbrFaceI] != currentZone)
+                        {
+                            FatalErrorIn
+                            (
+                                "PrimitivePatchExtra<Face, ListType, PointField>::markZone"
+                                "(const boolList&, const label, const label, labelList&) const"
+                            )
+                                << "Zones " << faceZone[nbrFaceI]
+                                << " at face " << nbrFaceI
+                                << " connects to zone " << currentZone
+                                << " at face " << faceI
+                                << abort(FatalError);
+                        }
+                    }
+                }
+            }
+        }
+
+        if (newChangedFaces.size() == 0)
+        {
+            break;
+        }
+
+        changedFaces.transfer(newChangedFaces.shrink());
+        newChangedFaces.clear();
+    }
+}
+
+
+// Finds areas delimited by borderEdge (or 'real' edges).
+// Fills faceZone accordingly
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+Foam::label Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+markZones
+(
+    const boolList& borderEdge,
+    labelList& faceZone
+) const
+{
+    const label numEdges = TemplateType::nEdges();
+    const label numFaces = TemplateType::size();
+
+    faceZone.setSize(numFaces);
+    faceZone = -1;
+
+    if (borderEdge.size() != numEdges)
+    {
+        FatalErrorIn
+        (
+            "PrimitivePatchExtra<Face, ListType, PointField>::markZones"
+            "(const boolList&, labelList&)"
+        )
+            << "borderEdge boolList not same size as number of edges" << endl
+            << "borderEdge:" << borderEdge.size() << endl
+            << "nEdges    :" << numEdges
+            << exit(FatalError);
+    }
+
+    label zoneI = 0;
+    label startFaceI = 0;
+
+    for (;;zoneI++)
+    {
+        // Find first non-coloured face
+        for (; startFaceI < numFaces; startFaceI++)
+        {
+            if (faceZone[startFaceI] == -1)
+            {
+                break;
+            }
+        }
+
+        if (startFaceI >= numFaces)
+        {
+            break;
+        }
+
+        faceZone[startFaceI] = zoneI;
+
+        markZone(borderEdge, startFaceI, zoneI, faceZone);
+    }
+
+    return zoneI;
+}
+
+
+
+// Finds areas delimited by borderEdge (or 'real' edges).
+// Fills faceZone accordingly
+template
+<
+    class Face,
+    template<class> class ListType,
+    class PointField,
+    class PointType
+>
+void Foam::PrimitivePatchExtra<Face, ListType, PointField, PointType>::
+subsetMap
+(
+    const boolList& include,
+    labelList& pointMap,
+    labelList& faceMap
+) const
+{
+    const List<FaceType>& locFaces = TemplateType::localFaces();
+    const label numPoints = TemplateType::nPoints();
+
+    label faceI = 0;
+    label pointI = 0;
+
+    faceMap.setSize(locFaces.size());
+    pointMap.setSize(numPoints);
+
+    boolList pointHad(numPoints, false);
+
+    forAll (include, oldFaceI)
+    {
+        if (include[oldFaceI])
+        {
+            // Store new faces compact
+            faceMap[faceI++] = oldFaceI;
+
+            // Renumber labels for face
+            const FaceType& f = locFaces[oldFaceI];
+
+            forAll (f, fp)
+            {
+                const label ptLabel = f[fp];
+                if (!pointHad[ptLabel])
+                {
+                    pointHad[ptLabel] = true;
+                    pointMap[pointI++] = ptLabel;
+                }
+            }
+        }
+    }
+
+    // Trim
+    faceMap.setSize(faceI);
+    pointMap.setSize(pointI);
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// ************************************************************************* //