C：使用查找表在有限的小整数值集上评估函数的最快方法？

Question

我目前正在做一个项目，我想通过调用 C 来优化 Python 中的一些数值计算。

简而言之，我需要为巨大数组x（通常有10^9 条目或更多）中的每个元素计算y[i] = f(x[i]) 的值。这里，x[i] 是介于 -10 和 10 之间的整数，f 是接受 x[i] 并返回双精度的函数。我的问题是 f 但需要 非常长 时间才能以数值稳定的方式进行评估。

为了加快速度，我想将所有2*10 + 1 可能的f(x[i]) 值硬编码到常量数组中，例如：

double table_of_values[] = {f(-10), ...., f(10)};

然后使用“查找表”方法评估f，如下所示：

for (i = 0; i < N; i++) {
    y[i] = table_of_values[x[i] + 11]; //instead of y[i] = f(x[i])
}

由于我不太擅长用 C 编写优化代码，我想知道：

1234563负指数（除了[x[i] + 10 + 1]）？

1. 假设 x[i] 不在 -10 和 10 之间，而是在 -20 和 20 之间。在这种情况下，我仍然可以使用相同的方法，但需要手动对查找表进行硬编码。有没有办法在代码中动态生成查找表，以便我使用相同的方法并允许x[i] 属于变量范围？

Answer 1

生成这样一个具有动态范围值的表格相当容易。

这是一个简单的单表方法：

#include <malloc.h>

#define VARIABLE_USED(_sym) \
    do { \
        if (1) \
            break; \
        if (!! _sym) \
            break; \
    } while (0)

double *table_of_values;
int table_bias;

// use the smallest of these that can contain the values the x array may have
#if 0
typedef int xval_t;
#endif
#if 0
typedef short xval_t;
#endif
#if 1
typedef char xval_t;
#endif

#define XLEN        (1 << 9)
xval_t *x;

// fslow -- your original function
double
fslow(int i)
{
    return 1;  // whatever
}

// ftablegen -- generate variable table
void
ftablegen(double (*f)(int),int lo,int hi)
{
    int len;

    table_bias = -lo;

    len = hi - lo;
    len += 1;

    // NOTE: you can do free(table_of_values) when no longer needed
    table_of_values = malloc(sizeof(double) * len);

    for (int i = lo;  i <= hi;  ++i)
        table_of_values[i + table_bias] = f(i);
}

// fcached -- retrieve cached table data
double
fcached(int i)
{
    return table_of_values[i + table_bias];
}

// fripper -- access x and table arrays
void
fripper(xval_t *x)
{
    double *tptr;
    int bias;
    double val;

    // ensure these go into registers to prevent needless extra memory fetches
    tptr = table_of_values;
    bias = table_bias;

    for (int i = 0;  i < XLEN;  ++i) {
        val = tptr[x[i] + bias];

        // do stuff with val
        VARIABLE_USED(val);
    }
}

int
main(void)
{

    ftablegen(fslow,-10,10);
    x = malloc(sizeof(xval_t) * XLEN);
    fripper(x);

    return 0;
}

这是一种稍微复杂的方法，可以生成许多类似的表：

#include <malloc.h>

#define VARIABLE_USED(_sym) \
    do { \
        if (1) \
            break; \
        if (!! _sym) \
            break; \
    } while (0)

// use the smallest of these that can contain the values the x array may have
#if 0
typedef int xval_t;
#endif
#if 1
typedef short xval_t;
#endif
#if 0
typedef char xval_t;
#endif

#define XLEN        (1 << 9)
xval_t *x;

struct table {
    int tbl_lo;                         // lowest index
    int tbl_hi;                         // highest index
    int tbl_bias;                       // bias for index
    double *tbl_data;                   // cached data
};

struct table ftable1;
struct table ftable2;

double
fslow(int i)
{
    return 1;  // whatever
}

double
f2(int i)
{
    return 2;  // whatever
}

// ftablegen -- generate variable table
void
ftablegen(double (*f)(int),int lo,int hi,struct table *tbl)
{
    int len;

    tbl->tbl_bias = -lo;

    len = hi - lo;
    len += 1;

    // NOTE: you can do free tbl_data when no longer needed
    tbl->tbl_data = malloc(sizeof(double) * len);

    for (int i = lo;  i <= hi;  ++i)
        tbl->tbl_data[i + tbl->tbl_bias] = fslow(i);
}

