tetrahedron: triangle: Improved barycentric handling on tets and tris

Updated the tetrahedron and triangle classes to use the barycentric
primitives. Removed duplicate code for generating random positions in
tets and tris, and fixed bug in tri random position.

applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControl/cellShapeControl.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -97,7 +97,7 @@ Foam::scalarField Foam::cellShapeControl::cellSize
 
 Foam::scalar Foam::cellShapeControl::cellSize(const point& pt) const
 {
-    FixedList<scalar, 4> bary;
+    barycentric bary;
     cellShapeControlMesh::Cell_handle ch;
 
     shapeControlMesh_.barycentricCoords(pt, bary, ch);
@@ -172,7 +172,7 @@ Foam::scalar Foam::cellShapeControl::cellSize(const point& pt) const
 
 Foam::tensor Foam::cellShapeControl::cellAlignment(const point& pt) const
 {
-    FixedList<scalar, 4> bary;
+    barycentric bary;
     cellShapeControlMesh::Cell_handle ch;
 
     shapeControlMesh_.barycentricCoords(pt, bary, ch);
@@ -244,7 +244,7 @@ void Foam::cellShapeControl::cellSizeAndAlignment
     tensor& alignment
 ) const
 {
-    FixedList<scalar, 4> bary;
+    barycentric bary;
     cellShapeControlMesh::Cell_handle ch;
 
     shapeControlMesh_.barycentricCoords(pt, bary, ch);

applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C

@@ -450,7 +450,7 @@ Foam::cellShapeControlMesh::~cellShapeControlMesh()
 void Foam::cellShapeControlMesh::barycentricCoords
 (
     const Foam::point& pt,
-    FixedList<scalar, 4>& bary,
+    barycentric& bary,
     Cell_handle& ch
 ) const
 {
@@ -470,7 +470,7 @@ void Foam::cellShapeControlMesh::barycentricCoords
             topoint(ch->vertex(3)->point())
         );
 
-        tet.barycentric(pt, bary);
+        bary = tet.pointToBarycentric(pt);
     }
 }
 

applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -122,7 +122,7 @@ public:
             void barycentricCoords
             (
                 const Foam::point& pt,
-                FixedList<scalar, 4>& bary,
+                barycentric& bary,
                 Cell_handle& ch
             ) const;
 

applications/utilities/mesh/generation/foamyMesh/conformalVoronoiMesh/cellSizeControlSurfaces/surfaceCellSizeFunction/nonUniformField/nonUniformField.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -118,8 +118,7 @@ Foam::scalar Foam::nonUniformField::interpolate
         pts[faceHitByPt[2]]
     );
 
-    FixedList<scalar, 3> bary;
-    tri.barycentric(pt, bary);
+    const barycentric2D bary = tri.pointToBarycentric(pt);
 
 //    return pointCellSize_[pMap[faceHitByPt[0]]]*bary[0]
 //         + pointCellSize_[pMap[faceHitByPt[1]]]*bary[1]

src/OpenFOAM/Make/files

@@ -134,6 +134,7 @@ primitives/polynomialEqns/cubicEqn/cubicEqn.C
 primitives/polynomialEqns/quadraticEqn/quadraticEqn.C
 
 primitives/Barycentric/barycentric/barycentric.C
+primitives/Barycentric2D/barycentric2D/barycentric2D.C
 
 containers/HashTables/HashTable/HashTableCore.C
 containers/HashTables/StaticHashTable/StaticHashTableCore.C

src/OpenFOAM/meshes/primitiveShapes/tetrahedron/tetrahedron.H

@@ -47,6 +47,7 @@ SourceFiles
 #include "UList.H"
 #include "triPointRef.H"
 #include "boundBox.H"
+#include "barycentric.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -169,13 +170,22 @@ public:
             //  uniform distribution
             inline Point randomPoint(Random& rndGen) const;
 
-            //- Calculate the barycentric coordinates of the given
-            //  point, in the same order as a, b, c, d.  Returns the
-            //  determinant of the solution.
-            inline scalar barycentric
+            //- Return a random point in the tetrahedron from a
+            //  uniform distribution
+            inline Point randomPoint(cachedRandom& rndGen) const;
+
+            //- Calculate the point from the given barycentric coordinates.
+            inline Point barycentricToPoint(const barycentric& bary) const;
+
+            //- Calculate the barycentric coordinates from the given point
+            inline barycentric pointToBarycentric(const point& pt) const;
+
+            //- Calculate the barycentric coordinates from the given point.
+            //  Returns the determinant.
+            inline scalar pointToBarycentric
             (
                 const point& pt,
-                FixedList<scalar, 4>& bary
+                barycentric& bary
             ) const;
 
             //- Return nearest point to p on tetrahedron. Is p itself

src/OpenFOAM/meshes/primitiveShapes/tetrahedron/tetrahedronI.H

