35 void Foam::blockMesh::createPoints()
const 41 prescaling_.x() * scaling_.x(),
42 prescaling_.y() * scaling_.y(),
43 prescaling_.z() * scaling_.z()
48 Info<<
"Creating points with scale " << scaleTot <<
endl;
51 points_.resize(nPoints_);
55 const pointField& blockPoints = blocks[blocki].points();
64 UIndirectList<point>(points_, pointAddr) = blockPoints;
68 Info<<
" Block " << blocki <<
" cell size :" <<
nl;
70 const label v0 = blocks[blocki].pointLabel(0, 0, 0);
73 const label nx = blocks[blocki].density().x();
74 const label v1 = blocks[blocki].pointLabel(1, 0, 0);
75 const label vn = blocks[blocki].pointLabel(nx, 0, 0);
76 const label vn1 = blocks[blocki].pointLabel(nx-1, 0, 0);
78 const scalar cwBeg =
mag(blockPoints[v1] - blockPoints[v0]);
79 const scalar cwEnd =
mag(blockPoints[vn] - blockPoints[vn1]);
82 << cwBeg*scaleTot.x() <<
" .. " 83 << cwEnd*scaleTot.x() <<
nl;
87 const label ny = blocks[blocki].density().y();
88 const label v1 = blocks[blocki].pointLabel(0, 1, 0);
89 const label vn = blocks[blocki].pointLabel(0, ny, 0);
90 const label vn1 = blocks[blocki].pointLabel(0, ny-1, 0);
92 const scalar cwBeg =
mag(blockPoints[v1] - blockPoints[v0]);
93 const scalar cwEnd =
mag(blockPoints[vn] - blockPoints[vn1]);
96 << cwBeg*scaleTot.y() <<
" .. " 97 << cwEnd*scaleTot.y() <<
nl;
101 const label nz = blocks[blocki].density().z();
102 const label v1 = blocks[blocki].pointLabel(0, 0, 1);
103 const label vn = blocks[blocki].pointLabel(0, 0, nz);
104 const label vn1 = blocks[blocki].pointLabel(0, 0, nz-1);
106 const scalar cwBeg =
mag(blockPoints[v1] - blockPoints[v0]);
107 const scalar cwEnd =
mag(blockPoints[vn] - blockPoints[vn1]);
110 << cwBeg*scaleTot.z() <<
" .. " 111 << cwEnd*scaleTot.z() <<
nl;
121 void Foam::blockMesh::createCells()
const 131 cells_.resize(nCells_);
139 for (
const hexCell& blockCell : blocks[blocki].
cells())
141 forAll(cellPoints, cellPointi)
143 cellPoints[cellPointi] =
146 blockCell[cellPointi]
147 + blockOffsets_[blocki]
152 cells_[celli].reset(
hex, cellPoints,
true);
161 const polyPatch& patchTopologyFaces
166 labelList blockLabels = patchTopologyFaces.polyPatch::faceCells();
170 forAll(patchTopologyFaces, patchTopologyFaceLabel)
172 const label blocki = blockLabels[patchTopologyFaceLabel];
174 faceList blockFaces = blocks[blocki].blockShape().faces();
176 forAll(blockFaces, blockFaceLabel)
180 blockFaces[blockFaceLabel]
181 == patchTopologyFaces[patchTopologyFaceLabel]
185 blocks[blocki].boundaryPatches()[blockFaceLabel].size();
195 forAll(patchTopologyFaces, patchTopologyFaceLabel)
197 const label blocki = blockLabels[patchTopologyFaceLabel];
199 faceList blockFaces = blocks[blocki].blockShape().faces();
201 forAll(blockFaces, blockFaceLabel)
205 blockFaces[blockFaceLabel]
206 == patchTopologyFaces[patchTopologyFaceLabel]
209 const List<FixedList<label, 4>>& blockPatchFaces =
210 blocks[blocki].boundaryPatches()[blockFaceLabel];
212 forAll(blockPatchFaces, blockFaceLabel)
220 blockPatchFaces[blockFaceLabel][0]
221 + blockOffsets_[blocki]
228 label facePointLabel = 1;
236 blockPatchFaces[blockFaceLabel][facePointLabel]
237 + blockOffsets_[blocki]
240 if (quad[nUnique] != quad[nUnique-1])
246 if (quad[nUnique-1] == quad[0])
253 patchFaces[faceLabel++] = quad;
255 else if (nUnique == 3)
257 patchFaces[faceLabel++] = face
268 patchFaces.resize(faceLabel);
274 void Foam::blockMesh::createPatches()
const 276 const polyPatchList& topoPatches = topology().boundaryMesh();
283 patches_.resize(topoPatches.size());
285 forAll(topoPatches, patchi)
287 patches_[patchi] = createPatchFaces(topoPatches[patchi]);
299 <<
"topology not allocated" 303 return *topologyPtr_;
310 const polyMesh& origTopo = topology();
312 if (applyTransform && hasPointTransforms())
317 newIO.registerObject(
false);
320 inplacePointTransforms(newPoints);
325 std::move(newPoints),
333 const polyBoundaryMesh& pbmOld = origTopo.
boundaryMesh();
340 newPatches.set(patchi, pbmOld[patchi].clone(pbmNew));
357 const blockMesh& blkMesh = *
this;
361 Info<<
nl <<
"Creating polyMesh from blockMesh" <<
endl;
370 blkMesh.patchNames(),
371 blkMesh.patchDicts(),
373 emptyPolyPatch::typeName
378 const label nZones = blkMesh.numZonedBlocks();
390 HashTable<label> zoneMap(2*nZones);
393 List<DynamicList<label>> zoneCells(nZones);
401 for (
const block&
b : blkMesh)
403 const word& zoneName =
b.zoneName();
404 const label nCellsInBlock =
b.cells().size();
408 const auto iter = zoneMap.cfind(zoneName);
410 label zonei = freeZonei;
418 zoneMap.insert(zoneName, zonei);
423 Info<<
" " << zonei <<
'\t' << zoneName <<
endl;
430 zoneCells[zonei].reserve
432 zoneCells[zonei].size() + nCellsInBlock
435 const label endOfFill = celli + nCellsInBlock;
437 for (; celli < endOfFill; ++celli)
439 zoneCells[zonei].append(celli);
444 celli += nCellsInBlock;
448 List<cellZone*> cz(zoneMap.size());
451 const word& zoneName = iter.key();
452 const label zonei = iter.val();
454 cz[zonei] =
new cellZone
457 zoneCells[zonei].shrink(),
463 pmesh.pointZones().clear();
464 pmesh.faceZones().clear();
465 pmesh.cellZones().clear();
466 pmesh.addZones(List<pointZone*>(), List<faceZone*>(), cz);
SubList< label > labelSubList
A SubList of labels.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
IOstream & hex(IOstream &io)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static refPtr< T > New(Args &&... args)
Construct refPtr with forwarding arguments.
constexpr char nl
The newline '\n' character (0x0a)
List< face > faceList
A List of faces.
PtrList< block > blockList
A PtrList of blocks.
Ostream & endl(Ostream &os)
Add newline and flush stream.
SubList< label > subList
Declare type of subList.
Ignore writing from objectRegistry::writeObject()
A class for managing references or pointers (no reference counting)
#define forAll(list, i)
Loop across all elements in list.
static const cellModel & ref(const modelType model)
Look up reference to cellModel by enumeration. Fatal on failure.
PtrList< polyPatch > polyPatchList
Store lists of polyPatch as a PtrList.
PtrList< block > blockList
The list of blocks is stored as a PtrList.
vectorField pointField
pointField is a vectorField.
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
Foam::autoPtr< Foam::dynamicFvMesh > meshPtr
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
const polyMesh & topology() const
The blockMesh topology as a polyMesh unscaled and non-transformed.
const polyMesh & topoMesh
autoPtr< polyMesh > mesh(const IOobject &io) const
Create polyMesh, with cell zones.
bool inplacePointTransforms(pointField &pts) const
Apply coordinate transforms and scaling.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Mesh consisting of general polyhedral cells.
List< label > labelList
A List of labels.
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
void addPatches(polyPatchList &plist, const bool validBoundary=true)
Add boundary patches.
forAllConstIters(mixture.phases(), phase)