diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
index 35d969714f8d3d7280a1487509ba158fdf98ab0a..0ab68ded742c8c3e17c2d0f93b3d8c14c9d4dafc 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
@@ -100,7 +100,8 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calcCoupled
     // ~~~~~~~~~~~~~~~~~~~~~~
     // Calculate phase change
     // ~~~~~~~~~~~~~~~~~~~~~~
-    calcPhaseChange(td, dt, T0, dMassMT);
+    scalarField X = td.cloud().composition().X(idLiquid, YLiquid_);
+    calcPhaseChange(td, dt, T0, X, dMassMT);
 
 
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~
@@ -271,7 +272,8 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calcUncoupled
     // ~~~~~~~~~~~~~~~~~~~~~~
     // Calculate phase change
     // ~~~~~~~~~~~~~~~~~~~~~~
-    calcPhaseChange(td, dt, T0, dMassMT);
+    scalarField X = td.cloud().composition().X(idLiquid, YLiquid_);
+    calcPhaseChange(td, dt, T0, YLiquid_, dMassMT);
 
 
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~
@@ -347,7 +349,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calcDevolatilisation
         (
             "void Foam::ReactingMultiphaseParcel<ParcelType>::"
             "calcDevolatilisation \n"
-            "("
+            "(\n"
             "    TrackData&,\n"
             "    const scalar,\n"
             "    const scalar,\n"
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
index daaa9ac66af65e9d40c6886d82228b0f4ff82ca2..266f6e25301424d1cf46a7cfee223abf78b2aa44 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
@@ -97,7 +97,8 @@ void Foam::ReactingParcel<ParcelType>::calcCoupled
     // ~~~~~~~~~~~~~~~~~~~~~~
     // Calculate phase change
     // ~~~~~~~~~~~~~~~~~~~~~~
-    calcPhaseChange(td, dt, T, dMassMT);
+    scalarField X = td.cliud().composition().X(0, YMixture_);
+    calcPhaseChange(td, dt, T, X, dMassMT);
 
 
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
@@ -229,7 +230,8 @@ void Foam::ReactingParcel<ParcelType>::calcUncoupled
     // ~~~~~~~~~~~~~~~~~~~~~~
     // Calculate phase change
     // ~~~~~~~~~~~~~~~~~~~~~~
-    calcPhaseChange(td, dt, T, dMassMT);
+    scalarField X = td.cloud().composition().X(0, YMixture_);
+    calcPhaseChange(td, dt, T, X, dMassMT);
 
 
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~
@@ -286,6 +288,7 @@ void Foam::ReactingParcel<ParcelType>::calcPhaseChange
     TrackData& td,
     const scalar dt,
     const scalar T,
+    scalarField& X,
     scalarList& dMassMT
 )
 {
@@ -296,35 +299,22 @@ void Foam::ReactingParcel<ParcelType>::calcPhaseChange
 
     // TODO: separate treatment for boiling
 
-    /*
-    // Determine mass to add to carrier phase
-    const scalar mass = this->mass();
-    const scalar dMassTot = td.cloud().devolatilisation().calculate
+    scalar dMassTot = td.cloud().phaseChange().calculate
     (
-        dt,
-        mass0_,
-        mass,
-        td.cloud().composition().YMixture0(),
-        YMixture_,
-        T0,
-        canCombust_
+        T,
+        this->d_,
+        X,
+        dMassMT,
+        this->U_ - this->Uc_,
+        this->Tc_,
+        pc_,
+        this->muc_/this->rhoc_,
+        dt
     );
 
-    // Update (total) mass fractions
-    YMixture_[0] = (YMixture_[0]*mass - dMassTot)/(mass - dMassTot);
-    YMixture_[1] = YMixture_[1]*mass/(mass - dMassTot);
-    YMixture_[2] = 1.0 - YMixture_[0] - YMixture_[1];
+    // TODO: Re-calculate mass fractions
 
-    // Add to cummulative mass transfer
-    forAll (YGas_, i)
-    {
-        label id = td.cloud().composition().gasGlobalIds()[i];
 
-        // Volatiles mass transfer
-        scalar volatileMass = YGas_[i]*dMassTot;
-        dMassMT[id] += volatileMass;
-    }
-        */
 }
 
 
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.H b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.H
index 763f1e2a6c1d0bf31e7c0b505c6961d3ff17639d..194ffc2d3864e36d693fba975f939ee6c467bf5f 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.H
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.H
@@ -175,6 +175,7 @@ protected:
             TrackData& td,
             const scalar dt,
             const scalar T,
