A coal mine microseismic positioning method based on riemannian manifold graph convolution network and physical constraint

Abstract

A method for microseismic location in coal mines based on Riemannian manifold graph convolutional networks and physical constraints, belonging to the field of coal mine microseismic prediction, comprises the following steps: Step 1) Creating a dataset; Step 2) Constructing a deep learning model based on a Riemannian manifold, including a front-end graph convolutional network (GCN) module and a back-end three-dimensional fully convolutional network (FCN) module; Step 3) Training the model: Training the obtained model using the Adam optimizer to obtain the final model; Step 4) Prediction: Inputting newly observed travel time difference datasets into the final model obtained in Step 3), propagating information through the graph convolutional network, constructing a three-dimensional spatial probability volume, and determining the most probable seismic source location from it. This invention improves the location accuracy through the above method, overcoming the long-standing contradiction between accuracy and efficiency.

Description

Technical Field

[0001] This invention belongs to the field of microseismic prediction in coal mining, specifically involving a coal mine microseismic location method based on Riemannian manifold convolutional networks and physical constraints. Background Technology

[0002] Common dynamic hazards during coal mining include mine tremors, gas explosions, and coal and gas outbursts. These hazards are mainly caused by factors such as stress concentration in the coal seam, rock fracturing, and gas escape. They not only pose a serious threat to the lives of miners but can also cause equipment damage and economic losses. Therefore, scientific prediction and monitoring, as well as effective preventive measures, are of paramount importance.

[0003] Connecting physical constraints with microseismic location in coal mines: Physical constraints directly embed the propagation laws of seismic waves into the deep learning model, enabling the model to more realistically reflect the geological structure of coal mines, effectively improving the model's location accuracy and providing more accurate and reliable location results.

[0004] Connecting Deep Learning with Coal Mine Microseismic Location: Traditional location methods struggle to balance computational efficiency and location accuracy, while deep learning, through end-to-end modeling and automated feature extraction, can quickly learn the spatial distribution patterns of coal mine microseismic events, achieving rapid and accurate location. This enables more effective identification of potential hazardous areas, providing strong support for mine safety monitoring and disaster early warning.

[0005] Existing technical solutions: Microseismic event localization methods mainly include traditional physical algorithms and artificial intelligence-based methods. Traditional physical algorithms typically rely on the accuracy of seismic wave first arrival time or velocity models, such as global grid search and nonlinear localization methods. Purely data-driven deep learning models built using artificial intelligence methods lack explicit physical constraints, and their performance significantly degrades when applied to unknown geological environments or velocity disturbance scenarios; examples include CNNs and FCNs. While traditional physical algorithms can achieve meter-level accuracy, each event localization takes several seconds to minutes. While artificial intelligence-based methods can achieve millisecond-level inference, they suffer from accuracy degradation, leading to a trade-off between accuracy and efficiency. Summary of the Invention

[0006] To overcome the shortcomings of existing technologies, this invention provides a coal mine microseismic location method based on Riemannian manifold convolutional networks and physical constraints.

[0007] The technical solution of this invention is: a coal mine microseismic location method based on Riemannian manifold graph convolutional network and physical constraints, the steps of which are as follows: Step 1) Create a dataset, the specific method is as follows: Step 1.1) Assume the study area is a three-dimensional geological body. Divide it into M×M×M cubic elements: (5) Each of them This represents a candidate earthquake source location, with the element size set to d, defining a three-dimensional spatial velocity model. Its value varies with spatial location Changes with the changes: (6) in, Represents background speed. This represents the velocity disturbance caused by geological heterogeneity; iterates through each cubic cell, takes the center of each cell as the simulation source, and uses the Fast March Method (FMM) to calculate the arrival time data of each station; Step 1.2) Calculate the propagation time difference between stations. As input data, a feature vector is constructed, and the location of the i-th candidate earthquake source is set as follows: The propagation times to the j-th and k-th stations are respectively and ,but: (7) Collect the travel time difference of all stations to form a correlation with the epicenter location. Corresponding feature vectors: (8) in, Indicates the number of stations; A three-dimensional probabilistic volume representation is used to model earthquake sources. For each training sample, the true source location is mapped to a three-dimensional Gaussian distribution centered at that location, as shown in the following formula: , (9) By extracting feature information from microseismic waveform data from multiple stations and fusing the spatial topological relationships between stations to construct training labels, the inversion process of microseismic source spatial coordinates is transformed into a conditional regression calculation problem. This is achieved by inputting observation features... Then, the neural network outputs the probability distribution corresponding to the earthquake source location. The coordinates of the earthquake source are used to characterize the spatial coordinates, and their calculation relationship is as follows: (10) in, This represents the probability of the location of the earthquake source given the observation data.

[0008] Step 2) Construct a deep learning model based on Riemannian manifolds: Establish a hybrid neural network architecture, adopt the inverse modeling method, train the neural network to establish an inverse mapping function, and convert the observed time difference into the location of the earthquake source; the hybrid neural network architecture includes a front-end graph convolutional GCN module and a back-end three-dimensional fully convolutional FCN module.

[0009] The aforementioned front-end graph convolutional GCN module is used to model the spatial relationships between stations and potential wave propagation paths, and to extract features that are sensitive to propagation.

[0010] The aforementioned three-dimensional fully convolutional FCN module is used to map the extracted graph features to a spatial probability volume representing the distribution of earthquake sources.

[0011] Step 2) includes: Step 2.1) Construct a graph convolutional network (GCN) based on Riemannian manifold embedding; for a given seismic event, construct input features by calculating the propagation time difference between each station and other stations to form a feature vector reflecting the relative arrival time relationship of each station; redistribute the features as the embedding representation of the corresponding station nodes to obtain the node feature matrix, and construct the inter-station graph structure based on the geodesic distance of the Riemannian manifold guided by the velocity model, and input the node features into the graph convolutional network for feature extraction.

