diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Kunz/Kunz.H b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Kunz/Kunz.H
index 74f8f2b931..39ccad5382 100644
--- a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Kunz/Kunz.H
+++ b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Kunz/Kunz.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -108,9 +108,8 @@ public:
         //  and a coefficient to multiply  alphal for the vaporisation rate
         virtual Pair<tmp<volScalarField> > mDotAlphal() const;
 
-        //- Return the mass condensation and vaporisation rates as an
-        //  explicit term for the condensation rate and a coefficient to
-        //  multiply (p - pSat) for the vaporisation rate
+        //- Return the mass condensation and vaporisation rates as coefficients
+        //  to multiply (p - pSat)
         virtual Pair<tmp<volScalarField> > mDotP() const;
 
         //- Correct the Kunz phaseChange model
diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Merkle/Merkle.H b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Merkle/Merkle.H
index 9351afd2de..7dd8e841c8 100644
--- a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Merkle/Merkle.H
+++ b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/Merkle/Merkle.H
@@ -1,8 +1,8 @@
 /*---------------------------------------------------------------------------*\
-  ========Merkle=                 |
+  =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -102,9 +102,8 @@ public:
         //  and a coefficient to multiply  alphal for the vaporisation rate
         virtual Pair<tmp<volScalarField> > mDotAlphal() const;
 
-        //- Return the mass condensation and vaporisation rates as an
-        //  explicit term for the condensation rate and a coefficient to
-        //  multiply (p - pSat) for the vaporisation rate
+        //- Return the mass condensation and vaporisation rates as coefficients
+        //  to multiply (p - pSat)
         virtual Pair<tmp<volScalarField> > mDotP() const;
 
         //- Correct the Merkle phaseChange model
diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/SchnerrSauer/SchnerrSauer.H b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/SchnerrSauer/SchnerrSauer.H
index dd676b8fd0..beef7fc7fc 100644
--- a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/SchnerrSauer/SchnerrSauer.H
+++ b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/SchnerrSauer/SchnerrSauer.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -115,9 +115,8 @@ public:
         //  and a coefficient to multiply  alphal for the vaporisation rate
         virtual Pair<tmp<volScalarField> > mDotAlphal() const;
 
-        //- Return the mass condensation and vaporisation rates as an
-        //  explicit term for the condensation rate and a coefficient to
-        //  multiply (p - pSat) for the vaporisation rate
+        //- Return the mass condensation and vaporisation rates as coefficients
+        //  to multiply (p - pSat)
         virtual Pair<tmp<volScalarField> > mDotP() const;
 
         //- Correct the SchnerrSauer phaseChange model
diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/phaseChangeTwoPhaseMixture/phaseChangeTwoPhaseMixture.H b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/phaseChangeTwoPhaseMixture/phaseChangeTwoPhaseMixture.H
index fc81d080d3..e05476390b 100644
--- a/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/phaseChangeTwoPhaseMixture/phaseChangeTwoPhaseMixture.H
+++ b/applications/solvers/multiphase/interPhaseChangeFoam/phaseChangeTwoPhaseMixtures/phaseChangeTwoPhaseMixture/phaseChangeTwoPhaseMixture.H
@@ -139,9 +139,8 @@ public:
         //  and a coefficient to multiply  alphal for the vaporisation rate
         virtual Pair<tmp<volScalarField> > mDotAlphal() const = 0;
 
-        //- Return the mass condensation and vaporisation rates as an
-        //  explicit term for the condensation rate and a coefficient to
-        //  multiply (p - pSat) for the vaporisation rate
+        //- Return the mass condensation and vaporisation rates as coefficients
+        //  to multiply (p - pSat)
         virtual Pair<tmp<volScalarField> > mDotP() const = 0;
 
         //- Return the volumetric condensation and vaporisation rates as a
@@ -149,9 +148,8 @@ public:
         //  and a coefficient to multiply  alphal for the vaporisation rate
         Pair<tmp<volScalarField> > vDotAlphal() const;
 
-        //- Return the volumetric condensation and vaporisation rates as an
-        //  explicit term for the condensation rate and a coefficient to
-        //  multiply (p - pSat) for the vaporisation rate
+        //- Return the volumetric condensation and vaporisation rates as
+        //  coefficients to multiply (p - pSat)
         Pair<tmp<volScalarField> > vDotP() const;
 
         //- Correct the phaseChange model
diff --git a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C
index 76b4dab8dd..a154565747 100644
--- a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C
+++ b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C
@@ -111,11 +111,6 @@ Foam::multiphaseMixture::multiphaseMixture
 {
     calcAlphas();
     alphas_.write();
-
-    forAllIter(PtrDictionary<phase>, phases_, iter)
-    {
-        alphaTable_.add(iter());
-    }
 }
 
 
diff --git a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.H b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.H
index 19b436ef2a..ce97221eea 100644
--- a/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.H
+++ b/applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -47,7 +47,6 @@ SourceFiles
 #include "PtrDictionary.H"
 #include "volFields.H"
 #include "surfaceFields.H"
-#include "multivariateSurfaceInterpolationScheme.H"
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -156,9 +155,6 @@ private:
         //- Conversion factor for degrees into radians
         static const scalar convertToRad;
 
-        //- Phase-fraction field table for multivariate discretisation
-        multivariateSurfaceInterpolationScheme<scalar>::fieldTable alphaTable_;
-
 
     // Private member functions
 
diff --git a/applications/utilities/preProcessing/viewFactorsGen/searchingEngine.H b/applications/utilities/preProcessing/viewFactorsGen/searchingEngine.H
index ed1d9031ec..ef5b264a54 100644
--- a/applications/utilities/preProcessing/viewFactorsGen/searchingEngine.H
+++ b/applications/utilities/preProcessing/viewFactorsGen/searchingEngine.H
@@ -26,7 +26,8 @@ dict.add("mergeDistance", SMALL);
 labelHashSet includePatches;
 forAll(patches, patchI)
 {
-    if (!isA<processorPolyPatch>(patches[patchI]))
+    const polyPatch& pp = patches[patchI];
+    if (!pp.coupled() && !isA<cyclicAMIPolyPatch>(pp))
     {
         includePatches.insert(patchI);
     }
diff --git a/applications/utilities/preProcessing/viewFactorsGen/viewFactorsGen.C b/applications/utilities/preProcessing/viewFactorsGen/viewFactorsGen.C
index c9aa219d13..8276528945 100644
--- a/applications/utilities/preProcessing/viewFactorsGen/viewFactorsGen.C
+++ b/applications/utilities/preProcessing/viewFactorsGen/viewFactorsGen.C
@@ -44,6 +44,7 @@ Description
 #include "volFields.H"
 #include "surfaceFields.H"
 #include "distributedTriSurfaceMesh.H"
+#include "cyclicAMIPolyPatch.H"
 #include "triSurfaceTools.H"
 #include "mapDistribute.H"
 
diff --git a/etc/controlDict b/etc/controlDict
index b9c41cb3a4..fba02069de 100644
--- a/etc/controlDict
+++ b/etc/controlDict
@@ -135,7 +135,7 @@ DebugSwitches
     FDIC                0;
     FaceCellWave        0;
     GAMG                0;
-    GAMGAgglomeration   0;
+    GAMGAgglomeration   1;
     GAMGInterface       0;
     GAMGInterfaceField  0;
     Gamma               0;
diff --git a/src/OpenFOAM/Make/files b/src/OpenFOAM/Make/files
index 19a8138fe5..6d346617c1 100644
--- a/src/OpenFOAM/Make/files
+++ b/src/OpenFOAM/Make/files
@@ -331,10 +331,8 @@ eagerGAMGProcAgglomeration = $(GAMGProcAgglomerations)/eagerGAMGProcAgglomeratio
 $(eagerGAMGProcAgglomeration)/eagerGAMGProcAgglomeration.C
 noneGAMGProcAgglomeration = $(GAMGProcAgglomerations)/noneGAMGProcAgglomeration
 $(noneGAMGProcAgglomeration)/noneGAMGProcAgglomeration.C
-/*
-cellFaceRatioGAMGProcAgglomeration = $(GAMGProcAgglomerations)/cellFaceRatioGAMGProcAgglomeration
-$(cellFaceRatioGAMGProcAgglomeration)/cellFaceRatioGAMGProcAgglomeration.C
-*/
+procFacesGAMGProcAgglomeration = $(GAMGProcAgglomerations)/procFacesGAMGProcAgglomeration
+$(procFacesGAMGProcAgglomeration)/procFacesGAMGProcAgglomeration.C
 
 
 meshes/lduMesh/lduMesh.C
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerateLduAddressing.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerateLduAddressing.C
index 511aeb70e2..066aa39d3b 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerateLduAddressing.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerateLduAddressing.C
@@ -351,7 +351,7 @@ void Foam::GAMGAgglomeration::agglomerateLduAddressing
     );
 
 
