diff options
Diffstat (limited to 'Numerical_Methods_by_E._Balaguruswamy/chapter11.ipynb')
| -rw-r--r-- | Numerical_Methods_by_E._Balaguruswamy/chapter11.ipynb | 428 |
1 files changed, 428 insertions, 0 deletions
diff --git a/Numerical_Methods_by_E._Balaguruswamy/chapter11.ipynb b/Numerical_Methods_by_E._Balaguruswamy/chapter11.ipynb new file mode 100644 index 00000000..c717bdbd --- /dev/null +++ b/Numerical_Methods_by_E._Balaguruswamy/chapter11.ipynb @@ -0,0 +1,428 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Chapter 11 - Numerical differentiation" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Example No. 11_01 Pg No. 348" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "h\ty\terr\n", + "0.2 \t2.2 \t0.2\n", + "0.1 \t2.1 \t0.1\n", + "0.05 \t2.05 \t0.05\n", + "0.01 \t2.01 \t0.01\n" + ] + } + ], + "source": [ + "from __future__ import division\n", + "from scipy.misc import derivative\n", + "\n", + "#First order forward difference\n", + "\n", + "def f(x):\n", + " F=x**2\n", + " return F\n", + " \n", + "\n", + "def df(x,h):\n", + " DF=(f(x+h)-f(x))/h\n", + " return DF\n", + "\n", + "dfactual = derivative(f,1)\n", + "h = [0.2, 0.1, 0.05, 0.01 ]\n", + "print 'h\\ty\\terr'\n", + "for hh in h:\n", + " y = df(1,hh)\n", + " err = y - dfactual\n", + " print hh,'\\t',y,'\\t',err" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Example No. 11_02 Pg No. 350" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "h\ty\terr\n", + "0.2 \t2.0 \t-4.4408920985e-16\n", + "0.1 \t2.0 \t4.4408920985e-16\n", + "0.05 \t2.0 \t4.4408920985e-16\n", + "0.01 \t2.0 \t1.7763568394e-15\n" + ] + } + ], + "source": [ + "from __future__ import division\n", + "from scipy.misc import derivative\n", + "\n", + "#Three-Point Formula\n", + "\n", + "def f(x):\n", + " F = x**2\n", + " return F\n", + "def df(x,h):\n", + " DF = (f(x+h)-f(x-h))/(2*h)\n", + " return DF\n", + "dfactual = derivative(f,1)\n", + "h = [0.2, 0.1, 0.05, 0.01 ]\n", + "print 'h\\ty\\terr'\n", + "for hh in h:\n", + " y = df(1,hh)\n", + " err = y - dfactual\n", + " print hh,'\\t',y,'\\t',err\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Example No. 11_03 Pg No. 353" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "h\ty\t\terr\n", + "0.01 \t0.898257285429 \t-0.00218981692372\n", + "0.015 \t0.897151155627 \t-0.00329594672573\n", + "0.02 \t0.896037563392 \t-0.00440953896077\n", + "0.025 \t0.894916522709 \t-0.00553057964367\n", + "0.03 \t0.893788047675 \t-0.00665905467759\n", + "0.035 \t0.892652152498 \t-0.00779494985432\n", + "0.04 \t0.891508851498 \t-0.00893825085448\n", + "\n", + "M2 = 2.061e-08 & hopt = 4.706e-01\n" + ] + } + ], + "source": [ + "from __future__ import division\n", + "from scipy.misc import derivative\n", + "from numpy import arange,sin,cos,sqrt\n", + "x = 0.45#\n", + "def f(x):\n", + " F = sin(x)\n", + " return F\n", + "def df(x,h):\n", + " DF = (f(x+h) - f(x))/h\n", + " return DF\n", + "dfactual = cos(x)#\n", + "h = arange(0.01,0.005+0.04,0.005)\n", + "n = len(h)#\n", + "print 'h\\ty\\t\\terr'\n", + "for hh in h:\n", + " y = df(x,hh)\n", + " err = y - dfactual \n", + " print hh,'\\t',y,'\\t',err\n", + "\n", + "#using 16 significant digits so the bound for roundoff error is 0.5*10**(-16)\n", + "e = 0.5*10**(-16)\n", + "M2 = max(sin(x+h))#\n", + "hopt = 2*sqrt(e/M2)#\n", + "print '\\nM2 = %.3e & hopt = %.3e'%(hopt,M2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Example No. 11_04 Pg No. 354" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + " d2fexact = -0.732 \n", + " err = -6.097e-06 \n", + " y = -0.732 \n" + ] + } + ], + "source": [ + "from __future__ import division\n", + "from numpy import arange,sin,cos,sqrt\n", + "\n", + "#Approximate second derivative\n", + "\n", + "x = 0.75#\n", + "h = 0.01#\n", + "def f(x):\n", + " F = cos(x)\n", + " return F\n", + "def d2f(x,h):\n", + " D2F = ( f(x+h) - 2*f(x) + f(x-h) )/h**2\n", + " return D2F\n", + "y = d2f(0.75,0.01)#\n", + "d2fexact = -cos(0.75)\n", + "err = d2fexact - y #\n", + "print ' d2fexact = %.3f \\n err = %.3e \\n y = %.3f '%(d2fexact,err,y)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Example No. 11_05 Pg No. 358" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "v(5) = 4.25\n", + "v(7) = 5.5\n", + "v(9) = 6.75\n" + ] + } + ], + "source": [ + "#Differentiation of tabulated data\n", + "\n", + "T = range(5,10)\n", + "s = [10, 14.5, 19.5, 25.5 , 32 ]\n", + "h = T[1]-T[0]\n", + "n = len(T)\n", + "def v(t):\n", + " if t in T and T.index(t) == 0:\n", + " V = (-3*s[T.index(t) ] + 4*s[T.index(t+h)] - s[T.index(t+2*h) ] ) /(2*h) #Three point forward difference formula\n", + " elif t in T and T.index(t) == n-1:\n", + " V = ( 3*s[T.index(t)] - 4*s[ T.index(t-h)] + s[T.index(t-2*h) ] )/(2*h) #Backward difference formula\n", + " else:\n", + " V = ( s[T.index(t+h)] - s[T.index(t-h)] )/(2*h) #Central difference formula\n", + " return V\n", + "\n", + "v_5 = v(5)\n", + "v_7 = v(7)\n", + "v_9 = v(9)\n", + "\n", + "print 'v(5) = ',v_5\n", + "print 'v(7) = ',v_7\n", + "print 'v(9) = ',v_9\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Example No. 11_06 Pg No. 359" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "a(7) = 0.5\n" + ] + } + ], + "source": [ + "T = range(5,10)\n", + "s = [10, 14.5 , 19.5 , 25.5, 32 ]\n", + "h = T[1]-T[0]\n", + "def a(t):\n", + " A = ( s[T.index(t+h) ] - 2*s[ T.index(t) ] + s[T.index(t-h) ] )/h**2\n", + " return A\n", + "a_7 = a(6)\n", + "\n", + "print 'a(7) = ',a_7" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Example No. 11_7 Pg No. 359" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "solving the above equations we get \n", + "\n", + " y1 = y(0.25) = -0.117562 \n", + " y2 = y(0.5) = -0.168475 \n", + " y3 = y(0.75) = -0.139088 \n", + " \n" + ] + } + ], + "source": [ + "from numpy import mat,divide,linalg\n", + "h = 0.25 #\n", + "#y''(x) = e**(x**2)\n", + "#y(0) = 0 , y(1) = 0\n", + "# y''(x) = y(x+h) - 2*y(x) + y(x-h)/h**2 = e**(x**2)\n", + "#( y(x + 0.25) - 2*y(x) + y(x-0.25) )/0.0625 = e**(x**2)\n", + "#y(x+0.25) - 2*y(x) + y(x - 0.25) = 0.0624*e**(x**2)\n", + "#y(0.5) - 2*y(0.25) + y(0) = 0.0665\n", + "#y(0.75) - 2*y(0.5) + y(0.25) = 0.0803\n", + "#y(1) - 2*y(0.75) + y(0.5) = 0.1097\n", + "#given y(0) = y(1) = 0\n", + "#\n", + "#0 + y2 - 2y1 = 0.06665\n", + "#y3 - 2*y2 + y1 = 0.0803\n", + "#-2*y3 + y2 + 0 = 0.1097\n", + "#Therefore\n", + "A = mat([[0, 1, -2],[1, -2, 1],[-2, 1, 0 ]])\n", + "B = mat([[ 0.06665],[0.0803] ,[0.1097 ]])\n", + "C = linalg.solve(A,B)\n", + "print 'solving the above equations we get \\n\\n y1 = y(0.25) = %f \\n y2 = y(0.5) = %f \\n y3 = y(0.75) = %f \\n '%(C[2],C[1],C[0])" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Example No. 11_08 Pg No. 362" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "richardsons technique -\n", + "\n", + "df(0.5) = 1.64850499031\n", + "D(rh) = D(0.25) = 1.71828182846\n", + "D(0.5) = 1.66594919985\n", + "Exact df(0.5) = 1.93757920531\n", + "The result by richardsons technique is much better than other results\n", + "for r = 2\n", + "df(x) = 1.64518270284\n", + "D(rh) = 1.93757920531\n" + ] + } + ], + "source": [ + "from numpy import arange,exp\n", + "from scipy.misc import derivative\n", + "\n", + "x = arange(-0.5,0.25+1.5,0.25)\n", + "h = 0.5 ;\n", + "r = 1.0/2 ;\n", + "\n", + "def f(x):\n", + " F = exp(x)\n", + " return F\n", + "def D(x,h):\n", + " D = (f(x + h) - f(x-h) )/(2*h) \n", + " return D\n", + "def df(x,h,r):\n", + " Df = (D(x,r*h) - r**2*D(x,h))/(1-r**2)\n", + " return Df\n", + "\n", + "df_05 = df(0.5,0.5,1.0/2);\n", + "print 'richardsons technique -\\n\\ndf(0.5) = ',df_05\n", + "print 'D(rh) = D(0.25) = ',D(0.5,0.5)\n", + "print 'D(0.5) = ',D(0.5,0.25)\n", + "\n", + "dfexact = derivative(f,0.5)\n", + "print 'Exact df(0.5) = ',dfexact\n", + "print 'The result by richardsons technique is much better than other results'\n", + "\n", + "#r = 2\n", + "print 'for r = 2'\n", + "print 'df(x) = ',df(0.5,0.5,2)\n", + "print 'D(rh) = ',D(0.5,2*0.5)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 2", + "language": "python", + "name": "python2" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 2 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython2", + "version": "2.7.9" + } + }, + "nbformat": 4, + "nbformat_minor": 0 +} |