diff --git a/src/lagrangian/coalCombustion/CoalCloud/CoalCloud.C b/src/lagrangian/coalCombustion/CoalCloud/CoalCloud.C
index 873b82d170ba37b1d33beb485a46533b3266a9f2..a4d249d265fe3f5fff3b62957e0a69eec4bd910c 100644
--- a/src/lagrangian/coalCombustion/CoalCloud/CoalCloud.C
+++ b/src/lagrangian/coalCombustion/CoalCloud/CoalCloud.C
@@ -35,7 +35,7 @@ Foam::CoalCloud<ThermoType>::CoalCloud
     const volScalarField& rho,
     const volVectorField& U,
     const dimensionedVector& g,
-    hCombustionThermo& thermo
+    basicThermo& thermo
 )
 :
     ReactingMultiphaseCloud<CoalParcel<ThermoType> >
diff --git a/src/lagrangian/coalCombustion/CoalCloud/CoalCloud.H b/src/lagrangian/coalCombustion/CoalCloud/CoalCloud.H
index 91a0289f0ffe3439d1f0900f7da70d073e2c96e1..36e3b4e03ee7e5d4e4875604c254c8fbf1b22416 100644
--- a/src/lagrangian/coalCombustion/CoalCloud/CoalCloud.H
+++ b/src/lagrangian/coalCombustion/CoalCloud/CoalCloud.H
@@ -76,7 +76,7 @@ public:
             const volScalarField& rho,
             const volVectorField& U,
             const dimensionedVector& g,
-            hCombustionThermo& thermo
+            basicThermo& thermo
         );
 
 
diff --git a/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationDiffusionLimitedRate/COxidationDiffusionLimitedRate.C b/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationDiffusionLimitedRate/COxidationDiffusionLimitedRate.C
index 6c44dc9da86e5cbe7ad406dc5799cfa32b07e361..2e252fbf596518b276d552225c0fc55353f48cc1 100644
--- a/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationDiffusionLimitedRate/COxidationDiffusionLimitedRate.C
+++ b/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationDiffusionLimitedRate/COxidationDiffusionLimitedRate.C
@@ -50,10 +50,10 @@ Foam::COxidationDiffusionLimitedRate<CloudType>::COxidationDiffusionLimitedRate
     CsLocalId_ = owner.composition().localId(idSolid, "C");
 
     // Set local copies of thermo properties
-    WO2_ = owner.composition().carrierSpecies()[O2GlobalId_].W();
-    scalar WCO2 = owner.composition().carrierSpecies()[CO2GlobalId_].W();
+    WO2_ = owner.mcCarrierThermo().speciesData()[O2GlobalId_].W();
+    scalar WCO2 = owner.mcCarrierThermo().speciesData()[CO2GlobalId_].W();
     WC_ = WCO2 - WO2_;
-    HcCO2_ = owner.composition().carrierSpecies()[CO2GlobalId_].Hc();
+    HcCO2_ = owner.mcCarrierThermo().speciesData()[CO2GlobalId_].Hc();
 
     if (Sb_ < 0)
     {
@@ -120,8 +120,7 @@ Foam::scalar Foam::COxidationDiffusionLimitedRate<CloudType>::calculate
     }
 
     // Local mass fraction of O2 in the carrier phase
-    const scalar YO2 =
-        this->owner().carrierThermo().composition().Y(O2GlobalId_)[cellI];
+    const scalar YO2 = this->owner().mcCarrierThermo().Y(O2GlobalId_)[cellI];
 
     // Change in C mass [kg]
     scalar dmC =
diff --git a/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationKineticDiffusionLimitedRate/COxidationKineticDiffusionLimitedRate.C b/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationKineticDiffusionLimitedRate/COxidationKineticDiffusionLimitedRate.C
index d21efd77e8b0b612b8a4c0d8bb22d3ce41a894f9..682bac3ec67c9bd4eb4634daa0e7315ea9d24140 100644
--- a/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationKineticDiffusionLimitedRate/COxidationKineticDiffusionLimitedRate.C
+++ b/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationKineticDiffusionLimitedRate/COxidationKineticDiffusionLimitedRate.C
@@ -58,10 +58,10 @@ COxidationKineticDiffusionLimitedRate
     CsLocalId_ = owner.composition().localId(idSolid, "C");
 
     // Set local copies of thermo properties
-    WO2_ = owner.composition().carrierSpecies()[O2GlobalId_].W();
-    scalar WCO2 = owner.composition().carrierSpecies()[CO2GlobalId_].W();
+    WO2_ = owner.mcCarrierThermo().speciesData()[O2GlobalId_].W();
+    scalar WCO2 = owner.mcCarrierThermo().speciesData()[CO2GlobalId_].W();
     WC_ = WCO2 - WO2_;
-    HcCO2_ = owner.composition().carrierSpecies()[CO2GlobalId_].Hc();
+    HcCO2_ = owner.mcCarrierThermo().speciesData()[CO2GlobalId_].Hc();
 
     if (Sb_ < 0)
     {
@@ -128,8 +128,7 @@ Foam::scalar Foam::COxidationKineticDiffusionLimitedRate<CloudType>::calculate
     }
 
     // Local mass fraction of O2 in the carrier phase
