imod.evaluate – Evaluate model output

imod.evaluate.facebudget(budgetzone, front=None, lower=None, right=None, netflow=True)[source]

Computes net face flow into a control volume, as defined by budgetzone.

Returns a three dimensional DataArray with in- and outgoing flows for every cell that is located on the edge of the control volume, thereby calculating the flow through the control surface of the control volume.

Front, lower, and right arguments refer to iMOD face flow budgets, in cubic meters per day. In terms of flow direction these are defined as:

  • front: positive with y (negative with row index)

  • lower: positive with layer (positive with layer index)

  • right: negative with x (negative with column index)

Only a single face flow has to be defined, for example, if only vertical fluxes (lower) are to be considered.

Note that you generally should not define a budgetzone that is only one cell wide. In that case, flow will both enter and leave the cell, resulting in a net flow of zero (given there are no boundary conditions present).

The output of this function defines ingoing flow as positive, and outgoing flow as negative. The output is a 3D array with net flows for every control surface cell. You can select specific parts of the control surface, for example only the east-facing side of the control volume. Please refer to the examples.

Parameters

  • budgetzone (xr.DataAray) – Binary array defining zone (1 or True) and outside of zone (0 or False). Dimensions must be exactly ("layer", "y", "x").

  • front (xr.DataArray of floats, optional) – Dimensions must be exactly ("layer", "y", "x").

  • lower (xr.DataArray of floats, optional) – Dimensions must be exactly ("layer", "y", "x").

  • right (xr.DataArray of floats, optional) – Dimensions must be exactly ("layer", "y", "x").

  • netflow (bool, optional) – Whether to split flows by direction (front, lower, right). True: sum all flows. False: return individual directions.

Returns

  • facebudget_front, facebudget_lower, face_budget_right (xr.DataArray of floats) – Only returned if netflow is False.

  • facebudget_net (xr.DataArray of floats) – Only returned if netflow is True.

Examples

Load the face flows, and select the last time using index selection:

>>> import imod
>>> lower = imod.idf.open("bdgflf*.idf").isel(time=-1)
>>> right = imod.idf.open("bdgfrf*.idf").isel(time=-1)
>>> front = imod.idf.open("bdgfff*.idf").isel(time=-1)

Define the zone of interest, e.g. via rasterizing a shapefile:

>>> import geopandas as gpd
>>> gdf = gpd.read_file("zone_of_interest.shp")
>>> zone2D = imod.prepare.rasterize(gdf, like=lower.isel(layer=0))

Broadcast it to three dimensions:

>>> zone = xr.full_like(flow, zone2D, dtype=np.bool)

Compute net flow through the (control) surface of the budget zone:

>>> flow = imod.evaluate.facebudget(
>>>     budgetzone=zone, front=front, lower=lower, right=right
>>> )

Or evaluate only horizontal fluxes:

>>> flow = imod.evaluate.facebudget(
>>>     budgetzone=zone, front=front, right=right
>>> )

Extract the net flow, only on the right side of the zone, for example as defined by x > 10000:

>>> netflow_right = flow.where(flow["x"] > 10_000.0).sum()

Extract the net flow for only the first four layers:

>>> netflow_layers = netflow_right.sel(layer=slice(1, 4)).sum()

Extract the net flow to the right of an arbitrary diagonal in x and y is simple using the equation for a straight line:

>>> netflow_right_of_diagonal = flow.where(
>>>    flow["y"] < (line_slope * flow["x"] + line_intercept)
>>> )

There are many ways to extract flows for a certain part of the zone of interest. The most flexible one with regards to the x and y dimensions is by drawing a vector file, rasterizing it, and using it to select with xarray.where().

To get the flows per direction, pass netflow=False.

>>> flowfront, flowlower, flowright = imod.evaluate.facebudget(
>>>    budgetzone=zone, front=front, lower=lower, right=right, netflow=False
>>> )
imod.evaluate.stability_constraint_advection(front, lower, right, top_bot, porosity=0.3, R=1.0)[source]

Computes advection stability constraint as applied in MT3D for adaptive timestepping (Zheng & Wang, 1999 p54). This function can be used to select which cells necessitate a small timestap, thereby slowing down calculations.

Front, lower, and right arguments refer to iMOD face flow budgets, in cubic meters per day. In terms of flow direction these are defined as:

  • front: positive with y (negative with row index)

  • lower: positive with layer (positive with layer index)

  • right: negative with x (negative with column index)

Returns the minimum timestep that is required to satisfy this constraint. The resulting dt xr.DataArray is the minimum timestep over all three directions, dt_xyz is an xr.Dataset containing minimum timesteps for the three directions separately.

Parameters

  • front (xr.DataArray of floats, optional) – Dimensions must be exactly ("layer", "y", "x").

  • lower (xr.DataArray of floats, optional) – Dimensions must be exactly ("layer", "y", "x").

  • right (xr.DataArray of floats, optional) – Dimensions must be exactly ("layer", "y", "x").

  • top_bot (xr.Dataset of floats, containing 'top', 'bot' and optionally) – ‘dz’ of layers. Dimensions must be exactly ("layer", "y", "x").

  • porosity (float or xr.DataArray of floats, optional (default 0.3)) – If xr.DataArray, dimensions must be exactly ("layer", "y", "x").

  • R (Retardation factor, optional (default)) – Only when sorption is a factor.

Returns

  • dt (xr.DataArray of floats)

  • dt_xyz (xr.Dataset of floats)

Examples

Load the face flows, and select the last time using index selection:

imod.evaluate.stability_constraint_wel(wel, top_bot, porosity=0.3, R=1.0)[source]

Computes sink/source stability constraint as applied in MT3D for adaptive timestepping (Zheng & Wang, 1999 p54). For the WEL package, a flux is known beforehand, so we can evaluate beforehand if a flux assigned to a cell will necessitate a small timestap, thereby slowing down calculations.

Returns a ipf DataFrame that includes a column for the specific discharge and resulting minimum timestep.

Parameters

  • wel (pd.DataFrame) – pd.DataFrame that defines a WEL package. Minimally includes x, y, layer and Q column.

  • top_bot (xr.Dataset of floats, containing 'top', 'bot' and optionally) – ‘dz’ of layers. Dimensions must be exactly ("layer", "y", "x").

  • porosity (float or xr.DataArray of floats, optional (default 0.3)) – If xr.DataArray, dimensions must be exactly ("layer", "y", "x").

  • R (Retardation factor, optional (default)) – Only when sorption is a factor.

Returns

wel – dt (minimum timestep) columns

Return type

pd.DataFrame containing addition qs (specific discharge) and