计算快速对数2的上限

43
有一个快速的方法可以计算输入整数 i 的 (long int) ceiling(log_2(i)) 值吗?输入和输出都是 64 位整数。支持有符号或无符号整数的解决方案都可行。我怀疑最好的方法将会是使用一种类似于这里找到的位操作方法,但与其自己尝试,我想使用已经经过充分测试的方法。通用的解决方案适用于所有正数。
例如,对于 2、3、4、5、6、7、8 这些值,应该返回 1、2、2、3、3、3、3。
编辑:目前看来最好的方法是使用任意数量的快速现有位运算或寄存器方法来计算整数/底数为 2 的对数(MSB 的位置),然后如果输入不是 2 的幂,则加 1。用于检查 2 的幂的快速位检查为(n&(n-1))
编辑 2:关于整数对数和前导零的方法,Henry S.Warren 在 Hacker's Delight 的第 5-3 和 11-4 节提供了很好的参考。这是我找到的最完整的参考资料。
编辑 3:这种技术看起来很有前途:https://dev59.com/UnA75IYBdhLWcg3weZDu#51351885 编辑 4:C23 显然正在添加 stdc_first_(leading/trailing)_(one/zero)。

必须准确地处理所有大于1且小于一个很大的数字(比如2^63或2^62)的值。 - kevinlawler
请看下面我的答案。我提供了说明和可以帮助您实现此操作的代码。 - Mahmoud Al-Qudsi
1
过去,我通常使用查找表和位操作的组合,类似于位操作页面链接中的内容。 - TechNeilogy
如果你只处理正值,一个简单的处理舍入的方法是找到 ((x << 1) - 1) 的最高有效位。你需要特殊处理 x == 0,并且如果最高位被设置,你会溢出,但这种方法可能比其他一些舍入技术更快。 - tomlogic
2
在C++20中,只需使用std::bit_ceil即可。不幸的是,这个问题是关于C语言的。 - phuclv
所以,你需要购买一个非常大的硬盘并创建一个大表格。 - Hunter Kohler
15个回答

23

如果您可以限制自己使用gcc,则有一组内置函数可返回前导零位数,并且可以在一些工作的帮助下用于实现您想要的功能:

int __builtin_clz (unsigned int x)
int __builtin_clzl (unsigned long)
int __builtin_clzll (unsigned long long)

23

这个算法已经发布过了,但下面的实现非常紧凑,并且应该能够优化为无分支代码。

