5 files changed, 99 insertions, 1 deletions
diff --git a/Eigen/src/Core/arch/AVX/PacketMath.h b/Eigen/src/Core/arch/AVX/PacketMath.h
index e3363d006..2e5f5e5bc 100644
--- a/Eigen/src/Core/arch/AVX/PacketMath.h
+++ b/Eigen/src/Core/arch/AVX/PacketMath.h
@@ -73,6 +73,8 @@ template<> struct packet_traits<float>  : default_packet_traits
     HasExpm1  = 1,
     HasExp  = 1,
     HasNdtri = 1,
+    HasI0e  = 1,
+    HasI1e  = 1,
     HasSqrt = 1,
     HasRsqrt = 1,
     HasTanh  = EIGEN_FAST_MATH,
diff --git a/Eigen/src/Core/arch/AVX512/PacketMath.h b/Eigen/src/Core/arch/AVX512/PacketMath.h
index 11c8dae02..67e667640 100644
--- a/Eigen/src/Core/arch/AVX512/PacketMath.h
+++ b/Eigen/src/Core/arch/AVX512/PacketMath.h
@@ -99,6 +99,8 @@ template<> struct packet_traits<float>  : default_packet_traits
     HasExpm1  = 1,
     HasNdtri = 1,
 #endif
+    HasI0e  = 1,
+    HasI1e  = 1,
     HasExp = 1,
     HasSqrt = EIGEN_FAST_MATH,
     HasRsqrt = EIGEN_FAST_MATH,
diff --git a/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h b/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h
index 518db2207..0a4b66089 100644
--- a/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h
+++ b/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h
@@ -570,6 +570,97 @@ struct ppolevl<Packet, 0> {
   }
 };
 
+/* chbevl (modified for Eigen)
+ *
+ *     Evaluate Chebyshev series
+ *
+ *
+ *
+ * SYNOPSIS:
+ *
+ * int N;
+ * Scalar x, y, coef[N], chebevl();
+ *
+ * y = chbevl( x, coef, N );
+ *
+ *
+ *
+ * DESCRIPTION:
+ *
+ * Evaluates the series
+ *
+ *        N-1
+ *         - '
+ *  y  =   >   coef[i] T (x/2)
+ *         -            i
+ *        i=0
+ *
+ * of Chebyshev polynomials Ti at argument x/2.
+ *
+ * Coefficients are stored in reverse order, i.e. the zero
+ * order term is last in the array.  Note N is the number of
+ * coefficients, not the order.
+ *
+ * If coefficients are for the interval a to b, x must
+ * have been transformed to x -> 2(2x - b - a)/(b-a) before
+ * entering the routine.  This maps x from (a, b) to (-1, 1),
+ * over which the Chebyshev polynomials are defined.
+ *
+ * If the coefficients are for the inverted interval, in
+ * which (a, b) is mapped to (1/b, 1/a), the transformation
+ * required is x -> 2(2ab/x - b - a)/(b-a).  If b is infinity,
+ * this becomes x -> 4a/x - 1.
+ *
+ *
+ *
+ * SPEED:
+ *
+ * Taking advantage of the recurrence properties of the
+ * Chebyshev polynomials, the routine requires one more
+ * addition per loop than evaluating a nested polynomial of
+ * the same degree.
+ *
+ */
+template <typename Packet, int N>
+struct generic_cheb_recurrence {
+  EIGEN_DEVICE_FUNC
+  static EIGEN_STRONG_INLINE Packet run(Packet x, const typename unpacket_traits<Packet>::type coef[]) {
+    EIGEN_STATIC_ASSERT((N > 2), YOU_MADE_A_PROGRAMMING_MISTAKE);
+    return pmadd(
+        generic_cheb_recurrence<Packet, N - 1>::run(x, coef), x,
+        psub(pset1<Packet>(coef[N - 1]), generic_cheb_recurrence<Packet, N -
+             2>::run(x, coef)));
+  }
+};
+
+template <typename Packet>
+struct generic_cheb_recurrence<Packet, 2> {
+  EIGEN_DEVICE_FUNC
+  static EIGEN_STRONG_INLINE Packet run(Packet x, const typename unpacket_traits<Packet>::type coef[]) {
+    return pmadd(pset1<Packet>(coef[0]), x, pset1<Packet>(coef[1]));
+  }
+};
+
+template <typename Packet>
+struct generic_cheb_recurrence<Packet, 1> {
+  EIGEN_DEVICE_FUNC
+  static EIGEN_STRONG_INLINE Packet run(Packet x, const typename unpacket_traits<Packet>::type coef[]) {
+    EIGEN_UNUSED_VARIABLE(x);
+    return pset1<Packet>(coef[0]);
+  }
+};
+
+template <typename Packet, int N>
+struct pchebevl {
+  EIGEN_DEVICE_FUNC
+  static EIGEN_STRONG_INLINE Packet run(Packet x, const typename unpacket_traits<Packet>::type coef[]) {
+    const Packet half = pset1<Packet>(0.5);
+    return pmul(half, psub(
+        generic_cheb_recurrence<Packet, N>::run(x, coef),
+        generic_cheb_recurrence<Packet, N - 2>::run(x, coef)));
+  }
+};
+
 } // end namespace internal
 } // end namespace Eigen
 
diff --git a/Eigen/src/Core/arch/SSE/PacketMath.h b/Eigen/src/Core/arch/SSE/PacketMath.h
index ddd2979af..0aadefab7 100755
--- a/Eigen/src/Core/arch/SSE/PacketMath.h
+++ b/Eigen/src/Core/arch/SSE/PacketMath.h
@@ -114,7 +114,8 @@ template<> struct packet_traits<float>  : default_packet_traits
     HasExpm1  = 1,
     HasNdtri = 1,
     HasExp  = 1,
-    HasNdtri = 1,
+    HasI0e  = 1,
+    HasI1e  = 1,
     HasSqrt = 1,
     HasRsqrt = 1,
     HasTanh  = EIGEN_FAST_MATH,
diff --git a/Eigen/src/Core/util/ForwardDeclarations.h b/Eigen/src/Core/util/ForwardDeclarations.h
index 18889fcf4..749945031 100644
--- a/Eigen/src/Core/util/ForwardDeclarations.h
+++ b/Eigen/src/Core/util/ForwardDeclarations.h
@@ -215,6 +215,8 @@ template<typename Scalar> struct scalar_digamma_op;
 template<typename Scalar> struct scalar_erf_op;
 template<typename Scalar> struct scalar_erfc_op;
 template<typename Scalar> struct scalar_ndtri_op;
+template<typename Scalar> struct scalar_i0e_op;
+template<typename Scalar> struct scalar_i1e_op;
 template<typename Scalar> struct scalar_igamma_op;
 template<typename Scalar> struct scalar_igammac_op;
 template<typename Scalar> struct scalar_zeta_op;