From 5d4030b7f2b5b2e165d193c498c73dd8d65417fb Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Wed, 21 Dec 2011 10:34:42 +0000
Subject: [PATCH 01/18] BUG: renumberMesh.C: reading bool from dictionary

---
 .../utilities/mesh/manipulation/renumberMesh/renumberMesh.C     | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/applications/utilities/mesh/manipulation/renumberMesh/renumberMesh.C b/applications/utilities/mesh/manipulation/renumberMesh/renumberMesh.C
index 7e090fa30e..924096a22f 100644
--- a/applications/utilities/mesh/manipulation/renumberMesh/renumberMesh.C
+++ b/applications/utilities/mesh/manipulation/renumberMesh/renumberMesh.C
@@ -562,7 +562,7 @@ int main(int argc, char *argv[])
                 << endl;
         }
 
-        writeMaps = readLabel(renumberDict.lookup("writeMaps"));
+        renumberDict.lookup("writeMaps") >> writeMaps;
         if (writeMaps)
         {
             Info<< "Writing renumber maps (new to old) to polyMesh." << nl

From 744162e74e1e74e890e66b2686f61f1648348ae2 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Wed, 21 Dec 2011 10:35:34 +0000
Subject: [PATCH 02/18] STYLE: tetIndices.H: improved comment

---
 .../meshes/polyMesh/polyMeshTetDecomposition/tetIndices.H   | 6 ++++--
 1 file changed, 4 insertions(+), 2 deletions(-)

diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshTetDecomposition/tetIndices.H b/src/OpenFOAM/meshes/polyMesh/polyMeshTetDecomposition/tetIndices.H
index afd31e642a..76f82fa78c 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshTetDecomposition/tetIndices.H
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshTetDecomposition/tetIndices.H
@@ -144,11 +144,13 @@ public:
             inline tetPointRef oldTet(const polyMesh& mesh) const;
 
             //- Return the geometry corresponding to the tri on the
-            //  mesh face for this tet from the supplied mesh
+            //  mesh face for this tet from the supplied mesh. Normal of
+            //  the tri points out of the cell.
             inline triPointRef faceTri(const polyMesh& mesh) const;
 
             //- Return the point indices corresponding to the tri on the mesh
-            //  face for this tet from the supplied mesh
+            //  face for this tet from the supplied mesh. Normal of
+            //  the tri points out of the cell.
             inline triFace faceTriIs(const polyMesh& mesh) const;
 
             //- Return the geometry corresponding to the tri on the

From 0504c7a4f6a79e4a4cc510fab8062eb5aa32facd Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Wed, 21 Dec 2011 11:22:02 +0000
Subject: [PATCH 03/18] STYLE: system/streamLines: documented new tracking
 option

---
 .../simpleFoam/motorBike/system/streamLines           | 11 +++++++++--
 1 file changed, 9 insertions(+), 2 deletions(-)

diff --git a/tutorials/incompressible/simpleFoam/motorBike/system/streamLines b/tutorials/incompressible/simpleFoam/motorBike/system/streamLines
index 8b82eccbb1..50b9eefdd1 100644
--- a/tutorials/incompressible/simpleFoam/motorBike/system/streamLines
+++ b/tutorials/incompressible/simpleFoam/motorBike/system/streamLines
@@ -28,8 +28,15 @@ streamLines
     // Steps particles can travel before being removed
     lifeTime        10000;
 
-    // Number of steps per cell (estimate). Set to 1 to disable subcycling.
-    nSubCycle 5;
+    //- Specify either absolute length of steps (trackLength) or a number
+    //  of subcycling steps per cell (nSubCycle)
+
+        // Size of single track segment [m]
+        //trackLength 1e-3;
+
+        // Number of steps per cell (estimate). Set to 1 to disable subcycling.
+        nSubCycle 5;
+
 
     // Cloud name to use
     cloudName       particleTracks;

From 704b33fa79397870c6ba22733e3cce7745b42122 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Wed, 21 Dec 2011 15:13:04 +0000
Subject: [PATCH 04/18] ENH: renumberMesh: get block ordering working

- added boundaryFirst renumbering method
- renumbering in parallel
- make api methods const
---
 .../manipulation/renumberMesh/Make/options    |   2 +-
 .../manipulation/renumberMesh/renumberMesh.C  | 125 +++++++----
 .../renumberMesh/renumberMeshDict             |  35 ++-
 .../CuthillMcKeeRenumber.C                    |  14 +-
 .../CuthillMcKeeRenumber.H                    |   6 +-
 src/renumberMethods/Make/files                |   1 +
 .../boundaryFirstRenumber.C                   | 201 ++++++++++++++++++
 .../boundaryFirstRenumber.H                   | 129 +++++++++++
 .../manualRenumber/manualRenumber.C           |   2 +-
 .../manualRenumber/manualRenumber.H           |   6 +-
 .../randomRenumber/randomRenumber.C           |   6 +-
 .../randomRenumber/randomRenumber.H           |   6 +-
 .../renumberMethod/renumberMethod.C           |   4 +-
 .../renumberMethod/renumberMethod.H           |   8 +-
 .../springRenumber/springRenumber.C           |   4 +-
 .../springRenumber/springRenumber.H           |   6 +-
 16 files changed, 473 insertions(+), 82 deletions(-)
 create mode 100644 src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
 create mode 100644 src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H

diff --git a/applications/utilities/mesh/manipulation/renumberMesh/Make/options b/applications/utilities/mesh/manipulation/renumberMesh/Make/options
index 2981cf30be..dc78db3436 100644
--- a/applications/utilities/mesh/manipulation/renumberMesh/Make/options
+++ b/applications/utilities/mesh/manipulation/renumberMesh/Make/options
@@ -11,4 +11,4 @@ EXE_LIBS = \
     -lfiniteVolume \
     -lgenericPatchFields \
     -lrenumberMethods \
-    -ldecompositionMethods
+    -ldecompositionMethods -L$(FOAM_LIBBIN)/dummy -lmetisDecomp -lscotchDecomp
diff --git a/applications/utilities/mesh/manipulation/renumberMesh/renumberMesh.C b/applications/utilities/mesh/manipulation/renumberMesh/renumberMesh.C
index 924096a22f..02e5a06567 100644
--- a/applications/utilities/mesh/manipulation/renumberMesh/renumberMesh.C
+++ b/applications/utilities/mesh/manipulation/renumberMesh/renumberMesh.C
@@ -29,7 +29,7 @@ Description
     renumbering all fields from all the time directories.
 
     By default uses bandCompression (CuthillMcKee) but will
-    read system/renumberMeshDict if present and use the method from there.
+    read system/renumberMeshDict if -dict option is present
 
 \*---------------------------------------------------------------------------*/
 
@@ -45,6 +45,7 @@ Description
 #include "renumberMethod.H"
 #include "zeroGradientFvPatchFields.H"
 #include "CuthillMcKeeRenumber.H"
+#include "fvMeshSubset.H"
 
 using namespace Foam;
 
@@ -455,6 +456,69 @@ autoPtr<mapPolyMesh> reorderMesh
 }
 
 
+// Return new to old cell numbering
+labelList regionRenumber
+(
+    const renumberMethod& method,
+    const fvMesh& mesh,
+    const labelList& cellToRegion
+)
+{
+    Info<< "Determining cell order:" << endl;
+
+    labelList cellOrder(cellToRegion.size());
+
+    label nRegions = max(cellToRegion)+1;
+
+    labelListList regionToCells(invertOneToMany(nRegions, cellToRegion));
+
+    label cellI = 0;
+
+    forAll(regionToCells, regionI)
+    {
+        Info<< "    region " << regionI << " starts at " << cellI << endl;
+
+        // Make sure no parallel comms
+        bool oldParRun = UPstream::parRun();
+        UPstream::parRun() = false;
+
+        // Per region do a reordering.
+        fvMeshSubset subsetter(mesh);
+        subsetter.setLargeCellSubset(cellToRegion, regionI);
+
+        const fvMesh& subMesh = subsetter.subMesh();
+
+        labelList subReverseCellOrder = method.renumber
+        (
+            subMesh,
+            subMesh.cellCentres()
+        );
+
+        labelList subCellOrder
+        (
+            invert
+            (
+                subMesh.nCells(),
+                subReverseCellOrder
+            )
+        );
+
+        // Restore state
+        UPstream::parRun() = oldParRun;
+
+        const labelList& cellMap = subsetter.cellMap();
+
+        forAll(subCellOrder, i)
+        {
+            cellOrder[cellI++] = cellMap[subCellOrder[i]];
+        }
+    }
+    Info<< endl;
+
+    return cellOrder;
+}
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 int main(int argc, char *argv[])
@@ -505,23 +569,28 @@ int main(int argc, char *argv[])
     label blockSize = 0;
 
     // Construct renumberMethod
+    autoPtr<IOdictionary> renumberDictPtr;
     autoPtr<renumberMethod> renumberPtr;
 
     if (readDict)
     {
         Info<< "Renumber according to renumberMeshDict." << nl << endl;
 
-        IOdictionary renumberDict
+        renumberDictPtr.reset
         (
-            IOobject
+            new IOdictionary
             (
-                "renumberMeshDict",
-                runTime.system(),
-                mesh,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
+                IOobject
+                (
+                    "renumberMeshDict",
+                    runTime.system(),
+                    mesh,
+                    IOobject::MUST_READ_IF_MODIFIED,
+                    IOobject::NO_WRITE
+                )
             )
         );
+        const IOdictionary renumberDict = renumberDictPtr();
 
         renumberPtr = renumberMethod::New(renumberDict);
 
@@ -683,22 +752,15 @@ int main(int argc, char *argv[])
         // fields is done correctly!
 
         label nBlocks = mesh.nCells() / blockSize;
-        Info<< "nBlocks = " << nBlocks << endl;
+        Info<< "nBlocks   = " << nBlocks << endl;
 
-        // Read decomposePar dictionary
-        IOdictionary decomposeDict
-        (
-            IOobject
-            (
-                "decomposeParDict",
-                runTime.system(),
-                mesh,
-                IOobject::MUST_READ_IF_MODIFIED,
-                IOobject::NO_WRITE
-            )
-        );
+        // Read decompositionMethod dictionary
+        dictionary decomposeDict(renumberDictPtr().subDict("blockCoeffs"));
         decomposeDict.set("numberOfSubdomains", nBlocks);
 
+        bool oldParRun = UPstream::parRun();
+        UPstream::parRun() = false;
+
         autoPtr<decompositionMethod> decomposePtr = decompositionMethod::New
         (
             decomposeDict
@@ -713,6 +775,9 @@ int main(int argc, char *argv[])
             )
         );
 
+        // Restore state
+        UPstream::parRun() = oldParRun;
+
         // For debugging: write out region
         createScalarField
         (
@@ -726,23 +791,7 @@ int main(int argc, char *argv[])
             << nl << endl;
 
 
-        // Find point per region
-        pointField regionPoints(nBlocks, vector::zero);
-        forAll(cellToRegion, cellI)
-        {
-            regionPoints[cellToRegion[cellI]] = mesh.cellCentres()[cellI];
-        }
-
-        // Use block based renumbering.
-        // Detemines old to new cell ordering
-        labelList reverseCellOrder = renumberPtr().renumber
-        (
-            mesh,
-            cellToRegion,
-            regionPoints
-        );
-
-        cellOrder = invert(mesh.nCells(), reverseCellOrder);
+        cellOrder = regionRenumber(renumberPtr(), mesh, cellToRegion);
 
         // Determine new to old face order with new cell numbering
         faceOrder = getRegionFaceOrder
diff --git a/applications/utilities/mesh/manipulation/renumberMesh/renumberMeshDict b/applications/utilities/mesh/manipulation/renumberMesh/renumberMeshDict
index a852fad7db..988379720b 100644
--- a/applications/utilities/mesh/manipulation/renumberMesh/renumberMeshDict
+++ b/applications/utilities/mesh/manipulation/renumberMesh/renumberMeshDict
@@ -22,10 +22,14 @@ writeMaps true;
 // e.g. nonBlockingGaussSeidel.
 sortCoupledFaceCells false;
 
-// Optional entry: renumber on a block-by-block basis. This can be used on
-// large cases to keep the blocks fitting in cache with all the the cache
-// missed bunched at the end.
-//blockSize 0;
+// Optional entry: renumber on a block-by-block basis. It uses a
+// blockCoeffs dictionary to construct a decompositionMethod to do
+// a block subdivision) and then applies the renumberMethod to each
+// block in turn. This can be used in large cases to keep the blocks
+// fitting in cache with all the the cache misses bunched at the end.
+// This number is the approximate size of the blocks - this gets converted
+// to a number of blocks that is the input to the decomposition method.
+//blockSize 1000;
 
 // Optional entry: sort points into internal and boundary points
 //orderPoints false;
@@ -37,11 +41,11 @@ method          CuthillMcKee;
 //method          random;
 //method          spring;
 
-CuthillMcKeeCoeffs
-{
-    // Reverse CuthillMcKee (RCM) or plain
-    reverse true;
-}
+//CuthillMcKeeCoeffs
+//{
+//    // Reverse CuthillMcKee (RCM) or plain
+//    reverse true;
+//}
 
 
 manualCoeffs
@@ -65,4 +69,17 @@ springCoeffs
 }
 
 
+blockCoeffs
+{
+    method          scotch;
+    //method hierarchical;
+    //hierarchicalCoeffs
+    //{
+    //    n           (1 2 1);
+    //    delta       0.001;
+    //    order       xyz;
+    //}
+}
+
+
 // ************************************************************************* //
diff --git a/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C b/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
index e1ec010335..2732484b36 100644
--- a/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
+++ b/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
@@ -63,7 +63,7 @@ Foam::labelList Foam::CuthillMcKeeRenumber::renumber
 (
     const polyMesh& mesh,
     const pointField& points
-)
+) const
 {
     CompactListList<label> cellCells;
     decompositionMethod::calcCellCells
@@ -79,10 +79,7 @@ Foam::labelList Foam::CuthillMcKeeRenumber::renumber
 
     if (reverse_)
     {
-        for (label i = 0; i < orderedToOld.size()/2; i++)
-        {
-            Swap(orderedToOld[i], orderedToOld[orderedToOld.size()-i-1]);
-        }
+        reverse(orderedToOld);
     }
 
     return invert(orderedToOld.size(), orderedToOld);
@@ -93,16 +90,13 @@ Foam::labelList Foam::CuthillMcKeeRenumber::renumber
 (
     const labelListList& cellCells,
     const pointField& points
-)
+) const
 {
     labelList orderedToOld = bandCompression(cellCells);
 
     if (reverse_)
     {
-        for (label i = 0; i < orderedToOld.size()/2; i++)
-        {
-            Swap(orderedToOld[i], orderedToOld[orderedToOld.size()-i-1]);
-        }
+        reverse(orderedToOld);
     }
 
     return invert(orderedToOld.size(), orderedToOld);
diff --git a/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H b/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
index 03148c8739..a4ec011b40 100644
--- a/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
+++ b/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
@@ -82,7 +82,7 @@ public:
 
         //- Return for every coordinate the wanted processor number.
         //  We need a polyMesh (to be able to load the file)
-        virtual labelList renumber(const pointField&)
+        virtual labelList renumber(const pointField&) const
         {
             notImplemented("CuthillMcKeeRenumber::renumber(const pointField&)");
             return labelList(0);
@@ -94,7 +94,7 @@ public:
         (
             const polyMesh& mesh,
             const pointField& cc
-        );
+        ) const;
 
         //- Return for every cell the new cell label.
         //  The connectivity is equal to mesh.cellCells() except
@@ -103,7 +103,7 @@ public:
         (
             const labelListList& cellCells,
             const pointField& cc
-        );
+        ) const;
 };
 
 
diff --git a/src/renumberMethods/Make/files b/src/renumberMethods/Make/files
index fb366f73c4..c19f1e6b6d 100644
--- a/src/renumberMethods/Make/files
+++ b/src/renumberMethods/Make/files
@@ -3,5 +3,6 @@ manualRenumber/manualRenumber.C
 CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
 randomRenumber/randomRenumber.C
 springRenumber/springRenumber.C
+boundaryFirstRenumber/boundaryFirstRenumber.C
 
 LIB = $(FOAM_LIBBIN)/librenumberMethods
