37 void Foam::triangleFuncs::setIntersection
39 const point& oppositeSidePt,
40 const scalar oppositeSign,
42 const point& thisSidePt,
43 const scalar thisSign,
50 const scalar denom = oppositeSign - thisSign;
59 pt = oppositeSidePt + oppositeSign/denom*(thisSidePt - oppositeSidePt);
73 const scalar maxLength,
82 label i1 = (i0 + 1) % 3;
83 label i2 = (i1 + 1) % 3;
93 scalar
det = v2*u1 - u2*v1;
99 if (localScale < VSMALL ||
Foam::mag(
det)/localScale < SMALL)
107 const point& P = origin[originI];
109 scalar
u0 = P[i1] - V0[i1];
110 scalar v0 = P[i2] - V0[i2];
138 scalar
s = (pInter - origin[originI])[i0];
140 if ((
s >= 0) && (
s <= maxLength))
153 const treeBoundBox& cubeBb
174 if (inter.hit() && inter.distance() <= 1)
187 cubeBb.intersects(tri.a(), tri.b(), pInter)
188 || cubeBb.intersects(tri.b(), tri.c(), pInter)
189 || cubeBb.intersects(tri.c(), tri.a(), pInter)
204 return intersectBb(tri, cubeBb);
223 scalar magArea =
mag(na);
226 if (
mag(na & normal) > (1 - SMALL))
232 const point va1 = va0 + va10;
233 const point va2 = va0 + va20;
236 scalar sign0 = (va0 - base) & normal;
237 scalar sign1 = (va1 - base) & normal;
238 scalar sign2 = (va2 - base) & normal;
240 label oppositeVertex = -1;
303 if (oppositeVertex == 0)
306 setIntersection(va0, sign0, va1, sign1, tol, pInter0);
307 setIntersection(va0, sign0, va2, sign2, tol, pInter1);
309 else if (oppositeVertex == 1)
312 setIntersection(va1, sign1, va0, sign0, tol, pInter0);
313 setIntersection(va1, sign1, va2, sign2, tol, pInter1);
318 setIntersection(va2, sign2, va0, sign0, tol, pInter0);
319 setIntersection(va2, sign2, va1, sign1, tol, pInter1);
350 if (!intersect(va0, va10, va20, vb0, nb, planeB0, planeB1))
358 if (!intersect(vb0, vb10, vb20, va0, na, planeA0, planeA1))
366 vector intersection(na ^ nb);
368 scalar coordB0 = planeB0 & intersection;
369 scalar coordB1 = planeB1 & intersection;
371 scalar coordA0 = planeA0 & intersection;
372 scalar coordA1 = planeA1 & intersection;
377 SortableList<scalar> sortCoords(4);
381 sortCoords[0] = coordB0;
385 sortCoords[1] = coordB1;
389 sortCoords[2] = coordA0;
393 sortCoords[3] = coordA1;
397 const labelList& indices = sortCoords.indices();
399 if (isFromB[indices[0]] == isFromB[indices[1]])
408 pInter0 = *
pts[indices[1]];
409 pInter1 = *
pts[indices[2]];
A triangle primitive used to calculate face normals and swept volumes. Uses referred points...
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar det(const dimensionedSphericalTensor &dt)
#define forAll(list, i)
Loop across all elements in list.
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
static const edgeList edges
Edge to point addressing, using octant corner points.
static bool intersect(const point &va0, const point &va10, const point &va20, const point &basePoint, const vector &normal, point &pInter0, point &pInter1)
Intersect triangle with plane.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
vector point
Point is a vector.
triangle< point, const point & > triPointRef
A triangle using referred points.
Standard boundBox with extra functionality for use in octree.
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
List< label > labelList
A List of labels.
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))
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
List< bool > boolList
A List of bools.
const volScalarField & p0
PointHit< point > pointHit
A PointHit with a 3D point.
static bool intersectAxesBundle(const point &V0, const point &V10, const point &V20, const label i0, const pointField &origin, const scalar maxLength, point &pInter)
Intersect triangle with parallel edges aligned with axis i0.
static bool intersectBb(const triPointRef &tri, const treeBoundBox &cubeBb)
Intersect triangle with bounding box.