Coverage for fixtures\package_instance_creation.py: 100%

70 statements  

« prev     ^ index     » next       coverage.py v7.5.1, created at 2024-05-08 14:15 +0200

1import geopandas as gpd 

2import numpy as np 

3import pandas as pd 

4import shapely 

5import xarray as xr 

6import xugrid as xu 

7 

8import imod 

9import imod.tests.fixtures.mf6_lake_package_fixture as mf_lake 

10from imod.mf6 import GWFGWF 

11 

12""" 

13This file is used to create instances of imod packages for testing purposes. 

14The main usage is importing ALL_PACKAGE_INSTANCES into a test- this list contains an instance of 

15each packages and boundary condition in mf6. 

16""" 

17 

18 

19def get_structured_grid_da(dtype, value=1): 

20 """ 

21 This function creates a dataarray with scalar values for a grid of 3 layers and 9 rows and columns. 

22 """ 

23 shape = nlay, nrow, ncol = 3, 9, 9 

24 dims = ("layer", "y", "x") 

25 

26 dx = 10.0 

27 dy = -10.0 

28 xmin = 0.0 

29 xmax = dx * ncol 

30 ymin = 0.0 

31 ymax = abs(dy) * nrow 

32 

33 layer = np.arange(1, nlay + 1) 

34 y = np.arange(ymax, ymin, dy) + 0.5 * dy 

35 x = np.arange(xmin, xmax, dx) + 0.5 * dx 

36 coords = {"layer": layer, "y": y, "x": x} 

37 

38 da = xr.DataArray(np.ones(shape, dtype=dtype) * value, coords=coords, dims=dims) 

39 return da 

40 

41 

42def get_unstructured_grid_da(dtype, value=1): 

43 """ 

44 This function creates an xugrid dataarray with scalar values for an unstructured grid 

45 """ 

46 grid = imod.data.circle() 

47 nface = grid.n_face 

48 nlayer = 2 

49 

50 dims = ("layer", grid.face_dimension) 

51 shape = (nlayer, nface) 

52 

53 uda = xu.UgridDataArray( 

54 xr.DataArray( 

55 np.ones(shape, dtype=dtype) * value, 

56 coords={"layer": [1, 2]}, 

57 dims=dims, 

58 ), 

59 grid=grid, 

60 ) 

61 return uda 

62 

63 

64def get_grid_da(is_unstructured, dtype, value=1): 

65 """ 

66 helper function for creating an xarray dataset of a given type 

67 Depending on the is_unstructured input parameter, it will create an array for a 

68 structured grid or for an unstructured grid. 

69 """ 

70 

71 if is_unstructured: 

72 return get_unstructured_grid_da(dtype, value) 

73 else: 

74 return get_structured_grid_da(dtype, value) 

75 

76 

77def create_lake_package(is_unstructured): 

78 is_lake1 = get_grid_da(is_unstructured, bool, False) 

79 times_of_numeric_timeseries = [ 

80 np.datetime64("2000-01-01"), 

81 np.datetime64("2000-03-01"), 

82 np.datetime64("2000-05-01"), 

83 ] 

84 numeric = xr.DataArray( 

85 np.full((len(times_of_numeric_timeseries)), 5.0), 

86 coords={"time": times_of_numeric_timeseries}, 

87 dims=["time"], 

88 ) 

89 if is_unstructured: 

90 is_lake1[1, 0] = True 

91 is_lake1[1, 1] = True 

92 lake1 = mf_lake.create_lake_data_unstructured( 

93 is_lake1, 11.0, "Naardermeer", stage=numeric 

94 ) 

95 else: 

96 is_lake1[1, 1, 0] = True 

97 is_lake1[1, 1, 1] = True 

98 lake1 = mf_lake.create_lake_data_structured( 

99 is_lake1, 11.0, "Naardermeer", stage=numeric 

100 ) 

101 outlet1 = imod.mf6.OutletManning("Naardermeer", "", 3.0, 2.0, 3.0, 4.0) 

102 return imod.mf6.Lake.from_lakes_and_outlets([lake1], [outlet1]) 

103 

104 

105def create_vertices_discretization(): 

106 """ 

107 return imod.mf6.VerticesDiscretization object 

108 """ 

109 idomain = get_grid_da(True, int, value=1) 

110 bottom = idomain * xr.DataArray([5.0, 0.0], dims=["layer"]) 

111 return imod.mf6.VerticesDiscretization(top=10.0, bottom=bottom, idomain=idomain) 

112 

113 

114def create_instance_packages(is_unstructured): 

115 """ 

116 creates instances of those modflow packages that are not boundary conditions. 

117 """ 

