Re: About diagonal matrices

[Jaroslav Hajek]

On Sat, Feb 21, 2009 at 10:46 PM, John W. Eaton <address@hidden> wrote:
> On 21-Feb-2009, Jaroslav Hajek wrote:
>
> | Mathematically, there's no such thing as a NaN.
>
> Substitute "numerically" for "mathematically" then, where
> "numerically" means according to accepted conventions of the IEEE
> standard for floating point arithmetic.
>

But the current behavior does not violate IEEE standard anywhere. It
just no longer makes certain multiplications that used to be performed
before, and assumes the result instead. That's the distinction between
assumed zero and numerical zero.


> | In my opinion, the current behaviour is the most useful and sane and
> | is OK if it's documented.
> | If you intend to work on "fixing" these "issues", the please make the
> | "mathematically wrong" behavior optional by a configuration variable,
> | because I will most definitely want it to be set.
>
> Yes, I understand that it is convenient for many uses for diagonal
> and sparse matrices to have the properties you want.  But I'm also
> don't like having things like
>
>  full_matrix == diag_matrix
>
> yet
>
>  diag_matrix * scalar != full_matrix * scalar
>
> for some values of scalar.
>
> jwe
>

Why? This behaviour has been around for ages for sparse matrices, and
nobody complained.
It's just the distinction between an assumed zero and numerical zero.
It's *standard* in numerical software.
It is usually both numerically (or practically) superior and more
effective, you just need to be aware of it in
certain cases.

Let me show you another example:
a = [3 1 1; 0 2 1; 0 0 4];
 a \ [NaN; 1; 1]
ans =

       NaN
   0.37500
   0.25000

Of course this is again "mathematically wrong" by David's definition.
The reason, obviously, is that here we use optimized treatment of the
triangular matrix.
Do you really want to insert checks for such stuff everywhere just to
deliver more slowly a result that is not numerically superior in any
reasonable sense? I think it's bad, but still OK if there's a runtime
or compile-time option that can switch off these "improvements". But
if there is not even such an option, and if I'm not allowed to add it,
then Octave will be a much less useful tool for me.
I know I can fork my own version if I'm not satisfied, but I think
this is how good projects die.
So, either of these options would be fine with me:

1. Leave the assumed zeros treatment. Document it (I'm working on this)
2. Implement the additional checks, but provide an option to turn them
off to get better performance and avoid memory explosions.

-- 
RNDr. Jaroslav Hajek
computing expert
Aeronautical Research and Test Institute (VZLU)
Prague, Czech Republic
url: www.highegg.matfyz.cz