From 0e2b526cfdabb9230c9eef142c204d67c6cbacc2 Mon Sep 17 00:00:00 2001
From: andy <a.heather@opencfd.co.uk>
Date: Fri, 6 Mar 2009 13:50:17 +0000
Subject: [PATCH] more updates

---
 .../Templates/KinematicCloud/KinematicCloud.H |   5 +
 .../KinematicCloud/KinematicCloudI.H          |   8 +
 .../Templates/ReactingCloud/ReactingCloud.H   |   4 +
 .../Templates/ReactingCloud/ReactingCloudI.H  |   8 +
 .../ReactingMultiphaseCloud.H                 |   5 +
 .../ReactingMultiphaseCloudI.H                |   8 +
 .../Templates/ThermoCloud/ThermoCloud.H       |   4 +
 .../Templates/ThermoCloud/ThermoCloudI.H      |   8 +
 .../KinematicParcel/KinematicParcel.C         |  19 +-
 .../ReactingMultiphaseParcel.C                | 196 +++++++++++-------
 .../ReactingMultiphaseParcel.H                |  18 +-
 .../ReactingMultiphaseParcelIO.C              |  71 +++++--
 .../Templates/ReactingParcel/ReactingParcel.C | 119 ++++++-----
 .../Templates/ReactingParcel/ReactingParcel.H |   8 +-
 .../ReactingParcel/ReactingParcelI.H          |  14 +-
 .../ReactingParcel/ReactingParcelIO.C         |  63 ++++--
 .../Templates/ThermoParcel/ThermoParcel.C     |  41 ++--
 .../phaseProperties/phaseProperties.C         |  63 +++++-
 .../phaseProperties/phaseProperties.H         |   9 +
 .../phaseProperties/phasePropertiesIO.C       |  15 +-
 .../phasePropertiesList/phasePropertiesList.C |  13 +-
 .../phasePropertiesList/phasePropertiesList.H |   6 +
 .../IO/DataEntry/DataEntry/NewDataEntry.C     |   2 +-
 .../DispersionModel/NewDispersionModel.C      |  12 +-
 .../DragModel/DragModel/NewDragModel.C        |  15 +-
 .../InjectionModel/InjectionModel.C           |  32 +--
 .../InjectionModel/NewInjectionModel.C        |  12 +-
 .../NewWallInteractionModel.C                 |  15 +-
 .../CompositionModel/CompositionModel.C       | 127 +++++++++---
 .../CompositionModel/CompositionModel.H       |  14 +-
 .../CompositionModel/NewCompositionModel.C    |  13 +-
 .../SingleMixtureFraction.C                   |  20 +-
 .../LiquidEvaporation/LiquidEvaporation.C     | 102 +++++----
 .../LiquidEvaporation/LiquidEvaporation.H     |   8 +-
 .../PhaseChangeModel/NewPhaseChangeModel.C    |  13 +-
 .../NewDevolatilisationModel.C                |  15 +-
 .../NoSurfaceReaction/NoSurfaceReaction.C     |   4 +-
 .../NewSurfaceReactionModel.C                 |  15 +-
 .../HeatTransferModel/NewHeatTransferModel.C  |  15 +-
 39 files changed, 732 insertions(+), 407 deletions(-)

diff --git a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/KinematicCloud.H b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/KinematicCloud.H
index 7d05141fca0..8a03a73cd7e 100644
--- a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/KinematicCloud.H
+++ b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/KinematicCloud.H
@@ -216,6 +216,11 @@ public:
                 //- Return particle properties dictionary
                 inline const IOdictionary& particleProperties() const;
 
+                //- Return the constant properties
+                inline const typename ParcelType::constantProperties&
+                    constProps() const;
+
+
             //- Return coupled flag
             inline const Switch coupled() const;
 
diff --git a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/KinematicCloudI.H b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/KinematicCloudI.H
index 948a6b5dfff..f129df9bdb0 100644
--- a/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/KinematicCloudI.H
+++ b/src/lagrangian/intermediate/clouds/Templates/KinematicCloud/KinematicCloudI.H
@@ -55,6 +55,14 @@ Foam::KinematicCloud<ParcelType>::particleProperties() const
 }
 
 
+template<class ParcelType>
+inline const typename ParcelType::constantProperties&
+Foam::KinematicCloud<ParcelType>::constProps() const
+{
+    return constProps_;
+}
+
+
 template<class ParcelType>
 inline const Foam::Switch Foam::KinematicCloud<ParcelType>::coupled() const
 {
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.H b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.H
index 3af70b1904e..6fe5bf68a65 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.H
@@ -143,6 +143,10 @@ public:
 
         // Access
 
+            //- Return the constant properties
+            inline const typename ParcelType::constantProperties&
+                constProps() const;
+
             //- Return const access to carrier phase thermo package
             inline const hCombustionThermo& carrierThermo() const;
 
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloudI.H b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloudI.H
index 4bf20a37a9c..1414b6b16d7 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloudI.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloudI.H
@@ -26,6 +26,14 @@ License
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class ParcelType>
+inline const typename ParcelType::constantProperties&
+Foam::ReactingCloud<ParcelType>::constProps() const
+{
+    return constProps_;
+}
+
+
 template<class ParcelType>
 inline const Foam::hCombustionThermo&
 Foam::ReactingCloud<ParcelType>::carrierThermo() const
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.H b/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.H
index b78a5c02426..ea3f585313a 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.H
@@ -137,6 +137,11 @@ public:
 
         // Access
 
+            //- Return the constant properties
+            inline const typename ParcelType::constantProperties&
+                constProps() const;
+
+
             // Sub-models
 
                 //- Return reference to devolatilisation model
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloudI.H b/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloudI.H
index 752f2b8ca32..cbf042f688f 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloudI.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloudI.H
@@ -26,6 +26,14 @@ License
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class ParcelType>
+inline const typename ParcelType::constantProperties&
+Foam::ReactingMultiphaseCloud<ParcelType>::constProps() const
+{
+    return constProps_;
+}
+
+
 template<class ParcelType>
 inline const Foam::DevolatilisationModel
 <
diff --git a/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.H b/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.H
index d845b725b19..9789bd0f5f0 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.H
@@ -142,6 +142,10 @@ public:
 
         // Access
 
+            //- Return the constant properties
+            inline const typename ParcelType::constantProperties&
+                constProps() const;
+
             //- Return const access to thermo package
             inline const basicThermo& carrierThermo() const;
 
diff --git a/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloudI.H b/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloudI.H
index 94f12b85a2c..2f3df96f41b 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloudI.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloudI.H
@@ -28,6 +28,14 @@ License
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class ParcelType>
+inline const typename ParcelType::constantProperties&
+Foam::ThermoCloud<ParcelType>::constProps() const
+{
+    return constProps_;
+}
+
+
 template<class ParcelType>
 inline const Foam::basicThermo&
 Foam::ThermoCloud<ParcelType>::carrierThermo() const
diff --git a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
index 32d40fc5923..96a777b1936 100644
--- a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
@@ -64,20 +64,19 @@ void Foam::KinematicParcel<ParcelType>::calc
     const label cellI
 )
 {
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Calculate velocity - update U, dUTrans
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 1. Calculate velocity - update U, dUTrans
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     scalar Cud = 0.0;
     vector dUTrans = vector::zero;
     const vector U1 = calcVelocity(td, dt, Cud, dUTrans);
 
 
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 2. Accumulate carrier phase source terms
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     if (td.cloud().coupled())
     {
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        // Accumulate carrier phase source terms
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
         // Update momentum transfer
         td.cloud().UTrans()[cellI] += nParticle_*dUTrans;
 
@@ -86,9 +85,9 @@ void Foam::KinematicParcel<ParcelType>::calc
     }
 
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Set new particle properties
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 3. Set new particle properties
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     this->U() = U1;
 }
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
index 4d96e786018..6bb568409cc 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
@@ -26,6 +26,40 @@ License
 
 #include "ReactingMultiphaseParcel.H"
 
+// * * * * * * * * * * * *  Private Member Functions * * * * * * * * * * * * //
+
+template<class ParcelType>
+template<class TrackData>
+Foam::scalar Foam::ReactingMultiphaseParcel<ParcelType>::cpEff
+(
+    TrackData& td,
+    const scalar p,
+    const scalar T
+) const
+{
+    return
+        this->Y_[0]*td.cloud().composition().cp(0, YGas_, p, T)
+      + this->Y_[1]*td.cloud().composition().cp(0, YLiquid_, p, T)
+      + this->Y_[2]*td.cloud().composition().cp(0, YSolid_, p, T);
+}
+
+
+template<class ParcelType>
+template<class TrackData>
+Foam::scalar Foam::ReactingMultiphaseParcel<ParcelType>::HEff
+(
+    TrackData& td,
+    const scalar p,
+    const scalar T
+) const
+{
+    return
+        this->Y_[0]*td.cloud().composition().H(0, YGas_, p, T)
+      + this->Y_[1]*td.cloud().composition().H(0, YLiquid_, p, T)
+      + this->Y_[2]*td.cloud().composition().H(0, YSolid_, p, T);
+}
+
+
 // * * * * * * * * * * *  Protected Member Functions * * * * * * * * * * * * //
 
 template<class ParcelType>
@@ -57,8 +91,10 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
     const scalar mass0 = this->mass();
     const scalar cp0 = this->cp_;
     const scalar np0 = this->nParticle_;
-    const scalar T0 = this->T_;
-    scalarField& YMix = this->Y_;
+    scalar T0 = this->T_;
+    const scalar pc = this->pc_;
+    scalarField& Y = this->Y_;
+    scalar mass1 = mass0;
 
     label idGas = td.cloud().composition().idGas();
     label idLiquid = td.cloud().composition().idLiquid();
