From b6f350f6f7e1f9db9d2d2896c0118e3e5dd6d2e9 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 24 Nov 2015 13:42:58 +0000
Subject: [PATCH 1/8] STYLE: viewFactorsGen: spelling

---
 applications/utilities/preProcessing/viewFactorsGen/shootRays.H | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/applications/utilities/preProcessing/viewFactorsGen/shootRays.H b/applications/utilities/preProcessing/viewFactorsGen/shootRays.H
index d20b204f03..8b5dcd6db7 100644
--- a/applications/utilities/preProcessing/viewFactorsGen/shootRays.H
+++ b/applications/utilities/preProcessing/viewFactorsGen/shootRays.H
@@ -89,7 +89,7 @@ for (label procI = 0; procI < Pstream::nProcs(); procI++)
 
     DynamicList<label> dRayIs;
 
-    // Collect the rays which has not abstacle in bettween in rayStartFace
+    // Collect the rays which has not hit obstacle inbetween rayStartFace
     // and rayEndFace. If the ray hit itself get stored in dRayIs
     forAll(hitInfo, rayI)
     {

From 1fcd27af2a17214a512fc582db7e03159b900b55 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 24 Nov 2015 13:42:58 +0000
Subject: [PATCH 2/8] STYLE: viewFactorsGen: spelling

---
 applications/utilities/preProcessing/viewFactorsGen/shootRays.H | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/applications/utilities/preProcessing/viewFactorsGen/shootRays.H b/applications/utilities/preProcessing/viewFactorsGen/shootRays.H
index d20b204f03..8b5dcd6db7 100644
--- a/applications/utilities/preProcessing/viewFactorsGen/shootRays.H
+++ b/applications/utilities/preProcessing/viewFactorsGen/shootRays.H
@@ -89,7 +89,7 @@ for (label procI = 0; procI < Pstream::nProcs(); procI++)
 
     DynamicList<label> dRayIs;
 
-    // Collect the rays which has not abstacle in bettween in rayStartFace
+    // Collect the rays which has not hit obstacle inbetween rayStartFace
     // and rayEndFace. If the ray hit itself get stored in dRayIs
     forAll(hitInfo, rayI)
     {

From f3787e5722006de6fb92c8ae524288a0bb5e5702 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 24 Nov 2015 13:46:48 +0000
Subject: [PATCH 3/8] ENH: offsetSurface: new extrusion model to extrude
 inbetween two surfaces

---
 .../extrude/extrudeMesh/extrudeMeshDict       |  14 ++
 .../offsetSurface/offsetSurface.C             | 165 ++++++++++++++++++
 .../offsetSurface/offsetSurface.H             | 114 ++++++++++++
 3 files changed, 293 insertions(+)
 create mode 100644 src/mesh/extrudeModel/offsetSurface/offsetSurface.C
 create mode 100644 src/mesh/extrudeModel/offsetSurface/offsetSurface.H

diff --git a/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict b/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
index 0728864488..d3018d209d 100644
--- a/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
+++ b/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
@@ -63,6 +63,10 @@ extrudeModel        wedge;
 //- Extrudes into sphere with grading according to pressure (atmospherics)
 //extrudeModel        sigmaRadial;
 
+//- Extrudes by interpolating along path inbetween two (topologically identical)
+//  surfaces (e.g. one is an offsetted version of the other)
+//extrudeModel        offsetSurface;
+
 nLayers             10;
 
 expansionRatio      1.0;
@@ -105,6 +109,16 @@ sigmaRadialCoeffs
     pStrat          1;
 }
 
+offsetSurfaceCoeffs
+{
+    // Surface that mesh has been meshed to
+    baseSurface "$FOAM_CASE/constant/triSurface/DTC-scaled-inflated.obj";
+
+    // Surface to fill in to
+    offsetSurface "$FOAM_CASE/constant/triSurface/DTC-scaled.obj";
+}
+
+
 // Do front and back need to be merged? Usually only makes sense for 360
 // degree wedges.
 mergeFaces false;
diff --git a/src/mesh/extrudeModel/offsetSurface/offsetSurface.C b/src/mesh/extrudeModel/offsetSurface/offsetSurface.C
new file mode 100644
index 0000000000..1b2872279f
--- /dev/null
+++ b/src/mesh/extrudeModel/offsetSurface/offsetSurface.C
@@ -0,0 +1,165 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2014 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 "offsetSurface.H"
+#include "addToRunTimeSelectionTable.H"
+#include "triSurface.H"
+#include "triSurfaceSearch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace extrudeModels
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(offsetSurface, 0);
+
+addToRunTimeSelectionTable(extrudeModel, offsetSurface, dictionary);
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+offsetSurface::offsetSurface(const dictionary& dict)
+:
+    extrudeModel(typeName, dict),
+    project_(coeffDict_.lookupOrDefault("project", false))
+{
+    // Read surface
+    fileName baseName(coeffDict_.lookup("baseSurface"));
+    baseName.expand();
+    baseSurfPtr_.reset(new triSurface(baseName));
+
+    // Construct search engine
+    baseSearchPtr_.reset(new triSurfaceSearch(baseSurfPtr_()));
+
+    // Read offsetted surface
+    fileName offsetName(coeffDict_.lookup("offsetSurface"));
+    offsetName.expand();
+    offsetSurfPtr_.reset(new triSurface(offsetName));
+
+    // Construct search engine
+    offsetSearchPtr_.reset(new triSurfaceSearch(offsetSurfPtr_()));
+
+
+    const triSurface& b = baseSurfPtr_();
+    const triSurface& o = offsetSurfPtr_();
+
+    if
+    (
+        b.size() != o.size()
+     || b.nPoints() != o.nPoints()
+     || b.nEdges() != o.nEdges()
+    )
+    {
+        FatalIOErrorIn("offsetSurface::offsetSurface(const dictionary&)", dict)
+            << "offsetSurface " << offsetName
+            << " should have exactly the same topology as the baseSurface "
+            << baseName << exit(FatalIOError);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+offsetSurface::~offsetSurface()
+{}
+
+
+// * * * * * * * * * * * * * * * * Operators * * * * * * * * * * * * * * * * //
+
+point offsetSurface::operator()
+(
+    const point& surfacePoint,
+    const vector& surfaceNormal,
+    const label layer
+) const
+{
+    if (layer == 0)
+    {
+        return surfacePoint;
+    }
+    else
+    {
+        pointField samples(1, surfacePoint);
+        scalarField nearestDistSqr(1, GREAT);
+        List<pointIndexHit> info;
+        baseSearchPtr_().findNearest(samples, nearestDistSqr, info);
+
+        label triI = info[0].index();
+
+
+        const triSurface& base = baseSurfPtr_();
+        const triPointRef baseTri(base[triI].tri(base.points()));
+
+        List<scalar> bary;
+        baseTri.barycentric(surfacePoint, bary);
+
+        const triSurface& offset = offsetSurfPtr_();
+        const triPointRef offsetTri(offset[triI].tri(offset.points()));
+
+        const point offsetPoint
+        (
+            bary[0]*offsetTri.a()
+           +bary[1]*offsetTri.b()
+           +bary[2]*offsetTri.c()
+        );
+
+        point interpolatedPoint
+        (
+            surfacePoint + sumThickness(layer)*(offsetPoint-surfacePoint)
+        );
+
+
+        //- Either return interpolatedPoint or re-project onto surface (since
+        //  snapping might not have do so exactly)
+
+        if (project_)
+        {
+            // Re-project onto surface
+            offsetSearchPtr_().findNearest
+            (
+                pointField(1, interpolatedPoint),
+                scalarField(1, GREAT),
+                info
+            );
+            return info[0].hitPoint();
+        }
+        else
+        {
+            return interpolatedPoint;
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace extrudeModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/mesh/extrudeModel/offsetSurface/offsetSurface.H b/src/mesh/extrudeModel/offsetSurface/offsetSurface.H
new file mode 100644
index 0000000000..f5df55cfd8
--- /dev/null
+++ b/src/mesh/extrudeModel/offsetSurface/offsetSurface.H
@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2014 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::extrudeModels::offsetSurface
+
+Description
+    Extrudes by interpolating points from one surface to the other. Surfaces
+    have to be topologically identical i.e. one has to be an offsetted version
+    of the other.
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef offsetSurface_H
+#define offsetSurface_H
+
+#include "point.H"
+#include "extrudeModel.H"
+#include "Switch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class triSurface;
+class triSurfaceSearch;
+
+
+namespace extrudeModels
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class offsetSurface Declaration
+\*---------------------------------------------------------------------------*/
+
+class offsetSurface
+:
+    public extrudeModel
+{
+    // Private data
+
+        //- surface
+        autoPtr<triSurface> baseSurfPtr_;
+
+        //- search engine
+        autoPtr<triSurfaceSearch> baseSearchPtr_;
+
+        //- offsets
+        autoPtr<triSurface> offsetSurfPtr_;
+
+        //- search engine
+        autoPtr<triSurfaceSearch> offsetSearchPtr_;
+
+        // Whether to re-project onto offsetted surface
+        const Switch project_;
+
+public:
+
+    //- Runtime type information
+    TypeName("offsetSurface");
+
+    // Constructors
+
+        //- Construct from dictionary
+        offsetSurface(const dictionary& dict);
+
+
+    //- Destructor
+    virtual ~offsetSurface();
+
+
+    // Member Operators
+
+        //- Return point
+        point operator()
+        (
+            const point& surfacePoint,
+            const vector& surfaceNormal,
+            const label layer
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace extrudeModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

From 5e4d40ba70f88f88136b03844ab0832cee167ad9 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 24 Nov 2015 13:46:48 +0000
Subject: [PATCH 4/8] ENH: offsetSurface: new extrusion model to extrude
 inbetween two surfaces

---
 .../extrude/extrudeMesh/extrudeMeshDict       |  14 ++
 .../offsetSurface/offsetSurface.C             | 165 ++++++++++++++++++
 .../offsetSurface/offsetSurface.H             | 114 ++++++++++++
 3 files changed, 293 insertions(+)
 create mode 100644 src/mesh/extrudeModel/offsetSurface/offsetSurface.C
 create mode 100644 src/mesh/extrudeModel/offsetSurface/offsetSurface.H

diff --git a/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict b/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
index 0728864488..d3018d209d 100644
--- a/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
+++ b/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
@@ -63,6 +63,10 @@ extrudeModel        wedge;
 //- Extrudes into sphere with grading according to pressure (atmospherics)
 //extrudeModel        sigmaRadial;
 
+//- Extrudes by interpolating along path inbetween two (topologically identical)
+//  surfaces (e.g. one is an offsetted version of the other)
+//extrudeModel        offsetSurface;
+
 nLayers             10;
 
 expansionRatio      1.0;
@@ -105,6 +109,16 @@ sigmaRadialCoeffs
     pStrat          1;
 }
 
+offsetSurfaceCoeffs
+{
+    // Surface that mesh has been meshed to
+    baseSurface "$FOAM_CASE/constant/triSurface/DTC-scaled-inflated.obj";
+
+    // Surface to fill in to
+    offsetSurface "$FOAM_CASE/constant/triSurface/DTC-scaled.obj";
+}
+
+
 // Do front and back need to be merged? Usually only makes sense for 360
 // degree wedges.
 mergeFaces false;
diff --git a/src/mesh/extrudeModel/offsetSurface/offsetSurface.C b/src/mesh/extrudeModel/offsetSurface/offsetSurface.C
new file mode 100644
index 0000000000..1b2872279f
--- /dev/null
+++ b/src/mesh/extrudeModel/offsetSurface/offsetSurface.C
@@ -0,0 +1,165 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2014 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 "offsetSurface.H"
+#include "addToRunTimeSelectionTable.H"
+#include "triSurface.H"
+#include "triSurfaceSearch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace extrudeModels
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(offsetSurface, 0);
+
+addToRunTimeSelectionTable(extrudeModel, offsetSurface, dictionary);
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+offsetSurface::offsetSurface(const dictionary& dict)
+:
+    extrudeModel(typeName, dict),
+    project_(coeffDict_.lookupOrDefault("project", false))
+{
+    // Read surface
+    fileName baseName(coeffDict_.lookup("baseSurface"));
+    baseName.expand();
+    baseSurfPtr_.reset(new triSurface(baseName));
+
+    // Construct search engine
+    baseSearchPtr_.reset(new triSurfaceSearch(baseSurfPtr_()));
+
+    // Read offsetted surface
+    fileName offsetName(coeffDict_.lookup("offsetSurface"));
+    offsetName.expand();
+    offsetSurfPtr_.reset(new triSurface(offsetName));
+
+    // Construct search engine
+    offsetSearchPtr_.reset(new triSurfaceSearch(offsetSurfPtr_()));
+
+
+    const triSurface& b = baseSurfPtr_();
+    const triSurface& o = offsetSurfPtr_();
+
+    if
+    (
+        b.size() != o.size()
+     || b.nPoints() != o.nPoints()
+     || b.nEdges() != o.nEdges()
+    )
+    {
+        FatalIOErrorIn("offsetSurface::offsetSurface(const dictionary&)", dict)
+            << "offsetSurface " << offsetName
+            << " should have exactly the same topology as the baseSurface "
+            << baseName << exit(FatalIOError);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+offsetSurface::~offsetSurface()
+{}
+
+
+// * * * * * * * * * * * * * * * * Operators * * * * * * * * * * * * * * * * //
+
+point offsetSurface::operator()
+(
+    const point& surfacePoint,
+    const vector& surfaceNormal,
+    const label layer
+) const
+{
+    if (layer == 0)
+    {
+        return surfacePoint;
+    }
+    else
+    {
+        pointField samples(1, surfacePoint);
+        scalarField nearestDistSqr(1, GREAT);
+        List<pointIndexHit> info;
+        baseSearchPtr_().findNearest(samples, nearestDistSqr, info);
+
+        label triI = info[0].index();
+
+
+        const triSurface& base = baseSurfPtr_();
+        const triPointRef baseTri(base[triI].tri(base.points()));
+
+        List<scalar> bary;
+        baseTri.barycentric(surfacePoint, bary);
+
+        const triSurface& offset = offsetSurfPtr_();
+        const triPointRef offsetTri(offset[triI].tri(offset.points()));
+
+        const point offsetPoint
+        (
+            bary[0]*offsetTri.a()
+           +bary[1]*offsetTri.b()
+           +bary[2]*offsetTri.c()
+        );
+
+        point interpolatedPoint
+        (
+            surfacePoint + sumThickness(layer)*(offsetPoint-surfacePoint)
+        );
+
+
+        //- Either return interpolatedPoint or re-project onto surface (since
+        //  snapping might not have do so exactly)
+
+        if (project_)
+        {
+            // Re-project onto surface
+            offsetSearchPtr_().findNearest
+            (
+                pointField(1, interpolatedPoint),
+                scalarField(1, GREAT),
+                info
+            );
+            return info[0].hitPoint();
+        }
+        else
+        {
+            return interpolatedPoint;
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace extrudeModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/mesh/extrudeModel/offsetSurface/offsetSurface.H b/src/mesh/extrudeModel/offsetSurface/offsetSurface.H
new file mode 100644
index 0000000000..f5df55cfd8
--- /dev/null
+++ b/src/mesh/extrudeModel/offsetSurface/offsetSurface.H
@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2014 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::extrudeModels::offsetSurface
+
+Description
+    Extrudes by interpolating points from one surface to the other. Surfaces
+    have to be topologically identical i.e. one has to be an offsetted version
+    of the other.
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef offsetSurface_H
+#define offsetSurface_H
+
+#include "point.H"
+#include "extrudeModel.H"
+#include "Switch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class triSurface;
+class triSurfaceSearch;
+
+
+namespace extrudeModels
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class offsetSurface Declaration
+\*---------------------------------------------------------------------------*/
+
+class offsetSurface
+:
+    public extrudeModel
+{
+    // Private data
+
+        //- surface
+        autoPtr<triSurface> baseSurfPtr_;
+
+        //- search engine
+        autoPtr<triSurfaceSearch> baseSearchPtr_;
+
+        //- offsets
+        autoPtr<triSurface> offsetSurfPtr_;
+
+        //- search engine
+        autoPtr<triSurfaceSearch> offsetSearchPtr_;
+
+        // Whether to re-project onto offsetted surface
+        const Switch project_;
+
+public:
+
+    //- Runtime type information
+    TypeName("offsetSurface");
+
+    // Constructors
+
+        //- Construct from dictionary
+        offsetSurface(const dictionary& dict);
+
+
+    //- Destructor
+    virtual ~offsetSurface();
+
+
+    // Member Operators
+
+        //- Return point
+        point operator()
+        (
+            const point& surfacePoint,
+            const vector& surfaceNormal,
+            const label layer
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace extrudeModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

From 5c459f7f61bab9760a5ef4ab41c7fe0202091c96 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 24 Nov 2015 13:46:48 +0000
Subject: [PATCH 5/8] ENH: offsetSurface: new extrusion model to extrude
 inbetween two surfaces

---
 .../extrude/extrudeMesh/extrudeMeshDict       |  14 ++
 src/mesh/extrudeModel/Make/files              |   1 +
 src/mesh/extrudeModel/Make/options            |   1 +
 .../offsetSurface/offsetSurface.C             | 165 ++++++++++++++++++
 .../offsetSurface/offsetSurface.H             | 114 ++++++++++++
 5 files changed, 295 insertions(+)
 create mode 100644 src/mesh/extrudeModel/offsetSurface/offsetSurface.C
 create mode 100644 src/mesh/extrudeModel/offsetSurface/offsetSurface.H

diff --git a/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict b/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
index 0728864488..d3018d209d 100644
--- a/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
+++ b/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMeshDict
@@ -63,6 +63,10 @@ extrudeModel        wedge;
 //- Extrudes into sphere with grading according to pressure (atmospherics)
 //extrudeModel        sigmaRadial;
 
+//- Extrudes by interpolating along path inbetween two (topologically identical)
+//  surfaces (e.g. one is an offsetted version of the other)
+//extrudeModel        offsetSurface;
+
 nLayers             10;
 
 expansionRatio      1.0;
@@ -105,6 +109,16 @@ sigmaRadialCoeffs
     pStrat          1;
 }
 
+offsetSurfaceCoeffs
+{
+    // Surface that mesh has been meshed to
+    baseSurface "$FOAM_CASE/constant/triSurface/DTC-scaled-inflated.obj";
+
+    // Surface to fill in to
+    offsetSurface "$FOAM_CASE/constant/triSurface/DTC-scaled.obj";
+}
+
+
 // Do front and back need to be merged? Usually only makes sense for 360
 // degree wedges.
 mergeFaces false;
diff --git a/src/mesh/extrudeModel/Make/files b/src/mesh/extrudeModel/Make/files
index f6c2f50b04..e18ad638e1 100644
--- a/src/mesh/extrudeModel/Make/files
+++ b/src/mesh/extrudeModel/Make/files
@@ -4,6 +4,7 @@ linearNormal/linearNormal.C
 planeExtrusion/planeExtrusion.C
 linearDirection/linearDirection.C
 linearRadial/linearRadial.C
+offsetSurface/offsetSurface.C
 radial/radial.C
 sigmaRadial/sigmaRadial.C
 sector/sector.C
diff --git a/src/mesh/extrudeModel/Make/options b/src/mesh/extrudeModel/Make/options
index eabd0b53a8..2dd02e7d86 100644
--- a/src/mesh/extrudeModel/Make/options
+++ b/src/mesh/extrudeModel/Make/options
@@ -1,5 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/dynamicMesh/lnInclude
 
 LIB_LIBS = \
diff --git a/src/mesh/extrudeModel/offsetSurface/offsetSurface.C b/src/mesh/extrudeModel/offsetSurface/offsetSurface.C
new file mode 100644
index 0000000000..1b2872279f
--- /dev/null
+++ b/src/mesh/extrudeModel/offsetSurface/offsetSurface.C
@@ -0,0 +1,165 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2014 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 "offsetSurface.H"
+#include "addToRunTimeSelectionTable.H"
+#include "triSurface.H"
+#include "triSurfaceSearch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace extrudeModels
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(offsetSurface, 0);
+
+addToRunTimeSelectionTable(extrudeModel, offsetSurface, dictionary);
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+offsetSurface::offsetSurface(const dictionary& dict)
+:
+    extrudeModel(typeName, dict),
+    project_(coeffDict_.lookupOrDefault("project", false))
+{
+    // Read surface
+    fileName baseName(coeffDict_.lookup("baseSurface"));
+    baseName.expand();
+    baseSurfPtr_.reset(new triSurface(baseName));
+
+    // Construct search engine
+    baseSearchPtr_.reset(new triSurfaceSearch(baseSurfPtr_()));
+
+    // Read offsetted surface
+    fileName offsetName(coeffDict_.lookup("offsetSurface"));
+    offsetName.expand();
+    offsetSurfPtr_.reset(new triSurface(offsetName));
+
+    // Construct search engine
+    offsetSearchPtr_.reset(new triSurfaceSearch(offsetSurfPtr_()));
+
+
+    const triSurface& b = baseSurfPtr_();
+    const triSurface& o = offsetSurfPtr_();
+
+    if
+    (
+        b.size() != o.size()
+     || b.nPoints() != o.nPoints()
+     || b.nEdges() != o.nEdges()
+    )
+    {
+        FatalIOErrorIn("offsetSurface::offsetSurface(const dictionary&)", dict)
+            << "offsetSurface " << offsetName
+            << " should have exactly the same topology as the baseSurface "
+            << baseName << exit(FatalIOError);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+offsetSurface::~offsetSurface()
+{}
+
+
+// * * * * * * * * * * * * * * * * Operators * * * * * * * * * * * * * * * * //
+
+point offsetSurface::operator()
+(
+    const point& surfacePoint,
+    const vector& surfaceNormal,
+    const label layer
+) const
+{
+    if (layer == 0)
+    {
+        return surfacePoint;
+    }
+    else
+    {
+        pointField samples(1, surfacePoint);
+        scalarField nearestDistSqr(1, GREAT);
+        List<pointIndexHit> info;
+        baseSearchPtr_().findNearest(samples, nearestDistSqr, info);
+
+        label triI = info[0].index();
+
+
+        const triSurface& base = baseSurfPtr_();
+        const triPointRef baseTri(base[triI].tri(base.points()));
+
+        List<scalar> bary;
+        baseTri.barycentric(surfacePoint, bary);
+
+        const triSurface& offset = offsetSurfPtr_();
+        const triPointRef offsetTri(offset[triI].tri(offset.points()));
+
+        const point offsetPoint
+        (
+            bary[0]*offsetTri.a()
+           +bary[1]*offsetTri.b()
+           +bary[2]*offsetTri.c()
+        );
+
+        point interpolatedPoint
+        (
+            surfacePoint + sumThickness(layer)*(offsetPoint-surfacePoint)
+        );
+
+
+        //- Either return interpolatedPoint or re-project onto surface (since
+        //  snapping might not have do so exactly)
+
+        if (project_)
+        {
+            // Re-project onto surface
+            offsetSearchPtr_().findNearest
+            (
+                pointField(1, interpolatedPoint),
+                scalarField(1, GREAT),
+                info
+            );
+            return info[0].hitPoint();
+        }
+        else
+        {
+            return interpolatedPoint;
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace extrudeModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/mesh/extrudeModel/offsetSurface/offsetSurface.H b/src/mesh/extrudeModel/offsetSurface/offsetSurface.H
new file mode 100644
index 0000000000..f5df55cfd8
--- /dev/null
+++ b/src/mesh/extrudeModel/offsetSurface/offsetSurface.H
@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2014 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::extrudeModels::offsetSurface
+
+Description
+    Extrudes by interpolating points from one surface to the other. Surfaces
+    have to be topologically identical i.e. one has to be an offsetted version
+    of the other.
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef offsetSurface_H
+#define offsetSurface_H
+
+#include "point.H"
+#include "extrudeModel.H"
+#include "Switch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class triSurface;
+class triSurfaceSearch;
+
+
+namespace extrudeModels
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class offsetSurface Declaration
+\*---------------------------------------------------------------------------*/
+
+class offsetSurface
+:
+    public extrudeModel
+{
+    // Private data
+
+        //- surface
+        autoPtr<triSurface> baseSurfPtr_;
+
+        //- search engine
+        autoPtr<triSurfaceSearch> baseSearchPtr_;
+
+        //- offsets
+        autoPtr<triSurface> offsetSurfPtr_;
+
+        //- search engine
+        autoPtr<triSurfaceSearch> offsetSearchPtr_;
+
+        // Whether to re-project onto offsetted surface
+        const Switch project_;
+
+public:
+
+    //- Runtime type information
+    TypeName("offsetSurface");
+
+    // Constructors
+
+        //- Construct from dictionary
+        offsetSurface(const dictionary& dict);
+
+
+    //- Destructor
+    virtual ~offsetSurface();
+
+
+    // Member Operators
+
+        //- Return point
+        point operator()
+        (
+            const point& surfacePoint,
+            const vector& surfaceNormal,
+            const label layer
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace extrudeModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

From 251839df0d850002dd37f54221bc879e0d1afcff Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 24 Nov 2015 13:46:48 +0000
Subject: [PATCH 6/8] ENH: offsetSurface: new extrusion model to extrude
 inbetween two surfaces

---
 src/mesh/extrudeModel/Make/files   | 1 +
 src/mesh/extrudeModel/Make/options | 1 +
 2 files changed, 2 insertions(+)

diff --git a/src/mesh/extrudeModel/Make/files b/src/mesh/extrudeModel/Make/files
index f6c2f50b04..e18ad638e1 100644
--- a/src/mesh/extrudeModel/Make/files
+++ b/src/mesh/extrudeModel/Make/files
@@ -4,6 +4,7 @@ linearNormal/linearNormal.C
 planeExtrusion/planeExtrusion.C
 linearDirection/linearDirection.C
 linearRadial/linearRadial.C
+offsetSurface/offsetSurface.C
 radial/radial.C
 sigmaRadial/sigmaRadial.C
 sector/sector.C
diff --git a/src/mesh/extrudeModel/Make/options b/src/mesh/extrudeModel/Make/options
index eabd0b53a8..2dd02e7d86 100644
--- a/src/mesh/extrudeModel/Make/options
+++ b/src/mesh/extrudeModel/Make/options
@@ -1,5 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/dynamicMesh/lnInclude
 
 LIB_LIBS = \

From c020702403704eea19546769a0ab7f054e174b3a Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 24 Nov 2015 17:35:18 +0000
Subject: [PATCH 7/8] ENH: parallel: added -decomposeParDict option to -
 foamyHexMesh - snappyHexMesh - decomposePar

---
 .../conformalVoronoiMesh/Make/options         |   2 +
 .../backgroundMeshDecomposition.C             |  31 +++--
 .../backgroundMeshDecomposition.H             |  13 +-
 .../backgroundMeshDecompositionI.H            |   7 --
 .../cellShapeControlMesh.C                    |  96 +++++++-------
 .../conformalVoronoiMesh.C                    |   6 +-
 .../conformalVoronoiMesh.H                    |  10 +-
 .../conformalVoronoiMeshCalcDualMesh.C        |   6 +-
 .../indexedCell/indexedCell.H                 |   4 +-
 .../foamyMesh/foamyHexMesh/Make/options       |   3 +-
 .../foamyMesh/foamyHexMesh/foamyHexMesh.C     |  15 ++-
 .../foamyMesh/foamyHexMesh/foamyHexMeshDict   |  17 ++-
 .../foamyHexMeshBackgroundMesh/Make/options   |   2 +
 .../foamyHexMeshBackgroundMesh.C              |  49 +++-----
 .../foamyMesh/foamyQuadMesh/foamyQuadMesh.C   |   2 +-
 .../generation/snappyHexMesh/Make/options     |   4 +
 .../generation/snappyHexMesh/snappyHexMesh.C  |  28 +++--
 .../decomposePar/decomposePar.C               |   4 +-
 .../preProcessing/mapFields/Make/options      |   3 +
 .../preProcessing/mapFields/mapFields.C       | 119 ++++++++++--------
 .../mapFieldsPar/MapLagrangianFields.H        |  41 ++++--
 .../mapFieldsPar/mapLagrangian.C              |   5 +-
 22 files changed, 265 insertions(+), 202 deletions(-)

diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/Make/options b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/Make/options
index 44c1d87758..bad58b9d58 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/Make/options
@@ -16,6 +16,7 @@ EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
     -I$(LIB_SRC)/edgeMesh/lnInclude \
     -I$(LIB_SRC)/fileFormats/lnInclude \
     -I$(LIB_SRC)/dynamicMesh/lnInclude \
@@ -29,6 +30,7 @@ EXE_INC = \
 LIB_LIBS = \
     ${CGAL_LIBS} \
     -lmeshTools \
+    -ldecompose \
     -ledgeMesh \
     -lfileFormats \
     -ltriSurface \
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
index 25e0f9a73e..b4428c3df8 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
@@ -30,6 +30,7 @@ License
 #include "Time.H"
 #include "Random.H"
 #include "pointConversion.H"
+#include "decompositionModel.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -149,7 +150,8 @@ void Foam::backgroundMeshDecomposition::initialRefinement()
 
     const conformationSurfaces& geometry = geometryToConformTo_;
 
-    decompositionMethod& decomposer = decomposerPtr_();
+    decompositionMethod& decomposer =
+        decompositionModel::New(mesh_).decomposer();
 
     volScalarField::InternalField& icellWeights = cellWeights.internalField();
 
@@ -782,7 +784,8 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
     const Time& runTime,
     Random& rndGen,
     const conformationSurfaces& geometryToConformTo,
-    const dictionary& coeffsDict
+    const dictionary& coeffsDict,
+    const fileName& decompDictFile
 )
 :
     runTime_(runTime),
@@ -810,18 +813,6 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
     bFTreePtr_(),
     allBackgroundMeshBounds_(Pstream::nProcs()),
     globalBackgroundBounds_(),
-    decomposeDict_
-    (
-        IOobject
-        (
-            "decomposeParDict",
-            runTime_.system(),
-            runTime_,
-            IOobject::MUST_READ_IF_MODIFIED,
-            IOobject::NO_WRITE
-        )
-    ),
-    decomposerPtr_(decompositionMethod::New(decomposeDict_)),
     mergeDist_(1e-6*mesh_.bounds().mag()),
     spanScale_(readScalar(coeffsDict.lookup("spanScale"))),
     minCellSizeLimit_
@@ -846,14 +837,17 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
             << exit(FatalError);
     }
 
-    if (!decomposerPtr_().parallelAware())
+    const decompositionMethod& decomposer =
+        decompositionModel::New(mesh_, decompDictFile).decomposer();
+
+    if (!decomposer.parallelAware())
     {
         FatalErrorIn
         (
             "void Foam::backgroundMeshDecomposition::initialRefinement() const"
         )
             << "You have selected decomposition method "
-            << decomposerPtr_().typeName
+            << decomposer.typeName
             << " which is not parallel aware." << endl
             << exit(FatalError);
     }
@@ -1008,7 +1002,10 @@ Foam::backgroundMeshDecomposition::distribute
             << endl;
     }
 
-    labelList newDecomp = decomposerPtr_().decompose
+    decompositionMethod& decomposer =
+        decompositionModel::New(mesh_).decomposer();
+
+    labelList newDecomp = decomposer.decompose
     (
         mesh_,
         mesh_.cellCentres(),
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
index 90dde14aa4..9d597732cb 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
@@ -127,13 +127,7 @@ class backgroundMeshDecomposition
         //  a point that is not found on any processor is in the domain at all
         treeBoundBox globalBackgroundBounds_;
 
-        //- Decomposition dictionary
-        IOdictionary decomposeDict_;
-
-        //- Decomposition method
-        autoPtr<decompositionMethod> decomposerPtr_;
-
-        //- Merge distance required by fvMeshDistribute
+        //- merge distance required by fvMeshDistribute
         scalar mergeDist_;
 
         //- Scale of a cell span vs cell size used to decide to refine a cell
@@ -204,7 +198,8 @@ public:
             const Time& runTime,
             Random& rndGen,
             const conformationSurfaces& geometryToConformTo,
-            const dictionary& coeffsDict
+            const dictionary& coeffsDict,
+            const fileName& decompDictFile = ""
         );
 
 
@@ -324,8 +319,6 @@ public:
             //- Return the point level of the underlying mesh
             inline const labelList& pointLevel() const;
 
-            //- Return the current decomposition method
-            inline const decompositionMethod& decomposer() const;
 };
 
 
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
index 199f4ae4ff..b8dcb5187f 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
@@ -57,11 +57,4 @@ const Foam::labelList& Foam::backgroundMeshDecomposition::pointLevel() const
 }
 
 
-const Foam::decompositionMethod&
-Foam::backgroundMeshDecomposition::decomposer() const
-{
-    return decomposerPtr_();
-}
-
-
 // ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
index e7f7f5816c..9b6f71ffac 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
@@ -392,61 +392,59 @@ Foam::cellShapeControlMesh::cellShapeControlMesh(const Time& runTime)
 
         if (mesh.nPoints() == this->vertexCount())
         {
-            pointScalarField sizes
+            IOobject io
             (
-                IOobject
-                (
-                    "sizes",
-                    runTime.timeName(),
-                    meshSubDir,
-                    runTime,
-                    IOobject::READ_IF_PRESENT,
-                    IOobject::NO_WRITE,
-                    false
-                ),
-                pointMesh::New(mesh)
+                "sizes",
+                runTime.timeName(),
+                meshSubDir,
+                runTime,
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE,
+                false
             );
 
-            triadIOField alignments
-            (
-                IOobject
-                (
-                    "alignments",
-                    mesh.time().timeName(),
-                    meshSubDir,
-                    mesh.time(),
-                    IOobject::READ_IF_PRESENT,
-                    IOobject::AUTO_WRITE,
-                    false
-                )
-            );
-
-            if
-            (
-                sizes.size() == this->vertexCount()
-             && alignments.size() == this->vertexCount()
-            )
+            if (io.headerOk())
             {
-                for
+                pointScalarField sizes(io, pointMesh::New(mesh));
+
+                triadIOField alignments
                 (
-                    Finite_vertices_iterator vit = finite_vertices_begin();
-                    vit != finite_vertices_end();
-                    ++vit
-                )
+                    IOobject
+                    (
+                        "alignments",
+                        mesh.time().timeName(),
+                        meshSubDir,
+                        mesh.time(),
+                        IOobject::MUST_READ,
+                        IOobject::NO_WRITE,
+                        false
+                    )
+                );
+
+                if (alignments.size() == this->vertexCount())
                 {
-                    vit->targetCellSize() = sizes[vit->index()];
-                    vit->alignment() = alignments[vit->index()];
+                    for
+                    (
+                        Finite_vertices_iterator vit = finite_vertices_begin();
+                        vit != finite_vertices_end();
+                        ++vit
+                    )
+                    {
+                        vit->targetCellSize() = sizes[vit->index()];
+                        vit->alignment() = alignments[vit->index()];
+                    }
+                }
+                else
+                {
+                    FatalErrorIn
+                    (
+                        "Foam::cellShapeControlMesh::cellShapeControlMesh"
+                        "(const Time&)"
+                    )   << "Cell alignments point field " << alignments.size()
+                        << " is not the same size as the number of vertices"
+                        << " in the mesh " << this->vertexCount()
+                        << abort(FatalError);
                 }
-            }
-            else
-            {
-                FatalErrorIn
-                (
-                    "Foam::cellShapeControlMesh::cellShapeControlMesh"
-                    "(const Time&)"
-                )   << "Cell size point field is not the same size as the "
-                    << "mesh."
-                    << abort(FatalError);
             }
         }
     }
@@ -672,7 +670,7 @@ void Foam::cellShapeControlMesh::write() const
             IOobject::AUTO_WRITE
         ),
         pointMesh::New(mesh),
-        scalar(0)
+        dimensionedScalar("zero", dimLength, scalar(0))
     );
 
     triadIOField alignments
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
index 9cc9862129..a10f19be22 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
@@ -837,7 +837,8 @@ bool Foam::conformalVoronoiMesh::ownerAndNeighbour
 Foam::conformalVoronoiMesh::conformalVoronoiMesh
 (
     const Time& runTime,
-    const dictionary& foamyHexMeshDict
+    const dictionary& foamyHexMeshDict,
+    const fileName& decompDictFile
 )
 :
     DistributedDelaunayMesh<Delaunay>(runTime),
@@ -876,7 +877,8 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
             foamyHexMeshControls().foamyHexMeshDict().subDict
             (
                 "backgroundMeshDecomposition"
-            )
+            ),
+            decompDictFile
         )
       : NULL
     ),
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
index f6d44489b7..59fbba90a9 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
@@ -163,10 +163,8 @@ private:
         //- Limiting bound box before infinity begins
         treeBoundBox limitBounds_;
 
-        //-
         mutable pointPairs<Delaunay> ptPairs_;
 
-        //-
         featurePointConformer ftPtConformer_;
 
         //- Search tree for edge point locations
@@ -546,7 +544,7 @@ private:
         ) const;
 
         //- Check if a location is in the exclusion range of an existing feature
-        //- Edge conformation location
+        //  edge conformation location
         bool nearFeatureEdgeLocation
         (
             const pointIndexHit& pHit,
@@ -730,8 +728,7 @@ private:
 
         label classifyBoundaryPoint(Cell_handle cit) const;
 
-        //- Index all of the the Delaunay cells and calculate their
-        //- Dual points
+        //- Index all of the the Delaunay cells and calculate their dual points
         void indexDualVertices
         (
             pointField& pts,
@@ -874,7 +871,8 @@ public:
         conformalVoronoiMesh
         (
             const Time& runTime,
-            const dictionary& foamyHexMeshDict
+            const dictionary& foamyHexMeshDict,
+            const fileName& decompDictFile = ""
         );
 
 
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
index e88ee94369..85c49a9ccd 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
@@ -1292,9 +1292,9 @@ void Foam::conformalVoronoiMesh::indexDualVertices
         }
     }
 
-    OBJstream snapping1("snapToSurface1.obj");
-    OBJstream snapping2("snapToSurface2.obj");
-    OFstream tetToSnapTo("tetsToSnapTo.obj");
+    //OBJstream snapping1("snapToSurface1.obj");
+    //OBJstream snapping2("snapToSurface2.obj");
+    //OFstream tetToSnapTo("tetsToSnapTo.obj");
 
     for
     (
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
index a2eeb21c99..7a19a25ba2 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
@@ -88,15 +88,13 @@ class indexedCell
     // Private data
 
         //- The index for this Delaunay tetrahedral cell.  Type information is
-        //- Also carried:
-
+        //  also carried:
         //  ctFar : the dual point of this cell does not form part of the
         //          internal or boundary of the dual mesh
         //  >=0   : the (local) index of an internal or boundary dual point,
         //           not on a processor face
         //  < 0 && > ctFar : the (global) index of a dual point on a processor
         //                   face
-
         Foam::label index_;
 
         //- The number of times that this Delaunay cell has been limited
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/Make/options b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/Make/options
index b3f2e7dd7c..7aa621f26e 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/Make/options
@@ -30,7 +30,8 @@ EXE_LIBS = \
     -lconformalVoronoiMesh \
     -lmeshTools \
     -ldecompositionMethods \
-    -L$(FOAM_LIBBIN)/dummy -lptscotchDecomp \
+    -ldecompose \
+    -L$(FOAM_LIBBIN)/dummy -lptscotchDecomp -lscotchDecomp \
     -ledgeMesh \
     -lfileFormats \
     -ltriSurface \
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMesh.C b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMesh.C
index 6e68f1c80e..c83800e031 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMesh.C
@@ -62,6 +62,17 @@ int main(int argc, char *argv[])
     const bool checkGeometry = args.optionFound("checkGeometry");
     const bool conformationOnly = args.optionFound("conformationOnly");
 
+    // Allow override of decomposeParDict location
+    fileName decompDictFile;
+    if (args.optionReadIfPresent("decomposeParDict", decompDictFile))
+    {
+        if (isDir(decompDictFile))
+        {
+            decompDictFile = decompDictFile / "decomposeParDict";
+        }
+    }
+
+
     IOdictionary foamyHexMeshDict
     (
         IOobject
@@ -114,7 +125,7 @@ int main(int argc, char *argv[])
 
     Info<< "Create mesh for time = " << runTime.timeName() << nl << endl;
 
-    conformalVoronoiMesh mesh(runTime, foamyHexMeshDict);
+    conformalVoronoiMesh mesh(runTime, foamyHexMeshDict, decompDictFile);
 
 
     if (conformationOnly)
@@ -145,7 +156,7 @@ int main(int argc, char *argv[])
     }
 
 
-    Info<< nl << "End" << nl << endl;
+    Info<< "\nEnd\n" << endl;
 
     return 0;
 }
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMeshDict b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMeshDict
index deca4e0edb..78f22b4e49 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMeshDict
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMeshDict
@@ -271,11 +271,18 @@ motionControl
             forceInitialPointInsertion  on;
             priority                    1;
             mode                        inside;
+
+            // Cell size at surface
             surfaceCellSizeFunction     uniformValue;
             uniformValueCoeffs
             {
                 surfaceCellSizeCoeff    0.5;
             }
+
+            // Cell size inside domain by having a region of thickness
+            // surfaceOffsetaround the surface with the surface cell size
+            // (so constant) and then down to distanceCellSize over a distance
+            // of linearDistance.
             cellSizeFunction            surfaceOffsetLinearDistance;
             surfaceOffsetLinearDistanceCoeffs
             {
@@ -375,9 +382,17 @@ polyMeshFiltering
     // Filter small and sliver faces
     filterFaces                             off;
 
-    // Write the underlying Delaunay tet mesh at output time
+    // Write the underlying Delaunay tet mesh (at output time)
     writeTetDualMesh                        false;
 
+    // Write the Delaunay tet mesh used for interpolating cell size and
+    // alignment (at output time)
+    writeCellShapeControlMesh           true;
+
+    // Write the hex/split-hex mesh used for parallel load balancing
+    // (at output time)
+    writeBackgroundMeshDecomposition    true;
+
     // Upper limit on the size of faces to be filtered.
     // fraction of the local target cell size
     filterSizeCoeff                         0.2;
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/Make/options b/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/Make/options
index 07f6800de6..00d889c6c6 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/Make/options
@@ -9,6 +9,7 @@ EXE_INC = \
     ${CGAL_INC} \
     -I../conformalVoronoiMesh/lnInclude \
     -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
     -I$(LIB_SRC)/edgeMesh/lnInclude \
     -I$(LIB_SRC)/mesh/autoMesh/lnInclude \
     -I$(LIB_SRC)/triSurface/lnInclude \
@@ -26,6 +27,7 @@ EXE_LIBS = \
     -lgmp \
     -lconformalVoronoiMesh \
     -ldecompositionMethods /* -L$(FOAM_LIBBIN)/dummy -lscotchDecomp */ \
+    -ldecompose \
     -ledgeMesh \
     -ltriSurface \
     -lmeshTools \
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C b/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C
index fc99f0fd05..62be3f1efe 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C
@@ -44,6 +44,7 @@ Description
 #include "isoSurfaceCell.H"
 #include "vtkSurfaceWriter.H"
 #include "syncTools.H"
+#include "decompositionModel.H"
 
 using namespace Foam;
 
@@ -467,7 +468,7 @@ int main(int argc, char *argv[])
 
         // Determine the number of cells in each direction.
         const vector span = bb.span();
-        vector nScalarCells = span/cellShapeControls().defaultCellSize();
+        vector nScalarCells = span/cellShapeControls.defaultCellSize();
 
         // Calculate initial cell size to be a little bit smaller than the
         // defaultCellSize to avoid initial refinement triggering.
@@ -521,28 +522,21 @@ int main(int argc, char *argv[])
         Info<< "Loaded mesh:" << endl;
         printMeshData(mesh);
 
-        // Allocate a decomposer
-        IOdictionary decompositionDict
-        (
-            IOobject
-            (
-                "decomposeParDict",
-                runTime.system(),
-                mesh,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
-        );
+        // Allow override of decomposeParDict location
+        fileName decompDictFile;
+        if (args.optionReadIfPresent("decomposeParDict", decompDictFile))
+        {
+            if (isDir(decompDictFile))
+            {
+                decompDictFile = decompDictFile / "decomposeParDict";
+            }
+        }
 
-        autoPtr<decompositionMethod> decomposer
+        labelList decomp = decompositionModel::New
         (
-            decompositionMethod::New
-            (
-                decompositionDict
-            )
-        );
-
-        labelList decomp = decomposer().decompose(mesh, mesh.cellCentres());
+            mesh,
+            decompDictFile
+        ).decomposer().decompose(mesh, mesh.cellCentres());
 
         // Global matching tolerance
         const scalar tolDim = getMergeDistance
@@ -574,18 +568,15 @@ int main(int argc, char *argv[])
     Info<< "Refining backgroud mesh according to cell size specification" << nl
         << endl;
 
+    const dictionary& backgroundMeshDict =
+        foamyHexMeshDict.subDict("backgroundMeshDecomposition");
+
     backgroundMeshDecomposition backgroundMesh
     (
-        1.0,    //spanScale,ratio of poly cell size v.s. hex cell size
-        0.0,    //minCellSizeLimit
-        0,      //minLevels
-        4,      //volRes, check multiple points per cell
-        20.0,   //maxCellWeightCoeff
         runTime,
-        geometryToConformTo,
-        cellShapeControls(),
         rndGen,
-        foamyHexMeshDict
+        geometryToConformTo,
+        backgroundMeshDict
     );
 
     if (writeMesh)
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyQuadMesh/foamyQuadMesh.C b/applications/utilities/mesh/generation/foamyMesh/foamyQuadMesh/foamyQuadMesh.C
index 066bc8205c..997726bba6 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyQuadMesh/foamyQuadMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyQuadMesh/foamyQuadMesh.C
@@ -221,7 +221,7 @@ int main(int argc, char *argv[])
     Info<< "Finished extruding in = "
         << runTime.cpuTimeIncrement() << " s." << endl;
 
-    Info<< nl << "End\n" << endl;
+    Info<< "\nEnd\n" << endl;
 
     return 0;
 }
diff --git a/applications/utilities/mesh/generation/snappyHexMesh/Make/options b/applications/utilities/mesh/generation/snappyHexMesh/Make/options
index 219e3d419c..d62b82a239 100644
--- a/applications/utilities/mesh/generation/snappyHexMesh/Make/options
+++ b/applications/utilities/mesh/generation/snappyHexMesh/Make/options
@@ -8,14 +8,18 @@ EXE_INC = \
     -I$(LIB_SRC)/surfMesh/lnInclude \
     -I$(LIB_SRC)/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/edgeMesh/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 
 EXE_LIBS = \
     -lfiniteVolume \
     -ldecompositionMethods \
     -L$(FOAM_LIBBIN)/dummy -lptscotchDecomp \
+    /* note: scotch < 6.0 does not like both scotch and ptscotch together */ \
+    -lscotchDecomp \
     -lmeshTools \
     -lsurfMesh \
     -lfileFormats \
     -ldynamicMesh \
+    -ldecompose \
     -lautoMesh
diff --git a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
index 7b4d8577c4..15bd5e965f 100644
--- a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
+++ b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+     \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -57,6 +57,7 @@ Description
 #include "MeshedSurface.H"
 #include "globalIndex.H"
 #include "IOmanip.H"
+#include "decompositionModel.H"
 
 using namespace Foam;
 
@@ -819,15 +820,28 @@ int main(int argc, char *argv[])
     {
         if (Pstream::parRun())
         {
+            fileName decompDictFile;
+            if (args.optionReadIfPresent("decomposeParDict", decompDictFile))
+            {
+                if (isDir(decompDictFile))
+                {
+                    decompDictFile = decompDictFile/"decomposeParDict";
+                }
+            }
+
             decomposeDict = IOdictionary
             (
-                IOobject
+                decompositionModel::selectIO
                 (
-                    "decomposeParDict",
-                    runTime.system(),
-                    mesh,
-                    IOobject::MUST_READ_IF_MODIFIED,
-                    IOobject::NO_WRITE
+                    IOobject
+                    (
+                        "decomposeParDict",
+                        runTime.system(),
+                        mesh,
+                        IOobject::MUST_READ_IF_MODIFIED,
+                        IOobject::NO_WRITE
+                    ),
+                    decompDictFile
                 )
             );
         }
diff --git a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
index 63303d9f8c..9c62a12837 100644
--- a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
+++ b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+     \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -143,7 +143,7 @@ int main(int argc, char *argv[])
     );
 
     argList::noParallel();
-    Foam::argList::addOption
+    argList::addOption
     (
         "decomposeParDict",
         "file",
diff --git a/applications/utilities/preProcessing/mapFields/Make/options b/applications/utilities/preProcessing/mapFields/Make/options
index 7bd964094e..383ac3c12c 100644
--- a/applications/utilities/preProcessing/mapFields/Make/options
+++ b/applications/utilities/preProcessing/mapFields/Make/options
@@ -2,6 +2,8 @@ EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/lagrangian/basic/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude
 
 EXE_LIBS = \
@@ -9,4 +11,5 @@ EXE_LIBS = \
     -lmeshTools \
     -llagrangian \
     -lfiniteVolume \
+    -ldecompose \
     -lgenericPatchFields
diff --git a/applications/utilities/preProcessing/mapFields/mapFields.C b/applications/utilities/preProcessing/mapFields/mapFields.C
index 929776d9df..6f72053a08 100644
--- a/applications/utilities/preProcessing/mapFields/mapFields.C
+++ b/applications/utilities/preProcessing/mapFields/mapFields.C
@@ -36,9 +36,48 @@ Description
 #include "meshToMesh0.H"
 #include "processorFvPatch.H"
 #include "MapMeshes.H"
+#include "decompositionModel.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+int readNumProcs
+(
+    const argList& args,
+    const word& optionName,
+    const Time& runTime
+)
+{
+    fileName dictFile;
+    if (args.optionReadIfPresent(optionName, dictFile))
+    {
+        if (isDir(dictFile))
+        {
+            dictFile = dictFile/"decomposeParDict";
+        }
+    }
+
+    return readInt
+    (
+        IOdictionary
+        (
+            decompositionModel::selectIO
+            (
+                IOobject
+                (
+                    "decomposeParDict",
+                    runTime.system(),
+                    runTime,
+                    IOobject::MUST_READ_IF_MODIFIED,
+                    IOobject::NO_WRITE,
+                    false
+                ),
+                dictFile
+            )
+        ).lookup("numberOfSubdomains")
+    );
+}
+
+
 void mapConsistentMesh
 (
     const fvMesh& meshSource,
@@ -225,6 +264,19 @@ int main(int argc, char *argv[])
         "subtract",
         "subtract mapped source from target"
     );
+    argList::addOption
+    (
+        "sourceDecomposeParDict",
+        "file",
+        "read decomposePar dictionary from specified location"
+    );
+    argList::addOption
+    (
+        "targetDecomposeParDict",
+        "file",
+        "read decomposePar dictionary from specified location"
+    );
+
 
     argList args(argc, argv);
 
@@ -320,20 +372,13 @@ int main(int argc, char *argv[])
 
     if (parallelSource && !parallelTarget)
     {
-        IOdictionary decompositionDict
+        int nProcs = readNumProcs
         (
-            IOobject
-            (
-                "decomposeParDict",
-                runTimeSource.system(),
-                runTimeSource,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
+            args,
+            "sourceDecomposeParDict",
+            runTimeSource
         );
 
-        int nProcs(readInt(decompositionDict.lookup("numberOfSubdomains")));
-
         Info<< "Create target mesh\n" << endl;
 
         fvMesh meshTarget
@@ -399,19 +444,13 @@ int main(int argc, char *argv[])
     }
     else if (!parallelSource && parallelTarget)
     {
-        IOdictionary decompositionDict
+        int nProcs = readNumProcs
         (
-            IOobject
-            (
-                "decomposeParDict",
-                runTimeTarget.system(),
-                runTimeTarget,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
+            args,
+            "targetDecomposeParDict",
+            runTimeTarget
         );
 
-        int nProcs(readInt(decompositionDict.lookup("numberOfSubdomains")));
 
         Info<< "Create source mesh\n" << endl;
 
@@ -478,39 +517,17 @@ int main(int argc, char *argv[])
     }
     else if (parallelSource && parallelTarget)
     {
-        IOdictionary decompositionDictSource
+        int nProcsSource = readNumProcs
         (
-            IOobject
-            (
-                "decomposeParDict",
-                runTimeSource.system(),
-                runTimeSource,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
+            args,
+            "sourceDecomposeParDict",
+            runTimeSource
         );
-
-        int nProcsSource
+        int nProcsTarget = readNumProcs
         (
-            readInt(decompositionDictSource.lookup("numberOfSubdomains"))
-        );
-
-
-        IOdictionary decompositionDictTarget
-        (
-            IOobject
-            (
-                "decomposeParDict",
-                runTimeTarget.system(),
-                runTimeTarget,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
-        );
-
-        int nProcsTarget
-        (
-            readInt(decompositionDictTarget.lookup("numberOfSubdomains"))
+            args,
+            "targetDecomposeParDict",
+            runTimeTarget
         );
 
         List<boundBox> bbsTarget(nProcsTarget);
diff --git a/applications/utilities/preProcessing/mapFieldsPar/MapLagrangianFields.H b/applications/utilities/preProcessing/mapFieldsPar/MapLagrangianFields.H
index e045522145..d1a010ce21 100644
--- a/applications/utilities/preProcessing/mapFieldsPar/MapLagrangianFields.H
+++ b/applications/utilities/preProcessing/mapFieldsPar/MapLagrangianFields.H
@@ -104,6 +104,8 @@ void MapLagrangianFields
             Info<< "    mapping lagrangian fieldField " << fieldName << endl;
 
             // Read field (does not need mesh)
+            // Note: some fieldFields are 0 size (e.g. collision records) if
+            //       not used
             IOField<Field<Type> > fieldSource(*fieldIter());
 
             // Map - use CompactIOField to automatically write in
@@ -120,12 +122,22 @@ void MapLagrangianFields
                     IOobject::NO_WRITE,
                     false
                 ),
-                addParticles.size()
+                min(fieldSource.size(), addParticles.size()) // handle 0 size
             );
 
-            forAll(addParticles, i)
+            if (fieldSource.size())
             {
-                fieldTarget[i] = fieldSource[addParticles[i]];
+                forAll(addParticles, i)
+                {
+                    fieldTarget[i] = fieldSource[addParticles[i]];
+                }
+            }
+            else if (cloud::debug)
+            {
+                Pout<< "Not mapping " << fieldName << " since source size "
+                    << fieldSource.size() << " different to"
+                    << " cloud size " << addParticles.size()
+                    << endl;
             }
 
             // Write field
@@ -139,8 +151,9 @@ void MapLagrangianFields
 
         forAllIter(IOobjectList, fieldFields, fieldIter)
         {
-            Info<< "    mapping lagrangian fieldField "
-                << fieldIter()->name() << endl;
+            const word& fieldName = fieldIter()->name();
+
+            Info<< "    mapping lagrangian fieldField " << fieldName << endl;
 
             // Read field (does not need mesh)
             CompactIOField<Field<Type>, Type> fieldSource(*fieldIter());
@@ -150,7 +163,7 @@ void MapLagrangianFields
             (
                 IOobject
                 (
-                    fieldIter()->name(),
+                    fieldName,
                     meshTarget.time().timeName(),
                     cloud::prefix/cloudName,
                     meshTarget,
@@ -158,12 +171,22 @@ void MapLagrangianFields
                     IOobject::NO_WRITE,
                     false
                 ),
-                addParticles.size()
+                min(fieldSource.size(), addParticles.size()) // handle 0 size
             );
 
-            forAll(addParticles, i)
+            if (fieldSource.size())
             {
-                fieldTarget[i] = fieldSource[addParticles[i]];
+                forAll(addParticles, i)
+                {
+                    fieldTarget[i] = fieldSource[addParticles[i]];
+                }
+            }
+            else if (cloud::debug)
+            {
+                Pout<< "Not mapping " << fieldName << " since source size "
+                    << fieldSource.size() << " different to"
+                    << " cloud size " << addParticles.size()
+                    << endl;
             }
 
             // Write field
diff --git a/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C b/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C
index d2ed7da591..6196b5e5e5 100644
--- a/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C
+++ b/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C
@@ -110,9 +110,10 @@ void mapLagrangian(const meshToMesh& interp)
             cloud::prefix/cloudDirs[cloudI]
         );
 
-        IOobject* positionsPtr = objects.lookup(word("positions"));
+        bool foundPositions =
+            returnReduce(objects.found("positions"), orOp<bool>());;
 
-        if (positionsPtr)
+        if (foundPositions)
         {
             Info<< nl << "    processing cloud " << cloudDirs[cloudI] << endl;
 

From f0bc2d9faa5125b59e778cfeb443defcfd83f976 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 24 Nov 2015 17:35:18 +0000
Subject: [PATCH 8/8] ENH: parallel: added -decomposeParDict option to -
 foamyHexMesh - snappyHexMesh - decomposePar

---
 .../conformalVoronoiMesh/Make/options         |   2 +
 .../backgroundMeshDecomposition.C             |  31 +++--
 .../backgroundMeshDecomposition.H             |  13 +-
 .../backgroundMeshDecompositionI.H            |   7 --
 .../cellShapeControlMesh.C                    |  96 +++++++-------
 .../conformalVoronoiMesh.C                    |   6 +-
 .../conformalVoronoiMesh.H                    |  10 +-
 .../conformalVoronoiMeshCalcDualMesh.C        |   6 +-
 .../indexedCell/indexedCell.H                 |   4 +-
 .../foamyMesh/foamyHexMesh/Make/options       |   3 +-
 .../foamyMesh/foamyHexMesh/foamyHexMesh.C     |  15 ++-
 .../foamyMesh/foamyHexMesh/foamyHexMeshDict   |  17 ++-
 .../foamyHexMeshBackgroundMesh/Make/options   |   2 +
 .../foamyHexMeshBackgroundMesh.C              |  49 +++-----
 .../foamyMesh/foamyQuadMesh/foamyQuadMesh.C   |   2 +-
 .../generation/snappyHexMesh/Make/options     |   4 +
 .../generation/snappyHexMesh/snappyHexMesh.C  |  28 +++--
 .../decomposePar/decomposePar.C               |   4 +-
 .../preProcessing/mapFields/Make/options      |   3 +
 .../preProcessing/mapFields/mapFields.C       | 119 ++++++++++--------
 .../mapFieldsPar/MapLagrangianFields.H        |  41 ++++--
 .../mapFieldsPar/mapLagrangian.C              |   5 +-
 22 files changed, 265 insertions(+), 202 deletions(-)

diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/Make/options b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/Make/options
index 44c1d87758..bad58b9d58 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/Make/options
@@ -16,6 +16,7 @@ EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
     -I$(LIB_SRC)/edgeMesh/lnInclude \
     -I$(LIB_SRC)/fileFormats/lnInclude \
     -I$(LIB_SRC)/dynamicMesh/lnInclude \
@@ -29,6 +30,7 @@ EXE_INC = \
 LIB_LIBS = \
     ${CGAL_LIBS} \
     -lmeshTools \
+    -ldecompose \
     -ledgeMesh \
     -lfileFormats \
     -ltriSurface \
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
index 25e0f9a73e..b4428c3df8 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.C
@@ -30,6 +30,7 @@ License
 #include "Time.H"
 #include "Random.H"
 #include "pointConversion.H"
+#include "decompositionModel.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -149,7 +150,8 @@ void Foam::backgroundMeshDecomposition::initialRefinement()
 
     const conformationSurfaces& geometry = geometryToConformTo_;
 
-    decompositionMethod& decomposer = decomposerPtr_();
+    decompositionMethod& decomposer =
+        decompositionModel::New(mesh_).decomposer();
 
     volScalarField::InternalField& icellWeights = cellWeights.internalField();
 
@@ -782,7 +784,8 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
     const Time& runTime,
     Random& rndGen,
     const conformationSurfaces& geometryToConformTo,
-    const dictionary& coeffsDict
+    const dictionary& coeffsDict,
+    const fileName& decompDictFile
 )
 :
     runTime_(runTime),
@@ -810,18 +813,6 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
     bFTreePtr_(),
     allBackgroundMeshBounds_(Pstream::nProcs()),
     globalBackgroundBounds_(),
-    decomposeDict_
-    (
-        IOobject
-        (
-            "decomposeParDict",
-            runTime_.system(),
-            runTime_,
-            IOobject::MUST_READ_IF_MODIFIED,
-            IOobject::NO_WRITE
-        )
-    ),
-    decomposerPtr_(decompositionMethod::New(decomposeDict_)),
     mergeDist_(1e-6*mesh_.bounds().mag()),
     spanScale_(readScalar(coeffsDict.lookup("spanScale"))),
     minCellSizeLimit_
@@ -846,14 +837,17 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
             << exit(FatalError);
     }
 
-    if (!decomposerPtr_().parallelAware())
+    const decompositionMethod& decomposer =
+        decompositionModel::New(mesh_, decompDictFile).decomposer();
+
+    if (!decomposer.parallelAware())
     {
         FatalErrorIn
         (
             "void Foam::backgroundMeshDecomposition::initialRefinement() const"
         )
             << "You have selected decomposition method "
-            << decomposerPtr_().typeName
+            << decomposer.typeName
             << " which is not parallel aware." << endl
             << exit(FatalError);
     }
@@ -1008,7 +1002,10 @@ Foam::backgroundMeshDecomposition::distribute
             << endl;
     }
 
-    labelList newDecomp = decomposerPtr_().decompose
+    decompositionMethod& decomposer =
+        decompositionModel::New(mesh_).decomposer();
+
+    labelList newDecomp = decomposer.decompose
     (
         mesh_,
         mesh_.cellCentres(),
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
index 90dde14aa4..9d597732cb 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecomposition.H
@@ -127,13 +127,7 @@ class backgroundMeshDecomposition
         //  a point that is not found on any processor is in the domain at all
         treeBoundBox globalBackgroundBounds_;
 
-        //- Decomposition dictionary
-        IOdictionary decomposeDict_;
-
-        //- Decomposition method
-        autoPtr<decompositionMethod> decomposerPtr_;
-
-        //- Merge distance required by fvMeshDistribute
+        //- merge distance required by fvMeshDistribute
         scalar mergeDist_;
 
         //- Scale of a cell span vs cell size used to decide to refine a cell
@@ -204,7 +198,8 @@ public:
             const Time& runTime,
             Random& rndGen,
             const conformationSurfaces& geometryToConformTo,
-            const dictionary& coeffsDict
+            const dictionary& coeffsDict,
+            const fileName& decompDictFile = ""
         );
 
 
@@ -324,8 +319,6 @@ public:
             //- Return the point level of the underlying mesh
             inline const labelList& pointLevel() const;
 
-            //- Return the current decomposition method
-            inline const decompositionMethod& decomposer() const;
 };
 
 
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
index 199f4ae4ff..b8dcb5187f 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/backgroundMeshDecomposition/backgroundMeshDecompositionI.H
@@ -57,11 +57,4 @@ const Foam::labelList& Foam::backgroundMeshDecomposition::pointLevel() const
 }
 
 
-const Foam::decompositionMethod&
-Foam::backgroundMeshDecomposition::decomposer() const
-{
-    return decomposerPtr_();
-}
-
-
 // ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
index e7f7f5816c..9b6f71ffac 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
@@ -392,61 +392,59 @@ Foam::cellShapeControlMesh::cellShapeControlMesh(const Time& runTime)
 
         if (mesh.nPoints() == this->vertexCount())
         {
-            pointScalarField sizes
+            IOobject io
             (
-                IOobject
-                (
-                    "sizes",
-                    runTime.timeName(),
-                    meshSubDir,
-                    runTime,
-                    IOobject::READ_IF_PRESENT,
-                    IOobject::NO_WRITE,
-                    false
-                ),
-                pointMesh::New(mesh)
+                "sizes",
+                runTime.timeName(),
+                meshSubDir,
+                runTime,
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE,
+                false
             );
 
-            triadIOField alignments
-            (
-                IOobject
-                (
-                    "alignments",
-                    mesh.time().timeName(),
-                    meshSubDir,
-                    mesh.time(),
-                    IOobject::READ_IF_PRESENT,
-                    IOobject::AUTO_WRITE,
-                    false
-                )
-            );
-
-            if
-            (
-                sizes.size() == this->vertexCount()
-             && alignments.size() == this->vertexCount()
-            )
+            if (io.headerOk())
             {
-                for
+                pointScalarField sizes(io, pointMesh::New(mesh));
+
+                triadIOField alignments
                 (
-                    Finite_vertices_iterator vit = finite_vertices_begin();
-                    vit != finite_vertices_end();
-                    ++vit
-                )
+                    IOobject
+                    (
+                        "alignments",
+                        mesh.time().timeName(),
+                        meshSubDir,
+                        mesh.time(),
+                        IOobject::MUST_READ,
+                        IOobject::NO_WRITE,
+                        false
+                    )
+                );
+
+                if (alignments.size() == this->vertexCount())
                 {
-                    vit->targetCellSize() = sizes[vit->index()];
-                    vit->alignment() = alignments[vit->index()];
+                    for
+                    (
+                        Finite_vertices_iterator vit = finite_vertices_begin();
+                        vit != finite_vertices_end();
+                        ++vit
+                    )
+                    {
+                        vit->targetCellSize() = sizes[vit->index()];
+                        vit->alignment() = alignments[vit->index()];
+                    }
+                }
+                else
+                {
+                    FatalErrorIn
+                    (
+                        "Foam::cellShapeControlMesh::cellShapeControlMesh"
+                        "(const Time&)"
+                    )   << "Cell alignments point field " << alignments.size()
+                        << " is not the same size as the number of vertices"
+                        << " in the mesh " << this->vertexCount()
+                        << abort(FatalError);
                 }
-            }
-            else
-            {
-                FatalErrorIn
-                (
-                    "Foam::cellShapeControlMesh::cellShapeControlMesh"
-                    "(const Time&)"
-                )   << "Cell size point field is not the same size as the "
-                    << "mesh."
-                    << abort(FatalError);
             }
         }
     }
@@ -672,7 +670,7 @@ void Foam::cellShapeControlMesh::write() const
             IOobject::AUTO_WRITE
         ),
         pointMesh::New(mesh),
-        scalar(0)
+        dimensionedScalar("zero", dimLength, scalar(0))
     );
 
     triadIOField alignments
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
index 9cc9862129..a10f19be22 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
@@ -837,7 +837,8 @@ bool Foam::conformalVoronoiMesh::ownerAndNeighbour
 Foam::conformalVoronoiMesh::conformalVoronoiMesh
 (
     const Time& runTime,
-    const dictionary& foamyHexMeshDict
+    const dictionary& foamyHexMeshDict,
+    const fileName& decompDictFile
 )
 :
     DistributedDelaunayMesh<Delaunay>(runTime),
@@ -876,7 +877,8 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
             foamyHexMeshControls().foamyHexMeshDict().subDict
             (
                 "backgroundMeshDecomposition"
-            )
+            ),
+            decompDictFile
         )
       : NULL
     ),
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
index f6d44489b7..59fbba90a9 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
@@ -163,10 +163,8 @@ private:
         //- Limiting bound box before infinity begins
         treeBoundBox limitBounds_;
 
