From 642b6da7ed320c32849ec02e54964336ac1e361d Mon Sep 17 00:00:00 2001
From: Joshua Larsen <jlarsen@usgs.gov>
Date: Fri, 4 Feb 2022 09:12:40 -0800
Subject: [PATCH 1/7] update(t074_test_geospatial_util.py): update tests for
 pyshp 2.2.0

* fix typo in AG package
* update Raster for slightly faster grid resampling
---
 autotest/t074_test_geospatial_util.py | 54 +++++++++++++--------------
 flopy/modflow/mfag.py                 |  8 ++--
 flopy/utils/rasters.py                | 24 +++++++-----
 3 files changed, 44 insertions(+), 42 deletions(-)

diff --git a/autotest/t074_test_geospatial_util.py b/autotest/t074_test_geospatial_util.py
index 90492d110f..8956e7fe42 100644
--- a/autotest/t074_test_geospatial_util.py
+++ b/autotest/t074_test_geospatial_util.py
@@ -3,12 +3,12 @@
     "coordinates": (
         (
             (-121.389308, 38.560816),
-            (-121.363391, 38.568835),
-            (-121.358641, 38.565972),
-            (-121.359327, 38.562767),
-            (-121.369932, 38.560575),
-            (-121.370609, 38.557232),
             (-121.385435, 38.555018),
+            (-121.370609, 38.557232),
+            (-121.369932, 38.560575),
+            (-121.359327, 38.562767),
+            (-121.358641, 38.565972),
+            (-121.363391, 38.568835),
             (-121.389308, 38.560816),
         ),
     ),
@@ -19,19 +19,19 @@
     "coordinates": (
         (
             (-121.383097, 38.565764),
-            (-121.342866, 38.579086),
-            (-121.342739, 38.578995),
-            (-121.323309, 38.578953),
-            (-121.358295, 38.561163),
-            (-121.379047, 38.559053),
             (-121.382318, 38.562934),
+            (-121.379047, 38.559053),
+            (-121.358295, 38.561163),
+            (-121.323309, 38.578953),
+            (-121.342739, 38.578995),
+            (-121.342866, 38.579086),
             (-121.383097, 38.565764),
         ),
         (
             (-121.367281, 38.567214),
-            (-121.362633, 38.562622),
-            (-121.345857, 38.570301),
             (-121.352168, 38.572258),
+            (-121.345857, 38.570301),
+            (-121.362633, 38.562622),
             (-121.367281, 38.567214),
         ),
     ),
@@ -43,17 +43,17 @@
         [
             (
                 (-121.433775, 38.544254),
-                (-121.424263, 38.547474),
                 (-121.422917, 38.540376),
+                (-121.424263, 38.547474),
                 (-121.433775, 38.544254),
             )
         ],
         [
             (
                 (-121.456113, 38.552220),
-                (-121.440092, 38.548303),
-                (-121.440053, 38.537820),
                 (-121.459991, 38.541350),
+                (-121.440053, 38.537820),
+                (-121.440092, 38.548303),
                 (-121.456113, 38.552220),
             )
         ],
@@ -113,16 +113,16 @@
 def test_import_geospatial_utils():
 
     from flopy.utils.geospatial_utils import (
-        GeoSpatialCollection,
         GeoSpatialUtil,
+        GeoSpatialCollection,
     )
 
     return
 
 
 def test_polygon():
-    from flopy.utils.geometry import Polygon, Shape
     from flopy.utils.geospatial_utils import GeoSpatialUtil
+    from flopy.utils.geometry import Shape, Polygon
 
     poly = Shape.from_geojson(polygon)
     gi1 = poly.__geo_interface__
@@ -153,8 +153,8 @@ def test_polygon():
 
 
 def test_polygon_with_hole():
-    from flopy.utils.geometry import Polygon, Shape
     from flopy.utils.geospatial_utils import GeoSpatialUtil
+    from flopy.utils.geometry import Shape, Polygon
 
     poly = Shape.from_geojson(poly_w_hole)
     gi1 = poly.__geo_interface__
@@ -185,8 +185,8 @@ def test_polygon_with_hole():
 
 
 def test_multipolygon():
-    from flopy.utils.geometry import MultiPolygon, Shape
     from flopy.utils.geospatial_utils import GeoSpatialUtil
+    from flopy.utils.geometry import Shape, MultiPolygon
 
     poly = Shape.from_geojson(multipolygon)
     gi1 = poly.__geo_interface__
@@ -219,8 +219,8 @@ def test_multipolygon():
 
 
 def test_point():
-    from flopy.utils.geometry import Point, Shape
     from flopy.utils.geospatial_utils import GeoSpatialUtil
+    from flopy.utils.geometry import Shape, Point
 
     pt = Shape.from_geojson(point)
     gi1 = pt.__geo_interface__
@@ -251,8 +251,8 @@ def test_point():
 
 
 def test_multipoint():
-    from flopy.utils.geometry import MultiPoint, Shape
     from flopy.utils.geospatial_utils import GeoSpatialUtil
+    from flopy.utils.geometry import Shape, MultiPoint
 
     mpt = Shape.from_geojson(multipoint)
     gi1 = mpt.__geo_interface__