@@ -97,16 +133,60 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
     scalar T1 = calcHeatTransfer(td, dt, cellI, htc, dhTrans);
 
 
-    // ~~~~~~~~~~~~~~~~~~~~~~
-    // Calculate phase change
-    // ~~~~~~~~~~~~~~~~~~~~~~
-    scalar dMassPC = calcPhaseChange(td, dt, cellI, T1, dMassMT);
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // Calculate phase change - update mass, Y, cp, T, dhTrans
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    const scalar dMassPC = calcPhaseChange(td, dt, cellI, T1, Y[idLiquid], dMassMT);
+
+    // Update particle mass
+    mass1 -= dMassPC;
+
+    // Update particle liquid component mass fractions
+    this->updateMassFraction(mass0, dMassMT, YLiquid_);
+
+    // New specific heat capacity of mixture
+    scalar cp1 = cpEff(td, pc, T1);
+
+    // Correct temperature due to evaporated components
+    // TODO: use hl function in liquidMixture???
+//    scalar dhPC = -dMassPCTot*td.cloud().composition().L(0, Y_, pc, T0);
+    scalar Lvap = td.cloud().composition().L(idLiquid, this->YLiquid_, pc, T0);
+    T1 -= Lvap*dMassPC/(0.5*(mass0 + mass1)*cp1);
+
+    // Correct dhTrans to account for the change in enthalpy due to the
+    // liquid phase change
+    dhTrans +=
+        dMassPC
+       *td.cloud().composition().H(idLiquid, YLiquid_, pc, 0.5*(T0 + T1));
+
+//????????????????????????????????????????????????????????????????????????????????
+    // Store temperature for the start of the next process
+    T0 = T1;
+//????????????????????????????????????????????????????????????????????????????????
 
 
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~
     // Calculate Devolatilisation
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~
-    calcDevolatilisation(td, dt, T0, T1, dMassMT);
+    const scalar dMassD = calcDevolatilisation(td, dt, T0, dMassMT);
+
+    // Update particle mass
+    mass1 -= dMassD;
+
+    // New specific heat capacity of mixture
+    cp1 = cpEff(td, pc, T1);
+
+    // Update gas and solid component mass fractions
+//    updateMassFraction(mass0, mass1, dMassMT, YGas_);
+//    updateMassFraction(mass0, mass1, dMassMT, YSolid_);
+
+    // Correct particle temperature to account for latent heat of
+    // devolatilisation
+    T1 -=
+        td.constProps().Ldevol()
+       *sum(dMassMT)
+       /(0.5*(mass0 + mass1)*cp1);
+
 
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
     // Calculate surface reactions
@@ -118,14 +198,8 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
     calcSurfaceReactions(td, dt, cellI, T0, T1, dMassMT, dMassSR, dhRet);
 
     // New total mass
-    const scalar mass1 = mass0 - sum(dMassMT) - sum(dMassSR);
+    mass1 -= sum(dMassSR);
 
-    // Correct particle temperature to account for latent heat of
-    // devolatilisation
-    T1 -=
-        td.constProps().Ldevol()
-       *sum(dMassMT)
-       /(0.5*(mass0 + mass1)*cp0);
 
     // Enthalpy retention divided between particle and carrier phase by the
     // fraction fh and (1 - fh)
@@ -135,36 +209,39 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
     // Correct dhTrans to account for enthalpy of evolved volatiles
     dhTrans +=
         sum(dMassMT)
-       *td.cloud().composition().H(idGas, YGas_, this->pc_, 0.5*(T0 + T1));
+       *td.cloud().composition().H(idGas, YGas_, pc, 0.5*(T0 + T1));
 
     // Correct dhTrans to account for enthalpy of consumed solids
     dhTrans +=
         sum(dMassSR)
-       *td.cloud().composition().H(idSolid, YSolid_, this->pc_, 0.5*(T0 + T1));
-
+       *td.cloud().composition().H(idSolid, YSolid_, pc, 0.5*(T0 + T1));
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~
-    // Accumulate source terms
-    // ~~~~~~~~~~~~~~~~~~~~~~~
 
-    // Transfer mass lost from particle to carrier mass source
-    forAll(dMassMT, i)
+    if (td.cloud().coupled())
     {
-        td.cloud().rhoTrans(i)[cellI] += np0*(dMassMT[i] + dMassSR[i]);
-    }
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        // Accumulate carrier phase source terms
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
-    // Update momentum transfer
-    td.cloud().UTrans()[cellI] += np0*dUTrans;
+        // Transfer mass lost from particle to carrier mass source
+        forAll(dMassMT, i)
+        {
+            td.cloud().rhoTrans(i)[cellI] += np0*(dMassMT[i] + dMassSR[i]);
+        }
 
-    // Accumulate coefficient to be applied in carrier phase momentum coupling
-    td.cloud().UCoeff()[cellI] += np0*mass0*Cud;
+        // Update momentum transfer
+        td.cloud().UTrans()[cellI] += np0*dUTrans;
 
-    // Update enthalpy transfer
-    // - enthalpy of lost solids already accounted for
-    td.cloud().hTrans()[cellI] += np0*dhTrans;
+        // Coefficient to be applied in carrier phase momentum coupling
+        td.cloud().UCoeff()[cellI] += np0*mass0*Cud;
 
-    // Accumulate coefficient to be applied in carrier phase enthalpy coupling
-    td.cloud().hCoeff()[cellI] += np0*htc*this->areaS();
+        // Update enthalpy transfer
+        // - enthalpy of lost solids already accounted for
+        td.cloud().hTrans()[cellI] += np0*dhTrans;
+
+        // Coefficient to be applied in carrier phase enthalpy coupling
+        td.cloud().hCoeff()[cellI] += np0*htc*this->areaS();
+    }
 
 
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
@@ -179,14 +256,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
         {
             td.cloud().rhoTrans(i)[cellI] += np0*dMassMT[i];
         }
-        td.cloud().hTrans()[cellI] +=
-            np0*mass1
-           *(
-                YMix[0]
-               *td.cloud().composition().H(idGas, YGas_, this->pc_, T1)
-              + YMix[2]
-               *td.cloud().composition().H(idSolid, YSolid_, this->pc_, T1)
-            );
+        td.cloud().hTrans()[cellI] += np0*mass1*HEff(td, pc, T1);
         td.cloud().UTrans()[cellI] += np0*mass1*U1;
     }
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
@@ -196,10 +266,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
     {
         this->U_ = U1;
         this->T_ = T1;
-        this->cp_ =
-            YMix[0]*td.cloud().composition().cp(idGas, YGas_, this->pc_, T1)
-          + YMix[1]*td.cloud().composition().cp(idLiquid, YLiquid_, this->pc_, T1)
-          + YMix[2]*td.cloud().composition().cp(idSolid, YSolid_, this->pc_, T1);
+        this->cp_ = cpEff(td, pc, T1);
 
         // Update particle density or diameter
         if (td.constProps().constantVolume())
@@ -216,33 +283,14 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
 
 template<class ParcelType>
 template<class TrackData>
-void Foam::ReactingMultiphaseParcel<ParcelType>::calcDevolatilisation
+Foam::scalar Foam::ReactingMultiphaseParcel<ParcelType>::calcDevolatilisation
 (
     TrackData& td,
     const scalar dt,
-    const scalar T0,
-    const scalar T1,
+    const scalar T,
     scalarList& dMassMT
 )
 {
-    if (td.cloud().composition().YMixture0()[1]>SMALL)
-    {
-        notImplemented
-        (
-            "void Foam::ReactingMultiphaseParcel<ParcelType>::"
-            "calcDevolatilisation \n"
-            "(\n"
-            "    TrackData&,\n"
-            "    const scalar,\n"
-            "    const scalar,\n"
-            "    const scalar,\n"
-            "    scalarList&\n"
-            ")\n"
-            "no treatment currently available for particles containing "
-            "liquid species"
-        );
-    }
-
     // Check that model is active, and that the parcel temperature is
     // within necessary limits for devolatilisation to occur
     if
@@ -252,27 +300,27 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calcDevolatilisation
      || this->T_<td.constProps().Tbp()
     )
     {
-        return;
+        return 0.0;
     }
 
     // Determine mass to add to carrier phase
     const scalar mass = this->mass();
-    scalarField& YMix = this->Y_;
+    scalarField& Y = this->Y_;
     const scalar dMassTot = td.cloud().devolatilisation().calculate
     (
         dt,
         this->mass0_,
         mass,
         td.cloud().composition().YMixture0(),
-        YMix,
-        T0,
+        Y,
+        T,
         canCombust_
     );
 
     // Update (total) mass fractions
-    YMix[0] = (YMix[0]*mass - dMassTot)/(mass - dMassTot);
-    YMix[1] = YMix[1]*mass/(mass - dMassTot);
-    YMix[2] = 1.0 - YMix[0] - YMix[1];
+    Y[0] = (Y[0]*mass - dMassTot)/(mass - dMassTot);
+    Y[1] = Y[1]*mass/(mass - dMassTot);
+    Y[2] = 1.0 - Y[0] - Y[1];
 
     // Add to cummulative mass transfer
     label idGas = td.cloud().composition().idGas();
@@ -286,6 +334,8 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calcDevolatilisation
     }
 
     td.cloud().addToMassDevolatilisation(dMassTot);
+
+    return dMassTot;
 }
 
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.H b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.H
index d53fc2a6d02..e70f0ae5c30 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.H
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.H
@@ -143,6 +143,19 @@ public:
     };
 
 
+private:
+
+    // Private member functions
+
+        //- Return the mixture effective specific heat capacity
+        template<class TrackData>
+        scalar cpEff(TrackData& td, const scalar p, const scalar T) const;
+
+        //- Return the mixture effective enthalpy
+        template<class TrackData>
+        scalar HEff(TrackData& td, const scalar p, const scalar T) const;
+
+
 protected:
 
     // Protected data
