From 68167d23734b29ac71bea32b9336c3c1037c9aa5 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 19 Jul 2013 11:24:50 +0100
Subject: [PATCH 01/11] ENH: searchableSurface: added boundingSphere member
 function

---
 .../searchableSurface/searchableBox.C         | 35 ++++++++++
 .../searchableSurface/searchableBox.H         | 10 ++-
 .../searchableSurface/searchableCylinder.C    | 17 +++++
 .../searchableSurface/searchableCylinder.H    | 10 ++-
 .../searchableSurface/searchablePlane.C       | 16 ++++-
 .../searchableSurface/searchablePlane.H       | 11 ++-
 .../searchableSurface/searchablePlate.C       | 67 ++++++++++++++++++-
 .../searchableSurface/searchablePlate.H       | 54 +++------------
 .../searchableSurface/searchableSphere.C      | 17 +++++
 .../searchableSurface/searchableSphere.H      | 12 +++-
 .../searchableSurface/searchableSurface.H     |  8 +++
 .../searchableSurfaceCollection.C             | 36 ++++++++++
 .../searchableSurfaceCollection.H             | 10 ++-
 .../searchableSurfaceWithGaps.H               | 11 +++
 .../searchableSurface/triSurfaceMesh.C        | 55 ++++++++++++---
 .../searchableSurface/triSurfaceMesh.H        |  8 +++
 16 files changed, 317 insertions(+), 60 deletions(-)

diff --git a/src/meshTools/searchableSurface/searchableBox.C b/src/meshTools/searchableSurface/searchableBox.C
index 7d639b99f8..5970e82575 100644
--- a/src/meshTools/searchableSurface/searchableBox.C
+++ b/src/meshTools/searchableSurface/searchableBox.C
@@ -249,6 +249,41 @@ Foam::tmp<Foam::pointField> Foam::searchableBox::coordinates() const
 }
 
 
+void Foam::searchableBox::boundingSpheres
+(
+    pointField& centres,
+    scalarField& radiusSqr
+) const
+{
+    centres.setSize(size());
+    radiusSqr.setSize(size());
+    radiusSqr = 0.0;
+
+    const pointField pts(treeBoundBox::points());
+    const faceList& fcs = treeBoundBox::faces;
+
+    forAll(fcs, i)
+    {
+        const face& f = fcs[i];
+
+        centres[i] = f.centre(pts);
+        forAll(f, fp)
+        {
+            const point& pt = pts[f[fp]];
+
+            radiusSqr[i] = Foam::max
+            (
+                radiusSqr[i],
+                Foam::magSqr(pt-centres[i])
+            );
+        }
+    }
+
+    // Add a bit to make sure all points are tested inside
+    radiusSqr += Foam::sqr(SMALL);
+}
+
+
 Foam::tmp<Foam::pointField> Foam::searchableBox::points() const
 {
     return treeBoundBox::points();
diff --git a/src/meshTools/searchableSurface/searchableBox.H b/src/meshTools/searchableSurface/searchableBox.H
index 40e5c5d711..d677c1b87d 100644
--- a/src/meshTools/searchableSurface/searchableBox.H
+++ b/src/meshTools/searchableSurface/searchableBox.H
@@ -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
@@ -129,6 +129,14 @@ public:
         //  Usually the element centres (should be of length size()).
         virtual tmp<pointField> coordinates() const;
 
+        //- Get bounding spheres (centre and radius squared), one per element.
+        //  Any point on element is guaranteed to be inside.
+        virtual void boundingSpheres
+        (
+            pointField& centres,
+            scalarField& radiusSqr
+        ) const;
+
         //- Get the points that define the surface.
         virtual tmp<pointField> points() const;
 
diff --git a/src/meshTools/searchableSurface/searchableCylinder.C b/src/meshTools/searchableSurface/searchableCylinder.C
index 70d674d521..7657e883c9 100644
--- a/src/meshTools/searchableSurface/searchableCylinder.C
+++ b/src/meshTools/searchableSurface/searchableCylinder.C
@@ -47,6 +47,23 @@ Foam::tmp<Foam::pointField> Foam::searchableCylinder::coordinates() const
 }
 
 
+void Foam::searchableCylinder::boundingSpheres
+(
+    pointField& centres,
+    scalarField& radiusSqr
+) const
+{
+    centres.setSize(1);
+    centres[0] = 0.5*(point1_ + point2_);
+
+    radiusSqr.setSize(1);
+    radiusSqr[0] = Foam::magSqr(point1_-centres[0]) + Foam::sqr(radius_);
+
+    // Add a bit to make sure all points are tested inside
+    radiusSqr += Foam::sqr(SMALL);
+}
+
+
 Foam::tmp<Foam::pointField> Foam::searchableCylinder::points() const
 {
     tmp<pointField> tPts(new pointField(2));
diff --git a/src/meshTools/searchableSurface/searchableCylinder.H b/src/meshTools/searchableSurface/searchableCylinder.H
index 0bb570e480..ccf6850dc1 100644
--- a/src/meshTools/searchableSurface/searchableCylinder.H
+++ b/src/meshTools/searchableSurface/searchableCylinder.H
@@ -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
@@ -154,6 +154,14 @@ public:
         //  Usually the element centres (should be of length size()).
         virtual tmp<pointField> coordinates() const;
 
+        //- Get bounding spheres (centre and radius squared), one per element.
+        //  Any point on element is guaranteed to be inside.
+        virtual void boundingSpheres
+        (
+            pointField& centres,
+            scalarField& radiusSqr
+        ) const;
+
         //- Get the points that define the surface.
         virtual tmp<pointField> points() const;
 
diff --git a/src/meshTools/searchableSurface/searchablePlane.C b/src/meshTools/searchableSurface/searchablePlane.C
index a76821d457..e6c9668b1f 100644
--- a/src/meshTools/searchableSurface/searchablePlane.C
+++ b/src/meshTools/searchableSurface/searchablePlane.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -134,6 +134,20 @@ const Foam::wordList& Foam::searchablePlane::regions() const
 }
 
 
+void Foam::searchablePlane::boundingSpheres
+(
+    pointField& centres,
+    scalarField& radiusSqr
+) const
+{
+    centres.setSize(1);
+    centres[0] = refPoint();
+
+    radiusSqr.setSize(1);
+    radiusSqr[0] = Foam::sqr(GREAT);
+}
+
+
 void Foam::searchablePlane::findNearest
 (
     const pointField& samples,
diff --git a/src/meshTools/searchableSurface/searchablePlane.H b/src/meshTools/searchableSurface/searchablePlane.H
index cc64844d0d..a79f93f902 100644
--- a/src/meshTools/searchableSurface/searchablePlane.H
+++ b/src/meshTools/searchableSurface/searchablePlane.H
@@ -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
@@ -130,6 +130,15 @@ public:
             return tCtrs;
         }
 
+        //- Get bounding spheres (centre and radius squared), one per element.
+        //  Any point on element is guaranteed to be inside.
+        //  Note: radius limited to sqr(GREAT)
+        virtual void boundingSpheres
+        (
+            pointField& centres,
+            scalarField& radiusSqr
+        ) const;
+
         //- Get the points that define the surface.
         virtual tmp<pointField> points() const
         {
diff --git a/src/meshTools/searchableSurface/searchablePlate.C b/src/meshTools/searchableSurface/searchablePlate.C
index 37404db051..0ed10bd675 100644
--- a/src/meshTools/searchableSurface/searchablePlate.C
+++ b/src/meshTools/searchableSurface/searchablePlate.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
@@ -274,6 +274,71 @@ const Foam::wordList& Foam::searchablePlate::regions() const
 }
 
 
