p13: Solve linear BVP - Spectral Methods

u_{xx} = \exp(4x), \qquad u(-1)=u(1)=0

(1)

#%config InlineBackend.figure_format='svg'
from chebPy import *
from numpy import dot,exp,zeros,sinh,cosh,max,linspace,inf
from numpy.linalg import norm
from scipy.linalg import solve
from scipy.interpolate import barycentric_interpolate
from matplotlib.pyplot import title,plot,legend,xlabel

N = 16
D,x = cheb(N)
D2 = dot(D,D)
D2 = D2[1:N,1:N]
f = exp(4.0*x[1:N])
u = solve(D2,f)
s = zeros(N+1)
s[1:N] = u

# Do Chebyshev interpolation to finer grid for plotting and error estimation
xx = linspace(-1.0,1.0,200)
uu = barycentric_interpolate(x,s,xx)
exact = (exp(4.0*xx) - sinh(4.0)*xx - cosh(4.0))/16.0
maxerr = norm(uu-exact,inf)

title('max err = %e' % maxerr)
plot(x,s,'o',label='Chebyshev')
plot(xx,exact,label='Exact')
legend(), xlabel('x');