p15: Solve eigenvalue BVP

p15: Solve eigenvalue BVP#

\[ u_{xx} = \lambda u, \qquad u(-1)=u(1)=0 \]
%config InlineBackend.figure_format='svg'
from pylab import *
from chebPy import cheb
from scipy.linalg import solve,eig
from scipy.interpolate import barycentric_interpolate
N = 36
D,x = cheb(N)
D2 = dot(D,D)
D2 = D2[1:N,1:N]

lam,V = eig(D2)   
ii = argsort(-lam)
lam = real(lam[ii])
V = V[:,ii]

fig = figure(figsize=(10,15))
for j in range(5,35,5):              
    lv = shape(V)[0]+2
    u = zeros(lv)
    u[1:lv-1] = V[:,int(j)]  
    subplot(6,1,j//5)
    plot(x,u,'bo')
    xx = linspace(-1.0,1.0,501)
    uu = barycentric_interpolate(x,u,xx)
    s = 'eig %d = %20.13f * $\\pi^2/4$' %(j,lam[j-1]*4/pi**2)
    s = s + '\t\t %4.1f ppw' % (4*N/(pi*j))
    title(s)
    plot(xx,uu,'b')
    axis('off')
_images/80d28fc75bef39af6668ebadef8d36e3150c5e9a98f65c31bae3e969a3febe4a.svg
lexact = - pi**2/4 * arange(1,N)**2
error = abs((lam - lexact)/lexact)
semilogy(error,'o')
grid(True), xlabel('n'), ylabel('Error in $\\lambda_n$');
_images/e2cfbf3f1cb8bc25f11377f04027c78fbf4c5f744ac8c8603a97439e0cae6761.svg