-        //-
         mutable pointPairs<Delaunay> ptPairs_;
 
-        //-
         featurePointConformer ftPtConformer_;
 
         //- Search tree for edge point locations
@@ -546,7 +544,7 @@ private:
         ) const;
 
         //- Check if a location is in the exclusion range of an existing feature
-        //- Edge conformation location
+        //  edge conformation location
         bool nearFeatureEdgeLocation
         (
             const pointIndexHit& pHit,
@@ -730,8 +728,7 @@ private:
 
         label classifyBoundaryPoint(Cell_handle cit) const;
 
-        //- Index all of the the Delaunay cells and calculate their
-        //- Dual points
+        //- Index all of the the Delaunay cells and calculate their dual points
         void indexDualVertices
         (
             pointField& pts,
@@ -874,7 +871,8 @@ public:
         conformalVoronoiMesh
         (
             const Time& runTime,
-            const dictionary& foamyHexMeshDict
+            const dictionary& foamyHexMeshDict,
+            const fileName& decompDictFile = ""
         );
 
 
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
index e88ee94369..85c49a9ccd 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
@@ -1292,9 +1292,9 @@ void Foam::conformalVoronoiMesh::indexDualVertices
         }
     }
 
-    OBJstream snapping1("snapToSurface1.obj");
-    OBJstream snapping2("snapToSurface2.obj");
-    OFstream tetToSnapTo("tetsToSnapTo.obj");
+    //OBJstream snapping1("snapToSurface1.obj");
+    //OBJstream snapping2("snapToSurface2.obj");
+    //OFstream tetToSnapTo("tetsToSnapTo.obj");
 
     for
     (
diff --git a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
index a2eeb21c99..7a19a25ba2 100644
--- a/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
+++ b/applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
@@ -88,15 +88,13 @@ class indexedCell
     // Private data
 
         //- The index for this Delaunay tetrahedral cell.  Type information is
-        //- Also carried:
-
+        //  also carried:
         //  ctFar : the dual point of this cell does not form part of the
         //          internal or boundary of the dual mesh
         //  >=0   : the (local) index of an internal or boundary dual point,
         //           not on a processor face
         //  < 0 && > ctFar : the (global) index of a dual point on a processor
         //                   face
-
         Foam::label index_;
 
         //- The number of times that this Delaunay cell has been limited
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/Make/options b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/Make/options
index b3f2e7dd7c..7aa621f26e 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/Make/options
@@ -30,7 +30,8 @@ EXE_LIBS = \
     -lconformalVoronoiMesh \
     -lmeshTools \
     -ldecompositionMethods \
-    -L$(FOAM_LIBBIN)/dummy -lptscotchDecomp \
+    -ldecompose \
+    -L$(FOAM_LIBBIN)/dummy -lptscotchDecomp -lscotchDecomp \
     -ledgeMesh \
     -lfileFormats \
     -ltriSurface \
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMesh.C b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMesh.C
index 6e68f1c80e..c83800e031 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMesh.C
@@ -62,6 +62,17 @@ int main(int argc, char *argv[])
     const bool checkGeometry = args.optionFound("checkGeometry");
     const bool conformationOnly = args.optionFound("conformationOnly");
 
+    // Allow override of decomposeParDict location
+    fileName decompDictFile;
+    if (args.optionReadIfPresent("decomposeParDict", decompDictFile))
+    {
+        if (isDir(decompDictFile))
+        {
+            decompDictFile = decompDictFile / "decomposeParDict";
+        }
+    }
+
+
     IOdictionary foamyHexMeshDict
     (
         IOobject
@@ -114,7 +125,7 @@ int main(int argc, char *argv[])
 
     Info<< "Create mesh for time = " << runTime.timeName() << nl << endl;
 
-    conformalVoronoiMesh mesh(runTime, foamyHexMeshDict);
+    conformalVoronoiMesh mesh(runTime, foamyHexMeshDict, decompDictFile);
 
 
     if (conformationOnly)
@@ -145,7 +156,7 @@ int main(int argc, char *argv[])
     }
 
 
-    Info<< nl << "End" << nl << endl;
+    Info<< "\nEnd\n" << endl;
 
     return 0;
 }
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMeshDict b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMeshDict
index deca4e0edb..78f22b4e49 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMeshDict
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMesh/foamyHexMeshDict
@@ -271,11 +271,18 @@ motionControl
             forceInitialPointInsertion  on;
             priority                    1;
             mode                        inside;
