{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Chapter 27 : Numerical Solution Of Ordinary Differential Equations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.1, page no. 686"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution through picards method \n",
      "The no of iterations required3\n",
      "y(0)=1 and y(x)=x+y \n",
      "Y =  x**4/24 + x**3/3 + x**2 + x + 1\n"
     ]
    }
   ],
   "source": [
    "import numpy,sympy,math\n",
    "\n",
    "x = sympy.Symbol('x')\n",
    "print \"Solution through picards method \"\n",
    "n = int(raw_input(\"The no of iterations required\")) \n",
    "print \"y(0)=1 and y(x)=x+y \"\n",
    "yo = 1\n",
    "yn = 1\n",
    "for i in range(1,n+1):\n",
    "    yn = yo+sympy.integrate(yn+x,(x,0,x))\n",
    "print \"Y = \",yn"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.2, page no. 687"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution through picards method \n",
      "The no of iterations required : 2\n",
      "y(0)=1 and y(x)=x+y \n",
      "Y =  -Integral(2*x/(2*log(x + 1) + 1), x) + Integral(2*x/(2*log(x + 1) + 1), (x, 0)) - Integral(-2*log(x + 1)/(2*log(x + 1) + 1), x) + Integral(-2*log(x + 1)/(2*log(x + 1) + 1), (x, 0)) - Integral(-1/(2*log(x + 1) + 1), x) + Integral(-1/(2*log(x + 1) + 1), (x, 0)) + 1.0\n"
     ]
    }
   ],
   "source": [
    "import numpy,sympy,math\n",
    "\n",
    "x = sympy.Symbol('x')\n",
    "print \"Solution through picards method \"\n",
    "n = int(raw_input(\"The no of iterations required : \")) \n",
    "print \"y(0)=1 and y(x)=x+y \"\n",
    "yo = 1\n",
    "y = 1\n",
    "for i in range(1,n+1):\n",
    "    f = (y-x)/(y+x) \n",
    "    y = yo+sympy.integrate(f,(x,0,x))\n",
    "x = 0.1\n",
    "print \"Y = \",y.evalf()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.5, page no. 690"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution using Eulers Method \n",
      "0.3\n",
      "1.362\n",
      "Input the number of iteration :− 4\n",
      "The value of y is :−  1.5282\n"
     ]
    }
   ],
   "source": [
    "print \"Solution using Eulers Method \"\n",
    "print x\n",
    "print y\n",
    "n = int(raw_input(\"Input the number of iteration :− \"))\n",
    "x = 0\n",
    "y = 1\n",
    "for i in range(1,n+1):\n",
    "    y1 = x+y\n",
    "    y = y+0.1*y1\n",
    "    x = x+0.1\n",
    "print \"The value of y is :− \",y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.6, page no. 691"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution using Eulers Method \n",
      "0.4\n",
      "1.5282\n",
      "Input the number of iteration :− 2\n",
      "1.02\n",
      "1.03642857143\n",
      "The value of y is :−  1.03642857143\n"
     ]
    }
   ],
   "source": [
    "print \"Solution using Eulers Method \"\n",
    "print x\n",
    "print y\n",
    "n = int(raw_input(\"Input the number of iteration :− \"))\n",
    "x = 0\n",
    "y = 1\n",
    "for i in range(1,n+1):\n",
    "    y1 = (y-x)/(y+x)\n",
    "    y = y+0.02*y1\n",
    "    x = x+0.1\n",
    "    print y\n",
    "print \"The value of y is :− \",y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.7, page no. 692"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution using Eulers Method \n",
      "0.2\n",
      "1.03642857143\n",
      "Input the number of iteration :− 3\n",
      "1.1\n",
      "1.11\n",
      "1.1105\n",
      "1.110525\n",
      "−−−−−−−−−−−−−−−−−−−−−−− \n",
      "1.2315775\n",
      "1.2315775\n",
      "1.242630125\n",
      "1.24318275625\n",
      "1.24321038781\n",
      "−−−−−−−−−−−−−−−−−−−−−−− \n",
      "1.38753142659\n",
      "1.38753142659\n",
      "1.39974747853\n",
      "1.40035828113\n",
      "1.40038882126\n",
      "−−−−−−−−−−−−−−−−−−−−−−− \n",
      "1.57042770339\n",
      "The value of y is :−  1.40038882126\n"
     ]
    }
   ],
   "source": [
    "print \"Solution using Eulers Method \"\n",
    "print x\n",
    "print y\n",
    "n = int(raw_input(\"Input the number of iteration :− \"))\n",
    "x = 0.1\n",
    "m = 1\n",
    "y = 1\n",
    "yn = 1\n",
    "y1 = 1\n",
    "k = 1\n",
    "for i in range(1,n+1):\n",
    "    yn = y \n",
    "    for i in range(1,5):\n",
    "        m = (k+y1)/2\n",
    "        yn = y+0.1*m\n",
    "        y1 = (yn+x)\n",
    "        print yn\n",
    "    print \"−−−−−−−−−−−−−−−−−−−−−−− \"\n",
    "    y = yn \n",
    "    m = y1\n",
    "    yn = yn+0.1*m \n",
    "    print yn\n",
    "    x = x+0.1\n",
    "    yn = y \n",
    "    k = m\n",
    "print \"The value of y is :− \",y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.8, page no. 694"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution using Eulers Method \n",
      "0.2\n",
      "2\n",
      "Input the number of iteration :− 3\n",
      "2.06020299957\n",
      "2.06551474469\n",
      "2.06561668931\n",
      "2.06561864352\n",
      "−−−−−−−−−−−−−−−−−−−−−−−\n",
      "2.13665600548\n",
      "2.13665600548\n",
      "2.14156348217\n",
      "2.14164742068\n",
      "2.14164885497\n",
      "−−−−−−−−−−−−−−−−−−−−−−−\n",
      "2.22267196493\n",
      "2.22267196493\n",
      "2.22722645093\n",
      "2.22729646948\n",
      "2.22729754504\n",
      "−−−−−−−−−−−−−−−−−−−−−−−\n",
      "2.31757184825\n",
      "The value of y is :−  2.22729754504\n"
     ]
    }
   ],
   "source": [
    "import math\n",
    "\n",
    "print \"Solution using Eulers Method \"\n",
    "print x\n",
    "print y\n",
    "n = int(raw_input(\"Input the number of iteration :− \"))\n",
    "x = 0.2\n",
    "m = 0.301\n",
    "y = 2\n",
    "yn = 2\n",
    "y1 = math.log10(2)\n",
    "k = 0.301\n",
    "for i in range(1,n+1):\n",
    "    yn = y\n",
    "    for i in range(1,5):\n",
    "        m = (k+y1)/2\n",
    "        yn = y+0.2*m\n",
    "        y1 = math.log10(yn+x)\n",
    "        print yn\n",
    "    print \"−−−−−−−−−−−−−−−−−−−−−−−\"\n",
    "    y = yn\n",
    "    m = y1\n",
    "    yn = yn+0.2*m \n",
    "    print yn\n",
    "    x = x+0.2\n",
    "    yn = y\n",
    "    k = m\n",
    "print \"The value of y is :− \",y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.9, page no. 696"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution using Eulers Method \n",
      "0.8\n",
      "2.22729754504\n",
      "Input the number of iteration :− 2\n",
      "1.2\n",
      "1.2295445115\n",
      "1.23088482816\n",
      "1.23094524903\n",
      "−−−−−−−−−−−−−−−−−−−−−−− \n",
      "1.49284119302\n",
      "1.49284119302\n",
      "1.52407510156\n",
      "1.52534665765\n",
      "1.52539814634\n",
      "−−−−−−−−−−−−−−−−−−−−−−− \n",
      "1.85241216588\n",
      "The value of y is :−  1.52539814634\n"
     ]
    }
   ],
   "source": [
    "print \"Solution using Eulers Method \"\n",
    "print x\n",
    "print y\n",
    "n = int(raw_input(\"Input the number of iteration :− \"))\n",
    "x = 0.2 \n",
    "m = 1\n",
    "y = 1\n",
    "yn = 1\n",
    "y1 = 1\n",
    "k = 1\n",
    "for i in range(1,n+1):\n",
    "    yn = y\n",
    "    for i in range(1,5):\n",
    "        m = (k+y1)/2\n",
    "        yn = y+0.2*m \n",
    "        y1 = math.sqrt(yn)+x\n",
    "        print yn\n",
    "    print \"−−−−−−−−−−−−−−−−−−−−−−− \"\n",
    "    y = yn\n",
    "    m = y1\n",
    "    yn = yn+0.2*m\n",
    "    print yn\n",
    "    x = x+0.2\n",
    "    yn = y\n",
    "    k = m\n",
    "print \"The value of y is :− \",y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## Example 27.10, page no. 697"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Runges method\n",
      "The required approximate value is :− \n",
      "1.0\n"
     ]
    }
   ],
   "source": [
    "print \"Runges method\"\n",
    "def f(x,y):\n",
    "    y = x+y\n",
    "    return y\n",
    "x = 0\n",
    "y = 1\n",
    "h = 0.2\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "kk = h*f(x+h,y+k1)\n",
    "k3 = h*f(x+h,y+kk)\n",
    "k = 1/6*(k1+4*k2+k3)\n",
    "print \"The required approximate value is :− \"\n",
    "y = y+k\n",
    "print y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.11, page no. 697"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Runges method\n",
      "The required approximate value is :− \n",
      "1.0\n"
     ]
    }
   ],
   "source": [
    "print \"Runges method\"\n",
    "def f(x,y):\n",
    "    y = x+y\n",
    "    return y\n",
    "x = 0\n",
    "y = 1\n",
    "h = 0.2\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "k = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y = y+k\n",
    "print y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.12, page no. 698"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Runga kutta method\n",
      "The required approximate value is :− \n",
      "1.0\n",
      "To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \n",
      "The required approximate value is :− \n",
      "1.0\n"
     ]
    }
   ],
   "source": [
    "print \"Runga kutta method\"\n",
    "def f(x,y):\n",
    "    y = (y**2-x**2)/(x**2+y**2)\n",
    "    return y\n",
    "x = 0\n",
    "y = 1\n",
    "h = 0.2\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "k = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y = y+k\n",
    "print y\n",
    "print \"To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \"\n",
    "x = 0.2\n",
    "h = 0.2\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "k = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y = y+k\n",
    "print y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.13, page no. 699"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Runga kutta method\n",
      "The required approximate value is :− \n",
      "1.0\n",
      "To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \n",
      "The required approximate value is :− \n",
      "1.0\n"
     ]
    }
   ],
   "source": [
    "print \"Runga kutta method\"\n",
    "def f(x,y):\n",
    "    yy = x+y**2\n",
    "    return yy\n",
    "x = 0\n",
    "y = 1\n",
    "h = 0.1\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2) \n",
    "k4 = h*f(x+h,y+k3)\n",
    "k = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y = y+k\n",
    "print y\n",
    "print \"To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \"\n",
    "x = 0.1\n",
    "h = 0.1\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "k = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y = y+k\n",
    "print y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.14, page no. 700"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "y1=  x**2/2\n",
      "y2=  -x**5/20 + x**2/2\n",
      "y3=  -x**5/20 + x**2/2\n",
      "Determining the initial values for milnes method using y3 \n",
      "x=0.0 y0=0.0 f0=0 \n",
      "x=0.2 y1=  x**2/2\n",
      "f1= -x**4/4 + 0.2\n",
      "x=0.4 y2=  -x**5/20 + x**2/2\n",
      "f2=  -(-x**5/20 + x**2/2)**2 + 0.4\n",
      "x=0.6 y3=  -x**5/20 + x**2/2\n",
      "f3= -(-x**5/20 + x**2/2)**2 + 0.6\n",
      "Using predictor method to find y4\n",
      "y4= -0.1*x**4 - 0.2*(-x**5/20 + x**2/2)**2 + 0.24\n",
      "f4=  -(-x**5/20 + x**2/2)**2 + 0.8\n",
      "Using predictor method to find y5 \n",
      "[-1.53353211362738, 2.59188514369719, -1.51825864909117 - 1.96576535664013*I, -1.51825864909117 + 1.96576535664013*I, -0.611593620665791 - 2.14551655586226*I, -0.611593620665791 + 2.14551655586226*I, 0.00495892544779777 - 0.780181689871346*I, 0.00495892544779777 + 0.780181689871346*I, 1.59571682927425 - 0.8271211188914*I, 1.59571682927425 + 0.8271211188914*I]\n",
      "f5 = [-1.28459666867223, 2.38248090853886, -1.33416685318757 - 1.68618279990398*I, -1.33416685318757 + 1.68618279990398*I, -1.00716479884721 - 2.10036734766602*I, -1.00716479884721 + 2.10036734766602*I, 0.142926398918135 - 1.33989714259572*I, 0.142926398918135 + 1.33989714259572*I, 1.64946313318333 - 0.385565643189579*I, 1.64946313318333 + 0.385565643189579*I]\n",
      "Hence y(1)=  [-1.53353211362738, 2.59188514369719, -1.51825864909117 - 1.96576535664013*I, -1.51825864909117 + 1.96576535664013*I, -0.611593620665791 - 2.14551655586226*I, -0.611593620665791 + 2.14551655586226*I, 0.00495892544779777 - 0.780181689871346*I, 0.00495892544779777 + 0.780181689871346*I, 1.59571682927425 - 0.8271211188914*I, 1.59571682927425 + 0.8271211188914*I]\n"
     ]
    }
   ],
   "source": [
    "import sympy\n",
    "\n",
    "x = sympy.Symbol('x')\n",
    "yo = 0\n",
    "y = 0\n",
    "h = 0.2\n",
    "f = x-y**2\n",
    "y = sympy.integrate(f,(x,0,x))\n",
    "y1 = yo+y\n",
    "print \"y1= \",y1\n",
    "f = x-y**2\n",
    "y = sympy.integrate(f,(x,0,x))\n",
    "y2 = yo+y\n",
    "print \"y2= \",y2\n",
    "y = x-y**2\n",
    "y = sympy.integrate(f,(x,0,x))\n",
    "y3 = yo+y\n",
    "print \"y3= \",y3\n",
    "print \"Determining the initial values for milnes method using y3 \"\n",
    "print \"x=0.0 y0=0.0 f0=0 \"\n",
    "print \"x=0.2 y1= \",\n",
    "x = 0.2\n",
    "print y1\n",
    "#y1 = sympy.solve(y1)\n",
    "print \"f1=\",\n",
    "f1 = x-y1**2\n",
    "print f1\n",
    "print \"x=0.4 y2= \",\n",
    "x = 0.4\n",
    "print y2\n",
    "print \"f2= \",\n",
    "f2 = x-y2**2\n",
    "print f2\n",
    "print \"x=0.6 y3= \",\n",
    "x = 0.6\n",
    "print y3\n",
    "print \"f3=\",\n",
    "f3 = x-y3**2\n",
    "print f3\n",
    "print \"Using predictor method to find y4\"\n",
    "x = 0.8\n",
    "y4 = yo+4/3*h*(2*f1-f2+2*f3)\n",
    "print \"y4=\",y4\n",
    "f4 = x-y**2\n",
    "print \"f4= \",f4\n",
    "print \"Using predictor method to find y5 \"\n",
    "x = 1.0\n",
    "y5 = sympy.solve(y1+4/3*h*(2*f2-f3+2*f4))\n",
    "print y5\n",
    "f5 = sympy.solve(x-y**2)\n",
    "print \"f5 =\",f5\n",
    "print \"Hence y(1)= \",y5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.15, page no. 704"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Runga kutta method\n",
      "The required approximate value is :− \n",
      "1.0\n",
      "To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \n",
      "The required approximate value is :− \n",
      "1.0\n",
      "To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \n",
      "The required approximate value is :− \n",
      "1.0\n",
      "y0 y1 y2 y3 are respectively: \n",
      "1.0 1.0 1.0 1\n",
      "f0 f1 f2 f3 are respectively: \n",
      "1.3 1.2 1.1 1\n",
      "finding y4 using predictors milne method x=0.4\n",
      "y4 =  1.48\n",
      "f4 =  2.7824\n",
      "Using corrector method : \n",
      "y4 =  -0.0333333333333333*(-x**5/20 + x**2/2)**2 + 1.24\n",
      "f4 =  -0.0133333333333333*(-x**5/20 + x**2/2)**2 + (-0.0333333333333333*(-x**5/20 + x**2/2)**2 + 1.24)**2 + 0.496\n"
     ]
    }
   ],
   "source": [
    "print \"Runga kutta method\"\n",
    "def f(x,y):\n",
    "    yy = x*y+y**2\n",
    "    return yy\n",
    "y0 = 1\n",
    "x = 0\n",
    "y = 1\n",
    "h = 0.1\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "ka = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y1 = y+ka \n",
    "y = y+ka \n",
    "print y\n",
    "print \"To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \"\n",
    "x = 0.1\n",
    "h = 0.1\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "kb = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y2 = y+kb \n",
    "y = y+kb\n",
    "print y\n",
    "print \"To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \"\n",
    "x = 0.2\n",
    "h = 0.1\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "kc = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y3 = y+kc\n",
    "y = y+kc \n",
    "print y\n",
    "f0 = f(0,y0)\n",
    "f1 = f(0.1,y1)\n",
    "f2 = f(0.2,y2)\n",
    "f3 = f(0.3,y3)\n",
    "print \"y0 y1 y2 y3 are respectively: \"\n",
    "print y3,y2,y1,y0\n",
    "print \"f0 f1 f2 f3 are respectively: \"\n",
    "print f3,f2,f1,f0\n",
    "print \"finding y4 using predictors milne method x=0.4\"\n",
    "h = 0.1\n",
    "y4 = y0+4*h/3*(2*f1-f2+2*f3)\n",
    "print \"y4 = \",y4\n",
    "print \"f4 = \",f(0.4,y4)\n",
    "print \"Using corrector method : \"\n",
    "y4 = y2+h/3*(f2+4*f3+f4) \n",
    "print \"y4 = \",y4\n",
    "print \"f4 = \",f(0.4,y4) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.16, page no. 705"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using predictor method \n",
      "y11 =  2.57229741667\n",
      "Using corrector method \n",
      "y11 =  2.57494727679\n",
      "f11 =  7.00689666251\n"
     ]
    }
   ],
   "source": [
    "def f(x,y):\n",
    "    yy = x**2*(1+y)\n",
    "    return yy\n",
    "y3 = 1\n",
    "y2 = 1.233\n",
    "y1 = 1.548\n",
    "y0 = 1.979\n",
    "f3 = f(1,y3)\n",
    "f2 = f(1.1,y2)\n",
    "f1 = f(1.2,y1)\n",
    "f0 = f(1.3,y0)\n",
    "print \"Using predictor method \"\n",
    "h = 0.1\n",
    "y11 = y0+h/24*(55*f0-59*f1+37*f2-9*f3)\n",
    "print \"y11 = \",y11\n",
    "x = 1.4\n",
    "f11 = f(1.4,y11)\n",
    "print \"Using corrector method \"\n",
    "y11 = y0+h/24*(9*f11+19*f0-5*f1+f2)\n",
    "print \"y11 = \",y11\n",
    "f11 = f(1.4,y11)\n",
    "print \"f11 = \",f11"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.17, page no. 706"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Runga kutta method\n",
      "The required approximate value is :− \n",
      "1.0\n",
      "To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \n",
      "The required approximate value is :− \n",
      "1.0\n",
      "To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \n",
      "The required approximate value is :− \n",
      "1.0\n",
      "y0 y1 y2 y3 are respectively :\n",
      "1.0 1.0 1.0 1\n",
      "f0 f1 f2 f3 are respectively :\n",
      "-0.7 -0.8 -0.9 -1\n",
      "Using adams method \n",
      "Using the predictor\n",
      "y4 =  0.935\n",
      "Using corrector method \n",
      "y4 =  0.9397165625\n",
      "f4 =  -0.483067217837\n"
     ]
    }
   ],
   "source": [
    "print \"Runga kutta method\"\n",
    "def f(x,y):\n",
    "    yy = x-y**2\n",
    "    return yy\n",
    "y0 = 1\n",
    "x = 0\n",
    "y = 1\n",
    "h = 0.1\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "ka = 1/6*(k1+2*k2+2*k3+k4) \n",
    "print \"The required approximate value is :− \"\n",
    "y1 = y+ka\n",
    "y = y+ka\n",
    "print y\n",
    "print \"To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \"\n",
    "x = 0.1\n",
    "h = 0.1\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "kb = 1/6*(k1+2*k2+2*k3+k4)\n",
    "print \"The required approximate value is :− \"\n",
    "y2 = y+kb\n",
    "y = y+kb\n",
    "print y\n",
    "print \"To find y(0.4) put x=0.2 y=above value ie 1.196 h=0.2 \"\n",
    "x = 0.2\n",
    "h = 0.1\n",
    "k1 = h*f(x,y)\n",
    "k2 = h*f(x+1/2*h,y+1/2*k1)\n",
    "k3 = h*f(x+1/2*h,y+1/2*k2)\n",
    "k4 = h*f(x+h,y+k3)\n",
    "kc = 1/6*(k1+2*k2+2*k3+k4) \n",
    "print \"The required approximate value is :− \"\n",
    "y3 = y+kc\n",
    "y = y+kc\n",
    "print y\n",
    "f0 = f(0,y0)\n",
    "f1 = f(0.1,y1)\n",
    "f2 = f(0.2,y2)\n",
    "f3 = f(0.3,y3) \n",
    "print \"y0 y1 y2 y3 are respectively :\"\n",
    "print y3,y2,y1,y0\n",
    "print \"f0 f1 f2 f3 are respectively :\"\n",
    "print f3,f2,f1,f0\n",
    "print \"Using adams method \"\n",
    "print \"Using the predictor\"\n",
    "h = 0.1\n",
    "y4 = y3+h/24*(55*f3-59*f2+37*f1-9*f0)\n",
    "x = 0.4\n",
    "f4 = f(0.4,y4)\n",
    "print \"y4 = \",y4\n",
    "print \"Using corrector method \"\n",
    "y4 = y3+h/24*(9*f4+19*f3-5*f2+f1) \n",
    "print \"y4 = \",y4\n",
    "f4 = f(0.4,y4)\n",
    "print \"f4 = \",f4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.18, page no. 708"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " Picards method\n",
      "First approximation\n",
      "y1 =  x**2/2 + x + 2\n",
      "z1 =  x**2/2 - 4*x + 1\n",
      "Second approximation\n",
      "y2 =  x**3/6 - 3*x**2/2 + x + 2\n",
      "z2 =  -x**5/20 - x**4/4 - x**3 - 3*x**2/2 - 4*x + 1\n",
      "Third approximation\n",
      "y3 =  -x**6/120 - x**5/20 - x**4/4 - x**3/2 - 3*x**2/2 + x + 2\n",
      "z3 =  -x**7/252 + x**6/12 - 31*x**5/60 + 7*x**4/12 + 5*x**3/3 - 3*x**2/2 - 4*x + 1\n",
      "y(0.1) =  -0.00833333333333333*x**6 - 0.05*x**5 - 0.25*x**4 - 0.5*x**3 - 1.5*x**2 + x + 2.0\n",
      "z(0.1) =  -0.00396825396825397*x**7 + 0.0833333333333333*x**6 - 0.516666666666667*x**5 + 0.583333333333333*x**4 + 1.66666666666667*x**3 - 1.5*x**2 - 4.0*x + 1.0\n"
     ]
    }
   ],
   "source": [
    "import sympy\n",
    "\n",
    "print \"Picards method\"\n",
    "x0 = 0\n",
    "y0 = 2\n",
    "z0 = 1\n",
    "x = sympy.Symbol('x')\n",
    "def f(x,y,z):\n",
    "    yy = x+z\n",
    "    return yy\n",
    "def g(x,y,z):\n",
    "    yy = x-y**2\n",
    "    return yy\n",
    "print \"First approximation\"\n",
    "y1 = y0+sympy.integrate(f(x,y0,z0),(x,x0,x))\n",
    "print \"y1 = \",y1\n",
    "z1 = z0+sympy.integrate(g(x,y0,z0),(x,x0,x))\n",
    "print \"z1 = \",z1\n",
    "print \"Second approximation\"\n",
    "y2 = y0+sympy.integrate(f(x,y1,z1),(x,x0,x))\n",
    "print \"y2 = \",y2\n",
    "z2 = z0+sympy.integrate(g(x,y1,z1),(x,x0,x))\n",
    "print \"z2 = \",z2\n",
    "print \"Third approximation\"\n",
    "y3 = y0+sympy.integrate(f(x,y2,z2),(x,x0,x))\n",
    "print \"y3 = \",y3\n",
    "z3 = z0+sympy.integrate(g(x,y2,z2),(x,x0,x))\n",
    "print \"z3 = \",z3\n",
    "x = 0.1\n",
    "print \"y(0.1) = \",y3.evalf()\n",
    "print \"z(0.1) = \",z3.evalf()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 27.19, page no. 710"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using k1 k2.. for f and l1 l2... for g runga kutta formula becomes \n",
      "y =  1.0\n",
      "y1 =  0.0\n"
     ]
    }
   ],
   "source": [
    "import sympy\n",
    "\n",
    "x = sympy.Symbol('x')\n",
    "def f(x,y,z):\n",
    "    yy = z \n",
    "    return yy\n",
    "def g(x,y,z):\n",
    "    yy = x*y**2-y**2\n",
    "    return yy\n",
    "x0 = 0\n",
    "y0 = 1\n",
    "z0 = 0\n",
    "h = 0.2\n",
    "print \"Using k1 k2.. for f and l1 l2... for g runga kutta formula becomes \"\n",
    "h = 0.2\n",
    "k1 = h*f(x0,y0,z0)\n",
    "l1 = h*g(x0,y0,z0)\n",
    "k2 = h*f(x0+1/2*h,y0+1/2*k1,z0+1/2*l1)\n",
    "l2 = h*g(x0+1/2*h,y0+1/2*k1,z0+1/2*l1)\n",
    "k3 = h*f(x0+1/2*h,y0+1/2*k2,z0+1/2*l2)\n",
    "l3 = h*g(x0+1/2*h,y0+1/2*k2,z0+1/2*l2)\n",
    "k4 = h*f(x0+h,y0+k3,z0+l3)\n",
    "l4 = h*g(x0+h,y0+k3,z0+l3)\n",
    "k = 1/6*(k1+2*k2+2*k3+k4)\n",
    "l = 1/6*(l1+2*l2+2*l3+2*l4)\n",
    "x = 0.2\n",
    "y = y0+k\n",
    "y1 = z0+l\n",
    "print \"y = \",y\n",
    "print \"y1 = \",y1"
   ]
  }
 ],
 "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.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}