Pull to refresh

Comments 47

И вот почему в Python 3 вместо решения одной проблемы — GIL — создали дополнительные в виде юникода, бинарей и др?
И вот почему некоторые не видят разницы между языком и его интерпретаторами?
И почему некоторые считают, что язык, у которого только одна реализация, и на особенности которой уже заточена вся инфраструктура, не равен его интерпретатору?

Когда-нибудь обязательно настанет счастье, но сейчас…
Основные реализации — CPython, PyPy, IronPython, Jython
Да, сейчас стало чуть легче, и IronPython более-менее исправляет ситуацию. Но неслабый пласт модулей прибиты к CPython гвоздями.
Тот же numpy, matplotlib, pygtk…
Люди, которые не понимают, зачем в Py3 unicode стал дефолтным строковым типом, раньше говнокодили страницы на cp1251 в пхп.
Не только «говнокодили страницы на cp1251 в пхп», но ещё писали под DOS, например (когда Python собственно и не было).

Реально Py3 и Py2 — это разные языки. И продолжение Py2 это скорее Pyston, а не Py3. А кому нужен дефолтный юникод, скорее будет уже выбирать между Py3 и Go.
Потому что, как говорит великодушный пожизненный диктатор, GIL, хоть и имеет слабые места, он просто не создан для того, чтобы делать многопроцессорные приложения, а скорее для работы с IO. multiprocessing для описанных в топике задач в помощь, если нужен чистый python.
Потоки от процессов отличаются наличием shared memory, на которую автор сделал упор. Между процессами такое придется делать руками, и желательно прямыми. Либо отказаться от общей памяти и делать общение через пайпы. Пайпы кстати автоматизированы для multiprocessing, например, очередями.
Ваша фантазия про отличия потоков от процессов имеет мало общего с нашей реальностью.
Мало того, что SM изначально является механизмом IPC и к потокам отношения по сути не имеет, различия между процессами и потоками совсем не в этом.
Извините, я конкретно про отличия питоновских модулей Threading и Multiprocessing и их функций из коробки.
Ничего страшного. Но в питоновских модулях с SM тоже всё достаточно просто — multiprocessing.sharedctypes. Ну или если это называется «руками», то я уже не знаю, что значит не-руками. =)
Сейчас расскажу, что значит не руками. Предположим, мы отправляем в новые потоки ссылку на некий объект. Так вот в случае с потоками, в них уйдут ссылки на объект, а в случае с процессами, (дальше я точно не интересовался, как это работает, но предположу) объект запиклится и уйдет копия в каждый подпроцесс. Это я имел в виду под словом «из коробки», возможно неправильно выразился, виноват.

Ну вот что-то типа того:
from threading import Thread, Lock
from multiprocessing import Process, Lock as PLock

class SharedObject(object):
    def __init__(self, property):
        self.property = property

def some_func(shared_obj, lock):
    lock.acquire()
    shared_obj.property += 1
    print shared_obj.property,
    lock.release()

if __name__ == '__main__':
    # работа с потоками:
    lock = Lock()
    shared_obj = SharedObject(0)
    for thread_number in range(5):
        thread = Thread(target=some_func, args=(shared_obj, lock))
        thread.start()
    # 1 2 3 4 5

    # работа с процессами:
    lock = PLock()
    shared_obj = SharedObject(0)
    for process_number in range(5):
        process = Process(target=some_func, args=(shared_obj, lock))
        process.start()
    # 1 1 1 1 1



Для работы с IO GIL отпускает структуры.
Я только не понимаю, почему нельзя использовать потокобезопасные структуры в современном мире?
Ох и переписывать там придется (если рассматривать CPython). Тут мне хочется надеяться что PyPy-stm будет более прогрессивным.
Ну так в том и суть, что python 3 это смена версии крупная, можно было и переписать, а то выходит это основной ограничитель развития языка на данный момент. Т.к. в других языках это уже решено.
Есть интересное видео David Beasley про анализ проблем GIL (youtube link). Там говорится, что есть надежда на улучшение GIL (не избавление, а именно исправление существующих недостатков реализации) в 3 версии Python.
Смотришь на этот код, и возникает ощущение, что писать на плюсах таки гораздо проще. %)
Автор не весь код на C++ выложил. Как минимум перегрузка * осталась за кадром.
Для плюсов тоже есть математические библиотеки, в которых всё нужное уже реализовано. Такого уровня велосипедостроение не нужно.
Я сам очень сильно люблю C++. Но вот в чем дело. Надо разделять две задачи.

