### Set Initial Condition (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Configures the initial condition for the rootzone pressure head. This example sets a uniform initial pF value of 2.2 across all cells.

```python
msw_model["ic"] = msw.InitialConditionsRootzonePressureHead(initial_pF=2.2)
```

--------------------------------

### Build Documentation with Quarto and Pixi

Source: https://context7.com/deltares/imod-documentation/llms.txt

This section outlines the process of building a documentation website using Quarto and Pixi, including automatic API documentation generation from Python docstrings. It covers installation, API building, local preview, PDF rendering, and checking Quarto installation. Dependencies are managed via Pixi.

```bash
# Install pixi and pre-commit hooks
pixi run install-pre-commit

# Build API documentation from Python docstrings
pixi run quartodoc-build

# Preview documentation website locally
pixi run docs

# Render PDF manual
quarto render docs --profile manual --to pdf

# Check Quarto installation and dependencies
pixi run quarto-check

# Render to HTML with execution
pixi run quarto-render
```

--------------------------------

### Create Time Discretization for Simulation

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Sets up the time discretization for the Modflow 6 simulation, defining the simulation period and frequency. This example models a period of 2 days.

```python
freq = "D"
times = pd.date_range(start="1/1/1971", end="1/3/1971", freq=freq)

simulation.create_time_discretization(additional_times=times)
```

--------------------------------

### Initialize Modflow 6 Simulation and Solver Settings

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Creates a Modflow 6 simulation object and attaches the configured groundwater flow model. It also defines the solver settings using a simple preset suitable for the example.

```python
simulation = mf6.Modflow6Simulation("test")
simulation["GWF_1"] = gwf_model

# Define solver settings, we'll use a preset that is sufficient for this example.
simulation["solver"] = mf6.SolutionPresetSimple(modelnames=["GWF_1"])
```

--------------------------------

### Clone and Build iMOD Documentation

Source: https://github.com/deltares/imod-documentation/blob/main/README.md

This snippet demonstrates the command-line steps to clone the iMOD documentation repository, install pre-commit hooks, and build the documentation locally using pixi and Quarto.

```bash
git clone  https://github.com/Deltares/iMOD-Documentation.git
cd iMOD-Documentation
pixi run install-pre-commit
pixi run docs
```

--------------------------------

### Define Scaling Factors (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Initializes scaling factors for soil moisture, hydraulic conductivity, and pressure head to 1.0. It also sets a uniform depth for the perched water table.

```python
msw_model["scaling"] = msw.ScalingFactors(
    scale_soil_moisture=xr.full_like(area, 1.0),
    scale_hydraulic_conductivity=xr.full_like(area, 1.0),
    scale_pressure_head=xr.full_like(area, 1.0),
    depth_perched_water_table=xr.full_like(msw_grid, 1.0),
)
```

--------------------------------

### Define Meteorology Grids (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Creates DataArrays for precipitation and evapotranspiration, both having a time dimension. Evapotranspiration is calculated as 1.5 times the precipitation.

```python
precipitation = msw_grid.expand_dims(time=times[:-1])
evapotranspiration = precipitation * 1.5

msw_model["meteo_grid"] = msw.MeteoGrid(precipitation, evapotranspiration)
msw_model["mapping_prec"] = msw.PrecipitationMapping(precipitation)
msw_model["mapping_evt"] = msw.EvapotranspirationMapping(evapotranspiration)
```

--------------------------------

### Configure Ponding Parameters (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Sets the parameters related to surface ponding, including ponding depth, runon resistance, and runoff resistance. All parameters are initialized to uniform values across the area grid.

```python
msw_model["ponding"] = msw.Ponding(
    ponding_depth=xr.full_like(area, 0.0),
    runon_resistance=xr.full_like(area, 1.0),
    runoff_resistance=xr.full_like(area, 1.0),
)
```

--------------------------------

### Initialize MetaSWAP Model

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Instantiates the MetaSwapModel from the imod.msw module. A critical step is providing the correct path to MetaSWAP's soil physical database, which is essential for its operation.

```python
msw_model = msw.MetaSwapModel(unsaturated_database="./path/to/unsaturated/database")
```

--------------------------------

### Set Rootzone Depth and Surface Elevation (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Initializes grids for rootzone depth and surface elevation. Rootzone depth is set uniformly across all subunits, while surface elevation is defined for each grid cell.

```python
rootzone_depth = xr.full_like(area, 0.5)
surface_elevation = xr.full_like(msw_grid, 2.0)

msw_model["grid"] = msw.GridData(
    area=area,
    landuse=landuse,
    rootzone_depth=rootzone_depth,
    surface_elevation=surface_elevation,
    soil_physical_unit=slt,
    active=active,
)
```

--------------------------------

### Create iMOD Recharge Package

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Initializes the iMOD Python Recharge package with the prepared recharge rate data. This package will be used in the MODFLOW 6 simulation setup.

```python
recharge_pkg_coarse = imod.mf6.Recharge(rate=recharge_ss_yearly)

recharge_pkg_coarse
```

--------------------------------

### Set Infiltration Parameters (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Configures infiltration parameters, setting infiltration capacity to 1.0 and resistances (downward, upward, bottom) to -9999.0 to ignore them. An extra storage coefficient is also set.