-    const scalar YO2 =
-        this->owner().carrierThermo().composition().Y(O2GlobalId_)[cellI];
+    const scalar YO2 = this->owner().mcCarrierThermo().Y(O2GlobalId_)[cellI];
 
     // Diffusion rate coefficient
     const scalar D0 = C1_/d*pow(0.5*(T + Tc), 0.75);
diff --git a/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationMurphyShaddix/COxidationMurphyShaddix.C b/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationMurphyShaddix/COxidationMurphyShaddix.C
index c7c8dd795582d9ce079240f1b65bf2bc69740080..e2a1c68700d84c7f2ad2470881d5ee946ee9d127 100644
--- a/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationMurphyShaddix/COxidationMurphyShaddix.C
+++ b/src/lagrangian/coalCombustion/submodels/surfaceReactionModel/COxidationMurphyShaddix/COxidationMurphyShaddix.C
@@ -69,8 +69,8 @@ Foam::COxidationMurphyShaddix<CloudType>::COxidationMurphyShaddix
     CsLocalId_ = owner.composition().localId(idSolid, "C");
 
     // Set local copies of thermo properties
-    WO2_ = owner.composition().carrierSpecies()[O2GlobalId_].W();
-    scalar WCO2 = owner.composition().carrierSpecies()[CO2GlobalId_].W();
+    WO2_ = owner.mcCarrierThermo().speciesData()[O2GlobalId_].W();
+    scalar WCO2 = owner.mcCarrierThermo().speciesData()[CO2GlobalId_].W();
     WC_ = WCO2 - WO2_;
 }
 
@@ -125,7 +125,7 @@ Foam::scalar Foam::COxidationMurphyShaddix<CloudType>::calculate
 
     // Cell carrier phase O2 species density [kg/m^3]
     const scalar rhoO2 =
-        rhoc*this->owner().carrierThermo().composition().Y(O2GlobalId_)[cellI];
+        rhoc*this->owner().mcCarrierThermo().Y(O2GlobalId_)[cellI];
 
     if (rhoO2 < SMALL)
     {
@@ -211,9 +211,9 @@ Foam::scalar Foam::COxidationMurphyShaddix<CloudType>::calculate
         this->owner().composition().solids().properties()[CsLocalId_].Hf()
       + this->owner().composition().solids().properties()[CsLocalId_].cp()*T;
     const scalar HCO2 =
-        this->owner().composition().carrierSpecies()[CO2GlobalId_].H(T);
+        this->owner().mcCarrierThermo().speciesData()[CO2GlobalId_].H(T);
     const scalar HO2 =
-        this->owner().composition().carrierSpecies()[O2GlobalId_].H(T);
+        this->owner().mcCarrierThermo().speciesData()[O2GlobalId_].H(T);
 
     // Heat of reaction
     return dOmega*(WC_*HC + WO2_*HO2 - (WC_ + WO2_)*HCO2);
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.C b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.C
index bebeab25a0d846335f8ad2d6e8299fd21ab753c8..d4f4f10c9abe7823f6c99a8fc6605fa8e17ed07b 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.C
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.C
@@ -28,7 +28,6 @@ License
 
 #include "CompositionModel.H"
 #include "PhaseChangeModel.H"
-#include "multiComponentMixture.H"
 
 // * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * * //
 
@@ -68,14 +67,16 @@ Foam::ReactingCloud<ParcelType>::ReactingCloud
     const volScalarField& rho,
     const volVectorField& U,
     const dimensionedVector& g,
-    hCombustionThermo& thermo
+    basicThermo& thermo
 )
 :
     ThermoCloud<ParcelType>(cloudName, rho, U, g, thermo),
     reactingCloud(),
     constProps_(this->particleProperties()),
-    carrierThermo_(thermo),
-    carrierSpecies_(thermo.composition().Y().size()),
+    mcCarrierThermo_
+    (
+        dynamic_cast<multiComponentMixture<thermoType>&>(thermo)
+    ),
     compositionModel_
     (
         CompositionModel<ReactingCloud<ParcelType> >::New
@@ -92,23 +93,9 @@ Foam::ReactingCloud<ParcelType>::ReactingCloud
             *this
         )
     ),
