From 95b2a4d39bed3805700af1bf9c8f30d43292cdd8 Mon Sep 17 00:00:00 2001
From: andy <a.heather@opencfd.co.uk>
Date: Tue, 24 Nov 2009 13:05:01 +0000
Subject: [PATCH] re-design to use RASModel

---
 .../applyBoundaryLayer/Make/options           |   6 +
 .../applyBoundaryLayer/applyBoundaryLayer.C   | 166 +++++++-----------
 .../applyBoundaryLayer/createFields.H         |  78 ++++++++
 3 files changed, 152 insertions(+), 98 deletions(-)
 create mode 100644 applications/utilities/preProcessing/applyBoundaryLayer/createFields.H

diff --git a/applications/utilities/preProcessing/applyBoundaryLayer/Make/options b/applications/utilities/preProcessing/applyBoundaryLayer/Make/options
index 06d42b6c3c9..31f1dd4d95e 100644
--- a/applications/utilities/preProcessing/applyBoundaryLayer/Make/options
+++ b/applications/utilities/preProcessing/applyBoundaryLayer/Make/options
@@ -1,8 +1,14 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/turbulenceModels \
+    -I$(LIB_SRC)/turbulenceModels/incompressible/RAS/RASModel \
+    -I$(LIB_SRC)/transportModels \
+    -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
     -I$(LIB_SRC)/meshTools/lnInclude
 
 EXE_LIBS = \
     -lfiniteVolume \
+    -lincompressibleRASModels \
+    -lincompressibleTransportModels \
     -lgenericPatchFields \
     -lmeshTools
diff --git a/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C b/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
index 5714cccbffc..52d79cec815 100644
--- a/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
+++ b/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
@@ -37,6 +37,8 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
+#include "singlePhaseTransportModel.H"
+#include "RASModel.H"
 #include "wallDist.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -52,54 +54,14 @@ int main(int argc, char *argv[])
     argList::validOptions.insert("Cbl", "scalar");
     argList::validOptions.insert("writenut", "");
 
-#   include "setRootCase.H"
-#   include "createTime.H"
-#   include "createMesh.H"
+    #include "setRootCase.H"
+    #include "createTime.H"
+    #include "createMesh.H"
+    #include "createFields.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-    Info<< "Reading field U\n" << endl;
-    volVectorField U
-    (
-        IOobject
-        (
-            "U",
-            runTime.timeName(),
-            mesh,
-            IOobject::MUST_READ,
-            IOobject::NO_WRITE
-        ),
-        mesh
-    );
-
-#   include "createPhi.H"
-
-    Info<< "Calculating wall distance field" << endl;
-    volScalarField y = wallDist(mesh).y();
-
-    // Set the mean boundary-layer thickness
-    dimensionedScalar ybl("ybl", dimLength, 0);
-
-    if (args.optionFound("ybl"))
-    {
-        // If the boundary-layer thickness is provided use it
-        ybl.value() = args.optionRead<scalar>("ybl");
-    }
-    else if (args.optionFound("Cbl"))
-    {
-        // Calculate boundary layer thickness as Cbl * mean distance to wall
-        ybl.value() = gAverage(y) * args.optionRead<scalar>("Cbl");
-    }
-    else
-    {
-        FatalErrorIn(args.executable())
-            << "Neither option 'ybl' or 'Cbl' have been provided to calculate"
-               " the boundary-layer thickness"
-            << exit(FatalError);
-    }
+    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-    Info<< "\nCreating boundary-layer for U of thickness "
-        << ybl.value() << " m" << nl << endl;
+    Info<< "Time = " << runTime.timeName() << nl << endl;
 
     // Modify velocity by applying a 1/7th power law boundary-layer
     // u/U0 = (y/ybl)^(1/7)
@@ -114,58 +76,88 @@ int main(int argc, char *argv[])
         }
     }
 
-    Info<< "Writing U" << endl;
+    Info<< "Writing U\n" << endl;
     U.write();
 
     // Update/re-write phi
     phi = fvc::interpolate(U) & mesh.Sf();
     phi.write();
 
-    // Read and modify turbulence fields if present
+    // Calculate nut
+    tmp<volScalarField> tnut = turbulence->nut();
+    volScalarField& nut = tnut();
+    volScalarField S = mag(dev(symm(fvc::grad(U))));
+    nut = sqr(kappa*min(y, ybl))*::sqrt(2)*S;
 
-    IOobject epsilonHeader
-    (
-        "epsilon",
-        runTime.timeName(),
-        mesh,
-        IOobject::MUST_READ
-    );
+    if (args.optionFound("writenut"))
+    {
+        Info<< "Writing nut" << endl;
+        nut.write();
+    }
+
+    // Create G field - used by RAS wall functions
+    volScalarField G("RASModel::G", nut*2*sqr(S));
+
+
+    //--- Read and modify turbulence fields
+
+    // Turbulence k
+    tmp<volScalarField> tk = turbulence->k();
+    volScalarField& k = tk();
+    scalar ck0 = pow025(Cmu)*kappa;
+    k = sqr(nut/(ck0*min(y, ybl)));
+    k.correctBoundaryConditions();
+
+    Info<< "Writing k\n" << endl;
+    k.write();
 
