{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Advection using TVD Flux Limiters\n", "## CH EN 6355 - Computational Fluid Dynamics\n", "**Prof. Tony Saad (www.tsaad.net)
\n", "slides at: www.tsaad.net
\n", "Department of Chemical Engineering
\n", "University of Utah**\n", "
" ] }, { "cell_type": "markdown", "metadata": { "variables": { "1 + k": {}, "1 - k": {} } }, "source": [ "Here, we will implement the k-scheme or kappa-schemes for advection. It is easiest to implement this scheme since for different values of k, we recover all sorts of high-order flux approximations. We will assume a positive advecting velocity for illustration purposes.\n", "\n", "We are solving the constant speed advection equation given by\n", "\\begin{equation}\n", "u_t = - c u_x = - F_x;\\quad F = cu\n", "\\end{equation}\n", "We will use a simple Forward Euler explicit method. Using a finite volume integration, we get\n", "\\begin{equation}\n", "u_i^{n+1} = u_i^n - \\frac{\\Delta t}{\\Delta x} (F_{i+\\tfrac{1}{2}}^n - F_{i-\\tfrac{1}{2}}^n)\n", "\\end{equation}\n", "Our TVD formulation is based on a finite volume discretization with face values for the dependent variable given by:\n", "\\begin{equation}\n", "{\\phi _f} = {\\phi _{\\rm{C}}} + {1 \\over 2}\\psi ({r_f})\\left( {{\\phi _{\\rm{D}}} - {\\phi _{\\rm{C}}}} \\right)\n", "\\end{equation}\n", "where $\\psi$ is the limiter function." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "%matplotlib inline\n", "%config InlineBackend.figure_format = 'svg'\n", "import matplotlib.pyplot as plt\n", "import matplotlib.animation as animation\n", "plt.rcParams['animation.html'] = 'html5'\n", "from matplotlib import cm" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "def step(x,x0):\n", " x0 = 0.6\n", " x1 = 0.8\n", " result = x - x0\n", " result[x-x1x1] = 0.0 \n", " return result\n", "\n", "def gaussian(x,x0):\n", " s = 0.08\n", " s = s*s\n", " result = np.exp( -(x-x0)**2/s)\n", " return result" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/svg+xml": [ "\n", "\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n" ], "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "L = 1.0\n", "n = 256 # cells\n", "dx = L/n # n intervals\n", "x = np.linspace(-3*dx/2, L + 3*dx/2, n+4) # include ghost cells - we will include 2 ghost cells on each side for high order schemes\n", "\n", "# create arrays\n", "phi = np.zeros(n+4) # cell centered quantity\n", "f = np.zeros(n+4+1) # flux\n", "u = np.ones(n+4+1) # velocity field - assumed to live on faces same as flux\n", "\n", "x0 = 0.3\n", "# u0 = np.zeros(N + 2)\n", "# u0[1:-1] = np.sin(2*np.pi*x)\n", "# u0 = np.zeros(N)\n", "# phi0 = np.sin(np.pi*x)\n", "phi0 = gaussian(x,x0) + step(x,x0)\n", "# u0 = triangle(x,0.5,0.75,1)\n", "# u0[0:N//2] = 1.0\n", "plt.plot(x,phi0)" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "dt= 0.001953125\n", "dx= 0.00390625\n" ] } ], "source": [ "cfl = 0.5\n", "c = 1.0\n", "dt = cfl*dx/abs(c)\n", "print('dt=',dt)\n", "print('dx=',dx)\n", "# μ = 0.001\n", "# dt = 0.02" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n" ], "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# finite volume implementation with arrays for fluxes\n", "t = 0\n", "tend = L/abs(c)\n", "sol = []\n", "sol.append(phi0)\n", "ims = []\n", "\n", "fig = plt.figure(figsize=[5,3],dpi=200)\n", "plt.rcParams[\"font.family\"] = \"serif\"\n", "plt.rcParams[\"font.size\"] = 10\n", "plt.rc('text', usetex=True)\n", "\n", "# plt.grid()\n", "plt.xlim([0.,L])\n", "plt.ylim([-0.25,1.25])\n", "plt.xlabel('$x$')\n", "plt.ylabel('$\\phi$')\n", "plt.tight_layout()\n", "# plot initial condition\n", "plt.plot(x,phi0,'darkred',animated=True)\n", "\n", "schemename='Bounded CD'\n", "\n", "i = 0\n", "while t < tend: \n", " phin = sol[-1]\n", "\n", " # impose periodic conditions\n", " phin[-2] = phin[2]\n", " phin[-1] = phin[3] \n", " phin[0] = phin[-4] \n", " phin[1] = phin[-3] \n", "\n", "\n", " if (i%2==0):\n", " shift =int(np.ceil(c*(t-dt)/dx))\n", " im = plt.plot(x[2:-2], np.roll(phin[2:-2], -shift) ,'k-o',markevery=2,markersize=3.5,markerfacecolor='deepskyblue',\n", " markeredgewidth=0.25, markeredgecolor='k',linewidth=0.45, animated=True)\n", " ims.append(im)\n", " \n", " phi = np.zeros_like(phi0)\n", " \n", " # predictor - take half a step and use upwind\n", " # du/dt = -c*du/dx\n", " if c >= 0:\n", " ϕc = phin[1:-2] # phi upwind\n", " else:\n", " ϕc = phin[2:-1] # phi upwind\n", " \n", " f[2:-2] = c*ϕc\n", " phi[2:-2] = phin[2:-2] - dt/2.0/dx*(f[3:-2] - f[2:-3])\n", " phi[-2] = phi[2]\n", " phi[-1] = phi[3] \n", " phi[0] = phi[-4] \n", " phi[1] = phi[-3] \n", "\n", " # du/dt = -c*du/dx\n", " if c >= 0:\n", " ϕc = phi[1:-2] # phi upwind\n", " ϕu = phi[:-3] # phi far upwind\n", " ϕd = phi[2:-1] # phi downwind\n", " else:\n", " ϕc = phi[2:-1] # phi upwind\n", " ϕu = phi[3:] # phi far upwind\n", " ϕd = phi[1:-2] # phi downwind\n", "\n", " # compute r - ratio of successive gradients\n", " numerator = ϕc - ϕu \n", " denominator = ϕd - ϕc\n", " # find the locations where the denominator is > 1e-12 - that's where we will do the divisions. otherwise, set phi to 1\n", " divideloc = np.where(np.abs(denominator) > 1e-6)\n", " i1 = np.where(numerator < 1e-12)\n", " r = np.ones_like(ϕc)\n", " r[divideloc] = numerator[divideloc]/denominator[divideloc]\n", "\n", " # compute face values - these actually live on the faces\n", " ϕf = np.zeros_like(ϕc)\n", " \n", "# # bounded CD \n", "# schemename='Bounded CD'\n", "# ϕftilde = 0.5 + 0.5*ϕctilde\n", " \n", "# # minmod\n", "# schemename='MinMod'\n", "# psi = np.maximum(0,np.minimum(1,r))\n", " \n", " # SUPERBEE\n", " schemename='Superbee'\n", " a1 = np.maximum(0, np.minimum(1.0,2.0*r))\n", " psi = np.maximum(a1, np.minimum(2.0,r))\n", "\n", "# # MUSCL\n", "# schemename='MUSCL'\n", "# psi = np.maximum(0.0, np.minimum(np.minimum(2*r, (r+1)/2.0), 2.0))\n", "\n", " # now retrieve phif from phiftilde\n", " ϕf = ϕc + 0.5*psi*(ϕd - ϕc)\n", " f[2:-2] = ϕf\n", " f = c*f # multiply the flux by the velocity\n", " # advect\n", " phi[2:-2] = phin[2:-2] - dt/dx*(f[3:-2] - f[2:-3]) #+ dt/dx/dx*diffusion\n", " t += dt \n", " i+=1\n", " sol.append(phi)\n", "\n", "\n", "plt.annotate(schemename, xy=(0.5, 0.8), xytext=(0.015, 0.9),fontsize=8)\n", "plt.legend(('exact','numerical'),loc='upper left',fontsize=7)\n", "ani = animation.ArtistAnimation(fig, ims, interval=100, blit=True,\n", " repeat_delay=1000)\n", "\n", "# ani.save('k-scheme-'+str(k)+'.mp4',dpi=300,fps=24)" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n" ], "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.plot(sol[0], label='initial condition')\n", "plt.plot(sol[-1], label='one residence time')\n", "plt.legend()\n", "plt.grid()" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ani" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Create Animation in Moving Reference Frame" ] }, { "cell_type": "code", "execution_count": 19, "metadata": { "scrolled": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "done!\n" ] }, { "data": { "image/svg+xml": [ "\n", "\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n" ], "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "\"\"\"\n", "Create Animation in Moving Reference Frame\n", "\"\"\"\n", "import matplotlib\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import matplotlib.animation as animation\n", "matplotlib.use(\"Agg\")\n", "fig, ax = plt.subplots(figsize=(4,3),dpi=150) \n", "ax.grid(True,color='0.9')\n", "f0 = sol[0]\n", "line0, = ax.plot(x[2:-2], f0[2:-2] ,'r-',linewidth=0.75, animated=True)\n", "line1, = ax.plot(x[2:-2], f0[2:-2] ,'k-o',markevery=2,markersize=3.5,markerfacecolor='deepskyblue',\n", " markeredgewidth=0.25, markeredgecolor='k',linewidth=0.45, animated=True)\n", "\n", "\n", "ann = ax.annotate('time ='+str(round(t,4))+' s.' + '\\t'+ schemename, xy=(0,0), xytext=(40, 200),xycoords='figure points')\n", "plt.tight_layout()\n", "# plt.locator_params(axis='x', nbins=15)\n", "\n", "# def animate(i):\n", "# # xt = x - i*c*dt/dx\n", "# # line.set_xdata(xt)\n", "# # ax.axes.set_xlim(xt[0],xt[-1])\n", "# # ax.grid()\n", "# # line.set_ydata(np.sin(xt)) # update the data\n", "# # line.set_ydata(sol[i])\n", "# print('time=',i*dt)\n", "# shift =int(np.ceil(i*c*dt/dx))\n", "# f = sol[i]\n", "# line.set_ydata(np.roll(f, -shift))\n", "# # im = plt.plot(x[2:-2], np.roll(phin[2:-2], -shift) ,'k-o',markevery=2,markersize=3.5,markerfacecolor='deepskyblue',\n", "# # markeredgewidth=0.25, markeredgecolor='k',linewidth=0.45, animated=True)\n", "# return line,\n", "\n", "def animate_moving(i):\n", " print('time=',i*dt)\n", " t = i*dt\n", " xt = x + i*1.1*c*dt\n", " line0.set_xdata(xt[2:-2])\n", " line1.set_xdata(xt[2:-2]) \n", " ax.axes.set_xlim(xt[0],0.0*dx + xt[-1])\n", " f = sol[i]\n", " ax.axes.set_ylim(-0.1,1.1*max(f))\n", " ann.set_text('time ='+str(round(t,4))+'s (' + str(i)+ ').\\t' + schemename)\n", " shift =int(np.ceil(i*c*dt/dx))\n", " line1.set_ydata(np.roll(f[2:-2], -shift))\n", "\n", " f0 = sol[0]\n", " line0.set_ydata(f0[2:-2])\n", " return line0,line1\n", "\n", "\n", "# Init only required for blitting to give a clean slate.\n", "def init():\n", " line0.set_ydata(np.ma.array(x[2:-2], mask=True))\n", " line1.set_ydata(np.ma.array(x[2:-2], mask=True)) \n", " return line0,line1\n", "\n", "ani = animation.FuncAnimation(fig, animate_moving, np.arange(0,len(sol),2*int(1.0/cfl)), init_func=init,\n", " interval=20, blit=False)\n", "print('done!')" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "time= 0.0\n", "time= 0.0078125\n", "time= 0.015625\n", "time= 0.0234375\n", "time= 0.03125\n", "time= 0.0390625\n", "time= 0.046875\n", "time= 0.0546875\n", "time= 0.0625\n", "time= 0.0703125\n", "time= 0.078125\n", "time= 0.0859375\n", "time= 0.09375\n", "time= 0.1015625\n", "time= 0.109375\n", "time= 0.1171875\n", "time= 0.125\n", "time= 0.1328125\n", "time= 0.140625\n", "time= 0.1484375\n", "time= 0.15625\n", "time= 0.1640625\n", "time= 0.171875\n", "time= 0.1796875\n", "time= 0.1875\n", "time= 0.1953125\n", "time= 0.203125\n", "time= 0.2109375\n", "time= 0.21875\n", "time= 0.2265625\n", "time= 0.234375\n", "time= 0.2421875\n", "time= 0.25\n", "time= 0.2578125\n", "time= 0.265625\n", "time= 0.2734375\n", "time= 0.28125\n", "time= 0.2890625\n", "time= 0.296875\n", "time= 0.3046875\n", "time= 0.3125\n", "time= 0.3203125\n", "time= 0.328125\n", "time= 0.3359375\n", "time= 0.34375\n", "time= 0.3515625\n", "time= 0.359375\n", "time= 0.3671875\n", "time= 0.375\n", "time= 0.3828125\n", "time= 0.390625\n", "time= 0.3984375\n", "time= 0.40625\n", "time= 0.4140625\n", "time= 0.421875\n", "time= 0.4296875\n", "time= 0.4375\n", "time= 0.4453125\n", "time= 0.453125\n", "time= 0.4609375\n", "time= 0.46875\n", "time= 0.4765625\n", "time= 0.484375\n", "time= 0.4921875\n", "time= 0.5\n", "time= 0.5078125\n", "time= 0.515625\n", "time= 0.5234375\n", "time= 0.53125\n", "time= 0.5390625\n", "time= 0.546875\n", "time= 0.5546875\n", "time= 0.5625\n", "time= 0.5703125\n", "time= 0.578125\n", "time= 0.5859375\n", "time= 0.59375\n", "time= 0.6015625\n", "time= 0.609375\n", "time= 0.6171875\n", "time= 0.625\n", "time= 0.6328125\n", "time= 0.640625\n", "time= 0.6484375\n", "time= 0.65625\n", "time= 0.6640625\n", "time= 0.671875\n", "time= 0.6796875\n", "time= 0.6875\n", "time= 0.6953125\n", "time= 0.703125\n", "time= 0.7109375\n", "time= 0.71875\n", "time= 0.7265625\n", "time= 0.734375\n", "time= 0.7421875\n", "time= 0.75\n", "time= 0.7578125\n", "time= 0.765625\n", "time= 0.7734375\n", "time= 0.78125\n", "time= 0.7890625\n", "time= 0.796875\n", "time= 0.8046875\n", "time= 0.8125\n", "time= 0.8203125\n", "time= 0.828125\n", "time= 0.8359375\n", "time= 0.84375\n", "time= 0.8515625\n", "time= 0.859375\n", "time= 0.8671875\n", "time= 0.875\n", "time= 0.8828125\n", "time= 0.890625\n", "time= 0.8984375\n", "time= 0.90625\n", "time= 0.9140625\n", "time= 0.921875\n", "time= 0.9296875\n", "time= 0.9375\n", "time= 0.9453125\n", "time= 0.953125\n", "time= 0.9609375\n", "time= 0.96875\n", "time= 0.9765625\n", "time= 0.984375\n", "time= 0.9921875\n", "time= 1.0\n" ] } ], "source": [ "ani.save(schemename+'_TVD.mp4',fps=24,dpi=200)" ] }, { "cell_type": "code", "execution_count": 21, "metadata": { "scrolled": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "time= 0.0\n", "time= 0.0078125\n", "time= 0.015625\n", "time= 0.0234375\n", "time= 0.03125\n", "time= 0.0390625\n", "time= 0.046875\n", "time= 0.0546875\n", "time= 0.0625\n", "time= 0.0703125\n", "time= 0.078125\n", "time= 0.0859375\n", "time= 0.09375\n", "time= 0.1015625\n", "time= 0.109375\n", "time= 0.1171875\n", "time= 0.125\n", "time= 0.1328125\n", "time= 0.140625\n", "time= 0.1484375\n", "time= 0.15625\n", "time= 0.1640625\n", "time= 0.171875\n", "time= 0.1796875\n", "time= 0.1875\n", "time= 0.1953125\n", "time= 0.203125\n", "time= 0.2109375\n", "time= 0.21875\n", "time= 0.2265625\n", "time= 0.234375\n", "time= 0.2421875\n", "time= 0.25\n", "time= 0.2578125\n", "time= 0.265625\n", "time= 0.2734375\n", "time= 0.28125\n", "time= 0.2890625\n", "time= 0.296875\n", "time= 0.3046875\n", "time= 0.3125\n", "time= 0.3203125\n", "time= 0.328125\n", "time= 0.3359375\n", "time= 0.34375\n", "time= 0.3515625\n", "time= 0.359375\n", "time= 0.3671875\n", "time= 0.375\n", "time= 0.3828125\n", "time= 0.390625\n", "time= 0.3984375\n", "time= 0.40625\n", "time= 0.4140625\n", "time= 0.421875\n", "time= 0.4296875\n", "time= 0.4375\n", "time= 0.4453125\n", "time= 0.453125\n", "time= 0.4609375\n", "time= 0.46875\n", "time= 0.4765625\n", "time= 0.484375\n", "time= 0.4921875\n", "time= 0.5\n", "time= 0.5078125\n", "time= 0.515625\n", "time= 0.5234375\n", "time= 0.53125\n", "time= 0.5390625\n", "time= 0.546875\n", "time= 0.5546875\n", "time= 0.5625\n", "time= 0.5703125\n", "time= 0.578125\n", "time= 0.5859375\n", "time= 0.59375\n", "time= 0.6015625\n", "time= 0.609375\n", "time= 0.6171875\n", "time= 0.625\n", "time= 0.6328125\n", "time= 0.640625\n", "time= 0.6484375\n", "time= 0.65625\n", "time= 0.6640625\n", "time= 0.671875\n", "time= 0.6796875\n", "time= 0.6875\n", "time= 0.6953125\n", "time= 0.703125\n", "time= 0.7109375\n", "time= 0.71875\n", "time= 0.7265625\n", "time= 0.734375\n", "time= 0.7421875\n", "time= 0.75\n", "time= 0.7578125\n", "time= 0.765625\n", "time= 0.7734375\n", "time= 0.78125\n", "time= 0.7890625\n", "time= 0.796875\n", "time= 0.8046875\n", "time= 0.8125\n", "time= 0.8203125\n", "time= 0.828125\n", "time= 0.8359375\n", "time= 0.84375\n", "time= 0.8515625\n", "time= 0.859375\n", "time= 0.8671875\n", "time= 0.875\n", "time= 0.8828125\n", "time= 0.890625\n", "time= 0.8984375\n", "time= 0.90625\n", "time= 0.9140625\n", "time= 0.921875\n", "time= 0.9296875\n", "time= 0.9375\n", "time= 0.9453125\n", "time= 0.953125\n", "time= 0.9609375\n", "time= 0.96875\n", "time= 0.9765625\n", "time= 0.984375\n", "time= 0.9921875\n", "time= 1.0\n" ] } ], "source": [ "ani.save(schemename+'_TVD.gif', writer='imagemagick',fps=24,dpi=200)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/html": [ "CSS style adapted from https://github.com/barbagroup/CFDPython. Copyright (c) Barba group\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "text/plain": [ "" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import urllib\n", "import requests\n", "from IPython.core.display import HTML\n", "def css_styling():\n", " styles = requests.get(\"https://raw.githubusercontent.com/saadtony/NumericalMethods/master/styles/custom.css\")\n", " return HTML(styles.text)\n", "css_styling()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "anaconda-cloud": {}, "kernelspec": { "display_name": "Python [default]", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.8" } }, "nbformat": 4, "nbformat_minor": 2 }