Rutkowska的位反转算法

6
我发现了一篇关于适用于原地FFT的比特反转算法的非常有趣的论文:“一个简单的比特反转置换算法”,作者是Urszula Rutkowska,发表于1990年(doi.org/10.1016/0165-1684(91)90008-7)。
然而,她的G1算法似乎不起作用,因为对于N1 = L << 1swap(a + 1, a + N1);,第一次迭代的结果超出了边界。我假设L表示输入向量的长度。
请问是否有任何关于该论文的勘误或如何修复算法的信息?
该论文的伪代码:
G1(L)
{int     i,j,L1
         N1,N2,a,b;
unsigned k;
j=0;    L1=L-1;
N1=L<<1;N2=N1+1;
for(i=0;i<L1;i++)
{if(i<j)
    { a=i<<1;
      b=j<<1;
      swap(a,b);
      swap(a+N2,b+N2);
      swap(a+1,b+N1);
      swap(b+1,a+N1);
    }
 else
    if(i==j)
    { a=i<<1;
      swap(a+1,a+N1);
    }
 k=L>>1;
 while(k<=j){ j=j-k;
              k=k>>1;
            }
 j+=k;
 }
 i<<=1;
 swap(i+1,i+N1);
}

论文截图:

在此输入图片描述


注:本段内容无需翻译,已是中文。

1
让我建议使用Elster的线性时间算法:https://www.idi.ntnu.no/~elster/pubs/elster-bit-rev-1989.pdf(参见算法3; n = 2 ^ t)。[免责声明:Elster是我的导师。] - Aasmund Eldhuset
1
这个内容是从 DSP 迁移过来的,但我看不出它为什么适合在这里而不是那里。如果迁移此内容的 moderator @jojek 能够解释原因那就太好了。首先,在这里放置代码图像是完全不被允许的。 - Cris Luengo
1
我在Meta DSP上发布了一个关于此迁移的问题:https://dsp.meta.stackexchange.com/questions/1602/why-was-this-question-migrated-to-so - Cris Luengo
可能移位运算符在错误的维度中。对数组大小的两倍进行操作没有太多意义。相反,更有可能将数组分成两个切片。 - Marcel
1个回答

1

说实话,它相当混乱。我不得不阅读论文以获取想法(运行Gold算法(G)的L/4,然后推导出L的交换),然后将代码按正确形式进行调整。这是我的最终结果。

#include <assert.h>
#include <stdbool.h>
#include <stdio.h>
#include <stdlib.h>

static bool is_power_of_two(int L) { return L > 0 && (L & (L - 1)) == 0; }

static void swap(int i, int j) { printf("swap %d,%d\n", i, j); }

static void G(int L) {
  assert(is_power_of_two(L));
  int j = 0;
  for (int i = 0; i < L - 1; i++) {
    if (i < j) {
      swap(i, j);
    }
    int k = L >> 1;
    while (j >= k) {
      j -= k;
      k >>= 1;
    }
    j += k;
  }
}

static void G1(int L) {
  assert(is_power_of_two(L));
  if (L < 4) {
    return;
  }
  int j = 0;
  int N1 = L >> 1;
  int N2 = N1 + 1;
  int L2 = L >> 2;
  for (int i = 0; i < L2 - 1; i++) {
    if (i < j) {
      int a = i << 1;
      int b = j << 1;
      swap(a, b);
      swap(a + N2, b + N2);
      swap(a + 1, b + N1);
      swap(b + 1, a + N1);
    } else if (i == j) {
      int a = i << 1;
      swap(a + 1, a + N1);
    }
    int k = L2 >> 1;
    while (j >= k) {
      j -= k;
      k >>= 1;
    }
    j += k;
  }
  int a = (L2 - 1) << 1;
  swap(a + 1, a + N1);
}

int main(int argc, char *argv[]) {
  assert(1 < argc);
  int L = atoi(argv[1]);
  G(L);
  putchar('\n');
  G1(L);
}

谢谢,我可以确认现在它可以工作了!PS:Karp在他的《单处理器上的位反转》一文中对几种位反转算法进行了基准测试,了解他是如何修复问题的将会很有趣。 - minmax
我的比特位逆序研究的官方引用为“单处理器上的比特逆序”,出自SIAM评论杂志,卷38,#1,1-26页,1996年3月。(我曾尝试在注释中加入此引用,但系统不允许。) - user1649736
1
在私人通信中,Alan Karp很友好地让我看到了他的实现,唯一的区别是这里的倒数第二行读作a = (L2 - 1) << 1,而在他的代码中类似于a = i << 1。当然,这会导致相同的输出结果。 - minmax
@minmax 是的,我希望循环变量在循环中具有作用域,即使这意味着要进行冗余算术运算。 - David Eisenstat

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接