[racket] matrix-solve and approximation errors

From: Neil Toronto (neil.toronto at gmail.com)
Date: Wed Apr 16 11:33:00 EDT 2014

FWIW, Jens Axel's "very close" needs quantifying: the entries are less 
than 2^-52 away from comprising a non-invertible matrix. :)

But approximation error is a big deal, and quantifying and responding to 
it currently has little support in `math/matrix'. We need two things to 
start with. One requires engineering, the other a design decision.

The engineering is implementing a function that computes condition 
numbers, which you can use to estimate error in the results of 
`matrix-solve' and other inverse-like functions. IIRC, it requires 
computing eigendecomposition or singular value decomposition, both of 
which are useful for other things (e.g. principal component analysis, 
different solving strategies) and are on the radar.

The design decision is to determine how much the solvers should automate 
detecting and responding to floating-point error. I have no idea what 
extent popular linear algebra libraries go to, nor whether we should try 
to do better.

Neil ⊥

On 04/16/2014 05:34 AM, Laurent wrote:
> Interesting; even more when using rounded numbers returns a correct
> solution!
>
> To protect against such errors, maybe `matrix-solve` could run a
> post-check to verify that M×X=B up to some epsilon, depending on an
> optional argument?
> Or maybe just mention to do this check systematically in the docs?
>
>
> On Wed, Apr 16, 2014 at 1:31 PM, Jens Axel Søgaard
> <jensaxel at soegaard.net <mailto:jensaxel at soegaard.net>> wrote:
>
>     2014-04-16 12:18 GMT+02:00 Jens Axel Søgaard <jensaxel at soegaard.net
>     <mailto:jensaxel at soegaard.net>>:
>      > Hi Laurent,
>      >
>      > I think the underlying problem is that the matrix is *very* close to
>      > an invertible one:
>
>     non-invertible
>
>
>      >> (matrix-determinant
>      >      (matrix [[ 1     0     9/10  1]
>      >                  [ 0     1     1/10  1]
>      >                  [ 9/10  1/10 82/100 1]
>      >                  [ 1     1      1   0]]))
>      > 0
>
>
>
>
> ____________________
>    Racket Users list:
>    http://lists.racket-lang.org/users
>


Posted on the users mailing list.