Please, help us to better know about our user community by answering the following short survey: https://forms.gle/wpyrxWi18ox9Z5ae9
Eigen  3.4.0
 
Loading...
Searching...
No Matches
AVX/PacketMath.h
1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 2014 Benoit Steiner (benoit.steiner.goog@gmail.com)
5//
6// This Source Code Form is subject to the terms of the Mozilla
7// Public License v. 2.0. If a copy of the MPL was not distributed
8// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9
10#ifndef EIGEN_PACKET_MATH_AVX_H
11#define EIGEN_PACKET_MATH_AVX_H
12
13namespace Eigen {
14
15namespace internal {
16
17#ifndef EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD
18#define EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD 8
19#endif
20
21#if !defined(EIGEN_VECTORIZE_AVX512) && !defined(EIGEN_ARCH_DEFAULT_NUMBER_OF_REGISTERS)
22#define EIGEN_ARCH_DEFAULT_NUMBER_OF_REGISTERS 16
23#endif
24
25#ifdef EIGEN_VECTORIZE_FMA
26#ifndef EIGEN_HAS_SINGLE_INSTRUCTION_MADD
27#define EIGEN_HAS_SINGLE_INSTRUCTION_MADD
28#endif
29#endif
30
31typedef __m256 Packet8f;
32typedef __m256i Packet8i;
33typedef __m256d Packet4d;
34typedef eigen_packet_wrapper<__m128i, 2> Packet8h;
35typedef eigen_packet_wrapper<__m128i, 3> Packet8bf;
36
37template<> struct is_arithmetic<__m256> { enum { value = true }; };
38template<> struct is_arithmetic<__m256i> { enum { value = true }; };
39template<> struct is_arithmetic<__m256d> { enum { value = true }; };
40template<> struct is_arithmetic<Packet8h> { enum { value = true }; };
41template<> struct is_arithmetic<Packet8bf> { enum { value = true }; };
42
43#define _EIGEN_DECLARE_CONST_Packet8f(NAME,X) \
44 const Packet8f p8f_##NAME = pset1<Packet8f>(X)
45
46#define _EIGEN_DECLARE_CONST_Packet4d(NAME,X) \
47 const Packet4d p4d_##NAME = pset1<Packet4d>(X)
48
49#define _EIGEN_DECLARE_CONST_Packet8f_FROM_INT(NAME,X) \
50 const Packet8f p8f_##NAME = _mm256_castsi256_ps(pset1<Packet8i>(X))
51
52#define _EIGEN_DECLARE_CONST_Packet8i(NAME,X) \
53 const Packet8i p8i_##NAME = pset1<Packet8i>(X)
54
55// Use the packet_traits defined in AVX512/PacketMath.h instead if we're going
56// to leverage AVX512 instructions.
57#ifndef EIGEN_VECTORIZE_AVX512
58template<> struct packet_traits<float> : default_packet_traits
59{
60 typedef Packet8f type;
61 typedef Packet4f half;
62 enum {
63 Vectorizable = 1,
64 AlignedOnScalar = 1,
65 size = 8,
66 HasHalfPacket = 1,
67
68 HasCmp = 1,
69 HasDiv = 1,
70 HasSin = EIGEN_FAST_MATH,
71 HasCos = EIGEN_FAST_MATH,
72 HasLog = 1,
73 HasLog1p = 1,
74 HasExpm1 = 1,
75 HasExp = 1,
76 HasNdtri = 1,
77 HasBessel = 1,
78 HasSqrt = 1,
79 HasRsqrt = 1,
80 HasTanh = EIGEN_FAST_MATH,
81 HasErf = EIGEN_FAST_MATH,
82 HasBlend = 1,
83 HasRound = 1,
84 HasFloor = 1,
85 HasCeil = 1,
86 HasRint = 1
87 };
88};
89template<> struct packet_traits<double> : default_packet_traits
90{
91 typedef Packet4d type;
92 typedef Packet2d half;
93 enum {
94 Vectorizable = 1,
95 AlignedOnScalar = 1,
96 size=4,
97 HasHalfPacket = 1,
98
99 HasCmp = 1,
100 HasDiv = 1,
101 HasLog = 1,
102 HasExp = 1,
103 HasSqrt = 1,
104 HasRsqrt = 1,
105 HasBlend = 1,
106 HasRound = 1,
107 HasFloor = 1,
108 HasCeil = 1,
109 HasRint = 1
110 };
111};
112
113template <>
114struct packet_traits<Eigen::half> : default_packet_traits {
115 typedef Packet8h type;
116 // There is no half-size packet for Packet8h.
117 typedef Packet8h half;
118 enum {
119 Vectorizable = 1,
120 AlignedOnScalar = 1,
121 size = 8,
122 HasHalfPacket = 0,
123
124 HasCmp = 1,
125 HasAdd = 1,
126 HasSub = 1,
127 HasMul = 1,
128 HasDiv = 1,
129 HasSin = EIGEN_FAST_MATH,
130 HasCos = EIGEN_FAST_MATH,
131 HasNegate = 1,
132 HasAbs = 1,
133 HasAbs2 = 0,
134 HasMin = 1,
135 HasMax = 1,
136 HasConj = 1,
137 HasSetLinear = 0,
138 HasLog = 1,
139 HasLog1p = 1,
140 HasExpm1 = 1,
141 HasExp = 1,
142 HasSqrt = 1,
143 HasRsqrt = 1,
144 HasTanh = EIGEN_FAST_MATH,
145 HasErf = EIGEN_FAST_MATH,
146 HasBlend = 0,
147 HasRound = 1,
148 HasFloor = 1,
149 HasCeil = 1,
150 HasRint = 1,
151 HasBessel = 1,
152 HasNdtri = 1
153 };
154};
155
156template <>
157struct packet_traits<bfloat16> : default_packet_traits {
158 typedef Packet8bf type;
159 // There is no half-size packet for current Packet8bf.
160 // TODO: support as SSE path.
161 typedef Packet8bf half;
162 enum {
163 Vectorizable = 1,
164 AlignedOnScalar = 1,
165 size = 8,
166 HasHalfPacket = 0,
167
168 HasCmp = 1,
169 HasAdd = 1,
170 HasSub = 1,
171 HasMul = 1,
172 HasDiv = 1,
173 HasSin = EIGEN_FAST_MATH,
174 HasCos = EIGEN_FAST_MATH,
175 HasNegate = 1,
176 HasAbs = 1,
177 HasAbs2 = 0,
178 HasMin = 1,
179 HasMax = 1,
180 HasConj = 1,
181 HasSetLinear = 0,
182 HasLog = 1,
183 HasLog1p = 1,
184 HasExpm1 = 1,
185 HasExp = 1,
186 HasSqrt = 1,
187 HasRsqrt = 1,
188 HasTanh = EIGEN_FAST_MATH,
189 HasErf = EIGEN_FAST_MATH,
190 HasBlend = 0,
191 HasRound = 1,
192 HasFloor = 1,
193 HasCeil = 1,
194 HasRint = 1,
195 HasBessel = 1,
196 HasNdtri = 1
197 };
198};
199#endif
200
201template<> struct scalar_div_cost<float,true> { enum { value = 14 }; };
202template<> struct scalar_div_cost<double,true> { enum { value = 16 }; };
203
204/* Proper support for integers is only provided by AVX2. In the meantime, we'll
205 use SSE instructions and packets to deal with integers.
206template<> struct packet_traits<int> : default_packet_traits
207{
208 typedef Packet8i type;
209 enum {
210 Vectorizable = 1,
211 AlignedOnScalar = 1,
212 size=8
213 };
214};
215*/
216
217template<> struct unpacket_traits<Packet8f> {
218 typedef float type;
219 typedef Packet4f half;
220 typedef Packet8i integer_packet;
221 typedef uint8_t mask_t;
222 enum {size=8, alignment=Aligned32, vectorizable=true, masked_load_available=true, masked_store_available=true};
223};
224template<> struct unpacket_traits<Packet4d> {
225 typedef double type;
226 typedef Packet2d half;
227 enum {size=4, alignment=Aligned32, vectorizable=true, masked_load_available=false, masked_store_available=false};
228};
229template<> struct unpacket_traits<Packet8i> { typedef int type; typedef Packet4i half; enum {size=8, alignment=Aligned32, vectorizable=false, masked_load_available=false, masked_store_available=false}; };
230template<> struct unpacket_traits<Packet8bf> { typedef bfloat16 type; typedef Packet8bf half; enum {size=8, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; };
231
232// Helper function for bit packing snippet of low precision comparison.
233// It packs the flags from 16x16 to 8x16.
234EIGEN_STRONG_INLINE __m128i Pack16To8(Packet8f rf) {
235 return _mm_packs_epi32(_mm256_extractf128_si256(_mm256_castps_si256(rf), 0),
236 _mm256_extractf128_si256(_mm256_castps_si256(rf), 1));
237}
238
239
240template<> EIGEN_STRONG_INLINE Packet8f pset1<Packet8f>(const float& from) { return _mm256_set1_ps(from); }
241template<> EIGEN_STRONG_INLINE Packet4d pset1<Packet4d>(const double& from) { return _mm256_set1_pd(from); }
242template<> EIGEN_STRONG_INLINE Packet8i pset1<Packet8i>(const int& from) { return _mm256_set1_epi32(from); }
243
244template<> EIGEN_STRONG_INLINE Packet8f pset1frombits<Packet8f>(unsigned int from) { return _mm256_castsi256_ps(pset1<Packet8i>(from)); }
245template<> EIGEN_STRONG_INLINE Packet4d pset1frombits<Packet4d>(uint64_t from) { return _mm256_castsi256_pd(_mm256_set1_epi64x(from)); }
246
247template<> EIGEN_STRONG_INLINE Packet8f pzero(const Packet8f& /*a*/) { return _mm256_setzero_ps(); }
248template<> EIGEN_STRONG_INLINE Packet4d pzero(const Packet4d& /*a*/) { return _mm256_setzero_pd(); }
249template<> EIGEN_STRONG_INLINE Packet8i pzero(const Packet8i& /*a*/) { return _mm256_setzero_si256(); }
250
251
252template<> EIGEN_STRONG_INLINE Packet8f peven_mask(const Packet8f& /*a*/) { return _mm256_castsi256_ps(_mm256_set_epi32(0, -1, 0, -1, 0, -1, 0, -1)); }
253template<> EIGEN_STRONG_INLINE Packet8i peven_mask(const Packet8i& /*a*/) { return _mm256_set_epi32(0, -1, 0, -1, 0, -1, 0, -1); }
254template<> EIGEN_STRONG_INLINE Packet4d peven_mask(const Packet4d& /*a*/) { return _mm256_castsi256_pd(_mm256_set_epi32(0, 0, -1, -1, 0, 0, -1, -1)); }
255
256template<> EIGEN_STRONG_INLINE Packet8f pload1<Packet8f>(const float* from) { return _mm256_broadcast_ss(from); }
257template<> EIGEN_STRONG_INLINE Packet4d pload1<Packet4d>(const double* from) { return _mm256_broadcast_sd(from); }
258
259template<> EIGEN_STRONG_INLINE Packet8f plset<Packet8f>(const float& a) { return _mm256_add_ps(_mm256_set1_ps(a), _mm256_set_ps(7.0,6.0,5.0,4.0,3.0,2.0,1.0,0.0)); }
260template<> EIGEN_STRONG_INLINE Packet4d plset<Packet4d>(const double& a) { return _mm256_add_pd(_mm256_set1_pd(a), _mm256_set_pd(3.0,2.0,1.0,0.0)); }
261
262template<> EIGEN_STRONG_INLINE Packet8f padd<Packet8f>(const Packet8f& a, const Packet8f& b) { return _mm256_add_ps(a,b); }
263template<> EIGEN_STRONG_INLINE Packet4d padd<Packet4d>(const Packet4d& a, const Packet4d& b) { return _mm256_add_pd(a,b); }
264template<> EIGEN_STRONG_INLINE Packet8i padd<Packet8i>(const Packet8i& a, const Packet8i& b) {
265#ifdef EIGEN_VECTORIZE_AVX2
266 return _mm256_add_epi32(a,b);
267#else
268 __m128i lo = _mm_add_epi32(_mm256_extractf128_si256(a, 0), _mm256_extractf128_si256(b, 0));
269 __m128i hi = _mm_add_epi32(_mm256_extractf128_si256(a, 1), _mm256_extractf128_si256(b, 1));
270 return _mm256_insertf128_si256(_mm256_castsi128_si256(lo), (hi), 1);
271#endif
272}
273
274template<> EIGEN_STRONG_INLINE Packet8f psub<Packet8f>(const Packet8f& a, const Packet8f& b) { return _mm256_sub_ps(a,b); }
275template<> EIGEN_STRONG_INLINE Packet4d psub<Packet4d>(const Packet4d& a, const Packet4d& b) { return _mm256_sub_pd(a,b); }
276template<> EIGEN_STRONG_INLINE Packet8i psub<Packet8i>(const Packet8i& a, const Packet8i& b) {
277#ifdef EIGEN_VECTORIZE_AVX2
278 return _mm256_sub_epi32(a,b);
279#else
280 __m128i lo = _mm_sub_epi32(_mm256_extractf128_si256(a, 0), _mm256_extractf128_si256(b, 0));
281 __m128i hi = _mm_sub_epi32(_mm256_extractf128_si256(a, 1), _mm256_extractf128_si256(b, 1));
282 return _mm256_insertf128_si256(_mm256_castsi128_si256(lo), (hi), 1);
283#endif
284}
285
286template<> EIGEN_STRONG_INLINE Packet8f pnegate(const Packet8f& a)
287{
288 return _mm256_sub_ps(_mm256_set1_ps(0.0),a);
289}
290template<> EIGEN_STRONG_INLINE Packet4d pnegate(const Packet4d& a)
291{
292 return _mm256_sub_pd(_mm256_set1_pd(0.0),a);
293}
294
295template<> EIGEN_STRONG_INLINE Packet8f pconj(const Packet8f& a) { return a; }
296template<> EIGEN_STRONG_INLINE Packet4d pconj(const Packet4d& a) { return a; }
297template<> EIGEN_STRONG_INLINE Packet8i pconj(const Packet8i& a) { return a; }
298
299template<> EIGEN_STRONG_INLINE Packet8f pmul<Packet8f>(const Packet8f& a, const Packet8f& b) { return _mm256_mul_ps(a,b); }
300template<> EIGEN_STRONG_INLINE Packet4d pmul<Packet4d>(const Packet4d& a, const Packet4d& b) { return _mm256_mul_pd(a,b); }
301template<> EIGEN_STRONG_INLINE Packet8i pmul<Packet8i>(const Packet8i& a, const Packet8i& b) {
302#ifdef EIGEN_VECTORIZE_AVX2
303 return _mm256_mullo_epi32(a,b);
304#else
305 const __m128i lo = _mm_mullo_epi32(_mm256_extractf128_si256(a, 0), _mm256_extractf128_si256(b, 0));
306 const __m128i hi = _mm_mullo_epi32(_mm256_extractf128_si256(a, 1), _mm256_extractf128_si256(b, 1));
307 return _mm256_insertf128_si256(_mm256_castsi128_si256(lo), (hi), 1);
308#endif
309}
310
311template<> EIGEN_STRONG_INLINE Packet8f pdiv<Packet8f>(const Packet8f& a, const Packet8f& b) { return _mm256_div_ps(a,b); }
312template<> EIGEN_STRONG_INLINE Packet4d pdiv<Packet4d>(const Packet4d& a, const Packet4d& b) { return _mm256_div_pd(a,b); }
313template<> EIGEN_STRONG_INLINE Packet8i pdiv<Packet8i>(const Packet8i& /*a*/, const Packet8i& /*b*/)
314{ eigen_assert(false && "packet integer division are not supported by AVX");
315 return pset1<Packet8i>(0);
316}
317
318#ifdef EIGEN_VECTORIZE_FMA
319template<> EIGEN_STRONG_INLINE Packet8f pmadd(const Packet8f& a, const Packet8f& b, const Packet8f& c) {
320#if ( (EIGEN_COMP_GNUC_STRICT && EIGEN_COMP_GNUC<80) || (EIGEN_COMP_CLANG) )
321 // Clang stupidly generates a vfmadd213ps instruction plus some vmovaps on registers,
322 // and even register spilling with clang>=6.0 (bug 1637).
323 // Gcc stupidly generates a vfmadd132ps instruction.
324 // So let's enforce it to generate a vfmadd231ps instruction since the most common use
325 // case is to accumulate the result of the product.
326 Packet8f res = c;
327 __asm__("vfmadd231ps %[a], %[b], %[c]" : [c] "+x" (res) : [a] "x" (a), [b] "x" (b));
328 return res;
329#else
330 return _mm256_fmadd_ps(a,b,c);
331#endif
332}
333template<> EIGEN_STRONG_INLINE Packet4d pmadd(const Packet4d& a, const Packet4d& b, const Packet4d& c) {
334#if ( (EIGEN_COMP_GNUC_STRICT && EIGEN_COMP_GNUC<80) || (EIGEN_COMP_CLANG) )
335 // see above
336 Packet4d res = c;
337 __asm__("vfmadd231pd %[a], %[b], %[c]" : [c] "+x" (res) : [a] "x" (a), [b] "x" (b));
338 return res;
339#else
340 return _mm256_fmadd_pd(a,b,c);
341#endif
342}
343#endif
344
345template<> EIGEN_STRONG_INLINE Packet8f pcmp_le(const Packet8f& a, const Packet8f& b) { return _mm256_cmp_ps(a,b,_CMP_LE_OQ); }
346template<> EIGEN_STRONG_INLINE Packet8f pcmp_lt(const Packet8f& a, const Packet8f& b) { return _mm256_cmp_ps(a,b,_CMP_LT_OQ); }
347template<> EIGEN_STRONG_INLINE Packet8f pcmp_lt_or_nan(const Packet8f& a, const Packet8f& b) { return _mm256_cmp_ps(a, b, _CMP_NGE_UQ); }
348template<> EIGEN_STRONG_INLINE Packet8f pcmp_eq(const Packet8f& a, const Packet8f& b) { return _mm256_cmp_ps(a,b,_CMP_EQ_OQ); }
349
350template<> EIGEN_STRONG_INLINE Packet4d pcmp_le(const Packet4d& a, const Packet4d& b) { return _mm256_cmp_pd(a,b,_CMP_LE_OQ); }
351template<> EIGEN_STRONG_INLINE Packet4d pcmp_lt(const Packet4d& a, const Packet4d& b) { return _mm256_cmp_pd(a,b,_CMP_LT_OQ); }
352template<> EIGEN_STRONG_INLINE Packet4d pcmp_lt_or_nan(const Packet4d& a, const Packet4d& b) { return _mm256_cmp_pd(a, b, _CMP_NGE_UQ); }
353template<> EIGEN_STRONG_INLINE Packet4d pcmp_eq(const Packet4d& a, const Packet4d& b) { return _mm256_cmp_pd(a,b,_CMP_EQ_OQ); }
354
355
356template<> EIGEN_STRONG_INLINE Packet8i pcmp_eq(const Packet8i& a, const Packet8i& b) {
357#ifdef EIGEN_VECTORIZE_AVX2
358 return _mm256_cmpeq_epi32(a,b);
359#else
360 __m128i lo = _mm_cmpeq_epi32(_mm256_extractf128_si256(a, 0), _mm256_extractf128_si256(b, 0));
361 __m128i hi = _mm_cmpeq_epi32(_mm256_extractf128_si256(a, 1), _mm256_extractf128_si256(b, 1));
362 return _mm256_insertf128_si256(_mm256_castsi128_si256(lo), (hi), 1);
363#endif
364}
365
366template<> EIGEN_STRONG_INLINE Packet8f pmin<Packet8f>(const Packet8f& a, const Packet8f& b) {
367#if EIGEN_COMP_GNUC && EIGEN_COMP_GNUC < 63
368 // There appears to be a bug in GCC, by which the optimizer may flip
369 // the argument order in calls to _mm_min_ps/_mm_max_ps, so we have to
370 // resort to inline ASM here. This is supposed to be fixed in gcc6.3,
371 // see also: https://gcc.gnu.org/bugzilla/show_bug.cgi?id=72867
372 Packet8f res;
373 asm("vminps %[a], %[b], %[res]" : [res] "=x" (res) : [a] "x" (a), [b] "x" (b));
374 return res;
375#else
376 // Arguments are swapped to match NaN propagation behavior of std::min.
377 return _mm256_min_ps(b,a);
378#endif
379}
380template<> EIGEN_STRONG_INLINE Packet4d pmin<Packet4d>(const Packet4d& a, const Packet4d& b) {
381#if EIGEN_COMP_GNUC && EIGEN_COMP_GNUC < 63
382 // See pmin above
383 Packet4d res;
384 asm("vminpd %[a], %[b], %[res]" : [res] "=x" (res) : [a] "x" (a), [b] "x" (b));
385 return res;
386#else
387 // Arguments are swapped to match NaN propagation behavior of std::min.
388 return _mm256_min_pd(b,a);
389#endif
390}
391
392template<> EIGEN_STRONG_INLINE Packet8f pmax<Packet8f>(const Packet8f& a, const Packet8f& b) {
393#if EIGEN_COMP_GNUC && EIGEN_COMP_GNUC < 63
394 // See pmin above
395 Packet8f res;
396 asm("vmaxps %[a], %[b], %[res]" : [res] "=x" (res) : [a] "x" (a), [b] "x" (b));
397 return res;
398#else
399 // Arguments are swapped to match NaN propagation behavior of std::max.
400 return _mm256_max_ps(b,a);
401#endif
402}
403template<> EIGEN_STRONG_INLINE Packet4d pmax<Packet4d>(const Packet4d& a, const Packet4d& b) {
404#if EIGEN_COMP_GNUC && EIGEN_COMP_GNUC < 63
405 // See pmin above
406 Packet4d res;
407 asm("vmaxpd %[a], %[b], %[res]" : [res] "=x" (res) : [a] "x" (a), [b] "x" (b));
408 return res;
409#else
410 // Arguments are swapped to match NaN propagation behavior of std::max.
411 return _mm256_max_pd(b,a);
412#endif
413}
414
415// Add specializations for min/max with prescribed NaN progation.
416template<>
417EIGEN_STRONG_INLINE Packet8f pmin<PropagateNumbers, Packet8f>(const Packet8f& a, const Packet8f& b) {
418 return pminmax_propagate_numbers(a, b, pmin<Packet8f>);
419}
420template<>
421EIGEN_STRONG_INLINE Packet4d pmin<PropagateNumbers, Packet4d>(const Packet4d& a, const Packet4d& b) {
422 return pminmax_propagate_numbers(a, b, pmin<Packet4d>);
423}
424template<>
425EIGEN_STRONG_INLINE Packet8f pmax<PropagateNumbers, Packet8f>(const Packet8f& a, const Packet8f& b) {
426 return pminmax_propagate_numbers(a, b, pmax<Packet8f>);
427}
428template<>
429EIGEN_STRONG_INLINE Packet4d pmax<PropagateNumbers, Packet4d>(const Packet4d& a, const Packet4d& b) {
430 return pminmax_propagate_numbers(a, b, pmax<Packet4d>);
431}
432template<>
433EIGEN_STRONG_INLINE Packet8f pmin<PropagateNaN, Packet8f>(const Packet8f& a, const Packet8f& b) {
434 return pminmax_propagate_nan(a, b, pmin<Packet8f>);
435}
436template<>
437EIGEN_STRONG_INLINE Packet4d pmin<PropagateNaN, Packet4d>(const Packet4d& a, const Packet4d& b) {
438 return pminmax_propagate_nan(a, b, pmin<Packet4d>);
439}
440template<>
441EIGEN_STRONG_INLINE Packet8f pmax<PropagateNaN, Packet8f>(const Packet8f& a, const Packet8f& b) {
442 return pminmax_propagate_nan(a, b, pmax<Packet8f>);
443}
444template<>
445EIGEN_STRONG_INLINE Packet4d pmax<PropagateNaN, Packet4d>(const Packet4d& a, const Packet4d& b) {
446 return pminmax_propagate_nan(a, b, pmax<Packet4d>);
447}
448
449template<> EIGEN_STRONG_INLINE Packet8f print<Packet8f>(const Packet8f& a) { return _mm256_round_ps(a, _MM_FROUND_CUR_DIRECTION); }
450template<> EIGEN_STRONG_INLINE Packet4d print<Packet4d>(const Packet4d& a) { return _mm256_round_pd(a, _MM_FROUND_CUR_DIRECTION); }
451
452template<> EIGEN_STRONG_INLINE Packet8f pceil<Packet8f>(const Packet8f& a) { return _mm256_ceil_ps(a); }
453template<> EIGEN_STRONG_INLINE Packet4d pceil<Packet4d>(const Packet4d& a) { return _mm256_ceil_pd(a); }
454
455template<> EIGEN_STRONG_INLINE Packet8f pfloor<Packet8f>(const Packet8f& a) { return _mm256_floor_ps(a); }
456template<> EIGEN_STRONG_INLINE Packet4d pfloor<Packet4d>(const Packet4d& a) { return _mm256_floor_pd(a); }
457
458
459template<> EIGEN_STRONG_INLINE Packet8i ptrue<Packet8i>(const Packet8i& a) {
460#ifdef EIGEN_VECTORIZE_AVX2
461 // vpcmpeqd has lower latency than the more general vcmpps
462 return _mm256_cmpeq_epi32(a,a);
463#else
464 const __m256 b = _mm256_castsi256_ps(a);
465 return _mm256_castps_si256(_mm256_cmp_ps(b,b,_CMP_TRUE_UQ));
466#endif
467}
468
469template<> EIGEN_STRONG_INLINE Packet8f ptrue<Packet8f>(const Packet8f& a) {
470#ifdef EIGEN_VECTORIZE_AVX2
471 // vpcmpeqd has lower latency than the more general vcmpps
472 const __m256i b = _mm256_castps_si256(a);
473 return _mm256_castsi256_ps(_mm256_cmpeq_epi32(b,b));
474#else
475 return _mm256_cmp_ps(a,a,_CMP_TRUE_UQ);
476#endif
477}
478
479template<> EIGEN_STRONG_INLINE Packet4d ptrue<Packet4d>(const Packet4d& a) {
480#ifdef EIGEN_VECTORIZE_AVX2
481 // vpcmpeqq has lower latency than the more general vcmppd
482 const __m256i b = _mm256_castpd_si256(a);
483 return _mm256_castsi256_pd(_mm256_cmpeq_epi64(b,b));
484#else
485 return _mm256_cmp_pd(a,a,_CMP_TRUE_UQ);
486#endif
487}
488
489template<> EIGEN_STRONG_INLINE Packet8f pand<Packet8f>(const Packet8f& a, const Packet8f& b) { return _mm256_and_ps(a,b); }
490template<> EIGEN_STRONG_INLINE Packet4d pand<Packet4d>(const Packet4d& a, const Packet4d& b) { return _mm256_and_pd(a,b); }
491template<> EIGEN_STRONG_INLINE Packet8i pand<Packet8i>(const Packet8i& a, const Packet8i& b) {
492#ifdef EIGEN_VECTORIZE_AVX2
493 return _mm256_and_si256(a,b);
494#else
495 return _mm256_castps_si256(_mm256_and_ps(_mm256_castsi256_ps(a),_mm256_castsi256_ps(b)));
496#endif
497}
498
499template<> EIGEN_STRONG_INLINE Packet8f por<Packet8f>(const Packet8f& a, const Packet8f& b) { return _mm256_or_ps(a,b); }
500template<> EIGEN_STRONG_INLINE Packet4d por<Packet4d>(const Packet4d& a, const Packet4d& b) { return _mm256_or_pd(a,b); }
501template<> EIGEN_STRONG_INLINE Packet8i por<Packet8i>(const Packet8i& a, const Packet8i& b) {
502#ifdef EIGEN_VECTORIZE_AVX2
503 return _mm256_or_si256(a,b);
504#else
505 return _mm256_castps_si256(_mm256_or_ps(_mm256_castsi256_ps(a),_mm256_castsi256_ps(b)));
506#endif
507}
508
509template<> EIGEN_STRONG_INLINE Packet8f pxor<Packet8f>(const Packet8f& a, const Packet8f& b) { return _mm256_xor_ps(a,b); }
510template<> EIGEN_STRONG_INLINE Packet4d pxor<Packet4d>(const Packet4d& a, const Packet4d& b) { return _mm256_xor_pd(a,b); }
511template<> EIGEN_STRONG_INLINE Packet8i pxor<Packet8i>(const Packet8i& a, const Packet8i& b) {
512#ifdef EIGEN_VECTORIZE_AVX2
513 return _mm256_xor_si256(a,b);
514#else
515 return _mm256_castps_si256(_mm256_xor_ps(_mm256_castsi256_ps(a),_mm256_castsi256_ps(b)));
516#endif
517}
518
519template<> EIGEN_STRONG_INLINE Packet8f pandnot<Packet8f>(const Packet8f& a, const Packet8f& b) { return _mm256_andnot_ps(b,a); }
520template<> EIGEN_STRONG_INLINE Packet4d pandnot<Packet4d>(const Packet4d& a, const Packet4d& b) { return _mm256_andnot_pd(b,a); }
521template<> EIGEN_STRONG_INLINE Packet8i pandnot<Packet8i>(const Packet8i& a, const Packet8i& b) {
522#ifdef EIGEN_VECTORIZE_AVX2
523 return _mm256_andnot_si256(b,a);
524#else
525 return _mm256_castps_si256(_mm256_andnot_ps(_mm256_castsi256_ps(b),_mm256_castsi256_ps(a)));
526#endif
527}
528
529template<> EIGEN_STRONG_INLINE Packet8f pround<Packet8f>(const Packet8f& a)
530{
531 const Packet8f mask = pset1frombits<Packet8f>(static_cast<numext::uint32_t>(0x80000000u));
532 const Packet8f prev0dot5 = pset1frombits<Packet8f>(static_cast<numext::uint32_t>(0x3EFFFFFFu));
533 return _mm256_round_ps(padd(por(pand(a, mask), prev0dot5), a), _MM_FROUND_TO_ZERO);
534}
535template<> EIGEN_STRONG_INLINE Packet4d pround<Packet4d>(const Packet4d& a)
536{
537 const Packet4d mask = pset1frombits<Packet4d>(static_cast<numext::uint64_t>(0x8000000000000000ull));
538 const Packet4d prev0dot5 = pset1frombits<Packet4d>(static_cast<numext::uint64_t>(0x3FDFFFFFFFFFFFFFull));
539 return _mm256_round_pd(padd(por(pand(a, mask), prev0dot5), a), _MM_FROUND_TO_ZERO);
540}
541
542template<> EIGEN_STRONG_INLINE Packet8f pselect<Packet8f>(const Packet8f& mask, const Packet8f& a, const Packet8f& b)
543{ return _mm256_blendv_ps(b,a,mask); }
544template<> EIGEN_STRONG_INLINE Packet4d pselect<Packet4d>(const Packet4d& mask, const Packet4d& a, const Packet4d& b)
545{ return _mm256_blendv_pd(b,a,mask); }
546
547template<int N> EIGEN_STRONG_INLINE Packet8i parithmetic_shift_right(Packet8i a) {
548#ifdef EIGEN_VECTORIZE_AVX2
549 return _mm256_srai_epi32(a, N);
550#else
551 __m128i lo = _mm_srai_epi32(_mm256_extractf128_si256(a, 0), N);
552 __m128i hi = _mm_srai_epi32(_mm256_extractf128_si256(a, 1), N);
553 return _mm256_insertf128_si256(_mm256_castsi128_si256(lo), (hi), 1);
554#endif
555}
556
557template<int N> EIGEN_STRONG_INLINE Packet8i plogical_shift_right(Packet8i a) {
558#ifdef EIGEN_VECTORIZE_AVX2
559 return _mm256_srli_epi32(a, N);
560#else
561 __m128i lo = _mm_srli_epi32(_mm256_extractf128_si256(a, 0), N);
562 __m128i hi = _mm_srli_epi32(_mm256_extractf128_si256(a, 1), N);
563 return _mm256_insertf128_si256(_mm256_castsi128_si256(lo), (hi), 1);
564#endif
565}
566
567template<int N> EIGEN_STRONG_INLINE Packet8i plogical_shift_left(Packet8i a) {
568#ifdef EIGEN_VECTORIZE_AVX2
569 return _mm256_slli_epi32(a, N);
570#else
571 __m128i lo = _mm_slli_epi32(_mm256_extractf128_si256(a, 0), N);
572 __m128i hi = _mm_slli_epi32(_mm256_extractf128_si256(a, 1), N);
573 return _mm256_insertf128_si256(_mm256_castsi128_si256(lo), (hi), 1);
574#endif
575}
576
577template<> EIGEN_STRONG_INLINE Packet8f pload<Packet8f>(const float* from) { EIGEN_DEBUG_ALIGNED_LOAD return _mm256_load_ps(from); }
578template<> EIGEN_STRONG_INLINE Packet4d pload<Packet4d>(const double* from) { EIGEN_DEBUG_ALIGNED_LOAD return _mm256_load_pd(from); }
579template<> EIGEN_STRONG_INLINE Packet8i pload<Packet8i>(const int* from) { EIGEN_DEBUG_ALIGNED_LOAD return _mm256_load_si256(reinterpret_cast<const __m256i*>(from)); }
580
581template<> EIGEN_STRONG_INLINE Packet8f ploadu<Packet8f>(const float* from) { EIGEN_DEBUG_UNALIGNED_LOAD return _mm256_loadu_ps(from); }
582template<> EIGEN_STRONG_INLINE Packet4d ploadu<Packet4d>(const double* from) { EIGEN_DEBUG_UNALIGNED_LOAD return _mm256_loadu_pd(from); }
583template<> EIGEN_STRONG_INLINE Packet8i ploadu<Packet8i>(const int* from) { EIGEN_DEBUG_UNALIGNED_LOAD return _mm256_loadu_si256(reinterpret_cast<const __m256i*>(from)); }
584
585template<> EIGEN_STRONG_INLINE Packet8f ploadu<Packet8f>(const float* from, uint8_t umask) {
586 Packet8i mask = _mm256_set1_epi8(static_cast<char>(umask));
587 const Packet8i bit_mask = _mm256_set_epi32(0xffffff7f, 0xffffffbf, 0xffffffdf, 0xffffffef, 0xfffffff7, 0xfffffffb, 0xfffffffd, 0xfffffffe);
588 mask = por<Packet8i>(mask, bit_mask);
589 mask = pcmp_eq<Packet8i>(mask, _mm256_set1_epi32(0xffffffff));
590 EIGEN_DEBUG_UNALIGNED_LOAD return _mm256_maskload_ps(from, mask);
591}
592
593// Loads 4 floats from memory a returns the packet {a0, a0 a1, a1, a2, a2, a3, a3}
594template<> EIGEN_STRONG_INLINE Packet8f ploaddup<Packet8f>(const float* from)
595{
596 // TODO try to find a way to avoid the need of a temporary register
597// Packet8f tmp = _mm256_castps128_ps256(_mm_loadu_ps(from));
598// tmp = _mm256_insertf128_ps(tmp, _mm_movehl_ps(_mm256_castps256_ps128(tmp),_mm256_castps256_ps128(tmp)), 1);
599// return _mm256_unpacklo_ps(tmp,tmp);
600
601 // _mm256_insertf128_ps is very slow on Haswell, thus:
602 Packet8f tmp = _mm256_broadcast_ps((const __m128*)(const void*)from);
603 // mimic an "inplace" permutation of the lower 128bits using a blend
604 tmp = _mm256_blend_ps(tmp,_mm256_castps128_ps256(_mm_permute_ps( _mm256_castps256_ps128(tmp), _MM_SHUFFLE(1,0,1,0))), 15);
605 // then we can perform a consistent permutation on the global register to get everything in shape:
606 return _mm256_permute_ps(tmp, _MM_SHUFFLE(3,3,2,2));
607}
608// Loads 2 doubles from memory a returns the packet {a0, a0 a1, a1}
609template<> EIGEN_STRONG_INLINE Packet4d ploaddup<Packet4d>(const double* from)
610{
611 Packet4d tmp = _mm256_broadcast_pd((const __m128d*)(const void*)from);
612 return _mm256_permute_pd(tmp, 3<<2);
613}
614
615// Loads 2 floats from memory a returns the packet {a0, a0 a0, a0, a1, a1, a1, a1}
616template<> EIGEN_STRONG_INLINE Packet8f ploadquad<Packet8f>(const float* from)
617{
618 Packet8f tmp = _mm256_castps128_ps256(_mm_broadcast_ss(from));
619 return _mm256_insertf128_ps(tmp, _mm_broadcast_ss(from+1), 1);
620}
621
622template<> EIGEN_STRONG_INLINE void pstore<float>(float* to, const Packet8f& from) { EIGEN_DEBUG_ALIGNED_STORE _mm256_store_ps(to, from); }
623template<> EIGEN_STRONG_INLINE void pstore<double>(double* to, const Packet4d& from) { EIGEN_DEBUG_ALIGNED_STORE _mm256_store_pd(to, from); }
624template<> EIGEN_STRONG_INLINE void pstore<int>(int* to, const Packet8i& from) { EIGEN_DEBUG_ALIGNED_STORE _mm256_storeu_si256(reinterpret_cast<__m256i*>(to), from); }
625
626template<> EIGEN_STRONG_INLINE void pstoreu<float>(float* to, const Packet8f& from) { EIGEN_DEBUG_UNALIGNED_STORE _mm256_storeu_ps(to, from); }
627template<> EIGEN_STRONG_INLINE void pstoreu<double>(double* to, const Packet4d& from) { EIGEN_DEBUG_UNALIGNED_STORE _mm256_storeu_pd(to, from); }
628template<> EIGEN_STRONG_INLINE void pstoreu<int>(int* to, const Packet8i& from) { EIGEN_DEBUG_UNALIGNED_STORE _mm256_storeu_si256(reinterpret_cast<__m256i*>(to), from); }
629
630template<> EIGEN_STRONG_INLINE void pstoreu<float>(float* to, const Packet8f& from, uint8_t umask) {
631 Packet8i mask = _mm256_set1_epi8(static_cast<char>(umask));
632 const Packet8i bit_mask = _mm256_set_epi32(0xffffff7f, 0xffffffbf, 0xffffffdf, 0xffffffef, 0xfffffff7, 0xfffffffb, 0xfffffffd, 0xfffffffe);
633 mask = por<Packet8i>(mask, bit_mask);
634 mask = pcmp_eq<Packet8i>(mask, _mm256_set1_epi32(0xffffffff));
635 EIGEN_DEBUG_UNALIGNED_STORE return _mm256_maskstore_ps(to, mask, from);
636}
637
638// NOTE: leverage _mm256_i32gather_ps and _mm256_i32gather_pd if AVX2 instructions are available
639// NOTE: for the record the following seems to be slower: return _mm256_i32gather_ps(from, _mm256_set1_epi32(stride), 4);
640template<> EIGEN_DEVICE_FUNC inline Packet8f pgather<float, Packet8f>(const float* from, Index stride)
641{
642 return _mm256_set_ps(from[7*stride], from[6*stride], from[5*stride], from[4*stride],
643 from[3*stride], from[2*stride], from[1*stride], from[0*stride]);
644}
645template<> EIGEN_DEVICE_FUNC inline Packet4d pgather<double, Packet4d>(const double* from, Index stride)
646{
647 return _mm256_set_pd(from[3*stride], from[2*stride], from[1*stride], from[0*stride]);
648}
649
650template<> EIGEN_DEVICE_FUNC inline void pscatter<float, Packet8f>(float* to, const Packet8f& from, Index stride)
651{
652 __m128 low = _mm256_extractf128_ps(from, 0);
653 to[stride*0] = _mm_cvtss_f32(low);
654 to[stride*1] = _mm_cvtss_f32(_mm_shuffle_ps(low, low, 1));
655 to[stride*2] = _mm_cvtss_f32(_mm_shuffle_ps(low, low, 2));
656 to[stride*3] = _mm_cvtss_f32(_mm_shuffle_ps(low, low, 3));
657
658 __m128 high = _mm256_extractf128_ps(from, 1);
659 to[stride*4] = _mm_cvtss_f32(high);
660 to[stride*5] = _mm_cvtss_f32(_mm_shuffle_ps(high, high, 1));
661 to[stride*6] = _mm_cvtss_f32(_mm_shuffle_ps(high, high, 2));
662 to[stride*7] = _mm_cvtss_f32(_mm_shuffle_ps(high, high, 3));
663}
664template<> EIGEN_DEVICE_FUNC inline void pscatter<double, Packet4d>(double* to, const Packet4d& from, Index stride)
665{
666 __m128d low = _mm256_extractf128_pd(from, 0);
667 to[stride*0] = _mm_cvtsd_f64(low);
668 to[stride*1] = _mm_cvtsd_f64(_mm_shuffle_pd(low, low, 1));
669 __m128d high = _mm256_extractf128_pd(from, 1);
670 to[stride*2] = _mm_cvtsd_f64(high);
671 to[stride*3] = _mm_cvtsd_f64(_mm_shuffle_pd(high, high, 1));
672}
673
674template<> EIGEN_STRONG_INLINE void pstore1<Packet8f>(float* to, const float& a)
675{
676 Packet8f pa = pset1<Packet8f>(a);
677 pstore(to, pa);
678}
679template<> EIGEN_STRONG_INLINE void pstore1<Packet4d>(double* to, const double& a)
680{
681 Packet4d pa = pset1<Packet4d>(a);
682 pstore(to, pa);
683}
684template<> EIGEN_STRONG_INLINE void pstore1<Packet8i>(int* to, const int& a)
685{
686 Packet8i pa = pset1<Packet8i>(a);
687 pstore(to, pa);
688}
689
690#ifndef EIGEN_VECTORIZE_AVX512
691template<> EIGEN_STRONG_INLINE void prefetch<float>(const float* addr) { _mm_prefetch((SsePrefetchPtrType)(addr), _MM_HINT_T0); }
692template<> EIGEN_STRONG_INLINE void prefetch<double>(const double* addr) { _mm_prefetch((SsePrefetchPtrType)(addr), _MM_HINT_T0); }
693template<> EIGEN_STRONG_INLINE void prefetch<int>(const int* addr) { _mm_prefetch((SsePrefetchPtrType)(addr), _MM_HINT_T0); }
694#endif
695
696template<> EIGEN_STRONG_INLINE float pfirst<Packet8f>(const Packet8f& a) {
697 return _mm_cvtss_f32(_mm256_castps256_ps128(a));
698}
699template<> EIGEN_STRONG_INLINE double pfirst<Packet4d>(const Packet4d& a) {
700 return _mm_cvtsd_f64(_mm256_castpd256_pd128(a));
701}
702template<> EIGEN_STRONG_INLINE int pfirst<Packet8i>(const Packet8i& a) {
703 return _mm_cvtsi128_si32(_mm256_castsi256_si128(a));
704}
705
706
707template<> EIGEN_STRONG_INLINE Packet8f preverse(const Packet8f& a)
708{
709 __m256 tmp = _mm256_shuffle_ps(a,a,0x1b);
710 return _mm256_permute2f128_ps(tmp, tmp, 1);
711}
712template<> EIGEN_STRONG_INLINE Packet4d preverse(const Packet4d& a)
713{
714 __m256d tmp = _mm256_shuffle_pd(a,a,5);
715 return _mm256_permute2f128_pd(tmp, tmp, 1);
716 #if 0
717 // This version is unlikely to be faster as _mm256_shuffle_ps and _mm256_permute_pd
718 // exhibit the same latency/throughput, but it is here for future reference/benchmarking...
719 __m256d swap_halves = _mm256_permute2f128_pd(a,a,1);
720 return _mm256_permute_pd(swap_halves,5);
721 #endif
722}
723
724// pabs should be ok
725template<> EIGEN_STRONG_INLINE Packet8f pabs(const Packet8f& a)
726{
727 const Packet8f mask = _mm256_castsi256_ps(_mm256_setr_epi32(0x7FFFFFFF,0x7FFFFFFF,0x7FFFFFFF,0x7FFFFFFF,0x7FFFFFFF,0x7FFFFFFF,0x7FFFFFFF,0x7FFFFFFF));
728 return _mm256_and_ps(a,mask);
729}
730template<> EIGEN_STRONG_INLINE Packet4d pabs(const Packet4d& a)
731{
732 const Packet4d mask = _mm256_castsi256_pd(_mm256_setr_epi32(0xFFFFFFFF,0x7FFFFFFF,0xFFFFFFFF,0x7FFFFFFF,0xFFFFFFFF,0x7FFFFFFF,0xFFFFFFFF,0x7FFFFFFF));
733 return _mm256_and_pd(a,mask);
734}
735
736template<> EIGEN_STRONG_INLINE Packet8f pfrexp<Packet8f>(const Packet8f& a, Packet8f& exponent) {
737 return pfrexp_generic(a,exponent);
738}
739
740// Extract exponent without existence of Packet4l.
741template<>
742EIGEN_STRONG_INLINE
743Packet4d pfrexp_generic_get_biased_exponent(const Packet4d& a) {
744 const Packet4d cst_exp_mask = pset1frombits<Packet4d>(static_cast<uint64_t>(0x7ff0000000000000ull));
745 __m256i a_expo = _mm256_castpd_si256(pand(a, cst_exp_mask));
746#ifdef EIGEN_VECTORIZE_AVX2
747 a_expo = _mm256_srli_epi64(a_expo, 52);
748 __m128i lo = _mm256_extractf128_si256(a_expo, 0);
749 __m128i hi = _mm256_extractf128_si256(a_expo, 1);
750#else
751 __m128i lo = _mm256_extractf128_si256(a_expo, 0);
752 __m128i hi = _mm256_extractf128_si256(a_expo, 1);
753 lo = _mm_srli_epi64(lo, 52);
754 hi = _mm_srli_epi64(hi, 52);
755#endif
756 Packet2d exponent_lo = _mm_cvtepi32_pd(vec4i_swizzle1(lo, 0, 2, 1, 3));
757 Packet2d exponent_hi = _mm_cvtepi32_pd(vec4i_swizzle1(hi, 0, 2, 1, 3));
758 Packet4d exponent = _mm256_insertf128_pd(_mm256_setzero_pd(), exponent_lo, 0);
759 exponent = _mm256_insertf128_pd(exponent, exponent_hi, 1);
760 return exponent;
761}
762
763
764template<> EIGEN_STRONG_INLINE Packet4d pfrexp<Packet4d>(const Packet4d& a, Packet4d& exponent) {
765 return pfrexp_generic(a, exponent);
766}
767
768template<> EIGEN_STRONG_INLINE Packet8f pldexp<Packet8f>(const Packet8f& a, const Packet8f& exponent) {
769 return pldexp_generic(a, exponent);
770}
771
772template<> EIGEN_STRONG_INLINE Packet4d pldexp<Packet4d>(const Packet4d& a, const Packet4d& exponent) {
773 // Clamp exponent to [-2099, 2099]
774 const Packet4d max_exponent = pset1<Packet4d>(2099.0);
775 const Packet4i e = _mm256_cvtpd_epi32(pmin(pmax(exponent, pnegate(max_exponent)), max_exponent));
776
777 // Split 2^e into four factors and multiply.
778 const Packet4i bias = pset1<Packet4i>(1023);
779 Packet4i b = parithmetic_shift_right<2>(e); // floor(e/4)
780
781 // 2^b
782 Packet4i hi = vec4i_swizzle1(padd(b, bias), 0, 2, 1, 3);
783 Packet4i lo = _mm_slli_epi64(hi, 52);
784 hi = _mm_slli_epi64(_mm_srli_epi64(hi, 32), 52);
785 Packet4d c = _mm256_castsi256_pd(_mm256_insertf128_si256(_mm256_castsi128_si256(lo), hi, 1));
786 Packet4d out = pmul(pmul(pmul(a, c), c), c); // a * 2^(3b)
787
788 // 2^(e - 3b)
789 b = psub(psub(psub(e, b), b), b); // e - 3b
790 hi = vec4i_swizzle1(padd(b, bias), 0, 2, 1, 3);
791 lo = _mm_slli_epi64(hi, 52);
792 hi = _mm_slli_epi64(_mm_srli_epi64(hi, 32), 52);
793 c = _mm256_castsi256_pd(_mm256_insertf128_si256(_mm256_castsi128_si256(lo), hi, 1));
794 out = pmul(out, c); // a * 2^e
795 return out;
796}
797
798template<> EIGEN_STRONG_INLINE float predux<Packet8f>(const Packet8f& a)
799{
800 return predux(Packet4f(_mm_add_ps(_mm256_castps256_ps128(a),_mm256_extractf128_ps(a,1))));
801}
802template<> EIGEN_STRONG_INLINE double predux<Packet4d>(const Packet4d& a)
803{
804 return predux(Packet2d(_mm_add_pd(_mm256_castpd256_pd128(a),_mm256_extractf128_pd(a,1))));
805}
806
807template<> EIGEN_STRONG_INLINE Packet4f predux_half_dowto4<Packet8f>(const Packet8f& a)
808{
809 return _mm_add_ps(_mm256_castps256_ps128(a),_mm256_extractf128_ps(a,1));
810}
811
812template<> EIGEN_STRONG_INLINE float predux_mul<Packet8f>(const Packet8f& a)
813{
814 Packet8f tmp;
815 tmp = _mm256_mul_ps(a, _mm256_permute2f128_ps(a,a,1));
816 tmp = _mm256_mul_ps(tmp, _mm256_shuffle_ps(tmp,tmp,_MM_SHUFFLE(1,0,3,2)));
817 return pfirst(_mm256_mul_ps(tmp, _mm256_shuffle_ps(tmp,tmp,1)));
818}
819template<> EIGEN_STRONG_INLINE double predux_mul<Packet4d>(const Packet4d& a)
820{
821 Packet4d tmp;
822 tmp = _mm256_mul_pd(a, _mm256_permute2f128_pd(a,a,1));
823 return pfirst(_mm256_mul_pd(tmp, _mm256_shuffle_pd(tmp,tmp,1)));
824}
825
826template<> EIGEN_STRONG_INLINE float predux_min<Packet8f>(const Packet8f& a)
827{
828 Packet8f tmp = _mm256_min_ps(a, _mm256_permute2f128_ps(a,a,1));
829 tmp = _mm256_min_ps(tmp, _mm256_shuffle_ps(tmp,tmp,_MM_SHUFFLE(1,0,3,2)));
830 return pfirst(_mm256_min_ps(tmp, _mm256_shuffle_ps(tmp,tmp,1)));
831}
832template<> EIGEN_STRONG_INLINE double predux_min<Packet4d>(const Packet4d& a)
833{
834 Packet4d tmp = _mm256_min_pd(a, _mm256_permute2f128_pd(a,a,1));
835 return pfirst(_mm256_min_pd(tmp, _mm256_shuffle_pd(tmp, tmp, 1)));
836}
837
838template<> EIGEN_STRONG_INLINE float predux_max<Packet8f>(const Packet8f& a)
839{
840 Packet8f tmp = _mm256_max_ps(a, _mm256_permute2f128_ps(a,a,1));
841 tmp = _mm256_max_ps(tmp, _mm256_shuffle_ps(tmp,tmp,_MM_SHUFFLE(1,0,3,2)));
842 return pfirst(_mm256_max_ps(tmp, _mm256_shuffle_ps(tmp,tmp,1)));
843}
844
845template<> EIGEN_STRONG_INLINE double predux_max<Packet4d>(const Packet4d& a)
846{
847 Packet4d tmp = _mm256_max_pd(a, _mm256_permute2f128_pd(a,a,1));
848 return pfirst(_mm256_max_pd(tmp, _mm256_shuffle_pd(tmp, tmp, 1)));
849}
850
851// not needed yet
852// template<> EIGEN_STRONG_INLINE bool predux_all(const Packet8f& x)
853// {
854// return _mm256_movemask_ps(x)==0xFF;
855// }
856
857template<> EIGEN_STRONG_INLINE bool predux_any(const Packet8f& x)
858{
859 return _mm256_movemask_ps(x)!=0;
860}
861
862EIGEN_DEVICE_FUNC inline void
863ptranspose(PacketBlock<Packet8f,8>& kernel) {
864 __m256 T0 = _mm256_unpacklo_ps(kernel.packet[0], kernel.packet[1]);
865 __m256 T1 = _mm256_unpackhi_ps(kernel.packet[0], kernel.packet[1]);
866 __m256 T2 = _mm256_unpacklo_ps(kernel.packet[2], kernel.packet[3]);
867 __m256 T3 = _mm256_unpackhi_ps(kernel.packet[2], kernel.packet[3]);
868 __m256 T4 = _mm256_unpacklo_ps(kernel.packet[4], kernel.packet[5]);
869 __m256 T5 = _mm256_unpackhi_ps(kernel.packet[4], kernel.packet[5]);
870 __m256 T6 = _mm256_unpacklo_ps(kernel.packet[6], kernel.packet[7]);
871 __m256 T7 = _mm256_unpackhi_ps(kernel.packet[6], kernel.packet[7]);
872 __m256 S0 = _mm256_shuffle_ps(T0,T2,_MM_SHUFFLE(1,0,1,0));
873 __m256 S1 = _mm256_shuffle_ps(T0,T2,_MM_SHUFFLE(3,2,3,2));
874 __m256 S2 = _mm256_shuffle_ps(T1,T3,_MM_SHUFFLE(1,0,1,0));
875 __m256 S3 = _mm256_shuffle_ps(T1,T3,_MM_SHUFFLE(3,2,3,2));
876 __m256 S4 = _mm256_shuffle_ps(T4,T6,_MM_SHUFFLE(1,0,1,0));
877 __m256 S5 = _mm256_shuffle_ps(T4,T6,_MM_SHUFFLE(3,2,3,2));
878 __m256 S6 = _mm256_shuffle_ps(T5,T7,_MM_SHUFFLE(1,0,1,0));
879 __m256 S7 = _mm256_shuffle_ps(T5,T7,_MM_SHUFFLE(3,2,3,2));
880 kernel.packet[0] = _mm256_permute2f128_ps(S0, S4, 0x20);
881 kernel.packet[1] = _mm256_permute2f128_ps(S1, S5, 0x20);
882 kernel.packet[2] = _mm256_permute2f128_ps(S2, S6, 0x20);
883 kernel.packet[3] = _mm256_permute2f128_ps(S3, S7, 0x20);
884 kernel.packet[4] = _mm256_permute2f128_ps(S0, S4, 0x31);
885 kernel.packet[5] = _mm256_permute2f128_ps(S1, S5, 0x31);
886 kernel.packet[6] = _mm256_permute2f128_ps(S2, S6, 0x31);
887 kernel.packet[7] = _mm256_permute2f128_ps(S3, S7, 0x31);
888}
889
890EIGEN_DEVICE_FUNC inline void
891ptranspose(PacketBlock<Packet8f,4>& kernel) {
892 __m256 T0 = _mm256_unpacklo_ps(kernel.packet[0], kernel.packet[1]);
893 __m256 T1 = _mm256_unpackhi_ps(kernel.packet[0], kernel.packet[1]);
894 __m256 T2 = _mm256_unpacklo_ps(kernel.packet[2], kernel.packet[3]);
895 __m256 T3 = _mm256_unpackhi_ps(kernel.packet[2], kernel.packet[3]);
896
897 __m256 S0 = _mm256_shuffle_ps(T0,T2,_MM_SHUFFLE(1,0,1,0));
898 __m256 S1 = _mm256_shuffle_ps(T0,T2,_MM_SHUFFLE(3,2,3,2));
899 __m256 S2 = _mm256_shuffle_ps(T1,T3,_MM_SHUFFLE(1,0,1,0));
900 __m256 S3 = _mm256_shuffle_ps(T1,T3,_MM_SHUFFLE(3,2,3,2));
901
902 kernel.packet[0] = _mm256_permute2f128_ps(S0, S1, 0x20);
903 kernel.packet[1] = _mm256_permute2f128_ps(S2, S3, 0x20);
904 kernel.packet[2] = _mm256_permute2f128_ps(S0, S1, 0x31);
905 kernel.packet[3] = _mm256_permute2f128_ps(S2, S3, 0x31);
906}
907
908EIGEN_DEVICE_FUNC inline void
909ptranspose(PacketBlock<Packet4d,4>& kernel) {
910 __m256d T0 = _mm256_shuffle_pd(kernel.packet[0], kernel.packet[1], 15);
911 __m256d T1 = _mm256_shuffle_pd(kernel.packet[0], kernel.packet[1], 0);
912 __m256d T2 = _mm256_shuffle_pd(kernel.packet[2], kernel.packet[3], 15);
913 __m256d T3 = _mm256_shuffle_pd(kernel.packet[2], kernel.packet[3], 0);
914
915 kernel.packet[1] = _mm256_permute2f128_pd(T0, T2, 32);
916 kernel.packet[3] = _mm256_permute2f128_pd(T0, T2, 49);
917 kernel.packet[0] = _mm256_permute2f128_pd(T1, T3, 32);
918 kernel.packet[2] = _mm256_permute2f128_pd(T1, T3, 49);
919}
920
921template<> EIGEN_STRONG_INLINE Packet8f pblend(const Selector<8>& ifPacket, const Packet8f& thenPacket, const Packet8f& elsePacket) {
922 const __m256 zero = _mm256_setzero_ps();
923 const __m256 select = _mm256_set_ps(ifPacket.select[7], ifPacket.select[6], ifPacket.select[5], ifPacket.select[4], ifPacket.select[3], ifPacket.select[2], ifPacket.select[1], ifPacket.select[0]);
924 __m256 false_mask = _mm256_cmp_ps(select, zero, _CMP_EQ_UQ);
925 return _mm256_blendv_ps(thenPacket, elsePacket, false_mask);
926}
927template<> EIGEN_STRONG_INLINE Packet4d pblend(const Selector<4>& ifPacket, const Packet4d& thenPacket, const Packet4d& elsePacket) {
928 const __m256d zero = _mm256_setzero_pd();
929 const __m256d select = _mm256_set_pd(ifPacket.select[3], ifPacket.select[2], ifPacket.select[1], ifPacket.select[0]);
930 __m256d false_mask = _mm256_cmp_pd(select, zero, _CMP_EQ_UQ);
931 return _mm256_blendv_pd(thenPacket, elsePacket, false_mask);
932}
933
934// Packet math for Eigen::half
935
936template<> 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; };
937
938template<> EIGEN_STRONG_INLINE Packet8h pset1<Packet8h>(const Eigen::half& from) {
939 return _mm_set1_epi16(numext::bit_cast<numext::uint16_t>(from));
940}
941
942template<> EIGEN_STRONG_INLINE Eigen::half pfirst<Packet8h>(const Packet8h& from) {
943 return numext::bit_cast<Eigen::half>(static_cast<numext::uint16_t>(_mm_extract_epi16(from, 0)));
944}
945
946template<> EIGEN_STRONG_INLINE Packet8h pload<Packet8h>(const Eigen::half* from) {
947 return _mm_load_si128(reinterpret_cast<const __m128i*>(from));
948}
949
950template<> EIGEN_STRONG_INLINE Packet8h ploadu<Packet8h>(const Eigen::half* from) {
951 return _mm_loadu_si128(reinterpret_cast<const __m128i*>(from));
952}
953
954template<> EIGEN_STRONG_INLINE void pstore<Eigen::half>(Eigen::half* to, const Packet8h& from) {
955 _mm_store_si128(reinterpret_cast<__m128i*>(to), from);
956}
957
958template<> EIGEN_STRONG_INLINE void pstoreu<Eigen::half>(Eigen::half* to, const Packet8h& from) {
959 _mm_storeu_si128(reinterpret_cast<__m128i*>(to), from);
960}
961
962template<> EIGEN_STRONG_INLINE Packet8h
963ploaddup<Packet8h>(const Eigen::half* from) {
964 const numext::uint16_t a = numext::bit_cast<numext::uint16_t>(from[0]);
965 const numext::uint16_t b = numext::bit_cast<numext::uint16_t>(from[1]);
966 const numext::uint16_t c = numext::bit_cast<numext::uint16_t>(from[2]);
967 const numext::uint16_t d = numext::bit_cast<numext::uint16_t>(from[3]);
968 return _mm_set_epi16(d, d, c, c, b, b, a, a);
969}
970
971template<> EIGEN_STRONG_INLINE Packet8h
972ploadquad<Packet8h>(const Eigen::half* from) {
973 const numext::uint16_t a = numext::bit_cast<numext::uint16_t>(from[0]);
974 const numext::uint16_t b = numext::bit_cast<numext::uint16_t>(from[1]);
975 return _mm_set_epi16(b, b, b, b, a, a, a, a);
976}
977
978template<> EIGEN_STRONG_INLINE Packet8h ptrue(const Packet8h& a) {
979 return _mm_cmpeq_epi32(a, a);
980}
981
982template <>
983EIGEN_STRONG_INLINE Packet8h pabs(const Packet8h& a) {
984 const __m128i sign_mask = _mm_set1_epi16(static_cast<numext::uint16_t>(0x8000));
985 return _mm_andnot_si128(sign_mask, a);
986}
987
988EIGEN_STRONG_INLINE Packet8f half2float(const Packet8h& a) {
989#ifdef EIGEN_HAS_FP16_C
990 return _mm256_cvtph_ps(a);
991#else
992 EIGEN_ALIGN32 Eigen::half aux[8];
993 pstore(aux, a);
994 float f0(aux[0]);
995 float f1(aux[1]);
996 float f2(aux[2]);
997 float f3(aux[3]);
998 float f4(aux[4]);
999 float f5(aux[5]);
1000 float f6(aux[6]);
1001 float f7(aux[7]);
1002
1003 return _mm256_set_ps(f7, f6, f5, f4, f3, f2, f1, f0);
1004#endif
1005}
1006
1007EIGEN_STRONG_INLINE Packet8h float2half(const Packet8f& a) {
1008#ifdef EIGEN_HAS_FP16_C
1009 return _mm256_cvtps_ph(a, _MM_FROUND_TO_NEAREST_INT|_MM_FROUND_NO_EXC);
1010#else
1011 EIGEN_ALIGN32 float aux[8];
1012 pstore(aux, a);
1013 const numext::uint16_t s0 = numext::bit_cast<numext::uint16_t>(Eigen::half(aux[0]));
1014 const numext::uint16_t s1 = numext::bit_cast<numext::uint16_t>(Eigen::half(aux[1]));
1015 const numext::uint16_t s2 = numext::bit_cast<numext::uint16_t>(Eigen::half(aux[2]));
1016 const numext::uint16_t s3 = numext::bit_cast<numext::uint16_t>(Eigen::half(aux[3]));
1017 const numext::uint16_t s4 = numext::bit_cast<numext::uint16_t>(Eigen::half(aux[4]));
1018 const numext::uint16_t s5 = numext::bit_cast<numext::uint16_t>(Eigen::half(aux[5]));
1019 const numext::uint16_t s6 = numext::bit_cast<numext::uint16_t>(Eigen::half(aux[6]));
1020 const numext::uint16_t s7 = numext::bit_cast<numext::uint16_t>(Eigen::half(aux[7]));
1021 return _mm_set_epi16(s7, s6, s5, s4, s3, s2, s1, s0);
1022#endif
1023}
1024
1025template <>
1026EIGEN_STRONG_INLINE Packet8h pmin<Packet8h>(const Packet8h& a,
1027 const Packet8h& b) {
1028 return float2half(pmin<Packet8f>(half2float(a), half2float(b)));
1029}
1030
1031template <>
1032EIGEN_STRONG_INLINE Packet8h pmax<Packet8h>(const Packet8h& a,
1033 const Packet8h& b) {
1034 return float2half(pmax<Packet8f>(half2float(a), half2float(b)));
1035}
1036
1037template <>
1038EIGEN_STRONG_INLINE Packet8h plset<Packet8h>(const half& a) {
1039 return float2half(plset<Packet8f>(static_cast<float>(a)));
1040}
1041
1042template<> EIGEN_STRONG_INLINE Packet8h por(const Packet8h& a,const Packet8h& b) {
1043 // in some cases Packet4i is a wrapper around __m128i, so we either need to
1044 // cast to Packet4i to directly call the intrinsics as below:
1045 return _mm_or_si128(a,b);
1046}
1047template<> EIGEN_STRONG_INLINE Packet8h pxor(const Packet8h& a,const Packet8h& b) {
1048 return _mm_xor_si128(a,b);
1049}
1050template<> EIGEN_STRONG_INLINE Packet8h pand(const Packet8h& a,const Packet8h& b) {
1051 return _mm_and_si128(a,b);
1052}
1053template<> EIGEN_STRONG_INLINE Packet8h pandnot(const Packet8h& a,const Packet8h& b) {
1054 return _mm_andnot_si128(b,a);
1055}
1056
1057template<> EIGEN_STRONG_INLINE Packet8h pselect(const Packet8h& mask, const Packet8h& a, const Packet8h& b) {
1058 return _mm_blendv_epi8(b, a, mask);
1059}
1060
1061template<> EIGEN_STRONG_INLINE Packet8h pround<Packet8h>(const Packet8h& a) {
1062 return float2half(pround<Packet8f>(half2float(a)));
1063}
1064
1065template<> EIGEN_STRONG_INLINE Packet8h print<Packet8h>(const Packet8h& a) {
1066 return float2half(print<Packet8f>(half2float(a)));
1067}
1068
1069template<> EIGEN_STRONG_INLINE Packet8h pceil<Packet8h>(const Packet8h& a) {
1070 return float2half(pceil<Packet8f>(half2float(a)));
1071}
1072
1073template<> EIGEN_STRONG_INLINE Packet8h pfloor<Packet8h>(const Packet8h& a) {
1074 return float2half(pfloor<Packet8f>(half2float(a)));
1075}
1076
1077template<> EIGEN_STRONG_INLINE Packet8h pcmp_eq(const Packet8h& a,const Packet8h& b) {
1078 return Pack16To8(pcmp_eq(half2float(a), half2float(b)));
1079}
1080
1081template<> EIGEN_STRONG_INLINE Packet8h pcmp_le(const Packet8h& a,const Packet8h& b) {
1082 return Pack16To8(pcmp_le(half2float(a), half2float(b)));
1083}
1084
1085template<> EIGEN_STRONG_INLINE Packet8h pcmp_lt(const Packet8h& a,const Packet8h& b) {
1086 return Pack16To8(pcmp_lt(half2float(a), half2float(b)));
1087}
1088
1089template<> EIGEN_STRONG_INLINE Packet8h pcmp_lt_or_nan(const Packet8h& a,const Packet8h& b) {
1090 return Pack16To8(pcmp_lt_or_nan(half2float(a), half2float(b)));
1091}
1092
1093template<> EIGEN_STRONG_INLINE Packet8h pconj(const Packet8h& a) { return a; }
1094
1095template<> EIGEN_STRONG_INLINE Packet8h pnegate(const Packet8h& a) {
1096 Packet8h sign_mask = _mm_set1_epi16(static_cast<numext::uint16_t>(0x8000));
1097 return _mm_xor_si128(a, sign_mask);
1098}
1099
1100template<> EIGEN_STRONG_INLINE Packet8h padd<Packet8h>(const Packet8h& a, const Packet8h& b) {
1101 Packet8f af = half2float(a);
1102 Packet8f bf = half2float(b);
1103 Packet8f rf = padd(af, bf);
1104 return float2half(rf);
1105}
1106
1107template<> EIGEN_STRONG_INLINE Packet8h psub<Packet8h>(const Packet8h& a, const Packet8h& b) {
1108 Packet8f af = half2float(a);
1109 Packet8f bf = half2float(b);
1110 Packet8f rf = psub(af, bf);
1111 return float2half(rf);
1112}
1113
1114template<> EIGEN_STRONG_INLINE Packet8h pmul<Packet8h>(const Packet8h& a, const Packet8h& b) {
1115 Packet8f af = half2float(a);
1116 Packet8f bf = half2float(b);
1117 Packet8f rf = pmul(af, bf);
1118 return float2half(rf);
1119}
1120
1121template<> EIGEN_STRONG_INLINE Packet8h pdiv<Packet8h>(const Packet8h& a, const Packet8h& b) {
1122 Packet8f af = half2float(a);
1123 Packet8f bf = half2float(b);
1124 Packet8f rf = pdiv(af, bf);
1125 return float2half(rf);
1126}
1127
1128template<> EIGEN_STRONG_INLINE Packet8h pgather<Eigen::half, Packet8h>(const Eigen::half* from, Index stride)
1129{
1130 const numext::uint16_t s0 = numext::bit_cast<numext::uint16_t>(from[0*stride]);
1131 const numext::uint16_t s1 = numext::bit_cast<numext::uint16_t>(from[1*stride]);
1132 const numext::uint16_t s2 = numext::bit_cast<numext::uint16_t>(from[2*stride]);
1133 const numext::uint16_t s3 = numext::bit_cast<numext::uint16_t>(from[3*stride]);
1134 const numext::uint16_t s4 = numext::bit_cast<numext::uint16_t>(from[4*stride]);
1135 const numext::uint16_t s5 = numext::bit_cast<numext::uint16_t>(from[5*stride]);
1136 const numext::uint16_t s6 = numext::bit_cast<numext::uint16_t>(from[6*stride]);
1137 const numext::uint16_t s7 = numext::bit_cast<numext::uint16_t>(from[7*stride]);
1138 return _mm_set_epi16(s7, s6, s5, s4, s3, s2, s1, s0);
1139}
1140
1141template<> EIGEN_STRONG_INLINE void pscatter<Eigen::half, Packet8h>(Eigen::half* to, const Packet8h& from, Index stride)
1142{
1143 EIGEN_ALIGN32 Eigen::half aux[8];
1144 pstore(aux, from);
1145 to[stride*0] = aux[0];
1146 to[stride*1] = aux[1];
1147 to[stride*2] = aux[2];
1148 to[stride*3] = aux[3];
1149 to[stride*4] = aux[4];
1150 to[stride*5] = aux[5];
1151 to[stride*6] = aux[6];
1152 to[stride*7] = aux[7];
1153}
1154
1155template<> EIGEN_STRONG_INLINE Eigen::half predux<Packet8h>(const Packet8h& a) {
1156 Packet8f af = half2float(a);
1157 float reduced = predux<Packet8f>(af);
1158 return Eigen::half(reduced);
1159}
1160
1161template<> EIGEN_STRONG_INLINE Eigen::half predux_max<Packet8h>(const Packet8h& a) {
1162 Packet8f af = half2float(a);
1163 float reduced = predux_max<Packet8f>(af);
1164 return Eigen::half(reduced);
1165}
1166
1167template<> EIGEN_STRONG_INLINE Eigen::half predux_min<Packet8h>(const Packet8h& a) {
1168 Packet8f af = half2float(a);
1169 float reduced = predux_min<Packet8f>(af);
1170 return Eigen::half(reduced);
1171}
1172
1173template<> EIGEN_STRONG_INLINE Eigen::half predux_mul<Packet8h>(const Packet8h& a) {
1174 Packet8f af = half2float(a);
1175 float reduced = predux_mul<Packet8f>(af);
1176 return Eigen::half(reduced);
1177}
1178
1179template<> EIGEN_STRONG_INLINE Packet8h preverse(const Packet8h& a)
1180{
1181 __m128i m = _mm_setr_epi8(14,15,12,13,10,11,8,9,6,7,4,5,2,3,0,1);
1182 return _mm_shuffle_epi8(a,m);
1183}
1184
1185EIGEN_STRONG_INLINE void
1186ptranspose(PacketBlock<Packet8h,8>& kernel) {
1187 __m128i a = kernel.packet[0];
1188 __m128i b = kernel.packet[1];
1189 __m128i c = kernel.packet[2];
1190 __m128i d = kernel.packet[3];
1191 __m128i e = kernel.packet[4];
1192 __m128i f = kernel.packet[5];
1193 __m128i g = kernel.packet[6];
1194 __m128i h = kernel.packet[7];
1195
1196 __m128i a03b03 = _mm_unpacklo_epi16(a, b);
1197 __m128i c03d03 = _mm_unpacklo_epi16(c, d);
1198 __m128i e03f03 = _mm_unpacklo_epi16(e, f);
1199 __m128i g03h03 = _mm_unpacklo_epi16(g, h);
1200 __m128i a47b47 = _mm_unpackhi_epi16(a, b);
1201 __m128i c47d47 = _mm_unpackhi_epi16(c, d);
1202 __m128i e47f47 = _mm_unpackhi_epi16(e, f);
1203 __m128i g47h47 = _mm_unpackhi_epi16(g, h);
1204
1205 __m128i a01b01c01d01 = _mm_unpacklo_epi32(a03b03, c03d03);
1206 __m128i a23b23c23d23 = _mm_unpackhi_epi32(a03b03, c03d03);
1207 __m128i e01f01g01h01 = _mm_unpacklo_epi32(e03f03, g03h03);
1208 __m128i e23f23g23h23 = _mm_unpackhi_epi32(e03f03, g03h03);
1209 __m128i a45b45c45d45 = _mm_unpacklo_epi32(a47b47, c47d47);
1210 __m128i a67b67c67d67 = _mm_unpackhi_epi32(a47b47, c47d47);
1211 __m128i e45f45g45h45 = _mm_unpacklo_epi32(e47f47, g47h47);
1212 __m128i e67f67g67h67 = _mm_unpackhi_epi32(e47f47, g47h47);
1213
1214 __m128i a0b0c0d0e0f0g0h0 = _mm_unpacklo_epi64(a01b01c01d01, e01f01g01h01);
1215 __m128i a1b1c1d1e1f1g1h1 = _mm_unpackhi_epi64(a01b01c01d01, e01f01g01h01);
1216 __m128i a2b2c2d2e2f2g2h2 = _mm_unpacklo_epi64(a23b23c23d23, e23f23g23h23);
1217 __m128i a3b3c3d3e3f3g3h3 = _mm_unpackhi_epi64(a23b23c23d23, e23f23g23h23);
1218 __m128i a4b4c4d4e4f4g4h4 = _mm_unpacklo_epi64(a45b45c45d45, e45f45g45h45);
1219 __m128i a5b5c5d5e5f5g5h5 = _mm_unpackhi_epi64(a45b45c45d45, e45f45g45h45);
1220 __m128i a6b6c6d6e6f6g6h6 = _mm_unpacklo_epi64(a67b67c67d67, e67f67g67h67);
1221 __m128i a7b7c7d7e7f7g7h7 = _mm_unpackhi_epi64(a67b67c67d67, e67f67g67h67);
1222
1223 kernel.packet[0] = a0b0c0d0e0f0g0h0;
1224 kernel.packet[1] = a1b1c1d1e1f1g1h1;
1225 kernel.packet[2] = a2b2c2d2e2f2g2h2;
1226 kernel.packet[3] = a3b3c3d3e3f3g3h3;
1227 kernel.packet[4] = a4b4c4d4e4f4g4h4;
1228 kernel.packet[5] = a5b5c5d5e5f5g5h5;
1229 kernel.packet[6] = a6b6c6d6e6f6g6h6;
1230 kernel.packet[7] = a7b7c7d7e7f7g7h7;
1231}
1232
1233EIGEN_STRONG_INLINE void
1234ptranspose(PacketBlock<Packet8h,4>& kernel) {
1235 EIGEN_ALIGN32 Eigen::half in[4][8];
1236 pstore<Eigen::half>(in[0], kernel.packet[0]);
1237 pstore<Eigen::half>(in[1], kernel.packet[1]);
1238 pstore<Eigen::half>(in[2], kernel.packet[2]);
1239 pstore<Eigen::half>(in[3], kernel.packet[3]);
1240
1241 EIGEN_ALIGN32 Eigen::half out[4][8];
1242
1243 for (int i = 0; i < 4; ++i) {
1244 for (int j = 0; j < 4; ++j) {
1245 out[i][j] = in[j][2*i];
1246 }
1247 for (int j = 0; j < 4; ++j) {
1248 out[i][j+4] = in[j][2*i+1];
1249 }
1250 }
1251
1252 kernel.packet[0] = pload<Packet8h>(out[0]);
1253 kernel.packet[1] = pload<Packet8h>(out[1]);
1254 kernel.packet[2] = pload<Packet8h>(out[2]);
1255 kernel.packet[3] = pload<Packet8h>(out[3]);
1256}
1257
1258// BFloat16 implementation.
1259
1260EIGEN_STRONG_INLINE Packet8f Bf16ToF32(const Packet8bf& a) {
1261#ifdef EIGEN_VECTORIZE_AVX2
1262 __m256i extend = _mm256_cvtepu16_epi32(a);
1263 return _mm256_castsi256_ps(_mm256_slli_epi32(extend, 16));
1264#else
1265 __m128i lo = _mm_cvtepu16_epi32(a);
1266 __m128i hi = _mm_cvtepu16_epi32(_mm_srli_si128(a, 8));
1267 __m128i lo_shift = _mm_slli_epi32(lo, 16);
1268 __m128i hi_shift = _mm_slli_epi32(hi, 16);
1269 return _mm256_castsi256_ps(_mm256_insertf128_si256(_mm256_castsi128_si256(lo_shift), hi_shift, 1));
1270#endif
1271}
1272
1273// Convert float to bfloat16 according to round-to-nearest-even/denormals algorithm.
1274EIGEN_STRONG_INLINE Packet8bf F32ToBf16(const Packet8f& a) {
1275 Packet8bf r;
1276
1277 __m256i input = _mm256_castps_si256(a);
1278
1279#ifdef EIGEN_VECTORIZE_AVX2
1280 // uint32_t lsb = (input >> 16);
1281 __m256i t = _mm256_srli_epi32(input, 16);
1282 // uint32_t lsb = lsb & 1;
1283 t = _mm256_and_si256(t, _mm256_set1_epi32(1));
1284 // uint32_t rounding_bias = 0x7fff + lsb;
1285 t = _mm256_add_epi32(t, _mm256_set1_epi32(0x7fff));
1286 // input += rounding_bias;
1287 t = _mm256_add_epi32(t, input);
1288 // input = input >> 16;
1289 t = _mm256_srli_epi32(t, 16);
1290 // Check NaN before converting back to bf16
1291 __m256 mask = _mm256_cmp_ps(a, a, _CMP_ORD_Q);
1292 __m256i nan = _mm256_set1_epi32(0x7fc0);
1293 t = _mm256_blendv_epi8(nan, t, _mm256_castps_si256(mask));
1294 // output = numext::bit_cast<uint16_t>(input);
1295 return _mm_packus_epi32(_mm256_extractf128_si256(t, 0),
1296 _mm256_extractf128_si256(t, 1));
1297#else
1298 // uint32_t lsb = (input >> 16);
1299 __m128i lo = _mm_srli_epi32(_mm256_extractf128_si256(input, 0), 16);
1300 __m128i hi = _mm_srli_epi32(_mm256_extractf128_si256(input, 1), 16);
1301 // uint32_t lsb = lsb & 1;
1302 lo = _mm_and_si128(lo, _mm_set1_epi32(1));
1303 hi = _mm_and_si128(hi, _mm_set1_epi32(1));
1304 // uint32_t rounding_bias = 0x7fff + lsb;
1305 lo = _mm_add_epi32(lo, _mm_set1_epi32(0x7fff));
1306 hi = _mm_add_epi32(hi, _mm_set1_epi32(0x7fff));
1307 // input += rounding_bias;
1308 lo = _mm_add_epi32(lo, _mm256_extractf128_si256(input, 0));
1309 hi = _mm_add_epi32(hi, _mm256_extractf128_si256(input, 1));
1310 // input = input >> 16;
1311 lo = _mm_srli_epi32(lo, 16);
1312 hi = _mm_srli_epi32(hi, 16);
1313 // Check NaN before converting back to bf16
1314 __m256 mask = _mm256_cmp_ps(a, a, _CMP_ORD_Q);
1315 __m128i nan = _mm_set1_epi32(0x7fc0);
1316 lo = _mm_blendv_epi8(nan, lo, _mm_castps_si128(_mm256_castps256_ps128(mask)));
1317 hi = _mm_blendv_epi8(nan, hi, _mm_castps_si128(_mm256_extractf128_ps(mask, 1)));
1318 // output = numext::bit_cast<uint16_t>(input);
1319 return _mm_packus_epi32(lo, hi);
1320#endif
1321}
1322
1323template<> EIGEN_STRONG_INLINE Packet8bf pset1<Packet8bf>(const bfloat16& from) {
1324 return _mm_set1_epi16(numext::bit_cast<numext::uint16_t>(from));
1325}
1326
1327template<> EIGEN_STRONG_INLINE bfloat16 pfirst<Packet8bf>(const Packet8bf& from) {
1328 return numext::bit_cast<bfloat16>(static_cast<numext::uint16_t>(_mm_extract_epi16(from, 0)));
1329}
1330
1331template<> EIGEN_STRONG_INLINE Packet8bf pload<Packet8bf>(const bfloat16* from) {
1332 return _mm_load_si128(reinterpret_cast<const __m128i*>(from));
1333}
1334
1335template<> EIGEN_STRONG_INLINE Packet8bf ploadu<Packet8bf>(const bfloat16* from) {
1336 return _mm_loadu_si128(reinterpret_cast<const __m128i*>(from));
1337}
1338
1339template<> EIGEN_STRONG_INLINE void pstore<bfloat16>(bfloat16* to, const Packet8bf& from) {
1340 _mm_store_si128(reinterpret_cast<__m128i*>(to), from);
1341}
1342
1343template<> EIGEN_STRONG_INLINE void pstoreu<bfloat16>(bfloat16* to, const Packet8bf& from) {
1344 _mm_storeu_si128(reinterpret_cast<__m128i*>(to), from);
1345}
1346
1347template<> EIGEN_STRONG_INLINE Packet8bf
1348ploaddup<Packet8bf>(const bfloat16* from) {
1349 const numext::uint16_t a = numext::bit_cast<numext::uint16_t>(from[0]);
1350 const numext::uint16_t b = numext::bit_cast<numext::uint16_t>(from[1]);
1351 const numext::uint16_t c = numext::bit_cast<numext::uint16_t>(from[2]);
1352 const numext::uint16_t d = numext::bit_cast<numext::uint16_t>(from[3]);
1353 return _mm_set_epi16(d, d, c, c, b, b, a, a);
1354}
1355
1356template<> EIGEN_STRONG_INLINE Packet8bf
1357ploadquad<Packet8bf>(const bfloat16* from) {
1358 const numext::uint16_t a = numext::bit_cast<numext::uint16_t>(from[0]);
1359 const numext::uint16_t b = numext::bit_cast<numext::uint16_t>(from[1]);
1360 return _mm_set_epi16(b, b, b, b, a, a, a, a);
1361}
1362
1363template<> EIGEN_STRONG_INLINE Packet8bf ptrue(const Packet8bf& a) {
1364 return _mm_cmpeq_epi32(a, a);
1365}
1366
1367template <>
1368EIGEN_STRONG_INLINE Packet8bf pabs(const Packet8bf& a) {
1369 const __m128i sign_mask = _mm_set1_epi16(static_cast<numext::uint16_t>(0x8000));
1370 return _mm_andnot_si128(sign_mask, a);
1371}
1372
1373template <>
1374EIGEN_STRONG_INLINE Packet8bf pmin<Packet8bf>(const Packet8bf& a,
1375 const Packet8bf& b) {
1376 return F32ToBf16(pmin<Packet8f>(Bf16ToF32(a), Bf16ToF32(b)));
1377}
1378
1379template <>
1380EIGEN_STRONG_INLINE Packet8bf pmax<Packet8bf>(const Packet8bf& a,
1381 const Packet8bf& b) {
1382 return F32ToBf16(pmax<Packet8f>(Bf16ToF32(a), Bf16ToF32(b)));
1383}
1384
1385template <>
1386EIGEN_STRONG_INLINE Packet8bf plset<Packet8bf>(const bfloat16& a) {
1387 return F32ToBf16(plset<Packet8f>(static_cast<float>(a)));
1388}
1389
1390template<> EIGEN_STRONG_INLINE Packet8bf por(const Packet8bf& a,const Packet8bf& b) {
1391 return _mm_or_si128(a,b);
1392}
1393template<> EIGEN_STRONG_INLINE Packet8bf pxor(const Packet8bf& a,const Packet8bf& b) {
1394 return _mm_xor_si128(a,b);
1395}
1396template<> EIGEN_STRONG_INLINE Packet8bf pand(const Packet8bf& a,const Packet8bf& b) {
1397 return _mm_and_si128(a,b);
1398}
1399template<> EIGEN_STRONG_INLINE Packet8bf pandnot(const Packet8bf& a,const Packet8bf& b) {
1400 return _mm_andnot_si128(b,a);
1401}
1402
1403template<> EIGEN_STRONG_INLINE Packet8bf pselect(const Packet8bf& mask, const Packet8bf& a, const Packet8bf& b) {
1404 return _mm_blendv_epi8(b, a, mask);
1405}
1406
1407template<> EIGEN_STRONG_INLINE Packet8bf pround<Packet8bf>(const Packet8bf& a)
1408{
1409 return F32ToBf16(pround<Packet8f>(Bf16ToF32(a)));
1410}
1411
1412template<> EIGEN_STRONG_INLINE Packet8bf print<Packet8bf>(const Packet8bf& a) {
1413 return F32ToBf16(print<Packet8f>(Bf16ToF32(a)));
1414}
1415
1416template<> EIGEN_STRONG_INLINE Packet8bf pceil<Packet8bf>(const Packet8bf& a) {
1417 return F32ToBf16(pceil<Packet8f>(Bf16ToF32(a)));
1418}
1419
1420template<> EIGEN_STRONG_INLINE Packet8bf pfloor<Packet8bf>(const Packet8bf& a) {
1421 return F32ToBf16(pfloor<Packet8f>(Bf16ToF32(a)));
1422}
1423
1424template<> EIGEN_STRONG_INLINE Packet8bf pcmp_eq(const Packet8bf& a,const Packet8bf& b) {
1425 return Pack16To8(pcmp_eq(Bf16ToF32(a), Bf16ToF32(b)));
1426}
1427
1428template<> EIGEN_STRONG_INLINE Packet8bf pcmp_le(const Packet8bf& a,const Packet8bf& b) {
1429 return Pack16To8(pcmp_le(Bf16ToF32(a), Bf16ToF32(b)));
1430}
1431
1432template<> EIGEN_STRONG_INLINE Packet8bf pcmp_lt(const Packet8bf& a,const Packet8bf& b) {
1433 return Pack16To8(pcmp_lt(Bf16ToF32(a), Bf16ToF32(b)));
1434}
1435
1436template<> EIGEN_STRONG_INLINE Packet8bf pcmp_lt_or_nan(const Packet8bf& a,const Packet8bf& b) {
1437 return Pack16To8(pcmp_lt_or_nan(Bf16ToF32(a), Bf16ToF32(b)));
1438}
1439
1440template<> EIGEN_STRONG_INLINE Packet8bf pconj(const Packet8bf& a) { return a; }
1441
1442template<> EIGEN_STRONG_INLINE Packet8bf pnegate(const Packet8bf& a) {
1443 Packet8bf sign_mask = _mm_set1_epi16(static_cast<numext::uint16_t>(0x8000));
1444 return _mm_xor_si128(a, sign_mask);
1445}
1446
1447template<> EIGEN_STRONG_INLINE Packet8bf padd<Packet8bf>(const Packet8bf& a, const Packet8bf& b) {
1448 return F32ToBf16(padd<Packet8f>(Bf16ToF32(a), Bf16ToF32(b)));
1449}
1450
1451template<> EIGEN_STRONG_INLINE Packet8bf psub<Packet8bf>(const Packet8bf& a, const Packet8bf& b) {
1452 return F32ToBf16(psub<Packet8f>(Bf16ToF32(a), Bf16ToF32(b)));
1453}
1454
1455template<> EIGEN_STRONG_INLINE Packet8bf pmul<Packet8bf>(const Packet8bf& a, const Packet8bf& b) {
1456 return F32ToBf16(pmul<Packet8f>(Bf16ToF32(a), Bf16ToF32(b)));
1457}
1458
1459template<> EIGEN_STRONG_INLINE Packet8bf pdiv<Packet8bf>(const Packet8bf& a, const Packet8bf& b) {
1460 return F32ToBf16(pdiv<Packet8f>(Bf16ToF32(a), Bf16ToF32(b)));
1461}
1462
1463
1464template<> EIGEN_STRONG_INLINE Packet8bf pgather<bfloat16, Packet8bf>(const bfloat16* from, Index stride)
1465{
1466 const numext::uint16_t s0 = numext::bit_cast<numext::uint16_t>(from[0*stride]);
1467 const numext::uint16_t s1 = numext::bit_cast<numext::uint16_t>(from[1*stride]);
1468 const numext::uint16_t s2 = numext::bit_cast<numext::uint16_t>(from[2*stride]);
1469 const numext::uint16_t s3 = numext::bit_cast<numext::uint16_t>(from[3*stride]);
1470 const numext::uint16_t s4 = numext::bit_cast<numext::uint16_t>(from[4*stride]);
1471 const numext::uint16_t s5 = numext::bit_cast<numext::uint16_t>(from[5*stride]);
1472 const numext::uint16_t s6 = numext::bit_cast<numext::uint16_t>(from[6*stride]);
1473 const numext::uint16_t s7 = numext::bit_cast<numext::uint16_t>(from[7*stride]);
1474 return _mm_set_epi16(s7, s6, s5, s4, s3, s2, s1, s0);
1475}
1476
1477template<> EIGEN_STRONG_INLINE void pscatter<bfloat16, Packet8bf>(bfloat16* to, const Packet8bf& from, Index stride)
1478{
1479 EIGEN_ALIGN32 bfloat16 aux[8];
1480 pstore(aux, from);
1481 to[stride*0] = aux[0];
1482 to[stride*1] = aux[1];
1483 to[stride*2] = aux[2];
1484 to[stride*3] = aux[3];
1485 to[stride*4] = aux[4];
1486 to[stride*5] = aux[5];
1487 to[stride*6] = aux[6];
1488 to[stride*7] = aux[7];
1489}
1490
1491template<> EIGEN_STRONG_INLINE bfloat16 predux<Packet8bf>(const Packet8bf& a) {
1492 return static_cast<bfloat16>(predux<Packet8f>(Bf16ToF32(a)));
1493}
1494
1495template<> EIGEN_STRONG_INLINE bfloat16 predux_max<Packet8bf>(const Packet8bf& a) {
1496 return static_cast<bfloat16>(predux_max<Packet8f>(Bf16ToF32(a)));
1497}
1498
1499template<> EIGEN_STRONG_INLINE bfloat16 predux_min<Packet8bf>(const Packet8bf& a) {
1500 return static_cast<bfloat16>(predux_min<Packet8f>(Bf16ToF32(a)));
1501}
1502
1503template<> EIGEN_STRONG_INLINE bfloat16 predux_mul<Packet8bf>(const Packet8bf& a) {
1504 return static_cast<bfloat16>(predux_mul<Packet8f>(Bf16ToF32(a)));
1505}
1506
1507template<> EIGEN_STRONG_INLINE Packet8bf preverse(const Packet8bf& a)
1508{
1509 __m128i m = _mm_setr_epi8(14,15,12,13,10,11,8,9,6,7,4,5,2,3,0,1);
1510 return _mm_shuffle_epi8(a,m);
1511}
1512
1513EIGEN_STRONG_INLINE void
1514ptranspose(PacketBlock<Packet8bf,8>& kernel) {
1515 __m128i a = kernel.packet[0];
1516 __m128i b = kernel.packet[1];
1517 __m128i c = kernel.packet[2];
1518 __m128i d = kernel.packet[3];
1519 __m128i e = kernel.packet[4];
1520 __m128i f = kernel.packet[5];
1521 __m128i g = kernel.packet[6];
1522 __m128i h = kernel.packet[7];
1523
1524 __m128i a03b03 = _mm_unpacklo_epi16(a, b);
1525 __m128i c03d03 = _mm_unpacklo_epi16(c, d);
1526 __m128i e03f03 = _mm_unpacklo_epi16(e, f);
1527 __m128i g03h03 = _mm_unpacklo_epi16(g, h);
1528 __m128i a47b47 = _mm_unpackhi_epi16(a, b);
1529 __m128i c47d47 = _mm_unpackhi_epi16(c, d);
1530 __m128i e47f47 = _mm_unpackhi_epi16(e, f);
1531 __m128i g47h47 = _mm_unpackhi_epi16(g, h);
1532
1533 __m128i a01b01c01d01 = _mm_unpacklo_epi32(a03b03, c03d03);
1534 __m128i a23b23c23d23 = _mm_unpackhi_epi32(a03b03, c03d03);
1535 __m128i e01f01g01h01 = _mm_unpacklo_epi32(e03f03, g03h03);
1536 __m128i e23f23g23h23 = _mm_unpackhi_epi32(e03f03, g03h03);
1537 __m128i a45b45c45d45 = _mm_unpacklo_epi32(a47b47, c47d47);
1538 __m128i a67b67c67d67 = _mm_unpackhi_epi32(a47b47, c47d47);
1539 __m128i e45f45g45h45 = _mm_unpacklo_epi32(e47f47, g47h47);
1540 __m128i e67f67g67h67 = _mm_unpackhi_epi32(e47f47, g47h47);
1541
1542 kernel.packet[0] = _mm_unpacklo_epi64(a01b01c01d01, e01f01g01h01);
1543 kernel.packet[1] = _mm_unpackhi_epi64(a01b01c01d01, e01f01g01h01);
1544 kernel.packet[2] = _mm_unpacklo_epi64(a23b23c23d23, e23f23g23h23);
1545 kernel.packet[3] = _mm_unpackhi_epi64(a23b23c23d23, e23f23g23h23);
1546 kernel.packet[4] = _mm_unpacklo_epi64(a45b45c45d45, e45f45g45h45);
1547 kernel.packet[5] = _mm_unpackhi_epi64(a45b45c45d45, e45f45g45h45);
1548 kernel.packet[6] = _mm_unpacklo_epi64(a67b67c67d67, e67f67g67h67);
1549 kernel.packet[7] = _mm_unpackhi_epi64(a67b67c67d67, e67f67g67h67);
1550}
1551
1552EIGEN_STRONG_INLINE void
1553ptranspose(PacketBlock<Packet8bf,4>& kernel) {
1554 __m128i a = kernel.packet[0];
1555 __m128i b = kernel.packet[1];
1556 __m128i c = kernel.packet[2];
1557 __m128i d = kernel.packet[3];
1558
1559 __m128i ab_03 = _mm_unpacklo_epi16(a, b);
1560 __m128i cd_03 = _mm_unpacklo_epi16(c, d);
1561 __m128i ab_47 = _mm_unpackhi_epi16(a, b);
1562 __m128i cd_47 = _mm_unpackhi_epi16(c, d);
1563
1564 kernel.packet[0] = _mm_unpacklo_epi32(ab_03, cd_03);
1565 kernel.packet[1] = _mm_unpackhi_epi32(ab_03, cd_03);
1566 kernel.packet[2] = _mm_unpacklo_epi32(ab_47, cd_47);
1567 kernel.packet[3] = _mm_unpackhi_epi32(ab_47, cd_47);
1568}
1569
1570} // end namespace internal
1571
1572} // end namespace Eigen
1573
1574#endif // EIGEN_PACKET_MATH_AVX_H
@ Aligned32
Definition: Constants.h:236
@ Aligned16
Definition: Constants.h:235
Namespace containing all symbols from the Eigen library.
Definition: Core:141
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74