From 8a3a7d0270bd55be321308534545fe7ea561db07 Mon Sep 17 00:00:00 2001 From: luitjan Date: Fri, 7 Jun 2024 09:46:14 +0200 Subject: [PATCH 1/5] work started --- imod/mf6/hfb.py | 22 +++++++++++++++++++++- imod/mf6/model_gwf.py | 6 ++++++ imod/tests/test_mf6/test_mf6_hfb.py | 16 ++++++++++++++++ imod/tests/test_mf6/test_mf6_simulation.py | 4 +++- 4 files changed, 46 insertions(+), 2 deletions(-) diff --git a/imod/mf6/hfb.py b/imod/mf6/hfb.py index 5163005fa..e5bbc8ab8 100644 --- a/imod/mf6/hfb.py +++ b/imod/mf6/hfb.py @@ -1080,5 +1080,25 @@ def _compute_barrier_values( edge_index, idomain, ) - return barrier_values + + @classmethod + def from_imod5_dataset(cls, imod5_data: dict[str, dict[str, GridDataArray]]): + imod5_keys = list(imod5_data.keys()) + hfb_keys = [key for key in imod5_keys if key[0:3] == "hfb"] + + first = True + for hfb_key in hfb_keys: + hfb_dict = imod5_data[hfb_key] + if not list(hfb_dict.keys()) == ["geodataframe", "layer"]: + raise ValueError + layer = hfb_dict["layer"] + geometry_layer = deepcopy(hfb_dict["geodataframe"]) + geometry_layer["layer"] = layer + if first: + compound_dataframe = deepcopy(geometry_layer) + first = False + else: + compound_dataframe = compound_dataframe._append(geometry_layer) + + return LayeredHorizontalFlowBarrierResistance(compound_dataframe) diff --git a/imod/mf6/model_gwf.py b/imod/mf6/model_gwf.py index 8736abd47..4e7d7e799 100644 --- a/imod/mf6/model_gwf.py +++ b/imod/mf6/model_gwf.py @@ -11,6 +11,7 @@ from imod.mf6.clipped_boundary_condition_creator import create_clipped_boundary from imod.mf6.dis import StructuredDiscretization from imod.mf6.drn import Drainage +from imod.mf6.hfb import LayeredHorizontalFlowBarrierResistance from imod.mf6.ic import InitialConditions from imod.mf6.model import Modflow6Model from imod.mf6.npf import NodePropertyFlow @@ -244,4 +245,9 @@ def from_imod5_data( ) result[drn_key] = drn_pkg + # import hfb + hfb_keys = [key for key in imod5_keys if key[0:3] == "hfb"] + if len(hfb_keys) != 0: + hfb = LayeredHorizontalFlowBarrierResistance.from_imod5_dataset(imod5_data) + result["hfb"] = hfb return result diff --git a/imod/tests/test_mf6/test_mf6_hfb.py b/imod/tests/test_mf6/test_mf6_hfb.py index 3e7c61109..4b1f31fc7 100644 --- a/imod/tests/test_mf6/test_mf6_hfb.py +++ b/imod/tests/test_mf6/test_mf6_hfb.py @@ -16,7 +16,9 @@ LayeredHorizontalFlowBarrierMultiplier, LayeredHorizontalFlowBarrierResistance, ) +from imod.mf6.dis import StructuredDiscretization from imod.mf6.hfb import to_connected_cells_dataset +from imod.mf6.npf import NodePropertyFlow from imod.mf6.utilities.regrid import RegridderWeightsCache from imod.tests.fixtures.flow_basic_fixture import BasicDisSettings from imod.typing.grid import ones_like @@ -440,3 +442,17 @@ def test_set_options(print_input, parameterizable_basic_dis): k = ones_like(top) mf6_package = hfb.to_mf6_pkg(idomain, top, bottom, k) assert mf6_package.dataset["print_input"].values[()] == print_input + + +@pytest.mark.usefixtures("imod5_dataset") +def test_hfb_from_imod5(imod5_dataset, tmp_path): + target_dis = StructuredDiscretization.from_imod5_data(imod5_dataset) + target_npf = NodePropertyFlow.from_imod5_data( + imod5_dataset, target_dis.dataset["idomain"] + ) + + hfb = LayeredHorizontalFlowBarrierResistance.from_imod5_dataset(imod5_dataset) + hfb_pack = hfb.to_mf6_pkg( + target_dis["idomain"], target_dis["top"], target_dis["bottom"], target_npf["k"] + ) + pass diff --git a/imod/tests/test_mf6/test_mf6_simulation.py b/imod/tests/test_mf6/test_mf6_simulation.py index fc2373fff..d997b1245 100644 --- a/imod/tests/test_mf6/test_mf6_simulation.py +++ b/imod/tests/test_mf6/test_mf6_simulation.py @@ -479,7 +479,9 @@ def test_import_from_imod5(imod5_dataset, tmp_path): default_simulation_distributing_options, ) - simulation["imported_model"] = OutputControl(save_head="last", save_budget="last") + simulation["imported_model"]["oc"] = OutputControl( + save_head="last", save_budget="last" + ) simulation.create_time_discretization(["01-01-2003", "02-01-2003"]) From 3b12ce583c2c11063bc26ef93495eddbdb7b11d4 Mon Sep 17 00:00:00 2001 From: luitjan Date: Fri, 7 Jun 2024 10:18:00 +0200 Subject: [PATCH 2/5] error messages --- imod/mf6/hfb.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/imod/mf6/hfb.py b/imod/mf6/hfb.py index e5bbc8ab8..4a4811051 100644 --- a/imod/mf6/hfb.py +++ b/imod/mf6/hfb.py @@ -1086,17 +1086,19 @@ def _compute_barrier_values( def from_imod5_dataset(cls, imod5_data: dict[str, dict[str, GridDataArray]]): imod5_keys = list(imod5_data.keys()) hfb_keys = [key for key in imod5_keys if key[0:3] == "hfb"] - + if len(hfb_keys) ==0: + raise ValueError("no hfb keys present.") + first = True for hfb_key in hfb_keys: hfb_dict = imod5_data[hfb_key] if not list(hfb_dict.keys()) == ["geodataframe", "layer"]: - raise ValueError + raise ValueError("hfb is not a LayeredHorizontalFlowBarrierResistance") layer = hfb_dict["layer"] geometry_layer = deepcopy(hfb_dict["geodataframe"]) geometry_layer["layer"] = layer if first: - compound_dataframe = deepcopy(geometry_layer) + compound_dataframe = geometry_layer first = False else: compound_dataframe = compound_dataframe._append(geometry_layer) From 3e2b8c5c937558055933ce9a14fc02ffeb7d5e20 Mon Sep 17 00:00:00 2001 From: luitjan Date: Fri, 7 Jun 2024 10:52:22 +0200 Subject: [PATCH 3/5] improved test for hfb's --- imod/mf6/hfb.py | 7 +++---- imod/tests/test_mf6/test_mf6_hfb.py | 4 ++-- 2 files changed, 5 insertions(+), 6 deletions(-) diff --git a/imod/mf6/hfb.py b/imod/mf6/hfb.py index 4a4811051..73b89e701 100644 --- a/imod/mf6/hfb.py +++ b/imod/mf6/hfb.py @@ -1089,17 +1089,16 @@ def from_imod5_dataset(cls, imod5_data: dict[str, dict[str, GridDataArray]]): if len(hfb_keys) ==0: raise ValueError("no hfb keys present.") - first = True + compound_dataframe = None for hfb_key in hfb_keys: hfb_dict = imod5_data[hfb_key] if not list(hfb_dict.keys()) == ["geodataframe", "layer"]: raise ValueError("hfb is not a LayeredHorizontalFlowBarrierResistance") layer = hfb_dict["layer"] - geometry_layer = deepcopy(hfb_dict["geodataframe"]) + geometry_layer = hfb_dict["geodataframe"] geometry_layer["layer"] = layer - if first: + if compound_dataframe is None: compound_dataframe = geometry_layer - first = False else: compound_dataframe = compound_dataframe._append(geometry_layer) diff --git a/imod/tests/test_mf6/test_mf6_hfb.py b/imod/tests/test_mf6/test_mf6_hfb.py index 4b1f31fc7..81a242d16 100644 --- a/imod/tests/test_mf6/test_mf6_hfb.py +++ b/imod/tests/test_mf6/test_mf6_hfb.py @@ -452,7 +452,7 @@ def test_hfb_from_imod5(imod5_dataset, tmp_path): ) hfb = LayeredHorizontalFlowBarrierResistance.from_imod5_dataset(imod5_dataset) - hfb_pack = hfb.to_mf6_pkg( + hfb_package = hfb.to_mf6_pkg( target_dis["idomain"], target_dis["top"], target_dis["bottom"], target_npf["k"] ) - pass + assert list(np.unique(hfb_package["layer"].values)) == [3, 21] From 2a78d13620003f4933e590ec4af82fab2badc4a1 Mon Sep 17 00:00:00 2001 From: luitjan Date: Fri, 7 Jun 2024 10:55:06 +0200 Subject: [PATCH 4/5] formatting --- imod/mf6/hfb.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/imod/mf6/hfb.py b/imod/mf6/hfb.py index 73b89e701..c2a46e6c0 100644 --- a/imod/mf6/hfb.py +++ b/imod/mf6/hfb.py @@ -1086,9 +1086,9 @@ def _compute_barrier_values( def from_imod5_dataset(cls, imod5_data: dict[str, dict[str, GridDataArray]]): imod5_keys = list(imod5_data.keys()) hfb_keys = [key for key in imod5_keys if key[0:3] == "hfb"] - if len(hfb_keys) ==0: + if len(hfb_keys) == 0: raise ValueError("no hfb keys present.") - + compound_dataframe = None for hfb_key in hfb_keys: hfb_dict = imod5_data[hfb_key] From db38be0985286ca12f04238cd2b110a34eb8c579 Mon Sep 17 00:00:00 2001 From: luitjan Date: Fri, 7 Jun 2024 14:24:30 +0200 Subject: [PATCH 5/5] added warning for layer=0 --- imod/mf6/hfb.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/imod/mf6/hfb.py b/imod/mf6/hfb.py index c2a46e6c0..c61146a35 100644 --- a/imod/mf6/hfb.py +++ b/imod/mf6/hfb.py @@ -1095,6 +1095,10 @@ def from_imod5_dataset(cls, imod5_data: dict[str, dict[str, GridDataArray]]): if not list(hfb_dict.keys()) == ["geodataframe", "layer"]: raise ValueError("hfb is not a LayeredHorizontalFlowBarrierResistance") layer = hfb_dict["layer"] + if layer == 0: + raise ValueError( + "assigning to layer 0 is not supported yet for import of HFB's" + ) geometry_layer = hfb_dict["geodataframe"] geometry_layer["layer"] = layer if compound_dataframe is None: