快速可分性测试(2,3,4,5,...,16)？

Question

什么是最快的可分性测试？比如说,给定一个小端架构和一个32位有符号整数:如何计算得非常快,一个数字可被2,3,4,5整除,......最多16？

警告:给定的代码仅为示例.每一行都是独立的!使用模运算的明显解决方案在许多处理器上都很慢,这些处理器没有DIV硬件(像许多ARM一样).有些编译器也无法进行这样的优化(例如,如果divisor是函数的参数或依赖于某些东西).

Divisible_by_1 = do();
Divisible_by_2 = if (!(number & 1)) do();
Divisible_by_3 = ?
Divisible_by_4 = ?
Divisible_by_5 = ?
Divisible_by_6 = ?
Divisible_by_7 = ?
Divisible_by_8 = ?
Divisible_by_9 = ?
Divisible_by_10 = ?
Divisible_by_11 = ?
Divisible_by_12 = ?
Divisible_by_13 = ?
Divisible_by_14 = ?
Divisible_by_15 = ?
Divisible_by_16 = if(!number & 0x0000000F) do();

和特殊情况:

Divisible_by_2k = if(number & (tk-1)) do();  //tk=2**k=(2*2*2*...) k times

Answer 1

在每种情况下(包括2整除):

if (number % n == 0) do();

使用低阶位掩码进行定义只是模糊处理,使用现代编译器并不比以可读方式编写代码更快.

如果你必须测试所有的情况,你可以通过将一些案例放在if另一个案例中来提高绩效:例如,如果2的可分性已经失败,那么它就没有必要测试可分性4.

Answer 2

在AT ALL中找出除法指令的替代品(包括x86/x64上的模数)并不是一个坏主意,因为它们非常慢.比大多数人意识到的更慢(甚至更慢).那些暗示"%n",其中n是变量的人给出了愚蠢的建议,因为它总是会导致使用除法指令.另一方面,"%c"(其中c是常量)将允许编译器确定其库中可用的最佳算法.有时这将是分裂指令,但很多时候它不会.

在本文档中, TorbjörnGranlund显示无符号32位多元化所需的时钟周期比率:在Sandybridge上为4:26(6.5x),在K10上为3:45(15x).对于64位,各自的比率是4:92(23x)和5:77(14.4x).

"L"列表示延迟."T"列表示吞吐量.这与处理器并行处理多个指令的能力有关.Sandybridge可以每隔一个周期发出一个32位乘法,或者每个周期发出一个64位乘法.对于K10,相应的吞吐量相反.对于分区,K10需要在它开始另一个序列之前完成整个序列.我怀疑Sandybridge是一样的.

使用K10作为示例意味着在32位除法(45)所需的周期期间,可以发出相同数量(45)的乘法,并且倒数第二个和最后一个将完成一个和两个分裂完成后的时钟周期.可以在45次乘法中完成大量工作.

值得注意的是,随着从K8-K9到K10的演变,div的效率降低:从32到64位的39到45和71到77个时钟周期.

格兰隆德在gmplib.org和斯德哥尔摩皇家理工学院的页面包含更多好东西,其中一些已被纳入gcc编译器.

Answer 3

正如@James所提到的,让编译器为您简化它.如果n是常量,任何下降编译器都能够识别模式并将其更改为更有效的等效项.

例如,代码

#include <stdio.h>

int main() {
    size_t x;
    scanf("%u\n", &x);
    __asm__ volatile ("nop;nop;nop;nop;nop;");
    const char* volatile foo = (x%3 == 0) ? "yes" : "no";
    __asm__ volatile ("nop;nop;nop;nop;nop;");
    printf("%s\n", foo);
    return 0;
}

用g ++ - 4.5 -O3编译,相关部分x%3 == 0将成为

