3 views (last 30 days)

w=0:(pi/100):pi;

w = w([17 66 93]);

x = exp(1i*w);

x./x

ans =

1.0000e+00 + 4.8645e-17i 1.0000e+00 - 4.9461e-17i 1.0000e+00 - 2.6884e-17i

I imagine this result has something to do with how Matlab implements division of complex numbers and I realize that the error is less than eps(1), But I thought it was an interesting result and it would never have occurred to me that this could happen. I checked the doc page for rdivide, but it doesn't say anything about its algorithm when both arguments are complex. Just curious if anyone can explain why this result happens (it only happened for these three particular values of the original w vector).

Matt J
on 14 Jun 2020

Edited: Matt J
on 14 Jun 2020

It depends what algorithm is used to perform complex division x./y. For example, if it is implemented as x./abs(y).^2.*conj(y), then there will be floating point errors introduced in these operations, even when x=y,

>> x./abs(x).^2.*conj(x) - 1

ans =

1.0e-15 *

0 0 0.2220

If it had been done, however, by conversion to complex exponentials, I don't think we should have seen residual errors,

>> [r,t]=deal(abs(x),angle(x));

>> (r./r).*exp(1i*(t-t))-1

ans =

0 0 0

James Tursa
on 16 Jun 2020

These are not supid questions IMO.

You would need to look at IEEE floating point standards for stuff like x/x = 1. I think IEEE requires that the calculation result be correct as if the calculation is done in infinite precision and rounded according to the lsb rounding rule in effect (e.g., typically "round to even"), so that would mean x/x = 1 would be required to be true for all possible real finite floating point bit patterns (except 0). It would not be hard to write a MATLAB program to test every possible single and double bit pattern to see if it is true. (But I would be shocked to find out that any such bit pattern failed this test ... this is implemented at the chip level)

For language standards regarding complex arithmetic, I am unaware of any particular requirements (e.g., z/z = 1 or avoiding overflow & underflow when possible). I think that Java has some floating point replicability requirements in their standards but I don't know the details.

Opportunities for recent engineering grads.

Apply Today
## 6 Comments

## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_897552

⋮## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_897552

## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_897654

⋮## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_897654

## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_897786

⋮## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_897786

## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_897852

⋮## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_897852

## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_898158

⋮## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_898158

## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_899268

⋮## Direct link to this comment

https://ww2.mathworks.cn/matlabcentral/answers/548010-why-is-a-number-divided-by-itself-not-equal-to-unity#comment_899268

Sign in to comment.