diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/Make/files b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/Make/files
index c6e553f64c..7e9d038441 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/Make/files
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/Make/files
@@ -29,6 +29,10 @@ $(cellSizeAndAlignmentControl)/searchableSurfaceControl/searchableSurfaceControl
 
 cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C
 
+cellShapeControl/smoothAlignmentSolver/smoothAlignmentSolver.C
+
+cellShapeControl/controlMeshRefinement/controlMeshRefinement.C
+
 /*cellSizeControlSurfaces/cellSizeControlSurfaces.C*/
 
 cellSizeFunctions = cellSizeControlSurfaces/cellSizeFunction
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.C
index 932b264f94..444a0a4bf0 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.C
@@ -31,213 +31,16 @@ License
 #include "searchableSurfaceControl.H"
 #include "cellSizeFunction.H"
 
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+defineTypeNameAndDebug(cellShapeControl, 0);
+}
+
+
 // * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //
 
-template <class Triangulation, class Type>
-Foam::tmp<Foam::Field<Type> > Foam::cellShapeControl::filterFarPoints
-(
-    const Triangulation& mesh,
-    const Field<Type>& field
-)
-{
-    tmp<Field<Type> > tNewField(new Field<Type>(field.size()));
-    Field<Type>& newField = tNewField();
-
-    label added = 0;
-    label count = 0;
-
-    for
-    (
-        typename Triangulation::Finite_vertices_iterator vit =
-            mesh.finite_vertices_begin();
-        vit != mesh.finite_vertices_end();
-        ++vit
-    )
-    {
-        if (vit->real())
-        {
-            newField[added++] = field[count];
-        }
-
-        count++;
-    }
-
-    newField.resize(added);
-
-    return tNewField;
-}
-
-
-template <class Triangulation>
-Foam::autoPtr<Foam::mapDistribute> Foam::cellShapeControl::buildReferredMap
-(
-    const Triangulation& mesh,
-    labelList& indices
-)
-{
-    globalIndex globalIndexing(mesh.vertexCount());
-
-    DynamicList<label> dynIndices(mesh.vertexCount()/10);
-
-    for
-    (
-        typename Triangulation::Finite_vertices_iterator vit =
-            mesh.finite_vertices_begin();
-        vit != mesh.finite_vertices_end();
-        ++vit
-    )
-    {
-        if (vit->referred())
-        {
-            dynIndices.append
-            (
-                globalIndexing.toGlobal(vit->procIndex(), vit->index())
-            );
-        }
-    }
-
-    indices.transfer(dynIndices);
-
-    List<Map<label> > compactMap;
-    return autoPtr<mapDistribute>
-    (
-        new mapDistribute
-        (
-            globalIndexing,
-            indices,
-            compactMap
-        )
-    );
-}
-
-
-template <class Triangulation>
-Foam::autoPtr<Foam::mapDistribute> Foam::cellShapeControl::buildMap
-(
-    const Triangulation& mesh,
-    labelListList& pointPoints
-)
-{
-    pointPoints.setSize(mesh.vertexCount());
-
-    globalIndex globalIndexing(mesh.vertexCount());
-
-    for
-    (
-        typename Triangulation::Finite_vertices_iterator vit =
-            mesh.finite_vertices_begin();
-        vit != mesh.finite_vertices_end();
-        ++vit
-    )
-    {
-        if (!vit->real())
-        {
-            continue;
-        }
-
-        std::list<typename Triangulation::Vertex_handle> adjVerts;
-        mesh.finite_adjacent_vertices(vit, std::back_inserter(adjVerts));
-
-        DynamicList<label> indices(adjVerts.size());
-
-        for
-        (
-            typename std::list<typename Triangulation::Vertex_handle>::
-                const_iterator adjVertI = adjVerts.begin();
-            adjVertI != adjVerts.end();
-            ++adjVertI
-        )
-        {
-            typename Triangulation::Vertex_handle vh = *adjVertI;
-
-            if (!vh->farPoint())
-            {
-                indices.append
-                (
-                    globalIndexing.toGlobal(vh->procIndex(), vh->index())
-                );
-            }
-        }
-
-        pointPoints[vit->index()].transfer(indices);
-    }
-
-    List<Map<label> > compactMap;
-    return autoPtr<mapDistribute>
-    (
-        new mapDistribute
-        (
-            globalIndexing,
-            pointPoints,
-            compactMap
-        )
-    );
-}
-
-
-template <class Triangulation>
-Foam::tmp<Foam::triadField> Foam::cellShapeControl::buildAlignmentField
-(
-    const Triangulation& mesh
-)
-{
-    tmp<triadField> tAlignments
-    (
-        new triadField(mesh.vertexCount(), triad::unset)
-    );
-    triadField& alignments = tAlignments();
-
-    for
-    (
-        typename Triangulation::Finite_vertices_iterator vit =
-            mesh.finite_vertices_begin();
-        vit != mesh.finite_vertices_end();
-        ++vit
-    )
-    {
-        if (!vit->real())
-        {
-            continue;
-        }
-
-        alignments[vit->index()] = vit->alignment();
-    }
-
-    return tAlignments;
-}
-
-
-template <class Triangulation>
-Foam::tmp<Foam::pointField> Foam::cellShapeControl::buildPointField
-(
-    const Triangulation& mesh
-)
-{
-    tmp<pointField> tPoints
-    (
-        new pointField(mesh.vertexCount(), point(GREAT, GREAT, GREAT))
-    );
-    pointField& points = tPoints();
-
-    for
-    (
-        typename Triangulation::Finite_vertices_iterator vit =
-            mesh.finite_vertices_begin();
-        vit != mesh.finite_vertices_end();
-        ++vit
-    )
-    {
-        if (!vit->real())
-        {
-            continue;
-        }
-
-        points[vit->index()] = topoint(vit->point());
-    }
-
-    return tPoints;
-}
-
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
@@ -261,7 +64,7 @@ Foam::cellShapeControl::cellShapeControl
     (
         runTime,
         motionDict.subDict("shapeControlFunctions"),
-        geometryToConformTo,
+        geometryToConformTo_,
         defaultCellSize_
     )
 {}
@@ -300,16 +103,7 @@ Foam::scalar Foam::cellShapeControl::cellSize(const point& pt) const
 
     scalar size = 0;
 
