FiniteVolumeGPU/OpenCL to CUDA.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Lets have matplotlib \"inline\"\n",
    "%matplotlib inline\n",
    "\n",
    "# Add line profiler\n",
    "%load_ext line_profiler\n",
    "\n",
    "#Import packages we need\n",
    "import numpy as np\n",
    "from matplotlib import animation, rc\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "import subprocess\n",
    "import os\n",
    "import gc\n",
    "import datetime\n",
    "import importlib\n",
    "\n",
    "import pycuda.driver as cuda\n",
    "\n",
    "try:\n",
    "    from StringIO import StringIO\n",
    "except ImportError:\n",
    "    from io import StringIO\n",
    "\n",
    "#Set large figure sizes\n",
    "#Note, this prevents nice figures for articles...\n",
    "rc('figure', figsize=(16.0, 12.0))\n",
    "rc('animation', html='html5')\n",
    "\n",
    "#Finally, import our simulator\n",
    "#from SWESimulators import FBL, CTCS, LxF, FORCE, HLL, HLL2, KP07, KP07_dimsplit, WAF, CDKLM16, DataOutput, PlotHelper"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CUDA version (9, 1, 0)\n",
      "Driver version 9010\n",
      "Using GeForce 840M\n",
      " => compute capability: (5, 0)\n",
      " => memory: 2048.0 MB\n"
     ]
    }
   ],
   "source": [
    "class CudaContext(object):\n",
    "    def __init__(self, verbose=True, blocking=False):\n",
    "        self.verbose = verbose\n",
    "        self.blocking = blocking\n",
    "        \n",
    "        cuda.init(flags=0)\n",
    "        \n",
    "        try:\n",
    "            cuda.Context.pop()\n",
    "            if (self.verbose):\n",
    "                print(\"=== WARNING ===\")\n",
    "                print(\"Popped existing context\")\n",
    "                print(\"=== WARNING ===\")\n",
    "        except:\n",
    "            pass\n",
    "        \n",
    "        if (self.verbose):\n",
    "            print(\"CUDA version \" + str(cuda.get_version()))\n",
    "            print(\"Driver version \" + str(cuda.get_driver_version()))\n",
    "\n",
    "        self.cuda_device = cuda.Device(0)\n",
    "        if (self.verbose):\n",
    "            print(\"Using \" + self.cuda_device.name())\n",
    "            print(\" => compute capability: \" + str(self.cuda_device.compute_capability()))\n",
    "            print(\" => memory: \" + str(self.cuda_device.total_memory() / (1024*1024)) + \" MB\")\n",
    "\n",
    "        if (self.blocking):\n",
    "            self.cuda_context = self.cuda_device.make_context(flags=cuda.ctx_flags.SCHED_BLOCKING_SYNC)\n",
    "            if (self.verbose):\n",
    "                print(\"=== WARNING ===\")\n",
    "                print(\"Using blocking context\")\n",
    "                print(\"=== WARNING ===\")\n",
    "        else:\n",
    "            self.cuda_context = self.cuda_device.make_context(flags=cuda.ctx_flags.SCHED_AUTO)\n",
    "            \n",
    "    \n",
    "    def __del__(self, *args):\n",
    "        if self.verbose:\n",
    "            print(\"Cleaning up CUDA context\")\n",
    "            \n",
    "        self.cuda_context.detach()\n",
    "        cuda.Context.pop()\n",
    "\n",
    "            \n",
    "my_context = CudaContext(verbose=True, blocking=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=> sleep 125.088930 ms\n",
      "=> elapsed time: 0.125089 s\n"
     ]
    }
   ],
   "source": [
    "import time\n",
    "class Timer(object):\n",
    "    def __init__(self, tag, verbose=True):\n",
    "        self.verbose = verbose\n",
    "        self.tag = tag\n",
    "        \n",
    "    def __enter__(self):\n",
    "        self.start = time.time()\n",
    "        return self\n",
    "    \n",
    "    def __exit__(self, *args):\n",
    "        self.end = time.time()\n",
    "        self.secs = self.end - self.start\n",
    "        self.msecs = self.secs * 1000 # millisecs\n",
    "        if self.verbose:\n",
    "            print(\"=> \" + self.tag + ' %f ms' % self.msecs)\n",
    "            \n",
    "with Timer(\"sleep\", verbose=True) as t:\n",
    "    time.sleep(0.125)\n",
    "    \n",
    "print(\"=> elapsed time: %f s\" % t.secs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def gen_test_data(nx, ny, num_ghost_cells):\n",
    "    width = 100.0\n",
    "    height = width\n",
    "    dx = width / float(nx)\n",
    "    dy = height / float(ny)\n",
    "\n",
    "    h = np.zeros((ny+2*num_ghost_cells, nx+2*num_ghost_cells), dtype=np.float32); \n",
    "    hu = np.zeros((ny+2*num_ghost_cells, nx+2*num_ghost_cells), dtype=np.float32);\n",
    "    hv = np.zeros((ny+2*num_ghost_cells, nx+2*num_ghost_cells), dtype=np.float32);\n",
    "\n",
    "    #Create a gaussian \"dam break\" that will not form shocks\n",
    "    x_center = dx*nx/2.0\n",
    "    y_center = dy*ny/2.0\n",
    "    size = width \n",
    "    for j in range(-num_ghost_cells, ny+num_ghost_cells):\n",
    "        y = dy*(j+0.5) - y_center\n",
    "        for i in range(-num_ghost_cells, nx+num_ghost_cells):\n",
    "            x = dx*(i+0.5) - x_center\n",
    "            \n",
    "            h[j+num_ghost_cells, i+num_ghost_cells] = 0.5 + 0.1*np.exp(-(x**2/size + y**2/size))\n",
    "            #hu[j+num_ghost_cells, i+num_ghost_cells] = 0.01*np.sin(x)*np.exp(-(x**2/size))\n",
    "            hu[j+num_ghost_cells, i+num_ghost_cells] = 0.1*np.exp(-(x**2/size + y**2/size))\n",
    "            hv[j+num_ghost_cells, i+num_ghost_cells] = 0.1*np.exp(-(x**2/size + y**2/size))\n",
    "    \n",
    "    return h, hu, hv, dx, dy, nx, ny"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "#%lprun -f gen_test_data gen_test_data(100, 150, 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=> upload (async) 8.011341 ms\n",
      "=> download (async) 20.012617 ms\n",
      "=> sync 0.000000 ms\n",
      "Sum of absolute difference:  0.0\n"
     ]
    }
   ],
   "source": [
    "from SWESimulators import Common\n",
    "importlib.reload(Common)\n",
    "\n",
    "nx = 1000\n",
    "ny = 1500\n",
    "nx_halo = 2\n",
    "ny_halo = 3\n",
    "a = np.random.rand(ny+2*ny_halo, nx+2*nx_halo).astype(np.float32)\n",
    "\n",
    "import pycuda.driver as cuda\n",
    "stream = cuda.Stream()\n",
    "\n",
    "with Timer(\"upload (async)\", verbose=True) as t:\n",
    "    a_gpu = Common.CUDAArray2D(stream, nx, ny, nx_halo, ny_halo, a)\n",
    "\n",
    "with Timer(\"download (async)\", verbose=True) as t:\n",
    "    b = a_gpu.download(stream, async=True)\n",
    "    \n",
    "with Timer(\"sync\", verbose=True) as t:\n",
    "    cuda.Context.synchronize()\n",
    "    \n",
    "print(\"Sum of absolute difference: \", np.sum(np.abs(a-b)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=> compile 6411.859989 ms\n"
     ]
    }
   ],
   "source": [
    "with Timer(\"compile\", verbose=True) as t:\n",
    "    module = Common.get_kernel(\"FORCE_kernel.cu\", 16, 16)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "#%lprun -f Common.get_kernel Common.get_kernel(\"FORCE_kernel.cu\", 16, 16)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=> construct 4934.113741 ms\n",
      "=> step 14.986992 ms\n",
      "=> download 2.002239 ms\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x94b81290f0>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZgAAADxCAYAAAD2t6FzAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3X+wXWV97/H3Jyc/kB8KIYr8qngrVy9jkWoK7TjtoBQIjA128EfwTpv0wqQ6ZWhnbkehP8AJbUesrdarV27ADGhbxOZe9HBFwwHlau9Vm0ABCT+alEE5JAVCUoRSAkk+94+1Nix29j77OXuvnb3Ofr6vmTXsvfaznrXW4bvyXb+e55FtQgghhLrNG/UGhBBCGE+RYEIIIQxFJJgQQghDEQkmhBDCUESCCSGEMBSRYEIIIQxFJJgQQghDEQkmhBDCUESCCSGEMBTzR70BIXRy9rsO8VM79yaVvfPe3RtsLxvyJoVQi9TYHoe4jgQTGumpnXv5hw0/k1R24ugtS4a8OSHUJjW2xyGuI8GERjKwj32j3owQapdTbEeCCY1kzItOu0UWwlySU2xHggmNlctZXshPLrEdCSY0kjF7YyiJMIZyiu1IMKGx9pHHQRjyk0tsR4IJjWRgbyYHYchLTrEdCSY0Vi5neSE/ucR2JJjQSAZezOQ+dchLTrEdCSY0knE2txFCXnKK7UgwoZkMe/M4BkNuMortSDChkYrWziGMn5xiOxJMaCixF416I0IYgnxiOxJMaKTiQWgeB2HIS06xHQkmNFLRViCPgzDkJafYjgQTGmtfJmd5IT+5xHaMaBkaqXWWlzKFMJekxnYKScskPSRpq6RLO/y+StKTku4up4sqv62UtKWcVlbmf0vSPZI2S7pa0kQ5/+OSHqvUdW6v7YsrmNBIRuyN858whuqK7fIf/s8DZwLTwEZJk7bvbyt6o+2L25ZdDFwBLKXIeXeWy+4CPmD7p5IErAfeD3ylXPTTtj+Vuo1xBIfG2mclTb30e5Yn6V2VeXdLel7Se4ewqyEzdcQ1cCqw1fbDtl+gSALnJW7C2cCU7Z1lUpkClgHY/mlZZj6wEPpvFRoJJjSSES94ImmaSeUs7xzgJOACSSd1KHqj7VPK6VoA299pzQPeDTwH3FrrjobspMZ2gmOBRyvfp8t57c6XdK+k9ZKOT1lW0gbgCeAZiquYlovLutZJOqLXBkaCCY1UNEablzT1MMhZXtX7gG/afq6PZUN4SWpsA0skbapMq9uq6nSZ0361cTNwgu2TgduA61OWtX02cDSwiOLkCuALwM8CpwDbgb/ota+RYEJjzeIh/0wH4iBneVUrgBtq27mQtcS43mF7aWVa21bNNFCN1eOAbdUCtp+yvbv8eg3wjlks+zwwSXlCZvtx23tt7yvrOrXXfjbyIf9CLfJBHDLqzQhD9Dz/xgve3fVGsy32Ovn8Z4ftpV1+Sz3Lu8H2bkkfpjjLa521Ielo4OeADakb1EnEdR6eYdcO26/t9vssY3smG4ETJb0ReIziJOhD1QKSjra9vfy6HHig/LwB+LPKba6zgMskHQocZnu7pPnAucD3OtT168B9vTZwoAQjaRnwV8AEcK3tT7T9vgj4EkXWfAr4oO1HetV7EIdwms4YZNNCw/3Qt/css6+eV5CTzvIqX68BrmqL7fuBm2y/2CrUT2xHXOfhNq//ca8ydcS27T2SLqZIFhPAOtubJa0BNtmeBC6RtBzYA+wEVpXL7pR0JUWSAlhTzjsKmCzjewL4NnB1WeaTkk6hOEF7BPjtXtvYd4JJfEXuQmCX7TdJWgFcBXyw33WGfBQPQmu5wO73LK8a27uA32mrN2I79KXG2Mb2LcAtbfMur3y+DLisy7LrgHVt8x4HfqFL+d+Y7fYNcp2W8vD0PF5+qLQeOKN8tzqEGdX1kN/2HqB1lvcA8NXWWV55ZgfFWd5mSfcAlwCfoYxt4BiKs79j2qqO2A59mcVD/jlvkDTa6eHpad3KlJdzTwNHAjsGWG/IxN6autOY7VmepPfxctw+Ut6GiNgOtakrtptukAST8vA0pUxRsHjzZzXAQRw8wGaFcTDilvy1xXbEdWiXUy8VgySYng9PK2WmyzcSXkPxoGk/5St4awFercWZjPcWZrKvnjdt+lFbbEdch05GGNsH1CB7+dLDU0kLKR6eTraVmQRanai9D/i27TjIQk9Fh4DzkqYhiNgOQ5Ma2+Og7yuYxFfkvgh8WdJWirO7FXVsdBh/RryY1l1G/euO2A5DNMrYPtAGelcu4eHp8xQ9cYYwKzZ1NUbrc/2ji+2f/N3PJZWT0i6Y5s1Lv7CaNy9ttPj5qeUm0sotmNibVi5xval/G4B5iWXn/+pPkuucyahj+0BqZEv+EEB1NbQMoWHyie1IMKGRTD5neSEvOcV2JJjQWOPyoDOEdrnEdiSY0EgmedClEOaUnGI7EkxoJAMv1tRfUwhNklNs57GXYQ56aUyMEMZMPrEdCSY0ksmntXPIS06xHQkmNFYuZ3khP7nEdiSY0Ei2sjnLC3nJKbYjwYRGKh6E5tGdRrvU1vSpo8+k1gcwkdiqfSKxd4CJxHWnrje1hX5qTwMA8zp38D40OcV2JJjQULWNWx5Cw+QT25FgQiMVD0LzuE8d8pJTbEeCCY2VS2vnkJ9cYjuPvQxzTqu1c8oUwlySGtspJC2T9JCkrZIu7fD7KklPSrq7nC6q/LZS0pZyWlmZ/y1J90jaLOlqSRPl/MWSpsryU5KO6LV9kWBCY+1jXtLUy4AH4c9IulXSA5Lul3RCrTsZslRTXE8AnwfOAU4CLpB0UoeiN9o+pZyuLZddDFwBnAacClxRSRgfsP024K3Aa3l5WIpLgdttnwjcXn6fUdwiC41kw4v7Bj//qRyEZ1IMc7xR0qTt+9uK3mj74g5VfAn4U9tTkg4F0l9PCqGDumKbIjFstf0wgKSvAOcB7bHdydnAlO2d5bJTwDLgBts/LcvMBxbCS6/ZnQecXn6+HrgD+NhMK4krmNBIxW2EeUlTDy8dhLZfAFoHYU/l2eB821MAtp+1/dwg+xVCamwDSyRtqkyr26o6Fni08n26nNfufEn3Slov6fiUZSVtAJ4AngHWl7OPsr0doPzv63rta98JRtLxkr5T3jrYLOl3O5Q5XdLTlVsPl3eqK4RO9pZ9NvWamPlA7OcgPFXSdyiGTH5beT/6HyX9eeV+dMR26FtiXO+wvbQyrW2rptODmvZGPTcDJ9g+GbiN4sqj57K2zwaOBhYB7+5jF4HBbpHtAf6r7bskHQbcKWmqw62H79l+zwDrCRma5aucO2wv7fJb6kF4g+3dkj4M/CVwCfAfgC8CrwaWAx8HVpXzIGI79KHG15SngeMr348Dtr1iXfZTla/XAFdVlj29bdk72pZ9XtIkxRX/FPC4pKNtb5d0NMUVzoz6voKxvd32XeXnZ4AH6HxmGEIfartFlnQQ2t5dfr0G+LkytqeBfwTuBV4PfA14ew07F7KWfIusl43AiZLeKGkhsAKYfMWaikTQspzi32kors7PknRE+XD/LGCDpENby0iaD5wLPFguMwm03jZbCXy91wbW8pC/fLPm54Efdvj5lyTdQ3FQ/77tzV3qWA2sBjiIg+vYrDDH1TRu+UsHIfAYxUH4oWqB1llZ+bV6EG6kuM/8JorY/itgU2XRnrHdT1xPTCR2FZNUKr27FoD5ietOrXPRxN5a61swL62+WXUVk9j9zO7eRZLVEdu290i6mCJZTADrbG+WtAbYZHsSuETScoo7TjsprsCxvVPSlRQxDrCmnHcUMClpUVnnt4GryzKfAL4q6ULgJ7z8dllXAyeY8s2a/wn8XuXtg5a7gDfYflbSuRRngCd2qqe8v7gW4NVafGA7BwqNU7xpM3h/TYMchMCrKI6RPcD/Be6kuMKBxNiOuA7t6ortoi7fAtzSNu/yyufLgMu6LLsOWNc273HgF7qUfwo4YzbbN1CCkbSAIrn8je3/1WGDflr5fIuk/y5pie0dg6w3jL86h5Xt5yCsxPbVtv+yQ50R26EvMWRyAkmieNj5QKcDsCzzeuBx25Z0KsUzn6c6lQ2hXU23yGYtYjsM26hi+0Ab5ArmncBvAD+SdHc57w+AnwGwfTXwPuAjkvYA/w6ssB23CUJPI+4QMGI7DE10dpnA9t/T4zmj7c8Bn+t3HSFvoxqUKWI7DFsMOBbCCNliTyYHYchLTrEdCSY0Vi63EUJ+contSDChkXK6Tx3yklNsR4IJjZXLQRjyk0tsR4IJjZRTW4F2i+antVZXYgv02bTkX5Daoj6xhX5qy/uFifXN19xvyZ9TbEeCCY2VS1uBkJ9cYjsSTGgkG/bUMyhTCI2SU2xHggmNlctthJCfXGI7EkxopJzuU4e85BTbkWBCYzmTgzDkJ5fYjgQTGiuXB6EhP7nEdiSY0Eh2PvepQ15yiu1IMKGhxN5M3rQJuckntiPBhMbK5T51yE8usR0JZsi0YGGt9fnFF2qtr6ly6q+p3SEL0/4fp7bkT22pDvW3vF84b09iubT6Fk2k1TeP9H1O/fvsSq5xZjnFdh7XaWHucXGvOmUKYU6pMa4lLZP0kKStki7t8PsqSU9KurucLqr8tlLSlnJaWc47WNI3JD0oabOkT6TU1c3ACUbSI5J+VK5wU4ffJemz5R/gXklvH3SdIQ/7UNLUSz8HYSWuLem5cv5kZZmI69C3muJ6Avg8cA5wEnCBpJM6FL3R9inldG257GLgCuA04FTgCklHlOU/ZfstwM8D75R0zkx1zaSuW2Tvsr2jy2/nACeW02nAF8r/htCVa3oQWjkIzwSmgY2SJm3f31b0RtsXV5b7I+BdwCO2D+1QdcR16EtdsU2RGLbafhhA0leA84D22O7kbGDK9s5y2Slgme0bgO8A2H5B0l3Acf1u4IG4RXYe8CUXfgAcLunoA7DeMMfVdIvspYPQ9gtA6yAcVMR16FtiXC+RtKkyrW6r5ljg0cr36XJeu/PLq+z1ko5PXVbS4cCvAbf3qKurOhKMgVsl3dnhDwDpf4QQXsFW0sTMB2K/B6GBW4FDyttlP5D03j7qDWE/iXG9w/bSyrS2rZpO99HaT7luBk6wfTJwG3B9yrKS5gM3AJ9tXSHNUFdXdSSYd9p+O8Utg9+R9Cttv6f8EZC0uvUPxIu1jbwQ5qriLC45wcx0IPZ7ELbi+mTgaeDTwGck/ews6o24DvtJje0E00D1KuI4YNsr1+WnbLcC7xrgHYnLrgW22P5MQl1dDZxgbG8r//sEcBPFLYmqnn+Ecvm1rX8gFrBo0M0KY2CflTT10NdBWInrH1HE9fHAHRQPPpPqLZePuA77qSGuATYCJ0p6o6SFwApgslqg7bbtcuCB8vMG4CxJR5QP988q5yHpT4DXAL+XWFdXAyUYSYdIOqz1udzI+9qKTQK/Wb5184vA07a3D7LekIeansH0cxA+JOmwtoPvJ8A7efkBasR16Fsdrynb3gNcTJEYHgC+anuzpDWSlpfFLilfN74HuARYVS67E7iS4vjYCKyxvVPSccAfUryVdlfb68gd65rJoG+RHQXcJKlV19/a/pakD5c7cTVwC3AusBV4DvitAdcZMmDEvhretLG9R1LrIJwA1rUOQmCT7UmKA2c5sAfYCfwB8PfAwRRXKbuAPwa+D/wKRZKJuA59qSu2AWzfQhGL1XmXVz5fBlzWZdl1wLq2edN0vv07Y13dDJRgyoc/b+sw/+rKZwO/M8h6mkiL0m537Hv7W2pd77y7Hkwq591z/35/XW0o+zwI94vrtuWHFteHLkz7f5faAj11HHtIH8s+ueV9Ykv+1Bb6qT0DLJjFPqf+HX+SXGNvubQPjq5iQjM5n/6aQmYyiu1IMKG5cjnNC/nJJLYjwYTGyuUsL+Qnl9iOBBMaycC+fXkchCEvOcV2JJjQTAYyOcsLmckotiPBhMaKrvjDuMoltiPBhObK5CAMGcoktiPBhIZK7o8phDkmn9iOBBOaK5OzvJChTGI7EkyfJl7/uqRy2z+WNr56quN+N229e378aO9CTWZwJm/atDt84b8nlZuntFb3CxLLAcyvu4V+zeUOmvdiUrnZteRP/fscllznjDKK7UgwocHyOAhDjvKI7UgwobkyuY0QMpRJbEeCCc2VyUEYMpRJbEeCCc2UUWO0kJmMYjsSTGisXBqjhfzkEtuRYEJzZfKmTchQJrHd97Bqkt5cDqfZmn4qqX0M59MlPV0pc3m3+kJoJ6dNta83YjsM2SjiehT6voKx/RBwCoCkCeAx4KYORb9n+z39ridkyozsQWjEdhiqEcb2gVbPwNBwBvDPtn9cU30heyoehKZMvWqSlkl6SNJWSZd2+H2VpCcrVyMXVX4+A3gE+H+SPlfb7oWM1RPXMFhsS1opaUs5rSznHSzpG5IelLRZ0icq5RdJurFc1w8lndBr++p6BrMCuKHLb78k6R5gG/D7tjfXtM6R2veaQ5LK3XNqtz9Lf855zYpa62u0Gs7yyiuQzwNnAtPARkmTtu9vK3qj7Ys7VLEC+Hfg/3T4bSixvXjhvyWVSx1Lfjat2lPLpra8P3hid1K5g5RW34Lkcun7PHGgW/LDyGNb0mLgCmBpuTV3SpoEdgOfsv0dSQuB2yWdY/ubwIXALttvkrQCuAr44EzbOPAVTLkRy4G/6/DzXcAbbL8N+G/A12aoZ7WkTZI2vUhaUIYxty9xmtmpwFbbD9t+AfgKcF7K6svY/nXgX4Fb235Oiu2I69DR4HENA8Q2cDYwZXun7V3AFLDM9nO2vwNQ1nkXcFy5zHnA9eXn9cAZkma81KrjFtk5wF22H2//wfZPbT9bfr4FWCBpSadKbK+1vdT20gUsqmGzwpzWaisw+K2EY4Fqx2zT5bx250u6V9J6SceX884tt+SS/TYvMbYjrsN+UmMblrROTsppdVtNg8R2z2UlHQ78GnB7+zK29wBPA0fOtKt1JJgL6HJ7TNLrWxlO0qnl+p6qYZ0hA7N4i2ymA7FTBmq/QXEzcILtk4HbePks7Y+B22zv13NoxHYYRGJc72idnJTT2vZqOlSdGtszLitpPsW/65+1/fAs1vcKAz2DkXQwxf2/367M+zCA7auB9wEfkbSH4j72CjuXJkZhYOmRssP20i6/TQPHV74fR/HM5OXV2NXEcA1wVRnbbwVeJ+kR4FDgEElvsr2MiO0wiHoipa/Yrix7etuyd1S+rwW22P5Mh/VNlwnoNcDOmTZwoARj+znaLpHKxNL6/Dkg3rwJo7QROFHSGyleN14BfKhaQNLRtreXX5cDD5SxvahSZhWwtPWwNGI7NEBfsV1+3gD8maQjyu9nAZeVy/wJRfKovk0JMAmsBL5PcYL17V4nVdGSPzRWHY3NbO+RdDHFATUBrLO9WdIaYJPtSeASScuBPRRnZKsGX3MI3Y06tm3vlHQlRZICWFPOOw74Q+BB4K7yLvDnbF8LfBH4sqStZV09X2mNBBOaydTWnUb5EP6WtnmXVz5fRnn2NkMd1wHX1bJBIW8NiW3b64B1bfOm6TJYje3ngffPZvsiwYTmiicaYVxlEtuRYEJjjUt/TCG0yyW2I8GE5srkIAwZyiS2I8H0ad7Tad15vO0fLqh1vcc9/WxSudTOLxotk4Ow3ZEL0mIrtYuT2XSbsmjei0nlDlJiuZrrS92X1C5lYDZdxRyTXGdPmcR2JJjQSOPUZXkIVTnFdiSY0FyZDMoUMpRJbEeCCY2Vy1leyE8usR0JJjRXJgdhyFAmsR0JJjRTRvepQ2Yyiu1IMKG5MjkIQ4Yyie1IMKGxkt8eDWGOySW26xgPJoQQQthPXMGE5srkNkLIUCaxHQmmT3v/5YmkckddtXgk653zMnoQ2m7JgmeSyk0k9tcwm5b8qS3g626hX3fPAAtJ3+d5B/p+VUaxHQkmNFcmB2HIUCaxnfQMRtI6SU9Iuq8yb7GkKUlbyv8e0WXZlWWZLZJW1rXhIQNOnPoUcR1GZohx3SSpD/mvA5a1zbsUuN32icDt5fdXkLQYuAI4DTgVuKLbARtClSjetEmZBnAdEdfhAEuN7XGQlGBsf5diiMyq84Dry8/XA+/tsOjZwJTtnbZ3AVPsf0CHsD+/3Clgr6nvVURch1EYclw3ySDPYI6yvR3A9nZJr+tQ5ljg0cr36XLefiStBlYDHMTBA2xWGBujOcgirsPwjUkC6WXY7WA6dRna8U9re63tpbaXLmDRkDcrzAk1PYORtEzSQ5K2StrvlhdwuKQnJd0NHCbponK5N5Tf75a0WdKHW1V22dr9Z0Zch05qegbTK7YlrWrFdjldVPmt43NESX8q6VFJz6bW1c0gCeZxSUeXKz4a6PT+7DRwfOX7ccC2AdYZMlLHLTJJE8DngXOAk4ALJJ3UoeiNtk8BtgLfKOcZeLicfxpwqaRjiLgOA6rjFtlsY7ucri2Xnek54s3lvE72q2smgySYSaCV9VYCX+9QZgNwlqQjyo0/q5wXQm/1XMGcCmy1/bDtF4CvUDxn6aYa1x8CvlZ+XsTLx0vEdRhMPVcws43tqq7PEW3/oHWbeFCprynfAHwfeLOkaUkXAp8AzpS0BTiz/I6kpZKuLTd0J3AlsLGc1pTzQpiZa3uLrOvzkkpcvx74iKRHgTcD76nE9Zck3Qs8Bjxqe1vEdRhIYmwDSyRtqkyr22pKfRZ4vqR7Ja2X1LryTn6OmFBXV0kP+W13G1j+jA5lNwEXVb6vA9alrGcu8e7dSeXmbXqg3vW++EKt9TVa+oPQJZI2Vb6vtb22/Nz1eUkrriUdCTxre3f5nOUD5WvKLSeXt8a+Juko248PM66PnHi2dyHSW6AvHEZL/rpb6CevN21fFsziPd+JUTxxT1vlDttLZ/g95VngzcANldi+Hnh34rLtutXVVXR2GRprFs9gdrQepJfT2ko1PZ+X2H7KduuM4RrgHe3bYnsbsBn45Rp3MWSqpteUB4ntWT9HTDlO2kWCCc1VzzOYjcCJkt4oaSGwguI5y0taL6uUlgMPlPOPk/Sq8vMRwDuBhwbapxCgrmcwfcc2fTxHnKGurqIvstBMs3hVc8Zq7D2SLqY4eCaAdbY3S1oDbLI9CVwiaTmwh6Lh5apy8f8E/IUkU9xS+JTtHw2+VSFrDYht2zsltZ4jQuU5oqRPUrzgcrCkaeBa2x/vVtdMIsGERhL1tWa2fQtwS9u8yyufLwMu67DcFHByPVsRQqEJsV3+1vE5ou2PAh/tML9rXd1EggmNNS7dZYTQLpfYjgQTmiuTgzBkKJPYjgQTmiuTgzBkKJPYjgQTmin9Vc0Q5paMYjsSTGiuTA7CkKFMYjsSzJBl1fK+ZuMy6NJsHT7xXFK5eaT9gSZm8YccVcv71HKLOrU/72AirVhRVomV1iiX2I4EExorl9sIIT+5xHYkmNBMNTVGC6FxMortSDChuTI5CEOGMontSDChkeps7RxCk+QU25FgQmNpXyZHYchOLrEdCSY0U0b3qUNmMortnt31S1on6QlJ91Xm/bmkB8uRzW6SdHiXZR+R9CNJd7cNCBVCT7MYD6a/+iO2w4gMM66bJGU8mOsox2qumALeavtk4J+YuYfNd9k+pcfIbCHsr57xYGZyHRHbYRSGG9eN0TPB2P4uRd//1Xm32m61tPoBxWhoIdRq2FcwEdthVHK5gqnjGcx/AW7s8puBW8sBm/5H21C2Icxs9AfZSGL74Hm7excCJhJb8i9IbCUP9be8PyjxX8qDElvTL+g4lHyHckofrHfeKAb2HX1sHxADJRhJf0gxutnfdCnyTtvbJL0OmJL0YHnW2Kmu1cBqgIM4eJDNCuPAo+1Oo67YjrgO+xlxbB9IfaduSSuB9wD/2XbHfGx7W/nfJ4CbgFO71Wd7re2ltpcuYFG/mxXGRKutwChuJdQZ2xHXoV1qbI+DvhKMpGXAx4Dltjv2zCfpEEmHtT4DZwH3dSobQkd22lSjiO1wQBzguB6VlNeUbwC+D7xZ0rSkC4HPAYdR3Bq4W9LVZdljJLXGhz4K+HtJ9wD/AHzD9reGshdhLNV1BSNpmaSHJG2VdGllfiu23yJpr6SfAF8HXk8R2w9J+hdJmyXdL+kfy0UjtsNA6rqC6Rbbld9XSXqy/Hf6bkkXVX5bKWlLOa2szP9TSY9KeratrkWSbizX9UNJJ/Tavp7PYGxf0GH2F7uU3QacW35+GHhbr/pD6KimVzUlTQCfB84EpoGNkiZt39+KbUmrgKW2L25b9j8Ctr1F0jHAnZIOj9gOAzkAsd1W9MYOsb0YuAJYWm7NneWyu4CbKS4itrTVcyGwy/abJK0ArgI+ONM2juD1iRDSaF/a1MOpwFbbD9t+AfgKcF7K+m3/k+0t5edtwBPAa/vfoxAKNcQ1DBDbwNnAlO2dZVKZomwTZvsHtrd3WOY84Pry83rgDGnm1/8iwYTGqinBHAs8Wvk+Xc5rd37Zen+9pOP32xbpVGAh8M997k4IL6kpwQwS26nLdlxf2VbsaeDImRaIBBOayczmIf8SSZsq0+pKTZ3OsNpvUNwMnFC23r+Nl8/Sigqko4EvA79lO5MXTMPQpMb2zHENg8V2yrLtZr1MdHYZGmsWr2rumKG7lmmgekVyHLCtWsD2U5Wv11DcWy62QXo18A3gj2z/IHmLQphBYmzPFNcwWGxPA6e3LXtHj+1prW9a0nzgNbT1hNEuEkxornre1NwInCjpjcBjwArgQ9UCko6u3HNeDjxQzl9I0cblS7b/rpatSXCQXkwqN5H4B1owi1Z9o2qhvyix5f0CJtLKKa0cwLzE3gFqNeLYBjYAfybpiPL7Wczc7x7AJLCS4s3L9wHf7tZOrCUSTGikugZlsr1H0sUUB9QEsM72ZklrgE22J4FLJC2naLm/E1hVLv4B4FeAI8s3zQBW2b578C0LuWpCbNveKelKiiQFsMb2TgBJn6RIVAdLmgautf1xireHvyxpa1nXil7bGAkmNJNd26BMtm8Bbmmbd3nl82V0OHuz/dfAX9eyESG0NCC2y9/WAes6zP8o8NEO858H3j+b7YsEE5prPBozh7C/TGI7EkxorHHpjymEdrnEdiSY0EwGMhm3PGSGy78KAAAImElEQVQmo9iOBBOaK49jMOQok9iOBBMaK5fbCCE/ucR2JJjQWHW9aRNC0+QS25FgQjPV1ONsCI2TUWxHggmNVDRGy+QoDFnJKbYjwYTmyrRbyYWJOz4v8Ub+gln8IVO7n0ntJTe13ERidy2pXcDMpquYkcgktlNGtFwn6QlJ91XmfVzSY5VR0s7tsuyMo62FMBPZSVPf9UdshxEZZlw3ScoJxnWUA9G0+bTtU8rplvYfK6OtnQOcBFwg6aRBNjZkxLOY+ncdEdvhQBt+XDdGzwRj+7v06JK5i0FGWwvZK/prSpn6XkPEdhiJ4cZ1kwwy4NjF5Shp6ypdPlf1M2JaCC9LH3CsbhHbYbhGE9cHXL8J5gvAzwKnANuBv+hQZlajn0la3Rq57UV297lZYWy4tiGTZ6vW2I64DvtJjO1x0FeCsf247b3l8LHXUNwyaNdztLW2OtfaXmp76QIW9bNZYdyM4Aqm7tiOuA4dxRVMd+UY5S2/DtzXodhLo62VIwOuoBgRLYQ0I3gYGrEdDohMHvL3bAcj6QaKsZuXlKObXQGcLukUij/DI8Bvl2WPoRj97Nxuo60NZS/CWNK+4d4niNgOozLs2G6KngnG9gUdZn+xS9ltwLmV7/uNthZCEjP0xmgR22EkDkBsN0W05A+NJMansdlsfeyNp416E8IQ5RTbg7ymHMJw1fSQv1ere0mrJD1Zab1/UeW3b0n6V0n/u+a9Czmr6SH/gLG9UtKWclpZmf8OST8q6/ysJJXzk3q5qIormNBcNZzlVVrdn0nx9tdGSZO2728reqPtiztU8efAwZTPYkKoxYhjW9JiimeOSylu2t1ZLruL4lX91cAPKG4DLwO+WS76adufSt3GuIIJzdS6T50yzWygVve2bweemeXWh9Bdamz3Nkhsnw1M2d5ZJpUpYFn5FuWrbX/ftoEvAe9N3bV2kWBCY2nfvqSph9RW9+eXrffXSzq+w+8h1KaGuIbBYrvbsseWn7vV2auXi1eIBBMaKvH5S3GrYUmrtXw5ra5UlNLq/mbgBNsnA7cB1w9nn0KA5NieOa5hsNjutuxMdab0cvEK8QwmNJOZzX3qHbaXdvmtZ6t7209Vvl4DXJW64hBmLT22Z4prGCy2pynagFWXvaOcf1ynOm0/3pop6Rqg54svcQUTmqueZzA9W923td5fDjxQw9aH0F09z2AGie0NwFmSjihvdZ0FbLC9HXhG0i+Wb4/9JvD1DnV16+XiFeIKJjRWHW0FurW6l7QG2GR7ErhE0nJgD0X3/ate2gbpe8BbgEPL1v4X2t4w8IaFrI06tm3vlHQlRZICWGO7NXTFRyjGSnoVxdtjrTfIPtmpl4uZRIIJzVVTY7ROre5tX175fBlwWZdlf7mWjQihqhmxvQ5Y12H+JuCtHeb/xmy3r5EJ5hl27bjN63/cNnsJsGMU21OzcdkPGGxf3jDjrzbsHa/+NMY8riH2pSW72O6mkQnG9mvb50na1OOB15wwLvsBB2Bfxqw7jXGOa4h9mZUxi+1uGplgQgCyOQhDhjKJ7UgwoZkMjMm45CG8QkaxPZcSzNpRb0BNxmU/YKj7YnAW96kjHpopYrsGcybB2B6L4B2X/YAh74vJ4kFoxEMzRWzXY84kmJChTO5ThwxlEtuRYEJzZXIQhgxlEtuN7yqm14A6c4mkR8qBfO6WtGnU2zMbZe+pT0i6rzJvsaSpcsCiqZTeVdPNqrPLOSliuxkaG9tjoNEJpjKgzjnAScAFkk4a7VYN7F22T5mD7QWuoxh4qOpS4HbbJwK3l9/rYWDfvrRpDorYbpTraGJsj4FGJxgGHCwq1Mf2dyn6Mqo6j5e7/76eAQYm6rLScT7Ti9huiMbG9hhoeoJJHVBnrjBwq6Q7O4ztMBcdVfa+Svnf19VXddmdRso0N0VsN9voY3sMNP0hf8qAOnPJO21vk/Q6YErSg+XZU2hn8Hi3FYjYztX4x/ZLmn4F03NAnbnEdmvgnieAmyhuk8xlj7fGiCj/+0Stte9z2jQ3RWw32+hjeww0PcH0HFBnrpB0iKTDWp8pBvjpOWBPw00CK8vPKykHJqrNeD+DidhuttHH9hho9C2ybgPqjHiz+nUUcFMxSBzzgb+1/a3RblI6STdQDLG6pBx46wrgE8BXJV0I/AR4f20rtMfmTZpOIrabI2J7eBqdYKDzgDpzke2HgbeNejv6ZfuCLj+dMcSVDq3qJojYboaI7eFpfIIJuTLeu3fUGxHCEOQT25FgQjNl1KV5yExGsd30h/whZ96XNvXQq0sWSaskPVl2c3K3pIsqv60suwvZImll+7Ih9KWGuIbhxLakd5Td/myV9FmVD9f66T4nEkxoJAPe56RpJrPokuXGspuTU2xfWy67mOKB72kUr91eUW+fVCFHqbHdyxBj+wvAauDEcmp1ozPr7nMiwYRmsuu6ghmkS5azgSnbO23vAqbYv8+qEGYnNbZ7qz22yzY/r7b9fdsGvsTL3eTMuvucSDChsbx3b9LUQ2qXLOdLulfSekmtBpDj1p1LaIga4hqGE9vHlp871Tnr7nPiIX9opGfYteE2r1+SWPygti7i11ZGJEzpkuVm4AbbuyV9mOLs7N2Jy4YwK7OI7ZniGoYT27XGfCSY0Ei267oV1bNLFttPVb5eA1xVWfb0tmXvqGm7QqYaHtvT5edOdT4u6Wjb21O7z4lbZGHc9eySpdXnVGk58ED5eQNwlqQjygegZ5XzQmiC2mO7vPX1jKRfLN8e+01e7iZn1t3nxBVMGGvdumSRtAbYZHsSuETScmAPxbggq8pld0q6kuJABlhju33ckBBGYoix/RGKQdheBXyznKCP7nPkTLosCCGEcGDFLbIQQghDEQkmhBDCUESCCSGEMBSRYEIIIQxFJJgQQghDEQkmhBDCUESCCSGEMBSRYEIIIQzF/wckIqabhj+HPAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 4 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from SWESimulators import LxF\n",
    "importlib.reload(LxF)\n",
    "\n",
    "nx = 10\n",
    "ny = 15\n",
    "num_ghost_cells = 1\n",
    "dt = 0.1\n",
    "g = 9.81\n",
    "\n",
    "h0, hu0, hv0, dx, dy, nx, ny = gen_test_data(nx, ny, num_ghost_cells)\n",
    "plt.figure()\n",
    "plt.subplot(121)\n",
    "plt.imshow(h0)\n",
    "plt.colorbar()\n",
    "\n",
    "with Timer(\"construct\") as t:\n",
    "    sim = LxF.LxF(h0, hu0, hv0, \\\n",
    "                    nx, ny, \\\n",
    "                    dx, dy, dt, \\\n",
    "                    g)\n",
    "\n",
    "with Timer(\"step\") as t:\n",
    "    t = sim.step(10.0)\n",
    "    \n",
    "with Timer(\"download\") as t:\n",
    "    h1, hu1, hv1 = sim.download()\n",
    "\n",
    "plt.subplot(122)\n",
    "plt.imshow(h1)\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=> construct 4916.617155 ms\n",
      "=> step 118.532658 ms\n",
      "=> download 0.992298 ms\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x94b88cda58>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZgAAADxCAYAAAD2t6FzAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3X/wXXV95/Hni19BUJQYQRAqbqV2qUKsWbDDtAMyhMDYaAfsQndq0sKkOMtQu7tjoW7BgboD1tYfiysTIAP0R4rNLhpWNHxBXe1WLIFGIPxoIoMSw4IQRNAaSL7v/eN8Ljm5uT/O995zc8+9n9dj5sz33nPPj8/NvE/e93w+n/P5KCIwMzOr2z7jLoCZmU0nJxgzMxsJJxgzMxsJJxgzMxsJJxgzMxsJJxgzMxsJJxgzMxsJJxgzMxsJJxgzMxuJ/cZdALNOzjj14Hh2285K2957//Z1EbFkxEUyq0XV2J6GuHaCsUZ6dttO/mndL1Tadt8jNi0YcXHMalM1tqchrp1grJECmGV23MUwq11Ose0EY40UBC9HtSoys0mSU2w7wVhj5fIrz/KTS2w7wVgjBcFOTyVhUyin2HaCscaaJY+L0PKTS2z7ORhrpAB2EpUWs0lSNbarkLRE0qOSNku6pMPnyyX9SNKGtFxQ+myZpE1pWVZa/1VJ35W0UdK1kvZN6+dLmknbz0g6tF/5nGCssWaJSovZpKkjrtN//J8DzgSOA86TdFyHTW+JiIVpuT7tOx+4HDgJOBG4vJQwfjsiTgDeDrwB+EBafwlwV0QcC9yV3vfkBGONFMDLEZUWs0lSNbYrOBHYHBGPRcRLwN8B76tYjDOAmYjYFhHPATPAEoCI+EnaZj/ggFRk0rFvSq9vAt7f7yROMNZIUbF6zFVkNmmqxjawQNL60rKi7VBvAp4ovd+S1rU7W9L9ktZIOrrKvpLWAU8DLwBr0urDI+JJgPT3sH7f1Y381kwBO507bBpVj+1nImJRj8/V+ei7uQ1YHRHbJV1Icefxnn77RsQZkg4E/iZtP1OpxG18B2ONVDztXG0xmyRVY7uCLcDRpfdHAVt3O1fEsxGxPb29DnjXHPb9ObCWXdVuT0k6AiD9fbpfAZ1grKHEzoqL2WSpLa7vAY6V9BZJBwDnUiSEXWdKCSFZCjycXq8DFks6NDXuLwbWSXp1KYnsB5wFPJL2WQu0epstA77Ur4CuIrNGKhpCnTxs+tQV2xGxQ9JFFMliX2BVRGyUdAWwPiLWAhdLWgrsALYBy9O+2yRdSZGkAK5I6w4H1kqal475NeDatM1VwBcknQ/8gF29y7pygrFGKp4VcIKx6VNnbEfE7cDtbesuK72+FLi0y76rgFVt654C/l2X7Z8FTptL+ZxgrLFmfQdjUyqX2HaCsUbyHYxNq5xi2wnGGikQO90HxaZQTrGdx7e0iTQbqrT0M+h4TZJOLa3bIOnnkvo+vWzWTx1xPQl8B2ONFIiXYt+hj1Mar+l0ir7/90haGxEPtW16S0RctFsZIr4OLEzHmQ9sBu4YulCWtbpiexI4wVgjFQ+j1XKD/cp4TQCSWuM1tSeYfs4BvhIRP6ujUJavGmO78fL4ljaR5vCgZa8xm4YZr6nsXGB1bV/OspbLA8SNvIM5QPPiQA4edzFshH7OT3kptne9iiLEzqj8+6fXmE3DjNdUHKB4svkdFA+0DcxxnYcXeO6ZiHhDt8/nGNsTbagEI2kJ8BmKJz6vj4ir2j6fB9xMMf7Ns8C/j4jH+x33QA7mJM3peR6bMN+Ju/puM1vPr7hK4zWV3l4HXN0W2w8Bt0bEy62NBoltx3Ue7ow13++3TU2x3XgDp9GKk92cDzwXEW8FPgVcPej5LC9FQ+h+lZY+Bh2vqRzbpwH/2HZcx7YNpGpsT4Nh7tOqTHZTnqBmDXCapDxStw2l1RBaZel5nIgdQGu8poeBL7TGa0pjNEExXtNGSd8FLgY+za6OAUdSjON0ZNuhHds2kKqxPQ2GSZOdGk9P6rZNGpjteeD1wDNDnNcysbOmZwHmOl6TpHPYFbePpwEFHdtWm7piu+mGSTBVGk+rbFNsWPT8WQFwIAcNUSybBmN+2rm22HZcW7ucnuQfJsH0bTwtbbMlzS3wWooho/cQESuBlQCHaL7nMjRmx9fTprbYdlxbJ2OM7b1qmG/Zt/GU3SeoOQf4WkT4IrO+igEB96m0jIBj20amamxPg4HvYCpOdnMD8FeSNlP8uju3jkLb9AvEy2MaTsOxbaM0ztje24bqC1eh8fTnVJj1zKxdBGN9GG2csb1u64ZRHNaGdMaRC2s5zrhje2+ajs7WNoWUzcNolpt8YtsJxhopyOdXnuUlp9h2grHGmpaGTrN2ucS2E4w1UjA9ky6ZleUU204w1kgBvDwl4zGZleUU23l8S5tA0zMnhtnu8oltJxhrpCCfp50tLznFthOMNVYuv/IsP7nEdh5p1CZOhJiNfSotZpOkamxXIWmJpEclbZZ0SYfPl0v6kaQNabmg9NkySZvSsiytO0jSlyU9kqawuKrKsbrxHYw1UtEQmsdwGpaXumK7NOnj6RSDr94jaW1EPNS26S0RcVHbvvOBy4FFqUj3SloLbAc+GRFfT+Pw3SXpzIj4Srdj9eIEYw2Vz7zllpvaYvuVSR8BJLUmfWxPMJ2cAcxExLa07wywJCJWA18HiIiXJN1HMZr4QHwFWyMVDaGqtJhNkqqxXUGnSR/f1GG7syXdL2mNpNY0FH33lfQ64DeBu/ocqysnGGusMQ7XbzZSFeN6gaT1pWVF22GqTHp3G3BMRBwP3Mmuab577pvmOFoNfLZ1h9TjWF25iswaKaennS0vc4jtZyJiUY/P+06MFxHPlt5eB1xd2veUtn2/UXq/EtgUEZ+ucKyu/PPPGmuWfSot/QzZ0+YXJN0h6WFJD0k6ptYvaVmqI66pMDGepCNKb5cCD6fX64DFkg6VdCiwOK1D0p9RzND64YrH6sp3MNZIEfDy7PC/f4bpaZPcDHw8ImYkvRqYHbpQlrW6YrvixHgXS1oK7KCYGG952nebpCspkhTAFWndUcBHgUeA+yQBXBMR13c7Vi9OMNZIRTXCeHvaSDoO2C8iZgAi4sU6CmR5qzG2q0yMdylwaZd9VwGr2tZtoXP7TM9jdTPwt5R0tKSvp6qDjZL+sMM2p0h6vlT1cFmnY5l1sjON2dRvoXdj6CA9bU6U9HWKX4YnSPqupH+W9OfpjsixbUOpGNcTb5g7mB3Af46I+yS9huJBnZkOVQ/fioj3DnEey1CrK2dFvRpDq/a0WR0R2yVdCPwlcDHwb4AbgEMo6pw/RlEtcEPaz7FtczbH2J5oA9/BRMSTEXFfev0CRYNPp1+GZgOobaiYSj1tImJ7ensd8I4U21uAfwbuB94IfBH41Rq+nGWtvqFimq6Wb5F61rwT+E6Hj38tVTF8RdKv9DjGilYVx8ts77aZZWQ2zV3eb+ljmJ429wCHAe+iiO33sHvbTd/YdlxbJzXE9UQYupE/9az5n8CHI+InbR/fB7w5Il6UdBbFL8BjOx0nIlZS9L3mEM1vr8KwzBQ9bYYfr2mYnjbAqyiukR3A/wXupbjDgYqx7bi2dnXF9iQYKsFI2p8iufxNRPyv9s/LCScibpf0PyQtiIhnhjmvTb86H7QcpKdNKbavjYi/7HBMx7YNJKeHiAdOMCo6SN8APNzpAkzbvBF4KiJC0okUVXLPdtrWrN24qgkc2zZq01IF1s8wdzAnA78LPCBpQ1r3J8AvAETEtcA5wIck7QD+FTg3IlxNYH2NuaeNY9tGJqdeZAMnmIj4B7o8kFPa5hrgmkHPYXkbV08ax7aN2rT0EuvHT/JbI0WIHZlchJaXnGLbCcYaK5dqBMtPLrHtBGONlFM9teUlp9h2grHGyuUitPzkEttOMNZIOT0rYHnJKbadYKyxcnlWwPKTS2w7wVgjRcCOGiZlMmuanGLbCcYaK5dqBMtPLrHtBGONlFM9teUlp9h2grHGikwuQstPLrHtBGONlUtDqOUnl9h2grFGisinntryklNsO8FYQ4mdmfS0sdzkE9tOMNZYudRTW35yiW0nmBHT/gfUerx4+aVaj9dUOY3XZHnJKbadYKyZoqirNps6GcX20BWBkh6X9ICkDZLWd/hckj4rabOk+yX96rDntDzMokpLP5KWSHo0xeAlHT5fLulHKYY3SLqgFNch6Wdp/drSPo5rG1gdcT0J6mppOjUiFkbEog6fnQkcm5YVwOdrOqdNsUgNoVWWXiTtC3yOIg6PA86TdFyHTW9JMbwwIq5P604FfhoRB6X1S0vbO65tIFVju4pBfjyVPlsmaVNalqV1B0n6sqRHJG2UdFVp+3mSbknn+o6kY/qVb290ZXgfcHMU7gZeJ+mIvXBem3AR1ZY+TgQ2R8RjEfES8HcUMTksx7UNrIa4HurHk6T5wOXASRTXyOWSDk3bfzIifhl4J3CypDPT+vOB5yLircCngKv7lbGOBBPAHZLulbSiw+dvAp4ovd+S1pn1FKFKC7BA0vrSUo7DqvF3dqrqWiPpaFJcAwen6rK7Jb1/gOOa7aFiXPczzI+nM4CZiNgWEc8BM8CSiPhZRHy9KGO8BNwHHJX2eR9wU3q9BjhNUs+C1tHIf3JEbJV0GDAj6ZGI+Gbp804F2CM/p/8UVgAcyEE1FMsmWfErrnI99DNdqmehWvzdBqyOiO2SLqS4iFpx/Q7gryl+sX1a0gMR8b2Kx3Vc2x7mENsL2tq1V0bEytL7Tj9yTupwnLMl/QbwL8AfRcQTXfbd7QeSpNcBvwl8pv18EbFD0vPA64Fnun2Boe9gImJr+vs0cCtFVi3bAhxden8UsLXDcVZGxKKIWLQ/84Ytlk2B2VClpY++8RcRz0bE9vT2OuBdpbh+gCKujwa+QVFtUOm4aX/Hte2hYlw/04qdtKxsO0zVH0/HRMTxwJ3sugPpua+k/YDVwGcj4rE5nG83QyUYSQdLek3rNbAYeLBts7XAB1Ovm3cDz0fEk8Oc1/JQUxvMPcCxkt4i6QDgXIqYfEVb28lS4FFJr5F0aKqXXgz8ADgZeCht57i2gdXRBsOAP54q7rsS2BQRn+50vpSAXgts61XAYavIDgduTdVw+wF/GxFfTdUMRMS1wO3AWcBm4GfA7w15TstAIGZrGE4j3cpfBKwD9gVWRcRGSVcA6yNiLXCxpKXADooL5k+AfwAOorigngP+FPg28BsUScZxbQOpK7Yp/XgCfkjx4+l3yhtIOqL0w2cp8HB6vQ74b6WG/cXApWmfP6NIHhewu7XAMorr4BzgaxG9U+FQCSbdOp3QYf21pdcB/MdhztNEmletumP2V3+51vPuc98jlbaL7dv7b9RwdT2LFhG3UySE8rrLSq8vJV1cJXvEddv+UxnXtnfUEdsD/nhanvbdJulKiiQFcEVadxTwUeAR4L5083BN6n12A/BXkjanY53br4x+kt+aaW6N/GaTo8bYHvDHU+uzVcCqtnVb6NzWQkT8HPjAXMrnBGPNlclwGpahTGLbCcYay3cwNq1yiW0nGGukAGZn87gILS85xbYTjDVTAJn8yrPMZBTbTjDWWLkMaW75ySW2nWCsuTK5CC1DmcS2E4w1VOUB/8wmTD6x7QRjzZXJrzzLUCax7QQzoH3feFil7Z7845dqPe9Rf1jtvDu+/0T/jZosIDLpaWOZySi2nWCswfK4CC1HecS2E4w1VybVCJahTGLbCcaaK5OL0DKUSWw7wVgzZfQwmmUmo9h2grHGyuVhNMtPLrHtBGPNlUlPG8tQJrE98LRqkt4maUNp+YmkD7dtc4qk50vbXNbteGbtFNWW2s/r2LYRG0dcj8PAdzAR8SiwEEDSvhRTdt7aYdNvRcR7Bz2PZSoYW0OoY9tGaoyxvbfVMjE0cBrwvYj4fk3Hs+ypaAitsvQ7krRE0qOSNku6pMPnyyX9qHQ3Up6L/DTgceAfJV1T29ezjNUT15OgrjaYc4HVXT77NUnfBbYC/yUiNtZ0zrGafe3Blbb77ond/lkGc+Zr+06DPT1q+JWX7kA+B5wObAHukbQ2Ih5q2/SWiLiowyHOBf4V+D8dPpvK2La9wHcw1Ug6AFgK/H2Hj+8D3hwRJwD/Hfhij+OskLRe0vqX2T5ssWwazFZcejsR2BwRj0XES8DfAe+rcvoU278F/Bi4o+3jSrHtuLaOho/riVBHFdmZwH0R8VT7BxHxk4h4Mb2+Hdhf0oJOB4mIlRGxKCIW7c+8GoplE631rMDwVQlvAsoDs21J69qdLel+SWskHZ3WnZVKcvEexasY245r20PV2J4CdSSY8+hSPSbpjZKUXp+YzvdsDee0DMyhF9mC1l1CWlaUD9Ph0O0VFLcBx0TE8cCdwE1p/Z8Cd0bEHiOHOrZtGO5FVoGkgyjqtv+gtO5CgIi4FjgH+JCkHRT12OdG5PKIkQ2teqQ8ExGLuny2BTi69P4oijaTXaeJKCeG64CrU2y/HThM0uPAq4GDJb01Ipbg2LZhZBIpQyWYiPgZ8Pq2ddeWXl8DuOeNjdM9wLGS3kLR3fhc4HfKG0g6IiKeTG+XAg+n2J5X2mY5sKjVEcCxbdafn+S3xqqjmiAidki6CFgH7AusioiNkq4A1kfEWuBiSUuBHcA2YPnwZzbrblqqwPqp6zkYs3oFxXAaVZZ+h4q4PSJ+KSJ+MSI+ntZdlpILEXFpRPxKRJwQEadGxCMdjnFjl27MZnNTNbYrGOYZL0nLJG1Ky7LS+o9LekLSi1WP1Y3vYKy5MvmVZxka8zNekuYDlwOLUmnuTfs+R9Hp5RpgU4fTdnterCPfwVhjjWssMrNRqymuB37GCzgDmImIbSmpzABLACLi7lKb5FCcYKy5ouJiNmnqiethnvGqum+VY3XlKrIB7fP8Tyttd8I/nVfreY96/sX+GzElDwI7edi0qhbbCyStL71fGRErS++rPuO1OiK2p0dIbgLeU3Hfdt2O1ZUTjDWSq79sWs0htns93wUDPuNV2veUtn2/0aswPY7VlavIrLlq6kVm1jj1xPUrz3ilcfPOBdaWN5B0ROntUuDh9HodsFjSoZIOBRandV31OFZXvoOxxvIdjE2rcT/jFRHbJF1JkaQAroiIbQCSPkHxMPJBkrYA10fEx7odqxcnGGsuJxibVjXFdhpo9fa2dZeVXl8KXNpl31XAqg7rPwJ8pMP6rsfqxgnGmsltMDatMoptJxhrrkwuQstQJrHtBGONpanoa222p1xi273IzMxsJHwHY82VSTWCZSiT2HaCGdDO//d0pe0Ov3r+WM478TJqCLXMZBTbTjDWXJlchJahTGK7UhuMpFWSnpb0YGndfEkzaS6BmfQ0aKd9O845YNbXiAe7dFzb2GQyiGvVRv4bSUM5l1wC3BURxwJ3pfe7Kc05cBLF0NKXd7tgzcpE0dOmyjKEG3Fc215WNbanQaUEExHfpBgaoOx9FKNpkv6+v8OuXeccMOup4lwww9RlO65tLDKa52iYNpjDW5PSRMSTkg7rsE3lOQckrQBWABzIQUMUy6bGeC4yx7WN3pQkkH5G/RxM5TkHImJlRCyKiEX7M2/ExbKJUFMbTL95y4HXteYaB17Tmmtc0pvT+w2SNqY5MMBxbcNyG0xfT7WGb05/O/Wf7TtfgVk3dVSRleYtPxM4DjhP0nEdNr0lIhYCm4Evp3UBPJbWnwRcIulIHNc2pFyqyIZJMGuBVu+ZZcCXOmwz5zkHzF5Rzx3MXOctL8f17wBfTK/nset6cVzbcHwHs4uk1cC3gbdJ2iLpfOAq4HRJm4DT03skLZJ0PRRzDgCtOQfuoTTngFlPUVsvsq7tJaW4fiPwIUlPAG8D3luK65sl3Q/8EHgiIrY6rm0oFWN7GlRq5I+IbhPLn9Zh2/XABaX3HeccmHSxfXul7fZZ33fSt7md9+WXaj1eo1X/Fddr7vKu7SWtuJb0euDF0lzjv526Kbccn6rGvijp8Ih4alrj2vaSKblD6cdP8ltjzaEeutfc5cPMW17eZqukjcCvA2sql8ysg2lpY+nHoylbc9XTBjPwvOWSjpL0qvT6UOBk4NGhvpMZZNMG4zsYa6aaLrJh5i0H/i3wF5KCoqrtkxHxwPClsqxNUQLpxwnGGknUV40w6LzlETEDHF9PKcwKdcZ20znBWGPlchFafnKJbScYa65MLkLLUCax7QRjzZXJRWgZyiS2nWCsmaZouAyz3WQU204w1lyZXISWoUxi2wlmxLJ68r5m0zJchlm7XGLbD1paY416wjGzcakrrvtNRSFpeWsqirRcUPqs47Tfkj4u6QlJL7Yda56kW9K5viPpmH7lc4KxZqr6FL8TjE2amuJ6rlNRpOX6tG+vab9vS+vanQ88FxFvBT5FhyGV2jnBWHM5wdi0qieu5zoVRVnXab8j4u7WrK5tytOJrwFOk9RpMNlXOMFYI7WednYVmU2bqrFdQdWpu8+WdL+kNZJaA79Wnva70/kiYgfwPPD6Xjs4wVhjaTYqLWaTpmJcL5C0vrSsaD9Mh0O3XxC3AcdExPHAney6A6k87fccz7cb9yKzZnL1l02r6rHdaxoKGG4qii3AKW37fqNPeVrn2yJpP+C1FIPDdtX3DkbSKklPS3qwtO7PJT2SbrtulfS6Lvs+LumB1HthfadtzLoZdRWZY9vGpaa4HngqCgab9rs8nfg5wNciomdJq1SR3Uhq/CmZAd6ebrv+hQ4j0Zacmnov9MrEZnsafSP/jTi2bRxqiOvUDtKaiuJh4AutqSjS9BNQTEWxUdJ3gYtJU1H0mvZb0ickbQEOkrRF0sfSsW4AXi9pM/CfgD26RbfrW0UWEd9s7+8cEXeU3t5Nkc3MajXqBnzHto3LuKeiSJ91nPY7Ij4CfKTD+p8DH5hL+epo5P994CtdPgvgDkn3dmigMutt/N2UHds2Gpl0vx+qkV/SRylmAfybLpucnOYyPwyYkfRIRHyzy7FWACsADuSgYYpl0yDGO5xGXbHtuLY9jDm296aB72DS0ALvBf5Dt4aeiNia/j4N3Ernp0Nb266MiEURsWh/5g1aLJsS43wOps7Ydlxbuxqfg2m8gRKMpCXAHwNLI+JnXbY5WNJrWq8peik82Glbs44iqi01cmzbXrGX43pc+laRSVpN0V96QepZcDlFo9E8iqoBgLsj4kJJRwLXR8RZwOHArenz/YC/jYivjuRb2FSq61dcShqfAfaliM+r0vpWbL9B0k7ghxRx+1OK2H4VRV//Zyl+eG6PiHfi2LYhTcsdSj9VepGd12H1DV223QqclV4/BpwwVOksXzU1dJYGBDyd4kGxeyStjYiHWrEtaTmwKCIuatv3l4CIiE3px9O9kl7n2LahTFEjfj9+kt8aq6aG0FcGBASQ1BoQ8KF+O0bEv5Reb5X0NPAG4Me1lMyy5UZ+szHTbLWlj2EGBNxVFulE4ADgewN+HbNX1BDXE8EJxpopmEsjf69BAYcZELA4QDHcxl8BvxcRU3Lp29hUje0p4Coya6w5NIT2GhRwmAEBkXQI8GXgv0bE3ZVLZNZDLo38voOx5qrnSf6BBwRM298K3BwRfz/ktzHbxU/ym41P62G0YUXEDkmtAQH3BVa1BgQE1kfEWooBAZdSPLm/jTQgIPDbwG9QDPDXWrc8IjYMXzLLVV2xPQmcYKyZor7JxAYdEDAi/hr461oKYdZSY2w3nROMNVce16DlKJPYdoKxxsqlGsHyk0tsO8FYMwWQSTWCZSaj2HaCsebK4xq0HGUS204w1li5VCNYfnKJbScYa6xcetpYfnKJbScYa6YpetjMbDcZxbYTjDVS8TBaJlehZSWn2HaCsebysJI2rTKJ7b5jkUlaJelpSQ+W1n1M0g8lbUjLWV32XSLpUUmbJV1SZ8Ft+imi0jLw8R3bNiajjOsmqTLY5Y3Akg7rPxURC9Nye/uHpZkEzwSOA86TdNwwhbWMVB3ocrjr8EYc27a3jT6uG6NvgomIb1IMADhXr8wkGBEvAa2ZBM0qKMZrqrIMfAbHto3FaOO6SYYZrv+iNAPgKkmHdvi86kyCZp1Vn3Csbo5tG61MJhwbNMF8HvhFYCHwJPAXHbapMpPgro2lFa0ZCV9m+4DFsqkRtU2ZPFe1xrbj2vZQMbar6NcWKGm5pB+V2hQvKH22TNKmtCwrrX+XpAfSMT8rSWl9pfbJsoF6kUXEU6XCXAf87w6b9Z1JsO2YK4GVAIdo/nSkbxvOGH7F1R3bjmvrqIbYLrUFnk4Rk/dIWhsRD7VtektEXNS273zgcmARxY+je9O+z1H8yFoB3E0xzcUS4Ctp109FxCerlnGgO5i2GQB/C3iww2Z9ZxI062kMjaGObdsr6onrYdoCzwBmImJbSiozwJIU/4dExLcjIoCbgfdX/l5tqnRTXg18G3ibpC2Szgc+kW6h7gdOBf4obXukpNuhmEkQaM0k+DDwhYjYOGhBLT+ana20DHx8x7aNScW4XtCqXk3LirbDVG0LPDu1Ka6R1Lrz7rbvm9Lrbsfs1z65m75VZBFxXofVN3TZditwVun9HjMJmlUSjPxhNMe2jUX12H4mIhb1+LxKW+BtwOqI2C7pQuAm4D099u11zM8DV6b3V1K0T/5+j/L5SX5rJjE9D5vN1RlHLhx3EWyEaoztvm2BEfFs6e11wNWlfU9p2/cbaf1RnY5ZsX1yN8N0UzYbrZq6KQ/Z0+arkn4sqe/FZFZZPd2U+7YFtrUpLqWo0oWienexpENTVddiYF1EPAm8IOndqffYB4EvdThWt/bJ3fgOxpprzD1tkj8HDgL+YOjCmLXUENsRsUNSqy1wX2BVRGyUdAWwPiLWAhdLWgrsoHioeHnad5ukKymSFMAVEdF66PhDFKNcvIqi91irB9knJC2kqCJ7nArXhBOMNVN9bTCv9LQBkNTqadOeYDoXI+IuSafUUhIzqLV9sVNbYERcVnp9KXBpl31XAas6rF8PvL3D+t+da/lcRWaNVVMvsmF62piNxCh7RzaJE4w1VMX2l6KqoVd3zqo9bY6JiOOBOyl62piNyNiGQNrrXEVmzRTM5SLr1Z1zmJ42ZvWbW2xPNN/BWHPQXanpAAAEsUlEQVTNVlx6G6anjdloDB/XE8F3MNZYdTwrMExPGwBJ3wJ+GXi1pC3A+RGxbuiCWdZyecbLCcaaq6aLcMieNr9eSyHMypxgxucFnnvmzljz/bbVC4BnxlGemk3L94Dhvsube34aATunpJ4gmfK4Bn+Xluxiu5tGJpiIeEP7Oknr+4zLMxGm5XvAXvguU/Yrb5rjGvxd5mTKYrubRiYYMyCbi9AylElsO8FYMwUwJfOSm+0mo9iepASzctwFqMm0fA8Y6XcJiCzqqR0PzeTYrsHEJJg09ezEm5bvASP+LkEWDaGOh2ZybNdjYhKMZSiTemrLUCax7QRjzZXJRWgZyiS2Gz9UTL/JoiaJpMfTfO8bJK0fd3nmIs3B/bSkB0vr5kuakbQp/e07R3d1cxrsciI5tpuhsbE9BRqdYEqTRZ0JHAecJ+m48ZZqaKdGxMIJfF7gRmBJ27pLgLsi4ljgrvS+HgHMzlZbJpBju1FupImxPQUanWAoTRYVES8BrcmibC+LiG9SjNNV9j52DW1/E/D+mk86zb/0HNsN0djYngJNTzBVJ4uaFAHcIenetjlLJtXhaQ5v0t/D6jt0Gk6jyjKZHNvNNv7YngJNb+SvMlnUJDk5IrZKOgyYkfRI+vVk7QJiup8VcGznavpj+xVNv4PpO1nUJImIrenv08CtFNUkk+yp1lwq6e/TtR59Nqotk8mx3Wzjj+0p0PQE03eyqEkh6WBJr2m9BhYDD/beq/HWAsvS62XAl2o9+nS3wTi2m238sT0FGl1F1m2yqDEXa1CHA7dKguLf/W8j4qvjLVJ1klYDpwAL0sRblwNXAV+QdD7wA+ADtZ0wYmp60nTi2G4Ox/boNDrBQOfJoiZRRDwGnDDucgwqIs7r8tFpIzzpyA7dBI7tZnBsj07jE4zlKoidO8ddCLMRyCe2nWCsmTIa0twyk1FsN72R33IWs9WWPvoNySJpuaQfpWFONki6oPTZsjRcyCZJy9r3NRtIDXENo4ltSe9Kw/5slvRZpca1QYbPcYKxRgogZqPS0sschmS5JQ1zsjAirk/7zqdo8D2Jotvt5fWOSWU5qhrb/Ywwtj8PrACOTUtrGJ05D5/jBGPNFFHXHcwwQ7KcAcxExLaIeA6YYc8xq8zmpmps91d7bKdnfg6JiG9HRAA3s2uYnDkPn+MEY40VO3dWWvqoOiTL2ZLul7RGUusByGkbzsUaooa4htHE9pvS607HnPPwOW7kt0Z6gefW3RlrFlTc/MC2IeJXlmYkrDIky23A6ojYLulCil9n76m4r9mczCG2e8U1jCa2a415JxhrpIioqyqq75AsEfFs6e11wNWlfU9p2/cbNZXLMtXw2N6SXnc65lOSjoiIJ6sOn+MqMpt2fYdkaY05lSwFHk6v1wGLJR2aGkAXp3VmTVB7bKeqrxckvTv1Hvsgu4bJmfPwOb6DsanWbUgWSVcA6yNiLXCxpKXADop5QZanfbdJupLiQga4IiLa5w0xG4sRxvaHKCZhexXwlbTAAMPnKDIZssDMzPYuV5GZmdlIOMGYmdlIOMGYmdlIOMGYmdlIOMGYmdlIOMGYmdlIOMGYmdlIOMGYmdlI/H9uO/ABFB/hSQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 4 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from SWESimulators import FORCE\n",
    "importlib.reload(FORCE)\n",
    "\n",
    "nx = 10\n",
    "ny = 15\n",
    "num_ghost_cells = 1\n",
    "dt = 0.01\n",
    "g = 9.81\n",
    "\n",
    "h0, hu0, hv0, dx, dy, nx, ny = gen_test_data(nx, ny, num_ghost_cells)\n",
    "plt.figure()\n",
    "plt.subplot(121)\n",
    "plt.imshow(h0)\n",
    "plt.colorbar()\n",
    "\n",
    "with Timer(\"construct\") as t:\n",
    "    sim = FORCE.FORCE(h0, hu0, hv0, \\\n",
    "                    nx, ny, \\\n",
    "                    dx, dy, dt, \\\n",
    "                    g)\n",
    "\n",
    "with Timer(\"step\") as t:\n",
    "    t = sim.step(10.0)\n",
    "    \n",
    "with Timer(\"download\") as t:\n",
    "    h1, hu1, hv1 = sim.download()\n",
    "\n",
    "plt.subplot(122)\n",
    "plt.imshow(h1)\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=> construct 4879.117727 ms\n",
      "=> step 109.936714 ms\n",
      "=> download 2.000093 ms\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x94ba01c3c8>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAADxCAYAAADhlTG6AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3XuQXOV55/HvT5fRIAFGIIMxYOMEYq/KFy4qSEKFAlOAoLzCWUhWeCuRErOKXcuSeJMqQy6QgDeFYzu2s/aaCKzFjmOM4yy2SAhCELM4duwgCPebBCWCEEFIwpiLkTQzz/5xTg9HPT3d73SfUZ8+8/tUnaL79HtOvz08rafPeW+KCMzMzLo1q98VMDOzweZEYmZmPXEiMTOznjiRmJlZT5xIzMysJ04kZmbWEycSMzPriROJmZn1xInEzMx6MqffFTBr5ezTF8SOnaNJZe95YNe6iFg6zVUyK0VqbA9SXDuRWCXt2DnKv6x7W1LZ2YdvXDTN1TErTWpsD1JcO5FYJQUwxli/q2FWujrGthOJVVIQ7Im0W1tmg6SOse1EYpVVt19tZg11i2332rJKCoLRSNvMBklqbKeQtFTS45I2Sbq0xesrJb0g6b58u6jw2q2Sfizp75qOeYekH0naKOlGSUOd6uFEYpU1RiRtZoOmjLiWNBv4InAOsBi4UNLiFkVvjIjj8u26wv5PAb/Wovwngc9GxLHAi8CHO9XFicQqKYBRImkzGySpsZ3gJGBTRDwVEbuBbwDnJdcj4g7g5eI+SQLeD3wr3/UV4IOdzuVEYpXlKxKrq5Li+gjgmcLzLfm+ZudLekDStyQd1eGchwA/joiRDufcixvbrZIC2OP2D6uhKcT2IkkbCs9XR8TqwnNNcvqim4EbImKXpI+QXWG8v817ppxzAicSq6TwbSurqSnE9vaIWNLm9S1A8QrjSGDrXu8VsaPw9Fqy9o+27wkcJGlOflUy4Zyt+NaWVVPAaOJmNlDKi+u7gWPzXlZDwHJgbbGApMMLT5cBj7atWkQA3wUuyHetAL7TqSK+IrFKykb/mtVPWbEdESOSLgbWAbOBNRHxsKQrgQ0RsRa4RNIyYATYCaxsHC/pe8C7gP0lbQE+HBHrgI8D35D0CeBfgS93qosTiVWUGG15u9Zs0JUX2xFxC3BL077LC48vAy6b5NhfmmT/U2Q9wpI5kVglZQ2STiRWP3WMbScSq6Ssr329vmxmUM/YdiKxyhqr2a82s4a6xbYTiVVSHX+1mUE9Y9uJxCopEKPunW41VMfYrtensVoZCyVtnXQ7Q6qk0wv77pP0uqSO8w6ZdVJGXFeJr0iskgKxO2b3fJ7CDKlnko0EvlvS2oh4pKnojRFx8V51iPgucFx+noOBTcBtPVfKZrSyYrtKnEiskrJBW6VcMI/PkAogqTFDanMi6eQC4B8i4rUyKmUzV4mxXRn1+jRWK6P5wK1OG/nkdoVtVeE0Zc2Quhy4obQPZzNaYlwPjEpekQxpXgyzoN/VsGn0Oq+yO3ZN+m2JEKOR/Dun3eR2Pc+Qms9X9B6yqSi6NndoQQzvt7CXU9gAeOUnz26PiDdP9voUY3sg9JRIJC0FPk82z8t1EXF10+vzgK8CJwI7gP8cEZs7nXeYBZysM3qpmlXcj+KOjmXGyvlV1tUMqU2x/QhwU0TsaRTqJraH91vICb/437v/JDYQ7rr10qc7lSkptiuj67SYuMzjh4EXI+IY4LN0nsLYDGg0SM5J2jrodobUYmyfAfyg6byObetKamwPkl6ur1KWeTyP7DYBZEs3npEv5WjWVqNBMmVre55sTYXGDKmPAt9szJCaz4oK2QypD0u6H7gE+BxvNNC/lWzm1Lc2ndqxbV1Jje1B0kvaa9WIefJkZfIpj18iW8pxew/vazPEaEl96ac6Q6qkC3gjbjfnU3U7tq00ZcV2VfSSSFIaMZOXbcx72qwCGGZ+D9WyOujz6N/SYrsY1/OGD+q9ZjbwPLJ9bx0bMYtlJM0B3kS2uMoEEbE6IpZExJK5zOuhWlYXYzEraZsGpcX2XnE95J6IlulTXE+bXmrbsREzf74if3wB8I/5Uo5mbWUT281K2qaBY9umTWpsD5Kub20lLvP4ZeCvJG0i+7W2vIxKW/0FYk+fppFwbNt06mdsT5ee+pglNGK+DvxKL+9hM1MEfR201c/YjlmJDbGJf57RofS/48hw2nuP7JdWbmxu2vtqJK3cnNfTLvrm/jT94nDWnrSys3aXsdJ6/2N7OgxWZ2WbQVS7QVtmmfrFthOJVVJQv19tZlDP2HYiscoatAZHs1R1i20nEqukYPAW9zFLUcfYdiKxSgpgz4DNN2SWoo6xXa/rK6uRtLVIBm3dBrPU2E46U5fLSOevrZC0Md9WFPZfKOnBfH2eWyUt6lSPeqVFq42AgRvda5airNjuZRnpfOnoK4AleZXukbQWeJls+YTFEbFd0p+RTXr6x+3q4m+qVZavSKyuSorrlBnYJ3M2sD4idkbEi8B6YCnZHHICFuSzWR/IxOmBJnAisUqKUD/n2jKbNqmxnaCXZaRbHpsv3vZR4EGyBLKYbBaHtnxryyopa5Cs1zQSyRIvslJHwE/lzzg2N3Fke+II+NHhtPfVaFq51L+NpjAIPfUfwVm708/ZzhRie5GkDYXnqyNideF5L8tItzxW0lyyRHI88BTwv8iWWPhEu4o6kVhF1W9da7NMcmxvj4glbV7vahnpwrGnNR17J3BcftyTAJK+CUxoxG/mb6pVUtYgqaTNbJCkxnaCbpeRhmxC0rMkLZS0EDgr3/cssFjSm/NyZxaOmZSvSKyy6jb616yhjNhOnKX6knxJ6RGyWapX5sfulHQVWTICuDIidgJI+hPgLkl7gKcbx7TjRGKVVMfRv2ZQbmxPdRnppnJrgDUt9l8DXDOVevgnn1XWGLOStk56HLT1Nkm3SXpU0iOSji71Q9qMVEZcV4mvSKySImDPWH8HbeW+CvzPiFgvaX+gnEUpbMYqK7arxInEKim7/C/lyzY+aAtAUmPQVnMimUDSYmBORKwHiIhXyqiQzWwlxnZldP1pJB0l6bv5Jf/Dkn67RZnTJL1UuGVweatzmbUyhZHtiyRtKGyrCqfpZtDWSZK+S9aI+T5J90v6V0mfyq9wHNvWk7rN2NDLFckI8LsRca+kA8jmalnf4pbB9yLiAz28j81AjS6Sidr1t+9m0NafA5cAP0M2qvdAsq6Tf0zWg6Ux0texbVM2xdgeCF1fkUTEcxFxb/74ZbK+xq1+6Zl1obQpUpIGbUXErvzptcB78tjeAvwr8ADwFuDbwAklfDib0UqbIqUySmkjyXuyHA/8qMXLvyDpfrIv7+9FxMOTnGMVsApgmPllVMsGXEnrWo8P2iIbbLUc+FCxgKTDI+K5/Glx0NbdwKHAMWSx/XmgOGVFx9guxvW84YOSKhyz0z536nQmo0Ppf8fRobRyI4lf0T0HNl/8tZb6A300cWqWqfw7nFp2zqvp5+zEa7Y3yXuy/C3wOxHxk6aX7wXeHhGvSDqX7Bfdsa3Ok88hsxrgQB2cFn1WW1nPlt7n2upl0BawH9l3ZAT4PnAP2RULJMZ2Ma4PeNORjmsrLbarpKdEkk/w9bfAX0fE/21+vZhYIuIWSf9b0qKI2N7L+1r99XvQViG2r4mIP29xTse2daWOg227TiT5XPVfBh5t9UXLy7wFeD4iQtJJZG0yO1qVNWvWr8t/x7ZNN9/aesMpwK8BD0q6L9/3+8DbYHyY/QXARyWNAD8FlkeEL++toz73bHFs27SpY6+trhNJRPwTHVYHiIgvAF/o9j1sZutXzxXHtk23QeuV1YlHtlslRYiRmn3ZzKCese1EYpVVt8t/s4a6xbYTiVVSHe8jm0E9Y9uJxCqrbl82s4a6xbYTiVVSHfvap0oesZ5YbmxO+t8x+dZ9YrmR/dLKjR4wmna+BWlvHLPS2yBS/z77lTRCqI6x7URilVW3vvZmDXWLbScSq6QIGKnZ4j9mUM/YdiKxyqrb5b9ZQ91i24nEKqmO95HNoJ6x7URilRU1+7KZNdQttp1IrLLq1iBp1lC32K5Xi4/VRkR2HzllMxskqbGdQtJSSY9L2iTp0havr5T0gqT78u2iwmsrJG3MtxWF/UOSVkt6QtJjks7vVA9fkVhFidGa9Wwxy5QT25JmA18EziRbFvpuSWsj4pGmojdGxMVNxx4MXAEsIRtsf09+7IvAHwDbIuLnJM0CDu5UFycSq6y63Uc2aygptk8CNkXEUwCSvgGcBzQnklbOBtZHxM782PXAUuAG4DeBd2X1jDGg41BMJ5JpprmJi2Anij27Sz1fVdVxPqJU6SPW0843lT/j3FfTllRR2kB0dh+QVm7ssJG09z0grX6vz5qX9sZAzEm8Ongy+ZTt34/SYvsI4JnC8y3AyS3KnS/pVOAJ4GMR8cwkxx4h6aD8+VWSTiP71BdHxPPtKuJ7B1ZNkd1LTtnMBkp6XC+StKGwrWo6U6ts1PyNuBk4OiLeC9wOfKXDsXOAI4HvR8QJwD8Dn+70kXpOJJI2S3owb8jZ0OJ1SfqLvDHoAUkn9PqeNjOMoaStk24aJAtxHZJey/evLRzjuLauJcb19ohYUthWN51mC3BU4fmRwNZigYjYERG78qfXAid2OHYH8BpwU77/b4COsV3Wra3TI2Ky+2jnAMfm28nAl2h9+WU2LvrcICnpD4HTgc0RsX+LUzuurStlxTZwN3CspHcAzwLLgQ8VC0g6PCKey58uAx7NH68D/lTSwvz5WcBlERGSbgZOA/4ROIOENpd90UZyHvDVfD3rH0o6qOnDmbVU0m2rXhok23FcW9fKiO2IGJF0MVlSmA2siYiHJV0JbIiItcAlkpYBI8BOYGV+7E5JV5ElI4ArGw3vwMeBv5L0OeAF4Dc61aWMRBLAbZIC+MsWl18tG3UAf+GsrSn0bFnUdFt1dSEOu2qQJI9rYIGkzcC/A1dHxLfbnNdxbUnK6pEYEbcAtzTtu7zw+DLgskmOXQOsabH/aeDUqdSjjERySkRslXQosF7SYxFxV+H1lAYh8oakVQDDzC+hWjbIsgbH5C/b9ohYMslrqQ2SN0TELkkfIWuQbMT1e4CvAZ8FPifpwYh4MvG8e8X1vOGDJhxgM88UY3sg9HyjLiK25v/dRtZAc1JTkY4NQvnxqxuNSnNJ77pn9VXSyPauGiQLcf0gWVwfBdwJHJ963vz4N+J6aEHKx7YZoG4zNvSUSCQtkHRA4zFZg81DTcXWAr+e93L5eeAl30e2FCV1/x1vkJQ0RNYgubZYQNLhhafLgMclHSBpYd4YeRbwb8ApvNG24ri2rtWtW3uvt7YOA26S1DjX1yPi1vz2ABFxDdn9u3OBTWTdyjo23JgFYqyEni1dNkj+PvBPwHyyq44XgT8i61N/KlkycVxbV8qK7SrpKZHkPWHe12L/NYXHAfy3Xt6nijQv7fbb2AnvKvV9Z937WFK52LWrc6GKK+tHWZcNkhPiuun4aYvr0aG02xqp66vP3pP+l1zwXNrMCXNe2ZNUTi17Tk+0/Yi0D3PwIS8lldup9M/8uoaTy5ZlwC44OvIUKVZNNWyQNANqGdtOJFZddfvZZtZQs9h2IrHKqtuvNrOGusW2E4lVUgBjY/X6splBPWPbicSqKZja/Odmg6KGse1EYpU1aH3pzVLVLbadSKy6avZlMxtXs9h2IrGKUu0aJM0y9YttJxKrrpr9ajMbV7PYdiLp0uy3HJpU7rmPl7vG+pG/nfa+I08/07lQlQVEzXq2pBodSiyY+OeZ83r6ew89mzZyfPSJtAXMD5p3XFK5HSemzRTxlv1fTio3NDtxUXngudHUOJudfM62ahjbTiRWYfX6spm9oV6x7URi1VWzy3+zcTWLbScSq66afdnMxtUstp1IrJpqOGjLDKhlbDuRWGXVbdCWWUPdYtuJxKqrZj1bzMbVLLa7XqZL0jsl3VfYfiLpd5rKnCbppUKZyyc7n1kzRdpW+vs6tm2a9SOup1PXVyQR8ThwHICk2cCzwE0tin4vIj7Q7fvYDBX0rUHSsW3Tqo+xPV3KWjj4DODJiHi6pPPZjKesQTJl63QmaamkxyVtknRpi9dXSnqhcHVxUeHlM4DNwA8kfaG0j2czWDlxXSVltZEsB26Y5LVfkHQ/sBX4vYh4uKT37KuxNy1IKnf/SZP9WbpzzpuWl3q+SivhV1t+RfFF4ExgC3C3pLUR8UhT0Rsj4uIWp1gO/BT4fy1em5bYHp2XuGZ74kBrpQ/yJobmphdOMTvxs8wfSSr3M/tvTyo3f86BSeUAXvpp6prt85PP2VFJVySSlgKfJxt2f11EXN30+krgU2RX1QBfiIjr8tdWAH+Y7/9ERHyl6di1wM9ExLs71aPnRCJpCFgGXNbi5XuBt0fEK5LOBb4NHDvJeVYBqwCGy/wfZoNrrJSznARsioinACR9AzgPaE4kE+Sx/cvAncBtwJLCy0mxXYzrecMH9fRBrEZKiO1efiRJOhi4giymA7gnP/bF/PX/BLySWpcybm2dA9wbEc83vxARP4mIV/LHtwBzJS1qdZKIWB0RSyJiyVzS5t2xGmv0te/9FsARQHHisS35vmbnS3pA0rckHZXvOzevySUTqpcY23vF9VDaVazVXGpsdzb+IykidgONH0kpzgbWR8TOPHmsB5YCSNof+B/AJ1I/UhmJ5EImua0l6S2SlD8+KX+/HSW8p80AU+i1tUjShsK2qniaFqduvrFwM3B0RLwXuB1oXOL/EXB7REyYAdOxbb0oIa6htx9J7Y69CvgM8Frq5+np1pak+WSXVb9V2PcRgIi4BrgA+KikEbL7zMsj6jYUx6ZNeqRsj4glk7y2BTiq8PxIsjaNN94mopgArgU+mcf2u4FDJW0G9gcWSDomIpbi2LZepEVKu7iG9B9JN0TErvzf5q8A75/sWEnHAcdExMckHZ1US3pMJBHxGnBI075rCo+/ALini/XT3cCxkt5B1uC4HPhQsYCkwyPiufzpMuDRPLbnFcqsBJY07jU7tq0CuvqRVDj2tKZj7wR+ATgx//E0h+yH1J0RUSw7gUe2W2WVMSgrIkYkXQysI+vZsiYiHpZ0JbAhItYCl0haBowAO4GVvb+z2eRKGnDY1Y+k/PE64E8lLcyfnwVcFhE7gS/lxx4N/F2nJAJOJFZVQWnTSOSN4bc07bu88PgyWvc6LJa/Hri+lArZzFZSbPfyIykidkq6iiwZAVyZJ5GuOJFYdbnFweqqpNju5UdSRKwB1rQ592aydsKOnEissgZtviGzVHWLbScSq66afdnMxtUstp1IujTrpVeTyr3vXy4s9X2PfCltsGk5g8L7rGZftlRjQ2nlRvZLK7dn//T78TtPWNi5EDB89ElJ5ba/N+2fmON/dmNSuRMWbE4q9/3RlhNotDSrH5cHNYttJxKrpEGcStssRR1j24nEqqtmi/+YjatZbDuRWGXV7VebWUPdYtuJxKqrZl82s3E1i20nEqumGt5HNgNqGdtOJFZdNfuymY2rWWw7kVhlqRZ9mM0mqltsl7Vmu5mZzVC+IrHqqtnlv9m4msW2E0mXRv99W1K5wz55cF/ed+DVsEEyVSTeJxgdTiu3+5DR5Pd+7e1pf3TNH0kq94vHPJFU7pLD1yeV2xOzk8rdNpI01yAAr72WtrR3ad/kGsa2E4lVV82+bGbjahbbSb99JK2RtE3SQ4V9B0taL2lj/t+Wk/RIWpGX2ShpRVkVtxkgErcuOa6tb6YxrvshtbH9emBp075LgTsi4ljgjvz5XiQdDFwBnAycBFwx2RfTrEhkPVtSth5cj+Pa9rHU2B4kSYkkIu4iW12r6DyyheTJ//vBFoeeDayPiJ0R8SKwnolfXLOJ4o3J7TptXb+F49r6YZrjuh96aSM5rLEWcEQ8J+nQFmWOAJ4pPN+S75tA0ipgFcAw83uoltVGf75M0xbX84YPKrmqNrAGLFF0Mt3jSFpNcdnyTxgRqyNiSUQsmUtaLwqruZLaSCQtlfS4pE2SJtyqAg6S9IKk+4ADJF2UH/f2/Pl9kh6W9JHGKSep7cSdxbgeWtC5sjYzzNA2klael3Q4QP7fVv1StwBHFZ4fCWzt4T1tBinj1pak2cAXgXOAxcCFkha3KHpjRBwHbAL+Pt8XwFP5/pOBSyW9Fce19ahut7Z6SSRrgUZvlRXAd1qUWQecJWlh3hh5Vr7PrLNyrkhOAjZFxFMRsRv4Blk7yGSKcf0h4Nv543m88X1xXFtvZuIViaQbgH8G3ilpi6QPA1cDZ0raCJyZP0fSEknXAUTETuAq4O58uzLfZ9ZelNZra9L2jEJcvwX4qKRngHcCHyjE9VclPQA8CzwTEVsd19aTxNgeJEmN7REx2cLjZ7QouwG4qPB8DbCmq9pVWOzalVRu1oZHy33fPbtLPV+lpf8qWyRpQ+H56ohYnT+etD2jEdeSDgFeiYhdeTvIr+bdfxvem9/S+rakwyLi+emM61l70sql3v6IBWmj0AGOeVvazAkfPPy+pHIfOuDxpHILZ6d1sPnbVw5MKrfpx4uSygGM7khtk02fIaCjAbvi6MSTNlplTaGNZHujQTvfVhdO07E9IyJ2RETjl8G1wInNdYmIrcDDwC+V+BFthiqrjaRTRxJJKxsdSfLtosJrEwbVSpov6e8lPZZ3MLk6pR5OJFZd5bSR3A0cK+kdkoaA5WTtIOManUZyy4BH8/1HStovf7wQOAVI+4lt1k45vRGn1JEk367Lj203qPbTEfEu4HjgFEnndKqL59qyaiqpwTEiRiRdTNYYPhtYExEPS7oS2BARa4FLJC0DRsgGKK7MD/8PwGckBdktsk9HxIO918pmtPIa08c7kgBIanQkeSTh2PFBtfmx64GlEXED8F2AiNgt6V6yq/i2nEiskkR5XSAj4hbglqZ9lxceXwZc1uK49cB7y6mFWabE2G7VkeTkFuXOl3Qq8ATwsYh4ZpJj9xpUK+kg4D8Cn+9UEd/assqa7ilSzPolMa4XSdpQ2FY1n6bFqZu/ETcDR0fEe4HbeWP6n7bHSpoD3AD8ReOKpx1fkVh1OUlYXaXF9vaIWNLm9aSOJIWn1wKfLBx7WtOxdxaerwY2RsTnUirqKxKrrnIa282qp5y47rojCW0G1Ur6BPAm4HdSP46vSKyafNvK6qqk2O6lI0lE7JTUGFQL+aBaSUcCfwA8BtwrCeALjd5ek3EisepyIrG66nNHkvy1CYNqI2ILrdtP2nIimWYzaiR6yQZtmoiyzPlp2r8yIz9N/L7vSlvnHGDe7LRR8McO/XtSuf1npY0af2D360nlvrPj1KRyz29JX2dsv+dT/z7ljWyvW2w7kVhl+daW1VXdYtuJxKrJDelWVzWMbScSq66afdnMxtUstp1IrJLKHNluViV1jG0nEqssjdXs22aWq1tsO5FYNdXwPrIZUMvY7jiyXdIaSdskPVTY96l8vvoHJN2UT+7V6tjNkh7M58Hf0KqM2WSme64tx7b1S93mkEuZIuV6YGnTvvXAu/OJwJ5gkgEvudPzefDbzRljNtH0T5FyPY5t64eaTf3TMZFExF1kQ+uL+26LiMbIpR+SMF+92VRN9xWJY9v6pW5XJGW0kfwmcOMkrwVwW74w0F82LYFq1l7/v0x9ie25r6V98LG5aefbvTN9ZPtT2w9JKrf+wHcnlXtmz3NJ5X7w0jFJ5b7/5M8mlRvekvjHAYa39yHQ+h/bpeopkUj6A7LJwP56kiKnRMRWSYcC6yU9lv8KbHWuVcAqgGHm91Itq4Po7zQSZcV2Ma7nDbdsbrGZps+xPR26nkY+Xyz+A8B/iYiW+TUitub/3QbcRLY0ZEsRsToilkTEkrmkzc9j9dXoa9+PWwBlxvZecT20oPzK2sBJje1B0lUikbQU+DiwLCJem6TMAkkHNB6TzXf/UKuyZi1FpG0lcmzbPrGP43q6dby1JekGspW0FknaAlxB1pNlHtklPcAPI+Ijkt4KXBcR5wKHATflr88Bvh4Rt07Lp7BaKutXWZ4cPk+2ZsN1EXF1vr8R22+WNAo8Sxa3r5LF9n5kC/zsIPshuSsijsexbT0atCuOTjomkoi4sMXuL09Sditwbv74KeB9PdXOZq6SukBKmg18ETiTbHnRuyWtjYhHGrEtaSWwJCIubjr254CIiI35j6R7JB3k2LaeDGD33k48st0qq6QGyZOATfk//kj6BnAe8EinAyPiicLjrZK2AW8GflxKzWzGcmO72T6isbStgyOAZwrPt+T7mp2fj2b/lqSjJtRFOgkYAp7s8uOYjSshrivFicSqKZhKY/siSRsK26rCmVotI9h8Y+Fm4Oh8NPvtwFeKL0o6HPgr4DciYsC+4lY5qbE9QHxryyprCg2S29tMU7IFKF5hHAlsLRaIiB2Fp9cCnxyvg3Qg8PfAH0bED5NrZNbGjGtsN+ubcr5sdwPHSnoHWa+s5cCHigUkHR4RjSHYy4BH8/1DZGNEvhoRf1NKbRLM3pX2wee+krZm+/D2xLXdgVfn759U7h9mLU4q973htJHo27YfmFRu7tNpY8z225YePPN+4pHtvXIisUoqa/GfiBiRdDGwjqz775qIeFjSlcCGiFgLXCJpGdlI9p3AyvzwXwVOBQ7Je3YBrIyI+3qvmc1UXtjKbF+JKG3xn4i4Bbilad/lhceX0WKW34j4GvC1Uiph1lBibFeFG9utulKm2q7X99FmipLiWtJSSY9L2iTp0havr5T0Qr5uzn2SLiq8tkLSxnxbUdh/Yr7WziZJf6F85G07viKxyqrb5b9ZQxmx3W6wbVPRG1sMtj2YbJaSJWRp65782BeBL5FNNPpDsiv5pcA/tKuLr0ismgIYi7TNbJCkxnZn44NtI2I30Bhsm+JsYH1E7MyTx3pgad7V/cCI+Od8wtKvAh/sdDInEqsu39qyuionrnsZbDvZsUfkjzudcy9OJFZZ/ZpG3my6JcZ1u4G20Ntg28mOTTnnBG4jscqqW88Ws4bE2G430BZ6G2y7hWzm6+Kxd+b7j2zav9c5W/EViVVT6m0t5xobNOXF9fhg23zw7HJgbbFA3ubRMD7quLgEAAAIWElEQVTYlmxc1VmSFkpaSLamzrp8YO7Lkn4+763168B3OlXEVyRWSdmgLWcJq5+yYruXwbYRsVPSVWTJCODKiNiZP/4ocD2wH1lvrbY9tsCJxKpshk6POGsk7R+ZOclTqaS/97wdaTcpXpubNpXKq3PnJ5WbszPtn6LhnWnTvcx9Nf0f6tQpaUpVUmx3O9g2f20NsKbF/g3Au6dSj45RI2mNpG2SHirs+2NJzxYGuZw7ybFtB8uYtaOIpK3r8zu2rU+mM677IeXnx/VkA1KafTYijsu3W5pfLAyWOQdYDFwoKW2mN7N900ZyPY5t29dq2PbXMZFExF1k99amqpfBMjbjZfMRpWxdv4Nj2/pieuO6H3rptXVxPshlTd7q3yx1sIxZa+kLW5XNsW3Tq2YLW3WbSL4E/CxwHPAc8JkWZaY0sEXSqsbAmz3s6rJaVhtR2lK7U1VqbO8V17tfLa+WNrgSY3uQdJVIIuL5iBjNlx29luxSv1nHwTJN51wdEUsiYslc0havsZrrwxVJ2bG9V1wPLSi1rjbAfEUyYZDLLwMPtSjWcbCMWVt9aJR0bNs+UbPG9o6dtyXdQDaUfpGkLWRTD58m6Tiyj7sZ+K287FuB6yLi3MkGy0zLp7Ba0tj0Xt87tq1fpju297WOiSQiLmyx+8uTlN0KnFt4PmGwjFmSYNoHJDq2rS/2QWzvax7ZbpUkBm9QVllmv572r0xquaGX0t/7wKfTy5ZrT7/eeJ+rY2x70karrpIa23tcjvRWST+W9HclfzqbyWrW2O4rEquuEr5MvSxHmvsUMJ+8rcSsFAOWKDrxFYlVU+M+csrWXk+j0CPiDuDlKdbebHKpsT1AnEissjQ2lrR10MtypGbTooS4rhQnEquoxPaR7BZBuyVJe1mO1Gwa9G3qn2njNhKrpmAqX6Z2S5L2shypWfmmFtsDwVckVl3ltJH0shyp2fSoWRuJr0issvq9HCmApO8B7wL2z0e/fzgi1vVcMZvR6jaOxInEqqukL1uPy5H+UimVMCtyIpl+L/Pi9tvjW81jbBcB2/tRn5LV5XNAb5/l7W1fjYDRAbu+7+CVnzy7/a5bL61rXIM/S8OMi+1KJpKIeHPzPkkb2jSoDoy6fA7YB5+lZr/a6hzX4M8yJTWL7UomEjOgdl82s3E1i20nEqumAAZs3WqzJDWM7UFKJKv7XYGS1OVzwLR+loCo133kSTgeqsmxPQUDM44kImoRpHX5HDDNnyXIGiRTtgHmeKimSsR2gk4zWxfKXSApJC3Jnw9J+j+SHpR0v6TTCmUvzPc/kM9+vahTPQYmkdgM1Ic12832iXKWR2jMbH0OsBi4UNLiFuUOAC4BflTY/V+zasR7yGbG/oykWZLmAJ8HTs+nDHoAaDUr9l6cSKy6nEisrsqJ69SZra8C/gx4vbBvMXBHVpXYBvwYWEI2N52ABZIEHEjTlEKtVD6RpF66DQJJm/NLxvskbeh3faZC0hpJ2yQ9VNh3sKT1kjbm/11Y3jtOadLGgeTYrobKxnb7yUghYWZrSccDR0VE88Js9wPnSZoj6R3AiXm5PcBHgQfJEshiJll+uqjSiST10m3AnB4Rxw1gf/vrgaVN+y4F7oiIY8l+3ZT3j2EAY2Np2wBybFfK9VQxtvPJSAtbc7tN25mtJc0CPgv8botya8gSzwbgc8APgBFJc8kSyfHAW8lubbWc9aGo0omEHhclsvJExF1k81AVnccbU65/BfhgyW9a5ysSx3ZFVDa2O+s0s/UBwLuBOyVtBn4eWCtpSUSMRMTH8sR/HnAQsBE4LqtePBkRAXwT+MVOFal6IkldlGhQBHCbpHtaXKYOosMi4jmA/L+HlnfqqHuvLcd2tfU/tjtrO7N1RLwUEYsi4uiIOBr4IbAsIjZImi9pAYCkM4GRfPnpZ4HFkhqzMJxJwmzYVR9HkrIo0SA5JSK2SjoUWC/psfzXkDULiJr1tW/i2J6pSortxJmtJ3MosE7SGFny+LX8nFsl/Qlwl6Q9wNMUZsOeTNUTScdFiQZJRGzN/7tN0k1ktzcG+cv2vKTDI+K5fE2PbaWevWajf5s4tqttIGK708zWTftPKzzeDLxzknLXANdMpR5Vv7XVcVGiQSFpQd6fm/yS8izgofZHVd5aYEX+eAXwnVLPXu82Esd2tfU/tgdIpa9IJrt063O1unUYcFPWNZs5wNcj4tb+VimdpBuA08i6JG4BrgCuBr4p6cPAvwG/UtobRgxsj6wUju3qcGz3rtKJBFpfug2iiHgKeF+/69GtiLhwkpfOmMY3nbZTV4Fjuxoc272rfCKxmSqI0dF+V8JsGtQvtp1IrJpqONW2GVDL2K56Y7vNZDGWtnXQaSoSSSslvZBP73GfpIsKr63Ip8nYKGlF87FmXSkhrqvEVyRWSQFECb/aClORnEnW5fZuSWvzwVdFN0bExU3HHkzW8Lokr9I9+bEv9lwxm7HKiu0q8RWJVVNEWVckvUxFcjawPiJ25sljPRPnZDKbmtTYHiC+IrHKKqlBstVUJCe3KHe+pFOBJ4CPRcQzkxw7yNOYWEW4sd1sH3iZF9fdHt/quDJbbrhp6vLVhZlSU6YiuRm4ISJ2SfoI2SR970881mxKphDb26e9MiVxIrFKioiybiF1nIokInYUnl4LfLJw7GlNx95ZUr1shioxtivDbSRWdx2nIsnnUmpYxhuzna4DzpK0MF/Y6Kx8n5kV+IrEai1xhtRLJC0DRsjWpViZH7tT0lVkyQjgyohoXrfCbMZT1GyovpmZ7Vu+tWVmZj1xIjEzs544kZiZWU+cSMzMrCdOJGZm1hMnEjMz64kTiZmZ9cSJxMzMevL/AcwMKth6J6IIAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 4 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from SWESimulators import HLL\n",
    "importlib.reload(HLL)\n",
    "\n",
    "nx = 10\n",
    "ny = 15\n",
    "num_ghost_cells = 1\n",
    "dt = 0.01\n",
    "g = 9.81\n",
    "\n",
    "h0, hu0, hv0, dx, dy, nx, ny = gen_test_data(nx, ny, num_ghost_cells)\n",
    "plt.figure()\n",
    "plt.subplot(121)\n",
    "plt.imshow(h0)\n",
    "plt.colorbar()\n",
    "\n",
    "with Timer(\"construct\") as t:\n",
    "    sim = HLL.HLL(h0, hu0, hv0, \\\n",
    "                    nx, ny, \\\n",
    "                    dx, dy, dt, \\\n",
    "                    g)\n",
    "\n",
    "with Timer(\"step\") as t:\n",
    "    t = sim.step(10.0)\n",
    "    \n",
    "with Timer(\"download\") as t:\n",
    "    h1, hu1, hv1 = sim.download()\n",
    "\n",
    "plt.subplot(122)\n",
    "plt.imshow(h1)\n",
    "plt.colorbar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "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.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}