## 2 Methodology
### Synthetic Waveforms Based on the Reciprocity Theorem
\[u(\mathbf{x}_{r};\mathbf{x}_{s})= G(\mathbf{x}_{r};\mathbf{x}_{s})*S(t) \tag{2}\]  
Assuming that the wavefield modeling is carried out at a preset highest frequency of 100 Hz, and the average P wave velocity is set as 4500 m/s, then the shortest wavelength of modeled wavefield will be 45 m. According to the sampling conditions of accurate wavefield simulation described in Section 2.3: for an accurate SEM numerical simulation, each wavelength should contain at least 4.5 grids [47], that is, the grid length should be less than 10 m. By referring to the mining zone of the Yongshaba mine (details in Section 4.1), the dimension is about \(3000\times 1000\times 500\) m\({}^{3}\), then the number of grids reaches \(300\times 100\times 50=\) 1,500,000. If the recording waveforms are fitted through forward modeling waveforms, the Green's function at each candidate source grid point will be modeled in the GSWI with relatively large computation. To solve this computational problem, it is necessary to store the Green's function in wavefield modeling as the database. Only a small number of wavefield-forward modeling is carried out based on the reciprocity theorem, that is, the number of wavefield modeling is directly proportional to the number of sensors instead of all searching grid points as candidate source. There are usually 10-50 sensors in a mine MS monitoring system, thus, it is easy to obtain the Green's function from each sensor as virtual source to each grid point. In accordance to the reciprocity theorem in scalar acoustic wave representation, Equation (2) can be rewritten as:  
\[u(\mathbf{x}_{r};\mathbf{x}_{s})= G(\mathbf{x}_{s};\mathbf{x}_{r})*S(t) \tag{3}\]  
where, \(G(\mathbf{x}_{s};\mathbf{x}_{r})\) denotes the Green's function from the sensor location \(\mathbf{x}_{r}\) to source location \(\mathbf{x}_{s}\). In a practical application, coarse grids are used in the regions far away from the mining area, and fine grids are used in the regions near the active mining area (where MS and blasting events occur frequently in practice). The combination of coarse and fine grids greatly reduce modeling computation and storage space of [PERSON]'s function database, which is called the multi-scale strategy. It is also conducive to improve the grid-search speed in the GSWI process. For the broadband waveform inversion-based location method, it needs to define a reasonable objective function to determine the fitting degree of the target observed waveform and theoretical modeling waveform. After that, inversion algorithm searches for the minimum value corresponding point as the source location. Due to the complex geological structures and velocity distribution in a mining area, the inversion problem will be highly nonlinear if the L2 norm misfit between the synthetic and observed waveforms is directly used as the objective function. This research utilized the L2 norm difference of cross-correlation travel time as the objective function of GSWI location study, which is more robust and linear for inversion problems. The cross-correlation function of two discrete waveform time series is defined as follows:  
\[f(\mathbf{x}_{r},\mathbf{u},\mathbf{d})=\sum_{t=-\infty}^{\infty}\mathbf{u}( \mathbf{x}_{r},t)\cdot\mathbf{d}(\mathbf{x}_{r},t+\tau) \tag{4}\]  
The cross-correlation travel time difference is determined by the maximum value of Equation (4) based on the first-order approximation, that is, \(\Delta\tau(\mathbf{x}_{r},\mathbf{u},\mathbf{d})=\frac{1}{N_{r}}\int_{0}^{T}w( t)\frac{\partial\mathbf{u}(\mathbf{x}_{r},t)}{\partial t}[\mathbf{d}(\mathbf{x}_{r},t)-\mathbf{u}(\mathbf{x}_{r},t)]dt\) where, \(N_{r}=\int_{0}^{T}w(t)\mathbf{u}(\mathbf{x}_{r},t)\frac{\partial\mathbf{u}( \mathbf{x}_{r},t)}{\partial t}dt\), \(w(t)\) is the time window function of specified seismic phases in the interval [0 T].