diff --git a/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C b/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
index 4fda7ff5a9..f4ff9245e8 100644
--- a/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
+++ b/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
@@ -24,7 +24,6 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "backgroundMeshDecomposition.H"
-#include "zeroGradientFvPatchFields.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -40,21 +39,6 @@ defineTypeNameAndDebug(backgroundMeshDecomposition, 0);
 
 void Foam::backgroundMeshDecomposition::initialRefinement()
 {
-    //Read decomposePar dictionary
-    IOdictionary decomposeDict
-    (
-        IOobject
-        (
-            "decomposeParDict",
-            cvMesh_.time().system(),
-            cvMesh_.time(),
-            IOobject::MUST_READ_IF_MODIFIED,
-            IOobject::NO_WRITE
-        )
-    );
-
-    scalar mergeDist = 1e-6*mesh_.bounds().mag();
-
     volScalarField cellWeights
     (
         IOobject
@@ -72,25 +56,7 @@ void Foam::backgroundMeshDecomposition::initialRefinement()
 
     const conformationSurfaces& geometry = cvMesh_.geometryToConformTo();
 
-    // Decomposition
-    autoPtr<decompositionMethod> decomposerPtr
-    (
-        decompositionMethod::New(decomposeDict)
-    );
-
-    decompositionMethod& decomposer = decomposerPtr();
-
-    if (!decomposer.parallelAware())
-    {
-        FatalErrorIn
-        (
-            "void Foam::backgroundMeshDecomposition::initialRefinement() const"
-        )
-            << "You have selected decomposition method "
-            << decomposer.typeName
-            << " which is not parallel aware." << endl
-            << exit(FatalError);
-    }
+    decompositionMethod& decomposer = decomposerPtr_();
 
     hexRef8 meshCutter
     (
@@ -369,7 +335,7 @@ void Foam::backgroundMeshDecomposition::initialRefinement()
                 cellWeights
             );
 
-            fvMeshDistribute distributor(mesh_, mergeDist);
+            fvMeshDistribute distributor(mesh_, mergeDist_);
 
             autoPtr<mapDistributePolyMesh> mapDist =
                 distributor.distribute(newDecomp);
@@ -625,6 +591,88 @@ Foam::labelList Foam::backgroundMeshDecomposition::selectRefinementCells
 }
 
 
+void Foam::backgroundMeshDecomposition::buildPatchAndTree()
+{
+    primitivePatch tmpBoundaryFaces
+    (
+        SubList<face>
+        (
+            mesh_.faces(),
+            mesh_.nFaces() - mesh_.nInternalFaces(),
+            mesh_.nInternalFaces()
+        ),
+        mesh_.points()
+    );
+
+    boundaryFacesPtr_.reset
+    (
+        new bPatch
+        (
+            tmpBoundaryFaces.localFaces(),
+            tmpBoundaryFaces.localPoints()
+        )
+    );
+
+    // Overall bb
+    treeBoundBox overallBb(boundaryFacesPtr_().localPoints());
+
+    Random& rnd = cvMesh_.rndGen();
+
+    bFTreePtr_.reset
+    (
+        new indexedOctree<treeDataBPatch>
+        (
+            treeDataBPatch(false, boundaryFacesPtr_()),
+            overallBb.extend(rnd, 1e-4),
+            10, // maxLevel
+            10, // leafSize
+            3.0 // duplicity
+        )
+    );
+
+    if (debug)
+    {
+        OFstream fStr
+        (
+            cvMesh_.time().path()
+           /"backgroundMeshDecomposition_proc_"
+          + name(Pstream::myProcNo())
+          + "_boundaryFaces.obj"
+        );
+
+        const faceList& faces = boundaryFacesPtr_().localFaces();
+        const pointField& points = boundaryFacesPtr_().localPoints();
+
+        Map<label> foamToObj(points.size());
+
+        label vertI = 0;
+
+        forAll(faces, i)
+        {
+            const face& f = faces[i];
+
+            forAll(f, fPI)
+            {
+                if (foamToObj.insert(f[fPI], vertI))
+                {
+                    meshTools::writeOBJ(fStr, points[f[fPI]]);
+                    vertI++;
+                }
+            }
+
+            fStr<< 'f';
+
+            forAll(f, fPI)
+            {
+                fStr<< ' ' << foamToObj[f[fPI]] + 1;
+            }
+
+            fStr<< nl;
+        }
+    }
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
@@ -645,6 +693,21 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
             IOobject::MUST_READ
         )
     ),
