Почему нет numpy.mean multithreaded?

Я искал способы легко многопоточить некоторые из моего простого кода анализа, так как я заметил numpy, что он использовал только одно ядро, несмотря на то, что он должен быть многопоточным.

Я знаю, что numpy настроен для нескольких ядер, так как я могу видеть тесты с использованием numpy.dot, используя все мои ядра, поэтому я просто переопределял среднее значение как точечный продукт, и он работает быстрее. Есть ли какая-то причина, которая не может справиться с этим быстро? Я нахожу подобное поведение для больших массивов, хотя отношение близко к 2, чем показанное в моем примере.

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

In [27]: data = numpy.random.rand(10,10)

In [28]: a = numpy.ones(10)

In [29]: %timeit numpy.dot(data,a)/10.0
100000 loops, best of 3: 4.8 us per loop

In [30]: %timeit numpy.mean(data,axis=1)
100000 loops, best of 3: 14.8 us per loop

In [31]: numpy.dot(data,a)/10.0 - numpy.mean(data,axis=1)
Out[31]: 
array([  0.00000000e+00,   0.00000000e+00,   0.00000000e+00,
         0.00000000e+00,   1.11022302e-16,   0.00000000e+00,
         0.00000000e+00,   0.00000000e+00,   0.00000000e+00,
        -1.11022302e-16])

Ответ 1

Я искал способы легко многопоточить некоторые из моего простого кода анализа, так как я заметил numpy, что он использовал только одно ядро, несмотря на то, что он должен быть многопоточным.

Кто сказал, что он должен быть многопоточным?

numpy в первую очередь предназначен для того, чтобы быть как можно быстрее на одном ядре и быть как можно более распараллеливаемым, если вам это нужно. Но вам все равно придется распараллелить его.

В частности, вы можете одновременно работать с независимыми субобъектами, а медленные операции освобождают GIL, когда это возможно, хотя "когда возможно" может быть недостаточно. Кроме того, объекты numpy предназначены для совместного использования или передачи между процессами как можно проще, чтобы облегчить использование multiprocessing.

Есть несколько специализированных методов, которые автоматически распараллеливаются, но большинство основных методов - нет. В частности, dot реализуется поверх BLAS, когда это возможно, и BLAS автоматически распараллеливается на большинстве платформ, но mean реализуется в виде простого кода C.

Подробнее см. Параллельное программирование с numpy и scipy.

Итак, как вы знаете, какие методы распараллеливаются, а какие нет? И, из тех, которые не являются, как вы знаете, какие из них могут быть хорошо обработаны вручную и которые требуют многопроцессорности?

Там нет хорошего ответа. Вы можете сделать обоснованные догадки (X кажется, что он, вероятно, реализован поверх ATLAS, и моя копия ATLAS неявно переименована), или вы можете прочитать источник.

Но обычно лучше всего попробовать и попробовать. Если код использует 100% одного ядра и 0% других, добавьте ручную резьбу. Если теперь он использует 100% одного ядра и 10% других и едва работает быстрее, измените многопоточность на многопроцессорную обработку. (К счастью, Python делает это довольно легко, особенно если вы используете классы Executor из concurrent.futures или классы пула из multiprocessing. Но вам по-прежнему часто нужно вдуматься в это и проверить относительные затраты на совместное использование vs. если у вас есть большие массивы.)

Кроме того, как указывает kwatford, просто потому, что какой-то метод не кажется неявно параллельным, не означает, что он не будет параллелен в следующей версии numpy или следующей версии BLAS или другой платформу или даже на машине с немного отличающимися на ней материалами. Поэтому будьте готовы к повторной проверке. И сделайте что-нибудь вроде my_mean = numpy.mean, а затем используйте my_mean всюду, поэтому вы можете просто изменить одну строку на my_mean = pool_threaded_mean.

Ответ 2

В принципе, поскольку библиотека BLAS имеет оптимизированный точечный продукт, который они могут легко вызвать для dot, который по своей сути параллелен. Они признают, что могут расширять numpy, чтобы распараллелить другие операции, но отказались от этого маршрута. Тем не менее, они дают несколько советов о том, как распараллелить ваш код numpy (в основном, чтобы разделить работу между N ядрами (например, N = 4), разделить ваш массив на N подмассивов и отправить задания для каждого подматрица в собственный поток и затем объедините результаты).

См. http://wiki.scipy.org/ParallelProgramming:

Использовать параллельные примитивы

Одна из сильных сторон numpy заключается в том, что вы можете выразить операции массива очень чисто. Например, чтобы вычислить произведение матрицы A и матрицы B, вы просто выполните:

>>> C = numpy.dot(A,B)

Это не просто и понятно для чтения и записи, поскольку numpy знает, что вы хотите сделать продукт с матричной точкой, он может использовать оптимизированную реализацию, полученную как часть "BLAS" (подпрограммы базовой линейной алгебры). Обычно это библиотека, тщательно настроенная для быстрого запуска на вашем аппаратном обеспечении, используя преимущества кеш-памяти и реализации ассемблера. Но у многих архитектур теперь есть BLAS, который также использует многоядерную машину. Если ваш numpy/scipy скомпилирован с использованием одного из них, то точка() будет вычисляться параллельно (если это выполняется быстрее), если вы ничего не делаете. Аналогично для других матричных операций, таких как инверсия, декомпозиция сингулярных значений, детерминант и т.д. Например, библиотека с открытым исходным кодом ATLAS позволяет подбирать время выбора уровня parallelism (количество потоков). Собственная библиотека MKL от Intel предлагает возможность выбрать уровень parallelism во время выполнения. Существует также библиотека GOTO, которая позволяет выбирать время parallelism во время выполнения. Это коммерческий продукт, но исходный код распространяется бесплатно для академического использования.

Наконец, scipy/numpy не распараллеливает такие операции, как

>>> A = B + C

>>> A = numpy.sin(B)

>>> A = scipy.stats.norm.isf(B)

Эти операции выполняются последовательно, не пользуясь преимуществами многоядерных машин (но см. ниже). В принципе, это можно было бы изменить без особых усилий. OpenMP - это расширение для языка C, которое позволяет компиляторам создавать параллелизующий код для соответствующих аннотированных циклов (и других вещей). Если кто-то сел и аннотировал несколько циклов ядра в numpy (и, возможно, в scipy), и если один из них скомпилировал numpy/scipy с включенным OpenMP, все три из них автоматически будут запускаться параллельно. Конечно, на самом деле можно было бы иметь некоторый контроль времени исполнения - например, можно было бы отключить автоматическую распараллеливание, если бы планировалось запустить несколько заданий на одной и той же машине с несколькими процессорами.