Search code examples
pythonmatplotlibplot3droot

How to plot color plot or density plot for complex numbers?


I have a res function that takes as input complex value of a certain omega variable and returns the absolute value of the residue (also complex) of a continued fraction.

I'd like to test the res function for various omega values and see which of those values yields a residual value closest to zero. So I would see which omega corresponds, which would be, in a sense, the root of the equation for the residue of the continued function. Therefore, I generated a list of omegas as follows:

N = 400

omega = [complex(np.random.uniform(0.001, 15), np.random.uniform(0, 2)) for i in range(N)]

And I imposed each omega[i] value on the res function. Truncating the lists, for the omega list

omega = [(8.89186304186803+1.8580361935398448j), (4.250554058118386+0.34692636983137204j),(13.77975011058525+0.4067498913305867j), (11.848491447761512+1.1353493824887473j), (5.368197310760776+1.603262808639006j), (7.953834937525846+0.31087387982821735j),
 (12.589975089023985+1.9010925629632527j),(0.8533967396680111+0.9083830823882479j),
 (11.966071152073017+1.7535358294965886j), (6.723597895069045+0.5943732982185757j)]

I got the following the absolute value of the corresponding residuals

list_res = [2.1663149855057022, 12.202467603015915, 13.105388388645999, 3.524702335556132, 5.007101216080126, 2.1989779961023457, 4.6753139694285775, 1.82871868615581, 3.8247209107872857, 2.977690850044372]

I want to generate a color plot or density plot with the omega values and the respective absolute values of the residual. The idea is to locate the local minimum points (the points where the residual is closest to zero)

I tried using plot_comples from the spb library. But apparently you can only plot symbolic functions.

I would like to make a plot similar to the plot below. Where the black dots are the values ​​where res is minimum.

Color plot


Solution

  • This is how spb creates that kind of visualization:

    1. create a grid of omega values in the complex plane.
    2. evaluate the function over that complex grid.
    3. compute the magnitude of the result, which varies from 0 to (possibly) infinity. The range in values is usually large, so it's not a good idea to plot the magnitude because you'd end up with an image difficult to understand...
    4. Instead, we compute the brightness with this formula: b = m / (m + 1) where m is the magnitude. Note that b is zero when m=0 and b tends to one when m tends to infinity.
    5. plot the brightness.

    Here is an example:

    import numpy as np
    import matplotlib.pyplot as plt
    
    func = lambda w: np.sin(w**2)
    n = 200j
    re_w, im_w = np.mgrid[-3:3:n, -3:3:n]
    w = re_w + 1j * im_w
    res = func(w)
    mag = np.absolute(res)
    brightness = mag / (mag + 1)
    
    fig, ax = plt.subplots()
    img = ax.imshow(brightness, cmap="gray", vmin=brightness.min(), vmax=brightness.max())
    cbar = fig.colorbar(img, label="brightness")
    

    enter image description here