神秘常量复出！用0x077CB531计算末尾0的个数 -- De Bruijn 序列 - IAmAProgrammer

大家或许还记得 Quake III 里面的一段有如天书般的代码，其中用到的神秘常量 0x5F3759DF 究竟是怎么一回事，着实让不少人伤透了脑筋。今天，我见到了一段同样诡异的代码。
下面这个位运算小技巧可以迅速给出一个数的二进制表达中末尾有多少个 0 。比如， 123 456 的二进制表达是 1 11100010 01000000 ，因此这个程序给出的结果就是 6 。

unsigned int v;  // find the number of trailing zeros in 32-bit v
int r;           // result goes here
static const int MultiplyDeBruijnBitPosition[32] =
{
  0, 1, 28, 2, 29, 14, 24, 3, 30, 22, 20, 15, 25, 17, 4, 8,
  31, 27, 13, 23, 21, 19, 16, 7, 26, 12, 18, 6, 11, 5, 10, 9
};
r = MultiplyDeBruijnBitPosition[((uint32_t)((v & -v) * 0x077CB531U)) >> 27];

熟悉位运算的朋友们可以认出， v & -v 的作用就是取出右起连续的 0 以及首次出现的 1 。

当 v = 123 456 时， v & -v 就等于 64 ，即二进制的 1000000 。

怪就怪在，这个 0x077CB531 是怎么回事？数组 MultiplyDeBruijnBitPosition 又是什么玩意儿呢？

这还得从 0x077CB531 本身的一个性质开始说起。把这个常数写成 32 位二进制，可以得到

00000111011111001011010100110001

这个 01 串有一个无比牛 B 的地方：

如果把它看作是循环的，它正好包含了全部 32 种可能的 5 位 01 串，既无重复，又无遗漏！

其实，这样的 01 串并不稀奇，因为构造这样的 01 串完全等价于寻找一个有向图中的 Euler 回路。

如下图，构造一个包含 16 个顶点的图，顶点分别命名为 0000, 0001, 0010, …, 1111 。

如果某个点的后 3 位，正好等于另一个点的前 3 位，就画一条从前者出发指向后者的箭头。

也就是说，只要两个顶点上的数满足 abcd 和 bcde 的关系（ a 、 b 、 c 、 d 、 e 可能代表相同的数字），

就从 abcd 出发，连一条到 bcde 的路，这条路就记作 abcde 。

注意，有些点之间是可以相互到达的（比如 1010 和 0101 ），有些点甚至有一条到达自己的路（比如 0000 ）。

构造一个字符串使其包含所有可能的 5 位 01 子串，其实就相当于沿着箭头在上图中游走的过程。

不妨假设字符串以 0000 开头。

如果下一个数字是 1 ，那么 00001 这个子串就被包含了，同时最新的 4 位数就变成了 0001 ；

但若下一个数字还是 0 ，那么 00000 就被包含了进来，最新的 4 个数仍然是 0000 。

从图上看，这无非是一个从 0000 点出发走了哪条路的问题：

你是选择了沿 00001 这条路走到了 0001 这个点，还是沿着 00000 这条路走回了 0000 这个点。

同理，每添加一个数字，就相当于沿着某条路走到了一个新的点，路上所写的 5 位数就是刚被考虑到的 5 位数。

我们的目的便是既无重复又无遗漏地遍历所有的路。

显然图中的每个顶点入度和出度都是 2 ，因此这个图一定存在 Euler 回路，

我们便能轻易构造出一个满足要求的 01 串了。这样的 01 串就叫做 De Bruijn 序列。

De Bruijn 序列在这里究竟有什么用呢？

它的用途其实很简单，就是为 32 种不同的情况提供了一个唯一索引。

比方说， 1000000 后面有 6 个 0 ，将 1000000 乘以 0x077CB531 ，就得到

   00000111011111001011010100110001
-> 11011111001011010100110001000000

相当于把 De Bruijn 序列左移了 6 位。

再把这个数右移 27 位，就相当于提取出了这个数的头 5 位：

  11011111001011010100110001000000
->00000000000000000000000000011011

由于 De Bruijn 序列的性质，因此当输入数字的末尾 0 个数不同时，最后得到的这个 5 位数也不同。

而数组 MultiplyDeBruijnBitPosition 则相当于一个字典的功能。

11011 转回十进制是 27 ，于是我们查一查 MultiplyDeBruijnBitPosition[27] ，程序即返回 6 。

注意到当输入数字的末尾 0 个数超过 27 个时，程序也是正确的，因为左移时低位正好是用 0 填充的。

