diff --git a/applications/solvers/incompressible/pimpleDyMFoam/createFields.H b/applications/solvers/incompressible/pimpleDyMFoam/createFields.H
index ff3f9c5169..bc798ba988 100644
--- a/applications/solvers/incompressible/pimpleDyMFoam/createFields.H
+++ b/applications/solvers/incompressible/pimpleDyMFoam/createFields.H
@@ -1,4 +1,3 @@
-
     Info<< "Reading field p\n" << endl;
     volScalarField p
     (
@@ -33,7 +32,7 @@
 
     label pRefCell = 0;
     scalar pRefValue = 0.0;
-    setRefCell(p, mesh.solutionDict().subDict("PISO"), pRefCell, pRefValue);
+    setRefCell(p, mesh.solutionDict().subDict("PIMPLE"), pRefCell, pRefValue);
 
     singlePhaseTransportModel laminarTransport(U, phi);
 
diff --git a/applications/solvers/incompressible/pimpleDyMFoam/pimpleDyMFoam.C b/applications/solvers/incompressible/pimpleDyMFoam/pimpleDyMFoam.C
index b705d2bbb7..66692cc4a8 100644
--- a/applications/solvers/incompressible/pimpleDyMFoam/pimpleDyMFoam.C
+++ b/applications/solvers/incompressible/pimpleDyMFoam/pimpleDyMFoam.C
@@ -23,7 +23,7 @@ License
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
 Application
-    turbDyMFoam
+    pimpleDyMFoam.C
 
 Description
     Transient solver for incompressible, flow of Newtonian fluids
@@ -47,7 +47,7 @@ int main(int argc, char *argv[])
 
 #   include "createTime.H"
 #   include "createDynamicFvMesh.H"
-#   include "readPISOControls.H"
+#   include "readPIMPLEControls.H"
 #   include "initContinuityErrs.H"
 #   include "createFields.H"
 #   include "readTimeControls.H"
diff --git a/applications/solvers/incompressible/pimpleDyMFoam/readControls.H b/applications/solvers/incompressible/pimpleDyMFoam/readControls.H
index 9336616a1b..3bd20c5c5c 100644
--- a/applications/solvers/incompressible/pimpleDyMFoam/readControls.H
+++ b/applications/solvers/incompressible/pimpleDyMFoam/readControls.H
@@ -1,14 +1,14 @@
 #   include "readTimeControls.H"
-#   include "readPISOControls.H"
+#   include "readPIMPLEControls.H"
 
     bool correctPhi = false;
-    if (piso.found("correctPhi"))
+    if (pimple.found("correctPhi"))
     {
-        correctPhi = Switch(piso.lookup("correctPhi"));
+        correctPhi = Switch(pimple.lookup("correctPhi"));
     }
 
     bool checkMeshCourantNo = false;
-    if (piso.found("checkMeshCourantNo"))
+    if (pimple.found("checkMeshCourantNo"))
     {
-        checkMeshCourantNo = Switch(piso.lookup("checkMeshCourantNo"));
+        checkMeshCourantNo = Switch(pimple.lookup("checkMeshCourantNo"));
     }
diff --git a/applications/utilities/mesh/manipulation/setsToZones/setsToZones.C b/applications/utilities/mesh/manipulation/setsToZones/setsToZones.C
index 35eae96a58..ac32c1ccc9 100644
--- a/applications/utilities/mesh/manipulation/setsToZones/setsToZones.C
+++ b/applications/utilities/mesh/manipulation/setsToZones/setsToZones.C
@@ -48,6 +48,7 @@ Description
 #include "OFstream.H"
 #include "IFstream.H"
 #include "IOobjectList.H"
+#include "SortableList.H"
 
 using namespace Foam;
 
@@ -107,6 +108,7 @@ int main(int argc, char *argv[])
     {
         // Not in memory. Load it.
         pointSet set(*iter());
+        SortableList<label> pointLabels(set.toc());
 
         label zoneID = mesh.pointZones().findZoneID(set.name());
         if (zoneID == -1)
@@ -120,7 +122,7 @@ int main(int argc, char *argv[])
                 new pointZone
                 (
                     set.name(),             //name
-                    set.toc(),              //addressing
+                    pointLabels,            //addressing
                     sz,                     //index
                     mesh.pointZones()       //pointZoneMesh
                 )
@@ -131,7 +133,7 @@ int main(int argc, char *argv[])
         {
             Info<< "Overwriting contents of existing pointZone " << zoneID
                 << " with that of set " << set.name() << "." << endl;
-            mesh.pointZones()[zoneID] = set.toc();
+            mesh.pointZones()[zoneID] = pointLabels;
             mesh.pointZones().writeOpt() = IOobject::AUTO_WRITE;
         }
     }
@@ -150,6 +152,7 @@ int main(int argc, char *argv[])
     {
         // Not in memory. Load it.
         cellSet set(*iter());
+        SortableList<label> cellLabels(set.toc());
 
         label zoneID = mesh.cellZones().findZoneID(set.name());
         if (zoneID == -1)
@@ -163,7 +166,7 @@ int main(int argc, char *argv[])
                 new cellZone
                 (
                     set.name(),             //name
-                    set.toc(),              //addressing
+                    cellLabels,             //addressing
                     sz,                     //index
                     mesh.cellZones()        //pointZoneMesh
                 )
@@ -174,7 +177,7 @@ int main(int argc, char *argv[])
         {
             Info<< "Overwriting contents of existing cellZone " << zoneID
                 << " with that of set " << set.name() << "." << endl;
-            mesh.cellZones()[zoneID] = set.toc();
+            mesh.cellZones()[zoneID] = cellLabels;
             mesh.cellZones().writeOpt() = IOobject::AUTO_WRITE;
         }
     }
@@ -193,6 +196,7 @@ int main(int argc, char *argv[])
     {
         // Not in memory. Load it.
         faceSet set(*iter());
+        SortableList<label> faceLabels(set.toc());
 
         DynamicList<label> addressing(set.size());
         DynamicList<bool> flipMap(set.size());
@@ -214,9 +218,9 @@ int main(int argc, char *argv[])
             // Load corresponding cells
             cellSet cells(mesh, setName);
 
-            forAllConstIter(faceSet, set, iter)
+            forAll(faceLabels, i)
             {
-                label faceI = iter.key();
+                label faceI = faceLabels[i];
 
                 bool flip = false;
 
@@ -273,9 +277,10 @@ int main(int argc, char *argv[])
         else
         {
             // No flip map.
-            forAllConstIter(faceSet, set, iter)
+            forAll(faceLabels, i)
             {
-                addressing.append(iter.key());
+                label faceI = faceLabels[i];
+                addressing.append(faceI);
                 flipMap.append(false);
             }
         }
diff --git a/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx b/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx
index 0d99f3e95b..47e919c9c3 100644
--- a/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx
+++ b/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx
@@ -50,6 +50,7 @@ vtkPV3FoamReader::vtkPV3FoamReader()
 
     output0_  = NULL;
 
+#ifdef VTKPV3FOAM_DUALPORT
     // Add second output for the Lagrangian
     this->SetNumberOfOutputPorts(2);
     vtkMultiBlockDataSet *lagrangian = vtkMultiBlockDataSet::New();
@@ -57,6 +58,7 @@ vtkPV3FoamReader::vtkPV3FoamReader()
 
     this->GetExecutive()->SetOutputData(1, lagrangian);
     lagrangian->Delete();
+#endif
 
     TimeStepRange[0] = 0;
     TimeStepRange[1] = 0;
@@ -319,15 +321,6 @@ int vtkPV3FoamReader::RequestData
         )
     );
 
-    vtkMultiBlockDataSet* lagrangianOutput = vtkMultiBlockDataSet::SafeDownCast
-    (
-        outputVector->GetInformationObject(1)->Get
-        (
-            vtkMultiBlockDataSet::DATA_OBJECT()
-        )
-    );
-
-
     if (Foam::vtkPV3Foam::debug)
     {
         cout<< "update output with "
@@ -376,7 +369,21 @@ int vtkPV3FoamReader::RequestData
 
 #else
 
-    foamData_->Update(output, lagrangianOutput);
+#ifdef VTKPV3FOAM_DUALPORT
+    foamData_->Update
+    (
+        output,
+        vtkMultiBlockDataSet::SafeDownCast
+        (
+            outputVector->GetInformationObject(1)->Get
+            (
+                vtkMultiBlockDataSet::DATA_OBJECT()
+            )
+        );
+    );
+#else
+    foamData_->Update(output, output);
+#endif
 
     if (ShowPatchNames)
     {
diff --git a/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C b/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C
index 624358b4a0..118e99fc91 100644
--- a/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C
+++ b/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C
@@ -465,7 +465,6 @@ void Foam::vtkPV3Foam::Update
         cout<< "<beg> Foam::vtkPV3Foam::Update - output with "
             << output->GetNumberOfBlocks() << " and "
             << lagrangianOutput->GetNumberOfBlocks() << " blocks\n";
-
         output->Print(cout);
         lagrangianOutput->Print(cout);
         printMemory();
@@ -504,8 +503,10 @@ void Foam::vtkPV3Foam::Update
         reader_->UpdateProgress(0.7);
     }
 
+#ifdef VTKPV3FOAM_DUALPORT
     // restart port1 at block=0
     blockNo = 0;
+#endif
     convertMeshLagrangian(lagrangianOutput, blockNo);
 
     reader_->UpdateProgress(0.8);
diff --git a/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.H b/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.H
index c6dd65af9e..4e203eacec 100644
--- a/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.H
+++ b/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.H
@@ -71,6 +71,8 @@ SourceFiles
 #include "PrimitivePatchInterpolation.H"
 #include "volPointInterpolation.H"
 
+#undef VTKPV3FOAM_DUALPORT
+
 // * * * * * * * * * * * * * Forward Declarations  * * * * * * * * * * * * * //
 
 class vtkDataArraySelection;
diff --git a/applications/utilities/postProcessing/miscellaneous/engineCompRatio/engineCompRatio.C b/applications/utilities/postProcessing/miscellaneous/engineCompRatio/engineCompRatio.C
index 34ebbc75dc..65356b382b 100644
--- a/applications/utilities/postProcessing/miscellaneous/engineCompRatio/engineCompRatio.C
+++ b/applications/utilities/postProcessing/miscellaneous/engineCompRatio/engineCompRatio.C
@@ -40,57 +40,50 @@ Description
 
 int main(int argc, char *argv[])
 {
+    #include "setRootCase.H"
+    #include "createEngineTime.H"
+    #include "createEngineMesh.H"
 
-#   include "setRootCase.H"
-#   include "createEngineTime.H"
-#   include "createEngineMesh.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
     scalar eps = 1.0e-10;
     scalar fullCycle = 360.0;
 
     scalar ca0 = -180.0;
     scalar ca1 = 0.0;
-    scalar ca = runTime.theta();
 
-    while (ca > ca0)
+    while (runTime.theta() > ca0)
     {
         ca0 += fullCycle;
         ca1 += fullCycle;
     }
 
-    if (mag(ca - ca0) > eps)
+    while (mag(runTime.theta() - ca0) > eps)
     {
-        while(mag(ca - ca0) > eps)
-        {
-            ca = runTime.theta();
-            scalar t0 = runTime.userTimeToTime(ca0 - ca);
-            runTime.setDeltaT(t0);
-            runTime++;
-            Info << "CA = " << ca << endl;
-            mesh.move();
-        }
+        scalar t0 = runTime.userTimeToTime(ca0 - runTime.theta());
+        runTime.setDeltaT(t0);
+        runTime++;
+        Info<< "CA = " << runTime.theta() << endl;
+        mesh.move();
     }
 
     scalar Vmax = sum(mesh.V().field());
 
-    while(mag(ca-ca1) > eps)
+    while (mag(runTime.theta()-ca1) > eps)
     {
-        ca = runTime.theta();
-        scalar t1 = runTime.userTimeToTime(ca1-ca);
+        scalar t1 = runTime.userTimeToTime(ca1-runTime.theta());
         runTime.setDeltaT(t1);
         runTime++;
-        Info << "CA = " << runTime.theta() << endl;
+        Info<< "CA = " << runTime.theta() << endl;
         mesh.move();
     }
 
     scalar Vmin = sum(mesh.V().field());
 
-    Info << "\nVmax = " << Vmax;
-    Info << ", Vmin = " << Vmin << endl;
-    Info << "Vmax/Vmin = " << Vmax/Vmin << endl;
-    Info << "\n end\n";
+    Info<< "\nVmax = " << Vmax;
+    Info<< ", Vmin = " << Vmin << endl;
+    Info<< "Vmax/Vmin = " << Vmax/Vmin << endl;
+    Info<< "\nEnd" << endl;
 
     return 0;
 }
diff --git a/src/OpenFOAM/db/dictionary/dictionary.C b/src/OpenFOAM/db/dictionary/dictionary.C
index 8c351510e7..fd89f0f608 100644
--- a/src/OpenFOAM/db/dictionary/dictionary.C
+++ b/src/OpenFOAM/db/dictionary/dictionary.C
@@ -529,7 +529,7 @@ bool Foam::dictionary::add(entry* entryPtr, bool mergeEntry)
             }
             else
             {
-                IOWarningIn("dictionary::add(entry*)", (*this))
+                IOWarningIn("dictionary::add(entry*, bool)", (*this))
                     << "problem replacing entry "<< entryPtr->keyword()
                     << " in dictionary " << name() << endl;
 
@@ -558,7 +558,7 @@ bool Foam::dictionary::add(entry* entryPtr, bool mergeEntry)
     }
     else
     {
-        IOWarningIn("dictionary::add(entry* entryPtr)", (*this))
+        IOWarningIn("dictionary::add(entry*, bool)", (*this))
             << "attempt to add entry "<< entryPtr->keyword()
             << " which already exists in dictionary " << name()
             << endl;
@@ -574,11 +574,13 @@ void Foam::dictionary::add(const entry& e, bool mergeEntry)
     add(e.clone(*this).ptr(), mergeEntry);
 }
 
+
 void Foam::dictionary::add(const keyType& k, const word& w, bool overwrite)
 {
     add(new primitiveEntry(k, token(w)), overwrite);
 }
 
+
 void Foam::dictionary::add
 (
     const keyType& k,
@@ -589,16 +591,19 @@ void Foam::dictionary::add
     add(new primitiveEntry(k, token(s)), overwrite);
 }
 
+
 void Foam::dictionary::add(const keyType& k, const label l, bool overwrite)
 {
     add(new primitiveEntry(k, token(l)), overwrite);
 }
 
+
 void Foam::dictionary::add(const keyType& k, const scalar s, bool overwrite)
 {
     add(new primitiveEntry(k, token(s)), overwrite);
 }
 
+
 void Foam::dictionary::add
 (
     const keyType& k,
@@ -628,6 +633,7 @@ void Foam::dictionary::set(const entry& e)
     set(e.clone(*this).ptr());
 }
 
+
 void Foam::dictionary::set(const keyType& k, const dictionary& d)
 {
     set(new dictionaryEntry(k, *this, d));
@@ -770,7 +776,7 @@ bool Foam::dictionary::merge(const dictionary& dict)
 
     bool changed = false;
 
-    forAllConstIter(IDLList<entry>, *this, iter)
+    forAllConstIter(IDLList<entry>, dict, iter)
     {
         HashTable<entry*>::iterator fnd = hashedEntries_.find(iter().keyword());
 
diff --git a/src/OpenFOAM/interpolations/interpolationTable/interpolationTable.C b/src/OpenFOAM/interpolations/interpolationTable/interpolationTable.C
index aa6d0d750f..9ea1307ee1 100644
--- a/src/OpenFOAM/interpolations/interpolationTable/interpolationTable.C
+++ b/src/OpenFOAM/interpolations/interpolationTable/interpolationTable.C
@@ -67,11 +67,25 @@ Foam::interpolationTable<Type>::interpolationTable()
 
 
 template<class Type>
-Foam::interpolationTable<Type>::interpolationTable(const fileName& fn)
+Foam::interpolationTable<Type>::interpolationTable
+(
+    const List<Tuple2<scalar, Type> >& values,
+    const boundsHandling bounds,
+    const fileName& fName
+)
+:
+    List<Tuple2<scalar, Type> >(values),
+    boundsHandling_(bounds),
+    fileName_(fName)
+{}
+
+
+template<class Type>
+Foam::interpolationTable<Type>::interpolationTable(const fileName& fName)
 :
     List<Tuple2<scalar, Type> >(),
     boundsHandling_(interpolationTable::WARN),
-    fileName_(fn)
+    fileName_(fName)
 {
     readTable();
 }
diff --git a/src/OpenFOAM/interpolations/interpolationTable/interpolationTable.H b/src/OpenFOAM/interpolations/interpolationTable/interpolationTable.H
index 2637ccf808..884c3cabd7 100644
--- a/src/OpenFOAM/interpolations/interpolationTable/interpolationTable.H
+++ b/src/OpenFOAM/interpolations/interpolationTable/interpolationTable.H
@@ -102,8 +102,16 @@ public:
         //- Construct null
         interpolationTable();
 
+        //- Construct from components
+        interpolationTable
+        (
+            const List<Tuple2<scalar, Type> >& values,
+            const boundsHandling bounds,
+            const fileName& fName
+        );
+
         //- Construct given the name of the file containing the table of data
-        interpolationTable(const fileName& fn);
+        interpolationTable(const fileName& fName);
 
         //- Construct by reading the fileName and boundsHandling from dictionary
         //  and read the table from that file.
diff --git a/src/OpenFOAM/meshes/pointMesh/pointPatches/constraint/processor/processorPointPatch.C b/src/OpenFOAM/meshes/pointMesh/pointPatches/constraint/processor/processorPointPatch.C
index 6ef8787d77..eee664b1fb 100644
--- a/src/OpenFOAM/meshes/pointMesh/pointPatches/constraint/processor/processorPointPatch.C
+++ b/src/OpenFOAM/meshes/pointMesh/pointPatches/constraint/processor/processorPointPatch.C
@@ -236,6 +236,7 @@ void processorPointPatch::initPatchPatchPoints()
     );
 
     toNeighbProc
+        << ppmp.size()              // number of points for checking
         << patchPatchPoints
         << patchPatchPointNormals;
 
@@ -257,13 +258,33 @@ void Foam::processorPointPatch::calcPatchPatchPoints()
         neighbProcNo()
     );
 
+    label nbrNPoints(readLabel(fromNeighbProc));
     labelListList patchPatchPoints(fromNeighbProc);
     List<List<vector> > patchPatchPointNormals(fromNeighbProc);
 
     pointBoundaryMesh& pbm = const_cast<pointBoundaryMesh&>(boundaryMesh());
-
     const labelList& ppmp = meshPoints();
 
+    // Simple check for the very rare situation when not the same number
+    // of points on both sides. This can happen with decomposed cyclics.
+    // If on one side the cyclic shares a point with proc faces coming from
+    // internal faces it will have a different number of points from
+    // the situation where the cyclic and the 'normal' proc faces are fully
+    // separate.
+    if (nbrNPoints != ppmp.size())
+    {
+        WarningIn("processorPointPatch::calcPatchPatchPoints()")
+            << "Processor patch " << name()
+            << " has " << ppmp.size() << " points; coupled patch has "
+            << nbrNPoints << " points." << endl
+            << "   (usually due to decomposed cyclics)."
+            << " This might give problems" << endl
+            << "    when using point fields (interpolation, mesh motion)."
+            << endl;
+    }
+
+
+
     // Loop over the patches looking for other patches that share points
     forAll(patchPatchPoints, patchi)
     {
diff --git a/src/dynamicFvMesh/dynamicFvMesh/dynamicFvMesh.C b/src/dynamicFvMesh/dynamicFvMesh/dynamicFvMesh.C
index aff80a2ac3..d0b81bc8c2 100644
--- a/src/dynamicFvMesh/dynamicFvMesh/dynamicFvMesh.C
+++ b/src/dynamicFvMesh/dynamicFvMesh/dynamicFvMesh.C
@@ -32,6 +32,7 @@ defineTypeNameAndDebug(Foam::dynamicFvMesh, 0);
 
 defineRunTimeSelectionTable(Foam::dynamicFvMesh, IOobject);
 
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::dynamicFvMesh::dynamicFvMesh(const IOobject& io)
@@ -45,4 +46,5 @@ Foam::dynamicFvMesh::dynamicFvMesh(const IOobject& io)
 Foam::dynamicFvMesh::~dynamicFvMesh()
 {}
 
+
 // ************************************************************************* //
diff --git a/src/dynamicFvMesh/dynamicFvMesh/dynamicFvMesh.H b/src/dynamicFvMesh/dynamicFvMesh/dynamicFvMesh.H
index dc437ea48e..04d7a69e59 100644
--- a/src/dynamicFvMesh/dynamicFvMesh/dynamicFvMesh.H
+++ b/src/dynamicFvMesh/dynamicFvMesh/dynamicFvMesh.H
@@ -26,7 +26,7 @@ Class
     Foam::dynamicFvMesh
 
 Description
-    The dynamicFvMesh is ..
+    Abstract base class for geometry and/or topology changing fvMesh.
 
 SourceFiles
     dynamicFvMesh.C
diff --git a/src/dynamicMesh/attachDetach/attachInterface.C b/src/dynamicMesh/attachDetach/attachInterface.C
index 21e4f882fd..4c4294f8f3 100644
--- a/src/dynamicMesh/attachDetach/attachInterface.C
+++ b/src/dynamicMesh/attachDetach/attachInterface.C
@@ -207,7 +207,20 @@ void Foam::attachDetach::attachInterface
                     mesh.faceZones()[modifiedFaceZone].whichFace(curFaceID)
                 ];
         }
