Release 6.2

.git-blame-ignore-revs

@@ -0,0 +1,2 @@
+# Migrate code style to black
+998b86fc128e35f1506e3d4facb169c558679bde

.github/workflows/ci.yml

@@ -5,9 +5,11 @@ on:
   push:
     branches:
       - master
+      - dev
   pull_request:
     branches:
       - master
+      - dev
 
 jobs:
   test:

CITATION.cff

@@ -0,0 +1,25 @@
+cff-version: 1.2.0
+message: "If you use this software, please cite it as below."
+authors:
+    - family-names: "Bakker"
+      given-names: "Mark"
+      orcid: "https://orcid.org/0000-0002-5629-2861"
+title: "TimML"
+url: "https://github.com/mbakker7/timml"
+
+preferred-citation:
+  type: article
+  authors:
+    - family-names: "Bakker"
+      given-names: "Mark"
+      orcid: "https://orcid.org/0000-0002-5629-2861"
+    - family-names: "Strack"  
+      given-names: "Otto D.L."
+  doi: "10.1016/S0022-1694(02)00319-0"
+  journal: "Journal of Hydrology"
+  title: "Analytic elements for multiaquifer flow"
+  volume: 271
+  number: 1
+  start: 119
+  end: 129
+  year: 2003

notebooks/benchmarking_besselaes.ipynb

@@ -24,7 +24,7 @@
     "import pandas as pd\n",
     "\n",
     "import matplotlib as mpl\n",
-    "import matplotlib.pyplot as plt\n"
+    "import matplotlib.pyplot as plt"
    ]
   },
   {
@@ -72,7 +72,7 @@
     "def ctxt(val):\n",
     "    bg = cmap(norm(val))\n",
     "    c = \"White\" if np.mean(bg[:3]) < 0.4 else \"Black\"\n",
-    "    return f\"color: {c}\"\n"
+    "    return f\"color: {c}\""
    ]
   },
   {
@@ -144,7 +144,7 @@
     "    ls2 = timml.HeadLineSinkString(ml, xy=xy2, hls=hls2, layers=0)\n",
     "    ls3 = timml.HeadLineSinkString(ml, xy=xy3, hls=hls3, layers=0)\n",
     "    ls4 = timml.HeadLineSinkString(ml, xy=xy4, hls=hls4, layers=0)\n",
-    "    \n",
+    "\n",
     "    return ml"
    ]
   },
@@ -179,31 +179,29 @@
     "print(timml.bessel.bessel)  # print module name to check if switch works\n",
     "\n",
     "for nproc in df_ls.index:\n",
-    "    \n",
     "    ml = model1()\n",
-    "    \n",
+    "\n",
     "    if nproc == 0:\n",
     "        t_n = %timeit -o ml.solve(silent=True)\n",
     "    else:\n",
     "        t_n = %timeit -o ml.solve_mp(nproc=nproc, silent=True)\n",
-    "    \n",
+    "\n",
     "    df_ls.loc[nproc, \"numba\"] = t_n.average\n",
     "\n",
     "# %%\n",
     "# Fortran\n",
-    "timml.bessel.set_bessel_method(method=\"fortran\")  \n",
+    "timml.bessel.set_bessel_method(method=\"fortran\")\n",
     "print(timml.bessel.bessel)  # print module name to check if switch works\n",
     "\n",
     "for nproc in df_ls.index:\n",
-    "\n",
     "    ml = model1()\n",
     "\n",
     "    if nproc == 0:\n",
     "        t_f = %timeit -o ml.solve(silent=True)\n",
     "    else:\n",
     "        t_f = %timeit -o ml.solve_mp(nproc=nproc, silent=True)\n",
