Head in leaky layer

docs/transient/02examples/head_in_vertical.ipynb

@@ -0,0 +1,193 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Head propagation in leaky layers"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import matplotlib.pyplot as plt\n",
+    "import numpy as np\n",
+    "\n",
+    "import timflow.transient as tft\n",
+    "\n",
+    "plt.rcParams[\"font.size\"] = 8.0"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Consider a three-aquifer system. Aquifer properties are given in Table 1. All aquifers have elastic storage. A well is located at $(x,y)=(0,0)$ and is screened in layer 1. The well starts pumping at time $t=0$ with a discharge $Q=1000$ m$^3$/d. The radius of the well is 0.2 m.\n",
+    "\n",
+    "**Table 1 - Aquifer properties for exercise 1**\n",
+    "|Layer          | $k$ (m/d) | $c$ (d) | $S_s$ (m$^{-1}$) | $z_t$ (m) | $z_b$ (m)|\n",
+    "|---------------| ---------:| -------:| -----:| ---------:| --------:|\n",
+    "|Aquifer 0      |      1    |    -    |0.0001 |   25      |        20|\n",
+    "|Leaky layer 1  |      -    |  1000   |0      |   20      |        18|\n",
+    "|Aquifer 1      |     20    |    -    |0.0001 |   18      |        10|\n",
+    "|Leaky layer 2  |      -    |  2000   |0      |   10      |         8|\n",
+    "|Aquifer 2      |      2    |    -    |0.0001 |    8      |         0|\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Confined top"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "ml = tft.ModelMaq(\n",
+    "    kaq=[5, 10, 20],\n",
+    "    z=[25, 22, 18, 13, 8, 0],\n",
+    "    c=[1000, 2000],\n",
+    "    Saq=[1e-4, 1e-4, 1e-4],\n",
+    "    Sll=[1e-3, 1e-3],\n",
+    "    topboundary=\"conf\",\n",
+    "    tmin=0.1,\n",
+    "    tmax=1000,\n",
+    ")\n",
+    "w1 = tft.Well(ml, xw=0, yw=0, rw=0.2, tsandQ=[(0, 1000)], layers=1)\n",
+    "w2 = tft.Well(ml, xw=20, yw=0, rw=0.2, tsandQ=[(0, 177)], layers=1)\n",
+    "ml.solve()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "z = np.linspace(0, 25, 100)\n",
+    "h = ml.headinvertical(10, 10, z, [0.1, 1])"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "z = np.linspace(0, 25, 100)\n",
+    "t = [0.1, 1, 10, 100]\n",
+    "h = ml.headinvertical(10, 10, z, t)\n",
+    "for i in range(len(t)):\n",
+    "    plt.plot(h[:, i], z)\n",
+    "for i in range(ml.aq.naq - 1):\n",
+    "    plt.axhspan(ml.aq.zaqtop[i + 1], ml.aq.zaqbot[i], color=[0.9, 0.9, 0.9])\n",
+    "plt.axhspan(25, 30, color=[0.6, 0.6, 0.6])\n",
+    "plt.ylim(z[0], 30);"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Semi-confined top"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "ml = tft.ModelMaq(\n",
+    "    kaq=[5, 10, 20],\n",
+    "    z=[30, 25, 22, 18, 13, 8, 0],\n",
+    "    c=[500, 1000, 2000],\n",
+    "    Saq=[1e-4, 1e-4, 1e-4],\n",
+    "    Sll=[1e-3, 1e-3, 1e-3],\n",
+    "    topboundary=\"semi\",\n",
+    "    tmin=0.1,\n",
+    "    tmax=1000,\n",
+    ")\n",
+    "w1 = tft.Well(ml, xw=0, yw=0, rw=0.2, tsandQ=[(0, 1000)], layers=1)\n",
+    "w2 = tft.Well(ml, xw=20, yw=0, rw=0.2, tsandQ=[(0, 177)], layers=1)\n",
+    "ml.solve()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "z = np.linspace(0, 30, 100)\n",
+    "t = [0.1, 1, 10, 100]\n",
+    "h = ml.headinvertical(10, 10, z, t)\n",
+    "for i in range(len(t)):\n",
+    "    plt.plot(h[:, i], z)\n",
+    "for i in range(len(ml.aq.z[::2]) - 1):\n",
+    "    plt.axhspan(ml.aq.z[2 * i + 1], ml.aq.z[2 * i], color=[0.9, 0.9, 0.9])\n",
+    "plt.ylim(z[0], z[-1]);"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Leaky layer top"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "ml = tft.ModelMaq(\n",
+    "    kaq=[5, 10, 20],\n",
+    "    z=[30, 25, 22, 18, 13, 8, 0],\n",
+    "    c=[500, 1000, 2000],\n",
+    "    Saq=[1e-4, 1e-4, 1e-4],\n",
+    "    Sll=[1e-3, 1e-3, 1e-3],\n",
+    "    topboundary=\"lea\",\n",
+    "    tmin=0.1,\n",
+    "    tmax=1000,\n",
+    ")\n",
+    "w1 = tft.Well(ml, xw=0, yw=0, rw=0.2, tsandQ=[(0, 1000)], layers=1)\n",
+    "w2 = tft.Well(ml, xw=20, yw=0, rw=0.2, tsandQ=[(0, 177)], layers=1)\n",
+    "ml.solve()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "z = np.linspace(0, 30, 100)\n",
+    "t = [0.1, 1, 10, 100]\n",
+    "h = ml.headinvertical(10, 10, z, t)\n",
+    "for i in range(len(t)):\n",
+    "    plt.plot(h[:, i], z)\n",
+    "for i in range(len(ml.aq.z[::2]) - 1):\n",
+    "    plt.axhspan(ml.aq.z[2 * i + 1], ml.aq.z[2 * i], color=[0.9, 0.9, 0.9])\n",
+    "plt.ylim(z[0], z[-1]);"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}

docs/transient/02examples/index.rst

@@ -73,4 +73,16 @@ Pathlines
     :caption: Pathlines
     :hidden:
 
-    pathline_trace
+    pathline_trace
+
+Head in leaky layers
+--------------------
+
+- :doc:`head_in_vertical`
+
+.. toctree::
+    :maxdepth: 2
+    :caption: Head in leaky layers
+    :hidden:
+
+    head_in_vertical