merge upstream/develop

.github/workflows/ci_examples.yaml

@@ -50,6 +50,23 @@ jobs:
           cache-environment: true
           cache-downloads: true
 
+      - name: Checkout MODFLOW 6
+        uses: actions/checkout@v4
+        with:
+          repository: MODFLOW-USGS/modflow6
+          ref: develop
+          path: modflow6
+
+      - name: Update flopy MODFLOW 6 classes
+        working-directory: modflow6/autotest
+        run: |
+          python update_flopy.py
+
+      - name: Install mf6 nightly build binaries
+        uses: modflowpy/install-modflow-action@v1
+        with:
+          repo: modflow6-nightly-build
+
       - name: Install error reporter
         run: |
           pip install pytest-github-actions-annotate-failures

autotest/test_mmr_to_mf6_dfw.py

@@ -15,14 +15,14 @@
 # Getting the gis files causes problems in parallel in CI. See regression test.
 
 # See below to check if these answers are still up-to-date with
-# mf6/autotest/test_swf_dfw.py
+# mf6/autotest/test_chf_dfw.py
 # Note the answers are from the binary FLW files in
-# mf6/autotest/test_swf_dfw.py, if you switch to text files, the line should be
+# mf6/autotest/test_chf_dfw.py, if you switch to text files, the line should be
 # flw_list = [
 #     (int(binary), 100),
 # ]  # one-based cell numbers here if binary, zero-based if text
 