[0012] The specific method is as follows: The Graph Convolutional Network (GCN) uses the standard spectral convolution formula, where the feature update rule for each layer is as follows: (11) in, It is the first The node feature matrix of the layer; It is the first Layer node feature matrix; adjacency matrix It is a station-to-station edge weight matrix constructed based on the geodesic distance of the Riemannian manifold guided by the velocity model; It is an adjacency matrix The degree matrix is ​​such that the diagonal elements are the sum of the edge weights between the node and all its neighbors; It is the first The trainable weight matrix of the layer is used to linearly map the feature dimension; function It is a non-linear activation function; A geodesic distance based on a Riemannian manifold, guided by a velocity model, is proposed to construct physically consistent inter-station edge weights. This method models the propagation path between stations as geodesics on the Riemannian manifold. The geodesic distance is obtained by solving the Eikonal equation and approximating the minimum propagation time using a fast-progression method. It is then mapped to edge weights through a kernel function, reflecting the influence of geological heterogeneity on the wave propagation path. These edge weights can be used in graph convolutional networks (GCNs) to integrate dynamic geometric information between stations and capture geological heterogeneity and wave propagation characteristics. Two edge weight construction methods are proposed: Euclidean distance connections and connections based on Riemannian manifold geodesic distances. Assume the geological body region is a three-dimensional spatial domain. The velocity model is The medium is modeled as a Riemannian manifold with a position-dependent metric. metric tensor Defined as: (12) in, The Kronecker delta (δ) is used to restrict the metric tensor to take non-zero values ​​only in the diagonal directions, indicating that the spatial directions are independent; the metric tensor Embedding velocity field information into the geometric structure makes the propagation path of the wave shorter in the high-velocity region, reflecting the variable characteristics of wave propagation. Based on the above metric tensor, the line elements on the corresponding Riemannian manifold can be obtained, and their squared distances are expressed as: (13) in, Representing a point in three-dimensional space Wave speed at that location; On a Riemannian manifold, the geodesic path between any two points refers to the shortest path given a metric. For any two stations... and Geodetic distance is defined as: (14) in, Indicates the integration path along the geodesic line; In the angle-size tensor method used, the above equation simplifies to: (15) in, Representing line elements in Euclidean space, a geodesic is equivalent to a physical path with the minimum propagation time; Theoretically, geodesics satisfy the following second-order ordinary differential equation, namely the geodesic equation: (16) in, These are Christopher notations, used to describe the curvature relations on Riemannian manifolds. Their calculation formula is: (17) In metric tensors Below, we obtain: (18) Further simplification of Christopher's symbol: (19) In a real three-dimensional heterogeneous velocity field, considering that only the first arrival of seismic waves is calculated, and that waves will continuously refract as they propagate in the non-uniform three-dimensional medium, an approximation strategy is adopted, introducing the Fast Mash Method (FMM) based on the travel time field for solution; this is achieved by constructing the minimum propagation time field from the source or station to various points in space. And its gradient norm describes the propagation characteristics of the wavefront in space, where Indicates the spatial location of the traveling field The gradient magnitude at point satisfies the following Eikonal equation: (20) From a certain station Starting from here, obtain the time field of the entire area. For the minimum propagation time at each point, the geodesic distance is approximately equal to the minimum cumulative propagation time between the two stations: (twenty one) Using the steepest descent path from Tracing back to To reconstruct the propagation path, the edge weights of the adjacency graph are constructed using geodesic distances, and defined in the form of a Gaussian kernel function: (twenty two) in, The typical propagation distance between stations is used to control weight decay; to further reflect the nonlinear effect of velocity changes on the propagation path, a velocity gradient term is introduced to enhance the sensitivity of the graph structure to geological heterogeneity. (twenty three) in, (twenty four) It is a hyperparameter that adjusts the strength of velocity change regularization, representing a weighted penalty for regions with significant velocity changes along the path, used to enhance the model's response to highly heterogeneous regions; Indicates the integral coefficient along the path; This represents the velocity gradient, i.e., the change of velocity with space. Let be the length of the infinitesimal line segment.

[0013] Step 2.2) Construct a three-dimensional fully convolutional FCN network; after completing graph-based feature extraction, map the node-level features output by GCN to a three-dimensional discrete space, and use the fully convolutional network as a decoding module to generate the source spatial probability volume, which is used to characterize the probability distribution of each candidate spatial location as a source.

[0014] The specific method is as follows: Given a set of propagation-related features extracted from GCN ,in It is the number of stations. It is the embedded dimension; First, aggregate these node features into a unified three-dimensional spatial representation; then design a spatial projection mechanism. ,in It is a discretized grid space. It is the number of feature channels; Each grid voxel corresponds to a candidate position. Its input features are constructed by interpolating or aggregating GCN feature information from nearby stations, and the weights are given by inverse geodesic distance: (25) in, It is a station and voxel center Geodetic distance between them By controlling the spatial influence scale, graph-based features can be smoothly embedded into the volumetric space. The weighting coefficients representing the characteristics of voxels and station nodes are calculated and normalized using geodesic distance; Indicates the first Feature vectors of each node; Subsequently, the fully convolutional network FCN obtains 3D tensors through convolutional layers: (26) in, Represents 3D convolution. and It is the first Learnable convolutional kernels and biases of layers, It is a non-linear activation function that ultimately outputs a feature tensor. The shape is M×M×M, which represents the predicted three-dimensional probability distribution in the discretized source space; Apply voxel-level softmax normalization to each voxel: (27) Observed propagation time difference characteristics Given the condition, generate a value at all candidate positions. probability volume on .

[0015] Step 2.3) Construct a composite loss function; For the source spatial probability volume output by the model, introduce distribution consistency loss, location accuracy loss and forward consistency loss based on wave propagation physics, weight and combine each loss term to form a total loss function, and use the total loss to backpropagate and update the network parameters to achieve model training and positioning accuracy optimization.

[0016] The specific method is as follows: A multi-component composite loss function was designed: (28) in, and These are weighting coefficients; Representing the KL divergence term, its objective is to make the predicted spatial distribution of earthquake sources... With real labels Consistency: (29) in, Represents the three-dimensional coordinates of the predicted earthquake source location; Three-dimensional coordinates representing the actual location of the earthquake source; Representing the mean squared error term, the point where the maximum value of the predicted probability volume is taken as the estimated source location, i.e.: (30) Calculate and estimate the location of the earthquake source Compared with the actual epicenter location The Euclidean distance between them is calculated, and the mean squared error loss term is constructed: (31) This allows for optimization of the loss term, directly constraining the positioning error and improving the accuracy of seismic source positioning. By introducing a forward modeling loss and using travel time residuals for inversion, the physical consistency of the model with the source location is enhanced. The theoretical travel time is obtained by forward modeling based on the velocity model by predicting the source location. It is a physical consistency term, characterized by its time difference with the input. By comparison, we obtain: (32) in, This represents the station calculated from the predicted location. and The time difference in propagation between them; Indicates the observed stations and The time difference in propagation between them; By reversing the predicted location to obtain the propagation path and matching it with the original input, reverse consistency and physical plausibility are enhanced.

[0017] Step 3) Train the model: Use the Adam optimizer to train the model obtained in step 2) to obtain the final model.

[0018] The specific method is as follows: During the training phase, forward modeling techniques were used to calculate the arrival time from each candidate source to all seismic stations. The Fast Mover Method (FMM) was used to solve the Eikonal equation to calculate the travel time difference between each station and all other stations, and feature vectors were constructed. All forward modeling results were pre-calculated and cached before training. At the same time, the Gaussian volume representation corresponding to the location of each candidate source was also pre-discretized. The Adam optimizer was used for training. The initial learning rate was set to 10⁻³, and an exponential decay strategy was adopted, with the learning rate decaying once every 10 epochs. The decay factor was 0.9, the batch size was set to 32, and a total of 200 epochs of training were performed. The overall training process includes the following steps: sample selection, feature construction through time difference, graph convolution processing, spatial mapping to generate 3D output, calculation of loss function, and subsequent optimization through backpropagation and parameter update.