diff --git a/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C b/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
new file mode 100644
index 0000000000..293e5ba19b
--- /dev/null
+++ b/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
@@ -0,0 +1,201 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "boundaryFirstRenumber.H"
+#include "addToRunTimeSelectionTable.H"
+#include "bandCompression.H"
+#include "decompositionMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(boundaryFirstRenumber, 0);
+
+    addToRunTimeSelectionTable
+    (
+        renumberMethod,
+        boundaryFirstRenumber,
+        dictionary
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::boundaryFirstRenumber::boundaryFirstRenumber
+(
+    const dictionary& renumberDict
+)
+:
+    renumberMethod(renumberDict),
+    reverse_
+    (
+        renumberDict.found(typeName + "Coeffs")
+      ? Switch(renumberDict.subDict(typeName + "Coeffs").lookup("reverse"))
+      : Switch(true)
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::labelList Foam::boundaryFirstRenumber::renumber
+(
+    const polyMesh& mesh,
+    const pointField& points
+) const
+{
+    const labelList& own = mesh.faceOwner();
+    const labelList& nei = mesh.faceNeighbour();
+
+    // Distance to boundary. Minimum of neighbours.
+    labelList distance(mesh.nCells(), -1);
+
+    // New order
+    labelList newOrder(mesh.nCells());
+    label cellInOrder = 0;
+
+
+    // Starting faces for walk. These are the zero distance faces.
+    DynamicList<label> frontFaces(mesh.nFaces()-mesh.nInternalFaces());
+
+    const polyBoundaryMesh& patches = mesh.boundaryMesh();
+
+    forAll(patches, patchI)
+    {
+        const polyPatch& pp = patches[patchI];
+
+        //- Note: cannot check for empty since these are introduced by
+        //  fvMeshSubset. Also most loops don't care about patch type.
+        //if (!isA<emptyPolyPatch>(pp))
+        {
+            forAll(pp, i)
+            {
+                frontFaces.append(pp.start()+i);
+            }
+        }
+    }
+
+    // Loop over all frontFaces
+    label currentDistance = 0;
+    while (frontFaces.size() > 0)
+    {
+        DynamicList<label> frontCells(frontFaces.size());
+
+        // Set all of frontFaces' neighbours to current distance
+
+        forAll(frontFaces, i)
+        {
+            label faceI = frontFaces[i];
+
+            if (mesh.isInternalFace(faceI))
+            {
+                label ownCellI = own[faceI];
+                if (distance[ownCellI] == -1)
+                {
+                    distance[ownCellI] = currentDistance;
+                    frontCells.append(ownCellI);
+                }
+                label neiCellI = nei[faceI];
+                if (distance[neiCellI] == -1)
+                {
+                    distance[neiCellI] = currentDistance;
+                    frontCells.append(neiCellI);
+                }
+            }
+            else
+            {
+                label ownCellI = own[faceI];
+                if (distance[ownCellI] == -1)
+                {
+                    distance[ownCellI] = currentDistance;
+                    frontCells.append(ownCellI);
+                }
+            }
+        }
+
+        //Pout<< "For distance:" << currentDistance
+        //    << " from " << frontFaces.size() << " faces to "
+        //    << frontCells.size() << " cells." << endl;
+
+
+        // TBD. Determine order within current shell (frontCells). For now
+        // just add them.
+        forAll(frontCells, i)
+        {
+            newOrder[cellInOrder] = frontCells[i];
+            cellInOrder++;
+        }
+
+        // From cells to faces
+        frontFaces.clear();
+        forAll(frontCells, i)
+        {
+            label cellI = frontCells[i];
+            const cell& cFaces = mesh.cells()[cellI];
+
+            forAll(cFaces, i)
+            {
+                label faceI = cFaces[i];
+                if (mesh.isInternalFace(faceI))
+                {
+                    label nbrCellI =
+                    (
+                        mesh.faceOwner()[faceI] == cellI
+                      ? mesh.faceNeighbour()[faceI]
+                      : mesh.faceOwner()[faceI]
+                    );
+                    if (distance[nbrCellI] == -1)
+                    {
+                        frontFaces.append(faceI);
+                    }
+                }
+            }
+        }
+
+        currentDistance++;
+    }
+
+    // Return furthest away cell first
+    if (reverse_)
+    {
+        reverse(newOrder);
+    }
+
+    //forAll(newOrder, i)
+    //{
+    //    label cellI = newOrder[i];
+    //
+    //    Pout<< "cell:" << cellI << endl;
+    //    Pout<< " at distance:" << distance[cellI]
+    //        << endl;
+    //}
+
+    return invert(newOrder.size(), newOrder);
+}
+
+
+// ************************************************************************* //
diff --git a/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H b/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H
new file mode 100644
index 0000000000..690af21455
--- /dev/null
+++ b/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H
@@ -0,0 +1,129 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::boundaryFirstRenumber
+
+Description
+    Cuthill-McKee renumbering
+
+SourceFiles
+    boundaryFirstRenumber.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef boundaryFirstRenumber_H
+#define boundaryFirstRenumber_H
+
+#include "renumberMethod.H"
+#include "Switch.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                   Class boundaryFirstRenumber Declaration
+\*---------------------------------------------------------------------------*/
+
+class boundaryFirstRenumber
+:
+    public renumberMethod
+{
+    // Private data
+
+        const Switch reverse_;
+
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct and assignment
+        void operator=(const boundaryFirstRenumber&);
+        boundaryFirstRenumber(const boundaryFirstRenumber&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("boundaryFirst");
+
+
+    // Constructors
+
+        //- Construct given the renumber dictionary
+        boundaryFirstRenumber(const dictionary& renumberDict);
+
+
+    //- Destructor
+    virtual ~boundaryFirstRenumber()
+    {}
+
+
+    // Member Functions
+
+        //- Return for every coordinate the wanted processor number.
+        //  We need a polyMesh (to be able to load the file)
+        virtual labelList renumber(const pointField&) const
+        {
+            notImplemented
+            (
+                "boundaryFirstRenumber::renumber(const pointField&)"
+            );
+            return labelList(0);
+        }
+
+        //- Return for every coordinate the wanted processor number. Use the
+        //  mesh connectivity (if needed)
+        virtual labelList renumber
+        (
+            const polyMesh& mesh,
+            const pointField& cc
+        ) const;
+
+        //- Return for every cell the new cell label.
+        //  The connectivity is equal to mesh.cellCells() except
+        //  - the connections are across coupled patches
+        virtual labelList renumber
+        (
+            const labelListList& cellCells,
+            const pointField& cc
+        ) const
+        {
+            notImplemented
+            (
+                "boundaryFirstRenumber::renumber"
+                "(const labelListList&, const pointField&)"
+            );
+            return labelList(0);
+        }
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/renumberMethods/manualRenumber/manualRenumber.C b/src/renumberMethods/manualRenumber/manualRenumber.C
index 40989264b5..92b5c8e00c 100644
--- a/src/renumberMethods/manualRenumber/manualRenumber.C
+++ b/src/renumberMethods/manualRenumber/manualRenumber.C
@@ -61,7 +61,7 @@ Foam::labelList Foam::manualRenumber::renumber
 (
     const polyMesh& mesh,
     const pointField& points
-)
+) const
 {
     labelIOList oldToNew
     (
diff --git a/src/renumberMethods/manualRenumber/manualRenumber.H b/src/renumberMethods/manualRenumber/manualRenumber.H
index 635e93c614..a4d4b2f46f 100644
--- a/src/renumberMethods/manualRenumber/manualRenumber.H
+++ b/src/renumberMethods/manualRenumber/manualRenumber.H
@@ -81,7 +81,7 @@ public:
 
         //- Return for every coordinate the wanted processor number.
         //  We need a polyMesh (to be able to load the file)
-        virtual labelList renumber(const pointField&)
+        virtual labelList renumber(const pointField&) const
         {
             notImplemented("manualRenumber::renumber(const pointField&)");
             return labelList(0);
@@ -93,7 +93,7 @@ public:
         (
             const polyMesh& mesh,
             const pointField& cc
-        );
+        ) const;
 
         //- Return for every cell the new cell label.
         //  The connectivity is equal to mesh.cellCells() except
@@ -102,7 +102,7 @@ public:
         (
             const labelListList& cellCells,
             const pointField& cc
-        )
+        ) const
         {
             notImplemented
             (
diff --git a/src/renumberMethods/randomRenumber/randomRenumber.C b/src/renumberMethods/randomRenumber/randomRenumber.C
index 27ddc29c3f..2dae28e438 100644
--- a/src/renumberMethods/randomRenumber/randomRenumber.C
+++ b/src/renumberMethods/randomRenumber/randomRenumber.C
@@ -55,7 +55,7 @@ Foam::randomRenumber::randomRenumber(const dictionary& renumberDict)
 Foam::labelList Foam::randomRenumber::renumber
 (
     const pointField& points
-)
+) const
 {
     Random rndGen(0);
 
@@ -77,7 +77,7 @@ Foam::labelList Foam::randomRenumber::renumber
 (
     const polyMesh& mesh,
     const pointField& points
-)
+) const
 {
     return renumber(points);
 }
@@ -87,7 +87,7 @@ Foam::labelList Foam::randomRenumber::renumber
 (
     const labelListList& cellCells,
     const pointField& points
-)
+) const
 {
     return renumber(points);
 }
diff --git a/src/renumberMethods/randomRenumber/randomRenumber.H b/src/renumberMethods/randomRenumber/randomRenumber.H
index be704d1097..0665d5159b 100644
--- a/src/renumberMethods/randomRenumber/randomRenumber.H
+++ b/src/renumberMethods/randomRenumber/randomRenumber.H
@@ -76,7 +76,7 @@ public:
 
         //- Return for every coordinate the wanted processor number.
         //  We need a polyMesh (to be able to load the file)
-        virtual labelList renumber(const pointField&);
+        virtual labelList renumber(const pointField&) const;
 
         //- Return for every coordinate the wanted processor number. Use the
         //  mesh connectivity (if needed)
@@ -84,7 +84,7 @@ public:
         (
             const polyMesh& mesh,
             const pointField& cc
-        );
+        ) const;
 
         //- Return for every cell the new cell label.
         //  The connectivity is equal to mesh.cellCells() except
@@ -93,7 +93,7 @@ public:
         (
             const labelListList& cellCells,
             const pointField& cc
-        );
+        ) const;
 };
 
 
diff --git a/src/renumberMethods/renumberMethod/renumberMethod.C b/src/renumberMethods/renumberMethod/renumberMethod.C
index 53138673b8..76ba1911e2 100644
--- a/src/renumberMethods/renumberMethod/renumberMethod.C
+++ b/src/renumberMethods/renumberMethod/renumberMethod.C
@@ -72,7 +72,7 @@ Foam::labelList Foam::renumberMethod::renumber
 (
     const polyMesh& mesh,
     const pointField& points
-)
+) const
 {
     CompactListList<label> cellCells;
     decompositionMethod::calcCellCells
@@ -94,7 +94,7 @@ Foam::labelList Foam::renumberMethod::renumber
     const polyMesh& mesh,
     const labelList& fineToCoarse,
     const pointField& coarsePoints
-)
+) const
 {
     CompactListList<label> coarseCellCells;
     decompositionMethod::calcCellCells
diff --git a/src/renumberMethods/renumberMethod/renumberMethod.H b/src/renumberMethods/renumberMethod/renumberMethod.H
index b4b6e89f22..ead92b1e41 100644
--- a/src/renumberMethods/renumberMethod/renumberMethod.H
+++ b/src/renumberMethods/renumberMethod/renumberMethod.H
@@ -111,7 +111,7 @@ public:
 
         //- Return for every cell the new cell label.
         //  This is only defined for geometric renumberMethods.
-        virtual labelList renumber(const pointField&)
+        virtual labelList renumber(const pointField&) const
         {
             notImplemented
             (
@@ -122,7 +122,7 @@ public:
 
         //- Return for every cell the new cell label. Use the
         //  mesh connectivity (if needed)
-        virtual labelList renumber(const polyMesh&, const pointField&);
+        virtual labelList renumber(const polyMesh&, const pointField&) const;
 
         //- Return for every cell the new cell label. Gets
         //  passed agglomeration map (from fine to coarse cells) and coarse
@@ -136,7 +136,7 @@ public:
             const polyMesh& mesh,
             const labelList& cellToRegion,
             const pointField& regionPoints
-        );
+        ) const;
 
         //- Return for every cell the new cell label.
         //  The connectivity is equal to mesh.cellCells() except
@@ -145,7 +145,7 @@ public:
         (
             const labelListList& cellCells,
             const pointField& cc
-        ) = 0;
+        ) const = 0;
 
 };
 
diff --git a/src/renumberMethods/springRenumber/springRenumber.C b/src/renumberMethods/springRenumber/springRenumber.C
index 5c6fa49350..64a15145ce 100644
--- a/src/renumberMethods/springRenumber/springRenumber.C
+++ b/src/renumberMethods/springRenumber/springRenumber.C
@@ -60,7 +60,7 @@ Foam::labelList Foam::springRenumber::renumber
 (
     const polyMesh& mesh,
     const pointField& points
-)
+) const
 {
     CompactListList<label> cellCells;
     decompositionMethod::calcCellCells
@@ -80,7 +80,7 @@ Foam::labelList Foam::springRenumber::renumber
 (
     const labelListList& cellCells,
     const pointField& points
-)
+) const
 {
     // Look at cell index as a 1D position parameter.
     // Move cells to the average 'position' of their neighbour.
diff --git a/src/renumberMethods/springRenumber/springRenumber.H b/src/renumberMethods/springRenumber/springRenumber.H
index 64342af1ec..30f7a3d838 100644
--- a/src/renumberMethods/springRenumber/springRenumber.H
+++ b/src/renumberMethods/springRenumber/springRenumber.H
@@ -96,7 +96,7 @@ public:
     // Member Functions
 
         //- Return for every coordinate the wanted processor number.
-        virtual labelList renumber(const pointField&)
+        virtual labelList renumber(const pointField&) const
         {
             notImplemented("springRenumber::renumber(const pointField&)");
             return labelList(0);
@@ -108,7 +108,7 @@ public:
         (
             const polyMesh& mesh,
             const pointField& cc
-        );
+        ) const;
 
         //- Return for every cell the new cell label.
         //  The connectivity is equal to mesh.cellCells() except
@@ -117,7 +117,7 @@ public:
         (
             const labelListList& cellCells,
             const pointField& cc
-        );
+        ) const;
 };
 
 

From 7e95308c6b886ebc37e3d7a15933644d320bb929 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Wed, 21 Dec 2011 16:16:54 +0000
Subject: [PATCH 05/18] ENH: decomposePar.C: add shortcircuit to avoid
 allocating point mappers

---
 .../parallelProcessing/decomposePar/decomposePar.C     | 10 +++++++++-
 1 file changed, 9 insertions(+), 1 deletion(-)

diff --git a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
index fd58570b34..02acd59236 100644
--- a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
+++ b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
@@ -528,7 +528,7 @@ int main(int argc, char *argv[])
                         << "Cell number should be between 0 and "
                         << mesh.nCells()-1 << nl
                         << "On this mesh the particle should be in cell "
-                        <<  mesh.findCell(iter().position())
+                        << mesh.findCell(iter().position())
                         << exit(FatalError);
                 }
 
@@ -789,6 +789,14 @@ int main(int argc, char *argv[])
 
 
         // Point fields
+        if
+        (
+            pointScalarFields.size()
+         || pointVectorFields.size()
+         || pointSphericalTensorFields.size()
+         || pointSymmTensorFields.size()
+         || pointTensorFields.size()
+        )
         {
             labelIOList pointProcAddressing
             (

From 0f2e8be91c6bea6de15e1bb7f828e7ec4eac683b Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Thu, 22 Dec 2011 11:39:54 +0000
Subject: [PATCH 06/18] ENH: mapFields: add -subtract option to subtract whilst
 mapping

---
 .../mapFields/MapConsistentVolFields.H        |  20 +-
 .../preProcessing/mapFields/MapMeshes.H       | 263 ++++++++++++++++++
 .../preProcessing/mapFields/MapVolFields.H    |  13 +-
 .../preProcessing/mapFields/mapFields.C       | 217 ++++++---------
 .../meshToMesh/meshToMesh.H                   |  45 +--
 .../meshToMesh/meshToMeshInterpolate.C        | 121 ++++----
 6 files changed, 476 insertions(+), 203 deletions(-)
 create mode 100644 applications/utilities/preProcessing/mapFields/MapMeshes.H

diff --git a/applications/utilities/preProcessing/mapFields/MapConsistentVolFields.H b/applications/utilities/preProcessing/mapFields/MapConsistentVolFields.H
index 2a203eac48..9015c473a2 100644
--- a/applications/utilities/preProcessing/mapFields/MapConsistentVolFields.H
+++ b/applications/utilities/preProcessing/mapFields/MapConsistentVolFields.H
@@ -35,12 +35,13 @@ License
 namespace Foam
 {
 
-template<class Type>
+template<class Type, class CombineOp>
 void MapConsistentVolFields
 (
     const IOobjectList& objects,
     const meshToMesh& meshToMeshInterp,
-    const meshToMesh::order& mapOrder
+    const meshToMesh::order& mapOrder,
+    const CombineOp& cop
 )
 {
     const fvMesh& meshSource = meshToMeshInterp.fromMesh();
@@ -84,7 +85,13 @@ void MapConsistentVolFields
             );
 
             // Interpolate field
-            meshToMeshInterp.interpolate(fieldTarget, fieldSource, mapOrder);
+            meshToMeshInterp.interpolate//<Type, eqOp<Type> >
+            (
+                fieldTarget,
+                fieldSource,
+                mapOrder,
+                cop
+            );
 
             // Write field
             fieldTarget.write();
@@ -97,7 +104,12 @@ void MapConsistentVolFields
             GeometricField<Type, fvPatchField, volMesh> fieldTarget
             (
                 fieldTargetIOobject,
-                meshToMeshInterp.interpolate(fieldSource, mapOrder)
+                meshToMeshInterp.interpolate//<Type, eqOp<Type> >
+                (
+                    fieldSource,
+                    mapOrder,
+                    cop
+                )
             );
 
             // Write field
diff --git a/applications/utilities/preProcessing/mapFields/MapMeshes.H b/applications/utilities/preProcessing/mapFields/MapMeshes.H
new file mode 100644
index 0000000000..e497eb0ccd
--- /dev/null
+++ b/applications/utilities/preProcessing/mapFields/MapMeshes.H
@@ -0,0 +1,263 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef MapMeshes_H
+#define MapMeshes_H
+
+#include "MapVolFields.H"
+#include "MapConsistentVolFields.H"
+#include "mapLagrangian.H"
+#include "UnMapped.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+template<template<class> class CombineOp>
+void MapConsistentMesh
+(
+    const fvMesh& meshSource,
+    const fvMesh& meshTarget,
+    const meshToMesh::order& mapOrder
+)
+{
+    // Create the interpolation scheme
+    meshToMesh meshToMeshInterp(meshSource, meshTarget);
+
+    Info<< nl
+        << "Consistently creating and mapping fields for time "
+        << meshSource.time().timeName() << nl << endl;
+
+    {
+        // Search for list of objects for this time
+        IOobjectList objects(meshSource, meshSource.time().timeName());
+
+        // Map volFields
+        // ~~~~~~~~~~~~~
+        MapConsistentVolFields<scalar>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<scalar>()
+        );
+        MapConsistentVolFields<vector>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<vector>()
+        );
+        MapConsistentVolFields<sphericalTensor>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<sphericalTensor>()
+        );
+        MapConsistentVolFields<symmTensor>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<symmTensor>()
+        );
+        MapConsistentVolFields<tensor>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<tensor>()
+        );
+    }
+
+    {
+        // Search for list of target objects for this time
+        IOobjectList objects(meshTarget, meshTarget.time().timeName());
+
+        // Mark surfaceFields as unmapped
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        UnMapped<surfaceScalarField>(objects);
+        UnMapped<surfaceVectorField>(objects);
+        UnMapped<surfaceSphericalTensorField>(objects);
+        UnMapped<surfaceSymmTensorField>(objects);
+        UnMapped<surfaceTensorField>(objects);
+
+        // Mark pointFields as unmapped
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        UnMapped<pointScalarField>(objects);
+        UnMapped<pointVectorField>(objects);
+        UnMapped<pointSphericalTensorField>(objects);
+        UnMapped<pointSymmTensorField>(objects);
+        UnMapped<pointTensorField>(objects);
+    }
+
+    mapLagrangian(meshToMeshInterp);
+}
+
+
+template<template<class> class CombineOp>
+void MapSubMesh
+(
+    const fvMesh& meshSource,
+    const fvMesh& meshTarget,
+    const HashTable<word>& patchMap,
+    const wordList& cuttingPatches,
+    const meshToMesh::order& mapOrder
+)
+{
+    // Create the interpolation scheme
+    meshToMesh meshToMeshInterp
+    (
+        meshSource,
+        meshTarget,
+        patchMap,
+        cuttingPatches
+    );
+
+    Info<< nl
+        << "Mapping fields for time " << meshSource.time().timeName()
+        << nl << endl;
+
+    {
+        // Search for list of source objects for this time
+        IOobjectList objects(meshSource, meshSource.time().timeName());
+
+        // Map volFields
+        // ~~~~~~~~~~~~~
+        MapVolFields<scalar>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<scalar>()
+        );
+        MapVolFields<vector>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<vector>()
+        );
+        MapVolFields<sphericalTensor>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<sphericalTensor>()
+        );
+        MapVolFields<symmTensor>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<symmTensor>()
+        );
+        MapVolFields<tensor>
+        (
+            objects,
+            meshToMeshInterp,
+            mapOrder,
+            CombineOp<tensor>()
+        );
+    }
+
+    {
+        // Search for list of target objects for this time
+        IOobjectList objects(meshTarget, meshTarget.time().timeName());
+
+        // Mark surfaceFields as unmapped
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        UnMapped<surfaceScalarField>(objects);
+        UnMapped<surfaceVectorField>(objects);
+        UnMapped<surfaceSphericalTensorField>(objects);
+        UnMapped<surfaceSymmTensorField>(objects);
+        UnMapped<surfaceTensorField>(objects);
+
+        // Mark pointFields as unmapped
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        UnMapped<pointScalarField>(objects);
+        UnMapped<pointVectorField>(objects);
+        UnMapped<pointSphericalTensorField>(objects);
+        UnMapped<pointSymmTensorField>(objects);
+        UnMapped<pointTensorField>(objects);
+    }
+
+    mapLagrangian(meshToMeshInterp);
+}
+
+
+template<template<class> class CombineOp>
+void MapConsistentSubMesh
+(
+    const fvMesh& meshSource,
+    const fvMesh& meshTarget,
+    const meshToMesh::order& mapOrder
+)
+{
+    HashTable<word> patchMap;
+    HashTable<label> cuttingPatchTable;
+
+    forAll(meshTarget.boundary(), patchi)
+    {
+        if (!isA<processorFvPatch>(meshTarget.boundary()[patchi]))
+        {
+            patchMap.insert
+            (
+                meshTarget.boundary()[patchi].name(),
+                meshTarget.boundary()[patchi].name()
+            );
+        }
+        else
+        {
+            cuttingPatchTable.insert
+            (
+                meshTarget.boundaryMesh()[patchi].name(),
+                -1
+            );
+        }
+    }
+
+    MapSubMesh<CombineOp>
+    (
+        meshSource,
+        meshTarget,
+        patchMap,
+        cuttingPatchTable.toc(),
+        mapOrder
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/applications/utilities/preProcessing/mapFields/MapVolFields.H b/applications/utilities/preProcessing/mapFields/MapVolFields.H
index 64547a54de..65eb44d12d 100644
--- a/applications/utilities/preProcessing/mapFields/MapVolFields.H
+++ b/applications/utilities/preProcessing/mapFields/MapVolFields.H
@@ -35,12 +35,13 @@ License
 namespace Foam
 {
 
-template<class Type>
+template<class Type, class CombineOp>
 void MapVolFields
 (
     const IOobjectList& objects,
     const meshToMesh& meshToMeshInterp,
-    const meshToMesh::order& mapOrder
+    const meshToMesh::order& mapOrder,
+    const CombineOp& cop
 )
 {
     const fvMesh& meshSource = meshToMeshInterp.fromMesh();
@@ -84,7 +85,13 @@ void MapVolFields
             );
 
             // Interpolate field
-            meshToMeshInterp.interpolate(fieldTarget, fieldSource, mapOrder);
+            meshToMeshInterp.interpolate//<Type, eqOp<Type> >
+            (
+                fieldTarget,
+                fieldSource,
+                mapOrder,
+                cop
+            );
 
             // Write field
             fieldTarget.write();
diff --git a/applications/utilities/preProcessing/mapFields/mapFields.C b/applications/utilities/preProcessing/mapFields/mapFields.C
index 8d70911cab..9571be6fb0 100644
--- a/applications/utilities/preProcessing/mapFields/mapFields.C
+++ b/applications/utilities/preProcessing/mapFields/mapFields.C
@@ -34,11 +34,8 @@ Description
 
 #include "fvCFD.H"
 #include "meshToMesh.H"
-#include "MapVolFields.H"
-#include "MapConsistentVolFields.H"
-#include "UnMapped.H"
 #include "processorFvPatch.H"
-#include "mapLagrangian.H"
+#include "MapMeshes.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -46,56 +43,28 @@ void mapConsistentMesh
 (
     const fvMesh& meshSource,
     const fvMesh& meshTarget,
-    const meshToMesh::order& mapOrder
+    const meshToMesh::order& mapOrder,
+    const bool subtract
 )
 {
-    // Create the interpolation scheme
-    meshToMesh meshToMeshInterp(meshSource, meshTarget);
-
-    Info<< nl
-        << "Consistently creating and mapping fields for time "
-        << meshSource.time().timeName() << nl << endl;
-
+    if (subtract)
     {
-        // Search for list of objects for this time
-        IOobjectList objects(meshSource, meshSource.time().timeName());
-
-        // Map volFields
-        // ~~~~~~~~~~~~~
-        MapConsistentVolFields<scalar>(objects, meshToMeshInterp, mapOrder);
-        MapConsistentVolFields<vector>(objects, meshToMeshInterp, mapOrder);
-        MapConsistentVolFields<sphericalTensor>
+        MapConsistentMesh<minusEqOp>
         (
-            objects,
-            meshToMeshInterp,
+            meshSource,
+            meshTarget,
             mapOrder
         );
-        MapConsistentVolFields<symmTensor>(objects, meshToMeshInterp, mapOrder);
-        MapConsistentVolFields<tensor>(objects, meshToMeshInterp, mapOrder);
     }
-
+    else
     {
-        // Search for list of target objects for this time
-        IOobjectList objects(meshTarget, meshTarget.time().timeName());
-
-        // Mark surfaceFields as unmapped
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        UnMapped<surfaceScalarField>(objects);
-        UnMapped<surfaceVectorField>(objects);
-        UnMapped<surfaceSphericalTensorField>(objects);
-        UnMapped<surfaceSymmTensorField>(objects);
-        UnMapped<surfaceTensorField>(objects);
-
-        // Mark pointFields as unmapped
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        UnMapped<pointScalarField>(objects);
-        UnMapped<pointVectorField>(objects);
-        UnMapped<pointSphericalTensorField>(objects);
-        UnMapped<pointSymmTensorField>(objects);
-        UnMapped<pointTensorField>(objects);
+        MapConsistentMesh<eqOp>
+        (
+            meshSource,
+            meshTarget,
+            mapOrder
+        );
     }
-
-    mapLagrangian(meshToMeshInterp);
 }
 
 
@@ -105,57 +74,32 @@ void mapSubMesh
     const fvMesh& meshTarget,
     const HashTable<word>& patchMap,
     const wordList& cuttingPatches,
-    const meshToMesh::order& mapOrder
+    const meshToMesh::order& mapOrder,
+    const bool subtract
 )
 {
-    // Create the interpolation scheme
-    meshToMesh meshToMeshInterp
-    (
-        meshSource,
-        meshTarget,
-        patchMap,
-        cuttingPatches
-    );
-
-    Info<< nl
-        << "Mapping fields for time " << meshSource.time().timeName()
-        << nl << endl;
-
+    if (subtract)
     {
-        // Search for list of source objects for this time
-        IOobjectList objects(meshSource, meshSource.time().timeName());
-
-        // Map volFields
-        // ~~~~~~~~~~~~~
-        MapVolFields<scalar>(objects, meshToMeshInterp, mapOrder);
-        MapVolFields<vector>(objects, meshToMeshInterp, mapOrder);
-        MapVolFields<sphericalTensor>(objects, meshToMeshInterp, mapOrder);
-        MapVolFields<symmTensor>(objects, meshToMeshInterp, mapOrder);
-        MapVolFields<tensor>(objects, meshToMeshInterp, mapOrder);
+        MapSubMesh<minusEqOp>
+        (
+            meshSource,
+            meshTarget,
+            patchMap,
+            cuttingPatches,
+            mapOrder
+        );
     }
-
+    else
     {
-        // Search for list of target objects for this time
-        IOobjectList objects(meshTarget, meshTarget.time().timeName());
-
-        // Mark surfaceFields as unmapped
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        UnMapped<surfaceScalarField>(objects);
-        UnMapped<surfaceVectorField>(objects);
-        UnMapped<surfaceSphericalTensorField>(objects);
-        UnMapped<surfaceSymmTensorField>(objects);
-        UnMapped<surfaceTensorField>(objects);
-
-        // Mark pointFields as unmapped
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        UnMapped<pointScalarField>(objects);
-        UnMapped<pointVectorField>(objects);
-        UnMapped<pointSphericalTensorField>(objects);
-        UnMapped<pointSymmTensorField>(objects);
-        UnMapped<pointTensorField>(objects);
+        MapSubMesh<eqOp>
+        (
+            meshSource,
+            meshTarget,
+            patchMap,
+            cuttingPatches,
+            mapOrder
+        );
     }
-
-    mapLagrangian(meshToMeshInterp);
 }
 
 
@@ -163,40 +107,28 @@ void mapConsistentSubMesh
 (
     const fvMesh& meshSource,
     const fvMesh& meshTarget,
-    const meshToMesh::order& mapOrder
+    const meshToMesh::order& mapOrder,
+    const bool subtract
 )
 {
-    HashTable<word> patchMap;
-    HashTable<label> cuttingPatchTable;
-
-    forAll(meshTarget.boundary(), patchi)
+    if (subtract)
     {
-        if (!isA<processorFvPatch>(meshTarget.boundary()[patchi]))
-        {
-            patchMap.insert
-            (
-                meshTarget.boundary()[patchi].name(),
-                meshTarget.boundary()[patchi].name()
-            );
-        }
-        else
-        {
-            cuttingPatchTable.insert
-            (
-                meshTarget.boundaryMesh()[patchi].name(),
-                -1
-            );
-        }
+        MapConsistentSubMesh<minusEqOp>
+        (
+            meshSource,
+            meshTarget,
+            mapOrder
+        );
+    }
+    else
+    {
+        MapConsistentSubMesh<eqOp>
+        (
+            meshSource,
+            meshTarget,
+            mapOrder
+        );
     }
-
-    mapSubMesh
-    (
-        meshSource,
-        meshTarget,
-        patchMap,
-        cuttingPatchTable.toc(),
-        mapOrder
-    );
 }
 
 
@@ -288,6 +220,11 @@ int main(int argc, char *argv[])
         "word",
         "specify the mapping method"
     );
+    argList::addBoolOption
+    (
+        "subtract",
+        "subtract mapped source from target"
+    );
 
     argList args(argc, argv);
 
@@ -350,6 +287,13 @@ int main(int argc, char *argv[])
         Info<< "Mapping method: " << mapMethod << endl;
     }
 
+    const bool subtract = args.optionFound("subtract");
+    if (subtract)
+    {
+        Info<< "Subtracting mapped source field from target" << endl;
+    }
+
+
     #include "createTimes.H"
 
     HashTable<word> patchMap;
@@ -431,7 +375,13 @@ int main(int argc, char *argv[])
 
             if (consistent)
             {
-                mapConsistentSubMesh(meshSource, meshTarget, mapOrder);
+                mapConsistentSubMesh
+                (
+                    meshSource,
+                    meshTarget,
+                    mapOrder,
+                    subtract
+                );
             }
             else
             {
@@ -441,7 +391,8 @@ int main(int argc, char *argv[])
                     meshTarget,
                     patchMap,
                     cuttingPatches,
-                    mapOrder
+                    mapOrder,
+                    subtract
                 );
             }
         }
@@ -503,7 +454,13 @@ int main(int argc, char *argv[])
 
             if (consistent)
             {
-                mapConsistentSubMesh(meshSource, meshTarget, mapOrder);
+                mapConsistentSubMesh
+                (
+                    meshSource,
+                    meshTarget,
+                    mapOrder,
+                    subtract
+                );
             }
             else
             {
@@ -513,7 +470,8 @@ int main(int argc, char *argv[])
                     meshTarget,
                     patchMap,
                     addProcessorPatches(meshTarget, cuttingPatches),
-                    mapOrder
+                    mapOrder,
+                    subtract
                 );
             }
         }
@@ -629,7 +587,8 @@ int main(int argc, char *argv[])
                             (
                                 meshSource,
                                 meshTarget,
-                                mapOrder
+                                mapOrder,
+                                subtract
                             );
                         }
                         else
@@ -640,7 +599,8 @@ int main(int argc, char *argv[])
                                 meshTarget,
                                 patchMap,
                                 addProcessorPatches(meshTarget, cuttingPatches),
-                                mapOrder
+                                mapOrder,
+                                subtract
                             );
                         }
                     }
@@ -679,7 +639,7 @@ int main(int argc, char *argv[])
 
         if (consistent)
         {
-            mapConsistentMesh(meshSource, meshTarget, mapOrder);
+            mapConsistentMesh(meshSource, meshTarget, mapOrder, subtract);
         }
         else
         {
@@ -689,7 +649,8 @@ int main(int argc, char *argv[])
                 meshTarget,
                 patchMap,
                 cuttingPatches,
-                mapOrder
+                mapOrder,
+                subtract
             );
         }
     }
diff --git a/src/sampling/meshToMeshInterpolation/meshToMesh/meshToMesh.H b/src/sampling/meshToMeshInterpolation/meshToMesh/meshToMesh.H
index 785829cc60..83d1ee902a 100644
--- a/src/sampling/meshToMeshInterpolation/meshToMesh/meshToMesh.H
+++ b/src/sampling/meshToMeshInterpolation/meshToMesh/meshToMesh.H
@@ -219,84 +219,93 @@ public:
         // Interpolation
 
             //- Map field
-            template<class Type>
+            template<class Type, class CombineOp>
             void mapField
             (
                 Field<Type>&,
                 const Field<Type>&,
-                const labelList& adr
+                const labelList& adr,
+                const CombineOp& cop
             ) const;
 
             //- Interpolate field using inverse-distance weights
-            template<class Type>
+            template<class Type, class CombineOp>
             void interpolateField
             (
                 Field<Type>&,
                 const GeometricField<Type, fvPatchField, volMesh>&,
                 const labelList& adr,
-                const scalarListList& weights
+                const scalarListList& weights,
+                const CombineOp& cop
             ) const;
 
             //- Interpolate field using cell-point interpolation
-            template<class Type>
+            template<class Type, class CombineOp>
             void interpolateField
             (
                 Field<Type>&,
                 const GeometricField<Type, fvPatchField, volMesh>&,
                 const labelList& adr,
-                const vectorField& centres
+                const vectorField& centres,
+                const CombineOp& cop
             ) const;
 
 
             //- Interpolate internal volume field
-            template<class Type>
+            template<class Type, class CombineOp>
             void interpolateInternalField
             (
                 Field<Type>&,
                 const GeometricField<Type, fvPatchField, volMesh>&,
-                order=INTERPOLATE
+                order=INTERPOLATE,
+                const CombineOp& cop = eqOp<Type>()
             ) const;
 
-            template<class Type>
+            template<class Type, class CombineOp>
             void interpolateInternalField
             (
                 Field<Type>&,
                 const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
-                order=INTERPOLATE
+                order=INTERPOLATE,
+                const CombineOp& cop = eqOp<Type>()
             ) const;
 
 
             //- Interpolate volume field
-            template<class Type>
+            template<class Type, class CombineOp>
             void interpolate
             (
                 GeometricField<Type, fvPatchField, volMesh>&,
                 const GeometricField<Type, fvPatchField, volMesh>&,
-                order=INTERPOLATE
+                order=INTERPOLATE,
+                const CombineOp& cop = eqOp<Type>()
             ) const;
 
-            template<class Type>
+            template<class Type, class CombineOp>
             void interpolate
             (
                 GeometricField<Type, fvPatchField, volMesh>&,
                 const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
-                order=INTERPOLATE
+                order=INTERPOLATE,
+                const CombineOp& cop = eqOp<Type>()
             ) const;
 
 
             //- Interpolate volume field
-            template<class Type>
+            template<class Type, class CombineOp>
             tmp<GeometricField<Type, fvPatchField, volMesh> > interpolate
             (
                 const GeometricField<Type, fvPatchField, volMesh>&,
-                order=INTERPOLATE
+                order=INTERPOLATE,
+                const CombineOp& cop = eqOp<Type>()
             ) const;
 
-            template<class Type>
+            template<class Type, class CombineOp>
             tmp<GeometricField<Type, fvPatchField, volMesh> > interpolate
             (
                 const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
-                order=INTERPOLATE
+                order=INTERPOLATE,
+                const CombineOp& cop = eqOp<Type>()
             ) const;
 };
 
diff --git a/src/sampling/meshToMeshInterpolation/meshToMesh/meshToMeshInterpolate.C b/src/sampling/meshToMeshInterpolation/meshToMesh/meshToMeshInterpolate.C
index 991ae31826..7ed543ed3c 100644
--- a/src/sampling/meshToMeshInterpolation/meshToMesh/meshToMeshInterpolate.C
+++ b/src/sampling/meshToMeshInterpolation/meshToMesh/meshToMeshInterpolate.C
@@ -31,12 +31,13 @@ License
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-template<class Type>
+template<class Type, class CombineOp>
 void Foam::meshToMesh::mapField
 (
     Field<Type>& toF,
     const Field<Type>& fromVf,
-    const labelList& adr
+    const labelList& adr,
+    const CombineOp& cop
 ) const
 {
     // Direct mapping of nearest-cell values
@@ -45,7 +46,7 @@ void Foam::meshToMesh::mapField
     {
         if (adr[celli] != -1)
         {
-            toF[celli] = fromVf[adr[celli]];
+            cop(toF[celli], fromVf[adr[celli]]);
         }
     }
 
@@ -53,13 +54,14 @@ void Foam::meshToMesh::mapField
 }
 
 
-template<class Type>
+template<class Type, class CombineOp>
 void Foam::meshToMesh::interpolateField
 (
     Field<Type>& toF,
     const GeometricField<Type, fvPatchField, volMesh>& fromVf,
     const labelList& adr,
-    const scalarListList& weights
+    const scalarListList& weights,
+    const CombineOp& cop
 ) const
 {
     // Inverse distance weighted interpolation
@@ -74,24 +76,27 @@ void Foam::meshToMesh::interpolateField
             const labelList& neighbours = cc[adr[celli]];
             const scalarList& w = weights[celli];
 
-            toF[celli] = fromVf[adr[celli]]*w[0];
+            Type f = fromVf[adr[celli]]*w[0];
 
             for (label ni = 1; ni < w.size(); ni++)
             {
-                toF[celli] += fromVf[neighbours[ni - 1]]*w[ni];
+                f += fromVf[neighbours[ni - 1]]*w[ni];
             }
+
+            cop(toF[celli], f);
         }
     }
 }
 
 
-template<class Type>
+template<class Type, class CombineOp>
 void Foam::meshToMesh::interpolateField
 (
     Field<Type>& toF,
     const GeometricField<Type, fvPatchField, volMesh>& fromVf,
     const labelList& adr,
-    const vectorField& centres
+    const vectorField& centres,
+    const CombineOp& cop
 ) const
 {
     // Cell-Point interpolation
@@ -101,31 +106,36 @@ void Foam::meshToMesh::interpolateField
     {
         if (adr[celli] != -1)
         {
-            toF[celli] = interpolator.interpolate
+            cop
             (
-                centres[celli],
-                adr[celli]
+                toF[celli],
+                interpolator.interpolate
+                (
+                    centres[celli],
+                    adr[celli]
+                )
             );
         }
     }
 }
 
 
-template<class Type>
+template<class Type, class CombineOp>
 void Foam::meshToMesh::interpolateInternalField
 (
     Field<Type>& toF,
     const GeometricField<Type, fvPatchField, volMesh>& fromVf,
-    meshToMesh::order ord
+    meshToMesh::order ord,
+    const CombineOp& cop
 ) const
 {
     if (fromVf.mesh() != fromMesh_)
     {
         FatalErrorIn
         (
-            "meshToMesh::interpolateInternalField(Field<Type>& toF, "
-            "const GeometricField<Type, fvPatchField, volMesh>& fromVf, "
-            "meshToMesh::order ord) const"
+            "meshToMesh::interpolateInternalField(Field<Type>&, "
+            "const GeometricField<Type, fvPatchField, volMesh>&, "
+            "meshToMesh::order, const CombineOp&) const"
         )   << "the argument field does not correspond to the right mesh. "
             << "Field size: " << fromVf.size()
             << " mesh size: " << fromMesh_.nCells()
@@ -136,9 +146,9 @@ void Foam::meshToMesh::interpolateInternalField
     {
         FatalErrorIn
         (
-            "meshToMesh::interpolateInternalField(Field<Type>& toF, "
-            "const GeometricField<Type, fvPatchField, volMesh>& fromVf, "
-            "meshToMesh::order ord) const"
+            "meshToMesh::interpolateInternalField(Field<Type>&, "
+            "const GeometricField<Type, fvPatchField, volMesh>&, "
+            "meshToMesh::order, const CombineOp&) const"
         )   << "the argument field does not correspond to the right mesh. "
             << "Field size: " << toF.size()
             << " mesh size: " << toMesh_.nCells()
@@ -148,7 +158,7 @@ void Foam::meshToMesh::interpolateInternalField
     switch(ord)
     {
         case MAP:
-            mapField(toF, fromVf, cellAddressing_);
+            mapField(toF, fromVf, cellAddressing_, cop);
         break;
 
         case INTERPOLATE:
@@ -157,7 +167,8 @@ void Foam::meshToMesh::interpolateInternalField
                 toF,
                 fromVf,
                 cellAddressing_,
-                inverseDistanceWeights()
+                inverseDistanceWeights(),
+                cop
             );
         break;
 
@@ -167,44 +178,47 @@ void Foam::meshToMesh::interpolateInternalField
                 toF,
                 fromVf,
                 cellAddressing_,
-                toMesh_.cellCentres()
+                toMesh_.cellCentres(),
+                cop
             );
         break;
 
         default:
             FatalErrorIn
             (
-                "meshToMesh::interpolateInternalField(Field<Type>& toF, "
-                "const GeometricField<Type, fvPatchField, volMesh>& fromVf, "
-                "meshToMesh::order ord) const"
+                "meshToMesh::interpolateInternalField(Field<Type>&, "
+                "const GeometricField<Type, fvPatchField, volMesh>&, "
+                "meshToMesh::order, const CombineOp&) const"
             )   << "unknown interpolation scheme " << ord
                 << exit(FatalError);
     }
 }
 
 
-template<class Type>
+template<class Type, class CombineOp>
 void Foam::meshToMesh::interpolateInternalField
 (
     Field<Type>& toF,
     const tmp<GeometricField<Type, fvPatchField, volMesh> >& tfromVf,
-    meshToMesh::order ord
+    meshToMesh::order ord,
+    const CombineOp& cop
 ) const
 {
-    interpolateInternalField(toF, tfromVf(), ord);
+    interpolateInternalField(toF, tfromVf(), ord, cop);
     tfromVf.clear();
 }
 
 
-template<class Type>
+template<class Type, class CombineOp>
 void Foam::meshToMesh::interpolate
 (
     GeometricField<Type, fvPatchField, volMesh>& toVf,
     const GeometricField<Type, fvPatchField, volMesh>& fromVf,
-    meshToMesh::order ord
+    meshToMesh::order ord,
+    const CombineOp& cop
 ) const
 {
-    interpolateInternalField(toVf, fromVf, ord);
+    interpolateInternalField(toVf, fromVf, ord, cop);
 
     forAll(toMesh_.boundaryMesh(), patchi)
     {
@@ -219,7 +233,8 @@ void Foam::meshToMesh::interpolate
                     (
                         toVf.boundaryField()[patchi],
                         fromVf,
-                        boundaryAddressing_[patchi]
+                        boundaryAddressing_[patchi],
+                        cop
                     );
                 break;
 
@@ -229,7 +244,8 @@ void Foam::meshToMesh::interpolate
                         toVf.boundaryField()[patchi],
                         fromVf,
                         boundaryAddressing_[patchi],
-                        toPatch.Cf()
+                        toPatch.Cf(),
+                        cop
                     );
                 break;
 
@@ -239,7 +255,8 @@ void Foam::meshToMesh::interpolate
                         toVf.boundaryField()[patchi],
                         fromVf,
                         boundaryAddressing_[patchi],
-                        toPatch.Cf()
+                        toPatch.Cf(),
+                        cop
                     );
                 break;
 
@@ -247,9 +264,9 @@ void Foam::meshToMesh::interpolate
                     FatalErrorIn
                     (
                         "meshToMesh::interpolate("
-                        "GeometricField<Type, fvPatchField, volMesh>& toVf, "
-                        "const GeometricField<Type, fvPatchField, volMesh>& "
-                        "fromVf, meshToMesh::order ord) const"
+                        "GeometricField<Type, fvPatchField, volMesh>&, "
+                        "const GeometricField<Type, fvPatchField, volMesh>&, "
+                        "meshToMesh::order, const CombineOp&) const"
                     )   << "unknown interpolation scheme " << ord
                         << exit(FatalError);
             }
@@ -286,37 +303,40 @@ void Foam::meshToMesh::interpolate
                 [
                     fromMeshPatches_.find(patchMap_.find(toPatch.name())())()
                 ],
-                boundaryAddressing_[patchi]
+                boundaryAddressing_[patchi],
+                cop
             );
         }
     }
 }
 
 
-template<class Type>
+template<class Type, class CombineOp>
 void Foam::meshToMesh::interpolate
 (
     GeometricField<Type, fvPatchField, volMesh>& toVf,
     const tmp<GeometricField<Type, fvPatchField, volMesh> >& tfromVf,
-    meshToMesh::order ord
+    meshToMesh::order ord,
+    const CombineOp& cop
 ) const
 {
-    interpolate(toVf, tfromVf(), ord);
+    interpolate(toVf, tfromVf(), ord, cop);
     tfromVf.clear();
 }
 
 
-template<class Type>
+template<class Type, class CombineOp>
 Foam::tmp< Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh> >
 Foam::meshToMesh::interpolate
 (
     const GeometricField<Type, fvPatchField, volMesh>& fromVf,
-    meshToMesh::order ord
+    meshToMesh::order ord,
+    const CombineOp& cop
 ) const
 {
     // Create and map the internal-field values
     Field<Type> internalField(toMesh_.nCells());
-    interpolateInternalField(internalField, fromVf, ord);
+    interpolateInternalField(internalField, fromVf, ord, cop);
 
     // check whether both meshes have got the same number
     // of boundary patches
@@ -325,8 +345,8 @@ Foam::meshToMesh::interpolate
         FatalErrorIn
         (
             "meshToMesh::interpolate"
-            "(const GeometricField<Type, fvPatchField, volMesh>& fromVf,"
-            "meshToMesh::order ord) const"
+            "(const GeometricField<Type, fvPatchField, volMesh>&,"
+            "meshToMesh::order, const CombineOp&) const"
         )   << "Incompatible meshes: different number of boundaries, "
                "only internal field may be interpolated"
             << exit(FatalError);
@@ -381,16 +401,17 @@ Foam::meshToMesh::interpolate
 }
 
 
-template<class Type>
+template<class Type, class CombineOp>
 Foam::tmp< Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh> >
 Foam::meshToMesh::interpolate
 (
     const tmp<GeometricField<Type, fvPatchField, volMesh> >& tfromVf,
-    meshToMesh::order ord
+    meshToMesh::order ord,
+    const CombineOp& cop
 ) const
 {
     tmp<GeometricField<Type, fvPatchField, volMesh> > tint =
-        interpolate(tfromVf(), ord);
+        interpolate(tfromVf(), ord, cop);
     tfromVf.clear();
 
     return tint;

From 06699987a8afe1465064c5484b519cb60191177b Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Thu, 22 Dec 2011 11:40:40 +0000
Subject: [PATCH 07/18] ENH: redistributePar: print decomposition stats upon
 entry

---
 .../parallelProcessing/redistributePar/redistributePar.C   | 7 +++----
 1 file changed, 3 insertions(+), 4 deletions(-)

diff --git a/applications/utilities/parallelProcessing/redistributePar/redistributePar.C b/applications/utilities/parallelProcessing/redistributePar/redistributePar.C
index 5eaff912ef..cde7f5b12d 100644
--- a/applications/utilities/parallelProcessing/redistributePar/redistributePar.C
+++ b/applications/utilities/parallelProcessing/redistributePar/redistributePar.C
@@ -795,9 +795,9 @@ int main(int argc, char *argv[])
     );
     fvMesh& mesh = meshPtr();
 
-    //Pout<< "Read mesh:" << endl;
-    //printMeshData(Pout, mesh);
-    //Pout<< endl;
+    // Print some statistics
+    Info<< "Before distribution:" << endl;
+    printMeshData(mesh);
 
 
 
