55 return (f0.edgeDirection(
e) > 0) ^ (f1.edgeDirection(
e) > 0);
69 labelList changedEdges(3*changedFaces.size());
72 for (
const label facei : changedFaces)
74 for (
const label edgei :
s.faceEdges()[facei])
76 changedEdges[changedI++] = edgei;
92 labelList changedFaces(2*changedEdges.size());
95 for (
const label edgeI : changedEdges)
97 const labelList& eFaces =
s.edgeFaces()[edgeI];
99 if (eFaces.size() < 2)
103 else if (eFaces.size() == 2)
105 const label face0 = eFaces[0];
106 const label face1 = eFaces[1];
111 if (flip[face0] == UNVISITED)
113 if (flip[face1] == UNVISITED)
124 flip[face0] = (flip[face1] == FLIP ? FLIP : NOFLIP);
129 flip[face0] = (flip[face1] == FLIP ? NOFLIP : FLIP);
131 changedFaces[changedI++] = face0;
136 if (flip[face1] == UNVISITED)
141 flip[face1] = (flip[face0] == FLIP ? FLIP : NOFLIP);
145 flip[face1] = (flip[face0] == FLIP ? NOFLIP : FLIP);
147 changedFaces[changedI++] = face1;
156 changedFaces.setSize(changedI);
162 void Foam::orientedSurface::walkSurface
165 const label startFacei,
176 changedEdges = faceToEdge(
s, changedFaces);
178 if (changedEdges.empty())
183 changedFaces = edgeToFace(
s, changedEdges, flipState);
185 if (changedFaces.empty())
193 void Foam::orientedSurface::propagateOrientation
196 const point& samplePoint,
197 const bool orientOutside,
198 const label nearestFacei,
199 const point& nearestPt,
217 flipState[nearestFacei] = NOFLIP;
223 flipState[nearestFacei] = NOFLIP;
228 flipState[nearestFacei] = FLIP;
235 Pout<<
"orientedSurface::propagateOrientation : starting face" 236 <<
" orientation:" <<
nl 237 <<
" for samplePoint:" << samplePoint <<
nl 238 <<
" starting from point:" << nearestPt <<
nl 239 <<
" on face:" << nearestFacei <<
nl 240 <<
" with normal:" <<
n <<
nl 241 <<
" decided side:" << label(side)
246 walkSurface(
s, nearestFacei, flipState);
252 void Foam::orientedSurface::findZoneSide
254 const triSurfaceSearch& surfSearches,
257 const point& outsidePoint,
262 const triSurface&
s = surfSearches.surface();
268 List<List<pointIndexHit>> hits(1, List<pointIndexHit>());
272 if (faceZone[facei] == zoneI)
274 const point& fc =
s.faceCentres()[facei];
275 const vector&
n =
s.faceNormals()[facei];
277 const vector d = fc - outsidePoint;
278 const scalar magD =
mag(d);
281 if (magD > SMALL && (
mag(
n & d/magD) > 1
e-6))
292 surfSearches.findLineAll(start,
end, hits);
294 label zoneIndex = -1;
297 if (hits[0][i].index() == facei)
308 if ((zoneIndex%2) == 0)
312 isOutside = ((
n & d) < 0);
316 isOutside = ((
n & d) > 0);
326 bool Foam::orientedSurface::flipSurface
332 bool hasFlipped =
false;
337 if (flipState[facei] == UNVISITED)
340 <<
"unvisited face " << facei
343 else if (flipState[facei] == FLIP)
345 labelledTri& tri =
s[facei];
364 bool Foam::orientedSurface::orientConsistent(triSurface&
s)
366 bool anyFlipped =
false;
382 label startFacei = -1;
383 while (facei <
s.size())
385 if (flipState[facei] == UNVISITED)
393 if (startFacei == -1)
398 flipState[startFacei] = NOFLIP;
399 walkSurface(
s, startFacei, flipState);
402 anyFlipped = flipSurface(
s, flipState);
421 const point& samplePoint,
422 const bool orientOutside
427 orient(*
this, samplePoint, orientOutside);
435 const bool orientOutside
443 point outsidePoint = bb.max() + bb.span();
445 orient(*
this, outsidePoint, orientOutside);
454 const point& samplePoint,
455 const bool orientOutside
461 bool topoFlipped = orientConsistent(
s);
475 scalar minDist = GREAT;
481 if (flipState[facei] == UNVISITED)
484 s[facei].nearestPoint(samplePoint,
s.points());
489 minPoint = curHit.
point();
515 bool geomFlipped = flipSurface(
s, flipState);
517 return topoFlipped || geomFlipped;
525 const point& samplePoint,
526 const bool orientOutside
532 bool topoFlipped = orientConsistent(
s);
538 if (
s.edgeFaces()[edgeI].size() != 2)
540 borderEdge[edgeI] =
true;
544 label nZones =
s.markZones(borderEdge, faceZone);
549 for (label zoneI = 0; zoneI < nZones; zoneI++)
551 label zoneFacei = -1;
564 if (isOutside == orientOutside)
566 flipState[zoneFacei] = NOFLIP;
570 flipState[zoneFacei] = FLIP;
572 walkSurface(
s, zoneFacei, flipState);
576 bool geomFlipped = flipSurface(
s, flipState);
578 return topoFlipped || geomFlipped;
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
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.
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
static bool orient(triSurface &, const point &, const bool orientOutside)
Flip faces such that normals are consistent with point:
orientedSurface()
Default construct.
const labelList & meshPoints() const
Return labelList of mesh points in patch.
#define forAll(list, i)
Loop across all elements in list.
Helper class to search on triSurface.
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
void setSize(const label n)
Alias for resize()
An edge is a list of two vertex labels. This can correspond to a directed graph edge or an edge on a ...
labelledTri face_type
The face type (same as the underlying PrimitivePatch)
Describes the interaction of a object and a (templated) point. It carries the info of a successful hi...
const Field< point_type > & points() const noexcept
Return reference to global points.
errorManip< error > abort(error &err)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
defineTypeNameAndDebug(combustionModel, 0)
vector point
Point is a vector.
scalar distance() const noexcept
Return distance to hit.
const point_type & point() const noexcept
Return the point, no checks.
Standard boundBox with extra functionality for use in octree.
List< label > labelList
A List of labels.
Triangulated surface description with patch information.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
static bool consistentEdge(const edge &e, const Face &f0, const Face &f1)