详解二进制数中1的个数

作者: | 更新日期:

之前记录了一些得到二进制中1的个数的算法, 但是没有讲解为什么, 今天就详细的讲解一下啊.

本文首发于公众号:天空的代码世界,微信号:tiankonguse

mit-hackmem-count

前言

之前写了一篇文章《二进制数中1的个数》.
然后有人说看不懂, 然后我自己再看看, 原来对每个代码我只是随意的说了几句, 并没有详细的讲解, 不是 acmer 的话还真是看不懂的.
于是今天再详细的记录一下.

突然想起一件事:小舟学长曾给我们讲过二进制数中1的个数这个专题, 而我只是轻视的说了一句:统计多少个1, 遍历一遍不就行了, 就32位. 现在想想, 只能说当时自己太年轻了.

更新细节:添加了64位的代码样例,代码解释请参考32位的。

更新细节: 纠正了64位代码考虑不全的bug, 感谢热心网友的提醒。

基础知识

在计算机中, 数字是以二进制储存的.
我们平常生活中用的是十进制.
假设对于十进制 10030, 我们想得到第二位的数字是多少的话需要这样:

(10030 / 10^1 ) % 10 

也就是需要先把需要查询的数字移到个位, 然后在取模即可.

而对于二进制 10010, 想得到第二位的数字就有很多方法了.

(10010 >> 1 ) & 1
(10010 & 10) >> 1

其中上面的第一个方法和十进制的思想相同, 第二个是想把其他位置0, 然后把第二位移到个位.

然后我们就可以判断这位是不是我们想要的数字了。

if((10010 >> 1 ) & 1 == 1){
    //do some thing
}

然后对于二进制,每一位只有0或1,所以当那一位是1时,和1判断是否相等得到true, 而数字1依然是true
于是我们可以简写为下面的形式

if((10010 >> 1) & 1){
    //do some thing
}

而现在我们要做的是统计一个32位整数的二进制有多少个1.

最简单的算法

实际上最简单的方法是下面的代码, 直接一位一位的判断即可.

unsigned __countbits(unsigned x){
    unsigned n=0;
    for(int i=0;i<32;i++){
        n += x&1;
        x>>=1;
    }
    return n;
}

相对暴力的算法

大家可以看到, 当 x 为0 的时候, 我们并没有判断还有没有1, 而是继续循环, 直到32位遍历结束.
什么意思呢, 假设传入的是 1, 我们进入循环之后, x就变成0了, 然后再循环是不会更新答案n的, 所以我们可以结束循环了.

于是下面的代码就出来了.
其中 使用了逗号表达式, ture 与 false 只与最后一个表达式有关系.

unsigned __countbits(unsigned x){
    unsigned n=0;
    while(n +=(x&1) , x >>= 1);
    return n;
}

与 1 的个数有关的算法

大家先来看一个实验.
假设我们有 二进制 101100, 这个二进制减一后是 101100 - 1 = 101011, 这两个数字的后三位不同, 且数字相反.
然后这两个数字进行与操作后得到 101100 & 101011 = 101000, 我们可以发现少了一个 1.
然后再来一次, 原始数字: 101000, 减一: 101000 - 1 = 100111, 与操作 : 101000 & 100111 = 100000,
最后再来一次, 原始数字: 100000, 减一: 100000 - 1 = 011111, 与操作 : 100000 & 011111 = 0.

我们发现1的个数与我们进行的那个运算次数相等.
那这是不是运气呢?
假设这个数字的二进制有1, 则我们找到最后一个1, 1后面可能若干个0. 这个数字减1时, 前面若干0那些位都不够减的, 需要向高位借1, 直到最后一个1, 借到了1, 于是1变成0.
后面的位向高位借到1时是2, 再借给低位一个, 又变成了1.
最后一位借到1是2, 减去我们需要减的1也变成1了.
原来这个不是运气, 也不是偶然, 所有数字都遵循这个规律的.
于是我们可以使用这个规律来计算1的个数了.

复杂度就是1的个数了.
需要注意的是对于0, 是没有1的, 我们需要特判, 不然会出错的, 我就踩过这个.

unsigned _countbits(unsigned x){
    unsigned n=0;
    while(x && (++n , x&=x-1));
    return n;
}

log 算法 分支法

用过二分查找的人应该很喜欢二分吧, 复杂度那是 log 级的, 用着太爽了.
而我们32位整数, log(32) = 5, 也就是我们只需要5步就可以得到我们想要的答案了.

但是怎么二分了.
首先将32位数看成相互独立的32个 0-1 数字.
然后两两相加, 这样我们就有16个数字了.
然后再两两相加, 就变成8个数字了, 接着再两两相加变成4个, 然后两个, 最后一个数字.
由于是32位数字的和, 所以答案刚好就是1的个数啦.

现在问题是怎么数字两两相加呢?
32位太多了, 我们看看4位的数字吧, 假设是 1101.
要想相邻两位相加, 那就需要右移一下, 这样就能对齐了.
于是我们得到下面的数据.

11 01
01 10

然后我们发现这个时候还不能相加, 因为对于两位数字来说, 高位还有垃圾数据, 我们需要清理掉.
需要清理的数据刚好都是高位, 那就把那些高位的位置置为0吧.

11 01 & 01 01 = 01 01
01 10 & 01 01 = 01 00

然后这个时候就可以相加了.

01 01 + 01 00 = 10 01

这两个数字接下来怎么相加呢?
应该是右移两位了.

1001
0010

发现高位还有垃圾数据, 需要清理掉.

1001 & 0011 = 0001
0010 & 0011 = 0010

然后相加

0001 + 0010 = 0011

0011 的十进制就是3, 这样果然可以得到二进制1的个数.

而且可以总结为

右移 1位. 与数字: 0101 => 0x5
右移 2位. 与数字: 0011 => 0x3

那对于32位整数需要操作哪五部呢?

右移 1位. 与数字: 01010101010101010101010101010101 => 0x55555555
右移 2位. 与数字: 00110011001100110011001100110011 => 0x33333333
右移 4位. 与数字: 00001111000011110000111100001111 => 0x0F0F0F0F
右移 8位. 与数字: 00000000111111110000000011111111 => 0x00FF00FF
右移16位. 与数字: 00000000000000001111111111111111 => 0x0000FFFF

由上面的方法, 我们可以快速的写出32位对应的代码了.

uint countbits(uint x) {
    uint mask[]= {0x55555555,0x33333333, 0x0F0F0F0F,
                  0x00FF00FF, 0x0000FFFF
                 };
    for(uint i=0,j=1; i<5; i++,j<<=1) {
        x=(x & mask[i]) + ((x>>j) & mask[i]);
    }
    return x;
}

打表法

分治法虽然是 log 级的, 但是有时候还是感觉复杂度太高, 想追求O(1)的做法.
这个时候就需要拿出空间来换取时间了.
什么意思呢?
我们想求x的二进制1的个数, 如果我们已经事先求出来了的话, 那是不是只需要取得那个值不就行了.
但是32位整数太多了, 我们存不下.
那16位就很小了, 就几万, 可以试试.

const int OneNumMax = 1<<16;
int oneNum[OneNumMax];
int init(){
    for(int i = 1; i < OneNumMax; ++i) {
        oneNum[i] = oneNum[i >> 1] + (i & 1);
    }
    return 0;
}

那么对于32位整数, 不就是两个16位整数吗?

对于高16位数字, 我们直接右移16位就得到了.
对于低16位, 我们可以进行与运算把高位置0得到.
于是代码就可以敲出来了.

unsigned countbits(unsigned x) {
    return oneNum[x >> 16] + oneNum[x & ((1 << 16) - 1)];
}

MIT HACKMEM

这个算法是 MIT 实验室的备忘录中的一个算法.
上面我们看到了二分的好处, 按我们如果三分会怎么样呢?

三分的第一步是得到每三位的数字和.

于是有人这样想了

master = 011111111111;
x = (x & master) + ((x >> 1) & master) + ((n >> 2) & master);

那有没有其他方法呢?

对于三位, 我们的目的是下面的映射

000 => 0 
001 => 1
010 => 1
011 => 2
100 => 1
101 => 2
110 => 2
111 => 3

然后我们发现一个奇妙的等价公式

321         321   32   3
000 => 0 <= 000 - 00 - 0
001 => 1 <= 001 - 00 - 0
010 => 1 <= 010 - 00 - 1
011 => 2 <= 011 - 01 - 0
100 => 1 <= 100 - 10 - 1
101 => 2 <= 101 - 10 - 1
110 => 2 <= 110 - 11 - 1
111 => 3 <= 111 - 11 - 1

竟然这个巧, 于是下面的式子也是成立的.

unsigned tmp = n - ((n >> 1) & 033333333333) - ((n >> 2) & 011111111111);

现在我们得到了每三位的数字和, 然后你会不会想我是该计算每六位的数字和了?
恭喜你, 打对了.
32位共有10个三位和一个两位.

每个三位最高值为 11, 相加是 110, 所以低3位就可以表示了, 高三位是没有用的, 因此可以计算结束后在置0.

