Adds a curved polygon class to primal

.mailmap

@@ -15,10 +15,12 @@ Benjamin Curtice Corbett <corbett5@llnl.gov>    Benjamin Corbett <corbett5@llnl.
 Benjamin Curtice Corbett <corbett5@llnl.gov>    Ben Corbett <32752943+corbett5@users.noreply.github.com>
 Brian Manh Hien Han <han12@llnl.gov>            Brian Han <han12@llnl.gov>          
 Brian Manh Hien Han <han12@llnl.gov>            Brian Manh Hien Han <han12@rzansel19.coral.llnl.gov>
+Brian T.N. Gunney <gunney1@llnl.gov>            Brian Gunney <45609916+gunney1@users.noreply.github.com>
 Chris White <white238@llnl.gov>                 Christopher A. White <white238@llnl.gov>
 Cyrus D. Harrison <harrison37@llnl.gov>         Cyrus Harrison <cyrush@llnl.gov>
 Cyrus D. Harrison <harrison37@llnl.gov>         Cyrus Harrison <harrison37@llnl.gov>
 Cyrus D. Harrison <harrison37@llnl.gov>         Cyrus <cyrush@llnl.gov>
+Daniel Taller <taller1@llnl.gov>                Danny Taller <66029857+dtaller@users.noreply.github.com>
 Esteban Pauli <pauli2@llnl.gov>                 Esteban Pauli <40901502+estebanpauli@users.noreply.github.com>
 Evan Taylor Desantola <desantola1@llnl.gov>     Evan Taylor DeSantola <desantola1@llnl.gov>
 George Zagaris <zagaris2@llnl.gov>              George Zagaris <george.zagaris@gmail.com>

RELEASE-NOTES.md

@@ -88,6 +88,10 @@ The Axom project release numbers follow [Semantic Versioning](http://semver.org/
 - Adds an overload to quest's `SignedDistance` query to return the closest point on the surface 
   to the query point and the surface normal at that point. Also exposes this functionality 
   in quest's signed_distance C API.
+- Adds utility function for linear interpolation (`lerp`) of two numbers
+- Adds utility function to compute binomial coefficients  
+- Adds a `CurvedPolygon` class to primal representing a polygon with `BezierCurves` as edges
+- Adds functions to compute the moments (area, centroid) of a `CurvedPolygon`
 
 ###  Changed
 - Axom now requires C++14 and will default to that if not specified via `BLT_CXX_STD`.
@@ -140,6 +144,7 @@ The Axom project release numbers follow [Semantic Versioning](http://semver.org/
   and added a new derived class `sidre::IndexedCollection`
 - Spin: `BVH::findPoints/Rays/BoundingBoxes()` candidate search methods now accept an `axom::ArrayView<IndexType>`
   for the `offsets` and `counts` output arrays, and return `candidates` as an `axom::Array<IndexType>`.
+- Renamed `primal::Polygon::centroid()` to `primal::Polygon::vertexMean()` because it was not actually computing the centroid.
 
 ###  Fixed
 - Fixed a bug relating to swap and assignment operations for multidimensional `axom::Array`s

src/axom/core/numerics/Matrix.hpp

@@ -117,6 +117,11 @@ template <typename T>
 class Matrix
 {
 public:
+  /*!
+   * \brief Default constructor
+   */
+  Matrix() : m_rows(0), m_cols(0), m_data(nullptr), m_usingExternal(false) { }
+
   /*!
    * \brief Constructor, creates a Matrix with the given rows and columns.
    *
@@ -527,12 +532,6 @@ class Matrix
   /// @}
 
 private:
-  /*!
-   * \brief Default constructor. Does nothing.
-   * \note Made private to prevent host-code from calling this.
-   */
-  Matrix() : m_rows(0), m_cols(0), m_data(nullptr) {};
-
   /// \name Private Helper Methods
   /// @{
 

src/axom/core/tests/utils_utilities.hpp

@@ -158,3 +158,139 @@ TEST(utils_utilities, floor_ceil)
     EXPECT_EQ(-5.0, axom::utilities::ceil(val));
   }
 }
+
+//------------------------------------------------------------------------------
+TEST(utils_utilities, binomial_coefficient)
+{
+  std::cout << "Testing binomial coefficient function." << std::endl;
+
+  // test n less than zero
+  {
+    const int n = -1;
+    const int exp = 0;
+    for(int k = -1; k < 10; ++k)
+    {
+      auto binom_k_n = axom::utilities::binomialCoefficient(n, k);
+      EXPECT_EQ(exp, binom_k_n);
+    }
+  }
+
+  // test n := 0
+  {
+    const int n = 0;
+
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 0));
+
+    EXPECT_EQ(0, axom::utilities::binomialCoefficient(n, -1));
+    EXPECT_EQ(0, axom::utilities::binomialCoefficient(n, 1));
+  }
+
+  // test n := 1
+  {
+    const int n = 1;
+
+    EXPECT_EQ(0, axom::utilities::binomialCoefficient(n, -1));
+
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 0));
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 1));
+
+    EXPECT_EQ(0, axom::utilities::binomialCoefficient(n, 2));
+  }
+
+  // test n := 2
+  {
+    const int n = 2;
+
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 0));
+    EXPECT_EQ(2, axom::utilities::binomialCoefficient(n, 1));
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 2));
+  }
+
+  // test n := 3
+  {
+    const int n = 3;
+
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 0));
+    EXPECT_EQ(3, axom::utilities::binomialCoefficient(n, 1));
+    EXPECT_EQ(3, axom::utilities::binomialCoefficient(n, 2));
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 3));
+  }
+
+  // test n := 4
+  {
+    const int n = 4;
+
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 0));
+    EXPECT_EQ(4, axom::utilities::binomialCoefficient(n, 1));
+    EXPECT_EQ(6, axom::utilities::binomialCoefficient(n, 2));
+    EXPECT_EQ(4, axom::utilities::binomialCoefficient(n, 3));
+    EXPECT_EQ(1, axom::utilities::binomialCoefficient(n, 4));
+  }
+
+  // test recurrence relation  nCk = (n-1)C(k-1) + (n-1)C(k)
+  {
+    for(int n = 1; n < 10; ++n)
+    {
+      for(int k = 1; k <= n; ++k)
+      {
+        auto binom_n_k = axom::utilities::binomialCoefficient(n, k);
+        auto binom_n1_k1 = axom::utilities::binomialCoefficient(n - 1, k - 1);
+        auto binom_n1_k = axom::utilities::binomialCoefficient(n - 1, k);
+
+        EXPECT_EQ(binom_n_k, binom_n1_k1 + binom_n1_k);
+      }
+    }
+  }
+}
+
+TEST(utils_utilities, lerp)
+{
+  std::cout << "Testing lerp function." << std::endl;
+
+  {
+    double A = 0.;
+    double B = 1.;
+
+    EXPECT_DOUBLE_EQ(0., axom::utilities::lerp(A, B, 0.));
+    EXPECT_DOUBLE_EQ(1., axom::utilities::lerp(A, B, 1.));
+    EXPECT_DOUBLE_EQ(.5, axom::utilities::lerp(A, B, .5));
+
+    for(double t = -2.0; t < 2.0; t += 0.05)
+    {
+      EXPECT_DOUBLE_EQ(t, axom::utilities::lerp(A, B, t));
+    }
+  }
+
+  // Test endpoint interpolation
+  {
+    const double lower = -.23;  // Arbitrary end points
+    const double upper = 5.73;
+    for(int i = 0; i < 100; ++i)
+    {
+      double A = axom::utilities::random_real(lower, upper);
+      double B = axom::utilities::random_real(lower, upper);
+
+      EXPECT_DOUBLE_EQ(A, axom::utilities::lerp(A, B, 0.));
+      EXPECT_DOUBLE_EQ(B, axom::utilities::lerp(B, A, 0.));
+      EXPECT_DOUBLE_EQ(B, axom::utilities::lerp(A, B, 1.));
+      EXPECT_DOUBLE_EQ(A, axom::utilities::lerp(B, A, 1.));
+    }
+  }
+
+  // Compute using different form
+  {
+    const double lower = -.23;  // Arbitrary end points
+    const double upper = 5.73;
+    for(int i = 0; i < 100; ++i)
+    {
+      double A = axom::utilities::random_real(lower, upper);
+      double B = axom::utilities::random_real(lower, upper);
+
+      // Test interpolation and also extrapolation beyond endpoints.
+      double t = axom::utilities::random_real(-1.5, 1.5);
+
+      double exp = A + (B - A) * t;
+      EXPECT_NEAR(exp, axom::utilities::lerp(A, B, t), 1e-12);
+    }
+  }
+}

