// This file is part of Eigen, a lightweight C++ template library // for linear algebra. // // Copyright (C) 2010 Gael Guennebaud // // This Source Code Form is subject to the terms of the Mozilla // Public License v. 2.0. If a copy of the MPL was not distributed // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. #ifndef EIGEN_COMPLEX_SSE_H #define EIGEN_COMPLEX_SSE_H // IWYU pragma: private #include "../../InternalHeaderCheck.h" namespace Eigen { namespace internal { //---------- float ---------- struct Packet2cf { EIGEN_STRONG_INLINE Packet2cf() {} EIGEN_STRONG_INLINE explicit Packet2cf(const __m128& a) : v(a) {} Packet4f v; }; // Use the packet_traits defined in AVX/PacketMath.h instead if we're going // to leverage AVX instructions. #ifndef EIGEN_VECTORIZE_AVX template <> struct packet_traits > : default_packet_traits { typedef Packet2cf type; typedef Packet2cf half; enum { Vectorizable = 1, AlignedOnScalar = 1, size = 2, HasAdd = 1, HasSub = 1, HasMul = 1, HasDiv = 1, HasNegate = 1, HasSqrt = 1, HasLog = 1, HasExp = 1, HasAbs = 0, HasAbs2 = 0, HasMin = 0, HasMax = 0, HasSetLinear = 0, HasBlend = 1 }; }; #endif template <> struct unpacket_traits { typedef std::complex type; typedef Packet2cf half; typedef Packet4f as_real; enum { size = 2, alignment = Aligned16, vectorizable = true, masked_load_available = false, masked_store_available = false }; }; template <> EIGEN_STRONG_INLINE Packet2cf padd(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(_mm_add_ps(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet2cf psub(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(_mm_sub_ps(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet2cf pnegate(const Packet2cf& a) { const __m128 mask = _mm_castsi128_ps(_mm_setr_epi32(0x80000000, 0x80000000, 0x80000000, 0x80000000)); return Packet2cf(_mm_xor_ps(a.v, mask)); } template <> EIGEN_STRONG_INLINE Packet2cf pconj(const Packet2cf& a) { const __m128 mask = _mm_castsi128_ps(_mm_setr_epi32(0x00000000, 0x80000000, 0x00000000, 0x80000000)); return Packet2cf(_mm_xor_ps(a.v, mask)); } template <> EIGEN_STRONG_INLINE Packet2cf pmul(const Packet2cf& a, const Packet2cf& b) { #ifdef EIGEN_VECTORIZE_SSE3 return Packet2cf(_mm_addsub_ps(_mm_mul_ps(_mm_moveldup_ps(a.v), b.v), _mm_mul_ps(_mm_movehdup_ps(a.v), vec4f_swizzle1(b.v, 1, 0, 3, 2)))); // return Packet2cf(_mm_addsub_ps(_mm_mul_ps(vec4f_swizzle1(a.v, 0, 0, 2, 2), b.v), // _mm_mul_ps(vec4f_swizzle1(a.v, 1, 1, 3, 3), // vec4f_swizzle1(b.v, 1, 0, 3, 2)))); #else const __m128 mask = _mm_castsi128_ps(_mm_setr_epi32(0x80000000, 0x00000000, 0x80000000, 0x00000000)); return Packet2cf( _mm_add_ps(_mm_mul_ps(vec4f_swizzle1(a.v, 0, 0, 2, 2), b.v), _mm_xor_ps(_mm_mul_ps(vec4f_swizzle1(a.v, 1, 1, 3, 3), vec4f_swizzle1(b.v, 1, 0, 3, 2)), mask))); #endif } template <> EIGEN_STRONG_INLINE Packet2cf ptrue(const Packet2cf& a) { return Packet2cf(ptrue(Packet4f(a.v))); } template <> EIGEN_STRONG_INLINE Packet2cf pand(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(_mm_and_ps(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet2cf por(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(_mm_or_ps(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet2cf pxor(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(_mm_xor_ps(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet2cf pandnot(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(_mm_andnot_ps(b.v, a.v)); } template <> EIGEN_STRONG_INLINE Packet2cf pload(const std::complex* from) { EIGEN_DEBUG_ALIGNED_LOAD return Packet2cf(pload(&numext::real_ref(*from))); } template <> EIGEN_STRONG_INLINE Packet2cf ploadu(const std::complex* from) { EIGEN_DEBUG_UNALIGNED_LOAD return Packet2cf(ploadu(&numext::real_ref(*from))); } template <> EIGEN_STRONG_INLINE Packet2cf pset1(const std::complex& from) { const float re = std::real(from); const float im = std::imag(from); return Packet2cf(_mm_set_ps(im, re, im, re)); } template <> EIGEN_STRONG_INLINE Packet2cf ploaddup(const std::complex* from) { return pset1(*from); } template <> EIGEN_STRONG_INLINE void pstore >(std::complex* to, const Packet2cf& from) { EIGEN_DEBUG_ALIGNED_STORE pstore(&numext::real_ref(*to), Packet4f(from.v)); } template <> EIGEN_STRONG_INLINE void pstoreu >(std::complex* to, const Packet2cf& from) { EIGEN_DEBUG_UNALIGNED_STORE pstoreu(&numext::real_ref(*to), Packet4f(from.v)); } template <> EIGEN_DEVICE_FUNC inline Packet2cf pgather, Packet2cf>(const std::complex* from, Index stride) { return Packet2cf(_mm_set_ps(std::imag(from[1 * stride]), std::real(from[1 * stride]), std::imag(from[0 * stride]), std::real(from[0 * stride]))); } template <> EIGEN_DEVICE_FUNC inline void pscatter, Packet2cf>(std::complex* to, const Packet2cf& from, Index stride) { to[stride * 0] = std::complex(_mm_cvtss_f32(_mm_shuffle_ps(from.v, from.v, 0)), _mm_cvtss_f32(_mm_shuffle_ps(from.v, from.v, 1))); to[stride * 1] = std::complex(_mm_cvtss_f32(_mm_shuffle_ps(from.v, from.v, 2)), _mm_cvtss_f32(_mm_shuffle_ps(from.v, from.v, 3))); } template <> EIGEN_STRONG_INLINE void prefetch >(const std::complex* addr) { _mm_prefetch((SsePrefetchPtrType)(addr), _MM_HINT_T0); } template <> EIGEN_STRONG_INLINE std::complex pfirst(const Packet2cf& a) { alignas(alignof(__m64)) std::complex res; _mm_storel_pi((__m64*)&res, a.v); return res; } template <> EIGEN_STRONG_INLINE Packet2cf preverse(const Packet2cf& a) { return Packet2cf(_mm_castpd_ps(preverse(Packet2d(_mm_castps_pd(a.v))))); } template <> EIGEN_STRONG_INLINE std::complex predux(const Packet2cf& a) { return pfirst(Packet2cf(_mm_add_ps(a.v, _mm_movehl_ps(a.v, a.v)))); } template <> EIGEN_STRONG_INLINE std::complex predux_mul(const Packet2cf& a) { return pfirst(pmul(a, Packet2cf(_mm_movehl_ps(a.v, a.v)))); } EIGEN_STRONG_INLINE Packet2cf pcplxflip /* */ (const Packet2cf& x) { return Packet2cf(vec4f_swizzle1(x.v, 1, 0, 3, 2)); } EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet2cf, Packet4f) template <> EIGEN_STRONG_INLINE Packet2cf pdiv(const Packet2cf& a, const Packet2cf& b) { return pdiv_complex(a, b); } //---------- double ---------- struct Packet1cd { EIGEN_STRONG_INLINE Packet1cd() {} EIGEN_STRONG_INLINE explicit Packet1cd(const __m128d& a) : v(a) {} Packet2d v; }; // Use the packet_traits defined in AVX/PacketMath.h instead if we're going // to leverage AVX instructions. #ifndef EIGEN_VECTORIZE_AVX template <> struct packet_traits > : default_packet_traits { typedef Packet1cd type; typedef Packet1cd half; enum { Vectorizable = 1, AlignedOnScalar = 0, size = 1, HasAdd = 1, HasSub = 1, HasMul = 1, HasDiv = 1, HasNegate = 1, HasSqrt = 1, HasLog = 1, HasAbs = 0, HasAbs2 = 0, HasMin = 0, HasMax = 0, HasSetLinear = 0 }; }; #endif template <> struct unpacket_traits { typedef std::complex type; typedef Packet1cd half; typedef Packet2d as_real; enum { size = 1, alignment = Aligned16, vectorizable = true, masked_load_available = false, masked_store_available = false }; }; template <> EIGEN_STRONG_INLINE Packet1cd padd(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(_mm_add_pd(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet1cd psub(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(_mm_sub_pd(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet1cd pnegate(const Packet1cd& a) { return Packet1cd(pnegate(Packet2d(a.v))); } template <> EIGEN_STRONG_INLINE Packet1cd pconj(const Packet1cd& a) { const __m128d mask = _mm_castsi128_pd(_mm_set_epi32(0x80000000, 0x0, 0x0, 0x0)); return Packet1cd(_mm_xor_pd(a.v, mask)); } template <> EIGEN_STRONG_INLINE Packet1cd pmul(const Packet1cd& a, const Packet1cd& b) { #ifdef EIGEN_VECTORIZE_SSE3 return Packet1cd(_mm_addsub_pd(_mm_mul_pd(_mm_movedup_pd(a.v), b.v), _mm_mul_pd(vec2d_swizzle1(a.v, 1, 1), vec2d_swizzle1(b.v, 1, 0)))); #else const __m128d mask = _mm_castsi128_pd(_mm_set_epi32(0x0, 0x0, 0x80000000, 0x0)); return Packet1cd(_mm_add_pd(_mm_mul_pd(vec2d_swizzle1(a.v, 0, 0), b.v), _mm_xor_pd(_mm_mul_pd(vec2d_swizzle1(a.v, 1, 1), vec2d_swizzle1(b.v, 1, 0)), mask))); #endif } template <> EIGEN_STRONG_INLINE Packet1cd ptrue(const Packet1cd& a) { return Packet1cd(ptrue(Packet2d(a.v))); } template <> EIGEN_STRONG_INLINE Packet1cd pand(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(_mm_and_pd(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet1cd por(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(_mm_or_pd(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet1cd pxor(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(_mm_xor_pd(a.v, b.v)); } template <> EIGEN_STRONG_INLINE Packet1cd pandnot(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(_mm_andnot_pd(b.v, a.v)); } // FIXME force unaligned load, this is a temporary fix template <> EIGEN_STRONG_INLINE Packet1cd pload(const std::complex* from) { EIGEN_DEBUG_ALIGNED_LOAD return Packet1cd(pload((const double*)from)); } template <> EIGEN_STRONG_INLINE Packet1cd ploadu(const std::complex* from) { EIGEN_DEBUG_UNALIGNED_LOAD return Packet1cd(ploadu((const double*)from)); } template <> EIGEN_STRONG_INLINE Packet1cd pset1(const std::complex& from) { /* here we really have to use unaligned loads :( */ return ploadu(&from); } template <> EIGEN_STRONG_INLINE Packet1cd ploaddup(const std::complex* from) { return pset1(*from); } // FIXME force unaligned store, this is a temporary fix template <> EIGEN_STRONG_INLINE void pstore >(std::complex* to, const Packet1cd& from) { EIGEN_DEBUG_ALIGNED_STORE pstore((double*)to, Packet2d(from.v)); } template <> EIGEN_STRONG_INLINE void pstoreu >(std::complex* to, const Packet1cd& from) { EIGEN_DEBUG_UNALIGNED_STORE pstoreu((double*)to, Packet2d(from.v)); } template <> EIGEN_STRONG_INLINE void prefetch >(const std::complex* addr) { _mm_prefetch((SsePrefetchPtrType)(addr), _MM_HINT_T0); } template <> EIGEN_STRONG_INLINE std::complex pfirst(const Packet1cd& a) { EIGEN_ALIGN16 double res[2]; _mm_store_pd(res, a.v); return std::complex(res[0], res[1]); } template <> EIGEN_STRONG_INLINE Packet1cd preverse(const Packet1cd& a) { return a; } template <> EIGEN_STRONG_INLINE std::complex predux(const Packet1cd& a) { return pfirst(a); } template <> EIGEN_STRONG_INLINE std::complex predux_mul(const Packet1cd& a) { return pfirst(a); } EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet1cd, Packet2d) template <> EIGEN_STRONG_INLINE Packet1cd pdiv(const Packet1cd& a, const Packet1cd& b) { return pdiv_complex(a, b); } EIGEN_STRONG_INLINE Packet1cd pcplxflip /* */ (const Packet1cd& x) { return Packet1cd(preverse(Packet2d(x.v))); } EIGEN_DEVICE_FUNC inline void ptranspose(PacketBlock& kernel) { __m128d w1 = _mm_castps_pd(kernel.packet[0].v); __m128d w2 = _mm_castps_pd(kernel.packet[1].v); __m128 tmp = _mm_castpd_ps(_mm_unpackhi_pd(w1, w2)); kernel.packet[0].v = _mm_castpd_ps(_mm_unpacklo_pd(w1, w2)); kernel.packet[1].v = tmp; } template <> EIGEN_STRONG_INLINE Packet2cf pcmp_eq(const Packet2cf& a, const Packet2cf& b) { __m128 eq = _mm_cmpeq_ps(a.v, b.v); return Packet2cf(pand(eq, vec4f_swizzle1(eq, 1, 0, 3, 2))); } template <> EIGEN_STRONG_INLINE Packet1cd pcmp_eq(const Packet1cd& a, const Packet1cd& b) { __m128d eq = _mm_cmpeq_pd(a.v, b.v); return Packet1cd(pand(eq, vec2d_swizzle1(eq, 1, 0))); } template <> EIGEN_STRONG_INLINE Packet2cf pblend(const Selector<2>& ifPacket, const Packet2cf& thenPacket, const Packet2cf& elsePacket) { __m128d result = pblend(ifPacket, _mm_castps_pd(thenPacket.v), _mm_castps_pd(elsePacket.v)); return Packet2cf(_mm_castpd_ps(result)); } template <> EIGEN_STRONG_INLINE Packet1cd psqrt(const Packet1cd& a) { return psqrt_complex(a); } template <> EIGEN_STRONG_INLINE Packet2cf psqrt(const Packet2cf& a) { return psqrt_complex(a); } template <> EIGEN_STRONG_INLINE Packet1cd plog(const Packet1cd& a) { return plog_complex(a); } template <> EIGEN_STRONG_INLINE Packet2cf plog(const Packet2cf& a) { return plog_complex(a); } template <> EIGEN_STRONG_INLINE Packet2cf pexp(const Packet2cf& a) { return pexp_complex(a); } } // end namespace internal } // end namespace Eigen #endif // EIGEN_COMPLEX_SSE_H