{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# p33: Solve linear BVP, Neumann bc\n",
"\n",
"$$\n",
"u_{xx} = \\exp(4x)\n",
"$$\n",
"\n",
"$$\n",
"u'(-1)=0, \\qquad u(1)=0\n",
"$$"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"%config InlineBackend.figure_format='svg'\n",
"from chebPy import *\n",
"from numpy import dot,exp,zeros,max,linspace,polyval,polyfit,inf\n",
"from numpy.linalg import norm\n",
"from scipy.linalg import solve\n",
"from matplotlib.pyplot import title,plot,xlabel,ylabel,grid"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n",
"\n"
],
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"N = 16\n",
"\n",
"# Build matrix\n",
"D,x = cheb(N)\n",
"D2 = dot(D,D)\n",
"D2[0,:] = D[0,:] # First eqn has neumann bc\n",
"D2 = D2[0:N,0:N] # Remove last row and column\n",
"\n",
"# RHS\n",
"f = zeros(N)\n",
"f[1:] = exp(4.0*x[1:N])\n",
"\n",
"# Solve\n",
"u = solve(D2,f)\n",
"s = zeros(N+1)\n",
"s[0:N] = u\n",
"\n",
"# Compute error\n",
"xx = linspace(-1.0,1.0,200)\n",
"uu = polyval(polyfit(x,s,N),xx) # interpolate grid data\n",
"exact = (exp(4.0*xx) - 4.0*exp(-4.0)*(xx-1.0) - exp(4.0))/16.0\n",
"maxerr = norm(uu-exact,inf)\n",
"\n",
"title('max err = %e' % maxerr)\n",
"plot(x,s,'o',xx,exact)\n",
"xlabel('x'); ylabel('u'); grid(True);"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.12.8"
}
},
"nbformat": 4,
"nbformat_minor": 4
}