Merge commit 'OpenCFD/master' into olesenm

applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict

@@ -316,6 +316,11 @@ meshQualityControls
     //must be >0 for Fluent compatibility
     minTriangleTwist -1;
 
+    //- if >0 : preserve single cells with all points on the surface if the
+    //  resulting volume after snapping is larger than minVolFraction times old
+    //  volume. If <0 : delete always.
+    minVolFraction 0.1;
+
 
     // Advanced
 

applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C

@@ -23,16 +23,13 @@ Foam::label Foam::checkGeometry(const polyMesh& mesh, const bool allGeometry)
     scalar minDistSqr = magSqr(1e-6 * globalBb.span());
 
     // Non-empty directions
-    const Vector<label> validDirs = (mesh.directions() + Vector<label>::one)/2;
+    const Vector<label> validDirs = (mesh.geometricD() + Vector<label>::one)/2;
+    Info<< "    Mesh (non-empty, non-wedge) directions " << validDirs << endl;
 
-    Info<< "    Mesh (non-empty) directions " << validDirs << endl;
+    const Vector<label> solDirs = (mesh.solutionD() + Vector<label>::one)/2;
+    Info<< "    Mesh (non-empty) directions " << solDirs << endl;
 
-    scalar nGeomDims = mesh.nGeometricD();
+    if (mesh.nGeometricD() < 3)
-
-    Info<< "    Mesh (non-empty, non-wedge) dimensions "
-        << nGeomDims << endl;
-
-    if (nGeomDims < 3)
     {
         pointSet nonAlignedPoints(mesh, "nonAlignedEdges", mesh.nPoints()/100);
 

applications/utilities/parallelProcessing/decomposePar/fvFieldDecomposer.C

@@ -59,7 +59,6 @@ processorVolPatchFieldDecomposer
     addressing_(addressingSlice.size()),
     weights_(addressingSlice.size())
 {
-    const scalarField& weights = mesh.weights().internalField();
     const labelList& own = mesh.faceOwner();
     const labelList& neighb = mesh.faceNeighbour();
 
@@ -72,15 +71,22 @@ processorVolPatchFieldDecomposer
         {
             // This is a regular face. it has been an internal face
             // of the original mesh and now it has become a face
-            // on the parallel boundary
+            // on the parallel boundary.
-            addressing_[i].setSize(2);
+            // Give face the value of the neighbour.
-            weights_[i].setSize(2);
 
+            addressing_[i].setSize(1);
+            weights_[i].setSize(1);
+            weights_[i][0] = 1.0;
+
+            if (addressingSlice[i] >= 0)
+            {
+                // I have the owner so use the neighbour value
+                addressing_[i][0] = neighb[ai];
+            }
+            else
+            {
                 addressing_[i][0] = own[ai];
-            addressing_[i][1] = neighb[ai];
+            }
-
-            weights_[i][0] = weights[ai];
-            weights_[i][1] = 1.0 - weights[ai];
         }
         else
         {
@@ -89,7 +95,7 @@ processorVolPatchFieldDecomposer
             // do the interpolation properly (I would need to look
             // up the different (face) list of data), so I will
             // just grab the value from the owner cell
-            // 
+
             addressing_[i].setSize(1);
             weights_[i].setSize(1);
 

applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/PV3FoamReader_SM.xml

@@ -30,16 +30,16 @@
 
     <!-- Global settings -->
 
-    <!-- Extrapolate Walls check-box -->
+    <!-- Extrapolate Patches check-box -->
     <IntVectorProperty
-      name="ExtrapolateWalls"
+      name="ExtrapolatePatches"
-      command="SetExtrapolateWalls"
+      command="SetExtrapolatePatches"
       number_of_elements="1"
       default_values="0"
       animateable="0">
       <BooleanDomain name="bool"/>
       <Documentation>
-        Extrapolate internalField to wall and empty patches
+        Extrapolate internalField to non-constraint patches
       </Documentation>
     </IntVectorProperty>
 

applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx

@@ -64,7 +64,7 @@ vtkPV3FoamReader::vtkPV3FoamReader()
 
     CacheMesh = 1;
 
-    ExtrapolateWalls = 0;
+    ExtrapolatePatches = 0;
     IncludeSets = 0;
     IncludeZones = 0;
     ShowPatchNames = 0;

applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.h

@@ -65,9 +65,9 @@ public:
     vtkGetMacro(CacheMesh, int);
 
     // Description:
-    // FOAM extrapolate internal values onto the walls
+    // FOAM extrapolate internal values onto the patches
-    vtkSetMacro(ExtrapolateWalls, int);
+    vtkSetMacro(ExtrapolatePatches, int);
-    vtkGetMacro(ExtrapolateWalls, int);
+    vtkGetMacro(ExtrapolatePatches, int);
 
     // FOAM read sets control
     vtkSetMacro(IncludeSets, int);
@@ -183,7 +183,7 @@ private:
     int TimeStepRange[2];
     int CacheMesh;
 
-    int ExtrapolateWalls;
+    int ExtrapolatePatches;
     int IncludeSets;
     int IncludeZones;
     int ShowPatchNames;

applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C

@@ -659,29 +659,55 @@ void Foam::vtkPV3Foam::addPatchNames(vtkRenderer* renderer)
         }
     }
 
+    // Count number of zones we're actually going to display. This is truncated
+    // to a max per patch
+
+    const label MAXPATCHZONES = 20;
+
+    label displayZoneI = 0;
+
+    forAll(pbMesh, patchI)
+    {
+        displayZoneI += min(MAXPATCHZONES, nZones[patchI]);
+    }
+
+
     zoneCentre.shrink();
 
     if (debug)
     {
         Info<< "patch zone centres = " << zoneCentre << nl
+            << "displayed zone centres = " << displayZoneI << nl
             << "zones per patch = " << nZones << endl;
     }
 
     // Set the size of the patch labels to max number of zones
-    patchTextActorsPtrs_.setSize(zoneCentre.size());
+    patchTextActorsPtrs_.setSize(displayZoneI);
 
     if (debug)
     {
         Info<< "constructing patch labels" << endl;
     }
 
+    // Actor index
+    displayZoneI = 0;
+
+    // Index in zone centres
     label globalZoneI = 0;
+
     forAll(pbMesh, patchI)
     {
         const polyPatch& pp = pbMesh[patchI];
 
         // Only selected patches will have a non-zero number of zones
-        for (label i=0; i<nZones[patchI]; i++)
+        label nDisplayZones = min(MAXPATCHZONES, nZones[patchI]);
+        label increment = 1;
+        if (nZones[patchI] >= MAXPATCHZONES)
+        {
+            increment = nZones[patchI]/MAXPATCHZONES;
+        }
+
+        for (label i = 0; i < nDisplayZones; i++)
         {
             if (debug)
             {
@@ -719,14 +745,15 @@ void Foam::vtkPV3Foam::addPatchNames(vtkRenderer* renderer)
 
             // Maintain a list of text labels added so that they can be
             // removed later
-            patchTextActorsPtrs_[globalZoneI] = txt;
+            patchTextActorsPtrs_[displayZoneI] = txt;
 
-            globalZoneI++;
+            globalZoneI += increment;
+            displayZoneI++;
         }
     }
 
     // Resize the patch names list to the actual number of patch names added
-    patchTextActorsPtrs_.setSize(globalZoneI);
+    patchTextActorsPtrs_.setSize(displayZoneI);
 
     if (debug)
     {

applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3FoamVolFields.H

@@ -132,8 +132,8 @@ void Foam::vtkPV3Foam::convertVolFields
                 isType<emptyFvPatchField<Type> >(ptf)
              ||
                 (
-                    typeid(patches[patchId]) == typeid(wallPolyPatch)
+                    reader_->GetExtrapolatePatches()
-                 && reader_->GetExtrapolateWalls()
+                && !polyPatch::constraintType(patches[patchId].type())
                 )
             )
             {

src/OpenFOAM/db/Time/findInstance.C

@@ -58,8 +58,9 @@ Foam::word Foam::Time::findInstance
     {
         if (debug)
         {
-            Info<< "Time::findInstance(const fileName&, const word&) : "
+            Info<< "Time::findInstance"
-                << "found \"" << name
+                "(const fileName&, const word&, const IOobject::readOption)"
+                << " : found \"" << name
                 << "\" in " << timeName()/dir
                 << endl;
         }
@@ -98,8 +99,8 @@ Foam::word Foam::Time::findInstance
             if (debug)
             {
                 Info<< "Time::findInstance"
-                    "(const fileName&,const word&) : "
+                    "(const fileName&, const word&, const IOobject::readOption)"
-                    << "found \"" << name
+                    << " : found \"" << name
                     << "\" in " << ts[instanceI].name()/dir
                     << endl;
             }
@@ -129,8 +130,8 @@ Foam::word Foam::Time::findInstance
         if (debug)
         {
             Info<< "Time::findInstance"
-                   "(const fileName&,const word&) : "
+                "(const fileName&, const word&, const IOobject::readOption)"
-                << "found \"" << name
+                << " : found \"" << name
                 << "\" in " << constant()/dir
                 << endl;
         }
@@ -142,9 +143,9 @@ Foam::word Foam::Time::findInstance
     {
         FatalErrorIn
         (
-                "Time::findInstance(const fileName&,const word&)"
+            "Time::findInstance"
-            )
+            "(const fileName&, const word&, const IOobject::readOption)"
-            << "Cannot find file \"" << name << "\" in directory "
+        )   << "Cannot find file \"" << name << "\" in directory "
             << constant()/dir
             << exit(FatalError);
     }

src/OpenFOAM/meshes/polyMesh/polyMesh.C

@@ -54,40 +54,79 @@ void Foam::polyMesh::calcDirections() const
 {
     for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
     {
-        directions_[cmpt] = 1;
+        solutionD_[cmpt] = 1;
     }
 
-    label nEmptyPatches = 0;
+    // Knock out empty and wedge directions. Note:they will be present on all
+    // domains.
 
-    vector dirVec = vector::zero;
+    label nEmptyPatches = 0;
+    label nWedgePatches = 0;
+
+    vector emptyDirVec = vector::zero;
+    vector wedgeDirVec = vector::zero;
 
     forAll(boundaryMesh(), patchi)
-    {
-        if (isA<emptyPolyPatch>(boundaryMesh()[patchi]))
     {
         if (boundaryMesh()[patchi].size())
+        {
+            if (isA<emptyPolyPatch>(boundaryMesh()[patchi]))
             {
                 nEmptyPatches++;
-                dirVec += sum(cmptMag(boundaryMesh()[patchi].faceAreas()));
+                emptyDirVec += sum(cmptMag(boundaryMesh()[patchi].faceAreas()));
+            }
+            else if (isA<wedgePolyPatch>(boundaryMesh()[patchi]))
+            {
+                const wedgePolyPatch& wpp = refCast<const wedgePolyPatch>
+                (
+                    boundaryMesh()[patchi]
+                );
+
+                nWedgePatches++;
+                wedgeDirVec += cmptMag(wpp.centreNormal());
             }
         }
     }
 
     if (nEmptyPatches)
     {
-        reduce(dirVec, sumOp<vector>());
+        reduce(emptyDirVec, sumOp<vector>());
 
-        dirVec /= mag(dirVec);
+        emptyDirVec /= mag(emptyDirVec);
 
         for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
         {
-            if (dirVec[cmpt] > 1e-6)
+            if (emptyDirVec[cmpt] > 1e-6)
             {
-                directions_[cmpt] = -1;
+                solutionD_[cmpt] = -1;
             }
             else
             {
-                directions_[cmpt] = 1;
+                solutionD_[cmpt] = 1;
+            }
+        }
+    }
+
+
+    // Knock out wedge directions
+
+    geometricD_ = solutionD_;
+
+    if (nWedgePatches)
+    {
+        reduce(wedgeDirVec, sumOp<vector>());
+
+        wedgeDirVec /= mag(wedgeDirVec);
+
+        for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+        {
+            if (wedgeDirVec[cmpt] > 1e-6)
+            {
+                geometricD_[cmpt] = -1;
+            }
+            else
+            {
+                geometricD_[cmpt] = 1;
             }
         }
     }
@@ -163,7 +202,8 @@ Foam::polyMesh::polyMesh(const IOobject& io)
         *this
     ),
     bounds_(points_),
-    directions_(Vector<label>::zero),
+    geometricD_(Vector<label>::zero),
+    solutionD_(Vector<label>::zero),
     pointZones_
     (
         IOobject
@@ -350,7 +390,8 @@ Foam::polyMesh::polyMesh
         0
     ),
     bounds_(points_, syncPar),
-    directions_(Vector<label>::zero),
+    geometricD_(Vector<label>::zero),
+    solutionD_(Vector<label>::zero),
     pointZones_
     (
         IOobject
@@ -505,7 +546,8 @@ Foam::polyMesh::polyMesh
         0
     ),
     bounds_(points_, syncPar),
-    directions_(Vector<label>::zero),
+    geometricD_(Vector<label>::zero),
+    solutionD_(Vector<label>::zero),
     pointZones_
     (
         IOobject
@@ -766,44 +808,37 @@ const Foam::fileName& Foam::polyMesh::facesInstance() const
 }
 
 
-const Foam::Vector<Foam::label>& Foam::polyMesh::directions() const
+const Foam::Vector<Foam::label>& Foam::polyMesh::geometricD() const
 {
-    if (directions_.x() == 0)
+    if (geometricD_.x() == 0)
     {
         calcDirections();
     }
 
-    return directions_;
+    return geometricD_;
 }
 
 
 Foam::label Foam::polyMesh::nGeometricD() const
 {
-    label nWedges = 0;
+    return cmptSum(geometricD() + Vector<label>::one)/2;
-
-    forAll(boundary_, patchi)
-    {
-        if (isA<wedgePolyPatch>(boundary_[patchi]))
-        {
-            nWedges++;
-        }
 }
 
-    if (nWedges != 0 && nWedges != 2 && nWedges != 4)
+
+const Foam::Vector<Foam::label>& Foam::polyMesh::solutionD() const
 {
-        FatalErrorIn("label polyMesh::nGeometricD() const")
+    if (solutionD_.x() == 0)
-            << "Number of wedge patches " << nWedges << " is incorrect, "
+    {
-               "should be 0, 2 or 4"
+        calcDirections();
-            << exit(FatalError);
     }
 
-    return nSolutionD() - nWedges/2;
+    return solutionD_;
 }
 
 
 Foam::label Foam::polyMesh::nSolutionD() const
 {
-    return cmptSum(directions() + Vector<label>::one)/2;
+    return cmptSum(solutionD() + Vector<label>::one)/2;
 }
 
 
@@ -823,6 +858,10 @@ void Foam::polyMesh::addPatches
             << abort(FatalError);
     }
 
+    // Reset valid directions
+    geometricD_ = Vector<label>::zero;
+    solutionD_ = Vector<label>::zero;
+
     boundary_.setSize(p.size());
 
     // Copy the patch pointers
@@ -1037,6 +1076,10 @@ Foam::tmp<Foam::scalarField> Foam::polyMesh::movePoints
     faceZones_.movePoints(points_);
     cellZones_.movePoints(points_);
 
+    // Reset valid directions (could change with rotation)
+    geometricD_ = Vector<label>::zero;
+    solutionD_ = Vector<label>::zero;
+
 
     // Hack until proper callbacks. Below are all the polyMeh MeshObjects with a
     // movePoints function.

src/OpenFOAM/meshes/polyMesh/polyMesh.H

@@ -120,9 +120,13 @@ private:
             //  Created from points on construction, updated when the mesh moves
             boundBox bounds_;
 
+            //- vector of non-constrained directions in mesh
+            //  defined according to the presence of empty and wedge patches
+            mutable Vector<label> geometricD_;
+
             //- vector of valid directions in mesh
             //  defined according to the presence of empty patches
-            mutable Vector<label> directions_;
+            mutable Vector<label> solutionD_;
 
 
         // Zoning information
@@ -309,17 +313,22 @@ public:
                 return bounds_;
             }
 
-            //- Return the vector of valid directions in mesh.
+            //- Return the vector of geometric directions in mesh.
-            //  Defined according to the presence of empty patches.
+            //  Defined according to the presence of empty and wedge patches.
-            //  1 indicates valid direction and -1 an invalid direction.
+            //  1 indicates unconstrained direction and -1 a constrained
-            const Vector<label>& directions() const;
+            //  direction.
+            const Vector<label>& geometricD() const;
 
             //- Return the number of valid geometric dimensions in the mesh
             label nGeometricD() const;
 
-            //- Return the number of valid solution dimensions in the mesh.
+            //- Return the vector of solved-for directions in mesh.
-            //  For wedge cases this includes the circumferential direction
+            //  Differs from geometricD in that it includes for wedge cases
-            //  in case of swirl.
+            //  the circumferential direction in case of swirl.
+            //  1 indicates valid direction and -1 an invalid direction.
+            const Vector<label>& solutionD() const;
+
+            //- Return the number of valid solved-for dimensions in the mesh
             label nSolutionD() const;
 
             //- Return point zone mesh

src/OpenFOAM/meshes/polyMesh/polyMeshClear.C

@@ -68,6 +68,10 @@ void Foam::polyMesh::clearGeom()
         boundary_[patchI].clearGeom();
     }
 
+    // Reset valid directions (could change with rotation)
+    geometricD_ = Vector<label>::zero;
+    solutionD_ = Vector<label>::zero;
+
     pointMesh::Delete(*this);
 }
 
@@ -87,6 +91,10 @@ void Foam::polyMesh::clearAddressing()
     // recalculation
     deleteDemandDrivenData(globalMeshDataPtr_);
 
+    // Reset valid directions
+    geometricD_ = Vector<label>::zero;
+    solutionD_ = Vector<label>::zero;
+
     pointMesh::Delete(*this);
 }
 

