使用 Intel 内在函数的位反向重新排序优化

Posted 2023-02-16

【中文标题】使用 Intel 内在函数的位反向重新排序优化

【英文标题】：Bit reverse reorder optimization using Intel intrinsics

【发布时间】：2016-11-10 14:28:56

【问题描述】：

前段时间，我使用 SSE 内在函数优化了一个 radix2 函数，我几乎接近 FFTW 性能，所以，我的下一步是优化我在原始代码中找到的位反向重新排序函数，老实说，我想优化它。原来的代码是这样的：

void bit_reverse_reorder(float *x,float *y, int N)

   int bits=0;
   int i, j, k;
   float tempr, tempi;
   //MAXPOW = 12

   for (i=0; i<MAXPOW; i++)
      if (pow_2[i]==N) bits=i;

   for (i=0; i<N; i++)
   
      j=0;
      for (k=0; k<bits; k++)
         if (i&pow_2[k]) j+=pow_2[bits-k-1];
      if (j>i)
      
        tempr=x[i];
        tempi=y[i];
        x[i]=x[j];
        y[i]=y[j];
        x[j]=tempr;
        y[j]=tempi;
      
   


int main()

   radix2(re,im,N);

   bit_reverse_reorder(re,im,N);

PS：pow_2[] 是一个预先计算好的数组，包含 2 的幂 (1,2,4,8,16,32,...)，N 是元素的数量 = 4096 , *x 和 *y 分别表示输入数据的每个元素的实部和虚部。

radix2 生成未排序的结果，因此所述函数重新排序结果。

首先，我并没有完全理解这个位反转是如何工作的！所以，我认为如果有人给我关于这个功能如何工作的提示会很好。

其次，我打算使用 SSE 内在函数来增强此函数的性能，那么是否有任何 2 条指令向量可用于交换循环？

【问题讨论】：

我可能是错的，但我认为 bit_reverse_reorder 函数正在执行 FFT 的 butterfly reordering 步骤。如果你对 FFT 的那一步不太了解，那么很难优化这个函数。

相关：In-place bit-reversed shuffle on an array

位反转重新排序很可能是您执行配置文件中相当微不足道的部分 - 您是否实际测量/分析它以查看是否值得进一步优化？请记住，FFT 本身是 O(n log n)，而您的位反转阶段只是 O(n) - 除非 n 非常小，否则它不应该成为性能瓶颈。

即使在 AES 解压缩等经常使用此类东西的情况下，这一步也可能不是您最大的打击。如果是的话，看看那里的各种 AES-128 解压缩器 (openSSL)，你将有一个如何在 Intel 中执行此操作的示例。

看看my answer 我上面链接的问题。它不使用 SIMD 指令，所以我决定将我的答案添加到更一般的问题中。无论如何，SIMD 可能对您没有帮助，除非如果您交错 x 和 y 数组，可能会一次性交换元素。如果MAXPOW 是静态的，或者如果您多次执行 FFT，最好的方法可能是预先计算要交换的索引。

【参考方案1】：

感谢@nwellnhof 的评论，我对交换功能进行了如下修改：

void bit_reverse_reorder(float *x,float *y, int N)

   unsigned i,j;
   for (i = 0, j = 0; i < N; i++) 
      if (i < j) 
      
         float tmpx = x[i];
         x[i] = x[j];
         x[j] = tmpx;

         float tmpy = y[i];
         y[i] = y[j];
         y[j] = tmpy;
      
      unsigned bit = ~i & (i + 1);

      unsigned rev = (N / 2) / bit;

      j ^= (N - 1) & ~(rev - 1);
   

现在我在函数内部的 for 循环中获得了 54 900 个循环的性能，这也很好:)

【讨论】：

bit 是 2 的幂，对吧？如果编译器没有弄清楚，它将使用实际的 DIV 指令，而不仅仅是右移。有关使用 TZCNT 或 GNU C __builtin_ctz 实现有效的 log2 函数的一些想法，请参见 ***.com/questions/40431599/…，该函数假定输入 2 的幂，为您提供移位计数。 （这个问题的主要部分是关于四舍五入的，但你只需要unsigned rev = N >> (1+lg2(bit));

好的，我会试一试...这个问题将用于 intel 和 ARM 架构。所以，我不知道__builtin_ctz 是否适用于所有架构:)

它可能无法在 ARM 上编译成任何美妙的东西，但它肯定受到支持，因为它不是特定于平台的 GNU 扩展。

我很好奇并检查了 (godbolt.org/g/ewhq2i)：ARM 有一个计数前导零指令和一个位反转指令。如果您知道（但编译器不知道）您有 2 的幂，请选择 32 - count leading zeros 而不是尾随零，这样编译器就不必对寄存器进行位反转。

【参考方案2】：

根据您的建议，我进行了一些修改以增强功能性能。

首先，我通过移位指令替换了 power_2[] 调用。

然后，我制作了一个交换函数，它使用 add/sub 操作进行交换，而不使用第三个变量。

void swap(float* a, float* b)

 *a = *a+*b;
 *b = *a-*b;
 *a = *a-*b;


void bit_reverse_reorder(float *x,float *y, int N)

   int bits=0;
   int i, j, k;
   unsigned pow_bits=0;

   for (i=0; i<MAXPOW; i++)
      if (1 << i==N) bits=i-1;
      for (i=1; i<N; i++)
      
          j=0;
          unsigned pow_2k=0,pow_2kj;

          for (k=0; k<bits; k++)
          
             pow_2k = 1<<k;
             pow_2kj = 1<<(bits-k);
             if (i&pow_2k) j+=pow_2kj;
          
          if (j>i)
          
             swap(&x[i],&x[j]);
             swap(&y[i],&y[j]);
          
       

周期数从近 500 000 个周期减少到大约 180 000 个周期。

【讨论】：

你绝对不想那样交换，尤其是。如果你在没有-ffast-math 的情况下编译。以这种方式编写它会强制编译器实际生成具有正确舍入的临时变量，而不仅仅是交换。不要试图避开临时工；它们很便宜。编译器知道如何使用寄存器。

查看godbolt上的asm输出，例如：godbolt.org/g/NR4Z7K。 swap() 编译成一堆数学运算，而不仅仅是两次加载和两次存储。

另外它还存在别名问题。