47void Foam::LESModels::smoothDelta::setChangedFaces
50 const volScalarField&
delta,
51 DynamicList<label>& changedFaces,
52 DynamicList<deltaData>& changedFacesInfo
63 if (ownDelta > maxDeltaRatio_ * neiDelta)
65 changedFaces.append(facei);
66 changedFacesInfo.append(deltaData(ownDelta));
68 else if (neiDelta > maxDeltaRatio_ * ownDelta)
70 changedFaces.append(facei);
71 changedFacesInfo.append(deltaData(neiDelta));
85 label meshFacei =
patch.start() + patchFacei;
89 changedFaces.append(meshFacei);
90 changedFacesInfo.append(deltaData(ownDelta));
95 changedFaces.shrink();
96 changedFacesInfo.shrink();
100void Foam::LESModels::smoothDelta::calcDelta()
102 const fvMesh&
mesh = turbulenceModel_.mesh();
107 DynamicList<label> changedFaces(
mesh.nFaces()/100 + 100);
108 DynamicList<deltaData> changedFacesInfo(changedFaces.size());
110 setChangedFaces(
mesh, geometricDelta, changedFaces, changedFacesInfo);
113 List<deltaData> cellDeltaData(
mesh.nCells());
115 forAll(geometricDelta, celli)
117 cellDeltaData[celli] = geometricDelta[celli];
121 List<deltaData> faceDeltaData(
mesh.nFaces());
125 FaceCellWave<deltaData, scalar> deltaCalc
132 mesh.globalData().nTotalCells()+1,
138 delta_[celli] = cellDeltaData[celli].delta();
142 delta_.correctBoundaryConditions();
162 dict.optionalSubDict(
type() +
"Coeffs")
167 dict.optionalSubDict(
type() +
"Coeffs").get<scalar>(
"maxDeltaRatio")
180 geometricDelta_().read(coeffsDict);
181 coeffsDict.
readEntry(
"maxDeltaRatio", maxDeltaRatio_);
188 geometricDelta_().correct();
190 if (turbulenceModel_.mesh().changing())
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Abstract base class for LES deltas.
virtual bool read()
Re-read model coefficients if they have changed.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
virtual const labelList & faceOwner() const
Return face owner.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const labelList & faceNeighbour() const
Return face neighbour.
label nInternalFaces() const noexcept
Number of internal faces.
Abstract base class for turbulence models (RAS, LES and laminar).
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
compressible::turbulenceModel & turbulence
const std::string patch
OpenFOAM patch number as a std::string.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
#define forAll(list, i)
Loop across all elements in list.