src/OpenFOAM/meshes/polyMesh/polyMeshFromShapeMesh.C

@@ -217,7 +217,8 @@ Foam::polyMesh::polyMesh
         boundaryFaces.size() + 1    // add room for a default patch
     ),
     bounds_(points_, syncPar),
-    directions_(Vector<label>::zero),
+    geometricD_(Vector<label>::zero),
+    solutionD_(Vector<label>::zero),
     pointZones_
     (
         IOobject

src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C

@@ -68,6 +68,11 @@ void Foam::polyMesh::updateMesh(const mapPolyMesh& mpm)
         newMotionPoints.map(oldMotionPoints, mpm.pointMap());
     }
 
+    // Reset valid directions (could change by faces put into empty patches)
+    geometricD_ = Vector<label>::zero;
+    solutionD_ = Vector<label>::zero;
+
+
     // Hack until proper callbacks. Below are all the polyMesh-MeshObjects.
 
     // pointMesh

src/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C

@@ -29,6 +29,7 @@ License
 #include "mathematicalConstants.H"
 #include "refinementSurfaces.H"
 #include "searchableSurfaces.H"
+#include "regExp.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -300,10 +301,44 @@ Foam::layerParameters::layerParameters
             //    readScalar(layerDict.lookup("minThickness"));
         }
     }
+
+
+    // Check whether layer specification matches any patches
+    const List<keyType> wildCards = layersDict.keys(true);
+
+    forAll(wildCards, i)
+    {
+        regExp re(wildCards[i]);
+
+        bool hasMatch = false;
+        forAll(boundaryMesh, patchI)
+        {
+            if (re.match(boundaryMesh[patchI].name()))
+            {
+                hasMatch = true;
+                break;
+            }
+        }
+        if (!hasMatch)
+        {
+            IOWarningIn("layerParameters::layerParameters(..)", layersDict)
+                << "Wildcard layer specification for " << wildCards[i]
+                << " does not match any patch." << endl;
+        }
+    }
+
+    const List<keyType> nonWildCards = layersDict.keys(false);
+
+    forAll(nonWildCards, i)
+    {
+        if (boundaryMesh.findPatchID(nonWildCards[i]) == -1)
+        {
+            IOWarningIn("layerParameters::layerParameters(..)", layersDict)
+                << "Layer specification for " << nonWildCards[i]
+                << " does not match any patch." << endl;
+        }
+    }
 }
 
 
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-
 // ************************************************************************* //

src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C

@@ -265,31 +265,47 @@ Foam::refinementSurfaces::refinementSurfaces
     zoneInside_(surfacesDict.size()),
     regionOffset_(surfacesDict.size())
 {
-    labelList globalMinLevel(surfacesDict.size(), 0);
+    // Wilcard specification : loop over all surface, all regions
-    labelList globalMaxLevel(surfacesDict.size(), 0);
+    // and try to find a match.
-    scalarField globalAngle(surfacesDict.size(), -GREAT);
-    List<Map<label> > regionMinLevel(surfacesDict.size());
-    List<Map<label> > regionMaxLevel(surfacesDict.size());
-    List<Map<scalar> > regionAngle(surfacesDict.size());
 
+    // Count number of surfaces.
     label surfI = 0;
-    forAllConstIter(dictionary, surfacesDict, iter)
+    forAll(allGeometry.names(), geomI)
     {
-        names_[surfI] = iter().keyword();
+        const word& geomName = allGeometry_.names()[geomI];
 
-        surfaces_[surfI] = allGeometry_.findSurfaceID(names_[surfI]);
+        if (surfacesDict.found(geomName))
-
-        if (surfaces_[surfI] == -1)
         {
-            FatalErrorIn
+            surfI++;
-            (
-                "refinementSurfaces::refinementSurfaces"
-                "(const searchableSurfaces&, const dictionary>&"
-            )   << "No surface called " << iter().keyword() << endl
-                << "Valid surfaces are " << allGeometry_.names()
-                << exit(FatalError);
         }
-        const dictionary& dict = surfacesDict.subDict(iter().keyword());
+    }
+
+    // Size lists
+    surfaces_.setSize(surfI);
+    names_.setSize(surfI);
+    faceZoneNames_.setSize(surfI);
+    cellZoneNames_.setSize(surfI);
+    zoneInside_.setSize(surfI);
+    regionOffset_.setSize(surfI);
+
+    labelList globalMinLevel(surfI, 0);
+    labelList globalMaxLevel(surfI, 0);
+    scalarField globalAngle(surfI, -GREAT);
+    List<Map<label> > regionMinLevel(surfI);
+    List<Map<label> > regionMaxLevel(surfI);
+    List<Map<scalar> > regionAngle(surfI);
+
+    surfI = 0;
+    forAll(allGeometry.names(), geomI)
+    {
+        const word& geomName = allGeometry_.names()[geomI];
+
+        if (surfacesDict.found(geomName))
+        {
+            const dictionary& dict = surfacesDict.subDict(geomName);
+
+            names_[surfI] = geomName;
+            surfaces_[surfI] = geomI;
 
             const labelPair refLevel(dict.lookup("level"));
             globalMinLevel[surfI] = refLevel[0];
@@ -306,7 +322,10 @@ Foam::refinementSurfaces::refinementSurfaces
             // Global perpendicular angle
             if (dict.found("perpendicularAngle"))
             {
-            globalAngle[surfI] = readScalar(dict.lookup("perpendicularAngle"));
+                globalAngle[surfI] = readScalar
+                (
+                    dict.lookup("perpendicularAngle")
+                );
             }
 
             if (dict.found("regions"))
@@ -315,27 +334,15 @@ Foam::refinementSurfaces::refinementSurfaces
                 const wordList& regionNames =
                     allGeometry_[surfaces_[surfI]].regions();
 
-            forAllConstIter(dictionary, regionsDict, iter)
+                forAll(regionNames, regionI)
                 {
-                const word& key = iter().keyword();
+                    if (regionsDict.found(regionNames[regionI]))
-
-                if (regionsDict.isDict(key))
                     {
-                    // Get the dictionary for region iter.keyword()
+                        // Get the dictionary for region 
-                    const dictionary& regionDict = regionsDict.subDict(key);
+                        const dictionary& regionDict = regionsDict.subDict
-
-                    label regionI = findIndex(regionNames, key);
-                    if (regionI == -1)
-                    {
-                        FatalErrorIn
                         (
-                            "refinementSurfaces::refinementSurfaces"
+                            regionNames[regionI]
-                            "(const searchableSurfaces&, const dictionary>&"
+                        );
-                        )   << "No region called " << key << " on surface "
-                            << allGeometry_[surfaces_[surfI]].name() << endl
-                            << "Valid regions are " << regionNames
-                            << exit(FatalError);
-                    }
 
                         const labelPair refLevel(regionDict.lookup("level"));
 
@@ -347,7 +354,10 @@ Foam::refinementSurfaces::refinementSurfaces
                             regionAngle[surfI].insert
                             (
                                 regionI,
-                            readScalar(regionDict.lookup("perpendicularAngle"))
+                                readScalar
+                                (
+                                    regionDict.lookup("perpendicularAngle")
+                                )
                             );
                         }
                     }
@@ -355,11 +365,12 @@ Foam::refinementSurfaces::refinementSurfaces
             }
             surfI++;
         }
+    }
 
     // Calculate local to global region offset
     label nRegions = 0;
 
-    forAll(surfacesDict, surfI)
+    forAll(surfaces_, surfI)
     {
         regionOffset_[surfI] = nRegions;
         nRegions += allGeometry_[surfaces_[surfI]].regions().size();

src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGradData.C

@@ -140,11 +140,11 @@ void Foam::quadraticFitSnGradData::findFaceDirs
     #ifndef SPHERICAL_GEOMETRY
         if (mesh.nGeometricD() <= 2) // find the normal direcion
         {
-            if (mesh.directions()[0] == -1)
+            if (mesh.geometricD()[0] == -1)
             {
                 kdir = vector(1, 0, 0);
             }
-            else if (mesh.directions()[1] == -1)
+            else if (mesh.geometricD()[1] == -1)
             {
                 kdir = vector(0, 1, 0);
             }
@@ -153,7 +153,7 @@ void Foam::quadraticFitSnGradData::findFaceDirs
                 kdir = vector(0, 0, 1);
             }
         }