-            
+
+
+        label patchID = mesh.boundaryMesh().whichPatch(curFaceID);
+        label neiCell;
+        if (patchID == -1)
+        {
+            neiCell = nei[curFaceID];
+        }
+        else
+        {
+            neiCell = -1;
+        }
+
+
         // Modify the face
         ref.setAction
         (
@@ -216,9 +229,9 @@ void Foam::attachDetach::attachInterface
                 newFace,                // modified face
                 curFaceID,              // label of face being modified
                 own[curFaceID],         // owner
-                nei[curFaceID],         // neighbour
+                neiCell,                // neighbour
                 false,                  // face flip
-                mesh.boundaryMesh().whichPatch(curFaceID),// patch for face
+                patchID,                // patch for face
                 false,                  // remove from zone
                 modifiedFaceZone,       // zone for face
                 modifiedFaceZoneFlip    // face flip in zone
diff --git a/src/dynamicMesh/attachDetach/detachInterface.C b/src/dynamicMesh/attachDetach/detachInterface.C
index 7e0ac0a933..6b57ccb685 100644
--- a/src/dynamicMesh/attachDetach/detachInterface.C
+++ b/src/dynamicMesh/attachDetach/detachInterface.C
@@ -73,6 +73,37 @@ void Foam::attachDetach::detachInterface
     const polyMesh& mesh = topoChanger().mesh();
     const faceZoneMesh& zoneMesh = mesh.faceZones();
 
+    // Check that zone is in increasing order (needed since adding faces
+    // in same order - otherwise polyTopoChange face ordering will mess up
+    // correspondence)
+    if (debug)
+    {
+        const labelList& faceLabels = zoneMesh[faceZoneID_.index()];
+        if (faceLabels.size() > 0)
+        {
+            for (label i = 1; i < faceLabels.size(); i++)
+            {
+                if (faceLabels[i] <= faceLabels[i-1])
+                {
+                    FatalErrorIn
+                    (
+                        "attachDetach::detachInterface"
+                        "(polyTopoChange&) const"
+                    )   << "faceZone " << zoneMesh[faceZoneID_.index()].name()
+                        << " does not have mesh face labels in"
+                        << " increasing order." << endl
+                        << "Face label " << faceLabels[i]
+                        << " at position " << i
+                        << " is smaller than the previous value "
+                        << faceLabels[i-1]
+                        << exit(FatalError);
+                }
+            }
+        }
+    }
+
+
+
     const primitiveFacePatch& masterFaceLayer = zoneMesh[faceZoneID_.index()]();
     const pointField& points = mesh.points();
     const labelListList& meshEdgeFaces = mesh.edgeFaces();
@@ -109,14 +140,11 @@ void Foam::attachDetach::detachInterface
 
         if (edgeIsInternal)
         {
-// Pout<< "Internal edge found: (" << mp[zoneLocalEdges[curEdgeID].start()] << " " << mp[zoneLocalEdges[curEdgeID].end()] << ")" << endl;
+            const edge& e = zoneLocalEdges[curEdgeID];
 
             // Reset the point creation
-            addedPoints[zoneLocalEdges[curEdgeID].start()] =
-                mp[zoneLocalEdges[curEdgeID].start()];
-
-            addedPoints[zoneLocalEdges[curEdgeID].end()] =
-                mp[zoneLocalEdges[curEdgeID].end()];
+            addedPoints[e.start()] = mp[e.start()];
+            addedPoints[e.end()] = mp[e.end()];
         }
     }
 // Pout << "addedPoints before point creation: " << addedPoints << endl;
