FreeFOAM The Cross-Platform CFD Toolkit
applyWallFunctionBoundaryConditions.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-2007 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  applyWallFunctionBounaryConditions
26 
27 Description
28  Updates OpenFOAM RAS cases to use the new (v1.6) wall function framework
29 
30  Attempts to determine whether case is compressible or incompressible, or
31  can be supplied with -compressible command line argument
32 
33 Usage
34  - applyWallFunctionBoundaryConditions [OPTION]
35 
36  @param -compressible
37  Operate in compressible mode.
38 
39  @param -noZero \n
40  Ignore timestep 0.
41 
42  @param -constant \n
43  Include the constant directory.
44 
45  @param -time <time>\n
46  Apply only to specific time.
47 
48  @param -latestTime \n
49  Only apply to latest time step.
50 
51  @param -case <dir> \n
52  Specify the case directory
53 
54  @param -parallel \n
55  Run the case in parallel
56 
57  @param -help \n
58  Display short usage message
59 
60  @param -doc \n
61  Display Doxygen documentation page
62 
63  @param -srcDoc \n
64  Display source code
65 
66 \*---------------------------------------------------------------------------*/
67 
68 #include <OpenFOAM/argList.H>
69 #include <finiteVolume/fvMesh.H>
70 #include <OpenFOAM/Time.H>
71 #include <finiteVolume/volFields.H>
73 
74 #include <OpenFOAM/wallPolyPatch.H>
75 
80 
85 
86 using namespace Foam;
87 
88 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
89 
90 bool caseIsCompressible(const fvMesh& mesh)
91 {
92  // Attempt flux field
93  IOobject phiHeader
94  (
95  "phi",
96  mesh.time().timeName(),
97  mesh,
100  );
101 
102  if (phiHeader.headerOk())
103  {
104  surfaceScalarField phi(phiHeader, mesh);
105  if (phi.dimensions() == dimDensity*dimVelocity*dimArea)
106  {
107  return true;
108  }
109  }
110 
111  // Attempt density field
112  IOobject rhoHeader
113  (
114  "rho",
115  mesh.time().timeName(),
116  mesh,
119  );
120 
121  if (rhoHeader.headerOk())
122  {
123  volScalarField rho(rhoHeader, mesh);
124  if (rho.dimensions() == dimDensity)
125  {
126  return true;
127  }
128  }
129 
130  // Attempt pressure field
131  IOobject pHeader
132  (
133  "p",
134  mesh.time().timeName(),
135  mesh,
138  );
139 
140  if (pHeader.headerOk())
141  {
142  volScalarField p(pHeader, mesh);
143  if (p.dimensions() == dimMass/sqr(dimTime)/dimLength)
144  {
145  return true;
146  }
147  }
148 
149  // If none of the above are true, assume that the case is incompressible
150  return false;
151 }
152 
153 
154 void createVolScalarField
155 (
156  const fvMesh& mesh,
157  const word& fieldName,
158  const dimensionSet& dims
159 )
160 {
161  IOobject fieldHeader
162  (
163  fieldName,
164  mesh.time().timeName(),
165  mesh,
168  );
169 
170  if (!fieldHeader.headerOk())
171  {
172  Info<< "Creating field " << fieldName << nl << endl;
173 
174  volScalarField field
175  (
176  IOobject
177  (
178  fieldName,
179  mesh.time().timeName(),
180  mesh,
183  ),
184  mesh,
185  dimensionedScalar("zero", dims, 0.0)
186  );
187 
188  field.write();
189  }
190 }
191 
192 
193 void replaceBoundaryType
194 (
195  const fvMesh& mesh,
196  const word& fieldName,
197  const word& boundaryType,
198  const string& boundaryValue
199 )
200 {
201  IOobject header
202  (
203  fieldName,
204  mesh.time().timeName(),
205  mesh,
208  );
209 
210  if (!header.headerOk())
211  {
212  return;
213  }
214 
215  Info<< "Updating boundary types for field " << header.name() << endl;
216 
217  const word oldTypeName = IOdictionary::typeName;
218  const_cast<word&>(IOdictionary::typeName) = word::null;
219 
220  IOdictionary dict(header);
221 
222  const_cast<word&>(IOdictionary::typeName) = oldTypeName;
223  const_cast<word&>(dict.type()) = dict.headerClassName();
224 
225  // Make a backup of the old field
226  word backupName(dict.name() + ".old");
227  Info<< " copying " << dict.name() << " to "
228  << backupName << endl;
229  IOdictionary dictOld = dict;
230  dictOld.rename(backupName);
231  dictOld.regIOobject::write();
232 
233  // Loop through boundary patches and update
234  const polyBoundaryMesh& bMesh = mesh.boundaryMesh();
235  dictionary& boundaryDict = dict.subDict("boundaryField");
236  forAll(bMesh, patchI)
237  {
238  if (isA<wallPolyPatch>(bMesh[patchI]))
239  {
240  word patchName = bMesh[patchI].name();
241  dictionary& oldPatch = boundaryDict.subDict(patchName);
242 
243  dictionary newPatch(dictionary::null);
244  newPatch.add("type", boundaryType);
245  newPatch.add("value", ("uniform " + boundaryValue).c_str());
246 
247  oldPatch = newPatch;
248  }
249  }
250 
251  Info<< " writing updated " << dict.name() << nl << endl;
252  dict.regIOobject::write();
253 }
254 
255 
256 void updateCompressibleCase(const fvMesh& mesh)
257 {
258  Info<< "Case treated as compressible" << nl << endl;
259  createVolScalarField
260  (
261  mesh,
262  "mut",
264  );
265  replaceBoundaryType
266  (
267  mesh,
268  "mut",
269  compressible::RASModels::mutWallFunctionFvPatchScalarField::typeName,
270  "0"
271  );
272  replaceBoundaryType
273  (
274  mesh,
275  "epsilon",
276  compressible::RASModels::epsilonWallFunctionFvPatchScalarField::
277  typeName,
278  "0"
279  );
280  replaceBoundaryType
281  (
282  mesh,
283  "omega",
284  compressible::RASModels::omegaWallFunctionFvPatchScalarField::typeName,
285  "0"
286  );
287  replaceBoundaryType
288  (
289  mesh,
290  "k",
291  compressible::RASModels::kqRWallFunctionFvPatchField<scalar>::typeName,
292  "0"
293  );
294  replaceBoundaryType
295  (
296  mesh,
297  "q",
298  compressible::RASModels::kqRWallFunctionFvPatchField<scalar>::typeName,
299  "0"
300  );
301  replaceBoundaryType
302  (
303  mesh,
304  "R",
305  compressible::RASModels::kqRWallFunctionFvPatchField<symmTensor>::
306  typeName,
307  "(0 0 0 0 0 0)"
308  );
309 }
310 
311 
312 void updateIncompressibleCase(const fvMesh& mesh)
313 {
314  Info<< "Case treated as incompressible" << nl << endl;
315  createVolScalarField(mesh, "nut", dimArea/dimTime);
316 
317  replaceBoundaryType
318  (
319  mesh,
320  "nut",
321  incompressible::RASModels::nutWallFunctionFvPatchScalarField::typeName,
322  "0"
323  );
324  replaceBoundaryType
325  (
326  mesh,
327  "epsilon",
329  typeName,
330  "0"
331  );
332  replaceBoundaryType
333  (
334  mesh,
335  "omega",
337  typeName,
338  "0"
339  );
340  replaceBoundaryType
341  (
342  mesh,
343  "k",
345  typeName,
346  "0"
347  );
348  replaceBoundaryType
349  (
350  mesh,
351  "q",
353  typeName,
354  "0"
355  );
356  replaceBoundaryType
357  (
358  mesh,
359  "R",
361  typeName,
362  "(0 0 0 0 0 0)"
363  );
364 }
365 
366 
367 int main(int argc, char *argv[])
368 {
369  #include <OpenFOAM/addTimeOptions.H>
370  argList::validOptions.insert("compressible", "");
371 
372  #include <OpenFOAM/setRootCase.H>
373  #include <OpenFOAM/createTime.H>
374  #include <OpenFOAM/createMesh.H>
375 
376  bool compressible = args.optionFound("compressible");
377 
378  Info<< "Updating turbulence fields to operate using new run time "
379  << "selectable" << nl << "wall functions"
380  << nl << endl;
381 
382  if (compressible || caseIsCompressible(mesh))
383  {
384  updateCompressibleCase(mesh);
385  }
386  else
387  {
388  updateIncompressibleCase(mesh);
389  }
390 
391  Info<< "End\n" << endl;
392 
393  return 0;
394 }
395 
396 
397 // ************************ vim: set sw=4 sts=4 et: ************************ //