euler_1d_jax.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "23b3ca67",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import jax.numpy as jnp\n",
    "from jax import jit\n",
    "import jax\n",
    "from jax.lib import xla_bridge"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "23579a70",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "WARNING:absl:No GPU/TPU found, falling back to CPU. (Set TF_CPP_MIN_LOG_LEVEL=0 and rerun for more info.)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "cpu\n"
     ]
    }
   ],
   "source": [
    "print(xla_bridge.get_backend().platform)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6b4d85ee",
   "metadata": {},
   "source": [
    "### Grid "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "38ee8260",
   "metadata": {},
   "outputs": [],
   "source": [
    "class Grid:\n",
    "    def __init__(self, xmin: float, xmax: float, nx: int, n_ghost: int = 1):\n",
    "        self.xmin = xmin\n",
    "        self.xmax = xmax\n",
    "        self.nx = nx\n",
    "        self.tot_nx = nx + 2 * n_ghost\n",
    "        self.dx = (xmax - xmin) / (nx - 1)\n",
    "        self.x = jnp.linspace(xmin, xmax, nx)\n",
    "        self.x = jnp.insert(self.x, 0, self.x[0])\n",
    "        self.x = jnp.insert(self.x, self.x.shape[0], self.x[-1])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a8b7a1d0",
   "metadata": {},
   "source": [
    "### Primitive variables to conserved variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "aa482c99",
   "metadata": {},
   "outputs": [],
   "source": [
    "@jit\n",
    "def prim_to_cons_vars(pv: jnp.ndarray) -> jnp.ndarray:\n",
    "    cv = jnp.zeros(pv.shape)\n",
    "    cv = cv.at[:,0].set(pv[:, 0])\n",
    "    cv = cv.at[:,1].set(pv[:, 0] * pv[:, 1])\n",
    "    cv = cv.at[:,2].set(pv[:, 0] * (pv[:, 2] / (pv[:, 0] * (1.4 - 1)) + 0.5 * pv[:, 1] ** 2))\n",
    "\n",
    "    return cv"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a2b6f0a0",
   "metadata": {},
   "source": [
    "### Conserved variables to primitive variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "4ee818b4",
   "metadata": {},
   "outputs": [],
   "source": [
    "@jit\n",
    "def cons_to_prim_vars(cv: jnp.ndarray) -> jnp.ndarray:\n",
    "    pv = jnp.zeros(cv.shape)\n",
    "    pv = pv.at[:, 0].set(cv[:, 0])\n",
    "    pv = pv.at[:, 1].set(cv[:, 1] / cv[:, 0])\n",
    "    pv = pv.at[:, 2].set((cv[:, 2] / cv[:, 0] - 0.5 * pv[:, 1] ** 2) * cv[:, 0] * (1.4 - 1))\n",
    "\n",
    "    return pv"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8004a4fa",
   "metadata": {},
   "source": [
    "### Conserved variables to flux"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "568be14c",
   "metadata": {},
   "outputs": [],
   "source": [
    "@jit\n",
    "def conv_vars_to_flux(cv: jnp.ndarray) -> jnp.ndarray:\n",
    "    pv = cons_to_prim_vars(cv)\n",
    "    f = jnp.zeros(cv.shape)\n",
    "    f = f.at[:, 0].set(cv[:, 1])\n",
    "    f = f.at[:, 1].set(pv[:, 2] + cv[:, 1] * pv[:, 1])\n",
    "    f = f.at[:, 2].set(pv[:, 2] * pv[:, 1] + cv[:, 2] * pv[:, 1])\n",
    "\n",
    "    return f"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "51c2177f",
   "metadata": {},
   "source": [
    "### To compute time step"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "46c7b070",
   "metadata": {},
   "outputs": [],
   "source": [
    "@jit\n",
    "def time_step(cv: jnp.ndarray, dx: float, cfl: float = 0.2) -> jnp.ndarray:\n",
    "    pv = cons_to_prim_vars(cv)\n",
    "    a = jnp.sqrt(1.4 * pv[:, 2] / pv[:, 0])\n",
    "    eig_val = jnp.abs(pv[:, 1] + a)\n",
    "    dt = cfl * dx / eig_val\n",
    "    \n",
    "    return jnp.min(dt)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bf5c8435",
   "metadata": {},
   "source": [
    "### To compute Rusanov (LLF) flux"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "857e55b6",
   "metadata": {},
   "outputs": [],
   "source": [
    "@jit\n",
    "def llf_flux(cv: jnp.ndarray) -> jnp.ndarray:\n",
    "    f = conv_vars_to_flux(cv)\n",
    "    pv = cons_to_prim_vars(cv)\n",
    "    a = jnp.sqrt(1.4 * pv[:, 2] / pv[:, 0])\n",
    "    eig_val_l = jnp.abs(pv[1:-1, 1]) + a[1:-1]\n",
    "    eig_val_r = jnp.abs(pv[2:, 1]) + a[2:]\n",
    "    max_eig_val = jnp.maximum(eig_val_l, eig_val_r)\n",
    "    max_eig_val = max_eig_val.reshape(max_eig_val.size, 1)\n",
    "    f_iph = 0.5 * (f[1:-1, :] + f[2:, :]) - 0.5 * max_eig_val * (cv[2:, :] - cv[1:-1, :])\n",
    "\n",
    "    eig_val_l = jnp.abs(pv[:-2, 1]) + a[:-2]\n",
    "    eig_val_r = jnp.abs(pv[1:-1, 1]) + a[1:-1]\n",
    "    max_eig_val = jnp.maximum(eig_val_l, eig_val_r)\n",
    "    max_eig_val = max_eig_val.reshape(max_eig_val.size, 1)\n",
    "    f_imh = 0.5 * (f[1:-1, :] + f[:-2, :]) - 0.5 * max_eig_val * (cv[1:-1, :] - cv[:-2, :])\n",
    "\n",
    "    return f_iph - f_imh"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "06694378",
   "metadata": {},
   "source": [
    "### Initial conditions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "0a1dee46",
   "metadata": {},
   "outputs": [],
   "source": [
    "def initialize(pv_l: jnp.ndarray, pv_r: jnp.ndarray, x: jnp.ndarray, pv: jnp.ndarray, xloc):\n",
    "    pv=pv.at[jnp.argwhere(x <= xloc), :].set(pv_l)\n",
    "    pv=pv.at[jnp.argwhere(x > xloc), :].set(pv_r)\n",
    "    \n",
    "    return pv"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f2a49f3b",
   "metadata": {},
   "source": [
    "### Boundary conditions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "ae0f9610",
   "metadata": {},
   "outputs": [],
   "source": [
    "@jit\n",
    "def bc(cv: jnp.ndarray) -> jnp.ndarray:\n",
    "    pv = cons_to_prim_vars(cv)\n",
    "    pv = jnp.insert(pv, 0, pv[0], axis=0)\n",
    "    pv = jnp.insert(pv, pv.shape[0], pv[-1], axis=0)\n",
    "    con_var = prim_to_cons_vars(pv)\n",
    "\n",
    "    return con_var"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ee1213b9",
   "metadata": {},
   "source": [
    "### Simulation parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "0b89b92e",
   "metadata": {},
   "outputs": [],
   "source": [
    "prim_var_l = jnp.array([1.0, 0.75, 1.0])\n",
    "prim_var_r = jnp.array([0.125, 0.0, 0.1])\n",
    "\n",
    "cfl_num = 0.9\n",
    "grid = Grid(0, 1, 1001)\n",
    "x_grid = grid.x\n",
    "x_loc = 0.3\n",
    "prim_var = jnp.zeros((x_grid.size, 3))\n",
    "prim_var = initialize(prim_var_l, prim_var_r, x_grid, prim_var, x_loc)\n",
    "# print(prim_var)\n",
    "cons_var = prim_to_cons_vars(prim_var)\n",
    "# print(cons_var)\n",
    "flux = conv_vars_to_flux(cons_var)\n",
    "prim_var = cons_to_prim_vars(cons_var)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b97101da",
   "metadata": {},
   "source": [
    "### Run the simulation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "db27f072",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0004655455\n",
      "0.00088304735\n",
      "0.0012766281\n",
      "0.001655027\n",
      "0.0020256073\n",
      "0.0023936369\n",
      "0.0027578631\n",
      "0.0031157434\n",
      "0.0034709622\n",
      "0.0038232978\n",
      "0.0041723833\n",
      "0.0045201397\n",
      "0.004867265\n",
      "0.005211972\n",
      "0.005555554\n",
      "0.005898653\n",
      "0.006239865\n",
      "0.006580139\n",
      "0.00692052\n",
      "0.0072592623\n",
      "0.0075971293\n",
      "0.007935035\n",
      "0.008271811\n",
      "0.008607763\n",
      "0.008943652\n",
      "0.009278902\n",
      "0.009613367\n",
      "0.009947683\n",
      "0.010281693\n",
      "0.01061497\n",
      "0.010948034\n",
      "0.011281036\n",
      "0.011613363\n",
      "0.011945435\n",
      "0.012277599\n",
      "0.012609157\n",
      "0.012940436\n",
      "0.013271902\n",
      "0.013602829\n",
      "0.013933472\n",
      "0.014264209\n",
      "0.014594615\n",
      "0.014924748\n",
      "0.015254904\n",
      "0.015584874\n",
      "0.015914591\n",
      "0.016244283\n",
      "0.016573887\n",
      "0.016903263\n",
      "0.017232586\n",
      "0.01756188\n",
      "0.017890973\n",
      "0.018219998\n",
      "0.01854903\n",
      "0.018877884\n",
      "0.019206667\n",
      "0.019535476\n",
      "0.01986413\n",
      "0.020192716\n",
      "0.020521339\n",
      "0.020849824\n",
      "0.021178244\n",
      "0.021506699\n",
      "0.02183504\n",
      "0.022163322\n",
      "0.022491625\n",
      "0.022819845\n",
      "0.023148011\n",
      "0.023476189\n",
      "0.023804309\n",
      "0.024132378\n",
      "0.024460452\n",
      "0.02478849\n",
      "0.025116477\n",
      "0.025444465\n",
      "0.025772434\n",
      "0.026100352\n",
      "0.026428267\n",
      "0.026756182\n",
      "0.027084045\n",
      "0.0274119\n",
      "0.027739769\n",
      "0.028067587\n",
      "0.028395392\n",
      "0.028723225\n",
      "0.02905101\n",
      "0.029378774\n",
      "0.029706558\n",
      "0.030034315\n",
      "0.03036205\n",
      "0.030689795\n",
      "0.03101753\n",
      "0.03134524\n",
      "0.031672955\n",
      "0.032000672\n",
      "0.032328364\n",
      "0.03265605\n",
      "0.032983758\n",
      "0.033311438\n",
      "0.033639107\n",
      "0.033966795\n",
      "0.034294464\n",
      "0.03462212\n",
      "0.034949787\n",
      "0.035277452\n",
      "0.035605103\n",
      "0.035932753\n",
      "0.03626041\n",
      "0.036588054\n",
      "0.036915693\n",
      "0.037243348\n",
      "0.03757099\n",
      "0.037898622\n",
      "0.03822627\n",
      "0.03855391\n",
      "0.03888154\n",
      "0.039209176\n",
      "0.039536815\n",
      "0.039864447\n",
      "0.040192075\n",
      "0.04051971\n",
      "0.040847342\n",
      "0.041174967\n",
      "0.041502602\n",
      "0.041830234\n",
      "0.04215786\n",
      "0.042485487\n",
      "0.04281312\n",
      "0.043140747\n",
      "0.04346837\n",
      "0.043796003\n",
      "0.04412363\n",
      "0.044451255\n",
      "0.044778887\n",
      "0.045106515\n",
      "0.045434143\n",
      "0.04576177\n",
      "0.0460894\n",
      "0.046417024\n",
      "0.046744652\n",
      "0.047072284\n",
      "0.04739991\n",
      "0.047727536\n",
      "0.04805517\n",
      "0.048382796\n",
      "0.048710424\n",
      "0.049038056\n",
      "0.049365684\n",
      "0.049693312\n",
      "0.05002094\n",
      "0.050348572\n",
      "0.0506762\n",
      "0.05100383\n",
      "0.051331464\n",
      "0.051659092\n",
      "0.05198672\n",
      "0.052314352\n",
      "0.052641984\n",
      "0.052969612\n",
      "0.053297244\n",
      "0.05362488\n",
      "0.05395251\n",
      "0.054280143\n",
      "0.05460778\n",
      "0.05493541\n",
      "0.055263042\n",
      "0.055590678\n",
      "0.055918314\n",
      "0.056245945\n",
      "0.05657358\n",
      "0.05690122\n",
      "0.057228856\n",
      "0.05755649\n",
      "0.05788413\n",
      "0.05821177\n",
      "0.058539405\n",
      "0.058867045\n",
      "0.059194684\n",
      "0.059522323\n",
      "0.059849963\n",
      "0.060177606\n",
      "0.060505245\n",
      "0.060832884\n",
      "0.061160527\n",
      "0.06148817\n",
      "0.06181581\n",
      "0.062143452\n",
      "0.0624711\n",
      "0.06279874\n",
      "0.06312638\n",
      "0.063454024\n",
      "0.06378167\n",
      "0.06410932\n",
      "0.064436965\n",
      "0.06476461\n",
      "0.06509226\n",
      "0.065419905\n",
      "0.06574755\n",
      "0.0660752\n",
      "0.066402845\n",
      "0.06673049\n",
      "0.06705814\n",
      "0.067385785\n",
      "0.06771343\n",
      "0.068041086\n",
      "0.06836873\n",
      "0.06869638\n",
      "0.069024034\n",
      "0.06935168\n",
      "0.06967933\n",
      "0.07000698\n",
      "0.070334636\n",
      "0.07066228\n",
      "0.07098994\n",
      "0.07131759\n",
      "0.071645245\n",
      "0.0719729\n",
      "0.07230055\n",
      "0.07262821\n",
      "0.07295586\n",
      "0.073283516\n",
      "0.07361117\n",
      "0.073938824\n",
      "0.07426648\n",
      "0.07459413\n",
      "0.07492179\n",
      "0.07524944\n",
      "0.075577095\n",
      "0.07590475\n",
      "0.0762324\n",
      "0.07656006\n",
      "0.07688771\n",
      "0.077215366\n",
      "0.07754302\n",
      "0.077870674\n",
      "0.07819833\n",
      "0.07852598\n",
      "0.07885364\n",
      "0.0791813\n",
      "0.07950895\n",
      "0.07983661\n",
      "0.08016427\n",
      "0.08049192\n",
      "0.08081958\n",
      "0.08114724\n",
      "0.0814749\n",
      "0.081802554\n",
      "0.082130216\n",
      "0.08245788\n",
      "0.08278553\n",
      "0.083113186\n",
      "0.08344085\n",
      "0.08376851\n",
      "0.08409616\n",
      "0.084423825\n",
      "0.08475149\n",
      "0.08507915\n",
      "0.08540681\n",
      "0.08573447\n",
      "0.08606213\n",
      "0.086389795\n",
      "0.08671746\n",
      "0.08704512\n",
      "0.08737278\n",
      "0.08770044\n",
      "0.0880281\n",
      "0.088355765\n",
      "0.08868343\n",
      "0.08901109\n",
      "0.08933875\n",
      "0.08966641\n",
      "0.08999407\n",
      "0.090321735\n",
      "0.090649396\n",
      "0.09097706\n",
      "0.09130472\n",
      "0.09163238\n",
      "0.09196004\n",
      "0.092287704\n",
      "0.092615366\n",
      "0.09294303\n",
      "0.09327069\n",
      "0.09359835\n",
      "0.09392601\n",
      "0.094253674\n",
      "0.094581336\n",
      "0.094909\n",
      "0.09523666\n",
      "0.09556432\n",
      "0.09589198\n",
      "0.096219644\n",
      "0.096547306\n",
      "0.09687497\n",
      "0.09720263\n",
      "0.09753029\n",
      "0.09785795\n",
      "0.098185614\n",
      "0.098513275\n",
      "0.09884094\n",
      "0.0991686\n",
      "0.09949626\n",
      "0.09982392\n",
      "0.10015158\n",
      "0.100479245\n",
      "0.10080691\n",
      "0.10113457\n",
      "0.10146223\n",
      "0.10178989\n",
      "0.10211755\n",
      "0.102445215\n",
      "0.10277288\n",
      "0.10310054\n",
      "0.1034282\n",
      "0.10375586\n",
      "0.10408352\n",
      "0.104411185\n",
      "0.10473885\n",
      "0.10506651\n",
      "0.10539417\n",
      "0.10572183\n",
      "0.10604949\n",
      "0.106377155\n",
      "0.106704816\n",
      "0.10703248\n",
      "0.10736014\n",
      "0.1076878\n",
      "0.10801546\n",
      "0.108343124\n",
      "0.108670786\n",
      "0.10899845\n",
      "0.10932611\n",
      "0.10965377\n",
      "0.10998143\n",
      "0.110309094\n",
      "0.110636756\n",
      "0.11096442\n",
      "0.11129208\n",
      "0.11161974\n",
      "0.1119474\n",
      "0.112275064\n",
      "0.112602726\n",
      "0.11293039\n",
      "0.11325805\n",
      "0.11358571\n",
      "0.11391337\n",
      "0.114241034\n",
      "0.114568695\n",
      "0.11489636\n",
      "0.11522402\n",
      "0.11555168\n",
      "0.11587934\n",
      "0.116207\n",
      "0.116534665\n",
      "0.11686233\n",
      "0.11718999\n",
      "0.11751765\n",
      "0.11784531\n",
      "0.11817297\n",
      "0.118500635\n",
      "0.1188283\n",
      "0.11915596\n",
      "0.11948362\n",
      "0.11981128\n",
      "0.12013894\n",
      "0.120466605\n",
      "0.12079427\n",
      "0.12112193\n",
      "0.12144959\n",
      "0.12177725\n",
      "0.12210491\n",
      "0.122432575\n",
      "0.122760236\n",
      "0.1230879\n",
      "0.12341556\n",
      "0.12374322\n",
      "0.124070875\n",
      "0.12439853\n",
      "0.12472618\n",
      "0.12505384\n",
      "0.1253815\n",
      "0.12570916\n",
      "0.12603682\n",
      "0.12636448\n",
      "0.12669215\n",
      "0.12701981\n",
      "0.12734747\n",
      "0.12767513\n",
      "0.1280028\n",
      "0.12833045\n",
      "0.12865812\n",
      "0.12898578\n",
      "0.12931344\n",
      "0.1296411\n",
      "0.12996876\n",
      "0.13029642\n",
      "0.13062409\n",
      "0.13095175\n",
      "0.13127941\n",
      "0.13160707\n",
      "0.13193473\n",
      "0.1322624\n",
      "0.13259006\n",
      "0.13291772\n",
      "0.13324538\n",
      "0.13357304\n",
      "0.1339007\n",
      "0.13422836\n",
      "0.13455603\n",
      "0.13488369\n",
      "0.13521135\n",
      "0.13553901\n",
      "0.13586667\n",
      "0.13619433\n",
      "0.136522\n",
      "0.13684966\n",
      "0.13717732\n",
      "0.13750498\n",
      "0.13783264\n",
      "0.1381603\n",
      "0.13848796\n",
      "0.13881563\n",
      "0.13914329\n",
      "0.13947095\n",
      "0.13979861\n",
      "0.14012627\n",
      "0.14045393\n",
      "0.1407816\n",
      "0.14110924\n",
      "0.14143689\n",
      "0.14176454\n",
      "0.14209218\n",
      "0.14241983\n",
      "0.14274748\n",
      "0.14307512\n",
      "0.14340277\n",
      "0.14373042\n",
      "0.14405806\n",
      "0.14438571\n",
      "0.14471336\n",
      "0.145041\n",
      "0.14536865\n",
      "0.1456963\n",
      "0.14602394\n",
      "0.14635159\n",
      "0.14667924\n",
      "0.14700688\n",
      "0.14733453\n",
      "0.14766218\n",
      "0.14798982\n",
      "0.14831747\n",
      "0.14864512\n",
      "0.14897276\n",
      "0.14930041\n",
      "0.14962806\n",
      "0.1499557\n",
      "0.15028335\n",
      "0.150611\n",
      "0.15093865\n",
      "0.15126629\n",
      "0.15159394\n",
      "0.15192159\n",
      "0.15224923\n",
      "0.15257688\n",
      "0.15290453\n",
      "0.15323217\n",
      "0.15355982\n",
      "0.15388747\n",
      "0.15421511\n",
      "0.15454276\n",
      "0.1548704\n",
      "0.15519805\n",
      "0.1555257\n",
      "0.15585335\n",
      "0.156181\n",
      "0.15650864\n",
      "0.15683629\n",
      "0.15716393\n",
      "0.15749158\n",
      "0.15781923\n",
      "0.15814687\n",
      "0.15847452\n",
      "0.15880217\n",
      "0.15912981\n",
      "0.15945746\n",
      "0.1597851\n",
      "0.16011275\n",
      "0.1604404\n",
      "0.16076805\n",
      "0.1610957\n",
      "0.16142334\n",
      "0.16175099\n",
      "0.16207863\n",
      "0.16240628\n",
      "0.16273393\n",
      "0.16306157\n",
      "0.16338922\n",
      "0.16371687\n",
      "0.16404451\n",
      "0.16437216\n",
      "0.16469981\n",
      "0.16502745\n",
      "0.1653551\n",
      "0.16568275\n",
      "0.1660104\n",
      "0.16633804\n",
      "0.16666569\n",
      "0.16699333\n",
      "0.16732098\n",
      "0.16764863\n",
      "0.16797628\n",
      "0.16830392\n",
      "0.16863157\n",
      "0.16895922\n",
      "0.16928686\n",
      "0.16961451\n",
      "0.16994216\n",
      "0.1702698\n",
      "0.17059745\n",
      "0.1709251\n",
      "0.17125274\n",
      "0.17158039\n",
      "0.17190804\n",
      "0.17223568\n",
      "0.17256333\n",
      "0.17289098\n",
      "0.17321862\n",
      "0.17354627\n",
      "0.17387392\n",
      "0.17420156\n",
      "0.17452921\n",
      "0.17485686\n",
      "0.1751845\n",
      "0.17551215\n",
      "0.1758398\n",
      "0.17616744\n",
      "0.17649509\n",
      "0.17682274\n",
      "0.17715038\n",
      "0.17747803\n",
      "0.17780568\n",
      "0.17813332\n",
      "0.17846097\n",
      "0.17878862\n",
      "0.17911626\n",
      "0.17944391\n",
      "0.17977156\n",
      "0.1800992\n",
      "0.18042685\n",
      "0.1807545\n",
      "0.18108214\n",
      "0.18140979\n",
      "0.18173744\n",
      "0.18206508\n",
      "0.18239273\n",
      "0.18272038\n",
      "0.18304802\n",
      "0.18337567\n",
      "0.18370332\n",
      "0.18403096\n",
      "0.18435861\n",
      "0.18468626\n",
      "0.1850139\n",
      "0.18534155\n",
      "0.1856692\n",
      "0.18599685\n",
      "0.18632449\n",
      "0.18665214\n",
      "0.18697979\n",
      "0.18730743\n",
      "0.18763508\n",
      "0.18796273\n",
      "0.18829037\n",
      "0.18861802\n",
      "0.18894567\n",
      "0.18927331\n",
      "0.18960096\n",
      "0.1899286\n",
      "0.19025625\n",
      "0.1905839\n",
      "0.19091155\n",
      "0.1912392\n",
      "0.19156684\n",
      "0.19189449\n",
      "0.19222213\n",
      "0.19254978\n",
      "0.19287743\n",
      "0.19320507\n",
      "0.19353272\n",
      "0.19386037\n",
      "0.19418801\n",
      "0.19451566\n",
      "0.1948433\n",
      "0.19517095\n",
      "0.1954986\n",
      "0.19582625\n",
      "0.1961539\n",
      "0.19648154\n",
      "0.19680919\n",
      "0.19713683\n",
      "0.19746448\n",
      "0.19779213\n",
      "0.19811977\n",
      "0.19844742\n",
      "0.19877507\n",
      "0.19910271\n",
      "0.19943036\n",
      "0.19975801\n",
      "0.20008565\n",
      "CPU times: user 2.13 s, sys: 21.9 ms, total: 2.15 s\n",
      "Wall time: 2.09 s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "time = 0\n",
    "while time <= 0.2:\n",
    "    t_dt = time_step(cons_var, grid.dx, cfl_num)\n",
    "    flux = llf_flux(cons_var)\n",
    "    cv_new = cons_var[1:-1, :] - (t_dt / grid.dx) * flux\n",
    "    cons_var = bc(cv_new)\n",
    "    time += t_dt\n",
    "    print(time)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "705d52c5",
   "metadata": {},
   "source": [
    "### Plot results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "aa868f03",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([0. , 0.2, 0.4, 0.6, 0.8, 1. , 1.2]),\n",
       " [Text(0, 0, ''),\n",
       "  Text(0, 0, ''),\n",
       "  Text(0, 0, ''),\n",
       "  Text(0, 0, ''),\n",
       "  Text(0, 0, ''),\n",
       "  Text(0, 0, ''),\n",
       "  Text(0, 0, '')])"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZEAAAEQCAYAAABxzUkqAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/NK7nSAAAACXBIWXMAAAsTAAALEwEAmpwYAAAuYklEQVR4nO3dd3xV9f3H8dfn3iwCJKwgIAJSUAEV0agIDsSBOLBV1FpXbQFX1Z/a4aqCpa22zurPQbU/EBW1LkCtdQCCAzUOkCGiLMVBgDADIeP7++PchJBBkjtyTm7ez4fXm3zP+Z77+ZLkfu53nHPMOYeIiEg0Qn4HICIiTZeSiIiIRE1JREREoqYkIiIiUVMSERGRqKX4HUBj69Chg+vRo4ffYYiINCkff/zxWudcTtXyZpdEevToQV5ent9hiIg0KWa2sqZyDWeJiEjUlERERCRqSiIiIhI1X5KImXU1s/vN7H0zKzQzZ2Y96lk3ZGY3mNkKM9tuZvPM7MwEhywiIjXwqyfSCzgbKADmNLDun4CxwAPAcGAu8G8zOzmeAYqISN38Wp012zm3B4CZjQJOrE8lM+sI/Ba43Tl3Z6R4ppn1Am4HXk1EsCIiUjNfkohzrizKqsOANOCJKuVPAP8ys72dc8tjCk4kSW3cuJG1a9eyY8cOv0ORgAiHw7Ru3Zp27dqRnp4e1TGa2nki/YAi4Ksq5Qsjz32BhCSRt597kPQNVV8WsPInq6m40vaG1dl1S6V9rJY6QJVNlY63c0M4ZGSmhclMSyE7I4WM1HAtr1SDUAr0/zm07V7/OhII27dv58cff6Rr1660aNGi2u+RND/OOYqLi9m0aROrVq2iW7duUSWSppZE2gEbXPWboKyvtL0aMxsDjAHo1q1bVC8c+mIahxXPjapukJVVSXK7f2txsGwm/Oq1xAYlcZefn09OTg6ZmZl+hyIBYWakpaXRoUMHANavX0/nzp0bfJymlkSi4pybAEwAyM3NjeouXEfdXP2NszyXlac0V7V8l2277rvzGNW313ScXepGUad836LiMvK3FJG/uYiv87cw/5uNzFmaz9YdpfTMacnVx/VmRP8uNX9SfeceeHMs5H8JOftU3y6BtX37djp16uR3GBJQWVlZrFixolkkkQKgjZlZld5IeQ9kfQ11Eqb8jbb6+22whwr2arfrp9HtxaW8vuhHHpz5FVc//RkvfLKav408kD2yMnateNB5MGM8fPo4nDi+ESOWWJWUlJCS0tT+3KWxpKamUlpaGlXdpnay4UIgHfhJlfK+kedFjRtOcshIDTOifxdeveooxo3ox0cr1jPigXeY/+2GXXds1RH2HQ6fTYESTc42NZoHkdrE8rvR1JLIa0AxcF6V8vOBBVqZFZtQyLhoUA+ev2wQKaEQ5zwylw+Wrdt1p4MvgsK1sESrqUXExyRiZiPNbCRwSKRoeKTsmEr7lJjZY+XfO+fWAHcDN5jZtWY2xMweAoYCNzRm/MmsT+csXrxiEF3aZHDxxI/4eGXBzo0/GQpZXeGTx/0LUEQCw8+eyL8jj0sj3z8Y+X5cpX3CkUdlNwHjgauB/wKDgbOdcy8nNNpmpmPrDKaMHkhO63TGPJ7HtwWF3oZQGAacD1/PgIIarwwt0igmTpyImfHVVzUsvQfGjh2LmVFSUlLrMX75y19iZjU+dldPdvJtps05V+cgXE37OOdK8ZKIZnYTrGNWBo9ddCg/e/Bdfj0xj+cvH0Sr9BQYcB68fQd89iQce6PfYYrEJCcnh2nTplUr10KE+tG/kuxWr46tePC8g7noXx8ydtpC7jyrP7TpBr2Og0+fgGP+4PVORJqotLQ0Bg4c6HcYTVZTm1gXHxzVO4crju3Fcx9/y/R533mFB18Im1bDV2/5G5yI+Eo9EamXq47rzTtfreXGFz9nQLc2dN1nOLTMgU8mwT71un6mBMy46QtZ9N0mX2Po2yWLW0/r52sMQLX5j1AoRCikz9j1oX8lqZfUcIj7zhlAWZnjphcX4MKp0P9c+PI12Pyj3+GJRG316tWkpqbu8rjlllv8DqvJUE9E6q1b+0x+O2xfxk1fxLR533H6wRfCe/+AeU/Bkdf4HZ40UBB6AEHQsWNHXnnllV3KunTp4lM0TY+SiDTIhUf04KVPV3Pb9EUcc90xtOk+GD6eBIOuBnX/pQlKTU0lNzfX7zCaLP3VS4OEQ8ZfzziQDduK+curi+GQi6FguXd1XxFpdpREpMH6dsli1JF782zet8zPOhoyO0Dev/wOS0R8oOEsicqVx/Xm+U9WM/bVr3h+wAXYe/fBxtWQvaffoUkz89prr1W7zH12dnbF1y+88EK1lVadO3dm8ODBjRJfslMSkai0Sk/h98P25ffPz+eN/idzorsXPp4IQ2/yOzRpZq688spqZf369WPkyJEAnHPOOdW2n3LKKbz8sq6UFA9W/SaByS03N9fl5eX5HUZSKCtzjPjfd1i3ZQfv7PUw4R/mwzULIJzqd2hSxeLFi+nTp4/fYUiA1fU7YmYfO+eqrUDQnIhELRQybjm1H99v3M701OGw5Qf44pW6K4pI0lASkZgctnc7TjmwMzcu2IPSrL0g77G6K4lI0lASkZj9fti+7Cg13sg8GZbP9u7BLiLNgpKIxKx7+5acd3g3blk1ABdK1XJfkWZESUTi4jdDe7MlpS15mUfBZ0/Bjq1+hyQijUBJROIip3U6o47qyR3rjoSijbDgeb9DEpFGoCQicTP6qL1Z1uIAVqX0wH30KDSz5eMizZGSiMRN64xUrjyuNxO2HYt9Pw++/cjvkEQkwZREJK5+cXg3Pmh9AlusJW7uQ36HIyIJpiQicZWeEuayE/vzVPEQ3KKp3vW0RCRpKYlI3J1+0J7MaftTnHOUffio3+FIkhs9ejRmxjXXJP7GaBMnTsTMWLFiRYPqzZo1i7Fjx1JWVrZL+YoVKzAzJk6cGL8gKzEzxo4dm5Bjl1MSkbgLh4wLTjqaN0sPpvjDf0HxNr9DkiS1bds2nn32WQCeeuqpavdKD4pZs2Yxbty4akmkc+fOvP/++5xyyik+RRY7JRFJiBP67sGc9iNJL95A8WfP+h2OJKmXXnqJTZs2cfLJJ7NmzRpee+01v0NqkPT0dAYOHEhOTo7foURNSUQSwsw4+dSRLC7rxqa379dyX0mISZMm0bZtWyZOnEiLFi2YNGnSLtvHjh2LmbF06VJOOeUUWrVqRffu3bntttt26RVs376da665hv33359WrVrRqVMnTjvtNL744ovdvv5pp53GgAEDqpUvX76cUCjEww8/zNixYxk3bhzg3YrXzDAzoPbhrLfffpsTTjiB7OxsWrZsSf/+/XnssZ3XpXv66acZOnQoOTk5tGrVigEDBlRre2PR/UQkYQb1yuHR9mcyquAeCr+cRea+x/odklT2n+vhh8/9jaHTATD89qiqfvfdd7z55puMHj2anJwcfvrTn/LCCy9QUFBA27Ztd9n3Zz/7GRdffDHXXHMN06dP59Zbb2Wvvfbi4osvBqCoqIjNmzdz880307lzZ9avX8+DDz7IEUccweLFi6vd9KrcZZddximnnMKHH37IYYcdVlE+YcIEWrZsyXnnncfGjRv59ttveeyxx3jnnXcIh8O7bdfUqVM588wzGTx4MI888ggdOnRg4cKFrFy5smKfZcuWMXLkSK6//npCoRCzZ89m1KhRbNu2jUsvvTSqf89oKYlIQh1++qWs/79/su6/99JbSUTi6IknnqC0tJQLL7wQgIsuuogpU6bwzDPPVHsjve666yoSxvHHH8+MGTOYMmVKRVl2djaPPrpzEUhpaSnDhg1jjz32YMqUKbVO2p900kn07NmTRx55pCKJFBcX83//93+cd955tG7dmtatW9O1a1cADj/8cFJSan/bdc5x9dVXc9BBBzFz5syKOzIef/zxu+x34403VnxdVlbGkCFD+P7773nooYeURCS5HNCjE/9pdzonrnuKgtVLabtnb79DknJR9gCCYtKkSfTu3ZsjjjgC8N5ou3TpwqRJk6q9kVaduN5///359NNPdyl79tlnueuuu1iyZAkbN26sKF+yZEmtMYRCIS655BLGjRvH3XffTXZ2Ni+99BI//vgjl1xySYPbtGTJElauXFnRw6jN0qVLueWWW5g9ezY//PBDxdBcenp6g18zVpoTkYTb77SrcRhfTLvL71AkSeTl5bFo0SLOOOMMNmzYwIYNG9i8eTNnnHEGc+fO5csvd70dQbt27Xb5Pj09ne3bt1d8P336dM455xz69OnDU089xQcffMBHH31ETk7OLvvV5Ne//jWlpaVMnjwZgIcffpjDDjusxrmSuqxbtw6goudSky1btnDCCScwb948br/9dubMmcNHH33Er371K4qKihr8mrFST0QSbu+e+/JZ9hD6/TCV79eMp3PHDn6HJE1c+STyHXfcwR133FFt++OPP8748ePrfbynn36aXr167TLBXVxczPr16+us2759e84++2weeeQRhg0bxsyZM3cZGmuIDh28v43Vq2s/Sff9999n5cqVzJkzhyOPPLKi3K/lzeqJSKPoMvw6sqyQD1/8h9+hSBO3Y8cOpkyZwuGHH87MmTOrPQ466CAmT56Ma8CKwMLCwmpzFZMnT6a0tLRe9S+//HIWLFjAqFGjyM7O5uc///ku28uHmbZt2/05U/vssw89evTg0UcfrTX+wsJCwFvpVa6goICpU6fWK9Z4U09EGkXHPkeyquWBDPjuKb7+8Vp+skcbv0OSJuqVV15h3bp13HXXXQwZMqTa9ksuuYTLLruMWbNm1fuYJ510Ei+99BLXXHMNp556Knl5edx///20adOmXvUHDhzIgAEDmD17NldeeSWZmZm7bO/bty8Ad911F8OHDyccDpObm1vtOGbGvffeyxlnnMHQoUO59NJLycnJYfHixaxZs4Zx48YxaNAgsrKyuOKKKxg3bhxbt25l/PjxdOjQYZe5nMbiS0/EzPYys+fMbKOZbTKzF8ysWz3rdjOzSWa2ysy2mdmXZjbezFomOm6JTZvjr6Wb5TPzRd35UKI3adIkWrduzVlnnVXj9nPPPbfGc0Z2Z/To0dx0000888wznHbaabz66qtMnz6d7Ozseh+jPJ6aJtRPPfVULr/88oplw4ceemitxzn99NN54403AG++ZcSIEUyYMIEePXoAkJOTw4svvkhpaSkjR47khhtuYNSoUZx//vn1jjWerCFdvri8oFkmMA8oAm4GHDAeyAQOdM7Veku8SKL4FEgFxgKrgEOBccA059w5db1+bm6uy8vLi7EVEpWyUgr+1p+VhamkjJnJ/l3b+B1Rs7F48WL69OnjdxhJbfDgwYRCIebMmeN3KFGp63fEzD52zlXrPvkxnDUa6Ans65z7CsDM5gNLgUuAu3dTdzDQGxjmnHs9UjbTzNoBvzWzTOdcYeJCl5iEwmQccxUH/fd3/Hnac+x/+Si/IxKJSVFREZ988glvvvkm7733nm/zEn7yYzhrBDC3PIEAOOeWA+8Cp9dRNy3yvKlK+Qa8tlicYpQEaZF7AdtS23DY908yd9k6v8MRicn333/PoEGDuPvuu7nxxhsZMWKE3yE1Oj+SSD9gQQ3lC4G+ddR9E6/HcoeZ9TWzVmY2FLgaeHh3Q2ESEKktSBl4CSeEP+HJl19v0AoakaDp0aMHzjkKCgr485//7Hc4vvAjibQDCmooXw+0raG8gnNuO3AkXtwLgc3AW8DLwG9qq2dmY8wsz8zy8vPzo41b4iR14BhKQukMWvMMby1e43c4IhKDJnWeiJllAM8AHYELgGOA3wHnAP9bWz3n3ATnXK5zLrcpX3I5abTsgA34BWemzOGx1+ZSVqbeSGNQr09qE8vvhh9JpICaexy19VAq+zUwBDjZOfeEc262c+5O4DrgUjPrH9dIJWHCg64klVKOWP8C0+Z953c4SS81NbXOE92k+dq2bVvU193yI4ksxJsXqaovsKiOugcABc65r6uUfxh51hrGpqL9T2C/U/hl6ls8+Pp8dpSU1V1HotaxY0dWr15NYWGheiQCeL2P8ku7fPvtt7Rv3z6q4/ixxHcacKeZ9XTOLQMwsx54y3evr6PuD0BbM+tVeXUXcHjkufYLzkjg2KCryPriZQZveoVn8vpywcDufoeUtLKysgDvHhzFxcU+RyNBkZKSQkZGBt26dSMjIyOqY/hxsmFLvJMNt7HzZMM/Aa3xTjbcEtmvO/A1cJtz7rZIWQ9gPl4y+TPeyYa5wB+BL4HDnHO7/Uirkw2Dxf3rJNZ9u5QRoQd46/cn0iJt9zfsERF/1HayYaMPZ0WW4Q7Fe9OfDDwJLAeGlieQCAPClWN0zq0ABgKf4Z3l/ireyYsTgBPqSiASPHb0b+lQtpbB22Yw8b0VfocjIg3kywUYnXOrgDPr2GcFNZw86JxbBJydmMik0f3kOOjcn+vyX2H4rKH84vBuZLdIrbueiARCk1riK0nIDI66jk4lqxm8410mzK66ZkJEgkxJRPy332nQYR9uaP0q/3pnOWs27/5OciISHEoi4r9QCI68lj2LvmZw2cf874yv6q4jIoGgJCLBcMBIyO7GLdmv8tSHK/lmvS7GLNIUKIlIMIRTYfBVdCtcyMDQYu5580u/IxKRelASkeAYcAG07MhtbV/jxU9X8+WPm/2OSETqoCQiwZGaAYN+w96bPmJg2gr+/t8lfkckInVQEpFgyf0VZLThL+3/wxuLfuT9r3XjKpEgUxKRYElv7fVG1s/huKzV3PbyIkp1qXiRwFISkeA57BJo0Za/tHuZxd9v4pmPvvE7IhGphZKIBE9GFhzxG/b44W3O7ZLPXa8vYdN2XXlWJIiURCSYDvd6Ize2nMr6wh3c/9ZSvyMSkRooiUgwpbeGQVfS+psZXNdnExPfW8Gy/C111xORRqUkIsF12Bho0Y4xZc+SnhLmz68s9jsiEalCSUSCK701DL6KtBUzuO2QQt76Yg1vf5nvd1QiUomSiATboaMhsz0/3fA4e3doya1TF7C9uNTvqEQkQklEgi29FQy6itCyGdw7aAcr1hXy8Nu654hIUCiJSPAdOgoyO9B/6QOM6N+FB2d+zfK1W/2OSkRQEpGmIL0VHHUtLH+bcQfkk54a4o8vLcA5ncku4jclEWkacn8NWV1p+/7t/O7EfXjnq7VMn/+931GJNHtKItI0pGbAkOth9cec1+ZzDuyazZ9eXqQz2UV8piQiTUf/c6HDPoRnjOfPI/qybksRf31V546I+ElJRJqOcAoMvRnWLuGAdf9h9NE9mfLhN8zWuSMivlESkaalzwjoMgBm/ZVrhnSnV8dWXP/8fA1rifhESUSaFjM47lbY+A0Z8yZx51n9+WHTdv78soa1RPygJCJNz0+Ohb2Phtl3clDHMGOO/gnP5H3DrCVr/I5MpNlREpGm6bixULgW5j7E/xzfm94dW/H75+azbkuR35GJNCtKItI0dT0E+pwG795HRtE67vv5ADZsK+a6f8+jTLfTFWk0SiLSdB03Fkq2w6y/0rdLFn88pQ+zluTz2DvL/Y5MpNlQEpGmq0Mv70z2jyfCmi84f2B3hvXbgzte+4JPVxX4HZ1Is6AkIk3bMX+AtNbwxh8xM/52Zn86ZWdw2ROfsGbTdr+jE0l6SiLStLVsD0dfB0tfh69nkp2Zyj8vzGXT9mLGTP5Y9x4RSTBfkoiZ7WVmz5nZRjPbZGYvmFm3BtTvY2b/NrO1ZrbNzJaY2dWJjFkC7LBLoE03eP2PUFZKn85Z3H12fz77ZgM3vvC5rvYrkkCNnkTMLBOYAewHXARcAPQGZppZy3rUzwU+ANKBUcDJwF1AOFExS8ClZsDxY+HHz2HeFABO2r8z156wDy98upq//ucLJRKRBEnx4TVHAz2BfZ1zXwGY2XxgKXAJcHdtFc0sBDwOvOWc+1mlTTMTF640Cf3OgPcfhLf+BP1+BmktuXJoL9ZuKWLC7GVkt0jlimN7+R2lSNLxYzhrBDC3PIEAOOeWA+8Cp9dRdwjQh90kGmmmzGDYX2DLD/DeA5EiY+xp/fjpQV34+3+XcM8bX6pHIhJnfiSRfsCCGsoXAn3rqHtk5DnDzOaaWbGZrTGzf5hZi7hGKU1Pt8Oh7+nw7r2w6TsAQiHjzrP6M/KQrtz31lJunbaQUp2MKBI3fiSRdkBNi/jXA23rqNsl8vwM8DpwAvA3vLmRp+IVoDRhx4+DslJ4c2xFUUo4xN9HHsiYo3vy+PsrueCxD8jfrMujiMRDU1viWx7vE865W5xzs5xzdwLjgJ+aWZ+aKpnZGDPLM7O8/HzdeyKptdsbBv0G5j8D33xYUWxm3HhyH/428kA+XlnA8PvmMH3edxreEomRH0mkgJp7HLX1UCpbF3l+o0r565HnATVVcs5NcM7lOudyc3Jy6h2oNFFHXgutO8N/fg9lZbtsOjt3L166YjCdstO5csqn/OKfH/D+1+uUTESi5EcSWYg3L1JVX2BRPeruTlkd26U5SG8FJ9wG330Knz1ZbXOfzllMveJIxo3ox9I1Wzj3n3MZft8c/nfmVyxYvZGSUv0aidSXH0t8pwF3mllP59wyADPrAQwGrq+j7n+AImAYML1S+UmR57z4hipN1gFnwUePwn9vgtV5gO2yOYx3ktJ5BziWr93KsrVbWftWEfPegoVhI6tFKq3SUshMD5MWDpESNlJDIUIhwwwMI2TeYUNEyswwvIViIaulLPL1znLvOOXHNPPmcOL66W7ABd5Vj0USwI8k8k/gN8BUM7sZcMCfgG+AR8p3MrPuwNfAbc652wCcc+vM7K/AH81sE95Ji7nALcCkysuGpZkzg1Pvhed+BV+8WutuKXhnuvYOQWmWo7jUUVJaRskOR1mRo3RjbMNcLvIAaMgFWMoTjBmkhIxwyEgNh0gNh6qkwzpszYeyYiURSZh6JREzOwvvHI5WwHLgRefc7Ghe0Dm31cyGAvcAk/H+Xt4C/sc5t6Xyy+J9YKz6oew2YDNwOfBb4Hvg73iJSGSnPfrCFXPrvXuYmi97UFJaxtYdpWwpKqGktIzSMkeZc5SUOe/rMigpK6PMOUrL2GV7WWSfkkhZaVmlh/O2V95WUurYUlTC5u0lbN5ezNotRSxfu5VV6wspc9ApK4NRR+3NLwf1ICVcj/7K3f1A8z2SQHUmETO7FbgV7417JXAUcLWZfQac65xb0tAXdc6tAs6sY58VVB2D8Mod3smGOuFQGkVKOER2ixDZLVJ9i2HT9mLeWbqWye+vZPwri3nps9VMuCCXLm3qOD3KQkoiklD1GXq9Avg30N45d6Bzrj1eItkCfFjbsloRiZ+sjFROPqAzU8YM5KHzDmbl2kLOfOg9Vq0r3H1FA5wWCkji1CeJZAOPOedKygucc+8CxwCf4J3sJyKNZPgBnXn6koEU7ihl9ON5bC0q2c3exs5ZGZH4q08S+RbYq2phZFjpfrzrWYlII+rXJZsHfjGApWs2c9v03ayM13CWJFh9kshDwK1mtmct23X7OBEfHNU7h9FH9eSZvG/IW7G+5p3MNJwlCVWfJHIP3kl+C8zsT2Z2ROSmUkOA8XhLdkXEB1cf35su2RmMm76o5rPuLYSGsySR6kwizrlS4DTgDmAM8A6wAm9ZbktghZkNMDM/zjkRadYy01L4n+P34fPVG5m1pKbrwqknIolVrxNjnXMlzrnbgU7AIOBa4Em8FVoP4p0pvtnMPqz9KCKSCD87eE/2bNOC+2csrb5RcyKSYA3qPUQm0z+IPICK290OwDtz/OC4RicidUoNhxhzdE9unbaQ+d9u4MCubXZu1JyIJFjMl+hxzhU65951zt3nnLsoHkGJSMP87OA9aZEaZsqHq6ps0RJfSaymdj8REalBVkYqp/XvzNTPvmPz9uKdGzScJQmmJCKSJM45tBuFO0p5feGPOwvNlEQkoZRERJLEwd3asGebFrw8/7udhabhLEksJRGRJGFmnNq/M3OWrmVD4Y7yUk2sS0IpiYgkkdMO7EJJmeO1BT94BZoTkQRTEhFJIv26ZNGjfSavfP69V6AlvpJgSiIiScTMOLFfJz5Ytp4tRSW67IkknJKISJI5br+O7Cgt452l+WhORBJNSUQkyRzSvS1ZGSm8tXiNlvhKwimJiCSZlHCIIft2ZOaSNTg0nCWJpSQikoSO69ORtVt2sHVHqYazJKGURESS0DH75BAOGQXbSjScJQmlJCKShNpkpnFI97ZsUBKRBFMSEUlSQ/fryJaiMopKSvwORZKYkohIkjpuv444YGPFJVBE4k9JRCRJ9erYirTUMBsLi/wORZKYkohIkjIz2rZMZ8v2YrYXl/odjiQpJRGRJNa2ZTrOlTF32Tq/Q5EkpSQiksSyM9MJGcz8Yo3foUiSUhIRSWLhUIis9DAzlqzBaamvJICSiEgysxDZGSl8s34bX+dv8TsaSUJKIiJJzWidEQZghoa0JAGURESSmYVIC8F+nVp7V/UViTMlEZFkFrkU/ND9OpK3soCN24r9jkiSjC9JxMz2MrPnzGyjmW0ysxfMrFsUx7nezJyZvZOIOEWaPgO8JFJa5pizNN/vgCTJNHoSMbNMYAawH3ARcAHQG5hpZi0bcJyewM2A+ugitYncY31At7a0b5nGawt+8DsiSTIpPrzmaKAnsK9z7isAM5sPLAUuAe6u53EeAp4E9sWfdogEX2Q4Kxwyhh/Qiec/Xk3hjhIy0/QnI/Hhx3DWCGBueQIBcM4tB94FTq/PAczsF8DBwA0JiVAkWdjOOxueemAXthWXaoJd4sqPJNIPWFBD+UKgb12VzawtcA/we+fc+jjHJpJkrOLOhof2aMceWelMn/edzzFJMvEjibQDCmooXw+0rUf9vwNfAhPr+4JmNsbM8swsLz9fE4vSjFio4qZU4ZBxygFdmLUkn4Ktujy8xEeTWuJrZkcBFwKXuQZcw8E5N8E5l+ucy83JyUlcgCJBY7bLPdbPPrQrO0rLeO7jb30MSpKJH0mkgJp7HLX1UCp7BHgM+NbM2phZG7xJ9XDk+/S4RirS1FWaEwHYr1MWud3b8uQHKykr07W0JHZ+JJGFePMiVfUFFtVRtw9wKV6yKX8MBgZGvr4sfmGKJINdeyIA5w/szop1hbz79VqfYpJk4kcSmQYMjJznAYCZ9cBLBtPqqHtsDY95eBP1xwLPJSBekaYrssS3suEHdKJDqzQmzF7mU1CSTPxIIv8EVgBTzex0MxsBTAW+wRuuAsDMuptZiZndUl7mnJtV9QFsADZGvtdAr0hlVYazANJTwow5uidzlq4lb4UWOEpsGj2JOOe2AkPxVlhNxjthcDkw1DlX+VrVBoT9iFEkeVTviYA3pNWhVRp3vPaF5kYkJr6ctuqcWwWcWcc+K/ASSV3HGhKfqESSUKUlvpVlpqXwu2H78ofnP+ffH3/DOYc2+NJ1IoA+5YskN6PaxHq5s3P34vC92zH+lcUs0w2rJEpKIiLJrIY5kYpNZtx5Vn9SwyFGPZ7H2i1FjRubJAUlEZGkVn2Jb2V7tcvk4fMP4bsN2xj50Hss/XFzI8YmyUBJRCSZ1bDEt6rD9m7Hk6MOZ9P2Ek7+xxz+9PIilq/d2kgBSlOn60GLJDMLQfE2WPDCbnc7BJg1vISXPl3NR3Pf5873IKd1Gnu2yaRDqzQy08JkpIVJDYcImWF4+ckMDKt4jvznrYgp/9q89THldSivXx5itbo7v05PCdGtZx/oekgi/nUkDpRERJJZZnvYsRmeu7jOXbPwLkx3Yfm7QhHwY+Tho9JZaYT/uCaSgSRolEREktmQG+CAs+oc0qrLjtJSCotK2VZcCnhT9WUOypzDOYdz3ks4HI7I15HXrCh3kSl+B2WRo5SXucj/vON6hSXOsXDqPfwy5XVvByWRQFISEUlmoTDk7BvzYdIijzYxH6lhPn/jCdhBZHGApnCDSD8VEQmsFump3he7WWEm/lISEZHAalFxL3hdmiWolEREJLBSwmHvC/VEAktJREQCK6wkEnhKIiISWOFQ5C0qxtVlkjhKIiISWKFweRJRTySolEREJLDKh7OckkhgKYmISGClRIazikuURIJKSUREAqu8J1JUUuxzJFIbJRERCazyifXiyOVWJHiUREQksMrPE9lRUuJzJFIbJRERCaxQxXCW5kSCSklERAKrfDirVD2RwFISEZHACoW8nkhJqeZEgkpJREQCKxzy7iFSWqbhrKBSEhGRwCqfEykuUU8kqJRERCSwKuZEStUTCSolEREJrPI5kdIy9USCSklERALLIj2REq3OCiwlEREJrIqJ9VJdCj6olEREJLDKJ9Y1nBVcSiIiElgVE+tKIoGlJCIigVV+FV+tzgouJRERCayK1Vk6Yz2wfEkiZraXmT1nZhvNbJOZvWBm3epRL9fMJpjZF2ZWaGarzOxJM9u7MeIWkcalM9aDr9GTiJllAjOA/YCLgAuA3sBMM2tZR/WfA/2AfwDDgeuBg4E8M9srYUGLiC8qJtbVEwmsFB9eczTQE9jXOfcVgJnNB5YClwB376buHc65/MoFZvYusDxy3FsSErGI+CKl4ox1JZGg8mM4awQwtzyBADjnlgPvAqfvrmLVBBIpWwnkA3vGOU4R8dnOJb4azgoqP5JIP2BBDeULgb4NPZiZ9QE6AotjjEtEAiYcmVgvUxIJLD+SSDugoIby9UDbhhzIzFKAh/F6Io/tZr8xZpZnZnn5+dU6MyISUOGwLsAYdE19ie8DwCDgfOdcTYkJAOfcBOdcrnMuNycnp/GiE5GYhCJzImWaEwksPybWC6i5x1FbD6VGZnY7MAa4yDn3epxiE5EAMdMZ60HnRxJZiDcvUlVfYFF9DmBmNwF/AK50zk2OY2wiEiSm80SCzo/hrGnAQDPrWV5gZj2AwZFtu2VmVwHjgZuccw8kKkgRCQDTnEjQ+ZFE/gmsAKaa2elmNgKYCnwDPFK+k5l1N7MSM7ulUtnPgXuB14AZZjaw0qPBK7tEJOi8nkiZhrMCq9GHs5xzW81sKHAPMBnvt+Qt4H+cc1sq7WpAmF0T3UmR8pMij8reBoYkKGwR8UOkJ6IlvsHlx5wIzrlVwJl17LOC8o8hO8t+CfwyUXGJSMBYeU9ESSSomvoSXxFJZpoTCTwlEREJsEhPxGlOJKiUREQkuEwnGwadkoiIBFfFxLrzORCpjZKIiASXTjYMPCUREQkuLfENPCUREQmwyMS65kQCS0lERIIrMpy1efsOnwOR2vhysqGISL1EhrN+s+0hVtz2uM/BNG3rrR39bphJeko4rsdVEhGR4MrZjy39fsHab76jRPMiMdkSzsZ2vQhIXCiJiEhwpWbQ6qyHONDvOKRWmhMREZGoKYmIiEjUlERERCRqSiIiIhI1JREREYmakoiIiERNSURERKKmJCIiIlEz55rXdfrNLB9YGWX1DsDaOIbTFKjNzUNza3Nzay/E3ubuzrmcqoXNLonEwszynHO5fsfRmNTm5qG5tbm5tRcS12YNZ4mISNSUREREJGpKIg0zwe8AfKA2Nw/Nrc3Nrb2QoDZrTkRERKKmnoiIiERNSURERKLW7JOIme1lZs+Z2UYz22RmL5hZt3rWzTCzv5vZ92a2zczeN7OjEx1zrKJts5nlmtkEM/vCzArNbJWZPWlmezdG3LGI5edc5TjXm5kzs3cSEWc8xdpmM+tjZv82s7WR3+8lZnZ1ImOOVYx/z93MbFLk93qbmX1pZuPNrGWi446WmXU1s/sj7z2Fkd/NHvWsGzKzG8xshZltN7N5ZnZmQ2No1knEzDKBGcB+wEXABUBvYGY9f3EeA0YDtwCnAt8D/zWzgxIScBzE2OafA/2AfwDDgeuBg4E8M9srYUHHKA4/5/Lj9ARuBtYkIs54irXNZpYLfACkA6OAk4G7gPjeoDuOYmlzZPubwNHAH/Ha+yhwHfCvBIYdq17A2UABMKeBdf8EjAUewPt7ngv828xObtBRnHPN9gFcDZQCvSqV7Q2UANfWUbc/4ICLK5WlAEuAaX63LUFtzqmhrDtQBtzmd9sS0eYqx/kv8AgwC3jH73Yl8OccAhYBL/rdjkZs84mRv+cTq5TfHqmf6Xf7avtZVfp6VKQNPepRryNQBIyrUv4WML8hMTTrnggwApjrnPuqvMA5txx4Fzi9HnWLgWcq1S0BngaGmVl6/MONi6jb7JzLr6FsJZAP7BnnOOMplp8zAGb2C7xe1w0JiTD+YmnzEKAPcHfCokuMWNqcFnneVKV8A15StTjFGFfOubIoqw7Da/MTVcqfAA5oyBB1c08i/YAFNZQvBPrWo+5y51xhDXXT8LqZQRRLm6sxsz54n2oWxxhXIsXUZjNrC9wD/N45tz7OsSVKLG0+MvKcYWZzzazYzNaY2T/MrEVco4yvWNr8JrAUuMPM+ppZKzMbite7edg5tzW+ofquH15P5Ksq5Qsjz/V+L2juSaQd3lhiVeuBtjHULd8eRLG0eRdmlgI8jNcTeSz20BIm1jb/HfgSmBjHmBItljZ3iTw/A7wOnAD8DW+45Kl4BZgAUbfZObcdL3mG8N5IN+MN7bwM/Ca+YQZCO2CDi4xhVdLg96+UuIUkzdEDwCDgFOdcTX+8TZ6ZHQVcCBxcwx9csir/cPmEc+6WyNezzCwM3G5mfZxzQe55NpiZZeAlzY54E/KrgMPwFs2UAJf5F12wNfckUkDNn1Bq+0RTtW73WurCzoweNLG0uYKZ3Q6MAS5yzr0ep9gSJZY2P4LXy/rWzNpEylKAcOT7bc65ojjFGU+xtHld5PmNKuWv4000DyCYw5extPnXeHNBvZxzX0fKZpvZRmCCmT3snJsXt0j9VwC0MTOr8uGowe9fzX04ayHe2GBVffFWp9RVd+/IssKqdXdQfawxKGJpMwBmdhPwB+Aq59zkOMaWKLG0uQ9wKd4fXfljMDAw8nVQP6HG+ru9O9FO5iZaLG0+ACiolEDKfRh57hNjbEGzEG/59k+qlJfPhdTrvQCURKYBAyPr/wGInKgzOLJtd6YDqcBZleqmAOcArwf00ynE1mbM7CpgPHCTc+6BRAUZZ7G0+dgaHvPwJnCPBZ5LQLzxEEub/4M36TqsSvlJkee8OMUYb7G0+QegrZlVXRBzeOR5dbyCDIjX8FaXnlel/HxgQWRVW/34vc7Z5zXWLfF6DJ/jLQEcgfcGsQxoVWm/7njjordUqf803qfRUcBxeG8o2/HGz31vX7zbjHeyYRnem8zAKo++frctUT/nGo43i+CfJxLr7/atkfK/AMfjnVi6DZjod9sS0WagB97y3i/xTlQ8FvhdpCyPSudjBO0BjIw8HsI7T+SyyPfHVNqnBHisSr3bI+9X1+IN5T0U+fs+tUGv7/c/gN8PoBvwfOSXZTPwElVO1on8gjlgbJXyFnhr6X+I/DA+AIb43aZEtRlvdZKr5THL73Yl6udcw7ECn0RibTPeeRHXRt6Ud+DdUvo2INXvdiWwzX2BZ4Fv8BLml8CdQFu/21VHm+v8m4x8P7FKvTDeFRhW4vU85wMjG/r6uhS8iIhErbnPiYiISAyUREREJGpKIiIiEjUlERERiZqSiIiIRE1JREREoqYkIiIiUVMSERGRqCmJiIhI1JRERHxiZi3N7Asz+9DMUiuVn2hmZWZ2hZ/xidSHLnsi4iMzGwDMBe5xzl1vZnvgXTTwA+dcve7/LuInJRERn5nZNXgX+hsG/Bbv3hb9nXNrfQ1MpB6URER8ZmYGvAIMBdKAE5xzb/kblUj9aE5ExGfO+yQ3Ge9Oc/OUQKQpURIR8ZmZdQLuAz4B+pvZ1T6HJFJvSiIiPooMZU3CuynQ8cC9wB1mdqCfcYnUl+ZERHxkZtcBfwOGOufeNrM0vNVa6UCuc26brwGK1EE9ERGfmNnBePcw/6tz7m0A59wO4Fy8W7je7V90IvWjnoiIiERNPREREYmakoiIiERNSURERKKmJCIiIlFTEhERkagpiYiISNSUREREJGpKIiIiErX/B/EczNyvQir2AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "prim_var = cons_to_prim_vars(cons_var)\n",
    "analytical = np.loadtxt('/mnt/external_data/maruthi/my_codes/compressible_euler_1d/output/shock_tube.dat')\n",
    "plt.plot(x_grid, prim_var[:, 0], label='LLF')\n",
    "plt.plot(analytical[:, 0], analytical[:, 1], label='Analytical')\n",
    "plt.xlabel('x', fontsize=16)\n",
    "plt.ylabel(r'$\\rho $', fontsize=16)\n",
    "plt.legend(fontsize=16)\n",
    "plt.xticks(fontsize=16)\n",
    "plt.yticks(fontsize=16)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "31ddb5f0",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}