@@ -137,7 +165,10 @@ void Foam::attachDetach::detachInterface
                         true                       // supports a cell
                     )
                 );
-// Pout << "Adding point " << points[mp[pointI]] << " for original point " << mp[pointI] << endl;
+            //Pout<< "Adding point " << addedPoints[pointI]
+            //    << " coord1:" << points[mp[pointI]]
+            //    << " coord2:" << masterFaceLayer.localPoints()[pointI]
+            //    << " for original point " << mp[pointI] << endl;
         }
     }
 
@@ -185,6 +216,7 @@ void Foam::attachDetach::detachInterface
             );
 
             // Add renumbered face into the slave patch
+            //label addedFaceI =
             ref.setAction
             (
                 polyAddFace
@@ -201,7 +233,15 @@ void Foam::attachDetach::detachInterface
                     false                           // zone flip
                 )
             );
-// Pout << "Flip.  Modifying face: " << faces[curFaceID].reverseFace() << " next to cell: " << nei[curFaceID] << " and adding face: " << newFace << " next to cell: " << own[curFaceID] << endl;
+            //{
+            //    pointField newPts(ref.points());
+            //Pout<< "Flip.  Modifying face: " << ref.faces()[curFaceID]
+            //    << " fc:" <<  ref.faces()[curFaceID].centre(newPts)
+            //    << " next to cell: " << nei[curFaceID]
+            //    << " and adding face: " << newFace
+            //    << " fc:" << ref.faces()[addedFaceI].centre(newPts)
+            //    << " next to cell: " << own[curFaceID] << endl;
+            //}
         }
         else
         {
@@ -223,6 +263,7 @@ void Foam::attachDetach::detachInterface
             );
 
             // Add renumbered face into the slave patch
+            //label addedFaceI =
             ref.setAction
             (
                 polyAddFace
@@ -239,7 +280,15 @@ void Foam::attachDetach::detachInterface
                     false                           // face flip in zone
                 )
             );
