{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# 01 - Introduction to seismic modelling\n",
    "\n",
    "This notebook is the first in a series of tutorials highlighting various aspects of seismic inversion based on Devito operators. In this first example we aim to highlight the core ideas behind seismic modelling, where we create a numerical model that captures the processes involved in a seismic survey. This forward model will then form the basis for further tutorials on the implementation of inversion processes using Devito operators.\n",
    "\n",
    "## Modelling workflow\n",
    "\n",
    "The core process we are aiming to model is a seismic survey, which consists of two main components:\n",
    "\n",
    "- **Source** - A source is positioned at a single or a few physical locations where artificial pressure is injected into the domain we want to model. In the case of land survey, it is usually dynamite blowing up at a given location, or a vibroseis (a vibrating engine generating continuous sound waves). For a marine survey, the source is an air gun sending a bubble of compressed air into the water that will expand and generate a seismic wave.\n",
    "- **Receiver** - A set of microphones or hydrophones are used to measure the resulting wave and create a set of measurements called a *Shot Record*. These measurements are recorded at multiple locations, and usually at the surface of the domain or at the bottom of the ocean in some marine cases.\n",
    "\n",
    "In order to create a numerical model of a seismic survey, we need to solve the wave equation and implement source and receiver interpolation to inject the source and record the seismic wave at sparse point locations in the grid.\n",
    "\n",
    "\n",
    "<img src='./survey-ship-diagram.png' width=400>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The acoustic seismic wave equation\n",
    "The acoustic wave equation for the square slowness $m$, defined as $m=\\frac{1}{c^2}$, where $c$ is the speed of sound in the given physical media, and a source $q$ is given by:\n",
    "\n",
    "\\begin{cases}\n",
    " &m \\frac{d^2 u(x,t)}{dt^2} - \\nabla^2 u(x,t) = q \\ \\text{in } \\Omega \\\\\n",
    " &u(.,t=0) = 0 \\\\\n",
    " &\\frac{d u(x,t)}{dt}|_{t=0} = 0 \n",
    "\\end{cases}\n",
    "\n",
    "with the zero initial conditions to guarantee unicity of the solution.\n",
    "The boundary conditions are Dirichlet conditions:\n",
    "\\begin{equation}\n",
    " u(x,t)|_\\delta\\Omega = 0\n",
    "\\end{equation}\n",
    "\n",
    "where $\\delta\\Omega$ is the surface of the boundary of the model $\\Omega$.\n",
    "\n",
    "\n",
    "# Finite domains\n",
    "\n",
    "The last piece of the puzzle is the computational limitation. In the field, the seismic wave propagates in every direction to an \"infinite\" distance. However, solving the wave equation in a mathematically/discrete infinite domain is not feasible. In order to compensate, Absorbing Boundary Conditions (ABC) or Perfectly Matched Layers (PML) are required to mimic an infinite domain. These two methods allow to approximate an infinite media by damping and absorbing the waves at the limit of the domain to avoid reflections.\n",
    "\n",
    "The simplest of these methods is the absorbing damping mask. The core idea is to extend the physical domain and to add a Sponge mask in this extension that will absorb the incident waves. The acoustic wave equation with this damping mask can be rewritten as:\n",
    "\n",
    "\\begin{cases} \n",
    " &m \\frac{d^2 u(x,t)}{dt^2} - \\nabla^2 u(x,t) + \\eta \\frac{d u(x,t)}{dt}=q  \\ \\text{in } \\Omega \\\\\n",
    " &u(.,0) = 0 \\\\\n",
    " &\\frac{d u(x,t)}{dt}|_{t=0} = 0 \n",
    "\\end{cases}\n",
    "\n",
    "where $\\eta$ is the damping mask equal to $0$ inside the physical domain and increasing inside the sponge layer. Multiple choice of profile can be chosen for $\\eta$ from linear to exponential."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Seismic modelling with devito\n",
    "\n",
    "We describe here a step by step setup of seismic modelling with Devito in a simple 2D case. We will create a physical model of our domain and define a single source and an according set of receivers to model for the forward model. But first, we initialize some basic utilities."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Define the physical problem\n",
    "\n",
    "The first step is to define the physical model:\n",
    "\n",
    "- What are the physical dimensions of interest\n",
    "- What is the velocity profile of this physical domain\n",
    "\n",
    "We will create a simple velocity model here by hand for demonstration purposes. This model essentially consists of two layers, each with a different velocity: $1.5km/s$ in the top layer and $2.5km/s$ in the bottom layer. We will use this simple model a lot in the following tutorials, so we will rely on a utility function to create it again later."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Operator `initdamp` run in 0.01 s\n",
      "Operator `padfunc` run in 0.01 s\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAc8AAAGDCAYAAABN4ps8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3Xm8XdP9//HXW0JN/bVUUTFEqiXRwZC2Kd8S+m1RSvs11LfUUJHOaKtKEUT6paWor7bEUIRvTdWipVJDorSGdEIQVIghKiTmCEk+vz/WPnKcnHvO3vcM99xz38/HYz/OPXvvtfY6Kzf3c9baa6+liMDMzMzyW6avC2BmZtbfOHiamZkV5OBpZmZWkIOnmZlZQQ6eZmZmBTl4mpmZFdT24ClpHUlXSHpB0ouSrpS0bs60y0s6SdJsSfMl/UXSVq0us5mZWbm2Bk9JKwI3ARsB+wJfAt4H3CxppRxZnAscCIwDdgJmA9dL2qQ1JTYzM1ua2jlJgqSDgVOADSPi4Wzf+sBDwGERcUqNtB8G/gF8OSJ+me0bDEwHZkTEzq0uv5mZGbS/23Zn4PZS4ASIiJnAbcAuOdK+AVxalnYhcAmwnaS3Nb+4ZmZmS2t38NwYuLfK/unAiBxpZ0bEq1XSLgds0HjxzMzM6mt38FwVmFdl/1xglQbSlo6bmZm13OC+LkCrSRoLjE3vlt0cVuvT8piZdZbniXhVABtIS3XtFTUbro+I7ZtQsI7W7uA5j+otzJ5alZVp1+shLSxpgb5FREwEJgJIa8WbcdTMzMj+PAIwH/hGg7kdNUBaKO3utp1OundZaQRwX46062ePu1SmfR14eOkkZmZmzdfu4Hk1MErSsNIOSUOBLbNjtVwDLAvsXpZ2MPAFYHJELGh2Yc3MBhKR/sg2sg0U7Q6eZwOPAldJ2kXSzsBVwOPAWaWTJK0naaGkcaV9EfF30mMqp0kaI+mTpMdU1geOaeNnMDPrSiLdy2tkGyja+lkj4hVJ2wKnApNI/1Y3AodExMtlpwoYxNLBfX/gh8AE4J3AP4HtI+JvrS67mVm3K7U8rb62f1GIiFnArnXOeZT071i5fz7wnWwzMzPrEwOplW1mZjWUum2tPteTmZkB7rYtwsHTzMwAtzyLcD2ZmRnglmcRbV8M28zMrL9zy9PMzAB32xbhejIzM8DdtkU4eJqZGeDgWYSDp5mZvclBIR8PGDIzMyvIXzLMzAxwt20RDp5mZgZ4tG0RriczMwPc8izC9zzNzMwKcsvTzMwAd9sW4XoyMzPA3bZFOHiamRnglmcRriczMwPc8izCA4bMzMwKcsvTzMwAd9sW4ZanmZkBS7ptG9nqXkPaTdKvJT0mab6kGZJOkPT2QmWVDpcUkm4tkq5Z/CXDzMyAtt3zPBSYBfwAeALYFDgW2EbSFhGxuF4GkoYBRwHPtLCcNTl4mpnZm9oQFD4bEXPK3k+VNBe4ABgN3JQjj18AFwMb0kdxzN22ZmbWNhWBs+Su7HVIvfSSvghsBhzRzHIV5ZanmZkBWbdto1FhYa9SbZ293l/rJEmrAKcCh0XEXEm9ulgzOHiamRkAEgxuPHiuJmla2Z6JETGx52tqCDAeuCEipvV0XuYk4EHg/AZL2TAHTzMzA1LwXHZQw9k8GxEj811PKwNXkdqr+9c59xPAPsBmERENl7JBDp5mZtZ2klYArgGGAVtHxBN1kpwFnAs8Iemd2b7BwKDs/fyIWNCyAldw8DQzM6BJ3ba5rqNlgSuAkcCnIuKeHMmGZ9tXqxybB3wbOK1phazDwdPMzIAmDRiqdw1pGdJjJtsCO0XE7TmTblNl32nAIOBbwMPNKWE+Dp5mZpaIFIpa62fA7sAPgVckjSo79kREPCFpPeBfwPiIGA8QEVOWKq70PDC42rFWc/A0M7OkPZPb7pC9Hplt5Y4jzTZUCuMdOxeBg6eZmbVNRAzNcc6jpABa77zRjZeodxw8zcws8bIqubmazMxsCUeFXFxNZmaWtGfAUFdw8DQzs8Tdtrl17EgmMzOzTuXvGGZmlrjlmZuryczMlvA9z1wcPM3MLHHLMzff8zQzMyvI3zHMzCxxyzM3V5OZmS3he565OHiamVnilmduriYzM0scPHPzgCEzM7OC/B3DzMwStzxzczWZmdkSHjCUi4OnmZklbnnm5moyM7PEwTM3DxgyMzMrqO3BU9I6kq6Q9IKkFyVdKWndHOlGSpoo6QFJr0qaJeliSeu3o9xmZl2vtBh2I9sA0dYGuqQVgZuABcC+QAATgJslfSgiXqmRfE9gY+B0YDowBDgamCZpk4h4vKWFNzPrdu62za3d1XQgMAzYMCIeBpB0N/AQ8BXglBppfxQRc8p3SLoNmJnlO64lJTYzG0gcPHNpd7ftzsDtpcAJEBEzgduAXWolrAyc2b7HgDmkVqiZmVlbtDt4bgzcW2X/dGBE0cwkDQdWB+5vsFxmZuZ7nrm1u4G+KjCvyv65wCpFMpI0GDiT1PI8t/GimZkNcL7nmVt/rqYzgC2AHSOiWkAGQNJYYGx69462FMzMrF9y8Myt3dU0j+otzJ5apFVJOpEUEPeNiMm1zo2IicDElG6tyF9UM7MByMEzl3ZX03TSfc9KI4D78mQg6Ujg+8C3ImJSE8tmZmaWS7sHDF0NjJI0rLRD0lBgy+xYTZIOIj0XemREnNGiMpqZDUweMJRbu4Pn2cCjwFWSdpG0M3AV8DhwVukkSetJWihpXNm+PYHTgD8AN0kaVbYVHqlrZmYVSvc8G9kGiLZ+1Ih4RdK2wKnAJNI/1Y3AIRHxctmppe8/5cF9+2z/9tlWbiowukXFNjMbGDxgKLe2V1NEzAJ2rXPOo6R/xvJ9+wH7tapcZmbGgOp6bYRXVTEzMyvIDXQzM0vcbZubq8nMzBIHz9xcTWZmljh45uZ7nmZmZgX5O4aZmS3h0ba5OHiamVnibtvcXE1mZpY4eObmajIzsyXcbZuLBwyZmZkV5JanmZkl7rbNzdVkZmaJg2duriYzM0tK61lZXb7naWZmVpBbnmZmlrjbNjdXk5mZLeGokIuryczMErc8c3M1mZlZ4gFDuXnAkJmZWUFueZqZWeJu29zc8jQzsyUGN7jVIWk3Sb+W9Jik+ZJmSDpB0tvrpBspaaKkByS9KmmWpIslrd/rz9oAf8cwM7OkPfc8DwVmAT8AngA2BY4FtpG0RUQs7iHdnsDGwOnAdGAIcDQwTdImEfF4qwtezsHTzMyS9nTbfjYi5pS9nyppLnABMBq4qYd0P6pIh6TbgJnAgcC4FpS1R+62NTOztqkMgJm7stchRdJFxGPAnFrpWsUtTzMzS/puwNDW2ev9RRJJGg6sXjRdMzh4mpnZEo3f81xN0rSy9xMjYmJPJ0saAowHboiIaT2dVyXdYOBMUsvz3N4WtrccPM3MLGlOy/PZiBiZ63LSysBVwEJg/4LXOQPYAtgxIuYVTNswB08zM2s7SSsA1wDDgK0j4okCaU8ExgL7RsTkFhWxJgdPMzNL2nTPU9KywBXASOBTEXFPgbRHAt8HvhURk1pUxLocPM3MLGlD8JS0DHAxsC2wU0TcXiDtQcAE4MiIOKNFRczFwdPMzJZo/SQJPwN2B34IvCJpVNmxJyLiCUnrAf8CxkfEeABJewKnAX8AbqpI92JE3Nfykpdx8DQzs6Q93bY7ZK9HZlu540izDZXmOiqfi2D7bP/22VZuKmmChbZx8DQzs7aJiKE5znmUFCjL9+0H7NeKMvWGg6eZmSVduqpK9kzo50gt1lHAWsAKwLPADFLL9dKIeDBvnl1YTWZm1mtdtBh29jjMd4CDgNWAB4G/AzcC84FVgfWB7wLHSpoC/CAi7qiXt4OnmZkl3dfyfBh4jjRC99KIeKbaSZIEbAXsTRqMdEhEnF0r4+6qJjMz673uC54HRcSv650UEUHqup0q6RhgvXppuquazMzMMnkCZ5U0TwFP1TvPwdPMzJLua3nWJOn9wHDgjoh4ukjaAVRNZmZWT3TRgKFykn4KLBsRX8/e7wJcToqDL0j6ZET8LW9+XgzbzMwACMGiwY1tHWxHoHwqwONJsxVtDvyNNEFDbg6eZmY2ELwHeBTeXEP0A8API+LvpGn/PlIks87+nmBmZu2jjm89NuI1YKXs562Bl4C7svcvAf+vSGbdW01mZlZICBYOarRDcnFTytICfwO+Lmkm8HXgjxFRKuxQYHaRzBw8zcwMgJBYNLjRsPB6U8rSAkcD1wLTgReBb5Yd+xxLWqG5OHiamdmbFg3qzuG2EXG7pKGkR1NmRMTzZYfPI03dl5uDp5mZdSVJ1wK/Aa6OiH9HxIvAUvPWRsTVRfMuFDwlrclbZ6OfGREd20Y3M7P8ArGom2aGTzMFHQ/8QtKdwG+B3xZZPaUndYOnpJHAGGA7YN2Kw69Lugv4FXBRRLzUaIHMzKxvBGJhFwXPiBiTTfq+JbALKZadIGkGKZD+JiIK3ess6TF4ZkHzZNJM8/cAvyMt5TKHty7l8jHgROBEST8GfhIRr/WmMGZm1rcWddndvGzS91uz7XuSPkAaILQLcLikp4CrSd27N0fEwjz51qqlqcDZwNci4v5amUhaPivIYaSJF47Pc3EzM+scXdhtu5SIuBe4F5iQTZbweVL8+j3wCrBKnnxqBc/35p0oN2tpXgpcKmmNPGnMzMz6UkQ8CZwBnCHpnaQp/HLpMXgWnWG+LN2/e5POzMz61kBoeQJIWh1YvnJ/RFycN49eTSUhaZnKrUDadSRdIekFSS9KulJS5UCkPPkcLikk3Vo0rZmZVbeIQQ1tnUrSqpIukjSfNJvQzCpbbrnuDEtaATgG2B1Yu0q6yJOXpBWBm4AFwL5ZugnAzZI+FBGv5CzPMOAo4Jk855uZWX3dNtq2wjnAfwITgQdocCqkvMOqfg7sBVwDXNLARQ8EhgEbRsTDAJLuBh4CvgKckjOfXwAXAxviiR7MzKy+bYGDI+KXzcgsb+DZGTg0Ik5v8Ho7A7eXAidARMyUdBtptFPd4Cnpi8BmwH8DVzZYHjMzy6R7nl3bHnke6NVYnmry3qtcANR8XCWnjUlDhCtNB0bUSyxpFeBU4LCImNuE8piZWZluvecJ/IzU+9kUeb9inA/sCfyxweutCsyrsn8u+Z6tOYk0ee/5eS8oaSwwNr17R95kZmYDTjePto2IkyT9RNI9wA0sHYsiInLPUZA3eB5NmhtwMnB9lYsSEeflvWhvSPoEsA+wWTZjRC4RMZF0gxhprdzpzMwGmoCuHTAkaXvgq6S52TeuckpQYIKfvMFzc9L9ytVJo5WqXTRP8JxH9RZmTy3ScmcB5wJPZA+zQir/oOz9/IhYkKMMZmY28JwK/IO0jmfbRtueCTxH6i9u5KLTqR7xRwD31Uk7PNu+WuXYPODbwGm9LJeZmXX3gKH1SKNt/96MzPLW0kbAbhFxbYPXuxo4WdKwiHgEIFucdEvg8Dppt6my7zRgEPAt4OEqx83MLKduvudJanW+p1mZ5Q2eM4CVmnC9s0lN5qskHcWSPubHSd2yAEhaD/gXMD4ixgNExJTKzCQ9DwyudszMzIrr4uB5CHCepAciYqkFsYvKGzwPB34s6c6IeKy3F4uIVyRtS+p7ngQIuBE4JCJeLjtVpBZlr6YPNDOz4rq85XkpaczNnyW9SPXRtu/Nm1ne4HkUabDQg5Ie7OGiW+fJKCJmAbvWOedRUgCtl9foPNc0M7MB7zZSb2dT5A2ei0gDhczMrEt189y2EbF3M/PLFTzdwjMzGxi6dbStpI0jYnqN47tFxBV588t1T1HS2nWO5+qyNTOzzlW659ml0/P9oadYJmlX0mIjueUdkHN92cQElRf9BPC7Ihc1MzNrs3uAydkc6W+S9HngV8D/Fsksb/B8Gfi9pLesvC3pP4BrSc9vmplZP9blLc/dgJdIsWwFAEm7kJbZ/HlEHFoks7zBc0fgXcDlkpbJLroFKXD+HmjqjVgzM+sbCxnU0NapIuJV3hrLPg9cBkyMiEOK5pd3wNCz2aS6twHnSpoIXEeaJH6vIhO1m5lZZ+ry9TxLsezTwJ+BK4AzI+Jbvckrdy1FxKOSdgCmAnsB1wB7RsSi3lzYzMw6S7dNkiBpXA+H/gxsDTxTdk5zliST9OUeDl0N7ABMBvaVVLpqS5ckMzMzK+jYOsePKfu5aUuSnVMn7S8qLurgaWbWz3VTyxNYtlUZ1wqe67fqomZm1nm6bYahVt5W7DF4NjIBvJmZ9T/dNmBI0tsiYkEr0nnVEjMze1OXPec5U9K3JL09z8mSPirpSuCweuf2GDwl/UPS51UaEVT/omtLOl1S3YuamZm1wSHAQcDTki6XdJCkrSWNkPReSSMl7SHpZEkzSI9jzqP+mJ+a7fMLSYtXnyHpMuBPwD+BOcAC0rpow4CPAp8lDfu9ETij1x/TzMz6TLc9qhIRl2UtyV2BA4Afs/QgIgFPktb7PCsiHsqTd617nqdIOhcYk130YJZeC02kQHoV8MmImJrnomZm1nm6LXgCRMRCUmC8NJtidjNgLWB54DnggYiYWTTfmneGI+IF4CfATyStC4yqvChwZ29uyJqZWefpptG2lSLiNdIECQ0rMsPQLGBWMy5qZmbWn3XPmGQzM2tItz2q0kp+VMXMzID2LEkmaTdJv5b0mKT5kmZIOiHP4ySSlpd0kqTZWdq/SNqqKR++IH/FMDOzN7VhwNChpFuAPwCeADYlzUG7jaQtImJxjbTnkpYV+x7wCPAN4HpJH4+If7S01BUcPM3MDGjb9HyfjYg5Ze+nSpoLXACMBm6qlkjSh4EvAl+OiF9m+6YC04HxwM6tLHQld9uamVnbVATOkruy1yE1ku4MvEF67KSU10LgEmA7SW+rdV1JX5a0QsHi9sjB08zMgCUDhhrZemnr7PX+GudsDMyMiFcr9k8HlgM2qHONs4GnJP1U0ojeFXOJ3J9U0jBgD2Bd0nOe5SIiDmi0MGZm1rfaPUmCpCGkbtcbImJajVNXJU2dV2lu2fFaNgLGAvsC35T0Z+BM4PKIeL1YqXMGT0mfAy4jtVSfIc0qVK5y5iEzM+tnmjTD0GqSyoPgxIiYWO1ESSuTZqhbCOzf6IVryabd+56kHwC7AV8FJgGnSTo/K2euqfkgf8vzeGAKsFcP/dX9wnuYzViO6+timJl1jMqo1oTg+WxEjKx3Unb/8RrSHOlbR8QTdZLMA9arsr/U4pxb5dhSIuIN4FfAryRtRGp9fgf4jqQpwI8j4vp6+eS95zkMOLk/B04zM+sMkpYFrgBGAp+JiHtyJJsOrC9pxYr9I4DXgYcLXH8lSWOBi4GtgHuAY4CVgWslHVMvj7zB8wHgXXkLZmZm/U/pUZVGtnokLUMKWtsCn4uI23MW7xrSiii7l+U1GPgCMDnPHOuSNpH0C+Ap4HRSbPtERGwSERMi4mOkntZv1csrb7ftYaR+4Tsi4pGcaczMrB9p0/R8PyMFwB8Cr0gaVXbsiYh4QtJ6wL+A8RExHiAi/i7pUlIsWhaYCXwNWB/Yq95FJd0JbA48DpwInNNDb+ofgHH18uuxliTdUrHrXcD9kh5i6b7liIitMTOzfq0No213yF6PzLZyx5FmGxIwiKV7R/cnBd0JwDtJa0xvHxF/y3HdZ4HPAb+LiFqDXP8GvK9eZrW+YizmraNoZ+QonJmZWY8iYmiOcx4lBdDK/fPJBvf04tITgH9WC5ySVgI+HBF/zh5b+Ve9zGothj26F4UzM7N+qhsXwy7zJ+DjwJ1Vjm2UHc/94XMNGJK0j6SqA4YkrSppn7wXNDOzztSOAUN9aKmWbJm3AYuKZJb3zvAvSRH7uSrH1s+OX1jkwmZm1nm6aT1PSesCQ8t2bSqpcoa8FYADSAOJcstbS7Ui9kqk2SHMzKwf68Ju2/1Jz29Gtv28yjkitTrrPp5SrtZo202Azcp2fVbSBypOWwHYE8g9pZGZmVmbXAjcSgqQk4GDWHry+QXAjKKTANVqee5CitiQInblkOKS50hNXjMz68e6reUZETNJz4Mi6VPAnRHxUjPyrhU8TwPOJ0XsR4D/Av5ecc4C4N91npkxM7N+opuCZ7mIuLGZ+dV6VOUF4AUASesDs3uzbIuZmfUPpdG23ULSg8BuEXF3NsFPrYZeRMSGefPONWAoIh7LCrINadTtEOBJ4C8RcXPei5mZmbXRHcBLZT83rZc073qeqwKXA9uQZh6aB6ySDulmYI+IyLUcjJmZdaY2zW3bNhHxpbKf925m3nlXVTkd+AiwN7BCRLybNNJ2n2z/T5tZKDMz6xuLGNTQNlDk/YrxWeCIiPi/0o5sQdGLs1bphFYUzszM2qfbRtuWk3QysHpELDUjnqQLSINfD8ubX96W5yJ6fpZzBgWnNTIzs87T5dPzfR64oYdjN2THc8sbPK8iLThazZ7Ab4tc1MzMrM2GALN6ODYrO55b3m7ba4BTJf2eNHDo38AawB7AxsDBkrYtnRwRNxUphJmZdYZuGjBU4XlgGDClyrENgFeKZJa3lq7IXtdhyUKm5X6dvYo0FLij2+5mZra0br7nCdwIHCnpmvKp+CStBhxBz126VeUNntsUydTMzPqfLg+eRwN3AQ9Juhp4gtRVuwvwBj1PQVtV3kkSphYspJmZ9UMdPuin1yLiEUkfIT0dsgOwKmlu9t8BR2fz4OZWqHM7a96OAt4FXBMRc7O10V6PiMVF8jIzM2uniHgE+GIz8so7w5CAH5PWO1uOdF/zI8Bc0kjcW4Hjm1EgMzPrG902w1BPJG1I1vKMiAd7k0feR1WOAL4JjAc+xlsXx74G2CnvBSWtI+kKSS9IelHSldlq33nTD5d0uaRnJc2XNEPSwXnTm5lZdaV7nt06w5Ck/SQ9CdxHavTdL+lJSfsWzSvvV4wxwPiIOEFSZe08DLw3TyaSVgRuIi1lti+pBTsBuFnShyKi5lBhSSOz9FOyMr0AvA9YOefnMDOzGjo9APaWpP8GzgOmAuOAp4E1gb2A8yS9FhGX5s0vb/AcAtzew7HXgZVy5nMg6TmbDSPiYQBJd5NmL/oKcEpPCSUtQ1oV/MaIKJ8Jwqu6mJlZPd8HfhURe1XsP1fSxcDhQO7gmbfb9kngAz0c+zDZSt057AzcXgqc8OZK37eRhgvXMhoYTo0Aa2Zmvdfl3bYbkhpg1UwCNiqSWd7geTkwTtKWZftC0vuB7wKX5MxnY+DeKvunAyPqpP2P7HV5SbdLekPSM5JOl7RCzuubmVkPArp5btuX6XkKvrWy47nlDZ7HAg8At7BkgvjLgXuy9yfmzGdV0lqgleaS1getZa3s9VJgMvAp0gjgMcD/9ZRI0lhJ0yRNezVnIc3MBqY02raRrYNdD/yPpI+X78ye/TweuK5IZnknSZgvaTTp+ZjtSIOEnssueHFELCxy0V4qBfqLImJc9vOUbADTiZKGR8T9lYkiYiIwEWAtqWmriJuZdZsun2HoMFID8FZJjwGzSQOGhgKPkO6J5pb7a0JELCL1C08qcoEK86jewuypRVruuez1jxX7J5NavpsCSwVPMzOziHhK0iak3spPkOLOP4CfAudFRKFu27yTJCwPjATeQ+oWnw38NSJeK3Ix0r3NjavsH0F67qZe2lo8w5GZWYO6uOVJFiBPy7aG1LznKeltkn5Kuic5lXS/8TJS0/c5SSdLWq7A9a4GRkkaVnaNocCW2bFariM9H7pdxf7ts9dpBcphZmYVunwx7Kaq1/L8HbAtaQq+a0kLhoq0NNlOwLdJrcbP5Lze2aSZiq6SdBSpFXs88DhwVukkSesB/yJNzDAeICKek3QCcLSkF0mTJYwkPex6QfnjL2ZmVly3Tc8n6SFSnMkjImLDvHn3WEuSdictRbZbRPymyinnSNoVuFTSf0XElTlK9kq2aPappHunIq2xdkhFf7NIa4JWtozHAy8BXwcOJXUfn4Tn1TUza4ou67a9g/zBsxBFVM9X0pXAaxFRcwZ6Sb8ClouIXVtQvqZaS4qxfV0IM7MOMhF4KkIAy438YKw2rd4dtNpma9hfI2JkM8rWyWrd89wU+H2OPH4HbNac4piZWV/p8hmGmqpW8Hw36R5nPbOA1ZtTHDMz6yuBWLR4UENbJ5P0IUmXSXpa0uuSNsv2T5D06SJ51QqeK5JGt9bzOrB8kYuamVkHCli4cFBDW6eStAXpHuiHgSvhLc3kZYCvFsmv3rCqIeWPlfRg7SIXNDMz6wM/Ig1Q3Zmlg+U00tJkudULnlfkyEO0aDSTmZm1T4RYtLB7HlWpsDmwa0QslqSKY88CaxTJrFYt7V+0ZGZm1n+l4Nm5Xa8NWgD0tALXmsALRTLrMXhGxAVFMjIzs34u6ObgeStwkKTflu0r9Zp+Gbi5SGZd2z43M7NiIsTCN7o2eI4jBdC/k5bUDGBvST8GRgEfLZJZ3vU8zczM+q2I+DswGnietEa1gENIT4tsU21Jy1rc8jQzs4xYvKh7w0JE3AVsLWlFYDVgXkS81Ju8ureWzMysmAC66J6npPOA8yPilvL9EfEq+SYB6pGDp5mZJaGuCp7AF4B9Jc0CLgQmNWsFLt/zNDOzJICFamzrLGsAY4BHgaOAGZJuk3SgpHc0krGDp5mZdaWIeDkifhkR2wBDgaOBVUjrR8+WdImkHSQVjoUOnmZmtsTCBrcOFRGPR8T/RMQI0qMp5wHbklYGe1LSyUXyc/A0M7Mk6NrgWS4i7oyIbwJDgFNJK4N9u0geHjBkZmZJKXh2OUkbAPsAe5O6c18ELiuSh4OnmZl1PUmrAHuSguZHSV8V/gj8APhtRLxWJD8HTzMzSwJ4o68L0TySlgV2IgXMHYDlgPuAw4GLImJ2b/N28DQzsySARX1diKb6N/AOYC4wEbggIv7ajIw9YMjMzJZow4AhSWtL+l9Jf5H0qqSQNDRn2ndJ+qmkRyTNlzRT0hmS3l3l9KnArsBaEXFQswInuOVpZmYl7RswtAGwB/BX4E/Ap/Mkyhaxvhp4P2mVlPuBEcB4YKSkj0dEaZkxIuLzTS73mxw8zcys3W6JiDUAJI0hZ/AE3gdsAXwlIiZm+6ZIWgz8ghRUZzS7sNU4eJqZWdKmlmdELO4VErRwAAAUjElEQVRl0uWy1xcr9j+fvbbtVqSDp5mZJZ3/nOd04BbgaEkPAw+Qum3HAdcVXZOzEQ6eZmaWNCd4riZpWtn7iWVdrA2JiJD0GWAScFfZod8DuzfjGnk5eJqZ2RKNB89nI2JkE0rSk7NJc9N+lTRgaDhwHHCFpM820CVciIOnmZn1C5J2BP4b+M+IuDHbfYukR4DJwGeBq9pRFj/naWZmSWmGoUa21vpg9npXxf47s9fhLS9BxsHTzMyS0gxDjWyt9XT2+tGK/R/LXp9seQky7rY1M7OkjaNtJe2W/bh59rqDpDnAnIiYmp2zkDSl3gHZOVcCPwQulHQ8abTtRsAxwOPAb9pTegdPMzPrG5dXvP959joVGJ39PCjbAIiIFyWNAo4FDgPeA8wGrgGOjYiXW1jet3DwNDOzpI0tz4hQb86JiMeBA6qc3lYOnmZmlnT+JAkdw8HTzMyWcPDMxcHTzMwStzxz86MqZmZmBbnlaWZmiVueuTl4mplZUpphyOpy8DQzs6Q0w5DV5eBpZmZLuNs2Fw8YMjMzK8gtTzMzSzxgKDcHTzMzSxw8c3PwNDOzxKNtc/M9TzMzs4Lc8jQzs8SPquTm4GlmZkv4nmcuDp5mZpZ4wFBuDp5mZpZ4wFBuHjBkZmZWkFueZmaWeMBQbm1veUpaR9IVkl6Q9KKkKyWtmzPtupIukDRL0nxJD0qaIGmlVpfbzKzrle55NrINEG1teUpaEbgJWADsS/qnmgDcLOlDEfFKjbQrATcAywJHA7OAjwDHAe8DvtDa0puZDQADKAA2ot3dtgcCw4ANI+JhAEl3Aw8BXwFOqZF2S1KQ3C4iJmf7bpa0KnCopBUj4tXWFd3MrMt5wFBu7e623Rm4vRQ4ASJiJnAbsEudtMtlry9W7H+e9DnUrEKamZnV0u7guTFwb5X904ERddLeQGqh/kjSCEkrS9oWOBg4s1aXr5mZ5VAaMNTINkC0u9t2VWBelf1zgVVqJYyI1yT9B/BrUrAtOQf4ZtNKaGY2UHmShNz6zaMqkpYHLgVWB75EGjD0UWAc6Z/7az2kGwuMBXhHW0pqZtZPOXjm1u7gOY/qLcyeWqTlDgBGAxtExL+yfbdIegGYKOnMiPhnZaKImAhMBFhLit4W3MzMrKTdwXM66b5npRHAfXXSfhCYVxY4S+7MXocDSwVPMzPLyaNtc2v3gKGrgVGShpV2SBpKegzl6jppnwZWkbRBxf6PZa9PNqmMZmYDlwcM5dLu4Hk28ChwlaRdJO0MXAU8DpxVOknSepIWShpXlvZ84CXgWkn7StpG0veAk4G/kh53MTOz3vIMQ7m1NXhmj5NsCzwITAIuBmYC20bEy2WnChhUXr6IeBQYBfyDNCvRtaRJFyYCn4qIxW34CGZm3cvBM7e2j7aNiFnArnXOeZQqkx5ExH3AHq0pmZmZWT795lEVMzNrMQ8Yys3B08zMEi9JlpuDp5mZLTGA7ls2ou3reZqZmfV3bnmamVni6flyc/A0M7PEA4Zyc/A0M7PEA4Zyc/A0M7PE3ba5ecCQmZlZQW55mpnZEm555uLgaWZmiQcM5ebgaWZmiQcM5ebgaWZmiQcM5eYBQ2ZmZgW55WlmZolbnrk5eJqZWeIBQ7k5eJqZ2RIeMJSL73mamZkV5JanmZktEX1dgP7BLU8zM7OCHDzNzKytJK0t6X8l/UXSq5JC0tAC6YdIOk/S05IWSJop6YTWlXhp7rY1M7N22wDYA/gr8Cfg03kTZkH2NmAmcBDwb2BolmfbOHiamVm73RIRawBIGkOB4AmcCTwJbBMRpQdrpja5fHU5eJqZWaY9D3pGxOLepJP0XmA7YJ+ywNknfM/TzMwypSmGGtlaasvsdb6kP2b3O+dJulDSu1p98XIOnmZmlim1PBvZWE3StLJtbBMLuFb2eh7wILAD8H1gR+B6SW2Lae62NTOzTFMmt302IkY2oTDVlILjlIj4RvbzTZJeAC4hdele16JrVy2ImZlZp3sue/1jxf7J2eum7SqIW55mZpbp+Jnhp9c53quBSL3hlqeZmWWacs+zlW4HniZ1z5bbPnu9q9UFKHHL08zMyrRnQU9Ju2U/bp697iBpDjAnIqZm5ywELoiIAwAiYqGkw4HzJZ0JXEmaHOGHwBTgprYUHgdPMzPrG5dXvP959joVGJ39PCjb3hQRF0haTBpluz8wF7gIOCIi2jatvYOnmZll2nfPMyLU23MiYhIwqemFKsDB08zMMk15VGVAcPA0M7NMx4+27RgOnmZmlnHLMy8/qmJmZlaQW55mZpZxt21eDp5mZpZxt21eDp5mZpZxyzMvB08zM8u45ZmXBwyZmZkV5JanmZll3G2bl4OnmZmVcbdtHg6eZmaWccszL9/zNDMzK8gtTzMzy7jlmZeDp5mZZfyoSl4OnmZmlnHLMy8HTzMzy7jlmZcHDJmZmRXklqeZmWXcbZtX21uektaW9L+S/iLpVUkhaWjOtMtIOkLSo5Jek/RPSbu2tsRmZgNFqdu2kW1g6Itu2w2APYB5wJ8Kpj0eOBY4A9gBuB24XNJnmllAM7OBqdTybGQbGPqi2/aWiFgDQNIY4NN5EklaHTgUODEiTs523yxpA+BE4NpWFNbMbODwgKG82t7yjIjFvUy6HbAccFHF/ouAD0pav6GCmZmZ5dSfBgxtDCwAHq7YPz17HQHMbGuJzMy6igcM5dWfgueqwPMRERX755YdNzOzXnO3bV79KXj2iqSxwNjs7YLj4N6+LE8HWA14tq8L0cdcB66DEtcDbLjkx9nXw7GrNZjfgKjP/hQ85wHvlKSK1mepxTm3ShoiYiIwEUDStIgY2dpidjbXgesAXAclrodUB6WfI2L7vixLf9KfZhiaDrwNeG/F/hHZ633tLY6ZmQ1U/Sl4/oF0J3uviv17A/dGhAcLmZlZW/RJt62k3bIfN89ed5A0B5gTEVOzcxYCF0TEAQAR8YykU4AjJL0E/A34ArAtsHPOS09s1mfox1wHrgNwHZS4HlwHvaKlB6+24aJSTxedGhGjy865ICL2K0s3CDgCOBBYE5gBjI+IK1paYDMzszJ9EjzNzMz6s/50z7MqSetIukLSC5JelHSlpHVzpl1e0kmSZkuan01Wv1Wry9xsva0DSSMlTZT0QDZJ/yxJF/fH2Zoa+T2oyOfwbLGCW1tRzlZrtB4kDZd0uaRns/8TMyQd3MoyN1uDfxPWlXRB9n9hvqQHJU2QtFKry91MXoCj9fp18JS0InATsBGwL/Al4H2kOW/z/LKfS+oCHgfsBMwGrpe0SWtK3HwN1sGepJmbTidNtH84sBkwTdI6LSt0kzXh96CUzzDgKOCZVpSz1RqtB0kjgTtIo9rHAJ8BfgIMalWZm62ROsiO3wBsBRxN+vznAN8FzmthsVvBC3C0WkT02w04GFgEbFC2b33SFBnfqZP2w6TpNPYv2zeYdB/16r7+bG2qg3dX2bcesJh0L7nPP1+r66Ain+uBs4ApwK19/bna/LuwDOlxr9/09efowzr4dPY34dMV+0/M0q/Y15+vQD0sU/bzmOxzDc2RbnXSNKjHVey/Ebi7rz9XJ239uuVJGmV7e0S8Od9tpEdWbgN2yZH2DeDSsrQLgUuA7SS9rfnFbYle10FEzKmy7zFgDjCkyeVspUZ+DwCQ9EVSq/uIlpSwPRqph9HAcOCUlpWuPRqpg+Wy1xcr9j9P+nKhZhWy1cILcLRcfw+eG1N9ur3pLJk8oVbamRHxapW0y5G6PfqDRupgKZKGk7593t9gudqpoTqQtApwKnBYRFSdqaqfaKQe/iN7XV7S7ZLekPSMpNMlrdDUUrZWI3VwA/AQ8CNJIyStLGlbUmv2zIh4pblF7Uh5FuAw+n/wXJXUp19pLrBKA2lLx/uDRurgLSQNBs4ktTzPbbxobdNoHZwEPAic38Qy9YVG6mGt7PVSYDLwKeDHpC6//2tWAdug13UQEa+RvkQsQwoWL5G6K38HfLO5xexYXoAjp/40t6213hnAFsCOEVHtD1DXkfQJYB9gsyp/MAaS0hfpiyJiXPbzlOzZ6hMlDY+I/tQbUZik5UlfHlYnDTSaBXyUNKBwIfC1viuddZr+HjznUf3bZE/fPivTrtdDWuhhovkO1EgdvEnSiaTVZ/aNiMlNKlu7NFIHZ5Fa2U9Ieme2bzAwKHs/PyIWNK2krdVIPTyXvf6xYv9k0oCZTekfXfmN1MEBpHu/G0TEv7J9t0h6AZgo6cyI+GfTStqZerUAx0DU37ttp5P66CuNoP5E8dOB9bOh7ZVpX2fpPv9O1UgdACDpSOD7wEERMamJZWuXRupgOPBV0h+N0rYlMCr7uT+1Nhr9/1BLbwegtFsjdfBBYF5Z4Cy5M3sd3mDZ+gMvwJFTfw+eVwOjsufzAMgeBN4yO1bLNcCywO5laQeT5sud3I9aG43UAZIOAiYAR0bEGS0qY6s1UgfbVNn+SRp0sg3Qn6Z+bKQeriMNFNmuYn9piapp9A+N1MHTwCqSKgcLfix7fbJJZexkXoAjr75+VqaRDViJ1EK8hzQMfWfSH75HgJXLzluPdM9iXEX6S0itizHAJ0l/KF8j3f/q88/X6jogTZKwmPSHc1TFNqKvP1u7fg+q5DeF/vmcZ6P/H47J9v8P8J+kSTPmA+f39WdrRx0AQ0mPqTxImmBhG+B72b5plD072R82YLds+wXpOc+vZe+3LjtnIXBuRboTs7+D3yF1Y/8i+zuxU19/pk7a+rwATfgFWRf4dfYL/hLwWyoeBs7+UwRwbMX+FUjPtT2d/bLcAYzu68/UrjogjS6NHrYpff252vV7UCWvfhk8G60H0nOM38mCz+vAY8B4YNm+/lxtrIMRwGXA46QvDg8CJwOr9PXn6kU91P2/nb0/vyLdINJMW4+ReiPuBnbr68/TaZsnhjczMyuov9/zNDMzazsHTzMzs4IcPM3MzApy8DQzMyvIwdPMzKwgB08zM7OCHDytI0i6VNJcSWtW7B8k6S5JD3XS0liShkoKSfuV7dtP0pernLtfdu7QNhaxdO1lJP1D0qFl+47NytOyua0lHSLpHkn+G2Ndyb/Y1im+RXpg++cV+w8FNgfGRMT8tpeqZ7OBjwO/L9u3H7BU8MzO+XiWpt32Bt7D0vXaamcB7ybN1GPWdRw8rSNExDPAt4HPS9odQNL7gWOBsyJiah8WbykRsSAibo+IOTnOnZOd2xfzJR8KXBhLL/reUtkXnQuz65t1HQdP6xgRcSFpYuozJK1GWipsDnBYvbRlXaNbSfqtpJclPSfpZ5XdvZLeI+lCSc9KWiDpbkl7V5yzpqQLJD2VnTNb0u8krZ4df0u3raQpwNbAltn+yPZV7baVtKykCZIelfR69jpB0rJl55Su8RVJ47MyPC/pGklr56iTj5FWCqm7mLWk7bM6OyPr6i1d+6uSTpD0tKSXJF0kaUVJG0i6PkvzsKRqLcxLgBGStqh3fbP+pr+v52nd5yukZZHuAIaRFuZ+qUD6i0hzk/6cJQsZr0TqUkXSSsBU0pqPPyDNYbo3MEnSihExMctnEmny8O9l56xBWjygcgm7kq9n1x6UfQZIc6v25AJgD9Ik7LeSFiE/MvvMX6w49wjgz6Qu4dWBn2TXGl0jf0grorxEmhi9R5L2Ac4BxkfEhGxf+bWnkLpfRwA/Jk0SvilwNmne168Bv5Q0LSLKlzb7R3b97bPym3WPvp5c15u3yg04gXT/89cF0uyXpTmzYv+RwCLg/dn7b2bnja447wbgGWBQ9v5l0vqmPV1vaJbPfmX7plBlQvmysg3N3n+A6pOSH5Xt/1DFNaZUnHdotn+tOnVyHXBblf3HZukHk1r1b5DuKVf7fDdV7L8y27932b5VSKtzHFPlWn8iLfHX579X3rw1c3O3rXUUSf8P+BLpD/RHJL29YBaXVby/hHR74qPZ+62AJyNiSsV5F5EGuJQW/b0L+J6kgyV9UGVNsSbYquyalWWA1P1b7tqK9/dkr+vWuc5apG7vnpwKHEdaMeOcHs65ruL9A9nr9aUdETGP9MVjnSrp52TlMOsqDp7WaU4itWR2JHVRnlAw/b97eD8ke12V6qNeny47DmlR9KtJLbO7gScljWvSoxela1SWo7IMJXMr3pcGHi1f5zrLl51bzX+TFv2+ocY58yrev15jf7XyzCct/WfWVRw8rWNIGg0cCBwVEdcBE4CvFRxwskYP75/MXucCa7K0NcuOExHPRMQ3ImIIsBFp7dPjWHI/sxGlYFhZjjUrjjfqOdIXkZ58ktR6vU7Syk26ZqVVgWdblLdZn3HwtI6QjYg9m9Rd+tNs949Ig4fOkbRczqz2qHi/J2mAyx3Z+6nA2pK2rDjvi6Sux/sqM4yIGRHxA1Jr6wM1rr2AfK2sW8rKVm6v7HVKjjzyeIA0AKkn00mDjt5H6wLo+sCMFuRr1qccPK1TjCeNbh0TEYsBIuINYAywIWngTx6fkXSSpE9JOhI4hvSc40PZ8fOBh4ArJY3JHtGYBHwKODoiFkl6Rzar0SHZ8U9KOp3Uiptc49r3AR+Q9AVJIyVtWO2kiLgX+BVwrKRjsrKOIw3k+VVE3FMtXS/cArxX0rt6OiEi7icF0PcC1/fiHnOPJL0TeD9LviyYdQ0/qmJ9TtJI0gQJ/1MZOCLiTkk/BQ6XdFm89VGIavYGvkt6fOJ1Umv2zQf1I+IVSVuTHrk4EXg7qWX0pYgoDdh5DfgbqQt5PVLLdQawV0RcVePaPyIF+nOAlUmt3NE9nLsf8Ajp8ZOjgKey9MfV+XxFXEX6LDuRHo2pKiJmZHVyMzBZ0nZNuv6OpH+D3zQpP7OOoYjo6zKYNSybrOCXwPsi4uE+Lk7HkHQ+sHZE/GcfXPs64NmI+FK7r23Wam55mnW344D7JY2MiGntuqikTYBtgY3bdU2zdvI9T7MuFhEzSV3Eq7f50muSJpBwL4B1JXfbmpmZFeSWp5mZWUEOnmZmZgU5eJqZmRXk4GlmZlaQg6eZmVlBDp5mZmYF/X/aZzJmhfncMAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 576x432 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#NBVAL_IGNORE_OUTPUT\n",
    "from examples.seismic import Model, plot_velocity\n",
    "\n",
    "# Define a physical size\n",
    "shape = (101, 101)  # Number of grid point (nx, nz)\n",
    "spacing = (10., 10.)  # Grid spacing in m. The domain size is now 1km by 1km\n",
    "origin = (0., 0.)  # What is the location of the top left corner. This is necessary to define\n",
    "# the absolute location of the source and receivers\n",
    "\n",
    "# Define a velocity profile. The velocity is in km/s\n",
    "v = np.empty(shape, dtype=np.float32)\n",
    "v[:, :51] = 1.5\n",
    "v[:, 51:] = 2.5\n",
    "\n",
    "# With the velocity and model size defined, we can create the seismic model that\n",
    "# encapsulates this properties. We also define the size of the absorbing layer as 10 grid points\n",
    "model = Model(vp=v, origin=origin, shape=shape, spacing=spacing,\n",
    "              space_order=2, nbl=10)\n",
    "\n",
    "plot_velocity(model)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Acquisition geometry\n",
    "\n",
    "To fully define our problem setup we also need to define the source that injects the wave to model and the set of receiver locations at which to sample the wavefield. The source time signature will be modelled using a Ricker wavelet defined as\n",
    "\n",
    "\\begin{equation}\n",
    "  q(t) = (1-2\\pi^2 f_0^2 (t - \\frac{1}{f_0})^2 )e^{- \\pi^2 f_0^2 (t - \\frac{1}{f_0})}\n",
    "\\end{equation}\n",
    "\n",
    "To fully define the source signature we first need to define the time duration for our model and the timestep size, which is dictated by the CFL condition and our grid spacing. Luckily, our `Model` utility provides us with the critical timestep size, so we can fully discretize our model time axis as an array:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from examples.seismic import TimeAxis\n",
    "\n",
    "t0 = 0.  # Simulation starts a t=0\n",
    "tn = 1000.  # Simulation last 1 second (1000 ms)\n",
    "dt = model.critical_dt  # Time step from model grid spacing\n",
    "\n",
    "time_range = TimeAxis(start=t0, stop=tn, step=dt)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The source is positioned at a $20m$ depth and at the middle of the $x$ axis ($x_{src}=500m$), with a peak wavelet frequency of $10Hz$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAhEAAAGBCAYAAADYEOPMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3XmYZHV97/H3t6t6nZmenhWYnQEEBgHFUYlLQBIDLgENLokrLuASE71eE/UxoiLXYDToVaKCyRUVY4wGAzFGkUWMCNFxARnWkRmGGQZm656ZXqq6lu/945zqrqmpqq6qPlXdp+vzep5+qvvUOad+dRi6Pv39LcfcHREREZF6dcx0A0RERCSeFCJERESkIQoRIiIi0hCFCBEREWmIQoSIiIg0RCFCREREGqIQISIiIg1RiBAREZGGKESIiIhIQxQiREREpCHJmW7AbLd06VJft27dTDdDRESkJX75y1/udfdlteyrEDGFdevWsWnTppluhoiISEuY2aO17qvuDBEREWmIQoSIiIg0RCFCREREGqIQISIiIg1RiBAREZGGKESIiIhIQ2IRIsxslZl93szuNLNRM3MzW1fjsR1m9kEz22ZmKTO728wubG6LRURE5r5YhAjgeOBVwCDw33Ue+3Hgo8BVwIuAu4Bvm9mLo2ygiIhIu4nLYlM/cfejAMzsrcAf1XKQmS0H3gdc4e6fDjffZmbHA1cA329GY0VERNpBLCoR7p5v8NBzgS7gupLt1wGnmtmx02qYiIhIG4tFiJiGU4A0sKVk++bwcUNrmyMiIjJ3zPUQsRgYcncv2b6/6HkJuTv37Bhix+DoTDdFRERiYK6HiIaY2SVmtsnMNu3Zs2emm9MyN933JOdfdQd/es1dM90UERGJgbkeIgaBATOzku2FCsR+ynD3a9x9o7tvXLaspruhzgm/enQQgB2DYwyOjM9wa0REZLab6yFiM9ANHFeyvTAW4r7WNmd2u/fxAxPf/3bngSp7ioiIzP0Q8QMgA7y2ZPvrgHvdfWvrmzQ7uTv37jzIS047BlCIEBGRqcVlnQjM7BXht88IH19kZnuAPe5+e7hPFviqu78FwN13m9mVwAfN7BDwK+DVwDnA+S19A7PcwbEsB8YyPH31AP/zyH4NrhQRkSnFJkQA3y75+Qvh4+3A2eH3ifCr2IeAYeDdwNHAg8Cr3P17zWlmPD15KAXAUf09HNXfzZMH0zPcIhERme1iEyLcvXRwZE37uHsOuDz8kgqePBiEiOULujmqv2fiZxERkUrm+pgIqdHusPKgSoSIiNRKIUKAye6M5f1BJWLfSJpMrtHVxkVEpB0oRAgQVCIWdCfp60pyVH8P7rDnkKoRIiJSmUKEALD7UIpl/d1AMC4i2KYQISIilSlECAD7R8ZZOi8IDwN9XQAMjWrVShERqUwhQgAYGs3Q39sJwEBf8HhgLDOTTRIRkVlOIUKAIDAUwsNAr0KEiIhMTSFCgCAwLAzDQ6EiMTSqECEiIpUpRAjj2Tyj47mJCkRnooP53UmFCBERqUohQia6LRaG3RkAC3s7GRrTwEoREalMIUI4EIaFQncGBIMrD6gSISIiVShEyGQlojREaGCliIhUoRAhE2MfCutDQKE7QyFCREQqU4iQiRBRXIlY2NupgZUiIlKVQoQwnM4C0N8zeWf4+d1JRsLtIiIi5ShEyESImNddHCI6GcvkyOV9ppolIiKznEKEMJzO0pkwupOT/xzmh1WJYVUjRESkAoUIYTiVZV53EjOb2Da/OxE8pxAhIiIVKEQII+ks84u6MiDozoAgYIiIiJSjECEcKhci1J0hIiJTUIgQhlPlKhEKESIiUp1ChDAynp2oPBRMhAh1Z4iISAUKETIxsLJYIVRorQgREalEIUIYTmdZUKE745BChIiIVKAQIQyXnZ2h7gwREalOIaLN5fLO6HjuiO6MRIfR25lgOK37Z4iISHkKEW1uZDyoNJRWIgDmdScYGc+1ukkiIhITChFtbjQdhITSSgRAb1eClEKEiIhUEIsQYWarzew7ZnbAzA6a2fVmtqbGY9eY2VfNbLuZjZnZQ2Z2uZnNa3a742A0rET0dSWOeK6vM8moQoSIiFRw5J+fs4yZ9QG3AmngjYADlwO3mdlp7j5S5dh5wM1AJ/BhYDvwTOBjwAnAq5vb+tmvEBJ6y4SInq4EoxmFCBERKW/WhwjgYmA9cKK7bwEws3uAh4G3AVdWOfa5BGHhXHe/Kdx2m5ktBt5nZn3uPtq8ps9+hRAxr+vIfwp9nQnGxjU7Q0REyotDd8b5wF2FAAHg7luBO4ALpji2K3w8WLJ9iOC9G22u0J1RrhLR15VgTJUIERGpIA4h4hTg3jLbNwMbpjj2ZoKKxSfNbIOZzTezc4B3A1+q1hXSLiYqEd0VujM0JkJERCqIQ4hYDAyW2b4fWFTtQHdPAc8jeJ+bgUPALcD3gHdF28x4KoSEvs5K3RkKESIiUl4cxkQ0zMx6gG8By4HXEwysfBZwKZAF3lHhuEuASwDWrKlpEkhsFcY89JWpRPSpEiEiIlXEIUQMUr7iUKlCUewtwNnA8e7+u3DbT8zsAHCNmX3J3e8uPcjdrwGuAdi4caM32vA4KCwmVW6KZ29XUmMiRESkojh0Z2wmGBdRagNw3xTHngoMFgWIgp+HjydPs22xV6g09CTLhIjOBOPZPLn8nM5RIiLSoDiEiBuBM81sfWGDma0jmL554xTHPgEsMrPjS7Y/O3zcGVEbY2s0naWvK0FHx5ETVQrViVFN8xQRkTLiECK+DGwDbjCzC8zsfOAG4DHg6sJOZrbWzLJmdmnRsdcSDKb8vpm90cxeYGZ/BXwa+CXBNNG2NprJle3KgMlpnxpcKSIi5cz6EBFOwzwHeAj4OvANYCtwjrsPF+1qQIKi9+Tu24Azgd8QrHL5fYLFq64BXuju+Ra8hVltbDxHX5mFpiDozgA0LkJERMqKw8BK3H07cOEU+2yjzOJR7n4f8KrmtCz+RsLujHImuzMUIkRE5EizvhIhzTWWyZVdrRImuzMUIkREpByFiDY3Np6b6LYoVejmSKk7Q0REylCIaHPpbJ6eCiGiEC5UiRARkXIUItpcKpOjp7P8P4NeTfEUEZEqFCLaXCqbK7vQFEwOrNQUTxERKUchos2lMnm6K1UiNMVTRESqUIhoc+lMju4KlQjNzhARkWoUItpcqsrAyu5kBx2m7gwRESlPIaKN5fPOeDZfcWClmdHbqduBi4hIeQoRbSydDVb9rtSdAboduIiIVKYQ0cbS2fA24BUqERDM0BjTFE8RESlDIaKNpTJBJaLSmAgIQoS6M0REpByFiDZWWM66O1n5n0FPZ0LdGSIiUpZCRBsrjImYqhKh2RkiIlKOQkQbK1QiphoToe4MEREpRyGijU2EiCqzM9SdISIilShEtLFUYYrnlLMzFCJERORIChFtLD0xsLLamIik7uIpIiJlKUS0sVQNAyu7Ozsm9hMRESmmENHGahlY2Z1MMJ7Nk897q5olIiIxoRDRxmrpzigEjLSqESIiUkIhoo1NrhNRZbGpMGCkNENDRERKKES0scnujOpTPAFSWYUIERE5nEJEG0tl8nQYJDus4j4T3RkZdWeIiMjhFCLaWCqTo6czgVm1EKFKhIiIlKcQ0cbS2XzVrgyYrESkVIkQEZESChFtLJXJ0VPlDp6ggZUiIlKZQkQbS2XzdE9RiSg8rxAhIiKlYhMizGy1mX3HzA6Y2UEzu97M1tRx/Mlm9m0z22tmY2b2oJm9u5ltnu3SmRzdU1Ui1J0hIiIVJGe6AbUwsz7gViANvBFw4HLgNjM7zd1Hpjh+Y3j8j4G3AgeAE4D5TWz2rJeqaUxE8HxaAytFRKRELEIEcDGwHjjR3bcAmNk9wMPA24ArKx1oZh3A14Bb3P3lRU/d1rzmxkMwO6N6JaJQqVB3hoiIlIpLd8b5wF2FAAHg7luBO4ALpjj2bOBkqgSNdhV0Z9RWiVB3hoiIlIpLiDgFuLfM9s3AhimOfV742GNmd5lZxsx2m9nnzKw30lbGTDDFc6oxERpYKSIi5cUlRCwGBsts3w8smuLYFeHjt4CbgBcCf0cwNuKfo2pgHBUWm6qmJ6mBlSIiUl5cxkRMRyEoXeful4bf/9jMEsAVZnayu99ffICZXQJcArBmTc0TQGInlclPOTsjmegg2WEaWCkiIkeISyVikPIVh0oVimL7wscflWy/KXx8eukB7n6Nu290943Lli2rq6FxkspOXYmAoEtDlQgRESkVlxCxmWBcRKkNwH01HFtN2346pjNTT/GEYK0I3TtDRERKxSVE3AicaWbrCxvMbB3w3PC5av6LYH2Jc0u2nxc+boqmifHi7kElYoruDIDuZEIDK0VE5AhxCRFfBrYBN5jZBWZ2PnAD8BhwdWEnM1trZlkzK4x9wN33AX8LvN3MPmFmf2hmHwAuBb5aPG20nYzn8rgz5bLXEFQidCtwEREpFYuBle4+YmbnAJ8Bvg4YcAvwHncfLtrVgARHhqPLgEPAO4H3AbuATwEfb3LTZ63CGIepBlZCYUyEKhEiInK4WIQIAHffDlw4xT7bCIJE6XYnWGxKC06FCrMtahkT0Z3UmAgRETlSXLozJGLpuisR6s4QEZHDKUS0qUL3RO1TPFWJEBGRwylEtKl0Nqgs1DzFUyFCRERKKES0qclKRA3dGcnEROgQEREpUIhoU5OzM2oYWKkxESIiUoZCRJuqqxLR2UFa3RkiIlJCIaJN1TcmIqEpniIicgSFiDY1UYmooTujJ5kgk3NyeW92s0REJEYUItpUobLQXWN3BqAZGiIichiFiDZVGChZUyUi7PJQiBARkWIKEW0qXUclorCqZUrTPEVEpIhCRJuq9wZcwTGqRIiIyCSFiDaVzuToTnZgdsT9yo6gMREiIlKOQkSbSmVyNU3vhGCxqeAYdWeIiMgkhYg2lc7ma1poCiYHX6a1VoSIiBRRiGhTqUyupiWvYbI7I61KhIiIFFGIaFOpTB2VCA2sFBGRMhQi2lQ6W8eYiHAGh+7kKSIixRQi2lQqk69poSlQJUJERMpTiGhTqWyupoWmQJUIEREpTyGiTaUy+ToGVqoSISIiR6o5RFjgfDP7tJl9xczWhtvPMrMVzWuiNEMwJkKVCBERaVyylp3MbBHwfeDZwCFgPvB54FHgYmA/8JdNaqM0QbqOSkQy0UGyw1SJEBGRw9RaifgUsBp4LrAEKF4r+WbgDyJulzRZsGJl7b1Z3ckOVSJEROQwNVUigAuA97n7nWZW+ufrdoKAITFSz7LXEIyLUCVCRESK1fqn6HxgZ4Xneji8MiExUM+y16BKhIiIHKnWT5EHgT+q8NxZwG+jaY60QjaXJ5v3msdEgCoRIiJypFq7M74AXGVmB4B/DrcNmNmbgHcBlzSjcdIcqbCiUFclojOhSoSIiBympk8Rd78GuBL4GLAl3Pwj4Brgs+7+jeY0b5KZrTaz75jZATM7aGbXm9maBs7zATNzM/tpM9oZB4WKQj1jIrqTHapEiIjIYWqtRODuHzCzLwIvBJYD+4AfufsjzWpcgZn1AbcCaeCNgAOXA7eZ2WnuPlLjedYDfwPsblZb46BQUah12WsIqhaqRIiISLGaQwSAuz8K/GOT2lLNxcB64ER33wJgZvcADwNvI6iS1OKLwDeAE6nzvc8lhYpCrcteA3QnEwyNjjerSSIiEkMVP0jr7Spw9+3Tb05F5wN3FQJE+HpbzewOgumnU4YIM3sNcAbwZ8D1zWpoHEyECFUiRERkGqr9Nb6NoNugVrV/ItXvFOCGMts3A6+c6uBwxc3PAH/t7vvN2ntGairTwMDKpGZniIjI4aqFiDczGSK6CcYSHAT+FXgSOBp4FbAA+HgT2wiwGBgss30/sKiG4z8FPARcG2GbYiudVSVCRESmr2KIcPdrC9+b2WeBXwEvd3cv2n4Z8O/Ahia2cVrM7PnAG4Azits+xTGXEE5bXbOm7gkgs15alQgREYlArZ8ifwZcXfohHP78JeA1UTesxCDlKw6VKhTFrgb+CdhhZgNmNkAQnhLhz92lB7j7Ne6+0d03Llu2bLptn3UameLZ09kx0Q0iIiICtc9QmA9U+jRdDsyLpjkVbSYYF1FqA3DfFMeeHH69vcxzg8D/Aj47rdbFzMQUz7rWiUiQzuZwd9p9TImIiARqDRE/Bj5hZve7+y8KG83sWcD/CZ9vphuBT5vZ+sK6FGa2juCuoh+Y4tgXlNn2WYKBoH/B5OJZbWNydkbt3Rk9nR3kHTI5pyupECEiIrWHiHcR3PL7LjN7jGBg5VEEd+/cGj7fTF8OX+MGM/sbggGfHwceI+iuAMDM1gK/Ay5z98sA3P3HpSczsyEgWe65dtDYipXBvulsjq46woeIiMxdtS57vRU4iaBL4BaC1SpvIVjo6WR339asBoavPwKcQzDD4usEC0ZtBc5x9+GiXY2gwqBPuSoauXdGYV+NixARkYJ6lr3OEFQEvty85lR9/e3AhVPss40abkvu7mdH06p4KszOqGeKZ3ElQkREBPQXe1tKZXN0JoxER+1jG7pViRARkRI1VSLMbCvVV690dz8umiZJs6UyubpuvgWqRIiIyJFq7c64nSNDxBLgOcAwwR02JSZSmTzddQyqBI2JEBGRI9UUItz9onLbw4WbfkAwc0NiIp3N1TW9E1SJEBGRI01rTIS7DxHcl+LSaJojrZDO5OuamQGTlYi0KhEiIhKKYmBlClgVwXmkRVKZXF1rRMDkmhKqRIiISEHNUzxLmVkSeCrwUYJlqSUmUtn6Q0Sh+0NjIkREpKDW2Rl5Ks/OOAi8JLIWSdOlM/m6x0SoEiEiIqVqrURcxpEhIgU8CvyXux+ItFXSVKlsjv7ezrqOUSVCRERK1To746NNboe0UKqhgZWqRIiIyOFq+iQxs1vN7KQKzz3FzLRORIyks40sNqVKhIiIHK7WP0fPBvorPLcAOCuS1khLBItN1VeJSCY6SHaYKhEiIjKhnk+SSgMrjyNYtVJiIpXJ1XXzrYLuZIcqESIiMqHimAgzexPwpvBHB64xs0Mlu/USTPO8pTnNk2YIFpuqP0T0dCZUiRARkQnVKhF5IBd+WcnPha99wBeBtzS3mRKVfN4Zz9U/xRNUiRARkcNVrES4+1eBrwKY2W3AO9z9gVY1TJojnQ1CQOOVCIUIEREJ1DrF8wXNboi0RioTdEfUO8UToCvZMXG8iIhItTERbwD+0933hd9X5e5fi7Rl0hSpbCFEqBIhIiLTU60ScS1wJsG4h2unOI8DChExULgLZ+NjIlSJEBGRQLUQcSywq+h7mQOmW4kYGstE3SQREYmpagMrHy33vcRbYXZFI2MiupMdpFWJEBGRUP2fJBJrEwMrG1hsSmMiRESkWLWBlVupvEplKXf346JpkjRTIQTUu+w1aEyEiIgcrtqYiNupPURITBRCQCPLXvd0JhQiRERkQrUxERe1sB3SIpPrRDQSIjrUnSEiIhM0JqLNTHRnNDTFM6hEuKtAJSIidYQIMzvBzL5qZg+Z2Uj4eK2ZHd/MBkq00tOsROQdsnmFCBERqXHZazM7G/g+MAb8J/AkcBTwx8Crzew8d7+9WY2U6ExvimciPEeOzoSKWCIi7a7WT4K/B34NrHX3N7j7X7n7G4B1wG/C55vKzFab2XfM7ICZHTSz681sTQ3HbTSza8zsATMbNbPtZvYNM2vLBbSmOyYC0LgIEREBag8RG4BPuvtw8UZ3PwR8Ejgl6oYVM7M+4FbgJOCNwOuBE4DbzGzeFIf/adi+zwEvAj4AnAFsMrPVTWv0LJXO5ukwSHZY3ccWVyJERERq6s4AdgBdFZ7rAnZG05yKLgbWAye6+xYAM7sHeBh4G3BllWM/6e57ijeY2R3A1vC8lzalxbNUKpOjpzOBWQMhQpUIEREpUmsl4pPAx8xsRfFGM1sJfAT4RNQNK3E+cFchQAC4+1bgDuCCageWBohw26PAHmBlxO2c9VLZXENdGaBKhIiIHK7WSsRZQD/wiJndxeTAyjPD788OB19CsHrlGyNu5ynADWW2bwZeWe/JzOxkYDlw/zTbFTupTJ6eBqZ3gsZEiIjI4WoNEc8DsgR39VwbfsHkXT6fX7RvM+b/LQYGy2zfDyyq50RmlgS+RFCJ+KfpNy1e0tk83apEiIhIBGoKEe4+l2YyXAU8B3iJu5cLJpjZJcAlAGvWTDkBJFZSmVxDC02BKhEiInK4uEz2H6R8xaFShaIsM7uCIBy82d1vqrSfu1/j7hvdfeOyZcvqbuxsVhhY2YhCJUK3AxcREai9OwMI1moAVgM9pc+5+61RNaqMzZSfRroBuK+WE5jZh4D3A3/h7l+PsG2xks7kVYkQEZFI1Lpi5XrgG8CzCpvCRw+/d6CxP29rcyPwaTNb7+6PhG1aBzyXYN2HqszsL4HLgQ+5+1VNbOesl87mGOirNFu3usJYCo2JEBERqL0S8Y/AGuA9wAPAeNNaVN6XgXcBN5jZ3xCElo8DjwFXF3Yys7XA74DL3P2ycNufAp8FfgDcamZnFp33oLvXVMmYK1KZfENLXgMTszpUiRAREag9RDwTuMjd/62ZjanE3UfM7BzgM8DXCaoftwDvKVlF0wgqIsWfkueF288Lv4rdDpzdpGbPStNaJ0KVCBERKVLPipWtrj4cxt23AxdOsc82JrtaCtsuAi5qVrviZlpjIgqViIwqESIiUvvsjE8A76/hPhUyy02nEpFMdJDsMFJZVSJERKT2dSK+bmYnAdvCFStLp1U2Y5VKaYLpTPEE6E52qBIhIiJA7bMzLgI+COQI7oBZ2rXRjFUqJWLuPq1lryG4hbgqESIiArWPifgY8F3gLe4+1MT2SBON54IKQqPLXoMqESIiMqnWP0mXAF9QgIi3VPjh3+jASihUIhQiRESk9hDxU+DkZjZEmq+wXPV0xkR0JTu07LWIiAC1d2e8G/hXMxskWLTpiPtVuLv+PJ3lVIkQEZEo1Roi7g8fv1Zln2Yuey0RSGenX4noViVCRERCtYaIy9AMjNgrVCKmEyJ6OhMMjWWiapKIiMRYretEfLTSc2Z2NvCGiNojTZSaqEQ03p2hSoSIiBQ09GliZseb2WVmtpXgHhavirZZ0gyFe150J6dXidANuEREBOoIEWa20MwuMbM7gAeBDxEMsHwnsKJJ7ZMIpSe6M6ZXidANuEREBKYIEWbWYWYvNrNvAbuALwFrgX8Id3mPu1/t7geb3E6JQCqCgZWqRIiISEHFMRFm9vfAa4DlQIpgxcqvAjcD/cC7WtFAic7EwMppdGeoEiEiIgXVBlb+L4IZGd8HLnL3fYUnzEwzNWKoMMWzexrdGapEiIhIQbVPk38CDgEvAR40s6vM7FmtaZY0Q1SViFzeyeQUJERE2l3FEOHuFwNHA68FNgFvA+40s/uB96N1I2JnYnbGNCsRgKoRIiJSfWClu6fc/Zvufh6whsnbgX8AMOAKM3udmfU0v6kyXemJKZ7TmJ0RBhCNixARkZo/Tdx9l7v/nbs/FXgWwQyNEwiWwt7VpPZJhFLZPD2dHZhZw+codIWoEiEiIg39Serum9z9LwjWh7gQ+HGUjZLmSGVy05reCapEiIjIpFrvnVGWu2cIpn5+N5rmSDONjefonW6IKFQiMqpEiIi0u8Y7xyV2gu6M6YWIwmqXhYWrRESkfSlEtJGx8Qi6M1SJEBGRkEJEG0lnc9O6bwaoEiEiIpMUItqIxkSIiEiUFCLaSCo7/e6MQiUirUqEiEjbU4hoI5FUIjpViRARkYBCRBtJZfLTWvIaoCepMREiIhKITYgws9Vm9h0zO2BmB83sejNbU+OxPWb2KTPbZWZjZnanmf1+s9s826QyqkSIiEh0YhEizKwPuBU4CXgj8HqCJbdvM7N5NZzin4CLgUuBlxIs0/1DM3tac1o8O0WxYuVEJUIrVoqItL1prVjZQhcD64ET3X0LgJndAzxMcHfRKysdaGanA68B3uzuXwm33Q5sBi4Dzm9u02cHd2csgkpEMtFBosN07wwREYlHJYLgg/6uQoAAcPetwB3ABTUcmwG+VXRsFvgX4Fwz646+ubNPJufknWmvEwFBNUKVCBERiUuIOAW4t8z2zcCGGo7d6u6jZY7tAo6ffvNmv7HwQ3+63RkQjItQJUJEROLSnbEYGCyzfT+waBrHFp5vid2HUqQzeY5Z2EMy0dr8lo4wRMxkJeJQKsOTB1OMpHOMjGcZG8+RzTvuQZeNA+6Qn/jeZ6SdIiKtsrC3k7NPXD4jrx2XENFSZnYJcAnAmjU1TQCpyXV3bedztzxMssNYu6SPZx27hD8+7Rh+77glmFlkr1NOXCsRQ6PjXP+rndz+0B5+u/MA+0fGW/K6IiJxseGYfoWIKQxSvuJQqcpQeuzaCsfCZEVigrtfA1wDsHHjxsj+lH3pacewYmEP2/eP8uATh/je3Y/zzZ9v56kr+/nEy0/ltFUDUb3UEVLhlMzpDqwE6G5BJWI8m+ean/yOz9+6hXQ2z/HL5/PCk49i3dJ5rBjoYX53kr6uJH1dCRIdRocZZmBG8D3B92A0OZ+JiMyorhZXtovFJURsJhjbUGoDcF8Nx77czPpKxkVsAMaBLeUPi95TjlrAU45aMPFzKpPjxrsf59M/fJCXf+FnXPmq07ngaSub8tqTlYjp/2Pr7UpMnK8ZRtJZ3vLVX3DXI/t50VOP5i/OOYENK/qb9noiItKYuAysvBE408zWFzaY2TrgueFz1fwH0Am8sujYJPBq4CZ3T0fd2Fr1dCZ41cbV/Oi9Z7Fx7SLe863f8L17Hm/KaxUqB1FUIno7E02rRGRyed70lV/wi22D/P0rT+eLr3uGAoSIyCwVlxDxZWAbcIOZXWBm5wM3AI8BVxd2MrO1ZpY1s0sL29z91wTTOz9rZm81sz8gmN55LPCRFr6Hihb2dnLtm57FM9Ys4q++fQ/b9o5E/hqFD/3uiEJEsyoRn7/lYX6+bT9//8rTufAZq5ryGiIiEo1YhAh3HwHOAR4Cvg58A9gKnOPuw0W7GpDgyPf1JuArwOXAfwKrgfPc/VdNbnrNersSXPWaM0jzHGmKAAAdVElEQVQmjA/9+28jn1UQaSWiK8HoePQh4rc7DnDVbVu48IxVvOzpzenWERGR6MRlTATuvh24cIp9thEEidLtY8B7w69Z6+iFPfzvFz6Fj/7Hfdz+0J5IR9sWBlZGMiaiM0GqCSHiMzc/RH9vJx85f6qlP0REZDaIRSWinbzm2WtZsbCHL//3I5Get9D90NsVTSUi6u6Me3ce4NYHdvOW5x5Lf09npOcWEZHmUIiYZbqSHbz2zLXcsWUfW3YPT31AjQrdGT3JaMZERN2dcfVPHmFBT5I3PnddpOcVEZHmUYiYhV61cTWdCeMb//NoZOeMuhKRzubJ56MZtzGcznLT5id42dNWqgohIhIjChGz0LIF3Zz31GP47q93ks1FszJkYUxEdzKaMREAqWw01YibNj9BOpvnZU9fEcn5RESkNRQiZqmXnHo0Q6MZNj061YKctUllcvR0dkSyvHahmhFVl8YNv3mcVYt6OWPNVLdBERGR2UQhYpZ6/gnL6Ep0cPN9T0ZyviBETL8rAyYrEWMRhIgDoxl+umUvf3z6iqbfP0RERKKlEDFLzetO8nvHLeHm+5+MZM2IsfFcJGtEwGQlIopVK+98ZC+5vHPOSTNz8xgREWmcQsQs9ocbjmLbvlF+t2f6K1imsvnIKxFRdGfcsWUf87oSPG11824+JiIizaEQMYv9/glLAbjrkX3TPtfYeITdGWElIoq1Iu743V6edexiOmfwLnQiItIY/eaexdYs7mPp/G5+GcHgynQ2F8lqlVA0JmKaIWLXgTEe2TPCc49fGkWzRESkxRQiZjEz45nrFrHp0f3TPlczxkRMd2Dlnb8LKizPOU4hQkQkjhQiZrlnrF3EY/vHePJgalrnSWVn3+yM3zw2xLyuBCcdvSCKZomISIspRMxyG9ctBmDTtul1aTSlEjHN7oy7dxzgqSsX0tGhqZ0iInGkEDHLnbKin57Ojml3aaQyebqjHhMxjUrEeDbP/bsOcrpmZYiIxJZCxCzXmejg5GP62fz4wWmdZ3Q8y7yuaO78HsXAyoeePMR4Ns+pKxdG0iYREWk9hYgYOPmYfh7YdXBai06NjOfoi+DmWwDJRAddiY5phYjf7jwAwGmrFCJEROJKISIGTj6mn4OpLI8faGxwZTaXZzybpy+iSgRAT2fHtLoz7tlxgP6eJGsW90XWJhERaS2FiBjYcEwwe+H+Brs0RsOKwbzuaCoREAyunE6IuHfnAU5dtVD3yxARiTGFiBg48eh+AO7f1ViIKHzY90bUnQHQ15VsuDsjn3e27B7mxKP6I2uPiIi0nkJEDMzvTrJ2SR/3P9FYiBhJZwEiG1gJ0NOZaPjeGTuHxhjL5DjhqPmRtUdERFpPISImTjp6AffvOtTQsYUP+6gGVgL0dnY0fBfPLbuHAThhuUKEiEicKUTExElH97Nt30hDH9yTISK6SsR0ujMe3h2EoeMVIkREYk0hIibWL5uHOzy6b7TuY0fGg+6MvggHVk6nO+PhJ4dZtqCbgb6uyNojIiKtpxARE8ctC/5qf2TPcN3Hjqab0J3RlWi4O+Oh3cM8ReMhRERiTyEiJtYtnQfAI3tH6j52dDz6gZW9Da4T4e5sefIQJyzXTbdEROJOISIm5ncnOaq/m0f2NBIioq9E9HUlJ8JJPXYdSDEynuM4jYcQEYk9hYgYWb90Po/sbaA7owkDK3s6E6Qy+bqP27YvCEHHhZUVERGJL4WIGDl22Ty2NtidYRYsVR2V3s4E47k82Vx9QaIwMHStQoSISOwpRMTI+qXzGBrNsH9kvK7jRtI55nUlI11iutA1Uu80z237RuhKdnBMf09kbRERkZkRixBhZh1m9kEz22ZmKTO728wurOG4fjO71Mx+Zmb7zGwo/P5lrWh31BqdoTGWyUa65DVMThcdSdcXIh7dO8rqRb10dOieGSIicReLEAF8HPgocBXwIuAu4Ntm9uIpjlsDvBO4HXgd8GrgIeC7ZvbnTWttkxRmaNTbpRFUIqINEfO7g/EVw+n6Blc+un+UdUvUlSEiMhdEN9KuScxsOfA+4Ap3/3S4+TYzOx64Avh+lcO3AuvdvXiFph+a2Wrg/cA/NKPNzbJioAcz2DE4Vtdxo+NZeiMcVAmTIWKkjhDh7jy6b4TfW78k0raIiMjMiEMl4lygC7iuZPt1wKlmdmylA919pCRAFGwCVkTXxNboTiY4ur+HxwbrW7VydDz6SsS8BioRe4bTjI7nWLukL9K2iIjIzIhDiDgFSANbSrZvDh83NHDO3wcemE6jZsrqRX3s2F9fJWJkPEdfd3MqEfWEiO2FmRkKESIic0IcQsRiYMjdvWT7/qLna2ZmlwBnAn9bbR8z22Rmm/bs2VNXY5tt1eLeuisRY+NZ+jqbMyainu6MbWGI0JgIEZG5oeUhwsz+0My8hq8fN+G1zwY+B3zN3b9RaT93v8bdN7r7xmXLlkXdjGlZvaiPJw6mSGdrnxUxks5FevMtmOzOqCdEbN83QofBykW9kbZFRERmxkwMrPwZcHIN+xX+3B4EBszMSqoRhQrEfmpgZs8EbgRuBd5aY1tnndWL+3CHXUOpidkaUxkdz0Z63wyYrEQcqiNE7Bgc45iFvXQm4lAAExGRqbQ8RIQDHesZj7AZ6AaO4/BxEYWxEPdNdQIzOxX4IfAb4EJ3z9Tx+rPKqvCv+McGR+sIEblI75sBweqXHVZfJWLH0JiqECIic0gc/iT8AZABXluy/XXAve6+tdrBZnYC8CPgEeCl7l7fqMRZZvXiYFDiYzUOrszm8qSz+UjvmwFgZszvTta12NTOwTFWDShEiIjMFbN+nQh3321mVwIfNLNDwK8IFo06Bzi/eF8zuwVY6+7Hhz8vJwgQXcBHgA0lSz//2t3TzX8X0Tm6v4fOhNU8uHI0E/0dPAvmdydrnp2RzeV54mBKlQgRkTlk1oeI0IeAYeDdwNHAg8Cr3P17JfslOPw9bQDWht+X7gtwLLAt0pY2WaLDWDHQy2P7awsRY4U7eEY8sBKCwZXDqdpCxK4DKXJ5n+iOERGR+ItFiHD3HHB5+FVtv7NLfv4xMOdu0rB6UR+P1bhqZWHMQtQDKyEIESPjtYWInUNBe1cOaI0IEZG5Ig5jIqTE6sW97KixEjEaViKivgEXwIKe2rszCkt1qztDRGTuUIiIoVWL+tg3Ms5oDVWAplYiupI1z87YGYaIFQO6BbiIyFyhEBFDhXEFtdyIq1ApmN/TpO6MGmdn7BgcZfmCbrqT0VdERERkZihExNDkNM+puzQKIWJBE0LEgp4kh1K1Lbmxc2hMgypFROYYhYgYWr0oCBG1VCIOppoXIuZ1JxgZz3HkbU2OtGNwjJWLNKhSRGQuUYiIoaXzu+hOdkzMeKimMAVzQXdn5O2Y150kl3fS2XzV/XJ5Z9eBMVZqoSkRkTlFISKGzIyVA70TgxWrOZTKkOwwejqj/09d6+3Adx9KkclpjQgRkblGISKmVi7qZUctlYh0lvk9SUpW6ozERIiYYsGpnZreKSIyJylExFTtlYhsU8ZDwOTtwKeqRBTGbqxWiBARmVMUImJq5UAve4fTpDLVp1geSmWZ34TxEDBZiZhqrYjC2I0VGhMhIjKnKETEVKFrYKrBlYdSmaZXIqZa+nrH4CiL53VFfidRERGZWQoRMVWY6TBVl8ZwOsuC7uZ8eBcqEYemGBOxY1BrRIiIzEUKETFVeyWieWMiJrszqnep7BzS9E4RkblIISKmju7vIdFhU1YiDqUyTVnyGqC/NzjvwSqrVro7OwcVIkRE5iKFiJhKJjo4ur+naiUin3cOjGUY6O1qSht6OxN0JowDY5VDxN7hcdLZvLozRETmIIWIGJtqmufweJa8w8Le5szOMDMW9nZWDRGFkKMlr0VE5h6FiBhbuai3aiXiwGjw4b6wrzkhAoKAUnidciYWmlJ3hojInKMQEWMrB3p54mCKbK78vSsKFYJmVSIK565WidgxGNxpVKtViojMPQoRMbZyUS+5vPPEwVTZ5wsf7gNNDBEDfV0MjY1XfH7n0BgLepJNDTIiIjIzFCJibKq1IoZa1Z1RbUyEZmaIiMxZChExNtVaEZOViObMzoAaxkQMjbFKgypFROYkhYgYm7ISEXYzNHtMxMFUllzej3jO3bVapYjIHKYQEWM9nQmWzu+qWonoSnbQ09m8/8yFgHKozIJTB8eyDKez6s4QEZmjFCJibuVA5WmeQyMZBno7MbOmvf7ieUFXyf6RIwdX7hgKZmaoEiEiMjcpRMTcykWVF5zaN5Jmyfzupr7+kvld4WsdGSIm1ohQiBARmZMUImJu5UAvO4bGyJcZk7B3eJyl85s3qBJgybwgpOwbTh/x3A4tNCUiMqcpRMTcqkV9jGfz7C3zIb53OM3SJlciCiFl73CZSsTQGL2diYkuDxERmVsUImJu7ZJg+uS2faNHPLdveJwlTf4AXxSef1+ZEPHovhHWLulr6pgMERGZObEIEWbWYWYfNLNtZpYys7vN7MIGzrPezEbNzM3s+Ga0tdXWL50PwLa9I4dtHx3PMpbJsXRBcysRnYkOBvo62TdyZCXkkb0jHLt0XlNfX0REZk4sQgTwceCjwFXAi4C7gG+b2YvrPM8XgAPRNm1mrVzUS2fCeKQkROw9FFQGml2JKLxGaSUim8uzfd+oQoSIyBw260OEmS0H3gdc4e6fdvfb3P1twG3AFXWc5zXA04FPNqelMyPRYaxZ3MfWvcOHbd8bVgaaPSYCYMn8bvaUjMnYMThGNu+sU4gQEZmzZn2IAM4FuoDrSrZfB5xqZsdOdQIzWwRcSRBGhiJv4Qw7dul8tu09fEzE7oPBh/qyJndnABzV38OTJTcB2xpWRtYrRIiIzFlxCBGnAGlgS8n2zeHjhhrO8XfAA+7+9SgbNlusXzaPrftGDpvmWViAqhXTK1cM9LBrKHXY6xe6V9SdISIyd8UhRCwGhty9dCGE/UXPV2RmzwfeALyzCW2bFdYtmcd4Ns/jByYXnXo8nF450MQ7eBasHOhlPJef6EIB2Lp3mP6epKZ3iojMYS0PEWb2h+HsiKm+fhzBa3UBVwOfcff76jjuEjPbZGab9uzZM91mNF3hr/3iLo3Hh8ZYMdDTkumVKxb2hq852aWxde8Ixy6br+mdIiJzWHIGXvNnwMk17Ff4RBwEBszMSqoRhQrEfip7D7AI+JyZDYTbCvelXmBmC9z9UOlB7n4NcA3Axo0bj1wKcpZZvywIEY/sHeZ5JywFgu6MlS26BfeKgUKIGONpq4PLvHXPCM9ev6Qlry8iIjOj5SHC3UeBB+o4ZDPQDRzH4eMiCmMhqlUYNgBHAzvLPPcr4G7gaXW0ZVZavqCbhb2d3L9rMg89PjTGKSv6W/L6hXEXOwaD3HdgNMPjB1Icv3x+S15fRERmxkxUIur1AyADvBb4WNH21wH3uvvWKsdeAVxbsu084P3h8Q9G18yZY2Y8dWU/mx8PlsAYHBln7/A465a0ZlDjwr5Ols7vYsvuYJppoR1PXbmwJa8vIiIzY9aHCHffbWZXAh80s0MEFYRXA+cA5xfva2a3AGvd/fjw2AcoqXqY2brw2/9x99IZH7H11BUL+cod28jk8jzwRFCROPmY1lQiAE48egEPhq/7251BiDhVIUJEZE6b9SEi9CFgGHg3QffEg8Cr3P17JfsliM97itSpqxYynstz784DPPDEQQBOOmZBy17/xKP6+ebPt5PPO7/ePsTKgV7NzBARmeNi8YHr7jng8vCr2n5n13CuazmyiyP2nnPcUszgJw/t5dH9IyyZ18WyFqxWWXDSMQsYy+S4/4mD3LFlLy857ZiWvbaIiMyMWIQImdrieV2ctmqA/7jncZ48mOKFJx/V0umVZz1lGWZw6Q2bOZTOcvaJy1r22iIiMjPisNiU1Oh1z17Dlt3DHEplOf9pK1r62kf19/CsdYv55aODHLOwh7Oesrylry8iIq2nSsQccuEZq3joyUOsGOjlrKe0vhLw6VeezudvfZg/OWMVvV2Jlr++iIi0lh25mrQU27hxo2/atGmmmyEiItISZvZLd99Yy77qzhAREZGGKESIiIhIQxQiREREpCEKESIiItIQhQgRERFpiEKEiIiINEQhQkRERBqiECEiIiINUYgQERGRhihEiIiISEMUIkRERKQhChEiIiLSEIUIERERaYju4jkFM9sDPBrhKZcCeyM8XzvSNZw+XcPp0zWcPl3DaER9Hde6+7JadlSIaDEz21TrLValPF3D6dM1nD5dw+nTNYzGTF5HdWeIiIhIQxQiREREpCEKEa13zUw3YA7QNZw+XcPp0zWcPl3DaMzYddSYCBEREWmIKhEiIiLSEIWIFjCz1Wb2HTM7YGYHzex6M1sz0+2aaWb2CjP7NzN71MzGzOxBM/tbM1tQst8iM/tHM9trZiNmdrOZnVrmfD1m9ikz2xWe704z+/3WvaPZwcx+YGZuZpeXbNd1rMLMXmxmPzGz4fD/001mdk7R87p+VZjZc83sJjPbbWaHzOxXZvbmkn1qujZm1mFmHzSzbWaWMrO7zezC1r2b5jOzVWb2+fAajIb/z64rs1/k18zMLjazB8wsHf7efXvDb8Td9dXEL6APeBi4F3gZcAHwW+B3wLyZbt8MX5u7gH8FXgucBbwHGAq3d4T7GPBTYAfwZ8B5wO0Ec6JXlZzvG+HxFwN/AFwPjAFPm+n32sJr+mfALsCBy4u26zpWv25vAzLAZ4AXAucC7wdequtX0/U7LXyPt4W/414IXB3+O3xHvdcG+D9AGngf8ILwXHngxTP9XiO8ZmcDTwLfB34YXqt1ZfaL9JqF58mH+78AuDz8+R0NvY+ZvpBz/Qt4N5ADji/adiyQBd470+2b4WuzrMy2N4T/M50T/nxB+PMLivZZCOwHPle07fRwvzcVbUsCDwI3zvR7bdH1XAQ8EX7IlYYIXcfK121d+Ev5PVX20fWrfg0/AYwD80u23wncWc+1AZaHH4YfKznXLcA9M/1eI7xmHUXfv7VciIj6moXH7ga+WrLf/yMIxJ31vg91ZzTf+cBd7r6lsMHdtwJ3EPxialvuvqfM5l+EjyvDx/OBx939tqLjDgD/weHX73yCvyS/VbRfFvgX4Fwz646w6bPVJ4F73f2bZZ7TdazszQR/iX2pyj66ftV1EbzvsZLtB5jsNq/12pwbnu+6knNdB5xqZsdG2/SZ4e75GnaL+pr9HrCszH5fB5YAz6vnPYDGRLTCKQRdGaU2Axta3JY4OCt8vD98rHb91pjZ/KL9trr7aJn9uoDjo27obGJmzyOo4vx5hV10HSt7HvAA8Kdm9jszy5rZFjMrvpa6ftVdGz5+zsxWmNmAmRXK758Jn6v12pxC8Ff1ljL7QXv93oz6mp0SPpb+W2742ipENN9iYLDM9v0E5WcJmdlK4DLgZnffFG6udv1g8hpOtd/iqNo525hZF0H/56fd/cEKu+k6VrYCOAH4FHAF8EfAj4CrzOzd4T66flW4+70EffwXADsJrsE/AG93938Jd6v12iwGhjyss1fZrx1Efc0Kj6XnbPjaJus9QKQZwr/kbiAYK/KmGW5O3Pw10EswUErq1wEsAC5y9+vDbbeGI+U/aGafm6mGxYWZnQD8G8FftG8n6Na4APiSmaXc/Rsz2T5pHoWI5hukfMWhUsJsO2bWS9C3vB44y913FD1d7foVni88rq2y3/4yz8VeOFX4QwQDs7pL+ty7zWwAOISuYzX7CCoRPyrZfhPBLIxj0PWbyicI+u5f6u6ZcNstZrYE+L9m9k1qvzaDwICZWclf1nP9GpYT9TUr/DtdRDCLq9J+NVN3RvNtZrIfqtgG4L4Wt2XWMbNO4DvARoKpSL8t2aXa9dvu7sNF+x1rZn1l9hvnyL7CuWI90EMwUGqw6AuCqV6DwKnoOlazeYrn8+j6TeVU4O6iAFHwc4IBe8up/dpsBrqB48rsB+31ezPqa1b4t176b7nha6sQ0Xw3Amea2frChrBM+tzwubZlZh0Ec6DPAV7m7neV2e1GYKWZnVV0XD/wxxx+/f4D6AReWbRfEng1cJO7p6N/B7PCbwjmepd+QRAsXkDwi0bXsbLvho/nlmw/D9jh7k+g6zeVJ4CnheNzij0bSBH8hVvrtfkBQVXjtSXneh3B7KOt0Td/1or6mt1JMJWz3H77CWYN1mem58rO9S9gHsEv8d8S9BGeD9wNPELJnOp2+wK+SLieAXBmydeqcJ8O4GfAY8CfEvyi/3H4D351yfn+heAv77cSjAr/DsEvsDNm+r3OwLUtXSdC17HytTLgVoJujbcTDKz8cngNL9L1q+kaviK8Xj8Mf8/9EXBVuO3Keq8NwQDXFPBeggGbXySoCL10pt9rE67bK4p+F74j/PmsZl2z8N94Pvy9ezbBYPY88OcNvYeZvojt8AWsIRh0dJCgf/rfKbMyWbt9AdvC/3HKfX20aL/FBIuh7AdGCRZQOb3M+XqBKwn+KkoB/wOcPdPvc4au7WEhQtdxyuvVTzCb4EmCMvE9wGt0/eq6hi8iCFZ7wt9zvwHeCSTqvTZAAvgb4FGCqYv3AK+Y6ffYhGtW6fffj5t5zQhWaH0o3O9h4J2NvgfdxVNEREQaojERIiIi0hCFCBEREWmIQoSIiIg0RCFCREREGqIQISIiIg1RiBAREZGGKESItCEz8xq+toX7Xlv4frYws8+Z2fda+HovM7Mni275LSKgdSJE2pGZnVmy6bsEK6l+tGhb2t1/bWbHAf3u/utWta+asD33A8/xyVvGN/s1Dfg1cIO7f6QVrykSBwoRIkJYafipu79uptsyFTP7PHCmuz+zxa/7TuDjwEp3T7XytUVmK3VniEhVpd0ZZrYu7O54u5n9rZk9YWaHzOw6M+szs+PN7IdmNmxmW8zsjWXOebqZ3Whmg2Y2ZmZ3mNnza2hLN8HNgv65ZPvZYZteZmZXm9l+Mxsys8+aWcLMnmlmPzWzETPbbGbnlhz/TDP7kZntC9vziJl9oeTl/xUYAP6k9qsnMrcpRIhIoz4IrADeCFxKcGfBLxF0jfwn8HKC9fu/YmYTtx42szMIbma1GLgYuJDg5lc3m9kzpnjNMwk+yP+7wvOfBUbCtnweeHe47WsE9734E4J7X1xvZkvD9swnuHFUDriI4B4QlwHJ4hO7+16CbpTzpmijSNtITr2LiEhZv3P3QpXhh2El4fXA6939OgAz20Rw59pXAJvDfT8FbAfOcffxcL8fAvcCHwZeVuU1zyS4QdE9FZ6/1d3fG37/IzN7CfAu4Pnu/tPwtXYRjP94CfBV4CRgEfDX7l583mvLnP/XYRtEBFUiRKRx/1Xy8wPh4w8LG9x9ENgNrAYws17gLODbQN7MkmaWJLgd983A70/xmiuAg4XwUWObRgoBoqSdq8PHh4Eh4Goze52ZraayPWEbRASFCBFp3GDJz+NVtveE3y8muGXxh4FMyde7gEVmVu33Ug/B7YvradNQ8YaiANIT/nwAeAHwOPAFYLuZ3WtmF5Y5/1jRexFpe+rOEJFWGgLywD8QjFM4grvnqxy/j2BMRKTc/TfAhWFVZCPBeI9/NbPT3f3eol0Xh20QERQiRKSF3H3EzP4bOB341RSBoZwHgC4zW+XuO5rQvixwl5l9mGAsx8kEYzUKjgUejPp1ReJKIUJEWu29wE8IBmP+E7ALWAqcASTc/QNVjv1J+PgsIJIQYWYvBS4B/h3YCswD/hI4BNxZtJ+Fr1s69VOkbWlMhIi0lLv/CngmQbfA54CbgP8LnMpkSKh07Dbg58AfR9ikhwnGOnyYYGDmV4As8MKSasdzCGZx/EuEry0Sa1qxUkRixcwuIggdx7j7aAtf94vAU919ykWxRNqFKhEiEjfXEcykeGerXtDMjiZYVOtDrXpNkThQiBCRWAkHP74JaFkVAlgH/G93r9rdItJu1J0hIiIiDVElQkRERBqiECEiIiINUYgQERGRhihEiIiISEMUIkRERKQhChEiIiLSkP8P8P9ftGVAXqMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 576x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#NBVAL_IGNORE_OUTPUT\n",
    "from examples.seismic import RickerSource\n",
    "\n",
    "f0 = 0.010  # Source peak frequency is 10Hz (0.010 kHz)\n",
    "src = RickerSource(name='src', grid=model.grid, f0=f0,\n",
    "                   npoint=1, time_range=time_range)\n",
    "\n",
    "# First, position source centrally in all dimensions, then set depth\n",
    "src.coordinates.data[0, :] = np.array(model.domain_size) * .5\n",
    "src.coordinates.data[0, -1] = 20.  # Depth is 20m\n",
    "\n",
    "# We can plot the time signature to see the wavelet\n",
    "src.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Similarly to our source object, we can now define our receiver geometry as a symbol of type `Receiver`. It is worth noting here that both utility classes, `RickerSource` and `Receiver` are thin wrappers around the Devito's `SparseTimeFunction` type, which encapsulates sparse point data and allows us to inject and interpolate values into and out of the computational grid. As we have already seen, both types provide a `.coordinates` property to define the position within the domain of all points encapsulated by that symbol. \n",
    "\n",
    "In this example we will position receivers at the same depth as the source, every $10m$ along the x axis. The `rec.data` property will be initialized, but left empty, as we will compute the receiver readings during the simulation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAc8AAAGDCAYAAABN4ps8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3XvcpVP9//HXew5y6iuSNE5jUhgd0IT4lqFvIaJ+DilEOXT4RuorEQZj+lJUkopxCMPXWaHIJGZEOUxUGgwyDsPIMOM8Ztwzn98f69pm27Pve1/XvY/3vt/Px+N67Huv61rrWvuae+7PXuta11qKCMzMzCy/Ie2ugJmZ2UDj4GlmZlaQg6eZmVlBDp5mZmYFOXiamZkV5OBpZmZWUMuDp6S1JF0h6QVJL0q6StLaOfMuK+lkSbMlzZf0F0kfa3adzczMyrU0eEpaHrgJ2ADYF9gHeA9ws6QVchRxDnAgMA7YCZgN3CBp4+bU2MzMbGlq5SQJkr4J/BhYPyIeztLWBR4CDo+IH/eR94PA34AvR8SvsrRhwHRgRkTs3Oz6m5mZQeu7bXcGbi8FToCImAncBuySI+/rwKVleXuAS4DtJL2l8dU1MzNbWquD50bAP6ukTwdG58g7MyJerZJ3GWC9+qtnZmZWW6uD5yrAvCrpc4GV68hb2m9mZtZ0w9pdgWaTdBBwUHo3/EOwalvrY2bWWZ4n4lUBrCct1bVX1Gy4ISK2b0DFOlqrg+c8qrcwe2tVVuZdp5e8sKQF+iYRMRGYCCCNiDfiqJmZkf15BGA+8N91lnb0IGmhtLrbdjrp3mWl0cB9OfKumz3uUpl3IfDw0lkqLP8KDOvJUc3MsB7Y7E7ncR7naVWeTq9fN+ZZ/pX8x9sbWh08rwG2kDSqlCBpJLBVtq8v1wLDgd3L8g4DPgdMjogFNc++0ouwz6R8v1zDetKx293gPM7TsDxDWMyOPMjRTGVHHmQIizumbm3P0+n169Y8K734RpJIf2Tr2QaLVgfPs4BHgasl7SJpZ+Bq4AngzNJBktaR1CNpXCktIu4hPaZyqqQDJH2c9JjKusCxuWsw4qnav1ylX6oRT8HQxc7jPA3JM4TF3MCFXMyVHMcULuZKbuBChgxd2Pa6tT1Pp9evm/OUETCszm2waGnwjIhXgG2BB4FJwEXATGDbiHi57FABQ6vU70vAr4AJwO+AtYDtI+Lu3JUY3tP3L1f5L9XwHudxnobl2WHIDDbnSd7KQoYCb2UhmzOLHbY5u+11a2ueTq/fYMiTccszv5bPbRsRj0fErhHxHxHx1oj4TEQ8WnHMoxGhiDiuIn1+RHw7IlaPiGUjYvOImFK4EsN7YM1ZsGmVmLvp3WlfxS+V8zhPvXk2edc9LM/CN+1antfZuOfZttetrXk6vX6DJY8VMjhXVXl9GMxaE+7edOl9d2+a9r0+zHmcp6F57pm9Ca+yzJt2vcpw/jZs1bbXra15Or1+gyUP7rYtYvAFz9eHwVMjYNI+0FPln7pnWNr31Iglv1zO4zwNyHP94vW5gzV4iWVYBLzEMtzBmlx/8wFtr1tb83R6/QZDnoy7bfMbfMGz1n9iePMv16IhzuM8DcmzmCFsx958nl05lm34PLuyHXuzeNEyba9b2/N0ev26OU8Ztzzza+mqKu2mFVYKFh5c+z9xybCedH/g7k2dx3mcpxV5Or1+3ZhnmZ8Rr7wggFFSnJAvZ6/2hr9GxJg6i+l4gyt4eoYhM7MKE4l4ysGzoMHUyjYzsz6Uum2tNl8nMzMDlgwYstocPM3MDHDwLMLB08zM3uCgkM/ge1TFzMysTv6SYWZmgLtti3DwNDMzwKNti/B1MjMzwC3PInzP08zMrCC3PM3MDHC3bRG+TmZmBrjbtggHTzMzA9zyLMLXyczMALc8i/CAITMzs4Lc8jQzM8DdtkW45WlmZsCSbtt6tprnkHaTdKWkxyTNlzRD0omS3lqortIRkkLSrUXyNYq/ZJiZGdCye56HAY8D3wNmAZsAxwHbSNoyIhbXKkDSKOBo4Jkm1rNPDp5mZvaGFgSFT0fEnLL3UyXNBc4HxgI35Sjjl8BFwPq0KY6529bMzFqmInCW3JW9rlErv6QvAJsCRzayXkW55WlmZkDWbVtvVOjpV66ts9f7+zpI0srAT4DDI2KupH6drBEcPM3MDAAJhtUfPFeVNK0sZWJETOz9nFoDGA/cGBHTejsuczLwIHBenbWsm4OnmZkBKXgOH1p3Mc9GxJh859OKwNWk9uqXahz7UeCLwKYREXXXsk4OnmZm1nKSlgOuBUYBW0fErBpZzgTOAWZJeluWNgwYmr2fHxELmlbhCg6eZmYGNKjbNtd5NBy4AhgDfCIi7s2RbcNs+2qVffOAbwGnNqySNTh4mpkZ0KABQ7XOIQ0hPWayLbBTRNyeM+s2VdJOBYYCBwMPN6aG+Th4mplZIlIoaq6fA7sD3wdekbRF2b5ZETFL0jrAv4DxETEeICKmLFVd6XlgWLV9zebgaWZmSWsmt90hez0q28odT5ptqBTGO3YuAgdPMzNrmYgYmeOYR0kBtNZxY+uvUf84eJqZWeJlVXLzZTIzsyUcFXLxZTIzs6Q1A4a6goOnmZkl7rbNrWNHMpmZmXUqf8cwM7PELc/cfJnMzGwJ3/PMxcHTzMwStzxz8z1PMzOzgvwdw8zMErc8c/NlMjOzJXzPMxcHTzMzS9zyzM2XyczMEgfP3DxgyMzMrCB/xzAzs8Qtz9x8mczMbAkPGMrFwdPMzBK3PHPzZTIzs8TBMzcPGDIzMyuo5cFT0lqSrpD0gqQXJV0lae0c+cZImijpAUmvSnpc0kWS1m1Fvc3Mul5pMex6tkGipQ10ScsDNwELgH2BACYAN0v6QES80kf2PYGNgNOA6cAawDHANEkbR8QTTa28mVm3c7dtbq2+TAcCo4D1I+JhAEn/AB4CvgL8uI+8P4iIOeUJkm4DZmbljmtKjc3MBhMHz1xa3W27M3B7KXACRMRM4DZgl74yVgbOLO0xYA6pFWpmZtYSrQ6eGwH/rJI+HRhdtDBJGwKrAffXWS8zM/M9z9xa3UBfBZhXJX0usHKRgiQNA84gtTzPqb9qZmaDnO955jaQL9PpwJbAjhFRLSADIOkg4KD0bqWWVMzMbEBy8Myt1ZdpHtVbmL21SKuSdBIpIO4bEZP7OjYiJgITU74Rkb+qZmaDkINnLq2+TNNJ9z0rjQbuy1OApKOA7wIHR8SkBtbNzMwsl1YPGLoG2ELSqFKCpJHAVtm+Pkk6hPRc6FERcXqT6mhmNjh5wFBurQ6eZwGPAldL2kXSzsDVwBPAmaWDJK0jqUfSuLK0PYFTgd8DN0naomwrPFLXzMwqlO551rMNEi39qBHxiqRtgZ8Ak0j/VH8EDo2Il8sOLX3/KQ/u22fp22dbuanA2CZV28xscPCAodxafpki4nFg1xrHPEr6ZyxP2w/Yr1n1MjMzBlXXaz28qoqZmVlBbqCbmVnibtvcfJnMzCxx8MzNl8nMzBIHz9x8z9PMzKwgf8cwM7MlPNo2FwdPMzNL3G2bmy+TmZklDp65+TKZmdkS7rbNxQOGzMzMCnLL08zMEnfb5ubLZGZmiYNnbr5MZmaWlNazspp8z9PMzKwgtzzNzCxxt21uvkxmZraEo0IuvkxmZpa45ZmbL5OZmSUeMJSbBwyZmZkV5JanmZkl7rbNzS1PMzNbYlidWw2SdpN0paTHJM2XNEPSiZLeWiPfGEkTJT0g6VVJj0u6SNK6/f6sdfB3DDMzS1pzz/Mw4HHge8AsYBPgOGAbSVtGxOJe8u0JbAScBkwH1gCOAaZJ2jginmh2xcs5eJqZWdKabttPR8ScsvdTJc0FzgfGAjf1ku8HFfmQdBswEzgQGNeEuvbK3bZmZtYylQEwc1f2ukaRfBHxGDCnr3zN4panmZkl7RswtHX2en+RTJI2BFYrmq8RHDzNzGyJ+u95rippWtn7iRExsbeDJa0BjAdujIhpvR1XJd8w4AxSy/Oc/la2vxw8zcwsaUzL89mIGJPrdNKKwNVAD/Clguc5HdgS2DEi5hXMWzcHTzMzazlJywHXAqOArSNiVoG8JwEHAftGxOQmVbFPDp5mZpa06J6npOHAFcAY4BMRcW+BvEcB3wUOjohJTapiTQ6eZmaWtCB4ShoCXARsC+wUEbcXyHsIMAE4KiJOb1IVc3HwNDOzJZo/ScLPgd2B7wOvSNqibN+siJglaR3gX8D4iBgPIGlP4FTg98BNFflejIj7ml7zMg6eZmaWtKbbdofs9ahsK3c8abah0lxH5XMRbJ+lb59t5aaSJlhoGQdPMzNrmYgYmeOYR0mBsjxtP2C/ZtSpPxw8zcws6dJVVbJnQj9DarFuAYwAlgOeBWaQWq6XRsSDecvswstkZmb91kWLYWePw3wbOARYFXgQuAf4IzAfWAVYF/gf4DhJU4DvRcQdtcp28DQzs6T7Wp4PA8+RRuheGhHPVDtIkoCPAXuTBiMdGhFn9VVwd10mMzPrv+4LnodExJW1DoqIIHXdTpV0LLBOrTzddZnMzMwyeQJnlTxPAU/VOs7B08zMku5refZJ0nuBDYE7IuLpInkH0WUyM7NaoosGDJWT9FNgeER8PXu/C3A5KQ6+IOnjEXF33vK8GLaZmQEQgkXD6ts62I5A+VSAJ5BmK/oQcDdpgobcHDzNzGwweBfwKLyxhuj7gO9HxD2kaf8+XKSwzv6eYGZmraOObz3W4zVgheznrYGXgLuy9y8B/1GksO69TGZmVkgIeobW2yG5uCF1aYK7ga9Lmgl8HfhDRJQqOxKYXaQwB08zMwMgJBYNqzcsLGxIXZrgGOA6YDrwIvCNsn2fYUkrNBcHTzMze8Oiod053DYibpc0kvRoyoyIeL5s97mkqftyc/A0M7OuJOk64NfANRHx74h4EVhq3tqIuKZo2YWCp6TVefNs9DMjomPb6GZmll8gFnXTzPBppqATgF9KuhP4DfCbIqun9KZm8JQ0BjgA2A5Yu2L3Qkl3ARcDF0bES/VWyMzM2iMQPV0UPCPigGzS962AXUix7ERJM0iB9NcRUeheZ0mvwTMLmqeQZpq/F/gtaSmXObx5KZfNgZOAkyT9EPhRRLzWn8qYmVl7Leqyu3nZpO+3Ztt3JL2PNEBoF+AISU8B15C6d2+OiJ485fZ1laYCZwFfi4j7+ypE0rJZRQ4nTbxwQp6Tm5lZ5+jCbtulRMQ/gX8CE7LJEj5Lil+/A14BVs5TTl/B8915J8rNWpqXApdKemeePGZmZu0UEU8CpwOnS3obaQq/XHoNnkVnmC/L9+/+5DMzs/YaDC1PAEmrActWpkfERXnL6NdUEpKGVG4F8q4l6QpJL0h6UdJVkioHIuUp5whJIenWonnNzKy6RQyta+tUklaRdKGk+aTZhGZW2XLLdWdY0nLAscDuwJpV8kWesiQtD9wELAD2zfJNAG6W9IGIeCVnfUYBRwPP5DnezMxq67bRthXOBv4LmAg8QJ1TIeUdVvULYC/gWuCSOk56IDAKWD8iHgaQ9A/gIeArwI9zlvNL4CJgfTzRg5mZ1bYt8M2I+FUjCssbeHYGDouI0+o8387A7aXACRARMyXdRhrtVDN4SvoCsCnweeCqOutjZmaZdM+za9sjzwP9GstTTd57lQuAPh9XyWkj0hDhStOB0bUyS1oZ+AlweETMbUB9zMysTLfe8wR+Tur9bIi8XzHOA/YE/lDn+VYB5lVJn0u+Z2tOJk3ee17eE0o6CDgovVspbzYzs0Gnm0fbRsTJkn4k6V7gRpaORRERuecoyBs8jyHNDTgZuKHKSYmIc/OetD8kfRT4IrBpNmNELhExkXSDGGlE7nxmZoNNQNcOGJK0PfBV0tzsG1U5JCgwwU/e4Pkh0v3K1UijlaqdNE/wnEf1FmZvLdJyZwLnALOyh1kh1X9o9n5+RCzIUQczMxt8fgL8jbSOZ8tG254BPEfqL67npNOpHvFHA/fVyLthtn21yr55wLeAU/tZLzMz6+4BQ+uQRtve04jC8l6lDYDdIuK6Os93DXCKpFER8QhAtjjpVsARNfJuUyXtVGAocDDwcJX9ZmaWUzff8yS1Ot/VqMLyBs8ZwAoNON9ZpCbz1ZKOZkkf8xOkblkAJK0D/AsYHxHjASJiSmVhkp4HhlXbZ2ZmxXVx8DwUOFfSAxGx1ILYReUNnkcAP5R0Z0Q81t+TRcQrkrYl9T1PAgT8ETg0Il4uO1SkFmW/pg80M7PiurzleSlpzM2fJb1I9dG2785bWN7geTRpsNCDkh7s5aRb5ykoIh4Hdq1xzKOkAFqrrLF5zmlmZoPebaTezobIGzwXkQYKmZlZl+rmuW0jYu9GlpcreLqFZ2Y2OHTraFtJG0XE9D727xYRV+QtL9c9RUlr1tifq8vWzMw6V+meZ5dOz/f73mKZpF1Ji43klndAzg1lExNUnvSjwG+LnNTMzKzF7gUmZ3Okv0HSZ4GLgZ8VKSxv8HwZ+J2kN628Lek/getIz2+amdkA1uUtz92Al0ixbDkASbuQltn8RUQcVqSwvMFzR+DtwOWShmQn3ZIUOH8HNPRGrJmZtUcPQ+vaOlVEvMqbY9lngcuAiRFxaNHy8g4YejabVPc24BxJE4HrSZPE71VkonYzM+tMXb6eZymWfRL4M3AFcEZEHNyfsnJfpYh4VNIOwFRgL+BaYM+IWNSfE5uZWWfptkkSJI3rZdefga2BZ8qOacySZJK+3Muua4AdgMnAvpJKZ23qkmRmZmYFHVdj/7FlPzdsSbKza+T9ZcVJHTzNzAa4bmp5AsObVXBfwXPdZp3UzMw6T7fNMNTM24q9Bs96JoA3M7OBp9sGDEl6S0QsaEY+r1piZmZv6LLnPGdKOljSW/McLGkzSVcBh9c6ttfgKelvkj6r0oig2iddU9Jpkmqe1MzMrAUOBQ4BnpZ0uaRDJG0tabSkd0saI2kPSadImkF6HHMetcf89Nk+v4C0ePXpki4D/gT8HZgDLCCtizYK2Az4NGnY7x+B0/v9Mc3MrG267VGViLgsa0nuCuwP/JClBxEJeJK03ueZEfFQnrL7uuf5Y0nnAAdkJ/0mS6+FJlIgvRr4eERMzXNSMzPrPN0WPAEioocUGC/NppjdFBgBLAs8BzwQETOLltvnneGIeAH4EfAjSWsDW1SeFLizPzdkzcys83TTaNtKEfEaaYKEuhWZYehx4PFGnNTMzGwg654xyWZmVpdue1SlmfyoipmZAa1ZkkzSbpKulPSYpPmSZkg6Mc/jJJKWlXSypNlZ3r9I+lhDPnxB/ophZmZvaMGAocNItwC/B8wCNiHNQbuNpC0jYnEfec8hLSv2HeAR4L+BGyR9JCL+1tRaV3DwNDMzoGXT8306IuaUvZ8qaS5wPjAWuKlaJkkfBL4AfDkifpWlTQWmA+OBnZtZ6UrutjUzs5apCJwld2Wva/SRdWfgddJjJ6WyeoBLgO0kvaWv80r6sqTlCla3Vw6eZmYGLBkwVM/WT1tnr/f3ccxGwMyIeLUifTqwDLBejXOcBTwl6aeSRvevmkvk/qSSRgF7AGuTnvMsFxGxf72VMTOz9mr1JAmS1iB1u94YEdP6OHQV0tR5leaW7e/LBsBBwL7ANyT9GTgDuDwiFhardc7gKekzwGWkluozpFmFylXOPGRmZgNMg2YYWlVSeRCcGBETqx0oaUXSDHU9wJfqPXFfsmn3viPpe8BuwFeBScCpks7L6plraj7I3/I8AZgC7NVLf/WA8C5mcxDHt7saZmYdozKqNSB4PhsRY2odlN1/vJY0R/rWETGrRpZ5wDpV0kstzrlV9i0lIl4HLgYulrQBqfX5beDbkqYAP4yIG2qVk/ee5yjglIEcOM3MrDNIGg5cAYwBPhUR9+bINh1YV9LyFemjgYXAwwXOv4Kkg4CLgI8B9wLHAisC10k6tlYZeYPnA8Db81bMzMwGntKjKvVstUgaQgpa2wKfiYjbc1bvWtKKKLuXlTUM+BwwOc8c65I2lvRL4CngNFJs+2hEbBwREyJic1JP68G1ysrbbXs4qV/4joh4JGceMzMbQFo0Pd/PSQHw+8ArkrYo2zcrImZJWgf4FzA+IsYDRMQ9ki4lxaLhwEzga8C6wF61TirpTuBDwBPAScDZvfSm/h4YV6u8Xq+SpFsqkt4O3C/pIZbuW46I2BozMxvQWjDadofs9ahsK3c8abYhAUNZunf0S6SgOwF4G2mN6e0j4u4c530W+Azw24joa5Dr3cB7ahXW11eMxbx5FO2MHJUzMzPrVUSMzHHMo6QAWpk+n2xwTz9OPQH4e7XAKWkF4IMR8efssZV/1Sqsr8Wwx/ajcmZmNkB142LYZf4EfAS4s8q+DbL9uT98rgFDkr4oqeqAIUmrSPpi3hOamVlnasWAoTZaqiVb5i3AoiKF5b0z/CtSxH6uyr51s/0XFDmxmZl1nm5az1PS2sDIsqRNJFXOkLccsD9pIFFuea9SXxF7BdLsEGZmNoB1Ybftl0jPb0a2/aLKMSK1Oms+nlKur9G2GwObliV9WtL7Kg5bDtgTyD2lkZmZWYtcANxKCpCTgUNYevL5BcCMopMA9dXy3IUUsSFF7MohxSXPkZq8ZmY2gHVbyzMiZpKeB0XSJ4A7I+KlRpTdV/A8FTiPFLEfAf4fcE/FMQuAf9d4ZsbMzAaIbgqe5SLij40sr69HVV4AXgCQtC4wuz/LtpiZ2cBQGm3bLSQ9COwWEf/IJvjpq6EXEbF+3rJzDRiKiMeyimxDGnW7BvAk8JeIuDnvyczMzFroDuClsp8b1kuadz3PVYDLgW1IMw/NA1ZOu3QzsEdE5FoOxszMOlOL5rZtmYjYp+znvRtZdt5VVU4DPgzsDSwXEe8gjbT9Ypb+00ZWyszM2mMRQ+vaBou8XzE+DRwZEf9XSsgWFL0oa5VOaEblzMysdbpttG05SacAq0XEUjPiSTqfNPj18Lzl5W15LqL3ZzlnUHBaIzMz6zxdPj3fZ4Ebe9l3Y7Y/t7zB82rSgqPV7An8pshJzczMWmwN4PFe9j2e7c8tb7fttcBPJP2ONHDo38A7gT2AjYBvStq2dHBE3FSkEmZm1hm6acBQheeBUcCUKvvWA14pUljeq3RF9roWSxYyLXdl9irSUOCObrubmdnSuvmeJ/BH4ChJ15ZPxSdpVeBIeu/SrSpv8NymSKFmZjbwdHnwPAa4C3hI0jXALFJX7S7A6/Q+BW1VeSdJmFqwkmZmNgB1+KCffouIRyR9mPR0yA7AKqS52X8LHJPNg5tboc7trHm7BfB24NqImJutjbYwIhYXKcvMzKyVIuIR4AuNKCvvDEMCfkha72wZ0n3NDwNzSSNxbwVOaESFzMysPbpthqHeSFqfrOUZEQ/2p4y8j6ocCXwDGA9szpsXx74W2CnvCSWtJekKSS9IelHSVdlq33nzbyjpcknPSpovaYakb+bNb2Zm1ZXueXbrDEOS9pP0JHAfqdF3v6QnJe1btKy8XzEOAMZHxImSKq/Ow8C78xQiaXngJtJSZvuSWrATgJslfSAi+hwqLGlMln9KVqcXgPcAK+b8HGZm1odOD4D9JenzwLnAVGAc8DSwOrAXcK6k1yLi0rzl5Q2eawC397JvIbBCznIOJD1ns35EPAwg6R+k2Yu+Avy4t4yShpBWBf9jRJTPBOFVXczMrJbvAhdHxF4V6edIugg4AsgdPPN22z4JvK+XfR8kW6k7h52B20uBE95Y6fs20nDhvowFNqSPAGtmZv3X5d2265MaYNVMAjYoUlje4Hk5ME7SVmVpIem9wP8Al+QsZyPgn1XSpwOja+T9z+x1WUm3S3pd0jOSTpO0XM7zm5lZLwK6eW7bl+l9Cr4R2f7c8gbP44AHgFtYMkH85cC92fuTcpazCmkt0EpzSeuD9mVE9nopMBn4BGkE8AHA//WWSdJBkqZJmvZqzkqamQ1OabRtPVsHuwH4X0kfKU/Mnv08Abi+SGF5J0mYL2ks6fmY7UiDhJ7LTnhRRPQUOWk/lQL9hRExLvt5SjaA6SRJG0bE/ZWZImIiMBFghNSwVcTNzLpNl88wdDipAXirpMeA2aQBQyOBR0j3RHPL/TUhIhaR+oUnFTlBhXlUb2H21iIt91z2+oeK9Mmklu8mwFLB08zMLCKekrQxqbfyo6S48zfgp8C5EVGo2zbvJAnLAmOAd5G6xWcDf42I14qcjHRvc6Mq6aNJz93UytsXz3BkZlanLm55kgXIU7OtLn3e85T0Fkk/Jd2TnEq633gZqen7nKRTJC1T4HzXAFtIGlV2jpHAVtm+vlxPej50u4r07bPXaQXqYWZmFbp8MeyGqtXy/C2wLWkKvutIC4aKtDTZTsC3SK3GT+U831mkmYqulnQ0qRV7AvAEcGbpIEnrAP8iTcwwHiAinpN0InCMpBdJkyWMIT3sen754y9mZlZct03PJ+khUpzJIyJi/bxl93qVJO1OWopst4j4dZVDzpa0K3CppP8XEVflqNkr2aLZPyHdOxVpjbVDK/qbRVoTtLJlPB54Cfg6cBip+/hkPK+umVlDdFm37R3kD56FKKJ6uZKuAl6LiD5noJd0MbBMROzahPo11AgpDmp3JczMOshE4KkIASwz5v2x6rRad9D6Nluj/hoRYxpRt07W1z3PTYDf5Sjjt8CmjamOmZm1S5fPMNRQfQXPd5DucdbyOLBaY6pjZmbtEohFi4fWtXUySR+QdJmkpyUtlLRplj5B0ieLlNVX8FyeNLq1loXAskVOamZmHSigp2doXVunkrQl6R7oB4Gr4E3N5CHAV4uUV2tY1Rrlj5X0Ys0iJzQzM2uDH5AGqO7M0sFyGmlpstxqBc8rcpQhmjSayczMWidCLOrpnkdVKnwI2DUiFktSxb5ngXcWKayvq/SlojUzM7OBKwXPzu16rdMCoLcVuFYHXihSWK/BMyLOL1KQmZkNcEE3B89bgUMk/aYsrdRr+mXg5iKFdW373MzMiokQPa93bfAcRwqg95CW1Axgb0k/BLYANitSWN71PM3MzAasiLgHGAs8T1qjWsChpKdFtqm2pGVf3PI0M7OMWLyoe8NCRNwFbC1peWB3b5IEAAAWVUlEQVRVYF5EvNSfsrr3KpmZWTEBdNE9T0nnAudFxC3l6RHxKvkmAeqVg6eZmSWhrgqewOeAfSU9DlwATGrUCly+52lmZkkAPapv6yzvBA4AHgWOBmZIuk3SgZJWqqdgB08zM+tKEfFyRPwqIrYBRgLHACuT1o+eLekSSTtIKhwLHTzNzGyJnjq3DhURT0TE/0bEaNKjKecC25JWBntS0ilFynPwNDOzJOja4FkuIu6MiG8AawA/Ia0M9q0iZXjAkJmZJaXg2eUkrQd8Edib1J37InBZkTIcPM3MrOtJWhnYkxQ0NyN9VfgD8D3gNxHxWpHyHDzNzCwJ4PV2V6JxJA0HdiIFzB2AZYD7gCOACyNidn/LdvA0M7MkgEXtrkRD/RtYCZgLTATOj4i/NqJgDxgyM7MlWjBgSNKakn4m6S+SXpUUkkbmzPt2ST+V9Iik+ZJmSjpd0juqHD4V2BUYERGHNCpwglueZmZW0roBQ+sBewB/Bf4EfDJPpmwR62uA95JWSbkfGA2MB8ZI+khElJYZIyI+2+B6v8HB08zMWu2WiHgngKQDyBk8gfcAWwJfiYiJWdoUSYuBX5KC6oxGV7YaB08zM0ta1PKMiMX9zLpM9vpiRfrz2WvLbkU6eJqZWdL5z3lOB24BjpH0MPAAqdt2HHB90TU56+HgaWZmSWOC56qSppW9n1jWxVqXiAhJnwImAXeV7fodsHsjzpGXg6eZmS1Rf/B8NiLGNKAmvTmLNDftV0kDhjYEjgeukPTpOrqEC3HwNDOzAUHSjsDngf+KiD9mybdIegSYDHwauLoVdfFznmZmlpRmGKpna673Z693VaTfmb1u2PQaZBw8zcwsKc0wVM/WXE9nr5tVpG+evT7Z9Bpk3G1rZmZJC0fbStot+/FD2esOkuYAcyJianZMD2lKvf2zY64Cvg9cIOkE0mjbDYBjgSeAX7em9g6eZmbWHpdXvP9F9joVGJv9PDTbAIiIFyVtARwHHA68C5gNXAscFxEvN7G+b+LgaWZmSQtbnhGh/hwTEU8A+1c5vKUcPM3MLOn8SRI6hoOnmZkt4eCZi4OnmZklbnnm5kdVzMzMCnLL08zMErc8c3PwNDOzpDTDkNXk4GlmZklphiGrycHTzMyWcLdtLh4wZGZmVpBbnmZmlnjAUG4OnmZmljh45ubgaWZmiUfb5uZ7nmZmZgW55WlmZokfVcnNwdPMzJbwPc9cHDzNzCzxgKHcHDzNzCzxgKHcPGDIzMysILc8zcws8YCh3Fre8pS0lqQrJL0g6UVJV0laO2fetSWdL+lxSfMlPShpgqQVml1vM7OuV7rnWc82SLS05SlpeeAmYAGwL+mfagJws6QPRMQrfeRdAbgRGA4cAzwOfBg4HngP8Lnm1t7MbBAYRAGwHq3utj0QGAWsHxEPA0j6B/AQ8BXgx33k3YoUJLeLiMlZ2s2SVgEOk7R8RLzavKqbmXU5DxjKrdXdtjsDt5cCJ0BEzARuA3apkXeZ7PXFivTnSZ9DjaqkmZlZX1odPDcC/lklfTowukbeG0kt1B9IGi1pRUnbAt8Ezuiry9fMzHIoDRiqZxskWt1tuwowr0r6XGDlvjJGxGuS/hO4khRsS84GvtGwGpqZDVaeJCG3AfOoiqRlgUuB1YB9SAOGNgPGkf65v9ZLvoOAgwBWaklNzcwGKAfP3FodPOdRvYXZW4u03P7AWGC9iPhXlnaLpBeAiZLOiIi/V2aKiInARIARUvS34mZmZiWtDp7TSfc9K40G7quR9/3AvLLAWXJn9rohsFTwNDOznDzaNrdWDxi6BthC0qhSgqSRpMdQrqmR92lgZUnrVaRvnr0+2aA6mpkNXh4wlEurg+dZwKPA1ZJ2kbQzcDXwBHBm6SBJ60jqkTSuLO95wEvAdZL2lbSNpO8ApwB/JT3uYmZm/eUZhnJrafDMHifZFngQmARcBMwEto2Il8sOFTC0vH4R8SiwBfA30qxE15EmXZgIfCIiFrfgI5iZdS8Hz9xaPto2Ih4Hdq1xzKNUmfQgIu4D9mhOzczMzPIZMI+qmJlZk3nAUG4OnmZmlnhJstwcPM3MbIlBdN+yHi1fz9PMzGygc8vTzMwST8+Xm4OnmZklHjCUm4OnmZklHjCUm4OnmZkl7rbNzQOGzMzMCnLL08zMlnDLMxcHTzMzSzxgKDcHTzMzSzxgKDcHTzMzSzxgKDcPGDIzMyvILU8zM0vc8szNwdPMzBIPGMrNwdPMzJbwgKFcfM/TzMysILc8zcxsiWh3BQYGtzzNzMwKcvA0M7OWkrSmpJ9J+oukVyWFpJEF8q8h6VxJT0taIGmmpBObV+OludvWzMxabT1gD+CvwJ+AT+bNmAXZ24CZwCHAv4GRWZkt4+BpZmatdktEvBNA0gEUCJ7AGcCTwDYRUXqwZmqD61eTg6eZmWVa86BnRCzuTz5J7wa2A75YFjjbwvc8zcwsU5piqJ6tqbbKXudL+kN2v3OepAskvb3ZJy/n4GlmZplSy7OejVUlTSvbDmpgBUdkr+cCDwI7AN8FdgRukNSymOZuWzMzyzRkcttnI2JMAypTTSk4TomI/85+vknSC8AlpC7d65t07qoVMTMz63TPZa9/qEifnL1u0qqKuOVpZmaZjp8ZfnqN/f0aiNQfbnmamVmmIfc8m+l24GlS92y57bPXu5pdgRK3PM3MrExrFvSUtFv244ey1x0kzQHmRMTU7Jge4PyI2B8gInokHQGcJ+kM4CrS5AjfB6YAN7Wk8jh4mplZe1xe8f4X2etUYGz289Bse0NEnC9pMWmU7ZeAucCFwJER0bJp7R08zcws07p7nhGh/h4TEZOASQ2vVAEOnmZmlmnIoyqDgoOnmZllOn60bcdw8DQzs4xbnnn5URUzM7OC3PI0M7OMu23zcvA0M7OMu23zcvA0M7OMW555OXiamVnGLc+8PGDIzMysILc8zcws427bvBw8zcysjLtt83DwNDOzjFueefmep5mZWUFueZqZWcYtz7wcPM3MLONHVfJy8DQzs4xbnnk5eJqZWcYtz7w8YMjMzKwgtzzNzCzjbtu8Wt7ylLSmpJ9J+oukVyWFpJE58w6RdKSkRyW9JunvknZtbo3NzAaLUrdtPdvg0I5u2/WAPYB5wJ8K5j0BOA44HdgBuB24XNKnGllBM7PBqdTyrGcbHNrRbXtLRLwTQNIBwCfzZJK0GnAYcFJEnJIl3yxpPeAk4LpmVNbMbPDwgKG8Wt7yjIjF/cy6HbAMcGFF+oXA+yWtW1fFzMzMchpIA4Y2AhYAD1ekT89eRwMzW1ojM7Ou4gFDeQ2k4LkK8HxEREX63LL9ZmbWb+62zWsgBc9+kXQQcFD2dsHx8M921qcDrAo82+5KtJmvga9Bia8DrL/kx9k3wHGr1lneoLieAyl4zgPeJkkVrc9Si3NulTxExERgIoCkaRExprnV7Gy+Br4G4GtQ4uuQrkHp54jYvp11GUgG0gxD04G3AO+uSB+dvd7X2uqYmdlgNZCC5+9Jd7L3qkjfG/hnRHiwkJmZtURbum0l7Zb9+KHsdQdJc4A5ETE1O6YHOD8i9geIiGck/Rg4UtJLwN3A54BtgZ1znnpioz7DAOZr4GsAvgYlvg6+Bv2ipQevtuCkUm8nnRoRY8uOOT8i9ivLNxQ4EjgQWB2YAYyPiCuaWmEzM7MybQmeZmZmA9lAuudZlaS1JF0h6QVJL0q6StLaOfMuK+lkSbMlzc8mq/9Ys+vcaP29BpLGSJoo6YFskv7HJV00EGdrquf3oKKcI7LFCm5tRj2brd7rIGlDSZdLejb7PzFD0jebWedGq/NvwtqSzs/+L8yX9KCkCZJWaHa9G8kLcDTfgA6ekpYHbgI2APYF9gHeQ5rzNs8v+zmkLuBxwE7AbOAGSRs3p8aNV+c12JM0c9NppIn2jwA2BaZJWqtplW6wBvwelMoZBRwNPNOMejZbvddB0hjgDtKo9gOATwE/AoY2q86NVs81yPbfCHwMOIb0+c8G/gc4t4nVbgYvwNFsETFgN+CbwCJgvbK0dUlTZHy7Rt4PkqbT+FJZ2jDSfdRr2v3ZWnQN3lElbR1gMelects/X7OvQUU5NwBnAlOAW9v9uVr8uzCE9LjXr9v9Odp4DT6Z/U34ZEX6SVn+5dv9+QpchyFlPx+Qfa6ROfKtRpoG9fiK9D8C/2j35+qkbUC3PEmjbG+PiDfmu430yMptwC458r4OXFqWtwe4BNhO0lsaX92m6Pc1iIg5VdIeA+YAazS4ns1Uz+8BAJK+QGp1H9mUGrZGPddhLLAh8OOm1a416rkGy2SvL1akP0/6cqFGVbLZwgtwNN1AD54bUX26veksmTyhr7wzI+LVKnmXIXV7DAT1XIOlSNqQ9O3z/jrr1Up1XQNJKwM/AQ6PiKozVQ0Q9VyH/8xel5V0u6TXJT0j6TRJyzW0ls1VzzW4EXgI+IGk0ZJWlLQtqTV7RkS80tiqdqQ8C3AYAz94rkLq0680F1i5jryl/QNBPdfgTSQNA84gtTzPqb9qLVPvNTgZeBA4r4F1aod6rsOI7PVSYDLwCeCHpC6//2tUBVug39cgIl4jfYkYQgoWL5G6K38LfKOx1exYXoAjp4E0t6013+nAlsCOEVHtD1DXkfRR4IvAplX+YAwmpS/SF0bEuOznKdmz1SdJ2jAiBlJvRGGSliV9eViNNNDocWAz0oDCHuBr7auddZqBHjznUf3bZG/fPivzrtNLXuhlovkOVM81eIOkk0irz+wbEZMbVLdWqecanElqZc+S9LYsbRgwNHs/PyIWNKymzVXPdXgue/1DRfpk0oCZTRgYXfn1XIP9Sfd+14uIf2Vpt0h6AZgo6YyI+HvDatqZ+rUAx2A00Lttp5P66CuNpvZE8dOBdbOh7ZV5F7J0n3+nqucaACDpKOC7wCERMamBdWuVeq7BhsBXSX80SttWwBbZzwOptVHv/4e+9HcASqvVcw3eD8wrC5wld2avG9ZZt4HAC3DkNNCD5zXAFtnzeQBkDwJvle3ry7XAcGD3srzDSPPlTh5ArY16rgGSDgEmAEdFxOlNqmOz1XMNtqmy/Z006GQbYCBN/VjPdbieNFBku4r00hJV0xgY6rkGTwMrS6ocLLh59vpkg+rYybwAR17tflamng1YgdRCvJc0DH1n0h++R4AVy45bh3TPYlxF/ktIrYsDgI+T/lC+Rrr/1fbP1+xrQJokYTHpD+cWFdvodn+2Vv0eVClvCgPzOc96/z8cm6X/L/BfpEkz5gPntfuzteIaACNJj6k8SJpgYRvgO1naNMqenRwIG7Bbtv2S9Jzn17L3W5cd0wOcU5HvpOzv4LdJ3di/zP5O7NTuz9RJW9sr0IBfkLWBK7Nf8JeA31DxMHD2nyKA4yrSlyM91/Z09styBzC23Z+pVdeANLo0etmmtPtzter3oEpZAzJ41nsdSM8xfjsLPguBx4DxwPB2f64WXoPRwGXAE6QvDg8CpwArt/tz9eM61Py/nb0/ryLfUNJMW4+ReiP+AezW7s/TaZsnhjczMytooN/zNDMzazkHTzMzs4IcPM3MzApy8DQzMyvIwdPMzKwgB08zM7OCHDytI0i6VNJcSatXpA+VdJekhzppaSxJIyWFpP3K0vaT9OUqx+6XHTuyhVUsnXuIpL9JOqws7bisPk2b21rSoZLuleS/MdaV/IttneJg0gPbv6hIPwz4EHBARMxvea16Nxv4CPC7srT9gKWCZ3bMR7I8rbY38C6Wvq7NdibwDtJMPWZdx8HTOkJEPAN8C/ispN0BJL0XOA44MyKmtrF6S4mIBRFxe0TMyXHsnOzYdsyXfBhwQSy96HtTZV90LsjOb9Z1HDytY0TEBaSJqU+XtCppqbA5wOG18pZ1jX5M0m8kvSzpOUk/r+zulfQuSRdIelbSAkn/kLR3xTGrSzpf0lPZMbMl/VbSatn+N3XbSpoCbA1slaVHlla121bScEkTJD0qaWH2OkHS8LJjSuf4iqTxWR2el3StpDVzXJPNSSuF1FzMWtL22TU7PevqLZ37q5JOlPS0pJckXShpeUnrSbohy/OwpGotzEuA0ZK2rHV+s4FmoK/nad3nK6Rlke4ARpEW5n6pQP4LSXOT/oIlCxmvQOpSRdIKwFTSmo/fI81hujcwSdLyETExK2cSafLw72THvJO0eEDlEnYlX8/OPTT7DJDmVu3N+cAepEnYbyUtQn5U9pm/UHHskcCfSV3CqwE/ys41to/yIa2I8hJpYvReSfoicDYwPiImZGnl555C6n4dDfyQNEn4JsBZpHlfvwb8StK0iChf2uxv2fm3z+pv1j3aPbmuN2+VG3Ai6f7nlQXy7JflOaMi/ShgEfDe7P03suPGVhx3I/AMMDR7/zJpfdPezjcyK2e/srQpVJlQvqxuI7P376P6pORHZ+kfqDjHlIrjDsvSR9S4JtcDt1VJPy7LP4zUqn+ddE+52ue7qSL9qix977K0lUmrcxxb5Vx/Ii3x1/bfK2/eGrm529Y6iqT/APYh/YH+sKS3Fizisor3l5BuT2yWvf8Y8GRETKk47kLSAJfSor93Ad+R9E1J71dZU6wBPlZ2zso6QOr+LXddxft7s9e1a5xnBKnbuzc/AY4nrZhxdi/HXF/x/oHs9YZSQkTMI33xWKtK/jlZPcy6ioOndZqTSS2ZHUldlCcWzP/vXt6vkb2uQvVRr0+X7Ye0KPo1pJbZP4AnJY1r0KMXpXNU1qOyDiVzK96XBh4tW+M8y5YdW83nSYt+39jHMfMq3i/sI71afeaTlv4z6yoOntYxJI0FDgSOjojrgQnA1woOOHlnL++fzF7nAquztNXL9hMRz0TEf0fEGsAGpLVPj2fJ/cx6lIJhZT1Wr9hfr+dIX0R683FS6/V6SSs26JyVVgGebVLZZm3j4GkdIRsRexapu/SnWfIPSIOHzpa0TM6i9qh4vydpgMsd2fupwJqStqo47gukrsf7KguMiBkR8T1Sa+t9fZx7AflaWbeU1a3cXtnrlBxl5PEAaQBSb6aTBh29h+YF0HWBGU0o16ytHDytU4wnjW49ICIWA0TE68ABwPqkgT95fErSyZI+Ieko4FjSc44PZfvPAx4CrpJ0QPaIxiTgE8AxEbFI0krZrEaHZvs/Luk0Uituch/nvg94n6TPSRojaf1qB0XEP4GLgeMkHZvVdRxpIM/FEXFvtXz9cAvwbklv7+2AiLifFEDfDdzQj3vMvZL0NuC9LPmyYNY1/KiKtZ2kMaQJEv63MnBExJ2SfgocIemyePOjENXsDfwP6fGJhaTW7BsP6kfEK5K2Jj1ycRLwVlLLaJ+IKA3YeQ24m9SFvA6p5ToD2Csiru7j3D8gBfqzgRVJrdyxvRy7H/AI6fGTo4GnsvzH1/h8RVxN+iw7kR6NqSoiZmTX5GZgsqTtGnT+HUn/Br9uUHlmHUMR0e46mNUtm6zgV8B7IuLhNlenY0g6D1gzIv6rDee+Hng2IvZp9bnNms0tT7Pudjxwv6QxETGtVSeVtDGwLbBRq85p1kq+52nWxSJiJqmLeLUWn3p10gQS7gWwruRuWzMzs4Lc8jQzMyvIwdPMzKwgB08zM7OCHDzNzMwKcvA0MzMryMHTzMysoP8Py6vbY+i+KC4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 576x432 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#NBVAL_IGNORE_OUTPUT\n",
    "from examples.seismic import Receiver\n",
    "\n",
    "# Create symbol for 101 receivers\n",
    "rec = Receiver(name='rec', grid=model.grid, npoint=101, time_range=time_range)\n",
    "\n",
    "# Prescribe even spacing for receivers along the x-axis\n",
    "rec.coordinates.data[:, 0] = np.linspace(0, model.domain_size[0], num=101)\n",
    "rec.coordinates.data[:, 1] = 20.  # Depth is 20m\n",
    "\n",
    "# We can now show the source and receivers within our domain:\n",
    "# Red dot: Source location\n",
    "# Green dots: Receiver locations (every 4th point)\n",
    "plot_velocity(model, source=src.coordinates.data,\n",
    "              receiver=rec.coordinates.data[::4, :])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Finite-difference discretization\n",
    "\n",
    "Devito is a finite-difference DSL that solves the discretized wave-equation on a Cartesian grid. The finite-difference approximation is derived from Taylor expansions of the continuous field after removing the error term.\n",
    "\n",
    "## Time discretization\n",
    "\n",
    "We only consider the second order time discretization for now. From the Taylor expansion, the second order discrete approximation of the second order time derivative is:\n",
    "\\begin{equation}\n",
    "\\begin{aligned}\n",
    " \\frac{d^2 u(x,t)}{dt^2} = \\frac{\\mathbf{u}(\\mathbf{x},\\mathbf{t+\\Delta t}) - 2 \\mathbf{u}(\\mathbf{x},\\mathbf{t}) + \\mathbf{u}(\\mathbf{x},\\mathbf{t-\\Delta t})}{\\mathbf{\\Delta t}^2} + O(\\mathbf{\\Delta t}^2).\n",
    "\\end{aligned}\n",
    "\\end{equation} \n",
    "\n",
    "where $\\mathbf{u}$ is the discrete wavefield, $\\mathbf{\\Delta t}$ is the discrete\n",
    "time-step (distance between two consecutive discrete time points) and $O(\\mathbf{\\Delta\n",
    "  t}^2)$ is the discretization error term. The discretized approximation of the\n",
    "second order time derivative is then given by dropping the error term. This derivative is represented in Devito by `u.dt2` where u is a `TimeFunction` object.\n",
    "\n",
    "## Spatial discretization \n",
    "\n",
    "We define the discrete Laplacian as the sum of the second order spatial\n",
    "derivatives in the three dimensions:\n",
    "\\begin{equation}\n",
    "\\begin{aligned}\n",
    "\\Delta \\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t})= \\sum_{j=1}^{j=\\frac{k}{2}} \\Bigg[\\alpha_j \\Bigg(&\n",
    "\\mathbf{u}(\\mathbf{x+jdx},\\mathbf{y},\\mathbf{z},\\mathbf{t})+\\mathbf{u}(\\mathbf{x-jdx},\\mathbf{y},\\mathbf{z},\\mathbf{t}) + \\\\\n",
    "&\\mathbf{u}(\\mathbf{x},\\mathbf{y+jdy},\\mathbf{z},\\mathbf{t})+\\mathbf{u}(\\mathbf{x},\\mathbf{y-jdy},\\mathbf{z}\\mathbf{t}) + \\\\\n",
    "&\\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z+jdz},\\mathbf{t})+\\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z-jdz},\\mathbf{t})\\Bigg) \\Bigg] + \\\\\n",
    "&3\\alpha_0 \\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t}).\n",
    "\\end{aligned}\n",
    "\\end{equation}\n",
    "\n",
    "This derivative is represented in Devito by `u.laplace` where u is a `TimeFunction` object.\n",
    "\n",
    "## Wave equation\n",
    "\n",
    "With the space and time discretization defined, we can fully discretize the wave-equation with the combination of time and space discretizations and obtain the following second order in time and $k^{th}$ order in space discrete stencil to update one grid point at position $\\mathbf{x}, \\mathbf{y},\\mathbf{z}$ at time $\\mathbf{t}$, i.e.\n",
    "\\begin{equation}\n",
    "\\begin{aligned}\n",
    "\\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t+\\Delta t}) = &2\\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t}) - \\mathbf{u}(\\mathbf{x},\\mathbf{y}, \\mathbf{z},\\mathbf{t-\\Delta t}) +\\\\\n",
    "& \\frac{\\mathbf{\\Delta t}^2}{\\mathbf{m(\\mathbf{x},\\mathbf{y},\\mathbf{z})}} \\Big(\\Delta \\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t}) + \\mathbf{q}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t}) \\Big). \n",
    "\\end{aligned}\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\operatorname{damp}{\\left(x,y \\right)} \\frac{\\partial}{\\partial t} u{\\left(t,x,y \\right)} - \\frac{\\partial^{2}}{\\partial x^{2}} u{\\left(t,x,y \\right)} - \\frac{\\partial^{2}}{\\partial y^{2}} u{\\left(t,x,y \\right)} + \\frac{\\frac{\\partial^{2}}{\\partial t^{2}} u{\\left(t,x,y \\right)}}{\\operatorname{vp}^{2}{\\left(x,y \\right)}}$"
      ],
      "text/plain": [
       "damp(x, y)*Derivative(u(t, x, y), t) - Derivative(u(t, x, y), (x, 2)) - Derivative(u(t, x, y), (y, 2)) + Derivative(u(t, x, y), (t, 2))/vp(x, y)**2"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# In order to represent the wavefield u and the square slowness we need symbolic objects \n",
    "# corresponding to time-space-varying field (u, TimeFunction) and \n",
    "# space-varying field (m, Function)\n",
    "from devito import TimeFunction\n",
    "\n",
    "# Define the wavefield with the size of the model and the time dimension\n",
    "u = TimeFunction(name=\"u\", grid=model.grid, time_order=2, space_order=2)\n",
    "\n",
    "# We can now write the PDE\n",
    "pde = model.m * u.dt2 - u.laplace + model.damp * u.dt\n",
    "\n",
    "# The PDE representation is as on paper\n",
    "pde"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# This discrete PDE can be solved in a time-marching way updating u(t+dt) from the previous time step\n",
    "# Devito as a shortcut for u(t+dt) which is u.forward. We can then rewrite the PDE as \n",
    "# a time marching updating equation known as a stencil using customized SymPy functions\n",
    "from devito import Eq, solve\n",
    "\n",
    "stencil = Eq(u.forward, solve(pde, u.forward))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Source injection and receiver interpolation\n",
    "\n",
    "With a numerical scheme to solve the homogenous wave equation, we need to add the source to introduce seismic waves and to implement the measurement operator, and interpolation operator. This operation is linked to the discrete scheme and needs to be done at the proper time step. The semi-discretized in time wave equation with a source reads:\n",
    "\n",
    "\\begin{equation}\n",
    "\\begin{aligned}\n",
    "\\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t+\\Delta t}) = &2\\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t}) - \\mathbf{u}(\\mathbf{x},\\mathbf{y}, \\mathbf{z},\\mathbf{t-\\Delta t}) +\\\\\n",
    "& \\frac{\\mathbf{\\Delta t}^2}{\\mathbf{m(\\mathbf{x},\\mathbf{y},\\mathbf{z})}} \\Big(\\Delta \\mathbf{u}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t}) + \\mathbf{q}(\\mathbf{x},\\mathbf{y},\\mathbf{z},\\mathbf{t}) \\Big). \n",
    "\\end{aligned}\n",
    "\\end{equation}\n",
    "\n",
    "It shows that in order to update $\\mathbf{u}$ at time $\\mathbf{t+\\Delta t}$ we have to inject the value of the source term $\\mathbf{q}$ of time $\\mathbf{t}$. In Devito, it corresponds the update of $u$ at index $t+1$ (t = time implicitly) with the source of time $t$.\n",
    "On the receiver side, the problem is either as it only requires to record the data at the given time step $t$ for the receiver at time $time=t$.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Finally we define the source injection and receiver read function to generate the corresponding code\n",
    "src_term = src.inject(field=u.forward, expr=src * dt**2 / model.m)\n",
    "\n",
    "# Create interpolation expression for receivers\n",
    "rec_term = rec.interpolate(expr=u.forward)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Devito operator and solve\n",
    "After constructing all the necessary expressions for updating the wavefield, injecting the source term and interpolating onto the receiver points, we can now create the Devito operator that will generate the C code at runtime. When creating the operator, Devito's two optimization engines will log which performance optimizations have been performed:\n",
    "* **DSE:** The Devito Symbolics Engine will attempt to reduce the number of operations required by the kernel.\n",
    "* **DLE:** The Devito Loop Engine will perform various loop-level optimizations to improve runtime performance.\n",
    "\n",
    "**Note**: The argument `subs=model.spacing_map` causes the operator to substitute values for our current grid spacing into the expressions before code generation. This reduces the number of floating point operations executed by the kernel by pre-evaluating certain coefficients."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "#NBVAL_IGNORE_OUTPUT\n",
    "from devito import Operator\n",
    "\n",
    "op = Operator([stencil] + src_term + rec_term, subs=model.spacing_map)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can execute the create operator for a number of timesteps. We specify the number of timesteps to compute with the keyword `time` and the timestep size with `dt`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Operator `Kernel` run in 0.01 s\n"
     ]
    }
   ],
   "source": [
    "#NBVAL_IGNORE_OUTPUT\n",
    "op(time=time_range.num-1, dt=model.critical_dt)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "After running our operator kernel, the data associated with the receiver symbol `rec.data` has now been populated due to the interpolation expression we inserted into the operator. This allows us the visualize the shot record:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAboAAAGDCAYAAABdgtXgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzsvWnYdUdVJnyvNwMhIJAQIgYyGoaEGSIgQQwgELtpEJEWtTX4AQFbFJqmBQWRUVFoQBo+4UUUPkDlE9GAV7fMYBM7QBgCJEyBDMzkTd5AyARJqn+cIfdZOauedapqD+ecdV/Xcz37nLNr79q1q+quNZaklBAIBAKBwKZi19AVCAQCgUCgSwTRBQKBQGCjEUQXCAQCgY1GEF0gEAgENhpBdIFAIBDYaATRBQKBQGCj0TvRicjhIvJ2EfmeiHxfRN4hIkc4yx4gIi8VkW+JyJUi8n9E5AFd1zkQCAQC64teiU5EDgTwAQB3BHAKgF8HcDsAHxSRmzgu8QYATwTwXAAPB/AtAO8Wkbt3U+NAIBAIrDukz4BxEXkqgJcDuENK6dzpd0cD+DKA30spvTxT9m4APg3g/0kp/fX0u30BnA3giymlR3Rd/0AgEAisH/pWXT4CwBkzkgOAlNJ5AE4H8EhH2R8BeBuVvQbA3wF4mIjcqH11A4FAILDu6Jvo7gTgc0u+PxvA8Y6y56WUrlhSdn8Ax9ZXLxAIBAKbhr6J7mAAe5d8fwmAgyrKzn4PBAKBQGAB+w5dga4hIqcCOBUA9ttvv3sdcsghA9coEAgExoNLL70UV1xxhQDAySefnPbs2VN1vU984hPvTimd3KRyjdA30e3FcsnNktZ02SONssD1kt0CUkq7AewGgMMOOyydeuqpvpoGAh3huuuumx/v2hWhrIFhsXv37vnxnj178PGPf7zqert27RqdNNH3KDsbE1ubxvEAznGUPXoaoqDL/hDAuTcsEgiMD9ddd938LxAIdI++ie6dAO4rIsfMvhCRowCcOP0th3cB2A/AY6jsvgB+GcB7UkpXt65sINAKKaX537XXXjv/4+8DgTGA+2TJ3xjRN9G9HsD5AE4TkUeKyCMAnAbgawBeNztJRI4UkWtE5Lmz71JKn8IktOCVIvIEEXkwJqEFRwP4ox6fIRBYGSzFXXPNNfO/kO4CY8MmEl2vNrqU0uUi8iAArwDwZgAC4P0AnpZS+gGdKgD2wQ2J+DcBvBjAiwDcAsBZAE5OKX2y67oHAq1w7bXXDl2FQGApxkxWNejd6zKldCGAR+9wzvmYkJ3+/koAT5/+BQKjhZ4smNyuvPLK+fEBBxwwP9aOKSI3GAKBQKAAGx9eEAgMAa2KvOaaa+bHl1122fz4x37sx+bH++yzz0IZ/TkQ6AMh0QUCARM8QWj15A9/+MP58cUXXzw/5rjOffddHI4s4YV0F+gLQXSBQMAETxBaomN15Te/+c358RFHXL9DFasx9fWC6AJ9IYguEAiYyBEdS3SXXnrp0u91mU2ccALjxyb2uyC6QKAClrryRz/60cJ53/ve9+bHF1544dLvb3GLWyyUYVVmqDEDgXIE0QUCFbCkuKuuumrhPLbLWd//+I//+MJvN7rR9TtPhRoz0AcivCAQCACwyY3VkJdffvlCGbbLWd8feeRiKtcb3/jG82OW7pjogvQCrRFEFwgEXOpKVkkCwAUXXLD0Wvz98ccvbsl485vffOl9Qo0Z6BJBdIHAFkIPfJbiOD6OPSu1qlLb7JZ9r8vc6la3mh+zRybH12miC+IL1CKILhAILBAdS1o51aUHugxfj+/D94+g8kBgZwTRBQJLkAsVsOxy3//+9+fHlk0uB12GbXY3u9nN5sf777+/WbdQawZqERJdILCF0FlOWF15xRVXzI/37r1+7+Dvfve7K99Hl+HrHXzwwfPjXH7M2Mg1UIPwugwENhyWFKcHvqWiLCG3HPh6t771refHN73pTefH++2330IZa/fykO4CXgTRBQJbAh7sLMEBi04nnOXk29/+dtM68PU4VRirMTnWDgibXSCwDEF0ga0GE5oVNsASHLAoxV1yySXzYya9FuDr8X1uectbzo8PPPDAhTJWNhVGSHeBHEKiCwQ2GJbqMrcTwQ9+8AP0Ab6P5Y0J2CrXILeAF0F0gcCGwRP8zQ4nwGIweGu7nAW+z21ve9v5MasxgUWPTCveLkgvYCGcUQKBDYAexJYUlyM6ViPycZew7slqTAC4yU1uMj9m0rNUtEF6AY0gukBgzaEHMUtx7HTCSZk5Pg4ALrrooo5q5wPfXyeC5h3L2VHF8sAMoguMASLyLwAeBuDFKaXntL5+EF1g42FJM4AvKTPHswH16srTTz99fnziiSeuXJ7vr+vGW/1YSaFz7RHEF+hbohORXwFwty7vEUQX2HhYdjhgUYq7+uqr58fs/KHVk3yeFxwHd7/73W/p997wBL6/rhvnx2TpjtWYufi6ILpAn0QnIgcBeAWA/wLgb7q6TxBdYKth2eWYTLSNrgQPfvCDd/z+rW9968rX1XXjevPz6FRhgYCFniW6PwXwuZTS34pIEF0gsApKQgWs+LgWnpWvec1rdvy+hOh03VhCZDUmx9txNhUdaxf2u+1Gn16XInJ/AL+BjtWWQBBdYEPh8aYEFp1OLrvssvkx2750wLgHj3rUoxY+895y1ve6zD/+4z/ueB9dN643P4+lxtSZVCJtWKAPiMj+AF4H4GUppS92fb8gusDGwCPFaWJg1R/Hx+3Zs2fl+zMx/PZv//bK5XWZf/qnf5ofe1fZXG/2yGTpjpNC61yZ4Z0ZaCDRHSIiZ9Ln3Sml3eqc3wNwYwAvrr2ZB0F0gY2BR4rTjiQs9bC6Uu8Q7sHDH/7w+bFlk8tBl+Hrvetd73Jdg+vNz8O7H1ixdsCihBfS3XaiAdHtSSmdYP0oIkcAeDaAJwC4kYhwwtYbicgtAFyWUrp26QUKEEQXWGtYrvJe1SUTH6sxS3CnO92pqnzuel6iY/DzeJ1UPO0ZpLfZ6MFGdwyAAwC8Zclvz5j+3QPAp1vdMIgusNawJmMOG2B1ZS6dV0mWkzvf+c7z41/7tV9buXwOfL1//ud/nh9/7nOfc5W3JFTLSQVYVGWydBdEF2iITwN44JLvP4gJ+b0BwLktbxhEF1gr5FJ4ebKcsKoSWNwhQGdA8eBhD3vY/JhJrwX4enwfL9Hx8/BzshqT97YDFrOpcJA5k1sEmW82upboUkqXAviQ/n7ajy5IKd3gt1oE0QXWCl6is+LgNJnpzCIe3Ote95of/8Iv/MLK5UvA9/nQhz40P/7EJz7hKs/PyW2gvUEtRxVLugOC6DYJkdQ5EBgIuZRV1o4DLMVxlhPtZOLZQ45TaQHASSedND++//3vv2P5FuD78P3POeechfN4U1gGPye3wUEHHbRwHqsyWbqzdkIAQq25aRiK6FJKnXWeILrAWiFHdJZjCU/smti0c8oy3O1ui/GsJ598sq+yHYHvz3kzAeCMM85YWoaf0yI9YHHbH490F9g8hEQXCAwAy4NSf7YSNPOx17OSJZPDDz984bfb3/72rmt0Bb6/rttHP/rR+bE1YXEb6LhC/my1rX4HEYYQGDuC6AKjhCdUAFhdXamdUSwwmTzgAQ9Y+O2II45wXaMr8P113T7zmc/Mj7/4xeUJJ7gNtETHqkx2VLGcVIAIQ9g0hEQXCPQEjzclYIcO8GTOzhe5dF4cSH3Pe95zfvzQhz7UW20XPv/5z8+PjzvuuKpr6br927/92/z461//+vyY83hyG2jnHG43VmOynTKyqWw2gugCgQ7hkeI00VlSHE/gXinuJ3/yJ+fHvJVOa1Ul57T8whe+UHUtXTeu99lnnz0/ZkmPodvG8sjkRQBLd4Bvr7sgvfVAeF0GAh3DkuJyuSrZy9CS4ixPRAC45S1vOT9mpxP2bGyB73znO/NjViny93q38BJwvT/2sY/Nj7/xjW/Mjy+++OL5sW4ba4FgSXeAvRtCSHfriSC6QKAhcjFxVtou7UzCKjmW6PiYiVKr3Y4++uj58b3vfe/5cevg7ze84Q07fv8Hf/AH1ffhevPzsLqUyUx7nVptyO2ss6lYuyHETuaBsSCILjAaWBOj10bHJGiFDWhphHfk7tLJxMpm4s1yUgJ+Hn5ObgPdTtaigttZvwN+PzlyC6wHNvG9BdEFekVuIrSkuFyuSpZOLIcLljgOO+ywhfJ3v/vd58c/8zM/s/MDOPHxj3984fP73ve+pefx97rMT/3UT1XVgZ+H4+vOO++8+bGWkK2NaLmdeW87wHZUsXZC0AjpblwIogsEGkKHCnhyVfLkCyyq15gEebByEmPtvHGf+9xnfqyzhNTgve9978Lniy66aOl5/L0uU0t0/Dz8nOykohNZczA5t6elxgRsRxV+n5roIuh8vAiiCwQKUBITx5IFO0zwhKs/cxmecH/iJ35ifnzXu951ofxP//RP7/wATnzyk5+cH1sZSnLQZfh6HO5QAn7OM8+8fk/MCy64YOE8bmtuT8t2ByxKeCzdsTdmbidzRkh3gS4QRBfoFdbO34AtxeUmWSvTCbvGH3vssfNjTsgMAIceeqin2i585CMfmR+z84cXugxfr5bo+Dm5DfQ92SOTg8mt96E/s6MKq4z1u2aiC+luPIjwgkBgBXikOD35WVIcq9C0jY6lQJbibn3rW8+PeQNT9kRsAbarffrT1+8Tye78XugyfD2+T61Kk9vgU5/61MJvLOFZzj36HfBnfm+cK9O7kzkjpLthEEQXCBTA402pP/PEysmadRwdX4+Jju1TTHraGaUWnH1kz5498+Nc7J4FXYavx/epJTpuA24bYLHdvvnNb86Pc/kxrd3LrVhIIALLx4wgukDACU/wt3Zt58nUkuL0JMuSAW8oyk4nWl1Zi/PPP39+zBKRtnfVgq/H97nHPe4xPz7qqKOq7qHbhh1VvvWtb82PObBdvwPrXeWyqbB3phVYHqQ3DILoAgEDueBvPmY7XC7LCXv28eSpHVh4MuVM/hw43ZroWI3IWU5YAmsBvh7fh+/fmujYAearX/3q/JjtcNouyu+H3xsnhdbxixFkHugTQXSBJsjtE+fZHBWwJQMuox0XOIXX7W53u/kxx8fVOjvofJSf/exn58cXXnjh/FjbrmrB1+P78P3vcpe7LJS54x3vuNI9dNtwu7Gjyne/+9358be//e2FMpb9jo91NhVrN4ScRBdE1w9CogsECLnVtyXFsT1H26SsSZKlOC0Z3OY2t5kfH3/88fNjVu/VgqUcADj33HPnxxyHlosb05LoMugyfD2+D99f121VotPgdjvrrLPmx6yu1Vv78Du13qF+155NXXN9KkivG4TXZSCQgTdvJZOeTiXFkgGX4UmNJ0hgMRickyJrm9CqYAmGj4HFoGpr2x+dU5OflQmMJ3a9zxuD78P313XjzyWhE9xu3J7czvodWJu18vvU75o/c5lIITY8NrHdg+gCK8Eb/G3Z4iyHE/2Zy7M9R3sJsrqSdx+oBUtKX/rSlxZ+42wmLM0wUWnSYrK2VHq5CYbvw/fXdeN6n3zyyeb1POD25Jyc7I0JLL5Tfm+5sBCWzJk4ud20hBthCP0giC4QIDC5tSA6SzriLWLY4QRYdDrhnQhKwF6OPLGzfQxYzPvIkwITspboWHKz1HbaBZ8lIr4P31/XjevNm7oeeeSRWBXcntzO7KQCLKpVLWlXv2t2IrLykurFAvexCDIPrILeiU5EDgfwCgAPASAA3gfgaSmlC3codwKAUwE8AMARAPYA+N8AnpNSOi9XNlCH2uBvJje21WhnFEuKYxUcS3DADZ0xasAB2l/+8pfnxzpPJT+bpXrMbU5qqS61eo+lOGvhoOvG9ebnKSE6BrezzqbCYQiWVK7ftdUPcosFluhCuusOIdFVQkQOBPABAFcDOAVAAvAiAB8UkbumlC7PFH8sgDsBeBWAswHcBsAfAjhTRO6eUvpap5UPALClOO92Lzmi40mKU3ixFMdSCrC4K3gJ2KPynHPOmR9zlhKd8oqfmydjJjdtx7J24eZn1kTHvzHp8f113bje/Dx3uMMd5sclDivczvod8G4Ie/funR+zdOclutwWQrwoCOmuOwTR1eOJAI4BcIeU0rkAICKfAfBlAE8C8PJM2T9NKS0sX0XkdADnTa/73E5qvIXwOpbkdv62VJT8vVZ3MlHw/mk8yXI6rxLoZ2NXfVbJcQwbkwywKE1Y5JYLkLZg7aGnwfXRdeN68/PwczLplUhA+h1wjB9Ld0xmmsSt/mFt+aM/e6Q7ICS8VRFel23wCABnzEgOAFJK500J65HIEJ0muel3F4jIRZhId4GO4EnhlbPR5TZOZfBExhnxOVaOCbAEOgaMVX/sNs8Eoge+pa7MOaN4pI5cLCJfj9tQl+F68/Pwc3Ib8M4OXuh3wO+H3xu/z5zXpdU/cts4eSTkQGCGvonuTgBOW/L92QAes+rFROQ4AIcCWD1VfMBESaiAV6LjMpoMOM8i25RY1cZu7iVgyQZYtGmxaz7XX0+elrrSOgYWiY6vx22bI0NrgZFLo8bPw8/JbVBCdPod8Pv5yle+Mj9mctX1tLwz2Qs1l02F30FkU2mLkOjqcTCAvUu+vwTASrteisi+AF4L4CIAb6iv2najZOfvXJYTK/M9X1tPZBw6wFvraJvQqmAbEtutAOBrX7vetMsSUI6Qraz8TG45r0sLORWcZRvV74rrzc/Dz8ltoJ17SjxX+f2wGpMTUV922WULZVitafUj3aesbCrWIgIIaa8EQXTjwqsB3A/Av08pLSNPAICInIqJt+aCg0NgETmis9J5sZrMS3Q8QekdvdnphO1IrBorAbvcs8QBLEodXGcmHU1all3O6zHoleis86ysM/o3fh5+Tm4DbhugjOj4/fB74/voPKCeBZPuU1Zb5zLSBNGtjiC6euzFcsnNkvSWQkReggl5nZJSek/u3JTSbgC7AeCwww7bvDdYAW/wt6Uqs0II9Ge+HsdPceYNYFGK4wmzBOwCz56VLGUAdkwcE5X2oGSis8gtFzDOkzG3TU6i84Z1WCTIz8ltoPN41krS/N44mF0HmbOEx4mgc0TH/Y37oZVCDLAXC0F6NoLo6nE2JnY6jeMBnLPk+xtARJ4N4JkAfiel9OaGdds6ePeJ85BbLiaOJ3227+jgb1ajlTidcN14uxnLm1KX4cnPssPpz16iy5HYDHqBYU3GuXdlOXnwc1qemcBiu7G3q9441QK/N36frD7WdbDiBXOLJ0+uTGCx3YPothd9E907AbxMRI5JKX0VAETkKAAnAnjWToVF5Hcxibt7dkrp1R3Wc2NhrdZydh9Prspc0mKe9Fmi004NtSpKdr7giZRtVdppxorH8npQWsfe9FU5eOLGdN0s25VFejpBM7cbt+dtb3vbleoOLL5P/a65H1gepTlPXk+uTP1bkN7OiPCCNng9gKcAOE1EnoNJwPgLAXwNwOtmJ4nIkQC+AuAFKaUXTL97LIBXAvgXAB8QkfvSdb+fUnJJhIHr4bX78MRoxXBprzqeVDiFF3v56WBvVpuVgNWVlvdfbisdS12pJTpLirPivJZ9XgY9sVtlvGpmS7rjNtCZVbjduD1LiI7fp37X7BzDgeV8rPuU1fe8TkBBbj4E0VUipXS5iDwIkxRgb8YkBdj7MUkBxikeBMA+AHiknzz9/uTpH+PDAE7qqNprD4+KMuey7iE6PeEyGfDK/ogjjpgfa2LzqscYPDGzvcmaSHP1tMjNS3SWNAX4JtacE1BJGUvq4e+5bYDFduP2POaYY+bH3mw03E76XbPKlCVHKxQFsPuelSsT8O91t9P324QgugaY5rR89A7nnI8JqfF3jwPwuK7qtS3wSHG5lbTlTamlDw4c5rCBkgmToYmKpQ5r8uSJUNfTIjrLww+wpThLjQmUEZ1nYs6p6jykp1W53G7cntzO2jPTI63qd839gJ1jmHi1WtXyzvRmpAnpzocgusBaIBfwbUlx3uBvKyZOk8HBBx88P2anE57wdBydB9pLkD37OM8je/XxxJ6LibPITT+blWE/Z6MrITorzVXuWh6iy+0JyO3G7cntzCQFLG54a0G/a+4H7KjCRKs3a/WEJGgPWetdeRYRO/0WWB8E0W0BVpXitLebNcHwJHDTm950oQxLcSwBlMRp8YTHAcnA4s7X7EhhSXF6IvSoK3MSnTV5tlZdeolu1Q1vc1slcXtyO+t3wO/Uu3ix+gSTq5boLPtdru96gsxDultESHSBtUBOovPEY+W83awAZ00g7IDCgeGaED24+OKLlx4Di/Fhls0wZzvz5KrU0pmViaO1RMewtvbRdjzP8+Tqye+a25PbWb8D/ux1WuF+wP2D+43uUxY5efNjenYy3/YUYuF1GRg1vCm8POpKnRHfyqrPq3cd98ZOJ7UborIEoeO+rPyUPCnyql5LHJ7MJjkHB4s0chOkFTCeO4+PPU4q+jyPGhOwM6tYtjtg8f2UeGdy/2A1pk7AzVsSsSTK/VP3XctRxSstb6O0F0QXGC28RGepK3Nbv1gxXJxSTScHPuqoo+bHTHpesDcl24fYKxBYVG9Z2e0tO1zuN+9u4TxJ5nIuWpOkJaXoz/xsXvuSda2c45FFglbeTGDx/fB79zobcf/g8vo+vJO5Jb3rvusJQ8hJuNtIdJuIILo1RssUXjwh5CQGzi5/yCGHzI91lpOSHa25bjx58ipfZzaxdhlg0mIpTkt01uSXy3JieVfmCKhEdcmfLYmuxF6X2zaJP1upufQ74PfD740JzLMfH7DYb9guCADf+c535seWdKf7rhWGYL3fnMp5W0gvJLrAqOAlOo8Ux9/n7D5sQ2GHE01smvg8OPfcc5ces/s524qAxWe1ki0zueVi4iwpTtv1chKA9X0J0VlqzRJCtaRDTXTcDzz2OmDx/fB7Y+nMmzeT+43uU3wflu6YwHLp6/gZvO96G7OpBNEFBkdtCi/LbpOT4iw1IJOe3hlCS0EesFPD3r3X5/jO7VRtqQ6tFXtOOrNsYt4sJ7mteEomRsvpxJp89WerbXLt4VHp6XfA74ffm3Za8YDro/uU5ajCZXRojKe/RwqxRQTRBUYFS6LTdhfLFduS4vREbu04wKvvEjuctsFYUhy7letJlic5S0VZks4rp84qiW+r9brk8jkHllXvk0sKbUl3+h3w+7GkuxKJX/cpLvOtb31rfsyxf7pultreOt52iS68LgODwOugkNu6xbLLeYO/OSEvO53wRMT2uhy4nnpvOE/gsCYdS0VpkVvOscQit1xMnCXdlUyEJR5/ufvw83C753ZmsIiOj3XsHb8ffm/8PvW7vs1tbjM/ttpQ9ynub7xIYpuhDjK3+rvljan7hzWuIsh8vRBEt2bw5K3MSXTWlig82bAEByyGDrD7eIkrOTsY8IofsIOFLW9KoJ7orETM3pg476RWmwKMkVNdWrCcQby5MpkYtHrQ8s7k96nfNZOWzrRiwep7HIag7YeWJGp5Y+qxY/WPTZbuQqIL9ALvbt+ezBfA4sRkBX/zYNdbqrAUx+ojvXGqBfbY45X9BRdcsHCeZ7dv7TXJXqAecstJdN6YuKEnNm+Mnge6f3G7eaQ7wPbO5Pep3zX3g8MOO2x+rMM/GJbanAlV2wW5PtYiL5cUms+zbKabJt0F0QV6h1d1aalo9GePFKf3hWM1Ex97wVIcBxvzBAUsrsaZxK3Ey4AtxXnyUerPFrm1mKy6mjxqSU9LyFZQtSUNAbZKkN+nftfcD5i07njHO+5YZ8Duk2y7Axbtd9Y4yJE4P4+l2t60bCpBdIFekJPorDCCnIOB5VFpZQ/hnQeARW+3khRe1n5jLAkAi5MPTxY8wXizlHhCAPR9atHlBFFbT48aFLDb0PsOLDLR79rqE15wP+T+qfsu92tLIst5H3vGW26MrjvpbQqC6EYCb0ycJ8tJzobCYGmI8w3qLCe8YuZJJQd2FuDVu+VNqevJKsZcwPeqUlzOscS7FY73t5YomTw92Vhy97Gku5zUY6nQ9bvmfsD9g/tazjOT+yGX4esCi+pTJlsOOM/t3OHxytVjdN2DzEOiC3SGkuBvjyu4vh4PXF4VH3roofNjrZ7kwHALeoJgdSUfc3aL3G7fVvB3bhsWTwyYN5NIDtZ5fUl0tZNnzmvUikVkotPvwMo+ktslgfsB94+cDdjamJf7p+677KjCZGsFmXsD6C0TgP68bkQX4QWB5vCQmzedV47oGDxJ8Z5xPFloiY6lPQt6Jc2u5fwbe+XpevJExpIbO5zk4uBqky0zcqop67xcfFstrGco2azVEx6hr83trN8bvx+rT+qFEPcD7h/cbzRpWd6ZOW0E92sOQ2DiZUkvR3R8zIuAnNZlHaW7ILpAU3i2CvEOPCsmDlgclCzFcZwSe755JDhgUf2jveouvPDC+TFPMNYOA8DqoQKA7SBQq4YsIbBaFad3Q9Ucah1qLBLkds7F3llhCNr2ZeXO5H6j+xRrHSxbse673K9Xle6AsiBzK1F3EN1wCKLrEbkO5EnTBdiTSi6FFxMFr34tKU57XVr45je/OT/WCXj5NysmTquiLFtcyU4CFrzOPd4ynu93+m2GEjd1/b31DCUE6EmvBtiSeM4TmH/j/pHrUyzh3f72t19aZ913uV9zf+cwBEu6A2zPZst+qT9bC4dN89QcO4LoBoRnwvRugmrlOAQWJyZWM7GHGq+Qc4POCg7Wu0HzhGHtSq5d2y3Vo7WL9051naFkher1qvPeZ4hVcsu2yZGj573pd21JR9xvdJ+y1N65/sH9mvs7jwMrKBzwaVr0GPX0jzETW9d9VUR+CcCvADgBwKEALgTwDgB/nFK6LFe2FEF0HaMkVCAX/G39ZgV/A4seapzlhI393hReHKfEqiVto+PkvjyRWQmiAVtdmZswV1VX5lSAHlXyKuetCm+uy5JJ0iPp5e5jqTEB21ElZ6OziI77je5T3N9YOst5Z3K/5v7O44B3QmCVKmAHmXvHqKVx8L7rIdDDouwZmJDbHwD4OoB7AHgegAeKyP1SSs2N3UF0PcK7q4B3Y0wrETOvVgHb6YSPdfwRg/MH8uTD9hTOcQgYbfEQAAAgAElEQVQsqoO4brW7fZd4Teay03tIyyvR9eWMkpMMLKKyjnWdPepffU/LMYPfpyY6/sxluN/oPsX9jdWYTGa6T3G/tvo+hyCw3RmwE5/nxqi1MBtztp0ZevK6/A8ppYvo84dF5BIAbwJwEoAPtL5hEF0H8E6YtUTHYALRsW68euUB7pXi2IjPk00uJs4y1ns3Qc3tHuBBSXBvyXvrC16JrOV9vNf2OK3od205rXC/ycXecT/kPn300Ueb9eT+zmUsJxXAl0Islx+Tx7V3J3PGEATYNdEpkpvh49P/q6deciCIrmPk1GYWueUGEV/DIhAdDsAqG/Zcs4K/cxtrWlnjdZyUtdt3zoPSioMr8Ub0elB6VI9DEJsXOYnMkvy84QXe5+brWU4r+l1zP2AC4X6k+xT3N+6H7FmpnVG4j/MxjwMeH3r3dGsLIG6b3GLUo8YE2qecW0P87PT/57u4eBBdI3ikAS/R5Xb7ZvBEwoNYS2o8kHkisCZFDuYFFicV9orj1a+uJ9fNChvIpZIqkVQ87b6JrtOrItcna7cdsqS7nIettdN9bt877ofcP3V4AasuuW48Dnh86L7P92TpzhuSYGVT8e57x+iLAPseIyJyGwAvAPC+lNKZXdwjiK5jeHcvto5zbsgWmegYI/6sV9Yz8CDmVSyw6O1meVNq8KC2JLWSdFw5lDiGlAzq1urCGpTE+JXYOUvi8CxJT3+2yCAXQ2p5Z+q+y/2ayZXHAY8PPXZ4XHE9c3GrJeO6RGXcFRoQ3SEiwoS1O6W0e9mJInJTAKcBuAbAb9be2EIQXSN4JAhvTJzlTQksTgTWjgOslgEWJTztqDIDq2y0txuvntlDjVfi3uBvyw637BoeeGxxtQO3r6TQXuRCSWqulUPJRGyFGgC2/S4Xe8f9jfsh90/dd9k+zd6ZPA54fOixw44qrNK3xqv+bO1+kPMettp6jSS6PSmlE3Y6SURuDOBdAI4B8LMppa/vUKQYQXSF8AZ/W2ED+rO1z5zu3DxB3PzmN58f5wYrn8cDh1e/HDbAhn5gUZ2jA2pn0Kopz95wOYnOQmtX/xY7gbcqU1Ln1momz/W8wc45ic6y31m2O2BROuN+yP1T911WZfI+i6zS5PGhxw5fm8mV7Yd6AetxMtNzgWWT3tQgcxHZD8DbMYmle0hK6bNd3i+IrhDeCdeb5cRSheiVH69EOWyAV67aIM8rZr4nS3G8XxivkIHF2CZrhapj4jy7Cnglg5bSWQ7e7COrlq+9P+BXKZYEs6+Kkhgw/a4t+12O6CxJifun7rvcr3mMWN6/euxwGR4vrC7VoRNWPlrLG1N/HjLero/wAhHZBeCtAB4E4OEppTM6vSGC6FZCrgNYElkuM7qlorQ8FoHFVaklxfEKFVhcLbL6hd2qWeWjPc8sKc5aievfcpknLNTGtDFqiarFDtK1El2XUpynDiX39KYX80h3wCJpcGwn90/dd7lfW9IdO3LpscPjitWYTK46yJylPSubip4LLEcVnj/0orerTCs9OKO8BsBjALwYwOUicl/67etdqDCD6ArhTfeU24nAsi9ZkwCwOPjZRmelNNL34QmCVZd8zOcAdu6+XP5DT6CsVxopQUtya0F0NXUphTUR9iUleAPbLU/NXF5TLs/9U/ddTx9nBxR9Tx5XPN4sJxX9PB7HFP3Zu9joSnXZA9H9/PT/s6d/jOdjkiWlKYLodkBJELGlkvQGf/PA4cEFLKorWaLjFaq2l/Gg5mwTrNZhW4SOo+NnszZEzcXEebeYYbR0LPESVa3qMgePJFsSt+Y9r6UXagtbokV0udg7j9OK7rvcr7m/Wzsh6HtaGhQehy2CzC1vZG6bvoLMewgYP6rTGyxBEN0KyMUfeXYfyMXEcSfOBX+zzYAHHg9WXTceiOx0wgPf8qbUdfMGf6+6qwDQHaG1ls5KPB09167NAOO9T1cEmDuvZJcE3aeYHLiP8rjSfZf7Nfd3VmNatm7A3taKz2M1JrCourRCcHI2OssMknvX6+yY0geC6JagJPjb8qD0Ep0V16OJzpLiWOWj1TdsW2BjPUt3rNbRz2Y5C3hDBfpS79USHWOIkIKS63ozm4x5IrTU4bmQBCuPpnYM4X7N/Z3HQW6fO+77lnTH4wtYdFRh4rXIWX+27PrencwZtfGom4Iguh1QYm+zjM65LBQ8cK2wAWDRK4zVmlxPvb0JS3F8bElxekBZUpzlTQl0t32OhidxcQ61mUCGjqPzSkqMHAF6pL2S91YSkqD7lLXgysW0WbF3PA44MwpLd8Di+LPiVvUYZQnPSrCQWyhbWxB5dzKv6ZN9eF0OgSA6lMVmebfPyUlxPHAtKU4PIiZBtmewXUCvMHlQ86qWM7XzoMnlJfQmXh5CiitRD66q6hsDvGEIHlf/nF1wCGkgl4fTShid2yWBiY77O48DHh86js6y3+UWo+z5yTZDrov21PTE2+n28OxkrtHXAnRsCKLbATkpbtXzvAOXPb20i7X28JrBGtDAovqG7Qc8iHLBvZbrc0nWhtpA7pLtTbqyqbUo40VXweS6bay+W0uuJWX095YtL5dDkq/B/Z3HAY8PPXZ4XFkhM3qM8vgt2ZGjZM7JLRBWRRDdlsB60bksJ5YB2Rv8zfp/VotoVYrlecZqGY6PA2yPSktFUhIT1+UkX6uSbB0eMISEZ92zdaiAR9orCWwvQY7oSmLvuL/zOODxoccOS2tWWjs9Rnn88rhkEtWSp2XiyMXhWp6aXcXXrTO2luhyon5J8Lelf+eOplWCnIaIBwsfa+M4X88ytLMqBgAuvvji+TE7qvC1LPuH/q3Em9KLWvd+D7mN2fZWgtpsKl63/1oVZwsCtCS6XN+1HEB4HPD40GPH2taKj/UYtcYy285zsarWFkBeTQujZJEXEt2GwpvxJOdBaREnd0adUJn1/DwgWLrLrVDZ6M0rUb0zsyXFWeoXTcjWgOoy5ZU3lZSnzKaRmxceNWYuaN+rDrMWdjnJojZ8xFJd5uzLVhwbjw89dnhcsaTGEl0uexGPa0u6A3xJ3b2emrX5MYPoNgCzl5jrQJYUl9snzso4boUNAItOJzwgWNLTHZOlOHY64QHJK1Rg0fDtkeJy+8QFua0nSjKjeDOrWJKfV3IsgeWdmUsubpEJjw89dnhceaQ7wNbU8HjXntEs4Vk2uty+d7ntrxg7tXt4XW4YSqS4XKfzSHE6jx4PAv6NyVEHnHqkOL0n16pSXG4Lkb68KUvc/rvMbLJJKHEs8UpnHtLTZUpgSZve2DsrY5EeOzyueLzx2GXpTt/T0troIHMOQ7AI2bsgt6Q7YHvtd1tHdMskOm94gXdDVCsHpFZxWCm0+FraaM3eYqz+8GRj0HXzpB3S9alFa+ePbRqsQ8IrnXXlHeq9Vs4F3zrOpekrGW88lq1993Le1FZS51zawZIk6C3jJMeOrSW6nAuv5XSSIzor+NtSYwCLunweBHwfnbuPV5h8nNs2xEqrZCWmbU104UyyPujKrud1Zqm13em+ay00rVyZeuzwuLLGns5eZNnyLNsdsJimzwoB8iapyKkxPR7MQXRrDtY/lxh5rSz+gJ2ImTu3HhCs2+cO6bUZcGAqqz70gPBkNvF6U7YkEC+hBrkNg5Z2vRbxehYs1Sng887MBXJbWwDxONQB4zz++Z483vVcwPMES45ct9y+d16nuWWOKjmpb1OwVUQHXD/AcvY2Sy+eG1A8cNjphDswHwOLqz2+ds4LjB1QeLVpbfCo62bZ4nLelF2RWwuVZJBbP2hJeoBf2lu1bhoe70weHzmHDx5vPA71GLViYi3pTn/2SHeAL/RJS3Qer+kgug3AMtVlzshrGdQ1mVhSXC4mzpLi2A1ZDyKW8Hi1yZ0zR3RWfsqSoGwvvPF2nvuuK7Ftklu317HES45dkp7HOzOXQozHJY83Hod6jFr5MXMe2Na2P3xPvTMD19XrqWl5h286tpbovEZehuVwAtgDJ7dBI4M7MXduvbs3Dzyr0+Z2EvBugtoV1pWoPOjSlrkuJOhFS2cU730sVWpuxwSLhC0C1J89KcSA1c0L+nkYOV8CnveWhSFEeMEGIKXkUl1aHYU7pw7+tvTv7IyS21+LV3G8QmS7ALCov+d65zas9OTb61KKW3fb29D3X6UOQ0xSQziwlNTNWgzy+NASneW0wuNQj1HLUSW3UTHPE1a8nSZUKwA+N7cxuc3OCxvdBmL2cnMGW0tdyasrvfM3x8vwsT6PYalCWP+v4214JWkFpud2++5qn7hNcCwZA6HVYmi7yxAOLCV1s0KA9NhhouN5gsehHqM8flmNyYthJjZgcZ6w5hLtgc118OTK1J+t9gyiW3OklOYv2vPCgcVBwCsynQ2B7XL8G5OjjrexjNu8QtSpgriu1qo0l9mkpRTX0va2ynktsQnk5kFfqsKSe/Zly/NId3rs8LhiNSQTix6jPH55XLMdLpdPlucPnld0NhUrli+3792yBX5IdBuI2cvNrXp4EFhSnM5ywp9ZrcmDizsmYEtxrMbUBmgrXi+n/++D3FoHgneFbSG2HIYmPe99SzKwlNTH2hUBsL0zeSzrMcrj15LutETHn61sSnrO4SwuLN3ltFVWMPqmY6uILqU0X4nlBoe1ZQ13Rt3p+DdLiuMBACwOAiY9a0NUwFaz5Hb7rvWuapmOqy8Mff91wdgcYEoysNSSXs6Ri8eVZcvTmhoev9ZiVocXWA4ouTnHst/lAuCXBZnrbCkh0TWAiBwO4BUAHgJAALwPwNNSSheueJ1nAfgTAKenlO6/aj1yMXGefa9ydjAekNzp9NYc3Dl5hZhbaXk2nBwihdeY1ZOB8aGlA0tfKcSs8aaJjscvj2se73ou4GtYnto5JzMrs1Eu16WFILpKiMiBAD4A4GoApwBIAF4E4IMicteU0uW58nSdYwA8B8B3dzqXwV6X/DK1m60V78IrKm2j4zLc0Vm9oLOcWFJcbg8qq3N7dy/2oK995lojSLQthggBqL1PrQOLLs/jyiIdTXSWdyaPdz0X8Nxi3UfPOZajCpOoxxklbHTt8UQAxwC4Q0rpXAAQkc8A+DKAJwF4ufM6fwHgrQDugBWfYfZyczFxltMJH2tvSh4QrLO39PWAvSEqI2cc5+Nab8out8LpEmOowzZgCLuetw65uq2q1tT9ycoTm1MPsr3M2uBVzwXsqGIFmes5x5qbrBRiQD7H7yajb6J7BIAzZiQHACml80TkdACPhIPoRORXAdwTwK8AeMeqFZgNBO7A2gOKpTirM+ky3GlYimMPLL2K4w7JXly5jSQ9MXE5eM4bc2quILbhUeJYUnLtFqnGciaKGXITfknsHY9lyztTzwU8T/A8k9v9xJqbWLrTuTtzTngzhERXjzsBOG3J92cDeMxOhUXkIEzse7+XUrqkxtMvt/O3R4rTUqAlxXGH5tRewA074Qw5vXxJTFwf5NYlAQW5jRtdSXstUo156paT+rqKvdNzAc8THulO/2bNWdrT20obxgiiq8fBAPYu+f4SAAct+V7jpQC+BOCN3huKyKkATgUmktqss+aCv63ATu5oupNYtjgrNyXgk+J05/ZIcbXxbUNLbV1fO9AdulRx9uHAkrPxeaQ7wN44lY/1XGARnSXd6fvyPMVl9Eayy5zedLsE0Q0IEfkZAL8B4J5phTeRUtoNYDcAHHrooWn2UpnocnYw/o07uhb7WaLjzpTboNFSq+TIbJsSsQYCY0XOs9Eav1bKLsCeM3he8cb75pJH8GdLo7T1RCci+2NiHzsMwI0B7AHwxZTS+c5L7MVyyc2S9BivA/AGAF8XkVkQyr4A9pl+vjKldLVZelL/eSdgdaX2ZuLVETumMLRKgFURvDrjFZXW5fNKyuqcueBvz/caIcUF+sLQ0l3uvNpsLLkgc4t0cg4slkYot/sJf+Z5ypLugOUJp7dh3O1IdCKyD4BHAXgCgJ8FsD8m8W8zJBH5BoC/BfB6djRZgrMxsdNpHA/gnB2qctz078lLftsL4L8AeGXuArt27ZobdLmT6CwF/Bt3Yl6R6bxzTHR8nNsQ1bMKy2Usb52Cq2+MtV6B9ujSgaUWtaSX28nEclrRRMfzhDWX5ILMeZ7KzW1MqDOP0JJYu3VDluhE5JcwCco+HMC7MYld+xSAiwBciYkkdjSA+2BChk8XkTcCeE5K6TtLLvlOAC8TkWNSSl+d3uMoACcCeNYOdX3gku9eCWAfAL8DIEewACYvdCbJcQfQKyVLiuPOqBO5eqQ43aEsyc1yOPGitXTWR8B5YLvQVcB3bU7NEkLWY9Qayzze9SaqPE9Y0p3elZx9CywtlJ7b+LzZfLb1RAfgVQD+DMAbU0qXGud8DMDbMCG5+wB4JibOHy9ccu7rATwFwGki8hxMAsZfCOBrmKgmAQAiciSArwB4QUrpBQCQUvqQvpiIXApg32W/LcOuXbvmL55fuF71sI7dcg/WXlMeKS4XKmDFxC17hhoM4UEZ5BbIoSvSy12v9p45T01P7J220Vn7UVrSHbAo4TG58fyl5zaew2bHQXTAMSml5RbLJUgpfRTAL4rIAcbvl4vIgzAJEXgzJirQ92OSAoxTgAsmklpTz4t99tln/uJ5paPjU3gQcMAnS3E6DoZXYVZmk5zTi2XAbqGeDHILrAta2/VWleJqwxM0SmLveP7IZVZiCY/ns9zctkytyXPUVnpdrkJy3nLTnJaP3qH8+Vi0A1rnnbRKvURkTjY5CYpXaLzSYtLL5aqz0ovp++TyU3oQZBII9IcSErZ2MtdzgeXRncuVy595nmJHu1wKwdlcuA3ziNvrUkRuD+AWKaWPTT/fGMBzAdwZwLtTSq/uportwBIdi/q54G92OmHVQW4TRIZlmAZsW1xJOq5Vz1nlvL6vFQgA7R1YWtrovGnHGLkgc4vQeF7JOcCxdJcLMl/mnannv62T6BReDeDTmNjkAODFmNjbPgvgFSKSUkqvaVy/pmAbHYv0WsduOZ1wx9KbLfIqzJMIVp9nYcwEMua6BTYPY8i3uQxecswFmfPimkmP5xU95/B85A0y53lvNheGjW4RdwPwGgAQkV2YBG8/M6X0ChH5I0wcUEZPdLMVDb9cnfiU93mypDhdhsFSHHes3D5xtZujjiEOLhBYF7T01MyVsaQ7PRfwPMH2Oks9CdjaJnZS0ekNmWCXzYVAEN3NAcysoffAJPD77dPPHwLwjHbV6gYcMM5SnE7Hw7kq+ZjP01kKuONaWVdyRLcu5BZEGRgDxhqekCvj3eDV2rhZzzk8H1lzlg4YX6Zh2oYxvQrRfQfAsQA+AuChAL6SUvra9LebArjGKjgWiMj8pTLRafsaZz3h33TsC8PaoLGWzAKBwOZBzwUWCfK8oomO5yOep6z5C1ie3WXMwfytsArRvRPAn4jInQE8DhT3BuAuAL7asF6dg1UC2sjLKyJrQ9Rc8LcVCK47VFdSXIQKBLYJ6xKHx2X0fazA8lyQubXBK89fem7juDqdzB7Y0vAChWcBOADAwzAhvRfTb48A8J6G9eoEKaV5Z+H4FB2IyTY6JkTuAN5E0DmXYg+GIJkgtsC6YgjppEU2FYYnVyawGIvH85TlYwAsqjJn144dxgkppcsx2SF82W/3a1ajDnHdddfNO4Sl0wZsKc7ajwqwbXF9xcfVklOQW2ATUSvt1UpuJfF2lr1fzzlWvJ0l3enPM4lOe51vNdFtAq699tq5JJcT7z1SnM44YG3n06XDSJDTemKsbvKBPLokPctpJTfnMNFZ0l3OLDOT7qydxjcJOyV1fieAP0opfcpzsWnqr/8M4IqU0msb1K8prr322rlYzy9cb07InYY7YC742zIgj9leFkRZjyEk6SDHMgzhqWmVyUGn5JohF2RubQGk5zae925+85sDuCHRbWL/2kmiOx/AGSLyaQBvxcTj8jMppblVVEQOA3BvAP8BwC8C+CaA3+yktpVIKc29kNgbSeu+GRaB5bbPCWwWxvZuc/XZxElqm8HvOpdC0EpEr+c2nvdmx1tvo0sp/a6I/DmApwF4HiaxdElEvg/gagC3wPX7031set5bUkqjlIWvueaaeaYTb/C3lc1A68u7kuLGLBFuMta1nUIt6sMQnpol18v5BVhJopnockHms7mQz99ar8uU0lcA/I6I/FcAP43J3nOHYeKBeTGALwD415TSBV1WtAWuvfbauejOBttc8DfrxXMbolpOJ7UG8BZY10m7D2xy21jPtokTWQ2G2DHBWz4n0Vn2OyYuHZKwzFElVJeElNIPAXx4+reWuPbaa+crGr1NBsMTE7fuu3tvK+J9hNTXJ1r2Nz3nlMTe8bw3mwu33hll03DdddfN0+bwy82FClhSXI7ohp48YjJfRLSHjW3IirEKxtYeuSBza69LPtYLerbZzebCCC/YMKSUltrjcnnnSmLihlgxx2S+iGiPMoS0t4gxj2VP7J2e25joZnPh1jujbBpSSnNJzrsJohUT12WqIC9iMg8E+sPQYzl3jjVn5bzDZ3NhEN0GYvYSc95MJbY4D7zE5A0sDSwi2qYtQrqzURII3iV4nuL5K+cdnktSv2nYOqKbdTzLkAuUhQq09HCLCdtGtM0wiNg9G2Pwkra8M/XcxvPeMoe8rQ0v2CSIyLwTePeJs6A7g9U5Y4KoR5DbuBGS3+roMiOOd9+7WcD4QImwDwfwCgAPwSQO+30AnpZSurCL+61EdCJyDwB/COABmASL3zul9EkR+WNMYun+pYM6NgNvvJqLietrYo0JwkaQ23oi+rSNIeaVXOydtfFq1+9NRA4E8AFMko6cAiABeBGAD4rIXacbCDSFm+hE5P6YsO5XAfwNgKfQz9cBeDKA0RPd7EXnvCn5RXucTFrVbad7BgKB9UKX5GbNU3ys57Zl3pkDSHRPBHAMgDuklM6d1uEzAL4M4EkAXt76hqt4WrwEwLsB3AnA09VvnwRwz1aV6goz1eU+++yDXbt2zf9kuvP47I8x01n3ST5WXTYd2/rcm4ptfZ9DPLc1T+m5jee92Vw4wPt5BIAzZiQ3rf95AE4H8MgubriK6vKeAH4xpZRERM/6ewDcql21uoGIzI2xJZugcieqdVLR1/OW3yRpb9smwG1G9GN/GW/blLThMm/zAd7NnQCctuT7swE8posbrkJ0VwE40PjtJwB8z/htNGCiy3U0zhRgxdG1cCkuUVeuu4ozyC0AbGc/rnVAybUT/6YznVh1WDYX9qS9OhjA3iXfXwLgoC5uuArRfQTA00SEmXjWIo/HxLg4esxWNN5Ox50mFzBei20kvUBgndAXuXnhJTdGzlFl2XULcYiInEmfd6eUdtdetAarEN0fYqJDPQvA2zEhuVNE5OUA7gXgp9pXry1mOmoNb6hAbRkvNo30QooL5LBpfber/l7SNrkys7mwA9XlnpTSCZnf92K55GZJetVwO6OklM7CJKzgOwCejUnsw8zz8mdTSl9sX7328BiJLcNuzjFlCKcVCznnmkAgYGNsY6dkzvHORQM+49mY2Ok0jgdwThc3XCmOLqX0SQAPFpEDMGHfS1NKV3RRsa5Ru4rMSXQtQxLW1Yg/hkkisH4Ys3RnoUVf99jiWu+V56lLR3gngJeJyDEppa9O63UUgBMBPKuLGxZlRkkpXQXgm43rMnqwHlzrt0s8MktQq9ZktOzQQWyB1uhykVfbX/vq77ln9trlWt6zEV6PiTbwNBF5DiZmsBcC+BqA13Vxw1UzoxwH4JcAHI7JDuOMlFI6pVXFusJML205mXjhlejGTAC1q+cxP1tg87At/bVWovM+p+Wv0DXRpZQuF5EHYZIC7M2YmMHej0kKsB9kCxdilcwovwHgrzBh3+8C0BlB10PP0BEscusys3lLNY83rmddJovAZmOIhWXLa3u3xlkX9e2qmOa0fHRf91vV6/I0AI9PKV3aUX06BRtfeTVTogLQZViVWWujG5vXZZBbYMwYG7mVjEWPja6FqnJZXPC6+gGsglWI7tYAnryuJNc1ulpVbkMnDATWGX3FyvWFTZxjViG60wEch4kudW2xrFNqXXXtyqlWjRmqwkBgfeFRQ3qznNRCz22epBfbTnRPAfAOEbkYwHuwJLAvpdSNG1BDeIijVq3JaG2jKwlj2MSOGwh0gRY5bD3w2uhKYKUt1J+D6Jbj6wA+BeAtxu9pxesNgmUvt8WLtfJj5u7dly1vHWOTAoG+0CW5lUhxJam9vOgyjeGYsQoxvR7ALwP4JwBfwA29LtcCy15u7oV3Kd0NnWosEAjUo3VqrhJYUtyqEt1Ysju1xipE90gA/y2l9OddVWYolMSdeEkvZ2/zdKiQ7gKB9hhCimttl/OQW4sY4U3AKkR3OTrKQ9YXrNxurV+spcYM6SwQCADDZDwBlpPgNtj0VyG6vwbwqwDe21FdBkMJ4bTw1Cyx0YVHZiAwPFrGx3nh9aDcVjtcDqsQ3QUAfkVE3gvgX7Dc6/KvWlWsa7QW6T1qTa9EV+upWeL0somruECAUTJ2PCjxoPSON8uxTcNLdCNJ6tw7ViG6v5j+PxLAg5f8njBJEbZ2aCFB8Xnezllro/PWLQgtELDhJbdaKa6W3LykFTa6G2IVoju6s1qMGEN3lC5zZVrX2sSOHthOdKW6a6mSzKHGg7IEW+91mVK6oMuK9I3cxN6XeN/S6zJsd4FAe9SGDnRJiCV2Oc+CdquJbtPQ2u2/JO2XB6235ojMKoFNxtjCBkrgUUmWXmtbx3KW6ETkqwAelVI6S0TOQ34rnpRS+smmtWsMFstrB0ROCmxJbi1UpyHhBQI+9KWp8cIrndXME12mIxsLdpLoPgzg+3S8eS2wBCVSjwfegNFa0owg88A2YWgprkUgeMtQgdZ5ODcBWaJLKf0mHT+u89r0gGUSXU4Cau2RNQTCfhcI2Bjz2GW0TtAcNropWHXZU316Q+vJv6WE15fO34MiIq4AACAASURBVDonpLvA2NGVdJNT43kluhIbf18B3zvVbVu9Lo8CcKMe6tEbVrXRWVhXI29Id4HA+oxXRt+kt0kIr8sdYK2uWjt8lASjBlEFAsOg1nbeQqLzqCu7zNqyTvCk8Gj61CJyuIi8XUS+JyLfF5F3iMgRK5Q/TkT+XkT2iMiVIvJFEXmqt/xMNOe/FpBpwmj9V1I+V9+dnmWnP6t87llWrXMg0AVW7ZO5sWR9XzuOcmOp5DlbjzHP/FfSHq3n09bwSHTPF5E9jvNSSumU3AkiciCADwC4GsApmJDoiwB8UETumlK6fIfyJ0zLfwjAEwB8D8DtANzUUb9cxfV9lv6WW2nVemeOtYMEAoE28Ep0nvIaLe39mzgXeYju7pgQ007wtM4TARwD4A4ppXMBQEQ+A+DLAJ4E4OVWQRHZBeD/A/D+lNKj6KcPOu57fSWXvERv8HcOnvNK1KUtXJdXvU9uEHoJfRMHS6BftBh7HgLwqvpKxq8XLZ1RvHPGNo1RD9H9QkrpY43u9wgAZ8xIDgBSSueJyOmYbOxqEh2AkwAchwkhNoXXy7AWJdJdCxIOBAJt4R2HfXlTMmrJbBMJ0Jdmvx3uBOBzS74/G8DxO5S9//T/ASJyhoj8SES+KyKvEpEbN60lYV300C1tdC2eNex3gRLU9psW/X3dx3wX19wGG11LHIwl+9gBuATAQTuUPWz6/20AXg3gWQBOAPACAIcDeNSyQiJyKoBTAeBmN7vZ/Pu+1BC1yNkPa69X+2xj7dSBzUdfXs4l9/dIcV6vyxJ4VbGe8puCdQovmEmfb0kpPXd6/CER2QfAS0TkuJTS53WhlNJuALsB4LDDDkuzTuR9mS0JsS87Vq0tsLW6NMgxkEOXpFVro+sKLWxvtfX2hCdsCnZKAdZatbkXyyU3S9JjXDz9/171/XsAvATAPQDcgOgslBiwW0/yrQ3NNdcK79DAOmFMtievdOaVAkvgJb1tNSX0LdGdjYmdTuN4AOc4yuZwnacCq0p0jNbqziEkHQ+Jt1aXeu8T2B50NeGO2eOwdmFbooYsqc8mjsu+nVHeCeC+InLM7AsROQrAidPfcvhfmIQ5PEx9f/L0/5meCqxq6G5pcM0FgFpBorm/kmeoOScQWDd01fdrx2vpWLbqXfMM3uuuszNK30T3egDnAzhNRB4pIo8AcBqArwF43ewkETlSRK4RkZktDimliwH8CYAni8gfi8jPicizADwXwJsShSx40KV3VyAQCLRG68W25x6bQnS9qi5TSpeLyIMAvALAmwEIgPcDeFpK6Qd0qgDYBzck4hcAuAzAfwbwDADfAvBSAC/01qFGdZlDiV2vS+/M2vMsdBlXONZBEqjHUHFjLZ1Oaj0oS+zgXY6JbVJd9u51mVK6EMCjdzjnfEzITn+fMAkqzwWWu1DijNLaRhcIBAI51NrociS8TVin8IJqeFWVtSuaWq/NvmKExiTd6ett4qpy2zCEFNelC36JROdFV1KcR9rcBs3KVhGdF54Jt0S6y3kzxiQfCGweStSDterWvhbK64StJbq+4sa4fG0cXos61JyTQ0h3AWB4Ka7WBb/E3tbXGPWiVtrcxPG2tUTnhXfCrSW0knsGAoHxYQgns221vXkRRNcBam0BJefl7tkHcebUsoFACfpa5JWMvS6Dv4dc3I45RKAGQXQYPhVXDl5DcYkXaavzd7p/LULCHTf6WtSUqC5XRYmqb6hFXVdq1U0cY0F0O6Bkku3L3jUm6c57z03AmCaCTW7bIWLIcue1lOK69A6txZj6dytsFdGxWO4N5Pba3oaQ8PqQ7rxl+ppw+3KFXpfBPsQCYyhy7UOK6ysEqAReaXNVR7taB551wFYR3RCoDUPYBJVgX88TWMS6tPvYbHG18IQkbQO5jAlbS3TeSaAvHfcQYQhdBcZrjHmSZcRkMzz6IoMug78ZXUqbtenFas5ZN2wt0Q2BMQeZ9xVGMTYpYxMH9TJsWrsPscgaIvi7b4TX5QbD6xof0l19+QgsX08MHQjeAusS/N1S2uzKJr9uCKIbECVOHhYJdiXd5a7XUrorvUZgPdBi8hzCKWobpDiNdannKgiiWwLPBD6EdFd631XLDDUp9aFe28RBvCr6UmOOgQz6kOJa1LllPaOP3xBBdGuMMdldQjoLtMbYJuyx1acrbOJzBtE1gkfFkZNgWgaJ1saaDRGHp+87tvjFVdF6lT80at9H6wBpRkvprERqKnk2r6f3EH1grGOqBkF0O6ClY8gYYtW6MmCX1sdTpi9Hma7K93nfWjV1y/4xBkeIdVEDjkF9OrtOEF1gZXhDBVpKd1121C5Tog0RNM/YhAFeIkH0cf/WEt0QNlyvXbFWigu0RxDdCiiRmryqyyGwLmEItdJhi/O6Kl+CLuMXW7Z1CRmUoHVQdR9oQXpd9b2h56UuEETXI7xehmNd7a4Cz6q2lvRa1qv0vCFQG7PYpf2zxJtyaCk0B8/z5J6tpdd2Xxhz3y9FEF0hvLaNEtVlYBGt4/Vqz+vrWi3tNt54sCGk/G3BGGznHozxfYrI7QH8NoAHAjgGwGUAPg7gD1NKZ+1UPohuzdCHhNell+PYBpGnPpvgjNKXFGd935fduLVTR1/xfmPBiJ1RHooJyb0JwCcB3ALA7wE4Q0Tun1L6RK5wEF0jeCaSoVbLLR07SjJFXHfddUu/37Vrl/saHpTUrcvz+oDX7uvVHtRmAvGW8fYJq24tVX9DeHoObYdbQ/wdgNckahAR+QCA8wE8FcBv5AoH0a0x1sWb0YI12QH5Ca8GXZJZ7nlqYbWHl7S8drAS0rJ+67I9atGXd+g6YozPmVLas+S774nIlwDcZqfyQXQdoIRMWnq4ee/TpU2qdpIrWfGXwNMGY5iwvXXg9mnpEFQinZWAr6XfdVexZrXeoa09PYcmmqHv74WIHAzgzgD+eqdzg+gGxBDOKCUOCmOCnlQ9E3uL+6wLckQxQ207jbltvM416xIWMgTGPP4V/gcAAfDKnU4MousYQ6gKa9VZ3uu1dFzwQrehZ9LlCb+vSXoMXperto23jEZXttXcfVrGqo5BPTkmcmlQl0NE5Ez6vDultJtPEJGfA/Bex7U+nFI6SX8pIr8P4FcBPD6ldO5OFwmi6xEtXLlz16vBEMbxLp0AGK3JbUxu3n2Q4U516cMTOHfPLgPbLYxp7I4Qe1JKJ+xwzr8BOM5xrSv0FyLyZAB/DOA5KaW/8lQoiG5AlATKDiFFtSxfu3ouWSwMMSn1JdH1tXgaIlTAiyE0E4ycja7lffpAX+EFKaUrAHxh1XIi8usA/l8A/z2l9GJvuSC6NUNtBpVaVap1zxZS0zpMBDm0rn9LqaXkntb3fUlwOXusdW2v5Fn7DOtibyvBWMehiDwKE8eTv0wpPWOVskF0HaBW/dGXs4DlyNFCovRIcS083Fb1LOxygmrtXNTHZNqCGLqaGFv3Q0bJGPGipcPXEIQ6RqITkQcA+FsAZwF4o4jcl36+OqX0qVz5ILpGGEq9VQPLQ89bryG870q86oaYLNZ1xb+q80af2IT+WlN+XftUIzwIwI0A3BPA6eq3CwAclSscRFeIoRwSuho4Y3YZz6HEzmlhbNk2hq7PGLwRLdT215YqfI2uCKm1PdZ7nzEgpfQ8AM8rLR9EtwKGdurIITfwParDdUk42xpdTRZjW32vYwLxkom9VsWaC1/ZZ599dizvvc9O9131PkMvisaOILol6MuFvqtVdk4KHHMnLklZZf1WkmPQ4+wwdnhsUq2fp6VUXXIf7zm1BNRlOE9XEuaq1+rL67JvBNGhTQceYhB6sSo5tk5bVlKmxJnEQ265tFItk3EPhVWdc1qrrLvsO4ySRZ4HXYY35EizdpHnrdsma2dy2Dqia/USa72+9G8l9bImqRZ1s9CS3LxedSXkxmBy8xJdrSq3L5T0r1wZS/LLlR96gVDiyev1hrTs2LU5V0skuhY2uk0kMQ+2juhqUEsgLcls7Bja9sUTUY7AatWdQ7t/l6S8YujyVrq0XHt0pboc86Rcu/OGV6JrsRhcFWNu91IE0e2AluTmXUmXoKUUV+JRVmL7aoFViSo3sXskwtx5Q0C/K4+Un7PX8XmWC39fTlVeQi2pW5dkUrvzxtCkF0S3JSiRyDxlWuQVrC1TYgsoITf2UGup2mrtjOKR/Fq0Ry28fdJDTrk+UKLubNlfuW7XXnut+VsJAVjXKqmnF94Qnpaq0C5DJ9YVW0t0rdWQHm8373366mgljgNWGe163ZJAvHWzUGKj80h6ud9akF7tgov7niWp6f7p6a8tiK2EUJn4urRjWWhJehol9j8Pua2qnQmvyw1BjeouNyBrJ4haoiuRdEokOj5mcishutyArvXeKyHXErue5/sWKFlwWYRmnaN/q7U1e1X1HslTQ0t7M+SehzE2Sccj+XnHS620u4nYOqJbBq90VquGrFV9Ar4B2lJSA2yJKCcp8WdLjVlCdC0JUH/2SnRjttF5bHGtVZIlhGgRnUVgGhYxeBMneFGyMPZKVBas/uVNcl1ybcYmEuDWEl3JIB6K3Gqx6sRc4oKfI7oSAunKrtfa65LRpdON1fe8CxSvqrDEvlOiAbGkzVwbWnXL9cNVA+VbjD2PGrGF56rVbrVOK0F0G4BlL9G7KvZet9YLzIuWDhfWdXPlLalNf7YkvxKiK0GJWpfhlTy7lPS8UkJfUpwHXkK1JLKcdGf1Q6/UY9WttcNYie3MQgsbn6cOQXQbhBJ7W1+SWq1nYU66qnUSsa6Vs9F51Z1jUgnmUCv5tURJnxpDPT0xerl+aJGj1+boJbeWYT9ex5ASM4bVNi3s4JuArSO6WWfxdvQSSa3laq3EMaSE6HLhAB5C9d6zNdHVeru2RK1bd+56Q6Bl7k/vJF1CdF7JkSVEjxpUf/baD0ukuJbvOufYUhLMvgnYKqJj19mxkVuJGrLWMaRWpZgrU0uoFrzvIKfOWhW5ycI7qbV8tlzdauEhtxbvzbOwy/WpEqKzyufCFjxSk34HXXmueu16DF1mmbSn6xtEtwFYNjHkBscQeuzWqksPuZUQXe6eud9awiKXlraW3GTBqH1Or+2t9BoztHZwqIXXMcWz8arXg9oLr+THqLXFeepSep+ctLfqtdYJW0d0M5QMjtYGZI+Th5e0uPy+++7rKtPS7V+jxBZYglXdv1tMhCXeal1JdC01Ccs+L7tPyXsr6Td6Irb2g/MuECw1ptcz05MHVP9m1bNFGESto0t4XW4wltnoWq8CS2wOHjLK/dbSA7KFBFY7sVvneR2HrAkup2Yqgde2WoKutAleSdxCyVZHLepWez2LGHKkZZ2XU3dav/H33sD2Lr22twlbR3QztO4YLR02chlHPORWch/vZOUdkC0nae/q17KbeImuqwVOKfqom1cd5rV5tlRTe6VNrxRmSXReu6BFVLmcnNZ9cqnXulQpeiTzTSTN3olORA4H8AoADwEgAN4H4GkppQsdZY8A8EIADwRwKwBfA/D/A/iTlNLlnvsve4kt1FctY828ROeR1Lx1q/Wk80paVnnvb16i8pDeTr9Z9yzxRuzLGcVTN2/cmdWeXieR2vYoSVyQu4/Vhrn2sOyCtQ4sJR6cJdqIErteEF0lRORAAB8AcDWAUwAkAC8C8EERuWuOrETkJpiQ4n4A/hDAhQB+CsDzAdwOwC/vdH/2KCrRY5eoIUtIy5tD0mOv8z5DrRSXG4QeMtG/dUV0Jd62O9VhGfqyP2p44tNyKrRVYy6XXcNTxjovRzqe65XYAnNSjofEvW3oJUc+z6vuLJH+LTNOEF09ngjgGAB3SCmdCwAi8hkAXwbwJAAvz5Q9ERNCe1hK6T3T7z4oIgcDeIaIHJhSusJbEa/6xvo+NyC83oy1oQLeuq2qUvN6oZY4AdSSnrduXomudrKwMBTRWXXwSBkl1wLqPRMtMtGLvJL305UNNaei9fQpryRd4sBS6/UZRFePRwA4Y0ZyAJBSOk9ETgfwSOSJbv/p/++r7y8FsAsTNWgRSiSgFtKZ57yc12VLw703O4QVdKvLW15tXtKxvveqSL3kWkJ0JR6IXTnneMnEK8l7JtYSZ5TcGMst7Dx1y40XT1vVqkFzNjpvu3tshvo+lrS3iURVi76J7k4ATlvy/dkAHrND2fdhIvn9qYj8Fiaqy3sDeCqA13psdCKytCOXDEIvabGrfwt7W63qcdVjwKcGzA3CWjVi7vtV1aq58kNkU+kLXjWkhxhyKkXP97o+1rF+VzlCs8qsSjolji26jGXXK7GjezU1NYvJXB/YFPRNdAcD2Lvk+0sAHJQrmFK6SkTuD+AfMCHGGf4SwFO8FZi91FyH9hBQLlbNQ3q6DiX2NgveydwyiHs9wloTXUs1S629bRMGuzUR5shkVQmopC66Dl6i48+1qn4Paep6W8decs3VyyrjXfRaY9G7SztjE/q+xtqEF4jIAQDeBuBQAL+O6yW65wK4BsBvGeVOBXAqANzsZjebv2ivO76XtGodS1pKaiVElSMgz3nXXHONWTeL9LzSWS1qs4qs08C3JjKvXa4rFZhX6uExklsIlRCdRwpsIWlZ7e51tPE+m6XWzJk0PAvddervXvRNdHuxXHKzJD3G4wGcBODYlNJXpt/9q4h8D8BuEXltSuksXSiltBvAbgA47LDD0oygStSQJURXoobMweNwoTuwh3Ss73PnlZBrzq3aGmAtHDs88JCE/m1oeG1KXtSqb0veIY+XHOlY76Rk0eol15a5YXN9qpboPJ6e+rNHu7Up6JvozsbETqdxPIBzdih7FwB7ieRm+Nj0/3EAbkB0GrNO4PWG9DqJlHRUq3yJl2GJjc2jXtTn1XpDWufkzsuhq0HqVU0NQXqtyW1VeB11cvD0d2//8MakeerSwoW/xDnHQi73pyclmjfOksuOaSHXCn0T3TsBvExEjkkpfRUAROQoTEIHnrVD2W8DOEhEjmWvTQD3mf7/xk43F5G5JJZTKVoOJCUZS0omHq+k5JHOvOflJDqPY4fXLlgykTFK2tPrkWpJM7XvsBQtbWS1bZCDlwBWvWduweZVCVqLVo8Hpz7PIx3qz16tjUdy1OO6JDxp2Zgfs8aiFfomutdj4jhymog8B0DCJNPJ1wC8bnaSiBwJ4CsAXpBSesH06zcCeDqA/ykiL8bERncCJsHjnwBw+k43Z6JjMsupIS3SK/GG9JKBl7Q8x/p6bEuzVJq5CcYaBGPzWCxZbFgToYY3+0ctPNf2agm8yJFGHyhJL5aTeqy+y/2dx3huYWn1D0uy0vCqOz02Pn1ffh4e47oM/2apLoPoKpFSulxEHoRJCrA3AxAA78ckBdgP6FQBsA8m8XGzsueLyH0BPA+TbCqHYEKQuwG8OKW046gUEey3334A8m7/s3P0byXekF6VotVRc2pEL9Gt6gGZq2eX8EgwXvdvbxnPBJNzc7cWK14nIi8ZeYnb04beicyTFaRPePqhV81sOd14g9S941ovopfdX9fTGzrhkeJ+9KMfmfcJousQaZLT8tE7nHM+cMMA8JTSOQD+Y+m9meiYzLzb2ng9rbzekBa55UjLc57XxlYb08bwTn7eNvRee1WvuFrSy9XTKw3VklvJ81jf15LeTnVY5Zyd4KlrSYo37yLPuzD0SngMz/vN2Q+9qkuPRLeJWJvwghawiC5nbyvxhrQGkZe0+Hvttu+R4rwekFb9vWitGqsdcC3JgLFOK1zPc5dIlCVOJjmUqEg9BO1VOXvjCj2puVo7VXn7sdXWq2qeQqLbQMxIrNbt3+t55o1PqyWtEtVjXx3aUsXk7Es5uwtjVXIskTxzGNrr0jv5WeeUEFiL98boyy5oEXxOLWv9lnsHnnReuXHtbbcSG+5OtsXwutwAeL0uPfDaznLSmfWbJenpzy1Ta+XQ0inCa7sqeSeBfuD1HvSWycXOMTx2Qm//tkgv59lolW8h0TFKnKes8p5xFBLdhkFElsbReeGNNfMQGLBoKObjnLqzq8wGtWpIr3pwKFf9lvcceiLw2g+HsL20DInQ8Eh+XkL1fK/v4/HA1J89Hssa3jb0StmMZarYsNFtIGYv2vtyPZ6RwCJR8W/W9/o3b5aSWnKrJTTvtUrIrcQD0kJLMipRU5fUJ6c2sxYVXi9Db71aLpJq1cS5Z/OEGpTex/otJ/l5oD0gPXUoeR9e7YrXu3MTsHVEN3uJuQHhyeeYU0N6SE9f2xv8baF2VVZiB+tyJTg0ueWuW0IgJfepVc+V3LMvtVtX1/KqPmvvw2hNel5Y6lNvNqaW5Dp2bC3R5b7zOImU2M5yk0pO598VaiW1dURrMqnNJFKiWcjZXTzX3oSJrCQUZAinl5LYO2v+yakua7UmjE3oHxpbRXQppaWrr9bxbVaZnKdVy85V66XYmsy6dGevbbcSFWkLaW2n+5dcq0WZ2pjJvp7HWx8LJVlwSuBVbVsZi7zkaGV6ye2YYt0jiG7NkVKadyivw0cJ0XnDCzwdqqWr8SrndYUuJ2lrIilxPGo9+XklEA9KXNG7bMNae2qXqO3vXkL0wEt6jNz8YXmOexNBz34bw3vqGltHdD/84Q8BlKkhLTLT55WoJGs9E4dGi8HS1XN78ydaz9BazdXyOXNJjPt6Hg/G3D+GWCR654/cwoHnI5bi+HvPFmSt7bRjxFYR3XXXXYerrroKQN6D0iORte4MLQfb0Kot73ljiN/patIfarEydHLtIfpHl0TVcrx51eHWebndwvm3XNjPMrVmiaZp3bBVRJdSwtVXXw3An+2/ZOLw5sf0Zlpoia5IsEsvQy+Gzry/aShR+Vroywu1NZnV2rS9Dh8W0dXOP7ltekKi21CwjS5HdJ7VVW3qJf25xIW/L3f6vu5fMmGNzSa0qRhDoH+JR+kQ6s6uHHK8hJQjymVzWxDdhoFtdLmtcDyBuroDW8bgnCt4LdExvC7Wnu+99xkC3kwgfUl3LaUeL/p6Hm+fHKJP1JJJyXhrqXXxOjvxed4Fec4ZhT/P1JhDj+k+sLVEl3P7Z1idRhOYl9w8qFW/tCawvgZC7eRj1bM16ZXEH1rnlUjPLZ+nhVdvy4w2JfAu8rqy0eXK1D53zi/AmsNy5Lis70R4wYYhpTTPSOB1LCmZCLyeltuOEtvIGFRoDG89PeflVEh9kUaJ3dhTzzFIgWOGx1Pb60uwatanUF1uGFK6PmC85R5Y+jevGqFlh6q9ljd/4hBo6VWn30FLNaCXGEpi6lr2Fe9+ZV1KqH31qZYqxlpJPLcA9oQnea/t7dNho9tgLOsEuYm91mtqp3vP0DowfNXyXU64XU5qQzuj1HrOjk1CtdAXmXXZv7rybPaSUclWWiWkl/uen3ubPJO3juiWoXWqIu9KyzN59KXy8U5kQ6/EawmkpP1yG/N60TIey3K6AVafvFp4D7dUsbbuX30RmnWedy4oMXd4zluXhVTXCKJbAo/nWU71YEE7qVhZz0ti73Log5y6VHe2VKHpa7V07Biz6tK7dyCjtt1beya2XGC0rJtXu8Pj3etY4vUlaOlxHES3AWi1c3Wuc3s7l6Wu9G6zwWg5sXYpRQ5h7/NKGSWhArWkUSJJt3yGMUhQ3kVJibTpQa0dPrfo9TqMeMitROJfdb4Lr8sNgIiY2QAseM/jzmkNNt25OR2PFx7X9hJvxk1DrZSRKzN0G5Z4QHrLD70QqT2vtQ27RA3pITeddtAqX0Juq7bnOnrEishjAfwtgG+klG670/lbRXTAzttULPu87PsWBJhLvmrB6tBdTsxjsiW29C7NXbtWRdyX6nJoR4wurrfqdbu0w3nmghzRWRmYvOrOWqIb2yKtBUTkFgBeCeDb3jJbRXSWRJdzErHUkN6Omut0TLrWBou6jOe8nG1jXTp6S0L1to0XLe2ktU5A3udpaavqi4D6uk+JZ2OOtCyiq3VG8drrc+rrZbbaNZTo/gzAWQC+BeDnPAW2iuiA5YMi57lW60WWW7F7vNW8oQ8tnQByqJ1wS1aoXdoI+2q3rtAlUZVgzMRZAg8BlUiBq/zmgdU2OY/h2jE6BETkRAD/CcBdATzHW26riE5Elm42mNvTqyvSy/1mSW25OnivPbSayYuWLus59DV5jtkFvxa10uLQz1PrgKLnD0uKK1FJ5uCR6HJapHWT6ERkPwC7Abw0pXTuKv1mq4gOWD6oWts5ald0uYnQGmC16rna5MQlmWY0+iI3CyX3LJmka8ltaKnYe0/P96ugjz7qHW8li8yWtrdln1cts6z8yL0unwngRgD+ZNWCW0d0y5AbhJZ01zqVlCf2BvDZ5byB6S0nJe1MYw2WXDtt0sTcIlTAU76kzNBtq1GbJDuHVR3LvLYzyw6XK8MoWVzn1JDeXVZ6wiEiciZ93p1S2s0niMjPAXiv41ofTimdJCLHAng2gEellK5atUJbR3SrvnirM+XUna1VFAyPRJdTRdSqMWsHztCbo+bq73222nfapc3Rc88hSE+jZDsgC7VmhNx4LQnetlCrOfKoIXcq05ONbk9K6YQdzvk3AMc5rnXF9P+rAHwAwBky8boEgP0ByPTz1SmlK62LbB3RrYqWdoYSrzrvNXKD1UO8uW2GxixBjGHS9mBV79Cxoa+6te5rHq9p79jxfL9K3azzvBJZV5JbH+MopXQFgC+sUOR4AEcC2Lvkt70A/hzA06zCQXRor6rz2s5KOpRHisupLj1qzFqX9dy1GV16u3aJEocgD0omwrGhpfdfl3Zwa7zmVJfepMwWSuxtJTtN1PaVsY23KR4L4AD13bMA3AvAYwB8PVd4a4muRFLzrvxqbWdeeMjVe14uYL1E3TlWoipx1GlxnxKMuW5docuJuSTNVkt1pdcXoKUacsyLolWQUjpDfycij8NEZfmhncpvLdHl4Ok03pVjCQGVwEuoluRXMgi96NAWsDZYd9Vla7SUUHN930NuXonOi1rpWO4m3wAAFN9JREFUzDvetnEhVIoguiXwuDF7vQdzKsUSKdCqg+UMk7uGl/Q8hFiiZhobvJKrx77TQk3tQbR7XuqyfrN2C/CGF+RQIp2VLCZr56llGHl4wQJSSo/znrt1ROfpHJ6O5lX15Zw8PCvMkjAGrzMK59rkBNO5QcgTRMlOEOsyiLpErefrJqDEC9QiKq9Ex/3du0DxoMTtf9U0XaUoiT3cxDG6dUS3DK1dzr0hCSXSkcc930t0Le2PGusygXfldt/Ce6/Wk7fkvDGhpA1ry3jbqdZJZJXrecozcgspz/OtY1/ZCVtFdCKytLO0IDqPN6O+Vi7V1ww5Cc67WmtpC6ytSy1yA7cvB5hal/OS+3j6V2v05eDgbU/P9jc6eLt2Q1OGNyjbs5+kN4SnjwXjuixKa7BVRAcsf6mtDb7eCZc7e4nnGMOr8vFcqzWswV6SNqwkfolR0k459LX67YrQW0xynmvkEizkzrN+82yFkytTIsVZ5KafxSK6Wg/KvhAS3QZg1olKJLWa+wF+Q7dXOvKSnufarUkv54TjQcnkU+LZOGZy89y/pRdsizLeuq0ayK3LWOSWs9F5nFm88X45qc1jl8sRXZfwjJeh+3cXCKJznt81vDk1GSUd1ZpgalWP+p4lufc8be1VXebqZl2vpD3HJtHVqsBK3oEH3jb0OpZ4VJI7/eZBCdF5yG0ootsJ6+R1uQqC6ArKroKS1bc3H6RnOyEvWsbx6c9dJpz1kFtrF/yWsVUl6Euiaxnf5nWQKiEtbxkPvPFtuTIl4QW563nQ0m4bRLfmsJxRlp3X8p4zeFffuU5rufrnsGoiZW8Gl9zq2yMplUhnXpRId9uILiWJEiku16dKynhQos0oCepuTW6e8h5txlikyS6xVUTHaL2qbWkfKhk4uYBxD3Ik49lI0qtizZ3nmRT6Wu2OwcmjFkNIE9b3JR6UJRJdrg6ePtVCOutqMd3aTho2ug1GV6uYWgnCW94iN2+QeUkcnnWtXDB8CSyy7tIhqKUHo3XdUvRRn9bjoMTeVhIqYN3TC28gt9eZxCrTl4mkZfkgug3AGFbaq6Bk9V3ruuxVOfG1cmrUEtVlrW2T4c032mXfKGn3EsmgViqu9cr1BmJ7HEtKNjTVWFWKKyGwdZtTZtgm9f7WEd0mwSPd6fO8ZRglE4xFfN6JpCSI13OONztNa4eTIVbZnsnYa5+q7Sve+DZLTZ67tleN2FcOysD4EETXAbp0U/esUEuM67WxbrnwAisDjH4Wz350tWTi3QOvtfqoK9tiiURXMnnXhgd4Jbpax5IcrJAAb6hASf8oCQvpS8uwDBFeEBg1vNKQZ/JrsVWJZ4DnJr8WNr8adGknaTmptSDUruDdFcDzvRc5iSx3nvd624AguoCJoT2YvOqXkiwp1mAvsdvk8nt64gJbT1AtYxFL7t+X5FgSy+jt0x4PSP29x5M3N3a8UqknHdcYMpZ4FoZD1GVTEERXiBadoeXqNefCv2p6MS1NeextJfaU3IRpeWC2lgI9NsshPOda3KelfckrnXlVl57yGh41ou4DHhWlV+3f0s5aq8puUR/vfTYBW0d0s5dYspJd9R6tr5uDd8LzkF6unp4dF/Q1vA4KfD3vprLWPa3r5rwuGSUOLK1Rcp/a/RY9pJMjOm+ogOc+Xvsjk1tJguUhHE5a2OhKYndXLb8p2Dqim6ErMmuN1is/D+mVkEmJGjNHOtaEWaKq88IrIVvocjeI2udheN9vbSB360VJiTNJLbkNbaMrsbd7y7csM3b0TnQiclsAzwRwAoC7AbgxgKNTSuc7yu6aln0SgFsD+CKAF6SU/qGmTl2+2CE6TY5ArIkkR0Aem06LlaNHBaYnsq7at+R5cmTkIcGWZNYaXvVxSSgKw2s/LPGGrO3HXaJErdk6BGZ2zSC6NjgWwH8E8AkA/xvAQ1co+0IAzwDw7Gn5xwL4exF5eErpf3ousOpLHPNLrx2UJROr1R7evJs5eCRMb4yghZz7uIUc8XvVey1JrHZR4XXh90p0ngVKDl7SslSUXomuq128W6O15+mq1x7znFeKIYjuX1NKPw4AIvIEOIlORA7FhOReklJ62fTrD4rIsQBeAmBHoluX1UrrVWUJAVjwtp9lyyuR7kpUl4ycva2r7YRyqI0LzLWnJ5TDS24lW+F436+nrUscS7yqy768aoeYb2oX8+swR66K3okupVRqyHgYgP0BvEV9/xYAfyUiR6eUzluhHoXVWI6+PPP6mGS9qksvOXodWKwyubqtipxHquW12YL0atvde/9Vbag5t38P6eXOy8FDbiWk1Vrt7kGJHW0MtrOh7Y99Yp2cUe4E4GoA56rvz57+Px7AjkTX1WqlZCKz0GUH9EhXOScABk9qJRJhbRxeC5Qkku7q/ZQQXa3jUC41V0kOSi88TiI5NbNFiEOMsRJVoXeBMoSKNCS6YXEwgEvTDd/CJfS7G13qwVter/V9amN5PKSXOy/X7h6pQ183lxtxGUomQi/p1bZtyX00VvWazEl0FtHlnFEstPCG9JDbEBKdhtUPWhNIV3NYEN0aQkROBXDq9OPVz3/+8z83ZH1GgEMA7Bm6EgMj2iDaYIZoB+AOdPxuTNqkBqNrz3Uiur0AbiEioqS6mSR3yZIySCntBrAbAETkzJTSCd1Wc9yINog2AKINZoh2mLTB7DildPKQdekKq/s8D4ezAdwIwE+q74+f/j+n3+oEAoFAYB2wTkT3LwB+BODX1Pf/CcDnVvG4DAQCgcD2YBDVpYj80vTwXtP/Py8iFwG4KKX04ek51wB4U0rp8QCQUvquiLwcwO+LyGUAPgnglwE8CMAjnLfe3eoZ1hjRBtEGQLTBDNEOW9AGMlCKKuumH04pnUTnvCml9Dgqtw+A3wfwRCymAHt7pxUOBAKBwNpiEKILBAKBQKAvrJONbilE5HARebuIfE9Evi8i7xCRI5xlDxCRl4rIt0TkShH5PyLygK7r3BqlbSAiJ4jIbhH5gohcISIXishbReToPurdEjX9QF3nWSKSROQjXdSza9S2g4gcJyJ/LyJ7pmPiiyLy1C7r3BqVc8IRIvKm6Vi4UkS+JCIvEpGbdF3vlhCR24rI/5jOaVdM+/RRzrK7ROT3ReR8EblKRM4SkUd3W+NusdZEJyIHAvgAgDsCOAXArwO4HSY5MD0d8w2YqEGfC+DhAL4F4N0icvduatwelW3wWEwyzrwKwM8DeBaAewI4U0QO76zSjdGgH8yucwyA5wD4bhf17Bq17SAiJwD4KCbezU8A8O8A/HcAq+9oOxBq2mD6+/sAPADAH2Ly/H8J4L8C+KsOq90FZsnz92KSPH8VvBDA8wC8GpN54QxMkuf/u5YV7BWzRMfr+AfgqQCuBXAsfXc0gGsAPH2HsncDkAD8Jn23LyZ2v3cO/Ww9tcGtlnx3JIDrMLF9Dv58XbeBus67AbwOwIcAfGTo5+q5L+zCJETnH4d+jgHb4KHTOeGh6vuXTMsfOPTzrdAOu+j4CdPnOspR7lBMUi0+X33/fgCfGfq5Sv/WWqLDxNvyjJTSPP9lmoQZnA7gkY6yPwLwNip7DYC/A/AwEblR++p2guI2SCldtOS7CwBcBOA2jevZJWr6AQBARH4VE2n29zupYT+oaYeTABwH4OWd1a4f1LTB/tP/31ffX4rJQmBtsiCnbpLn32UdzRrAmqsuMVG7LUvpdTauDyTPlT0vpXTFkrL7YyL6rwNq2uAGEJHjMFnVfb6yXn2iqg1E5CAArwDweymlpRl21gQ17XD/6f8DROQMEfmRiHxXRF4lIjduWstuUdMG7wPwZQB/KiLHi8hNReRBmEiJr00pXd62qqOEJ3n+2mHdie5gTHTQGpcAOKii7Oz3dUBNGyxARPYF8FpMJLo31FetN9S2wUsBfAnAGxvWaQjUtMNh0/9vA/AeAA8B8GeYqL3+plUFe0BxG6SUrsKE8HdhMrFfhonK7p8BPKVtNUeLpsnzx4J1ynUZ6B6vBnA/AP8+pbRsstg4iMjPAPgNAPdcMri3CbNF71tSSs+dHn9oGrv6EhE5LqW0TlL+yhCRAzAh+kMxcWK5EMC9MXFWuwbAbw1Xu0AN1p3o9mL5Ks1a1emyRxplASNJ9AhR0wZziMhLMNnl4ZSU0nsa1a0v1LTB6zCRXr8uIreYfrcvgH2mn69MKV3drKbdoqYdLp7+f6/6/j2YOGPcA+uhzq5pg8djYqs8NqX0lel3/yoi3wOwW0Rem1I6q1lNx4mi5Pljx7qrLs/GRKescTx2TvJ8NoCjp+7IuuwPcUMd9VhR0wYAABF5NoBnAvjdlNKbG9atL9S0wXEAnozJAJ/9nQjgvtPjdVrF146HHEqdG/pGTRvcBcBeIrkZPjb9f1xl3dYBG5k8f92J7p0A7juNfwIATIMiT5z+lsO7AOwH4DFUdl9M8me+Z41W8TVtABH5XQAvAvDslNKrO6pj16hpgwcu+TsLE4eGBwJYp/RyNe3wvzBxQniY+n62bcuZWA/UtMG3ARwkItoR7T7T/99oVMcxYzOT5w8d31DzB+AmmEhen8XEdfgRmExSXwVwUzrvSEx07M9V5f8Ok1X7EwA8GJNJ7SpM7DWDP1/XbYBJwPh1mExy91V/xw/9bH31gyXX+xDWM46udjz80fT7Pwbwc5gkELgSwBuHfrY+2gDAUZiEFnwJk2DzBwL4b9PvzgTFpq3DH4Bfmv79BSZxdL81/fyzdM41AN6gyr1kOg8+HRNV7l9M54mHD/1MxW0xdAUavMwjAPzDtDNeBuCfoAIjpx04AXie+v7GmMQNfXv6Yj8K4KT/297dx9hVlHEc//6sNYQXFdSyQNUqUkytxsSqAZLSiEQi/mMiGGFrNqYG8SWiFqO00pc0VoLGYJCINKG0TcAmoPi2KSG6xZeIEENaBWqJSgwvtqWYFAItgcc/nrn0cvbcfenu3nN7/X2Sm5uZO+fMnNPtPjtz58w0fU3dugfkLMPo8Bpp+rq69XNQc66jMtBN9T6Qz4l9tQSKQ8CjwFpgdtPX1cV7sADYCvybDPJ/B74LnNj0dR3BfRj3/3ZJb6wcN4tcIehRspe/A/hE09czlZcXdTYzs752tH9HZ2ZmNiYHOjMz62sOdGZm1tcc6MzMrK850JmZWV9zoDMzs77mQGc9QdJPJO2XNFDJnyXpPkm7e2m7GEnzJIWkoba8IUmfqSk7VMrO62ITW3W/StIDkpa35a0u7ZmxtW4lXSFppyT/jrHG+YfQesWXyIdXb6jkLwfeByyLiOe63qrOngDOAn7VljcEjAp0pcxZ5ZhuGwROYfR9nWk3Am8iVxgxa5QDnfWEiNgDfAX4uKSLACTNB1YDN0bE9gabN0pEHIyIP0XNLu01ZfeWsk2sn7oc2BSjNxieUeWPkk2lfrNGOdBZz4iITeSistdLeiO5fc5e4OvjHds2PLhY0s8kPSPpKUk/rA55SjpF0iZJ+yQdlLRD0mClzICkWyQ9Xso8IemXkuaUz18xdClpBDgXOKfkR8mrHbqUNFvSOkn/knSovK+TNLutTKuOyyStLW34r6RfSJo7gXvyQXJF/nE3TpV0Qbln15fhzlbdn5O0XtKTkg5I2iLpWEnvkLStHPOIpLqe223AAklnj1e/2Uw62vejs/5zGblVyL3A28lNYA9M4vgt5FqFN3B408zjyGFFJB0HbCf3LLuKXNNwENgs6diI+HE5z2Zy4d8rS5mTyYW/q9s6tXy+1D2rXAPkWoud3AJcTC6g/Htyw9sV5ZovqZT9JvBHclh0DvC9UteSMc4PufPAAXJR444kfRrYAKyNiHUlr73uEXIIcgG56/hL5P50N5HrQF4O3Czp/oho3+7ngVL/BaX9Zs1oerFNv/yqvoD15Pd1t0/imKFyzI8q+SuAF4H5Jf3FUm5JpdzdwB5gVkk/Q+7P16m+eeU8Q215I9QsBt3WtnklvZD6BYVXlvz3VOoYqZRbXvJPHeeeDAN/qMlfXY5/NdlbfoH8DrTu+n5Tyb+j5A+25Z1IroK/qqau35HbXjX+c+XX/+/LQ5fWUyS9FlhK/jJ9v6QTJnmKrZX0beQQ/QdKejHwWESMVMptISdPtDaYvA+4UtKXJb1bbV2cabC4rc5qGyCHQNv9upLeWd7fMk49p5JDv518H1hDrky/oUOZ4Ur64fK+rZUREU+TfyS8ueb4vaUdZo1xoLNecy3ZQ7iQHKZbP8nj/9MhfVp5P4n62Y9Ptn0OuQHvz8kezw7gMUlXT9N0+VYd1XZU29Cyv5JuTWo5Zpx6jmkrW+dT5Aazd49R5ulK+tAY+XXteY7cDsusMQ501jMkLQE+C6yMiGFy5/PLJzmZ4eQO6dbu0PuBAUYbaPuciNgTEV+IiNOAd5J7963h8PdvU9EKXNV2DFQ+n6qnyD8aOjmP7BUOSzp+muqsOgnYN0PnNpsQBzrrCWVm5E3kkOF1JfsacmLKBkmvmeCpLq6kW7uo31vS24G5ks6plLuEHH57sHrCiNgVEVeRvZiFY9R9kIn1Xu5pa1u7S8v7yATOMREPk5NbOvkbOaHlDGYu2L0N2DUD5zWbMAc66xVryVmOyyLiJYCIeAFYBpxJTiqZiI9KulbS+ZJWAKvI58h2l883AruBOyQtK9PqNwPnA9+KiBclva6sxnJF+fw8ST8ge0d3jVH3g8BCSZ+UtEjSmXWFIuKvwK3AakmrSluvJieJ3BoRO+uOOwL3AKdLekOnAhHxEBnsTge2HcF3oh1Jej0wn8OB3awRfrzAGidpEfmw+Lerv+Qj4s+SrgO+IWlrvHL6ep1B4GvklPdDZC/x5YeWI+JZSeeS0+S/A5xA9jiWRkRrMsjzwF/IYdS3kj3CXcClEXHnGHVfQwblDcDxZO9xSYeyQ8A/yEcGVgKPl+PXjHN9k3EneS0fIx9nqBURu8o9+S1wl6SPTFP9F5L/Bj+dpvOZHRFFRNNtMJuy8uD2zcAZEfFIw83pGZI2AnMj4sMN1D0M7IuIpd2u26yde3Rm/W0N8JCkRRFxf7cqlfRe4EPAu7pVp1kn/o7OrI9FxD/JYdI5Xa56gHyY3r1ra5yHLs3MrK+5R2dmZn3Ngc7MzPqaA52ZmfU1BzozM+trDnRmZtbXHOjMzKyv/Q/c5aSi2RxsvwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 576x432 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#NBVAL_IGNORE_OUTPUT\n",
    "from examples.seismic import plot_shotrecord\n",
    "\n",
    "plot_shotrecord(rec.data, model, t0, tn)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "assert np.isclose(np.linalg.norm(rec.data), 450, rtol=10)"
   ]
  }
 ],
 "metadata": {
  "hide_input": false,
  "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.7.4"
  },
  "latex_envs": {
   "LaTeX_envs_menu_present": true,
   "autoclose": false,
   "autocomplete": true,
   "bibliofile": "biblio.bib",
   "cite_by": "apalike",
   "current_citInitial": 1,
   "eqLabelWithNumbers": true,
   "eqNumInitial": 1,
   "hotkeys": {
    "equation": "Ctrl-E",
    "itemize": "Ctrl-I"
   },
   "labels_anchors": false,
   "latex_user_defs": false,
   "report_style_numbering": false,
   "user_envs_cfg": false
  },
  "widgets": {
   "state": {},
   "version": "1.1.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}