118 return [ 

119 imod.mf6.Dispersion( 

120 diffusion_coefficient=get_grid_da(is_unstructured, np.float32, 1e-4), 

121 longitudinal_horizontal=get_grid_da(is_unstructured, np.float32, 10), 

122 transversal_horizontal1=get_grid_da(is_unstructured, np.float32, 10), 

123 longitudinal_vertical=get_grid_da(is_unstructured, np.float32, 5), 

124 transversal_horizontal2=get_grid_da(is_unstructured, np.float32, 2), 

125 transversal_vertical=get_grid_da(is_unstructured, np.float32, 4), 

126 ), 

127 imod.mf6.InitialConditions(start=get_grid_da(is_unstructured, np.float32)), 

128 imod.mf6.MobileStorageTransfer( 

129 porosity=get_grid_da(is_unstructured, np.float32, 0.35), 

130 decay=get_grid_da(is_unstructured, np.float32, 0.01), 

131 decay_sorbed=get_grid_da(is_unstructured, np.float32, 0.02), 

132 bulk_density=get_grid_da(is_unstructured, np.float32, 1300), 

133 distcoef=get_grid_da(is_unstructured, np.float32, 0.1), 

134 ), 

135 imod.mf6.NodePropertyFlow( 

136 get_grid_da(is_unstructured, np.int32), 3.0, True, 32.0, 34.0, 7 

137 ), 

138 imod.mf6.SpecificStorage( 

139 0.001, 0.1, True, get_grid_da(is_unstructured, np.int32) 

140 ), 

141 imod.mf6.StorageCoefficient( 

142 0.001, 0.1, True, get_grid_da(is_unstructured, np.int32) 

143 ), 

144 imod.mf6.ApiPackage( 

145 maxbound=33, print_input=True, print_flows=True, save_flows=True 

146 ), 

147 ] 

148 

149 

150def create_instance_boundary_condition_packages(is_unstructured): 

151 """ 

152 creates instances of those modflow packages that are boundary conditions. 

153 """ 

154 return [ 

155 imod.mf6.ConstantConcentration( 

156 get_grid_da(is_unstructured, np.float32, 2), 

157 print_input=True, 

158 print_flows=True, 

159 save_flows=True, 

160 ), 

161 imod.mf6.ConstantHead( 

162 get_grid_da(is_unstructured, np.float32, 2), 

163 print_input=True, 

164 print_flows=True, 

165 save_flows=True, 

166 ), 

167 imod.mf6.Drainage( 

168 elevation=get_grid_da(is_unstructured, np.float64, 4), 

169 conductance=get_grid_da(is_unstructured, np.float64, 1e-3), 

170 ), 

171 imod.mf6.Evapotranspiration( 

172 surface=get_grid_da(is_unstructured, np.float64, 3), 

173 rate=get_grid_da(is_unstructured, np.float64, 2), 

174 depth=get_grid_da(is_unstructured, np.float64, 1), 

175 proportion_rate=get_grid_da(is_unstructured, np.float64, 0.2), 

176 proportion_depth=get_grid_da(is_unstructured, np.float64, 0.2), 

177 fixed_cell=True, 

178 ), 

179 imod.mf6.GeneralHeadBoundary( 

180 head=get_grid_da(is_unstructured, np.float64, 3), 

181 conductance=get_grid_da(is_unstructured, np.float64, 0.33), 

182 ), 

183 imod.mf6.Recharge( 

184 rate=get_grid_da(is_unstructured, np.float64, 0.33), 

185 ), 

186 imod.mf6.River( 

187 stage=get_grid_da(is_unstructured, np.float64, 0.33), 

188 conductance=get_grid_da(is_unstructured, np.float64, 0.33), 

189 bottom_elevation=get_grid_da(is_unstructured, np.float64, 0.33), 

190 ), 

191 imod.mf6.SourceSinkMixing( 

192 package_names=["a", "b"], 

193 concentration_boundary_type=["a", "b"], 

194 auxiliary_variable_name=["a", "b"], 

195 ), 

196 create_lake_package(is_unstructured), 

197 ] 

198 

199 

200STRUCTURED_GRID_PACKAGES = [ 

201 imod.mf6.StructuredDiscretization( 

202 2.0, get_grid_da(False, np.float32), get_grid_da(False, np.int32) 

203 ), 

204 imod.mf6.WellDisStructured( 

205 layer=[3, 2, 1], 

206 row=[1, 2, 3], 

207 column=[2, 2, 2], 

208 rate=[-5.0] * 3, 

209 ), 

210 imod.mf6.MassSourceLoading( 

211 rate=get_grid_da(False, np.float64, 0.33), 

212 print_input=True, 

213 print_flows=False, 

214 save_flows=False, 

215 ), 

216] + [ 

217 *create_instance_packages(is_unstructured=False), 

218 *create_instance_boundary_condition_packages(False), 

219] 

