[SciPy-User] Better way to multiply block diagonal object array?
Thu Oct 1 18:52:14 CDT 2009
On Thu, Oct 1, 2009 at 7:46 PM, Skipper Seabold <firstname.lastname@example.org> wrote:
> 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?
Err, I must have missed something when I was doing this earlier
XX = scipy.linalg.block_diag(x_1, x_2, x_3)
XTX2 = np.dot(XX.T,XX)
Seems to work fine.
More information about the SciPy-User