使用AVX/AVX2转换8x8浮点

我已经回答了这个问题 Fast memory transpose with SSE, AVX, and OpenMP .

让我重复使用AVX转换8x8浮点矩阵的解决方案。让我知道这是否比使用4x4块和 _MM_TRANSPOSE4_PS 我将其用于更大矩阵转置中的内核,这是内存受限的,因此这可能不是一个公平的测试。

inline void transpose8_ps(__m256 &row0, __m256 &row1, __m256 &row2, __m256 &row3, __m256 &row4, __m256 &row5, __m256 &row6, __m256 &row7) {
__m256 __t0, __t1, __t2, __t3, __t4, __t5, __t6, __t7;
__m256 __tt0, __tt1, __tt2, __tt3, __tt4, __tt5, __tt6, __tt7;
__t0 = _mm256_unpacklo_ps(row0, row1);
__t1 = _mm256_unpackhi_ps(row0, row1);
__t2 = _mm256_unpacklo_ps(row2, row3);
__t3 = _mm256_unpackhi_ps(row2, row3);
__t4 = _mm256_unpacklo_ps(row4, row5);
__t5 = _mm256_unpackhi_ps(row4, row5);
__t6 = _mm256_unpacklo_ps(row6, row7);
__t7 = _mm256_unpackhi_ps(row6, row7);
__tt0 = _mm256_shuffle_ps(__t0,__t2,_MM_SHUFFLE(1,0,1,0));
__tt1 = _mm256_shuffle_ps(__t0,__t2,_MM_SHUFFLE(3,2,3,2));
__tt2 = _mm256_shuffle_ps(__t1,__t3,_MM_SHUFFLE(1,0,1,0));
__tt3 = _mm256_shuffle_ps(__t1,__t3,_MM_SHUFFLE(3,2,3,2));
__tt4 = _mm256_shuffle_ps(__t4,__t6,_MM_SHUFFLE(1,0,1,0));
__tt5 = _mm256_shuffle_ps(__t4,__t6,_MM_SHUFFLE(3,2,3,2));
__tt6 = _mm256_shuffle_ps(__t5,__t7,_MM_SHUFFLE(1,0,1,0));
__tt7 = _mm256_shuffle_ps(__t5,__t7,_MM_SHUFFLE(3,2,3,2));
row0 = _mm256_permute2f128_ps(__tt0, __tt4, 0x20);
row1 = _mm256_permute2f128_ps(__tt1, __tt5, 0x20);
row2 = _mm256_permute2f128_ps(__tt2, __tt6, 0x20);
row3 = _mm256_permute2f128_ps(__tt3, __tt7, 0x20);
row4 = _mm256_permute2f128_ps(__tt0, __tt4, 0x31);
row5 = _mm256_permute2f128_ps(__tt1, __tt5, 0x31);
row6 = _mm256_permute2f128_ps(__tt2, __tt6, 0x31);
row7 = _mm256_permute2f128_ps(__tt3, __tt7, 0x31);
}

基于 this comment 我了解到有更有效的方法来进行8x8转置。参见 Intel optimization manual 根据“11.11处理端口5压力”一节。实施例11-19使用相同数量的指令,但通过使用也到端口0的混合来减少端口5上的压力。我可能会在某个时候用intrinsic实现这一点,但现在我不需要这样做。

我在上面提到的英特尔手册中更仔细地查看了示例11-19和11-20。事实证明,示例11-19使用了4个多余的洗牌操作。它有8个解包、12个洗牌和8个128位排列。我的方法少用4次洗牌。他们用混合物替换了8个洗牌。因此,4次洗牌和8次混合。我怀疑这比我只有八次洗牌的方法更好。

然而,如果需要从内存加载矩阵,则示例11-20是一个改进。这使用8个解包、8个插入、8个洗牌、8个128位加载和8个存储。128位负载降低了端口压力。我继续使用内部函数实现了这一点。

//Example 11-20. 8x8 Matrix Transpose Using VINSERTF128 loads
void tran(float* mat, float* matT) {
  __m256  r0, r1, r2, r3, r4, r5, r6, r7;
  __m256  t0, t1, t2, t3, t4, t5, t6, t7;

  r0 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[0*8+0])), _mm_load_ps(&mat[4*8+0]), 1);
  r1 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[1*8+0])), _mm_load_ps(&mat[5*8+0]), 1);
  r2 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[2*8+0])), _mm_load_ps(&mat[6*8+0]), 1);
  r3 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[3*8+0])), _mm_load_ps(&mat[7*8+0]), 1);
  r4 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[0*8+4])), _mm_load_ps(&mat[4*8+4]), 1);
  r5 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[1*8+4])), _mm_load_ps(&mat[5*8+4]), 1);
  r6 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[2*8+4])), _mm_load_ps(&mat[6*8+4]), 1);
  r7 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[3*8+4])), _mm_load_ps(&mat[7*8+4]), 1);

  t0 = _mm256_unpacklo_ps(r0,r1);
  t1 = _mm256_unpackhi_ps(r0,r1);
  t2 = _mm256_unpacklo_ps(r2,r3);
  t3 = _mm256_unpackhi_ps(r2,r3);
  t4 = _mm256_unpacklo_ps(r4,r5);
  t5 = _mm256_unpackhi_ps(r4,r5);
  t6 = _mm256_unpacklo_ps(r6,r7);
  t7 = _mm256_unpackhi_ps(r6,r7);

  r0 = _mm256_shuffle_ps(t0,t2, 0x44);
  r1 = _mm256_shuffle_ps(t0,t2, 0xEE);
  r2 = _mm256_shuffle_ps(t1,t3, 0x44);
  r3 = _mm256_shuffle_ps(t1,t3, 0xEE);
  r4 = _mm256_shuffle_ps(t4,t6, 0x44);
  r5 = _mm256_shuffle_ps(t4,t6, 0xEE);
  r6 = _mm256_shuffle_ps(t5,t7, 0x44);
  r7 = _mm256_shuffle_ps(t5,t7, 0xEE);

  _mm256_store_ps(&matT[0*8], r0);
  _mm256_store_ps(&matT[1*8], r1);
  _mm256_store_ps(&matT[2*8], r2);
  _mm256_store_ps(&matT[3*8], r3);
  _mm256_store_ps(&matT[4*8], r4);
  _mm256_store_ps(&matT[5*8], r5);
  _mm256_store_ps(&matT[6*8], r6);
  _mm256_store_ps(&matT[7*8], r7);
}

所以我再次研究了示例11-19。据我所知,基本的想法是两个混洗指令(shufps)可以被一个混洗和两个混合指令所取代。例如

r0 = _mm256_shuffle_ps(t0,t2, 0x44);
r1 = _mm256_shuffle_ps(t0,t2, 0xEE);

可以替换为

v = _mm256_shuffle_ps(t0,t2, 0x4E);
r0 = _mm256_blend_ps(t0, v, 0xCC);
r1 = _mm256_blend_ps(t2, v, 0x33);

这解释了为什么我的原始代码使用了8次洗牌,而示例11-19使用了4次洗牌和8次混合。

混合对吞吐量有好处,因为混洗只会到达一个端口(在混洗端口上造成瓶颈),但混合可以在多个端口上运行,因此不会竞争。但什么更好:8次洗牌或4次洗牌和8次混合?

这必须经过测试,并且可能取决于周围的代码。如果您使用 大量 对于循环中不需要端口5的其他uop,您可能会选择纯shuffle版本。理想情况下,您应该在存储转置数据之前对其进行一些计算,因为它已经在寄存器中。看见 https://agner.org/optimize/ 和其他性能链接 the x86 tag wiki .

然而,我看不出用混合物代替解包说明的方法。

这里是结合了示例11-19将2个混洗转换为1个混洗和两个混合的完整代码,以及示例11-20 vinsertf128 负载(Intel Haswell/Skylake CPU上的负载是2 uops:任何端口一个ALU,一个内存。不幸的是,它们没有微型熔断器。 葡萄酒128 对于英特尔上的随机端口,所有寄存器操作数都是1uop,所以这是很好的,因为编译器会将负载折叠到内存操作数中 葡萄酒128 .)这样做的优点是,只需要源数据16字节对齐即可获得最大性能,避免任何缓存线拆分。

#include <stdio.h>
#include <x86intrin.h>
#include <omp.h>   

/*
void tran(float* mat, float* matT) {
  __m256  r0, r1, r2, r3, r4, r5, r6, r7;
  __m256  t0, t1, t2, t3, t4, t5, t6, t7;

  r0 = _mm256_load_ps(&mat[0*8]);
  r1 = _mm256_load_ps(&mat[1*8]);
  r2 = _mm256_load_ps(&mat[2*8]);
  r3 = _mm256_load_ps(&mat[3*8]);
  r4 = _mm256_load_ps(&mat[4*8]);
  r5 = _mm256_load_ps(&mat[5*8]);
  r6 = _mm256_load_ps(&mat[6*8]);
  r7 = _mm256_load_ps(&mat[7*8]);

  t0 = _mm256_unpacklo_ps(r0, r1);
  t1 = _mm256_unpackhi_ps(r0, r1);
  t2 = _mm256_unpacklo_ps(r2, r3);
  t3 = _mm256_unpackhi_ps(r2, r3);
  t4 = _mm256_unpacklo_ps(r4, r5);
  t5 = _mm256_unpackhi_ps(r4, r5);
  t6 = _mm256_unpacklo_ps(r6, r7);
  t7 = _mm256_unpackhi_ps(r6, r7);

  r0 = _mm256_shuffle_ps(t0,t2,_MM_SHUFFLE(1,0,1,0));  
  r1 = _mm256_shuffle_ps(t0,t2,_MM_SHUFFLE(3,2,3,2));
  r2 = _mm256_shuffle_ps(t1,t3,_MM_SHUFFLE(1,0,1,0));
  r3 = _mm256_shuffle_ps(t1,t3,_MM_SHUFFLE(3,2,3,2));
  r4 = _mm256_shuffle_ps(t4,t6,_MM_SHUFFLE(1,0,1,0));
  r5 = _mm256_shuffle_ps(t4,t6,_MM_SHUFFLE(3,2,3,2));
  r6 = _mm256_shuffle_ps(t5,t7,_MM_SHUFFLE(1,0,1,0));
  r7 = _mm256_shuffle_ps(t5,t7,_MM_SHUFFLE(3,2,3,2));

  t0 = _mm256_permute2f128_ps(r0, r4, 0x20);
  t1 = _mm256_permute2f128_ps(r1, r5, 0x20);
  t2 = _mm256_permute2f128_ps(r2, r6, 0x20);
  t3 = _mm256_permute2f128_ps(r3, r7, 0x20);
  t4 = _mm256_permute2f128_ps(r0, r4, 0x31);
  t5 = _mm256_permute2f128_ps(r1, r5, 0x31);
  t6 = _mm256_permute2f128_ps(r2, r6, 0x31);
  t7 = _mm256_permute2f128_ps(r3, r7, 0x31);

  _mm256_store_ps(&matT[0*8], t0);
  _mm256_store_ps(&matT[1*8], t1);
  _mm256_store_ps(&matT[2*8], t2);
  _mm256_store_ps(&matT[3*8], t3);
  _mm256_store_ps(&matT[4*8], t4);
  _mm256_store_ps(&matT[5*8], t5);
  _mm256_store_ps(&matT[6*8], t6);
  _mm256_store_ps(&matT[7*8], t7);
}
*/

void tran(float* mat, float* matT) {
  __m256  r0, r1, r2, r3, r4, r5, r6, r7;
  __m256  t0, t1, t2, t3, t4, t5, t6, t7;

  r0 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[0*8+0])), _mm_load_ps(&mat[4*8+0]), 1);
  r1 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[1*8+0])), _mm_load_ps(&mat[5*8+0]), 1);
  r2 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[2*8+0])), _mm_load_ps(&mat[6*8+0]), 1);
  r3 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[3*8+0])), _mm_load_ps(&mat[7*8+0]), 1);
  r4 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[0*8+4])), _mm_load_ps(&mat[4*8+4]), 1);
  r5 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[1*8+4])), _mm_load_ps(&mat[5*8+4]), 1);
  r6 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[2*8+4])), _mm_load_ps(&mat[6*8+4]), 1);
  r7 = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm_load_ps(&mat[3*8+4])), _mm_load_ps(&mat[7*8+4]), 1);

  t0 = _mm256_unpacklo_ps(r0,r1);
  t1 = _mm256_unpackhi_ps(r0,r1);
  t2 = _mm256_unpacklo_ps(r2,r3);
  t3 = _mm256_unpackhi_ps(r2,r3);
  t4 = _mm256_unpacklo_ps(r4,r5);
  t5 = _mm256_unpackhi_ps(r4,r5);
  t6 = _mm256_unpacklo_ps(r6,r7);
  t7 = _mm256_unpackhi_ps(r6,r7);

  __m256 v;

  //r0 = _mm256_shuffle_ps(t0,t2, 0x44);
  //r1 = _mm256_shuffle_ps(t0,t2, 0xEE);  
  v = _mm256_shuffle_ps(t0,t2, 0x4E);
  r0 = _mm256_blend_ps(t0, v, 0xCC);
  r1 = _mm256_blend_ps(t2, v, 0x33);

  //r2 = _mm256_shuffle_ps(t1,t3, 0x44);
  //r3 = _mm256_shuffle_ps(t1,t3, 0xEE);
  v = _mm256_shuffle_ps(t1,t3, 0x4E);
  r2 = _mm256_blend_ps(t1, v, 0xCC);
  r3 = _mm256_blend_ps(t3, v, 0x33);

  //r4 = _mm256_shuffle_ps(t4,t6, 0x44);
  //r5 = _mm256_shuffle_ps(t4,t6, 0xEE);
  v = _mm256_shuffle_ps(t4,t6, 0x4E);
  r4 = _mm256_blend_ps(t4, v, 0xCC);
  r5 = _mm256_blend_ps(t6, v, 0x33);

  //r6 = _mm256_shuffle_ps(t5,t7, 0x44);
  //r7 = _mm256_shuffle_ps(t5,t7, 0xEE);
  v = _mm256_shuffle_ps(t5,t7, 0x4E);
  r6 = _mm256_blend_ps(t5, v, 0xCC);
  r7 = _mm256_blend_ps(t7, v, 0x33);

  _mm256_store_ps(&matT[0*8], r0);
  _mm256_store_ps(&matT[1*8], r1);
  _mm256_store_ps(&matT[2*8], r2);
  _mm256_store_ps(&matT[3*8], r3);
  _mm256_store_ps(&matT[4*8], r4);
  _mm256_store_ps(&matT[5*8], r5);
  _mm256_store_ps(&matT[6*8], r6);
  _mm256_store_ps(&matT[7*8], r7);
}