[0019] Step 4) Prediction: Input the newly observed travel time difference dataset into the final model obtained in Step 3), propagate information through a graph convolutional network, construct a three-dimensional spatial probability volume, and determine the location of the most likely earthquake source from it.

[0020] The beneficial effects of this invention are as follows: It tightly integrates known knowledge of earthquake physics into a deep neural network. By combining Graph Convolutional Networks (GCNs) with Fully Convolutional Networks (FCNs) and introducing three layers of physical constraints (structural constraints, training constraints, and loss constraints), the seismic wave propagation mechanism is embedded into the deep learning framework, achieving accurate modeling of the source probability distribution. The high positioning accuracy and fast inference capability of RiemannianGF have been verified through numerical simulations and real-world data testing, demonstrating superior positioning speed and accuracy compared to traditional positioning algorithms. With sufficient training samples, RiemannianGF has the potential to achieve high-precision classification, with direct positioning errors reaching the meter level, and exhibits strong robustness against velocity perturbations and complex geological structures. Attached Figure Description

[0021] Figure 1 This is a schematic diagram illustrating the implementation process of the method of the present invention; Figure 2 This is a diagram of the Riemannian GF model architecture. Figure 3 This is a schematic diagram of two edge weight construction methods; Figure 4 This is a schematic diagram of the propagation of seismic waves in a three-dimensional heterogeneous medium. Figure 5 This is a schematic diagram of the numerical simulation model; where red dots represent stations. Figure 6 A schematic diagram of the model architecture for the numerical simulation experiment of the Marmousi model; Figure 7 A graph showing the changes in loss function and accuracy; Figure 8 This is a histogram of the positioning error distribution. Figure 9 A scatter plot comparing the predicted results with the actual results; Figure 10 This is a distribution map of the microseismic dataset; Among them, (a) is the distribution map of the stations of Dongtan Coal Mine, (b) is the distribution map of microseismic events of the 6306 working face of Dongtan Coal Mine, and (c) is the geological information map of the 6306 working face of Dongtan Coal Mine. Figure 11 A comparison of training loss and validation accuracy between RiemannianGF and LGF; Figure 12 Map showing the location error distribution at Dongtan Coal Mine; Figure 13 A schematic diagram of the predicted probability volume and slices for representative synthetic events; Figure 14 This diagram illustrates the location error of the blasting event using three different methods. Detailed Implementation

[0022] A method for locating microseismic events in coal mines based on Riemannian graph convolutional networks and physical constraints is proposed, which includes a deep learning model based on a Riemannian manifold, named RiemannianGF, aimed at improving the accuracy of coal mine microseismic event location. We first embed detailed geological structure data and construct a refined three-dimensional coal mine velocity model that can capture faults, coal seam fluctuations, and lithological changes. Then, we propose RiemannianGF, a hybrid graph convolutional network (GCN) – fully convolutional network (FCN) model that embeds the coal mine velocity field as a Riemannian manifold. The model integrates three physical constraints at each stage of learning: (1) structural constraints, where graph edges are geological propagation time distances that preserve the physical wave propagation path; (2) training constraints, where the supervision signal comes from synthetic propagation time pairs generated through forward modeling; and (3) loss constraints, where the propagation time residual term is regularized and optimized. By aligning the architecture, data, and targets with seismic physics, RiemannianGF improves positioning accuracy, interpretability, and robustness while maintaining millisecond-level inference speed, thus overcoming the long-standing trade-off between accuracy and efficiency.

[0023] In the problem of microseismic location, seismic waves originate from the source point. Propagation occurs through a medium, reaching multiple receiving stations (rj). The propagation process is determined by a geological velocity model. Control. Given a velocity model, the source location uniquely determines the propagation time at each station. This further uniquely determined the propagation time difference between the stations: (1) This process is essentially a forward mapping determined by the laws of physics: (2) in, Indicates the source location The corresponding propagation time difference vector, function Subject to velocity model The control reflects the influence of wave propagation path and geological heterogeneity on propagation time. Traditional grid search methods are based on this relationship: by traversing all candidate locations, the corresponding theoretical propagation time difference is calculated. The difference between the observed propagation time and the actual propagation time. Compare the results and select the location with the smallest error as the estimated source location: (3) Although this method is accurate, it is computationally intensive, especially in complex three-dimensional velocity fields, requiring multiple forward modeling iterations.

[0024] In contrast, the deep learning model proposed in this study employs a reverse modeling approach, training a neural network to establish an inverse mapping function that converts the observed time difference into the location of the earthquake source: (4) Specifically, it bypasses the traditional traversal and matching process, predicting the source location directly from the input features. This modeling paradigm not only improves computational efficiency but also opens up new avenues for the introduction of physical constraints.

[0025] Our approach achieves physically consistent modeling through three mechanisms: (1) Structural constraints: Graph construction on Riemann manifolds.

[0026] We propose a non-Euclidean graph architecture that uses edge weights of a geodesic distance-defined graph convolutional network (GCN) induced by a velocity model to dynamically capture wave propagation paths. This structure simultaneously reflects the distribution of space stations and rigorously characterizes seismic wave propagation dominated by the velocity field.

[0027] (2) Training constraints: Physics-based synthetic training data.

[0028] All training samples were generated by solving the Eikonal equations, ensuring that the time difference characteristics of the observations inherently incorporate wave propagation physics. These graph edges are explicitly derived from geodesic paths reconstructed from the forward-modeled wavefront.

[0029] (3) Loss constraint: Physical regularization loss function.

[0030] In addition to the conventional KL divergence and position error supervised loss, we also introduce a penalty term based on the temporal difference mismatch residual of the source position prediction. This forces the model not only to learn the position prediction but also to grasp the physical plausibility through inverse wavefield reconstruction.

[0031] Through this multi-physics embedding scheme, we developed a physics-driven neural model for microseismic source localization, achieving high accuracy and enhanced robustness under complex geological heterogeneity.

[0032] The specific method for microseismic location of coal mines using the above model is as follows: Figure 1 As shown, it is divided into dataset creation, model building, model training and prediction, and experimental verification.

[0033] Step 1: Dataset creation.

[0034] Assume the study area is a three-dimensional geological body Divide it into M×M×M cubic elements: (5) Each of them This represents a candidate seismic source location. The element size is set to d (e.g., 10 meters) to balance resolution and computational load. A three-dimensional spatial velocity model is defined. Its value varies with spatial location Changes with the changes: (6) in, Represents background speed. This represents velocity disturbances caused by geological heterogeneity. We traverse each cubic cell, taking the center of each cell as the simulation source, and use the Fast March Method (FMM) to calculate the arrival time data for each station.

[0035] Then, by calculating the propagation time difference between stations As input data, a feature vector is constructed. Let the location of the i-th candidate earthquake source be... The propagation times to the j-th and k-th stations are respectively and ,but: (7) Collect the travel time difference of all stations to form a correlation with the epicenter location. Corresponding feature vectors: (8) in, Indicates the number of stations.