-    if (debug)
+    if (debug & 2)
     {
         Pout<< "GAMGAgglomeration :"
             << " agglomerated level " << fineLevelIndex
@@ -488,7 +488,10 @@ void Foam::GAMGAgglomeration::procAgglomerateRestrictAddressing
         comm,
         procIDs,
         restrictAddressing_[levelIndex],
-        procRestrictAddressing
+        procRestrictAddressing,
+
+        UPstream::msgType(),
+        Pstream::nonBlocking    //Pstream::scheduled
     );
 
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C
index d1551eb399..a53a5fe9ab 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C
@@ -27,9 +27,9 @@ License
 #include "lduMesh.H"
 #include "lduMatrix.H"
 #include "Time.H"
-#include "dlLibraryTable.H"
 #include "GAMGInterface.H"
 #include "GAMGProcAgglomeration.H"
+#include "IOmanip.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -67,23 +67,94 @@ void Foam::GAMGAgglomeration::compactLevels(const label nCreatedLevels)
         procBoundaryFaceMap_.setSize(nCreatedLevels);
 
         procAgglomeratorPtr_().agglomerate();
-    }
 
-    if (debug)
-    {
-        for (label levelI = 0; levelI <= size(); levelI++)
+
+        if (debug)
         {
-            if (hasMeshLevel(levelI))
+
+            Info<< "GAMGAgglomeration:" << nl
+                << "    local agglomerator     : " << type() << nl
+                << "    processor agglomerator : "
+                << procAgglomeratorPtr_().type() << nl
+                << nl;
+
+            Info<< setw(40) << "nCells"
+                << setw(24) << "nInterfaces"
+                << setw(24) << "Ratio" << nl
+                << setw(8) << "Level"
+                << setw(8) << "nProcs"
+                << "        "
+                << setw(8) << "avg"
+                << setw(8) << "max"
+                << "        "
+                << setw(8) << "avg"
+                << setw(8) << "max"
+                << "            " << setw(4) << "avg"
+                << "    " << setw(4) << "max"
+                << nl
+                << setw(8) << "-----"
+                << setw(8) << "------"
+                << "        "
+                << setw(8) << "---"
+                << setw(8) << "---"
+                << "        "
+                << setw(8) << "---"
+                << setw(8) << "---"
+                << "            " << setw(4) << "---"
+                << "    " << setw(4) << "---"
+                << nl;
+
+            for (label levelI = 0; levelI <= size(); levelI++)
             {
-                const lduMesh& fineMesh = meshLevel(levelI);
-                Pout<< "Level " << levelI << " fine mesh:"<< nl;
-                Pout<< fineMesh.info() << endl;
-            }
-            else
-            {
-                Pout<< "Level " << levelI << " has no fine mesh:" << nl
-                    << endl;
+                label nProcs = 0;
+                label nCells = 0;
+                label nInterfaces = 0;
+                label nIntFaces = 0;
+                scalar ratio = 0.0;
+
+                if (hasMeshLevel(levelI))
+                {
+                    nProcs = 1;
+
+                    const lduMesh& fineMesh = meshLevel(levelI);
+                    nCells = fineMesh.lduAddr().size();
+
+                    const lduInterfacePtrsList interfaces =
+                        fineMesh.interfaces();
+                    forAll(interfaces, i)
+                    {
+                        if (interfaces.set(i))
+                        {
+                            nInterfaces++;
+                            nIntFaces += interfaces[i].faceCells().size();
+                        }
+                    }
+                    ratio = scalar(nIntFaces)/nCells;
+                }
+
+                label totNprocs = returnReduce(nProcs, sumOp<label>());
+
+                label maxNCells = returnReduce(nCells, maxOp<label>());
+                label totNCells = returnReduce(nCells, sumOp<label>());
+
+                label maxNInt = returnReduce(nInterfaces, maxOp<label>());
+                label totNInt = returnReduce(nInterfaces, sumOp<label>());
+
+                scalar maxRatio = returnReduce(ratio, maxOp<scalar>());
+                scalar totRatio = returnReduce(ratio, sumOp<scalar>());
+
+                Info<< setw(8) << levelI
+                    << setw(8) << totNprocs
+                    << setw(16) << totNCells/totNprocs
+                    << setw(8) << maxNCells
+                    << setw(16) << totNInt/totNprocs
+                    << setw(8) << maxNInt
+                    << "        "
+                    << setw(8) << setprecision(4) << totRatio/totNprocs
+                    << setw(8) << setprecision(4) << maxRatio
+                    << nl;
             }
+            Info<< endl;
         }
     }
 }
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerationTemplates.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerationTemplates.C
index ce4fa0073a..39d69458e2 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerationTemplates.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerationTemplates.C
@@ -126,7 +126,15 @@ void Foam::GAMGAgglomeration::restrictField
         const List<int>& procIDs = agglomProcIDs(coarseLevelIndex);
         const labelList& offsets = cellOffsets(coarseLevelIndex);
 
-        globalIndex::gather(offsets, fineComm, procIDs, cf);
+        globalIndex::gather
+        (
+            offsets,
+            fineComm,
+            procIDs,
+            cf,
+            UPstream::msgType(),
+            Pstream::nonBlocking    //Pstream::scheduled
+        );
     }
 }
 
@@ -194,7 +202,16 @@ void Foam::GAMGAgglomeration::prolongField
         label localSize = nCells_[levelIndex];
 
         Field<Type> allCf(localSize);
-        globalIndex::scatter(offsets, coarseComm, procIDs, cf, allCf);
+        globalIndex::scatter
+        (
+            offsets,
+            coarseComm,
+            procIDs,
+            cf,
+            allCf,
+            UPstream::msgType(),
+            Pstream::nonBlocking    //Pstream::scheduled
+        );
 
         forAll(fineToCoarse, i)
         {
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomeration.H b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomeration.H
index 2867cac2f6..59de4b1e4d 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomeration.H
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomeration.H
@@ -61,14 +61,6 @@ protected:
 
     // Protected Member Functions
 
-        //- Calculate and return agglomeration of given level
-        tmp<labelField> agglomerate
-        (
-            label& nCoarseCells,
-            const lduAddressing& fineMatrixAddressing,
-            const scalarField& faceWeights
-        );
-
         //- Agglomerate all levels starting from the given face weights
         void agglomerate
         (
@@ -97,6 +89,15 @@ public:
             const lduMesh& mesh,
             const dictionary& controlDict
         );
+
+        //- Calculate and return agglomeration
+        static tmp<labelField> agglomerate
+        (
+            label& nCoarseCells,
+            const lduAddressing& fineMatrixAddressing,
+            const scalarField& faceWeights
+        );
+
 };
 
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/manualGAMGProcAgglomeration/manualGAMGProcAgglomeration.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/manualGAMGProcAgglomeration/manualGAMGProcAgglomeration.C
index 1fea00057c..c67ce2b18b 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/manualGAMGProcAgglomeration/manualGAMGProcAgglomeration.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/manualGAMGProcAgglomeration/manualGAMGProcAgglomeration.C
@@ -109,7 +109,7 @@ bool Foam::manualGAMGProcAgglomeration::agglomerate()
                     // My processor id
                     const label myProcID = Pstream::myProcNo(levelMesh.comm());
 
-                    const List<clusterAndMaster>& clusters =
+                    const List<labelList>& clusters =
                         procAgglomMaps_[i].second();
 
                     // Coarse to fine master processor
