CDMATH/tests/doc/1DTransportEquation/UnstableSchemes/TransportEquation1D_UnstableSchemes.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Examples of unstable finite volume approximations of the 1D transport equation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The transport equation with periodic boundary conditions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We are interested in the finite volume approximation of the following partial differential equation  \n",
    "$$\n",
    "\\partial_t u + c \\partial_x u=0,\n",
    "$$  \n",
    "on the 1D domain $[0,1]$ with periodic boundary condition at $x=0$ and $x=1$ and initial data  \n",
    "$$\n",
    "u_0(x)=\\frac{1}{2} (1+sin(\\pi(4x-0.5)) 1_{[0,0.5]} + 1_{[0.6,0.85]}.\n",
    "$$  \n",
    "\n",
    "The initial data consists in a smooth part ($x\\in [0,0.5]$) and a stiff part ($x\\in [0.5,1]$).  \n",
    "\n",
    "The exact solution is given by \n",
    "$$\n",
    "u(x,t)=u_0(x-ct).\n",
    "$$  \n",
    "\n",
    "Since we sused periodic boundary condition, and whatever the initial data, the exact solution is periodic with period $T=\\frac{1}{c}$ and therefore \n",
    "$$\n",
    "u(x,T)=u_0(x).\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Finite volume approximations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In 1D finite volume approximations, the domain $\\Omega=[0,1]$ is decomposed into $N$ intervals $C_i=[x_i,x_{i+1}], i=1,\\dots N$, and we seek the average values  \n",
    "\n",
    "$$\n",
    "u_i(t) = \\frac{1}{x_{i+1} - x_i}\\int_{x_i}^{x_{i+1}} u(x,t) dx\n",
    "$$\n",
    "of the exact solution $u(x,t)$ in each cell $C_i$.  \n",
    "Similarly we decompose the time domain $\\mathbb{R}_+$ into finite length intervals $[t_n,t_{n+1}]$.\n",
    "Denotig $\\Delta t_n=t_{n+1}-t_n$ the time step and $\\Delta x_i=x_{i+1} - x_i$ the space step, the double integration  \n",
    "\n",
    "$$\n",
    "\\begin{array}{ccc}\\displaystyle\n",
    "\\frac{1}{\\Delta x_i}\\int_{x_i}^{x_{i+1}}\\frac{1}{\\Delta t_n}\\int_{t_n}^{t_{n+1}}(\\partial_t u + c \\partial_x u)\\, dt\\, dx\n",
    "&=&\\displaystyle\n",
    "\\frac{1}{\\Delta x_i}\\int_{x_i}^{x_{i+1}}\\frac{1}{\\Delta t_n}\\int_{t_n}^{t_{n+1}}\\partial_t u(x,t)\\, dt\\, dx\n",
    "+ c \n",
    "\\frac{1}{\\Delta t_n}\\int_{t_n}^{t_{n+1}}\\frac{1}{\\Delta x_i}\\int_{x_i}^{x_{i+1}}\\partial_x u(x,t)\\, dx\\, dt\\\\[1.5ex]\n",
    " &=&\\displaystyle \\frac{1}{\\Delta x_i}\\int_{x_i}^{x_{i+1}} \\frac{u(x,t_{n+1}) - u(x,t_n)}{\\Delta t_n}\\, dx\n",
    " + c\n",
    " \\frac{1}{\\Delta t_n}\\int_{t_n}^{t_{n+1}}\\frac{u(x_{i+1},t) - u(x_i,t)}{\\Delta x_i}\\, dt,\n",
    "\\end{array}\n",
    "$$  \n",
    "\n",
    "yields\n",
    "$$\n",
    "\\frac{u_i^{n+1}-u_i^n}{\\Delta t_n} + c \\frac{u_{i+\\frac{1}{2}}^{n,n+1}-u_{i-\\frac{1}{2}}^{n,n+1}}{\\Delta x_i}=0,\n",
    "$$  \n",
    "\n",
    "where\n",
    "$\n",
    "u_i^n = \\frac{1}{\\Delta x_i}\\int_{x_i}^{x_{i+1}} u(x,t_n)\\, dx\n",
    "$\n",
    "is the average of the exact solution on the space interval $[x_i,x_{i+1}]$ at time $t=t_n$,  \n",
    "and $u_{i+\\frac{1}{2}}^{n,n+1}$ is the average of the exact solution on the time interval $[t_n,t_{n+1}]$ at position $x=x_{i+1}$.  \n",
    "In order to close the systemne needs to express $u_{i+\\frac{1}{2}^{n,n+1}}$ using the unknown variables $u_i^n, i=1,\\dots N, n=0,\\infty$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###  The upwind explicit scheme with cfl>1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the upwind explicit scheme we take \n",
    "$$\n",
    "\\begin{array}{ccc}\n",
    "u_{i+\\frac{1}{2}^{n,n+1}} &=& u_i^n     \\textrm{ if } c\\geq 0\\\\\n",
    "u_{i+\\frac{1}{2}^{n,n+1}} &=& u_{i-1}^n \\textrm{ if } c\\leq 0\n",
    "\\end{array}\n",
    "$$\n",
    "and assuming $c>0$, we obtain the discrete equation\n",
    "$$\n",
    "\\frac{u_i^{n+1}-u_i^n}{\\Delta t_n} + c \\frac{u_{i}^{n}-u_{i-1}^{n}}{\\Delta x_i}=0.\n",
    "$$  \n",
    "The scheme is stable provided $\\Delta t_n \\leq\\Delta t_{stab}=\\frac{\\min_i\\Delta x_i}{c}$, and is able to capture the exact solution without numerical diffusion.  \n",
    "We however choose a large time step $\\Delta t_n =2\\Delta t_{stab}=2\\frac{\\min_i\\Delta x_i}{c}$ and observe the apparition of instabilities."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We choose an initial data consisting in a smooth and a stiff part.\n",
    "![\"The initial data for the unstable (cfl>1) explicit scheme\"](TransportEquation_UpwindExplicit_50Cells_CFL2_ResultField_0.png)\n",
    "We quickly observe from the first time step instabilities emerging from the stiff part."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![\"Strong instabilities obtained with the upwind explicit scheme with cfl>1\"](TransportEquation_UpwindExplicit_50Cells_CFL2_ResultField_2.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The instabilities grow and propagate to the smooth part.\n",
    "![\"Strong instabilities obtained with the upwind explicit scheme with cfl>1\"](TransportEquation_UpwindExplicit_50Cells_CFL2_ResultField_11.png)\n",
    "The numerical solution diverges."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "scrolled": true
   },
   "source": [
    "<div align=\"middle\">\n",
    "<video width=\"80%\" controls>\n",
    "      <source src=\"1DTransportEquation_UpwindExplicit_100Cells_CFL2.0.mp4\" type=\"video/mp4\">\n",
    "</video></div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[\"Click here for a video animation of the unstable (cfl>1) upwind explicit scheme\"](1DTransportEquation_UpwindExplicit_100Cells_CFL2.0.mp4)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###  The centered explicit scheme"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the centered explicit scheme we take \n",
    "$$\n",
    "\\begin{array}{ccc}\n",
    "u_{i+\\frac{1}{2}^{n,n+1}} &=& \\frac{u_{i+1}^{n}+u_i^{n}}{2}\n",
    "\\end{array}\n",
    "$$\n",
    "and we obtain the discrete equation\n",
    "$$\n",
    "\\frac{u_i^{n+1}-u_i^n}{\\Delta t_n} + c \\frac{u_{i+1}^{n}-u_{i-1}^{n}}{2\\Delta x_i}=0.\n",
    "$$  \n",
    "The scheme is unstable for any $\\Delta t_n > 0$. It introduces a numerical antidiffusion. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Stiff initial data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We start with a stiff initial data\n",
    "![\"Stiff initial data\"](TransportEquation_CenteredExplicit_50Cells_Stiff_CFL0.99_ResultField_0.png)\n",
    "We quickly (10th time step) see a violent instability"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![\"Strong instabilities obtained with the centered explicit scheme (stiff initial data)\"](TransportEquation_CenteredExplicit_50Cells_Stiff_CFL0.99_ResultField_10.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We observe strong spurious oscillations that emerge from the stiffest parts of the initial data as can be seen on the picture above. The numerical solution quickly diverges."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "scrolled": true
   },
   "source": [
    "<div align=\"middle\">\n",
    "<video width=\"80%\" controls>\n",
    "      <source src=\"1DTransportEquation_CenteredExplicit_50Cells_Stiff_CFL0.99.mp4\" type=\"video/mp4\">\n",
    "</video></div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[\"Click here for a video animation of the unstable centered explicit scheme with stiff intial data\"](1DTransportEquation_CenteredExplicit_50Cells_Stiff_CFL0.99.mp4)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Smooth initial data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We start with a smooth initial data\n",
    "![\"Stiff initial data\"](TransportEquation_CenteredExplicit_50Cells_Smooth_CFL0.99_ResultField_0.png)\n",
    "After a while (110th time step), numerical instabilities eventually emerge"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![\"Instabilities obtained with the centered explicit scheme (smooth initial data)\"](TransportEquation_CenteredExplicit_50Cells_Smooth_CFL0.99_ResultField_110.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We observe growing spurious oscillations and the numerical solution starts diverging."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "scrolled": true
   },
   "source": [
    "<div align=\"middle\">\n",
    "<video width=\"80%\" controls>\n",
    "      <source src=\"1DTransportEquation_CenteredExplicit_50Cells_Smooth_CFL0.99.mp4\" type=\"video/mp4\">\n",
    "</video></div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[\"Click here for a video animation of the unstable centered explicit scheme with smooth intial data\"](1DTransportEquation_CenteredExplicit_50Cells_Smooth_CFL0.99.mp4)"
   ]
  }
 ],
 "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.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}