diff --git a/applications/utilities/preProcessing/applyBoundaryLayer/Make/options b/applications/utilities/preProcessing/applyBoundaryLayer/Make/options
index 9778b15c09da197c25fef02b86cd13e6b2dcecde..6b1672b25db2f2cdb2abc7d42ef7d556564c4ba4 100644
--- a/applications/utilities/preProcessing/applyBoundaryLayer/Make/options
+++ b/applications/utilities/preProcessing/applyBoundaryLayer/Make/options
@@ -1,7 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
+ -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
+ -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/transportModels \
+ -I$(LIB_SRC)/transportModels/compressible/lnInclude \
-I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
@@ -9,7 +12,9 @@ EXE_INC = \
EXE_LIBS = \
-lturbulenceModels \
-lincompressibleTurbulenceModels \
+ -lcompressibleTurbulenceModels \
-lincompressibleTransportModels \
+ -lcompressibleTransportModels \
-lgenericPatchFields \
-lfiniteVolume \
-lmeshTools
diff --git a/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C b/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
index 4f2a30fb5f40539e52e5a4298213dd60168ad95c..e00f12b48ebfccf6f14fbf25ef276358deaee5ae 100644
--- a/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
+++ b/applications/utilities/preProcessing/applyBoundaryLayer/applyBoundaryLayer.C
@@ -38,7 +38,9 @@ Description
#include "fvCFD.H"
#include "singlePhaseTransportModel.H"
#include "turbulentTransportModel.H"
+#include "turbulentFluidThermoModel.H"
#include "wallDist.H"
+#include "processorFvPatchField.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -46,6 +48,266 @@ Description
static const scalar Cmu(0.09);
static const scalar kappa(0.41);
+void correctProcessorPatches(volScalarField& vf)
+{
+ if (!Pstream::parRun())
+ {
+ return;
+ }
+
+ // Not possible to use correctBoundaryConditions on fields as they may
+ // use local info as opposed to the constraint values employed here,
+ // but still need to update processor patches
+
+ volScalarField::GeometricBoundaryField& bf = vf.boundaryField();
+
+ forAll(bf, patchI)
+ {
+ if (isA<processorFvPatchField<scalar> >(bf[patchI]))
+ {
+ bf[patchI].initEvaluate();
+ }
+ }
+
+ forAll(bf, patchI)
+ {
+ if (isA<processorFvPatchField<scalar> >(bf[patchI]))
+ {
+ bf[patchI].evaluate();
+ }
+ }
+}
+
+
+template<class TurbulenceModel>
+Foam::tmp<Foam::volScalarField> calcK
+(
+ TurbulenceModel& turbulence,
+ const volScalarField& mask,
+ const volScalarField& nut,
+ const volScalarField& y,
+ const dimensionedScalar& ybl,
+ const scalar Cmu,
+ const scalar kappa
+)
+{
+ // Turbulence k
+ tmp<volScalarField> tk = turbulence->k();
+ volScalarField& k = tk();
+ scalar ck0 = pow025(Cmu)*kappa;
+ k = (1 - mask)*k + mask*sqr(nut/(ck0*min(y, ybl)));
+
+ // Do not correct BC
+ // - operation may use inconsistent fields wrt these local manipulations
+ // k.correctBoundaryConditions();
+ correctProcessorPatches(k);
+
+ Info<< "Writing k\n" << endl;
+ k.write();
+
+ return tk;
+}
+
+
+template<class TurbulenceModel>
+Foam::tmp<Foam::volScalarField> calcEpsilon
+(
+ TurbulenceModel& turbulence,
+ const volScalarField& mask,
+ const volScalarField& k,
+ const volScalarField& y,
+ const dimensionedScalar& ybl,
+ const scalar Cmu,
+ const scalar kappa
+)
+{
+ // Turbulence epsilon
+ tmp<volScalarField> tepsilon = turbulence->epsilon();
+ volScalarField& epsilon = tepsilon();
+ scalar ce0 = ::pow(Cmu, 0.75)/kappa;
+ epsilon = (1 - mask)*epsilon + mask*ce0*k*sqrt(k)/min(y, ybl);
+ epsilon.max(SMALL);
+
+ // Do not correct BC
+ // - operation may use inconsistent fields wrt these local manipulations
+ // epsilon.correctBoundaryConditions();
+ correctProcessorPatches(epsilon);
+
+ Info<< "Writing epsilon\n" << endl;
+ epsilon.write();
+
+ return tepsilon;
+}
+
+
+void calcOmega
+(
+ const fvMesh& mesh,
+ const volScalarField& mask,
+ const volScalarField& k,
+ const volScalarField& epsilon
+)
+{
+ // Turbulence omega
+ IOobject omegaHeader
+ (
+ "omega",
+ mesh.time().timeName(),
+ mesh,
+ IOobject::MUST_READ,
+ IOobject::NO_WRITE,
+ false
+ );
+
+ if (omegaHeader.headerOk())
+ {
+ volScalarField omega(omegaHeader, mesh);
+ dimensionedScalar k0("SMALL", k.dimensions(), SMALL);
+
+ omega = (1 - mask)*omega + mask*epsilon/(Cmu*k + k0);
+ omega.max(SMALL);
+
+ // Do not correct BC
+ // - operation may use inconsistent fields wrt these local
+ // manipulations
+ // omega.correctBoundaryConditions();
+ correctProcessorPatches(omega);
+
+ Info<< "Writing omega\n" << endl;
+ omega.write();
+ }
+}
+
+
+void setField
+(
+ const fvMesh& mesh,
+ const word& fieldName,
+ const volScalarField& value
+)
+{
+ IOobject fldHeader
+ (
+ fieldName,
+ mesh.time().timeName(),
+ mesh,
+ IOobject::MUST_READ,
+ IOobject::NO_WRITE,
+ false
+ );
+
+ if (fldHeader.headerOk())
+ {
+ volScalarField fld(fldHeader, mesh);
+ fld = value;
+
+ // Do not correct BC
+ // - operation may use inconsistent fields wrt these local
+ // manipulations
+ // fld.correctBoundaryConditions();
+ correctProcessorPatches(fld);
+
+ Info<< "Writing " << fieldName << nl << endl;
+ fld.write();
+ }
+}
+
+
+void calcCompressible
+(
+ const fvMesh& mesh,
+ const volScalarField& mask,
+ const volVectorField& U,
+ const volScalarField& y,
+ const dimensionedScalar& ybl
+)
+{
+ const Time& runTime = mesh.time();
+
+ autoPtr<fluidThermo> pThermo(fluidThermo::New(mesh));
+ fluidThermo& thermo = pThermo();
+ volScalarField rho(thermo.rho());
+
+ // Update/re-write phi
+ #include "compressibleCreatePhi.H"
+ phi.write();
+
+ autoPtr<compressible::turbulenceModel> turbulence
+ (
+ compressible::turbulenceModel::New
+ (
+ rho,
+ U,
+ phi,
+ thermo
+ )
+ );
+
+ // Calculate nut - reference nut is calculated by the turbulence model
+ // on its construction
+ tmp<volScalarField> tnut = turbulence->nut();
+
+ volScalarField& nut = tnut();
+ volScalarField S(mag(dev(symm(fvc::grad(U)))));
+ nut = (1 - mask)*nut + mask*sqr(kappa*min(y, ybl))*::sqrt(2)*S;
+
+ // Do not correct BC - wall functions will 'undo' manipulation above
+ // by using nut from turbulence model
+ correctProcessorPatches(nut);
+ nut.write();
+
+ tmp<volScalarField> k =
+ calcK(turbulence, mask, nut, y, ybl, Cmu, kappa);
+ tmp<volScalarField> epsilon =
+ calcEpsilon(turbulence, mask, k, y, ybl, Cmu, kappa);
+ calcOmega(mesh, mask, k, epsilon);
+ setField(mesh, "nuTilda", nut);
+}
+
+
+void calcIncompressible
+(
+ const fvMesh& mesh,
+ const volScalarField& mask,
+ const volVectorField& U,
+ const volScalarField& y,
+ const dimensionedScalar& ybl
+)
+{
+ const Time& runTime = mesh.time();
+
+ // Update/re-write phi
+ #include "createPhi.H"
+ phi.write();
+
+ singlePhaseTransportModel laminarTransport(U, phi);
+
+ autoPtr<incompressible::turbulenceModel> turbulence
+ (
+ incompressible::turbulenceModel::New(U, phi, laminarTransport)
+ );
+
+ tmp<volScalarField> tnut = turbulence->nut();
+
+ // Calculate nut - reference nut is calculated by the turbulence model
+ // on its construction
+ volScalarField& nut = tnut();
+ volScalarField S(mag(dev(symm(fvc::grad(U)))));
+ nut = (1 - mask)*nut + mask*sqr(kappa*min(y, ybl))*::sqrt(2)*S;
+
+ // Do not correct BC - wall functions will 'undo' manipulation above
+ // by using nut from turbulence model
+ correctProcessorPatches(nut);
+ nut.write();
+
+ tmp<volScalarField> k =
+ calcK(turbulence, mask, nut, y, ybl, Cmu, kappa);
+ tmp<volScalarField> epsilon =
+ calcEpsilon(turbulence, mask, k, y, ybl, Cmu, kappa);
+ calcOmega(mesh, mask, k, epsilon);
+ setField(mesh, "nuTilda", nut);
+}
+
int main(int argc, char *argv[])
{
@@ -55,6 +317,8 @@ int main(int argc, char *argv[])
"turbulence fields based on the 1/7th power-law."
);
+ #include "addRegionOption.H"
+
argList::addOption
(
"ybl",
@@ -69,8 +333,8 @@ int main(int argc, char *argv[])
);
argList::addBoolOption
(
- "writenut",
- "write nut field"
+ "compressible",
+ "apply to compressible case"
);
#include "setRootCase.H"
@@ -93,9 +357,11 @@ int main(int argc, char *argv[])
}
#include "createTime.H"
- #include "createMesh.H"
+ #include "createNamedMesh.H"
#include "createFields.H"
+ const bool compressible = args.optionFound("compressible");
+
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Modify velocity by applying a 1/7th power law boundary-layer
@@ -117,113 +383,16 @@ int main(int argc, char *argv[])
Info<< "Writing U\n" << endl;
U.write();
- // Update/re-write phi
- #include "createPhi.H"
- phi.write();
-
- singlePhaseTransportModel laminarTransport(U, phi);
-
- autoPtr<incompressible::turbulenceModel> turbulence
- (
- incompressible::turbulenceModel::New(U, phi, laminarTransport)
- );
- if (isA<incompressible::RASModel>(turbulence()))
+ if (compressible)
{
- // Calculate nut - reference nut is calculated by the turbulence model
- // on its construction
- tmp<volScalarField> tnut = turbulence->nut();
- volScalarField& nut = tnut();
- volScalarField S(mag(dev(symm(fvc::grad(U)))));
- nut = (1 - mask)*nut + mask*sqr(kappa*min(y, ybl))*::sqrt(2)*S;
-
- // do not correct BC - wall functions will 'undo' manipulation above
- // by using nut from turbulence model
-
- if (args.optionFound("writenut"))
- {
- Info<< "Writing nut" << endl;
- nut.write();
- }
-
-
- //--- Read and modify turbulence fields
-
- // Turbulence k
- tmp<volScalarField> tk = turbulence->k();
- volScalarField& k = tk();
- scalar ck0 = pow025(Cmu)*kappa;
- k = (1 - mask)*k + mask*sqr(nut/(ck0*min(y, ybl)));
-
- // do not correct BC - operation may use inconsistent fields wrt these
- // local manipulations
- // k.correctBoundaryConditions();
-
- Info<< "Writing k\n" << endl;
- k.write();
-
-
- // Turbulence epsilon
- tmp<volScalarField> tepsilon = turbulence->epsilon();
- volScalarField& epsilon = tepsilon();
- scalar ce0 = ::pow(Cmu, 0.75)/kappa;
- epsilon = (1 - mask)*epsilon + mask*ce0*k*sqrt(k)/min(y, ybl);
-
- // do not correct BC - wall functions will use non-updated k from
- // turbulence model
- // epsilon.correctBoundaryConditions();
-
- Info<< "Writing epsilon\n" << endl;
- epsilon.write();
-
- // Turbulence omega
- IOobject omegaHeader
- (
- "omega",
- runTime.timeName(),
- mesh,
- IOobject::MUST_READ,
- IOobject::NO_WRITE,
- false
- );
-
- if (omegaHeader.headerOk())
- {
- volScalarField omega(omegaHeader, mesh);
- dimensionedScalar k0("VSMALL", k.dimensions(), VSMALL);
- omega = (1 - mask)*omega + mask*epsilon/(Cmu*k + k0);
-
- // do not correct BC - wall functions will use non-updated k from
- // turbulence model
- // omega.correctBoundaryConditions();
-
- Info<< "Writing omega\n" << endl;
- omega.write();
- }
-
- // Turbulence nuTilda
- IOobject nuTildaHeader
- (
- "nuTilda",
- runTime.timeName(),
- mesh,
- IOobject::MUST_READ,
- IOobject::NO_WRITE,
- false
- );
-
- if (nuTildaHeader.headerOk())
- {
- volScalarField nuTilda(nuTildaHeader, mesh);
- nuTilda = nut;
-
- // do not correct BC
- // nuTilda.correctBoundaryConditions();
-
- Info<< "Writing nuTilda\n" << endl;
- nuTilda.write();
- }
+ calcCompressible(mesh, mask, U, y, ybl);
}
+ else
+ {
+ calcIncompressible(mesh, mask, U, y, ybl);
+ }
+
Info<< nl << "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"