Pull to refresh

Python, scipy.weave и openMP — разгоняем код

Reading time 3 min
Views 13K
Здравствуйте %username%, данная статья посвящена проблеме увеличения скорости математических вычислений на основе языка python с использованием scipy.weave и openMP.

Многие могут задаться вопросом: «Зачем вообще использовать python для математических вычислений?», но мы не будем отвечать на «вечные» вопросы, как и не будем рассматривать множество других решений данной проблемы, таких как, например, psyco.

Инструменты


Как описано выше, наш инструмент — это библиотека scipy.weave, а также библиотека openMP.
scipy — набор библиотек для вычислений в прикладной математике и науке. openMP — открытый стандарт для распараллеливания программ на языках Си, Си++ и Фортран.

Установка пакетов


В Debian-подобных Linux системах необходимо выполнить:
apt-get install python-scipy
apt-get install libgomp1

Метод


Для увеличения скорости вычислений необходимо реализовать «узкую» часть python кода (обычно это цикл, в котором происходят какие-либо действия с матрицей) на C и добавить директивы openMP для распараллеливания.

Пример


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

Python реализация


В python c использование numpy данная задача, не учитывая различных подготовительных операций, как генерация матрицы и прочего, решает в пару строк кода:
  1. # цикл по строкам матрицы, где i - номер строки
  2. # с - целое число, randRow - случайный вектор
  3. for i in xrange(N):
  4.     matrix[i,:] -= c*randRow
  5.  
Генерация случайной матрицы x на y, в нашем случае x=y:
  1. # генерация случайной матрицы x на y
  2. # элементы матрицы - случайные числа от 0 до 99 включ.
  3. def randMat(x,y):  
  4.     randRaw = lambda a:  [ randint(0100) for i in xrange(0,a) ]
  5.     randConst = lambda x, y: [ randRaw(x) for e in xrange(0,y)]
  6.     return array(randConst(x, y))
  7.  

Реализация scipy.weave без openMP


scipy.weave — часть библиотеки scipy, которая позволяет использовать C/C++ код внутри кода python.
Происходит это следующим образом:
  1. # C код
  2. codeC = 
  3. """
  4. int i = 0;
  5.  
  6. for(i = 0; i < N*M; i++) {
  7.     matrix[0,i] = matrix[0,i] - (c * randRow[i%M]);
  8. }
  9. """
  10. weave.inline(codeC, ['matrix','c''randRow','N''M'], compiler = 'gcc')

т.е. сам С код хранится в виде multiline string, а переменные python кода, передаются в С списком, где элементы — одноименные текстовые константы. Также numpy массивы передаются в C не в виде матрицы, а в виде вектора, именно поэтому в коде один цикл, а не два.

Кстати, получившийся С код можно поискать в /tmp/%user%/python2x_intermediate/compiler_x

Реализация scipy.weave с openMP


Теперь уже к добавленной версии необходимо добавить директивы openMP и в вызове inline добавить недостающие параметры, а именно:
  1. # C и openMP код
  2. codeOpenMP = 
  3. """
  4. int i = 0;
  5.  
  6. omp_set_num_threads(2);
  7. #pragma omp parallel shared(matrix, randRow, c) private(i)
  8. {
  9. #pragma omp for 
  10. for(i = 0; i < N*M; i++) {
  11.     matrix[0,i] = matrix[0,i] - (c * randRow[i%M]);
  12. }
  13. }
  14. """
  15.  
  16. ...
  17.  
  18. weave.inline(codeOpenMP, ['matrix','c''randRow','N''M'],
  19.     extra_compile_args =['-O3 -fopenmp'],
  20.     compiler = 'gcc',
  21.     libraries=['gomp'],
  22.     headers=['<omp.h>'])
Полный исходный код со всеми реализациями можно скачать здесь

Сравнение результатов


Выше изложенный исходный код можно запустить и убедиться в том, что scipy.weave действительно даёт прирост в скорости:
Test on size: 100x100
	Pure python:  0.0725984573364
	Pure C: 0.303888320923
	C plus OpenMP: 0.109100341797
	Test - ok

Test on size: 1000x1000
	Pure python:  1.00839138031
	Pure C: 0.506997108459
	C plus OpenMP: 0.333213806152
	Test - ok

Test on size: 2000x2000
	Pure python:  3.24151515961
	Pure C: 2.10800170898
	C plus OpenMP: 1.17690563202
	Test - ok

Test on size: 3000x3000
	Pure python:  5.54490089417
	Pure C: 4.61800098419
	C plus OpenMP: 2.56960391998
	Test - ok

Литература


В написании кода использовались следующие ресурсы:
Tags:
Hubs:
+32
Comments 13
Comments Comments 13

Articles