yPlus.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2013-2016 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "yPlus.H"
27 #include "turbulenceModel.H"
29 #include "wallFvPatch.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace functionObjects
37 {
38  defineTypeNameAndDebug(yPlus, 0);
39 
41  (
42  functionObject,
43  yPlus,
44  dictionary
45  );
46 }
47 }
48 
49 
50 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51 
52 void Foam::functionObjects::yPlus::writeFileHeader(const label i)
53 {
54  writeHeader(file(), "y+ ()");
55 
56  writeCommented(file(), "Time");
57  writeTabbed(file(), "patch");
58  writeTabbed(file(), "min");
59  writeTabbed(file(), "max");
60  writeTabbed(file(), "average");
61  file() << endl;
62 }
63 
64 
65 void Foam::functionObjects::yPlus::calcYPlus
66 (
67  const turbulenceModel& turbModel,
68  const fvMesh& mesh,
69  volScalarField& yPlus
70 )
71 {
72  volScalarField::Boundary d = nearWallDist(mesh).y();
73 
74  const volScalarField::Boundary nutBf =
75  turbModel.nut()().boundaryField();
76 
77  const volScalarField::Boundary nuEffBf =
78  turbModel.nuEff()().boundaryField();
79 
80  const volScalarField::Boundary nuBf =
81  turbModel.nu()().boundaryField();
82 
83  const fvPatchList& patches = mesh.boundary();
84 
85  volScalarField::Boundary& yPlusBf = yPlus.boundaryFieldRef();
86 
87  forAll(patches, patchi)
88  {
89  const fvPatch& patch = patches[patchi];
90 
91  if (isA<nutWallFunctionFvPatchScalarField>(nutBf[patchi]))
92  {
93  const nutWallFunctionFvPatchScalarField& nutPf =
94  dynamic_cast<const nutWallFunctionFvPatchScalarField&>
95  (
96  nutBf[patchi]
97  );
98 
99  yPlusBf[patchi] = nutPf.yPlus();
100  }
101  else if (isA<wallFvPatch>(patch))
102  {
103  yPlusBf[patchi] =
104  d[patchi]
105  *sqrt
106  (
107  nuEffBf[patchi]
108  *mag(turbModel.U().boundaryField()[patchi].snGrad())
109  )/nuBf[patchi];
110  }
111  }
112 }
113 
114 
115 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
116 
117 Foam::functionObjects::yPlus::yPlus
118 (
119  const word& name,
120  const Time& runTime,
121  const dictionary& dict
122 )
123 :
124  writeFiles(name, runTime, dict, name)
125 {
126  if (!isA<fvMesh>(obr_))
127  {
129  << "objectRegistry is not an fvMesh" << exit(FatalError);
130  }
131 
132  const fvMesh& mesh = refCast<const fvMesh>(obr_);
133 
134  volScalarField* yPlusPtr
135  (
136  new volScalarField
137  (
138  IOobject
139  (
140  type(),
141  mesh.time().timeName(),
142  mesh,
145  ),
146  mesh,
147  dimensionedScalar("0", dimless, 0.0)
148  )
149  );
150 
151  mesh.objectRegistry::store(yPlusPtr);
152 
153  resetName(typeName);
154 }
155 
156 
157 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
158 
160 {}
161 
162 
163 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
164 
166 {
167  writeFiles::read(dict);
168 
169  return true;
170 }
171 
172 
174 {
175  const fvMesh& mesh = refCast<const fvMesh>(obr_);
176 
177  volScalarField& yPlus =
178  const_cast<volScalarField&>
179  (
181  );
182 
184  {
185  const turbulenceModel& model =
187 
188  calcYPlus(model, mesh, yPlus);
189  }
190  else
191  {
193  << "Unable to find turbulence model in the "
194  << "database" << exit(FatalError);
195  }
196 
197  return true;
198 }
199 
200 
202 {
203  const volScalarField& yPlus =
204  obr_.lookupObject<volScalarField>(type());
205 
206  Log << type() << " " << name() << " write:" << nl
207  << " writing field " << yPlus.name() << endl;
208 
209  yPlus.write();
210 
212 
213  const volScalarField::Boundary& yPlusBf = yPlus.boundaryField();
214 
215  const fvMesh& mesh = refCast<const fvMesh>(obr_);
216  const fvPatchList& patches = mesh.boundary();
217 
218  forAll(patches, patchi)
219  {
220  const fvPatch& patch = patches[patchi];
221 
222  if (isA<wallFvPatch>(patch))
223  {
224  const scalarField& yPlusp = yPlusBf[patchi];
225 
226  const scalar minYplus = gMin(yPlusp);
227  const scalar maxYplus = gMax(yPlusp);
228  const scalar avgYplus = gAverage(yPlusp);
229 
230  if (Pstream::master())
231  {
232  Log << " patch " << patch.name()
233  << " y+ : min = " << minYplus << ", max = " << maxYplus
234  << ", average = " << avgYplus << nl;
235 
236  writeTime(file());
237  file()
238  << token::TAB << patch.name()
239  << token::TAB << minYplus
240  << token::TAB << maxYplus
241  << token::TAB << avgYplus
242  << endl;
243  }
244  }
245  }
246 
247  return true;
248 }
249 
250 
251 // ************************************************************************* //
virtual bool read(const dictionary &)
Read the yPlus data.
Definition: yPlus.C:165
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
virtual bool write()
Write the yPlus field.
Definition: yPlus.C:201
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
bool foundObject(const word &name) const
Is the named Type found?
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual bool execute()
Calculate the yPlus field.
Definition: yPlus.C:173
Type gMin(const FieldField< Field, Type > &f)
functionObject base class for writing files
Definition: writeFiles.H:56
const word & name() const
Return name.
Definition: fvPatch.H:149
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:411
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
Abstract base class for turbulence models (RAS, LES and laminar).
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
virtual bool read(const dictionary &)
Read optional controls.
dynamicFvMesh & mesh
const Boundary & boundaryField() const
Return const-reference to the boundary field.
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
static const char nl
Definition: Ostream.H:262
Type gMax(const FieldField< Field, Type > &f)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
PtrList< fvPatch > fvPatchList
container classes for fvPatch
Definition: fvPatchList.H:45
label patchi
virtual ~yPlus()
Destructor.
Definition: yPlus.C:159
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Type gAverage(const FieldField< Field, Type > &f)
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
addToRunTimeSelectionTable(functionObject, blendingFactor, dictionary)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual bool write()
Write function.
Definition: writeFiles.C:197
virtual bool write() const
Write using setting from DB.
#define Log
Report write to Foam::Info if the local log switch is true.
dimensioned< scalar > mag(const dimensioned< Type > &)
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:545
const word & name() const
Return name.
Definition: IOobject.H:260
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243