37 const Foam::scalar Foam::attachDetach::positionDifference_ = 1
e-8;
41 void Foam::attachDetach::attachInterface
61 Pout<<
"void attachDetach::attachInterface(" 62 <<
"polyTopoChange& ref) const " 63 <<
" for object " <<
name() <<
" : " 64 <<
"Attaching interface" <<
endl;
68 const faceList& faces = mesh.faces();
70 const labelList& nei = mesh.faceNeighbour();
72 const polyPatch& masterPatch = mesh.boundaryMesh()[masterPatchID_.
index()];
73 const polyPatch& slavePatch = mesh.boundaryMesh()[slavePatchID_.
index()];
75 const label masterPatchStart = masterPatch.start();
76 const label slavePatchStart = slavePatch.start();
78 const labelList& slaveMeshPoints = slavePatch.meshPoints();
80 const Map<label>& removedPointMap = pointMatchMap();
82 const labelList removedPoints = removedPointMap.toc();
84 forAll(removedPoints, pointi)
90 ref.setAction(polyRemovePoint(removedPoints[pointi]));
106 ref.setAction(polyRemoveFace(i + slavePatchStart));
110 const labelList& masterFaceCells = masterPatch.faceCells();
111 const labelList& slaveFaceCells = slavePatch.faceCells();
113 const boolList& mfFlip = mesh.faceZones()[faceZoneID_.
index()].flipMap();
115 forAll(masterFaceCells, facei)
119 if (masterFaceCells[facei] < slaveFaceCells[facei])
125 faces[masterPatchStart + facei],
126 masterPatchStart + facei,
127 masterFaceCells[facei],
128 slaveFaceCells[facei],
144 faces[masterPatchStart + facei].reverseFace(),
145 masterPatchStart + facei,
146 slaveFaceCells[facei],
147 masterFaceCells[facei],
170 forAll(slaveMeshPoints, pointi)
172 const labelList& curFaces = pf[slaveMeshPoints[pointi]];
176 if (!ref.faceRemoved(curFaces[facei]))
178 facesToModifyMap.insert(curFaces[facei]);
184 const labelList ftm = facesToModifyMap.toc();
191 label curFaceID = ftm[facei];
193 face newFace(faces[curFaceID]);
198 removedPointMap.find(newFace[pointi]);
200 if (rpmIter != removedPointMap.end())
203 newFace[pointi] = rpmIter();
213 label modifiedFaceZone = mesh.faceZones().whichZone(curFaceID);
214 bool modifiedFaceZoneFlip =
false;
216 if (modifiedFaceZone >= 0)
218 modifiedFaceZoneFlip =
219 mesh.faceZones()[modifiedFaceZone].flipMap()
221 mesh.faceZones()[modifiedFaceZone].whichFace(curFaceID)
226 label patchID = mesh.boundaryMesh().whichPatch(curFaceID);
230 neiCell = nei[curFaceID];
258 Pout<<
"void attachDetach::attachInterface(" 259 <<
"polyTopoChange& ref) const " 260 <<
" for object " <<
name() <<
" : " 261 <<
"Finished attaching interface" <<
endl;
271 const Map<label>& removedPointMap = pointMatchMap();
277 Pout<<
"void attachDetach::modifyMotionPoints(" 278 <<
"pointField& motionPoints) const " 279 <<
" for object " <<
name() <<
" : " 280 <<
"Adjusting motion points." <<
endl;
283 scalar pointDiff = 0;
285 forAll(removedPoints, pointi)
290 motionPoints[removedPoints[pointi]]
291 - motionPoints[removedPointMap.
find(removedPoints[pointi])()]
295 if (pointDiff > removedPoints.
size()*positionDifference_)
297 Pout<<
"Point motion difference = " << pointDiff <<
endl;
302 forAll(removedPoints, pointi)
304 motionPoints[removedPoints[pointi]] =
305 motionPoints[removedPointMap.
find(removedPoints[pointi])()];
List< labelList > labelListList
A List of labelList.
#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.
void size(const label)
Override size to be inconsistent with allocated storage.
Ostream & endl(Ostream &os)
Add newline and flush stream.
HashTable< label, label, Hash< label > >::const_iterator const_iterator
static const unsigned facesPerPoint_
Estimated number of faces per point.
List< bool > boolList
Bool container classes.
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
const dimensionedScalar & e
Elementary charge.
List< label > labelList
A List of labels.
label index() const
Return index of first matching zone.
const polyTopoChanger & topoChanger() const
Return reference to morph engine.
prefixOSstream Pout(cout, "Pout")
const word & name() const
Return name of this modifier.
virtual void modifyMotionPoints(pointField &motionPoints) const
Modify motion points to comply with the topological change.
dimensioned< scalar > mag(const dimensioned< Type > &)
List< Key > toc() const
Return the table of contents.
const polyMesh & mesh() const
Return the mesh reference.