Объяснения кода NumPy на C

Фанатизм заключается в удвоении усилий, когда вы забыли свою цель. — Джордж Сантаяна

Авторитет — это человек, который может рассказать вам больше о чем-то, чем вам действительно хочется знать. — Неизвестно

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

Модель памяти

Одним из фундаментальных аспектов ndarray является то, что массив рассматривается как «блок» памяти, начинающийся в некотором месте. Интерпретация этой памяти зависит от информации о шагах. Для каждой размерности в \(N\)-мерном массиве целое число (шаг) диктует, сколько байтов нужно пропустить, чтобы перейти к следующему элементу в этой размерности. Если у вас не один сегмент массива, эта информация о шагах должна учитываться при прохождении по массиву. Написание кода, принимающего шаги, несложно, вам просто нужно использовать указатели (char *), потому что шаги выражены в байтах. Также следует учитывать, что шаги не обязательно должны быть кратными размеру элемента. Кроме того, помните, что если количество измерений массива равно 0 (иногда называемый массивом ранга 0), то переменные strides и dimensions являются NULL.

Помимо структурной информации, содержащейся в членах strides и dimensions объекта PyArrayObject, флаги содержат важную информацию о том, как можно получить доступ к данным. В частности, флаг NPY_ARRAY_ALIGNED устанавливается, когда память находится на подходящей границе в соответствии с массивом типа данных. Даже если у вас есть непрерывный блок памяти, вы не можете просто предположить, что безопасно ссылаться на указатель типа данных на элемент. Только если установлен флаг NPY_ARRAY_ALIGNED, эта операция является безопасной (на некоторых платформах она будет работать, но на других, таких как Solaris, она вызовет ошибку шины). Флаг NPY_ARRAY_WRITEABLE также должен быть обеспечен, если вы планируете записывать в область памяти массива. Также можно получить указатель на область памяти, которая недоступна для записи. Иногда запись в область памяти, когда флаг NPY_ARRAY_WRITEABLE не установлен, просто невежливо. В других случаях это может привести к сбою программы (например, область данных — это только для чтения, отображаемый в памяти файл).

Инкапсуляция типа данных

Тип данных — важная абстракция ndarray. Операции обратятся к типу данных, чтобы получить ключевую функциональность, необходимую для работы с массивом. Эта функциональность предоставляется в списке указателей на функции, на которые указывает член ‘f’ структуры PyArray_Descr. Таким образом, количество типов данных можно расширить, просто предоставив структуру PyArray_Descr с соответствующими указателями на функции в члене ‘f’. Для встроенных типов существуют некоторые оптимизации, которые обходят этот механизм, но смысл абстракции типа данных заключается в возможности добавления новых типов данных.

Один из встроенных типов данных, тип данных void, позволяет использовать произвольные структурированные типы, содержащие 1 или более полей в качестве элементов массива. Поле — это просто другой объект типа данных вместе со смещением в текущем структурированном типе. Для поддержки произвольно вложенных полей реализованы несколько рекурсивных реализаций доступа к типам данных для типа void. Общим приёмом является цикл по элементам словаря и выполнение определённой операции на основе объекта типа данных, хранящегося по заданному смещению. Эти смещения могут быть произвольными числами. Поэтому необходимо учитывать возможность столкновения с невыровненными данными, если это необходимо.

N-мерные итераторы

Очень распространённой операцией во многом коде NumPy является необходимость итерироваться по всем элементам общего, с шагами, N-мерного массива. Эта операция общего назначения N-мерного цикла абстрагируется понятием объекта-итератора. Чтобы написать N-мерный цикл, вам нужно только создать объект-итератор из ndarray, работать с членом dataptr структуры объекта-итератора и вызвать макрос PyArray_ITER_NEXT (it) на объекте-итератора, чтобы перейти к следующему элементу. «Следующий» элемент всегда находится в C-непрерывном порядке.

Макрос работает, сначала обработав случаи C-непрерывного, 1-D и 2-D, которые работают очень просто.