```python
msw_model["infiltration"] = msw.Infiltration(
    infiltration_capacity=xr.full_like(area, 1.0),
    downward_resistance=xr.full_like(msw_grid, -9999.0),
    upward_resistance=xr.full_like(msw_grid, -9999.0),
    bottom_resistance=xr.full_like(msw_grid, -9999.0),
    extra_storage_coefficient=xr.full_like(msw_grid, 0.1),
)
```

--------------------------------

### Assign Initial Conditions

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Netherlands_mesh.ipynb

Assigns the InitialConditions package to the groundwater flow model. This sets the starting head values for the simulation, often based on a drainage level.

```python
gwf_model["ic"] = imod.mf6.InitialConditions(start=drainage_lvl)
```

--------------------------------

### Configure Specific Storage and Initial Conditions

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Sets the SpecificStorage package, defining aquifer storage properties and marking cells as non-convertible for MetaSWAP coupling. It also sets the InitialConditions package for the Modflow 6 model.

```python
gwf_model["sto"] = mf6.SpecificStorage(
    specific_storage=1e-3, specific_yield=0.0, transient=True, convertible=0
)

gwf_model["ic"] = mf6.InitialConditions(start=0.5)
```

--------------------------------

### Import necessary libraries with imod

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Imports essential Python libraries for numerical operations, data handling, and the imod package for groundwater and unsaturated flow modeling. These libraries are fundamental for setting up and running the coupled models.

```python
import numpy as np
import pandas as pd
import xarray as xr

import primod
import imod
from imod import mf6, msw
```

--------------------------------

### Define Solver Settings

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Netherlands_mesh.ipynb

Specifies the solver settings for the MODFLOW 6 simulation using a preset. This example uses the 'Moderate' preset for the 'gwf' model.

```python
simulation["solver"] = imod.mf6.SolutionPresetModerate(modelnames=["gwf"])
```

--------------------------------

### Initialize Modflow 6 Groundwater Flow Model

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Creates an instance of the Modflow 6 GroundwaterFlowModel. This is the base object for defining the groundwater flow simulation, to which various packages and settings will be added.

```python
gwf_model = mf6.GroundwaterFlowModel()
```

--------------------------------

### Define Soil Cover for Crops

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This Python code defines the soil cover for different crops over a year using xarray. It starts with a grid of zeros and then applies a step function to represent periods of soil cover. The resulting DataArray can be used to calculate leaf area index.

```python
soil_cover = xr.DataArray(
    data=np.zeros(day_of_year.shape + vegetation_index.shape),
    coords=coords,
    dims=("day_of_year", "vegetation_index"),
)

soil_cover[132:254, :] = 1.0

soil_cover.sel(vegetation_index=1).plot()
```

--------------------------------

### Map MODFLOW 6 Discretization and Wells for MetaSWAP

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This snippet demonstrates how to map the MODFLOW 6 structured discretization and well package to the MetaSWAP model for internal use. It utilizes `mf6.StructuredDiscretization` and `mf6.Well` which are provided to `mf6.CouplerMapping`.

```python
msw_model["mod2svat"] = msw.CouplerMapping(
    modflow_dis=gwf_model["dis"],
    well=gwf_model["wells_msw"]
)
```

--------------------------------

### Configure IMOD Output Controls

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This Python snippet demonstrates the configuration of output controls for the IMOD model. It sets up an `IdfMapping` for area and a default value, a `VariableOutputControl`, and a `TimeOutputControl` specifying the output times. These settings determine what data is written and when during the simulation.

```python
msw_model["oc_idf"] = msw.IdfMapping(area, -9999.0)
msw_model["oc_var"] = msw.VariableOutputControl()
msw_model["oc_time"] = msw.TimeOutputControl(time=times)
```

--------------------------------

### Initialize Crop Growth Coordinates

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This Python snippet sets up the coordinate arrays for crop growth calculations. It defines the range of days in a year and the vegetation indices, then combines them into a dictionary for use in xarray DataArrays. This establishes the dimensions for crop-specific data.

```python
day_of_year = np.arange(1, 367)
vegetation_index = np.arange(1, 4)

coords = {"day_of_year": day_of_year, "vegetation_index": vegetation_index}
```

--------------------------------

### Define Dummy Recharge Package for MetaSWAP Coupling

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Creates a dummy Recharge package required by the iMOD Coupler when MetaSWAP is enabled. This package is configured with no recharge at the ditches' locations, ensuring proper matrix allocation during coupled simulation.

```python
recharge = xr.zeros_like(idomain.sel(layer=1), dtype=float)
recharge[:, 0] = np.nan
recharge[:, -1] = np.nan

gwf_model["rch_msw"] = mf6.Recharge(recharge)
```

--------------------------------

### Import iMOD Python, xarray, and xugrid Packages

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Netherlands_mesh.ipynb

Imports the necessary Python libraries for working with iMOD, xarray datasets, and unstructured grids. These are essential for all subsequent operations in the example.

```python
import imod
import xarray as xr
import xugrid as xu
```

--------------------------------

### Define Landuse Grid (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Creates a DataArray representing landuse classes for each cell and subunit. It initializes all cells with landuse class 1 and then assigns landuse class 2 to the cells corresponding to the second subunit.

```python
landuse = xr.full_like(area, 1, dtype=np.int16)
landuse[1, :, :] = 2

landuse
```

