Python
Comments 38
0
Я не знаю C++ вообще, но честное слово, код на C++ выглядит проще всех этих извращений.
А если вместо векторов использовать массивы, то думаю, что C++ будет работать ещё быстрее.
0
Тут вот в чем дело. Python очень красивый язык с большим количеством красивых решений. Но вот иногда в рамках красивого и чистого Python проекта возникает вычислительная задача, где нужно порой нужно пройтись по 2D или 3D циклу и что-то вычислить. Кроме цикла сложно что-то придумать, и чистый Python для этого как раз работает плохо.

И тут как раз возникает ряд возможных решений с той иной степенью красоты:
  • отказаться от итераций и работать только с матричными операциями (NumPy)
  • свести задачу к набору стандартных базовых операций (SciPy)
  • спуститься к уровню C/C++ (и только здесь возникает С++ код)
  • добавить информацию о типах для JIT компилятора

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

И по всей статье ни одного (даже плохого примера) с параллельностью. Есть и multiprocessing, и concurrent futures…
Да и тема «чистого» питона не раскрыта, можно генератор написать, можно map функцию использовать.
В общем посыл не понятен, что Вы хотели этой статьей показать, что питон медленный и не «параллельный»? Так мы это и так знаем)
0
Хм… А вы точно читали статью, например, фразу:
"… заменить 'range' на 'prange' для параллельного выполнения внешнего цикла"?
0
Я читал статью, и имел ввиду «встроенную» параллельность, на встроенных модулях которые я перечислил.
0
А можно ли на чистом Python писать на самом деле быстрый и параллельный код?

Вот с этого места.
не увидел «чистого» питона в статье, уж простите.
Возможно у нас просто разные взгляды на «чистый питон», тогда я уточню. Для меня это пайтон, который ставится из коробки и используется только функционал его встроенных модулей.
0
Ну хорошо. А давайте попробуем быть чуть ближе к практике? Вы пишите на чистом Python очень красивый web сервис. Но вот проблема, вам в рамках это сервиса нужно решить задачу, приведенную в статье. Это очень простая задача. И достаточно распространенная.

Что для вас наиболее чистое решение этой задачи?
+2
Если касаться именно веб сервисов, то нужно разделять быстрые задачи и медленные. Вторые обычно кладутся с очередь задач и там исполняются, по окончанию пользователь может забрать результат.
Если же требуется решение в реальном времени, то задачи с недетерминированным временем в рантайм выносить нельзя. Сегодня у Вас 5000 точек, завтра будет 10к, после -1млн. Очевидно что человек, выполнивший такую операцию повесит сервер для всех.
Что касается реализации, о которой вы спросили, то лично я бы таки выносил это в очередь, и делал на concurrent futures, с очередью, что касается чистого пайтона. Если бы скорость меня не устроила, я бы написал расширение DLL/SO на С++/С и просто бы его использовал и пайтона. Возможно это не самый простой способ, ноя не люблю полагаться на сторонний код, если реализация подобного занимает немного времени.
Второй момент, я не придирался к вашей реализации, возможно в жизни она действительно хороша и эффективна, я даже не сомневаюсь в этом. Единственное что меня смутило, что Вы обозвали питон не параллельным из коробки, указали примеры на сторонних библиотеках, при этом не привели решения, который таки дает «коробочный» питон.
+3
В этой задаче немаловажным должен быть порядок хранения матриц. Почему для q в питоновском коде используется транспонированная форма?
0

Думаю, с транспонированием q что-то не так. Зачем там транспонирование — не понятно. cdist принимает на вход данные в виде такой матрицы:


[(x1, y1, z1, ..., n1),
 (x2, y2, z2, ..., n2),
 ...
 (xK, yK, zK, ..., nK)
]

Поэтому, если у автора данные хранятся именно в такой форме и если транспонировать один из наборов данных, то будет ошибка:


ValueError: XA and XB must have the same number of columns (i.e. feature dimension.)
+1
Да, тут вы правы, спасибо за замечание! Тест вырос из NumPy кода, а там я использую такое хранение для того, чтобы избежать итераций и сразу получить матрицу за счет использования numpy broadcasing rules, например:

import numpy as np

N = 5000

p = np.random.rand(N, 1)
q = np.random.rand(N, 1)

# no broadcasting
print((p - q).shape)

>>> (5000, 1)

# with broadcasting
print((p - q.T).shape)

>>> (5000, 5000)


Сейчас я попробую сделать тест без транспонирования и погляжу на результаты.
+1
После тестов время осталось примерно такое же. Но вообще, вы правы! Порядок очень важен, поскольку от него зависит попадание в кеш процессора.
0
А еще я бы в с++ коде попробовал бы использовать Point4D вместо Point3D. В этом случае, мне кажется, выравнивание может еще ускорить код. Насколько я знаю, стандартные библиотеки для работы с матрицами стараются выравнивать строки.
+1
Не ожидал такого результата от Numba! Поразительно.

А не пробовали запустить тест на PyPy? Интересно посмотреть, что он может продемонстрировать в сравнении с другими участниками тестирования.

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

Спасибо за тест!
0

Как насчёт того, чтобы попробовать сравнить еще и Intel Distribution for Python?