mov    rcx,QWORD PTR [rbp-0x8]   # rbp-0x8 = &x
mov    rdx,0xaaaaaaaaaaaaaaab
mov    rax,rcx
mul    rdx
lea    rax,"yes"
shr    rdx,1
lea    rdx,[rdx+rdx*2]
cmp    rcx,rdx
lea    rdx,"no"
cmovne rax,rdx
mov    QWORD PTR [rbp-0x10],rax

其中,翻译成C代码,意思是

(hi64bit(x * 0xaaaaaaaaaaaaaaab) / 2) * 3 == x ? "yes" : "no"
// equivalatent to:                 x % 3 == 0 ? "yes" : "no"

这里没有分工.(注意0xaaaaaaaaaaaaaaab == 0x20000000000000001L/3)

编辑:

• 魔术常数0xaaaaaaaaaaaaaa可以在http://www.hackersdelight.org/magic.htm中计算出来.
• 对于形式为2 ⁿ - 1的除数,请查看http://graphics.stanford.edu/~seander/bithacks.html#ModulusDivision

Answer 4

脸上有点舌头,但假设你得到了其余的答案:

Divisible_by_6  = Divisible_by_3 && Divisible_by_2;
Divisible_by_10 = Divisible_by_5 && Divisible_by_2;
Divisible_by_12 = Divisible_by_4 && Divisible_by_3;
Divisible_by_14 = Divisible_by_7 && Divisible_by_2;
Divisible_by_15 = Divisible_by_5 && Divisible_by_3;

Answer 5

假设number是unsigned(32位).然后以下是非常快速的方法来计算高达16的可分性.(我没有测量,但汇编代码表示如此.)

bool divisible_by_2 = number % 2 == 0;
bool divisible_by_3 = number * 2863311531u <= 1431655765u;
bool divisible_by_4 = number % 4 == 0;
bool divisible_by_5 = number * 3435973837u <= 858993459u;
bool divisible_by_6 = divisible_by_2 && divisible_by_3;
bool divisible_by_7 = number * 3067833783u <= 613566756u;
bool divisible_by_8 = number % 8 == 0;
bool divisible_by_9 = number * 954437177u <= 477218588u;
bool divisible_by_10 = divisible_by_2 && divisible_by_5;
bool divisible_by_11 = number * 3123612579u <= 390451572u;
bool divisible_by_12 = divisible_by_3 && divisible_by_4;
bool divisible_by_13 = number * 3303820997u <= 330382099u;
bool divisible_by_14 = divisible_by_2 && divisible_by_7;
bool divisible_by_15 = number * 4008636143u <= 286331153u;
bool divisible_by_16 = number % 16 == 0;

关于按d以下规则划分的可分性:

• 当d是2的幂:

  正如指出的由詹姆斯·甘孜,你可以使用is_divisible_by_d = (number % d == 0).编译器非常聪明,可以实现这一点,因为(number & (d - 1)) == 0它非常有效但是模糊不清.

但是,当d不是2的幂时,看起来上面显示的混淆比当前编译器更有效.(稍后会详细介绍).

• 什么时候d奇怪:

  该技术采用以下形式is_divisible_by_d = number * a <= b,其中a与b被巧妙地求出的常数.请注意,我们需要的只是1次乘法和1次比较:

• 什么时候d是偶数但不是2的幂:

  然后,写d = p * q在那里p是2的幂和q是奇数,并使用"在颊舌"所建议unpythonic,即,is_divisible_by_d = is_divisible_by_p && is_divisible_by_q.同样,仅执行1次乘法(在计算中is_divisible_by_q).

许多编译器(我已使用godbolt测试过clang 5.0.0,gcc 7.3,icc 18和msvc 19 )替换number % d == 0为(number / d) * d == number.他们使用一种聪明的技术(参见Olof Forshell的答案中的参考文献)用乘法和位移代替除法.他们最终做了2次乘法.相反,上述技术仅执行1次乘法.

2018年10月1日更新

