diff --git a/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/Make/options b/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/Make/options
index 7782cc45ccb80be0afa7df6874b97d62551a200a..495defbeec48e847eb02e3d684c8b96de37139f0 100644
--- a/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/Make/options
+++ b/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/Make/options
@@ -9,7 +9,7 @@ EXE_INC = \
     -I$(LIB_SRC)/thermophysicalModels/basicSolidThermo/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/solidChemistryModel/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/combustionModels/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/properties/solidProperties/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/properties/solidMixtureProperties/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/properties/liquidProperties/lnInclude \
diff --git a/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/YhsEqn.H b/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/YhsEqn.H
index 7606d9ce727e04c4e69d949ea2f425fa99c8384f..f985a650b9bcc3316f881294058dacb4884e7ab8 100644
--- a/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/YhsEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/YhsEqn.H
@@ -9,6 +9,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
     )
 );
 {
+    radiation->correct();
     combustion->correct();
     dQ = combustion->dQ();
     label inertIndex = -1;
@@ -65,7 +66,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 
     thermo.correct();
 
-    radiation->correct();
+    //radiation->correct();
 
     Info<< "min/max(T) = " << min(T).value() << ", " << max(T).value() << endl;
 }
diff --git a/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/reactingParcelFilmPyrolysisFoam.C b/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/reactingParcelFilmPyrolysisFoam.C
index ac029bea8642464f79f80b037ee791504bc34a92..178ceadcd62a02f137e75fd554930776c25b335c 100644
--- a/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/reactingParcelFilmPyrolysisFoam.C
+++ b/applications/solvers/lagrangian/reactingParcelFilmPyrolysisFoam/reactingParcelFilmPyrolysisFoam.C
@@ -30,6 +30,7 @@ Description
 
 \*---------------------------------------------------------------------------*/
 
+#include "mapDistribute.H"
 #include "fvCFD.H"
 #include "turbulenceModel.H"
 #include "basicReactingCloud.H"
@@ -54,9 +55,9 @@ int main(int argc, char *argv[])
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "createClouds.H"
-    #include "createRadiationModel.H"
     #include "createSurfaceFilmModel.H"
     #include "createPyrolysisModel.H"
+    #include "createRadiationModel.H"
     #include "initContinuityErrs.H"
     #include "readTimeControls.H"
     #include "compressibleCourantNo.H"
@@ -111,10 +112,8 @@ int main(int argc, char *argv[])
 
             rho = thermo.rho();
         }
-        else
-        {
-            runTime.write();
-        }
+
+        runTime.write();
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
             << "  ClockTime = " << runTime.elapsedClockTime() << " s"
diff --git a/src/thermophysicalModels/radiationModels/derivedFvPatchFields/greyDiffusiveRadiation/greyDiffusiveRadiationMixedFvPatchScalarField.C b/src/thermophysicalModels/radiationModels/derivedFvPatchFields/greyDiffusiveRadiation/greyDiffusiveRadiationMixedFvPatchScalarField.C
index e231144b9208021b21bebff47cc170a4b256e0cf..90727e381bd1e118b279af75a8ad46cab3e01957 100644
--- a/src/thermophysicalModels/radiationModels/derivedFvPatchFields/greyDiffusiveRadiation/greyDiffusiveRadiationMixedFvPatchScalarField.C
+++ b/src/thermophysicalModels/radiationModels/derivedFvPatchFields/greyDiffusiveRadiation/greyDiffusiveRadiationMixedFvPatchScalarField.C
@@ -189,6 +189,8 @@ updateCoeffs()
 
     ray.Qr().boundaryField()[patchI] += Iw*(n & ray.dAve());
 
+    scalarList temissivity = emissivity();
+
     forAll(Iw, faceI)
     {
         scalar Ir = 0.0;
@@ -217,8 +219,8 @@ updateCoeffs()
             valueFraction()[faceI] = 1.0;
             refValue()[faceI] =
                 (
-                    Ir*(scalar(1.0) - emissivity()()[faceI])
-                  + emissivity()()[faceI]*physicoChemical::sigma.value()
+                    Ir*(scalar(1.0) - temissivity[faceI])
+                  + temissivity[faceI]*physicoChemical::sigma.value()
                   * pow4(Tp[faceI])
                 )/pi;
 
@@ -238,7 +240,6 @@ updateCoeffs()
                 Iw[faceI]*(n[faceI] & ray.dAve());
         }
     }
