add notebook example for specific parameters

Author: mjreno

.docs/Notebooks/mf6_netcdf01_tutorial.py

@@ -2131,7 +2131,7 @@ def test_grid_array(function_tmpdir):
 
     assert not wel.has_stress_period_data
     q_nan = np.where(wel.q.array == DNODATA, np.nan, wel.q.array)
-    val_q = np.nansum(q_nan, axis=(1, 2, 3, 4))
+    val_q = np.nansum(q_nan, axis=(1, 2, 3))
     assert val_q[0] == 0.0
     assert val_q[1] == 0.25
     assert val_q[2] == 0.1
@@ -2189,7 +2189,6 @@ def test_grid_array(function_tmpdir):
     assert len(wel.aux.get_data()) == 4
     assert wel.q.get_data()[0] is None
     assert wel.q.get_data(0) is None
-    wel_q_array = wel.q.array
     assert np.allclose(wel.q.get_data()[1], wel.q.get_data(1))
     assert np.allclose(wel.q.get_data()[2], wel.q.get_data(2))
     assert np.allclose(wel.q.array[1], wel.q.get_data(1))
@@ -2199,12 +2198,15 @@ def test_grid_array(function_tmpdir):
     assert np.allclose(wel.aux.array[1][0], wel.aux.get_data(1)[0])
     assert np.allclose(wel.aux.array[2][0], wel.aux.get_data(2)[0])
     assert not wel.has_stress_period_data
+    wel_q_array = wel.q.array
+    print(wel_q_array)
     q_nan = np.where(wel_q_array == DNODATA, np.nan, wel_q_array)
-    val_q = np.nansum(q_nan, axis=(1, 2, 3, 4))
+    val_q = np.nansum(q_nan, axis=(1, 2, 3))
+    print(val_q)
     assert val_q[0] == 0.0
     assert val_q[1] == 0.25
     assert val_q[2] == 0.1
-    assert val_q[3] == 0.1
+    assert val_q[3] == 0.0
     val_q_2 = wel.q.get_data()
     assert val_q_2[0] is None
     assert val_q_2[1][0, 0, 0] == 0.25
@@ -2416,9 +2418,9 @@ def test_grid_array(function_tmpdir):
 
     wel_q_array = wel.q.array
     assert np.all(wel_q_array[0][0] == 0.0)
-    assert wel_q_array[1][0][0, 0, 0] == 0.25
-    assert wel_q_array[2][0][0, 0, 0] == 0.1
-    assert wel_q_array[3][0][0, 0, 0] == 0.1
+    assert wel_q_array[1][0, 0, 0] == 0.25
+    assert wel_q_array[2][0, 0, 0] == 0.1
+    assert wel_q_array[3][0, 0, 0] == 0.1
     welg_q_per = wel.q.get_data()
     assert welg_q_per[0] is None
     assert welg_q_per[1][0, 0, 0] == 0.25
@@ -2437,20 +2439,20 @@ def test_grid_array(function_tmpdir):
 
     # first axis is nper, second is naux, then grid
     drn_elev_array = drn.elev.array
-    assert drn_elev_array[0][0][0, 0, 0] == 60.0
-    assert drn_elev_array[1][0][0, 0, 0] == 60.0
-    assert drn_elev_array[2][0][0, 0, 0] == DNODATA
-    assert drn_elev_array[3][0][0, 0, 0] == 55.0
+    assert drn_elev_array[0][0, 0, 0] == 60.0
+    assert drn_elev_array[1][0, 0, 0] == 60.0
+    assert drn_elev_array[2][0, 0, 0] == DNODATA
+    assert drn_elev_array[3][0, 0, 0] == 55.0
     assert np.allclose(drn_elev_array[0], drn.elev.get_data(0))
     assert drn.elev.get_data(1) is None
     assert np.allclose(drn_elev_array[2], drn.elev.get_data(2))
     assert np.allclose(drn_elev_array[3], drn.elev.get_data(3))
 
     ghb_bhead_array = ghb.bhead.array
-    assert ghb_bhead_array[0][0][0, 1, 1] == 60.0
-    assert ghb_bhead_array[1][0][0, 1, 1] == 60.0
-    assert ghb_bhead_array[2][0][0, 1, 1] == 60.0
-    assert ghb_bhead_array[3][0][0, 1, 1] == 60.0
+    assert ghb_bhead_array[0][0, 1, 1] == 60.0
+    assert ghb_bhead_array[1][0, 1, 1] == 60.0
+    assert ghb_bhead_array[2][0, 1, 1] == 60.0
+    assert ghb_bhead_array[3][0, 1, 1] == 60.0
     assert np.allclose(ghb_bhead_array[0], ghb.bhead.get_data(0))
     assert ghb.bhead.get_data(1) is None
     assert ghb.bhead.get_data(2) is None

autotest/test_mf6.py

@@ -1210,7 +1210,7 @@ def read_usgs_model_reference_file(self, reffile="usgs.model.reference"):
                                 elif info[0] == "proj4":
                                     self.crs = data
                                 else:
-                                    print("   ->warn: update not applied.")
+                                    print(" -> warn: update not applied.")
 
             # model must be rotated first, before setting xoff and yoff
             # when xul and yul are provided.

flopy/discretization/grid.py

@@ -831,7 +831,7 @@ def read_usgs_model_reference_file(self, reffile="usgs.model.reference"):
                                     else:
                                         start_date_time = data
                                 else:
-                                    print("   ->warn: update not applied.")
+                                    print(" -> warn: update not applied.")
 
             if len(start_date_time) > 0:
                 self.start_datetime = start_date_time

flopy/discretization/modeltime.py

@@ -1835,7 +1835,7 @@ def dataset(self, modeltime=None, mesh=None, configuration=None):
     def _layered_mesh_dataset(self, ds, modeltime=None, configuration=None):
         FILLNA_INT32 = np.int32(-2147483647)
         FILLNA_DBL = 9.96920996838687e36
-        lenunits = {0: "m", 1: "ft", 2: "m", 3: "m"}
+        lenunits = {0: "u", 1: "ft", 2: "m", 3: "cm"}
 
         # create dataset coordinate vars
         var_d = {
@@ -1958,7 +1958,7 @@ def _layered_mesh_dataset(self, ds, modeltime=None, configuration=None):
         return ds
 
     def _structured_dataset(self, ds, modeltime=None, configuration=None):
-        lenunits = {0: "m", 1: "ft", 2: "m", 3: "m"}
+        lenunits = {0: "u", 1: "ft", 2: "m", 3: "cm"}
 
         xc = self.xoffset + self.xycenters[0]
         yc = self.yoffset + self.xycenters[1]
@@ -2024,6 +2024,8 @@ def _structured_dataset(self, ds, modeltime=None, configuration=None):
             and configuration["longitude"] is not None
         )
 
+        lats = None
+        lons = None
         if latlon_cfg:
             lats = configuration["latitude"]
             lons = configuration["longitude"]

flopy/discretization/structuredgrid.py

@@ -614,7 +614,7 @@ def dataset(self, modeltime=None, mesh=None, configuration=None):
 
         FILLNA_INT32 = np.int32(-2147483647)
         FILLNA_DBL = 9.96920996838687e36
-        lenunits = {0: "m", 1: "ft", 2: "m", 3: "m"}
+        lenunits = {0: "u", 1: "ft", 2: "m", 3: "cm"}
 
         if mesh is None or mesh.upper() != "LAYERED":
             raise ValueError("Vextex grid only supports layered mesh datasets")

flopy/discretization/vertexgrid.py

@@ -735,6 +735,7 @@ def _get_data(self, layer=None, apply_mult=False, **kwargs):
                     and isinstance(self, MFTransientArray)
                     and data is not []  # noqa: F632
                     and not self._is_grid_aux()
+                    and not "nodes" in self.structure.shape
                 ):
                     data = np.expand_dims(data, 0)
                 return data

flopy/mf6/data/mfdataarray.py

@@ -2411,7 +2411,6 @@ def _build_full_data(self, apply_multiplier=False):
             layers_to_process = [0]
         else:
             layers_to_process = self.layer_storage.indexes()
-            layers_to_process = self.layer_storage.indexes()
         for layer in layers_to_process:
             if (
                 self.layer_storage[layer].factor is not None

flopy/mf6/data/mfdatastorage.py

@@ -1402,8 +1402,9 @@ def write(
             pp._set_netcdf_storage(reset=True)
 
         # write netcdf file
-        if write_netcdf and netcdf != "":
-            self._write_netcdf(mesh=netcdf)
+        if write_netcdf:
+            if netcdf != "":
+                self._write_netcdf(mesh=netcdf)
             if nc_fname is not None:
                 self.name_file.nc_filerecord = None