[0036] Traditional positioning methods provide a single-point source location estimate. However, in situations with data uncertainty or unclear distribution, such as low-precision time-of-flight measurements or sparse station coverage, single-point regression can lead to positioning errors with skipped or discontinuous source locations.

[0037] To address this issue, we employ a three-dimensional probabilistic volume representation to model earthquake sources. Specifically, for each training sample, the true source location is mapped to a three-dimensional Gaussian distribution centered at that location. The formula is as follows: , (9) This tag structure design is essentially equivalent to defining a posterior distribution estimation problem for microseismic location. Given observation characteristics... We hope to output via the network The true posterior distribution that can approximate the epicenter location: (10) in, This represents the probability of the source location given the observation data. This modeling approach allows the network to more strongly express spatial uncertainty and is more consistent with the Bayesian source location concept used in physical inversion.

[0038] Step 2: Model building.

[0039] To effectively integrate the spatial configuration of stations, the characteristics of wave propagation paths, and the spatial distribution of seismic sources, we designed a hybrid Graph Convolutional Network (GCN) – Fully Convolutional Network (FCN) neural network architecture. This framework consists of two key components: The front-end graph convolutional module (GCN) is used to model the spatial relationships between stations and potential wave propagation paths, and extract propagation-sensitive features; The backend 3D fully convolutional module (FCN) is used to map the extracted graph features to a spatial probability volume representing the distribution of earthquake sources.

[0040] This architecture combines the relational modeling capabilities of graph neural networks with the spatial continuity learning advantages of convolutional networks, representing a structure-aware design specifically tailored for microseismic localization tasks in complex underground environments. A schematic diagram of the model architecture is shown below. Figure 2 As shown.

[0041] Step 2.1) Graph Convolutional Network (GCN) based on Riemannian manifold embedding

[0042] In the graph structure, nodes represent seismic stations. For a given seismic event, the input features are constructed by calculating the propagation time difference between each station and all other stations, forming a feature vector that captures the relative arrival times between stations in that event. To align with the graph representation, these features are reassigned as individual node embeddings, resulting in a node feature matrix.

[0043] Graph Convolutional Networks (GCNs) employ the standard spectral convolution formula, where the feature update rule for each layer is as follows: (11) Among them, the adjacency matrix It is based on geodesic distances defined on Riemannian manifolds guided by velocity models.

[0044] Early studies on GCN-style localization used Euclidean distance to construct the graph structure, with edge weights reflecting geometric proximity. However, in geologically heterogeneous media, seismic wave propagation is strongly influenced by spatial velocity variations, leading to nonlinear curvature of the propagation path, which Euclidean distance cannot accurately capture.

[0045] To address this limitation, we propose a geodesic distance based on a Riemannian manifold, guided by a velocity model, to construct a physically consistent inter-station graph structure. This method models the propagation paths between stations as geodesics on the Riemannian manifold, enabling the GCN to integrate dynamic geometry and reflect the influence of geological heterogeneity and wave propagation characteristics, rather than relying solely on static spatial distances. Figure 3 Two schematic diagrams illustrating edge weight construction methods are shown. (a) is a graph connected by Euclidean distance; (b) is a graph connected by geodesic distance based on Riemannian manifolds.

[0046] Assume the geological body region is a three-dimensional spatial domain. The velocity model is We model the medium as a Riemannian manifold with a position-dependent metric. metric tensor Defined as: (12) in, The Kroneck δ (i.e., the identity matrix) indicates that measurements in each direction are independent. This metric tensor embeds velocity field information into the geometry, resulting in shorter propagation paths for waves in high-velocity regions and naturally reflecting the variable characteristics of wave propagation.

[0047] The derived line element (the square of the distance) is: (13) In high-speed regions, the physical distance per unit length is shorter, which is consistent with the physical intuition that seismic waves travel at a relatively fast speed.

[0048] On a Riemannian manifold, the geodesic path between any two points refers to the shortest path given a metric. For any two stations... and Geodetic distance is defined as: (14) In the angle-size tensor method used in this study, the above equation simplifies to: (15) in, A geodesic represents a line element in Euclidean space. This definition is essentially equivalent to Fermat's principle (that waves travel along the path with the shortest propagation time), therefore a geodesic is equivalent to a physical path with the shortest propagation time.

[0049] Theoretically, geodesics satisfy the following second-order ordinary differential equation (geodesic equation): (16) in, It is the Christopher symbol, and the calculation formula is: (17) The metric tensor we use The following simplified form was obtained: (18) Further simplification of Christopher's symbol: (19) Although the equations are complete in form, numerical solutions are extremely challenging in real three-dimensional heterogeneous velocity fields. Considering that only the first arrival of seismic waves is calculated, the waves undergo continuous refraction as they propagate in the inhomogeneous three-dimensional medium, such as... Figure 4 As shown.

[0050] Therefore, we adopt an approximation strategy, using the Fast Progress Method (FMM) to solve the Eikonal equations: (20) From a certain station Starting from here, obtain the time field of the entire area. This corresponds to the minimum propagation time at each point. The geodesic distance is approximately equal to the minimum cumulative propagation time between two stations. (twenty one) In addition, the steepest descent path can be used from Tracing back to The transmission path can be roughly reconstructed.

[0051] Finally, the edge weights of the adjacency graph are constructed using geodesic distances and defined in the form of a Gaussian kernel function: (twenty two) in, The decay scale for controlling the weights is set to the typical propagation distance between stations. To further reflect the nonlinear effect of velocity variations on the propagation path, we introduce a velocity gradient term to enhance the sensitivity of the graph structure to geological heterogeneity: (twenty three) in, (twenty four) λ is a hyperparameter that adjusts the strength of the velocity variation regularization. This term represents a weighted penalty for regions with significant velocity variations along the path, enhancing the model's response to highly heterogeneous regions such as faults and interfaces.

[0052] Step 2.2) Fully Convolutional Network (FCN)

[0053] Following graph-based feature extraction, a fully convolutional network (FCN) is used as a decoding module to map the node-level features learned from the GCN to a dense spatial probability volume, representing the probability of the microseismic source location.

[0054] Given a set of propagation-related features extracted from GCN ,in It is the number of stations. This involves embedding the dimension. First, these node features are aggregated into a unified three-dimensional spatial representation. To achieve this, we designed a spatial projection mechanism. ,in It is a discretized grid space. It is the number of feature channels.

[0055] Each grid voxel corresponds to a candidate position. Its input features are constructed by interpolating or aggregating GCN feature information from nearby stations, and the weights are given by inverse geodesic distance: (25) in, It is a station and voxel center Geodetic distance between them Controlling the spatial influence of scale. This step enables the smooth embedding of graph-based features into the volumetric space.

