diff --git a/src/thermophysicalModels/reactionThermo/mixtures/basicMultiComponentMixture/basicMultiComponentMixture.H b/src/thermophysicalModels/reactionThermo/mixtures/basicMultiComponentMixture/basicMultiComponentMixture.H
index 0a5abd06..fdfae416 100644
--- a/src/thermophysicalModels/reactionThermo/mixtures/basicMultiComponentMixture/basicMultiComponentMixture.H
+++ b/src/thermophysicalModels/reactionThermo/mixtures/basicMultiComponentMixture/basicMultiComponentMixture.H
@@ -63,8 +63,7 @@ protected:
         PtrList<volScalarField> Y_;
 
         //- Mixture-averaged species diffusivities
-        PtrList<volScalarField> D_;
-
+        mutable PtrList<volScalarField> D_;
 
 public:
 
diff --git a/src/thermophysicalModels/reactionThermo/mixtures/multiComponentMixture/multiComponentMixture.C b/src/thermophysicalModels/reactionThermo/mixtures/multiComponentMixture/multiComponentMixture.C
index f10a1ce4..c7f4b6d8 100644
--- a/src/thermophysicalModels/reactionThermo/mixtures/multiComponentMixture/multiComponentMixture.C
+++ b/src/thermophysicalModels/reactionThermo/mixtures/multiComponentMixture/multiComponentMixture.C
@@ -122,6 +122,55 @@ Foam::multiComponentMixture<ThermoType>::multiComponentMixture
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class ThermoType>
+void Foam::multiComponentMixture<ThermoType>::
+calculateDiffusivities
+(
+    const volScalarField &p,
+    const volScalarField &T
+) const
+{
+    const volScalarField& Yi = Y_[0];
+    const scalarField& YiCells = Yi.internalField();
+
+    forAll(YiCells, celli)
+    {
+        const ThermoType &mixture = cellMixture(celli);
+        mixture.calculateDiffusivities(p[celli], T[celli]);
+
+        forAll(speciesData_, i)
+        {
+            volScalarField& Di = D_[i];
+            scalarField& DiCells = Di.internalField();
+            DiCells[celli] = mixture.D(i);
+        }
+    }
+
+    forAll(Yi.boundaryField(), patchi)
+    {
+        const scalarField& Yip = Yi.boundaryField()[patchi];
+        const scalarField& pp = p.boundaryField()[patchi];
+        const scalarField& Tp = T.boundaryField()[patchi];
+
+        forAll(Yip, facei)
+        {
+            const ThermoType &mixture = patchFaceMixture(patchi, facei);
+            mixture.calculateDiffusivities(pp[facei], Tp[facei]);
+
+
+            forAll(speciesData_, i)
+            {
+                volScalarField& Di = D_[i];
+                scalarField& Dip = Di.boundaryField()[patchi];
+                Dip[facei] = mixture.D(i);
+            }
+        }
+    }
+
+    return;
+}
+
+
 template<class ThermoType>
 const ThermoType& Foam::multiComponentMixture<ThermoType>::cellMixture
 (
diff --git a/src/thermophysicalModels/reactionThermo/mixtures/multiComponentMixture/multiComponentMixture.H b/src/thermophysicalModels/reactionThermo/mixtures/multiComponentMixture/multiComponentMixture.H
index c7b50a77..f51bf14e 100644
--- a/src/thermophysicalModels/reactionThermo/mixtures/multiComponentMixture/multiComponentMixture.H
+++ b/src/thermophysicalModels/reactionThermo/mixtures/multiComponentMixture/multiComponentMixture.H
@@ -106,6 +106,9 @@ public:
 
     // Member functions
 
+        //- Calculate mixture-aveaged diffusivities
+        virtual void calculateDiffusivities(const volScalarField& p, const volScalarField& T) const;
+
         const ThermoType& cellMixture(const label celli) const;
 
         const ThermoType& patchFaceMixture