@@ -283,8 +283,8 @@ def test_multipoint():
 
 
 def test_linestring():
-    from flopy.utils.geometry import LineString, Shape
     from flopy.utils.geospatial_utils import GeoSpatialUtil
+    from flopy.utils.geometry import Shape, LineString
 
     lstr = Shape.from_geojson(linestring)
     gi1 = lstr.__geo_interface__
@@ -315,8 +315,8 @@ def test_linestring():
 
 
 def test_multilinestring():
-    from flopy.utils.geometry import MultiLineString, Shape
     from flopy.utils.geospatial_utils import GeoSpatialUtil
+    from flopy.utils.geometry import Shape, MultiLineString
 
     mlstr = Shape.from_geojson(multilinestring)
     gi1 = mlstr.__geo_interface__
@@ -349,8 +349,8 @@ def test_multilinestring():
 
 
 def test_polygon_collection():
-    from flopy.utils.geometry import Collection, Shape
     from flopy.utils.geospatial_utils import GeoSpatialCollection
+    from flopy.utils.geometry import Shape, Collection
 
     col = [
         Shape.from_geojson(polygon),
@@ -387,8 +387,8 @@ def test_polygon_collection():
 
 
 def test_point_collection():
-    from flopy.utils.geometry import Collection, Shape
     from flopy.utils.geospatial_utils import GeoSpatialCollection
+    from flopy.utils.geometry import Shape, Collection
 
     col = [Shape.from_geojson(point), Shape.from_geojson(multipoint)]
 
@@ -421,8 +421,8 @@ def test_point_collection():
 
 
 def test_linestring_collection():
-    from flopy.utils.geometry import Collection, Shape
     from flopy.utils.geospatial_utils import GeoSpatialCollection
+    from flopy.utils.geometry import Shape, Collection
 
     col = [Shape.from_geojson(linestring), Shape.from_geojson(multilinestring)]
 
@@ -455,8 +455,8 @@ def test_linestring_collection():
 
 
 def test_mixed_collection():
-    from flopy.utils.geometry import Collection, Shape
     from flopy.utils.geospatial_utils import GeoSpatialCollection
+    from flopy.utils.geometry import Shape, Collection
 
     col = [
         Shape.from_geojson(polygon),
diff --git a/flopy/modflow/mfag.py b/flopy/modflow/mfag.py
index 5eae021880..ebb193bca7 100644
--- a/flopy/modflow/mfag.py
+++ b/flopy/modflow/mfag.py
@@ -9,13 +9,11 @@
 """
 
 import os
-
 import numpy as np
-
-from ..pakbase import Package
 from ..utils.flopy_io import multi_line_strip
-from ..utils.optionblock import OptionBlock
+from ..pakbase import Package
 from ..utils.recarray_utils import create_empty_recarray
+from ..utils.optionblock import OptionBlock
 
 
 class ModflowAg(Package):
@@ -360,7 +358,7 @@ def write_file(self, check=False):
 
             # check if item 12 exists and write item 12 - 14
             if self.segment_list is not None:
-                foo.write("# segment list for irriagation diversions\n")
+                foo.write("# segment list for irrigation diversions\n")
                 foo.write("SEGMENT LIST\n")
                 for iseg in self.segment_list:
                     foo.write(f"{iseg}\n")
diff --git a/flopy/utils/rasters.py b/flopy/utils/rasters.py
index f67b072f00..3971b25cb8 100644
--- a/flopy/utils/rasters.py
+++ b/flopy/utils/rasters.py
@@ -265,7 +265,7 @@ def sample_point(self, *point, band=1):
 
         return value
 
-    def sample_polygon(self, polygon, band, invert=False):
+    def sample_polygon(self, polygon, band, invert=False, **kwargs):
         """
         Method to get an unordered list of raster values that are located
         within a arbitrary polygon
@@ -309,7 +309,7 @@ def sample_polygon(self, polygon, band, invert=False):
                     arr_dict[b] = t
 
         else:
-            mask = self._intersection(polygon, invert)
+            mask = self._intersection(polygon, invert, **kwargs)
 
             arr_dict = {}
             for b, arr in self.__arr_dict.items():
@@ -455,7 +455,7 @@ def resample_to_grid(
             else:
                 for node in range(ncpl):
                     verts = modelgrid.get_cell_vertices(node)
-                    rstr_data = self.sample_polygon(verts, band).astype(float)
+                    rstr_data = self.sample_polygon(verts, band, convert=False).astype(float)
                     msk = np.in1d(rstr_data, self.nodatavals)
                     rstr_data[msk] = np.nan
 
@@ -540,7 +540,7 @@ def __threaded_resampling(
         """
         container.acquire()
         verts = modelgrid.get_cell_vertices(node)
-        rstr_data = self.sample_polygon(verts, band).astype(float)
+        rstr_data = self.sample_polygon(verts, band, convert=False).astype(float)
         msk = np.in1d(rstr_data, self.nodatavals)
         rstr_data[msk] = np.nan
 
@@ -715,7 +715,7 @@ def _sample_rio_dataset(self, polygon, invert):
 
         return arr_dict, rstr_crp_meta
 
-    def _intersection(self, polygon, invert):
+    def _intersection(self, polygon, invert, **kwargs):
         """
         Internal method to create an intersection mask, used for cropping
         arrays and sampling arrays.
@@ -740,14 +740,18 @@ def _intersection(self, polygon, invert):
             mask : np.ndarray (dtype = bool)
 
         """
-        from .geospatial_utils import GeoSpatialUtil
+        # the convert kwarg is to speed up the resample_to_grid method
+        #  which already provides the proper datatype to _intersect()
+        convert = kwargs.pop("convert", True)
+        if convert:
+            from .geospatial_utils import GeoSpatialUtil
 
-        if isinstance(polygon, (list, tuple, np.ndarray)):
-            polygon = [polygon]
+            if isinstance(polygon, (list, tuple, np.ndarray)):
+                polygon = [polygon]
 
-        geom = GeoSpatialUtil(polygon, shapetype="Polygon")
+            geom = GeoSpatialUtil(polygon, shapetype="Polygon")
 
-        polygon = geom.points[0]
+            polygon = geom.points[0]
 
         # step 2: create a grid of centoids
         xc = self.xcenters

From 050b4b65ed7ed9e5ee487442ce02b45d8ff3c48a Mon Sep 17 00:00:00 2001
From: Joshua Larsen <jlarsen@usgs.gov>
Date: Fri, 4 Feb 2022 11:18:23 -0800
Subject: [PATCH 2/7] update(geometry.py): added __reversed__ to geometery
 object as patch for pyshp 2.2.0 vs 2.1.0 differences

* updated t074_test_geospatial_util.py
* linted repository
---
 autotest/t074_test_geospatial_util.py | 59 ++++++++++++++++++++-------
 flopy/discretization/modeltime.py     |  2 +-
 flopy/modflow/mfdis.py                |  2 +-
 flopy/modflow/mfmlt.py                |  2 +-
 flopy/modflow/mfsfr2.py               |  2 +-
 flopy/plot/crosssection.py            | 10 ++---
 flopy/plot/map.py                     | 10 ++---
 flopy/utils/binaryfile.py             |  2 +-
 flopy/utils/flopy_io.py               |  2 +-
 flopy/utils/geometry.py               | 25 ++++++++++++
 flopy/utils/rasters.py                |  8 +++-
 11 files changed, 91 insertions(+), 33 deletions(-)

diff --git a/autotest/t074_test_geospatial_util.py b/autotest/t074_test_geospatial_util.py
index 8956e7fe42..0da761bcfe 100644
--- a/autotest/t074_test_geospatial_util.py
+++ b/autotest/t074_test_geospatial_util.py
@@ -149,7 +149,11 @@ def test_polygon():
         is_equal = gi1 == gi2
 
         if not is_equal:
-            raise AssertionError("GeoSpatialUtil polygon conversion error")
+            # pyshp < 2.2.0 sorts coordinates in opposite direction
+            gi2['coordinates'] = (gi2['coordinates'][0][::-1],)
+            is_equal = gi1 == gi2
+            if not is_equal:
+                raise AssertionError("GeoSpatialUtil polygon conversion error")
 
 
 def test_polygon_with_hole():
@@ -181,7 +185,12 @@ def test_polygon_with_hole():
         is_equal = gi1 == gi2
 
         if not is_equal:
-            raise AssertionError("GeoSpatialUtil polygon conversion error")
+            # pyshp < 2.2.0 sorts coordinates in opposite direction
+            t = reversed(t)
+            gi2 = t.__geo_interface__
+            is_equal = gi1 == gi2
+            if not is_equal:
+                raise AssertionError("GeoSpatialUtil polygon conversion error")
 
 
 def test_multipolygon():
@@ -213,9 +222,14 @@ def test_multipolygon():
         is_equal = gi1 == gi2
 
         if not is_equal:
-            raise AssertionError(
-                "GeoSpatialUtil multipolygon conversion error"
-            )
+            # pyshp < 2.2.0 sorts coordinates in opposite direction
+            t = reversed(t)
+            gi2 = t.__geo_interface__
+            is_equal = gi1 == gi2
+            if not is_equal:
+                raise AssertionError(
+                    "GeoSpatialUtil multipolygon conversion error"
+                )
 
 
 def test_point():
@@ -375,15 +389,21 @@ def test_polygon_collection():
             continue
 
         gc2 = GeoSpatialCollection(col, shapetype)
-        gi2 = [i.flopy_geometry.__geo_interface__ for i in gc2]
 
-        for ix, gi in enumerate(gi2):
-            is_equal = gi == gi1[ix]
+        for ix, gi in enumerate(gc2):
+            t = gi.flopy_geometry
+            gi2 = t.__geo_interface__
+            is_equal = gi2 == gi1[ix]
 
             if not is_equal:
-                raise AssertionError(
-                    "GeoSpatialCollection Polygon conversion error"
-                )
+                # pyshp < 2.2.0 sorts coordinates in opposite direction
+                t = reversed(t)
+                gi2 = t.__geo_interface__
+                is_equal = gi2 == gi1[ix]
+                if not is_equal:
+                    raise AssertionError(
+                        "GeoSpatialCollection Polygon conversion error"
+                    )
 
 
 def test_point_collection():
@@ -485,13 +505,22 @@ def test_mixed_collection():
             continue
 
         gc2 = GeoSpatialCollection(col, shapetype)
-        gi2 = [i.flopy_geometry.__geo_interface__ for i in gc2]
 
-        for ix, gi in enumerate(gi2):
-            is_equal = gi == gi1[ix]
+        for ix, gi in enumerate(gc2):
+            t = gi.flopy_geometry
+            gi2 = t.__geo_interface__
+
+            is_equal = gi2 == gi1[ix]
 
             if not is_equal:
-                raise AssertionError("GeoSpatialCollection conversion error")
+                t = reversed(t)
+                gi2 = t.__geo_interface__
+                is_equal = gi2 == gi1[ix]
+
+                if not is_equal:
+                    raise AssertionError(
+                        "GeoSpatialCollection conversion error"
+                    )
 
 
 if __name__ == "__main__":
diff --git a/flopy/discretization/modeltime.py b/flopy/discretization/modeltime.py
index da7014931e..83994ea3ed 100644
--- a/flopy/discretization/modeltime.py
+++ b/flopy/discretization/modeltime.py
@@ -65,7 +65,7 @@ def totim(self):
             for stp in range(nstp):
                 if stp == 0:
                     if tsmult != 1.0:
-                        dt = perlen * (tsmult - 1) / ((tsmult**nstp) - 1)
+                        dt = perlen * (tsmult - 1) / ((tsmult ** nstp) - 1)
                     else:
                         dt = perlen / nstp
                 else:
diff --git a/flopy/modflow/mfdis.py b/flopy/modflow/mfdis.py
index 879337bc00..836d0fc313 100644
--- a/flopy/modflow/mfdis.py
+++ b/flopy/modflow/mfdis.py
@@ -300,7 +300,7 @@ def get_totim(self, use_cached=False):
                 p = float(nstp[kper])
                 dt = perlen[kper]
                 if m > 1:
-                    dt *= (m - 1.0) / (m**p - 1.0)
+                    dt *= (m - 1.0) / (m ** p - 1.0)
                 else:
                     dt = dt / p
                 for kstp in range(nstp[kper]):
diff --git a/flopy/modflow/mfmlt.py b/flopy/modflow/mfmlt.py
index f827b97224..eaea8b6d80 100644
--- a/flopy/modflow/mfmlt.py
+++ b/flopy/modflow/mfmlt.py
@@ -243,7 +243,7 @@ def mult_function(mult_dict, line):
             elif op == "/":
                 multarray = multarray / atemp
             elif op == "^":
-                multarray = multarray**atemp
+                multarray = multarray ** atemp
             else:
                 raise Exception(f"Invalid MULT operation {op}")
         return multarray
diff --git a/flopy/modflow/mfsfr2.py b/flopy/modflow/mfsfr2.py
index 9da7cfb6ab..50a613057d 100644
--- a/flopy/modflow/mfsfr2.py
+++ b/flopy/modflow/mfsfr2.py
@@ -2499,7 +2499,7 @@ def routing(self):
 
             dx = delr[rd.j]
             dy = delc[rd.i]
-            hyp = np.sqrt(dx**2 + dy**2)
+            hyp = np.sqrt(dx ** 2 + dy ** 2)
 
             # breaks are when the connection distance is greater than
             # max node with * a tolerance
diff --git a/flopy/plot/crosssection.py b/flopy/plot/crosssection.py
index 4049684944..e84950249c 100644
--- a/flopy/plot/crosssection.py
+++ b/flopy/plot/crosssection.py
@@ -537,13 +537,13 @@ def contour_array(self, a, masked_values=None, head=None, **kwargs):
             kwargs["levels"] = levels
 
         # workaround for tri-contour nan issue
-        plotarray[np.isnan(plotarray)] = -(2**31)
+        plotarray[np.isnan(plotarray)] = -(2 ** 31)
         if masked_values is None:
-            masked_values = [-(2**31)]
+            masked_values = [-(2 ** 31)]
         else:
             masked_values = list(masked_values)
-            if -(2**31) not in masked_values:
-                masked_values.append(-(2**31))
+            if -(2 ** 31) not in masked_values:
+                masked_values.append(-(2 ** 31))
 
         ismasked = None
         if masked_values is not None:
@@ -965,7 +965,7 @@ def plot_vector(
 
         # normalize
         if normalize:
-            vmag = np.sqrt(u**2.0 + v**2.0)
+            vmag = np.sqrt(u ** 2.0 + v ** 2.0)
             idx = vmag > 0.0
             u[idx] /= vmag[idx]
             v[idx] /= vmag[idx]
diff --git a/flopy/plot/map.py b/flopy/plot/map.py
index fafaf6eb4c..e0ef87bfe3 100644
--- a/flopy/plot/map.py
+++ b/flopy/plot/map.py
@@ -182,13 +182,13 @@ def contour_array(self, a, masked_values=None, **kwargs):
 
         # workaround for tri-contour nan issue
         # use -2**31 to allow for 32 bit int arrays
-        plotarray[np.isnan(plotarray)] = -(2**31)
+        plotarray[np.isnan(plotarray)] = -(2 ** 31)
         if masked_values is None:
-            masked_values = [-(2**31)]
+            masked_values = [-(2 ** 31)]
         else:
             masked_values = list(masked_values)
-            if -(2**31) not in masked_values:
-                masked_values.append(-(2**31))
+            if -(2 ** 31) not in masked_values:
+                masked_values.append(-(2 ** 31))
 
         ismasked = None
         if masked_values is not None:
@@ -684,7 +684,7 @@ def plot_vector(
 
         # normalize
         if normalize:
-            vmag = np.sqrt(u**2.0 + v**2.0)
+            vmag = np.sqrt(u ** 2.0 + v ** 2.0)
             idx = vmag > 0.0
             u[idx] /= vmag[idx]
             v[idx] /= vmag[idx]
diff --git a/flopy/utils/binaryfile.py b/flopy/utils/binaryfile.py
index 75c43fe13d..e9be3346c8 100755
--- a/flopy/utils/binaryfile.py
+++ b/flopy/utils/binaryfile.py
@@ -781,7 +781,7 @@ def _totim_from_kstpkper(self, kstpkper):
         if tsmult == 1:
             dt1 = this_perlen / float(nstp)
         else:
-            dt1 = this_perlen * (tsmult - 1.0) / ((tsmult**nstp) - 1.0)
+            dt1 = this_perlen * (tsmult - 1.0) / ((tsmult ** nstp) - 1.0)
         kstp_len = [dt1]
         for i in range(kstp + 1):
             kstp_len.append(kstp_len[-1] * tsmult)
diff --git a/flopy/utils/flopy_io.py b/flopy/utils/flopy_io.py
index 6061a223bc..084c9146df 100755
--- a/flopy/utils/flopy_io.py
+++ b/flopy/utils/flopy_io.py
@@ -184,7 +184,7 @@ def write_fixed_var(v, length=10, ipos=None, free=False, comment=None):
             width = ipos[n]
             if isinstance(v[n], (float, np.float32, np.float64)):
                 decimal = width - 6
-                vmin, vmax = 10**-decimal, 10**decimal
+                vmin, vmax = 10 ** -decimal, 10 ** decimal
                 if abs(v[n]) < vmin or abs(v[n]) > vmax:
                     ctype = "g"  # default precision is 6 if not specified
                 else:
diff --git a/flopy/utils/geometry.py b/flopy/utils/geometry.py
index e095a51897..490fbd82d6 100644
--- a/flopy/utils/geometry.py
+++ b/flopy/utils/geometry.py
@@ -199,6 +199,11 @@ def bounds(self):
 
         return xmin, ymin, xmax, ymax
 
+    def __reversed__(self):
+        for shp in self:
+            reversed(shp)
+        return self
+
     def plot(self, ax=None, **kwargs):
         """
         Plotting method for collection
@@ -409,6 +414,19 @@ def pyshp_parts(self):
     def patch(self):
         return self.get_patch()
 
+    def __reversed__(self):
+        # method to reverse the sorting on polygon points, patch for
+        # differences in pyshp 2.2.0 vs pyshp < 2.2.0
+        if self.exterior:
+            self.exterior = self.exterior[::-1]
+        if self.interiors:
+            interiors = []
+            for i in self.interiors:
+                interiors.append(i[::-1])
+            self.interiors = interiors
+
+        return self
+
     def get_patch(self, **kwargs):
         descartes = import_optional_dependency("descartes")
         from descartes import PolygonPatch
@@ -515,6 +533,10 @@ def bounds(self):
     def pyshp_parts(self):
         return [self.coords]
 
+    def __reversed__(self):
+        self.coords = self.coords[::-1]
+        return self
+
     def plot(self, ax=None, **kwargs):
         import matplotlib.pyplot as plt
 
@@ -611,6 +633,9 @@ def bounds(self):
     def pyshp_parts(self):
         return self.coords
 
+    def __reversed__(self):
+        return self
+
     def plot(self, ax=None, **kwargs):
         import matplotlib.pyplot as plt
 
diff --git a/flopy/utils/rasters.py b/flopy/utils/rasters.py
index 3971b25cb8..df56397d45 100644
--- a/flopy/utils/rasters.py
+++ b/flopy/utils/rasters.py
@@ -455,7 +455,9 @@ def resample_to_grid(
             else:
                 for node in range(ncpl):
                     verts = modelgrid.get_cell_vertices(node)
-                    rstr_data = self.sample_polygon(verts, band, convert=False).astype(float)
+                    rstr_data = self.sample_polygon(
+                        verts, band, convert=False
+                    ).astype(float)
                     msk = np.in1d(rstr_data, self.nodatavals)
                     rstr_data[msk] = np.nan
 
@@ -540,7 +542,9 @@ def __threaded_resampling(
         """
         container.acquire()
         verts = modelgrid.get_cell_vertices(node)
-        rstr_data = self.sample_polygon(verts, band, convert=False).astype(float)
+        rstr_data = self.sample_polygon(verts, band, convert=False).astype(
+            float
+        )
         msk = np.in1d(rstr_data, self.nodatavals)
         rstr_data[msk] = np.nan
 

From c4168bbb3eaa46aa21d89af0f1f0d1a533aef259 Mon Sep 17 00:00:00 2001
From: Joshua Larsen <jlarsen@usgs.gov>
Date: Fri, 4 Feb 2022 11:31:44 -0800
Subject: [PATCH 3/7] isort mfag.py

---
 flopy/modflow/mfag.py | 6 ++++--
 1 file changed, 4 insertions(+), 2 deletions(-)

diff --git a/flopy/modflow/mfag.py b/flopy/modflow/mfag.py
index ebb193bca7..786e2ec37c 100644
--- a/flopy/modflow/mfag.py
+++ b/flopy/modflow/mfag.py
@@ -9,11 +9,13 @@
 """
 
 import os
+
 import numpy as np
-from ..utils.flopy_io import multi_line_strip
+
 from ..pakbase import Package
-from ..utils.recarray_utils import create_empty_recarray
+from ..utils.flopy_io import multi_line_strip
 from ..utils.optionblock import OptionBlock
+from ..utils.recarray_utils import create_empty_recarray
 
 
 class ModflowAg(Package):

From f3ec13a471c3482eae0af47c7e5c0c2513025124 Mon Sep 17 00:00:00 2001
From: Joshua Larsen <jlarsen@usgs.gov>
Date: Fri, 4 Feb 2022 11:34:28 -0800
Subject: [PATCH 4/7] update black and perform linting

---
 flopy/discretization/modeltime.py |  2 +-
 flopy/modflow/mfdis.py            |  2 +-
 flopy/modflow/mfmlt.py            |  2 +-
 flopy/modflow/mfsfr2.py           |  2 +-
 flopy/plot/crosssection.py        | 10 +++++-----
 flopy/plot/map.py                 | 10 +++++-----
 flopy/utils/binaryfile.py         |  2 +-
 flopy/utils/flopy_io.py           |  2 +-
 8 files changed, 16 insertions(+), 16 deletions(-)

diff --git a/flopy/discretization/modeltime.py b/flopy/discretization/modeltime.py
index 83994ea3ed..da7014931e 100644
--- a/flopy/discretization/modeltime.py
+++ b/flopy/discretization/modeltime.py
@@ -65,7 +65,7 @@ def totim(self):
             for stp in range(nstp):
                 if stp == 0:
                     if tsmult != 1.0:
-                        dt = perlen * (tsmult - 1) / ((tsmult ** nstp) - 1)
+                        dt = perlen * (tsmult - 1) / ((tsmult**nstp) - 1)
                     else:
                         dt = perlen / nstp
                 else:
diff --git a/flopy/modflow/mfdis.py b/flopy/modflow/mfdis.py
index 836d0fc313..879337bc00 100644
--- a/flopy/modflow/mfdis.py
+++ b/flopy/modflow/mfdis.py
@@ -300,7 +300,7 @@ def get_totim(self, use_cached=False):
                 p = float(nstp[kper])
                 dt = perlen[kper]
                 if m > 1:
-                    dt *= (m - 1.0) / (m ** p - 1.0)
+                    dt *= (m - 1.0) / (m**p - 1.0)
                 else:
                     dt = dt / p
                 for kstp in range(nstp[kper]):
diff --git a/flopy/modflow/mfmlt.py b/flopy/modflow/mfmlt.py
index eaea8b6d80..f827b97224 100644
--- a/flopy/modflow/mfmlt.py
+++ b/flopy/modflow/mfmlt.py
@@ -243,7 +243,7 @@ def mult_function(mult_dict, line):
             elif op == "/":
                 multarray = multarray / atemp
             elif op == "^":
-                multarray = multarray ** atemp
+                multarray = multarray**atemp
             else:
                 raise Exception(f"Invalid MULT operation {op}")
         return multarray
diff --git a/flopy/modflow/mfsfr2.py b/flopy/modflow/mfsfr2.py
index 50a613057d..9da7cfb6ab 100644
--- a/flopy/modflow/mfsfr2.py
+++ b/flopy/modflow/mfsfr2.py
@@ -2499,7 +2499,7 @@ def routing(self):
 
             dx = delr[rd.j]
             dy = delc[rd.i]
-            hyp = np.sqrt(dx ** 2 + dy ** 2)
+            hyp = np.sqrt(dx**2 + dy**2)
 
             # breaks are when the connection distance is greater than
             # max node with * a tolerance
diff --git a/flopy/plot/crosssection.py b/flopy/plot/crosssection.py
index e84950249c..4049684944 100644
--- a/flopy/plot/crosssection.py
+++ b/flopy/plot/crosssection.py
@@ -537,13 +537,13 @@ def contour_array(self, a, masked_values=None, head=None, **kwargs):
             kwargs["levels"] = levels
 
         # workaround for tri-contour nan issue
-        plotarray[np.isnan(plotarray)] = -(2 ** 31)
+        plotarray[np.isnan(plotarray)] = -(2**31)
         if masked_values is None:
-            masked_values = [-(2 ** 31)]
+            masked_values = [-(2**31)]
         else:
             masked_values = list(masked_values)
-            if -(2 ** 31) not in masked_values:
-                masked_values.append(-(2 ** 31))
+            if -(2**31) not in masked_values:
+                masked_values.append(-(2**31))
 
         ismasked = None
         if masked_values is not None:
@@ -965,7 +965,7 @@ def plot_vector(
 
         # normalize
         if normalize:
-            vmag = np.sqrt(u ** 2.0 + v ** 2.0)
+            vmag = np.sqrt(u**2.0 + v**2.0)
             idx = vmag > 0.0
             u[idx] /= vmag[idx]
             v[idx] /= vmag[idx]
diff --git a/flopy/plot/map.py b/flopy/plot/map.py
index e0ef87bfe3..fafaf6eb4c 100644
--- a/flopy/plot/map.py
+++ b/flopy/plot/map.py
@@ -182,13 +182,13 @@ def contour_array(self, a, masked_values=None, **kwargs):
 
         # workaround for tri-contour nan issue
         # use -2**31 to allow for 32 bit int arrays
-        plotarray[np.isnan(plotarray)] = -(2 ** 31)
+        plotarray[np.isnan(plotarray)] = -(2**31)
         if masked_values is None:
-            masked_values = [-(2 ** 31)]
+            masked_values = [-(2**31)]
         else:
             masked_values = list(masked_values)
-            if -(2 ** 31) not in masked_values:
-                masked_values.append(-(2 ** 31))
+            if -(2**31) not in masked_values:
+                masked_values.append(-(2**31))
 
         ismasked = None
         if masked_values is not None:
@@ -684,7 +684,7 @@ def plot_vector(
 
         # normalize
         if normalize:
-            vmag = np.sqrt(u ** 2.0 + v ** 2.0)
+            vmag = np.sqrt(u**2.0 + v**2.0)
             idx = vmag > 0.0
             u[idx] /= vmag[idx]
             v[idx] /= vmag[idx]
diff --git a/flopy/utils/binaryfile.py b/flopy/utils/binaryfile.py
index e9be3346c8..75c43fe13d 100755
--- a/flopy/utils/binaryfile.py
+++ b/flopy/utils/binaryfile.py
@@ -781,7 +781,7 @@ def _totim_from_kstpkper(self, kstpkper):
         if tsmult == 1:
             dt1 = this_perlen / float(nstp)
         else:
-            dt1 = this_perlen * (tsmult - 1.0) / ((tsmult ** nstp) - 1.0)
+            dt1 = this_perlen * (tsmult - 1.0) / ((tsmult**nstp) - 1.0)
         kstp_len = [dt1]
         for i in range(kstp + 1):
             kstp_len.append(kstp_len[-1] * tsmult)
diff --git a/flopy/utils/flopy_io.py b/flopy/utils/flopy_io.py
index 084c9146df..6061a223bc 100755
--- a/flopy/utils/flopy_io.py
+++ b/flopy/utils/flopy_io.py
@@ -184,7 +184,7 @@ def write_fixed_var(v, length=10, ipos=None, free=False, comment=None):
             width = ipos[n]
             if isinstance(v[n], (float, np.float32, np.float64)):
                 decimal = width - 6
-                vmin, vmax = 10 ** -decimal, 10 ** decimal
+                vmin, vmax = 10**-decimal, 10**decimal
                 if abs(v[n]) < vmin or abs(v[n]) > vmax:
                     ctype = "g"  # default precision is 6 if not specified
                 else:

From a744058109403b793456f2583eef487de4987092 Mon Sep 17 00:00:00 2001
From: Joshua Larsen <jlarsen@usgs.gov>
Date: Fri, 4 Feb 2022 12:15:11 -0800
Subject: [PATCH 5/7] update(requirements): merge upstream pyshp requirements
 with fork

---
 binder/environment.yml           | 2 +-
 etc/environment.yml              | 2 +-
 etc/requirements.full.pip.txt    | 2 +-
 etc/requirements.pip.txt         | 2 +-
 etc/requirements.windows.pip.txt | 2 +-
 requirements.txt                 | 2 +-
 6 files changed, 6 insertions(+), 6 deletions(-)

diff --git a/binder/environment.yml b/binder/environment.yml
index 1889ff19d2..0d9e1f248a 100644
--- a/binder/environment.yml
+++ b/binder/environment.yml
@@ -9,7 +9,7 @@ dependencies:
   - conda-forge::scipy
   - conda-forge::pandas
   - conda-forge::netcdf4
-  - conda-forge::pyshp
+  - conda-forge::pyshp<2.2.0
   - conda-forge::rasterio
   - conda-forge::fiona
   - conda-forge::descartes
diff --git a/etc/environment.yml b/etc/environment.yml
index 2f8f314353..ca9de03292 100644
--- a/etc/environment.yml
+++ b/etc/environment.yml
@@ -24,7 +24,7 @@ dependencies:
   - scipy
   - pandas
   - netcdf4
-  - pyshp>=2.0.0
+  - pyshp<2.2.0
   - rasterio
   - fiona
   - descartes
diff --git a/etc/requirements.full.pip.txt b/etc/requirements.full.pip.txt
index 39413e8dce..5f0e12a4d7 100644
--- a/etc/requirements.full.pip.txt
+++ b/etc/requirements.full.pip.txt
@@ -18,7 +18,7 @@ affine
 scipy
 pandas
 netcdf4
-pyshp>=2.0.0
+pyshp<2.2.0
 rasterio
 fiona
 descartes
diff --git a/etc/requirements.pip.txt b/etc/requirements.pip.txt
index 31feb1653d..bf0b434c73 100644
--- a/etc/requirements.pip.txt
+++ b/etc/requirements.pip.txt
@@ -7,6 +7,6 @@ requests
 appdirs
 numpy
 matplotlib
-pytest
+pytest<2.2.0
 pytest-xdist
 pytest-cov
diff --git a/etc/requirements.windows.pip.txt b/etc/requirements.windows.pip.txt
index 2eb4e382a4..ddee12acbb 100644
--- a/etc/requirements.windows.pip.txt
+++ b/etc/requirements.windows.pip.txt
@@ -18,7 +18,7 @@ affine
 scipy
 pandas
 netcdf4
-pyshp>=2.0.0
+pyshp<2.2.0
 descartes
 pyproj
 shapely
diff --git a/requirements.txt b/requirements.txt
index 3f16414286..139ee51355 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -5,7 +5,7 @@ netcdf4
 fiona ; python_version < '3.8'
 descartes
 pyproj
-pyshp
+pyshp<2.2.0
 pandas
 scipy
 affine

From d05e2b25ad5e2b3d14db07ffc9cf8c94f3cd3bfc Mon Sep 17 00:00:00 2001
From: Joshua Larsen <jlarsen@usgs.gov>
Date: Fri, 4 Feb 2022 12:17:54 -0800
Subject: [PATCH 6/7] update(requirements): update pyshp requirements to pyshp
 >= 2.0.0

---
 binder/environment.yml           | 2 +-
 etc/environment.yml              | 2 +-
 etc/requirements.full.pip.txt    | 2 +-
 etc/requirements.pip.txt         | 2 +-
 etc/requirements.windows.pip.txt | 2 +-
 requirements.txt                 | 2 +-
 6 files changed, 6 insertions(+), 6 deletions(-)

diff --git a/binder/environment.yml b/binder/environment.yml
index 0d9e1f248a..1889ff19d2 100644
--- a/binder/environment.yml
+++ b/binder/environment.yml
@@ -9,7 +9,7 @@ dependencies:
   - conda-forge::scipy
   - conda-forge::pandas
   - conda-forge::netcdf4
-  - conda-forge::pyshp<2.2.0
+  - conda-forge::pyshp
   - conda-forge::rasterio
   - conda-forge::fiona
   - conda-forge::descartes
diff --git a/etc/environment.yml b/etc/environment.yml
index ca9de03292..2f8f314353 100644
--- a/etc/environment.yml
+++ b/etc/environment.yml
@@ -24,7 +24,7 @@ dependencies:
   - scipy
   - pandas
   - netcdf4
-  - pyshp<2.2.0
+  - pyshp>=2.0.0
   - rasterio
   - fiona
   - descartes
diff --git a/etc/requirements.full.pip.txt b/etc/requirements.full.pip.txt
index 5f0e12a4d7..39413e8dce 100644
--- a/etc/requirements.full.pip.txt
+++ b/etc/requirements.full.pip.txt
@@ -18,7 +18,7 @@ affine
 scipy
 pandas
 netcdf4
-pyshp<2.2.0
+pyshp>=2.0.0
 rasterio
 fiona
 descartes
diff --git a/etc/requirements.pip.txt b/etc/requirements.pip.txt
index bf0b434c73..31feb1653d 100644
--- a/etc/requirements.pip.txt
+++ b/etc/requirements.pip.txt
@@ -7,6 +7,6 @@ requests
 appdirs
 numpy
 matplotlib
-pytest<2.2.0
+pytest
 pytest-xdist
 pytest-cov
diff --git a/etc/requirements.windows.pip.txt b/etc/requirements.windows.pip.txt
index ddee12acbb..2eb4e382a4 100644
--- a/etc/requirements.windows.pip.txt
+++ b/etc/requirements.windows.pip.txt
@@ -18,7 +18,7 @@ affine
 scipy
 pandas
 netcdf4
-pyshp<2.2.0
+pyshp>=2.0.0
 descartes
 pyproj
 shapely
diff --git a/requirements.txt b/requirements.txt
index 139ee51355..e3c5c34e2a 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -5,7 +5,7 @@ netcdf4
 fiona ; python_version < '3.8'
 descartes
 pyproj
-pyshp<2.2.0
+pyshp>=2.0.0
 pandas
 scipy
 affine

From c3e7db5d5e55bc1f18f60d59e4a6780909af6653 Mon Sep 17 00:00:00 2001
From: Joshua Larsen <jlarsen@usgs.gov>
Date: Fri, 4 Feb 2022 12:21:06 -0800
Subject: [PATCH 7/7] update(requirements.txt): resolve conflict

---
 requirements.txt | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/requirements.txt b/requirements.txt
index e3c5c34e2a..3f16414286 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -5,7 +5,7 @@ netcdf4
 fiona ; python_version < '3.8'
 descartes
 pyproj
-pyshp>=2.0.0
+pyshp
 pandas
 scipy
 affine