## 2 Materials and Methods
### Description of the Numerical Model
#### 2.2.2 Characteristics of the Model  
The numerical model is a squared domain with a thickness of 180 m and a side of 2200 m (Figures 2 and 3). The simulated mine is placed just in the center of the model. Thus, the distance between the underground reservoir and the outer boundaries is of 1000 m. This distance allows to minimize the effects of the outer BCs on the results. The underground reservoir is represented by the nine underground chambers as described above (CH1 to CH9) that are hydraulically connected by galleries. The operation shaft is modeled using a rectangular prism adjacent to CH1 that connects the mine to the surface (Figures 2 and 3). The chambers and the operation shaft are modelled as independent domains that behave as linear reservoirs. Hydrogeologically, the groundwater flow behavior in the Martelange site is governed by fractures. However, it is known that a fracture-based groundwater flow is particularly difficult to implement in classical groundwater numerical models. In these types of cases, the underground medium is modeled using an equivalent porous medium (EPM) approach [42]. EPM approaches have been proven to be suitable to estimate the global groundwater behavior by numerous authors, including [48; 49], particularly in sites with the presence of a high density of fractures as in the study site. Concerning the hydraulic parameters, the \(K\) was modified in four different scenarios to observe its influence on the system behavior (note that the given EPM approach is considered, and \(K\) refers to the equivalent hydraulic conductivity). The values adopted for \(K\) in the different scenarios are specified in Section 2.2.3. The other hydraulic parameters were the same in all scenarios: the saturated and residual water contents were 0.05 and 0.01, respectively, while the specific storage coefficient was \(10^{-4}\) m\({}^{-1}\). The chosen parameters are typical of slate mines [22; 37]. The spatial and temporal discretizations were as follows: The domain was divided vertically in 29 layers and horizontally discretized using 3D prismatic elements. The maximum horizontal size of the elements was 150 m at the outer boundaries, and the mesh was refined around the mine where the horizontal size was about 5 m (Figure 3). The model was composed of 64,844 nodes and 38,680 elements. The total simulation length was one year and was divided into constant time steps of 15 minutes. Larger time steps would induce convergence problems and errors. Three different types of BCs were implemented. First, Dirichlet BCs were used to prescribe the piezometric head on the W and E outer boundaries of the model. Two different hypotheses (scenarios TOP and MIDDLE) were considered to assess the influence of the elevation of the piezometric head on the system performance. In the TOP scenarios, the piezometric head was fixed at an elevation (i.e., with respect to the bottom of the model) of 121 and 120 m on the upgradient (W) and downgradient (E) sides, respectively. Consequently, the mine was totally full of water in the natural conditions. In the MIDDLE scenarios, groundwater head was prescribed at 76 and 75 m on the upgradient (W) and downgradient (E) sides, respectively. In this case, only half of the volume of the mine was full of water under natural conditions. As a result of the prescribed heads, the hydraulic gradient was 4.6 \(\times\)\(10^{-4}\) for both hypotheses, and the groundwater flowed from W to E (Figure 2). The BCs implemented at the outer boundaries did not change through the simulations. Secondly, no-flow BCs were assigned to the top and the bottom of the modeled domain and to the N and S boundaries. Finally, pumping and discharge operations were simulated by use of a Neuman BC prescribing discharged or pumped water from the underground reservoir through the operation shaft.