-    rhoTrans_(thermo.composition().Y().size()),
+    rhoTrans_(mcCarrierThermo_.species().size()),
     dMassPhaseChange_(0.0)
 {
-    // Create the carrier species
-    forAll(carrierSpecies_, specieI)
-    {
-        carrierSpecies_.set
-        (
-            specieI,
-            new thermoType
-            (
-                dynamic_cast<const multiComponentMixture<thermoType>&>
-                    (thermo).speciesData()[specieI]
-            )
-        );
-    }
-
     // Set storage for mass source fields and initialise to zero
     forAll(rhoTrans_, i)
     {
@@ -119,9 +106,7 @@ Foam::ReactingCloud<ParcelType>::ReactingCloud
             (
                 IOobject
                 (
-                    this->name()
-                      + "rhoTrans_"
-                      + thermo.composition().Y()[i].name(),
+                    this->name() + "rhoTrans_" + mcCarrierThermo_.species()[i],
                     this->db().time().timeName(),
                     this->db(),
                     IOobject::NO_READ,
@@ -193,9 +178,9 @@ void Foam::ReactingCloud<ParcelType>::resetSourceTerms()
 template<class ParcelType>
 void Foam::ReactingCloud<ParcelType>::evolve()
 {
-    const volScalarField& T = carrierThermo_.T();
-    const volScalarField cp = carrierThermo_.Cp();
-    const volScalarField& p = carrierThermo_.p();
+    const volScalarField& T = this->carrierThermo().T();
+    const volScalarField cp = this->carrierThermo().Cp();
+    const volScalarField& p = this->carrierThermo().p();
 
     autoPtr<interpolation<scalar> > rhoInterp = interpolation<scalar>::New
     (
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.H b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.H
index da346d64ff9c1ac6e56522eac5eab11d7eec9ffe..4a79c03b6f42ca3ec614c5f782bffcff2bd26c44 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloud.H
@@ -41,11 +41,9 @@ SourceFiles
 #ifndef ReactingCloud_H
 #define ReactingCloud_H
 
-#include "autoPtr.H"
-#include "hCombustionThermo.H"
-
 #include "ThermoCloud.H"
 #include "reactingCloud.H"
+#include "multiComponentMixture.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -94,11 +92,8 @@ protected:
         //- Parcel constant properties
         typename ParcelType::constantProperties constProps_;
 
-        //- Thermodynamics package (combustion)
-        hCombustionThermo& carrierThermo_;
-
-        //- Gas phase properties
-        PtrList<thermoType> carrierSpecies_;
+        //- Multi-component carrier phase thermo
+        multiComponentMixture<thermoType>& mcCarrierThermo_;
 
 
         // References to the cloud sub-models
@@ -148,7 +143,7 @@ public:
             const volScalarField& rho,
             const volVectorField& U,
             const dimensionedVector& g,
-            hCombustionThermo& thermo
+            basicThermo& thermo
         );
 
 
@@ -168,14 +163,12 @@ public:
             inline const typename ParcelType::constantProperties&
                 constProps() const;
 
-            //- Return const access to carrier phase thermo package
-            inline const hCombustionThermo& carrierThermo() const;
-
-            //- Return access to carrier phase thermo package
-            inline hCombustionThermo& carrierThermo();
+            //- Return const access to multi-component carrier phase thermo
+            inline const multiComponentMixture<thermoType>&
+                mcCarrierThermo() const;
 
-            //- Gas phase properties
-            inline const PtrList<thermoType>& carrierSpecies() const;
+            //- Return access to multi-component carrier phase thermo
+            inline multiComponentMixture<thermoType>& mcCarrierThermo();
 
 
             // Sub-models
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloudI.H b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloudI.H
index 8dbaf52bbcf81675641d903cc1be6b8e171bde5a..b59686482706c47710e3041ad8b4f85e13b693c4 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloudI.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingCloud/ReactingCloudI.H
@@ -35,26 +35,18 @@ Foam::ReactingCloud<ParcelType>::constProps() const
 
 
 template<class ParcelType>
-inline const Foam::hCombustionThermo&
-Foam::ReactingCloud<ParcelType>::carrierThermo() const
+inline const Foam::multiComponentMixture<typename ParcelType::thermoType>&
+Foam::ReactingCloud<ParcelType>::mcCarrierThermo() const
 {
-    return carrierThermo_;
+    return mcCarrierThermo_;
 }
 
 
 template<class ParcelType>
-inline Foam::hCombustionThermo&
-Foam::ReactingCloud<ParcelType>::carrierThermo()
+inline Foam::multiComponentMixture<typename ParcelType::thermoType>&
+Foam::ReactingCloud<ParcelType>::mcCarrierThermo()
 {
-    return carrierThermo_;
-}
-
-
-template<class ParcelType>
-inline const Foam::PtrList<typename ParcelType::thermoType>&
-Foam::ReactingCloud<ParcelType>::carrierSpecies() const
-{
-    return carrierSpecies_;
+    return mcCarrierThermo_;
 }
 
 
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.C b/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.C
index 36521105b5d6b0da2675ccebaf82843227e83a7c..d07db7ed622a5b6efb51f27c2436da58b6f9e478 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.C
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.C
@@ -38,7 +38,7 @@ Foam::ReactingMultiphaseCloud<ParcelType>::ReactingMultiphaseCloud
     const volScalarField& rho,
     const volVectorField& U,
     const dimensionedVector& g,
-    hCombustionThermo& thermo
+    basicThermo& thermo
 )
 :
     ReactingCloud<ParcelType>(cloudName, rho, U, g, thermo),
diff --git a/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.H b/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.H
index e99375cf2316be171998ad65508ec8035d02e16d..525b03983cabf044fbd7ec87ceb5d43d4e151005 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ReactingMultiphaseCloud/ReactingMultiphaseCloud.H
@@ -123,7 +123,7 @@ public:
             const volScalarField& rho,
             const volVectorField& U,
             const dimensionedVector& g,
-            hCombustionThermo& thermo
+            basicThermo& thermo
         );
 
 
diff --git a/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.H b/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.H
index aaebdff80ce2859f9c6fe96e86080fb7d78bfa07..cc8f0fadc8e1d313cf5e1b74fef6eb5ab97ad490 100644
--- a/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.H
+++ b/src/lagrangian/intermediate/clouds/Templates/ThermoCloud/ThermoCloud.H
@@ -40,11 +40,9 @@ SourceFiles
 #ifndef ThermoCloud_H
 #define ThermoCloud_H
 
-#include "autoPtr.H"
-#include "hCombustionThermo.H"
-
 #include "KinematicCloud.H"
 #include "thermoCloud.H"
+#include "basicThermo.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
diff --git a/src/lagrangian/intermediate/clouds/derived/BasicReactingCloud/BasicReactingCloud.C b/src/lagrangian/intermediate/clouds/derived/BasicReactingCloud/BasicReactingCloud.C
index 1a639b5820b0d008932fe9454bc1df9e62bc0cee..e59ca18ecf7dcf98670722534fc12059f9cd6665 100644
--- a/src/lagrangian/intermediate/clouds/derived/BasicReactingCloud/BasicReactingCloud.C
+++ b/src/lagrangian/intermediate/clouds/derived/BasicReactingCloud/BasicReactingCloud.C
@@ -35,7 +35,7 @@ Foam::BasicReactingCloud<ThermoType>::BasicReactingCloud
     const volScalarField& rho,
     const volVectorField& U,
     const dimensionedVector& g,
-    hCombustionThermo& thermo
+    basicThermo& thermo
 )
 :
     ReactingCloud<BasicReactingParcel<ThermoType> >
diff --git a/src/lagrangian/intermediate/clouds/derived/BasicReactingCloud/BasicReactingCloud.H b/src/lagrangian/intermediate/clouds/derived/BasicReactingCloud/BasicReactingCloud.H
index aa3048e6663f3db16cbad66a4392a13f397ac2a6..49a51a230e11e72d0e413366a0a0759e30ab2a1d 100644
--- a/src/lagrangian/intermediate/clouds/derived/BasicReactingCloud/BasicReactingCloud.H
+++ b/src/lagrangian/intermediate/clouds/derived/BasicReactingCloud/BasicReactingCloud.H
@@ -81,7 +81,7 @@ public:
             const volScalarField& rho,
             const volVectorField& U,
             const dimensionedVector& g,
-            hCombustionThermo& thermo
+            basicThermo& thermo
         );
 
 
diff --git a/src/lagrangian/intermediate/clouds/derived/BasicReactingMultiphaseCloud/BasicReactingMultiphaseCloud.C b/src/lagrangian/intermediate/clouds/derived/BasicReactingMultiphaseCloud/BasicReactingMultiphaseCloud.C
index 07cb83587ecd620b23e26fb2d874afdbf7a66ee2..6f08db20167346b0b1fc41e9c669661c38d39861 100644
--- a/src/lagrangian/intermediate/clouds/derived/BasicReactingMultiphaseCloud/BasicReactingMultiphaseCloud.C
+++ b/src/lagrangian/intermediate/clouds/derived/BasicReactingMultiphaseCloud/BasicReactingMultiphaseCloud.C
@@ -35,8 +35,7 @@ Foam::BasicReactingMultiphaseCloud<ThermoType>::BasicReactingMultiphaseCloud
     const volScalarField& rho,
     const volVectorField& U,
     const dimensionedVector& g,
-    hCombustionThermo& thermo,
-    PtrList<ThermoType>& carrierSpecies
+    basicThermo& thermo
 )
 :
     ReactingMultiphaseCloud<BasicReactingMultiphaseParcel<ThermoType> >