int verify(float *mat) {
    int i,j;
    int error = 0;
    for(i=0; i<8; i++) {
      for(j=0; j<8; j++) {
        if(mat[j*8+i] != 1.0f*i*8+j) error++;
      }
    }
    return error;
}

void print_mat(float *mat) {
    int i,j;
    for(i=0; i<8; i++) {
      for(j=0; j<8; j++) printf("%2.0f ", mat[i*8+j]);
      puts("");
    }
    puts("");
}

int main(void) {
    int i,j, rep;
    float mat[64] __attribute__((aligned(64)));
    float matT[64] __attribute__((aligned(64)));
    double dtime;

    rep = 10000000;
    for(i=0; i<64; i++) mat[i] = i;
    print_mat(mat);

    tran(mat,matT);
    //dtime = -omp_get_wtime();
    //tran(mat, matT, rep);
    //dtime += omp_get_wtime();
    printf("errors %d\n", verify(matT));
    //printf("dtime %f\n", dtime);
    print_mat(matT);
}

这是一个适用于8 x 8 32位整数的AVX2解决方案。当然,如果您想调换8 x 8个浮点值,可以将浮点向量转换为int和back。也有可能只针对浮点运算执行AVX版本(即不需要AVX2),但我还没有尝试过。

//
// tranpose_8_8_avx2.c
//

#include <stdio.h>

#include <immintrin.h>

#define V_ELEMS 8

static inline void _mm256_merge_epi32(const __m256i v0, const __m256i v1, __m256i *vl, __m256i *vh)
{
    __m256i va = _mm256_permute4x64_epi64(v0, _MM_SHUFFLE(3, 1, 2, 0));
    __m256i vb = _mm256_permute4x64_epi64(v1, _MM_SHUFFLE(3, 1, 2, 0));
    *vl = _mm256_unpacklo_epi32(va, vb);
    *vh = _mm256_unpackhi_epi32(va, vb);
}

static inline void _mm256_merge_epi64(const __m256i v0, const __m256i v1, __m256i *vl, __m256i *vh)
{
    __m256i va = _mm256_permute4x64_epi64(v0, _MM_SHUFFLE(3, 1, 2, 0));
    __m256i vb = _mm256_permute4x64_epi64(v1, _MM_SHUFFLE(3, 1, 2, 0));
    *vl = _mm256_unpacklo_epi64(va, vb);
    *vh = _mm256_unpackhi_epi64(va, vb);
}

static inline void _mm256_merge_si128(const __m256i v0, const __m256i v1, __m256i *vl, __m256i *vh)
{
    *vl = _mm256_permute2x128_si256(v0, v1, _MM_SHUFFLE(0, 2, 0, 0));
    *vh = _mm256_permute2x128_si256(v0, v1, _MM_SHUFFLE(0, 3, 0, 1));
}

//
// Transpose_8_8
//
// in place transpose of 8 x 8 int array
//

static void Transpose_8_8(
    __m256i *v0,
    __m256i *v1,
    __m256i *v2,
    __m256i *v3,
    __m256i *v4,
    __m256i *v5,
    __m256i *v6,
    __m256i *v7)
{
    __m256i w0, w1, w2, w3, w4, w5, w6, w7;
    __m256i x0, x1, x2, x3, x4, x5, x6, x7;

    _mm256_merge_epi32(*v0, *v1, &w0, &w1);
    _mm256_merge_epi32(*v2, *v3, &w2, &w3);
    _mm256_merge_epi32(*v4, *v5, &w4, &w5);
    _mm256_merge_epi32(*v6, *v7, &w6, &w7);

    _mm256_merge_epi64(w0, w2, &x0, &x1);
    _mm256_merge_epi64(w1, w3, &x2, &x3);
    _mm256_merge_epi64(w4, w6, &x4, &x5);
    _mm256_merge_epi64(w5, w7, &x6, &x7);

    _mm256_merge_si128(x0, x4, v0, v1);
    _mm256_merge_si128(x1, x5, v2, v3);
    _mm256_merge_si128(x2, x6, v4, v5);
    _mm256_merge_si128(x3, x7, v6, v7);
}

