fix issue and add test

imod_coupler/drivers/ribametamod/exchange.py

@@ -187,17 +187,18 @@ def flux_to_ribasim(self, delt_gw: float, delt_sw: float) -> None:
         demand_per_subtimestep = self.get_demand_flux_sec(delt_gw, delt_sw)
 
         # exchange to Ribasim; negative demand in exchange class means infiltration from Ribasim
-        coupled_index = self.mapping.coupled_mod2rib
-        self.ribasim_infiltration[coupled_index] = np.where(
+        self.ribasim_infiltration[self.mapping.coupled_index] = np.where(
             demand_per_subtimestep < 0, -demand_per_subtimestep, 0
-        )[coupled_index]
-        self.ribasim_drainage[coupled_index] = np.where(
+        )[self.mapping.coupled_index]
+        self.ribasim_drainage[self.mapping.coupled_index] = np.where(
             demand_per_subtimestep > 0, demand_per_subtimestep, 0
-        )[coupled_index]
+        )[self.mapping.coupled_index]
 
     def flux_to_modflow(
         self, realised_volume: NDArray[np.float64], delt_gw: float
     ) -> None:
+        if not self.mf6_river_packages:
+            return  # no active coupling
         super().compute_realised(realised_volume)
         for key in self.mf6_active_river_api_packages.keys():
             realised_fraction = np.where(

imod_coupler/drivers/ribametamod/mapping.py

@@ -32,10 +32,12 @@ def __init__(
         self.coupling = coupling
         if has_ribasim:
             self.set_ribasim_modflow_mapping(packages)
+            self.coupled_index = self.coupled_mod2rib
         if has_metaswap and mod2svat is not None:
             self.set_metaswap_modflow_mapping(packages, mod2svat)
         if has_ribasim and has_metaswap and mod2svat is not None:
             self.set_metaswap_ribasim_mapping(packages, mod2svat)
+            self.coupled_index |= self.coupled_msw2rib
 
     def set_ribasim_modflow_mapping(self, packages: ChainMap[str, Any]) -> None:
         self.map_mod2rib = {}
@@ -200,6 +202,9 @@ def set_metaswap_ribasim_mapping(
         self, packages: ChainMap[str, Any], mod2svat: Path
     ) -> None:
         table_node2svat: NDArray[np.int32]
+        self.coupled_msw2rib: NDArray[np.bool_] = np.full(
+            packages["ribasim_nbasin"], False
+        )
         self.msw2rib = {}
 
         # surface water ponding mapping
@@ -222,6 +227,7 @@ def set_metaswap_ribasim_mapping(
                 packages["ribasim_nbasin"],
                 "sum",
             )
+            self.coupled_msw2rib |= self.msw2rib["sw_ponding"].getnnz(axis=1) > 0
 
         # surface water sprinkling mapping
         if self.coupling.rib_msw_sprinkling_map_surface_water is not None:

imod_coupler/drivers/ribametamod/ribametamod.py

@@ -218,12 +218,12 @@ def couple_ribasim(self, mf6_flowmodel_key: str) -> ChainMap[str, Any]:
                     self.mf6_drainage_packages,
                     {
                         "ribasim_nbasin": len(self.ribasim_level),
-                        "ribasim_nuser": (
-                            len(self.ribasim_user_realized)
-                            if self.ribasim_user_realized.ndim > 0
-                            else 0
-                        ),
-                        "ribasim_nsubgrid": len(self.subgrid_level),
+                        "ribasim_nuser": len(self.ribasim_user_realized)
+                        if self.ribasim_user_realized.ndim > 0
+                        else 0,
+                        "ribasim_nsubgrid": len(self.subgrid_level)
+                        if self.subgrid_level.ndim > 0
+                        else 0,
                     },
                 )
             )

tests/fixtures/fixture_ribasim.py

@@ -36,6 +36,13 @@ def ribasim_bucket_model() -> ribasim.Model:
     return add_subgrid(bucket)
 
 
+@pytest_cases.fixture(scope="function")
+def ribasim_bucket_model_no_subgrid() -> ribasim.Model:
+    bucket = ribasim_testmodels.bucket_model()
+    bucket.endtime = datetime(2023, 1, 1, 0, 0)
+    return bucket
+
+
 @pytest_cases.fixture(scope="function")
 def ribasim_backwater_model() -> ribasim.Model:
     return add_subgrid(ribasim_testmodels.backwater_model())

tests/test_imod_coupler/test_ribametamod.py