-// Pout << "No flip.  Modifying face: " << faces[curFaceID] << " next to cell: " << own[curFaceID] << " and adding face: " << newFace << " next to cell: " << nei[curFaceID] << endl;
+            //{
+            //    pointField newPts(ref.points());
+            //Pout<< "No flip.  Modifying face: " << ref.faces()[curFaceID]
+            //    << " fc:" <<  ref.faces()[curFaceID].centre(newPts)
+            //    << " next to cell: " << own[curFaceID]
+            //    << " and adding face: " << newFace
+            //    << " fc:" << ref.faces()[addedFaceI].centre(newPts)
+            //    << " next to cell: " << nei[curFaceID] << endl;
+            //}
         }
     }
 
diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C
index f317573d8a..d2728f8b9d 100644
--- a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C
+++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C
@@ -1464,6 +1464,11 @@ void Foam::polyTopoChange::resetZones
 
             addressing[zoneI][nPoints[zoneI]++] = iter.key();
         }
+        // Sort the addressing
+        forAll(addressing, zoneI)
+        {
+            stableSort(addressing[zoneI]);
+        }
 
         // So now we both have old zones and the new addressing.
         // Invert the addressing to get pointZoneMap.
@@ -1551,6 +1556,28 @@ void Foam::polyTopoChange::resetZones
             addressing[zoneI][index] = faceI;
             flipMode[zoneI][index] = faceZoneFlip_[faceI];
         }
+        // Sort the addressing
+        forAll(addressing, zoneI)
+        {
+            labelList newToOld;
+            sortedOrder(addressing[zoneI], newToOld);
+            {
+                labelList newAddressing(addressing[zoneI].size());
+                forAll(newAddressing, i)
+                {
+                    newAddressing[i] = addressing[zoneI][newToOld[i]];
+                }
+                addressing[zoneI].transfer(newAddressing);
+            }
+            {
+                boolList newFlipMode(flipMode[zoneI].size());
+                forAll(newFlipMode, i)
+                {
+                    newFlipMode[i] = flipMode[zoneI][newToOld[i]];
+                }
+                flipMode[zoneI].transfer(newFlipMode);
+            }
+        }
 
         // So now we both have old zones and the new addressing.
         // Invert the addressing to get faceZoneFaceMap.
@@ -1644,6 +1671,11 @@ void Foam::polyTopoChange::resetZones
                 addressing[zoneI][nCells[zoneI]++] = cellI;
             }
         }
+        // Sort the addressing
+        forAll(addressing, zoneI)
+        {
+            stableSort(addressing[zoneI]);
+        }
 
         // So now we both have old zones and the new addressing.
         // Invert the addressing to get cellZoneMap.
diff --git a/src/fvMotionSolver/fvMotionSolvers/displacement/displacementFvMotionSolver/displacementFvMotionSolver.C b/src/fvMotionSolver/fvMotionSolvers/displacement/displacementFvMotionSolver/displacementFvMotionSolver.C
index aee52c717c..aa3d51bd4f 100644
--- a/src/fvMotionSolver/fvMotionSolvers/displacement/displacementFvMotionSolver/displacementFvMotionSolver.C
+++ b/src/fvMotionSolver/fvMotionSolvers/displacement/displacementFvMotionSolver/displacementFvMotionSolver.C
@@ -32,23 +32,13 @@ License
 
 namespace Foam
 {
-//    defineTypeNameAndDebug(displacementFvMotionSolver, 0);
-//
-//    addToRunTimeSelectionTable
-//    (
-//        fvMotionSolver,
-//        displacementFvMotionSolver,
-//        dictionary
-//    );
+    defineTypeNameAndDebug(displacementFvMotionSolver, 0);
 }
 
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::displacementFvMotionSolver::
-displacementFvMotionSolver
+Foam::displacementFvMotionSolver::displacementFvMotionSolver
 (
     const polyMesh& mesh,
     Istream&
diff --git a/src/thermophysicalModels/radiation/derivedFvPatchFields/greyDiffusiveRadiation/greyDiffusiveRadiationMixedFvPatchScalarField.C b/src/thermophysicalModels/radiation/derivedFvPatchFields/greyDiffusiveRadiation/greyDiffusiveRadiationMixedFvPatchScalarField.C
index e7daa58aaf..1f9a6ac63c 100644
--- a/src/thermophysicalModels/radiation/derivedFvPatchFields/greyDiffusiveRadiation/greyDiffusiveRadiationMixedFvPatchScalarField.C
+++ b/src/thermophysicalModels/radiation/derivedFvPatchFields/greyDiffusiveRadiation/greyDiffusiveRadiationMixedFvPatchScalarField.C
@@ -151,7 +151,10 @@ updateCoeffs()
     const scalarField& Tp =
         patch().lookupPatchField<volScalarField, scalar>(TName_);
 
-    const fvDOM& dom = db().lookupObject<fvDOM>("radiationProperties");
+    const radiationModel& radiation =
+            db().lookupObject<radiationModel>("radiationProperties");
+
+    const fvDOM& dom(refCast<const fvDOM>(radiation));
 
     const label patchI = patch().index();
 
diff --git a/src/thermophysicalModels/radiation/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.C b/src/thermophysicalModels/radiation/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.C
index eb65337aca..c2d89834cd 100644
--- a/src/thermophysicalModels/radiation/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.C
+++ b/src/thermophysicalModels/radiation/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.C
@@ -147,7 +147,10 @@ updateCoeffs()
         return;
     }
 
-    const fvDOM& dom = db().lookupObject<fvDOM>("radiationProperties");
+    const radiationModel& radiation =
+        db().lookupObject<radiationModel>("radiationProperties");
+
+    const fvDOM& dom(refCast<const fvDOM>(radiation));
 
     const label patchI = patch().index();
 
diff --git a/src/thermophysicalModels/radiation/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.H b/src/thermophysicalModels/radiation/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.H
index 3ff3990786..9690b0cf5d 100644
--- a/src/thermophysicalModels/radiation/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.H
+++ b/src/thermophysicalModels/radiation/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.H
@@ -66,9 +66,6 @@ class wideBandDiffusiveRadiationMixedFvPatchScalarField
         //- Wavelength index
         label lambdaId_;
 
-        //- Radiative heat flux on walls.
-        scalarField qr_;
-
 
 public:
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmission.C b/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmission.C
index 39d02831bb..c78041bd71 100644
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmission.C
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmission.C
@@ -27,19 +27,130 @@ License
 #include "blackBodyEmission.H"
 #include "dimensionedConstants.H"
 
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+const Foam::List<Foam::Tuple2<Foam::scalar, Foam::scalar> >
+Foam::radiation::blackBodyEmission::emissivePowerTable
+(
+    IStringStream
+    (
+        "( 1000 0.00032)"
+        "( 1100 0.00091)"
+        "( 1200 0.00213)"
+        "( 1300 0.00432)"
+        "( 1400 0.00779)"
+        "( 1500 0.01280)"
+        "( 1600 0.01972)"
+        "( 1700 0.02853)"
+        "( 1800 0.03934)"
+        "( 1900 0.05210)"
+        "( 2000 0.06672)"
+        "( 2100 0.08305)"
+        "( 2200 0.10088)"
+        "( 2300 0.12002)"
+        "( 2400 0.14025)"
+        "( 2500 0.16135)"
+        "( 2600 0.18311)"
+        "( 2700 0.20535)"
+        "( 2800 0.22788)"
+        "( 2900 0.25055)"
+        "( 3000 0.27322)"
+        "( 3100 0.29576)"
+        "( 3200 0.31809)"
+        "( 3300 0.34009)"
+        "( 3400 0.36172)"
+        "( 3500 0.38290)"
+        "( 3600 0.40359)"
+        "( 3700 0.42375)"
+        "( 3800 0.44336)"
+        "( 3900 0.46240)"
+        "( 4000 0.48085)"
+        "( 4100 0.49872)"
+        "( 4200 0.51599)"
+        "( 4300 0.53267)"
+        "( 4400 0.54877)"
+        "( 4500 0.56429)"
+        "( 4600 0.57925)"
+        "( 4700 0.59366)"
+        "( 4800 0.60753)"
+        "( 4900 0.62088)"
+        "( 5000 0.63372)"
+        "( 5100 0.64606)"
+        "( 5200 0.65794)"
+        "( 5300 0.66935)"
+        "( 5400 0.68033)"
+        "( 5500 0.69087)"
+        "( 5600 0.70101)"
+        "( 5700 0.71076)"
+        "( 5800 0.72012)"
+        "( 5900 0.72913)"
+        "( 6000 0.73778)"
+        "( 6100 0.74610)"
+        "( 6200 0.75410)"
+        "( 6300 0.76180)"
+        "( 6400 0.76920)"
+        "( 6500 0.77631)"
+        "( 6600 0.78316)"
+        "( 6700 0.78975)"
+        "( 6800 0.79609)"
+        "( 6900 0.80219)"
+        "( 7000 0.80807)"
+        "( 7100 0.81373)"
+        "( 7200 0.81918)"
+        "( 7300 0.82443)"
+        "( 7400 0.82949)"
+        "( 7500 0.83436)"
+        "( 7600 0.83906)"
+        "( 7700 0.84359)"
+        "( 7800 0.84796)"
+        "( 7900 0.85218)"
+        "( 8000 0.85625)"
+        "( 8100 0.86017)"
+        "( 8200 0.86396)"
+        "( 8300 0.86762)"
+        "( 8400 0.87115)"
+        "( 8500 0.87456)"
+        "( 8600 0.87786)"
+        "( 8700 0.88105)"
+        "( 8800 0.88413)"
+        "( 8900 0.88711)"
+        "( 9000 0.88999)"
+        "( 9100 0.89277)"
+        "( 9200 0.89547)"
+        "( 9300 0.89807)"
+        "( 9400 0.90060)"
+        "( 9500 0.90304)"
+        "( 9600 0.90541)"
+        "( 9700 0.90770)"
+        "( 9800 0.90992)"
+        "( 9900 0.91207)"
+        "(10000 0.91415)"
+        "(12000 0.94505)"
+        "(15000 0.96893)"
+        "(20000 0.98555)"
+        "(30000 0.99529)"
+        "(40000 0.99792)"
+        "(50000 0.99890)"
+    )()
+);
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::radiation::blackBodyEmission::blackBodyEmission
 (
-    const fileName& name,
-    const word& instance,
     const label nLambda,
     const volScalarField& T
 )
 :
-    blackBodyEmissiveTable_(name, instance, T.mesh()),
-    C1_("C1",dimensionSet(1, 4, 3, 0, 0, 0, 0), 3.7419e-16),
-    C2_("C2",dimensionSet(0, 1, 0, 1, 0, 0, 0), 14.388e-6),
+    table_
+    (
+        emissivePowerTable,
+        interpolationTable<scalar>::CLAMP,
+        "blackBodyEmissivePower"
+    ),
+    C1_("C1", dimensionSet(1, 4, 3, 0, 0, 0, 0), 3.7419e-16),
+    C2_("C2", dimensionSet(0, 1, 0, 1, 0, 0, 0), 14.388e-6),
     bLambda_(nLambda),
     T_(T)
 {
@@ -79,7 +190,7 @@ Foam::scalar Foam::radiation::blackBodyEmission::fLambdaT
     const scalar lambdaT
 ) const
 {
-    return  blackBodyEmissiveTable_.lookUp(lambdaT*1.0e6)[1];
+    return  table_(lambdaT*1.0e6);
 }
 
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmission.H b/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmission.H
index d21582938f..01719e5076 100644
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmission.H
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmission.H
@@ -28,6 +28,9 @@ Class
 Description
     Class black body emission
 
+    Table of black body emissive power taken from:
+        Modest, "Radiative Heat Transfer", pp.775-777, 1993
+
 SourceFiles
     blackBodyEmission.C
 
@@ -40,7 +43,7 @@ SourceFiles
 #include "dimensionedScalar.H"
 #include "mathematicalConstants.H"
 #include "radiationConstants.H"
-#include "interpolationLookUpTable.H"
+#include "interpolationTable.H"
 #include "Vector2D.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -56,9 +59,18 @@ namespace radiation
 
 class blackBodyEmission
 {
+public:
+
+    //- Static table of black body emissive power
+    static const List<Tuple2<scalar, scalar> > emissivePowerTable;
+
+
+private:
+
     // Private data
 
-        mutable interpolationLookUpTable<scalar> blackBodyEmissiveTable_;
+        //- Interpolation table of black body emissive power
+        mutable interpolationTable<scalar> table_;
 
         //- Constant C1
         const dimensionedScalar C1_;
@@ -85,8 +97,6 @@ public:
         //- Construct from components
         blackBodyEmission
         (
-            const fileName& name,
-            const word& instance,
             const label nLambda,
             const volScalarField& T
         );
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmissivePower b/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmissivePower
deleted file mode 100644
index c0151b436a..0000000000
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/blackBodyEmission/blackBodyEmissivePower
+++ /dev/null
@@ -1,174 +0,0 @@
-171
-(
-1000    0.00032
-1100    0.00091
-1200    0.00213
-1300    0.00432
-1400    0.00779
-1500    0.01285
-1600    0.01972
-1700    0.02853
-1800    0.03934
-1900    0.05210
-2000    0.06672
-2100    0.08305
-2200    0.10088
-2300    0.12002
-2400    0.14025
-2500    0.16135
-2600    0.18311
-2700    0.20535
-2800    0.22788
-2900    0.25055
-3000    0.27322
-3100    0.29576
-3200    0.31809
-3300    0.34009
-3400    0.36172
-3500    0.38290
-3600    0.40359
-3700    0.42375
-3800    0.44336
-3900    0.46240
-4000    0.48085
-4100    0.49872
-4200    0.51599
-4300    0.53267
-4400    0.54877
-4500    0.56429
-4600    0.57925
-4700    0.59366
-4800    0.60753
-4900    0.62088
-5000    0.63372
-5100    0.64606
-5200    0.65794
-5300    0.66935
-5400    0.68033
-5500    0.69087
-5600    0.70101
-5700    0.71076
-5800    0.72012
-5900    0.72913
-6000    0.73778
-6100    0.74610
-6200    0.75410
-6300    0.76180
-6400    0.76920
-6500    0.77631
-6600    0.78316
-6700    0.78975
-6800    0.79609
-6900    0.80219
-7000    0.80807
-7100    0.81373
-7200    0.81918
-7300    0.82443
-7400    0.82949
-7500    0.83436
-7600    0.83906
-7700    0.84359
-7800    0.84796
-7900    0.85218
-8000    0.85625
-8100    0.86017
-8200    0.86396
-8300    0.86762
-8400    0.87115
-8500    0.87456
-8600    0.87786
-8700    0.88105
-8800    0.88413
-8900    0.88711
-9000    0.88999
-9100   0.89277
-9200   0.89547
-9300   0.89807
-9400   0.90060
-9500   0.90304
-9600   0.90541
-9700   0.90770
-9800   0.90992
-9900   0.91207
-10000   0.91415
-10200   0.91813
-10400   0.92188
-10600   0.92540
-10800   0.92872
-11000   0.93184
-11200   0.93479
-11400   0.93758
-11600   0.94021
-11800   0.94270
-12000   0.94505
-12200   0.94728
-12400   0.94939
-12600   0.95139
-12800   0.95329
-13000   0.95509
-13200   0.95680
-13400   0.95843
-13600   0.95998
-13800   0.96145
-14000   0.96285
-14200   0.96418
-14400   0.96546
-14600   0.96667
-14800   0.96783
-15000   0.96893
-15200   0.96999
-15400   0.97100
-15600   0.97196
-15800   0.97288
-16000   0.97377
-16200   0.97461
-16400   0.97542
-16600   0.97620
-16800   0.97694
-17000   0.97765
-17200   0.97834
-17400   0.97899
-17600   0.97962
-17800   0.98023
-18000   0.98080
-18200   0.98137
-18400   0.98191
-18600   0.98243
-18900   0.98293
-19000   0.98340
-19200   0.98387
-19400   0.98431
-19600   0.98474
-19800   0.98515
-20000   0.98555
-21000   0.98735
-22000   0.98886
-23000   0.99014
-24000   0.99123
-25000   0.99217
-26000   0.99297
-27000   0.99367
-28000   0.99429
-29000   0.99482
-30000   0.99529
-31000   0.99571
-32000   0.99607
-33000   0.99640
-34000   0.99669
-35000   0.99695
-35000   0.99719
-36000   0.99740
-37000   0.99759
-38000   0.99776
-39000   0.99792
-40000   0.99806
-41000   0.99819
-42000   0.99831
-43000   0.99842
-44000   0.99851
-45000   0.99861
-46000   0.99869
-47000   0.99877
-48000   0.99884
-49000   0.99890
-)
\ No newline at end of file
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.C b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.C
index af3b14ebb0..c8fb606b4c 100644
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.C
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.C
@@ -126,20 +126,14 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
     nRay_(0),
     nLambda_(absorptionEmission_->nBands()),
     aLambda_(nLambda_),
-    blackBody_
-    (
-        fileName("blackBodyEmissivePower"),
-        mesh_.time().constant(),
-        nLambda_,
-        T
-    ),
+    blackBody_(nLambda_, T),
     IRay_(0),
     convergence_(coeffs_.lookupOrDefault<scalar>("convergence", 0.0))
 {
     if (mesh_.nSolutionD() == 3)    //3D
     {
-        IRay_.setSize(4*nPhi_*nTheta_);
         nRay_ = 4*nPhi_*nTheta_;
+        IRay_.setSize(nRay_);
         scalar deltaPhi = mathematicalConstant::pi/(2.0*nPhi_);
         scalar deltaTheta = mathematicalConstant::pi/nTheta_;
         label i = 0;
@@ -175,8 +169,8 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
         {
             scalar thetai = mathematicalConstant::pi/2.0;
             scalar deltaTheta = mathematicalConstant::pi;
-            IRay_.setSize(4*nPhi_);
             nRay_ = 4*nPhi_;
+            IRay_.setSize(nRay_);
             scalar deltaPhi = mathematicalConstant::pi/(2.0*nPhi_);
             label i = 0;
             for (label m = 1; m <= 4*nPhi_; m++)
@@ -205,8 +199,8 @@ Foam::radiation::fvDOM::fvDOM(const volScalarField& T)
         {
             scalar thetai = mathematicalConstant::pi/2.0;
             scalar deltaTheta = mathematicalConstant::pi;
-            IRay_.setSize(2);
             nRay_ = 2;
+            IRay_.setSize(nRay_);
             scalar deltaPhi = mathematicalConstant::pi;
             label i = 0;
             for (label m = 1; m <= 2; m++)
@@ -302,7 +296,7 @@ void Foam::radiation::fvDOM::calculate()
             maxResidual = max(maxBandResidual, maxResidual);
         }
 
-        Info << "Radiation solver iter: " <<  radIter << endl;
+        Info << "Radiation solver iter: " << radIter << endl;
 
     } while(maxResidual > convergence_);
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.H b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.H
index 3cbd921406..1d2662e307 100644
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.H
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/fvDOM/fvDOM.H
@@ -37,19 +37,19 @@ Description
     i.e. dictionary
     fvDOMCoeffs
     {
-        Nphi    1;          // azimuthal angles in PI/2 on X-Y.(from Y to X)
-        Ntheta  2;          // polar angles in P1 (from Z to X-Y plane)
+        nPhi    1;          // azimuthal angles in PI/2 on X-Y.(from Y to X)
+        nTheta  2;          // polar angles in PI (from Z to X-Y plane)
         convergence 1e-4;   // convergence criteria for radiation iteration
     }
 
-    nFlowIterPerRadIter   1; // Number of flow iterations per radiation
+    solverFreq   1; // Number of flow iterations per radiation
                                 iteration
 
-    The total number of solid angles is  4*Nphi*Ntheta.
+    The total number of solid angles is  4*nPhi*nTheta.
 
-    In 1D the direction of the rays is X (Nphi and Ntheta are ignored)
-    In 2D the direction of the rays is on X-Y plane (only Nphi is considered)
-    In 3D (Nphi and Ntheta are considered)
+    In 1D the direction of the rays is X (nPhi and nTheta are ignored)
+    In 2D the direction of the rays is on X-Y plane (only nPhi is considered)
+    In 3D (nPhi and nTheta are considered)
 
 SourceFiles
     fvDOM.C
@@ -127,9 +127,6 @@ class fvDOM
         //- Disallow default bitwise assignment
         void operator=(const fvDOM&);
 
-        //- Update Absorption Coefficients
-//        void updateAbsorptionCoeffs(void);
-
         //- Update nlack body emission
         void updateBlackBodyEmission();
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/interpolationLookUpTable/interpolationLookUpTable.C b/src/thermophysicalModels/radiation/radiationModel/fvDOM/interpolationLookUpTable/interpolationLookUpTable.C
index e9e28658c4..f000e43fa3 100755
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/interpolationLookUpTable/interpolationLookUpTable.C
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/interpolationLookUpTable/interpolationLookUpTable.C
@@ -29,23 +29,23 @@ License
 // * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //
 
 template <class Type>
-Foam::label Foam::interpolationLookUpTable <Type>::index
+Foam::label Foam::interpolationLookUpTable<Type>::index
 (
     const List<scalar>& indices,
     const bool lastDim
 ) const
 {
-    label totalindex = 0;
+    label totalIndex = 0;
 
-    for (int i = 0; i < dim_.size() - 1; i++)
+    forAll(dim_, i)
     {
         label dim = 1;
         for (int j = i + 1; j < dim_.size(); j++)
         {
-            dim *=(dim_[j]+1);
+            dim *= dim_[j] + 1;
         }
 
-        totalindex +=
+        totalIndex +=
             dim
            *min
             (
@@ -57,7 +57,7 @@ Foam::label Foam::interpolationLookUpTable <Type>::index
     if (lastDim)
     {
         label iLastdim = dim_.size() - 1;
-        totalindex += Foam::min
+        totalIndex += Foam::min
         (
             max
             (
@@ -68,12 +68,12 @@ Foam::label Foam::interpolationLookUpTable <Type>::index
         );
     }
 
-    return totalindex;
+    return totalIndex;
 }
 
 
 template <class Type>
-Foam::label Foam::interpolationLookUpTable <Type>::index
+Foam::label Foam::interpolationLookUpTable<Type>::index
 (
     const scalar indice
 ) const
@@ -171,15 +171,15 @@ void Foam::interpolationLookUpTable<Type>::dimensionTable()
         max_[i] = readScalar(entries_[i].lookup("max"));
         min_[i] = readScalar(entries_[i].lookup("min"));
         delta_[i] = (max_[i] - min_[i])/dim_[i];
-        tableDim *= (dim_[i] + 1);
-        fieldIndices_.insert(entries_[i].lookup("name"),index);
+        tableDim *= dim_[i] + 1;
+        fieldIndices_.insert(entries_[i].lookup("name"), index);
         entryIndices_[i] = index;
         index++;
     }
 
     forAll(output_,i)
     {
-        fieldIndices_.insert(output_[i].lookup("name"),index);
+        fieldIndices_.insert(output_[i].lookup("name"), index);
         outputIndices_[i] = index;
         index++;
     }
@@ -229,8 +229,7 @@ void Foam::interpolationLookUpTable<Type>::readTable
         FatalErrorIn
         (
             "Foam::interpolationLookUpTable<Type>::readTable()"
-        )   << "table is empty" << nl
-            << exit(FatalError);
+        )   << "table is empty" << nl << exit(FatalError);
     }
 }
 
@@ -240,7 +239,7 @@ void Foam::interpolationLookUpTable<Type>::readTable
 template<class Type>
 Foam::interpolationLookUpTable<Type>::interpolationLookUpTable()
 :
-    List<scalarField >(),
+    List<scalarField>(),
     fileName_("fileNameIsUndefined")
 {}
 
@@ -251,7 +250,7 @@ Foam::interpolationLookUpTable<Type>::interpolationLookUpTable
     const fileName& fn, const word& instance, const fvMesh& mesh
 )
 :
-    List<scalarField >(),
+    List<scalarField>(),
     fileName_(fn),
     dim_(0),
     min_(0),
@@ -273,7 +272,7 @@ Foam::interpolationLookUpTable<Type>::interpolationLookUpTable
      const interpolationLookUpTable& interpTable
 )
 :
