Lei Yu
Lei Yu

Reputation: 367

Python - Matrix outer product

Given two matrices

A: m * r
B: n * r

I want to generate another matrix C: m * n, with each entry C_ij being a matrix calculated by the outer product of A_i and B_j.

For example,

A: [[1, 2],
    [3, 4]]

B: [[3, 1],
    [1, 2]]

gives

C: [[[3, 1],  [[1 ,2],
     [6, 2]],  [2 ,4]],
     [9, 3],  [[3, 6],
     [12,4]],  [4, 8]]]

I can do it using for loops, like

    for i in range (A.shape(0)):
      for j in range (B.shape(0)):
         C_ij = np.outer(A_i, B_j)

I wonder If there is a vectorised way of doing this calculation to speed it up?

Upvotes: 18

Views: 27250

Answers (5)

user3187724
user3187724

Reputation: 218

You can use

C = numpy.tensordot(A, B, axes=0)

tensordot does exactly what you want. The axes parameter is used to in addition perform sums over certain axes (for >2d tensors, the default value is 2 and it will sum away 2 axes of each of the arrays), but by setting it to 0 it simply doesn't reduce the dimensions and keeps the whole outer product.

Upvotes: 1

CypherX
CypherX

Reputation: 7353

Simple Solution with Numpy Array Broadcasting

Since, you want C_ij = A_i * B_j, this can be achieved simply by numpy broadcasting on element-wise-product of column-vector-A and row-vector-B, as shown below:

# import numpy as np
# A = [[1, 2], [3, 4]]
# B = [[3, 1], [1, 2]]
A, B = np.array(A), np.array(B)
C = A.reshape(-1,1) * B.reshape(1,-1)
# same as: 
# C = np.einsum('i,j->ij', A.flatten(), B.flatten())
print(C)

Output:

array([[ 3,  1,  1,  2],
       [ 6,  2,  2,  4],
       [ 9,  3,  3,  6],
       [12,  4,  4,  8]])

You could then get your desired four sub-matrices by using numpy.dsplit() or numpy.array_split() as follows:

np.dsplit(C.reshape(2, 2, 4), 2)
# same as:
# np.array_split(C.reshape(2,2,4), 2, axis=2)

Output:

[array([[[ 3,  1],
         [ 6,  2]],

        [[ 9,  3],
         [12,  4]]]), 
array([[[1, 2],
         [2, 4]],

        [[3, 6],
         [4, 8]]])]

Upvotes: 1

hpaulj
hpaulj

Reputation: 231540

The Einstein notation expresses this problem nicely

In [85]: np.einsum('ac,bd->abcd',A,B)
Out[85]: 
array([[[[ 3,  1],
         [ 6,  2]],

        [[ 1,  2],
         [ 2,  4]]],


       [[[ 9,  3],
         [12,  4]],

        [[ 3,  6],
         [ 4,  8]]]])

Upvotes: 19

user2357112
user2357112

Reputation: 281594

temp = numpy.multiply.outer(A, B)
C = numpy.swapaxes(temp, 1, 2)

NumPy ufuncs, such as multiply, have an outer method that almost does what you want. The following:

temp = numpy.multiply.outer(A, B)

produces a result such that temp[a, b, c, d] == A[a, b] * B[c, d]. You want C[a, b, c, d] == A[a, c] * B[b, d]. The swapaxes call rearranges temp to put it in the order you want.

Upvotes: 9

Roland Smith
Roland Smith

Reputation: 43533

Use numpy;

In [1]: import numpy as np

In [2]: A = np.array([[1, 2], [3, 4]])

In [3]: B = np.array([[3, 1], [1, 2]])

In [4]: C = np.outer(A, B)

In [5]: C
Out[5]: 
array([[ 3,  1,  1,  2],
       [ 6,  2,  2,  4],
       [ 9,  3,  3,  6],
       [12,  4,  4,  8]])

Once you have the desired result, you can use numpy.reshape() to mold it in almost any shape you want;

In [6]: C.reshape([4,2,2])
Out[6]: 
array([[[ 3,  1],
        [ 1,  2]],

       [[ 6,  2],
        [ 2,  4]],

       [[ 9,  3],
        [ 3,  6]],

       [[12,  4],
        [ 4,  8]]])

Upvotes: 0

Related Questions