[0056] Subsequently, the 3D tensor obtained by the fully convolutional network (FCN) through a series of convolutional layers is: (26) in, Represents 3D convolution. and It is the first Learnable convolutional kernels and biases of layers, It is a non-linear activation function (such as ReLU). The final output feature tensor The shape is M×M×M, which represents the predicted three-dimensional probability distribution in the discretized source space.

[0057] To ensure probabilistic interpretability, we apply voxel-level softmax normalization to each voxel: (27) This will generate a list of all candidate locations. probability volume on And based on the observed propagation time difference characteristics This is a condition. This probability formula allows the model to naturally express spatial uncertainty and provides a theoretical basis for source location estimation under data ambiguity.

[0058] In summary, the FCN module decodes high-level relational features into a spatially dense posterior distribution, capturing the physical propagation properties learned by GCN and the continuity of the source location field. Its fully convolutional structure also enables efficient inference over large-scale 3D regions without the need for manual post-processing.

[0059] Step 2.3) Loss Function

[0060] To improve the accuracy, stability, and physical consistency of the model, this study designed a multi-component composite loss function: (28) in, and It is the weighting coefficient.

[0061] Representing the KL divergence term, its objective is to make the predicted spatial distribution of earthquake sources... With real labels Consistency: (29) This represents the mean squared error term, designed to improve the accuracy of the hypocenter location. The point of maximum predicted probability is selected as the estimated hypocenter location. (30) And apply Euclidean distance penalty based on the actual epicenter location: (31) It directly optimizes the positioning error and improves the accuracy of seismic source positioning.

[0062] By introducing a forward modeling loss and using travel time residuals for inversion, the physical consistency of the model with respect to the source location is enhanced. In other words, by predicting the source location, the theoretical travel time is obtained through forward modeling based on the velocity model. It is a physical consistency term. This is achieved through the time difference feature with the input. By comparison, we obtain: (32) The propagation path is obtained by reversing the predicted position and matching it with the original input to enhance reverse consistency and physical plausibility. The table below summarizes the objectives and effects of various loss functions.

[0063] Table 1: Objectives and Effects of Various Loss Functions

[0064] Step 3: Model Training and Prediction

[0065] During the training phase, we used forward modeling techniques to calculate the arrival time from each candidate source to all seismic stations. The Fast Streaming Model (FMM) was employed to solve the Eikonal equations to calculate the travel time difference between each station and all other stations, constructing feature vectors. To improve computational efficiency, all forward modeling results were pre-calculated and cached before training, and the Gaussian volume representation corresponding to each candidate source location was also pre-discretized.

[0066] We used the Adam optimizer for training, with an initial learning rate of 10⁻³ and an exponential decay strategy, where the learning rate decreases every 10 epochs with a decay factor of 0.9. The batch size was set to 32, and training lasted for 200 epochs.

[0067] The overall training process includes the following steps: sample selection, feature construction through time difference, graph convolution processing, spatial mapping to generate 3D output, calculation of loss function, and subsequent optimization through backpropagation and parameter update.

[0068] After model training is complete, the prediction phase begins. The goal of the prediction phase is to efficiently predict the source location of microseismic events based on the observed travel time differences. For newly observed travel time difference datasets, the model propagates information through a graph convolutional network to construct a three-dimensional spatial probability volume and determine the most probable source location from it. Table 2 below summarizes the model training and prediction process in detail.

[0069] Table 2: Training and Prediction Flowchart

[0070] Step 4: Experimental Verification

[0071] I. Numerical Simulation

[0072] To quantitatively evaluate the performance of the proposed RiemannianGF model, we designed a numerical localization experiment based on the three-dimensional Marmousi velocity model (Versteeg, 1994; Martin et al., 2006). An event distribution area of ​​700 m × 800 m × 200 m was divided into a 70 × 80 × 40 grid, with each grid cell measuring 10 m (north-south), 10 m (east-west), and 5 m (depth). Eight seismic stations were placed at the vertices of the cuboid to ensure maximum spatial coverage, as shown in Figure 5.

[0073] A total of 224,000 synthetic microseismic events were uniformly distributed across each node of a 70×80×40 Marmousi grid. For each event, the accurate first arrival time to all stations was calculated using the Fast March Method (FMM) without any added perturbations, ensuring that the location error only reflects the algorithm itself. For the grid search method, each grid node in the 3D volume was considered a candidate source. The pre-calculated FMM time difference table for that node was compared with the actual FMM time of the active event, and the candidate point with the smallest root mean square residual was selected as the estimate for the grid search. For RiemannianGF, the time differences between stations from the same FMM results were extracted and used as input features. The graph structure was constructed based on geodesic distances guided by the velocity model, consistent with the model training process. A detailed schematic diagram of the model architecture is shown in Figure 6.

[0074] LGF SeismoLocator is a purely data-driven network that directly takes three-component waveform data as input; while RiemannianGF utilizes the differential propagation time between stations, whereas traditional grid search inversion relies on the absolute arrival time picking for each station. We first created a general, noise-free synthetic dataset using the SW4 wave propagation solver. Complete seismic waveforms were simulated at each grid node, with identical source and receiver geometry and velocity models. These waveform data were then processed through the LGF preprocessing pipeline, including normalization, time window selection, and spectrogram generation, enabling all three methods to be benchmarked in physically consistent scenarios while preserving their original input formats.

[0075] The training curves in Figure 7 show that RiemannianGF converges faster and reaches a lower final loss than LGF, indicating that embedding geodesic constraints based on the velocity model into the network architecture improves optimization efficiency. At the end of training (200 epochs), the loss of RiemannianGF drops to 0.07 with a validation accuracy of 0.95, while LGFSeismoLocator stagnates at a higher loss of 0.15 with an accuracy of 0.89. The validation accuracy curve also shows that RiemannianGF consistently outperforms LGF throughout the training process, reflecting its stronger generalization ability and capacity to adapt to unseen synthetic events. The error distribution in Figure 8 further shows that the error distribution of RiemannianGF is more concentrated, with fewer large outliers, demonstrating better robustness and stability. Finally, the scatter plot in Figure 9 confirms that the predictions of RiemannianGF are more closely aligned with the actual source locations in all three spatial dimensions, while LGF shows a larger deviation, especially in the depth direction.

[0076] Table 3 summarizes the localization experimental results of all methods. RiemannianGF has significantly lower mean error and maximum error compared to LGF, while maintaining millisecond-level inference time comparable to LGF and being several orders of magnitude faster than grid search.

[0077] Table 3. Results of the Positioning Experiment

[0078] Comprehensive numerical simulation results show that RiemannianGF effectively narrows the accuracy gap with the ideal location, surpasses LGF, and provides a high-precision and efficient source location solution.

[0079] Example 1: To verify the performance of the Riemannian GF model in practical applications, we selected the 6306 working face of the Dongtan Coal Mine in Shandong Province, China, as the field test area. After ignoring phase picking errors in the numerical simulation experiment, we evaluated 241 marked microseismic events (-0.26 ≤ MW ≤ 1.73) on this working face.

