| // This file is part of Eigen, a lightweight C++ template library |
| // for linear algebra. |
| // |
| // Copyright (C) 2008-2009 Gael Guennebaud <gael.guennebaud@inria.fr> |
| // |
| // 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_GENERAL_BLOCK_PANEL_H |
| #define EIGEN_GENERAL_BLOCK_PANEL_H |
| |
| // IWYU pragma: private |
| #include "../InternalHeaderCheck.h" |
| |
| namespace Eigen { |
| |
| namespace internal { |
| |
| enum GEBPPacketSizeType { GEBPPacketFull = 0, GEBPPacketHalf, GEBPPacketQuarter }; |
| |
| template <typename LhsScalar_, typename RhsScalar_, bool ConjLhs_ = false, bool ConjRhs_ = false, |
| int Arch = Architecture::Target, int PacketSize_ = GEBPPacketFull> |
| class gebp_traits; |
| |
| /** \internal \returns b if a<=0, and returns a otherwise. */ |
| inline std::ptrdiff_t manage_caching_sizes_helper(std::ptrdiff_t a, std::ptrdiff_t b) { return a <= 0 ? b : a; } |
| |
| #if defined(EIGEN_DEFAULT_L1_CACHE_SIZE) |
| #define EIGEN_SET_DEFAULT_L1_CACHE_SIZE(val) EIGEN_DEFAULT_L1_CACHE_SIZE |
| #else |
| #define EIGEN_SET_DEFAULT_L1_CACHE_SIZE(val) val |
| #endif // defined(EIGEN_DEFAULT_L1_CACHE_SIZE) |
| |
| #if defined(EIGEN_DEFAULT_L2_CACHE_SIZE) |
| #define EIGEN_SET_DEFAULT_L2_CACHE_SIZE(val) EIGEN_DEFAULT_L2_CACHE_SIZE |
| #else |
| #define EIGEN_SET_DEFAULT_L2_CACHE_SIZE(val) val |
| #endif // defined(EIGEN_DEFAULT_L2_CACHE_SIZE) |
| |
| #if defined(EIGEN_DEFAULT_L3_CACHE_SIZE) |
| #define EIGEN_SET_DEFAULT_L3_CACHE_SIZE(val) EIGEN_DEFAULT_L3_CACHE_SIZE |
| #else |
| #define EIGEN_SET_DEFAULT_L3_CACHE_SIZE(val) val |
| #endif // defined(EIGEN_DEFAULT_L3_CACHE_SIZE) |
| |
| #if EIGEN_ARCH_i386_OR_x86_64 |
| const std::ptrdiff_t defaultL1CacheSize = EIGEN_SET_DEFAULT_L1_CACHE_SIZE(32 * 1024); |
| const std::ptrdiff_t defaultL2CacheSize = EIGEN_SET_DEFAULT_L2_CACHE_SIZE(256 * 1024); |
| const std::ptrdiff_t defaultL3CacheSize = EIGEN_SET_DEFAULT_L3_CACHE_SIZE(2 * 1024 * 1024); |
| #elif EIGEN_ARCH_PPC |
| const std::ptrdiff_t defaultL1CacheSize = EIGEN_SET_DEFAULT_L1_CACHE_SIZE(64 * 1024); |
| #ifdef _ARCH_PWR10 |
| const std::ptrdiff_t defaultL2CacheSize = EIGEN_SET_DEFAULT_L2_CACHE_SIZE(2 * 1024 * 1024); |
| const std::ptrdiff_t defaultL3CacheSize = EIGEN_SET_DEFAULT_L3_CACHE_SIZE(8 * 1024 * 1024); |
| #else |
| const std::ptrdiff_t defaultL2CacheSize = EIGEN_SET_DEFAULT_L2_CACHE_SIZE(512 * 1024); |
| const std::ptrdiff_t defaultL3CacheSize = EIGEN_SET_DEFAULT_L3_CACHE_SIZE(4 * 1024 * 1024); |
| #endif |
| #else |
| const std::ptrdiff_t defaultL1CacheSize = EIGEN_SET_DEFAULT_L1_CACHE_SIZE(16 * 1024); |
| const std::ptrdiff_t defaultL2CacheSize = EIGEN_SET_DEFAULT_L2_CACHE_SIZE(512 * 1024); |
| const std::ptrdiff_t defaultL3CacheSize = EIGEN_SET_DEFAULT_L3_CACHE_SIZE(512 * 1024); |
| #endif |
| |
| #undef EIGEN_SET_DEFAULT_L1_CACHE_SIZE |
| #undef EIGEN_SET_DEFAULT_L2_CACHE_SIZE |
| #undef EIGEN_SET_DEFAULT_L3_CACHE_SIZE |
| |
| /** \internal */ |
| struct CacheSizes { |
| CacheSizes() : m_l1(-1), m_l2(-1), m_l3(-1) { |
| int l1CacheSize, l2CacheSize, l3CacheSize; |
| queryCacheSizes(l1CacheSize, l2CacheSize, l3CacheSize); |
| m_l1 = manage_caching_sizes_helper(l1CacheSize, defaultL1CacheSize); |
| m_l2 = manage_caching_sizes_helper(l2CacheSize, defaultL2CacheSize); |
| m_l3 = manage_caching_sizes_helper(l3CacheSize, defaultL3CacheSize); |
| } |
| |
| std::ptrdiff_t m_l1; |
| std::ptrdiff_t m_l2; |
| std::ptrdiff_t m_l3; |
| }; |
| |
| /** \internal */ |
| inline void manage_caching_sizes(Action action, std::ptrdiff_t* l1, std::ptrdiff_t* l2, std::ptrdiff_t* l3) { |
| static CacheSizes m_cacheSizes; |
| |
| if (action == SetAction) { |
| // set the cpu cache size and cache all block sizes from a global cache size in byte |
| eigen_internal_assert(l1 != 0 && l2 != 0); |
| m_cacheSizes.m_l1 = *l1; |
| m_cacheSizes.m_l2 = *l2; |
| m_cacheSizes.m_l3 = *l3; |
| } else if (action == GetAction) { |
| eigen_internal_assert(l1 != 0 && l2 != 0); |
| *l1 = m_cacheSizes.m_l1; |
| *l2 = m_cacheSizes.m_l2; |
| *l3 = m_cacheSizes.m_l3; |
| } else { |
| eigen_internal_assert(false); |
| } |
| } |
| |
| /* Helper for computeProductBlockingSizes. |
| * |
| * Given a m x k times k x n matrix product of scalar types \c LhsScalar and \c RhsScalar, |
| * this function computes the blocking size parameters along the respective dimensions |
| * for matrix products and related algorithms. The blocking sizes depends on various |
| * parameters: |
| * - the L1 and L2 cache sizes, |
| * - the register level blocking sizes defined by gebp_traits, |
| * - the number of scalars that fit into a packet (when vectorization is enabled). |
| * |
| * \sa setCpuCacheSizes */ |
| |
| template <typename LhsScalar, typename RhsScalar, int KcFactor, typename Index> |
| void evaluateProductBlockingSizesHeuristic(Index& k, Index& m, Index& n, Index num_threads = 1) { |
| typedef gebp_traits<LhsScalar, RhsScalar> Traits; |
| |
| // Explanations: |
| // Let's recall that the product algorithms form mc x kc vertical panels A' on the lhs and |
| // kc x nc blocks B' on the rhs. B' has to fit into L2/L3 cache. Moreover, A' is processed |
| // per mr x kc horizontal small panels where mr is the blocking size along the m dimension |
| // at the register level. This small horizontal panel has to stay within L1 cache. |
| std::ptrdiff_t l1, l2, l3; |
| manage_caching_sizes(GetAction, &l1, &l2, &l3); |
| #ifdef EIGEN_VECTORIZE_AVX512 |
| // We need to find a rationale for that, but without this adjustment, |
| // performance with AVX512 is pretty bad, like -20% slower. |
| // One reason is that with increasing packet-size, the blocking size k |
| // has to become pretty small if we want that 1 lhs panel fit within L1. |
| // For instance, with the 3pX4 kernel and double, the size of the lhs+rhs panels are: |
| // k*(3*64 + 4*8) Bytes, with l1=32kBytes, and k%8=0, we have k=144. |
| // This is quite small for a good reuse of the accumulation registers. |
| l1 *= 4; |
| #endif |
| |
| if (num_threads > 1) { |
| typedef typename Traits::ResScalar ResScalar; |
| enum { |
| kdiv = KcFactor * (Traits::mr * sizeof(LhsScalar) + Traits::nr * sizeof(RhsScalar)), |
| ksub = Traits::mr * (Traits::nr * sizeof(ResScalar)), |
| kr = 8, |
| mr = Traits::mr, |
| nr = Traits::nr |
| }; |
| // Increasing k gives us more time to prefetch the content of the "C" |
| // registers. However once the latency is hidden there is no point in |
| // increasing the value of k, so we'll cap it at 320 (value determined |
| // experimentally). |
| // To avoid that k vanishes, we make k_cache at least as big as kr |
| const Index k_cache = numext::maxi<Index>(kr, (numext::mini<Index>)((l1 - ksub) / kdiv, 320)); |
| if (k_cache < k) { |
| k = k_cache - (k_cache % kr); |
| eigen_internal_assert(k > 0); |
| } |
| |
| const Index n_cache = (l2 - l1) / (nr * sizeof(RhsScalar) * k); |
| const Index n_per_thread = numext::div_ceil(n, num_threads); |
| if (n_cache <= n_per_thread) { |
| // Don't exceed the capacity of the l2 cache. |
| eigen_internal_assert(n_cache >= static_cast<Index>(nr)); |
| n = n_cache - (n_cache % nr); |
| eigen_internal_assert(n > 0); |
| } else { |
| n = (numext::mini<Index>)(n, (n_per_thread + nr - 1) - ((n_per_thread + nr - 1) % nr)); |
| } |
| |
| if (l3 > l2) { |
| // l3 is shared between all cores, so we'll give each thread its own chunk of l3. |
| const Index m_cache = (l3 - l2) / (sizeof(LhsScalar) * k * num_threads); |
| const Index m_per_thread = numext::div_ceil(m, num_threads); |
| if (m_cache < m_per_thread && m_cache >= static_cast<Index>(mr)) { |
| m = m_cache - (m_cache % mr); |
| eigen_internal_assert(m > 0); |
| } else { |
| m = (numext::mini<Index>)(m, (m_per_thread + mr - 1) - ((m_per_thread + mr - 1) % mr)); |
| } |
| } |
| } else { |
| // In unit tests we do not want to use extra large matrices, |
| // so we reduce the cache size to check the blocking strategy is not flawed |
| #ifdef EIGEN_DEBUG_SMALL_PRODUCT_BLOCKS |
| l1 = 9 * 1024; |
| l2 = 32 * 1024; |
| l3 = 512 * 1024; |
| #endif |
| |
| // Early return for small problems because the computation below are time consuming for small problems. |
| // Perhaps it would make more sense to consider k*n*m?? |
| // Note that for very tiny problem, this function should be bypassed anyway |
| // because we use the coefficient-based implementation for them. |
| if ((numext::maxi)(k, (numext::maxi)(m, n)) < 48) return; |
| |
| typedef typename Traits::ResScalar ResScalar; |
| enum { |
| k_peeling = 8, |
| k_div = KcFactor * (Traits::mr * sizeof(LhsScalar) + Traits::nr * sizeof(RhsScalar)), |
| k_sub = Traits::mr * (Traits::nr * sizeof(ResScalar)) |
| }; |
| |
| // ---- 1st level of blocking on L1, yields kc ---- |
| |
| // Blocking on the third dimension (i.e., k) is chosen so that an horizontal panel |
| // of size mr x kc of the lhs plus a vertical panel of kc x nr of the rhs both fits within L1 cache. |
| // We also include a register-level block of the result (mx x nr). |
| // (In an ideal world only the lhs panel would stay in L1) |
| // Moreover, kc has to be a multiple of 8 to be compatible with loop peeling, leading to a maximum blocking size of: |
| const Index max_kc = numext::maxi<Index>(((l1 - k_sub) / k_div) & (~(k_peeling - 1)), 1); |
| const Index old_k = k; |
| if (k > max_kc) { |
| // We are really blocking on the third dimension: |
| // -> reduce blocking size to make sure the last block is as large as possible |
| // while keeping the same number of sweeps over the result. |
| k = (k % max_kc) == 0 ? max_kc |
| : max_kc - k_peeling * ((max_kc - 1 - (k % max_kc)) / (k_peeling * (k / max_kc + 1))); |
| |
| eigen_internal_assert(((old_k / k) == (old_k / max_kc)) && "the number of sweeps has to remain the same"); |
| } |
| |
| // ---- 2nd level of blocking on max(L2,L3), yields nc ---- |
| |
| // TODO find a reliable way to get the actual amount of cache per core to use for 2nd level blocking, that is: |
| // actual_l2 = max(l2, l3/nb_core_sharing_l3) |
| // The number below is quite conservative: it is better to underestimate the cache size rather than overestimating it) |
| // For instance, it corresponds to 6MB of L3 shared among 4 cores. |
| #ifdef EIGEN_DEBUG_SMALL_PRODUCT_BLOCKS |
| const Index actual_l2 = l3; |
| #else |
| const Index actual_l2 = 1572864; // == 1.5 MB |
| #endif |
| |
| // Here, nc is chosen such that a block of kc x nc of the rhs fit within half of L2. |
| // The second half is implicitly reserved to access the result and lhs coefficients. |
| // When k<max_kc, then nc can arbitrarily growth. In practice, it seems to be fruitful |
| // to limit this growth: we bound nc to growth by a factor x1.5. |
| // However, if the entire lhs block fit within L1, then we are not going to block on the rows at all, |
| // and it becomes fruitful to keep the packed rhs blocks in L1 if there is enough remaining space. |
| Index max_nc; |
| const Index lhs_bytes = m * k * sizeof(LhsScalar); |
| const Index remaining_l1 = l1 - k_sub - lhs_bytes; |
| if (remaining_l1 >= Index(Traits::nr * sizeof(RhsScalar)) * k) { |
| // L1 blocking |
| max_nc = remaining_l1 / (k * sizeof(RhsScalar)); |
| } else { |
| // L2 blocking |
| max_nc = (3 * actual_l2) / (2 * 2 * max_kc * sizeof(RhsScalar)); |
| } |
| // WARNING Below, we assume that Traits::nr is a power of two. |
| Index nc = numext::mini<Index>(actual_l2 / (2 * k * sizeof(RhsScalar)), max_nc) & (~(Traits::nr - 1)); |
| if (n > nc) { |
| // We are really blocking over the columns: |
| // -> reduce blocking size to make sure the last block is as large as possible |
| // while keeping the same number of sweeps over the packed lhs. |
| // Here we allow one more sweep if this gives us a perfect match, thus the commented "-1" |
| n = (n % nc) == 0 ? nc : (nc - Traits::nr * ((nc /*-1*/ - (n % nc)) / (Traits::nr * (n / nc + 1)))); |
| } else if (old_k == k) { |
| // So far, no blocking at all, i.e., kc==k, and nc==n. |
| // In this case, let's perform a blocking over the rows such that the packed lhs data is kept in cache L1/L2 |
| // TODO: part of this blocking strategy is now implemented within the kernel itself, so the L1-based heuristic |
| // here should be obsolete. |
| Index problem_size = k * n * sizeof(LhsScalar); |
| Index actual_lm = actual_l2; |
| Index max_mc = m; |
| if (problem_size <= 1024) { |
| // problem is small enough to keep in L1 |
| // Let's choose m such that lhs's block fit in 1/3 of L1 |
| actual_lm = l1; |
| } else if (l3 != 0 && problem_size <= 32768) { |
| // we have both L2 and L3, and problem is small enough to be kept in L2 |
| // Let's choose m such that lhs's block fit in 1/3 of L2 |
| actual_lm = l2; |
| max_mc = (numext::mini<Index>)(576, max_mc); |
| } |
| Index mc = (numext::mini<Index>)(actual_lm / (3 * k * sizeof(LhsScalar)), max_mc); |
| if (mc > Traits::mr) |
| mc -= mc % Traits::mr; |
| else if (mc == 0) |
| return; |
| m = (m % mc) == 0 ? mc : (mc - Traits::mr * ((mc /*-1*/ - (m % mc)) / (Traits::mr * (m / mc + 1)))); |
| } |
| } |
| } |
| |
| template <typename Index> |
| inline bool useSpecificBlockingSizes(Index& k, Index& m, Index& n) { |
| #ifdef EIGEN_TEST_SPECIFIC_BLOCKING_SIZES |
| if (EIGEN_TEST_SPECIFIC_BLOCKING_SIZES) { |
| k = numext::mini<Index>(k, EIGEN_TEST_SPECIFIC_BLOCKING_SIZE_K); |
| m = numext::mini<Index>(m, EIGEN_TEST_SPECIFIC_BLOCKING_SIZE_M); |
| n = numext::mini<Index>(n, EIGEN_TEST_SPECIFIC_BLOCKING_SIZE_N); |
| return true; |
| } |
| #else |
| EIGEN_UNUSED_VARIABLE(k) |
| EIGEN_UNUSED_VARIABLE(m) |
| EIGEN_UNUSED_VARIABLE(n) |
| #endif |
| return false; |
| } |
| |
| /** \brief Computes the blocking parameters for a m x k times k x n matrix product |
| * |
| * \param[in,out] k Input: the third dimension of the product. Output: the blocking size along the same dimension. |
| * \param[in,out] m Input: the number of rows of the left hand side. Output: the blocking size along the same dimension. |
| * \param[in,out] n Input: the number of columns of the right hand side. Output: the blocking size along the same |
| * dimension. |
| * |
| * Given a m x k times k x n matrix product of scalar types \c LhsScalar and \c RhsScalar, |
| * this function computes the blocking size parameters along the respective dimensions |
| * for matrix products and related algorithms. |
| * |
| * The blocking size parameters may be evaluated: |
| * - either by a heuristic based on cache sizes; |
| * - or using fixed prescribed values (for testing purposes). |
| * |
| * \sa setCpuCacheSizes */ |
| |
| template <typename LhsScalar, typename RhsScalar, int KcFactor, typename Index> |
| void computeProductBlockingSizes(Index& k, Index& m, Index& n, Index num_threads = 1) { |
| if (!useSpecificBlockingSizes(k, m, n)) { |
| evaluateProductBlockingSizesHeuristic<LhsScalar, RhsScalar, KcFactor, Index>(k, m, n, num_threads); |
| } |
| } |
| |
| template <typename LhsScalar, typename RhsScalar, typename Index> |
| inline void computeProductBlockingSizes(Index& k, Index& m, Index& n, Index num_threads = 1) { |
| computeProductBlockingSizes<LhsScalar, RhsScalar, 1, Index>(k, m, n, num_threads); |
| } |
| |
| template <typename RhsPacket, typename RhsPacketx4, int registers_taken> |
| struct RhsPanelHelper { |
| private: |
| static constexpr int remaining_registers = |
| (std::max)(int(EIGEN_ARCH_DEFAULT_NUMBER_OF_REGISTERS) - registers_taken, 0); |
| |
| public: |
| typedef std::conditional_t<remaining_registers >= 4, RhsPacketx4, RhsPacket> type; |
| }; |
| |
| template <typename Packet> |
| struct QuadPacket { |
| Packet B_0, B1, B2, B3; |
| const Packet& get(const FixedInt<0>&) const { return B_0; } |
| const Packet& get(const FixedInt<1>&) const { return B1; } |
| const Packet& get(const FixedInt<2>&) const { return B2; } |
| const Packet& get(const FixedInt<3>&) const { return B3; } |
| }; |
| |
| template <int N, typename T1, typename T2, typename T3> |
| struct packet_conditional { |
| typedef T3 type; |
| }; |
| |
| template <typename T1, typename T2, typename T3> |
| struct packet_conditional<GEBPPacketFull, T1, T2, T3> { |
| typedef T1 type; |
| }; |
| |
| template <typename T1, typename T2, typename T3> |
| struct packet_conditional<GEBPPacketHalf, T1, T2, T3> { |
| typedef T2 type; |
| }; |
| |
| #define PACKET_DECL_COND_POSTFIX(postfix, name, packet_size) \ |
| typedef typename packet_conditional< \ |
| packet_size, typename packet_traits<name##Scalar>::type, typename packet_traits<name##Scalar>::half, \ |
| typename unpacket_traits<typename packet_traits<name##Scalar>::half>::half>::type name##Packet##postfix |
| |
| #define PACKET_DECL_COND(name, packet_size) \ |
| typedef typename packet_conditional< \ |
| packet_size, typename packet_traits<name##Scalar>::type, typename packet_traits<name##Scalar>::half, \ |
| typename unpacket_traits<typename packet_traits<name##Scalar>::half>::half>::type name##Packet |
| |
| #define PACKET_DECL_COND_SCALAR_POSTFIX(postfix, packet_size) \ |
| typedef typename packet_conditional< \ |
| packet_size, typename packet_traits<Scalar>::type, typename packet_traits<Scalar>::half, \ |
| typename unpacket_traits<typename packet_traits<Scalar>::half>::half>::type ScalarPacket##postfix |
| |
| #define PACKET_DECL_COND_SCALAR(packet_size) \ |
| typedef typename packet_conditional< \ |
| packet_size, typename packet_traits<Scalar>::type, typename packet_traits<Scalar>::half, \ |
| typename unpacket_traits<typename packet_traits<Scalar>::half>::half>::type ScalarPacket |
| |
| /* Vectorization logic |
| * real*real: unpack rhs to constant packets, ... |
| * |
| * cd*cd : unpack rhs to (b_r,b_r), (b_i,b_i), mul to get (a_r b_r,a_i b_r) (a_r b_i,a_i b_i), |
| * storing each res packet into two packets (2x2), |
| * at the end combine them: swap the second and addsub them |
| * cf*cf : same but with 2x4 blocks |
| * cplx*real : unpack rhs to constant packets, ... |
| * real*cplx : load lhs as (a0,a0,a1,a1), and mul as usual |
| */ |
| template <typename LhsScalar_, typename RhsScalar_, bool ConjLhs_, bool ConjRhs_, int Arch, int PacketSize_> |
| class gebp_traits { |
| public: |
| typedef LhsScalar_ LhsScalar; |
| typedef RhsScalar_ RhsScalar; |
| typedef typename ScalarBinaryOpTraits<LhsScalar, RhsScalar>::ReturnType ResScalar; |
| |
| PACKET_DECL_COND_POSTFIX(_, Lhs, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Rhs, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Res, PacketSize_); |
| |
| enum { |
| ConjLhs = ConjLhs_, |
| ConjRhs = ConjRhs_, |
| Vectorizable = unpacket_traits<LhsPacket_>::vectorizable && unpacket_traits<RhsPacket_>::vectorizable, |
| LhsPacketSize = Vectorizable ? unpacket_traits<LhsPacket_>::size : 1, |
| RhsPacketSize = Vectorizable ? unpacket_traits<RhsPacket_>::size : 1, |
| ResPacketSize = Vectorizable ? unpacket_traits<ResPacket_>::size : 1, |
| |
| NumberOfRegisters = EIGEN_ARCH_DEFAULT_NUMBER_OF_REGISTERS, |
| |
| // register block size along the N direction must be 1 or 4 |
| nr = 4, |
| |
| // register block size along the M direction (currently, this one cannot be modified) |
| default_mr = (plain_enum_min(16, NumberOfRegisters) / 2 / nr) * LhsPacketSize, |
| #if defined(EIGEN_HAS_SINGLE_INSTRUCTION_MADD) && !defined(EIGEN_VECTORIZE_ALTIVEC) && \ |
| !defined(EIGEN_VECTORIZE_VSX) && ((!EIGEN_COMP_MSVC) || (EIGEN_COMP_MSVC >= 1914)) |
| // we assume 16 registers or more |
| // See bug 992, if the scalar type is not vectorizable but that EIGEN_HAS_SINGLE_INSTRUCTION_MADD is defined, |
| // then using 3*LhsPacketSize triggers non-implemented paths in syrk. |
| // Bug 1515: MSVC prior to v19.14 yields to register spilling. |
| mr = Vectorizable ? 3 * LhsPacketSize : default_mr, |
| #else |
| mr = default_mr, |
| #endif |
| |
| LhsProgress = LhsPacketSize, |
| RhsProgress = 1 |
| }; |
| |
| typedef std::conditional_t<Vectorizable, LhsPacket_, LhsScalar> LhsPacket; |
| typedef std::conditional_t<Vectorizable, RhsPacket_, RhsScalar> RhsPacket; |
| typedef std::conditional_t<Vectorizable, ResPacket_, ResScalar> ResPacket; |
| typedef LhsPacket LhsPacket4Packing; |
| |
| typedef QuadPacket<RhsPacket> RhsPacketx4; |
| typedef ResPacket AccPacket; |
| |
| EIGEN_STRONG_INLINE void initAcc(AccPacket& p) { p = pset1<ResPacket>(ResScalar(0)); } |
| |
| template <typename RhsPacketType> |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, RhsPacketType& dest) const { |
| dest = pset1<RhsPacketType>(*b); |
| } |
| |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, RhsPacketx4& dest) const { |
| pbroadcast4(b, dest.B_0, dest.B1, dest.B2, dest.B3); |
| } |
| |
| template <typename RhsPacketType> |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar* b, RhsPacketType& dest) const { |
| loadRhs(b, dest); |
| } |
| |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar*, RhsPacketx4&) const {} |
| |
| EIGEN_STRONG_INLINE void loadRhsQuad(const RhsScalar* b, RhsPacket& dest) const { dest = ploadquad<RhsPacket>(b); } |
| |
| template <typename LhsPacketType> |
| EIGEN_STRONG_INLINE void loadLhs(const LhsScalar* a, LhsPacketType& dest) const { |
| dest = pload<LhsPacketType>(a); |
| } |
| |
| template <typename LhsPacketType> |
| EIGEN_STRONG_INLINE void loadLhsUnaligned(const LhsScalar* a, LhsPacketType& dest) const { |
| dest = ploadu<LhsPacketType>(a); |
| } |
| |
| template <typename LhsPacketType, typename RhsPacketType, typename AccPacketType, typename LaneIdType> |
| EIGEN_STRONG_INLINE void madd(const LhsPacketType& a, const RhsPacketType& b, AccPacketType& c, RhsPacketType& tmp, |
| const LaneIdType&) const { |
| conj_helper<LhsPacketType, RhsPacketType, ConjLhs, ConjRhs> cj; |
| // It would be a lot cleaner to call pmadd all the time. Unfortunately if we |
| // let gcc allocate the register in which to store the result of the pmul |
| // (in the case where there is no FMA) gcc fails to figure out how to avoid |
| // spilling register. |
| #ifdef EIGEN_HAS_SINGLE_INSTRUCTION_MADD |
| EIGEN_UNUSED_VARIABLE(tmp); |
| c = cj.pmadd(a, b, c); |
| #else |
| tmp = b; |
| tmp = cj.pmul(a, tmp); |
| c = padd(c, tmp); |
| #endif |
| } |
| |
| template <typename LhsPacketType, typename AccPacketType, typename LaneIdType> |
| EIGEN_STRONG_INLINE void madd(const LhsPacketType& a, const RhsPacketx4& b, AccPacketType& c, RhsPacket& tmp, |
| const LaneIdType& lane) const { |
| madd(a, b.get(lane), c, tmp, lane); |
| } |
| |
| EIGEN_STRONG_INLINE void acc(const AccPacket& c, const ResPacket& alpha, ResPacket& r) const { |
| r = pmadd(c, alpha, r); |
| } |
| |
| template <typename ResPacketHalf> |
| EIGEN_STRONG_INLINE void acc(const ResPacketHalf& c, const ResPacketHalf& alpha, ResPacketHalf& r) const { |
| r = pmadd(c, alpha, r); |
| } |
| }; |
| |
| template <typename RealScalar, bool ConjLhs_, int Arch, int PacketSize_> |
| class gebp_traits<std::complex<RealScalar>, RealScalar, ConjLhs_, false, Arch, PacketSize_> { |
| public: |
| typedef std::complex<RealScalar> LhsScalar; |
| typedef RealScalar RhsScalar; |
| typedef typename ScalarBinaryOpTraits<LhsScalar, RhsScalar>::ReturnType ResScalar; |
| |
| PACKET_DECL_COND_POSTFIX(_, Lhs, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Rhs, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Res, PacketSize_); |
| |
| enum { |
| ConjLhs = ConjLhs_, |
| ConjRhs = false, |
| Vectorizable = unpacket_traits<LhsPacket_>::vectorizable && unpacket_traits<RhsPacket_>::vectorizable, |
| LhsPacketSize = Vectorizable ? unpacket_traits<LhsPacket_>::size : 1, |
| RhsPacketSize = Vectorizable ? unpacket_traits<RhsPacket_>::size : 1, |
| ResPacketSize = Vectorizable ? unpacket_traits<ResPacket_>::size : 1, |
| |
| NumberOfRegisters = EIGEN_ARCH_DEFAULT_NUMBER_OF_REGISTERS, |
| nr = 4, |
| #if defined(EIGEN_HAS_SINGLE_INSTRUCTION_MADD) && !defined(EIGEN_VECTORIZE_ALTIVEC) && !defined(EIGEN_VECTORIZE_VSX) |
| // we assume 16 registers |
| mr = 3 * LhsPacketSize, |
| #else |
| mr = (plain_enum_min(16, NumberOfRegisters) / 2 / nr) * LhsPacketSize, |
| #endif |
| |
| LhsProgress = LhsPacketSize, |
| RhsProgress = 1 |
| }; |
| |
| typedef std::conditional_t<Vectorizable, LhsPacket_, LhsScalar> LhsPacket; |
| typedef std::conditional_t<Vectorizable, RhsPacket_, RhsScalar> RhsPacket; |
| typedef std::conditional_t<Vectorizable, ResPacket_, ResScalar> ResPacket; |
| typedef LhsPacket LhsPacket4Packing; |
| |
| typedef QuadPacket<RhsPacket> RhsPacketx4; |
| |
| typedef ResPacket AccPacket; |
| |
| EIGEN_STRONG_INLINE void initAcc(AccPacket& p) { p = pset1<ResPacket>(ResScalar(0)); } |
| |
| template <typename RhsPacketType> |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, RhsPacketType& dest) const { |
| dest = pset1<RhsPacketType>(*b); |
| } |
| |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, RhsPacketx4& dest) const { |
| pbroadcast4(b, dest.B_0, dest.B1, dest.B2, dest.B3); |
| } |
| |
| template <typename RhsPacketType> |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar* b, RhsPacketType& dest) const { |
| loadRhs(b, dest); |
| } |
| |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar*, RhsPacketx4&) const {} |
| |
| EIGEN_STRONG_INLINE void loadRhsQuad(const RhsScalar* b, RhsPacket& dest) const { |
| loadRhsQuad_impl(b, dest, std::conditional_t<RhsPacketSize == 16, true_type, false_type>()); |
| } |
| |
| EIGEN_STRONG_INLINE void loadRhsQuad_impl(const RhsScalar* b, RhsPacket& dest, const true_type&) const { |
| // FIXME we can do better! |
| // what we want here is a ploadheight |
| RhsScalar tmp[4] = {b[0], b[0], b[1], b[1]}; |
| dest = ploadquad<RhsPacket>(tmp); |
| } |
| |
| EIGEN_STRONG_INLINE void loadRhsQuad_impl(const RhsScalar* b, RhsPacket& dest, const false_type&) const { |
| eigen_internal_assert(RhsPacketSize <= 8); |
| dest = pset1<RhsPacket>(*b); |
| } |
| |
| EIGEN_STRONG_INLINE void loadLhs(const LhsScalar* a, LhsPacket& dest) const { dest = pload<LhsPacket>(a); } |
| |
| template <typename LhsPacketType> |
| EIGEN_STRONG_INLINE void loadLhsUnaligned(const LhsScalar* a, LhsPacketType& dest) const { |
| dest = ploadu<LhsPacketType>(a); |
| } |
| |
| template <typename LhsPacketType, typename RhsPacketType, typename AccPacketType, typename LaneIdType> |
| EIGEN_STRONG_INLINE void madd(const LhsPacketType& a, const RhsPacketType& b, AccPacketType& c, RhsPacketType& tmp, |
| const LaneIdType&) const { |
| madd_impl(a, b, c, tmp, std::conditional_t<Vectorizable, true_type, false_type>()); |
| } |
| |
| template <typename LhsPacketType, typename RhsPacketType, typename AccPacketType> |
| EIGEN_STRONG_INLINE void madd_impl(const LhsPacketType& a, const RhsPacketType& b, AccPacketType& c, |
| RhsPacketType& tmp, const true_type&) const { |
| #ifdef EIGEN_HAS_SINGLE_INSTRUCTION_MADD |
| EIGEN_UNUSED_VARIABLE(tmp); |
| c.v = pmadd(a.v, b, c.v); |
| #else |
| tmp = b; |
| tmp = pmul(a.v, tmp); |
| c.v = padd(c.v, tmp); |
| #endif |
| } |
| |
| EIGEN_STRONG_INLINE void madd_impl(const LhsScalar& a, const RhsScalar& b, ResScalar& c, RhsScalar& /*tmp*/, |
| const false_type&) const { |
| c += a * b; |
| } |
| |
| template <typename LhsPacketType, typename AccPacketType, typename LaneIdType> |
| EIGEN_STRONG_INLINE void madd(const LhsPacketType& a, const RhsPacketx4& b, AccPacketType& c, RhsPacket& tmp, |
| const LaneIdType& lane) const { |
| madd(a, b.get(lane), c, tmp, lane); |
| } |
| |
| template <typename ResPacketType, typename AccPacketType> |
| EIGEN_STRONG_INLINE void acc(const AccPacketType& c, const ResPacketType& alpha, ResPacketType& r) const { |
| conj_helper<ResPacketType, ResPacketType, ConjLhs, false> cj; |
| r = cj.pmadd(c, alpha, r); |
| } |
| |
| protected: |
| }; |
| |
| template <typename Packet> |
| struct DoublePacket { |
| Packet first; |
| Packet second; |
| }; |
| |
| template <typename Packet> |
| DoublePacket<Packet> padd(const DoublePacket<Packet>& a, const DoublePacket<Packet>& b) { |
| DoublePacket<Packet> res; |
| res.first = padd(a.first, b.first); |
| res.second = padd(a.second, b.second); |
| return res; |
| } |
| |
| // note that for DoublePacket<RealPacket> the "4" in "downto4" |
| // corresponds to the number of complexes, so it means "8" |
| // it terms of real coefficients. |
| |
| template <typename Packet> |
| const DoublePacket<Packet>& predux_half_dowto4(const DoublePacket<Packet>& a, |
| std::enable_if_t<unpacket_traits<Packet>::size <= 8>* = 0) { |
| return a; |
| } |
| |
| template <typename Packet> |
| DoublePacket<typename unpacket_traits<Packet>::half> predux_half_dowto4( |
| const DoublePacket<Packet>& a, std::enable_if_t<unpacket_traits<Packet>::size == 16>* = 0) { |
| // yes, that's pretty hackish :( |
| DoublePacket<typename unpacket_traits<Packet>::half> res; |
| typedef std::complex<typename unpacket_traits<Packet>::type> Cplx; |
| typedef typename packet_traits<Cplx>::type CplxPacket; |
| res.first = predux_half_dowto4(CplxPacket(a.first)).v; |
| res.second = predux_half_dowto4(CplxPacket(a.second)).v; |
| return res; |
| } |
| |
| // same here, "quad" actually means "8" in terms of real coefficients |
| template <typename Scalar, typename RealPacket> |
| void loadQuadToDoublePacket(const Scalar* b, DoublePacket<RealPacket>& dest, |
| std::enable_if_t<unpacket_traits<RealPacket>::size <= 8>* = 0) { |
| dest.first = pset1<RealPacket>(numext::real(*b)); |
| dest.second = pset1<RealPacket>(numext::imag(*b)); |
| } |
| |
| template <typename Scalar, typename RealPacket> |
| void loadQuadToDoublePacket(const Scalar* b, DoublePacket<RealPacket>& dest, |
| std::enable_if_t<unpacket_traits<RealPacket>::size == 16>* = 0) { |
| // yes, that's pretty hackish too :( |
| typedef typename NumTraits<Scalar>::Real RealScalar; |
| RealScalar r[4] = {numext::real(b[0]), numext::real(b[0]), numext::real(b[1]), numext::real(b[1])}; |
| RealScalar i[4] = {numext::imag(b[0]), numext::imag(b[0]), numext::imag(b[1]), numext::imag(b[1])}; |
| dest.first = ploadquad<RealPacket>(r); |
| dest.second = ploadquad<RealPacket>(i); |
| } |
| |
| template <typename Packet> |
| struct unpacket_traits<DoublePacket<Packet> > { |
| typedef DoublePacket<typename unpacket_traits<Packet>::half> half; |
| enum { size = 2 * unpacket_traits<Packet>::size }; |
| }; |
| // template<typename Packet> |
| // DoublePacket<Packet> pmadd(const DoublePacket<Packet> &a, const DoublePacket<Packet> &b) |
| // { |
| // DoublePacket<Packet> res; |
| // res.first = padd(a.first, b.first); |
| // res.second = padd(a.second,b.second); |
| // return res; |
| // } |
| |
| template <typename RealScalar, bool ConjLhs_, bool ConjRhs_, int Arch, int PacketSize_> |
| class gebp_traits<std::complex<RealScalar>, std::complex<RealScalar>, ConjLhs_, ConjRhs_, Arch, PacketSize_> { |
| public: |
| typedef std::complex<RealScalar> Scalar; |
| typedef std::complex<RealScalar> LhsScalar; |
| typedef std::complex<RealScalar> RhsScalar; |
| typedef std::complex<RealScalar> ResScalar; |
| |
| PACKET_DECL_COND_POSTFIX(_, Lhs, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Rhs, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Res, PacketSize_); |
| PACKET_DECL_COND(Real, PacketSize_); |
| PACKET_DECL_COND_SCALAR(PacketSize_); |
| |
| enum { |
| ConjLhs = ConjLhs_, |
| ConjRhs = ConjRhs_, |
| Vectorizable = unpacket_traits<RealPacket>::vectorizable && unpacket_traits<ScalarPacket>::vectorizable, |
| ResPacketSize = Vectorizable ? unpacket_traits<ResPacket_>::size : 1, |
| LhsPacketSize = Vectorizable ? unpacket_traits<LhsPacket_>::size : 1, |
| RhsPacketSize = Vectorizable ? unpacket_traits<RhsScalar>::size : 1, |
| RealPacketSize = Vectorizable ? unpacket_traits<RealPacket>::size : 1, |
| |
| // FIXME: should depend on NumberOfRegisters |
| nr = 4, |
| mr = ResPacketSize, |
| |
| LhsProgress = ResPacketSize, |
| RhsProgress = 1 |
| }; |
| |
| typedef DoublePacket<RealPacket> DoublePacketType; |
| |
| typedef std::conditional_t<Vectorizable, ScalarPacket, Scalar> LhsPacket4Packing; |
| typedef std::conditional_t<Vectorizable, RealPacket, Scalar> LhsPacket; |
| typedef std::conditional_t<Vectorizable, DoublePacketType, Scalar> RhsPacket; |
| typedef std::conditional_t<Vectorizable, ScalarPacket, Scalar> ResPacket; |
| typedef std::conditional_t<Vectorizable, DoublePacketType, Scalar> AccPacket; |
| |
| // this actually holds 8 packets! |
| typedef QuadPacket<RhsPacket> RhsPacketx4; |
| |
| EIGEN_STRONG_INLINE void initAcc(Scalar& p) { p = Scalar(0); } |
| |
| EIGEN_STRONG_INLINE void initAcc(DoublePacketType& p) { |
| p.first = pset1<RealPacket>(RealScalar(0)); |
| p.second = pset1<RealPacket>(RealScalar(0)); |
| } |
| |
| // Scalar path |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, ScalarPacket& dest) const { dest = pset1<ScalarPacket>(*b); } |
| |
| // Vectorized path |
| template <typename RealPacketType> |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, DoublePacket<RealPacketType>& dest) const { |
| dest.first = pset1<RealPacketType>(numext::real(*b)); |
| dest.second = pset1<RealPacketType>(numext::imag(*b)); |
| } |
| |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, RhsPacketx4& dest) const { |
| loadRhs(b, dest.B_0); |
| loadRhs(b + 1, dest.B1); |
| loadRhs(b + 2, dest.B2); |
| loadRhs(b + 3, dest.B3); |
| } |
| |
| // Scalar path |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar* b, ScalarPacket& dest) const { loadRhs(b, dest); } |
| |
| // Vectorized path |
| template <typename RealPacketType> |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar* b, DoublePacket<RealPacketType>& dest) const { |
| loadRhs(b, dest); |
| } |
| |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar*, RhsPacketx4&) const {} |
| |
| EIGEN_STRONG_INLINE void loadRhsQuad(const RhsScalar* b, ResPacket& dest) const { loadRhs(b, dest); } |
| EIGEN_STRONG_INLINE void loadRhsQuad(const RhsScalar* b, DoublePacketType& dest) const { |
| loadQuadToDoublePacket(b, dest); |
| } |
| |
| // nothing special here |
| EIGEN_STRONG_INLINE void loadLhs(const LhsScalar* a, LhsPacket& dest) const { |
| dest = pload<LhsPacket>((const typename unpacket_traits<LhsPacket>::type*)(a)); |
| } |
| |
| template <typename LhsPacketType> |
| EIGEN_STRONG_INLINE void loadLhsUnaligned(const LhsScalar* a, LhsPacketType& dest) const { |
| dest = ploadu<LhsPacketType>((const typename unpacket_traits<LhsPacketType>::type*)(a)); |
| } |
| |
| template <typename LhsPacketType, typename RhsPacketType, typename ResPacketType, typename TmpType, |
| typename LaneIdType> |
| EIGEN_STRONG_INLINE std::enable_if_t<!is_same<RhsPacketType, RhsPacketx4>::value> madd(const LhsPacketType& a, |
| const RhsPacketType& b, |
| DoublePacket<ResPacketType>& c, |
| TmpType& /*tmp*/, |
| const LaneIdType&) const { |
| c.first = padd(pmul(a, b.first), c.first); |
| c.second = padd(pmul(a, b.second), c.second); |
| } |
| |
| template <typename LaneIdType> |
| EIGEN_STRONG_INLINE void madd(const LhsPacket& a, const RhsPacket& b, ResPacket& c, RhsPacket& /*tmp*/, |
| const LaneIdType&) const { |
| c = cj.pmadd(a, b, c); |
| } |
| |
| template <typename LhsPacketType, typename AccPacketType, typename LaneIdType> |
| EIGEN_STRONG_INLINE void madd(const LhsPacketType& a, const RhsPacketx4& b, AccPacketType& c, RhsPacket& tmp, |
| const LaneIdType& lane) const { |
| madd(a, b.get(lane), c, tmp, lane); |
| } |
| |
| EIGEN_STRONG_INLINE void acc(const Scalar& c, const Scalar& alpha, Scalar& r) const { r += alpha * c; } |
| |
| template <typename RealPacketType, typename ResPacketType> |
| EIGEN_STRONG_INLINE void acc(const DoublePacket<RealPacketType>& c, const ResPacketType& alpha, |
| ResPacketType& r) const { |
| // assemble c |
| ResPacketType tmp; |
| if ((!ConjLhs) && (!ConjRhs)) { |
| tmp = pcplxflip(pconj(ResPacketType(c.second))); |
| tmp = padd(ResPacketType(c.first), tmp); |
| } else if ((!ConjLhs) && (ConjRhs)) { |
| tmp = pconj(pcplxflip(ResPacketType(c.second))); |
| tmp = padd(ResPacketType(c.first), tmp); |
| } else if ((ConjLhs) && (!ConjRhs)) { |
| tmp = pcplxflip(ResPacketType(c.second)); |
| tmp = padd(pconj(ResPacketType(c.first)), tmp); |
| } else if ((ConjLhs) && (ConjRhs)) { |
| tmp = pcplxflip(ResPacketType(c.second)); |
| tmp = psub(pconj(ResPacketType(c.first)), tmp); |
| } |
| |
| r = pmadd(tmp, alpha, r); |
| } |
| |
| protected: |
| conj_helper<LhsScalar, RhsScalar, ConjLhs, ConjRhs> cj; |
| }; |
| |
| template <typename RealScalar, bool ConjRhs_, int Arch, int PacketSize_> |
| class gebp_traits<RealScalar, std::complex<RealScalar>, false, ConjRhs_, Arch, PacketSize_> { |
| public: |
| typedef std::complex<RealScalar> Scalar; |
| typedef RealScalar LhsScalar; |
| typedef Scalar RhsScalar; |
| typedef Scalar ResScalar; |
| |
| PACKET_DECL_COND_POSTFIX(_, Lhs, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Rhs, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Res, PacketSize_); |
| PACKET_DECL_COND_POSTFIX(_, Real, PacketSize_); |
| PACKET_DECL_COND_SCALAR_POSTFIX(_, PacketSize_); |
| |
| #undef PACKET_DECL_COND_SCALAR_POSTFIX |
| #undef PACKET_DECL_COND_POSTFIX |
| #undef PACKET_DECL_COND_SCALAR |
| #undef PACKET_DECL_COND |
| |
| enum { |
| ConjLhs = false, |
| ConjRhs = ConjRhs_, |
| Vectorizable = unpacket_traits<RealPacket_>::vectorizable && unpacket_traits<ScalarPacket_>::vectorizable, |
| LhsPacketSize = Vectorizable ? unpacket_traits<LhsPacket_>::size : 1, |
| RhsPacketSize = Vectorizable ? unpacket_traits<RhsPacket_>::size : 1, |
| ResPacketSize = Vectorizable ? unpacket_traits<ResPacket_>::size : 1, |
| |
| NumberOfRegisters = EIGEN_ARCH_DEFAULT_NUMBER_OF_REGISTERS, |
| // FIXME: should depend on NumberOfRegisters |
| nr = 4, |
| mr = (plain_enum_min(16, NumberOfRegisters) / 2 / nr) * ResPacketSize, |
| |
| LhsProgress = ResPacketSize, |
| RhsProgress = 1 |
| }; |
| |
| typedef std::conditional_t<Vectorizable, LhsPacket_, LhsScalar> LhsPacket; |
| typedef std::conditional_t<Vectorizable, RhsPacket_, RhsScalar> RhsPacket; |
| typedef std::conditional_t<Vectorizable, ResPacket_, ResScalar> ResPacket; |
| typedef LhsPacket LhsPacket4Packing; |
| typedef QuadPacket<RhsPacket> RhsPacketx4; |
| typedef ResPacket AccPacket; |
| |
| EIGEN_STRONG_INLINE void initAcc(AccPacket& p) { p = pset1<ResPacket>(ResScalar(0)); } |
| |
| template <typename RhsPacketType> |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, RhsPacketType& dest) const { |
| dest = pset1<RhsPacketType>(*b); |
| } |
| |
| EIGEN_STRONG_INLINE void loadRhs(const RhsScalar* b, RhsPacketx4& dest) const { |
| pbroadcast4(b, dest.B_0, dest.B1, dest.B2, dest.B3); |
| } |
| |
| template <typename RhsPacketType> |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar* b, RhsPacketType& dest) const { |
| loadRhs(b, dest); |
| } |
| |
| EIGEN_STRONG_INLINE void updateRhs(const RhsScalar*, RhsPacketx4&) const {} |
| |
| EIGEN_STRONG_INLINE void loadLhs(const LhsScalar* a, LhsPacket& dest) const { dest = ploaddup<LhsPacket>(a); } |
| |
| EIGEN_STRONG_INLINE void loadRhsQuad(const RhsScalar* b, RhsPacket& dest) const { dest = ploadquad<RhsPacket>(b); } |
| |
| template <typename LhsPacketType> |
| EIGEN_STRONG_INLINE void loadLhsUnaligned(const LhsScalar* a, LhsPacketType& dest) const { |
| dest = ploaddup<LhsPacketType>(a); |
| } |
| |
| template <typename LhsPacketType, typename RhsPacketType, typename AccPacketType, typename LaneIdType> |
| EIGEN_STRONG_INLINE void madd(const LhsPacketType& a, const RhsPacketType& b, AccPacketType& c, RhsPacketType& tmp, |
| const LaneIdType&) const { |
| madd_impl(a, b, c, tmp, std::conditional_t<Vectorizable, true_type, false_type>()); |
| } |
| |
| template <typename LhsPacketType, typename RhsPacketType, typename AccPacketType> |
| EIGEN_STRONG_INLINE void madd_impl(const LhsPacketType& a, const RhsPacketType& b, AccPacketType& c, |
| RhsPacketType& tmp, const true_type&) const { |
| #ifdef EIGEN_HAS_SINGLE_INSTRUCTION_MADD |
| EIGEN_UNUSED_VARIABLE(tmp); |
| c.v = pmadd(a, b.v, c.v); |
| #else |
| tmp = b; |
| tmp.v = pmul(a, tmp.v); |
| c = padd(c, tmp); |
| #endif |
| } |
| |
| EIGEN_STRONG_INLINE void madd_impl(const LhsScalar& a, const RhsScalar& b, ResScalar& c, RhsScalar& /*tmp*/, |
| const false_type&) const { |
| c += a * b; |
| } |
| |
| template <typename LhsPacketType, typename AccPacketType, typename LaneIdType> |
| EIGEN_STRONG_INLINE void madd(const LhsPacketType& a, const RhsPacketx4& b, AccPacketType& c, RhsPacket& tmp, |
| const LaneIdType& lane) const { |
| madd(a, b.get(lane), c, tmp, lane); |
| } |
| |
| template <typename ResPacketType, typename AccPacketType> |
| EIGEN_STRONG_INLINE void acc(const AccPacketType& c, const ResPacketType& alpha, ResPacketType& r) const { |
| conj_helper<ResPacketType, ResPacketType, false, ConjRhs> cj; |
| r = cj.pmadd(alpha, c, r); |
| } |
| |
| protected: |
| }; |
| |
| /* optimized General packed Block * packed Panel product kernel |
| * |
| * Mixing type logic: C += A * B |
| * | A | B | comments |
| * |real |cplx | no vectorization yet, would require to pack A with duplication |
| * |cplx |real | easy vectorization |
| */ |
| template <typename LhsScalar, typename RhsScalar, typename Index, typename DataMapper, int mr, int nr, |
| bool ConjugateLhs, bool ConjugateRhs> |
| struct gebp_kernel { |
| typedef gebp_traits<LhsScalar, RhsScalar, ConjugateLhs, ConjugateRhs, Architecture::Target> Traits; |
| typedef gebp_traits<LhsScalar, RhsScalar, ConjugateLhs, ConjugateRhs, Architecture::Target, GEBPPacketHalf> |
| HalfTraits; |
| typedef gebp_traits<LhsScalar, RhsScalar, ConjugateLhs, ConjugateRhs, Architecture::Target, GEBPPacketQuarter> |
| QuarterTraits; |
| |
| typedef typename Traits::ResScalar ResScalar; |
| typedef typename Traits::LhsPacket LhsPacket; |
| typedef typename Traits::RhsPacket RhsPacket; |
| typedef typename Traits::ResPacket ResPacket; |
| typedef typename Traits::AccPacket AccPacket; |
| typedef typename Traits::RhsPacketx4 RhsPacketx4; |
| |
| typedef typename RhsPanelHelper<RhsPacket, RhsPacketx4, 15>::type RhsPanel15; |
| typedef typename RhsPanelHelper<RhsPacket, RhsPacketx4, 27>::type RhsPanel27; |
| |
| typedef gebp_traits<RhsScalar, LhsScalar, ConjugateRhs, ConjugateLhs, Architecture::Target> SwappedTraits; |
| |
| typedef typename SwappedTraits::ResScalar SResScalar; |
| typedef typename SwappedTraits::LhsPacket SLhsPacket; |
| typedef typename SwappedTraits::RhsPacket SRhsPacket; |
| typedef typename SwappedTraits::ResPacket SResPacket; |
| typedef typename SwappedTraits::AccPacket SAccPacket; |
| |
| typedef typename HalfTraits::LhsPacket LhsPacketHalf; |
| typedef typename HalfTraits::RhsPacket RhsPacketHalf; |
| typedef typename HalfTraits::ResPacket ResPacketHalf; |
| typedef typename HalfTraits::AccPacket AccPacketHalf; |
| |
| typedef typename QuarterTraits::LhsPacket LhsPacketQuarter; |
| typedef typename QuarterTraits::RhsPacket RhsPacketQuarter; |
| typedef typename QuarterTraits::ResPacket ResPacketQuarter; |
| typedef typename QuarterTraits::AccPacket AccPacketQuarter; |
| |
| typedef typename DataMapper::LinearMapper LinearMapper; |
| |
| enum { |
| Vectorizable = Traits::Vectorizable, |
| LhsProgress = Traits::LhsProgress, |
| LhsProgressHalf = HalfTraits::LhsProgress, |
| LhsProgressQuarter = QuarterTraits::LhsProgress, |
| RhsProgress = Traits::RhsProgress, |
| RhsProgressHalf = HalfTraits::RhsProgress, |
| RhsProgressQuarter = QuarterTraits::RhsProgress, |
| ResPacketSize = Traits::ResPacketSize |
| }; |
| |
| EIGEN_DONT_INLINE void operator()(const DataMapper& res, const LhsScalar* blockA, const RhsScalar* blockB, Index rows, |
| Index depth, Index cols, ResScalar alpha, Index strideA = -1, Index strideB = -1, |
| Index offsetA = 0, Index offsetB = 0); |
| }; |
| |
| template <typename LhsScalar, typename RhsScalar, typename Index, typename DataMapper, int mr, int nr, |
| bool ConjugateLhs, bool ConjugateRhs, |
| int SwappedLhsProgress = |
| gebp_traits<RhsScalar, LhsScalar, ConjugateRhs, ConjugateLhs, Architecture::Target>::LhsProgress> |
| struct last_row_process_16_packets { |
| typedef gebp_traits<LhsScalar, RhsScalar, ConjugateLhs, ConjugateRhs, Architecture::Target> Traits; |
| typedef gebp_traits<RhsScalar, LhsScalar, ConjugateRhs, ConjugateLhs, Architecture::Target> SwappedTraits; |
| |
| typedef typename Traits::ResScalar ResScalar; |
| typedef typename SwappedTraits::LhsPacket SLhsPacket; |
| typedef typename SwappedTraits::RhsPacket SRhsPacket; |
| typedef typename SwappedTraits::ResPacket SResPacket; |
| typedef typename SwappedTraits::AccPacket SAccPacket; |
| |
| EIGEN_STRONG_INLINE void operator()(const DataMapper& res, SwappedTraits& straits, const LhsScalar* blA, |
| const RhsScalar* blB, Index depth, const Index endk, Index i, Index j2, |
| ResScalar alpha, SAccPacket& C0) { |
| EIGEN_UNUSED_VARIABLE(res); |
| EIGEN_UNUSED_VARIABLE(straits); |
| EIGEN_UNUSED_VARIABLE(blA); |
| EIGEN_UNUSED_VARIABLE(blB); |
| EIGEN_UNUSED_VARIABLE(depth); |
| EIGEN_UNUSED_VARIABLE(endk); |
| EIGEN_UNUSED_VARIABLE(i); |
| EIGEN_UNUSED_VARIABLE(j2); |
| EIGEN_UNUSED_VARIABLE(alpha); |
| EIGEN_UNUSED_VARIABLE(C0); |
| } |
| }; |
| |
| template <typename LhsScalar, typename RhsScalar, typename Index, typename DataMapper, int mr, int nr, |
| bool ConjugateLhs, bool ConjugateRhs> |
| struct last_row_process_16_packets<LhsScalar, RhsScalar, Index, DataMapper, mr, nr, ConjugateLhs, ConjugateRhs, 16> { |
| typedef gebp_traits<LhsScalar, RhsScalar, ConjugateLhs, ConjugateRhs, Architecture::Target> Traits; |
| typedef gebp_traits<RhsScalar, LhsScalar, ConjugateRhs, ConjugateLhs, Architecture::Target> SwappedTraits; |
| |
| typedef typename Traits::ResScalar ResScalar; |
| typedef typename SwappedTraits::LhsPacket SLhsPacket; |
| typedef typename SwappedTraits::RhsPacket SRhsPacket; |
| typedef typename SwappedTraits::ResPacket SResPacket; |
| typedef typename SwappedTraits::AccPacket SAccPacket; |
| |
| EIGEN_STRONG_INLINE void operator()(const DataMapper& res, SwappedTraits& straits, const LhsScalar* blA, |
| const RhsScalar* blB, Index depth, const Index endk, Index i, Index j2, |
| ResScalar alpha, SAccPacket& C0) { |
| typedef typename unpacket_traits<typename unpacket_traits<SResPacket>::half>::half SResPacketQuarter; |
| typedef typename unpacket_traits<typename unpacket_traits<SLhsPacket>::half>::half SLhsPacketQuarter; |
| typedef typename unpacket_traits<typename unpacket_traits<SRhsPacket>::half>::half SRhsPacketQuarter; |
| typedef typename unpacket_traits<typename unpacket_traits<SAccPacket>::half>::half SAccPacketQuarter; |
| |
| SResPacketQuarter R = res.template gatherPacket<SResPacketQuarter>(i, j2); |
| SResPacketQuarter alphav = pset1<SResPacketQuarter>(alpha); |
| |
| if (depth - endk > 0) { |
| // We have to handle the last row(s) of the rhs, which |
| // correspond to a half-packet |
| SAccPacketQuarter c0 = predux_half_dowto4(predux_half_dowto4(C0)); |
| |
| for (Index kk = endk; kk < depth; kk++) { |
| SLhsPacketQuarter a0; |
| SRhsPacketQuarter b0; |
| straits.loadLhsUnaligned(blB, a0); |
| straits.loadRhs(blA, b0); |
| straits.madd(a0, b0, c0, b0, fix<0>); |
| blB += SwappedTraits::LhsProgress / 4; |
| blA += 1; |
| } |
| straits.acc(c0, alphav, R); |
| } else { |
| straits.acc(predux_half_dowto4(predux_half_dowto4(C0)), alphav, R); |
| } |
| res.scatterPacket(i, j2, R); |
| } |
| }; |
| |
| template <int nr, Index LhsProgress, Index RhsProgress, typename LhsScalar, typename RhsScalar, typename ResScalar, |
| typename AccPacket, typename LhsPacket, typename RhsPacket, typename ResPacket, typename GEBPTraits, |
| typename LinearMapper, typename DataMapper> |
| struct lhs_process_one_packet { |
| typedef typename GEBPTraits::RhsPacketx4 RhsPacketx4; |
| |
| EIGEN_STRONG_INLINE void peeled_kc_onestep(Index K, const LhsScalar* blA, const RhsScalar* blB, GEBPTraits traits, |
| LhsPacket* A0, RhsPacketx4* rhs_panel, RhsPacket* T0, AccPacket* C0, |
| AccPacket* C1, AccPacket* C2, AccPacket* C3) { |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 1X4"); |
| EIGEN_ASM_COMMENT("Note: these asm comments work around bug 935!"); |
| traits.loadLhs(&blA[(0 + 1 * K) * LhsProgress], *A0); |
| traits.loadRhs(&blB[(0 + 4 * K) * RhsProgress], *rhs_panel); |
| traits.madd(*A0, *rhs_panel, *C0, *T0, fix<0>); |
| traits.madd(*A0, *rhs_panel, *C1, *T0, fix<1>); |
| traits.madd(*A0, *rhs_panel, *C2, *T0, fix<2>); |
| traits.madd(*A0, *rhs_panel, *C3, *T0, fix<3>); |
| #if EIGEN_GNUC_STRICT_AT_LEAST(6, 0, 0) && defined(EIGEN_VECTORIZE_SSE) && !(EIGEN_COMP_LCC) |
| __asm__("" : "+x,m"(*A0)); |
| #endif |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 1X4"); |
| } |
| |
| EIGEN_STRONG_INLINE void operator()(const DataMapper& res, const LhsScalar* blockA, const RhsScalar* blockB, |
| ResScalar alpha, Index peelStart, Index peelEnd, Index strideA, Index strideB, |
| Index offsetA, Index offsetB, int prefetch_res_offset, Index peeled_kc, Index pk, |
| Index cols, Index depth, Index packet_cols4) { |
| GEBPTraits traits; |
| Index packet_cols8 = nr >= 8 ? (cols / 8) * 8 : 0; |
| // loops on each largest micro horizontal panel of lhs |
| // (LhsProgress x depth) |
| for (Index i = peelStart; i < peelEnd; i += LhsProgress) { |
| #if EIGEN_ARCH_ARM64 |
| EIGEN_IF_CONSTEXPR(nr >= 8) { |
| for (Index j2 = 0; j2 < packet_cols8; j2 += 8) { |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (LhsProgress)]; |
| prefetch(&blA[0]); |
| |
| // gets res block as register |
| AccPacket C0, C1, C2, C3, C4, C5, C6, C7; |
| traits.initAcc(C0); |
| traits.initAcc(C1); |
| traits.initAcc(C2); |
| traits.initAcc(C3); |
| traits.initAcc(C4); |
| traits.initAcc(C5); |
| traits.initAcc(C6); |
| traits.initAcc(C7); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2 + 0); |
| LinearMapper r1 = res.getLinearMapper(i, j2 + 1); |
| LinearMapper r2 = res.getLinearMapper(i, j2 + 2); |
| LinearMapper r3 = res.getLinearMapper(i, j2 + 3); |
| LinearMapper r4 = res.getLinearMapper(i, j2 + 4); |
| LinearMapper r5 = res.getLinearMapper(i, j2 + 5); |
| LinearMapper r6 = res.getLinearMapper(i, j2 + 6); |
| LinearMapper r7 = res.getLinearMapper(i, j2 + 7); |
| r0.prefetch(prefetch_res_offset); |
| r1.prefetch(prefetch_res_offset); |
| r2.prefetch(prefetch_res_offset); |
| r3.prefetch(prefetch_res_offset); |
| r4.prefetch(prefetch_res_offset); |
| r5.prefetch(prefetch_res_offset); |
| r6.prefetch(prefetch_res_offset); |
| r7.prefetch(prefetch_res_offset); |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB * 8]; |
| prefetch(&blB[0]); |
| |
| LhsPacket A0; |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| RhsPacketx4 rhs_panel; |
| RhsPacket T0; |
| #define EIGEN_GEBGP_ONESTEP(K) \ |
| do { \ |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 1pX8"); \ |
| traits.loadLhs(&blA[(0 + 1 * K) * LhsProgress], A0); \ |
| traits.loadRhs(&blB[(0 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C0, T0, fix<0>); \ |
| traits.updateRhs(&blB[(1 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C1, T0, fix<1>); \ |
| traits.updateRhs(&blB[(2 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C2, T0, fix<2>); \ |
| traits.updateRhs(&blB[(3 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C3, T0, fix<3>); \ |
| traits.loadRhs(&blB[(4 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C4, T0, fix<0>); \ |
| traits.updateRhs(&blB[(5 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C5, T0, fix<1>); \ |
| traits.updateRhs(&blB[(6 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C6, T0, fix<2>); \ |
| traits.updateRhs(&blB[(7 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C7, T0, fix<3>); \ |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 1pX8"); \ |
| } while (false) |
| |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 1pX8"); |
| |
| EIGEN_GEBGP_ONESTEP(0); |
| EIGEN_GEBGP_ONESTEP(1); |
| EIGEN_GEBGP_ONESTEP(2); |
| EIGEN_GEBGP_ONESTEP(3); |
| EIGEN_GEBGP_ONESTEP(4); |
| EIGEN_GEBGP_ONESTEP(5); |
| EIGEN_GEBGP_ONESTEP(6); |
| EIGEN_GEBGP_ONESTEP(7); |
| |
| blB += pk * 8 * RhsProgress; |
| blA += pk * (1 * LhsProgress); |
| |
| EIGEN_ASM_COMMENT("end gebp micro kernel 1pX8"); |
| } |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPacketx4 rhs_panel; |
| RhsPacket T0; |
| EIGEN_GEBGP_ONESTEP(0); |
| blB += 8 * RhsProgress; |
| blA += 1 * LhsProgress; |
| } |
| |
| #undef EIGEN_GEBGP_ONESTEP |
| |
| ResPacket R0, R1; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| |
| R0 = r0.template loadPacket<ResPacket>(0); |
| R1 = r1.template loadPacket<ResPacket>(0); |
| traits.acc(C0, alphav, R0); |
| traits.acc(C1, alphav, R1); |
| r0.storePacket(0, R0); |
| r1.storePacket(0, R1); |
| |
| R0 = r2.template loadPacket<ResPacket>(0); |
| R1 = r3.template loadPacket<ResPacket>(0); |
| traits.acc(C2, alphav, R0); |
| traits.acc(C3, alphav, R1); |
| r2.storePacket(0, R0); |
| r3.storePacket(0, R1); |
| |
| R0 = r4.template loadPacket<ResPacket>(0); |
| R1 = r5.template loadPacket<ResPacket>(0); |
| traits.acc(C4, alphav, R0); |
| traits.acc(C5, alphav, R1); |
| r4.storePacket(0, R0); |
| r5.storePacket(0, R1); |
| |
| R0 = r6.template loadPacket<ResPacket>(0); |
| R1 = r7.template loadPacket<ResPacket>(0); |
| traits.acc(C6, alphav, R0); |
| traits.acc(C7, alphav, R1); |
| r6.storePacket(0, R0); |
| r7.storePacket(0, R1); |
| } |
| } |
| #endif |
| |
| // loops on each largest micro vertical panel of rhs (depth * nr) |
| for (Index j2 = packet_cols8; j2 < packet_cols4; j2 += 4) { |
| // We select a LhsProgress x nr micro block of res |
| // which is entirely stored into 1 x nr registers. |
| |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (LhsProgress)]; |
| prefetch(&blA[0]); |
| |
| // gets res block as register |
| AccPacket C0, C1, C2, C3; |
| traits.initAcc(C0); |
| traits.initAcc(C1); |
| traits.initAcc(C2); |
| traits.initAcc(C3); |
| // To improve instruction pipelining, let's double the accumulation registers: |
| // even k will accumulate in C*, while odd k will accumulate in D*. |
| // This trick is crutial to get good performance with FMA, otherwise it is |
| // actually faster to perform separated MUL+ADD because of a naturally |
| // better instruction-level parallelism. |
| AccPacket D0, D1, D2, D3; |
| traits.initAcc(D0); |
| traits.initAcc(D1); |
| traits.initAcc(D2); |
| traits.initAcc(D3); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2 + 0); |
| LinearMapper r1 = res.getLinearMapper(i, j2 + 1); |
| LinearMapper r2 = res.getLinearMapper(i, j2 + 2); |
| LinearMapper r3 = res.getLinearMapper(i, j2 + 3); |
| |
| r0.prefetch(prefetch_res_offset); |
| r1.prefetch(prefetch_res_offset); |
| r2.prefetch(prefetch_res_offset); |
| r3.prefetch(prefetch_res_offset); |
| |
| // performs "inner" products |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB * 4]; |
| prefetch(&blB[0]); |
| LhsPacket A0, A1; |
| |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 1/half/quarterX4"); |
| RhsPacketx4 rhs_panel; |
| RhsPacket T0; |
| |
| internal::prefetch(blB + (48 + 0)); |
| peeled_kc_onestep(0, blA, blB, traits, &A0, &rhs_panel, &T0, &C0, &C1, &C2, &C3); |
| peeled_kc_onestep(1, blA, blB, traits, &A1, &rhs_panel, &T0, &D0, &D1, &D2, &D3); |
| peeled_kc_onestep(2, blA, blB, traits, &A0, &rhs_panel, &T0, &C0, &C1, &C2, &C3); |
| peeled_kc_onestep(3, blA, blB, traits, &A1, &rhs_panel, &T0, &D0, &D1, &D2, &D3); |
| internal::prefetch(blB + (48 + 16)); |
| peeled_kc_onestep(4, blA, blB, traits, &A0, &rhs_panel, &T0, &C0, &C1, &C2, &C3); |
| peeled_kc_onestep(5, blA, blB, traits, &A1, &rhs_panel, &T0, &D0, &D1, &D2, &D3); |
| peeled_kc_onestep(6, blA, blB, traits, &A0, &rhs_panel, &T0, &C0, &C1, &C2, &C3); |
| peeled_kc_onestep(7, blA, blB, traits, &A1, &rhs_panel, &T0, &D0, &D1, &D2, &D3); |
| |
| blB += pk * 4 * RhsProgress; |
| blA += pk * LhsProgress; |
| |
| EIGEN_ASM_COMMENT("end gebp micro kernel 1/half/quarterX4"); |
| } |
| C0 = padd(C0, D0); |
| C1 = padd(C1, D1); |
| C2 = padd(C2, D2); |
| C3 = padd(C3, D3); |
| |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPacketx4 rhs_panel; |
| RhsPacket T0; |
| peeled_kc_onestep(0, blA, blB, traits, &A0, &rhs_panel, &T0, &C0, &C1, &C2, &C3); |
| blB += 4 * RhsProgress; |
| blA += LhsProgress; |
| } |
| |
| ResPacket R0, R1; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| |
| R0 = r0.template loadPacket<ResPacket>(0); |
| R1 = r1.template loadPacket<ResPacket>(0); |
| traits.acc(C0, alphav, R0); |
| traits.acc(C1, alphav, R1); |
| r0.storePacket(0, R0); |
| r1.storePacket(0, R1); |
| |
| R0 = r2.template loadPacket<ResPacket>(0); |
| R1 = r3.template loadPacket<ResPacket>(0); |
| traits.acc(C2, alphav, R0); |
| traits.acc(C3, alphav, R1); |
| r2.storePacket(0, R0); |
| r3.storePacket(0, R1); |
| } |
| |
| // Deal with remaining columns of the rhs |
| for (Index j2 = packet_cols4; j2 < cols; j2++) { |
| // One column at a time |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (LhsProgress)]; |
| prefetch(&blA[0]); |
| |
| // gets res block as register |
| AccPacket C0; |
| traits.initAcc(C0); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2); |
| |
| // performs "inner" products |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB]; |
| LhsPacket A0; |
| |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 1/half/quarterX1"); |
| RhsPacket B_0; |
| |
| #define EIGEN_GEBGP_ONESTEP(K) \ |
| do { \ |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 1/half/quarterX1"); \ |
| EIGEN_ASM_COMMENT("Note: these asm comments work around bug 935!"); \ |
| /* FIXME: why unaligned???? */ \ |
| traits.loadLhsUnaligned(&blA[(0 + 1 * K) * LhsProgress], A0); \ |
| traits.loadRhs(&blB[(0 + K) * RhsProgress], B_0); \ |
| traits.madd(A0, B_0, C0, B_0, fix<0>); \ |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 1/half/quarterX1"); \ |
| } while (false); |
| |
| EIGEN_GEBGP_ONESTEP(0); |
| EIGEN_GEBGP_ONESTEP(1); |
| EIGEN_GEBGP_ONESTEP(2); |
| EIGEN_GEBGP_ONESTEP(3); |
| EIGEN_GEBGP_ONESTEP(4); |
| EIGEN_GEBGP_ONESTEP(5); |
| EIGEN_GEBGP_ONESTEP(6); |
| EIGEN_GEBGP_ONESTEP(7); |
| |
| blB += pk * RhsProgress; |
| blA += pk * LhsProgress; |
| |
| EIGEN_ASM_COMMENT("end gebp micro kernel 1/half/quarterX1"); |
| } |
| |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPacket B_0; |
| EIGEN_GEBGP_ONESTEP(0); |
| blB += RhsProgress; |
| blA += LhsProgress; |
| } |
| #undef EIGEN_GEBGP_ONESTEP |
| ResPacket R0; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| R0 = r0.template loadPacket<ResPacket>(0); |
| traits.acc(C0, alphav, R0); |
| r0.storePacket(0, R0); |
| } |
| } |
| } |
| }; |
| |
| template <int nr, Index LhsProgress, Index RhsProgress, typename LhsScalar, typename RhsScalar, typename ResScalar, |
| typename AccPacket, typename LhsPacket, typename RhsPacket, typename ResPacket, typename GEBPTraits, |
| typename LinearMapper, typename DataMapper> |
| struct lhs_process_fraction_of_packet |
| : lhs_process_one_packet<nr, LhsProgress, RhsProgress, LhsScalar, RhsScalar, ResScalar, AccPacket, LhsPacket, |
| RhsPacket, ResPacket, GEBPTraits, LinearMapper, DataMapper> { |
| EIGEN_STRONG_INLINE void peeled_kc_onestep(Index K, const LhsScalar* blA, const RhsScalar* blB, GEBPTraits traits, |
| LhsPacket* A0, RhsPacket* B_0, RhsPacket* B1, RhsPacket* B2, RhsPacket* B3, |
| AccPacket* C0, AccPacket* C1, AccPacket* C2, AccPacket* C3) { |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 1X4"); |
| EIGEN_ASM_COMMENT("Note: these asm comments work around bug 935!"); |
| traits.loadLhsUnaligned(&blA[(0 + 1 * K) * (LhsProgress)], *A0); |
| traits.broadcastRhs(&blB[(0 + 4 * K) * RhsProgress], *B_0, *B1, *B2, *B3); |
| traits.madd(*A0, *B_0, *C0, *B_0); |
| traits.madd(*A0, *B1, *C1, *B1); |
| traits.madd(*A0, *B2, *C2, *B2); |
| traits.madd(*A0, *B3, *C3, *B3); |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 1X4"); |
| } |
| }; |
| |
| template <typename LhsScalar, typename RhsScalar, typename Index, typename DataMapper, int mr, int nr, |
| bool ConjugateLhs, bool ConjugateRhs> |
| EIGEN_DONT_INLINE void gebp_kernel<LhsScalar, RhsScalar, Index, DataMapper, mr, nr, ConjugateLhs, |
| ConjugateRhs>::operator()(const DataMapper& res, const LhsScalar* blockA, |
| const RhsScalar* blockB, Index rows, Index depth, |
| Index cols, ResScalar alpha, Index strideA, Index strideB, |
| Index offsetA, Index offsetB) { |
| Traits traits; |
| SwappedTraits straits; |
| |
| if (strideA == -1) strideA = depth; |
| if (strideB == -1) strideB = depth; |
| conj_helper<LhsScalar, RhsScalar, ConjugateLhs, ConjugateRhs> cj; |
| Index packet_cols4 = nr >= 4 ? (cols / 4) * 4 : 0; |
| Index packet_cols8 = nr >= 8 ? (cols / 8) * 8 : 0; |
| const Index peeled_mc3 = mr >= 3 * Traits::LhsProgress ? (rows / (3 * LhsProgress)) * (3 * LhsProgress) : 0; |
| const Index peeled_mc2 = |
| mr >= 2 * Traits::LhsProgress ? peeled_mc3 + ((rows - peeled_mc3) / (2 * LhsProgress)) * (2 * LhsProgress) : 0; |
| const Index peeled_mc1 = |
| mr >= 1 * Traits::LhsProgress ? peeled_mc2 + ((rows - peeled_mc2) / (1 * LhsProgress)) * (1 * LhsProgress) : 0; |
| const Index peeled_mc_half = |
| mr >= LhsProgressHalf ? peeled_mc1 + ((rows - peeled_mc1) / (LhsProgressHalf)) * (LhsProgressHalf) : 0; |
| const Index peeled_mc_quarter = |
| mr >= LhsProgressQuarter |
| ? peeled_mc_half + ((rows - peeled_mc_half) / (LhsProgressQuarter)) * (LhsProgressQuarter) |
| : 0; |
| enum { pk = 8 }; // NOTE Such a large peeling factor is important for large matrices (~ +5% when >1000 on Haswell) |
| const Index peeled_kc = depth & ~(pk - 1); |
| const int prefetch_res_offset = 32 / sizeof(ResScalar); |
| // const Index depth2 = depth & ~1; |
| |
| //---------- Process 3 * LhsProgress rows at once ---------- |
| // This corresponds to 3*LhsProgress x nr register blocks. |
| // Usually, make sense only with FMA |
| if (mr >= 3 * Traits::LhsProgress) { |
| // Here, the general idea is to loop on each largest micro horizontal panel of the lhs (3*Traits::LhsProgress x |
| // depth) and on each largest micro vertical panel of the rhs (depth * nr). Blocking sizes, i.e., 'depth' has been |
| // computed so that the micro horizontal panel of the lhs fit in L1. However, if depth is too small, we can extend |
| // the number of rows of these horizontal panels. This actual number of rows is computed as follow: |
| const Index l1 = defaultL1CacheSize; // in Bytes, TODO, l1 should be passed to this function. |
| // The max(1, ...) here is needed because we may be using blocking params larger than what our known l1 cache size |
| // suggests we should be using: either because our known l1 cache size is inaccurate (e.g. on Android, we can only |
| // guess), or because we are testing specific blocking sizes. |
| const Index actual_panel_rows = |
| (3 * LhsProgress) * std::max<Index>(1, ((l1 - sizeof(ResScalar) * mr * nr - depth * nr * sizeof(RhsScalar)) / |
| (depth * sizeof(LhsScalar) * 3 * LhsProgress))); |
| for (Index i1 = 0; i1 < peeled_mc3; i1 += actual_panel_rows) { |
| const Index actual_panel_end = (std::min)(i1 + actual_panel_rows, peeled_mc3); |
| #if EIGEN_ARCH_ARM64 |
| EIGEN_IF_CONSTEXPR(nr >= 8) { |
| for (Index j2 = 0; j2 < packet_cols8; j2 += 8) { |
| for (Index i = i1; i < actual_panel_end; i += 3 * LhsProgress) { |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (3 * LhsProgress)]; |
| prefetch(&blA[0]); |
| // gets res block as register |
| AccPacket C0, C1, C2, C3, C4, C5, C6, C7, C8, C9, C10, C11, C12, C13, C14, C15, C16, C17, C18, C19, C20, |
| C21, C22, C23; |
| traits.initAcc(C0); |
| traits.initAcc(C1); |
| traits.initAcc(C2); |
| traits.initAcc(C3); |
| traits.initAcc(C4); |
| traits.initAcc(C5); |
| traits.initAcc(C6); |
| traits.initAcc(C7); |
| traits.initAcc(C8); |
| traits.initAcc(C9); |
| traits.initAcc(C10); |
| traits.initAcc(C11); |
| traits.initAcc(C12); |
| traits.initAcc(C13); |
| traits.initAcc(C14); |
| traits.initAcc(C15); |
| traits.initAcc(C16); |
| traits.initAcc(C17); |
| traits.initAcc(C18); |
| traits.initAcc(C19); |
| traits.initAcc(C20); |
| traits.initAcc(C21); |
| traits.initAcc(C22); |
| traits.initAcc(C23); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2 + 0); |
| LinearMapper r1 = res.getLinearMapper(i, j2 + 1); |
| LinearMapper r2 = res.getLinearMapper(i, j2 + 2); |
| LinearMapper r3 = res.getLinearMapper(i, j2 + 3); |
| LinearMapper r4 = res.getLinearMapper(i, j2 + 4); |
| LinearMapper r5 = res.getLinearMapper(i, j2 + 5); |
| LinearMapper r6 = res.getLinearMapper(i, j2 + 6); |
| LinearMapper r7 = res.getLinearMapper(i, j2 + 7); |
| |
| r0.prefetch(0); |
| r1.prefetch(0); |
| r2.prefetch(0); |
| r3.prefetch(0); |
| r4.prefetch(0); |
| r5.prefetch(0); |
| r6.prefetch(0); |
| r7.prefetch(0); |
| |
| // performs "inner" products |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB * 8]; |
| prefetch(&blB[0]); |
| LhsPacket A0, A1; |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 3pX8"); |
| // 27 registers are taken (24 for acc, 3 for lhs). |
| RhsPanel27 rhs_panel; |
| RhsPacket T0; |
| LhsPacket A2; |
| #if EIGEN_ARCH_ARM64 && defined(EIGEN_VECTORIZE_NEON) && EIGEN_GNUC_STRICT_LESS_THAN(9, 0, 0) |
| // see http://eigen.tuxfamily.org/bz/show_bug.cgi?id=1633 |
| // without this workaround A0, A1, and A2 are loaded in the same register, |
| // which is not good for pipelining |
| #define EIGEN_GEBP_3Px8_REGISTER_ALLOC_WORKAROUND __asm__("" : "+w,m"(A0), "+w,m"(A1), "+w,m"(A2)); |
| #else |
| #define EIGEN_GEBP_3Px8_REGISTER_ALLOC_WORKAROUND |
| #endif |
| |
| #define EIGEN_GEBP_ONESTEP(K) \ |
| do { \ |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 3pX8"); \ |
| traits.loadLhs(&blA[(0 + 3 * K) * LhsProgress], A0); \ |
| traits.loadLhs(&blA[(1 + 3 * K) * LhsProgress], A1); \ |
| traits.loadLhs(&blA[(2 + 3 * K) * LhsProgress], A2); \ |
| EIGEN_GEBP_3Px8_REGISTER_ALLOC_WORKAROUND traits.loadRhs(blB + (0 + 8 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C0, T0, fix<0>); \ |
| traits.madd(A1, rhs_panel, C8, T0, fix<0>); \ |
| traits.madd(A2, rhs_panel, C16, T0, fix<0>); \ |
| traits.updateRhs(blB + (1 + 8 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C1, T0, fix<1>); \ |
| traits.madd(A1, rhs_panel, C9, T0, fix<1>); \ |
| traits.madd(A2, rhs_panel, C17, T0, fix<1>); \ |
| traits.updateRhs(blB + (2 + 8 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C2, T0, fix<2>); \ |
| traits.madd(A1, rhs_panel, C10, T0, fix<2>); \ |
| traits.madd(A2, rhs_panel, C18, T0, fix<2>); \ |
| traits.updateRhs(blB + (3 + 8 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C3, T0, fix<3>); \ |
| traits.madd(A1, rhs_panel, C11, T0, fix<3>); \ |
| traits.madd(A2, rhs_panel, C19, T0, fix<3>); \ |
| traits.loadRhs(blB + (4 + 8 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C4, T0, fix<0>); \ |
| traits.madd(A1, rhs_panel, C12, T0, fix<0>); \ |
| traits.madd(A2, rhs_panel, C20, T0, fix<0>); \ |
| traits.updateRhs(blB + (5 + 8 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C5, T0, fix<1>); \ |
| traits.madd(A1, rhs_panel, C13, T0, fix<1>); \ |
| traits.madd(A2, rhs_panel, C21, T0, fix<1>); \ |
| traits.updateRhs(blB + (6 + 8 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C6, T0, fix<2>); \ |
| traits.madd(A1, rhs_panel, C14, T0, fix<2>); \ |
| traits.madd(A2, rhs_panel, C22, T0, fix<2>); \ |
| traits.updateRhs(blB + (7 + 8 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C7, T0, fix<3>); \ |
| traits.madd(A1, rhs_panel, C15, T0, fix<3>); \ |
| traits.madd(A2, rhs_panel, C23, T0, fix<3>); \ |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 3pX8"); \ |
| } while (false) |
| |
| EIGEN_GEBP_ONESTEP(0); |
| EIGEN_GEBP_ONESTEP(1); |
| EIGEN_GEBP_ONESTEP(2); |
| EIGEN_GEBP_ONESTEP(3); |
| EIGEN_GEBP_ONESTEP(4); |
| EIGEN_GEBP_ONESTEP(5); |
| EIGEN_GEBP_ONESTEP(6); |
| EIGEN_GEBP_ONESTEP(7); |
| |
| blB += pk * 8 * RhsProgress; |
| blA += pk * 3 * Traits::LhsProgress; |
| EIGEN_ASM_COMMENT("end gebp micro kernel 3pX8"); |
| } |
| |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPanel27 rhs_panel; |
| RhsPacket T0; |
| LhsPacket A2; |
| EIGEN_GEBP_ONESTEP(0); |
| blB += 8 * RhsProgress; |
| blA += 3 * Traits::LhsProgress; |
| } |
| |
| #undef EIGEN_GEBP_ONESTEP |
| |
| ResPacket R0, R1, R2; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| |
| R0 = r0.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r0.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r0.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C0, alphav, R0); |
| traits.acc(C8, alphav, R1); |
| traits.acc(C16, alphav, R2); |
| r0.storePacket(0 * Traits::ResPacketSize, R0); |
| r0.storePacket(1 * Traits::ResPacketSize, R1); |
| r0.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r1.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r1.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r1.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C1, alphav, R0); |
| traits.acc(C9, alphav, R1); |
| traits.acc(C17, alphav, R2); |
| r1.storePacket(0 * Traits::ResPacketSize, R0); |
| r1.storePacket(1 * Traits::ResPacketSize, R1); |
| r1.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r2.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r2.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r2.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C2, alphav, R0); |
| traits.acc(C10, alphav, R1); |
| traits.acc(C18, alphav, R2); |
| r2.storePacket(0 * Traits::ResPacketSize, R0); |
| r2.storePacket(1 * Traits::ResPacketSize, R1); |
| r2.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r3.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r3.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r3.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C3, alphav, R0); |
| traits.acc(C11, alphav, R1); |
| traits.acc(C19, alphav, R2); |
| r3.storePacket(0 * Traits::ResPacketSize, R0); |
| r3.storePacket(1 * Traits::ResPacketSize, R1); |
| r3.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r4.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r4.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r4.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C4, alphav, R0); |
| traits.acc(C12, alphav, R1); |
| traits.acc(C20, alphav, R2); |
| r4.storePacket(0 * Traits::ResPacketSize, R0); |
| r4.storePacket(1 * Traits::ResPacketSize, R1); |
| r4.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r5.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r5.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r5.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C5, alphav, R0); |
| traits.acc(C13, alphav, R1); |
| traits.acc(C21, alphav, R2); |
| r5.storePacket(0 * Traits::ResPacketSize, R0); |
| r5.storePacket(1 * Traits::ResPacketSize, R1); |
| r5.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r6.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r6.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r6.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C6, alphav, R0); |
| traits.acc(C14, alphav, R1); |
| traits.acc(C22, alphav, R2); |
| r6.storePacket(0 * Traits::ResPacketSize, R0); |
| r6.storePacket(1 * Traits::ResPacketSize, R1); |
| r6.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r7.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r7.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r7.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C7, alphav, R0); |
| traits.acc(C15, alphav, R1); |
| traits.acc(C23, alphav, R2); |
| r7.storePacket(0 * Traits::ResPacketSize, R0); |
| r7.storePacket(1 * Traits::ResPacketSize, R1); |
| r7.storePacket(2 * Traits::ResPacketSize, R2); |
| } |
| } |
| } |
| #endif |
| for (Index j2 = packet_cols8; j2 < packet_cols4; j2 += 4) { |
| for (Index i = i1; i < actual_panel_end; i += 3 * LhsProgress) { |
| // We selected a 3*Traits::LhsProgress x nr micro block of res which is entirely |
| // stored into 3 x nr registers. |
| |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (3 * LhsProgress)]; |
| prefetch(&blA[0]); |
| |
| // gets res block as register |
| AccPacket C0, C1, C2, C3, C4, C5, C6, C7, C8, C9, C10, C11; |
| traits.initAcc(C0); |
| traits.initAcc(C1); |
| traits.initAcc(C2); |
| traits.initAcc(C3); |
| traits.initAcc(C4); |
| traits.initAcc(C5); |
| traits.initAcc(C6); |
| traits.initAcc(C7); |
| traits.initAcc(C8); |
| traits.initAcc(C9); |
| traits.initAcc(C10); |
| traits.initAcc(C11); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2 + 0); |
| LinearMapper r1 = res.getLinearMapper(i, j2 + 1); |
| LinearMapper r2 = res.getLinearMapper(i, j2 + 2); |
| LinearMapper r3 = res.getLinearMapper(i, j2 + 3); |
| |
| r0.prefetch(0); |
| r1.prefetch(0); |
| r2.prefetch(0); |
| r3.prefetch(0); |
| |
| // performs "inner" products |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB * 4]; |
| prefetch(&blB[0]); |
| LhsPacket A0, A1; |
| |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 3pX4"); |
| // 15 registers are taken (12 for acc, 3 for lhs). |
| RhsPanel15 rhs_panel; |
| RhsPacket T0; |
| LhsPacket A2; |
| #if EIGEN_ARCH_ARM64 && defined(EIGEN_VECTORIZE_NEON) && EIGEN_GNUC_STRICT_LESS_THAN(9, 0, 0) |
| // see http://eigen.tuxfamily.org/bz/show_bug.cgi?id=1633 |
| // without this workaround A0, A1, and A2 are loaded in the same register, |
| // which is not good for pipelining |
| #define EIGEN_GEBP_3PX4_REGISTER_ALLOC_WORKAROUND __asm__("" : "+w,m"(A0), "+w,m"(A1), "+w,m"(A2)); |
| #else |
| #define EIGEN_GEBP_3PX4_REGISTER_ALLOC_WORKAROUND |
| #endif |
| #define EIGEN_GEBP_ONESTEP(K) \ |
| do { \ |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 3pX4"); \ |
| EIGEN_ASM_COMMENT("Note: these asm comments work around bug 935!"); \ |
| internal::prefetch(blA + (3 * K + 16) * LhsProgress); \ |
| if (EIGEN_ARCH_ARM || EIGEN_ARCH_MIPS) { \ |
| internal::prefetch(blB + (4 * K + 16) * RhsProgress); \ |
| } /* Bug 953 */ \ |
| traits.loadLhs(&blA[(0 + 3 * K) * LhsProgress], A0); \ |
| traits.loadLhs(&blA[(1 + 3 * K) * LhsProgress], A1); \ |
| traits.loadLhs(&blA[(2 + 3 * K) * LhsProgress], A2); \ |
| EIGEN_GEBP_3PX4_REGISTER_ALLOC_WORKAROUND \ |
| traits.loadRhs(blB + (0 + 4 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C0, T0, fix<0>); \ |
| traits.madd(A1, rhs_panel, C4, T0, fix<0>); \ |
| traits.madd(A2, rhs_panel, C8, T0, fix<0>); \ |
| traits.updateRhs(blB + (1 + 4 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C1, T0, fix<1>); \ |
| traits.madd(A1, rhs_panel, C5, T0, fix<1>); \ |
| traits.madd(A2, rhs_panel, C9, T0, fix<1>); \ |
| traits.updateRhs(blB + (2 + 4 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C2, T0, fix<2>); \ |
| traits.madd(A1, rhs_panel, C6, T0, fix<2>); \ |
| traits.madd(A2, rhs_panel, C10, T0, fix<2>); \ |
| traits.updateRhs(blB + (3 + 4 * K) * Traits::RhsProgress, rhs_panel); \ |
| traits.madd(A0, rhs_panel, C3, T0, fix<3>); \ |
| traits.madd(A1, rhs_panel, C7, T0, fix<3>); \ |
| traits.madd(A2, rhs_panel, C11, T0, fix<3>); \ |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 3pX4"); \ |
| } while (false) |
| |
| internal::prefetch(blB); |
| EIGEN_GEBP_ONESTEP(0); |
| EIGEN_GEBP_ONESTEP(1); |
| EIGEN_GEBP_ONESTEP(2); |
| EIGEN_GEBP_ONESTEP(3); |
| EIGEN_GEBP_ONESTEP(4); |
| EIGEN_GEBP_ONESTEP(5); |
| EIGEN_GEBP_ONESTEP(6); |
| EIGEN_GEBP_ONESTEP(7); |
| |
| blB += pk * 4 * RhsProgress; |
| blA += pk * 3 * Traits::LhsProgress; |
| |
| EIGEN_ASM_COMMENT("end gebp micro kernel 3pX4"); |
| } |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPanel15 rhs_panel; |
| RhsPacket T0; |
| LhsPacket A2; |
| EIGEN_GEBP_ONESTEP(0); |
| blB += 4 * RhsProgress; |
| blA += 3 * Traits::LhsProgress; |
| } |
| |
| #undef EIGEN_GEBP_ONESTEP |
| |
| ResPacket R0, R1, R2; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| |
| R0 = r0.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r0.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r0.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C0, alphav, R0); |
| traits.acc(C4, alphav, R1); |
| traits.acc(C8, alphav, R2); |
| r0.storePacket(0 * Traits::ResPacketSize, R0); |
| r0.storePacket(1 * Traits::ResPacketSize, R1); |
| r0.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r1.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r1.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r1.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C1, alphav, R0); |
| traits.acc(C5, alphav, R1); |
| traits.acc(C9, alphav, R2); |
| r1.storePacket(0 * Traits::ResPacketSize, R0); |
| r1.storePacket(1 * Traits::ResPacketSize, R1); |
| r1.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r2.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r2.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r2.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C2, alphav, R0); |
| traits.acc(C6, alphav, R1); |
| traits.acc(C10, alphav, R2); |
| r2.storePacket(0 * Traits::ResPacketSize, R0); |
| r2.storePacket(1 * Traits::ResPacketSize, R1); |
| r2.storePacket(2 * Traits::ResPacketSize, R2); |
| |
| R0 = r3.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r3.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r3.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C3, alphav, R0); |
| traits.acc(C7, alphav, R1); |
| traits.acc(C11, alphav, R2); |
| r3.storePacket(0 * Traits::ResPacketSize, R0); |
| r3.storePacket(1 * Traits::ResPacketSize, R1); |
| r3.storePacket(2 * Traits::ResPacketSize, R2); |
| } |
| } |
| |
| // Deal with remaining columns of the rhs |
| for (Index j2 = packet_cols4; j2 < cols; j2++) { |
| for (Index i = i1; i < actual_panel_end; i += 3 * LhsProgress) { |
| // One column at a time |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (3 * Traits::LhsProgress)]; |
| prefetch(&blA[0]); |
| |
| // gets res block as register |
| AccPacket C0, C4, C8; |
| traits.initAcc(C0); |
| traits.initAcc(C4); |
| traits.initAcc(C8); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2); |
| r0.prefetch(0); |
| |
| // performs "inner" products |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB]; |
| LhsPacket A0, A1, A2; |
| |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 3pX1"); |
| RhsPacket B_0; |
| #define EIGEN_GEBGP_ONESTEP(K) \ |
| do { \ |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 3pX1"); \ |
| EIGEN_ASM_COMMENT("Note: these asm comments work around bug 935!"); \ |
| traits.loadLhs(&blA[(0 + 3 * K) * LhsProgress], A0); \ |
| traits.loadLhs(&blA[(1 + 3 * K) * LhsProgress], A1); \ |
| traits.loadLhs(&blA[(2 + 3 * K) * LhsProgress], A2); \ |
| traits.loadRhs(&blB[(0 + K) * RhsProgress], B_0); \ |
| traits.madd(A0, B_0, C0, B_0, fix<0>); \ |
| traits.madd(A1, B_0, C4, B_0, fix<0>); \ |
| traits.madd(A2, B_0, C8, B_0, fix<0>); \ |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 3pX1"); \ |
| } while (false) |
| |
| EIGEN_GEBGP_ONESTEP(0); |
| EIGEN_GEBGP_ONESTEP(1); |
| EIGEN_GEBGP_ONESTEP(2); |
| EIGEN_GEBGP_ONESTEP(3); |
| EIGEN_GEBGP_ONESTEP(4); |
| EIGEN_GEBGP_ONESTEP(5); |
| EIGEN_GEBGP_ONESTEP(6); |
| EIGEN_GEBGP_ONESTEP(7); |
| |
| blB += int(pk) * int(RhsProgress); |
| blA += int(pk) * 3 * int(Traits::LhsProgress); |
| |
| EIGEN_ASM_COMMENT("end gebp micro kernel 3pX1"); |
| } |
| |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPacket B_0; |
| EIGEN_GEBGP_ONESTEP(0); |
| blB += RhsProgress; |
| blA += 3 * Traits::LhsProgress; |
| } |
| #undef EIGEN_GEBGP_ONESTEP |
| ResPacket R0, R1, R2; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| |
| R0 = r0.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r0.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r0.template loadPacket<ResPacket>(2 * Traits::ResPacketSize); |
| traits.acc(C0, alphav, R0); |
| traits.acc(C4, alphav, R1); |
| traits.acc(C8, alphav, R2); |
| r0.storePacket(0 * Traits::ResPacketSize, R0); |
| r0.storePacket(1 * Traits::ResPacketSize, R1); |
| r0.storePacket(2 * Traits::ResPacketSize, R2); |
| } |
| } |
| } |
| } |
| |
| //---------- Process 2 * LhsProgress rows at once ---------- |
| if (mr >= 2 * Traits::LhsProgress) { |
| const Index l1 = defaultL1CacheSize; // in Bytes, TODO, l1 should be passed to this function. |
| // The max(1, ...) here is needed because we may be using blocking params larger than what our known l1 cache size |
| // suggests we should be using: either because our known l1 cache size is inaccurate (e.g. on Android, we can only |
| // guess), or because we are testing specific blocking sizes. |
| Index actual_panel_rows = |
| (2 * LhsProgress) * std::max<Index>(1, ((l1 - sizeof(ResScalar) * mr * nr - depth * nr * sizeof(RhsScalar)) / |
| (depth * sizeof(LhsScalar) * 2 * LhsProgress))); |
| |
| for (Index i1 = peeled_mc3; i1 < peeled_mc2; i1 += actual_panel_rows) { |
| Index actual_panel_end = (std::min)(i1 + actual_panel_rows, peeled_mc2); |
| #if EIGEN_ARCH_ARM64 |
| EIGEN_IF_CONSTEXPR(nr >= 8) { |
| for (Index j2 = 0; j2 < packet_cols8; j2 += 8) { |
| for (Index i = i1; i < actual_panel_end; i += 2 * LhsProgress) { |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (2 * Traits::LhsProgress)]; |
| prefetch(&blA[0]); |
| |
| AccPacket C0, C1, C2, C3, C4, C5, C6, C7, C8, C9, C10, C11, C12, C13, C14, C15; |
| traits.initAcc(C0); |
| traits.initAcc(C1); |
| traits.initAcc(C2); |
| traits.initAcc(C3); |
| traits.initAcc(C4); |
| traits.initAcc(C5); |
| traits.initAcc(C6); |
| traits.initAcc(C7); |
| traits.initAcc(C8); |
| traits.initAcc(C9); |
| traits.initAcc(C10); |
| traits.initAcc(C11); |
| traits.initAcc(C12); |
| traits.initAcc(C13); |
| traits.initAcc(C14); |
| traits.initAcc(C15); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2 + 0); |
| LinearMapper r1 = res.getLinearMapper(i, j2 + 1); |
| LinearMapper r2 = res.getLinearMapper(i, j2 + 2); |
| LinearMapper r3 = res.getLinearMapper(i, j2 + 3); |
| LinearMapper r4 = res.getLinearMapper(i, j2 + 4); |
| LinearMapper r5 = res.getLinearMapper(i, j2 + 5); |
| LinearMapper r6 = res.getLinearMapper(i, j2 + 6); |
| LinearMapper r7 = res.getLinearMapper(i, j2 + 7); |
| r0.prefetch(prefetch_res_offset); |
| r1.prefetch(prefetch_res_offset); |
| r2.prefetch(prefetch_res_offset); |
| r3.prefetch(prefetch_res_offset); |
| r4.prefetch(prefetch_res_offset); |
| r5.prefetch(prefetch_res_offset); |
| r6.prefetch(prefetch_res_offset); |
| r7.prefetch(prefetch_res_offset); |
| |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB * 8]; |
| prefetch(&blB[0]); |
| LhsPacket A0, A1; |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| RhsPacketx4 rhs_panel; |
| RhsPacket T0; |
| // NOTE: the begin/end asm comments below work around bug 935! |
| // but they are not enough for gcc>=6 without FMA (bug 1637) |
| #if EIGEN_GNUC_STRICT_AT_LEAST(6, 0, 0) && defined(EIGEN_VECTORIZE_SSE) |
| #define EIGEN_GEBP_2Px8_SPILLING_WORKAROUND __asm__("" : [a0] "+x,m"(A0), [a1] "+x,m"(A1)); |
| #else |
| #define EIGEN_GEBP_2Px8_SPILLING_WORKAROUND |
| #endif |
| #define EIGEN_GEBGP_ONESTEP(K) \ |
| do { \ |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 2pX8"); \ |
| traits.loadLhs(&blA[(0 + 2 * K) * LhsProgress], A0); \ |
| traits.loadLhs(&blA[(1 + 2 * K) * LhsProgress], A1); \ |
| traits.loadRhs(&blB[(0 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C0, T0, fix<0>); \ |
| traits.madd(A1, rhs_panel, C8, T0, fix<0>); \ |
| traits.updateRhs(&blB[(1 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C1, T0, fix<1>); \ |
| traits.madd(A1, rhs_panel, C9, T0, fix<1>); \ |
| traits.updateRhs(&blB[(2 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C2, T0, fix<2>); \ |
| traits.madd(A1, rhs_panel, C10, T0, fix<2>); \ |
| traits.updateRhs(&blB[(3 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C3, T0, fix<3>); \ |
| traits.madd(A1, rhs_panel, C11, T0, fix<3>); \ |
| traits.loadRhs(&blB[(4 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C4, T0, fix<0>); \ |
| traits.madd(A1, rhs_panel, C12, T0, fix<0>); \ |
| traits.updateRhs(&blB[(5 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C5, T0, fix<1>); \ |
| traits.madd(A1, rhs_panel, C13, T0, fix<1>); \ |
| traits.updateRhs(&blB[(6 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C6, T0, fix<2>); \ |
| traits.madd(A1, rhs_panel, C14, T0, fix<2>); \ |
| traits.updateRhs(&blB[(7 + 8 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C7, T0, fix<3>); \ |
| traits.madd(A1, rhs_panel, C15, T0, fix<3>); \ |
| EIGEN_GEBP_2Px8_SPILLING_WORKAROUND EIGEN_ASM_COMMENT("end step of gebp micro kernel 2pX8"); \ |
| } while (false) |
| |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 2pX8"); |
| |
| EIGEN_GEBGP_ONESTEP(0); |
| EIGEN_GEBGP_ONESTEP(1); |
| EIGEN_GEBGP_ONESTEP(2); |
| EIGEN_GEBGP_ONESTEP(3); |
| EIGEN_GEBGP_ONESTEP(4); |
| EIGEN_GEBGP_ONESTEP(5); |
| EIGEN_GEBGP_ONESTEP(6); |
| EIGEN_GEBGP_ONESTEP(7); |
| |
| blB += pk * 8 * RhsProgress; |
| blA += pk * (2 * Traits::LhsProgress); |
| |
| EIGEN_ASM_COMMENT("end gebp micro kernel 2pX8"); |
| } |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPacketx4 rhs_panel; |
| RhsPacket T0; |
| EIGEN_GEBGP_ONESTEP(0); |
| blB += 8 * RhsProgress; |
| blA += 2 * Traits::LhsProgress; |
| } |
| |
| #undef EIGEN_GEBGP_ONESTEP |
| |
| ResPacket R0, R1, R2, R3; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| |
| R0 = r0.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r0.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r1.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R3 = r1.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| traits.acc(C0, alphav, R0); |
| traits.acc(C8, alphav, R1); |
| traits.acc(C1, alphav, R2); |
| traits.acc(C9, alphav, R3); |
| r0.storePacket(0 * Traits::ResPacketSize, R0); |
| r0.storePacket(1 * Traits::ResPacketSize, R1); |
| r1.storePacket(0 * Traits::ResPacketSize, R2); |
| r1.storePacket(1 * Traits::ResPacketSize, R3); |
| |
| R0 = r2.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r2.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r3.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R3 = r3.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| traits.acc(C2, alphav, R0); |
| traits.acc(C10, alphav, R1); |
| traits.acc(C3, alphav, R2); |
| traits.acc(C11, alphav, R3); |
| r2.storePacket(0 * Traits::ResPacketSize, R0); |
| r2.storePacket(1 * Traits::ResPacketSize, R1); |
| r3.storePacket(0 * Traits::ResPacketSize, R2); |
| r3.storePacket(1 * Traits::ResPacketSize, R3); |
| |
| R0 = r4.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r4.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r5.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R3 = r5.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| traits.acc(C4, alphav, R0); |
| traits.acc(C12, alphav, R1); |
| traits.acc(C5, alphav, R2); |
| traits.acc(C13, alphav, R3); |
| r4.storePacket(0 * Traits::ResPacketSize, R0); |
| r4.storePacket(1 * Traits::ResPacketSize, R1); |
| r5.storePacket(0 * Traits::ResPacketSize, R2); |
| r5.storePacket(1 * Traits::ResPacketSize, R3); |
| |
| R0 = r6.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r6.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r7.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R3 = r7.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| traits.acc(C6, alphav, R0); |
| traits.acc(C14, alphav, R1); |
| traits.acc(C7, alphav, R2); |
| traits.acc(C15, alphav, R3); |
| r6.storePacket(0 * Traits::ResPacketSize, R0); |
| r6.storePacket(1 * Traits::ResPacketSize, R1); |
| r7.storePacket(0 * Traits::ResPacketSize, R2); |
| r7.storePacket(1 * Traits::ResPacketSize, R3); |
| } |
| } |
| } |
| #endif |
| for (Index j2 = packet_cols8; j2 < packet_cols4; j2 += 4) { |
| for (Index i = i1; i < actual_panel_end; i += 2 * LhsProgress) { |
| // We selected a 2*Traits::LhsProgress x nr micro block of res which is entirely |
| // stored into 2 x nr registers. |
| |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (2 * Traits::LhsProgress)]; |
| prefetch(&blA[0]); |
| |
| // gets res block as register |
| AccPacket C0, C1, C2, C3, C4, C5, C6, C7; |
| traits.initAcc(C0); |
| traits.initAcc(C1); |
| traits.initAcc(C2); |
| traits.initAcc(C3); |
| traits.initAcc(C4); |
| traits.initAcc(C5); |
| traits.initAcc(C6); |
| traits.initAcc(C7); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2 + 0); |
| LinearMapper r1 = res.getLinearMapper(i, j2 + 1); |
| LinearMapper r2 = res.getLinearMapper(i, j2 + 2); |
| LinearMapper r3 = res.getLinearMapper(i, j2 + 3); |
| |
| r0.prefetch(prefetch_res_offset); |
| r1.prefetch(prefetch_res_offset); |
| r2.prefetch(prefetch_res_offset); |
| r3.prefetch(prefetch_res_offset); |
| |
| // performs "inner" products |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB * 4]; |
| prefetch(&blB[0]); |
| LhsPacket A0, A1; |
| |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 2pX4"); |
| RhsPacketx4 rhs_panel; |
| RhsPacket T0; |
| |
| // NOTE: the begin/end asm comments below work around bug 935! |
| // but they are not enough for gcc>=6 without FMA (bug 1637) |
| #if EIGEN_GNUC_STRICT_AT_LEAST(6, 0, 0) && defined(EIGEN_VECTORIZE_SSE) && !(EIGEN_COMP_LCC) |
| #define EIGEN_GEBP_2PX4_SPILLING_WORKAROUND __asm__("" : [a0] "+x,m"(A0), [a1] "+x,m"(A1)); |
| #else |
| #define EIGEN_GEBP_2PX4_SPILLING_WORKAROUND |
| #endif |
| #define EIGEN_GEBGP_ONESTEP(K) \ |
| do { \ |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 2pX4"); \ |
| traits.loadLhs(&blA[(0 + 2 * K) * LhsProgress], A0); \ |
| traits.loadLhs(&blA[(1 + 2 * K) * LhsProgress], A1); \ |
| traits.loadRhs(&blB[(0 + 4 * K) * RhsProgress], rhs_panel); \ |
| traits.madd(A0, rhs_panel, C0, T0, fix<0>); \ |
| traits.madd(A1, rhs_panel, C4, T0, fix<0>); \ |
| traits.madd(A0, rhs_panel, C1, T0, fix<1>); \ |
| traits.madd(A1, rhs_panel, C5, T0, fix<1>); \ |
| traits.madd(A0, rhs_panel, C2, T0, fix<2>); \ |
| traits.madd(A1, rhs_panel, C6, T0, fix<2>); \ |
| traits.madd(A0, rhs_panel, C3, T0, fix<3>); \ |
| traits.madd(A1, rhs_panel, C7, T0, fix<3>); \ |
| EIGEN_GEBP_2PX4_SPILLING_WORKAROUND \ |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 2pX4"); \ |
| } while (false) |
| |
| internal::prefetch(blB + (48 + 0)); |
| EIGEN_GEBGP_ONESTEP(0); |
| EIGEN_GEBGP_ONESTEP(1); |
| EIGEN_GEBGP_ONESTEP(2); |
| EIGEN_GEBGP_ONESTEP(3); |
| internal::prefetch(blB + (48 + 16)); |
| EIGEN_GEBGP_ONESTEP(4); |
| EIGEN_GEBGP_ONESTEP(5); |
| EIGEN_GEBGP_ONESTEP(6); |
| EIGEN_GEBGP_ONESTEP(7); |
| |
| blB += pk * 4 * RhsProgress; |
| blA += pk * (2 * Traits::LhsProgress); |
| |
| EIGEN_ASM_COMMENT("end gebp micro kernel 2pX4"); |
| } |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPacketx4 rhs_panel; |
| RhsPacket T0; |
| EIGEN_GEBGP_ONESTEP(0); |
| blB += 4 * RhsProgress; |
| blA += 2 * Traits::LhsProgress; |
| } |
| #undef EIGEN_GEBGP_ONESTEP |
| |
| ResPacket R0, R1, R2, R3; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| |
| R0 = r0.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r0.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r1.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R3 = r1.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| traits.acc(C0, alphav, R0); |
| traits.acc(C4, alphav, R1); |
| traits.acc(C1, alphav, R2); |
| traits.acc(C5, alphav, R3); |
| r0.storePacket(0 * Traits::ResPacketSize, R0); |
| r0.storePacket(1 * Traits::ResPacketSize, R1); |
| r1.storePacket(0 * Traits::ResPacketSize, R2); |
| r1.storePacket(1 * Traits::ResPacketSize, R3); |
| |
| R0 = r2.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r2.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| R2 = r3.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R3 = r3.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| traits.acc(C2, alphav, R0); |
| traits.acc(C6, alphav, R1); |
| traits.acc(C3, alphav, R2); |
| traits.acc(C7, alphav, R3); |
| r2.storePacket(0 * Traits::ResPacketSize, R0); |
| r2.storePacket(1 * Traits::ResPacketSize, R1); |
| r3.storePacket(0 * Traits::ResPacketSize, R2); |
| r3.storePacket(1 * Traits::ResPacketSize, R3); |
| } |
| } |
| |
| // Deal with remaining columns of the rhs |
| for (Index j2 = packet_cols4; j2 < cols; j2++) { |
| for (Index i = i1; i < actual_panel_end; i += 2 * LhsProgress) { |
| // One column at a time |
| const LhsScalar* blA = &blockA[i * strideA + offsetA * (2 * Traits::LhsProgress)]; |
| prefetch(&blA[0]); |
| |
| // gets res block as register |
| AccPacket C0, C4; |
| traits.initAcc(C0); |
| traits.initAcc(C4); |
| |
| LinearMapper r0 = res.getLinearMapper(i, j2); |
| r0.prefetch(prefetch_res_offset); |
| |
| // performs "inner" products |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB]; |
| LhsPacket A0, A1; |
| |
| for (Index k = 0; k < peeled_kc; k += pk) { |
| EIGEN_ASM_COMMENT("begin gebp micro kernel 2pX1"); |
| RhsPacket B_0, B1; |
| |
| #define EIGEN_GEBGP_ONESTEP(K) \ |
| do { \ |
| EIGEN_ASM_COMMENT("begin step of gebp micro kernel 2pX1"); \ |
| EIGEN_ASM_COMMENT("Note: these asm comments work around bug 935!"); \ |
| traits.loadLhs(&blA[(0 + 2 * K) * LhsProgress], A0); \ |
| traits.loadLhs(&blA[(1 + 2 * K) * LhsProgress], A1); \ |
| traits.loadRhs(&blB[(0 + K) * RhsProgress], B_0); \ |
| traits.madd(A0, B_0, C0, B1, fix<0>); \ |
| traits.madd(A1, B_0, C4, B_0, fix<0>); \ |
| EIGEN_ASM_COMMENT("end step of gebp micro kernel 2pX1"); \ |
| } while (false) |
| |
| EIGEN_GEBGP_ONESTEP(0); |
| EIGEN_GEBGP_ONESTEP(1); |
| EIGEN_GEBGP_ONESTEP(2); |
| EIGEN_GEBGP_ONESTEP(3); |
| EIGEN_GEBGP_ONESTEP(4); |
| EIGEN_GEBGP_ONESTEP(5); |
| EIGEN_GEBGP_ONESTEP(6); |
| EIGEN_GEBGP_ONESTEP(7); |
| |
| blB += int(pk) * int(RhsProgress); |
| blA += int(pk) * 2 * int(Traits::LhsProgress); |
| |
| EIGEN_ASM_COMMENT("end gebp micro kernel 2pX1"); |
| } |
| |
| // process remaining peeled loop |
| for (Index k = peeled_kc; k < depth; k++) { |
| RhsPacket B_0, B1; |
| EIGEN_GEBGP_ONESTEP(0); |
| blB += RhsProgress; |
| blA += 2 * Traits::LhsProgress; |
| } |
| #undef EIGEN_GEBGP_ONESTEP |
| ResPacket R0, R1; |
| ResPacket alphav = pset1<ResPacket>(alpha); |
| |
| R0 = r0.template loadPacket<ResPacket>(0 * Traits::ResPacketSize); |
| R1 = r0.template loadPacket<ResPacket>(1 * Traits::ResPacketSize); |
| traits.acc(C0, alphav, R0); |
| traits.acc(C4, alphav, R1); |
| r0.storePacket(0 * Traits::ResPacketSize, R0); |
| r0.storePacket(1 * Traits::ResPacketSize, R1); |
| } |
| } |
| } |
| } |
| //---------- Process 1 * LhsProgress rows at once ---------- |
| if (mr >= 1 * Traits::LhsProgress) { |
| lhs_process_one_packet<nr, LhsProgress, RhsProgress, LhsScalar, RhsScalar, ResScalar, AccPacket, LhsPacket, |
| RhsPacket, ResPacket, Traits, LinearMapper, DataMapper> |
| p; |
| p(res, blockA, blockB, alpha, peeled_mc2, peeled_mc1, strideA, strideB, offsetA, offsetB, prefetch_res_offset, |
| peeled_kc, pk, cols, depth, packet_cols4); |
| } |
| //---------- Process LhsProgressHalf rows at once ---------- |
| if ((LhsProgressHalf < LhsProgress) && mr >= LhsProgressHalf) { |
| lhs_process_fraction_of_packet<nr, LhsProgressHalf, RhsProgressHalf, LhsScalar, RhsScalar, ResScalar, AccPacketHalf, |
| LhsPacketHalf, RhsPacketHalf, ResPacketHalf, HalfTraits, LinearMapper, DataMapper> |
| p; |
| p(res, blockA, blockB, alpha, peeled_mc1, peeled_mc_half, strideA, strideB, offsetA, offsetB, prefetch_res_offset, |
| peeled_kc, pk, cols, depth, packet_cols4); |
| } |
| //---------- Process LhsProgressQuarter rows at once ---------- |
| if ((LhsProgressQuarter < LhsProgressHalf) && mr >= LhsProgressQuarter) { |
| lhs_process_fraction_of_packet<nr, LhsProgressQuarter, RhsProgressQuarter, LhsScalar, RhsScalar, ResScalar, |
| AccPacketQuarter, LhsPacketQuarter, RhsPacketQuarter, ResPacketQuarter, |
| QuarterTraits, LinearMapper, DataMapper> |
| p; |
| p(res, blockA, blockB, alpha, peeled_mc_half, peeled_mc_quarter, strideA, strideB, offsetA, offsetB, |
| prefetch_res_offset, peeled_kc, pk, cols, depth, packet_cols4); |
| } |
| //---------- Process remaining rows, 1 at once ---------- |
| if (peeled_mc_quarter < rows) { |
| #if EIGEN_ARCH_ARM64 |
| EIGEN_IF_CONSTEXPR(nr >= 8) { |
| // loop on each panel of the rhs |
| for (Index j2 = 0; j2 < packet_cols8; j2 += 8) { |
| // loop on each row of the lhs (1*LhsProgress x depth) |
| for (Index i = peeled_mc_quarter; i < rows; i += 1) { |
| const LhsScalar* blA = &blockA[i * strideA + offsetA]; |
| prefetch(&blA[0]); |
| // gets a 1 x 1 res block as registers |
| ResScalar C0(0), C1(0), C2(0), C3(0), C4(0), C5(0), C6(0), C7(0); |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB * 8]; |
| for (Index k = 0; k < depth; k++) { |
| LhsScalar A0 = blA[k]; |
| RhsScalar B_0; |
| |
| B_0 = blB[0]; |
| C0 = cj.pmadd(A0, B_0, C0); |
| |
| B_0 = blB[1]; |
| C1 = cj.pmadd(A0, B_0, C1); |
| |
| B_0 = blB[2]; |
| C2 = cj.pmadd(A0, B_0, C2); |
| |
| B_0 = blB[3]; |
| C3 = cj.pmadd(A0, B_0, C3); |
| |
| B_0 = blB[4]; |
| C4 = cj.pmadd(A0, B_0, C4); |
| |
| B_0 = blB[5]; |
| C5 = cj.pmadd(A0, B_0, C5); |
| |
| B_0 = blB[6]; |
| C6 = cj.pmadd(A0, B_0, C6); |
| |
| B_0 = blB[7]; |
| C7 = cj.pmadd(A0, B_0, C7); |
| |
| blB += 8; |
| } |
| res(i, j2 + 0) += alpha * C0; |
| res(i, j2 + 1) += alpha * C1; |
| res(i, j2 + 2) += alpha * C2; |
| res(i, j2 + 3) += alpha * C3; |
| res(i, j2 + 4) += alpha * C4; |
| res(i, j2 + 5) += alpha * C5; |
| res(i, j2 + 6) += alpha * C6; |
| res(i, j2 + 7) += alpha * C7; |
| } |
| } |
| } |
| #endif |
| |
| for (Index j2 = packet_cols8; j2 < packet_cols4; j2 += 4) { |
| // loop on each row of the lhs (1*LhsProgress x depth) |
| for (Index i = peeled_mc_quarter; i < rows; i += 1) { |
| const LhsScalar* blA = &blockA[i * strideA + offsetA]; |
| prefetch(&blA[0]); |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB * 4]; |
| |
| // If LhsProgress is 8 or 16, it assumes that there is a |
| // half or quarter packet, respectively, of the same size as |
| // nr (which is currently 4) for the return type. |
| const int SResPacketHalfSize = unpacket_traits<typename unpacket_traits<SResPacket>::half>::size; |
| const int SResPacketQuarterSize = |
| unpacket_traits<typename unpacket_traits<typename unpacket_traits<SResPacket>::half>::half>::size; |
| // The following code assumes we can load SRhsPacket in such a way that |
| // it multiplies blocks of 4 elements in SLhsPacket. This is not the |
| // case for some customized kernels (i.e. NEON fp16). If the assumption |
| // fails, drop down to the scalar path. |
| constexpr bool kCanLoadSRhsQuad = |
| (unpacket_traits<SLhsPacket>::size < 4) || |
| (unpacket_traits<SRhsPacket>::size % ((std::max<int>)(unpacket_traits<SLhsPacket>::size, 4) / 4)) == 0; |
| if (kCanLoadSRhsQuad && (SwappedTraits::LhsProgress % 4) == 0 && (SwappedTraits::LhsProgress <= 16) && |
| (SwappedTraits::LhsProgress != 8 || SResPacketHalfSize == nr) && |
| (SwappedTraits::LhsProgress != 16 || SResPacketQuarterSize == nr)) { |
| SAccPacket C0, C1, C2, C3; |
| straits.initAcc(C0); |
| straits.initAcc(C1); |
| straits.initAcc(C2); |
| straits.initAcc(C3); |
| |
| const Index spk = (std::max)(1, SwappedTraits::LhsProgress / 4); |
| const Index endk = (depth / spk) * spk; |
| const Index endk4 = (depth / (spk * 4)) * (spk * 4); |
| |
| Index k = 0; |
| for (; k < endk4; k += 4 * spk) { |
| SLhsPacket A0, A1; |
| SRhsPacket B_0, B_1; |
| |
| straits.loadLhsUnaligned(blB + 0 * SwappedTraits::LhsProgress, A0); |
| straits.loadLhsUnaligned(blB + 1 * SwappedTraits::LhsProgress, A1); |
| |
| straits.loadRhsQuad(blA + 0 * spk, B_0); |
| straits.loadRhsQuad(blA + 1 * spk, B_1); |
| straits.madd(A0, B_0, C0, B_0, fix<0>); |
| straits.madd(A1, B_1, C1, B_1, fix<0>); |
| |
| straits.loadLhsUnaligned(blB + 2 * SwappedTraits::LhsProgress, A0); |
| straits.loadLhsUnaligned(blB + 3 * SwappedTraits::LhsProgress, A1); |
| straits.loadRhsQuad(blA + 2 * spk, B_0); |
| straits.loadRhsQuad(blA + 3 * spk, B_1); |
| straits.madd(A0, B_0, C2, B_0, fix<0>); |
| straits.madd(A1, B_1, C3, B_1, fix<0>); |
| |
| blB += 4 * SwappedTraits::LhsProgress; |
| blA += 4 * spk; |
| } |
| C0 = padd(padd(C0, C1), padd(C2, C3)); |
| for (; k < endk; k += spk) { |
| SLhsPacket A0; |
| SRhsPacket B_0; |
| |
| straits.loadLhsUnaligned(blB, A0); |
| straits.loadRhsQuad(blA, B_0); |
| straits.madd(A0, B_0, C0, B_0, fix<0>); |
| |
| blB += SwappedTraits::LhsProgress; |
| blA += spk; |
| } |
| if (SwappedTraits::LhsProgress == 8) { |
| // Special case where we have to first reduce the accumulation register C0 |
| typedef std::conditional_t<SwappedTraits::LhsProgress >= 8, typename unpacket_traits<SResPacket>::half, |
| SResPacket> |
| SResPacketHalf; |
| typedef std::conditional_t<SwappedTraits::LhsProgress >= 8, typename unpacket_traits<SLhsPacket>::half, |
| SLhsPacket> |
| SLhsPacketHalf; |
| typedef std::conditional_t<SwappedTraits::LhsProgress >= 8, typename unpacket_traits<SRhsPacket>::half, |
| SRhsPacket> |
| SRhsPacketHalf; |
| typedef std::conditional_t<SwappedTraits::LhsProgress >= 8, typename unpacket_traits<SAccPacket>::half, |
| SAccPacket> |
| SAccPacketHalf; |
| |
| SResPacketHalf R = res.template gatherPacket<SResPacketHalf>(i, j2); |
| SResPacketHalf alphav = pset1<SResPacketHalf>(alpha); |
| |
| if (depth - endk > 0) { |
| // We have to handle the last row of the rhs which corresponds to a half-packet |
| SLhsPacketHalf a0; |
| SRhsPacketHalf b0; |
| straits.loadLhsUnaligned(blB, a0); |
| straits.loadRhs(blA, b0); |
| SAccPacketHalf c0 = predux_half_dowto4(C0); |
| straits.madd(a0, b0, c0, b0, fix<0>); |
| straits.acc(c0, alphav, R); |
| } else { |
| straits.acc(predux_half_dowto4(C0), alphav, R); |
| } |
| res.scatterPacket(i, j2, R); |
| } else if (SwappedTraits::LhsProgress == 16) { |
| // Special case where we have to first reduce the |
| // accumulation register C0. We specialize the block in |
| // template form, so that LhsProgress < 16 paths don't |
| // fail to compile |
| last_row_process_16_packets<LhsScalar, RhsScalar, Index, DataMapper, mr, nr, ConjugateLhs, ConjugateRhs> p; |
| p(res, straits, blA, blB, depth, endk, i, j2, alpha, C0); |
| } else { |
| SResPacket R = res.template gatherPacket<SResPacket>(i, j2); |
| SResPacket alphav = pset1<SResPacket>(alpha); |
| straits.acc(C0, alphav, R); |
| res.scatterPacket(i, j2, R); |
| } |
| } else // scalar path |
| { |
| // get a 1 x 4 res block as registers |
| ResScalar C0(0), C1(0), C2(0), C3(0); |
| |
| for (Index k = 0; k < depth; k++) { |
| LhsScalar A0; |
| RhsScalar B_0, B_1; |
| |
| A0 = blA[k]; |
| |
| B_0 = blB[0]; |
| B_1 = blB[1]; |
| C0 = cj.pmadd(A0, B_0, C0); |
| C1 = cj.pmadd(A0, B_1, C1); |
| |
| B_0 = blB[2]; |
| B_1 = blB[3]; |
| C2 = cj.pmadd(A0, B_0, C2); |
| C3 = cj.pmadd(A0, B_1, C3); |
| |
| blB += 4; |
| } |
| res(i, j2 + 0) += alpha * C0; |
| res(i, j2 + 1) += alpha * C1; |
| res(i, j2 + 2) += alpha * C2; |
| res(i, j2 + 3) += alpha * C3; |
| } |
| } |
| } |
| // remaining columns |
| for (Index j2 = packet_cols4; j2 < cols; j2++) { |
| // loop on each row of the lhs (1*LhsProgress x depth) |
| for (Index i = peeled_mc_quarter; i < rows; i += 1) { |
| const LhsScalar* blA = &blockA[i * strideA + offsetA]; |
| prefetch(&blA[0]); |
| // gets a 1 x 1 res block as registers |
| ResScalar C0(0); |
| const RhsScalar* blB = &blockB[j2 * strideB + offsetB]; |
| for (Index k = 0; k < depth; k++) { |
| LhsScalar A0 = blA[k]; |
| RhsScalar B_0 = blB[k]; |
| C0 = cj.pmadd(A0, B_0, C0); |
| } |
| res(i, j2) += alpha * C0; |
| } |
| } |
| } |
| } |
| |
| // pack a block of the lhs |
| // The traversal is as follow (mr==4): |
| // 0 4 8 12 ... |
| // 1 5 9 13 ... |
| // 2 6 10 14 ... |
| // 3 7 11 15 ... |
| // |
| // 16 20 24 28 ... |
| // 17 21 25 29 ... |
| // 18 22 26 30 ... |
| // 19 23 27 31 ... |
| // |
| // 32 33 34 35 ... |
| // 36 36 38 39 ... |
| template <typename Scalar, typename Index, typename DataMapper, int Pack1, int Pack2, typename Packet, bool Conjugate, |
| bool PanelMode> |
| struct gemm_pack_lhs<Scalar, Index, DataMapper, Pack1, Pack2, Packet, ColMajor, Conjugate, PanelMode> { |
| typedef typename DataMapper::LinearMapper LinearMapper; |
| EIGEN_DONT_INLINE void operator()(Scalar* blockA, const DataMapper& lhs, Index depth, Index rows, Index stride = 0, |
| Index offset = 0); |
| }; |
| |
| template <typename Scalar, typename Index, typename DataMapper, int Pack1, int Pack2, typename Packet, bool Conjugate, |
| bool PanelMode> |
| EIGEN_DONT_INLINE void gemm_pack_lhs<Scalar, Index, DataMapper, Pack1, Pack2, Packet, ColMajor, Conjugate, |
| PanelMode>::operator()(Scalar* blockA, const DataMapper& lhs, Index depth, |
| Index rows, Index stride, Index offset) { |
| typedef typename unpacket_traits<Packet>::half HalfPacket; |
| typedef typename unpacket_traits<typename unpacket_traits<Packet>::half>::half QuarterPacket; |
| enum { |
| PacketSize = unpacket_traits<Packet>::size, |
| HalfPacketSize = unpacket_traits<HalfPacket>::size, |
| QuarterPacketSize = unpacket_traits<QuarterPacket>::size, |
| HasHalf = (int)HalfPacketSize < (int)PacketSize, |
| HasQuarter = (int)QuarterPacketSize < (int)HalfPacketSize |
| }; |
| |
| EIGEN_ASM_COMMENT("EIGEN PRODUCT PACK LHS"); |
| EIGEN_UNUSED_VARIABLE(stride); |
| EIGEN_UNUSED_VARIABLE(offset); |
| eigen_assert(((!PanelMode) && stride == 0 && offset == 0) || (PanelMode && stride >= depth && offset <= stride)); |
| eigen_assert(((Pack1 % PacketSize) == 0 && Pack1 <= 4 * PacketSize) || (Pack1 <= 4)); |
| conj_if<NumTraits<Scalar>::IsComplex && Conjugate> cj; |
| Index count = 0; |
| |
| const Index peeled_mc3 = Pack1 >= 3 * PacketSize ? (rows / (3 * PacketSize)) * (3 * PacketSize) : 0; |
| const Index peeled_mc2 = |
| Pack1 >= 2 * PacketSize ? peeled_mc3 + ((rows - peeled_mc3) / (2 * PacketSize)) * (2 * PacketSize) : 0; |
| const Index peeled_mc1 = |
| Pack1 >= 1 * PacketSize ? peeled_mc2 + ((rows - peeled_mc2) / (1 * PacketSize)) * (1 * PacketSize) : 0; |
| const Index peeled_mc_half = |
| Pack1 >= HalfPacketSize ? peeled_mc1 + ((rows - peeled_mc1) / (HalfPacketSize)) * (HalfPacketSize) : 0; |
| const Index peeled_mc_quarter = Pack1 >= QuarterPacketSize ? (rows / (QuarterPacketSize)) * (QuarterPacketSize) : 0; |
| const Index last_lhs_progress = rows > peeled_mc_quarter ? (rows - peeled_mc_quarter) & ~1 : 0; |
| const Index peeled_mc0 = Pack2 >= PacketSize ? peeled_mc_quarter |
| : Pack2 > 1 && last_lhs_progress ? (rows / last_lhs_progress) * last_lhs_progress |
| : 0; |
| |
| Index i = 0; |
| |
| // Pack 3 packets |
| if (Pack1 >= 3 * PacketSize) { |
| for (; i < peeled_mc3; i += 3 * PacketSize) { |
| if (PanelMode) count += (3 * PacketSize) * offset; |
| |
| for (Index k = 0; k < depth; k++) { |
| Packet A, B, C; |
| A = lhs.template loadPacket<Packet>(i + 0 * PacketSize, k); |
| B = lhs.template loadPacket<Packet>(i + 1 * PacketSize, k); |
| C = lhs.template loadPacket<Packet>(i + 2 * PacketSize, k); |
| pstore(blockA + count, cj.pconj(A)); |
| count += PacketSize; |
| pstore(blockA + count, cj.pconj(B)); |
| count += PacketSize; |
| pstore(blockA + count, cj.pconj(C)); |
| count += PacketSize; |
| } |
| if (PanelMode) count += (3 * PacketSize) * (stride - offset - depth); |
| } |
| } |
| // Pack 2 packets |
| if (Pack1 >= 2 * PacketSize) { |
| for (; i < peeled_mc2; i += 2 * PacketSize) { |
| if (PanelMode) count += (2 * PacketSize) * offset; |
| |
| for (Index k = 0; k < depth; k++) { |
| Packet A, B; |
| A = lhs.template loadPacket<Packet>(i + 0 * PacketSize, k); |
| B = lhs.template loadPacket<Packet>(i + 1 * PacketSize, k); |
| pstore(blockA + count, cj.pconj(A)); |
| count += PacketSize; |
| pstore(blockA + count, cj.pconj(B)); |
| count += PacketSize; |
| } |
| if (PanelMode) count += (2 * PacketSize) * (stride - offset - depth); |
| } |
| } |
| // Pack 1 packets |
| if (Pack1 >= 1 * PacketSize) { |
| for (; i < peeled_mc1; i += 1 * PacketSize) { |
| if (PanelMode) count += (1 * PacketSize) * offset; |
| |
| for (Index k = 0; k < depth; k++) { |
| Packet A; |
| A = lhs.template loadPacket<Packet>(i + 0 * PacketSize, k); |
| pstore(blockA + count, cj.pconj(A)); |
| count += PacketSize; |
| } |
| if (PanelMode) count += (1 * PacketSize) * (stride - offset - depth); |
| } |
| } |
| // Pack half packets |
| if (HasHalf && Pack1 >= HalfPacketSize) { |
| for (; i < peeled_mc_half; i += HalfPacketSize) { |
| if (PanelMode) count += (HalfPacketSize)*offset; |
| |
| for (Index k = 0; k < depth; k++) { |
| HalfPacket A; |
| A = lhs.template loadPacket<HalfPacket>(i + 0 * (HalfPacketSize), k); |
| pstoreu(blockA + count, cj.pconj(A)); |
| count += HalfPacketSize; |
| } |
| if (PanelMode) count += (HalfPacketSize) * (stride - offset - depth); |
| } |
| } |
| // Pack quarter packets |
| if (HasQuarter && Pack1 >= QuarterPacketSize) { |
| for (; i < peeled_mc_quarter; i += QuarterPacketSize) { |
| if (PanelMode) count += (QuarterPacketSize)*offset; |
| |
| for (Index k = 0; k < depth; k++) { |
| QuarterPacket A; |
| A = lhs.template loadPacket<QuarterPacket>(i + 0 * (QuarterPacketSize), k); |
| pstoreu(blockA + count, cj.pconj(A)); |
| count += QuarterPacketSize; |
| } |
| if (PanelMode) count += (QuarterPacketSize) * (stride - offset - depth); |
| } |
| } |
| // Pack2 may be *smaller* than PacketSize—that happens for |
| // products like real * complex, where we have to go half the |
| // progress on the lhs in order to duplicate those operands to |
| // address both real & imaginary parts on the rhs. This portion will |
| // pack those half ones until they match the number expected on the |
| // last peeling loop at this point (for the rhs). |
| if (Pack2 < PacketSize && Pack2 > 1) { |
| for (; i < peeled_mc0; i += last_lhs_progress) { |
| if (PanelMode) count += last_lhs_progress * offset; |
| |
| for (Index k = 0; k < depth; k++) |
| for (Index w = 0; w < last_lhs_progress; w++) blockA[count++] = cj(lhs(i + w, k)); |
| |
| if (PanelMode) count += last_lhs_progress * (stride - offset - depth); |
| } |
| } |
| // Pack scalars |
| for (; i < rows; i++) { |
| if (PanelMode) count += offset; |
| for (Index k = 0; k < depth; k++) blockA[count++] = cj(lhs(i, k)); |
| if (PanelMode) count += (stride - offset - depth); |
| } |
| } |
| |
| template <typename Scalar, typename Index, typename DataMapper, int Pack1, int Pack2, typename Packet, bool Conjugate, |
| bool PanelMode> |
| struct gemm_pack_lhs<Scalar, Index, DataMapper, Pack1, Pack2, Packet, RowMajor, Conjugate, PanelMode> { |
| typedef typename DataMapper::LinearMapper LinearMapper; |
| EIGEN_DONT_INLINE void operator()(Scalar* blockA, const DataMapper& lhs, Index depth, Index rows, Index stride = 0, |
| Index offset = 0); |
| }; |
| |
| template <typename Scalar, typename Index, typename DataMapper, int Pack1, int Pack2, typename Packet, bool Conjugate, |
| bool PanelMode> |
| EIGEN_DONT_INLINE void gemm_pack_lhs<Scalar, Index, DataMapper, Pack1, Pack2, Packet, RowMajor, Conjugate, |
| PanelMode>::operator()(Scalar* blockA, const DataMapper& lhs, Index depth, |
| Index rows, Index stride, Index offset) { |
| typedef typename unpacket_traits<Packet>::half HalfPacket; |
| typedef typename unpacket_traits<typename unpacket_traits<Packet>::half>::half QuarterPacket; |
| enum { |
| PacketSize = unpacket_traits<Packet>::size, |
| HalfPacketSize = unpacket_traits<HalfPacket>::size, |
| QuarterPacketSize = unpacket_traits<QuarterPacket>::size, |
| HasHalf = (int)HalfPacketSize < (int)PacketSize, |
| HasQuarter = (int)QuarterPacketSize < (int)HalfPacketSize |
| }; |
| |
| EIGEN_ASM_COMMENT("EIGEN PRODUCT PACK LHS"); |
| EIGEN_UNUSED_VARIABLE(stride); |
| EIGEN_UNUSED_VARIABLE(offset); |
| eigen_assert(((!PanelMode) && stride == 0 && offset == 0) || (PanelMode && stride >= depth && offset <= stride)); |
| conj_if<NumTraits<Scalar>::IsComplex && Conjugate> cj; |
| Index count = 0; |
| bool gone_half = false, gone_quarter = false, gone_last = false; |
| |
| Index i = 0; |
| Index pack = Pack1; |
| Index psize = PacketSize; |
| while (pack > 0) { |
| Index remaining_rows = rows - i; |
| Index peeled_mc = gone_last ? Pack2 > 1 ? (rows / pack) * pack : 0 : i + (remaining_rows / pack) * pack; |
| Index starting_pos = i; |
| for (; i < peeled_mc; i += pack) { |
| if (PanelMode) count += pack * offset; |
| |
| Index k = 0; |
| if (pack >= psize && psize >= QuarterPacketSize) { |
| const Index peeled_k = (depth / psize) * psize; |
| for (; k < peeled_k; k += psize) { |
| for (Index m = 0; m < pack; m += psize) { |
| if (psize == PacketSize) { |
| PacketBlock<Packet> kernel; |
| for (Index p = 0; p < psize; ++p) kernel.packet[p] = lhs.template loadPacket<Packet>(i + p + m, k); |
| ptranspose(kernel); |
| for (Index p = 0; p < psize; ++p) pstore(blockA + count + m + (pack)*p, cj.pconj(kernel.packet[p])); |
| } else if (HasHalf && psize == HalfPacketSize) { |
| gone_half = true; |
| PacketBlock<HalfPacket> kernel_half; |
| for (Index p = 0; p < psize; ++p) |
| kernel_half.packet[p] = lhs.template loadPacket<HalfPacket>(i + p + m, k); |
| ptranspose(kernel_half); |
| for (Index p = 0; p < psize; ++p) pstore(blockA + count + m + (pack)*p, cj.pconj(kernel_half.packet[p])); |
| } else if (HasQuarter && psize == QuarterPacketSize) { |
| gone_quarter = true; |
| PacketBlock<QuarterPacket> kernel_quarter; |
| for (Index p = 0; p < psize; ++p) |
| kernel_quarter.packet[p] = lhs.template loadPacket<QuarterPacket>(i + p + m, k); |
| ptranspose(kernel_quarter); |
| for (Index p = 0; p < psize; ++p) |
| pstore(blockA + count + m + (pack)*p, cj.pconj(kernel_quarter.packet[p])); |
| } |
| } |
| count += psize * pack; |
| } |
| } |
| |
| for (; k < depth; k++) { |
| Index w = 0; |
| for (; w < pack - 3; w += 4) { |
| Scalar a(cj(lhs(i + w + 0, k))), b(cj(lhs(i + w + 1, k))), c(cj(lhs(i + w + 2, k))), d(cj(lhs(i + w + 3, k))); |
| blockA[count++] = a; |
| blockA[count++] = b; |
| blockA[count++] = c; |
| blockA[count++] = d; |
| } |
| if (pack % 4) |
| for (; w < pack; ++w) blockA[count++] = cj(lhs(i + w, k)); |
| } |
| |
| if (PanelMode) count += pack * (stride - offset - depth); |
| } |
| |
| pack -= psize; |
| Index left = rows - i; |
| if (pack <= 0) { |
| if (!gone_last && (starting_pos == i || left >= psize / 2 || left >= psize / 4) && |
| ((psize / 2 == HalfPacketSize && HasHalf && !gone_half) || |
| (psize / 2 == QuarterPacketSize && HasQuarter && !gone_quarter))) { |
| psize /= 2; |
| pack = psize; |
| continue; |
| } |
| // Pack2 may be *smaller* than PacketSize—that happens for |
| // products like real * complex, where we have to go half the |
| // progress on the lhs in order to duplicate those operands to |
| // address both real & imaginary parts on the rhs. This portion will |
| // pack those half ones until they match the number expected on the |
| // last peeling loop at this point (for the rhs). |
| if (Pack2 < PacketSize && !gone_last) { |
| gone_last = true; |
| psize = pack = left & ~1; |
| } |
| } |
| } |
| |
| for (; i < rows; i++) { |
| if (PanelMode) count += offset; |
| for (Index k = 0; k < depth; k++) blockA[count++] = cj(lhs(i, k)); |
| if (PanelMode) count += (stride - offset - depth); |
| } |
| } |
| |
| // copy a complete panel of the rhs |
| // this version is optimized for column major matrices |
| // The traversal order is as follow: (nr==4): |
| // 0 1 2 3 12 13 14 15 24 27 |
| // 4 5 6 7 16 17 18 19 25 28 |
| // 8 9 10 11 20 21 22 23 26 29 |
| // . . . . . . . . . . |
| template <typename Scalar, typename Index, typename DataMapper, int nr, bool Conjugate, bool PanelMode> |
| struct gemm_pack_rhs<Scalar, Index, DataMapper, nr, ColMajor, Conjugate, PanelMode> { |
| typedef typename packet_traits<Scalar>::type Packet; |
| typedef typename DataMapper::LinearMapper LinearMapper; |
| enum { PacketSize = packet_traits<Scalar>::size }; |
| EIGEN_DONT_INLINE void operator()(Scalar* blockB, const DataMapper& rhs, Index depth, Index cols, Index stride = 0, |
| Index offset = 0); |
| }; |
| |
| template <typename Scalar, typename Index, typename DataMapper, int nr, bool Conjugate, bool PanelMode> |
| EIGEN_DONT_INLINE void gemm_pack_rhs<Scalar, Index, DataMapper, nr, ColMajor, Conjugate, PanelMode>::operator()( |
| Scalar* blockB, const DataMapper& rhs, Index depth, Index cols, Index stride, Index offset) { |
| EIGEN_ASM_COMMENT("EIGEN PRODUCT PACK RHS COLMAJOR"); |
| EIGEN_UNUSED_VARIABLE(stride); |
| EIGEN_UNUSED_VARIABLE(offset); |
| eigen_assert(((!PanelMode) && stride == 0 && offset == 0) || (PanelMode && stride >= depth && offset <= stride)); |
| conj_if<NumTraits<Scalar>::IsComplex && Conjugate> cj; |
| Index packet_cols8 = nr >= 8 ? (cols / 8) * 8 : 0; |
| Index packet_cols4 = nr >= 4 ? (cols / 4) * 4 : 0; |
| Index count = 0; |
| const Index peeled_k = (depth / PacketSize) * PacketSize; |
| |
| #if EIGEN_ARCH_ARM64 |
| EIGEN_IF_CONSTEXPR(nr >= 8) { |
| for (Index j2 = 0; j2 < packet_cols8; j2 += 8) { |
| // skip what we have before |
| if (PanelMode) count += 8 * offset; |
| const LinearMapper dm0 = rhs.getLinearMapper(0, j2 + 0); |
| const LinearMapper dm1 = rhs.getLinearMapper(0, j2 + 1); |
| const LinearMapper dm2 = rhs.getLinearMapper(0, j2 + 2); |
| const LinearMapper dm3 = rhs.getLinearMapper(0, j2 + 3); |
| const LinearMapper dm4 = rhs.getLinearMapper(0, j2 + 4); |
| const LinearMapper dm5 = rhs.getLinearMapper(0, j2 + 5); |
| const LinearMapper dm6 = rhs.getLinearMapper(0, j2 + 6); |
| const LinearMapper dm7 = rhs.getLinearMapper(0, j2 + 7); |
| Index k = 0; |
| if (PacketSize % 2 == 0 && PacketSize <= 8) // 2 4 8 |
| { |
| for (; k < peeled_k; k += PacketSize) { |
| if (PacketSize == 2) { |
| PacketBlock<Packet, PacketSize == 2 ? 2 : PacketSize> kernel0, kernel1, kernel2, kernel3; |
| kernel0.packet[0 % PacketSize] = dm0.template loadPacket<Packet>(k); |
| kernel0.packet[1 % PacketSize] = dm1.template loadPacket<Packet>(k); |
| kernel1.packet[0 % PacketSize] = dm2.template loadPacket<Packet>(k); |
| kernel1.packet[1 % PacketSize] = dm3.template loadPacket<Packet>(k); |
| kernel2.packet[0 % PacketSize] = dm4.template loadPacket<Packet>(k); |
| kernel2.packet[1 % PacketSize] = dm5.template loadPacket<Packet>(k); |
| kernel3.packet[0 % PacketSize] = dm6.template loadPacket<Packet>(k); |
| kernel3.packet[1 % PacketSize] = dm7.template loadPacket<Packet>(k); |
| ptranspose(kernel0); |
| ptranspose(kernel1); |
| ptranspose(kernel2); |
| ptranspose(kernel3); |
| |
| pstoreu(blockB + count + 0 * PacketSize, cj.pconj(kernel0.packet[0 % PacketSize])); |
| pstoreu(blockB + count + 1 * PacketSize, cj.pconj(kernel1.packet[0 % PacketSize])); |
| pstoreu(blockB + count + 2 * PacketSize, cj.pconj(kernel2.packet[0 % PacketSize])); |
| pstoreu(blockB + count + 3 * PacketSize, cj.pconj(kernel3.packet[0 % PacketSize])); |
| |
| pstoreu(blockB + count + 4 * PacketSize, cj.pconj(kernel0.packet[1 % PacketSize])); |
| pstoreu(blockB + count + 5 * PacketSize, cj.pconj(kernel1.packet[1 % PacketSize])); |
| pstoreu(blockB + count + 6 * PacketSize, cj.pconj(kernel2.packet[1 % PacketSize])); |
| pstoreu(blockB + count + 7 * PacketSize, cj.pconj(kernel3.packet[1 % PacketSize])); |
| count += 8 * PacketSize; |
| } else if (PacketSize == 4) { |
| PacketBlock<Packet, PacketSize == 4 ? 4 : PacketSize> kernel0, kernel1; |
| |
| kernel0.packet[0 % PacketSize] = dm0.template loadPacket<Packet>(k); |
| kernel0.packet[1 % PacketSize] = dm1.template loadPacket<Packet>(k); |
| kernel0.packet[2 % PacketSize] = dm2.template loadPacket<Packet>(k); |
| kernel0.packet[3 % PacketSize] = dm3.template loadPacket<Packet>(k); |
| kernel1.packet[0 % PacketSize] = dm4.template loadPacket<Packet>(k); |
| kernel1.packet[1 % PacketSize] = dm5.template loadPacket<Packet>(k); |
| kernel1.packet[2 % PacketSize] = dm6.template loadPacket<Packet>(k); |
| kernel1.packet[3 % PacketSize] = dm7.template loadPacket<Packet>(k); |
| ptranspose(kernel0); |
| ptranspose(kernel1); |
| |
| pstoreu(blockB + count + 0 * PacketSize, cj.pconj(kernel0.packet[0 % PacketSize])); |
| pstoreu(blockB + count + 1 * PacketSize, cj.pconj(kernel1.packet[0 % PacketSize])); |
| pstoreu(blockB + count + 2 * PacketSize, cj.pconj(kernel0.packet[1 % PacketSize])); |
| pstoreu(blockB + count + 3 * PacketSize, cj.pconj(kernel1.packet[1 % PacketSize])); |
| pstoreu(blockB + count + 4 * PacketSize, cj.pconj(kernel0.packet[2 % PacketSize])); |
| pstoreu(blockB + count + 5 * PacketSize, cj.pconj(kernel1.packet[2 % PacketSize])); |
| pstoreu(blockB + count + 6 * PacketSize, cj.pconj(kernel0.packet[3 % PacketSize])); |
| pstoreu(blockB + count + 7 * PacketSize, cj.pconj(kernel1.packet[3 % PacketSize])); |
| count += 8 * PacketSize; |
| } else if (PacketSize == 8) { |
| PacketBlock<Packet, PacketSize == 8 ? 8 : PacketSize> kernel0; |
| |
| kernel0.packet[0 % PacketSize] = dm0.template loadPacket<Packet>(k); |
| kernel0.packet[1 % PacketSize] = dm1.template loadPacket<Packet>(k); |
| kernel0.packet[2 % PacketSize] = dm2.template loadPacket<Packet>(k); |
| kernel0.packet[3 % PacketSize] = dm3.template loadPacket<Packet>(k); |
| kernel0.packet[4 % PacketSize] = dm4.template loadPacket<Packet>(k); |
| kernel0.packet[5 % PacketSize] = dm5.template loadPacket<Packet>(k); |
| kernel0.packet[6 % PacketSize] = dm6.template loadPacket<Packet>(k); |
| kernel0.packet[7 % PacketSize] = dm7.template loadPacket<Packet>(k); |
| ptranspose(kernel0); |
| |
| pstoreu(blockB + count + 0 * PacketSize, cj.pconj(kernel0.packet[0 % PacketSize])); |
| pstoreu(blockB + count + 1 * PacketSize, cj.pconj(kernel0.packet[1 % PacketSize])); |
| pstoreu(blockB + count + 2 * PacketSize, cj.pconj(kernel0.packet[2 % PacketSize])); |
| pstoreu(blockB + count + 3 * PacketSize, cj.pconj(kernel0.packet[3 % PacketSize])); |
| pstoreu(blockB + count + 4 * PacketSize, cj.pconj(kernel0.packet[4 % PacketSize])); |
| pstoreu(blockB + count + 5 * PacketSize, cj.pconj(kernel0.packet[5 % PacketSize])); |
| pstoreu(blockB + count + 6 * PacketSize, cj.pconj(kernel0.packet[6 % PacketSize])); |
| pstoreu(blockB + count + 7 * PacketSize, cj.pconj(kernel0.packet[7 % PacketSize])); |
| count += 8 * PacketSize; |
| } |
| } |
| } |
| |
| for (; k < depth; k++) { |
| blockB[count + 0] = cj(dm0(k)); |
| blockB[count + 1] = cj(dm1(k)); |
| blockB[count + 2] = cj(dm2(k)); |
| blockB[count + 3] = cj(dm3(k)); |
| blockB[count + 4] = cj(dm4(k)); |
| blockB[count + 5] = cj(dm5(k)); |
| blockB[count + 6] = cj(dm6(k)); |
| blockB[count + 7] = cj(dm7(k)); |
| count += 8; |
| } |
| // skip what we have after |
| if (PanelMode) count += 8 * (stride - offset - depth); |
| } |
| } |
| #endif |
| |
| EIGEN_IF_CONSTEXPR(nr >= 4) { |
| for (Index j2 = packet_cols8; j2 < packet_cols4; j2 += 4) { |
| // skip what we have before |
| if (PanelMode) count += 4 * offset; |
| const LinearMapper dm0 = rhs.getLinearMapper(0, j2 + 0); |
| const LinearMapper dm1 = rhs.getLinearMapper(0, j2 + 1); |
| const LinearMapper dm2 = rhs.getLinearMapper(0, j2 + 2); |
| const LinearMapper dm3 = rhs.getLinearMapper(0, j2 + 3); |
| |
| Index k = 0; |
| if ((PacketSize % 4) == 0) // TODO enable vectorized transposition for PacketSize==2 ?? |
| { |
| for (; k < peeled_k; k += PacketSize) { |
| PacketBlock<Packet, (PacketSize % 4) == 0 ? 4 : PacketSize> kernel; |
| kernel.packet[0] = dm0.template loadPacket<Packet>(k); |
| kernel.packet[1 % PacketSize] = dm1.template loadPacket<Packet>(k); |
| kernel.packet[2 % PacketSize] = dm2.template loadPacket<Packet>(k); |
| kernel.packet[3 % PacketSize] = dm3.template loadPacket<Packet>(k); |
| ptranspose(kernel); |
| pstoreu(blockB + count + 0 * PacketSize, cj.pconj(kernel.packet[0])); |
| pstoreu(blockB + count + 1 * PacketSize, cj.pconj(kernel.packet[1 % PacketSize])); |
| pstoreu(blockB + count + 2 * PacketSize, cj.pconj(kernel.packet[2 % PacketSize])); |
| pstoreu(blockB + count + 3 * PacketSize, cj.pconj(kernel.packet[3 % PacketSize])); |
| count += 4 * PacketSize; |
| } |
| } |
| for (; k < depth; k++) { |
| blockB[count + 0] = cj(dm0(k)); |
| blockB[count + 1] = cj(dm1(k)); |
| blockB[count + 2] = cj(dm2(k)); |
| blockB[count + 3] = cj(dm3(k)); |
| count += 4; |
| } |
| // skip what we have after |
| if (PanelMode) count += 4 * (stride - offset - depth); |
| } |
| } |
| |
| // copy the remaining columns one at a time (nr==1) |
| for (Index j2 = packet_cols4; j2 < cols; ++j2) { |
| if (PanelMode) count += offset; |
| const LinearMapper dm0 = rhs.getLinearMapper(0, j2); |
| for (Index k = 0; k < depth; k++) { |
| blockB[count] = cj(dm0(k)); |
| count += 1; |
| } |
| if (PanelMode) count += (stride - offset - depth); |
| } |
| } |
| |
| // this version is optimized for row major matrices |
| template <typename Scalar, typename Index, typename DataMapper, int nr, bool Conjugate, bool PanelMode> |
| struct gemm_pack_rhs<Scalar, Index, DataMapper, nr, RowMajor, Conjugate, PanelMode> { |
| typedef typename packet_traits<Scalar>::type Packet; |
| typedef typename unpacket_traits<Packet>::half HalfPacket; |
| typedef typename unpacket_traits<typename unpacket_traits<Packet>::half>::half QuarterPacket; |
| typedef typename DataMapper::LinearMapper LinearMapper; |
| enum { |
| PacketSize = packet_traits<Scalar>::size, |
| HalfPacketSize = unpacket_traits<HalfPacket>::size, |
| QuarterPacketSize = unpacket_traits<QuarterPacket>::size |
| }; |
| EIGEN_DONT_INLINE void operator()(Scalar* blockB, const DataMapper& rhs, Index depth, Index cols, Index stride = 0, |
| Index offset = 0) { |
| EIGEN_ASM_COMMENT("EIGEN PRODUCT PACK RHS ROWMAJOR"); |
| EIGEN_UNUSED_VARIABLE(stride); |
| EIGEN_UNUSED_VARIABLE(offset); |
| eigen_assert(((!PanelMode) && stride == 0 && offset == 0) || (PanelMode && stride >= depth && offset <= stride)); |
| const bool HasHalf = (int)HalfPacketSize < (int)PacketSize; |
| const bool HasQuarter = (int)QuarterPacketSize < (int)HalfPacketSize; |
| conj_if<NumTraits<Scalar>::IsComplex && Conjugate> cj; |
| Index packet_cols8 = nr >= 8 ? (cols / 8) * 8 : 0; |
| Index packet_cols4 = nr >= 4 ? (cols / 4) * 4 : 0; |
| Index count = 0; |
| |
| #if EIGEN_ARCH_ARM64 |
| EIGEN_IF_CONSTEXPR(nr >= 8) { |
| for (Index j2 = 0; j2 < packet_cols8; j2 += 8) { |
| // skip what we have before |
| if (PanelMode) count += 8 * offset; |
| for (Index k = 0; k < depth; k++) { |
| if (PacketSize == 8) { |
| Packet A = rhs.template loadPacket<Packet>(k, j2); |
| pstoreu(blockB + count, cj.pconj(A)); |
| count += PacketSize; |
| } else if (PacketSize == 4) { |
| Packet A = rhs.template loadPacket<Packet>(k, j2); |
| Packet B = rhs.template loadPacket<Packet>(k, j2 + 4); |
| pstoreu(blockB + count, cj.pconj(A)); |
| pstoreu(blockB + count + PacketSize, cj.pconj(B)); |
| count += 2 * PacketSize; |
| } else { |
| const LinearMapper dm0 = rhs.getLinearMapper(k, j2); |
| blockB[count + 0] = cj(dm0(0)); |
| blockB[count + 1] = cj(dm0(1)); |
| blockB[count + 2] = cj(dm0(2)); |
| blockB[count + 3] = cj(dm0(3)); |
| blockB[count + 4] = cj(dm0(4)); |
| blockB[count + 5] = cj(dm0(5)); |
| blockB[count + 6] = cj(dm0(6)); |
| blockB[count + 7] = cj(dm0(7)); |
| count += 8; |
| } |
| } |
| // skip what we have after |
| if (PanelMode) count += 8 * (stride - offset - depth); |
| } |
| } |
| #endif |
| |
| if (nr >= 4) { |
| for (Index j2 = packet_cols8; j2 < packet_cols4; j2 += 4) { |
| // skip what we have before |
| if (PanelMode) count += 4 * offset; |
| for (Index k = 0; k < depth; k++) { |
| if (PacketSize == 4) { |
| Packet A = rhs.template loadPacket<Packet>(k, j2); |
| pstoreu(blockB + count, cj.pconj(A)); |
| count += PacketSize; |
| } else if (HasHalf && HalfPacketSize == 4) { |
| HalfPacket A = rhs.template loadPacket<HalfPacket>(k, j2); |
| pstoreu(blockB + count, cj.pconj(A)); |
| count += HalfPacketSize; |
| } else if (HasQuarter && QuarterPacketSize == 4) { |
| QuarterPacket A = rhs.template loadPacket<QuarterPacket>(k, j2); |
| pstoreu(blockB + count, cj.pconj(A)); |
| count += QuarterPacketSize; |
| } else { |
| const LinearMapper dm0 = rhs.getLinearMapper(k, j2); |
| blockB[count + 0] = cj(dm0(0)); |
| blockB[count + 1] = cj(dm0(1)); |
| blockB[count + 2] = cj(dm0(2)); |
| blockB[count + 3] = cj(dm0(3)); |
| count += 4; |
| } |
| } |
| // skip what we have after |
| if (PanelMode) count += 4 * (stride - offset - depth); |
| } |
| } |
| // copy the remaining columns one at a time (nr==1) |
| for (Index j2 = packet_cols4; j2 < cols; ++j2) { |
| if (PanelMode) count += offset; |
| for (Index k = 0; k < depth; k++) { |
| blockB[count] = cj(rhs(k, j2)); |
| count += 1; |
| } |
| if (PanelMode) count += stride - offset - depth; |
| } |
| } |
| }; |
| |
| } // end namespace internal |
| |
| /** \returns the currently set level 1 cpu cache size (in bytes) used to estimate the ideal blocking size parameters. |
| * \sa setCpuCacheSize */ |
| inline std::ptrdiff_t l1CacheSize() { |
| std::ptrdiff_t l1, l2, l3; |
| internal::manage_caching_sizes(GetAction, &l1, &l2, &l3); |
| return l1; |
| } |
| |
| /** \returns the currently set level 2 cpu cache size (in bytes) used to estimate the ideal blocking size parameters. |
| * \sa setCpuCacheSize */ |
| inline std::ptrdiff_t l2CacheSize() { |
| std::ptrdiff_t l1, l2, l3; |
| internal::manage_caching_sizes(GetAction, &l1, &l2, &l3); |
| return l2; |
| } |
| |
| /** \returns the currently set level 3 cpu cache size (in bytes) used to estimate the ideal blocking size paramete\ |
| rs. |
| * \sa setCpuCacheSize */ |
| inline std::ptrdiff_t l3CacheSize() { |
| std::ptrdiff_t l1, l2, l3; |
| internal::manage_caching_sizes(GetAction, &l1, &l2, &l3); |
| return l3; |
| } |
| |
| /** Set the cpu L1 and L2 cache sizes (in bytes). |
| * These values are use to adjust the size of the blocks |
| * for the algorithms working per blocks. |
| * |
| * \sa computeProductBlockingSizes */ |
| inline void setCpuCacheSizes(std::ptrdiff_t l1, std::ptrdiff_t l2, std::ptrdiff_t l3) { |
| internal::manage_caching_sizes(SetAction, &l1, &l2, &l3); |
| } |
| |
| } // end namespace Eigen |
| |
| #endif // EIGEN_GENERAL_BLOCK_PANEL_H |