-    "    \n",
-    "    df_ls.loc[nproc, \"fortran\"] = t_f.average\n"
+    "\n",
+    "    df_ls.loc[nproc, \"fortran\"] = t_f.average"
    ]
   },
   {
@@ -308,7 +306,7 @@
     "norm = mpl.colors.LogNorm(vmin=50, vmax=100)\n",
     "\n",
     "df_ls.index.name = \"nproc\"\n",
-    "(df_ls * 1e3).style.applymap(cbg).applymap(ctxt)\n"
+    "(df_ls * 1e3).style.applymap(cbg).applymap(ctxt)"
    ]
   },
   {
@@ -471,15 +469,45 @@
     "    layers = np.arange(np.sum(z_dw <= ml.aq.zaqbot))\n",
     "    last_lay_dw = layers[-1]\n",
     "\n",
-    "    inhom = timml.BuildingPit(ml, xy, kaq=kh_arr, z=z[1:], topboundary=\"conf\", \n",
-    "                            c=c[1:], order=4, ndeg=3, layers=layers)\n",
+    "    inhom = timml.BuildingPit(\n",
+    "        ml,\n",
+    "        xy,\n",
+    "        kaq=kh_arr,\n",
+    "        z=z[1:],\n",
+    "        topboundary=\"conf\",\n",
+    "        c=c[1:],\n",
+    "        order=4,\n",
+    "        ndeg=3,\n",
+    "        layers=layers,\n",
+    "    )\n",
     "\n",
-    "    timml.HeadLineSink(ml, x1=-l/2+offset, y1=b/2-offset, x2=l/2-offset, y2=b/2-offset, hls=h_bem, \n",
-    "                    layers=np.arange(last_lay_dw+1))\n",
-    "    timml.HeadLineSink(ml, x1=-l/2+offset, y1=0, x2=l/2-offset, y2=0, hls=h_bem, \n",
-    "                    layers=np.arange(last_lay_dw+1))\n",
-    "    timml.HeadLineSink(ml, x1=-l/2+offset, y1=-b/2+offset, x2=l/2-offset, y2=-b/2+offset, hls=h_bem, \n",
-    "                    layers=np.arange(last_lay_dw+1))\n",
+    "    timml.HeadLineSink(\n",
+    "        ml,\n",
+    "        x1=-l / 2 + offset,\n",
+    "        y1=b / 2 - offset,\n",
+    "        x2=l / 2 - offset,\n",
+    "        y2=b / 2 - offset,\n",
+    "        hls=h_bem,\n",
+    "        layers=np.arange(last_lay_dw + 1),\n",
+    "    )\n",
+    "    timml.HeadLineSink(\n",
+    "        ml,\n",
+    "        x1=-l / 2 + offset,\n",
+    "        y1=0,\n",
+    "        x2=l / 2 - offset,\n",
+    "        y2=0,\n",
+    "        hls=h_bem,\n",
+    "        layers=np.arange(last_lay_dw + 1),\n",
+    "    )\n",
+    "    timml.HeadLineSink(\n",
+    "        ml,\n",
+    "        x1=-l / 2 + offset,\n",
+    "        y1=-b / 2 + offset,\n",
+    "        x2=l / 2 - offset,\n",
+    "        y2=-b / 2 + offset,\n",
+    "        hls=h_bem,\n",
+    "        layers=np.arange(last_lay_dw + 1),\n",
+    "    )\n",
     "\n",
     "    return ml"
    ]