[0080] Figure 10(a) shows an overview of the Dongtan Coal Mine and the layout of the KJ874(A) microseismic monitoring system. This system integrates surface three-component stations (TXL, FPZ), deep-hole three-component stations (XFJ, CCH), underground three-component stations (HDS, SCJD, BDS), and single-component seismographs (WDZ-3, 4, 5, 6, 8, 10). These sensors record the event catalog used in our study.

[0081] The corresponding epicenter locations are shown in Figure 10(b), distributed within a 700 m × 800 m × 200 m region centered on the 6306 working face, with a range of Y = [20489900, 20490700] m and X = [3921100, 3921800] m, and a depth of –500 m to –700 m. This region was discretized into 10 m × 10 m × 5 m voxels, resulting in a 70 × 80 × 40 grid (224,000 elements), providing the same number of training samples for RiemannianGF.

[0082] For reference, Figure 10(c) shows the local stratigraphic structure of the 6306 working face, highlighting faults, coal seam undulations, and lithological stratification that affect wave propagation. All events were initially located using a grid search algorithm. The grid search method had a mean error of 6.20 meters, but required approximately 12.6 seconds per event, which was insufficient for real-time applications. LGF reduced the inference time to approximately 0.001 seconds per event, but its accuracy decreased to 7.40 meters, revealing a trade-off between speed and accuracy. To address this issue, we applied the RiemannianGF method to the Dongtan Coal Mine dataset.

[0083] After establishing the field geometry, velocity model, and training protocol for working face 6306, we applied the calibrated model to a catalog of 241 events recorded between September 15 and November 15, 2023. The mesh resolution and hyperparameters remained consistent with the numerical simulation experiments; only the edge weights of the updated velocity model and the actual station geometry were recalculated. Figures 11-12 visualize the training and prediction results, and Table 4 summarizes the localization errors.

[0084] The training dynamics of RiemannianGF are as follows: Figure 11 As shown, its convergence behavior over 200 training epochs is illustrated. RiemannianGF not only achieves a lower final training loss (0.02), but also exhibits smoother and faster convergence than LGFSeismoLocator, whose loss stagnates around 0.06. In terms of validation accuracy, RiemannianGF achieves 0.94, consistently outperforming LGFSeismoLocator, which peaks at 0.88 during training, indicating that RiemannianGF has stronger generalization ability. These results demonstrate the benefits of directly embedding three-layer physical constraints, map structure, synthetic forward model training supervision, and time-delayed residual regularization into the model optimization process.

[0085] Figure 12 shows a comparative positioning error distribution of two methods (LGF SeismoLocator and RiemannianGF) based on the actual microseismic dataset from the Dongtan Coal Mine. Figures (a), (b), and (c) in Figure 12 show histograms of the horizontal, depth, and absolute positioning errors for 230 test data points, respectively. These figures demonstrate that RiemannianGF consistently outperforms LGF, particularly in terms of error distribution. The error distribution of RiemannianGF exhibits more concentrated peaks at lower error values, indicating higher positioning accuracy and stability. In contrast, LGF shows a wider error distribution, indicating poorer accuracy. Figure 12 (d) shows a box plot of the absolute positioning error, highlighting the larger range of RiemannianGF across the interquartile range, while the median position of RiemannianGF is more stable compared to LGF. These observations suggest that RiemannianGF provides better robustness in real-world applications and is less sensitive to outliers, which is particularly important in practical applications where accuracy and consistency are critical.

[0086] Table 4: Positioning Performance Table of Dongtan Coal Mine

[0087] Table 4 summarizes the localization performance of the two methods on the Dongtan Coal Mine dataset. The LGF SeismoLocator has a mean error of 8.68 meters, a maximum error of 22.36 meters, and an inference time of 0.0010 seconds per event. RiemannianGF is slightly more accurate, with a mean error of 7.85 meters and a maximum error of 20.62 meters, and is also faster, with a processing time of 0.0008 seconds per event. These results demonstrate that RiemannianGF outperforms LGF in both efficiency and accuracy, making it a highly efficient choice for practical coal mine safety monitoring.

[0088] Despite RiemannianGF's strong localization performance in both synthetic and field scenarios, its error in the Dongtan coal mine experiment was significantly greater than that in the Marmousi simulation. This increase is primarily due to the lack of an accurate velocity model in the field, which only provides a simplified hierarchical approximation. This inaccuracy introduces a systematic bias in the arrival estimation. Furthermore, arrival picking errors (whether manual or automatic) further degrade the quality of the input features, increasing the uncertainty in localization. Compared to LGF SeismoLocator, RiemannianGF achieved lower mean and maximum errors on the Dongtan dataset (Table 4). LGF is trained on historical locations generated by a grid search based on the velocity model, thus implicitly embedding its labels into the velocity field; however, LGF does not re-evaluate the velocity model during inference. Both methods output Gaussian probability volumes, but LGF primarily derives its estimation through an LSTM-based temporal sequence representation, which is more sensitive to waveform noise. RiemannianGF obtains a Gaussian posterior with stronger spatial consistency and higher accuracy by explicitly embedding the velocity model at the edge of the map and adding time-delay residual loss, while maintaining robustness under velocity perturbations and noisy inputs.

[0089] RiemannianGF addresses the long-standing trade-off between accuracy and speed in microseismic location. On the Dongtan Coal Mine dataset, it reduces the mean error from 8.68 meters to 7.85 meters, the maximum error from 22.36 meters to 20.62 meters, and the inference time from 0.0010 seconds per event to 0.0008 seconds. By providing higher accuracy and faster computation speed than LGFSeismoLocator, far surpassing slower grid search methods, RiemannianGF sets a new technical standard for real-time monitoring in coal mine environments.

[0090] III. Explosion Experiment

[0091] To quantitatively evaluate the source location accuracy of RiemannianGF under actual mining conditions, we designed a rooftop deep-hole blasting experiment. In the 6306 working face of the Dongtan Coal Mine, the model's location performance was verified by comparing the source location with known blasting locations.

[0092] From November 23 to November 28, 2023, we conducted eight deep-hole blasting experiments on the roof of the key layer of the transport roadway in the 6306 working face of the Dongtan Coal Mine. The blasting locations are listed in Table 5. Table 5: Information Table of Blasting Incidents at Dongtan Coal Mine

[0093] We applied the RiemannianGF model trained in the previous stage to the blasting events at these 8 known locations for validation. Figure 13 The image shows a slice view and 3D Gaussian volume result of the first event obtained using the RiemannianGF and LGF methods. Figure 14 This demonstrates the absolute positioning errors of the three methods for eight blasting events.

