31 void Foam::patchToPoly2DMesh::flipFaceOrder()
37 Info<<
"Flipping face order if necessary." <<
endl;
40 const edge&
e = edges[edgeI];
44 label edgeOwner = owner_[edgeI];
46 const face&
f = localFaces[edgeOwner];
50 if (
f.nextLabel(fp) !=
e[1])
52 Info<<
"Flipping face " << faces_[edgeI] <<
endl;
53 faces_[edgeI][0] = meshPoints[
e[1]];
54 faces_[edgeI][1] = meshPoints[
e[0]];
58 faces_[edgeI][0] = meshPoints[
e[0]];
59 faces_[edgeI][1] = meshPoints[
e[1]];
65 void Foam::patchToPoly2DMesh::createNeighbours()
67 const edgeList& edges = patch_.edges();
70 Info<<
"Calculating neighbours." <<
endl;
73 const labelList& eFaces = edgeFaces[edgeI];
74 if (eFaces.size() == 2)
76 if (owner_[edgeI] == eFaces[0])
78 neighbour_[edgeI] = eFaces[1];
82 neighbour_[edgeI] = eFaces[0];
85 else if (eFaces.size() == 1)
108 const labelList& fEdges = faceEdges[facei];
110 SortableList<label> nbr(fEdges.size(), -1);
114 if (fEdges[feI] < neighbour_.size())
118 label nbrFacei = neighbour_[fEdges[feI]];
120 if (nbrFacei == facei)
122 nbrFacei = owner_[fEdges[feI]];
125 if (facei < nbrFacei)
139 oldToNew[fEdges[nbr.indices()[i]]] = newFacei++;
148 void Foam::patchToPoly2DMesh::addPatchFacesToFaces()
150 const labelList& meshPoints = patch_.meshPoints();
161 f[0] = meshPoints[eIter.key().start()];
162 f[1] = meshPoints[eIter.key().end()];
172 void Foam::patchToPoly2DMesh::addPatchFacesToOwner()
174 const label nInternalEdges = patch_.nInternalEdges();
175 const faceList& faces = patch_.faces();
176 const label nExternalEdges = patch_.edges().size() - nInternalEdges;
177 const labelList& meshPoints = patch_.meshPoints();
186 label bFacei = nInternalEdges;
187 bFacei < faces_.size();
191 const face&
e = faces_[bFacei];
193 bool matched =
false;
197 label bEdgeI = nInternalEdges;
198 bEdgeI < faces_.size();
204 e[0] == meshPoints[patch_.edges()[bEdgeI][0]]
205 &&
e[1] == meshPoints[patch_.edges()[bEdgeI][1]]
208 const face&
f = faces[owner_[bEdgeI]];
212 newOwner[bFacei] = owner_[bEdgeI];
214 if (
f.nextLabel(fp) !=
e[1])
218 faces_[bFacei][0] =
e[1];
219 faces_[bFacei][1] =
e[0];
228 e[0] == meshPoints[patch_.edges()[bEdgeI][1]]
229 &&
e[1] == meshPoints[patch_.edges()[bEdgeI][0]]
232 Info<<
"Warning: Wrong orientation." <<
endl;
239 Info<<
"No match for edge." <<
endl;
243 if (nMatched != nExternalEdges)
245 Info<<
"Number of matched edges, " << nMatched
246 <<
", does not match number of external edges, "
247 << nExternalEdges <<
endl;
250 owner_ = move(newOwner);
254 void Foam::patchToPoly2DMesh::createPolyMeshComponents()
261 labelList oldToNew = internalFaceOrder();
268 addPatchFacesToFaces();
270 addPatchFacesToOwner();
278 const MeshedSurface<face>& patch,
281 const EdgeMap<label>& mapEdgesRegion
286 patchSizes_(patchSizes),
288 mapEdgesRegion_(mapEdgesRegion),
290 faces_(patch.nEdges()),
291 owner_(PatchTools::edgeOwner(patch)),
292 neighbour_(patch.nInternalEdges())
306 for (
label edgeI = 0; edgeI < patch_.nInternalEdges(); edgeI++)
308 if (patch_.edgeFaces()[edgeI].size() != 2)
311 <<
"internal edge:" << edgeI
312 <<
" patch.edgeFaces()[edgeI]:" << patch_.edgeFaces()[edgeI]
319 label edgeI = patch_.nInternalEdges();
320 edgeI < patch_.nEdges();
324 if (patch_.edgeFaces()[edgeI].size() != 1)
327 <<
"boundary edge:" << edgeI
328 <<
" patch.edgeFaces()[edgeI]:" << patch_.edgeFaces()[edgeI]
333 createPolyMeshComponents();
335 label startFace = patch_.nInternalEdges();
338 patchStarts_[
patchi] = startFace;
339 startFace += patchSizes_[
patchi];
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
void clear()
Clear the list, i.e. set size to zero.
void setSize(const label)
Reset size of List.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
const List< FaceType > & localFaces() const
Return patch faces addressing into local point list.
void createMesh()
Create the mesh.
~patchToPoly2DMesh()
Destructor.
patchToPoly2DMesh(const MeshedSurface< face > &patch, const wordList &patchNames, const labelList &patchSizes, const EdgeMap< label > &mapEdgesRegion)
Construct from a primitivePatch.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
List< word > wordList
A List of words.
List< label > labelList
A List of labels.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
errorManip< error > abort(error &err)
List< labelList > labelListList
A List of labelList.
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
void offset(label &lst, const label o)
void inplaceReorder(const labelUList &oldToNew, ListType &)
Inplace reorder the elements of a list.
wordList patchNames(nPatches)