FreeFOAM The Cross-Platform CFD Toolkit
execFlowFunctionObjects.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 Global
25  execFlowFunctionObjects
26 
27 Application
28  execFlowFunctionObjects
29 
30 Description
31  Execute the set of functionObjects specified in the selected dictionary
32  (which defaults to system/controlDict) for the selected set of times.
33 
34  The flow (p-U) and optionally turbulence fields are available for the
35  function objects to operate on allowing forces and other related properties
36  to be calculated in addition to cutting planes etc.
37 
38 Usage
39 
40  - execFlowFunctionObjects [OPTIONS]
41 
42  @param -noWrite \n
43  Suppress output to files.
44 
45  @param -dict <dictionary name>\n
46  Use named dictionary instead of system/controlDict.
47 
48  @param -noZero \n
49  Ignore timestep 0.
50 
51  @param -constant \n
52  Include the constant directory.
53 
54  @param -time <time>\n
55  Apply only to specific time.
56 
57  @param -latestTime \n
58  Only apply to latest time step.
59 
60  @param -case <dir>\n
61  Case directory.
62 
63  @param -parallel \n
64  Run in parallel.
65 
66  @param -help \n
67  Display help message.
68 
69  @param -doc \n
70  Display Doxygen API documentation page for this application.
71 
72  @param -srcDoc \n
73  Display Doxygen source documentation page for this application.
74 
75 \*---------------------------------------------------------------------------*/
76 
77 #include <postCalc/calc.H>
78 
80 
83 
87 
88 
89 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
90 
91 namespace Foam
92 {
93  void execFlowFunctionObjects(const argList& args, const Time& runTime)
94  {
95  if (args.optionFound("dict"))
96  {
97  IOdictionary dict
98  (
99  IOobject
100  (
101  args.option("dict"),
102  runTime.system(),
103  runTime,
105  )
106  );
107 
108  functionObjectList fol(runTime, dict);
109  fol.start();
110  fol.execute();
111  }
112  else
113  {
114  functionObjectList fol(runTime);
115  fol.start();
116  fol.execute();
117  }
118  }
119 }
120 
121 
122 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
123 
124 void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh)
125 {
126  Info<< " Reading phi" << endl;
128  (
129  IOobject
130  (
131  "phi",
132  runTime.timeName(),
133  mesh,
134  IOobject::MUST_READ
135  ),
136  mesh
137  );
138 
139  Info<< " Reading U" << endl;
141  (
142  IOobject
143  (
144  "U",
145  runTime.timeName(),
146  mesh,
147  IOobject::MUST_READ
148  ),
149  mesh
150  );
151 
152  Info<< " Reading p" << endl;
154  (
155  IOobject
156  (
157  "p",
158  runTime.timeName(),
159  mesh,
160  IOobject::MUST_READ
161  ),
162  mesh
163  );
164 
165  if (phi.dimensions() == dimensionSet(0, 3, -1, 0, 0))
166  {
167  IOobject RASPropertiesHeader
168  (
169  "RASProperties",
170  runTime.constant(),
171  mesh,
172  IOobject::MUST_READ,
173  IOobject::NO_WRITE,
174  false
175  );
176 
177  IOobject LESPropertiesHeader
178  (
179  "LESProperties",
180  runTime.constant(),
181  mesh,
182  IOobject::MUST_READ,
183  IOobject::NO_WRITE,
184  false
185  );
186 
187  singlePhaseTransportModel laminarTransport(U, phi);
188 
189  if (RASPropertiesHeader.headerOk())
190  {
191  IOdictionary RASProperties(RASPropertiesHeader);
192 
193  autoPtr<incompressible::RASModel> RASModel
194  (
195  incompressible::RASModel::New
196  (
197  U,
198  phi,
200  )
201  );
202  execFlowFunctionObjects(args, runTime);
203  }
204  else if (LESPropertiesHeader.headerOk())
205  {
206  IOdictionary LESProperties(LESPropertiesHeader);
207 
208  autoPtr<incompressible::LESModel> sgsModel
209  (
210  incompressible::LESModel::New(U, phi, laminarTransport)
211  );
212 
213  execFlowFunctionObjects(args, runTime);
214  }
215  else
216  {
217  IOdictionary transportProperties
218  (
219  IOobject
220  (
221  "transportProperties",
222  runTime.constant(),
223  mesh,
224  IOobject::MUST_READ,
225  IOobject::NO_WRITE
226  )
227  );
228 
230 
231  execFlowFunctionObjects(args, runTime);
232  }
233  }
234  else if (phi.dimensions() == dimensionSet(1, 0, -1, 0, 0))
235  {
236  autoPtr<basicPsiThermo> thermo(basicPsiThermo::New(mesh));
237 
239  (
240  IOobject
241  (
242  "rho",
243  runTime.timeName(),
244  mesh
245  ),
246  thermo->rho()
247  );
248 
249  IOobject RASPropertiesHeader
250  (
251  "RASProperties",
252  runTime.constant(),
253  mesh,
254  IOobject::MUST_READ,
255  IOobject::NO_WRITE,
256  false
257  );
258 
259  IOobject LESPropertiesHeader
260  (
261  "LESProperties",
262  runTime.constant(),
263  mesh,
264  IOobject::MUST_READ,
265  IOobject::NO_WRITE,
266  false
267  );
268 
269  if (RASPropertiesHeader.headerOk())
270  {
271  IOdictionary RASProperties(RASPropertiesHeader);
272 
273  autoPtr<compressible::RASModel> RASModel
274  (
275  compressible::RASModel::New
276  (
277  rho,
278  U,
279  phi,
280  thermo()
281  )
282  );
283 
284  execFlowFunctionObjects(args, runTime);
285  }
286  else if (LESPropertiesHeader.headerOk())
287  {
288  IOdictionary LESProperties(LESPropertiesHeader);
289 
290  autoPtr<compressible::LESModel> sgsModel
291  (
292  compressible::LESModel::New(rho, U, phi, thermo())
293  );
294 
295  execFlowFunctionObjects(args, runTime);
296  }
297  else
298  {
299  IOdictionary transportProperties
300  (
301  IOobject
302  (
303  "transportProperties",
304  runTime.constant(),
305  mesh,
306  IOobject::MUST_READ,
307  IOobject::NO_WRITE
308  )
309  );
310 
312 
313  execFlowFunctionObjects(args, runTime);
314  }
315  }
316  else
317  {
318  FatalErrorIn(args.executable())
319  << "Incorrect dimensions of phi: " << phi.dimensions()
320  << nl << exit(FatalError);
321  }
322 }
323 
324 
325 // ************************ vim: set sw=4 sts=4 et: ************************ //