Лучший способ улучшить суммирование элементов комплексной матрицы в Numpy по MATLAB

Я пытаюсь переписать приведенный ниже код MATLAB на Python и обнаружил, что мой код Python (2,7 секунды) медленнее, чем MATLAB (1,2 секунды). Я пробовал много разных способов, включая модуль numba, но пока не повезло. Как я могу сделать код Python быстрее?

Код MATLAB:

szA=[1024,1280]; HfszA=[512,640];
[aPx,aPy]=meshgrid(-HfszA(2):HfszA(2)-1,-HfszA(1):HfszA(1)-1);
img=randi(255,1024,1280);
fx=rand(); fy=rand();
tic
for i=1:20
    F=abs(sum(sum(img.*exp(-1i*2*pi*(fx*aPx+fy*aPy)))));
end
toc

Код Python:

import numpy as np
import time

szA=[1024,1280]; HfszA=[512,640]
aPx,aPy=np.meshgrid(np.arange(-HfszA[1],HfszA[1]),np.arange(-HfszA[0],HfszA[0]))
img=np.array(np.random.randint(256,size=(1024,1280)))

fx=np.random.rand()
fy=np.random.rand()

start = time.time()
for i in range(20):
    F=abs(np.sum(img*np.exp(-1j*2*np.pi*(fx*aPx+fy*aPy))))
end = time.time()
print("Elapsed (after compilation) = %s" % (end - start))
print(F)

Вы перебираете переменную i, но не используете эту переменную в цикле. Вы просто вычисляете F 20 раз, перезаписывая предыдущий результат. Этот цикл предназначен только для определения времени выполнения вычислений?

Cris Luengo 26.09.2018 18:35

Я задал вышеупомянутый вопрос, потому что (1) ваш вопрос может содержать некоторые пояснения относительно того, что происходит, и (2) вы должны использовать timeit в MATLAB и timeit.timeit в Python для определения времени.

Cris Luengo 26.09.2018 18:42
Почему в Python есть оператор "pass"?
Почему в Python есть оператор "pass"?
Оператор pass в Python - это простая концепция, которую могут быстро освоить даже новички без опыта программирования.
Некоторые методы, о которых вы не знали, что они существуют в Python
Некоторые методы, о которых вы не знали, что они существуют в Python
Python - самый известный и самый простой в изучении язык в наши дни. Имея широкий спектр применения в области машинного обучения, Data Science,...
Основы Python Часть I
Основы Python Часть I
Вы когда-нибудь задумывались, почему в программах на Python вы видите приведенный ниже код?
LeetCode - 1579. Удаление максимального числа ребер для сохранения полной проходимости графа
LeetCode - 1579. Удаление максимального числа ребер для сохранения полной проходимости графа
Алиса и Боб имеют неориентированный граф из n узлов и трех типов ребер:
Оптимизация кода с помощью тернарного оператора Python
Оптимизация кода с помощью тернарного оператора Python
И последнее, что мы хотели бы показать вам, прежде чем двигаться дальше, это
Советы по эффективной веб-разработке с помощью Python
Советы по эффективной веб-разработке с помощью Python
Как веб-разработчик, Python может стать мощным инструментом для создания эффективных и масштабируемых веб-приложений.
1
4
70
2
Перейти к ответу Данный вопрос помечен как решенный

Ответы 2

Я считаю, что вы можете попробовать даже более простую задачу

MATLAB

function test
    szA=[1024,1280]; HfszA=[512,640];
    [aPx,aPy]=meshgrid(-HfszA(2):HfszA(2)-1,-HfszA(1):HfszA(1)-1);
    fx=1.0; fy=2.0;
    tic
    for i=1:2000
        F=sum(sum(fx*aPx+fy*aPy));
    end
    toc
    disp(F)

выходы

Elapsed time is 9.566274 seconds.

-1966080

Python

import numpy as np
import time

szA=[1024,1280]; HfszA=[512,640]
aPx,aPy=np.meshgrid(np.arange(-HfszA[1],HfszA[1]),np.arange(-HfszA[0],HfszA[0]))

fx=1.0
fy=2.0

start = time.time()
for i in range(2000):
    F = np.sum(np.sum(fx*aPx+fy*aPy, axis=0))
end = time.time()
print("Elapsed (after compilation) = %s" % (end - start))
print(F)

выходы

Elapsed (after compilation) = 33.3840000629

-1966080.0

Я считаю, что разница между ними двоякая:

  • Matlab использует OpenMP при выполнении многих операций с массивами
  • использует набор инструкций AVX.

Причина, по которой я пришел к такому выводу, заключается в том, что я смотрю на использование ЦП на моем Core i3 с четырьмя ядрами (вы можете проверить количество ядер на своем ЦП), со сценарием python - 30%, с Matlab - 100%.

Что касается набора инструкций AVX, я просто однажды сравнил операцию MATLAB matmul с операцией Eigen (http://eigen.tuxfamily.org/index.php?title=Main_Page), и для обеспечения соответствия производительности мне пришлось скомпилировать Eigen с -openmp и -axAVX.

Чтобы наконец ответить на ваш вопрос, я не думаю, что вы можете сделать код Python быстрее, если вы не можете скомпилировать numpy базовую библиотеку с openmp, директивой AVX.

Вот учебник https://docs.scipy.org/doc/scipy/reference/building/linux.html

Удачи, дайте нам знать, как все прошло.

Ответ принят как подходящий

Пожалуйста, всегда публикуйте то, что вы пробовали до сих пор. Что касается вашей версии Numba, я думаю, что вы сделали что-то, что привело к плохой работе.

Пример

import numpy as np
import numba as nb
import time

@nb.njit(fastmath=True)
def your_function(fx,fy,aPx,aPy,img):
  pi=np.pi
  sum=0.
  for i in range(aPx.shape[0]):
    for j in range(aPx.shape[1]):
      sum+=img[i,j]*np.exp(-1j*2*pi*(fx*aPx[i,j]+fy*aPy[i,j]))
  return np.abs(sum)

@nb.njit(fastmath=True,parallel=True)
def your_function_p(fx,fy,aPx,aPy,img):
  pi=np.pi
  sum=0.
  for i in nb.prange(aPx.shape[0]):
    for j in range(aPx.shape[1]):
      sum+=img[i,j]*np.exp(-1j*2*pi*(fx*aPx[i,j]+fy*aPy[i,j]))
  return np.abs(sum)

#The function gets compiled at the first call
#you may also use cache=True, which only works in single threaded code
F=your_function(fx,fy,aPx,aPy,img)
start = time.time()
for i in range(20):
    F=your_function(fx,fy,aPx,aPy,img)

end = time.time()
print("Elapsed (after compilation) = %s" % (end - start))
print(F)

F=your_function_p(fx,fy,aPx,aPy,img)
start = time.time()
for i in range(20):
    F=your_function_p(fx,fy,aPx,aPy,img)

end = time.time()
print("Elapsed (after compilation) = %s" % (end - start))
print(F)

Тайминги (4C / 8T)

your_version: 2.45s
Numba single threaded: 0.17s
Numba parallel: 0.07s

Большое спасибо, max9111. Ваш ответ не только очень помогает мне в работе, но и исправляет мои возможные неправильные представления о самом языке Python. Еще раз спасибо.

S. Lee 27.09.2018 10:21

Другие вопросы по теме