diff --git a/docs/transient/04pumpingtests/confined1_oude_korendijk.ipynb b/docs/transient/04pumpingtests/confined1_oude_korendijk.ipynb index b5ce880..acb713c 100644 --- a/docs/transient/04pumpingtests/confined1_oude_korendijk.ipynb +++ b/docs/transient/04pumpingtests/confined1_oude_korendijk.ipynb @@ -2,7 +2,13 @@ "cells": [ { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "# 1. Confined Aquifer Test - Oude Korendijk" ] @@ -49,7 +55,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "\n" + "" ] }, { @@ -129,7 +135,7 @@ "cal.set_parameter(name=\"Saq\", initial=1e-4, layers=0)\n", "cal.series(name=\"obs1\", x=ro1, y=0, t=to1, h=ho1, layer=0) # Adding well 1\n", "cal.series(name=\"obs2\", x=ro2, y=0, t=to2, h=ho2, layer=0) # Adding well 2\n", - "cal.fit(report=True)" + "cal.fit(report=False)" ] }, { @@ -138,14 +144,20 @@ "metadata": {}, "outputs": [], "source": [ - "display(cal.parameters)\n", - "print(\"RMSE:\", cal.rmse())" + "display(cal.parameters.loc[:, [\"optimal\"]])\n", + "print(f\"RMSE: {cal.rmse():.3f} m\")" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "hm1 = ml.head(ro1, 0, to1)\n", @@ -155,8 +167,8 @@ "plt.semilogx(to2, ho2, \"C1.\", label=\"obs at 90m\")\n", "plt.semilogx(to2, hm2[0], \"C1\", label=\"timflow at 90m\")\n", "plt.title(\"Model Results\")\n", - "plt.xlabel(\"time(d)\")\n", - "plt.ylabel(\"drawdown(m)\")\n", + "plt.xlabel(\"time (d)\")\n", + "plt.ylabel(\"head change (m)\")\n", "plt.legend()\n", "plt.grid()" ] @@ -172,13 +184,21 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "The performance of `timflow` was evaluated by comparison with AQTESOLV (Duffield, 2007), MLU (Carlson and Randall, 2012), and Kruseman and de Ridder (1970), here abbreviated to K&dR. Results from `timflow` and AQTESOLV are identical, those from MLU are very similar, while the solution of Kruseman and de Ridder shows small deviations." + "The performance of `timflow` was compared with AQTESOLV (Duffield, 2007), MLU (Carlson and Randall, 2012), and Kruseman and de Ridder (1970), here abbreviated to K&dR. Results from `timflow` and AQTESOLV are identical, those from MLU are very similar, while the solution of Kruseman and de Ridder shows small differences." ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, "outputs": [], "source": [ "t = pd.DataFrame(\n", @@ -195,7 +215,7 @@ " {\n", " \"k [m/d]\": \"{:.2f}\",\n", " \"Ss [1/m]\": \"{:.2e}\",\n", - " \"RMSE [m]\": lambda x: \"-\" if x == \"-\" else f\"{float(x):.2f}\",\n", + " \"RMSE [m]\": lambda x: \"-\" if x == \"-\" else f\"{float(x):.3f}\",\n", " }\n", ")\n", "t_formatted" diff --git a/docs/transient/04pumpingtests/confined2_grindley.ipynb b/docs/transient/04pumpingtests/confined2_grindley.ipynb index cd3cf73..db6c4ec 100644 --- a/docs/transient/04pumpingtests/confined2_grindley.ipynb +++ b/docs/transient/04pumpingtests/confined2_grindley.ipynb @@ -43,7 +43,13 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Introduction and Conceptual Model" ] @@ -58,14 +64,22 @@ "tags": [] }, "source": [ - "This example is a pumping test conducted in 1953 in Gridley, Illinois, US. It was reported by Walton (1962). The aquifer is an 18 ft thick sand and gravel layer under confined conditions. The pumping well fully penetrates the formation, and pumping was conducted for 8 hours at a rate of 220 gallons per minute. The effect of pumping was observed at observation well 1, located 824 ft away from the well.\n", + "This example is a pumping test conducted in 1953 in Gridley, Illinois, US. It was reported by Walton (1962). \n", + "\n", + "The aquifer is an 18 ft thick sand and gravel layer under confined conditions. The pumping well fully penetrates the formation, and pumping was conducted for 8 hours at a rate of 220 gallons per minute. The effect of pumping was observed at observation well 1, located 824 ft away from the well.\n", "\n", "The time-drawdown data for the observation well were obtained from the AQTESOLV documentation (Duffield, 2007), while data from the pumping well were obtained from the original paper from Walton (1962). Following AQTESOLV documentation (Duffield, 2007), radii of 0.5 ft were used for both the pumping and observation wells. " ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "" ] @@ -201,7 +215,13 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "# unknown parameters: kaq, Saq\n", @@ -212,23 +232,35 @@ "cal.set_parameter(name=\"Saq\", initial=1e-4, pmin=1e-7, layers=0)\n", "cal.series(name=\"obs1\", x=r, y=0, t=to1, h=ho1, layer=0) # setting the observation data\n", "cal.seriesinwell(name=\"obs2\", element=w, t=to2, h=ho2) # setting the observation data\n", - "cal.fit(report=True)" + "cal.fit(report=False)" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ - "display(cal.parameters)\n", - "print(\"RMSE:\", cal.rmse())" + "display(cal.parameters.loc[:, [\"optimal\"]])\n", + "print(f\"RMSE: {cal.rmse():.3f} m\")" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "hm1 = ml.head(r, 0, to1)\n", @@ -239,7 +271,7 @@ "plt.semilogx(to2, ho2, \"C1.\", label=\"obs well\")\n", "plt.title(\"Model Results\")\n", "plt.xlabel(\"time [d]\")\n", - "plt.ylabel(\"head [m]\")\n", + "plt.ylabel(\"head change [m]\")\n", "plt.legend()\n", "plt.grid()" ] @@ -267,7 +299,7 @@ "tags": [] }, "source": [ - "The performance of `timflow` with two datasets simultaneously was evaluated by comparison with AQTESOLV (Duffield, 2007), and MLU (Carlson and Randall, 2012). Results from `timflow` with added wellbore storage and MLU are identical, while those from AQTESOLV show small deviations." + "The performance of `timflow` with two datasets simultaneously was compared with AQTESOLV (Duffield, 2007) and MLU (Carlson and Randall, 2012). Results from `timflow` with added wellbore storage and MLU are identical, while those from AQTESOLV show small differences." ] }, { @@ -293,7 +325,7 @@ "t.loc[\"AQTESOLV\"] = [37.803, 1.356e-06, 0.270]\n", "\n", "t_formatted = t.style.format(\n", - " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"RMSE [m]\": \"{:.2f}\"}\n", + " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"RMSE [m]\": \"{:.3f}\"}\n", ")\n", "t_formatted" ] diff --git a/docs/transient/04pumpingtests/confined3_sioux.ipynb b/docs/transient/04pumpingtests/confined3_sioux.ipynb new file mode 100644 index 0000000..153e37f --- /dev/null +++ b/docs/transient/04pumpingtests/confined3_sioux.ipynb @@ -0,0 +1,357 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "# 3. Confined Aquifer Test - Sioux Flats" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Import packages" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "import pandas as pd\n", + "\n", + "import timflow.transient as tft\n", + "\n", + "plt.rcParams[\"figure.figsize\"] = (5, 3) # default figure size" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Introduction and Conceptual Model" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "This example is a pumping test done in Sioux Flats, South Dakota, USA. The data comes from the AQTESOLV documentation (Duffield, 2007). \n", + "\n", + "The aquifer is 50 ft thick and is bounded by impermeable layers. The test was conducted for 2045 minutes (~34 hours), with a constant pumping rate of 2.7 ft$^3$/s. Drawdown data has been collected at three piezometers located 100, 200 and 400 ft away, respectively. The well radius is 0.5 ft. " + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Load data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# time and drawdown of piezometer 100ft away from pumping well\n", + "data1 = np.loadtxt(\"data/sioux100.txt\")\n", + "t1 = data1[:, 0]\n", + "h1 = data1[:, 1]\n", + "\n", + "# time and drawdown of piezometer 200ft away from pumping well\n", + "data2 = np.loadtxt(\"data/sioux200.txt\")\n", + "t2 = data2[:, 0]\n", + "h2 = data2[:, 1]\n", + "\n", + "# time and drawdown of piezometer 300ft away from pumping well\n", + "data3 = np.loadtxt(\"data/sioux400.txt\")\n", + "t3 = data3[:, 0]\n", + "h3 = data3[:, 1]" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Parameters and model" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# known parameters\n", + "Q = (\n", + " (2.7 * 0.3048**3) * 60 * 60 * 24\n", + ") # constant discharge in m^3/d (2.7 ft^3/s = 6605.754 m^3/d)\n", + "b = 50 * 0.3048 # aquifer thickness in m (50 ft = 15.24 m)\n", + "rw = 0.5 * 0.3048 # well radius in m (0.5 ft = 0.1524 m)\n", + "r1 = 100 * 0.3048 # distance between obs1 to pumping well (100 ft = 30.48 m)\n", + "r2 = 200 * 0.3048 # distance between obs2 to pumping well (200 ft = 60.96 m)\n", + "r3 = 400 * 0.3048 # distance between obs3 to pumping well (400 ft = 121.92 m)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# timflow model\n", + "ml = tft.ModelMaq(kaq=10, z=[0, -b], Saq=0.001, tmin=0.001, tmax=10, topboundary=\"conf\")\n", + "w = tft.Well(ml, xw=0, yw=0, rw=rw, tsandQ=[(0, Q)], layers=0)\n", + "ml.solve()" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Estimate aquifer parameters" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# unknown parameters: k, Saq\n", + "cal = tft.Calibrate(ml)\n", + "cal.set_parameter(name=\"kaq\", initial=10, layers=0)\n", + "cal.set_parameter(name=\"Saq\", initial=1e-4, layers=0)\n", + "cal.series(name=\"obs1\", x=r1, y=0, t=t1, h=h1, layer=0) # Adding well 1\n", + "cal.series(name=\"obs2\", x=r2, y=0, t=t2, h=h2, layer=0) # Adding well 2\n", + "cal.series(name=\"obs3\", x=r3, y=0, t=t3, h=h3, layer=0) # Adding well 3\n", + "cal.fit(report=False)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "display(cal.parameters.loc[:, [\"optimal\"]])\n", + "print(f\"RMSE: {cal.rmse():.3f} m\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "hm1 = ml.head(r1, 0, t1)\n", + "hm2 = ml.head(r2, 0, t2)\n", + "hm3 = ml.head(r3, 0, t3)\n", + "plt.semilogx(t1, h1, \".\", label=\"obs1\")\n", + "plt.semilogx(t1, hm1[0], label=\"timflow result 1\")\n", + "plt.semilogx(t2, h2, \".\", label=\"obs2\")\n", + "plt.semilogx(t2, hm2[0], label=\"timflow result 2\")\n", + "plt.semilogx(t3, h3, \".\", label=\"obs3\")\n", + "plt.semilogx(t3, hm3[0], label=\"timflow result 3\")\n", + "plt.title(\"Model Results\")\n", + "plt.xlabel(\"time [d]\")\n", + "plt.ylabel(\"head change [m]\")\n", + "plt.legend()\n", + "plt.grid()" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Comparison of results\n", + "The performance of `timflow` was compared with AQTESOLV (Duffield, 2007) and MLU (Carlson and Randall, 2012).`timflow` achieved similar results as the other software. " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "t = pd.DataFrame(\n", + " columns=[\"k [m/d]\", \"Ss [1/m]\", \"RMSE [m]\"], index=[\"timflow\", \"AQTESOLV\", \"MLU\"]\n", + ")\n", + "\n", + "t.loc[\"timflow\"] = np.append(cal.parameters[\"optimal\"].values, cal.rmse())\n", + "t.loc[\"AQTESOLV\"] = [282.659, 4.211e-03, 0.003925]\n", + "t.loc[\"MLU\"] = [282.684, 4.209e-03, 0.003897]\n", + "\n", + "t_formatted = t.style.format(\n", + " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"RMSE [m]\": \"{:.3f}\"}\n", + ")\n", + "t_formatted" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "## References" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "* Carlson, F. and Randall, J. (2012), MLU: a Windows application for the analysis of aquifer tests and the design of well fields in layered systems, Ground Water 50(4):504–510\n", + "* Duffield, G.M. (2007), AQTESOLV for Windows Version 4.5 User's Guide, HydroSOLVE, Inc., Reston, VA.\n", + "* Newville, M., Stensitzki, T., Allen, D.B., Ingargiola, A. (2014), LMFIT: Non Linear Least-Squares Minimization and Curve Fitting for Python, https://dx.doi.org/10.5281/zenodo.11813, https://lmfit.github.io/lmfit-py/intro.html (last access: August,2021)." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "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.13.5" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/docs/transient/04pumpingtests/figs/Dalem.png b/docs/transient/04pumpingtests/figs/Dalem.png index 05ba9e5..d28b157 100644 Binary files a/docs/transient/04pumpingtests/figs/Dalem.png and b/docs/transient/04pumpingtests/figs/Dalem.png differ diff --git a/docs/transient/04pumpingtests/figs/Dalem_model1.jpeg b/docs/transient/04pumpingtests/figs/Dalem_model1.jpeg deleted file mode 100644 index 7b4c892..0000000 Binary files a/docs/transient/04pumpingtests/figs/Dalem_model1.jpeg and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/Dalem_model1.png b/docs/transient/04pumpingtests/figs/Dalem_model1.png deleted file mode 100644 index d28b157..0000000 Binary files a/docs/transient/04pumpingtests/figs/Dalem_model1.png and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/Dalem_model2.png b/docs/transient/04pumpingtests/figs/Dalem_model2.png deleted file mode 100644 index 2b7ac1e..0000000 Binary files a/docs/transient/04pumpingtests/figs/Dalem_model2.png and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/Hardinxveld.jpeg b/docs/transient/04pumpingtests/figs/Hardinxveld.jpeg deleted file mode 100644 index a5e38cf..0000000 Binary files a/docs/transient/04pumpingtests/figs/Hardinxveld.jpeg and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/Multi_well.png b/docs/transient/04pumpingtests/figs/Multi_well.png new file mode 100644 index 0000000..2668c51 Binary files /dev/null and b/docs/transient/04pumpingtests/figs/Multi_well.png differ diff --git a/docs/transient/04pumpingtests/figs/Nevada.png b/docs/transient/04pumpingtests/figs/Nevada.png deleted file mode 100644 index fc333b5..0000000 Binary files a/docs/transient/04pumpingtests/figs/Nevada.png and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/NevadaTTim.png b/docs/transient/04pumpingtests/figs/NevadaTTim.png deleted file mode 100644 index cc4353e..0000000 Binary files a/docs/transient/04pumpingtests/figs/NevadaTTim.png and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/Schroth.jpeg b/docs/transient/04pumpingtests/figs/Schroth.jpeg deleted file mode 100644 index 3f8c6e6..0000000 Binary files a/docs/transient/04pumpingtests/figs/Schroth.jpeg and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/Schroth.png b/docs/transient/04pumpingtests/figs/Schroth.png deleted file mode 100644 index 075ecac..0000000 Binary files a/docs/transient/04pumpingtests/figs/Schroth.png and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/Schroth_New.png b/docs/transient/04pumpingtests/figs/Schroth_New.png deleted file mode 100644 index 8823895..0000000 Binary files a/docs/transient/04pumpingtests/figs/Schroth_New.png and /dev/null differ diff --git a/docs/transient/04pumpingtests/figs/TexasHill.png b/docs/transient/04pumpingtests/figs/TexasHill.png deleted file mode 100644 index f88f5ae..0000000 Binary files a/docs/transient/04pumpingtests/figs/TexasHill.png and /dev/null differ diff --git a/docs/transient/04pumpingtests/index.rst b/docs/transient/04pumpingtests/index.rst index 2fcb1ec..1726a26 100644 --- a/docs/transient/04pumpingtests/index.rst +++ b/docs/transient/04pumpingtests/index.rst @@ -12,6 +12,7 @@ problems such as pumping tests and slug tests. More to come soon. confined1_oude_korendijk confined2_grindley + confined3_sioux .. toctree:: :maxdepth: 1 @@ -20,6 +21,7 @@ problems such as pumping tests and slug tests. More to come soon. leaky1_dalem leaky2_hardixveld + leaky3_texashill .. toctree:: :maxdepth: 1 @@ -35,15 +37,17 @@ Confined Pumping Tests 1. :doc:`Oude Korendijk ` - One layer confined pumping test with two observation wells 2. :doc:`Grindley ` - One layer confined pumping test with data from both observation well and pumping well. +3. :doc:`Sioux ` - One layer confined pumping test with three observation wells. Leaky Pumping Tests (Semi-confined) ----------------------------------- 1. :doc:`Dalem ` - One layer, semi-confined aquifer test with four observation wells. 2. :doc:`Hardixveld ` - Four layers (Aquitard - Aquifer - Aquitard - Aquifer) semi-confined pumping test with data from the pumping well. +3. :doc:`Texas Hill ` - Two layers, semi-confined aquifer test with three observation wells. Slug Tests ---------- -1. :doc:`Pratt County ` - One layer partially penetrated slug test with data from the test well. -2. :doc:`Falling Head ` - One layer partially penetrated slug test with data from the test well. +1. :doc:`Pratt County ` - Slug test with data from the partially penetrating test well. +2. :doc:`Multi well ` - Slug test with data from both fully penetrating test well and observation well. diff --git a/docs/transient/04pumpingtests/leaky1_dalem.ipynb b/docs/transient/04pumpingtests/leaky1_dalem.ipynb index c97544c..30e1aa6 100644 --- a/docs/transient/04pumpingtests/leaky1_dalem.ipynb +++ b/docs/transient/04pumpingtests/leaky1_dalem.ipynb @@ -2,14 +2,26 @@ "cells": [ { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "# 1. Leaky Aquifer Test - Dalem" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Import packages" ] @@ -37,15 +49,16 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Introduction and Conceptual Model\n", - "This example is a pumping test from Dalem (Kruseman et al., 1970), the Netherlands. The hydrogeological cross-section is composed of the following elements: \n", - "* an initial 8 m deep aquitard layer,\n", - "* followed by an aquifer from 8 m to 45 m depth,\n", - "* the layer underlying the aquifer is considered an aquiclude.\n", - "\n", - "The pumping well is placed at the aquifer, and drawdown is recorded at four different piezometers, 30, 60, 90 and 120 m away from the well. The pumping lasted 8 hours in total at a rate of 761 m$^3$/d. There is a river 1500 m away from the well. The tide affects both river and well levels. Data has been previously corrected for the tide effect." + "This example is a pumping test from Dalem, the Netherlands (Kruseman and de Ridder, 1970). An aquifer of 37 m thickness is overlain by a semi-confining layer of 8 m thick. Water is pumped from the aquifer and drawdown is recorded at four different piezometers, 30, 60, 90 and 120 m away from the well. The pumping lasted 8 hours in total at a rate of 761 m$^3$/d. There is a river 1500 m away from the well. The tide affects both river and well levels. Data was corrected for the tide effect." ] }, { @@ -60,12 +73,18 @@ ] }, "source": [ - "\n" + "" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Load Data" ] @@ -73,13 +92,19 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "# data of observation well 30 m away from pumping well\n", "data1 = np.loadtxt(\"data/dalem_p30.txt\", skiprows=1)\n", - "t1 = data1[:, 0] # days\n", - "h1 = data1[:, 1] # m\n", + "t1 = data1[:, 0]\n", + "h1 = data1[:, 1]\n", "\n", "# data of observation well 60 m away from pumping well\n", "data2 = np.loadtxt(\"data/dalem_p60.txt\", skiprows=1)\n", @@ -99,7 +124,13 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Parameters and model" ] @@ -107,9 +138,16 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ + "# known parameters\n", "H = 37 # aquifer thickness in m\n", "zt = -8 # top boundary of aquifer in m\n", "zb = zt - H # bottom boundary of the aquifer in m\n", @@ -124,36 +162,61 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "In many situations, we cannot ignore the leakage potential of overlying and underlying formations to aquifers, and we cannot conceptualize them as confined. `timflow` is capable of modelling and adjusting parameters to leaky, semi-confined aquifers." + "The `timflow` model consists of a single semi-confined aquifer and a pumping well." ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ - "# unknown parameters: kaq, Saq, c\n", + "# unkonwn parameters: kaq, Saq, c\n", "ml = tft.ModelMaq(\n", - " kaq=10, z=[0, zt, zb], c=500, Saq=0.001, topboundary=\"semi\", tmin=0.01, tmax=1\n", + " kaq=10, z=[0, zt, zb], c=500, Saq=1e-4, topboundary=\"semi\", tmin=0.01, tmax=1\n", ")\n", - "w = tft.Well(ml, xw=0, yw=0, tsandQ=[(0, Q), (0.34, 0)])\n", + "w = tft.Well(ml, xw=0, yw=0, tsandQ=[(0, Q)])\n", "ml.solve(silent=\"True\")" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "### Estimate aquifer parameters" + "### Estimate aquifer parameters\n", + "The hydraulic condivity and specific storage coefficient of the aquifer are estimated and the resistance of the leaky semi-confining layer. " ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "cal = tft.Calibrate(ml)\n", @@ -165,23 +228,35 @@ "cal.series(name=\"obs2\", x=60, y=0, t=t2, h=h2, layer=0)\n", "cal.series(name=\"obs3\", x=90, y=0, t=t3, h=h3, layer=0)\n", "cal.series(name=\"obs4\", x=120, y=0, t=t4, h=h4, layer=0)\n", - "cal.fit(report=True)" + "cal.fit(report=False)" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ - "display(cal.parameters)\n", - "print(\"RMSE:\", cal.rmse())" + "display(cal.parameters.loc[:, [\"optimal\"]])\n", + "print(f\"RMSE: {cal.rmse():.5f} m\")" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "hm_1 = ml.head(r1, 0, t1)\n", @@ -200,13 +275,55 @@ "plt.semilogx(t4, h4, \".\", label=\"obs at 120 m\")\n", "plt.semilogx(t4, hm_4[0], label=\"timflow at 120 m\")\n", "\n", - "plt.xlabel(\"time(d)\")\n", - "plt.ylabel(\"drawdown(m)\")\n", + "plt.xlabel(\"time (d)\")\n", + "plt.ylabel(\"head change (m)\")\n", "plt.title(\"Model Results\")\n", "plt.legend(bbox_to_anchor=(1.05, 1))\n", "plt.grid();" ] }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Estimate aquifer parameters using weights\n", + "The head change is smaller in the observation wells that are farther away from the pumping well. It may possibly be expected that the difference between the modeled and measured values is smaller when the total drawdown is smaller. This may be taken into account by giving the measurements weights. In the calibration below, the weight for each observation well is equal to one divided by the maximum absolute head change in that overservation well. It is up to the modeler to decide whether such a weighing is desirable. Note that the root mean squared error (computed without weights) is a bit larger than for the case without the weights, but the increase is small. The fitted values of $k$ and $S_s$ are similar to the case without weights, but the resistance $c$ is quite a bit larger. " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "cal2 = tft.Calibrate(ml)\n", + "cal2.set_parameter(name=\"kaq\", initial=10, layers=0)\n", + "cal2.set_parameter(name=\"Saq\", initial=1e-4, layers=0)\n", + "cal2.set_parameter(name=\"c\", initial=500, pmin=0, layers=0)\n", + "\n", + "cal2.series(name=\"obs1\", x=30, y=0, t=t1, h=h1, layer=0, weights=1 / np.max(np.abs(h1)))\n", + "cal2.series(name=\"obs2\", x=60, y=0, t=t2, h=h2, layer=0, weights=1 / np.max(np.abs(h2)))\n", + "cal2.series(name=\"obs3\", x=90, y=0, t=t3, h=h3, layer=0, weights=1 / np.max(np.abs(h3)))\n", + "cal2.series(name=\"obs4\", x=120, y=0, t=t4, h=h4, layer=0, weights=1 / np.max(np.abs(h4)))\n", + "cal2.fit(report=True)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "display(cal2.parameters[[\"layers\", \"optimal\"]])\n", + "print(f\"RMSE: {cal2.rmse(weighted=False):.5f} m\")" + ] + }, { "cell_type": "markdown", "metadata": {}, @@ -218,29 +335,38 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "The performance of `timflow` was evaluated by comparison with other numerical models and with parameter values reported in Kruseman and de Ridder (1970). The published values were obtained by graphical matching to the Hantush family of type curves (Hantush, 1955). In addition to the `timflow` results and the literature values, results from AQTESOLV (Duffield, 2007) and MLU (Carlson & Randall, 2012), as reported by Yang (2020), are included for comparison.\n", + "The performance of `timflow` is compared with other models and with parameter values reported in Kruseman and de Ridder (1970). The published values were obtained by graphical matching to a Hantush type curve (Hantush, 1955). In addition to the `timflow` results and the literature values, results from AQTESOLV (Duffield, 2007) and MLU (Carlson & Randall, 2012) are included for comparison. Both AQTESOLV and MLU were applied to minimize the sum of the squared log-transformed drawdowns. \n", "\n", - "Overall, all models yield similar estimates of the aquifer hydraulic conductivity. However, the `timflow` results differ substantially from the MLU and AQTESOLV solutions with respect to aquitard storage and hydraulic resistance." + "Overall, all models yield similar estimates of the hydraulic conductivity and specific storage coefficient. The fitted value of the resistance differs. However, the `timflow` results differ substantially from the MLU and AQTESOLV solutions with respect to aquitard storage and hydraulic resistance." ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "t = pd.DataFrame(\n", " columns=[\"k [m/d]\", \"Ss [1/m]\", \"c [d]\", \"RMSE [m]\"],\n", - " index=[\"timflow\", \"AQTESOLV\", \"MLU\", \"Hantush\"],\n", + " index=[\"timflow\", \"timflow weights\", \"AQTESOLV\", \"MLU\", \"Hantush\"],\n", ")\n", "\n", "t.loc[\"timflow\"] = np.append(cal.parameters[\"optimal\"].values, cal.rmse())\n", + "t.loc[\"timflow weights\"] = np.append(\n", + " cal2.parameters[\"optimal\"].values, cal2.rmse(weighted=False)\n", + ")\n", "t.loc[\"AQTESOLV\"] = [49.286, 4.559e-05, 745.156, 0.007245]\n", "t.loc[\"MLU\"] = [45.186, 3.941e-05, 769.200, 0.005941]\n", "t.loc[\"Hantush\"] = [45.332, 4.762e-5, 331.141, 0.005917]\n", "\n", "t_formatted = t.style.format(\n", - " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"c [d]\": \"{:.0f}\", \"RMSE [m]\": \"{:.2f}\"}\n", + " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"c [d]\": \"{:.0f}\", \"RMSE [m]\": \"{:.4f}\"}\n", ")\n", "t_formatted" ] diff --git a/docs/transient/04pumpingtests/leaky2_hardixveld.ipynb b/docs/transient/04pumpingtests/leaky2_hardixveld.ipynb index 981be67..930504f 100644 --- a/docs/transient/04pumpingtests/leaky2_hardixveld.ipynb +++ b/docs/transient/04pumpingtests/leaky2_hardixveld.ipynb @@ -2,14 +2,26 @@ "cells": [ { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "# 2. Leaky Aquifer Test - Hardinxveld " ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Import packages" ] @@ -17,7 +29,13 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", @@ -32,26 +50,27 @@ { "cell_type": "markdown", "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, "tags": [] }, "source": [ "### Introduction and Conceptual Model\n", "\n", - "This example, taken from MLU examples (Carlson and Randall, 2012), is a pumping test from Hardinxveld-Giessendam, Netherlands, in 1981 to quantify the head-loss at each pumping well by clogging and to assess the transmissivity variation in the area.\n", - "\n", - "The hydrogeological conceptualization can be described as the following:\n", - "* The first ten meters depth is an aquitard\n", - "* Followed by the first aquifer from 10 to 37 m depth, this is also the test aquifer.\n", - "* A new aquitard is present from 37 m depth to 68 m depth\n", - "* A final aquifer is from 68 to 88 m depth.\n", - "* Below 88 m depth the formations are considered an aquiclude\n", - "\n", - "Five pumping wells are screened in the first aquifer. The drawdown of one of them is available in the MLU documentation (Carlson & Randall, 2012). The provided pumping well was operated for 20 minutes at 1848 m3/d. Drawdown was recorded for a total of 50 minutes during and after pumping. The radius of the pumped well is 0.155 m." + "Consider a pumping test conducted near the town of Hardinxveld-Giessendam in the Netherlands in 1981. The objective was to quantify the aquifer parameters and the head-loss at the pumping well caused by clogging.\n", + "The subsurface consists of a 10 meter thick aquitard, followed by a 27 m aquifer. Below the aquifer is another aquitard of 31 m thickness and a 20 m thick second aquifer. \n", + "Five pumping wells are screened in the first aquifer. Drawdown measurements at of one of wells is obtained from the MLU documentation (Carlson & Randall, 2012). The pumping well was operated for 20 minutes at 1848 m$^3$/d. Drawdown was measured inside the pumping well during pumping and during 30 minutes of recovery. The radius of the pumped well is 0.155 m. Two `timflow` models are calibrated. The first one has an impermeable top and bottom. For the second model, the effect of the semi-confining top and the aquitard and aquifer below the pumped aquifer are combined into one semi-confining layer. By trial and error, the resistance of this semi-confining layer was set to $c=2000$ d. For both models, three parameters are calibrated: the hydraulic conductivity and specific storage of the aquifer and the skin resistance of the well. " ] }, { "cell_type": "markdown", "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, "tags": [] }, "source": [ @@ -60,7 +79,13 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Load data" ] @@ -68,7 +93,13 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "data = np.loadtxt(\"data/recovery.txt\", skiprows=1)\n", @@ -78,7 +109,13 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Parameters and model" ] @@ -86,79 +123,167 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ + "# known parameters\n", "H = 27 # aquifer thickness, m\n", "zt = -10 # upper boundary of aquifer, m\n", "zb = zt - H # lower boundary of the aquifer, m\n", "rw = 0.155 # well screen radius, m\n", "Q = 1848 # constant discharge rate, m^3/d\n", - "t0 = 0.013889 # time stop pumping, d" + "tstop = 0.013889 # time at which pumping stops, d" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ + "# model with impermeable top and bottom\n", "ml = tft.ModelMaq(\n", - " kaq=[50, 40],\n", - " z=[0, zt, zb, -68, -88],\n", - " c=[1000, 1000],\n", - " Saq=[1e-4, 5e-5],\n", - " topboundary=\"semi\",\n", + " kaq=[50],\n", + " z=[zt, zb],\n", + " Saq=[1e-4],\n", + " topboundary=\"conf\",\n", " tmin=1e-4,\n", - " tmax=0.04,\n", + " tmax=1,\n", ")\n", - "w = tft.Well(ml, xw=0, yw=0, rw=rw, res=1, tsandQ=[(0, Q), (t0, 0)], layers=0)\n", + "w = tft.Well(ml, xw=0, yw=0, rw=rw, res=1, tsandQ=[(0, Q), (tstop, 0)], layers=0)\n", "ml.solve()" ] }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Estimate aquifer parameters\n", - "The parameters to be calibrated are the hydraulic conductivity and specific storage of the first layer, and the skin resistance of the well. The parameters of the aquitards and the second aquifer are kept fixed. " - ] - }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "cal = tft.Calibrate(ml)\n", "cal.set_parameter(name=\"kaq\", initial=50, pmin=0, layers=0)\n", - "cal.set_parameter(name=\"Saq\", initial=1e-4, pmin=0, layers=0)\n", + "cal.set_parameter(name=\"Saq\", initial=1e-4, pmin=0, pmax=1e-3, layers=0)\n", "\n", "cal.set_parameter_by_reference(name=\"res\", parameter=w.res[:], initial=1, pmin=0)\n", "cal.seriesinwell(name=\"obs\", element=w, t=to, h=ho)\n", - "cal.fit(report=True)" + "cal.fit(report=False)" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ - "display(cal.parameters)\n", - "print(\"RMSE:\", cal.rmse())" + "display(cal.parameters.loc[:, [\"optimal\"]])\n", + "print(f\"RMSE: {cal.rmse():.5f} m\")" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# model with semi-confining top\n", + "ml1 = tft.ModelMaq(\n", + " kaq=[50],\n", + " z=[0, zt, zb],\n", + " Saq=[1e-4],\n", + " c=2000,\n", + " topboundary=\"semi\",\n", + " tmin=1e-4,\n", + " tmax=1,\n", + ")\n", + "w1 = tft.Well(ml1, xw=0, yw=0, rw=rw, res=1, tsandQ=[(0, Q), (tstop, 0)], layers=0)\n", + "ml1.solve()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "cal1 = tft.Calibrate(ml1)\n", + "cal1.set_parameter(name=\"kaq\", initial=50, pmin=0, layers=0)\n", + "cal1.set_parameter(name=\"Saq\", initial=1e-4, pmin=0, pmax=1e-3, layers=0)\n", + "\n", + "cal1.set_parameter_by_reference(name=\"res\", parameter=w1.res[:], initial=1, pmin=0)\n", + "cal1.seriesinwell(name=\"obs\", element=w1, t=to, h=ho)\n", + "cal1.fit(report=False)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "display(cal1.parameters[\"optimal\"])\n", + "print(f\"RMSE: {cal1.rmse():.5f} m\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "tm = np.logspace(np.log10(to[0]), np.log10(to[-1]), 100)\n", - "hm = w.headinside(tm)\n", - "plt.loglog(to, -ho, \".\", label=\"obs - pumping well\")\n", - "plt.loglog(tm, -hm[0], label=\"timflow results\")\n", + "h = w.headinside(tm)\n", + "h1 = w1.headinside(tm)\n", + "plt.loglog(to, -ho, \".\", label=\"observed\")\n", + "plt.loglog(tm, -h[0], label=\"timflow confined\")\n", + "plt.loglog(tm, -h1[0], label=\"timflow semi-confined\")\n", "plt.xlabel(\"time [d]\")\n", "plt.ylabel(\"drawdown [m]\")\n", "plt.legend()\n", @@ -168,34 +293,55 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Comparison of results" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "The performance of `timflow` was evaluated by comparison with MLU (Carlson and Randall, 2012). The MLU parameters are higher than the ones in `timflow`. " + "The performance of `timflow` is compared to MLU (Carlson and Randall, 2012). The MLU parameters are similar to the `timflow` parameters for the confined case. The semi-confined case gives somewhat lower RMSE and a slightly better fit. It is noted that this model is exceedingly sensitive to the specified moment that the well is turned off. Even a change of 1 second will produce significantly different results. " ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, "outputs": [], "source": [ "t = pd.DataFrame(\n", - " columns=[\"k [m/d]\", \"Ss [1/m]\", \"res [d]\", \"RMSE [m]\"],\n", - " index=[\"timflow\", \"MLU\"],\n", + " columns=[\"k [m/d]\", \"Ss [1/m]\", \"res\", \"RMSE [m]\"],\n", + " index=[\"timflow confined\", \"timflow semi-confined\", \"MLU\"],\n", ")\n", "\n", - "t.loc[\"timflow\"] = np.append(cal.parameters[\"optimal\"].values, cal.rmse())\n", + "t.loc[\"timflow confined\"] = np.append(cal.parameters[\"optimal\"].values, cal.rmse())\n", + "t.loc[\"timflow semi-confined\"] = np.append(cal1.parameters[\"optimal\"].values, cal1.rmse())\n", "t.loc[\"MLU\"] = [51.530, 8.16e-04, 0.022, 0.00756]\n", "\n", "t_formatted = t.style.format(\n", - " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"res [d]\": \"{:.2f}\", \"RMSE [m]\": \"{:.2f}\"}\n", + " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"res\": \"{:.3f}\", \"RMSE [m]\": \"{:.4f}\"}\n", ")\n", "t_formatted" ] diff --git a/docs/transient/04pumpingtests/leaky3_texashill.ipynb b/docs/transient/04pumpingtests/leaky3_texashill.ipynb new file mode 100644 index 0000000..7ed25dc --- /dev/null +++ b/docs/transient/04pumpingtests/leaky3_texashill.ipynb @@ -0,0 +1,359 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "# 3. Leaky Aquifer Test - Texas Hill" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Import packages" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "import pandas as pd\n", + "\n", + "import timflow.transient as tft\n", + "\n", + "plt.rcParams[\"figure.figsize\"] = [5, 3]" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Introduction and Conceptual Model\n", + "\n", + "This example concerns a pumping test conducted in the location of Texas Hill and is taken from the AQTESOLV examples (Duffield, 2007). A pumping well is screened in an aquifer located between 20 ft and 70 ft depths. The aquifer is overlain by an aquitard. The formation at the base of the aquifer is considered impermeable.\n", + "\n", + "Three observation wells are located at 40, 80, and 160 ft distance. They are called OW1, OW2, and OW3, respectively. Pumping lasted for 420 minutes at a rate of 4488 gallons per minute. " + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Load data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# data of OW1\n", + "data1 = np.loadtxt(\"data/texas40.txt\")\n", + "t1 = data1[:, 0] # days\n", + "h1 = data1[:, 1] # meters\n", + "r1 = 40 * 0.3048 # distance between obs1 to pumping well in m\n", + "\n", + "# data of OW2\n", + "data2 = np.loadtxt(\"data/texas80.txt\")\n", + "t2 = data2[:, 0]\n", + "h2 = data2[:, 1]\n", + "r2 = 80 * 0.3048 # distance between obs2 to pumping well in m\n", + "\n", + "# data of OW3\n", + "data3 = np.loadtxt(\"data/texas160.txt\")\n", + "t3 = data3[:, 0]\n", + "h3 = data3[:, 1]\n", + "r3 = 160 * 0.3048 # distance between obs3 to pumping well in m" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Parameters and model" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# known parameters\n", + "Q = (4488 * 0.00378541) * 60 * 24 # discharge, m^3/d\n", + "b1 = 20 * 0.3048 # overlying aquitard thickness, m\n", + "b2 = 50 * 0.3048 # aquifer thickness, m\n", + "zt = -b1 # top of aquifer, m\n", + "zb = -b1 - b2 # bottom of aquifer, m\n", + "rw = 0.5 * 0.3048 # well radius, m" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "The overlying layer is modeld as an aquitard without storage (```Sll```, the storage parameter, is set to zero in ModelMaq)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "ml = tft.ModelMaq(\n", + " kaq=10,\n", + " z=[0, zt, zb],\n", + " Saq=0.001,\n", + " Sll=0,\n", + " c=10,\n", + " tmin=0.001,\n", + " tmax=1,\n", + " topboundary=\"semi\",\n", + ")\n", + "w = tft.Well(ml, xw=0, yw=0, rw=rw, tsandQ=[(0, Q)], layers=0)\n", + "ml.solve()" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Estimate aquifer parameters\n", + "The hydraulic parameters of the aquifer (`kaq` and `Saq`) and the aquitard (`c`) are estimated using observations in all three observation wells. " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# unknown parameters: kaq, Saq, c\n", + "cal = tft.Calibrate(ml)\n", + "cal.set_parameter(name=\"kaq\", layers=0, initial=10)\n", + "cal.set_parameter(name=\"Saq\", layers=0, initial=1e-4)\n", + "cal.set_parameter(name=\"c\", layers=0, initial=100)\n", + "cal.series(name=\"OW1\", x=r1, y=0, t=t1, h=h1, layer=0)\n", + "cal.series(name=\"OW2\", x=r2, y=0, t=t2, h=h2, layer=0)\n", + "cal.series(name=\"OW3\", x=r3, y=0, t=t3, h=h3, layer=0)\n", + "cal.fit(report=False)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "display(cal.parameters.loc[:, [\"optimal\"]])\n", + "print(f\"RMSE: {cal.rmse():.3f} m\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# plot the fitted curves\n", + "hm1_1 = ml.head(r1, 0, t1)\n", + "hm2_1 = ml.head(r2, 0, t2)\n", + "hm3_1 = ml.head(r3, 0, t3)\n", + "plt.semilogx(t1, h1, \".\", label=\"OW1\")\n", + "plt.semilogx(t1, hm1_1[0], label=\"timflow OW1\")\n", + "plt.semilogx(t2, h2, \".\", label=\"OW2\")\n", + "plt.semilogx(t2, hm2_1[0], label=\"timflow OW2\")\n", + "plt.semilogx(t3, h3, \".\", label=\"OW3\")\n", + "plt.semilogx(t3, hm3_1[0], label=\"timflow OW3\")\n", + "plt.xlabel(\"time (d)\")\n", + "plt.ylabel(\"head change (m)\")\n", + "plt.legend(bbox_to_anchor=(1.05, 1))\n", + "plt.grid()" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Comparison of Results " + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "The aquifer response obtained with `timflow` is compared to the results based on Hantush’s analytical solution (Hantush, 1955), implemented in the software AQTESOLV (Duffield, 2007). The results are almost identical." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "t = pd.DataFrame(\n", + " columns=[\"k [m/d]\", \"Ss [1/m]\", \"c [d]\", \"RMSE [m]\"],\n", + " index=[\"timflow\", \"AQTESOLV\"],\n", + ")\n", + "\n", + "t.loc[\"AQTESOLV\"] = [224.726, 2.125e-4, 43.964, 0.059627]\n", + "t.loc[\"timflow\"] = np.append(cal.parameters[\"optimal\"].values, cal.rmse())\n", + "\n", + "t_formatted = t.style.format(\n", + " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"c [d]\": \"{:.2f}\", \"RMSE [m]\": \"{:.3f}\"}\n", + ")\n", + "t_formatted" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## References\n", + "\n", + "* Carlson F. and Randall J. (2012), MLU: a Windows application for the analysis of aquifer tests and the design of well fields in layered systems, Ground Water 50(4):504–510.\n", + "* Duffield, G.M. (2007), AQTESOLV for Windows Version 4.5 User's Guide, HydroSOLVE, Inc., Reston, VA.\n", + "* Newville, M.,Stensitzki, T., Allen, D.B. and Ingargiola, A. (2014), LMFIT: Non Linear Least-Squares Minimization and Curve Fitting for Python, https://dx.doi.org/10.5281/zenodo.11813, https://lmfit.github.io/lmfit-py/intro.html (last access: August,2021)." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "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.13.5" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/docs/transient/04pumpingtests/slug1_pratt_county.ipynb b/docs/transient/04pumpingtests/slug1_pratt_county.ipynb index bb4150f..92a1eb6 100644 --- a/docs/transient/04pumpingtests/slug1_pratt_county.ipynb +++ b/docs/transient/04pumpingtests/slug1_pratt_county.ipynb @@ -2,14 +2,26 @@ "cells": [ { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "# 1. Slug Test - Pratt County" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Import packages" ] @@ -37,30 +49,53 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Introduction and Conceptual Model" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "This slug test, taken from AQTESOLV examples (Duffield, 2007), was conducted in Pratt County Monitoring Site, US, and reported by Butler (1998). \n", - "\n", - "A partially penetrating well is screened in unconsolidated alluvial deposits, consisting of sand and gravel interbedded by clay. The total thickness of the aquifer is 47.87 m. The screen is located at 16.77 m depth and has a screen length 1.52 m the well radius is 0.125 m, and the casing radius 0.064 m. The slug displacement is 0.671 m. Head change has been recorded at the slug well." + "Consider the slug test conducted in Pratt County Monitoring Site, US, and reported by Butler (1998). \n", + "A partially penetrating well is screened in unconsolidated alluvial deposits consisting of sand and gravel interbedded by clay. The total thickness of the aquifer is 47.87 m. The screen is located at 16.77 m depth and has a screen length of 1.52 m. The well radius is 0.125 m and the casing radius 0.064 m. The slug displacement is 0.671 m. The head change was recorded inside the well." ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "\n" + "" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Load data" ] @@ -68,7 +103,13 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "data = np.loadtxt(\"data/slug.txt\", skiprows=1)\n", @@ -78,7 +119,13 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Parameters and model" ] @@ -86,62 +133,96 @@ { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ - "rw = 0.125 # well radius in m\n", - "rc = 0.064 # well casing radius in m\n", - "L = 1.52 # screen length in m\n", - "b = -47.87 # aquifer thickness in m\n", - "zt = -16.77 # depth to top of screen in m\n", - "H0 = 0.671 # initial displacement in the well in m\n", + "rw = 0.125 # well radius, m\n", + "rc = 0.064 # well casing radius, m\n", + "L = 1.52 # screen length, m\n", + "b = 47.87 # aquifer thickness, m\n", + "zt = -16.77 # elevation to top of screen, m\n", + "H0 = 0.671 # head displacement in the well, m\n", "zb = zt - L # bottom of screen in m" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "The input for the slug test in `timflow` is the added or removed volume. Therefore, we must first convert our measured displacement into volume." + "The total volume if the slut is computed, as this is used as the instantaneous discharge for the well in `timflow`. " ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ - "Q = np.pi * rc**2 * H0\n", - "print(\"slug:\", round(Q, 5), \"m^3\")" + "Q = np.pi * rc**2 * H0 # instantaneous discharge\n", + "print(f\"volume of slug: {Q:.5f} m^3\")" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "We conceptualize the aquifer as a three-layer model, one layer above the screen, one layer at the screen top and bottom and another layer just below it.\n", - "\n", - "The setting of the slug well is slightly different from the pumping well. We detail the differences below:\n", - "* the ```tsandQ``` argument in the ```Well``` object has a different meaning. Instead of meaning the time of start or shutdown and the pumping rate of the pumping well, it means the time a volume is instantaneously added or removed from the well. In our case, we defined it as ```[(0, -Q)]``` where ```0``` is the moment in time when we added the slug and ```Q``` is the added volume. A negative sign means a volume is added. Otherwise, it would mean an extracted volume.\n", - "* the ```wbstype``` argument is set to ```' slug'```, so `timflow` knows the ```tsandQ``` argument means time and instant volumes, instead of pumping rates." + "The aquifer is represented by a three-layer model, one layer above the screen, one layer at interval of the screen top, and one layer below the screen.\n", + "A slug test is simulated by specify the instantaneous volume that is added to the well and by defining the well type `wbstype` as `\"slug\"`. " ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ - "ml = tft.Model3D(kaq=10, z=[0, zt, zb, b], Saq=1e-4, kzoverkh=1, tmin=1e-6, tmax=0.01)\n", + "ml = tft.Model3D(kaq=10, z=[0, zt, zb, -b], Saq=1e-4, kzoverkh=1, tmin=1e-6, tmax=0.01)\n", "w = tft.Well(ml, xw=0, yw=0, rw=rw, rc=rc, tsandQ=[(0, -Q)], layers=1, wbstype=\"slug\")\n", "ml.solve()" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "### Estimate aquifer parameters " + "### Estimate aquifer parameters \n", + "The hydraulic conductivity and specific storage coeffient are calibrated. They are the same for all three layers. " ] }, { @@ -160,23 +241,35 @@ "cal.set_parameter(name=\"kaq\", layers=[0, 1, 2], initial=10)\n", "cal.set_parameter(name=\"Saq\", layers=[0, 1, 2], initial=1e-4)\n", "cal.seriesinwell(name=\"obs\", element=w, t=to, h=ho)\n", - "cal.fit(report=True)" + "cal.fit(report=False)" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ - "display(cal.parameters)\n", - "print(\"RMSE:\", cal.rmse())" + "display(cal.parameters.loc[:, [\"optimal\"]])\n", + "print(f\"RMSE: {cal.rmse():.5f} m\")" ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "outputs": [], "source": [ "tm = np.logspace(np.log10(to[0]), np.log10(to[-1]), 100)\n", @@ -193,22 +286,42 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "### Comparison of results" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ - "Here, we reproduce the work of Yang (2020) to check the `timflow` performance in analysing slug tests. We compare the solution in `timflow` with the KGS analytical model (Hyder et al. 1994) implemented in AQTESOLV (Duffield, 2007). Both models show similarly low RMSE values. However, the estimated hydraulic conductivity and specific storage parameters differ substantially." + "The solution in `timflow` is compared with the KGS analytical model (Hyder et al. 1994) implemented in AQTESOLV (Duffield, 2007). Both models show similarly low RMSE values. However, the estimated hydraulic conductivity and specific storage parameters differ substantially." ] }, { "cell_type": "code", "execution_count": null, - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, "outputs": [], "source": [ "t = pd.DataFrame(\n", @@ -220,7 +333,7 @@ "t.loc[\"AQTESOLV\"] = [4.034, 3.834e-04, 0.002976]\n", "\n", "t_formatted = t.style.format(\n", - " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"RMSE [m]\": \"{:.3f}\"}\n", + " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"RMSE [m]\": \"{:.4f}\"}\n", ")\n", "t_formatted" ] @@ -238,7 +351,7 @@ "## References\n", "* Butler Jr., J.J. (1998), The Design, Performance, and Analysis of Slug Tests, Lewis Publishers, Boca Raton, Florida, 252p.\n", "* Duffield, G.M. (2007), AQTESOLV for Windows Version 4.5 User's Guide, HydroSOLVE, Inc., Reston, VA.\n", - "* Hyder, Z., Butler Jr, J.J., McElwee, C.D., Liu, W. (1994), Slug tests in partially penetrating wells, Water Resources Research 30, 2945–2957." + "* Hyder, Z., Butler Jr, J.J., McElwee, C.D. and Liu, W. (1994), Slug tests in partially penetrating wells, Water Resources Research 30, 2945–2957." ] } ], diff --git a/docs/transient/04pumpingtests/slug2_falling_head.ipynb b/docs/transient/04pumpingtests/slug2_falling_head.ipynb deleted file mode 100644 index dd936cb..0000000 --- a/docs/transient/04pumpingtests/slug2_falling_head.ipynb +++ /dev/null @@ -1,279 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "metadata": { - "editable": true, - "slideshow": { - "slide_type": "" - }, - "tags": [] - }, - "source": [ - "# 2. Slug Test - Falling Head" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Import packages" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "editable": true, - "slideshow": { - "slide_type": "" - }, - "tags": [] - }, - "outputs": [], - "source": [ - "import matplotlib.pyplot as plt\n", - "import numpy as np\n", - "import pandas as pd\n", - "\n", - "import timflow.transient as tft\n", - "\n", - "plt.rcParams[\"figure.figsize\"] = [5, 3]" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Introduction and Conceptual Model\n", - "\n", - "This slug test, taken from the AQTESOLV examples (Duffield, 2007), was reported in Batu (1998). \n", - "\n", - "A well partially penetrates a sandy unconfined aquifer that has a saturated depth of 32.57 ft. The top of the screen is located 0.47 ft below the water table and has 13.8 ft in length. The well and casing radii are 5 and 2 inches, respectively. The slug displacement is 1.48 ft. Head change has been recorded at the slug well." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Load data" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "data = np.loadtxt(\"data/falling_head.txt\", skiprows=2)\n", - "to = data[:, 0] / 60 / 60 / 24 # convert time from seconds to days\n", - "ho = (10 - data[:, 1]) * 0.3048 # convert drawdown from ft to meters" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Parameters and model" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "rw = 5 * 0.0254 # well radius in m (5 inch = 0.127 m)\n", - "rc = 2 * 0.0254 # well casing radius in m (2 inch = 0.0508 m)\n", - "L = 13.8 * 0.3048 # screen length in m (13.8 ft = 4.2 m)\n", - "b = -32.57 * 0.3048 # aquifer thickness in m (-32.57 ft = -9.92 m)\n", - "zt = -0.47 * 0.3048 # depth to top of the screen in m (-0.47 ft = -0.14 m)\n", - "H0 = 1.48 * 0.3048 # initial displacement in the well in m (1.48 ft = 0.4511 m)\n", - "zb = zt - L # bottom of the screen in m" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "convert measured displacement into volume" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "Q = np.pi * rc**2 * H0\n", - "print(f\"slug: {Q:.5f} m^3\")" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "We will create a multi-layer model. For this, we divide the second and third layers into 0.5 m thick layers. " - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# The thickness of each layer is set to be 0.5 m\n", - "z0 = np.arange(zt, zb, -0.5)\n", - "z1 = np.arange(zb, b, -0.5)\n", - "zlay = np.append(z0, z1)\n", - "zlay = np.append(zlay, b)\n", - "zlay = np.insert(zlay, 0, 0)\n", - "nlay = len(zlay) - 1 # number of layers\n", - "Saq = 1e-4 * np.ones(nlay)\n", - "Saq[0] = 0.1" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "ml = tft.Model3D(\n", - " kaq=10, z=zlay, Saq=Saq, kzoverkh=1, tmin=1e-5, tmax=0.01, phreatictop=True\n", - ")\n", - "w = tft.Well(\n", - " ml,\n", - " xw=0,\n", - " yw=0,\n", - " rw=rw,\n", - " tsandQ=[(0, -Q)],\n", - " layers=[1, 2, 3, 4, 5, 6, 7, 8],\n", - " rc=rc,\n", - " wbstype=\"slug\",\n", - ")\n", - "ml.solve()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Estimate aquifer parameters" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "cal = tft.Calibrate(ml)\n", - "cal.set_parameter(name=\"kaq\", layers=list(range(nlay)), initial=10, pmin=0)\n", - "cal.set_parameter(name=\"Saq\", layers=list(range(nlay)), initial=1e-4, pmin=0)\n", - "cal.seriesinwell(name=\"obs\", element=w, t=to, h=ho)\n", - "cal.fit(report=True)" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "display(cal.parameters)\n", - "print(\"RMSE:\", cal.rmse())" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "tm = np.logspace(np.log10(to[0]), np.log10(to[-1]), 100)\n", - "hm = w.headinside(tm)\n", - "plt.semilogx(to, ho / H0, \".\", label=\"obs\")\n", - "plt.semilogx(tm, hm[0] / H0, label=\"timflow\")\n", - "plt.xlabel(\"time [d]\")\n", - "plt.ylabel(\"Normalized head (h/H0)\")\n", - "plt.title(\"Model results - multi-layer model\")\n", - "plt.legend()\n", - "plt.grid()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Comparison of results\n", - "Here, we reproduce the work of Yang (2020) to check the `timflow` performance in analysing slug tests. We compare the solution in `timflow` with the KGS analytical model (Hyder et al. 1994) implemented in AQTESOLV (Duffield, 2007). AQTESOLV parameters are quite different from the set parameters in `timflow`. Furthermore, AQTESOLV also has a better RMSE performance." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "t = pd.DataFrame(\n", - " columns=[\"k [m/d]\", \"Ss [1/m]\", \"RMSE [m]\"],\n", - " index=[\"timflow-multi\", \"AQTESOLV\"],\n", - ")\n", - "\n", - "t.loc[\"timflow-multi\"] = np.append(cal.parameters[\"optimal\"].values, cal.rmse())\n", - "t.loc[\"AQTESOLV\"] = [2.616, 7.894e-5, 0.001197]\n", - "\n", - "t_formatted = t.style.format(\n", - " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"RMSE [m]\": \"{:.3f}\"}\n", - ")\n", - "t_formatted" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "editable": true, - "slideshow": { - "slide_type": "" - }, - "tags": [] - }, - "source": [ - "## References\n", - "\n", - "* Batu, V. (1998), Aquifer hydraulics: a comprehensive guide to hydrogeologic data analysis, John Wiley & Sons\n", - "* Duffield, G.M. (2007), AQTESOLV for Windows Version 4.5 User's Guide, HydroSOLVE, Inc., Reston, VA.\n", - "* Hyder, Z., Butler Jr, J.J., McElwee, C.D. and Liu, W. (1994), Slug tests in partially penetrating wells, Water Resources Research 30, 2945–2957." - ] - } - ], - "metadata": { - "kernelspec": { - "display_name": "Python 3 (ipykernel)", - "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.13.5" - } - }, - "nbformat": 4, - "nbformat_minor": 4 -} diff --git a/docs/transient/04pumpingtests/slug2_multiwell.ipynb b/docs/transient/04pumpingtests/slug2_multiwell.ipynb new file mode 100644 index 0000000..90fee63 --- /dev/null +++ b/docs/transient/04pumpingtests/slug2_multiwell.ipynb @@ -0,0 +1,352 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "# 2. Slug Test for Confined Aquifer - Multi-well" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Import packages" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "import pandas as pd\n", + "\n", + "import timflow.transient as tft\n", + "\n", + "plt.rcParams[\"figure.figsize\"] = [5, 3]" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "## Introduction and Conceptual Model\n", + "\n", + "A well (Ln-2) fully penetrates a sandy confined aquifer of 6.1 m thickness. Additionally, a fully penetrating observation well (Ln-3) is located 6.45 m away from the test well. The slug displacement is 2.798 m. Head change was recorded at the slug well and the observation well. The well and casing radii of the slug well are 0.102 and 0.051 m, respectively. For the observation well, they are 0.051 and 0.025 m, respectively." + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Load data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "data1 = np.loadtxt(\"data/ln-2.txt\")\n", + "t1 = data1[:, 0] / 60 / 60 / 24 # convert time from seconds to days\n", + "h1 = data1[:, 1]\n", + "\n", + "data2 = np.loadtxt(\"data/ln-3.txt\")\n", + "t2 = data2[:, 0] / 60 / 60 / 24\n", + "h2 = data2[:, 1]" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Parameters and model" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# known parameters\n", + "H0 = 2.798 # initial displacement, m\n", + "b = 6.1 # aquifer thickness, m\n", + "rw1 = 0.102 # well radius of Ln-2 Well, m\n", + "rw2 = 0.071 # well radius of observation Ln-3 Well, m\n", + "rc1 = 0.051 # casing radius of Ln-2 Well, m\n", + "r = 6.45 # distance from observation well to test well, m" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "Converting slug displacement into volume" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "Q = np.pi * rc1**2 * H0\n", + "print(f\"Slug: {Q:.5f} m^3\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "ml = tft.ModelMaq(kaq=10, z=[0, -b], Saq=1e-4, tmin=1e-5, tmax=0.01)\n", + "w = tft.Well(ml, xw=0, yw=0, rw=rw1, rc=rc1, tsandQ=[(0, -Q)], layers=0, wbstype=\"slug\")\n", + "ml.solve()" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Estimate aquifer parameters\n", + "The hydraulic conductivity and specific storage are calibrated." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "# unknown parameters: kaq, Saq\n", + "cal = tft.Calibrate(ml)\n", + "cal.set_parameter(name=\"kaq\", layers=0, initial=10)\n", + "cal.set_parameter(name=\"Saq\", layers=0, initial=1e-4)\n", + "cal.seriesinwell(name=\"Ln-2\", element=w, t=t1, h=h1)\n", + "cal.series(name=\"Ln-3\", x=r, y=0, layer=0, t=t2, h=h2)\n", + "cal.fit()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "display(cal.parameters.loc[:, [\"optimal\"]])\n", + "print(f\"RMSE: {cal.rmse():.5f} m\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "outputs": [], + "source": [ + "hm1 = w.headinside(t1)\n", + "hm2 = ml.head(r, 0, t2, layers=0)\n", + "plt.semilogx(t1, h1 / H0, \".\", label=\"obs ln-2\")\n", + "plt.semilogx(t1, hm1[0] / H0, label=\"timflow ln-2\")\n", + "plt.semilogx(t2, h2 / H0, \".\", label=\"obs ln-3\")\n", + "plt.semilogx(t2, hm2[0] / H0, label=\"timflow ln-3\")\n", + "plt.xlabel(\"time [d]\")\n", + "plt.ylabel(\"Normalized Head: h/H0\")\n", + "plt.title(\"Model Results - Single layer model\")\n", + "plt.legend()\n", + "plt.grid()" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### Comparison of results\n", + "\n", + "The values in `timflow` are compared with the KGS analytical model (Hyder et al. 1994) implemented in AQTESOLV (Duffield, 2007), and to the MLU model (Carlson & Randall, 2012). The parameters of all three models are very similar." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "t = pd.DataFrame(\n", + " columns=[\"k [m/d]\", \"Ss [1/m]\", \"RMSE [m]\"],\n", + " index=[\"timflow\", \"AQTESOLV\", \"MLU\"],\n", + ")\n", + "\n", + "t.loc[\"timflow\"] = np.append(cal.parameters[\"optimal\"].values, cal.rmse())\n", + "t.loc[\"AQTESOLV\"] = [1.166, 9.368e-06, 0.009151]\n", + "t.loc[\"MLU\"] = [1.311, 8.197e-06, 0.010373]\n", + "\n", + "t_formatted = t.style.format(\n", + " {\"k [m/d]\": \"{:.2f}\", \"Ss [1/m]\": \"{:.2e}\", \"RMSE [m]\": \"{:.3f}\"}\n", + ")\n", + "t_formatted" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "### References\n", + "\n", + "* Butler Jr., J.J. (1998), The Design, Performance, and Analysis of Slug Tests, Lewis Publishers, Boca Raton, Florida, 252p.\n", + "* Carlson F. and Randall J. (2012), MLU: a Windows application for the analysis of aquifer tests and the design of well fields in layered systems, Ground Water 50(4):504–510.\n", + "* Duffield, G.M. (2007), AQTESOLV for Windows Version 4.5 User's Guide, HydroSOLVE, Inc., Reston, VA.\n", + "* Hyder, Z., Butler Jr, J.J., McElwee, C.D. and Liu, W. (1994), Slug tests in partially penetrating wells, Water Resources Research 30, 2945–2957." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "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.13.5" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +}