46 static const scalar Cmu(0.09);
47 static const scalar
kappa(0.41);
50 int main(
int argc,
char *argv[])
54 "apply a simplified boundary-layer model to the velocity and\n" 55 "turbulence fields based on the 1/7th power-law." 62 "specify the boundary-layer thickness" 68 "boundary-layer thickness as Cbl * mean distance to wall" 70 argList::addBoolOption
81 <<
"Neither option 'ybl' or 'Cbl' have been provided to calculate " 82 <<
"the boundary-layer thickness.\n" 83 <<
"Please choose either 'ybl' OR 'Cbl'." 89 <<
"Both 'ybl' and 'Cbl' have been provided to calculate " 90 <<
"the boundary-layer thickness.\n" 91 <<
"Please choose either 'ybl' OR 'Cbl'." 97 #include "createFields.H" 105 Info<<
"Setting boundary layer velocity" <<
nl <<
endl;
106 scalar yblv = ybl.value();
109 if (
y[celli] <= yblv)
112 U[celli] *=
::pow(
y[celli]/yblv, (1.0/7.0));
115 mask.correctBoundaryConditions();
121 #include "createPhi.H" 126 autoPtr<incompressible::turbulenceModel>
turbulence 131 if (isA<incompressible::RASModel>(
turbulence()))
156 k = (1 - mask)*k + mask*
sqr(nut/(ck0*
min(
y, ybl)));
167 tmp<volScalarField> tepsilon =
turbulence->epsilon();
170 epsilon = (1 - mask)*epsilon + mask*ce0*k*
sqrt(k)/
min(
y, ybl);
194 omega = (1 - mask)*omega + mask*epsilon/(Cmu*k + k0);
205 IOobject nuTildaHeader
228 Info<<
nl <<
"ExecutionTime = " << runTime.elapsedCpuTime() <<
" s" 229 <<
" ClockTime = " << runTime.elapsedClockTime() <<
" s"
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
#define forAll(list, i)
Loop across all elements in list.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar pow025(const dimensionedScalar &ds)
bool optionFound(const word &opt) const
Return true if the named option is found.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
singlePhaseTransportModel laminarTransport(U, phi)
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual Ostream & write(const token &)=0
Write next token to stream.
Foam::argList args(argc, argv)