@@ -167,12 +180,11 @@ protected:
 
         //- Calculate Devolatilisation
         template<class TrackData>
-        void calcDevolatilisation
+        scalar calcDevolatilisation
         (
             TrackData& td,
             const scalar dt,
-            const scalar T0,
-            const scalar T1,
+            const scalar T,
             scalarList& dMassMT
         );
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcelIO.C b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcelIO.C
index e94d3b9177a..61513ae7830 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcelIO.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcelIO.C
@@ -77,11 +77,11 @@ Foam::ReactingMultiphaseParcel<ParcelType>::ReactingMultiphaseParcel
     // Check state of Istream
     is.check
     (
-        "ReactingMultiphaseParcel<ParcelType>::ReactingMultiphaseParcel\n"
-        "(\n"
-        "    const Cloud<ParcelType>&,\n"
-        "    Istream&,\n"
-        "    bool\n"
+        "ReactingMultiphaseParcel<ParcelType>::ReactingMultiphaseParcel"
+        "("
+            "const Cloud<ParcelType>&, "
+            "Istream&, "
+            "bool"
         ")"
     );
 }
@@ -102,11 +102,12 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::readFields
 
     // Get names and sizes for each Y...
     const label idGas = c.composition().idGas();
-    const wordList gasNames = c.composition().componentNames(idGas);
+    const wordList& gasNames = c.composition().componentNames(idGas);
     const label idLiquid = c.composition().idLiquid();
-    const wordList liquidNames = c.composition().componentNames(idLiquid);
+    const wordList& liquidNames = c.composition().componentNames(idLiquid);
     const label idSolid = c.composition().idSolid();
-    const wordList solidNames = c.composition().componentNames(idSolid);
+    const wordList& solidNames = c.composition().componentNames(idSolid);
+    const wordList& stateLabels = c.composition().stateLabels();
 
     // Set storage for each Y... for each parcel
     forAllIter(typename Cloud<ParcelType>, c, iter)