-    label nFarPoints = 0;
-    for (label pI = 0; pI < 4; ++pI)
-    {
-        if (ch->vertex(pI)->farPoint())
-        {
-            nFarPoints++;
-        }
-    }
-
-    if (shapeControlMesh_.is_infinite(ch))
+    if (shapeControlMesh_.dimension() < 3 || shapeControlMesh_.is_infinite(ch))
     {
 //        if (nFarPoints != 0)
 //        {
@@ -323,22 +117,25 @@ Foam::scalar Foam::cellShapeControl::cellSize(const point& pt) const
 //            }
 //        }
 
-        // Look up nearest point
-        cellShapeControlMesh::Vertex_handle nearV =
-            shapeControlMesh_.nearest_vertex
-            (
-                toPoint<cellShapeControlMesh::Point>(pt)
-            );
-
-        size = nearV->targetCellSize();
+        // Find nearest surface
+        size = sizeAndAlignment_.cellSize(pt);
     }
     else
     {
+        label nFarPoints = 0;
+        for (label pI = 0; pI < 4; ++pI)
+        {
+            if (ch->vertex(pI)->farPoint())
+            {
+                nFarPoints++;
+            }
+        }
+
         if (nFarPoints != 0)
         {
             for (label pI = 0; pI < 4; ++pI)
             {
-                if (!ch->vertex(pI)->farPoint())
+                if (!ch->vertex(pI)->uninitialised())
                 {
                     size = ch->vertex(pI)->targetCellSize();
                     return size;
@@ -368,54 +165,58 @@ Foam::tensor Foam::cellShapeControl::cellAlignment(const point& pt) const
 
     tensor alignment = tensor::zero;
 
-    label nFarPoints = 0;
-    for (label pI = 0; pI < 4; ++pI)
+    if (shapeControlMesh_.dimension() < 3 || shapeControlMesh_.is_infinite(ch))
     {
-        if (ch->vertex(pI)->farPoint())
-        {
-            nFarPoints++;
-        }
-    }
-
-    if (shapeControlMesh_.is_infinite(ch) || nFarPoints == 4)
-    {
-        Pout<< "At Infinite vertex" << endl;
-
-        if (nFarPoints != 0)
-        {
-            for (label pI = 0; pI < 4; ++pI)
-            {
-                if (!ch->vertex(pI)->farPoint())
-                {
-                    alignment = ch->vertex(pI)->alignment();
-                    return alignment;
-                }
-            }
-        }
-
-//        cellShapeControlMesh::Vertex_handle nearV =
-//            shapeControlMesh_.nearest_vertex
-//            (
-//                toPoint<cellShapeControlMesh::Point>(pt)
-//            );
-//
-//        alignment = nearV->alignment();
+        alignment = tensor::I;
     }
     else
     {
-//        forAll(bary, pI)
+        label nFarPoints = 0;
+        for (label pI = 0; pI < 4; ++pI)
+        {
+            if (ch->vertex(pI)->farPoint())
+            {
+                nFarPoints++;
+            }
+        }
+
+//        if (nFarPoints != 0)
 //        {
-//            alignment += bary[pI]*ch->vertex(pI)->alignment();
+//            for (label pI = 0; pI < 4; ++pI)
+//            {
+//                if (!ch->vertex(pI)->farPoint())
+//                {
+//                    alignment = ch->vertex(pI)->alignment();
+//                }
+//            }
 //        }
+//        else
+        {
+            triad tri;
 
-        cellShapeControlMesh::Vertex_handle nearV =
-            shapeControlMesh_.nearest_vertex_in_cell
-            (
-                toPoint<cellShapeControlMesh::Point>(pt),
-                ch
-            );
+            for (label pI = 0; pI < 4; ++pI)
+            {
+                if (bary[pI] > SMALL)
+                {
+                    tri += triad(bary[pI]*ch->vertex(pI)->alignment());
+                }
+            }
 
-        alignment = nearV->alignment();
+            tri.normalize();
+            tri.orthogonalize();
+            //            tri = tri.sortxyz();
+
+            alignment = tri;
+        }
+
+//        cellShapeControlMesh::Vertex_handle nearV =
+//            shapeControlMesh_.nearest_vertex_in_cell
+//            (
+//                toPoint<cellShapeControlMesh::Point>(pt),
+//                ch
+//            );
+//
+//        alignment = nearV->alignment();
     }
 
     return alignment;
@@ -437,491 +238,81 @@ void Foam::cellShapeControl::cellSizeAndAlignment
     alignment = tensor::zero;
     size = 0;
 
-    label nFarPoints = 0;
-    for (label pI = 0; pI < 4; ++pI)
+    if (shapeControlMesh_.dimension() < 3 || shapeControlMesh_.is_infinite(ch))
     {
-        if (ch->vertex(pI)->farPoint())
-        {
-            nFarPoints++;
-        }
-    }
-
-    if (shapeControlMesh_.is_infinite(ch))
-    {
-        if (nFarPoints != 0)
-        {
-            for (label pI = 0; pI < 4; ++pI)
-            {
-                if (!ch->vertex(pI)->farPoint())
-                {
-                    size = ch->vertex(pI)->targetCellSize();
-                    alignment = ch->vertex(pI)->alignment();
-                    return;
-                }
-            }
-        }
-
-//        cellShapeControlMesh::Vertex_handle nearV =
-//            shapeControlMesh_.nearest_vertex
-//            (
-//                toPoint<cellShapeControlMesh::Point>(pt)
-//            );
-//
-//        size = nearV->targetCellSize();
-//        alignment = nearV->alignment();
+        // Find nearest surface
+        size = sizeAndAlignment_.cellSize(pt);
+        alignment = tensor::I;
     }
     else
     {
+        label nFarPoints = 0;
+        for (label pI = 0; pI < 4; ++pI)
+        {
+            if (ch->vertex(pI)->farPoint())
+            {
+                nFarPoints++;
+            }
+        }
+
         if (nFarPoints != 0)
         {
             for (label pI = 0; pI < 4; ++pI)
             {
-                if (!ch->vertex(pI)->farPoint())
+                if (!ch->vertex(pI)->uninitialised())
                 {
                     size = ch->vertex(pI)->targetCellSize();
                     alignment = ch->vertex(pI)->alignment();
-                    return;
                 }
             }
         }
         else
         {
-//            triad tri;
+            triad tri;
 
-            forAll(bary, pI)
+            for (label pI = 0; pI < 4; ++pI)
             {
                 size += bary[pI]*ch->vertex(pI)->targetCellSize();
 
-//                triad triTmp2 = ch->vertex(pI)->alignment();
-//                tri += triTmp2*bary[pI];
+                if (bary[pI] > SMALL)
+                {
+                    tri += triad(bary[pI]*ch->vertex(pI)->alignment());
+                }
             }
 
-//            tri.normalize();
-//            tri.orthogonalize();
+            tri.normalize();
+            tri.orthogonalize();
 //            tri = tri.sortxyz();
+
+            alignment = tri;
+
+//            cellShapeControlMesh::Vertex_handle nearV =
+//                shapeControlMesh_.nearest_vertex
+//                (
+//                    toPoint<cellShapeControlMesh::Point>(pt)
+//                );
 //
-//            alignment = tri;
-
-            cellShapeControlMesh::Vertex_handle nearV =
-                shapeControlMesh_.nearest_vertex
-                (
-                    toPoint<cellShapeControlMesh::Point>(pt)
-                );
-
-            alignment = nearV->alignment();
+//            alignment = nearV->alignment();
         }
     }
-}
 
-
-void Foam::cellShapeControl::initialMeshPopulation
-(
-    const autoPtr<backgroundMeshDecomposition>& decomposition
-)
-{
-    // Need to pass in the background mesh decomposition so that can test if
-    // a point to insert is on the processor.
-    if (Pstream::parRun())
+    for (label dir = 0; dir < 3; dir++)
     {
-        shapeControlMesh_.insertBoundingPoints(decomposition().procBounds());
-    }
-    else
-    {
-        shapeControlMesh_.insertBoundingPoints(allGeometry_.bounds());
-    }
+        const triad v = alignment;
 
-    const PtrList<cellSizeAndAlignmentControl>& controlFunctions =
-        sizeAndAlignment_.controlFunctions();
-
-    forAll(controlFunctions, fI)
-    {
-        const cellSizeAndAlignmentControl& controlFunction =
-            controlFunctions[fI];
-
-        Info<< "Inserting points from " << controlFunction.name()
-            << " (" << controlFunction.type() << ")" << endl;
-
-        pointField pts;
-        scalarField sizes;
-        triadField alignments;
-
-        controlFunction.initialVertices(pts, sizes, alignments);
-
-        List<Vb> vertices(pts.size());
-
-        // Clip the minimum size
-        forAll(vertices, vI)
+        if (!v.set(dir) || size == 0)
         {
-            vertices[vI] = Vb(pts[vI], Vb::vtInternal);
-            vertices[vI].targetCellSize() = max(sizes[vI], minimumCellSize_);
-            vertices[vI].alignment() = alignments[vI];
-        }
-
-        pts.clear();
-        sizes.clear();
-        alignments.clear();
-
-        label nRejected = 0;
-
-        PackedBoolList keepVertex(vertices.size(), true);
-
-        if (Pstream::parRun())
-        {
-            forAll(vertices, vI)
-            {
-                const bool onProc = decomposition().positionOnThisProcessor
-                (
-                    topoint(vertices[vI].point())
-                );
-
-                if (!onProc)
-                {
-                    keepVertex[vI] = false;
-                }
-            }
-        }
-
-        inplaceSubset(keepVertex, vertices);
-
-        const label preInsertedSize = shapeControlMesh_.number_of_vertices();
-
-        shapeControlMesh_.rangeInsertWithInfo(vertices.begin(), vertices.end());
-
-        Info<< "    Inserted "
-            << returnReduce
-               (
-                   shapeControlMesh_.number_of_vertices()
-                 - preInsertedSize, sumOp<label>()
-               )
-            << "/" << vertices.size()
-            << endl;
-    }
-}
-
-
-Foam::label Foam::cellShapeControl::refineMesh
-(
-    const autoPtr<backgroundMeshDecomposition>& decomposition
-)
-{
-    const pointField cellCentres(shapeControlMesh_.cellCentres());
-
-    Info<< "    Created cell centres" << endl;
-
-    const PtrList<cellSizeAndAlignmentControl>& controlFunctions =
-        sizeAndAlignment_.controlFunctions();
-
-    DynamicList<Vb> verts(shapeControlMesh_.number_of_vertices());
-
-    forAll(cellCentres, cellI)
-    {
-        const Foam::point& pt = cellCentres[cellI];
-
-        if (!geometryToConformTo_.inside(pt))
-        {
-            continue;
-        }
-
-        scalarList bary;
-        cellShapeControlMesh::Cell_handle ch;
-
-        shapeControlMesh_.barycentricCoords(pt, bary, ch);
-
-        if (shapeControlMesh_.is_infinite(ch))
-        {
-            continue;
-        }
-
-        scalar interpolatedSize = 0;
-        forAll(bary, pI)
-        {
-            interpolatedSize += bary[pI]*ch->vertex(pI)->targetCellSize();
-        }
-
-        label lastPriority = labelMax;
-        scalar lastCellSize = GREAT;
-        forAll(controlFunctions, fI)
-        {
-            const cellSizeAndAlignmentControl& controlFunction =
-                controlFunctions[fI];
-
-            if (controlFunction.priority() > lastPriority)
-            {
-                continue;
-            }
-
-            if (isA<searchableSurfaceControl>(controlFunction))
-            {
-                const cellSizeFunction& sizeFunction =
-                    dynamicCast<const searchableSurfaceControl>
-                    (
-                        controlFunction
-                    ).sizeFunction();
-
-                scalar cellSize = 0;
-                if (sizeFunction.cellSize(pt, cellSize))
-                {
-                    if (controlFunction.priority() == lastPriority)
-                    {
-                        if (cellSize < lastCellSize)
-                        {
-                            lastCellSize = cellSize;
-                        }
-                    }
-                    else
-                    {
-                        lastCellSize = cellSize;
-                    }
-
-                    lastPriority = controlFunction.priority();
-                }
-            }
-        }
-
-        if
-        (
-            lastCellSize < GREAT
-         //&& mag(interpolatedSize - lastCellSize)/lastCellSize > 0.2
-        )
-        {
-            if (Pstream::parRun())
-            {
-                if (!decomposition().positionOnThisProcessor(pt))
-                {
-                    continue;
-                }
-            }
-
-            verts.append
+            FatalErrorIn
             (
-                Vb
-                (
-                    toPoint<cellShapeControlMesh::Point>(pt),
-                    Vb::vtInternal
-                )
-            );
-            verts.last().targetCellSize() = lastCellSize;
-            verts.last().alignment() = triad::unset;
-        }
-    }
-
-    shapeControlMesh_.insertPoints(verts);
-
-    return verts.size();
-}
-
-
-void Foam::cellShapeControl::smoothMesh()
-{
-    label maxSmoothingIterations = readLabel(lookup("maxSmoothingIterations"));
-    scalar minResidual = 0;
-
-    labelListList pointPoints;
-    autoPtr<mapDistribute> meshDistributor = buildMap
-    (
-        shapeControlMesh_,
-        pointPoints
-    );
-
-    triadField alignments(buildAlignmentField(shapeControlMesh_));
-    pointField points(buildPointField(shapeControlMesh_));
-    // Setup the sizes and alignments on each point
-    triadField fixedAlignments(shapeControlMesh_.vertexCount(), triad::unset);
-
-    for
-    (
-        CellSizeDelaunay::Finite_vertices_iterator vit =
-            shapeControlMesh_.finite_vertices_begin();
-        vit != shapeControlMesh_.finite_vertices_end();
-        ++vit
-    )
-    {
-        if (vit->real())
-        {
-            fixedAlignments[vit->index()] = vit->alignment();
-        }
-    }
-
-    Info<< nl << "Smoothing alignments" << endl;
-
-    for (label iter = 0; iter < maxSmoothingIterations; iter++)
-    {
-        Info<< "Iteration " << iter;
-
-        meshDistributor().distribute(points);
-        meshDistributor().distribute(alignments);
-
-        scalar residual = 0;
-
-        triadField triadAv(alignments.size(), triad::unset);
-
-        forAll(pointPoints, pI)
-        {
-            const labelList& pPoints = pointPoints[pI];
-
-            if (pPoints.empty())
-            {
-                continue;
-            }
-
-            triad& newTriad = triadAv[pI];
-
-            forAll(pPoints, adjPointI)
-            {
-                const label adjPointIndex = pPoints[adjPointI];
-
-                scalar dist = mag(points[pI] - points[adjPointIndex]);
-
-                triad tmpTriad = alignments[adjPointIndex];
-
-                for (direction dir = 0; dir < 3; dir++)
-                {
-                    if (tmpTriad.set(dir))
-                    {
-                        tmpTriad[dir] *= (1.0/dist);
-                    }
-                }
-
-                newTriad += tmpTriad;
-            }
-
-            newTriad.normalize();
-            newTriad.orthogonalize();
-            newTriad = newTriad.sortxyz();
-
-            // Enforce the boundary conditions
-            const triad& fixedAlignment = fixedAlignments[pI];
-
-            label nFixed = 0;
-
-            forAll(fixedAlignment, dirI)
-            {
-                if (fixedAlignment.set(dirI))
-                {
-                    nFixed++;
-                }
-            }
-
-            if (nFixed == 1)
-            {
-                forAll(fixedAlignment, dirI)
-                {
-                    if (fixedAlignment.set(dirI))
-                    {
-                        newTriad.align(fixedAlignment[dirI]);
-                    }
-                }
-            }
-            else if (nFixed == 2)
-            {
-                forAll(fixedAlignment, dirI)
-                {
-                    if (fixedAlignment.set(dirI))
-                    {
-                        newTriad[dirI] = fixedAlignment[dirI];
-                    }
-                    else
-                    {
-                        newTriad[dirI] = triad::unset[dirI];
-                    }
-                }
-
-                newTriad.orthogonalize();
-            }
-            else if (nFixed == 3)
-            {
-                forAll(fixedAlignment, dirI)
-                {
-                    if (fixedAlignment.set(dirI))
-                    {
-                        newTriad[dirI] = fixedAlignment[dirI];
-                    }
-                }
-            }
-
-            const triad& oldTriad = alignments[pI];
-
-            for (direction dir = 0; dir < 3; ++dir)
-            {
-                if
-                (
-                    newTriad.set(dir)
-                 && oldTriad.set(dir)
-                 //&& !fixedAlignment.set(dir)
-                )
-                {
-                    scalar dotProd = (oldTriad[dir] & newTriad[dir]);
-                    scalar diff = mag(dotProd) - 1.0;
-
-                    residual += mag(diff);
-                }
-            }
-
-//            if (iter == 198 || iter == 199)
-//            {
-//                Info<< "Triad" << nl
-//                    << "    Fixed (" << nFixed << ") = " << fixedAlignment
-//                    << nl
-//                    << "    Old      = " << oldTriad << nl
-//                    << "    Pre-Align= " << triadAv[pI] << nl
-//                    << "    New      = " << newTriad << nl
-//                    << "    Residual = " << residual << endl;
-//            }
-        }
-
-        forAll(alignments, pI)
-        {
-            alignments[pI] = triadAv[pI].sortxyz();
-        }
-
-        reduce(residual, sumOp<scalar>());
-
-        Info<< ", Residual = " << residual << endl;
-
-        if (iter > 0 && residual <= minResidual)
-        {
-            break;
-        }
-    }
-
-    meshDistributor().distribute(alignments);
-
-    for
-    (
-        CellSizeDelaunay::Finite_vertices_iterator vit =
-            shapeControlMesh_.finite_vertices_begin();
-        vit != shapeControlMesh_.finite_vertices_end();
-        ++vit
-    )
-    {
-        if (vit->real())
-        {
-            vit->alignment() = alignments[vit->index()];
-        }
-    }
-
-    labelList referredPoints;
-    autoPtr<mapDistribute> referredDistributor = buildReferredMap
-    (
-        shapeControlMesh_,
-        referredPoints
-    );
-
-    alignments.setSize(shapeControlMesh_.vertexCount());
-    referredDistributor().distribute(alignments);
-
-    label referredI = 0;
-    for
-    (
-        CellSizeDelaunay::Finite_vertices_iterator vit =
-            shapeControlMesh_.finite_vertices_begin();
-        vit != shapeControlMesh_.finite_vertices_end();
-        ++vit
-    )
-    {
-        if (vit->referred())
-        {
-            vit->alignment() = alignments[referredPoints[referredI++]];
+                "Foam::conformalVoronoiMesh::setVertexSizeAndAlignment()"
+            )   << "Point has bad alignment! "
+                << pt << " " << size << " " << alignment << nl
+                << "Bary Coords = " << bary <<  nl
+                << ch->vertex(0)->info() << nl
+                << ch->vertex(1)->info() << nl
+                << ch->vertex(2)->info() << nl
+                << ch->vertex(3)->info()
+                << abort(FatalError);
         }
     }
 }
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.H
index 040ff6f2e0..5c59f648e5 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.H
@@ -84,33 +84,6 @@ class cellShapeControl
 
     // Private Member Functions
 
-        template <class Triangulation, class Type>
-        tmp<Field<Type> > filterFarPoints
-        (
-            const Triangulation& mesh,
-            const Field<Type>& field
-        );
-
-        template <class Triangulation>
-        autoPtr<mapDistribute> buildMap
-        (
-            const Triangulation& mesh,
-            labelListList& pointPoints
-        );
-
-        template <class Triangulation>
-        autoPtr<mapDistribute> buildReferredMap
-        (
-            const Triangulation& mesh,
-            labelList& indices
-        );
-
-        template <class Triangulation>
-        tmp<triadField> buildAlignmentField(const Triangulation& mesh);
-
-        template <class Triangulation>
-        tmp<pointField> buildPointField(const Triangulation& mesh);
-
         //- Disallow default bitwise copy construct
         cellShapeControl(const cellShapeControl&);
 
@@ -155,6 +128,8 @@ public:
 
             inline const cellSizeAndAlignmentControls& sizeAndAlignment() const;
 
+            inline const scalar& minimumCellSize() const;
+
 
         // Query
 
@@ -172,29 +147,6 @@ public:
                 scalar& size,
                 tensor& alignment
             ) const;
-
-
-        // Edit
-
-            void initialMeshPopulation
-            (
-                const autoPtr<backgroundMeshDecomposition>& decomposition
-            );
-
-            label refineMesh
-            (
-                const autoPtr<backgroundMeshDecomposition>& decomposition
-            );
-
-            void smoothMesh();
-
-            //- Add a control point with a specified size and alignment
-//            virtual void addControlPoint
-//            (
-//                const point& pt,
-//                const scalar& size,
-//                const tensor& alignment
-//            );
 };
 
 
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControlI.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControlI.H
index f7fdebc931..2e20866a80 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControlI.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControlI.H
@@ -52,4 +52,17 @@ Foam::cellShapeControl::aspectRatio() const
 }
 
 
