From d2bb8ce3373a91baeaa7d288f015c34b9d343acf Mon Sep 17 00:00:00 2001
From: mattijs <m.janssens@opencfd.co.uk>
Date: Tue, 3 Jun 2008 23:52:11 +0100
Subject: [PATCH] changed dictionary

---
 .../generation/snappyHexMesh/snappyHexMesh.C  |  65 +-
 .../snappyHexMesh/snappyHexMeshDict           | 314 ++++++++
 .../autoHexMeshDriver/autoHexMeshDriver.C     | 701 +++++++++++++-----
 .../autoHexMeshDriver/autoHexMeshDriver.H     |  79 +-
 .../autoHexMeshDriverLayers.C                 | 134 +---
 .../autoHexMeshDriver/autoHexMeshDriverSnap.C |  14 +-
 .../meshRefinement/meshRefinementMerge.C      |   2 +-
 .../refinementSurfaces/refinementSurfaces.C   | 185 ++++-
 .../refinementSurfaces/refinementSurfaces.H   |  30 +-
 9 files changed, 1160 insertions(+), 364 deletions(-)
 create mode 100644 applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict

diff --git a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
index afdba507673..5fc097dcc70 100644
--- a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
+++ b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
@@ -34,9 +34,6 @@ Description
 #include "Time.H"
 #include "fvMesh.H"
 #include "autoHexMeshDriver.H"
-#include "pointMesh.H"
-#include "motionSmoother.H"
-#include "mapDistributePolyMesh.H"
 
 using namespace Foam;
 
@@ -52,6 +49,18 @@ int main(int argc, char *argv[])
     Info<< "Read mesh in = "
         << runTime.cpuTimeIncrement() << " s" << endl;
 
+    // Read decomposePar dictionary
+    IOdictionary decomposeDict
+    (
+        IOobject
+        (
+            "decomposeParDict",
+            runTime.system(),
+            mesh,
+            IOobject::MUST_READ,
+            IOobject::NO_WRITE
+        )
+    );
 
     // Read meshing dictionary
     IOdictionary meshDict
@@ -66,24 +75,46 @@ int main(int argc, char *argv[])
        )
     );
 
-    // Read decomposePar dictionary
-    IOdictionary decomposeDict
+    // refinement parameters
+    const dictionary& refineDict = meshDict.subDict("refineDict");
+
+    // snap-to-surface parameters
+    const dictionary& snapDict = meshDict.subDict("snapDict");
+
+    // mesh motion and mesh quality parameters
+    const dictionary& motionDict = meshDict.subDict("motionDict");
+
+    // layer addition parameters
+    const dictionary& layerDict = meshDict.subDict("layerDict");
+
+
+    // Main meshing driver. Read surfaces. Determine initial intersections.
+    autoHexMeshDriver meshEngine
     (
-        IOobject
-        (
-            "decomposeParDict",
-            runTime.system(),
-            mesh,
-            IOobject::MUST_READ,
-            IOobject::NO_WRITE
-        )
+        mesh,
+        meshDict,       // global control parameters
+        refineDict,     // refinement parameters
+        decomposeDict
     );
 
-    // Main meshing driver. Read surfaces. Determine intersections.
-    autoHexMeshDriver meshEngine(mesh, meshDict, decomposeDict);
+    Switch wantRefine(meshDict.lookup("doRefine"));
+    Switch wantSnap(meshDict.lookup("doSnap"));
+    Switch wantLayers(meshDict.lookup("doLayers"));
+
+    if (wantRefine)
+    {
+        meshEngine.doRefine(refineDict, wantSnap);
+    }
+
+    if (wantSnap)
+    {
+        meshEngine.doSnap(snapDict, motionDict);
+    }
 
-    // Do all: refine, snap, add layers
-    meshEngine.doMesh();
+    if (wantLayers)
+    {
+        meshEngine.doLayers(layerDict, motionDict);
+    }
 
     Info<< "Finished meshing in = "
         << runTime.elapsedCpuTime() << " s." << endl;