tmp = (tmp + (tmp >> 3) ) & 030707070707

这个时候, 我们得到的是每六位的数字和了.
共有5个六位和一个两位.
其中5个六位的数字都在低三位, 高三位确定为0.
数字可以表示为

tmp = a*64^5 + b*64^4 + c*64^3 + d*64^2 + e*64^1 + f

我们的目标是求 a+b+c+d+e+f 的答案. 然后又有人发现一个公式.
前面的数字都是64的幂数倍, 对于这个规律, 可以使用模 64-1 来抵消的.

  64^k % 63 
= (64 % 63)^k % 63 
= 1

所以

a * 64^k % 63 = a

也就是说

  tmp % 63
= (a*64^5 + b*64^4 + c*64^3 + d*64^2 + e*64^1 + f) % 63
= ((a*64^5 % 63) + (b*64^4 % 63) + (c*64^3 % 63) + (d*64^2 % 63) + (e*64^1 % 63) + (f % 63)) % 63
= (a % 63 + b%63 + c%63 + d%63 + e%63 + f%63) % 65
= (a + b + c + d + e + f) % 63

而我们32位的答案刚好最多是32, 所以可以模上63, 因此下面的公式就诞生了.

  tmp % 63
= a + b + c + d + e + f

于是我们的二进制中1的个数也可以用这个方法做了.

ans = tmp % 63

代码综合起来就是

int bitCount ( unsigned n ) {
    unsigned tmp = n - ((n >> 1) & 033333333333) - ((n >> 2) & 011111111111);
    return ( (tmp + (tmp >> 3) ) & 030707070707) % 63;
}

注:最后一个算法时我在看 sphinx 源代码的时候看到的.

64位中二进制1的个数

有位同学读了我的记录,然后很热情的反馈了很多问题。
最后他说他想要的其实是一个计算64位数字二进制1的个数,于是我便附加了一个小节。

windows 一个坑

在开始之前我们定义一些类型吧。
按照 acm 的习惯,我喜欢把 long long 定义为 LL 。

为什么呢? 这个你就该去查查 windows 的手册了,虽然现在 windows 也慢慢的开源了,但是它在程序员心中的丑陋印象还是很难改变的。

简单的说就是在 windows 下, 你使用 long long 的话, 需要使用 __int64 代替,否则肯定会出错。
至于为什么,那是因为在 windows 下根本就没有 long long 这个关键字。
大家可以看看windows下的C++ Keywords
当然,这个是 Microsoft C++ 编译器下的关键字,有人会说我安装的是 GNU 的 c++ 编译器。
虽然 windows 在官网文档上一直强调 long long 与 __int64 等价,但是事实上并不等价的。
好了,不说废话了,我的所有程序都会加上下面这个宏的。

#ifdef __int64
    typedef __int64 LL;
    typedef  unsigned __int64  ULL;
#else
    typedef long long LL;
    typedef unsigned long long ULL;
#endif

暴力求64位1的个数

暴力程序虽然暴力,但是有时候是最可靠的程序。
比如用来验证我们的其他程序是否正确。

int _countbits(ULL x) {
    int n=0;
    while(n +=(x&1) , x >>= 1);
    return n;
}

与1的个数有关复杂度

这个程序完全可以代替上面的暴力程序,而且只是还快点。

int _countbits(ULL x) {
    int n=0;
    while(x && (++n , x&=x-1));
    return n;
}

打表法

这个打表的复杂度看起来就稍微多了一些。

const int OneNumMax = 1<<16;
ULL oneNum[OneNumMax];

int init() {
    for(int i = 1; i < OneNumMax; ++i) {
        oneNum[i] = oneNum[i >> 1] + (i & 1);
    }
    return 0;
}

int countbits(ULL x) {
    static int a = init();
    static  ULL mask = (1ULL << 16) - 1;
    ULL ans = 0;
    for(int i=0;i<4;i++,x>>=16){
        ans += oneNum[x & mask];
    }
    return ans;
}

错误的MIT HACKMEM

关于测试代码可以参考这里

ULL mitHeakArray[3];
bool okMitHeak = false;
int initMitHeakArray() {
    if(okMitHeak) {
        return true;
    }
    ULL a = 011ULL, b = 033ULL, c = 07ULL;
    for(int i=0; i<11; i++) {
        a = (a << 6) | 011ULL;
        b = (b << 6) | 033ULL;
        c = (c << 6) | 07ULL;
    }
    mitHeakArray[0] = a;
    mitHeakArray[1] = b;
    mitHeakArray[2] = c;
    okMitHeak = true;
    return 0;

}
int bitCount ( ULL n ) {
    static int a = initMitHeakArray();
    ULL tmp = n - ((n >> 1) & mitHeakArray[1]) - ((n >> 2) & mitHeakArray[0]);
    return ( (tmp + (tmp >> 3) ) & mitHeakArray[2]) % 63;
}

