最快的4x4字节矩阵转置方法

Question

我有一块 4x4 字节的区块，希望使用通用硬件进行转置。换句话说，对于字节 A-P，我正在寻找最有效（指指令数量）的方法来转置。

A B C D
E F G H
I J K L
M N O P

to

A E I M
B F J N
C G K O
D H L P

我们假设我在内存中有指向 A、E、I 和 M 的有效指针（这样从 A 中读取 32 位将得到包含字节 ABCD 的整数）。

这不是该问题的重复，因为存在大小和数据类型方面的限制。我的矩阵的每一行都可以适合于一个 32 位整数，并且我正在寻找可以使用通用硬件快速执行转置的答案，类似于 SSE 宏 _MM_TRANSPOSE4_PS 的实现。

- Mokosha

---

了解为什么需要这样做将非常有用 - 我所知道的大多数矩阵处理库很少实际上对矩阵进行转置，它们只是标记它为转置，并将 (行，列) 访问为 (列，行)。不需要移动任何内存。 - Roger Rowland

当转置一个16byte矩阵时，你可能会发现任何转置标志在重复访问时的成本都比转置本身更高。 - Kuba hasn't forgotten Monica

@KubaOber 是的，确实如此，这就是为什么我很想了解 OP 提出这个问题背后的原因 - 对于一个大矩阵，转置可以通过引用局部性来提高性能，但正如你所说，一个非常小的矩阵经常通过交换索引访问可能会在其他方面遭受开销。我只是想要更多的背景信息... - Roger Rowland

@RogerRowland，一个有用的字节转置的例子是将RGBARGBARGBARGBA中的四个像素转置为RRRRGGGGBBBBAAAA。 - Z boson

@Zboson 很酷，你因为那个回答又获得了一个赞 :-) - Roger Rowland

Answer 1

你想要可移植性和效率，但两者不可兼得。你说你想用最少的指令来完成这个任务。使用x86指令集中的SSE3中的pshufb指令（见下文）可以仅使用一个指令完成。

也许ARM Neon有类似的东西。如果你想要效率（并且确定你需要它），那就学习硬件。

SSE在字节级别上等价于_MM_TRANSPOSE4_PS的方法是使用具有掩码的_mm_shuffle_epi8（pshufb的内联函数）。在主循环外定义掩码。

//use -msse3 with GCC or /arch:SSE2 with MSVC
#include <stdio.h>
#include <tmmintrin.h> //SSSE3
int main() {
    char x[] = {0,1,2,3, 4,5,6,7, 8,9,10,11, 12,13,15,16};
    __m128i mask = _mm_setr_epi8(0x0,0x04,0x08,0x0c, 0x01,0x05,0x09,0x0d, 0x02,0x06,0x0a,0x0e, 0x03,0x07,0x0b,0x0f);

    __m128i v = _mm_loadu_si128((__m128i*)x);
    v = _mm_shuffle_epi8(v,mask);
    _mm_storeu_si128((__m128i*)x,v);
    for(int i=0; i<16; i++) printf("%d ", x[i]); printf("\n");
    //output: 0 4 8 12 1 5 9 13 2 6 10 15 3 7 11 16   
}

- Z boson

---

在某些情况下，使用 _mm_lddqu_si128（） 可以提高性能，但不会使性能变差。如果 x[] 对齐到 16 字节，则使用 _mm_load_si128 和 _mm_store_si128 可能会更快。 - St0fF

@St0fF，谢谢！我之前从未听说过 _mm_lddqu_si128()。很有趣。 - Z boson

Answer 2

让我重新阐述您的问题：您正在寻求一种仅限于C或C ++的可移植解决方案。那么：

void transpose(uint32_t const in[4], uint32_t out[4]) {
  // A B C D    A E I M
  // E F G H    B F J N
  // I J K L    C G K O
  // M N O P    D H L P

  out[0] = in[0] & 0xFF000000U; // A . . .
  out[1] = in[1] & 0x00FF0000U; // . F . .
  out[2] = in[2] & 0x0000FF00U; // . . K .
  out[3] = in[3] & 0x000000FFU; // . . . P

  out[1] |= (in[0] <<  8) & 0xFF000000U; // B F . .
  out[2] |= (in[0] << 16) & 0xFF000000U; // C . K .
  out[3] |= (in[0] << 24);               // D . . P

  out[0] |= (in[1] >>  8) & 0x00FF0000U; // A E . .
  out[2] |= (in[1] <<  8) & 0x00FF0000U; // C G K .
  out[3] |= (in[1] << 16) & 0x00FF0000U; // D H . P

  out[0] |= (in[2] >> 16) & 0x0000FF00U; // A E I .
  out[1] |= (in[2] >>  8) & 0x0000FF00U; // B F J .
  out[3] |= (in[2] <<  8) & 0x0000FF00U; // D H L P

  out[0] |= (in[3] >> 24);               // A E I M
  out[1] |= (in[3] >>  8) & 0x000000FFU; // B F J N
  out[2] |= (in[3] <<  8) & 0x000000FFU; // C G K O
}