-        else // 3D so find a direction in the place of the face
+        else // 3D so find a direction in the plane of the face
         {
             const face& f = mesh.faces()[faci];
             kdir = mesh.points()[f[0]] - mesh.points()[f[1]];

src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C

@@ -713,7 +713,7 @@ Foam::fvMatrix<Type>::H() const
     (
         pow
         (
-            psi_.mesh().directions(),
+            psi_.mesh().solutionD(),
             pTraits<typename powProduct<Vector<label>, Type::rank>::type>::zero
         )
     );

src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C

@@ -82,7 +82,7 @@ Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::solve
     (
         pow
         (
-            psi_.mesh().directions(),
+            psi_.mesh().solutionD(),
             pTraits<typename powProduct<Vector<label>, Type::rank>::type>::zero
         )
     );

src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C

@@ -83,11 +83,11 @@ void Foam::FitData<FitDataType, ExtendedStencil, Polynomial>::findFaceDirs
 #   ifndef SPHERICAL_GEOMETRY
     if (mesh.nGeometricD() <= 2) // find the normal direction
     {
-        if (mesh.directions()[0] == -1)
+        if (mesh.geometricD()[0] == -1)
         {
             kdir = vector(1, 0, 0);
         }
-        else if (mesh.directions()[1] == -1)
+        else if (mesh.geometricD()[1] == -1)
         {
             kdir = vector(0, 1, 0);
         }
@@ -115,7 +115,7 @@ void Foam::FitData<FitDataType, ExtendedStencil, Polynomial>::findFaceDirs
 
         if (magk < SMALL)
         {
-            FatalErrorIn("findFaceDirs") << " calculated kdir = zero"
+            FatalErrorIn("findFaceDirs(..)") << " calculated kdir = zero"
                 << exit(FatalError);
         }
         else

src/meshTools/Make/files

@@ -60,6 +60,7 @@ $(searchableSurface)/searchablePlane.C
 $(searchableSurface)/searchablePlate.C
 $(searchableSurface)/searchableSphere.C
 $(searchableSurface)/searchableSurface.C
+$(searchableSurface)/searchableSurfaceCollection.C
 $(searchableSurface)/searchableSurfaces.C
 $(searchableSurface)/searchableSurfacesQueries.C
 $(searchableSurface)/searchableSurfaceWithGaps.C

src/meshTools/meshTools/meshTools.C

@@ -25,7 +25,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "meshTools.H"
-#include "primitiveMesh.H"
+#include "polyMesh.H"
 #include "hexMatcher.H"
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
@@ -635,6 +635,103 @@ Foam::label Foam::meshTools::walkFace
 }
 
 
+void Foam::meshTools::constrainToMeshCentre(const polyMesh& mesh, point& pt)
+{
+    const Vector<label>& dirs = mesh.geometricD();
+    const point& min = mesh.bounds().min();
+    const point& max = mesh.bounds().max();
+
+    for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+    {
+        if (dirs[cmpt] == -1)
+        {
+            pt[cmpt] = 0.5*(min[cmpt]+max[cmpt]);
+        }
+    }
+}
+
+
+void Foam::meshTools::constrainToMeshCentre
+(
+    const polyMesh& mesh,
+    pointField& pts
+)
+{
+    const Vector<label>& dirs = mesh.geometricD();
+    const point& min = mesh.bounds().min();
+    const point& max = mesh.bounds().max();
+
+    bool isConstrained = false;
+    for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+    {
+        if (dirs[cmpt] == -1)
+        {
+            isConstrained = true;
+            break;
+        }
+    }
+
+    if (isConstrained)
+    {
+        forAll(pts, i)
+        {
+            for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+            {
+                if (dirs[cmpt] == -1)
+                {
+                    pts[i][cmpt] = 0.5*(min[cmpt]+max[cmpt]);
+                }
+            }
+        }
+    }
+}
+
+
+//- Set the constrained components of directions/velocity to zero
+void Foam::meshTools::constrainDirection(const polyMesh& mesh, vector& d)
+{
+    const Vector<label>& dirs = mesh.geometricD();
+
+    for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+    {
+        if (dirs[cmpt] == -1)
+        {
+            d[cmpt] = 0.0;
+        }
+    }
+}
+
+
+void Foam::meshTools::constrainDirection(const polyMesh& mesh, vectorField& d)
+{
+    const Vector<label>& dirs = mesh.geometricD();
+
+    bool isConstrained = false;
+    for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+    {
+        if (dirs[cmpt] == -1)
+        {
+            isConstrained = true;
+            break;
+        }
+    }
+
+    if (isConstrained)
+    {
+        forAll(d, i)
+        {
+            for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+            {
+                if (dirs[cmpt] == -1)
+                {
+                    d[i][cmpt] = 0.0;
+                }
+            }
+        }
+    }
+}
+
+
 void Foam::meshTools::getParallelEdges
 (
     const primitiveMesh& mesh,

src/meshTools/meshTools/meshTools.H

@@ -50,7 +50,7 @@ namespace Foam
 {
 
 class primitiveMesh;
-
+class polyMesh;
 
 /*---------------------------------------------------------------------------*\
                         Namespace meshTools Declaration
@@ -81,6 +81,8 @@ namespace meshTools
     static const label pXpYpZMask = 1 << pXpYpZ;
 
 
+    // Normal handling
+
         //- Check if n is in same direction as normals of all faceLabels
         bool visNormal
         (
@@ -96,6 +98,9 @@ namespace meshTools
         //- Normalized edge vector
         vector normEdgeVec(const primitiveMesh&, const label edgeI);
 
+
+    // OBJ writing
+
         //- Write obj representation of point
         void writeOBJ
         (
@@ -103,7 +108,6 @@ namespace meshTools
             const point& pt
         );
 
-
         //- Write obj representation of faces subset
         void writeOBJ
         (
@@ -131,6 +135,9 @@ namespace meshTools
             const labelList& cellLabels
         );
 
+
+    // Cell/face/edge walking
+
         //- Is edge used by cell
         bool edgeOnCell
         (
@@ -239,9 +246,18 @@ namespace meshTools
         );
 
 
-    //
+    // Constraints on position
+
+        //- Set the constrained components of position to mesh centre
+        void constrainToMeshCentre(const polyMesh&, point&);
+        void constrainToMeshCentre(const polyMesh&, pointField&);
+
+        //- Set the constrained components of directions/velocity to zero
+        void constrainDirection(const polyMesh&, vector&);
+        void constrainDirection(const polyMesh&, vectorField&);
+
+    
     // Hex only functionality.
-    //
 
         //- Given edge on hex find other 'parallel', non-connected edges.
         void getParallelEdges

src/meshTools/searchableSurface/distributedTriSurfaceMesh.C

@@ -1347,19 +1347,6 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh
 Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh(const IOobject& io)
 :
     triSurfaceMesh(io),
-//    triSurfaceMesh
-//    (
-//        IOobject
-//        (
-//            io.name(),
-//            io.db().time().findInstanceDir(io.local()),
-//            io.local(),
-//            io.db(),
-//            io.readOpt(),
-//            io.writeOpt(),
-//            io.registerObject()
-//        )
-//    ),
     dict_
     (
         IOobject
@@ -1385,20 +1372,6 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh
 )
 :
     triSurfaceMesh(io, dict),
-//    triSurfaceMesh
-//    (
-//        IOobject
-//        (
-//            io.name(),
-//            io.db().time().findInstanceDir(io.local()),
-//            io.local(),
-//            io.db(),
-//            io.readOpt(),
-//            io.writeOpt(),
-//            io.registerObject()
-//        ),
-//        dict
-//    ),
     dict_
     (
         IOobject

src/meshTools/searchableSurface/searchableSurfaceCollection.C

@@ -0,0 +1,524 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+     \\/     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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "searchableSurfaceCollection.H"
+#include "addToRunTimeSelectionTable.H"
+#include "SortableList.H"
+#include "Time.H"
+#include "ListOps.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+defineTypeNameAndDebug(searchableSurfaceCollection, 0);
+addToRunTimeSelectionTable(searchableSurface, searchableSurfaceCollection, dict);
+
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+void Foam::searchableSurfaceCollection::findNearest
+(
+    const pointField& samples,
+    scalarField& minDistSqr,
+    List<pointIndexHit>& nearestInfo,
+    labelList& nearestSurf
+) const
+{
+    // Initialise
+    nearestInfo.setSize(samples.size());
+    nearestInfo = pointIndexHit();
+    nearestSurf.setSize(samples.size());
+    nearestSurf = -1;
+
+    List<pointIndexHit> hitInfo(samples.size());
+
+    forAll(subGeom_, surfI)
+    {
+        // Transform then divide
+        tmp<pointField> localSamples = cmptDivide
+        (
+            transform_[surfI].localPosition
+            (
+                samples
+            ),
+            scale_[surfI]
+        );
+        
+        subGeom_[surfI].findNearest(localSamples, minDistSqr, hitInfo);
+
+        forAll(hitInfo, pointI)
+        {
+            if (hitInfo[pointI].hit())
+            {
+                minDistSqr[pointI] = magSqr
+                (   
+                    hitInfo[pointI].hitPoint()
+                  - localSamples()[pointI]
+                );
+
+                // Rework back into global coordinate sys. Multiply then
+                // transform
+                nearestInfo[pointI] = hitInfo[pointI];
+                nearestInfo[pointI].rawPoint() = 
+                    transform_[surfI].globalPosition
+                    (
+                        cmptMultiply
+                        (
+                            nearestInfo[pointI].rawPoint(),
+                            scale_[surfI]
+                        )
+                    );
+                nearestSurf[pointI] = surfI;
+            }
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::searchableSurfaceCollection::searchableSurfaceCollection
+(
+    const IOobject& io,
+    const dictionary& dict
+)
+:
+    searchableSurface(io),
+    instance_(dict.size()),
+    scale_(dict.size()),
+    transform_(dict.size()),
+    subGeom_(dict.size())
+{
+    Info<< "SearchableCollection : " << name() << endl;
+
+    label surfI = 0;
+    forAllConstIter(dictionary, dict, iter)
+    {
+        if (dict.isDict(iter().keyword()))
+        {
+            instance_[surfI] = iter().keyword();
+
+            const dictionary& subDict = dict.subDict(instance_[surfI]);
+
+            scale_[surfI] = subDict.lookup("scale");
+            transform_.set
+            (
+                surfI,
+                coordinateSystem::New
+                (
+                    "",
+                    subDict.subDict("transform")
+                )
+            );
+
+            const word subGeomName(subDict.lookup("surface"));
+            //Pout<< "Trying to find " << subGeomName << endl;
+
+            const searchableSurface& s =
+                io.db().lookupObject<searchableSurface>(subGeomName);
+
+            subGeom_.set(surfI, &const_cast<searchableSurface&>(s));
+
+            Info<< "    instance : " << instance_[surfI] << endl;
+            Info<< "    surface  : " << s.name() << endl;
+            Info<< "    scale    : " << scale_[surfI] << endl;
+            Info<< "    coordsys : " << transform_[surfI] << endl;
+
+            surfI++;
+        }
+    }
+    instance_.setSize(surfI);
+    scale_.setSize(surfI);
+    transform_.setSize(surfI);
+    subGeom_.setSize(surfI);
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::searchableSurfaceCollection::~searchableSurfaceCollection()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+const Foam::wordList& Foam::searchableSurfaceCollection::regions() const
+{
+    if (regions_.size() == 0)
+    {
+        regionOffset_.setSize(subGeom_.size());
+
+        DynamicList<word> allRegions;
+        forAll(subGeom_, surfI)
+        {
+            regionOffset_[surfI] = allRegions.size();
+
+            const wordList& subRegions = subGeom_[surfI].regions();
+
+            forAll(subRegions, i)
+            {
+                //allRegions.append(subRegions[i] + "_" + Foam::name(surfI));
+                allRegions.append(instance_[surfI] + "_" + subRegions[i]);
+            }
+        }
+        regions_.transfer(allRegions.shrink());
+    }
+    return regions_;
+}
+
+
+Foam::label Foam::searchableSurfaceCollection::size() const
+{
+    label n = 0;
+    forAll(subGeom_, surfI)
+    {
+        n += subGeom_[surfI].size();
+    }
+    return n;
+}
+
+
+void Foam::searchableSurfaceCollection::findNearest
+(
+    const pointField& samples,
+    const scalarField& nearestDistSqr,
+    List<pointIndexHit>& nearestInfo
+) const
+{
+    // How to scale distance?
+    scalarField minDistSqr(nearestDistSqr);
+
+    labelList nearestSurf;
+    findNearest
+    (
+        samples,
+        minDistSqr,
+        nearestInfo,
+        nearestSurf
+    );
+}
+
+
+void Foam::searchableSurfaceCollection::findLine
+(
+    const pointField& start,
+    const pointField& end,
+    List<pointIndexHit>& info
+) const
+{
+    info.setSize(start.size());
+    info = pointIndexHit();
+
+    // Current nearest (to start) intersection
+    pointField nearest(end);
+
+    List<pointIndexHit> hitInfo(start.size());
+
+    forAll(subGeom_, surfI)
+    {
+        // Starting point
+        tmp<pointField> e0 = cmptDivide
+        (
+            transform_[surfI].localPosition
+            (
+                start
+            ),
+            scale_[surfI]
+        );
+
+        // Current best end point
+        tmp<pointField> e1 = cmptDivide
+        (
+            transform_[surfI].localPosition
+            (
+                nearest
+            ),
+            scale_[surfI]
+        );
+
+        subGeom_[surfI].findLine(e0, e1, hitInfo);
+
+        forAll(hitInfo, pointI)
+        {
+            if (hitInfo[pointI].hit())
+            {
+                // Transform back to global coordinate sys.
+                nearest[pointI] = transform_[surfI].globalPosition
+                (
+                    cmptMultiply
+                    (
+                        hitInfo[pointI].rawPoint(),
+                        scale_[surfI]
+                    )
+                );
+                info[pointI] = hitInfo[pointI];
+                info[pointI].rawPoint() = nearest[pointI];
+            }
+        }
+    }
+
+
+    // Debug check
+    if (false)
+    {
+        forAll(info, pointI)
+        {
+            if (info[pointI].hit())
+            {
+                vector n(end[pointI] - start[pointI]);
+                scalar magN = mag(n);
+
+                if (magN > SMALL)
+                {
+                    n /= mag(n);
+
+                    scalar s = ((info[pointI].rawPoint()-start[pointI])&n);
+
+                    if (s < 0 || s > 1)
+                    {
+                        FatalErrorIn
+                        (
+                            "searchableSurfaceCollection::findLine(..)"
+                        )   << "point:" << info[pointI]
+                            << " s:" << s
+                            << " outside vector "
+                            << " start:" << start[pointI]
+                            << " end:" << end[pointI]
+                            << abort(FatalError);
+                    }
+                }
+            }
+        }
+    }
+}
+
+
+void Foam::searchableSurfaceCollection::findLineAny
+(
+    const pointField& start,
+    const pointField& end,
+    List<pointIndexHit>& info
+) const
+{
+    // To be done ...
+    findLine(start, end, info);
+}
+
+
+void Foam::searchableSurfaceCollection::findLineAll
+(
+    const pointField& start,
+    const pointField& end,
+    List<List<pointIndexHit> >& info
+) const
+{
+    // To be done. Assume for now only one intersection.
+    List<pointIndexHit> nearestInfo;
+    findLine(start, end, nearestInfo);
+
+    info.setSize(start.size());
+    forAll(info, pointI)
+    {
+        if (nearestInfo[pointI].hit())
+        {
+            info[pointI].setSize(1);
+            info[pointI][0] = nearestInfo[pointI];
+        }
+        else
+        {
+            info[pointI].clear();
+        }
+    }
+}
+
+
+void Foam::searchableSurfaceCollection::getRegion
+(
+    const List<pointIndexHit>& info,
+    labelList& region
+) const
+{
+    if (subGeom_.size() == 0)
+    {}
+    else if (subGeom_.size() == 1)
+    {
+        subGeom_[0].getRegion(info, region);
+    }
+    else
+    {
+        region.setSize(info.size());
+        region = -1;
+
+        // Which region did point come from. Retest for now to see which
+        // surface it originates from - crap solution! Should use global indices
+        // in index inside pointIndexHit to do this better.
+
+        pointField samples(info.size());
+        forAll(info, pointI)
+        {
+            if (info[pointI].hit())
+            {
+                samples[pointI] = info[pointI].hitPoint();
+            }
+            else
+            {
+                samples[pointI] = vector::zero;
+            }
+        }
+        //scalarField minDistSqr(info.size(), SMALL);
+        scalarField minDistSqr(info.size(), GREAT);
+
+        labelList nearestSurf;
+        List<pointIndexHit> nearestInfo;
+        findNearest
+        (
+            samples,
+            minDistSqr,
+            nearestInfo,
+            nearestSurf
+        );
+
+        // Check
+        {
+            forAll(info, pointI)
+            {
+                if (info[pointI].hit() && nearestSurf[pointI] == -1)
+                {
+                    FatalErrorIn
+                    (
+                        "searchableSurfaceCollection::getRegion(..)"
+                    )   << "pointI:" << pointI
+                        << " sample:" << samples[pointI]
+                        << " nearest:" << nearestInfo[pointI]
+                        << " nearestsurf:" << nearestSurf[pointI]
+                        << abort(FatalError);
+                }
+            }
+        }
+
+        forAll(subGeom_, surfI)
+        {
+            // Collect points from my surface
+            labelList indices(findIndices(nearestSurf, surfI));
+
+            labelList surfRegion;
+            subGeom_[surfI].getRegion
+            (
+                IndirectList<pointIndexHit>(info, indices),
+                surfRegion
+            );
+            forAll(indices, i)
+            {
+                region[indices[i]] = regionOffset_[surfI] + surfRegion[i];
+            }
+        }
+    }
+}
+
+
+void Foam::searchableSurfaceCollection::getNormal
+(
+    const List<pointIndexHit>& info,
+    vectorField& normal
+) const
+{
+    if (subGeom_.size() == 0)
+    {}
+    else if (subGeom_.size() == 1)
+    {
+        subGeom_[0].getNormal(info, normal);
+    }
+    else
+    {
+        normal.setSize(info.size());
+
+        // See above - crap retest to find surface point originates from.
+        pointField samples(info.size());
+        forAll(info, pointI)
+        {
+            if (info[pointI].hit())
+            {
+                samples[pointI] = info[pointI].hitPoint();
+            }
+            else
+            {
+                samples[pointI] = vector::zero;
+            }
+        }
+        //scalarField minDistSqr(info.size(), SMALL);
+        scalarField minDistSqr(info.size(), GREAT);
+
+        labelList nearestSurf;
+        List<pointIndexHit> nearestInfo;
+        findNearest
+        (
+            samples,
+            minDistSqr,
+            nearestInfo,
+            nearestSurf
+        );
+
+
+        forAll(subGeom_, surfI)
+        {
+            // Collect points from my surface
+            labelList indices(findIndices(nearestSurf, surfI));
+
+            vectorField surfNormal;
+            subGeom_[surfI].getNormal
+            (
+                IndirectList<pointIndexHit>(info, indices),
+                surfNormal
+            );
+            forAll(indices, i)
+            {
+                normal[indices[i]] = surfNormal[i];
+            }
+        }
+    }
+}
+
+
+void Foam::searchableSurfaceCollection::getVolumeType
+(
+    const pointField& points,
+    List<volumeType>& volType
+) const
+{
+    FatalErrorIn
+    (
+        "searchableSurfaceCollection::getVolumeType(const pointField&"
+        ", List<volumeType>&) const"
+    )   << "Volume type not supported for collection."
+        << exit(FatalError);
+}
+
+
+// ************************************************************************* //

src/meshTools/searchableSurface/searchableSurfaceCollection.H

@@ -0,0 +1,212 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+     \\/     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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    Foam::searchableSurfaceCollection
+
+Description
+    Union of transformed searchableSurfaces
+
+SourceFiles
+    searchableSurfaceCollection.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef searchableSurfaceCollection_H
+#define searchableSurfaceCollection_H
+
+#include "searchableSurface.H"
+#include "treeBoundBox.H"
+#include "coordinateSystem.H"
+#include "UPtrList.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward declaration of classes
+
+/*---------------------------------------------------------------------------*\
+                           Class searchableSurfaceCollection Declaration
+\*---------------------------------------------------------------------------*/
+
+class searchableSurfaceCollection
+:
+    public searchableSurface
+{
+private:
+
+    // Private Member Data
+
+        // Per instance data
+
+            //- instance name
+            wordList instance_;
+
+            //- scaling vector
+            vectorField scale_;
+
+            //- transformation
+            PtrList<coordinateSystem> transform_;
+
+            UPtrList<searchableSurface> subGeom_;
+
+        //- Region names
+        mutable wordList regions_;
+        //- From individual regions to collection regions
+        mutable labelList regionOffset_;
+
+
+    // Private Member Functions
+
+        //- Find point nearest to sample. Updates minDistSqr. Sets nearestInfo
+        //  and surface index
+        void findNearest
+        (
+            const pointField& samples,
+            scalarField& minDistSqr,
+            List<pointIndexHit>& nearestInfo,
+            labelList& nearestSurf
+        ) const;
+
+
+        //- Disallow default bitwise copy construct
+        searchableSurfaceCollection(const searchableSurfaceCollection&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const searchableSurfaceCollection&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("searchableSurfaceCollection");
+
+
+    // Constructors
+
+        //- Construct from dictionary (used by searchableSurface)
+        searchableSurfaceCollection
+        (
+            const IOobject& io,
+            const dictionary& dict
+        );
+
+    // Destructor
+
+        virtual ~searchableSurfaceCollection();
+
+
+    // Member Functions
+
+        virtual const wordList& regions() const;
+
+        //- Whether supports volume type below
+        virtual bool hasVolumeType() const
+        {
+            return false;
+        }
+
+        //- Range of local indices that can be returned.
+        virtual label size() const;
+
+
+        // Multiple point queries.
+
+            virtual void findNearest
+            (
+                const pointField& sample,
+                const scalarField& nearestDistSqr,
+                List<pointIndexHit>&
+            ) const;
+
+            virtual void findLine
+            (
+                const pointField& start,
+                const pointField& end,
+                List<pointIndexHit>&
+            ) const;
+
+            virtual void findLineAny
+            (
+                const pointField& start,
+                const pointField& end,
+                List<pointIndexHit>&
+            ) const;
+
+            //- Get all intersections in order from start to end.
+            virtual void findLineAll
+            (
+                const pointField& start,
+                const pointField& end,
+                List<List<pointIndexHit> >&
+            ) const;
+
+            //- From a set of points and indices get the region
+            virtual void getRegion
+            (
+                const List<pointIndexHit>&,
+                labelList& region
+            ) const;
+
+            //- From a set of points and indices get the normal
+            virtual void getNormal
+            (
+                const List<pointIndexHit>&,
+                vectorField& normal
+            ) const;
+
+            //- Determine type (inside/outside/mixed) for point. unknown if
+            //  cannot be determined (e.g. non-manifold surface)
+            virtual void getVolumeType
+            (
+                const pointField&,
+                List<volumeType>&
+            ) const;
+
+
+        // regIOobject implementation
+
+            bool writeData(Ostream&) const
+            {
+                notImplemented
+                (
+                    "searchableSurfaceCollection::writeData(Ostream&) const"
+                );
+                return false;
+            }
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/searchableSurface/searchableSurfaces.C

@@ -39,8 +39,6 @@ defineTypeNameAndDebug(searchableSurfaces, 0);
 }
 
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Construct with length.
@@ -187,14 +185,6 @@ Foam::searchableSurfaces::searchableSurfaces
         // their object name. Maybe have stlTriSurfaceMesh which appends .stl
         // when reading/writing?
         namedIO().rename(key);  // names_[surfI]
-        if (namedIO().local() != word::null)
-        {
-            namedIO().instance() = namedIO().time().findInstance
-            (
-                namedIO().local(),
-                namedIO().name()
-            );
-        }
 
         // Create and hook surface
         set

src/meshTools/searchableSurface/searchableSurfaces.H

@@ -69,6 +69,8 @@ class searchableSurfaces
         labelList allSurfaces_;
 
 
+    // Private Member Functions
+
         //- Disallow default bitwise copy construct
         searchableSurfaces(const searchableSurfaces&);
 

src/meshTools/searchableSurface/triSurfaceMesh.C

@@ -43,6 +43,64 @@ addToRunTimeSelectionTable(searchableSurface, triSurfaceMesh, dict);
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+//// Special version of Time::findInstance that does not check headerOk
+//// to search for instances of raw files
+//Foam::word Foam::triSurfaceMesh::findRawInstance
+//(
+//    const Time& runTime,
+//    const fileName& dir,
+//    const word& name
+//)
+//{
+//    // Check current time first
+//    if (isFile(runTime.path()/runTime.timeName()/dir/name))
+//    {
+//        return runTime.timeName();
+//    }
+//    instantList ts = runTime.times();
+//    label instanceI;
+//
+//    for (instanceI = ts.size()-1; instanceI >= 0; --instanceI)
+//    {
+//        if (ts[instanceI].value() <= runTime.timeOutputValue())
+//        {
+//            break;
+//        }
+//    }
+//
+//    // continue searching from here
+//    for (; instanceI >= 0; --instanceI)
+//    {
+//        if (isFile(runTime.path()/ts[instanceI].name()/dir/name))
+//        {
+//            return ts[instanceI].name();
+//        }
+//    }
+//
+//
+//    // not in any of the time directories, try constant
+//
+//    // Note. This needs to be a hard-coded constant, rather than the
+//    // constant function of the time, because the latter points to
+//    // the case constant directory in parallel cases
+//
+//    if (isFile(runTime.path()/runTime.constant()/dir/name))
+//    {
+//        return runTime.constant();
+//    }
+//
+//    FatalErrorIn
+//    (
+//        "searchableSurfaces::findRawInstance"
+//        "(const Time&, const fileName&, const word&)"
+//    )   << "Cannot find file \"" << name << "\" in directory "
+//        << runTime.constant()/dir
+//        << exit(FatalError);
+//
+//    return runTime.constant();
+//}
+
+
 //- Check file existence
 const Foam::fileName& Foam::triSurfaceMesh::checkFile
 (
@@ -149,7 +207,19 @@ bool Foam::triSurfaceMesh::isSurfaceClosed() const
 Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
 :
     searchableSurface(io),
-    objectRegistry(io),
+    objectRegistry
+    (
+        IOobject
+        (
+            io.name(),
+            io.instance(),
+            io.local(),
+            io.db(),
+            io.readOpt(),
+            io.writeOpt(),
+            false       // searchableSurface already registered under name
+        )
+    ),
     triSurface(s),
     surfaceClosed_(-1)
 {}
@@ -157,8 +227,34 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
 
 Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io)
 :
-    searchableSurface(io),
+    // Find instance for triSurfaceMesh
-    objectRegistry(io),
+    searchableSurface
+    (
+        IOobject
+        (
+            io.name(),
+            io.time().findInstance(io.local(), word::null),
+            io.local(),
+            io.db(),
+            io.readOpt(),
+            io.writeOpt(),
+            io.registerObject()
+        )
+    ),
+    // Reused found instance in objectRegistry
+    objectRegistry
+    (
+        IOobject
+        (
+            io.name(),
+            static_cast<const searchableSurface&>(*this).instance(),
+            io.local(),
+            io.db(),
+            io.readOpt(),
+            io.writeOpt(),
+            false       // searchableSurface already registered under name
+        )
+    ),
     triSurface
     (
         checkFile
@@ -177,8 +273,33 @@ Foam::triSurfaceMesh::triSurfaceMesh
     const dictionary& dict
 )
 :
-    searchableSurface(io),
+    searchableSurface
-    objectRegistry(io),
+    (
+        IOobject
+        (
+            io.name(),
+            io.time().findInstance(io.local(), word::null),
+            io.local(),
+            io.db(),
+            io.readOpt(),
+            io.writeOpt(),
+            io.registerObject()
+        )
+    ),
+    // Reused found instance in objectRegistry
+    objectRegistry
+    (
+        IOobject
+        (
+            io.name(),
+            static_cast<const searchableSurface&>(*this).instance(),
+            io.local(),
+            io.db(),
+            io.readOpt(),
+            io.writeOpt(),
+            false       // searchableSurface already registered under name
+        )
+    ),
     triSurface
     (
         checkFile

src/meshTools/searchableSurface/triSurfaceMesh.H

@@ -76,6 +76,14 @@ private:
 
     // Private Member Functions
 
+        ////- Helper: find instance of files without header
+        //static word findRawInstance
+        //(
+        //    const Time&,
+        //    const fileName&,
+        //    const word&
+        //);
+
         //- Check file existence
         static const fileName& checkFile
         (
@@ -105,10 +113,11 @@ public:
         //- Construct from triSurface
         triSurfaceMesh(const IOobject&, const triSurface&);
 
-        //- Construct read
+        //- Construct read. Does timeInstance search.
         triSurfaceMesh(const IOobject& io);
 
-        //- Construct from dictionary (used by searchableSurface)
+        //- Construct from IO and dictionary (used by searchableSurface).
+        //  Does timeInstance search.
         //  Dictionary may contain a 'scale' entry (eg, 0.001: mm -> m)
         triSurfaceMesh
         (

src/sampling/sampledSurface/isoSurface/isoSurface.C

@@ -29,6 +29,8 @@ License
 #include "mergePoints.H"
 #include "syncTools.H"
 #include "addToRunTimeSelectionTable.H"
+#include "slicedVolFields.H"
+#include "surfaceFields.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -83,9 +85,74 @@ bool Foam::isoSurface::isEdgeOfFaceCut
 }
 
 
+// Get neighbour value and position.
+void Foam::isoSurface::getNeighbour
+(
+    const labelList& boundaryRegion,
+    const volScalarField& cVals,
+    const label cellI,
+    const label faceI,
+    scalar& nbrValue,
+    point& nbrPoint
+) const
+{
+    const labelList& own = mesh_.faceOwner();
+    const labelList& nei = mesh_.faceNeighbour();
+    const surfaceScalarField& weights = mesh_.weights();
+
+    if (mesh_.isInternalFace(faceI))
+    {
+        label nbr = (own[faceI] == cellI ? nei[faceI] : own[faceI]);
+        nbrValue = cVals[nbr];
+        nbrPoint = mesh_.cellCentres()[nbr];
+    }
+    else
+    {
+        label bFaceI = faceI-mesh_.nInternalFaces();
+        label patchI = boundaryRegion[bFaceI];
+        const polyPatch& pp = mesh_.boundaryMesh()[patchI];
+        label patchFaceI = faceI-pp.start();
+
+        if (isA<emptyPolyPatch>(pp))
+        {
+            // Assume zero gradient
+            nbrValue = cVals[own[faceI]];
+            nbrPoint = mesh_.faceCentres()[faceI];
+        }
+        else if (pp.coupled())
+        {
+            if (!refCast<const coupledPolyPatch>(pp).separated())
+            {
+                // Behave as internal face:
+                // other side value
+                nbrValue = cVals.boundaryField()[patchI][patchFaceI];
+                // other side cell centre
+                nbrPoint = mesh_.C().boundaryField()[patchI][patchFaceI];
+            }
+            else
+            {
+                // Do some interpolation for now
+                const scalarField& w = weights.boundaryField()[patchI];
+
+                nbrPoint = mesh_.faceCentres()[faceI];
+                nbrValue =
+                    (1-w[patchFaceI])*cVals[own[faceI]]
+                  + w[patchFaceI]*cVals.boundaryField()[patchI][patchFaceI];
+            }
+        }
+        else
+        {
+            nbrValue = cVals.boundaryField()[patchI][patchFaceI];
+            nbrPoint = mesh_.faceCentres()[faceI];
+        }
+    }
+}
+
+
 // Determine for every face/cell whether it (possibly) generates triangles.
 void Foam::isoSurface::calcCutTypes
 (
+    const labelList& boundaryRegion,
     const volScalarField& cVals,
     const scalarField& pVals
 )
@@ -101,7 +168,20 @@ void Foam::isoSurface::calcCutTypes
     {
         // CC edge.
         bool ownLower = (cVals[own[faceI]] < iso_);
-        bool neiLower = (cVals[nei[faceI]] < iso_);
+
+        scalar nbrValue;
+        point nbrPoint;
+        getNeighbour
+        (
+            boundaryRegion,
+            cVals,
+            own[faceI],
+            faceI,
+            nbrValue,
+            nbrPoint
+        );
+
+        bool neiLower = (nbrValue < iso_);
 
         if (ownLower != neiLower)
         {
@@ -127,15 +207,29 @@ void Foam::isoSurface::calcCutTypes
 
         if (isA<emptyPolyPatch>(pp))
         {
-            // Assume zero gradient so owner and neighbour/boundary value equal
+            // Still needs special treatment?
 
             forAll(pp, i)
             {
                 bool ownLower = (cVals[own[faceI]] < iso_);
 
+                scalar nbrValue;
+                point nbrPoint;
+                getNeighbour
+                (
+                    boundaryRegion,
+                    cVals,
+                    own[faceI],
+                    faceI,
+                    nbrValue,
+                    nbrPoint
+                );
+
+                bool neiLower = (nbrValue < iso_);
+
                 const face f = mesh_.faces()[faceI];
 
-                if (isEdgeOfFaceCut(pVals, f, ownLower, ownLower))
+                if (isEdgeOfFaceCut(pVals, f, ownLower, neiLower))
                 {
                     faceCutType_[faceI] = CUT;
                 }
@@ -148,7 +242,20 @@ void Foam::isoSurface::calcCutTypes
             forAll(pp, i)
             {
                 bool ownLower = (cVals[own[faceI]] < iso_);
-                bool neiLower = (cVals.boundaryField()[patchI][i] < iso_);
+
+                scalar nbrValue;
+                point nbrPoint;
+                getNeighbour
+                (
+                    boundaryRegion,
+                    cVals,
+                    own[faceI],
+                    faceI,
+                    nbrValue,
+                    nbrPoint
+                );
+
+                bool neiLower = (nbrValue < iso_);
 
                 if (ownLower != neiLower)
                 {
@@ -164,6 +271,7 @@ void Foam::isoSurface::calcCutTypes
                         faceCutType_[faceI] = CUT;
                     }
                 }
+
                 faceI++;
             }
         }
@@ -329,47 +437,6 @@ Foam::pointIndexHit Foam::isoSurface::collapseSurface
 }
 
 
-// Get neighbour value and position.
-void Foam::isoSurface::getNeighbour
-(
-    const labelList& boundaryRegion,
-    const volScalarField& cVals,
-    const label cellI,
-    const label faceI,
-    scalar& nbrValue,
-    point& nbrPoint
-) const
-{
-    const labelList& own = mesh_.faceOwner();
-    const labelList& nei = mesh_.faceNeighbour();
-
-    if (mesh_.isInternalFace(faceI))
-    {
-        label nbr = (own[faceI] == cellI ? nei[faceI] : own[faceI]);
-        nbrValue = cVals[nbr];
-        nbrPoint = mesh_.C()[nbr];
-    }
-    else
-    {
-        label bFaceI = faceI-mesh_.nInternalFaces();
-        label patchI = boundaryRegion[bFaceI];
-        label patchFaceI = faceI-mesh_.boundaryMesh()[patchI].start();
-
-        if (isA<emptyPolyPatch>(mesh_.boundaryMesh()[patchI]))
-        {
-            // Assume zero gradient
-            nbrValue = cVals[own[faceI]];
-            nbrPoint = mesh_.C().boundaryField()[patchI][patchFaceI];
-        }
-        else
-        {
-            nbrValue = cVals.boundaryField()[patchI][patchFaceI];
-            nbrPoint = mesh_.C().boundaryField()[patchI][patchFaceI];
-        }
-    }
-}
-
-
 // Determine per cell centre whether all the intersections get collapsed
 // to a single point
 void Foam::isoSurface::calcSnappedCc
@@ -383,6 +450,7 @@ void Foam::isoSurface::calcSnappedCc
 ) const
 {
     const pointField& pts = mesh_.points();
+    const pointField& cc = mesh_.cellCentres();
 
     snappedCc.setSize(mesh_.nCells());
     snappedCc = -1;
@@ -427,7 +495,7 @@ void Foam::isoSurface::calcSnappedCc
 
                 // From cc to neighbour cc.
                 s[2] = isoFraction(cVal, nbrValue);
-                pt[2] = (1.0-s[2])*mesh_.C()[cellI] + s[2]*nbrPoint;
+                pt[2] = (1.0-s[2])*cc[cellI] + s[2]*nbrPoint;
 
                 const face& f = mesh_.faces()[cFaces[cFaceI]];
 
@@ -436,12 +504,12 @@ void Foam::isoSurface::calcSnappedCc
                     // From cc to fp
                     label p0 = f[fp];
                     s[0] = isoFraction(cVal, pVals[p0]);
-                    pt[0] = (1.0-s[0])*mesh_.C()[cellI] + s[0]*pts[p0];
+                    pt[0] = (1.0-s[0])*cc[cellI] + s[0]*pts[p0];
 
                     // From cc to fp+1
                     label p1 = f[f.fcIndex(fp)];
                     s[1] = isoFraction(cVal, pVals[p1]);
-                    pt[1] = (1.0-s[1])*mesh_.C()[cellI] + s[1]*pts[p1];
+                    pt[1] = (1.0-s[1])*cc[cellI] + s[1]*pts[p1];
 
                     if
                     (
@@ -548,6 +616,7 @@ void Foam::isoSurface::calcSnappedPoint
 ) const
 {
     const pointField& pts = mesh_.points();
+    const pointField& cc = mesh_.cellCentres();
 
 
     const point greatPoint(VGREAT, VGREAT, VGREAT);
@@ -591,13 +660,14 @@ void Foam::isoSurface::calcSnappedPoint
 
         forAll(pFaces, pFaceI)
         {
+            // Create points for all intersections close to point
+            // (i.e. from pyramid edges)
+
             label faceI = pFaces[pFaceI];
             const face& f = mesh_.faces()[faceI];
             label own = mesh_.faceOwner()[faceI];
 
-            // Create points for all intersections close to point
+            // Get neighbour value
-            // (i.e. from pyramid edges)
-
             scalar nbrValue;
             point nbrPoint;
             getNeighbour
@@ -616,7 +686,7 @@ void Foam::isoSurface::calcSnappedPoint
 
             label fp = findIndex(f, pointI);
             s[0] = isoFraction(pVals[pointI], cVals[own]);
-            pt[0] = (1.0-s[0])*pts[pointI] + s[0]*mesh_.C()[own];
+            pt[0] = (1.0-s[0])*pts[pointI] + s[0]*cc[own];
 
             s[1] = isoFraction(pVals[pointI], nbrValue);
             pt[1] = (1.0-s[1])*pts[pointI] + s[1]*nbrPoint;
@@ -815,12 +885,6 @@ Foam::triSurface Foam::isoSurface::stitchTriPoints
         tris.transfer(dynTris);
     }
 
-    if (debug)
-    {
-        Pout<< "isoSurface : merged from " << nTris
-            << " down to " << tris.size() << " triangles." << endl;
-    }
-
 
 
     // Determine 'flat hole' situation (see RMT paper).
@@ -838,22 +902,29 @@ Foam::triSurface Foam::isoSurface::stitchTriPoints
         forAll(tris, triI)
         {
             const labelledTri& tri = tris[triI];
-
             const labelList& pFaces = pointFaces[tri[0]];
 
-            // Find the minimum of any duplicates
+            // Find the maximum of any duplicates. Maximum since the tris
+            // below triI
+            // get overwritten so we cannot use them in a comparison.
             label dupTriI = -1;
             forAll(pFaces, i)
             {
-                if (pFaces[i] < triI && tris[pFaces[i]] == tri)
+                label nbrTriI = pFaces[i];
+
+                if (nbrTriI > triI && (tris[nbrTriI] == tri))
                 {
-                    dupTriI = pFaces[i];
+                    //Pout<< "Duplicate : " << triI << " verts:" << tri
+                    //    << " to " << nbrTriI << " verts:" << tris[nbrTriI]
+                    //    << endl;
+                    dupTriI = nbrTriI;
+                    break;
                 }
             }
 
             if (dupTriI == -1)
             {
-                // There is no lower triangle
+                // There is no (higher numbered) duplicate triangle
                 label newTriI = newToOldTri.size();
                 newToOldTri.append(triI);
                 tris[newTriI] = tris[triI];
@@ -862,6 +933,49 @@ Foam::triSurface Foam::isoSurface::stitchTriPoints
 
         triMap.transfer(newToOldTri);
         tris.setSize(triMap.size());
+
+        if (debug)
+        {
+            Pout<< "isoSurface : merged from " << nTris
+                << " down to " << tris.size() << " unique triangles." << endl;
+        }
+
+
+        if (debug)
+        {
+            triSurface surf(tris, geometricSurfacePatchList(0), newPoints);
+
+            forAll(surf, faceI)
+            {
+                const labelledTri& f = surf[faceI];
+                const labelList& fFaces = surf.faceFaces()[faceI];
+
+                forAll(fFaces, i)
+                {
+                    label nbrFaceI = fFaces[i];
+
+                    if (nbrFaceI <= faceI)
+                    {
+                        // lower numbered faces already checked
+                        continue;
+                    }
+
+                    const labelledTri& nbrF = surf[nbrFaceI];
+
+                    if (f == nbrF)
+                    {
+                        FatalErrorIn("validTri(const triSurface&, const label)")
+                            << "Check : "
+                            << " triangle " << faceI << " vertices " << f
+                            << " fc:" << f.centre(surf.points())
+                            << " has the same vertices as triangle " << nbrFaceI
+                            << " vertices " << nbrF
+                            << " fc:" << nbrF.centre(surf.points())
+                            << abort(FatalError);
+                    }
+                }
+            }
+        }
     }
 
     return triSurface(tris, geometricSurfacePatchList(0), newPoints, true);
@@ -1461,6 +1575,12 @@ Foam::isoSurface::isoSurface
     {
         Pout<< "isoSurface:" << nl
             << "    isoField      : " << cVals.name() << nl
+            << "    cell min/max  : "
+            << min(cVals.internalField()) << " / "
+            << max(cVals.internalField()) << nl
+            << "    point min/max : "
+            << min(pVals) << " / "
+            << max(pVals) << nl
             << "    isoValue      : " << iso << nl
             << "    regularise    : " << regularise << nl
             << "    mergeTol      : " << mergeTol << nl
@@ -1489,23 +1609,13 @@ Foam::isoSurface::isoSurface
         }
     }
 
-
+    // Pre-calculate patch-per-face to avoid whichPatch call.
-    // Determine if any cut through face/cell
-    calcCutTypes(cVals, pVals);
-
-
-    // Determine if point is on boundary. Points on boundaries are never
-    // snapped. Coupled boundaries are handled explicitly so not marked here.
-    PackedBoolList isBoundaryPoint(mesh_.nPoints());
-
     labelList boundaryRegion(mesh_.nFaces()-mesh_.nInternalFaces());
 
     forAll(patches, patchI)
     {
         const polyPatch& pp = patches[patchI];
 
-        if (pp.coupled())
-        {
         label faceI = pp.start();
 
         forAll(pp, i)
@@ -1514,24 +1624,10 @@ Foam::isoSurface::isoSurface
             faceI++;
         }
     }
-        else
-        {
-            label faceI = pp.start();
 
-            forAll(pp, i)
-            {
-                boundaryRegion[faceI-mesh_.nInternalFaces()] = patchI;
 
-                const face& f = mesh_.faces()[faceI];
+    // Determine if any cut through face/cell
-
+    calcCutTypes(boundaryRegion, cVals, pVals);
-                forAll(f, fp)
-                {
-                    isBoundaryPoint.set(f[fp], 1);
-                }
-                faceI++;
-            }
-        }
-    }
 
 
     DynamicList<point> snappedPoints(nCutCells_);
@@ -1566,6 +1662,39 @@ Foam::isoSurface::isoSurface
 
     label nCellSnaps = snappedPoints.size();
 
+
+    // Determine if point is on boundary. Points on boundaries are never
+    // snapped. Coupled boundaries are handled explicitly so not marked here.
+    PackedBoolList isBoundaryPoint(mesh_.nPoints());
+
+    forAll(patches, patchI)
+    {
+        // Mark all boundary points that are not physically coupled (so anything
+        // but collocated coupled patches)
+        const polyPatch& pp = patches[patchI];
+
+        if
+        (
+           !pp.coupled()
+        || refCast<const coupledPolyPatch>(pp).separated()
+        )
+        {
+            label faceI = pp.start();
+
+            forAll(pp, i)
+            {
+                const face& f = mesh_.faces()[faceI];
+
+                forAll(f, fp)
+                {
+                    isBoundaryPoint.set(f[fp], 1);
+                }
+                faceI++;
+            }
+        }
+    }
+
+
     // Per point -1 or a point inside snappedPoints.
     labelList snappedPoint;
     if (regularise)
@@ -1594,6 +1723,48 @@ Foam::isoSurface::isoSurface
     }
 
 
+    // Generate field to interpolate. This is identical to the mesh.C()
+    // except on separated coupled patches.
+    slicedVolVectorField meshC
+    (
+        IOobject
+        (
+            "C",
+            mesh_.pointsInstance(),
+            mesh_.meshSubDir,
+            mesh_,
+            IOobject::NO_READ,
+            IOobject::NO_WRITE,
+            false
+        ),
+        mesh_,
+        dimLength,
+        mesh_.cellCentres(),
+        mesh_.faceCentres()
+    );
+    {
+        const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+        forAll(patches, patchI)
+        {
+            if
+            (
+                patches[patchI].coupled()
+             && refCast<const coupledPolyPatch>(patches[patchI]).separated()
+            )
+            {
+                fvPatchVectorField& pfld = const_cast<fvPatchVectorField&>
+                (
+                    meshC.boundaryField()[patchI]
+                );
+                pfld.operator==
+                (
+                    patches[patchI].patchSlice(mesh_.faceCentres())
+                );
+            }
+        }
+    }
+
+
     DynamicList<point> triPoints(nCutCells_);
     DynamicList<label> triMeshCells(nCutCells_);
 
@@ -1602,7 +1773,7 @@ Foam::isoSurface::isoSurface
         cVals,
         pVals,
 
-        mesh_.C(),
+        meshC,
         mesh_.points(),
 
         snappedPoints,
@@ -1616,7 +1787,8 @@ Foam::isoSurface::isoSurface
     if (debug)
     {
         Pout<< "isoSurface : generated " << triMeshCells.size()
-            << " unmerged triangles." << endl;
+            << " unmerged triangles from " << triPoints.size()
+            << " unmerged points." << endl;
     }
 
 

src/sampling/sampledSurface/isoSurface/isoSurface.H

@@ -135,9 +135,20 @@ class isoSurface
             const bool neiLower
         ) const;
 
+        void getNeighbour
+        (
+            const labelList& boundaryRegion,
+            const volScalarField& cVals,
+            const label cellI,
+            const label faceI,
+            scalar& nbrValue,
+            point& nbrPoint
+        ) const;
+
         //- Set faceCutType,cellCutType.
         void calcCutTypes
         (
+            const labelList& boundaryRegion,
             const volScalarField& cVals,
             const scalarField& pVals
         );
@@ -156,16 +167,6 @@ class isoSurface
             DynamicList<labelledTri, 64>& localTris
         );
 
-        void getNeighbour
-        (
-            const labelList& boundaryRegion,
-            const volScalarField& cVals,
-            const label cellI,
-            const label faceI,
-            scalar& nbrValue,
-            point& nbrPoint
-        ) const;
-
         //- Determine per cc whether all near cuts can be snapped to single
         //  point.
         void calcSnappedCc
@@ -193,37 +194,39 @@ class isoSurface
         template<class Type>
         Type generatePoint
         (
-            const DynamicList<Type>& snappedPoints,
-
             const scalar s0,
             const Type& p0,
-            const label p0Index,
+            const bool hasSnap0,
+            const Type& snapP0,
 
             const scalar s1,
             const Type& p1,
-            const label p1Index
+            const bool hasSnap1,
+            const Type& snapP1
         ) const;
 
         template<class Type>
         void generateTriPoints
         (
-            const DynamicList<Type>& snapped,
-
             const scalar s0,
             const Type& p0,
-            const label p0Index,
+            const bool hasSnap0,
+            const Type& snapP0,
 
             const scalar s1,
             const Type& p1,
-            const label p1Index,
+            const bool hasSnap1,
+            const Type& snapP1,
 
             const scalar s2,
             const Type& p2,
-            const label p2Index,
+            const bool hasSnap2,
+            const Type& snapP2,
 
             const scalar s3,
             const Type& p3,
-            const label p3Index,
+            const bool hasSnap3,
+            const Type& snapP3,
 
             DynamicList<Type>& points
         ) const;
@@ -244,7 +247,8 @@ class isoSurface
 
             const scalar neiVal,
             const Type& neiPt,
-            const label neiSnap,
+            const bool hasNeiSnap,
+            const Type& neiSnapPt,
 
             DynamicList<Type>& triPoints,
             DynamicList<label>& triMeshCells

src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C

@@ -33,15 +33,15 @@ License
 template<class Type>
 Type Foam::isoSurface::generatePoint
 (
-    const DynamicList<Type>& snappedPoints,
-
     const scalar s0,
     const Type& p0,
-    const label p0Index,
+    const bool hasSnap0,
+    const Type& snapP0,
 
     const scalar s1,
     const Type& p1,
-    const label p1Index
+    const bool hasSnap1,
+    const Type& snapP1
 ) const
 {
     scalar d = s1-s0;
@@ -50,13 +50,13 @@ Type Foam::isoSurface::generatePoint
     {
         scalar s = (iso_-s0)/d;
 
-        if (s >= 0.5 && s <= 1 && p1Index != -1)
+        if (hasSnap1 && s >= 0.5 && s <= 1)
         {
-            return snappedPoints[p1Index];
+            return snapP1;
         }
-        else if (s >= 0.0 && s <= 0.5 && p0Index != -1)
+        else if (hasSnap0 && s >= 0.0 && s <= 0.5)
         {
-            return snappedPoints[p0Index];
+            return snapP0;
         }
         else
         {
@@ -75,23 +75,25 @@ Type Foam::isoSurface::generatePoint
 template<class Type>
 void Foam::isoSurface::generateTriPoints
 (
-    const DynamicList<Type>& snapped,
-
     const scalar s0,
     const Type& p0,
-    const label p0Index,
+    const bool hasSnap0,
+    const Type& snapP0,
 
     const scalar s1,
     const Type& p1,
-    const label p1Index,
+    const bool hasSnap1,
+    const Type& snapP1,
 
     const scalar s2,
     const Type& p2,
-    const label p2Index,
+    const bool hasSnap2,
+    const Type& snapP2,
 
     const scalar s3,
     const Type& p3,
-    const label p3Index,
+    const bool hasSnap3,
+    const Type& snapP3,
 
     DynamicList<Type>& points
 ) const
@@ -123,29 +125,55 @@ void Foam::isoSurface::generateTriPoints
 
         case 0x0E:
         case 0x01:
-            points.append(generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index));
+            points.append
-            points.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
+            (
-            points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
+                generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1)
+            );
+            points.append
+            (
+                generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2)
+            );
+            points.append
+            (
+                generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
+            );
         break;
 
         case 0x0D:
         case 0x02:
-            points.append(generatePoint(snapped,s1,p1,p1Index,s0,p0,p0Index));
+            points.append
-            points.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
+            (
-            points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
+                generatePoint(s1,p1,hasSnap1,snapP1,s0,p0,hasSnap0,snapP0)
+            );
+            points.append
+            (
+                generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3)
+            );
+            points.append
+            (
+                generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
+            );
         break;
 
         case 0x0C:
         case 0x03:
         {
-            Type tp1 = generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index);
+            Type tp1 =
-            Type tp2 = generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index);
+                generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2);
+            Type tp2 =
+                generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3);
 