-answers_swf_dfw = {
+answers_chf_dfw = {
     "ia": np.array([0, 2, 5, 7]),
     "ja": np.array([0, 1, 1, 0, 2, 2, 1], dtype=np.int32),
     "stage": np.array([[[[1.00196123, 1.00003366, 1.0]]]]),
@@ -52,9 +52,9 @@
 @pytest.mark.skipif(mf6_bin_unavailable, reason="mf6 binary not available")
 @pytest.mark.domainless
 @pytest.mark.parametrize("binary_flw", [True, False])
-def test_mmr_to_mf6_swf_dfw(tmp_path, binary_flw):
+def test_mmr_to_mf6_chf_dfw(tmp_path, binary_flw):
     # The point of this test is to reproduce the
-    # modflow6/autotest/test_swf_dfw.py
+    # modflow6/autotest/test_chf_dfw.py
 
     # Here we supply "seg_mid_elevation" in the parameter data and
     # "stress_period_data" in chd options. This bypasses the calculation of
@@ -68,8 +68,8 @@ def test_mmr_to_mf6_swf_dfw(tmp_path, binary_flw):
     # stream segments, the other is about the units of volume/flow being in
     # cubicfeet.
 
-    name = "swf-dfw01"
-    output_dir = tmp_path / "test_swf_dfw01"
+    name = "chf-dfw01"
+    output_dir = tmp_path / "test_chf_dfw01"
     save_flows = True
     print_flows = True
 
@@ -123,16 +123,16 @@ def test_mmr_to_mf6_swf_dfw(tmp_path, binary_flw):
     # vertices could also be supplied by shapefiles
     vertices = []
     vertices = [[j, j * dx, 0.0] for j in range(nreach + 1)]
-    cell2d = []
+    cell1d = []
     for j in range(nreach):
-        cell2d.append([j, 0.5, 2, j, j + 1])
-    nodes = len(cell2d)
+        cell1d.append([j, 0.5, 2, j, j + 1])
+    nodes = len(cell1d)
     nvert = len(vertices)
     disv1d_options = {
         "nodes": nodes,
         "nvert": nvert,
         "vertices": vertices,
-        "cell2d": cell2d,
+        "cell1d": cell1d,
     }
 
     # dfw
@@ -246,26 +246,30 @@ def test_mmr_to_mf6_swf_dfw(tmp_path, binary_flw):
     # Checks
     assert success
 
-    # one can verify the answers match the current mf6 results for test_swf_dfw
-    # by placing the path to its output here
+    # one can verify the answers match the current mf6 results for test_chf_dfw
+    # by placing the path to its output here. Run the following lines to
+    # generate the mf6 results:
+    #   cd modflow6_for_pws_ci/autotest  # or your mf6 repo location
+    #   pytest -s -vv test_chf_dfw.py --keep=keepers
     # output_dir = pl.Path(
-    #     "../../modflow6/autotest/keepers/test_mf6model[0-swf-dfw01]0"
+    #     "../../modflow6_for_pws_ci/autotest/keepers/"
+    #     "test_mf6model[0-chf-dfw01]0"
     # )
 
     # check binary grid file
     grb = flopy.mf6.utils.MfGrdFile(output_dir / f"{name}.disv1d.grb")
     ia = grb.ia
     ja = grb.ja
-    assert (answers_swf_dfw["ia"] == ia).all()
-    assert (answers_swf_dfw["ja"] == ja).all()
+    assert (answers_chf_dfw["ia"] == ia).all()
+    assert (answers_chf_dfw["ja"] == ja).all()
     assert ia.shape[0] == grb.nodes + 1, "ia in grb file is not correct size"
 
     # check stage file
     qobj = flopy.utils.HeadFile(
         output_dir / f"{name}.stage", precision="double", text="STAGE"
     )
     stage = qobj.get_alldata()
-    assert ((answers_swf_dfw["stage"] - stage) < 1.0e-7).all()
+    assert ((answers_chf_dfw["stage"] - stage) < 1.0e-7).all()
 
     # read the budget file
     budobj = flopy.utils.binaryfile.CellBudgetFile(output_dir / f"{name}.bud")
@@ -275,17 +279,17 @@ def test_mmr_to_mf6_swf_dfw(tmp_path, binary_flw):
     qchd = np.array(budobj.get_data(text="CHD")[0].tolist()[0])
     qresidual = np.zeros(grb.nodes)[0]
 
-    assert (answers_swf_dfw["flowja"] - flowja < 1.0e-7).all()
-    assert (answers_swf_dfw["qstorage"] - qstorage < 1.0e-7).all()
-    assert (answers_swf_dfw["qflw"] - qflw < 1.0e-7).all()
-    assert (answers_swf_dfw["qchd"] - qchd < 1.0e-7).all()
-    assert (answers_swf_dfw["qresidual"] - qresidual < 1.0e-7).all()
+    assert (answers_chf_dfw["flowja"] - flowja < 1.0e-7).all()
+    assert (answers_chf_dfw["qstorage"] - qstorage < 1.0e-7).all()
+    assert (answers_chf_dfw["qflw"] - qflw < 1.0e-7).all()
+    assert (answers_chf_dfw["qchd"] - qchd < 1.0e-7).all()
+    assert (answers_chf_dfw["qresidual"] - qresidual < 1.0e-7).all()
 
 
 # <
 answers_regression_means = {
-    "stage_all": 1.03667372881148,
-    "flow_all": 44.685014111989425,
+    "stage_all": 1.0372047024253908,
+    "flow_all": 44.70210250910929,
 }
 
 
@@ -468,4 +472,5 @@ def get_outflow(itime):
 
     for kk, vv in answers_regression_means.items():
         abs_diff = abs(locals()[kk].mean() - vv)
+
         assert abs_diff < 1e-5, f"results for {kk} are not close"

autotest/test_obsin_flow_node.py

@@ -1,6 +1,6 @@
 import numpy as np
 import pytest
-from pyPRMS.DataFile import DataFile as PRMSStreamflowData
+from pyPRMS import DataFile as PRMSStreamflowData
 
 from pywatershed import PRMSChannel
 from pywatershed.base.adapter import Adapter, AdapterNetcdf, adapter_factory

examples/06_flow_graph_starfit.ipynb

@@ -48,7 +48,11 @@
     "plot_height = 600\n",
     "plot_width = 1000\n",
     "\n",
-    "jupyter_black.load()"
+    "jupyter_black.load()\n",
+    "\n",
+    "# to remove:\n",
+    "import warnings\n",
+    "warnings.filterwarnings(\"ignore\")"
    ]
   },
   {
@@ -292,26 +296,8 @@
     "    .sel(nhm_seg=44438)\n",
     "    .rename(\"modeled\")\n",
     "    .to_dataframe()[\"modeled\"]\n",
-    ")"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "id": "3b385029-6554-4f4d-9fd9-4c255d70e40d",
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "obs_all = pyPRMS.DataFile(domain_dir / \"sf_data\").data_by_variable(\"runoff\")"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "id": "54058459-fed5-4cc3-bed2-8a3edfb22a13",
-   "metadata": {},
-   "outputs": [],
-   "source": [
+    ")\n",
+    "obs_all = pyPRMS.DataFile(domain_dir / \"sf_data\").data_by_variable(\"runoff\")\n",
     "wh_poi_obs = np.where(params.parameters[\"poi_gage_segment\"] == 184)\n",
     "gage_id = params.parameters[\"poi_gage_id\"][wh_poi_obs][0]\n",
     "obs = obs_all[f\"runoff_{gage_id}\"]\n",

examples/mmr_to_mf6_dfw.ipynb → examples/07_mmr_to_mf6_chf_dfw.ipynb

@@ -7,9 +7,9 @@
    "source": [
     "# MMR To MF6 DFW\n",
     "\n",
-    "This notebook performs MF6 diffussive wave routing (DFW) simulaton based on PRMS MMR (Muskingum-Mann Routing) inputs and boundary flows from the PRMS model. We'll compare and contrast the MF6 DFW and PRMS MMR simulations and delve into their differences at the largest gaged flows in the domain (that are not affected by tides). \n",
+    "This notebook performs 1-D diffussive wave (DFW) routing in MODFLOW 6 using the CHF (channel flow) model. The static and time-varying boundary conditions, come from PRMS. Specifically, PRMS's MMR (Muskingum-Mann Routing) parameters and its channel inflows calculated during a pywatershed run. We'll compare and contrast the MF6 CHF-DFW and PRMS MMR simulations and delve into their differences at the largest gaged flows in the domain (that are not affected by tides). \n",
     "\n",
-    "This notebook requries the develop branch of MF6 and a flopy upto date with this branch by running `python update_flopy.py` in modflow6/autotest. You will also need to add the MF6 executable to your path below. You will also need to have run \n",
+    "This notebook requries the develop branch of MF6 and a flopy upto date with this branch by running `python update_flopy.py` in modflow6/autotest. You will also need to add the MF6 executable to your path below. \n",
     "\n",
     "## User configuration"
    ]
@@ -21,16 +21,20 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "# Set YOUR path to MF6 here\n",
+    "# Set YOUR path to MF6 in this block\n",
     "import pathlib as pl\n",
     "\n",
     "mf6_bin = pl.Path(\"../../modflow6/bin/mf6\")\n",
+    "# double check\n",
     "msg = \"A build of mf6/develop branch required for DFW simulation\"\n",
     "assert mf6_bin.exists, msg\n",
     "\n",
     "# Rerun MF6 model below or use an existing run?\n",
     "rerun_mf6 = True\n",
-    "rerun_prms = True"
+    "rerun_prms = True\n",
+    "\n",
+    "# Perform the full, 2 year run period or just the first 45 days?\n",
+    "full_run_period = False"
    ]
   },
   {
@@ -96,7 +100,7 @@
    "id": "e6a36353-3a61-4b57-8db5-e48293ffc0d1",
    "metadata": {},
    "source": [
-    "## Run PRMS\n",
+    "## Run PRMS NHM configuration using pywatershed\n",
     "Running PRMS gives the boundary conditions for the MF6 DFW run and it also produces its own streamflow simulation with its Muskingum-Mann routing method. "
    ]
   },
