为什么有些时候 Python 中乘法比位运算更快

开发 后端
某天,一个技术群里老哥提出了这样一个问题,为什么在一些情况下,Python 中的简单乘/除法比位运算要慢。

我本来以为我不再会写水文了,但是突然发现自己现在也只能勉强写写水文才能维持生活这样子。那就继续写水文吧!

某天,一个技术群里老哥提出了这样一个问题,为什么在一些情况下,Python 中的简单乘/除法比位运算要慢。

首先秉持着实事求是的精神,我们先来验证一下:

  1. In [33]: %timeit 1073741825*2                                                                                                                                                                                                                                                                            
  2. 7.47 ns ± 0.0843 ns per loop (mean ± std. dev. of 7 runs, 100000000 loops each) 
  3.  
  4. In [34]: %timeit 1073741825<<1                                                                                                                                                                                                                                                                           
  5. 7.43 ns ± 0.0451 ns per loop (mean ± std. dev. of 7 runs, 100000000 loops each) 
  6.  
  7. In [35]: %timeit 1073741823<<1                                                                                                                                                                                                                                                                           
  8. 7.48 ns ± 0.0621 ns per loop (mean ± std. dev. of 7 runs, 100000000 loops each) 
  9.  
  10. In [37]: %timeit 1073741823*2                                                                                                                                                                                                                                                                            
  11. 7.47 ns ± 0.0564 ns per loop (mean ± std. dev. of 7 runs, 100000000 loops each) 

我们发现几个很有趣的现象:

  • 在值 x<=2^30 时,乘法比直接位运算要快
  • 在值 x>2^32 时,乘法显著慢于位运算

这个现象很有趣,那么这个现象的 root cause 是什么?实际上这和 Python 底层的实现有关。

简单聊聊

1. PyLongObject 的实现

在 Python 2.x 时期,Python 中将整型分为两类,一类是 long, 一类是 int 。在 Python3 中这两者进行了合并。目前在 Python3 中这两者做了合并,仅剩一个 long。

首先来看看 long 这样一个数据结构底层的实现:

  1. struct _longobject { 
  2.     PyObject_VAR_HEAD 
  3.     digit ob_digit[1]; 
  4. }; 

在这里不用关心,PyObject_VAR_HEAD 的含义,我们只需要关心 ob_digit 即可。

在这里,ob_digit 是使用了 C99 中的“柔性数组”来实现任意长度的整数的存储。这里我们可以看一下官方代码中的文档:

Long integer representation.The absolute value of a number is equal to SUM(for i=0 through abs(ob_size)-1) ob_digit[i] * 2**(SHIFT*i) Negative numbers are represented with ob_size < 0; zero is represented by ob_size == 0. In a normalized number, ob_digit[abs(ob_size)-1] (the most significant digit) is never zero. Also, in all cases, for all valid i,0 <= ob_digit[i] <= MASK. The allocation function takes care of allocating extra memory so that ob_digit[0] ... ob_digit[abs(ob_size)-1] are actually available. CAUTION: Generic code manipulating subtypes of PyVarObject has to aware that ints abuse ob_size's sign bit.

简而言之,Python 是将一个十进制数转为 2^(SHIFT) 进制数来进行存储。这里可能不太好了理解。我来举个例子,在我的电脑上,SHIFT 为 30 ,假设现在有整数 1152921506754330628 ,那么将其转为 2^30 进制表示则为: 4*(2^30)^0+2*(2^30)^1+1*(2^30)^2 。那么此时 ob_digit 是一个含有三个元素的数组,其值为 [4,2,1]。

OK,在明白了这样一些基础知识后,我们回过头去看看 Python 中的乘法运算。

2. Python 中的乘法运算

Python 中的乘法运算,分为两部分,其中关于大数的乘法,Python 使用了 Karatsuba 算法1,具体实现如下:

  1. static PyLongObject * 
  2. k_mul(PyLongObject *a, PyLongObject *b) 
  3.     Py_ssize_t asize = Py_ABS(Py_SIZE(a)); 
  4.     Py_ssize_t bsize = Py_ABS(Py_SIZE(b)); 
  5.     PyLongObject *ah = NULL
  6.     PyLongObject *al = NULL
  7.     PyLongObject *bh = NULL
  8.     PyLongObject *bl = NULL
  9.     PyLongObject *ret = NULL
  10.     PyLongObject *t1, *t2, *t3; 
  11.     Py_ssize_t shift;           /* the number of digits we split off */ 
  12.     Py_ssize_t i; 
  13.  
  14.     /* (ah*X+al)(bh*X+bl) = ah*bh*X*X + (ah*bl + al*bh)*X + al*bl 
  15.      * Let k = (ah+al)*(bh+bl) = ah*bl + al*bh  + ah*bh + al*bl 
  16.      * Then the original product is 
  17.      *     ah*bh*X*X + (k - ah*bh - al*bl)*X + al*bl 
  18.      * By picking X to be a power of 2, "*X" is just shifting, and it's 
  19.      * been reduced to 3 multiplies on numbers half the size. 
  20.      */ 
  21.  
  22.     /* We want to split based on the larger number; fiddle so that b 
  23.      * is largest. 
  24.      */ 
  25.     if (asize > bsize) { 
  26.         t1 = a
  27.         a = b
  28.         b = t1
  29.  
  30.         i = asize
  31.         asize = bsize
  32.         bsize = i; 
  33.     } 
  34.  
  35.     /* Use gradeschool math when either number is too small. */ 
  36.     i = a == b ? KARATSUBA_SQUARE_CUTOFF : KARATSUBA_CUTOFF; 
  37.     if (asize <= i) { 
  38.         if (asize == 0) 
  39.             return (PyLongObject *)PyLong_FromLong(0); 
  40.         else 
  41.             return x_mul(a, b); 
  42.     } 
  43.  
  44.     /* If a is small compared to b, splitting on b gives a degenerate 
  45.      * case with ah==0, and Karatsuba may be (even much) less efficient 
  46.      * than "grade school" then.  However, we can still win, by viewing 
  47.      * b as a string of "big digits", each of width a->ob_size.  That 
  48.      * leads to a sequence of balanced calls to k_mul. 
  49.      */ 
  50.     if (2 * asize <= bsize) 
  51.         return k_lopsided_mul(a, b); 
  52.  
  53.     /* Split a & b into hi & lo pieces. */ 
  54.     shift = bsize >> 1; 
  55.     if (kmul_split(a, shift, &ah, &al) < 0) goto fail; 
  56.     assert(Py_SIZE(ah) > 0);            /* the split isn't degenerate */ 
  57.  
  58.     if (a == b) { 
  59.         bh = ah
  60.         bl = al
  61.         Py_INCREF(bh); 
  62.         Py_INCREF(bl); 
  63.     } 
  64.     else if (kmul_split(b, shift, &bh, &bl) < 0) goto fail; 
  65.  
  66.     /* The plan: 
  67.      * 1. Allocate result space (asize + bsize digits:  that's always 
  68.      *    enough). 
  69.      * 2. Compute ah*bh, and copy into result at 2*shift. 
  70.      * 3. Compute al*bl, and copy into result at 0.  Note that this 
  71.      *    can't overlap with #2. 
  72.      * 4. Subtract al*bl from the result, starting at shift.  This may 
  73.      *    underflow (borrow out of the high digit), but we don't care: 
  74.      *    we're effectively doing unsigned arithmetic mod 
  75.      *    BASE**(sizea + sizeb), and so long as the *final* result fits, 
  76.      *    borrows and carries out of the high digit can be ignored. 
  77.      * 5. Subtract ah*bh from the result, starting at shift. 
  78.      * 6. Compute (ah+al)*(bh+bl), and add it into the result starting 
  79.      *    at shift. 
  80.      */ 
  81.  
  82.     /* 1. Allocate result space. */ 
  83.     ret = _PyLong_New(asize + bsize); 
  84.     if (ret == NULL) goto fail; 
  85. #ifdef Py_DEBUG 
  86.     /* Fill with trash, to catch reference to uninitialized digits. */ 
  87.     memset(ret->ob_digit, 0xDF, Py_SIZE(ret) * sizeof(digit)); 
  88. #endif 
  89.  
  90.     /* 2. t1 <- ah*bh, and copy into high digits of result. */ 
  91.     if ((t1 = k_mul(ah, bh)) == NULL) goto fail; 
  92.     assert(Py_SIZE(t1) >= 0); 
  93.     assert(2*shift + Py_SIZE(t1) <= Py_SIZE(ret)); 
  94.     memcpy(ret->ob_digit + 2*shift, t1->ob_digit, 
  95.            Py_SIZE(t1) * sizeof(digit)); 
  96.  
  97.     /* Zero-out the digits higher than the ah*bh copy. */ 
  98.     i = Py_SIZE(ret) - 2*shift - Py_SIZE(t1); 
  99.     if (i) 
  100.         memset(ret->ob_digit + 2*shift + Py_SIZE(t1), 0, 
  101.                i * sizeof(digit)); 
  102.  
  103.     /* 3. t2 <- al*bl, and copy into the low digits. */ 
  104.     if ((t2 = k_mul(al, bl)) == NULL) { 
  105.         Py_DECREF(t1); 
  106.         goto fail; 
  107.     } 
  108.     assert(Py_SIZE(t2) >= 0); 
  109.     assert(Py_SIZE(t2) <= 2*shift); /* no overlap with high digits */ 
  110.     memcpy(ret->ob_digit, t2->ob_digit, Py_SIZE(t2) * sizeof(digit)); 
  111.  
  112.     /* Zero out remaining digits. */ 
  113.     i = 2*shift - Py_SIZE(t2);          /* number of uninitialized digits */ 
  114.     if (i) 
  115.         memset(ret->ob_digit + Py_SIZE(t2), 0, i * sizeof(digit)); 
  116.  
  117.     /* 4 & 5. Subtract ah*bh (t1) and al*bl (t2).  We do al*bl first 
  118.      * because it's fresher in cache. 
  119.      */ 
  120.     i = Py_SIZE(ret) - shift;  /* # digits after shift */ 
  121.     (void)v_isub(ret->ob_digit + shift, i, t2->ob_digit, Py_SIZE(t2)); 
  122.     Py_DECREF(t2); 
  123.  
  124.     (void)v_isub(ret->ob_digit + shift, i, t1->ob_digit, Py_SIZE(t1)); 
  125.     Py_DECREF(t1); 
  126.  
  127.     /* 6. t3 <- (ah+al)(bh+bl), and add into result. */ 
  128.     if ((t1 = x_add(ah, al)) == NULL) goto fail; 
  129.     Py_DECREF(ah); 
  130.     Py_DECREF(al); 
  131.     ah = al = NULL; 
  132.  
  133.     if (a == b) { 
  134.         t2 = t1
  135.         Py_INCREF(t2); 
  136.     } 
  137.     else if ((t2 = x_add(bh, bl)) == NULL) { 
  138.         Py_DECREF(t1); 
  139.         goto fail; 
  140.     } 
  141.     Py_DECREF(bh); 
  142.     Py_DECREF(bl); 
  143.     bh = bl = NULL; 
  144.  
  145.     t3 = k_mul(t1, t2); 
  146.     Py_DECREF(t1); 
  147.     Py_DECREF(t2); 
  148.     if (t3 == NULL) goto fail; 
  149.     assert(Py_SIZE(t3) >= 0); 
  150.  
  151.     /* Add t3.  It's not obvious why we can't run out of room here. 
  152.      * See the (*) comment after this function. 
  153.      */ 
  154.     (void)v_iadd(ret->ob_digit + shift, i, t3->ob_digit, Py_SIZE(t3)); 
  155.     Py_DECREF(t3); 
  156.  
  157.     return long_normalize(ret); 
  158.  
  159.   fail: 
  160.     Py_XDECREF(ret); 
  161.     Py_XDECREF(ah); 
  162.     Py_XDECREF(al); 
  163.     Py_XDECREF(bh); 
  164.     Py_XDECREF(bl); 
  165.     return NULL; 