+    boundaryFacesPtr_(),
+    bFTreePtr_(),
+    decomposeDict_
+    (
+        IOobject
+        (
+            "decomposeParDict",
+            cvMesh_.time().system(),
+            cvMesh_.time(),
+            IOobject::MUST_READ_IF_MODIFIED,
+            IOobject::NO_WRITE
+        )
+    ),
+    decomposerPtr_(decompositionMethod::New(decomposeDict_)),
+    mergeDist_(1e-6*mesh_.bounds().mag()),
     spanScale_(readScalar(coeffsDict_.lookup("spanScale"))),
     minCellSizeLimit_
     (
@@ -667,9 +730,23 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
             << exit(FatalError);
     }
 
+    if (!decomposerPtr_().parallelAware())
+    {
+        FatalErrorIn
+        (
+            "void Foam::backgroundMeshDecomposition::initialRefinement() const"
+        )
+            << "You have selected decomposition method "
+            << decomposerPtr_().typeName
+            << " which is not parallel aware." << endl
+            << exit(FatalError);
+    }
+
     Info<< nl << "Building initial background mesh decomposition" << endl;
 
     initialRefinement();
+
+    buildPatchAndTree();
 }
 
 
@@ -681,17 +758,54 @@ Foam::backgroundMeshDecomposition::~backgroundMeshDecomposition()
 
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
-bool Foam::conformationSurfaces::positionOnThisProc(const point& pt) const
+Foam::autoPtr<Foam::mapDistributePolyMesh>
+Foam::backgroundMeshDecomposition::distribute
+(
+    volScalarField& cellWeights
+)
 {
+    if (debug)
+    {
+        const_cast<Time&>(mesh_.time())++;
+        cellWeights.write();
+        mesh_.write();
+    }
 
+    labelList newDecomp = decomposerPtr_().decompose
+    (
+        mesh_,
+        mesh_.cellCentres(),
+        cellWeights
+    );
+
+    fvMeshDistribute distributor(mesh_, mergeDist_);
+
+    autoPtr<mapDistributePolyMesh> mapDist =
+        distributor.distribute(newDecomp);
+
+    if (debug)
+    {
+        printMeshData(mesh_);
+
+        const_cast<Time&>(mesh_.time())++;
+        mesh_.write();
+        cellWeights.write();
+    }
+
+    buildPatchAndTree();
+
+    return mapDist;
 }
 
 
-Foam::label Foam::conformationSurfaces::positionProc(const point& pt) const
+bool Foam::backgroundMeshDecomposition::positionInThisDomain
+(
+    const point& pt
+) const
 {
-
-
-    return -1;
+    return
+        bFTreePtr_().getVolumeType(pt)
+     == indexedOctree<treeDataBPatch>::INSIDE;
 }
 
 
diff --git a/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H b/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
index 879dcbfeb5..c6d77b89bd 100644
--- a/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
+++ b/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
@@ -66,13 +66,23 @@ SourceFiles
 #include "mapDistributePolyMesh.H"
 #include "globalIndex.H"
 #include "treeBoundBox.H"
-
+#include "primitivePatch.H"
+#include "face.H"
+#include "List.H"
+#include "pointField.H"
+#include "indexedOctree.H"
+#include "treeDataPrimitivePatch.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {
 
+typedef treeDataPrimitivePatch<face, List, const pointField, point>
+    treeDataBPatch;
+
+typedef PrimitivePatch<face, List, const pointField, point> bPatch;
+
 /*---------------------------------------------------------------------------*\
                   Class backgroundMeshDecomposition Declaration
 \*---------------------------------------------------------------------------*/
@@ -91,6 +101,22 @@ class backgroundMeshDecomposition
         //  is responsible for.
         fvMesh mesh_;
 
+        //- Patch containing an independent representation of the surface of the
+        //  mesh of this processor
+        autoPtr<bPatch> boundaryFacesPtr_;
+
+        //- Search tree for the boundaryFaces_ patch
+        autoPtr<indexedOctree<treeDataBPatch> > bFTreePtr_;
+
+        //- Decomposition dictionary
+        IOdictionary decomposeDict_;
+
+        //- Decomposition method
+        autoPtr<decompositionMethod> decomposerPtr_;
+
+        //- merge distance required by fvMeshDistribute
+        scalar mergeDist_;
+
         //- Scale of a cell span vs cell size used to decide to refine a cell
         scalar spanScale_;
 
