• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2016 Konstantinos Margaritis <markos@freevec.org>
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_ZVECTOR_H
11 #define EIGEN_PACKET_MATH_ZVECTOR_H
12 
13 #include <stdint.h>
14 
15 namespace Eigen {
16 
17 namespace internal {
18 
19 #ifndef EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD
20 #define EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD 4
21 #endif
22 
23 #ifndef EIGEN_HAS_SINGLE_INSTRUCTION_MADD
24 #define EIGEN_HAS_SINGLE_INSTRUCTION_MADD
25 #endif
26 
27 #ifndef EIGEN_HAS_SINGLE_INSTRUCTION_CJMADD
28 #define EIGEN_HAS_SINGLE_INSTRUCTION_CJMADD
29 #endif
30 
31 #ifndef EIGEN_ARCH_DEFAULT_NUMBER_OF_REGISTERS
32 #define EIGEN_ARCH_DEFAULT_NUMBER_OF_REGISTERS  16
33 #endif
34 
35 typedef __vector int                 Packet4i;
36 typedef __vector unsigned int        Packet4ui;
37 typedef __vector __bool int          Packet4bi;
38 typedef __vector short int           Packet8i;
39 typedef __vector unsigned char       Packet16uc;
40 typedef __vector double              Packet2d;
41 typedef __vector unsigned long long  Packet2ul;
42 typedef __vector long long           Packet2l;
43 
44 typedef struct {
45 	Packet2d  v4f[2];
46 } Packet4f;
47 
48 typedef union {
49   int32_t   i[4];
50   uint32_t ui[4];
51   int64_t   l[2];
52   uint64_t ul[2];
53   double    d[2];
54   Packet4i  v4i;
55   Packet4ui v4ui;
56   Packet2l  v2l;
57   Packet2ul v2ul;
58   Packet2d  v2d;
59 } Packet;
60 
61 // We don't want to write the same code all the time, but we need to reuse the constants
62 // and it doesn't really work to declare them global, so we define macros instead
63 
64 #define _EIGEN_DECLARE_CONST_FAST_Packet4i(NAME,X) \
65   Packet4i p4i_##NAME = reinterpret_cast<Packet4i>(vec_splat_s32(X))
66 
67 #define _EIGEN_DECLARE_CONST_FAST_Packet2d(NAME,X) \
68   Packet2d p2d_##NAME = reinterpret_cast<Packet2d>(vec_splat_s64(X))
69 
70 #define _EIGEN_DECLARE_CONST_FAST_Packet2l(NAME,X) \
71   Packet2l p2l_##NAME = reinterpret_cast<Packet2l>(vec_splat_s64(X))
72 
73 #define _EIGEN_DECLARE_CONST_Packet4i(NAME,X) \
74   Packet4i p4i_##NAME = pset1<Packet4i>(X)
75 
76 #define _EIGEN_DECLARE_CONST_Packet2d(NAME,X) \
77   Packet2d p2d_##NAME = pset1<Packet2d>(X)
78 
79 #define _EIGEN_DECLARE_CONST_Packet2l(NAME,X) \
80   Packet2l p2l_##NAME = pset1<Packet2l>(X)
81 
82 // These constants are endian-agnostic
83 //static _EIGEN_DECLARE_CONST_FAST_Packet4i(ZERO, 0); //{ 0, 0, 0, 0,}
84 static _EIGEN_DECLARE_CONST_FAST_Packet4i(ONE, 1); //{ 1, 1, 1, 1}
85 
86 static _EIGEN_DECLARE_CONST_FAST_Packet2d(ZERO, 0);
87 static _EIGEN_DECLARE_CONST_FAST_Packet2l(ZERO, 0);
88 static _EIGEN_DECLARE_CONST_FAST_Packet2l(ONE, 1);
89 
90 static Packet2d p2d_ONE = { 1.0, 1.0 };
91 static Packet2d p2d_ZERO_ = { -0.0, -0.0 };
92 
93 static Packet4i p4i_COUNTDOWN = { 0, 1, 2, 3 };
94 static Packet4f p4f_COUNTDOWN = { 0.0, 1.0, 2.0, 3.0 };
95 static Packet2d p2d_COUNTDOWN = reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet16uc>(p2d_ZERO), reinterpret_cast<Packet16uc>(p2d_ONE), 8));
96 
97 static Packet16uc p16uc_PSET64_HI = { 0,1,2,3, 4,5,6,7, 0,1,2,3, 4,5,6,7 };
98 static Packet16uc p16uc_DUPLICATE32_HI = { 0,1,2,3, 0,1,2,3, 4,5,6,7, 4,5,6,7 };
99 
100 // Mask alignment
101 #define _EIGEN_MASK_ALIGNMENT	0xfffffffffffffff0
102 
103 #define _EIGEN_ALIGNED_PTR(x)	((std::ptrdiff_t)(x) & _EIGEN_MASK_ALIGNMENT)
104 
105 // Handle endianness properly while loading constants
106 // Define global static constants:
107 
108 static Packet16uc p16uc_FORWARD =   { 0,1,2,3, 4,5,6,7, 8,9,10,11, 12,13,14,15 };
109 static Packet16uc p16uc_REVERSE32 = { 12,13,14,15, 8,9,10,11, 4,5,6,7, 0,1,2,3 };
110 static Packet16uc p16uc_REVERSE64 = { 8,9,10,11, 12,13,14,15, 0,1,2,3, 4,5,6,7 };
111 
112 static Packet16uc p16uc_PSET32_WODD   = vec_sld((Packet16uc) vec_splat((Packet4ui)p16uc_FORWARD, 0), (Packet16uc) vec_splat((Packet4ui)p16uc_FORWARD, 2), 8);//{ 0,1,2,3, 0,1,2,3, 8,9,10,11, 8,9,10,11 };
113 static Packet16uc p16uc_PSET32_WEVEN  = vec_sld(p16uc_DUPLICATE32_HI, (Packet16uc) vec_splat((Packet4ui)p16uc_FORWARD, 3), 8);//{ 4,5,6,7, 4,5,6,7, 12,13,14,15, 12,13,14,15 };
114 /*static Packet16uc p16uc_HALF64_0_16 = vec_sld((Packet16uc)p4i_ZERO, vec_splat((Packet16uc) vec_abs(p4i_MINUS16), 3), 8);      //{ 0,0,0,0, 0,0,0,0, 16,16,16,16, 16,16,16,16};
115 
116 static Packet16uc p16uc_PSET64_HI = (Packet16uc) vec_mergeh((Packet4ui)p16uc_PSET32_WODD, (Packet4ui)p16uc_PSET32_WEVEN);     //{ 0,1,2,3, 4,5,6,7, 0,1,2,3, 4,5,6,7 };*/
117 static Packet16uc p16uc_PSET64_LO = (Packet16uc) vec_mergel((Packet4ui)p16uc_PSET32_WODD, (Packet4ui)p16uc_PSET32_WEVEN);     //{ 8,9,10,11, 12,13,14,15, 8,9,10,11, 12,13,14,15 };
118 /*static Packet16uc p16uc_TRANSPOSE64_HI = vec_add(p16uc_PSET64_HI, p16uc_HALF64_0_16);                                         //{ 0,1,2,3, 4,5,6,7, 16,17,18,19, 20,21,22,23};
119 static Packet16uc p16uc_TRANSPOSE64_LO = vec_add(p16uc_PSET64_LO, p16uc_HALF64_0_16);                                         //{ 8,9,10,11, 12,13,14,15, 24,25,26,27, 28,29,30,31};*/
120 static Packet16uc p16uc_TRANSPOSE64_HI = { 0,1,2,3, 4,5,6,7, 16,17,18,19, 20,21,22,23};
121 static Packet16uc p16uc_TRANSPOSE64_LO = { 8,9,10,11, 12,13,14,15, 24,25,26,27, 28,29,30,31};
122 
123 //static Packet16uc p16uc_COMPLEX32_REV = vec_sld(p16uc_REVERSE32, p16uc_REVERSE32, 8);                                         //{ 4,5,6,7, 0,1,2,3, 12,13,14,15, 8,9,10,11 };
124 
125 //static Packet16uc p16uc_COMPLEX32_REV2 = vec_sld(p16uc_FORWARD, p16uc_FORWARD, 8);                                            //{ 8,9,10,11, 12,13,14,15, 0,1,2,3, 4,5,6,7 };
126 
127 
128 #if EIGEN_HAS_BUILTIN(__builtin_prefetch) || EIGEN_COMP_GNUC
129   #define EIGEN_ZVECTOR_PREFETCH(ADDR) __builtin_prefetch(ADDR);
130 #else
131   #define EIGEN_ZVECTOR_PREFETCH(ADDR) asm( "   pfd [%[addr]]\n" :: [addr] "r" (ADDR) : "cc" );
132 #endif
133 
134 template<> struct packet_traits<int>    : default_packet_traits
135 {
136   typedef Packet4i type;
137   typedef Packet4i half;
138   enum {
139     Vectorizable = 1,
140     AlignedOnScalar = 1,
141     size = 4,
142     HasHalfPacket = 0,
143 
144     HasAdd  = 1,
145     HasSub  = 1,
146     HasMul  = 1,
147     HasDiv  = 1,
148     HasBlend = 1
149   };
150 };
151 
152 template<> struct packet_traits<float> : default_packet_traits
153 {
154   typedef Packet4f type;
155   typedef Packet4f half;
156   enum {
157     Vectorizable = 1,
158     AlignedOnScalar = 1,
159     size=4,
160     HasHalfPacket = 0,
161 
162     HasAdd  = 1,
163     HasSub  = 1,
164     HasMul  = 1,
165     HasDiv  = 1,
166     HasMin  = 1,
167     HasMax  = 1,
168     HasAbs  = 1,
169     HasSin  = 0,
170     HasCos  = 0,
171     HasLog  = 0,
172     HasExp  = 1,
173     HasSqrt = 1,
174     HasRsqrt = 1,
175     HasRound = 1,
176     HasFloor = 1,
177     HasCeil = 1,
178     HasNegate = 1,
179     HasBlend = 1
180   };
181 };
182 
183 template<> struct packet_traits<double> : default_packet_traits
184 {
185   typedef Packet2d type;
186   typedef Packet2d half;
187   enum {
188     Vectorizable = 1,
189     AlignedOnScalar = 1,
190     size=2,
191     HasHalfPacket = 1,
192 
193     HasAdd  = 1,
194     HasSub  = 1,
195     HasMul  = 1,
196     HasDiv  = 1,
197     HasMin  = 1,
198     HasMax  = 1,
199     HasAbs  = 1,
200     HasSin  = 0,
201     HasCos  = 0,
202     HasLog  = 0,
203     HasExp  = 1,
204     HasSqrt = 1,
205     HasRsqrt = 1,
206     HasRound = 1,
207     HasFloor = 1,
208     HasCeil = 1,
209     HasNegate = 1,
210     HasBlend = 1
211   };
212 };
213 
214 template<> struct unpacket_traits<Packet4i> { typedef int    type; enum {size=4, alignment=Aligned16}; typedef Packet4i half; };
215 template<> struct unpacket_traits<Packet4f> { typedef float  type; enum {size=4, alignment=Aligned16}; typedef Packet4f half; };
216 template<> struct unpacket_traits<Packet2d> { typedef double type; enum {size=2, alignment=Aligned16}; typedef Packet2d half; };
217 
218 /* Forward declaration */
219 EIGEN_DEVICE_FUNC inline void ptranspose(PacketBlock<Packet4f,4>& kernel);
220 
221 inline std::ostream & operator <<(std::ostream & s, const Packet4i & v)
222 {
223   Packet vt;
224   vt.v4i = v;
225   s << vt.i[0] << ", " << vt.i[1] << ", " << vt.i[2] << ", " << vt.i[3];
226   return s;
227 }
228 
229 inline std::ostream & operator <<(std::ostream & s, const Packet4ui & v)
230 {
231   Packet vt;
232   vt.v4ui = v;
233   s << vt.ui[0] << ", " << vt.ui[1] << ", " << vt.ui[2] << ", " << vt.ui[3];
234   return s;
235 }
236 
237 inline std::ostream & operator <<(std::ostream & s, const Packet2l & v)
238 {
239   Packet vt;
240   vt.v2l = v;
241   s << vt.l[0] << ", " << vt.l[1];
242   return s;
243 }
244 
245 inline std::ostream & operator <<(std::ostream & s, const Packet2ul & v)
246 {
247   Packet vt;
248   vt.v2ul = v;
249   s << vt.ul[0] << ", " << vt.ul[1] ;
250   return s;
251 }
252 
253 inline std::ostream & operator <<(std::ostream & s, const Packet2d & v)
254 {
255   Packet vt;
256   vt.v2d = v;
257   s << vt.d[0] << ", " << vt.d[1];
258   return s;
259 }
260 
261 /* Helper function to simulate a vec_splat_packet4f
262  */
263 template<int element> EIGEN_STRONG_INLINE Packet4f vec_splat_packet4f(const Packet4f&   from)
264 {
265   Packet4f splat;
266   switch (element) {
267   case 0:
268     splat.v4f[0] = vec_splat(from.v4f[0], 0);
269     splat.v4f[1] = splat.v4f[0];
270     break;
271   case 1:
272     splat.v4f[0] = vec_splat(from.v4f[0], 1);
273     splat.v4f[1] = splat.v4f[0];
274     break;
275   case 2:
276     splat.v4f[0] = vec_splat(from.v4f[1], 0);
277     splat.v4f[1] = splat.v4f[0];
278     break;
279   case 3:
280     splat.v4f[0] = vec_splat(from.v4f[1], 1);
281     splat.v4f[1] = splat.v4f[0];
282     break;
283   }
284   return splat;
285 }
286 
287 template<int Offset>
288 struct palign_impl<Offset,Packet4i>
289 {
290   static EIGEN_STRONG_INLINE void run(Packet4i& first, const Packet4i& second)
291   {
292     switch (Offset % 4) {
293     case 1:
294       first = vec_sld(first, second, 4); break;
295     case 2:
296       first = vec_sld(first, second, 8); break;
297     case 3:
298       first = vec_sld(first, second, 12); break;
299     }
300   }
301 };
302 
303 /* This is a tricky one, we have to translate float alignment to vector elements of sizeof double
304  */
305 template<int Offset>
306 struct palign_impl<Offset,Packet4f>
307 {
308   static EIGEN_STRONG_INLINE void run(Packet4f& first, const Packet4f& second)
309   {
310     switch (Offset % 4) {
311     case 1:
312       first.v4f[0] = vec_sld(first.v4f[0], first.v4f[1], 8);
313       first.v4f[1] = vec_sld(first.v4f[1], second.v4f[0], 8);
314       break;
315     case 2:
316       first.v4f[0] = first.v4f[1];
317       first.v4f[1] = second.v4f[0];
318       break;
319     case 3:
320       first.v4f[0] = vec_sld(first.v4f[1],  second.v4f[0], 8);
321       first.v4f[1] = vec_sld(second.v4f[0], second.v4f[1], 8);
322       break;
323     }
324   }
325 };
326 
327 
328 template<int Offset>
329 struct palign_impl<Offset,Packet2d>
330 {
331   static EIGEN_STRONG_INLINE void run(Packet2d& first, const Packet2d& second)
332   {
333     if (Offset == 1)
334       first = reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4i>(first), reinterpret_cast<Packet4i>(second), 8));
335   }
336 };
337 
338 template<> EIGEN_STRONG_INLINE Packet4i pload<Packet4i>(const int*     from)
339 {
340   // FIXME: No intrinsic yet
341   EIGEN_DEBUG_ALIGNED_LOAD
342   Packet *vfrom;
343   vfrom = (Packet *) from;
344   return vfrom->v4i;
345 }
346 
347 template<> EIGEN_STRONG_INLINE Packet4f pload<Packet4f>(const float*   from)
348 {
349   // FIXME: No intrinsic yet
350   EIGEN_DEBUG_ALIGNED_LOAD
351   Packet4f vfrom;
352   vfrom.v4f[0] = vec_ld2f(&from[0]);
353   vfrom.v4f[1] = vec_ld2f(&from[2]);
354   return vfrom;
355 }
356 
357 template<> EIGEN_STRONG_INLINE Packet2d pload<Packet2d>(const double* from)
358 {
359   // FIXME: No intrinsic yet
360   EIGEN_DEBUG_ALIGNED_LOAD
361   Packet *vfrom;
362   vfrom = (Packet *) from;
363   return vfrom->v2d;
364 }
365 
366 template<> EIGEN_STRONG_INLINE void pstore<int>(int*       to, const Packet4i& from)
367 {
368   // FIXME: No intrinsic yet
369   EIGEN_DEBUG_ALIGNED_STORE
370   Packet *vto;
371   vto = (Packet *) to;
372   vto->v4i = from;
373 }
374 
375 template<> EIGEN_STRONG_INLINE void pstore<float>(float*   to, const Packet4f& from)
376 {
377   // FIXME: No intrinsic yet
378   EIGEN_DEBUG_ALIGNED_STORE
379   vec_st2f(from.v4f[0], &to[0]);
380   vec_st2f(from.v4f[1], &to[2]);
381 }
382 
383 
384 template<> EIGEN_STRONG_INLINE void pstore<double>(double*   to, const Packet2d& from)
385 {
386   // FIXME: No intrinsic yet
387   EIGEN_DEBUG_ALIGNED_STORE
388   Packet *vto;
389   vto = (Packet *) to;
390   vto->v2d = from;
391 }
392 
393 template<> EIGEN_STRONG_INLINE Packet4i pset1<Packet4i>(const int&    from)
394 {
395   return vec_splats(from);
396 }
397 template<> EIGEN_STRONG_INLINE Packet2d pset1<Packet2d>(const double& from) {
398   return vec_splats(from);
399 }
400 template<> EIGEN_STRONG_INLINE Packet4f pset1<Packet4f>(const float&    from)
401 {
402   Packet4f to;
403   to.v4f[0] = pset1<Packet2d>(static_cast<const double&>(from));
404   to.v4f[1] = to.v4f[0];
405   return to;
406 }
407 
408 template<> EIGEN_STRONG_INLINE void
409 pbroadcast4<Packet4i>(const int *a,
410                       Packet4i& a0, Packet4i& a1, Packet4i& a2, Packet4i& a3)
411 {
412   a3 = pload<Packet4i>(a);
413   a0 = vec_splat(a3, 0);
414   a1 = vec_splat(a3, 1);
415   a2 = vec_splat(a3, 2);
416   a3 = vec_splat(a3, 3);
417 }
418 
419 template<> EIGEN_STRONG_INLINE void
420 pbroadcast4<Packet4f>(const float *a,
421                       Packet4f& a0, Packet4f& a1, Packet4f& a2, Packet4f& a3)
422 {
423   a3 = pload<Packet4f>(a);
424   a0 = vec_splat_packet4f<0>(a3);
425   a1 = vec_splat_packet4f<1>(a3);
426   a2 = vec_splat_packet4f<2>(a3);
427   a3 = vec_splat_packet4f<3>(a3);
428 }
429 
430 template<> EIGEN_STRONG_INLINE void
431 pbroadcast4<Packet2d>(const double *a,
432                       Packet2d& a0, Packet2d& a1, Packet2d& a2, Packet2d& a3)
433 {
434   a1 = pload<Packet2d>(a);
435   a0 = vec_splat(a1, 0);
436   a1 = vec_splat(a1, 1);
437   a3 = pload<Packet2d>(a+2);
438   a2 = vec_splat(a3, 0);
439   a3 = vec_splat(a3, 1);
440 }
441 
442 template<> EIGEN_DEVICE_FUNC inline Packet4i pgather<int, Packet4i>(const int* from, Index stride)
443 {
444   int EIGEN_ALIGN16 ai[4];
445   ai[0] = from[0*stride];
446   ai[1] = from[1*stride];
447   ai[2] = from[2*stride];
448   ai[3] = from[3*stride];
449  return pload<Packet4i>(ai);
450 }
451 
452 template<> EIGEN_DEVICE_FUNC inline Packet4f pgather<float, Packet4f>(const float* from, Index stride)
453 {
454   float EIGEN_ALIGN16 ai[4];
455   ai[0] = from[0*stride];
456   ai[1] = from[1*stride];
457   ai[2] = from[2*stride];
458   ai[3] = from[3*stride];
459  return pload<Packet4f>(ai);
460 }
461 
462 template<> EIGEN_DEVICE_FUNC inline Packet2d pgather<double, Packet2d>(const double* from, Index stride)
463 {
464   double EIGEN_ALIGN16 af[2];
465   af[0] = from[0*stride];
466   af[1] = from[1*stride];
467  return pload<Packet2d>(af);
468 }
469 
470 template<> EIGEN_DEVICE_FUNC inline void pscatter<int, Packet4i>(int* to, const Packet4i& from, Index stride)
471 {
472   int EIGEN_ALIGN16 ai[4];
473   pstore<int>((int *)ai, from);
474   to[0*stride] = ai[0];
475   to[1*stride] = ai[1];
476   to[2*stride] = ai[2];
477   to[3*stride] = ai[3];
478 }
479 
480 template<> EIGEN_DEVICE_FUNC inline void pscatter<float, Packet4f>(float* to, const Packet4f& from, Index stride)
481 {
482   float EIGEN_ALIGN16 ai[4];
483   pstore<float>((float *)ai, from);
484   to[0*stride] = ai[0];
485   to[1*stride] = ai[1];
486   to[2*stride] = ai[2];
487   to[3*stride] = ai[3];
488 }
489 
490 template<> EIGEN_DEVICE_FUNC inline void pscatter<double, Packet2d>(double* to, const Packet2d& from, Index stride)
491 {
492   double EIGEN_ALIGN16 af[2];
493   pstore<double>(af, from);
494   to[0*stride] = af[0];
495   to[1*stride] = af[1];
496 }
497 
498 template<> EIGEN_STRONG_INLINE Packet4i padd<Packet4i>(const Packet4i& a, const Packet4i& b) { return (a + b); }
499 template<> EIGEN_STRONG_INLINE Packet4f padd<Packet4f>(const Packet4f& a, const Packet4f& b)
500 {
501   Packet4f c;
502   c.v4f[0] = a.v4f[0] + b.v4f[0];
503   c.v4f[1] = a.v4f[1] + b.v4f[1];
504   return c;
505 }
506 template<> EIGEN_STRONG_INLINE Packet2d padd<Packet2d>(const Packet2d& a, const Packet2d& b) { return (a + b); }
507 
508 template<> EIGEN_STRONG_INLINE Packet4i psub<Packet4i>(const Packet4i& a, const Packet4i& b) { return (a - b); }
509 template<> EIGEN_STRONG_INLINE Packet4f psub<Packet4f>(const Packet4f& a, const Packet4f& b)
510 {
511   Packet4f c;
512   c.v4f[0] = a.v4f[0] - b.v4f[0];
513   c.v4f[1] = a.v4f[1] - b.v4f[1];
514   return c;
515 }
516 template<> EIGEN_STRONG_INLINE Packet2d psub<Packet2d>(const Packet2d& a, const Packet2d& b) { return (a - b); }
517 
518 template<> EIGEN_STRONG_INLINE Packet4i pmul<Packet4i>(const Packet4i& a, const Packet4i& b) { return (a * b); }
519 template<> EIGEN_STRONG_INLINE Packet4f pmul<Packet4f>(const Packet4f& a, const Packet4f& b)
520 {
521   Packet4f c;
522   c.v4f[0] = a.v4f[0] * b.v4f[0];
523   c.v4f[1] = a.v4f[1] * b.v4f[1];
524   return c;
525 }
526 template<> EIGEN_STRONG_INLINE Packet2d pmul<Packet2d>(const Packet2d& a, const Packet2d& b) { return (a * b); }
527 
528 template<> EIGEN_STRONG_INLINE Packet4i pdiv<Packet4i>(const Packet4i& a, const Packet4i& b) { return (a / b); }
529 template<> EIGEN_STRONG_INLINE Packet4f pdiv<Packet4f>(const Packet4f& a, const Packet4f& b)
530 {
531   Packet4f c;
532   c.v4f[0] = a.v4f[0] / b.v4f[0];
533   c.v4f[1] = a.v4f[1] / b.v4f[1];
534   return c;
535 }
536 template<> EIGEN_STRONG_INLINE Packet2d pdiv<Packet2d>(const Packet2d& a, const Packet2d& b) { return (a / b); }
537 
538 template<> EIGEN_STRONG_INLINE Packet4i pnegate(const Packet4i& a) { return (-a); }
539 template<> EIGEN_STRONG_INLINE Packet4f pnegate(const Packet4f& a)
540 {
541   Packet4f c;
542   c.v4f[0] = -a.v4f[0];
543   c.v4f[1] = -a.v4f[1];
544   return c;
545 }
546 template<> EIGEN_STRONG_INLINE Packet2d pnegate(const Packet2d& a) { return (-a); }
547 
548 template<> EIGEN_STRONG_INLINE Packet4i pconj(const Packet4i& a) { return a; }
549 template<> EIGEN_STRONG_INLINE Packet4f pconj(const Packet4f& a) { return a; }
550 template<> EIGEN_STRONG_INLINE Packet2d pconj(const Packet2d& a) { return a; }
551 
552 template<> EIGEN_STRONG_INLINE Packet4i pmadd(const Packet4i& a, const Packet4i& b, const Packet4i& c) { return padd<Packet4i>(pmul<Packet4i>(a, b), c); }
553 template<> EIGEN_STRONG_INLINE Packet4f pmadd(const Packet4f& a, const Packet4f& b, const Packet4f& c)
554 {
555   Packet4f res;
556   res.v4f[0] = vec_madd(a.v4f[0], b.v4f[0], c.v4f[0]);
557   res.v4f[1] = vec_madd(a.v4f[1], b.v4f[1], c.v4f[1]);
558   return res;
559 }
560 template<> EIGEN_STRONG_INLINE Packet2d pmadd(const Packet2d& a, const Packet2d& b, const Packet2d& c) { return vec_madd(a, b, c); }
561 
562 template<> EIGEN_STRONG_INLINE Packet4i plset<Packet4i>(const int& a)    { return padd<Packet4i>(pset1<Packet4i>(a), p4i_COUNTDOWN); }
563 template<> EIGEN_STRONG_INLINE Packet4f plset<Packet4f>(const float& a)  { return padd<Packet4f>(pset1<Packet4f>(a), p4f_COUNTDOWN); }
564 template<> EIGEN_STRONG_INLINE Packet2d plset<Packet2d>(const double& a) { return padd<Packet2d>(pset1<Packet2d>(a), p2d_COUNTDOWN); }
565 
566 template<> EIGEN_STRONG_INLINE Packet4i pmin<Packet4i>(const Packet4i& a, const Packet4i& b) { return vec_min(a, b); }
567 template<> EIGEN_STRONG_INLINE Packet2d pmin<Packet2d>(const Packet2d& a, const Packet2d& b) { return vec_min(a, b); }
568 template<> EIGEN_STRONG_INLINE Packet4f pmin<Packet4f>(const Packet4f& a, const Packet4f& b)
569 {
570   Packet4f res;
571   res.v4f[0] = pmin(a.v4f[0], b.v4f[0]);
572   res.v4f[1] = pmin(a.v4f[1], b.v4f[1]);
573   return res;
574 }
575 
576 template<> EIGEN_STRONG_INLINE Packet4i pmax<Packet4i>(const Packet4i& a, const Packet4i& b) { return vec_max(a, b); }
577 template<> EIGEN_STRONG_INLINE Packet2d pmax<Packet2d>(const Packet2d& a, const Packet2d& b) { return vec_max(a, b); }
578 template<> EIGEN_STRONG_INLINE Packet4f pmax<Packet4f>(const Packet4f& a, const Packet4f& b)
579 {
580   Packet4f res;
581   res.v4f[0] = pmax(a.v4f[0], b.v4f[0]);
582   res.v4f[1] = pmax(a.v4f[1], b.v4f[1]);
583   return res;
584 }
585 
586 template<> EIGEN_STRONG_INLINE Packet4i pand<Packet4i>(const Packet4i& a, const Packet4i& b) { return vec_and(a, b); }
587 template<> EIGEN_STRONG_INLINE Packet2d pand<Packet2d>(const Packet2d& a, const Packet2d& b) { return vec_and(a, b); }
588 template<> EIGEN_STRONG_INLINE Packet4f pand<Packet4f>(const Packet4f& a, const Packet4f& b)
589 {
590   Packet4f res;
591   res.v4f[0] = pand(a.v4f[0], b.v4f[0]);
592   res.v4f[1] = pand(a.v4f[1], b.v4f[1]);
593   return res;
594 }
595 
596 template<> EIGEN_STRONG_INLINE Packet4i por<Packet4i>(const Packet4i& a, const Packet4i& b) { return vec_or(a, b); }
597 template<> EIGEN_STRONG_INLINE Packet2d por<Packet2d>(const Packet2d& a, const Packet2d& b) { return vec_or(a, b); }
598 template<> EIGEN_STRONG_INLINE Packet4f por<Packet4f>(const Packet4f& a, const Packet4f& b)
599 {
600   Packet4f res;
601   res.v4f[0] = pand(a.v4f[0], b.v4f[0]);
602   res.v4f[1] = pand(a.v4f[1], b.v4f[1]);
603   return res;
604 }
605 
606 template<> EIGEN_STRONG_INLINE Packet4i pxor<Packet4i>(const Packet4i& a, const Packet4i& b) { return vec_xor(a, b); }
607 template<> EIGEN_STRONG_INLINE Packet2d pxor<Packet2d>(const Packet2d& a, const Packet2d& b) { return vec_xor(a, b); }
608 template<> EIGEN_STRONG_INLINE Packet4f pxor<Packet4f>(const Packet4f& a, const Packet4f& b)
609 {
610   Packet4f res;
611   res.v4f[0] = pand(a.v4f[0], b.v4f[0]);
612   res.v4f[1] = pand(a.v4f[1], b.v4f[1]);
613   return res;
614 }
615 
616 template<> EIGEN_STRONG_INLINE Packet4i pandnot<Packet4i>(const Packet4i& a, const Packet4i& b) { return pand<Packet4i>(a, vec_nor(b, b)); }
617 template<> EIGEN_STRONG_INLINE Packet2d pandnot<Packet2d>(const Packet2d& a, const Packet2d& b) { return vec_and(a, vec_nor(b, b)); }
618 template<> EIGEN_STRONG_INLINE Packet4f pandnot<Packet4f>(const Packet4f& a, const Packet4f& b)
619 {
620   Packet4f res;
621   res.v4f[0] = pandnot(a.v4f[0], b.v4f[0]);
622   res.v4f[1] = pandnot(a.v4f[1], b.v4f[1]);
623   return res;
624 }
625 
626 template<> EIGEN_STRONG_INLINE Packet4f pround<Packet4f>(const Packet4f& a)
627 {
628   Packet4f res;
629   res.v4f[0] = vec_round(a.v4f[0]);
630   res.v4f[1] = vec_round(a.v4f[1]);
631   return res;
632 }
633 template<> EIGEN_STRONG_INLINE Packet2d pround<Packet2d>(const Packet2d& a) { return vec_round(a); }
634 template<> EIGEN_STRONG_INLINE Packet4f pceil<Packet4f>(const  Packet4f& a)
635 {
636   Packet4f res;
637   res.v4f[0] = vec_ceil(a.v4f[0]);
638   res.v4f[1] = vec_ceil(a.v4f[1]);
639   return res;
640 }
641 template<> EIGEN_STRONG_INLINE Packet2d pceil<Packet2d>(const  Packet2d& a) { return vec_ceil(a); }
642 template<> EIGEN_STRONG_INLINE Packet4f pfloor<Packet4f>(const Packet4f& a)
643 {
644   Packet4f res;
645   res.v4f[0] = vec_floor(a.v4f[0]);
646   res.v4f[1] = vec_floor(a.v4f[1]);
647   return res;
648 }
649 template<> EIGEN_STRONG_INLINE Packet2d pfloor<Packet2d>(const Packet2d& a) { return vec_floor(a); }
650 
651 template<> EIGEN_STRONG_INLINE Packet4i ploadu<Packet4i>(const int*       from) { return pload<Packet4i>(from); }
652 template<> EIGEN_STRONG_INLINE Packet4f ploadu<Packet4f>(const float*     from) { return pload<Packet4f>(from); }
653 template<> EIGEN_STRONG_INLINE Packet2d ploadu<Packet2d>(const double*    from) { return pload<Packet2d>(from); }
654 
655 
656 template<> EIGEN_STRONG_INLINE Packet4i ploaddup<Packet4i>(const int*     from)
657 {
658   Packet4i p = pload<Packet4i>(from);
659   return vec_perm(p, p, p16uc_DUPLICATE32_HI);
660 }
661 
662 template<> EIGEN_STRONG_INLINE Packet4f ploaddup<Packet4f>(const float*    from)
663 {
664   Packet4f p = pload<Packet4f>(from);
665   p.v4f[1] = vec_splat(p.v4f[0], 1);
666   p.v4f[0] = vec_splat(p.v4f[0], 0);
667   return p;
668 }
669 
670 template<> EIGEN_STRONG_INLINE Packet2d ploaddup<Packet2d>(const double*   from)
671 {
672   Packet2d p = pload<Packet2d>(from);
673   return vec_perm(p, p, p16uc_PSET64_HI);
674 }
675 
676 template<> EIGEN_STRONG_INLINE void pstoreu<int>(int*        to, const Packet4i& from) { pstore<int>(to, from); }
677 template<> EIGEN_STRONG_INLINE void pstoreu<float>(float*    to, const Packet4f& from) { pstore<float>(to, from); }
678 template<> EIGEN_STRONG_INLINE void pstoreu<double>(double*  to, const Packet2d& from) { pstore<double>(to, from); }
679 
680 template<> EIGEN_STRONG_INLINE void prefetch<int>(const int*       addr) { EIGEN_ZVECTOR_PREFETCH(addr); }
681 template<> EIGEN_STRONG_INLINE void prefetch<float>(const float*   addr) { EIGEN_ZVECTOR_PREFETCH(addr); }
682 template<> EIGEN_STRONG_INLINE void prefetch<double>(const double* addr) { EIGEN_ZVECTOR_PREFETCH(addr); }
683 
684 template<> EIGEN_STRONG_INLINE int    pfirst<Packet4i>(const Packet4i& a) { int    EIGEN_ALIGN16 x[4]; pstore(x, a); return x[0]; }
685 template<> EIGEN_STRONG_INLINE float  pfirst<Packet4f>(const Packet4f& a) { float  EIGEN_ALIGN16 x[2]; vec_st2f(a.v4f[0], &x[0]); return x[0]; }
686 template<> EIGEN_STRONG_INLINE double pfirst<Packet2d>(const Packet2d& a) { double EIGEN_ALIGN16 x[2]; pstore(x, a); return x[0]; }
687 
688 template<> EIGEN_STRONG_INLINE Packet4i preverse(const Packet4i& a)
689 {
690   return reinterpret_cast<Packet4i>(vec_perm(reinterpret_cast<Packet16uc>(a), reinterpret_cast<Packet16uc>(a), p16uc_REVERSE32));
691 }
692 
693 template<> EIGEN_STRONG_INLINE Packet2d preverse(const Packet2d& a)
694 {
695   return reinterpret_cast<Packet2d>(vec_perm(reinterpret_cast<Packet16uc>(a), reinterpret_cast<Packet16uc>(a), p16uc_REVERSE64));
696 }
697 
698 template<> EIGEN_STRONG_INLINE Packet4f preverse(const Packet4f& a)
699 {
700   Packet4f rev;
701   rev.v4f[0] = preverse<Packet2d>(a.v4f[1]);
702   rev.v4f[1] = preverse<Packet2d>(a.v4f[0]);
703   return rev;
704 }
705 
706 template<> EIGEN_STRONG_INLINE Packet4i pabs<Packet4i>(const Packet4i& a) { return vec_abs(a); }
707 template<> EIGEN_STRONG_INLINE Packet2d pabs<Packet2d>(const Packet2d& a) { return vec_abs(a); }
708 template<> EIGEN_STRONG_INLINE Packet4f pabs<Packet4f>(const Packet4f& a)
709 {
710   Packet4f res;
711   res.v4f[0] = pabs(a.v4f[0]);
712   res.v4f[1] = pabs(a.v4f[1]);
713   return res;
714 }
715 
716 template<> EIGEN_STRONG_INLINE int predux<Packet4i>(const Packet4i& a)
717 {
718   Packet4i b, sum;
719   b   = vec_sld(a, a, 8);
720   sum = padd<Packet4i>(a, b);
721   b   = vec_sld(sum, sum, 4);
722   sum = padd<Packet4i>(sum, b);
723   return pfirst(sum);
724 }
725 
726 template<> EIGEN_STRONG_INLINE double predux<Packet2d>(const Packet2d& a)
727 {
728   Packet2d b, sum;
729   b   = reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4i>(a), reinterpret_cast<Packet4i>(a), 8));
730   sum = padd<Packet2d>(a, b);
731   return pfirst(sum);
732 }
733 template<> EIGEN_STRONG_INLINE float predux<Packet4f>(const Packet4f& a)
734 {
735   Packet2d sum;
736   sum = padd<Packet2d>(a.v4f[0], a.v4f[1]);
737   double first = predux<Packet2d>(sum);
738   return static_cast<float>(first);
739 }
740 
741 template<> EIGEN_STRONG_INLINE Packet4i preduxp<Packet4i>(const Packet4i* vecs)
742 {
743   Packet4i v[4], sum[4];
744 
745   // It's easier and faster to transpose then add as columns
746   // Check: http://www.freevec.org/function/matrix_4x4_transpose_floats for explanation
747   // Do the transpose, first set of moves
748   v[0] = vec_mergeh(vecs[0], vecs[2]);
749   v[1] = vec_mergel(vecs[0], vecs[2]);
750   v[2] = vec_mergeh(vecs[1], vecs[3]);
751   v[3] = vec_mergel(vecs[1], vecs[3]);
752   // Get the resulting vectors
753   sum[0] = vec_mergeh(v[0], v[2]);
754   sum[1] = vec_mergel(v[0], v[2]);
755   sum[2] = vec_mergeh(v[1], v[3]);
756   sum[3] = vec_mergel(v[1], v[3]);
757 
758   // Now do the summation:
759   // Lines 0+1
760   sum[0] = padd<Packet4i>(sum[0], sum[1]);
761   // Lines 2+3
762   sum[1] = padd<Packet4i>(sum[2], sum[3]);
763   // Add the results
764   sum[0] = padd<Packet4i>(sum[0], sum[1]);
765 
766   return sum[0];
767 }
768 
769 template<> EIGEN_STRONG_INLINE Packet2d preduxp<Packet2d>(const Packet2d* vecs)
770 {
771   Packet2d v[2], sum;
772   v[0] = padd<Packet2d>(vecs[0], reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4ui>(vecs[0]), reinterpret_cast<Packet4ui>(vecs[0]), 8)));
773   v[1] = padd<Packet2d>(vecs[1], reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4ui>(vecs[1]), reinterpret_cast<Packet4ui>(vecs[1]), 8)));
774 
775   sum = reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4ui>(v[0]), reinterpret_cast<Packet4ui>(v[1]), 8));
776 
777   return sum;
778 }
779 
780 template<> EIGEN_STRONG_INLINE Packet4f preduxp<Packet4f>(const Packet4f* vecs)
781 {
782   PacketBlock<Packet4f,4> transpose;
783   transpose.packet[0] = vecs[0];
784   transpose.packet[1] = vecs[1];
785   transpose.packet[2] = vecs[2];
786   transpose.packet[3] = vecs[3];
787   ptranspose(transpose);
788 
789   Packet4f sum = padd(transpose.packet[0], transpose.packet[1]);
790   sum = padd(sum, transpose.packet[2]);
791   sum = padd(sum, transpose.packet[3]);
792   return sum;
793 }
794 
795 // Other reduction functions:
796 // mul
797 template<> EIGEN_STRONG_INLINE int predux_mul<Packet4i>(const Packet4i& a)
798 {
799   EIGEN_ALIGN16 int aux[4];
800   pstore(aux, a);
801   return aux[0] * aux[1] * aux[2] * aux[3];
802 }
803 
804 template<> EIGEN_STRONG_INLINE double predux_mul<Packet2d>(const Packet2d& a)
805 {
806   return pfirst(pmul(a, reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4i>(a), reinterpret_cast<Packet4i>(a), 8))));
807 }
808 
809 template<> EIGEN_STRONG_INLINE float predux_mul<Packet4f>(const Packet4f& a)
810 {
811   // Return predux_mul<Packet2d> of the subvectors product
812   return static_cast<float>(pfirst(predux_mul(pmul(a.v4f[0], a.v4f[1]))));
813 }
814 
815 // min
816 template<> EIGEN_STRONG_INLINE int predux_min<Packet4i>(const Packet4i& a)
817 {
818   Packet4i b, res;
819   b   = pmin<Packet4i>(a, vec_sld(a, a, 8));
820   res = pmin<Packet4i>(b, vec_sld(b, b, 4));
821   return pfirst(res);
822 }
823 
824 template<> EIGEN_STRONG_INLINE double predux_min<Packet2d>(const Packet2d& a)
825 {
826   return pfirst(pmin<Packet2d>(a, reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4i>(a), reinterpret_cast<Packet4i>(a), 8))));
827 }
828 
829 template<> EIGEN_STRONG_INLINE float predux_min<Packet4f>(const Packet4f& a)
830 {
831   Packet2d b, res;
832   b   = pmin<Packet2d>(a.v4f[0], a.v4f[1]);
833   res = pmin<Packet2d>(b, reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4i>(b), reinterpret_cast<Packet4i>(b), 8)));
834   return static_cast<float>(pfirst(res));
835 }
836 
837 // max
838 template<> EIGEN_STRONG_INLINE int predux_max<Packet4i>(const Packet4i& a)
839 {
840   Packet4i b, res;
841   b = pmax<Packet4i>(a, vec_sld(a, a, 8));
842   res = pmax<Packet4i>(b, vec_sld(b, b, 4));
843   return pfirst(res);
844 }
845 
846 // max
847 template<> EIGEN_STRONG_INLINE double predux_max<Packet2d>(const Packet2d& a)
848 {
849   return pfirst(pmax<Packet2d>(a, reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4i>(a), reinterpret_cast<Packet4i>(a), 8))));
850 }
851 
852 template<> EIGEN_STRONG_INLINE float predux_max<Packet4f>(const Packet4f& a)
853 {
854   Packet2d b, res;
855   b   = pmax<Packet2d>(a.v4f[0], a.v4f[1]);
856   res = pmax<Packet2d>(b, reinterpret_cast<Packet2d>(vec_sld(reinterpret_cast<Packet4i>(b), reinterpret_cast<Packet4i>(b), 8)));
857   return static_cast<float>(pfirst(res));
858 }
859 
860 EIGEN_DEVICE_FUNC inline void
861 ptranspose(PacketBlock<Packet4i,4>& kernel) {
862   Packet4i t0 = vec_mergeh(kernel.packet[0], kernel.packet[2]);
863   Packet4i t1 = vec_mergel(kernel.packet[0], kernel.packet[2]);
864   Packet4i t2 = vec_mergeh(kernel.packet[1], kernel.packet[3]);
865   Packet4i t3 = vec_mergel(kernel.packet[1], kernel.packet[3]);
866   kernel.packet[0] = vec_mergeh(t0, t2);
867   kernel.packet[1] = vec_mergel(t0, t2);
868   kernel.packet[2] = vec_mergeh(t1, t3);
869   kernel.packet[3] = vec_mergel(t1, t3);
870 }
871 
872 EIGEN_DEVICE_FUNC inline void
873 ptranspose(PacketBlock<Packet2d,2>& kernel) {
874   Packet2d t0 = vec_perm(kernel.packet[0], kernel.packet[1], p16uc_TRANSPOSE64_HI);
875   Packet2d t1 = vec_perm(kernel.packet[0], kernel.packet[1], p16uc_TRANSPOSE64_LO);
876   kernel.packet[0] = t0;
877   kernel.packet[1] = t1;
878 }
879 
880 /* Split the Packet4f PacketBlock into 4 Packet2d PacketBlocks and transpose each one
881  */
882 EIGEN_DEVICE_FUNC inline void
883 ptranspose(PacketBlock<Packet4f,4>& kernel) {
884   PacketBlock<Packet2d,2> t0,t1,t2,t3;
885   // copy top-left 2x2 Packet2d block
886   t0.packet[0] = kernel.packet[0].v4f[0];
887   t0.packet[1] = kernel.packet[1].v4f[0];
888 
889   // copy top-right 2x2 Packet2d block
890   t1.packet[0] = kernel.packet[0].v4f[1];
891   t1.packet[1] = kernel.packet[1].v4f[1];
892 
893   // copy bottom-left 2x2 Packet2d block
894   t2.packet[0] = kernel.packet[2].v4f[0];
895   t2.packet[1] = kernel.packet[3].v4f[0];
896 
897   // copy bottom-right 2x2 Packet2d block
898   t3.packet[0] = kernel.packet[2].v4f[1];
899   t3.packet[1] = kernel.packet[3].v4f[1];
900 
901   // Transpose all 2x2 blocks
902   ptranspose(t0);
903   ptranspose(t1);
904   ptranspose(t2);
905   ptranspose(t3);
906 
907   // Copy back transposed blocks, but exchange t1 and t2 due to transposition
908   kernel.packet[0].v4f[0] = t0.packet[0];
909   kernel.packet[0].v4f[1] = t2.packet[0];
910   kernel.packet[1].v4f[0] = t0.packet[1];
911   kernel.packet[1].v4f[1] = t2.packet[1];
912   kernel.packet[2].v4f[0] = t1.packet[0];
913   kernel.packet[2].v4f[1] = t3.packet[0];
914   kernel.packet[3].v4f[0] = t1.packet[1];
915   kernel.packet[3].v4f[1] = t3.packet[1];
916 }
917 
918 template<> EIGEN_STRONG_INLINE Packet4i pblend(const Selector<4>& ifPacket, const Packet4i& thenPacket, const Packet4i& elsePacket) {
919   Packet4ui select = { ifPacket.select[0], ifPacket.select[1], ifPacket.select[2], ifPacket.select[3] };
920   Packet4ui mask = vec_cmpeq(select, reinterpret_cast<Packet4ui>(p4i_ONE));
921   return vec_sel(elsePacket, thenPacket, mask);
922 }
923 
924 template<> EIGEN_STRONG_INLINE Packet4f pblend(const Selector<4>& ifPacket, const Packet4f& thenPacket, const Packet4f& elsePacket) {
925   Packet2ul select_hi = { ifPacket.select[0], ifPacket.select[1] };
926   Packet2ul select_lo = { ifPacket.select[2], ifPacket.select[3] };
927   Packet2ul mask_hi = vec_cmpeq(select_hi, reinterpret_cast<Packet2ul>(p2l_ONE));
928   Packet2ul mask_lo = vec_cmpeq(select_lo, reinterpret_cast<Packet2ul>(p2l_ONE));
929   Packet4f result;
930   result.v4f[0] = vec_sel(elsePacket.v4f[0], thenPacket.v4f[0], mask_hi);
931   result.v4f[1] = vec_sel(elsePacket.v4f[1], thenPacket.v4f[1], mask_lo);
932   return result;
933 }
934 
935 template<> EIGEN_STRONG_INLINE Packet2d pblend(const Selector<2>& ifPacket, const Packet2d& thenPacket, const Packet2d& elsePacket) {
936   Packet2ul select = { ifPacket.select[0], ifPacket.select[1] };
937   Packet2ul mask = vec_cmpeq(select, reinterpret_cast<Packet2ul>(p2l_ONE));
938   return vec_sel(elsePacket, thenPacket, mask);
939 }
940 
941 } // end namespace internal
942 
943 } // end namespace Eigen
944 
945 #endif // EIGEN_PACKET_MATH_ZVECTOR_H
946