Multidimensional matrix exponentiation 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.
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.