+Foam::tmp<Foam::pointField> Foam::searchablePlate::coordinates() const
+{
+    return tmp<pointField>(new pointField(1, origin_ + 0.5*span_));
+}
+
+
+void Foam::searchablePlate::boundingSpheres
+(
+    pointField& centres,
+    scalarField& radiusSqr
+) const
+{
+    centres.setSize(1);
+    centres[0] = origin_ + 0.5*span_;
+
+    radiusSqr.setSize(1);
+    radiusSqr[0] = Foam::magSqr(0.5*span_);
+
+    // Add a bit to make sure all points are tested inside
+    radiusSqr += Foam::sqr(SMALL);
+}
+
+
+Foam::tmp<Foam::pointField> Foam::searchablePlate::points() const
+{
+    tmp<pointField> tPts(new pointField(4));
+    pointField& pts = tPts();
+
+    pts[0] = origin_;
+    pts[2] = origin_ + span_;
+
+    if (span_.x() < span_.y() && span_.x() < span_.z())
+    {
+        pts[1] = origin_ + point(0, span_.y(), 0);
+        pts[3] = origin_ + point(0, 0, span_.z());
+    }
+    else if (span_.y() < span_.z())
+    {
+        pts[1] = origin_ + point(span_.x(), 0, 0);
+        pts[3] = origin_ + point(0, 0, span_.z());
+    }
+    else
+    {
+        pts[1] = origin_ + point(span_.x(), 0, 0);
+        pts[3] = origin_ + point(0, span_.y(), 0);
+    }
+
+    return tPts;
+}
+
+
+bool Foam::searchablePlate::overlaps(const boundBox& bb) const
+{
+    return
+    (
+        (origin_.x() + span_.x()) >= bb.min().x()
+     && origin_.x() <= bb.max().x()
+     && (origin_.y() + span_.y()) >= bb.min().y()
+     && origin_.y() <= bb.max().y()
+     && (origin_.z() + span_.z()) >= bb.min().z()
+     && origin_.z() <= bb.max().z()
+    );
+}
+
+
 void Foam::searchablePlate::findNearest
 (
     const pointField& samples,
diff --git a/src/meshTools/searchableSurface/searchablePlate.H b/src/meshTools/searchableSurface/searchablePlate.H
index 86b0785077..af2d732ff6 100644
--- a/src/meshTools/searchableSurface/searchablePlate.H
+++ b/src/meshTools/searchableSurface/searchablePlate.H
@@ -145,53 +145,21 @@ public:
 
         //- Get representative set of element coordinates
         //  Usually the element centres (should be of length size()).
-        virtual tmp<pointField> coordinates() const
-        {
-            tmp<pointField> tCtrs(new pointField(1, origin_ + 0.5*span_));
-            return tCtrs;
-        }
+        virtual tmp<pointField> coordinates() const;
+
+        //- Get bounding spheres (centre and radius squared), one per element.
+        //  Any point on element is guaranteed to be inside.
+        virtual void boundingSpheres
+        (
+            pointField& centres,
+            scalarField& radiusSqr
+        ) const;
 
         //- Get the points that define the surface.
-        virtual tmp<pointField> points() const
-        {
-            tmp<pointField> tPts(new pointField(4));
-            pointField& pts = tPts();
-
-            pts[0] = origin_;
-            pts[2] = origin_ + span_;
-
-            if (span_.x() < span_.y() && span_.x() < span_.z())
-            {
-                pts[1] = origin_ + point(0, span_.y(), 0);
-                pts[3] = origin_ + point(0, 0, span_.z());
-            }
-            else if (span_.y() < span_.z())
-            {
-                pts[1] = origin_ + point(span_.x(), 0, 0);
-                pts[3] = origin_ + point(0, 0, span_.z());
-            }
-            else
-            {
-                pts[1] = origin_ + point(span_.x(), 0, 0);
-                pts[3] = origin_ + point(0, span_.y(), 0);
-            }
-
-            return tPts;
-        }
+        virtual tmp<pointField> points() const;
 
         //- Does any part of the surface overlap the supplied bound box?
-        virtual bool overlaps(const boundBox& bb) const
-        {
-            return
-            (
-                (origin_.x() + span_.x()) >= bb.min().x()
-             && origin_.x() <= bb.max().x()
-             && (origin_.y() + span_.y()) >= bb.min().y()
-             && origin_.y() <= bb.max().y()
-             && (origin_.z() + span_.z()) >= bb.min().z()
-             && origin_.z() <= bb.max().z()
-            );
-        }
+        virtual bool overlaps(const boundBox& bb) const;
 
 
         // Multiple point queries.
diff --git a/src/meshTools/searchableSurface/searchableSphere.C b/src/meshTools/searchableSurface/searchableSphere.C
index c367c6281b..e4eae4ab27 100644
--- a/src/meshTools/searchableSurface/searchableSphere.C
+++ b/src/meshTools/searchableSurface/searchableSphere.C
@@ -185,6 +185,23 @@ const Foam::wordList& Foam::searchableSphere::regions() const
 
 
 
+void Foam::searchableSphere::boundingSpheres
+(
+    pointField& centres,
+    scalarField& radiusSqr
+) const
+{
+    centres.setSize(1);
+    centres[0] = centre_;
+
+    radiusSqr.setSize(1);
+    radiusSqr[0] = Foam::sqr(radius_);
+
+    // Add a bit to make sure all points are tested inside
+    radiusSqr += Foam::sqr(SMALL);
+}
+
+
 void Foam::searchableSphere::findNearest
 (
     const pointField& samples,
diff --git a/src/meshTools/searchableSurface/searchableSphere.H b/src/meshTools/searchableSurface/searchableSphere.H
index 5d46f8ad4b..50900f9361 100644
--- a/src/meshTools/searchableSurface/searchableSphere.H
+++ b/src/meshTools/searchableSurface/searchableSphere.H
@@ -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
@@ -60,7 +60,7 @@ private:
         //- Centre point
         const point centre_;
 
-        //- Radius squared
+        //- Radius
         const scalar radius_;
 
         //- Names of regions
@@ -139,6 +139,14 @@ public:
             return tCtrs;
         }
 
+        //- Get bounding spheres (centre and radius squared), one per element.
+        //  Any point on element is guaranteed to be inside.
+        virtual void boundingSpheres
+        (
+            pointField& centres,
+            scalarField& radiusSqr
+        ) const;
+
         //- Get the points that define the surface.
         virtual tmp<pointField> points() const
         {
diff --git a/src/meshTools/searchableSurface/searchableSurface.H b/src/meshTools/searchableSurface/searchableSurface.H
index 90224ae9ce..6a2adcc946 100644
--- a/src/meshTools/searchableSurface/searchableSurface.H
+++ b/src/meshTools/searchableSurface/searchableSurface.H
@@ -190,6 +190,14 @@ public:
         //  Usually the element centres (should be of length size()).
         virtual tmp<pointField> coordinates() const = 0;
 
+        //- Get bounding spheres (centre and radius squared), one per element.
+        //  Any point on element is guaranteed to be inside.
+        virtual void boundingSpheres
+        (
+            pointField& centres,
+            scalarField& radiusSqr
+        ) const = 0;
+
         //- Get the points that define the surface.
         virtual tmp<pointField> points() const = 0;
 
diff --git a/src/meshTools/searchableSurface/searchableSurfaceCollection.C b/src/meshTools/searchableSurface/searchableSurfaceCollection.C
index 3a02066afa..3aad85f72d 100644
--- a/src/meshTools/searchableSurface/searchableSurfaceCollection.C
+++ b/src/meshTools/searchableSurface/searchableSurfaceCollection.C
@@ -354,6 +354,42 @@ Foam::searchableSurfaceCollection::coordinates() const
 }
 
 
+void Foam::searchableSurfaceCollection::boundingSpheres
+(
+    pointField& centres,
+    scalarField& radiusSqr
+) const
+{
+    centres.setSize(size());
+    radiusSqr.setSize(centres.size());
+
+    // Append individual coordinates
+    label coordI = 0;
+
+    forAll(subGeom_, surfI)
+    {
+        scalar maxScale = cmptMax(scale_[surfI]);
+
+        pointField subCentres;
+        scalarField subRadiusSqr;
+        subGeom_[surfI].boundingSpheres(subCentres, subRadiusSqr);
+
+        forAll(subCentres, i)
+        {
+            centres[coordI++] = transform_[surfI].globalPosition
+            (
+                cmptMultiply
+                (
+                    subCentres[i],
+                    scale_[surfI]
+                )
+            );
+            radiusSqr[coordI++] = maxScale*subRadiusSqr[i];
+        }
+    }
+}
+
+
 Foam::tmp<Foam::pointField>
 Foam::searchableSurfaceCollection::points() const
 {
diff --git a/src/meshTools/searchableSurface/searchableSurfaceCollection.H b/src/meshTools/searchableSurface/searchableSurfaceCollection.H
index 53e677eaa6..dd8a21664b 100644
--- a/src/meshTools/searchableSurface/searchableSurfaceCollection.H
+++ b/src/meshTools/searchableSurface/searchableSurfaceCollection.H
@@ -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
@@ -176,6 +176,14 @@ public:
         //  Usually the element centres (should be of length size()).
         virtual tmp<pointField> coordinates() const;
 
+        //- Get bounding spheres (centre and radius squared), one per element.
+        //  Any point on element is guaranteed to be inside.
+        virtual void boundingSpheres
+        (
+            pointField& centres,
+            scalarField& radiusSqr
+        ) const;
+
         //- Get the points that define the surface.
         virtual tmp<pointField> points() const;
 
diff --git a/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H b/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H
index 1344cda2d9..d07d331c8a 100644
--- a/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H
+++ b/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H
@@ -156,6 +156,17 @@ public:
             return surface().coordinates();
         }
 
+        //- Get bounding spheres (centre and radius squared), one per element.
+        //  Any point on element is guaranteed to be inside.
+        virtual void boundingSpheres
+        (
+            pointField& centres,
+            scalarField& radiusSqr
+        ) const
+        {
+            surface().boundingSpheres(centres, radiusSqr);
+        }
+
         //- Get the points that define the surface.
         virtual tmp<pointField> points() const
         {
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.C b/src/meshTools/searchableSurface/triSurfaceMesh.C
index 125450603b..7b6b8ad6ad 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.C
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.C
@@ -145,10 +145,12 @@ bool Foam::triSurfaceMesh::addFaceToEdge
 
 bool Foam::triSurfaceMesh::isSurfaceClosed() const
 {
+    const pointField& pts = triSurface::points();
+
     // Construct pointFaces. Let's hope surface has compact point
     // numbering ...
     labelListList pointFaces;
-    invertManyToMany(points()().size(), *this, pointFaces);
+    invertManyToMany(pts.size(), *this, pointFaces);
 
     // Loop over all faces surrounding point. Count edges emanating from point.
     // Every edge should be used by two faces exactly.
@@ -241,7 +243,9 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
     minQuality_(-1),
     surfaceClosed_(-1)
 {
-    bounds() = boundBox(points());
+    const pointField& pts = triSurface::points();
+
+    bounds() = boundBox(pts);
 }
 
 
@@ -287,7 +291,9 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io)
     minQuality_(-1),
     surfaceClosed_(-1)
 {
-    bounds() = boundBox(points());
+    const pointField& pts = triSurface::points();
+
+    bounds() = boundBox(pts);
 }
 
 
@@ -347,7 +353,9 @@ Foam::triSurfaceMesh::triSurfaceMesh
         triSurface::scalePoints(scaleFactor);
     }
 
-    bounds() = boundBox(points());
+    const pointField& pts = triSurface::points();
+
+    bounds() = boundBox(pts);
 
     // Have optional minimum quality for normal calculation
     if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0)