@@ -122,7 +123,11 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::readFields
     {
         IOField<scalar> YGas
         (
-            c.fieldIOobject("Y" + gasNames[j], IOobject::MUST_READ)
+            c.fieldIOobject
+            (
+                "Y" + gasNames[j] + stateLabels[idGas],
+                IOobject::MUST_READ
+            )
         );
 
         label i = 0;
@@ -137,7 +142,11 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::readFields
     {
         IOField<scalar> YLiquid
         (
-            c.fieldIOobject("Y" + liquidNames[j], IOobject::MUST_READ)
+            c.fieldIOobject
+            (
+                "Y" + liquidNames[j] + stateLabels[idLiquid],
+                 IOobject::MUST_READ
+            )
         );
 
         label i = 0;
@@ -152,7 +161,11 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::readFields
     {
         IOField<scalar> YSolid
         (
-            c.fieldIOobject("Y" + solidNames[j], IOobject::MUST_READ)
+            c.fieldIOobject
+            (
+                "Y" + solidNames[j] + stateLabels[idSolid],
+                IOobject::MUST_READ
+            )
         );
 
         label i = 0;
@@ -178,13 +191,19 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::writeFields
     // Write the composition fractions
     if (np > 0)
     {
+        const wordList& stateLabels = c.composition().stateLabels();
+
         const label idGas = c.composition().idGas();
-        const wordList gasNames = c.composition().componentNames(idGas);
+        const wordList& gasNames = c.composition().componentNames(idGas);
         forAll(gasNames, j)
         {
             IOField<scalar> YGas
             (
-                c.fieldIOobject("Y" + gasNames[j], IOobject::NO_READ),
+                c.fieldIOobject
+                (
+                    "Y" + gasNames[j] + stateLabels[idGas],
+                    IOobject::NO_READ
+                ),
                 np
             );
 
@@ -197,13 +216,18 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::writeFields
 
             YGas.write();
         }
+
         const label idLiquid = c.composition().idLiquid();
-        const wordList liquidNames = c.composition().componentNames(idLiquid);
+        const wordList& liquidNames = c.composition().componentNames(idLiquid);
         forAll(liquidNames, j)
         {
             IOField<scalar> YLiquid
             (
-                c.fieldIOobject("Y" + liquidNames[j], IOobject::NO_READ),
+                c.fieldIOobject
+                (
+                    "Y" + liquidNames[j] + stateLabels[idLiquid],
+                    IOobject::NO_READ
+                ),
                 np
             );
 
@@ -216,13 +240,18 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::writeFields
 
             YLiquid.write();
         }
+
         const label idSolid = c.composition().idSolid();
-        const wordList solidNames = c.composition().componentNames(idSolid);
+        const wordList& solidNames = c.composition().componentNames(idSolid);
         forAll(solidNames, j)
         {
             IOField<scalar> YSolid
             (
-                c.fieldIOobject("Y" + solidNames[j], IOobject::NO_READ),
+                c.fieldIOobject
+                (
+                    "Y" + solidNames[j] + stateLabels[idSolid],
+                    IOobject::NO_READ
+                ),
                 np
             );
 
@@ -267,10 +296,10 @@ Foam::Ostream& Foam::operator<<
     // Check state of Ostream
     os.check
     (
-        "Ostream& operator<<\n"
-        "(\n"
-        "    Ostream&,\n"
-        "    const ReactingMultiphaseParcel<ParcelType>&\n"
+        "Ostream& operator<<"
+        "("
+            "Ostream&, "
+            "const ReactingMultiphaseParcel<ParcelType>&"
         ")"
     );
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
index 6240d503bae..d605841d406 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
@@ -47,11 +47,12 @@ template<class ParcelType>
 void Foam::ReactingParcel<ParcelType>::updateMassFraction
 (
     const scalar mass0,
-    const scalar mass1,
     const scalarList& dMass,
     scalarField& Y
 )
 {
+    scalar mass1 = mass0 + sum(dMass);
+
     forAll(Y, i)
     {
         Y[i] = (Y[i]*mass0 - dMass[i])/mass1;
@@ -68,76 +69,76 @@ void Foam::ReactingParcel<ParcelType>::calc
     const label cellI
 )
 {
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Define local properties at beginning of timestep
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    const vector U0 = this->U_;
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // Define local properties at beginning of time step
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     const scalar mass0 = this->mass();
-    const scalar cp0 = this->cp_;
     const scalar np0 = this->nParticle_;
     const scalar T0 = this->T_;
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Initialise transfer terms
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~
+    // 1. Calculate velocity
+    // ~~~~~~~~~~~~~~~~~~~~~
 
-    // Momentum transfer from the particle to the carrier phase
+    scalar Cud = 0.0;
     vector dUTrans = vector::zero;
+    const vector U1 = calcVelocity(td, dt, Cud, dUTrans);
 
-    // Enthalpy transfer from the particle to the carrier phase
-    scalar dhTrans = 0.0;
 
-    // Mass transfer from particle to carrier phase
-    scalarList dMassMT(td.cloud().gases().size(), 0.0);
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 2. Calculate heat transfer
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~
+    scalar htc = 0.0;
+    scalar dhHT = 0.0;
+//    TODO: T1 no longer used - return dhHT instead??????????????????????????????
+//    scalar T1 = calcHeatTransfer(td, dt, cellI, htc, dhHT);
+    calcHeatTransfer(td, dt, cellI, htc, dhHT);
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Calculate velocity - update U
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    scalar Cud = 0.0;
-    const vector U1 = calcVelocity(td, dt, Cud, dUTrans);
 
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 3. Calculate phase change
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~
+    // Mass transfer from particle to carrier phase
+    scalarList dMassPC(td.cloud().gases().size(), 0.0);
+    const scalar dMassPCTot = calcPhaseChange(td, dt, cellI, T0, 1.0, dMassPC);
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Calculate heat transfer - update T
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    scalar htc = 0.0;
-    scalar T1 = calcHeatTransfer(td, dt, cellI, htc, dhTrans);
+    // Update particle component mass fractions
+    updateMassFraction(mass0, dMassPC, Y_);
 
+    // Enthalpy change due to change in particle composition (sink)
+    scalar dhPC = -dMassPCTot*td.cloud().composition().L(0, Y_, pc, T0);
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Calculate phase change - update mass, Y, cp, T, dhTrans
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    const scalar dMassPC = calcPhaseChange(td, dt, cellI, T1, dMassMT);
+    // Enthalpy change due to species released into the carrier (source)
+    scalar HEff = td.cloud().composition().H(0, Y_, pc, T0);
+    dhPC += dMassPCTot*HEff;
 
-    // Update particle mass
-    scalar mass1 = mass0 - dMassPC;
 
-    // Update particle component mass fractions
-    updateMassFraction(mass0, mass1, dMassMT, Y_);
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 4. Collect contributions to determine new particle thermo properties
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // New mass
+    scalar mass1 = mass0 - dMassPCTot;
 
-    // New specific heat capacity of
-    scalar cp1 = td.cloud().composition().cp(0, Y_, pc_, T1);
+    // Total enthalpy transfer from the particle to the carrier phase
+    scalar dhTrans = dhHT + dhPC;
 
-    // Correct temperature due to evaporated components
-    // TODO: use hl function in liquidMixture???
-    T1 -= td.constProps().Lvap()*dMassPC/(0.5*(mass0 + mass1)*cp1);
+    // New specific heat capacity of mixture - using old temperature
+    scalar cp1 = td.cloud().composition().cp(0, Y_, pc_, T0);
 
-    // Correct dhTrans to account for the change in enthalpy due to the
-    // phase change
-    scalar HMix = td.cloud().composition().H(0, Y_, pc, 0.5*(T0 + T1));
-    dhTrans += dMassPC*HMix;
+    // New particle temperature - using average mass over the time interval
+    scalar T1 = T0 + dhTrans/(0.5*(mass0 + mass1)*cp1);
 
 
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 5. Accumulate carrier phase source terms
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     if (td.cloud().coupled())
     {
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        // Accumulate carrier phase source terms
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
         // Transfer mass lost from particle to carrier mass source
-        forAll(dMassMT, i)
+        forAll(dMassPC, i)
         {
-            td.cloud().rhoTrans(i)[cellI] += np0*dMassMT[i];
+            td.cloud().rhoTrans(i)[cellI] += np0*dMassPC[i];
         }
 
         // Update momentum transfer
@@ -154,9 +155,9 @@ void Foam::ReactingParcel<ParcelType>::calc
     }
 
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Remove the particle when mass falls below minimum threshold
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 6a. Remove the particle when mass falls below minimum threshold
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     if (mass1 < td.constProps().minParticleMass())
     {
         td.keepParticle = false;
@@ -164,17 +165,17 @@ void Foam::ReactingParcel<ParcelType>::calc
         if (td.cloud().coupled())
         {
             // Absorb particle(s) into carrier phase
-            forAll(dMassMT, i)
+            forAll(Y_, i)
             {
                 td.cloud().rhoTrans(i)[cellI] += np0*mass1*Y_[i];
             }
             td.cloud().UTrans()[cellI] += np0*mass1*U1;
-            td.cloud().hTrans()[cellI] += np0*mass1*HMix;
+            td.cloud().hTrans()[cellI] += np0*mass1*HEff;
         }
     }
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Set new particle properties
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 6b. Set new particle properties
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     else
     {
         this->U_ = U1;
@@ -202,10 +203,16 @@ Foam::scalar Foam::ReactingParcel<ParcelType>::calcPhaseChange
     const scalar dt,
     const label cellI,
     const scalar T,
+    const scalar YPhase,
     scalarList& dMassMT
 )
 {
-    if (!td.cloud().phaseChange().active())
+    if
+    (
+        !td.cloud().phaseChange().active()
+     || T < td.constProps().Tvap()
+     || YPhase > SMALL
+    )
     {
         return 0.0;
     }
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.H b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.H
index 6848a58d24e..d1c41acc402 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.H
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.H
@@ -89,9 +89,6 @@ public:
             //- Vapourisation temperature [K]
             const scalar Tvap_;
 
-            //- Latent heat of vaporisation [J/kg]
-            const scalar Lvap_;
-
 
     public:
 
@@ -108,9 +105,6 @@ public:
 
             //- Return const access to the vapourisation temperature
             inline scalar Tvap() const;
-
-            //- Return const access to the latent heat of vaporisation
-            inline scalar Lvap() const;
     };
 
 
@@ -194,6 +188,7 @@ protected:
             const scalar dt,
             const label cellI,
             const scalar T,
+            const scalar YPhase,
             scalarList& dMassMT
         );
 
@@ -201,7 +196,6 @@ protected:
         void updateMassFraction
         (
             const scalar mass0,
-            const scalar mass1,
             const scalarList& dMass,
             scalarField& Y
         );
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcelI.H b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcelI.H
index b3d6fd0b4a6..8c0d09e7578 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcelI.H
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcelI.H
@@ -35,8 +35,7 @@ inline Foam::ReactingParcel<ParcelType>::constantProperties::constantProperties
     ThermoParcel<ParcelType>::constantProperties(dict),
     constantVolume_(dict.lookup("constantVolume")),
     Tbp_(dimensionedScalar(dict.lookup("Tbp")).value()),
-    Tvap_(dimensionedScalar(dict.lookup("Tvap")).value()),
-    Lvap_(dimensionedScalar(dict.lookup("Lvap")).value())
+    Tvap_(dimensionedScalar(dict.lookup("Tvap")).value())
 {}
 
 
@@ -131,14 +130,6 @@ Foam::ReactingParcel<ParcelType>::constantProperties::Tvap() const
 }
 
 
-template<class ParcelType>
-inline Foam::scalar
-Foam::ReactingParcel<ParcelType>::constantProperties::Lvap() const
-{
-    return Lvap_;
-}
-
-
 // * * * * * * * * * * * trackData Member Functions  * * * * * * * * * * * * //
 
 template<class ParcelType>
@@ -175,8 +166,7 @@ inline Foam::scalar Foam::ReactingParcel<ParcelType>::mass0() const
 
 
 template<class ParcelType>
-inline const Foam::scalarField& Foam::ReactingParcel<ParcelType>::
-Y() const
+inline const Foam::scalarField& Foam::ReactingParcel<ParcelType>::Y() const
 {
     return Y_;
 }
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcelIO.C b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcelIO.C
index 585f688eef2..eacf49b1126 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcelIO.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcelIO.C
@@ -68,12 +68,12 @@ Foam::ReactingParcel<ParcelType>::ReactingParcel
     // Check state of Istream
     is.check
     (
-        "ReactingParcel<ParcelType>::ReactingParcel\n"
-        "(\n"
-        "    const Cloud<ParcelType>&,\n"
-        "    Istream&,\n"
-        "    bool\n"
-        ")\n"
+        "ReactingParcel<ParcelType>::ReactingParcel"
+        "("
+            "const Cloud<ParcelType>&, "
+            "Istream&, "
+            "bool"
+        ")"
     );
 }
 
@@ -102,8 +102,14 @@ void Foam::ReactingParcel<ParcelType>::readFields
     }
 
     // Get names and sizes for each Y...
-    const wordList phaseTypes = c.composition().phaseTypes();
+    const wordList& phaseTypes = c.composition().phaseTypes();
     const label nPhases = phaseTypes.size();
+    wordList stateLabels(nPhases, "");
+    if (c.composition().nPhase() == 1)
+    {
+        stateLabels = c.composition().stateLabels();
+    }
+
 
     // Set storage for each Y... for each parcel
     forAllIter(typename Cloud<ParcelType>, c, iter)
@@ -117,7 +123,11 @@ void Foam::ReactingParcel<ParcelType>::readFields
     {
         IOField<scalar> Y
         (
-            c.fieldIOobject("Y" + phaseTypes[j], IOobject::MUST_READ)
+            c.fieldIOobject
+            (
+                "Y" + phaseTypes[j] + stateLabels[j],
+                 IOobject::MUST_READ
+            )
         );
 
         label i = 0;
@@ -138,28 +148,37 @@ void Foam::ReactingParcel<ParcelType>::writeFields
 {
     ThermoParcel<ParcelType>::writeFields(c);
 
-    label np =  c.size();
+    const label np =  c.size();
 
-    IOField<scalar> mass0(c.fieldIOobject("mass0", IOobject::NO_READ), np);
-
-    label i = 0;
-    forAllConstIter(typename Cloud<ParcelType>, c, iter)
+    if (np > 0)
     {
-        const ReactingParcel<ParcelType>& p = iter();
-        mass0[i++] = p.mass0_;
-    }
+        IOField<scalar> mass0(c.fieldIOobject("mass0", IOobject::NO_READ), np);
 
-    mass0.write();
+        label i = 0;
+        forAllConstIter(typename Cloud<ParcelType>, c, iter)
+        {
+            const ReactingParcel<ParcelType>& p = iter();
+            mass0[i++] = p.mass0_;
+        }
+        mass0.write();
+
+        // Write the composition fractions
+        const wordList& phaseTypes = c.composition().phaseTypes();
+        wordList stateLabels(phaseTypes.size(), "");
+        if (c.composition().nPhase() == 1)
+        {
+            stateLabels = c.composition().stateLabels();
+        }
 
-    // Write the composition fractions
-    if (np > 0)
-    {
-        const wordList phaseTypes = c.composition().phaseTypes();
         forAll(phaseTypes, j)
         {
             IOField<scalar> Y
             (
-                c.fieldIOobject("Y" + phaseTypes[j], IOobject::NO_READ),
+                c.fieldIOobject
+                (
+                    "Y" + phaseTypes[j] + stateLabels[j],
+                    IOobject::NO_READ
+                ),
                 np
             );
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
index 99bf172920f..a158fc8e52b 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
@@ -53,17 +53,17 @@ void Foam::ThermoParcel<ParcelType>::calc
     const label cellI
 )
 {
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Define local properties at beginning of timestep
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // Define local properties at beginning of time step
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     const vector U0 = this->U_;
     const scalar mass0 = this->mass();
     const scalar np0 = this->nParticle_;
 
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Initialise transfer terms
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 1. Initialise transfer terms
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
     // Momentum transfer from the particle to the carrier phase
     vector dUTrans = vector::zero;
@@ -72,26 +72,25 @@ void Foam::ThermoParcel<ParcelType>::calc
     scalar dhTrans = 0.0;
 
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Calculate velocity - update U
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 2. Calculate velocity - update U
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     scalar Cud = 0.0;
     const vector U1 = calcVelocity(td, dt, Cud, dUTrans);
 
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Calculate heat transfer - update T
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 3. Calculate heat transfer - update T
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     scalar htc = 0.0;
     const scalar T1 = calcHeatTransfer(td, dt, cellI, htc, dhTrans);
 
 
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 4. Accumulate carrier phase source terms
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     if (td.cloud().coupled())
     {
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        // Accumulate carrier phase source terms
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
         // Update momentum transfer
         td.cloud().UTrans()[cellI] += np0*dUTrans;
 
@@ -105,9 +104,9 @@ void Foam::ThermoParcel<ParcelType>::calc
         td.cloud().hCoeff()[cellI] += np0*htc*this->areaS();
     }
 
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Set new particle properties
-    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // 5. Set new particle properties
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     this->U() = U1;
     this->T() = T1;
 }
@@ -164,10 +163,6 @@ Foam::scalar Foam::ThermoParcel<ParcelType>::calcHeatTransfer
     bp *= 6.0/(this->rho_*this->d_*cp_);
 
 
-    //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    // Set new particle temperature
-    //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
     // Integrate to find the new parcel temperature
     IntegrationScheme<scalar>::integrationResult Tres =
         td.cloud().TIntegrator().integrate(T_, dt, ap, bp);
diff --git a/src/lagrangian/intermediate/phaseProperties/phaseProperties/phaseProperties.C b/src/lagrangian/intermediate/phaseProperties/phaseProperties/phaseProperties.C
index afec9962107..9df83bbbaf4 100644
--- a/src/lagrangian/intermediate/phaseProperties/phaseProperties/phaseProperties.C
+++ b/src/lagrangian/intermediate/phaseProperties/phaseProperties/phaseProperties.C
@@ -72,9 +72,9 @@ void Foam::phaseProperties::setGlobalGasIds
 
             FatalErrorIn
             (
-                "void phaseProperties::setGlobalGasIds\n"
-                "(\n"
-                "    const PtrList<volScalarField>&\n"
+                "void phaseProperties::setGlobalGasIds"
+                "("
+                    "const PtrList<volScalarField>&"
                 ")"
             )   << "Could not find gas species " << names_[i]
                 << " in species list" <<  nl
@@ -131,11 +131,47 @@ void Foam::phaseProperties::checkTotalMassFraction() const
 }
 
 
+Foam::word Foam::phaseProperties::phaseToStateLabel(phaseType pt)
+{
+    word state = "(unknown)";
+    switch (pt)
+    {
+        case GAS:
+        {
+            state = "(g)";
+            break;
+        }
+        case LIQUID:
+        {
+            state = "(l)";
+            break;
+        }
+        case SOLID:
+        {
+            state = "(s)";
+            break;
+        }
+        default:
+        {
+            FatalErrorIn
+            (
+                "Foam::phaseProperties::phaseToStateLabel(phaseType pt)"
+            )   << "Invalid phase: " << phaseTypeNames_[pt] << nl
+                << "    phase must be gas, liquid or solid" << nl
+                << exit(FatalError);
+        }
+    }
+
+    return state;
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::phaseProperties::phaseProperties()
 :
     phase_(UNKNOWN),
+    stateLabel_("(unknown)"),
     names_(0),
     Y_(0),
     globalIds_(0)
@@ -145,6 +181,7 @@ Foam::phaseProperties::phaseProperties()
 Foam::phaseProperties::phaseProperties(const phaseProperties& pp)
 :
     phase_(pp.phase_),
+    stateLabel_(pp.stateLabel_),
     names_(pp.names_),
     Y_(pp.Y_),
     globalIds_(pp.globalIds_)
@@ -189,11 +226,11 @@ void Foam::phaseProperties::initialiseGlobalIds
         {
             FatalErrorIn
             (
-                "Foam::phaseProperties::setGlobalIds\n"
-                "(\n"
-                "    const PtrList<volScalarField>&,\n"
-                "    const wordList&,\n"
-                "    const wordList&\n"
+                "Foam::phaseProperties::setGlobalIds"
+                "("
+                    "const PtrList<volScalarField>&, "
+                    "const wordList&, "
+                    "const wordList&"
                 ")"
             )   << "Invalid phase: " << phaseTypeNames_[phase_] << nl
                 << "    phase must be gas, liquid or solid" << nl
@@ -209,6 +246,12 @@ Foam::phaseProperties::phaseType Foam::phaseProperties::phase() const
 }
 
 
+const Foam::word& Foam::phaseProperties::stateLabel() const
+{
+    return stateLabel_;
+}
+
+
 Foam::word Foam::phaseProperties::phaseTypeName() const
 {
     return phaseTypeNames_[phase_];
@@ -229,7 +272,7 @@ const Foam::word& Foam::phaseProperties::name(const label cmptI) const
         (
             "const Foam::word& Foam::phaseProperties::name"
             "("
-                "const label cmptI"
+                "const label"
             ") const"
         )   << "Requested component " << cmptI << "out of range" << nl
             << "Available phase components:" << nl << names_ << nl
@@ -254,7 +297,7 @@ Foam::scalar& Foam::phaseProperties::Y(const label cmptI)
         (
             "const Foam::scalar& Foam::phaseProperties::Y"
             "("
-                "const label cmptI"
+                "const label"
             ") const"
         )   << "Requested component " << cmptI << "out of range" << nl
             << "Available phase components:" << nl << names_ << nl
diff --git a/src/lagrangian/intermediate/phaseProperties/phaseProperties/phaseProperties.H b/src/lagrangian/intermediate/phaseProperties/phaseProperties/phaseProperties.H
index 258d0069aed..2998a9831a5 100644
--- a/src/lagrangian/intermediate/phaseProperties/phaseProperties/phaseProperties.H
+++ b/src/lagrangian/intermediate/phaseProperties/phaseProperties/phaseProperties.H
@@ -76,6 +76,9 @@ private:
         //- Phase type
         phaseType phase_;
 
+        //- State label (s), (l), (g) etc.
+        word stateLabel_;
+
         //- List of component names
         List<word> names_;
 
@@ -97,6 +100,9 @@ private:
         //- Check the total mass fraction
         void checkTotalMassFraction() const;
 
+        //- Set the state label
+        word phaseToStateLabel(phaseType pt);
+
 
 public:
 
@@ -132,6 +138,9 @@ public:
             //- Return const access to the phase type
             phaseType phase() const;
 
+            //- Return const access to the phase state label
+            const word& stateLabel() const;
+
             //- Return word representation of the phase type
             word phaseTypeName() const;
 
diff --git a/src/lagrangian/intermediate/phaseProperties/phaseProperties/phasePropertiesIO.C b/src/lagrangian/intermediate/phaseProperties/phaseProperties/phasePropertiesIO.C
index 192136ed8d0..b0c3937ea4b 100644
--- a/src/lagrangian/intermediate/phaseProperties/phaseProperties/phasePropertiesIO.C
+++ b/src/lagrangian/intermediate/phaseProperties/phaseProperties/phasePropertiesIO.C
@@ -32,6 +32,7 @@ License
 Foam::phaseProperties::phaseProperties(Istream& is)
 :
     phase_(UNKNOWN),
+    stateLabel_("(unknown)"),
     names_(0),
     Y_(0),
     globalIds_(0)
@@ -53,6 +54,8 @@ Foam::phaseProperties::phaseProperties(Istream& is)
 
     phase_ = phaseTypeNames_[phaseInfo.keyword()];
 
+    stateLabel_ = phaseToStateLabel(phase_);
+
     label cmptI = 0;
     forAllConstIter(IDLList<entry>, phaseInfo, iter)
     {
@@ -95,6 +98,8 @@ Foam::Istream& Foam::operator>>(Istream& is, phaseProperties& pp)
 
     pp.phase_ = pp.phaseTypeNames_[phaseInfo.keyword()];
 
+    pp.stateLabel_ = pp.phaseToStateLabel(pp.phase_);
+
     label cmptI = 0;
     forAllConstIter(IDLList<entry>, phaseInfo, iter)
     {
@@ -117,7 +122,10 @@ Foam::Ostream& Foam::operator<<(Ostream& os, const phaseProperties& pp)
     os.check
     (
         "Foam::Ostream& Foam::operator<<"
-        "(Foam::Ostream&, const Foam::phaseProperties&)"
+        "("
+            "Foam::Ostream&, "
+            "const Foam::phaseProperties&"
+        ")"
     );
 
     os  << pp.phaseTypeNames_[pp.phase_] << nl << token::BEGIN_BLOCK << nl
@@ -134,7 +142,10 @@ Foam::Ostream& Foam::operator<<(Ostream& os, const phaseProperties& pp)
     os.check
     (
         "Foam::Ostream& Foam::operator<<"
-        "(Foam::Ostream&, const Foam::phaseProperties&)"
+        "("
+            "Foam::Ostream&, "
+            "const Foam::phaseProperties&"
+        ")"
     );
 
     return os;
diff --git a/src/lagrangian/intermediate/phaseProperties/phasePropertiesList/phasePropertiesList.C b/src/lagrangian/intermediate/phaseProperties/phasePropertiesList/phasePropertiesList.C
index 9e5aeb69c3b..f1146e0245c 100644
--- a/src/lagrangian/intermediate/phaseProperties/phasePropertiesList/phasePropertiesList.C
+++ b/src/lagrangian/intermediate/phaseProperties/phasePropertiesList/phasePropertiesList.C
@@ -37,7 +37,8 @@ Foam::phasePropertiesList::phasePropertiesList
 )
 :
     props_(is),
-    phaseTypeNames_()
+    phaseTypeNames_(),
+    stateLabels_()
 {
     forAll(props_, i)
     {
@@ -45,9 +46,11 @@ Foam::phasePropertiesList::phasePropertiesList
     }
 
     phaseTypeNames_.setSize(props_.size());
-    forAll(phaseTypeNames_, i)
+    stateLabels_.setSize(props_.size());
+    forAll(props_, i)
     {
         phaseTypeNames_[i] = props_[i].phaseTypeName();
+        stateLabels_[i] = props_[i].stateLabel();
     }
 }
 
@@ -73,6 +76,12 @@ const Foam::wordList& Foam::phasePropertiesList::phaseTypes() const
 }
 
 
+const Foam::wordList& Foam::phasePropertiesList::stateLabels() const
+{
+    return stateLabels_;
+}
+
+
 Foam::label Foam::phasePropertiesList::size() const
 {
     return props_.size();
diff --git a/src/lagrangian/intermediate/phaseProperties/phasePropertiesList/phasePropertiesList.H b/src/lagrangian/intermediate/phaseProperties/phasePropertiesList/phasePropertiesList.H
index 80e8f5c859d..4247f2c4bad 100644
--- a/src/lagrangian/intermediate/phaseProperties/phasePropertiesList/phasePropertiesList.H
+++ b/src/lagrangian/intermediate/phaseProperties/phasePropertiesList/phasePropertiesList.H
@@ -60,6 +60,9 @@ class phasePropertiesList
         //- List of word representation of phase types
         wordList phaseTypeNames_;
 
+        //- List of state labels
+        wordList stateLabels_;
+
 
 public:
 
@@ -84,6 +87,9 @@ public:
         //- Return the list of word representation of phase types
         const wordList& phaseTypes() const;
 
+        //- Return the list of state labels
+        const wordList& stateLabels() const;
+
         //- Return the size (number of phases)
         label size() const;
 
diff --git a/src/lagrangian/intermediate/submodels/IO/DataEntry/DataEntry/NewDataEntry.C b/src/lagrangian/intermediate/submodels/IO/DataEntry/DataEntry/NewDataEntry.C
index 425dcfe6f77..3bed47c540e 100644
--- a/src/lagrangian/intermediate/submodels/IO/DataEntry/DataEntry/NewDataEntry.C
+++ b/src/lagrangian/intermediate/submodels/IO/DataEntry/DataEntry/NewDataEntry.C
@@ -47,7 +47,7 @@ Foam::autoPtr<Foam::DataEntry<Type> > Foam::DataEntry<Type>::New
         FatalErrorIn("DataEntry<Type>::New(Istream&)")
             << "Unknown DataEntry type " << DataEntryType << " for DataEntry "
             << entryName << ". Constructor not in hash table" << nl << nl
-            << "    Valid DataEntry types are :" << nl
+            << "    Valid DataEntry types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << nl
             << exit(FatalError);
     }
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/DispersionModel/DispersionModel/NewDispersionModel.C b/src/lagrangian/intermediate/submodels/Kinematic/DispersionModel/DispersionModel/NewDispersionModel.C
index 4db4d791943..09717b8e864 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/DispersionModel/DispersionModel/NewDispersionModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/DispersionModel/DispersionModel/NewDispersionModel.C
@@ -38,10 +38,7 @@ Foam::DispersionModel<CloudType>::New
     CloudType& owner
 )
 {
-    word DispersionModelType
-    (
-        dict.lookup("DispersionModel")
-    );
+    word DispersionModelType(dict.lookup("DispersionModel"));
 
     Info<< "Selecting DispersionModel " << DispersionModelType << endl;
 
@@ -53,11 +50,14 @@ Foam::DispersionModel<CloudType>::New
         FatalErrorIn
         (
             "DispersionModel<CloudType>::New"
-            "(const dictionary&, CloudType&)"
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
         )   << "Unknown DispersionModelType type "
             << DispersionModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid DispersionModel types are :" << nl
+            << "    Valid DispersionModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << exit(FatalError);
     }
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/DragModel/DragModel/NewDragModel.C b/src/lagrangian/intermediate/submodels/Kinematic/DragModel/DragModel/NewDragModel.C
index 3cee2ec548d..6e1e80d3cda 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/DragModel/DragModel/NewDragModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/DragModel/DragModel/NewDragModel.C
@@ -35,10 +35,7 @@ Foam::autoPtr<Foam::DragModel<CloudType> > Foam::DragModel<CloudType>::New
     CloudType& owner
 )
 {
-    word DragModelType
-    (
-        dict.lookup("DragModel")
-    );
+    word DragModelType(dict.lookup("DragModel"));
 
     Info<< "Selecting DragModel " << DragModelType << endl;
 
@@ -50,12 +47,14 @@ Foam::autoPtr<Foam::DragModel<CloudType> > Foam::DragModel<CloudType>::New
         FatalErrorIn
         (
             "DragModel<CloudType>::New"
-            "(const dictionary&, CloudType&)"
-        )
-            << "Unknown DragModelType type "
+            "("
+                "const dictionary&,"
+                "CloudType&"
+            ")"
+        )   << "Unknown DragModelType type "
             << DragModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid DragModel types are :" << nl
+            << "    Valid DragModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc()
             << exit(FatalError);
     }
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
index a21580845e4..adf2559e985 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/InjectionModel.C
@@ -170,10 +170,10 @@ void Foam::InjectionModel<CloudType>::findCellAtPosition
     {
         FatalErrorIn
         (
-            "Foam::InjectionModel<CloudType>::findCellAtPosition\n"
-            "(\n"
-            "    label&,\n"
-            "    vector&\n"
+            "Foam::InjectionModel<CloudType>::findCellAtPosition"
+            "("
+                "label&, "
+                "vector&"
             ")"
         )<< "Cannot find parcel injection cell. "
          << "Parcel position = " << p0 << nl
@@ -212,13 +212,13 @@ Foam::scalar Foam::InjectionModel<CloudType>::setNumberOfParticles
             FatalErrorIn
             (
                 "Foam::scalar "
-                "Foam::InjectionModel<CloudType>::setNumberOfParticles\n"
-                "(\n"
-                "    const label,\n"
-                "    const scalar,\n"
-                "    const scalar,\n"
-                "    const scalar,\n"
-                "    const scalar\n"
+                "Foam::InjectionModel<CloudType>::setNumberOfParticles"
+                "("
+                "    const label, "
+                "    const scalar, "
+                "    const scalar, "
+                "    const scalar, "
+                "    const scalar"
                 ")"
             )<< "Unknown parcelBasis type" << nl
              << exit(FatalError);
@@ -314,11 +314,11 @@ Foam::InjectionModel<CloudType>::InjectionModel
     {
         FatalErrorIn
         (
-            "Foam::InjectionModel<CloudType>::InjectionModel\n"
-            "(\n"
-            "    const dictionary&,\n"
-            "    CloudType&,\n"
-            "    const word&\n"
+            "Foam::InjectionModel<CloudType>::InjectionModel"
+            "("
+                "const dictionary&, "
+                "CloudType&, "
+                "const word&"
             ")"
         )<< "parcelBasisType must be either 'number' or 'mass'" << nl
          << exit(FatalError);
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/NewInjectionModel.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/NewInjectionModel.C
index 3501e367aad..c0926272c62 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/NewInjectionModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/InjectionModel/NewInjectionModel.C
@@ -36,10 +36,7 @@ Foam::InjectionModel<CloudType>::New
     CloudType& owner
 )
 {
-    word InjectionModelType
-    (
-        dict.lookup("InjectionModel")
-    );
+    word InjectionModelType(dict.lookup("InjectionModel"));
 
     Info<< "Selecting InjectionModel " << InjectionModelType << endl;
 
@@ -51,11 +48,14 @@ Foam::InjectionModel<CloudType>::New
         FatalErrorIn
         (
             "InjectionModel<CloudType>::New"
-            "(const dictionary&, CloudType&)"
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
         )   << "Unknown InjectionModelType type "
             << InjectionModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid InjectionModel types are :" << nl
+            << "    Valid InjectionModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << exit(FatalError);
     }
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/WallInteractionModel/WallInteractionModel/NewWallInteractionModel.C b/src/lagrangian/intermediate/submodels/Kinematic/WallInteractionModel/WallInteractionModel/NewWallInteractionModel.C
index ee1fd8a169e..522432191af 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/WallInteractionModel/WallInteractionModel/NewWallInteractionModel.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/WallInteractionModel/WallInteractionModel/NewWallInteractionModel.C
@@ -36,10 +36,7 @@ Foam::WallInteractionModel<CloudType>::New
     CloudType& owner
 )
 {
-    word WallInteractionModelType
-    (
-        dict.lookup("WallInteractionModel")
-    );
+    word WallInteractionModelType(dict.lookup("WallInteractionModel"));
 
     Info<< "Selecting WallInteractionModel " << WallInteractionModelType
         << endl;
@@ -52,12 +49,14 @@ Foam::WallInteractionModel<CloudType>::New
         FatalErrorIn
         (
             "WallInteractionModel<CloudType>::New"
-            "(const dictionary&, CloudType&)"
-        )
-            << "Unknown WallInteractionModelType type "
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
+        )   << "Unknown WallInteractionModelType type "
             << WallInteractionModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid WallInteractionModel types are :" << nl
+            << "    Valid WallInteractionModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << exit(FatalError);
     }
 
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
index df898464fda..5408ad656e8 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
@@ -161,6 +161,14 @@ Foam::CompositionModel<CloudType>::phaseTypes() const
 }
 
 
+template<class CloudType>
+const Foam::wordList&
+Foam::CompositionModel<CloudType>::stateLabels() const
+{
+    return phaseProps_.stateLabels();
+}
+
+
 template<class CloudType>
 const Foam::wordList&
 Foam::CompositionModel<CloudType>::componentNames(const label phaseI) const
@@ -182,10 +190,10 @@ Foam::label Foam::CompositionModel<CloudType>::globalId
     {
         FatalErrorIn
         (
-            "Foam::label Foam::CompositionModel<CloudType>::globalId\n"
-            "(\n"
-            "    const label phaseI,\n"
-            "    const word& cmptName\n"
+            "Foam::label Foam::CompositionModel<CloudType>::globalId"
+            "("
+                "const label, "
+                "const word&"
             ") const"
         )   << "Unable to determine global id for requested component "
             << cmptName << nl << abort(FatalError);
@@ -218,10 +226,10 @@ Foam::label Foam::CompositionModel<CloudType>::localId
     {
         FatalErrorIn
         (
-            "Foam::label Foam::CompositionModel<CloudType>::localId\n"
-            "(\n"
-            "    const label phaseI,\n"
-            "    const word& cmptName\n"
+            "Foam::label Foam::CompositionModel<CloudType>::localId"
+            "("
+                "const label, "
+                "const word&"
             ") const"
         )   << "Unable to determine local id for component " << cmptName
             << nl << abort(FatalError);
@@ -277,10 +285,10 @@ Foam::scalarField Foam::CompositionModel<CloudType>::X
         {
             FatalErrorIn
             (
-                "Foam::scalarField Foam::CompositionModel<CloudType>::X\n"
-                "(\n"
-                "    const label phaseI,\n"
-                "    const scalarField Y\n"
+                "Foam::scalarField Foam::CompositionModel<CloudType>::X"
+                "("
+                    "const label, "
+                    "const scalarField&"
                 ") const"
             )   << "Only possible to convert gas and liquid mass fractions"
                 << nl << abort(FatalError);
@@ -340,12 +348,12 @@ Foam::scalar Foam::CompositionModel<CloudType>::H
         {
             FatalErrorIn
             (
-                "Foam::scalar Foam::CompositionModel<CloudType>::H\n"
-                "(\n"
-                "    const label phaseI,\n"
-                "    const scalarField& Y,"
-                "    const scalar p,\n"
-                "    const scalar T\n"
+                "Foam::scalar Foam::CompositionModel<CloudType>::H"
+                "("
+                "    const label, "
+                "    const scalarField&, "
+                "    const scalar, "
+                "    const scalar"
                 ") const"
             )   << "Unknown phase enumeration" << nl << abort(FatalError);
         }
@@ -399,12 +407,12 @@ Foam::scalar Foam::CompositionModel<CloudType>::cp
         {
             FatalErrorIn
             (
-                "Foam::scalar Foam::CompositionModel<CloudType>::cp\n"
-                "(\n"
-                "    const label phaseI,\n"
-                "    const scalarField& Y,"
-                "    const scalar p,\n"
-                "    const scalar T\n"
+                "Foam::scalar Foam::CompositionModel<CloudType>::cp"
+                "("
+                    "const label, "
+                    "const scalarField&, "
+                    "const scalar, "
+                    "const scalar"
                 ") const"
             )   << "Unknown phase enumeration" << nl << abort(FatalError);
         }
@@ -414,6 +422,77 @@ Foam::scalar Foam::CompositionModel<CloudType>::cp
 }
 
 
+template<class CloudType>
+Foam::scalar Foam::CompositionModel<CloudType>::L
+(
+    const label phaseI,
+    const scalarField& Y,
+    const scalar p,
+    const scalar T
+) const
+{
+    const phaseProperties& props = phaseProps_[phaseI];
+    scalar LMixture = 0.0;
+    switch (props.phase())
+    {
+        case phaseProperties::GAS:
+        {
+            notImplemented
+            (
+                "Foam::scalar Foam::CompositionModel<CloudType>::L"
+                "("
+                    "const label, "
+                    "const scalarField&, "
+                    "const scalar, "
+                    "const scalar"
+                ") const\n"
+                "*** No support for gaseous components"
+            );
+            break;
+        }
+        case phaseProperties::LIQUID:
+        {
+            forAll(Y, i)
+            {
+                label gid = props.globalIds()[i];
+                LMixture += Y[i]*this->liquids().properties()[gid].hl(p, T);
+            }
+            break;
+        }
+        case phaseProperties::SOLID:
+        {
+            notImplemented
+            (
+                "Foam::scalar Foam::CompositionModel<CloudType>::L"
+                "("
+                    "const label, "
+                    "const scalarField&, "
+                    "const scalar, "
+                    "const scalar"
+                ") const\n"
+                "*** No support for solid components"
+            );
+            break;
+        }
+        default:
+        {
+            FatalErrorIn
+            (
+                "Foam::scalar Foam::CompositionModel<CloudType>::L"
+                "("
+                    "const label, "
+                    "const scalarField&, "
+                    "const scalar, "
+                    "const scalar"
+                ") const"
+            )   << "Unknown phase enumeration" << nl << abort(FatalError);
+        }
+    }
+
+    return LMixture;
+}
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 #include "NewCompositionModel.C"
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
index 5be0e2c1b92..c5382a4f9f8 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
@@ -172,6 +172,9 @@ public:
                 //  If only 1 phase, return the component names of that phase
                 const wordList& phaseTypes() const;
 
+                //- Return the list of state labels (s), (l), (g) etc.
+                const wordList& stateLabels() const;
+
                 //- Return the list of component names for phaseI
                 const wordList& componentNames(const label phaseI) const;
 
@@ -226,7 +229,7 @@ public:
                 const scalar T
             ) const;
 
-            //- Return cp for the phase phaseI
+            //- Return specific heat caoacity for the phase phaseI
             virtual scalar cp
             (
                 const label phaseI,
@@ -234,6 +237,15 @@ public:
                 const scalar p,
                 const scalar T
             ) const;
+
+            //- Return latent heat for the phase phaseI
+            virtual scalar L
+            (
+                const label phaseI,
+                const scalarField& Y,
+                const scalar p,
+                const scalar T
+            ) const;
 };
 
 
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/NewCompositionModel.C b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/NewCompositionModel.C
index 707c158cc6a..4a9f8d7b745 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/NewCompositionModel.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/NewCompositionModel.C
@@ -47,16 +47,15 @@ Foam::CompositionModel<CloudType>::New
     {
         FatalErrorIn
         (
-            "CompositionModel<CloudType>::New\n"
-            "(\n"
-            "    const dictionary&,\n"
-            "    CloudType&\n"
+            "CompositionModel<CloudType>::New"
+            "("
+                "const dictionary&, "
+                "CloudType&"
             ")"
-        )
-            << "Unknown CompositionModelType type "
+        )   << "Unknown CompositionModelType type "
             << CompositionModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid CompositionModel types are :" << nl
+            << "    Valid CompositionModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << nl
             << exit(FatalError);
     }
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/SingleMixtureFraction/SingleMixtureFraction.C b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/SingleMixtureFraction/SingleMixtureFraction.C
index 9f5eed2c7c9..fb00ddc7633 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/SingleMixtureFraction/SingleMixtureFraction.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/SingleMixtureFraction/SingleMixtureFraction.C
@@ -106,11 +106,11 @@ Foam::SingleMixtureFraction<CloudType>::SingleMixtureFraction
         FatalErrorIn
         (
             "Foam::SingleMixtureFraction<CloudType>::"
-            "SingleMixtureFraction\n"
-            "(\n"
-            "    const dictionary&,\n"
-            "    CloudType&\n"
-            ")\n"
+            "SingleMixtureFraction"
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
         )   << "Incorrect numebr of phases: " << nl
             << "    Please specify 1 gas, 1 liquid and 1 solid" << nl
             << exit(FatalError);
@@ -125,11 +125,11 @@ Foam::SingleMixtureFraction<CloudType>::SingleMixtureFraction
         FatalErrorIn
         (
             "Foam::SingleMixtureFraction<CloudType>::"
-            "SingleMixtureFraction\n"
-            "(\n"
-            "    const dictionary&,\n"
-            "    CloudType&\n"
-            ")\n"
+            "SingleMixtureFraction"
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
         )   << "Sum of phases should be 1. Phase fractions:" << nl
             << YMixture0_ << nl << exit(FatalError);
     }