int main(void)
{
    int32_t buff[V_ELEMS][V_ELEMS] __attribute__ ((aligned(32)));
    int i, j;
    int k = 0;

    // init buff
    for (i = 0; i < V_ELEMS; ++i)
    {
        for (j = 0; j < V_ELEMS; ++j)
        {
            buff[i][j] = k++;
        }
    }

    // print buff
    printf("\nBEFORE:\n");
    for (i = 0; i < V_ELEMS; ++i)
    {
        for (j = 0; j < V_ELEMS; ++j)
        {
            printf("%4d", buff[i][j]);
        }
        printf("\n");
    }

    // transpose
    Transpose_8_8((__m256i *)buff[0], (__m256i *)buff[1], (__m256i *)buff[2], (__m256i *)buff[3], (__m256i *)buff[4], (__m256i *)buff[5], (__m256i *)buff[6], (__m256i *)buff[7]);

    // print buff
    printf("\nAFTER:\n");
    for (i = 0; i < V_ELEMS; ++i)
    {
        for (j = 0; j < V_ELEMS; ++j)
        {
            printf("%4d", buff[i][j]);
        }
        printf("\n");
    }

    // transpose
    Transpose_8_8((__m256i *)buff[0], (__m256i *)buff[1], (__m256i *)buff[2], (__m256i *)buff[3], (__m256i *)buff[4], (__m256i *)buff[5], (__m256i *)buff[6], (__m256i *)buff[7]);

    // print buff
    printf("\nAFTER x2:\n");
    for (i = 0; i < V_ELEMS; ++i)
    {
        for (j = 0; j < V_ELEMS; ++j)
        {
            printf("%4d", buff[i][j]);
        }
        printf("\n");
    }

    return 0;
}

Transpose_8_8 编译到大约56条指令时发出叮当声,包括加载和存储。

编译和测试:

$ gcc -Wall -mavx2 -O3 transpose_8_8_avx2.c && ./a.out

BEFORE:
   0   1   2   3   4   5   6   7
   8   9  10  11  12  13  14  15
  16  17  18  19  20  21  22  23
  24  25  26  27  28  29  30  31
  32  33  34  35  36  37  38  39
  40  41  42  43  44  45  46  47
  48  49  50  51  52  53  54  55
  56  57  58  59  60  61  62  63

AFTER:
   0   8  16  24  32  40  48  56
   1   9  17  25  33  41  49  57
   2  10  18  26  34  42  50  58
   3  11  19  27  35  43  51  59
   4  12  20  28  36  44  52  60
   5  13  21  29  37  45  53  61
   6  14  22  30  38  46  54  62
   7  15  23  31  39  47  55  63

AFTER x2:
   0   1   2   3   4   5   6   7
   8   9  10  11  12  13  14  15
  16  17  18  19  20  21  22  23
  24  25  26  27  28  29  30  31
  32  33  34  35  36  37  38  39
  40  41  42  43  44  45  46  47
  48  49  50  51  52  53  54  55
  56  57  58  59  60  61  62  63

$ 