-    List<scalarField >(interpTable),
+    List<scalarField>(interpTable),
     fileName_(interpTable.fileName_),
     entryIndices_(interpTable.entryIndices_),
     outputIndices_(interpTable.outputIndices_),
@@ -293,7 +292,7 @@ Foam::interpolationLookUpTable<Type>::interpolationLookUpTable
     const dictionary& dict
 )
 :
-    List<scalarField >(),
+    List<scalarField>(),
     fileName_(fileName(dict.lookup("fileName")).expand()),
     dim_(0),
     min_(0.0),
@@ -315,20 +314,19 @@ Foam::interpolationLookUpTable<Type>::interpolationLookUpTable
 template<class Type>
 void Foam::interpolationLookUpTable<Type>::check() const
 {
-//  check order in the first dimension.
-
-    scalar prevValue = List<scalarField >::operator[](0).operator[](0);
-    label dim = 1 ;
+    // check order in the first dimension.
+    scalar prevValue = List<scalarField>::operator[](0).operator[](0);
+    label dim = 1;
     for (int j = 1; j < dim_.size(); j++)
     {
-        dim *=(dim_[j]+1);
+        dim *= dim_[j] + 1;
     }
 
     for (label i = 1; i < dim_[0]; i++)
     {
         label index = i*dim;
         const scalar currValue =
-            List<scalarField >::operator[](0).operator[](index);
+            List<scalarField>::operator[](0).operator[](index);
 
         // avoid duplicate values (divide-by-zero error)
         if (currValue <= prevValue)
@@ -336,9 +334,8 @@ void Foam::interpolationLookUpTable<Type>::check() const
             FatalErrorIn
             (
                 "Foam::interpolationLookUpTable<Type>::checkOrder() const"
-            )   << "out-of-order value: "
-                << currValue << " at index " << index << nl
-                << exit(FatalError);
+            )   << "out-of-order value: " << currValue
+                << " at index " << index << nl << exit(FatalError);
         }
         prevValue = currValue;
     }
@@ -379,8 +376,7 @@ void Foam::interpolationLookUpTable<Type>::write
         FatalErrorIn
         (
             "Foam::interpolationTable<Type>::write()"
-        )   << "table is empty" << nl
-            << exit(FatalError);
+        )   << "table is empty" << nl << exit(FatalError);
     }
     os.writeKeyword("values");
     os << *this << token::END_STATEMENT << nl;