@@ -393,6 +401,34 @@ Foam::tmp<Foam::pointField> Foam::triSurfaceMesh::coordinates() const
 }
 
 
+void Foam::triSurfaceMesh::boundingSpheres
+(
+    pointField& centres,
+    scalarField& radiusSqr
+) const
+{
+    centres = coordinates();
+    radiusSqr.setSize(size());
+    radiusSqr = 0.0;
+
+    const pointField& pts = triSurface::points();
+
+    forAll(*this, faceI)
+    {
+        const labelledTri& f = triSurface::operator[](faceI);
+        const point& fc = centres[faceI];
+        forAll(f, fp)
+        {
+            const point& pt = pts[f[fp]];
+            radiusSqr[faceI] = max(radiusSqr[faceI], Foam::magSqr(fc-pt));
+        }
+    }
+
+    // Add a bit to make sure all points are tested inside
+    radiusSqr += Foam::sqr(SMALL);
+}
+
+
 Foam::tmp<Foam::pointField> Foam::triSurfaceMesh::points() const
 {
     return triSurface::points();
@@ -606,6 +642,7 @@ void Foam::triSurfaceMesh::getNormal
 ) const
 {
     const triSurface& s = static_cast<const triSurface&>(*this);
+    const pointField& pts = s.points();
 
     normal.setSize(info.size());
 
@@ -621,9 +658,9 @@ void Foam::triSurfaceMesh::getNormal
             if (info[i].hit())
             {
                 label faceI = info[i].index();
-                normal[i] = s[faceI].normal(points());
+                normal[i] = s[faceI].normal(pts);
 
-                scalar qual = s[faceI].tri(points()).quality();
+                scalar qual = s[faceI].tri(pts).quality();
 
                 if (qual < minQuality_)
                 {
@@ -633,11 +670,11 @@ void Foam::triSurfaceMesh::getNormal
                     forAll(fFaces, j)
                     {
                         label nbrI = fFaces[j];
-                        scalar nbrQual = s[nbrI].tri(points()).quality();
+                        scalar nbrQual = s[nbrI].tri(pts).quality();
                         if (nbrQual > qual)
                         {
                             qual = nbrQual;
-                            normal[i] = s[nbrI].normal(points());
+                            normal[i] = s[nbrI].normal(pts);
                         }
                     }
                 }
@@ -662,7 +699,7 @@ void Foam::triSurfaceMesh::getNormal
                 //normal[i] = faceNormals()[faceI];
 
                 //- Uncached
-                normal[i] = s[faceI].normal(points());
+                normal[i] = s[faceI].normal(pts);
                 normal[i] /= mag(normal[i]) + VSMALL;
             }
             else
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.H b/src/meshTools/searchableSurface/triSurfaceMesh.H
index fa7e9e511c..e6ecd39594 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.H
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.H
@@ -187,6 +187,14 @@ public:
             //  Usually the element centres (should be of length size()).
             virtual tmp<pointField> coordinates() const;
 
+            //- Get bounding spheres (centre and radius squared). Any point
+            //  on surface is guaranteed to be inside.
+            virtual void boundingSpheres
+            (
+                pointField& centres,
+                scalarField& radiusSqr
+            ) const;
+
             //- Get the points that define the surface.
             virtual tmp<pointField> points() const;
 

From 25804eaad915a11f7502c784739ef6e7f53a0989 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 19 Jul 2013 11:43:57 +0100
Subject: [PATCH 02/11] ENH: refinementSurface: search bounding sphere only

---
 .../refinementSurfaces/refinementSurfaces.C       | 15 ++++++---------
 1 file changed, 6 insertions(+), 9 deletions(-)

diff --git a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
index 98b16e9f6b..f0277b2990 100644
--- a/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
+++ b/src/mesh/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
@@ -605,8 +605,10 @@ void Foam::refinementSurfaces::setMinLevelFields
         // searchableBox (size=6)
         if (geom.regions().size() > 1 && geom.globalSize() > 10)
         {
-            // Representative local coordinates
-            const pointField ctrs(geom.coordinates());
+            // Representative local coordinates and bounding sphere
+            pointField ctrs;
+            scalarField radiusSqr;
+            geom.boundingSpheres(ctrs, radiusSqr);
 
             labelList minLevelField(ctrs.size(), -1);
             {
@@ -614,12 +616,7 @@ void Foam::refinementSurfaces::setMinLevelFields
                 // distributed surface where local indices differ from global
                 // ones (needed for getRegion call)
                 List<pointIndexHit> info;
-                geom.findNearest
-                (
-                    ctrs,
-                    scalarField(ctrs.size(), sqr(GREAT)),
-                    info
-                );
+                geom.findNearest(ctrs, radiusSqr, info);
 
                 // Get per element the region
                 labelList region;
@@ -628,7 +625,7 @@ void Foam::refinementSurfaces::setMinLevelFields
                 // From the region get the surface-wise refinement level
                 forAll(minLevelField, i)
                 {
-                    if (info[i].hit())
+                    if (info[i].hit()) //Note: should not be necessary
                     {
                         minLevelField[i] = minLevel(surfI, region[i]);
                     }

From 53c8a1266866f0e13c2a759ebe3b9f4aaed60c47 Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 15:48:33 +0100
Subject: [PATCH 03/11] ENH: Pair: Add extra comparison operators. Make member
 functions non-member

---
 src/OpenFOAM/primitives/Pair/Pair.H           | 93 ++++++++++++-------
 .../meshRefinement/meshRefinementBaffles.C    |  2 +-
 2 files changed, 60 insertions(+), 35 deletions(-)

diff --git a/src/OpenFOAM/primitives/Pair/Pair.H b/src/OpenFOAM/primitives/Pair/Pair.H
index 12c16fc6c3..6715bffd7e 100644
--- a/src/OpenFOAM/primitives/Pair/Pair.H
+++ b/src/OpenFOAM/primitives/Pair/Pair.H
@@ -108,12 +108,6 @@ public:
             return this->operator[](1);
         }
 
-        //- Return reverse pair
-        inline Pair<Type> reversePair() const
-        {
-            return Pair<Type>(second(), first());
-        }
-
         //- Return other
         inline const Type& other(const Type& a) const
         {
@@ -147,11 +141,11 @@ public:
         //  - -1: same pair, but reversed order
         static inline int compare(const Pair<Type>& a, const Pair<Type>& b)
         {
-            if (a[0] == b[0] && a[1] == b[1])
+            if (a == b)
             {
                 return 1;
             }
-            else if (a[0] == b[1] && a[1] == b[0])
+            else if (a == reverse(b))
             {
                 return -1;
             }
@@ -160,35 +154,66 @@ public:
                 return 0;
             }
         }
-
-
-    // Friend Operators
-
-        friend bool operator==(const Pair<Type>& a, const Pair<Type>& b)
-        {
-            return (a.first() == b.first() && a.second() == b.second());
-        }
-
-        friend bool operator!=(const Pair<Type>& a, const Pair<Type>& b)
-        {
-            return !(a == b);
-        }
-
-        friend bool operator<(const Pair<Type>& a, const Pair<Type>& b)
-        {
-            return
-            (
-                a.first() < b.first()
-             ||
-                (
-                    !(b.first() < a.first())
-                 && a.second() < b.second()
-                )
-            );
-        }
 };
 
 
+template<class Type>
+Pair<Type> reverse(const Pair<Type>& p)
+{
+    return Pair<Type>(p.second(), p.first());
+}
+
+
+template<class Type>
+bool operator==(const Pair<Type>& a, const Pair<Type>& b)
+{
+    return (a.first() == b.first() && a.second() == b.second());
+}
+
+
+template<class Type>
+bool operator!=(const Pair<Type>& a, const Pair<Type>& b)
+{
+    return !(a == b);
+}
+
+
+template<class Type>
+bool operator<(const Pair<Type>& a, const Pair<Type>& b)
+{
+    return
+    (
+        a.first() < b.first()
+     ||
+        (
+            !(b.first() < a.first())
+         && a.second() < b.second()
+        )
+    );
+}
+
+
+template<class Type>
+bool operator<=(const Pair<Type>& a, const Pair<Type>& b)
+{
+    return !(b < a);
+}
+
+
+template<class Type>
+bool operator>(const Pair<Type>& a, const Pair<Type>& b)
+{
+    return (b < a);
+}
+
+
+template<class Type>
+bool operator>=(const Pair<Type>& a, const Pair<Type>& b)
+{
+    return !(a < b);
+}
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam
diff --git a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C
index 6b9d1074c5..16528d3ea7 100644
--- a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C
+++ b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementBaffles.C
@@ -450,7 +450,7 @@ Foam::Map<Foam::labelPair>  Foam::meshRefinement::getZoneBafflePatches
                     labelPair patches = zPatches;
                     if (fZone.flipMap()[i])
                     {
-                       patches = patches.reversePair();
+                       patches = reverse(patches);
                     }
 
                     if (!bafflePatch.insert(faceI, patches))

From 686b8f19d1750d8144f765036662689a9baa5167 Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 15:50:48 +0100
Subject: [PATCH 04/11] ENH: Move longestEdge function to face.H as a global
 function

---
 src/OpenFOAM/meshes/meshShapes/face/face.C    | 24 ++++++++++++++++-
 src/OpenFOAM/meshes/meshShapes/face/face.H    |  8 +++++-
 .../polyTopoChange/edgeCollapser.C            | 26 -------------------
 .../polyTopoChange/edgeCollapser.H            |  3 ---
 4 files changed, 30 insertions(+), 31 deletions(-)

diff --git a/src/OpenFOAM/meshes/meshShapes/face/face.C b/src/OpenFOAM/meshes/meshShapes/face/face.C
index 21576535e9..41650aad6c 100644
--- a/src/OpenFOAM/meshes/meshShapes/face/face.C
+++ b/src/OpenFOAM/meshes/meshShapes/face/face.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
@@ -829,4 +829,26 @@ Foam::label Foam::face::trianglesQuads
 }
 
 
+Foam::label Foam::longestEdge(const face& f, const pointField& pts)
+{
+    const edgeList& eds = f.edges();
+
+    label longestEdgeI = -1;
+    scalar longestEdgeLength = -SMALL;
+
+    forAll(eds, edI)
+    {
+        scalar edgeLength = eds[edI].mag(pts);
+
+        if (edgeLength > longestEdgeLength)
+        {
+            longestEdgeI = edI;
+            longestEdgeLength = edgeLength;
+        }
+    }
+
+    return longestEdgeI;
+}
+
+
 // ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/meshShapes/face/face.H b/src/OpenFOAM/meshes/meshShapes/face/face.H
index c481258a3c..adedf991bd 100644
--- a/src/OpenFOAM/meshes/meshShapes/face/face.H
+++ b/src/OpenFOAM/meshes/meshShapes/face/face.H
@@ -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
@@ -408,6 +408,12 @@ public:
 };
 
 
