是否可以优化和更快地将数组乘以 C 中的数字?

Is it possible to optimize and faster multiplying an array by a number in C?

我在 C 中有一个代码,它将数组的每个元素乘以一个数字 (0-9),得到一系列以 10 为基数的数字。

我的问题是此函数花费的时间比我预期的要长 运行。我需要它更快。 我知道在优化我的功能时我的问题是对进位的依赖。如何修改此代码以解决此问题并使代码更快? 解决方案可以使用内在函数或其他专门技术。

到目前为止我最快的版本是:

void ConstMult( uint8_t *V, size_t N, uint8_t digit )
{
  uint8_t CARRY = 0;
  for ( size_t i=0; i< N; ++i )
  {
    V[i] = V[i] * digit + CARRY;
    CARRY = ((uint32_t)V[i] * (uint32_t)0xCCCD) >> 19;
    V[i] -= (CARRY << 3) + (CARRY << 1);
  }
}

但我也尝试了这些较慢的方法:

uint8_t ConstMult( uint8_t *V, size_t N, uint8_t digit )
{
  uint8_t CARRY = 0;
  for ( int i=0; i< N; i++ ) 
  {
    char R = V[i] * digit + CARRY;
    CARRY = R / 10;
    R = R - CARRY*10;
    V[i] = R;
  }
  return CARRY; // may be from 0 to 9
}
uint8_t ConstMult(uint8_t *V, size_t N, uint8_t digit)
{
  uint8_t CARRY = 0;
  uint8_t ja = 0;
  for (size_t i = 0; i < N; ++i) {
    uint8_t aux = V[i] * digit;
    uint8_t R = aux + CARRY;
    CARRY = ((u_int32_t)R*(u_int32_t)0xCCCD) >> 19;
    ja = (CARRY << 3) + 2*CARRY;
    R -= ja;
    V[i] = R;
  }
  return CARRY;
}

这是一个函数,它使用辅助函数 table:

一次处理 2 个字节的块而不除法
uint8_t ConstMult3(uint8_t *V, size_t N, uint8_t digit) {
#define TABLE_SIZE  ((9 * 256 + 9) * 9 + 9 + 1)
    static uint32_t table[TABLE_SIZE];
    if (!table[1]) {
        for (uint32_t x = 0; x < TABLE_SIZE; x++) {
            uint32_t u = x % 256 % 10;
            uint32_t d = (x / 256 + x % 256 / 10) % 10;
            uint32_t c = (x / 256 + x % 256 / 10) / 10;
            //table[x] = u | (d << 8) | (c << 16);
            // modified following Jerome Richard's comment
            table[x] = c | (u << 8) | (d << 16);
        }
    }
    if (N == 0 || digit <= 1) {
        if (digit == 0)
            memset(V, 0, N);
        return 0;
    } else {
        size_t CARRY = 0;

        if ((uintptr_t)V & 1) {  // V is misaligned
            int R = V[0] * digit + (uint8_t)CARRY;
            CARRY = (uint8_t)(R / 10);
            V[0] = (uint8_t)(R - CARRY * 10);
            V++;
            N--;
        }
        {   // handle aligned block 2 bytes at a time
            uint16_t *V2 = (uint16_t *)(void *)V;
            size_t N2 = N / 2;
            for (size_t i = 0; i < N2; i++) {
                uint32_t x = table[V2[i] * digit + CARRY];
                //V2[i] = (uint16_t)x;
                //CARRY = x >> 16;
                // modified following Jerome Richard's comment
                V2[i] = (uint16_t)(x >> 8);
                CARRY = (uint8_t)x;
            }
        }
        if (N & 1) {    // handle last byte
            int R = V[N - 1] * digit + (uint8_t)CARRY;
            CARRY = (uint8_t)(R / 10);
            V[N - 1] = (uint8_t)(R - CARRY * 10);
        }
        return (uint8_t)CARRY;
    }
#undef TABLE_SIZE
}

在我的慢速笔记本电脑上,在 64 位模式下使用 clang 9.0,我用 ConstMult0ConstMult1ConstMult2 得到这些时间是问题中发布的函数:

ConstMult0(1000000): 15.816ms sum0=4495507, sum=4501418
ConstMult1(1000000): 16.464ms sum0=4495507, sum=4501418
ConstMult2(1000000): 16.483ms sum0=4495507, sum=4501418
ConstMult3(1000000): 9.644ms sum0=4495507, sum=4501418

编辑: 根据 Jérôme Richard 的评论,table 内容中的一个小改动可以额外提高 11% 的性能:

ConstMult0(1000000): 15.837ms sum0=4500384, sum=4495487
ConstMult1(1000000): 16.494ms sum0=4500384, sum=4495487
ConstMult2(1000000): 16.482ms sum0=4500384, sum=4495487
ConstMult3(1000000): 8.537ms sum0=4500384, sum=4495487

这是另一个实现(比其他实现快得多):

void ConstMult4(uint8_t *V, size_t N, uint8_t digit)
{
    uint8_t CARRY = 0;

    const uint32_t coef7  = digit * 10000000;
    const uint32_t coef6  = digit * 1000000;
    const uint32_t coef5  = digit * 100000;
    const uint32_t coef4  = digit * 10000;
    const uint32_t coef3  = digit * 1000;
    const uint32_t coef2  = digit * 100;
    const uint32_t coef1  = digit * 10;
    const uint32_t coef0  = digit;

    static uint8_t table[10000][4];
    static int init = 1;

    if(init)
    {
        for(int i=0 ; i<10000 ; ++i)
        {
            table[i][0] = (i / 1) % 10;
            table[i][1] = (i / 10) % 10;
            table[i][2] = (i / 100) % 10;
            table[i][3] = (i / 1000) % 10;
        }

        init = 0;
    }

    for(size_t i=0 ; i<N/8*8 ; i+=8)
    {
        const uint32_t val = V[i+7]*coef7 + V[i+6]*coef6 + V[i+5]*coef5 + V[i+4]*coef4 + V[i+3]*coef3 + V[i+2]*coef2 + V[i+1]*coef1 + V[i+0]*coef0 + CARRY;

        CARRY = val / 100000000;

        const uint32_t loVal = val % 10000;
        const uint32_t hiVal = val / 10000 - CARRY * 10000;
        const uint8_t* loTablePtr = &table[loVal][0];
        const uint8_t* hiTablePtr = &table[hiVal][0];

        // Assume the compiler optimize the 2 following calls
        // (otherwise the performance could be quite bad).
        // memcpy is used to prevent performance issue due to pointer aliasing. 
        memcpy(V+i, loTablePtr, 4);
        memcpy(V+i+4, hiTablePtr, 4);
    }

    for(size_t i=N/8*8 ; i<N ; ++i)
    {
        V[i] = V[i] * digit + CARRY;
        CARRY = V[i] / 10;
        V[i] -= CARRY * 10;
    }
}

此实现假设 Vdigit 中的计算数字实际上是数字。 它比其他方法快得多:

  • 按照@phuclv 的建议在内部使用更大的基础(它减少了关键路径并引入了更多的并行性);
  • 使用@chqrlieforyellowblockquotes 提出的查找 table(它可以非常快速地计算 division/modulus 操作)。

甚至可以使用 SSE 4.1 内在函数(SIMD 指令)改进此代码。但是以更少的 portable 代码为代价(尽管它可以在大多数现代 x86_64-based 处理器上运行)。这是实现:

void ConstMult5(uint8_t *V, size_t N, uint8_t digit)
{
    uint8_t CARRY = 0;

    static uint8_t table[10000][4];
    static int init = 1;

    if(init)
    {
        for(int i=0 ; i<10000 ; ++i)
        {
            table[i][0] = (i / 1) % 10;
            table[i][1] = (i / 10) % 10;
            table[i][2] = (i / 100) % 10;
            table[i][3] = (i / 1000) % 10;
        }

        init = 0;
    }

    __m128i coefs1 = _mm_set_epi16(1000, 100, 10, 1, 1000, 100, 10, 1);
    __m128i coefs2 = _mm_set_epi32(10000*digit, 10000*digit, digit, digit);

    for(size_t i=0 ; i<N/16*16 ; i+=8)
    {
        // Require SSE 4.1 (thus smmintrin.h need to be included)
        const __m128i vBlock = _mm_loadu_si128((const __m128i*)&V[i]); // load 16 x uint8_t values (only half is used)
        const __m128i v = _mm_cvtepu8_epi16(vBlock); // Convert the block to 8 x int16_t values
        const __m128i tmp1 = _mm_madd_epi16(v, coefs1); // Compute the sum of adjacent pairs of v * coefs1 and put this in 4 x int32_t values
        const __m128i tmp2 = _mm_add_epi32(tmp1, _mm_shuffle_epi32(tmp1, 0b10110001)); // Horizontal partial sum of 4 x int32_t values
        const __m128i tmp3 = _mm_mul_epu32(tmp2, coefs2); // Compute tmp2 * coefs2 and put this in 2 x int64_t values
        const uint32_t val = _mm_extract_epi64(tmp3, 1) + _mm_extract_epi64(tmp3, 0) + CARRY; // Final horizontal sum with CARRY

        CARRY = val / 100000000;

        const uint32_t loVal = val % 10000;
        const uint32_t hiVal = val / 10000 - CARRY * 10000;
        const uint8_t* loTablePtr = &table[loVal][0];
        const uint8_t* hiTablePtr = &table[hiVal][0];

        // See the memcpy remark in the code above (alternative version).
        memcpy(V+i, loTablePtr, 4);
        memcpy(V+i+4, hiTablePtr, 4);
    }

    for(size_t i=N/16*16 ; i<N ; ++i)
    {
        V[i] = V[i] * digit + CARRY;
        CARRY = V[i] / 10;
        V[i] -= CARRY * 10;
    }
}

以下是在我的机器(配备 i7-9700KF 处理器)上的性能结果(使用随机输入在 1000 运行 上重复并取平均值):

ConstMult0(10000): 11.702 us
ConstMult3(10000): 6.768 us (last optimized version)
ConstMult4(10000): 3.569 us
ConstMult5(10000): 2.552 us

基于 SSE 的最终版本比您的原始实现快 4.6 倍!