diff --git a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
index 21ed8c53a8af1e8b6646348400f5de70f3a1d9aa..c00679e5be652b4f455d4078a321f9aca9ede983 100644
--- a/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
+++ b/applications/solvers/multiphase/reactingTwoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
@@ -68,16 +68,16 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::WenYu::CdRe() const
         max(scalar(1) - pair_.dispersed(), pair_.continuous().residualAlpha())
     );
 
-    volScalarField Re(pair_.Re());
-    volScalarField CdsRe
+    volScalarField Res(alpha2*pair_.Re());
+    volScalarField CdsRes
     (
-        neg(Re - 1000)*24.0*(1.0 + 0.15*pow(Re, 0.687))
-      + pos(Re - 1000)*0.44*max(Re, residualRe_)
+        neg(Res - 1000)*24.0*(1.0 + 0.15*pow(Res, 0.687))
+      + pos(Res - 1000)*0.44*max(Res, residualRe_)
     );
 
     return
-        CdsRe
-       *pow(alpha2, -2.65)
+        CdsRes
+       *pow(alpha2, -3.65)
        *max(pair_.continuous(), pair_.continuous().residualAlpha());
 }
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
index 21ed8c53a8af1e8b6646348400f5de70f3a1d9aa..c00679e5be652b4f455d4078a321f9aca9ede983 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
@@ -68,16 +68,16 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::WenYu::CdRe() const
         max(scalar(1) - pair_.dispersed(), pair_.continuous().residualAlpha())
     );
 
-    volScalarField Re(pair_.Re());
-    volScalarField CdsRe
+    volScalarField Res(alpha2*pair_.Re());
+    volScalarField CdsRes
     (
-        neg(Re - 1000)*24.0*(1.0 + 0.15*pow(Re, 0.687))
-      + pos(Re - 1000)*0.44*max(Re, residualRe_)
+        neg(Res - 1000)*24.0*(1.0 + 0.15*pow(Res, 0.687))
+      + pos(Res - 1000)*0.44*max(Res, residualRe_)
     );
 
     return
-        CdsRe
-       *pow(alpha2, -2.65)
+        CdsRes
+       *pow(alpha2, -3.65)
        *max(pair_.continuous(), pair_.continuous().residualAlpha());
 }