-
     mixedFvPatchScalarField::updateCoeffs();
 }
 
diff --git a/src/thermophysicalModels/radiationModels/derivedFvPatchFields/radiationCoupledBase/radiationCoupledBase.C b/src/thermophysicalModels/radiationModels/derivedFvPatchFields/radiationCoupledBase/radiationCoupledBase.C
index f8d034ca34963952df97d50bd540e15c05e87379..580342d92a50ae650ec6bb653fcc5c9b0ad44d2b 100644
--- a/src/thermophysicalModels/radiationModels/derivedFvPatchFields/radiationCoupledBase/radiationCoupledBase.C
+++ b/src/thermophysicalModels/radiationModels/derivedFvPatchFields/radiationCoupledBase/radiationCoupledBase.C
@@ -127,7 +127,7 @@ Foam::radiationCoupledBase::radiationCoupledBase
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::scalarField> Foam::radiationCoupledBase::emissivity() const
+Foam::scalarField Foam::radiationCoupledBase::emissivity() const
 {
     switch (method_)
     {
@@ -142,6 +142,9 @@ Foam::tmp<Foam::scalarField> Foam::radiationCoupledBase::emissivity() const
 
             const polyMesh& nbrMesh = mpp.sampleMesh();
 
+            // Force recalculation of mapping and schedule
+            const mapDistribute& distMap = mpp.map();
+
             const fvPatch& nbrPatch = refCast<const fvMesh>
             (
                 nbrMesh
@@ -160,12 +163,11 @@ Foam::tmp<Foam::scalarField> Foam::radiationCoupledBase::emissivity() const
                     )
                 );
 
-                scalarField& emissivity = temissivity();
-
+                scalarField emissivity(temissivity);
                 // Use direct map mapping to exchange data
-                mpp.map().distribute(emissivity);
-
-                return temissivity;
+                distMap.distribute(emissivity);
+                //Pout << emissivity << endl;
+                return emissivity;
             }
             else
             {
diff --git a/src/thermophysicalModels/radiationModels/derivedFvPatchFields/radiationCoupledBase/radiationCoupledBase.H b/src/thermophysicalModels/radiationModels/derivedFvPatchFields/radiationCoupledBase/radiationCoupledBase.H
index 405d0f75a54ef6723d169f997bf1cdaaadf9d706..0b9fe86c74a618fd28976be156856cd6389e8c9a 100644
--- a/src/thermophysicalModels/radiationModels/derivedFvPatchFields/radiationCoupledBase/radiationCoupledBase.H
+++ b/src/thermophysicalModels/radiationModels/derivedFvPatchFields/radiationCoupledBase/radiationCoupledBase.H
@@ -114,7 +114,7 @@ public:
 
 
         //- Calculate corresponding emissivity field
-        tmp<scalarField> emissivity() const;
+        scalarField emissivity() const;
 
         //- Write
         void write(Ostream&) const;
diff --git a/src/thermophysicalModels/radiationModels/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.C b/src/thermophysicalModels/radiationModels/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.C
index 45ab196b5c0780f5bae069d5aa06d97e2e2393dc..c985ea5e3a708b0126073f4a15415ba15c03dcfd 100644
--- a/src/thermophysicalModels/radiationModels/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.C
+++ b/src/thermophysicalModels/radiationModels/derivedFvPatchFields/wideBandDiffusiveRadiation/wideBandDiffusiveRadiationMixedFvPatchScalarField.C
@@ -189,6 +189,8 @@ updateCoeffs()
         dom.blackBody().bLambda(lambdaId).boundaryField()[patchI]
     );
 
+    scalarList temissivity = emissivity();
+
     forAll(Iw, faceI)
     {
         scalar Ir = 0.0;
@@ -215,8 +217,8 @@ updateCoeffs()
             valueFraction()[faceI] = 1.0;
             refValue()[faceI] =
                 (
-                    Ir*(1.0 - emissivity()()[faceI])
-                  + emissivity()()[faceI]*Eb[faceI]
+                    Ir*(1.0 - temissivity[faceI])
+                  + temissivity[faceI]*Eb[faceI]
                 )/pi;
         }
         else
