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();
59 <<
"Can not find edge " << mesh_.
edges()[edgeI]
68 Foam::label Foam::cellFeatures::nextEdge
70 const Map<label>& toSuperFace,
71 const label superFacei,
72 const label thisEdgeI,
76 const labelList& pEdges = mesh_.pointEdges()[thisVertI];
80 label edgeI = pEdges[pEdgeI];
82 if ((edgeI != thisEdgeI) && featureEdge_.found(edgeI))
86 const labelList& eFaces = mesh_.edgeFaces()[edgeI];
90 label facei = eFaces[eFacei];
95 && (toSuperFace[facei] == superFacei)
105 <<
"Can not find edge in " << featureEdge_ <<
" connected to edge " 106 << thisEdgeI <<
" at vertex " << thisVertI <<
endl 107 <<
"This might mean that the externalEdges do not form a closed loop" 115 bool Foam::cellFeatures::isCellFeatureEdge
132 scalar cosAngle = n0 & n1;
134 const edge&
e = mesh_.edges()[edgeI];
136 const face& f0 = mesh_.faces()[face0];
138 label face0Start = f0.find(
e.start());
139 label face0End = f0.fcIndex(face0Start);
141 const face& f1 = mesh_.faces()[face1];
143 label face1Start = f1.find(
e.start());
144 label face1End = f1.fcIndex(face1Start);
149 (f0[face0End] ==
e.end())
150 && (f1[face1End] !=
e.end())
153 (f0[face0End] !=
e.end())
154 && (f1[face1End] ==
e.end())
161 cosAngle = -cosAngle;
164 if (cosAngle < minCos)
175 void Foam::cellFeatures::walkSuperFace
178 const label superFacei,
179 Map<label>& toSuperFace
182 if (toSuperFace.insert(facei, superFacei))
184 const labelList& fEdges = mesh_.faceEdges()[facei];
186 for (
const label edgeI : fEdges)
188 if (!featureEdge_.found(edgeI))
211 void Foam::cellFeatures::calcSuperFaces()
const 215 const labelList& cFaces = mesh_.cells()[celli_];
220 Map<label> toSuperFace(10*cFaces.size());
222 label superFacei = 0;
226 label facei = cFaces[cFacei];
228 if (!toSuperFace.found(facei))
242 faceMap_.setSize(superFacei);
246 label facei = cFaces[cFacei];
248 faceMap_[toSuperFace[facei]].append(facei);
253 faceMap_[superI].shrink();
259 facesPtr_.reset(
new faceList(superFacei));
260 auto& faces = *facesPtr_;
264 label facei = cFaces[cFacei];
266 label superFacei = toSuperFace[facei];
268 if (faces[superFacei].empty())
273 label startEdgeI = -1;
275 const labelList& fEdges = mesh_.faceEdges()[facei];
279 label edgeI = fEdges[fEdgeI];
281 if (featureEdge_.found(edgeI))
290 if (startEdgeI != -1)
294 DynamicList<label> superFace(10*mesh_.faces()[facei].size());
296 const edge&
e = mesh_.edges()[startEdgeI];
300 bool flipOrientation =
301 (mesh_.faceOwner()[facei] == celli_)
302 ^ (faceAlignedEdge(facei, startEdgeI));
304 label startVertI = -1;
308 startVertI =
e.end();
312 startVertI =
e.start();
315 label edgeI = startEdgeI;
317 label vertI =
e.otherVertex(startVertI);
321 label newEdgeI = nextEdge
330 if (isFeaturePoint(edgeI, newEdgeI))
332 superFace.append(vertI);
337 if (vertI == startVertI)
342 vertI = mesh_.edges()[edgeI].otherVertex(vertI);
346 if (superFace.size() <= 2)
349 <<
" Can not collapse faces " << faceMap_[superFacei]
350 <<
" into one big face on cell " << celli_ <<
endl 351 <<
"Try decreasing minCos:" << minCos_ <<
endl;
355 faces[superFacei].transfer(superFace);
366 Foam::cellFeatures::cellFeatures
376 featureEdge_(10*
mesh.cellEdges()[celli].size()),
384 label edgeI = cEdges[cEdgeI];
386 if (isCellFeatureEdge(minCos_, edgeI))
388 featureEdge_.
insert(edgeI);
408 || (edge0 >= mesh_.nEdges())
410 || (edge1 >= mesh_.nEdges())
414 <<
"Illegal edge labels : edge0:" << edge0 <<
" edge1:" << edge1
418 const edge& e0 = mesh_.edges()[edge0];
420 const vector e0Vec = e0.unitVec(mesh_.points());
422 const edge& e1 = mesh_.edges()[edge1];
424 const vector e1Vec = e1.unitVec(mesh_.points());
430 (e0.start() == e1.end())
431 || (e0.end() == e1.start())
435 cosAngle = e0Vec & e1Vec;
439 (e0.start() == e1.start())
440 || (e0.end() == e1.end())
444 cosAngle = - e0Vec & e1Vec;
451 <<
"Edges do not share common vertex. e0:" << e0
455 if (cosAngle < minCos_)
474 || (facei >= mesh_.nFaces())
476 || (vertI >= mesh_.nPoints())
480 <<
"Illegal face " << facei <<
" or vertex " << vertI
484 const labelList& pEdges = mesh_.pointEdges()[vertI];
491 label edgeI = pEdges[pEdgeI];
512 <<
"Did not find two edges sharing vertex " << vertI
513 <<
" on face " << facei <<
" vertices:" << mesh_.faces()[facei]
517 return isFeaturePoint(edge0, edge1);
const labelListList & cellEdges() const
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Cell-face mesh analysis engine.
bool isFeatureVertex(const label facei, const label vertI) const
Is vertexI on facei used by two edges that form feature.
Ostream & endl(Ostream &os)
Add newline and flush stream.
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Various functions to operate on Lists.
#define forAll(list, i)
Loop across all elements in list.
List< face > faceList
List of faces.
const dimensionedScalar e
Elementary charge.
~cellFeatures()
Destructor.
errorManip< error > abort(error &err)
bool isFeaturePoint(const label edge0, const label edge1) const
Are two edges connected at feature point?
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
#define WarningInFunction
Report a warning using Foam::Warning.
virtual const faceList & faces() const =0
Return faces.
List< label > labelList
A List of labels.