23

下限除法是当结果总是下限(朝向 -∞),而不是朝向 0:

部门类型

是否可以在 C/C++ 中有效地实现下限或欧几里得整数除法?

(显而易见的解决方案是检查股息的符号)

4

6 回答 6

8

我编写了一个测试程序来对这里提出的想法进行基准测试:

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

#define N 10000000
#define M 100

int dividends[N], divisors[N], results[N];

__forceinline int floordiv_signcheck(int a, int b)
{
    return (a<0 ? a-(b-1) : a) / b;
}

__forceinline int floordiv_signcheck2(int a, int b)
{
    return (a - (a<0 ? b-1 : 0)) / b;
}

__forceinline int floordiv_signmultiply(int a, int b)
{
    return (a + (a>>(sizeof(a)*8-1))*(b-1)) / b;
}

__forceinline int floordiv_floatingpoint(int a, int b)
{
    // I imagine that the call to floor can be replaced to a cast
    // if you can get FPU rounding control to work (I couldn't).
    return floor((double)a / b);
}

void main()
{
    for (int i=0; i<N; i++)
    {
        dividends[i] = rand();
        do
            divisors[i] = rand();
        while (divisors[i]==0);
    }

    LARGE_INTEGER t0, t1;

    QueryPerformanceCounter(&t0);
    for (int j=0; j<M; j++)
        for (int i=0; i<N; i++)
            results[i] = floordiv_signcheck(dividends[i], divisors[i]);
    QueryPerformanceCounter(&t1);
    printf("signcheck    : %9llu\n", t1.QuadPart-t0.QuadPart);

    QueryPerformanceCounter(&t0);
    for (int j=0; j<M; j++)
        for (int i=0; i<N; i++)
            results[i] = floordiv_signcheck2(dividends[i], divisors[i]);
    QueryPerformanceCounter(&t1);
    printf("signcheck2   : %9llu\n", t1.QuadPart-t0.QuadPart);

    QueryPerformanceCounter(&t0);
    for (int j=0; j<M; j++)
        for (int i=0; i<N; i++)
            results[i] = floordiv_signmultiply(dividends[i], divisors[i]);
    QueryPerformanceCounter(&t1);
    printf("signmultiply : %9llu\n", t1.QuadPart-t0.QuadPart);

    QueryPerformanceCounter(&t0);
    for (int j=0; j<M; j++)
        for (int i=0; i<N; i++)
            results[i] = floordiv_floatingpoint(dividends[i], divisors[i]);
    QueryPerformanceCounter(&t1);
    printf("floatingpoint: %9llu\n", t1.QuadPart-t0.QuadPart);
}

结果:

signcheck    :  61458768
signcheck2   :  61284370
signmultiply :  61625076
floatingpoint: 287315364

所以,根据我的结果,检查标志是最快的:

(a - (a<0 ? b-1 : 0)) / b
于 2010-11-05T22:25:53.177 回答
4

五年后我将重新审视这个问题,因为这也与我有关。我对 x86-64 的两个纯 C 版本和两个内联汇编版本进行了一些性能测量,结果可能很有趣。

测试的地板分割变体是:

  1. 我已经使用了一段时间的实现;
  2. 上述仅使用一种划分的轻微变体;
  3. 前一个,但在 inline-assembly 中手动实现;和
  4. CMOV在汇编中实现的版本。

以下是我的基准程序:

#include <stdio.h>
#include <stdlib.h>
#include <sys/time.h>

#ifndef VARIANT
#define VARIANT 3
#endif

#if VARIANT == 0
#define floordiv(a, b) (((a) < 0)?((((a) + 1) / (b)) - 1):((a) / (b)))
#elif VARIANT == 1
#define floordiv(a, b) ((((a) < 0)?((a) - ((b) - 1)):(a)) / (b))
#elif VARIANT == 2
#define floordiv(a, b) ({                                   \
    int result;                                             \
    asm("test %%eax, %%eax; jns 1f; sub %1, %%eax;"         \
        "add $1, %%eax; 1: cltd; idivl %1;"                 \
        : "=a" (result)                                     \
        : "r" (b),                                          \
          "0" (a)                                           \
        : "rdx");                                           \
    result;})
#elif VARIANT == 3
#define floordiv(a, b) ({                                           \
    int result;                                                     \
    asm("mov %%eax, %%edx; sub %1, %%edx; add $1, %%edx;"           \
        "test %%eax, %%eax; cmovs %%edx, %%eax; cltd;"              \
        "idivl %1;"                                                 \
        : "=a" (result)                                             \
        : "r" (b),                                                  \
          "0" (a)                                                   \
        : "rdx");                                                   \
    result;})
#endif

double ntime(void)
{
    struct timeval tv;

    gettimeofday(&tv, NULL);
    return(tv.tv_sec + (((double)tv.tv_usec) / 1000000.0));
}

void timediv(int n, int *p, int *q, int *r)
{
    int i;

    for(i = 0; i < n; i++)
        r[i] = floordiv(p[i], q[i]);
}

int main(int argc, char **argv)
{
    int n, i, *q, *p, *r;
    double st;

    n = 10000000;
    p = malloc(sizeof(*p) * n);
    q = malloc(sizeof(*q) * n);
    r = malloc(sizeof(*r) * n);
    for(i = 0; i < n; i++) {
        p[i] = (rand() % 1000000) - 500000;
        q[i] = (rand() % 1000000) + 1;
    }

    st = ntime();
    for(i = 0; i < 100; i++)
        timediv(n, p, q, r);
    printf("%g\n", ntime() - st);
    return(0);
}

gcc -march=native -Ofast使用 GCC 4.9.2 编译了它,在我的 Core i5-2400 上,结果如下。每次运行的结果都是相当可重复的——至少它们总是以相同的顺序降落。

  • 变体 0:7.21 秒
  • 变体 1:7.26 秒
  • 变体 2:6.73 秒
  • 变体 3:4.32 秒

因此CMOV,至少,该实施将其他人从水中吹了出来。令我惊讶的是,变体 2 以相当大的优势超越了它的纯 C 版本(变体 1)。我原以为编译器应该能够发出至少和我一样高效的代码。

以下是一些其他平台,用于比较:

AMD 速龙 64 X2 4200+、GCC 4.7.2:

  • 变体 0:26.33 秒
  • 变体 1:25.38 秒
  • 变体 2:25.19 秒
  • 变体 3:22.39 秒

至强 E3-1271 v3,GCC 4.9.2:

  • 变体 0:5.95 秒
  • 变体 1:5.62 秒
  • 变体 2:5.40 秒
  • 变体 3:3.44 秒

作为最后一点,我或许应该警告不要CMOV过于认真地对待版本的明显性能优势,因为在现实世界中,其他版本中的分支可能不会像这个基准测试那样完全随机,如果分支预测器可以做一个合理的工作,分支版本可能会变得更好。然而,实际情况在很大程度上取决于实际使用的数据,因此尝试进行任何通用基准测试可能毫无意义。

于 2016-04-22T12:40:08.283 回答
2

提出一些免费的分支来根据符号更正结果可能会更有效,因为分支很昂贵。

请参阅Hacker's Delight第 2 章第 20页,了解如何访问该标志。

于 2010-11-05T06:55:00.963 回答
1

请注意:sar当涉及到 2 的幂时,x86 指令执行地板除法。

于 2012-01-14T05:38:18.143 回答
0

是否可以在 C/C++ 中有效地实现下限或欧几里得整数除法?

是的。

(显而易见的解决方案是检查股息的符号)

我完全同意,并且很难相信存在明显更快的替代方案。

于 2010-11-05T02:53:31.180 回答
0

由于 IEEE-754 将向 -inf 的舍入指定为所需的舍入模式之一,我想您的问题的答案是肯定的。但也许您可以解释一下,如果您正在编写编译器,您是想知道如何实现该过程,还是想知道如何使用特定的编译器来执行该操作?

于 2010-11-04T23:48:03.360 回答