+inline const Foam::cellSizeAndAlignmentControls&
+Foam::cellShapeControl::sizeAndAlignment() const
+{
+    return sizeAndAlignment_;
+}
+
+
+inline const Foam::scalar& Foam::cellShapeControl::minimumCellSize() const
+{
+    return minimumCellSize_;
+}
+
+
 // ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
index ba1b46d4c7..842fbe8b1f 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
@@ -24,6 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "cellShapeControlMesh.H"
+#include "cellSizeAndAlignmentControls.H"
 #include "pointIOField.H"
 #include "scalarIOField.H"
 #include "tensorIOField.H"
@@ -36,9 +37,7 @@ License
 
 namespace Foam
 {
-
 defineTypeNameAndDebug(cellShapeControlMesh, 0);
-
 }
 
 
@@ -537,12 +536,9 @@ void Foam::cellShapeControlMesh::barycentricCoords
 ) const
 {
     // Use the previous cell handle as a hint on where to start searching
-    ch = locate
-    (
-        Point(pt.x(), pt.y(), pt.z())
-    );
+    ch = locate(toPoint<Point>(pt));
 
-    if (!is_infinite(ch))
+    if (dimension() > 2 && !is_infinite(ch))
     {
         oldCellHandle_ = ch;
 
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H
index e6cb37a2c3..aeaef8baaa 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H
@@ -53,6 +53,8 @@ SourceFiles
 namespace Foam
 {
 
+class cellSizeAndAlignmentControls;
+
 /*---------------------------------------------------------------------------*\
                       Class cellShapeControlMesh Declaration
 \*---------------------------------------------------------------------------*/
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl.C
index 31d5d61406..8cbc08d773 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl.C
@@ -48,7 +48,7 @@ Foam::cellSizeAndAlignmentControl::cellSizeAndAlignmentControl
     const Time& runTime,
     const word& name,
     const dictionary& controlFunctionDict,
-    const conformationSurfaces& allGeometry
+    const conformationSurfaces& geometryToConformTo
 )
 :
     runTime_(runTime),
@@ -64,7 +64,7 @@ Foam::cellSizeAndAlignmentControl::New
     const Time& runTime,
     const word& name,
     const dictionary& controlFunctionDict,
-    const conformationSurfaces& allGeometry
+    const conformationSurfaces& geometryToConformTo
 )
 {
     word cellSizeAndAlignmentControlTypeName
@@ -72,7 +72,7 @@ Foam::cellSizeAndAlignmentControl::New
         controlFunctionDict.lookup("type")
     );
 
-    Info<< nl << "Selecting cellSizeAndAlignmentControl "
+    Info<< "    Selecting cellSizeAndAlignmentControl "
         << cellSizeAndAlignmentControlTypeName << endl;
 
     dictionaryConstructorTable::iterator cstrIter =
@@ -101,7 +101,7 @@ Foam::cellSizeAndAlignmentControl::New
             runTime,
             name,
             controlFunctionDict,
-            allGeometry
+            geometryToConformTo
         )
     );
 }
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl.H
index 243211929e..1f8fd878bc 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControl.H
@@ -91,9 +91,9 @@ public:
                 const Time& runTime,
                 const word& name,
                 const dictionary& controlFunctionDict,
-                const conformationSurfaces& allGeometry
+                const conformationSurfaces& geometryToConformTo
             ),
-            (runTime, name, controlFunctionDict, allGeometry)
+            (runTime, name, controlFunctionDict, geometryToConformTo)
         );
 
 
@@ -106,7 +106,7 @@ public:
             const Time& runTime,
             const word& name,
             const dictionary& controlFunctionDict,
-            const conformationSurfaces& allGeometry
+            const conformationSurfaces& geometryToConformTo
         );
 
 
@@ -118,7 +118,7 @@ public:
             const Time& runTime,
             const word& name,
             const dictionary& controlFunctionDict,
-            const conformationSurfaces& allGeometry
+            const conformationSurfaces& geometryToConformTo
         );
 
 
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControls.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControls.C
index 515482689d..d793bba878 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControls.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControls.C
@@ -30,9 +30,7 @@ License
 
 namespace Foam
 {
-
 defineTypeNameAndDebug(cellSizeAndAlignmentControls, 0);
-
 }
 
 