[0094] Figure 13 In the table, (a) to (d) show the localization results for Riemannian GF, and (e) to (h) show the localization results for LGF. Figure 13 As can be seen, the RiemannianGF slices exhibit a tightly concentrated, almost symmetrical Gaussian peak centered at the true event location, reflecting high localization accuracy and isotropic spatial uncertainty. The corresponding 3D Gaussian volume is compact and centered, indicating low variance across the three dimensions. In contrast, the LGF SeismoLocator results show a wider probability distribution, significantly stretched in some planes, indicating anisotropic uncertainty and lower spatial resolution. The LGF's 3D Gaussian volume is larger and stretched, indicating greater overall localization uncertainty. These differences highlight the accuracy and robustness of RiemannianGF in capturing true event locations, with smaller errors and a more uniform spatial distribution.

[0095] Figure 14 In the diagram, the solid gray line represents the positioning error of the grid search method, the dashed gray line represents the computation time of the grid search method, the solid blue line represents the positioning error of the Riemannian GF method, the dashed blue line represents the computation time of the Riemannian GF method, the solid red line represents the positioning error of the LGF method, and the dashed red line represents the computation time of the LGF method.

[0096] In the blasting calibration experiment at the Dongtan Coal Mine, the average positioning error of RiemannianGF was 6.80 m, which is roughly equivalent to that of the grid search method (6.20 m), but far exceeds that of LGF SeismoLocator (7.40 m). Furthermore, the inference time per event was only 0.0008 seconds, demonstrating a significant improvement in computational efficiency. This method also exhibits higher stability than both LGF and the grid search method. In events 2 and 3, the accuracy of RiemannianGF surpassed that of the two comparative methods, while in events 7 and 8, the error of LGF deviated significantly from the results of the grid search method, and the results of RiemannianGF closely aligned with the positioning error of the grid search method. These results indicate that RiemannianGF combines near-optimal accuracy with excellent robustness and efficiency.

[0097] In the East Beach coal mine dataset, RiemannianGF reduced the average error from 8.68 m to 7.85 m and the maximum error from 22.36 m to 20.62 m, while simultaneously reducing inference time from 0.0010 seconds to 0.0008 seconds per event, demonstrating higher accuracy and faster computation speed. By providing higher accuracy and faster computation speed than LGF SeismoLocator, RiemannianGF surpasses traditional grid search methods and sets a new standard for real-time monitoring.

Claims

1. A method for microseismic localization of coal mines based on Riemannian manifold graph convolutional networks and physical constraints, characterized in that, The steps are as follows: Step 1) Create a dataset; Step 2) Construct a deep learning model based on Riemannian manifolds: Establish a hybrid neural network architecture, adopt the inverse modeling method, train the neural network to establish an inverse mapping function, and convert the observed time difference into the location of the earthquake source; the hybrid neural network architecture includes a front-end graph convolutional GCN module and a back-end three-dimensional fully convolutional FCN module. The aforementioned front-end graph convolutional GCN module is used to model the spatial relationships between stations and potential wave propagation paths, and to extract propagation-sensitive features; The aforementioned three-dimensional fully convolutional FCN module is used to map the extracted graph features to a spatial probability volume representing the distribution of earthquake sources; Step 3) Train the model: Use the Adam optimizer to train the model obtained in Step 2) to obtain the final model; Step 4) Prediction: Input the newly observed travel time difference dataset into the final model obtained in Step 3), propagate information through a graph convolutional network, construct a three-dimensional spatial probability volume, and determine the location of the most likely earthquake source from it.

2. The coal mine microseismic location method based on Riemannian manifold graph convolutional network and physical constraints according to claim 1, characterized in that, The specific method in step 1) is as follows: Step 1.1) Assume the study area is a three-dimensional geological body. Divide it into M×M×M cubic elements: (5) Each of them This represents a candidate earthquake source location, with the element size set to d, defining a three-dimensional spatial velocity model. Its value varies with spatial location Changes with the changes: (6) in, Represents background speed. This represents the velocity disturbance caused by geological heterogeneity; iterates through each cubic cell, takes the center of each cell as the simulation source, and uses the Fast March Method (FMM) to calculate the arrival time data of each station; Step 1.2) Calculate the propagation time difference between stations. As input data, a feature vector is constructed, and the location of the i-th candidate earthquake source is set as follows: The propagation times to the j-th and k-th stations are respectively and ,but: (7) Collect the travel time difference of all stations to form a correlation with the epicenter location. Corresponding feature vectors: (8) in, Indicates the number of stations; A three-dimensional probabilistic volume representation is used to model earthquake sources. For each training sample, the true source location is mapped to a three-dimensional Gaussian distribution centered at that location, as shown in the following formula: , (9) By extracting feature information from microseismic waveform data from multiple stations and fusing the spatial topological relationships between stations to construct training labels, the inversion process of microseismic source spatial coordinates is transformed into a conditional regression calculation problem. This is achieved by inputting observation features... Then, the neural network outputs the probability distribution corresponding to the earthquake source location. The coordinates of the earthquake source are used to characterize the spatial coordinates, and their calculation relationship is as follows: (10) in, This represents the probability of the location of the earthquake source given the observation data.

3. The coal mine microseismic location method based on Riemannian manifold graph convolutional network and physical constraints according to claim 1, characterized in that, Step 2) includes: Step 2.1) Construct a graph convolutional network (GCN) based on Riemannian manifold embedding; for a given seismic event, construct input features by calculating the propagation time difference between each station and other stations to form a feature vector reflecting the relative arrival time sequence of each station; redistribute the features as the embedding representation of the corresponding station nodes to obtain the node feature matrix, and construct the inter-station graph structure based on the geodesic distance of the Riemannian manifold guided by the velocity model, and input the node features into the graph convolutional network for feature extraction; Step 2.2) Construct a three-dimensional fully convolutional FCN network; After completing graph-based feature extraction, map the node-level features output by GCN to a three-dimensional discrete space, and use the fully convolutional network as a decoding module to generate the source spatial probability volume, which is used to characterize the probability distribution of each candidate spatial location as a source. Step 2.3) Construct a composite loss function; For the source spatial probability volume output by the model, introduce distribution consistency loss, location accuracy loss and forward consistency loss based on wave propagation physics, weight and combine each loss term to form a total loss function, and use the total loss to backpropagate and update the network parameters to achieve model training and positioning accuracy optimization.