@@ -400,31 +396,25 @@ Foam::interpolationLookUpTable<Type>::operator[](const label i)
     {
         FatalErrorIn
         (
-            "Foam::interpolationLookUpTable<Type>::operator[]"
-            "(const label) const"
-        )   << "table has (" << n << ") columns" << nl
-            << exit(FatalError);
+            "Foam::interpolationLookUpTable<Type>::operator[](const label)"
+        )   << "table has (" << n << ") columns" << nl << exit(FatalError);
     }
     else if (ii < 0)
     {
         FatalErrorIn
         (
-            "Foam::interpolationLookUpTable<Type>::operator[]"
-            "(const label) const"
-        )   << "index (" << ii << ") underflow" << nl
-            << exit(FatalError);
+            "Foam::interpolationLookUpTable<Type>::operator[](const label)"
+        )   << "index (" << ii << ") underflow" << nl << exit(FatalError);
     }
     else if (ii > n)
     {
         FatalErrorIn
         (
-            "Foam::interpolationLookUpTable<Type>::operator[]"
-            "(const label) const"
-        )   << "index (" << ii << ") overflow" << nl
-            << exit(FatalError);
+            "Foam::interpolationLookUpTable<Type>::operator[](const label)"
+        )   << "index (" << ii << ") overflow" << nl << exit(FatalError);
     }
 
-    return List<scalarField >::operator[](ii);
+    return List<scalarField>::operator[](ii);
 }
 
 
@@ -441,8 +431,7 @@ Foam::interpolationLookUpTable<Type>::operator[](const label i) const
         (
             "Foam::interpolationLookUpTable<Type>::operator[]"
             "(const label) const"
-        )   << "table has (" << n << ") columns" << nl
-            << exit(FatalError);
+        )   << "table has (" << n << ") columns" << nl << exit(FatalError);
     }
     else if (ii < 0)
     {
@@ -450,8 +439,7 @@ Foam::interpolationLookUpTable<Type>::operator[](const label i) const
         (
             "Foam::interpolationLookUpTable<Type>::operator[]"
             "(const label) const"
-        )   << "index (" << ii << ") underflow" << nl
-            << exit(FatalError);
+        )   << "index (" << ii << ") underflow" << nl << exit(FatalError);
     }
 
     else if (ii > n)
@@ -464,7 +452,7 @@ Foam::interpolationLookUpTable<Type>::operator[](const label i) const
             << exit(FatalError);
     }
 
-    return List<scalarField >::operator[](ii);
+    return List<scalarField>::operator[](ii);
 }
 
 
