Why does numpy.angle(0+0j) answer pi instead of zero

Question

I am creating an array of complex numbers and compute the angle.

When I set the magnitude to zero, I expect the angle to always be zero. However, I get pi if the phase is in the range of pi/2 to pi.

Here is my example code:

import numpy as np
import matplotlib.pyplot as plt

def main():
    # s = np.full(360, 1)    # Use this instead of np.zeros to get expected results
    s = np.zeros(360)
    f = 1 / s.size
    t = np.arange(s.size)
    s = s * np.exp(1j * 2*np.pi * f * t)

    fig, axs = plt.subplots(2, 2)

    axs[0,0].set_title("Magnitude")
    axs[0,0].grid()
    axs[0,0].plot(t, abs(s))

    axs[1,0].set_title("Phase")
    axs[1,0].grid()
    axs[1,0].plot(t, np.angle(s))
    # axs[1,0].plot(t, np.arctan2(np.imag(s), np.real(s)))

    axs[0,1].set_title("Real")
    axs[0,1].grid()
    axs[0,1].plot(t, np.real(s))

    axs[1,1].set_title("Imag")
    axs[1,1].grid()
    axs[1,1].plot(t, np.imag(s))

    plt.show()

if __name__ == "__main__":
    main()

Here is the output of the test program

I have also tried arctan2 and get the same behaviour.

Answer 1

Why does numpy.angle(0+0j) answer pi instead of zero

I did check numpy.angle docs, they claim

Although the angle of the complex number 0 is undefined, numpy.angle(0) returns the value 0.

so I did checked that

import numpy as np
assert np.angle(0+0j) == 0

and it does work as described in docs, I have used NumPy version 1.18.1, please check your version of NumPy. Maybe you are using one such ancient that described extension to numpy.angle was not implemented?

Answer 2

As @harold points out, not all elements are the same:

In [154]:     s = np.zeros(360)
     ...:     f = 1 / s.size
     ...:     t = np.arange(s.size)
     ...:     s = s * np.exp(1j * 2*np.pi * f * t)
     ...: 
In [155]: s
Out[155]: 
array([ 0.+0.j,  0.+0.j,  0.+0.j,  0.+0.j,  0.+0.j,  0.+0.j,  0.+0.j,
        0.+0.j,  0.+0.j,  0.+0.j,  0.+0.j,  0.+0.j,  0.+0.j,  0.+0.j,
       ...
       -0.+0.j, -0.+0.j, -0.+0.j, -0.+0.j, -0.+0.j, -0.+0.j, -0.+0.j,
       -0.+0.j, -0.+0.j, -0.+0.j, -0.+0.j, -0.+0.j, -0.+0.j, -0.+0.j,
        ...])
In [156]: np.angle(s)
Out[156]: 
array([ 0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
        0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
        0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
        ...
        0.        ,  3.14159265,  3.14159265,  3.14159265,  3.14159265,
        3.14159265,  3.14159265,  3.14159265,  3.14159265,  3.14159265,
        3.14159265,  3.14159265,  3.14159265,  3.14159265,  3.14159265,
        ...])

You had a question about the results of np.angle, so the logical debugging step is to examine the argument to that function. Next step is to figure out why some elements are -0.0+0.0j.

Answer 3

Answering to my own question:

The root cause seems to be the behaviour of math.atan2, as in the following code:

print(math.atan2(0.0, 0.0))
print(math.atan2(0.0, -0.0))
print(math.atan2(-0.0, 0.0))
print(math.atan2(-0.0, -0.0))

which prints on my system:

0.0
3.141592653589793
-0.0
-3.141592653589793

It is the same results as if the real part was non-zero, so that seems to be the rationale:

print(math.atan2(0.0, 1))
print(math.atan2(0.0, -1))
print(math.atan2(-0.0, 1))
print(math.atan2(-0.0, -1))

So coming back to numpy, the documentation of numpy.angle says that:

Although the angle of the complex number 0 is undefined, numpy.angle(0) returns the value 0.

This seems to be incorrect if the real part of the complex number is -0.0.