Commit 6d73c32b authored by Mark Olesen's avatar Mark Olesen
Browse files

ENH: new triFace methods to match those in face. FaceType typedef for surface meshes

parent 0a17e071
......@@ -618,19 +618,17 @@ Foam::face Foam::face::reverseFace() const
Foam::label Foam::face::which(const label globalIndex) const
{
label pointInFace = -1;
const labelList& f = *this;
forAll(f, i)
forAll(f, localIdx)
{
if (f[i] == globalIndex)
if (f[localIdx] == globalIndex)
{
pointInFace = i;
break;
return localIdx;
}
}
return pointInFace;
return -1;
}
......@@ -654,9 +652,7 @@ Foam::scalar Foam::face::sweptVol
point nextOldPoint = centreOldPoint;
point nextNewPoint = centreNewPoint;
register label pI;
for (pI = 0; pI < nPoints; pI++)
for (register label pI = 0; pI < nPoints; ++pI)
{
if (pI < nPoints - 1)
{
......@@ -708,20 +704,18 @@ Foam::tensor Foam::face::inertia
).inertia(refPt, density);
}
point c = centre(p);
const point ctr = centre(p);
tensor J = tensor::zero;
forAll(*this, i)
{
triPointRef t
J += triPointRef
(
p[operator[](i)],
p[operator[](fcIndex(i))],
c
);
J += t.inertia(refPt, density);
ctr
).inertia(refPt, density);
}
return J;
......@@ -734,9 +728,7 @@ Foam::edgeList Foam::face::edges() const
edgeList e(points.size());
label pointI;
for (pointI = 0; pointI < points.size() - 1; pointI++)
for (label pointI = 0; pointI < points.size() - 1; ++pointI)
{
e[pointI] = edge(points[pointI], points[pointI + 1]);
}
......@@ -839,9 +831,4 @@ Foam::label Foam::face::trianglesQuads
}
// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// ************************************************************************* //
......@@ -27,12 +27,16 @@ Class
Description
A face is a list of labels corresponding to mesh vertices.
SeeAlso
Foam::triFace
SourceFiles
faceI.H
face.C
faceIntersection.C
faceContactSphere.C
faceAreaInContact.C
faceTemplates.C
\*---------------------------------------------------------------------------*/
......@@ -191,6 +195,7 @@ public:
//- Navigation through face vertices
//- Which vertex on face (face index given a global index)
// returns -1 if not found
label which(const label globalIndex) const;
//- Next vertex on face
......@@ -289,8 +294,8 @@ public:
//- Return number of edges
inline label nEdges() const;
//- Return edges in face point ordering, i.e. edges()[0] is edge
// between [0] and [1]
//- Return edges in face point ordering,
// i.e. edges()[0] is edge between [0] and [1]
edgeList edges() const;
//- Return n-th face edge
......
......@@ -35,7 +35,7 @@ inline Foam::label Foam::face::right(const label i) const
// Edge to the left of face vertex i
inline Foam::label Foam::face::left(const label i) const
{
return i ? i-1 : size()-1;
return rcIndex(i);
}
......
......@@ -21,9 +21,6 @@ License
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Description
Return intersection of a line with the face
\*---------------------------------------------------------------------------*/
#include "face.H"
......@@ -51,6 +48,17 @@ Foam::pointHit Foam::face::ray
const intersection::direction dir
) const
{
// If the face is a triangle, do a direct calculation
if (size() == 3)
{
return triPointRef
(
meshPoints[operator[](0)],
meshPoints[operator[](1)],
meshPoints[operator[](2)]
).ray(p, n, alg, dir);
}
point ctr = Foam::average(points(meshPoints));
scalar nearestHitDist = GREAT;
......@@ -139,6 +147,17 @@ Foam::pointHit Foam::face::intersection
const scalar tol
) const
{
// If the face is a triangle, do a direct calculation
if (size() == 3)
{
return triPointRef
(
meshPoints[operator[](0)],
meshPoints[operator[](1)],
meshPoints[operator[](2)]
).intersection(p, q, alg, tol);
}
scalar nearestHitDist = VGREAT;
// Initialize to miss, distance = GREAT
......@@ -198,6 +217,17 @@ Foam::pointHit Foam::face::nearestPointClassify
label& nearLabel
) const
{
// If the face is a triangle, do a direct calculation
if (size() == 3)
{
return triPointRef
(
meshPoints[operator[](0)],
meshPoints[operator[](1)],
meshPoints[operator[](2)]
).nearestPointClassify(p, nearType, nearLabel);
}
const face& f = *this;
point ctr = centre(meshPoints);
......
......@@ -50,7 +50,7 @@ template<class Type>
Type Foam::face::average
(
const pointField& meshPoints,
const Field<Type>& f
const Field<Type>& fld
) const
{
// Calculate the average by breaking the face into triangles and
......@@ -62,9 +62,9 @@ Type Foam::face::average
return
(1.0/3.0)
*(
f[operator[](0)]
+ f[operator[](1)]
+ f[operator[](2)]
fld[operator[](0)]
+ fld[operator[](1)]
+ fld[operator[](2)]
);
}
......@@ -76,7 +76,7 @@ Type Foam::face::average
for (register label pI=0; pI<nPoints; pI++)
{
centrePoint += meshPoints[operator[](pI)];
cf += f[operator[](pI)];
cf += fld[operator[](pI)];
}
centrePoint /= nPoints;
......@@ -90,8 +90,8 @@ Type Foam::face::average
// Calculate 3*triangle centre field value
Type ttcf =
(
f[operator[](pI)]
+ f[operator[]((pI + 1) % nPoints)]
fld[operator[](pI)]
+ fld[operator[]((pI + 1) % nPoints)]
+ cf
);
......@@ -116,4 +116,5 @@ Type Foam::face::average
}
}
// ************************************************************************* //
......@@ -25,7 +25,7 @@ Class
Foam::oppositeFace
Description
Class containing opposite face for a prismatic cell with addresing
Class containing opposite face for a prismatic cell with addressing
and a possibility of failure.
SourceFiles
......
......@@ -25,10 +25,15 @@ Class
Foam::triFace
Description
A triangle face primitive using a FixedList.
A triangular face using a FixedList of labels corresponding to mesh
vertices.
SeeAlso
Foam::face, Foam::triangle
SourceFiles
triFaceI.H
triFaceTemplates.C
\*---------------------------------------------------------------------------*/
......@@ -47,17 +52,17 @@ SourceFiles
namespace Foam
{
class face;
// Forward declaration of friend functions and operators
class face;
class triFace;
inline bool operator==(const triFace&, const triFace&);
inline bool operator!=(const triFace&, const triFace&);
/*---------------------------------------------------------------------------*\
class triFace Declaration
Class triFace Declaration
\*---------------------------------------------------------------------------*/
class triFace
......@@ -93,62 +98,108 @@ public:
// return the collapsed size, set collapsed point labels to -1
inline label collapse();
//- Return the edge direction on the face
// - +1: forward (counter-clockwise) on the face
// - -1: reverse (clockwise) on the face
// - 0: edge not found on the face
inline int edgeDirection(const edge&) const;
//- Return the points corresponding to this face
inline pointField points(const pointField& meshPoints) const;
//- Return triangle as a face
inline face triFaceFace() const;
// Properties
//- Return the triangle
inline triPointRef tri(const pointField&) const;
//- Return the points corresponding to this face
inline pointField points(const pointField& points) const;
//- Return centre (centroid)
inline point centre(const pointField&) const;
//- Return triangle as a face
inline face triFaceFace() const;
//- Calculate average value at centroid of face
template<class Type>
Type average(const pointField&, const Field<Type>&) const;
//- Return number of edges
inline label nEdges() const;
//- Return scalar magnitude
inline scalar mag(const pointField&) const;
//- Return edges
inline edgeList edges() const;
//- Return vector normal
inline vector normal(const pointField&) const;
//- Return centre (centroid)
inline point centre(const pointField&) const;
//- Number of triangles after splitting
inline label nTriangles() const;
//- Return scalar magnitude
inline scalar mag(const pointField&) const;
//- Return face with reverse direction
inline triFace reverseFace() const;
//- Return vector normal
inline vector normal(const pointField&) const;
//- Return swept-volume
inline scalar sweptVol
(
const pointField& oldPoints,
const pointField& newPoints
) const;
//- Number of triangles after splitting
inline label nTriangles() const;
//- Return the inertia tensor, with optional reference
// point and density specification
inline tensor inertia
(
const pointField&,
const point& refPt = vector::zero,
scalar density = 1.0
) const;
//- Return point intersection with a ray starting at p,
// with direction q.
inline pointHit ray
(
const point& p,
const vector& q,
const pointField& points,
const intersection::algorithm = intersection::FULL_RAY,
const intersection::direction dir = intersection::VECTOR
) const;
//- Fast intersection with a ray.
inline pointHit intersection
(
const point& p,
const vector& q,
const pointField& points,
const intersection::algorithm alg,
const scalar tol = 0.0
) const;
//- Return nearest point to face
inline pointHit nearestPoint
(
const point& p,
const pointField& points
) const;
//- Return nearest point to face and classify it:
// + near point (nearType=POINT, nearLabel=0, 1, 2)
// + near edge (nearType=EDGE, nearLabel=0, 1, 2)
// Note: edges are counted from starting vertex so
// e.g. edge n is from f[n] to f[0], where the face has n + 1
// points
inline pointHit nearestPointClassify
(
const point& p,
const pointField& points,
label& nearType,
label& nearLabel
) const;
//- Return face with reverse direction
inline triFace reverseFace() const;
//- Return number of edges
inline label nEdges() const;
//- Return swept-volume
inline scalar sweptVol
(
const pointField& oldPoints,
const pointField& newPoints
) const;
//- Return edges in face point ordering,
// i.e. edges()[0] is edge between [0] and [1]
inline edgeList edges() const;
//- Return point intersection with a ray starting at p, with
// direction n.
inline pointHit ray
(
const point& p,
const vector& q,
const pointField& points,
const intersection::algorithm = intersection::FULL_RAY,
const intersection::direction dir = intersection::VECTOR
) const;
//- Return n-th face edge
inline edge faceEdge(const label n) const;
//- Return the triangle
inline triPointRef tri(const pointField&) const;
//- Return the edge direction on the face
// - +1: forward (counter-clockwise) on the face
// - -1: reverse (clockwise) on the face
// - 0: edge not found on the face
inline int edgeDirection(const edge&) const;
//- compare triFaces
// - 0: different
......@@ -199,6 +250,10 @@ inline bool contiguous<triFace>() {return true;}
#include "triFaceI.H"
#ifdef NoRepository
# include "triFaceTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
......
......@@ -142,57 +142,14 @@ inline Foam::face Foam::triFace::triFaceFace() const
}
inline Foam::label Foam::triFace::nEdges() const
{
return 3;
}
inline Foam::edgeList Foam::triFace::edges() const
{
edgeList e(3);
e[0].start() = operator[](0);
e[0].end() = operator[](1);
e[1].start() = operator[](1);
e[1].end() = operator[](2);
e[2].start() = operator[](2);
e[2].end() = operator[](0);
return e;
}
// return
// - +1: forward (counter-clockwise) on the face
// - -1: reverse (clockwise) on the face
// - 0: edge not found on the face
inline int Foam::triFace::edgeDirection(const edge& e) const
inline Foam::triPointRef Foam::triFace::tri(const pointField& points) const
{
if
(
(operator[](0) == e.start() && operator[](1) == e.end())
|| (operator[](1) == e.start() && operator[](2) == e.end())
|| (operator[](2) == e.start() && operator[](0) == e.end())
)
{
return 1;
}
else if
return triPointRef
(
(operator[](0) == e.end() && operator[](1) == e.start())
|| (operator[](1) == e.end() && operator[](2) == e.start())
|| (operator[](2) == e.end() && operator[](0) == e.start())
)
{
return -1;
}
else
{
return 0;
}
points[operator[](0)],
points[operator[](1)],
points[operator[](2)]
);
}
......@@ -269,6 +226,18 @@ inline Foam::scalar Foam::triFace::sweptVol
}
Foam::tensor Foam::triFace::inertia
(
const pointField& points,
const point& refPt,
scalar density
) const
{
// a triangle, do a direct calculation
return this->tri(points).inertia(refPt, density);
}
inline Foam::pointHit Foam::triFace::ray
(
const point& p,
......@@ -278,23 +247,103 @@ inline Foam::pointHit Foam::triFace::ray
const intersection::direction dir
) const
{
return triPointRef
(
points[operator[](0)],
points[operator[](1)],
points[operator[](2)]
).ray(p, q, alg, dir);
return this->tri(points).ray(p, q, alg, dir);
}
inline Foam::triPointRef Foam::triFace::tri(const pointField& points) const
inline Foam::pointHit Foam::triFace::intersection
(
const point& p,
const vector& q,
const pointField& points,
const intersection::algorithm alg,
const scalar tol
) const
{
return triPointRef
return this->tri(points).intersection(p, q, alg, tol);
}
inline Foam::pointHit Foam::triFace::nearestPoint
(
const point& p,
const pointField& points
) const
{
return this->tri(points).nearestPoint(p);
}
inline Foam::pointHit Foam::triFace::nearestPointClassify
(
const point& p,
const pointField& points,
label& nearType,
label& nearLabel
) const
{
return this->tri(points).nearestPointClassify(p, nearType, nearLabel);
}
inline Foam::label Foam::triFace::nEdges() const
{
return 3;
}
inline Foam::edgeList Foam::triFace::edges() const
{
edgeList e(3);
e[0].start() = operator[](0);
e[0].end() = operator[](1);
e[1].start() = operator[](1);
e[1].end() = operator[](2);
e[2].start() = operator[](2);
e[2].end() = operator[](0);
return e;
}
inline Foam::edge Foam::triFace::faceEdge(const label n) const
{
return edge(operator[](n), operator[](fcIndex(n)));
}