Многоточечный точечный продукт очень медленный, используя ints

Извините за столько вопросов. Я запускаю Mac OSX 10.6 на Intel Core 2 Duo. Я использую некоторые тесты для своих исследований, и я столкнулся с другой проблемой, которая меня озадачивает.

Если я запустил

python -mtimeit -s 'import numpy as np; a = np.random.randn(1e3,1e3)' 'np.dot(a,a)'

Я получаю следующий вывод: 10 loops, best of 3: 142 msec per loop

Однако, если я запустил

python -mtimeit -s 'import numpy as np; a = np.random.randint(10,size=1e6).reshape(1e3,1e3)' 'np.dot(a,a)'

Я получаю следующий вывод: 10 loops, best of 3: 7.57 sec per loop

Затем я запустил

python -mtimeit -s 'import numpy as np; a = np.random.randn(1e3,1e3)' 'a*a' И затем

python -mtimeit -s 'import numpy as np; a = np.random.randint(10,size=1e6).reshape(1e3,1e3)' 'a*a'

Оба выполнялись примерно в 7,6 мс за цикл, поэтому это не умножение. Добавление тоже схожих скоростей, поэтому ни одно из них не должно влиять на точечный продукт, верно? Итак, почему он более чем в 50 раз медленнее вычисляет точечный продукт с помощью ints, чем с помощью float?

Ответ 1

Очень интересно, мне было любопытно посмотреть, как это было реализовано, поэтому я сделал:

>>> import inspect
>>> import numpy as np
>>> inspect.getmodule(np.dot)
<module 'numpy.core._dotblas' from '/Library/Python/2.6/site-packages/numpy-1.6.1-py2.6-macosx-10.6-universal.egg/numpy/core/_dotblas.so'>
>>>

Итак, это похоже на использование библиотеки BLAS.

так:

>>> help(np.core._dotblas)

из которого я нашел это:

Когда Numpy построен с ускоренным BLAS, таким как ATLAS, эти функции заменяются, чтобы использовать более быстрые реализации. Быстрее реализации влияют только на float32, float64, complex64 и complex128 массивы. Кроме того, API BLAS включает только матричную матрицу, матрично-векторных и вектор-векторных произведений. Продукты массивов с большими мерности используют встроенные функции и не ускоряются.

Таким образом, похоже, что ATLAS точно настраивает определенные функции, но применим только к определенным типам данных, очень интересно.

Так что да, похоже, я буду чаще использовать поплавки...

Ответ 2

Использование типов данных int vs float приводит к выполнению различных путей кода:

Трассировка стека для float выглядит следующим образом:

(gdb) backtr
#0  0x007865a0 in dgemm_ () from /usr/lib/libblas.so.3gf
#1  0x007559d5 in cblas_dgemm () from /usr/lib/libblas.so.3gf
#2  0x00744108 in dotblas_matrixproduct (__NPY_UNUSED_TAGGEDdummy=0x0, args=(<numpy.ndarray at remote 0x85d9090>, <numpy.ndarray at remote 0x85d9090>), 
kwargs=0x0) at numpy/core/blasdot/_dotblas.c:798
#3  0x08088ba1 in PyEval_EvalFrameEx ()
...

.. в то время как трассировка стека для int выглядит так:

(gdb) backtr
#0  LONG_dot (ip1=0xb700a280 "\t", is1=4, ip2=0xb737dc64 "\a", is2=4000, op=0xb6496fc4 "", n=1000, __NPY_UNUSED_TAGGEDignore=0x85fa960)
at numpy/core/src/multiarray/arraytypes.c.src:3076
#1  0x00659d9d in PyArray_MatrixProduct2 (op1=<numpy.ndarray at remote 0x85dd628>, op2=<numpy.ndarray at remote 0x85dd628>, out=0x0)
at numpy/core/src/multiarray/multiarraymodule.c:847
#2  0x00742b93 in dotblas_matrixproduct (__NPY_UNUSED_TAGGEDdummy=0x0, args=(<numpy.ndarray at remote 0x85dd628>, <numpy.ndarray at remote 0x85dd628>), 
kwargs=0x0) at numpy/core/blasdot/_dotblas.c:254
#3  0x08088ba1 in PyEval_EvalFrameEx ()
...

Оба вызова приводят к методу dotblas_matrixproduct, но похоже, что поплавочный вызов остается в библиотеке BLAS (возможно, доступ к некоторому хорошо оптимизированному коду), в то время как int-вызов отбрасывается обратно в numpy PyArray_MatrixProduct2.

Так что это либо ошибка, либо BLAS просто не поддерживает целые типы в matrixproduct (что кажется маловероятным).

Здесь легко и недорого обходное решение:

af = a.astype(float)
np.dot(af, af).astype(int)