我不认为还有其他答案，因为这样的话你就得依赖于特定的编译器以特定的方式进行编译等。

当然，如果那些操作本身能够以某种方式简化，那么就会有所帮助。所以这是进一步追求的唯一途径。目前还没有什么突出的地方，但对我来说，已经是漫长的一天了。

到目前为止，成本是12个移位，12个或运算，16个与运算。如果编译器和平台都很好，可以在9个32位寄存器中完成。

如果编译器非常糟糕，或者平台没有桶形移位器，则一些类型转换可以帮助强调移位和掩码只是字节提取：

void transpose(uint8_t const in[16], uint8_t out[16]) {
  // A B C D    A E I M
  // E F G H    B F J N
  // I J K L    C G K O
  // M N O P    D H L P

  out[0]  = in[0];  // A . . .
  out[1]  = in[4];  // A E . .
  out[2]  = in[8];  // A E I .
  out[3]  = in[12]; // A E I M
  out[4]  = in[1];  // B . . .
  out[5]  = in[5];  // B F . .
  out[6]  = in[9];  // B F J .
  out[7]  = in[13]; // B F J N
  out[8]  = in[2];  // C . . .
  out[9]  = in[6];  // C G . .
  out[10] = in[10]; // C G K .
  out[11] = in[14]; // C G K O
  out[12] = in[3];  // D . . .
  out[13] = in[7];  // D H . .
  out[14] = in[11]; // D H L .
  out[15] = in[15]; // D H L P
}

如果您真的想要在原地打乱它，那么以下方法可以实现。

void transpose(uint8_t m[16]) {
  std::swap(m[1], m[4]);
  std::swap(m[2], m[8]);
  std::swap(m[3], m[12]);
  std::swap(m[6], m[9]);
  std::swap(m[7], m[13]);
  std::swap(m[11], m[14]);
}

字节定向版本在现代平台上可能会生成更差的代码。只有基准测试才能说明问题。

- Kuba hasn't forgotten Monica

---

+1 教会我一个新词。我以前从没听说过“移位器”（太沉迷于高级语言了）。 - Roger Rowland

这可以通过一个 SSSE3 指令完成：pshufb（忽略加载和存储）。请参见我的答案。 - Z boson

@Zboson 这个问题没有提供平台信息，所以对我来说，任何明确绑定到一个平台的内容都不在回答的范围内。了解这些内容仍然很有用，但我特别不想调用任何特定于平台的内部函数。 - Kuba hasn't forgotten Monica

@KubaOber，我明白。虽然我不确定我的回答是否超出了OP问题的范围。OP正在询问两件事情，他不能同时拥有它们。你回答了其中之一（可移植性），而我想表明另一个（效率）无法在不违反第一个条件的情况下实现（12个移位，12个OR和16个AND比一个洗牌要糟糕得多）。 - Z boson

你第一个代码示例的最后两行是不正确的。你应该分别进行右移16位和8位。 - Chris_F

Answer 3

不确定速度，但这些还可以。

template<typename T, std::size_t Size>
void Transpose(T (&Data)[Size][Size])
{
    for (int I = 0; I < Size; ++I)
    {
        for (int J = 0; J < I; ++J)
        {
            std::swap(Data[I][J], Data[J][I]);
        }
    }
}

template<typename T, std::size_t Size>
void Transpose(T (&Data)[Size * Size])
{
    for (int I = 0; I < Size; ++I)
    {
        for (int J = 0; J < I; ++J)
        {
            std::swap(Data[I * Size + J], Data[J * Size + I]);
        }
    }
}

- Brandon

Answer 4

如果您接受的话，64位机器上可以提供高效的解决方案。首先将32位整数常量分别向左移动（0,）1、2和3个字节 [3次移位]。然后屏蔽掉不需要的位并执行逻辑或运算 [12次与常量相与，12次逻辑或运算]。最后，向右移回32位 [3次移位] 并读取32位。

ABCD
EFGH
IJKL
MNOP

ABCD
 EFGH
  IJKL
   MNOP

A---
 E---
  I---
   MNOP
=======
AEIMNOP
AEIM

AB--
 -F--
  -J--
   -NOP
=======
ABFJNOP
BFJN

ABC-
 --G-
  --K-
   --OP
=======
ABCGKOP
CGKO

ABCD
 ---H
  ---L
   ---P
=======
ABCDHLP
DHLP

- user1196549

Answer 5

我之前为SSE 这里的同一个问题提供了一个答案。

唯一需要添加的是向量化的加载/存储操作。

这个答案类似于Z boson对这个问题的回答。可以在那里看到加载/存储的例子。这个答案不同之处在于除了SSE3实现外，还有一个SSE2实现，它保证可以运行在任何x64处理器上。

值得注意的是，这两个解决方案都假定整个矩阵在内存中是按行主序存储的，但OP的问题说明每行可能有自己的指针，这意味着数组可能是分散的。

- Apriori