@@ -133,6 +159,9 @@ class backgroundMeshDecomposition
             volScalarField& cellWeights
         ) const;
 
+        //- Build the surface patch and search tree
+        void buildPatchAndTree();
+
         //- Disallow default bitwise copy construct
         backgroundMeshDecomposition(const backgroundMeshDecomposition&);
 
@@ -162,23 +191,20 @@ public:
 
     // Member Functions
 
-    void distribute() const;
-
-    //- Check if the point is in the domain handled by this processor
-    bool positionOnThisProc(const point& pt) const;
-
-    //- Which processor's domain handles this point
-    label positionProc(const point& pt) const;
+        //- Redistribute the background mesh based on a supplied weight field,
+        //  returning a map to use to redistribute vertices.
+        autoPtr<mapDistributePolyMesh> distribute
+        (
+            volScalarField& cellWeights
+        );
 
+        //- Is the given position inside the domain of this decomposition
+        bool positionInThisDomain(const point& pt) const;
 
         // Access
 
-        // Check
-
-        // Edit
-
-        // Write
-
+            //- Return access to the underlying mesh
+            inline const fvMesh& mesh() const;
 };
 
 
diff --git a/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H b/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
index 718179dd65..a857ede285 100644
--- a/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
+++ b/src/mesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
@@ -23,36 +23,11 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * Destructors * * * * * * * * * * * * * * * //
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Ostream Operator  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
+const Foam::fvMesh& Foam::backgroundMeshDecomposition::mesh() const
+{
+    return mesh_;
+}
 
 // ************************************************************************* //