@@ -245,41 +245,37 @@ inline Point Foam::tetrahedron<Point, PointRef>::randomPoint
     Random& rndGen
 ) const
 {
-    // Adapted from
-    // http://vcg.isti.cnr.it/activities/geometryegraphics/pointintetraedro.html
+    return barycentricToPoint(barycentric01(rndGen));
+}
 
-    scalar s = rndGen.sample01<scalar>();
-    scalar t = rndGen.sample01<scalar>();
-    scalar u = rndGen.sample01<scalar>();
 
-    if (s + t > 1.0)
-    {
-        s = 1.0 - s;
-        t = 1.0 - t;
-    }
+template<class Point, class PointRef>
+inline Point Foam::tetrahedron<Point, PointRef>::barycentricToPoint
+(
+    const barycentric& bary
+) const
+{
+    return bary[0]*a_ + bary[1]*b_ + bary[2]*c_ + bary[3]*d_;
+}
 
-    if (t + u > 1.0)
-    {
-        scalar tmp = u;
-        u = 1.0 - s - t;
-        t = 1.0 - tmp;
-    }
-    else if (s + t + u > 1.0)
-    {
-        scalar tmp = u;
-        u = s + t + u - 1.0;
-        s = 1.0 - t - tmp;
-    }
 
-    return (1 - s - t - u)*a_ + s*b_ + t*c_ + u*d_;
+template<class Point, class PointRef>
+inline Foam::barycentric Foam::tetrahedron<Point, PointRef>::pointToBarycentric
+(
+    const point& pt
+) const
+{
+    barycentric bary;
+    pointToBarycentric(pt, bary);
+    return bary;
 }
 
 
 template<class Point, class PointRef>
-Foam::scalar Foam::tetrahedron<Point, PointRef>::barycentric
+inline Foam::scalar Foam::tetrahedron<Point, PointRef>::pointToBarycentric
 (
     const point& pt,
-    FixedList<scalar, 4>& bary
+    barycentric& bary
 ) const
 {
     // Reference:
@@ -300,9 +296,9 @@ Foam::scalar Foam::tetrahedron<Point, PointRef>::barycentric
 
     if (Foam::mag(detT) < SMALL)
     {
-        // Degenerate tetrahedron, returning 1/4 barycentric coordinates.
+        // Degenerate tetrahedron, returning 1/4 barycentric coordinates
 
-        bary = List<scalar>(4, 0.25);
+        bary = barycentric(0.25, 0.25, 0.25, 0.25);
 
         return detT;
     }
@@ -312,7 +308,7 @@ Foam::scalar Foam::tetrahedron<Point, PointRef>::barycentric
     bary[0] = res.x();
     bary[1] = res.y();
     bary[2] = res.z();
-    bary[3] = (1.0 - res.x() - res.y() - res.z());
+    bary[3] = 1 - cmptSum(res);
 
     return detT;
 }

src/OpenFOAM/meshes/primitiveShapes/triangle/triangle.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,6 +43,7 @@ SourceFiles
 #include "FixedList.H"
 #include "UList.H"
 #include "linePointRef.H"
+#include "barycentric2D.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -237,13 +238,18 @@ public:
             //  distribution
             inline Point randomPoint(Random& rndGen) const;
 
-            //- Calculate the barycentric coordinates of the given
-            //  point, in the same order as a, b, c.  Returns the
-            //  determinant of the solution.
-            inline scalar barycentric
+            //- Calculate the point from the given barycentric coordinates.
+            inline Point barycentricToPoint(const barycentric2D& bary) const;
+
+            //- Calculate the barycentric coordinates from the given point
+            inline barycentric2D pointToBarycentric(const point& pt) const;
+
+            //- Calculate the barycentric coordinates from the given point.
+            //  Returns the determinant.
+            inline scalar pointToBarycentric
             (
                 const point& pt,
-                FixedList<scalar, 3>& bary
+                barycentric2D& bary
             ) const;
 
             //- Return point intersection with a ray.

src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -244,24 +244,37 @@ inline Foam::tensor Foam::triangle<Point, PointRef>::inertia
 template<class Point, class PointRef>
 inline Point Foam::triangle<Point, PointRef>::randomPoint(Random& rndGen) const
 {
-    // Generating Random Points in Triangles
-    // by Greg Turk
-    // from "Graphics Gems", Academic Press, 1990
-    // http://tog.acm.org/GraphicsGems/gems/TriPoints.c
+    return barycentricToPoint(barycentric2D01(rndGen));
+}
 
-    scalar s = rndGen.sample01<scalar>();
 
-    scalar t = sqrt(rndGen.sample01<scalar>());
+template<class Point, class PointRef>
+inline Point Foam::triangle<Point, PointRef>::barycentricToPoint
+(
+    const barycentric2D& bary
+) const
+{
+    return bary[0]*a_ + bary[1]*b_ + bary[2]*c_;
+}
 