+1
В таблице не хватает только теста на фортране с параллелизацией )))
+1
У меня на интеловском компиляторе код ниже(переписал исходники автора) с O3 оптимизацией выходит 4.43E-2 секунды в среднем. Навряд ли у кого-то с 8 ядрами выйдет сильно быстрее.
    USE IFPORT
    REAL(4) p(5000,3),q(5000,3),R(5000,5000)
    integer i,o
    CALL SRAND(100)
    do j=1,1000
       N = 5000
       do i=1,N
           do o=1,3
               p(i,o)=RAND()
               q(i,o)=RAND()
           end do
       end do
        CALL CPU_time(t1)
        DO CONCURRENT  (i=1:size(p,1))
            DO CONCURRENT  (o=1:size(q,1))
                rx=p(i,1)-q(o,1)
                ry=p(i,2)-q(o,2)
                rz=p(i,3)-q(o,3)
                R(i,o)=1 / (1 + sqrt(rx * rx + ry * ry + rz * rz))
            end do
        end do
        CALL CPU_time(t2)
        print *, t1, t2, t2-t1
        dt=dt+(t2-t1)
        print *, R(1,1) !чтобы компилятор не халтурил и выполнял всю работу
    end do
    print *, 'Avg time:', dt/1000.
    read(*,*)
0
Не придраться ради, а просто из спортивного интереса:
1) real(8) будет много шустрее работать
2) print обязательно выносить из-под цикла
Можно еще поиграть с опциями, я как-то «оптимизировал» код методом подбора параметров компиляции с 80+ секунд до 29…
0
  • real(8) — плавающая точка с двойной точностью, такие операции выполняются медленнее. У меня вышло 96 мс.
  • принт в данном случае выносить из цикла не надо, так как затраты на вывод информации не учитываются (считается время исполнения кода между двумя вызовами cpu_time)

А оптимизацией я заниматься не буду, так как оптимизируя опции компилятора вкупе с использованием старого кода, досконально не разбираясь в каждом привнесенном или унесенном эффекте, легко что-нибудь сломать (например, изменив опции qsave и qzero можно получить совсем неадекватные результаты, если в коде все переменные «инициализируются» implicit'ом), а в моем случае скорость работы кода не критична, но все должно считаться правильно)
0
Про print принимается, в бенчмарк он не попадает, это точно)
А вот про замедление на real(8) это неожиданно, у меня он самый быстрый…
0
Он не может быть самым быстрым, так как всегда медленнее real 4. А если взять real(16), то разница в скорости будет уже на порядки по сравнению с real(8).
Надо смотреть условия, код, настройки компилятора, при которых real(8) самый быстрый, там точно что-то не так и это значит, что не реал8 быстрый, а реал4 слишком медленный
0
Если SSE задействовано то реал(8) будет в 2 раза медленнее реал(4) так как обрабатывается в 2 раза больше данных за такт + еще больше вероятность попадания в кеш, так что больше 2х раз.
0
Чистый С на одном ядре выдаёт 127мс без трюков с оптимизацией и -О3 флагом. Ирония в том, что количество расчётов можно уменьшить в ~два раза (N(N-1)/2), но время вычисления только возрастает, где-то процентов на 80 =\
Upd. все расчёты велись для double, float выдаёт 041мс
0
Чёрт, я про корень забыл, с ним время сразу скачет до Плеяд. Тогда получаются совсем другие результаты: float — 494мс без оптимизации, 250 — с ней самой, как и ожидалось. Смена float на double увеличивает эти значения на 10-15мс.
0
новые, озвученные вами цифры характерны для одноядерных расчетов, вы видимо использовали чистый чистый си на одном ядре) У фортрана в таких условиях аналогичная производительность была.
0
Да, чистый С, одно ядро :) Распараллелить такую задачу — не проблема, а вот если что-то более изощрённое, вроде barnes-hut или cell-linked list, то там беда.
0
На мой взгляд не хватает теста на Matlab. Все-таки его можно в данной области считать стандартом. К тому-же он достаточно быстрый. Ну вообще было-бы интересно.
0
Тут обновилась Numba, реализация на которой свелась к красивому декоратору, а производительность значительно выросла.
0

Я конечно тоже считаю numpy неотъемлемой частью питона и не буду городить списки списков для двумерных массивов чисел, но всё-таки код на чистом питоне будет выглядеть как-то так:


def get_R_pure_python(p, q):
    R = [[0 for x in range(len(q))] for y in range(len(p))]

    for i in range(len(p)):
        for j in range(len(q)):
            rx = p[i][0] - q[j][0]
            ry = p[i][1] - q[j][1]
            rz = p[i][2] - q[j][2]
            R[i][j] = 1 / (1 + math.sqrt(rx * rx + ry * ry + rz * rz))

    return R

Кроме того, этот код быстрее чем numpy без broadcasting и его уже намного проще подсунуть pypy. В итоге на i5-2500 я получил следующие результаты:


| numpy 1.13.1 (без broadcasting) | 42532 мс | 477.89x |
| python 3.5.2 | 22247 мс | 249.97x |
| python 2.7.12 | 17773 мс | 199.70x |
| pypy 5.1.2 | 621 мс | 6.98x |
| numpy 1.13.1 (с broadcasting) | 511 мс | 5.74x |
| scipy 0.19.1 | 280 мс | 3.15x |
| numba 0.35.0rc1 | 89 мс | 1x |


Для python 2.7.12 функцию range я заменил на xrange там где это надо.
P.S. У кого-нибудь есть возможность поправить комментарий, чтобы заработал Markdown, потому что у меня, к сожалению, не работает.

0
Не знал, что Markdown начинает работать только через какое-то время.
0
Начнем тест скорости с кода на чистом Python:

    R = np.empty((p.shape[0], q.shape[1]))

это разве не Numpy методы?
т.е. это уже не совсем чистый python?
/занудамодеоф
Only those users with full accounts are able to leave comments., please.