@@ -45,8 +44,7 @@ Foam::BasicReactingMultiphaseCloud<ThermoType>::BasicReactingMultiphaseCloud
         rho,
         U,
         g,
-        thermo,
-        carrierSpecies
+        thermo
     )
 {
     BasicReactingMultiphaseParcel<ThermoType>::readFields(*this);
diff --git a/src/lagrangian/intermediate/clouds/derived/BasicReactingMultiphaseCloud/BasicReactingMultiphaseCloud.H b/src/lagrangian/intermediate/clouds/derived/BasicReactingMultiphaseCloud/BasicReactingMultiphaseCloud.H
index 7b093c529c03d7e2c9c8871b493cc93ae46f3293..ac0853823fce7f97dc2ca15c2afa5bd7520b44af 100644
--- a/src/lagrangian/intermediate/clouds/derived/BasicReactingMultiphaseCloud/BasicReactingMultiphaseCloud.H
+++ b/src/lagrangian/intermediate/clouds/derived/BasicReactingMultiphaseCloud/BasicReactingMultiphaseCloud.H
@@ -81,8 +81,7 @@ public:
             const volScalarField& rho,
             const volVectorField& U,
             const dimensionedVector& g,
-            hCombustionThermo& thermo,
-            PtrList<ThermoType>& carrierSpecies
+            basicThermo& thermo
         );
 
 
