diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseIncompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseIncompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C
index bb2658a3ef..1fc59d91fb 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseIncompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseIncompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.C
@@ -55,8 +55,6 @@ Foam::RASModels::kineticTheoryModel::kineticTheoryModel
 
     phase_(phase),
 
-    draga_(phase.fluid().drag1()),
-
     viscosityModel_
     (
         kineticTheoryModels::viscosityModel::New
@@ -401,7 +399,11 @@ void Foam::RASModels::kineticTheoryModel::correct()
             (
                 alpha*(1.0 - alpha),
                 phase_.fluid().residualPhaseFraction()
-            )*draga_.K(magUr + phase_.fluid().residualSlip())/rho
+            )
+           *phase_.fluid().drag(phase_).K
+            (
+                magUr + phase_.fluid().residualSlip()
+            )/rho
         );
 
         // Eq. 3.25, p. 50 Js = J1 - J2
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseIncompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.H b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseIncompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.H
index 776b523f56..ecf03b105e 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/phaseIncompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/phaseIncompressibleTurbulenceModels/kineticTheoryModels/kineticTheoryModel/kineticTheoryModel.H
@@ -82,8 +82,8 @@ class kineticTheoryModel
             const phaseModel& phase_;
 
 
-        //- Drag model
-        const dragModel& draga_;
+        ////- Drag model
+        //const dragModel& draga_;
 
         // Sub-models
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.H b/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.H
index 967ffc2226..cc781bcbc1 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseSystem/twoPhaseSystem.H
@@ -82,7 +82,7 @@ class twoPhaseSystem
         //- Total volumetric flux
         surfaceScalarField phi_;
 
-        //- 
+        //-  Dilatation term
         volScalarField dgdt_;
 
         //- Surface tension coefficient
@@ -224,7 +224,7 @@ public:
                 return phi_;
             }
 
-            //- Return
+            //- Return the dilatation term
             const volScalarField& dgdt() const
             {
                 return dgdt_;