|
From: | Daniel J Sebald |
Subject: | Re: About diagonal matrices |
Date: | Fri, 20 Feb 2009 22:41:49 -0600 |
User-agent: | Mozilla/5.0 (X11; U; Linux i686; en-US; rv:1.7.3) Gecko/20041020 |
John W. Eaton wrote:
On 20-Feb-2009, Jaroslav Hajek wrote: | On Fri, Feb 20, 2009 at 5:45 PM, John W. Eaton <address@hidden> wrote:| | > But I think there are also some bugs. Shouldn't the following return| > full matrices with the zero elements replaced by NaN (or -0, in the | > case of dividing by -Inf)? | > | > diag ([1,2,3]) / 0 | > diag ([1,2,3]) / NaN | > diag ([1,2,3]) / -Inf | > diag ([1,2,3]) * NaN | > diag ([1,2,3]) * Inf | >| | I think it's better in the current manner. I don't like the idea that| the memory can suddenly explode just because the multiplier happened | to be Inf. Right now, scalar * diag gives invariantly diag. This is | somewhat analogous to how sparse matrices behave. I see octave:3> speye (3)*Inf ans = Compressed Column Sparse (rows = 3, cols = 3, nnz = 3 [33%]) ) (1, 1) -> Inf (2, 2) -> Inf (3, 3) -> Inf octave:4> speye (3)/0 warning: division by zero ans = Inf NaN NaN NaN Inf NaN NaN NaN Inf So there seems to be some disagreement here, even within Octave. I'm not sure what is best here, but I don't like it that we produce numerically incorrect or incompatible answers just because of a storage scheme.
It depends on how one wants to interpret 1/0, 1/Inf I suppose. (And perhaps the solution is to have some type of option specifying how to treat it.) If I'm a user who wants to work with just rational numbers (notice I didn't say "real"), then 1/0 is NaN, i.e., undefined. Also, the "number" Inf then means the CPU has reached the limits of its number system, i.e., 1/Inf is meaningless (or NaN) in this context as well. If I'm a user who wants to think in terms of the extended real number system, then limits are of interest, i.e., Inf is the limit as series x_i becomes positively unbounded (but not negatively unbounded). 1/0 is fraught with problems though because it doesn't indicate how the limit is approached. That is, how do we know that 1/0 is +Inf and not -Inf? 1/+0 could mean approach zero from the left, so 1/+0 = +Inf. Likewise 1/-0 = -Inf. Octave agrees with that train of thought, e.g., octave:20> 1/(+0) warning: division by zero ans = Inf octave:21> 1/(-0) warning: division by zero ans = -Inf But here is where the trouble shows up: octave:22> 1/(+0-0) warning: division by zero ans = Inf octave:23> 1/(-0+0) warning: division by zero ans = Inf The problem is trying to treat 0 as both a number in the number system and the notion of a one-sided limit. (Perhaps what some mathematician should have done long ago was to define the extended real number system as R union {+Inf, -Inf, +Zed, -Zed}. Anyway...) And how about this one, irrespective of the previous comments? octave:42> sqrt(-0) ans = -0 Shouldn't that be complex 0 + 0i? That is, octave:44> sqrt(-2) ans = 0.00000 + 1.41421i octave:45> sqrt(-1) ans = 0 + 1i octave:46> sqrt(-0.5) ans = 0.00000 + 0.70711i ... sqrt(lim x-> 0 from left) = 0 + 0i. Dan
[Prev in Thread] | Current Thread | [Next in Thread] |