В общем случае итерация работает путём отслеживания списка счётчиков координат в объекте-итератора. На каждой итерации последний счётчик координат увеличивается (начиная с 0). Если этот счётчик меньше, чем на единицу меньше размера массива в этой размерности (предвычисленное и сохранённое значение), то счётчик увеличивается, а член dataptr увеличивается на шаги в этой размерности, и макрос завершается. Если конец размерности достигнут, счётчик последней размерности сбрасывается в ноль, и dataptr перемещается обратно в начало этой размерности, вычитая значение strides, умноженное на значение на единицу меньше числа элементов в этой размерности (это также предварительно вычислено и сохранено в члене backstrides объекта-итератора). В этом случае макрос не завершается, а счётчик локальной размерности уменьшается, так что предпоследняя размерность занимает роль последней размерности, и описанные выше проверки выполняются снова для предпоследней размерности. Таким образом, dataptr корректируется должным образом для произвольных шагов.

Член coordinates структуры PyArrayIterObject сохраняет текущий N-мерный счётчик, если только базовый массив не является C-непрерывным, в этом случае подсчёт координат пропускается. Член index структуры PyArrayIterObject отслеживает текущий плоский индекс итератора. Он обновляется макросом PyArray_ITER_NEXT.

Распространение

В Numeric, предке Numpy, распространение реализовывалось в нескольких строках кода, глубоко спрятанных в ufuncobject.c. В NumPy понятие распространения было абстрагировано, чтобы его можно было выполнять в нескольких местах. Распространение обрабатывается функцией PyArray_Broadcast. Эта функция требует передачи объекта PyArrayMultiIterObject (или чего-то, что является двоичным эквивалентом).

PyArrayMultiIterObject отслеживает число размерностей распространения и размер в каждой размерности вместе с общим размером результата распространения. Он также отслеживает число массивов, которые распространяются, и указатель на итератор для каждого из распространяемых массивов.

Функция PyArray_Broadcast принимает уже определенные итераторы и использует их для определения формы распространения в каждой размерности (для создания итераторов в то же время, что и распространение, используйте функцию PyArray_MultiIterNew). Затем итераторы корректируются таким образом, чтобы каждый итератор думал, что итерируется по массиву с размером распространения. Это делается путём корректировки числа размерностей итератора и формы в каждой размерности. Это работает, потому что шаги итератора также корректируются. Распространение корректирует (или добавляет) размерности длины 1. Для этих размерностей переменная strides просто устанавливается в 0, чтобы указатель данных для итератора по этому массиву не перемещался во время операции распространения по расширенной размерности.

Распространение всегда реализовывалось в Numeric с помощью шагов со значением 0 для расширенных размерностей. Это делается точно так же в NumPy. Основное отличие заключается в том, что сейчас массив шагов отслеживается в объекте PyArrayIterObject, итераторы, участвующие в результате распространения, отслеживаются в PyArrayMultiIterObject, а вызов PyArray_Broadcast реализует правила распространения.

Массивовые скаляры

Массивовые скаляры предлагают иерархию типов Python, которые позволяют установить взаимно однозначное соответствие между типом данных, хранящимся в массиве, и типом Python, который возвращается при извлечении элемента из массива. Исключение из этого правила было сделано для массивов объектов. Массивы объектов — это гетерогенные коллекции произвольных объектов Python. Когда вы выбираете элемент из массива объектов, вы получаете исходный объект Python (а не скаляр массива объектов, который существует, но редко используется в практических целях).

Массивовые скаляры также предлагают те же методы и атрибуты, что и массивы, с намерением, что один и тот же код может использоваться для поддержки произвольных размерностей (включая размерность 0). Массивовые скаляры являются только для чтения (неизменяемые), за исключением скаляра void, который также может быть изменён, чтобы работа с полями структурированных массивов была более естественной (a[0]['f1'] = value).

Индексирование

Все операции индексирования Python arr[index] организованы путём предварительной подготовки индекса и определения типа индекса. Поддерживаемые типы индексов:

• целое число
• newaxis
• срез
• многоточие
• массивы/массивоподобные объекты с целыми числами (fancy)
• булевы (один булев массив); если в качестве индекса используется более одного булева массива или форма не соответствует точно, булев массив будет преобразован в массив целых чисел вместо этого.
• булевы 0-мерные (и также целые); 0-мерные булевы массивы — это особый случай, который нужно обрабатывать в коде расширенного индексирования. Они сигнализируют о том, что 0-мерный булев массив должен был быть интерпретирован как массив целых чисел.

А также особый случай скалярного массива, сигнализирующий о том, что целочисленный массив был интерпретирован как целочисленный индекс, что важно, поскольку целочисленный индекс массива приводит к копированию, но игнорируется, если возвращается скаляр (полный целочисленный индекс). Подготовленный индекс гарантируется как корректный, за исключением значений, выходящих за пределы границ, и ошибок трансляции для расширенного индексирования. Это включает добавление многоточия для неполных индексов, например, когда двухмерный массив индексируется одним целым числом.

