Итерация по массивам
Объект-итератор nditer, представленный в NumPy 1.6, предоставляет множество гибких способов посещения всех элементов одного или нескольких массивов систематическим образом. Эта страница знакомит с некоторыми базовыми способами использования объекта для вычислений с массивами в Python, а затем завершается демонстрацией ускорения внутреннего цикла в Cython. Поскольку Python-интерфейс nditer является относительно простым отображением API C-итератора массивов, эти идеи также помогут при работе с итерацией по массивам из C или C++.
Итерация по одному массиву
Самая базовая задача, которую можно выполнить с помощью nditer, заключается в посещении каждого элемента массива. Каждый элемент предоставляется по одному с использованием стандартного интерфейса Python-итератора.
Пример
>>> a = np.arange(6).reshape(2,3) >>> for x in np.nditer(a): ... print(x, end=' ') ... 0 1 2 3 4 5
Важно помнить, что порядок выбирается для соответствия структуре памяти массива, а не с использованием стандартного C или Fortran порядка. Это делается для повышения эффективности доступа, отражая идею, что по умолчанию требуется просто посетить каждый элемент без учета конкретного порядка. Мы можем увидеть это, итерируя по транспонированному массиву по сравнению с копированием этого транспонированного массива в порядке C.
Пример
>>> a = np.arange(6).reshape(2,3) >>> for x in np.nditer(a.T): ... print(x, end=' ') ... 0 1 2 3 4 5
>>> for x in np.nditer(a.T.copy(order='C')): ... print(x, end=' ') ... 0 3 1 4 2 5
Элементы как a , так и a.T просматриваются в одном и том же порядке, а именно в порядке их хранения в памяти, тогда как элементы a.T.copy(order=’C’) посещаются в другом порядке, потому что они были размещены в другой структуре памяти.
Управление порядком итерации
Бывают случаи, когда важно посещать элементы массива в определенном порядке, независимо от расположения элементов в памяти. Объект nditer предоставляет параметр order для управления этим аспектом итерации. По умолчанию поведение, описанное выше, имеет order=’K’ для сохранения существующего порядка. Это можно переопределить с order=’C’ для порядка C и order=’F’ для порядка Fortran.
Пример
>>> a = np.arange(6).reshape(2,3) >>> for x in np.nditer(a, order='F'): ... print(x, end=' ') ... 0 3 1 4 2 5 >>> for x in np.nditer(a.T, order='C'): ... print(x, end=' ') ... 0 3 1 4 2 5
Изменение значений массива
По умолчанию nditer обрабатывает входной операнд как объект только для чтения. Чтобы иметь возможность изменять элементы массива, необходимо указать режим чтения-записи или только записи, используя флаги ‘readwrite’ или ‘writeonly’ для каждого операнда.
nditer тогда вернет буферизованные массивы для записи, которые можно изменить. Однако, поскольку nditer должен скопировать эти данные буфера обратно в исходный массив после завершения итерации, необходимо указать, когда итерация завершена, двумя способами. Можно либо:
- использовать nditer как контекстный менеджер, используя инструкцию
with, и временные данные будут записаны обратно при выходе из контекста. - вызвать метод
closeитератора после завершения итерации, что запустит запись обратно.
После вызова close или выхода из контекста nditer больше нельзя итерировать.
Пример
>>> a = np.arange(6).reshape(2,3)
>>> a
array([[0, 1, 2],
[3, 4, 5]])
>>> with np.nditer(a, op_flags=['readwrite']) as it:
... for x in it:
... x[...] = 2 * x
...
>>> a
array([[ 0, 2, 4],
[ 6, 8, 10]])
Если вы пишете код, который должен поддерживать более старые версии numpy, обратите внимание, что до версии 1.15 nditer не был контекстным менеджером и не имел метода close. Вместо этого он полагался на деструктор для запуска записи буфера.
Использование внешнего цикла
Во всех приведенных примерах элементы a предоставляются итератором по одному, поскольку вся логика цикла находится внутри итератора. Хотя это просто и удобно, это не очень эффективно. Более эффективный подход заключается в перемещении одномерного внутреннего цикла в свой код, внешнего по отношению к итератору. Таким образом, можно использовать векторные операции NumPy для обработки больших блоков посещаемых элементов.
nditer попытается предоставить куски максимально возможного размера для внутреннего цикла. Принудительный порядок ‘C’ и ‘F’ приводит к разным размерам внешнего цикла. Этот режим включается путем задания флага итератора.
Обратите внимание, что по умолчанию при сохранении родного порядка памяти итератор может предоставить один одномерный блок, в то время как принудительном порядке Fortran он должен предоставить три блока по два элемента каждый.
Пример
>>> a = np.arange(6).reshape(2,3) >>> for x in np.nditer(a, flags=['external_loop']): ... print(x, end=' ') ... [0 1 2 3 4 5]
>>> for x in np.nditer(a, flags=['external_loop'], order='F'): ... print(x, end=' ') ... [0 3] [1 4] [2 5]
Отслеживание индекса или многомерного индекса
Во время итерации может потребоваться использовать индекс текущего элемента в вычислениях. Например, может потребоваться посещать элементы массива в порядке памяти, но использовать индекс в порядке C, Fortran или многомерный индекс для поиска значений в другом массиве.
Индекс отслеживается самим объектом итератора и доступен через свойства index или multi_index, в зависимости от запрошенного значения. В примерах ниже показаны выводы, демонстрирующие ход индекса:
Пример
>>> a = np.arange(6).reshape(2,3)
>>> it = np.nditer(a, flags=['f_index'])
>>> for x in it:
... print("%d <%d>" % (x, it.index), end=' ')
...
0 <0> 1 <2> 2 <4> 3 <1> 4 <3> 5 <5>
>>> it = np.nditer(a, flags=['multi_index'])
>>> for x in it:
... print("%d <%s>" % (x, it.multi_index), end=' ')
...
0 <(0, 0)> 1 <(0, 1)> 2 <(0, 2)> 3 <(1, 0)> 4 <(1, 1)> 5 <(1, 2)>
>>> with np.nditer(a, flags=['multi_index'], op_flags=['writeonly']) as it:
... for x in it:
... x[...] = it.multi_index[1] - it.multi_index[0]
...
>>> a
array([[ 0, 1, 2],
[-1, 0, 1]])
Отслеживание индекса или многомерного индекса несовместимо с использованием внешнего цикла, поскольку требует разных значений индексов для каждого элемента. Если вы попытаетесь объединить эти флаги, объект nditer генерирует исключение.
Пример
>>> a = np.zeros((2,3)) >>> it = np.nditer(a, flags=['c_index', 'external_loop']) Traceback (most recent call last): File "<stdin>", line 1, in <module> ValueError: Iterator flag EXTERNAL_LOOP cannot be used if an index or multi-index is being tracked
Альтернативная итерация и доступ к элементам
Для более удобного доступа к свойствам во время итерации nditer имеет альтернативную синтаксическую конструкцию для итерации, которая работает явно с самим объектом итератора. С помощью этой конструкции цикла текущее значение доступно путем индексирования итератора. Другие свойства, такие как отслеживаемые индексы, остаются прежними. Примеры ниже дают такие же результаты, как и в предыдущем разделе.
Пример
>>> a = np.arange(6).reshape(2,3)
>>> it = np.nditer(a, flags=['f_index'])
>>> while not it.finished:
... print("%d <%d>" % (it[0], it.index), end=' ')
... is_not_finished = it.iternext()
...
0 <0> 1 <2> 2 <4> 3 <1> 4 <3> 5 <5>
>>> it = np.nditer(a, flags=['multi_index'])
>>> while not it.finished:
... print("%d <%s>" % (it[0], it.multi_index), end=' ')
... is_not_finished = it.iternext()
...
0 <(0, 0)> 1 <(0, 1)> 2 <(0, 2)> 3 <(1, 0)> 4 <(1, 1)> 5 <(1, 2)>
>>> with np.nditer(a, flags=['multi_index'], op_flags=['writeonly']) as it:
... while not it.finished:
... it[0] = it.multi_index[1] - it.multi_index[0]
... is_not_finished = it.iternext()
...
>>> a
array([[ 0, 1, 2],
[-1, 0, 1]])
Буферизация элементов массива
Принудительном порядке итерации мы заметили, что внешняя петля может предоставлять элементы меньшими блоками, поскольку элементы не могут быть посещены в нужном порядке с постоянным шагом. При написании кода на C это обычно приемлемо, но в чистом Python-коде это может привести к значительному снижению производительности.
Включив буферизацию, можно увеличить блоки, предоставляемые итератором внутреннему циклу, значительно уменьшив нагрузку интерпретатора Python. В примере с принудительным порядком итерации Fortran внутренний цикл получает все элементы сразу при включенной буферизации.
Пример
>>> a = np.arange(6).reshape(2,3) >>> for x in np.nditer(a, flags=['external_loop'], order='F'): ... print(x, end=' ') ... [0 3] [1 4] [2 5]
>>> for x in np.nditer(a, flags=['external_loop','buffered'], order='F'): ... print(x, end=' ') ... [0 3 1 4 2 5]
Итерация как определенный тип данных
Иногда необходимо рассматривать массив как тип данных, отличный от того, в котором он хранится. Например, можно выполнить все вычисления с 64-битными числами с плавающей запятой, даже если обрабатываемые массивы имеют 32-битные числа с плавающей запятой. За исключением написания низкоуровневого C-кода, лучше позволить итератору обрабатывать копирование или буферизацию, а не самому преобразовывать тип данных во внутреннем цикле.
Существуют два механизма, которые позволяют это сделать: временные копии и режим буферизации. При временных копиях создается копия всего массива с новым типом данных, а затем итерация выполняется по копии. Доступ к записи разрешен в режиме, который обновляет исходный массив после завершения всей итерации. Главный недостаток временных копий заключается в том, что временная копия может занимать много памяти, особенно если тип итерируемых данных имеет больший размер элемента, чем исходный.
Режим буферизации смягчает проблему с использованием памяти и является более эффективным с точки зрения кэша, чем создание временных копий. За исключением особых случаев, когда весь массив нужен сразу вне итератора, буферизация предпочтительнее временного копирования. В NumPy буферизация используется функциями ufunc и другими функциями для поддержки гибких входных данных с минимальной нагрузкой на память.
В наших примерах мы будем рассматривать входной массив со сложным типом данных, чтобы можно было извлекать квадратные корни из отрицательных чисел. Без включения режима копий или буферизации итератор генерирует исключение, если тип данных не совпадает точно.
Пример
>>> a = np.arange(6).reshape(2,3) - 3 >>> for x in np.nditer(a, op_dtypes=['complex128']): ... print(np.sqrt(x), end=' ') ... Traceback (most recent call last): File "<stdin>", line 1, in <module> TypeError: Iterator operand required copying or buffering, but neither copying nor buffering was enabled
В режиме копирования флаг ‘copy’ задается для каждого операнда. Это делается для управления на уровне каждого операнда. Режим буферизации задается как флаг итератора.
Пример
>>> a = np.arange(6).reshape(2,3) - 3 >>> for x in np.nditer(a, op_flags=['readonly','copy'], ... op_dtypes=['complex128']): ... print(np.sqrt(x), end=' ') ... 1.7320508075688772j 1.4142135623730951j 1j 0j (1+0j) (1.4142135623730951+0j)
>>> for x in np.nditer(a, flags=['buffered'], op_dtypes=['complex128']): ... print(np.sqrt(x), end=' ') ... 1.7320508075688772j 1.4142135623730951j 1j 0j (1+0j) (1.4142135623730951+0j)
Итератор использует правила преобразования NumPy, чтобы определить, разрешено ли конкретное преобразование. По умолчанию он принуждает «безопасное» преобразование. Это означает, например, что он генерирует исключение, если вы попытаетесь рассматривать массив с 64-битными числами с плавающей запятой как массив с 32-битными числами с плавающей запятой. Во многих случаях правило «same_kind» является наиболее разумным, поскольку оно позволит преобразование с 64 на 32-битные числа с плавающей запятой, но не с числа с плавающей запятой на целое число или со сложного числа на число с плавающей запятой.
Пример
>>> a = np.arange(6.)
>>> for x in np.nditer(a, flags=['buffered'], op_dtypes=['float32']):
... print(x, end=' ')
...
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
TypeError: Iterator operand 0 dtype could not be cast from dtype('float64') to dtype('float32') according to the rule 'safe'
>>> for x in np.nditer(a, flags=['buffered'], op_dtypes=['float32'], ... casting='same_kind'): ... print(x, end=' ') ... 0.0 1.0 2.0 3.0 4.0 5.0
>>> for x in np.nditer(a, flags=['buffered'], op_dtypes=['int32'], casting='same_kind'):
... print(x, end=' ')
...
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
TypeError: Iterator operand 0 dtype could not be cast from dtype('float64') to dtype('int32') according to the rule 'same_kind'
Следует учитывать преобразования обратно к исходному типу данных при использовании операнда с режимом чтения-записи или только записи. Общий случай — реализация внутреннего цикла в терминах 64-битных чисел с плавающей запятой и использование преобразования «same_kind», чтобы позволить обрабатывать и другие типы чисел с плавающей запятой. В режиме только чтения можно предоставить целочисленный массив, но в режиме чтения-записи будет сгенерировано исключение, так как преобразование обратно в массив нарушит правило преобразования.
Пример
>>> a = np.arange(6)
>>> for x in np.nditer(a, flags=['buffered'], op_flags=['readwrite'],
... op_dtypes=['float64'], casting='same_kind'):
... x[...] = x / 2.0
...
Traceback (most recent call last):
File "<stdin>", line 2, in <module>
TypeError: Iterator requested dtype could not be cast from dtype('float64') to dtype('int64'), the operand 0 dtype, according to the rule 'same_kind'
Итерация с трансляцией массивов
NumPy имеет набор правил для работы с массивами, имеющими разные формы, которые применяются всякий раз, когда функции принимают несколько операндов, которые объединяют элементы поэлементно. Это называется трансляцией. Объект nditer может применять эти правила за вас, когда вам нужно написать такую функцию.
В качестве примера выведем результат трансляции одномерного и двумерного массивов вместе.
Пример
>>> a = np.arange(3)
>>> b = np.arange(6).reshape(2,3)
>>> for x, y in np.nditer([a,b]):
... print("%d:%d" % (x,y), end=' ')
...
0:0 1:1 2:2 0:3 1:4 2:5
При возникновении ошибки трансляции итератор генерирует исключение, которое включает формы входных данных, чтобы помочь в диагностике проблемы.
Пример
>>> a = np.arange(2)
>>> b = np.arange(6).reshape(2,3)
>>> for x, y in np.nditer([a,b]):
... print("%d:%d" % (x,y), end=' ')
...
Traceback (most recent call last):
...
ValueError: operands could not be broadcast together with shapes (2,) (2,3)
Итератор с выделенными выходными массивами
Частым случаем в функциях NumPy является выделение вывода на основе трансляции входных данных и дополнительный необязательный параметр «out», в котором результат будет размещен при его предоставлении. Объект nditer предоставляет удобный idiom, который делает этот механизм очень простым.
Покажем, как это работает, создав функцию square, которая возводит свой вход в квадрат. Начнем с минимального определения функции без поддержки параметра «out».
Пример
>>> def square(a): ... with np.nditer([a, None]) as it: ... for x, y in it: ... y[...] = x*x ... return it.operands[1] ... >>> square([1,2,3]) array([1, 4, 9])
По умолчанию nditer использует флаги «allocate» и «writeonly» для операндов, переданных как None. Это означает, что мы смогли предоставить только два операнда итератору, и он справился с остальным.
При добавлении параметра «out» мы должны явно предоставить эти флаги, потому что если кто-то передаст массив как «out», итератор по умолчанию будет «readonly», и наш внутренний цикл завершится ошибкой. Причина, по которой «readonly» является значением по умолчанию для входных массивов, заключается в предотвращении путаницы относительно непреднамеренного запуска операции сокращения. Если бы значение по умолчанию было «readwrite», любая операция трансляции также запускала бы сокращение, тема, которая рассматривается позже в этом документе.
Заодно добавим флаг «no_broadcast», который предотвратит трансляцию результата. Это важно, потому что мы хотим только одно входное значение для каждого выхода. Агрегация более одного входного значения является операцией сокращения, которая требует специальной обработки. Это уже вызвало бы ошибку, потому что сокращения должны быть явно включены во флаге итератора, но сообщение об ошибке, полученное при отключении трансляции, гораздо понятнее для конечных пользователей. Чтобы узнать, как обобщить функцию квадратного корня на сокращение, обратитесь к функции суммы квадратов в разделе о Cython.
Для полноты добавим также флаги «external_loop» и «buffered», поскольку они обычно необходимы по соображениям производительности.
Пример
>>> def square(a, out=None): ... it = np.nditer([a, out], ... flags = ['external_loop', 'buffered'], ... op_flags = [['readonly'], ... ['writeonly', 'allocate', 'no_broadcast']]) ... with it: ... for x, y in it: ... y[...] = x*x ... return it.operands[1] ...
>>> square([1,2,3]) array([1, 4, 9])
>>> b = np.zeros((3,)) >>> square([1,2,3], out=b) array([ 1., 4., 9.]) >>> b array([ 1., 4., 9.])
>>> square(np.arange(6).reshape(2,3), out=b) Traceback (most recent call last): ... ValueError: non-broadcastable output operand with shape (3,) doesn't match the broadcast shape (2,3)
Итерация внешнего произведения
Любая бинарная операция может быть расширена до операции над массивом во внешнем произведении, как в outer, а объект nditer предоставляет способ достижения этого путем явного сопоставления осей операндов. Это также возможно с помощью индексации newaxis, но мы покажем, как напрямую использовать параметр nditer op_axes для достижения этого без промежуточных представлений.
Мы выполним простое внешнее произведение, разместив измерения первого операнда перед измерениями второго операнда. Параметр op_axes требует одного списка осей для каждого операнда и предоставляет отображение осей итератора на оси операнда.
Предположим, что первый операнд одномерный, а второй операнд двумерный. Итератор будет иметь три измерения, поэтому op_axes будет иметь два списка по 3 элемента. Первый список выбирает одну ось первого операнда и имеет значение -1 для остальных осей итератора, с конечным результатом [0, -1, -1]. Второй список выбирает две оси второго операнда, но не должен перекрываться с осями, выбранными в первом операнде. Его список - [-1, 0, 1]. Операнд вывода отображается на оси итератора стандартным образом, поэтому мы можем предоставить None вместо создания другого списка.
Операция во внутреннем цикле - это простое умножение. Все, что связано с внешним произведением, обрабатывается настройкой итератора.
Пример
>>> a = np.arange(3)
>>> b = np.arange(8).reshape(2,4)
>>> it = np.nditer([a, b, None], flags=['external_loop'],
... op_axes=[[0, -1, -1], [-1, 0, 1], None])
>>> with it:
... for x, y, z in it:
... z[...] = x*y
... result = it.operands[2] # same as z
...
>>> result
array([[[ 0, 0, 0, 0],
[ 0, 0, 0, 0]],
[[ 0, 1, 2, 3],
[ 4, 5, 6, 7]],
[[ 0, 2, 4, 6],
[ 8, 10, 12, 14]]])
Обратите внимание, что после закрытия итератора мы не можем получить доступ к operands и должны использовать ссылку, созданную внутри контекстного менеджера.
Итерация сокращения
Всякий раз, когда имеющий запись операнд имеет меньше элементов, чем весь пространство итерации, этот операнд подвергается сокращению. Объект nditer требует, чтобы любой операнд сокращения был помечен как читаемый и записываемый, и разрешает сокращения только тогда, когда «reduce_ok» предоставляется в качестве флага итератора.
В качестве простого примера рассмотрим подсчет суммы всех элементов в массиве.
Пример
>>> a = np.arange(24).reshape(2,3,4) >>> b = np.array(0) >>> with np.nditer([a, b], flags=['reduce_ok'], ... op_flags=[['readonly'], ['readwrite']]) as it: ... for x,y in it: ... y[...] += x ... >>> b array(276) >>> np.sum(a) 276
Дело немного сложнее, когда комбинируются сокращение и выделяемые операнды. Перед началом итерации любой операнд сокращения должен быть инициализирован своими начальными значениями. Вот как мы можем сделать это, суммируя по последней оси a.
Пример
>>> a = np.arange(24).reshape(2,3,4)
>>> it = np.nditer([a, None], flags=['reduce_ok'],
... op_flags=[['readonly'], ['readwrite', 'allocate']],
... op_axes=[None, [0,1,-1]])
>>> with it:
... it.operands[1][...] = 0
... for x, y in it:
... y[...] += x
... result = it.operands[1]
...
>>> result
array([[ 6, 22, 38],
[54, 70, 86]])
>>> np.sum(a, axis=2)
array([[ 6, 22, 38],
[54, 70, 86]])
Для буферизованного сокращения требуется еще одна корректировка во время настройки. Обычно построение итератора включает копирование первого буфера данных из читаемых массивов в буфер. Любой операнд сокращения читаем, поэтому он может быть считан в буфер. К сожалению, инициализация операнда после завершения этой операции буферизации не отразится в буфере, с которым начинается итерация, и будут получены ложные результаты.
Флаг итератора «delay_bufalloc» существует для того, чтобы разрешить выделенным итератором операндам сокращения существовать вместе с буферизацией. При установке этого флага итератор оставит свои буферы неинициализированными до тех пор, пока не получит сброс, после чего он будет готов к стандартной итерации. Вот как выглядит предыдущий пример, если мы также включим буферизацию.
Пример
>>> a = np.arange(24).reshape(2,3,4)
>>> it = np.nditer([a, None], flags=['reduce_ok',
... 'buffered', 'delay_bufalloc'],
... op_flags=[['readonly'], ['readwrite', 'allocate']],
... op_axes=[None, [0,1,-1]])
>>> with it:
... it.operands[1][...] = 0
... it.reset()
... for x, y in it:
... y[...] += x
... result = it.operands[1]
...
>>> result
array([[ 6, 22, 38],
[54, 70, 86]])
Размещение внутреннего цикла в Cython
Те, кто хочет получить действительно высокую производительность от своих низкоуровневых операций, должны серьезно рассмотреть вопрос о непосредственном использовании API итерации, предоставленного на C, но для тех, кто не знаком с C или C++, Cython является хорошей промежуточной точкой с разумными компромиссами в плане производительности. Для объекта nditer это означает, что итератор позаботится о трансляции, преобразовании типов и буферизации, а внутренний цикл будет передан Cython.
В нашем примере мы создадим функцию суммы квадратов. Начнем с реализации этой функции в обычном Python. Мы хотим поддерживать параметр «axis», аналогичный функции numpy sum, поэтому нам нужно будет создать список для параметра op_axes. Вот как это выглядит.
Пример
>>> def axis_to_axeslist(axis, ndim): ... if axis is None: ... return [-1] * ndim ... else: ... if type(axis) is not tuple: ... axis = (axis,) ... axeslist = [1] * ndim ... for i in axis: ... axeslist[i] = -1 ... ax = 0 ... for i in range(ndim): ... if axeslist[i] != -1: ... axeslist[i] = ax ... ax += 1 ... return axeslist ... >>> def sum_squares_py(arr, axis=None, out=None): ... axeslist = axis_to_axeslist(axis, arr.ndim) ... it = np.nditer([arr, out], flags=['reduce_ok', ... 'buffered', 'delay_bufalloc'], ... op_flags=[['readonly'], ['readwrite', 'allocate']], ... op_axes=[None, axeslist], ... op_dtypes=['float64', 'float64']) ... with it: ... it.operands[1][...] = 0 ... it.reset() ... for x, y in it: ... y[...] += x*x ... return it.operands[1] ... >>> a = np.arange(6).reshape(2,3) >>> sum_squares_py(a) array(55.0) >>> sum_squares_py(a, axis=-1) array([ 5., 50.])
Для преобразования этой функции в Cython мы заменяем внутренний цикл (y[…] += x*x) кодом Cython, специализированным для типа float64. При включенном флаге «external_loop» массивы, предоставляемые внутреннему циклу, всегда будут одномерными, поэтому проверки почти не требуется.
Вот список sum_squares.pyx:
import numpy as np
cimport numpy as np
cimport cython
def axis_to_axeslist(axis, ndim):
if axis is None:
return [-1] * ndim
else:
if type(axis) is not tuple:
axis = (axis,)
axeslist = [1] * ndim
for i in axis:
axeslist[i] = -1
ax = 0
for i in range(ndim):
if axeslist[i] != -1:
axeslist[i] = ax
ax += 1
return axeslist
@cython.boundscheck(False)
def sum_squares_cy(arr, axis=None, out=None):
cdef np.ndarray[double] x
cdef np.ndarray[double] y
cdef int size
cdef double value
axeslist = axis_to_axeslist(axis, arr.ndim)
it = np.nditer([arr, out], flags=['reduce_ok', 'external_loop',
'buffered', 'delay_bufalloc'],
op_flags=[['readonly'], ['readwrite', 'allocate']],
op_axes=[None, axeslist],
op_dtypes=['float64', 'float64'])
with it:
it.operands[1][...] = 0
it.reset()
for xarr, yarr in it:
x = xarr
y = yarr
size = x.shape[0]
for i in range(size):
value = x[i]
y[i] = y[i] + value * value
return it.operands[1]
На этом компьютере создание файла .pyx в модуль выглядело следующим образом, но вам, возможно, придется найти некоторые учебники по Cython, чтобы узнать особенности вашей конфигурации системы.:
$ cython sum_squares.pyx $ gcc -shared -pthread -fPIC -fwrapv -O2 -Wall -I/usr/include/python2.7 -fno-strict-aliasing -o sum_squares.so sum_squares.c
Запуск этого из интерпретатора Python дает те же результаты, что и наш исходный код Python/NumPy.
Пример
>>> from sum_squares import sum_squares_cy >>> a = np.arange(6).reshape(2,3) >>> sum_squares_cy(a) array(55.0) >>> sum_squares_cy(a, axis=-1) array([ 5., 50.])
Небольшое измерение времени в IPython показывает, что снижение накладных расходов и выделение памяти внутреннего цикла Cython обеспечивают значительное ускорение по сравнению как со стандартным кодом Python, так и с выражением, использующим встроенную функцию суммы NumPy.:
>>> a = np.random.rand(1000,1000) >>> timeit sum_squares_py(a, axis=-1) 10 loops, best of 3: 37.1 ms per loop >>> timeit np.sum(a*a, axis=-1) 10 loops, best of 3: 20.9 ms per loop >>> timeit sum_squares_cy(a, axis=-1) 100 loops, best of 3: 11.8 ms per loop >>> np.all(sum_squares_cy(a, axis=-1) == np.sum(a*a, axis=-1)) True >>> np.all(sum_squares_py(a, axis=-1) == np.sum(a*a, axis=-1)) True
© 2005–2020 NumPy Developers
Licensed under the 3-clause BSD License.
https://numpy.org/doc/1.19/reference/arrays.nditer.html