--------------------------------

### Define Landuse Names and Vegetation Indices

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This Python snippet initializes lists for vegetation indices and their corresponding names, then creates xarray DataArrays for landuse names and vegetation indices. These DataArrays are structured with coordinates, forming the basis for landuse parameter mapping.

```python
vegetation_index = [1, 2, 3]
names = ["grassland", "maize", "potatoes"]

landuse_index = [1, 2, 3]
coords = {"landuse_index": landuse_index}

landuse_names = xr.DataArray(data=names, coords=coords, dims=("landuse_index",))
vegetation_index_da = xr.DataArray(
    data=vegetation_index, coords=coords, dims=("landuse_index",)
)
```

--------------------------------

### Configure LanduseOptions in IMOD

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This Python code configures the `LanduseOptions` class in the `msw` library. It initializes a DataArray of ones for broadcasting parameters and then assigns various crop-specific stress and seasonal parameters, including Jarvis stress, Feddes parameters, and sprinkling thresholds.

```python
lu = xr.ones_like(vegetation_index_da, dtype=float)

msw_model["landuse_options"] = msw.LanduseOptions(
    landuse_name=landuse_names,
    vegetation_index=vegetation_index_da,
    jarvis_o2_stress=xr.ones_like(lu),
    jarvis_drought_stress=xr.ones_like(lu),
    feddes_p1=xr.full_like(lu, 99.0),
    feddes_p2=xr.full_like(lu, 99.0),
    feddes_p3h=lu * [-2.0, -4.0, -3.0],
    feddes_p3l=lu * [-8.0, -5.0, -5.0],
    feddes_p4=lu * [-80.0, -100.0, -100.0],
    feddes_t3h=xr.full_like(lu, 5.0),
    feddes_t3l=xr.full_like(lu, 1.0),
    threshold_sprinkling=lu * [-8.0, -5.0, -5.0],
    fraction_evaporated_sprinkling=xr.full_like(lu, 0.05),
    gift=xr.full_like(lu, 20.0),
    gift_duration=xr.full_like(lu, 0.25),
    rotational_period=lu * [10, 7, 7],
    start_sprinkling_season=lu * [120, 180, 150],
    end_sprinkling_season=lu * [230, 230, 240],
    interception_option=xr.ones_like(lu, dtype=int),
    interception_capacity_per_LAI=xr.zeros_like(lu),
    interception_intercept=xr.ones_like(lu),
)
```

--------------------------------

### Write Coupled Models using iMOD Coupler

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This code demonstrates how to write the coupled MetaSWAP and MODFLOW 6 models using the iMOD Coupler. It requires specifying the output directory and the paths to the MODFLOW 6 library, MetaSWAP library, and the directory containing MetaSWAP's dependent libraries.

```python
metamod_dir = imod.util.temporary_directory()
mf6_dll = "./path/to/mf6.dll"
metaswap_dll = "./path/to/metaswap.dll"
metaswap_dll_dependency_dir = "./path/to/metaswap/dll/dependency/directory"

metamod.write(metamod_dir, mf6_dll, metaswap_dll, metaswap_dll_dependency_dir)
```

--------------------------------

### Configure Node Property Flow Package (Hydraulic Conductivity)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Assigns the NodePropertyFlow package to the Modflow 6 model, specifying hydraulic conductivity (k and k33) for each layer. This configuration sets up the hydrogeological properties, including an aquitard in the middle layer.

```python
k = xr.DataArray([10.0, 0.1, 10.0], {"layer": layer}, ("layer",))
k33 = xr.DataArray([1.0, 0.01, 1.0], {"layer": layer}, ("layer",))
gwf_model["npf"] = mf6.NodePropertyFlow(
    icelltype=0,
    k=k,
    k33=k33,
    save_flows=True,
)
```

--------------------------------

### Initialize MetaMod for Coupled Model Simulation

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Initializes the `MetaMod` class, which manages the coupling between the MetaSWAP model and the MODFLOW 6 simulation. It takes the MetaSWAP model object, the MODFLOW 6 simulation object, and a list of coupling configurations as input.

```python
metamod = primod.MetaMod(
    msw_model=msw_model,
    mf6_simulation=simulation,
    coupling_list=[driver_coupling]
)
```

--------------------------------

### Define Dummy Well Package for MetaSWAP Coupling

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Sets up a dummy Well package, necessary for MetaSWAP's sprinkling functionality when coupled with Modflow 6. The input data is structured as long tables, specifying well layer, row, column, and rate.

```python
wel_layer = 3

ix = np.tile(np.arange(ncol) + 1, nrow)
iy = np.repeat(np.arange(nrow) + 1, ncol)
rate = np.zeros(ix.shape)
layer = np.full_like(ix, wel_layer)

gwf_model["wells_msw"] = mf6.WellDisStructured(
    layer=layer, row=iy, column=ix, rate=rate
)
```

--------------------------------

### Visualize Time Series in QGIS using IPF Files

Source: https://context7.com/deltares/imod-documentation/llms.txt

This guide explains how to load and visualize IPF (iMOD Point File) time series data within the QGIS iMOD plugin. It details the steps for selecting layers, choosing ID columns and variables, enabling interactive updates, and using selection tools. Data can be exported in various formats like PNG, JPG, SVG, and CSV.

