Coverage for C:\src\imod-python\imod\flow\cap.py: 100%
50 statements
« prev ^ index » next coverage.py v7.4.4, created at 2024-04-08 13:27 +0200
« prev ^ index » next coverage.py v7.4.4, created at 2024-04-08 13:27 +0200
1import os
2import pathlib
4import jinja2
6from imod.flow.pkgbase import Package
9class MetaSwap(Package):
10 """
11 The MetaSWAP package (CAP), provides the input to be converted to a
12 MetaSWAP model, which is an external model code used to simulate the
13 unsaturated zone.
15 Note that only two-dimensional DataArrays with dimensions ``("y", "x")``
16 should be supplied to this package. In the current implementation
17 time-related files are provided as external files ("extra files"). Similar
18 to the iMODFLOW implementation of the projectfile. For now these need to be
19 provided as a path.
21 MetaSWAP is developed by Alterra, Wageningen as part of the SIMGRO model
22 code. The SIMGRO framework is intended for regions with an undulating
23 topography and unconsolidated sediments in the (shallow) subsoil. Both
24 shallow and deep groundwater levels can be modelled by MetaSWAP. This
25 model is based on a simplification of ‘straight Richards’, meaning that no
26 special processes like hysteresis, preferential flow and bypass flow are
27 modelled. Snow is not modelled, and neither the influence of frost on the
28 soil water conductivity. A perched watertable can be present in the SVAT
29 column model, but interflow is not modelled. Processes that are typical for
30 steep slopes are not included. The code contains several parameterized
31 water management schemes, including irrigation and water level management.
33 References:
35 * Van Walsum, P. E. V., 2017a. SIMGRO V7.3.3.2, Input and Output reference
36 manual. Tech. Rep. Alterra-Report 913.3, Alterra, Wageningen. 98 pp.
38 * Van Walsum, P. E. V., 2017b. SIMGRO V7.3.3.2, Users Guide. Tech. Rep.
39 Alterra-Report 913.2, Alterra, Wageningen. 111 pp.
41 * Van Walsum, P. E. V. and P. Groenendijk, 2008. "Quasi Steady-State
42 Simulation on of the Unsaturated Zone in Groundwater Modeling of Lowland
43 Regions." Vadose Zone Journal 7: 769-778.
45 * Van Walsum, P. E. V., A. A. Veldhuizen and P. Groenendijk, 2016. SIMGRO
46 V7.2.27, Theory and model implementation. Tech. Rep. Alterra-Report 913.1,
47 Alterra, Wageningen. 93 pp 491.
49 Parameters
50 ----------
51 boundary : int or xr.DataArray of ints
52 2D boundary, used to specify active MetaSWAP elements, similar to
53 ibound in the Boundary package
55 landuse : int or xr.DataArray of ints
56 Landuse codes, referred to in the lookup table file luse_mswp.inp
58 rootzone_thickness : float or xr.DataArray of floats
59 Rootzone thickness in cm (min. value is 10 centimeter).
61 soil_physical_unit : int or xr.DataArray of ints
62 Soil Physical Unit, referred to in the lookup table file fact_mswp.inp.
64 meteo_station_number : float or xr.DataArray of ints
65 Meteo station number, referred to by mete_mswp.inp.
67 surface_elevation : float or xr.DataArray of floats
68 Surface Elevation (m+MSL)
70 sprinkling_type : int or xr.DataArray of ints
71 Sprinkling type ("Artificial Recharge Type" in iMOD manual):
73 * 0 = no occurrence
74 * 1 = from groundwater
75 * 2 = from surface water
77 sprinkling_layer : int or xr.DataArray of ints
78 Number of modellayer from which water is extracted ("Artificial
79 Recharge Location" in iMOD manual)
81 sprinkling_capacity : float or xr.DataArray of floats
82 Sprinkling capacity (mm/d) sets the maximum amount extracted for
83 sprinkling ("Artificial Recharge Capacity" in iMOD manual)
85 wetted_area : float or xr.DataArray of floats
86 Total area (m2) occupied by surface water elements. Values will be
87 truncated by maximum cellsize.
89 urban_area : float or xr.DataArray of floats
90 Total area (m2) occupied by urban area. Values will be truncated by
91 maximum cellsize.
93 ponding_depth_urban : float or xr.DataArray of floats
94 Ponding Depth Urban Area (m), specifying the acceptable depth of the
95 ponding of water on the surface in the urban area before surface runoff
96 occurs.
98 ponding_depth_rural : float or xr.DataArray of floats
99 Ponding Depth Rural Area (m), specifying the acceptable depth of the
100 ponding of water on the surface in the rural area before surface runoff
101 occurs.
103 runoff_resistance_urban : float or xr.DataArray of floats
104 Runoff Resistance Urban Area (day), specifying the resistance surface
105 flow encounters in the urban area. The minimum value is equal to the
106 model time period.
108 runoff_resistance_rural : float or xr.DataArray of floats
109 Runoff Resistance Rural Area (day), specifying the resistance surface
110 flow encounters in the rural area. The minimum value is equal to the
111 model time period.
113 runon_resistance_urban : float or xr.DataArray of floats
114 Runon Resistance Urban Area (day), specifying the resistance surface
115 flow encounters to a model cell from an adjacent cell in the urban
116 area. The minimum value is equal to the model time period.
118 runon_resistance_rural : float or xr.DataArray of floats
119 Runon Resistance Rural Area (day), specifying the resistance surface
120 flow encounters to a model cell from an adjacent cell in the rural
121 area. The minimum value is equal to the model time period.
123 infiltration_capacity_urban : float or xr.DataArray of floats
124 the infiltration capacity (m/d) of the soil surface in the urban area.
125 The range is 0-1000 m/d. The NoDataValue -9999 indicates unlimited
126 infiltration is possible.
128 infiltration_capacity_rural : float or xr.DataArray of floats
129 the infiltration capacity (m/d) of the soil surface in the urban area.
130 The range is 0-1000 m/d. The NoDataValue -9999 indicates unlimited
131 infiltration is possible.
133 perched_water_table : float or xr.DataArray of floats
134 Depth of the perched water table level (m)
136 soil_moisture_factor : float
137 Soil Moisture Factor to adjust the soil moisture coefficient. This
138 factor may be used during calibration. Default value is 1.0.
140 conductivity_factor : float
141 Conductivity Factor to adjust the vertical conductivity. This factor
142 may be used during calibration. Default value is 1.0.
144 lookup_and_forcing_files : list of pathlib.Path or str
145 List of paths to files required by MetaSWAP. This a list of
146 lookup tables and meteorological information that is required by
147 MetaSwap. Note that MetaSwap looks for files with a specific name, so
148 calling "luse_svat.inp" something else will result in errors. To view
149 the files required, you can call: ``print(MetaSwap()._required_extra)``
151 """
153 _pkg_id = "cap"
154 _variable_order = [
155 "boundary",
156 "landuse",
157 "rootzone_thickness",
158 "soil_physical_unit",
159 "meteo_station_number",
160 "surface_elevation",
161 "sprinkling_type",
162 "sprinkling_layer",
163 "sprinkling_capacity",
164 "wetted_area",
165 "urban_area",
166 "ponding_depth_urban",
167 "ponding_depth_rural",
168 "runoff_resistance_urban",
169 "runoff_resistance_rural",
170 "runon_resistance_urban",
171 "runon_resistance_rural",
172 "infiltration_capacity_urban",
173 "infiltration_capacity_rural",
174 "perched_water_table",
175 "soil_moisture_factor",
176 "conductivity_factor",
177 ]
179 _template_projectfile = jinja2.Template(
180 "0001, ({{pkg_id}}), 1, {{name}}, {{variable_order}}\n"
181 '{{"{:03d}".format(variable_order|length)}}, {{"{:03d}".format(n_entry)}}\n'
182 "{%- for variable in variable_order%}\n" # Preserve variable order
183 "{%- for layer, value in package_data[variable].items()%}\n"
184 "{%- if value is string %}\n"
185 # If string then assume path:
186 # 1 indicates the layer is activated
187 # 2 indicates the second element of the final two elements should be read
188 # 1.000 is the multiplication factor
189 # 0.000 is the addition factor
190 # -9999 indicates there is no data, following iMOD usual practice
191 '1, 2, {{"{:03d}".format(layer)}}, 1.000, 0.000, -9999., {{value}}\n'
192 "{%- else %}\n"
193 # Else assume a constant value is provided
194 '1, 1, {{"{:03d}".format(layer)}}, 1.000, 0.000, {{value}}, ""\n'
195 "{%- endif %}\n"
196 "{%- endfor %}\n"
197 "{%- endfor %}\n"
198 # Section for EXTRA FILES comes below
199 '{{"{:03d}".format(lookup_and_forcing_files|length)}},extra files\n'
200 "{%- for file in lookup_and_forcing_files %}\n"
201 "{{file}}\n"
202 "{%- endfor %}\n"
203 )
205 # TODO: Check which of these actually are required.
206 _required_extra = [
207 "fact_svat.inp",
208 "luse_svat.inp",
209 "mete_grid.inp",
210 "para_sim.inp",
211 "tiop_sim.inp",
212 "init_svat.inp",
213 "comp_post.inp",
214 "sel_key_svat_per.inp",
215 ]
217 def __init__(
218 self,
219 boundary,
220 landuse,
221 rootzone_thickness,
222 soil_physical_unit,
223 meteo_station_number,
224 surface_elevation,
225 sprinkling_type,
226 sprinkling_layer,
227 sprinkling_capacity,
228 wetted_area,
229 urban_area,
230 ponding_depth_urban,
231 ponding_depth_rural,
232 runoff_resistance_urban,
233 runoff_resistance_rural,
234 runon_resistance_urban,
235 runon_resistance_rural,
236 infiltration_capacity_urban,
237 infiltration_capacity_rural,
238 perched_water_table,
239 lookup_and_forcing_files,
240 soil_moisture_factor=1.0,
241 conductivity_factor=1.0,
242 ):
243 super().__init__()
244 self.dataset["boundary"] = boundary
245 self.dataset["landuse"] = landuse
246 self.dataset["rootzone_thickness"] = rootzone_thickness
247 self.dataset["soil_physical_unit"] = soil_physical_unit
248 self.dataset["meteo_station_number"] = meteo_station_number
249 self.dataset["surface_elevation"] = surface_elevation
250 self.dataset["sprinkling_type"] = sprinkling_type
251 self.dataset["sprinkling_layer"] = sprinkling_layer
252 self.dataset["sprinkling_capacity"] = sprinkling_capacity
253 self.dataset["wetted_area"] = wetted_area
254 self.dataset["urban_area"] = urban_area
255 self.dataset["ponding_depth_urban"] = ponding_depth_urban
256 self.dataset["ponding_depth_rural"] = ponding_depth_rural
257 self.dataset["runoff_resistance_urban"] = runoff_resistance_urban
258 self.dataset["runoff_resistance_rural"] = runoff_resistance_rural
259 self.dataset["runon_resistance_urban"] = runon_resistance_urban
260 self.dataset["runon_resistance_rural"] = runon_resistance_rural
261 self.dataset["infiltration_capacity_urban"] = infiltration_capacity_urban
262 self.dataset["infiltration_capacity_rural"] = infiltration_capacity_rural
263 self.dataset["perched_water_table"] = perched_water_table
264 self.dataset["soil_moisture_factor"] = soil_moisture_factor
265 self.dataset["conductivity_factor"] = conductivity_factor
266 self.lookup_and_forcing_files = lookup_and_forcing_files
268 def _force_absolute_path(self, f):
269 """Force absolute path, because projectfile cannot handle relative paths"""
270 return str(pathlib.Path(f).resolve())
272 def _render_projectfile(self, **kwargs):
273 """
274 Returns
275 -------
276 rendered : str
277 The rendered projfectfile part,
278 if part of PkgGroup: for a single boundary condition system.
279 """
280 lookup_and_forcing_files = [
281 self._force_absolute_path(file) for file in self.lookup_and_forcing_files
282 ]
283 kwargs["lookup_and_forcing_files"] = lookup_and_forcing_files
284 return self._template_projectfile.render(**kwargs)
286 def check_lookup_and_forcing_files(self):
287 """Check for presence of required MetaSWAP input files."""
288 filenames = set([os.path.basename(f) for f in self.lookup_and_forcing_files])
289 missing_files = set(self._required_extra) - filenames
290 if len(missing_files) > 0:
291 raise ValueError(f"Missing MetaSWAP input files {missing_files}")
293 def _pkgcheck(self, active_cells=None):
294 # Dataset.dims does not return a tuple, like DataArray does.
295 # http://xarray.pydata.org/en/stable/generated/xarray.Dataset.dims.html
296 # Frozen(SortedKeysDict).keys() does not preserve ordering in keys
297 # Therefore convert to set and check for symmetric difference,
298 # to see if any dimensions are missing or unexpected dimensions are included.
299 dims = set(self.dataset.dims.keys())
301 if len(dims.symmetric_difference(("y", "x"))) > 0:
302 raise ValueError(f'Dataset dims not ("y", "x"), instead got {dims}')
303 self.check_lookup_and_forcing_files()