Следующий шаг зависит от типа найденного индекса. Если все измерения индексируются целым числом, возвращается или устанавливается скаляр. Один массив булевых индексов вызовет специализированные булевы функции. Индексы, содержащие многоточие или срез, но не расширенное индексирование, всегда создадут представление в старом массиве, вычислив новые шаги и смещение памяти. Это представление может быть возвращено или, для присваиваний, заполнено с помощью PyArray_CopyObject. Обратите внимание, что PyArray_CopyObject также может вызываться для временных массивов в других ветвях, чтобы поддерживать сложные присваивания, когда массив имеет тип объекта.

Расширенное индексирование

Наиболее сложным случаем является расширенное индексирование, которое может или не может быть комбинировано с типичным индексированием на основе представления. Здесь целочисленные индексы интерпретируются как основанные на представлении. Прежде чем пытаться понять это, вы можете ознакомиться с его тонкостями. Код расширенного индексирования имеет три разные ветви и один специальный случай:

• Существует один массив индексов, и он, а также массив присваивания, могут быть тривиально итерированы. Например, они могут быть смежными. Кроме того, массив индексов должен быть типа intp, а массив значений в присваиваниях должен иметь правильный тип. Это чисто быстрый путь.
• Существуют только целочисленные индексы массивов, так что ни один подмассив не существует.
• Индексирование на основе представления и расширенное индексирование смешаны. В этом случае индексирование на основе представления определяет набор подмассивов, которые комбинируются расширенным индексированием. Например, arr[[1, 2, 3], :] создается путем вертикального объединения подмассивов arr[1, :], arr[2,:], и arr[3, :].
• Существует подмассив, но у него ровно один элемент. Этот случай может быть обработан так, как будто подмассива нет, но требует некоторой осторожности при настройке.

Определение применяемого случая, проверка трансляции и определение необходимого типа транспонирования выполняются в PyArray_MapIterNew. После настройки есть два случая. Если нет подмассива или он содержит только один элемент, итерация по подмассивам не требуется, и подготавливается итератор, который перебирает все массивы индексов и массив результатов или значений. Если подмассив существует, подготовлены три итератора. Один для массивов индексов, один для массива результатов или значений (без его подмассива) и один для подмассивов исходного массива и массива результатов/присваивания. Первые два итератора дают (или позволяют вычислить) указатели в начало подмассива, что затем позволяет перезапустить итерацию по подмассиву.

Когда расширенные индексы расположены рядом, может потребоваться транспонирование. Все необходимые транспонирования обрабатываются PyArray_MapIterSwapAxes и должны обрабатываться вызывающим объектом, если не нужно, чтобы PyArray_MapIterNew выделял результат.

После подготовки получение и установка относительно просты, хотя необходимо учитывать различные режимы итерации. Если во время получения элемента существует только один массив индексов, проверка корректности индексов выполняется предварительно. В противном случае это обрабатывается внутри цикла для оптимизации.

Универсальные функции

Универсальные функции — это вызываемые объекты, которые принимают \(N\) входных данных и производят \(M\) выходных данных, обертывая базовые циклы по 1-му измерению, работающие с элементом за элементом, в полноценные простые в использовании функции, которые бесшовно реализуют трансляцию, проверку типа и буферизованное приведение типов, а также обработку аргументов вывода. Новые универсальные функции обычно создаются на C, хотя существует механизм для создания ufunc из функций Python (frompyfunc). Пользователь должен предоставить цикл по 1-му измерению, который реализует базовую функцию, принимающую скалярные входные значения и помещающую полученные скаляры в соответствующие слоты вывода, как описано в реализации.

Настройка

Каждый расчет ufunc включает определённую накладную стоимость, связанную с настройкой расчёта. Практическое значение этой накладной стоимости заключается в том, что, хотя фактический расчёт ufunc очень быстрый, вы сможете написать код, специфичный для массивов и типов, который будет работать быстрее для небольших массивов, чем ufunc. В частности, использование ufunc для выполнения многих вычислений над 0-мерными массивами будет медленнее, чем другие решения на основе Python (существует молчаливо импортируемый модуль scalarmath, специально предназначенный для обеспечения массива скаляров внешнего вида расчётов на основе ufunc с значительно сниженной накладной стоимостью).