```python
# After loading IPF file via QGIS plugin:
# 1. Select IPF layer in Time Series widget dropdown
# 2. Choose ID column for legend labels
# 3. Select variable(s) to plot
# 4. Enable "update on selection" checkbox
# 5. Draw selection box on map (left-click and drag)
# 6. SHIFT + click to add more points
# 7. CTRL + ALT + A to deselect all

# Export time series data
# Click "Export" button and choose format:
# - .png or .jpg for images
# - .svg for vector graphics
# - .csv for further analysis

```

--------------------------------

### Access River Package Data from Groundwater Model

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Illustrates how to retrieve a specific model package, such as the river package ('riv'), from a `GroundwaterFlowModel` object. The package object itself is accessed like a dictionary entry. To get the actual data, the `.dataset` attribute of the package object must be used, which returns an xarray Dataset.

```python
riv_pkg = gwf_model["riv"]
```

--------------------------------

### Load MODFLOW 6 Simulation from TOML

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Loads a MODFLOW 6 simulation configuration from a TOML file using iMOD Python. This function creates a `Modflow6Simulation` object, which acts like a Python dictionary. It utilizes `imod.data.hondsrug_simulation` to load example data and `imod.util.temporary_directory` for temporary storage. The code also demonstrates disabling the XT3D option for the NPF package.

```python
tmpdir = imod.util.temporary_directory()

simulation = imod.data.hondsrug_simulation(tmpdir / "hondsrug_saved")

# Ensure xt3d is not used
simulation["GWF"]["npf"].set_xt3d_option(is_xt3d_used=False, is_rhs=False)
```

--------------------------------

### Set Up MetaMod Driver Coupling for MODFLOW 6 and MetaSWAP

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This snippet defines the coupling between a MODFLOW 6 simulation and MetaSWAP. It specifies the MODFLOW 6 model name, the recharge package key, and the well package key, which are essential for iMOD Python to identify the connection points.

```python
driver_coupling = primod.MetaModDriverCoupling(
    mf6_model="GWF_1",
    mf6_recharge_package="rch_msw",
    mf6_wel_package="wells_msw"
)
```

--------------------------------

### Define Modflow 6 Structured Discretization

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Sets up the spatial discretization for the Modflow 6 model, defining the grid dimensions, cell sizes, and layer elevations. This includes creating coordinate arrays and an idomain DataArray to represent the active cells.

```python
shape = nlay, nrow, ncol = 3, 9, 9

dx = 10.0
dy = -10.0
dz = np.array([1.0, 2.0, 10.0])
xmin = 0.0
xmax = dx * ncol
ymin = 0.0
ymax = abs(dy) * nrow
dims = ("layer", "y", "x")

layer = np.arange(1, nlay + 1)
y = np.arange(ymax, ymin, dy) + 0.5 * dy
x = np.arange(xmin, xmax, dx) + 0.5 * dx
coords = {"layer": layer, "y": y, "x": x}

idomain = xr.DataArray(np.ones(shape, dtype=int), coords=coords, dims=dims)

top = 0.0
bottom = top - xr.DataArray(
    np.cumsum(layer * dz), coords={"layer": layer}, dims="layer"
)

gwf_model["dis"] = mf6.StructuredDiscretization(idomain=idomain, top=top, bottom=bottom)
```

--------------------------------

### Read and Write IPF Files with iMOD Python

Source: https://context7.com/deltares/imod-documentation/llms.txt

This example demonstrates how to read and write IPF (iMOD Point File) format for point data, including time series. It uses the imod.formats.ipf module to read data into a GeoDataFrame with associated time series and to write this data back to an IPF file. Dependencies include imod.

```python
import imod

# Read IPF file (flexible parser for various formats)
ipf_data = imod.formats.ipf.read("timeseries.ipf")

# Inspect data structure
print(ipf_data)  # Returns geodataframe with associated time series

# Write IPF file (standardized format for QGIS compatibility)
imod.formats.ipf.write(
    "output_timeseries.ipf",
    ipf_data,
    delimiter=","
)

```

--------------------------------

### Configure MetaSWAP Sprinkling Package with MODFLOW 6 Wells

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This code configures the MetaSWAP sprinkling package, which requires the MODFLOW 6 wells. It sets maximum abstraction values for groundwater and surface water, linking the `well` parameter to the MODFLOW 6 well package.

```python
msw_model["sprinkling"] = msw.Sprinkling(
    max_abstraction_groundwater=xr.full_like(msw_grid, 100.0),
    max_abstraction_surfacewater=xr.full_like(msw_grid, 100.0),
    well=gwf_model["wells_msw"]
)
```

--------------------------------

### Select River Stage Data for a Specific Layer

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Shows how to extract specific data variables and subset them by dimension from an xarray Dataset representing a model package. In this example, it selects the 'stage' variable from the river package and then filters it to include only data for layer 1 using the `.sel()` method. The result is a 2D xarray Dataset.

```python
riv_stage_layer_1 = riv_pkg["stage"].sel(layer=1)
```

--------------------------------

### Define Modflow 6 Constant Head Boundary Conditions

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Implements Constant Head boundary conditions at the left and right edges of the model grid, simulating ditches with fixed water levels. The boundary conditions are defined using an xarray DataArray.

```python
head = xr.full_like(idomain, np.nan, dtype=float)
head[0, :, 0] = -1.0
head[0, :, -1] = -1.0

gwf_model["chd"] = mf6.ConstantHead(
    head, print_input=True, print_flows=True, save_flows=True
)
```

--------------------------------

### Configure AnnualCropFactors in IMOD

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This Python code configures the `AnnualCropFactors` package within the `msw` library. It assigns the previously defined `soil_cover`, `leaf_area_index`, and `vegetation_factor` DataArrays, along with default values for interception capacity and other factors. This setup is crucial for simulating crop-specific water use.

```python
msw_model["crop_factors"] = msw.AnnualCropFactors(
    soil_cover=soil_cover,
    leaf_area_index=leaf_area_index,
    interception_capacity=xr.zeros_like(soil_cover),
    vegetation_factor=vegetation_factor,
    interception_factor=xr.ones_like(soil_cover),
    bare_soil_factor=xr.ones_like(soil_cover),
    ponding_factor=xr.ones_like(soil_cover),
)
```

--------------------------------

### Clone Repository and Create Branch

Source: https://github.com/deltares/imod-documentation/blob/main/CONTRIBUTING.md

This batch script demonstrates how to clone the iMOD-Documentation repository, create a new branch for your changes, add modified files, commit them with a message, and push the branch to the origin. Remember to replace placeholders like '<my_branch_name>' and '<my commit message>' with your specific details.

```batch
git clone https://github.com/Deltares/iMOD-Documentation.git
git checkout -b <my_branch_name>

rem Change some files here 

git add <files>
git commit -m "<my commit message>"

git push origin <my_branch_name>
```

--------------------------------

### Render iMOD Documentation to PDF

Source: https://github.com/deltares/imod-documentation/blob/main/README.md

This command renders the iMOD documentation to a PDF format using Quarto. It utilizes the 'manual' profile and ensures any previously built documentation is deleted before rendering.

```bash
quarto render docs --profile manual --to pdf
```

--------------------------------

### Write and Run iMOD Simulation

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Writes the simulation configuration to a specified directory and then runs the simulation. This is a standard procedure for executing an iMOD model after modifications.

```python
modeldir = prj_dir / "new_recharge"
simulation.write(modeldir)
simulation.run(mf6_path)
```

--------------------------------

### Run iMOD Simulation with Extraction Well

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Writes the simulation files to a specified directory and then runs the iMOD simulation. This step includes configuring the output directory and executing the simulation using the mf6_path. It's a crucial step before analyzing results with the added well.

```python
# Step 1: run the model
modeldir = prj_dir / "extraction"
simulation.write(modeldir)
simulation.run(mf6_path)

```

--------------------------------

### Running Partitioned Simulations with iMOD

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

This snippet executes the partitioned simulations. It first specifies the directory where the partitioned models will be written using `modeldir`. Then, `simulation_split.write()` saves the partitioned simulation files, and `simulation_split.run()` executes each partition using the specified MODFLOW 6 path. The output is the execution of multiple groundwater flow models.

```python
modeldir = prj_dir / "partitioning"
simulation_split.write(modeldir)
simulation_split.run(mf6_path)
```

--------------------------------

### Define Model Idomain for Discretization

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Netherlands_mesh.ipynb

Prepares the model's idomain, specifying active [1], inactive [0], or vertical pass-through [-1] cells. This example sets all cells to be active.

```python
idomain = template.copy()
```

--------------------------------

### Import iMOD Python Library

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

This code snippet demonstrates how to import the iMOD Python library. After importing, all iMOD Python functions and classes become available for use. This is typically the first step in any iMOD Python script.

```python
import imod
```

--------------------------------

### Select Heads at Final Timestep using isel

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Selects data based on integer index. This example extracts head data for the final time step (index -1) from an xarray Dataset. It uses the `.isel()` method, which is suitable for index-based selection.

```python
heads_end = calculated_heads.isel(time=-1)

heads_end

```

--------------------------------

### Partitioning a Simulation with iMOD

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

This code partitions an existing simulation object based on provided submodel labels. The `split` function creates a new simulation object where the original simulation is divided according to the spatial distribution of the submodel labels. The output is a partitioned simulation object ready for execution.

```python
simulation_split = simulation.split(submodel_labels)

simulation_split
```

--------------------------------

### Plot River Stage with iMOD Visualize (Viridis Colormap)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Demonstrates plotting river stage data with a different colormap ('viridis'). It utilizes the iMOD plot_map function with predefined levels and the specified colormap.

```python
colors = "viridis"
imod.visualize.plot_map(riv_stage_layer_1, colors, levels)
```

--------------------------------

### Opening and Merging Partitioned Heads with iMOD

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

This code opens the head results from the partitioned simulations and automatically merges them into a single grid. The `open_head` function handles the merging process, allowing for direct comparison with the original simulation grid. The output is a DataArray containing the merged head data.

```python
head_split = simulation_split.open_head()["head"]

head_split
```

--------------------------------

### Configure Modflow 6 Output Control

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Sets the OutputControl package for the Modflow 6 model to manage the saving of simulation results, specifically the last time step's head and budget data. This is crucial for post-processing and analysis.

```python
gwf_model["oc"] = mf6.OutputControl(save_head="last", save_budget="last")
```