+
+            // Cell size at surface
             surfaceCellSizeFunction     uniformValue;
             uniformValueCoeffs
             {
                 surfaceCellSizeCoeff    0.5;
             }
+
+            // Cell size inside domain by having a region of thickness
+            // surfaceOffsetaround the surface with the surface cell size
+            // (so constant) and then down to distanceCellSize over a distance
+            // of linearDistance.
             cellSizeFunction            surfaceOffsetLinearDistance;
             surfaceOffsetLinearDistanceCoeffs
             {
@@ -375,9 +382,17 @@ polyMeshFiltering
     // Filter small and sliver faces
     filterFaces                             off;
 
-    // Write the underlying Delaunay tet mesh at output time
+    // Write the underlying Delaunay tet mesh (at output time)
     writeTetDualMesh                        false;
 
+    // Write the Delaunay tet mesh used for interpolating cell size and
+    // alignment (at output time)
+    writeCellShapeControlMesh           true;
+
+    // Write the hex/split-hex mesh used for parallel load balancing
+    // (at output time)
+    writeBackgroundMeshDecomposition    true;
+
     // Upper limit on the size of faces to be filtered.
     // fraction of the local target cell size
     filterSizeCoeff                         0.2;
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/Make/options b/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/Make/options
index 07f6800de6..00d889c6c6 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/Make/options
@@ -9,6 +9,7 @@ EXE_INC = \
     ${CGAL_INC} \
     -I../conformalVoronoiMesh/lnInclude \
     -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
     -I$(LIB_SRC)/edgeMesh/lnInclude \
     -I$(LIB_SRC)/mesh/autoMesh/lnInclude \
     -I$(LIB_SRC)/triSurface/lnInclude \
@@ -26,6 +27,7 @@ EXE_LIBS = \
     -lgmp \
     -lconformalVoronoiMesh \
     -ldecompositionMethods /* -L$(FOAM_LIBBIN)/dummy -lscotchDecomp */ \
+    -ldecompose \
     -ledgeMesh \
     -ltriSurface \
     -lmeshTools \
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C b/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C
index fc99f0fd05..62be3f1efe 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C
@@ -44,6 +44,7 @@ Description
 #include "isoSurfaceCell.H"
 #include "vtkSurfaceWriter.H"
 #include "syncTools.H"
+#include "decompositionModel.H"
 
 using namespace Foam;
 
@@ -467,7 +468,7 @@ int main(int argc, char *argv[])
 
         // Determine the number of cells in each direction.
         const vector span = bb.span();
-        vector nScalarCells = span/cellShapeControls().defaultCellSize();
+        vector nScalarCells = span/cellShapeControls.defaultCellSize();
 
         // Calculate initial cell size to be a little bit smaller than the
         // defaultCellSize to avoid initial refinement triggering.
@@ -521,28 +522,21 @@ int main(int argc, char *argv[])
         Info<< "Loaded mesh:" << endl;
         printMeshData(mesh);
 
-        // Allocate a decomposer
-        IOdictionary decompositionDict
-        (
-            IOobject
-            (
-                "decomposeParDict",
-                runTime.system(),
-                mesh,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
-        );
+        // Allow override of decomposeParDict location
+        fileName decompDictFile;
+        if (args.optionReadIfPresent("decomposeParDict", decompDictFile))
+        {
+            if (isDir(decompDictFile))
+            {
+                decompDictFile = decompDictFile / "decomposeParDict";
+            }
+        }
 
-        autoPtr<decompositionMethod> decomposer
+        labelList decomp = decompositionModel::New
         (
-            decompositionMethod::New
-            (
-                decompositionDict
-            )
-        );
-
-        labelList decomp = decomposer().decompose(mesh, mesh.cellCentres());
+            mesh,
+            decompDictFile
+        ).decomposer().decompose(mesh, mesh.cellCentres());
 
         // Global matching tolerance
         const scalar tolDim = getMergeDistance
@@ -574,18 +568,15 @@ int main(int argc, char *argv[])
     Info<< "Refining backgroud mesh according to cell size specification" << nl
         << endl;
 
+    const dictionary& backgroundMeshDict =
+        foamyHexMeshDict.subDict("backgroundMeshDecomposition");
+
     backgroundMeshDecomposition backgroundMesh
     (
-        1.0,    //spanScale,ratio of poly cell size v.s. hex cell size
-        0.0,    //minCellSizeLimit
-        0,      //minLevels
-        4,      //volRes, check multiple points per cell
-        20.0,   //maxCellWeightCoeff
         runTime,
-        geometryToConformTo,
-        cellShapeControls(),
         rndGen,
-        foamyHexMeshDict
+        geometryToConformTo,
+        backgroundMeshDict
     );
 
     if (writeMesh)
diff --git a/applications/utilities/mesh/generation/foamyMesh/foamyQuadMesh/foamyQuadMesh.C b/applications/utilities/mesh/generation/foamyMesh/foamyQuadMesh/foamyQuadMesh.C
index 066bc8205c..997726bba6 100644
--- a/applications/utilities/mesh/generation/foamyMesh/foamyQuadMesh/foamyQuadMesh.C
+++ b/applications/utilities/mesh/generation/foamyMesh/foamyQuadMesh/foamyQuadMesh.C
@@ -221,7 +221,7 @@ int main(int argc, char *argv[])
     Info<< "Finished extruding in = "
         << runTime.cpuTimeIncrement() << " s." << endl;
 
-    Info<< nl << "End\n" << endl;
+    Info<< "\nEnd\n" << endl;
 
     return 0;
 }
