aboutsummaryrefslogtreecommitdiffhomepage
path: root/Eigen/src
diff options
context:
space:
mode:
authorGravatar Srinivas Vasudevan <srvasude@google.com>2019-09-03 15:34:47 -0400
committerGravatar Srinivas Vasudevan <srvasude@google.com>2019-09-03 15:34:47 -0400
commit99036a3615a57315564ab86f1d8754bc6d77c8f3 (patch)
treeef0a22c09ac900224ce2243561b019c66752f372 /Eigen/src
parent18ceb3413d09afc4f143014f89552f941321209b (diff)
parenta8d264fa9c56e42f77e2129d4e504f5c854821c2 (diff)
Merging from eigen/eigen.
Diffstat (limited to 'Eigen/src')
-rw-r--r--Eigen/src/Core/MathFunctions.h8
-rw-r--r--Eigen/src/Core/arch/AVX/MathFunctions.h10
-rw-r--r--Eigen/src/Core/arch/AVX/PacketMath.h364
-rw-r--r--Eigen/src/Core/arch/AVX/TypeCasting.h30
-rw-r--r--Eigen/src/Core/arch/AVX512/MathFunctions.h12
-rw-r--r--Eigen/src/Core/arch/AVX512/PacketMath.h503
-rw-r--r--Eigen/src/Core/arch/AVX512/TypeCasting.h47
-rwxr-xr-xEigen/src/Core/arch/AltiVec/PacketMath.h79
-rw-r--r--Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h45
-rw-r--r--Eigen/src/Core/arch/Default/Half.h (renamed from Eigen/src/Core/arch/GPU/Half.h)18
-rw-r--r--Eigen/src/Core/arch/Default/TypeCasting.h77
-rw-r--r--Eigen/src/Core/arch/GPU/PacketMath.h540
-rw-r--r--Eigen/src/Core/arch/GPU/PacketMathHalf.h1630
-rw-r--r--Eigen/src/Core/arch/GPU/TypeCasting.h161
-rw-r--r--Eigen/src/Core/arch/SSE/MathFunctions.h13
-rwxr-xr-xEigen/src/Core/arch/SSE/PacketMath.h212
-rw-r--r--Eigen/src/Core/arch/SSE/TypeCasting.h51
-rwxr-xr-xEigen/src/Core/util/DisableStupidWarnings.h5
-rw-r--r--Eigen/src/OrderingMethods/Eigen_Colamd.h568
-rw-r--r--Eigen/src/OrderingMethods/Ordering.h10
-rw-r--r--Eigen/src/PardisoSupport/PardisoSupport.h13
21 files changed, 2255 insertions, 2141 deletions
diff --git a/Eigen/src/Core/MathFunctions.h b/Eigen/src/Core/MathFunctions.h
index 8bef59354..fbec39d83 100644
--- a/Eigen/src/Core/MathFunctions.h
+++ b/Eigen/src/Core/MathFunctions.h
@@ -501,7 +501,8 @@ namespace std_fallback {
}
EIGEN_USING_STD_MATH(log);
- return (u - RealScalar(1)) * x / log(u);
+ Scalar logu = log(u);
+ return numext::equal_strict(u, logu) ? u : (u - RealScalar(1)) * x / logu;
}
}
@@ -548,7 +549,10 @@ namespace std_fallback {
typedef typename NumTraits<Scalar>::Real RealScalar;
EIGEN_USING_STD_MATH(log);
Scalar x1p = RealScalar(1) + x;
- return numext::equal_strict(x1p, Scalar(1)) ? x : x * ( log(x1p) / (x1p - RealScalar(1)) );
+ Scalar log_1p = log(x1p);
+ const bool is_small = numext::equal_strict(x1p, Scalar(1));
+ const bool is_inf = numext::equal_strict(x1p, log_1p);
+ return (is_small || is_inf) ? x : x * (log_1p / (x1p - RealScalar(1)));
}
}
diff --git a/Eigen/src/Core/arch/AVX/MathFunctions.h b/Eigen/src/Core/arch/AVX/MathFunctions.h
index 9f375ed98..c6d3cf6a0 100644
--- a/Eigen/src/Core/arch/AVX/MathFunctions.h
+++ b/Eigen/src/Core/arch/AVX/MathFunctions.h
@@ -36,6 +36,16 @@ plog<Packet8f>(const Packet8f& _x) {
return plog_float(_x);
}
+template<> EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED
+Packet8f plog1p<Packet8f>(const Packet8f& _x) {
+ return generic_plog1p(_x);
+}
+
+template<> EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED
+Packet8f pexpm1<Packet8f>(const Packet8f& _x) {
+ return generic_expm1(_x);
+}
+
// Exponential function. Works by writing "x = m*log(2) + r" where
// "m = floor(x/log(2)+1/2)" and "r" is the remainder. The result is then
// "exp(x) = 2^m*exp(r)" where exp(r) is in the range [-1,1).
diff --git a/Eigen/src/Core/arch/AVX/PacketMath.h b/Eigen/src/Core/arch/AVX/PacketMath.h
index 9feb96f8b..e3363d006 100644
--- a/Eigen/src/Core/arch/AVX/PacketMath.h
+++ b/Eigen/src/Core/arch/AVX/PacketMath.h
@@ -31,10 +31,14 @@ namespace internal {
typedef __m256 Packet8f;
typedef __m256i Packet8i;
typedef __m256d Packet4d;
+typedef struct {
+ __m128i x;
+} Packet8h;
template<> struct is_arithmetic<__m256> { enum { value = true }; };
template<> struct is_arithmetic<__m256i> { enum { value = true }; };
template<> struct is_arithmetic<__m256d> { enum { value = true }; };
+template<> struct is_arithmetic<Packet8h> { enum { value = true }; };
#define _EIGEN_DECLARE_CONST_Packet8f(NAME,X) \
const Packet8f p8f_##NAME = pset1<Packet8f>(X)
@@ -65,6 +69,8 @@ template<> struct packet_traits<float> : default_packet_traits
HasSin = EIGEN_FAST_MATH,
HasCos = EIGEN_FAST_MATH,
HasLog = 1,
+ HasLog1p = 1,
+ HasExpm1 = 1,
HasExp = 1,
HasNdtri = 1,
HasSqrt = 1,
@@ -96,6 +102,35 @@ template<> struct packet_traits<double> : default_packet_traits
HasCeil = 1
};
};
+
+template <>
+struct packet_traits<Eigen::half> : default_packet_traits {
+ typedef Packet8h type;
+ // There is no half-size packet for Packet8h.
+ typedef Packet8h half;
+ enum {
+ Vectorizable = 1,
+ AlignedOnScalar = 1,
+ size = 8,
+ HasHalfPacket = 0,
+ HasAdd = 1,
+ HasSub = 1,
+ HasMul = 1,
+ HasDiv = 1,
+ HasNegate = 1,
+ HasAbs = 0,
+ HasAbs2 = 0,
+ HasMin = 0,
+ HasMax = 0,
+ HasConj = 0,
+ HasSetLinear = 0,
+ HasSqrt = 0,
+ HasRsqrt = 0,
+ HasExp = 0,
+ HasLog = 0,
+ HasBlend = 0
+ };
+};
#endif
template<> struct scalar_div_cost<float,true> { enum { value = 14 }; };
@@ -846,6 +881,335 @@ template<> EIGEN_STRONG_INLINE Packet4d pinsertlast(const Packet4d& a, double b)
return _mm256_blend_pd(a,pset1<Packet4d>(b),(1<<3));
}
+
+// Packet math for Eigen::half
+template<> struct unpacket_traits<Packet8h> { typedef Eigen::half type; enum {size=8, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef Packet8h half; };
+
+template<> EIGEN_STRONG_INLINE Packet8h pset1<Packet8h>(const Eigen::half& from) {
+ Packet8h result;
+ result.x = _mm_set1_epi16(from.x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Eigen::half pfirst<Packet8h>(const Packet8h& from) {
+ return half_impl::raw_uint16_to_half(static_cast<unsigned short>(_mm_extract_epi16(from.x, 0)));
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pload<Packet8h>(const Eigen::half* from) {
+ Packet8h result;
+ result.x = _mm_load_si128(reinterpret_cast<const __m128i*>(from));
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h ploadu<Packet8h>(const Eigen::half* from) {
+ Packet8h result;
+ result.x = _mm_loadu_si128(reinterpret_cast<const __m128i*>(from));
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE void pstore<Eigen::half>(Eigen::half* to, const Packet8h& from) {
+ _mm_store_si128(reinterpret_cast<__m128i*>(to), from.x);
+}
+
+template<> EIGEN_STRONG_INLINE void pstoreu<Eigen::half>(Eigen::half* to, const Packet8h& from) {
+ _mm_storeu_si128(reinterpret_cast<__m128i*>(to), from.x);
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h
+ploaddup<Packet8h>(const Eigen::half* from) {
+ Packet8h result;
+ unsigned short a = from[0].x;
+ unsigned short b = from[1].x;
+ unsigned short c = from[2].x;
+ unsigned short d = from[3].x;
+ result.x = _mm_set_epi16(d, d, c, c, b, b, a, a);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h
+ploadquad<Packet8h>(const Eigen::half* from) {
+ Packet8h result;
+ unsigned short a = from[0].x;
+ unsigned short b = from[1].x;
+ result.x = _mm_set_epi16(b, b, b, b, a, a, a, a);
+ return result;
+}
+
+EIGEN_STRONG_INLINE Packet8f half2float(const Packet8h& a) {
+#ifdef EIGEN_HAS_FP16_C
+ return _mm256_cvtph_ps(a.x);
+#else
+ EIGEN_ALIGN32 Eigen::half aux[8];
+ pstore(aux, a);
+ float f0(aux[0]);
+ float f1(aux[1]);
+ float f2(aux[2]);
+ float f3(aux[3]);
+ float f4(aux[4]);
+ float f5(aux[5]);
+ float f6(aux[6]);
+ float f7(aux[7]);
+
+ return _mm256_set_ps(f7, f6, f5, f4, f3, f2, f1, f0);
+#endif
+}
+
+EIGEN_STRONG_INLINE Packet8h float2half(const Packet8f& a) {
+#ifdef EIGEN_HAS_FP16_C
+ Packet8h result;
+ result.x = _mm256_cvtps_ph(a, _MM_FROUND_TO_NEAREST_INT|_MM_FROUND_NO_EXC);
+ return result;
+#else
+ EIGEN_ALIGN32 float aux[8];
+ pstore(aux, a);
+ Eigen::half h0(aux[0]);
+ Eigen::half h1(aux[1]);
+ Eigen::half h2(aux[2]);
+ Eigen::half h3(aux[3]);
+ Eigen::half h4(aux[4]);
+ Eigen::half h5(aux[5]);
+ Eigen::half h6(aux[6]);
+ Eigen::half h7(aux[7]);
+
+ Packet8h result;
+ result.x = _mm_set_epi16(h7.x, h6.x, h5.x, h4.x, h3.x, h2.x, h1.x, h0.x);
+ return result;
+#endif
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h ptrue(const Packet8h& a) {
+ Packet8h r; r.x = _mm_cmpeq_epi32(a.x, a.x); return r;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h por(const Packet8h& a,const Packet8h& b) {
+ // in some cases Packet4i is a wrapper around __m128i, so we either need to
+ // cast to Packet4i to directly call the intrinsics as below:
+ Packet8h r; r.x = _mm_or_si128(a.x,b.x); return r;
+}
+template<> EIGEN_STRONG_INLINE Packet8h pxor(const Packet8h& a,const Packet8h& b) {
+ Packet8h r; r.x = _mm_xor_si128(a.x,b.x); return r;
+}
+template<> EIGEN_STRONG_INLINE Packet8h pand(const Packet8h& a,const Packet8h& b) {
+ Packet8h r; r.x = _mm_and_si128(a.x,b.x); return r;
+}
+template<> EIGEN_STRONG_INLINE Packet8h pandnot(const Packet8h& a,const Packet8h& b) {
+ Packet8h r; r.x = _mm_andnot_si128(b.x,a.x); return r;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pselect(const Packet8h& mask, const Packet8h& a, const Packet8h& b) {
+ Packet8h r; r.x = _mm_blendv_epi8(b.x, a.x, mask.x); return r;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pcmp_eq(const Packet8h& a,const Packet8h& b) {
+ Packet8f af = half2float(a);
+ Packet8f bf = half2float(b);
+ Packet8f rf = pcmp_eq(af, bf);
+ // Pack the 32-bit flags into 16-bits flags.
+ Packet8h result; result.x = _mm_packs_epi32(_mm256_extractf128_si256(_mm256_castps_si256(rf), 0),
+ _mm256_extractf128_si256(_mm256_castps_si256(rf), 1));
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pconj(const Packet8h& a) { return a; }
+
+template<> EIGEN_STRONG_INLINE Packet8h pnegate(const Packet8h& a) {
+ Packet8h sign_mask; sign_mask.x = _mm_set1_epi16(static_cast<unsigned short>(0x8000));
+ Packet8h result; result.x = _mm_xor_si128(a.x, sign_mask.x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h padd<Packet8h>(const Packet8h& a, const Packet8h& b) {
+ Packet8f af = half2float(a);
+ Packet8f bf = half2float(b);
+ Packet8f rf = padd(af, bf);
+ return float2half(rf);
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h psub<Packet8h>(const Packet8h& a, const Packet8h& b) {
+ Packet8f af = half2float(a);
+ Packet8f bf = half2float(b);
+ Packet8f rf = psub(af, bf);
+ return float2half(rf);
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pmul<Packet8h>(const Packet8h& a, const Packet8h& b) {
+ Packet8f af = half2float(a);
+ Packet8f bf = half2float(b);
+ Packet8f rf = pmul(af, bf);
+ return float2half(rf);
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pdiv<Packet8h>(const Packet8h& a, const Packet8h& b) {
+ Packet8f af = half2float(a);
+ Packet8f bf = half2float(b);
+ Packet8f rf = pdiv(af, bf);
+ return float2half(rf);
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pgather<Eigen::half, Packet8h>(const Eigen::half* from, Index stride)
+{
+ Packet8h result;
+ result.x = _mm_set_epi16(from[7*stride].x, from[6*stride].x, from[5*stride].x, from[4*stride].x, from[3*stride].x, from[2*stride].x, from[1*stride].x, from[0*stride].x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE void pscatter<Eigen::half, Packet8h>(Eigen::half* to, const Packet8h& from, Index stride)
+{
+ EIGEN_ALIGN32 Eigen::half aux[8];
+ pstore(aux, from);
+ to[stride*0].x = aux[0].x;
+ to[stride*1].x = aux[1].x;
+ to[stride*2].x = aux[2].x;
+ to[stride*3].x = aux[3].x;
+ to[stride*4].x = aux[4].x;
+ to[stride*5].x = aux[5].x;
+ to[stride*6].x = aux[6].x;
+ to[stride*7].x = aux[7].x;
+}
+
+template<> EIGEN_STRONG_INLINE Eigen::half predux<Packet8h>(const Packet8h& a) {
+ Packet8f af = half2float(a);
+ float reduced = predux<Packet8f>(af);
+ return Eigen::half(reduced);
+}
+
+template<> EIGEN_STRONG_INLINE Eigen::half predux_max<Packet8h>(const Packet8h& a) {
+ Packet8f af = half2float(a);
+ float reduced = predux_max<Packet8f>(af);
+ return Eigen::half(reduced);
+}
+
+template<> EIGEN_STRONG_INLINE Eigen::half predux_min<Packet8h>(const Packet8h& a) {
+ Packet8f af = half2float(a);
+ float reduced = predux_min<Packet8f>(af);
+ return Eigen::half(reduced);
+}
+
+template<> EIGEN_STRONG_INLINE Eigen::half predux_mul<Packet8h>(const Packet8h& a) {
+ Packet8f af = half2float(a);
+ float reduced = predux_mul<Packet8f>(af);
+ return Eigen::half(reduced);
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h preduxp<Packet8h>(const Packet8h* p) {
+ Packet8f pf[8];
+ pf[0] = half2float(p[0]);
+ pf[1] = half2float(p[1]);
+ pf[2] = half2float(p[2]);
+ pf[3] = half2float(p[3]);
+ pf[4] = half2float(p[4]);
+ pf[5] = half2float(p[5]);
+ pf[6] = half2float(p[6]);
+ pf[7] = half2float(p[7]);
+ Packet8f reduced = preduxp<Packet8f>(pf);
+ return float2half(reduced);
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h preverse(const Packet8h& a)
+{
+ __m128i m = _mm_setr_epi8(14,15,12,13,10,11,8,9,6,7,4,5,2,3,0,1);
+ Packet8h res;
+ res.x = _mm_shuffle_epi8(a.x,m);
+ return res;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pinsertfirst(const Packet8h& a, Eigen::half b)
+{
+ Packet8h res;
+ res.x = _mm_insert_epi16(a.x,int(b.x),0);
+ return res;
+}
+
+template<> EIGEN_STRONG_INLINE Packet8h pinsertlast(const Packet8h& a, Eigen::half b)
+{
+ Packet8h res;
+ res.x = _mm_insert_epi16(a.x,int(b.x),7);
+ return res;
+}
+
+template<int Offset>
+struct palign_impl<Offset,Packet8h>
+{
+ static EIGEN_STRONG_INLINE void run(Packet8h& first, const Packet8h& second)
+ {
+ if (Offset!=0)
+ first.x = _mm_alignr_epi8(second.x,first.x, Offset*2);
+ }
+};
+
+EIGEN_STRONG_INLINE void
+ptranspose(PacketBlock<Packet8h,8>& kernel) {
+ __m128i a = kernel.packet[0].x;
+ __m128i b = kernel.packet[1].x;
+ __m128i c = kernel.packet[2].x;
+ __m128i d = kernel.packet[3].x;
+ __m128i e = kernel.packet[4].x;
+ __m128i f = kernel.packet[5].x;
+ __m128i g = kernel.packet[6].x;
+ __m128i h = kernel.packet[7].x;
+
+ __m128i a03b03 = _mm_unpacklo_epi16(a, b);
+ __m128i c03d03 = _mm_unpacklo_epi16(c, d);
+ __m128i e03f03 = _mm_unpacklo_epi16(e, f);
+ __m128i g03h03 = _mm_unpacklo_epi16(g, h);
+ __m128i a47b47 = _mm_unpackhi_epi16(a, b);
+ __m128i c47d47 = _mm_unpackhi_epi16(c, d);
+ __m128i e47f47 = _mm_unpackhi_epi16(e, f);
+ __m128i g47h47 = _mm_unpackhi_epi16(g, h);
+
+ __m128i a01b01c01d01 = _mm_unpacklo_epi32(a03b03, c03d03);
+ __m128i a23b23c23d23 = _mm_unpackhi_epi32(a03b03, c03d03);
+ __m128i e01f01g01h01 = _mm_unpacklo_epi32(e03f03, g03h03);
+ __m128i e23f23g23h23 = _mm_unpackhi_epi32(e03f03, g03h03);
+ __m128i a45b45c45d45 = _mm_unpacklo_epi32(a47b47, c47d47);
+ __m128i a67b67c67d67 = _mm_unpackhi_epi32(a47b47, c47d47);
+ __m128i e45f45g45h45 = _mm_unpacklo_epi32(e47f47, g47h47);
+ __m128i e67f67g67h67 = _mm_unpackhi_epi32(e47f47, g47h47);
+
+ __m128i a0b0c0d0e0f0g0h0 = _mm_unpacklo_epi64(a01b01c01d01, e01f01g01h01);
+ __m128i a1b1c1d1e1f1g1h1 = _mm_unpackhi_epi64(a01b01c01d01, e01f01g01h01);
+ __m128i a2b2c2d2e2f2g2h2 = _mm_unpacklo_epi64(a23b23c23d23, e23f23g23h23);
+ __m128i a3b3c3d3e3f3g3h3 = _mm_unpackhi_epi64(a23b23c23d23, e23f23g23h23);
+ __m128i a4b4c4d4e4f4g4h4 = _mm_unpacklo_epi64(a45b45c45d45, e45f45g45h45);
+ __m128i a5b5c5d5e5f5g5h5 = _mm_unpackhi_epi64(a45b45c45d45, e45f45g45h45);
+ __m128i a6b6c6d6e6f6g6h6 = _mm_unpacklo_epi64(a67b67c67d67, e67f67g67h67);
+ __m128i a7b7c7d7e7f7g7h7 = _mm_unpackhi_epi64(a67b67c67d67, e67f67g67h67);
+
+ kernel.packet[0].x = a0b0c0d0e0f0g0h0;
+ kernel.packet[1].x = a1b1c1d1e1f1g1h1;
+ kernel.packet[2].x = a2b2c2d2e2f2g2h2;
+ kernel.packet[3].x = a3b3c3d3e3f3g3h3;
+ kernel.packet[4].x = a4b4c4d4e4f4g4h4;
+ kernel.packet[5].x = a5b5c5d5e5f5g5h5;
+ kernel.packet[6].x = a6b6c6d6e6f6g6h6;
+ kernel.packet[7].x = a7b7c7d7e7f7g7h7;
+}
+
+EIGEN_STRONG_INLINE void
+ptranspose(PacketBlock<Packet8h,4>& kernel) {
+ EIGEN_ALIGN32 Eigen::half in[4][8];
+ pstore<Eigen::half>(in[0], kernel.packet[0]);
+ pstore<Eigen::half>(in[1], kernel.packet[1]);
+ pstore<Eigen::half>(in[2], kernel.packet[2]);
+ pstore<Eigen::half>(in[3], kernel.packet[3]);
+
+ EIGEN_ALIGN32 Eigen::half out[4][8];
+
+ for (int i = 0; i < 4; ++i) {
+ for (int j = 0; j < 4; ++j) {
+ out[i][j] = in[j][2*i];
+ }
+ for (int j = 0; j < 4; ++j) {
+ out[i][j+4] = in[j][2*i+1];
+ }
+ }
+
+ kernel.packet[0] = pload<Packet8h>(out[0]);
+ kernel.packet[1] = pload<Packet8h>(out[1]);
+ kernel.packet[2] = pload<Packet8h>(out[2]);
+ kernel.packet[3] = pload<Packet8h>(out[3]);
+}
+
} // end namespace internal
} // end namespace Eigen
diff --git a/Eigen/src/Core/arch/AVX/TypeCasting.h b/Eigen/src/Core/arch/AVX/TypeCasting.h
index 7d2e1e67f..181043588 100644
--- a/Eigen/src/Core/arch/AVX/TypeCasting.h
+++ b/Eigen/src/Core/arch/AVX/TypeCasting.h
@@ -52,6 +52,36 @@ template<> EIGEN_STRONG_INLINE Packet8f preinterpret<Packet8f,Packet8i>(const Pa
return _mm256_castsi256_ps(a);
}
+#ifndef EIGEN_VECTORIZE_AVX512
+
+template <>
+struct type_casting_traits<Eigen::half, float> {
+ enum {
+ VectorizedCast = 1,
+ SrcCoeffRatio = 1,
+ TgtCoeffRatio = 1
+ };
+};
+
+template<> EIGEN_STRONG_INLINE Packet8f pcast<Packet8h, Packet8f>(const Packet8h& a) {
+ return half2float(a);
+}
+
+template <>
+struct type_casting_traits<float, Eigen::half> {
+ enum {
+ VectorizedCast = 1,
+ SrcCoeffRatio = 1,
+ TgtCoeffRatio = 1
+ };
+};
+
+#endif // EIGEN_VECTORIZE_AVX512
+
+template<> EIGEN_STRONG_INLINE Packet8h pcast<Packet8f, Packet8h>(const Packet8f& a) {
+ return float2half(a);
+}
+
} // end namespace internal
} // end namespace Eigen
diff --git a/Eigen/src/Core/arch/AVX512/MathFunctions.h b/Eigen/src/Core/arch/AVX512/MathFunctions.h
index c2158c538..9e37a720b 100644
--- a/Eigen/src/Core/arch/AVX512/MathFunctions.h
+++ b/Eigen/src/Core/arch/AVX512/MathFunctions.h
@@ -393,6 +393,18 @@ pcos<Packet16f>(const Packet16f& _x) {
return pcos_float(_x);
}
+#if defined(EIGEN_VECTORIZE_AVX512DQ)
+template<> EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED
+Packet16f plog1p<Packet16f>(const Packet16f& _x) {
+ return generic_plog1p(_x);
+}
+
+template<> EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED
+Packet16f pexpm1<Packet16f>(const Packet16f& _x) {
+ return generic_expm1(_x);
+}
+#endif
+
} // end namespace internal
} // end namespace Eigen
diff --git a/Eigen/src/Core/arch/AVX512/PacketMath.h b/Eigen/src/Core/arch/AVX512/PacketMath.h
index 1312829d8..744d7c4e4 100644
--- a/Eigen/src/Core/arch/AVX512/PacketMath.h
+++ b/Eigen/src/Core/arch/AVX512/PacketMath.h
@@ -44,6 +44,41 @@ template <>
struct is_arithmetic<__m512d> {
enum { value = true };
};
+typedef struct {
+ __m256i x;
+} Packet16h;
+
+
+template<> struct is_arithmetic<Packet16h> { enum { value = true }; };
+
+template <>
+struct packet_traits<half> : default_packet_traits {
+ typedef Packet16h type;
+ // There is no half-size packet for Packet16h.
+ typedef Packet16h half;
+ enum {
+ Vectorizable = 1,
+ AlignedOnScalar = 1,
+ size = 16,
+ HasHalfPacket = 0,
+ HasAdd = 1,
+ HasSub = 1,
+ HasMul = 1,
+ HasDiv = 1,
+ HasNegate = 1,
+ HasAbs = 0,
+ HasAbs2 = 0,
+ HasMin = 0,
+ HasMax = 0,
+ HasConj = 0,
+ HasSetLinear = 0,
+ HasSqrt = 0,
+ HasRsqrt = 0,
+ HasExp = 0,
+ HasLog = 0,
+ HasBlend = 0
+ };
+};
template<> struct packet_traits<float> : default_packet_traits
{
@@ -121,6 +156,13 @@ struct unpacket_traits<Packet16i> {
enum { size = 16, alignment=Aligned64, vectorizable=false, masked_load_available=false, masked_store_available=false };
};
+template<>
+struct unpacket_traits<Packet16h> {
+ typedef Eigen::half type;
+ typedef Packet16h half;
+ enum {size=16, alignment=Aligned32, vectorizable=true, masked_load_available=false, masked_store_available=false};
+};
+
template <>
EIGEN_STRONG_INLINE Packet16f pset1<Packet16f>(const float& from) {
return _mm512_set1_ps(from);
@@ -1398,6 +1440,467 @@ template<> EIGEN_STRONG_INLINE Packet16f preinterpret<Packet16f,Packet16i>(const
return _mm512_castsi512_ps(a);
}
+
+// Packet math for Eigen::half
+template<> EIGEN_STRONG_INLINE Packet16h pset1<Packet16h>(const Eigen::half& from) {
+ Packet16h result;
+ result.x = _mm256_set1_epi16(from.x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Eigen::half pfirst<Packet16h>(const Packet16h& from) {
+ return half_impl::raw_uint16_to_half(static_cast<unsigned short>(_mm256_extract_epi16(from.x, 0)));
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pload<Packet16h>(const Eigen::half* from) {
+ Packet16h result;
+ result.x = _mm256_load_si256(reinterpret_cast<const __m256i*>(from));
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h ploadu<Packet16h>(const Eigen::half* from) {
+ Packet16h result;
+ result.x = _mm256_loadu_si256(reinterpret_cast<const __m256i*>(from));
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE void pstore<half>(Eigen::half* to, const Packet16h& from) {
+ // (void*) -> workaround clang warning:
+ // cast from 'Eigen::half *' to '__m256i *' increases required alignment from 2 to 32
+ _mm256_store_si256((__m256i*)(void*)to, from.x);
+}
+
+template<> EIGEN_STRONG_INLINE void pstoreu<half>(Eigen::half* to, const Packet16h& from) {
+ // (void*) -> workaround clang warning:
+ // cast from 'Eigen::half *' to '__m256i *' increases required alignment from 2 to 32
+ _mm256_storeu_si256((__m256i*)(void*)to, from.x);
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h
+ploaddup<Packet16h>(const Eigen::half* from) {
+ Packet16h result;
+ unsigned short a = from[0].x;
+ unsigned short b = from[1].x;
+ unsigned short c = from[2].x;
+ unsigned short d = from[3].x;
+ unsigned short e = from[4].x;
+ unsigned short f = from[5].x;
+ unsigned short g = from[6].x;
+ unsigned short h = from[7].x;
+ result.x = _mm256_set_epi16(h, h, g, g, f, f, e, e, d, d, c, c, b, b, a, a);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h
+ploadquad(const Eigen::half* from) {
+ Packet16h result;
+ unsigned short a = from[0].x;
+ unsigned short b = from[1].x;
+ unsigned short c = from[2].x;
+ unsigned short d = from[3].x;
+ result.x = _mm256_set_epi16(d, d, d, d, c, c, c, c, b, b, b, b, a, a, a, a);
+ return result;
+}
+
+EIGEN_STRONG_INLINE Packet16f half2float(const Packet16h& a) {
+#ifdef EIGEN_HAS_FP16_C
+ return _mm512_cvtph_ps(a.x);
+#else
+ EIGEN_ALIGN64 half aux[16];
+ pstore(aux, a);
+ float f0(aux[0]);
+ float f1(aux[1]);
+ float f2(aux[2]);
+ float f3(aux[3]);
+ float f4(aux[4]);
+ float f5(aux[5]);
+ float f6(aux[6]);
+ float f7(aux[7]);
+ float f8(aux[8]);
+ float f9(aux[9]);
+ float fa(aux[10]);
+ float fb(aux[11]);
+ float fc(aux[12]);
+ float fd(aux[13]);
+ float fe(aux[14]);
+ float ff(aux[15]);
+
+ return _mm512_set_ps(
+ ff, fe, fd, fc, fb, fa, f9, f8, f7, f6, f5, f4, f3, f2, f1, f0);
+#endif
+}
+
+EIGEN_STRONG_INLINE Packet16h float2half(const Packet16f& a) {
+#ifdef EIGEN_HAS_FP16_C
+ Packet16h result;
+ result.x = _mm512_cvtps_ph(a, _MM_FROUND_TO_NEAREST_INT|_MM_FROUND_NO_EXC);
+ return result;
+#else
+ EIGEN_ALIGN64 float aux[16];
+ pstore(aux, a);
+ half h0(aux[0]);
+ half h1(aux[1]);
+ half h2(aux[2]);
+ half h3(aux[3]);
+ half h4(aux[4]);
+ half h5(aux[5]);
+ half h6(aux[6]);
+ half h7(aux[7]);
+ half h8(aux[8]);
+ half h9(aux[9]);
+ half ha(aux[10]);
+ half hb(aux[11]);
+ half hc(aux[12]);
+ half hd(aux[13]);
+ half he(aux[14]);
+ half hf(aux[15]);
+
+ Packet16h result;
+ result.x = _mm256_set_epi16(
+ hf.x, he.x, hd.x, hc.x, hb.x, ha.x, h9.x, h8.x,
+ h7.x, h6.x, h5.x, h4.x, h3.x, h2.x, h1.x, h0.x);
+ return result;
+#endif
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pnot(const Packet16h& a) {
+ Packet16h r; r.x = _mm256_xor_si256(a.x, pcmp_eq(a.x, a.x)); return r;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h ptrue(const Packet16h& a) {
+ Packet16h r; r.x = Packet8i(ptrue(a.x)); return r;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h por(const Packet16h& a,const Packet16h& b) {
+ // in some cases Packet8i is a wrapper around __m256i, so we need to
+ // cast to Packet8i to call the correct overload.
+ Packet16h r; r.x = por(Packet8i(a.x),Packet8i(b.x)); return r;
+}
+template<> EIGEN_STRONG_INLINE Packet16h pxor(const Packet16h& a,const Packet16h& b) {
+ Packet16h r; r.x = pxor(Packet8i(a.x),Packet8i(b.x)); return r;
+}
+template<> EIGEN_STRONG_INLINE Packet16h pand(const Packet16h& a,const Packet16h& b) {
+ Packet16h r; r.x = pand(Packet8i(a.x),Packet8i(b.x)); return r;
+}
+template<> EIGEN_STRONG_INLINE Packet16h pandnot(const Packet16h& a,const Packet16h& b) {
+ Packet16h r; r.x = pandnot(Packet8i(a.x),Packet8i(b.x)); return r;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pselect(const Packet16h& mask, const Packet16h& a, const Packet16h& b) {
+ Packet16h r; r.x = _mm256_blendv_epi8(b.x, a.x, mask.x); return r;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pcmp_eq(const Packet16h& a,const Packet16h& b) {
+ Packet16f af = half2float(a);
+ Packet16f bf = half2float(b);
+ Packet16f rf = pcmp_eq(af, bf);
+ // Pack the 32-bit flags into 16-bits flags.
+ __m256i lo = _mm256_castps_si256(extract256<0>(rf));
+ __m256i hi = _mm256_castps_si256(extract256<1>(rf));
+ __m128i result_lo = _mm_packs_epi32(_mm256_extractf128_si256(lo, 0),
+ _mm256_extractf128_si256(lo, 1));
+ __m128i result_hi = _mm_packs_epi32(_mm256_extractf128_si256(hi, 0),
+ _mm256_extractf128_si256(hi, 1));
+ Packet16h result; result.x = _mm256_insertf128_si256(_mm256_castsi128_si256(result_lo), result_hi, 1);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pnegate(const Packet16h& a) {
+ Packet16h sign_mask; sign_mask.x = _mm256_set1_epi16(static_cast<unsigned short>(0x8000));
+ Packet16h result; result.x = _mm256_xor_si256(a.x, sign_mask.x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h padd<Packet16h>(const Packet16h& a, const Packet16h& b) {
+ Packet16f af = half2float(a);
+ Packet16f bf = half2float(b);
+ Packet16f rf = padd(af, bf);
+ return float2half(rf);
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h psub<Packet16h>(const Packet16h& a, const Packet16h& b) {
+ Packet16f af = half2float(a);
+ Packet16f bf = half2float(b);
+ Packet16f rf = psub(af, bf);
+ return float2half(rf);
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pmul<Packet16h>(const Packet16h& a, const Packet16h& b) {
+ Packet16f af = half2float(a);
+ Packet16f bf = half2float(b);
+ Packet16f rf = pmul(af, bf);
+ return float2half(rf);
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pdiv<Packet16h>(const Packet16h& a, const Packet16h& b) {
+ Packet16f af = half2float(a);
+ Packet16f bf = half2float(b);
+ Packet16f rf = pdiv(af, bf);
+ return float2half(rf);
+}
+
+template<> EIGEN_STRONG_INLINE half predux<Packet16h>(const Packet16h& from) {
+ Packet16f from_float = half2float(from);
+ return half(predux(from_float));
+}
+
+template<> EIGEN_STRONG_INLINE half predux_mul<Packet16h>(const Packet16h& from) {
+ Packet16f from_float = half2float(from);
+ return half(predux_mul(from_float));
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h preduxp<Packet16h>(const Packet16h* p) {
+ Packet16f pf[16];
+ pf[0] = half2float(p[0]);
+ pf[1] = half2float(p[1]);
+ pf[2] = half2float(p[2]);
+ pf[3] = half2float(p[3]);
+ pf[4] = half2float(p[4]);
+ pf[5] = half2float(p[5]);
+ pf[6] = half2float(p[6]);
+ pf[7] = half2float(p[7]);
+ pf[8] = half2float(p[8]);
+ pf[9] = half2float(p[9]);
+ pf[10] = half2float(p[10]);
+ pf[11] = half2float(p[11]);
+ pf[12] = half2float(p[12]);
+ pf[13] = half2float(p[13]);
+ pf[14] = half2float(p[14]);
+ pf[15] = half2float(p[15]);
+ Packet16f reduced = preduxp<Packet16f>(pf);
+ return float2half(reduced);
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h preverse(const Packet16h& a)
+{
+ __m128i m = _mm_setr_epi8(14,15,12,13,10,11,8,9,6,7,4,5,2,3,0,1);
+ Packet16h res;
+ res.x = _mm256_insertf128_si256(
+ _mm256_castsi128_si256(_mm_shuffle_epi8(_mm256_extractf128_si256(a.x,1),m)),
+ _mm_shuffle_epi8(_mm256_extractf128_si256(a.x,0),m), 1);
+ return res;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pinsertfirst(const Packet16h& a, Eigen::half b)
+{
+ Packet16h res;
+ res.x = _mm256_insert_epi16(a.x,b.x,0);
+ return res;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pinsertlast(const Packet16h& a, Eigen::half b)
+{
+ Packet16h res;
+ res.x = _mm256_insert_epi16(a.x,b.x,15);
+ return res;
+}
+
+template<> EIGEN_STRONG_INLINE Packet16h pgather<Eigen::half, Packet16h>(const Eigen::half* from, Index stride)
+{
+ Packet16h result;
+ result.x = _mm256_set_epi16(
+ from[15*stride].x, from[14*stride].x, from[13*stride].x, from[12*stride].x,
+ from[11*stride].x, from[10*stride].x, from[9*stride].x, from[8*stride].x,
+ from[7*stride].x, from[6*stride].x, from[5*stride].x, from[4*stride].x,
+ from[3*stride].x, from[2*stride].x, from[1*stride].x, from[0*stride].x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE void pscatter<half, Packet16h>(half* to, const Packet16h& from, Index stride)
+{
+ EIGEN_ALIGN64 half aux[16];
+ pstore(aux, from);
+ to[stride*0].x = aux[0].x;
+ to[stride*1].x = aux[1].x;
+ to[stride*2].x = aux[2].x;
+ to[stride*3].x = aux[3].x;
+ to[stride*4].x = aux[4].x;
+ to[stride*5].x = aux[5].x;
+ to[stride*6].x = aux[6].x;
+ to[stride*7].x = aux[7].x;
+ to[stride*8].x = aux[8].x;
+ to[stride*9].x = aux[9].x;
+ to[stride*10].x = aux[10].x;
+ to[stride*11].x = aux[11].x;
+ to[stride*12].x = aux[12].x;
+ to[stride*13].x = aux[13].x;
+ to[stride*14].x = aux[14].x;
+ to[stride*15].x = aux[15].x;
+}
+
+EIGEN_STRONG_INLINE void
+ptranspose(PacketBlock<Packet16h,16>& kernel) {
+ __m256i a = kernel.packet[0].x;
+ __m256i b = kernel.packet[1].x;
+ __m256i c = kernel.packet[2].x;
+ __m256i d = kernel.packet[3].x;
+ __m256i e = kernel.packet[4].x;
+ __m256i f = kernel.packet[5].x;
+ __m256i g = kernel.packet[6].x;
+ __m256i h = kernel.packet[7].x;
+ __m256i i = kernel.packet[8].x;
+ __m256i j = kernel.packet[9].x;
+ __m256i k = kernel.packet[10].x;
+ __m256i l = kernel.packet[11].x;
+ __m256i m = kernel.packet[12].x;
+ __m256i n = kernel.packet[13].x;
+ __m256i o = kernel.packet[14].x;
+ __m256i p = kernel.packet[15].x;
+
+ __m256i ab_07 = _mm256_unpacklo_epi16(a, b);
+ __m256i cd_07 = _mm256_unpacklo_epi16(c, d);
+ __m256i ef_07 = _mm256_unpacklo_epi16(e, f);
+ __m256i gh_07 = _mm256_unpacklo_epi16(g, h);
+ __m256i ij_07 = _mm256_unpacklo_epi16(i, j);
+ __m256i kl_07 = _mm256_unpacklo_epi16(k, l);
+ __m256i mn_07 = _mm256_unpacklo_epi16(m, n);
+ __m256i op_07 = _mm256_unpacklo_epi16(o, p);
+
+ __m256i ab_8f = _mm256_unpackhi_epi16(a, b);
+ __m256i cd_8f = _mm256_unpackhi_epi16(c, d);
+ __m256i ef_8f = _mm256_unpackhi_epi16(e, f);
+ __m256i gh_8f = _mm256_unpackhi_epi16(g, h);
+ __m256i ij_8f = _mm256_unpackhi_epi16(i, j);
+ __m256i kl_8f = _mm256_unpackhi_epi16(k, l);
+ __m256i mn_8f = _mm256_unpackhi_epi16(m, n);
+ __m256i op_8f = _mm256_unpackhi_epi16(o, p);
+
+ __m256i abcd_03 = _mm256_unpacklo_epi32(ab_07, cd_07);
+ __m256i abcd_47 = _mm256_unpackhi_epi32(ab_07, cd_07);
+ __m256i efgh_03 = _mm256_unpacklo_epi32(ef_07, gh_07);
+ __m256i efgh_47 = _mm256_unpackhi_epi32(ef_07, gh_07);
+ __m256i ijkl_03 = _mm256_unpacklo_epi32(ij_07, kl_07);
+ __m256i ijkl_47 = _mm256_unpackhi_epi32(ij_07, kl_07);
+ __m256i mnop_03 = _mm256_unpacklo_epi32(mn_07, op_07);
+ __m256i mnop_47 = _mm256_unpackhi_epi32(mn_07, op_07);
+
+ __m256i abcd_8b = _mm256_unpacklo_epi32(ab_8f, cd_8f);
+ __m256i abcd_cf = _mm256_unpackhi_epi32(ab_8f, cd_8f);
+ __m256i efgh_8b = _mm256_unpacklo_epi32(ef_8f, gh_8f);
+ __m256i efgh_cf = _mm256_unpackhi_epi32(ef_8f, gh_8f);
+ __m256i ijkl_8b = _mm256_unpacklo_epi32(ij_8f, kl_8f);
+ __m256i ijkl_cf = _mm256_unpackhi_epi32(ij_8f, kl_8f);
+ __m256i mnop_8b = _mm256_unpacklo_epi32(mn_8f, op_8f);
+ __m256i mnop_cf = _mm256_unpackhi_epi32(mn_8f, op_8f);
+
+ __m256i abcdefgh_01 = _mm256_unpacklo_epi64(abcd_03, efgh_03);
+ __m256i abcdefgh_23 = _mm256_unpackhi_epi64(abcd_03, efgh_03);
+ __m256i ijklmnop_01 = _mm256_unpacklo_epi64(ijkl_03, mnop_03);
+ __m256i ijklmnop_23 = _mm256_unpackhi_epi64(ijkl_03, mnop_03);
+ __m256i abcdefgh_45 = _mm256_unpacklo_epi64(abcd_47, efgh_47);
+ __m256i abcdefgh_67 = _mm256_unpackhi_epi64(abcd_47, efgh_47);
+ __m256i ijklmnop_45 = _mm256_unpacklo_epi64(ijkl_47, mnop_47);
+ __m256i ijklmnop_67 = _mm256_unpackhi_epi64(ijkl_47, mnop_47);
+ __m256i abcdefgh_89 = _mm256_unpacklo_epi64(abcd_8b, efgh_8b);
+ __m256i abcdefgh_ab = _mm256_unpackhi_epi64(abcd_8b, efgh_8b);
+ __m256i ijklmnop_89 = _mm256_unpacklo_epi64(ijkl_8b, mnop_8b);
+ __m256i ijklmnop_ab = _mm256_unpackhi_epi64(ijkl_8b, mnop_8b);
+ __m256i abcdefgh_cd = _mm256_unpacklo_epi64(abcd_cf, efgh_cf);
+ __m256i abcdefgh_ef = _mm256_unpackhi_epi64(abcd_cf, efgh_cf);
+ __m256i ijklmnop_cd = _mm256_unpacklo_epi64(ijkl_cf, mnop_cf);
+ __m256i ijklmnop_ef = _mm256_unpackhi_epi64(ijkl_cf, mnop_cf);
+
+ // NOTE: no unpacklo/hi instr in this case, so using permute instr.
+ __m256i a_p_0 = _mm256_permute2x128_si256(abcdefgh_01, ijklmnop_01, 0x20);
+ __m256i a_p_1 = _mm256_permute2x128_si256(abcdefgh_23, ijklmnop_23, 0x20);
+ __m256i a_p_2 = _mm256_permute2x128_si256(abcdefgh_45, ijklmnop_45, 0x20);
+ __m256i a_p_3 = _mm256_permute2x128_si256(abcdefgh_67, ijklmnop_67, 0x20);
+ __m256i a_p_4 = _mm256_permute2x128_si256(abcdefgh_89, ijklmnop_89, 0x20);
+ __m256i a_p_5 = _mm256_permute2x128_si256(abcdefgh_ab, ijklmnop_ab, 0x20);
+ __m256i a_p_6 = _mm256_permute2x128_si256(abcdefgh_cd, ijklmnop_cd, 0x20);
+ __m256i a_p_7 = _mm256_permute2x128_si256(abcdefgh_ef, ijklmnop_ef, 0x20);
+ __m256i a_p_8 = _mm256_permute2x128_si256(abcdefgh_01, ijklmnop_01, 0x31);
+ __m256i a_p_9 = _mm256_permute2x128_si256(abcdefgh_23, ijklmnop_23, 0x31);
+ __m256i a_p_a = _mm256_permute2x128_si256(abcdefgh_45, ijklmnop_45, 0x31);
+ __m256i a_p_b = _mm256_permute2x128_si256(abcdefgh_67, ijklmnop_67, 0x31);
+ __m256i a_p_c = _mm256_permute2x128_si256(abcdefgh_89, ijklmnop_89, 0x31);
+ __m256i a_p_d = _mm256_permute2x128_si256(abcdefgh_ab, ijklmnop_ab, 0x31);
+ __m256i a_p_e = _mm256_permute2x128_si256(abcdefgh_cd, ijklmnop_cd, 0x31);
+ __m256i a_p_f = _mm256_permute2x128_si256(abcdefgh_ef, ijklmnop_ef, 0x31);
+
+ kernel.packet[0].x = a_p_0;
+ kernel.packet[1].x = a_p_1;
+ kernel.packet[2].x = a_p_2;
+ kernel.packet[3].x = a_p_3;
+ kernel.packet[4].x = a_p_4;
+ kernel.packet[5].x = a_p_5;
+ kernel.packet[6].x = a_p_6;
+ kernel.packet[7].x = a_p_7;
+ kernel.packet[8].x = a_p_8;
+ kernel.packet[9].x = a_p_9;
+ kernel.packet[10].x = a_p_a;
+ kernel.packet[11].x = a_p_b;
+ kernel.packet[12].x = a_p_c;
+ kernel.packet[13].x = a_p_d;
+ kernel.packet[14].x = a_p_e;
+ kernel.packet[15].x = a_p_f;
+}
+
+EIGEN_STRONG_INLINE void
+ptranspose(PacketBlock<Packet16h,8>& kernel) {
+ EIGEN_ALIGN64 half in[8][16];
+ pstore<half>(in[0], kernel.packet[0]);
+ pstore<half>(in[1], kernel.packet[1]);
+ pstore<half>(in[2], kernel.packet[2]);
+ pstore<half>(in[3], kernel.packet[3]);
+ pstore<half>(in[4], kernel.packet[4]);
+ pstore<half>(in[5], kernel.packet[5]);
+ pstore<half>(in[6], kernel.packet[6]);
+ pstore<half>(in[7], kernel.packet[7]);
+
+ EIGEN_ALIGN64 half out[8][16];
+
+ for (int i = 0; i < 8; ++i) {
+ for (int j = 0; j < 8; ++j) {
+ out[i][j] = in[j][2*i];
+ }
+ for (int j = 0; j < 8; ++j) {
+ out[i][j+8] = in[j][2*i+1];
+ }
+ }
+
+ kernel.packet[0] = pload<Packet16h>(out[0]);
+ kernel.packet[1] = pload<Packet16h>(out[1]);
+ kernel.packet[2] = pload<Packet16h>(out[2]);
+ kernel.packet[3] = pload<Packet16h>(out[3]);
+ kernel.packet[4] = pload<Packet16h>(out[4]);
+ kernel.packet[5] = pload<Packet16h>(out[5]);
+ kernel.packet[6] = pload<Packet16h>(out[6]);
+ kernel.packet[7] = pload<Packet16h>(out[7]);
+}
+
+EIGEN_STRONG_INLINE void
+ptranspose(PacketBlock<Packet16h,4>& kernel) {
+ EIGEN_ALIGN64 half in[4][16];
+ pstore<half>(in[0], kernel.packet[0]);
+ pstore<half>(in[1], kernel.packet[1]);
+ pstore<half>(in[2], kernel.packet[2]);
+ pstore<half>(in[3], kernel.packet[3]);
+
+ EIGEN_ALIGN64 half out[4][16];
+
+ for (int i = 0; i < 4; ++i) {
+ for (int j = 0; j < 4; ++j) {
+ out[i][j] = in[j][4*i];
+ }
+ for (int j = 0; j < 4; ++j) {
+ out[i][j+4] = in[j][4*i+1];
+ }
+ for (int j = 0; j < 4; ++j) {
+ out[i][j+8] = in[j][4*i+2];
+ }
+ for (int j = 0; j < 4; ++j) {
+ out[i][j+12] = in[j][4*i+3];
+ }
+ }
+
+ kernel.packet[0] = pload<Packet16h>(out[0]);
+ kernel.packet[1] = pload<Packet16h>(out[1]);
+ kernel.packet[2] = pload<Packet16h>(out[2]);
+ kernel.packet[3] = pload<Packet16h>(out[3]);
+}
+
+
} // end namespace internal
} // end namespace Eigen
diff --git a/Eigen/src/Core/arch/AVX512/TypeCasting.h b/Eigen/src/Core/arch/AVX512/TypeCasting.h
new file mode 100644
index 000000000..a82176941
--- /dev/null
+++ b/Eigen/src/Core/arch/AVX512/TypeCasting.h
@@ -0,0 +1,47 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra.
+//
+// Copyright (C) 2019 Rasmus Munk Larsen <rmlarsen@google.com>
+//
+// This Source Code Form is subject to the terms of the Mozilla
+// Public License v. 2.0. If a copy of the MPL was not distributed
+// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
+
+#ifndef EIGEN_TYPE_CASTING_AVX512_H
+#define EIGEN_TYPE_CASTING_AVX512_H
+
+namespace Eigen {
+
+namespace internal {
+
+template <>
+struct type_casting_traits<half, float> {
+ enum {
+ VectorizedCast = 1,
+ SrcCoeffRatio = 1,
+ TgtCoeffRatio = 1
+ };
+};
+
+template<> EIGEN_STRONG_INLINE Packet16f pcast<Packet16h, Packet16f>(const Packet16h& a) {
+ return half2float(a);
+}
+
+template <>
+struct type_casting_traits<float, half> {
+ enum {
+ VectorizedCast = 1,
+ SrcCoeffRatio = 1,
+ TgtCoeffRatio = 1
+ };
+};
+
+template<> EIGEN_STRONG_INLINE Packet16h pcast<Packet16f, Packet16h>(const Packet16f& a) {
+ return float2half(a);
+}
+
+} // end namespace internal
+
+} // end namespace Eigen
+
+#endif // EIGEN_TYPE_CASTING_AVX512_H
diff --git a/Eigen/src/Core/arch/AltiVec/PacketMath.h b/Eigen/src/Core/arch/AltiVec/PacketMath.h
index 4b770d036..7ee290a29 100755
--- a/Eigen/src/Core/arch/AltiVec/PacketMath.h
+++ b/Eigen/src/Core/arch/AltiVec/PacketMath.h
@@ -83,15 +83,6 @@ static Packet4i p4i_COUNTDOWN = { 0, 1, 2, 3 };
static Packet16uc p16uc_REVERSE32 = { 12,13,14,15, 8,9,10,11, 4,5,6,7, 0,1,2,3 };
static Packet16uc p16uc_DUPLICATE32_HI = { 0,1,2,3, 0,1,2,3, 4,5,6,7, 4,5,6,7 };
-// Mask alignment
-#ifdef __PPC64__
-#define _EIGEN_MASK_ALIGNMENT 0xfffffffffffffff0
-#else
-#define _EIGEN_MASK_ALIGNMENT 0xfffffff0
-#endif
-
-#define _EIGEN_ALIGNED_PTR(x) ((std::ptrdiff_t)(x) & _EIGEN_MASK_ALIGNMENT)
-
// Handle endianness properly while loading constants
// Define global static constants:
#ifdef _BIG_ENDIAN
@@ -249,42 +240,38 @@ inline std::ostream & operator <<(std::ostream & s, const Packet4ui & v)
// Need to define them first or we get specialization after instantiation errors
template<> EIGEN_STRONG_INLINE Packet4f pload<Packet4f>(const float* from)
{
+ // some versions of GCC throw "unused-but-set-parameter".
+ // ignoring these warnings for now.
+ EIGEN_UNUSED_VARIABLE(from);
EIGEN_DEBUG_ALIGNED_LOAD
-#ifdef __VSX__
- return vec_vsx_ld(0, from);
-#else
return vec_ld(0, from);
-#endif
}
template<> EIGEN_STRONG_INLINE Packet4i pload<Packet4i>(const int* from)
{
+ // some versions of GCC throw "unused-but-set-parameter".
+ // ignoring these warnings for now.
+ EIGEN_UNUSED_VARIABLE(from);
EIGEN_DEBUG_ALIGNED_LOAD
-#ifdef __VSX__
- return vec_vsx_ld(0, from);
-#else
return vec_ld(0, from);
-#endif
}
template<> EIGEN_STRONG_INLINE void pstore<float>(float* to, const Packet4f& from)
{
+ // some versions of GCC throw "unused-but-set-parameter" (float *to).
+ // ignoring these warnings for now.
+ EIGEN_UNUSED_VARIABLE(to);
EIGEN_DEBUG_ALIGNED_STORE
-#ifdef __VSX__
- vec_vsx_st(from, 0, to);
-#else
vec_st(from, 0, to);
-#endif
}
template<> EIGEN_STRONG_INLINE void pstore<int>(int* to, const Packet4i& from)
{
+ // some versions of GCC throw "unused-but-set-parameter" (float *to).
+ // ignoring these warnings for now.
+ EIGEN_UNUSED_VARIABLE(to);
EIGEN_DEBUG_ALIGNED_STORE
-#ifdef __VSX__
- vec_vsx_st(from, 0, to);
-#else
vec_st(from, 0, to);
-#endif
}
template<> EIGEN_STRONG_INLINE Packet4f pset1<Packet4f>(const float& from) {
@@ -452,7 +439,7 @@ template<> EIGEN_STRONG_INLINE Packet4f pandnot<Packet4f>(const Packet4f& a, con
template<> EIGEN_STRONG_INLINE Packet4i pandnot<Packet4i>(const Packet4i& a, const Packet4i& b) { return vec_and(a, vec_nor(b, b)); }
template<> EIGEN_STRONG_INLINE Packet4f pselect(const Packet4f& mask, const Packet4f& a, const Packet4f& b) {
- return vec_sel(b, a, mask);
+ return vec_sel(b, a, reinterpret_cast<Packet4ui>(mask));
}
template<> EIGEN_STRONG_INLINE Packet4f pround<Packet4f>(const Packet4f& a) { return vec_round(a); }
@@ -487,12 +474,12 @@ template<> EIGEN_STRONG_INLINE Packet4i ploadu<Packet4i>(const int* from)
template<> EIGEN_STRONG_INLINE Packet4i ploadu<Packet4i>(const int* from)
{
EIGEN_DEBUG_UNALIGNED_LOAD
- return (Packet4i) vec_vsx_ld((long)from & 15, (const int*) _EIGEN_ALIGNED_PTR(from));
+ return vec_vsx_ld(0, from);
}
template<> EIGEN_STRONG_INLINE Packet4f ploadu<Packet4f>(const float* from)
{
EIGEN_DEBUG_UNALIGNED_LOAD
- return (Packet4f) vec_vsx_ld((long)from & 15, (const float*) _EIGEN_ALIGNED_PTR(from));
+ return vec_vsx_ld(0, from);
}
#endif
@@ -552,13 +539,13 @@ template<> EIGEN_STRONG_INLINE void pstoreu<int>(int* to, const Packet4i& f
// We also need to redefine little endian loading of Packet4i/Packet4f using VSX
template<> EIGEN_STRONG_INLINE void pstoreu<int>(int* to, const Packet4i& from)
{
- EIGEN_DEBUG_ALIGNED_STORE
- vec_vsx_st(from, (long)to & 15, (int*) _EIGEN_ALIGNED_PTR(to));
+ EIGEN_DEBUG_UNALIGNED_STORE
+ vec_vsx_st(from, 0, to);
}
template<> EIGEN_STRONG_INLINE void pstoreu<float>(float* to, const Packet4f& from)
{
- EIGEN_DEBUG_ALIGNED_STORE
- vec_vsx_st(from, (long)to & 15, (float*) _EIGEN_ALIGNED_PTR(to));
+ EIGEN_DEBUG_UNALIGNED_STORE
+ vec_vsx_st(from, 0, to);
}
#endif
@@ -949,21 +936,13 @@ inline std::ostream & operator <<(std::ostream & s, const Packet2d & v)
template<> EIGEN_STRONG_INLINE Packet2d pload<Packet2d>(const double* from)
{
EIGEN_DEBUG_ALIGNED_LOAD
-#ifdef __VSX__
- return vec_vsx_ld(0, from);
-#else
- return vec_ld(0, from);
-#endif
+ return vec_xl(0, const_cast<double *>(from)); // cast needed by Clang
}
template<> EIGEN_STRONG_INLINE void pstore<double>(double* to, const Packet2d& from)
{
EIGEN_DEBUG_ALIGNED_STORE
-#ifdef __VSX__
- vec_vsx_st(from, 0, to);
-#else
- vec_st(from, 0, to);
-#endif
+ vec_xst(from, 0, to);
}
template<> EIGEN_STRONG_INLINE Packet2d pset1<Packet2d>(const double& from) {
@@ -1030,6 +1009,14 @@ template<> EIGEN_STRONG_INLINE Packet2d pmax<Packet2d>(const Packet2d& a, const
return ret;
}
+template<> EIGEN_STRONG_INLINE Packet2d pcmp_le(const Packet2d& a, const Packet2d& b) { return reinterpret_cast<Packet2d>(vec_cmple(a,b)); }
+template<> EIGEN_STRONG_INLINE Packet2d pcmp_lt(const Packet2d& a, const Packet2d& b) { return reinterpret_cast<Packet2d>(vec_cmplt(a,b)); }
+template<> EIGEN_STRONG_INLINE Packet2d pcmp_eq(const Packet2d& a, const Packet2d& b) { return reinterpret_cast<Packet2d>(vec_cmpeq(a,b)); }
+template<> EIGEN_STRONG_INLINE Packet2d pcmp_lt_or_nan(const Packet2d& a, const Packet2d& b) {
+ Packet2d c = reinterpret_cast<Packet2d>(vec_cmpge(a,b));
+ return vec_nor(c,c);
+}
+
template<> EIGEN_STRONG_INLINE Packet2d pand<Packet2d>(const Packet2d& a, const Packet2d& b) { return vec_and(a, b); }
template<> EIGEN_STRONG_INLINE Packet2d por<Packet2d>(const Packet2d& a, const Packet2d& b) { return vec_or(a, b); }
@@ -1044,8 +1031,8 @@ template<> EIGEN_STRONG_INLINE Packet2d pfloor<Packet2d>(const Packet2d& a) { re
template<> EIGEN_STRONG_INLINE Packet2d ploadu<Packet2d>(const double* from)
{
- EIGEN_DEBUG_ALIGNED_LOAD
- return (Packet2d) vec_vsx_ld((long)from & 15, (const double*) _EIGEN_ALIGNED_PTR(from));
+ EIGEN_DEBUG_UNALIGNED_LOAD
+ return vec_vsx_ld(0, from);
}
template<> EIGEN_STRONG_INLINE Packet2d ploaddup<Packet2d>(const double* from)
@@ -1058,8 +1045,8 @@ template<> EIGEN_STRONG_INLINE Packet2d ploaddup<Packet2d>(const double* from)
template<> EIGEN_STRONG_INLINE void pstoreu<double>(double* to, const Packet2d& from)
{
- EIGEN_DEBUG_ALIGNED_STORE
- vec_vsx_st((Packet4f)from, (long)to & 15, (float*) _EIGEN_ALIGNED_PTR(to));
+ EIGEN_DEBUG_UNALIGNED_STORE
+ vec_vsx_st(from, 0, to);
}
template<> EIGEN_STRONG_INLINE void prefetch<double>(const double* addr) { EIGEN_PPC_PREFETCH(addr); }
diff --git a/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h b/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h
index b021bd0b7..13351d5ec 100644
--- a/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h
+++ b/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h
@@ -126,6 +126,51 @@ Packet plog_float(const Packet _x)
por(pselect(pos_inf_mask,cst_pos_inf,x), invalid_mask));
}
+/** \internal \returns log(1 + x) computed using W. Kahan's formula.
+ See: http://www.plunk.org/~hatch/rightway.php
+ */
+template<typename Packet>
+Packet generic_plog1p(const Packet& x)
+{
+ typedef typename unpacket_traits<Packet>::type ScalarType;
+ const Packet one = pset1<Packet>(ScalarType(1));
+ Packet xp1 = padd(x, one);
+ Packet small_mask = pcmp_eq(xp1, one);
+ Packet log1 = plog(xp1);
+ Packet inf_mask = pcmp_eq(xp1, log1);
+ Packet log_large = pmul(x, pdiv(log1, psub(xp1, one)));
+ return pselect(por(small_mask, inf_mask), x, log_large);
+}
+
+/** \internal \returns exp(x)-1 computed using W. Kahan's formula.
+ See: http://www.plunk.org/~hatch/rightway.php
+ */
+template<typename Packet>
+Packet generic_expm1(const Packet& x)
+{
+ typedef typename unpacket_traits<Packet>::type ScalarType;
+ const Packet one = pset1<Packet>(ScalarType(1));
+ const Packet neg_one = pset1<Packet>(ScalarType(-1));
+ Packet u = pexp(x);
+ Packet one_mask = pcmp_eq(u, one);
+ Packet u_minus_one = psub(u, one);
+ Packet neg_one_mask = pcmp_eq(u_minus_one, neg_one);
+ Packet logu = plog(u);
+ // The following comparison is to catch the case where
+ // exp(x) = +inf. It is written in this way to avoid having
+ // to form the constant +inf, which depends on the packet
+ // type.
+ Packet pos_inf_mask = pcmp_eq(logu, u);
+ Packet expm1 = pmul(u_minus_one, pdiv(x, logu));
+ expm1 = pselect(pos_inf_mask, u, expm1);
+ return pselect(one_mask,
+ x,
+ pselect(neg_one_mask,
+ neg_one,
+ expm1));
+}
+
+
// Exponential function. Works by writing "x = m*log(2) + r" where
// "m = floor(x/log(2)+1/2)" and "r" is the remainder. The result is then
// "exp(x) = 2^m*exp(r)" where exp(r) is in the range [-1,1).
diff --git a/Eigen/src/Core/arch/GPU/Half.h b/Eigen/src/Core/arch/Default/Half.h
index 655dc20d5..56782b340 100644
--- a/Eigen/src/Core/arch/GPU/Half.h
+++ b/Eigen/src/Core/arch/Default/Half.h
@@ -33,8 +33,8 @@
// to disk and the likes), but fast on GPUs.
-#ifndef EIGEN_HALF_GPU_H
-#define EIGEN_HALF_GPU_H
+#ifndef EIGEN_HALF_H
+#define EIGEN_HALF_H
#if __cplusplus > 199711L
#define EIGEN_EXPLICIT_CAST(tgt_type) explicit operator tgt_type()
@@ -76,7 +76,6 @@ EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC float half_to_float(__half_raw h);
struct half_base : public __half_raw {
EIGEN_DEVICE_FUNC half_base() {}
- EIGEN_DEVICE_FUNC half_base(const half_base& h) : __half_raw(h) {}
EIGEN_DEVICE_FUNC half_base(const __half_raw& h) : __half_raw(h) {}
#if defined(EIGEN_HAS_GPU_FP16)
@@ -114,8 +113,7 @@ struct half : public half_impl::half_base {
EIGEN_DEVICE_FUNC half() {}
EIGEN_DEVICE_FUNC half(const __half_raw& h) : half_impl::half_base(h) {}
- EIGEN_DEVICE_FUNC half(const half& h) : half_impl::half_base(h) {}
-
+
#if defined(EIGEN_HAS_GPU_FP16)
#if defined(EIGEN_HAS_HIP_FP16)
EIGEN_DEVICE_FUNC half(const __half& h) : half_impl::half_base(h) {}
@@ -125,7 +123,7 @@ struct half : public half_impl::half_base {
#endif
#endif
#endif
-
+
explicit EIGEN_DEVICE_FUNC half(bool b)
: half_impl::half_base(half_impl::raw_uint16_to_half(b ? 0x3c00 : 0)) {}
@@ -175,12 +173,6 @@ struct half : public half_impl::half_base {
EIGEN_DEVICE_FUNC EIGEN_EXPLICIT_CAST(double) const {
return static_cast<double>(half_impl::half_to_float(*this));
}
-
- EIGEN_DEVICE_FUNC half& operator=(const half& other) {
- x = other.x;
- return *this;
- }
-
};
} // end namespace Eigen
@@ -761,4 +753,4 @@ bool (isfinite)(const Eigen::half& h) {
} // namespace numext
#endif
-#endif // EIGEN_HALF_GPU_H
+#endif // EIGEN_HALF_H
diff --git a/Eigen/src/Core/arch/Default/TypeCasting.h b/Eigen/src/Core/arch/Default/TypeCasting.h
new file mode 100644
index 000000000..b6df98468
--- /dev/null
+++ b/Eigen/src/Core/arch/Default/TypeCasting.h
@@ -0,0 +1,77 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra.
+//
+// Copyright (C) 2016 Benoit Steiner <benoit.steiner.goog@gmail.com>
+// Copyright (C) 2019 Rasmus Munk Larsen <rmlarsen@google.com>
+//
+// This Source Code Form is subject to the terms of the Mozilla
+// Public License v. 2.0. If a copy of the MPL was not distributed
+// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
+
+#ifndef EIGEN_GENERIC_TYPE_CASTING_H
+#define EIGEN_GENERIC_TYPE_CASTING_H
+
+namespace Eigen {
+
+namespace internal {
+
+template<>
+struct scalar_cast_op<float, Eigen::half> {
+ EIGEN_EMPTY_STRUCT_CTOR(scalar_cast_op)
+ typedef Eigen::half result_type;
+ EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half operator() (const float& a) const {
+ #if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
+ (defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIP_DEVICE_COMPILE))
+ return __float2half(a);
+ #else
+ return Eigen::half(a);
+ #endif
+ }
+};
+
+template<>
+struct functor_traits<scalar_cast_op<float, Eigen::half> >
+{ enum { Cost = NumTraits<float>::AddCost, PacketAccess = false }; };
+
+
+template<>
+struct scalar_cast_op<int, Eigen::half> {
+ EIGEN_EMPTY_STRUCT_CTOR(scalar_cast_op)
+ typedef Eigen::half result_type;
+ EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half operator() (const int& a) const {
+ #if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
+ (defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIP_DEVICE_COMPILE))
+ return __float2half(static_cast<float>(a));
+ #else
+ return Eigen::half(static_cast<float>(a));
+ #endif
+ }
+};
+
+template<>
+struct functor_traits<scalar_cast_op<int, Eigen::half> >
+{ enum { Cost = NumTraits<float>::AddCost, PacketAccess = false }; };
+
+
+template<>
+struct scalar_cast_op<Eigen::half, float> {
+ EIGEN_EMPTY_STRUCT_CTOR(scalar_cast_op)
+ typedef float result_type;
+ EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE float operator() (const Eigen::half& a) const {
+ #if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
+ (defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIP_DEVICE_COMPILE))
+ return __half2float(a);
+ #else
+ return static_cast<float>(a);
+ #endif
+ }
+};
+
+template<>
+struct functor_traits<scalar_cast_op<Eigen::half, float> >
+{ enum { Cost = NumTraits<float>::AddCost, PacketAccess = false }; };
+
+}
+}
+
+#endif // EIGEN_GENERIC_TYPE_CASTING_H
diff --git a/Eigen/src/Core/arch/GPU/PacketMath.h b/Eigen/src/Core/arch/GPU/PacketMath.h
index 6ba2990d1..bdbaa5362 100644
--- a/Eigen/src/Core/arch/GPU/PacketMath.h
+++ b/Eigen/src/Core/arch/GPU/PacketMath.h
@@ -460,6 +460,546 @@ ptranspose(PacketBlock<double2,2>& kernel) {
#endif
+// Packet math for Eigen::half
+// Most of the following operations require arch >= 3.0
+#if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDACC) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
+ (defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIPCC) && defined(EIGEN_HIP_DEVICE_COMPILE)) || \
+ (defined(EIGEN_HAS_CUDA_FP16) && defined(__clang__) && defined(__CUDA__))
+
+template<> struct is_arithmetic<half2> { enum { value = true }; };
+
+template<> struct packet_traits<Eigen::half> : default_packet_traits
+{
+ typedef half2 type;
+ typedef half2 half;
+ enum {
+ Vectorizable = 1,
+ AlignedOnScalar = 1,
+ size=2,
+ HasHalfPacket = 0,
+ HasAdd = 1,
+ HasSub = 1,
+ HasMul = 1,
+ HasDiv = 1,
+ HasSqrt = 1,
+ HasRsqrt = 1,
+ HasExp = 1,
+ HasExpm1 = 1,
+ HasLog = 1,
+ HasLog1p = 1
+ };
+};
+
+template<> struct unpacket_traits<half2> { typedef Eigen::half type; enum {size=2, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef half2 half; };
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pset1<half2>(const Eigen::half& from) {
+#if !defined(EIGEN_CUDA_ARCH) && !defined(EIGEN_HIP_DEVICE_COMPILE)
+ half2 r;
+ r.x = from;
+ r.y = from;
+ return r;
+#else
+ return __half2half2(from);
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pload<half2>(const Eigen::half* from) {
+ return *reinterpret_cast<const half2*>(from);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 ploadu<half2>(const Eigen::half* from) {
+ return __halves2half2(from[0], from[1]);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 ploaddup<half2>(const Eigen::half* from) {
+ return __halves2half2(from[0], from[0]);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void pstore<Eigen::half>(Eigen::half* to, const half2& from) {
+ *reinterpret_cast<half2*>(to) = from;
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void pstoreu<Eigen::half>(Eigen::half* to, const half2& from) {
+#if !defined(EIGEN_CUDA_ARCH) && !defined(EIGEN_HIP_DEVICE_COMPILE)
+ to[0] = from.x;
+ to[1] = from.y;
+#else
+ to[0] = __low2half(from);
+ to[1] = __high2half(from);
+#endif
+}
+
+template<>
+ EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE half2 ploadt_ro<half2, Aligned>(const Eigen::half* from) {
+
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __ldg((const half2*)from);
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 350
+ return __ldg((const half2*)from);
+#else
+ return __halves2half2(*(from+0), *(from+1));
+#endif
+
+#endif
+}
+
+template<>
+EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE half2 ploadt_ro<half2, Unaligned>(const Eigen::half* from) {
+
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __halves2half2(__ldg(from+0), __ldg(from+1));
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 350
+ return __halves2half2(__ldg(from+0), __ldg(from+1));
+#else
+ return __halves2half2(*(from+0), *(from+1));
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pgather<Eigen::half, half2>(const Eigen::half* from, Index stride) {
+ return __halves2half2(from[0*stride], from[1*stride]);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void pscatter<Eigen::half, half2>(Eigen::half* to, const half2& from, Index stride) {
+ to[stride*0] = __low2half(from);
+ to[stride*1] = __high2half(from);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half pfirst<half2>(const half2& a) {
+ return __low2half(a);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pabs<half2>(const half2& a) {
+ half a1 = __low2half(a);
+ half a2 = __high2half(a);
+ half result1 = half_impl::raw_uint16_to_half(a1.x & 0x7FFF);
+ half result2 = half_impl::raw_uint16_to_half(a2.x & 0x7FFF);
+ return __halves2half2(result1, result2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 ptrue<half2>(const half2& a) {
+ half true_half = half_impl::raw_uint16_to_half(0xffffu);
+ return pset1<half2>(true_half);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pzero<half2>(const half2& a) {
+ half false_half = half_impl::raw_uint16_to_half(0x0000u);
+ return pset1<half2>(false_half);
+}
+
+EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void
+ptranspose(PacketBlock<half2,2>& kernel) {
+ __half a1 = __low2half(kernel.packet[0]);
+ __half a2 = __high2half(kernel.packet[0]);
+ __half b1 = __low2half(kernel.packet[1]);
+ __half b2 = __high2half(kernel.packet[1]);
+ kernel.packet[0] = __halves2half2(a1, b1);
+ kernel.packet[1] = __halves2half2(a2, b2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 plset<half2>(const Eigen::half& a) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __halves2half2(a, __hadd(a, __float2half(1.0f)));
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ return __halves2half2(a, __hadd(a, __float2half(1.0f)));
+#else
+ float f = __half2float(a) + 1.0f;
+ return __halves2half2(a, __float2half(f));
+#endif
+
+#endif
+}
+
+template <>
+EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pselect<half2>(const half2& mask,
+ const half2& a,
+ const half2& b) {
+ half mask_low = __low2half(mask);
+ half mask_high = __high2half(mask);
+ half result_low = mask_low == half(0) ? __low2half(b) : __low2half(a);
+ half result_high = mask_high == half(0) ? __high2half(b) : __high2half(a);
+ return __halves2half2(result_low, result_high);
+}
+
+template <>
+EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pcmp_eq<half2>(const half2& a,
+ const half2& b) {
+ half true_half = half_impl::raw_uint16_to_half(0xffffu);
+ half false_half = half_impl::raw_uint16_to_half(0x0000u);
+ half a1 = __low2half(a);
+ half a2 = __high2half(a);
+ half b1 = __low2half(b);
+ half b2 = __high2half(b);
+ half eq1 = __half2float(a1) == __half2float(b1) ? true_half : false_half;
+ half eq2 = __half2float(a2) == __half2float(b2) ? true_half : false_half;
+ return __halves2half2(eq1, eq2);
+}
+
+template <>
+EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pand<half2>(const half2& a,
+ const half2& b) {
+ half a1 = __low2half(a);
+ half a2 = __high2half(a);
+ half b1 = __low2half(b);
+ half b2 = __high2half(b);
+ half result1 = half_impl::raw_uint16_to_half(a1.x & b1.x);
+ half result2 = half_impl::raw_uint16_to_half(a2.x & b2.x);
+ return __halves2half2(result1, result2);
+}
+
+template <>
+EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 por<half2>(const half2& a,
+ const half2& b) {
+ half a1 = __low2half(a);
+ half a2 = __high2half(a);
+ half b1 = __low2half(b);
+ half b2 = __high2half(b);
+ half result1 = half_impl::raw_uint16_to_half(a1.x | b1.x);
+ half result2 = half_impl::raw_uint16_to_half(a2.x | b2.x);
+ return __halves2half2(result1, result2);
+}
+
+template <>
+EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pxor<half2>(const half2& a,
+ const half2& b) {
+ half a1 = __low2half(a);
+ half a2 = __high2half(a);
+ half b1 = __low2half(b);
+ half b2 = __high2half(b);
+ half result1 = half_impl::raw_uint16_to_half(a1.x ^ b1.x);
+ half result2 = half_impl::raw_uint16_to_half(a2.x ^ b2.x);
+ return __halves2half2(result1, result2);
+}
+
+template <>
+EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pandnot<half2>(const half2& a,
+ const half2& b) {
+ half a1 = __low2half(a);
+ half a2 = __high2half(a);
+ half b1 = __low2half(b);
+ half b2 = __high2half(b);
+ half result1 = half_impl::raw_uint16_to_half(a1.x & ~b1.x);
+ half result2 = half_impl::raw_uint16_to_half(a2.x & ~b2.x);
+ return __halves2half2(result1, result2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 padd<half2>(const half2& a, const half2& b) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __hadd2(a, b);
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ return __hadd2(a, b);
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float b1 = __low2float(b);
+ float b2 = __high2float(b);
+ float r1 = a1 + b1;
+ float r2 = a2 + b2;
+ return __floats2half2_rn(r1, r2);
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 psub<half2>(const half2& a, const half2& b) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __hsub2(a, b);
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ return __hsub2(a, b);
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float b1 = __low2float(b);
+ float b2 = __high2float(b);
+ float r1 = a1 - b1;
+ float r2 = a2 - b2;
+ return __floats2half2_rn(r1, r2);
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pnegate(const half2& a) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __hneg2(a);
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ return __hneg2(a);
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ return __floats2half2_rn(-a1, -a2);
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pconj(const half2& a) { return a; }
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pmul<half2>(const half2& a, const half2& b) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __hmul2(a, b);
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ return __hmul2(a, b);
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float b1 = __low2float(b);
+ float b2 = __high2float(b);
+ float r1 = a1 * b1;
+ float r2 = a2 * b2;
+ return __floats2half2_rn(r1, r2);
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pmadd<half2>(const half2& a, const half2& b, const half2& c) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __hfma2(a, b, c);
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ return __hfma2(a, b, c);
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float b1 = __low2float(b);
+ float b2 = __high2float(b);
+ float c1 = __low2float(c);
+ float c2 = __high2float(c);
+ float r1 = a1 * b1 + c1;
+ float r2 = a2 * b2 + c2;
+ return __floats2half2_rn(r1, r2);
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pdiv<half2>(const half2& a, const half2& b) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __h2div(a, b);
+
+#else // EIGEN_CUDA_ARCH
+
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float b1 = __low2float(b);
+ float b2 = __high2float(b);
+ float r1 = a1 / b1;
+ float r2 = a2 / b2;
+ return __floats2half2_rn(r1, r2);
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pmin<half2>(const half2& a, const half2& b) {
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float b1 = __low2float(b);
+ float b2 = __high2float(b);
+ __half r1 = a1 < b1 ? __low2half(a) : __low2half(b);
+ __half r2 = a2 < b2 ? __high2half(a) : __high2half(b);
+ return __halves2half2(r1, r2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pmax<half2>(const half2& a, const half2& b) {
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float b1 = __low2float(b);
+ float b2 = __high2float(b);
+ __half r1 = a1 > b1 ? __low2half(a) : __low2half(b);
+ __half r2 = a2 > b2 ? __high2half(a) : __high2half(b);
+ return __halves2half2(r1, r2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half predux<half2>(const half2& a) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __hadd(__low2half(a), __high2half(a));
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ return __hadd(__low2half(a), __high2half(a));
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ return Eigen::half(__float2half(a1 + a2));
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half predux_max<half2>(const half2& a) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ __half first = __low2half(a);
+ __half second = __high2half(a);
+ return __hgt(first, second) ? first : second;
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ __half first = __low2half(a);
+ __half second = __high2half(a);
+ return __hgt(first, second) ? first : second;
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ return a1 > a2 ? __low2half(a) : __high2half(a);
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half predux_min<half2>(const half2& a) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ __half first = __low2half(a);
+ __half second = __high2half(a);
+ return __hlt(first, second) ? first : second;
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ __half first = __low2half(a);
+ __half second = __high2half(a);
+ return __hlt(first, second) ? first : second;
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ return a1 < a2 ? __low2half(a) : __high2half(a);
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half predux_mul<half2>(const half2& a) {
+#if defined(EIGEN_HIP_DEVICE_COMPILE)
+
+ return __hmul(__low2half(a), __high2half(a));
+
+#else // EIGEN_CUDA_ARCH
+
+#if EIGEN_CUDA_ARCH >= 530
+ return __hmul(__low2half(a), __high2half(a));
+#else
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ return Eigen::half(__float2half(a1 * a2));
+#endif
+
+#endif
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 plog1p<half2>(const half2& a) {
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float r1 = log1pf(a1);
+ float r2 = log1pf(a2);
+ return __floats2half2_rn(r1, r2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pexpm1<half2>(const half2& a) {
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float r1 = expm1f(a1);
+ float r2 = expm1f(a2);
+ return __floats2half2_rn(r1, r2);
+}
+
+#if (EIGEN_CUDA_SDK_VER >= 80000 && defined EIGEN_CUDA_ARCH && EIGEN_CUDA_ARCH >= 530) || \
+ defined(EIGEN_HIP_DEVICE_COMPILE)
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
+half2 plog<half2>(const half2& a) {
+ return h2log(a);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
+half2 pexp<half2>(const half2& a) {
+ return h2exp(a);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
+half2 psqrt<half2>(const half2& a) {
+ return h2sqrt(a);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
+half2 prsqrt<half2>(const half2& a) {
+ return h2rsqrt(a);
+}
+
+#else
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 plog<half2>(const half2& a) {
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float r1 = logf(a1);
+ float r2 = logf(a2);
+ return __floats2half2_rn(r1, r2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pexp<half2>(const half2& a) {
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float r1 = expf(a1);
+ float r2 = expf(a2);
+ return __floats2half2_rn(r1, r2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 psqrt<half2>(const half2& a) {
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float r1 = sqrtf(a1);
+ float r2 = sqrtf(a2);
+ return __floats2half2_rn(r1, r2);
+}
+
+template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 prsqrt<half2>(const half2& a) {
+ float a1 = __low2float(a);
+ float a2 = __high2float(a);
+ float r1 = rsqrtf(a1);
+ float r2 = rsqrtf(a2);
+ return __floats2half2_rn(r1, r2);
+}
+#endif
+
+#endif
+
} // end namespace internal
} // end namespace Eigen
diff --git a/Eigen/src/Core/arch/GPU/PacketMathHalf.h b/Eigen/src/Core/arch/GPU/PacketMathHalf.h
deleted file mode 100644
index 5e143e65a..000000000
--- a/Eigen/src/Core/arch/GPU/PacketMathHalf.h
+++ /dev/null
@@ -1,1630 +0,0 @@
-// This file is part of Eigen, a lightweight C++ template library
-// for linear algebra.
-//
-// Copyright (C) 2016 Benoit Steiner <benoit.steiner.goog@gmail.com>
-//
-// This Source Code Form is subject to the terms of the Mozilla
-// Public License v. 2.0. If a copy of the MPL was not distributed
-// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
-
-#ifndef EIGEN_PACKET_MATH_HALF_GPU_H
-#define EIGEN_PACKET_MATH_HALF_GPU_H
-
-
-namespace Eigen {
-namespace internal {
-
-// Most of the following operations require arch >= 3.0
-#if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDACC) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
- (defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIPCC) && defined(EIGEN_HIP_DEVICE_COMPILE)) || \
- (defined(EIGEN_HAS_CUDA_FP16) && defined(__clang__) && defined(__CUDA__))
-
-template<> struct is_arithmetic<half2> { enum { value = true }; };
-
-template<> struct packet_traits<Eigen::half> : default_packet_traits
-{
- typedef half2 type;
- typedef half2 half;
- enum {
- Vectorizable = 1,
- AlignedOnScalar = 1,
- size=2,
- HasHalfPacket = 0,
- HasAdd = 1,
- HasSub = 1,
- HasMul = 1,
- HasDiv = 1,
- HasSqrt = 1,
- HasRsqrt = 1,
- HasExp = 1,
- HasExpm1 = 1,
- HasLog = 1,
- HasLog1p = 1
- };
-};
-
-template<> struct unpacket_traits<half2> { typedef Eigen::half type; enum {size=2, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef half2 half; };
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pset1<half2>(const Eigen::half& from) {
-#if !defined(EIGEN_CUDA_ARCH) && !defined(EIGEN_HIP_DEVICE_COMPILE)
- half2 r;
- r.x = from;
- r.y = from;
- return r;
-#else
- return __half2half2(from);
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pload<half2>(const Eigen::half* from) {
- return *reinterpret_cast<const half2*>(from);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 ploadu<half2>(const Eigen::half* from) {
- return __halves2half2(from[0], from[1]);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 ploaddup<half2>(const Eigen::half* from) {
- return __halves2half2(from[0], from[0]);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void pstore<Eigen::half>(Eigen::half* to, const half2& from) {
- *reinterpret_cast<half2*>(to) = from;
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void pstoreu<Eigen::half>(Eigen::half* to, const half2& from) {
-#if !defined(EIGEN_CUDA_ARCH) && !defined(EIGEN_HIP_DEVICE_COMPILE)
- to[0] = from.x;
- to[1] = from.y;
-#else
- to[0] = __low2half(from);
- to[1] = __high2half(from);
-#endif
-}
-
-template<>
- EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE half2 ploadt_ro<half2, Aligned>(const Eigen::half* from) {
-
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __ldg((const half2*)from);
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 350
- return __ldg((const half2*)from);
-#else
- return __halves2half2(*(from+0), *(from+1));
-#endif
-
-#endif
-}
-
-template<>
-EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE half2 ploadt_ro<half2, Unaligned>(const Eigen::half* from) {
-
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __halves2half2(__ldg(from+0), __ldg(from+1));
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 350
- return __halves2half2(__ldg(from+0), __ldg(from+1));
-#else
- return __halves2half2(*(from+0), *(from+1));
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pgather<Eigen::half, half2>(const Eigen::half* from, Index stride) {
- return __halves2half2(from[0*stride], from[1*stride]);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void pscatter<Eigen::half, half2>(Eigen::half* to, const half2& from, Index stride) {
- to[stride*0] = __low2half(from);
- to[stride*1] = __high2half(from);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half pfirst<half2>(const half2& a) {
- return __low2half(a);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pabs<half2>(const half2& a) {
- half a1 = __low2half(a);
- half a2 = __high2half(a);
- half result1 = half_impl::raw_uint16_to_half(a1.x & 0x7FFF);
- half result2 = half_impl::raw_uint16_to_half(a2.x & 0x7FFF);
- return __halves2half2(result1, result2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 ptrue<half2>(const half2& a) {
- half true_half = half_impl::raw_uint16_to_half(0xffffu);
- return pset1<half2>(true_half);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pzero<half2>(const half2& a) {
- half false_half = half_impl::raw_uint16_to_half(0x0000u);
- return pset1<half2>(false_half);
-}
-
-EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void
-ptranspose(PacketBlock<half2,2>& kernel) {
- __half a1 = __low2half(kernel.packet[0]);
- __half a2 = __high2half(kernel.packet[0]);
- __half b1 = __low2half(kernel.packet[1]);
- __half b2 = __high2half(kernel.packet[1]);
- kernel.packet[0] = __halves2half2(a1, b1);
- kernel.packet[1] = __halves2half2(a2, b2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 plset<half2>(const Eigen::half& a) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __halves2half2(a, __hadd(a, __float2half(1.0f)));
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- return __halves2half2(a, __hadd(a, __float2half(1.0f)));
-#else
- float f = __half2float(a) + 1.0f;
- return __halves2half2(a, __float2half(f));
-#endif
-
-#endif
-}
-
-template <>
-EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pselect<half2>(const half2& mask,
- const half2& a,
- const half2& b) {
- half mask_low = __low2half(mask);
- half mask_high = __high2half(mask);
- half result_low = mask_low == half(0) ? __low2half(b) : __low2half(a);
- half result_high = mask_high == half(0) ? __high2half(b) : __high2half(a);
- return __halves2half2(result_low, result_high);
-}
-
-template <>
-EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pcmp_eq<half2>(const half2& a,
- const half2& b) {
- half true_half = half_impl::raw_uint16_to_half(0xffffu);
- half false_half = half_impl::raw_uint16_to_half(0x0000u);
- half a1 = __low2half(a);
- half a2 = __high2half(a);
- half b1 = __low2half(b);
- half b2 = __high2half(b);
- half eq1 = __half2float(a1) == __half2float(b1) ? true_half : false_half;
- half eq2 = __half2float(a2) == __half2float(b2) ? true_half : false_half;
- return __halves2half2(eq1, eq2);
-}
-
-template <>
-EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pand<half2>(const half2& a,
- const half2& b) {
- half a1 = __low2half(a);
- half a2 = __high2half(a);
- half b1 = __low2half(b);
- half b2 = __high2half(b);
- half result1 = half_impl::raw_uint16_to_half(a1.x & b1.x);
- half result2 = half_impl::raw_uint16_to_half(a2.x & b2.x);
- return __halves2half2(result1, result2);
-}
-
-template <>
-EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 por<half2>(const half2& a,
- const half2& b) {
- half a1 = __low2half(a);
- half a2 = __high2half(a);
- half b1 = __low2half(b);
- half b2 = __high2half(b);
- half result1 = half_impl::raw_uint16_to_half(a1.x | b1.x);
- half result2 = half_impl::raw_uint16_to_half(a2.x | b2.x);
- return __halves2half2(result1, result2);
-}
-
-template <>
-EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pxor<half2>(const half2& a,
- const half2& b) {
- half a1 = __low2half(a);
- half a2 = __high2half(a);
- half b1 = __low2half(b);
- half b2 = __high2half(b);
- half result1 = half_impl::raw_uint16_to_half(a1.x ^ b1.x);
- half result2 = half_impl::raw_uint16_to_half(a2.x ^ b2.x);
- return __halves2half2(result1, result2);
-}
-
-template <>
-EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pandnot<half2>(const half2& a,
- const half2& b) {
- half a1 = __low2half(a);
- half a2 = __high2half(a);
- half b1 = __low2half(b);
- half b2 = __high2half(b);
- half result1 = half_impl::raw_uint16_to_half(a1.x & ~b1.x);
- half result2 = half_impl::raw_uint16_to_half(a2.x & ~b2.x);
- return __halves2half2(result1, result2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 padd<half2>(const half2& a, const half2& b) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __hadd2(a, b);
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- return __hadd2(a, b);
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float b1 = __low2float(b);
- float b2 = __high2float(b);
- float r1 = a1 + b1;
- float r2 = a2 + b2;
- return __floats2half2_rn(r1, r2);
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 psub<half2>(const half2& a, const half2& b) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __hsub2(a, b);
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- return __hsub2(a, b);
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float b1 = __low2float(b);
- float b2 = __high2float(b);
- float r1 = a1 - b1;
- float r2 = a2 - b2;
- return __floats2half2_rn(r1, r2);
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pnegate(const half2& a) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __hneg2(a);
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- return __hneg2(a);
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- return __floats2half2_rn(-a1, -a2);
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pconj(const half2& a) { return a; }
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pmul<half2>(const half2& a, const half2& b) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __hmul2(a, b);
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- return __hmul2(a, b);
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float b1 = __low2float(b);
- float b2 = __high2float(b);
- float r1 = a1 * b1;
- float r2 = a2 * b2;
- return __floats2half2_rn(r1, r2);
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pmadd<half2>(const half2& a, const half2& b, const half2& c) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __hfma2(a, b, c);
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- return __hfma2(a, b, c);
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float b1 = __low2float(b);
- float b2 = __high2float(b);
- float c1 = __low2float(c);
- float c2 = __high2float(c);
- float r1 = a1 * b1 + c1;
- float r2 = a2 * b2 + c2;
- return __floats2half2_rn(r1, r2);
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pdiv<half2>(const half2& a, const half2& b) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __h2div(a, b);
-
-#else // EIGEN_CUDA_ARCH
-
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float b1 = __low2float(b);
- float b2 = __high2float(b);
- float r1 = a1 / b1;
- float r2 = a2 / b2;
- return __floats2half2_rn(r1, r2);
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pmin<half2>(const half2& a, const half2& b) {
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float b1 = __low2float(b);
- float b2 = __high2float(b);
- __half r1 = a1 < b1 ? __low2half(a) : __low2half(b);
- __half r2 = a2 < b2 ? __high2half(a) : __high2half(b);
- return __halves2half2(r1, r2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pmax<half2>(const half2& a, const half2& b) {
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float b1 = __low2float(b);
- float b2 = __high2float(b);
- __half r1 = a1 > b1 ? __low2half(a) : __low2half(b);
- __half r2 = a2 > b2 ? __high2half(a) : __high2half(b);
- return __halves2half2(r1, r2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half predux<half2>(const half2& a) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __hadd(__low2half(a), __high2half(a));
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- return __hadd(__low2half(a), __high2half(a));
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- return Eigen::half(__float2half(a1 + a2));
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half predux_max<half2>(const half2& a) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- __half first = __low2half(a);
- __half second = __high2half(a);
- return __hgt(first, second) ? first : second;
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- __half first = __low2half(a);
- __half second = __high2half(a);
- return __hgt(first, second) ? first : second;
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- return a1 > a2 ? __low2half(a) : __high2half(a);
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half predux_min<half2>(const half2& a) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- __half first = __low2half(a);
- __half second = __high2half(a);
- return __hlt(first, second) ? first : second;
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- __half first = __low2half(a);
- __half second = __high2half(a);
- return __hlt(first, second) ? first : second;
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- return a1 < a2 ? __low2half(a) : __high2half(a);
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half predux_mul<half2>(const half2& a) {
-#if defined(EIGEN_HIP_DEVICE_COMPILE)
-
- return __hmul(__low2half(a), __high2half(a));
-
-#else // EIGEN_CUDA_ARCH
-
-#if EIGEN_CUDA_ARCH >= 530
- return __hmul(__low2half(a), __high2half(a));
-#else
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- return Eigen::half(__float2half(a1 * a2));
-#endif
-
-#endif
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 plog1p<half2>(const half2& a) {
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float r1 = log1pf(a1);
- float r2 = log1pf(a2);
- return __floats2half2_rn(r1, r2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pexpm1<half2>(const half2& a) {
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float r1 = expm1f(a1);
- float r2 = expm1f(a2);
- return __floats2half2_rn(r1, r2);
-}
-
-#if (EIGEN_CUDA_SDK_VER >= 80000 && defined EIGEN_CUDA_ARCH && EIGEN_CUDA_ARCH >= 530) || \
- defined(EIGEN_HIP_DEVICE_COMPILE)
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
-half2 plog<half2>(const half2& a) {
- return h2log(a);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
-half2 pexp<half2>(const half2& a) {
- return h2exp(a);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
-half2 psqrt<half2>(const half2& a) {
- return h2sqrt(a);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
-half2 prsqrt<half2>(const half2& a) {
- return h2rsqrt(a);
-}
-
-#else
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 plog<half2>(const half2& a) {
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float r1 = logf(a1);
- float r2 = logf(a2);
- return __floats2half2_rn(r1, r2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pexp<half2>(const half2& a) {
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float r1 = expf(a1);
- float r2 = expf(a2);
- return __floats2half2_rn(r1, r2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 psqrt<half2>(const half2& a) {
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float r1 = sqrtf(a1);
- float r2 = sqrtf(a2);
- return __floats2half2_rn(r1, r2);
-}
-
-template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 prsqrt<half2>(const half2& a) {
- float a1 = __low2float(a);
- float a2 = __high2float(a);
- float r1 = rsqrtf(a1);
- float r2 = rsqrtf(a2);
- return __floats2half2_rn(r1, r2);
-}
-
-#endif
-
-#elif defined EIGEN_VECTORIZE_AVX512
-
-typedef struct {
- __m256i x;
-} Packet16h;
-
-
-template<> struct is_arithmetic<Packet16h> { enum { value = true }; };
-
-template <>
-struct packet_traits<half> : default_packet_traits {
- typedef Packet16h type;
- // There is no half-size packet for Packet16h.
- typedef Packet16h half;
- enum {
- Vectorizable = 1,
- AlignedOnScalar = 1,
- size = 16,
- HasHalfPacket = 0,
- HasAdd = 1,
- HasSub = 1,
- HasMul = 1,
- HasDiv = 1,
- HasNegate = 1,
- HasAbs = 0,
- HasAbs2 = 0,
- HasMin = 0,
- HasMax = 0,
- HasConj = 0,
- HasSetLinear = 0,
- HasSqrt = 0,
- HasRsqrt = 0,
- HasExp = 0,
- HasLog = 0,
- HasBlend = 0
- };
-};
-
-
-template<> struct unpacket_traits<Packet16h> { typedef Eigen::half type; enum {size=16, alignment=Aligned32, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef Packet16h half; };
-
-template<> EIGEN_STRONG_INLINE Packet16h pset1<Packet16h>(const Eigen::half& from) {
- Packet16h result;
- result.x = _mm256_set1_epi16(from.x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Eigen::half pfirst<Packet16h>(const Packet16h& from) {
- return half_impl::raw_uint16_to_half(static_cast<unsigned short>(_mm256_extract_epi16(from.x, 0)));
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pload<Packet16h>(const Eigen::half* from) {
- Packet16h result;
- result.x = _mm256_load_si256(reinterpret_cast<const __m256i*>(from));
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h ploadu<Packet16h>(const Eigen::half* from) {
- Packet16h result;
- result.x = _mm256_loadu_si256(reinterpret_cast<const __m256i*>(from));
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE void pstore<half>(Eigen::half* to, const Packet16h& from) {
- // (void*) -> workaround clang warning:
- // cast from 'Eigen::half *' to '__m256i *' increases required alignment from 2 to 32
- _mm256_store_si256((__m256i*)(void*)to, from.x);
-}
-
-template<> EIGEN_STRONG_INLINE void pstoreu<half>(Eigen::half* to, const Packet16h& from) {
- // (void*) -> workaround clang warning:
- // cast from 'Eigen::half *' to '__m256i *' increases required alignment from 2 to 32
- _mm256_storeu_si256((__m256i*)(void*)to, from.x);
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h
-ploaddup<Packet16h>(const Eigen::half* from) {
- Packet16h result;
- unsigned short a = from[0].x;
- unsigned short b = from[1].x;
- unsigned short c = from[2].x;
- unsigned short d = from[3].x;
- unsigned short e = from[4].x;
- unsigned short f = from[5].x;
- unsigned short g = from[6].x;
- unsigned short h = from[7].x;
- result.x = _mm256_set_epi16(h, h, g, g, f, f, e, e, d, d, c, c, b, b, a, a);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h
-ploadquad(const Eigen::half* from) {
- Packet16h result;
- unsigned short a = from[0].x;
- unsigned short b = from[1].x;
- unsigned short c = from[2].x;
- unsigned short d = from[3].x;
- result.x = _mm256_set_epi16(d, d, d, d, c, c, c, c, b, b, b, b, a, a, a, a);
- return result;
-}
-
-EIGEN_STRONG_INLINE Packet16f half2float(const Packet16h& a) {
-#ifdef EIGEN_HAS_FP16_C
- return _mm512_cvtph_ps(a.x);
-#else
- EIGEN_ALIGN64 half aux[16];
- pstore(aux, a);
- float f0(aux[0]);
- float f1(aux[1]);
- float f2(aux[2]);
- float f3(aux[3]);
- float f4(aux[4]);
- float f5(aux[5]);
- float f6(aux[6]);
- float f7(aux[7]);
- float f8(aux[8]);
- float f9(aux[9]);
- float fa(aux[10]);
- float fb(aux[11]);
- float fc(aux[12]);
- float fd(aux[13]);
- float fe(aux[14]);
- float ff(aux[15]);
-
- return _mm512_set_ps(
- ff, fe, fd, fc, fb, fa, f9, f8, f7, f6, f5, f4, f3, f2, f1, f0);
-#endif
-}
-
-EIGEN_STRONG_INLINE Packet16h float2half(const Packet16f& a) {
-#ifdef EIGEN_HAS_FP16_C
- Packet16h result;
- result.x = _mm512_cvtps_ph(a, _MM_FROUND_TO_NEAREST_INT|_MM_FROUND_NO_EXC);
- return result;
-#else
- EIGEN_ALIGN64 float aux[16];
- pstore(aux, a);
- half h0(aux[0]);
- half h1(aux[1]);
- half h2(aux[2]);
- half h3(aux[3]);
- half h4(aux[4]);
- half h5(aux[5]);
- half h6(aux[6]);
- half h7(aux[7]);
- half h8(aux[8]);
- half h9(aux[9]);
- half ha(aux[10]);
- half hb(aux[11]);
- half hc(aux[12]);
- half hd(aux[13]);
- half he(aux[14]);
- half hf(aux[15]);
-
- Packet16h result;
- result.x = _mm256_set_epi16(
- hf.x, he.x, hd.x, hc.x, hb.x, ha.x, h9.x, h8.x,
- h7.x, h6.x, h5.x, h4.x, h3.x, h2.x, h1.x, h0.x);
- return result;
-#endif
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pnot(const Packet16h& a) {
- Packet16h r; r.x = _mm256_xor_si256(a.x, pcmp_eq(a.x, a.x)); return r;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h ptrue(const Packet16h& a) {
- Packet16h r; r.x = Packet8i(ptrue(a.x)); return r;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h por(const Packet16h& a,const Packet16h& b) {
- // in some cases Packet8i is a wrapper around __m256i, so we need to
- // cast to Packet8i to call the correct overload.
- Packet16h r; r.x = por(Packet8i(a.x),Packet8i(b.x)); return r;
-}
-template<> EIGEN_STRONG_INLINE Packet16h pxor(const Packet16h& a,const Packet16h& b) {
- Packet16h r; r.x = pxor(Packet8i(a.x),Packet8i(b.x)); return r;
-}
-template<> EIGEN_STRONG_INLINE Packet16h pand(const Packet16h& a,const Packet16h& b) {
- Packet16h r; r.x = pand(Packet8i(a.x),Packet8i(b.x)); return r;
-}
-template<> EIGEN_STRONG_INLINE Packet16h pandnot(const Packet16h& a,const Packet16h& b) {
- Packet16h r; r.x = pandnot(Packet8i(a.x),Packet8i(b.x)); return r;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pselect(const Packet16h& mask, const Packet16h& a, const Packet16h& b) {
- Packet16h r; r.x = _mm256_blendv_epi8(b.x, a.x, mask.x); return r;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pcmp_eq(const Packet16h& a,const Packet16h& b) {
- Packet16f af = half2float(a);
- Packet16f bf = half2float(b);
- Packet16f rf = pcmp_eq(af, bf);
- // Pack the 32-bit flags into 16-bits flags.
- __m256i lo = _mm256_castps_si256(extract256<0>(rf));
- __m256i hi = _mm256_castps_si256(extract256<1>(rf));
- __m128i result_lo = _mm_packs_epi32(_mm256_extractf128_si256(lo, 0),
- _mm256_extractf128_si256(lo, 1));
- __m128i result_hi = _mm_packs_epi32(_mm256_extractf128_si256(hi, 0),
- _mm256_extractf128_si256(hi, 1));
- Packet16h result; result.x = _mm256_insertf128_si256(_mm256_castsi128_si256(result_lo), result_hi, 1);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pnegate(const Packet16h& a) {
- Packet16h sign_mask; sign_mask.x = _mm256_set1_epi16(static_cast<unsigned short>(0x8000));
- Packet16h result; result.x = _mm256_xor_si256(a.x, sign_mask.x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h padd<Packet16h>(const Packet16h& a, const Packet16h& b) {
- Packet16f af = half2float(a);
- Packet16f bf = half2float(b);
- Packet16f rf = padd(af, bf);
- return float2half(rf);
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h psub<Packet16h>(const Packet16h& a, const Packet16h& b) {
- Packet16f af = half2float(a);
- Packet16f bf = half2float(b);
- Packet16f rf = psub(af, bf);
- return float2half(rf);
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pmul<Packet16h>(const Packet16h& a, const Packet16h& b) {
- Packet16f af = half2float(a);
- Packet16f bf = half2float(b);
- Packet16f rf = pmul(af, bf);
- return float2half(rf);
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pdiv<Packet16h>(const Packet16h& a, const Packet16h& b) {
- Packet16f af = half2float(a);
- Packet16f bf = half2float(b);
- Packet16f rf = pdiv(af, bf);
- return float2half(rf);
-}
-
-template<> EIGEN_STRONG_INLINE half predux<Packet16h>(const Packet16h& from) {
- Packet16f from_float = half2float(from);
- return half(predux(from_float));
-}
-
-template<> EIGEN_STRONG_INLINE half predux_mul<Packet16h>(const Packet16h& from) {
- Packet16f from_float = half2float(from);
- return half(predux_mul(from_float));
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h preduxp<Packet16h>(const Packet16h* p) {
- Packet16f pf[16];
- pf[0] = half2float(p[0]);
- pf[1] = half2float(p[1]);
- pf[2] = half2float(p[2]);
- pf[3] = half2float(p[3]);
- pf[4] = half2float(p[4]);
- pf[5] = half2float(p[5]);
- pf[6] = half2float(p[6]);
- pf[7] = half2float(p[7]);
- pf[8] = half2float(p[8]);
- pf[9] = half2float(p[9]);
- pf[10] = half2float(p[10]);
- pf[11] = half2float(p[11]);
- pf[12] = half2float(p[12]);
- pf[13] = half2float(p[13]);
- pf[14] = half2float(p[14]);
- pf[15] = half2float(p[15]);
- Packet16f reduced = preduxp<Packet16f>(pf);
- return float2half(reduced);
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h preverse(const Packet16h& a)
-{
- __m128i m = _mm_setr_epi8(14,15,12,13,10,11,8,9,6,7,4,5,2,3,0,1);
- Packet16h res;
- res.x = _mm256_insertf128_si256(
- _mm256_castsi128_si256(_mm_shuffle_epi8(_mm256_extractf128_si256(a.x,1),m)),
- _mm_shuffle_epi8(_mm256_extractf128_si256(a.x,0),m), 1);
- return res;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pinsertfirst(const Packet16h& a, Eigen::half b)
-{
- Packet16h res;
- res.x = _mm256_insert_epi16(a.x,b.x,0);
- return res;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pinsertlast(const Packet16h& a, Eigen::half b)
-{
- Packet16h res;
- res.x = _mm256_insert_epi16(a.x,b.x,15);
- return res;
-}
-
-template<> EIGEN_STRONG_INLINE Packet16h pgather<Eigen::half, Packet16h>(const Eigen::half* from, Index stride)
-{
- Packet16h result;
- result.x = _mm256_set_epi16(
- from[15*stride].x, from[14*stride].x, from[13*stride].x, from[12*stride].x,
- from[11*stride].x, from[10*stride].x, from[9*stride].x, from[8*stride].x,
- from[7*stride].x, from[6*stride].x, from[5*stride].x, from[4*stride].x,
- from[3*stride].x, from[2*stride].x, from[1*stride].x, from[0*stride].x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE void pscatter<half, Packet16h>(half* to, const Packet16h& from, Index stride)
-{
- EIGEN_ALIGN64 half aux[16];
- pstore(aux, from);
- to[stride*0].x = aux[0].x;
- to[stride*1].x = aux[1].x;
- to[stride*2].x = aux[2].x;
- to[stride*3].x = aux[3].x;
- to[stride*4].x = aux[4].x;
- to[stride*5].x = aux[5].x;
- to[stride*6].x = aux[6].x;
- to[stride*7].x = aux[7].x;
- to[stride*8].x = aux[8].x;
- to[stride*9].x = aux[9].x;
- to[stride*10].x = aux[10].x;
- to[stride*11].x = aux[11].x;
- to[stride*12].x = aux[12].x;
- to[stride*13].x = aux[13].x;
- to[stride*14].x = aux[14].x;
- to[stride*15].x = aux[15].x;
-}
-
-EIGEN_STRONG_INLINE void
-ptranspose(PacketBlock<Packet16h,16>& kernel) {
- __m256i a = kernel.packet[0].x;
- __m256i b = kernel.packet[1].x;
- __m256i c = kernel.packet[2].x;
- __m256i d = kernel.packet[3].x;
- __m256i e = kernel.packet[4].x;
- __m256i f = kernel.packet[5].x;
- __m256i g = kernel.packet[6].x;
- __m256i h = kernel.packet[7].x;
- __m256i i = kernel.packet[8].x;
- __m256i j = kernel.packet[9].x;
- __m256i k = kernel.packet[10].x;
- __m256i l = kernel.packet[11].x;
- __m256i m = kernel.packet[12].x;
- __m256i n = kernel.packet[13].x;
- __m256i o = kernel.packet[14].x;
- __m256i p = kernel.packet[15].x;
-
- __m256i ab_07 = _mm256_unpacklo_epi16(a, b);
- __m256i cd_07 = _mm256_unpacklo_epi16(c, d);
- __m256i ef_07 = _mm256_unpacklo_epi16(e, f);
- __m256i gh_07 = _mm256_unpacklo_epi16(g, h);
- __m256i ij_07 = _mm256_unpacklo_epi16(i, j);
- __m256i kl_07 = _mm256_unpacklo_epi16(k, l);
- __m256i mn_07 = _mm256_unpacklo_epi16(m, n);
- __m256i op_07 = _mm256_unpacklo_epi16(o, p);
-
- __m256i ab_8f = _mm256_unpackhi_epi16(a, b);
- __m256i cd_8f = _mm256_unpackhi_epi16(c, d);
- __m256i ef_8f = _mm256_unpackhi_epi16(e, f);
- __m256i gh_8f = _mm256_unpackhi_epi16(g, h);
- __m256i ij_8f = _mm256_unpackhi_epi16(i, j);
- __m256i kl_8f = _mm256_unpackhi_epi16(k, l);
- __m256i mn_8f = _mm256_unpackhi_epi16(m, n);
- __m256i op_8f = _mm256_unpackhi_epi16(o, p);
-
- __m256i abcd_03 = _mm256_unpacklo_epi32(ab_07, cd_07);
- __m256i abcd_47 = _mm256_unpackhi_epi32(ab_07, cd_07);
- __m256i efgh_03 = _mm256_unpacklo_epi32(ef_07, gh_07);
- __m256i efgh_47 = _mm256_unpackhi_epi32(ef_07, gh_07);
- __m256i ijkl_03 = _mm256_unpacklo_epi32(ij_07, kl_07);
- __m256i ijkl_47 = _mm256_unpackhi_epi32(ij_07, kl_07);
- __m256i mnop_03 = _mm256_unpacklo_epi32(mn_07, op_07);
- __m256i mnop_47 = _mm256_unpackhi_epi32(mn_07, op_07);
-
- __m256i abcd_8b = _mm256_unpacklo_epi32(ab_8f, cd_8f);
- __m256i abcd_cf = _mm256_unpackhi_epi32(ab_8f, cd_8f);
- __m256i efgh_8b = _mm256_unpacklo_epi32(ef_8f, gh_8f);
- __m256i efgh_cf = _mm256_unpackhi_epi32(ef_8f, gh_8f);
- __m256i ijkl_8b = _mm256_unpacklo_epi32(ij_8f, kl_8f);
- __m256i ijkl_cf = _mm256_unpackhi_epi32(ij_8f, kl_8f);
- __m256i mnop_8b = _mm256_unpacklo_epi32(mn_8f, op_8f);
- __m256i mnop_cf = _mm256_unpackhi_epi32(mn_8f, op_8f);
-
- __m256i abcdefgh_01 = _mm256_unpacklo_epi64(abcd_03, efgh_03);
- __m256i abcdefgh_23 = _mm256_unpackhi_epi64(abcd_03, efgh_03);
- __m256i ijklmnop_01 = _mm256_unpacklo_epi64(ijkl_03, mnop_03);
- __m256i ijklmnop_23 = _mm256_unpackhi_epi64(ijkl_03, mnop_03);
- __m256i abcdefgh_45 = _mm256_unpacklo_epi64(abcd_47, efgh_47);
- __m256i abcdefgh_67 = _mm256_unpackhi_epi64(abcd_47, efgh_47);
- __m256i ijklmnop_45 = _mm256_unpacklo_epi64(ijkl_47, mnop_47);
- __m256i ijklmnop_67 = _mm256_unpackhi_epi64(ijkl_47, mnop_47);
- __m256i abcdefgh_89 = _mm256_unpacklo_epi64(abcd_8b, efgh_8b);
- __m256i abcdefgh_ab = _mm256_unpackhi_epi64(abcd_8b, efgh_8b);
- __m256i ijklmnop_89 = _mm256_unpacklo_epi64(ijkl_8b, mnop_8b);
- __m256i ijklmnop_ab = _mm256_unpackhi_epi64(ijkl_8b, mnop_8b);
- __m256i abcdefgh_cd = _mm256_unpacklo_epi64(abcd_cf, efgh_cf);
- __m256i abcdefgh_ef = _mm256_unpackhi_epi64(abcd_cf, efgh_cf);
- __m256i ijklmnop_cd = _mm256_unpacklo_epi64(ijkl_cf, mnop_cf);
- __m256i ijklmnop_ef = _mm256_unpackhi_epi64(ijkl_cf, mnop_cf);
-
- // NOTE: no unpacklo/hi instr in this case, so using permute instr.
- __m256i a_p_0 = _mm256_permute2x128_si256(abcdefgh_01, ijklmnop_01, 0x20);
- __m256i a_p_1 = _mm256_permute2x128_si256(abcdefgh_23, ijklmnop_23, 0x20);
- __m256i a_p_2 = _mm256_permute2x128_si256(abcdefgh_45, ijklmnop_45, 0x20);
- __m256i a_p_3 = _mm256_permute2x128_si256(abcdefgh_67, ijklmnop_67, 0x20);
- __m256i a_p_4 = _mm256_permute2x128_si256(abcdefgh_89, ijklmnop_89, 0x20);
- __m256i a_p_5 = _mm256_permute2x128_si256(abcdefgh_ab, ijklmnop_ab, 0x20);
- __m256i a_p_6 = _mm256_permute2x128_si256(abcdefgh_cd, ijklmnop_cd, 0x20);
- __m256i a_p_7 = _mm256_permute2x128_si256(abcdefgh_ef, ijklmnop_ef, 0x20);
- __m256i a_p_8 = _mm256_permute2x128_si256(abcdefgh_01, ijklmnop_01, 0x31);
- __m256i a_p_9 = _mm256_permute2x128_si256(abcdefgh_23, ijklmnop_23, 0x31);
- __m256i a_p_a = _mm256_permute2x128_si256(abcdefgh_45, ijklmnop_45, 0x31);
- __m256i a_p_b = _mm256_permute2x128_si256(abcdefgh_67, ijklmnop_67, 0x31);
- __m256i a_p_c = _mm256_permute2x128_si256(abcdefgh_89, ijklmnop_89, 0x31);
- __m256i a_p_d = _mm256_permute2x128_si256(abcdefgh_ab, ijklmnop_ab, 0x31);
- __m256i a_p_e = _mm256_permute2x128_si256(abcdefgh_cd, ijklmnop_cd, 0x31);
- __m256i a_p_f = _mm256_permute2x128_si256(abcdefgh_ef, ijklmnop_ef, 0x31);
-
- kernel.packet[0].x = a_p_0;
- kernel.packet[1].x = a_p_1;
- kernel.packet[2].x = a_p_2;
- kernel.packet[3].x = a_p_3;
- kernel.packet[4].x = a_p_4;
- kernel.packet[5].x = a_p_5;
- kernel.packet[6].x = a_p_6;
- kernel.packet[7].x = a_p_7;
- kernel.packet[8].x = a_p_8;
- kernel.packet[9].x = a_p_9;
- kernel.packet[10].x = a_p_a;
- kernel.packet[11].x = a_p_b;
- kernel.packet[12].x = a_p_c;
- kernel.packet[13].x = a_p_d;
- kernel.packet[14].x = a_p_e;
- kernel.packet[15].x = a_p_f;
-}
-
-EIGEN_STRONG_INLINE void
-ptranspose(PacketBlock<Packet16h,8>& kernel) {
- EIGEN_ALIGN64 half in[8][16];
- pstore<half>(in[0], kernel.packet[0]);
- pstore<half>(in[1], kernel.packet[1]);
- pstore<half>(in[2], kernel.packet[2]);
- pstore<half>(in[3], kernel.packet[3]);
- pstore<half>(in[4], kernel.packet[4]);
- pstore<half>(in[5], kernel.packet[5]);
- pstore<half>(in[6], kernel.packet[6]);
- pstore<half>(in[7], kernel.packet[7]);
-
- EIGEN_ALIGN64 half out[8][16];
-
- for (int i = 0; i < 8; ++i) {
- for (int j = 0; j < 8; ++j) {
- out[i][j] = in[j][2*i];
- }
- for (int j = 0; j < 8; ++j) {
- out[i][j+8] = in[j][2*i+1];
- }
- }
-
- kernel.packet[0] = pload<Packet16h>(out[0]);
- kernel.packet[1] = pload<Packet16h>(out[1]);
- kernel.packet[2] = pload<Packet16h>(out[2]);
- kernel.packet[3] = pload<Packet16h>(out[3]);
- kernel.packet[4] = pload<Packet16h>(out[4]);
- kernel.packet[5] = pload<Packet16h>(out[5]);
- kernel.packet[6] = pload<Packet16h>(out[6]);
- kernel.packet[7] = pload<Packet16h>(out[7]);
-}
-
-EIGEN_STRONG_INLINE void
-ptranspose(PacketBlock<Packet16h,4>& kernel) {
- EIGEN_ALIGN64 half in[4][16];
- pstore<half>(in[0], kernel.packet[0]);
- pstore<half>(in[1], kernel.packet[1]);
- pstore<half>(in[2], kernel.packet[2]);
- pstore<half>(in[3], kernel.packet[3]);
-
- EIGEN_ALIGN64 half out[4][16];
-
- for (int i = 0; i < 4; ++i) {
- for (int j = 0; j < 4; ++j) {
- out[i][j] = in[j][4*i];
- }
- for (int j = 0; j < 4; ++j) {
- out[i][j+4] = in[j][4*i+1];
- }
- for (int j = 0; j < 4; ++j) {
- out[i][j+8] = in[j][4*i+2];
- }
- for (int j = 0; j < 4; ++j) {
- out[i][j+12] = in[j][4*i+3];
- }
- }
-
- kernel.packet[0] = pload<Packet16h>(out[0]);
- kernel.packet[1] = pload<Packet16h>(out[1]);
- kernel.packet[2] = pload<Packet16h>(out[2]);
- kernel.packet[3] = pload<Packet16h>(out[3]);
-}
-
-
-#elif defined EIGEN_VECTORIZE_AVX
-
-typedef struct {
- __m128i x;
-} Packet8h;
-
-
-template<> struct is_arithmetic<Packet8h> { enum { value = true }; };
-
-template <>
-struct packet_traits<Eigen::half> : default_packet_traits {
- typedef Packet8h type;
- // There is no half-size packet for Packet8h.
- typedef Packet8h half;
- enum {
- Vectorizable = 1,
- AlignedOnScalar = 1,
- size = 8,
- HasHalfPacket = 0,
- HasAdd = 1,
- HasSub = 1,
- HasMul = 1,
- HasDiv = 1,
- HasNegate = 1,
- HasAbs = 0,
- HasAbs2 = 0,
- HasMin = 0,
- HasMax = 0,
- HasConj = 0,
- HasSetLinear = 0,
- HasSqrt = 0,
- HasRsqrt = 0,
- HasExp = 0,
- HasLog = 0,
- HasBlend = 0
- };
-};
-
-
-template<> struct unpacket_traits<Packet8h> { typedef Eigen::half type; enum {size=8, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef Packet8h half; };
-
-template<> EIGEN_STRONG_INLINE Packet8h pset1<Packet8h>(const Eigen::half& from) {
- Packet8h result;
- result.x = _mm_set1_epi16(from.x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Eigen::half pfirst<Packet8h>(const Packet8h& from) {
- return half_impl::raw_uint16_to_half(static_cast<unsigned short>(_mm_extract_epi16(from.x, 0)));
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pload<Packet8h>(const Eigen::half* from) {
- Packet8h result;
- result.x = _mm_load_si128(reinterpret_cast<const __m128i*>(from));
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h ploadu<Packet8h>(const Eigen::half* from) {
- Packet8h result;
- result.x = _mm_loadu_si128(reinterpret_cast<const __m128i*>(from));
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE void pstore<Eigen::half>(Eigen::half* to, const Packet8h& from) {
- _mm_store_si128(reinterpret_cast<__m128i*>(to), from.x);
-}
-
-template<> EIGEN_STRONG_INLINE void pstoreu<Eigen::half>(Eigen::half* to, const Packet8h& from) {
- _mm_storeu_si128(reinterpret_cast<__m128i*>(to), from.x);
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h
-ploaddup<Packet8h>(const Eigen::half* from) {
- Packet8h result;
- unsigned short a = from[0].x;
- unsigned short b = from[1].x;
- unsigned short c = from[2].x;
- unsigned short d = from[3].x;
- result.x = _mm_set_epi16(d, d, c, c, b, b, a, a);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h
-ploadquad<Packet8h>(const Eigen::half* from) {
- Packet8h result;
- unsigned short a = from[0].x;
- unsigned short b = from[1].x;
- result.x = _mm_set_epi16(b, b, b, b, a, a, a, a);
- return result;
-}
-
-EIGEN_STRONG_INLINE Packet8f half2float(const Packet8h& a) {
-#ifdef EIGEN_HAS_FP16_C
- return _mm256_cvtph_ps(a.x);
-#else
- EIGEN_ALIGN32 Eigen::half aux[8];
- pstore(aux, a);
- float f0(aux[0]);
- float f1(aux[1]);
- float f2(aux[2]);
- float f3(aux[3]);
- float f4(aux[4]);
- float f5(aux[5]);
- float f6(aux[6]);
- float f7(aux[7]);
-
- return _mm256_set_ps(f7, f6, f5, f4, f3, f2, f1, f0);
-#endif
-}
-
-EIGEN_STRONG_INLINE Packet8h float2half(const Packet8f& a) {
-#ifdef EIGEN_HAS_FP16_C
- Packet8h result;
- result.x = _mm256_cvtps_ph(a, _MM_FROUND_TO_NEAREST_INT|_MM_FROUND_NO_EXC);
- return result;
-#else
- EIGEN_ALIGN32 float aux[8];
- pstore(aux, a);
- Eigen::half h0(aux[0]);
- Eigen::half h1(aux[1]);
- Eigen::half h2(aux[2]);
- Eigen::half h3(aux[3]);
- Eigen::half h4(aux[4]);
- Eigen::half h5(aux[5]);
- Eigen::half h6(aux[6]);
- Eigen::half h7(aux[7]);
-
- Packet8h result;
- result.x = _mm_set_epi16(h7.x, h6.x, h5.x, h4.x, h3.x, h2.x, h1.x, h0.x);
- return result;
-#endif
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h ptrue(const Packet8h& a) {
- Packet8h r; r.x = _mm_cmpeq_epi32(a.x, a.x); return r;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h por(const Packet8h& a,const Packet8h& b) {
- // in some cases Packet4i is a wrapper around __m128i, so we either need to
- // cast to Packet4i to directly call the intrinsics as below:
- Packet8h r; r.x = _mm_or_si128(a.x,b.x); return r;
-}
-template<> EIGEN_STRONG_INLINE Packet8h pxor(const Packet8h& a,const Packet8h& b) {
- Packet8h r; r.x = _mm_xor_si128(a.x,b.x); return r;
-}
-template<> EIGEN_STRONG_INLINE Packet8h pand(const Packet8h& a,const Packet8h& b) {
- Packet8h r; r.x = _mm_and_si128(a.x,b.x); return r;
-}
-template<> EIGEN_STRONG_INLINE Packet8h pandnot(const Packet8h& a,const Packet8h& b) {
- Packet8h r; r.x = _mm_andnot_si128(b.x,a.x); return r;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pselect(const Packet8h& mask, const Packet8h& a, const Packet8h& b) {
- Packet8h r; r.x = _mm_blendv_epi8(b.x, a.x, mask.x); return r;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pcmp_eq(const Packet8h& a,const Packet8h& b) {
- Packet8f af = half2float(a);
- Packet8f bf = half2float(b);
- Packet8f rf = pcmp_eq(af, bf);
- // Pack the 32-bit flags into 16-bits flags.
- Packet8h result; result.x = _mm_packs_epi32(_mm256_extractf128_si256(_mm256_castps_si256(rf), 0),
- _mm256_extractf128_si256(_mm256_castps_si256(rf), 1));
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pconj(const Packet8h& a) { return a; }
-
-template<> EIGEN_STRONG_INLINE Packet8h pnegate(const Packet8h& a) {
- Packet8h sign_mask; sign_mask.x = _mm_set1_epi16(static_cast<unsigned short>(0x8000));
- Packet8h result; result.x = _mm_xor_si128(a.x, sign_mask.x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h padd<Packet8h>(const Packet8h& a, const Packet8h& b) {
- Packet8f af = half2float(a);
- Packet8f bf = half2float(b);
- Packet8f rf = padd(af, bf);
- return float2half(rf);
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h psub<Packet8h>(const Packet8h& a, const Packet8h& b) {
- Packet8f af = half2float(a);
- Packet8f bf = half2float(b);
- Packet8f rf = psub(af, bf);
- return float2half(rf);
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pmul<Packet8h>(const Packet8h& a, const Packet8h& b) {
- Packet8f af = half2float(a);
- Packet8f bf = half2float(b);
- Packet8f rf = pmul(af, bf);
- return float2half(rf);
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pdiv<Packet8h>(const Packet8h& a, const Packet8h& b) {
- Packet8f af = half2float(a);
- Packet8f bf = half2float(b);
- Packet8f rf = pdiv(af, bf);
- return float2half(rf);
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pgather<Eigen::half, Packet8h>(const Eigen::half* from, Index stride)
-{
- Packet8h result;
- result.x = _mm_set_epi16(from[7*stride].x, from[6*stride].x, from[5*stride].x, from[4*stride].x, from[3*stride].x, from[2*stride].x, from[1*stride].x, from[0*stride].x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE void pscatter<Eigen::half, Packet8h>(Eigen::half* to, const Packet8h& from, Index stride)
-{
- EIGEN_ALIGN32 Eigen::half aux[8];
- pstore(aux, from);
- to[stride*0].x = aux[0].x;
- to[stride*1].x = aux[1].x;
- to[stride*2].x = aux[2].x;
- to[stride*3].x = aux[3].x;
- to[stride*4].x = aux[4].x;
- to[stride*5].x = aux[5].x;
- to[stride*6].x = aux[6].x;
- to[stride*7].x = aux[7].x;
-}
-
-template<> EIGEN_STRONG_INLINE Eigen::half predux<Packet8h>(const Packet8h& a) {
- Packet8f af = half2float(a);
- float reduced = predux<Packet8f>(af);
- return Eigen::half(reduced);
-}
-
-template<> EIGEN_STRONG_INLINE Eigen::half predux_max<Packet8h>(const Packet8h& a) {
- Packet8f af = half2float(a);
- float reduced = predux_max<Packet8f>(af);
- return Eigen::half(reduced);
-}
-
-template<> EIGEN_STRONG_INLINE Eigen::half predux_min<Packet8h>(const Packet8h& a) {
- Packet8f af = half2float(a);
- float reduced = predux_min<Packet8f>(af);
- return Eigen::half(reduced);
-}
-
-template<> EIGEN_STRONG_INLINE Eigen::half predux_mul<Packet8h>(const Packet8h& a) {
- Packet8f af = half2float(a);
- float reduced = predux_mul<Packet8f>(af);
- return Eigen::half(reduced);
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h preduxp<Packet8h>(const Packet8h* p) {
- Packet8f pf[8];
- pf[0] = half2float(p[0]);
- pf[1] = half2float(p[1]);
- pf[2] = half2float(p[2]);
- pf[3] = half2float(p[3]);
- pf[4] = half2float(p[4]);
- pf[5] = half2float(p[5]);
- pf[6] = half2float(p[6]);
- pf[7] = half2float(p[7]);
- Packet8f reduced = preduxp<Packet8f>(pf);
- return float2half(reduced);
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h preverse(const Packet8h& a)
-{
- __m128i m = _mm_setr_epi8(14,15,12,13,10,11,8,9,6,7,4,5,2,3,0,1);
- Packet8h res;
- res.x = _mm_shuffle_epi8(a.x,m);
- return res;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pinsertfirst(const Packet8h& a, Eigen::half b)
-{
- Packet8h res;
- res.x = _mm_insert_epi16(a.x,int(b.x),0);
- return res;
-}
-
-template<> EIGEN_STRONG_INLINE Packet8h pinsertlast(const Packet8h& a, Eigen::half b)
-{
- Packet8h res;
- res.x = _mm_insert_epi16(a.x,int(b.x),7);
- return res;
-}
-
-template<int Offset>
-struct palign_impl<Offset,Packet8h>
-{
- static EIGEN_STRONG_INLINE void run(Packet8h& first, const Packet8h& second)
- {
- if (Offset!=0)
- first.x = _mm_alignr_epi8(second.x,first.x, Offset*2);
- }
-};
-
-EIGEN_STRONG_INLINE void
-ptranspose(PacketBlock<Packet8h,8>& kernel) {
- __m128i a = kernel.packet[0].x;
- __m128i b = kernel.packet[1].x;
- __m128i c = kernel.packet[2].x;
- __m128i d = kernel.packet[3].x;
- __m128i e = kernel.packet[4].x;
- __m128i f = kernel.packet[5].x;
- __m128i g = kernel.packet[6].x;
- __m128i h = kernel.packet[7].x;
-
- __m128i a03b03 = _mm_unpacklo_epi16(a, b);
- __m128i c03d03 = _mm_unpacklo_epi16(c, d);
- __m128i e03f03 = _mm_unpacklo_epi16(e, f);
- __m128i g03h03 = _mm_unpacklo_epi16(g, h);
- __m128i a47b47 = _mm_unpackhi_epi16(a, b);
- __m128i c47d47 = _mm_unpackhi_epi16(c, d);
- __m128i e47f47 = _mm_unpackhi_epi16(e, f);
- __m128i g47h47 = _mm_unpackhi_epi16(g, h);
-
- __m128i a01b01c01d01 = _mm_unpacklo_epi32(a03b03, c03d03);
- __m128i a23b23c23d23 = _mm_unpackhi_epi32(a03b03, c03d03);
- __m128i e01f01g01h01 = _mm_unpacklo_epi32(e03f03, g03h03);
- __m128i e23f23g23h23 = _mm_unpackhi_epi32(e03f03, g03h03);
- __m128i a45b45c45d45 = _mm_unpacklo_epi32(a47b47, c47d47);
- __m128i a67b67c67d67 = _mm_unpackhi_epi32(a47b47, c47d47);
- __m128i e45f45g45h45 = _mm_unpacklo_epi32(e47f47, g47h47);
- __m128i e67f67g67h67 = _mm_unpackhi_epi32(e47f47, g47h47);
-
- __m128i a0b0c0d0e0f0g0h0 = _mm_unpacklo_epi64(a01b01c01d01, e01f01g01h01);
- __m128i a1b1c1d1e1f1g1h1 = _mm_unpackhi_epi64(a01b01c01d01, e01f01g01h01);
- __m128i a2b2c2d2e2f2g2h2 = _mm_unpacklo_epi64(a23b23c23d23, e23f23g23h23);
- __m128i a3b3c3d3e3f3g3h3 = _mm_unpackhi_epi64(a23b23c23d23, e23f23g23h23);
- __m128i a4b4c4d4e4f4g4h4 = _mm_unpacklo_epi64(a45b45c45d45, e45f45g45h45);
- __m128i a5b5c5d5e5f5g5h5 = _mm_unpackhi_epi64(a45b45c45d45, e45f45g45h45);
- __m128i a6b6c6d6e6f6g6h6 = _mm_unpacklo_epi64(a67b67c67d67, e67f67g67h67);
- __m128i a7b7c7d7e7f7g7h7 = _mm_unpackhi_epi64(a67b67c67d67, e67f67g67h67);
-
- kernel.packet[0].x = a0b0c0d0e0f0g0h0;
- kernel.packet[1].x = a1b1c1d1e1f1g1h1;
- kernel.packet[2].x = a2b2c2d2e2f2g2h2;
- kernel.packet[3].x = a3b3c3d3e3f3g3h3;
- kernel.packet[4].x = a4b4c4d4e4f4g4h4;
- kernel.packet[5].x = a5b5c5d5e5f5g5h5;
- kernel.packet[6].x = a6b6c6d6e6f6g6h6;
- kernel.packet[7].x = a7b7c7d7e7f7g7h7;
-}
-
-EIGEN_STRONG_INLINE void
-ptranspose(PacketBlock<Packet8h,4>& kernel) {
- EIGEN_ALIGN32 Eigen::half in[4][8];
- pstore<Eigen::half>(in[0], kernel.packet[0]);
- pstore<Eigen::half>(in[1], kernel.packet[1]);
- pstore<Eigen::half>(in[2], kernel.packet[2]);
- pstore<Eigen::half>(in[3], kernel.packet[3]);
-
- EIGEN_ALIGN32 Eigen::half out[4][8];
-
- for (int i = 0; i < 4; ++i) {
- for (int j = 0; j < 4; ++j) {
- out[i][j] = in[j][2*i];
- }
- for (int j = 0; j < 4; ++j) {
- out[i][j+4] = in[j][2*i+1];
- }
- }
-
- kernel.packet[0] = pload<Packet8h>(out[0]);
- kernel.packet[1] = pload<Packet8h>(out[1]);
- kernel.packet[2] = pload<Packet8h>(out[2]);
- kernel.packet[3] = pload<Packet8h>(out[3]);
-}
-
-
-// Disable the following code since it's broken on too many platforms / compilers.
-//#elif defined(EIGEN_VECTORIZE_SSE) && (!EIGEN_ARCH_x86_64) && (!EIGEN_COMP_MSVC)
-#elif 0
-
-typedef struct {
- __m64 x;
-} Packet4h;
-
-
-template<> struct is_arithmetic<Packet4h> { enum { value = true }; };
-
-template <>
-struct packet_traits<Eigen::half> : default_packet_traits {
- typedef Packet4h type;
- // There is no half-size packet for Packet4h.
- typedef Packet4h half;
- enum {
- Vectorizable = 1,
- AlignedOnScalar = 1,
- size = 4,
- HasHalfPacket = 0,
- HasAdd = 1,
- HasSub = 1,
- HasMul = 1,
- HasDiv = 1,
- HasNegate = 0,
- HasAbs = 0,
- HasAbs2 = 0,
- HasMin = 0,
- HasMax = 0,
- HasConj = 0,
- HasSetLinear = 0,
- HasSqrt = 0,
- HasRsqrt = 0,
- HasExp = 0,
- HasLog = 0,
- HasBlend = 0
- };
-};
-
-
-template<> struct unpacket_traits<Packet4h> { typedef Eigen::half type; enum {size=4, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef Packet4h half; };
-
-template<> EIGEN_STRONG_INLINE Packet4h pset1<Packet4h>(const Eigen::half& from) {
- Packet4h result;
- result.x = _mm_set1_pi16(from.x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Eigen::half pfirst<Packet4h>(const Packet4h& from) {
- return half_impl::raw_uint16_to_half(static_cast<unsigned short>(_mm_cvtsi64_si32(from.x)));
-}
-
-template<> EIGEN_STRONG_INLINE Packet4h pconj(const Packet4h& a) { return a; }
-
-template<> EIGEN_STRONG_INLINE Packet4h padd<Packet4h>(const Packet4h& a, const Packet4h& b) {
- __int64_t a64 = _mm_cvtm64_si64(a.x);
- __int64_t b64 = _mm_cvtm64_si64(b.x);
-
- Eigen::half h[4];
-
- Eigen::half ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64));
- Eigen::half hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64));
- h[0] = ha + hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 16));
- h[1] = ha + hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 32));
- h[2] = ha + hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 48));
- h[3] = ha + hb;
- Packet4h result;
- result.x = _mm_set_pi16(h[3].x, h[2].x, h[1].x, h[0].x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet4h psub<Packet4h>(const Packet4h& a, const Packet4h& b) {
- __int64_t a64 = _mm_cvtm64_si64(a.x);
- __int64_t b64 = _mm_cvtm64_si64(b.x);
-
- Eigen::half h[4];
-
- Eigen::half ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64));
- Eigen::half hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64));
- h[0] = ha - hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 16));
- h[1] = ha - hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 32));
- h[2] = ha - hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 48));
- h[3] = ha - hb;
- Packet4h result;
- result.x = _mm_set_pi16(h[3].x, h[2].x, h[1].x, h[0].x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet4h pmul<Packet4h>(const Packet4h& a, const Packet4h& b) {
- __int64_t a64 = _mm_cvtm64_si64(a.x);
- __int64_t b64 = _mm_cvtm64_si64(b.x);
-
- Eigen::half h[4];
-
- Eigen::half ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64));
- Eigen::half hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64));
- h[0] = ha * hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 16));
- h[1] = ha * hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 32));
- h[2] = ha * hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 48));
- h[3] = ha * hb;
- Packet4h result;
- result.x = _mm_set_pi16(h[3].x, h[2].x, h[1].x, h[0].x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet4h pdiv<Packet4h>(const Packet4h& a, const Packet4h& b) {
- __int64_t a64 = _mm_cvtm64_si64(a.x);
- __int64_t b64 = _mm_cvtm64_si64(b.x);
-
- Eigen::half h[4];
-
- Eigen::half ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64));
- Eigen::half hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64));
- h[0] = ha / hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 16));
- h[1] = ha / hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 32));
- h[2] = ha / hb;
- ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
- hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 48));
- h[3] = ha / hb;
- Packet4h result;
- result.x = _mm_set_pi16(h[3].x, h[2].x, h[1].x, h[0].x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet4h pload<Packet4h>(const Eigen::half* from) {
- Packet4h result;
- result.x = _mm_cvtsi64_m64(*reinterpret_cast<const __int64_t*>(from));
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE Packet4h ploadu<Packet4h>(const Eigen::half* from) {
- Packet4h result;
- result.x = _mm_cvtsi64_m64(*reinterpret_cast<const __int64_t*>(from));
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE void pstore<Eigen::half>(Eigen::half* to, const Packet4h& from) {
- __int64_t r = _mm_cvtm64_si64(from.x);
- *(reinterpret_cast<__int64_t*>(to)) = r;
-}
-
-template<> EIGEN_STRONG_INLINE void pstoreu<Eigen::half>(Eigen::half* to, const Packet4h& from) {
- __int64_t r = _mm_cvtm64_si64(from.x);
- *(reinterpret_cast<__int64_t*>(to)) = r;
-}
-
-template<> EIGEN_STRONG_INLINE Packet4h
-ploadquad<Packet4h>(const Eigen::half* from) {
- return pset1<Packet4h>(*from);
-}
-
-template<> EIGEN_STRONG_INLINE Packet4h pgather<Eigen::half, Packet4h>(const Eigen::half* from, Index stride)
-{
- Packet4h result;
- result.x = _mm_set_pi16(from[3*stride].x, from[2*stride].x, from[1*stride].x, from[0*stride].x);
- return result;
-}
-
-template<> EIGEN_STRONG_INLINE void pscatter<Eigen::half, Packet4h>(Eigen::half* to, const Packet4h& from, Index stride)
-{
- __int64_t a = _mm_cvtm64_si64(from.x);
- to[stride*0].x = static_cast<unsigned short>(a);
- to[stride*1].x = static_cast<unsigned short>(a >> 16);
- to[stride*2].x = static_cast<unsigned short>(a >> 32);
- to[stride*3].x = static_cast<unsigned short>(a >> 48);
-}
-
-EIGEN_STRONG_INLINE void
-ptranspose(PacketBlock<Packet4h,4>& kernel) {
- __m64 T0 = _mm_unpacklo_pi16(kernel.packet[0].x, kernel.packet[1].x);
- __m64 T1 = _mm_unpacklo_pi16(kernel.packet[2].x, kernel.packet[3].x);
- __m64 T2 = _mm_unpackhi_pi16(kernel.packet[0].x, kernel.packet[1].x);
- __m64 T3 = _mm_unpackhi_pi16(kernel.packet[2].x, kernel.packet[3].x);
-
- kernel.packet[0].x = _mm_unpacklo_pi32(T0, T1);
- kernel.packet[1].x = _mm_unpackhi_pi32(T0, T1);
- kernel.packet[2].x = _mm_unpacklo_pi32(T2, T3);
- kernel.packet[3].x = _mm_unpackhi_pi32(T2, T3);
-}
-
-#endif
-
-}
-}
-
-#endif // EIGEN_PACKET_MATH_HALF_GPU_H
diff --git a/Eigen/src/Core/arch/GPU/TypeCasting.h b/Eigen/src/Core/arch/GPU/TypeCasting.h
index 57a55d08b..c278f3fe8 100644
--- a/Eigen/src/Core/arch/GPU/TypeCasting.h
+++ b/Eigen/src/Core/arch/GPU/TypeCasting.h
@@ -14,64 +14,6 @@ namespace Eigen {
namespace internal {
-template<>
-struct scalar_cast_op<float, Eigen::half> {
- EIGEN_EMPTY_STRUCT_CTOR(scalar_cast_op)
- typedef Eigen::half result_type;
- EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half operator() (const float& a) const {
- #if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
- (defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIP_DEVICE_COMPILE))
- return __float2half(a);
- #else
- return Eigen::half(a);
- #endif
- }
-};
-
-template<>
-struct functor_traits<scalar_cast_op<float, Eigen::half> >
-{ enum { Cost = NumTraits<float>::AddCost, PacketAccess = false }; };
-
-
-template<>
-struct scalar_cast_op<int, Eigen::half> {
- EIGEN_EMPTY_STRUCT_CTOR(scalar_cast_op)
- typedef Eigen::half result_type;
- EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Eigen::half operator() (const int& a) const {
- #if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
- (defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIP_DEVICE_COMPILE))
- return __float2half(static_cast<float>(a));
- #else
- return Eigen::half(static_cast<float>(a));
- #endif
- }
-};
-
-template<>
-struct functor_traits<scalar_cast_op<int, Eigen::half> >
-{ enum { Cost = NumTraits<float>::AddCost, PacketAccess = false }; };
-
-
-template<>
-struct scalar_cast_op<Eigen::half, float> {
- EIGEN_EMPTY_STRUCT_CTOR(scalar_cast_op)
- typedef float result_type;
- EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE float operator() (const Eigen::half& a) const {
- #if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
- (defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIP_DEVICE_COMPILE))
- return __half2float(a);
- #else
- return static_cast<float>(a);
- #endif
- }
-};
-
-template<>
-struct functor_traits<scalar_cast_op<Eigen::half, float> >
-{ enum { Cost = NumTraits<float>::AddCost, PacketAccess = false }; };
-
-
-
#if (defined(EIGEN_HAS_CUDA_FP16) && defined(EIGEN_CUDA_ARCH) && EIGEN_CUDA_ARCH >= 300) || \
(defined(EIGEN_HAS_HIP_FP16) && defined(EIGEN_HIP_DEVICE_COMPILE))
@@ -104,109 +46,6 @@ template<> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE half2 pcast<float4, half2>(cons
return __floats2half2_rn(a.x, a.y);
}
-#elif defined EIGEN_VECTORIZE_AVX512
-template <>
-struct type_casting_traits<half, float> {
- enum {
- VectorizedCast = 1,
- SrcCoeffRatio = 1,
- TgtCoeffRatio = 1
- };
-};
-
-template<> EIGEN_STRONG_INLINE Packet16f pcast<Packet16h, Packet16f>(const Packet16h& a) {
- return half2float(a);
-}
-
-template <>
-struct type_casting_traits<float, half> {
- enum {
- VectorizedCast = 1,
- SrcCoeffRatio = 1,
- TgtCoeffRatio = 1
- };
-};
-
-template<> EIGEN_STRONG_INLINE Packet16h pcast<Packet16f, Packet16h>(const Packet16f& a) {
- return float2half(a);
-}
-
-#elif defined EIGEN_VECTORIZE_AVX
-
-template <>
-struct type_casting_traits<Eigen::half, float> {
- enum {
- VectorizedCast = 1,
- SrcCoeffRatio = 1,
- TgtCoeffRatio = 1
- };
-};
-
-template<> EIGEN_STRONG_INLINE Packet8f pcast<Packet8h, Packet8f>(const Packet8h& a) {
- return half2float(a);
-}
-
-template <>
-struct type_casting_traits<float, Eigen::half> {
- enum {
- VectorizedCast = 1,
- SrcCoeffRatio = 1,
- TgtCoeffRatio = 1
- };
-};
-
-template<> EIGEN_STRONG_INLINE Packet8h pcast<Packet8f, Packet8h>(const Packet8f& a) {
- return float2half(a);
-}
-
-// Disable the following code since it's broken on too many platforms / compilers.
-//#elif defined(EIGEN_VECTORIZE_SSE) && (!EIGEN_ARCH_x86_64) && (!EIGEN_COMP_MSVC)
-#elif 0
-
-template <>
-struct type_casting_traits<Eigen::half, float> {
- enum {
- VectorizedCast = 1,
- SrcCoeffRatio = 1,
- TgtCoeffRatio = 1
- };
-};
-
-template<> EIGEN_STRONG_INLINE Packet4f pcast<Packet4h, Packet4f>(const Packet4h& a) {
- __int64_t a64 = _mm_cvtm64_si64(a.x);
- Eigen::half h = raw_uint16_to_half(static_cast<unsigned short>(a64));
- float f1 = static_cast<float>(h);
- h = raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
- float f2 = static_cast<float>(h);
- h = raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
- float f3 = static_cast<float>(h);
- h = raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
- float f4 = static_cast<float>(h);
- return _mm_set_ps(f4, f3, f2, f1);
-}
-
-template <>
-struct type_casting_traits<float, Eigen::half> {
- enum {
- VectorizedCast = 1,
- SrcCoeffRatio = 1,
- TgtCoeffRatio = 1
- };
-};
-
-template<> EIGEN_STRONG_INLINE Packet4h pcast<Packet4f, Packet4h>(const Packet4f& a) {
- EIGEN_ALIGN16 float aux[4];
- pstore(aux, a);
- Eigen::half h0(aux[0]);
- Eigen::half h1(aux[1]);
- Eigen::half h2(aux[2]);
- Eigen::half h3(aux[3]);
-
- Packet4h result;
- result.x = _mm_set_pi16(h3.x, h2.x, h1.x, h0.x);
- return result;
-}
-
#endif
} // end namespace internal
diff --git a/Eigen/src/Core/arch/SSE/MathFunctions.h b/Eigen/src/Core/arch/SSE/MathFunctions.h
index 0d491ab88..02c8f3c2f 100644
--- a/Eigen/src/Core/arch/SSE/MathFunctions.h
+++ b/Eigen/src/Core/arch/SSE/MathFunctions.h
@@ -22,12 +22,21 @@ namespace Eigen {
namespace internal {
template<> EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED
-Packet4f plog<Packet4f>(const Packet4f& _x)
-{
+Packet4f plog<Packet4f>(const Packet4f& _x) {
return plog_float(_x);
}
template<> EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED
+Packet4f plog1p<Packet4f>(const Packet4f& _x) {
+ return generic_plog1p(_x);
+}
+
+template<> EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED
+Packet4f pexpm1<Packet4f>(const Packet4f& _x) {
+ return generic_expm1(_x);
+}
+
+template<> EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED
Packet4f pexp<Packet4f>(const Packet4f& _x)
{
return pexp_float(_x);
diff --git a/Eigen/src/Core/arch/SSE/PacketMath.h b/Eigen/src/Core/arch/SSE/PacketMath.h
index 69daea8f7..5da8ff5f4 100755
--- a/Eigen/src/Core/arch/SSE/PacketMath.h
+++ b/Eigen/src/Core/arch/SSE/PacketMath.h
@@ -110,9 +110,9 @@ template<> struct packet_traits<float> : default_packet_traits
HasSin = EIGEN_FAST_MATH,
HasCos = EIGEN_FAST_MATH,
HasLog = 1,
+ HasLog1p = 1,
+ HasExpm1 = 1,
HasExp = 1,
- HasLog1p = 1,
- HasExpm1 = 1,
HasNdtri = 1,
HasSqrt = 1,
HasRsqrt = 1,
@@ -1056,6 +1056,214 @@ template<> EIGEN_STRONG_INLINE double pmadd(const double& a, const double& b, co
}
#endif
+
+// Packet math for Eigen::half
+// Disable the following code since it's broken on too many platforms / compilers.
+//#elif defined(EIGEN_VECTORIZE_SSE) && (!EIGEN_ARCH_x86_64) && (!EIGEN_COMP_MSVC)
+#if 0
+
+typedef struct {
+ __m64 x;
+} Packet4h;
+
+
+template<> struct is_arithmetic<Packet4h> { enum { value = true }; };
+
+template <>
+struct packet_traits<Eigen::half> : default_packet_traits {
+ typedef Packet4h type;
+ // There is no half-size packet for Packet4h.
+ typedef Packet4h half;
+ enum {
+ Vectorizable = 1,
+ AlignedOnScalar = 1,
+ size = 4,
+ HasHalfPacket = 0,
+ HasAdd = 1,
+ HasSub = 1,
+ HasMul = 1,
+ HasDiv = 1,
+ HasNegate = 0,
+ HasAbs = 0,
+ HasAbs2 = 0,
+ HasMin = 0,
+ HasMax = 0,
+ HasConj = 0,
+ HasSetLinear = 0,
+ HasSqrt = 0,
+ HasRsqrt = 0,
+ HasExp = 0,
+ HasLog = 0,
+ HasBlend = 0
+ };
+};
+
+
+template<> struct unpacket_traits<Packet4h> { typedef Eigen::half type; enum {size=4, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef Packet4h half; };
+
+template<> EIGEN_STRONG_INLINE Packet4h pset1<Packet4h>(const Eigen::half& from) {
+ Packet4h result;
+ result.x = _mm_set1_pi16(from.x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Eigen::half pfirst<Packet4h>(const Packet4h& from) {
+ return half_impl::raw_uint16_to_half(static_cast<unsigned short>(_mm_cvtsi64_si32(from.x)));
+}
+
+template<> EIGEN_STRONG_INLINE Packet4h pconj(const Packet4h& a) { return a; }
+
+template<> EIGEN_STRONG_INLINE Packet4h padd<Packet4h>(const Packet4h& a, const Packet4h& b) {
+ __int64_t a64 = _mm_cvtm64_si64(a.x);
+ __int64_t b64 = _mm_cvtm64_si64(b.x);
+
+ Eigen::half h[4];
+
+ Eigen::half ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64));
+ Eigen::half hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64));
+ h[0] = ha + hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 16));
+ h[1] = ha + hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 32));
+ h[2] = ha + hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 48));
+ h[3] = ha + hb;
+ Packet4h result;
+ result.x = _mm_set_pi16(h[3].x, h[2].x, h[1].x, h[0].x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet4h psub<Packet4h>(const Packet4h& a, const Packet4h& b) {
+ __int64_t a64 = _mm_cvtm64_si64(a.x);
+ __int64_t b64 = _mm_cvtm64_si64(b.x);
+
+ Eigen::half h[4];
+
+ Eigen::half ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64));
+ Eigen::half hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64));
+ h[0] = ha - hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 16));
+ h[1] = ha - hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 32));
+ h[2] = ha - hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 48));
+ h[3] = ha - hb;
+ Packet4h result;
+ result.x = _mm_set_pi16(h[3].x, h[2].x, h[1].x, h[0].x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet4h pmul<Packet4h>(const Packet4h& a, const Packet4h& b) {
+ __int64_t a64 = _mm_cvtm64_si64(a.x);
+ __int64_t b64 = _mm_cvtm64_si64(b.x);
+
+ Eigen::half h[4];
+
+ Eigen::half ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64));
+ Eigen::half hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64));
+ h[0] = ha * hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 16));
+ h[1] = ha * hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 32));
+ h[2] = ha * hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 48));
+ h[3] = ha * hb;
+ Packet4h result;
+ result.x = _mm_set_pi16(h[3].x, h[2].x, h[1].x, h[0].x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet4h pdiv<Packet4h>(const Packet4h& a, const Packet4h& b) {
+ __int64_t a64 = _mm_cvtm64_si64(a.x);
+ __int64_t b64 = _mm_cvtm64_si64(b.x);
+
+ Eigen::half h[4];
+
+ Eigen::half ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64));
+ Eigen::half hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64));
+ h[0] = ha / hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 16));
+ h[1] = ha / hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 32));
+ h[2] = ha / hb;
+ ha = half_impl::raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
+ hb = half_impl::raw_uint16_to_half(static_cast<unsigned short>(b64 >> 48));
+ h[3] = ha / hb;
+ Packet4h result;
+ result.x = _mm_set_pi16(h[3].x, h[2].x, h[1].x, h[0].x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet4h pload<Packet4h>(const Eigen::half* from) {
+ Packet4h result;
+ result.x = _mm_cvtsi64_m64(*reinterpret_cast<const __int64_t*>(from));
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE Packet4h ploadu<Packet4h>(const Eigen::half* from) {
+ Packet4h result;
+ result.x = _mm_cvtsi64_m64(*reinterpret_cast<const __int64_t*>(from));
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE void pstore<Eigen::half>(Eigen::half* to, const Packet4h& from) {
+ __int64_t r = _mm_cvtm64_si64(from.x);
+ *(reinterpret_cast<__int64_t*>(to)) = r;
+}
+
+template<> EIGEN_STRONG_INLINE void pstoreu<Eigen::half>(Eigen::half* to, const Packet4h& from) {
+ __int64_t r = _mm_cvtm64_si64(from.x);
+ *(reinterpret_cast<__int64_t*>(to)) = r;
+}
+
+template<> EIGEN_STRONG_INLINE Packet4h
+ploadquad<Packet4h>(const Eigen::half* from) {
+ return pset1<Packet4h>(*from);
+}
+
+template<> EIGEN_STRONG_INLINE Packet4h pgather<Eigen::half, Packet4h>(const Eigen::half* from, Index stride)
+{
+ Packet4h result;
+ result.x = _mm_set_pi16(from[3*stride].x, from[2*stride].x, from[1*stride].x, from[0*stride].x);
+ return result;
+}
+
+template<> EIGEN_STRONG_INLINE void pscatter<Eigen::half, Packet4h>(Eigen::half* to, const Packet4h& from, Index stride)
+{
+ __int64_t a = _mm_cvtm64_si64(from.x);
+ to[stride*0].x = static_cast<unsigned short>(a);
+ to[stride*1].x = static_cast<unsigned short>(a >> 16);
+ to[stride*2].x = static_cast<unsigned short>(a >> 32);
+ to[stride*3].x = static_cast<unsigned short>(a >> 48);
+}
+
+EIGEN_STRONG_INLINE void
+ptranspose(PacketBlock<Packet4h,4>& kernel) {
+ __m64 T0 = _mm_unpacklo_pi16(kernel.packet[0].x, kernel.packet[1].x);
+ __m64 T1 = _mm_unpacklo_pi16(kernel.packet[2].x, kernel.packet[3].x);
+ __m64 T2 = _mm_unpackhi_pi16(kernel.packet[0].x, kernel.packet[1].x);
+ __m64 T3 = _mm_unpackhi_pi16(kernel.packet[2].x, kernel.packet[3].x);
+
+ kernel.packet[0].x = _mm_unpacklo_pi32(T0, T1);
+ kernel.packet[1].x = _mm_unpackhi_pi32(T0, T1);
+ kernel.packet[2].x = _mm_unpacklo_pi32(T2, T3);
+ kernel.packet[3].x = _mm_unpackhi_pi32(T2, T3);
+}
+
+#endif
+
+
} // end namespace internal
} // end namespace Eigen
diff --git a/Eigen/src/Core/arch/SSE/TypeCasting.h b/Eigen/src/Core/arch/SSE/TypeCasting.h
index f607366f0..1b8e9a550 100644
--- a/Eigen/src/Core/arch/SSE/TypeCasting.h
+++ b/Eigen/src/Core/arch/SSE/TypeCasting.h
@@ -77,6 +77,57 @@ template<> EIGEN_STRONG_INLINE Packet4f preinterpret<Packet4f,Packet4i>(const Pa
return _mm_castsi128_ps(a);
}
+
+// Disable the following code since it's broken on too many platforms / compilers.
+//#elif defined(EIGEN_VECTORIZE_SSE) && (!EIGEN_ARCH_x86_64) && (!EIGEN_COMP_MSVC)
+#if 0
+
+template <>
+struct type_casting_traits<Eigen::half, float> {
+ enum {
+ VectorizedCast = 1,
+ SrcCoeffRatio = 1,
+ TgtCoeffRatio = 1
+ };
+};
+
+template<> EIGEN_STRONG_INLINE Packet4f pcast<Packet4h, Packet4f>(const Packet4h& a) {
+ __int64_t a64 = _mm_cvtm64_si64(a.x);
+ Eigen::half h = raw_uint16_to_half(static_cast<unsigned short>(a64));
+ float f1 = static_cast<float>(h);
+ h = raw_uint16_to_half(static_cast<unsigned short>(a64 >> 16));
+ float f2 = static_cast<float>(h);
+ h = raw_uint16_to_half(static_cast<unsigned short>(a64 >> 32));
+ float f3 = static_cast<float>(h);
+ h = raw_uint16_to_half(static_cast<unsigned short>(a64 >> 48));
+ float f4 = static_cast<float>(h);
+ return _mm_set_ps(f4, f3, f2, f1);
+}
+
+template <>
+struct type_casting_traits<float, Eigen::half> {
+ enum {
+ VectorizedCast = 1,
+ SrcCoeffRatio = 1,
+ TgtCoeffRatio = 1
+ };
+};
+
+template<> EIGEN_STRONG_INLINE Packet4h pcast<Packet4f, Packet4h>(const Packet4f& a) {
+ EIGEN_ALIGN16 float aux[4];
+ pstore(aux, a);
+ Eigen::half h0(aux[0]);
+ Eigen::half h1(aux[1]);
+ Eigen::half h2(aux[2]);
+ Eigen::half h3(aux[3]);
+
+ Packet4h result;
+ result.x = _mm_set_pi16(h3.x, h2.x, h1.x, h0.x);
+ return result;
+}
+
+#endif
+
} // end namespace internal
} // end namespace Eigen
diff --git a/Eigen/src/Core/util/DisableStupidWarnings.h b/Eigen/src/Core/util/DisableStupidWarnings.h
index 6c7c2d655..4501d3248 100755
--- a/Eigen/src/Core/util/DisableStupidWarnings.h
+++ b/Eigen/src/Core/util/DisableStupidWarnings.h
@@ -44,6 +44,11 @@
#if __clang_major__ >= 3 && __clang_minor__ >= 5
#pragma clang diagnostic ignored "-Wabsolute-value"
#endif
+ #if ( defined(__ALTIVEC__) || defined(__VSX__) ) && __cplusplus < 201103L
+ // warning: generic selections are a C11-specific feature
+ // ignoring warnings thrown at vec_ctf in Altivec/PacketMath.h
+ #pragma clang diagnostic ignored "-Wc11-extensions"
+ #endif
#elif defined __GNUC__
diff --git a/Eigen/src/OrderingMethods/Eigen_Colamd.h b/Eigen/src/OrderingMethods/Eigen_Colamd.h
index 67fcad3f7..8e339a704 100644
--- a/Eigen/src/OrderingMethods/Eigen_Colamd.h
+++ b/Eigen/src/OrderingMethods/Eigen_Colamd.h
@@ -13,115 +13,119 @@
// Davis (davis@cise.ufl.edu), University of Florida. The algorithm was
// developed in collaboration with John Gilbert, Xerox PARC, and Esmond
// Ng, Oak Ridge National Laboratory.
-//
+//
// Date:
-//
+//
// September 8, 2003. Version 2.3.
-//
+//
// Acknowledgements:
-//
+//
// This work was supported by the National Science Foundation, under
// grants DMS-9504974 and DMS-9803599.
-//
+//
// Notice:
-//
+//
// Copyright (c) 1998-2003 by the University of Florida.
// All Rights Reserved.
-//
+//
// THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY
// EXPRESSED OR IMPLIED. ANY USE IS AT YOUR OWN RISK.
-//
+//
// Permission is hereby granted to use, copy, modify, and/or distribute
// this program, provided that the Copyright, this License, and the
// Availability of the original version is retained on all copies and made
// accessible to the end-user of any code or package that includes COLAMD
-// or any modified version of COLAMD.
-//
+// or any modified version of COLAMD.
+//
// Availability:
-//
+//
// The colamd/symamd library is available at
-//
+//
// http://www.suitesparse.com
-
+
#ifndef EIGEN_COLAMD_H
#define EIGEN_COLAMD_H
namespace internal {
+
+namespace Colamd {
+
/* Ensure that debugging is turned off: */
#ifndef COLAMD_NDEBUG
#define COLAMD_NDEBUG
#endif /* NDEBUG */
+
+
/* ========================================================================== */
/* === Knob and statistics definitions ====================================== */
/* ========================================================================== */
/* size of the knobs [ ] array. Only knobs [0..1] are currently used. */
-#define COLAMD_KNOBS 20
+const int NKnobs = 20;
/* number of output statistics. Only stats [0..6] are currently used. */
-#define COLAMD_STATS 20
+const int NStats = 20;
-/* knobs [0] and stats [0]: dense row knob and output statistic. */
-#define COLAMD_DENSE_ROW 0
+/* Indices into knobs and stats array. */
+enum KnobsStatsIndex {
+ /* knobs [0] and stats [0]: dense row knob and output statistic. */
+ DenseRow = 0,
-/* knobs [1] and stats [1]: dense column knob and output statistic. */
-#define COLAMD_DENSE_COL 1
+ /* knobs [1] and stats [1]: dense column knob and output statistic. */
+ DenseCol = 1,
-/* stats [2]: memory defragmentation count output statistic */
-#define COLAMD_DEFRAG_COUNT 2
+ /* stats [2]: memory defragmentation count output statistic */
+ DefragCount = 2,
-/* stats [3]: colamd status: zero OK, > 0 warning or notice, < 0 error */
-#define COLAMD_STATUS 3
+ /* stats [3]: colamd status: zero OK, > 0 warning or notice, < 0 error */
+ Status = 3,
-/* stats [4..6]: error info, or info on jumbled columns */
-#define COLAMD_INFO1 4
-#define COLAMD_INFO2 5
-#define COLAMD_INFO3 6
+ /* stats [4..6]: error info, or info on jumbled columns */
+ Info1 = 4,
+ Info2 = 5,
+ Info3 = 6
+};
/* error codes returned in stats [3]: */
-#define COLAMD_OK (0)
-#define COLAMD_OK_BUT_JUMBLED (1)
-#define COLAMD_ERROR_A_not_present (-1)
-#define COLAMD_ERROR_p_not_present (-2)
-#define COLAMD_ERROR_nrow_negative (-3)
-#define COLAMD_ERROR_ncol_negative (-4)
-#define COLAMD_ERROR_nnz_negative (-5)
-#define COLAMD_ERROR_p0_nonzero (-6)
-#define COLAMD_ERROR_A_too_small (-7)
-#define COLAMD_ERROR_col_length_negative (-8)
-#define COLAMD_ERROR_row_index_out_of_bounds (-9)
-#define COLAMD_ERROR_out_of_memory (-10)
-#define COLAMD_ERROR_internal_error (-999)
-
+enum Status {
+ Ok = 0,
+ OkButJumbled = 1,
+ ErrorANotPresent = -1,
+ ErrorPNotPresent = -2,
+ ErrorNrowNegative = -3,
+ ErrorNcolNegative = -4,
+ ErrorNnzNegative = -5,
+ ErrorP0Nonzero = -6,
+ ErrorATooSmall = -7,
+ ErrorColLengthNegative = -8,
+ ErrorRowIndexOutOfBounds = -9,
+ ErrorOutOfMemory = -10,
+ ErrorInternalError = -999
+};
/* ========================================================================== */
/* === Definitions ========================================================== */
/* ========================================================================== */
-#define ONES_COMPLEMENT(r) (-(r)-1)
+template <typename IndexType>
+IndexType ones_complement(const IndexType r) {
+ return (-(r)-1);
+}
/* -------------------------------------------------------------------------- */
-
-#define COLAMD_EMPTY (-1)
+const int Empty = -1;
/* Row and column status */
-#define ALIVE (0)
-#define DEAD (-1)
+enum RowColumnStatus {
+ Alive = 0,
+ Dead = -1
+};
/* Column status */
-#define DEAD_PRINCIPAL (-1)
-#define DEAD_NON_PRINCIPAL (-2)
-
-/* Macros for row and column status update and checking. */
-#define ROW_IS_DEAD(r) ROW_IS_MARKED_DEAD (Row[r].shared2.mark)
-#define ROW_IS_MARKED_DEAD(row_mark) (row_mark < ALIVE)
-#define ROW_IS_ALIVE(r) (Row [r].shared2.mark >= ALIVE)
-#define COL_IS_DEAD(c) (Col [c].start < ALIVE)
-#define COL_IS_ALIVE(c) (Col [c].start >= ALIVE)
-#define COL_IS_DEAD_PRINCIPAL(c) (Col [c].start == DEAD_PRINCIPAL)
-#define KILL_ROW(r) { Row [r].shared2.mark = DEAD ; }
-#define KILL_PRINCIPAL_COL(c) { Col [c].start = DEAD_PRINCIPAL ; }
-#define KILL_NON_PRINCIPAL_COL(c) { Col [c].start = DEAD_NON_PRINCIPAL ; }
+enum ColumnStatus {
+ DeadPrincipal = -1,
+ DeadNonPrincipal = -2
+};
/* ========================================================================== */
/* === Colamd reporting mechanism =========================================== */
@@ -129,9 +133,9 @@ namespace internal {
// == Row and Column structures ==
template <typename IndexType>
-struct colamd_col
+struct ColStructure
{
- IndexType start ; /* index for A of first row in this column, or DEAD */
+ IndexType start ; /* index for A of first row in this column, or Dead */
/* if column is dead */
IndexType length ; /* number of rows in this column */
union
@@ -159,11 +163,21 @@ struct colamd_col
IndexType degree_next ; /* next column, if col is in a degree list */
IndexType hash_next ; /* next column, if col is in a hash list */
} shared4 ;
-
+
+ inline bool is_dead() const { return start < Alive; }
+
+ inline bool is_alive() const { return start >= Alive; }
+
+ inline bool is_dead_principal() const { return start == DeadPrincipal; }
+
+ inline void kill_principal() { start = DeadPrincipal; }
+
+ inline void kill_non_principal() { start = DeadNonPrincipal; }
+
};
-
+
template <typename IndexType>
-struct Colamd_Row
+struct RowStructure
{
IndexType start ; /* index for A of first col in this row */
IndexType length ; /* number of principal columns in this row */
@@ -177,13 +191,19 @@ struct Colamd_Row
IndexType mark ; /* for computing set differences and marking dead rows*/
IndexType first_column ;/* first column in row (used in garbage collection) */
} shared2 ;
-
+
+ inline bool is_dead() const { return shared2.mark < Alive; }
+
+ inline bool is_alive() const { return shared2.mark >= Alive; }
+
+ inline void kill() { shared2.mark = Dead; }
+
};
-
+
/* ========================================================================== */
/* === Colamd recommended memory size ======================================= */
/* ========================================================================== */
-
+
/*
The recommended length Alen of the array A passed to colamd is given by
the COLAMD_RECOMMENDED (nnz, n_row, n_col) macro. It returns -1 if any
@@ -192,41 +212,41 @@ struct Colamd_Row
required for the Col and Row arrays, respectively, which are internal to
colamd. An additional n_col space is the minimal amount of "elbow room",
and nnz/5 more space is recommended for run time efficiency.
-
+
This macro is not needed when using symamd.
-
+
Explicit typecast to IndexType added Sept. 23, 2002, COLAMD version 2.2, to avoid
gcc -pedantic warning messages.
*/
template <typename IndexType>
-inline IndexType colamd_c(IndexType n_col)
-{ return IndexType( ((n_col) + 1) * sizeof (colamd_col<IndexType>) / sizeof (IndexType) ) ; }
+inline IndexType colamd_c(IndexType n_col)
+{ return IndexType( ((n_col) + 1) * sizeof (ColStructure<IndexType>) / sizeof (IndexType) ) ; }
template <typename IndexType>
inline IndexType colamd_r(IndexType n_row)
-{ return IndexType(((n_row) + 1) * sizeof (Colamd_Row<IndexType>) / sizeof (IndexType)); }
+{ return IndexType(((n_row) + 1) * sizeof (RowStructure<IndexType>) / sizeof (IndexType)); }
// Prototypes of non-user callable routines
template <typename IndexType>
-static IndexType init_rows_cols (IndexType n_row, IndexType n_col, Colamd_Row<IndexType> Row [], colamd_col<IndexType> col [], IndexType A [], IndexType p [], IndexType stats[COLAMD_STATS] );
+static IndexType init_rows_cols (IndexType n_row, IndexType n_col, RowStructure<IndexType> Row [], ColStructure<IndexType> col [], IndexType A [], IndexType p [], IndexType stats[NStats] );
template <typename IndexType>
-static void init_scoring (IndexType n_row, IndexType n_col, Colamd_Row<IndexType> Row [], colamd_col<IndexType> Col [], IndexType A [], IndexType head [], double knobs[COLAMD_KNOBS], IndexType *p_n_row2, IndexType *p_n_col2, IndexType *p_max_deg);
+static void init_scoring (IndexType n_row, IndexType n_col, RowStructure<IndexType> Row [], ColStructure<IndexType> Col [], IndexType A [], IndexType head [], double knobs[NKnobs], IndexType *p_n_row2, IndexType *p_n_col2, IndexType *p_max_deg);
template <typename IndexType>
-static IndexType find_ordering (IndexType n_row, IndexType n_col, IndexType Alen, Colamd_Row<IndexType> Row [], colamd_col<IndexType> Col [], IndexType A [], IndexType head [], IndexType n_col2, IndexType max_deg, IndexType pfree);
+static IndexType find_ordering (IndexType n_row, IndexType n_col, IndexType Alen, RowStructure<IndexType> Row [], ColStructure<IndexType> Col [], IndexType A [], IndexType head [], IndexType n_col2, IndexType max_deg, IndexType pfree);
template <typename IndexType>
-static void order_children (IndexType n_col, colamd_col<IndexType> Col [], IndexType p []);
+static void order_children (IndexType n_col, ColStructure<IndexType> Col [], IndexType p []);
template <typename IndexType>
-static void detect_super_cols (colamd_col<IndexType> Col [], IndexType A [], IndexType head [], IndexType row_start, IndexType row_length ) ;
+static void detect_super_cols (ColStructure<IndexType> Col [], IndexType A [], IndexType head [], IndexType row_start, IndexType row_length ) ;
template <typename IndexType>
-static IndexType garbage_collection (IndexType n_row, IndexType n_col, Colamd_Row<IndexType> Row [], colamd_col<IndexType> Col [], IndexType A [], IndexType *pfree) ;
+static IndexType garbage_collection (IndexType n_row, IndexType n_col, RowStructure<IndexType> Row [], ColStructure<IndexType> Col [], IndexType A [], IndexType *pfree) ;
template <typename IndexType>
-static inline IndexType clear_mark (IndexType n_row, Colamd_Row<IndexType> Row [] ) ;
+static inline IndexType clear_mark (IndexType n_row, RowStructure<IndexType> Row [] ) ;
/* === No debugging ========================================================= */
@@ -240,37 +260,37 @@ static inline IndexType clear_mark (IndexType n_row, Colamd_Row<IndexType> Row
/**
- * \brief Returns the recommended value of Alen
- *
- * Returns recommended value of Alen for use by colamd.
- * Returns -1 if any input argument is negative.
- * The use of this routine or macro is optional.
- * Note that the macro uses its arguments more than once,
- * so be careful for side effects, if you pass expressions as arguments to COLAMD_RECOMMENDED.
- *
+ * \brief Returns the recommended value of Alen
+ *
+ * Returns recommended value of Alen for use by colamd.
+ * Returns -1 if any input argument is negative.
+ * The use of this routine or macro is optional.
+ * Note that the macro uses its arguments more than once,
+ * so be careful for side effects, if you pass expressions as arguments to COLAMD_RECOMMENDED.
+ *
* \param nnz nonzeros in A
* \param n_row number of rows in A
* \param n_col number of columns in A
* \return recommended value of Alen for use by colamd
*/
template <typename IndexType>
-inline IndexType colamd_recommended ( IndexType nnz, IndexType n_row, IndexType n_col)
+inline IndexType recommended ( IndexType nnz, IndexType n_row, IndexType n_col)
{
if ((nnz) < 0 || (n_row) < 0 || (n_col) < 0)
return (-1);
else
- return (2 * (nnz) + colamd_c (n_col) + colamd_r (n_row) + (n_col) + ((nnz) / 5));
+ return (2 * (nnz) + colamd_c (n_col) + colamd_r (n_row) + (n_col) + ((nnz) / 5));
}
/**
* \brief set default parameters The use of this routine is optional.
- *
- * Colamd: rows with more than (knobs [COLAMD_DENSE_ROW] * n_col)
+ *
+ * Colamd: rows with more than (knobs [DenseRow] * n_col)
* entries are removed prior to ordering. Columns with more than
- * (knobs [COLAMD_DENSE_COL] * n_row) entries are removed prior to
- * ordering, and placed last in the output column ordering.
+ * (knobs [DenseCol] * n_row) entries are removed prior to
+ * ordering, and placed last in the output column ordering.
*
- * COLAMD_DENSE_ROW and COLAMD_DENSE_COL are defined as 0 and 1,
+ * DenseRow and DenseCol are defined as 0 and 1,
* respectively, in colamd.h. Default values of these two knobs
* are both 0.5. Currently, only knobs [0] and knobs [1] are
* used, but future versions may use more knobs. If so, they will
@@ -279,37 +299,37 @@ inline IndexType colamd_recommended ( IndexType nnz, IndexType n_row, IndexType
* not need to change, assuming that you either use
* colamd_set_defaults, or pass a (double *) NULL pointer as the
* knobs array to colamd or symamd.
- *
+ *
* \param knobs parameter settings for colamd
*/
-static inline void colamd_set_defaults(double knobs[COLAMD_KNOBS])
+static inline void set_defaults(double knobs[NKnobs])
{
/* === Local variables ================================================== */
-
+
int i ;
if (!knobs)
{
return ; /* no knobs to initialize */
}
- for (i = 0 ; i < COLAMD_KNOBS ; i++)
+ for (i = 0 ; i < NKnobs ; i++)
{
knobs [i] = 0 ;
}
- knobs [COLAMD_DENSE_ROW] = 0.5 ; /* ignore rows over 50% dense */
- knobs [COLAMD_DENSE_COL] = 0.5 ; /* ignore columns over 50% dense */
+ knobs [Colamd::DenseRow] = 0.5 ; /* ignore rows over 50% dense */
+ knobs [Colamd::DenseCol] = 0.5 ; /* ignore columns over 50% dense */
}
-/**
+/**
* \brief Computes a column ordering using the column approximate minimum degree ordering
- *
+ *
* Computes a column ordering (Q) of A such that P(AQ)=LU or
* (AQ)'AQ=LL' have less fill-in and require fewer floating point
* operations than factorizing the unpermuted matrix A or A'A,
* respectively.
- *
- *
+ *
+ *
* \param n_row number of rows in A
* \param n_col number of columns in A
* \param Alen, size of the array A
@@ -319,143 +339,143 @@ static inline void colamd_set_defaults(double knobs[COLAMD_KNOBS])
* \param stats colamd output statistics and error codes
*/
template <typename IndexType>
-static bool colamd(IndexType n_row, IndexType n_col, IndexType Alen, IndexType *A, IndexType *p, double knobs[COLAMD_KNOBS], IndexType stats[COLAMD_STATS])
+static bool compute_ordering(IndexType n_row, IndexType n_col, IndexType Alen, IndexType *A, IndexType *p, double knobs[NKnobs], IndexType stats[NStats])
{
/* === Local variables ================================================== */
-
+
IndexType i ; /* loop index */
IndexType nnz ; /* nonzeros in A */
IndexType Row_size ; /* size of Row [], in integers */
IndexType Col_size ; /* size of Col [], in integers */
IndexType need ; /* minimum required length of A */
- Colamd_Row<IndexType> *Row ; /* pointer into A of Row [0..n_row] array */
- colamd_col<IndexType> *Col ; /* pointer into A of Col [0..n_col] array */
+ Colamd::RowStructure<IndexType> *Row ; /* pointer into A of Row [0..n_row] array */
+ Colamd::ColStructure<IndexType> *Col ; /* pointer into A of Col [0..n_col] array */
IndexType n_col2 ; /* number of non-dense, non-empty columns */
IndexType n_row2 ; /* number of non-dense, non-empty rows */
IndexType ngarbage ; /* number of garbage collections performed */
IndexType max_deg ; /* maximum row degree */
- double default_knobs [COLAMD_KNOBS] ; /* default knobs array */
-
-
+ double default_knobs [NKnobs] ; /* default knobs array */
+
+
/* === Check the input arguments ======================================== */
-
+
if (!stats)
{
COLAMD_DEBUG0 (("colamd: stats not present\n")) ;
return (false) ;
}
- for (i = 0 ; i < COLAMD_STATS ; i++)
+ for (i = 0 ; i < NStats ; i++)
{
stats [i] = 0 ;
}
- stats [COLAMD_STATUS] = COLAMD_OK ;
- stats [COLAMD_INFO1] = -1 ;
- stats [COLAMD_INFO2] = -1 ;
-
+ stats [Colamd::Status] = Colamd::Ok ;
+ stats [Colamd::Info1] = -1 ;
+ stats [Colamd::Info2] = -1 ;
+
if (!A) /* A is not present */
{
- stats [COLAMD_STATUS] = COLAMD_ERROR_A_not_present ;
+ stats [Colamd::Status] = Colamd::ErrorANotPresent ;
COLAMD_DEBUG0 (("colamd: A not present\n")) ;
return (false) ;
}
-
+
if (!p) /* p is not present */
{
- stats [COLAMD_STATUS] = COLAMD_ERROR_p_not_present ;
+ stats [Colamd::Status] = Colamd::ErrorPNotPresent ;
COLAMD_DEBUG0 (("colamd: p not present\n")) ;
return (false) ;
}
-
+
if (n_row < 0) /* n_row must be >= 0 */
{
- stats [COLAMD_STATUS] = COLAMD_ERROR_nrow_negative ;
- stats [COLAMD_INFO1] = n_row ;
+ stats [Colamd::Status] = Colamd::ErrorNrowNegative ;
+ stats [Colamd::Info1] = n_row ;
COLAMD_DEBUG0 (("colamd: nrow negative %d\n", n_row)) ;
return (false) ;
}
-
+
if (n_col < 0) /* n_col must be >= 0 */
{
- stats [COLAMD_STATUS] = COLAMD_ERROR_ncol_negative ;
- stats [COLAMD_INFO1] = n_col ;
+ stats [Colamd::Status] = Colamd::ErrorNcolNegative ;
+ stats [Colamd::Info1] = n_col ;
COLAMD_DEBUG0 (("colamd: ncol negative %d\n", n_col)) ;
return (false) ;
}
-
+
nnz = p [n_col] ;
if (nnz < 0) /* nnz must be >= 0 */
{
- stats [COLAMD_STATUS] = COLAMD_ERROR_nnz_negative ;
- stats [COLAMD_INFO1] = nnz ;
+ stats [Colamd::Status] = Colamd::ErrorNnzNegative ;
+ stats [Colamd::Info1] = nnz ;
COLAMD_DEBUG0 (("colamd: number of entries negative %d\n", nnz)) ;
return (false) ;
}
-
+
if (p [0] != 0)
{
- stats [COLAMD_STATUS] = COLAMD_ERROR_p0_nonzero ;
- stats [COLAMD_INFO1] = p [0] ;
+ stats [Colamd::Status] = Colamd::ErrorP0Nonzero ;
+ stats [Colamd::Info1] = p [0] ;
COLAMD_DEBUG0 (("colamd: p[0] not zero %d\n", p [0])) ;
return (false) ;
}
-
+
/* === If no knobs, set default knobs =================================== */
-
+
if (!knobs)
{
- colamd_set_defaults (default_knobs) ;
+ set_defaults (default_knobs) ;
knobs = default_knobs ;
}
-
+
/* === Allocate the Row and Col arrays from array A ===================== */
-
+
Col_size = colamd_c (n_col) ;
Row_size = colamd_r (n_row) ;
need = 2*nnz + n_col + Col_size + Row_size ;
-
+
if (need > Alen)
{
/* not enough space in array A to perform the ordering */
- stats [COLAMD_STATUS] = COLAMD_ERROR_A_too_small ;
- stats [COLAMD_INFO1] = need ;
- stats [COLAMD_INFO2] = Alen ;
+ stats [Colamd::Status] = Colamd::ErrorATooSmall ;
+ stats [Colamd::Info1] = need ;
+ stats [Colamd::Info2] = Alen ;
COLAMD_DEBUG0 (("colamd: Need Alen >= %d, given only Alen = %d\n", need,Alen));
return (false) ;
}
-
+
Alen -= Col_size + Row_size ;
- Col = (colamd_col<IndexType> *) &A [Alen] ;
- Row = (Colamd_Row<IndexType> *) &A [Alen + Col_size] ;
+ Col = (ColStructure<IndexType> *) &A [Alen] ;
+ Row = (RowStructure<IndexType> *) &A [Alen + Col_size] ;
/* === Construct the row and column data structures ===================== */
-
- if (!Eigen::internal::init_rows_cols (n_row, n_col, Row, Col, A, p, stats))
+
+ if (!Colamd::init_rows_cols (n_row, n_col, Row, Col, A, p, stats))
{
/* input matrix is invalid */
COLAMD_DEBUG0 (("colamd: Matrix invalid\n")) ;
return (false) ;
}
-
+
/* === Initialize scores, kill dense rows/columns ======================= */
- Eigen::internal::init_scoring (n_row, n_col, Row, Col, A, p, knobs,
+ Colamd::init_scoring (n_row, n_col, Row, Col, A, p, knobs,
&n_row2, &n_col2, &max_deg) ;
-
+
/* === Order the supercolumns =========================================== */
-
- ngarbage = Eigen::internal::find_ordering (n_row, n_col, Alen, Row, Col, A, p,
+
+ ngarbage = Colamd::find_ordering (n_row, n_col, Alen, Row, Col, A, p,
n_col2, max_deg, 2*nnz) ;
-
+
/* === Order the non-principal columns ================================== */
-
- Eigen::internal::order_children (n_col, Col, p) ;
-
+
+ Colamd::order_children (n_col, Col, p) ;
+
/* === Return statistics in stats ======================================= */
-
- stats [COLAMD_DENSE_ROW] = n_row - n_row2 ;
- stats [COLAMD_DENSE_COL] = n_col - n_col2 ;
- stats [COLAMD_DEFRAG_COUNT] = ngarbage ;
- COLAMD_DEBUG0 (("colamd: done.\n")) ;
+
+ stats [Colamd::DenseRow] = n_row - n_row2 ;
+ stats [Colamd::DenseCol] = n_col - n_col2 ;
+ stats [Colamd::DefragCount] = ngarbage ;
+ COLAMD_DEBUG0 (("colamd: done.\n")) ;
return (true) ;
}
@@ -465,7 +485,6 @@ static bool colamd(IndexType n_row, IndexType n_col, IndexType Alen, IndexType *
/* There are no user-callable routines beyond this point in the file */
-
/* ========================================================================== */
/* === init_rows_cols ======================================================= */
/* ========================================================================== */
@@ -485,11 +504,11 @@ static IndexType init_rows_cols /* returns true if OK, or false otherwise */
IndexType n_row, /* number of rows of A */
IndexType n_col, /* number of columns of A */
- Colamd_Row<IndexType> Row [], /* of size n_row+1 */
- colamd_col<IndexType> Col [], /* of size n_col+1 */
+ RowStructure<IndexType> Row [], /* of size n_row+1 */
+ ColStructure<IndexType> Col [], /* of size n_col+1 */
IndexType A [], /* row indices of A, of size Alen */
IndexType p [], /* pointers to columns in A, of size n_col+1 */
- IndexType stats [COLAMD_STATS] /* colamd statistics */
+ IndexType stats [NStats] /* colamd statistics */
)
{
/* === Local variables ================================================== */
@@ -512,24 +531,24 @@ static IndexType init_rows_cols /* returns true if OK, or false otherwise */
if ((Col [col].length) < 0) // extra parentheses to work-around gcc bug 10200
{
/* column pointers must be non-decreasing */
- stats [COLAMD_STATUS] = COLAMD_ERROR_col_length_negative ;
- stats [COLAMD_INFO1] = col ;
- stats [COLAMD_INFO2] = Col [col].length ;
+ stats [Colamd::Status] = Colamd::ErrorColLengthNegative ;
+ stats [Colamd::Info1] = col ;
+ stats [Colamd::Info2] = Col [col].length ;
COLAMD_DEBUG0 (("colamd: col %d length %d < 0\n", col, Col [col].length)) ;
return (false) ;
}
Col [col].shared1.thickness = 1 ;
Col [col].shared2.score = 0 ;
- Col [col].shared3.prev = COLAMD_EMPTY ;
- Col [col].shared4.degree_next = COLAMD_EMPTY ;
+ Col [col].shared3.prev = Empty ;
+ Col [col].shared4.degree_next = Empty ;
}
/* p [0..n_col] no longer needed, used as "head" in subsequent routines */
/* === Scan columns, compute row degrees, and check row indices ========= */
- stats [COLAMD_INFO3] = 0 ; /* number of duplicate or unsorted row indices*/
+ stats [Info3] = 0 ; /* number of duplicate or unsorted row indices*/
for (row = 0 ; row < n_row ; row++)
{
@@ -551,10 +570,10 @@ static IndexType init_rows_cols /* returns true if OK, or false otherwise */
/* make sure row indices within range */
if (row < 0 || row >= n_row)
{
- stats [COLAMD_STATUS] = COLAMD_ERROR_row_index_out_of_bounds ;
- stats [COLAMD_INFO1] = col ;
- stats [COLAMD_INFO2] = row ;
- stats [COLAMD_INFO3] = n_row ;
+ stats [Colamd::Status] = Colamd::ErrorRowIndexOutOfBounds ;
+ stats [Colamd::Info1] = col ;
+ stats [Colamd::Info2] = row ;
+ stats [Colamd::Info3] = n_row ;
COLAMD_DEBUG0 (("colamd: row %d col %d out of bounds\n", row, col)) ;
return (false) ;
}
@@ -563,10 +582,10 @@ static IndexType init_rows_cols /* returns true if OK, or false otherwise */
{
/* row index are unsorted or repeated (or both), thus col */
/* is jumbled. This is a notice, not an error condition. */
- stats [COLAMD_STATUS] = COLAMD_OK_BUT_JUMBLED ;
- stats [COLAMD_INFO1] = col ;
- stats [COLAMD_INFO2] = row ;
- (stats [COLAMD_INFO3]) ++ ;
+ stats [Colamd::Status] = Colamd::OkButJumbled ;
+ stats [Colamd::Info1] = col ;
+ stats [Colamd::Info2] = row ;
+ (stats [Colamd::Info3]) ++ ;
COLAMD_DEBUG1 (("colamd: row %d col %d unsorted/duplicate\n",row,col));
}
@@ -604,7 +623,7 @@ static IndexType init_rows_cols /* returns true if OK, or false otherwise */
/* === Create row form ================================================== */
- if (stats [COLAMD_STATUS] == COLAMD_OK_BUT_JUMBLED)
+ if (stats [Status] == OkButJumbled)
{
/* if cols jumbled, watch for repeated row indices */
for (col = 0 ; col < n_col ; col++)
@@ -646,7 +665,7 @@ static IndexType init_rows_cols /* returns true if OK, or false otherwise */
/* === See if we need to re-create columns ============================== */
- if (stats [COLAMD_STATUS] == COLAMD_OK_BUT_JUMBLED)
+ if (stats [Status] == OkButJumbled)
{
COLAMD_DEBUG0 (("colamd: reconstructing column form, matrix jumbled\n")) ;
@@ -701,11 +720,11 @@ static void init_scoring
IndexType n_row, /* number of rows of A */
IndexType n_col, /* number of columns of A */
- Colamd_Row<IndexType> Row [], /* of size n_row+1 */
- colamd_col<IndexType> Col [], /* of size n_col+1 */
+ RowStructure<IndexType> Row [], /* of size n_row+1 */
+ ColStructure<IndexType> Col [], /* of size n_col+1 */
IndexType A [], /* column form and row form of A */
IndexType head [], /* of size n_col+1 */
- double knobs [COLAMD_KNOBS],/* parameters */
+ double knobs [NKnobs],/* parameters */
IndexType *p_n_row2, /* number of non-dense, non-empty rows */
IndexType *p_n_col2, /* number of non-dense, non-empty columns */
IndexType *p_max_deg /* maximum row degree */
@@ -732,8 +751,8 @@ static void init_scoring
/* === Extract knobs ==================================================== */
- dense_row_count = numext::maxi(IndexType(0), numext::mini(IndexType(knobs [COLAMD_DENSE_ROW] * n_col), n_col)) ;
- dense_col_count = numext::maxi(IndexType(0), numext::mini(IndexType(knobs [COLAMD_DENSE_COL] * n_row), n_row)) ;
+ dense_row_count = numext::maxi(IndexType(0), numext::mini(IndexType(knobs [Colamd::DenseRow] * n_col), n_col)) ;
+ dense_col_count = numext::maxi(IndexType(0), numext::mini(IndexType(knobs [Colamd::DenseCol] * n_row), n_row)) ;
COLAMD_DEBUG1 (("colamd: densecount: %d %d\n", dense_row_count, dense_col_count)) ;
max_deg = 0 ;
n_col2 = n_col ;
@@ -750,7 +769,7 @@ static void init_scoring
{
/* this is a empty column, kill and order it last */
Col [c].shared2.order = --n_col2 ;
- KILL_PRINCIPAL_COL (c) ;
+ Col[c].kill_principal() ;
}
}
COLAMD_DEBUG1 (("colamd: null columns killed: %d\n", n_col - n_col2)) ;
@@ -761,7 +780,7 @@ static void init_scoring
for (c = n_col-1 ; c >= 0 ; c--)
{
/* skip any dead columns */
- if (COL_IS_DEAD (c))
+ if (Col[c].is_dead())
{
continue ;
}
@@ -777,7 +796,7 @@ static void init_scoring
{
Row [*cp++].shared1.degree-- ;
}
- KILL_PRINCIPAL_COL (c) ;
+ Col[c].kill_principal() ;
}
}
COLAMD_DEBUG1 (("colamd: Dense and null columns killed: %d\n", n_col - n_col2)) ;
@@ -791,7 +810,7 @@ static void init_scoring
if (deg > dense_row_count || deg == 0)
{
/* kill a dense or empty row */
- KILL_ROW (r) ;
+ Row[r].kill() ;
--n_row2 ;
}
else
@@ -813,7 +832,7 @@ static void init_scoring
for (c = n_col-1 ; c >= 0 ; c--)
{
/* skip dead column */
- if (COL_IS_DEAD (c))
+ if (Col[c].is_dead())
{
continue ;
}
@@ -826,7 +845,7 @@ static void init_scoring
/* get a row */
row = *cp++ ;
/* skip if dead */
- if (ROW_IS_DEAD (row))
+ if (Row[row].is_dead())
{
continue ;
}
@@ -845,7 +864,7 @@ static void init_scoring
/* and have already been killed) */
COLAMD_DEBUG2 (("Newly null killed: %d\n", c)) ;
Col [c].shared2.order = --n_col2 ;
- KILL_PRINCIPAL_COL (c) ;
+ Col[c].kill_principal() ;
}
else
{
@@ -870,7 +889,7 @@ static void init_scoring
/* clear the hash buckets */
for (c = 0 ; c <= n_col ; c++)
{
- head [c] = COLAMD_EMPTY ;
+ head [c] = Empty ;
}
min_score = n_col ;
/* place in reverse order, so low column indices are at the front */
@@ -878,7 +897,7 @@ static void init_scoring
for (c = n_col-1 ; c >= 0 ; c--)
{
/* only add principal columns to degree lists */
- if (COL_IS_ALIVE (c))
+ if (Col[c].is_alive())
{
COLAMD_DEBUG4 (("place %d score %d minscore %d ncol %d\n",
c, Col [c].shared2.score, min_score, n_col)) ;
@@ -891,16 +910,16 @@ static void init_scoring
COLAMD_ASSERT (min_score <= n_col) ;
COLAMD_ASSERT (score >= 0) ;
COLAMD_ASSERT (score <= n_col) ;
- COLAMD_ASSERT (head [score] >= COLAMD_EMPTY) ;
+ COLAMD_ASSERT (head [score] >= Empty) ;
/* now add this column to dList at proper score location */
next_col = head [score] ;
- Col [c].shared3.prev = COLAMD_EMPTY ;
+ Col [c].shared3.prev = Empty ;
Col [c].shared4.degree_next = next_col ;
/* if there already was a column with the same score, set its */
/* previous pointer to this new column */
- if (next_col != COLAMD_EMPTY)
+ if (next_col != Empty)
{
Col [next_col].shared3.prev = c ;
}
@@ -939,8 +958,8 @@ static IndexType find_ordering /* return the number of garbage collections */
IndexType n_row, /* number of rows of A */
IndexType n_col, /* number of columns of A */
IndexType Alen, /* size of A, 2*nnz + n_col or larger */
- Colamd_Row<IndexType> Row [], /* of size n_row+1 */
- colamd_col<IndexType> Col [], /* of size n_col+1 */
+ RowStructure<IndexType> Row [], /* of size n_row+1 */
+ ColStructure<IndexType> Col [], /* of size n_col+1 */
IndexType A [], /* column form and row form of A */
IndexType head [], /* of size n_col+1 */
IndexType n_col2, /* Remaining columns to order */
@@ -986,7 +1005,7 @@ static IndexType find_ordering /* return the number of garbage collections */
/* === Initialization and clear mark ==================================== */
max_mark = INT_MAX - n_col ; /* INT_MAX defined in <limits.h> */
- tag_mark = Eigen::internal::clear_mark (n_row, Row) ;
+ tag_mark = Colamd::clear_mark (n_row, Row) ;
min_score = 0 ;
ngarbage = 0 ;
COLAMD_DEBUG1 (("colamd: Ordering, n_col2=%d\n", n_col2)) ;
@@ -1001,10 +1020,10 @@ static IndexType find_ordering /* return the number of garbage collections */
/* make sure degree list isn't empty */
COLAMD_ASSERT (min_score >= 0) ;
COLAMD_ASSERT (min_score <= n_col) ;
- COLAMD_ASSERT (head [min_score] >= COLAMD_EMPTY) ;
+ COLAMD_ASSERT (head [min_score] >= Empty) ;
/* get pivot column from head of minimum degree list */
- while (min_score < n_col && head [min_score] == COLAMD_EMPTY)
+ while (min_score < n_col && head [min_score] == Empty)
{
min_score++ ;
}
@@ -1012,12 +1031,12 @@ static IndexType find_ordering /* return the number of garbage collections */
COLAMD_ASSERT (pivot_col >= 0 && pivot_col <= n_col) ;
next_col = Col [pivot_col].shared4.degree_next ;
head [min_score] = next_col ;
- if (next_col != COLAMD_EMPTY)
+ if (next_col != Empty)
{
- Col [next_col].shared3.prev = COLAMD_EMPTY ;
+ Col [next_col].shared3.prev = Empty ;
}
- COLAMD_ASSERT (COL_IS_ALIVE (pivot_col)) ;
+ COLAMD_ASSERT (Col[pivot_col].is_alive()) ;
COLAMD_DEBUG3 (("Pivot col: %d\n", pivot_col)) ;
/* remember score for defrag check */
@@ -1036,12 +1055,12 @@ static IndexType find_ordering /* return the number of garbage collections */
needed_memory = numext::mini(pivot_col_score, n_col - k) ;
if (pfree + needed_memory >= Alen)
{
- pfree = Eigen::internal::garbage_collection (n_row, n_col, Row, Col, A, &A [pfree]) ;
+ pfree = Colamd::garbage_collection (n_row, n_col, Row, Col, A, &A [pfree]) ;
ngarbage++ ;
/* after garbage collection we will have enough */
COLAMD_ASSERT (pfree + needed_memory < Alen) ;
/* garbage collection has wiped out the Row[].shared2.mark array */
- tag_mark = Eigen::internal::clear_mark (n_row, Row) ;
+ tag_mark = Colamd::clear_mark (n_row, Row) ;
}
@@ -1064,9 +1083,9 @@ static IndexType find_ordering /* return the number of garbage collections */
{
/* get a row */
row = *cp++ ;
- COLAMD_DEBUG4 (("Pivot col pattern %d %d\n", ROW_IS_ALIVE (row), row)) ;
+ COLAMD_DEBUG4 (("Pivot col pattern %d %d\n", Row[row].is_alive(), row)) ;
/* skip if row is dead */
- if (ROW_IS_DEAD (row))
+ if (Row[row].is_dead())
{
continue ;
}
@@ -1078,7 +1097,7 @@ static IndexType find_ordering /* return the number of garbage collections */
col = *rp++ ;
/* add the column, if alive and untagged */
col_thickness = Col [col].shared1.thickness ;
- if (col_thickness > 0 && COL_IS_ALIVE (col))
+ if (col_thickness > 0 && Col[col].is_alive())
{
/* tag column in pivot row */
Col [col].shared1.thickness = -col_thickness ;
@@ -1105,7 +1124,7 @@ static IndexType find_ordering /* return the number of garbage collections */
/* may be killing an already dead row */
row = *cp++ ;
COLAMD_DEBUG3 (("Kill row in pivot col: %d\n", row)) ;
- KILL_ROW (row) ;
+ Row[row].kill() ;
}
/* === Select a row index to use as the new pivot row =============== */
@@ -1120,7 +1139,7 @@ static IndexType find_ordering /* return the number of garbage collections */
else
{
/* there is no pivot row, since it is of zero length */
- pivot_row = COLAMD_EMPTY ;
+ pivot_row = Empty ;
COLAMD_ASSERT (pivot_row_length == 0) ;
}
COLAMD_ASSERT (Col [pivot_col].length > 0 || pivot_row_length == 0) ;
@@ -1157,7 +1176,7 @@ static IndexType find_ordering /* return the number of garbage collections */
while (rp < rp_end)
{
col = *rp++ ;
- COLAMD_ASSERT (COL_IS_ALIVE (col) && col != pivot_col) ;
+ COLAMD_ASSERT (Col[col].is_alive() && col != pivot_col) ;
COLAMD_DEBUG3 (("Col: %d\n", col)) ;
/* clear tags used to construct pivot row pattern */
@@ -1172,8 +1191,8 @@ static IndexType find_ordering /* return the number of garbage collections */
next_col = Col [col].shared4.degree_next ;
COLAMD_ASSERT (cur_score >= 0) ;
COLAMD_ASSERT (cur_score <= n_col) ;
- COLAMD_ASSERT (cur_score >= COLAMD_EMPTY) ;
- if (prev_col == COLAMD_EMPTY)
+ COLAMD_ASSERT (cur_score >= Empty) ;
+ if (prev_col == Empty)
{
head [cur_score] = next_col ;
}
@@ -1181,7 +1200,7 @@ static IndexType find_ordering /* return the number of garbage collections */
{
Col [prev_col].shared4.degree_next = next_col ;
}
- if (next_col != COLAMD_EMPTY)
+ if (next_col != Empty)
{
Col [next_col].shared3.prev = prev_col ;
}
@@ -1194,12 +1213,12 @@ static IndexType find_ordering /* return the number of garbage collections */
{
/* get a row */
row = *cp++ ;
- row_mark = Row [row].shared2.mark ;
/* skip if dead */
- if (ROW_IS_MARKED_DEAD (row_mark))
+ if (Row[row].is_dead())
{
continue ;
}
+ row_mark = Row [row].shared2.mark ;
COLAMD_ASSERT (row != pivot_row) ;
set_difference = row_mark - tag_mark ;
/* check if the row has been seen yet */
@@ -1215,7 +1234,7 @@ static IndexType find_ordering /* return the number of garbage collections */
if (set_difference == 0)
{
COLAMD_DEBUG3 (("aggressive absorption. Row: %d\n", row)) ;
- KILL_ROW (row) ;
+ Row[row].kill() ;
}
else
{
@@ -1237,7 +1256,7 @@ static IndexType find_ordering /* return the number of garbage collections */
{
/* get a column */
col = *rp++ ;
- COLAMD_ASSERT (COL_IS_ALIVE (col) && col != pivot_col) ;
+ COLAMD_ASSERT (Col[col].is_alive() && col != pivot_col) ;
hash = 0 ;
cur_score = 0 ;
cp = &A [Col [col].start] ;
@@ -1252,12 +1271,12 @@ static IndexType find_ordering /* return the number of garbage collections */
/* get a row */
row = *cp++ ;
COLAMD_ASSERT(row >= 0 && row < n_row) ;
- row_mark = Row [row].shared2.mark ;
/* skip if dead */
- if (ROW_IS_MARKED_DEAD (row_mark))
+ if (Row [row].is_dead())
{
continue ;
}
+ row_mark = Row [row].shared2.mark ;
COLAMD_ASSERT (row_mark > tag_mark) ;
/* compact the column */
*new_cp++ = row ;
@@ -1278,7 +1297,7 @@ static IndexType find_ordering /* return the number of garbage collections */
{
COLAMD_DEBUG4 (("further mass elimination. Col: %d\n", col)) ;
/* nothing left but the pivot row in this column */
- KILL_PRINCIPAL_COL (col) ;
+ Col[col].kill_principal() ;
pivot_row_degree -= Col [col].shared1.thickness ;
COLAMD_ASSERT (pivot_row_degree >= 0) ;
/* order it */
@@ -1302,7 +1321,7 @@ static IndexType find_ordering /* return the number of garbage collections */
COLAMD_ASSERT (hash <= n_col) ;
head_column = head [hash] ;
- if (head_column > COLAMD_EMPTY)
+ if (head_column > Empty)
{
/* degree list "hash" is non-empty, use prev (shared3) of */
/* first column in degree list as head of hash bucket */
@@ -1319,7 +1338,7 @@ static IndexType find_ordering /* return the number of garbage collections */
/* save hash function in Col [col].shared3.hash */
Col [col].shared3.hash = (IndexType) hash ;
- COLAMD_ASSERT (COL_IS_ALIVE (col)) ;
+ COLAMD_ASSERT (Col[col].is_alive()) ;
}
}
@@ -1329,11 +1348,11 @@ static IndexType find_ordering /* return the number of garbage collections */
COLAMD_DEBUG3 (("** Supercolumn detection phase. **\n")) ;
- Eigen::internal::detect_super_cols (Col, A, head, pivot_row_start, pivot_row_length) ;
+ Colamd::detect_super_cols (Col, A, head, pivot_row_start, pivot_row_length) ;
/* === Kill the pivotal column ====================================== */
- KILL_PRINCIPAL_COL (pivot_col) ;
+ Col[pivot_col].kill_principal() ;
/* === Clear mark =================================================== */
@@ -1341,7 +1360,7 @@ static IndexType find_ordering /* return the number of garbage collections */
if (tag_mark >= max_mark)
{
COLAMD_DEBUG2 (("clearing tag_mark\n")) ;
- tag_mark = Eigen::internal::clear_mark (n_row, Row) ;
+ tag_mark = Colamd::clear_mark (n_row, Row) ;
}
/* === Finalize the new pivot row, and column scores ================ */
@@ -1357,7 +1376,7 @@ static IndexType find_ordering /* return the number of garbage collections */
{
col = *rp++ ;
/* skip dead columns */
- if (COL_IS_DEAD (col))
+ if (Col[col].is_dead())
{
continue ;
}
@@ -1391,11 +1410,11 @@ static IndexType find_ordering /* return the number of garbage collections */
COLAMD_ASSERT (min_score <= n_col) ;
COLAMD_ASSERT (cur_score >= 0) ;
COLAMD_ASSERT (cur_score <= n_col) ;
- COLAMD_ASSERT (head [cur_score] >= COLAMD_EMPTY) ;
+ COLAMD_ASSERT (head [cur_score] >= Empty) ;
next_col = head [cur_score] ;
Col [col].shared4.degree_next = next_col ;
- Col [col].shared3.prev = COLAMD_EMPTY ;
- if (next_col != COLAMD_EMPTY)
+ Col [col].shared3.prev = Empty ;
+ if (next_col != Empty)
{
Col [next_col].shared3.prev = col ;
}
@@ -1448,7 +1467,7 @@ static inline void order_children
/* === Parameters ======================================================= */
IndexType n_col, /* number of columns of A */
- colamd_col<IndexType> Col [], /* of size n_col+1 */
+ ColStructure<IndexType> Col [], /* of size n_col+1 */
IndexType p [] /* p [0 ... n_col-1] is the column permutation*/
)
{
@@ -1464,15 +1483,15 @@ static inline void order_children
for (i = 0 ; i < n_col ; i++)
{
/* find an un-ordered non-principal column */
- COLAMD_ASSERT (COL_IS_DEAD (i)) ;
- if (!COL_IS_DEAD_PRINCIPAL (i) && Col [i].shared2.order == COLAMD_EMPTY)
+ COLAMD_ASSERT (col_is_dead(Col, i)) ;
+ if (!Col[i].is_dead_principal() && Col [i].shared2.order == Empty)
{
parent = i ;
/* once found, find its principal parent */
do
{
parent = Col [parent].shared1.parent ;
- } while (!COL_IS_DEAD_PRINCIPAL (parent)) ;
+ } while (!Col[parent].is_dead_principal()) ;
/* now, order all un-ordered non-principal columns along path */
/* to this parent. collapse tree at the same time */
@@ -1482,7 +1501,7 @@ static inline void order_children
do
{
- COLAMD_ASSERT (Col [c].shared2.order == COLAMD_EMPTY) ;
+ COLAMD_ASSERT (Col [c].shared2.order == Empty) ;
/* order this column */
Col [c].shared2.order = order++ ;
@@ -1495,7 +1514,7 @@ static inline void order_children
/* continue until we hit an ordered column. There are */
/* guaranteed not to be anymore unordered columns */
/* above an ordered column */
- } while (Col [c].shared2.order == COLAMD_EMPTY) ;
+ } while (Col [c].shared2.order == Empty) ;
/* re-order the super_col parent to largest order for this group */
Col [parent].shared2.order = order ;
@@ -1547,8 +1566,8 @@ template <typename IndexType>
static void detect_super_cols
(
/* === Parameters ======================================================= */
-
- colamd_col<IndexType> Col [], /* of size n_col+1 */
+
+ ColStructure<IndexType> Col [], /* of size n_col+1 */
IndexType A [], /* row indices of A */
IndexType head [], /* head of degree lists and hash buckets */
IndexType row_start, /* pointer to set of columns to check */
@@ -1578,7 +1597,7 @@ static void detect_super_cols
while (rp < rp_end)
{
col = *rp++ ;
- if (COL_IS_DEAD (col))
+ if (Col[col].is_dead())
{
continue ;
}
@@ -1590,7 +1609,7 @@ static void detect_super_cols
/* === Get the first column in this hash bucket ===================== */
head_column = head [hash] ;
- if (head_column > COLAMD_EMPTY)
+ if (head_column > Empty)
{
first_col = Col [head_column].shared3.headhash ;
}
@@ -1601,10 +1620,10 @@ static void detect_super_cols
/* === Consider each column in the hash bucket ====================== */
- for (super_c = first_col ; super_c != COLAMD_EMPTY ;
+ for (super_c = first_col ; super_c != Empty ;
super_c = Col [super_c].shared4.hash_next)
{
- COLAMD_ASSERT (COL_IS_ALIVE (super_c)) ;
+ COLAMD_ASSERT (Col [super_c].is_alive()) ;
COLAMD_ASSERT (Col [super_c].shared3.hash == hash) ;
length = Col [super_c].length ;
@@ -1614,10 +1633,10 @@ static void detect_super_cols
/* === Compare super_c with all columns after it ================ */
for (c = Col [super_c].shared4.hash_next ;
- c != COLAMD_EMPTY ; c = Col [c].shared4.hash_next)
+ c != Empty ; c = Col [c].shared4.hash_next)
{
COLAMD_ASSERT (c != super_c) ;
- COLAMD_ASSERT (COL_IS_ALIVE (c)) ;
+ COLAMD_ASSERT (Col[c].is_alive()) ;
COLAMD_ASSERT (Col [c].shared3.hash == hash) ;
/* not identical if lengths or scores are different */
@@ -1635,8 +1654,8 @@ static void detect_super_cols
for (i = 0 ; i < length ; i++)
{
/* the columns are "clean" (no dead rows) */
- COLAMD_ASSERT (ROW_IS_ALIVE (*cp1)) ;
- COLAMD_ASSERT (ROW_IS_ALIVE (*cp2)) ;
+ COLAMD_ASSERT ( cp1->is_alive() );
+ COLAMD_ASSERT ( cp2->is_alive() );
/* row indices will same order for both supercols, */
/* no gather scatter necessary */
if (*cp1++ != *cp2++)
@@ -1658,9 +1677,9 @@ static void detect_super_cols
Col [super_c].shared1.thickness += Col [c].shared1.thickness ;
Col [c].shared1.parent = super_c ;
- KILL_NON_PRINCIPAL_COL (c) ;
+ Col[c].kill_non_principal() ;
/* order c later, in order_children() */
- Col [c].shared2.order = COLAMD_EMPTY ;
+ Col [c].shared2.order = Empty ;
/* remove c from hash bucket */
Col [prev_c].shared4.hash_next = Col [c].shared4.hash_next ;
}
@@ -1668,15 +1687,15 @@ static void detect_super_cols
/* === Empty this hash bucket ======================================= */
- if (head_column > COLAMD_EMPTY)
+ if (head_column > Empty)
{
/* corresponding degree list "hash" is not empty */
- Col [head_column].shared3.headhash = COLAMD_EMPTY ;
+ Col [head_column].shared3.headhash = Empty ;
}
else
{
/* corresponding degree list "hash" is empty */
- head [hash] = COLAMD_EMPTY ;
+ head [hash] = Empty ;
}
}
}
@@ -1698,11 +1717,11 @@ template <typename IndexType>
static IndexType garbage_collection /* returns the new value of pfree */
(
/* === Parameters ======================================================= */
-
+
IndexType n_row, /* number of rows */
IndexType n_col, /* number of columns */
- Colamd_Row<IndexType> Row [], /* row info */
- colamd_col<IndexType> Col [], /* column info */
+ RowStructure<IndexType> Row [], /* row info */
+ ColStructure<IndexType> Col [], /* column info */
IndexType A [], /* A [0 ... Alen-1] holds the matrix */
IndexType *pfree /* &A [0] ... pfree is in use */
)
@@ -1721,7 +1740,7 @@ static IndexType garbage_collection /* returns the new value of pfree */
pdest = &A[0] ;
for (c = 0 ; c < n_col ; c++)
{
- if (COL_IS_ALIVE (c))
+ if (Col[c].is_alive())
{
psrc = &A [Col [c].start] ;
@@ -1732,7 +1751,7 @@ static IndexType garbage_collection /* returns the new value of pfree */
for (j = 0 ; j < length ; j++)
{
r = *psrc++ ;
- if (ROW_IS_ALIVE (r))
+ if (Row[r].is_alive())
{
*pdest++ = r ;
}
@@ -1745,22 +1764,22 @@ static IndexType garbage_collection /* returns the new value of pfree */
for (r = 0 ; r < n_row ; r++)
{
- if (ROW_IS_ALIVE (r))
+ if (Row[r].is_alive())
{
if (Row [r].length == 0)
{
- /* this row is of zero length. cannot compact it, so kill it */
- COLAMD_DEBUG3 (("Defrag row kill\n")) ;
- KILL_ROW (r) ;
+ /* this row is of zero length. cannot compact it, so kill it */
+ COLAMD_DEBUG3 (("Defrag row kill\n")) ;
+ Row[r].kill() ;
}
else
{
- /* save first column index in Row [r].shared2.first_column */
- psrc = &A [Row [r].start] ;
- Row [r].shared2.first_column = *psrc ;
- COLAMD_ASSERT (ROW_IS_ALIVE (r)) ;
- /* flag the start of the row with the one's complement of row */
- *psrc = ONES_COMPLEMENT (r) ;
+ /* save first column index in Row [r].shared2.first_column */
+ psrc = &A [Row [r].start] ;
+ Row [r].shared2.first_column = *psrc ;
+ COLAMD_ASSERT (Row[r].is_alive()) ;
+ /* flag the start of the row with the one's complement of row */
+ *psrc = ones_complement(r) ;
}
}
@@ -1776,11 +1795,11 @@ static IndexType garbage_collection /* returns the new value of pfree */
{
psrc-- ;
/* get the row index */
- r = ONES_COMPLEMENT (*psrc) ;
+ r = ones_complement(*psrc) ;
COLAMD_ASSERT (r >= 0 && r < n_row) ;
/* restore first column index */
*psrc = Row [r].shared2.first_column ;
- COLAMD_ASSERT (ROW_IS_ALIVE (r)) ;
+ COLAMD_ASSERT (Row[r].is_alive()) ;
/* move and compact the row */
COLAMD_ASSERT (pdest <= psrc) ;
@@ -1789,7 +1808,7 @@ static IndexType garbage_collection /* returns the new value of pfree */
for (j = 0 ; j < length ; j++)
{
c = *psrc++ ;
- if (COL_IS_ALIVE (c))
+ if (Col[c].is_alive())
{
*pdest++ = c ;
}
@@ -1821,7 +1840,7 @@ static inline IndexType clear_mark /* return the new value for tag_mark */
/* === Parameters ======================================================= */
IndexType n_row, /* number of rows in A */
- Colamd_Row<IndexType> Row [] /* Row [0 ... n_row-1].shared2.mark is set to zero */
+ RowStructure<IndexType> Row [] /* Row [0 ... n_row-1].shared2.mark is set to zero */
)
{
/* === Local variables ================================================== */
@@ -1830,7 +1849,7 @@ static inline IndexType clear_mark /* return the new value for tag_mark */
for (r = 0 ; r < n_row ; r++)
{
- if (ROW_IS_ALIVE (r))
+ if (Row[r].is_alive())
{
Row [r].shared2.mark = 0 ;
}
@@ -1838,6 +1857,7 @@ static inline IndexType clear_mark /* return the new value for tag_mark */
return (1) ;
}
+} // namespace Colamd
-} // namespace internal
+} // namespace internal
#endif
diff --git a/Eigen/src/OrderingMethods/Ordering.h b/Eigen/src/OrderingMethods/Ordering.h
index 8791158be..c57897014 100644
--- a/Eigen/src/OrderingMethods/Ordering.h
+++ b/Eigen/src/OrderingMethods/Ordering.h
@@ -129,17 +129,17 @@ class COLAMDOrdering
StorageIndex n = StorageIndex(mat.cols());
StorageIndex nnz = StorageIndex(mat.nonZeros());
// Get the recommended value of Alen to be used by colamd
- StorageIndex Alen = internal::colamd_recommended(nnz, m, n);
+ StorageIndex Alen = internal::Colamd::recommended(nnz, m, n);
// Set the default parameters
- double knobs [COLAMD_KNOBS];
- StorageIndex stats [COLAMD_STATS];
- internal::colamd_set_defaults(knobs);
+ double knobs [internal::Colamd::NKnobs];
+ StorageIndex stats [internal::Colamd::NStats];
+ internal::Colamd::set_defaults(knobs);
IndexVector p(n+1), A(Alen);
for(StorageIndex i=0; i <= n; i++) p(i) = mat.outerIndexPtr()[i];
for(StorageIndex i=0; i < nnz; i++) A(i) = mat.innerIndexPtr()[i];
// Call Colamd routine to compute the ordering
- StorageIndex info = internal::colamd(m, n, Alen, A.data(), p.data(), knobs, stats);
+ StorageIndex info = internal::Colamd::compute_ordering(m, n, Alen, A.data(), p.data(), knobs, stats);
EIGEN_UNUSED_VARIABLE(info);
eigen_assert( info && "COLAMD failed " );
diff --git a/Eigen/src/PardisoSupport/PardisoSupport.h b/Eigen/src/PardisoSupport/PardisoSupport.h
index 07006b5c4..f89b79bd5 100644
--- a/Eigen/src/PardisoSupport/PardisoSupport.h
+++ b/Eigen/src/PardisoSupport/PardisoSupport.h
@@ -386,14 +386,15 @@ class PardisoLU : public PardisoImpl< PardisoLU<MatrixType> >
{
protected:
typedef PardisoImpl<PardisoLU> Base;
- typedef typename Base::Scalar Scalar;
- typedef typename Base::RealScalar RealScalar;
using Base::pardisoInit;
using Base::m_matrix;
friend class PardisoImpl< PardisoLU<MatrixType> >;
public:
+ typedef typename Base::Scalar Scalar;
+ typedef typename Base::RealScalar RealScalar;
+
using Base::compute;
using Base::solve;
@@ -441,14 +442,14 @@ class PardisoLLT : public PardisoImpl< PardisoLLT<MatrixType,_UpLo> >
{
protected:
typedef PardisoImpl< PardisoLLT<MatrixType,_UpLo> > Base;
- typedef typename Base::Scalar Scalar;
- typedef typename Base::RealScalar RealScalar;
using Base::pardisoInit;
using Base::m_matrix;
friend class PardisoImpl< PardisoLLT<MatrixType,_UpLo> >;
public:
+ typedef typename Base::Scalar Scalar;
+ typedef typename Base::RealScalar RealScalar;
typedef typename Base::StorageIndex StorageIndex;
enum { UpLo = _UpLo };
using Base::compute;
@@ -504,14 +505,14 @@ class PardisoLDLT : public PardisoImpl< PardisoLDLT<MatrixType,Options> >
{
protected:
typedef PardisoImpl< PardisoLDLT<MatrixType,Options> > Base;
- typedef typename Base::Scalar Scalar;
- typedef typename Base::RealScalar RealScalar;
using Base::pardisoInit;
using Base::m_matrix;
friend class PardisoImpl< PardisoLDLT<MatrixType,Options> >;
public:
+ typedef typename Base::Scalar Scalar;
+ typedef typename Base::RealScalar RealScalar;
typedef typename Base::StorageIndex StorageIndex;
using Base::compute;
enum { UpLo = Options&(Upper|Lower) };