pi = np.pi
alpha = np.linspace(0,2*pi,201)
def fdlambda(alpha): # eigenvalues for finite-difference scheme
term = np.array( 1 - sigma**2*(1-np.cos(alpha)), dtype=complex )
lplus = term + np.sqrt( term**2 - 1 )
lminus = term - np.sqrt( term**2 - 1 )
return lplus,lminus
fig,ax = plt.subplots(2,3,figsize=(20,15))
ms = 4
for m,sigma in enumerate([.5,1,1.5]):#, 1, 1.1]:
lambdaplus,lambdaminus = fdlambda(alpha)
ax[1,m].plot(alpha,np.ones_like(alpha),ms=ms,label='PDE')
ax[1,m].plot(alpha,np.abs(lambdaplus),'.',ms=ms,label='FD')
ax[1,m].plot(alpha,np.abs(lambdaminus),'.',ms=ms)
ax[0,m].plot(alpha,np.angle(lambdaplus),'.',ms=ms,label='FD')
ax[0,m].plot(alpha,np.angle(lambdaminus),'.',ms=ms)
ax[0,m].plot(alpha, (alpha*sigma),label='PDE ') # arg eigenvalues for PDE
ax[0,m].plot(alpha,-(alpha*sigma))
ax[0,m].set_ylabel('arg of eigenvalues',fontsize=20)
ax[1,m].set_ylabel('modulus of eigenvalues',fontsize=20)
ax[0,m].set_xlabel('$\\alpha$',fontsize=20)
ax[1,m].set_xlabel('$\\alpha$',fontsize=20)
ax[0,m].legend();
ax[1,m].legend();
ax[0,m].set_xticks( [0,pi,2*pi], ['0','$\pi$','$2\pi$'],fontsize=20 )
ax[1,m].set_xticks( [0,pi,2*pi], ['0','$\pi$','$2\pi$'],fontsize=20 )
ax[0,m].set_title('$\\sigma =$ '+str(sigma),fontsize=20);