+// Global functions
+
+//- Find the longest edge on a face. Face point labels index into pts.
+label longestEdge(const face& f, const pointField& pts);
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C
index 70c2e7c558..6ea1e748b0 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C
@@ -42,32 +42,6 @@ defineTypeNameAndDebug(edgeCollapser, 0);
 }
 
 
-Foam::label Foam::edgeCollapser::longestEdge
-(
-    const face& f,
-    const pointField& pts
-)
-{
-    const edgeList& eds = f.edges();
-
-    label longestEdgeI = -1;
-    scalar longestEdgeLength = -SMALL;
-
-    forAll(eds, edI)
-    {
-        scalar edgeLength = eds[edI].mag(pts);
-
-        if (edgeLength > longestEdgeLength)
-        {
-            longestEdgeI = edI;
-            longestEdgeLength = edgeLength;
-        }
-    }
-
-    return longestEdgeI;
-}
-
-
 Foam::HashSet<Foam::label> Foam::edgeCollapser::checkBadFaces
 (
     const polyMesh& mesh,
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.H b/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.H
index 2393458cf7..dab8bd298b 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.H
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.H
@@ -262,9 +262,6 @@ public:
 
         // Check
 
-            //- Find the longest edge in a face
-            static label longestEdge(const face& f, const pointField& pts);
-
             //- Calls motionSmoother::checkMesh and returns a set of bad faces
             static HashSet<label> checkBadFaces
             (

From 0efac64b9199662e31b913789408306848c2830c Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 15:51:46 +0100
Subject: [PATCH 05/11] BUG: pointZone: Check point in same (or no) zone on all
 processors.

---
 .../meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C    |  4 ++--
 .../meshes/polyMesh/zones/pointZone/pointZone.C  | 16 ++++++++++++----
 2 files changed, 14 insertions(+), 6 deletions(-)

diff --git a/src/OpenFOAM/meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C b/src/OpenFOAM/meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C
index 33e17117c4..eeec0974d6 100644
--- a/src/OpenFOAM/meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C
+++ b/src/OpenFOAM/meshes/polyMesh/zones/ZoneMesh/ZoneMesh.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -492,7 +492,7 @@ bool Foam::ZoneMesh<ZoneType, MeshType>::checkParallelSync
                         << " of type " << zones[zoneI].type()
                         << " is not correctly synchronised"
                         << " across coupled boundaries."
-                        << " (coupled faces are either not both "
+                        << " (coupled faces are either not both"
                         << " present in set or have same flipmap)" << endl;
                 }
             }
diff --git a/src/OpenFOAM/meshes/polyMesh/zones/pointZone/pointZone.C b/src/OpenFOAM/meshes/polyMesh/zones/pointZone/pointZone.C
index ff15dcea3b..30d1161bb7 100644
--- a/src/OpenFOAM/meshes/polyMesh/zones/pointZone/pointZone.C
+++ b/src/OpenFOAM/meshes/polyMesh/zones/pointZone/pointZone.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -154,8 +154,15 @@ bool Foam::pointZone::checkParallelSync(const bool report) const
 
     forAll(maxZone, pointI)
     {
-        // Check point in zone on both sides
-        if (maxZone[pointI] != minZone[pointI])
+        // Check point in same (or no) zone on all processors
+        if
+        (
+            (
+                maxZone[pointI] != -1
+             || minZone[pointI] != labelMax
+            )
+         && (maxZone[pointI] != minZone[pointI])
+        )
         {
             if (report && !error)
             {
@@ -167,7 +174,8 @@ bool Foam::pointZone::checkParallelSync(const bool report) const
                     << (minZone[pointI] == labelMax ? -1 : minZone[pointI])
                     << " on some processors and in zone "
                     << maxZone[pointI]
-                    << " on some other processors."
+                    << " on some other processors." << nl
+                    << "(suppressing further warnings)"
                     << endl;
             }
             error = true;

From 97bac47725fbed54465a588ad16b4439ac7b5a20 Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 15:53:08 +0100
Subject: [PATCH 06/11] ENH: Tuple2: Remove comparison operators as friends and
 make non-members.

---
 src/OpenFOAM/primitives/Tuple2/Tuple2.H | 45 +++++--------------------
 1 file changed, 9 insertions(+), 36 deletions(-)

diff --git a/src/OpenFOAM/primitives/Tuple2/Tuple2.H b/src/OpenFOAM/primitives/Tuple2/Tuple2.H
index b7aae7598c..26ca000575 100644
--- a/src/OpenFOAM/primitives/Tuple2/Tuple2.H
+++ b/src/OpenFOAM/primitives/Tuple2/Tuple2.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -47,20 +47,6 @@ namespace Foam
 template<class Type1, class Type2>
 class Tuple2;
 
-template<class Type1, class Type2>
-inline bool operator==
-(
-    const Tuple2<Type1, Type2>&,
-    const Tuple2<Type1, Type2>&
-);
-
-template<class Type1, class Type2>
-inline bool operator!=
-(
-    const Tuple2<Type1, Type2>&,
-    const Tuple2<Type1, Type2>&
-);
-
 template<class Type1, class Type2>
 inline Istream& operator>>(Istream&, Tuple2<Type1, Type2>&);
 
@@ -129,27 +115,6 @@ public:
             return s_;
         }
 
-        //- Return reverse pair
-        inline Tuple2<Type2, Type1> reverseTuple2() const
-        {
-            return Tuple2<Type2, Type1>(second(), first());
-        }
-
-
-    // Friend Operators
-
-        friend bool operator== <Type1, Type2>
-        (
-            const Tuple2<Type1, Type2>& a,
-            const Tuple2<Type1, Type2>& b
-        );
-
-        friend bool operator!= <Type1, Type2>
-        (
-            const Tuple2<Type1, Type2>& a,
-            const Tuple2<Type1, Type2>& b
-        );
-
 
     // IOstream operators
 
@@ -169,6 +134,14 @@ public:
 };
 
 
+//- Return reverse of a tuple2
+template<class Type1, class Type2>
+inline Tuple2<Type2, Type1> reverse(const Tuple2<Type1, Type2>& t)
+{
+    return Tuple2<Type2, Type1>(t.second(), t.first());
+}
+
+
 template<class Type1, class Type2>
 inline bool operator==
 (

From bc6ed9b70abba4d9ec43836a8bad16ddd0330211 Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 15:53:59 +0100
Subject: [PATCH 07/11] ENH: triad: Do not sort unless all components are set.

---
 src/OpenFOAM/primitives/triad/triad.C | 5 +++++
 1 file changed, 5 insertions(+)

diff --git a/src/OpenFOAM/primitives/triad/triad.C b/src/OpenFOAM/primitives/triad/triad.C
index ac997d8bd2..a642463f42 100644
--- a/src/OpenFOAM/primitives/triad/triad.C
+++ b/src/OpenFOAM/primitives/triad/triad.C
@@ -282,6 +282,11 @@ void Foam::triad::align(const vector& v)
 
 Foam::triad Foam::triad::sortxyz() const
 {
+    if (!this->set())
+    {
+        return *this;
+    }
+
     triad t;
 
     if

From 1c55166d36b56643e79f11da4f55692530de89fd Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 16:42:16 +0100
Subject: [PATCH 08/11] BUG: DelaunayMesh: Do not allow vertexCount to be
 externally updated

---
 .../conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C        | 2 +-
 .../conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H        | 1 -
 .../conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H       | 6 ------
 3 files changed, 1 insertion(+), 8 deletions(-)

diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C
index d78011c9db..7eb8216b00 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C
@@ -122,7 +122,7 @@ Foam::DelaunayMesh<Triangulation>::DelaunayMesh
             if (indices.headerOk())
             {
                 vh->index() = indices[ptI];
-                vertexCount()++;
+                vertexCount_++;
             }
             else
             {
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H
index 1fbf19c1d6..63d227fcca 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H
@@ -201,7 +201,6 @@ public:
         inline void resetCellCount();
 
         inline label vertexCount() const;
-        inline label& vertexCount();
 
         inline void resetVertexCount();
 
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H
index bf2277e356..bad25b2b83 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMeshI.H
@@ -124,12 +124,6 @@ Foam::label Foam::DelaunayMesh<Triangulation>::vertexCount() const
     return vertexCount_;
 }
 
-template<class Triangulation>
-Foam::label& Foam::DelaunayMesh<Triangulation>::vertexCount()
-{
-    return vertexCount_;
-}
-
 
 template<class Triangulation>
 void Foam::DelaunayMesh<Triangulation>::resetVertexCount()

From f008250b03686f2e2826ba33235ae45129998efd Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 16:43:07 +0100
Subject: [PATCH 09/11] ENH: foamyHexMesh: Normalise the aspect ratio direction

---
 .../cellAspectRatioControl.C                  | 47 ++++++++++---------
 .../cellAspectRatioControl.H                  |  7 +--
 2 files changed, 26 insertions(+), 28 deletions(-)

diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C
index 37125bb76b..a033e5076d 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.C
@@ -40,11 +40,15 @@ Foam::cellAspectRatioControl::cellAspectRatioControl
         aspectRatioDict_.lookupOrDefault<vector>
         (
             "aspectRatioDirection",
-            vector(0, 0, 0)
+            vector::zero
         )
     )
 {
-    Info<< nl << "Cell Aspect Ratio Control" << nl
+    // Normalise the direction
+    aspectRatioDirection_ /= mag(aspectRatioDirection_) + SMALL;
+
+    Info<< nl
+        << "Cell Aspect Ratio Control" << nl
         << "    Ratio     : " << aspectRatio_ << nl
         << "    Direction : " << aspectRatioDirection_
         << endl;
@@ -66,22 +70,20 @@ void Foam::cellAspectRatioControl::updateCellSizeAndFaceArea
     scalar& targetCellSize
 ) const
 {
-    const scalar cosAngle = mag
-    (
-        vectorTools::cosPhi(alignmentDir, aspectRatioDirection_)
-    );
+    const scalar cosAngle =
+        mag(vectorTools::cosPhi(alignmentDir, aspectRatioDirection_));
 
     // Change target face area based on aspect ratio
-    targetFaceArea
-        += targetFaceArea
-          *(aspectRatio_ - 1.0)
-          *(1.0 - cosAngle);
+    targetFaceArea +=
+        targetFaceArea
+       *(aspectRatio_ - 1.0)
+       *(1.0 - cosAngle);
 
     // Change target cell size based on aspect ratio
-    targetCellSize
-        += targetCellSize
-          *(aspectRatio_ - 1.0)
-          *cosAngle;
+    targetCellSize +=
+        targetCellSize
+       *(aspectRatio_ - 1.0)
+       *cosAngle;
 
     alignmentDir *= 0.5*targetCellSize;
 }
@@ -95,16 +97,15 @@ void Foam::cellAspectRatioControl::updateDeltaVector
     vector& delta
 ) const
 {
-    const scalar cosAngle = mag
-    (
-        vectorTools::cosPhi(alignmentDir, aspectRatioDirection_)
-    );
+    const scalar cosAngle =
+        mag(vectorTools::cosPhi(alignmentDir, aspectRatioDirection_));
 
-    delta += 0.5
-            *delta
-            *cosAngle
-            *(targetCellSize/rABMag)
-            *(aspectRatio_ - 1.0);
+    delta +=
+        0.5
+       *delta
+       *cosAngle
+       *(targetCellSize/rABMag)
+       *(aspectRatio_ - 1.0);
 }
 
 
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.H
index cdc7510f23..b3964432b1 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellAspectRatioControl/cellAspectRatioControl.H
@@ -56,7 +56,7 @@ class cellAspectRatioControl
 
         const scalar aspectRatio_;
 