@@ -107,11 +111,9 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "%time\n",
-    "\n",
-    "prms_run_dir = repo_root_dir / \"examples/mmr_to_mf6_dfw/prms_run\"\n",
+    "prms_run_dir = repo_root_dir / \"examples/07_mmr_to_mf6_chf_dfw/prms_run\"\n",
     "\n",
-    "if rerun_prms:\n",
+    "if rerun_prms and prms_run_dir.exists():\n",
     "    shutil.rmtree(prms_run_dir)\n",
     "\n",
     "if not prms_run_dir.exists():\n",
@@ -186,12 +188,10 @@
     "control_file = domain_dir / \"nhm.control\"\n",
     "control = pws.Control.load_prms(control_file)\n",
     "ndays_run = control.n_times\n",
-    "# Could shorten the run duration\n",
-    "# Subtract one from ndays_run becase end day/time would be included in the PRMS run\n",
-    "# ndays_run = 45\n",
-    "# control.edit_end_time(\n",
-    "#     control.start_time + ((ndays_run - 1) * control.time_step)\n",
-    "# )"
+    "\n",
+    "if not full_run_period:\n",
+    "    ndays_run = 45\n",
+    "    control.edit_n_time_steps(ndays_run)"
    ]
   },
   {
@@ -482,7 +482,10 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "tt = 3000\n",
+    "if ndays_run < 3000:\n",
+    "    tt = ndays_run\n",
+    "else:\n",
+    "    tt = 3000\n",
     "zoom = 1.0\n",
     "figsize = (7 * zoom, 10 * zoom)\n",
     "dt = tdis_perlen / tdis_nstp\n",