умножение трехмерной матрицы в numpy

Кажется, вчера я задал неправильный вопрос. На самом деле я хочу умножить две матрицы 2x2xN A и B, чтобы

C[:,:,i] = dot(A[:,:,i], B[:,:,i])

Например, если у меня есть матрица

A = np.arange(12).reshape(2, 2, 3)

Как я могу получить C = A x A с определением, описанным выше? Есть ли встроенная функция для этого?


Кроме того, если я умножу A (shape 2x2xN) на B (shape 2x2x1, instead of N), я хочу получить

C[:,:,i] = dot(A[:,:,i], B[:,:,1])

person LWZ    schedule 20.03.2013    source источник


Ответы (1)


Попробуйте использовать numpy.einsum, у него есть небольшая кривая обучения, но он должен дать вам то, что вы хотите. Вот пример для начала.

import numpy as np

A = np.random.random((2, 2, 3))
B = np.random.random((2, 2, 3))

C1 = np.empty((2, 2, 3))
for i in range(3):
    C1[:, :, i] = np.dot(A[:, :, i], B[:, :, i])

C2 = np.einsum('ijn,jkn->ikn', A, B)
np.allclose(C1, C2)
person Bi Rico    schedule 20.03.2013
comment
Это потрясающе. Это работает для обоих случаев в вопросе. Это даже работает, когда я меняю местами A и B (разной формы, во втором вопросе). - person LWZ; 21.03.2013
comment
+1 Когда я впервые увидел einsum, у него был тот неотличимый от магии привкус, отличная функция и сверхбыстрый. - person Jaime; 21.03.2013
comment
@Bi Rico, @Jaime, ладно, кроме крутизны, может кто-нибудь немного объяснить, что происходит с 'ijn,jkn->ikn'? Мне было трудно связать документ numpy с этим конкретным случаем. - person LWZ; 21.03.2013
comment
@LWZ Во-первых, вам нужно понять, что скалярное произведение матрицы будет записано как ('ij,jk->ik') (см. mathworld.wolfram. com/MatrixMultiplication.html). Затем вы хотите сделать то же самое для всех матриц в последнем измерении, что и делает добавление n в конце. Возможно, это поможет вам увидеть это более четко, если вы замените переменную цикла в ответе выше с i на n. - person jorgeca; 21.03.2013
comment
@jorgeca Спасибо! Теперь я понимаю эту часть. Но что касается последней части моего вопроса, когда последнее измерение двух входных данных отличается, N против 1, почему это также работает? Не вызовет ли это путаницу, каково последнее измерение вывода? - person LWZ; 21.03.2013
comment
@LWZ Это действительно интересная проблема. Я не знаю, почему это все еще работает. Если бы мы написали ('ij...,jk...->ik...'), я бы не удивился, что это тоже сработало. - person jorgeca; 21.03.2013
comment
Причина, по которой это работает, когда последнее измерение отличается, заключается в том, что numpy передает, scipy.org/EricsBroadcastingDoc, массивы, чтобы иметь совместимые формы. Вот сообщение об ошибке, которое вы получаете, когда последние измерения равны 3 и 4: ValueError: operands could not be broadcast together with remapped shapes [original->remapped]: (2,2,3)->(2,newaxis,3,2) (2,2,4)->(2,4,2). - person Bi Rico; 22.03.2013
comment
Но einsum не разрешает вещание, если вы не добавите многоточие (...). В документах четко сказано einsum does not allow broadcasting by default. To enable and control broadcasting, use an ellipsis. Это неправильно или это ошибка? - person jorgeca; 23.03.2013