[Numpy-discussion] Toward release 1.0 of NumPy

[Charles R Harris]

On 4/13/06, Tim Hochberg <tim.hochberg at cox.net> wrote:
>
> Charles R Harris wrote:
>
> > Tim,
> >
> > On 4/13/06, *Tim Hochberg* <tim.hochberg at cox.net
> > <mailto:tim.hochberg at cox.net>> wrote:
> >
> >     Alan G Isaac wrote:
> >
> >     >On Thu, 13 Apr 2006, Charles R Harris apparently wrote:
> >     >
> >     >
> >     >>The Kronecker product (aka Tensor product) of two
> >     >>matrices isn't a matrix.
> >     >>
> >     >>
> >     >
> >     >That is an unusual way to describe things in
> >     >the world of econometrics.  Here is a more
> >     >common way:
> >     >http://planetmath.org/encyclopedia/KroneckerProduct.html
> >     <http://planetmath.org/encyclopedia/KroneckerProduct.html>
> >     >I share Sven's expectation.
> >     >
> >     >
> >     mathworld also agrees with you. As does the documentation (as best
> >     as I
> >     can tell) and the actual output of kron. I think Charles must be
> >     thinking of the tensor product instead.
> >
> >
> > It *is* the tensor product, A \tensor B, but it is not the most
> > general tensor with four indices just as a bivector is not the most
> > general tensor with two indices. Numerically, kron chooses to
> > represent the tensor product of two vector spaces a, b with dimensions
> > n,m respectively as the direct sum of n copies of b, and the  tensor
> > product of two operators takes the given form. More generally, the B
> > matrix in each spot could be replaced with an arbitrary matrix of the
> > correct dimensions and you would recover the general tensor with four
> > indices.
> >
> > Anyway, it sounds like you are proposing that the tensor (outer)
> > product of two matrices be reshaped to run over two indices. It seems
> > that likewise the tensor (outer) product of two vectors should be
> > reshaped to run over one index ( i.e. flat). That would do the trick.
>
> I'm not proposing anything. I don't care at all what kron does. I just
> want to fix the return type if that's feasible so that people stop
> complaining about it. As far as I can tell, kron already returns a
> flattened tensor product of some sort. I believe the general tensor
> product that you are talking about is already covered by multiply.outer,
> but I'm not sure so correct me if I'm wrong. Here's what kron does as
> present:
>
> >>> a
> array([[1, 1],
>        [1, 1]])
> >>> kron(a,a) # => 4x4 matrix
> array([[1, 1, 1, 1],
>        [1, 1, 1, 1],
>        [1, 1, 1, 1],
>        [1, 1, 1, 1]])


Good at first look. Lets see a simpler version... Nevermind, seems numpy
isn't working on this machine (X86_64, fc5 64 bit) at the moment, maybe I
need to check out a clean version.

>>> kron(a,a[0]) => 8x1
> array([1, 1, 1, 1, 1, 1, 1, 1])


Looks broken. a[0] should be an operator (matrix), so either it should be
(2,1) or (1,2). In the first case, the return should have shape (4,2), in
the latter (2,4). Should probably raise an error as the result strikes me as
ambiguous. But I have to admit I am not sure what the point of this
particular construction is.

>>> kron(a[0], a[0])
> Traceback (most recent call last):
>   File "<stdin>", line 1, in ?
>   File "C:\Python24\Lib\site-packages\numpy\lib\shape_base.py", line
> 577, in kron
>     result = concatenate(concatenate(o, axis=1), axis=1)
> ValueError: 0-d arrays can't be concatenated


See above. this could be (1,4) or (4,1), depending.

>>>  b.shape
> (2, 2, 2)
> >>> kron(b,b).shape
> (4, 4, 2, 2)


I think this is doing transpose(outer(b,b), axis=(0,2,1,3)) and reshaping
the first 4 indices into 2. Again, I am not sure what the point is for these
operators. Now another way to get all this functionality is to have a
contraction function or method with a list of axis. For instance, consider
the matrices A(i,j) and B(k,l) operating on x(j) and y(l) like A(i,j)x(j)
and B(k,l)y(l), then the outer product of all of these is

A(i,j)B(k,l)x(j)y(l)

with the summation convention on the indices j and l. The result should be
the same as kron(A,B)*kron(x,y) up to a permutation of rows and columes. It
is just a question of which basis is used and how the elements are indexed.

So, it looks like the 2d x 2d product obeys Alan's definition. The other
> products are probably all broken.
>
>
Chuck
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20060413/4983be21/attachment.html>