I am attempting to learn Julia and am experimenting with the matrix
function inv to calculate the inverse of some matrices. Some are singular
and some are nonsingular.
One case (below) gives me the wrong result but it is a simple example so I
must be misusing Julia (I found a similar issue with large numbers in the
basic factorial function in the documentation returning 0 which came down
to me not realizing the overflow behavior immediately)
These attempts works as expected
julia> A = [1 2 3; 0 2 2; 1 2 3]; inv(A)
ERROR: SingularException(3)
in inv at linalg/lu.jl:149
in inv at linalg/dense.jl:328
julia> A = [1 1 1; 0 2 3; 5 5 1]; inv(A)
3x3 Array{Float64,2}:
1.625 -0.5 -0.125
-1.875 0.5 0.375
1.25 0.0 -0.25
julia> inv(A) * A
3x3 Array{Float64,2}:
1.0 0.0 0.0
0.0 1.0 0.0
0.0 0.0 1.0
This does not
jjulia> A = [2 1 -1; 1 -2 -3; -3 -1 2]; inv(A)
3x3 Array{Float64,2}:
-4.5036e15 -6.43371e14 -3.21686e15
4.5036e15 6.43371e14 3.21686e15
-4.5036e15 -6.43371e14 -3.21686e15
julia> A * inv(A)
3x3 Array{Float64,2}:
0.0 -0.125 -0.5
0.0 0.75 0.0
0.0 0.0 0.0
Double checking manually, in julia with rref, and with wolframalpha.com
shows the expected results on this singular matrix
julia> A = [2 1 -1; 1 -2 -3; -3 -1 2]; rref([A eye(3)])
3x6 Array{Float64,2}:
1.0 0.0 -1.0 0.0 0.142857 -0.285714
0.0 1.0 1.0 0.0 -0.428571 -0.142857
* 0.0 0.0 0.0 * 1.0 0.142857 0.714286
inv {{2,1,-1}, {1,-2,-3}, {-3,-1,2}}
http://www.wolframalpha.com/input/?i=inv+%7B%7B2%2C1%2C-1%7D%2C+%7B1%2C-2%2C-3%7D%2C+%7B-3%2C-1%2C2%7D%7D
Thanks for any help,
Ben