--------------------------------

### Define Active MetaSWAP Grid (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Creates a boolean grid to specify active MetaSWAP cells, excluding cells occupied by ditches in Modflow 6. It initializes the grid based on the msw_grid and sets the first and last columns to False.

```python
active = msw_grid.astype(bool)
active[..., 0] = False
active[..., -1] = False

active
```

--------------------------------

### Open and Load Simulation Results

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Netherlands_mesh.ipynb

Opens the simulation results, specifically the head, loads them, and selects the last time step and first layer for analysis and plotting.

```python
head = simulation.open_head().load().isel(time=-1, layer=0)
# calculate ground water below surface level.
groundwater_depth = elevation - head
# Plot
groundwater_depth.ugrid.plot()
```

--------------------------------

### Couple MODFLOW 6 with MetaSWAP using iMOD Coupler and Xmipy

Source: https://context7.com/deltares/imod-documentation/llms.txt

This snippet shows how to configure and run coupled hydrological models (MODFLOW 6 and MetaSWAP) using iMOD Coupler and the BMI interface via xmipy. It covers initializing MODFLOW 6, updating the model during runtime, exchanging data, and finalizing the simulation. It requires downloading iMOD Coupler and having a MODFLOW 6 model file.

```python
# Installation of iMOD Coupler
# Download from: https://github.com/Deltares/imod_coupler/releases/latest/

# Using xmipy for runtime control
import xmipy

# Initialize MODFLOW 6 model with BMI
mf6_bmi = xmipy.XmiWrapper(lib_path="path/to/mf6.dll")
mf6_bmi.initialize("mfsim.nam")

# Update model during runtime
timestep = 1.0
while mf6_bmi.get_current_time() < mf6_bmi.get_end_time():
    mf6_bmi.update()

    # Get head values
    heads = mf6_bmi.get_value("X")

    # Exchange data with other models (e.g., MetaSWAP)
    # Modify values if needed
    mf6_bmi.set_value("X", modified_heads)

# Finalizemf6_bmi.finalize()

```

--------------------------------

### Define Cell Area per Subunit (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Calculates and assigns the area for each subunit within each model cell. It determines the total cell area and distributes it equally among the defined subunits, creating an 'area' DataArray that can be broadcast across subunits.

```python
subunit = [0, 1]

total_cell_area = abs(dx * dy)
equal_area_per_subunit = total_cell_area / len(subunit)

total_cell_area

# Create a full grid equal to the msw_grid. And expand_dims() to broadcast this
# grid along a new dimension, named "subunit"
area = (
    xr.full_like(msw_grid, equal_area_per_subunit, dtype=float)
    .expand_dims(subunit=subunit)
    .copy()  # expand_dims creates a view, so copy it to allow setting values.
)

# To the left we only have subunit 0
area[0, :, :3] = total_cell_area
area[1, :, :3] = np.nan

# To the right we only have subunit 1
area[0, :, -3:] = np.nan
area[1, :, -3:] = total_cell_area

area
```

--------------------------------

### Partition iMOD Simulation for Parallel Computation

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Partitions an iMOD simulation into multiple sub-models for parallel computation. It utilizes the 'get_label_array' function to generate subdomain labels based on the simulation grid and a specified number of partitions.

```python
from imod.mf6.multimodel.partition_generator import get_label_array

npartitions = 4
submodel_labels = get_label_array(simulation, npartitions=npartitions)
```

--------------------------------

### Access Groundwater Flow Model from Simulation Object

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Demonstrates how to access a specific groundwater flow model (e.g., 'GWF') from a loaded `Modflow6Simulation` object. The simulation object is treated as a dictionary, allowing direct access to its constituent models by name. The returned `GroundwaterFlowModel` object also behaves like a dictionary, mapping package names to their respective objects.

```python
gwf_model = simulation["GWF"]
```

--------------------------------

### Resample Recharge to Yearly Timestep

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Resamples the daily recharge data to a yearly timestep using xarray's resample function. The 'A' alias denotes annual frequency, and 'label="left"' sets the year label to its starting point, aligning with iMOD Python's stress period definition.

```python
recharge_yearly = recharge_daily.resample(time="A", label="left").mean()

recharge_yearly
```

--------------------------------

### Write and Run Partitioned Simulation

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Netherlands_mesh.ipynb

Writes the split simulation to a directory and then runs it. iMOD Python handles the merging of results from the partitioned models after execution.

```python
modeldir = tmpdir / "partitioned"
simulation_split.write(modeldir)
simulation_split.run(mf6_path)
```

--------------------------------

### Create and Plot Custom iMOD Cross-Section

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Demonstrates creating a cross-section plot using a custom-defined line geometry from Shapely. This allows for flexible cross-section slicing through the model domain.

```python
from shapely.geometry import LineString

geometry = LineString([[244000,563000],[244000 ,559000]])
cross_k = imod.select.cross_section_linestring(geology, geometry)
imod.visualize.cross_section(cross_k, colors, levels)
```

--------------------------------

### Open Calculated Heads from iMOD Simulation

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Opens the calculated heads from an iMOD simulation. This method automatically locates and loads the necessary output files, returning an xarray Dataset containing the head data. Dependencies include the simulation object and potentially the iMOD Python library.

```python
calculated_heads = simulation.open_head()

```