@@ -125,8 +125,8 @@ bool Foam::manualGAMGProcAgglomeration::agglomerate()
 
                     forAll(clusters, coarseI)
                     {
-                        const labelList& cluster = clusters[coarseI].first();
-                        coarseToMaster[coarseI] = clusters[coarseI].second();
+                        const labelList& cluster = clusters[coarseI];
+                        coarseToMaster[coarseI] = cluster[0];
 
                         forAll(cluster, i)
                         {
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/manualGAMGProcAgglomeration/manualGAMGProcAgglomeration.H b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/manualGAMGProcAgglomeration/manualGAMGProcAgglomeration.H
index c7efec92f0..1c93711dd3 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/manualGAMGProcAgglomeration/manualGAMGProcAgglomeration.H
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/manualGAMGProcAgglomeration/manualGAMGProcAgglomeration.H
@@ -30,19 +30,23 @@ Description
     In the GAMG control dictionary:
 
         processorAgglomerator manual;
+        // List of level+procagglomeration where
+        // procagglomeration is a set of labelLists. Each labelList is
+        // a cluster of processor which gets combined onto the first element
+        // in the list.
         processorAgglomeration
         (
             (
                 3           //at level 3
                 (
-                    ((0 1) 0)   //coarse 0 from 0,1 (and moved onto 0)
-                    ((2 3) 3)   //coarse 1 from 2,3 (and moved onto 3)
+                    (0 1)   //coarse 0 from 0,1 (and moved onto 0)
+                    (3 2)   //coarse 1 from 2,3 (and moved onto 3)
                 )
             )
             (
                 6           //at level6
                 (
-                    ((0 1) 0)   //coarse 0 from 0,1 (and moved onto 0)
+                    (0 1)   //coarse 0 from 0,1 (and moved onto 0)
                 )
             )
         );
@@ -76,10 +80,8 @@ class manualGAMGProcAgglomeration
 {
     // Private data
 
-        typedef Tuple2<labelList, label> clusterAndMaster;
-
         //- Per level the agglomeration map
-        const List<Tuple2<label, List<clusterAndMaster> > > procAgglomMaps_;
+        const List<Tuple2<label, List<labelList> > > procAgglomMaps_;
 
         //- Any allocated communicators
         DynamicList<label> comms_;
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/procFacesGAMGProcAgglomeration/procFacesGAMGProcAgglomeration.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/procFacesGAMGProcAgglomeration/procFacesGAMGProcAgglomeration.C
new file mode 100644
index 0000000000..6e13d52c99
--- /dev/null
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/procFacesGAMGProcAgglomeration/procFacesGAMGProcAgglomeration.C
@@ -0,0 +1,333 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 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 "procFacesGAMGProcAgglomeration.H"
+#include "addToRunTimeSelectionTable.H"
+#include "GAMGAgglomeration.H"
+#include "Random.H"
+#include "lduMesh.H"
+#include "processorLduInterface.H"
+#include "processorGAMGInterface.H"
+#include "pairGAMGAgglomeration.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(procFacesGAMGProcAgglomeration, 0);
+
+    addToRunTimeSelectionTable
+    (
+        GAMGProcAgglomeration,
+        procFacesGAMGProcAgglomeration,
+        GAMGAgglomeration
+    );
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+// Create single cell mesh
+Foam::autoPtr<Foam::lduPrimitiveMesh>
+Foam::procFacesGAMGProcAgglomeration::singleCellMesh
+(
+    const label singleCellMeshComm,
+    const lduMesh& mesh,
+    scalarField& faceWeights
+) const
+{
+    // Count number of faces per processor
+    List<Map<label> > procFaces(UPstream::nProcs(mesh.comm()));
+    Map<label>& myNeighbours = procFaces[UPstream::myProcNo(mesh.comm())];
+
+    {
+        const lduInterfacePtrsList interfaces(mesh.interfaces());
+        forAll(interfaces, intI)
+        {
+            if (interfaces.set(intI))
+            {
+                const processorLduInterface& pp =
+                    refCast<const processorLduInterface>
+                    (
+                        interfaces[intI]
+                    );
+                label size = interfaces[intI].faceCells().size();
+                myNeighbours.insert(pp.neighbProcNo(), size);
+            }
+        }
+    }
+
+    Pstream::gatherList(procFaces, Pstream::msgType(), mesh.comm());
+    Pstream::scatterList(procFaces, Pstream::msgType(), mesh.comm());
+
+    autoPtr<lduPrimitiveMesh> singleCellMeshPtr;
+
+    if (Pstream::master(mesh.comm()))
+    {
+        // I am master
+        label nCells = Pstream::nProcs(mesh.comm());
+
+        DynamicList<label> l(3*nCells);
+        DynamicList<label> u(3*nCells);
+        DynamicList<scalar> weight(3*nCells);
+
+        DynamicList<label> nbrs;
+        DynamicList<scalar> weights;
+
+        forAll(procFaces, procI)
+        {
+            const Map<label>& neighbours = procFaces[procI];
+
+            // Add all the higher processors
+            nbrs.clear();
+            weights.clear();
+            forAllConstIter(Map<label>, neighbours, iter)
+            {
+                if (iter.key() > procI)
+                {
+                    nbrs.append(iter.key());
+                    weights.append(iter());
+                }
+                sort(nbrs);
+                forAll(nbrs, i)
+                {
+                    l.append(procI);
+                    u.append(nbrs[i]);
+                    weight.append(weights[i]);
+                }
+            }
+        }
+
+        faceWeights.transfer(weight);
+
+        PtrList<const lduInterface> primitiveInterfaces(0);
+        const lduSchedule ps(0);
+
+        singleCellMeshPtr.reset
+        (
+            new lduPrimitiveMesh
+            (
+                nCells,
+                l,
+                u,
+                primitiveInterfaces,
+                ps,
+                singleCellMeshComm
+            )
+        );
+    }
+    return singleCellMeshPtr;
+}
+
+
+Foam::tmp<Foam::labelField>
+Foam::procFacesGAMGProcAgglomeration::processorAgglomeration
+(
+    const lduMesh& mesh
+) const
+{
+    label singleCellMeshComm = UPstream::allocateCommunicator
+    (
+        mesh.comm(),
+        labelList(1, 0)            // only processor 0
+    );
+
+    scalarField faceWeights;
+    autoPtr<lduPrimitiveMesh> singleCellMeshPtr
+    (
+        singleCellMesh
+        (
+            singleCellMeshComm,
+            mesh,
+            faceWeights
+        )
+    );
+
+    tmp<labelField> tfineToCoarse(new labelField(0));
+    labelField& fineToCoarse = tfineToCoarse();
+
+    if (singleCellMeshPtr.valid())
+    {
+        // On master call the agglomerator
+        const lduPrimitiveMesh& singleCellMesh = singleCellMeshPtr();
+
+        label nCoarseProcs;
+        fineToCoarse = pairGAMGAgglomeration::agglomerate
+        (
+            nCoarseProcs,
+            singleCellMesh,
+            faceWeights
+        );
+
+        labelList coarseToMaster(nCoarseProcs, labelMax);
+        forAll(fineToCoarse, cellI)
+        {
+            label coarseI = fineToCoarse[cellI];
+            coarseToMaster[coarseI] = min(coarseToMaster[coarseI], cellI);
+        }
+
+        // Sort according to master and redo restriction
+        labelList newToOld;
+        sortedOrder(coarseToMaster, newToOld);
+        labelList oldToNew(invert(newToOld.size(), newToOld));
+
+        fineToCoarse = UIndirectList<label>(oldToNew, fineToCoarse)();
+    }
+
+    Pstream::scatter(fineToCoarse, Pstream::msgType(), mesh.comm());
+    UPstream::freeCommunicator(singleCellMeshComm);
+
+    return tfineToCoarse;
+}
+
+
+bool Foam::procFacesGAMGProcAgglomeration::doProcessorAgglomeration
+(
+    const lduMesh& mesh
+) const
+{
+    // Check the need for further agglomeration on all processors
+    bool doAgg = mesh.lduAddr().size() < nAgglomeratingCells_;
+    mesh.reduce(doAgg, orOp<bool>());
+    return doAgg;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::procFacesGAMGProcAgglomeration::procFacesGAMGProcAgglomeration
+(
+    GAMGAgglomeration& agglom,
+    const dictionary& controlDict
+)
+:
+    GAMGProcAgglomeration(agglom, controlDict),
+    nAgglomeratingCells_(readLabel(controlDict.lookup("nAgglomeratingCells")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::procFacesGAMGProcAgglomeration::~procFacesGAMGProcAgglomeration()
+{
+    forAllReverse(comms_, i)
+    {
+        if (comms_[i] != -1)
+        {
+            UPstream::freeCommunicator(comms_[i]);
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+bool Foam::procFacesGAMGProcAgglomeration::agglomerate()
+{
+    if (debug)
+    {
+        Pout<< nl << "Starting mesh overview" << endl;
+        printStats(Pout, agglom_);
+    }
+
+    if (agglom_.size() >= 1)
+    {
+        Random rndGen(0);
+
+        for
+        (
+            label fineLevelIndex = 2;
+            fineLevelIndex < agglom_.size();
+            fineLevelIndex++
+        )
+        {
+            if (agglom_.hasMeshLevel(fineLevelIndex))
+            {
+                // Get the fine mesh
+                const lduMesh& levelMesh = agglom_.meshLevel(fineLevelIndex);
+
+                label levelComm = levelMesh.comm();
+                label nProcs = UPstream::nProcs(levelComm);
+
+                if (nProcs > 1 && doProcessorAgglomeration(levelMesh))
+                {
+                    tmp<labelField> tprocAgglomMap
+                    (
+                        processorAgglomeration(levelMesh)
+                    );
+                    const labelField& procAgglomMap = tprocAgglomMap();
+
+                    // Master processor
+                    labelList masterProcs;
+                    // Local processors that agglomerate. agglomProcIDs[0] is in
+                    // masterProc.
+                    List<int> agglomProcIDs;
+                    GAMGAgglomeration::calculateRegionMaster
+                    (
+                        levelComm,
+                        procAgglomMap,
+                        masterProcs,
+                        agglomProcIDs
+                    );
+
+                    // Allocate a communicator for the processor-agglomerated
+                    // matrix
+                    comms_.append
+                    (
+                        UPstream::allocateCommunicator
+                        (
+                            levelComm,
+                            masterProcs
+                        )
+                    );
+
+
+                    // Use procesor agglomeration maps to do the actual
+                    // collecting.
+                    GAMGProcAgglomeration::agglomerate
+                    (
+                        fineLevelIndex,
+                        procAgglomMap,
+                        masterProcs,
+                        agglomProcIDs,
+                        comms_.last()
+                    );
+                }
+            }
+        }
+    }
+
+    // Print a bit
+    if (debug)
+    {
+        Pout<< nl << "Agglomerated mesh overview" << endl;
+        printStats(Pout, agglom_);
+    }
+
+    return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/procFacesGAMGProcAgglomeration/procFacesGAMGProcAgglomeration.H b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/procFacesGAMGProcAgglomeration/procFacesGAMGProcAgglomeration.H
new file mode 100644
index 0000000000..5c574c4158
--- /dev/null
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/procFacesGAMGProcAgglomeration/procFacesGAMGProcAgglomeration.H
@@ -0,0 +1,136 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 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::procFacesGAMGProcAgglomeration
+
+Description
+    Processor agglomeration of GAMGAgglomerations. Needs nAgglomeratingCells
+    which is when to start agglomerating processors. Processors get agglomerated
+    by constructing a single cell mesh for each processor with each
+    processor interface a face. This then gets agglomerated using
+    the pairGAMGAgglomeration algorithm with the number of faces
+    on the original processor interface as face weight.
+
+SourceFiles
+    procFacesGAMGProcAgglomeration.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef procFacesGAMGProcAgglomeration_H
+#define procFacesGAMGProcAgglomeration_H
+
+#include "GAMGProcAgglomeration.H"
+#include "DynamicList.H"
+#include "labelField.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class GAMGAgglomeration;
+class lduMesh;
+class lduPrimitiveMesh;
+
+/*---------------------------------------------------------------------------*\
+                    Class procFacesGAMGProcAgglomeration Declaration
+\*---------------------------------------------------------------------------*/
+
+class procFacesGAMGProcAgglomeration
+:
+    public GAMGProcAgglomeration
+{
+    // Private data
+
+        //- When to processor agglomerate
+        const label nAgglomeratingCells_;
+
+        //- Allocated communicators
+        DynamicList<label> comms_;
+
+    // Private Member Functions
+
+        //- Return (on master) all single-cell meshes collected. single-cell
+        //  meshes are just one cell with all proc faces intact.
+        autoPtr<lduPrimitiveMesh> singleCellMesh
+        (
+            const label singleCellMeshComm,
+            const lduMesh& mesh,
+            scalarField& faceWeights
+        ) const;
+
+        //- Construct processor agglomeration: for every processor the
+        //  coarse processor-cluster it agglomerates onto
+        tmp<labelField> processorAgglomeration(const lduMesh&) const;
+
+        //- Do we need to agglomerate across processors?
+        bool doProcessorAgglomeration(const lduMesh&) const;
+
+        //- Disallow default bitwise copy construct
+        procFacesGAMGProcAgglomeration
+        (
+            const procFacesGAMGProcAgglomeration&
+        );
+
+        //- Disallow default bitwise assignment
+        void operator=(const procFacesGAMGProcAgglomeration&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("procFaces");
+
+
+    // Constructors
+
+        //- Construct given agglomerator and controls
+        procFacesGAMGProcAgglomeration
+        (
+            GAMGAgglomeration& agglom,
+            const dictionary& controlDict
+        );
+
+
+    //- Destructor
+    virtual ~procFacesGAMGProcAgglomeration();
+
+
+    // Member Functions
+
+       //- Modify agglomeration. Return true if modified
+        virtual bool agglomerate();
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
index ae2cfa21ad..9677a17a9b 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
@@ -114,13 +114,6 @@ Foam::solverPerformance Foam::GAMGSolver::solve
                 cmpt
             );
 
-
-            //Pout<< "finestCorrection:" << finestCorrection << endl;
-            //Pout<< "finestResidual:" << finestResidual << endl;
-            //Pout<< "psi:" << psi << endl;
-            //Pout<< "Apsi:" << Apsi << endl;
-
-
             // Calculate finest level residual field
             matrix_.Amul(Apsi, psi, interfaceBouCoeffs_, interfaces_, cmpt);
             finestResidual = source;
@@ -205,7 +198,6 @@ void Foam::GAMGSolver::Vcycle
                     scratch1,
                     coarseCorrFields[leveli].size()
                 );
-                //scalarField ACf(coarseCorrFields[leveli].size(), VGREAT);
 
                 // Scale coarse-grid correction field
                 // but not on the coarsest level because it evaluates to 1
@@ -218,8 +210,6 @@ void Foam::GAMGSolver::Vcycle
                         (
                             ACf.operator const scalarField&()
                         ),
-                        //ACf,
-
                         matrixLevels_[leveli],
                         interfaceLevelsBouCoeffs_[leveli],
                         interfaceLevels_[leveli],
@@ -235,7 +225,6 @@ void Foam::GAMGSolver::Vcycle
                     (
                         ACf.operator const scalarField&()
                     ),
-                    //ACf,
                     coarseCorrFields[leveli],
                     interfaceLevelsBouCoeffs_[leveli],
                     interfaceLevels_[leveli],
@@ -294,11 +283,6 @@ void Foam::GAMGSolver::Vcycle
                 scratch2,
                 coarseCorrFields[leveli].size()
             );
-            //scalarField preSmoothedCoarseCorrField
-            //(
-            //    coarseCorrFields[leveli].size(),
-            //    VGREAT
-            //);
 
             // Only store the preSmoothedCoarseCorrField if pre-smoothing is
             // used
@@ -328,12 +312,6 @@ void Foam::GAMGSolver::Vcycle
             );
             scalarField& ACfRef =
                 const_cast<scalarField&>(ACf.operator const scalarField&());
-            //scalarField ACfRef
-            //(
-            //    coarseCorrFields[leveli].size(),
-            //    VGREAT
-            //);
-
 
             if (interpolateCorrection_)
             {
diff --git a/src/OpenFOAM/meshes/MeshObject/meshObject.C b/src/OpenFOAM/meshes/MeshObject/meshObject.C
index 7130f38b5d..a616f777b0 100644
--- a/src/OpenFOAM/meshes/MeshObject/meshObject.C
+++ b/src/OpenFOAM/meshes/MeshObject/meshObject.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
diff --git a/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.C b/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.C
index 070edd4b69..507fc425ef 100644
--- a/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.C
+++ b/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.C
@@ -254,7 +254,7 @@ Foam::lduPrimitiveMesh::lduPrimitiveMesh
     const label nCells,
     labelList& l,
     labelList& u,
-    const Xfer<PtrList<const lduInterface> >& primitiveInterfaces,
+    PtrList<const lduInterface>& primitiveInterfaces,
     const lduSchedule& ps,
     const label comm
 )
@@ -262,10 +262,12 @@ Foam::lduPrimitiveMesh::lduPrimitiveMesh
     lduAddressing(nCells),
     lowerAddr_(l, true),
     upperAddr_(u, true),
-    primitiveInterfaces_(primitiveInterfaces),
+    primitiveInterfaces_(0),
     patchSchedule_(ps),
     comm_(comm)
 {
+    primitiveInterfaces_.transfer(primitiveInterfaces);
+
     // Create interfaces
     interfaces_.setSize(primitiveInterfaces_.size());
     forAll(primitiveInterfaces_, i)
diff --git a/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.H b/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.H
index 66c26dc37e..5bd7901c8b 100644
--- a/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.H
+++ b/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.H
@@ -133,7 +133,7 @@ public:
             const label nCells,
             labelList& l,
             labelList& u,
-            const Xfer<PtrList<const lduInterface> >& primitiveInterfaces,
+            PtrList<const lduInterface>& primitiveInterfaces,
             const lduSchedule& ps,
             const label comm
         );
diff --git a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndex.H b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndex.H
index 80263b91d3..1aa39c1754 100644
--- a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndex.H
+++ b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndex.H
@@ -160,7 +160,8 @@ public:
                 const labelList& procIDs,
                 const UList<Type>& fld,
                 List<Type>& allFld,
-                const int tag = UPstream::msgType()
+                const int tag = UPstream::msgType(),
+                const Pstream::commsTypes commsType=Pstream::nonBlocking
             );
 
             //- Collect data in processor order on master (== procIDs[0]).
@@ -172,10 +173,11 @@ public:
                 const labelList& procIDs,
                 const UList<Type>& fld,
                 List<Type>& allFld,
-                const int tag = UPstream::msgType()
+                const int tag = UPstream::msgType(),
+                const Pstream::commsTypes commsType=Pstream::nonBlocking
             ) const
             {
-                gather(offsets_, comm, procIDs, fld, allFld, tag);
+                gather(offsets_, comm, procIDs, fld, allFld, tag, commsType);
             }
 
             //- Inplace collect data in processor order on master
@@ -187,7 +189,8 @@ public:
                 const label comm,
                 const labelList& procIDs,
                 List<Type>& fld,
-                const int tag = UPstream::msgType()
+                const int tag = UPstream::msgType(),
+                const Pstream::commsTypes commsType=Pstream::nonBlocking
             );
 
             //- Inplace collect data in processor order on master
@@ -198,10 +201,11 @@ public:
                 const label comm,
                 const labelList& procIDs,
                 List<Type>& fld,
-                const int tag = UPstream::msgType()
+                const int tag = UPstream::msgType(),
+                const Pstream::commsTypes commsType=Pstream::nonBlocking
             ) const
             {
-                gather(offsets_, comm, procIDs, fld, tag);
+                gather(offsets_, comm, procIDs, fld, tag, commsType);
             }
 
             //- Distribute data in processor order. Requires fld to be sized!
@@ -213,7 +217,8 @@ public:
                 const labelList& procIDs,
                 const UList<Type>& allFld,
                 UList<Type>& fld,
-                const int tag = UPstream::msgType()
+                const int tag = UPstream::msgType(),
+                const Pstream::commsTypes commsType=Pstream::nonBlocking
             );
 
             //- Distribute data in processor order. Requires fld to be sized!
@@ -224,10 +229,11 @@ public:
                 const labelList& procIDs,
                 const UList<Type>& allFld,
                 UList<Type>& fld,
-                const int tag = UPstream::msgType()
+                const int tag = UPstream::msgType(),
+                const Pstream::commsTypes commsType=Pstream::nonBlocking
             ) const
             {
-                scatter(offsets_, comm, procIDs, allFld, fld, tag);
+                scatter(offsets_, comm, procIDs, allFld, fld, tag, commsType);
             }
 
 
diff --git a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndexTemplates.C b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndexTemplates.C
index 73d2ea25a8..08435991ee 100644
--- a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndexTemplates.C
+++ b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndexTemplates.C
@@ -30,30 +30,75 @@ License
 template<class Type>
 void Foam::globalIndex::gather
 (
-    const labelUList& offsets,
+    const labelUList& off,
     const label comm,
     const labelList& procIDs,
     const UList<Type>& fld,
     List<Type>& allFld,
-    const int tag
+    const int tag,
+    const Pstream::commsTypes commsType
 )
 {
     if (Pstream::myProcNo(comm) == procIDs[0])
     {
-        allFld.setSize(offsets.last());
+        allFld.setSize(off.last());
 
         // Assign my local data
         SubList<Type>(allFld, fld.size(), 0).assign(fld);
 
-        for (label i = 1; i < procIDs.size(); i++)
+        if (commsType == Pstream::scheduled || commsType == Pstream::blocking)
         {
-            SubList<Type> procSlot(allFld, offsets[i+1]-offsets[i], offsets[i]);
-
-            if (contiguous<Type>())
+            for (label i = 1; i < procIDs.size(); i++)
             {
+                SubList<Type> procSlot(allFld, off[i+1]-off[i], off[i]);
+
+                if (contiguous<Type>())
+                {
+                    IPstream::read
+                    (
+                        commsType,
+                        procIDs[i],
+                        reinterpret_cast<char*>(procSlot.begin()),
+                        procSlot.byteSize(),
+                        tag,
+                        comm
+                    );
+                }
+                else
+                {
+                    IPstream fromSlave
+                    (
+                        commsType,
+                        procIDs[i],
+                        0,
+                        tag,
+                        comm
+                    );
+                    fromSlave >> procSlot;
+                }
+            }
+        }
+        else
+        {
+            // nonBlocking
+
+            if (!contiguous<Type>())
+            {
+                FatalErrorIn("globalIndex::gather(..)")
+                    << "nonBlocking not supported for non-contiguous data"
+                    << exit(FatalError);
+            }
+
+            label startOfRequests = Pstream::nRequests();
+
+            // Set up reads
+            for (label i = 1; i < procIDs.size(); i++)
+            {
+                SubList<Type> procSlot(allFld, off[i+1]-off[i], off[i]);
+
                 IPstream::read
                 (
-                    Pstream::scheduled,
+                    commsType,
                     procIDs[i],
                     reinterpret_cast<char*>(procSlot.begin()),
                     procSlot.byteSize(),
@@ -61,45 +106,66 @@ void Foam::globalIndex::gather
                     comm
                 );
             }
-            else
-            {
-                IPstream fromSlave
-                (
-                    Pstream::scheduled,
-                    procIDs[i],
-                    0,
-                    tag,
-                    comm
-                );
-                fromSlave >> procSlot;
-            }
+
+            // Wait for all to finish
+            Pstream::waitRequests(startOfRequests);
         }
     }
     else
     {
-        if (contiguous<Type>())
+        if (commsType == Pstream::scheduled || commsType == Pstream::blocking)
         {
+            if (contiguous<Type>())
+            {
+                OPstream::write
+                (
+                    commsType,
+                    procIDs[0],
+                    reinterpret_cast<const char*>(fld.begin()),
+                    fld.byteSize(),
+                    tag,
+                    comm
+                );
+            }
+            else
+            {
+                OPstream toMaster
+                (
+                    commsType,
+                    procIDs[0],
+                    0,
+                    tag,
+                    comm
+                );
+                toMaster << fld;
+            }
+        }
+        else
+        {
+            // nonBlocking
+
+            if (!contiguous<Type>())
+            {
+                FatalErrorIn("globalIndex::gather(..)")
+                    << "nonBlocking not supported for non-contiguous data"
+                    << exit(FatalError);
+            }
+
+            label startOfRequests = Pstream::nRequests();
+
+            // Set up write
             OPstream::write
             (
-                Pstream::scheduled,
+                commsType,
                 procIDs[0],
                 reinterpret_cast<const char*>(fld.begin()),
                 fld.byteSize(),
                 tag,
                 comm
             );
-        }
-        else
-        {
-            OPstream toMaster
-            (
-                Pstream::scheduled,
-                procIDs[0],
-                0,
-                tag,
-                comm
-            );
-            toMaster << fld;
+
+            // Wait for all to finish
+            Pstream::waitRequests(startOfRequests);
         }
     }
 }
@@ -108,111 +174,104 @@ void Foam::globalIndex::gather
 template<class Type>
 void Foam::globalIndex::gather
 (
-    const labelUList& offsets,
+    const labelUList& off,
     const label comm,
     const labelList& procIDs,
     List<Type>& fld,
-    const int tag
+    const int tag,
+    const Pstream::commsTypes commsType
 )
 {
+    List<Type> allFld;
+
+    gather(off, comm, procIDs, fld, allFld, tag, commsType);
+
     if (Pstream::myProcNo(comm) == procIDs[0])
     {
-        List<Type> allFld(offsets.last());
-
-        // Assign my local data
-        SubList<Type>(allFld, fld.size(), 0).assign(fld);
-
-        for (label i = 1; i < procIDs.size(); i++)
-        {
-            SubList<Type> procSlot(allFld, offsets[i+1]-offsets[i], offsets[i]);
-
-            if (contiguous<Type>())
-            {
-                IPstream::read
-                (
-                    Pstream::scheduled,
-                    procIDs[i],
-                    reinterpret_cast<char*>(procSlot.begin()),
-                    procSlot.byteSize(),
-                    tag,
-                    comm
-                );
-            }
-            else
-            {
-                IPstream fromSlave
-                (
-                    Pstream::scheduled,
-                    procIDs[i],
-                    0,
-                    tag,
-                    comm
-                );
-                fromSlave >> procSlot;
-            }
-        }
-
         fld.transfer(allFld);
     }
-    else
-    {
-        if (contiguous<Type>())
-        {
-            OPstream::write
-            (
-                Pstream::scheduled,
-                procIDs[0],
-                reinterpret_cast<const char*>(fld.begin()),
-                fld.byteSize(),
-                tag,
-                comm
-            );
-        }
-        else
-        {
-            OPstream toMaster
-            (
-                Pstream::scheduled,
-                procIDs[0],
-                0,
-                tag,
-                comm
-            );
-            toMaster << fld;
-        }
-    }
 }
 
 
 template<class Type>
 void Foam::globalIndex::scatter
 (
-    const labelUList& offsets,
+    const labelUList& off,
     const label comm,
     const labelList& procIDs,
     const UList<Type>& allFld,
     UList<Type>& fld,
-    const int tag
+    const int tag,
+    const Pstream::commsTypes commsType
 )
 {
     if (Pstream::myProcNo(comm) == procIDs[0])
     {
-        fld.assign(SubList<Type>(allFld, offsets[1]-offsets[0]));
+        fld.assign(SubList<Type>(allFld, off[1]-off[0]));
 
-        for (label i = 1; i < procIDs.size(); i++)
+        if (commsType == Pstream::scheduled || commsType == Pstream::blocking)
         {
-            const SubList<Type> procSlot
-            (
-                allFld,
-                offsets[i+1]-offsets[i],
-                offsets[i]
-            );
-
-            if (contiguous<Type>())
+            for (label i = 1; i < procIDs.size(); i++)
             {
+                const SubList<Type> procSlot
+                (
+                    allFld,
+                    off[i+1]-off[i],
+                    off[i]
+                );
+
+                if (contiguous<Type>())
+                {
+                    OPstream::write
+                    (
+                        commsType,
+                        procIDs[i],
+                        reinterpret_cast<const char*>(procSlot.begin()),
+                        procSlot.byteSize(),
+                        tag,
+                        comm
+                    );
+                }
+                else
+                {
+                    OPstream toSlave
+                    (
+                        commsType,
+                        procIDs[i],
+                        0,
+                        tag,
+                        comm
+                    );
+                    toSlave << procSlot;
+                }
+            }
+        }
+        else
+        {
+            // nonBlocking
+
+            if (!contiguous<Type>())
+            {
+                FatalErrorIn("globalIndex::scatter(..)")
+                    << "nonBlocking not supported for non-contiguous data"
+                    << exit(FatalError);
+            }
+
+            label startOfRequests = Pstream::nRequests();
+
+            // Set up writes
+            for (label i = 1; i < procIDs.size(); i++)
+            {
+                const SubList<Type> procSlot
+                (
+                    allFld,
+                    off[i+1]-off[i],
+                    off[i]
+                );
+
                 OPstream::write
                 (
-                    Pstream::scheduled,
+                    commsType,
                     procIDs[i],
                     reinterpret_cast<const char*>(procSlot.begin()),
                     procSlot.byteSize(),
@@ -220,45 +279,66 @@ void Foam::globalIndex::scatter
                     comm
                 );
             }
-            else
-            {
-                OPstream toSlave
-                (
-                    Pstream::scheduled,
-                    procIDs[i],
-                    0,
-                    tag,
-                    comm
-                );
-                toSlave << procSlot;
-            }
+
+            // Wait for all to finish
+            Pstream::waitRequests(startOfRequests);
         }
     }
     else
     {
-        if (contiguous<Type>())
+        if (commsType == Pstream::scheduled || commsType == Pstream::blocking)
         {
+            if (contiguous<Type>())
+            {
+                IPstream::read
+                (
+                    commsType,
+                    procIDs[0],
+                    reinterpret_cast<char*>(fld.begin()),
+                    fld.byteSize(),
+                    tag,
+                    comm
+                );
+            }
+            else
+            {
+                IPstream fromMaster
+                (
+                    commsType,
+                    procIDs[0],
+                    0,
+                    tag,
+                    comm
+                );
+                fromMaster >> fld;
+            }
+        }
+        else
+        {
+            // nonBlocking
+
+            if (!contiguous<Type>())
+            {
+                FatalErrorIn("globalIndex::scatter(..)")
+                    << "nonBlocking not supported for non-contiguous data"
+                    << exit(FatalError);
+            }
+
+            label startOfRequests = Pstream::nRequests();
+
+            // Set up read
             IPstream::read
             (
-                Pstream::scheduled,
+                commsType,
                 procIDs[0],
                 reinterpret_cast<char*>(fld.begin()),
                 fld.byteSize(),
                 tag,
                 comm
             );
-        }
-        else
-        {
-            IPstream fromMaster
-            (
-                Pstream::scheduled,
-                procIDs[0],
-                0,
-                tag,
-                comm
-            );
-            fromMaster >> fld;
+
+            // Wait for all to finish
+            Pstream::waitRequests(startOfRequests);
         }
     }
 }
diff --git a/src/OpenFOAM/primitives/strings/word/wordIOList.C b/src/OpenFOAM/primitives/strings/word/wordIOList.C
index ede987af53..58456cb0ac 100644
--- a/src/OpenFOAM/primitives/strings/word/wordIOList.C
+++ b/src/OpenFOAM/primitives/strings/word/wordIOList.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -38,12 +38,17 @@ namespace Foam
 }
 
 
-void Foam::printTable(const List<wordList>& wll, Ostream& os)
+void Foam::printTable
+(
+    const List<wordList>& wll,
+    List<string::size_type>& columnWidth,
+    Ostream& os
+)
 {
     if (!wll.size()) return;
 
     // Find the maximum word length for each column
-    List<string::size_type> columnWidth(wll[0].size(), string::size_type(0));
+    columnWidth.setSize(wll[0].size(), string::size_type(0));
     forAll(columnWidth, j)
     {
         forAll(wll, i)
@@ -75,4 +80,11 @@ void Foam::printTable(const List<wordList>& wll, Ostream& os)
 }
 
 
+void Foam::printTable(const List<wordList>& wll, Ostream& os)
+{
+    List<string::size_type> columnWidth;
+    printTable(wll, columnWidth, os);
+}
+
+
 // ************************************************************************* //
diff --git a/src/OpenFOAM/primitives/strings/word/wordIOList.H b/src/OpenFOAM/primitives/strings/word/wordIOList.H
index c560a33a95..da9d5aa426 100644
--- a/src/OpenFOAM/primitives/strings/word/wordIOList.H
+++ b/src/OpenFOAM/primitives/strings/word/wordIOList.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,6 +43,7 @@ namespace Foam
     typedef IOList<wordList> wordListIOList;
 
     // Print word list list as a table
+    void printTable(const List<wordList>&, List<string::size_type>&, Ostream&);
     void printTable(const List<wordList>&, Ostream&);
 }
 
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSup.C b/src/finiteVolume/finiteVolume/fvc/fvcSup.C
index d1af6ad007..67b0e11df7 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcSup.C
+++ b/src/finiteVolume/finiteVolume/fvc/fvcSup.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,7 +43,7 @@ tmp<GeometricField<Type, fvPatchField, volMesh> >
 Su
 (
     const GeometricField<Type, fvPatchField, volMesh>& su,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return su;
@@ -54,7 +54,7 @@ tmp<GeometricField<Type, fvPatchField, volMesh> >
 Su
 (
     const tmp<GeometricField<Type, fvPatchField, volMesh> >& tsu,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return tsu;
@@ -66,7 +66,7 @@ tmp<GeometricField<Type, fvPatchField, volMesh> >
 Sp
 (
     const volScalarField& sp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return sp*vf;
@@ -77,7 +77,7 @@ tmp<GeometricField<Type, fvPatchField, volMesh> >
 Sp
 (
     const tmp<volScalarField>& tsp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return tsp*vf;
@@ -89,7 +89,7 @@ tmp<GeometricField<Type, fvPatchField, volMesh> >
 Sp
 (
     const dimensionedScalar& sp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return sp*vf;
@@ -101,7 +101,7 @@ tmp<GeometricField<Type, fvPatchField, volMesh> >
 SuSp
 (
     const volScalarField& sp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return sp*vf;
@@ -112,7 +112,7 @@ tmp<GeometricField<Type, fvPatchField, volMesh> >
 SuSp
 (
     const tmp<volScalarField>& tsp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return tsp*vf;
diff --git a/src/finiteVolume/finiteVolume/fvc/fvcSup.H b/src/finiteVolume/finiteVolume/fvc/fvcSup.H
index 8ac75a2150..2958d676b9 100644
--- a/src/finiteVolume/finiteVolume/fvc/fvcSup.H
+++ b/src/finiteVolume/finiteVolume/fvc/fvcSup.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -55,14 +55,14 @@ namespace fvc
         tmp<GeometricField<Type, fvPatchField, volMesh> > Su
         (
             const GeometricField<Type, fvPatchField, volMesh>&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<GeometricField<Type, fvPatchField, volMesh> > Su
         (
             const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
 
@@ -72,14 +72,14 @@ namespace fvc
         tmp<GeometricField<Type, fvPatchField, volMesh> > Sp
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<GeometricField<Type, fvPatchField, volMesh> > Sp
         (
             const tmp<volScalarField>&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
 
@@ -87,7 +87,7 @@ namespace fvc
         tmp<GeometricField<Type, fvPatchField, volMesh> > Sp
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
 
@@ -97,14 +97,14 @@ namespace fvc
         tmp<GeometricField<Type, fvPatchField, volMesh> > SuSp
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<GeometricField<Type, fvPatchField, volMesh> > SuSp
         (
             const tmp<volScalarField>&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 }
 
diff --git a/src/finiteVolume/fvMatrices/solvers/GAMGSymSolver/GAMGAgglomerations/faceAreaPairGAMGAgglomeration/faceAreaPairGAMGAgglomeration.C b/src/finiteVolume/fvMatrices/solvers/GAMGSymSolver/GAMGAgglomerations/faceAreaPairGAMGAgglomeration/faceAreaPairGAMGAgglomeration.C
index 9bfcda020c..42a35f5c1b 100644
--- a/src/finiteVolume/fvMatrices/solvers/GAMGSymSolver/GAMGAgglomerations/faceAreaPairGAMGAgglomeration/faceAreaPairGAMGAgglomeration.C
+++ b/src/finiteVolume/fvMatrices/solvers/GAMGSymSolver/GAMGAgglomerations/faceAreaPairGAMGAgglomeration/faceAreaPairGAMGAgglomeration.C
@@ -90,8 +90,6 @@ Foam::faceAreaPairGAMGAgglomeration::faceAreaPairGAMGAgglomeration
 :
     pairGAMGAgglomeration(mesh, controlDict)
 {
-    vectorField n(faceAreas/mag(faceAreas));
-
     //agglomerate(mesh, sqrt(mag(faceAreas)));
     agglomerate
     (
@@ -100,7 +98,8 @@ Foam::faceAreaPairGAMGAgglomeration::faceAreaPairGAMGAgglomeration
         (
             cmptMultiply
             (
-                n,
+                faceAreas
+               /sqrt(mag(faceAreas)),
                 vector(1, 1.01, 1.02)
                 //vector::one
             )
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.C b/src/meshTools/searchableSurface/triSurfaceMesh.C
index 2222981793..6009ec6670 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.C
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.C
@@ -433,24 +433,28 @@ Foam::triSurfaceMesh::edgeTree() const
           + nInternalEdges()
         );
 
-        treeBoundBox bb;
-        label nPoints;
-        PatchTools::calcBounds
-        (
-            static_cast<const triSurface&>(*this),
-            bb,
-            nPoints
-        );
+        treeBoundBox bb(vector::zero, vector::zero);
 
-        // Random number generator. Bit dodgy since not exactly random ;-)
-        Random rndGen(65431);
+        if (bEdges.size())
+        {
+            label nPoints;
+            PatchTools::calcBounds
+            (
+                static_cast<const triSurface&>(*this),
+                bb,
+                nPoints
+            );
 
-        // Slightly extended bb. Slightly off-centred just so on symmetric
-        // geometry there are less face/edge aligned items.
+            // Random number generator. Bit dodgy since not exactly random ;-)
+            Random rndGen(65431);
 
-        bb = bb.extend(rndGen, 1e-4);
-        bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
-        bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+            // Slightly extended bb. Slightly off-centred just so on symmetric
+            // geometry there are less face/edge aligned items.
+
+            bb = bb.extend(rndGen, 1e-4);
+            bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+            bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+        }
 
         scalar oldTol = indexedOctree<treeDataEdge>::perturbTol();
         indexedOctree<treeDataEdge>::perturbTol() = tolerance();
diff --git a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.C b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.C
index 738cb8dcdd..34dc89f920 100644
--- a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.C
+++ b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceRegionSearch.C
@@ -121,33 +121,40 @@ Foam::triSurfaceRegionSearch::treeByRegion() const
             );
 
             // Calculate bb without constructing local point numbering.
-            treeBoundBox bb;
-            label nPoints;
-            PatchTools::calcBounds
-            (
-                indirectRegionPatches_[regionI],
-                bb,
-                nPoints
-            );
+            treeBoundBox bb(vector::zero, vector::zero);
 
-//            if (nPoints != surface().points().size())
-//            {
-//                WarningIn("triSurfaceRegionSearch::treeByRegion() const")
-//                    << "Surface does not have compact point numbering. "
-//                    << "Of " << surface().points().size()
-//                    << " only " << nPoints
-//                    << " are used. This might give problems in some routines."
-//                    << endl;
-//            }
+            if (indirectRegionPatches_[regionI].size())
+            {
+                label nPoints;
+                PatchTools::calcBounds
+                (
+                    indirectRegionPatches_[regionI],
+                    bb,
+                    nPoints
+                );
 
-            // Random number generator. Bit dodgy since not exactly random ;-)
-            Random rndGen(65431);
+    //            if (nPoints != surface().points().size())
+    //            {
+    //                WarningIn("triSurfaceRegionSearch::treeByRegion() const")
+    //                    << "Surface does not have compact point numbering. "
+    //                    << "Of " << surface().points().size()
+    //                    << " only " << nPoints
+    //                    << " are used."
+    //                    << " This might give problems in some routines."
+    //                    << endl;
+    //            }
 
-            // Slightly extended bb. Slightly off-centred just so on symmetric
-            // geometry there are fewer face/edge aligned items.
-            bb = bb.extend(rndGen, 1e-4);
-            bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
-            bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+                // Random number generator. Bit dodgy since not exactly
+                // random ;-)
+                Random rndGen(65431);
+
+                // Slightly extended bb. Slightly off-centred just so
+                // on symmetric geometry there are fewer face/edge
+                // aligned items.
+                bb = bb.extend(rndGen, 1e-4);
+                bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+                bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+            }
 
             treeByRegion_.set
             (
diff --git a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.C b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.C
index fbff791040..43877e3393 100644
--- a/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.C
+++ b/src/meshTools/triSurface/triSurfaceSearch/triSurfaceSearch.C
@@ -200,28 +200,33 @@ Foam::triSurfaceSearch::tree() const
     if (treePtr_.empty())
     {
         // Calculate bb without constructing local point numbering.
-        treeBoundBox bb;
-        label nPoints;
-        PatchTools::calcBounds(surface(), bb, nPoints);
+        treeBoundBox bb(vector::zero, vector::zero);
 
-        if (nPoints != surface().points().size())
+        if (surface().size())
         {
-            WarningIn("triSurfaceSearch::tree() const")
-                << "Surface does not have compact point numbering."
-                << " Of " << surface().points().size() << " only " << nPoints
-                << " are used. This might give problems in some routines."
-                << endl;
+            label nPoints;
+            PatchTools::calcBounds(surface(), bb, nPoints);
+
+            if (nPoints != surface().points().size())
+            {
+                WarningIn("triSurfaceSearch::tree() const")
+                    << "Surface does not have compact point numbering."
+                    << " Of " << surface().points().size()
+                    << " only " << nPoints
+                    << " are used. This might give problems in some routines."
+                    << endl;
+            }
+
+            // Random number generator. Bit dodgy since not exactly random ;-)
+            Random rndGen(65431);
+
+            // Slightly extended bb. Slightly off-centred just so on symmetric
+            // geometry there are less face/edge aligned items.
+            bb = bb.extend(rndGen, 1e-4);
+            bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
+            bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
         }
 
-        // Random number generator. Bit dodgy since not exactly random ;-)
-        Random rndGen(65431);
-
-        // Slightly extended bb. Slightly off-centred just so on symmetric
-        // geometry there are less face/edge aligned items.
-        bb = bb.extend(rndGen, 1e-4);
-        bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
-        bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
-
         scalar oldTol = indexedOctree<treeDataTriSurface>::perturbTol();
         indexedOctree<treeDataTriSurface>::perturbTol() = tolerance_;
 
diff --git a/src/parallel/decompose/decompositionMethods/noDecomp/noDecomp.H b/src/parallel/decompose/decompositionMethods/noDecomp/noDecomp.H
index 13368c9970..154838d35c 100644
--- a/src/parallel/decompose/decompositionMethods/noDecomp/noDecomp.H
+++ b/src/parallel/decompose/decompositionMethods/noDecomp/noDecomp.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -78,7 +78,7 @@ public:
         //  up to the user.
         virtual bool parallelAware() const
         {
-            return false;
+            return true;
         }
 
         //- Return for every coordinate the wanted processor number. Use the
@@ -90,12 +90,7 @@ public:
             const scalarField& cWeights
         )
         {
-            notImplemented
-            (
-                "decompose(const polyMesh&, const pointField&"
-                ", const scalarField&)"
-            );
-            return labelList(0);
+            return labelList(cc.size(), Pstream::myProcNo());
         }
 
         //- Return for every coordinate the wanted processor number. Explicitly
@@ -112,12 +107,7 @@ public:
             const scalarField& cWeights
         )
         {
-            notImplemented
-            (
-                "decompose(const labelListList&, const pointField&"
-                ", const scalarField&)"
-            );
-            return labelList(0);
+            return labelList(globalCellCells.size(), Pstream::myProcNo());
         }
 
 };
diff --git a/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.C b/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.C
index 80308fac5e..7c80bd1a56 100644
--- a/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.C
+++ b/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.C
@@ -2360,10 +2360,28 @@ bool Foam::distributedTriSurfaceMesh::writeObject
     // Make sure dictionary goes to same directory as surface
     const_cast<fileName&>(dict_.instance()) = searchableSurface::instance();
 
+    // Copy of triSurfaceMesh::writeObject except for the sorting of
+    // triangles by region. This is done so we preserve region names,
+    // even if locally we have zero triangles.
+    {
+        fileName fullPath(searchableSurface::objectPath());
+
+        if (!mkDir(fullPath.path()))
+        {
+            return false;
+        }
+
+        // Important: preserve any zero-sized patches
+        triSurface::write(fullPath, true);
+
+        if (!isFile(fullPath))
+        {
+            return false;
+        }
+    }
+
     // Dictionary needs to be written in ascii - binary output not supported.
-    return
-        triSurfaceMesh::writeObject(fmt, ver, cmp)
-     && dict_.writeObject(IOstream::ASCII, ver, cmp);
+    return dict_.writeObject(IOstream::ASCII, ver, cmp);
 }
 
 
diff --git a/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.H b/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.H
index 32a062a573..7e4d27ef0c 100644
--- a/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.H
+++ b/src/parallel/distributed/distributedTriSurfaceMesh/distributedTriSurfaceMesh.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -446,6 +446,10 @@ public:
         // regIOobject implementation
 
             //- Write using given format, version and compression
+            //  Do not use the triSurfaceMesh::writeObject since it
+            //  would filter out empty regions. These need to be preserved
+            //  in case we want to make decisions based on the number of
+            //  regions.
             virtual bool writeObject
             (
                 IOstream::streamFormat fmt,
diff --git a/src/surfMesh/surfaceFormats/obj/OBJstream.C b/src/surfMesh/surfaceFormats/obj/OBJstream.C
index ce0eb2721b..46a6b7e0c0 100644
--- a/src/surfMesh/surfaceFormats/obj/OBJstream.C
+++ b/src/surfMesh/surfaceFormats/obj/OBJstream.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -241,6 +241,20 @@ Foam::Ostream& Foam::OBJstream::write(const linePointRef& ln)
 }
 
 
+Foam::Ostream& Foam::OBJstream::write
+(
+    const linePointRef& ln,
+    const vector& n0,
+    const vector& n1
+)
+{
+    write(ln.start(), n0);
+    write(ln.end(), n1);
+    write("l ") << nVertices_-1 << ' ' << nVertices_ << nl;
+    return *this;
+}
+
+
 Foam::Ostream& Foam::OBJstream::write
 (
     const triPointRef& f,
diff --git a/src/surfMesh/surfaceFormats/obj/OBJstream.H b/src/surfMesh/surfaceFormats/obj/OBJstream.H
index ec7e93d24e..ee6ae58143 100644
--- a/src/surfMesh/surfaceFormats/obj/OBJstream.H
+++ b/src/surfMesh/surfaceFormats/obj/OBJstream.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -92,7 +92,7 @@ public:
 
         // Access
 
-            //- Return the name of the stream
+            //- Return the number of vertices written
             label nVertices() const
             {
                 return nVertices_;
@@ -135,6 +135,14 @@ public:
             //- Write line
             Ostream& write(const linePointRef&);
 
+            //- Write line with points and vector normals ('vn')
+            Ostream& write
+            (
+                const linePointRef&,
+                const vector& n0,
+                const vector& n1
+            );
+
             //- Write triangle as points with lines or filled polygon
             Ostream& write(const triPointRef&, const bool lines = true);
 
diff --git a/src/turbulenceModels/derivedFvPatchFields/porousBafflePressure/porousBafflePressureFvPatchField.H b/src/turbulenceModels/derivedFvPatchFields/porousBafflePressure/porousBafflePressureFvPatchField.H
index a52e8e5b34..ed22575785 100644
--- a/src/turbulenceModels/derivedFvPatchFields/porousBafflePressure/porousBafflePressureFvPatchField.H
+++ b/src/turbulenceModels/derivedFvPatchFields/porousBafflePressure/porousBafflePressureFvPatchField.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -34,7 +34,7 @@ Description
     The porous baffle introduces a pressure jump defined by:
 
         \f[
-            \Delta p = -(I \mu U + 0.5 D \rho |U|^2 L)
+            \Delta p = -(I \mu U + 0.5 D \rho |U|^2 )L
         \f]
 
     where
diff --git a/src/turbulenceModels/derivedFvPatchFields/porousBafflePressure/porousBafflePressureFvPatchFields.C b/src/turbulenceModels/derivedFvPatchFields/porousBafflePressure/porousBafflePressureFvPatchFields.C
index 08855fa7fb..32317c0ad1 100644
--- a/src/turbulenceModels/derivedFvPatchFields/porousBafflePressure/porousBafflePressureFvPatchFields.C
+++ b/src/turbulenceModels/derivedFvPatchFields/porousBafflePressure/porousBafflePressureFvPatchFields.C
@@ -78,7 +78,7 @@ void Foam::porousBafflePressureFvPatchField<Foam::scalar>::updateCoeffs()
 
         const scalarField nuEffw = turbModel.nuEff()().boundaryField()[patchI];
 
-        jump_ = -sign(Un)*(I_*nuEffw + D_*0.5*magUn*length_)*magUn;
+        jump_ = -sign(Un)*(I_*nuEffw + D_*0.5*magUn)*magUn*length_;
     }
     else
     {
@@ -95,7 +95,7 @@ void Foam::porousBafflePressureFvPatchField<Foam::scalar>::updateCoeffs()
 
         Un /= rhow;
 
-        jump_ = -sign(Un)*(I_*muEffw + D_*0.5*rhow*magUn*length_)*magUn;
+        jump_ = -sign(Un)*(I_*muEffw + D_*0.5*rhow*magUn)*magUn*length_;
     }
 
     if (debug)
diff --git a/wmake/wmake b/wmake/wmake
index ced70188cc..ebcd38cb6e 100755
--- a/wmake/wmake
+++ b/wmake/wmake
@@ -3,7 +3,7 @@
 # =========                 |
 # \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
 #  \\    /   O peration     |
-#   \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+#   \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
 #    \\/     M anipulation  |
 #-------------------------------------------------------------------------------
 # License
@@ -202,9 +202,9 @@ fi
     exit 1
 }
 
-# transform "all" option to "libso" if that looks appropriate or remove it
+# transform "all" or no option to "libso" if that looks appropriate or remove it
 # so that the call to make builds the application
-if [ "$makeType" = all ]
+if [ "$makeType" = all -o "$makeType" = "" ]
 then
     unset makeType
     if grep -e '^ *LIB *=' "$MakeDir/files" >/dev/null 2>&1