Who would have thought what a long story this would be?:
https://github.com/JuliaLang/julia/issues/4774
Turns out the linear algebra world and the multilinear algebra world don't
play nicely together, and the dictates of performance can really splist
the needs of the human and the computer.
some day I hope to write a paper on this subject
On Saturday, December 6, 2014 8:40:59 AM UTC-5, Andrew wrote:
>
> (Edit: I answered my question somewhat and posted it at the end)
>
> I recently ran into a bug in my code in which Julia gave an error when
> computing the matrix inverse of what I thought was a 1x1 matrix. I looked
> into it, and it looks like my sequence of matrix multiplications resulted
> in a 5x5 matrix being converted into a 1 dimensional array. I'm not really
> sure why this happens. Here's an example:
>
> A = rand(2,2)
> B = rand(2,1)
> B'*A*B
>
> 1x1 Array{Float64,2}:
> 0.836342
>
>
> This is fine. But,
> R = [0;1]
> R'*A*R
>
> 1-element Array{Float64,1}:
> 0.921516
>
>
> I think the problem is that B is explicitly defined as a 2 dimensional
> array. On the other hand, R is 2-element Array{Int64,1}: , a 1
> dimensional array. This is a problem because
> (R'*A*R)^(-1)
>
> DomainError
>
>
> It looks like, in Julia, I have to be very careful about combining
> 1-dimensional vectors and 2 dimensional matrices. Is there a better way to
> do this? Also, when constructing a matrix(or vector) by hand, as in R =
> [0;1], can I force it to be 2 dimensional?
>
>
> Edit: I just sort of answered my own question, but I'll post this anyway
> in case anyone else has this question, or if there are any comments.
> If I define R as
> R = [0 1]'
>
> 2x1 Array{Int64,2}:
> 0
> 1
>
>
> then there is no issue.
>
