From e2c6985ea51298556927d90e11dc90bd75c9dcc1 Mon Sep 17 00:00:00 2001
From: Andrew Heather <a.heather@opencfd.co.uk>
Date: Thu, 17 May 2018 11:08:00 +0100
Subject: [PATCH] ENH: turbulenceFields FO - added turbulence intensity to list
 of available fields

---
 .../field/turbulenceFields/turbulenceFields.C | 14 ++++-
 .../field/turbulenceFields/turbulenceFields.H | 11 +++-
 .../turbulenceFieldsTemplates.C               | 51 ++++++++++++-------
 3 files changed, 54 insertions(+), 22 deletions(-)

diff --git a/src/functionObjects/field/turbulenceFields/turbulenceFields.C b/src/functionObjects/field/turbulenceFields/turbulenceFields.C
index 054ad322cb..0fc33a74e7 100644
--- a/src/functionObjects/field/turbulenceFields/turbulenceFields.C
+++ b/src/functionObjects/field/turbulenceFields/turbulenceFields.C
@@ -61,7 +61,8 @@ Foam::functionObjects::turbulenceFields::compressibleFieldNames_
     { compressibleField::cfAlphaEff, "alphaEff" },
     { compressibleField::cfR, "R" },
     { compressibleField::cfDevRhoReff, "devRhoReff" },
-    { compressibleField::cfL, "L" }
+    { compressibleField::cfL, "L" },
+    { compressibleField::cfI, "I" }
 };
 
 
@@ -80,6 +81,7 @@ Foam::functionObjects::turbulenceFields::incompressibleFieldNames_
     { incompressibleField::ifR, "R" },
     { incompressibleField::ifDevReff, "devReff" },
     { incompressibleField::ifL, "L" },
+    { incompressibleField::ifI, "I" }
 };
 
 
@@ -236,6 +238,11 @@ bool Foam::functionObjects::turbulenceFields::execute()
                     processField<scalar>(f, L(model));
                     break;
                 }
+                case cfI:
+                {
+                    processField<scalar>(f, I(model));
+                    break;
+                }
                 default:
                 {
                     FatalErrorInFunction
@@ -298,6 +305,11 @@ bool Foam::functionObjects::turbulenceFields::execute()
                     processField<scalar>(f, L(model));
                     break;
                 }
+                case ifI:
+                {
+                    processField<scalar>(f, I(model));
+                    break;
+                }
                 default:
                 {
                     FatalErrorInFunction
diff --git a/src/functionObjects/field/turbulenceFields/turbulenceFields.H b/src/functionObjects/field/turbulenceFields/turbulenceFields.H
index a4ae31d199..881a772f3b 100644
--- a/src/functionObjects/field/turbulenceFields/turbulenceFields.H
+++ b/src/functionObjects/field/turbulenceFields/turbulenceFields.H
@@ -76,6 +76,7 @@ Usage
         devReff     | Deviatoric part of the effective Reynolds stress
         devRhoReff  | Divergence of the Reynolds stress
         L           | turbulence length scale
+        I           | turbulence intensity
     \endplaintable
 
 See also
@@ -125,7 +126,8 @@ public:
         cfAlphaEff,
         cfR,
         cfDevRhoReff,
-        cfL
+        cfL,
+        cfI
     };
     static const Enum<compressibleField> compressibleFieldNames_;
 
@@ -139,7 +141,8 @@ public:
         ifNuEff,
         ifR,
         ifDevReff,
-        ifL
+        ifL,
+        ifI
     };
     static const Enum<incompressibleField> incompressibleFieldNames_;
 
@@ -179,6 +182,10 @@ protected:
         template<class Model>
         tmp<volScalarField> L(const Model& model) const;
 
+        //- Return I calculated from k and U
+        template<class Model>
+        tmp<volScalarField> I(const Model& model) const;
+
 
 private:
 
diff --git a/src/functionObjects/field/turbulenceFields/turbulenceFieldsTemplates.C b/src/functionObjects/field/turbulenceFields/turbulenceFieldsTemplates.C
index ec6b70d92d..4ded70c6be 100644
--- a/src/functionObjects/field/turbulenceFields/turbulenceFieldsTemplates.C
+++ b/src/functionObjects/field/turbulenceFields/turbulenceFieldsTemplates.C
@@ -85,19 +85,16 @@ Foam::functionObjects::turbulenceFields::omega
     const volScalarField k(model.k());
     const volScalarField epsilon(model.epsilon());
 
-    return tmp<volScalarField>
+    return tmp<volScalarField>::New
     (
-        new volScalarField
+        IOobject
         (
-            IOobject
-            (
-                "omega.tmp",
-                k.mesh().time().timeName(),
-                k.mesh()
-            ),
-            epsilon/(Cmu*k),
-            epsilon.boundaryField().types()
-        )
+            "omega.tmp",
+            k.mesh().time().timeName(),
+            k.mesh()
+        ),
+        epsilon/(Cmu*k),
+        epsilon.boundaryField().types()
     );
 }
 
@@ -109,9 +106,10 @@ Foam::functionObjects::turbulenceFields::nuTilda
     const Model& model
 ) const
 {
-    return tmp<volScalarField>
+    return tmp<volScalarField>::New
     (
-        new volScalarField("nuTilda.tmp", model.k()/omega(model))
+        "nuTilda.tmp",
+        model.k()/omega(model)
     );
 }
 
@@ -130,15 +128,30 @@ Foam::functionObjects::turbulenceFields::L
     const volScalarField epsilon(model.epsilon());
     const dimensionedScalar eps0("eps0", epsilon.dimensions(), SMALL);
 
-    return tmp<volScalarField>
+    return tmp<volScalarField>::New
     (
-        new volScalarField
-        (
-            "L.tmp",
-            pow(Cmu, 0.75)*pow(k, 1.5)/(epsilon + eps0)
-        )
+        "L.tmp",
+        pow(Cmu, 0.75)*pow(k, 1.5)/(epsilon + eps0)
     );
 }
 
 
+template<class Model>
+Foam::tmp<Foam::volScalarField>
+Foam::functionObjects::turbulenceFields::I
+(
+    const Model& model
+) const
+{
+    // Assume k is available
+    const volScalarField uPrime(sqrt((2.0/3.0)*model.k()));
+    const dimensionedScalar U0("U0", dimVelocity, SMALL);
+
+    return tmp<volScalarField>::New
+    (
+        "I.tmp",
+        uPrime/max(max(uPrime, mag(model.U())), U0)
+    );
+}
+
 // ************************************************************************* //