diff --git a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict
new file mode 100644
index 00000000000..e9aa8f189d9
--- /dev/null
+++ b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict
@@ -0,0 +1,314 @@
+/*---------------------------------------------------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  1.0                                   |
+|   \\  /    A nd           | Web:      http://www.openfoam.org               |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+
+FoamFile
+{
+    version         2.0;
+    format          ascii;
+
+    root            "/home/penfold/mattijs/foam/mattijs2.1/run/icoFoam";
+    case            "cavity";
+    instance        "system";
+    local           "";
+
+    class           dictionary;
+    object          autoHexMeshDict;
+}
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// Which phases to run.
+doRefine true;
+doSnap   true;
+doLayers true; // includes autoMergeFaces
+
+
+// Whether to dump intermediate meshes and print lots
+// 1 : write mesh
+// 2 : write volScalarField with cellLevel for postprocessing
+// 4 : write current intersections as .obj files
+debug 0;
+
+refineDict
+{
+    // Which part to keep.
+    // NOTE: This point should never be on a face, always inside a cell, even
+    // after refinement.
+    keepPoints ((3 0.28 0.43));
+
+    // Whether to remove/split cells likely to give problems when snapping
+    handleSnapProblems on;
+
+    // Merge tolerance. Is fraction of overall bounding box of initial mesh
+    mergeTolerance 1E-6;
+
+    // While refining maximum number of cells per processor. This is basically
+    // the number of cells that fit on a processor. If you choose this too small
+    // it will do just more refinement iterations to obtain a similar mesh.
+    procCellLimit 1000000;
+
+
+    // Overall cell limit (approximately). Refinement will stop immediately
+    // upon reaching this number so a refinement level might not complete.
+    // Note that this is the number of cells before removing the part which
+    // is not 'visible' from the keepPoint. The final number of cells might actually
+    // be a lot less.
+    cellLimit 2000000;
+
+
+    // The surface refinement loop might spend lots of iterations refining just a
+    // few cells. This setting will cause refinement to stop if <= minimumRefine
+    // are selected for refinement. Note: it will at least do one iteration
+    // (unless the number of cells to refine is 0)
+    minimumRefine 0;
+
+
+
+
+    // Number of buffer layers between different levels.
+    // 1 means normal 2:1 refinement restriction, larger means slower
+    // refinement.
+    nBufferLayers 1;
+
+
+    // Feature Edge Refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~~
+
+    // External feature file. Read from constant/triSurface for now.
+    // Limitations:
+    // - either start or edge of any feature line has to be inside domain
+    //   and be visible from keepPoint on starting mesh.
+    // - refining cells containing features is separate phase before refining
+    //   based on surface.
+    // - no load balancing, no check for cellLimit is done while doing this.
+    features
+    (
+    //    {
+    //        file "someLine.eMesh";
+    //        level 2;
+    //    }
+    );
+
+
+    // Internal Mesh Refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // Specifies the areas where the refinement has to be a certain level.
+    // These surfaces have to be closed. Refinement is either inside
+    // (refineInside = true) or outside. (note:insideness or outsideness
+    // is with respect to far away)
+    // Note:that even using these the transition can never be faster than
+    // nBufferLayers!
+    // Note:the refinement level can never be higher than any of the surfaces
+    // they overlap with. See below for the surface refinement level specification.
+    refinementShells
+    (
+        {
+            //type triSurfaceMesh;
+            //name cube1x1x1.stl;
+            type searchableBox;
+            name box1x1x1;
+            min (2.5 -0.5 -0.5);
+            max (3.5 0.5 0.5);
+
+            level 4;
+            refineInside true;
+        }
+    );
+
+
+    // Surface based refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // Curvature. Cosine of angle between two neighbouring surface intersections.
+    // Only used if cell level > minLevel and < maxLevel.
+    curvature 0.5;
+
+
+    // Overall the refinement is according to the following rules:
+    // - if the cell-cell vector intersects a surface any cell that
+    //   is less refined than the minRefinementLevel of the surface gets refined.
+    // - else if the refinement level of the cell is between the
+    //   minRefinementLevel and maxRefinementLevel the cell gets refined if
+    //      - the normal of neighbouring surface intersections differ by more
+    //        than above curvature
+    //      - or if neighbouring surface intersections are on different surfaces or
+    //        different surface regions.
+
+
+    // surfaces
+    surfaces
+    (
+        {
+            name    sphere;
+            file    "sphere.stl";
+
+            // Surface wide refinement level
+            minRefinementLevel 1;
+            maxRefinementLevel 1;
+
+            // Layers
+            surfaceLayers 1;
+
+            // Region specific refinement level
+            regions
+            (
+                {
+                    name firstSolid;
+                    minRefinementLevel 3;
+                    maxRefinementLevel 3;
+                    surfaceLayers 2;
+                }
+                {
+                    name secondSolid;
+                    minRefinementLevel 1;
+                    maxRefinementLevel 1;
+                    surfaceLayers 1;
+                }
+            );
+        }
+    );
+}
+
+// For snapping
+snapDict
+{
+    //- Number of patch smoothing iterations before finding correspondence
+    //  to surface
+    nSmoothPatch 3;
+
+    //- Relative distance for points to be attracted by surface feature point
+    //  or edge. True distance is this factor times local
+    //  maximum edge length.
+    snapTol 4.0;
+
+    //- Whether to move internal mesh as well as boundary
+    smoothMesh true;
+
+    //- Number of mesh displacement smoothing iterations.
+    nSmoothDispl 30;
+
+    //- Maximum number of snapping relaxation iterations. Should stop
+    //  before upon reaching a correct mesh.
+    nSnap 5;
+}
+
+
+// For cell layers
+layerDict
+{
+    //- When not to extrude surface. 0 is flat surface, 90 is when two faces
+    //  make straight angle.
+    featureAngle 60;
+
+    //- Maximum number of snapping relaxation iterations. Should stop
+    //  before upon reaching a correct mesh.
+    nSnap 5;
+
+
+    //- Minimum thickness of cell layer. If for any reason layer cannot be
+    //  above minThickness do not add layer if thickness below minThickNess.
+    //  Relative to undistorted cell size
+    minThickness 0.25;
+
+    //- If points get not extruded do nGrow layers of connected faces that are
+    //  not grown. Is used to not do layers at all close to features.
+    nGrow 1;
+
+    // Expansion factor for layer mesh
+    expansionRatio 1.3;	
+
+    // Ratio of cell size in final added cell layer to cell size 
+    // outside layer 
+    finalLayerRatio 0.3;
+
+    // Number of smoothing iterations of surface normals 
+    nSmoothSurfaceNormals 1;
+
+    // Number of smoothing iterations of interior mesh movement direction  
+    nSmoothNormals 3;
+
+    // Smooth layer thickness over surface patches
+    nSmoothThickness 10;
+
+    // Stop layer growth on highly warped cells 
+    maxFaceThicknessRatio 0.5;
+
+    // Reduce layer growth where ratio thickness to medial 
+    // distance is large 
+    maxThicknessToMedialRatio 0.3;
+
+    // Angle used to pick up medial axis points
+    minMedianAxisAngle 130;
+
+    // Create buffer region for new layer terminations
+    nBufferCellsNoExtrude 0;
+
+    thickness 0.5;
+    nSmoothDispl 4;
+}
+
+
+// For mesh motion
+motionDict
+{
+    //
+    // Mesh Quality Parameters. Decide when mesh is good enough to stop
+    // smoothing.
+    //
+
+    //- Maximum non-orthogonality allowed. Set to 180 to disable.
+    maxNonOrtho 65;
+
+    //- Max skewness allowed. Set to <0 to disable.
+    maxBoundarySkewness 20;
+    maxInternalSkewness 4;
+
+    //- Max concaveness allowed. Is angle (in degrees) below which concavity
+    //  is allowed. 0 is straight face, <0 would be convex face.
+    //  Set to 180 to disable.
+    maxConcave 80;
+
+    //- Minimum projected area v.s. actual area. Set to -1 to disable.
+    minFlatness 0.5;
+
+    //- Minimum pyramid volume. Is absolute volume of cell pyramid.
+    //  Set to very negative number (e.g. -1E30) to disable.
+    minVol 1e-13;
+
+    //- Minimum face area. Set to <0 to disable.
+    minArea -1;
+
+    //- Minimum face twist. Set to <-1 to disable. dot product of face normal
+    //- and face centre triangles normal
+    minTwist 0.05;
+
+    //- minimum normalised cell determinant
+    //- 1 = hex, <= 0 = folded or flattened illegal cell
+    minDeterminant 0.001;
+
+    //- minFaceWeight (0 -> 0.5) 
+    minFaceWeight 0.05;
+
+    //- minVolRatio (0 -> 1)
+    minVolRatio 0.01;
+
+
+    //must be >0 for Fluent compatibility
+    minTriangleTwist -1;
+
+    // Advanced
+
+    //- Number of error distribution iterations
+    nSmoothScale 4;
+    //- amount to scale back displacement at error points
+    errorReduction 0.75;
+}
+
+
+// ************************************************************************* //
diff --git a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.C b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.C
index f01b31d70d3..6435d491d9b 100644
--- a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.C
+++ b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.C
@@ -52,9 +52,9 @@ defineTypeNameAndDebug(autoHexMeshDriver, 0);
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Check writing tolerance before doing any serious work
-Foam::scalar Foam::autoHexMeshDriver::getMergeDistance() const
+Foam::scalar Foam::autoHexMeshDriver::getMergeDistance(const scalar mergeTol)
+ const
 {
-    const scalar mergeTol = readScalar(dict_.lookup("mergeTolerance"));
     const boundBox& meshBb = mesh_.bounds();
     scalar mergeDist = mergeTol*mag(meshBb.max() - meshBb.min());
     scalar writeTol = std::pow
@@ -70,7 +70,7 @@ Foam::scalar Foam::autoHexMeshDriver::getMergeDistance() const
 
     if (mesh_.time().writeFormat() == IOstream::ASCII && mergeTol < writeTol)
     {
-        FatalErrorIn("autoHexMeshDriver::getMergeDistance() const")
+        FatalErrorIn("autoHexMeshDriver::getMergeDistance(const scalar) const")
             << "Your current settings specify ASCII writing with "
             << IOstream::defaultPrecision() << " digits precision." << endl
             << "Your merging tolerance (" << mergeTol << ") is finer than this."
@@ -306,9 +306,8 @@ Foam::autoHexMeshDriver::autoHexMeshDriver
     curvature_(readScalar(dict_.lookup("curvature"))),
     nBufferLayers_(readLabel(dict_.lookup("nBufferLayers"))),
     keepPoints_(dict_.lookup("keepPoints")),
-    mergeDist_(getMergeDistance())
+    mergeDist_(getMergeDistance(readScalar(dict_.lookup("mergeTolerance"))))
 {
-
     if (debug_ > 0)
     {
         meshRefinement::debug = debug_;
@@ -471,17 +470,34 @@ Foam::autoHexMeshDriver::autoHexMeshDriver
             {
                 if (nTrisPerRegion[i] > 0)
                 {
+                    label globalRegionI = surfaces().globalRegion(surfI, i);
+
+                    // Use optionally specified patch type and name
+                    word patchType = surfaces().patchType()[globalRegionI];
+                    if (patchType == "")
+                    {
+                        patchType = wallPolyPatch::typeName;
+                    }
+
+                    word patchName = surfaces().patchName()[globalRegionI];
+                    if (patchName == "")
+                    {
+                        patchName =
+                            surfaces().names()[surfI]
+                          + '_'
+                          + regions[i].name();
+                    }
+
                     label patchI = meshRefinement::addPatch
                     (
                         mesh,
-                        //s.searchableSurface::name() + '_' + regions[i].name(),
-                        surfaces().names()[surfI] + '_' + regions[i].name(),
-                        wallPolyPatch::typeName
+                        patchName,
+                        patchType
                     );
 
                     Info<< patchI << '\t' << regions[i].name() << nl;
 
-                    globalToPatch_[surfaces().globalRegion(surfI, i)] = patchI;
+                    globalToPatch_[globalRegionI] = patchI;
                 }
             }
 
@@ -593,15 +609,283 @@ Foam::autoHexMeshDriver::autoHexMeshDriver
 }
 
 
+// Construct from separate dictionaries.
+Foam::autoHexMeshDriver::autoHexMeshDriver
+(
+    fvMesh& mesh,
+    const dictionary& controlDict,
+    const dictionary& refineDict,
+    const dictionary& decomposeDict
+)
+:
+    mesh_(mesh),
+    dict_(controlDict),
+    debug_(readLabel(controlDict.lookup("debug"))),
+    maxGlobalCells_(readLabel(refineDict.lookup("cellLimit"))),
+    maxLocalCells_(readLabel(refineDict.lookup("procCellLimit"))),
+    minRefineCells_(readLabel(refineDict.lookup("minimumRefine"))),
+    curvature_(readScalar(refineDict.lookup("curvature"))),
+    nBufferLayers_(readLabel(refineDict.lookup("nBufferLayers"))),
+    keepPoints_(refineDict.lookup("keepPoints")),
+    mergeDist_
+    (
+        getMergeDistance(readScalar(refineDict.lookup("mergeTolerance")))
+    )
+{
+    if (debug_ > 0)
+    {
+        meshRefinement::debug = debug_;
+        autoHexMeshDriver::debug = debug_;
+    }
+
+    Info<< "Overall cell limit                         : " << maxGlobalCells_
+        << endl;
+    Info<< "Per processor cell limit                   : " << maxLocalCells_
+        << endl;
+    Info<< "Minimum number of cells to refine          : " << minRefineCells_
+        << endl;
+    Info<< "Curvature                                  : " << curvature_
+        << nl << endl;
+    Info<< "Layers between different refinement levels : " << nBufferLayers_
+        << endl;
+
+    // Check keepPoints are sensible
+    findCells(keepPoints_);
+
+
+    // Read refinement shells
+    // ~~~~~~~~~~~~~~~~~~~~~~
+
+    {
+        Info<< "Reading refinement shells." << endl;
+
+        PtrList<dictionary> shellDicts(refineDict.lookup("refinementShells"));
+
+        shells_.setSize(shellDicts.size());
+        shellLevels_.setSize(shellDicts.size());
+        shellRefineInside_.setSize(shellDicts.size());
+
+        forAll(shellDicts, i)
+        {
+            const dictionary& dict = shellDicts[i];
+
+            shells_.set
+            (
+                i,
+                searchableSurface::New
+                (
+                    dict.lookup("type"),
+                    dict.lookup("name"),
+                    mesh_.time(),
+                    dict
+                )
+            );
+            shellLevels_[i] = readLabel(dict.lookup("level"));
+            shellRefineInside_[i] = Switch(dict.lookup("refineInside"));
+
+            if (shellRefineInside_[i])
+            {
+                Info<< "Refinement level " << shellLevels_[i]
+                    << " for all cells inside " << shells_[i].name() << endl;
+            }
+            else
+            {
+                Info<< "Refinement level " << shellLevels_[i]
+                    << " for all cells outside " << shells_[i].name() << endl;
+            }
+        }
+
+        Info<< "Read refinement shells in = "
+            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+
+        // Orient shell surfaces before any searching is done.
+        Info<< "Orienting triSurface shells so point far away is outside."
+            << endl;
+        orientOutside(shells_);
+        Info<< "Oriented shells in = "
+            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+    }
+
+
+    // Read refinement surfaces
+    // ~~~~~~~~~~~~~~~~~~~~~~~~
+
+    {
+        Info<< "Reading surfaces and constructing search trees." << endl;
+
+        surfacesPtr_.reset
+        (
+            new refinementSurfaces
+            (
+                IOobject
+                (
+                    "",                                 // dummy name
+                    mesh_.time().constant(),            // directory
+                    "triSurface",                       // instance
+                    mesh_.time(),                       // registry
+                    IOobject::MUST_READ,
+                    IOobject::NO_WRITE
+                ),
+                refineDict.lookup("surfaces")
+            )
+        );
+        Info<< "Read surfaces in = "
+            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+
+        // Orient surfaces (if they're closed) before any searching is done.
+        Info<< "Orienting (closed) surfaces so keepPoint is outside." << endl;
+        forAll(surfaces(), i)
+        {
+            if (refinementSurfaces::isSurfaceClosed(surfaces()[i]))
+            {
+                refinementSurfaces::orientSurface
+                (
+                    keepPoints_[0],
+                    surfacesPtr_()[i]
+                );
+            }
+        }
+        Info<< "Oriented closed surfaces in = "
+            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+
+        Info<< "Setting refinement level of surface to be consistent"
+            << " with shells." << endl;
+        surfacesPtr_().setMinLevelFields
+        (
+            shells_,
+            shellLevels_,
+            shellRefineInside_
+        );
+        Info<< "Checked shell refinement in = "
+            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+    }
+
+    // Check faceZones are synchronised
+    checkCoupledFaceZones();
+
+
+    // Add all the surface regions as patches
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    {
+        Info<< nl
+            << "Adding patches for surface regions" << nl
+            << "----------------------------------" << nl
+            << endl;
+
+        // From global region number to mesh patch.
+        globalToPatch_.setSize(surfaces().nRegions(), -1);
+
+        Info<< "Patch\tRegion" << nl
+            << "-----\t------"
+            << endl;
+
+        forAll(surfaces(), surfI)
+        {
+            const triSurfaceMesh& s = surfaces()[surfI];
+
+            Info<< surfaces().names()[surfI] << ':' << nl << nl;
+
+            const geometricSurfacePatchList& regions = s.patches();
+
+            labelList nTrisPerRegion(surfaces().countRegions(s));
+
+            forAll(regions, i)
+            {
+                if (nTrisPerRegion[i] > 0)
+                {
+                    label patchI = meshRefinement::addPatch
+                    (
+                        mesh,
+                        //s.searchableSurface::name() + '_' + regions[i].name(),
+                        surfaces().names()[surfI] + '_' + regions[i].name(),
+                        wallPolyPatch::typeName
+                    );
+
+                    Info<< patchI << '\t' << regions[i].name() << nl;
+
+                    globalToPatch_[surfaces().globalRegion(surfI, i)] = patchI;
+                }
+            }
+
+            Info<< nl;
+        }
+        Info<< "Added patches in = "
+            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+    }
+
+    // Parallel
+    // ~~~~~~~~
+
+    {
+        // Decomposition
+        decomposerPtr_ = decompositionMethod::New
+        (
+            decomposeDict,
+            mesh_
+        );
+        decompositionMethod& decomposer = decomposerPtr_();
+
+
+        if (Pstream::parRun() && !decomposer.parallelAware())
+        {
+            FatalErrorIn("autoHexMeshDriver::autoHexMeshDriver(const IOobject&, fvMesh&)")
+                << "You have selected decomposition method "
+                << decomposer.typeName
+                << " which is not parallel aware." << endl
+                << "Please select one that is (parMetis, hierarchical)"
+                << exit(FatalError);
+        }
+
+        // Mesh distribution engine (uses tolerance to reconstruct meshes)
+        distributorPtr_.reset(new fvMeshDistribute(mesh_, mergeDist_));
+    }
+
+
+    // Refinement engine
+    // ~~~~~~~~~~~~~~~~~
+
+    {
+        Info<< nl
+            << "Determining initial surface intersections" << nl
+            << "-----------------------------------------" << nl
+            << endl;
+
+        // Main refinement engine
+        meshRefinerPtr_.reset
+        (
+            new meshRefinement
+            (
+                mesh,
+                mergeDist_,         // tolerance used in sorting coordinates
+                surfaces()
+            )
+        );
+        Info<< "Calculated surface intersections in = "
+            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+
+        // Some stats
+        meshRefinerPtr_().printMeshInfo(debug_, "Initial mesh");
+
+        meshRefinerPtr_().write
+        (
+            debug_&meshRefinement::OBJINTERSECTIONS,
+            mesh_.time().path()/mesh_.time().timeName()
+        );
+    }
+}
+
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 // Read explicit feature edges
-Foam::label Foam::autoHexMeshDriver::readFeatureEdges()
+Foam::label Foam::autoHexMeshDriver::readFeatureEdges
+(
+    const PtrList<dictionary>& featDicts
+)
 {
     Info<< "Reading external feature lines." << endl;
 
-    PtrList<dictionary> featDicts(dict_.lookup("features"));
-
     featureMeshes_.setSize(featDicts.size());
     featureLevels_.setSize(featDicts.size());
 
@@ -647,13 +931,13 @@ Foam::label Foam::autoHexMeshDriver::readFeatureEdges()
 
 Foam::label Foam::autoHexMeshDriver::featureEdgeRefine
 (
+    const PtrList<dictionary>& featDicts,
     const label maxIter,
     const label minRefine
 )
 {
     // Read explicit feature edges
-    readFeatureEdges();
-
+    readFeatureEdges(featDicts);
 
     meshRefinement& meshRefiner = meshRefinerPtr_();
 
@@ -1359,122 +1643,232 @@ void Foam::autoHexMeshDriver::writeMesh(const string& msg) const
 }
 
 
-void Foam::autoHexMeshDriver::doMesh()
-{
-    Switch doRefine(dict_.lookup("doRefine"));
-    Switch doSnap(dict_.lookup("doSnap"));
-    Switch doLayers(dict_.lookup("doLayers"));
 
-    Info<< "Do refinement : " << doRefine << nl
-        << "Do snapping   : " << doSnap << nl
-        << "Do layers     : " << doLayers << nl
+
+void Foam::autoHexMeshDriver::doRefine
+(
+    const dictionary& refineDict,
+    const bool prepareForSnapping
+)
+{
+    Info<< nl
+        << "Refinement phase" << nl
+        << "----------------" << nl
         << endl;
 
-    if (doRefine)
+    const_cast<Time&>(mesh_.time())++;
+
+    PtrList<dictionary> featDicts(refineDict.lookup("features"));
+
+    // Refine around feature edges
+    featureEdgeRefine
+    (
+        featDicts,
+        100,    // maxIter
+        0       // min cells to refine
+    );
+
+    // Refine based on surface
+    surfaceOnlyRefine
+    (
+        100     // maxIter
+    );
+
+    // Remove cells (a certain distance) beyond surface intersections
+    removeInsideCells
+    (
+        1       // nBufferLayers
+    );
+
+    // Internal mesh refinement
+    shellRefine
+    (
+        100    // maxIter
+    );
+
+    // Introduce baffles at surface intersections
+    baffleAndSplitMesh(prepareForSnapping);
+
+    // Mesh is at its finest. Do optional zoning.
+    zonify();
+
+    // Pull baffles apart
+    splitAndMergeBaffles(prepareForSnapping);
+
+    // Do something about cells with refined faces on the boundary
+    if (prepareForSnapping)
     {
-        Info<< nl
-            << "Refinement phase" << nl
-            << "----------------" << nl
-            << endl;
+        mergePatchFaces();
+    }
 
-        const_cast<Time&>(mesh_.time())++;
+    // Do final balancing. Keep zoned faces on one processor.
+    balance(true, false);
 
-        // Refine around feature edges
-        featureEdgeRefine
-        (
-            100,    // maxIter
-            0       // min cells to refine
-        );
+    // Write mesh
+    writeMesh("Refined mesh");
+}
 
-        // Refine based on surface
-        surfaceOnlyRefine
+
+void Foam::autoHexMeshDriver::doSnap
+(
+    const dictionary& snapDict,
+    const dictionary& motionDict
+)
+{
+    Info<< nl
+        << "Morphing phase" << nl
+        << "--------------" << nl
+        << endl;
+
+    const_cast<Time&>(mesh_.time())++;
+
+    // Get the labels of added patches.
+    labelList adaptPatchIDs(getSurfacePatches());
+
+    // Create baffles (pairs of faces that share the same points)
+    // Baffles stored as owner and neighbour face that have been created.
+    List<labelPair> baffles;
+    createZoneBaffles(baffles);
+
+    {
+        autoPtr<indirectPrimitivePatch> ppPtr
         (
-            100     // maxIter
+            meshRefinement::makePatch
+            (
+                mesh_,
+                adaptPatchIDs
+            )
         );
+        indirectPrimitivePatch& pp = ppPtr();
+
+        // Distance to attact to nearest feature on surface
+        const scalarField snapDist(calcSnapDistance(snapDict, pp));
+
 
-        // Remove cells (a certain distance) beyond surface intersections
-        removeInsideCells
+        // Construct iterative mesh mover.
+        Info<< "Constructing mesh displacer ..." << endl;
+        Info<< "Using mesh parameters " << motionDict << nl << endl;
+
+        pointMesh pMesh(mesh_);
+
+        motionSmoother meshMover
         (
-            1       // nBufferLayers
+            mesh_,
+            pp,
+            adaptPatchIDs,
+            meshRefinement::makeDisplacementField(pMesh, adaptPatchIDs),
+            motionDict
         );
 
-        // Internal mesh refinement
-        shellRefine
+
+        // Check initial mesh
+        Info<< "Checking initial mesh ..." << endl;
+        labelHashSet wrongFaces(mesh_.nFaces()/100);
+        motionSmoother::checkMesh(false, mesh_, motionDict, wrongFaces);
+        const label nInitErrors = returnReduce
         (
-            100    // maxIter
+            wrongFaces.size(),
+            sumOp<label>()
         );
 
-        // Introduce baffles at surface intersections
-        baffleAndSplitMesh(doSnap);
+        Info<< "Detected " << nInitErrors << " illegal faces"
+            << " (concave, zero area or negative cell pyramid volume)"
+            << endl;
 
-        // Mesh is at its finest. Do optional zoning.
-        zonify();
 
-        // Pull baffles apart
-        splitAndMergeBaffles(doSnap);
+        Info<< "Checked initial mesh in = "
+            << mesh_.time().cpuTimeIncrement() << " s\n" << nl << endl;
 
-        // Do something about cells with refined faces on the boundary
-        if (doSnap)
-        {
-            mergePatchFaces();
-        }
+        // Pre-smooth patch vertices (so before determining nearest)
+        preSmoothPatch(snapDict, nInitErrors, baffles, meshMover);
+
+        // Calculate displacement at every patch point. Insert into
+        // meshMover.
+        calcNearestSurface(snapDist, meshMover);
 
-        // Do final balancing. Keep zoned faces on one processor.
-        balance(true, false);
+        // Get smoothly varying internal displacement field.
+        smoothDisplacement(snapDict, meshMover);
 
-        // Write mesh
-        writeMesh("Refined mesh");
+        // Apply internal displacement to mesh.
+        scaleMesh(snapDict, nInitErrors, baffles, meshMover);
     }
 
+    // Merge any introduced baffles.
+    mergeZoneBaffles(baffles);
 
-    if (doSnap)
-    {
-        Info<< nl
-            << "Morphing phase" << nl
-            << "--------------" << nl
-            << endl;
+    // Write mesh.
+    writeMesh("Snapped mesh");
+}
 
-        const_cast<Time&>(mesh_.time())++;
 
-        // Get the labels of added patches.
-        labelList adaptPatchIDs(getSurfacePatches());
+void Foam::autoHexMeshDriver::doLayers
+(
+    const dictionary& shrinkDict,
+    const dictionary& motionDict
+)
+{
+    Info<< nl
+        << "Shrinking and layer addition phase" << nl
+        << "----------------------------------" << nl
+        << endl;
 
-        // Create baffles (pairs of faces that share the same points)
-        // Baffles stored as owner and neighbour face that have been created.
-        List<labelPair> baffles;
-        createZoneBaffles(baffles);
+    const_cast<Time&>(mesh_.time())++;
 
-        {
-            autoPtr<indirectPrimitivePatch> ppPtr
-            (
-                meshRefinement::makePatch
-                (
-                    mesh_,
-                    adaptPatchIDs
-                )
-            );
-            indirectPrimitivePatch& pp = ppPtr();
+    Info<< "Using mesh parameters " << motionDict << nl << endl;
+
+    // Merge coplanar boundary faces
+    mergePatchFacesUndo(shrinkDict, motionDict);
 
-            // Distance to attact to nearest feature on surface
-            const scalarField snapDist(calcSnapDistance(pp));
+    // Per global region the number of layers (0 if no layer)
+    const labelList& numLayers = surfaces().numLayers();
 
+    // Patches that need to get a layer
+    DynamicList<label> patchIDs(numLayers.size());
+    label nFacesWithLayers = 0;
+    forAll(numLayers, region)
+    {
+        if (numLayers[region] > 0 && globalToPatch()[region] != -1)
+        {
+            label patchI = globalToPatch()[region];
+            patchIDs.append(patchI);
+            nFacesWithLayers += mesh_.boundaryMesh()[patchI].size();
+        }
+    }
+    patchIDs.shrink();
+
+    if (returnReduce(nFacesWithLayers, sumOp<label>()) == 0)
+    {
+        Info<< nl << "No layers to generate ..." << endl;
+    }
+    else
+    {
+        autoPtr<indirectPrimitivePatch> ppPtr
+        (
+            meshRefinement::makePatch
+            (
+                mesh_,
+                patchIDs
+            )
+        );
+        indirectPrimitivePatch& pp = ppPtr();
 
-            // Construct iterative mesh mover.
-            Info<< "Constructing mesh displacer ..." << endl;
-            const dictionary& motionDict = dict_.subDict("motionDict");
-            Info<< "Using mesh parameters " << motionDict << nl << endl;
+        // Construct iterative mesh mover.
+        Info<< "Constructing mesh displacer ..." << endl;
 
+        {
             pointMesh pMesh(mesh_);
 
             motionSmoother meshMover
             (
                 mesh_,
                 pp,
-                adaptPatchIDs,
-                meshRefinement::makeDisplacementField(pMesh, adaptPatchIDs),
+                patchIDs,
+                meshRefinement::makeDisplacementField(pMesh, patchIDs),
                 motionDict
             );
 
+            // Check that outside of mesh is not multiply connected.
+            checkMeshManifold();
 
             // Check initial mesh
             Info<< "Checking initial mesh ..." << endl;
@@ -1490,119 +1884,46 @@ void Foam::autoHexMeshDriver::doMesh()
                 << " (concave, zero area or negative cell pyramid volume)"
                 << endl;
 
-
-            Info<< "Checked initial mesh in = "
-                << mesh_.time().cpuTimeIncrement() << " s\n" << nl << endl;
-
-            // Pre-smooth patch vertices (so before determining nearest)
-            preSmoothPatch(nInitErrors, baffles, meshMover);
-
-            // Calculate displacement at every patch point. Insert into
-            // meshMover.
-            calcNearestSurface(snapDist, meshMover);
-
-            // Get smoothly varying internal displacement field.
-            smoothDisplacement(meshMover);
-
-            // Apply internal displacement to mesh.
-            scaleMesh(nInitErrors, baffles, meshMover);
+            // Do all topo changes
+            addLayers(shrinkDict, motionDict, nInitErrors, meshMover);
         }
 
-        // Merge any introduced baffles.
-        mergeZoneBaffles(baffles);
-
         // Write mesh.
-        writeMesh("Snapped mesh");
+        writeMesh("Layer mesh");
     }
+}
 
 
-    if (doLayers)
-    {
-        Info<< nl
-            << "Shrinking and layer addition phase" << nl
-            << "----------------------------------" << nl
-            << endl;
-
-        const_cast<Time&>(mesh_.time())++;
-
-        // Merge coplanar boundary faces
-        mergePatchFacesUndo();
-
-        // Per global region the number of layers (0 if no layer)
-        labelList nLayers(readNumLayers());
-
-        // Patches that need to get a layer
-        DynamicList<label> patchIDs(nLayers.size());
-        label nFacesWithLayers = 0;
-        forAll(nLayers, region)
-        {
-            if (nLayers[region] > 0 && globalToPatch()[region] != -1)
-            {
-                label patchI = globalToPatch()[region];
-                patchIDs.append(patchI);
-                nFacesWithLayers += mesh_.boundaryMesh()[patchI].size();
-            }
-        }
-        patchIDs.shrink();
-
-        if (returnReduce(nFacesWithLayers, sumOp<label>()) == 0)
-        {
-            Info<< nl << "No layers to generate ..." << endl;
-        }
-        else
-        {
-            autoPtr<indirectPrimitivePatch> ppPtr
-            (
-                meshRefinement::makePatch
-                (
-                    mesh_,
-                    patchIDs
-                )
-            );
-            indirectPrimitivePatch& pp = ppPtr();
-
-            // Construct iterative mesh mover.
-            Info<< "Constructing mesh displacer ..." << endl;
-            const dictionary& motionDict = dict_.subDict("motionDict");
-            Info<< "Using mesh parameters " << motionDict << nl << endl;
-
-            {
-                pointMesh pMesh(mesh_);
-
-                motionSmoother meshMover
-                (
-                    mesh_,
-                    pp,
-                    patchIDs,
-                    meshRefinement::makeDisplacementField(pMesh, patchIDs),
-                    motionDict
-                );
+void Foam::autoHexMeshDriver::doMesh()
+{
+    Switch wantRefine(dict_.lookup("doRefine"));
+    Switch wantSnap(dict_.lookup("doSnap"));
+    Switch wantLayers(dict_.lookup("doLayers"));
 
-                // Check that outside of mesh is not multiply connected.
-                checkMeshManifold();
+    Info<< "Do refinement : " << wantRefine << nl
+        << "Do snapping   : " << wantSnap << nl
+        << "Do layers     : " << wantLayers << nl
+        << endl;
 
-                // Check initial mesh
-                Info<< "Checking initial mesh ..." << endl;
-                labelHashSet wrongFaces(mesh_.nFaces()/100);
-                motionSmoother::checkMesh(false, mesh_, motionDict, wrongFaces);
-                const label nInitErrors = returnReduce
-                (
-                    wrongFaces.size(),
-                    sumOp<label>()
-                );
+    if (wantRefine)
+    {
+        doRefine(dict_, wantSnap);
+    }
 
-                Info<< "Detected " << nInitErrors << " illegal faces"
-                    << " (concave, zero area or negative cell pyramid volume)"
-                    << endl;
+    if (wantSnap)
+    {
+        const dictionary& snapDict = dict_.subDict("snapDict");
+        const dictionary& motionDict = dict_.subDict("motionDict");
 
+        doSnap(snapDict, motionDict);
+    }
 
-                // Do all topo changes
-                addLayers(nInitErrors, meshMover);
-            }
+    if (wantLayers)
+    {
+        const dictionary& motionDict = dict_.subDict("motionDict");
+        const dictionary& shrinkDict = dict_.subDict("shrinkDict");
 
-            // Write mesh.
-            writeMesh("Layer mesh");
-        }
+        doLayers(shrinkDict, motionDict);
     }
 }
 
diff --git a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.H b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.H
index 3ff2513f2c8..fccb96a337d 100644
--- a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.H
+++ b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.H
@@ -135,7 +135,7 @@ class autoHexMeshDriver
         //- Reference to mesh
         fvMesh& mesh_;
 
-        //- Input dictionarys
+        //- Input dictionary
         const dictionary dict_;
 
         //- Debug level
@@ -164,7 +164,6 @@ class autoHexMeshDriver
 
 
 
-
         //- Explicit features
         PtrList<featureEdgeMesh> featureMeshes_;
         //- Per feature the refinement level
@@ -183,10 +182,6 @@ class autoHexMeshDriver
         //- Per refinement surface region the patch
         labelList globalToPatch_;
 
-        ////- Per refinement surface with a valid faceZone the cyclicpatch
-        ////  (or -1)
-        //labelList surfaceToCyclicPatch_;
-
         //- Mesh refinement engine
         autoPtr<meshRefinement> meshRefinerPtr_;
 
@@ -197,12 +192,13 @@ class autoHexMeshDriver
         autoPtr<fvMeshDistribute> distributorPtr_;
 
 
+
     // Private Member Functions
 
         // Refinement
 
-            //- Get merge tolerance. Check against writing tolerance.
-            scalar getMergeDistance() const;
+            //- Calculate merge distance. Check against writing tolerance.
+            scalar getMergeDistance(const scalar mergeTol) const;
 
             // Return per keeppoint -1 or the local cell label the point is in.
             // Guaranteed to be only on one processor.
@@ -211,7 +207,7 @@ class autoHexMeshDriver
             static void orientOutside(PtrList<searchableSurface>&);
 
             //- Read feature edges
-            label readFeatureEdges();
+            label readFeatureEdges(const PtrList<dictionary>&);
 
         // Snapping
 
@@ -375,7 +371,6 @@ class autoHexMeshDriver
                 //  layers per surface.
                 void setNumLayers
                 (
-                    const labelList& nLayers,
                     const labelList& patchIDs,
                     const indirectPrimitivePatch& pp,
                     pointField& patchDisp,
@@ -651,6 +646,16 @@ public:
             const dictionary& decomposeDict
         );
 
+        //- Alternative constructor from top-level controldictionary and
+        //  refinement specific dictionary
+        autoHexMeshDriver
+        (
+            fvMesh& mesh,
+            const dictionary& controlDict,
+            const dictionary& refineDict,
+            const dictionary& decomposeDict
+        );
+
 
     // Member Functions
 
@@ -684,6 +689,7 @@ public:
             //- Refine around explicit feature edges
             label featureEdgeRefine
             (
+                const PtrList<dictionary>& featDicts,
                 const label maxIter,
                 const label minRefine
             );
@@ -734,7 +740,11 @@ public:
             autoPtr<mapPolyMesh> mergeZoneBaffles(const List<labelPair>&);
 
             //- Calculate edge length per patch point.
-            scalarField calcSnapDistance(const indirectPrimitivePatch&) const;
+            scalarField calcSnapDistance
+            (
+                const dictionary& snapDict,
+                const indirectPrimitivePatch&
+            ) const;
 
             //- Get patches generated for surfaces.
             labelList getSurfacePatches() const;
@@ -743,6 +753,7 @@ public:
             //  of surface points (on castellated mesh) w.r.t. surface.
             void preSmoothPatch
             (
+                const dictionary& snapDict,
                 const label nInitErrors,
                 const List<labelPair>& baffles,
                 motionSmoother&
@@ -758,12 +769,17 @@ public:
             ) const;
 
             //- Smooth the displacement field to the internal.
-            void smoothDisplacement(motionSmoother&) const;
+            void smoothDisplacement
+            (
+                const dictionary& snapDict,
+                motionSmoother&
+            ) const;
 
             //- Do the hard work: move the mesh according to displacement,
             //  locally relax the displacement.
             void scaleMesh
             (
+                const dictionary& snapDict,
                 const label nInitErrors,
                 const List<labelPair>& baffles,
                 motionSmoother&
@@ -773,23 +789,50 @@ public:
         // Layers
 
             //- Merge patch faces on same cell.
-            void mergePatchFacesUndo();
-
-            //- Read layer per region
-            labelList readNumLayers() const;
+            void mergePatchFacesUndo
+            (
+                const dictionary& shrinkDict,
+                const dictionary& motionDict
+            );
 
             //- Check that mesh outside is not multiply connected.
             void checkMeshManifold() const;
 
             //- Add cell layers
-            void addLayers(const scalar nAllowableErrors, motionSmoother&);
+            void addLayers
+            (
+                const dictionary& shrinkDict,
+                const dictionary& motionDict,
+                const scalar nAllowableErrors,
+                motionSmoother&
+            );
 
 
         // Other
 
+            //- Do all refinement.
+            void doRefine
+            (
+                const dictionary& refineDict,
+                const bool prepareForSnapping
+            );
+
+            //- Do all snapping.
+            void doSnap
+            (
+                const dictionary& snapDict,
+                const dictionary& motionDict
+            );
+
+            //- Do alllayer addition.
+            void doLayers
+            (
+                const dictionary& shrinkDict,
+                const dictionary& motionDict
+            );
+
             //- Do all : refine, snap, layers
             void doMesh();
-
 };
 
 
diff --git a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriverLayers.C b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriverLayers.C
index 8ae25d7b655..b75439cb570 100644
--- a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriverLayers.C
+++ b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriverLayers.C
@@ -61,6 +61,30 @@ Foam::label Foam::autoHexMeshDriver::mergePatchFacesUndo
     // Patch face merging engine
     combineFaces faceCombiner(mesh_, true);
 
+    // Pick up all candidate cells on boundary
+    labelHashSet boundaryCells(mesh_.nFaces()-mesh_.nInternalFaces());
+
+    {
+        labelList patchIDs(meshRefinement::addedPatches(globalToPatch_));
+
+        const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+
+        forAll(patchIDs, i)
+        {
+            label patchI = patchIDs[i];
+
+            const polyPatch& patch = patches[patchI];
+
+            if (!patch.coupled())
+            {
+                forAll(patch, i)
+                {
+                    boundaryCells.insert(mesh_.faceOwner()[patch.start()+i]);
+                }
+            }
+        }
+    }
+
     // Get all sets of faces that can be merged
     labelListList allFaceSets
     (
@@ -68,7 +92,7 @@ Foam::label Foam::autoHexMeshDriver::mergePatchFacesUndo
         (
             minCos,
             concaveCos,
-            meshRefinement::addedPatches(globalToPatch_)
+            boundaryCells
         )
     );
 
@@ -1164,91 +1188,9 @@ void Foam::autoHexMeshDriver::handleWarpedFaces
 //}
 
 
-Foam::labelList Foam::autoHexMeshDriver::readNumLayers() const
-{
-    // Read dictionary information
-    // (could already be extracted in refinementSurfaces?)
-
-    const PtrList<dictionary>& surfaceDicts = dict_.lookup("surfaces");
-
-    labelList globalSurfLayers(surfaceDicts.size());
-    List<HashTable<label> > regionSurfLayers(surfaceDicts.size());
-
-    forAll(surfaceDicts, surfI)
-    {
-        const dictionary& dict = surfaceDicts[surfI];
-
-        globalSurfLayers[surfI] = readLabel(dict.lookup("surfaceLayers"));
-
-        if (dict.found("regions"))
-        {
-            // Per-region layer information
-
-            PtrList<dictionary> regionDicts(dict.lookup("regions"));
-
-            forAll(regionDicts, dictI)
-            {
-                const dictionary& regionDict = regionDicts[dictI];
-
-                const word regionName(regionDict.lookup("name"));
-
-                label nLayers = readLabel(regionDict.lookup("surfaceLayers"));
-
-                Info<< "    region " << regionName << ':'<< nl
-                    << "        surface layers:" << nLayers << nl;
-
-                regionSurfLayers[surfI].insert(regionName, nLayers);
-            }
-        }
-    }
-
-
-    // Transfer per surface/region information into global region info
-
-    labelList nLayers(surfaces().minLevel().size(), 0);
-
-    forAll(surfaces(), surfI)
-    {
-        const geometricSurfacePatchList& regions = surfaces()[surfI].patches();
-
-        forAll(regions, regionI)
-        {
-            label global = surfaces().globalRegion(surfI, regionI);
-
-            // Initialise to surface-wise layers
-            nLayers[global] = globalSurfLayers[surfI];
-
-            // Override with region specific data if available
-            HashTable<label>::const_iterator iter =
-                regionSurfLayers[surfI].find(regions[regionI].name());
-
-            if (iter != regionSurfLayers[surfI].end())
-            {
-                nLayers[global] = iter();
-            }
-
-            // Check
-            if (nLayers[global] < 0)
-            {
-                FatalErrorIn
-                (
-                    "autoHexMeshDriver::readNumLayers()"
-                )   << "Illegal number of layers " << nLayers[global]
-                    << " for surface "
-                    << surfaces().names()[surfI]
-                    << " region " << regions[regionI].name() << endl
-                    << exit(FatalError);
-            }
-        }
-    }
-    return nLayers;
-}
-
-
 // No extrusion on faces with differing number of layers for points
 void Foam::autoHexMeshDriver::setNumLayers
 (
-    const labelList& nLayers,
     const labelList& patchIDs,
     const indirectPrimitivePatch& pp,
     pointField& patchDisp,
@@ -1259,13 +1201,15 @@ void Foam::autoHexMeshDriver::setNumLayers
     Info<< nl << "Handling points with inconsistent layer specification ..."
         << endl;
 
+    const labelList& nSurfLayers = surfaces().numLayers();
+
     // Build map from patch to layers
-    Map<label> patchToNLayers(nLayers.size());
-    forAll(nLayers, region)
+    Map<label> patchToNLayers(nSurfLayers.size());
+    forAll(nSurfLayers, region)
     {
         if (globalToPatch_[region] != -1)
         {
-            patchToNLayers.insert(globalToPatch_[region], nLayers[region]);
+            patchToNLayers.insert(globalToPatch_[region], nSurfLayers[region]);
         }
     }
 
@@ -2363,10 +2307,12 @@ void Foam::autoHexMeshDriver::getLayerCellsFaces
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-void Foam::autoHexMeshDriver::mergePatchFacesUndo()
+void Foam::autoHexMeshDriver::mergePatchFacesUndo
+(
+    const dictionary& shrinkDict,
+    const dictionary& motionDict
+)
 {
-    const dictionary& shrinkDict = dict_.subDict("shrinkDict");
-
     const scalar featureAngle(readScalar(shrinkDict.lookup("featureAngle")));
     scalar minCos = Foam::cos(featureAngle*mathematicalConstant::pi/180.0);
 
@@ -2390,9 +2336,6 @@ void Foam::autoHexMeshDriver::mergePatchFacesUndo()
         << concaveAngle << " degrees (0=straight, 180=fully concave)" << nl
         << endl;
 
-
-    const dictionary& motionDict = dict_.subDict("motionDict");
-
     label nChanged = mergePatchFacesUndo(minCos, concaveCos, motionDict);
 
     nChanged += mergeEdgesUndo(minCos, motionDict);
@@ -2401,6 +2344,8 @@ void Foam::autoHexMeshDriver::mergePatchFacesUndo()
 
 void Foam::autoHexMeshDriver::addLayers
 (
+    const dictionary& shrinkDict,
+    const dictionary& motionDict,
     const scalar nAllowableErrors,
     motionSmoother& meshMover
 )
@@ -2408,8 +2353,6 @@ void Foam::autoHexMeshDriver::addLayers
     // Read some more dictionary settings
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
-    const dictionary& shrinkDict = dict_.subDict("shrinkDict");
-
     // Min thickness per cell
     const scalar relMinThickness(readScalar(shrinkDict.lookup("minThickness")));
 
@@ -2506,7 +2449,6 @@ void Foam::autoHexMeshDriver::addLayers
 
     setNumLayers
     (
-        readNumLayers(),
         meshMover.adaptPatchIDs(),
         pp,
 
@@ -2923,7 +2865,7 @@ void Foam::autoHexMeshDriver::addLayers
         label nTotChanged = checkAndUnmark
         (
             addLayer,
-            dict_.subDict("motionDict"),
+            motionDict,
             pp,
             newMesh,
 
diff --git a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriverSnap.C b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriverSnap.C
index c8b08b22d90..97abfc7ea7c 100644
--- a/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriverSnap.C
+++ b/src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriverSnap.C
@@ -755,10 +755,10 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::autoHexMeshDriver::mergeZoneBaffles
 
 Foam::scalarField Foam::autoHexMeshDriver::calcSnapDistance
 (
+    const dictionary& snapDict,
     const indirectPrimitivePatch& pp
 ) const
 {
-    const dictionary& snapDict = dict_.subDict("snapDict");
     // When to snap
     scalar snapTol(readScalar(snapDict.lookup("snapTol")));
 
@@ -826,12 +826,12 @@ Foam::labelList Foam::autoHexMeshDriver::getSurfacePatches() const
 
 void Foam::autoHexMeshDriver::preSmoothPatch
 (
+    const dictionary& snapDict,
     const label nInitErrors,
     const List<labelPair>& baffles,
     motionSmoother& meshMover
 ) const
 {
-    const dictionary& snapDict = dict_.subDict("snapDict");
     // Smoothing iterations
     label nSmoothPatch(readLabel(snapDict.lookup("nSmoothPatch")));
     // Snapping iterations
@@ -1076,15 +1076,17 @@ Foam::vectorField Foam::autoHexMeshDriver::calcNearestSurface
 }
 
 
-void Foam::autoHexMeshDriver::smoothDisplacement(motionSmoother& meshMover)
- const
+void Foam::autoHexMeshDriver::smoothDisplacement
+(
+    const dictionary& snapDict,
+    motionSmoother& meshMover
+) const
 {
     const pointMesh& pMesh = meshMover.pMesh();
     const indirectPrimitivePatch& pp = meshMover.patch();
 
     Info<< "Smoothing displacement ..." << endl;
 
-    const dictionary& snapDict = dict_.subDict("snapDict");
     // Smoothing iterations
     label nSmoothDisp(readLabel(snapDict.lookup("nSmoothDispl")));
 
@@ -1138,12 +1140,12 @@ void Foam::autoHexMeshDriver::smoothDisplacement(motionSmoother& meshMover)
 
 void Foam::autoHexMeshDriver::scaleMesh
 (
+    const dictionary& snapDict,
     const label nInitErrors,
     const List<labelPair>& baffles,
     motionSmoother& meshMover
 )
 {
-    const dictionary& snapDict = dict_.subDict("snapDict");
     // Snapping iterations
     label nSnap(readLabel(snapDict.lookup("nSnap")));
 
diff --git a/src/autoMesh/autoHexMesh/meshRefinement/meshRefinementMerge.C b/src/autoMesh/autoHexMesh/meshRefinement/meshRefinementMerge.C
index 0f489029af3..74ee5b4bf3b 100644
--- a/src/autoMesh/autoHexMesh/meshRefinement/meshRefinementMerge.C
+++ b/src/autoMesh/autoHexMesh/meshRefinement/meshRefinementMerge.C
@@ -57,7 +57,7 @@ Foam::label Foam::meshRefinement::mergePatchFaces
 
     const polyBoundaryMesh& patches = mesh_.boundaryMesh();
 
-    // Pick up all cells on boundary
+    // Pick up all candidate cells on boundary
     labelHashSet boundaryCells(mesh_.nFaces()-mesh_.nInternalFaces());
 
     forAll(patchIDs, i)
diff --git a/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C b/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
index 2ab52e0c4ee..6dbccc7c142 100644
--- a/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
+++ b/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
@@ -73,38 +73,41 @@ Foam::refinementSurfaces::refinementSurfaces
     List<HashTable<label> > regionMinLevel(surfaceDicts.size());
     List<HashTable<label> > regionMaxLevel(surfaceDicts.size());
 
-    forAll(surfaceDicts, i)
+    labelList globalSurfLayers(surfaceDicts.size());
+    List<HashTable<label> > regionSurfLayers(surfaceDicts.size());
+
+    wordList globalPatchType(surfaceDicts.size());
+    List<HashTable<word> > regionPatchType(surfaceDicts.size());
+    List<HashTable<word> > regionPatchName(surfaceDicts.size());
+
+    forAll(surfaceDicts, surfI)
     {
-        const dictionary& dict = surfaceDicts[i];
+        const dictionary& dict = surfaceDicts[surfI];
 
-        names_[i] = word(dict.lookup("name"));
+        names_[surfI] = word(dict.lookup("name"));
 
-        globalMinLevel[i] = readLabel(dict.lookup("minRefinementLevel"));
-        globalMaxLevel[i] = readLabel(dict.lookup("maxRefinementLevel"));
+        // Global refinement level
+        globalMinLevel[surfI] = readLabel(dict.lookup("minRefinementLevel"));
+        globalMaxLevel[surfI] = readLabel(dict.lookup("maxRefinementLevel"));
 
-        if
-        (
-            globalMinLevel[i] < 0
-         || globalMaxLevel[i] < globalMinLevel[i]
-        )
+        // Global number of layers
+        globalSurfLayers[surfI] = readLabel(dict.lookup("surfaceLayers"));
+
+        // Global zone names per surface
+        if (dict.found("faceZone"))
         {
-            FatalErrorIn
-            (
-                "refinementSurfaces::refinementSurfaces"
-                "(const IOobject&, const PtrList<dictionary>&)"
-            )   << "Illegal level specification for surface " << names_[i]
-                << " : minLevel:" << globalMinLevel[i]
-                << " maxLevel:" << globalMaxLevel[i]
-                << exit(FatalError);
+            dict.lookup("faceZone") >> faceZoneNames_[surfI];
+            dict.lookup("cellZone") >> cellZoneNames_[surfI];
+            dict.lookup("zoneInside") >> zoneInside_[surfI];
         }
 
-        if (dict.found("faceZone"))
+        // Global patch name per surface
+        if (dict.found("patchType"))
         {
-            dict.lookup("faceZone") >> faceZoneNames_[i];
-            dict.lookup("cellZone") >> cellZoneNames_[i];
-            dict.lookup("zoneInside") >> zoneInside_[i];
+            dict.lookup("patchType") >> globalPatchType[surfI];
         }
 
+
         if (dict.found("regions"))
         {
             PtrList<dictionary> regionDicts(dict.lookup("regions"));
@@ -117,12 +120,72 @@ Foam::refinementSurfaces::refinementSurfaces
                 label min = readLabel(regionDict.lookup("minRefinementLevel"));
                 label max = readLabel(regionDict.lookup("maxRefinementLevel"));
 
-                regionMinLevel[i].insert(regionName, min);
-                regionMaxLevel[i].insert(regionName, max);
+                regionMinLevel[surfI].insert(regionName, min);
+                regionMaxLevel[surfI].insert(regionName, max);
+
+                label nLayers = readLabel(regionDict.lookup("surfaceLayers"));
+                regionSurfLayers[surfI].insert(regionName, nLayers);
+
+                if (regionDict.found("patchType"))
+                {
+                    regionPatchType[surfI].insert
+                    (
+                        regionName,
+                        regionDict.lookup("patchType")
+                    );
+                    regionPatchName[surfI].insert
+                    (
+                        regionName,
+                        regionDict.lookup("patchName")
+                    );
+                }
             }
         }
     }
 
+
+    // Check for duplicate surface or region names
+    {
+        HashTable<label> surfaceNames(names_.size());
+
+        forAll(names_, surfI)
+        {
+            if (!surfaceNames.insert(names_[surfI], surfI))
+            {
+                FatalErrorIn
+                (
+                    "refinementSurfaces::refinementSurfaces"
+                    "(const IOobject&, const PtrList<dictionary>&)"
+                )   << "Duplicate surface name " << names_[surfI] << endl
+                    << "Previous occurrence of name at surface "
+                    << surfaceNames[names_[surfI]]
+                    << exit(FatalError);
+            }
+
+            // Check for duplicate region names
+            const geometricSurfacePatchList& patches =
+                operator[](surfI).patches();
+
+            HashTable<label> regionNames(patches.size());
+            forAll(patches, i)
+            {
+                if (!regionNames.insert(patches[i].name(), i))
+                {
+                    FatalErrorIn
+                    (
+                        "refinementSurfaces::refinementSurfaces"
+                        "(const IOobject&, const PtrList<dictionary>&)"
+                    )   << "Duplicate region name " << patches[i].name()
+                        << " on surface " << names_[surfI] << endl
+                        << "Previous occurrence of region at index "
+                        << regionNames[patches[i].name()]
+                        << exit(FatalError);
+                }
+            }
+        }
+    }
+
+
     // Calculate closedness
     forAll(closed_, surfI)
     {
@@ -151,38 +214,90 @@ Foam::refinementSurfaces::refinementSurfaces
         nRegions += operator[](surfI).patches().size();
     }
 
-    // From global region number to refinement level
+    // Rework surface specific information into information per global region
     minLevel_.setSize(nRegions);
     minLevel_ = 0;
     maxLevel_.setSize(nRegions);
     maxLevel_ = 0;
+    numLayers_.setSize(nRegions);
+    numLayers_ = 0;
+    patchName_.setSize(nRegions);
+    patchType_.setSize(nRegions);
 
     forAll(surfaceDicts, surfI)
     {
         const geometricSurfacePatchList& regions = operator[](surfI).patches();
 
+        // Initialise to global (i.e. per surface)
         forAll(regions, i)
         {
             minLevel_[regionOffset_[surfI] + i] = globalMinLevel[surfI];
             maxLevel_[regionOffset_[surfI] + i] = globalMaxLevel[surfI];
+            numLayers_[regionOffset_[surfI] + i] = globalSurfLayers[surfI];
+            patchType_[regionOffset_[surfI] + i] = globalPatchType[surfI];
+        }
+
+        // Get the region names
+        wordList regionNames(regions.size());
+        forAll(regions, regionI)
+        {
+            regionNames[regionI] = regions[regionI].name();
         }
 
+        // Overwrite with region specific information
         forAllConstIter(HashTable<label>, regionMinLevel[surfI], iter)
         {
-            // Find the patch
-            forAll(regions, regionI)
+            // Find the local region number.
+            label regionI = findIndex(regionNames, iter.key());
+
+            if (regionI == -1)
             {
-                if (regions[regionI].name() == iter.key())
-                {
-                    label globalRegionI = regionOffset_[surfI] + regionI;
+                FatalErrorIn
+                (
+                    "refinementSurfaces::refinementSurfaces"
+                    "(const IOobject&, const PtrList<dictionary>&)"
+                )   << "Cannot find region " << iter.key()
+                    << " in surface " << names_[surfI]
+                    << " which has regions " << regionNames
+                    << abort(FatalError);
+            }
+            
+            label globalRegionI = regionOffset_[surfI] + regionI;
 
-                    minLevel_[globalRegionI] = iter();
-                    maxLevel_[globalRegionI] =
-                        regionMaxLevel[surfI][iter.key()];
-                    break;
-                }
+            minLevel_[globalRegionI] = iter();
+            maxLevel_[globalRegionI] = regionMaxLevel[surfI][iter.key()];
+            numLayers_[globalRegionI] = regionSurfLayers[surfI][iter.key()];
+
+            // Check validity
+            if
+            (
+                minLevel_[globalRegionI] < 0
+             || maxLevel_[globalRegionI] < minLevel_[globalRegionI]
+             || numLayers_[globalRegionI] < 0
+            )
+            {
+                FatalErrorIn
+                (
+                    "refinementSurfaces::refinementSurfaces"
+                    "(const IOobject&, const PtrList<dictionary>&)"
+                )   << "Illegal level or layer specification for surface "
+                    << names_[surfI]
+                    << " : minLevel:" << minLevel_[globalRegionI]
+                    << " maxLevel:" << maxLevel_[globalRegionI]
+                    << " numLayers:" << numLayers_[globalRegionI]
+                    << exit(FatalError);
             }
         }
+
+        // Optional patch names and patch types
+        forAllConstIter(HashTable<word>, regionPatchName[surfI], iter)
+        {
+            label regionI = findIndex(regionNames, iter.key());
+            label globalRegionI = regionOffset_[surfI] + regionI;
+
+            patchName_[globalRegionI] = iter();
+            patchType_[globalRegionI] = regionPatchType[surfI][iter.key()];
+        }
     }
 }
 
diff --git a/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.H b/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.H
index 9a366d44018..9e0b720d680 100644
--- a/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.H
+++ b/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.H
@@ -60,7 +60,7 @@ class refinementSurfaces
 {
     // Private data
 
-        //- Name (word)
+        //- Surface name (word)
         wordList names_;
 
         //- Per surface whether is closed
@@ -86,6 +86,16 @@ class refinementSurfaces
         //- From global region number to refinement level
         labelList maxLevel_;
 
+        //- From global region number to added layers
+        labelList numLayers_;
+
+        //- From global region number to patch name
+        wordList patchName_;
+
+        //- From global region number to patch name
+        wordList patchType_;
+
+
         //- Per surface refinement level adapted for shells.
         PtrList<triSurfaceLabelField> minLevelFields_;
 
@@ -169,6 +179,24 @@ public:
                 return maxLevel_;
             }
 
+            //- From global region number to added layers
+            const labelList& numLayers() const
+            {
+                return numLayers_;
+            }
+
+            //- From global region number to patch name. Patchnames can be empty
+            const wordList& patchName() const
+            {
+                return patchName_;
+            }
+
+            //- From global region number to patch name
+            const wordList& patchType() const
+            {
+                return patchType_;
+            }
+
 
         // Helper
 
-- 
GitLab