int ceil_log2(unsigned long long x)
{
  static const unsigned long long t[6] = {
    0xFFFFFFFF00000000ull,
    0x00000000FFFF0000ull,
    0x000000000000FF00ull,
    0x00000000000000F0ull,
    0x000000000000000Cull,
    0x0000000000000002ull
  };

  int y = (((x & (x - 1)) == 0) ? 0 : 1);
  int j = 32;
  int i;

  for (i = 0; i < 6; i++) {
    int k = (((x & t[i]) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;
  }

  return y;
}


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

int main(int argc, char *argv[])
{
  printf("%d\n", ceil_log2(atol(argv[1])));

  return 0;
}

1
我已经将此标记为最佳答案。对于跟进讨论的任何人来说,目前这是最佳答案的原因是汇编语言解决方案是特定于平台的。 - kevinlawler
5
能否重写一下,使用有意义的变量名,这样就可以理解它是如何工作的了? - anon
为什么你使用int而不是unsigned - sebastian
@sebastian 答案总是在0到64之间,因此任何整数类型都可以使用。由于可以使用任何整数类型,我选择了默认类型“int”。 - dgobbi

13

如果你在Windows上为64位处理器进行编译,我认为这应该能够正常工作。_BitScanReverse64是一个内置函数。

#include <intrin.h>
__int64 log2ceil( __int64 x )
{
  unsigned long index;
  if ( !_BitScanReverse64( &index, x ) )
     return -1LL; //dummy return value for x==0

  // add 1 if x is NOT a power of 2 (to do the ceil)
  return index + (x&(x-1)?1:0);
}

对于32位系统,你可以通过 1 或者 2 次调用 _BitScanReverse 来模拟 _BitScanReverse64。 检查 x 的高32位 ((long*)&x)[1],然后如果需要的话再检查低32位 ((long*)&x)[0]。


我认为这个问题不仅限于Windows,但是到目前为止你的回答最接近正确答案了。简单来说,我们可以使用一种快速位运算方法来检查一个数是否为2的幂(如上所示)。我们可以将这种方法与寄存器方法结合使用,以确定MSB的位置并检索答案。 - kevinlawler
1
如果您只想在x86处理器上运行代码,您可以通过一些汇编语言自行进行位扫描操作(具体来说是BSR指令)。 @highperformance +1 - casablanca
@casablanca:最好使用编译器理解的内置函数,而不是实际的内联汇编。编译器可以通过 _BitScanReverse64 / __builtin_clzll 进行常量传播,但不能通过内联汇编语句进行常量传播。此外,MSVC 内联汇编对于包装单个指令来说是完全垃圾的。(GNU C 内联汇编做得还可以,但仍然没有常量传播)。此外,__builtin_clzll 是可移植的 GNU C,并且将编译为任何目标机器上的优秀代码。 - Peter Cordes

8
我所知道的最快方法是使用一个快速的 log2 函数进行向下取整,结合输入值的无条件调整来处理向上取整的情况,就像下面展示的 lg_down() 函数一样。
/* base-2 logarithm, rounding down */
static inline uint64_t lg_down(uint64_t x) {
  return 63U - __builtin_clzl(x);
}

/* base-2 logarithm, rounding up */
static inline uint64_t lg_up(uint64_t x) {
  return lg_down(x - 1) + 1;
}

基本上将结果向下取整后加1对于除了2的幂次方数以外的所有值都是正确的(因为在这种情况下,floor和ceil方法应该返回相同的答案),因此只需从输入值中减去1来处理这种情况(它不会改变其他情况的答案)并将结果加1。

这通常比通过显式检查是否为2的幂次方数来调整值的方法稍微快一些(例如,添加一个!!(x&(x-1))项)。 它避免了任何比较和条件操作或分支,更容易进行内联,更适合矢量化等。

这依赖于大多数CPU使用clang / icc / gcc内置的__builtin_clzl提供的“计算前导位数”的功能,但其他平台也提供类似的功能(例如Visual Studio中的BitScanReverse内置函数)。

很不幸,这会对log(1)返回错误的答案,因为它会导致基于gcc文档而言的未定义行为__builtin_clzl(0)。当然,总体上,“计算前导零”函数在零处有完全定义的行为,但gcc内置的方式是这样定义的,因为在x86平台的BMI ISA扩展之前,它将使用bsr指令,而该指令本身具有未定义的行为。 如果您知道自己拥有lzcnt指令,则可以直接使用lzcnt内部函数来解决此问题。除了x86之外的大多数平台一开始就没有经历过bsr错误,并且如果它们有“计算前导零”指令,那么可能也提供访问它们的方法。

这看起来很有前途,因为我已经编辑了主要问题以表明。您是否在相关案例和角落上进行了测试,例如0、1和2,完整的32位(无符号)整数空间,对数输出中每个步骤更改周围的三个整数以及接近uint64_t最大值的值? - kevinlawler
@kevinlawler - 我使用的是32位版本,并测试了所有32位值。然而,关于__builtin_clzl有一个警告 - 传递零是未定义的,这意味着log(1)是未定义的。在实践中,当clz编译为x86上的lzcnt指令或ARM等效指令时,它可以工作,但如果它最终成为bsr,则不行。也许有一些方法可以使其定义但仍然快速,考虑到GCC内部函数的这种行为。如果您知道正在编译具有lzcnt的平台,则可以使用该内置函数。我将更新答案以包含这些警告。 - BeeOnRope
@BeeOnRope,你的尝试很不错,但遗憾的是,我的答案在x86-64上比你的好33%(三个相当复杂的汇编指令与你的四个相比),可以在编译器浏览器中查看差异这里。我不确定为什么GCC和Clang(顺便说一下,后者对此生成了更优秀的结果)会惩罚你的解决方案。 - TheCppZoo
@TheCppZoo - 是的,你可以找到很多类似的性能错误,内在定义会扩展到一些固定的配方,但太晚了,在优化过程列表中删除不必要的操作,所以你只剩下表面上冗余的东西。在GCC上,有时取决于你使用什么“-march”标志,即使选择相同的指令也是如此。关于我的log(1)解决方案失败,是的,这很遗憾。原则上,在汇编级别上使用lzcnt是可行的,但不幸的是必须使用英特尔特定的指令来可靠地实现它。 - BeeOnRope
GCC仍然没有修复这个错误,所以在x86_64上返回无符号值时可能会生成次优代码。这可以通过(目前)未记录的特定于目标的内嵌函数或内联汇编来解决。尽管如此,Clang并没有遭受这个错误的困扰。 - FrankHB
显示剩余3条评论

6
#include "stdafx.h"
#include "assert.h"

int getpos(unsigned __int64 value)
{
    if (!value)
    {
      return -1; // no bits set
    }
    int pos = 0;
    if (value & (value - 1ULL))
    {
      pos = 1;
    }
    if (value & 0xFFFFFFFF00000000ULL)
    {
      pos += 32;
      value = value >> 32;
    }
    if (value & 0x00000000FFFF0000ULL)
    {
      pos += 16;
      value = value >> 16;
    }
    if (value & 0x000000000000FF00ULL)
    {
      pos += 8;
      value = value >> 8;
    }
    if (value & 0x00000000000000F0ULL)
    {
      pos += 4;
      value = value >> 4;
    }
    if (value & 0x000000000000000CULL)
    {
      pos += 2;
      value = value >> 2;
    }
    if (value & 0x0000000000000002ULL)
    {
      pos += 1;
      value = value >> 1;
    }
    return pos;
}

int _tmain(int argc, _TCHAR* argv[])
{    
    assert(getpos(0ULL) == -1); // None bits set, return -1.
    assert(getpos(1ULL) == 0);
    assert(getpos(2ULL) == 1);
    assert(getpos(3ULL) == 2);
    assert(getpos(4ULL) == 2);
    for (int k = 0; k < 64; ++k)
    {
        int pos = getpos(1ULL << k);
        assert(pos == k);
    }
    for (int k = 0; k < 64; ++k)
    {
        int pos = getpos( (1ULL << k) - 1);
        assert(pos == (k < 2 ? k - 1 : k) );
    }
    for (int k = 0; k < 64; ++k)
    {
        int pos = getpos( (1ULL << k) | 1);
        assert(pos == (k < 1 ? k : k + 1) );
    }
    for (int k = 0; k < 64; ++k)
    {
        int pos = getpos( (1ULL << k) + 1);
        assert(pos == k + 1);
    }
    return 0;
}

查找表可以消除最后几个 if 语句,具体取决于可用内存的数量(如果内存越大,可以消除的 if 语句就越多)。 - kevinlawler
1
要使用16位查找表:声明全局变量short getpos_lookup [1 << 16];。 预先填充:int i; for(i=1;i<1<<16;i++) getpos_lookup[i]=log2(i); 然后在16个案例下注释掉8/4/2/1个案例,并放置pos += getpos_lookup[v]; - kevinlawler
实际上,我甚至不确定我的版本是否真的很快。使用 for 循环可能比其他任何方法都更快。 - rwong
1
我进行了一些简单的测试。你的方法比循环(例如while(value>>=1)pos++;)明显更快。将你的方法修改为使用查找表稍微更快,但我不会称之为明显更快。对于我的目的来说,你的方法已经足够快了。然而,如果有人想继续改进它,我会考虑以下几点:1. 用寄存器调用替换MSB检测(可能使用#ifdef语句实现可移植性)。2. 使用一些启发式算法来利用输入的已知分布(例如90%的输入数字都在1000以下)。3. 使用查找表。 - kevinlawler
1
通过切换到常量,如 0xFFFF0000FFFF0000、0xFF00FF00FF00FF00 等,您可以消除移位操作。 - BCS
如果您这样做,您需要屏蔽输入值而不是移位。否则,您无法确定哪个1的集合正在让该值通过。例如,值“0x0000000100100000”通过第一个掩码和您建议的第二个掩码“0xFFFF0000FFFF0000”,即使第32位是最高设置的位。 - plasmoidia

6

使用@egosys提到的gcc内置函数可以构建一些有用的宏。 要快速粗略地计算floor(log2(x)),您可以使用以下方法:

#define FAST_LOG2(x) (sizeof(unsigned long)*8 - 1 - __builtin_clzl((unsigned long)(x)))

如果需要类似于ceil(log2(x))的函数,请使用:

#define FAST_LOG2_UP(x) (((x) - (1 << FAST_LOG2(x))) ? FAST_LOG2(x) + 1 : FAST_LOG2(x))

后者可以进一步优化,使用更多gcc独特的技巧来避免双重调用内置函数,但我不确定您在这里是否需要。

谢谢。我担心的是内置函数的调用。它可以在gcc和clang上编译,这应该涵盖大多数情况。如果我知道它可以在icc上编译,我可能会选择使用它。跨平台兼容性是一个问题。我也不介意按照你建议的方式清理双重调用。(你不能使用(x&(x-1))吗?) - kevinlawler
1
你可以使用#define FAST_LOG2_UP(x) ({ unsigned long log = FAST_LOG2(x); ((x) - (1 << log)) ? log + 1 : log; })来避免多次调用内置函数。请注意,这仅适用于gcc。 - Alexander Amelkin
1
另一个改进可以是消除分支: log + !!(x ^ (1 << log)) - alveko
1
在我的测试中,AND-MINUS方法的两个变体((x&(x-1))?1:0和!!(x&(x-1)))之间的运行时间相差不到1%:可能是相同的。XOR-SHIFT方法MSB+!!(x^(1<<MSB))的速度明显较慢,大约慢了50%。 - kevinlawler

5
以下代码片段是一种安全且可移植的方法,用于扩展纯C方法(例如@dgobbi)以在使用支持编译器(Clang)编译时使用编译器内置函数。将其放置在方法顶部将使方法在可用时使用内置函数。当内置函数不可用时,该方法将回退到标准C代码。
#ifndef __has_builtin
#define __has_builtin(x) 0
#endif

#if __has_builtin(__builtin_clzll) //use compiler if possible
  return ((sizeof(unsigned long long) * 8 - 1) - __builtin_clzll(x)) + (!!(x & (x - 1)));
#endif

4

最快的解决方案:

一个包含63个条目的二叉搜索树。这些条目是从0到63的2的幂。生成函数只需生成一次即可创建整个树。叶子节点表示各幂次数的以2为底的对数(基本上是数字1-63)。

要找到答案,只需将数字输入树中,并导航到大于该项的叶节点。如果叶子节点恰好相等,则结果为叶子值。否则,结果为叶子值+1。

复杂度固定为O(6)。


1
当然,您不一定需要一棵树,只需将二分法用于寻找根来每次缩小间隔即可。 - President James K. Polk
Dave:谢谢。 Greg:当然。对于未接触过编程的人来说,树更容易形象化。 - Mahmoud Al-Qudsi
9
实际上,我认为这种分支在CPU上会比较慢(有个基准测试会很有趣!)。更快的方法是可能的,因为我们可以使用能够并行操作64位的指令。总之,你的方法就像链接中“在O(lg(N))个操作中找到一个N位整数的以2为底的对数”的展开版本。而且... O(6)=O(1)=O(99999)。 - Tom Sirgedas
5
二分查找不需要任何分支,完全可以通过符号位/进位标志的算术运算完成。 - R.. GitHub STOP HELPING ICE

3
我将为您提供x86-64的最快方法,并提供一种通用技术,如果您有一个快速的底部,那么对于<2^63的参数有效,如果您关心全部范围,请参见下文。
我对其他答案的质量感到惊讶,因为它们告诉您如何获得底部,但使用非常昂贵的方式(使用条件语句等!)将底部转换为顶部。
如果您可以快速获取对数底,例如使用__builtin_clzll,则可以很容易地获得底部,如下所示:
unsigned long long log2Floor(unsigned long long x) {
    return 63 - __builtin_clzll(x);
}

unsigned long long log2Ceiling(unsigned long long x) {
    return log2Floor(2*x - 1);
}

它之所以有效,是因为它将结果加1,除非x恰好是2的幂。请参见编译器浏览器中的x86-64汇编器差异,了解另一种类似于此的向上取整实现。
auto log2CeilingDumb(unsigned long x) {
    return log2Floor(x) + (!!(x & (x - 1)));
}

给:

log2Floor(unsigned long): # @log2Floor(unsigned long)
  bsr rax, rdi
  ret
log2CeilingDumb(unsigned long): # @log2CeilingDumb(unsigned long)
  bsr rax, rdi
  lea rcx, [rdi - 1]
  and rcx, rdi
  cmp rcx, 1
  sbb eax, -1
  ret
log2Ceiling(unsigned long): # @log2Ceiling(unsigned long)
  lea rax, [rdi + rdi]
  add rax, -1
  bsr rax, rax
  ret

针对完整的范围,可以参考之前的一个答案:return log2Floor(x - 1) + 1,但这种方法较慢,因为它在x86-64中需要四个指令,而比上面的三个指令多。


3
我认为你对其他答案“因为它们告诉你如何得到地板却将地板转换成非常昂贵的东西”的“低质量”描述是错误的,因为其他答案包括我的答案,它与你的答案基本相同(你只是用* 2替换了一个加1操作)。额外的指令或多或少只是编译器的一个奇怪特性。如果使用gcc,它需要6个指令 - BeeOnRope
2
通常情况下,我期望生成的代码在各种编译器和硬件上都非常相似,除了这个版本对于一半的输入范围返回错误答案,而我的版本在输入“1”时表现出未定义行为(但在使用lzcnt或类似功能时实际上返回正确答案)。 - BeeOnRope

3

我已经对几个 64 位“最高位”的实现进行了基准测试。最不需要分支的代码并不是最快的。

这是我的highest-bit.c源文件:

int highest_bit_unrolled(unsigned long long n)
{
  if (n & 0xFFFFFFFF00000000) {
    if (n & 0xFFFF000000000000) {
      if (n & 0xFF00000000000000) {
        if (n & 0xF000000000000000) {
          if (n & 0xC000000000000000)
            return (n & 0x8000000000000000) ? 64 : 63;
          else
            return (n & 0x2000000000000000) ? 62 : 61;
        } else {
          if (n & 0x0C00000000000000)
            return (n & 0x0800000000000000) ? 60 : 59;
          else
            return (n & 0x0200000000000000) ? 58 : 57;
        }
      } else {
        if (n & 0x00F0000000000000) {
          if (n & 0x00C0000000000000)
            return (n & 0x0080000000000000) ? 56 : 55;
          else
            return (n & 0x0020000000000000) ? 54 : 53;
        } else {
          if (n & 0x000C000000000000)
            return (n & 0x0008000000000000) ? 52 : 51;
          else
            return (n & 0x0002000000000000) ? 50 : 49;
        }
      }
    } else {
      if (n & 0x0000FF0000000000) {
        if (n & 0x0000F00000000000) {
          if (n & 0x0000C00000000000)
            return (n & 0x0000800000000000) ? 48 : 47;
          else
            return (n & 0x0000200000000000) ? 46 : 45;
        } else {
          if (n & 0x00000C0000000000)
            return (n & 0x0000080000000000) ? 44 : 43;
          else
            return (n & 0x0000020000000000) ? 42 : 41;
        }
      } else {
        if (n & 0x000000F000000000) {
          if (n & 0x000000C000000000)
            return (n & 0x0000008000000000) ? 40 : 39;
          else
            return (n & 0x0000002000000000) ? 38 : 37;
        } else {
          if (n & 0x0000000C00000000)
            return (n & 0x0000000800000000) ? 36 : 35;
          else
            return (n & 0x0000000200000000) ? 34 : 33;
        }
      }
    }
  } else {
    if (n & 0x00000000FFFF0000) {
      if (n & 0x00000000FF000000) {
        if (n & 0x00000000F0000000) {
          if (n & 0x00000000C0000000)
            return (n & 0x0000000080000000) ? 32 : 31;
          else
            return (n & 0x0000000020000000) ? 30 : 29;
        } else {
          if (n & 0x000000000C000000)
            return (n & 0x0000000008000000) ? 28 : 27;
          else
            return (n & 0x0000000002000000) ? 26 : 25;
        }
      } else {
        if (n & 0x0000000000F00000) {
          if (n & 0x0000000000C00000)
            return (n & 0x0000000000800000) ? 24 : 23;
          else
            return (n & 0x0000000000200000) ? 22 : 21;
        } else {
          if (n & 0x00000000000C0000)
            return (n & 0x0000000000080000) ? 20 : 19;
          else
            return (n & 0x0000000000020000) ? 18 : 17;
        }
      }
    } else {
      if (n & 0x000000000000FF00) {
        if (n & 0x000000000000F000) {
          if (n & 0x000000000000C000)
            return (n & 0x0000000000008000) ? 16 : 15;
          else
            return (n & 0x0000000000002000) ? 14 : 13;
        } else {
          if (n & 0x0000000000000C00)
            return (n & 0x0000000000000800) ? 12 : 11;
          else
            return (n & 0x0000000000000200) ? 10 : 9;
        }
      } else {
        if (n & 0x00000000000000F0) {
          if (n & 0x00000000000000C0)
            return (n & 0x0000000000000080) ? 8 : 7;
          else
            return (n & 0x0000000000000020) ? 6 : 5;
        } else {
          if (n & 0x000000000000000C)
            return (n & 0x0000000000000008) ? 4 : 3;
          else
            return (n & 0x0000000000000002) ? 2 : (n ? 1 : 0);
        }
      }
    }
  }
}

int highest_bit_bs(unsigned long long n)
{
  const unsigned long long mask[] = {
    0x000000007FFFFFFF,
    0x000000000000FFFF,
    0x00000000000000FF,
    0x000000000000000F,
    0x0000000000000003,
    0x0000000000000001
  };
  int hi = 64;
  int lo = 0;
  int i = 0;

  if (n == 0)
    return 0;

  for (i = 0; i < sizeof mask / sizeof mask[0]; i++) {
    int mi = lo + (hi - lo) / 2;

    if ((n >> mi) != 0)
      lo = mi;
    else if ((n & (mask[i] << lo)) != 0)
      hi = mi;
  }

  return lo + 1;
}

int highest_bit_shift(unsigned long long n)
{
  int i = 0;
  for (; n; n >>= 1, i++)
    ; /* empty */
  return i;
}

static int count_ones(unsigned long long d)
{
  d = ((d & 0xAAAAAAAAAAAAAAAA) >>  1) + (d & 0x5555555555555555);
  d = ((d & 0xCCCCCCCCCCCCCCCC) >>  2) + (d & 0x3333333333333333);
  d = ((d & 0xF0F0F0F0F0F0F0F0) >>  4) + (d & 0x0F0F0F0F0F0F0F0F);
  d = ((d & 0xFF00FF00FF00FF00) >>  8) + (d & 0x00FF00FF00FF00FF);
  d = ((d & 0xFFFF0000FFFF0000) >> 16) + (d & 0x0000FFFF0000FFFF);
  d = ((d & 0xFFFFFFFF00000000) >> 32) + (d & 0x00000000FFFFFFFF);
  return d;
}

int highest_bit_parallel(unsigned long long n)
{
  n |= n >> 1;
  n |= n >> 2;
  n |= n >> 4;
  n |= n >> 8;
  n |= n >> 16;
  n |= n >> 32;
  return count_ones(n);
}

int highest_bit_so(unsigned long long x)
{
  static const unsigned long long t[6] = {
    0xFFFFFFFF00000000ull,
    0x00000000FFFF0000ull,
    0x000000000000FF00ull,
    0x00000000000000F0ull,
    0x000000000000000Cull,
    0x0000000000000002ull
  };

  int y = (((x & (x - 1)) == 0) ? 0 : 1);
  int j = 32;
  int i;

  for (i = 0; i < 6; i++) {
    int k = (((x & t[i]) == 0) ? 0 : j);
    y += k;
    x >>= k;
    j >>= 1;
  }

  return y;
}

int highest_bit_so2(unsigned long long value)
{
  int pos = 0;
  if (value & (value - 1ULL))
  {
    pos = 1;
  }
  if (value & 0xFFFFFFFF00000000ULL)
  {
    pos += 32;
    value = value >> 32;
  }
  if (value & 0x00000000FFFF0000ULL)
  {
    pos += 16;
    value = value >> 16;
  }
  if (value & 0x000000000000FF00ULL)
  {
    pos += 8;
    value = value >> 8;
  }
  if (value & 0x00000000000000F0ULL)
  {
    pos += 4;
    value = value >> 4;
  }
  if (value & 0x000000000000000CULL)
  {
    pos += 2;
    value = value >> 2;
  }
  if (value & 0x0000000000000002ULL)
  {
    pos += 1;
    value = value >> 1;
  }
  return pos;
}

这是highest-bit.h文件:

int highest_bit_unrolled(unsigned long long n);
int highest_bit_bs(unsigned long long n);
int highest_bit_shift(unsigned long long n);
int highest_bit_parallel(unsigned long long n);
int highest_bit_so(unsigned long long n);
int highest_bit_so2(unsigned long long n);

以下是主程序(抱歉,需要复制粘贴):

#include <stdlib.h>
#include <stdio.h>
#include <time.h>
#include "highest-bit.h"

double timedelta(clock_t start, clock_t end)
{
  return (end - start)*1.0/CLOCKS_PER_SEC;
}

int main(int argc, char **argv)
{
  int i;
  volatile unsigned long long v;
  clock_t start, end;

  start = clock();

  for (i = 0; i < 10000000; i++) {
    for (v = 0x8000000000000000; v; v >>= 1)
      highest_bit_unrolled(v);
  }

  end = clock();

  printf("highest_bit_unrolled = %6.3fs\n", timedelta(start, end));

  start = clock();

  for (i = 0; i < 10000000; i++) {
    for (v = 0x8000000000000000; v; v >>= 1)
      highest_bit_parallel(v);
  }

  end = clock();

  printf("highest_bit_parallel = %6.3fs\n", timedelta(start, end));

  start = clock();

  for (i = 0; i < 10000000; i++) {
    for (v = 0x8000000000000000; v; v >>= 1)
      highest_bit_bs(v);
  }

  end = clock();

  printf("highest_bit_bs = %6.3fs\n", timedelta(start, end));

  start = clock();

  for (i = 0; i < 10000000; i++) {
    for (v = 0x8000000000000000; v; v >>= 1)
      highest_bit_shift(v);
  }

  end = clock();

  printf("highest_bit_shift = %6.3fs\n", timedelta(start, end));

  start = clock();

  for (i = 0; i < 10000000; i++) {
    for (v = 0x8000000000000000; v; v >>= 1)
      highest_bit_so(v);
  }

  end = clock();

  printf("highest_bit_so = %6.3fs\n", timedelta(start, end));

  start = clock();

  for (i = 0; i < 10000000; i++) {
    for (v = 0x8000000000000000; v; v >>= 1)
      highest_bit_so2(v);
  }

  end = clock();

  printf("highest_bit_so2 = %6.3fs\n", timedelta(start, end));

  return 0;
}

我已尝试了各种不同的Intel x86盒子,无论是新的还是旧的。

highest_bit_unrolled(展开二进制搜索)始终比highest_bit_parallel(无分支比特运算)快得多。这比highest_bit_bs(二进制搜索循环)更快,而highest_bit_bs又比highest_bit_shift(朴素移位和计数循环)更快。

highest_bit_unrolled也比被接受的SO答案中的解法highest_bit_so以及另一个答案中给出的解法highest_bit_so2都要快。

基准测试循环遍历了覆盖连续位的一位掩码。这样做是为了尝试击败展开二进制搜索中的分支预测,这是现实中的情况:在真实世界的程序中,输入用例不太可能表现出位位置的局部性。

以下是在旧的Intel(R) Core(TM)2 Duo CPU E4500 @ 2.20GHz上的结果:

$ ./highest-bit
highest_bit_unrolled =  6.090s
highest_bit_parallel =  9.260s
highest_bit_bs = 19.910s
highest_bit_shift = 21.130s
highest_bit_so =  8.230s
highest_bit_so2 =  6.960s

在较新的型号 Intel(R) Core(TM) i7-6700K CPU @ 4.00GHz 上:

highest_bit_unrolled =  1.555s
highest_bit_parallel =  3.420s
highest_bit_bs =  6.486s
highest_bit_shift =  9.505s
highest_bit_so =  4.127s
highest_bit_so2 =  1.645s

在更新的硬件上,highest_bit_so2 更接近于 highest_bit_unrolled。顺序不完全相同;现在 highest_bit_so 真正落后,比 highest_bit_parallel 慢。

最快的一个,highest_bit_unrolled 包含了最多的代码和最多的分支。每个不同条件下到达的返回值都有自己专用的代码。

“避免所有分支”的直觉(由于担心分支错误预测)并不总是正确的。现代(甚至不那么现代的)处理器中包含了相当多的技巧,以便不受分支的影响。


P.S. highest_bit_unrolled 是在 2011 年 12 月 在 TXR 语言中引入的(存在错误,已经调试)。

最近,我开始想知道是否可以使用一些更简洁的没有分支的代码来提高速度。

结果让我有些惊讶。

可以说,该代码应该真正地使用#ifdef 来适配 GNU C 并使用一些编译器原语,但就可移植性而言,该版本是保留的。


这是很好的工作。我认为你在这个问题上进入了新领域。尽管如此,重新评估目标可能是值得的。我猜测,在没有证据的情况下,这里的加速取决于基准测试中增量for循环的可预测性,虽然它有其应用,但更现实的模型可能是随机抽样,这会增加产生合法时间的障碍。你可能已经打破了这个问题,因此使其恢复意味着以有意义和有用的方式重新构建它,这可能是可能的,也可能不可能。 - kevinlawler
不错。我怀疑highest_bit_unrolled最快的原因是它小心地发出了不止一个CPU 存储周期。只有取指令周期,CPU才能真正地前进。虽然我不是专家,但这似乎表明总的来说,总线一致性延迟比分支预测惩罚昂贵得多(这符合我的粗略 intuitions)。 - Glenn Slayden
@GlennSlayden 我们需要检查机器代码才能确定,但是对本地变量的存储应该只会转换为寄存器操作。然而,存在数据流依赖关系,这些依赖关系无法并行化。例如,在 highest_bit_so2 函数中,每次将输出累加器相加都取决于检查输入值,而输入值正在被有条件地移位。因此,这个寄存器不能轻易地映射到一组寄存器以并行化代码;并且会出现管道停顿。 - Kaz

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