diff --git a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.C b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.C
index e3b9383ecd4..a3fde009a16 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.C
@@ -118,9 +118,9 @@ Foam::LiquidEvaporation<CloudType>::LiquidEvaporation
             )
         )
     ),
-    Tvap_(readScalar(this->coeffDict().lookup("Tvap"))),
     activeLiquids_(this->coeffDict().lookup("activeLiquids")),
-    liqToGasMap_(activeLiquids_.size())
+    liqToGasMap_(activeLiquids_.size(), -1),
+    liqToLiqMap_(activeLiquids_.size(), -1)
 {
     if (activeLiquids_.size() == 0)
     {
@@ -135,11 +135,37 @@ Foam::LiquidEvaporation<CloudType>::LiquidEvaporation
             << nl << endl;
     }
 
-    // Calculate mapping between liquid and carrier phase species
+    // Determine mapping between liquid and carrier phase species
     forAll(activeLiquids_, i)
     {
         liqToGasMap_[i] = carrierSpecieId(activeLiquids_[i]);
     }
+
+    // Determine mapping between local and global liquids
+    forAll(activeLiquids_, i)
+    {
+        forAll(liquids_->components(), j)
+        {
+            if (liquids_->components()[j] == activeLiquids_[i])
+            {
+                liqToLiqMap_[i] = j;
+            }
+        }
+
+        if (liqToLiqMap_[i] < 0)
+        {
+            FatalErrorIn
+            (
+                "Foam::LiquidEvaporation<CloudType>::LiquidEvaporation"
+                "("
+                    "const dictionary& dict, "
+                    "CloudType& owner"
+                ")"
+            )   << "Unable to find liquid species " << activeLiquids_[i]
+                << " in global liquids list. Avaliable liquids:" << nl
+                << liquids_->components() << nl << exit(FatalError);
+        }
+    }
 }
 
 