+            scalarField& X,
             scalarList& dMassMT
         );
 
diff --git a/src/lagrangian/intermediate/phaseProperties/phaseProperties/phasePropertiesIO.C b/src/lagrangian/intermediate/phaseProperties/phaseProperties/phasePropertiesIO.C
index 34641fe4d4bca52f4b807b23d7943443d39c7e30..192136ed8d070071da361b74768c8fababe341be 100644
--- a/src/lagrangian/intermediate/phaseProperties/phaseProperties/phasePropertiesIO.C
+++ b/src/lagrangian/intermediate/phaseProperties/phaseProperties/phasePropertiesIO.C
@@ -29,7 +29,6 @@ License
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-
 Foam::phaseProperties::phaseProperties(Istream& is)
 :
     phase_(UNKNOWN),
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
index 5d1de160e8db82a3644f8f72f0a8ce371b94e4a9..df898464fda582dc1505293fcc4aa0d78735c2ae 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.C
@@ -241,6 +241,56 @@ const Foam::scalarField& Foam::CompositionModel<CloudType>::Y0
 }
 
 
+template<class CloudType>
+Foam::scalarField Foam::CompositionModel<CloudType>::X
+(
+    const label phaseI,
+    const scalarField& Y
+) const
+{
+    const phaseProperties& props = phaseProps_[phaseI];
+    scalarField X(Y.size(), 0.0);
+    scalar WInv = 0.0;
+    switch (props.phase())
+    {
+        case phaseProperties::GAS:
+        {
+            forAll(Y, i)
+            {
+                label gid = props.globalIds()[i];
+                WInv += Y[i]/this->gases()[gid].W();
+                X[i] = Y[i]/this->gases()[gid].W();
+            }
+            break;
+        }
+        case phaseProperties::LIQUID:
+        {
+            forAll(Y, i)
+            {
+                label gid = props.globalIds()[i];
+                WInv += Y[i]/this->liquids().properties()[gid].W();
+                X[i] += Y[i]/this->liquids().properties()[gid].W();
+            }
+            break;
+        }
+        default:
+        {
+            FatalErrorIn
+            (
+                "Foam::scalarField Foam::CompositionModel<CloudType>::X\n"
+                "(\n"
+                "    const label phaseI,\n"
+                "    const scalarField Y\n"
+                ") const"
+            )   << "Only possible to convert gas and liquid mass fractions"
+                << nl << abort(FatalError);
+        }
+    }
+
+    return X/WInv;
+}
+
+
 template<class CloudType>
 Foam::scalar Foam::CompositionModel<CloudType>::H
 (
diff --git a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
index 9dd21529c899f3985df171f45f11d9dce5820670..5be0e2c1b9242bc75477c188718c541f04e2bf01 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/CompositionModel/CompositionModel/CompositionModel.H
@@ -187,6 +187,14 @@ public:
                 //- Return the list of phase phaseI mass fractions
                 const scalarField& Y0(const label phaseI) const;
 
+                //- Return the list of phase phaseI volume fractions fractions
+                //  based on supplied mass fractions Y
+                scalarField X
+                (
+                    const label phaseI,
+                    const scalarField& Y
+                ) const;
+
 
             // Mixture properties
 
diff --git a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/NoPhaseChange/NoPhaseChange.C b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/NoPhaseChange/NoPhaseChange.C
index 791512789c47c11e742f889bf2dc5a5d69f30cdf..4029aa0570885ab33509d7804fa2f5369deea657 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/NoPhaseChange/NoPhaseChange.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/NoPhaseChange/NoPhaseChange.C
@@ -61,8 +61,8 @@ Foam::scalar Foam::NoPhaseChange<CloudType>::calculate
     const scalar,
     const scalar,
     const scalarField&,
-    const scalarField&,
-    const vector,
+    const scalarList&,
+    const vector&,
     const scalar,
     const scalar,
     const scalar,
diff --git a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/NoPhaseChange/NoPhaseChange.H b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/NoPhaseChange/NoPhaseChange.H
index 07341719db59b5734e01c5bf215ba87e9c8d6ceb..74b065b75a8123ce6f07bf948a3b6dfa55763932 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/NoPhaseChange/NoPhaseChange.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/NoPhaseChange/NoPhaseChange.H
@@ -79,8 +79,8 @@ public:
             const scalar T,
             const scalar d,
             const scalarField& Xc,
-            const scalarField& dMassMT,
-            const vector Ur,
+            const scalarList& dMassMT,
+            const vector& Ur,
             const scalar Tc,
             const scalar pc,
             const scalar nuc,
diff --git a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/PhaseChangeModel/PhaseChangeModel.H b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/PhaseChangeModel/PhaseChangeModel.H
index b06506d8281b5424fa5efba0bfd74643b44ba3c7..1a97370529e83f32a3f7f857851a4f30d8950ac6 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/PhaseChangeModel/PhaseChangeModel.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/PhaseChangeModel/PhaseChangeModel.H
@@ -138,8 +138,8 @@ public:
             const scalar T,
             const scalar d,
             const scalarField& Xc,
-            const scalarField& dMassMT,
-            const vector Ur,
+            const scalarList& dMassMT,
+            const vector& Ur,
             const scalar Tc,
             const scalar pc,
             const scalar nuc,
diff --git a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/liquidEvaporation/liquidEvaporation.C b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/liquidEvaporation/liquidEvaporation.C
index 792ae89c64078aed4d1d559810bef2a157552960..8bba9d784d1c7ae20051f887d5a2e8124ecae772 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/liquidEvaporation/liquidEvaporation.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/liquidEvaporation/liquidEvaporation.C
@@ -146,14 +146,16 @@ Foam::scalar Foam::liquidEvaporation<CloudType>::calculate
     const scalar T,
     const scalar d,
     const scalarField& Xc,
-    const scalarField& dMassMT,
-    const vector Ur
+    const scalarList& dMassMT,
+    const vector& Ur
     const scalar Tc,
     const scalar pc,
     const scalar nuc
     const scalar dt,
 ) const
 {
+    scalar dMassTot = 0.0;
+
     if (T < Tvap_)
     {
         // not reached temperature threshold
@@ -192,9 +194,13 @@ Foam::scalar Foam::liquidEvaporation<CloudType>::calculate
             scalar Ni = max(kc*(Cs - Cinf), 0.0);
 
             // mass transfer
-            dMassMT[i] -= Ni*A*liquids_.properies()[i].W()*dt;
+            scalar dm = Ni*A*liquids_.properies()[i].W()*dt;
+            dMassMT[i] -= 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 efc97fecabfe7e782b921a052ccf4a133ad1af06..a7e2325f17b5a9aa042797afe3c834796f540cfb 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/liquidEvaporation/liquidEvaporation.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/PhaseChangeModel/liquidEvaporation/liquidEvaporation.H
@@ -111,8 +111,8 @@ public:
             const scalar T,
             const scalar d,
             const scalarField& Xc,
-            const scalarField& dMassMT,
-            const vector Ur,
+            const scalarList& dMassMT,
+            const vector& Ur,
             const scalar Tc,
             const scalar pc,
             const scalar nuc,