diff --git a/src/lagrangian/intermediate/clouds/derived/BasicTrackedReactingCloud/BasicTrackedReactingCloud.C b/src/lagrangian/intermediate/clouds/derived/BasicTrackedReactingCloud/BasicTrackedReactingCloud.C
index 4b94788812116369124865291873fe467239e525..15b7774a0b86e74e817de607189194721533ff6d 100644
--- a/src/lagrangian/intermediate/clouds/derived/BasicTrackedReactingCloud/BasicTrackedReactingCloud.C
+++ b/src/lagrangian/intermediate/clouds/derived/BasicTrackedReactingCloud/BasicTrackedReactingCloud.C
@@ -35,7 +35,7 @@ Foam::BasicTrackedReactingCloud<ThermoType>::BasicTrackedReactingCloud
     const volScalarField& rho,
     const volVectorField& U,
     const dimensionedVector& g,
-    hCombustionThermo& thermo
+    basicThermo& thermo
 )
 :
     ReactingCloud<BasicTrackedReactingParcel<ThermoType> >
diff --git a/src/lagrangian/intermediate/clouds/derived/BasicTrackedReactingCloud/BasicTrackedReactingCloud.H b/src/lagrangian/intermediate/clouds/derived/BasicTrackedReactingCloud/BasicTrackedReactingCloud.H
index 5c3cd8148f2654b2588fb4cd3282768268cdb8a5..bb2b0bde1efa0f3b07b908ed92893ae298ada56b 100644
--- a/src/lagrangian/intermediate/clouds/derived/BasicTrackedReactingCloud/BasicTrackedReactingCloud.H
+++ b/src/lagrangian/intermediate/clouds/derived/BasicTrackedReactingCloud/BasicTrackedReactingCloud.H
@@ -81,7 +81,7 @@ public:
             const volScalarField& rho,
             const volVectorField& U,
             const dimensionedVector& g,
-            hCombustionThermo& thermo
+            basicThermo& thermo
         );
 
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
index 695ff95a8ef96e42020adfe63d7da93834d4ac93..682cad335ad2c0c9ad36aa89c11f5799adfc47dc 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
@@ -168,7 +168,8 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::cellValueSourceCorrection
         forAll(td.cloud().rhoTrans(), i)
         {
             scalar Y = td.cloud().rhoTrans(i)[cellI]/addedMass;
-            cpEff += Y*td.cloud().carrierSpecies()[i].Cp(this->Tc_);
+            cpEff +=
+                Y*td.cloud().mcCarrierThermo().speciesData()[i].Cp(this->Tc_);
         }
     }
     const scalar cpc = td.cpInterp().psi()[cellI];
@@ -248,7 +249,12 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
     scalarField dMassSRGas(YGas_.size(), 0.0);
     scalarField dMassSRLiquid(YLiquid_.size(), 0.0);
     scalarField dMassSRSolid(YSolid_.size(), 0.0);
-    scalarField dMassSRCarrier(td.cloud().carrierSpecies().size(), 0.0);
+    scalarField
+        dMassSRCarrier
+        (
+            td.cloud().mcCarrierThermo().species().size(),
+            0.0
+        );
 
     // Clac mass and enthalpy transfer due to surface reactions
     calcSurfaceReactions
@@ -340,7 +346,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
             td.cloud().hcTrans()[cellI] +=
                 np0
                *dMassGas[i]
-               *td.cloud().composition().carrierSpecies()[gid].H(T0);
+               *td.cloud().mcCarrierThermo().speciesData()[gid].H(T0);
         }
         forAll(YLiquid_, i)
         {
@@ -349,7 +355,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
             td.cloud().hcTrans()[cellI] +=
                 np0
                *dMassLiquid[i]
-               *td.cloud().composition().carrierSpecies()[gid].H(T0);
+               *td.cloud().mcCarrierThermo().speciesData()[gid].H(T0);
         }
 /*
         // No mapping between solid components and carrier phase
@@ -360,7 +366,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
             td.cloud().hcTrans()[cellI] +=
                 np0
                *dMassSolid[i]
-               *td.cloud().composition().carrierSpecies()[gid].H(T0);
+               *td.cloud().mcCarrierThermo().speciesData()[gid].H(T0);
         }
 */
         forAll(dMassSRCarrier, i)
@@ -369,7 +375,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
             td.cloud().hcTrans()[cellI] +=
                 np0
                *dMassSRCarrier[i]
-               *td.cloud().composition().carrierSpecies()[i].H(T0);
+               *td.cloud().mcCarrierThermo().speciesData()[i].H(T0);
         }
 
         // Update momentum transfer
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
index b52bfeef29f4c81dc18d8f038fc42a6198b197c3..e4137e90dc6eb1f46e70df5414c1fcb1a218ee4d 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
@@ -86,7 +86,8 @@ void Foam::ReactingParcel<ParcelType>::cellValueSourceCorrection
         forAll(td.cloud().rhoTrans(), i)
         {
             scalar Y = td.cloud().rhoTrans(i)[cellI]/addedMass;
-            cpEff += Y*td.cloud().carrierSpecies()[i].Cp(this->Tc_);
+            cpEff +=
+                Y*td.cloud().mcCarrierThermo().speciesData()[i].Cp(this->Tc_);
         }
     }
     const scalar cpc = td.cpInterp().psi()[cellI];