4. The coal mine microseismic location method based on Riemannian manifold graph convolutional network and physical constraints according to claim 3, characterized in that, In step 2.1), the specific method is as follows: The Graph Convolutional Network (GCN) uses the standard spectral convolution formula, where the feature update rule for each layer is as follows: (11) in, It is the first The node feature matrix of the layer; It is the first Layer node feature matrix; adjacency matrix It is a station-to-station edge weight matrix constructed based on the geodesic distance of the Riemannian manifold guided by the velocity model; It is an adjacency matrix The degree matrix is ​​such that the diagonal elements are the sum of the edge weights between the node and all its neighbors; It is the first The trainable weight matrix of the layer is used to linearly map the feature dimension; function It is a non-linear activation function; A geodesic distance based on a Riemannian manifold, guided by a velocity model, is proposed to construct physically consistent inter-station edge weights. This method models the propagation path between stations as geodesics on the Riemannian manifold. The geodesic distance is obtained by solving the Eikonal equation and approximating the minimum propagation time using a fast-progression method. It is then mapped to edge weights through a kernel function, reflecting the influence of geological heterogeneity on the wave propagation path. These edge weights can be used in graph convolutional networks (GCNs) to integrate dynamic geometric information between stations and capture geological heterogeneity and wave propagation characteristics. Two edge weight construction methods are proposed: Euclidean distance connections and connections based on Riemannian manifold geodesic distances. Assume the geological body region is a three-dimensional spatial domain. The velocity model is The medium is modeled as a Riemannian manifold with a position-dependent metric. metric tensor Defined as: (12) in, The Kronecker delta (δ) is used to restrict the metric tensor to take non-zero values ​​only in the diagonal directions, indicating that the spatial directions are independent; the metric tensor Embedding velocity field information into the geometric structure makes the propagation path of the wave shorter in the high-velocity region, reflecting the variable characteristics of wave propagation. Based on the above metric tensor, the line elements on the corresponding Riemannian manifold can be obtained, and their squared distances are expressed as: (13) in, Representing a point in three-dimensional space Wave speed at that location; On a Riemannian manifold, the geodesic path between any two points refers to the shortest path given a metric. For any two stations... and Geodetic distance is defined as: (14) in, Indicates the integration path along the geodesic line; In the angle-size tensor method used, the above equation simplifies to: (15) in, Representing line elements in Euclidean space, a geodesic is equivalent to a physical path with the minimum propagation time; Theoretically, geodesics satisfy the following second-order ordinary differential equation, namely the geodesic equation: (16) in, These are Christopher notations, used to describe the curvature relations on Riemannian manifolds. Their calculation formula is: (17) In metric tensors Below, we obtain: (18) Further simplification of Christopher's symbol: (19) In a real three-dimensional heterogeneous velocity field, considering that only the first arrival of seismic waves is calculated, and that waves will continuously refract as they propagate in the non-uniform three-dimensional medium, an approximation strategy is adopted, introducing the Fast Mash Method (FMM) based on the travel time field for solution; this is achieved by constructing the minimum propagation time field from the source or station to various points in space. And its gradient norm describes the propagation characteristics of the wavefront in space, where Indicates the spatial location of the traveling field The gradient magnitude at a given point satisfies the following Eikonal equation: (20) From a certain station Starting from here, obtain the time field of the entire area. For the minimum propagation time at each point, the geodesic distance is approximately equal to the minimum cumulative propagation time between the two stations: (21) Using the steepest descent path from Tracing back to To reconstruct the propagation path, the edge weights of the adjacency graph are constructed using geodesic distances, and defined in the form of a Gaussian kernel function: (22) in, The typical propagation distance between stations is used to control weight decay; to further reflect the nonlinear effect of velocity changes on the propagation path, a velocity gradient term is introduced to enhance the sensitivity of the graph structure to geological heterogeneity. (23) in, (24) It is a hyperparameter that adjusts the strength of velocity change regularization, representing a weighted penalty for regions with significant velocity changes along the path, used to enhance the model's response to highly heterogeneous regions; Indicates the integral coefficient along the path; This represents the velocity gradient, i.e., the change of velocity with space. Let be the length of the infinitesimal line segment.

5. The coal mine microseismic location method based on Riemannian manifold graph convolutional network and physical constraints according to claim 3, characterized in that, In step 2.2), the specific method is as follows: Given a set of propagation-related features extracted from GCN ,in It is the number of stations. It is the embedded dimension; First, aggregate these node features into a unified three-dimensional spatial representation; then design a spatial projection mechanism. ,in It is a discretized grid space. It is the number of feature channels; Each grid voxel corresponds to a candidate position. Its input features are constructed by interpolating or aggregating GCN feature information from nearby stations, and the weights are given by inverse geodesic distance: (25) in, It is a station and voxel center Geodetic distance between them By controlling the spatial influence scale, graph-based features can be smoothly embedded into the volumetric space. The weighting coefficients representing the characteristics of voxels and station nodes are calculated and normalized using geodesic distance; Indicates the first Feature vectors of each node; Subsequently, the fully convolutional network FCN obtains 3D tensors through convolutional layers: (26) in, Represents 3D convolution. and It is the first Learnable convolutional kernels and biases of layers, It is a non-linear activation function that ultimately outputs a feature tensor. The shape is M×M×M, which represents the predicted three-dimensional probability distribution in the discretized source space; Apply voxel-level softmax normalization to each voxel: (27) Observed propagation time difference characteristics Given the condition, generate a value at all candidate positions. probability volume on .

6. The coal mine microseismic location method based on Riemannian manifold graph convolutional network and physical constraints according to claim 3, characterized in that, In step 2.3), the specific method is as follows: A multi-component composite loss function was designed: (28) in, and These are weighting coefficients; Representing the KL divergence term, its objective is to make the predicted spatial distribution of earthquake sources... With real labels Consistency: (29) in, Represents the three-dimensional coordinates of the predicted earthquake source location; Three-dimensional coordinates representing the actual location of the earthquake source; Representing the mean squared error term, the point where the maximum value of the predicted probability volume is taken as the estimated source location, i.e.: (30) Calculate and estimate the location of the earthquake source Compared with the actual epicenter location The Euclidean distance between them is calculated, and the mean squared error loss term is constructed: (31) This allows for optimization of the loss term, directly constraining the positioning error and improving the accuracy of seismic source positioning. By introducing a forward modeling loss and using travel time residuals for inversion, the physical consistency of the model with the source location is enhanced. The theoretical travel time is obtained by forward modeling based on the velocity model by predicting the source location. It is a physical consistency term, characterized by its time difference with the input. By comparison, we obtain: (32) in, This represents the station calculated from the predicted location. and The time difference in propagation between them; Indicates the observed stations and The time difference in propagation between them; By reversing the predicted location to obtain the propagation path and matching it with the original input, reverse consistency and physical plausibility are enhanced.

7. The coal mine microseismic location method based on Riemannian manifold graph convolutional network and physical constraints according to claim 1, characterized in that, In step 3), the specific method is as follows: During the training phase, forward modeling techniques were used to calculate the arrival time from each candidate source to all seismic stations. The Fast Mover Method (FMM) was used to solve the Eikonal equation to calculate the travel time difference between each station and all other stations, and feature vectors were constructed. All forward modeling results were pre-calculated and cached before training. At the same time, the Gaussian volume representation corresponding to the location of each candidate source was also pre-discretized. The Adam optimizer was used for training. The initial learning rate was set to 10⁻³, and an exponential decay strategy was adopted, with the learning rate decaying once every 10 epochs. The decay factor was 0.9, the batch size was set to 32, and a total of 200 epochs of training were performed. The overall training process includes the following steps: sample selection, feature construction through time difference, graph convolution processing, spatial mapping to generate 3D output, calculation of loss function, and subsequent optimization through backpropagation and parameter update.

Images