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

Skipper Seabold jsseabold@gmail....
Thu Oct 1 18:52:14 CDT 2009


On Thu, Oct 1, 2009 at 7:46 PM, Skipper Seabold <jsseabold@gmail.com> 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?
>
> Skipper
>

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.

Skipper


More information about the SciPy-User mailing list