59static const scalar Cmu(0.09);
60static const scalar
kappa(0.41);
64void correctProcessorPatches
66 GeometricField<Type, fvPatchField, volMesh>& vf
69 if (!Pstream::parRun())
77 auto& bf = vf.boundaryFieldRef();
81 if (isA<processorFvPatchField<Type>>(bf[patchi]))
83 bf[patchi].initEvaluate();
89 if (isA<processorFvPatchField<Type>>(bf[patchi]))
91 bf[patchi].evaluate();
99 const word& fieldName,
101 const scalarField& mask,
102 const scalarField& boundaryLayerField
108 mesh.time().timeName(),
115 if (fieldHeader.typeHeaderOk<volScalarField>(
true))
119 pf = (1 - mask)*pf + mask*boundaryLayerField;
126 correctProcessorPatches<scalar>(
fld);
137 const scalarField& mask,
138 const scalarField& kBL,
139 const scalarField& epsilonBL
146 mesh.time().timeName(),
153 if (omegaHeader.typeHeaderOk<volScalarField>(
true))
158 pf = (1 - mask)*pf + mask*epsilonBL/(Cmu*kBL + SMALL);
165 correctProcessorPatches<scalar>(omega);
176 const word& fieldName,
177 const volScalarField& value
183 mesh.time().timeName(),
190 if (fldHeader.typeHeaderOk<volScalarField>(
true))
199 correctProcessorPatches<scalar>(
fld);
207tmp<volScalarField> calcNut
210 const volVectorField&
U
219 basicThermo::dictName,
222 ).typeHeaderOk<IOdictionary>(
true)
226 autoPtr<fluidThermo>
pThermo(fluidThermo::New(
mesh));
234 autoPtr<compressible::turbulenceModel>
turbulence
236 compressible::turbulenceModel::New
248 return tmp<volScalarField>::New(
turbulence->nut());
255 #include "createPhi.H"
260 autoPtr<incompressible::turbulenceModel>
turbulence
268 return tmp<volScalarField>::New(
turbulence->nut());
273int main(
int argc,
char *argv[])
277 "Apply a simplified boundary-layer model to the velocity and"
278 " turbulence fields based on the 1/7th power-law."
287 "Specify the boundary-layer thickness"
293 "Boundary-layer thickness as Cbl * mean distance to wall"
295 argList::addBoolOption
297 "writeTurbulenceFields",
298 "Write the turbulence fields"
300 argList::addOptionCompat
302 "writeTurbulenceFields", {
"write-nut", 1906}
310 <<
"Neither option 'ybl' or 'Cbl' have been provided to calculate "
311 <<
"the boundary-layer thickness.\n"
312 <<
"Please choose either 'ybl' OR 'Cbl'."
318 <<
"Both 'ybl' and 'Cbl' have been provided to calculate "
319 <<
"the boundary-layer thickness.\n"
320 <<
"Please choose either 'ybl' OR 'Cbl'."
324 const bool writeTurbulenceFields =
args.
found(
"writeTurbulenceFields");
328 #include "createFields.H"
335 Info<<
"Setting boundary layer velocity" <<
nl <<
endl;
336 const scalar yblv = ybl.value();
339 if ((
y[celli] > 0) && (
y[celli] <= yblv))
342 U[celli] *= ::pow(
y[celli]/yblv, (1.0/7.0));
345 mask.correctBoundaryConditions();
346 correctProcessorPatches<vector>(
U);
348 if (writeTurbulenceFields)
355 nut = (1 - mask)*
nut + mask*
sqr(kappa*
min(
y, ybl))*::sqrt(2)*S;
359 correctProcessorPatches<scalar>(
nut);
369 scalar ce0 = ::pow(Cmu, 0.75)/
kappa;
373 blendField(
"k",
mesh, mask, kBL);
374 blendField(
"epsilon",
mesh, mask, epsilonBL);
375 calcOmegaField(
mesh, mask, kBL, epsilonBL);
384 runTime.printExecutionTime(Info);
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
bool found(const word &optName) const
Return true if the named option is found.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Creates and initialises the face-flux field phi.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
compressible::turbulenceModel & turbulence
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar pow025(const dimensionedScalar &ds)
constexpr char nl
The newline '\n' character (0x0a)
surfacesMesh setField(triSurfaceToAgglom)
Foam::argList args(argc, argv)
Info<< "Reading thermophysical properties\n"<< endl;autoPtr< psiReactionThermo > pThermo(psiReactionThermo::New(mesh))
singlePhaseTransportModel laminarTransport(U, phi)
#define forAll(list, i)
Loop across all elements in list.