diff --git a/applications/utilities/mesh/generation/snappyHexMesh/Make/options b/applications/utilities/mesh/generation/snappyHexMesh/Make/options
index 219e3d419c..d62b82a239 100644
--- a/applications/utilities/mesh/generation/snappyHexMesh/Make/options
+++ b/applications/utilities/mesh/generation/snappyHexMesh/Make/options
@@ -8,14 +8,18 @@ EXE_INC = \
     -I$(LIB_SRC)/surfMesh/lnInclude \
     -I$(LIB_SRC)/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/edgeMesh/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 
 EXE_LIBS = \
     -lfiniteVolume \
     -ldecompositionMethods \
     -L$(FOAM_LIBBIN)/dummy -lptscotchDecomp \
+    /* note: scotch < 6.0 does not like both scotch and ptscotch together */ \
+    -lscotchDecomp \
     -lmeshTools \
     -lsurfMesh \
     -lfileFormats \
     -ldynamicMesh \
+    -ldecompose \
     -lautoMesh
diff --git a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
index 7b4d8577c4..15bd5e965f 100644
--- a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
+++ b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+     \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -57,6 +57,7 @@ Description
 #include "MeshedSurface.H"
 #include "globalIndex.H"
 #include "IOmanip.H"
+#include "decompositionModel.H"
 
 using namespace Foam;
 
