From b8f338497f58cb796dc94159b4c3e68a9374d3a5 Mon Sep 17 00:00:00 2001
From: andy <andy>
Date: Thu, 5 Jul 2012 15:15:10 +0100
Subject: [PATCH] BUG: avoid calculation of rho field for every cell in PaSR
 combustion model

---
 src/combustionModels/PaSR/PaSR.C | 20 +++++++-------------
 1 file changed, 7 insertions(+), 13 deletions(-)

diff --git a/src/combustionModels/PaSR/PaSR.C b/src/combustionModels/PaSR/PaSR.C
index 286b9f3f..3e001612 100644
--- a/src/combustionModels/PaSR/PaSR.C
+++ b/src/combustionModels/PaSR/PaSR.C
@@ -94,35 +94,29 @@ void Foam::combustionModels::PaSR<CombThermoType>::correct()
 
         if (turbulentReaction_)
         {
+            tmp<volScalarField> trho(this->rho());
+            const volScalarField& rho = trho();
             tmp<volScalarField> tepsilon(this->turbulence().epsilon());
             const volScalarField& epsilon = tepsilon();
             tmp<volScalarField> tmuEff(this->turbulence().muEff());
             const volScalarField& muEff = tmuEff();
             tmp<volScalarField> ttc(tc());
             const volScalarField& tc = ttc();
+
+            const scalar dt = this->mesh().time().deltaTValue();
+
             forAll(epsilon, i)
             {
                 if (epsilon[i] > 0)
                 {
-                    const dimensionedScalar e0
-                    (
-                        "e0",
-                        sqr(dimLength)/pow3(dimTime), SMALL
-                    );
-
                     scalar tk =
                         Cmix_.value()
-                       *Foam::sqrt
-                       (
-                            muEff[i]/this->rho()()[i]/(epsilon[i] + e0.value())
-                       );
+                       *Foam::sqrt(muEff[i]/rho[i]/(epsilon[i] + SMALL));
 
                     // Chalmers PaSR model
                     if (!useReactionRate_)
                     {
-                        kappa_[i] =
-                            ( this->mesh().time().deltaTValue() + tc[i])
-                        /( this->mesh().time().deltaTValue() + tc[i] + tk);
+                        kappa_[i] = (dt + tc[i])/(dt + tc[i] + tk);
                     }
                     else
                     {