52 const word& pType =
p.type();
53 const word& pName =
p.
name();
55 bool nameFound =
false;
57 forAll(patchNames_, patchi)
59 if (patchNames_[patchi] == pName)
61 if (patchDicts_[patchi].get<word>(
"type") == pType)
78 patchDicts_.append(dictionary(IStringStream(
os.str())()));
89 Info<<
"label patchIndex(const polyPatch& p) : "
90 <<
"Patch " <<
p.index() <<
" named "
92 <<
" already exists, but patch types "
93 <<
" do not match.\nCreating a composite name as "
101 return patchNames_.
size() - 1;
105Foam::label Foam::mergePolyMesh::zoneIndex
107 DynamicList<word>& names,
113 if (names[zonei] == curName)
126void Foam::mergePolyMesh::sortProcessorPatches()
128 Info<<
"Reordering processor patches last" <<
endl;
137 const polyBoundaryMesh& oldPatches = boundaryMesh();
139 DynamicList<polyPatch*> newPatches(oldPatches.size());
141 labelList oldToSorted(oldPatches.size());
143 forAll(oldPatches, patchi)
145 const polyPatch& pp = oldPatches[patchi];
147 if (!isA<processorPolyPatch>(pp))
149 oldToSorted[patchi] = newPatches.size();
162 forAll(oldPatches, patchi)
164 const polyPatch& pp = oldPatches[patchi];
166 if (isA<processorPolyPatch>(pp))
168 oldToSorted[patchi] = newPatches.size();
184 addPatches(newPatches);
188 DynamicList<label>&
patchID =
const_cast<DynamicList<label>&
>
210 patchNames_(2*boundaryMesh().size()),
211 patchDicts_(2*boundaryMesh().size()),
217 wordList curPatchNames = boundaryMesh().names();
219 forAll(boundaryMesh(), patchi)
221 patchNames_.append(boundaryMesh()[patchi].
name());
224 boundaryMesh()[patchi].write(
os);
225 patchDicts_.append(dictionary(IStringStream(
os.str())()));
231 wordList curPointZoneNames = pointZones().names();
232 if (curPointZoneNames.size())
234 pointZoneNames_.setCapacity(2*curPointZoneNames.size());
235 pointZoneNames_.
append(curPointZoneNames);
239 wordList curFaceZoneNames = faceZones().names();
240 if (curFaceZoneNames.size())
242 faceZoneNames_.setCapacity(2*curFaceZoneNames.size());
243 faceZoneNames_.
append(curFaceZoneNames);
247 wordList curCellZoneNames = cellZones().names();
248 if (curCellZoneNames.size())
250 cellZoneNames_.setCapacity(2*curCellZoneNames.size());
251 cellZoneNames_.
append(curCellZoneNames);
274 pointZoneIndices[zoneI] = zoneIndex(pointZoneNames_, pz[zoneI].
name());
280 zoneID = pz.whichZone(pointi);
288 renumberPoints[pointi] =
311 cellZoneIndices[zoneI] = zoneIndex(cellZoneNames_, cz[zoneI].
name());
317 zoneID = cz.whichZone(celli);
325 renumberCells[celli] =
340 const polyBoundaryMesh& bm = m.boundaryMesh();
345 forAll(patchIndices, patchi)
347 patchIndices[patchi] = patchIndex(bm[patchi]);
361 faceZoneIndices[zoneI] = zoneIndex(faceZoneNames_, fz[zoneI].
name());
368 const labelList& nei = m.faceNeighbour();
370 label newOwn, newNei, newPatch, newZone;
375 const face& curFace =
f[facei];
377 face newFace(curFace.size());
381 newFace[pointi] = renumberPoints[curFace[pointi]];
387 if (
min(newFace) < 0)
390 <<
"Error in point mapping for face " << facei
391 <<
". Old face: " << curFace <<
" New face: " << newFace
396 if (facei < m.nInternalFaces() || facei >= m.nFaces())
402 newPatch = patchIndices[bm.whichPatch(facei)];
406 if (newOwn > -1) newOwn = renumberCells[newOwn];
415 newNei = renumberCells[newNei];
419 newZone = fz.whichZone(facei);
424 newZoneFlip = fz[newZone].flipMap()[fz[newZone].whichFace(facei)];
427 newZone = faceZoneIndices[newZone];
430 renumberFaces[facei] =
453 Info<<
"patch names: " << patchNames_ <<
nl
454 <<
"patch dicts: " << patchDicts_ <<
nl
455 <<
"point zone names: " << pointZoneNames_ <<
nl
456 <<
"face zone names: " << faceZoneNames_ <<
nl
457 <<
"cell zone names: " << cellZoneNames_ <<
endl;
460 if (patchNames_.size() != boundaryMesh().size())
462 Info<<
"Copying old patches" <<
endl;
464 List<polyPatch*> newPatches(patchNames_.size());
466 const polyBoundaryMesh& oldPatches = boundaryMesh();
471 for (patchi = 0; patchi < oldPatches.size(); patchi++)
473 newPatches[patchi] = oldPatches[patchi].clone(oldPatches).ptr();
476 Info<<
"Adding new patches. " <<
endl;
478 label endOfLastPatch =
481 : oldPatches[patchi - 1].start() + oldPatches[patchi - 1].size();
483 for (; patchi < patchNames_.size(); patchi++)
486 dictionary
dict(patchDicts_[patchi]);
488 dict.
set(
"startFace", endOfLastPatch);
503 addPatches(newPatches);
507 if (pointZoneNames_.size() > pointZones().size())
509 Info<<
"Adding new pointZones." <<
endl;
511 label zonei = pointZones().size();
513 const label
nZones = pointZoneNames_.size();
515 pointZones().setSize(
nZones);
517 for (; zonei <
nZones; ++zonei)
522 new pointZone(pointZoneNames_[zonei], zonei, pointZones())
526 if (cellZoneNames_.size() > cellZones().size())
528 Info<<
"Adding new cellZones." <<
endl;
530 label zonei = cellZones().size();
532 const label
nZones = cellZoneNames_.size();
534 cellZones().setSize(cellZoneNames_.size());
536 for (; zonei <
nZones; ++zonei)
541 new cellZone(cellZoneNames_[zonei], zonei, cellZones())
545 if (faceZoneNames_.size() > faceZones().size())
547 Info<<
"Adding new faceZones." <<
endl;
549 label zonei = faceZones().size();
551 const label
nZones = faceZoneNames_.size();
553 faceZones().setSize(
nZones);
555 for (; zonei <
nZones; ++zonei)
560 new faceZone(faceZoneNames_[zonei], zonei, faceZones())
567 sortProcessorPatches();
570 meshMod_.changeMesh(*
this,
false);
const Mesh & mesh() const
Return mesh.
void append(const T &val)
Copy append an element to the end of this list.
const fileName & caseName() const
Return the Time::caseName()
const word & name() const noexcept
Return the object name.
void append(const T &val)
Append an element at the end of the list.
const fileName & caseName() const
Return case name.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
void size(const label n)
Older name for setAddressableSize.
T & last()
Return the last element of the list.
entry * set(entry *entryPtr)
Assign a new entry, overwriting any existing entry.
const Time & time() const
Return reference to time.
Add a given mesh to the original mesh to create a single new mesh.
void addMesh(const polyMesh &m)
Add a mesh.
void merge()
Merge meshes.
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
void setNumPatches(const label nPatches)
label patchIndex() const
Return access to the patch index.
virtual bool write(const bool valid=true) const
Write using setting from DB.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
const labelIOList & zoneID
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
List< word > names(const UPtrList< T > &list, const UnaryMatchPredicate &matcher)
const dimensionedScalar c
Speed of light in a vacuum.
List< word > wordList
A List of words.
List< label > labelList
A List of labels.
ZoneMesh< cellZone, polyMesh > cellZoneMesh
A ZoneMesh with the type cellZone.
List< cell > cellList
A List of cells.
ZoneMesh< pointZone, polyMesh > pointZoneMesh
A ZoneMesh with the type pointZone.
messageStream Info
Information stream (stdout output on master, null elsewhere)
vectorField pointField
pointField is a vectorField.
ZoneMesh< faceZone, polyMesh > faceZoneMesh
A ZoneMesh with the type faceZone.
Ostream & endl(Ostream &os)
Add newline and flush stream.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManip< error > abort(error &err)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
List< face > faceList
A List of faces.
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.