@@ -173,58 +199,52 @@ Foam::scalar Foam::LiquidEvaporation<CloudType>::calculate
     scalarList& dMassMT
 ) const
 {
+    // initialise total mass transferred from the particle to carrier phase
     scalar dMassTot = 0.0;
 
-    if (T < Tvap_)
-    {
-        // not reached model activation temperature
-        return dMassTot;
-    }
-    else
-    {
-        // Construct carrier phase species volume fractions
-        scalarField Xc = calcXc(cellI);
+    // construct carrier phase species volume fractions for cell, cellI
+    scalarField Xc = calcXc(cellI);
 
-        // droplet surface area
-        scalar A = mathematicalConstant::pi*sqr(d);
+    // droplet surface area
+    scalar A = mathematicalConstant::pi*sqr(d);
 
-        // Reynolds number
-        scalar Re = mag(Ur)*d/(nuc + ROOTVSMALL);
+    // Reynolds number
+    scalar Re = mag(Ur)*d/(nuc + ROOTVSMALL);
 
-        // Calculate mass transfer of each specie in liquid
-        forAll(activeLiquids_, i)
-        {
-            label gid = liqToGasMap_[i];
+    // calculate mass transfer of each specie in liquid
+    forAll(activeLiquids_, i)
+    {
+        label gid = liqToGasMap_[i];
+        label lid = liqToLiqMap_[i];
 
-            // Vapour diffusivity [m2/s]
-            scalar Dab = liquids_->properties()[gid].D(pc, T);
+        // vapour diffusivity [m2/s]
+        scalar Dab = liquids_->properties()[lid].D(pc, T);
 
-            // Saturation pressure for species i [pa]
-            scalar pSat = liquids_->properties()[gid].pv(pc, T);
+        // saturation pressure for species i [pa]
+        scalar pSat = liquids_->properties()[lid].pv(pc, T);
 
-            // Schmidt number
-            scalar Sc = nuc/(Dab + ROOTVSMALL);
+        // Schmidt number
+        scalar Sc = nuc/(Dab + ROOTVSMALL);
 
-            // Sherwood number
-            scalar Sh = this->Sh(Re, Sc);
+        // Sherwood number
+        scalar Sh = this->Sh(Re, Sc);
 
-            // mass transfer coefficient [m/s]
-            scalar kc = Sh*Dab/(d + ROOTVSMALL);
+        // mass transfer coefficient [m/s]
+        scalar kc = Sh*Dab/(d + ROOTVSMALL);
 
-            // vapour concentration at droplet surface [kgmol/m3]
-            scalar Cs = pSat/(specie::RR*T);
+        // vapour concentration at droplet surface [kgmol/m3]
+        scalar Cs = pSat/(specie::RR*T);
 
-            // vapour concentration in bulk gas [kgmol/m3]
-            scalar Cinf = Xc[gid]*pc/(specie::RR*Tc);
+        // vapour concentration in bulk gas [kgmol/m3]
+        scalar Cinf = Xc[gid]*pc/(specie::RR*Tc);
 
-            // molar flux of vapour [kgmol/m2/s]
-            scalar Ni = max(kc*(Cs - Cinf), 0.0);
+        // molar flux of vapour [kgmol/m2/s]
+        scalar Ni = max(kc*(Cs - Cinf), 0.0);
 
-            // mass transfer [kg]
-            scalar dm = Ni*A*liquids_->properties()[gid].W()*dt;
-            dMassMT[gid] -= dm;
-            dMassTot += dm;
-        }
+        // mass transfer [kg]
+        scalar dm = Ni*A*liquids_->properties()[lid].W()*dt;
+        dMassMT[gid] += dm;
+        dMassTot += dm;
     }
 
     return dMassTot;
diff --git a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.H b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.H
index 85d3714df34..a9321af4a91 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.H
@@ -56,17 +56,15 @@ protected:
         //- Global liquid properties data
         autoPtr<liquidMixture> liquids_;
 
-        //- Vaporisation temperature [K]
-        //  Local droplet temperature must exceed Tvap before evaporation is
-        //  allowed
-        scalar Tvap_;
-
         //- List of active liquid names
         List<word> activeLiquids_;
 
         //- Mapping between liquid and carrier gas species
         List<label> liqToGasMap_;
 
+        //- Mapping between local and global liquid species
+        List<label> liqToLiqMap_;
+
 
     // Protected member functions
 
diff --git a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/PhaseChangeModel/NewPhaseChangeModel.C b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/PhaseChangeModel/NewPhaseChangeModel.C
index 091fd8b810b..efb611522a8 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/PhaseChangeModel/NewPhaseChangeModel.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/PhaseChangeModel/NewPhaseChangeModel.C
@@ -47,16 +47,15 @@ Foam::PhaseChangeModel<CloudType>::New
     {
         FatalErrorIn
         (
-            "PhaseChangeModel<CloudType>::New\n"
-            "(\n"
-            "    const dictionary&,\n"
-            "    CloudType&\n"
+            "PhaseChangeModel<CloudType>::New"
+            "("
+                "const dictionary&, "
+                "CloudType&"
             ")"
-        )
-            << "Unknown PhaseChangeModelType type "
+        )   << "Unknown PhaseChangeModelType type "
             << PhaseChangeModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid PhaseChangeModel types are :" << nl
+            << "    Valid PhaseChangeModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << exit(FatalError);
     }
 
diff --git a/src/lagrangian/intermediate/submodels/ReactingMultiphase/DevolatilisationModel/DevolatilisationModel/NewDevolatilisationModel.C b/src/lagrangian/intermediate/submodels/ReactingMultiphase/DevolatilisationModel/DevolatilisationModel/NewDevolatilisationModel.C
index f51c4568c01..838bdf75e21 100644
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/DevolatilisationModel/DevolatilisationModel/NewDevolatilisationModel.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/DevolatilisationModel/DevolatilisationModel/NewDevolatilisationModel.C
@@ -36,10 +36,7 @@ Foam::DevolatilisationModel<CloudType>::New
     CloudType& owner
 )
 {
-    word DevolatilisationModelType
-    (
-        dict.lookup("DevolatilisationModel")
-    );
+    word DevolatilisationModelType(dict.lookup("DevolatilisationModel"));
 
     Info<< "Selecting DevolatilisationModel " << DevolatilisationModelType
         << endl;
@@ -52,12 +49,14 @@ Foam::DevolatilisationModel<CloudType>::New
         FatalErrorIn
         (
             "DevolatilisationModel<CloudType>::New"
-            "(const dictionary&, CloudType&)"
-        )
-            << "Unknown DevolatilisationModelType type "
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
+        )   << "Unknown DevolatilisationModelType type "
             << DevolatilisationModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid DevolatilisationModel types are :" << nl