看起来上面的算法很快就会进入GCC(已经在主干中):

https://gcc.gnu.org/bugzilla/show_bug.cgi?id=82853

海湾合作委员会的实施似乎更有效率.实际上,上面的实现有三个部分:1)除数的偶数部分的可分性; 2)除数的奇数部分的可除性; 3)&&连接前两个步骤的结果.通过使用在标准C++(ror)中无法有效使用的汇编程序指令,GCC将这三个部分组合成一个部分,这与奇数部分的可分性非常相似.好东西!有了这个实现,最好(为了清晰度和性能)可以回溯到%所有时间.

Answer 6

这些数字的LCM似乎是720720.它非常小,因此您可以执行单模运算并将余数用作预计算LUT中的索引.

Answer 7

首先,我提醒你,二进制形式的bn ... b2b1b0中的数字有值:

number = bn*2^n+...+b2*4+b1*2+b0

现在,当你说数字%3时,你有:

number%3 =3= bn*(2^n % 3)+...+b2*1+b1*2+b0

(我用= 3 =表示同余模3).另请注意b1*2 =3= -b1*1

现在我将使用+和 - 编写所有16个除法并且可能是乘法(注意乘法可以写为移位或相同值的和移位到不同的位置.例如,5*x意味着x+(x<<2)你x只计算一次)

让我们调用数字n,让我们说Divisible_by_i是一个布尔值.作为一个中间值,想象Congruence_by_i是一个与n模数一致的值i.

另外,假设n0n表示零n1位,表示位1等,即

ni = (n >> i) & 1;

Congruence_by_1 = 0
Congruence_by_2 = n&0x1
Congruence_by_3 = n0-n1+n2-n3+n4-n5+n6-n7+n8-n9+n10-n11+n12-n13+n14-n15+n16-n17+n18-n19+n20-n21+n22-n23+n24-n25+n26-n27+n28-n29+n30-n31
Congruence_by_4 = n&0x3
Congruence_by_5 = n0+2*n1-n2-2*n3+n4+2*n5-n6-2*n7+n8+2*n9-n10-2*n11+n12+2*n13-n14-2*n15+n16+2*n17-n18-2*n19+n20+2*n21-n22-2*n23+n24+2*n25-n26-2*n27+n28+2*n29-n30-2*n31
Congruence_by_7 = n0+2*n1+4*n2+n3+2*n4+4*n5+n6+2*n7+4*n8+n9+2*n10+4*n11+n12+2*n13+4*n14+n15+2*n16+4*n17+n18+2*n19+4*n20+n21+2*n22+4*n23+n24+2*n25+4*n26+n27+2*n28+4*n29+n30+2*n31
Congruence_by_8 = n&0x7
Congruence_by_9 = n0+2*n1+4*n2-n3-2*n4-4*n5+n6+2*n7+4*n8-n9-2*n10-4*n11+n12+2*n13+4*n14-n15-2*n16-4*n17+n18+2*n19+4*n20-n21-2*n22-4*n23+n24+2*n25+4*n26-n27-2*n28-4*n29+n30+2*n31
Congruence_by_11 = n0+2*n1+4*n2+8*n3+5*n4-n5-2*n6-4*n7-8*n8-5*n9+n10+2*n11+4*n12+8*n13+5*n14-n15-2*n16-4*n17-8*n18-5*n19+n20+2*n21+4*n22+8*n23+5*n24-n25-2*n26-4*n27-8*n28-5*n29+n30+2*n31
Congruence_by_13 = n0+2*n1+4*n2+8*n3+3*n4+6*n5-n6-2*n7-4*n8-8*n9-3*n10-6*n11+n12+2*n13+4*n14+8*n15+3*n16+6*n17-n18-2*n19-4*n20-8*n21-3*n22-6*n3+n24+2*n25+4*n26+8*n27+3*n28+6*n29-n30-2*n31
Congruence_by_16 = n&0xF