--------------------------------

### Define Soil Type Grid (Python)

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

Defines the soil type grid using xarray's `where()` method for efficient conditional assignment. All cells are initially set to soil type 1, and then cells in the right half of the grid are updated to soil type 2.

```python
slt = xr.full_like(msw_grid, 1, dtype=np.int16)
# Set all cells on the right half to 2.
slt = slt.where((slt.x < (xmax / 2)), 2)

slt
```

--------------------------------

### Create Time Discretization from Boundary Conditions (Python)

Source: https://context7.com/deltares/imod-documentation/llms.txt

Automatically generates a consistent time discretization for a MODFLOW 6 simulation based on the times assigned to boundary conditions. It allows specifying additional simulation end times. Dependency is 'imod'.

```python
import imod

# Create simulation and add packages
simulation = imod.mf6.Modflow6Simulation("example")
simulation["gwf"] = imod.mf6.GroundwaterFlowModel()
simulation["gwf"]["dis"] = dis_pkg
simulation["gwf"]["wel"] = wel_pkg

# Create time discretization from all package times
# Specify simulation end time as additional_times
simulation.create_time_discretization(additional_times=["2000-01-07"])

# Inspect generated time discretization
print(simulation["time_discretization"].dataset)
```

--------------------------------

### Add Specific Storage Package

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Netherlands_mesh.ipynb

Adds the SpecificStorage package to the groundwater flow model. This package defines parameters related to storage, such as specific storage, specific yield, and transient flow conditions.

```python
gwf_model["sto"] = imod.mf6.SpecificStorage(
    specific_storage=1e-5,
    specific_yield=0.01,
    transient=False,
    convertible=1
)
```

--------------------------------

### Calculate Leaf Area Index

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This Python snippet calculates the leaf area index (LAI) for crops by tripling the soil cover values. LAI is a crucial parameter in crop growth models, influencing evapotranspiration and other physiological processes. The result is an xarray DataArray with dimensions for day of year and vegetation index.

```python
leaf_area_index = soil_cover * 3
```

--------------------------------

### Define Vegetation Factors for Crops

Source: https://github.com/deltares/imod-documentation/blob/main/docs/coupler_metamod_example.ipynb

This Python code defines vegetation factors that adjust potential evapotranspiration based on crop type and day of the year. It initializes an xarray DataArray with zeros and then assigns specific factor values for different growth stages of grass, maize, and potatoes. The grass vegetation factor is fixed to 1.0 as it serves as the reference.

```python
vegetation_names = ["grass", "maize", "potatoes"]

vegetation_factor = xr.zeros_like(soil_cover)

vegetation_factor[120:132, :] = [1.0, 0.5, 0.0]
vegetation_factor[132:142, :] = [1.0, 0.7, 0.7]
vegetation_factor[142:152, :] = [1.0, 0.8, 0.9]
vegetation_factor[152:162, :] = [1.0, 0.9, 1.0]
vegetation_factor[162:172, :] = [1.0, 1.0, 1.2]
vegetation_factor[172:182, :] = [1.0, 1.2, 1.2]
vegetation_factor[182:192, :] = [1.0, 1.3, 1.2]
vegetation_factor[192:244, :] = [1.0, 1.2, 1.1]
vegetation_factor[244:254, :] = [1.0, 1.2, 0.7]
vegetation_factor[254:283, :] = [1.0, 1.2, 0.0]

# Since grass is the reference crop, force all grass to 1.0
vegetation_factor[:, 0] = 1.0


# Assign vegetation names for the plot
vegetation_factor.assign_coords(
    vegetation_names=("vegetation_index", vegetation_names)
).plot.line(x="day_of_year", hue="vegetation_names")
```

--------------------------------

### Plot River Stage with iMOD Visualize and Basemap

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Enhances river stage plotting by adding a background map using the 'contextily' package. It requires importing 'contextily', defining a basemap provider, and passing it to the plot_map function.

```python
import contextily as ctx

background_map = ctx.providers["OpenStreetMap"]["Mapnik"]
imod.visualize.plot_map(riv_stage_layer_1, colors, levels, basemap=background_map)
```

--------------------------------

### Write MODFLOW 6 Simulation Input Files

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Writes the MODFLOW 6 simulation model object into native MODFLOW 6 input files, including the mfsim.nam and other package files, to a specified directory.

```python
# Write the MODFLOW 6 simulation files from the model object
simulation.write(modeldir)
```

--------------------------------

### Add Output Control Package

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Netherlands_mesh.ipynb

Adds the OutputControl package to the groundwater flow model, specifying which outputs (e.g., head, budget) to save and when (e.g., at the last time step).

```python
gwf_model["oc"] = imod.mf6.OutputControl(save_head="last", save_budget="last")
```

--------------------------------

### Create Geology Dataset for iMOD Cross-Section

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Creates a geology dataset required for plotting cross-sections in iMOD Python. It extracts permeability ('k') from the 'npf' package and layer top/bottom information from the 'dis' package. Utilizes xarray's 'shift' for efficient layer top calculation.

```python
geology = gwf_model["npf"]["k"]

DEM = gwf_model["dis"]["top"]
geology["top"] = geology["bottom"].shift(layer=1).fillna(DEM)
```

--------------------------------

