diff --git a/docs/03examples/watex_example.ipynb b/docs/03examples/watex_example.ipynb
index e10271e..b21b1ac 100644
--- a/docs/03examples/watex_example.ipynb
+++ b/docs/03examples/watex_example.ipynb
@@ -323,7 +323,9 @@
     "    where=\"post\",\n",
     "    label=\"shift -1\",\n",
     ")\n",
-    "mid = (data_norm.loc[:, \"rivier\"] + data_norm.loc[:, \"rivier\"].shift(-1).values).divide(2)\n",
+    "mid = (data_norm.loc[:, \"rivier\"] + data_norm.loc[:, \"rivier\"].shift(-1).values).divide(\n",
+    "    2\n",
+    ")\n",
     "plt.step(\n",
     "    mid.loc[:1].index, mid.loc[:1], where=\"post\", label=\"midpoint\", color=\"k\", lw=1.0\n",
     ")\n",
diff --git a/ttim/besselnumba.py b/ttim/besselnumba.py
index bb37dd9..557129d 100644
--- a/ttim/besselnumba.py
+++ b/ttim/besselnumba.py
@@ -364,7 +364,7 @@ def bessells_gauss(x, y, z1, z2, lab):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsuni(x, y, z1, z2, lab):
-    """bessellsuni.
+    """Bessellsuni.
 
     # Uniform strength
     implicit none
@@ -402,7 +402,7 @@ def bessellsuni(x, y, z1, z2, lab):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsuniv(x, y, z1, z2, lab, rzero):
-    """bessellsuniv.
+    """Bessellsuniv.
 
     # Uniform strength
     implicit none
@@ -513,7 +513,7 @@ def find_d1d2(z1, z2, zc, R):
 
 @numba.njit(nogil=True, cache=True)
 def bessells(x, y, z1, z2, lab, order, d1in, d2in):
-    """bessells.
+    """Bessells.
 
     implicit none
     integer, intent(in) :: order
@@ -783,7 +783,7 @@ def bessellsqxqyv2(x, y, z1, z2, lab, order, R):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsqxqy(x, y, z1, z2, lab, order, d1in, d2in):
-    """bessellsqxqy.
+    """Bessellsqxqy.
 
     implicit none
     integer, intent(in) :: order
@@ -1240,7 +1240,7 @@ def besselldv2(x, y, z1, z2, lab, order, R):
 
 @numba.njit(nogil=True, cache=True)
 def besselld(x, y, z1, z2, lab, order, d1in, d2in):
-    """besselld.
+    """Besselld.
 
     implicit none
     integer, intent(in) :: order
@@ -1495,7 +1495,7 @@ def besselldqxqyv2(x, y, z1, z2, lab, order, R):
 
 @numba.njit(nogil=True, cache=True)
 def besselldqxqy(x, y, z1, z2, lab, order, d1in, d2in):
-    """besselldqxqy.
+    """Besselldqxqy.
 
     implicit none
     integer, intent(in) :: order
@@ -1816,7 +1816,7 @@ def besselld_gauss_ho_qxqy(x, y, z1, z2, lab, order):
 
 @numba.njit(nogil=True, cache=True)
 def besselldpart(x, y, z1, z2, lab, order, d1, d2):
-    """besselldpart.
+    """Besselldpart.
 
     implicit none
     integer, intent(in) :: order
diff --git a/ttim/besselnumba_total.py b/ttim/besselnumba_total.py
index bee0402..9b4c970 100644
--- a/ttim/besselnumba_total.py
+++ b/ttim/besselnumba_total.py
@@ -533,7 +533,7 @@ def lapls_int_ho(x, y, z1, z2, order):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsreal(x, y, x1, y1, x2, y2, lab):
-    """bessellsreal.
+    """Bessellsreal.
 
     implicit none
     real(kind=8), intent(in) :: x,y,x1,y1,x2,y2,lab
@@ -596,7 +596,7 @@ def bessellsreal(x, y, x1, y1, x2, y2, lab):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsrealho(x, y, x1, y1, x2, y2, lab, order):
-    """bessellsrealho.
+    """Bessellsrealho.
 
     implicit none
     real(kind=8), intent(in) :: x,y,x1,y1,x2,y2,lab
@@ -1183,7 +1183,7 @@ def bessells_gauss_ho_qxqy_d1d2(x, y, z1, z2, lab, order, d1, d2):
 
 @numba.njit(nogil=True, cache=True)
 def bessells(x, y, z1, z2, lab, order, d1in, d2in):
-    """bessells.
+    """Bessells.
 
     implicit none
     integer, intent(in) :: order
@@ -1229,7 +1229,7 @@ def bessells(x, y, z1, z2, lab, order, d1in, d2in):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsv(x, y, z1, z2, lab, order, R, nlab):
-    """bessellsv.
+    """Bessellsv.
 
     implicit none
     integer, intent(in) :: order
@@ -1278,7 +1278,7 @@ def bessellsv2(x, y, z1, z2, lab, order, R, nlab):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsqxqy(x, y, z1, z2, lab, order, d1in, d2in):
-    """bessellsqxqy.
+    """Bessellsqxqy.
 
     implicit none
     integer, intent(in) :: order
@@ -1328,7 +1328,7 @@ def bessellsqxqy(x, y, z1, z2, lab, order, d1in, d2in):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsqxqyv(x, y, z1, z2, lab, order, R, nlab):
-    """bessellsqxqyv.
+    """Bessellsqxqyv.
 
     implicit none
     integer, intent(in) :: order
@@ -1382,7 +1382,7 @@ def bessellsqxqyv2(x, y, z1, z2, lab, order, R, nlab):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsuni(x, y, z1, z2, lab):
-    """bessellsuni.
+    """Bessellsuni.
 
     # Uniform strength
     implicit none
@@ -1420,7 +1420,7 @@ def bessellsuni(x, y, z1, z2, lab):
 
 @numba.njit(nogil=True, cache=True)
 def bessellsuniv(x, y, z1, z2, lab, nlab):
-    """bessellsuniv.
+    """Bessellsuniv.
 
     # Uniform strength
     implicit none
@@ -1777,7 +1777,7 @@ def besselld_gauss_ho_d1d2(x, y, z1, z2, lab, order, d1, d2):
 
 @numba.njit(nogil=True, cache=True)
 def besselld(x, y, z1, z2, lab, order, d1in, d2in):
-    """besselld.
+    """Besselld.
 
     implicit none
     integer, intent(in) :: order
@@ -1824,7 +1824,7 @@ def besselld(x, y, z1, z2, lab, order, d1in, d2in):
 
 @numba.njit(nogil=True, cache=True)
 def besselldv(x, y, z1, z2, lab, order, R, nlab):
-    """besselldv.
+    """Besselldv.
 
     implicit none
     integer, intent(in) :: order
@@ -1875,7 +1875,7 @@ def besselldv2(x, y, z1, z2, lab, order, R, nlab):
 
 # @numba.njit(nogil=True, cache=True)
 def besselldpart(x, y, z1, z2, lab, order, d1, d2):
-    """besselldpart.
+    """Besselldpart.
 
     implicit none
     integer, intent(in) :: order
@@ -2256,7 +2256,7 @@ def besselld_gauss_ho_qxqy_d1d2(x, y, z1, z2, lab, order, d1, d2):
 
 # @numba.njit(nogil=True, cache=True)
 def besselldqxqy(x, y, z1, z2, lab, order, d1in, d2in):
-    """besselldqxqy.
+    """Besselldqxqy.
 
     implicit none
     integer, intent(in) :: order
@@ -2308,7 +2308,7 @@ def besselldqxqy(x, y, z1, z2, lab, order, d1in, d2in):
 
 # @numba.njit(nogil=True, cache=True)
 def besselldqxqyv(x, y, z1, z2, lab, order, R, nlab):
-    """besselldqxqyv.
+    """Besselldqxqyv.
 
     implicit none
     integer, intent(in) :: order
diff --git a/ttim/circinhom.py b/ttim/circinhom.py
index 1c01f44..c5b3d27 100644
--- a/ttim/circinhom.py
+++ b/ttim/circinhom.py
@@ -198,9 +198,9 @@ def initialize(self):
         self.aqin = self.model.aq.findAquiferData(
             self.x0 + (1 - 1e-8) * self.R, self.y0
         )
-        assert (
-            self.aqin.R == self.R
-        ), "Radius of CircInhom and CircInhomData must be equal"
+        assert self.aqin.R == self.R, (
+            "Radius of CircInhom and CircInhomData must be equal"
+        )
         self.aqout = self.model.aq.find_aquifer_data(
             self.x0 + (1 + 1e-8) * self.R, self.y0
         )
diff --git a/ttim/fit.py b/ttim/fit.py
index b80c97e..dff9352 100755
--- a/ttim/fit.py
+++ b/ttim/fit.py
@@ -168,9 +168,9 @@ def set_parameter_by_reference(
         """
         assert isinstance(name, str), "Error: name must be string"
         if parameter is not None:
-            assert isinstance(
-                parameter, np.ndarray
-            ), "Error: parameter needs to be numpy array"
+            assert isinstance(parameter, np.ndarray), (
+                "Error: parameter needs to be numpy array"
+            )
             p = parameter
         self.parameters.loc[name] = {
             "optimal": initial,
diff --git a/ttim/plots.py b/ttim/plots.py
index a555472..b5a5059 100644
--- a/ttim/plots.py
+++ b/ttim/plots.py
@@ -153,8 +153,7 @@ def xsection(
             if params and self._ml.aq.ltype[i] == "a":
                 if aqi == 0 and self._ml.aq.phreatictop:
                     paramtxt = (
-                        f"$k_h$ = {self._ml.aq.kaq[aqi]}, "
-                        f"$S$ = {self._ml.aq.Saq[aqi]}"
+                        f"$k_h$ = {self._ml.aq.kaq[aqi]}, $S$ = {self._ml.aq.Saq[aqi]}"
                     )
                 else:
                     paramtxt = (