2

我有一個整數,這是256和131072 之間256多,我想用一個整數分才1和1024生成一個查找表由10位整數通過乘法來劃分

之間這是一個我的代碼內部循環中的熱點和評論它顯着加速我的應用程序。

我可以製作一個1024位的查找表,它可以在比x86_64 cpu上的實際分割更少的時間內將分割轉換爲「乘加加移位」嗎?

有人可以幫助想出代碼來生成查找表,允許由這1024個可能除數之一進行有效除法嗎?

我很想看到一個模板元編程的方式來生成相關表作爲constexpr。

+0

http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.1.2556這似乎與我的問題有關 – hellcatv

回答

1

由於股息和除數的選擇非常有限,因此可以使用比您找到的TorbjörnGranlund和Peter Montgomery的paper更簡單的方法。我不熟悉C++模板元編程,但我可以演示生成和使用適當縮放倒數表的方法。

首先我們注意到股利是256的倍數,所以它們可以通過簡單的預先移位到8位右邊來減少到1 ... 0x200。由於我們不希望在減小的被除數與比例倒數相乘期間溢出無符號的32位整數,所以倒數理想地被縮放到範圍0x00200000 < rcp <= 0x00400000

如果快速計數前導零指令可用,則可用於在表預計算期間根據除數的對數基數2將倒數倒數計算到此範圍內,然後在運行時時間使用相同的動態計算的比例因子來縮減下降減少的股利和縮放倒數乘以相同因子的乘積。當按比例放大時,我們需要將向上調整結果到下一個整數以補償經右移的縮小的截斷性質。下面代碼中的變體0使用這種方法。

當沒有快速計數前導零指令可用時,我們該怎麼辦?我們需要存儲具有足夠位數的縮放倒數,以保持計算的準確性。事實證明,由於對除數範圍的嚴格限制,我們很幸運,並且可以在運行時只使用兩個不同的比例因子,這些比例因子可以很容易地從除數計算出來:除數< = 32,另一個因子(32,1024)中的除數。這用於代碼的變體1中,其中兩個比例因子計算爲2 和2 。最後,我們可能不想動態地計算比例因子,而是通過使用每個表項的最高有效位來存儲它們,而將較低有效位用於互惠本身。一個缺點是需要額外的操作來從表格條目中提取縮放比例和比例因子,這使得這種方法比其他兩種方法更不適合。這在下面的代碼變體2中顯示。

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

#define VARIANT 0 // variants available: 0 ... 2 

extern int __clz (uint32_t x); // intrinsic for access to CLZ instruction 

uint32_t rcp [1025]; // table of scaled reciprocals (& possibly scale factors) 

#define PRE_SCALE  (8) // downscaling for dividends, since multiples of 256 
#define POST_SCALE_1 (14) // upscale factor #1 for reciprocals 
#define POST_SCALE_2 (19) // upscale factor #2 for reciprocals 
#define RECIP_BITS (24) // bits used in table entry for scaled reciprocal 

// helper function: logarithm base-2 of an 32-bit integer  
int ilog2 (uint32_t x) 
{ 
    return 31 - __clz (x); 
} 

// division for dividends n*256, n in [1,512], divisor in [1,1024]  
uint32_t fast_div (uint32_t dividend, uint32_t divisor) 
{ 
#if VARIANT == 0 
    uint32_t scale = POST_SCALE_1 + ilog2 (divisor); 
    return ((dividend >> PRE_SCALE) * rcp[divisor]) >> scale; 
#elif VARIANT == 1 
    uint32_t scale = (divisor > 0x20) ? POST_SCALE_2 : POST_SCALE_1; 
    return ((dividend >> PRE_SCALE) * rcp[divisor]) >> scale; 
#elif VARIANT == 2 
    uint32_t scale = rcp[divisor] >> RECIP_BITS; 
    return ((dividend >> PRE_SCALE) * (rcp[divisor] & ((1 << RECIP_BITS) - 1))) >> scale; 
#else 
#error non-existing VARIANT 
#endif 
} 

int main (void) 
{ 
    uint32_t dividend, divisor, res, ref; 
    int i; 

    // precompute table of recprocals 
    for (i = 1; i < 1025; i++) { 
#if VARIANT == 0 
     uint32_t scale = POST_SCALE_1 + ilog2 (i); 
     rcp[i] = ((uint32_t)(pow (2.0, PRE_SCALE + scale)/i + 0.99999)); 
#elif VARIANT == 1 
     uint32_t scale = (i > 0x20) ? POST_SCALE_2 : POST_SCALE_1; 
     rcp[i] = ((uint32_t)(pow (2.0, PRE_SCALE + scale)/i + 0.99999)); 
#elif VARIANT == 2 
     uint32_t scale = (i > 0x20) ? POST_SCALE_2 : POST_SCALE_1; 
     rcp[i] = ((uint32_t)(pow (2.0, PRE_SCALE + scale)/i + 0.99999) + 
        (scale << RECIP_BITS)); 
#else 
#error non-existing VARIANT 
#endif 
    } 

    // test all supported dividens and divisors exhaustively 
    divisor = 1; 
    while (divisor <= 1024) { 
     dividend = 256; 
     while (dividend <= 131072) { 
      res = fast_div (dividend, divisor); 
      ref = dividend/divisor; 
      if (res != ref) { 
       printf ("n=%08x d=%08x res=%08x ref=%08x rcp=%08x\n", 
         dividend, divisor, res, ref, rcp[divisor]); 
       return EXIT_FAILURE; 
      } 
      dividend += 256; 
     } 
     divisor++; 
    } 
    printf ("division test passed\n"); 
    return EXIT_SUCCESS; 
}