@@ -511,33 +539,31 @@
    "source": [
     "df_bp = pd.DataFrame(index=nproc_list, columns=[\"numba\", \"fortran\"])\n",
     "\n",
-    "timml.bessel.set_bessel_method(method=\"numba\")  \n",
+    "timml.bessel.set_bessel_method(method=\"numba\")\n",
     "print(timml.bessel.bessel)\n",
     "\n",
     "for nproc in df_bp.index:\n",
-    "\n",
     "    ml = model2()\n",
     "\n",
     "    if nproc == 0:\n",
     "        t_n = %timeit -o ml.solve(silent=True)\n",
     "    else:\n",
     "        t_n = %timeit -o ml.solve_mp(nproc=nproc, silent=True)\n",
-    "    \n",
+    "\n",
     "    df_bp.loc[nproc, \"numba\"] = t_n.average\n",
     "\n",
     "timml.bessel.set_bessel_method(method=\"fortran\")\n",
     "print(timml.bessel.bessel)\n",
     "\n",
     "for nproc in df_bp.index:\n",
-    "\n",
     "    ml = model2()\n",
-    "    \n",
+    "\n",
     "    if nproc == 0:\n",
     "        t_f = %timeit -o ml.solve(silent=True)\n",
     "    else:\n",
     "        t_f = %timeit -o ml.solve_mp(nproc=nproc, silent=True)\n",
-    "    \n",
-    "    df_bp.loc[nproc, \"fortran\"] = t_f.average\n"
+    "\n",
+    "    df_bp.loc[nproc, \"fortran\"] = t_f.average"
    ]
   },
   {
@@ -770,15 +796,48 @@
    "source": [
     "def model3():\n",
     "    ml = timml.ModelMaq(kaq=[10, 20], z=[20, 0, -10, -30], c=[4000])\n",
-    "    xy1 = [(0, 600), (-100, 400), (-100, 200), (100, 100), (300, 100), (500, 100),\n",
-    "        (700, 300), (700, 500), (600, 700), (400, 700), (200, 600)]\n",
-    "    p1 = timml.PolygonInhomMaq(ml, xy=xy1, \n",
-    "                        kaq=[2, 80], z=[20, 0, -10, -30], c=[500], \n",
-    "                        topboundary='conf', order=4, ndeg=2)\n",
-    "    xy2 = [(0, 600), (200, 600), (400, 700), (400, 900), (200, 1100), (0, 1000), (-100, 800)]\n",
-    "    p2 = timml.PolygonInhomMaq(ml, xy=xy2, \n",
-    "                        kaq=[2, 8], z=[20, 0, -10, -30], c=[50], \n",
-    "                        topboundary='conf', order=4, ndeg=2)\n",
+    "    xy1 = [\n",
+    "        (0, 600),\n",
+    "        (-100, 400),\n",
+    "        (-100, 200),\n",
+    "        (100, 100),\n",
+    "        (300, 100),\n",
+    "        (500, 100),\n",
+    "        (700, 300),\n",
+    "        (700, 500),\n",
+    "        (600, 700),\n",
+    "        (400, 700),\n",
+    "        (200, 600),\n",
+    "    ]\n",
+    "    p1 = timml.PolygonInhomMaq(\n",
+    "        ml,\n",
+    "        xy=xy1,\n",
+    "        kaq=[2, 80],\n",
+    "        z=[20, 0, -10, -30],\n",
+    "        c=[500],\n",
+    "        topboundary=\"conf\",\n",
+    "        order=4,\n",
+    "        ndeg=2,\n",
+    "    )\n",
+    "    xy2 = [\n",
+    "        (0, 600),\n",
+    "        (200, 600),\n",
+    "        (400, 700),\n",
+    "        (400, 900),\n",
+    "        (200, 1100),\n",
+    "        (0, 1000),\n",
+    "        (-100, 800),\n",
+    "    ]\n",
+    "    p2 = timml.PolygonInhomMaq(\n",
+    "        ml,\n",
+    "        xy=xy2,\n",
+    "        kaq=[2, 8],\n",
+    "        z=[20, 0, -10, -30],\n",
+    "        c=[50],\n",
+    "        topboundary=\"conf\",\n",
+    "        order=4,\n",
+    "        ndeg=2,\n",
+    "    )\n",
     "    rf = timml.Constant(ml, xr=1000, yr=0, hr=40)\n",
     "    uf = timml.Uflow(ml, slope=0.002, angle=-45)\n",
     "    w = timml.Well(ml, xw=400, yw=400, Qw=500, rw=0.2, layers=0)\n",
@@ -812,33 +871,31 @@
    "source": [
     "df = pd.DataFrame(index=nproc_list, columns=[\"numba\", \"fortran\"])\n",
     "\n",
-    "timml.bessel.set_bessel_method(method=\"numba\")  \n",
+    "timml.bessel.set_bessel_method(method=\"numba\")\n",
     "print(timml.bessel.bessel)\n",
     "\n",
     "for nproc in df.index:\n",
-    "\n",
     "    ml = model3()\n",
     "\n",
     "    if nproc == 0:\n",
     "        t_n = %timeit -o ml.solve(silent=True)\n",
     "    else:\n",
     "        t_n = %timeit -o ml.solve_mp(nproc=nproc, silent=True)\n",
-    "    \n",
+    "\n",
     "    df.loc[nproc, \"numba\"] = t_n.average\n",
     "\n",
-    "timml.bessel.set_bessel_method(method=\"fortran\")  \n",
+    "timml.bessel.set_bessel_method(method=\"fortran\")\n",
     "print(timml.bessel.bessel)\n",
     "\n",
     "for nproc in df.index:\n",
-    "\n",
     "    ml = model3()\n",
-    "    \n",
+    "\n",
     "    if nproc == 0:\n",
     "        t_f = %timeit -o ml.solve(silent=True)\n",
     "    else:\n",
     "        t_f = %timeit -o ml.solve_mp(nproc=nproc, silent=True)\n",
-    "    \n",
-    "    df.loc[nproc, \"fortran\"] = t_f.average\n"
+    "\n",
+    "    df.loc[nproc, \"fortran\"] = t_f.average"
    ]
   },
   {

notebooks/circareasink_example.ipynb

@@ -67,9 +67,9 @@
     "for i in range(len(x)):\n",
     "    qx[i], qy = ml.disvec(x[i], 1e-6)\n",
     "plt.plot(x, qx)\n",
-    "qxb = N * np.pi * R ** 2 / (2 * np.pi * R)\n",
-    "plt.axhline(qxb, color='r', ls='--')\n",
-    "plt.axhline(-qxb, color='r', ls='--');"
+    "qxb = N * np.pi * R**2 / (2 * np.pi * R)\n",
+    "plt.axhline(qxb, color=\"r\", ls=\"--\")\n",
+    "plt.axhline(-qxb, color=\"r\", ls=\"--\");"
    ]
   },
   {
@@ -101,7 +101,7 @@
    "source": [
     "N = 0.001\n",
     "R = 100\n",
-    "Q = N * np.pi * R ** 2\n",
+    "Q = N * np.pi * R**2\n",
     "ml = tml.ModelMaq(kaq=5, z=[10, 0])\n",
     "ca = tml.CircAreaSink(ml, xc=-200, yc=0, R=100, N=0.001)\n",
     "w = tml.Well(ml, 200, 0, Qw=Q, rw=0.1)\n",
@@ -140,7 +140,7 @@
    "source": [
     "N = 0.001\n",
     "R = 100\n",
-    "Q = N * np.pi * R ** 2\n",
+    "Q = N * np.pi * R**2\n",
     "ml = tml.ModelMaq(kaq=5, z=[10, 0])\n",
     "ca = tml.CircAreaSink(ml, xc=-200, yc=0, R=100, N=0.001)\n",
     "w = tml.Well(ml, 200, 0, Qw=Q, rw=0.1)\n",
@@ -184,15 +184,15 @@
    "source": [
     "N = 0.001\n",
     "R = 100\n",
-    "Q = N * np.pi * R ** 2\n",
+    "Q = N * np.pi * R**2\n",
     "ml = tml.ModelMaq(kaq=[5, 20], z=[20, 12, 10, 0], c=[1000])\n",
     "ca = tml.CircAreaSink(ml, xc=0, yc=0, R=100, N=0.001)\n",
     "w = tml.Well(ml, 0, 0, Qw=Q, rw=0.1, layers=1)\n",
     "rf = tml.Constant(ml, 1000, 0, 0)\n",
     "ml.solve()\n",
     "x = np.linspace(-200, 200, 100)\n",
     "h = ml.headalongline(x, 0)\n",
-    "plt.plot(x, h[0]);\n",
+    "plt.plot(x, h[0])\n",
     "plt.plot(x, h[1]);"
    ]
   },
@@ -215,7 +215,7 @@
    "source": [
     "x = np.linspace(-1000, 1000, 101)\n",
     "h = ml.headalongline(x, 0)\n",
-    "plt.plot(x, h[0]);\n",
+    "plt.plot(x, h[0])\n",
     "plt.plot(x, h[1]);"
    ]
   },

notebooks/circular_buildingpit.ipynb

@@ -39,7 +39,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "R = 100.0 # radius of building pit"
+    "R = 100.0  # radius of building pit"
    ]
   },
   {
@@ -346,16 +346,16 @@
     "Qr = disrvec(r) * 2 * np.pi * r\n",
     "plt.figure(figsize=(10, 3))\n",
     "plt.subplot(121)\n",
-    "plt.plot(r, h, label='head with wall')\n",
-    "plt.plot(r, head_nowall(r), '--', label='head no wall')\n",
-    "plt.xlabel('r (m)')\n",
-    "plt.ylabel('head (m)')\n",
+    "plt.plot(r, h, label=\"head with wall\")\n",
+    "plt.plot(r, head_nowall(r), \"--\", label=\"head no wall\")\n",
+    "plt.xlabel(\"r (m)\")\n",
+    "plt.ylabel(\"head (m)\")\n",
     "plt.grid()\n",
     "plt.legend()\n",
     "plt.subplot(122)\n",
     "plt.plot(r, Qr)\n",
-    "plt.xlabel('r (m)')\n",
-    "plt.ylabel('$Q_r$ (m$^2$/d)')\n",
+    "plt.xlabel(\"r (m)\")\n",
+    "plt.ylabel(\"$Q_r$ (m$^2$/d)\")\n",
     "plt.grid()\n",
     "plt.tight_layout()"
    ]