33 bool Foam::blockMesh::readPatches
35 const dictionary& meshDescription,
42 bool topologyOK =
true;
44 ITstream& patchStream(meshDescription.lookup(
"patches"));
49 token firstToken(patchStream);
51 if (firstToken.isLabel())
53 nPatches = firstToken.labelToken();
55 tmpBlocksPatches.setSize(nPatches);
56 patchNames.setSize(nPatches);
57 patchTypes.setSize(nPatches);
58 nbrPatchNames.setSize(nPatches);
62 patchStream.putBack(firstToken);
66 patchStream.readBegin(
"patches");
70 token lastToken(patchStream);
74 lastToken.isPunctuation()
79 if (tmpBlocksPatches.size() <=
nPatches)
81 tmpBlocksPatches.setSize(nPatches + 1);
82 patchNames.setSize(nPatches + 1);
83 patchTypes.setSize(nPatches + 1);
84 nbrPatchNames.setSize(nPatches + 1);
87 patchStream.putBack(lastToken);
94 patchStream >> tmpBlocksPatches[
nPatches];
100 if (patchNames[nPatches] == patchNames[i])
103 <<
"Duplicate patch " << patchNames[
nPatches]
104 <<
" at line " << patchStream.lineNumber()
105 <<
". Exiting !" <<
nl 110 topologyOK = topologyOK && patchLabelsOK
114 tmpBlocksPatches[nPatches]
122 if (patchTypes[nPatches-1] == cyclicPolyPatch::typeName)
124 word halfA = patchNames[nPatches-1] +
"_half0";
125 word halfB = patchNames[nPatches-1] +
"_half1";
128 <<
"Old-style cyclic definition." 129 <<
" Splitting patch " 130 << patchNames[nPatches-1] <<
" into two halves " 131 << halfA <<
" and " << halfB << endl
132 <<
" Alternatively use new 'boundary' dictionary syntax." 136 if (tmpBlocksPatches.size() <=
nPatches)
138 tmpBlocksPatches.setSize(nPatches + 1);
139 patchNames.setSize(nPatches + 1);
140 patchTypes.setSize(nPatches + 1);
141 nbrPatchNames.setSize(nPatches + 1);
145 patchNames[nPatches-1] = halfA;
146 nbrPatchNames[nPatches-1] = halfB;
148 patchTypes[
nPatches] = patchTypes[nPatches-1];
153 if ((tmpBlocksPatches[nPatches-1].
size() % 2) != 0)
156 <<
"Size of cyclic faces is not a multiple of 2 :" 157 << tmpBlocksPatches[nPatches-1]
160 label sz = tmpBlocksPatches[nPatches-1].size()/2;
161 faceList unsplitFaces(tmpBlocksPatches[nPatches-1],
true);
162 tmpBlocksPatches[nPatches-1] =
faceList 164 SubList<face>(unsplitFaces, sz)
168 SubList<face>(unsplitFaces, sz, sz)
174 patchStream >> lastToken;
176 patchStream.putBack(lastToken);
179 patchStream.readEnd(
"patches");
185 bool Foam::blockMesh::readBoundary
187 const dictionary& meshDescription,
190 PtrList<dictionary>& patchDicts
193 bool topologyOK =
true;
196 const PtrList<entry> patchesInfo
198 meshDescription.lookup(
"boundary")
201 patchNames.setSize(patchesInfo.size());
202 tmpBlocksPatches.setSize(patchesInfo.size());
203 patchDicts.setSize(patchesInfo.size());
205 forAll(tmpBlocksPatches, patchi)
207 const entry& patchInfo = patchesInfo[
patchi];
209 if (!patchInfo.isDict())
212 <<
"Entry " << patchInfo <<
" in boundary section is not a" 216 patchNames[
patchi] = patchInfo.keyword();
218 patchDicts.set(patchi,
new dictionary(patchInfo.dict()));
220 patchDicts[
patchi].lookup(
"faces") >> tmpBlocksPatches[
patchi];
222 topologyOK = topologyOK && patchLabelsOK
226 tmpBlocksPatches[patchi]
234 void Foam::blockMesh::createCellShapes
239 const blockMesh& blocks = *
this;
241 tmpBlockCells.setSize(blocks.size());
244 tmpBlockCells[blockI] = cellShape(blocks[blockI].blockShape());
246 if (tmpBlockCells[blockI].
mag(blockPointField_) < 0.0)
249 <<
"negative volume block : " << blockI
250 <<
", probably defined inside-out" <<
endl;
260 const IOdictionary& meshDescription,
261 const word& regionName
264 bool topologyOK =
true;
268 word defaultPatchName =
"defaultFaces";
269 word defaultPatchType = emptyPolyPatch::typeName;
274 if (
const dictionary* dictPtr = meshDescription.subDictPtr(
"defaultPatch"))
276 dictPtr->readIfPresent(
"name", defaultPatchName);
277 dictPtr->readIfPresent(
"type", defaultPatchType);
281 if (!meshDescription.readIfPresent(
"convertToMeters", scaleFactor_))
283 meshDescription.readIfPresent(
"scale", scaleFactor_);
290 if (meshDescription.found(
"edges"))
294 Info<<
"Creating curved edges" <<
endl;
297 ITstream& is(meshDescription.lookup(
"edges"));
302 token firstToken(is);
304 if (firstToken.isLabel())
306 nEdges = firstToken.labelToken();
311 is.putBack(firstToken);
315 is.readBegin(
"edges");
323 lastToken.isPunctuation()
328 if (edges_.
size() <= nEdges)
333 is.putBack(lastToken);
345 is.putBack(lastToken);
350 else if (verboseOutput)
352 Info<<
"No non-linear edges defined" <<
endl;
361 Info<<
"Creating topology blocks" <<
endl;
365 ITstream& is(meshDescription.lookup(
"blocks"));
370 token firstToken(is);
372 if (firstToken.isLabel())
374 nBlocks = firstToken.labelToken();
375 blocks.setSize(nBlocks);
379 is.putBack(firstToken);
383 is.readBegin(
"blocks");
391 lastToken.isPunctuation()
396 if (blocks.size() <= nBlocks)
398 blocks.setSize(nBlocks + 1);
401 is.putBack(lastToken);
414 topologyOK = topologyOK && blockLabelsOK
418 blocks[nBlocks].blockShape()
425 is.putBack(lastToken);
428 is.readEnd(
"blocks");
432 polyMesh* blockMeshPtr = NULL;
439 Info<<
"Creating topology patches" <<
endl;
442 if (meshDescription.found(
"patches"))
444 Info<<
nl <<
"Reading patches section" <<
endl;
451 topologyOK = topologyOK && readPatches
463 <<
"Cannot create mesh due to errors in topology, exiting !" 467 Info<<
nl <<
"Creating block mesh topology" <<
endl;
470 createCellShapes(tmpBlockCells);
473 Info<<
nl <<
"Reading physicalType from existing boundary file" <<
endl;
475 PtrList<dictionary>
patchDicts(patchNames.size());
480 meshDescription.time(),
481 meshDescription.time().constant(),
491 forAll(patchDicts, patchi)
493 if (!patchDicts.set(patchi))
495 patchDicts.set(patchi,
new dictionary());
498 dictionary& dict = patchDicts[
patchi];
501 if (!dict.found(
"type"))
503 dict.add(
"type", patchTypes[patchi],
false);
505 else if (word(dict.lookup(
"type")) != patchTypes[patchi])
510 ) <<
"For patch " << patchNames[
patchi]
511 <<
" overriding type '" << patchTypes[
patchi]
512 <<
"' with '" << word(dict.lookup(
"type"))
513 <<
"' (read from boundary file)" 520 dict.set(
"neighbourPatch", nbrPatchNames[patchi]);
525 blockMeshPtr =
new polyMesh
530 meshDescription.time().constant(),
531 meshDescription.time(),
545 else if (meshDescription.found(
"boundary"))
551 topologyOK = topologyOK && readBoundary
562 <<
"Cannot create mesh due to errors in topology, exiting !" 567 Info<<
nl <<
"Creating block mesh topology" <<
endl;
570 createCellShapes(tmpBlockCells);
574 blockMeshPtr =
new polyMesh
579 meshDescription.time().constant(),
580 meshDescription.time(),
595 checkBlockMesh(*blockMeshPtr);
List< faceList > faceListList
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
wordList patchTypes(nPatches)
PtrList< block > blockList
A PtrList of blocks.
Ostream & endl(Ostream &os)
Add newline and flush stream.
wordList patchNames() const
Return patch names.
Xfer< T > xferCopy(const T &)
Construct by copying the contents of the arg.
PtrList< dictionary > patchDicts() const
Get patch information from the topology mesh.
List< cellShape > cellShapeList
List of cellShapes and PtrList of List of cellShape.
bool set(const label) const
Is element set.
static const word null
An empty word.
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
static autoPtr< curvedEdge > New(const pointField &, Istream &)
New function which constructs and returns pointer to a curvedEdge.
void preservePatchTypes(const objectRegistry &obr, const word &meshInstance, const fileName &meshDir, const wordList &patchNames, PtrList< dictionary > &patchDicts, const word &defaultFacesName, word &defaultFacesType)
Preserve patch types.
List< word > wordList
A List of words.
#define WarningInFunction
Report a warning using Foam::Warning.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
dimensioned< scalar > mag(const dimensioned< Type > &)
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
Mesh consisting of general polyhedral cells.
label size() const
Return the number of elements in the UPtrList.