imod.flow - Create iMODFLOW model¶
Create iMODFLOW model.
With this module full iMODFLOW models can be written to disk.
Models in imod.flow
are created in a similar way to those created in imod.wq
.
The ImodflowModel object by default writes iMODFLOW projectfiles now
instead of the older iMODFLOW runfiles.
- class imod.flow.Bottom(bottom)[source]¶
Bases:
imod.flow.pkgbase.Package
The bottom of the aquifers
- Parameters
bottom (xr.DataArray of floats) – is the bottom elevation of model layers or Quasi-3d confining beds, with dimensions
("layer", "y", "x")
. The DataArray should at least include the layer dimension.
- dataset¶
- class imod.flow.Boundary(ibound)[source]¶
Bases:
imod.flow.pkgbase.Package
Specify the locations of active, inactive, and specified head in cells
- Parameters
ibound (xr.DataArray of ints) –
is the boundary variable with dimensions
("layer", "y", "x")
.If IBOUND(J,I,K) < 0, cell J,I,K has a constant head.
If IBOUND(J,I,K) = 0, cell J,I,K is inactive.
If IBOUND(J,I,K) > 0, cell J,I,K is active.
- dataset¶
- class imod.flow.ConstantHead(head)[source]¶
Bases:
imod.flow.pkgbase.BoundaryCondition
The Constant Head package. The Time-Variant Specified-Head package is used to simulate specified head boundaries that can change between stress periods.
- Parameters
head (xr.DataArray of floats) – is the head at the boundary
- dataset¶
- class imod.flow.Drain(conductance, elevation)[source]¶
Bases:
imod.flow.pkgbase.BoundaryCondition
The Drain package is used to simulate head-dependent flux boundaries. In the Drain package if the head in the cell falls below a certain threshold, the flux from the drain to the model cell drops to zero.
- Parameters
elevation (float or xr.DataArray of floats) – elevation of the drain, dims
("layer", "y", "x")
.conductance (float or xr.DataArray of floats) – is the conductance of the drain, dims
("layer", "y", "x")
.
- dataset¶
- class imod.flow.EvapoTranspiration(rate, top_elevation, extinction_depth)[source]¶
Bases:
imod.flow.pkgbase.TopBoundaryCondition
Recharge provides a fixed flux boundary condition to the top layer of the groundwater system. Note that unlike in iMOD-WQ, there is only the option in iMODFLOW to apply the recharge package to the top layer.
- Parameters
rate (float or xr.DataArray of floats) – evaporation rate in mm/day (NOTA BENE!), dims
("time", "y", "x")
.top_elevation (floats or xr.DataArray of floats) – Top elevation in m+MSL for maximal evapotranspiration strength.
extinction_depth (float or xr.Datarray of floats) – Depth [m] in which evapotranspiration strength reduced to zero.
- dataset¶
- class imod.flow.GeneralHeadBoundary(conductance, head)[source]¶
Bases:
imod.flow.pkgbase.BoundaryCondition
The General-Head Boundary package is used to simulate head-dependent flux boundaries. In the General-Head Boundary package the flux is always proportional to the difference in head.
- Parameters
head (float or xr.DataArray of floats) – head value for the GHB (BHEAD), dims
("layer", "y", "x")
.conductance (float or xr.DataArray of floats) – the conductance of the GHB (COND), dims
("layer", "y", "x")
.
- dataset¶
- class imod.flow.HorizontalAnisotropy(anisotropy_factor, anisotropy_angle)[source]¶
Bases:
imod.flow.pkgbase.Package
Horizontal anisotropy is a phenomenon in which the horizontal hydraulic conductivity is not equal along the x and y Cartesian axes. iMODFLOW can calculate this anisotropy based on a anisotropy factor and an anisotropy angle. iMODFLOW also accounts for the cross-terms in the horizontal hydraulic conductivity tensor.
See also section 12.14 “ANI Horizontal anisotropy module” in the iMOD v5.2 manual for further explanation.
- Parameters
anisotropy_factor (xr.DataArray) – The anisotropy factor is defined perpendicular to the main principal axis. The factor is between 0.0 (full anisotropic) and 1.0 (full isotropic)
anisotropy_angle (xr.DataArray) – The anisotropy angle denotes the angle along the main principal axis (highest permeability k) measured in degrees from north (0°), east (90°), south (180°) and west (270°).
- dataset¶
- class imod.flow.HorizontalFlowBarrier(id_name, geometry, resistance, layer=None)[source]¶
Bases:
imod.flow.pkgbase.Package
Horizontal barriers obstructing flow such as semi- or impermeable fault zone or a sheet pile wall are defined for each model layer by a *.GEN line file.
- Parameters
id_name (str or list of str) – name of the barrier
geometry (object array of shapely LineStrings) – geometry of barriers, should be lines
resistance (float or list of floats) – resistance of the barrier (d).
layer (Optional, int) – layer where barrier is located. Defaults to None.
- dataset¶
- class imod.flow.HorizontalHydraulicConductivity(k_horizontal=None)[source]¶
Bases:
imod.flow.pkgbase.Package
Horizontal hydraulic conductivity [L/T] of the aquifers, between TOP and BOT.
This variable behaves somewhat similar to the horizontal hydraulic conductivity in MODFLOW 2005’s “Layer Property Flow” schematization.
Note however that this does not hold for the vertical hydraulic conductivity: iMODFLOW uses the vertical hydraulic conductivity to specify the hydraulic conductivity of aquitards (between BOT and TOP)
- Parameters
k_horizontal (xr.DataArray) – Horizontal hydraulic conductivity, dims
("layer", "y", "x")
.
- dataset¶
- class imod.flow.ImodflowModel(modelname, check='defer')[source]¶
Bases:
imod.flow.model.Model
- modelname¶
- Type
str
- check¶
When to perform model checks {None, “defer”, “eager”}. Defaults to “defer”.
- Type
str, optional
Examples
>>> m = Imodflow("example") >>> m["riv"] = River(...) >>> # ...etc. >>> m.time_discretization(endtime) >>> m.write()
- render(directory, render_projectfile=True)[source]¶
Render the runfile as a string, package by package.
- time_discretization(times)[source]¶
Collect all unique times from model packages and additional given times. These unique times are used as stress periods in the model. All stress packages must have the same starting time.
The time discretization in imod-python works as follows:
The datetimes of all packages you send in are always respected
Subsequently, the input data you use is always included fully as well
All times are treated as starting times for the stress: a stress is always applied until the next specified date
For this reason, a final time is required to determine the length of the last stress period
Additional times can be provided to force shorter stress periods & more detailed output
Every stress has to be defined on the first stress period (this is a modflow requirement)
Or visually (every letter a date in the time axes):
recharge a - b - c - d - e - f river g - - - - h - - - - j times - - - - - - - - - - - i model a - b - c h d - e - f i
with the stress periods defined between these dates. I.e. the model times are the set of all times you include in the model.
- Parameters
times (str, datetime; or iterable of str, datetimes.) – Times to add to the time discretization. At least one single time should be given, which will be used as the ending time of the simulation.
Examples
Add a single time:
>>> m.time_discretization("2001-01-01")
Add a daterange:
>>> m.time_discretization(pd.daterange("2000-01-01", "2001-01-01"))
Add a list of times:
>>> m.time_discretization(["2000-01-01", "2001-01-01"])
- write(directory=PosixPath('.'), result_dir=None, resultdir_is_workdir=False, render_projectfile=True)[source]¶
Writes model input files.
- Parameters
directory (str, pathlib.Path) – Directory into which the model input will be written. The model input will be written into a directory called modelname.
result_dir (str, pathlib.Path) –
Path to directory in which output will be written when running the model. Is written as the value of the
result_dir
key in the runfile.See the examples.
resultdir_is_workdir (boolean, optional) – Wether the set all input paths in the runfile relative to the output directory. Because iMOD-wq generates a number of files in its working directory, it may be advantageous to set the working directory to a different path than the runfile location.
- Returns
- Return type
None
Examples
Say we wish to write the model input to a file called input, and we desire that when running the model, the results end up in a directory called output. We may run:
>>> model.write(directory="input", result_dir="output")
And in the runfile, a value of
../../output
will be written for result_dir.
- class imod.flow.MetaSwap(boundary, landuse, rootzone_thickness, soil_physical_unit, meteo_station_number, surface_elevation, sprinkling_type, sprinkling_layer, sprinkling_capacity, wetted_area, urban_area, ponding_depth_urban, ponding_depth_rural, runoff_resistance_urban, runoff_resistance_rural, runon_resistance_urban, runon_resistance_rural, infiltration_capacity_urban, infiltration_capacity_rural, perched_water_table, lookup_and_forcing_files, soil_moisture_factor=1.0, conductivity_factor=1.0)[source]¶
Bases:
imod.flow.pkgbase.Package
The MetaSWAP package (CAP), provides the input to be converted to a MetaSWAP model, which is an external model code used to simulate the unsaturated zone.
Note that only two-dimensional DataArrays with dimensions
("y", "x")
should be supplied to this package. In the current implementation time-related files are provided as external files (“extra files”). Similar to the iMODFLOW implementation of the projectfile. For now these need to be provided as a path.MetaSWAP is developed by Alterra, Wageningen as part of the SIMGRO model code. The SIMGRO framework is intended for regions with an undulating topography and unconsolidated sediments in the (shallow) subsoil. Both shallow and deep groundwater levels can be modelled by MetaSWAP. This model is based on a simplification of ‘straight Richards’, meaning that no special processes like hysteresis, preferential flow and bypass flow are modelled. Snow is not modelled, and neither the influence of frost on the soil water conductivity. A perched watertable can be present in the SVAT column model, but interflow is not modelled. Processes that are typical for steep slopes are not included. The code contains several parameterized water management schemes, including irrigation and water level management.
References:
Van Walsum, P. E. V., 2017a. SIMGRO V7.3.3.2, Input and Output reference manual. Tech. Rep. Alterra-Report 913.3, Alterra, Wageningen. 98 pp.
Van Walsum, P. E. V., 2017b. SIMGRO V7.3.3.2, Users Guide. Tech. Rep. Alterra-Report 913.2, Alterra, Wageningen. 111 pp.
Van Walsum, P. E. V. and P. Groenendijk, 2008. “Quasi Steady-State Simulation on of the Unsaturated Zone in Groundwater Modeling of Lowland Regions.” Vadose Zone Journal 7: 769-778.
Van Walsum, P. E. V., A. A. Veldhuizen and P. Groenendijk, 2016. SIMGRO V7.2.27, Theory and model implementation. Tech. Rep. Alterra-Report 913.1, Alterra, Wageningen. 93 pp 491.
- Parameters
boundary (int or xr.DataArray of ints) – 2D boundary, used to specify active MetaSWAP elements, similar to ibound in the Boundary package
landuse (int or xr.DataArray of ints) – Landuse codes, referred to in the lookup table file luse_mswp.inp
rootzone_thickness (float or xr.DataArray of floats) – Rootzone thickness in cm (min. value is 10 centimeter).
soil_physical_unit (int or xr.DataArray of ints) – Soil Physical Unit, referred to in the lookup table file fact_mswp.inp.
meteo_station_number (float or xr.DataArray of ints) – Meteo station number, referred to by mete_mswp.inp.
surface_elevation (float or xr.DataArray of floats) – Surface Elevation (m+MSL)
sprinkling_type (int or xr.DataArray of ints) –
Sprinkling type (“Artificial Recharge Type” in iMOD manual):
0 = no occurrence
1 = from groundwater
2 = from surface water
sprinkling_layer (int or xr.DataArray of ints) –
- Number of modellayer from which water is extracted (“Artificial
Recharge Location” in iMOD manual)
sprinkling_capacity (float or xr.DataArray of floats) – Sprinkling capacity (mm/d) sets the maximum amount extracted for sprinkling (“Artificial Recharge Capacity” in iMOD manual)
wetted_area (float or xr.DataArray of floats) – Total area (m2) occupied by surface water elements. Values will be truncated by maximum cellsize.
urban_area (float or xr.DataArray of floats) – Total area (m2) occupied by urban area. Values will be truncated by maximum cellsize.
ponding_depth_urban (float or xr.DataArray of floats) – Ponding Depth Urban Area (m), specifying the acceptable depth of the ponding of water on the surface in the urban area before surface runoff occurs.
ponding_depth_rural (float or xr.DataArray of floats) – Ponding Depth Rural Area (m), specifying the acceptable depth of the ponding of water on the surface in the rural area before surface runoff occurs.
runoff_resistance_urban (float or xr.DataArray of floats) – Runoff Resistance Urban Area (day), specifying the resistance surface flow encounters in the urban area. The minimum value is equal to the model time period.
runoff_resistance_rural (float or xr.DataArray of floats) – Runoff Resistance Rural Area (day), specifying the resistance surface flow encounters in the rural area. The minimum value is equal to the model time period.
runon_resistance_urban (float or xr.DataArray of floats) – Runon Resistance Urban Area (day), specifying the resistance surface flow encounters to a model cell from an adjacent cell in the urban area. The minimum value is equal to the model time period.
runon_resistance_rural (float or xr.DataArray of floats) – Runon Resistance Rural Area (day), specifying the resistance surface flow encounters to a model cell from an adjacent cell in the rural area. The minimum value is equal to the model time period.
infiltration_capacity_urban (float or xr.DataArray of floats) – the infiltration capacity (m/d) of the soil surface in the urban area. The range is 0-1000 m/d. The NoDataValue -9999 indicates unlimited infiltration is possible.
infiltration_capacity_rural (float or xr.DataArray of floats) – the infiltration capacity (m/d) of the soil surface in the urban area. The range is 0-1000 m/d. The NoDataValue -9999 indicates unlimited infiltration is possible.
perched_water_table (float or xr.DataArray of floats) – Depth of the perched water table level (m)
soil_moisture_factor (float) – Soil Moisture Factor to adjust the soil moisture coefficient. This factor may be used during calibration. Default value is 1.0.
conductivity_factor (float) – Conductivity Factor to adjust the vertical conductivity. This factor may be used during calibration. Default value is 1.0.
lookup_and_forcing_files (list of pathlib.Path or str) – List of paths to files required by MetaSWAP. This a list of lookup tables and meteorological information that is required by MetaSwap. Note that MetaSwap looks for files with a specific name, so calling “luse_svat.inp” something else will result in errors. To view the files required, you can call:
print(MetaSwap()._required_extra)
- dataset¶
- class imod.flow.PreconditionedConjugateGradientSolver(max_iter=150, inner_iter=100, rclose=1.0, hclose=0.0001, relax=0.98, damp=1.0, damp_transient=1.0, matrix_conditioning_method=1, printout_interval=0, print_convergence_info=1)[source]¶
Bases:
imod.flow.pkgbase.Package
The Preconditioned Conjugate Gradient Solver is used to solve the finite difference equations in each step of a MODFLOW stress period.
- Parameters
max_iter (int) – is the maximum number of outer iterations - that is, calss to the solutions routine (MXITER). For a linear problem max_iter should be 1, unless more than 50 inner iterations are required, when max_iter could be as large as 10. A larger number (generally less than 100) is required for a nonlinear problem.
inner_iter (int) – is the number of inner iterations (iter1). For nonlinear problems, inner_iter usually ranges from 10 to 30; a value of 30 will be sufficient for most linear problems.
rclose (float) –
is the residual criterion for convergence, in units of cubic length per time. The units for length and time are the same as established for all model data. When the maximum absolute value of the residual at all nodes during an iteration is less than or equal to RCLOSE, and the criterion for HCLOSE is also satisfied (see below), iteration stops.
Default value: 100.0. Nota bene: this is aimed at regional modelling. For detailed studies (e.g. lab experiments) much smaller values can be required. Very general rule of thumb: should be less than 10% of smallest cell volume.
hclose (float) – is the head change criterion for convergence, in units of length. When the maximum absolute value of head change from all nodes during an iteration is less than or equal to HCLOSE, and the criterion for RCLOSE is also satisfied (see above), iteration stops. Default value: 1.0e-4. Nota bene: This is aimed at regional modelling, ` for detailed studies (e.g. lab experiments) much smaller values can be required.
relax (float, optional) – is the relaxation parameter used. Usually, RELAX = 1.0, but for some problems a value of 0.99, 0.98, or 0.97 will reduce the number of iterations required for convergence. Default value: 0.98.
matrix_conditioning_method (int, optional) –
- the flag used to select the matrix conditioning method
1 is for Modified Incomplete Cholesky (for use on scalar computers) 2 is for Polynomial (for use on vector computers or to conserve computer memory)
damp (float, optional) – the damping factor. It is typically set equal to one, which indicates no damping. A value less than 1 and greater than 0 causes damping. DAMP does not affect inner iterations; instead, it affects outer iterations. Default value: 1.0.
damp_transient (float, optional) – the damping factor for transient stress periods. it is read only when damp is specified as a negative value. If damp_transient is not read, then the single damping factor, damp, is used for both transient and steady-state stress periods.
printout_interval (int, optional) – is the printout interval for PCG. If equal to zero, it is changed to 999. The maximum head change (positive or negative) and residual change are printed for each iteration of a time step whenever the time step is an even multiple of printout_interval. This printout also occurs at the end of each stress period regardless of the value of printout_interval.
print_convergence_info (int, optional) –
- a flag that controls printing of convergence information from the solver:
0 is for printing tables of maximum head change and residual each iteration 1 is for printing only the total number of iterations 2 is for no printing 3 is for printing only if convergence fails
- dataset¶
- class imod.flow.Recharge(rate)[source]¶
Bases:
imod.flow.pkgbase.TopBoundaryCondition
Recharge provides a fixed flux boundary condition to the top layer of the groundwater system. Note that unlike in iMOD-WQ, there is only the option in iMODFLOW to apply the recharge package to the top layer.
- Parameters
rate (float or xr.DataArray of floats) – recharge rate in mm/day (NOTA BENE!), dims
("time", "y", "x")
.
- dataset¶
- class imod.flow.River(conductance, stage, bottom_elevation, infiltration_factor)[source]¶
Bases:
imod.flow.pkgbase.BoundaryCondition
The River package is used to simulate head-dependent flux boundaries. In the River package if the head in the cell falls below a certain threshold, the flux from the river to the model cell is set to a specified lower bound.
- Parameters
stage (float or xr.DataArray of floats) – is the head in the river (STAGE), dims
("layer", "y", "x")
.bottom_elevation (float or xr.DataArray of floats) – is the bottom of the riverbed (RBOT), dims
("layer", "y", "x")
.conductance (float or xr.DataArray of floats) – is the conductance of the river, dims
("layer", "y", "x")
.infiltration_factor (float or xr.DataArray of floats) –
is the infiltration factor, dims
("layer", "y", "x")
. This factor defines the reduces the conductance for infiltrating water compared to exfiltrating water:cond = A / (c * iff)
where
A
[L2] is the area where surface water and groundwater interact,c
[L] is the resistance, andiff
is the infiltration factor.The infiltration factor is thus equal or larger than 1.
- dataset¶
- class imod.flow.SpecificStorage(specific_storage)[source]¶
Bases:
imod.flow.pkgbase.Package
Specific storage [L-1]. Be careful, this is not the same as the storage coefficient.
From wikipedia (https://en.wikipedia.org/wiki/Specific_storage):
The specific storage is the amount of water that a portion of an aquifer releases from storage, per unit mass or volume of aquifer, per unit change in hydraulic head, while remaining fully saturated.
- Parameters
specific_storage (float or xr.DataArray) – Specific storage, dims
("layer", "y", "x")
.
- dataset¶
- class imod.flow.StartingHead(starting_head)[source]¶
Bases:
imod.flow.pkgbase.Package
The initial head in all cells
- Parameters
starting_head (float or xr.DataArray of floats) – is initial (starting) head—that is, head at the beginning of the simulation (SHD). starting_head must be specified for all simulations, including steady-state simulations. One value is read for every model cell. Usually, these values are read a layer at a time.
- dataset¶
- class imod.flow.StorageCoefficient(storage_coefficient)[source]¶
Bases:
imod.flow.pkgbase.Package
Storage coefficient [-]. Be careful, this is not the same as the specific storage.
From wikipedia (https://en.wikipedia.org/wiki/Specific_storage):
Storativity or the storage coefficient is the volume of water released from storage per unit decline in hydraulic head in the aquifer, per unit area of the aquifer. Storativity is a dimensionless quantity, and is always greater than 0.
Under confined conditions:
S = Ss * b, where S is the storage coefficient, Ss the specific storage, and b the aquifer thickness.
Under unconfined conditions:
S = Sy, where Sy is the specific yield
- Parameters
storage_coefficient (float or xr.DataArray) – Storage coefficient, dims = (“layer”, “y”, “x”).
- dataset¶
- class imod.flow.TimeDiscretization(timestep_duration, endtime, n_timesteps=1, transient=True, timestep_multiplier=1.0)[source]¶
Bases:
imod.flow.pkgbase.Package
Time discretisation package class.
- Parameters
timestep_duration (xr.DataArray) – is the length of the current stress period (PERLEN). If the flow solution is transient, timestep_duration specified here must be equal to that specified for the flow model. If the flow solution is steady-state, timestep_duration can be set to any desired length.
n_timesteps (int, optional) – is the number of time steps for the transient flow solution in the current stress period (NSTP). If the flow solution is steady-state, n_timestep=1. Default value is 1.
transient (bool, optional) – Flag indicating wether the flow simulation is transient (True) or False (Steady State). Default is True.
timestep_multiplier (float, optional) – is the multiplier for the length of successive time steps used in the transient flow solution (TSMULT); it is used only if n_timesteps>1. timestep_multiplier>0, the length of each flow time step within the current stress period is calculated using the geometric progression as in MODFLOW. Note that both n_timesteps and timestep_multiplier specified here must be identical to those specified in the flow model if the flow model is transient. If timestep_multiplier ≤ 0, the length of each flow time step within the current stress period is read from the record TSLNGH. This option is needed in case the length of time steps for the flow solution is not based on a geometric progression in a flow model, unlike MODFLOW. Default is 1.0.
- dataset¶
- class imod.flow.Top(top)[source]¶
Bases:
imod.flow.pkgbase.Package
The top of the aquifers
- Parameters
top (xr.DataArray of floats) – is the top elevation with dimensions
("layer", "y", "x")
. For the common situation in which the top layer represents a water-table aquifer, it may be reasonable to set`top` equal to land-surface elevation. The DataArray should at least include the layer dimension.
- dataset¶
- class imod.flow.VerticalAnisotropy(vertical_anisotropy)[source]¶
Bases:
imod.flow.pkgbase.Package
Vertical anisotropy for aquifers [-], defined as the horizontal hydraulic conductivity over the vertical hydraulic conductivity.
Use this package in combination with HorizontalHydraulicConductivity to specify the vertical hydraulic conductivity.
- Parameters
vertical_anisotropy (xr.DataArray) – Vertical anisotropy factor (Kv/Kh), dims
("layer", "y", "x")
.
- dataset¶
- class imod.flow.VerticalHydraulicConductivity(k_vertical=None)[source]¶
Bases:
imod.flow.pkgbase.Package
Vertical hydraulic conductivity [L/T] for aquitards (between BOT and TOP).
Note that this is different from MODFLOW 2005’s “Layer Property Flow” schematization. To specify the vertical hydraulic conductivity for aquifers, use VerticalAnisotropy in combination with HorizontalHydraulicConductivity.
- Parameters
k_vertical (xr.DataArray) – Vertical hydraulic conductivity, dims
("layer", "y", "x")
.
- dataset¶
- class imod.flow.Well(id_name, x, y, rate, layer, time=None)[source]¶
Bases:
imod.flow.pkgbase.BoundaryCondition
The Well package is used to simulate a specified flux to individual cells and specified in units of length3/time.
- Parameters
id_name (str or list of str) – name of the well(s).
x (float or list of floats) – x coordinate of the well(s).
y (float or list of floats) – y coordinate of the well(s).
rate (float or list of floats.) – pumping rate in the well(s).
layer ("None" or int, optional) – layer from which the pumping takes place.
time ("None" or listlike of np.datetime64, datetime.datetime, pd.Timestamp,) – cftime.datetime time during which the pumping takes place. Only need to specify if model is transient.
- dataset¶