Issue #1041 import simulation (#1059)

Fixes #1041

# Description
Adds a function that imports a simulation from an imod5 project file. It
creates a simulation, containing 1 GroundwaterFlowModel, containing the
packages for which we support importing currently. This includes the
mandatory flow packages like NPF and STO, and some boundary conditions
(like DRN and RCH) but not yet CHD or WEL.

# Checklist
- [X] Links to correct issue
- [ ] Update changelog, if changes affect users
- [X] PR title starts with ``Issue #nr``, e.g. ``Issue #737``
- [X] Unit tests were added
- [ ] **If feature added**: Added/extended example

---------

Co-authored-by: Joeri van Engelen <joerivanengelen@hotmail.com>

Author: luitjansl

imod/mf6/model_gwf.py

@@ -6,9 +6,21 @@
 import numpy as np
 
 from imod.logging import init_log_decorator
+from imod.logging.logging_decorators import standard_log_decorator
 from imod.mf6 import ConstantHead
 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.ic import InitialConditions
 from imod.mf6.model import Modflow6Model
+from imod.mf6.npf import NodePropertyFlow
+from imod.mf6.rch import Recharge
+from imod.mf6.sto import StorageCoefficient
+from imod.mf6.utilities.regridding_types import RegridderType
+from imod.prepare.topsystem.default_allocation_methods import (
+    SimulationAllocationOptions,
+    SimulationDistributingOptions,
+)
 from imod.typing import GridDataArray
 
 
@@ -153,3 +165,83 @@ def update_buoyancy_package(self, transport_models_per_flow_model) -> None:
         transport_models_old = buoyancy_package.get_transport_model_names()
         if len(transport_models_old) == len(transport_models_per_flow_model):
             buoyancy_package.update_transport_models(transport_models_per_flow_model)
+
+    @classmethod
+    @standard_log_decorator()
+    def from_imod5_data(
+        cls,
+        imod5_data: dict[str, dict[str, GridDataArray]],
+        allocation_options: SimulationAllocationOptions,
+        distributing_options: SimulationDistributingOptions,
+        regridder_types: Optional[dict[str, tuple[RegridderType, str]]] = None,
+    ) -> "GroundwaterFlowModel":
+        """
+        Imports a GroundwaterFlowModel (GWF) from the data in an IMOD5 project file.
+        It adds the packages for which import from imod5 is supported.
+        Some packages (like OC) must be added manually later.
+
+
+        Parameters
+        ----------
+        imod5_data: dict[str, dict[str, GridDataArray]]
+            dictionary containing the arrays mentioned in the project file as xarray datasets,
+            under the key of the package type to which it belongs
+        allocation_options: SimulationAllocationOptions
+            object containing the allocation options per package type.
+            If you want a package to have a different allocation option,
+            then it should be imported separately
+        distributing_options: SimulationDistributingOptions
+            object containing the conductivity distribution options per package type.
+            If you want a package to have a different allocation option,
+            then it should be imported separately
+        regridder_types: Optional[dict[str, dict[str, tuple[RegridderType, str]]]]
+            the first key is the package name. The second key is the array name, and the value is
+            the RegridderType tuple (method + function)
+
+        Returns
+        -------
+        A GWF model containing the packages that could be imported form IMOD5. Users must still
+        add the OC package to the model.
+
+        """
+        # first import the singleton packages
+        # import discretization
+        dis_pkg = StructuredDiscretization.from_imod5_data(imod5_data, regridder_types)
+        grid = dis_pkg.dataset["idomain"]
+
+        # import npf
+        npf_pkg = NodePropertyFlow.from_imod5_data(imod5_data, grid, regridder_types)
+
+        # import sto
+        sto_pkg = StorageCoefficient.from_imod5_data(imod5_data, grid, regridder_types)
+
+        # import initial conditions
+        ic_pkg = InitialConditions.from_imod5_data(imod5_data, grid, regridder_types)
+
+        # import recharge
+        rch_pkg = Recharge.from_imod5_data(imod5_data, dis_pkg, regridder_types)
+
+        result = GroundwaterFlowModel()
+        result["dis"] = dis_pkg
+        result["npf"] = npf_pkg
+        result["sto"] = sto_pkg
+        result["ic"] = ic_pkg
+        result["rch"] = rch_pkg
+
+        # now import the non-singleton packages
+        # import drainage
+        imod5_keys = list(imod5_data.keys())
+        drainage_keys = [key for key in imod5_keys if key[0:3] == "drn"]
+        for drn_key in drainage_keys:
+            drn_pkg = Drainage.from_imod5_data(
+                drn_key,
+                imod5_data,
+                dis_pkg,
+                npf_pkg,
+                allocation_options.drn,
+                distributing_option=distributing_options.drn,
+                regridder_types=regridder_types,
+            )
+            result[drn_key] = drn_pkg
+
+        return result