+            << "    Valid DevolatilisationModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << exit(FatalError);
     }
 
diff --git a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/NoSurfaceReaction/NoSurfaceReaction.C b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/NoSurfaceReaction/NoSurfaceReaction.C
index 3923455a921..8556282b9df 100644
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/NoSurfaceReaction/NoSurfaceReaction.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/NoSurfaceReaction/NoSurfaceReaction.C
@@ -31,11 +31,11 @@ License
 template <class CloudType>
 Foam::NoSurfaceReaction<CloudType>::NoSurfaceReaction
 (
-    const dictionary& dict,
+    const dictionary&,
     CloudType& owner
 )
 :
-    SurfaceReactionModel<CloudType>(dict, owner, typeName)
+    SurfaceReactionModel<CloudType>(owner)
 {}
 
 
diff --git a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/SurfaceReactionModel/NewSurfaceReactionModel.C b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/SurfaceReactionModel/NewSurfaceReactionModel.C
index 2303e58f21c..348bfa8a8c7 100644
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/SurfaceReactionModel/NewSurfaceReactionModel.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/SurfaceReactionModel/NewSurfaceReactionModel.C
@@ -36,10 +36,7 @@ Foam::SurfaceReactionModel<CloudType>::New
     CloudType& owner
 )
 {
-    word SurfaceReactionModelType
-    (
-        dict.lookup("SurfaceReactionModel")
-    );
+    word SurfaceReactionModelType(dict.lookup("SurfaceReactionModel"));
 
     Info<< "Selecting SurfaceReactionModel " << SurfaceReactionModelType
         << endl;
@@ -52,12 +49,14 @@ Foam::SurfaceReactionModel<CloudType>::New
         FatalErrorIn
         (
             "SurfaceReactionModel<CloudType>::New"
-            "(const dictionary&, CloudType&)"
-        )
-            << "Unknown SurfaceReactionModelType type "
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
+        )   << "Unknown SurfaceReactionModelType type "
             << SurfaceReactionModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid SurfaceReactionModel types are :" << nl