diff --git a/src/mesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeControlSurfaces.C b/src/mesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeControlSurfaces.C
index 54090d8769..84e472c5db 100644
--- a/src/mesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeControlSurfaces.C
+++ b/src/mesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeControlSurfaces.C
@@ -286,7 +286,7 @@ Foam::scalar Foam::cellSizeControlSurfaces::cellSize
 {
     if (isSurfacePoint)
     {
-        Pout<< "WARNING: isSurfacePoint is broken!" << endl;
+        // Pout<< "WARNING: isSurfacePoint is broken!" << endl;
     }
 
     scalar size = defaultCellSize_;
@@ -326,17 +326,17 @@ Foam::scalar Foam::cellSizeControlSurfaces::cellSize
             }
             else
             {
-                FatalErrorIn
-                (
-                    "Foam::scalar Foam::cellSizeControlSurfaces::cellSize"
-                    "("
-                        "const point& pt, "
-                        "bool isSurfacePoint"
-                    ") const"
-                )
-                    << "Point " << pt
-                    << " Cannot use isSurfacePoint here, or at all!"
-                    << nl << exit(FatalError) << endl;
+                // FatalErrorIn
+                // (
+                //     "Foam::scalar Foam::cellSizeControlSurfaces::cellSize"
+                //     "("
+                //         "const point& pt, "
+                //         "bool isSurfacePoint"
+                //     ") const"
+                // )
+                //     << "Point " << pt
+                //     << " Cannot use isSurfacePoint here, or at all!"
+                //     << nl << exit(FatalError) << endl;
 
                 // Evaluating the cell size at the nearest surface
                 evalCellSizeFunctions(surfHit.hitPoint(), size, true);
diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
index ce3cbd9f12..08f2f0d596 100644
--- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
+++ b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
@@ -28,6 +28,7 @@ License
 #include "relaxationModel.H"
 #include "faceAreaWeightModel.H"
 #include "backgroundMeshDecomposition.H"
+#include "meshSearch.H"
 
 // * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //
 
@@ -576,7 +577,7 @@ void Foam::conformalVoronoiMesh::insertConvexFeaturePoints()
         {
             const Foam::point& featPt = feMesh.points()[ptI];
 
-            if (!geometryToConformTo_.positionOnThisProc(featPt))
+            if (!positionOnThisProc(featPt))
             {
                 continue;
             }
@@ -628,7 +629,7 @@ void Foam::conformalVoronoiMesh::insertConcaveFeaturePoints()
         {
             const Foam::point& featPt = feMesh.points()[ptI];
 
-            if (!geometryToConformTo_.positionOnThisProc(featPt))
+            if (!positionOnThisProc(featPt))
             {
                 continue;
             }
@@ -728,7 +729,7 @@ void Foam::conformalVoronoiMesh::insertMixedFeaturePoints()
 
                 const Foam::point& pt(feMesh.points()[ptI]);
 
-                if (!geometryToConformTo_.positionOnThisProc(pt))
+                if (!positionOnThisProc(pt))
                 {
                     continue;
                 }
@@ -856,6 +857,116 @@ void Foam::conformalVoronoiMesh::insertInitialPoints()
 }
 
 
+void Foam::conformalVoronoiMesh::distribute()
+{
+    Info<< nl << "Redistributing points" << endl;
+
+    timeCheck("Before distribute");
+
+    const fvMesh& bMesh = decomposition_().mesh();
+
+    volScalarField cellWeights
+    (
+        IOobject
+        (
+            "cellWeights",
+            bMesh.time().timeName(),
+            bMesh,
+            IOobject::NO_READ,
+            IOobject::NO_WRITE
+        ),
+        bMesh,
+        dimensionedScalar("weight", dimless, 0),
+        zeroGradientFvPatchScalarField::typeName
+    );
+
+    meshSearch cellSearch(bMesh);
+
+    List<DynamicList<Foam::point> > cellVertices(bMesh.nCells());
+
+    for
+    (
+        Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
+        vit != finite_vertices_end();
+        vit++
+    )
+    {
+        if (vit->real())
+        {
+            Foam::point v = topoint(vit->point());
+
+            label cellI = cellSearch.findCell(v);
+
+            if (cellI == -1)
+            {
+                Pout<< "findCell "
+                    << vit->type() << " "
+                    << vit->index() << " "
+                    << v << " "
+                    << cellI << endl;
+            }
+            else
+            {
+                cellVertices[cellI].append(topoint(vit->point()));
+            }
+        }
+    }
+
+    forAll(cellVertices, cI)
+    {
+        cellWeights.internalField()[cI] = min
+        (
+            max(cellVertices[cI].size(), 1e-3),
+            100
+        );
+    }
+
+    autoPtr<mapDistributePolyMesh> mapDist = decomposition_().distribute
+    (
+        cellWeights
+    );
+
+    Pout<< "number_of_vertices() before distribution "
+        << label(number_of_vertices()) << endl;
+
+    mapDist().distributeCellData(cellVertices);
+
+    // Remove the entire tessellation
+    this->clear();
+
+    // Info<< "NEED TO MAP FEATURE POINTS" << endl;
+    // reinsertFeaturePoints();
+
+    startOfInternalPoints_ = number_of_vertices();
+
+    timeCheck("Distribution performed");
+
+    Info<< nl << "Inserting distributed tessellation" << endl;
+
+    std::vector<Point> pointsToInsert;
+
+    forAll(cellVertices, cI)
+    {
+        forAll(cellVertices[cI], cVPI)
+        {
+            pointsToInsert.push_back(toPoint(cellVertices[cI][cVPI]));
+        }
+    }
+
+    insertPoints(pointsToInsert);
+
+    Pout<< "number_of_vertices() after distribution "
+        << label(number_of_vertices()) << endl;
+
+    if(cvMeshControls().objOutput())
+    {
+        writePoints("distributedPoints.obj", true);
+    }
+
+    timeCheck("After distribute");
+}
+
+
 void Foam::conformalVoronoiMesh::storeSizesAndAlignments
 (
     const std::vector<Point>& storePts
@@ -1266,7 +1377,7 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
         );
     }
 
-    createFeaturePoints();
+    // createFeaturePoints();
 
     if (cvMeshControls().objOutput())
     {
@@ -1275,6 +1386,10 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
 
     insertInitialPoints();
 
+    distribute();
+
+    createFeaturePoints();
+
     buildSurfaceConformation(rmCoarse);
 
     if(cvMeshControls().objOutput())
@@ -1675,4 +1790,36 @@ void Foam::conformalVoronoiMesh::move()
 }
 
 
+bool Foam::conformalVoronoiMesh::positionOnThisProc
+(
+    const Foam::point& pt
+) const
+{
+    if (Pstream::parRun())
+    {
+        return decomposition_().positionInThisDomain(pt);
+    }
+
+    return true;
+}
+
+
+Foam::label Foam::conformalVoronoiMesh::positionProc
+(
+    const Foam::point& pt
+) const
+{
+    Pout<< "Dummy positionProc" << endl;
+
+    if (!Pstream::parRun())
+    {
+        return -1;
+    }
+    else
+    {
+        return Pstream::myProcNo();
+    }
+}
+
+
 // ************************************************************************* //
diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
index 4ade1c61cd..4a18ddfc23 100644
--- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
+++ b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
@@ -67,6 +67,7 @@ SourceFiles
 #include "cellSet.H"
 #include "wallPolyPatch.H"
 #include "processorPolyPatch.H"
+#include "zeroGradientFvPatchFields.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -82,7 +83,6 @@ class faceAreaWeightModel;
 
 class backgroundMeshDecomposition;
 
-
 /*---------------------------------------------------------------------------*\
                     Class conformalVoronoiMesh Declaration
 \*---------------------------------------------------------------------------*/
@@ -372,6 +372,10 @@ private:
         //  initialPointsMethod
         void insertInitialPoints();
 
+        //- In parallel redistribute the backgroundMeshDecomposition and
+        //- vertices to balance the number of vertices on each processor
+        void distribute();
+
         //- Store data for sizeAndAlignmentLocations_, storedSizes_ and
         //  storedAlignments_ and initialise the sizeAndAlignmentTreePtr_
         void storeSizesAndAlignments(const std::vector<Point>& storePts);
@@ -786,6 +790,12 @@ public:
         //  surface as required
         void move();
 
+        //- Check if the point is in the domain handled by this processor
+        bool positionOnThisProc(const Foam::point& pt) const;
+
+        //- Which processor's domain handles this point
+        label positionProc(const Foam::point& pt) const;
+
         //- Conversion functions between point (OpenFOAM) and Point (CGAL)
 
 #       ifdef CGAL_INEXACT
diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
index 6c83bff8b4..d0772251d7 100644
--- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
+++ b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
@@ -446,43 +446,43 @@ bool Foam::conformalVoronoiMesh::dualCellSurfaceAnyIntersection
         Foam::point dE0 = topoint(dual(fit->first));
         Foam::point dE1 = topoint(dual(fit->first->neighbor(fit->second)));
 
-        if (Pstream::parRun())
-        {
-            const treeBoundBoxList& procBbs =
-                geometryToConformTo_.processorDomains()[Pstream::myProcNo()];
+        // if (Pstream::parRun())
+        // {
+        //     const treeBoundBoxList& procBbs =
+        //         geometryToConformTo_.processorDomains()[Pstream::myProcNo()];
 
-            forAll(procBbs, pBI)
-            {
-                const treeBoundBox& procBb = procBbs[pBI];
+        //     forAll(procBbs, pBI)
+        //     {
+        //         const treeBoundBox& procBb = procBbs[pBI];
 
-                Foam::point a = dE0;
-                Foam::point b = dE1;
+        //         Foam::point a = dE0;
+        //         Foam::point b = dE1;
 
-                bool inBox = clipLineToBox(a, b, procBb);
+        //         bool inBox = clipLineToBox(a, b, procBb);
 
-                // Check for the edge passing through a surface
-                if
-                (
-                    inBox
-                 && geometryToConformTo_.findSurfaceAnyIntersection(a, b)
-                )
-                {
-                    // Pout<< "# findSurfaceAnyIntersection" << endl;
-                    // meshTools::writeOBJ(Pout, a);
-                    // meshTools::writeOBJ(Pout, b);
-                    // Pout<< "l cr0 cr1" << endl;
+        //         // Check for the edge passing through a surface
+        //         if
+        //         (
+        //             inBox
+        //          && geometryToConformTo_.findSurfaceAnyIntersection(a, b)
+        //         )
+        //         {
+        //             // Pout<< "# findSurfaceAnyIntersection" << endl;
+        //             // meshTools::writeOBJ(Pout, a);
+        //             // meshTools::writeOBJ(Pout, b);
+        //             // Pout<< "l cr0 cr1" << endl;
 
-                    return true;
-                }
-            }
-        }
-        else
-        {
+        //             return true;
+        //         }
+        //     }
+        // }
+        // else
+        // {
             if (geometryToConformTo_.findSurfaceAnyIntersection(dE0, dE1))
             {
                 return true;
             }
-        }
+        // }
     }
 
     return false;