-        const vector aspectRatioDirection_;
+        vector aspectRatioDirection_;
 
 
     // Private Member Functions
@@ -73,10 +73,7 @@ public:
     // Constructors
 
         //- Construct from dictionary
-        cellAspectRatioControl
-        (
-            const dictionary& motionDict
-        );
+        cellAspectRatioControl(const dictionary& motionDict);
 
 
     //- Destructor

From 1a9329ed4883f41954778284361decae2c414c78 Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 16:46:24 +0100
Subject: [PATCH 10/11] STYLE: foamyHexMesh: Remove old code related mainly to
 the deprecated cell sizing methods

---
 .../conformalVoronoiMesh.C                    | 199 ------------------
 .../conformalVoronoiMesh.H                    |  70 ------
 .../conformalVoronoiMeshCalcDualMesh.C        |   6 -
 .../conformalVoronoiMeshConformToSurface.C    |  36 ----
 .../conformalVoronoiMeshI.H                   |  95 ---------
 .../indexedCell/indexedCellI.H                |   9 +-
 .../pointFeatureEdgesTypes.H                  |   3 +-
 7 files changed, 5 insertions(+), 413 deletions(-)

diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
index 62d0c02afe..3ba5b8c9e7 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
@@ -570,159 +570,8 @@ void Foam::conformalVoronoiMesh::buildCellSizeAndAlignmentMesh()
 }
 
 
-void Foam::conformalVoronoiMesh::storeSizesAndAlignments()
-{
-    DynamicList<Point> storePts(number_of_vertices());
-
-    for
-    (
-        Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
-        vit != finite_vertices_end();
-        vit++
-    )
-    {
-        if (vit->internalPoint())
-        {
-            storePts.append(vit->point());
-        }
-    }
-
-    storePts.shrink();
-
-    storeSizesAndAlignments(storePts);
-}
-
-
-void Foam::conformalVoronoiMesh::storeSizesAndAlignments
-(
-    const List<Point>& storePts
-)
-{
-//    timeCheck("Start of storeSizesAndAlignments");
-//
-//    Info << nl << "Store size and alignment" << endl;
-//
-//    sizeAndAlignmentLocations_.setSize(storePts.size());
-//
-//    storedSizes_.setSize(sizeAndAlignmentLocations_.size());
-//
-//    storedAlignments_.setSize(sizeAndAlignmentLocations_.size());
-//
-//    label i = 0;
-//
-//    //checkCellSizing();
-//
-//    for
-//    (
-//        List<Point>::const_iterator pit = storePts.begin();
-//        pit != storePts.end();
-//        ++pit
-//    )
-//    {
-//        pointFromPoint pt = topoint(*pit);
-//
-////        storedAlignments_[i] = requiredAlignment(pt);
-////
-////        storedSizes_[i] = cellShapeControls().cellSize(pt);
-//
-//        cellShapeControls().cellSizeAndAlignment
-//        (
-//            pt,
-//            storedSizes_[i],
-//            storedAlignments_[i]
-//        );
-//
-//        i++;
-//    }
-//
-//    timeCheck("Sizes and alignments calculated, build tree");
-//
-//    buildSizeAndAlignmentTree();
-//
-//    timeCheck("Size and alignment tree built");
-}
-
-
-void Foam::conformalVoronoiMesh::updateSizesAndAlignments
-(
-    const List<Point>& storePts
-)
-{
-    // This function is only used in serial, the background redistribution
-    // triggers this when unbalance is detected in parallel.
-
-    if
-    (
-        !Pstream::parRun()
-     && runTime_.run()
-     && runTime_.timeIndex()
-    )
-    {
-        storeSizesAndAlignments(storePts);
-
-        timeCheck("Updated sizes and alignments");
-    }
-}
-
-
-const Foam::indexedOctree<Foam::treeDataPoint>&
-Foam::conformalVoronoiMesh::sizeAndAlignmentTree() const
-{
-    if (sizeAndAlignmentTreePtr_.empty())
-    {
-        buildSizeAndAlignmentTree();
-    }
-
-    return sizeAndAlignmentTreePtr_();
-}
-
-
 void Foam::conformalVoronoiMesh::setVertexSizeAndAlignment()
 {
-//    Info<< nl << "Looking up target cell alignment and size" << endl;
-//
-//    const indexedOctree<treeDataPoint>& tree = sizeAndAlignmentTree();
-//
-//    for
-//    (
-//        Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
-//        vit != finite_vertices_end();
-//        vit++
-//    )
-//    {
-//        if
-//        (
-//            vit->internalOrBoundaryPoint()
-//         || vit->referredInternalOrBoundaryPoint()
-//        )
-//        {
-//            pointFromPoint pt = topoint(vit->point());
-//
-//            pointIndexHit info = tree.findNearest(pt, sqr(GREAT));
-//
-//            if (info.hit())
-//            {
-//                vit->alignment() = storedAlignments_[info.index()];
-//
-//                vit->targetCellSize() = storedSizes_[info.index()];
-//            }
-//            else
-//            {
-//                WarningIn
-//                (
-//                    "void "
-//                    "Foam::conformalVoronoiMesh::setVertexSizeAndAlignment()"
-//                )
-//                    << "Point " << pt << " did not find a nearest point "
-//                    << " for alignment and size lookup." << endl;
-//
-//                vit->alignment() = cellShapeControls().cellAlignment(pt);
-//
-//                vit->targetCellSize() = cellShapeControls().cellSize(pt);
-//            }
-//        }
-//    }
-
     Info<< nl << "Calculating target cell alignment and size" << endl;
 
     for
@@ -744,44 +593,6 @@ void Foam::conformalVoronoiMesh::setVertexSizeAndAlignment()
             );
         }
     }