这里不对 Karatsuba 算法1 的实现做单独解释,有兴趣的朋友可以参考文末的 reference 去了解具体的详情。

在普通情况下,普通乘法的时间复杂度为 n^2 (n 为位数),而 K 算法的时间复杂度为 3n^(log3) ≈ 3n^1.585 ,看起来 K 算法的性能要优于普通乘法,那么为什么 Python 不全部使用 K 算法呢?

很简单,K 算法的优势实际上要在当 n 足够大的时候,才会对普通乘法形成优势。同时考虑到内存访问等因素,当 n 不够大时,实际上采用 K 算法的性能将差于直接进行乘法。

所以我们来看看 Python 中乘法的实现:

  1. static PyObject * 
  2. long_mul(PyLongObject *a, PyLongObject *b) 
  3.     PyLongObject *z; 
  4.  
  5.     CHECK_BINOP(a, b); 
  6.  
  7.     /* fast path for single-digit multiplication */ 
  8.     if (Py_ABS(Py_SIZE(a)) <= 1 && Py_ABS(Py_SIZE(b)) <= 1) { 
  9.         stwodigits v = (stwodigits)(MEDIUM_VALUE(a)) * MEDIUM_VALUE(b); 
  10.         return PyLong_FromLongLong((long long)v); 
  11.     } 
  12.  
  13.     z = k_mul(a, b); 
  14.     /* Negate if exactly one of the inputs is negative. */ 
  15.     if (((Py_SIZE(a) ^ Py_SIZE(b)) < 0) && z) { 
  16.         _PyLong_Negate(&z); 
  17.         if (z == NULL) 
  18.             return NULL; 
  19.     } 
  20.     return (PyObject *)z; 

