Python, numpy, einsum умножают стек матриц

По соображениям производительности,

Мне любопытно, есть ли способ умножить стек стека матриц . У меня есть 4-мерный массив (500, 201, 2, 2). Его в основном 500-разрядный стек (201,2,2) матриц, где для каждого из 500 я хочу умножить смежные матрицы с использованием einsum и получить другую (201,2,2) матрицу.

Я только делаю умножение матрицы на [2×2] матрицах в конце. Поскольку мое объяснение уже направлено на рельсы, я просто покажу, что я делаю сейчас, а также эквивалент «уменьшить» и почему его не полезно (потому что его скорость равна вычислительной мощности). Предпочтительно, это будет несколько однострочный, но я не знаю, что это такое, или даже если это возможно.

Код:

Arr = rand(500,201,2,2) def loopMult(Arr): ArrMult = Arr[0] for i in range(1,len(Arr)): ArrMult = np.einsum('fij,fjk->fik', ArrMult, Arr[i]) return ArrMult def myeinsum(A1, A2): return np.einsum('fij,fjk->fik', A1, A2) A1 = loopMult(Arr) A2 = reduce(myeinsum, Arr) print np.all(A1 == A2) print shape(A1); print shape(A2) %timeit loopMult(Arr) %timeit reduce(myeinsum, Arr) 

Возвращает:

 True (201, 2, 2) (201, 2, 2) 10 loops, best of 3: 34.8 ms per loop 10 loops, best of 3: 35.2 ms per loop 

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

  • Найти разность фаз между двумя (ангармоническими) волнами
  • Выполнение программ Java через Python
  • Инвертировать кортежи в списке кортежей (Python)
  • Как выполнить скрипт Python на стороне сервера, используя jQuery?
  • python как «пакетный» скрипт (т. е. запуск команд из python)
  • Как использовать библиотеку FastFM Python (машины факторизации) для задач рекомендации?
  • возможно ли слияние пустых нечетких совпадений с python pandas?
  • Ошибка атрибута при попытке запустить быстрый запуск API Gmail в Python
  • One Solution collect form web for “Python, numpy, einsum умножают стек матриц”

    Я не думаю, что это можно сделать эффективно с помощью numpy (решение cumprod было элегантным, хотя). Это такая ситуация, когда я буду использовать f2py . Это самый простой способ вызова более быстрого языка, который я знаю, и требует только одного дополнительного файла.

    fortran.f90:

     subroutine multimul(a, b) implicit none real(8), intent(in) :: a(:,:,:,:) real(8), intent(out) :: b(size(a,1),size(a,2),size(a,3)) real(8) :: work(size(a,1),size(a,2)) integer i, j, k, l, m !$omp parallel do private(work,i,j) do i = 1, size(b,3) b(:,:,i) = a(:,:,i,size(a,4)) do j = size(a,4)-1, 1, -1 work = matmul(b(:,:,i),a(:,:,i,j)) b(:,:,i) = work end do end do end subroutine 

    Скомпилируйте с f2py -c -m fortran fortran.f90 (или F90FLAGS="-fopenmp" f2py -c -m fortran fortran.f90 -lgomp чтобы включить ускорение OpenMP). Затем вы будете использовать его в своем скрипте как

     import numpy as np, fmuls Arr = np.random.standard_normal([500,201,2,2]) def loopMult(Arr): ArrMult = Arr[0] for i in range(1,len(Arr)): ArrMult = np.einsum('fij,fjk->fik', ArrMult, Arr[i]) return ArrMult def myeinsum(A1, A2): return np.einsum('fij,fjk->fik', A1, A2) A1 = loopMult(Arr) A2 = reduce(myeinsum, Arr) A3 = fmuls.multimul(Arr.T).T print np.allclose(A1,A2) print np.allclose(A1,A3) %timeit loopMult(Arr) %timeit reduce(myeinsum, Arr) %timeit fmuls.multimul(Arr.T).T 

    Какие результаты

     True True 10 loops, best of 3: 48.4 ms per loop 10 loops, best of 3: 48.8 ms per loop 100 loops, best of 3: 5.82 ms per loop 

    Так что это фактор 8 ускорения. Причиной для всех транспозиций является то, что f2py неявно переносит все массивы, и нам нужно их вручную транспонировать, чтобы сказать, что наш фортран-код ожидает, что вещи будут транспонированы. Это позволяет избежать операции копирования. Стоимость заключается в том, что каждая из наших матриц 2×2 транспонируется, поэтому, чтобы избежать неправильной операции, мы должны выполнить обратное преобразование.

    Большее ускорение, чем 8, должно быть возможным – я не тратил времени, пытаясь оптимизировать это.

    Interesting Posts

    Python Imaging: проблемы с YCbCr

    Это способ проверки полей модели Django?

    Как элегантно проверить существование объекта / экземпляра / переменной и одновременно назначить его переменной, если она существует в python?

    Может ли sys.argv обрабатывать дополнительные аргументы?

    Гистограмма выравнивания изображений в оттенках серого с помощью NumPy

    Как можно красиво печатать таблицы ASCII с Python?

    Перемещение по виджетам в макете PyQt

    Как установить «авто» для верхнего предела, но сохранить фиксированный нижний предел с помощью matplotlib.pyplot

    Как автоматически регистрировать фид Atom с помощью Python?

    Вычисление стандартного отклонения в потоке

    Каковы правила определения контекста списка в классе Python?

    Python ctypes: объект файла Python <-> C FILE *

    Как реализовать свойство () с динамическим именем (в python)

    Python – Загрузить изображения из Google Поиск изображений?

    Использование векторов freebase с gensim

    Python - лучший язык программирования в мире.