FreeFOAM The Cross-Platform CFD Toolkit
yPlusRAS.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) 1991-2010 OpenCFD Ltd.
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 Application
25  yPlusRAS
26 
27 Description
28  Calculates and reports yPlus for all wall patches, for the specified times
29  when using RAS turbulence models.
30 
31  Default behaviour assumes operating in incompressible mode. To apply to
32  compressible RAS cases, use the -compressible option.
33 
34 Usage
35 
36  - yPlusRAS [OPTIONS]
37 
38  @param -compressible
39  Operate in compressible mode.
40 
41  @param -noZero \n
42  Ignore timestep 0.
43 
44  @param -constant \n
45  Include the constant directory.
46 
47  @param -time <time>\n
48  Apply only to specific time.
49 
50  @param -latestTime \n
51  Only apply to latest time step.
52 
53  @param -case <dir>\n
54  Case directory.
55 
56  @param -parallel \n
57  Run in parallel.
58 
59  @param -help \n
60  Display help message.
61 
62  @param -doc \n
63  Display Doxygen API documentation page for this application.
64 
65  @param -srcDoc \n
66  Display Doxygen source documentation page for this application.
67 
68 \*---------------------------------------------------------------------------*/
69 
70 #include <finiteVolume/fvCFD.H>
71 
75 
79 
80 #include <finiteVolume/wallDist.H>
81 
82 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
83 
84 void calcIncompressibleYPlus
85 (
86  const fvMesh& mesh,
87  const Time& runTime,
88  const volVectorField& U,
89  volScalarField& yPlus
90 )
91 {
93  wallFunctionPatchField;
94 
95  #include <finiteVolume/createPhi.H>
96 
98 
100  (
101  incompressible::RASModel::New(U, phi, laminarTransport)
102  );
103 
104  const volScalarField::GeometricBoundaryField nutPatches =
105  RASModel->nut()().boundaryField();
106 
107  bool foundNutPatch = false;
108  forAll(nutPatches, patchi)
109  {
110  if (isA<wallFunctionPatchField>(nutPatches[patchi]))
111  {
112  foundNutPatch = true;
113 
114  const wallFunctionPatchField& nutPw =
115  dynamic_cast<const wallFunctionPatchField&>
116  (nutPatches[patchi]);
117 
118  yPlus.boundaryField()[patchi] = nutPw.yPlus();
119  const scalarField& Yp = yPlus.boundaryField()[patchi];
120 
121  Info<< "Patch " << patchi
122  << " named " << nutPw.patch().name()
123  << " y+ : min: " << min(Yp) << " max: " << max(Yp)
124  << " average: " << average(Yp) << nl << endl;
125  }
126  }
127 
128  if (!foundNutPatch)
129  {
130  Info<< " no " << wallFunctionPatchField::typeName << " patches"
131  << endl;
132  }
133 }
134 
135 
136 void calcCompressibleYPlus
137 (
138  const fvMesh& mesh,
139  const Time& runTime,
140  const volVectorField& U,
141  volScalarField& yPlus
142 )
143 {
144  typedef compressible::RASModels::mutWallFunctionFvPatchScalarField
145  wallFunctionPatchField;
146 
147  IOobject rhoHeader
148  (
149  "rho",
150  runTime.timeName(),
151  mesh,
152  IOobject::MUST_READ,
153  IOobject::NO_WRITE
154  );
155 
156  if (!rhoHeader.headerOk())
157  {
158  Info<< " no rho field" << endl;
159  return;
160  }
161 
162  Info << "Reading field rho\n" << endl;
163  volScalarField rho(rhoHeader, mesh);
164 
166 
168  (
169  basicPsiThermo::New(mesh)
170  );
172 
174  (
175  compressible::RASModel::New
176  (
177  rho,
178  U,
179  phi,
180  thermo
181  )
182  );
183 
184  const volScalarField::GeometricBoundaryField mutPatches =
185  RASModel->mut()().boundaryField();
186 
187  bool foundMutPatch = false;
188  forAll(mutPatches, patchi)
189  {
190  if (isA<wallFunctionPatchField>(mutPatches[patchi]))
191  {
192  foundMutPatch = true;
193 
194  const wallFunctionPatchField& mutPw =
195  dynamic_cast<const wallFunctionPatchField&>
196  (mutPatches[patchi]);
197 
198  yPlus.boundaryField()[patchi] = mutPw.yPlus();
199  const scalarField& Yp = yPlus.boundaryField()[patchi];
200 
201  Info<< "Patch " << patchi
202  << " named " << mutPw.patch().name()
203  << " y+ : min: " << min(Yp) << " max: " << max(Yp)
204  << " average: " << average(Yp) << nl << endl;
205  }
206  }
207 
208  if (!foundMutPatch)
209  {
210  Info<< " no " << wallFunctionPatchField::typeName << " patches"
211  << endl;
212  }
213 }
214 
215 
216 int main(int argc, char *argv[])
217 {
218  timeSelector::addOptions();
219 
220  #include <OpenFOAM/addRegionOption.H>
221 
222  argList::validOptions.insert("compressible","");
223 
224  #include <OpenFOAM/setRootCase.H>
225  #include <OpenFOAM/createTime.H>
226  instantList timeDirs = timeSelector::select0(runTime, args);
227  #include <OpenFOAM/createNamedMesh.H>
228 
229  bool compressible = args.optionFound("compressible");
230 
231  forAll(timeDirs, timeI)
232  {
233  runTime.setTime(timeDirs[timeI], timeI);
234  Info<< "Time = " << runTime.timeName() << endl;
235  fvMesh::readUpdateState state = mesh.readUpdate();
236 
237  // Wall distance
238  if (timeI == 0 || state != fvMesh::UNCHANGED)
239  {
240  Info<< "Calculating wall distance\n" << endl;
241  wallDist y(mesh, true);
242  Info<< "Writing wall distance to field "
243  << y.name() << nl << endl;
244  y.write();
245  }
246 
247  volScalarField yPlus
248  (
249  IOobject
250  (
251  "yPlus",
252  runTime.timeName(),
253  mesh,
254  IOobject::NO_READ,
255  IOobject::NO_WRITE
256  ),
257  mesh,
258  dimensionedScalar("yPlus", dimless, 0.0)
259  );
260 
261  IOobject UHeader
262  (
263  "U",
264  runTime.timeName(),
265  mesh,
266  IOobject::MUST_READ,
267  IOobject::NO_WRITE
268  );
269 
270  if (UHeader.headerOk())
271  {
272  Info << "Reading field U\n" << endl;
273  volVectorField U(UHeader, mesh);
274 
275  if (compressible)
276  {
277  calcCompressibleYPlus(mesh, runTime, U, yPlus);
278  }
279  else
280  {
281  calcIncompressibleYPlus(mesh, runTime, U, yPlus);
282  }
283  }
284  else
285  {
286  Info<< " no U field" << endl;
287  }
288 
289  Info<< "Writing yPlus to field " << yPlus.name() << nl << endl;
290 
291  yPlus.write();
292  }
293 
294  Info<< "End\n" << endl;
295 
296  return 0;
297 }
298 
299 // ************************ vim: set sw=4 sts=4 et: ************************ //