{ "metadata": { "name": "", "signature": "sha256:5ad505722b4750f74beb79bb35accdae672d552c4292ed1c194aafb6e3675d9a" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "Conside the ODE\n", "$$\n", "y' = -y + 2 \\exp(-t) \\cos(2t)\n", "$$\n", "with initial condition\n", "$$\n", "y(0) = 0\n", "$$\n", "The exact solution is\n", "$$\n", "y(t) = \\exp(-t) \\sin(2t)\n", "$$" ] }, { "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": "markdown", "metadata": {}, "source": [ "Right hand side function" ] }, { "cell_type": "code", "collapsed": false, "input": [ "def f(t,y):\n", " return -y + 2.0*np.exp(-t)*np.cos(2.0*t)" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 2 }, { "cell_type": "markdown", "metadata": {}, "source": [ "Exact solution" ] }, { "cell_type": "code", "collapsed": false, "input": [ "def yexact(t):\n", " return np.exp(-t)*np.sin(2.0*t)" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 3 }, { "cell_type": "markdown", "metadata": {}, "source": [ "This implements Euler method\n", "$$\n", "y_n = y_{n-1} + h f(t_{n-1},y_{n-1})\n", "$$" ] }, { "cell_type": "code", "collapsed": false, "input": [ "def euler(t0,T,y0,h):\n", " N = int((T-t0)/h)\n", " y = np.zeros(N)\n", " t = np.zeros(N)\n", " y[0] = y0\n", " t[0] = t0\n", " for n in range(1,N):\n", " y[n] = y[n-1] + h*f(t[n-1],y[n-1])\n", " t[n] = t[n-1] + h\n", " return t, y" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 4 }, { "cell_type": "code", "collapsed": false, "input": [ "t0 = 0\n", "T = 10\n", "h = 1.0/20.0\n", "t,y = euler(t0,T,0,h)\n", "te = np.linspace(t0,T,100)\n", "ye = yexact(te)\n", "plt.plot(t,y,te,ye,'--')\n", "plt.legend(('Numerical','Exact'))\n", "plt.xlabel('t')\n", "plt.ylabel('y')\n", "plt.title('Step size = ' + str(h))" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 5, "text": [ "" ] }, { "metadata": {}, "output_type": "display_data", "svg": [ "\n", "\n", "\n", "\n" ], "text": [ "" ] } ], "prompt_number": 5 }, { "cell_type": "markdown", "metadata": {}, "source": [ "Study the effect of decreasing step size. The error is plotted in log scale." ] }, { "cell_type": "code", "collapsed": false, "input": [ "hh = [1.0/2.0, 1.0/10.0, 1.0/50.0]\n", "for h in hh:\n", " t,y = euler(t0,T,0,h)\n", " ye = yexact(t)\n", " plt.semilogy(t,np.abs(y-ye))\n", " plt.legend(hh)\n", " plt.xlabel('t')\n", " plt.ylabel('log(error)')" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "display_data", "svg": [ "\n", "\n", "\n", "\n" ], "text": [ "" ] } ], "prompt_number": 6 }, { "cell_type": "code", "collapsed": false, "input": [], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 6 } ], "metadata": {} } ] }