在这里我们看到,当两个数皆小于 2^30-1 时,Python 将直接使用普通乘法并返回,否则将使用 K 算法进行计算

这个时候,我们来看一下位运算的实现,以右移为例:

  1. static PyObject * 
  2. long_rshift(PyObject *a, PyObject *b) 
  3.     Py_ssize_t wordshift; 
  4.     digit remshift; 
  5.  
  6.     CHECK_BINOP(a, b); 
  7.  
  8.     if (Py_SIZE(b) < 0) { 
  9.         PyErr_SetString(PyExc_ValueError, "negative shift count"); 
  10.         return NULL; 
  11.     } 
  12.     if (Py_SIZE(a) == 0) { 
  13.         return PyLong_FromLong(0); 
  14.     } 
  15.     if (divmod_shift(b, &wordshift, &remshift) < 0
  16.         return NULL; 
  17.     return long_rshift1((PyLongObject *)a, wordshift, remshift); 
  18.  
  19. static PyObject * 
  20. long_rshift1(PyLongObject *a, Py_ssize_t wordshift, digit remshift) 
  21.     PyLongObject *z = NULL
  22.     Py_ssize_t newsize, hishift, i, j; 
  23.     digit lomask, himask; 
  24.  
  25.     if (Py_SIZE(a) < 0) { 
  26.         /* Right shifting negative numbers is harder */ 
  27.         PyLongObject *a1, *a2; 
  28.         a1 = (PyLongObject *) long_invert(a); 
  29.         if (a1 == NULL) 
  30.             return NULL; 
  31.         a2 = (PyLongObject *) long_rshift1(a1, wordshift, remshift); 
  32.         Py_DECREF(a1); 
  33.         if (a2 == NULL) 
  34.             return NULL; 
  35.         z = (PyLongObject *) long_invert(a2); 
  36.         Py_DECREF(a2); 
  37.     } 
  38.     else { 
  39.         newsize = Py_SIZE(a) - wordshift; 
  40.         if (newsize <= 0) 
  41.             return PyLong_FromLong(0); 
  42.         hishift = PyLong_SHIFT - remshift; 
  43.         lomask = ((digit)1 << hishift) - 1; 
  44.         himask = PyLong_MASK ^ lomask; 
  45.         z = _PyLong_New(newsize); 
  46.         if (z == NULL) 
  47.             return NULL; 
  48.         for (i = 0j = wordshift; i < newsize; i++, j++) { 
  49.             z->ob_digit[i] = (a->ob_digit[j] >> remshift) & lomask; 
  50.             if (i+1 < newsize
  51.                 z->ob_digit[i] |= (a->ob_digit[j+1] << hishift) & himask; 
  52.         } 
  53.         z = maybe_small_long(long_normalize(z)); 
  54.     } 
  55.     return (PyObject *)z; 

在这里我们能看到,在两侧都是小数的情况下,位移动算法将比普通乘法,存在更多的内存分配等操作。这样也会回答了我们文初所提到的一个问题,“为什么一些时候乘法比位运算更快”。

总结本文差不多就到这里了,实际上通过这次分析我们能得到一些很有趣但是也很冷门的知识。实际上我们目前看到这样一个结果,是 Python 对于我们常见且高频的操作所做的一个特定的设计。而这也提醒我们,Python 实际上对于很多操作都存在自己内建的设计哲学,在日常使用的时候,其余语言的经验,可能无法复用。

 

责任编辑:赵宁宁 来源: Manjusaka的编程屋
相关推荐

2021-01-13 10:51:08

PromissetTimeout(函数

2023-09-14 15:48:53

排序测试

2023-09-20 00:06:30

Python代码函数

2013-07-24 09:47:52

语言语速环境语言

2012-08-23 09:28:01

编程编程语言

2021-01-30 10:51:07

Python编程语言开发

2021-04-23 11:11:59

加密货币硬币数字货币

2021-07-23 16:30:36

PythonC++代码

2014-08-29 09:56:47

排序数组编程技巧

2024-02-05 22:51:49

AGIRustPython

2020-09-15 09:55:30

类比Python开发

2015-08-06 12:50:47

技术人员博客

2022-11-10 15:32:29

2020-07-17 19:31:19

PythonR编程

2020-02-14 13:53:33

Python 开发编程语言

2019-09-16 12:00:03

constC编程语言

2015-07-31 16:29:15

DockerJavaLinux

2021-12-27 07:10:26

ClassmethodStaticmetho函数

2019-04-24 08:00:00

HTTPSHTTP前端

2016-12-14 12:02:01

StormHadoop大数据
点赞
收藏

51CTO技术栈公众号