@@ -500,7 +488,7 @@ void Foam::interpolationLookUpTable<Type>::findHi
 
         forAll(entryIndices_,i)
         {
-            if (checkRange(retvals,entryIndices_[i]))
+            if (checkRange(retvals, entryIndices_[i]))
             {
                 label dim = 1;
 
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/interpolationLookUpTable/interpolationLookUpTable.H b/src/thermophysicalModels/radiation/radiationModel/fvDOM/interpolationLookUpTable/interpolationLookUpTable.H
index d2518d4e02..6140ae4bc3 100755
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/interpolationLookUpTable/interpolationLookUpTable.H
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/interpolationLookUpTable/interpolationLookUpTable.H
@@ -54,7 +54,7 @@ namespace Foam
 {
 
 /*---------------------------------------------------------------------------*\
-                    Class interpolationLookUpTable Declaration
+                  Class interpolationLookUpTable Declaration
 \*---------------------------------------------------------------------------*/
 
 template<class Type>
@@ -64,22 +64,22 @@ class interpolationLookUpTable
 {
 private:
 
-    // Private data
+    // Privsate data
 
         //- File name
         fileName fileName_;
 
         //- Table dimensions
-        List<label>  dim_;
+        List<label> dim_;
 
         //- Min on each dimension
-        List<scalar>  min_;
+        List<scalar> min_;
 
         //- Deltas on each dimension
-        List<scalar>  delta_;
+        List<scalar> delta_;
 
         //- Maximum on each dimension
-        List<scalar>  max_;
+        List<scalar> max_;
 
         //- Dictionary entries
         List<dictionary> entries_;
diff --git a/src/thermophysicalModels/radiation/radiationModel/fvDOM/radiativeIntensityRay/radiativeIntensityRay.C b/src/thermophysicalModels/radiation/radiationModel/fvDOM/radiativeIntensityRay/radiativeIntensityRay.C
index 96090a3ef2..cba8dd855a 100644
--- a/src/thermophysicalModels/radiation/radiationModel/fvDOM/radiativeIntensityRay/radiativeIntensityRay.C
+++ b/src/thermophysicalModels/radiation/radiationModel/fvDOM/radiativeIntensityRay/radiativeIntensityRay.C
@@ -113,11 +113,11 @@ Foam::radiation::radiativeIntensityRay::radiativeIntensityRay
         0.5*deltaPhi*Foam::sin(2.0*theta)*Foam::sin(deltaTheta)
     );
 
-    forAll(ILambda_, i)
+    forAll(ILambda_, lambdaI)
     {
         IOobject IHeader
         (
-            "ILambda_" + name(rayId) + "_" + name(i),
+            "ILambda_" + name(rayId) + "_" + name(lambdaI),
             mesh_.time().timeName(),
             mesh_,
             IOobject::MUST_READ,
@@ -129,7 +129,7 @@ Foam::radiation::radiativeIntensityRay::radiativeIntensityRay
         {
             ILambda_.set
             (
-                i,
+                lambdaI,
                 new volScalarField(IHeader, mesh_)
             );
         }
@@ -150,7 +150,7 @@ Foam::radiation::radiativeIntensityRay::radiativeIntensityRay
 
             ILambda_.set
             (
-                i,
+                lambdaI,
                 new volScalarField(IHeader, IDefault)
             );
         }
@@ -182,7 +182,7 @@ Foam::scalar Foam::radiation::radiativeIntensityRay::correct()
 
         fvScalarMatrix IiEq
         (
-            fvm::div(Ji, ILambda_[lambdaI], " div(Ji,Ii_h)")
+            fvm::div(Ji, ILambda_[lambdaI], "div(Ji,Ii_h)")
           + fvm::Sp(k*omega_, ILambda_[lambdaI])
          ==
             1.0/Foam::mathematicalConstant::pi
@@ -201,7 +201,6 @@ Foam::scalar Foam::radiation::radiativeIntensityRay::correct()
         ).initialResidual();
 
         maxResidual = max(eqnResidual, maxResidual);
-
     }
 
     return maxResidual;
diff --git a/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.C b/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.C
index 131221c308..92c6397fb5 100644
--- a/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.C
+++ b/src/thermophysicalModels/radiation/radiationModel/radiationModel/radiationModel.C
@@ -51,7 +51,7 @@ Foam::radiation::radiationModel::radiationModel(const volScalarField& T)
         (
             "radiationProperties",
             T.time().constant(),
-            T.mesh().objectRegistry::db(),
+            T.mesh(),
             IOobject::MUST_READ,
             IOobject::NO_WRITE
         )
@@ -79,7 +79,7 @@ Foam::radiation::radiationModel::radiationModel
         (
             "radiationProperties",
             T.time().constant(),
-            T.mesh().objectRegistry::db(),
+            T.mesh(),
             IOobject::MUST_READ,
             IOobject::NO_WRITE
         )
diff --git a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.C b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.C
index 96eb031fdd..dccbbbad25 100644
--- a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.C
+++ b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.C
@@ -203,7 +203,7 @@ Foam::radiation::greyMeanAbsorptionEmission::aCont(const label bandI) const
             {
                 Ti = 1./T[i];
             }
-            a[i]+=
+            a[i] +=
                 Yipi
                *(
                     ((((b[5]*Ti + b[4])*Ti + b[3])*Ti + b[2])*Ti + b[1])*Ti
diff --git a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.H b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.H
index 7eaf936ec5..2c24af42c4 100644
--- a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.H
+++ b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/greyMeanAbsorptionEmission/greyMeanAbsorptionEmission.H
@@ -100,7 +100,7 @@ namespace radiation
 {
 
 /*---------------------------------------------------------------------------*\
-                  Class constantAbsorptionEmission Declaration
+                 Class greyMeanAbsorptionEmission Declaration
 \*---------------------------------------------------------------------------*/
 
 class greyMeanAbsorptionEmission
diff --git a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
index 6a3f290706..a96deeb0b4 100644
--- a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
+++ b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
@@ -80,7 +80,7 @@ Foam::radiation::wideBandAbsorptionEmission::wideBandAbsorptionEmission
         const dictionary& dict = iter().dict();
         dict.lookup("bandLimits") >> iBands_[nBand];
         dict.lookup("EhrrCoeff") >> iEhrrCoeffs_[nBand];
-        totalWaveLength_ += (iBands_[nBand][1] - iBands_[nBand][0]);
+        totalWaveLength_ += iBands_[nBand][1] - iBands_[nBand][0];
 
         label nSpec = 0;
         const dictionary& specDicts = dict.subDict("species");
@@ -111,7 +111,7 @@ Foam::radiation::wideBandAbsorptionEmission::wideBandAbsorptionEmission
     nBands_ = nBand;
 
     // Check that all the species on the dictionary are present in the
-    // look-up table  and save the corresponding indexes of the look-up table
+    // look-up table  and save the corresponding indices of the look-up table
 
     label j = 0;
     forAllConstIter(HashTable<label>, speciesNames_, iter)
@@ -180,7 +180,7 @@ Foam::radiation::wideBandAbsorptionEmission::aCont(const label bandI) const
                 IOobject::NO_WRITE
             ),
             mesh(),
-            dimensionedScalar("a",dimless/dimLength, 0.0)
+            dimensionedScalar("a", dimless/dimLength, 0.0)
         )
     );
 
@@ -193,7 +193,7 @@ Foam::radiation::wideBandAbsorptionEmission::aCont(const label bandI) const
         for (label n=0; n<nSpecies; n++)
         {
             label l = 0;
-            scalar Yipi = 0;
+            scalar Yipi = 0.0;
             if (specieIndex_[n] != 0)
             {
                 // moles x pressure [atm]
@@ -213,7 +213,7 @@ Foam::radiation::wideBandAbsorptionEmission::aCont(const label bandI) const
 
             if (coeffs_[n][bandI].invTemp())
             {
-                Ti = 1./T[i];
+                Ti = 1.0/T[i];
             }
 
             a[i]+=
@@ -245,7 +245,7 @@ Foam::radiation::wideBandAbsorptionEmission::eCont(const label bandI) const
                 IOobject::NO_WRITE
             ),
             mesh(),
-            dimensionedScalar("e",dimless/dimLength, 0.0)
+            dimensionedScalar("e", dimless/dimLength, 0.0)
         )
     );
 
@@ -275,7 +275,7 @@ Foam::radiation::wideBandAbsorptionEmission::ECont(const label bandI) const
 
     if (mesh().foundObject<volScalarField>("hrr"))
     {
-        const volScalarField& hrr =  mesh().lookupObject<volScalarField>("hrr");
+        const volScalarField& hrr = mesh().lookupObject<volScalarField>("hrr");
         E().internalField() =
             iEhrrCoeffs_[bandI]
            *hrr.internalField()
@@ -305,7 +305,7 @@ void Foam::radiation::wideBandAbsorptionEmission::correct
 {
     a = dimensionedScalar("zero", dimless/dimLength, 0.0);
 
-    for (label j = 0; j < nBands_; j++)
+    for (label j=0; j<nBands_; j++)
     {
         Info<< "Calculating absorption in band: " << j << endl;
         aLambda[j].internalField() = this->a(j);
@@ -317,4 +317,6 @@ void Foam::radiation::wideBandAbsorptionEmission::correct
     }
 
 }
+
+
 // ************************************************************************* //
diff --git a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.H b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.H
index d0763595ea..21e10c400b 100644
--- a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.H
+++ b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.H
@@ -36,9 +36,6 @@ Description
     There is no check of continuity of the bands. They should not ovelap or
     have gaps.
 
-    The black body emission power table(constant/blackBodyEmissivePower) range
-    of lambda * T = [1000; 10000] x 10E-6 (90% of the total emission).
-
     The emission constant proportionality is specified per band (EhrrCoeff).
 
     The coefficients for the species in the LookUpTable have to be specified
@@ -153,13 +150,13 @@ private:
         HashTable<label> speciesNames_;
 
         //- Indices of species in the look-up table
-        FixedList<label,nSpecies_> specieIndex_;
+        FixedList<label, nSpecies_> specieIndex_;
 
         //- Bands
-        FixedList<Vector2D<scalar>, maxBands_ > iBands_;
+        FixedList<Vector2D<scalar>, maxBands_> iBands_;
 
         //- Proportion of the heat released rate emitted
-        FixedList<scalar, maxBands_ > iEhrrCoeffs_;
+        FixedList<scalar, maxBands_> iEhrrCoeffs_;
 
         //- Look-up table of species related to ft
         mutable interpolationLookUpTable<scalar> lookUpTable_;
@@ -168,7 +165,7 @@ private:
         const basicThermo& thermo_;
 
         //- Bands
-        label nBands_ ;
+        label nBands_;
 
         //- Pointer list of species being solved involved in the absorption
         UPtrList<volScalarField> Yj_;
diff --git a/src/topoChangerFvMesh/movingConeTopoFvMesh/movingConeTopoFvMesh.C b/src/topoChangerFvMesh/movingConeTopoFvMesh/movingConeTopoFvMesh.C
index 2f3804afde..dedf725e24 100644
--- a/src/topoChangerFvMesh/movingConeTopoFvMesh/movingConeTopoFvMesh.C
+++ b/src/topoChangerFvMesh/movingConeTopoFvMesh/movingConeTopoFvMesh.C
@@ -343,7 +343,7 @@ bool Foam::movingConeTopoFvMesh::update()
             //    << endl;
 
             {
-                OFstream str(thisDb().path()/"meshPoints.obj");
+                OFstream str(time().timePath()/"meshPoints.obj");
                 Pout<< "Writing mesh with meshPoints to " << str.name()
                     << endl;
 
@@ -361,7 +361,7 @@ bool Foam::movingConeTopoFvMesh::update()
                 }
             }
             {
-                OFstream str(thisDb().path()/"preMotionPoints.obj");
+                OFstream str(time().timePath()/"preMotionPoints.obj");
                 Pout<< "Writing mesh with preMotionPoints to " << str.name()
                     << endl;
 
diff --git a/tutorials/incompressible/pimpleFoam/t-junction/0/epsilon b/tutorials/incompressible/pimpleFoam/t-junction/0/epsilon
index 2c4cfab1c0..5b446a2d8a 100644
--- a/tutorials/incompressible/pimpleFoam/t-junction/0/epsilon
+++ b/tutorials/incompressible/pimpleFoam/t-junction/0/epsilon
@@ -24,7 +24,7 @@ boundaryField
     inlet
     {
         type            turbulentMixingLengthDissipationRateInlet;
-        mixingLength    0.01;
+        mixingLength    0.01;       // 1cm - half channel height
         value           uniform 1;
     }
 
diff --git a/tutorials/incompressible/pimpleFoam/t-junction/0/k b/tutorials/incompressible/pimpleFoam/t-junction/0/k
index 51c0ab7cd3..ad217bc9e4 100644
--- a/tutorials/incompressible/pimpleFoam/t-junction/0/k
+++ b/tutorials/incompressible/pimpleFoam/t-junction/0/k
@@ -24,7 +24,7 @@ boundaryField
     inlet
     {
         type            turbulentIntensityKineticEnergyInlet;
-        intensity       0.05;
+        intensity       0.05;       // 5% turbulent intensity
         value           uniform 1;
     }
 
diff --git a/tutorials/incompressible/pimpleFoam/t-junction/0/p b/tutorials/incompressible/pimpleFoam/t-junction/0/p
index b3d4646058..2446eb2f7c 100644
--- a/tutorials/incompressible/pimpleFoam/t-junction/0/p
+++ b/tutorials/incompressible/pimpleFoam/t-junction/0/p
@@ -22,8 +22,14 @@ boundaryField
 {
     inlet      
     {
-        type            totalPressure;
-        p0              uniform 100040;
+        //type            totalPressure;
+        //p0              uniform 100040;
+
+        type            timeVaryingTotalPressure;
+        p0              100040; // only used for restarts
+        outOfBounds     clamp;
+        fileName        "$FOAM_CASE/constant/p0vsTime";
+
         U               U;
         phi             phi;
         rho             none;
diff --git a/tutorials/incompressible/pimpleFoam/t-junction/constant/p0vsTime b/tutorials/incompressible/pimpleFoam/t-junction/constant/p0vsTime
new file mode 100644
index 0000000000..a9ee6bb3e7
--- /dev/null
+++ b/tutorials/incompressible/pimpleFoam/t-junction/constant/p0vsTime
@@ -0,0 +1,4 @@
+(
+    (0 100010)
+    (1 100040)
+)
diff --git a/tutorials/incompressible/pimpleFoam/t-junction/constant/turbulenceProperties b/tutorials/incompressible/pimpleFoam/t-junction/constant/turbulenceProperties
new file mode 100644
index 0000000000..a18cc53b5b
--- /dev/null
+++ b/tutorials/incompressible/pimpleFoam/t-junction/constant/turbulenceProperties
@@ -0,0 +1,20 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  1.5                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "constant";
+    object      turbulenceProperties;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+simulationType RASModel;
+
+// ************************************************************************* //
diff --git a/tutorials/incompressible/pimpleFoam/t-junction/system/controlDict b/tutorials/incompressible/pimpleFoam/t-junction/system/controlDict
index a5ba827eec..2d7aaa81da 100644
--- a/tutorials/incompressible/pimpleFoam/t-junction/system/controlDict
+++ b/tutorials/incompressible/pimpleFoam/t-junction/system/controlDict
@@ -52,16 +52,24 @@ functions
     probes
     {
         type        probes;
-        name        probes;
+        // Where to load it from
         functionObjectLibs ( "libsampling.so" );
+        // Name of the directory for probe data
+        name        probes;
         probeLocations
         (
-            ( 1e-06 0 0.01 )
-            ( 0.21 -0.20999 0.01 )
-            ( 0.21 0.20999 0.01 )
-            ( 0.21 0 0.01 )
+            ( 1e-06 0 0.01 )        // at inlet
+            ( 0.21 -0.20999 0.01 )  // at outlet1
+            ( 0.21 0.20999 0.01 )   // at outlet2
+            ( 0.21 0 0.01 )         // at central block
         );
+
+        // Fields to be probed
         fields ( p U );
+
+        // Write at same frequency as fields
+        outputControl   outputTime;
+        outputInterval  1;
     }
 }