[SciPy-User] Better way to multiply block diagonal object array?

Skipper Seabold jsseabold@gmail....
Thu Oct 1 18:46:24 CDT 2009


This is somewhat related to the block diagonal matrix thread from a
few months back.  I'm wondering if there's a better (faster, cleaner)
way to do what I'm trying to do.

Say I have a "block diagonal" matrix X = diag(x_1, x_2, x_3) where the
x_# are 2d arrays that all have the same number of rows (but don't
have to be square).  I would like to do something like dot(X.T,X), but
putting these all into a single array (eg., scipy.linalg.block_diag)
doesn't do what I want.  I've been messing with kron and the sparse
functions, but I haven't hit on what I'm looking for yet.

x_1 = np.arange(4).reshape(2,-1)
x_2 = np.arange(4,8).reshape(2,-1)
x_3 = np.arange(8,16).reshape(2,-1)
X = np.zeros((3,3), dtype=object)
X[0,0] = x_1
X[1,1] = x_2
X[2,2] = x_3

# The resulting matrix will be 3 x 3 with position i,j = np.dot(x_i.T,x_j)

XTX = np.zeros((3,3), dtype=object)
for i in range(3):
    for j in range(3):
        XTX[i,j] = np.dot(X[i,i].T,X[j,j])

I could even work with this if I could take a view on it, so that I
have an 8x8 float.  Any ideas?

Skipper


More information about the SciPy-User mailing list