Первая задача — это промышленная реализация известного алгоритма в рамках системы с уже выработанным дизайном. И тут я бы тоже выбрал C++.

Вторая задача — это активный поиск нужного алгоритма, его прототипирование, анализ и доработка до рабочего состояния. И здесь системы научных расчетов типа MATLAB или Scientific Python просто незаменимы. Часто за день нужно по нескольким научных статьям попробовать алгоритмы, использующие совсем разную математику. И системы расчетов предоставляют тебе готовую интерактивную среду с очень удобными абстракциями, позволяющие мговенно визуализировать результаты.
На плюсах писать проще, когда есть сформулированная задача, есть время продумать архитектуру, классы, вот это всё. А когда идёт научный поиск, то первоначальные абстракции могут быть совершенно неправильными, так как понимание задачи меняется в процессе решения. И вот питоне, в силу гибкости, можно довольно быстро поменять всё в программе или применить какой-нибудь dirty hack и проверить идею о которой в начале даже не задумывались (например, приватные переменные класса подкрутить). А вот когда ужё всё отлажено можно и посадить студента на C++ переписать.
А numexpr не пробовали? Если заменить строку в вашем numpy-варианте на R = numexpr.evaluate('1 / (1 + sqrt(Rx * Rx + Ry * Ry + Rz * Rz))'), то у меня сильно ускоряется (стаёт 0.6 секунд против 1.4 с numpy). Тут и многопоточность используется (процессор двухядерный).
Кстати, если то, что приведено в посте — не пример, а действительно именно такие вычисления нужны, то есть готовая реализация: scipy.spatial.distance.cdist. А именно, просто cdist(p, q) (считает просто ||p-q||) работает у меня 0.35 с, а 1/(1+cdist(p, q)) (считает искомую матрицу) — 0.6 с. Ну и, numexpr.evaluate('1/(1+d)', local_dict={'d': cdist(p, q)}) работает менее 0.5 c :)
Спасибо за комментарии! Numexpr действительно еще один способ параллельных вычислений на Python. Насколько я помню, numexpr использует собственную виртуальную машину, написанную на C, и собственный JIT компилятор. Как-то получилось, что я numexp редко использую. Часто было, что сложные выражение он не смог распарсить. И часто часть кода, которую хотелось распараллелить, было проще выделять в крупные блоки. Но действительно, numexpr для простых выражений очень удобный способ.

По поводу cdist да, все верно. Но для поста я такую функцию выбрал для примера, чтобы было не слишком сложно, но вместе с тем и полезно.

Добавлю в PS к посту эти моменты. Спасибо!
Если вы хотите рассмотреть всевозможные пути улучшения производительности, то попробуйте ещё Intel MKL или IPP (смотря что подходит в конкретных задачах). Использовать их из Cython очень просто, и получается весьма значительный прирост в скорости — у меня от использования только BLAS Level 1 (операции с векторами) из Intel MKL вычисления ускорились в пару раз по сравнению с наивным C/Cython кодом. Про умножение матриц даже говорить нечего.
Мораль: в питоне быстро работает то, что написано на C, а не на питоне.
Но Питон и сам написан на Си…
На таких операциях (число операций линейно по объёму данных), и тем более таких размерах данных, странно ожидать прироста скорости от видеокарты. Конечно, кроме случая, когда до/после операции данные и так нужны на видеокарте.
Я недавно делал доклад на похожую тему — тоже с разбором полетов кода на чистом питоне, с постепенной заменой кусков на плюсы, и там же сбоку numpy для сравнения — только там задачкой было множество Мандельброта. И в числе прочих вариантов на плюсах там был вариант на C++AMP. Так вот, прирост производительности там был очень серьезный, где-то на порядок относительно многопоточной версии на PPL.
На днях довелось пообщаться с людьми из коммьюнити scientific Python, в т.ч. Оливье Гризелем и Брайаном Грейнджером. Так вот, в дискуссии на аналогичную тему, они сделали такое интересное замечание: производительность — это хорошо, но в рамках научных исследований важно, чтобы читатели всего этого дела могли потом легко разобраться в коде и воспроизвести те же подсчеты. Поэтому зачастую лучше использовать более медленный numpy вместо всякой экзотики, просто потому, что его знают и понимают все в данной области.
Для cuda, кстати, я приглядываюсь к NumbaPro (http://docs.continuum.io/numbapro/), это расширенный вариант Numba. Синтаксис там такой же простой, одной лишь аннотацией задача комплируется под CUDA:

from numbapro import cuda

@cuda.jit('void(float32[:], float32[:], float32[:])')
def sum(a, b, result):
    i = cuda.grid(1)   # equals to threadIdx.x + blockIdx.x * blockDim.x
    result[i] = a[i] + b[i]

Но тестов пока не делал.
UFO just landed and posted this here
В новой версии питона появился ещё модуль statistics. Много не щупал, но то, что нужно было — работало хорошо и интуитивно по интерфейсам.
Интересно было бы сделать то же самое на Java (а может даже на Jython), и сравнить с результатами автора. Судя по тому, что Java хорошо умеет много потоков вокруг shared memory для потоков. Вдруг это оказалось бы эффективно, особенно если использовать готовые Java-библиотеки для расчётов (которые внутри, конечно же, сильно параллельные).
Для научных дел Java это ещё хуже чем C++. Да и вряд ли она будет считать быстрее чем си.
Раньше не слышал о Numba. Но раз уж о нем (ней?) зашла речь, то интересно было бы добавить в сравнение Nuitka. Понятно, что статью вы переписывать не будете, но буду благодарен, если в коментах поделитесь хотя бы какими-нибудь субъективными впечатлениями/ощущениями.
Большое спасибо за статью! Я сам физик, использую scipy/numpy в своей работе постоянно. Тоже прошёл путь MATLAB->C++/openMP->Python. Ваши примеры очень интересны, я давно хотел разобраться с cython, к сожалению, по нему очень мало документации для начинающих.

Пробовал запустить ваши примеры, но ничего не вышло. Было бы здорово, чтобы вы написали какую-то минимальную инструкцию как это запускать. С помощью build.bat я модуль собрал, но дальше он запускаться оказывается:
$python test_mp_python_ext.py 
Traceback (most recent call last):
  File "test_mp_python_ext.py", line 20, in <module>
    import mp_cython_ext as mpcy
  File "-py/tests/mp_cython.pyx", line 16, in init tests.mp_cython_ext (tests/mp_cython.c:16853)
    import os, sys
SystemError: Parent module '' not loaded, cannot perform relative import

Ещё было бы хорошо, если бы вы положили в репозиторий тестовые данные, а то повторить ваши измерения невозможно. Спасибо.
Вы случайно не Python 3 используете? Насколько я знаю, в Python 3 были изменения в синтаксисе относительных импортов. Я пока еще использую Python 2.7.6. На Python 3 перейти пока не могу из-за библиотеки трехмерной визуализации Mayavi2, она пока сидит на Python 2.

Тестовые данные генерируются при запуске теста на c++. Ниже краткая инструкция по запуску тестов.

Окружение:
Windows 8.1 x64
Python 2.7.6
MinGW x64-4.8.1-posix-seh-rev5

Настройка рабочей среды:
1. Scientific Python. Скачать и установить Anaconda 2.1.0 для Windows x64 с сайта continuum.io/downloads.
2. MinGW compiler. Скачать MinGW x64-4.8.1-posix-seh-rev5 и установить с сайта sourceforge.net/projects/mingwbuilds/. При установке важно точно указать архитектуру x86_64 и posix threads.
3. Добавить путь к компилятору в переменные среды. Дать команду gcc -v и убедиться, что используется именно установленный компилятор.

Запуск на C++:
1. Скачать ветку open-nuance с github.com/alec-kalinin/open-nuance
2. Перейти в каталог src-cpp и в файле SConstruct в переменной CXX указать путь к установленному компилятору.
3. Дать команду scons для сборки тестов. Если пакет scons не установлен в системе, дать команду conda install scons.
4. Запустить benchmark-omp.exe. При этом выполнится тест на C++ и сгенерируются данные для python тестов.

Запуск тестов:
1. Перейти в каталог src-py и дать команду build.bat
2. Перейти в каталог test и дать команды python test_mp_python.py, python test_mp_python_ext.py.
Спасибо, да у меня python3. Я попробую на втором это запустить.
Попробовал применить это дело в реальном проекте. Результат неутешителен. Внутри параллельного цикла openmp нельзя вызвать никакие питоновские функции. А как как сам scipy не предоставляет интерфейс cython, то никакие функции из этой библиотеки вызвать нельзя.

Например, вам хочется посчитать значения функции Бесселя нулевого порядка в точках 1, 2, 3… и так далее.
@cython.boundscheck(False)
@cython.wraparound(False)
def par(N):
    pR = np.empty(N)    
    
    cdef int i
    cdef int nN = N
    cdef double[:] R = pR

    with nogil, parallel():
        for i in prange(nN):
                R[i] = scipy.special.j0(i)
                
    return R


Это не скомпилится, так как scipy.special.j0 это чисто питоновская функция. Конечно, можно вытащить Бесселя, например, из GSL и импортировать. Но зачем тогда нужна вся эта питоновская машинерия, когда уж можно на чистом C написать. Так что приведенные выше примеры работают только для простейших операций, для чего-то более серьёзного это уже не годится. Так что проблему GIL не преодолели, а замели под ковер.
Да, согласен. Более-менее работающая идея преодоления GIL заключается в отказе от Python объектов и в переходе на уровень С. Да, при этом действительно нельзя использовать Python методы. Но я не согласен, что это такая уж сильная проблема. Вот ведь в чем дело.

Сами по себе NumPy и SciPy работают достаточно быстро и вполне параллельно, особенно в части линейной алгебры. Проблемы возникают тогда, когда нам приходится делать что-то нестандартное, написанное нами самостоятельно. И почти всегда в нашем самостоятельном коде можно выделить небольшую часть кода — узкое место, оптимизировав которое мы добьемся значительного прироста в производительности. Например, в моем коде это специальная процедура вычисления сингулярных интегралов.

Очень многие качественные Python пакеты работают именно по такому принципу. Хороший пример — это пакет конечноэлементных расчетов SfePy (http://sfepy.org/) на Python. Там лишь маленькая часть всего кода оптимизирована при помощи Cython, но и этого хватает для очень даже неплохой производительности.
Да и вообще, а какой есть выход? На C++ экспериментировать с математикой слишком долго и непроизводительно. На MATLAB мне очень не хвататет возможностей языка и окружения, слишком все завязано на базовую линейную алгебру, шаг в сторону слишком тяжелый. SciPy идеальный вариант, но есть GIL. Но все-таки всегда можно так построить архитектуру своего приложения, что выделить небольшой низкоуровненый кусок кода и оптимизировать его при помощи Cython. И это работающий вариант.
Спасибо за разъяснения. Но выходов на самом деле несколько. Можно смотреть в сторону pypy, но пока numpy для него не полностью готов. Можно смотреть на новые языки типа rust и julia, там нет GIL и нативные треды. Но что меня сильно держит на python это matplotlib, уж больно приятные графики получаются без усилий.

Сами по себе NumPy и SciPy работают достаточно быстро и вполне параллельно, особенно в части линейной алгебры.

А как это сделать-то? Вот есть у меня СЛАУ, пишу я a=numpy.linalg.solve(M, b). Считает всё-равно в один тред, так что приходится использовать multiprocessing, что ещё те костыли, ни лямбд, ни динамических объектов.

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

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

F(k) = \int f(x,k)dx


и нам надо посчитать значения F(k) для разных k. Идеальная параллельная задача, так как f(x,k) — чистая функция. И сделать это с помощью cython не представляется возможным, потому что:
1. Внутри nogil ни о каком scipy.integrate.quad не может быть и речи — надо писать свой интегратор или брать из GSL
2. f(x,k) должна быть полностью реализована на cython. А это можно быть значительный кусок программы.

Либо я не понимаю чего-то :)
А как это сделать-то? Вот есть у меня СЛАУ, пишу я a=numpy.linalg.solve(M, b). Считает всё-равно в один тред, так что приходится использовать multiprocessing, что ещё те костыли, ни лямбд, ни динамических объектов.

Скорее всего, у вас Numpy не слинкован с какой-нибудь быстрой BLAS библиотекой (OpenBLAS, Intel MKL, ...) и поэтому использует свою реализацию, которая значительно медленнее (и, видимо, не параллельная). Сейчас проверил — linalg.solve работает вполне себе параллельно и грузит все ядра (у меня MKL стоит).

2. f(x,k) должна быть полностью реализована на cython. А это можно быть значительный кусок программы.

Если f() большая и сложная, вполне возможно, что в ней есть пара небольших кусков, которые занимают большую часть времени — их и нужно оптимизировать на Cython. Остальная же часть функции может быть также написана «на cython», но без оптимизаций и указания типов, к примеру — т.е. почти на обычном питоне.
Чтобы параллельно исполнять f(x,k) она должна быть объявлена как nogil, так что «почти как на обычном питоне» уже не выйдет. Тем более cython похож на python, но не равен. Особенно большие различия в классах, просто так приписать к классу cdef даже не исправляя функции не получается, особенно если есть что-то сложное типа декораторов, например.
Sign up to leave a comment.

Articles