Медленное деление на cython

Ответ 1

Во-первых, вам нужно вызывать функции много ( > 1000) раз и брать среднее время, затраченное в каждом, чтобы получить точное представление о том, насколько они отличаются. Вызов каждой функции один раз не будет достаточно точным.

Во-вторых, время, затрачиваемое на функцию, будет зависеть от других вещей, а не только от цикла с делениями. Вызов def то есть функция Python, подобная этой, включает некоторые накладные расходы при передаче и возврате аргументов. Кроме того, создание массива numpy в функции потребует времени, поэтому любые различия в циклах в двух функциях будут менее очевидными.

Наконец, см. здесь (https://github.com/cython/cython/wiki/enhancements-compilerdirectives), установка директивы c-division False имеет штраф в размере 35%. Я думаю, этого недостаточно, чтобы появиться в вашем примере, учитывая другие накладные расходы. Я проверил вывод кода C с помощью Cython, а код для примера2 явно отличается и содержит дополнительную проверку нулевого деления, но когда я его просматриваю, разница в run- время незначительно.

Чтобы проиллюстрировать это, я запустил код ниже, где я взял ваш код и сделал функции def функциями cdef, то есть Cython, а не Python. Это значительно уменьшает накладные расходы при передаче и возвращении аргументов. Я также изменил example1 и example2, чтобы просто вычислить сумму по значениям в массивах numpy, а не создавать новый массив и заполнять его. Это означает, что почти все время, проведенное в каждой функции, теперь находится в цикле, поэтому должно быть легче увидеть какие-либо различия. Я также выполнял каждую функцию много раз и делал D больше.

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
@cython.cdivision(True) 
@cython.profile(True)
cdef double example1(double[:] xi, double[:] a, double[:] b, int D):

    cdef int k
    cdef double theSum = 0.0

    for k in range(D):
        theSum += (xi[k] - a[k]) / (b[k] - a[k])

    return theSum

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
@cython.profile(True)
@cython.cdivision(False)
cdef double example2(double[:] xi, double[:] a, double[:] b, int D):

    cdef int k
    cdef double theSum = 0.0

    for k in range(D):
        theSum += (xi[k] - a[k]) / (b[k] - a[k])

    return theSum


def testExamples():
    D = 100000
    x = np.random.rand(D)
    a = np.zeros(D)
    b = np.random.rand(D) + 1

    for i in xrange(10000):
        example1(x, a, b, D)
        example2(x, a, b,D) 

Я запустил этот код через профайлер (python -m cProfile -s кумулятивный), а соответствующий вывод ниже:

ncalls  tottime  percall  cumtime  percall filename:lineno(function)
10000    1.546    0.000    1.546    0.000 test.pyx:26(example2)
10000    0.002    0.000    0.002    0.000 test.pyx:11(example1)

который показывает, что пример2 намного медленнее. Если я включу c-деление в примере2, тогда затраченное время будет идентичным для примера1 и example2, поэтому это явно имеет значительный эффект.

Ответ 2

Моя проблема в том, что я вижу все, что происходит в Assembly. Попытка использовать один язык для того, чтобы говорить на другом языке о том, что именно я хочу, чтобы извлечь производительность, кажется более сложной и трудной, чем это должно быть. Например, Сеймур Крей всегда переделывал такое разделение. C=A/B стал:

R = reciprocalApprox(B);
R = reciprocalImprove(R);   //M-Add performed 3x to get an accurate 1/B
C = A * R;

Однажды Крей спросил, почему он использовал этот подход Ньютона-Рафсона, и он объяснил, что он может получить больше операций по трубопроводу, чем с аппаратным разделением. Аналогичным образом использовался AMD 3DNow, хотя с 32-битными поплавками. Для SSE с использованием gcc попробуйте флаг -mrecip, а также -funsafe-math-optimizations, -finite-math-only и -fno-trapping-math. Эта печально известная опция -ffast-math также выполняет эту функцию. Вы теряете 2 ulps или единицы на последнем месте.

http://gcc.gnu.org/onlinedocs/gcc/i386-and-x86_002d64-Options.html

Возможно, вам захочется попробовать libdivide.h(на libdivide.com). Он очень опирался на память и использует серию "дешевых" сдвигов и умножений, и заканчивается примерно в десять раз быстрее, чем целочисленное разделение. Он также генерирует код C или С++ для компилятора.