-
-//    OFstream str(runTime_.path()/"alignments_internal.obj");
-//
-//    for
-//    (
-//        Finite_vertices_iterator vit = finite_vertices_begin();
-//        vit != finite_vertices_end();
-//        ++vit
-//    )
-//    {
-//        if (!vit->farPoint())
-//        {
-//            // Write alignments
-//            const tensor& alignment = vit->alignment();
-//            pointFromPoint pt = topoint(vit->point());
-//
-//            if
-//            (
-//                alignment.x() == triad::unset[0]
-//             || alignment.y() == triad::unset[0]
-//             || alignment.z() == triad::unset[0]
-//            )
-//            {
-//                Info<< "Bad alignment = " << vit->info();
-//
-//                vit->alignment() = tensor::I;
-//
-//                Info<< "New alignment = " << vit->info();
-//
-//                continue;
-//            }
-//
-//            meshTools::writeOBJ(str, pt, alignment.x() + pt);
-//            meshTools::writeOBJ(str, pt, alignment.y() + pt);
-//            meshTools::writeOBJ(str, pt, alignment.z() + pt);
-//        }
-//    }
-
 }
 
 
@@ -1039,10 +850,6 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
     featurePointLocations_(),
     edgeLocationTreePtr_(),
     surfacePtLocationTreePtr_(),
-    sizeAndAlignmentLocations_(),
-    storedSizes_(),
-    storedAlignments_(),
-    sizeAndAlignmentTreePtr_(),
     surfaceConformationVertices_(),
     initialPointsMethod_
     (
@@ -1133,10 +940,6 @@ void Foam::conformalVoronoiMesh::initialiseForMotion()
     // (potentially) redistributed.
     storeSurfaceConformation();
 
-    // Use storeSizesAndAlignments with no feed points because all background
-    // points may have been distributed.
-    storeSizesAndAlignments();
-
     // Report any Delaunay vertices that do not think that they are in the
     // domain the processor they are on.
     // reportProcessorOccupancy();
@@ -1836,8 +1639,6 @@ void Foam::conformalVoronoiMesh::move()
         writeMesh(time().timeName());
     }
 
-    updateSizesAndAlignments(pointsToInsert);
-
     Info<< nl
         << "Total displacement = " << totalDisp << nl
         << "Total distance = " << totalDist << nl
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
index ac11d1a8fa..777f277d58 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
@@ -177,19 +177,6 @@ private:
 
         mutable DynamicList<Foam::point> existingSurfacePtLocations_;
 
-        //- Store locations where the cell size and alignments will be
-        //  pre-calculated and looked up
-        pointField sizeAndAlignmentLocations_;
-
-        //- Stored cell size at sizeAndAlignmentLocations_
-        scalarField storedSizes_;
-
-        //- Stored alignments at sizeAndAlignmentLocations_
-        tensorField storedAlignments_;
-
-        //- Search tree for size and alignment lookup points
-        mutable autoPtr<indexedOctree<treeDataPoint> > sizeAndAlignmentTreePtr_;
-
         //- Store the surface and feature edge conformation locations to be
         //  reinserted
         List<Vb> surfaceConformationVertices_;
@@ -271,40 +258,6 @@ private:
         //- Return the local maximum surface protrusion distance
         inline scalar maxSurfaceProtrusion(const Foam::point& pt) const;
 
-        //- Insert Point and return its auto-generated index
-        inline bool insertPoint
-        (
-            const Point& P,
-            const indexedVertexEnum::vertexType type
-        );
-
-        //- Insert Foam::point with specified index and type
-        inline bool insertPoint
-        (
-            const Foam::point& p,
-            const indexedVertexEnum::vertexType type
-        );
-
-        //- Insert Point with specified index, type and original processor
-        inline bool insertReferredPoint
-        (
-            const Point& P,
-            const label index,
-            const indexedVertexEnum::vertexType type,
-            const label processor
-        );
-
-        inline bool insertReferredPoint(const Vb& P);
-
-        //- Insert Foam::point with specified index, type and original processor
-        inline bool insertReferredPoint
-        (
-            const Foam::point& p,
-            const label index,
-            const indexedVertexEnum::vertexType type,
-            const label processor
-        );
-
         //- Insert Delaunay vertices using the CGAL range insertion method,
         //  optionally check processor occupancy and distribute to other
         //  processors
@@ -539,22 +492,6 @@ private:
 
         void buildCellSizeAndAlignmentMesh();
 
-        //- Store data for sizeAndAlignmentLocations_, storedSizes_ and
-        //  storedAlignments_ and initialise the sizeAndAlignmentTreePtr_,
-        //  determining the appropriate sizeAndAlignmentLocations_
-        //  automatically
-        void storeSizesAndAlignments();
-
-        //- Store data for sizeAndAlignmentLocations_, storedSizes_ and
-        //  storedAlignments_ and initialise the sizeAndAlignmentTreePtr_
-        void storeSizesAndAlignments(const List<Point>& storePts);
-
-        //- Restore the sizes and alignments if required
-        void updateSizesAndAlignments(const List<Point>& storePts);
-
-        //- Demand driven construction of octree for and alignment points
-        const indexedOctree<treeDataPoint>& sizeAndAlignmentTree() const;
-
         //- Set the size and alignment data for each vertex
         void setVertexSizeAndAlignment();
 
@@ -750,9 +687,6 @@ private:
             const DynamicList<Foam::point>& existingSurfacePtLocations
         ) const;
 
-        //- Build or rebuild the sizeAndAlignmentTree
-        void buildSizeAndAlignmentTree() const;
-
         //- Process the surface conformation locations to decide which surface
         //  and edge conformation locations to add
         void addSurfaceAndEdgeHits
@@ -1063,10 +997,6 @@ public:
 //            const List<scalar>& radiusSqrs
 //        ) const;
 
-        typedef K::Vector_3 CGALVector;
-
-        inline CGALVector toCGALVector(const Foam::vector& v) const;
-
 
         // Access
 
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
index 7c07d12e67..5a45b8ab4f 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
@@ -211,11 +211,6 @@ void Foam::conformalVoronoiMesh::checkDuals()
 
     List<Point> duals(number_of_finite_cells());
 
-    typedef CGAL::Exact_predicates_exact_constructions_kernel       EK2;
-    typedef CGAL::Regular_triangulation_euclidean_traits_3<EK2>     EK;
-    typedef CGAL::Cartesian_converter<baseK::Kernel, EK2>  To_exact;
-    typedef CGAL::Cartesian_converter<EK2, baseK::Kernel>  Back_from_exact;
-
 //    PackedBoolList bPoints(number_of_finite_cells());
 
 //    indexDualVertices(duals, bPoints);
@@ -1787,7 +1782,6 @@ void Foam::conformalVoronoiMesh::indexDualVertices
     OBJstream snapping1("snapToSurface1.obj");
     OBJstream snapping2("snapToSurface2.obj");
     OFstream tetToSnapTo("tetsToSnapTo.obj");
-    label offset = 0;
 
     for
     (
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
index 8be496816f..43feda3654 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
@@ -76,10 +76,6 @@ void Foam::conformalVoronoiMesh::conformToSurface()
             {
                 sync(decomposition_().procBounds());
             }
-
-            // Use storeSizesAndAlignments with no feed points because all
-            // background points may have been distributed.
-            storeSizesAndAlignments();
         }
 
         // Do not store the surface conformation until after it has been
@@ -2055,38 +2051,6 @@ void Foam::conformalVoronoiMesh::buildSurfacePtLocationTree
 }
 
 