При вызове ufunc необходимо сделать много вещей. Информация, собранная из этих операций настройки, хранится в объекте цикла. Этот объект цикла — структура C (которая может стать объектом Python, но не инициализируется как таковой, потому что используется только во внутреннем использовании). Этот объект цикла имеет структуру, необходимую для использования с PyArray_Broadcast, чтобы обработка трансляции выполнялась так же, как и в других разделах кода.

Сначала выполняется поиск в глобальном словаре потока текущих значений для размера буфера, маски ошибок и связанного объекта ошибки. Состояние маски ошибок контролирует, что происходит при обнаружении условия ошибки. Следует отметить, что проверка флагов ошибок оборудования выполняется только после выполнения каждого цикла по 1-му измерению. Это означает, что если входные и выходные массивы смежные и имеют правильный тип, так что выполняется один цикл по 1-му измерению, то флаги могут не проверяться до тех пор, пока не будут вычислены все элементы массива. Поиск этих значений в словаре потока занимает время, которое легко игнорируется для всех, кроме очень малых массивов.

После проверки потоковых глобальных переменных входные данные оцениваются, чтобы определить, как ufunc должен продолжить работу, и входные и выходные массивы строятся при необходимости. Любые входные данные, которые не являются массивами, преобразуются в массивы (если необходимо, используя контекст). Замечается, какие из входных данных являются скалярами (и поэтому преобразуются в 0-мерные массивы).

Далее, из доступных ufunc выбирается подходящий цикл по 1-му измерению на основе типов входных массивов. Этот цикл по 1-му измерению выбирается путем попытки сопоставления сигнатуры типов данных входных данных с доступными сигнатурами. Сигнатуры, соответствующие встроенным типам, хранятся в члене types структуры ufunc. Сигнатуры, соответствующие пользовательским типам, хранятся в связанном списке информации о функциях, а головной элемент хранится в виде CObject в словаре userloops с ключом по номеру типа данных (используется первый пользовательский тип в списке аргументов в качестве ключа). Сигнатуры ищутся, пока не будет найдена сигнатура, в которую можно безопасно преобразовать все входные массивы (игнорируя любые скалярные аргументы, которые не могут определить тип результата). Следствием этой процедуры поиска является то, что «менее важные типы» должны располагаться ниже «более важных типов» при хранении сигнатур. Если цикл по 1-му измерению не найден, сообщается об ошибке. В противном случае список аргументов обновляется сохранённой сигнатурой — в случае необходимости преобразования и для исправления типов вывода, предполагаемых циклом по 1-му измерению.

Если ufunc имеет 2 входа и 1 выход, и второй вход — массив объектов, выполняется специальная проверка, чтобы при возвращении NotImplemented второй вход не был ndarray, не имел атрибута __array_priority__ и не имел специального метода __r{op}__. Таким образом, Python получает сигнал, чтобы дать другому объекту возможность завершить операцию, вместо использования общих расчётов массивов объектов. Это позволяет (например) разреженным матрицам переопределять оператор умножения в цикле по 1-му измерению.

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

Затем обрабатываются аргументы вывода (если таковые имеются), и строятся любые отсутствующие выходные массивы. Если какой-либо предоставленный выходной массив не имеет правильного типа (или не выровнен) и меньше размера буфера, то создаётся новый выходной массив со специальным флагом NPY_ARRAY_WRITEBACKIFCOPY. В конце функции вызывается PyArray_ResolveWritebackIfCopy, чтобы его содержимое было скопировано обратно в выходной массив. Затем обрабатываются итераторы для аргументов вывода.

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

Вызов функции

Этот раздел описывает, как настраивается и выполняется базовый цикл вычисления универсальной функции для каждого из трёх типов выполнения. Если NPY_ALLOW_THREADS определён во время компиляции, то при отсутствии массивов объектов блокировка глобального интерпретатора Python (GIL) освобождается перед вызовом циклов. Она снова приобретается при необходимости для обработки условий ошибок. Флаги ошибок оборудования проверяются только после завершения цикла по 1-му измерению.

Один цикл

