Coverage for fixtures\flow_basic_fixture.py: 100%
139 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
1from typing import Callable, NamedTuple
3import numpy as np
4import pandas as pd
5import pytest
6import xarray as xr
9@pytest.fixture(scope="module")
10def basic_dis():
11 """Basic model discretization"""
13 shape = nlay, nrow, ncol = 3, 9, 9
15 dx = 10.0
16 dy = -10.0
17 dz = np.array([5, 30, 100])
18 xmin = 0.0
19 xmax = dx * ncol
20 ymin = 0.0
21 ymax = abs(dy) * nrow
22 dims = ("layer", "y", "x")
24 layer = np.arange(1, nlay + 1)
25 y = np.arange(ymax, ymin, dy) + 0.5 * dy
26 x = np.arange(xmin, xmax, dx) + 0.5 * dx
27 coords = {"layer": layer, "y": y, "x": x}
29 ibound = xr.DataArray(np.ones(shape, dtype=np.int32), coords=coords, dims=dims)
31 surface = 0.0
32 interfaces = np.insert((surface - np.cumsum(dz)), 0, surface)
34 bottom = xr.DataArray(interfaces[1:], coords={"layer": layer}, dims="layer")
35 top = xr.DataArray(interfaces[:-1], coords={"layer": layer}, dims="layer")
37 return ibound, top, bottom
40class BasicDisSettings(NamedTuple):
41 nlay: int = 1
42 nrow: int = 1
43 ncol: int = 1
44 zstart: float = 0.0
45 zstop: float = -1.0
46 ystart: float = 0.0
47 ystop: float = 1.0
48 xstart: float = 0.0
49 xstop: float = 1.0
50 dx: float = 1.0
51 dy: float = 1.0
52 dz: float = 1.0
53 space_generator: Callable = np.linspace
56@pytest.fixture
57def parameterizable_basic_dis(request):
58 settings = request.param
59 shape = (settings.nlay, settings.nrow, settings.ncol)
61 x = (
62 settings.space_generator(
63 settings.xstart + 1, settings.xstop + 1, settings.ncol + 1, endpoint=True
64 )
65 - 1
66 )
67 y = (
68 settings.space_generator(
69 settings.ystop + 1, settings.ystart + 1, settings.nrow + 1, endpoint=True
70 )
71 - 1
72 )
73 z = (
74 settings.space_generator(
75 settings.zstart - 1, settings.zstop - 1, settings.nlay + 1, endpoint=True
76 )
77 + 1
78 )
80 xc = (x[:-1] + x[1:]) / 2
81 yc = (y[:-1] + y[1:]) / 2
83 dx = x[1:] - x[:-1]
84 dy = y[1:] - y[:-1]
86 layers = np.arange(settings.nlay) + 1
88 idomain = xr.DataArray(
89 np.ones(shape, dtype=np.int32),
90 coords={"layer": layers, "y": yc, "x": xc},
91 name="idomain",
92 )
94 # Assign dx and dy coordinates. They are needed for certain methods like 'coord_reference'
95 if np.all(np.isclose(dx, dx[0])):
96 idomain = idomain.assign_coords({"dx": dx[0]})
97 else:
98 idomain = idomain.assign_coords({"dx": ("x", dx)})
99 if np.all(np.isclose(dy, dy[0])):
100 idomain = idomain.assign_coords({"dy": dy[0]})
101 else:
102 idomain = idomain.assign_coords({"dy": ("y", dy)})
104 top = xr.DataArray(z[:-1], coords={"layer": layers})
105 bottom = xr.DataArray(z[1:], coords={"layer": layers})
107 return idomain, top, bottom
110@pytest.fixture(scope="module")
111def three_days():
112 """Simple time discretization of three days"""
113 return pd.date_range(start="2018-01-01", end="2018-01-03", freq="D")
116@pytest.fixture(scope="module")
117def two_days():
118 """Simple time discretization of two days"""
119 return pd.date_range(start="2018-01-01", end="2018-01-02", freq="D")
122@pytest.fixture(scope="module")
123def well_df(three_days):
124 nrow = 9
125 ncol = 9
126 dx = 10.0
127 dy = -10.0
128 xmin = 0.0
129 xmax = dx * ncol
130 ymin = 0.0
131 ymax = abs(dy) * nrow
133 y = np.arange(ymax, ymin, dy) + 0.5 * dy
134 x = np.arange(xmin, xmax, dx) + 0.5 * dx
136 times = three_days
137 n_repeats = int(len(x) / len(times)) + 1
139 df = pd.DataFrame()
140 df["id_name"] = np.arange(len(x)).astype(str)
141 df["x"] = x
142 df["y"] = y
143 df["rate"] = dx * dy * -1 * 0.5
144 df["time"] = np.tile(times, n_repeats)[: len(x)]
145 df["layer"] = 2
146 return df
149@pytest.fixture(scope="module")
150def get_render_dict():
151 """
152 Helper function to return dict to render.
154 Fixture returns local helper function, so that the helper function
155 is only evaluated when called.
156 See: https://stackoverflow.com/a/51389067
157 """
159 def _get_render_dict(package, directory, globaltimes, nlayer, system_index=1):
160 composition = package.compose(
161 directory,
162 globaltimes,
163 nlayer,
164 system_index=system_index,
165 )
167 return dict(
168 pkg_id=package._pkg_id,
169 name=package.__class__.__name__,
170 variable_order=package._variable_order,
171 package_data=composition[package._pkg_id],
172 )
174 return _get_render_dict
177@pytest.fixture(scope="module")
178def metaswap_dict(basic_dis):
179 ibound, _, _ = basic_dis
181 active = ibound.isel(layer=0)
183 d = {}
184 d["boundary"] = active
185 d["landuse"] = active
186 d["rootzone_thickness"] = 1.2
187 d["soil_physical_unit"] = active
188 d["meteo_station_number"] = active
189 d["surface_elevation"] = 0.0
190 d["sprinkling_type"] = active
191 d["sprinkling_layer"] = active
192 d["sprinkling_capacity"] = 1000.0
193 d["wetted_area"] = 30.0
194 d["urban_area"] = 30.0
195 d["ponding_depth_urban"] = 0.02
196 d["ponding_depth_rural"] = 0.005
197 d["runoff_resistance_urban"] = 1.5
198 d["runoff_resistance_rural"] = 1.5
199 d["runon_resistance_urban"] = 1.5
200 d["runon_resistance_rural"] = 1.5
201 d["infiltration_capacity_urban"] = 10.0
202 d["infiltration_capacity_rural"] = 2.0
203 d["perched_water_table"] = 0.5
204 d["soil_moisture_factor"] = 1.0
205 d["conductivity_factor"] = 1.0
207 d["lookup_and_forcing_files"] = [
208 "fact_svat.inp",
209 "luse_svat.inp",
210 "mete_grid.inp",
211 "para_sim.inp",
212 "tiop_sim.inp",
213 "init_svat.inp",
214 "comp_post.inp",
215 "sel_key_svat_per.inp",
216 ]
218 return d
221@pytest.fixture(scope="module")
222def horizontal_flow_barrier_gdf(basic_dis):
223 """GeoDataframe that can be used to initiate HorizontalFlowBarriers"""
224 import geopandas as gpd
225 from shapely.geometry import LineString
227 ibound, _, _ = basic_dis
229 x = ibound.x.values
230 y = ibound.y.values
232 line1 = LineString([(x[1], y[1]), (x[1], y[-2])])
233 line2 = LineString([(x[4], y[1]), (x[4], y[-2])])
235 lines = np.array([line1, line2, line1, line2], dtype="object")
236 hfb_layers = np.array([1, 1, 3, 3])
237 id_name = ["left_upper", "right_upper", "left_lower", "right_lower"]
239 hfb_gdf = gpd.GeoDataFrame(
240 geometry=lines, data=dict(layer=hfb_layers, resistance=100.0, id_name=id_name)
241 )
243 return hfb_gdf