From d7b7643ff4c394a2ea18dffcd7ffe0b6bfafe081 Mon Sep 17 00:00:00 2001
From: andy <andy>
Date: Tue, 12 Nov 2013 15:17:44 +0000
Subject: [PATCH] ENH: spray - parcel code clean-up

---
 .../Templates/SprayParcel/SprayParcel.C       | 37 +++++++------------
 .../BreakupModel/ReitzDiwakar/ReitzDiwakar.C  |  6 +--
 2 files changed, 14 insertions(+), 29 deletions(-)

diff --git a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
index 15d5cb1b622..16af6adc2af 100644
--- a/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
+++ b/src/lagrangian/spray/parcels/Templates/SprayParcel/SprayParcel.C
@@ -93,13 +93,11 @@ void Foam::SprayParcel<ParcelType>::calc
     td.cloud().constProps().setTMax(TMax);
 
     // store the parcel properties
-    const scalarField& Y(this->Y());
-    scalarField X(composition.liquids().X(Y));
-
-    this->Cp() = composition.liquids().Cp(this->pc_, T0, X);
-    sigma_ = composition.liquids().sigma(this->pc_, T0, X);
-    scalar rho0 = composition.liquids().rho(this->pc_, T0, X);
+    this->Cp() = composition.liquids().Cp(pc0, T0, X0);
+    sigma_ = composition.liquids().sigma(pc0, T0, X0);
+    scalar rho0 = composition.liquids().rho(pc0, T0, X0);
     this->rho() = rho0;
+    mu_ = composition.liquids().mu(pc0, T0, X0);
 
     ParcelType::calc(td, dt, cellI);
 
@@ -113,11 +111,13 @@ void Foam::SprayParcel<ParcelType>::calc
 
         this->Cp() = composition.liquids().Cp(this->pc_, T1, X1);
 
-        sigma_ = composition.liquids().sigma(this->pc_, T1, X);
+        sigma_ = composition.liquids().sigma(this->pc_, T1, X1);
 
         scalar rho1 = composition.liquids().rho(this->pc_, T1, X1);
         this->rho() = rho1;
 
+        mu_ = composition.liquids().mu(this->pc_, T1, X1);
+
         scalar d1 = this->d()*cbrt(rho0/rho1);
         this->d() = d1;
 
@@ -360,10 +360,6 @@ void Foam::SprayParcel<ParcelType>::solveTABEq
     const scalar dt
 )
 {
-    typedef typename TrackData::cloudType::reactingCloudType reactingCloudType;
-    const CompositionModel<reactingCloudType>& composition =
-        td.cloud().composition();
-
     const scalar& TABCmu = td.cloud().breakup().TABCmu();
     const scalar& TABWeCrit = td.cloud().breakup().TABWeCrit();
     const scalar& TABComega = td.cloud().breakup().TABComega();
@@ -372,26 +368,19 @@ void Foam::SprayParcel<ParcelType>::solveTABEq
     scalar r2 = r*r;
     scalar r3 = r*r2;
 
-    const scalarField& Y(this->Y());
-    scalarField X(composition.liquids().X(Y));
-
-    scalar rho = composition.liquids().rho(this->pc(), this->T(), X);
-    scalar mu = composition.liquids().mu(this->pc(), this->T(), X);
-    scalar sigma = composition.liquids().sigma(this->pc(), this->T(), X);
-
     // inverse of characteristic viscous damping time
-    scalar rtd = 0.5*TABCmu*mu/(rho*r2);
+    scalar rtd = 0.5*TABCmu*mu_/(this->rho()*r2);
 
     // oscillation frequency (squared)
-    scalar omega2 = TABComega*sigma/(rho*r3) - rtd*rtd;
+    scalar omega2 = TABComega*sigma_/(this->rho()*r3) - rtd*rtd;
 
     if (omega2 > 0)
     {
         scalar omega = sqrt(omega2);
         scalar rhoc = this->rhoc();
-        scalar Wetmp = this->We(this->U(), r, rhoc, sigma)/TABWeCrit;
+        scalar We = this->We(this->U(), r, rhoc, sigma_)/TABWeCrit;
 
-        scalar y1 = this->y() - Wetmp;
+        scalar y1 = this->y() - We;
         scalar y2 = this->yDot()/omega;
 
         // update distortion parameters
@@ -400,7 +389,7 @@ void Foam::SprayParcel<ParcelType>::solveTABEq
         scalar e = exp(-rtd*dt);
         y2 = (this->yDot() + y1*rtd)/omega;
 
-        this->y() = Wetmp + e*(y1*c + y2*s);
+        this->y() = We + e*(y1*c + y2*s);
         if (this->y() < 0)
         {
             this->y() = 0.0;
@@ -408,7 +397,7 @@ void Foam::SprayParcel<ParcelType>::solveTABEq
         }
         else
         {
-            this->yDot() = (Wetmp - this->y())*rtd + e*omega*(y2*c - y1*s);
+            this->yDot() = (We - this->y())*rtd + e*omega*(y2*c - y1*s);
         }
     }
     else
diff --git a/src/lagrangian/spray/submodels/BreakupModel/ReitzDiwakar/ReitzDiwakar.C b/src/lagrangian/spray/submodels/BreakupModel/ReitzDiwakar/ReitzDiwakar.C
index 79b26186eb1..504b681dc63 100644
--- a/src/lagrangian/spray/submodels/BreakupModel/ReitzDiwakar/ReitzDiwakar.C
+++ b/src/lagrangian/spray/submodels/BreakupModel/ReitzDiwakar/ReitzDiwakar.C
@@ -101,11 +101,9 @@ bool Foam::ReitzDiwakar<CloudType>::update
     scalar We = 0.5*rhoc*sqr(Urmag)*d/sigma;
     scalar Re = Urmag*d/nuc;
 
-    scalar sqRey = sqrt(Re);
-
     if (We > Cbag_)
     {
-        if (We > Cstrip_*sqRey)
+        if (We > Cstrip_*sqrt(Re))
         {
             scalar dStrip = sqr(2.0*Cstrip_*sigma)/(rhoc*pow3(Urmag)*muc);
             scalar tauStrip = Cs_*d*sqrt(rho/rhoc)/Urmag;
@@ -117,9 +115,7 @@ bool Foam::ReitzDiwakar<CloudType>::update
         else
         {
             scalar dBag = 2.0*Cbag_*sigma/(rhoc*sqr(Urmag));
-
             scalar tauBag = Cb_*d*sqrt(rho*d/sigma);
-
             scalar fraction = dt/tauBag;
 
             // new droplet diameter, implicit calculation
-- 
GitLab