@@ -1022,7 +1022,6 @@ int main(int argc, char *argv[])
     //map().distributeFaceData(faceCc);
 
 
-    // Print a bit
     // Print some statistics
     Info<< "After distribution:" << endl;
     printMeshData(mesh);

From b752c60cd1850b5f4b07d83cf2de7895353844af Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Thu, 22 Dec 2011 11:46:14 +0000
Subject: [PATCH 08/18] ENH: streamLine: added preliminary support for surface
 following

---
 .../field/streamLine/streamLine.C             | 239 +++++++++++++++++-
 .../field/streamLine/streamLine.H             |  11 +
 .../field/streamLine/streamLineParticle.C     |  78 +++++-
 .../field/streamLine/streamLineParticle.H     |  29 ++-
 4 files changed, 336 insertions(+), 21 deletions(-)

diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLine.C b/src/postProcessing/functionObjects/field/streamLine/streamLine.C
index 30d8052513..ab815a8df8 100644
--- a/src/postProcessing/functionObjects/field/streamLine/streamLine.C
+++ b/src/postProcessing/functionObjects/field/streamLine/streamLine.C
@@ -34,6 +34,7 @@ License
 #include "globalIndex.H"
 #include "mapDistribute.H"
 #include "interpolationCellPoint.H"
+#include "PatchTools.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -41,6 +42,57 @@ defineTypeNameAndDebug(Foam::streamLine, 0);
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+Foam::autoPtr<Foam::indirectPrimitivePatch>
+Foam::streamLine::wallPatch() const
+{
+    const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
+
+    const polyBoundaryMesh& patches = mesh.boundaryMesh();
+
+    label nFaces = 0;
+
+    forAll(patches, patchI)
+    {
+        //if (!polyPatch::constraintType(patches[patchI].type()))
+        if (isA<wallPolyPatch>(patches[patchI]))
+        {
+            nFaces += patches[patchI].size();
+        }
+    }
+
+    labelList addressing(nFaces);
+
+    nFaces = 0;
+
+    forAll(patches, patchI)
+    {
+        //if (!polyPatch::constraintType(patches[patchI].type()))
+        if (isA<wallPolyPatch>(patches[patchI]))
+        {
+            const polyPatch& pp = patches[patchI];
+
+            forAll(pp, i)
+            {
+                addressing[nFaces++] = pp.start()+i;
+            }
+        }
+    }
+
+    return autoPtr<indirectPrimitivePatch>
+    (
+        new indirectPrimitivePatch
+        (
+            IndirectList<face>
+            (
+                mesh.faces(),
+                addressing
+            ),
+            mesh.points()
+        )
+    );
+}
+
+
 void Foam::streamLine::track()
 {
     const Time& runTime = obr_.time();
@@ -72,7 +124,7 @@ void Foam::streamLine::track()
 
     label nSeeds = returnReduce(particles.size(), sumOp<label>());
 
-    Info<< "Seeded " << nSeeds << " particles." << endl;
+    Info<< type() << " : seeded " << nSeeds << " particles." << endl;
 
     // Read or lookup fields
     PtrList<volScalarField> vsFlds;
@@ -163,7 +215,6 @@ void Foam::streamLine::track()
                 vsInterp.set
                 (
                     nScalar++,
-                    //new interpolationCellPoint<scalar>(f)
                     interpolation<scalar>::New
                     (
                         interpolationScheme_,
@@ -186,7 +237,6 @@ void Foam::streamLine::track()
                 vvInterp.set
                 (
                     nVector++,
-                    //new interpolationCellPoint<vector>(f)
                     interpolation<vector>::New
                     (
                         interpolationScheme_,
@@ -241,6 +291,120 @@ void Foam::streamLine::track()
         allVectors_[i].setCapacity(nSeeds);
     }
 
+
+    //- For surface following: per face the normal to constrain the velocity
+    //  with.
+    pointField patchNormals;
+
+
+    if (followSurface_)
+    {
+        // Determine geometric data on the non-constraint patches
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        // - point normals
+        // - neighbouring cells
+        autoPtr<indirectPrimitivePatch> boundaryPatch(wallPatch());
+
+        // Calculate parallel consistent normal (on boundaryPatch points only)
+        pointField boundaryNormals
+        (
+            PatchTools::pointNormals
+            (
+                mesh,
+                boundaryPatch(),
+                boundaryPatch().addressing()
+            )
+        );
+        // Make sure all points know about point normal. There might be points
+        // which are not on boundaryPatch() but coupled to points that are.
+        pointField pointNormals(mesh.nPoints(), vector::zero);
+        UIndirectList<point>
+        (
+            pointNormals,
+            boundaryPatch().meshPoints()
+        ) = boundaryNormals;
+        boundaryNormals.clear();
+
+        mesh.globalData().syncPointData
+        (
+            pointNormals,
+            maxMagSqrEqOp<point>(),
+            mapDistribute::transform()
+        );
+
+        // Calculate per-face a single constraint normal:
+        // - faces on patch: face-normal
+        // - faces on cell with point on patch: point-normal
+        // - rest: vector::zero
+        //
+        // Note: the disadvantage is that we can only handle faces with one
+        //       point on the patch - otherwise we might pick up the wrong
+        //       velocity.
+
+
+
+        patchNormals.setSize(mesh.nFaces(), vector::zero);
+        label nPointNormals = 0;
+
+        // 1. faces on patch
+        UIndirectList<point>(patchNormals, boundaryPatch().addressing()) =
+             boundaryPatch().faceNormals();
+
+        // 2. faces with a point on the patch
+        forAll(mesh.faces(), faceI)
+        {
+            if (patchNormals[faceI] == vector::zero)
+            {
+                const face& f = mesh.faces()[faceI];
+                forAll(f, fp)
+                {
+                    label pointI = f[fp];
+                    if (pointNormals[pointI] != vector::zero)
+                    {
+                        patchNormals[faceI] = pointNormals[pointI];
+                        nPointNormals++;
+                    }
+                }
+            }
+        }
+
+        // 3. cells on faces with a point on the patch
+        forAll(mesh.points(), pointI)
+        {
+            if (pointNormals[pointI] != vector::zero)
+            {
+                const labelList& pCells = mesh.pointCells(pointI);
+                forAll(pCells, i)
+                {
+                    const cell& cFaces = mesh.cells()[pCells[i]];
+                    forAll(cFaces, j)
+                    {
+                        label faceI = cFaces[j];
+
+                        if (patchNormals[faceI] == vector::zero)
+                        {
+                            patchNormals[faceI] = pointNormals[pointI];
+                            nPointNormals++;
+                        }
+                    }
+                }
+            }
+        }
+
+        Info<< type() << " : calculated constrained-tracking normals" << nl
+            << "Number of faces in mesh                     : "
+            << returnReduce(mesh.nFaces(), sumOp<label>()) << nl
+            << "    of which on walls                       : "
+            <<  returnReduce
+                (
+                    boundaryPatch().addressing().size(),
+                    sumOp<label>()
+                ) << nl
+            << "    of which on cell with point on boundary : "
+            << returnReduce(nPointNormals, sumOp<label>()) << nl
+            << endl;
+    }
+
     // additional particle info
     streamLineParticle::trackingData td
     (
@@ -248,8 +412,16 @@ void Foam::streamLine::track()
         vsInterp,
         vvInterp,
         UIndex,         // index of U in vvInterp
-        trackForward_,  // track in +u direction
-        nSubCycle_,     // step through cells in steps?
+        trackForward_,  // track in +u direction?
+        nSubCycle_,     // automatic track control:step through cells in steps?
+        trackLength_,   // fixed track length
+
+        //- For surface-constrained streamlines
+        followSurface_, // stay on surface?
+        patchNormals,   // per face normal to constrain tracking to
+        //pointNormals,   // per point          ,,
+        trackHeight_,   // track height
+
         allTracks_,
         allScalars_,
         allVectors_
@@ -279,7 +451,8 @@ Foam::streamLine::streamLine
     name_(name),
     obr_(obr),
     loadFromFiles_(loadFromFiles),
-    active_(true)
+    active_(true),
+    nSubCycle_(0)
 {
     // Only active if a fvMesh is available
     if (isA<fvMesh>(obr_))
@@ -334,6 +507,16 @@ void Foam::streamLine::read(const dictionary& dict)
                 dict.lookup("U") >> UName_;
             }
         }
+
+        if (findIndex(fields_, UName_) == -1)
+        {
+            FatalIOErrorIn("streamLine::read(const dictionary&)", dict)
+                << "Velocity field for tracking " << UName_
+                << " should be present in the list of fields " << fields_
+                << exit(FatalIOError);
+        }
+
+
         dict.lookup("trackForward") >> trackForward_;
         dict.lookup("lifeTime") >> lifeTime_;
         if (lifeTime_ < 1)
@@ -343,10 +526,48 @@ void Foam::streamLine::read(const dictionary& dict)
                 << exit(FatalError);
         }
 
-        dict.lookup("nSubCycle") >> nSubCycle_;
-        if (nSubCycle_ < 1)
+
+        bool subCycling = dict.found("nSubCycle");
+        bool fixedLength = dict.found("trackLength");
+
+        if (subCycling && fixedLength)
         {
-            nSubCycle_ = 1;
+            FatalIOErrorIn("streamLine::read(const dictionary&)", dict)
+                << "Cannot both specify automatic time stepping (through '"
+                << "nSubCycle' specification) and fixed track length (through '"
+                << "trackLength')"
+                << exit(FatalIOError);
+        }
+
+
+        nSubCycle_ = 1;
+        if (dict.readIfPresent("nSubCycle", nSubCycle_))
+        {
+            trackLength_ = VGREAT;
+            if (nSubCycle_ < 1)
+            {
+                nSubCycle_ = 1;
+            }
+            Info<< type() << " : automatic track length specified through"
+                << " number of sub cycles : " << nSubCycle_ << nl << endl;
+        }
+        else
+        {
+            dict.lookup("trackLength") >> trackLength_;
+
+            Info<< type() << " : fixed track length specified : "
+                << trackLength_ << nl << endl;
+        }
+
+
+        followSurface_ = false;
+        trackHeight_ = VGREAT;
+        if (dict.readIfPresent("followSurface", followSurface_))
+        {
+            dict.lookup("trackHeight") >> trackHeight_;
+
+            Info<< type() << " : surface-constrained streamline at "
+                << trackHeight_ << " distance above the surface" << nl << endl;
         }
 
         interpolationScheme_ = dict.lookupOrDefault
diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLine.H b/src/postProcessing/functionObjects/field/streamLine/streamLine.H
index 9fc981034c..77984ed2c8 100644
--- a/src/postProcessing/functionObjects/field/streamLine/streamLine.H
+++ b/src/postProcessing/functionObjects/field/streamLine/streamLine.H
@@ -43,6 +43,7 @@ SourceFiles
 #include "vectorList.H"
 #include "polyMesh.H"
 #include "writer.H"
+#include "indirectPrimitivePatch.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -95,9 +96,17 @@ class streamLine
         //- Maximum lifetime (= number of cells) of particle
         label lifeTime_;
 
+        //- Surface-constrained?
+        bool followSurface_;
+
+        scalar trackHeight_;
+
         //- Number of subcycling steps
         label nSubCycle_;
 
+        //- Track length
+        scalar trackLength_;
+
         //- Optional specified name of particles
         word cloudName_;
 
@@ -141,6 +150,8 @@ class streamLine
             List<DynamicList<vectorList> > allVectors_;
 
 
+        //- Construct patch out of all wall patch faces
+        autoPtr<indirectPrimitivePatch> wallPatch() const;
 
         //- Do all seeding and tracking
         void track();
diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.C b/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.C
index c7808b9900..279e17ef80 100644
--- a/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.C
+++ b/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.C
@@ -43,7 +43,7 @@ Foam::scalar Foam::streamLineParticle::calcSubCycleDeltaT
     const vector& U
 ) const
 {
-    streamLineParticle testParticle(*this);
+    particle testParticle(*this);
 
     bool oldKeepParticle = td.keepParticle;
     bool oldSwitchProcessor = td.switchProcessor;
@@ -59,7 +59,8 @@ Foam::vector Foam::streamLineParticle::interpolateFields
 (
     const trackingData& td,
     const point& position,
-    const label cellI
+    const label cellI,
+    const label faceI
 )
 {
     if (cellI == -1)
@@ -76,7 +77,8 @@ Foam::vector Foam::streamLineParticle::interpolateFields
             td.vsInterp_[scalarI].interpolate
             (
                 position,
-                cellI
+                cellI,
+                faceI
             )
         );
     }
@@ -89,7 +91,8 @@ Foam::vector Foam::streamLineParticle::interpolateFields
             td.vvInterp_[vectorI].interpolate
             (
                 position,
-                cellI
+                cellI,
+                faceI
             )
         );
     }
@@ -201,13 +204,24 @@ bool Foam::streamLineParticle::move
 
             // Store current position and sampled velocity.
             sampledPositions_.append(position());
-            vector U = interpolateFields(td, position(), cell());
+            vector U = interpolateFields(td, position(), cell(), face());
 
             if (!td.trackForward_)
             {
                 U = -U;
             }
 
+            if (td.followSurface_)
+            {
+                //- Simple: remove wall-normal component. Should keep particle
+                //  at same height. Does not work well on warped faces.
+                const vector& n = td.patchNormals_[tetFaceI_];
+                if (n != vector::zero)
+                {
+                    U -= (U&n)*n;
+                }
+            }
+
             scalar magU = mag(U);
 
             if (magU < SMALL)
@@ -219,7 +233,13 @@ bool Foam::streamLineParticle::move
 
             U /= magU;
 
-            if (subIter == 0 && td.nSubCycle_ > 1)
+            if (td.trackLength_ < GREAT)
+            {
+                dt = td.trackLength_;
+                //Pout<< "    subiteration " << subIter
+                //    << " : fixed length: updated dt:" << dt << endl;
+            }
+            else if (subIter == 0 && td.nSubCycle_ > 1)
             {
                 // Adapt dt to cross cell in a few steps
                 dt = calcSubCycleDeltaT(td, dt, U);
@@ -264,7 +284,8 @@ bool Foam::streamLineParticle::move
             if (debug)
             {
                 Pout<< "streamLineParticle : Removing stagnant particle:"
-                    << p << " sampled positions:" << sampledPositions_.size()
+                    << p.position()
+                    << " sampled positions:" << sampledPositions_.size()
                     << endl;
             }
             td.keepParticle = false;
@@ -273,12 +294,13 @@ bool Foam::streamLineParticle::move
         {
             // Normal exit. Store last position and fields
             sampledPositions_.append(position());
-            interpolateFields(td, position(), cell());
+            interpolateFields(td, position(), cell(), face());
 
             if (debug)
             {
                 Pout<< "streamLineParticle : Removing particle:"
-                    << p << " sampled positions:" << sampledPositions_.size()
+                    << p.position()
+                    << " sampled positions:" << sampledPositions_.size()
                     << endl;
             }
         }
@@ -374,8 +396,42 @@ void Foam::streamLineParticle::hitWallPatch
     const tetIndices&
 )
 {
-    // Remove particle
-    td.keepParticle = false;
+    if (td.followSurface_)
+    {
+        // Push slightly in. Since we're pushing to tet centre we should still
+        // be in the same tet.
+        tetPointRef tet = currentTet();
+
+        const point oldPosition(position());
+
+        position() += trackingCorrectionTol*(tet.centre() - position());
+
+        if (debug)
+        {
+            label patchID = patch(face());
+
+            Pout<< "hitWallPatch : on wall "
+                << mesh().boundaryMesh()[patchID].name()
+                << " moved from " << oldPosition
+                << " to " << position()
+                << " due to centre " << tet.centre() << endl;
+        }
+
+        // Check that still in tet (should always be the case
+        if (debug && !tet.inside(position()))
+        {
+            FatalErrorIn("streamLineParticle::hitWallPatch()")
+                << "Pushing particle " << *this
+                << " away from wall " << wpp.name()
+                << " moved it outside of tet."
+                << exit(FatalError);
+        }
+    }
+    else
+    {
+        // Remove particle
+        td.keepParticle = false;
+    }
 }
 
 
diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.H b/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.H
index b1ee9acf88..c5feff4285 100644
--- a/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.H
+++ b/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.H
@@ -73,6 +73,13 @@ public:
         const label UIndex_;
         const bool trackForward_;
         const label nSubCycle_;
+        const scalar trackLength_;
+
+        //- For surface-constrained tracking
+        const bool followSurface_;
+        const pointField& patchNormals_;
+        const scalar trackHeight_;
+
 
         DynamicList<vectorList>& allPositions_;
         List<DynamicList<scalarList> >& allScalars_;
@@ -89,6 +96,12 @@ public:
                 const label UIndex,
                 const bool trackForward,
                 const label nSubCycle,
+                const scalar trackLength,
+                const bool followSurface,
+
+                const pointField& patchNormals,
+                const scalar trackHeight,
+
                 DynamicList<List<point> >& allPositions,
                 List<DynamicList<scalarList> >& allScalars,
                 List<DynamicList<vectorList> >& allVectors
@@ -100,6 +113,12 @@ public:
                 UIndex_(UIndex),
                 trackForward_(trackForward),
                 nSubCycle_(nSubCycle),
+                trackLength_(trackLength),
+                followSurface_(followSurface),
+
+                patchNormals_(patchNormals),
+                trackHeight_(trackHeight),
+
                 allPositions_(allPositions),
                 allScalars_(allScalars),
                 allVectors_(allVectors)
@@ -135,12 +154,20 @@ private:
             const vector& U
         ) const;
 
+        void constrainVelocity
+        (
+            trackingData& td,
+            const scalar dt,
+            vector& U
+        );
+
         //- Interpolate all quantities; return interpolated velocity.
         vector interpolateFields
         (
             const trackingData&,
             const point&,
-            const label cellI
+            const label cellI,
+            const label faceI
         );
 
 

From 2d7e76951e12aad43bf7820d65ea05682b1f2864 Mon Sep 17 00:00:00 2001
From: andy <andy>
Date: Thu, 22 Dec 2011 12:26:33 +0000
Subject: [PATCH 09/18] ENH: Provide message if ParaView is not installed
 instead of swallowing the error

---
 .../utilities/postProcessing/graphics/PV3Readers/Allwmake     | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake b/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake
index 6851210bd7..1568f8afb8 100755
--- a/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake
+++ b/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake
@@ -1,6 +1,6 @@
 #!/bin/sh
 cd ${0%/*} || exit 1    # run from this directory
-set -x
+#set -x
 
 if [ -d "$ParaView_DIR" -a -r "$ParaView_DIR" ]
 then
@@ -15,6 +15,8 @@ then
     wmake libso vtkPV3Readers
     PV3blockMeshReader/Allwmake
     PV3FoamReader/Allwmake
+else
+    echo "ERROR: ParaView not found in $ParaView_DIR"
 fi
 
 # ----------------------------------------------------------------- end-of-file

From ddfd2bfb382670ce1a666dd74559228a002175a9 Mon Sep 17 00:00:00 2001
From: andy <andy>
Date: Thu, 22 Dec 2011 12:27:03 +0000
Subject: [PATCH 10/18] ENH: Updated wallShearStress utility to include
 compressible RAS flows

---
 .../wall/wallShearStress/Make/options         |   7 +-
 .../wall/wallShearStress/createFields.H       |  22 ---
 .../wall/wallShearStress/wallShearStress.C    | 153 ++++++++++++++++--
 3 files changed, 144 insertions(+), 38 deletions(-)
 delete mode 100644 applications/utilities/postProcessing/wall/wallShearStress/createFields.H

diff --git a/applications/utilities/postProcessing/wall/wallShearStress/Make/options b/applications/utilities/postProcessing/wall/wallShearStress/Make/options
index 8862565863..b386fa4540 100644
--- a/applications/utilities/postProcessing/wall/wallShearStress/Make/options
+++ b/applications/utilities/postProcessing/wall/wallShearStress/Make/options
@@ -1,11 +1,14 @@
 EXE_INC = \
     -I$(LIB_SRC)/transportModels \
     -I$(LIB_SRC)/turbulenceModels \
-    -I$(LIB_SRC)/turbulenceModels/incompressible/RAS/RASModel \
+    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 
 EXE_LIBS = \
-    -lincompressibleRASModels \
     -lincompressibleTransportModels \
+    -lincompressibleRASModels \
+    -lbasicThermophysicalModels \
+    -lspecie \
+    -lcompressibleRASModels \
     -lfiniteVolume \
     -lgenericPatchFields
diff --git a/applications/utilities/postProcessing/wall/wallShearStress/createFields.H b/applications/utilities/postProcessing/wall/wallShearStress/createFields.H
deleted file mode 100644
index d069c4545d..0000000000
--- a/applications/utilities/postProcessing/wall/wallShearStress/createFields.H
+++ /dev/null
@@ -1,22 +0,0 @@
-    Info<< "Reading field U\n" << endl;
-    volVectorField U
-    (
-        IOobject
-        (
-            "U",
-            runTime.timeName(),
-            mesh,
-            IOobject::MUST_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh
-    );
-
-#   include "createPhi.H"
-
-    singlePhaseTransportModel laminarTransport(U, phi);
-
-    autoPtr<incompressible::RASModel> RASModel
-    (
-        incompressible::RASModel::New(U, phi, laminarTransport)
-    );
diff --git a/applications/utilities/postProcessing/wall/wallShearStress/wallShearStress.C b/applications/utilities/postProcessing/wall/wallShearStress/wallShearStress.C
index 7d440417df..6069837b95 100644
--- a/applications/utilities/postProcessing/wall/wallShearStress/wallShearStress.C
+++ b/applications/utilities/postProcessing/wall/wallShearStress/wallShearStress.C
@@ -25,23 +25,130 @@ Application
     wallShearStress
 
 Description
-    Calculates and writes the wall shear stress, for the specified times.
+    Calculates and reports wall shear stress for all patches, for the
+    specified times when using RAS turbulence models.
+
+    Default behaviour assumes operating in incompressible mode.
+    Use the -compressible option for compressible RAS cases.
 
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
+
 #include "incompressible/singlePhaseTransportModel/singlePhaseTransportModel.H"
-#include "RASModel.H"
+#include "incompressible/RAS/RASModel/RASModel.H"
+
+#include "basicPsiThermo.H"
+#include "compressible/RAS/RASModel/RASModel.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+void calcIncompressible
+(
+    const fvMesh& mesh,
+    const Time& runTime,
+    const volVectorField& U,
+    volVectorField& wallShearStress
+)
+{
+    #include "createPhi.H"
+
+    singlePhaseTransportModel laminarTransport(U, phi);
+
+    autoPtr<incompressible::RASModel> model
+    (
+        incompressible::RASModel::New(U, phi, laminarTransport)
+    );
+
+    const volSymmTensorField Reff(model->devReff());
+
+    forAll(wallShearStress.boundaryField(), patchI)
+    {
+        wallShearStress.boundaryField()[patchI] =
+        (
+           -mesh.Sf().boundaryField()[patchI]
+           /mesh.magSf().boundaryField()[patchI]
+        ) & Reff.boundaryField()[patchI];
+    }
+}
+
+
+void calcCompressible
+(
+    const fvMesh& mesh,
+    const Time& runTime,
+    const volVectorField& U,
+    volVectorField& wallShearStress
+)
+{
+    IOobject rhoHeader
+    (
+        "rho",
+        runTime.timeName(),
+        mesh,
+        IOobject::MUST_READ,
+        IOobject::NO_WRITE
+    );
+
+    if (!rhoHeader.headerOk())
+    {
+        Info<< "    no rho field" << endl;
+        return;
+    }
+
+    Info<< "Reading field rho\n" << endl;
+    volScalarField rho(rhoHeader, mesh);
+
+    #include "compressibleCreatePhi.H"
+
+    autoPtr<basicPsiThermo> pThermo
+    (
+        basicPsiThermo::New(mesh)
+    );
+    basicPsiThermo& thermo = pThermo();
+
+    autoPtr<compressible::RASModel> model
+    (
+        compressible::RASModel::New
+        (
+            rho,
+            U,
+            phi,
+            thermo
+        )
+    );
+
+    const volSymmTensorField Reff(model->devRhoReff());
+
+    forAll(wallShearStress.boundaryField(), patchI)
+    {
+        wallShearStress.boundaryField()[patchI] =
+        (
+           -mesh.Sf().boundaryField()[patchI]
+           /mesh.magSf().boundaryField()[patchI]
+        ) & Reff.boundaryField()[patchI];
+    }
+}
+
+
 int main(int argc, char *argv[])
 {
     timeSelector::addOptions();
+
+    #include "addRegionOption.H"
+
+    argList::addBoolOption
+    (
+        "compressible",
+        "calculate compressible wall shear stress"
+    );
+
     #include "setRootCase.H"
     #include "createTime.H"
     instantList timeDirs = timeSelector::select0(runTime, args);
-    #include "createMesh.H"
+    #include "createNamedMesh.H"
+
+    const bool compressible = args.optionFound("compressible");
 
     forAll(timeDirs, timeI)
     {
@@ -49,10 +156,6 @@ int main(int argc, char *argv[])
         Info<< "Time = " << runTime.timeName() << endl;
         mesh.readUpdate();
 
-        #include "createFields.H"
-
-        volSymmTensorField Reff(RASModel->devReff());
-
         volVectorField wallShearStress
         (
             IOobject
@@ -67,19 +170,41 @@ int main(int argc, char *argv[])
             dimensionedVector
             (
                 "wallShearStress",
-                Reff.dimensions(),
+                sqr(dimLength)/sqr(dimTime),
                 vector::zero
             )
         );
 
-        forAll(wallShearStress.boundaryField(), patchi)
+        IOobject UHeader
+        (
+            "U",
+            runTime.timeName(),
+            mesh,
+            IOobject::MUST_READ,
+            IOobject::NO_WRITE
+        );
+
+        if (UHeader.headerOk())
         {
-            wallShearStress.boundaryField()[patchi] =
-            (
-                -mesh.Sf().boundaryField()[patchi]
-                /mesh.magSf().boundaryField()[patchi]
-            ) & Reff.boundaryField()[patchi];
+            Info<< "Reading field U\n" << endl;
+            volVectorField U(UHeader, mesh);
+
+            if (compressible)
+            {
+                calcCompressible(mesh, runTime, U, wallShearStress);
+            }
+            else
+            {
+                calcIncompressible(mesh, runTime, U, wallShearStress);
+            }
         }
+        else
+        {
+            Info<< "    no U field" << endl;
+        }
+
+        Info<< "Writing wall shear stress to field " << wallShearStress.name()
+            << nl << endl;
 
         wallShearStress.write();
     }

From 9ed246733b308d0bfdbf90fcd899a6d101e36bec Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Thu, 22 Dec 2011 14:54:58 +0000
Subject: [PATCH 11/18] COMP: renumberMethods: move into subdirectory

---
 src/Allwmake                                  |   2 +-
 src/renumber/Allwmake                         |  30 +++
 .../CuthillMcKeeRenumber.C                    | 106 +++++++++
 .../CuthillMcKeeRenumber.H                    | 118 ++++++++++
 src/renumber/renumberMethods/Make/files       |   8 +
 src/renumber/renumberMethods/Make/options     |   7 +
 .../boundaryFirstRenumber.C                   | 201 ++++++++++++++++++
 .../boundaryFirstRenumber.H                   | 129 +++++++++++
 .../manualRenumber/manualRenumber.C           | 137 ++++++++++++
 .../manualRenumber/manualRenumber.H           | 125 +++++++++++
 .../randomRenumber/randomRenumber.C           |  96 +++++++++
 .../randomRenumber/randomRenumber.H           | 108 ++++++++++
 .../renumberMethod/renumberMethod.C           | 131 ++++++++++++
 .../renumberMethod/renumberMethod.H           | 161 ++++++++++++++
 .../springRenumber/springRenumber.C           | 172 +++++++++++++++
 .../springRenumber/springRenumber.H           | 132 ++++++++++++
 16 files changed, 1662 insertions(+), 1 deletion(-)
 create mode 100755 src/renumber/Allwmake
 create mode 100644 src/renumber/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
 create mode 100644 src/renumber/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
 create mode 100644 src/renumber/renumberMethods/Make/files
 create mode 100644 src/renumber/renumberMethods/Make/options
 create mode 100644 src/renumber/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
 create mode 100644 src/renumber/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H
 create mode 100644 src/renumber/renumberMethods/manualRenumber/manualRenumber.C
 create mode 100644 src/renumber/renumberMethods/manualRenumber/manualRenumber.H
 create mode 100644 src/renumber/renumberMethods/randomRenumber/randomRenumber.C
 create mode 100644 src/renumber/renumberMethods/randomRenumber/randomRenumber.H
 create mode 100644 src/renumber/renumberMethods/renumberMethod/renumberMethod.C
 create mode 100644 src/renumber/renumberMethods/renumberMethod/renumberMethod.H
 create mode 100644 src/renumber/renumberMethods/springRenumber/springRenumber.C
 create mode 100644 src/renumber/renumberMethods/springRenumber/springRenumber.H

diff --git a/src/Allwmake b/src/Allwmake
index a0622b30f7..97023aa700 100755
--- a/src/Allwmake
+++ b/src/Allwmake
@@ -46,7 +46,7 @@ wmake $makeType genericPatchFields
 # Build the proper scotchDecomp, metisDecomp etc.
 parallel/Allwmake $*
 
-wmake $makeType renumberMethods
+renumber/Allwmake $*
 
 wmake $makeType conversion
 
diff --git a/src/renumber/Allwmake b/src/renumber/Allwmake
new file mode 100755
index 0000000000..36c057dafe
--- /dev/null
+++ b/src/renumber/Allwmake
@@ -0,0 +1,30 @@
+#!/bin/sh
+cd ${0%/*} || exit 1    # run from this directory
+makeType=${1:-libso}
+
+## get ZOLTAN_ARCH_PATH
+#if settings=`$WM_PROJECT_DIR/bin/foamEtcFile config/zoltan.sh`
+#then
+#    . $settings
+#    echo "using ZOLTAN_ARCH_PATH=$ZOLTAN_ARCH_PATH"
+#else
+#    echo
+#    echo "Error: no config/zoltan.sh settings"
+#    echo
+#fi
+
+
+set -x
+
+wmake $makeType renumberMethods
+
+#if [ -n "$ZOLTAN_ARCH_PATH" ]
+#then
+#    wmake $makeType zoltanRenumber
+#else
+#    echo
+#    echo "Skipping zoltanRenumber"
+#    echo
+#fi
+
+# ----------------------------------------------------------------- end-of-file
diff --git a/src/renumber/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C b/src/renumber/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
new file mode 100644
index 0000000000..2732484b36
--- /dev/null
+++ b/src/renumber/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
@@ -0,0 +1,106 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "CuthillMcKeeRenumber.H"
+#include "addToRunTimeSelectionTable.H"
+#include "bandCompression.H"
+#include "decompositionMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(CuthillMcKeeRenumber, 0);
+
+    addToRunTimeSelectionTable
+    (
+        renumberMethod,
+        CuthillMcKeeRenumber,
+        dictionary
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::CuthillMcKeeRenumber::CuthillMcKeeRenumber(const dictionary& renumberDict)
+:
+    renumberMethod(renumberDict),
+    reverse_
+    (
+        renumberDict.found(typeName + "Coeffs")
+      ? Switch(renumberDict.subDict(typeName + "Coeffs").lookup("reverse"))
+      : Switch(false)
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::labelList Foam::CuthillMcKeeRenumber::renumber
+(
+    const polyMesh& mesh,
+    const pointField& points
+) const
+{
+    CompactListList<label> cellCells;
+    decompositionMethod::calcCellCells
+    (
+        mesh,
+        identity(mesh.nCells()),
+        mesh.nCells(),
+        false,                      // local only
+        cellCells
+    );
+
+    labelList orderedToOld = bandCompression(cellCells());
+
+    if (reverse_)
+    {
+        reverse(orderedToOld);
+    }
+
+    return invert(orderedToOld.size(), orderedToOld);
+}
+
+
+Foam::labelList Foam::CuthillMcKeeRenumber::renumber
+(
+    const labelListList& cellCells,
+    const pointField& points
+) const
+{
+    labelList orderedToOld = bandCompression(cellCells);
+
+    if (reverse_)
+    {
+        reverse(orderedToOld);
+    }
+
+    return invert(orderedToOld.size(), orderedToOld);
+}
+
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H b/src/renumber/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
new file mode 100644
index 0000000000..a4ec011b40
--- /dev/null
+++ b/src/renumber/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
@@ -0,0 +1,118 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::CuthillMcKeeRenumber
+
+Description
+    Cuthill-McKee renumbering
+
+SourceFiles
+    CuthillMcKeeRenumber.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef CuthillMcKeeRenumber_H
+#define CuthillMcKeeRenumber_H
+
+#include "renumberMethod.H"
+#include "Switch.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                   Class CuthillMcKeeRenumber Declaration
+\*---------------------------------------------------------------------------*/
+
+class CuthillMcKeeRenumber
+:
+    public renumberMethod
+{
+    // Private data
+
+        const Switch reverse_;
+
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct and assignment
+        void operator=(const CuthillMcKeeRenumber&);
+        CuthillMcKeeRenumber(const CuthillMcKeeRenumber&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("CuthillMcKee");
+
+
+    // Constructors
+
+        //- Construct given the renumber dictionary
+        CuthillMcKeeRenumber(const dictionary& renumberDict);
+
+
+    //- Destructor
+    virtual ~CuthillMcKeeRenumber()
+    {}
+
+
+    // Member Functions
+
+        //- Return for every coordinate the wanted processor number.
+        //  We need a polyMesh (to be able to load the file)
+        virtual labelList renumber(const pointField&) const
+        {
+            notImplemented("CuthillMcKeeRenumber::renumber(const pointField&)");
+            return labelList(0);
+        }
+
+        //- Return for every coordinate the wanted processor number. Use the
+        //  mesh connectivity (if needed)
+        virtual labelList renumber
+        (
+            const polyMesh& mesh,
+            const pointField& cc
+        ) const;
+
+        //- Return for every cell the new cell label.
+        //  The connectivity is equal to mesh.cellCells() except
+        //  - the connections are across coupled patches
+        virtual labelList renumber
+        (
+            const labelListList& cellCells,
+            const pointField& cc
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/Make/files b/src/renumber/renumberMethods/Make/files
new file mode 100644
index 0000000000..c19f1e6b6d
--- /dev/null
+++ b/src/renumber/renumberMethods/Make/files
@@ -0,0 +1,8 @@
+renumberMethod/renumberMethod.C
+manualRenumber/manualRenumber.C
+CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
+randomRenumber/randomRenumber.C
+springRenumber/springRenumber.C
+boundaryFirstRenumber/boundaryFirstRenumber.C
+
+LIB = $(FOAM_LIBBIN)/librenumberMethods
diff --git a/src/renumber/renumberMethods/Make/options b/src/renumber/renumberMethods/Make/options
new file mode 100644
index 0000000000..03cb68d946
--- /dev/null
+++ b/src/renumber/renumberMethods/Make/options
@@ -0,0 +1,7 @@
+EXE_INC = \
+    -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
+    -I$(LIB_SRC)/meshTools/lnInclude
+
+LIB_LIBS = \
+    -ldecompositionMethods \
+    -lmeshTools
diff --git a/src/renumber/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C b/src/renumber/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
new file mode 100644
index 0000000000..293e5ba19b
--- /dev/null
+++ b/src/renumber/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
@@ -0,0 +1,201 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "boundaryFirstRenumber.H"
+#include "addToRunTimeSelectionTable.H"
+#include "bandCompression.H"
+#include "decompositionMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(boundaryFirstRenumber, 0);
+
+    addToRunTimeSelectionTable
+    (
+        renumberMethod,
+        boundaryFirstRenumber,
+        dictionary
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::boundaryFirstRenumber::boundaryFirstRenumber
+(
+    const dictionary& renumberDict
+)
+:
+    renumberMethod(renumberDict),
+    reverse_
+    (
+        renumberDict.found(typeName + "Coeffs")
+      ? Switch(renumberDict.subDict(typeName + "Coeffs").lookup("reverse"))
+      : Switch(true)
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::labelList Foam::boundaryFirstRenumber::renumber
+(
+    const polyMesh& mesh,
+    const pointField& points
+) const
+{
+    const labelList& own = mesh.faceOwner();
+    const labelList& nei = mesh.faceNeighbour();
+
+    // Distance to boundary. Minimum of neighbours.
+    labelList distance(mesh.nCells(), -1);
+
+    // New order
+    labelList newOrder(mesh.nCells());
+    label cellInOrder = 0;
+
+
+    // Starting faces for walk. These are the zero distance faces.
+    DynamicList<label> frontFaces(mesh.nFaces()-mesh.nInternalFaces());
+
+    const polyBoundaryMesh& patches = mesh.boundaryMesh();
+
+    forAll(patches, patchI)
+    {
+        const polyPatch& pp = patches[patchI];
+
+        //- Note: cannot check for empty since these are introduced by
+        //  fvMeshSubset. Also most loops don't care about patch type.
+        //if (!isA<emptyPolyPatch>(pp))
+        {
+            forAll(pp, i)
+            {
+                frontFaces.append(pp.start()+i);
+            }
+        }
+    }
+
+    // Loop over all frontFaces
+    label currentDistance = 0;
+    while (frontFaces.size() > 0)
+    {
+        DynamicList<label> frontCells(frontFaces.size());
+
+        // Set all of frontFaces' neighbours to current distance
+
+        forAll(frontFaces, i)
+        {
+            label faceI = frontFaces[i];
+
+            if (mesh.isInternalFace(faceI))
+            {
+                label ownCellI = own[faceI];
+                if (distance[ownCellI] == -1)
+                {
+                    distance[ownCellI] = currentDistance;
+                    frontCells.append(ownCellI);
+                }
+                label neiCellI = nei[faceI];
+                if (distance[neiCellI] == -1)
+                {
+                    distance[neiCellI] = currentDistance;
+                    frontCells.append(neiCellI);
+                }
+            }
+            else
+            {
+                label ownCellI = own[faceI];
+                if (distance[ownCellI] == -1)
+                {
+                    distance[ownCellI] = currentDistance;
+                    frontCells.append(ownCellI);
+                }
+            }
+        }
+
+        //Pout<< "For distance:" << currentDistance
+        //    << " from " << frontFaces.size() << " faces to "
+        //    << frontCells.size() << " cells." << endl;
+
+
+        // TBD. Determine order within current shell (frontCells). For now
+        // just add them.
+        forAll(frontCells, i)
+        {
+            newOrder[cellInOrder] = frontCells[i];
+            cellInOrder++;
+        }
+
+        // From cells to faces
+        frontFaces.clear();
+        forAll(frontCells, i)
+        {
+            label cellI = frontCells[i];
+            const cell& cFaces = mesh.cells()[cellI];
+
+            forAll(cFaces, i)
+            {
+                label faceI = cFaces[i];
+                if (mesh.isInternalFace(faceI))
+                {
+                    label nbrCellI =
+                    (
+                        mesh.faceOwner()[faceI] == cellI
+                      ? mesh.faceNeighbour()[faceI]
+                      : mesh.faceOwner()[faceI]
+                    );
+                    if (distance[nbrCellI] == -1)
+                    {
+                        frontFaces.append(faceI);
+                    }
+                }
+            }
+        }
+
+        currentDistance++;
+    }
+
+    // Return furthest away cell first
+    if (reverse_)
+    {
+        reverse(newOrder);
+    }
+
+    //forAll(newOrder, i)
+    //{
+    //    label cellI = newOrder[i];
+    //
+    //    Pout<< "cell:" << cellI << endl;
+    //    Pout<< " at distance:" << distance[cellI]
+    //        << endl;
+    //}
+
+    return invert(newOrder.size(), newOrder);
+}
+
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H b/src/renumber/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H
new file mode 100644
index 0000000000..690af21455
--- /dev/null
+++ b/src/renumber/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H
@@ -0,0 +1,129 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::boundaryFirstRenumber
+
+Description
+    Cuthill-McKee renumbering
+
+SourceFiles
+    boundaryFirstRenumber.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef boundaryFirstRenumber_H
+#define boundaryFirstRenumber_H
+
+#include "renumberMethod.H"
+#include "Switch.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                   Class boundaryFirstRenumber Declaration
+\*---------------------------------------------------------------------------*/
+
+class boundaryFirstRenumber
+:
+    public renumberMethod
+{
+    // Private data
+
+        const Switch reverse_;
+
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct and assignment
+        void operator=(const boundaryFirstRenumber&);
+        boundaryFirstRenumber(const boundaryFirstRenumber&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("boundaryFirst");
+
+
+    // Constructors
+
+        //- Construct given the renumber dictionary
+        boundaryFirstRenumber(const dictionary& renumberDict);
+
+
+    //- Destructor
+    virtual ~boundaryFirstRenumber()
+    {}
+
+
+    // Member Functions
+
+        //- Return for every coordinate the wanted processor number.
+        //  We need a polyMesh (to be able to load the file)
+        virtual labelList renumber(const pointField&) const
+        {
+            notImplemented
+            (
+                "boundaryFirstRenumber::renumber(const pointField&)"
+            );
+            return labelList(0);
+        }
+
+        //- Return for every coordinate the wanted processor number. Use the
+        //  mesh connectivity (if needed)
+        virtual labelList renumber
+        (
+            const polyMesh& mesh,
+            const pointField& cc
+        ) const;
+
+        //- Return for every cell the new cell label.
+        //  The connectivity is equal to mesh.cellCells() except
+        //  - the connections are across coupled patches
+        virtual labelList renumber
+        (
+            const labelListList& cellCells,
+            const pointField& cc
+        ) const
+        {
+            notImplemented
+            (
+                "boundaryFirstRenumber::renumber"
+                "(const labelListList&, const pointField&)"
+            );
+            return labelList(0);
+        }
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/manualRenumber/manualRenumber.C b/src/renumber/renumberMethods/manualRenumber/manualRenumber.C
new file mode 100644
index 0000000000..92b5c8e00c
--- /dev/null
+++ b/src/renumber/renumberMethods/manualRenumber/manualRenumber.C
@@ -0,0 +1,137 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "manualRenumber.H"
+#include "addToRunTimeSelectionTable.H"
+#include "IFstream.H"
+#include "labelIOList.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(manualRenumber, 0);
+
+    addToRunTimeSelectionTable
+    (
+        renumberMethod,
+        manualRenumber,
+        dictionary
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::manualRenumber::manualRenumber(const dictionary& renumberDict)
+:
+    renumberMethod(renumberDict),
+    dataFile_
+    (
+        renumberDict.subDict(typeName+"Coeffs").lookup("dataFile")
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::labelList Foam::manualRenumber::renumber
+(
+    const polyMesh& mesh,
+    const pointField& points
+) const
+{
+    labelIOList oldToNew
+    (
+        IOobject
+        (
+            dataFile_,
+            mesh.facesInstance(),
+            mesh,
+            IOobject::MUST_READ,
+            IOobject::AUTO_WRITE,
+            false
+        )
+    );
+
+    // check if the final renumbering is OK
+
+    if (oldToNew.size() != points.size())
+    {
+        FatalErrorIn
+        (
+            "manualRenumber::renumber(const pointField&, const scalarField&)"
+        )   << "Size of renumber list does not correspond "
+            << "to the number of points.  Size: "
+            << oldToNew.size() << " Number of points: "
+            << points.size()
+            << ".\n" << "Manual renumbering data read from file "
+            << dataFile_ << "." << endl
+            << exit(FatalError);
+    }
+
+    // Invert to see if one to one
+    labelList newToOld(points.size(), -1);
+    forAll(oldToNew, i)
+    {
+        label newI = oldToNew[i];
+
+        if (newI < 0 || newI >= oldToNew.size())
+        {
+            FatalErrorIn
+            (
+                "manualRenumber::renumber(const pointField&"
+                ", const scalarField&)"
+            )   << "Renumbering is not one-to-one. Index "
+                << i << " maps onto " << newI
+                << ".\n" << "Manual renumbering data read from file "
+                << dataFile_ << "." << endl
+                << exit(FatalError);
+        }
+
+        if (newToOld[newI] == -1)
+        {
+            newToOld[newI] = i;
+        }
+        else
+        {
+            FatalErrorIn
+            (
+                "manualRenumber::renumber(const pointField&"
+                ", const scalarField&)"
+            )   << "Renumbering is not one-to-one. Both index "
+                << newToOld[newI]
+                << " and " << i << " map onto " << newI
+                << ".\n" << "Manual renumbering data read from file "
+                << dataFile_ << "." << endl
+                << exit(FatalError);
+        }
+    }
+
+    return oldToNew;
+}
+
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/manualRenumber/manualRenumber.H b/src/renumber/renumberMethods/manualRenumber/manualRenumber.H
new file mode 100644
index 0000000000..a4d4b2f46f
--- /dev/null
+++ b/src/renumber/renumberMethods/manualRenumber/manualRenumber.H
@@ -0,0 +1,125 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::manualRenumber
+
+Description
+    Renumber given a cell-to-new cell association in a file
+
+SourceFiles
+    manualRenumber.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef manualRenumber_H
+#define manualRenumber_H
+
+#include "renumberMethod.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                      Class manualRenumber Declaration
+\*---------------------------------------------------------------------------*/
+
+class manualRenumber
+:
+    public renumberMethod
+{
+    // Private data
+
+        const fileName dataFile_;
+
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct and assignment
+        void operator=(const manualRenumber&);
+        manualRenumber(const manualRenumber&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("manual");
+
+
+    // Constructors
+
+        //- Construct given the renumber dictionary
+        manualRenumber(const dictionary& renumberDict);
+
+
+    //- Destructor
+    virtual ~manualRenumber()
+    {}
+
+
+    // Member Functions
+
+        //- Return for every coordinate the wanted processor number.
+        //  We need a polyMesh (to be able to load the file)
+        virtual labelList renumber(const pointField&) const
+        {
+            notImplemented("manualRenumber::renumber(const pointField&)");
+            return labelList(0);
+        }
+
+        //- Return for every coordinate the wanted processor number. Use the
+        //  mesh connectivity (if needed)
+        virtual labelList renumber
+        (
+            const polyMesh& mesh,
+            const pointField& cc
+        ) const;
+
+        //- Return for every cell the new cell label.
+        //  The connectivity is equal to mesh.cellCells() except
+        //  - the connections are across coupled patches
+        virtual labelList renumber
+        (
+            const labelListList& cellCells,
+            const pointField& cc
+        ) const
+        {
+            notImplemented
+            (
+                "manualRenumber::renumber"
+                "(const labelListList&, const pointField&)"
+            );
+            return labelList(0);
+        }
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/randomRenumber/randomRenumber.C b/src/renumber/renumberMethods/randomRenumber/randomRenumber.C
new file mode 100644
index 0000000000..2dae28e438
--- /dev/null
+++ b/src/renumber/renumberMethods/randomRenumber/randomRenumber.C
@@ -0,0 +1,96 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "randomRenumber.H"
+#include "addToRunTimeSelectionTable.H"
+#include "Random.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(randomRenumber, 0);
+
+    addToRunTimeSelectionTable
+    (
+        renumberMethod,
+        randomRenumber,
+        dictionary
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::randomRenumber::randomRenumber(const dictionary& renumberDict)
+:
+    renumberMethod(renumberDict)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::labelList Foam::randomRenumber::renumber
+(
+    const pointField& points
+) const
+{
+    Random rndGen(0);
+
+    labelList oldToNew(identity(points.size()));
+
+    for (label iter = 0; iter < 10; iter++)
+    {
+        forAll(oldToNew, i)
+        {
+            label j = rndGen.integer(0, oldToNew.size()-1);
+            Swap(oldToNew[i], oldToNew[j]);
+        }
+    }
+    return oldToNew;
+}
+
+
+Foam::labelList Foam::randomRenumber::renumber
+(
+    const polyMesh& mesh,
+    const pointField& points
+) const
+{
+    return renumber(points);
+}
+
+
+Foam::labelList Foam::randomRenumber::renumber
+(
+    const labelListList& cellCells,
+    const pointField& points
+) const
+{
+    return renumber(points);
+}
+
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/randomRenumber/randomRenumber.H b/src/renumber/renumberMethods/randomRenumber/randomRenumber.H
new file mode 100644
index 0000000000..0665d5159b
--- /dev/null
+++ b/src/renumber/renumberMethods/randomRenumber/randomRenumber.H
@@ -0,0 +1,108 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::randomRenumber
+
+Description
+    Random renumber. Just to see effect of renumbering.
+
+SourceFiles
+    randomRenumber.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef randomRenumber_H
+#define randomRenumber_H
+
+#include "renumberMethod.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                      Class randomRenumber Declaration
+\*---------------------------------------------------------------------------*/
+
+class randomRenumber
+:
+    public renumberMethod
+{
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct and assignment
+        void operator=(const randomRenumber&);
+        randomRenumber(const randomRenumber&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("random");
+
+
+    // Constructors
+
+        //- Construct given the renumber dictionary
+        randomRenumber(const dictionary& renumberDict);
+
+
+    //- Destructor
+    virtual ~randomRenumber()
+    {}
+
+
+    // Member Functions
+
+        //- Return for every coordinate the wanted processor number.
+        //  We need a polyMesh (to be able to load the file)
+        virtual labelList renumber(const pointField&) const;
+
+        //- Return for every coordinate the wanted processor number. Use the
+        //  mesh connectivity (if needed)
+        virtual labelList renumber
+        (
+            const polyMesh& mesh,
+            const pointField& cc
+        ) const;
+
+        //- Return for every cell the new cell label.
+        //  The connectivity is equal to mesh.cellCells() except
+        //  - the connections are across coupled patches
+        virtual labelList renumber
+        (
+            const labelListList& cellCells,
+            const pointField& cc
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/renumberMethod/renumberMethod.C b/src/renumber/renumberMethods/renumberMethod/renumberMethod.C
new file mode 100644
index 0000000000..76ba1911e2
--- /dev/null
+++ b/src/renumber/renumberMethods/renumberMethod/renumberMethod.C
@@ -0,0 +1,131 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+InClass
+    renumberMethod
+
+\*---------------------------------------------------------------------------*/
+
+#include "renumberMethod.H"
+#include "decompositionMethod.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(renumberMethod, 0);
+    defineRunTimeSelectionTable(renumberMethod, dictionary);
+}
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::renumberMethod> Foam::renumberMethod::New
+(
+    const dictionary& renumberDict
+)
+{
+    const word methodType(renumberDict.lookup("method"));
+
+    //Info<< "Selecting renumberMethod " << methodType << endl;
+
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(methodType);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalErrorIn
+        (
+            "renumberMethod::New"
+            "(const dictionary& renumberDict)"
+        )   << "Unknown renumberMethod "
+            << methodType << nl << nl
+            << "Valid renumberMethods are : " << endl
+            << dictionaryConstructorTablePtr_->sortedToc()
+            << exit(FatalError);
+    }
+
+    return autoPtr<renumberMethod>(cstrIter()(renumberDict));
+}
+
+
+Foam::labelList Foam::renumberMethod::renumber
+(
+    const polyMesh& mesh,
+    const pointField& points
+) const
+{
+    CompactListList<label> cellCells;
+    decompositionMethod::calcCellCells
+    (
+        mesh,
+        identity(mesh.nCells()),
+        mesh.nCells(),
+        false,                      // local only
+        cellCells
+    );
+
+    // Renumber based on agglomerated points
+    return renumber(cellCells(), points);
+}
+
+
+Foam::labelList Foam::renumberMethod::renumber
+(
+    const polyMesh& mesh,
+    const labelList& fineToCoarse,
+    const pointField& coarsePoints
+) const
+{
+    CompactListList<label> coarseCellCells;
+    decompositionMethod::calcCellCells
+    (
+        mesh,
+        fineToCoarse,
+        coarsePoints.size(),
+        false,                      // local only
+        coarseCellCells
+    );
+
+    // Renumber based on agglomerated points
+    labelList coarseDistribution
+    (
+        renumber
+        (
+            coarseCellCells(),
+            coarsePoints
+        )
+    );
+
+    // Rework back into renumbering for original mesh_
+    labelList fineDistribution(fineToCoarse.size());
+
+    forAll(fineDistribution, i)
+    {
+        fineDistribution[i] = coarseDistribution[fineToCoarse[i]];
+    }
+
+    return fineDistribution;
+}
+
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/renumberMethod/renumberMethod.H b/src/renumber/renumberMethods/renumberMethod/renumberMethod.H
new file mode 100644
index 0000000000..ead92b1e41
--- /dev/null
+++ b/src/renumber/renumberMethods/renumberMethod/renumberMethod.H
@@ -0,0 +1,161 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::renumberMethod
+
+Description
+    Abstract base class for renumbering
+
+SourceFiles
+    renumberMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef renumberMethod_H
+#define renumberMethod_H
+
+#include "polyMesh.H"
+#include "pointField.H"
+#include "CompactListList.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class renumberMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class renumberMethod
+{
+
+protected:
+
+    // Protected data
+
+        const dictionary& renumberDict_;
+
+private:
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct and assignment
+        renumberMethod(const renumberMethod&);
+        void operator=(const renumberMethod&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("renumberMethod");
+
+
+    // Declare run-time constructor selection tables
+
+        declareRunTimeSelectionTable
+        (
+            autoPtr,
+            renumberMethod,
+            dictionary,
+            (
+                const dictionary& renumberDict
+            ),
+            (renumberDict)
+        );
+
+
+    // Selectors
+
+        //- Return a reference to the selected renumbering method
+        static autoPtr<renumberMethod> New
+        (
+            const dictionary& renumberDict
+        );
+
+
+    // Constructors
+
+        //- Construct given the renumber dictionary
+        renumberMethod(const dictionary& renumberDict)
+        :
+            renumberDict_(renumberDict)
+        {}
+
+
+    //- Destructor
+    virtual ~renumberMethod()
+    {}
+
+
+    // Member Functions
+
+        //- Return for every cell the new cell label.
+        //  This is only defined for geometric renumberMethods.
+        virtual labelList renumber(const pointField&) const
+        {
+            notImplemented
+            (
+                "renumberMethod:renumber(const pointField&)"
+            );
+            return labelList(0);
+        }
+
+        //- Return for every cell the new cell label. Use the
+        //  mesh connectivity (if needed)
+        virtual labelList renumber(const polyMesh&, const pointField&) const;
+
+        //- Return for every cell the new cell label. Gets
+        //  passed agglomeration map (from fine to coarse cells) and coarse
+        //  cell
+        //  location. Can be overridden by renumberMethods that provide this
+        //  functionality natively. Coarse cells are local to the processor
+        //  (if in parallel). If you want to have coarse cells spanning
+        //  processors use the globalCellCells instead.
+        virtual labelList renumber
+        (
+            const polyMesh& mesh,
+            const labelList& cellToRegion,
+            const pointField& regionPoints
+        ) const;
+
+        //- Return for every cell the new cell label.
+        //  The connectivity is equal to mesh.cellCells() except
+        //  - the connections are across coupled patches
+        virtual labelList renumber
+        (
+            const labelListList& cellCells,
+            const pointField& cc
+        ) const = 0;
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/springRenumber/springRenumber.C b/src/renumber/renumberMethods/springRenumber/springRenumber.C
new file mode 100644
index 0000000000..64a15145ce
--- /dev/null
+++ b/src/renumber/renumberMethods/springRenumber/springRenumber.C
@@ -0,0 +1,172 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "springRenumber.H"
+#include "addToRunTimeSelectionTable.H"
+#include "decompositionMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(springRenumber, 0);
+
+    addToRunTimeSelectionTable
+    (
+        renumberMethod,
+        springRenumber,
+        dictionary
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::springRenumber::springRenumber(const dictionary& renumberDict)
+:
+    renumberMethod(renumberDict),
+    dict_(renumberDict.subDict(typeName+"Coeffs")),
+    maxCo_(readScalar(dict_.lookup("maxCo"))),
+    maxIter_(readLabel(dict_.lookup("maxIter"))),
+    freezeFraction_(readScalar(dict_.lookup("freezeFraction")))
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::labelList Foam::springRenumber::renumber
+(
+    const polyMesh& mesh,
+    const pointField& points
+) const
+{
+    CompactListList<label> cellCells;
+    decompositionMethod::calcCellCells
+    (
+        mesh,
+        identity(mesh.nCells()),
+        mesh.nCells(),
+        false,                      // local only
+        cellCells
+    );
+
+    return renumber(cellCells(), points);
+}
+
+
+Foam::labelList Foam::springRenumber::renumber
+(
+    const labelListList& cellCells,
+    const pointField& points
+) const
+{
+    // Look at cell index as a 1D position parameter.
+    // Move cells to the average 'position' of their neighbour.
+
+    scalarField position(cellCells.size());
+    forAll(position, cellI)
+    {
+        position[cellI] = cellI;
+    }
+
+    labelList oldToNew(identity(cellCells.size()));
+
+    scalar maxCo = maxCo_ * cellCells.size();
+
+    for (label iter = 0; iter < maxIter_; iter++)
+    {
+        //Pout<< "Iteration : " << iter << nl
+        //    << "------------"
+        //    << endl;
+
+        //Pout<< "Position :" << nl
+        //    << "    min : " << min(position) << nl
+        //    << "    max : " << max(position) << nl
+        //    << "    avg : " << average(position) << nl
+        //    << endl;
+
+        // Sum force per cell.
+        scalarField sumForce(cellCells.size(), 0.0);
+        forAll(cellCells, oldCellI)
+        {
+            const labelList& cCells = cellCells[oldCellI];
+            label cellI = oldToNew[oldCellI];
+
+            forAll(cCells, i)
+            {
+                label nbrCellI = oldToNew[cCells[i]];
+
+                sumForce[cellI] += (position[nbrCellI]-position[cellI]);
+            }
+        }
+
+        //Pout<< "Force :" << nl
+        //    << "    min    : " << min(sumForce) << nl
+        //    << "    max    : " << max(sumForce) << nl
+        //    << "    avgMag : " << average(mag(sumForce)) << nl
+        //    << "DeltaT : " << deltaT << nl
+        //    << endl;
+
+        // Limit displacement
+        scalar deltaT = maxCo / max(mag(sumForce));
+
+        Info<< "Iter:" << iter
+            << "  maxCo:" << maxCo
+            << "  deltaT:" << deltaT
+            << "  average force:" << average(mag(sumForce)) << endl;
+
+        // Determine displacement.
+        scalarField displacement = deltaT*sumForce;
+
+        //Pout<< "Displacement :" << nl
+        //    << "    min    : " << min(displacement) << nl
+        //    << "    max    : " << max(displacement) << nl
+        //    << "    avgMag : " << average(mag(displacement)) << nl
+        //    << endl;
+
+        // Calculate new position and scale to be within original range
+        // (0..nCells-1) for ease of postprocessing.
+        position += displacement;
+        position -= min(position);
+        position *= (position.size()-1)/max(position);
+
+        // Slowly freeze.
+        maxCo *= freezeFraction_;
+    }
+
+    //writeOBJ("endPosition.obj", cellCells, position);
+
+    // Move cells to new position
+    labelList shuffle;
+    sortedOrder(position, shuffle);
+
+    // Reorder oldToNew
+    inplaceReorder(shuffle, oldToNew);
+
+    return oldToNew;
+}
+
+
+// ************************************************************************* //
diff --git a/src/renumber/renumberMethods/springRenumber/springRenumber.H b/src/renumber/renumberMethods/springRenumber/springRenumber.H
new file mode 100644
index 0000000000..30f7a3d838
--- /dev/null
+++ b/src/renumber/renumberMethods/springRenumber/springRenumber.H
@@ -0,0 +1,132 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::springRenumber
+
+Description
+    Use spring analogy - attract neighbouring cells according to the distance
+    of their cell indices.
+
+    // Maximum jump of cell indices. Is fraction of number of cells
+    maxCo 0.1;
+
+    // Limit the amount of movement; the fraction maxCo gets decreased
+    // with every iteration.
+    freezeFraction 0.9;
+
+    // Maximum number of iterations
+    maxIter 1000;
+
+SourceFiles
+    springRenumber.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef springRenumber_H
+#define springRenumber_H
+
+#include "renumberMethod.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                       Class springRenumber Declaration
+\*---------------------------------------------------------------------------*/
+
+class springRenumber
+:
+    public renumberMethod
+{
+    // Private data
+
+        const dictionary& dict_;
+
+        const scalar maxCo_;
+
+        const label maxIter_;
+
+        const scalar freezeFraction_;
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct and assignment
+        void operator=(const springRenumber&);
+        springRenumber(const springRenumber&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("spring");
+
+
+    // Constructors
+
+        //- Construct given the renumber dictionary
+        springRenumber(const dictionary& renumberDict);
+
+
+    //- Destructor
+    virtual ~springRenumber()
+    {}
+
+
+    // Member Functions
+
+        //- Return for every coordinate the wanted processor number.
+        virtual labelList renumber(const pointField&) const
+        {
+            notImplemented("springRenumber::renumber(const pointField&)");
+            return labelList(0);
+        }
+
+        //- Return for every coordinate the wanted processor number. Use the
+        //  mesh connectivity (if needed)
+        virtual labelList renumber
+        (
+            const polyMesh& mesh,
+            const pointField& cc
+        ) const;
+
+        //- Return for every cell the new cell label.
+        //  The connectivity is equal to mesh.cellCells() except
+        //  - the connections are across coupled patches
+        virtual labelList renumber
+        (
+            const labelListList& cellCells,
+            const pointField& cc
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

From efadefb7092554809382e8a5b8c39f03ba5b6940 Mon Sep 17 00:00:00 2001
From: andy <andy>
Date: Fri, 23 Dec 2011 09:19:36 +0000
Subject: [PATCH 12/18] BUG: Pimple corrector loop counter not re-initialised
 on criteria satisfied - mantis #371

---
 .../general/solutionControl/pimpleControl/pimpleControl.C      | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.C b/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.C
index 98a742e16b..7ff31db0ec 100644
--- a/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.C
+++ b/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.C
@@ -196,9 +196,10 @@ bool Foam::pimpleControl::loop()
     bool completed = false;
     if (criteriaSatisfied())
     {
-        Info<< algorithmName_ << ": converged in " << corr_ << " iterations"
+        Info<< algorithmName_ << ": converged in " << corr_ - 1 << " iterations"
             << endl;
         completed = true;
+        corr_ = 0;
     }
     else
     {

From d7d0a843926b8d07399fd535f90b22fddb55642b Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 23 Dec 2011 11:38:35 +0000
Subject: [PATCH 13/18] BUG: meshRefinementProblemCells: conversion of angle
 into cos

---
 .../meshRefinementProblemCells.C              | 26 ++++++++++++-------
 1 file changed, 16 insertions(+), 10 deletions(-)

diff --git a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementProblemCells.C b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementProblemCells.C
index 296029b7d5..48eff18670 100644
--- a/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementProblemCells.C
+++ b/src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementProblemCells.C
@@ -286,6 +286,10 @@ bool Foam::meshRefinement::isCollapsedFace
     const label faceI
 ) const
 {
+    // Severe nonorthogonality threshold
+    const scalar severeNonorthogonalityThreshold =
+        ::cos(degToRad(maxNonOrtho));
+
     vector s = mesh_.faces()[faceI].normal(points);
     scalar magS = mag(s);
 
@@ -301,11 +305,11 @@ bool Foam::meshRefinement::isCollapsedFace
     if (mesh_.isInternalFace(faceI))
     {
         label nei = mesh_.faceNeighbour()[faceI];
-        vector d = ownCc - mesh_.cellCentres()[nei];
+        vector d = mesh_.cellCentres()[nei] - ownCc;
 
         scalar dDotS = (d & s)/(mag(d)*magS + VSMALL);
 
-        if (dDotS < maxNonOrtho)
+        if (dDotS < severeNonorthogonalityThreshold)
         {
             return true;
         }
@@ -320,11 +324,11 @@ bool Foam::meshRefinement::isCollapsedFace
 
         if (mesh_.boundaryMesh()[patchI].coupled())
         {
-            vector d = ownCc - neiCc[faceI-mesh_.nInternalFaces()];
+            vector d = neiCc[faceI-mesh_.nInternalFaces()] - ownCc;
 
             scalar dDotS = (d & s)/(mag(d)*magS + VSMALL);
 
-            if (dDotS < maxNonOrtho)
+            if (dDotS < severeNonorthogonalityThreshold)
             {
                 return true;
             }
@@ -337,7 +341,7 @@ bool Foam::meshRefinement::isCollapsedFace
         {
             // Collapsing normal boundary face does not cause problems with
             // non-orthogonality
-            return true;
+            return false;
         }
     }
 }
@@ -663,7 +667,8 @@ Foam::labelList Foam::meshRefinement::markFacesOnProblemCells
                 )
                 {
                     nPrevented++;
-                    //Pout<< "Preventing collapse of 8 anchor point cell "
+                    //Pout<< "Preventing baffling/removal of 8 anchor point"
+                    //    << " cell "
                     //    << cellI << " at " << mesh_.cellCentres()[cellI]
                     //    << " since new volume "
                     //    << mesh_.cells()[cellI].mag(newPoints, mesh_.faces())
@@ -747,7 +752,7 @@ Foam::labelList Foam::meshRefinement::markFacesOnProblemCells
                     )
                     {
                         nPrevented++;
-                        //Pout<< "Preventing collapse of 7 anchor cell "
+                        //Pout<< "Preventing baffling of 7 anchor cell "
                         //    << cellI
                         //    << " at " << mesh_.cellCentres()[cellI]
                         //    << " since new volume "
@@ -853,7 +858,8 @@ Foam::labelList Foam::meshRefinement::markFacesOnProblemCells
                 )
                 {
                     nPrevented++;
-                    //Pout<< "Preventing collapse of face " << faceI
+                    //Pout<< "Preventing baffling (to avoid collapse) of face "
+                    //    << faceI
                     //    << " with all boundary edges "
                     //    << " at " << mesh_.faceCentres()[faceI]
                     //    << endl;
@@ -909,7 +915,7 @@ Foam::labelList Foam::meshRefinement::markFacesOnProblemCells
                         )
                         {
                             nPrevented++;
-                            //Pout<< "Preventing collapse of coupled face "
+                            //Pout<< "Preventing baffling of coupled face "
                             //    << faceI
                             //    << " with all boundary edges "
                             //    << " at " << mesh_.faceCentres()[faceI]
@@ -943,7 +949,7 @@ Foam::labelList Foam::meshRefinement::markFacesOnProblemCells
     {
         Info<< "markFacesOnProblemCells : prevented "
             << returnReduce(nPrevented, sumOp<label>())
-            << " internal faces fom getting converted into baffles."
+            << " internal faces from getting converted into baffles."
             << endl;
     }
 

From e4fe63a696a525aebcc90644951978c0a52381d1 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 23 Dec 2011 12:11:37 +0000
Subject: [PATCH 14/18] BUG: decompositionMethod: numbering across coupled
 boundaries

---
 .../decompositionMethod/decompositionMethod.C | 19 ++++++++-------
 .../decompose/ptscotchDecomp/ptscotchDecomp.C | 23 -------------------
 2 files changed, 9 insertions(+), 33 deletions(-)

diff --git a/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C b/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C
index 363000bad6..e0e96412d3 100644
--- a/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C
+++ b/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C
@@ -236,8 +236,8 @@ void Foam::decompositionMethod::calcCellCells
                 label globalNei = globalNeighbour[bFaceI];
                 if
                 (
-                    globalAgglom.isLocal(globalNei)
-                 && globalAgglom.toLocal(globalNei) != own
+                   !globalAgglom.isLocal(globalNei)
+                 || globalAgglom.toLocal(globalNei) != own
                 )
                 {
                     nFacesPerCell[own]++;
@@ -285,10 +285,11 @@ void Foam::decompositionMethod::calcCellCells
                 label own = agglom[faceOwner[faceI]];
 
                 label globalNei = globalNeighbour[bFaceI];
+
                 if
                 (
-                    globalAgglom.isLocal(globalNei)
-                 && globalAgglom.toLocal(globalNei) != own
+                   !globalAgglom.isLocal(globalNei)
+                 || globalAgglom.toLocal(globalNei) != own
                 )
                 {
                     m[offsets[own] + nFacesPerCell[own]++] = globalNei;
@@ -304,15 +305,15 @@ void Foam::decompositionMethod::calcCellCells
     // Check for duplicates connections between cells
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     // Done as postprocessing step since we now have cellCells.
-    label startIndex = 0;
     label newIndex = 0;
     labelHashSet nbrCells;
     forAll(cellCells, cellI)
     {
         nbrCells.clear();
-        nbrCells.insert(cellI);
+        nbrCells.insert(globalAgglom.toGlobal(cellI));
 
-        label& endIndex = cellCells.offsets()[cellI+1];
+        label startIndex = cellCells.offsets()[cellI];
+        label endIndex = cellCells.offsets()[cellI+1];
 
         for (label i = startIndex; i < endIndex; i++)
         {
@@ -321,9 +322,7 @@ void Foam::decompositionMethod::calcCellCells
                 cellCells.m()[newIndex++] = cellCells.m()[i];
             }
         }
-
-        startIndex = endIndex;
-        endIndex = newIndex;
+        cellCells.offsets()[cellI+1] = newIndex;
     }
 
     cellCells.m().setSize(newIndex);
diff --git a/src/parallel/decompose/ptscotchDecomp/ptscotchDecomp.C b/src/parallel/decompose/ptscotchDecomp/ptscotchDecomp.C
index 54be6ad834..e7c4995a36 100644
--- a/src/parallel/decompose/ptscotchDecomp/ptscotchDecomp.C
+++ b/src/parallel/decompose/ptscotchDecomp/ptscotchDecomp.C
@@ -371,29 +371,6 @@ Foam::label Foam::ptscotchDecomp::decompose
         Pout<< "ptscotchDecomp : entering with xadj:" << xadj.size() << endl;
     }
 
-
-if (debug)
-{
-    Pout<< "nProcessors_:" << nProcessors_ << endl;
-
-    globalIndex globalCells(xadj.size()-1);
-
-    Pout<< "Xadj:" << endl;
-    for (label cellI = 0; cellI < xadj.size()-1; cellI++)
-    {
-        Pout<< "cell:" << cellI
-            << "  global:" << globalCells.toGlobal(cellI)
-            << " connected to:" << endl;
-        label start = xadj[cellI];
-        label end = xadj[cellI+1];
-        for (label i = start; i < end; i++)
-        {
-            Pout<< "    cell:" << adjncy[i] << endl;
-        }
-    }
-    Pout<< endl;
-}
-
     // Dump graph
     if (decompositionDict_.found("ptscotchCoeffs"))
     {

From 72520e6c9554b90bac4ee2351f03763858fa2ed2 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 23 Dec 2011 12:12:50 +0000
Subject: [PATCH 15/18] ENH: snappyHexMeshDict: fixed (unused) setting

---
 .../simpleFoam/turbineSiting/system/snappyHexMeshDict           | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/tutorials/incompressible/simpleFoam/turbineSiting/system/snappyHexMeshDict b/tutorials/incompressible/simpleFoam/turbineSiting/system/snappyHexMeshDict
index 1e6f477fc4..8054858f2a 100644
--- a/tutorials/incompressible/simpleFoam/turbineSiting/system/snappyHexMeshDict
+++ b/tutorials/incompressible/simpleFoam/turbineSiting/system/snappyHexMeshDict
@@ -173,7 +173,7 @@ castellatedMeshControls
         }
     }
 
-    resolveFeatureAngle 2;
+    resolveFeatureAngle 30;
 
 
     // Region-wise refinement

From 208cfcfc3163724171728c1dd04dc177385d8e5c Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 23 Dec 2011 12:13:09 +0000
Subject: [PATCH 16/18] GIT: renumberMethods: moved to subdirectory

---
 .../CuthillMcKeeRenumber.C                    | 106 ---------
 .../CuthillMcKeeRenumber.H                    | 118 ----------
 src/renumberMethods/Make/files                |   8 -
 src/renumberMethods/Make/options              |   7 -
 .../boundaryFirstRenumber.C                   | 201 ------------------
 .../boundaryFirstRenumber.H                   | 129 -----------
 .../manualRenumber/manualRenumber.C           | 137 ------------
 .../manualRenumber/manualRenumber.H           | 125 -----------
 .../randomRenumber/randomRenumber.C           |  96 ---------
 .../randomRenumber/randomRenumber.H           | 108 ----------
 .../renumberMethod/renumberMethod.C           | 131 ------------
 .../renumberMethod/renumberMethod.H           | 161 --------------
 .../springRenumber/springRenumber.C           | 172 ---------------
 .../springRenumber/springRenumber.H           | 132 ------------
 14 files changed, 1631 deletions(-)
 delete mode 100644 src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
 delete mode 100644 src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
 delete mode 100644 src/renumberMethods/Make/files
 delete mode 100644 src/renumberMethods/Make/options
 delete mode 100644 src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
 delete mode 100644 src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H
 delete mode 100644 src/renumberMethods/manualRenumber/manualRenumber.C
 delete mode 100644 src/renumberMethods/manualRenumber/manualRenumber.H
 delete mode 100644 src/renumberMethods/randomRenumber/randomRenumber.C
 delete mode 100644 src/renumberMethods/randomRenumber/randomRenumber.H
 delete mode 100644 src/renumberMethods/renumberMethod/renumberMethod.C
 delete mode 100644 src/renumberMethods/renumberMethod/renumberMethod.H
 delete mode 100644 src/renumberMethods/springRenumber/springRenumber.C
 delete mode 100644 src/renumberMethods/springRenumber/springRenumber.H

diff --git a/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C b/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
deleted file mode 100644
index 2732484b36..0000000000
--- a/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
+++ /dev/null
@@ -1,106 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "CuthillMcKeeRenumber.H"
-#include "addToRunTimeSelectionTable.H"
-#include "bandCompression.H"
-#include "decompositionMethod.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(CuthillMcKeeRenumber, 0);
-
-    addToRunTimeSelectionTable
-    (
-        renumberMethod,
-        CuthillMcKeeRenumber,
-        dictionary
-    );
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::CuthillMcKeeRenumber::CuthillMcKeeRenumber(const dictionary& renumberDict)
-:
-    renumberMethod(renumberDict),
-    reverse_
-    (
-        renumberDict.found(typeName + "Coeffs")
-      ? Switch(renumberDict.subDict(typeName + "Coeffs").lookup("reverse"))
-      : Switch(false)
-    )
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::labelList Foam::CuthillMcKeeRenumber::renumber
-(
-    const polyMesh& mesh,
-    const pointField& points
-) const
-{
-    CompactListList<label> cellCells;
-    decompositionMethod::calcCellCells
-    (
-        mesh,
-        identity(mesh.nCells()),
-        mesh.nCells(),
-        false,                      // local only
-        cellCells
-    );
-
-    labelList orderedToOld = bandCompression(cellCells());
-
-    if (reverse_)
-    {
-        reverse(orderedToOld);
-    }
-
-    return invert(orderedToOld.size(), orderedToOld);
-}
-
-
-Foam::labelList Foam::CuthillMcKeeRenumber::renumber
-(
-    const labelListList& cellCells,
-    const pointField& points
-) const
-{
-    labelList orderedToOld = bandCompression(cellCells);
-
-    if (reverse_)
-    {
-        reverse(orderedToOld);
-    }
-
-    return invert(orderedToOld.size(), orderedToOld);
-}
-
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H b/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
deleted file mode 100644
index a4ec011b40..0000000000
--- a/src/renumberMethods/CuthillMcKeeRenumber/CuthillMcKeeRenumber.H
+++ /dev/null
@@ -1,118 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::CuthillMcKeeRenumber
-
-Description
-    Cuthill-McKee renumbering
-
-SourceFiles
-    CuthillMcKeeRenumber.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef CuthillMcKeeRenumber_H
-#define CuthillMcKeeRenumber_H
-
-#include "renumberMethod.H"
-#include "Switch.H"
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
-                   Class CuthillMcKeeRenumber Declaration
-\*---------------------------------------------------------------------------*/
-
-class CuthillMcKeeRenumber
-:
-    public renumberMethod
-{
-    // Private data
-
-        const Switch reverse_;
-
-
-    // Private Member Functions
-
-        //- Disallow default bitwise copy construct and assignment
-        void operator=(const CuthillMcKeeRenumber&);
-        CuthillMcKeeRenumber(const CuthillMcKeeRenumber&);
-
-
-public:
-
-    //- Runtime type information
-    TypeName("CuthillMcKee");
-
-
-    // Constructors
-
-        //- Construct given the renumber dictionary
-        CuthillMcKeeRenumber(const dictionary& renumberDict);
-
-
-    //- Destructor
-    virtual ~CuthillMcKeeRenumber()
-    {}
-
-
-    // Member Functions
-
-        //- Return for every coordinate the wanted processor number.
-        //  We need a polyMesh (to be able to load the file)
-        virtual labelList renumber(const pointField&) const
-        {
-            notImplemented("CuthillMcKeeRenumber::renumber(const pointField&)");
-            return labelList(0);
-        }
-
-        //- Return for every coordinate the wanted processor number. Use the
-        //  mesh connectivity (if needed)
-        virtual labelList renumber
-        (
-            const polyMesh& mesh,
-            const pointField& cc
-        ) const;
-
-        //- Return for every cell the new cell label.
-        //  The connectivity is equal to mesh.cellCells() except
-        //  - the connections are across coupled patches
-        virtual labelList renumber
-        (
-            const labelListList& cellCells,
-            const pointField& cc
-        ) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/Make/files b/src/renumberMethods/Make/files
deleted file mode 100644
index c19f1e6b6d..0000000000
--- a/src/renumberMethods/Make/files
+++ /dev/null
@@ -1,8 +0,0 @@
-renumberMethod/renumberMethod.C
-manualRenumber/manualRenumber.C
-CuthillMcKeeRenumber/CuthillMcKeeRenumber.C
-randomRenumber/randomRenumber.C
-springRenumber/springRenumber.C
-boundaryFirstRenumber/boundaryFirstRenumber.C
-
-LIB = $(FOAM_LIBBIN)/librenumberMethods
diff --git a/src/renumberMethods/Make/options b/src/renumberMethods/Make/options
deleted file mode 100644
index 03cb68d946..0000000000
--- a/src/renumberMethods/Make/options
+++ /dev/null
@@ -1,7 +0,0 @@
-EXE_INC = \
-    -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
-    -I$(LIB_SRC)/meshTools/lnInclude
-
-LIB_LIBS = \
-    -ldecompositionMethods \
-    -lmeshTools
diff --git a/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C b/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
deleted file mode 100644
index 293e5ba19b..0000000000
--- a/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.C
+++ /dev/null
@@ -1,201 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "boundaryFirstRenumber.H"
-#include "addToRunTimeSelectionTable.H"
-#include "bandCompression.H"
-#include "decompositionMethod.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(boundaryFirstRenumber, 0);
-
-    addToRunTimeSelectionTable
-    (
-        renumberMethod,
-        boundaryFirstRenumber,
-        dictionary
-    );
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::boundaryFirstRenumber::boundaryFirstRenumber
-(
-    const dictionary& renumberDict
-)
-:
-    renumberMethod(renumberDict),
-    reverse_
-    (
-        renumberDict.found(typeName + "Coeffs")
-      ? Switch(renumberDict.subDict(typeName + "Coeffs").lookup("reverse"))
-      : Switch(true)
-    )
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::labelList Foam::boundaryFirstRenumber::renumber
-(
-    const polyMesh& mesh,
-    const pointField& points
-) const
-{
-    const labelList& own = mesh.faceOwner();
-    const labelList& nei = mesh.faceNeighbour();
-
-    // Distance to boundary. Minimum of neighbours.
-    labelList distance(mesh.nCells(), -1);
-
-    // New order
-    labelList newOrder(mesh.nCells());
-    label cellInOrder = 0;
-
-
-    // Starting faces for walk. These are the zero distance faces.
-    DynamicList<label> frontFaces(mesh.nFaces()-mesh.nInternalFaces());
-
-    const polyBoundaryMesh& patches = mesh.boundaryMesh();
-
-    forAll(patches, patchI)
-    {
-        const polyPatch& pp = patches[patchI];
-
-        //- Note: cannot check for empty since these are introduced by
-        //  fvMeshSubset. Also most loops don't care about patch type.
-        //if (!isA<emptyPolyPatch>(pp))
-        {
-            forAll(pp, i)
-            {
-                frontFaces.append(pp.start()+i);
-            }
-        }
-    }
-
-    // Loop over all frontFaces
-    label currentDistance = 0;
-    while (frontFaces.size() > 0)
-    {
-        DynamicList<label> frontCells(frontFaces.size());
-
-        // Set all of frontFaces' neighbours to current distance
-
-        forAll(frontFaces, i)
-        {
-            label faceI = frontFaces[i];
-
-            if (mesh.isInternalFace(faceI))
-            {
-                label ownCellI = own[faceI];
-                if (distance[ownCellI] == -1)
-                {
-                    distance[ownCellI] = currentDistance;
-                    frontCells.append(ownCellI);
-                }
-                label neiCellI = nei[faceI];
-                if (distance[neiCellI] == -1)
-                {
-                    distance[neiCellI] = currentDistance;
-                    frontCells.append(neiCellI);
-                }
-            }
-            else
-            {
-                label ownCellI = own[faceI];
-                if (distance[ownCellI] == -1)
-                {
-                    distance[ownCellI] = currentDistance;
-                    frontCells.append(ownCellI);
-                }
-            }
-        }
-
-        //Pout<< "For distance:" << currentDistance
-        //    << " from " << frontFaces.size() << " faces to "
-        //    << frontCells.size() << " cells." << endl;
-
-
-        // TBD. Determine order within current shell (frontCells). For now
-        // just add them.
-        forAll(frontCells, i)
-        {
-            newOrder[cellInOrder] = frontCells[i];
-            cellInOrder++;
-        }
-
-        // From cells to faces
-        frontFaces.clear();
-        forAll(frontCells, i)
-        {
-            label cellI = frontCells[i];
-            const cell& cFaces = mesh.cells()[cellI];
-
-            forAll(cFaces, i)
-            {
-                label faceI = cFaces[i];
-                if (mesh.isInternalFace(faceI))
-                {
-                    label nbrCellI =
-                    (
-                        mesh.faceOwner()[faceI] == cellI
-                      ? mesh.faceNeighbour()[faceI]
-                      : mesh.faceOwner()[faceI]
-                    );
-                    if (distance[nbrCellI] == -1)
-                    {
-                        frontFaces.append(faceI);
-                    }
-                }
-            }
-        }
-
-        currentDistance++;
-    }
-
-    // Return furthest away cell first
-    if (reverse_)
-    {
-        reverse(newOrder);
-    }
-
-    //forAll(newOrder, i)
-    //{
-    //    label cellI = newOrder[i];
-    //
-    //    Pout<< "cell:" << cellI << endl;
-    //    Pout<< " at distance:" << distance[cellI]
-    //        << endl;
-    //}
-
-    return invert(newOrder.size(), newOrder);
-}
-
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H b/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H
deleted file mode 100644
index 690af21455..0000000000
--- a/src/renumberMethods/boundaryFirstRenumber/boundaryFirstRenumber.H
+++ /dev/null
@@ -1,129 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::boundaryFirstRenumber
-
-Description
-    Cuthill-McKee renumbering
-
-SourceFiles
-    boundaryFirstRenumber.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef boundaryFirstRenumber_H
-#define boundaryFirstRenumber_H
-
-#include "renumberMethod.H"
-#include "Switch.H"
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
-                   Class boundaryFirstRenumber Declaration
-\*---------------------------------------------------------------------------*/
-
-class boundaryFirstRenumber
-:
-    public renumberMethod
-{
-    // Private data
-
-        const Switch reverse_;
-
-
-    // Private Member Functions
-
-        //- Disallow default bitwise copy construct and assignment
-        void operator=(const boundaryFirstRenumber&);
-        boundaryFirstRenumber(const boundaryFirstRenumber&);
-
-
-public:
-
-    //- Runtime type information
-    TypeName("boundaryFirst");
-
-
-    // Constructors
-
-        //- Construct given the renumber dictionary
-        boundaryFirstRenumber(const dictionary& renumberDict);
-
-
-    //- Destructor
-    virtual ~boundaryFirstRenumber()
-    {}
-
-
-    // Member Functions
-
-        //- Return for every coordinate the wanted processor number.
-        //  We need a polyMesh (to be able to load the file)
-        virtual labelList renumber(const pointField&) const
-        {
-            notImplemented
-            (
-                "boundaryFirstRenumber::renumber(const pointField&)"
-            );
-            return labelList(0);
-        }
-
-        //- Return for every coordinate the wanted processor number. Use the
-        //  mesh connectivity (if needed)
-        virtual labelList renumber
-        (
-            const polyMesh& mesh,
-            const pointField& cc
-        ) const;
-
-        //- Return for every cell the new cell label.
-        //  The connectivity is equal to mesh.cellCells() except
-        //  - the connections are across coupled patches
-        virtual labelList renumber
-        (
-            const labelListList& cellCells,
-            const pointField& cc
-        ) const
-        {
-            notImplemented
-            (
-                "boundaryFirstRenumber::renumber"
-                "(const labelListList&, const pointField&)"
-            );
-            return labelList(0);
-        }
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/manualRenumber/manualRenumber.C b/src/renumberMethods/manualRenumber/manualRenumber.C
deleted file mode 100644
index 92b5c8e00c..0000000000
--- a/src/renumberMethods/manualRenumber/manualRenumber.C
+++ /dev/null
@@ -1,137 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "manualRenumber.H"
-#include "addToRunTimeSelectionTable.H"
-#include "IFstream.H"
-#include "labelIOList.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(manualRenumber, 0);
-
-    addToRunTimeSelectionTable
-    (
-        renumberMethod,
-        manualRenumber,
-        dictionary
-    );
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::manualRenumber::manualRenumber(const dictionary& renumberDict)
-:
-    renumberMethod(renumberDict),
-    dataFile_
-    (
-        renumberDict.subDict(typeName+"Coeffs").lookup("dataFile")
-    )
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::labelList Foam::manualRenumber::renumber
-(
-    const polyMesh& mesh,
-    const pointField& points
-) const
-{
-    labelIOList oldToNew
-    (
-        IOobject
-        (
-            dataFile_,
-            mesh.facesInstance(),
-            mesh,
-            IOobject::MUST_READ,
-            IOobject::AUTO_WRITE,
-            false
-        )
-    );
-
-    // check if the final renumbering is OK
-
-    if (oldToNew.size() != points.size())
-    {
-        FatalErrorIn
-        (
-            "manualRenumber::renumber(const pointField&, const scalarField&)"
-        )   << "Size of renumber list does not correspond "
-            << "to the number of points.  Size: "
-            << oldToNew.size() << " Number of points: "
-            << points.size()
-            << ".\n" << "Manual renumbering data read from file "
-            << dataFile_ << "." << endl
-            << exit(FatalError);
-    }
-
-    // Invert to see if one to one
-    labelList newToOld(points.size(), -1);
-    forAll(oldToNew, i)
-    {
-        label newI = oldToNew[i];
-
-        if (newI < 0 || newI >= oldToNew.size())
-        {
-            FatalErrorIn
-            (
-                "manualRenumber::renumber(const pointField&"
-                ", const scalarField&)"
-            )   << "Renumbering is not one-to-one. Index "
-                << i << " maps onto " << newI
-                << ".\n" << "Manual renumbering data read from file "
-                << dataFile_ << "." << endl
-                << exit(FatalError);
-        }
-
-        if (newToOld[newI] == -1)
-        {
-            newToOld[newI] = i;
-        }
-        else
-        {
-            FatalErrorIn
-            (
-                "manualRenumber::renumber(const pointField&"
-                ", const scalarField&)"
-            )   << "Renumbering is not one-to-one. Both index "
-                << newToOld[newI]
-                << " and " << i << " map onto " << newI
-                << ".\n" << "Manual renumbering data read from file "
-                << dataFile_ << "." << endl
-                << exit(FatalError);
-        }
-    }
-
-    return oldToNew;
-}
-
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/manualRenumber/manualRenumber.H b/src/renumberMethods/manualRenumber/manualRenumber.H
deleted file mode 100644
index a4d4b2f46f..0000000000
--- a/src/renumberMethods/manualRenumber/manualRenumber.H
+++ /dev/null
@@ -1,125 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::manualRenumber
-
-Description
-    Renumber given a cell-to-new cell association in a file
-
-SourceFiles
-    manualRenumber.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef manualRenumber_H
-#define manualRenumber_H
-
-#include "renumberMethod.H"
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
-                      Class manualRenumber Declaration
-\*---------------------------------------------------------------------------*/
-
-class manualRenumber
-:
-    public renumberMethod
-{
-    // Private data
-
-        const fileName dataFile_;
-
-
-    // Private Member Functions
-
-        //- Disallow default bitwise copy construct and assignment
-        void operator=(const manualRenumber&);
-        manualRenumber(const manualRenumber&);
-
-
-public:
-
-    //- Runtime type information
-    TypeName("manual");
-
-
-    // Constructors
-
-        //- Construct given the renumber dictionary
-        manualRenumber(const dictionary& renumberDict);
-
-
-    //- Destructor
-    virtual ~manualRenumber()
-    {}
-
-
-    // Member Functions
-
-        //- Return for every coordinate the wanted processor number.
-        //  We need a polyMesh (to be able to load the file)
-        virtual labelList renumber(const pointField&) const
-        {
-            notImplemented("manualRenumber::renumber(const pointField&)");
-            return labelList(0);
-        }
-
-        //- Return for every coordinate the wanted processor number. Use the
-        //  mesh connectivity (if needed)
-        virtual labelList renumber
-        (
-            const polyMesh& mesh,
-            const pointField& cc
-        ) const;
-
-        //- Return for every cell the new cell label.
-        //  The connectivity is equal to mesh.cellCells() except
-        //  - the connections are across coupled patches
-        virtual labelList renumber
-        (
-            const labelListList& cellCells,
-            const pointField& cc
-        ) const
-        {
-            notImplemented
-            (
-                "manualRenumber::renumber"
-                "(const labelListList&, const pointField&)"
-            );
-            return labelList(0);
-        }
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/randomRenumber/randomRenumber.C b/src/renumberMethods/randomRenumber/randomRenumber.C
deleted file mode 100644
index 2dae28e438..0000000000
--- a/src/renumberMethods/randomRenumber/randomRenumber.C
+++ /dev/null
@@ -1,96 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "randomRenumber.H"
-#include "addToRunTimeSelectionTable.H"
-#include "Random.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(randomRenumber, 0);
-
-    addToRunTimeSelectionTable
-    (
-        renumberMethod,
-        randomRenumber,
-        dictionary
-    );
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::randomRenumber::randomRenumber(const dictionary& renumberDict)
-:
-    renumberMethod(renumberDict)
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::labelList Foam::randomRenumber::renumber
-(
-    const pointField& points
-) const
-{
-    Random rndGen(0);
-
-    labelList oldToNew(identity(points.size()));
-
-    for (label iter = 0; iter < 10; iter++)
-    {
-        forAll(oldToNew, i)
-        {
-            label j = rndGen.integer(0, oldToNew.size()-1);
-            Swap(oldToNew[i], oldToNew[j]);
-        }
-    }
-    return oldToNew;
-}
-
-
-Foam::labelList Foam::randomRenumber::renumber
-(
-    const polyMesh& mesh,
-    const pointField& points
-) const
-{
-    return renumber(points);
-}
-
-
-Foam::labelList Foam::randomRenumber::renumber
-(
-    const labelListList& cellCells,
-    const pointField& points
-) const
-{
-    return renumber(points);
-}
-
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/randomRenumber/randomRenumber.H b/src/renumberMethods/randomRenumber/randomRenumber.H
deleted file mode 100644
index 0665d5159b..0000000000
--- a/src/renumberMethods/randomRenumber/randomRenumber.H
+++ /dev/null
@@ -1,108 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::randomRenumber
-
-Description
-    Random renumber. Just to see effect of renumbering.
-
-SourceFiles
-    randomRenumber.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef randomRenumber_H
-#define randomRenumber_H
-
-#include "renumberMethod.H"
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
-                      Class randomRenumber Declaration
-\*---------------------------------------------------------------------------*/
-
-class randomRenumber
-:
-    public renumberMethod
-{
-    // Private Member Functions
-
-        //- Disallow default bitwise copy construct and assignment
-        void operator=(const randomRenumber&);
-        randomRenumber(const randomRenumber&);
-
-
-public:
-
-    //- Runtime type information
-    TypeName("random");
-
-
-    // Constructors
-
-        //- Construct given the renumber dictionary
-        randomRenumber(const dictionary& renumberDict);
-
-
-    //- Destructor
-    virtual ~randomRenumber()
-    {}
-
-
-    // Member Functions
-
-        //- Return for every coordinate the wanted processor number.
-        //  We need a polyMesh (to be able to load the file)
-        virtual labelList renumber(const pointField&) const;
-
-        //- Return for every coordinate the wanted processor number. Use the
-        //  mesh connectivity (if needed)
-        virtual labelList renumber
-        (
-            const polyMesh& mesh,
-            const pointField& cc
-        ) const;
-
-        //- Return for every cell the new cell label.
-        //  The connectivity is equal to mesh.cellCells() except
-        //  - the connections are across coupled patches
-        virtual labelList renumber
-        (
-            const labelListList& cellCells,
-            const pointField& cc
-        ) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/renumberMethod/renumberMethod.C b/src/renumberMethods/renumberMethod/renumberMethod.C
deleted file mode 100644
index 76ba1911e2..0000000000
--- a/src/renumberMethods/renumberMethod/renumberMethod.C
+++ /dev/null
@@ -1,131 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-InClass
-    renumberMethod
-
-\*---------------------------------------------------------------------------*/
-
-#include "renumberMethod.H"
-#include "decompositionMethod.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(renumberMethod, 0);
-    defineRunTimeSelectionTable(renumberMethod, dictionary);
-}
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::autoPtr<Foam::renumberMethod> Foam::renumberMethod::New
-(
-    const dictionary& renumberDict
-)
-{
-    const word methodType(renumberDict.lookup("method"));
-
-    //Info<< "Selecting renumberMethod " << methodType << endl;
-
-    dictionaryConstructorTable::iterator cstrIter =
-        dictionaryConstructorTablePtr_->find(methodType);
-
-    if (cstrIter == dictionaryConstructorTablePtr_->end())
-    {
-        FatalErrorIn
-        (
-            "renumberMethod::New"
-            "(const dictionary& renumberDict)"
-        )   << "Unknown renumberMethod "
-            << methodType << nl << nl
-            << "Valid renumberMethods are : " << endl
-            << dictionaryConstructorTablePtr_->sortedToc()
-            << exit(FatalError);
-    }
-
-    return autoPtr<renumberMethod>(cstrIter()(renumberDict));
-}
-
-
-Foam::labelList Foam::renumberMethod::renumber
-(
-    const polyMesh& mesh,
-    const pointField& points
-) const
-{
-    CompactListList<label> cellCells;
-    decompositionMethod::calcCellCells
-    (
-        mesh,
-        identity(mesh.nCells()),
-        mesh.nCells(),
-        false,                      // local only
-        cellCells
-    );
-
-    // Renumber based on agglomerated points
-    return renumber(cellCells(), points);
-}
-
-
-Foam::labelList Foam::renumberMethod::renumber
-(
-    const polyMesh& mesh,
-    const labelList& fineToCoarse,
-    const pointField& coarsePoints
-) const
-{
-    CompactListList<label> coarseCellCells;
-    decompositionMethod::calcCellCells
-    (
-        mesh,
-        fineToCoarse,
-        coarsePoints.size(),
-        false,                      // local only
-        coarseCellCells
-    );
-
-    // Renumber based on agglomerated points
-    labelList coarseDistribution
-    (
-        renumber
-        (
-            coarseCellCells(),
-            coarsePoints
-        )
-    );
-
-    // Rework back into renumbering for original mesh_
-    labelList fineDistribution(fineToCoarse.size());
-
-    forAll(fineDistribution, i)
-    {
-        fineDistribution[i] = coarseDistribution[fineToCoarse[i]];
-    }
-
-    return fineDistribution;
-}
-
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/renumberMethod/renumberMethod.H b/src/renumberMethods/renumberMethod/renumberMethod.H
deleted file mode 100644
index ead92b1e41..0000000000
--- a/src/renumberMethods/renumberMethod/renumberMethod.H
+++ /dev/null
@@ -1,161 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::renumberMethod
-
-Description
-    Abstract base class for renumbering
-
-SourceFiles
-    renumberMethod.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef renumberMethod_H
-#define renumberMethod_H
-
-#include "polyMesh.H"
-#include "pointField.H"
-#include "CompactListList.H"
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
-                           Class renumberMethod Declaration
-\*---------------------------------------------------------------------------*/
-
-class renumberMethod
-{
-
-protected:
-
-    // Protected data
-
-        const dictionary& renumberDict_;
-
-private:
-
-    // Private Member Functions
-
-        //- Disallow default bitwise copy construct and assignment
-        renumberMethod(const renumberMethod&);
-        void operator=(const renumberMethod&);
-
-
-public:
-
-    //- Runtime type information
-    TypeName("renumberMethod");
-
-
-    // Declare run-time constructor selection tables
-
-        declareRunTimeSelectionTable
-        (
-            autoPtr,
-            renumberMethod,
-            dictionary,
-            (
-                const dictionary& renumberDict
-            ),
-            (renumberDict)
-        );
-
-
-    // Selectors
-
-        //- Return a reference to the selected renumbering method
-        static autoPtr<renumberMethod> New
-        (
-            const dictionary& renumberDict
-        );
-
-
-    // Constructors
-
-        //- Construct given the renumber dictionary
-        renumberMethod(const dictionary& renumberDict)
-        :
-            renumberDict_(renumberDict)
-        {}
-
-
-    //- Destructor
-    virtual ~renumberMethod()
-    {}
-
-
-    // Member Functions
-
-        //- Return for every cell the new cell label.
-        //  This is only defined for geometric renumberMethods.
-        virtual labelList renumber(const pointField&) const
-        {
-            notImplemented
-            (
-                "renumberMethod:renumber(const pointField&)"
-            );
-            return labelList(0);
-        }
-
-        //- Return for every cell the new cell label. Use the
-        //  mesh connectivity (if needed)
-        virtual labelList renumber(const polyMesh&, const pointField&) const;
-
-        //- Return for every cell the new cell label. Gets
-        //  passed agglomeration map (from fine to coarse cells) and coarse
-        //  cell
-        //  location. Can be overridden by renumberMethods that provide this
-        //  functionality natively. Coarse cells are local to the processor
-        //  (if in parallel). If you want to have coarse cells spanning
-        //  processors use the globalCellCells instead.
-        virtual labelList renumber
-        (
-            const polyMesh& mesh,
-            const labelList& cellToRegion,
-            const pointField& regionPoints
-        ) const;
-
-        //- Return for every cell the new cell label.
-        //  The connectivity is equal to mesh.cellCells() except
-        //  - the connections are across coupled patches
-        virtual labelList renumber
-        (
-            const labelListList& cellCells,
-            const pointField& cc
-        ) const = 0;
-
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/springRenumber/springRenumber.C b/src/renumberMethods/springRenumber/springRenumber.C
deleted file mode 100644
index 64a15145ce..0000000000
--- a/src/renumberMethods/springRenumber/springRenumber.C
+++ /dev/null
@@ -1,172 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "springRenumber.H"
-#include "addToRunTimeSelectionTable.H"
-#include "decompositionMethod.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(springRenumber, 0);
-
-    addToRunTimeSelectionTable
-    (
-        renumberMethod,
-        springRenumber,
-        dictionary
-    );
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::springRenumber::springRenumber(const dictionary& renumberDict)
-:
-    renumberMethod(renumberDict),
-    dict_(renumberDict.subDict(typeName+"Coeffs")),
-    maxCo_(readScalar(dict_.lookup("maxCo"))),
-    maxIter_(readLabel(dict_.lookup("maxIter"))),
-    freezeFraction_(readScalar(dict_.lookup("freezeFraction")))
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::labelList Foam::springRenumber::renumber
-(
-    const polyMesh& mesh,
-    const pointField& points
-) const
-{
-    CompactListList<label> cellCells;
-    decompositionMethod::calcCellCells
-    (
-        mesh,
-        identity(mesh.nCells()),
-        mesh.nCells(),
-        false,                      // local only
-        cellCells
-    );
-
-    return renumber(cellCells(), points);
-}
-
-
-Foam::labelList Foam::springRenumber::renumber
-(
-    const labelListList& cellCells,
-    const pointField& points
-) const
-{
-    // Look at cell index as a 1D position parameter.
-    // Move cells to the average 'position' of their neighbour.
-
-    scalarField position(cellCells.size());
-    forAll(position, cellI)
-    {
-        position[cellI] = cellI;
-    }
-
-    labelList oldToNew(identity(cellCells.size()));
-
-    scalar maxCo = maxCo_ * cellCells.size();
-
-    for (label iter = 0; iter < maxIter_; iter++)
-    {
-        //Pout<< "Iteration : " << iter << nl
-        //    << "------------"
-        //    << endl;
-
-        //Pout<< "Position :" << nl
-        //    << "    min : " << min(position) << nl
-        //    << "    max : " << max(position) << nl
-        //    << "    avg : " << average(position) << nl
-        //    << endl;
-
-        // Sum force per cell.
-        scalarField sumForce(cellCells.size(), 0.0);
-        forAll(cellCells, oldCellI)
-        {
-            const labelList& cCells = cellCells[oldCellI];
-            label cellI = oldToNew[oldCellI];
-
-            forAll(cCells, i)
-            {
-                label nbrCellI = oldToNew[cCells[i]];
-
-                sumForce[cellI] += (position[nbrCellI]-position[cellI]);
-            }
-        }
-
-        //Pout<< "Force :" << nl
-        //    << "    min    : " << min(sumForce) << nl
-        //    << "    max    : " << max(sumForce) << nl
-        //    << "    avgMag : " << average(mag(sumForce)) << nl
-        //    << "DeltaT : " << deltaT << nl
-        //    << endl;
-
-        // Limit displacement
-        scalar deltaT = maxCo / max(mag(sumForce));
-
-        Info<< "Iter:" << iter
-            << "  maxCo:" << maxCo
-            << "  deltaT:" << deltaT
-            << "  average force:" << average(mag(sumForce)) << endl;
-
-        // Determine displacement.
-        scalarField displacement = deltaT*sumForce;
-
-        //Pout<< "Displacement :" << nl
-        //    << "    min    : " << min(displacement) << nl
-        //    << "    max    : " << max(displacement) << nl
-        //    << "    avgMag : " << average(mag(displacement)) << nl
-        //    << endl;
-
-        // Calculate new position and scale to be within original range
-        // (0..nCells-1) for ease of postprocessing.
-        position += displacement;
-        position -= min(position);
-        position *= (position.size()-1)/max(position);
-
-        // Slowly freeze.
-        maxCo *= freezeFraction_;
-    }
-
-    //writeOBJ("endPosition.obj", cellCells, position);
-
-    // Move cells to new position
-    labelList shuffle;
-    sortedOrder(position, shuffle);
-
-    // Reorder oldToNew
-    inplaceReorder(shuffle, oldToNew);
-
-    return oldToNew;
-}
-
-
-// ************************************************************************* //
diff --git a/src/renumberMethods/springRenumber/springRenumber.H b/src/renumberMethods/springRenumber/springRenumber.H
deleted file mode 100644
index 30f7a3d838..0000000000
--- a/src/renumberMethods/springRenumber/springRenumber.H
+++ /dev/null
@@ -1,132 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::springRenumber
-
-Description
-    Use spring analogy - attract neighbouring cells according to the distance
-    of their cell indices.
-
-    // Maximum jump of cell indices. Is fraction of number of cells
-    maxCo 0.1;
-
-    // Limit the amount of movement; the fraction maxCo gets decreased
-    // with every iteration.
-    freezeFraction 0.9;
-
-    // Maximum number of iterations
-    maxIter 1000;
-
-SourceFiles
-    springRenumber.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef springRenumber_H
-#define springRenumber_H
-
-#include "renumberMethod.H"
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
-                       Class springRenumber Declaration
-\*---------------------------------------------------------------------------*/
-
-class springRenumber
-:
-    public renumberMethod
-{
-    // Private data
-
-        const dictionary& dict_;
-
-        const scalar maxCo_;
-
-        const label maxIter_;
-
-        const scalar freezeFraction_;
-
-    // Private Member Functions
-
-        //- Disallow default bitwise copy construct and assignment
-        void operator=(const springRenumber&);
-        springRenumber(const springRenumber&);
-
-
-public:
-
-    //- Runtime type information
-    TypeName("spring");
-
-
-    // Constructors
-
-        //- Construct given the renumber dictionary
-        springRenumber(const dictionary& renumberDict);
-
-
-    //- Destructor
-    virtual ~springRenumber()
-    {}
-
-
-    // Member Functions
-
-        //- Return for every coordinate the wanted processor number.
-        virtual labelList renumber(const pointField&) const
-        {
-            notImplemented("springRenumber::renumber(const pointField&)");
-            return labelList(0);
-        }
-
-        //- Return for every coordinate the wanted processor number. Use the
-        //  mesh connectivity (if needed)
-        virtual labelList renumber
-        (
-            const polyMesh& mesh,
-            const pointField& cc
-        ) const;
-
-        //- Return for every cell the new cell label.
-        //  The connectivity is equal to mesh.cellCells() except
-        //  - the connections are across coupled patches
-        virtual labelList renumber
-        (
-            const labelListList& cellCells,
-            const pointField& cc
-        ) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

From f4eeef024e1252e1a248add41b792fe4767916aa Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 23 Dec 2011 15:36:43 +0000
Subject: [PATCH 17/18] STYLE: duplicate ;;

---
 .../utilities/mesh/conversion/gambitToFoam/gambitToFoam.L | 2 +-
 .../USERD_get_number_of_files_in_dataset.H                | 2 +-
 .../cellPointWeight/cellPointWeight.C                     | 8 ++++----
 3 files changed, 6 insertions(+), 6 deletions(-)

diff --git a/applications/utilities/mesh/conversion/gambitToFoam/gambitToFoam.L b/applications/utilities/mesh/conversion/gambitToFoam/gambitToFoam.L
index 7706ccb272..dc4531fef6 100644
--- a/applications/utilities/mesh/conversion/gambitToFoam/gambitToFoam.L
+++ b/applications/utilities/mesh/conversion/gambitToFoam/gambitToFoam.L
@@ -433,7 +433,7 @@ mtype                 {space}"MTYPE:"{space}
         }
         else
         {
-            curGroupID = readLabel(groupStream);;
+            curGroupID = readLabel(groupStream);
         }
 
         BEGIN(cellStreams);
diff --git a/applications/utilities/postProcessing/graphics/ensightFoamReader/USERD_get_number_of_files_in_dataset.H b/applications/utilities/postProcessing/graphics/ensightFoamReader/USERD_get_number_of_files_in_dataset.H
index e2246d27da..ed936ee082 100644
--- a/applications/utilities/postProcessing/graphics/ensightFoamReader/USERD_get_number_of_files_in_dataset.H
+++ b/applications/utilities/postProcessing/graphics/ensightFoamReader/USERD_get_number_of_files_in_dataset.H
@@ -9,5 +9,5 @@ int USERD_get_number_of_files_in_dataset(void)
 
     // use 1 insted of 0 which gives an un-necessary warning.
     Num_dataset_files = 1;
-    return Num_dataset_files;;
+    return Num_dataset_files;
 }
diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/cellPointWeight/cellPointWeight.C b/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/cellPointWeight/cellPointWeight.C
index 3ef207b702..13828ff7cf 100644
--- a/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/cellPointWeight/cellPointWeight.C
+++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/cellPointWeight/cellPointWeight.C
@@ -83,8 +83,8 @@ void Foam::cellPointWeight::findTetrahedron
             )
             {
                 faceVertices_[0] = f[tetIs.faceBasePt()];
-                faceVertices_[1] = f[tetIs.facePtA()];;
-                faceVertices_[2] = f[tetIs.facePtB()];;
+                faceVertices_[1] = f[tetIs.facePtA()];
+                faceVertices_[2] = f[tetIs.facePtB()];
 
                 return;
             }
@@ -191,8 +191,8 @@ void Foam::cellPointWeight::findTriangle
                 weights_[3] = triWeights[2];
 
                 faceVertices_[0] = f[tetIs.faceBasePt()];
-                faceVertices_[1] = f[tetIs.facePtA()];;
-                faceVertices_[2] = f[tetIs.facePtB()];;
+                faceVertices_[1] = f[tetIs.facePtA()];
+                faceVertices_[2] = f[tetIs.facePtB()];
 
                 return;
             }

From b51dbec8d042fc4563014d4bc45537ffd7eb076e Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 23 Dec 2011 15:41:50 +0000
Subject: [PATCH 18/18] BUG: triangeI.H: barycentric coordinates error

---
 src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H b/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
index 3867915943..496442012d 100644
--- a/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
+++ b/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
@@ -277,9 +277,9 @@ Foam::scalar Foam::triangle<Point, PointRef>::barycentric
 
     bary.setSize(3);
 
-    bary[0] = (d11*d20 - d01*d21)/denom;
-    bary[1] = (d00*d21 - d01*d20)/denom;
-    bary[2] = 1.0 - bary[0] - bary[1];
+    bary[1] = (d11*d20 - d01*d21)/denom;
+    bary[2] = (d00*d21 - d01*d20)/denom;
+    bary[0] = 1.0 - bary[1] - bary[2];
 
     return denom;
 }