| // 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(
|