Ohm
Ohm

Reputation: 2442

Multiply each block of a block matrix with different coefficient in NumPy

Suppose I have a (4*n)*(4*n) block matrix, and I would like to multiply each block of (n*n) with a different coefficient from a corresponding (4*4) matrix - what is the way to do it in NumPy?

For example, I have:

>>> mat
matrix([[1, 1, 1, 1],
        [1, 1, 1, 1],
        [1, 1, 1, 1],
        [1, 1, 1, 1]])
>>> x
array([[1, 2],
       [3, 4]])

And the requested result should be like this:

>>> result
array([[1, 1, 2, 2],
       [1, 1, 2, 2],
       [3, 3, 4, 4],
       [3, 3, 4, 4]])

Upvotes: 3

Views: 917

Answers (2)

Divakar
Divakar

Reputation: 221614

Here's a vectorized approach using NumPy broadcasting -

m1,n1 = x.shape
m2,n2 = mat.shape
out = (mat.reshape(m1,m2//m1,n1,n2//n1)*x[:,None,:,None]).reshape(m2,n2)

Sample run -

In [41]: mat
Out[41]: 
array([[8, 8, 7, 2, 3, 4],
       [2, 4, 7, 5, 4, 8],
       [7, 2, 4, 5, 6, 5],
       [4, 3, 3, 6, 5, 3],
       [7, 3, 5, 8, 7, 7],
       [2, 5, 2, 4, 2, 7]])

In [42]: x
Out[42]: 
array([[1, 2],
       [3, 4]])

In [43]: out
Out[43]: 
array([[ 8,  8,  7,  4,  6,  8],
       [ 2,  4,  7, 10,  8, 16],
       [ 7,  2,  4, 10, 12, 10],
       [12,  9,  9, 24, 20, 12],
       [21,  9, 15, 32, 28, 28],
       [ 6, 15,  6, 16,  8, 28]])

I think np.kron as suggested in @Tim Fuchs's post would be the safest option here without fiddling around with the size of input arrays. For performance, here's some timings on a decent sized input arrays comparing kron and broadcasting-based approaches -

In [56]: mat = np.random.randint(2,9,(100,100))

In [57]: x = np.random.randint(2,9,(50,50))

In [58]: n = 2

In [59]: m1,n1 = x.shape # Ignoring timings from these as negligible
    ...: m2,n2 = mat.shape
    ...: 

In [60]: %timeit np.multiply(mat, np.kron(x, np.ones((n,n))))
1000 loops, best of 3: 312 µs per loop

In [61]: %timeit (mat.reshape(m1,m2//m1,n1,n2//n1)*x[:,None,:,None]).reshape(m2,n2)
10000 loops, best of 3: 83.7 µs per loop

Upvotes: 2

Tim Fuchs
Tim Fuchs

Reputation: 1191

One way would be to construct a (4*n)*(4*n) block matrix from x:

result = np.multiply(mat, np.kron(x, np.ones((n,n))))

Upvotes: 3

Related Questions