220 

221 

222UNSTRUCTURED_GRID_PACKAGES = ( 

223 [ 

224 imod.mf6.WellDisVertices( 

225 layer=[1, 2, 1], 

226 cell2d=[3, 12, 23], 

227 rate=[-0.1, 0.2, 0.3], 

228 print_input=False, 

229 print_flows=False, 

230 save_flows=False, 

231 ) 

232 ] 

233 + [create_vertices_discretization()] 

234 + [ 

235 *create_instance_packages(is_unstructured=True), 

236 *create_instance_boundary_condition_packages(True), 

237 ] 

238) 

239 

240GRIDLESS_PACKAGES = [ 

241 imod.mf6.adv.Advection("upstream"), 

242 imod.mf6.Buoyancy( 

243 reference_density=1000.0, 

244 reference_concentration=[4.0, 25.0], 

245 density_concentration_slope=[0.7, -0.375], 

246 modelname=["gwt-1", "gwt-2"], 

247 species=["salinity", "temperature"], 

248 ), 

249 imod.mf6.OutputControl(), 

250 imod.mf6.SolutionPresetSimple(modelnames=["gwf-1"]), 

251 imod.mf6.TimeDiscretization( 

252 xr.DataArray( 

253 data=[0.001, 7.0, 365.0], 

254 coords={"time": pd.date_range("2000-01-01", "2000-01-03")}, 

255 dims=["time"], 

256 ), 

257 23, 

258 1.02, 

259 ), 

260 imod.mf6.Well( 

261 x=[1.0, 6002.0], 

262 y=[3.0, 5004.0], 

263 screen_top=[0.0, 0.0], 

264 screen_bottom=[-10.0, -10.0], 

265 rate=[1.0, 3.0], 

266 print_flows=True, 

267 validate=True, 

268 ), 

269 imod.mf6.HorizontalFlowBarrierHydraulicCharacteristic( 

270 geometry=gpd.GeoDataFrame( 

271 geometry=[ 

272 shapely.linestrings([-1000.0, 1000.0], [500.0, 500.0]), 

273 ], 

274 data={ 

275 "hydraulic_characteristic": [1e-3], 

276 "ztop": [10.0], 

277 "zbottom": [0.0], 

278 }, 

279 ), 

280 ), 

281 imod.mf6.HorizontalFlowBarrierMultiplier( 

282 geometry=gpd.GeoDataFrame( 

283 geometry=[ 

284 shapely.linestrings([-1000.0, 1000.0], [500.0, 500.0]), 

285 ], 

286 data={ 

287 "multiplier": [1.5], 

288 "ztop": [10.0], 

289 "zbottom": [0.0], 

290 }, 

291 ), 

292 ), 

293 imod.mf6.HorizontalFlowBarrierResistance( 

294 geometry=gpd.GeoDataFrame( 

295 geometry=[ 

296 shapely.linestrings([-1000.0, 1000.0], [500.0, 500.0]), 

297 ], 

298 data={ 

299 "resistance": [1e3], 

300 "ztop": [10.0], 

301 "zbottom": [0.0], 

302 }, 

303 ), 

304 ), 

305 imod.mf6.GWFGWT("model_name1", "model_name2"), 

306] 

307 

308 

309def create_exchange_package() -> list[GWFGWF]: 

310 cell_id1 = xr.DataArray([[1, 1], [2, 1], [3, 1]]) 

311 cell_id2 = xr.DataArray([[1, 2], [2, 2], [3, 2]]) 

312 layer = np.array([12, 13, 14]) 

313 cl1 = np.ones(len(cell_id1)) 

314 cl2 = np.ones(len(cell_id1)) 

315 hwva = cl1 + cl2 

316 

317 return [ 

318 imod.mf6.GWFGWF( 

319 "submodel_1", 

320 "submodel_2", 

321 cell_id1=cell_id1, 

322 cell_id2=cell_id2, 

323 layer=layer, 

324 cl1=cl1, 

325 cl2=cl2, 

326 hwva=hwva, 

327 ) 

328 ] 

329 

330 

331ALL_PACKAGE_INSTANCES = ( 

332 GRIDLESS_PACKAGES 

333 + STRUCTURED_GRID_PACKAGES 

334 + UNSTRUCTURED_GRID_PACKAGES 

335 + create_exchange_package() 

336)