-            points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
+            points.append
+            (
+                generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
+            );
             points.append(tp1);
             points.append(tp2);
             points.append(tp2);
-            points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
+            points.append
+            (
+                generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
+            );
             points.append(tp1);
         }
         break;
@@ -153,23 +181,40 @@ void Foam::isoSurface::generateTriPoints
         case 0x0B:
         case 0x04:
         {
-            points.append(generatePoint(snapped,s2,p2,p2Index,s0,p0,p0Index));
+            points.append
-            points.append(generatePoint(snapped,s2,p2,p2Index,s1,p1,p1Index));
+            (
-            points.append(generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index));
+                generatePoint(s2,p2,hasSnap2,snapP2,s0,p0,hasSnap0,snapP0)
+            );
+            points.append
+            (
+                generatePoint(s2,p2,hasSnap2,snapP2,s1,p1,hasSnap1,snapP1)
+            );
+            points.append
+            (
+                generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3)
+            );
         }
         break;
 
         case 0x0A:
         case 0x05:
         {
-            Type tp0 = generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
+            Type tp0 =
-            Type tp1 = generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
+                generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1);
+            Type tp1 =
+                generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3);
 
             points.append(tp0);
             points.append(tp1);
-            points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
+            points.append
+            (
+                generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
+            );
             points.append(tp0);
