Я использую cython в первый раз, чтобы получить некоторую скорость для функции. Функция принимает квадратную матрицу A
чисел с плавающей запятой и выводит одно число с плавающей запятой. Функция его вычисления является permanent of a matrixМожно ли оптимизировать этот код cython?
Когда А 30 по 30 мой код занимает около 60 секунд, в настоящее время на моем компьютере.
В приведенном ниже коде реализована формула Balasubramanian-Bax/Franklin-Glynn для постоянной страницы wiki. Я назвал матрицу M.
Одна сложная часть кода - это массив f, который используется для удерживания указателя следующей позиции для перевода в массив d. Массив d содержит значения, равные + -1. Манипуляция f и j в цикле - это просто умный способ быстро обновить код серого.
from __future__ import division
import numpy as np
cimport numpy as np
cimport cython
DTYPE_int = np.int
ctypedef np.int_t DTYPE_int_t
DTYPE_float = np.float64
ctypedef np.float64_t DTYPE_float_t
@cython.boundscheck(False) # turn off bounds-checking for entire function
@cython.wraparound(False) # turn off negative index wrapping for entire function
def permfunc(np.ndarray [DTYPE_float_t, ndim =2, mode='c'] M):
cdef int n = M.shape[0]
cdef np.ndarray[DTYPE_float_t, ndim =1, mode='c' ] d = np.ones(n, dtype=DTYPE_float)
cdef int j = 0
cdef int s = 1
cdef np.ndarray [DTYPE_int_t, ndim =1, mode='c'] f = np.arange(n, dtype=DTYPE_int)
cdef np.ndarray [DTYPE_float_t, ndim =1, mode='c'] v = M.sum(axis=0)
cdef DTYPE_float_t p = 1
cdef int i
cdef DTYPE_float_t prod
for i in range(n):
p *= v[i]
while (j < n-1):
for i in range(n):
v[i] -= 2*d[j]*M[j, i]
d[j] = -d[j]
s = -s
prod = 1
for i in range(n):
prod *= v[i]
p += s*prod
f[0] = 0
f[j] = f[j+1]
f[j+1] = j+1
j = f[0]
return p/2**(n-1)
Я использовал все простые оптимизации, которые я нашел в учебнике cython. Некоторые аспекты, которые я должен признать, я не совсем понимаю. Например, если я делаю массив d
ints, так как значения всегда + -1, код работает примерно на 10% медленнее, поэтому я оставил его как float64.
Есть ли что-нибудь еще, что я могу сделать для ускорения кода?
Это является результатом Cython -a. Как вы можете видеть, все в цикле скомпилировано на C, поэтому основные оптимизации уже сработали.
Вот та же функция в NumPy, которая составляет более 100 раз медленнее, чем моя текущая версия Cython.
def npperm(M):
n = M.shape[0]
d = np.ones(n)
j = 0
s = 1
f = np.arange(n)
v = M.sum(axis=0)
p = np.prod(v)
while (j < n-1):
v -= 2*d[j]*M[j]
d[j] = -d[j]
s = -s
prod = np.prod(v)
p += s*prod
f[0] = 0
f[j] = f[j+1]
f[j+1] = j+1
j = f[0]
return p/2**(n-1)
Timings обновленный
Вот тайминги (с использованием IPython) моей версии Cython, версия Numpy и улучшение romeric к коду Cython. Я установил семя для воспроизводимости.
from scipy.stats import ortho_group
import pyximport; pyximport.install()
import permlib # This loads in the functions from permlib.pyx
import numpy as np; np.random.seed(7)
M = ortho_group.rvs(23) #Creates a random orthogonal matrix
%timeit permlib.npperm(M) # The numpy version
1 loop, best of 3: 44.5 s per loop
%timeit permlib.permfunc(M) # The cython version
1 loop, best of 3: 273 ms per loop
%timeit permlib.permfunc_modified(M) #romeric's improvement
10 loops, best of 3: 198 ms per loop
M = ortho_group.rvs(28)
%timeit permlib.permfunc(M) # The cython version run on a 28x28 matrix
1 loop, best of 3: 15.8 s per loop
%timeit permlib.permfunc_modified(M) # romeric's improvement run on a 28x28 matrix
1 loop, best of 3: 12.4 s per loop
Может код Cython быть ускорен на всех?
Я использую GCC и CPU является FX AMD 8350.
Да вы могли бы спросить это на [codereview.se]. – usr2564301
@RadLexus Спасибо. Однако, похоже, вопросы цитона встречаются редко. Было 30 когда-либо! – eleanora
@eleanora: с такими рассуждениями, что это число остается низким. – Mat