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 | Copyright (C) 2012-2016 OpenFOAM Foundation
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 \*---------------------------------------------------------------------------*/
25 
26 #include "loadOrCreateMesh.H"
27 #include "processorPolyPatch.H"
29 #include "Time.H"
30 #include "IOPtrList.H"
31 
32 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTemplateTypeNameAndDebug(IOPtrList<entry>, 0);
37 }
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 // Read mesh if available. Otherwise create empty mesh with same non-proc
42 // patches as proc0 mesh. Requires all processors to have all patches
43 // (and in same order).
45 (
46  const IOobject& io
47 )
48 {
49  fileName meshSubDir;
50 
51  if (io.name() == polyMesh::defaultRegion)
52  {
53  meshSubDir = polyMesh::meshSubDir;
54  }
55  else
56  {
57  meshSubDir = io.name()/polyMesh::meshSubDir;
58  }
59 
60 
61  // Scatter master patches
62  PtrList<entry> patchEntries;
63  if (Pstream::master())
64  {
65  // Read PtrList of dictionary as dictionary.
66  const word oldTypeName = IOPtrList<entry>::typeName;
67  const_cast<word&>(IOPtrList<entry>::typeName) = word::null;
68  IOPtrList<entry> dictList
69  (
70  IOobject
71  (
72  "boundary",
73  io.time().findInstance
74  (
75  meshSubDir,
76  "boundary",
78  ),
79  meshSubDir,
80  io.db(),
83  false
84  )
85  );
86  const_cast<word&>(IOPtrList<entry>::typeName) = oldTypeName;
87  // Fake type back to what was in field
88  const_cast<word&>(dictList.type()) = dictList.headerClassName();
89 
90  patchEntries.transfer(dictList);
91 
92  // Send patches
93  for
94  (
95  int slave=Pstream::firstSlave();
96  slave<=Pstream::lastSlave();
97  slave++
98  )
99  {
100  OPstream toSlave(Pstream::scheduled, slave);
101  toSlave << patchEntries;
102  }
103  }
104  else
105  {
106  // Receive patches
107  IPstream fromMaster(Pstream::scheduled, Pstream::masterNo());
108  fromMaster >> patchEntries;
109  }
110 
111 
112 
113  // Check who has a mesh
114  const bool haveMesh = isDir(io.time().path()/io.instance()/meshSubDir);
115 
116  if (!haveMesh)
117  {
118  bool oldParRun = Pstream::parRun();
119  Pstream::parRun() = false;
120 
121 
122  // Create dummy mesh. Only used on procs that don't have mesh.
123  IOobject noReadIO(io);
124  noReadIO.readOpt() = IOobject::NO_READ;
125  fvMesh dummyMesh
126  (
127  noReadIO,
128  xferCopy(pointField()),
129  xferCopy(faceList()),
130  xferCopy(labelList()),
131  xferCopy(labelList()),
132  false
133  );
134 
135  // Add patches
136  List<polyPatch*> patches(patchEntries.size());
137  label nPatches = 0;
138 
139  forAll(patchEntries, patchi)
140  {
141  const entry& e = patchEntries[patchi];
142  const word type(e.dict().lookup("type"));
143  const word& name = e.keyword();
144 
145  if
146  (
147  type != processorPolyPatch::typeName
148  && type != processorCyclicPolyPatch::typeName
149  )
150  {
151  dictionary patchDict(e.dict());
152  patchDict.set("nFaces", 0);
153  patchDict.set("startFace", 0);
154 
155  patches[patchi] = polyPatch::New
156  (
157  name,
158  patchDict,
159  nPatches++,
160  dummyMesh.boundaryMesh()
161  ).ptr();
162  }
163  }
164  patches.setSize(nPatches);
165  dummyMesh.addFvPatches(patches, false); // no parallel comms
166 
167  // Add some dummy zones so upon reading it does not read them
168  // from the undecomposed case. Should be done as extra argument to
169  // regIOobject::readStream?
170  List<pointZone*> pz
171  (
172  1,
173  new pointZone
174  (
175  "dummyPointZone",
176  labelList(0),
177  0,
178  dummyMesh.pointZones()
179  )
180  );
181  List<faceZone*> fz
182  (
183  1,
184  new faceZone
185  (
186  "dummyFaceZone",
187  labelList(0),
188  boolList(0),
189  0,
190  dummyMesh.faceZones()
191  )
192  );
193  List<cellZone*> cz
194  (
195  1,
196  new cellZone
197  (
198  "dummyCellZone",
199  labelList(0),
200  0,
201  dummyMesh.cellZones()
202  )
203  );
204  dummyMesh.addZones(pz, fz, cz);
205  //Pout<< "Writing dummy mesh to " << dummyMesh.polyMesh::objectPath()
206  // << endl;
207  dummyMesh.write();
208 
209  Pstream::parRun() = oldParRun;
210  }
211 
212  //Pout<< "Reading mesh from " << io.objectPath() << endl;
213  autoPtr<fvMesh> meshPtr(new fvMesh(io));
214  fvMesh& mesh = meshPtr();
215 
216 
217  // Sync patches
218  // ~~~~~~~~~~~~
219 
220  if (!Pstream::master() && haveMesh)
221  {
222  // Check master names against mine
223 
224  const polyBoundaryMesh& patches = mesh.boundaryMesh();
225 
226  forAll(patchEntries, patchi)
227  {
228  const entry& e = patchEntries[patchi];
229  const word type(e.dict().lookup("type"));
230  const word& name = e.keyword();
231 
232  if (type == processorPolyPatch::typeName)
233  {
234  break;
235  }
236 
237  if (patchi >= patches.size())
238  {
240  << "Non-processor patches not synchronised."
241  << endl
242  << "Processor " << Pstream::myProcNo()
243  << " has only " << patches.size()
244  << " patches, master has "
245  << patchi
246  << exit(FatalError);
247  }
248 
249  if
250  (
251  type != patches[patchi].type()
252  || name != patches[patchi].name()
253  )
254  {
256  << "Non-processor patches not synchronised."
257  << endl
258  << "Master patch " << patchi
259  << " name:" << type
260  << " type:" << type << endl
261  << "Processor " << Pstream::myProcNo()
262  << " patch " << patchi
263  << " has name:" << patches[patchi].name()
264  << " type:" << patches[patchi].type()
265  << exit(FatalError);
266  }
267  }
268  }
269 
270 
271  // Determine zones
272  // ~~~~~~~~~~~~~~~
273 
274  wordList pointZoneNames(mesh.pointZones().names());
275  Pstream::scatter(pointZoneNames);
276  wordList faceZoneNames(mesh.faceZones().names());
277  Pstream::scatter(faceZoneNames);
278  wordList cellZoneNames(mesh.cellZones().names());
279  Pstream::scatter(cellZoneNames);
280 
281  if (!haveMesh)
282  {
283  // Add the zones. Make sure to remove the old dummy ones first
284  mesh.pointZones().clear();
285  mesh.faceZones().clear();
286  mesh.cellZones().clear();
287 
288  List<pointZone*> pz(pointZoneNames.size());
289  forAll(pointZoneNames, i)
290  {
291  pz[i] = new pointZone
292  (
293  pointZoneNames[i],
294  labelList(0),
295  i,
296  mesh.pointZones()
297  );
298  }
299  List<faceZone*> fz(faceZoneNames.size());
300  forAll(faceZoneNames, i)
301  {
302  fz[i] = new faceZone
303  (
304  faceZoneNames[i],
305  labelList(0),
306  boolList(0),
307  i,
308  mesh.faceZones()
309  );
310  }
311  List<cellZone*> cz(cellZoneNames.size());
312  forAll(cellZoneNames, i)
313  {
314  cz[i] = new cellZone
315  (
316  cellZoneNames[i],
317  labelList(0),
318  i,
319  mesh.cellZones()
320  );
321  }
322  mesh.addZones(pz, fz, cz);
323  }
324 
325 
326  if (!haveMesh)
327  {
328  // We created a dummy mesh file above. Delete it.
329  const fileName meshFiles = io.time().path()/io.instance()/meshSubDir;
330  //Pout<< "Removing dummy mesh " << meshFiles << endl;
331  mesh.removeFiles();
332  rmDir(meshFiles);
333  }
334 
335  // Force recreation of globalMeshData.
336  mesh.clearOut();
337  mesh.globalData();
338 
339 
340  // Do some checks.
341 
342  // Check if the boundary definition is unique
343  mesh.boundaryMesh().checkDefinition(true);
344  // Check if the boundary processor patches are correct
345  mesh.boundaryMesh().checkParallelSync(true);
346  // Check names of zones are equal
347  mesh.cellZones().checkDefinition(true);
348  mesh.cellZones().checkParallelSync(true);
349  mesh.faceZones().checkDefinition(true);
350  mesh.faceZones().checkParallelSync(true);
351  mesh.pointZones().checkDefinition(true);
352  mesh.pointZones().checkParallelSync(true);
353 
354  return meshPtr;
355 }
356 
357 
358 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
static int masterNo()
Process index of the master.
Definition: UPstream.H:405
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static int firstSlave()
Process index of first slave.
Definition: UPstream.H:434
static fvMesh * meshPtr
Definition: globalFoam.H:52
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:309
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:417
List< face > faceList
Definition: faceListFwd.H:43
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:306
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:411
Xfer< T > xferCopy(const T &)
Construct by copying the contents of the arg.
bool isDir(const fileName &)
Does the name exist as a DIRECTORY in the file system?
Definition: POSIX.C:486
patches[0]
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Load or create (0 size) a mesh. Used in distributing meshes to a larger number of processors...
static const word null
An empty word.
Definition: word.H:77
List< label > labelList
A List of labels.
Definition: labelList.H:56
defineTemplateTypeNameAndDebug(IOPtrList< ensightPart >, 0)
static void scatter(const List< commsStruct > &comms, T &Value, const int tag, const label comm)
Scatter data. Distribute without modification. Reverse of gather.
bool rmDir(const fileName &)
Remove a dirctory and its contents.
Definition: POSIX.C:839
autoPtr< fvMesh > loadOrCreateMesh(const IOobject &io)
Load (if it exists) or create zero cell mesh given an IOobject:
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
List< word > wordList
A List of words.
Definition: fileName.H:54
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:393
label patchi
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
Namespace for OpenFOAM.
static int lastSlave(const label communicator=0)
Process index of last slave.
Definition: UPstream.H:440