Numerically drawing the region of absolute stability

In [4]:
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

def eigenvalue(hlre,hlim):
    hlambda = hlre+hlim*1j # form the complex hlambda from real and imaginary parts
    return abs( 1 + hlambda + hlambda**2/2 )  # midpoint method

r = 3  # plot box radius
x = np.arange(-r,r, 0.1)
y = np.arange(-r,r, 0.1)
X, Y = np.meshgrid(x, y)  # grid of points to evaluate eigenvalue

Z = eigenvalue(X,Y)

fig, ax = plt.subplots(figsize=(6,6))
ax.contourf(X,Y,Z,levels=[0,1])
ax.contour(X,Y,Z,levels=[1])
plt.grid()
plt.title('region of absolute stability in $\mathbb{C}$ for midpoint method')
plt.show()
In [ ]: