FreeFOAM The Cross-Platform CFD Toolkit
reconstructPar.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  reconstructPar
26 
27 Description
28  Reconstructs a mesh and fields of a case that is decomposed for parallel
29  execution of OpenFOAM.
30 
31 Usage
32 
33  - reconstructPar [OPTIONS]
34 
35  @param -region <name>\n
36  Only apply to named mesh region.
37 
38  @param -fields (field names)
39  Only apply to named fields.
40 
41  @param -noLagrangian
42  Exclude Lagrangian fields.
43 
44  @param -noZero \n
45  Ignore timestep 0.
46 
47  @param -constant \n
48  Include the constant directory.
49 
50  @param -time <time>\n
51  Apply only to specific time.
52 
53  @param -latestTime \n
54  Only apply to latest time step.
55 
56  @param -case <dir>\n
57  Case directory.
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 <OpenFOAM/argList.H>
71 #include <OpenFOAM/timeSelector.H>
72 
73 #include <finiteVolume/fvCFD.H>
74 #include <OpenFOAM/IOobjectList.H>
75 #include "processorMeshes.H"
76 #include "fvFieldReconstructor.H"
78 #include "reconstructLagrangian.H"
79 
80 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
81 
82 int main(int argc, char *argv[])
83 {
84  // enable -constant ... if someone really wants it
85  // enable -zeroTime to prevent accidentally trashing the initial fields
86  timeSelector::addOptions(true, true);
87  argList::noParallel();
89  argList::validOptions.insert("fields", "\"(list of fields)\"");
90  argList::validOptions.insert("noLagrangian", "");
91 
92 # include <OpenFOAM/setRootCase.H>
93 # include <OpenFOAM/createTime.H>
94 
95  HashSet<word> selectedFields;
96  if (args.optionFound("fields"))
97  {
98  args.optionLookup("fields")() >> selectedFields;
99  }
100 
101  bool noLagrangian = args.optionFound("noLagrangian");
102 
103  // determine the processor count directly
104  label nProcs = 0;
105  while (isDir(args.path()/(word("processor") + name(nProcs))))
106  {
107  ++nProcs;
108  }
109 
110  if (!nProcs)
111  {
113  << "No processor* directories found"
114  << exit(FatalError);
115  }
116 
117  // Create the processor databases
118  PtrList<Time> databases(nProcs);
119 
120  forAll (databases, procI)
121  {
122  databases.set
123  (
124  procI,
125  new Time
126  (
127  Time::controlDictName,
128  args.rootPath(),
129  args.caseName()/fileName(word("processor") + name(procI))
130  )
131  );
132  }
133 
134  // use the times list from the master processor
135  // and select a subset based on the command-line options
136  instantList timeDirs = timeSelector::select
137  (
138  databases[0].times(),
139  args
140  );
141 
142  if (timeDirs.empty())
143  {
145  << "No times selected"
146  << exit(FatalError);
147  }
148 
149 # include <OpenFOAM/createNamedMesh.H>
150  fileName regionPrefix = "";
151  if (regionName != fvMesh::defaultRegion)
152  {
153  regionPrefix = regionName;
154  }
155 
156  // Set all times on processor meshes equal to reconstructed mesh
157  forAll (databases, procI)
158  {
159  databases[procI].setTime(runTime.timeName(), runTime.timeIndex());
160  }
161 
162  // Read all meshes and addressing to reconstructed mesh
163  processorMeshes procMeshes(databases, regionName);
164 
165 
166  // check face addressing for meshes that have been decomposed
167  // with a very old foam version
168 # include "checkFaceAddressingComp.H"
169 
170  // Loop over all times
171  forAll (timeDirs, timeI)
172  {
173  // Set time for global database
174  runTime.setTime(timeDirs[timeI], timeI);
175 
176  Info << "Time = " << runTime.timeName() << endl << endl;
177 
178  // Set time for all databases
179  forAll (databases, procI)
180  {
181  databases[procI].setTime(timeDirs[timeI], timeI);
182  }
183 
184  // Check if any new meshes need to be read.
185  fvMesh::readUpdateState meshStat = mesh.readUpdate();
186 
187  fvMesh::readUpdateState procStat = procMeshes.readUpdate();
188 
189  if (procStat == fvMesh::POINTS_MOVED)
190  {
191  // Reconstruct the points for moving mesh cases and write them out
192  procMeshes.reconstructPoints(mesh);
193  }
194  else if (meshStat != procStat)
195  {
197  << "readUpdate for the reconstructed mesh:" << meshStat << nl
198  << "readUpdate for the processor meshes :" << procStat << nl
199  << "These should be equal or your addressing"
200  << " might be incorrect."
201  << " Please check your time directories for any "
202  << "mesh directories." << endl;
203  }
204 
205 
206  // Get list of objects from processor0 database
207  IOobjectList objects(procMeshes.meshes()[0], databases[0].timeName());
208 
209 
210  // If there are any FV fields, reconstruct them
211 
212  if
213  (
214  objects.lookupClass(volScalarField::typeName).size()
215  || objects.lookupClass(volVectorField::typeName).size()
216  || objects.lookupClass(volSphericalTensorField::typeName).size()
217  || objects.lookupClass(volSymmTensorField::typeName).size()
218  || objects.lookupClass(volTensorField::typeName).size()
219  || objects.lookupClass(surfaceScalarField::typeName).size()
220  || objects.lookupClass(surfaceVectorField::typeName).size()
221  || objects.lookupClass(surfaceSphericalTensorField::typeName).size()
222  || objects.lookupClass(surfaceSymmTensorField::typeName).size()
223  || objects.lookupClass(surfaceTensorField::typeName).size()
224  )
225  {
226  Info << "Reconstructing FV fields" << nl << endl;
227 
228  fvFieldReconstructor fvReconstructor
229  (
230  mesh,
231  procMeshes.meshes(),
232  procMeshes.faceProcAddressing(),
233  procMeshes.cellProcAddressing(),
234  procMeshes.boundaryProcAddressing()
235  );
236 
237  fvReconstructor.reconstructFvVolumeFields<scalar>
238  (
239  objects,
240  selectedFields
241  );
242  fvReconstructor.reconstructFvVolumeFields<vector>
243  (
244  objects,
245  selectedFields
246  );
247  fvReconstructor.reconstructFvVolumeFields<sphericalTensor>
248  (
249  objects,
250  selectedFields
251  );
252  fvReconstructor.reconstructFvVolumeFields<symmTensor>
253  (
254  objects,
255  selectedFields
256  );
257  fvReconstructor.reconstructFvVolumeFields<tensor>
258  (
259  objects,
260  selectedFields
261  );
262 
263  fvReconstructor.reconstructFvSurfaceFields<scalar>
264  (
265  objects,
266  selectedFields
267  );
268  fvReconstructor.reconstructFvSurfaceFields<vector>
269  (
270  objects,
271  selectedFields
272  );
273  fvReconstructor.reconstructFvSurfaceFields<sphericalTensor>
274  (
275  objects,
276  selectedFields
277  );
278  fvReconstructor.reconstructFvSurfaceFields<symmTensor>
279  (
280  objects,
281  selectedFields
282  );
283  fvReconstructor.reconstructFvSurfaceFields<tensor>
284  (
285  objects,
286  selectedFields
287  );
288  }
289  else
290  {
291  Info << "No FV fields" << nl << endl;
292  }
293 
294 
295  // If there are any point fields, reconstruct them
296  if
297  (
298  objects.lookupClass(pointScalarField::typeName).size()
299  || objects.lookupClass(pointVectorField::typeName).size()
300  || objects.lookupClass(pointSphericalTensorField::typeName).size()
301  || objects.lookupClass(pointSymmTensorField::typeName).size()
302  || objects.lookupClass(pointTensorField::typeName).size()
303  )
304  {
305  Info << "Reconstructing point fields" << nl << endl;
306 
307  pointMesh pMesh(mesh);
308  PtrList<pointMesh> pMeshes(procMeshes.meshes().size());
309 
310  forAll (pMeshes, procI)
311  {
312  pMeshes.set(procI, new pointMesh(procMeshes.meshes()[procI]));
313  }
314 
315  pointFieldReconstructor pointReconstructor
316  (
317  pMesh,
318  pMeshes,
319  procMeshes.pointProcAddressing(),
320  procMeshes.boundaryProcAddressing()
321  );
322 
323  pointReconstructor.reconstructFields<scalar>(objects);
324  pointReconstructor.reconstructFields<vector>(objects);
325  pointReconstructor.reconstructFields<sphericalTensor>(objects);
326  pointReconstructor.reconstructFields<symmTensor>(objects);
327  pointReconstructor.reconstructFields<tensor>(objects);
328  }
329  else
330  {
331  Info << "No point fields" << nl << endl;
332  }
333 
334 
335  // If there are any clouds, reconstruct them.
336  // The problem is that a cloud of size zero will not get written so
337  // in pass 1 we determine the cloud names and per cloud name the
338  // fields. Note that the fields are stored as IOobjectList from
339  // the first processor that has them. They are in pass2 only used
340  // for name and type (scalar, vector etc).
341 
342  if (!noLagrangian)
343  {
344  HashTable<IOobjectList> cloudObjects;
345 
346  forAll (databases, procI)
347  {
348  fileNameList cloudDirs
349  (
350  readDir
351  (
352  databases[procI].timePath()/regionPrefix/cloud::prefix,
353  fileName::DIRECTORY
354  )
355  );
356 
357  forAll (cloudDirs, i)
358  {
359  // Check if we already have cloud objects for this cloudname
361  cloudObjects.find(cloudDirs[i]);
362 
363  if (iter == cloudObjects.end())
364  {
365  // Do local scan for valid cloud objects
366  IOobjectList sprayObjs
367  (
368  procMeshes.meshes()[procI],
369  databases[procI].timeName(),
370  cloud::prefix/cloudDirs[i]
371  );
372 
373  IOobject* positionsPtr = sprayObjs.lookup("positions");
374 
375  if (positionsPtr)
376  {
377  cloudObjects.insert(cloudDirs[i], sprayObjs);
378  }
379  }
380  }
381  }
382 
383 
384  if (cloudObjects.size())
385  {
386  // Pass2: reconstruct the cloud
387  forAllConstIter(HashTable<IOobjectList>, cloudObjects, iter)
388  {
389  const word cloudName = string::validate<word>(iter.key());
390 
391  // Objects (on arbitrary processor)
392  const IOobjectList& sprayObjs = iter();
393 
394  Info<< "Reconstructing lagrangian fields for cloud "
395  << cloudName << nl << endl;
396 
398  (
399  mesh,
400  cloudName,
401  procMeshes.meshes(),
402  procMeshes.faceProcAddressing(),
403  procMeshes.cellProcAddressing()
404  );
405  reconstructLagrangianFields<label>
406  (
407  cloudName,
408  mesh,
409  procMeshes.meshes(),
410  sprayObjs
411  );
412  reconstructLagrangianFields<scalar>
413  (
414  cloudName,
415  mesh,
416  procMeshes.meshes(),
417  sprayObjs
418  );
419  reconstructLagrangianFields<vector>
420  (
421  cloudName,
422  mesh,
423  procMeshes.meshes(),
424  sprayObjs
425  );
426  reconstructLagrangianFields<sphericalTensor>
427  (
428  cloudName,
429  mesh,
430  procMeshes.meshes(),
431  sprayObjs
432  );
433  reconstructLagrangianFields<symmTensor>
434  (
435  cloudName,
436  mesh,
437  procMeshes.meshes(),
438  sprayObjs
439  );
440  reconstructLagrangianFields<tensor>
441  (
442  cloudName,
443  mesh,
444  procMeshes.meshes(),
445  sprayObjs
446  );
447  }
448  }
449  else
450  {
451  Info << "No lagrangian fields" << nl << endl;
452  }
453  }
454 
455  // If there are any "uniform" directories copy them from
456  // the master processor
457 
458  fileName uniformDir0 = databases[0].timePath()/"uniform";
459  if (isDir(uniformDir0))
460  {
461  cp(uniformDir0, runTime.timePath());
462  }
463  }
464 
465  Info<< "End.\n" << endl;
466 
467  return 0;
468 }
469 
470 
471 // ************************ vim: set sw=4 sts=4 et: ************************ //