search

Daniel-Diaz / matrix

A Haskell native implementation of matrices and their operations.
BSD 3-Clause "New" or "Revised" License
35 stars 31 forks source link

Matrix Inversion Should be fixed for nearly singular matricies #45

Closed octopuscabbage closed 6 years ago

octopuscabbage commented 7 years ago

Still need tests, Waiting on test case from @LATBauerdick Copied from his issue:

Hi -- hope this is the right way to report bugs for the Data.Matrix package!

There are two issues showing up when matrix inversion runs into problems (e.g. for nearly singular matrices):

in function ref:

goodRow = case listToMaybe (filter significantRow [1..ncols mtx]) of
            Nothing -> error "Attempt to invert a non-invertible matrix"
            Just x -> x

ncol needs to be nrows, current code crashes with access error instead of returning error.

as ref is called from inverse (and rref), which promises to return Either String (Matrix a), calling error in case of non-invertible matrices is not a very nice thing to do. ref should be modified to also return Either String (Matrix a), so inverse can return Left "can't invert". It's important that user code can handle singular or numerically unstable matrices themselves.

Thanks for your consideration (and thanks for providing Data.Matrix, which is otherwise a joy to use).

Best regards, LatB
octopuscabbage commented 7 years ago

Should fix #45

LATBauerdick commented 7 years ago

Here is a (trivial) test case:

λ> let m = Data.Matrix.fromList 3 3 [0,0,0,0,0,0,0,0,0::Double] λ> Data.Matrix.inverse m Right ( *** Exception: getElem: Trying to get the (4,1) element from a 3x6 matrix. CallStack (from HasCallStack): error, called at ./Data/Matrix.hs:427:16 in matrix-0.3.5.0-Vu1MfM9DVkGI45qY2HPrU:Data.Matrix

octopuscabbage commented 7 years ago

Alright, test should be there.