我决定在一个苹果对苹果的比较中,对3种不同的程序进行全面测试。

  // transpose.cpp : 
  /*
  Transpose8x8Shuff 100,000,000 times
  (0.750000 seconds).
  Transpose8x8Permute 100,000,000 times
  (0.749000 seconds).
  Transpose8x8Insert 100,000,000 times
  (0.858000 seconds).
  */

  #include <stdio.h>
  #include <time.h>
  #include <thread>
  #include <chrono>
  #include <xmmintrin.h>
  #include <emmintrin.h>
  #include <tmmintrin.h>
  #include <immintrin.h>

  inline void Transpose8x8Shuff(unsigned long *in)
  {
     __m256 *inI = reinterpret_cast<__m256 *>(in);
     __m256 rI[8];
     rI[0] = _mm256_unpacklo_ps(inI[0], inI[1]); 
     rI[1] = _mm256_unpackhi_ps(inI[0], inI[1]); 
     rI[2] = _mm256_unpacklo_ps(inI[2], inI[3]); 
     rI[3] = _mm256_unpackhi_ps(inI[2], inI[3]); 
     rI[4] = _mm256_unpacklo_ps(inI[4], inI[5]); 
     rI[5] = _mm256_unpackhi_ps(inI[4], inI[5]); 
     rI[6] = _mm256_unpacklo_ps(inI[6], inI[7]); 
     rI[7] = _mm256_unpackhi_ps(inI[6], inI[7]); 

     __m256 rrF[8];
     __m256 *rF = reinterpret_cast<__m256 *>(rI);
     rrF[0] = _mm256_shuffle_ps(rF[0], rF[2], _MM_SHUFFLE(1,0,1,0));
     rrF[1] = _mm256_shuffle_ps(rF[0], rF[2], _MM_SHUFFLE(3,2,3,2));
     rrF[2] = _mm256_shuffle_ps(rF[1], rF[3], _MM_SHUFFLE(1,0,1,0)); 
     rrF[3] = _mm256_shuffle_ps(rF[1], rF[3], _MM_SHUFFLE(3,2,3,2));
     rrF[4] = _mm256_shuffle_ps(rF[4], rF[6], _MM_SHUFFLE(1,0,1,0));
     rrF[5] = _mm256_shuffle_ps(rF[4], rF[6], _MM_SHUFFLE(3,2,3,2));
     rrF[6] = _mm256_shuffle_ps(rF[5], rF[7], _MM_SHUFFLE(1,0,1,0));
     rrF[7] = _mm256_shuffle_ps(rF[5], rF[7], _MM_SHUFFLE(3,2,3,2));

     rF = reinterpret_cast<__m256 *>(in);
     rF[0] = _mm256_permute2f128_ps(rrF[0], rrF[4], 0x20);
     rF[1] = _mm256_permute2f128_ps(rrF[1], rrF[5], 0x20);
     rF[2] = _mm256_permute2f128_ps(rrF[2], rrF[6], 0x20);
     rF[3] = _mm256_permute2f128_ps(rrF[3], rrF[7], 0x20);
     rF[4] = _mm256_permute2f128_ps(rrF[0], rrF[4], 0x31);
     rF[5] = _mm256_permute2f128_ps(rrF[1], rrF[5], 0x31);
     rF[6] = _mm256_permute2f128_ps(rrF[2], rrF[6], 0x31);
     rF[7] = _mm256_permute2f128_ps(rrF[3], rrF[7], 0x31);
  }

  inline void Transpose8x8Permute(unsigned long *in)
  {
     __m256i *inI = reinterpret_cast<__m256i *>(in);
     __m256i rI[8];
     rI[0] = _mm256_permute2f128_si256(inI[0], inI[4], 0x20); 
     rI[1] = _mm256_permute2f128_si256(inI[0], inI[4], 0x31); 
     rI[2] = _mm256_permute2f128_si256(inI[1], inI[5], 0x20); 
     rI[3] = _mm256_permute2f128_si256(inI[1], inI[5], 0x31); 
     rI[4] = _mm256_permute2f128_si256(inI[2], inI[6], 0x20); 
     rI[5] = _mm256_permute2f128_si256(inI[2], inI[6], 0x31); 
     rI[6] = _mm256_permute2f128_si256(inI[3], inI[7], 0x20); 
     rI[7] = _mm256_permute2f128_si256(inI[3], inI[7], 0x31); 

     __m256 rrF[8];
     __m256 *rF = reinterpret_cast<__m256 *>(rI);
     rrF[0] = _mm256_unpacklo_ps(rF[0], rF[4]); 
     rrF[1] = _mm256_unpackhi_ps(rF[0], rF[4]); 
     rrF[2] = _mm256_unpacklo_ps(rF[1], rF[5]); 
     rrF[3] = _mm256_unpackhi_ps(rF[1], rF[5]); 
     rrF[4] = _mm256_unpacklo_ps(rF[2], rF[6]); 
     rrF[5] = _mm256_unpackhi_ps(rF[2], rF[6]); 
     rrF[6] = _mm256_unpacklo_ps(rF[3], rF[7]); 
     rrF[7] = _mm256_unpackhi_ps(rF[3], rF[7]); 

     rF = reinterpret_cast<__m256 *>(in);
     rF[0] = _mm256_unpacklo_ps(rrF[0], rrF[4]); 
     rF[1] = _mm256_unpackhi_ps(rrF[0], rrF[4]); 
     rF[2] = _mm256_unpacklo_ps(rrF[1], rrF[5]); 
     rF[3] = _mm256_unpackhi_ps(rrF[1], rrF[5]); 
     rF[4] = _mm256_unpacklo_ps(rrF[2], rrF[6]); 
     rF[5] = _mm256_unpackhi_ps(rrF[2], rrF[6]); 
     rF[6] = _mm256_unpacklo_ps(rrF[3], rrF[7]); 
     rF[7] = _mm256_unpackhi_ps(rrF[3], rrF[7]); 
  }

  inline void Transpose8x8Insert(unsigned long *in)
  {
     __m256i *inIa = reinterpret_cast<__m256i *>(in);
     __m256i *inIb = reinterpret_cast<__m256i *>(&reinterpret_cast<__m128i *>(in)[1]);
     __m128i *inI128 = reinterpret_cast<__m128i *>(in);
     __m256i rI[8];
     rI[0] = _mm256_insertf128_si256(inIa[0], inI128[8], 1); 
     rI[1] = _mm256_insertf128_si256(inIb[0], inI128[9], 1); 
     rI[2] = _mm256_insertf128_si256(inIa[1], inI128[10], 1); 
     rI[3] = _mm256_insertf128_si256(inIb[1], inI128[11], 1); 
     rI[4] = _mm256_insertf128_si256(inIa[2], inI128[12], 1); 
     rI[5] = _mm256_insertf128_si256(inIb[2], inI128[13], 1); 
     rI[6] = _mm256_insertf128_si256(inIa[3], inI128[14], 1); 
     rI[7] = _mm256_insertf128_si256(inIb[3], inI128[15], 1); 

     __m256 rrF[8];
     __m256 *rF = reinterpret_cast<__m256 *>(rI);
     rrF[0] = _mm256_unpacklo_ps(rF[0], rF[4]); 
     rrF[1] = _mm256_unpackhi_ps(rF[0], rF[4]); 
     rrF[2] = _mm256_unpacklo_ps(rF[1], rF[5]); 
     rrF[3] = _mm256_unpackhi_ps(rF[1], rF[5]); 
     rrF[4] = _mm256_unpacklo_ps(rF[2], rF[6]); 
     rrF[5] = _mm256_unpackhi_ps(rF[2], rF[6]); 
     rrF[6] = _mm256_unpacklo_ps(rF[3], rF[7]); 
     rrF[7] = _mm256_unpackhi_ps(rF[3], rF[7]); 

     rF = reinterpret_cast<__m256 *>(in);
     rF[0] = _mm256_unpacklo_ps(rrF[0], rrF[4]); 
     rF[1] = _mm256_unpackhi_ps(rrF[0], rrF[4]); 
     rF[2] = _mm256_unpacklo_ps(rrF[1], rrF[5]); 
     rF[3] = _mm256_unpackhi_ps(rrF[1], rrF[5]); 
     rF[4] = _mm256_unpacklo_ps(rrF[2], rrF[6]); 
     rF[5] = _mm256_unpackhi_ps(rrF[2], rrF[6]); 
     rF[6] = _mm256_unpacklo_ps(rrF[3], rrF[7]); 
     rF[7] = _mm256_unpackhi_ps(rrF[3], rrF[7]); 
  }

  int main()
  {
  #define dwordCount 64
     alignas(32) unsigned long mat[dwordCount];
     for (int i = 0; i < dwordCount; i++) {
        mat[i] = i;
     }
     clock_t t;
     printf ("Transpose8x8Shuff 100,000,000 times\n");
     Transpose8x8Shuff(mat);
     t = clock();
     int i = 100000000;
     do {
        Transpose8x8Shuff(mat);
     } while (--i >= 0);
     t = clock() - t;
     volatile int dummy = mat[2];
     printf ("(%f seconds).\n",((float)t)/CLOCKS_PER_SEC);
     printf ("Transpose8x8Permute 100,000,000 times\n");
     Transpose8x8Permute(mat);
     t = clock();
     i = 100000000;
     do {
        Transpose8x8Permute(mat);
     } while (--i >= 0);
     t = clock() - t;
     volatile int dummy = mat[2];
     printf ("(%f seconds).\n",((float)t)/CLOCKS_PER_SEC);
     printf ("Transpose8x8Insert 100,000,000 times\n");
     Transpose8x8Insert(mat);
     t = clock();
     i = 100000000;
     do {
        Transpose8x8Insert(mat);
     } while (--i >= 0);
     t = clock() - t;
     volatile int dummy = mat[2];
     printf ("(%f seconds).\n",((float)t)/CLOCKS_PER_SEC);
     char c = getchar(); 
     return 0;
  }