-    return (1 - t)*a_ + (1 - s)*t*b_ + s*t*c_;
+
+template<class Point, class PointRef>
+inline Foam::barycentric2D Foam::triangle<Point, PointRef>::pointToBarycentric
+(
+    const point& pt
+) const
+{
+    barycentric2D bary;
+    pointToBarycentric(pt, bary);
+    return bary;
 }
 
 
 template<class Point, class PointRef>
-Foam::scalar Foam::triangle<Point, PointRef>::barycentric
+inline Foam::scalar Foam::triangle<Point, PointRef>::pointToBarycentric
 (
     const point& pt,
-    FixedList<scalar, 3>& bary
+    barycentric2D& bary
 ) const
 {
     // Reference:
@@ -283,7 +296,7 @@ Foam::scalar Foam::triangle<Point, PointRef>::barycentric
     {
         // Degenerate triangle, returning 1/3 barycentric coordinates.
 
-        bary = FixedList<scalar, 3>(1.0/3.0);
+        bary = barycentric2D(1.0/3.0, 1.0/3.0, 1.0/3.0);
 
         return denom;
     }

src/OpenFOAM/primitives/Barycentric2D/Barycentric2D.H

+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::Barycentric2D
+
+Description
+    Templated 2D Barycentric derived from VectorSpace. Has 3 components, one of
+    which is redundant.
+
+SourceFiles
+    Barycentric2DI.H
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Barycentric2D_H
+#define Barycentric2D_H
+
+#include "VectorSpace.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                         Class Barycentric2D Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Cmpt>
+class Barycentric2D
+:
+    public VectorSpace<Barycentric2D<Cmpt>, Cmpt, 3>
+{
+public:
+
+    //- Equivalent type of labels used for valid component indexing
+    typedef Barycentric2D<label> labelType;
+
+
+    // Member constants
+
+        //- Rank of Barycentric2D is 1
+        static const direction rank = 1;
+
+
+    //- Component labeling enumeration
+    enum components { A, B, C };
+
+
+    // Constructors
+
+        //- Construct null
+        inline Barycentric2D();
+
+        //- Construct initialized to zero
+        inline Barycentric2D(const Foam::zero);
+
+        //- Construct given four components
+        inline Barycentric2D
+        (
+            const Cmpt& va,
+            const Cmpt& vb,
+            const Cmpt& vc
+        );
+
+
+    // Member Functions
+
+        // Access
+
+            inline const Cmpt& a() const;
+            inline const Cmpt& b() const;
+            inline const Cmpt& c() const;
+
+            inline Cmpt& a();
+            inline Cmpt& b();
+            inline Cmpt& c();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "Barycentric2DI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/OpenFOAM/primitives/Barycentric2D/Barycentric2DI.H

+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANB WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class Cmpt>
+inline Foam::Barycentric2D<Cmpt>::Barycentric2D()
+{}
+
+
+template<class Cmpt>
+inline Foam::Barycentric2D<Cmpt>::Barycentric2D(const Foam::zero)
+:
+    Barycentric2D::vsType(Zero)
+{}
+
+
+template<class Cmpt>
+inline Foam::Barycentric2D<Cmpt>::Barycentric2D
+(
+    const Cmpt& va,
+    const Cmpt& vb,
+    const Cmpt& vc
+)
+{
+    this->v_[A] = va;
+    this->v_[B] = vb;
+    this->v_[C] = vc;
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class Cmpt>
+inline const Cmpt& Foam::Barycentric2D<Cmpt>::a() const
+{
+    return this->v_[A];
+}
+
+
+template<class Cmpt>
+inline const Cmpt& Foam::Barycentric2D<Cmpt>::b() const
+{
+    return this->v_[B];
+}
+
+
+template<class Cmpt>
+inline const Cmpt& Foam::Barycentric2D<Cmpt>::c() const
+{
+    return this->v_[C];
+}
+
+
+template<class Cmpt>
+inline Cmpt& Foam::Barycentric2D<Cmpt>::a()
+{
+    return this->v_[A];
+}
+
+
+template<class Cmpt>
+inline Cmpt& Foam::Barycentric2D<Cmpt>::b()
+{
+    return this->v_[B];
+}
+
+
+template<class Cmpt>
+inline Cmpt& Foam::Barycentric2D<Cmpt>::c()
+{
+    return this->v_[C];
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * Global Operators  * * * * * * * * * * * * * //
+
+template<class Cmpt>
+inline Cmpt operator&
+(
+    const Barycentric2D<Cmpt>& b1,
+    const Barycentric2D<Cmpt>& b2
+)
+{
+    return b1.a()*b2.a() + b1.b()*b2.b() + b1.c()*b2.c();
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/OpenFOAM/primitives/Barycentric2D/barycentric2D/barycentric2D.C

+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License