-    IOobject kHeader
+
+    // Turbulence epsilon
+    tmp<volScalarField> tepsilon = turbulence->epsilon();
+    volScalarField& epsilon = tepsilon();
+    scalar ce0 = ::pow(Cmu, 0.75)/kappa;
+    epsilon = ce0*k*sqrt(k)/min(y, ybl);
+    epsilon.correctBoundaryConditions();
+
+    Info<< "Writing epsilon\n" << endl;
+    epsilon.write();
+
+
+    // Turbulence omega
+    IOobject omegaHeader
     (
-        "k",
+        "omega",
         runTime.timeName(),
         mesh,
-        IOobject::MUST_READ
+        IOobject::MUST_READ,
+        IOobject::NO_WRITE,
+        false
     );
 
+    if (omegaHeader.headerOk())
+    {
+        volScalarField omega(omegaHeader, mesh);
+        omega = epsilon/(Cmu*k);
+        omega.correctBoundaryConditions();
+
+        Info<< "Writing omega\n" << endl;
+        omega.write();
+    }
+
+
+    // Turbulence nuTilda
     IOobject nuTildaHeader
     (
         "nuTilda",
         runTime.timeName(),
         mesh,
-        IOobject::MUST_READ
+        IOobject::MUST_READ,
+        IOobject::NO_WRITE,
+        false
     );
 
-    // First calculate nut
-    volScalarField nut
-    (
-        "nut",
-        sqr(kappa*min(y, ybl))*::sqrt(2)*mag(dev(symm(fvc::grad(U))))
-    );
-
-    if (args.optionFound("writenut"))
-    {
-        Info<< "Writing nut" << endl;
-        nut.write();
-    }
-
-
-    // Read and modify turbulence fields if present
-
     if (nuTildaHeader.headerOk())
     {
-        Info<< "Reading field nuTilda\n" << endl;
         volScalarField nuTilda(nuTildaHeader, mesh);
         nuTilda = nut;
         nuTilda.correctBoundaryConditions();
@@ -174,28 +166,6 @@ int main(int argc, char *argv[])
         nuTilda.write();
     }
 
-    if (kHeader.headerOk() && epsilonHeader.headerOk())
-    {
-        Info<< "Reading field k\n" << endl;
-        volScalarField k(kHeader, mesh);
-
-        Info<< "Reading field epsilon\n" << endl;
-        volScalarField epsilon(epsilonHeader, mesh);
-
-        scalar ck0 = ::pow(Cmu, 0.25)*kappa;
-        k = sqr(nut/(ck0*min(y, ybl)));
-        k.correctBoundaryConditions();
-
-        scalar ce0 = ::pow(Cmu, 0.75)/kappa;
-        epsilon = ce0*k*sqrt(k)/min(y, ybl);
-        epsilon.correctBoundaryConditions();
-
-        Info<< "Writing k\n" << endl;
-        k.write();
-
-        Info<< "Writing epsilon\n" << endl;
-        epsilon.write();
-    }
 
     Info<< nl << "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
         << "  ClockTime = " << runTime.elapsedClockTime() << " s"
diff --git a/applications/utilities/preProcessing/applyBoundaryLayer/createFields.H b/applications/utilities/preProcessing/applyBoundaryLayer/createFields.H
new file mode 100644
index 00000000000..e168ed978d5
--- /dev/null
+++ b/applications/utilities/preProcessing/applyBoundaryLayer/createFields.H
@@ -0,0 +1,78 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+    Info<< "Reading field U\n" << endl;
+    volVectorField U
+    (
+        IOobject
+        (
+            "U",
+            runTime.timeName(),
+            mesh,
+            IOobject::MUST_READ,
+            IOobject::NO_WRITE
+        ),
+        mesh
+    );
+
+    #include "createPhi.H"
+
+    singlePhaseTransportModel laminarTransport(U, phi);
+
+    autoPtr<incompressible::RASModel> turbulence
+    (
+        incompressible::RASModel::New(U, phi, laminarTransport)
+    );
+
+    Info<< "Calculating wall distance field" << endl;
+    volScalarField y = wallDist(mesh).y();
+
+    // Set the mean boundary-layer thickness
+    dimensionedScalar ybl("ybl", dimLength, 0);
+
+    if (args.optionFound("ybl"))
+    {
+        // If the boundary-layer thickness is provided use it
+        ybl.value() = args.optionRead<scalar>("ybl");
+    }
+    else if (args.optionFound("Cbl"))
+    {
+        // Calculate boundary layer thickness as Cbl * mean distance to wall
+        ybl.value() = gAverage(y)*args.optionRead<scalar>("Cbl");
+    }
+    else
+    {
+        FatalErrorIn(args.executable())
+            << "Neither option 'ybl' or 'Cbl' have been provided to calculate "
+            << "the boundary-layer thickness"
+            << exit(FatalError);
+    }
+
+    Info<< "\nCreating boundary-layer for U of thickness "
+        << ybl.value() << " m" << nl << endl;
+
+
+// ************************************************************************* //
-- 
GitLab