或者当分解时:

Congruence_by_1 = 0
Congruence_by_2 = n&0x1
Congruence_by_3 = (n0+n2+n4+n6+n8+n10+n12+n14+n16+n18+n20+n22+n24+n26+n28+n30)-(n1+n3+n5+n7+n9+n11+n13+n15+n17+n19+n21+n23+n25+n27+n29+n31)
Congruence_by_4 = n&0x3
Congruence_by_5 = n0+n4+n8+n12+n16+n20+n24+n28-(n2+n6+n10+n14+n18+n22+n26+n30)+2*(n1+n5+n9+n13+n17+n21+n25+n29-(n3+n7+n11+n15+n19+n23+n27+n31))
Congruence_by_7 = n0+n3+n6+n9+n12+n15+n18+n21+n24+n27+n30+2*(n1+n4+n7+n10+n13+n16+n19+n22+n25+n28+n31)+4*(n2+n5+n8+n11+n14+n17+n20+n23+n26+n29)
Congruence_by_8 = n&0x7
Congruence_by_9 = n0+n6+n12+n18+n24+n30-(n3+n9+n15+n21+n27)+2*(n1+n7+n13+n19+n25+n31-(n4+n10+n16+n22+n28))+4*(n2+n8+n14+n20+n26-(n5+n11+n17+n23+n29))
// and so on

如果这些值最终为负值,请将其添加,i直到它们变为正值.

现在你应该做的是通过我们刚才做的同一过程递归地提供这些值,直到Congruence_by_i变得小于i(显然>= 0).这与我们想要找到3或9的数字的余数时的情况类似,还记得吗？总结数字,如果它有多个数字,那么一些数字会再次增加结果的数字,直到你只得到一个数字.

现在i = 1, 2, 3, 4, 5, 7, 8, 9, 11, 13, 16:

Divisible_by_i = (Congruence_by_i == 0);

其余的:

Divisible_by_6 = Divisible_by_3 && Divisible_by_2;
Divisible_by_10 = Divisible_by_5 && Divisible_by_2;
Divisible_by_12 = Divisible_by_4 && Divisible_by_3;
Divisible_by_14 = Divisible_by_7 && Divisible_by_2;
Divisible_by_15 = Divisible_by_5 && Divisible_by_3;

编辑:请注意,从一开始就可以避免一些添加.例如n0+2*n1+4*n2是一样的n&0x7,同样n3+2*n4+4*n5是(n>>3)&0x7从而与每个公式,你不必让每一位单独,我写的一样,在操作的清晰度和相似的缘故.要优化每个公式,您应该自己处理; 组操作数和分解操作.

Answer 8

您应该使用(i%N)== 0作为您的测试.

我的编译器(gcc的一个相当旧的版本)为我尝试的所有情况生成了良好的代码.比特测试是合适的,它就是这样做的.在N是常数的情况下,它并没有为任何情况产生明显的"鸿沟",它总是使用一些"技巧".

只需让编译器为您生成代码,它几乎肯定会比您更了解机器的体系结构:)这些都是简单的优化,您不太可能想到比编译器更好的东西.

这是一个有趣的问题.我不能列出编译器为每个常量使用的技巧,因为我必须在另一台计算机上编译.但是如果没有人打败我,我会稍后更新这个回复:)

Answer 9

这可能对代码没有帮助,但在某些情况下,有一个巧妙的技巧可以帮助你做到这一点:

除以3:对于以十进制表示的数字,您可以对所有数字求和,并检查总和是否可被3整除.

示例:12345 => 1+2+3+4+5 = 15 => 1+5 = 6,可被3整除(3 x 4115 = 12345).

更有趣的是,相同的技术适用于X-1的所有因子,其中X是表示数字的基础.因此对于十进制数,您可以检查除以3或9.对于十六进制,您可以检查除以3,5或15.对于八进制数,您可以检查除以7.