### Create an iMOD Extraction Well Object

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Defines an extraction well with specific properties such as screen top, screen bottom, location (x, y), and pumping rate. This creates an imod.mf6.Well object that can be added to the simulation. Coordinates and rates are provided as lists.

```python
x_well = 246000.0
y_well = 560000.0
# create a well object
well = imod.mf6.Well(
    screen_top=[-13.0],
    screen_bottom=[-50.0],
    y=[y_well],
    x=[x_well],
    rate=[-60000.0],
)

```

--------------------------------

### Plot Bottom Elevation for Layer 3 with iMOD Visualize

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Placeholder for plotting the bottom elevation data for layer 3 using iMOD Python. Users should replace the comment with specific commands to select and plot the data.

```python
#
# type here you python commands to select and plot bottom data for layer 3
```

--------------------------------

### Plotting Partitioned Results with Overlays using iMOD

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

This snippet visualizes the difference between split and calculated heads, overlaid with submodel partitions. It first calculates the difference (`diff_split`) and selects a specific slice for plotting (`diff_split_end_first_aqf`). The `imod.visualize.plot_map` function is used with an `overlays` argument, specifying polygon geometry, edge color, and face color to distinguish partitions.

```python
overlays = [{"gdf": submodel_polygons, "edgecolor": "black", "facecolor": "None"}]

diff_split = head_split - calculated_heads

diff_split_end_first_aqf = (
    diff_split.isel(time=-1).sel(layer=slice(3, 5)).mean(dim="layer")
)

levels = np.arange(-10.0, 0.0, 1.0)

imod.visualize.plot_map(
    diff_split_end_first_aqf, "hot", levels, basemap=background_map, overlays=overlays
)
```

--------------------------------

### Open Calculated Heads After Adding Well

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Loads the head results from a simulation that includes an extraction well. Similar to the initial head opening, this retrieves the output data as an xarray Dataset, allowing for subsequent analysis of the well's impact.

```python
# Step 2: export the calculated heads

heads_with_well = simulation.open_head()

heads_with_well

```

--------------------------------

### Create MODFLOW 6 Result Directory

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Creates a directory for storing MODFLOW 6 simulation results using the 'pathlib' package. It defines the project directory and the specific model directory name.

```python
# Define the project directory where you want to write and run the Modflow 6
# simulation.
prj_dir = tmpdir / "results"

modeldir = prj_dir / "reference"
```

--------------------------------

### Add Well Package to Groundwater Model

Source: https://github.com/deltares/imod-documentation/blob/main/docs/tutorial_Hondsrug.ipynb

Adds a well package to an existing groundwater model ('GWF') within the iMOD simulation. This is done by assigning the created well object to a new key (e.g., 'wel') within the simulation's model dictionary.

```python
simulation["GWF"]["wel"] = well

```

--------------------------------

### Creating Time Discretization from Boundary Conditions

Source: https://context7.com/deltares/imod-documentation/llms.txt

Automatically generates consistent stress periods based on the times assigned to various boundary conditions within the simulation.

```APIDOC
## Creating Time Discretization from Boundary Conditions

### Description
Automatically generate consistent stress periods based on times assigned to boundary conditions.

### Method
(Implicitly demonstrated through Python code)

### Endpoint
(Not applicable, function-based)

### Parameters
(Implicitly handled by function arguments)

### Request Example
```python
import imod

# Assuming dis_pkg and wel_pkg are already defined
# simulation = imod.mf6.Modflow6Simulation("example")
# simulation["gwf"] = imod.mf6.GroundwaterFlowModel()
# simulation["gwf"]["dis"] = dis_pkg
# simulation["gwf"]["wel"] = wel_pkg

# Create time discretization from all package times
# Specify simulation end time as additional_times
# simulation.create_time_discretization(additional_times=["2000-01-07"])

# Inspect generated time discretization
# print(simulation["time_discretization"].dataset)
```

### Response
#### Success Response (200)
(Implicitly modifies the simulation object with a `time_discretization` dataset)

#### Response Example
```
# Example output from print(simulation["time_discretization"].dataset):
# <xarray.Dataset>
# Dimensions:        (time: 3)
# Coordinates:
#   * time            (time) object 2000-01-01 2000-01-03 2000-01-07
# Data variables:
#     stress_period   (time) int64 0 1 2
#     layer_group     (time) int64 0 0 0
```
```

--------------------------------

### Create Cross-Sections in QGIS with iMOD Plugin

Source: https://context7.com/deltares/imod-documentation/llms.txt

This section details the workflow for generating cross-sections from mesh and raster data using the QGIS iMOD plugin. It covers selecting line locations, styling datasets with customizable variables and layers, setting buffers for borehole selection, and plotting the cross-section. Export options include PNG, JPG, SVG, and CSV.

```python
# Cross-section workflow in QGIS plugin:
# 1. Click "Select location" button
# 2. Left-click to start drawing line, left-click again to end
# 3. Right-click to finish selection
# 4. Add datasets to styling table:
#    - Select dataset from dropdown
#    - Choose variable to plot
#    - Select layers (or keep default "all")
#    - Enable "As line(s)" for layer boundaries
# 5. Set search buffer for borehole selection
# 6. Click "Plot" to render cross-section

# Export cross-section
# Click "Export" button for .png, .jpg, .svg, or .csv output

```