这是基准测试延迟,而不是吞吐量(因为一个转置的输出是下一个转座的输入),但无论如何,这可能会对混洗吞吐量造成瓶颈。

在Skylake i7-6700k@3.9GHz上对上述代码的修改版本的结果(参见 on the Godbolt compiler explorer ),修复了以下错误:

• printf 在定时区域之外,启动 clock()
• volatile dummy = in[2] 最后,所有转置都不会优化掉(gcc实际上是这样做的)。
• 便携式C++11,不需要MSVC( alignas(32) 而不是 __declspec ,不包括 stdafx.h .)
• 睡眠被移除,所以CPU不会在测试之间降低到空闲速度。

我没有解决不必要的混合 __m256i* / __m256* ,我没有检查这是否导致gcc或clang的代码生成更差。我也没有使用 std::chrono 高可靠性时钟,因为 时钟() 对于如此多的重复来说足够准确。

克++7.3 -O3 -march=native 在ArchLinux:ZBoson的版本是最快的

Transpose8x8Shuff 100,000,000 times
(0.637479 seconds).
Transpose8x8Permute 100,000,000 times
(0.792658 seconds).
Transpose8x8Insert 100,000,000 times
(0.846590 seconds).

叮当声++5.0.1 -O3-三月=本地 :8x8Permute优化到比gcc更快的速度,但8x8Insert令人讨厌。

Transpose8x8Shuff 100,000,000 times
(0.642185 seconds).
Transpose8x8Permute 100,000,000 times
(0.622157 seconds).
Transpose8x8Insert 100,000,000 times
(2.149958 seconds).

从源代码生成的asm指令与内部函数不完全匹配:尤其是clang有一个shuffle优化器 编译 洗牌的方式与优化标量代码的方式相同,如 + 关于整数。 Transpose8x8Insert 应该不会慢那么多,所以clang肯定选得不好。

除了前面的答案之外,在这种情况下,shuffleps的使用非常过分,因为我们可以通过解包/解包的方式获得结果。洗牌&解包指令具有相同的延迟/吞吐量,但是shuffle在机器代码op中生成一个额外的字节(即,5字节用于shuffle,4字节用于解包)。