// fcached -- retrieve cached table data
double
fcached(struct table *tbl,int i)
{
    return tbl->tbl_data[i + tbl->tbl_bias];
}

// fripper -- access x and table arrays
void
fripper(xval_t *x,struct table *tbl)
{
    double *tptr;
    int bias;
    double val;

    // ensure these go into registers to prevent needless extra memory fetches
    tptr = tbl->tbl_data;
    bias = tbl->tbl_bias;

    for (int i = 0;  i < XLEN;  ++i) {
        val = tptr[x[i] + bias];

        // do stuff with val
        VARIABLE_USED(val);
    }
}

int
main(void)
{

    x = malloc(sizeof(xval_t) * XLEN);

    // NOTE: we could use 'char' for xval_t ...
    ftablegen(fslow,-37,62,&ftable1);
    fripper(x,&ftable1);

    // ... but, this forces us to use a 'short' for xval_t
    ftablegen(f2,-99,307,&ftable2);

    return 0;
}

注意事项：

fcached 可以/应该是一个inline 函数以提高速度。请注意，一旦表格被计算一次，fcached(x[i]) 就会非常快。您提到的索引偏移问题[由“偏差”解决]在计算时间上非常小。

虽然x 可能是一个大数组，但f() 值的缓存数组相当小（例如-10 到10）。即使它是（例如）-100 到 100，这仍然是大约 200 个元素。这个小缓存数组[可能] 会保留在硬件内存缓存中，因此访问速度将保持相当快。

因此，对x 进行排序以优化查找表的硬件缓存性能几乎没有[可衡量的] 影响。

x 的访问模式是独立的。如果您以线性方式访问x（例如for (i = 0; i < 999999999; ++i) x[i]），您将获得最佳性能。如果您以半随机方式访问它，它将对硬件缓存逻辑及其保持所需/想要的x 值“缓存热”的能力造成压力

即使使用线性访问，因为x 非常大，当你到达最后时，第一个元素将被从硬件缓存中逐出（例如，大多数 CPU 缓存大约为几个兆字节）

但是，如果 x 仅具有有限范围内的值，则将类型从 int x[...] 更改为 short x[...] 甚至 char x[...] 会将 大小 减少 2 倍 [或4x]。而且， 可以显着提高性能。

更新：我添加了一个 fripper 函数来显示 [据我所知] 在循环中访问表和 x 数组的最快方式。我还添加了一个名为 xval_t 的 typedef 以允许 x 数组消耗更少的空间（即具有更好的硬件缓存性能）。

---

更新 #2：

根据你的 cmets ...

fcached [主要] 用于说明简单/单一访问。但是，它并没有在最后一个例子中使用。

多年来，内联的确切要求各不相同（例如，是外部内联）。现在最好使用：static inline。但是，如果使用c++，它可能会再次不同。有整页专门讨论这个问题。原因是因为在不同的.c 文件中编译，优化打开或关闭时会发生什么。此外，请考虑使用 gcc 扩展名。所以，一直强制内联：

__attribute__((__always_inline__)) static inline

fripper 是最快的，因为它避免在每次循环迭代时重新获取全局变量 table_of_values 和 table_bias。在 fripper 中，编译器优化器将确保它们保留在寄存器中。请参阅我的回答：Is accessing statically or dynamically allocated memory faster? 至于为什么。

然而，我编写了一个使用fcached 的fripper 变体，反汇编代码是相同的[并且最佳]。所以，我们可以无视这一点……或者，我们可以吗？有时，反汇编代码是一种很好的交叉检查，也是唯一确定的方法。在创建完全优化的 C 代码时只是一个额外的项目。关于代码生成，编译器有很多选择，所以有时只是反复试验。

因为基准测试很重要，所以我加入了时间戳记（仅供参考，[AFAIK] 底层的 clock_gettime 调用是 python 的 time.clock() 的基础）。

所以，这是更新的版本：

#include <malloc.h>
#include <time.h>

typedef long long s64;

#define SUPER_INLINE \
    __attribute__((__always_inline__)) static inline

#define VARIABLE_USED(_sym) \
    do { \
        if (1) \
            break; \
        if (!! _sym) \
            break; \
    } while (0)

#define TVSEC           1000000000LL    // nanoseconds in a second
#define TVSECF          1e9             // nanoseconds in a second

