cart-elc

Source code for CART-ELC
git clone git://git.laack.co/cart-elc.git
Log | Files | Refs | README | LICENSE

Complex.h (17160B)


      1 // This file is part of Eigen, a lightweight C++ template library
      2 // for linear algebra.
      3 //
      4 // Copyright (C) 2018 Gael Guennebaud <gael.guennebaud@inria.fr>
      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_COMPLEX_AVX512_H
     11 #define EIGEN_COMPLEX_AVX512_H
     12 
     13 namespace Eigen {
     14 
     15 namespace internal {
     16 
     17 //---------- float ----------
     18 struct Packet8cf
     19 {
     20   EIGEN_STRONG_INLINE Packet8cf() {}
     21   EIGEN_STRONG_INLINE explicit Packet8cf(const __m512& a) : v(a) {}
     22   __m512  v;
     23 };
     24 
     25 template<> struct packet_traits<std::complex<float> >  : default_packet_traits
     26 {
     27   typedef Packet8cf type;
     28   typedef Packet4cf half;
     29   enum {
     30     Vectorizable = 1,
     31     AlignedOnScalar = 1,
     32     size = 8,
     33     HasHalfPacket = 1,
     34 
     35     HasAdd    = 1,
     36     HasSub    = 1,
     37     HasMul    = 1,
     38     HasDiv    = 1,
     39     HasNegate = 1,
     40     HasSqrt   = 1,
     41     HasAbs    = 0,
     42     HasAbs2   = 0,
     43     HasMin    = 0,
     44     HasMax    = 0,
     45     HasSetLinear = 0
     46   };
     47 };
     48 
     49 template<> struct unpacket_traits<Packet8cf> {
     50   typedef std::complex<float> type;
     51   typedef Packet4cf half;
     52   typedef Packet16f as_real;
     53   enum {
     54     size = 8,
     55     alignment=unpacket_traits<Packet16f>::alignment,
     56     vectorizable=true,
     57     masked_load_available=false,
     58     masked_store_available=false
     59   };
     60 };
     61 
     62 template<> EIGEN_STRONG_INLINE Packet8cf ptrue<Packet8cf>(const Packet8cf& a) { return Packet8cf(ptrue(Packet16f(a.v))); }
     63 template<> EIGEN_STRONG_INLINE Packet8cf padd<Packet8cf>(const Packet8cf& a, const Packet8cf& b) { return Packet8cf(_mm512_add_ps(a.v,b.v)); }
     64 template<> EIGEN_STRONG_INLINE Packet8cf psub<Packet8cf>(const Packet8cf& a, const Packet8cf& b) { return Packet8cf(_mm512_sub_ps(a.v,b.v)); }
     65 template<> EIGEN_STRONG_INLINE Packet8cf pnegate(const Packet8cf& a)
     66 {
     67   return Packet8cf(pnegate(a.v));
     68 }
     69 template<> EIGEN_STRONG_INLINE Packet8cf pconj(const Packet8cf& a)
     70 {
     71   const __m512 mask = _mm512_castsi512_ps(_mm512_setr_epi32(
     72     0x00000000,0x80000000,0x00000000,0x80000000,0x00000000,0x80000000,0x00000000,0x80000000,
     73     0x00000000,0x80000000,0x00000000,0x80000000,0x00000000,0x80000000,0x00000000,0x80000000));
     74   return Packet8cf(pxor(a.v,mask));
     75 }
     76 
     77 template<> EIGEN_STRONG_INLINE Packet8cf pmul<Packet8cf>(const Packet8cf& a, const Packet8cf& b)
     78 {
     79   __m512 tmp2 = _mm512_mul_ps(_mm512_movehdup_ps(a.v), _mm512_permute_ps(b.v, _MM_SHUFFLE(2,3,0,1)));
     80   return Packet8cf(_mm512_fmaddsub_ps(_mm512_moveldup_ps(a.v), b.v, tmp2));
     81 }
     82 
     83 template<> EIGEN_STRONG_INLINE Packet8cf pand   <Packet8cf>(const Packet8cf& a, const Packet8cf& b) { return Packet8cf(pand(a.v,b.v)); }
     84 template<> EIGEN_STRONG_INLINE Packet8cf por    <Packet8cf>(const Packet8cf& a, const Packet8cf& b) { return Packet8cf(por(a.v,b.v)); }
     85 template<> EIGEN_STRONG_INLINE Packet8cf pxor   <Packet8cf>(const Packet8cf& a, const Packet8cf& b) { return Packet8cf(pxor(a.v,b.v)); }
     86 template<> EIGEN_STRONG_INLINE Packet8cf pandnot<Packet8cf>(const Packet8cf& a, const Packet8cf& b) { return Packet8cf(pandnot(a.v,b.v)); }
     87 
     88 template <>
     89 EIGEN_STRONG_INLINE Packet8cf pcmp_eq(const Packet8cf& a, const Packet8cf& b) {
     90   __m512 eq = pcmp_eq<Packet16f>(a.v, b.v);
     91   return Packet8cf(pand(eq, _mm512_permute_ps(eq, 0xB1)));
     92 }
     93 
     94 template<> EIGEN_STRONG_INLINE Packet8cf pload <Packet8cf>(const std::complex<float>* from) { EIGEN_DEBUG_ALIGNED_LOAD return Packet8cf(pload<Packet16f>(&numext::real_ref(*from))); }
     95 template<> EIGEN_STRONG_INLINE Packet8cf ploadu<Packet8cf>(const std::complex<float>* from) { EIGEN_DEBUG_UNALIGNED_LOAD return Packet8cf(ploadu<Packet16f>(&numext::real_ref(*from))); }
     96 
     97 
     98 template<> EIGEN_STRONG_INLINE Packet8cf pset1<Packet8cf>(const std::complex<float>& from)
     99 {
    100   return Packet8cf(_mm512_castpd_ps(pload1<Packet8d>((const double*)(const void*)&from)));
    101 }
    102 
    103 template<> EIGEN_STRONG_INLINE Packet8cf ploaddup<Packet8cf>(const std::complex<float>* from)
    104 {
    105   return Packet8cf( _mm512_castpd_ps( ploaddup<Packet8d>((const double*)(const void*)from )) );
    106 }
    107 template<> EIGEN_STRONG_INLINE Packet8cf ploadquad<Packet8cf>(const std::complex<float>* from)
    108 {
    109   return Packet8cf( _mm512_castpd_ps( ploadquad<Packet8d>((const double*)(const void*)from )) );
    110 }
    111 
    112 template<> EIGEN_STRONG_INLINE void pstore <std::complex<float> >(std::complex<float>* to, const Packet8cf& from) { EIGEN_DEBUG_ALIGNED_STORE pstore(&numext::real_ref(*to), from.v); }
    113 template<> EIGEN_STRONG_INLINE void pstoreu<std::complex<float> >(std::complex<float>* to, const Packet8cf& from) { EIGEN_DEBUG_UNALIGNED_STORE pstoreu(&numext::real_ref(*to), from.v); }
    114 
    115 template<> EIGEN_DEVICE_FUNC inline Packet8cf pgather<std::complex<float>, Packet8cf>(const std::complex<float>* from, Index stride)
    116 {
    117   return Packet8cf(_mm512_castpd_ps(pgather<double,Packet8d>((const double*)(const void*)from, stride)));
    118 }
    119 
    120 template<> EIGEN_DEVICE_FUNC inline void pscatter<std::complex<float>, Packet8cf>(std::complex<float>* to, const Packet8cf& from, Index stride)
    121 {
    122   pscatter((double*)(void*)to, _mm512_castps_pd(from.v), stride);
    123 }
    124 
    125 template<> EIGEN_STRONG_INLINE std::complex<float>  pfirst<Packet8cf>(const Packet8cf& a)
    126 {
    127   return pfirst(Packet2cf(_mm512_castps512_ps128(a.v)));
    128 }
    129 
    130 template<> EIGEN_STRONG_INLINE Packet8cf preverse(const Packet8cf& a) {
    131   return Packet8cf(_mm512_castsi512_ps(
    132             _mm512_permutexvar_epi64( _mm512_set_epi32(0, 0, 0, 1, 0, 2, 0, 3, 0, 4, 0, 5, 0, 6, 0, 7),
    133                                       _mm512_castps_si512(a.v))));
    134 }
    135 
    136 template<> EIGEN_STRONG_INLINE std::complex<float> predux<Packet8cf>(const Packet8cf& a)
    137 {
    138   return predux(padd(Packet4cf(extract256<0>(a.v)),
    139                      Packet4cf(extract256<1>(a.v))));
    140 }
    141 
    142 template<> EIGEN_STRONG_INLINE std::complex<float> predux_mul<Packet8cf>(const Packet8cf& a)
    143 {
    144   return predux_mul(pmul(Packet4cf(extract256<0>(a.v)),
    145                          Packet4cf(extract256<1>(a.v))));
    146 }
    147 
    148 template <>
    149 EIGEN_STRONG_INLINE Packet4cf predux_half_dowto4<Packet8cf>(const Packet8cf& a) {
    150   __m256 lane0 = extract256<0>(a.v);
    151   __m256 lane1 = extract256<1>(a.v);
    152   __m256 res = _mm256_add_ps(lane0, lane1);
    153   return Packet4cf(res);
    154 }
    155 
    156 EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet8cf,Packet16f)
    157 
    158 template<> EIGEN_STRONG_INLINE Packet8cf pdiv<Packet8cf>(const Packet8cf& a, const Packet8cf& b)
    159 {
    160   Packet8cf num = pmul(a, pconj(b));
    161   __m512 tmp = _mm512_mul_ps(b.v, b.v);
    162   __m512 tmp2    = _mm512_shuffle_ps(tmp,tmp,0xB1);
    163   __m512 denom = _mm512_add_ps(tmp, tmp2);
    164   return Packet8cf(_mm512_div_ps(num.v, denom));
    165 }
    166 
    167 template<> EIGEN_STRONG_INLINE Packet8cf pcplxflip<Packet8cf>(const Packet8cf& x)
    168 {
    169   return Packet8cf(_mm512_shuffle_ps(x.v, x.v, _MM_SHUFFLE(2, 3, 0 ,1)));
    170 }
    171 
    172 //---------- double ----------
    173 struct Packet4cd
    174 {
    175   EIGEN_STRONG_INLINE Packet4cd() {}
    176   EIGEN_STRONG_INLINE explicit Packet4cd(const __m512d& a) : v(a) {}
    177   __m512d  v;
    178 };
    179 
    180 template<> struct packet_traits<std::complex<double> >  : default_packet_traits
    181 {
    182   typedef Packet4cd type;
    183   typedef Packet2cd half;
    184   enum {
    185     Vectorizable = 1,
    186     AlignedOnScalar = 0,
    187     size = 4,
    188     HasHalfPacket = 1,
    189 
    190     HasAdd    = 1,
    191     HasSub    = 1,
    192     HasMul    = 1,
    193     HasDiv    = 1,
    194     HasNegate = 1,
    195     HasSqrt   = 1,
    196     HasAbs    = 0,
    197     HasAbs2   = 0,
    198     HasMin    = 0,
    199     HasMax    = 0,
    200     HasSetLinear = 0
    201   };
    202 };
    203 
    204 template<> struct unpacket_traits<Packet4cd> {
    205   typedef std::complex<double> type;
    206   typedef Packet2cd half;
    207   typedef Packet8d as_real;
    208   enum {
    209     size = 4,
    210     alignment = unpacket_traits<Packet8d>::alignment,
    211     vectorizable=true,
    212     masked_load_available=false,
    213     masked_store_available=false
    214   };
    215 };
    216 
    217 template<> EIGEN_STRONG_INLINE Packet4cd padd<Packet4cd>(const Packet4cd& a, const Packet4cd& b) { return Packet4cd(_mm512_add_pd(a.v,b.v)); }
    218 template<> EIGEN_STRONG_INLINE Packet4cd psub<Packet4cd>(const Packet4cd& a, const Packet4cd& b) { return Packet4cd(_mm512_sub_pd(a.v,b.v)); }
    219 template<> EIGEN_STRONG_INLINE Packet4cd pnegate(const Packet4cd& a) { return Packet4cd(pnegate(a.v)); }
    220 template<> EIGEN_STRONG_INLINE Packet4cd pconj(const Packet4cd& a)
    221 {
    222   const __m512d mask = _mm512_castsi512_pd(
    223           _mm512_set_epi32(0x80000000,0x0,0x0,0x0,0x80000000,0x0,0x0,0x0,
    224                            0x80000000,0x0,0x0,0x0,0x80000000,0x0,0x0,0x0));
    225   return Packet4cd(pxor(a.v,mask));
    226 }
    227 
    228 template<> EIGEN_STRONG_INLINE Packet4cd pmul<Packet4cd>(const Packet4cd& a, const Packet4cd& b)
    229 {
    230   __m512d tmp1 = _mm512_shuffle_pd(a.v,a.v,0x0);
    231   __m512d tmp2 = _mm512_shuffle_pd(a.v,a.v,0xFF);
    232   __m512d tmp3 = _mm512_shuffle_pd(b.v,b.v,0x55);
    233   __m512d odd  = _mm512_mul_pd(tmp2, tmp3);
    234   return Packet4cd(_mm512_fmaddsub_pd(tmp1, b.v, odd));
    235 }
    236 
    237 template<> EIGEN_STRONG_INLINE Packet4cd ptrue<Packet4cd>(const Packet4cd& a) { return Packet4cd(ptrue(Packet8d(a.v))); }
    238 template<> EIGEN_STRONG_INLINE Packet4cd pand   <Packet4cd>(const Packet4cd& a, const Packet4cd& b) { return Packet4cd(pand(a.v,b.v)); }
    239 template<> EIGEN_STRONG_INLINE Packet4cd por    <Packet4cd>(const Packet4cd& a, const Packet4cd& b) { return Packet4cd(por(a.v,b.v)); }
    240 template<> EIGEN_STRONG_INLINE Packet4cd pxor   <Packet4cd>(const Packet4cd& a, const Packet4cd& b) { return Packet4cd(pxor(a.v,b.v)); }
    241 template<> EIGEN_STRONG_INLINE Packet4cd pandnot<Packet4cd>(const Packet4cd& a, const Packet4cd& b) { return Packet4cd(pandnot(a.v,b.v)); }
    242 
    243 template <>
    244 EIGEN_STRONG_INLINE Packet4cd pcmp_eq(const Packet4cd& a, const Packet4cd& b) {
    245   __m512d eq = pcmp_eq<Packet8d>(a.v, b.v);
    246   return Packet4cd(pand(eq, _mm512_permute_pd(eq, 0x55)));
    247 }
    248 
    249 template<> EIGEN_STRONG_INLINE Packet4cd pload <Packet4cd>(const std::complex<double>* from)
    250 { EIGEN_DEBUG_ALIGNED_LOAD return Packet4cd(pload<Packet8d>((const double*)from)); }
    251 template<> EIGEN_STRONG_INLINE Packet4cd ploadu<Packet4cd>(const std::complex<double>* from)
    252 { EIGEN_DEBUG_UNALIGNED_LOAD return Packet4cd(ploadu<Packet8d>((const double*)from)); }
    253 
    254 template<> EIGEN_STRONG_INLINE Packet4cd pset1<Packet4cd>(const std::complex<double>& from)
    255 {
    256   #ifdef EIGEN_VECTORIZE_AVX512DQ
    257   return Packet4cd(_mm512_broadcast_f64x2(pset1<Packet1cd>(from).v));
    258   #else
    259   return Packet4cd(_mm512_castps_pd(_mm512_broadcast_f32x4( _mm_castpd_ps(pset1<Packet1cd>(from).v))));
    260   #endif
    261 }
    262 
    263 template<> EIGEN_STRONG_INLINE Packet4cd ploaddup<Packet4cd>(const std::complex<double>* from) {
    264   return Packet4cd(_mm512_insertf64x4(
    265           _mm512_castpd256_pd512(ploaddup<Packet2cd>(from).v), ploaddup<Packet2cd>(from+1).v, 1));
    266 }
    267 
    268 template<> EIGEN_STRONG_INLINE void pstore <std::complex<double> >(std::complex<double> *   to, const Packet4cd& from) { EIGEN_DEBUG_ALIGNED_STORE pstore((double*)to, from.v); }
    269 template<> EIGEN_STRONG_INLINE void pstoreu<std::complex<double> >(std::complex<double> *   to, const Packet4cd& from) { EIGEN_DEBUG_UNALIGNED_STORE pstoreu((double*)to, from.v); }
    270 
    271 template<> EIGEN_DEVICE_FUNC inline Packet4cd pgather<std::complex<double>, Packet4cd>(const std::complex<double>* from, Index stride)
    272 {
    273   return Packet4cd(_mm512_insertf64x4(_mm512_castpd256_pd512(
    274             _mm256_insertf128_pd(_mm256_castpd128_pd256(ploadu<Packet1cd>(from+0*stride).v), ploadu<Packet1cd>(from+1*stride).v,1)),
    275             _mm256_insertf128_pd(_mm256_castpd128_pd256(ploadu<Packet1cd>(from+2*stride).v), ploadu<Packet1cd>(from+3*stride).v,1), 1));
    276 }
    277 
    278 template<> EIGEN_DEVICE_FUNC inline void pscatter<std::complex<double>, Packet4cd>(std::complex<double>* to, const Packet4cd& from, Index stride)
    279 {
    280   __m512i fromi = _mm512_castpd_si512(from.v);
    281   double* tod = (double*)(void*)to;
    282   _mm_storeu_pd(tod+0*stride, _mm_castsi128_pd(_mm512_extracti32x4_epi32(fromi,0)) );
    283   _mm_storeu_pd(tod+2*stride, _mm_castsi128_pd(_mm512_extracti32x4_epi32(fromi,1)) );
    284   _mm_storeu_pd(tod+4*stride, _mm_castsi128_pd(_mm512_extracti32x4_epi32(fromi,2)) );
    285   _mm_storeu_pd(tod+6*stride, _mm_castsi128_pd(_mm512_extracti32x4_epi32(fromi,3)) );
    286 }
    287 
    288 template<> EIGEN_STRONG_INLINE std::complex<double> pfirst<Packet4cd>(const Packet4cd& a)
    289 {
    290   __m128d low = extract128<0>(a.v);
    291   EIGEN_ALIGN16 double res[2];
    292   _mm_store_pd(res, low);
    293   return std::complex<double>(res[0],res[1]);
    294 }
    295 
    296 template<> EIGEN_STRONG_INLINE Packet4cd preverse(const Packet4cd& a) {
    297   return Packet4cd(_mm512_shuffle_f64x2(a.v, a.v, (shuffle_mask<3,2,1,0>::mask)));
    298 }
    299 
    300 template<> EIGEN_STRONG_INLINE std::complex<double> predux<Packet4cd>(const Packet4cd& a)
    301 {
    302   return predux(padd(Packet2cd(_mm512_extractf64x4_pd(a.v,0)),
    303                      Packet2cd(_mm512_extractf64x4_pd(a.v,1))));
    304 }
    305 
    306 template<> EIGEN_STRONG_INLINE std::complex<double> predux_mul<Packet4cd>(const Packet4cd& a)
    307 {
    308   return predux_mul(pmul(Packet2cd(_mm512_extractf64x4_pd(a.v,0)),
    309                          Packet2cd(_mm512_extractf64x4_pd(a.v,1))));
    310 }
    311 
    312 template<> struct conj_helper<Packet4cd, Packet4cd, false,true>
    313 {
    314   EIGEN_STRONG_INLINE Packet4cd pmadd(const Packet4cd& x, const Packet4cd& y, const Packet4cd& c) const
    315   { return padd(pmul(x,y),c); }
    316 
    317   EIGEN_STRONG_INLINE Packet4cd pmul(const Packet4cd& a, const Packet4cd& b) const
    318   {
    319     return internal::pmul(a, pconj(b));
    320   }
    321 };
    322 
    323 template<> struct conj_helper<Packet4cd, Packet4cd, true,false>
    324 {
    325   EIGEN_STRONG_INLINE Packet4cd pmadd(const Packet4cd& x, const Packet4cd& y, const Packet4cd& c) const
    326   { return padd(pmul(x,y),c); }
    327 
    328   EIGEN_STRONG_INLINE Packet4cd pmul(const Packet4cd& a, const Packet4cd& b) const
    329   {
    330     return internal::pmul(pconj(a), b);
    331   }
    332 };
    333 
    334 template<> struct conj_helper<Packet4cd, Packet4cd, true,true>
    335 {
    336   EIGEN_STRONG_INLINE Packet4cd pmadd(const Packet4cd& x, const Packet4cd& y, const Packet4cd& c) const
    337   { return padd(pmul(x,y),c); }
    338 
    339   EIGEN_STRONG_INLINE Packet4cd pmul(const Packet4cd& a, const Packet4cd& b) const
    340   {
    341     return pconj(internal::pmul(a, b));
    342   }
    343 };
    344 
    345 EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet4cd,Packet8d)
    346 
    347 template<> EIGEN_STRONG_INLINE Packet4cd pdiv<Packet4cd>(const Packet4cd& a, const Packet4cd& b)
    348 {
    349   Packet4cd num = pmul(a, pconj(b));
    350   __m512d tmp = _mm512_mul_pd(b.v, b.v);
    351   __m512d denom =  padd(_mm512_permute_pd(tmp,0x55), tmp);
    352   return Packet4cd(_mm512_div_pd(num.v, denom));
    353 }
    354 
    355 template<> EIGEN_STRONG_INLINE Packet4cd pcplxflip<Packet4cd>(const Packet4cd& x)
    356 {
    357   return Packet4cd(_mm512_permute_pd(x.v,0x55));
    358 }
    359 
    360 EIGEN_DEVICE_FUNC inline void
    361 ptranspose(PacketBlock<Packet8cf,4>& kernel) {
    362   PacketBlock<Packet8d,4> pb;
    363   
    364   pb.packet[0] = _mm512_castps_pd(kernel.packet[0].v);
    365   pb.packet[1] = _mm512_castps_pd(kernel.packet[1].v);
    366   pb.packet[2] = _mm512_castps_pd(kernel.packet[2].v);
    367   pb.packet[3] = _mm512_castps_pd(kernel.packet[3].v);
    368   ptranspose(pb);
    369   kernel.packet[0].v = _mm512_castpd_ps(pb.packet[0]);
    370   kernel.packet[1].v = _mm512_castpd_ps(pb.packet[1]);
    371   kernel.packet[2].v = _mm512_castpd_ps(pb.packet[2]);
    372   kernel.packet[3].v = _mm512_castpd_ps(pb.packet[3]);
    373 }
    374 
    375 EIGEN_DEVICE_FUNC inline void
    376 ptranspose(PacketBlock<Packet8cf,8>& kernel) {
    377   PacketBlock<Packet8d,8> pb;
    378   
    379   pb.packet[0] = _mm512_castps_pd(kernel.packet[0].v);
    380   pb.packet[1] = _mm512_castps_pd(kernel.packet[1].v);
    381   pb.packet[2] = _mm512_castps_pd(kernel.packet[2].v);
    382   pb.packet[3] = _mm512_castps_pd(kernel.packet[3].v);
    383   pb.packet[4] = _mm512_castps_pd(kernel.packet[4].v);
    384   pb.packet[5] = _mm512_castps_pd(kernel.packet[5].v);
    385   pb.packet[6] = _mm512_castps_pd(kernel.packet[6].v);
    386   pb.packet[7] = _mm512_castps_pd(kernel.packet[7].v);
    387   ptranspose(pb);
    388   kernel.packet[0].v = _mm512_castpd_ps(pb.packet[0]);
    389   kernel.packet[1].v = _mm512_castpd_ps(pb.packet[1]);
    390   kernel.packet[2].v = _mm512_castpd_ps(pb.packet[2]);
    391   kernel.packet[3].v = _mm512_castpd_ps(pb.packet[3]);
    392   kernel.packet[4].v = _mm512_castpd_ps(pb.packet[4]);
    393   kernel.packet[5].v = _mm512_castpd_ps(pb.packet[5]);
    394   kernel.packet[6].v = _mm512_castpd_ps(pb.packet[6]);
    395   kernel.packet[7].v = _mm512_castpd_ps(pb.packet[7]);
    396 }
    397 
    398 EIGEN_DEVICE_FUNC inline void
    399 ptranspose(PacketBlock<Packet4cd,4>& kernel) {
    400   __m512d T0 = _mm512_shuffle_f64x2(kernel.packet[0].v, kernel.packet[1].v, (shuffle_mask<0,1,0,1>::mask)); // [a0 a1 b0 b1]
    401   __m512d T1 = _mm512_shuffle_f64x2(kernel.packet[0].v, kernel.packet[1].v, (shuffle_mask<2,3,2,3>::mask)); // [a2 a3 b2 b3]
    402   __m512d T2 = _mm512_shuffle_f64x2(kernel.packet[2].v, kernel.packet[3].v, (shuffle_mask<0,1,0,1>::mask)); // [c0 c1 d0 d1]
    403   __m512d T3 = _mm512_shuffle_f64x2(kernel.packet[2].v, kernel.packet[3].v, (shuffle_mask<2,3,2,3>::mask)); // [c2 c3 d2 d3]
    404 
    405   kernel.packet[3] = Packet4cd(_mm512_shuffle_f64x2(T1, T3, (shuffle_mask<1,3,1,3>::mask))); // [a3 b3 c3 d3]
    406   kernel.packet[2] = Packet4cd(_mm512_shuffle_f64x2(T1, T3, (shuffle_mask<0,2,0,2>::mask))); // [a2 b2 c2 d2]
    407   kernel.packet[1] = Packet4cd(_mm512_shuffle_f64x2(T0, T2, (shuffle_mask<1,3,1,3>::mask))); // [a1 b1 c1 d1]
    408   kernel.packet[0] = Packet4cd(_mm512_shuffle_f64x2(T0, T2, (shuffle_mask<0,2,0,2>::mask))); // [a0 b0 c0 d0]
    409 }
    410 
    411 template<> EIGEN_STRONG_INLINE Packet4cd psqrt<Packet4cd>(const Packet4cd& a) {
    412   return psqrt_complex<Packet4cd>(a);
    413 }
    414 
    415 template<> EIGEN_STRONG_INLINE Packet8cf psqrt<Packet8cf>(const Packet8cf& a) {
    416   return psqrt_complex<Packet8cf>(a);
    417 }
    418 
    419 } // end namespace internal
    420 } // end namespace Eigen
    421 
    422 #endif // EIGEN_COMPLEX_AVX512_H