[JuliaLang/LinearAlgebra.jl] Bug in UMFPack Interface (#149)
Christoph Ortner
I had some unexplained \ behaviour that is, according to Ivar Nesje, a bug in /base/linalg/umfpack.jl#L18
Here is a short code snippet that works ok. I get the correct solution from it. Primarily the last line is relevant. ( I can produce the complete code on request. )
# construct a finite element grid
X, T = square(10)
# assemble the stiffness matrix and rhs (for - \Delta u = 1)
A, F = simple_stiffmat(X, T)
# find the interior nodes
Ifree = find( (minimum(X,1) .> 0) & (maximum(X, 1) .< 1) )
# solve the linear system and plot the solution
U = zeros(size(X,2))
U[Ifree] = cholfact(A[Ifree,Ifree]) \ F[Ifree]
If I replace the last line with
U[Ifree] = A[Ifree,Ifree] \ F[Ifree]
then I get the following error message:
umferror has no method matching umferror(::Int32)
while loading In[112], in expression starting on line 3
in umfpack_symbolic! at linalg/umfpack.jl:155
in umfpack_numeric! at linalg/umfpack.jl:176
in lufact at linalg/umfpack.jl:118
in A_ldiv_B! at linalg/sparse.jl:210
in \ at linalg/generic.jl:233
Same happens with lufact. I also tried copying B = copy(A[Ifree,Ifree]) first. Same result.
I have used UMFPack successfully in other cases. E.g.,
U = (A+speye(size(X,2)) \ F
works fine.
I copy Ivar Nesje's comment from Google Groups:
This is at least one bug.
In /base/linalg/umfpack.jl#L18, we have a convenience function to lookup the status code, seems to be restricted to Int64, even though it is called with Int32 on some occations. If you compile yourself, or know how to regenerate the system image, it is pretty easy to change Int to Integer on that line. A real fix for this would require someone to actually check the return types of the C functions called, and include some more tests to ensure that all code paths are followed in the test suite.
—
Reply to this email directly, view it on GitHub, or unsubscribe.
You are receiving this because you were mentioned.