Synchronized with upstream 5.0.1 release (via git tag in 5.0 branch)

Author: eddelbuettel

inst/include/Eigen/CholmodSupport

@@ -12,7 +12,7 @@
 
 #include "src/Core/util/DisableStupidWarnings.h"
 
-#include <RcppEigenCholmod.h>
+#include <cholmod.h>
 
 /** \ingroup Support_modules
  * \defgroup CholmodSupport_Module CholmodSupport module

inst/include/Eigen/Core

@@ -106,6 +106,11 @@
 #include <thread>
 #endif
 
+// for __cpp_lib feature test macros
+#if defined(__has_include) && __has_include(<version>)
+#include <version>
+#endif
+
 // for std::bit_cast()
 #if defined(__cpp_lib_bit_cast) && __cpp_lib_bit_cast >= 201806L
 #include <bit>

inst/include/Eigen/Version

@@ -6,9 +6,9 @@
 // As of Eigen3 5.0.0, we have moved to Semantic Versioning (semver.org).
 #define EIGEN_MAJOR_VERSION 5
 #define EIGEN_MINOR_VERSION 0
-#define EIGEN_PATCH_VERSION 0
+#define EIGEN_PATCH_VERSION 1
 #define EIGEN_PRERELEASE_VERSION ""
 #define EIGEN_BUILD_VERSION ""
-#define EIGEN_VERSION_STRING "5.0.0"
+#define EIGEN_VERSION_STRING "5.0.1"
 
 #endif  // EIGEN_VERSION_H

inst/include/Eigen/src/CholmodSupport/CholmodSupport.h

@@ -13,10 +13,6 @@
 // IWYU pragma: private
 #include "./InternalHeaderCheck.h"
 
-#ifndef R_MATRIX_CHOLMOD
-# define R_MATRIX_CHOLMOD(_NAME_) cholmod_ ## _NAME_
-#endif
-
 namespace Eigen {
 
 namespace internal {
@@ -88,8 +84,8 @@ cholmod_sparse viewAsCholmod(Ref<SparseMatrix<Scalar_, Options_, StorageIndex_>
 
   if (internal::is_same<StorageIndex_, int>::value) {
     res.itype = CHOLMOD_INT;
-  // } else if (internal::is_same<StorageIndex_, SuiteSparse_long>::value) {
-  //   res.itype = CHOLMOD_LONG;
+  } else if (internal::is_same<StorageIndex_, SuiteSparse_long>::value) {
+    res.itype = CHOLMOD_LONG;
   } else {
     eigen_assert(false && "Index type not supported yet");
   }
@@ -176,14 +172,22 @@ namespace internal {
 #define EIGEN_CHOLMOD_SPECIALIZE0(ret, name)                        \
   template <typename StorageIndex_>                                 \
   inline ret cm_##name(cholmod_common& Common) {                    \
-    return R_MATRIX_CHOLMOD(name)(&Common);			    \
-  }                                                                 
+    return cholmod_##name(&Common);                                 \
+  }                                                                 \
+  template <>                                                       \
+  inline ret cm_##name<SuiteSparse_long>(cholmod_common & Common) { \
+    return cholmod_l_##name(&Common);                               \
+  }
 
-#define EIGEN_CHOLMOD_SPECIALIZE1(ret, name, t1, a1)                    \
-  template <typename StorageIndex_>                                     \
-  inline ret cm_##name(t1& a1, cholmod_common& Common) {                \
-    return R_MATRIX_CHOLMOD(name) (&a1, &Common);			\
-  }                                                                          
+#define EIGEN_CHOLMOD_SPECIALIZE1(ret, name, t1, a1)                         \
+  template <typename StorageIndex_>                                          \
+  inline ret cm_##name(t1& a1, cholmod_common& Common) {                     \
+    return cholmod_##name(&a1, &Common);                                     \
+  }                                                                          \
+  template <>                                                                \
+  inline ret cm_##name<SuiteSparse_long>(t1 & a1, cholmod_common & Common) { \
+    return cholmod_l_##name(&a1, &Common);                                   \
+  }
 
 EIGEN_CHOLMOD_SPECIALIZE0(int, start)
 EIGEN_CHOLMOD_SPECIALIZE0(int, finish)
@@ -197,33 +201,33 @@ EIGEN_CHOLMOD_SPECIALIZE1(cholmod_sparse*, factor_to_sparse, cholmod_factor, L)
 
 template <typename StorageIndex_>
 inline cholmod_dense* cm_solve(int sys, cholmod_factor& L, cholmod_dense& B, cholmod_common& Common) {
-  return R_MATRIX_CHOLMOD(solve) (sys, &L, &B, &Common);
+  return cholmod_solve(sys, &L, &B, &Common);
+}
+template <>
+inline cholmod_dense* cm_solve<SuiteSparse_long>(int sys, cholmod_factor& L, cholmod_dense& B, cholmod_common& Common) {
+  return cholmod_l_solve(sys, &L, &B, &Common);
 }
-// template <>
-// inline cholmod_dense* cm_solve<SuiteSparse_long>(int sys, cholmod_factor& L, cholmod_dense& B, cholmod_common& Common) {
-//   return cholmod_l_solve(sys, &L, &B, &Common);
-// }
 
 template <typename StorageIndex_>
 inline cholmod_sparse* cm_spsolve(int sys, cholmod_factor& L, cholmod_sparse& B, cholmod_common& Common) {
-  return R_MATRIX_CHOLMOD(spsolve) (sys, &L, &B, &Common);
+  return cholmod_spsolve(sys, &L, &B, &Common);
+}
+template <>
+inline cholmod_sparse* cm_spsolve<SuiteSparse_long>(int sys, cholmod_factor& L, cholmod_sparse& B,
+                                                    cholmod_common& Common) {
+  return cholmod_l_spsolve(sys, &L, &B, &Common);
 }
-// template <>
-// inline cholmod_sparse* cm_spsolve<SuiteSparse_long>(int sys, cholmod_factor& L, cholmod_sparse& B,
-//                                                     cholmod_common& Common) {
-//   return cholmod_l_spsolve(sys, &L, &B, &Common);
-// }
 
 template <typename StorageIndex_>
 inline int cm_factorize_p(cholmod_sparse* A, double beta[2], StorageIndex_* fset, std::size_t fsize, cholmod_factor* L,
                           cholmod_common& Common) {
-  return R_MATRIX_CHOLMOD(factorize_p) (A, beta, fset, fsize, L, &Common);
+  return cholmod_factorize_p(A, beta, fset, fsize, L, &Common);
+}
+template <>
+inline int cm_factorize_p<SuiteSparse_long>(cholmod_sparse* A, double beta[2], SuiteSparse_long* fset,
+                                            std::size_t fsize, cholmod_factor* L, cholmod_common& Common) {
+  return cholmod_l_factorize_p(A, beta, fset, fsize, L, &Common);
 }
-// template <>
-// inline int cm_factorize_p<SuiteSparse_long>(cholmod_sparse* A, double beta[2], SuiteSparse_long* fset,
-//                                             std::size_t fsize, cholmod_factor* L, cholmod_common& Common) {
-//   return cholmod_l_factorize_p(A, beta, fset, fsize, L, &Common);
-// }
 
 #undef EIGEN_CHOLMOD_SPECIALIZE0
 #undef EIGEN_CHOLMOD_SPECIALIZE1

inst/include/Eigen/src/Core/Assign_MKL.h

@@ -56,11 +56,11 @@ class vml_assign_traits {
                                                    : int(Dst::MaxRowsAtCompileTime),
     MaxSizeAtCompileTime = Dst::SizeAtCompileTime,
 
-    MightEnableVml = StorageOrdersAgree && DstHasDirectAccess && SrcHasDirectAccess &&
+    MightEnableVml = bool(StorageOrdersAgree) && bool(DstHasDirectAccess) && bool(SrcHasDirectAccess) &&
                      Src::InnerStrideAtCompileTime == 1 && Dst::InnerStrideAtCompileTime == 1,
-    MightLinearize = MightEnableVml && (int(Dst::Flags) & int(Src::Flags) & LinearAccessBit),
-    VmlSize = MightLinearize ? MaxSizeAtCompileTime : InnerMaxSize,
-    LargeEnough = VmlSize == Dynamic || VmlSize >= EIGEN_MKL_VML_THRESHOLD
+    MightLinearize = bool(MightEnableVml) && (int(Dst::Flags) & int(Src::Flags) & LinearAccessBit),
+    VmlSize = bool(MightLinearize) ? MaxSizeAtCompileTime : InnerMaxSize,
+    LargeEnough = (VmlSize == Dynamic) || VmlSize >= EIGEN_MKL_VML_THRESHOLD
   };
 
  public:

inst/include/Eigen/src/Core/CwiseNullaryOp.h

@@ -94,7 +94,7 @@ class CwiseNullaryOp : public internal::dense_xpr_base<CwiseNullaryOp<NullaryOp,
  * the returned matrix. Must be compatible with this MatrixBase type.
  *
  * This variant is meant to be used for dynamic-size matrix types. For fixed-size types,
- * it is redundant to pass \a rows and \a cols as arguments, so Zero() should be used
+ * it is redundant to pass \a rows and \a cols as arguments, so NullaryExpr(const CustomNullaryOp&) should be used
  * instead.
  *
  * The template parameter \a CustomNullaryOp is the type of the functor.
@@ -121,7 +121,7 @@ EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
  * \only_for_vectors
  *
  * This variant is meant to be used for dynamic-size vector types. For fixed-size types,
- * it is redundant to pass \a size as argument, so Zero() should be used
+ * it is redundant to pass \a size as argument, so NullaryExpr(const CustomNullaryOp&) should be used
  * instead.
  *
  * The template parameter \a CustomNullaryOp is the type of the functor.
@@ -174,7 +174,7 @@ EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
  * the returned matrix. Must be compatible with this DenseBase type.
  *
  * This variant is meant to be used for dynamic-size matrix types. For fixed-size types,
- * it is redundant to pass \a rows and \a cols as arguments, so Zero() should be used
+ * it is redundant to pass \a rows and \a cols as arguments, so Constant(const Scalar&) should be used
  * instead.
  *
  * The template parameter \a CustomNullaryOp is the type of the functor.
@@ -195,7 +195,7 @@ DenseBase<Derived>::Constant(Index rows, Index cols, const Scalar& value) {
  * \only_for_vectors
  *
  * This variant is meant to be used for dynamic-size vector types. For fixed-size types,
- * it is redundant to pass \a size as argument, so Zero() should be used
+ * it is redundant to pass \a size as argument, so Constant(const Scalar&) should be used
  * instead.
  *
  * The template parameter \a CustomNullaryOp is the type of the functor.

inst/include/Eigen/src/Core/DenseStorage.h

@@ -65,20 +65,15 @@ struct plain_array {
 
 template <typename T, int Size, int MatrixOrArrayOptions>
 struct plain_array<T, Size, MatrixOrArrayOptions, 0> {
-  T array[Size];
+  // on some 32-bit platforms, stack-allocated arrays are aligned to 4 bytes, not the preferred alignment of T
+  EIGEN_ALIGN_TO_BOUNDARY(alignof(T)) T array[Size];
 #if defined(EIGEN_NO_DEBUG) || defined(EIGEN_TESTING_PLAINOBJECT_CTOR)
   EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE constexpr plain_array() = default;
 #else
   EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE constexpr plain_array() { EIGEN_MAKE_STACK_ALLOCATION_ASSERT(Size * sizeof(T)) }
 #endif
 };
 
-template <typename T, int MatrixOrArrayOptions, int Alignment>
-struct plain_array<T, 0, MatrixOrArrayOptions, Alignment> {
-  T array[1];
-  EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE constexpr plain_array() = default;
-};
-
 template <typename T, int Size, int Options, int Alignment>
 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE constexpr void swap_plain_array(plain_array<T, Size, Options, Alignment>& a,
                                                                       plain_array<T, Size, Options, Alignment>& b,

inst/include/Eigen/src/Core/Fill.h

@@ -115,17 +115,15 @@ struct eigen_zero_impl<Xpr, /*use_memset*/ false> {
 template <typename Xpr>
 struct eigen_zero_impl<Xpr, /*use_memset*/ true> {
   using Scalar = typename Xpr::Scalar;
-  static constexpr size_t max_bytes = (std::numeric_limits<std::ptrdiff_t>::max)();
   static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Xpr& dst) {
-    const size_t num_bytes = dst.size() * sizeof(Scalar);
-    if (num_bytes == 0) return;
+    const std::ptrdiff_t num_bytes = dst.size() * static_cast<std::ptrdiff_t>(sizeof(Scalar));
+    if (num_bytes <= 0) return;
     void* dst_ptr = static_cast<void*>(dst.data());
 #ifndef EIGEN_NO_DEBUG
-    if (num_bytes > max_bytes) throw_std_bad_alloc();
     eigen_assert((dst_ptr != nullptr) && "null pointer dereference error!");
 #endif
     EIGEN_USING_STD(memset);
-    memset(dst_ptr, 0, num_bytes);
+    memset(dst_ptr, 0, static_cast<std::size_t>(num_bytes));
   }
   template <typename SrcXpr>
   static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Xpr& dst, const SrcXpr& src) {

inst/include/Eigen/src/Core/InnerProduct.h

@@ -211,8 +211,14 @@ struct scalar_inner_product_op {
   static constexpr bool PacketAccess = false;
 };
 
+// Partial specialization for packet access if and only if
+// LhsScalar == RhsScalar == ScalarBinaryOpTraits<LhsScalar, RhsScalar>::ReturnType.
 template <typename Scalar, bool Conj>
-struct scalar_inner_product_op<Scalar, Scalar, Conj> {
+struct scalar_inner_product_op<
+    Scalar,
+    typename std::enable_if<internal::is_same<typename ScalarBinaryOpTraits<Scalar, Scalar>::ReturnType, Scalar>::value,
+                            Scalar>::type,
+    Conj> {
   using result_type = Scalar;
   using conj_helper = conditional_conj<Scalar, Conj>;
   EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar coeff(const Scalar& a, const Scalar& b) const {

inst/include/Eigen/src/Core/MathFunctions.h

@@ -1004,8 +1004,7 @@ struct madd_impl {
   }
 };
 
-// Use FMA if there is a single CPU instruction.
-#ifdef EIGEN_VECTORIZE_FMA
+#if EIGEN_SCALAR_MADD_USE_FMA
 template <typename Scalar>
 struct madd_impl<Scalar, std::enable_if_t<has_fma<Scalar>::value>> {
   static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar run(const Scalar& x, const Scalar& y, const Scalar& z) {
@@ -1927,7 +1926,6 @@ EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar arithmetic_shift_right(const Scalar
   return bit_cast<Scalar, SignedScalar>(bit_cast<SignedScalar, Scalar>(a) >> n);
 }
 
-// Otherwise, rely on template implementation.
 template <typename Scalar>
 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar fma(const Scalar& x, const Scalar& y, const Scalar& z) {
   return internal::fma_impl<Scalar>::run(x, y, z);