在 一些 点,我们需要8次穿越车道。这是一个较慢的操作(延迟为3个周期),因此我们希望尽早启动这些操作。假设转置8f方法被内联(它应该这样做!),那么a->h参数应融合到开箱说明中。

您可能面临的唯一一个小问题是,因为您在这里使用了8个以上的寄存器,您可能会溢出到YMM9及以上。这会导致VEX2操作被生成为VEX3,这将为每个操作添加一个字节。

结果,稍微晃动一下,你会得到这样的结果:

typedef __m256 f256;
#define unpacklo8f _mm256_unpacklo_ps
#define unpackhi8f _mm256_unpackhi_ps

template<uint8_t X, uint8_t Y>
inline f256 permute128f(const f256 a, const f256 b)
{
  return _mm256_permute2f128_ps(a, b, X | (Y << 4)); 
}

inline void transpose8f(
  const f256 a, const f256 b, const f256 c, const f256 d, 
  const f256 e, const f256 f, const f256 g, const f256 h, 
  f256& s, f256& t, f256& u, f256& v,
  f256& x, f256& y, f256& z, f256& w)
{
  const f256 t00 = unpacklo8f(a, c);
  const f256 t02 = unpacklo8f(b, d);
  const f256 t20 = unpacklo8f(e, g);
  const f256 t22 = unpacklo8f(f, h);

  const f256 t10 = unpacklo8f(t00, t02);
  const f256 t30 = unpacklo8f(t20, t22);
  const f256 t11 = unpackhi8f(t00, t02);
  s = permute128f<0, 2>(t10, t30);
  const f256 t31 = unpackhi8f(t20, t22);
  x = permute128f<1, 3>(t10, t30);
  const f256 t01 = unpackhi8f(a, c);
  t = permute128f<0, 2>(t11, t31);
  const f256 t21 = unpackhi8f(e, g);
  y = permute128f<1, 3>(t11, t31);
  const f256 t03 = unpackhi8f(b, d);
  const f256 t23 = unpackhi8f(f, h);

  const f256 t12 = unpacklo8f(t01, t03);
  const f256 t13 = unpackhi8f(t01, t03);
  const f256 t32 = unpacklo8f(t21, t23);
  const f256 t33 = unpackhi8f(t21, t23);

  u = permute128f<0, 2>(t12, t32);
  z = permute128f<1, 3>(t12, t32);
  v = permute128f<0, 2>(t13, t33);
  w = permute128f<1, 3>(t13, t33);
}

你不会改进的 (您可以先执行128位置换,然后再执行解包,但它们最终会完全相同) .

这是我的解决方案,指令更少,性能非常好,大约快8倍。我已经在Fedora中使用ICC、GCC和Clang进行了测试。

#include <stdio.h>
#include <x86intrin.h>
#define MAX1    128
#define MAX2 MAX1

float __attribute__(( aligned(32))) a_tra[MAX2][MAX1], __attribute__(( aligned(32))) a[MAX1][MAX2] ;
int main()
{

    int i,j;//, ii=0,jj=0;
    // variables for vector section
    int vindexm [8]={0, MAX1, MAX1*2, MAX1*3, MAX1*4, MAX1*5, MAX1*6, MAX1*7 };
    __m256i vindex = _mm256_load_si256((__m256i *) &vindexm[0]);
    __m256 vec1, vec2, vec3, vec4, vec5, vec6, vec7, vec8;

        for(i=0; i<MAX1;  i+=8){            
            for(j=0; j<MAX2;  j+=8){
                //loading from columns
                vec1 = _mm256_i32gather_ps (&a[i][j+0],vindex,4);
                vec2 = _mm256_i32gather_ps (&a[i][j+1],vindex,4);
                vec3 = _mm256_i32gather_ps (&a[i][j+2],vindex,4);
                vec4 = _mm256_i32gather_ps (&a[i][j+3],vindex,4);
                vec5 = _mm256_i32gather_ps (&a[i][j+4],vindex,4);
                vec6 = _mm256_i32gather_ps (&a[i][j+5],vindex,4);
                vec7 = _mm256_i32gather_ps (&a[i][j+6],vindex,4);
                vec8 = _mm256_i32gather_ps (&a[i][j+7],vindex,4);

                //storing to the rows
                _mm256_store_ps(&a_tra[j+0][i], vec1);
                _mm256_store_ps(&a_tra[j+1][i], vec2);
                _mm256_store_ps(&a_tra[j+2][i], vec3);
                _mm256_store_ps(&a_tra[j+3][i], vec4);
                _mm256_store_ps(&a_tra[j+4][i], vec5);
                _mm256_store_ps(&a_tra[j+5][i], vec6);
                _mm256_store_ps(&a_tra[j+6][i], vec7);
                _mm256_store_ps(&a_tra[j+7][i], vec8);  
            }
        }
    return 0;
}