Chapter18.ipynb

{
 "metadata": {
  "name": "",
  "signature": "sha256:c191052e4e716dc57a851e0779df4a7b41e07f15fcc8db18467801e675da7cb3"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Partial Differential Equations\n",
      "\n",
      "Many physical phenomena are a function of more than one independent variable and must be represented by partial-differential equations (PDE's). Most scientific problems have mathematical models that are second-order equations, with the highest order of derivative being the second. If $u$ is a function of the two independent variables $x$ and $y$, there are three second-order partial derivatives:\n",
      "$$\n",
      "\\frac{\\partial ^2 u}{\\partial x^2}, \\hspace{5pt} \\frac{\\partial^2}{\\partial x\\partial y}, \\hspace{5pt} \\frac{\\partial^2u}{\\partial y^2}\n",
      "$$\n",
      "\n",
      "In general an equation involving these partial derivatives might have the form:\n",
      "\n",
      "$$\n",
      "A\\frac{\\partial ^2 u}{\\partial x^2} + B\\frac{\\partial^2}{\\partial x\\partial y} + C\\frac{\\partial^2u}{\\partial y^2} + D\\left(x,y,u,\\frac{\\partial u}{\\partial x}, \\frac{\\partial u}{\\partial y}\\right) = 0\n",
      "$$\n",
      "\n",
      "Depending on the values of the coefficients of the second-derivative terms, partial-differential equations are classified as elliptic, parabolic, or hyperbolic according to the following:\n",
      "\n",
      "$$\n",
      "\\begin{matrix}\n",
      "B^2 - 4AC < 0 \\rightarrow &\\text{Elliptic}\\\\\n",
      "B^2 - 4AC = 0 \\rightarrow &\\text{Parabolic}\\\\\n",
      "B^2 - 4AC > 0 \\rightarrow &\\text{Hyperbolic}\\end{matrix}\n",
      "$$\n",
      "\n",
      "It is very important to understand the characteristic of the PDE because this will influence our solution technique. Notice that if the coefficients are functions of either independent variable, it's possible for the characteristic of the PDE's to change throughout the domain or during the solution time. We will begin our discussion of solution techniques with a very common hyperbolic PDE."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## The Wave Equation\n",
      "\n",
      "The wave equation has the form:\n",
      "\n",
      "$$\n",
      "\\frac{\\partial^2u}{\\partial t^2} - c^2 \\frac{\\partial u^2}{\\partial x^2} = 0\n",
      "$$\n",
      "\n",
      "We can use *python* to solve this equation analytically:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "simplify[Dsolve[{D[u[x,t],{t,2}]-c**2D[u[x,t],{x,2}] ==0},u[x,t],{x,t}],assumptions->{c>0}]//first"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "ename": "SyntaxError",
       "evalue": "invalid syntax (<ipython-input-2-558699e15459>, line 1)",
       "output_type": "pyerr",
       "traceback": [
        "\u001b[0;36m  File \u001b[0;32m\"<ipython-input-2-558699e15459>\"\u001b[0;36m, line \u001b[0;32m1\u001b[0m\n\u001b[0;31m    simplify[Dsolve[{D[u[x,t],{t,2}]-c**2D[u[x,t],{x,2}] ==0},u[x,t],{x,t}],assumptions->{c>0}]//first\u001b[0m\n\u001b[0m                                         ^\u001b[0m\n\u001b[0;31mSyntaxError\u001b[0m\u001b[0;31m:\u001b[0m invalid syntax\n"
       ]
      }
     ],
     "prompt_number": 2
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Here the $C_1$ and $C_2$ terms are functions of right-moving and left-moving waves. Any sum of two functions of this form that satisfies the boundary conditions are found to be solutions to the wave equation. More commonly you will see the solution to the wave equations written as:\n",
      "\n",
      "$$\n",
      "u[x,t] = f(x+ct)+g(x-ct)\n",
      "$$\n",
      "\n",
      "One method for solving the wave equation numerically is with a finite-difference scheme, just as we did for boundary value problems. Here we have two independent variables to keep track of so we will introduce some new notation: we will use superscripts to denote time, and subscripts for position. Substituting in the central difference approximations to the derivatives into the wave equations we have:\n",
      "\n",
      "$$\n",
      "\\frac{u_i^{j+1}-2u_i^j+u_i^{j-1}}{(\\Delta t)^2} - c^2 \\left(\\frac{u_{i+1}^j - 2y_i^j + y_{i-1}^j}{(\\Delta x)^2}\\right) = 0\n",
      "$$\n",
      "\n",
      "Solving for the displacement at the end of the current interval, at $t = t_{j+1}$ we get \n",
      "\n",
      "$$\n",
      "u_i^{j+1} = \\left(c\\frac{\\Delta t}{\\Delta x}\\right)^2\\left(u_{i+1}^j + u_{i-1}^j\\right)-u_i^{j-1}+2\\left(1-\\left(c\\frac{\\Delta t}{\\Delta x}\\right)^2\\right)u_i^j\n",
      "$$\n",
      "\n",
      "Note that if we let the ratio $\\left(c \\frac{\\Delta t}{\\Delta x}\\right)$ equal unity, we get some simplification:\n",
      "\n",
      "$$\n",
      "u_i^{j+1} = u_{i+1}^j + u_{i-1}^j - u_i^{j-1}, \\hspace{5pt} \\Delta t = \\frac{\\Delta x}{c}\n",
      "$$\n",
      "\n",
      "This is the most common numerical method for solving the wave equation. We will now investigate our choice of letting $\\left(c\\frac{\\Delta t}{\\Delta x}\\right) = 1$."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Comparison to the d'Alembert Solution\n",
      "\n",
      "The method of d'Alembert lets us find the analytical solution to the wave equation in one dimension. By direct substitution, it is readily seen that, for any arbitrary function $f$ and $g$ the wave equation is solved by:\n",
      "\n",
      "$$\n",
      "u[x,t] = f(x+ct)+g(x-ct)\n",
      "$$\n",
      "\n",
      "This is demonstrated easily:\n",
      "\n",
      "$$\n",
      "\\frac{\\partial u}{\\partial t} = f' \\frac{\\partial(x+ct)}{\\partial t} + g' \\frac{\\partial(x-ct)}{\\partial t} = cf' -cg'\\\\\n",
      "\\frac{\\partial^2u}{\\partial t^2}=c^2f''-c(-c)g''=c^2f''+c^2g''\\\\\n",
      "\\frac{\\partial u}{\\partial x} = f' \\frac{\\partial(x+ct)}{\\partial x} + g' \\frac{\\partial(x-ct)}{\\partial x} = f' + g'\\\\\n",
      "\\frac{\\partial^2u}{\\partial x^2} = f'' + g''\n",
      "$$\n",
      "\n",
      "Substituting these into the wave equation we have:\n",
      "\n",
      "$$\n",
      "\\frac{\\partial^2u}{\\partial t^2} - c^2 \\frac{\\partial u^2}{\\partial x^2} = c^2(f''+g'')-c^2(f''+g'')=0\n",
      "$$\n",
      "\n",
      "The solution to the problem is found, then, if we can find a pair of functions of the form whose sum matches the initial conditions and the boundary conditions for the problem. If the initial conditions are:\n",
      "\n",
      "$$\n",
      "y(x,0)=f(x)\\\\\n",
      "\\frac{\\partial u}{\\partial t}(x,0) = g(x)\n",
      "$$\n",
      "\n",
      "then $u[x,t]$ is given by\n",
      "\n",
      "$$\n",
      "u(x,t) = \\frac{1}{2}[f(x+ct)+f(x-ct)]+\\frac{1}{2c}\\int_{x-ct}^{x+ct}g(v)dv\n",
      "$$\n",
      "\n",
      "We are now in a position to verify that our numerical procedure gives the *exact* solution (except for roundoff), provided that two lines of correct $u$-values are known. We show that the solution ot the difference equation is a solution to the differential equation. Consider the function\n",
      "\n",
      "$$\n",
      "u_i^j = f(x_i + ct_j) + g(x_i - ct_j) \\hspace{5pt} (*)\n",
      "$$\n",
      "\n",
      "Because we chose $\\left(c\\frac{\\Delta t}{\\Delta x}\\right) = 1$, then $\\Delta x = c\\Delta t$. If we write $x_i$ and $t_j$ in terms of starting values $x_0$, and $t_0$:\n",
      "\n",
      "$$\n",
      "x_i = x_0 + i\\Delta x \\\\\n",
      "ct_j = c(t_0 + j\\Delta t) = cj \\Delta t = j \\Delta x\n",
      "$$\n",
      "\n",
      "substituting back into $(*)$ gives\n",
      "\n",
      "$$\n",
      "u_i^j = f(x_0 + i\\Delta x + j\\Delta x) + g(x_0 + i\\Delta x - j\\Delta x)\\\\\n",
      "= f(x_0 + (i+j)\\Delta x) + g(x_0 + (i-j)\\Delta x)\n",
      "$$\n",
      "\n",
      "now substitute back into the difference equation that defines our numerical procedure we have:\n",
      "\n",
      "$$\n",
      "\\begin{matrix} u_i^{j+1} = u_{i-1}^j + u_{i+1}^j - u_i^{j-1} & = & f(x_0 + (i-1 + j)\\Delta x)\\\\ & + & g(x_0 + (i-1-j)\\Delta x)\\\\ & + & f(x_0+(i+1+j)\\Delta x) \\\\ & + & g(x_0+(i+1-j)\\Delta x)\\\\ & - & f(x_0+(i+j-1)\\Delta x)\\\\ & - & g(x_0+(i-j+1)\\Delta x)\\end{matrix}\\\\\n",
      "= f(x_0 + (i+j+1)\\Delta x)+g(x_0+(i-j-1)\\Delta x)\n",
      "$$\n",
      "\n",
      "What this result illustrates is that if two previous lines of the numerical solution are of the form of $(*)$, and hence (exact solutions of the wave equation), it follows that the values on the next line are exact solutions also. Seems like our choice of $\\left(c\\frac{\\Delta t}{\\Delta x}\\right)=1$ was pretty good.\n",
      "\n",
      "We have seen that if two lines in the numerical solution are exact, then all subsequent lines can be exact as well (neglecting roundoff, of course). The first line will always be exact because it is the initially given values. From what we learned previously we know a very good way to get the second line of $u$ values. Using the solution that was found with the given velocity boundary condition and employing it at $t=\\Delta t$, we have \n",
      "\n",
      "$$\n",
      "u[x_i,t_i]=y_i^1 = \\frac{1}{2}[f(x_i+\\Delta x)+f(x_i - \\Delta x)]+\\frac{1}{2c}\\int_{x_i-\\Delta x}^{x_i+\\Delta x} g(v) dv\\\\\n",
      "= \\frac{1}{2} [y_{i+1}^0 + y_{i-1}^0]+\\frac{1}{2c}\\int_{x_i-1}^{x_i+1}g(v)dv\n",
      "$$\n",
      "\n",
      "If the initial velocity is zero, then $g(x)$ is zero everywhere and the second term above becomes irrelevant. We will always get exact answers using our numerical procedure in these cases."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": []
    }
   ],
   "metadata": {}
  }
 ]
}