10

I hate having to ask this because I assume the answer must be simple, but I cannot for the life of me seem to track down the source. While trying to rewrite a function I ran across this problem:

a = -j
x = real(a)  
y = imag(a)
y/x

Which spits out Inf, unexpectedly for me. However...

a = 0
b = -1
b/a

returns -Inf, like I would expect. Inquiring further, a == x, b == y. Clearly that isn't true however. I finally tracked down the problem to this after a lot of frustration. If the original input for a is instead 0-j (vs. -j) then there is no problem.

Both real(-j) and real(0-j) return zero and test as zero, but obviously seem to retain some metadata relating to their origin that I absolutely cannot discover. What precisely am I missing here? It will feel downright wrong if I have to solve this with something like if (x == 0) then x = 0;

Stunt
  • 103
  • 4

2 Answers2

10

Not metadata, just the sign bit of the double precision float.

>> a = 0-j;
>> b = -j;
>> ra = real(a)
ra =
     0
>> rb = real(b)
rb =
     0
>> ra==0
ans =
     1
>> isequal(ra,rb)
ans =
     1

Looks the same so far. However, the difference is that with b, we set the sign bit for both the real and imaginary parts when we do -j = -complex(0,1) vs. 0-j = complex(0,-1) (see Creating Complex Numbers). Looking deeper with typecast, which does no conversion of the underlying data:

>> dec2bin(typecast(ra,'uint64'),64)
ans =
0000000000000000000000000000000000000000000000000000000000000000
>> dec2bin(typecast(rb,'uint64'),64)
ans =
1000000000000000000000000000000000000000000000000000000000000000

That 1 is bit 63 (of 0) in the IEEE 754 double precision floating point representation:

enter image description here

Voila! -0 exists in MATLAB too!

chappjc
  • 30,359
  • 6
  • 75
  • 132
  • 2
    Another way to put it is that `-j` is computed as `-(0+j) = -0-j` while `0-j` is `+0-j`. Comparing `-0` to `+0` expectedly yields equality, but when you're reaching the realm of infinity through zero division, this sign actually does make a difference. – scenia Mar 01 '14 at 08:24
  • 1
    @scenia That's one way to view it. Another is `-j = -complex(0,1)` while `0-j = complex(0,-1)` (non-functional constructor for [complex](http://www.mathworks.com/help/matlab/matlab_prog/complex-numbers.html#f2-59290) number). Although, I'm surprised, perhaps naïvely, that `isequal` declares that the two are equal. I would expect that it would perform a bit-wise AND on _all_ 64-bits. – chappjc Mar 01 '14 at 09:11
  • Or rather an XNOR, if such an instruction exists. – chappjc Mar 01 '14 at 18:12
  • 1
    Yup, this makes total sense. Didn't even cross my mind that it would a signed zero and the sign being distributed to it during the construction of the complex argument. Thank you so much, this was driving me insane. Now that I know that matlab handles signed zeros in this way it wont be mystifying behavior. Will feel a little bit less insane to have (if x == 0, x = 0) going on. – Stunt Mar 01 '14 at 18:21
  • @Stunt Great, I'm glad that made sense! This is so odd, that I wonder if this is expected behavior, or if there is is a bug somewhere. The equality tests for one should at least detect this difference. Anyway, feel free to accept (checkbox indicating the best answer for you) and upvote my answer. Thanks! – chappjc Mar 01 '14 at 18:45
  • Yeah, it seems weird to me. I mean, ignoring signed zero makes sense, but if you're going to have it matter in some cases then you should be able to see it and detect a difference. Maddening. Incredibly glad the reason occurred to you, and glad I decided to register an account here to ask. – Stunt Mar 02 '14 at 04:19
  • I think this answer does not contain enough the analysis aspect of -0/+0 which can be reached. See here http://stackoverflow.com/a/7946995/54964 – Léo Léopold Hertz 준영 Aug 05 '16 at 13:56
  • 1
    @Masi Interesting paper. Thanks for the link. – chappjc Aug 05 '16 at 15:19
1

When using IEEE 754 floating point numbers there is a convention to have a number approaching zero that cannot be represented by the smallest possible float called an underflow where the precision of the number is being lost with each step below the smallest possible float. Some operating systems will consider the underflow to be equal to zero.

I was surprised to be testing some software and find that a threshold test of zero actually went below zero almost as far as the smallest possible negative float.

Perhaps this is why your getting a negative infinity instead of a divide by zero error which I am assuming is the problem your referring to.

Jason
  • 286
  • 2
  • 9
  • Your guess may be right... I tried the following things: `1. a=0-j, imag(a)/real(a) --> gives -Inf`, `2. a=-j, imag(a)/real(a) --> gives +Inf`, `3. a=-1e-1000, b=-1 --> b/a --> gives Inf` – Autonomous Mar 01 '14 at 03:04
  • Does Matlab complain about a = 1/0 or does it simply assign a +Inf to a ? – Jason Mar 01 '14 at 03:08
  • 2
    @Parag It's not that it's really really small. It's zero. The problem is that with floats, there is a _negative zero_. – chappjc Mar 01 '14 at 04:55