Это самый простой случай. Ufunc выполняется путём вызова основного цикла по 1-му измерению ровно один раз. Это возможно только тогда, когда у нас есть выровненные данные правильного типа (включая порядок байтов) как для входных, так и для выходных данных, и все массивы имеют равные шаги (смежные, 0-мерные или 1-мерные). В этом случае цикл вычислений по 1-му измерению вызывается один раз для вычисления расчёта для всего массива. Обратите внимание, что флаги ошибок оборудования проверяются только после завершения всего расчёта.

Цикл со смещением

Когда входные и выходные массивы выровнены и имеют правильный тип, но шаги не однородны (не смежные, а 2-мерные или более), для расчёта используется вторая структура цикла. Этот подход преобразует все итераторы для входных и выходных аргументов в итерацию по всем измерениям, кроме наибольшего. Внутренний цикл затем обрабатывается основным циклом вычислений по 1-му измерению. Внешний цикл — это стандартный цикл итерации по преобразованным итераторам. Флаги ошибок оборудования проверяются после завершения каждого цикла по 1-му измерению.

Буферизованный цикл

Вот код, который обрабатывает ситуации, когда входные и/или выходные массивы не выровнены или имеют неправильный тип данных (включая перестановку байтов) по сравнению с ожидаемым подлежащим 1-мерным циклом. Также предполагается, что массивы не являются непрерывными. Код очень похож на код с шагом, за исключением того, что внутренний 1-мерный цикл модифицирован таким образом, что предварительная обработка выполняется над входными данными, а пост-обработка — над выходными данными в блоках по bufsize (где bufsize — параметр, настраиваемый пользователем). Подлежащий 1-мерный вычислительный цикл вызывается с данными, которые копируются (если это необходимо). Код настройки и цикла в этом случае значительно сложнее, потому что он должен обрабатывать:

• выделение памяти для временных буферов
• определение, следует ли использовать буферы для входных и выходных данных (невыровненные и/или неправильный тип данных)
• копирование и, возможно, приведение типов данных для всех входных или выходных данных, для которых необходимы буферы.
• специальное обращение с массивами типа Object, чтобы правильно обрабатывать счетчики ссылок при необходимости копирования и/или приведения типов.
• разбиение внутреннего 1-мерного цикла на блоки по bufsize (с возможным остатком).

Опять же, флаги ошибок аппаратного обеспечения проверяются в конце каждого 1-мерного цикла.

Конечная обработка выходных данных

Ufunc позволяют беспрепятственно передавать другие массивоподобные классы через интерфейс, так как входные данные определенного класса приведут к тому, что выходные данные будут того же класса. Механизм, с помощью которого это работает, следующий. Если какие-либо из входных данных не являются ndarrays и определяют метод __array_wrap__, то класс с наибольшим атрибутом __array_priority__ определяет тип всех выходных данных (за исключением любых выходных массивов, которые были переданы). Метод __array_wrap__ входного массива будет вызван с ndarray, возвращаемым ufunc, в качестве входного аргумента. Поддерживаются два стиля вызова функции __array_wrap__. Первый принимает ndarray в качестве первого аргумента и кортеж «контекста» в качестве второго аргумента. Контекст — это (ufunc, аргументы, номер выходного аргумента). Это первый вызов, который будет выполнен. Если возникает ошибка TypeError, функция вызывается только с ndarray в качестве первого аргумента.

Методы

Существует три метода ufunc, которые требуют вычислений, аналогичных универсальным ufunc. Это reduce, accumulate и reduceat. Каждый из этих методов требует команды настройки, за которой следует цикл. Возможны четыре стиля циклов для методов: без элементов, один элемент, цикл с шагом и цикл с буфером. Это те же основные стили циклов, что и при реализации вызова универсальной функции, за исключением случаев без элементов и одного элемента, которые являются особыми случаями, возникающими, когда входные массивы имеют 0 и 1 элемент соответственно.

Настройка

