{
"metadata": {
"name": "",
"signature": "sha256:603c1d9f62c1fe2e9596bb9d4bada2c456d10b7d8c6fd0ac3894ca7b266609a5"
},
"nbformat": 3,
"nbformat_minor": 0,
"worksheets": [
{
"cells": [
{
"cell_type": "heading",
"level": 1,
"metadata": {},
"source": [
"BVP using finite difference method"
]
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"import numpy as np\n",
"from matplotlib import pyplot as plt"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 1
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"def yexact(x):\n",
" return 1.0/(np.exp(x)+np.exp(-x))\n",
"\n",
"def f(y,dy):\n",
" return -y + 2.0*dy**2/y\n",
"\n",
"def phi(y):\n",
" n = len(y) - 1\n",
" h = 2.0/n\n",
" res = np.zeros(n-1)\n",
" for i in range(1,n):\n",
" dy = (y[i+1] - y[i-1])/(2.0*h)\n",
" res[i-1] = (y[i-1]-2.0*y[i]+y[i+1])/h**2 - f(y[i],dy)\n",
" return res\n",
" \n",
"def dphi(y):\n",
" n = len(y) - 1\n",
" h = 2.0/n\n",
" res = np.zeros((n-1,n-1))\n",
" for i in range(1,n):\n",
" dy = (y[i+1] - y[i-1])/(2.0*h)\n",
" if i > 1:\n",
" res[i-1,i-2] = 1.0/h**2 + 4.0*dy/y[i]\n",
" res[i-1,i-1] =-2.0/h**2 + 1.0 - 2.0*(dy/y[i])**2\n",
" if i < n-1:\n",
" res[i-1,i ] = 1.0/h**2 - 4.0*dy/y[i]\n",
" return res"
],
"language": "python",
"metadata": {},
"outputs": [],
"prompt_number": 2
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"n = 100\n",
"\n",
"# Initial guess for y\n",
"y = np.zeros(n+1)\n",
"y[:] = 1.0/(exp(1) + exp(-1))\n",
"\n",
"maxiter = 100\n",
"TOL = 1.0e-6\n",
"it = 0\n",
"while it < maxiter:\n",
" b = phi(y)\n",
" if np.linalg.norm(b) < TOL:\n",
" break\n",
" A = dphi(y)\n",
" v = np.linalg.solve(A,b)\n",
" y[1:n] = y[1:n] - v\n",
" it += 1\n",
"print \"Number of iterations = %d\" % it\n",
"x = np.linspace(-1.0,1.0,n+1)\n",
"plt.plot(x,y,'o',x,yexact(x))\n",
"plt.legend((\"Numerical\",\"Exact\"));"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"Number of iterations = 13\n"
]
},
{
"metadata": {},
"output_type": "display_data",
"svg": [
"\n",
"\n",
"\n",
"\n"
],
"text": [
""
]
}
],
"prompt_number": 3
}
],
"metadata": {}
}
]
}