40const Foam::scalar Foam::attachDetach::positionDifference_ = 1
e-8;
44void Foam::attachDetach::attachInterface
64 Pout<<
"void attachDetach::attachInterface("
65 <<
"polyTopoChange& ref) const "
66 <<
" for object " <<
name() <<
" : "
67 <<
"Attaching interface" <<
endl;
78 const label masterPatchStart = masterPatch.
start();
79 const label slavePatchStart = slavePatch.start();
81 const labelList& slaveMeshPoints = slavePatch.meshPoints();
83 const Map<label>& removedPointMap = pointMatchMap();
85 const labelList removedPoints = removedPointMap.toc();
87 forAll(removedPoints, pointi)
93 ref.setAction(polyRemovePoint(removedPoints[pointi]));
109 ref.setAction(polyRemoveFace(i + slavePatchStart));
113 const labelList& masterFaceCells = masterPatch.faceCells();
114 const labelList& slaveFaceCells = slavePatch.faceCells();
121 forAll(masterFaceCells, facei)
125 if (masterFaceCells[facei] < slaveFaceCells[facei])
131 faces[masterPatchStart + facei],
132 masterPatchStart + facei,
133 masterFaceCells[facei],
134 slaveFaceCells[facei],
150 faces[masterPatchStart + facei].reverseFace(),
151 masterPatchStart + facei,
152 slaveFaceCells[facei],
153 masterFaceCells[facei],
162 faceModified[masterPatchStart + facei] =
true;
177 forAll(slaveMeshPoints, pointi)
179 const labelList& curFaces = pf[slaveMeshPoints[pointi]];
185 !
ref.faceRemoved(curFaces[facei])
186 && !faceModified[curFaces[facei]]
189 facesToModifyMap.insert(curFaces[facei]);
195 const labelList ftm = facesToModifyMap.toc();
202 label curFaceID = ftm[facei];
204 face newFace(faces[curFaceID]);
208 const auto rpmIter = removedPointMap.cfind(newFace[pointi]);
213 newFace[pointi] = rpmIter();
224 bool modifiedFaceZoneFlip =
false;
226 if (modifiedFaceZone >= 0)
228 modifiedFaceZoneFlip =
240 neiCell = nei[curFaceID];
268 Pout<<
"void attachDetach::attachInterface("
269 <<
"polyTopoChange& ref) const "
270 <<
" for object " <<
name() <<
" : "
271 <<
"Finished attaching interface" <<
endl;
281 const Map<label>& removedPointMap = pointMatchMap();
287 Pout<<
"void attachDetach::modifyMotionPoints("
288 <<
"pointField& motionPoints) const "
289 <<
" for object " <<
name() <<
" : "
290 <<
"Adjusting motion points." <<
endl;
293 scalar pointDiff = 0;
295 forAll(removedPoints, pointi)
300 motionPoints[removedPoints[pointi]]
301 - motionPoints[removedPointMap.
find(removedPoints[pointi])()]
305 if (pointDiff > removedPoints.
size()*positionDifference_)
307 Pout<<
"Point motion difference = " << pointDiff <<
endl;
312 forAll(removedPoints, pointi)
314 motionPoints[removedPoints[pointi]] =
315 motionPoints[removedPointMap.
find(removedPoints[pointi])()];
label index() const
The index of the first matching items, -1 if no matches.
List< Key > toc() const
The table of contents (the keys) in unsorted order.
iterator find(const Key &key)
Find and return an iterator set at the hashed entry.
A HashTable to objects of type <T> with a label key.
void size(const label n)
Older name for setAddressableSize.
label whichZone(const label objectIndex) const
Given a global object index, return the zone it is in.
virtual void modifyMotionPoints(pointField &motionPoints) const
Modify motion points to comply with the topological change.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
label start() const noexcept
The start label of boundary faces in the polyMesh face list.
const word & name() const
Return name of this modifier.
const polyTopoChanger & topoChanger() const
Return reference to morph engine.
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
virtual const faceList & faces() const
Return raw faces.
virtual const labelList & faceOwner() const
Return face owner.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const labelList & faceNeighbour() const
Return face neighbour.
const polyMesh & mesh() const
Return the mesh reference.
static const unsigned facesPerPoint_
Estimated number of faces per point.
label nFaces() const noexcept
Number of mesh faces.
const labelListList & pointFaces() const
List< label > labelList
A List of labels.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
List< labelList > labelListList
A List of labelList.
List< bool > boolList
A List of bools.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
List< face > faceList
A List of faces.
#define forAll(list, i)
Loop across all elements in list.