41 bool Foam::cellFeatures::faceAlignedEdge(
const label faceI,
const label edgeI)
44 const edge&
e = mesh_.
edges()[edgeI];
46 const face&
f = mesh_.
faces()[faceI];
50 if (f[fp] == e.start())
52 label fp1 = f.fcIndex(fp);
54 return f[fp1] == e.end();
60 "cellFeatures::faceAlignedEdge(const label, const label)"
61 ) <<
"Can not find edge " << mesh_.
edges()[edgeI]
70 Foam::label Foam::cellFeatures::nextEdge
72 const Map<label>& toSuperFace,
73 const label superFaceI,
74 const label thisEdgeI,
78 const labelList& pEdges = mesh_.pointEdges()[thisVertI];
82 label edgeI = pEdges[pEdgeI];
84 if ((edgeI != thisEdgeI) && featureEdge_.found(edgeI))
88 const labelList& eFaces = mesh_.edgeFaces()[edgeI];
92 label faceI = eFaces[eFaceI];
97 && (toSuperFace[faceI] == superFaceI)
108 "cellFeatures::nextEdge(const label, const Map<label>"
109 ", const labelHashSet&, const label, const label, const label)"
110 ) <<
"Can not find edge in " << featureEdge_ <<
" connected to edge "
111 << thisEdgeI <<
" at vertex " << thisVertI <<
endl
112 <<
"This might mean that the externalEdges do not form a closed loop"
120 bool Foam::cellFeatures::isCellFeatureEdge
134 vector n0 = mesh_.faceAreas()[face0];
137 vector n1 = mesh_.faceAreas()[face1];
140 scalar cosAngle = n0 & n1;
143 const edge& e = mesh_.edges()[edgeI];
145 const face& f0 = mesh_.faces()[face0];
147 label face0Start =
findIndex(f0, e.start());
148 label face0End = f0.fcIndex(face0Start);
150 const face& f1 = mesh_.faces()[face1];
152 label face1Start =
findIndex(f1, e.start());
153 label face1End = f1.fcIndex(face1Start);
158 (f0[face0End] == e.end())
159 && (f1[face1End] != e.end())
162 (f0[face0End] != e.end())
163 && (f1[face1End] == e.end())
170 cosAngle = -cosAngle;
173 if (cosAngle < minCos)
186 void Foam::cellFeatures::walkSuperFace
189 const label superFaceI,
190 Map<label>& toSuperFace
193 if (!toSuperFace.found(faceI))
195 toSuperFace.insert(faceI, superFaceI);
197 const labelList& fEdges = mesh_.faceEdges()[faceI];
201 label edgeI = fEdges[fEdgeI];
203 if (!featureEdge_.found(edgeI))
226 void Foam::cellFeatures::calcSuperFaces()
const
230 const labelList& cFaces = mesh_.cells()[cellI_];
235 Map<label> toSuperFace(10*cFaces.size());
237 label superFaceI = 0;
241 label faceI = cFaces[cFaceI];
243 if (!toSuperFace.found(faceI))
257 faceMap_.setSize(superFaceI);
261 label faceI = cFaces[cFaceI];
263 faceMap_[toSuperFace[faceI]].append(faceI);
268 faceMap_[superI].shrink();
274 facesPtr_ =
new faceList(superFaceI);
280 label faceI = cFaces[cFaceI];
282 label superFaceI = toSuperFace[faceI];
284 if (faces[superFaceI].empty())
289 label startEdgeI = -1;
291 const labelList& fEdges = mesh_.faceEdges()[faceI];
295 label edgeI = fEdges[fEdgeI];
297 if (featureEdge_.found(edgeI))
306 if (startEdgeI != -1)
310 DynamicList<label> superFace(10*mesh_.faces()[faceI].size());
312 const edge& e = mesh_.edges()[startEdgeI];
316 bool flipOrientation =
317 (mesh_.faceOwner()[faceI] == cellI_)
318 ^ (faceAlignedEdge(faceI, startEdgeI));
320 label startVertI = -1;
324 startVertI = e.end();
328 startVertI = e.start();
331 label edgeI = startEdgeI;
333 label vertI = e.otherVertex(startVertI);
337 label newEdgeI = nextEdge
346 if (isFeaturePoint(edgeI, newEdgeI))
348 superFace.append(vertI);
353 if (vertI == startVertI)
358 vertI = mesh_.edges()[edgeI].otherVertex(vertI);
362 if (superFace.size() <= 2)
364 WarningIn(
"cellFeatures::calcSuperFaces")
365 <<
" Can not collapse faces " << faceMap_[superFaceI]
366 <<
" into one big face on cell " << cellI_ <<
endl
367 <<
"Try decreasing minCos:" << minCos_ <<
endl;
371 faces[superFaceI].transfer(superFace);
382 Foam::cellFeatures::cellFeatures
396 const labelList& cEdges = mesh_.cellEdges()[cellI_];
400 label edgeI = cEdges[cEdgeI];
402 if (isCellFeatureEdge(minCos_, edgeI))
404 featureEdge_.insert(edgeI);
427 || (edge0 >= mesh_.nEdges())
429 || (edge1 >= mesh_.nEdges())
434 "cellFeatures::isFeatureVertex(const label, const label)"
435 ) <<
"Illegal edge labels : edge0:" << edge0 <<
" edge1:" << edge1
439 const edge& e0 = mesh_.edges()[edge0];
444 const edge& e1 = mesh_.edges()[edge1];
458 cosAngle = e0Vec & e1Vec;
467 cosAngle = - e0Vec & e1Vec;
475 "cellFeatures::isFeaturePoint(const label, const label"
477 ) <<
"Edges do not share common vertex. e0:" << e0
481 if (cosAngle < minCos_)
499 || (faceI >= mesh_.nFaces())
501 || (vertI >= mesh_.nPoints())
506 "cellFeatures::isFeatureVertex(const label, const label)"
507 ) <<
"Illegal face " << faceI <<
" or vertex " << vertI
511 const labelList& pEdges = mesh_.pointEdges()[vertI];
518 label edgeI = pEdges[pEdgeI];
540 "cellFeatures::isFeatureVertex(const label, const label)"
541 ) <<
"Did not find two edges sharing vertex " << vertI
542 <<
" on face " << faceI <<
" vertices:" << mesh_.faces()[faceI]
546 return isFeaturePoint(edge0, edge1);