-            points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
+            points.append
+            (
+                generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
+            );
             points.append(tp1);
         }
         break;
@@ -177,23 +222,40 @@ void Foam::isoSurface::generateTriPoints
         case 0x09:
         case 0x06:
         {
-            Type tp0 = generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
+            Type tp0 =
-            Type tp1 = generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
+                generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1);
+            Type tp1 =
+                generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3);
 
             points.append(tp0);
-            points.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
+            points.append
+            (
+                generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3)
+            );
             points.append(tp1);
             points.append(tp0);
-            points.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
+            points.append
+            (
+                generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2)
+            );
             points.append(tp1);
         }
         break;
 
         case 0x07:
         case 0x08:
-            points.append(generatePoint(snapped,s3,p3,p3Index,s0,p0,p0Index));
+            points.append
-            points.append(generatePoint(snapped,s3,p3,p3Index,s2,p2,p2Index));
+            (
-            points.append(generatePoint(snapped,s3,p3,p3Index,s1,p1,p1Index));
+                generatePoint(s3,p3,hasSnap3,snapP3,s0,p0,hasSnap0,snapP0)
+            );
+            points.append
+            (
+                generatePoint(s3,p3,hasSnap3,snapP3,s2,p2,hasSnap2,snapP2)
+            );
+            points.append
+            (
+                generatePoint(s3,p3,hasSnap3,snapP3,s1,p1,hasSnap1,snapP1)
+            );
         break;
     }
 }