Функция настройки для всех трех методов — construct_reduce. Эта функция создает объект цикла сокращения и заполняет его параметрами, необходимыми для выполнения цикла. Все методы работают только с ufunc, принимающими 2 входных значения и возвращающими 1 выходное значение. Таким образом, подлежащий 1-мерный цикл выбирается при условии подписи [ otype, otype, otype ], где otype — запрашиваемый тип данных сокращения. Размер буфера и обработка ошибок затем извлекаются из глобального хранилища (по потоку). Для небольших массивов, которые не выровнены или имеют неправильный тип данных, создается копия, чтобы использовался небуферизованный фрагмент кода. Затем выбирается стратегия циклирования. Если в массиве 1 элемент или 0 элементов, выбирается простой метод циклирования. Если массив не невыровнен и имеет правильный тип данных, выбирается циклирование с шагом. В противном случае необходимо выполнить буферизованное циклирование. Затем устанавливаются параметры циклирования, и создается возвращаемый массив. Форма выходного массива отличается в зависимости от того, является ли метод reduce, accumulate или reduceat. Если выходной массив уже предоставлен, проверяется его форма. Если выходной массив не является C-непрерывным, невыровненным и не имеет правильного типа данных, создается временная копия с установленным флагом WRITEBACKIFCOPY. Таким образом, методы смогут работать с корректным выходным массивом, но результат будет скопирован обратно в фактический выходной массив, когда PyArray_ResolveWritebackIfCopy вызывается по завершении функции. Наконец, инициализируются итераторы для перебора соответствующей оси (в зависимости от значения axis, переданного методу), и процедура настройки возвращается к фактической вычислительной процедуре.

Reduce

Все методы ufunc используют те же подлежащие 1-мерные вычислительные циклы со входными и выходными аргументами, скорректированными для выполнения соответствующего сокращения. Например, ключ к работе reduce заключается в том, что 1-мерный цикл вызывается с выходными данными и вторым входом, указывающими на одно и то же положение в памяти, и оба имеют размер шага 0. Первый вход указывает на входной массив с шагом, заданным соответствующим шагом для выбранной оси. Таким образом, выполняемая операция

где \(N+1\) — количество элементов во входных данных, \(i\), \(o\) — выходные данные, а \(i[k]\) — \(k^{\textrm{th}}\) элемент \(i\) вдоль выбранной оси. Эта базовая операция повторяется для массивов с более чем одной размерностью, чтобы сокращение выполнялось для каждого 1-мерного подмассива вдоль выбранной оси. Итератор с удаленной выбранной размерностью обрабатывает это циклирование.

Для буферизованных циклов необходимо позаботиться о копировании и приведении типов данных перед вызовом функции цикла, потому что подлежащий цикл ожидает выровненные данные правильного типа данных (включая порядок байтов). Буферизованный цикл должен обрабатывать это копирование и приведение типов перед вызовом функции цикла на блоках, не превышающих заданного пользователем параметра bufsize.

Accumulate

Функция accumulate очень похожа на функцию reduce, так как выходные данные и второй вход оба указывают на выходные данные. Разница в том, что второй вход указывает на память на один шаг позади текущего указателя на выходные данные. Таким образом, выполняемая операция

Выходные данные имеют ту же форму, что и входные данные, и каждый 1-мерный цикл работает с \(N\) элементами, когда форма на выбранной оси — \(N+1\). Опять же, буферизованные циклы следят за копированием и приведением типов данных перед вызовом подлежащего 1-мерного вычислительного цикла.

Reduceat

Функция reduceat — это обобщение функций reduce и accumulate. Она реализует сокращение по диапазонам входного массива, заданных индексами. Дополнительный аргумент indices проверяется, чтобы убедиться, что каждый вход не слишком большой для входного массива вдоль выбранного измерения перед вычислениями цикла. Реализация цикла обрабатывается кодом, очень похожим на код reduce, повторяемый столько раз, сколько элементов в входных данных. В частности: первый указатель на вход, переданный подлежащему 1-мерному вычислительному циклу, указывает на входной массив в правильной позиции, указанной массивом индексов. Кроме того, указатель на выходные данные и второй указатель на входные данные, переданные подлежащему 1-мерному циклу, указывают на одно и то же положение в памяти. Размер 1-мерного вычислительного цикла фиксируется как разность между текущим индексом и следующим индексом (если текущий индекс — последний индекс, следующий индекс предполагается равным длине массива вдоль выбранного измерения). Таким образом, 1-мерный цикл реализует сокращение по указанным индексам.

Обработка невыровненных данных или типа данных цикла, не соответствующего типу данных входных и/или выходных данных, выполняется с помощью буферизованного кода, в котором данные копируются во временной буфер и приводятся к нужному типу данных, если это необходимо, перед вызовом подлежащей 1-мерной функции. Временные буферы создаются размером (элементов), не превышающим значения параметра буфера, устанавливаемого пользователем. Таким образом, цикл должен быть достаточно гибким, чтобы вызывать подлежащий 1-мерный вычислительный цикл достаточно раз, чтобы завершить общее вычисление частями, не превышающими размер буфера.