{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Advection using the k-Scheme\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": {}, "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", "For constant grid spacing, the k-Scheme is given by\n", "\\begin{equation}\n", "{\\phi _f} = {\\phi _{\\rm{C}}} + \\frac{{1 - k}}{4}({\\phi _{\\rm{C}}} - {\\phi _{\\rm{U}}}) + \\frac{{1 + k}}{4}({\\phi _{\\rm{D}}} - {\\phi _{\\rm{C}}})\n", "\\end{equation}\n", "which, for a positive advecting velocity, gives us\n", "\\begin{equation}\n", "F_{i + {\\textstyle{1 \\over 2}}}^n = c\\phi _{i + {\\textstyle{1 \\over 2}}}^n = c{\\phi _i} + c\\frac{{1 - k}}{4}({\\phi _i} - {\\phi _{i - 1}}) + c\\frac{{1 + k}}{4}({\\phi _{i + 1}} - {\\phi _i})\n", "\\end{equation}" ] }, { "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": 14, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 14, "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 = 128 # 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": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "dt= 0.00390625\n", "dx= 0.0078125\n" ] } ], "source": [ "cfl =0.5\n", "c = 1.0 # use a negative value for left traveling waves\n", "dt = cfl*dx/abs(c)\n", "print('dt=',dt)\n", "print('dx=',dx)" ] }, { "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" ], "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# the k scheme\n", "k = 0.5\n", "# finite volume implementation with arrays for fluxes\n", "t = 0\n", "tend= L/abs(c)\n", "\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", "i = 0\n", "while t < tend: \n", " phin = sol[-1]\n", "\n", "# if (i%16==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", " # 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", " 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", " f[2:-2] = ϕc + (1-k)/4.0*(ϕc - ϕu) + (1+k)/4.0*(ϕd - ϕc)\n", " f = c*f # multiply the flux by the velocity\n", " # advect\n", " phi[2:-2] = phin[2:-2] - c * 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('k = '+ str(k), 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": 24, "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", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \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": "markdown", "metadata": {}, "source": [ "# Create Animation in Moving Reference Frame" ] }, { "cell_type": "code", "execution_count": 27, "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", " \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)\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", "ann = ax.annotate('time ='+str(round(t,3))+' s', xy=(2, 1), xytext=(40, 200),xycoords='figure points')\n", "ax.annotate('k ='+str(k) + ' (k-scheme)', xy=(2, 1), xytext=(40, 190),xycoords='figure points')\n", "plt.tight_layout()\n", "\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(1.1*min(f) - 0.1,1.1*max(f))\n", " ann.set_text('time ='+str(round(t,4))+'s (' + str(i)+ ').')\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/cfl)), init_func=init,\n", " interval=20, blit=False)\n", "print('done!')" ] }, { "cell_type": "code", "execution_count": 29, "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" ] }, { "ename": "KeyboardInterrupt", "evalue": "", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mKeyboardInterrupt\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mani\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0msave\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'__k-scheme_'\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0mstr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mk\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0;34m'.mp4'\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mfps\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m24\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mdpi\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m300\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/animation.py\u001b[0m in \u001b[0;36msave\u001b[0;34m(self, filename, writer, fps, dpi, codec, bitrate, extra_args, metadata, extra_anim, savefig_kwargs)\u001b[0m\n\u001b[1;32m 1172\u001b[0m \u001b[0;31m# TODO: See if turning off blit is really necessary\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1173\u001b[0m \u001b[0manim\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_draw_next_frame\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0md\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mblit\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mFalse\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1174\u001b[0;31m \u001b[0mwriter\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrab_frame\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m**\u001b[0m\u001b[0msavefig_kwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 1175\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1176\u001b[0m \u001b[0;31m# Reconnect signal for first draw if necessary\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/animation.py\u001b[0m in \u001b[0;36mgrab_frame\u001b[0;34m(self, **savefig_kwargs)\u001b[0m\n\u001b[1;32m 374\u001b[0m \u001b[0;31m# frame format and dpi.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 375\u001b[0m self.fig.savefig(self._frame_sink(), format=self.frame_format,\n\u001b[0;32m--> 376\u001b[0;31m dpi=self.dpi, **savefig_kwargs)\n\u001b[0m\u001b[1;32m 377\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mRuntimeError\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mIOError\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mas\u001b[0m \u001b[0me\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 378\u001b[0m \u001b[0mout\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0merr\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_proc\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcommunicate\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/figure.py\u001b[0m in \u001b[0;36msavefig\u001b[0;34m(self, fname, frameon, transparent, **kwargs)\u001b[0m\n\u001b[1;32m 2092\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mset_frameon\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mframeon\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2093\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2094\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcanvas\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mprint_figure\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfname\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 2095\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2096\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mframeon\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/backend_bases.py\u001b[0m in \u001b[0;36mprint_figure\u001b[0;34m(self, filename, dpi, facecolor, edgecolor, orientation, format, bbox_inches, **kwargs)\u001b[0m\n\u001b[1;32m 2073\u001b[0m \u001b[0morientation\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0morientation\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2074\u001b[0m \u001b[0mbbox_inches_restore\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0m_bbox_inches_restore\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2075\u001b[0;31m **kwargs)\n\u001b[0m\u001b[1;32m 2076\u001b[0m \u001b[0;32mfinally\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2077\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mbbox_inches\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0mrestore_bbox\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/backends/backend_agg.py\u001b[0m in \u001b[0;36mprint_raw\u001b[0;34m(self, filename_or_obj, *args, **kwargs)\u001b[0m\n\u001b[1;32m 459\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 460\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mprint_raw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfilename_or_obj\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 461\u001b[0;31m \u001b[0mFigureCanvasAgg\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 462\u001b[0m \u001b[0mrenderer\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_renderer\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 463\u001b[0m \u001b[0;32mwith\u001b[0m \u001b[0mcbook\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_setattr_cm\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrenderer\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdpi\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfigure\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;31m \u001b[0m\u001b[0;31m\\\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/backends/backend_agg.py\u001b[0m in \u001b[0;36mdraw\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 400\u001b[0m \u001b[0mtoolbar\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtoolbar\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 401\u001b[0m \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 402\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfigure\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mrenderer\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 403\u001b[0m \u001b[0;31m# A GUI class may be need to update a window using this draw, so\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 404\u001b[0m \u001b[0;31m# don't forget to call the superclass.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/artist.py\u001b[0m in \u001b[0;36mdraw_wrapper\u001b[0;34m(artist, renderer, *args, **kwargs)\u001b[0m\n\u001b[1;32m 48\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstart_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 49\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 50\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0martist\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 51\u001b[0m \u001b[0;32mfinally\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 52\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0martist\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_agg_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/figure.py\u001b[0m in \u001b[0;36mdraw\u001b[0;34m(self, renderer)\u001b[0m\n\u001b[1;32m 1647\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1648\u001b[0m mimage._draw_list_compositing_images(\n\u001b[0;32m-> 1649\u001b[0;31m renderer, self, artists, self.suppressComposite)\n\u001b[0m\u001b[1;32m 1650\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1651\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mclose_group\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'figure'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/image.py\u001b[0m in \u001b[0;36m_draw_list_compositing_images\u001b[0;34m(renderer, parent, artists, suppress_composite)\u001b[0m\n\u001b[1;32m 136\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mnot_composite\u001b[0m \u001b[0;32mor\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0mhas_images\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 137\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0ma\u001b[0m \u001b[0;32min\u001b[0m \u001b[0martists\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 138\u001b[0;31m \u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrenderer\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 139\u001b[0m \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 140\u001b[0m \u001b[0;31m# Composite any adjacent images together\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/artist.py\u001b[0m in \u001b[0;36mdraw_wrapper\u001b[0;34m(artist, renderer, *args, **kwargs)\u001b[0m\n\u001b[1;32m 48\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstart_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 49\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 50\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0martist\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 51\u001b[0m \u001b[0;32mfinally\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 52\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0martist\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_agg_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/axes/_base.py\u001b[0m in \u001b[0;36mdraw\u001b[0;34m(self, renderer, inframe)\u001b[0m\n\u001b[1;32m 2626\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstop_rasterizing\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2627\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2628\u001b[0;31m \u001b[0mmimage\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_draw_list_compositing_images\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrenderer\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0martists\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 2629\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2630\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mclose_group\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'axes'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/image.py\u001b[0m in \u001b[0;36m_draw_list_compositing_images\u001b[0;34m(renderer, parent, artists, suppress_composite)\u001b[0m\n\u001b[1;32m 136\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mnot_composite\u001b[0m \u001b[0;32mor\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0mhas_images\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 137\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0ma\u001b[0m \u001b[0;32min\u001b[0m \u001b[0martists\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 138\u001b[0;31m \u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrenderer\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 139\u001b[0m \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 140\u001b[0m \u001b[0;31m# Composite any adjacent images together\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/artist.py\u001b[0m in \u001b[0;36mdraw_wrapper\u001b[0;34m(artist, renderer, *args, **kwargs)\u001b[0m\n\u001b[1;32m 48\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstart_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 49\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 50\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0martist\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 51\u001b[0m \u001b[0;32mfinally\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 52\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0martist\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_agg_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/text.py\u001b[0m in \u001b[0;36mdraw\u001b[0;34m(self, renderer)\u001b[0m\n\u001b[1;32m 2391\u001b[0m \u001b[0;31m# Draw text, including FancyBboxPatch, after FancyArrowPatch.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2392\u001b[0m \u001b[0;31m# Otherwise, a wedge arrowstyle can land partly on top of the Bbox.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2393\u001b[0;31m \u001b[0mText\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 2394\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2395\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mget_window_extent\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/artist.py\u001b[0m in \u001b[0;36mdraw_wrapper\u001b[0;34m(artist, renderer, *args, **kwargs)\u001b[0m\n\u001b[1;32m 48\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstart_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 49\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 50\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mdraw\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0martist\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrenderer\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 51\u001b[0m \u001b[0;32mfinally\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 52\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0martist\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_agg_filter\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/text.py\u001b[0m in \u001b[0;36mdraw\u001b[0;34m(self, renderer)\u001b[0m\n\u001b[1;32m 752\u001b[0m textrenderer.draw_tex(gc, x, y, clean_line,\n\u001b[1;32m 753\u001b[0m \u001b[0mtextobj\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_fontproperties\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mangle\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 754\u001b[0;31m mtext=mtext)\n\u001b[0m\u001b[1;32m 755\u001b[0m \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 756\u001b[0m textrenderer.draw_text(gc, x, y, clean_line,\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/backends/backend_agg.py\u001b[0m in \u001b[0;36mdraw_tex\u001b[0;34m(self, gc, x, y, s, prop, angle, ismath, mtext)\u001b[0m\n\u001b[1;32m 231\u001b[0m \u001b[0mtexmanager\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_texmanager\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 232\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 233\u001b[0;31m \u001b[0mZ\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mtexmanager\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_grey\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msize\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 234\u001b[0m \u001b[0mZ\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0marray\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mZ\u001b[0m \u001b[0;34m*\u001b[0m \u001b[0;36m255.0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0muint8\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 235\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/texmanager.py\u001b[0m in \u001b[0;36mget_grey\u001b[0;34m(self, tex, fontsize, dpi)\u001b[0m\n\u001b[1;32m 418\u001b[0m \u001b[0malpha\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrey_arrayd\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 419\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0malpha\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 420\u001b[0;31m \u001b[0mpngfile\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmake_png\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mtex\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfontsize\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 421\u001b[0m \u001b[0mX\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0m_png\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mread_png\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mos\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpath\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mjoin\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtexcache\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpngfile\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 422\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgrey_arrayd\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0malpha\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mX\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m:\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m-\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/texmanager.py\u001b[0m in \u001b[0;36mmake_png\u001b[0;34m(self, tex, fontsize, dpi)\u001b[0m\n\u001b[1;32m 383\u001b[0m self._run_checked_subprocess(\n\u001b[1;32m 384\u001b[0m [\"dvipng\", \"-bg\", \"Transparent\", \"-D\", str(dpi),\n\u001b[0;32m--> 385\u001b[0;31m \"-T\", \"tight\", \"-o\", pngfile, dvifile], tex)\n\u001b[0m\u001b[1;32m 386\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mpngfile\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 387\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/site-packages/matplotlib/texmanager.py\u001b[0m in \u001b[0;36m_run_checked_subprocess\u001b[0;34m(self, command, tex)\u001b[0m\n\u001b[1;32m 296\u001b[0m report = subprocess.check_output(command,\n\u001b[1;32m 297\u001b[0m \u001b[0mcwd\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtexcache\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 298\u001b[0;31m stderr=subprocess.STDOUT)\n\u001b[0m\u001b[1;32m 299\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0msubprocess\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mCalledProcessError\u001b[0m \u001b[0;32mas\u001b[0m \u001b[0mexc\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 300\u001b[0m raise RuntimeError(\n", "\u001b[0;32m/anaconda3/lib/python3.6/subprocess.py\u001b[0m in \u001b[0;36mcheck_output\u001b[0;34m(timeout, *popenargs, **kwargs)\u001b[0m\n\u001b[1;32m 354\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 355\u001b[0m return run(*popenargs, stdout=PIPE, timeout=timeout, check=True,\n\u001b[0;32m--> 356\u001b[0;31m **kwargs).stdout\n\u001b[0m\u001b[1;32m 357\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 358\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/subprocess.py\u001b[0m in \u001b[0;36mrun\u001b[0;34m(input, timeout, check, *popenargs, **kwargs)\u001b[0m\n\u001b[1;32m 423\u001b[0m \u001b[0;32mwith\u001b[0m \u001b[0mPopen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mpopenargs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mas\u001b[0m \u001b[0mprocess\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 424\u001b[0m \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 425\u001b[0;31m \u001b[0mstdout\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mstderr\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mprocess\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcommunicate\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minput\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtimeout\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mtimeout\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 426\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0mTimeoutExpired\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 427\u001b[0m \u001b[0mprocess\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mkill\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m/anaconda3/lib/python3.6/subprocess.py\u001b[0m in \u001b[0;36mcommunicate\u001b[0;34m(self, input, timeout)\u001b[0m\n\u001b[1;32m 848\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_stdin_write\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minput\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 849\u001b[0m \u001b[0;32melif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstdout\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 850\u001b[0;31m \u001b[0mstdout\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstdout\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mread\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 851\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstdout\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mclose\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 852\u001b[0m \u001b[0;32melif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstderr\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;31mKeyboardInterrupt\u001b[0m: " ] } ], "source": [ "ani.save('__k-scheme_' + str(k)+'.mp4',fps=24,dpi=300)" ] }, { "cell_type": "code", "execution_count": 26, "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": 26, "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 3", "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 }