Is there a more efficient way of exponentiating a multidimensional matrix? This is how I am currently doing it:

import numpy as np 
from scipy import linalg

A = np.random.normal(0,0.01, (3, 2,2))
P = np.zeros((3,2,2))
for i in range(3):
     P[i] = linalg.expm( np.multiply(A[i], 0.5))

Note that my real array is much larger and any effort to speed up the process is greatly appreciated.


1 Answer 1


Dig through the source code for expm, which leads ultimately to a rather extensive bit of code in scipy.sparse.linalg.matfuncs.py. This uses a well known Pade approximation. It is specifically designed to work with on a square nxn matrix. So generalizing to a (m, n, n) matrix without iterating on m does not work.

Years ago I sped this up by using f2py to calculate an analytical (n,3,3) inverse. But I recall a simpler expm implementation. The current one appears to have a lot of helper functions and special cases, presumably hoping to gain speed.

If your matrix is always (m, 2,2) there might be a special case.

This might be better as a Stackoverflow question, since it has less to do with programming style or code organization, and more with calculation alternatives in the numpy and scipy libraries. But a SO search for expm and numpy/scipy turns up very few questions. The closest suggests looking at the source code just as I suggested.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.