写完上面的代码, 有人提醒我有问题哦。
64位的整数模63怎么可能得到答案呢?
此时我大呼:呀,只顾着套算法模板,把这个遗漏了。

什么情况呢?

应该还记得这段话吧

而我们32位的答案刚好最多是32, 所以可以模上63, 因此下面的公式就诞生了.

从而我们得到下面的公式

tmp % 63 = a + b + c + d + e + f

但是我们这里是 64位整数,答案最多是64位。
这怎么可以模上 63 呢?

这个时候有两个选择:

  1. 找个更宽的位数,比如7位或八位。
  2. 对63位和64位特殊判断。

现在我们分别来看看吧。

特殊判断 MIT HACKMEM

怎么特殊判断呢?
当 n 的答案是64时,就是全1了。
所以这个好判断,直接判断 ~n 是不是 0 啦。

比如下面的代码

if(~n == 0){
    return 64;
}

那对于63个1的情况,会怎么样呢?
首先没有进位,只是 答案应该是 63, 而我们模上 63, 得到一个答案0了。
所以我们判断这个数字到底是不是0就可以了。

if(ans == 0 && n != 0){
    ans = 63;
}

这样我们就可以得到答案了。
完整代码如下:

int bitCount ( ULL n ) {
    if(~n == 0){
        return 64;
    }
    static int a = initMitHeakArray();
    ULL tmp = n - ((n >> 1) & mitHeakArray[1]) - ((n >> 2) & mitHeakArray[0]);
    int ans = ( (tmp + (tmp >> 3) ) & mitHeakArray[2]) % 63;
    if(ans == 0 && n != 0){
        ans = 63;
    }
    return ans;
}

加宽的 MIT HACKMEM

64位,至少需要7位。
我们为了方便分治,只好采用8位了。

8位不刚好是二分的节奏嘛。
于是我们可以把二分的代码搬过来。

我们要做的是什么呢?
将64二进制数据按每八个一组,这8位并储存这八位1的个数。
做到这个很简单的,只需要第一次计算二位的,然后四位的,最后8位的即可。

for(uint i=0,j=1; i<3; i++,j<<=1) {
    x=(x & mitHeakMask[i]) + ((x>>j) & mitHeakMask[i]);
}

这样进行三次就得到每8位的1的个数了。

然后我们可以得到这个公式

x = a*256^8 + b*256^7 + c*256^6 + d*256^5 + e*256^4 + f*256^3 + g*256^2 + h*256^1 + i*256;

我们需要得到的是 a + b + c + d + e + f + g + h + i 的和。
这个和远远小于 255, 于是我们可以取模 255 了。

  x % 255
= (a*256^8 + b*256^7 + c*256^6 + d*256^5 + e*256^4 + f*256^3 + g*256^2 + h*256^1 + i*256) % 255
= ((a*256^8)  % 255 + (b*256^7)  % 255 + (c*256^6)  % 255 + (d*256^5)  % 255 + (e*256^4)  % 255 + (f*256^3)  % 255 + (g*256^2)  % 255 + (h*256^1)  % 255 + i*256) % 255
= (a + b + c + d + e + f + g + h + i) % 255
= a + b + c + d + e + f + g + h + i

所以答案就是 ans = x % 255

完整代码如下:

关于测试代码可以参考这里

ULL mitHeakMask[3];
bool okMitHeak = false;
int initMitHeakMask() {
    if(okMitHeak) {
        return true;
    }
    ULL a = 0x55ULL, b = 0x33ULL, c = 0x0FULL;
    for(int i=0; i<8; i++) {
        a = (a << 8) | 0x55ULL;
        b = (b << 8) | 0x33ULL;
        c = (c << 8) | 0x0FULL;
    }
    mitHeakMask[0] = a;
    mitHeakMask[1] = b;
    mitHeakMask[2] = c;
    okMitHeak = true;
    return 0;
}
uint countbits(ULL x) {
    int a = initMitHeakMask();
    for(uint i=0,j=1; i<3; i++,j<<=1) {
        x=(x & mitHeakMask[i]) + ((x>>j) & mitHeakMask[i]);
    }
    return ans%255;
}

《完》

本文首发于公众号:天空的代码世界,微信号:tiankonguse
如果你想留言,可以在微信里面关注公众号进行留言。

关注公众号,接收最新消息

tiankonguse +
穿越