Coverage for C:\src\imod-python\imod\mf6\model_gwf.py: 100%

47 statements  

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

1from __future__ import annotations 

2 

3from typing import Optional 

4 

5import cftime 

6import numpy as np 

7 

8from imod.logging import init_log_decorator 

9from imod.mf6 import ConstantHead 

10from imod.mf6.clipped_boundary_condition_creator import create_clipped_boundary 

11from imod.mf6.model import Modflow6Model 

12from imod.typing import GridDataArray 

13 

14 

15class GroundwaterFlowModel(Modflow6Model): 

16 """ 

17 The GroundwaterFlowModel (GWF) simulates flow of (liquid) groundwater. 

18 More information can be found here: 

19 https://water.usgs.gov/water-resources/software/MODFLOW-6/mf6io_6.4.2.pdf#page=27 

20 

21 Parameters 

22 ---------- 

23 

24 listing_file: Optional[str] = None 

25 name of the listing file to create for this GWF model. If not specified, 

26 then the name of the list file will be the basename of the GWF model 

27 name file and the 'lst' extension. 

28 print_input: bool = False 

29 keyword to indicate that the list of all model stress package 

30 information will be written to the listing file immediately after it is 

31 read. 

32 print_flows: bool = False 

33 keyword to indicate that the list of all model package flow rates will 

34 be printed to the listing file for every stress period time step in 

35 which "BUDGET PRINT" is specified in Output Control. 

36 save_flows: bool = False 

37 indicate that all model package flow terms will be written to the file 

38 specified with "BUDGET FILEOUT" in Output Control. 

39 newton: bool = False 

40 activates the Newton-Raphson formulation for groundwater flow between 

41 connected, convertible groundwater cells and stress packages that 

42 support calculation of Newton-Raphson terms for groundwater exchanges. 

43 under_relaxation: bool = False, 

44 indicates whether the groundwater head in a cell will be under-relaxed when 

45 water levels fall below the bottom of the model below any given cell. By 

46 default, Newton-Raphson UNDER_RELAXATION is not applied. 

47 """ 

48 

49 _mandatory_packages = ("npf", "ic", "oc", "sto") 

50 _model_id = "gwf6" 

51 _template = Modflow6Model._initialize_template("gwf-nam.j2") 

52 

53 @init_log_decorator() 

54 def __init__( 

55 self, 

56 listing_file: Optional[str] = None, 

57 print_input: bool = False, 

58 print_flows: bool = False, 

59 save_flows: bool = False, 

60 newton: bool = False, 

61 under_relaxation: bool = False, 

62 ): 

63 super().__init__() 

64 self._options = { 

65 "listing_file": listing_file, 

66 "print_input": print_input, 

67 "print_flows": print_flows, 

68 "save_flows": save_flows, 

69 "newton": newton, 

70 "under_relaxation": under_relaxation, 

71 } 

72 

73 def clip_box( 

74 self, 

75 time_min: Optional[cftime.datetime | np.datetime64 | str] = None, 

76 time_max: Optional[cftime.datetime | np.datetime64 | str] = None, 

77 layer_min: Optional[int] = None, 

78 layer_max: Optional[int] = None, 

79 x_min: Optional[float] = None, 

80 x_max: Optional[float] = None, 

81 y_min: Optional[float] = None, 

82 y_max: Optional[float] = None, 

83 state_for_boundary: Optional[GridDataArray] = None, 

84 ): 

85 clipped = super().clip_box( 

86 time_min, time_max, layer_min, layer_max, x_min, x_max, y_min, y_max 

87 ) 

88 

89 clipped_boundary_condition = self.__create_boundary_condition_clipped_boundary( 

90 self, clipped, state_for_boundary 

91 ) 

92 if clipped_boundary_condition is not None: 

93 clipped["chd_clipped"] = clipped_boundary_condition 

94 

95 clipped.purge_empty_packages() 

96 

97 return clipped 

98 

99 def __create_boundary_condition_clipped_boundary( 

100 self, 

101 original_model: Modflow6Model, 

102 clipped_model: Modflow6Model, 

103 state_for_boundary: Optional[GridDataArray], 

104 ): 

105 unassigned_boundary_original_domain = ( 

106 self.__create_boundary_condition_for_unassigned_boundary( 

107 original_model, state_for_boundary 

108 ) 

109 ) 

110 

111 return self.__create_boundary_condition_for_unassigned_boundary( 

112 clipped_model, state_for_boundary, [unassigned_boundary_original_domain] 

113 ) 

114 

115 @staticmethod 

116 def __create_boundary_condition_for_unassigned_boundary( 

117 model: Modflow6Model, 

118 state_for_boundary: Optional[GridDataArray], 

119 additional_boundaries: Optional[list[ConstantHead]] = None, 

120 ): 

121 if state_for_boundary is None: 

122 return None 

123 

124 constant_head_packages = [ 

125 pkg for name, pkg in model.items() if isinstance(pkg, ConstantHead) 

126 ] 

127 

128 additional_boundaries = [ 

129 item for item in additional_boundaries or [] if item is not None 

130 ] 

131 

132 constant_head_packages.extend(additional_boundaries) 

133 

134 return create_clipped_boundary( 

135 model.domain, state_for_boundary, constant_head_packages 

136 ) 

137 

138 def is_use_newton(self): 

139 return self._options["newton"] 

140 

141 def set_newton(self, is_newton: bool) -> None: 

142 self._options["newton"] = is_newton 

143 

144 def update_buoyancy_package(self, transport_models_per_flow_model) -> None: 

145 """ 

146 If the simulation is partitioned, then the buoyancy package, if present, 

147 must be updated for the renamed transport models. 

148 """ 

149 buoyancy_key = self._get_pkgkey("buy") 

150 if buoyancy_key is None: 

151 return 

152 buoyancy_package = self[buoyancy_key] 

153 transport_models_old = buoyancy_package.get_transport_model_names() 

154 if len(transport_models_old) == len(transport_models_per_flow_model): 

155 buoyancy_package.update_transport_models(transport_models_per_flow_model)