@@ -116,12 +114,12 @@ Foam::cellSizeAndAlignmentControls::cellSizeAndAlignmentControls
 (
     const Time& runTime,
     const dictionary& shapeControlDict,
-    const conformationSurfaces& allGeometry,
+    const conformationSurfaces& geometryToConformTo,
     const scalar defaultCellSize
 )
 :
     shapeControlDict_(shapeControlDict),
-    allGeometry_(allGeometry),
+    geometryToConformTo_(geometryToConformTo),
     controlFunctions_(shapeControlDict_.size()),
     defaultCellSize_(defaultCellSize)
 {
@@ -136,6 +134,8 @@ Foam::cellSizeAndAlignmentControls::cellSizeAndAlignmentControls
             shapeControlDict_.subDict(shapeControlEntryName)
         );
 
+        Info<< nl << "Shape Control : " << shapeControlEntryName << endl;
+
         controlFunctions_.set
         (
             functionI,
@@ -144,7 +144,7 @@ Foam::cellSizeAndAlignmentControls::cellSizeAndAlignmentControls
                 runTime,
                 shapeControlEntryName,
                 controlFunctionDict,
-                allGeometry
+                geometryToConformTo
             )
         );
 
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControls.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControls.H
index 76e6c801dd..7afa3ff82a 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControls.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/cellSizeAndAlignmentControls.H
@@ -52,7 +52,7 @@ class cellSizeAndAlignmentControls
 
         const dictionary& shapeControlDict_;
 
-        const conformationSurfaces& allGeometry_;
+        const conformationSurfaces& geometryToConformTo_;
 
         PtrList<cellSizeAndAlignmentControl> controlFunctions_;
 
@@ -83,7 +83,7 @@ public:
         (
             const Time& runTime,
             const dictionary& shapeControlDict,
-            const conformationSurfaces& allGeometry,
+            const conformationSurfaces& geometryToConformTo,
             const scalar defaultCellSize
         );
 
@@ -96,11 +96,17 @@ public:
 
         // Access
 
-            const PtrList<cellSizeAndAlignmentControl>& controlFunctions() const
+            inline const PtrList<cellSizeAndAlignmentControl>&
+            controlFunctions() const
             {
                 return controlFunctions_;
             }
 
+            inline const conformationSurfaces& geometryToConformTo() const
+            {
+                return geometryToConformTo_;
+            }
+
 
         // Query
 
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/fileControl/fileControl.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/fileControl/fileControl.C
index 9b86dbf2ba..63ac5960a1 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/fileControl/fileControl.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/fileControl/fileControl.C
@@ -58,7 +58,7 @@ Foam::fileControl::fileControl
     const Time& runTime,
     const word& name,
     const dictionary& controlFunctionDict,
-    const conformationSurfaces& allGeometry
+    const conformationSurfaces& geometryToConformTo
 )
 :
     cellSizeAndAlignmentControl
@@ -66,7 +66,7 @@ Foam::fileControl::fileControl
         runTime,
         name,
         controlFunctionDict,
-        allGeometry
+        geometryToConformTo
     ),
     pointsFile_(controlFunctionDict.lookup("pointsFile")),
     sizesFile_(controlFunctionDict.lookup("sizesFile")),
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/fileControl/fileControl.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/fileControl/fileControl.H
index f06892680d..93c631e90e 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/fileControl/fileControl.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/fileControl/fileControl.H
@@ -85,7 +85,7 @@ public:
             const Time& runTime,
             const word& name,
             const dictionary& controlFunctionDict,
-            const conformationSurfaces& allGeometry
+            const conformationSurfaces& geometryToConformTo
         );
 
     //- Destructor
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.C
index 0f2a4f6d8c..72251e88ce 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.C
@@ -158,7 +158,7 @@ Foam::searchableSurfaceControl::searchableSurfaceControl
     const Time& runTime,
     const word& name,
     const dictionary& controlFunctionDict,
-    const conformationSurfaces& allGeometry
+    const conformationSurfaces& geometryToConformTo
 )
 :
     cellSizeAndAlignmentControl
@@ -166,11 +166,11 @@ Foam::searchableSurfaceControl::searchableSurfaceControl
         runTime,
         name,
         controlFunctionDict,
-        allGeometry
+        geometryToConformTo
     ),
     surfaceName_(controlFunctionDict.lookupOrDefault<word>("surface", name)),