@@ -215,7 +277,8 @@ Foam::label Foam::isoSurface::generateFaceTriPoints
 
     const scalar neiVal,
     const Type& neiPt,
-    const label neiSnap,
+    const bool hasNeiSnap,
+    const Type& neiSnapPt,
 
     DynamicList<Type>& triPoints,
     DynamicList<label>& triMeshCells
@@ -234,23 +297,37 @@ Foam::label Foam::isoSurface::generateFaceTriPoints
 
         generateTriPoints
         (
-            snappedPoints,
-
             pVals[pointI],
             pCoords[pointI],
-            snappedPoint[pointI],
+            snappedPoint[pointI] != -1,
+            (
+                snappedPoint[pointI] != -1
+              ? snappedPoints[snappedPoint[pointI]]
+              : pTraits<Type>::zero
+            ),
 
             pVals[nextPointI],
             pCoords[nextPointI],
-            snappedPoint[nextPointI],
+            snappedPoint[nextPointI] != -1,
+            (
+                snappedPoint[nextPointI] != -1
+              ? snappedPoints[snappedPoint[nextPointI]]
+              : pTraits<Type>::zero
+            ),
 
             cVals[own],
             cCoords[own],
-            snappedCc[own],
+            snappedCc[own] != -1,
+            (
+                snappedCc[own] != -1
+              ? snappedPoints[snappedCc[own]]
+              : pTraits<Type>::zero
+            ),
 
             neiVal,
             neiPt,
-            neiSnap,
+            hasNeiSnap,
+            neiSnapPt,
 
             triPoints
         );
@@ -311,25 +388,6 @@ void Foam::isoSurface::generateTriPoints
             << abort(FatalError);
     }
 
-    // Determine neighbouring snap status
-    labelList neiSnappedCc(mesh_.nFaces()-mesh_.nInternalFaces(), -1);
-    forAll(patches, patchI)
-    {
-        const polyPatch& pp = patches[patchI];
-
-        if (pp.coupled())
-        {
-            label faceI = pp.start();
-            forAll(pp, i)
-            {
-                neiSnappedCc[faceI-mesh_.nInternalFaces()] =
-                    snappedCc[own[faceI]];
-                faceI++;
-            }
-        }
-    }
-    syncTools::swapBoundaryFaceList(mesh_, neiSnappedCc, false);
-
 
 
     // Generate triangle points
@@ -356,7 +414,12 @@ void Foam::isoSurface::generateTriPoints
 
                 cVals[nei[faceI]],
                 cCoords[nei[faceI]],
-                snappedCc[nei[faceI]],
+                snappedCc[nei[faceI]] != -1,
+                (
+                    snappedCc[nei[faceI]] != -1
+                  ? snappedPoints[snappedCc[nei[faceI]]]
+                  : pTraits<Type>::zero
+                ),
 
                 triPoints,
                 triMeshCells
@@ -365,13 +428,45 @@ void Foam::isoSurface::generateTriPoints
     }
 
 
+    // Determine neighbouring snap status
+    boolList neiSnapped(mesh_.nFaces()-mesh_.nInternalFaces(), false);
+    List<Type> neiSnappedPoint(neiSnapped.size(), pTraits<Type>::zero);
     forAll(patches, patchI)
     {
         const polyPatch& pp = patches[patchI];
 
-        if (isA<processorPolyPatch>(pp))
+        if (pp.coupled())
         {
-            if (refCast<const processorPolyPatch>(pp).owner())
+            label faceI = pp.start();
+            forAll(pp, i)
+            {
+                label bFaceI = faceI-mesh_.nInternalFaces();
+                label snappedIndex = snappedCc[own[faceI]];
+
+                if (snappedIndex != -1)
+                {
+                    neiSnapped[bFaceI] = true;
+                    neiSnappedPoint[bFaceI] = snappedPoints[snappedIndex];
+                }
+                faceI++;
+            }
+        }
+    }
+    syncTools::swapBoundaryFaceList(mesh_, neiSnapped, false);
+    syncTools::swapBoundaryFaceList(mesh_, neiSnappedPoint, false);
+
+
+    forAll(patches, patchI)
+    {
+        const polyPatch& pp = patches[patchI];
+
+        if
+        (
+            isA<processorPolyPatch>(pp)
+         && !refCast<const processorPolyPatch>(pp).separated()
+        )
+        {
+            //if (refCast<const processorPolyPatch>(pp).owner())
             {
                 label faceI = pp.start();
 
@@ -379,6 +474,18 @@ void Foam::isoSurface::generateTriPoints
                 {
                     if (faceCutType_[faceI] != NOTCUT)
                     {
+                        label bFaceI = faceI-mesh_.nInternalFaces();
+                        if
+                        (
+                            neiSnapped[bFaceI]
+                         && (neiSnappedPoint[bFaceI]==pTraits<Type>::zero)
+                        )
+                        {
+                            FatalErrorIn("isoSurface::generateTriPoints(..)")
+                                << "problem:" << abort(FatalError);
+                        }
+
+
                         generateFaceTriPoints
                         (
                             cVals,
@@ -394,7 +501,8 @@ void Foam::isoSurface::generateTriPoints
 
                             cVals.boundaryField()[patchI][i],
                             cCoords.boundaryField()[patchI][i],
-                            neiSnappedCc[faceI-mesh_.nInternalFaces()],
+                            neiSnapped[faceI-mesh_.nInternalFaces()],
+                            neiSnappedPoint[faceI-mesh_.nInternalFaces()],
 
                             triPoints,
                             triMeshCells
@@ -428,7 +536,8 @@ void Foam::isoSurface::generateTriPoints
 
                         cVals[own[faceI]],
                         cCoords.boundaryField()[patchI][i],
-                        -1, // fc not snapped
+                        false,  // fc not snapped
+                        pTraits<Type>::zero,
 
                         triPoints,
                         triMeshCells
@@ -460,7 +569,8 @@ void Foam::isoSurface::generateTriPoints
 
                         cVals.boundaryField()[patchI][i],
                         cCoords.boundaryField()[patchI][i],
-                        -1, // fc not snapped
+                        false,  // fc not snapped
+                        pTraits<Type>::zero,
 
                         triPoints,
                         triMeshCells
@@ -482,8 +592,8 @@ void Foam::isoSurface::generateTriPoints
 //{
 //    return tmp<Field<Type> >(new Field<Type>(vField, meshCells()));
 //}
-//
+
-//
+
 template <class Type>
 Foam::tmp<Foam::Field<Type> >
 Foam::isoSurface::interpolate

src/thermophysicalModels/liquids/liquid/liquid.H

@@ -41,7 +41,6 @@ SourceFiles
 #include "typeInfo.H"
 #include "autoPtr.H"
 #include "runTimeSelectionTables.H"
-//#include "eos.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -52,8 +51,6 @@ namespace Foam
                            Class liquid Declaration
 \*---------------------------------------------------------------------------*/
 
-    //class eos;
-
 class liquid
 {
     // Private data
@@ -85,13 +82,12 @@ class liquid
         //- Dipole moment []
         scalar dipm_;
 
-        //- Pitzer's acentric factor []
+        //- Pitzer's accentric factor []
         scalar omega_;
 
         //- Solubility parameter [(J/m^3)^0.5]
         scalar delta_;
 
-    //eos state_;
 
 public:
 
@@ -174,7 +170,7 @@ public:
 
     // Member Functions
 
-        // Phisical constants which define the specie
+        // Physical constants which define the specie
 
             //- Molecular weight [kg/kmol]
             scalar W() const

src/turbulenceModels/LES/LESdeltas/cubeRootVolDelta/cubeRootVolDelta.C

@@ -42,8 +42,7 @@ addToRunTimeSelectionTable(LESdelta, cubeRootVolDelta, dictionary);
 
 void cubeRootVolDelta::calcDelta()
 {
-    const Vector<label>& directions = mesh().directions();
+    label nD = mesh().nGeometricD();
-    label nD = (directions.nComponents + cmptSum(directions))/2;
 
     if (nD == 3)
     {
@@ -55,14 +54,15 @@ void cubeRootVolDelta::calcDelta()
             << "Case is 2D, LES is not strictly applicable\n"
             << endl;
 
+        const Vector<label>& directions = mesh().geometricD();
+
         scalar thickness = 0.0;
         for (direction dir=0; dir<directions.nComponents; dir++)
         {
             if (directions[dir] == -1)
             {
-                boundBox bb(mesh().points(), false);
+                thickness = mesh().bounds().span()[dir];
-
+                break;
-                thickness = bb.span()[dir];
             }
         }
 

src/turbulenceModels/incompressible/LES/SpalartAllmarasIDDES/IDDESDelta/IDDESDelta.C

@@ -41,8 +41,7 @@ namespace Foam
 
 void Foam::IDDESDelta::calcDelta()
 {
-    const Vector<label>& directions = mesh().directions();
+    label nD = mesh().nGeometricD();
-    label nD = (directions.nComponents + cmptSum(directions))/2;
 
     // initialise hwn as wall distance
     volScalarField hwn = wallDist(mesh()).y();

tutorials/combustion/XiFoam/Allrun

@@ -3,7 +3,7 @@
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
 # Get application name from directory
-application=${PWD##*/}
+application="XiFoam"
 
 setControlDict () {
     controlDict="system/controlDict"
@@ -20,6 +20,7 @@ setControlDict () {
 
 # Clone case
 cloneCase moriyoshiHomogeneous moriyoshiHomogeneousPart2
+
 # Modify and execute
 cd moriyoshiHomogeneousPart2
     cp -r ../moriyoshiHomogeneous/0.005 .

tutorials/combustion/XiFoam/ras/moriyoshiHomogeneous/0/epsilon

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      epsilon;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -46,4 +47,5 @@ boundaryField
     }
 }
 
+
 // ************************************************************************* //

tutorials/combustion/XiFoam/ras/moriyoshiHomogeneous/0/k

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      k;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -46,4 +47,5 @@ boundaryField
     }
 }
 
+
 // ************************************************************************* //

tutorials/combustion/XiFoam/ras/moriyoshiHomogeneous/0/mut

@@ -0,0 +1,51 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      mut;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [ 1 -1 -1 0 0 0 0 ];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    left
+    {
+        type            symmetryPlane;
+    }
+
+    right
+    {
+        type            symmetryPlane;
+    }
+
+    top
+    {
+        type            symmetryPlane;
+    }
+
+    bottom
+    {
+        type            symmetryPlane;
+    }
+
+    frontAndBack
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoPisoFoam/ras/cavity/0/R

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volSymmTensorField;
+    location    "0";
     object      R;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -22,12 +23,14 @@ boundaryField
 {
     movingWall
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform ( 0 0 0 0 0 0 );
     }
 
     fixedWalls
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform ( 0 0 0 0 0 0 );
     }
 
     frontAndBack
@@ -36,4 +39,5 @@ boundaryField
     }
 }
 
+
 // ************************************************************************* //

tutorials/compressible/rhoPisoFoam/ras/cavity/0/epsilon

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      epsilon;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -22,12 +23,14 @@ boundaryField
 {
     movingWall
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     fixedWalls
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     frontAndBack
@@ -36,4 +39,5 @@ boundaryField
     }
 }
 
+
 // ************************************************************************* //

tutorials/compressible/rhoPisoFoam/ras/cavity/0/k

@@ -24,13 +24,15 @@ boundaryField
     movingWall
     {
         type            kQRWallFunction;
-        value           uniform 0.00325;
+        value           uniform 0;
     }
+
     fixedWalls
     {
         type            kQRWallFunction;
-        value           uniform 0.00325;
+        value           uniform 0;
     }
+
     frontAndBack
     {
         type            empty;

tutorials/compressible/rhoPisoFoam/ras/cavity/0/mut

@@ -26,11 +26,13 @@ boundaryField
         type            mutWallFunction;
         value           uniform 0;
     }
+
     fixedWalls
     {
         type            mutWallFunction;
         value           uniform 0;
     }
+
     frontAndBack
     {
         type            empty;

tutorials/compressible/rhoPisoFoam/ras/cavity/0/omega

@@ -24,13 +24,15 @@ boundaryField
     movingWall
     {
         type            omegaWallFunction;
-        value           uniform 2.6;
+        value           uniform 0;
     }
+
     fixedWalls
     {
         type            omegaWallFunction;
-        value           uniform 2.6;
+        value           uniform 0;
     }
+
     frontAndBack
     {
         type            empty;

tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/0/epsilon

@@ -1,55 +1,64 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
     version     2.0;
-    format      ascii;
+    format      binary;
     class       volScalarField;
+    location    "0";
     object      epsilon;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
 dimensions      [ 0 2 -3 0 0 0 0 ];
 
 internalField   uniform 200;
 
 boundaryField
 {
-
     front
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
+
     back
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
+
     wall
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
+
     porosityWall
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     inlet
     {
         type            turbulentMixingLengthDissipationRateInlet;
         mixingLength    0.005;
-        value           $internalField;
+        value           uniform 200;
     }
 
     outlet
     {
         type            inletOutlet;
-        inletValue      $internalField;
+        inletValue      uniform 200;
-        value           $internalField;
+        value           uniform 200;
+    }
 }
 
-}
+
 // ************************************************************************* //

tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/0/k

@@ -1,55 +1,64 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
     version     2.0;
-    format      ascii;
+    format      binary;
     class       volScalarField;
+    location    "0";
     object      k;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
 dimensions      [ 0 2 -2 0 0 0 0 ];
 
 internalField   uniform 1;
 
 boundaryField
 {
-
     front
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
+
     back
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
+
     wall
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
+
     porosityWall
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 
     inlet
     {
         type            turbulentIntensityKineticEnergyInlet;
         intensity       0.05;
-        value           $internalField;
+        value           uniform 1;
     }
 
     outlet
     {
         type            inletOutlet;
-        inletValue      $internalField;
+        inletValue      uniform 1;
-        value           $internalField;
+        value           uniform 1;
+    }
 }
 
-}
+
 // ************************************************************************* //

tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/0/mut

@@ -0,0 +1,62 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      binary;
+    class       volScalarField;
+    location    "0";
+    object      mut;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [ 1 -1 -1 0 0 0 0 ];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    front
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    back
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    wall
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    porosityWall
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    inlet
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+
+    outlet
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/heatTransfer/buoyantPisoFoam/hotRoom/0/epsilon

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      epsilon;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -22,18 +23,22 @@ boundaryField
 {
     floor
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     ceiling
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     fixedWalls
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 }
 
+
 // ************************************************************************* //

tutorials/heatTransfer/buoyantPisoFoam/hotRoom/0/k

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      k;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -22,18 +23,22 @@ boundaryField
 {
     floor
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 
     ceiling
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 
     fixedWalls
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 }
 
+
 // ************************************************************************* //

tutorials/heatTransfer/buoyantPisoFoam/hotRoom/0/mut

@@ -0,0 +1,44 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      mut;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [ 1 -1 -1 0 0 0 0 ];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    floor
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    ceiling
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    fixedWalls
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/heatTransfer/buoyantPisoFoam/hotRoom/Allrun

@@ -3,7 +3,7 @@
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
 # Get application name from directory
-application="buoyantFoam"
+application="buoyantPisoFoam"
 
 compileApplication setHotRoom
 runApplication blockMesh

tutorials/heatTransfer/buoyantSimpleRadiationFoam/hotRadiationRoom/0/epsilon

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      epsilon;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -22,23 +23,28 @@ boundaryField
 {
     floor
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     ceiling
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     fixedWalls
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     box
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 }
 
+
 // ************************************************************************* //

tutorials/heatTransfer/buoyantSimpleRadiationFoam/hotRadiationRoom/0/k

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      k;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -22,23 +23,28 @@ boundaryField
 {
     floor
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 
     ceiling
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 
     fixedWalls
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 
     box
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 }
 
+
 // ************************************************************************* //

tutorials/heatTransfer/buoyantSimpleRadiationFoam/hotRadiationRoom/0/mut

@@ -0,0 +1,50 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      mut;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [ 1 -1 -1 0 0 0 0 ];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    box
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    floor
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    ceiling
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+
+    fixedWalls
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/0/alphat

@@ -0,0 +1,67 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0.001";
+    object      alphat;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [1 -1 -1 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    maxY
+    {
+        type            alphatWallFunction;
+        value           uniform 0;
+    }
+    minX
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    maxX
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    minZ
+    {
+        type            alphatWallFunction;
+        value           uniform 0;
+    }
+    maxZ
+    {
+        type            alphatWallFunction;
+        value           uniform 0;
+    }
+    topAir_to_leftSolid
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    topAir_to_heater
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    topAir_to_rightSolid
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/0/epsilon

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0.001";
     object      epsilon;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -20,6 +21,11 @@ internalField   uniform 0.01;
 
 boundaryField
 {
+    maxY
+    {
+        type            epsilonWallFunction;
+        value           uniform 0.01;
+    }
     minX
     {
         type            fixedValue;
@@ -29,24 +35,29 @@ boundaryField
     {
         type            zeroGradient;
     }
-
-    minY
-    {
-        type            zeroGradient;
-    }
-    maxY
-    {
-        type            zeroGradient;
-    }
-
     minZ
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0.01;
     }
     maxZ
+    {
+        type            epsilonWallFunction;
+        value           uniform 0.01;
+    }
+    topAir_to_leftSolid
+    {
+        type            zeroGradient;
+    }
+    topAir_to_heater
+    {
+        type            zeroGradient;
+    }
+    topAir_to_rightSolid
     {
         type            zeroGradient;
     }
 }
 
+
 // ************************************************************************* //

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/0/k

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0.001";
     object      k;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -20,6 +21,11 @@ internalField   uniform 0.1;
 
 boundaryField
 {
+    maxY
+    {
+        type            kQRWallFunction;
+        value           uniform 0.1;
+    }
     minX
     {
         type            fixedValue;
@@ -29,24 +35,29 @@ boundaryField
     {
         type            zeroGradient;
     }
-
-    minY
-    {
-        type            zeroGradient;
-    }
-    maxY
-    {
-        type            zeroGradient;
-    }
-
     minZ
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0.1;
     }
     maxZ
+    {
+        type            kQRWallFunction;
+        value           uniform 0.1;
+    }
+    topAir_to_leftSolid
+    {
+        type            zeroGradient;
+    }
+    topAir_to_heater
+    {
+        type            zeroGradient;
+    }
+    topAir_to_rightSolid
     {
         type            zeroGradient;
     }
 }
 
+
 // ************************************************************************* //

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/0/mut

@@ -0,0 +1,67 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0.001";
+    object      mut;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [1 -1 -1 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    maxY
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+    minX
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    maxX
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    minZ
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+    maxZ
+    {
+        type            mutWallFunction;
+        value           uniform 0;
+    }
+    topAir_to_leftSolid
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    topAir_to_heater
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    topAir_to_rightSolid
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/Allrun

@@ -1,4 +1,5 @@
 #!/bin/sh
+
 # Source tutorial run functions
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
@@ -20,7 +21,7 @@ done
 # remove fluid fields from solid regions (important for post-processing)
 for i in heater leftSolid rightSolid
 do
-   rm -f 0*/$i/{epsilon,k,p,pd,U}
+   rm -f 0*/$i/{mut,alphat,epsilon,k,p,pd,U}
 done
 
 # remove solid fields from fluid regions (important for post-processing)

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/constant/topAir/turbulenceProperties

@@ -0,0 +1 @@
+../bottomAir/turbulenceProperties

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/system/README

@@ -0,0 +1,3 @@
+fvSchemes and fvSolution are not used here
+- added so that pre-processing activities can proceed
+

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/system/fvSchemes

@@ -9,31 +9,10 @@ FoamFile
 {
     version     2.0;
     format      ascii;
-    class       volScalarField;
+    class       dictionary;
-    location    "0";
+    object      fvSchemes;
-    object      k.old;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-dimensions      [0 2 -2 0 0 0 0];
-
-internalField   uniform 0.1;
-
-boundaryField
-{
-    floor
-    {
-        type            zeroGradient;
-    }
-    ceiling
-    {
-        type            zeroGradient;
-    }
-    fixedWalls
-    {
-        type            zeroGradient;
-    }
-}
-
 
 // ************************************************************************* //

tutorials/heatTransfer/chtMultiRegionFoam/multiRegionHeater/system/fvSolution

@@ -9,31 +9,9 @@ FoamFile
 {
     version     2.0;
     format      ascii;
-    class       volScalarField;
+    class       dictionary;
-    location    "0";
+    object      fvSolution;
-    object      epsilon.old;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-dimensions      [0 2 -3 0 0 0 0];
-
-internalField   uniform 0.01;
-
-boundaryField
-{
-    floor
-    {
-        type            zeroGradient;
-    }
-    ceiling
-    {
-        type            zeroGradient;
-    }
-    fixedWalls
-    {
-        type            zeroGradient;
-    }
-}
-
-
 // ************************************************************************* //

tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/0/epsilon

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      epsilon;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -22,12 +23,14 @@ boundaryField
 {
     rotor
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     stator
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 
     front
@@ -41,4 +44,5 @@ boundaryField
     }
 }
 
+
 // ************************************************************************* //

tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/0/k

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      k;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -22,12 +23,14 @@ boundaryField
 {
     rotor
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 
     stator
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 
     front
@@ -41,4 +44,5 @@ boundaryField
     }
 }
 
+
 // ************************************************************************* //

tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/0/nut

@@ -0,0 +1,48 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      nut;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [ 0 2 -1 0 0 0 0 ];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    rotor
+    {
+        type            nutWallFunction;
+        value           uniform 0;
+    }
+
+    stator
+    {
+        type            nutWallFunction;
+        value           uniform 0;
+    }
+
+    front
+    {
+        type            empty;
+    }
+
+    back
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/incompressible/channelFoam/channel395/Allrun

@@ -3,8 +3,7 @@
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
 # Get application directory
-parentDir=`dirname $PWD`
+application="channelFoam"
-application=`basename $parentDir`
 
 runApplication blockMesh
 runApplication $application

tutorials/incompressible/channelFoam/channel395/constant/postChannelDict

@@ -15,10 +15,14 @@ FoamFile
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+// Seed patches to start layering from
 patches         ( bottomWall );
 
+// Direction in which the layers are
 component       y;
 
+// Is the mesh symmetric? If so average(symmetric fields) or
+// subtract(asymmetric) contributions from both halves
 symmetric       true;
 
 

tutorials/incompressible/channelFoam/channel395/system/controlDict

@@ -70,5 +70,4 @@ functions
     }
 }
 
-
 // ************************************************************************* //

tutorials/incompressible/pimpleFoam/t-junction/0/epsilon

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      epsilon;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -23,7 +24,7 @@ boundaryField
     inlet
     {
         type            turbulentMixingLengthDissipationRateInlet;
-        mixingLength    0.01;         // 1cm - half channel height
+        mixingLength    0.01;
         value           uniform 1;
     }
 
@@ -41,7 +42,8 @@ boundaryField
 
     defaultFaces
     {
-        type            zeroGradient;
+        type            epsilonWallFunction;
+        value           uniform 0;
     }
 }
 

tutorials/incompressible/pimpleFoam/t-junction/0/k

@@ -1,8 +1,8 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5                                   |
+|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      k;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -23,7 +24,7 @@ boundaryField
     inlet
     {
         type            turbulentIntensityKineticEnergyInlet;
-        intensity       0.05;           // 5% turbulence
+        intensity       0.05;
         value           uniform 1;
     }
 
@@ -41,8 +42,10 @@ boundaryField
 
     defaultFaces
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
+        value           uniform 0;
     }
 }
 
+
 // ************************************************************************* //

tutorials/incompressible/pimpleFoam/t-junction/0/nut

@@ -0,0 +1,50 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      nut;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [ 0 2 -1 0 0 0 0 ];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    inlet
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+
+    outlet1
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+
+    outlet2
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+
+    defaultFaces
+    {
+        type            nutWallFunction;
+        value           uniform 0;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/incompressible/pisoFoam/ras/cavity/0/R

@@ -22,12 +22,12 @@ boundaryField
 {
     movingWall
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
     }
 
     fixedWalls
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
     }
 
     frontAndBack

tutorials/incompressible/simpleFoam/pitzDaily/0/R

@@ -33,12 +33,12 @@ boundaryField
 
     upperWall
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
     }
 
     lowerWall
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
     }
 
     frontAndBack

tutorials/incompressible/simpleFoam/pitzDailyExptInlet/0/R

@@ -33,12 +33,12 @@ boundaryField
 
     upperWall
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
     }
 
     lowerWall
     {
-        type            zeroGradient;
+        type            kQRWallFunction;
     }
 
     frontAndBack

tutorials/mesh/snappyHexMesh/iglooWithFridges/system/snappyHexMeshDict

@@ -15,12 +15,18 @@ FoamFile
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+// Which of the steps to run
 castellatedMesh true;
-
 snap            true;
-
 addLayers       true;
 
+
+// Geometry. Definition of all surfaces. All surfaces are of class
+// searchableSurface.
+// Surfaces are used
+// - to specify refinement for any mesh cell intersecting it
+// - to specify refinement for any mesh cell inside/outside/near
+// - to 'snap' the mesh boundary to the surface
 geometry
 {
     fridgeA
@@ -45,17 +51,68 @@ geometry
     }
 }
 
+
+
+// Settings for the castellatedMesh generation.
 castellatedMeshControls
 {
+
+    // Refinement parameters
+    // ~~~~~~~~~~~~~~~~~~~~~
+
+    // 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.
     maxLocalCells   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.
     maxGlobalCells  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)
     minRefinementCells 0;
+
+    // Number of buffer layers between different levels.
+    // 1 means normal 2:1 refinement restriction, larger means slower
+    // refinement.
     nCellsBetweenLevels 1;
-    features        ( );
+
+
+
+    // Explicit feature edge refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // Specifies a level for any cell intersected by its edges.
+    // This is a featureEdgeMesh, read from constant/triSurface for now.
+    features
+    (
+        //{
+        //    file "someLine.eMesh";
+        //    level 2;
+        //}
+    );
+
+
+
+    // Surface based refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // Specifies two levels for every surface. The first is the minimum level,
+    // every cell intersecting a surface gets refined up to the minimum level.
+    // The second level is the maximum level. Cells that 'see' multiple
+    // intersections where the intersections make an
+    // angle > resolveFeatureAngle get refined up to the maximum level.
     refinementSurfaces
     {
         fridgeA
         {
+            // Surface-wise min and max refinement level
             level           ( 2 2 );
         }
 
@@ -71,23 +128,64 @@ castellatedMeshControls
     }
 
     resolveFeatureAngle 60;
+
+    // Region-wise refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~
+
+    // Specifies refinement level for cells in relation to a surface. One of
+    // three modes
+    // - distance. 'levels' specifies per distance to the surface the
+    //   wanted refinement level. The distances need to be specified in
+    //   descending order.
+    // - inside. 'levels' is only one entry and only the level is used. All
+    //   cells inside the surface get refined up to the level. The surface
+    //   needs to be closed for this to be possible.
+    // - outside. Same but cells outside.
+
     refinementRegions
     {
     }
 
+
+    // Mesh selection
+    // ~~~~~~~~~~~~~~
+
+    // After refinement patches get added for all refinementSurfaces and
+    // all cells intersecting the surfaces get put into these patches. The
+    // section reachable from the locationInMesh is kept.
+    // NOTE: This point should never be on a face, always inside a cell, even
+    // after refinement.
     locationInMesh  ( 3 0.28 0.43 );
 }
 
+
+
+// Settings for the snapping.
 snapControls
 {
+    //- 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.
     tolerance       4;
+
+    //- Number of mesh displacement relaxation iterations.
     nSolveIter      30;
+
+    //- Maximum number of snapping relaxation iterations. Should stop
+    //  before upon reaching a correct mesh.
     nRelaxIter      5;
 }
 
+
+
+// Settings for the layer addition.
 addLayersControls
 {
+    // Per final patch (so not geometry!) the layer information
     layers
     {
         fridgeA_region0
@@ -106,41 +204,129 @@ addLayersControls
         }
     }
 
+
+    // Expansion factor for layer mesh
     expansionRatio  1;
+
+    //- Wanted thickness of final added cell layer. If multiple layers
+    //  is the
+    //  thickness of the layer furthest away from the wall.
+    //  Relative to undistorted size of cell outside layer.
     finalLayerRatio 0.5;
+
+    //- Minimum thickness of cell layer. If for any reason layer
+    //  cannot be above minThickness do not add layer.
+    //  Relative to undistorted size of cell outside layer.
     minThickness    0.25;
+
+    //- If points get not extruded do nGrow layers of connected faces that are
+    //  also not grown. This helps convergence of the layer addition process
+    //  close to features.
     nGrow           0;
+
+
+    // Advanced settings
+
+    //- 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.
     nRelaxIter      5;
+
+    // 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;
 }
 
+
+
+// Generic mesh quality settings. At any undoable phase these determine
+// where to undo.
 meshQualityControls
 {
+    //- 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 a sensible fraction of the smallest cell volume expected.
+    //  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;
 }
 
+
+
+// Advanced
+
+// Flags for optional output
+// 0 : only write final meshes
+// 1 : write intermediate meshes
+// 2 : write volScalarField with cellLevel for postprocessing
+// 4 : write current intersections as .obj files
 debug           0;
 
+
+// Merge tolerance. Is fraction of overall bounding box of initial mesh.
+// Note: the write tolerance needs to be higher than this.
 mergeTolerance  1e-06;
 
 

tutorials/mesh/snappyHexMesh/motorBike/system/snappyHexMeshDict

@@ -13,14 +13,21 @@ FoamFile
     location    "system";
     object      snappyHexMeshDict;
 }
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+// Which of the steps to run
 castellatedMesh true;
-
 snap            true;
-
 addLayers       true;
 
+
+// Geometry. Definition of all surfaces. All surfaces are of class
+// searchableSurface.
+// Surfaces are used
+// - to specify refinement for any mesh cell intersecting it
+// - to specify refinement for any mesh cell inside/outside/near
+// - to 'snap' the mesh boundary to the surface
 geometry
 {
     motorBike.stl
@@ -37,22 +44,89 @@ geometry
     }
 }
 
+
+// Settings for the castellatedMesh generation.
 castellatedMeshControls
 {
+
+    // Refinement parameters
+    // ~~~~~~~~~~~~~~~~~~~~~
+
+    // 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.
     maxLocalCells   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.
     maxGlobalCells  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)
     minRefinementCells 10;
+
+    // Number of buffer layers between different levels.
+    // 1 means normal 2:1 refinement restriction, larger means slower
+    // refinement.
     nCellsBetweenLevels 2;
-    features        ( );
+
+
+
+    // Explicit feature edge refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // Specifies a level for any cell intersected by its edges.
+    // This is a featureEdgeMesh, read from constant/triSurface for now.
+    features
+    (
+        //{
+        //    file "someLine.eMesh";
+        //    level 2;
+        //}
+    );
+
+
+
+    // Surface based refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // Specifies two levels for every surface. The first is the minimum level,
+    // every cell intersecting a surface gets refined up to the minimum level.
+    // The second level is the maximum level. Cells that 'see' multiple
+    // intersections where the intersections make an
+    // angle > resolveFeatureAngle get refined up to the maximum level.
+
     refinementSurfaces
     {
         motorBike
         {
+            // Surface-wise min and max refinement level
             level           ( 5 6 );
         }
     }
 
+    // Resolve sharp angles
     resolveFeatureAngle 30;
+
+
+    // Region-wise refinement
+    // ~~~~~~~~~~~~~~~~~~~~~~
+
+    // Specifies refinement level for cells in relation to a surface. One of
+    // three modes
+    // - distance. 'levels' specifies per distance to the surface the
+    //   wanted refinement level. The distances need to be specified in
+    //   descending order.
+    // - inside. 'levels' is only one entry and only the level is used. All
+    //   cells inside the surface get refined up to the level. The surface
+    //   needs to be closed for this to be possible.
+    // - outside. Same but cells outside.
     refinementRegions
     {
         refinementBox
@@ -62,19 +136,45 @@ castellatedMeshControls
         }
     }
 
+
+
+    // Mesh selection
+    // ~~~~~~~~~~~~~~
+
+    // After refinement patches get added for all refinementSurfaces and
+    // all cells intersecting the surfaces get put into these patches. The
+    // section reachable from the locationInMesh is kept.
+    // NOTE: This point should never be on a face, always inside a cell, even
+    // after refinement.
     locationInMesh  ( 3 3 0.43 );
 }
 
+
+// Settings for the snapping.
 snapControls
 {
+    //- 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.
     tolerance       4;
+
+    //- Number of mesh displacement relaxation iterations.
     nSolveIter      30;
+
+    //- Maximum number of snapping relaxation iterations. Should stop
+    //  before upon reaching a correct mesh.
     nRelaxIter      5;
 }
 
+
+// Settings for the layer addition.
 addLayersControls
 {
+    // Per final patch (so not geometry!) the layer information
     layers
     {
         minZ
@@ -418,41 +518,126 @@ addLayersControls
         }
     }
 
+    // Expansion factor for layer mesh
     expansionRatio  1;
+
+    //- Wanted thickness of final added cell layer. If multiple layers
+    //  is the
+    //  thickness of the layer furthest away from the wall.
+    //  Relative to undistorted size of cell outside layer.
     finalLayerRatio 0.3;
+
+    //- Minimum thickness of cell layer. If for any reason layer
+    //  cannot be above minThickness do not add layer.
+    //  Relative to undistorted size of cell outside layer.
     minThickness    0.1;
+
+    //- If points get not extruded do nGrow layers of connected faces that are
+    //  also not grown. This helps convergence of the layer addition process
+    //  close to features.
     nGrow           1;
+
+
+    // Advanced settings
+
+    //- When not to extrude surface. 0 is flat surface, 90 is when two faces
+    //  make straight angle.
     featureAngle    30;
+
+    //- Maximum number of snapping relaxation iterations. Should stop
+    //  before upon reaching a correct mesh.
     nRelaxIter      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;
 }
 
+
+// Generic mesh quality settings. At any undoable phase these determine
+// where to undo.
 meshQualityControls
 {
+    //- 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.02;
+
+    //- minimum normalised cell determinant
+    //- 1 = hex, <= 0 = folded or flattened illegal cell
     minDeterminant  0.001;
+
+    //- minFaceWeight (0 -> 0.5)
     minFaceWeight   0.02;
+
+    //- 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;
 }
 
+
+
+// Advanced
+
+// Flags for optional output
+// 0 : only write final meshes
+// 1 : write intermediate meshes
+// 2 : write volScalarField with cellLevel for postprocessing
+// 4 : write current intersections as .obj files
 debug           0;
 
+
+// Merge tolerance. Is fraction of overall bounding box of initial mesh.
+// Note: the write tolerance needs to be higher than this.
 mergeTolerance  1e-06;
 
 

tutorials/multiphase/cavitatingFoam/les/Allrun

@@ -3,7 +3,7 @@
 # Source tutorial run functions
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
-application="lesCavitatingFoam"
+application="cavitatingFoam"
 
 refineMeshByCellSet()
 {

tutorials/multiphase/cavitatingFoam/ras/throttle/Allrun

@@ -3,7 +3,7 @@
 # Source tutorial run functions
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
-application="rasCavitatingFoam"
+application="cavitatingFoam"
 
 refineMeshByCellSet()
 {

tutorials/multiphase/compressibleInterFoam/les/depthCharge2D/Allrun

@@ -3,7 +3,7 @@
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
 # Get application name from directory
-application="compressibleLesInterFoam"
+application="compressibleInterFoam"
 
 runApplication blockMesh
 cp 0/alpha1.org 0/alpha1

tutorials/multiphase/compressibleInterFoam/les/depthCharge3D/Allrun

@@ -3,7 +3,7 @@
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
 # Get application name from directory
-application="compressibleLesInterFoam"
+application="compressibleInterFoam"
 
 runApplication blockMesh
 cp 0/alpha1.org 0/alpha1