// tvget -- get high resolution time of day
// RETURNS: absolute nanoseconds
s64
tvget(void)
{
    struct timespec ts;
    s64 nsec;

    clock_gettime(CLOCK_REALTIME,&ts);

    nsec = ts.tv_sec;
    nsec *= TVSEC;
    nsec += ts.tv_nsec;

    return nsec;
)

// tvgetf -- get high resolution time of day
// RETURNS: fractional seconds
double
tvgetf(void)
{
    struct timespec ts;
    double sec;

    clock_gettime(CLOCK_REALTIME,&ts);

    sec = ts.tv_nsec;
    sec /= TVSECF;
    sec += ts.tv_sec;

    return sec;
)

double *table_of_values;
int table_bias;

double *dummyptr;

// use the smallest of these that can contain the values the x array may have
#if 0
typedef int xval_t;
#endif
#if 0
typedef short xval_t;
#endif
#if 1
typedef char xval_t;
#endif

#define XLEN        (1 << 9)
xval_t *x;

// fslow -- your original function
double
fslow(int i)
{
    return 1;  // whatever
}

// ftablegen -- generate variable table
void
ftablegen(double (*f)(int),int lo,int hi)
{
    int len;

    table_bias = -lo;

    len = hi - lo;
    len += 1;

    // NOTE: you can do free(table_of_values) when no longer needed
    table_of_values = malloc(sizeof(double) * len);

    for (int i = lo;  i <= hi;  ++i)
        table_of_values[i + table_bias] = f(i);
}

// fcached -- retrieve cached table data
SUPER_INLINE double
fcached(int i)
{
    return table_of_values[i + table_bias];
}

// fripper_fcached -- access x and table arrays
void
fripper_fcached(xval_t *x)
{
    double val;
    double *dptr;

    dptr = dummyptr;

    for (int i = 0;  i < XLEN;  ++i) {
        val = fcached(x[i]);

        // do stuff with val
        dptr[i] = val;
    }
}

// fripper -- access x and table arrays
void
fripper(xval_t *x)
{
    double *tptr;
    int bias;
    double val;
    double *dptr;

    // ensure these go into registers to prevent needless extra memory fetches
    tptr = table_of_values;
    bias = table_bias;
    dptr = dummyptr;

    for (int i = 0;  i < XLEN;  ++i) {
        val = tptr[x[i] + bias];

        // do stuff with val
        dptr[i] = val;
    }
}

int
main(void)
{

    ftablegen(fslow,-10,10);
    x = malloc(sizeof(xval_t) * XLEN);
    dummyptr = malloc(sizeof(double) * XLEN);
    fripper(x);
    fripper_fcached(x);

    return 0;
}

Answer 2

您的数组中可以有负索引。 （我不确定这是否在规范中。）如果您有以下代码：

int arr[] = {1, 2 ,3, 4, 5};
int* lookupTable = arr + 3;
printf("%i", lookupTable[-2]);

它将打印出2。

---

这是因为 c 中的数组被定义为指针。如果指针不指向数组的开头，则可以访问指针之前的项。

请记住，如果您必须 malloc() 为 arr 分配内存，您可能无法使用 free(lookupTable) 来释放它。

Answer 3

我真的认为 Craig Estey 在以自动方式构建您的桌子方面走在了正确的轨道上。我只是想添加一个用于查找表格的注释。

如果您知道您将在 Haswell 机器（使用 AVX2）上运行代码，您应该确保您的代码使用 VGATHERDPD，您可以使用 _mm256_i32gather_pd 内在函数。如果你这样做，你的表查找就会飞起来！ （您甚至可以使用 cpuid() 即时检测 avx2，但这是另一回事）

编辑：
让我用一些代码来详细说明：

#include <stdint.h> 
#include <stdio.h> 
#include <immintrin.h> 
/* I'm not sure if you need the alignment */
double table[8]  __attribute__((aligned(16)))= { 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8 };

int main()
{
    int32_t i[4] = { 0,2,4,6 };
    __m128i index = _mm_load_si128( (__m128i*) i );
    __m256d result = _mm256_i32gather_pd( table, index, 8 );

    double* f = (double*)&result;
    printf("%f %f %f %f\n", f[0], f[1], f[2], f[3]);
    return 0;
}

编译运行：

$ gcc --std=gnu99 -mavx2 gathertest.c -o gathertest && ./gathertest
0.100000 0.300000 0.500000 0.700000

这很快！