43 Foam::scalar Foam::faceAreaIntersect::tol = 1
e-6;
47 void Foam::faceAreaIntersect::triSliceWithPlane
69 if (
mag(d[i]) < tol*len)
95 || ((nPos == 2) && (nCoPlanar == 1))
96 || ((nPos == 1) && (nCoPlanar == 2))
111 else if ((nPos == 2) && (nCoPlanar == 0))
134 point p01 = planeIntersection(d, tri, i0, i1);
135 point p02 = planeIntersection(d, tri, i0, i2);
139 setTriPoints(tri[i1], tri[i2], p02, nTris, tris);
140 setTriPoints(tri[i1], p02, p01, nTris, tris);
167 point p01 = planeIntersection(d, tri, i1, i0);
168 point p02 = planeIntersection(d, tri, i2, i0);
171 setTriPoints(tri[i0], p01, p02, nTris, tris);
194 point p01 = planeIntersection(d, tri, i1, i0);
199 setTriPoints(tri[i0], p01, tri[i2], nTris, tris);
203 setTriPoints(tri[i0], tri[i2], p01, nTris, tris);
222 Foam::scalar Foam::faceAreaIntersect::triangleIntersect
224 const triPoints& src,
225 const triPoints& tgt,
231 label nWorkTris1 = 0;
234 label nWorkTris2 = 0;
239 scalar t =
sqrt(triArea(src));
246 scalar
s =
mag(tgt[1] - tgt[0]);
247 plane pl0(tgt[0], tgt[1], tgt[1] + s*n);
248 triSliceWithPlane(src, pl0, workTris1, nWorkTris1, t);
261 scalar
s =
mag(tgt[2] - tgt[1]);
262 plane pl1(tgt[1], tgt[2], tgt[2] + s*n);
266 for (
label i = 0; i < nWorkTris1; i++)
268 triSliceWithPlane(workTris1[i], pl1, workTris2, nWorkTris2, t);
282 scalar
s =
mag(tgt[2] - tgt[0]);
283 plane pl2(tgt[2], tgt[0], tgt[0] + s*n);
287 for (
label i = 0; i < nWorkTris2; i++)
289 triSliceWithPlane(workTris2[i], pl2, workTris1, nWorkTris1, t);
300 for (
label i = 0; i < nWorkTris1; i++)
302 area += triArea(workTris1[i]);
344 faceTris[i] =
face(3);
345 faceTris[i][0] = f[0];
346 faceTris[i][1] = f[i + 1];
347 faceTris[i][2] = f[i + 2];
356 label nFaceTris1 = 0;
357 const label nFaceTris2 = f.
triangles(points, nFaceTris1, faceTris);
359 if (nFaceTris != nFaceTris1 || nFaceTris != nFaceTris2)
362 <<
"The numbers of reported triangles in the face do not " 363 <<
"match that generated by the triangulation" 381 triangulate(faceA, pointsA_, triMode, trisA);
382 triangulate(faceB, pointsB_, triMode, trisB);
385 scalar totalArea = 0.0;
388 triPoints tpA = getTriPoints(pointsA_, trisA[tA],
false);
394 triPoints tpB = getTriPoints(pointsB_, trisB[tB], !reverseB_);
398 totalArea += triangleIntersect(tpA, tpB, n);
#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.
errorManipArg< error, int > exit(error &err, const int errNo=1)
A face is a list of labels corresponding to mesh vertices.
A 1D vector of objects of type <T> with a fixed size <Size>.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Geometric class that creates a 2D plane and can return the intersection point between a line and the ...
static const NamedEnum< triangulationMode, 2 > triangulationModeNames_
void resize(const label)
Alias for setSize(const label)
Initialise the NamedEnum HashTable from the static list of names.
label triangles(const pointField &points, label &triI, faceList &triFaces) const
Split into triangles using existing points.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().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))
faceAreaIntersect(const pointField &pointsA, const pointField &pointsB, const bool reverseB=false)
Construct from components.
const point & refPoint() const
Return or return plane base point.
scalar calc(const face &faceA, const face &faceB, const vector &n, const triangulationMode &triMode)
Return area of intersection of faceA with faceB.
label fcIndex(const label i) const
Return the forward circular index, i.e. the next index.
const vector & normal() const
Return plane normal.
dimensioned< scalar > mag(const dimensioned< Type > &)
static void triangulate(const face &f, const pointField &points, const triangulationMode &triMode, faceList &faceTris)
Triangulate a face using the given triangulation mode.
const doubleScalar e
Elementary charge.
label nTriangles() const
Number of triangles after splitting.