src/axom/core/utilities/Utilities.cpp

@@ -4,8 +4,7 @@
 // SPDX-License-Identifier: (BSD-3-Clause)
 
 /*!
- *
- * \file
+ * \file    Utilities.cpp
  *
  * \brief   Implementation file for utility functions.
  *
@@ -39,5 +38,29 @@ void processAbort()
 #endif
 }
 
+int binomialCoefficient(int n, int k)
+{
+  if(k > n || k < 0)  // check if out-of-bounds
+  {
+    return 0;
+  }
+  if(k == n || k == 0)  // early return
+  {
+    return 1;
+  }
+  if(k > n - k)  // exploit symmetry to reduce work
+  {
+    k = n - k;
+  }
+
+  int val = 1;
+  for(int i = 1; i <= k; ++i)
+  {
+    val *= (n - k + i);
+    val /= i;
+  }
+  return val;
+}
+
 }  // end namespace utilities
 }  // end namespace axom

src/axom/core/utilities/Utilities.hpp

@@ -109,6 +109,15 @@ inline AXOM_HOST_DEVICE void swap(T& a, T& b)
   b = tmp;
 }
 
+/*! 
+ * \brief returns the linear interpolation of \a A and \a B at \a t. i.e. (1-t)A+tB
+ */
+template <typename T>
+inline AXOM_HOST_DEVICE T lerp(T A, T B, T t)
+{
+  return (1 - t) * A + t * B;
+}
+
 /*!
  * \brief Returns the base 2 logarithm of the input.
  * \param [in] val The input value
@@ -163,6 +172,14 @@ inline AXOM_HOST_DEVICE T clampLower(T val, T lower)
   return val < lower ? lower : val;
 }
 
+/*!
+ * \brief Computes the binomial coefficient `n choose k`
+ *
+ * \return \f$ {n\choose k}  = n! / (k! * (n-k)!)\f$
+ * when \f$ n \ge k \ge 0 \f$, 0 otherwise.
+ */
+int binomialCoefficient(int n, int k);
+
 /*!
  * \brief Returns a random real number within the specified interval
  *

src/axom/primal/CMakeLists.txt

@@ -20,6 +20,7 @@ set( primal_headers
      ## geometry
      geometry/BezierCurve.hpp
      geometry/BoundingBox.hpp
+     geometry/CurvedPolygon.hpp
      geometry/OrientedBoundingBox.hpp
      geometry/OrientationResult.hpp
      geometry/NumericArray.hpp
@@ -42,14 +43,16 @@ set( primal_headers
      operators/orientation.hpp
      operators/squared_distance.hpp
      operators/compute_bounding_box.hpp
+     operators/compute_moments.hpp
      operators/in_sphere.hpp
      operators/split.hpp
 
      operators/detail/clip_impl.hpp
+     operators/detail/compute_moments_impl.hpp
      operators/detail/intersect_bezier_impl.hpp
-     operators/detail/intersect_ray_impl.hpp
      operators/detail/intersect_bounding_box_impl.hpp
      operators/detail/intersect_impl.hpp
+     operators/detail/intersect_ray_impl.hpp
 
      ## utils
      utils/ZipIndexable.hpp