Numpy dot product очень медленный, используя целые числа

извините за столько вопросов. Я использую 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 раз медленнее, чем с использованием чисел с плавающей запятой?


person Nino    schedule 08.08.2012    source источник
comment
То же самое для меня в Linux - я получаю около 3 секунд для float64 и 10 секунд для int32 (это более старая машина). Не в 50 раз, но все равно очень странно.   -  person Luke    schedule 08.08.2012


Ответы (2)


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

>>> 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.

so:

>>> help(np.core._dotblas)

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

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

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

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

person Samy Vilar    schedule 08.08.2012
comment
@Luke спасибо +1 за то, что ты сделал обратную трассировку, это еще один полезный метод. - person Samy Vilar; 08.08.2012
comment
Хорошо знать. Когда я пойду на работу, я буду использовать тот же метод, чтобы проверить, верно ли это и для MKL. Спасибо за помощь. Люку больше нравится этот ответ, поэтому вы принимаете его. - person Nino; 09.08.2012
comment
Мне нравится, когда ответы не только отвечают на вопрос, но и показывают, как был найден ответ, попутно знакомя с техниками. +1. Нужно больше таких ответов здесь. - person naught101; 11.07.2016

Использование типов данных 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, но похоже, что вызов float остается в библиотеке BLAS (вероятно, доступ к какому-то хорошо оптимизированному коду), а вызов int возвращается обратно в PyArray_MatrixProduct2 numpy.

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

Вот простой и недорогой обходной путь:

af = a.astype(float)
np.dot(af, af).astype(int)
person Luke    schedule 08.08.2012
comment
Стоит отметить, что этот обходной путь может привести к ошибкам, если ваши данные имеют очень большие значения, и, вероятно, потребует копирования всей матрицы, поэтому это дорого, если матрица очень большая. - person Danica; 08.08.2012
comment
Спасибо, Люк. Этот обходной путь действительно копирует матрицы и оказывается довольно неприятным (для проблем с памятью), но с точки зрения времени он в тысячи раз быстрее для матриц 1e4x1e4! Чем больше, тем медленнее проверять умножение с помощью целых чисел. @Dougal Это верно только для чисел больше 2 ^ 52 с использованием 64-битного числа с плавающей запятой, верно? Числа не будут больше этого, и я хотел бы воспользоваться этим ускорением, если это возможно. - person Nino; 09.08.2012
comment
@ Нино Да, где-то там. Очень жаль, что библиотеки BLAS не работают с целочисленными типами, а встроенный в numpy dot намного медленнее. Если проблемы с памятью являются проблемой, вы можете написать небольшой интерфейс ctypes, который выполняет умножение в Eigen или подобном , что должно быть быстрее. - person Danica; 09.08.2012
comment
Фундаментальная матрица-матрица-умножение в BLAS существует в 4 версиях: sgemm.f, dgemm.f, cgemm.f, fzgemm.f для одинарного, двойного, комплексного, двойного комплексного соответственно. Линейная алгебра, которая включает в себя обратные матрицы и решения линейных уравнений, малопригодна для чисто целочисленных вычислений. netlib.org/blas - person hpaulj; 26.08.2013