diff --git a/src/turbulenceModels/compressible/turbulenceModel/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.C b/src/turbulenceModels/compressible/turbulenceModel/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.C
index 262b9da61f0256a1efae6baa70f5835c63a7bf6d..b5f81301b593e428c6bc56d5523323252fe20ecd 100644
--- a/src/turbulenceModels/compressible/turbulenceModel/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.C
+++ b/src/turbulenceModels/compressible/turbulenceModel/derivedFvPatchFields/turbulentTemperatureRadCoupledMixed/turbulentTemperatureRadCoupledMixedFvPatchScalarField.C
@@ -158,6 +158,9 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::updateCoeffs()
     const fvPatch& nbrPatch =
         refCast<const fvMesh>(nbrMesh).boundary()[samplePatchI];
 
+    // Force recalculation of mapping and schedule
+    const mapDistribute& distMap = mpp.map();
+
     scalarField Tc(patchInternalField());
     scalarField& Tp = *this;
 
@@ -170,13 +173,13 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::updateCoeffs()
 
     // Swap to obtain full local values of neighbour internal field
     scalarField TcNbr(nbrField.patchInternalField());
+    distMap.distribute(TcNbr);
 
-    mpp.map().distribute(TcNbr);
 
     // Swap to obtain full local values of neighbour K*delta
-    scalarField KDeltaNbr(nbrField.K(TcNbr)*nbrPatch.deltaCoeffs());
+    scalarField KDeltaNbr(nbrField.K(nbrField)*nbrPatch.deltaCoeffs());
+    distMap.distribute(KDeltaNbr);
 
-    mpp.map().distribute(KDeltaNbr);
 
     scalarField KDelta(K(*this)*patch().deltaCoeffs());
 
@@ -190,7 +193,7 @@ void turbulentTemperatureRadCoupledMixedFvPatchScalarField::updateCoeffs()
     if (QrNbrName_ != "none")
     {
         QrNbr = nbrPatch.lookupPatchField<volScalarField, scalar>(QrNbrName_);
-        mpp.map().distribute(QrNbr);
+        distMap.distribute(QrNbr);
     }
 
     scalarField alpha(KDeltaNbr - (Qr + QrNbr)/Tp);
diff --git a/tutorials/heatTransfer/chtMultiRegionFoam/snappyMultiRegionHeater/constant/polyMesh/boundary b/tutorials/heatTransfer/chtMultiRegionFoam/snappyMultiRegionHeater/constant/polyMesh/boundary
deleted file mode 100644
index 39c75e0d212e04c31d992b8a6bc1e7be7caa6a45..0000000000000000000000000000000000000000
--- a/tutorials/heatTransfer/chtMultiRegionFoam/snappyMultiRegionHeater/constant/polyMesh/boundary
+++ /dev/null
@@ -1,58 +0,0 @@
-/*--------------------------------*- C++ -*----------------------------------*\
-| =========                 |                                                 |
-| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      www.OpenFOAM.com                      |
-|    \\/     M anipulation  |                                                 |
-\*---------------------------------------------------------------------------*/
-FoamFile
-{
-    version     2.0;
-    format      ascii;
-    class       polyBoundaryMesh;
-    location    "constant/polyMesh";
-    object      boundary;
-}
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-6
-(
-    maxY
-    {
-        type            wall;
-        nFaces          300;
-        startFace       8300;
-    }
-    minX
-    {
-        type            patch;
-        nFaces          100;
-        startFace       8600;
-    }
-    maxX
-    {
-        type            patch;
-        nFaces          100;
-        startFace       8700;
-    }
-    minY
-    {
-        type            wall;
-        nFaces          300;
-        startFace       8800;
-    }
-    minZ
-    {
-        type            wall;
-        nFaces          300;
-        startFace       9100;
-    }
-    maxZ
-    {
-        type            wall;
-        nFaces          300;
-        startFace       9400;
-    }
-)
-
-// ************************************************************************* //