loadOrCreateMesh.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2012-2017 OpenFOAM Foundation
9  Copyright (C) 2015-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "loadOrCreateMesh.H"
30 #include "processorPolyPatch.H"
32 #include "Time.H"
33 //#include "IOPtrList.H"
35 #include "IOobjectList.H"
36 #include "pointSet.H"
37 #include "faceSet.H"
38 #include "cellSet.H"
39 #include "basicFvGeometryScheme.H"
40 
41 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
42 
43 //namespace Foam
44 //{
45 // defineTemplateTypeNameAndDebug(IOPtrList<entry>, 0);
46 //}
47 
48 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
49 
50 // Read mesh if available. Otherwise create empty mesh with same non-proc
51 // patches as proc0 mesh. Requires all processors to have all patches
52 // (and in same order).
54 (
55  const IOobject& io
56 )
57 {
58  // Region name
59  // ~~~~~~~~~~~
60 
61  fileName meshSubDir;
62 
63  if (io.name() == polyMesh::defaultRegion)
64  {
65  meshSubDir = polyMesh::meshSubDir;
66  }
67  else
68  {
69  meshSubDir = io.name()/polyMesh::meshSubDir;
70  }
71 
72 
73  // Patch types
74  // ~~~~~~~~~~~
75  // Read and scatter master patches (without reading master mesh!)
76 
77  PtrList<entry> patchEntries;
78  if (Pstream::master())
79  {
81  //const word oldTypeName = IOPtrList<entry>::typeName;
82  //const_cast<word&>(IOPtrList<entry>::typeName) = word::null;
83  //IOPtrList<entry> dictList
84  //(
85  // IOobject
86  // (
87  // "boundary",
88  // io.time().findInstance
89  // (
90  // meshSubDir,
91  // "boundary",
92  // IOobject::MUST_READ
93  // ),
94  // meshSubDir,
95  // io.db(),
96  // IOobject::MUST_READ,
97  // IOobject::NO_WRITE,
98  // false
99  // )
100  //);
101  //const_cast<word&>(IOPtrList<entry>::typeName) = oldTypeName;
103  //const_cast<word&>(dictList.type()) = dictList.headerClassName();
104  //
105  //patchEntries.transfer(dictList);
106  const fileName facesInstance = io.time().findInstance
107  (
108  meshSubDir,
109  "faces",
110  IOobject::MUST_READ
111  );
112 
113  patchEntries = polyBoundaryMeshEntries
114  (
115  IOobject
116  (
117  "boundary",
118  facesInstance,
119  meshSubDir,
120  io.db(),
121  IOobject::MUST_READ,
122  IOobject::NO_WRITE,
123  false
124  )
125  );
126 
127  // Send patches
128  for (const int slave : Pstream::subProcs())
129  {
130  OPstream toSlave(Pstream::commsTypes::scheduled, slave);
131  toSlave << patchEntries;
132  }
133  }
134  else
135  {
136  // Receive patches
137  IPstream fromMaster
138  (
139  Pstream::commsTypes::scheduled,
140  Pstream::masterNo()
141  );
142  fromMaster >> patchEntries;
143  }
144 
145 
146 
147  // Dummy meshes
148  // ~~~~~~~~~~~~
149 
150  // Check who has a mesh
151  const bool haveMesh = fileHandler().isFile
152  (
153  fileHandler().filePath
154  (
155  io.time().path()/io.instance()/meshSubDir/"faces"
156  )
157  );
158 
159  if (!haveMesh)
160  {
161  const bool oldParRun = Pstream::parRun(false);
162 
163 
164  // Create dummy mesh. Only used on procs that don't have mesh.
165  IOobject noReadIO(io);
166  noReadIO.readOpt() = IOobject::NO_READ;
167  noReadIO.writeOpt() = IOobject::AUTO_WRITE;
168  fvMesh dummyMesh(noReadIO, Zero, false);
169 
170  // Add patches
171  List<polyPatch*> patches(patchEntries.size());
172  label nPatches = 0;
173 
174  forAll(patchEntries, patchi)
175  {
176  const entry& e = patchEntries[patchi];
177  const word type(e.dict().get<word>("type"));
178  const word& name = e.keyword();
179 
180  if
181  (
182  type != processorPolyPatch::typeName
183  && type != processorCyclicPolyPatch::typeName
184  )
185  {
186  dictionary patchDict(e.dict());
187  patchDict.set("nFaces", 0);
188  patchDict.set("startFace", 0);
189 
190  patches[patchi] = polyPatch::New
191  (
192  name,
193  patchDict,
194  nPatches++,
195  dummyMesh.boundaryMesh()
196  ).ptr();
197  }
198  }
199  patches.setSize(nPatches);
200  dummyMesh.addFvPatches(patches, false); // no parallel comms
201 
202 
203  // Add some dummy zones so upon reading it does not read them
204  // from the undecomposed case. Should be done as extra argument to
205  // regIOobject::readStream?
206  List<pointZone*> pz
207  (
208  1,
209  new pointZone
210  (
211  "dummyPointZone",
212  0,
213  dummyMesh.pointZones()
214  )
215  );
216  List<faceZone*> fz
217  (
218  1,
219  new faceZone
220  (
221  "dummyFaceZone",
222  0,
223  dummyMesh.faceZones()
224  )
225  );
226  List<cellZone*> cz
227  (
228  1,
229  new cellZone
230  (
231  "dummyCellZone",
232  0,
233  dummyMesh.cellZones()
234  )
235  );
236  dummyMesh.addZones(pz, fz, cz);
237  dummyMesh.pointZones().clear();
238  dummyMesh.faceZones().clear();
239  dummyMesh.cellZones().clear();
240  //Pout<< "Writing dummy mesh to " << dummyMesh.polyMesh::objectPath()
241  // << endl;
242  dummyMesh.write();
243 
244  Pstream::parRun(oldParRun); // Restore parallel state
245  }
246 
247 
248 
249  // Read mesh
250  // ~~~~~~~~~
251  // Now all processors have a (possibly zero size) mesh so read in
252  // parallel
253 
254  //Pout<< "Reading mesh from " << io.objectPath() << endl;
255  auto meshPtr = autoPtr<fvMesh>::New(io);
256  fvMesh& mesh = *meshPtr;
257 
258  // Make sure to use a non-parallel geometry calculation method
259  {
260  tmp<fvGeometryScheme> basicGeometry
261  (
263  (
264  mesh,
265  dictionary(),
266  basicFvGeometryScheme::typeName
267  )
268  );
269  mesh.geometry(basicGeometry);
270  }
271 
272 
273  // Sync patches
274  // ~~~~~~~~~~~~
275 
276  if (!Pstream::master() && haveMesh)
277  {
278  // Check master names against mine
279 
280  const polyBoundaryMesh& patches = mesh.boundaryMesh();
281 
282  forAll(patchEntries, patchi)
283  {
284  const entry& e = patchEntries[patchi];
285  const word type(e.dict().get<word>("type"));
286  const word& name = e.keyword();
287 
288  if (type == processorPolyPatch::typeName)
289  {
290  break;
291  }
292 
293  if (patchi >= patches.size())
294  {
296  << "Non-processor patches not synchronised."
297  << endl
298  << "Processor " << Pstream::myProcNo()
299  << " has only " << patches.size()
300  << " patches, master has "
301  << patchi
302  << exit(FatalError);
303  }
304 
305  if
306  (
307  type != patches[patchi].type()
308  || name != patches[patchi].name()
309  )
310  {
312  << "Non-processor patches not synchronised."
313  << endl
314  << "Master patch " << patchi
315  << " name:" << type
316  << " type:" << type << endl
317  << "Processor " << Pstream::myProcNo()
318  << " patch " << patchi
319  << " has name:" << patches[patchi].name()
320  << " type:" << patches[patchi].type()
321  << exit(FatalError);
322  }
323  }
324  }
325 
326 
327  // Determine zones
328  // ~~~~~~~~~~~~~~~
329 
330  wordList pointZoneNames(mesh.pointZones().names());
331  Pstream::scatter(pointZoneNames);
332  wordList faceZoneNames(mesh.faceZones().names());
333  Pstream::scatter(faceZoneNames);
334  wordList cellZoneNames(mesh.cellZones().names());
335  Pstream::scatter(cellZoneNames);
336 
337  if (!haveMesh)
338  {
339  // Add the zones. Make sure to remove the old dummy ones first
340  mesh.pointZones().clear();
341  mesh.faceZones().clear();
342  mesh.cellZones().clear();
343 
344  List<pointZone*> pz(pointZoneNames.size());
345  forAll(pointZoneNames, i)
346  {
347  pz[i] = new pointZone
348  (
349  pointZoneNames[i],
350  i,
351  mesh.pointZones()
352  );
353  }
354  List<faceZone*> fz(faceZoneNames.size());
355  forAll(faceZoneNames, i)
356  {
357  fz[i] = new faceZone
358  (
359  faceZoneNames[i],
360  i,
361  mesh.faceZones()
362  );
363  }
364  List<cellZone*> cz(cellZoneNames.size());
365  forAll(cellZoneNames, i)
366  {
367  cz[i] = new cellZone
368  (
369  cellZoneNames[i],
370  i,
371  mesh.cellZones()
372  );
373  }
374  mesh.addZones(pz, fz, cz);
375  }
376 
377 
378  // Determine sets
379  // ~~~~~~~~~~~~~~
380 
381  wordList pointSetNames;
382  wordList faceSetNames;
383  wordList cellSetNames;
384  if (Pstream::master())
385  {
386  // Read sets
387  IOobjectList objects(mesh, mesh.facesInstance(), "polyMesh/sets");
388  pointSetNames = objects.sortedNames(pointSet::typeName);
389  faceSetNames = objects.sortedNames(faceSet::typeName);
390  cellSetNames = objects.sortedNames(cellSet::typeName);
391  }
392  Pstream::scatter(pointSetNames);
393  Pstream::scatter(faceSetNames);
394  Pstream::scatter(cellSetNames);
395 
396  if (!haveMesh)
397  {
398  for (const word& setName : pointSetNames)
399  {
400  pointSet(mesh, setName, 0).write();
401  }
402  for (const word& setName : faceSetNames)
403  {
404  faceSet(mesh, setName, 0).write();
405  }
406  for (const word& setName : cellSetNames)
407  {
408  cellSet(mesh, setName, 0).write();
409  }
410  }
411 
412 
413 // if (!haveMesh)
414 // {
415 // // We created a dummy mesh file above. Delete it.
416 // const fileName meshFiles = io.time().path()/io.instance()/meshSubDir;
417 // //Pout<< "Removing dummy mesh " << meshFiles << endl;
418 // mesh.removeFiles();
419 // rmDir(meshFiles);
420 // }
421 //
422  // Force recreation of globalMeshData.
423 // mesh.clearOut();
424  mesh.globalData();
425 
426 
427  // Do some checks.
428 
429  // Check if the boundary definition is unique
430  mesh.boundaryMesh().checkDefinition(true);
431  // Check if the boundary processor patches are correct
432  mesh.boundaryMesh().checkParallelSync(true);
433  // Check names of zones are equal
434  mesh.cellZones().checkDefinition(true);
435  mesh.cellZones().checkParallelSync(true);
436  mesh.faceZones().checkDefinition(true);
437  mesh.faceZones().checkParallelSync(true);
438  mesh.pointZones().checkDefinition(true);
439  mesh.pointZones().checkParallelSync(true);
440 
441  return meshPtr;
442 }
443 
444 
445 // ************************************************************************* //
meshPtr
Foam::autoPtr< Foam::fvMesh > meshPtr(nullptr)
nPatches
label nPatches
Definition: readKivaGrid.H:396
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::fileHandler
const fileOperation & fileHandler()
Get current file handler.
Definition: fileOperation.C:1354
IOobjectList.H
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:350
polyBoundaryMeshEntries.H
processorCyclicPolyPatch.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::wordList
List< word > wordList
A List of words.
Definition: fileName.H:59
Foam::name
word name(const complex &c)
Return string representation of complex.
Definition: complex.C:76
faceSet.H
Foam::FatalError
error FatalError
processorPolyPatch.H
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Time.H
Foam::autoPtr< Foam::fvMesh >
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:381
Foam::fileOperation::isFile
virtual bool isFile(const fileName &, const bool checkGzip=true, const bool followLink=true) const =0
Does the name exist as a FILE in the file system?
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:65
Foam::roots::type
type
Types of root.
Definition: Roots.H:54
cellSet.H
basicFvGeometryScheme.H
Foam::loadOrCreateMesh
autoPtr< fvMesh > loadOrCreateMesh(const IOobject &io)
Load (if it exists) or create zero cell mesh given an IOobject:
loadOrCreateMesh.H
Load or create (0 size) a mesh. Used in distributing meshes to a larger number of processors.
pointSet.H