@@ -545,7 +545,6 @@ def test_ribametamod_backwater(
 
 @pytest.mark.xdist_group(name="ribasim")
 @parametrize_with_cases("ribametamod_model", glob="bucket_model")
-# @pytest.mark.skip("potentially hangs, to be solved in model")
 def test_ribametamod_bucket(
     tmp_path_dev: Path,
     ribametamod_model: RibaMetaMod,
@@ -577,6 +576,34 @@ def test_ribametamod_bucket(
     assert_results(tmp_path_dev, ribametamod_model, results)
 
 
+@pytest.mark.xdist_group(name="ribasim")
+@parametrize_with_cases("ribametamod_model", glob="bucket_model_no_subgrid")
+def test_ribametamod_bucket_no_subgrid(
+    tmp_path_dev: Path,
+    ribametamod_model: RibaMetaMod,
+    modflow_dll_devel: Path,
+    ribasim_dll_devel: Path,
+    ribasim_dll_dep_dir_devel: Path,
+    metaswap_dll_devel: Path,
+    metaswap_dll_dep_dir_devel: Path,
+    run_coupler_function: Callable[[Path], None],
+) -> None:
+    """
+    Test if the bucket model runs without a subgrid in the Ribasim model
+    """
+    ribametamod_model.write(
+        tmp_path_dev,
+        modflow6_dll=modflow_dll_devel,
+        ribasim_dll=ribasim_dll_devel,
+        ribasim_dll_dependency=ribasim_dll_dep_dir_devel,
+        metaswap_dll=metaswap_dll_devel,
+        metaswap_dll_dependency=metaswap_dll_dep_dir_devel,
+        modflow6_write_kwargs={"binary": False},
+    )
+
+    run_coupler_function(tmp_path_dev / ribametamod_model._toml_name)
+
+
 @pytest.mark.xdist_group(name="ribasim")
 @parametrize_with_cases("ribametamod_model", glob="two_basin_model")
 def test_ribametamod_two_basin(

tests/test_imod_coupler/test_ribametamod_cases.py

@@ -7,6 +7,7 @@
 from fixtures.common import create_wells_max_layer
 from imod import msw
 from imod.mf6 import (
+    Drainage,
     Modflow6Simulation,
     Recharge,
     StorageCoefficient,
@@ -17,6 +18,7 @@
     RibaMetaDriverCoupling,
     RibaMetaMod,
     RibaModActiveDriverCoupling,
+    RibaModPassiveDriverCoupling,
 )
 from ribasim.geometry import NodeTable as RibaNodeTbl
 from ribasim.nodes import level_demand, user_demand
@@ -326,6 +328,59 @@ def case_bucket_model(
     )
 
 
+def case_bucket_model_no_subgrid(
+    mf6_bucket_model: Modflow6Simulation,
+    msw_bucket_model: MetaSwapModel,
+    ribasim_bucket_model_no_subgrid: ribasim.Model,
+) -> RibaMetaMod:
+    mf6_bucket_model = set_confined_storage_formulation(mf6_bucket_model)
+
+    # add rch and well-package for coupling MetaMod
+    mf6_bucket_model = add_rch_package(mf6_bucket_model)
+    mf6_bucket_model = add_well_package(mf6_bucket_model)
+
+    # get variables
+    mf6_modelname, mf6_model = get_mf6_gwf_modelnames(mf6_bucket_model)[0]
+    mf6_active_river_packages = get_mf6_river_packagenames(mf6_model)
+    basin_definition = create_basin_definition(
+        ribasim_bucket_model_no_subgrid.basin.node, buffersize=100.0
+    )
+
+    # riv to drn package
+    conductance = (
+        mf6_bucket_model["GWF_1"][mf6_active_river_packages[0]].dataset["conductance"]
+        / 400
+    )
+    stage = mf6_bucket_model["GWF_1"][mf6_active_river_packages[0]].dataset[
+        "conductance"
+    ]
+    mf6_bucket_model["GWF_1"].pop(mf6_active_river_packages[0])
+    mf6_bucket_model["GWF_1"][mf6_active_river_packages[0]] = Drainage(
+        elevation=stage, conductance=conductance
+    )
+
+    metamod_coupling = MetaModDriverCoupling(
+        mf6_model=mf6_modelname,
+        mf6_recharge_package="rch_msw",
+        mf6_wel_package="well_msw",
+    )
+    ribamod_coupling = RibaModPassiveDriverCoupling(
+        mf6_model=mf6_modelname,
+        mf6_packages=mf6_active_river_packages,
+        ribasim_basin_definition=basin_definition,
+    )
+    ribameta_coupling = RibaMetaDriverCoupling(
+        ribasim_basin_definition=basin_definition,
+    )
+
+    return RibaMetaMod(
+        ribasim_model=ribasim_bucket_model_no_subgrid,
+        msw_model=msw_bucket_model,
+        mf6_simulation=mf6_bucket_model,
+        coupling_list=[metamod_coupling, ribamod_coupling, ribameta_coupling],
+    )
+
+
 def case_backwater_model(
     mf6_backwater_model: Modflow6Simulation,
     msw_backwater_model: MetaSwapModel,