@@ -190,9 +191,6 @@ void Foam::ReactingParcel<ParcelType>::calc
     // Motion
     // ~~~~~~
 
-    // No additional forces
-    vector Fx = vector::zero;
-
     // Calculate new particle velocity
     vector U1 = calcVelocity(td, dt, cellI, d0, U0, rho0, mass0, Su, dUTrans);
 
@@ -210,7 +208,7 @@ void Foam::ReactingParcel<ParcelType>::calc
             td.cloud().hcTrans()[cellI] +=
                 np0
                *dMassPC[i]
-               *td.cloud().composition().carrierSpecies()[gid].H(T0);
+               *td.cloud().mcCarrierThermo().speciesData()[gid].H(T0);
         }
 
         // Update momentum transfer
@@ -317,16 +315,15 @@ void Foam::ReactingParcel<ParcelType>::calcPhaseChange
     td.cloud().addToMassPhaseChange(this->nParticle_*dMassTot);
 
     // Enthalphy transfer to carrier phase
+    label id;
     forAll(YComponents, i)
     {
-        label gid;
-
-        gid = td.cloud().composition().localToGlobalCarrierId(idPhase, i);
-        const scalar hv = td.cloud().composition().carrierSpecies()[gid].H(T);
+        id = td.cloud().composition().localToGlobalCarrierId(idPhase, i);
+        const scalar hv = td.cloud().mcCarrierThermo().speciesData()[id].H(T);
 
-        gid = td.cloud().composition().globalIds(idPhase)[i];
+        id = td.cloud().composition().globalIds(idPhase)[i];
         const scalar hl =
-            td.cloud().composition().liquids().properties()[gid].h(pc_, T);
+            td.cloud().composition().liquids().properties()[id].h(pc_, T);
 
         Sh += dMassPC[i]*(hl - hv)/dt;
     }
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
index 9ae011160404a5c051b39796812c3a8bfce51cf4..6ae3dabb6fd1ab3ac3efac2f5526adee30e0f3d1 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
@@ -39,15 +39,14 @@ Foam::CompositionModel<CloudType>::CompositionModel
     dict_(dict),
     owner_(owner),
     coeffDict_(dict.subDict(type + "Coeffs")),
-    carrierThermo_(owner.carrierThermo()),
-    carrierSpecies_(owner.carrierSpecies()),
+    mcCarrierThermo_(owner.mcCarrierThermo()),
     liquids_
     (
         liquidMixture::New
         (
             owner.mesh().objectRegistry::lookupObject<dictionary>
             (
-                carrierThermo_.name()
+                owner.carrierThermo().name()
             )
         )
     ),
@@ -57,14 +56,14 @@ Foam::CompositionModel<CloudType>::CompositionModel
         (
             owner.mesh().objectRegistry::lookupObject<dictionary>
             (
-                carrierThermo_.name()
+                owner.carrierThermo().name()
             )
         )
     ),
     phaseProps_
     (
         coeffDict_.lookup("phases"),
-        carrierThermo_.composition().species(),
+        mcCarrierThermo_.species(),
         liquids_().components(),
         solids_().components()
     )
@@ -102,18 +101,10 @@ const Foam::dictionary& Foam::CompositionModel<CloudType>::coeffDict() const
 
 
 template<class CloudType>
-const Foam::hCombustionThermo&
-Foam::CompositionModel<CloudType>::carrierThermo() const
+const Foam::multiComponentMixture<typename CloudType::thermoType>&
+Foam::CompositionModel<CloudType>::mcCarrierThermo() const
 {
-    return carrierThermo_;
-}
-
-
-template<class CloudType>
-const Foam::PtrList<typename CloudType::thermoType>&
-Foam::CompositionModel<CloudType>::carrierSpecies() const
-{
-    return carrierSpecies_;
+    return mcCarrierThermo_;
 }
 
 
@@ -184,11 +175,9 @@ Foam::label Foam::CompositionModel<CloudType>::globalCarrierId
     const word& cmptName
 ) const
 {
-    forAll(carrierThermo_.composition().species(), i)
+    forAll(mcCarrierThermo_.species(), i)
     {
-        const word& carrierSpecieName =
-            carrierThermo_.composition().species()[i];
-        if (cmptName == carrierSpecieName)
+        if (cmptName == mcCarrierThermo_.species()[i])
         {
             return i;
         }
@@ -324,8 +313,8 @@ Foam::scalarField Foam::CompositionModel<CloudType>::X
             forAll(Y, i)
             {
                 label gid = props.globalIds()[i];
-                WInv += Y[i]/this->carrierSpecies()[gid].W();
-                X[i] = Y[i]/this->carrierSpecies()[gid].W();
+                WInv += Y[i]/mcCarrierThermo_.speciesData()[gid].W();
+                X[i] = Y[i]/mcCarrierThermo_.speciesData()[gid].W();
             }
             break;
         }
@@ -375,7 +364,7 @@ Foam::scalar Foam::CompositionModel<CloudType>::H
             forAll(Y, i)
             {
                 label gid = props.globalIds()[i];
-                HMixture += Y[i]*this->carrierSpecies()[gid].H(T);
+                HMixture += Y[i]*mcCarrierThermo_.speciesData()[gid].H(T);
             }
             break;
         }
