diff --git a/notebooks/case-studies/DCIP/DC_over_TKC.ipynb b/notebooks/case-studies/DCIP/DC_over_TKC.ipynb
new file mode 100644
index 0000000..37086b0
--- /dev/null
+++ b/notebooks/case-studies/DCIP/DC_over_TKC.ipynb
@@ -0,0 +1,872 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# DC resistivity over TKC\n",
+    "\n",
+    "\n",
+    "## Purpose\n",
+    "\n",
+    "Understand basic setup and physics of a direct current (DC) resistivity\n",
+    "survey within the context of a kimberlite exploration. Run DC forward\n",
+    "modelling and inversion using SimPEG-Static package.\n",
+    "\n",
+    "<img src=\"../../../images/tkc/dc.png\" alt=\"Scipy_2016_PF_Thumbnail\" width=\"400px\">\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Set-up\n",
+    "\n",
+    "The physical behavior of DC resistivity survey is governed by the steady-state\n",
+    "Maxwell's equations:\n",
+    "\n",
+    "$$\n",
+    "\\vec{j} = \\sigma \\vec{e}\n",
+    "$$\n",
+    "\n",
+    "$$\n",
+    "\\vec{e} = -\\nabla \\phi\n",
+    "$$\n",
+    "\n",
+    "$$\n",
+    "\\nabla \\cdot \\vec{j} = -\\nabla \\cdot \\vec{j}_s = I_0 (\\delta(\\vec{r}-\\vec{r}_+)-\\delta(\\vec{r}-\\vec{r}_-))\n",
+    "$$\n",
+    "\n",
+    "$$\n",
+    "\\vec{j} \\cdot \\hat{n} \\ \\Big|_{\\partial \\Omega} = 0\n",
+    "$$\n",
+    "\n",
+    "where:\n",
+    "- $\\vec{j}$: Current density (A/m $^2$)\n",
+    "\n",
+    "- $\\vec{e}$: Electric field (V/m)\n",
+    "\n",
+    "- $I_0$: Current (A)\n",
+    "\n",
+    "- $\\delta$: Volumetric delta function ($m^{-3}$)\n",
+    "\n",
+    "Consider a simple gradient array having a pair of A (+) and B (+) current\n",
+    "electrodes (Tx) with multiple M (+) and N (-) potential electrodes (Rx). Using\n",
+    "giant battery (?), we setup a significant potential difference allowing\n",
+    "electrical currents to flow between the A to B electrodes. If the earth\n",
+    "includes conductors or resistors, these will distorts current flow, and\n",
+    "measured potential differences on the surface electrodes (MN) will be\n",
+    "reflective of those distortions. Typically kimberlitic pipes (including those\n",
+    "containing diamonds!) will be more conductive than the background rock\n",
+    "(granitic), hence, the measured potential difference will be low. That is,\n",
+    "contrasts in electrical conductivity between different rocks induce anomalous\n",
+    "voltages. From the observed voltages, we want to estimate conductivity\n",
+    "distribution of the earth. We use a geophysical inversion technique to do this\n",
+    "procedure.\n",
+    "\n",
+    "We work through each step of geophysical inversion using [SimPEG-Static\n",
+    "package](http://docs.simpeg.xyz/content/dc/index.html) under SimPEG's\n",
+    "framework having two main items: a) Forward simulation and b) Inversion.\n",
+    "\n",
+    "<img src=\"../../../images/SimPEGFramework.png\" alt=\"SimPEGframework\" width=\"400px\">"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Forward simulation\n",
+    "\n",
+    "\n",
+    "A forward simulation of a DC experiment requires [Survey](\n",
+    "http://docs.simpeg.xyz/content/api_core/api_ForwardProblem.html?#survey)\n",
+    "and [Problem](http://docs.simpeg.xyz/content/api_core/api_ForwardProblem.html\n",
+    "?#problem-class) classes. We need to the pass current and potential\n",
+    "electrode locations to a DC survey class. The physical behavior of DC fields\n",
+    "and fluxes are governed by the static Maxwell's equations. To numerically\n",
+    "work with these equations, we use the DC problem class, which handles this by solving a\n",
+    "corresponding partial different equation in a discrete space. For this, the\n",
+    "earth needs to be discretized to solve corresponding partial\n",
+    "differential equation. The Problem class computes fields in full discretized\n",
+    "domain, and the Survey class evaluates data at potential electrodes using the\n",
+    "fields. The Survey and Problem classes need to share information hence, we\n",
+    "pair them.\n",
+    "\n",
+    "\n",
+    "### Mesh\n",
+    "\n",
+    "We use a 3D tensor mesh to discretize the earth having 25x25x25 m core cell\n",
+    "size. Smaller vertical size of the cell (dz) is used close to the topographic\n",
+    "surface (12.5 m), and padding cells are used to satisfies the natural boundary\n",
+    "condition imposed."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 20,
+   "metadata": {
+    "collapsed": true
+   },
+   "outputs": [],
+   "source": [
+    "import numpy as np\n",
+    "import matplotlib.pyplot as plt\n",
+    "\n",
+    "from SimPEG import (\n",
+    "    Mesh, Maps, Utils, DataMisfit, Regularization, Optimization,\n",
+    "    InvProblem, Inversion, Directives\n",
+    ")\n",
+    "from SimPEG.EM.Static import DC\n",
+    "\n",
+    "try:\n",
+    "    from pymatsolver import Pardiso as Solver\n",
+    "except:\n",
+    "    from SimPEG import SolverLU as Solver\n",
+    "\n",
+    "%matplotlib inline"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaIAAAEWCAYAAAAkUJMMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHW5JREFUeJzt3X+wJXV55/H3RwgzDhod0BB+WDDoRlwxmazDoAnZZaEc\n0GxQjCQQowMxha4xxjVEnfxY/EFMTMU1G1MroM46mSRAZCVijM6yKyyVH7IMyggTJCAOAmsMMCDI\nKII8+8fpC8db557T9849p8+deb+qTnHvt59+vn26Ljx099PdqSokSerKk7reAEnS3s1CJEnqlIVI\nktQpC5EkqVMWIklSpyxEkqROWYikjiU5M8nf9v3+rSRHdrlN0iRZiKQJSHJckr9P8s0kO5P8XZJj\nBsVW1VOq6rZFnv9NSbYmeTjJxwYsX5HkvyW5p9nGqxdzfmmYfbveAGlPl+QHgb8G/iPwl8B+wE8B\nD09wM/4fcB5wEvDkAcsvpPffg+cBO4HVk9s07e08IpLG70cAquqiqvpeVX27qv5nVX1pUHCSSvKc\n5ucnJ3l/ktubI5W/TfLkZtmLmqOs+5NsS3L8XBtQVZ+oqr8C7h0w31HAKcDZVXV3s43X7f7Xltqx\nEEnj90/A95JsSvLSJCvnse4fAi8EfgI4AHgb8FiSQ4FP0zvKOQA4B/gfSZ65gO1bC9wOvKs5NXdD\nkp9dQB5pQSxE0phV1QPAcUABHwbuTnJ5koOGrZfkScAvAb9WVXc1Ryp/X1UPA78I/E1V/U1VPVZV\nVwBbgZctYBMPA44GvgkcArwJ2JTkeQvIJc2bhUiagKq6qarOrKqZ/+gfAvzRiNWeASwHvjJg2eHA\nac1pufuT3E+v2B28gM37NvAIcF5Vfbeq/g9wJbBuAbmkebMQSRNWVV8GPkavIA1zD/Ad4NkDlt0B\nbK6qp/d99q+q31/AJg26VuVj+TUxFiJpzJIcleTXkxzW/P4s4Azg88PWq6rHgI3Af0lySJJ9krw4\nyTLgz4CfSXJSM748yfEzcwzYhn2TLAf2AWbiZ7pmrwa+Bmxo4n4S+PfAlkX4+tJIFiJp/B4EjgWu\nSfIQvQJ0I/DrLdY9B7gBuJZeW/X7gCdV1R3Ay4HfBO6md4T0G8z97/Rv0zsF9w5615e+3YxRVY80\nuV5G7zrRh4HXNkdu0tjFF+NJkrrkEZEkqVMWIklSpyxEkqROWYgkSZ3yoactrFixog46aOhN8PP2\n8MO9510uW7Zs4jnarjcqbq7ls8fn8/s4li2lOZbStu4tc8w1Nmx81LLdid2ddcaRY1jur3/96/dU\n1cjHTlmIWjjooIP46le/uqg5P/axjwFw5plnTjxH2/VGxc21fPb4fH4fx7KlNMdS2ta9ZY65xoaN\nj1q2O7G7s844cgzLfdZZZ93eJtZTc5KkTlmIJEmdshBJkjplIZIkdcpCJEnqlIVIktQpC5EkqVMW\nIklSpyxEkqROWYgkSZ2yEEmSOmUhkiR1yleFt3DIIYfUe9/73kXNuWPHDgCOOOKIiedou96ouLmW\nzx6fz+/jWLaU5lhK27q3zDHX2LDxUct2J3Z31hlHjmG53/Wud11XVWtGxXpEJEnqlK+BaGHZsmWL\n/ph0XwOx570uYJryOIevgZhUjlG52/CISJLUKQuRJKlTFiJJUqcsRJKkTlmIJEmdGlshSvLcJNf3\nfR5I8pYkpyXZnuSxJGv64tf2xW5LcuqAnJcnubHv9zOT3N233i/3LVuf5Jbms75vfFWSa5LcmuSS\nJPuNax9IkkYbW/t2Vd0MrAZIsg9wF3AZsAJ4JXDBrFVuBNZU1aNJDga2JflUVT3a5Hgl8K0BU11S\nVW/qH0hyAHAusAYo4Lokl1fVfcD7gA9U1cVJzgdeB3xoUb60JGneJnVq7kTgK1V1e1Xd1BSp71NV\nu2aKDrCcXgEBIMlTgLcC57Wc7yTgiqra2RSfK4CTkwQ4Abi0idsEvGJB30iStCgm8oifJBuBL1TV\nn/SNXQWcU1Vb+8aOBTYChwOvqarLmvEPAFcDXwT+uqqObsbPBH4PuAe4GfhPVXVHknOA5VV1XhP3\nO8C3gY8Bn6+q5zTjzwI+M5NvLj7iZ3Ccj/iZ3jzO4SN+JpVjWO6pecRPcw3mFODjo2Kr6pqqej5w\nDLAhyfIkq4FnzxSlWT4FHFFVL6B31LNpEbf77CRbk2zdtWvXYqWVJM0yiUf8vJTe0dA32q5QVTcl\n+RZwNL2itCbJDnrb+0NJrqqq46vq3r7VPgL8QfPzXcDxfcsOA64C7gWenmTf5jTgYU3soG24ELgQ\nYNWqVeUjftov9xE/3edxDh/xM6kco3K3MYlrRGcAF40KarrZ9m1+Phw4CthRVR+qqkOq6gjgOOCf\nqur4Ju7gvhSnADc1P28B1iVZmWQlsA7YUr3zkFcCr2ri1gOf3M3vJ0naDWMtREn2B14CfKJv7NQk\ndwIvBj6dZEuz6Dh6nXLX0+uue2NV3TNiijc3reDbgDcDZwJU1U7gPcC1zefdzRjA24G3JrkVOBD4\n6O5/U0nSQo311FxVPUTvP/b9Y5fRKzSzYzcDm0fk20HvdN3M7xuADXPEbqTX+DB7/DZg7eitlyRN\ngk9WkCR1ykIkSeqUhUiS1KmJ3NC61HlD6+A4b2id3jzO4Q2tk8oxLPfU3NAqSdIwk7ihdclbtmzZ\not/w5Q2te96Nj9OUxzm8oXVSOUblbsMjIklSpyxEkqROWYgkSZ2ya64Fu+YGx9k1N715nMOuuUnl\nGJbbrjlJ0pJg11wLds0NjrNrbnrzOIddc5PKMSp3Gx4RSZI6ZSGSJHXKZoUWbFYYHGezwvTmcQ6b\nFSaVY1humxUkSUuCzQot2KwwOM5mhenN4xw2K0wqx6jcbXhEJEnqlNeIWvAa0eA4rxFNbx7n8BrR\npHIMy+01IknSkuA1oha8RjQ4zmtE05vHObxGNKkco3K34am5Fjw1NzjOU3PTm8c5PDU3qRzDcntq\nTpK0JHhqrgVPzQ2O89Tc9OZxDk/NTSrHqNxteGquBU/NDY7z1Nz05nEOT81NKsew3J6akyQtCRYi\nSVKnLESSpE5ZiCRJnbIQSZI6Zft2C7ZvD46zfXt68ziH7duTyjEqdxu2b7dg+/bgONu3pzePc9i+\nPakcw3Lbvi1JWhIsRJKkTlmIJEmdshBJkjplIZIkdcr27RZs3x4cZ/v29OZxDtu3J5VjVO42xta+\nneS5wCV9Q0cC/xm4C3gn8DxgbVVtbeLXAhfOrA68s6oum5XzcuDIqjq6+X0Z8KfAC4F7gZ+vqh3N\nsvXAbzernldVm5rxVcDFwIHAdcBrquq7w76L7duD42zfnt48zmH79qRyDMvdeft2Vd1cVaurajW9\nQrELuAy4EXglcPWsVW4E1jTxJwMXJHn8iC3JK4FvzVrndcB9VfUc4APA+5rYA4BzgWOBtcC5SVY2\n67wP+ECzzn1NDklSRyZ1jehE4CtVdXtV3VRVN88OqKpdVfVo8+ty4PFDtSRPAd4KnDdrtZcDm5qf\nLwVOTBLgJOCKqtpZVfcBVwAnN8tOaGJp1n3FonxDSdKCTKoQnQ5cNCooybFJtgM3AG/oK0zvAd5P\n76iq36HAHQBN7DfpnXJ7fLxxZzN2IHB/X96Z8UHbcnaSrUm27to1e1pJ0mIZeyFKsh9wCvDxUbFV\ndU1VPR84BtiQZHmS1cCzZ18vGrequrCq1lTVmhUrVkxyaknaq0ziiOilwBeq6httV6iqm+hdDzoa\neDGwJskO4G+BH0lyVRN6F/AsgOZ60tPoNS08Pt44rBm7F3h637WnmXFJUkcm0b59Bu1Oy60C7qiq\nR5McDhwF7Gi66j7UxBwB/HVVHd+sdjmwHvgH4FXA56qqkmwB3tvXoLAO2NAsu7KJvbhZ95Ojts32\n7cFxtm9Pbx7nsH17UjlG5W5jrE/fTrI/8DV6LdffbMZOBT4IPBO4H7i+qk5K8hrgHcAjwGPAu6vq\nr2blO4JeIZpp314ObAZ+HNgJnF5VtzXLfgn4zWbV362q/96MH0mvCB0AfBH4xap6eNj3sH17cJzt\n29Obxzls355UjmG527Zvj/WIqKoeotcg0D92Gb027tmxm+kVlWH5dtA7XTfz+3eA0+aI3QhsHDB+\nG72WbknSFPARP5KkTlmIJEmdshBJkjplIZIkdcqnb7dg+/bgONu3pzePc9i+Pakco3K3Mdb27T2F\n7duD42zfnt48zmH79qRyDMvd+dO3JUlqw0IkSeqUhUiS1CkLkSSpUxYiSVKnbN9uwfbtwXG2b09v\nHuewfXtSOUblbsP27RZs3x4cZ/v29OZxDtu3J5VjWG7btyVJS4KFSJLUKQuRJKlTFiJJUqcsRJKk\nTlmIJEmdshBJkjrlDa0teEPr4DhvaJ3ePM7hDa2TyjEqdxve0NqCN7QOjvOG1unN4xze0DqpHMNy\ne0OrJGlJsBBJkjplIZIkdcpCJEnqlIVIktQp27dbsH17cJzt29Obxzls355UjlG527B9uwXbtwfH\n2b49vXmcw/btSeUYltv2bUnSkjCyECX51SQrJ7ExkqS9T5sjooOAa5P8ZZKTk2TcGyVJ2nuMLERV\n9dvAvwI+CpwJ3JLkvUmePeZtkyTtBVpdI6peR8M/N59HgZXApUn+YIzbJknaC4xs307ya8BrgXuA\njwC/UVWPJHkScAvwtvFuYvds3x4cZ/v29OZxDtu3J5VjVO42RrZvJ3kXsLGqbh+w7HlVddN8N3Cp\nsX17cJzt29Obxzls355UjmG527Zvjzwiqqpzhyzb44uQJGm8xnYfUZLnJrm+7/NAkrckOS3J9iSP\nJVnTF7+2L3ZbklP7ln22Gdue5Pwk+zTjZya5u2+9X+5bZ32SW5rP+r7xVUmuSXJrkkuS7DeufSBJ\nGm1shaiqbq6q1VW1GnghsAu4DLgReCVw9axVbgTWNPEnAxckmTli+7mq+jHgaOCZwGl9610yM09V\nfQQgyQHAucCxwFrg3L57od4HfKCqngPcB7xuUb+4JGleJvVkhROBr1TV7VV1U1XdPDugqnZV1aPN\nr8uB6lv2QPPjvsB+/cvmcBJwRVXtrKr7gCuAmXugTgAubeI2Aa9Y6JeSJO2+SRWi04GLRgUlOTbJ\nduAG4A19hYkkW4B/AR7kiUIC8LNJbkhyaZJnNWOHAnf0xdzZjB0I3N+Xd2Z80LacnWRrkq27du1q\n9SUlSfM39qdvN9dgTgE2jIqtqmuA5yd5HrApyWeq6jvNspOSLAf+nN5RzRXAp4CLqurhJK+nd4Rz\nwmJsd1VdCFwIsGrVqrJ9u/1y27e7z+Mctm9PKseo3G2M/enbSV4O/EpVrZs1fhVwTlVtnWO9zwFv\nm708yWuBtVX1plnj+wA7q+ppSc4Ajq+q1zfLLgCuAi4G7gZ+uKoeTfJi4J1VddKw72D79uA427en\nN49z2L49qRzDck/T07fPoN1puVUzzQlJDgeOAnYkeUqSg5vxfYGfBr7c/H5wX4pTgJl28i3AuiQr\nmyaFdcCW5gkRVwKvauLWA5/cze8nSdoNYz01l2R/4CXA6/vGTgU+SK/77dNJrm+OSI4D3pHkEeAx\n4I1VdU+Sg4DLkyyjVzivBM5v0r05ySn0Hju0k96z8KiqnUneA1zbxL27qnY2P78duDjJecAX6T1D\nT5LUkbEWoqp6iF6DQP/YZfTauGfHbgY2Dxj/BnDMHPk3MMe1p6raCGwcMH4bvZZuSdIU8MV4kqRO\nWYgkSZ0ae/v2nsCnbw+Os317evM4h+3bk8oxKncbY2/f3hPYvj04zvbt6c3jHLZvTyrHsNzT1L4t\nSdKcLESSpE5ZiCRJnbIQSZI6ZSGSJHXK9u0WbN8eHGf79vTmcQ7btyeVY1TuNmzfbsH27cFxtm9P\nbx7nsH17UjmG5bZ9W5K0JFiIJEmdshBJkjplIZIkdcpCJEnqlO3bLdi+PTjO9u3pzeMctm9PKseo\n3G3Yvt2C7duD42zfnt48zmH79qRyDMtt+7YkaUmwEEmSOmUhkiR1ykIkSeqUhUiS1Cnbt1uwfXtw\nnO3b05vHOWzfnlSOUbnbsH27Bdu3B8fZvj29eZzD9u1J5RiW2/ZtSdKSYCGSJHXKQiRJ6pSFSJLU\nKQuRJKlTFiJJUqcsRJKkTnlDawve0Do4zhtapzePc3hD66RyjMrdhje0tuANrYPjvKF1evM4hze0\nTirHsNze0CpJWhIsRJKkTlmIJEmdshBJkjo1tkKU5LlJru/7PJDkLUlOS7I9yWNJ1vTFr+2L3Zbk\n1L5ln23Gtic5P8k+zfiyJJckuTXJNUmO6FtnfZJbms/6vvFVTeytzbr7jWsfSJJGm0jXXFM47gKO\nBVYAjwEXAOdU1dYmZgXw3ap6NMnBwDbgkOb3H6yqB5IEuBT4eFVdnOSNwI9W1RuSnA6cWlU/n+QA\nYCuwBijgOuCFVXVfkr8EPtGsfz6wrao+NGz7V61aVV/96lcXdZ/Yvr3ntflOUx7nsH17UjmG5T7r\nrLNadc1NqhCtA86tqp/sG7uKvkI0K34V8Hng0Kp6tG/8B4BPAH9WVZck2QK8s6r+Icm+wD8DzwRO\nB46vqtc3610AXAVcDNwN/HBT4F7crH/SsO23fXtwnO3b05vHOWzfnlSOYbmnrX37dOCiUUFJjk2y\nHbgBeMOsIrQF+BfgQXpHRQCHAncANLHfBA7sH2/c2YwdCNzfl3dmfNC2nJ1ka5Ktu3btavs9JUnz\nNPYnKzTXYE4BNoyKraprgOcneR6wKclnquo7zbKTkiwH/hw4AbhijJtNVV0IXAi9U3OLfejqqbk9\n7xTONOVxDk/NTSrHqNxtjP3UXJKXA79SVetmjV/FHKfmmuWfA942e3mS1wJrq+pNnppbWA5PzU3v\nHEtpW/eWOeYaGzY+atnuxO7OOuPIMSz3NJ2aO4N2p+VWNcWEJIcDRwE7kjylaV6gWf7TwJeb1S4H\nZjriXgV8rnqVdQuwLsnKJCuBdcCWZtmVTSzNup9chO8oSVqgsZ6aS7I/8BLg9X1jpwIfpHfk8ukk\n1zdHJMcB70jyCL2uujdW1T1JDgIuT7KMXuG8Eji/SfdRYHOSW4Gd9I6EqKqdSd4DXNvEvbuqdjY/\nvx24OMl5wBebHEP50NPBcZ6am948zuGpuUnlGJW7jbEWoqp6iF6DQP/YZcBlA2I3A5sHjH8DOGaO\n/N8BTptj2UZg44Dx24C1LTZfkjQBPn27Ba8RDY7zGtH05nEOrxFNKsew3NN0jUiSpDn5YrwWvEY0\nOM5rRNObxzm8RjSpHKNyt+ERkSSpUxYiSVKnbFZowWaFwXE2K0xvHuewWWFSOYbltllBkrQk2KzQ\ngs0Kg+NsVpjePM5hs8KkcozK3YZHRJKkTlmIJEmdshBJkjpl11wLds0NjrNrbnrzOIddc5PKMSy3\nXXOSpCXBrrkW7JobHGfX3PTmcQ675iaVY1TuNjwikiR1ykIkSeqUhUiS1CkLkSSpUxYiSVKnLESS\npE55Q2sL3tA6OM4bWqc3j3N4Q+ukcgzL7Q2tkqQlwRtaW/CG1sFx3tA6vXmcwxtaJ5VjVO42PCKS\nJHXKQiRJ6pSFSJLUKQuRJKlTFiJJUqcsRJKkTlmIJEmdshBJkjrlI35a8BE/g+N8xM/05nEOH/Ez\nqRzDcvuIH0nSkuAjflrwET+D43zEz/TmcQ4f8TOpHKNyt+ERkSSpUxYiSVKnLESSpE6NrRAleW6S\n6/s+DyR5S5LTkmxP8liSNX3xa/tityU5tRlfkeTTSb7crPf7feucmeTuvvV+uW/Z+iS3NJ/1feOr\nklyT5NYklyTZb1z7QJI02tgKUVXdXFWrq2o18EJgF3AZcCPwSuDqWavcCKxp4k8GLkgy00zxh1V1\nFPDjwE8meWnfepfMzFNVHwFIcgBwLnAssBY4N8nKJv59wAeq6jnAfcDrFvebS5LmY1Kn5k4EvlJV\nt1fVTVV18+yAqtpVVY82vy4Hqm/8yubn7wJfAA4bMd9JwBVVtbOq7gOuAE5OEuAE4NImbhPwit38\nbpKk3TCpQnQ6cNGooCTHJtkO3AC8oa8wzSx/OvAzwP/uG/7ZJDckuTTJs5qxQ4E7+mLubMYOBO7v\nyzszPmhbzk6yNcnWBx98cPQ3lCQtyNgLUXMN5hTg46Niq+qaqno+cAywIcnyvjz70itmf1xVtzXD\nnwKOqKoX0Dvq2bRY211VF1bVmqpa89SnPnWx0kqSZpnEEdFLgS9U1TfarlBVNwHfAo7uG74QuKWq\n/qgv7t6qerj59SP0rkUB3AU8q2/dw5qxe4Gn9117mhmXJHVkEoXoDNqdlls1UyCSHA4cBexofj8P\neBrwllnrHNz36ynATc3PW4B1SVY2TQrrgC3Ve7DelcCrmrj1wCcX9rUkSYthrIUoyf7AS4BP9I2d\nmuRO4MXAp5NsaRYdB2xLcj297ro3VtU9SQ4Dfgv418AXZrVpv7lp6d4GvBk4E6CqdgLvAa5tPu9u\nxgDeDrw1ya30rhl9dExfX5LUwlifNVdVD9H7j33/2GX0Cs3s2M3A5gHjdwKZI/8GYMMcyzYCGweM\n30avpVuSNAV8DUQLvgZicJyvgZjePM7hayAmlWNY7ravgbAQtZDkbuD2rrdjkT0DuKfrjZhS7pvh\n3D9zc998v8Or6pmjgixEe6kkW9v8n8reyH0znPtnbu6bhfGhp5KkTlmIJEmdshDtvS7segOmmPtm\nOPfP3Nw3C+A1IklSpzwikiR1ykIkSeqUhWgJSrKjefXF9Um2NmPvTHJX39tqX9YX/6NJ/qF5HNIN\nM081T3JVkpv71vmhZnxZ8/baW5u32R7Rl2vgm2+nxXz2TZJXz3qL8GNJVjfLXtjkuTXJHzfvsnLf\nsGf+3cC8988PJNnUxN+UZENfnj3ub2fsqsrPEvvQexjsM2aNvRM4Z0DsvsCXgB9rfj8Q2Kf5+Sp6\nb8Wdvc4bgfObn0+n9xZcgAOA25p/rmx+Xtn1/ljovpkV8wJ6L2+c+f3/Ai+i93ipzwAvdd88/vse\n93cz3/0D/AJwcfPzimbdI/bUv51xfzwi2vOtA75UVdvg8VdnfG/EOi/niXc7XQqc2Pxf3cA3345p\nuyftDOBiePyp7j9YVZ+v3n8p/pQn3uS7V++bEfamfVPA/um9MeDJwHeBB/zbWRgL0dJUwP9Kcl2S\ns/vGfzXJl5JsTO/1FwA/AlSSLUm+kORts3Jtak45/M7MKQT63nBbvbfZfpPekdRcb76dJvPZN/1+\nnideV3Iove82o/977u37Zsae9ncD89s/lwIPAV8Hvgb8YfWe8L+n/u2MlYVoaTquqlbTe+ngryT5\nt8CHgCOB1fT+5Xh/E7svvVdsvLr556lJTmyWvbp6b8T9qebzmsl9hbGZz74Beq+oB3ZV1Y2T3tgJ\nW6x9syf+3cD89s9a4HvAIcAq4NeTHDn5Td4zWIiWoKq6q/nnv9B7pcbaqvpGVX2vqh4DPswTr7q4\nE7i6qu6pql3A3wD/ZlaeB4G/6Fvn8TfcNqcenkbv7bZzvfl2asxz38w4ne//P/676H23Gf3fc2/f\nN3vk3w3Me//8AvDZqnqkif87YA176N/OuFmIlpgk+yd56szP9K4B3Zjvf1vtqcDM/8FuAV6QZEXz\nx//vgH9Msm+SZzR5fgD4D33rXE7v7bXQe5vt55rz3QPffDuu7zpfC9g3JHkS8HP0XQOpqq/TO9//\noua002t54k2+e/W+2RP/bmBB++drwAl98S8Cvrwn/u1MRNfdEn7m96F3mmBb89kO/FYzvhm4gV6H\n3OXAwX3r/GITeyPwB83Y/sB1Tfx24L/yRDfdcuDjwK30OoCO7Mv1S834rcBZXe+PRdg3xwOfH5Br\nTbO/vgL8CU88hWSv3jd74t/NQvYP8JTmu24H/hH4jT31b2cSHx/xI0nqlKfmJEmdshBJkjplIZIk\ndcpCJEnqlIVIktQpC5EkqVMWIklSpyxE0hKU5JjmQZzLm6cCbE9ydNfbJS2EN7RKS1SS8+jdrf9k\n4M6q+r2ON0laEAuRtEQl2Q+4FvgO8BM1+j1T0lTy1Jy0dB1I75lnT6V3ZCQtSR4RSUtUksvpPRl7\nFb2Hcb6p402SFmTfrjdA0vwleS3wSFX9RZJ9gL9PckJVfa7rbZPmyyMiSVKnvEYkSeqUhUiS1CkL\nkSSpUxYiSVKnLESSpE5ZiCRJnbIQSZI69f8Bk84YVIkvImIAAAAASUVORK5CYII=\n",
+      "text/plain": [
+       "<matplotlib.figure.Figure at 0x10d3c4710>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZEAAADhCAYAAAAJdIgiAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAFu1JREFUeJzt3X/wZXV93/HnS9AFfyAoUdld6n5pYTqAiRMUaaqNCVbR\nmmBmrN2oEdTIVDDKjEZFpmOYEWqTaI2xmtlUBogxSNWWbQslgmNTYgDXjAqLP7rKKuysvyHQkCxh\nefePe1YPX8733rtnv/d79977fMx8h/v5nM/5vM/3zGU+e97vc843VYUkSX08atoHIEmaXS4ikqTe\nXEQkSb25iEiSenMRkST15iIiSerNRUTqKcnZSW5stf9fkuOmeUzSWnMRkYZI8twkn0/yN0l+nOQv\nkzy7a2xVPb6qvrXK8T+W5LtJ7k3yjSS/2dr2mCSfTLIzSSV5/mrGlsbhIiKtIMkRwP8A/hB4ErAB\nuAjYs4aH8V7guKo6AvhV4D1JTmltvxF4NfDdNTwm6SdcRKSVnQBQVX9WVXur6u+q6s+r6itdg5ur\ngX/SfD48yfuSfLu5irkxyeHNttOaq5t7knx52BVEVd1WVffvazY//7jZ9kBVfaCqbgT2rtpvLe0H\nFxFpZd8A9ia5PMmLkxy1H/v+PnAK8AsMrmLeDjyUZAPwP4H3NP1vAz6V5GdWmijJh5PcD3wN2A1c\n0+u3kSbARURaQVXdCzyXwb/+/xj4QZKtSZ46bL8kjwJeB7ylqnY1VzGfr6o9DFJP11TVNVX1UFV9\nBtgGvGTIcZwLPAF4HvBp1jadJg3lIiINUVVfraqzq2ojcDKwHvjAiN2OBg4Dvtmx7enAv25SWfck\nuYfBQnXMiOPY26StNgJv3N/fQ5oUFxFpTFX1NeAyBovJMD8E/p6mdrHMncCfVNWRrZ/HVdV7xzyM\nQ1eYV5oKFxFpBUn+aZK3JtnYtI8Ffh24adh+VfUQcCnw/iTrkxyS5J8lWQd8DPiVJC9q+g9L8vx9\nMZbFf0qSzUke34x9URP/htaYdUkOa5qPaebLqpwAaQwuItLK7gOeA9yc5G8ZLB63AW8dY9+3AbcC\nXwB+DPwH4FFVdSdwJvAu4AcMrkx+m+7/F4tB6uou4G4Gxfrzq2pra8zXgb9jcPvxdc3np+/Xbykd\ngPhHqSRJfXklIknqzUVEktSbi4gkqTcXEUlSby4ikqTeDp32ASQ5hMFrH3ZV1UuTPAn4BLAJ2Am8\noqrubsZeALyewcvm3lxV142a/7GPfWw99alD31Kx3/bsGbx1Yt26dWs+x7j7jRq30vbl/fvTnsQ2\nYxhjteZZru+2Axl7IPtMYo5hc+/evfuHVbXiO932mfoiArwF+CpwRNN+J3BDVb03yTub9juSnAhs\nBk5i8OqJ65OcUFVD31565JFH8u53v3tVD3jnzp0AbNq0ac3nGHe/UeNW2r68f3/ak9hmDGOs1jzL\n9d12IGMPZJ9JzDFs7osuuujb44ydajqreUr3XwH/udV9JnB58/ly4GWt/iurak9V3QHsAE5dq2OV\nJD3StGsiH6B5RXar76lVtbv5/F1gXy5qA4One/e5q+mTJE3J1J5YT/JS4CVVdW7zR3ne1tRE7qmq\nI1vj7q6qo5J8CLipqj7W9H8UuLaqPtkx9znAOQBPfOITTzn//PPX4leSpLlx0UUXfbGqnjVq3DRr\nIv8c+NUkL2Hw2uwjknwM+F6SY6pqd5JjgO8343cBx7b239j0PUJVbQG2AKxfv75WO2doTWT28uPG\nWLwYXe02ayKj5x7H1NJZVXVBVW2sqk0MCuafrapXA1uBs5phZwFXN5+3Apubt5YuAccDt6zxYUuS\nWg6Gu7OWey9wVZLXA98GXgFQVduTXAXcDjwInDfqzixJ0mQdFItIVX0O+Fzz+UfA6SuMuxi4eM0O\nTJI01LTvzpIkzTAXEUlSby4ikqTeDoqayCStW7eOs88+e1XnvOyyywAOaN6+c4y736hxK21f3r8/\n7UlsM4YxVmue5fpuO5CxB7LPJOYYNfc45v7P465fv74uueSSVZ3T50Rm75kBYyxejK52m8+JDJ97\n3IcNTWdJknozndWD6azZS20YY/FidLXbTGeNnnscXolIknqzJtKDNZHZy48bY/FidLXbrIkMn9ua\niCRp4qyJ9GBNZPby48ZYvBhd7TZrIqPnHofprB5MZ81easMYixejq91mOmv43KazJEkTZzqrB9NZ\ns5faMMbixehqt5nOGj33OLwSkST1Zk2kB2sis5cfN8bixehqt1kTGT63NRFJ0sRZE+nBmsjs5ceN\nsXgxutpt1kRGzz0Or0QkSb1NrSaS5FjgCuCpQAFbquoPkjwJ+ASwCdgJvKKq7m72uQB4PbAXeHNV\nXTcqjjWR7nHWRIwx7zG62m3WRIbPPQs1kQeBt1bVicBpwHlJTgTeCdxQVccDNzRtmm2bgZOAM4AP\nJzlkKkcuSQKmWBOpqt3A7ubzfUm+CmwAzgSe3wy7HPgc8I6m/8qq2gPckWQHcCrwV8PiWBPpHmdN\nxBjzHqOr3WZNZPTc4zgobvFNsgn4C+Bk4DtVdWTTH+DuqjoyyYeAm6rqY822jwLXVtUnh81tOqt7\nnOksY8x7jK52m+ms4XPPQjoLgCSPBz4FnF9V97a31WCF2+9VLsk5SbYl2Xb//fev0pFKkpab6i2+\nSR7NYAH506r6dNP9vSTHVNXuJMcA32/6dwHHtnbf2PQ9QlVtAbYALC0tlems8bebzjLGvMToareZ\nzho99zimdiXSpKo+Cny1qt7f2rQVOKv5fBZwdat/c5J1SZaA44Fb1up4JUmPNM1bfJ8L/B/gVuCh\npvtdwM3AVcA/Ar7N4BbfHzf7XAi8jsGdXedX1bWj4lgT6R5nTcQY8x6jq91mTWT43OPWRKZ5d9aN\nQFbYfPoK+1wMXDyxg5Ik7Rdfe9KDNZHZy48bY/FidLXbrImMnnscB8UtvpNkOqt7nOksY8x7jK52\nm+ms4XPPzC2+kqTZZTqrB9NZs5faMMbixehqt5nOGj33OLwSkST1Zk2kB2sis5cfN8bixehqt1kT\nGT63NRFJ0sRZE+nBmsjs5ceNsXgxutpt1kRGzz0O01k9mM6avdSGMRYvRle7zXTW8LlNZ0mSJs50\nVg+ms2YvtWGMxYvR1W4znTV67nF4JSJJ6s2aSA/WRGYvP26MxYvR1W6zJjJ8bmsikqSJsybSgzWR\n2cuPG2PxYnS126yJjJ57HF6JSJJ6sybSgzWR2cuPG2PxYnS126yJDJ/bmogkaeKsifRgTWT28uPG\nWLwYXe02ayKj5x6H6aweTGfNXmrDGIsXo6vdZjpr+Nxzm85KckaSryfZkeSd0z4eSVpkM5XOSnII\n8J+AfwncBXwhydaqun2lfUxndY8znWWMeY/R1W4znTV67nHM2pXIqcCOqvpWVT0AXAmcOeVjkqSF\nNVM1kSQvB86oqt9s2r8BPKeq3rTSPtZEusdZEzHGvMfoardZExk+99zWRMaR5Jwk25Jsu//++6d9\nOJI0t2aqJgLsAo5ttTc2fQ9TVVuALQBLS0tlTWT87dZEjDEvMbrabdZERs89jllLZx0KfAM4ncHi\n8QXglVW1faV9TGd1jzOdZYx5j9HVbjOdNXzucdNZM3UlUlUPJnkTcB1wCHDpsAVEkjRZM7WIAFTV\nNcA14473Ft/ucaazjDHvMbrabaazRs89jrksrEuS1sZM1UT6sCbSPc6aiDHmPUZXu82ayPC5F/oW\nX0nS2pi5msj+sibSPc6aiDHmPUZXu82ayOi5x2E6qwfTWbOX2jDG4sXoareZzho+t+ksSdLEmc7q\nwXTW7KU2jLF4MbrabaazRs89Dq9EJEm9WRPpwZrI7OXHjbF4MbrabdZEhs9tTUSSNHHWRHqwJjJ7\n+XFjLF6MrnabNZHRc4/DKxFJUm/WRHqwJjJ7+XFjLF6MrnabNZHhc1sTkSRNnDWRHqyJzF5+3BiL\nF6Or3WZNZPTc4xgrnZXkBuB9zd/y2Ne3parO6XOAa8l0Vvc401nGmPcYXe0201nD517tdNYS8I4k\n7271jZxckjTfxk1n3cPg75p/MMl/B149uUNaXaazuseZzjLGvMfoareZzho99zjGvRJJVT1YVecC\nnwJuBJ6y/4cmSZon416J/NG+D1V1WZJbgfMmc0iSpFkxledEkvwe8CvAA8A3gddW1T3NtguA1wN7\ngTdX1XVN/ynAZcDhwDXAW2qMg7ew3j3Owrox5j1GV7vNwvrwuQ/250Q+A5xcVT8LfAO4ACDJicBm\n4CTgDODDSQ5p9vkI8Abg+ObnjLU+aEnSw03lOZGq+vNW8ybg5c3nM4Erq2oPcEeSHcCpSXYCR1TV\nTQBJrgBeBlw7KpaF9e5xFtaNMe8xutptFtZHzz2Og+GJ9dfx08VgA3Bna9tdTd+G5vPy/k5Jzkmy\nLcm2++67b5UPV5K0z8SuRJJcDzytY9OFVXV1M+ZC4EHgT1czdlVtAbYALC0tzffLwSRpiia2iFTV\nC4ZtT3I28FLg9FaBfBdwbGvYxqZvV/N5eb8kaYqmdXfWGcD7gV+sqh+0+k8CPg6cCqwHbgCOr6q9\nSW4B3gzczODurD9sv4ZlJd6d1T3Ou7OMMe8xutpt3p01fO5x786a1gsYPwSsAz6TBOCmqvq3VbU9\nyVXA7QzSXOdV1d5mn3P56S2+1zJGUV2SNFlz//dElpaW6o477ljVOb07a/bu1DHG4sXoard5d9bw\nuV/72tce1M+JSJLmgIuIJKk3FxFJUm8uIpKk3lxEJEm9zf3dWT4n0j3O50SMMe8xutptPicyfO6D\n/S2+kqQ5MK2HDdeMb/HtHudzIsaY9xhd7TafExk99zi8EpEk9eYiIknqzUVEktSbi4gkqTcXEUlS\nby4ikqTeXEQkSb35xHoPPrE+e08vG2PxYnS123xiffjcPrEuSZo4n1jvwSfWZ+/pZWMsXoyudptP\nrI+eexxeiUiSepvqIpLkrUkqydGtvguS7Ejy9SQvavWfkuTWZtsHk2Q6Ry1J2mdqi0iSY4EXAt9p\n9Z0IbAZOAs4APpzkkGbzR4A3AMc3P2es6QFLkh5hmlci/xF4O9C+PexM4Mqq2lNVdwA7gFOTHAMc\nUVU31eB2siuAl635EUuSHmYqi0iSM4FdVfXlZZs2AHe22nc1fRuaz8v7JUlTNLG7s5JcDzytY9OF\nwLsYpLImFfsc4ByAJz/5yZMKI0kLb2KLSFW9oKs/yTOAJeDLTW18I/DXSU4FdgHHtoZvbPp2NZ+X\n968UewuwBWBpaWm+n6aUpCla83RWVd1aVU+pqk1VtYlBaurnq+q7wFZgc5J1SZYYFNBvqardwL1J\nTmvuynoNcPVaH7sk6eEOqocNq2p7kquA24EHgfOqam+z+VzgMuBw4NrmR5I0RVNfRJqrkXb7YuDi\njnHbgJPX6LAkSWPwBYw9+ALG2XsZnzEWL0ZXu80XMA6f2xcwSpImburprEnzBYzd43wBozHmPUZX\nu80XMI6eexxeiUiSenMRkST15iIiSerNRUSS1JuLiCSpNxcRSVJvLiKSpN5cRCRJvbmISJJ6cxGR\nJPXmIiJJ6s1FRJLUm4uIJKk3FxFJUm/+Uaoe/KNUs/cHioyxeDG62m3+Uarhc/tHqSRJE+cfperB\nP0o1e3+gyBiLF6Or3eYfpRo99zimdiWS5LeSfC3J9iS/2+q/IMmOJF9P8qJW/ylJbm22fTBJpnPk\nkqR9pnIlkuSXgDOBn6uqPUme0vSfCGwGTgLWA9cnOaGq9gIfAd4A3AxcA5wBXDuN45ckDUzrSuSN\nwHurag9AVX2/6T8TuLKq9lTVHcAO4NQkxwBHVNVNNbgT4ArgZdM4cEnST01rETkBeF6Sm5P87yTP\nbvo3AHe2xt3V9G1oPi/v75TknCTbkmy77777VvnQJUn7TCydleR64Gkdmy5s4j4JOA14NnBVkuNW\nK3ZVbQG2ACwtLc33PcySNEUTW0Sq6gUrbUvyRuDTTWrqliQPAUcDu4BjW0M3Nn27ms/L+yVJUzSt\ndNZ/A34JIMkJwGOAHwJbgc1J1iVZAo4Hbqmq3cC9SU5r7sp6DXD1dA5dkrTPtJ4TuRS4NMltwAPA\nWc1VyfYkVwG3Aw8C5zV3ZgGcC1wGHM7grizvzJKkKZvKIlJVDwCvXmHbxcDFHf3bgJMnfGiSpP3g\na08kSb25iEiSenMRkST15iIiSerNRUSS1JuLiCSpNxcRSVJvLiKSpN5cRCRJvbmISJJ6cxGRJPWW\nwXsP59f69evrkksuWdU5d+7cCcCmTZvWfI5x9xs1bqXty/v3pz2JbcYwxmrNs1zfbQcy9kD2mcQc\nw+a+6KKLvlhVzxo1du4XkSQ/AL497eNYZUczeHW+HslzM5znZ2Wem4d7elX9zKhBc7+IzKMk28b5\nF8Ii8twM5/lZmeemH2sikqTeXEQkSb25iMymLdM+gIOY52Y4z8/KPDc9WBORJPXmlYgkqTcXkTWU\nZGeSW5N8Kcm2pu93kuxq+r6U5CWt8T+b5K+SbG/2O6zp/1ySr7f2eUrTvy7JJ5LsSHJzkk2tuc5K\n8n+bn7PW9jcfz/6cnySvavV9KclDSZ7ZbDulmWdHkg8mSdM/s+dnFc+N353k0Ukub8Z/NckFrXnm\n7rszcVXlzxr9ADuBo5f1/Q7wto6xhwJfAX6uaT8ZOKT5/DngWR37nAv8UfN5M/CJ5vOTgG81/z2q\n+XzUtM/HgZyfZWOeAXyz1b4FOA0IcC3w4lk/P6t4bhb+uwO8Eriy+fzYZt9N8/rdmfSPVyIHrxcC\nX6mqLwNU1Y+qau+Ifc4ELm8+fxI4vfmX1IuAz1TVj6vqbuAzwBkTOu5p+HXgSoAkxwBHVNVNNfi/\n/ArgZc24RTw/Pzk3IyzSuSngcUkOBQ4HHgDu9bvTj4vI2irg+iRfTHJOq/+3knwlyaVJjmr6TgAq\nyXVJ/jrJ25fNdXlzif7v9l1yAxuAOwGq6kHgbxhcwfykv3FX03ew2Z/z0/ZvgD9rPm9g8Pvt0/5d\nZ/n8rMa52WfRvzufBP4W2A18B/j9qvox8/vdmSgXkbX13Kp6JvBi4Lwk/wL4CHAc8EwGX+r3NWMP\nBZ4LvKr5768lOb3Z9qqqOgl4XvPzG2v3K0zU/pwfAJI8B7i/qm5b64NdY6t1bvzuwKnAXmA9sAS8\nNclxa3/I88FFZA1V1a7mv98H/itwalV9r6r2VtVDwB8z+ILD4F80f1FVP6yq+4FrgJ9fNs99wMdb\n++wCjgVoLtWfCPyo3d/Y2PQdVPbz/OyzmYf/S3sXg99vn/bvOrPnZ5XOjd+dgVcC/6uq/qEZ/5fA\ns5jT786kuYiskSSPS/KEfZ8Z1Dxua/Kw+/wasO9fjdcBz0jy2OZL+4vA7UkOTXJ0M8+jgZe29tkK\n7Ls75OXAZ5vc7nXAC5Mc1VzSv7DpO2j0OD8keRTwClo5/6razSC/fVqTqnkNcHWzeSbPz2qdG787\nP/ldvwP8cmv8acDX5vG7syamXdlflB8Gl9Vfbn62Axc2/X8C3MrgTqytwDGtfV7djL0N+N2m73HA\nF5vx24E/4Kd3bR0G/BdgB4O7TI5rzfW6pn8H8Nppn49VOj/PB27qmOtZzTn7JvAhfvpQ7Uyen9U6\nN353BucHeHzzu24Hbgd+e16/O2vx4xPrkqTeTGdJknpzEZEk9eYiIknqzUVEktSbi4gkqTcXEUlS\nby4ikqTeXESkNZTk2c0LAQ9rnrTenuTkaR+X1JcPG0prLMl7GDwBfThwV1X9+ykfktSbi4i0xpI8\nBvgC8PfAL9TovxMjHbRMZ0lr78kM3t/0BAZXJNLM8kpEWmNJtjJ4u+4Sg5cCvmnKhyT1dui0D0Ba\nJEleA/xDVX08ySHA55P8clV9dtrHJvXhlYgkqTdrIpKk3lxEJEm9uYhIknpzEZEk9eYiIknqzUVE\nktSbi4gkqTcXEUlSb/8fbTASWJQTSoQAAAAASUVORK5CYII=\n",
+      "text/plain": [
+       "<matplotlib.figure.Figure at 0x10d883390>"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "# Core cell sizes in x, y, and z\n",
+    "csx, csy, csz = 25., 25., 25.\n",
+    "\n",
+    "# Number of core cells in each directiPon s\n",
+    "ncx, ncy, ncz = 48, 48, 20\n",
+    "\n",
+    "# Number of padding cells to add in each direction\n",
+    "npad = 7\n",
+    "\n",
+    "# Vectors of cell lengthts in each direction\n",
+    "hx = [(csx,npad, -1.3),(csx,ncx),(csx,npad, 1.3)]\n",
+    "hy = [(csy,npad, -1.3),(csy,ncy),(csy,npad, 1.3)]\n",
+    "hz = [(csz,npad, -1.3),(csz,ncz), (csz/2.,6)]\n",
+    "\n",
+    "# Create mesh\n",
+    "mesh = Mesh.TensorMesh([hx, hy, hz],x0=\"CCN\")\n",
+    "\n",
+    "# Map mesh coordinates from local to UTM coordiantes\n",
+    "xc = 300+5.57e5\n",
+    "yc = 600+7.133e6\n",
+    "zc = 425.\n",
+    "mesh._x0 = mesh.x0 + np.r_[xc, yc, zc]\n",
+    "\n",
+    "# plot the mesh\n",
+    "mesh.plotSlice(np.ones(mesh.nC)*np.nan, grid=True)\n",
+    "mesh.plotSlice(np.ones(mesh.nC)*np.nan, grid=True, normal=\"Y\")\n",
+    "plt.gca().set_aspect('equal')\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Survey\n",
+    "\n",
+    "We use a simple gradient array having a pair of current electrodes (AB), and\n",
+    "multiple potential electrodes (MN).  The lengths of AB and MN electrodes are\n",
+    "1200 and 25 m, respectively.\n",
+    "\n",
+    "![Gradient array](../../../images/dc/GradientArray.png)\n",
+    "\n",
+    "\n",
+    "Once we have obtained locations of AB (Src) and MN (Rx) electrodes, we can\n",
+    "generate **Survey** class:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [],
+   "source": [
+    "# Create Src and Rx classes for DC problem\n",
+    "Aloc1_x = np.r_[-600., 0, 0.] + np.r_[xc, yc, zc]\n",
+    "Bloc1_x = np.r_[600., 0, 0.] + np.r_[xc, yc, zc]\n",
+    "\n",
+    "# Rx locations (M-N electrodes, x-direction)\n",
+    "x = mesh.vectorCCx[np.logical_and(mesh.vectorCCx>-300.+ xc, mesh.vectorCCx<300.+ xc)]\n",
+    "y = mesh.vectorCCy[np.logical_and(mesh.vectorCCy>-300.+ yc, mesh.vectorCCy<300.+ yc)]\n",
+    "# Grid selected cell centres to get M and N Rx electrode locations\n",
+    "Mx = Utils.ndgrid(x[:-1], y, np.r_[-12.5/2.])\n",
+    "Nx = Utils.ndgrid(x[1:], y, np.r_[-12.5/2.])\n",
+    "\n",
+    "rx = DC.Rx.Dipole(Mx, Nx)\n",
+    "src = DC.Src.Dipole([rx], Aloc1_x, Bloc1_x)\n",
+    "\n",
+    "# Form survey object using Srcs and Rxs that we have generated\n",
+    "survey = DC.Survey([src])"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Fields and Data\n",
+    "\n",
+    "By solving the DC equations, we compute electrical potential ($\\phi$) at\n",
+    "every cell. The **Problem** class does this, but it still requires survey\n",
+    "information hence we pair it to the **Survey** class:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 9,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [],
+   "source": [
+    "# Define problem and set solver\n",
+    "problem = DC.Problem3D_CC(mesh)\n",
+    "\n",
+    "problem.Solver = Solver\n",
+    "# Pair problem and survey\n",
+    "problem.pair(survey)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Here, we used ``DC.Problem3D_CC``, which means 3D space and $\\phi$ is\n",
+    "defined at the cell center. Now, we are ready to run DC forward modelling! For\n",
+    "this modelling, inside of the code, there are two steps:\n",
+    "\n",
+    "1. Compute fields ($\\phi$ at every cells)\n",
+    "2. Evaluate at Rx location (potential difference at MN electrodes)\n",
+    "\n",
+    "Consider two conductivity models:\n",
+    "\n",
+    "- Homogeneous background below topographic surface: ``sigma0``\n",
+    "  ($10^{-4}$ S/m)\n",
+    "\n",
+    "- Includes diamond pipes: ``sigma`` (S/m)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 10,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "IOError",
+     "evalue": "[Errno 2] No such file or directory: 'VTKout_DC.dat'",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mIOError\u001b[0m                                   Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-10-e6ad92636190>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m# Read pre-generated conductivity model in UBC format\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0msigma\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mmesh\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mreadModelUBC\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"VTKout_DC.dat\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      3\u001b[0m \u001b[0;31m# Identify air cells in the model\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0mairind\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0msigma\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m1e-8\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0;31m# Generate background model (constant conductiivty below topography)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+      "\u001b[0;32m/Users/lindseyjh/git/python_symlinks/discretize/MeshIO.pyc\u001b[0m in \u001b[0;36mreadModelUBC\u001b[0;34m(mesh, fileName)\u001b[0m\n\u001b[1;32m    205\u001b[0m         \u001b[0;34m:\u001b[0m\u001b[0;32mreturn\u001b[0m\u001b[0;34m:\u001b[0m \u001b[0mmodel\u001b[0m \u001b[0;32mwith\u001b[0m \u001b[0mTensorMesh\u001b[0m \u001b[0mordered\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    206\u001b[0m         \"\"\"\n\u001b[0;32m--> 207\u001b[0;31m         \u001b[0mf\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mopen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfileName\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'r'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    208\u001b[0m         \u001b[0mmodel\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0marray\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mlist\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmap\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfloat\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mf\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mreadlines\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    209\u001b[0m         \u001b[0mf\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mclose\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+      "\u001b[0;31mIOError\u001b[0m: [Errno 2] No such file or directory: 'VTKout_DC.dat'"
+     ]
+    }
+   ],
+   "source": [
+    "# Read pre-generated conductivity model in UBC format\n",
+    "sigma = mesh.readModelUBC(\"VTKout_DC.dat\")\n",
+    "# Identify air cells in the model\n",
+    "airind = sigma == 1e-8\n",
+    "# Generate background model (constant conductiivty below topography)\n",
+    "sigma0 = np.ones_like(sigma)*1e-4\n",
+    "sigma0[airind] = 1e-8"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "collapsed": true
+   },
+   "outputs": [],
+   "source": [
+    "# Forward model fields due to the reference model and true model\n",
+    "f0 = problem.fields(sigma0)\n",
+    "f = problem.fields(sigma)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Now ``f`` and ``f0`` are **Field** objects including computed $\\phi$\n",
+    "everywhere. However, this **Field** object know how to compute both\n",
+    "$\\vec{e}$, $\\vec{j}$, and electrical charge, $\\int_V \\rho_v\n",
+    "dV$ ($\\rho_v$ is volumetric charge density). Note that if we know\n",
+    "$\\phi$, all of them can be computed for a corresponding source:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 12,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'f' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-12-434e54a3c544>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mphi\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mf\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0msrc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'phi'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0me\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mf\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0msrc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'e'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0mj\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mf\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0msrc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'j'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0mcharge\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mf\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0msrc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'charge'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+      "\u001b[0;31mNameError\u001b[0m: name 'f' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "phi = f[src, 'phi']\n",
+    "e = f[src, 'e']\n",
+    "j = f[src, 'j']\n",
+    "charge = f[src, 'charge']"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Since the field object for the background model is generic, we can obtain\n",
+    "secondary potential:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "collapsed": true
+   },
+   "outputs": [],
+   "source": [
+    "# Secondary potential\n",
+    "phi0 = f0[src, 'phi']\n",
+    "phi_sec = phi - phi0"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "We present plan and section views of currents, charges, and secondary\n",
+    "potentials in the figure below (Left, middle, and right panels show currents, charges, and secondary potentials)\n",
+    "\n",
+    "![DCfields](../../../images/dc/DCfields.png)\n",
+    "\n",
+    "Current flows from A (+) to B (-) electrode (left to right). Kimberlite pipe\n",
+    "should be more conductive than the background considering more currents are\n",
+    "flowing through the pipe (See distortions of the current path in the left\n",
+    "panel).\n",
+    "\n",
+    "The distribution of electrical charges (the middle panel) supports that the\n",
+    "pipe is conductive since left and right side of the pipe has negative and\n",
+    "positive charges, respectively. In addition, charges only built on the\n",
+    "boundary of the conductive pipe.\n",
+    "\n",
+    "The secondary potential (the right panel) is important since it shows response\n",
+    "from the kimberlite pipe, which often called \"Anomalous potential\". Usually,\n",
+    "removing background response is a good way to see how much anomalous response\n",
+    "could be obtained for the target.\n",
+    "\n",
+    "On the other hand, we cannot measure those fields everywhere but measure\n",
+    "potential differences at MN electrodes (Rx) hence we need to evaluate them\n",
+    "from the fields:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 13,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'sigma' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-13-374816ad2ac7>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m# Get observed data\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mdobs\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0msurvey\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdpred\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msigma\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mf\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mf\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m: name 'sigma' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# Get observed data\n",
+    "dobs = survey.dpred(sigma, f=f)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "If the field has not been computed then we do:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 14,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'sigma' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-14-c9d725a72e5d>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m# Get observed data\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mdobs\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0msurvey\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdpred\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msigma\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m: name 'sigma' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# Get observed data\n",
+    "dobs = survey.dpred(sigma)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "This will compute the field inside of the code then evaluate for data at Rx\n",
+    "locations. Below image shows the computed DC data. Smaller potentials are\n",
+    "obtained at the center locations, which implies the existence of conductive\n",
+    "materials. Current easily flows with conductive materials, which means less\n",
+    "potential is required to path through them, hence for resistive materials we\n",
+    "get greater potential difference measured on the surface. The measured\n",
+    "potential provides some idea of the earth; however, this is not enough, we\n",
+    "want a 3D distribution of the conductivity!\n",
+    "\n",
+    "![DCdata](../../../images/dc/DCdata.png)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Inversion Elements\n",
+    "\n",
+    "Our goal here is finding a 3D conductivity model which explains the observed\n",
+    "data shown above. Inversion elements (red box in the\n",
+    "[SimPEG Framework](#Set-up)) will handle this task with an ability to simulate\n",
+    "forward problem. We go through each element and briefly explain.\n",
+    "\n",
+    "### Mapping\n",
+    "\n",
+    "For the simulation, we used a 3D conductivity model, with a value defined in\n",
+    "every cell center location. However, for the inversion, we may not want to\n",
+    "estimate conductivity at every cell. For instance, our domain include some air\n",
+    "cells, and we already know well about the conductivity of the air\n",
+    "($10^{-8} \\approx 0$) hence, those air cell should be excluded from the\n",
+    "inversion model, $m$. Accordingly, a mapping is required moving from the\n",
+    "inversion model to conductivity model defined at whole discrete domain:\n",
+    "\n",
+    "$$\n",
+    "    \\sigma  = \\mathcal{M}(m)\n",
+    "$$\n",
+    "\n",
+    "In addition, conductivity is strictly positive and varies logarithmically, so\n",
+    "we often use log conductivity as our inversion model ($m = log\n",
+    "(\\sigma)$). Our inversion model is log conductivity only defined below the\n",
+    "subsurface cells, and this can be expressed as\n",
+    "\n",
+    "$$\n",
+    "    \\sigma = \\mathcal{M}_{exp}\\Big(\\mathcal{M}_{act} (m)\\Big),\n",
+    "$$\n",
+    "\n",
+    "where $\\mathcal{M}_{act}(\\cdot)$ is a `InjectActiveCells` map, which\n",
+    "takes subsurface cell and surject to full domain including air cells, and\n",
+    "$\\mathcal{M}_{exp}(\\cdot)$ is an `ExpMap` map takes log conductivity\n",
+    "to conductivity. Combination of two maps are required to get $\\sigma$\n",
+    "from $m$, which can be codified as"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 15,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'Maps' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-15-7d7489ccdc1f>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m# from log conductivity to conductivity\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mexpmap\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mMaps\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mExpMap\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmesh\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      3\u001b[0m \u001b[0;31m# from subsurface cells to full 3D cells\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0mactmap\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mMaps\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mInjectActiveCells\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmesh\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m~\u001b[0m\u001b[0mairind\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mlog\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1e-8\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0mmapping\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mexpmap\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mactmap\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+      "\u001b[0;31mNameError\u001b[0m: name 'Maps' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# from log conductivity to conductivity\n",
+    "expmap = Maps.ExpMap(mesh)\n",
+    "# from subsurface cells to full 3D cells\n",
+    "actmap = Maps.InjectActiveCells(mesh, ~airind, np.log(1e-8))\n",
+    "mapping = expmap*actmap"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Generated mapping should be passed to **Problem** class:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 16,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'mapping' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-16-a60757074025>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m# Generate problem with mapping\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mproblem\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mDC\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mProblem3D_CC\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmesh\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msigmaMap\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mmapping\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m: name 'mapping' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# Generate problem with mapping\n",
+    "problem = DC.Problem3D_CC(mesh, sigmaMap=mapping)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Data Misfit\n",
+    "\n",
+    "Finding a model explaining the observed data requires a measure between\n",
+    "observed ($\\mathbf{d}^{obs}$) and predicted data\n",
+    "($\\mathbf{d}^{dpred}$):\n",
+    "\n",
+    "$$\n",
+    "    \\phi_d = 0.5\\| \\mathbf{W}_d (\\mathbf{d}^{pred}-\\mathbf{d}^{obs})\\|^2_2,\n",
+    "$$\n",
+    "\n",
+    "where $\\mathbf{W}_d = \\mathbf{diag}( \\frac{1}{\\% | \\mathbf{d}^{obs}\n",
+    "|+\\epsilon} )$ is the data weighting matrix. Uncertainty in the observed data\n",
+    "is approximated as $\\% | \\mathbf{d}^{obs} |+\\epsilon$."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 17,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'dobs' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-17-448954ed1fc5>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0msurvey\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstd\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mstd\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0msurvey\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0meps\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0meps\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 5\u001b[0;31m \u001b[0msurvey\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdobs\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mdobs\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      6\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      7\u001b[0m \u001b[0;31m# Define datamisfit portion of objective function\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+      "\u001b[0;31mNameError\u001b[0m: name 'dobs' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# percentage and floor for uncertainty in the observed data\n",
+    "std, eps = 0.05, 1e-3\n",
+    "survey.std = std\n",
+    "survey.eps = eps\n",
+    "survey.dobs = dobs\n",
+    "\n",
+    "# Define datamisfit portion of objective function\n",
+    "dmisfit = DataMisfit.l2_DataMisfit(survey)\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Regularization\n",
+    "\n",
+    "The objective function includes both data misfit and regularization terms,\n",
+    "$\\phi_m$ :\n",
+    "\n",
+    "$$\n",
+    "    \\phi = \\phi_d + \\beta \\phi_m\n",
+    "$$\n",
+    "\n",
+    "We use Tikhonov-style regularization including both smoothness and smallness\n",
+    "terms. For further details of this See XXX.\n",
+    "\n",
+    "In addition, considering the geometry of the gradient array: a single source\n",
+    "and distributed receivers, this specific DC survey may not have much depth\n",
+    "resolution similar to magnetic and gravity data. Depth weighting\n",
+    "($\\frac{1}{(z-z_0)^3}$) is often used to handle this. And with this\n",
+    "weight we form **Regularization** class:"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 18,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'Regularization' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-18-0093d24eefda>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0;31m# Define regulatization (model objective function)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 6\u001b[0;31m \u001b[0mreg\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mRegularization\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mSimple\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmesh\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmapping\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mregmap\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mindActive\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m~\u001b[0m\u001b[0mairind\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      7\u001b[0m \u001b[0mreg\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcell_weights\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mdepth\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m~\u001b[0m\u001b[0mairind\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      8\u001b[0m \u001b[0mreg\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0malpha_s\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;36m1e-1\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+      "\u001b[0;31mNameError\u001b[0m: name 'Regularization' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# Depth weighting\n",
+    "depth = 1./(abs(mesh.gridCC[:,2]-zc))**1.5\n",
+    "depth = depth/depth.max()\n",
+    "\n",
+    "# Define regulatization (model objective function)\n",
+    "reg = Regularization.Simple(mesh, mapping=regmap, indActive=~airind)\n",
+    "reg.cell_weights = depth[~airind]\n",
+    "reg.alpha_s = 1e-1\n",
+    "reg.alpha_x = 1.\n",
+    "reg.alpha_y = 1.\n",
+    "reg.alpha_z = 1."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Optimization\n",
+    "\n",
+    "To minimize the objective function, an optimization scheme is required. The\n",
+    "**Optimization** class handles this, and we use Inexact Gauss Newton Scheme\n",
+    "CITExxx."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 19,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'Optimization' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-19-e4631c7d7d37>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mopt\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mOptimization\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mInexactGaussNewton\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmaxIter\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;36m20\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m: name 'Optimization' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "opt = Optimization.InexactGaussNewton(maxIter = 20)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### InvProblem\n",
+    "\n",
+    "Both **DatamMisfit** and **Regularization** classes are created, and an\n",
+    "**Optimization** is chosen. To pose the inverse problem, they need to be\n",
+    "declared as an optimization problem:\n",
+    "\n",
+    "$$\n",
+    "    \\text{minimize} \\ \\phi_d + \\beta \\phi_m \\ \\\\\n",
+    "    s.t. \\ \\text{some constratins}\n",
+    "$$\n",
+    "\n",
+    "The **InvProblem** class can be set with **DatamMisfit**,  **Regularization** and **Optimiztion**\n",
+    "classes."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "collapsed": true
+   },
+   "outputs": [],
+   "source": [
+    "invProb = InvProblem.BaseInvProblem(dmisfit, reg, opt)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Inversion\n",
+    "\n",
+    "We have stated our inverse problem, but a conductor is required, who directs\n",
+    "our inverse problem. **Directives** conducts our **Inversion**. For instance,\n",
+    "the trade-off parameter, $\\beta$ needs to be estimated, and sometimes\n",
+    "cooled in the inversion iterations. A target misfit is need to be set usually\n",
+    "upon discrepancy principle ($\\phi_d^\\ast = 0.5 N_d$, where $N_d$\n",
+    "is the number of data).\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 21,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'invProb' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-21-e4ca6bc5751c>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      6\u001b[0m \u001b[0;31m# Define Inversion class\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 7\u001b[0;31m \u001b[0minv\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mInversion\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mBaseInversion\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minvProb\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdirectiveList\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mbeta\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mbetaest\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtarget\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m: name 'invProb' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# Define Directives\n",
+    "betaest = Directives.BetaEstimate_ByEig(beta0_ratio=1e0)\n",
+    "beta = Directives.BetaSchedule(coolingFactor=5, coolingRate=2)\n",
+    "target = Directives.TargetMisfit()\n",
+    "\n",
+    "# Define Inversion class\n",
+    "inv = Inversion.BaseInversion(invProb, directiveList=[beta, betaest, target])"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Run\n",
+    "\n",
+    "Now we all set. The initial model is assumed to be homogeneous."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 23,
+   "metadata": {
+    "collapsed": false
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'airind' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "\u001b[0;32m<ipython-input-23-2b3c69dc59ea>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m# Create inital and reference model (1e-4 S/m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mm0\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mones\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmesh\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnC\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m~\u001b[0m\u001b[0mairind\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mlog\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1e-4\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      3\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0;31m# Run inversion\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0mmopt\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0minv\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mrun\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mm0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
+      "\u001b[0;31mNameError\u001b[0m: name 'airind' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# Create inital and reference model (1e-4 S/m)\n",
+    "m0 = np.ones(mesh.nC)[~airind]*np.log(1e-4)\n",
+    "\n",
+    "# Run inversion\n",
+    "mopt = inv.run(m0)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "\n",
+    "The inversion runs and reaches to the target misfit, hence we fit the observed\n",
+    "data.\n",
+    "\n",
+    "![DCObsPred](../../../images/dc/DCObsPred.png)\n",
+    "\n",
+    "A 3D conductivity model is recovered and compared with the true conductivity\n",
+    "model. A conductive pipe at depth is recovered!\n",
+    "\n",
+    "![DCObsPred](../../../images/dc/Cond3D.png)\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "collapsed": true
+   },
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 2",
+   "language": "python",
+   "name": "python2"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 2
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython2",
+   "version": "2.7.13"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}