+            << "    Valid SurfaceReactionModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << exit(FatalError);
     }
 
diff --git a/src/lagrangian/intermediate/submodels/Thermodynamic/HeatTransferModel/HeatTransferModel/NewHeatTransferModel.C b/src/lagrangian/intermediate/submodels/Thermodynamic/HeatTransferModel/HeatTransferModel/NewHeatTransferModel.C
index fb0362bbc3a..d5d4748c4db 100644
--- a/src/lagrangian/intermediate/submodels/Thermodynamic/HeatTransferModel/HeatTransferModel/NewHeatTransferModel.C
+++ b/src/lagrangian/intermediate/submodels/Thermodynamic/HeatTransferModel/HeatTransferModel/NewHeatTransferModel.C
@@ -36,10 +36,7 @@ Foam::HeatTransferModel<CloudType>::New
     CloudType& owner
 )
 {
-    word HeatTransferModelType
-    (
-        dict.lookup("HeatTransferModel")
-    );
+    word HeatTransferModelType(dict.lookup("HeatTransferModel"));
 
     Info<< "Selecting HeatTransferModel " << HeatTransferModelType << endl;
 
@@ -51,12 +48,14 @@ Foam::HeatTransferModel<CloudType>::New
         FatalErrorIn
         (
             "HeatTransferModel<CloudType>::New"
-            "(const dictionary&, CloudType&)"
-        )
-            << "Unknown HeatTransferModelType type "
+            "("
+                "const dictionary&, "
+                "CloudType&"
+            ")"
+        )   << "Unknown HeatTransferModelType type "
             << HeatTransferModelType
             << ", constructor not in hash table" << nl << nl
-            << "    Valid HeatTransferModel types are :" << nl
+            << "    Valid HeatTransferModel types are:" << nl
             << dictionaryConstructorTablePtr_->toc() << exit(FatalError);
     }
 
-- 
GitLab