@@ -635,7 +635,7 @@ void Foam::conformalVoronoiMesh::buildParallelInterface
 
     // Hard coded switch, can be turned on for debugging purposes and all
     // vertices will be referred to all processors.
-    bool allPointReferral = false;
+    bool allPointReferral = true;
 
     if (allPointReferral)
     {
@@ -1121,48 +1121,48 @@ void Foam::conformalVoronoiMesh::parallelInterfaceIntersection
         label closestHitProc = -1;
         scalar closestHitDistSqr = GREAT;
 
-        forAll(geometryToConformTo_.processorDomains(), procI)
-        {
-            if (procI == Pstream::myProcNo())
-            {
-                continue;
-            }
+        // forAll(geometryToConformTo_.processorDomains(), procI)
+        // {
+        //     if (procI == Pstream::myProcNo())
+        //     {
+        //         continue;
+        //     }
 
-            const treeBoundBoxList& procBbs =
-                geometryToConformTo_.processorDomains()[procI];
+        //     const treeBoundBoxList& procBbs =
+        //         geometryToConformTo_.processorDomains()[procI];
 
-            forAll(procBbs, pBI)
-            {
-                const treeBoundBox& procBb = procBbs[pBI];
+        //     forAll(procBbs, pBI)
+        //     {
+        //         const treeBoundBox& procBb = procBbs[pBI];
 
-                intersects = procBb.intersects(dE0, dE1, boxPt);
+        //         intersects = procBb.intersects(dE0, dE1, boxPt);
 
-                if (intersects)
-                {
-                    hitDistSqr = magSqr(vert - boxPt);
+        //         if (intersects)
+        //         {
+        //             hitDistSqr = magSqr(vert - boxPt);
 
-                    if (hitDistSqr < closestHitDistSqr)
-                    {
-                        closestHitProc = procI;
-                        closestHitDistSqr = hitDistSqr;
-                    }
-                }
+        //             if (hitDistSqr < closestHitDistSqr)
+        //             {
+        //                 closestHitProc = procI;
+        //                 closestHitDistSqr = hitDistSqr;
+        //             }
+        //         }
 
-                // Perform the query in the opposite direction
-                intersects = procBb.intersects(dE1, dE0, boxPt);
+        //         // Perform the query in the opposite direction
+        //         intersects = procBb.intersects(dE1, dE0, boxPt);
 
-                if (intersects)
-                {
-                    hitDistSqr = magSqr(vert - boxPt);
+        //         if (intersects)
+        //         {
+        //             hitDistSqr = magSqr(vert - boxPt);
 
-                    if (hitDistSqr < closestHitDistSqr)
-                    {
-                        closestHitProc = procI;
-                        closestHitDistSqr = hitDistSqr;
-                    }
-                }
-            }
-        }
+        //             if (hitDistSqr < closestHitDistSqr)
+        //             {
+        //                 closestHitProc = procI;
+        //                 closestHitDistSqr = hitDistSqr;
+        //             }
+        //         }
+        //     }
+        // }
 
         if (closestHitProc >= 0)
         {
@@ -1204,21 +1204,21 @@ void Foam::conformalVoronoiMesh::parallelInterfaceInfluence
         // Pout<< nl << "# circumradius " << sqrt(circumradiusSqr) << endl;
         // drawDelaunayCell(Pout, cit);
 
-        forAll(geometryToConformTo_.processorMeshBounds(), procI)
-        {
-            if (procI == Pstream::myProcNo())
-            {
-                continue;
-            }
+        // forAll(geometryToConformTo_.processorMeshBounds(), procI)
+        // {
+        //     if (procI == Pstream::myProcNo())
+        //     {
+        //         continue;
+        //     }
 
-            const treeBoundBox& procBb =
-                geometryToConformTo_.processorMeshBounds()[procI];
+        //     const treeBoundBox& procBb =
+        //         geometryToConformTo_.processorMeshBounds()[procI];
 
-            if (procBb.overlaps(circumcentre, circumradiusSqr))
-            {
-                toProc[procI] = true;
-            }
-        }
+        //     if (procBb.overlaps(circumcentre, circumradiusSqr))
+        //     {
+        //         toProc[procI] = true;
+        //     }
+        // }
     }
 }
 
