35 Foam::epsilonWallFunctionFvPatchScalarField::calculate
42 const nutWallFunctionFvPatchScalarField& nutw =
47 const tmp<scalarField> tnuw = mtm.nu(
patchi);
50 const tmp<volScalarField> tk = mtm.k();
57 const scalar Cmu25 =
pow025(nutw.Cmu());
58 const scalar Cmu75 =
pow(nutw.Cmu(), 0.75);
61 Pair<tmp<scalarField>> GandEpsilon
74 const scalar
yPlus = Cmu25*
y[facei]*
sqrt(
k[celli])/nuw[facei];
76 if (
yPlus > nutw.yPlusLam())
79 (nutw[facei] + nuw[facei])
82 /(nutw.kappa()*
y[facei]);
85 Cmu75*
k[celli]*
sqrt(
k[celli])/(nutw.kappa()*
y[facei]);
99 void Foam::epsilonWallFunctionFvPatchScalarField::updateCoeffsMaster()
101 if (patch().index() != masterPatchIndex())
111 VolInternalField<scalar>&
G =
112 db().lookupObjectRef<VolInternalField<scalar>>(mtm.GName());
128 PtrList<scalarField> Gpfs(epsilonBf.size());
129 PtrList<scalarField> epsilonPfs(epsilonBf.size());
132 if (isA<epsilonWallFunctionFvPatchScalarField>(epsilonBf[
patchi]))
134 const epsilonWallFunctionFvPatchScalarField& epsilonPf =
135 refCast<const epsilonWallFunctionFvPatchScalarField>
138 Pair<tmp<scalarField>> GandEpsilon = epsilonPf.calculate(mtm);
140 Gpfs.set(
patchi, GandEpsilon.first().ptr());
141 epsilonPfs.set(
patchi, GandEpsilon.second().ptr());
146 wallCellGPtr_.reset(patchFieldsToWallCellField(Gpfs).ptr());
147 wallCellEpsilonPtr_.reset(patchFieldsToWallCellField(epsilonPfs).ptr());
150 UIndirectList<scalar>(
G, wallCells()) =
152 + wallCellFraction()*wallCellGPtr_();
153 UIndirectList<scalar>(
epsilon, wallCells()) =
155 + wallCellFraction()*wallCellEpsilonPtr_();
159 void Foam::epsilonWallFunctionFvPatchScalarField::manipulateMatrixMaster
161 fvMatrix<scalar>& matrix
164 if (patch().index() != masterPatchIndex())
169 matrix.setValues(wallCells(), wallCellEpsilonPtr_(), wallCellFraction());
184 wallCellGPtr_(nullptr),
185 wallCellEpsilonPtr_(nullptr)
199 wallCellGPtr_(nullptr),
200 wallCellEpsilonPtr_(nullptr)
212 wallCellGPtr_(nullptr),
213 wallCellEpsilonPtr_(nullptr)
226 wallCellGPtr_.clear();
227 wallCellEpsilonPtr_.clear();
237 wallCellGPtr_.clear();
238 wallCellEpsilonPtr_.clear();
251 updateCoeffsMaster();
264 if (manipulatedMatrix())
269 if (masterPatchIndex() == -1)
272 <<
"updateCoeffs must be called before manipulateMatrix"
276 manipulateMatrixMaster(matrix);
#define forAll(list, i)
Loop across all elements in list.
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
GeometricBoundaryField< Type, PatchField, GeoMesh > Boundary
Type of the boundary field.
label size() const
Return the number of elements in the UList.
An ordered pair of two objects of type <T> with first() and second() elements.
A list of keyword definitions, which are a keyword followed by any number of values (e....
This boundary condition provides a turbulence dissipation wall constraint for low- and high-Reynolds ...
virtual void reset(const fvPatchScalarField &)
Reset the fvPatchField to the given fvPatchField.
virtual void manipulateMatrix(fvMatrix< scalar > &matrix)
Manipulate matrix.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void map(const fvPatchScalarField &, const fieldMapper &)
Map the given fvPatchField onto this fvPatchField.
epsilonWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
Abstract base class for field mapping.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
const fvPatch & patch() const
Return patch.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
label index() const
Return the index of this patch in the fvBoundaryMesh.
virtual const labelUList & faceCells() const
Return faceCells.
static const nutWallFunctionFvPatchScalarField & nutw(const momentumTransportModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
Base class for wall functions that modify cell values.
virtual void reset(const fvPatchScalarField &)
Reset the fvPatchField to the given fvPatchField.
virtual void map(const fvPatchScalarField &, const fieldMapper &)
Map the given fvPatchField onto this fvPatchField.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
compressibleMomentumTransportModel momentumTransportModel
const dimensionedScalar G
Newtonian constant of gravitation.
errorManipArg< error, int > exit(error &err, const int errNo=1)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
VolField< scalar > volScalarField
dimensioned< scalar > mag(const dimensioned< Type > &)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
fvPatchField< vector > fvPatchVectorField
dimensionedScalar pow025(const dimensionedScalar &ds)