@@ -819,15 +820,28 @@ int main(int argc, char *argv[])
     {
         if (Pstream::parRun())
         {
+            fileName decompDictFile;
+            if (args.optionReadIfPresent("decomposeParDict", decompDictFile))
+            {
+                if (isDir(decompDictFile))
+                {
+                    decompDictFile = decompDictFile/"decomposeParDict";
+                }
+            }
+
             decomposeDict = IOdictionary
             (
-                IOobject
+                decompositionModel::selectIO
                 (
-                    "decomposeParDict",
-                    runTime.system(),
-                    mesh,
-                    IOobject::MUST_READ_IF_MODIFIED,
-                    IOobject::NO_WRITE
+                    IOobject
+                    (
+                        "decomposeParDict",
+                        runTime.system(),
+                        mesh,
+                        IOobject::MUST_READ_IF_MODIFIED,
+                        IOobject::NO_WRITE
+                    ),
+                    decompDictFile
                 )
             );
         }
diff --git a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
index 63303d9f8c..9c62a12837 100644
--- a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
+++ b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+     \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -143,7 +143,7 @@ int main(int argc, char *argv[])
     );
 
     argList::noParallel();
-    Foam::argList::addOption
+    argList::addOption
     (
         "decomposeParDict",
         "file",
diff --git a/applications/utilities/preProcessing/mapFields/Make/options b/applications/utilities/preProcessing/mapFields/Make/options
index 7bd964094e..383ac3c12c 100644
--- a/applications/utilities/preProcessing/mapFields/Make/options
+++ b/applications/utilities/preProcessing/mapFields/Make/options
@@ -2,6 +2,8 @@ EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/lagrangian/basic/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude
 
 EXE_LIBS = \
@@ -9,4 +11,5 @@ EXE_LIBS = \
     -lmeshTools \
     -llagrangian \
     -lfiniteVolume \
+    -ldecompose \
     -lgenericPatchFields
diff --git a/applications/utilities/preProcessing/mapFields/mapFields.C b/applications/utilities/preProcessing/mapFields/mapFields.C
index 929776d9df..6f72053a08 100644
--- a/applications/utilities/preProcessing/mapFields/mapFields.C
+++ b/applications/utilities/preProcessing/mapFields/mapFields.C
@@ -36,9 +36,48 @@ Description
 #include "meshToMesh0.H"
 #include "processorFvPatch.H"
 #include "MapMeshes.H"
+#include "decompositionModel.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+int readNumProcs
+(
+    const argList& args,
+    const word& optionName,
+    const Time& runTime
+)
+{
+    fileName dictFile;
+    if (args.optionReadIfPresent(optionName, dictFile))
+    {
+        if (isDir(dictFile))
+        {
+            dictFile = dictFile/"decomposeParDict";
+        }
+    }
+
+    return readInt
+    (
+        IOdictionary
+        (
+            decompositionModel::selectIO
+            (
+                IOobject
+                (
+                    "decomposeParDict",
+                    runTime.system(),
+                    runTime,
+                    IOobject::MUST_READ_IF_MODIFIED,
+                    IOobject::NO_WRITE,
+                    false
+                ),
+                dictFile
+            )
+        ).lookup("numberOfSubdomains")
+    );
+}
+
+
 void mapConsistentMesh
 (
     const fvMesh& meshSource,
@@ -225,6 +264,19 @@ int main(int argc, char *argv[])
         "subtract",
         "subtract mapped source from target"
     );
+    argList::addOption
+    (
+        "sourceDecomposeParDict",
+        "file",
+        "read decomposePar dictionary from specified location"
+    );
+    argList::addOption
+    (
+        "targetDecomposeParDict",
+        "file",
+        "read decomposePar dictionary from specified location"
+    );
+
 
     argList args(argc, argv);
 
@@ -320,20 +372,13 @@ int main(int argc, char *argv[])
 
     if (parallelSource && !parallelTarget)
     {
-        IOdictionary decompositionDict
+        int nProcs = readNumProcs
         (
-            IOobject
-            (
-                "decomposeParDict",
-                runTimeSource.system(),
-                runTimeSource,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
+            args,
+            "sourceDecomposeParDict",
+            runTimeSource
         );
 
-        int nProcs(readInt(decompositionDict.lookup("numberOfSubdomains")));
-
         Info<< "Create target mesh\n" << endl;
 
         fvMesh meshTarget
@@ -399,19 +444,13 @@ int main(int argc, char *argv[])
     }
     else if (!parallelSource && parallelTarget)
     {
-        IOdictionary decompositionDict
+        int nProcs = readNumProcs
         (
-            IOobject
-            (
-                "decomposeParDict",
-                runTimeTarget.system(),
-                runTimeTarget,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
+            args,
+            "targetDecomposeParDict",
+            runTimeTarget
         );
 
-        int nProcs(readInt(decompositionDict.lookup("numberOfSubdomains")));
 
         Info<< "Create source mesh\n" << endl;
 
@@ -478,39 +517,17 @@ int main(int argc, char *argv[])
     }
     else if (parallelSource && parallelTarget)
     {
-        IOdictionary decompositionDictSource
+        int nProcsSource = readNumProcs
         (
-            IOobject
-            (
-                "decomposeParDict",
-                runTimeSource.system(),
-                runTimeSource,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
+            args,
+            "sourceDecomposeParDict",
+            runTimeSource
         );
-
-        int nProcsSource
+        int nProcsTarget = readNumProcs
         (
-            readInt(decompositionDictSource.lookup("numberOfSubdomains"))
-        );
-
-
-        IOdictionary decompositionDictTarget
-        (
-            IOobject
-            (
-                "decomposeParDict",
-                runTimeTarget.system(),
-                runTimeTarget,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
-        );
-
-        int nProcsTarget
-        (
-            readInt(decompositionDictTarget.lookup("numberOfSubdomains"))
+            args,
+            "targetDecomposeParDict",
+            runTimeTarget
         );
 
         List<boundBox> bbsTarget(nProcsTarget);
diff --git a/applications/utilities/preProcessing/mapFieldsPar/MapLagrangianFields.H b/applications/utilities/preProcessing/mapFieldsPar/MapLagrangianFields.H
index e045522145..d1a010ce21 100644
--- a/applications/utilities/preProcessing/mapFieldsPar/MapLagrangianFields.H
+++ b/applications/utilities/preProcessing/mapFieldsPar/MapLagrangianFields.H
@@ -104,6 +104,8 @@ void MapLagrangianFields
             Info<< "    mapping lagrangian fieldField " << fieldName << endl;
 
             // Read field (does not need mesh)
+            // Note: some fieldFields are 0 size (e.g. collision records) if
+            //       not used
             IOField<Field<Type> > fieldSource(*fieldIter());
 
             // Map - use CompactIOField to automatically write in
@@ -120,12 +122,22 @@ void MapLagrangianFields
                     IOobject::NO_WRITE,
                     false
                 ),
-                addParticles.size()
+                min(fieldSource.size(), addParticles.size()) // handle 0 size
             );
 
-            forAll(addParticles, i)
+            if (fieldSource.size())
             {
-                fieldTarget[i] = fieldSource[addParticles[i]];
+                forAll(addParticles, i)
+                {
+                    fieldTarget[i] = fieldSource[addParticles[i]];
+                }
+            }
+            else if (cloud::debug)
+            {
+                Pout<< "Not mapping " << fieldName << " since source size "
+                    << fieldSource.size() << " different to"
+                    << " cloud size " << addParticles.size()
+                    << endl;
             }
 
             // Write field
@@ -139,8 +151,9 @@ void MapLagrangianFields
 
         forAllIter(IOobjectList, fieldFields, fieldIter)
         {
-            Info<< "    mapping lagrangian fieldField "
-                << fieldIter()->name() << endl;
+            const word& fieldName = fieldIter()->name();
+
+            Info<< "    mapping lagrangian fieldField " << fieldName << endl;
 
             // Read field (does not need mesh)
             CompactIOField<Field<Type>, Type> fieldSource(*fieldIter());
@@ -150,7 +163,7 @@ void MapLagrangianFields
             (
                 IOobject
                 (
-                    fieldIter()->name(),
+                    fieldName,
                     meshTarget.time().timeName(),
                     cloud::prefix/cloudName,
                     meshTarget,
@@ -158,12 +171,22 @@ void MapLagrangianFields
                     IOobject::NO_WRITE,
                     false
                 ),
-                addParticles.size()
+                min(fieldSource.size(), addParticles.size()) // handle 0 size
             );
 
-            forAll(addParticles, i)
+            if (fieldSource.size())
             {
-                fieldTarget[i] = fieldSource[addParticles[i]];
+                forAll(addParticles, i)
+                {
+                    fieldTarget[i] = fieldSource[addParticles[i]];
+                }
+            }
+            else if (cloud::debug)
+            {
+                Pout<< "Not mapping " << fieldName << " since source size "
+                    << fieldSource.size() << " different to"
+                    << " cloud size " << addParticles.size()
+                    << endl;
             }
 
             // Write field
diff --git a/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C b/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C
index d2ed7da591..6196b5e5e5 100644
--- a/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C
+++ b/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C
@@ -110,9 +110,10 @@ void mapLagrangian(const meshToMesh& interp)
             cloud::prefix/cloudDirs[cloudI]
         );
 
-        IOobject* positionsPtr = objects.lookup(word("positions"));
+        bool foundPositions =
+            returnReduce(objects.found("positions"), orOp<bool>());;
 
-        if (positionsPtr)
+        if (foundPositions)
         {
             Info<< nl << "    processing cloud " << cloudDirs[cloudI] << endl;