-void Foam::conformalVoronoiMesh::buildSizeAndAlignmentTree() const
-{
-    if (sizeAndAlignmentLocations_.empty())
-    {
-        FatalErrorIn("buildSizeAndAlignmentTree()")
-            << "sizeAndAlignmentLocations empty, must be populated before "
-            << "sizeAndAlignmentTree can be built."
-            << exit(FatalError);
-    }
-
-    treeBoundBox overallBb
-    (
-        geometryToConformTo_.globalBounds().extend(rndGen_, 1e-4)
-    );
-
-    overallBb.min() -= Foam::point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
-    overallBb.max() += Foam::point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
-
-    sizeAndAlignmentTreePtr_.reset
-    (
-        new indexedOctree<treeDataPoint>
-        (
-            treeDataPoint(sizeAndAlignmentLocations_),
-            overallBb,  // overall search domain
-            10,         // max levels
-            20.0,       // maximum ratio of cubes v.s. cells
-            100.0       // max. duplicity; n/a since no bounding boxes.
-        )
-    );
-}
-
-
 void Foam::conformalVoronoiMesh::addSurfaceAndEdgeHits
 (
     const Foam::point& vit,
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
index a912515fb3..9fa34b82d3 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
@@ -236,94 +236,6 @@ inline Foam::scalar Foam::conformalVoronoiMesh::maxSurfaceProtrusion
 }
 
 
-inline bool Foam::conformalVoronoiMesh::insertPoint
-(
-    const Foam::point& p,
-    const indexedVertexEnum::vertexType type
-)
-{
-    return insertPoint(toPoint<Point>(p), type);
-}
-
-
-inline bool Foam::conformalVoronoiMesh::insertPoint
-(
-    const Point& P,
-    const indexedVertexEnum::vertexType type
-)
-{
-    uint nVert = number_of_vertices();
-
-    Vertex_handle vh = insert(P);
-
-    bool pointInserted = true;
-
-    if (nVert == number_of_vertices())
-    {
-        Pout<< "Failed to insert point : " << topoint(P)
-            << " of type " << type << endl;
-        pointInserted = false;
-    }
-    else
-    {
-        vh->index() = getNewVertexIndex();
-        vh->type() = type;
-    }
-
-    return pointInserted;
-}
-
-
-inline bool Foam::conformalVoronoiMesh::insertReferredPoint(const Vb& P)
-{
-    return insertReferredPoint(P.point(), P.index(), P.type(), P.procIndex());
-}
-
-
-inline bool Foam::conformalVoronoiMesh::insertReferredPoint
-(
-    const Foam::point& p,
-    const label index,
-    const indexedVertexEnum::vertexType type,
-    const label processor
-)
-{
-    return insertReferredPoint(toPoint<Point>(p), index, type, processor);
-}
-
-
-inline bool Foam::conformalVoronoiMesh::insertReferredPoint
-(
-    const Point& P,
-    const label index,
-    const indexedVertexEnum::vertexType type,
-    const label processor
-)
-{
-    uint nVert = number_of_vertices();
-
-    Vertex_handle vh = insert(P);
-
-    bool pointInserted = true;
-
-    if (nVert == number_of_vertices())
-    {
-        Pout<< "Failed to insert point " << topoint(P)
-            << " type: " << type << " index: " << index
-            << " proc: " << processor << endl;
-        pointInserted = false;
-    }
-    else
-    {
-        vh->index() = index;
-        vh->type() = type;
-        vh->procIndex() = processor;
-    }
-
-    return pointInserted;
-}
-
-
 inline void Foam::conformalVoronoiMesh::createPointPair
 (
     const scalar ppDist,
@@ -663,13 +575,6 @@ inline bool Foam::conformalVoronoiMesh::isProcBoundaryEdge
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline Foam::conformalVoronoiMesh::CGALVector
-Foam::conformalVoronoiMesh::toCGALVector(const Foam::vector& v) const
-{
-    return CGALVector(v.x(), v.y(), v.z());
-}
-
-
 inline const Foam::Time& Foam::conformalVoronoiMesh::time() const
 {
     return runTime_;
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H
index 7875d67d39..a812fe4664 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H
@@ -108,17 +108,15 @@ int CGAL::indexedCell<Gt, Cb>::cellIndex() const
 
 
 #ifdef CGAL_INEXACT
+
     template<class Gt, class Cb>
     const Foam::point& CGAL::indexedCell<Gt, Cb>::dual()
     {
-    //    if (Foam::foamyHexMeshChecks::coplanarTet(*this, 1e-20) == 0)
-    //    {
-            // Do exact calc
-    //    }
-
         return reinterpret_cast<const Foam::point&>(this->circumcenter());
     }
+
 #else
+
     template<class Gt, class Cb>
     const Foam::point CGAL::indexedCell<Gt, Cb>::dual()
     {
@@ -131,6 +129,7 @@ int CGAL::indexedCell<Gt, Cb>::cellIndex() const
             CGAL::to_double(P.z())
         );
     }
+
 #endif
 
 
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/pointFeatureEdgesTypes.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/pointFeatureEdgesTypes.H
index 3dd3761c59..4277f0d3b7 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/pointFeatureEdgesTypes.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/pointFeatureEdgesTypes.H
@@ -25,8 +25,7 @@ Class
     pointFeatureEdgesTypes
 
 Description
-    struct for holding information on the types of feature edges attached to
-    feature points
+    Holds information on the types of feature edges attached to feature points.
 
 \*---------------------------------------------------------------------------*/
 

From ad2ef8e6fb9a784238bef4108b0fcff6d621a901 Mon Sep 17 00:00:00 2001
From: laurence <laurence>
Date: Fri, 19 Jul 2013 16:48:11 +0100
Subject: [PATCH 11/11] ENH: foamyHexMeshDict: Remove locationInMesh as a
 default and use default meshQualityDict

---
 etc/caseDicts/foamyHexMeshDict | 41 +++++++++++++++++++---------------
 1 file changed, 23 insertions(+), 18 deletions(-)

diff --git a/etc/caseDicts/foamyHexMeshDict b/etc/caseDicts/foamyHexMeshDict
index c866c832ad..cb8ac1517b 100644
--- a/etc/caseDicts/foamyHexMeshDict
+++ b/etc/caseDicts/foamyHexMeshDict
@@ -5,12 +5,19 @@
 |   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      foamyHexMeshDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 #inputMode merge;
 
 surfaceConformation
 {
-    locationInMesh                      (0 0 0);
     pointPairDistanceCoeff              0.1;
     mixedFeaturePointPPDistanceCoeff    5.0;
     featurePointExclusionDistanceCoeff  0.65;
@@ -22,11 +29,11 @@ surfaceConformation
 
     featurePointControls
     {
-        specialiseFeaturePoints on;
-        edgeAiming              on;
-        guardFeaturePoints      off;
-        snapFeaturePoints       off;
-        circulateEdges          off;
+        specialiseFeaturePoints         on;
+        edgeAiming                      on;
+        guardFeaturePoints              off;
+        snapFeaturePoints               off;
+        circulateEdges                  off;
     }
 
     conformationControls
@@ -115,20 +122,18 @@ polyMeshFiltering
 }
 
 
+backgroundMeshDecomposition
+{
+    minLevels           1;
+    sampleResolution    4;
+    spanScale           20;
+    maxCellWeightCoeff  20;
+}
+
+
 meshQualityControls
 {
-    maxNonOrtho         65;
-    maxBoundarySkewness 50;
-    maxInternalSkewness 10;
-    maxConcave          80;
-    minVol              -1E30;
-    minTetQuality       1e-30;
-    minArea             -1;
-    minTwist            0.02;
-    minDeterminant      0.001;
-    minFaceWeight       0.02;
-    minVolRatio         0.01;
-    minTriangleTwist    -1;
+    #include "meshQualityDict"
 }