Подтвердить что ты не робот

Скорость вращения шпинделя против Китона

У меня есть код анализа, который выполняет некоторые тяжелые числовые операции с использованием numpy. Просто для любопытства, попытался скомпилировать его с помощью cython с небольшими изменениями, а затем я переписал его с использованием циклов для части numpy.

К моему удивлению, код, основанный на циклах, был намного быстрее (8x). Я не могу опубликовать полный код, но я собрал очень простое несвязанное вычисление, которое показывает подобное поведение (хотя разница во времени не такая большая):

Версия 1 (без cython)

import numpy as np

def _process(array):

    rows = array.shape[0]
    cols = array.shape[1]

    out = np.zeros((rows, cols))

    for row in range(0, rows):
        out[row, :] = np.sum(array - array[row, :], axis=0)

    return out

def main():
    data = np.load('data.npy')
    out = _process(data)
    np.save('vianumpy.npy', out)

Версия 2 (создание модуля с cython)

import cython
cimport cython

import numpy as np
cimport numpy as np

DTYPE = np.float64
ctypedef np.float64_t DTYPE_t

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
cdef _process(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.zeros((rows, cols))

    for row in range(0, rows):
        out[row, :] = np.sum(array - array[row, :], axis=0)

    return out

def main():
    cdef np.ndarray[DTYPE_t, ndim=2] data
    cdef np.ndarray[DTYPE_t, ndim=2] out
    data = np.load('data.npy')
    out = _process(data)
    np.save('viacynpy.npy', out)

Версия 3 (создание модуля с помощью cython)

import cython
cimport cython

import numpy as np
cimport numpy as np

DTYPE = np.float64
ctypedef np.float64_t DTYPE_t

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
cdef _process(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.zeros((rows, cols))

    for row in range(0, rows):
        for col in range(0, cols):
            for row2 in range(0, rows):
                out[row, col] += array[row2, col] - array[row, col]

    return out

def main():
    cdef np.ndarray[DTYPE_t, ndim=2] data
    cdef np.ndarray[DTYPE_t, ndim=2] out
    data = np.load('data.npy')
    out = _process(data)
    np.save('vialoop.npy', out)

С матрицей 10000x10, сохраненной в data.npy, время:

$ python -m timeit -c "from version1 import main;main()"
10 loops, best of 3: 4.56 sec per loop

$ python -m timeit -c "from version2 import main;main()"
10 loops, best of 3: 4.57 sec per loop

$ python -m timeit -c "from version3 import main;main()"
10 loops, best of 3: 2.96 sec per loop

Ожидается ли это или есть оптимизация, которую мне не хватает? Тот факт, что версии 1 и 2 дает тот же результат, как-то ожидалось, но почему версия 3 быстрее?

Ps.- Это НЕ расчет, который мне нужно сделать, просто простой пример, который показывает одно и то же.

4b9b3361

ОТВЕТЫ

Ответ 1

Как упоминалось в других ответах, версия 2 по существу такая же, как и версия 1, поскольку cython не может выполнить поиск в операторе доступа к массиву, чтобы оптимизировать его. Для этого есть 2 причины

  • Во-первых, в каждом вызове функция numpy имеет определенный объем накладных расходов по сравнению с оптимизированным кодом C. Однако эти накладные расходы станут менее значительными, если каждая операция связана с большими массивами

  • Во-вторых, существует создание промежуточных массивов. Это более понятно, если вы рассматриваете более сложную операцию, например out[row, :] = A[row, :] + B[row, :]*C[row, :]. В этом случае в память должен быть создан целый массив B*C, а затем добавлен в A Это означает, что кэширование ЦП прерывается, поскольку данные считываются и записываются в память, а не сохраняются в ЦП и используются сразу. Важно отметить, что эта проблема ухудшается, если вы имеете дело с большими массивами.

В частности, поскольку вы заявляете, что ваш реальный код более сложный, чем ваш пример, и он показывает гораздо большее ускорение, я подозреваю, что вторая причина, вероятно, будет основным фактором в вашем случае.

В стороне, если ваши вычисления достаточно просты, вы можете преодолеть этот эффект, используя numexpr, хотя, конечно, cython полезен во многих других ситуациях, поэтому это может быть лучший подход для вас.

Ответ 2

С небольшой модификацией версия 3 становится в два раза быстрее:

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
def process2(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row, col, row2
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.empty((rows, cols))

    for row in range(rows):
        for row2 in range(rows):
            for col in range(cols):
                out[row, col] += array[row2, col] - array[row, col]

    return out

Узким местом в вашем расчете является доступ к памяти. Ваш входной массив C упорядочен, что означает, что движение по последней оси делает наименьший скачок в памяти. Поэтому ваш внутренний цикл должен быть вдоль оси 1, а не оси 0. При этом это изменение сокращает время выполнения в два раза.

Если вам нужно использовать эту функцию на небольших входных массивах, вы можете уменьшить накладные расходы, используя np.empty вместо np.ones. Чтобы уменьшить накладные расходы, используйте PyArray_EMPTY из API-интерфейса numpy C.

Если вы используете эту функцию на очень больших входных массивах (2 ** 31), то целые числа, используемые для индексирования (и в функции range), будут переполняться. Для безопасного использования:

cdef Py_ssize_t rows = array.shape[0]
cdef Py_ssize_t cols = array.shape[1]
cdef Py_ssize_t row, col, row2

вместо

cdef unsigned int rows = array.shape[0]
cdef unsigned int cols = array.shape[1]
cdef unsigned int row, col, row2

Timing:

In [2]: a = np.random.rand(10000, 10)
In [3]: timeit process(a)
1 loops, best of 3: 3.53 s per loop
In [4]: timeit process2(a)
1 loops, best of 3: 1.84 s per loop

где process - ваша версия 3.

Ответ 3

Я бы порекомендовал использовать флаг -a, чтобы cython сгенерировал html файл, который показывает, что происходит в чистом c, и вызывает API-интерфейс python:

http://docs.cython.org/src/quickstart/cythonize.html

Версия 2 дает почти тот же результат, что и в версии 1, потому что весь тяжелый подъем выполняется с помощью API Python (через numpy), а cython ничего не делает для вас. Фактически на моей машине numpy построен против MKL, поэтому, когда я компилирую cython сгенерированный c-код с помощью gcc, версия 3 на самом деле немного медленнее, чем две другие.

Cython светит, когда вы выполняете манипуляцию массивом, которую numpy не может выполнять "векторизованным" способом, или когда вы делаете что-то интенсивное, что позволяет избежать создания большого временного массива. Я получил 115-кратное ускорение, используя cython vs numpy для моего собственного кода:

https://github.com/synapticarbors/pylangevin-integrator

Часть этого вызова вызывала каталог randomkit на уровне кода c вместо того, чтобы вызывать его через numpy.random, но большая часть из них была cython, переводящая вычислительно интенсивную для циклов в чистую c без вызовов на python.

Ответ 4

Разница может быть связана с тем, что версии 1 и 2 выполняют вызов уровня на уровне Python для np.sum() для каждой строки, а версия 3, вероятно, компилируется в жесткий, чистый цикл C.

Изучая разницу между версиями 2 и 3, источником C-генерации, полученным Cython, следует просвещать.

Ответ 5

Я предполагаю, что основные накладные расходы, которые вы сохраняете, - это созданные временные массивы. Вы создаете большой массив array - array[row, :], а затем уменьшаете его до меньшего массива с помощью sum. Но создание этого большого временного массива не будет бесплатным, особенно если вам нужно выделить память.