Numpy multiplication of a matrix by an array of scalars without for loop

Viewed 946

Let's say I have the following code:

# define a 3x3 array
A = np.array([[-2, 1, 1], [1,-2,1], [1,1,-2]])
# define a 1D array of scalars to multiply A with
test = np.arange(10)

for i in range(len(test)):
    matrix = scipy.linalg.expm(A*test[i])

I want to see if there is a way to do this multiplication without using a for loop. I am not trying to multiply the two arrays using matrix multiplication. I am treating the test array as a bank of scalar values that I want to multiply A with one by one. There has got to be some kind of sneaky numpy way of doing this. Any ideas?

3 Answers

Here's an answer to the question "how to scale a matrix by multiple scalars without a for-loop". This disregards the matrix exponentiation, since I don't think there's a way to do that without the for-loop, but the question doesn't seem to be asking about the exponentiation step.

You can use 3D arrays to use numpy's vectorized multiplication.

import numpy as np

A = np.array([[[-2, 1, 1], [1,-2,1], [1,1,-2]]])
test = np.arange(10).reshape(-1,1,1)

result = test*A

scipy.linalg.expm cannot be applied with a not squared matrix:

ValueError: expected a square matrix

so I think the easiest thing is to do:

list_result = list(map(lambda i: scipy.linalg.expm(A*i), range(10)))

An then if you want only one array:

np.concatenate(list_result) #for 2D

or

np.stack(list_result) #for 3D

You can do the trick building a special matrix by stacking the identity matrix multiplied by the values you want to use.

I will use a 3x3 matrix as an example, but the idea is the same for an arbitrary matrix:

Imagine we want to multiply the matrix:

A = np.array([
   [1,0,1],
   [0,1,0],
   [0,0,1] 
])

By [2, 3, and 5]. Without loops.

We can build a special matrix contatenating the identity to do this:

[ 2 * np.eye(3), 3 * np.eye(3), 5 * np.eye(3)]

[ 0 0 2 ] [ 0 0 3] [ 0 0 5 ]
[ 0 2 0 ] [ 0 3 0] [ 0 5 0 ]
[ 2 0 0 ] [ 3 0 0] [ 5 0 0 ]

We can do this with a list comprehension:


multiply_values = [2,3,5]

special_matrix = np.concatenate( [ x * eye(3) for x in multiply_values ], axis=1)

# special_matrix is:
# array([[2., 0., 0., 3., 0., 0., 5., 0., 0.],
#        [0., 2., 0., 0., 3., 0., 0., 5., 0.],
#        [0., 0., 2., 0., 0., 3., 0., 0., 5.]])

Now we can multiply the matrix in one step:

result = np.dot( , special_matrix)

# result is
array([[2., 0., 0., 3., 0., 0., 5., 0., 0.],
       [0., 2., 0., 0., 3., 0., 0., 5., 0.],
       [0., 0., 2., 0., 0., 3., 0., 0., 5.]])

Probably this large wide matrix is not what we need.

We can get the partial results slicing the resulting product:

result[:,0:3]

# array([[2., 0., 2.],
#       [0., 2., 0.],
#       [0., 0., 2.]])

We can process the slices with another comprehension:

[ special_matrix[:,x:x+3] for x in [0,3,6] ]

# or in a more general way

[ special_matrix[:,x:x+A.shape[0]] for x in list(range(0,A.shape[0]**2,A.shape[0] ]
Related