这段神一般的代码取自 http://graphics.stanford.edu/~seander/bithacks.html ，欢迎大家前去围观。

http://blog.sina.com.cn/s/blog_000199c301018qai.html

二进制中的神秘数字，位运算

我在coding的时候，比较喜欢bitmap，而bitmap上的操作离不开位运算。

经常会遇到的一个问题就是，给定一个数字，向上取整到2的幂次，然后得到位数。

向上取整到2的幂次，非常简单：

v--;
v|=v>>1;
v|=v>>2;
v|=v>>4;
v|=v>>8;
v|=v>>16;
v++;

但是要得到2^n，这个n的值，就不是那么容易了。也就是求LOG2(v)

之前的时候，实现LOG2，我使用了一个65536个数的数组，表65536内的所有数字的log2都存起来。

这种方法浪费了内存，而且也不是非常高效，如果求32位数字的log2，需要做2次运算。64位的话，需要4次运算

昨天突然发现了一段代码：

unsigned int v;  // 整数
int r;           // 运算结果
static const int MultiplyDeBruijnBitPosition[32] =
{
  0, 1, 28, 2, 29, 14, 24, 3, 30, 22, 20, 15, 25, 17, 4, 8,
  31, 27, 13, 23, 21, 19, 16, 7, 26, 12, 18, 6, 11, 5, 10, 9
};
r = MultiplyDeBruijnBitPosition[((uint32_t)((v & -v) * 0x077CB531U)) >> 27];

只需要构造一个32数字的哈希表，一次运算就可以得到结果。
仔细研究了一下这个算法的奥妙，看似神秘的代码，其实并不难理解。

要求一个数字的log2，其实一共就只有32种情况，

假如我们能构造一个数字，使得这个数字的二进制表示中，任意取5个数，都能保证不重复，就可以了。

因为不重复，那么这个数字每左移一位，就会产生一个32内的数字，移动32次就遍历完了32内的所有数字。

然后，再用该数字的对应序列做一次映射，就得到了真正的结果。

使用这个算法，如果求64位数字的log2，就需要做2次运算，还不是很高效。

所以，自己就写程序来寻找64位的那个神秘数字，代码很简单，其实就是寻找有向图中的欧拉回路问题。

这个问题的解不是唯一的，有很多数字都满足这样的条件，只要找到一个就可以。

64位的如下：神秘数字：0x3f79d71b4cb0a89UL
64个数字的哈希表为：
const char g_math_64[] ={0,1,48,2,57,49,28,3,61,58,50,42,38,29,17,4,62,55,59,36,53,51,
43,22,45,39,33,30,24,18,12,5,63,47,56,27,60,41,37,16,54,35,52,21,44,32,23,11,46,26,40,
15,34,20,31,10,25,14,19,9,13,8,7,6};
求log2的算法：
#define LOG2_64(v) g_math_64[(uint64_t)((v) * 0x3f79d71b4cb0a89UL) >> 58]

32位，我自己找的一个数字：0x7dcd629
const char g_math_32[] ={0,1,23,2,29,24,14,3,30,27,25,18,20,15,10,4,31,22,28,13,26,17,19,9,21,12,16,8,11,7,6,5};
#define LOG2(v) g_math_32[(uint32_t)((v) * 0x7dcd629) >> 27]

128位的为：0x01fdf3d78edd3970d9ab464c582a5091
const char g_math_128[] ={0,1,101,2,116,102,60,3,124,117,103,94,82,61,33,4,125,121,118,87,111,104,95,53,90,83,
69,62,48,34,20,5,126,114,122,80,119,109,88,46,112,107,105,73,96,75,54,26,98,91,84,66,
77,70,63,39,56,49,42,35,28,21,14,6,127,100,115,59,123,93,81,32,120,86,110,52,89,68,47,
19,113,79,108,45,106,72,74,25,97,65,76,38,55,41,27,13,99,58,92,31,85,51,67,18,78,44,71,
24,64,37,40,12,57,30,50,17,43,23,36,11,29,16,22,10,15,9,8,7};
求log2的算法：
#define LOG2_128(v) g_math_128[(__uint128_t)((v) * 0x01fdf3d78edd3970d9ab464c582a5091) >> 121]

128位的这一个，gcc似乎并没有提供如何构造一个128位的数字常量。

如果直接写：0x01fdf3d78edd3970d9ab464c582a5091 并不能转成128位数字。所以这个算法不能直接用。
只是先把这个发出来，记录一下吧。
以后如果编译器完全支持128位数字了，就可以拿来用了。