@@ -1303,7 +1303,7 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceProtrusion
     if
     (
         surfHitLargest.hit()
-     && !geometryToConformTo_.positionOnThisProc(surfHitLargest.hitPoint())
+     && !positionOnThisProc(surfHitLargest.hitPoint())
     )
     {
         // A protrusion was identified, but not penetrating on this processor,
@@ -1400,7 +1400,7 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceIncursion
     if
     (
         surfHitLargest.hit()
-     && !geometryToConformTo_.positionOnThisProc(surfHitLargest.hitPoint())
+     && !positionOnThisProc(surfHitLargest.hitPoint())
     )
     {
         // A protrusion was identified, but not penetrating on this processor,
@@ -1757,7 +1757,7 @@ void Foam::conformalVoronoiMesh::addSurfaceAndEdgeHits
 
         if (edHit.hit())
         {
-            if(!geometryToConformTo_.positionOnThisProc(edHit.hitPoint()))
+            if(!positionOnThisProc(edHit.hitPoint()))
             {
                 continue;
             }
diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePointSpecialisations.C b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePointSpecialisations.C
index 928b964c87..8a4bc64b5d 100644
--- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePointSpecialisations.C
+++ b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePointSpecialisations.C
@@ -94,7 +94,7 @@ bool Foam::conformalVoronoiMesh::insertSpecialisedFeaturePoint
 
     if (nExternal == 2 && nInternal == 1)
     {
-        // if (!geometryToConformTo_.positionOnThisProc(pt))
+        // if (!positionOnThisProc(pt))
         // {
         //     continue;
         // }
diff --git a/src/mesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C b/src/mesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C
index 93e1e5041d..4e9c94ae48 100644
--- a/src/mesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C
+++ b/src/mesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C
@@ -771,7 +771,7 @@ Foam::label Foam::conformationSurfaces::findPatch(const point& pt) const
     findSurfaceNearest
     (
         pt,
-        spanMagSqr_,
+        sqr(GREAT),
         surfHit,
         hitSurface
     );
diff --git a/src/mesh/conformalVoronoiMesh/initialPointsMethod/bodyCentredCubic/bodyCentredCubic.C b/src/mesh/conformalVoronoiMesh/initialPointsMethod/bodyCentredCubic/bodyCentredCubic.C
index 1f2f10a068..d32e06d158 100644
--- a/src/mesh/conformalVoronoiMesh/initialPointsMethod/bodyCentredCubic/bodyCentredCubic.C
+++ b/src/mesh/conformalVoronoiMesh/initialPointsMethod/bodyCentredCubic/bodyCentredCubic.C
@@ -58,7 +58,7 @@ bodyCentredCubic::bodyCentredCubic
 
 std::vector<Vb::Point> bodyCentredCubic::initialPoints() const
 {
-    const boundBox& bb = cvMesh_.geometryToConformTo().bounds();
+    const boundBox& bb = cvMesh_.geometryToConformTo().globalBounds();
 
     scalar x0 = bb.min().x();
     scalar xR = bb.max().x() - x0;
diff --git a/src/mesh/conformalVoronoiMesh/initialPointsMethod/densityWeightedStochastic/densityWeightedStochastic.C b/src/mesh/conformalVoronoiMesh/initialPointsMethod/densityWeightedStochastic/densityWeightedStochastic.C
index 977c182d3b..662516a7f2 100644
--- a/src/mesh/conformalVoronoiMesh/initialPointsMethod/densityWeightedStochastic/densityWeightedStochastic.C
+++ b/src/mesh/conformalVoronoiMesh/initialPointsMethod/densityWeightedStochastic/densityWeightedStochastic.C
@@ -66,7 +66,7 @@ densityWeightedStochastic::densityWeightedStochastic
 
 std::vector<Vb::Point> densityWeightedStochastic::initialPoints() const
 {
-    const boundBox& bb = cvMesh_.geometryToConformTo().bounds();
+    const boundBox& bb = cvMesh_.geometryToConformTo().globalBounds();
 
     Random& rndGen = cvMesh_.rndGen();
 
diff --git a/src/mesh/conformalVoronoiMesh/initialPointsMethod/faceCentredCubic/faceCentredCubic.C b/src/mesh/conformalVoronoiMesh/initialPointsMethod/faceCentredCubic/faceCentredCubic.C
index 8c4fbf5399..49d61dbc83 100644
--- a/src/mesh/conformalVoronoiMesh/initialPointsMethod/faceCentredCubic/faceCentredCubic.C
+++ b/src/mesh/conformalVoronoiMesh/initialPointsMethod/faceCentredCubic/faceCentredCubic.C
@@ -58,7 +58,7 @@ faceCentredCubic::faceCentredCubic
 
 std::vector<Vb::Point> faceCentredCubic::initialPoints() const
 {
-    const boundBox& bb = cvMesh_.geometryToConformTo().bounds();
+    const boundBox& bb = cvMesh_.geometryToConformTo().globalBounds();
 
     scalar x0 = bb.min().x();
     scalar xR = bb.max().x() - x0;
diff --git a/src/mesh/conformalVoronoiMesh/initialPointsMethod/hierarchicalDensityWeightedStochastic/hierarchicalDensityWeightedStochastic.C b/src/mesh/conformalVoronoiMesh/initialPointsMethod/hierarchicalDensityWeightedStochastic/hierarchicalDensityWeightedStochastic.C
index 92fa817b0f..d6e441b8ba 100644
--- a/src/mesh/conformalVoronoiMesh/initialPointsMethod/hierarchicalDensityWeightedStochastic/hierarchicalDensityWeightedStochastic.C
+++ b/src/mesh/conformalVoronoiMesh/initialPointsMethod/hierarchicalDensityWeightedStochastic/hierarchicalDensityWeightedStochastic.C
@@ -758,7 +758,10 @@ hierarchicalDensityWeightedStochastic::initialPoints() const
 {
     const conformationSurfaces& geometry = cvMesh_.geometryToConformTo();
 
-    treeBoundBox hierBB = geometry.bounds().extend(cvMesh_.rndGen(), 1e-6);
+    treeBoundBox hierBB = geometry.globalBounds().extend
+    (
+        cvMesh_.rndGen(), 1e-6
+    );
 
     std::vector<Vb::Point> initialPoints;
 
diff --git a/src/mesh/conformalVoronoiMesh/initialPointsMethod/pointFile/pointFile.C b/src/mesh/conformalVoronoiMesh/initialPointsMethod/pointFile/pointFile.C
index 75fb9e6b1a..9bc15eed18 100644
--- a/src/mesh/conformalVoronoiMesh/initialPointsMethod/pointFile/pointFile.C
+++ b/src/mesh/conformalVoronoiMesh/initialPointsMethod/pointFile/pointFile.C
@@ -79,10 +79,7 @@ std::vector<Vb::Point> pointFile::initialPoints() const
 
     forAll(points, ptI)
     {
-        procPt[ptI] = cvMesh_.geometryToConformTo().positionOnThisProc
-        (
-            points[ptI]
-        );
+        procPt[ptI] = cvMesh_.positionOnThisProc(points[ptI]);
     }
 
     if (Pstream::parRun())
@@ -95,32 +92,6 @@ std::vector<Vb::Point> pointFile::initialPoints() const
 
         Pstream::scatterList(allProcPt);
 
-        // forAll(procPt, ptI)
-        // {
-        //     label nProcsTrue = 0;
-
-        //     forAll(allProcPt, procI)
-        //     {
-        //         if (allProcPt[procI][ptI])
-        //         {
-        //             nProcsTrue++;
-        //         }
-        //     }
-
-        //     if (nProcsTrue != 1)
-        //     {
-        //         Info<< points[ptI] << " " << nProcsTrue << " true" << endl;
-
-        //         if (nProcsTrue > 1)
-        //         {
-        //             forAll(allProcPt, procI)
-        //             {
-        //                 allProcPt[procI][ptI] = false;
-        //             }
-        //         }
-        //     }
-        // }
-
         forAll(procPt, ptI)
         {
             bool foundAlready = false;
diff --git a/src/mesh/conformalVoronoiMesh/initialPointsMethod/uniformGrid/uniformGrid.C b/src/mesh/conformalVoronoiMesh/initialPointsMethod/uniformGrid/uniformGrid.C
index d5895dc955..877b58635f 100644
--- a/src/mesh/conformalVoronoiMesh/initialPointsMethod/uniformGrid/uniformGrid.C
+++ b/src/mesh/conformalVoronoiMesh/initialPointsMethod/uniformGrid/uniformGrid.C
@@ -58,7 +58,7 @@ uniformGrid::uniformGrid
 
 std::vector<Vb::Point> uniformGrid::initialPoints() const
 {
-    const boundBox& bb = cvMesh_.geometryToConformTo().bounds();
+    const boundBox& bb = cvMesh_.geometryToConformTo().globalBounds();
 
     scalar x0 = bb.min().x();
     scalar xR = bb.max().x() - x0;