[racket] matrix-solve and approximation errors

Hi Laurent,

I think the underlying problem is that the matrix is *very* close to
an invertible one:

> (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

/Jens Axel




2014-04-16 12:02 GMT+02:00 Laurent <laurent.orseau at gmail.com>:
> Forgot to mention that with the true value of n=0.82, this of course returns
> the correct solution:
> (let ([n 0.82 #;(+ (* .9 .9)(* .1 .1))])
>
>   (matrix-solve
>    (matrix [[ 1  0 .9 1]
>             [ 0  1 .1 1]
>             [.9 .1  n 1]
>             [ 1  1  1 0]])
>    (col-matrix [0 0 0 1])))
> ; -> (array #[#[0.38] #[0.4866666666666667] #[0.13333333333333333] #[-0.5]])
>
>
> On Wed, Apr 16, 2014 at 11:10 AM, Laurent <laurent.orseau at gmail.com> wrote:
>>
>> I've just been bitten by a bad case of floating-point error with
>> `matrix-solve` (and a bad CPU that has some floating-point issues):
>>
>> (let ([n 0.8200000000000001 #;(+ (* .9 .9)(* .1 .1))])
>>   (matrix-solve
>>    (matrix [[ 1  0 .9 1]
>>             [ 0  1 .1 1]
>>             [.9 .1  n 1]
>>             [ 1  1  1 0]])
>>    (col-matrix [0 0 0 1])))
>> ; -> (array #[#[0.0] #[0.5] #[0.0] #[-0.5]])
>>
>> But clearly here M×X≠B, as is easily seen on the last row.
>> I've seen other situations where the approximation leads to an approximate
>> solution (which is okay of course), but this is the first case I see where
>> the result is completely off.
>>
>> I have no idea if anything can be done about it, though (apart from
>> throwing my computer through the window and buy a better one).
>>
>> Laurent
>
>
>
> ____________________
>   Racket Users list:
>   http://lists.racket-lang.org/users
>



-- 
--
Jens Axel Søgaard