imod/mf6/simulation.py

@@ -24,7 +24,7 @@
 from imod.mf6.gwfgwf import GWFGWF
 from imod.mf6.gwfgwt import GWFGWT
 from imod.mf6.gwtgwt import GWTGWT
-from imod.mf6.ims import Solution
+from imod.mf6.ims import Solution, SolutionPresetModerate
 from imod.mf6.interfaces.imodel import IModel
 from imod.mf6.interfaces.isimulation import ISimulation
 from imod.mf6.model import Modflow6Model
@@ -41,7 +41,12 @@
 from imod.mf6.statusinfo import NestedStatusInfo
 from imod.mf6.utilities.mask import _mask_all_models
 from imod.mf6.utilities.regrid import _regrid_like
+from imod.mf6.utilities.regridding_types import RegridderType
 from imod.mf6.write_context import WriteContext
+from imod.prepare.topsystem.default_allocation_methods import (
+    SimulationAllocationOptions,
+    SimulationDistributingOptions,
+)
 from imod.schemata import ValidationError
 from imod.typing import GridDataArray, GridDataset
 from imod.typing.grid import (
@@ -1308,3 +1313,62 @@ def mask_all_models(
             -1 sets cells to vertical passthrough
         """
         _mask_all_models(self, mask)
+
+    @classmethod
+    @standard_log_decorator()
+    def from_imod5_data(
+        cls,
+        imod5_data: dict[str, dict[str, GridDataArray]],
+        allocation_options: SimulationAllocationOptions,
+        distributing_options: SimulationDistributingOptions,
+        regridder_types: Optional[dict[str, tuple[RegridderType, str]]] = None,
+    ) -> "Modflow6Simulation":
+        """
+        Imports a GroundwaterFlowModel (GWF) from the data in an IMOD5 project file.
+        It adds the packages for which import from imod5 is supported.
+        Some packages (like OC) must be added manually later.
+
+
+        Parameters
+        ----------
+        imod5_data: dict[str, dict[str, GridDataArray]]
+            dictionary containing the arrays mentioned in the project file as xarray datasets,
+            under the key of the package type to which it belongs
+        allocation_options: SimulationAllocationOptions
+            object containing the allocation options per package type.
+            If you want a package to have a different allocation option,
+            then it should be imported separately
+        distributing_options: SimulationDistributingOptions
+            object containing the conductivity distribution options per package type.
+            If you want a package to have a different allocation option,
+            then it should be imported separately
+        regridder_types: Optional[dict[str, dict[str, tuple[RegridderType, str]]]]
+            the first key is the package name. The second key is the array name, and the value is
+            the RegridderType tuple (method + function)
+
+        Returns
+        -------
+        """
+        simulation = Modflow6Simulation("imported_simulation")
+
+        # import GWF model,
+        groundwaterFlowModel = GroundwaterFlowModel.from_imod5_data(
+            imod5_data,
+            allocation_options,
+            distributing_options,
+            regridder_types,
+        )
+        simulation["imported_model"] = groundwaterFlowModel
+
+        # generate ims package
+        solution = SolutionPresetModerate(
+            ["imported_model"],
+            print_option="all",
+        )
+        simulation["ims"] = solution
+
+        # cleanup packages for validation
+        idomain = groundwaterFlowModel.domain
+        simulation.mask_all_models(idomain)
+
+        return simulation

imod/mf6/sto.py

@@ -290,10 +290,6 @@ class StorageCoefficient(StorageBase):
             AllValueSchema(">=", 0.0),
             IdentityNoDataSchema(other="idomain", is_other_notnull=(">", 0)),
         ),
-        "convertible": (
-            IdentityNoDataSchema(other="idomain", is_other_notnull=(">", 0)),
-            # No need to check coords: dataset ensures they align with idomain.
-        ),
     }
 
     _template = Package._initialize_template(_pkg_id)
@@ -372,9 +368,7 @@ def from_imod5_data(
             data, target_grid, regridder_settings, regrid_context, {}
         )
 
-        new_package_data["convertible"] = zeros_like(
-            new_package_data["storage_coefficient"], dtype=int
-        )
+        new_package_data["convertible"] = zeros_like(target_grid, dtype=int)
         new_package_data["transient"] = np.any(
             new_package_data["storage_coefficient"].values > 0
         )

imod/prepare/topsystem/default_allocation_methods.py

@@ -0,0 +1,40 @@
+from dataclasses import dataclass
+
+from imod.prepare.topsystem.allocation import ALLOCATION_OPTION
+from imod.prepare.topsystem.conductance import DISTRIBUTING_OPTION
+
+
+@dataclass()
+class SimulationAllocationOptions:
+    """
+    Object containing allocation otpions, specified per packages type on
+    importing fron imod5. Can be used to set defaults when importing a
+    simulation or a GroundwaterFlowModel from imod5.
+
+    Parameters
+    ----------
+    drn: allocation option to be used for drainage packages
+    riv: allocation option to be used for river packages
+
+    """
+
+    drn: ALLOCATION_OPTION = ALLOCATION_OPTION.first_active_to_elevation
+    riv: ALLOCATION_OPTION = ALLOCATION_OPTION.stage_to_riv_bot
+
+
+@dataclass()
+class SimulationDistributingOptions:
+    """
+    Object containing conductivity distribution methods, specified per packages
+    type. Can be used to set defaults when importing a simulation or a
+    GroundwaterFlowModel from imod5.
+
+    Parameters
+    ----------
+    drn: distribution option to be used for drainage packages
+    riv: distribution option to be used for river packages
+
+    """
+
+    drn: DISTRIBUTING_OPTION = DISTRIBUTING_OPTION.by_corrected_transmissivity
+    riv: DISTRIBUTING_OPTION = DISTRIBUTING_OPTION.by_corrected_transmissivity

imod/tests/test_mf6/test_mf6_simulation.py

@@ -18,7 +18,13 @@
 import imod
 from imod.mf6.model import Modflow6Model
 from imod.mf6.multimodel.modelsplitter import PartitionInfo
+from imod.mf6.oc import OutputControl
+from imod.mf6.simulation import Modflow6Simulation
 from imod.mf6.statusinfo import NestedStatusInfo, StatusInfo
+from imod.prepare.topsystem.default_allocation_methods import (
+    SimulationAllocationOptions,
+    SimulationDistributingOptions,
+)
 from imod.schemata import ValidationError
 from imod.tests.fixtures.mf6_small_models_fixture import (
     grid_data_structured,
@@ -461,3 +467,21 @@ def compare_submodel_partition_info(first: PartitionInfo, second: PartitionInfo)
     return (first.id == second.id) and np.array_equal(
         first.active_domain, second.active_domain
     )
+
+
+@pytest.mark.usefixtures("imod5_dataset")
+def test_import_from_imod5(imod5_dataset, tmp_path):
+    default_simulation_allocation_options = SimulationAllocationOptions
+    default_simulation_distributing_options = SimulationDistributingOptions
+    simulation = Modflow6Simulation.from_imod5_data(
+        imod5_dataset,
+        default_simulation_allocation_options,
+        default_simulation_distributing_options,
+    )
+
+    simulation["imported_model"] = OutputControl(save_head="last", save_budget="last")
+
+    simulation.create_time_discretization(["01-01-2003", "02-01-2003"])
+
+    # write and validate the simulation
+    simulation.write(tmp_path, False, True, False)