@@ -439,7 +428,7 @@ Foam::scalar Foam::CompositionModel<CloudType>::cp
             forAll(Y, i)
             {
                 label gid = props.globalIds()[i];
-                cpMixture += Y[i]*this->carrierSpecies()[gid].Cp(T);
+                cpMixture += Y[i]*mcCarrierThermo_.speciesData()[gid].Cp(T);
             }
             break;
         }
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
index 1e0b800b45709680d18fc421aba7f9b5805e02df..b3f0949defcc84b06f7f8c4b4b0a9ff900d1b079 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
@@ -44,7 +44,7 @@ SourceFiles
 #include "runTimeSelectionTables.H"
 
 #include "PtrList.H"
-#include "hCombustionThermo.H"
+#include "multiComponentMixture.H"
 
 #include "liquidMixture.H"
 #include "solidMixture.H"
@@ -74,11 +74,8 @@ class CompositionModel
         //- The coefficients dictionary
         const dictionary& coeffDict_;
 
-        //- Reference to the carrier phase thermo package
-        hCombustionThermo& carrierThermo_;
-
-        //- Reference to the carrier phase species
-        const PtrList<typename CloudType::thermoType>& carrierSpecies_;
+        //- Reference to the multi-component carrier phase thermo
+        multiComponentMixture<typename CloudType::thermoType>& mcCarrierThermo_;
 
         //- Global (additional) liquid properties data
         autoPtr<liquidMixture> liquids_;
@@ -146,15 +143,12 @@ public:
             const dictionary& coeffDict() const;
 
             //- Return the carrier phase thermo package
-            const hCombustionThermo& carrierThermo() const;
+            const multiComponentMixture<typename CloudType::thermoType>&
+                mcCarrierThermo() const;
 
 
             // Composition lists
 
-                //- Return the carrier species
-                const PtrList<typename CloudType::thermoType>&
-                    carrierSpecies() const;
-
                 //- Return the global (additional) liquids
                 const liquidMixture& liquids() const;
 
diff --git a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.C b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.C
index b7caef3652d1f71f8b2406c6f4100adb5bd1f905..c486ccaa00a85d3c4d8b139da2824172f71aed98 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/LiquidEvaporation/LiquidEvaporation.C
@@ -35,14 +35,14 @@ Foam::scalarField Foam::LiquidEvaporation<CloudType>::calcXc
     const label cellI
 ) const
 {
-    scalarField Xc(this->owner().carrierThermo().composition().Y().size());
+    scalarField Xc(this->owner().mcCarrierThermo().Y().size());
 
     scalar Winv = 0.0;
     forAll(Xc, i)
     {
-        scalar Y = this->owner().carrierThermo().composition().Y()[i][cellI];
-        Winv += Y/this->owner().carrierSpecies()[i].W();
-        Xc[i] = Y/this->owner().carrierSpecies()[i].W();
+        scalar Y = this->owner().mcCarrierThermo().Y()[i][cellI];
+        Winv += Y/this->owner().mcCarrierThermo().speciesData()[i].W();
+        Xc[i] = Y/this->owner().mcCarrierThermo().speciesData()[i].W();
     }
 
     return Xc/Winv;
@@ -76,7 +76,7 @@ Foam::LiquidEvaporation<CloudType>::LiquidEvaporation
         (
             owner.mesh().objectRegistry::lookupObject<dictionary>
             (
-                owner.carrierThermo().name()
+                "thermophysicalProperties"
             )
         )
     ),
@@ -104,7 +104,7 @@ Foam::LiquidEvaporation<CloudType>::LiquidEvaporation
             owner.composition().globalCarrierId(activeLiquids_[i]);
     }
 
-    // Determine mapping between local and global liquids
+    // Determine mapping between model active liquids and global liquids
     label idLiquid = owner.composition().idLiquid();
     forAll(activeLiquids_, i)
     {
diff --git a/src/thermophysicalModels/reactionThermo/combustionThermo/hCombustionThermo/newhCombustionThermo.C b/src/thermophysicalModels/reactionThermo/combustionThermo/hCombustionThermo/newhCombustionThermo.C
index 175c7d4756bb432ad09067c01b481d56f72ae100..788a298d0d841b2d09d258cb2c35a904dd254bd7 100644
--- a/src/thermophysicalModels/reactionThermo/combustionThermo/hCombustionThermo/newhCombustionThermo.C
+++ b/src/thermophysicalModels/reactionThermo/combustionThermo/hCombustionThermo/newhCombustionThermo.C
@@ -103,6 +103,16 @@ Foam::autoPtr<Foam::hCombustionThermo> Foam::hCombustionThermo::NewType
 
         if (hCombustionThermoTypeName.find(thermoType) == string::npos)
         {
+            wordList allModels = fvMeshConstructorTablePtr_->toc();
+            DynamicList<word> validModels;
+            forAll(allModels, i)
+            {
+                if (allModels[i].find(thermoType) != string::npos)
+                {
+                    validModels.append(allModels[i]);
+                }
+            }
+
             FatalErrorIn
             (
                 "autoPtr<hCombustionThermo> hCombustionThermo::NewType"
@@ -112,7 +122,8 @@ Foam::autoPtr<Foam::hCombustionThermo> Foam::hCombustionThermo::NewType
                 ")"
             )   << "Inconsistent thermo package selected:" << nl << nl
                 << hCombustionThermoTypeName << nl << nl << "Please select a "
-                << "thermo package based on " << thermoType << nl << nl
+                << "thermo package based on " << thermoType
+                << ". Valid options include:" << nl << validModels << nl
                 << exit(FatalError);
         }
     }
diff --git a/src/thermophysicalModels/reactionThermo/reactionThermo/hReactionThermo/newhReactionThermo.C b/src/thermophysicalModels/reactionThermo/reactionThermo/hReactionThermo/newhReactionThermo.C
index ceeba8922558b19de31d2c65191e3702d5019c35..688b37b6d2a726a82771119fd896c3649bcf59df 100644
--- a/src/thermophysicalModels/reactionThermo/reactionThermo/hReactionThermo/newhReactionThermo.C
+++ b/src/thermophysicalModels/reactionThermo/reactionThermo/hReactionThermo/newhReactionThermo.C
@@ -103,6 +103,16 @@ Foam::autoPtr<Foam::hReactionThermo> Foam::hReactionThermo::NewType
 
         if (hReactionThermoTypeName.find(thermoType) == string::npos)
         {
+            wordList allModels = fvMeshConstructorTablePtr_->toc();
+            DynamicList<word> validModels;
+            forAll(allModels, i)
+            {
+                if (allModels[i].find(thermoType) != string::npos)
+                {
+                    validModels.append(allModels[i]);
+                }
+            }
+
             FatalErrorIn
             (
                 "autoPtr<hReactionThermo> hReactionThermo::NewType"
@@ -112,7 +122,8 @@ Foam::autoPtr<Foam::hReactionThermo> Foam::hReactionThermo::NewType
                 ")"
             )   << "Inconsistent thermo package selected:" << nl << nl
                 << hReactionThermoTypeName << nl << nl << "Please select a "
-                << "thermo package based on " << thermoType << nl << nl
+                << "thermo package based on " << thermoType
+                << ". Valid options include:" << nl << validModels << nl
                 << exit(FatalError);
         }
     }
diff --git a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/0/epsilon b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/0/epsilon
index e5ed472a8b8b551fb4ab7245a9318c590c271030..9fcf588c0c292537815ea7ef2c72563a85a14851 100644
--- a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/0/epsilon
+++ b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/0/epsilon
@@ -1,7 +1,7 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5.x                                 |
+|  \\    /   O peration     | Version:  dev                                   |
 |   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
diff --git a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/0/k b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/0/k
index 1d9c0d007fb294979895e61912032922c75639cd..9abb57fdb013a1d9dd1ac7859f9554dc0c642b71 100644
--- a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/0/k
+++ b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/0/k
@@ -1,7 +1,7 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
-|  \\    /   O peration     | Version:  1.5.x                                 |
+|  \\    /   O peration     | Version:  dev                                   |
 |   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
diff --git a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/chemistryProperties b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/chemistryProperties
index 53a5860d33ad658f415a68c944c3a89007fe22d3..cbffc12e84975d2b2e7be95277458e5042f53ea7 100644
--- a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/chemistryProperties
+++ b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/chemistryProperties
@@ -15,11 +15,13 @@ FoamFile
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+psiChemistryModel ODEChemistryModel<gasThermoPhysics>;
+
 chemistry       off;
 
 turbulentReaction off;
 
-chemistrySolver ODE;
+chemistrySolver ode;
 
 initialChemicalTimeStep 1e-07;
 
@@ -34,7 +36,7 @@ EulerImplicitCoeffs
     equilibriumRateLimiter off;
 }
 
-ODECoeffs
+odeCoeffs
 {
     ODESolver       RK;
     eps             0.05;
diff --git a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/thermophysicalProperties b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/thermophysicalProperties
index 04b13b12cf071d823be29a6c41dc37ec5ab58875..8b80f059a3230223aa3de5ae7ffb3c76ead945e0 100644
--- a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/thermophysicalProperties
+++ b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/thermophysicalProperties
@@ -15,7 +15,7 @@ FoamFile
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-thermoType      hMixtureThermo<reactingMixture>;
+thermoType      hPsiMixtureThermo<reactingMixture<gasThermoPhysics>>;
 
 chemistryReader foamChemistryReader;