Ускорение элементарного умножения массива в python
Я играл с numba и numexpr, пытаясь ускорить простое умножение матричных элементов. Я не смог получить лучшие результаты, они оба в основном (по скорости) эквивалентны функции numpys multiply. Кому-нибудь повезло в этой области? Я использую numba и numexpr неправильно (я совершенно новичок в этом), или это совсем плохой подход, чтобы попытаться ускорить это. Вот воспроизводимый код, спасибо вам за продвинутый:
import numpy as np
from numba import autojit
import numexpr as ne
a=np.random.rand(10,5000000)
# numpy
multiplication1 = np.multiply(a,a)
# numba
def multiplix(X,Y):
M = X.shape[0]
N = X.shape[1]
D = np.empty((M, N), dtype=np.float)
for i in range(M):
for j in range(N):
D[i,j] = X[i, j] * Y[i, j]
return D
mul = autojit(multiplix)
multiplication2 = mul(a,a)
# numexpr
def numexprmult(X,Y):
M = X.shape[0]
N = X.shape[1]
return ne.evaluate("X * Y")
multiplication3 = numexprmult(a,a)
Ответы
Ответ 1
Как насчет использования fortran и ctypes?
elementwise.F90:
subroutine elementwise( a, b, c, M, N ) bind(c, name='elementwise')
use iso_c_binding, only: c_float, c_int
integer(c_int),intent(in) :: M, N
real(c_float), intent(in) :: a(M, N), b(M, N)
real(c_float), intent(out):: c(M, N)
integer :: i,j
forall (i=1:M,j=1:N)
c(i,j) = a(i,j) * b(i,j)
end forall
end subroutine
elementwise.py:
from ctypes import CDLL, POINTER, c_int, c_float
import numpy as np
import time
fortran = CDLL('./elementwise.so')
fortran.elementwise.argtypes = [ POINTER(c_float),
POINTER(c_float),
POINTER(c_float),
POINTER(c_int),
POINTER(c_int) ]
# Setup
M=10
N=5000000
a = np.empty((M,N), dtype=c_float)
b = np.empty((M,N), dtype=c_float)
c = np.empty((M,N), dtype=c_float)
a[:] = np.random.rand(M,N)
b[:] = np.random.rand(M,N)
# Fortran call
start = time.time()
fortran.elementwise( a.ctypes.data_as(POINTER(c_float)),
b.ctypes.data_as(POINTER(c_float)),
c.ctypes.data_as(POINTER(c_float)),
c_int(M), c_int(N) )
stop = time.time()
print 'Fortran took ',stop - start,'seconds'
# Numpy
start = time.time()
c = np.multiply(a,b)
stop = time.time()
print 'Numpy took ',stop - start,'seconds'
Я скомпилировал файл Fortran с помощью
gfortran -O3 -funroll-loops -ffast-math -floop-strip-mine -shared -fPIC \
-o elementwise.so elementwise.F90
Выход дает ускорение ~ 10%:
$ python elementwise.py
Fortran took 0.213667869568 seconds
Numpy took 0.230120897293 seconds
$ python elementwise.py
Fortran took 0.209784984589 seconds
Numpy took 0.231616973877 seconds
$ python elementwise.py
Fortran took 0.214708089828 seconds
Numpy took 0.25369310379 seconds
Ответ 2
Как вы делаете свои тайминги?
Создание вашего случайного массива занимает верхнюю часть вашего расчета, и если вы включите его в свое время, вы вряд ли увидите какую-либо реальную разницу в результатах,
однако, если вы создадите его спереди, вы можете фактически сравнить методы.
Вот мои результаты, и я постоянно вижу, что вы видите. numpy и numba дают примерно одинаковые результаты (с numba немного быстрее).
(у меня нет numexpr)
In [1]: import numpy as np
In [2]: from numba import autojit
In [3]: a=np.random.rand(10,5000000)
In [4]: %timeit multiplication1 = np.multiply(a,a)
10 loops, best of 3: 90 ms per loop
In [5]: # numba
In [6]: def multiplix(X,Y):
...: M = X.shape[0]
...: N = X.shape[1]
...: D = np.empty((M, N), dtype=np.float)
...: for i in range(M):
...: for j in range(N):
...: D[i,j] = X[i, j] * Y[i, j]
...: return D
...:
In [7]: mul = autojit(multiplix)
In [26]: %timeit multiplication1 = np.multiply(a,a)
10 loops, best of 3: 182 ms per loop
In [27]: %timeit multiplication1 = np.multiply(a,a)
10 loops, best of 3: 185 ms per loop
In [28]: %timeit multiplication1 = np.multiply(a,a)
10 loops, best of 3: 181 ms per loop
In [29]: %timeit multiplication2 = mul(a,a)
10 loops, best of 3: 179 ms per loop
In [30]: %timeit multiplication2 = mul(a,a)
10 loops, best of 3: 180 ms per loop
In [31]: %timeit multiplication2 = mul(a,a)
10 loops, best of 3: 178 ms per loop
Обновление:
Я использовал последнюю версию numba, просто скомпилировал ее из источника: '0.11.0-3-gea20d11-dirty'
Я тестировал это с по умолчанию numpy в Fedora 19, '1.7.1'
и numpy '1.6.1', скомпилированный из источника, связанный с:
Update3
Мои ранние результаты были, конечно, неправильными, я вернул D во внутренний цикл, так что пропустил 90% вычислений.
Это дает больше доказательств предположения ali_m о том, что действительно трудно сделать лучше, чем уже оптимизированный код c.
Однако, если вы пытаетесь сделать что-то более сложное, например,
np.sqrt(((X[:, None, :] - X) ** 2).sum(-1))
Я могу воспроизвести фигуры Джейка Вандерпласа get's:
In [14]: %timeit pairwise_numba(X)
10000 loops, best of 3: 92.6 us per loop
In [15]: %timeit pairwise_numpy(X)
1000 loops, best of 3: 662 us per loop
Итак, кажется, что вы делаете то, что до сих пор оптимизировано numpy, трудно сделать лучше.
Ответ 3
Изменить: не верьте этому ответу, я ошибаюсь (см. комментарий ниже).
Я боюсь, будет очень и очень сложно иметь более быстрое умножение матрицы в python, чем при использовании numpy. NumPy обычно использует внутренние библиотеки fortran, такие как ATLAS/LAPACK, которые очень хорошо оптимизированы.
Чтобы проверить, была ли построена ваша версия NumPy с поддержкой LAPACK: откройте терминал, зайдите в каталог установки Python и введите:
for f in `find lib/python2.7/site-packages/numpy/* -name \*.so`; do echo $f; ldd $f;echo "\n";done | grep lapack
Обратите внимание, что путь может варьироваться в зависимости от вашей версии python.
Если вы напечатаете несколько строк, у вас наверняка будет поддержка LAPACK... так что добиться более быстрого умножения матрицы на одно ядро будет очень сложно.
Теперь я не знаю, как использовать несколько ядер для выполнения умножения на матрицу, поэтому вы можете посмотреть на это (см. комментарий ali_m).
Ответ 4
используйте графический процессор. используйте следующий пакет.
gnumpy