-    searchableSurface_(allGeometry.geometry()[surfaceName_]),
-    allGeometry_(allGeometry),
+    searchableSurface_(geometryToConformTo.geometry()[surfaceName_]),
+    geometryToConformTo_(geometryToConformTo),
     cellSizeFunction_
     (
         cellSizeFunction::New(controlFunctionDict, searchableSurface_)
@@ -579,7 +579,7 @@ void Foam::searchableSurfaceControl::initialVertices
 
         pointIndexHit info;
         label infoFeature;
-        allGeometry_.findFeaturePointNearest
+        geometryToConformTo_.findFeaturePointNearest
         (
             pts[pI],
             nearFeatDistSqrCoeff,
@@ -592,7 +592,7 @@ void Foam::searchableSurfaceControl::initialVertices
         if (info.hit())
         {
             const extendedFeatureEdgeMesh& features =
-                allGeometry_.features()[infoFeature];
+                geometryToConformTo_.features()[infoFeature];
 
             vectorField norms = features.featurePointNormals(info.index());
 
@@ -608,7 +608,7 @@ void Foam::searchableSurfaceControl::initialVertices
         }
         else
         {
-            allGeometry_.findEdgeNearest
+            geometryToConformTo_.findEdgeNearest
             (
                 pts[pI],
                 nearFeatDistSqrCoeff,
@@ -619,7 +619,7 @@ void Foam::searchableSurfaceControl::initialVertices
             if (info.hit())
             {
                 const extendedFeatureEdgeMesh& features =
-                    allGeometry_.features()[infoFeature];
+                    geometryToConformTo_.features()[infoFeature];
 
                 vectorField norms = features.edgeNormals(info.index());
 
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.H
index 7edc26cdc4..a6f9b0f24e 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.H
@@ -59,7 +59,7 @@ class searchableSurfaceControl
         //- Reference to the searchableSurface object holding the geometry data
         const searchableSurface& searchableSurface_;
 
-        const conformationSurfaces& allGeometry_;
+        const conformationSurfaces& geometryToConformTo_;
 
         autoPtr<cellSizeFunction> cellSizeFunction_;
 
@@ -106,7 +106,7 @@ public:
             const Time& runTime,
             const word& name,
             const dictionary& controlFunctionDict,
-            const conformationSurfaces& allGeometry
+            const conformationSurfaces& geometryToConformTo
         );
 
     //- Destructor
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C
new file mode 100644
index 0000000000..e6978a1ae1
--- /dev/null
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C
@@ -0,0 +1,683 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "controlMeshRefinement.H"
+#include "cellSizeAndAlignmentControl.H"
+#include "OFstream.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+defineTypeNameAndDebug(controlMeshRefinement, 0);
+}
+
+// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+bool Foam::controlMeshRefinement::detectEdge
+(
+    point a,
+    point b,
+    DynamicList<Vb>& pointsFound
+) const
+{
+    const Foam::point midpoint = 0.5*(a + b);
+    const scalar h = mag(a - b);
+
+    scalar cellSizeA = sizeControls_.cellSize(a);
+    scalar cellSizeMid = sizeControls_.cellSize(midpoint);
+    scalar cellSizeB = sizeControls_.cellSize(b);
+
+    scalar firstDerivative = (cellSizeA - cellSizeB)/h;
+
+    scalar secondDerivative =
+        (
+            cellSizeB
+          - 2*cellSizeMid
+          + cellSizeA
+        )
+       /magSqr(h/2.0);
+
+//    Info<< a << " " << midpoint << " " << b << endl;
+//    Info<< " length = " << h << endl;
+//    Info<< "size(a) = " << cellSizeA << endl;
+//    Info<< "size(m) = " << cellSizeMid << endl;
+//    Info<< "size(b) = " << cellSizeB << endl;
+//    Info<< "    f'  = " << firstDerivative << endl;
+//    Info<< "    f'' = " << secondDerivative << endl;
+
+    if (mag(secondDerivative) < 1e-1)
+    {
+        return false;
+    }
+
+    if (mag(a - b) < 1e-1*cellSizeMid)
+    {
+        if (geometryToConformTo_.outside(midpoint))
+        {
+            return false;
+        }
+        else
+        {
+//            Info<< "    Point Added" << endl;
+
+            pointsFound.append
+            (
+                Vb
+                (
+                    toPoint<cellShapeControlMesh::Point>(midpoint),
+                    Vb::vtInternal
+                )
+            );
+            pointsFound.last().targetCellSize() =
+                sizeControls_.cellSize(midpoint);
+            pointsFound.last().alignment() = triad::unset;
+
+            return true;
+        }
+    }
+
+    detectEdge(a, midpoint, pointsFound);
+    detectEdge(midpoint, b, pointsFound);
+}
+
+
+Foam::DynamicList<Vb>
+Foam::controlMeshRefinement::findDiscontinuities(const linePointRef& l) const
+{
+    DynamicList<Vb> verts(2);
+
+    // Divide up the line into a reasonable number of intervals
+
+
+    // Call detectEdge on each interval
+    detectEdge(l.start(), l.end(), verts);
+
+    if (verts.size())
+    {
+        Info<< l << " has " << verts.size() << " discontinuities" << endl;
+
+        OFstream str("plotSizes_" + name(verts.size()));
+
+        forAll(verts, vI)
+        {
+            scalar x =
+                mag(topoint(verts[vI].point()) - l.start())
+               /mag(l.end() - l.start());
+
+            str << x << " " << verts[vI].targetCellSize() << endl;
+        }
+
+//        forAll(verts, vI)
+//        {
+//            Info<< verts[vI].info();
+//        }
+    }
+
+    return verts;
+}
+
+
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::controlMeshRefinement::controlMeshRefinement
+(
+    cellShapeControl& shapeController
+)
+:
+    shapeController_(shapeController),
+    mesh_(shapeController.shapeControlMesh()),
+    sizeControls_(shapeController.sizeAndAlignment()),
+    geometryToConformTo_(sizeControls_.geometryToConformTo())
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::controlMeshRefinement::~controlMeshRefinement()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+void Foam::controlMeshRefinement::initialMeshPopulation
+(
+    const autoPtr<backgroundMeshDecomposition>& decomposition
+)
+{
+    autoPtr<boundBox> overallBoundBox;
+
+    // Need to pass in the background mesh decomposition so that can test if
+    // a point to insert is on the processor.
+    if (Pstream::parRun())
+    {
+        overallBoundBox.set(new boundBox(decomposition().procBounds()));
+    }
+    else
+    {
+        overallBoundBox.set
+        (
+            new boundBox(geometryToConformTo_.geometry().bounds())
+        );
+    }
+
+    mesh_.insertBoundingPoints
+    (
+        overallBoundBox(),
+        sizeControls_
+    );
+
+    const PtrList<cellSizeAndAlignmentControl>& controlFunctions =
+        sizeControls_.controlFunctions();
+
+    forAll(controlFunctions, fI)
+    {
+        const cellSizeAndAlignmentControl& controlFunction =
+            controlFunctions[fI];
+
+        Info<< "Inserting points from " << controlFunction.name()
+            << " (" << controlFunction.type() << ")" << endl;
+
+        pointField pts;
+        scalarField sizes;
+        triadField alignments;
+
+        controlFunction.initialVertices(pts, sizes, alignments);
+
+        List<Vb> vertices(pts.size());
+
+        // Clip the minimum size
+        forAll(vertices, vI)
+        {
+            vertices[vI] = Vb(pts[vI], Vb::vtInternalNearBoundary);
+            vertices[vI].targetCellSize() = max
+            (
+                sizes[vI],
+                shapeController_.minimumCellSize()
+            );
+            vertices[vI].alignment() = alignments[vI];
+        }
+
+        pts.clear();
+        sizes.clear();
+        alignments.clear();
+
+        label nRejected = 0;
+
+        PackedBoolList keepVertex(vertices.size(), true);
+
+        if (Pstream::parRun())
+        {
+            forAll(vertices, vI)
+            {
+                const bool onProc = decomposition().positionOnThisProcessor
+                (
+                    topoint(vertices[vI].point())
+                );
+
+                if (!onProc)
+                {
+                    keepVertex[vI] = false;
+                }
+            }
+        }
+
+        inplaceSubset(keepVertex, vertices);
+
+        const label preInsertedSize = mesh_.number_of_vertices();
+
+        forAll(vertices, vI)
+        {
+            bool insertPoint = false;
+
+            pointFromPoint pt(topoint(vertices[vI].point()));
+
+            if
+            (
+                mesh_.dimension() < 3
+             || mesh_.is_infinite
+                (
+                    mesh_.locate(vertices[vI].point())
+                )
+            )
+            {
+                insertPoint = true;
+            }
+
+            const scalar interpolatedCellSize = shapeController_.cellSize(pt);
+            const triad interpolatedAlignment =
+                shapeController_.cellAlignment(pt);
+            const scalar calculatedCellSize = vertices[vI].targetCellSize();
+            const triad calculatedAlignment = vertices[vI].alignment();
+
+            if (debug)
+            {
+                Info<< "Point = " << pt << nl
+                    << "  Size(interp) = " << interpolatedCellSize << nl
+                    << "    Size(calc) = " << calculatedCellSize << nl
+                    << " Align(interp) = " << interpolatedAlignment << nl
+                    << "   Align(calc) = " << calculatedAlignment << nl
+                    << endl;
+            }
+
+            const scalar sizeDiff =
+                mag(interpolatedCellSize - calculatedCellSize);
+            const scalar alignmentDiff =
+                diff(interpolatedAlignment, calculatedAlignment);
+
+            if (debug)
+            {
+                Info<< "    size difference = " << sizeDiff
+                    << ", alignment difference = " << alignmentDiff << endl;
+            }
+
+            // @todo Also need to base it on the alignments
+            if
+            (
+                sizeDiff/interpolatedCellSize > 0.1
+             || alignmentDiff > 0.15
+            )
+            {
+                insertPoint = true;
+            }
+
+            if (insertPoint)
+            {
+                mesh_.insert
+                (
+                    pt,
+                    calculatedCellSize,
+                    vertices[vI].alignment(),
+                    Vb::vtInternalNearBoundary
+                );
+            }
+        }
+
+        //mesh_.rangeInsertWithInfo(vertices.begin(), vertices.end());
+
+        Info<< "    Inserted "
+            << returnReduce
+               (
+                   label(mesh_.number_of_vertices()) - preInsertedSize,
+                   sumOp<label>()
+               )
+            << "/" << returnReduce(vertices.size(), sumOp<label>())
+            << endl;
+    }
+}
+
+
+Foam::label Foam::controlMeshRefinement::refineMesh
+(
+    const autoPtr<backgroundMeshDecomposition>& decomposition
+)
+{
+    Info<< "Iterate over cell size mesh edges" << endl;
+
+    DynamicList<Vb> verts(mesh_.number_of_vertices());
+
+    for
+    (
+        CellSizeDelaunay::Finite_edges_iterator eit =
+            mesh_.finite_edges_begin();
+        eit != mesh_.finite_edges_end();
+        ++eit
+    )
+    {
+        CellSizeDelaunay::Cell_handle c = eit->first;
+        CellSizeDelaunay::Vertex_handle vA = c->vertex(eit->second);
+        CellSizeDelaunay::Vertex_handle vB = c->vertex(eit->third);
+
+        if
+        (
+            mesh_.is_infinite(vA)
+         || mesh_.is_infinite(vB)
+         || (vA->referred() && vB->referred())
+         || (vA->referred() && vA->procIndex() > vB->procIndex())
+         || (vB->referred() && vB->procIndex() > vA->procIndex())
+        )
+        {
+            continue;
+        }
+
+        pointFromPoint ptA(topoint(vA->point()));
+        pointFromPoint ptB(topoint(vB->point()));
+
+        linePointRef l(ptA, ptB);
+
+//        Info<< "Edge " << l << endl;
+//
+//        Info<< vA->info() << vB->info();
+
+        const Foam::point midPoint = l.centre();
+        const scalar length = l.mag();
+
+        scalar evaluatedSize = sizeControls_.cellSize(midPoint);
+        scalar interpolatedSize =
+            (vA->targetCellSize() + vB->targetCellSize())/2;
+
+        const scalar diff = mag(interpolatedSize - evaluatedSize);
+
+//        Info<< "       Evaluated size = " << evaluatedSize << nl
+//            << "    Interpolated size = " << interpolatedSize << nl
+//            << "         vA cell size = " << vA->targetCellSize() << nl
+//            << "         vB cell size = " << vB->targetCellSize() << endl;
+
+//        verts.append(findDiscontinuities(l));
+
+        if (diff/interpolatedSize > 0.1)
+        {
+            // Create a new point
+            // Walk along the edge, placing points where the gradient in cell
+            // size changes
+            if (vA->targetCellSize() <= vB->targetCellSize())
+            {
+                // Walk from vA
+                const label nPoints = length/vA->targetCellSize();
+
+                if (nPoints == 0)
+                {
+                    continue;
+                }
+
+                scalar prevSize = vA->targetCellSize();
+
+                for (label pI = 1; pI < nPoints; ++pI)
+                {
+                    const Foam::point testPt =
+                        topoint(vA->point()) + pI*l.vec()/nPoints;
+
+//                    if (geometryToConformTo_.outside(testPt))
+//                    {
+//                        continue;
+//                    }
+
+                    scalar testPtSize = sizeControls_.cellSize(testPt);
+
+                    if (mag(testPtSize - prevSize) < 1e-3*testPtSize)
+                    {
+//                        Info<< "Adding " << testPt << " " << testPtSize
+//                            << endl;
+                        verts.append
+                        (
+                            Vb
+                            (
+                                toPoint<cellShapeControlMesh::Point>(testPt),
+                                Vb::vtInternal
+                            )
+                        );
+                        verts.last().targetCellSize() = testPtSize;
+                        verts.last().alignment() = triad::unset;
+
+                        break;
+                    }
+
+                    prevSize = testPtSize;
+                }
+            }
+            else
+            {
+                // Walk from vB
+                const label nPoints = length/vB->targetCellSize();
+
+                if (nPoints == 0)
+                {
+                    continue;
+                }
+
+                scalar prevSize = vA->targetCellSize();
+
+                for (label pI = 1; pI < nPoints; ++pI)
+                {
+                    const Foam::point testPt =
+                        topoint(vB->point()) - pI*l.vec()/nPoints;
+
+//                    if (geometryToConformTo_.outside(testPt))
+//                    {
+//                        continue;
+//                    }
+
+                    scalar testPtSize = sizeControls_.cellSize(testPt);
+
+                    if (mag(testPtSize - prevSize) < 1e-3*testPtSize)
+                    //if (testPtSize > prevSize)// || testPtSize < prevSize)
+                    {
+//                        Info<< "Adding " << testPt << " " << testPtSize
+//                            << endl;
+                        verts.append
+                        (
+                            Vb
+                            (
+                                toPoint<cellShapeControlMesh::Point>(testPt),
+                                Vb::vtInternal
+                            )
+                        );
+                        verts.last().targetCellSize() = testPtSize;
+                        verts.last().alignment() = triad::unset;
+
+                        break;
+                    }
+
+                    prevSize = testPtSize;
+                }
+            }
+        }
+    }
+
+//    const pointField cellCentres(mesh_.cellCentres());
+//
+//    Info<< "    Created cell centres" << endl;
+//
+//    const PtrList<cellSizeAndAlignmentControl>& controlFunctions =
+//        sizeControls_.controlFunctions();
+//
+//    DynamicList<Vb> verts(mesh_.number_of_vertices());
+//
+//    forAll(cellCentres, cellI)
+//    {
+//        Foam::point pt = cellCentres[cellI];
+//
+//        if (!geometryToConformTo_.inside(pt))
+//        {
+//            continue;
+//        }
+//
+//        scalarList bary;
+//        cellShapeControlMesh::Cell_handle ch;
+//
+//        mesh_.barycentricCoords(pt, bary, ch);
+//
+//        if (mesh_.is_infinite(ch))
+//        {
+//            continue;
+//        }
+//
+//        scalar interpolatedSize = 0;
+//        forAll(bary, pI)
+//        {
+//            interpolatedSize += bary[pI]*ch->vertex(pI)->targetCellSize();
+//        }
+//
+//        label lastPriority = labelMin;
+//        scalar lastCellSize = GREAT;
+//        forAll(controlFunctions, fI)
+//        {
+//            const cellSizeAndAlignmentControl& controlFunction =
+//                controlFunctions[fI];
+//
+//            if (controlFunction.priority() < lastPriority)
+//            {
+//                continue;
+//            }
+//
+//            if (isA<searchableSurfaceControl>(controlFunction))
+//            {
+//                const cellSizeFunction& sizeFunction =
+//                    dynamicCast<const searchableSurfaceControl>
+//                    (
+//                        controlFunction
+//                    ).sizeFunction();
+//
+//                scalar cellSize = 0;
+//                if (sizeFunction.cellSize(pt, cellSize))
+//                {
+//                    if (controlFunction.priority() == lastPriority)
+//                    {
+//                        if (cellSize < lastCellSize)
+//                        {
+//                            lastCellSize = cellSize;
+//                        }
+//                    }
+//                    else
+//                    {
+//                        lastCellSize = cellSize;
+//                    }
+//
+//                    lastPriority = controlFunction.priority();
+//                }
+//            }
+//        }
+//
+//        const scalar diff = mag(interpolatedSize - lastCellSize);
+//
+//        if (diff/interpolatedSize > 0.1)
+//        {
+//            bool pointFound = false;
+//
+//            // Travel along lines to each tet point
+//            for (label vI = 0; vI < 4; ++vI)
+//            {
+//                if (pointFound)
+//                {
+//                    break;
+//                }
+//
+//                pointFromPoint tetVertex = topoint(ch->vertex(vI)->point());
+//                scalar tetVertexSize = sizeControls_.cellSize(tetVertex);
+//
+//                if (tetVertexSize < interpolatedSize)
+//                {
+//                    linePointRef l(tetVertex, pt);
+//
+//                    label nTestPoints = l.mag()/tetVertexSize;
+//
+//                    for (label i = 1; i < nTestPoints; ++i)
+//                    {
+//                        const Foam::point testPt
+//                        (
+//                            tetVertex
+//                          + l.vec()*i/nTestPoints
+//                        );
+//
+//                        scalar testPtSize = sizeControls_.cellSize(testPt);
+//
+//                        if (testPtSize > tetVertexSize)
+//                        {
+//                            pt = testPt;
+//                            lastCellSize = testPtSize;
+//                            pointFound = true;
+//                            break;
+//                        }
+//                    }
+//                }
+//                else
+//                {
+//                    linePointRef l(pt, tetVertex);
+//
+//                    label nTestPoints = l.mag()/tetVertexSize;
+//
+//                    for (label i = 1; i < nTestPoints; ++i)
+//                    {
+//                        const Foam::point testPt
+//                        (
+//                            tetVertex
+//                          + l.vec()*i/nTestPoints
+//                        );
+//
+//                        scalar testPtSize = sizeControls_.cellSize(testPt);
+//
+//                        if (testPtSize < interpolatedSize)
+//                        {
+//                            pt = testPt;
+//                            lastCellSize = testPtSize;
+//                            pointFound = true;
+//                            break;
+//                        }
+//                    }
+//                }
+//            }
+//        }
+//
+//        // Decide whether to insert or not
+//        if (lastCellSize < GREAT)
+//        {
+//            if (!geometryToConformTo_.inside(pt))
+//            {
+//                continue;
+//            }
+//
+//            if (Pstream::parRun())
+//            {
+//                if (!decomposition().positionOnThisProcessor(pt))
+//                {
+//                    continue;
+//                }
+//            }
+//
+//            verts.append
+//            (
+//                Vb
+//                (
+//                    toPoint<cellShapeControlMesh::Point>(pt),
+//                    Vb::vtInternal
+//                )
+//            );
+//            verts.last().targetCellSize() = lastCellSize;
+//            verts.last().alignment() = triad::unset;
+//        }
+//    }
+
+    mesh_.insertPoints(verts);
+
+    return verts.size();
+}
+
+
+// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * //
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.H
new file mode 100644
index 0000000000..5f42fcf051
--- /dev/null
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.H
@@ -0,0 +1,133 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::controlMeshRefinement
+
+Description
+
+SourceFiles
+    controlMeshRefinementI.H
+    controlMeshRefinement.C
+    controlMeshRefinementIO.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef controlMeshRefinement_H
+#define controlMeshRefinement_H
+
+#include "cellShapeControl.H"
+#include "cellShapeControlMesh.H"
+#include "cellSizeAndAlignmentControls.H"
+#include "conformationSurfaces.H"
+#include "backgroundMeshDecomposition.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class controlMeshRefinement Declaration
+\*---------------------------------------------------------------------------*/
+
+class controlMeshRefinement
+{
+    // Private data
+
+        const cellShapeControl& shapeController_;
+
+        cellShapeControlMesh& mesh_;
+
+        const cellSizeAndAlignmentControls& sizeControls_;
+
+        const conformationSurfaces& geometryToConformTo_;
+
+
+    // Private Member Functions
+
+        bool detectEdge
+        (
+            point a,
+            point b,
+            DynamicList<Vb>& pointsFound
+        ) const;
+
+        DynamicList<Vb> findDiscontinuities(const linePointRef& l) const;
+
+        //- Disallow default bitwise copy construct
+        controlMeshRefinement(const controlMeshRefinement&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const controlMeshRefinement&);
+
+
+public:
+
+        //- Runtime type information
+        ClassName("controlMeshRefinement");
+
+
+    // Constructors
+
+        //- Construct null
+        controlMeshRefinement(cellShapeControl& shapeController);
+
+
+    //- Destructor
+    ~controlMeshRefinement();
+
+
+    // Member Functions
+
+        // Access
+
+        // Check
+
+        // Edit
+
+            void initialMeshPopulation
+            (
+                const autoPtr<backgroundMeshDecomposition>& decomposition
+            );
+
+            label refineMesh
+            (
+                const autoPtr<backgroundMeshDecomposition>& decomposition
+            );
+
+
+        // Write
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/smoothAlignmentSolver/smoothAlignmentSolver.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/smoothAlignmentSolver/smoothAlignmentSolver.C
new file mode 100644
index 0000000000..39eaf8b581
--- /dev/null
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/smoothAlignmentSolver/smoothAlignmentSolver.C
@@ -0,0 +1,488 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "smoothAlignmentSolver.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template <class Triangulation, class Type>
+Foam::tmp<Foam::Field<Type> > Foam::smoothAlignmentSolver::filterFarPoints
+(
+    const Triangulation& mesh,
+    const Field<Type>& field
+)
+{
+    tmp<Field<Type> > tNewField(new Field<Type>(field.size()));
+    Field<Type>& newField = tNewField();
+
+    label added = 0;
+    label count = 0;
+
+    for
+    (
+        typename Triangulation::Finite_vertices_iterator vit =
+            mesh.finite_vertices_begin();
+        vit != mesh.finite_vertices_end();
+        ++vit
+    )
+    {
+        if (vit->real())
+        {
+            newField[added++] = field[count];
+        }
+
+        count++;
+    }
+
+    newField.resize(added);
+
+    return tNewField;
+}
+
+
+template <class Triangulation>
+Foam::autoPtr<Foam::mapDistribute> Foam::smoothAlignmentSolver::buildReferredMap
+(
+    const Triangulation& mesh,
+    labelList& indices
+)
+{
+    globalIndex globalIndexing(mesh.vertexCount());
+
+    DynamicList<label> dynIndices(mesh.vertexCount()/10);
+
+    for
+    (
+        typename Triangulation::Finite_vertices_iterator vit =
+            mesh.finite_vertices_begin();
+        vit != mesh.finite_vertices_end();
+        ++vit
+    )
+    {
+        if (vit->referred())
+        {
+            dynIndices.append
+            (
+                globalIndexing.toGlobal(vit->procIndex(), vit->index())
+            );
+        }
+    }
+
+    indices.transfer(dynIndices);
+
+    List<Map<label> > compactMap;
+    return autoPtr<mapDistribute>
+    (
+        new mapDistribute
+        (
+            globalIndexing,
+            indices,
+            compactMap
+        )
+    );
+}
+
+
+template <class Triangulation>
+Foam::autoPtr<Foam::mapDistribute> Foam::smoothAlignmentSolver::buildMap
+(
+    const Triangulation& mesh,
+    labelListList& pointPoints
+)
+{
+    pointPoints.setSize(mesh.vertexCount());
+
+    globalIndex globalIndexing(mesh.vertexCount());
+
+    for
+    (
+        typename Triangulation::Finite_vertices_iterator vit =
+            mesh.finite_vertices_begin();
+        vit != mesh.finite_vertices_end();
+        ++vit
+    )
+    {
+        if (!vit->real())
+        {
+            continue;
+        }
+
+        std::list<typename Triangulation::Vertex_handle> adjVerts;
+        mesh.finite_adjacent_vertices(vit, std::back_inserter(adjVerts));
+
+        DynamicList<label> indices(adjVerts.size());
+
+        for
+        (
+            typename std::list<typename Triangulation::Vertex_handle>::
+                const_iterator adjVertI = adjVerts.begin();
+            adjVertI != adjVerts.end();
+            ++adjVertI
+        )
+        {
+            typename Triangulation::Vertex_handle vh = *adjVertI;
+
+            if (!vh->farPoint())
+            {
+                indices.append
+                (
+                    globalIndexing.toGlobal(vh->procIndex(), vh->index())
+                );
+            }
+        }
+
+        pointPoints[vit->index()].transfer(indices);
+    }
+
+    List<Map<label> > compactMap;
+    return autoPtr<mapDistribute>
+    (
+        new mapDistribute
+        (
+            globalIndexing,
+            pointPoints,
+            compactMap
+        )
+    );
+}
+
+
+template <class Triangulation>
+Foam::tmp<Foam::triadField> Foam::smoothAlignmentSolver::buildAlignmentField
+(
+    const Triangulation& mesh
+)
+{
+    tmp<triadField> tAlignments
+    (
+        new triadField(mesh.vertexCount(), triad::unset)
+    );
+    triadField& alignments = tAlignments();
+
+    for
+    (
+        typename Triangulation::Finite_vertices_iterator vit =
+            mesh.finite_vertices_begin();
+        vit != mesh.finite_vertices_end();
+        ++vit
+    )
+    {
+        if (!vit->real())
+        {
+            continue;
+        }
+
+        alignments[vit->index()] = vit->alignment();
+    }
+
+    return tAlignments;
+}
+
+
+template <class Triangulation>
+Foam::tmp<Foam::pointField> Foam::smoothAlignmentSolver::buildPointField
+(
+    const Triangulation& mesh
+)
+{
+    tmp<pointField> tPoints
+    (
+        new pointField(mesh.vertexCount(), point(GREAT, GREAT, GREAT))
+    );
+    pointField& points = tPoints();
+
+    for
+    (
+        typename Triangulation::Finite_vertices_iterator vit =
+            mesh.finite_vertices_begin();
+        vit != mesh.finite_vertices_end();
+        ++vit
+    )
+    {
+        if (!vit->real())
+        {
+            continue;
+        }
+
+        points[vit->index()] = topoint(vit->point());
+    }
+
+    return tPoints;
+}
+
+
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::smoothAlignmentSolver::smoothAlignmentSolver(cellShapeControlMesh& mesh)
+:
+    mesh_(mesh)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::smoothAlignmentSolver::~smoothAlignmentSolver()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+void Foam::smoothAlignmentSolver::smoothAlignments
+(
+    const label maxSmoothingIterations
+)
+{
+    scalar minResidual = 0;
+
+    labelListList pointPoints;
+    autoPtr<mapDistribute> meshDistributor = buildMap
+    (
+        mesh_,
+        pointPoints
+    );
+
+    triadField alignments(buildAlignmentField(mesh_));
+    pointField points(buildPointField(mesh_));
+    // Setup the sizes and alignments on each point
+    triadField fixedAlignments(mesh_.vertexCount(), triad::unset);
+
+    for
+    (
+        CellSizeDelaunay::Finite_vertices_iterator vit =
+            mesh_.finite_vertices_begin();
+        vit != mesh_.finite_vertices_end();
+        ++vit
+    )
+    {
+        if (vit->real())
+        {
+            fixedAlignments[vit->index()] = vit->alignment();
+        }
+    }
+
+    Info<< nl << "Smoothing alignments" << endl;
+
+    for (label iter = 0; iter < maxSmoothingIterations; iter++)
+    {
+        Info<< "Iteration " << iter;
+
+        meshDistributor().distribute(points);
+        meshDistributor().distribute(alignments);
+
+        scalar residual = 0;
+
+        triadField triadAv(alignments.size(), triad::unset);
+
+        forAll(pointPoints, pI)
+        {
+            const labelList& pPoints = pointPoints[pI];
+
+            if (pPoints.empty())
+            {
+                continue;
+            }
+
+            triad& newTriad = triadAv[pI];
+
+            forAll(pPoints, adjPointI)
+            {
+                const label adjPointIndex = pPoints[adjPointI];
+
+                scalar dist = mag(points[pI] - points[adjPointIndex]);
+
+                triad tmpTriad = alignments[adjPointIndex];
+
+                for (direction dir = 0; dir < 3; dir++)
+                {
+                    if (tmpTriad.set(dir))
+                    {
+                        tmpTriad[dir] *= 1.0/(dist + SMALL);
+                    }
+                }
+
+                newTriad += tmpTriad;
+            }
+
+            newTriad.normalize();
+            newTriad.orthogonalize();
+//            newTriad = newTriad.sortxyz();
+
+            // Enforce the boundary conditions
+            const triad& fixedAlignment = fixedAlignments[pI];
+
+            label nFixed = 0;
+
+            forAll(fixedAlignment, dirI)
+            {
+                if (fixedAlignment.set(dirI))
+                {
+                    nFixed++;
+                }
+            }
+
+            if (nFixed == 1)
+            {
+                forAll(fixedAlignment, dirI)
+                {
+                    if (fixedAlignment.set(dirI))
+                    {
+                        newTriad.align(fixedAlignment[dirI]);
+                    }
+                }
+            }
+            else if (nFixed == 2)
+            {
+                forAll(fixedAlignment, dirI)
+                {
+                    if (fixedAlignment.set(dirI))
+                    {
+                        newTriad[dirI] = fixedAlignment[dirI];
+                    }
+                    else
+                    {
+                        newTriad[dirI] = triad::unset[dirI];
+                    }
+                }
+
+                newTriad.orthogonalize();
+            }
+            else if (nFixed == 3)
+            {
+                forAll(fixedAlignment, dirI)
+                {
+                    if (fixedAlignment.set(dirI))
+                    {
+                        newTriad[dirI] = fixedAlignment[dirI];
+                    }
+                }
+            }
+
+            const triad& oldTriad = alignments[pI];
+
+            for (direction dir = 0; dir < 3; ++dir)
+            {
+                if
+                (
+                    newTriad.set(dir)
+                 && oldTriad.set(dir)
+                 && !fixedAlignment.set(dir)
+                )
+                {
+                    residual += diff(oldTriad, newTriad);
+                }
+            }
+
+//            if (iter == 198 || iter == 199)
+//            {
+//                Info<< "Triad" << nl
+//                    << "    Fixed (" << nFixed << ") = " << fixedAlignment
+//                    << nl
+//                    << "    Old      = " << oldTriad << nl
+//                    << "    Pre-Align= " << triadAv[pI] << nl
+//                    << "    New      = " << newTriad << nl
+//                    << "    Residual = " << residual << endl;
+//            }
+        }
+
+        forAll(alignments, pI)
+        {
+            alignments[pI] = triadAv[pI].sortxyz();
+        }
+
+        reduce(residual, sumOp<scalar>());
+
+        Info<< ", Residual = "
+            << residual
+              /returnReduce(points.size(), sumOp<label>())
+            << endl;
+
+        if (iter > 0 && residual <= minResidual)
+        {
+            break;
+        }
+    }
+
+    meshDistributor().distribute(alignments);
+
+    for
+    (
+        CellSizeDelaunay::Finite_vertices_iterator vit =
+            mesh_.finite_vertices_begin();
+        vit != mesh_.finite_vertices_end();
+        ++vit
+    )
+    {
+        if (vit->real())
+        {
+            vit->alignment() = alignments[vit->index()];
+        }
+    }
+
+    labelList referredPoints;
+    autoPtr<mapDistribute> referredDistributor = buildReferredMap
+    (
+        mesh_,
+        referredPoints
+    );
+
+    alignments.setSize(mesh_.vertexCount());
+    referredDistributor().distribute(alignments);
+
+    label referredI = 0;
+    for
+    (
+        CellSizeDelaunay::Finite_vertices_iterator vit =
+            mesh_.finite_vertices_begin();
+        vit != mesh_.finite_vertices_end();
+        ++vit
+    )
+    {
+        if (vit->referred())
+        {
+            vit->alignment() = alignments[referredPoints[referredI++]];
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * //
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/smoothAlignmentSolver/smoothAlignmentSolver.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/smoothAlignmentSolver/smoothAlignmentSolver.H
new file mode 100644
index 0000000000..fa158605c2
--- /dev/null
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellShapeControl/smoothAlignmentSolver/smoothAlignmentSolver.H
@@ -0,0 +1,139 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::smoothAlignmentSolver
+
+Description
+
+SourceFiles
+    smoothAlignmentSolverI.H
+    smoothAlignmentSolver.C
+    smoothAlignmentSolverIO.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef smoothAlignmentSolver_H
+#define smoothAlignmentSolver_H
+
+#include "cellShapeControlMesh.H"
+#include "triadField.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class smoothAlignmentSolver Declaration
+\*---------------------------------------------------------------------------*/
+
+class smoothAlignmentSolver
+{
+    // Private data
+
+        cellShapeControlMesh& mesh_;
+
+
+    // Private Member Functions
+
+        template <class Triangulation, class Type>
+        tmp<Field<Type> > filterFarPoints
+        (
+            const Triangulation& mesh,
+            const Field<Type>& field
+        );
+
+        template <class Triangulation>
+        autoPtr<mapDistribute> buildMap
+        (
+            const Triangulation& mesh,
+            labelListList& pointPoints
+        );
+
+        template <class Triangulation>
+        autoPtr<mapDistribute> buildReferredMap
+        (
+            const Triangulation& mesh,
+            labelList& indices
+        );
+
+        template <class Triangulation>
+        tmp<triadField> buildAlignmentField(const Triangulation& mesh);
+
+        template <class Triangulation>
+        tmp<pointField> buildPointField(const Triangulation& mesh);
+
+        //- Disallow default bitwise copy construct
+        smoothAlignmentSolver(const smoothAlignmentSolver&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const smoothAlignmentSolver&);
+
+
+public:
+
+    // Constructors
+
+        //- Construct null
+        smoothAlignmentSolver(cellShapeControlMesh& mesh);
+
+
+    //- Destructor
+    ~smoothAlignmentSolver();
+
+
+    // Member Functions
+
+        // Access
+
+        // Check
+
+        // Edit
+
+            //- Smooth the alignments on the mesh
+            void smoothAlignments(const label maxSmoothingIterations);
+
+
+        // Write
+
+
+    // Member Operators
+
+    // Friend Functions
+
+    // Friend Operators
+
+    // IOstream Operators
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
index ec6ca47901..165d0b41e8 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
@@ -31,15 +31,15 @@ License
 #include "vectorTools.H"
 #include "IOmanip.H"
 #include "indexedCellChecks.H"
+#include "controlMeshRefinement.H"
+#include "smoothAlignmentSolver.H"
 
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
 namespace Foam
 {
-
 defineTypeNameAndDebug(conformalVoronoiMesh, 0);
-
 }
 
 const Foam::scalar Foam::conformalVoronoiMesh::tolParallel = 1e-3;
@@ -634,160 +634,6 @@ void Foam::conformalVoronoiMesh::insertInitialPoints()
 }
 
 
-Foam::scalar Foam::conformalVoronoiMesh::calculateLoadUnbalance() const
-{
-    label nRealVertices = 0;
-
-    for
-    (
-        Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
-        vit != finite_vertices_end();
-        ++vit
-    )
-    {
-        // Only store real vertices that are not feature vertices
-        if (vit->real() && !vit->featurePoint())
-        {
-            nRealVertices++;
-        }
-    }
-
-    scalar globalNRealVertices = returnReduce
-    (
-        nRealVertices,
-        sumOp<label>()
-    );
-
-    scalar unbalance = returnReduce
-    (
-        mag(1.0 - nRealVertices/(globalNRealVertices/Pstream::nProcs())),
-        maxOp<scalar>()
-    );
-
-    Info<< "    Processor unbalance " << unbalance << endl;
-
-    return unbalance;
-}
-
-
-bool Foam::conformalVoronoiMesh::distributeBackground()
-{
-    if (!Pstream::parRun())
-    {
-        return false;
-    }
-
-    Info<< nl << "Redistributing points" << endl;
-
-    timeCheck("Before distribute");
-
-    label iteration = 0;
-
-    scalar previousLoadUnbalance = 0;
-
-    while (true)
-    {
-        scalar maxLoadUnbalance = calculateLoadUnbalance();
-
-        if
-        (
-            maxLoadUnbalance <= cvMeshControls().maxLoadUnbalance()
-         || maxLoadUnbalance <= previousLoadUnbalance
-        )
-        {
-            // If this is the first iteration, return false, if it was a
-            // subsequent one, return true;
-            return iteration != 0;
-        }
-
-        previousLoadUnbalance = maxLoadUnbalance;
-
-        Info<< "    Total number of vertices before redistribution "
-            << returnReduce(label(number_of_vertices()), sumOp<label>())
-            << endl;
-
-        const fvMesh& bMesh = decomposition_().mesh();
-
-        volScalarField cellWeights
-        (
-            IOobject
-            (
-                "cellWeights",
-                bMesh.time().timeName(),
-                bMesh,
-                IOobject::NO_READ,
-                IOobject::NO_WRITE
-            ),
-            bMesh,
-            dimensionedScalar("weight", dimless, 1e-2),
-            zeroGradientFvPatchScalarField::typeName
-        );
-
-        meshSearch cellSearch(bMesh, polyMesh::FACEPLANES);
-
-        labelList cellVertices(bMesh.nCells(), 0);
-
-        for
-        (
-            Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
-            vit != finite_vertices_end();
-            ++vit
-        )
-        {
-            // Only store real vertices that are not feature vertices
-            if (vit->real() && !vit->featurePoint())
-            {
-                pointFromPoint v = topoint(vit->point());
-
-                label cellI = cellSearch.findCell(v);
-
-                if (cellI == -1)
-                {
-//                     Pout<< "findCell conformalVoronoiMesh::distribute "
-//                         << "findCell "
-//                         << vit->type() << " "
-//                         << vit->index() << " "
-//                         << v << " "
-//                         << cellI
-//                         << " find nearest cellI ";
-
-                    cellI = cellSearch.findNearestCell(v);
-                }
-
-                cellVertices[cellI]++;
-            }
-        }
-
-        forAll(cellVertices, cI)
-        {
-            // Give a small but finite weight for empty cells.  Some
-            // decomposition methods have difficulty with integer overflows in
-            // the sum of the normalised weight field.
-            cellWeights.internalField()[cI] = max
-            (
-                cellVertices[cI],
-                1e-2
-            );
-        }
-
-        autoPtr<mapDistributePolyMesh> mapDist = decomposition_().distribute
-        (
-            cellWeights
-        );
-
-        cellShapeControl_.shapeControlMesh().distribute(decomposition_);
-
-        distribute();
-
-        timeCheck("After distribute");
-
-        iteration++;
-    }
-
-    return true;
-}
-
-
 void Foam::conformalVoronoiMesh::distribute()
 {
     if (!Pstream::parRun())
@@ -871,13 +717,24 @@ void Foam::conformalVoronoiMesh::distribute()
 
 void Foam::conformalVoronoiMesh::buildCellSizeAndAlignmentMesh()
 {
-    cellShapeControl_.initialMeshPopulation(decomposition_);
+    controlMeshRefinement meshRefinement
+    (
+        cellShapeControl_
+    );
+
+    smoothAlignmentSolver meshAlignmentSmoother
+    (
+        cellShapeControl_.shapeControlMesh()
+    );
+
+    meshRefinement.initialMeshPopulation(decomposition_);
 
     cellShapeControlMesh& cellSizeMesh = cellShapeControl_.shapeControlMesh();
 
     if (Pstream::parRun())
     {
-        cellSizeMesh.distribute(decomposition_);
+        distributeBackground(cellSizeMesh);
+//        cellSizeMesh.distribute(decomposition_);
     }
 
     label nMaxIter = readLabel
@@ -892,40 +749,48 @@ void Foam::conformalVoronoiMesh::buildCellSizeAndAlignmentMesh()
 
     for (label i = 0; i < nMaxIter; ++i)
     {
-//        label nRemoved = cellSizeMesh.removePoints();
-        label nRemoved = 0;
-        reduce(nRemoved, sumOp<label>());
-
-        label nAdded = cellShapeControl_.refineMesh(decomposition_);
-//        label nAdded = 0;
+        label nAdded = meshRefinement.refineMesh(decomposition_);
+        //cellShapeControl_.refineMesh(decomposition_);
         reduce(nAdded, sumOp<label>());
 
         if (Pstream::parRun())
         {
-            cellSizeMesh.distribute(decomposition_);
-        }
-
-        if (nRemoved + nAdded == 0)
-        {
-            break;
+            distributeBackground(cellSizeMesh);
         }
 
         Info<< "    Iteration " << i
             << " Added = " << nAdded << " points"
-            << ", Removed = " << nRemoved << " points"
             << endl;
+
+        if (nAdded == 0)
+        {
+            break;
+        }
     }
 
-    cellShapeControl_.smoothMesh();
+    label maxSmoothingIterations = readLabel
+    (
+        cvMeshControls().cvMeshDict().subDict("motionControl").lookup
+        (
+            "maxSmoothingIterations"
+        )
+    );
+    meshAlignmentSmoother.smoothAlignments(maxSmoothingIterations);
 
     Info<< "Background cell size and alignment mesh:" << endl;
     cellSizeMesh.printInfo(Info);
 
-    if (cvMeshControls().objOutput())
+    Info<< "Triangulation is "
+        << (cellSizeMesh.is_valid() ? "valid" : "not valid!" )
+        << endl;
+
+    if (!Pstream::parRun())
     {
-        cellSizeMesh.writeTriangulation();
+        //cellSizeMesh.writeTriangulation();
         cellSizeMesh.write();
     }
+
+    cellSizeMesh.printVertexInfo(Info);
 }
 
 
@@ -1102,11 +967,6 @@ void Foam::conformalVoronoiMesh::setVertexSizeAndAlignment()
                 vit->targetCellSize(),
                 vit->alignment()
             );
-
-            //vit->alignment() = tensor(1,0,0,0,1,0,0,0,1);
-            //vit->alignment() = requiredAlignment(pt);
-
-            //vit->targetCellSize() = cellShapeControls().cellSize(pt);
         }
     }
 }
@@ -1424,13 +1284,13 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
     // Improve the guess that the backgroundMeshDecomposition makes with the
     // initial positions.  Use before building the surface conformation to
     // better balance the surface conformation load.
-    distributeBackground();
+    distributeBackground(*this);
 
     buildSurfaceConformation();
 
     // The introduction of the surface conformation may have distorted the
     // balance of vertices, distribute if necessary.
-    distributeBackground();
+    distributeBackground(*this);
 
     if (Pstream::parRun())
     {
@@ -1701,22 +1561,22 @@ void Foam::conformalVoronoiMesh::move()
                         delta
                     );
 
-                    if (targetFaceArea == 0)
-                    {
-                        Pout<< vA->info() << vB->info();
-
-                        Cell_handle ch = locate(vA->point());
-                        if (is_infinite(ch))
-                        {
-                            Pout<< "vA " << vA->targetCellSize() << endl;
-                        }
-
-                        ch = locate(vB->point());
-                        if (is_infinite(ch))
-                        {
-                            Pout<< "vB " << vB->targetCellSize() << endl;
-                        }
-                    }
+//                    if (targetFaceArea == 0)
+//                    {
+//                        Pout<< vA->info() << vB->info();
+//
+//                        Cell_handle ch = locate(vA->point());
+//                        if (is_infinite(ch))
+//                        {
+//                            Pout<< "vA " << vA->targetCellSize() << endl;
+//                        }
+//
+//                        ch = locate(vB->point());
+//                        if (is_infinite(ch))
+//                        {
+//                            Pout<< "vB " << vB->targetCellSize() << endl;
+//                        }
+//                    }
 
                     delta *= faceAreaWeightModel_->faceAreaWeight
                     (
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
index 912a426dd1..e0ab96a35a 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
@@ -34,6 +34,7 @@ SourceFiles
     conformalVoronoiMeshFeaturePoints.C
     conformalVoronoiMeshFeaturePointSpecialisations.C
     conformalVoronoiMeshCalcDualMesh.C
+    conformalVoronoiMeshTemplates.C
 
 \*---------------------------------------------------------------------------*/
 
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
index 14c5428c81..e277405937 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
@@ -69,7 +69,7 @@ void Foam::conformalVoronoiMesh::conformToSurface()
         // Rebuild, insert and store new surface conformation
         buildSurfaceConformation();
 
-        if (distributeBackground())
+        if (distributeBackground(*this))
         {
             if (Pstream::parRun())
             {
diff --git a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C
index 028dd5a0df..7abf6f23b7 100644
--- a/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C
+++ b/applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C
@@ -613,17 +613,9 @@ void Foam::conformalVoronoiMesh::insertFeaturePoints()
     }
 
     label nFeatureVertices = number_of_vertices() - preFeaturePointSize;
+    reduce(nFeatureVertices, sumOp<label>());
 
-    if (Pstream::parRun())
-    {
-        reduce(nFeatureVertices, sumOp<label>());
-    }
-
-    if (nFeatureVertices > 0)
-    {
-        Info<< "    Inserted " << nFeatureVertices
-            << " feature vertices" << endl;
-    }
+    Info<< "    Inserted " << nFeatureVertices << " feature vertices" << endl;
 
     featureVertices_.clear();
     featureVertices_.setSize(pts.size());
diff --git a/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtractDict b/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtractDict
index a52ebdfa55..0137468804 100644
--- a/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtractDict
+++ b/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtractDict
@@ -92,6 +92,7 @@ surface2.nas
 
     // Output the proximity of feature points and edges to each other
     featureProximity        no;
+
     // The maximum search distance to use when looking for other feature
     // points and edges
     maxFeatureProximity     1;
diff --git a/etc/bashrc b/etc/bashrc
index f077631ac8..4fb5d511df 100644
--- a/etc/bashrc
+++ b/etc/bashrc
@@ -2,7 +2,7 @@
 # =========                 |
 # \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
 #  \\    /   O peration     |
-#   \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+#   \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
 #    \\/     M anipulation  |
 #------------------------------------------------------------------------------
 # License
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchPointAddressing.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchPointAddressing.C
index 36f78d8b1e..b8f78906f3 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchPointAddressing.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchPointAddressing.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
diff --git a/src/OpenFOAM/primitives/triad/triad.C b/src/OpenFOAM/primitives/triad/triad.C
index e5bc4dcfff..a5c0d38b42 100644
--- a/src/OpenFOAM/primitives/triad/triad.C
+++ b/src/OpenFOAM/primitives/triad/triad.C
@@ -75,8 +75,6 @@ const triad triad::unset
 
 }
 
-}
-
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //