Tunnel seismic travel time tomography method and system

Abstract

The invention provides a tunnel seismic travel time tomography method and system, and the method comprises the steps: constructing a target geologic model, and forming a tight connection diagram; executing a shortest path algorithm on the close connection diagram, determining a shortest path between any two points, and determining a propagation path of seismic waves in combination with the Fermat theorem; based on the propagation path, constructing a calculation graph of seismic wave propagation, obtaining travel time, estimating travel time deviation through a loss function, obtaining gradients of all variables on the calculation graph through deep learning back propagation, iteratively updating speed distribution of the model and possibly existing geological interface positions; and optimizing the iterative updating process by using a stochastic gradient descent optimization method, and in the optimization process, applying a linear smooth loss function constraint to nodes of the fault zone at the same place, and smoothing the wave velocity distribution in the fault under the condition of ensuring that the model conforms to a seismic wave observation result. According to the invention, the calculation speed and prediction accuracy of seismic travel time imaging are improved.

Description

Technical Field

[0001] The invention belongs to the field of geophysical exploration, and in particular relates to a tunnel seismic wave travel time tomography method and system. Background Art

[0002] The statements in this section merely provide background information related to the present invention and do not necessarily constitute prior art.

[0003] A tunnel is an underground engineering facility and a way for humans to utilize underground space. It is of great benefit in shortening travel time, improving transportation efficiency, alleviating traffic pressure, and increasing the convenience of travel for residents. It has great significance in the fields of economy, people's livelihood, regional integration and development.

[0004] During tunnel construction, we often face complex geological conditions, such as fault zones, caves, broken zones, water-rich zones and other adverse geological phenomena. In order to ensure the safety of the project and avoid damage to personnel and property, it is necessary to predict the adverse geology ahead of the tunnel.

[0005] At present, the commonly used methods for tunnel advance geological prediction include high-density resistivity method, geological radar, seismic wave method, etc. The detection technology based on seismic waves has the advantages of high interface sensitivity and long detection distance, and is widely popular in tunnel advance geological prediction. However, the currently commonly used seismic wave travel time imaging method often has problems such as high model calculation complexity and unstable calculation results, which are subject to certain limitations in practical applications. Summary of the invention

[0006] In order to solve the above problems, the present invention proposes a tunnel seismic wave traveltime tomography method and system. The present invention refers to the deep learning optimization paradigm, first constructs a forward calculation graph for seismic wave traveltime calculation, then combines the efficient deep learning back propagation method to calculate the gradient, and finally uses the efficient optimizer of deep learning to update the model parameters, thereby avoiding the instability of traditional methods for solving large sparse linear equations, improving the calculation speed and prediction accuracy of seismic wave traveltime imaging, and enhancing the versatility of seismic wave traveltime imaging technology.

[0007] According to some embodiments, the present invention adopts the following technical solutions:

[0008] A tunnel seismic wave travel time tomography method comprises the following steps:

[0009] Construct the target geological model to form a tightly connected map;

[0010] Executing the shortest path algorithm on the densely connected graph to determine the shortest path between any two points, and combining Fermat's theorem to determine the propagation path of the seismic wave;

[0011] Based on the propagation path, a calculation graph of seismic wave propagation is constructed to obtain travel time, and the travel time deviation is estimated through a loss function. The gradients of all variables on the calculation graph are obtained through deep learning back propagation, and the velocity distribution of the model and the possible geological interface positions are iteratively updated;

[0012] The stochastic gradient descent optimization method is used to optimize the above iterative update process. During the optimization process, a linear smoothing loss function constraint is imposed on the nodes in the same fault zone to smooth the wave velocity distribution within the fault while ensuring that the model is consistent with the seismic wave observation results.

[0013] As an optional implementation, the process of constructing a target geological model and forming a tightly connected graph includes: for an unknown geological model, using a finite element analysis method, dividing the continuous bottom layer into a plurality of regular cells, each cell being a node of the tightly connected graph, and recording the seismic wave velocity of the corresponding geological block;

[0014] Construct the edges connected to each node, and for each edge, record the length and weight of the edge;

[0015] A densely connected graph is constructed using nodes and edges.

[0016] As an optional implementation, the shortest path algorithm is executed on the densely connected graph, and the process of determining the shortest path between any two points includes:

[0017] Execute the shortest path algorithm on the densely connected graph to obtain the edge set of the shortest path between any two points and the corresponding edge length set and edge weight set;

[0018] The shortest path edge set is further divided into two segments according to the reflection surface: from the transmitter to the reflection point on the interface; from the reflection point on the interface to the receiver;

[0019] In addition, a correction item is added to the edge set of the shortest path to compensate for the error between the current interface and the real interface, and the correction item will share the slowness parameter of the node at the current position.

[0020] As an optional implementation, in combination with Fermat's theorem, the process of determining the propagation path of seismic waves includes determining the shortest path between each pair of transmitters and receivers based on the shortest path edge set, and performing forward propagation of seismic wave travel time.

[0021] As an optional implementation, the loss function is a mean square error loss function.

[0022] As an optional implementation, the process of obtaining the gradients of all variables on the calculation graph through deep learning back propagation and iteratively updating the wave velocity model parameters and the possible geological interface positions includes: recording the calculation graph during forward calculation, and calculating the gradient of each parameter through back propagation, and back propagation applies the chain rule followed in derivation to the calculation graph to obtain the derivative of the combination of differentiable functions of velocity distribution and interface position;

[0023] The gradients generated from the loss function propagate along the computational graph to the velocity distribution and interface positions of the model.

[0024] As a further limited implementation, a gradient is set for the interface update, and the gradient limits the interface to update only one unit length each time.

[0025] As an optional implementation, the process of applying a linear smooth loss function constraint to nodes in the same fault zone includes:

[0026] For each model node, its neighborhood is intercepted by calculating the distance to other nodes, and this area does not cross the fracture surface of the current model;

[0027] Define the domain U(p i ) for each node p j The relative slowness difference is: where d ij For p i to p j The Euclidean distance, s i 、s j is the corresponding current slowness;

[0028] The constraint of the entire neighborhood is represented by using the average of the maximum speed difference and the minimum speed difference.

[0029] As a further limited implementation method, it also includes adding a borehole prior constraint in the optimization process. If there is a borehole that samples part of the tunnel front address, the wave velocity at the corresponding position is estimated by density measurement and other means. Then we get the estimated slowness A priori constraints are placed on the estimated slowness and the model slowness.

[0030] A tunnel seismic wave travel time tomography system, comprising:

[0031] a model building module configured to build a target geological model to form a tightly connected graph;

[0032] a propagation path determination module, configured to execute a shortest path algorithm on the densely connected graph to determine the shortest path between any two points, and determine the propagation path of the seismic wave in combination with Fermat's theorem;

[0033] A back propagation optimization module is configured to construct a computational graph of seismic wave propagation based on the propagation path, obtain travel time, estimate travel time deviation through a loss function, obtain the gradient of all variables on the computational graph through deep learning back propagation, and iteratively update the velocity distribution of the model and the possible geological interface position;

[0034] The prior information optimization module is configured to optimize the above-mentioned iterative update process using a stochastic gradient descent optimization method, and during the optimization process, a linear smoothing loss function constraint is imposed on the nodes of the same fault zone to smooth the wave velocity distribution within the fault while ensuring that the model conforms to the seismic wave observation results.

[0035] Compared with the prior art, the present invention has the following beneficial effects:

[0036] The present invention utilizes the shortest path search on the compact graph, which greatly saves the calculation time of the seismic wave path and improves the execution speed of the method.

[0037] The present invention calculates gradients and updates wave velocity model parameters through deep learning back propagation and optimizer, avoids solving large sparse linear equations, and improves solution stability and accuracy.

[0038] The present invention achieves synchronous updating of wave velocity and interface by adding geological interface position prediction to seismic wave travel time prediction, improves calculation efficiency and enhances the versatility of the method.

[0039] The present invention specifies known geological information through preset slowness uniformity constraints, drilling prior constraints and other methods, and participates in calculations during the optimization process, so that the final result is more stable and practical.

[0040] The present invention provides a practical tunnel advanced geological prediction imaging technology, which can be widely used in the fields of tunnel exploration, construction engineering, etc.

[0041] In order to make the above-mentioned objects, features and advantages of the present invention more obvious and easy to understand, preferred embodiments are given below and described in detail with reference to the accompanying drawings. BRIEF DESCRIPTION OF THE DRAWINGS

[0042] The accompanying drawings in the specification, which constitute a part of the present invention, are used to provide a further understanding of the present invention. The exemplary embodiments of the present invention and their descriptions are used to explain the present invention and do not constitute improper limitations on the present invention.

[0043] Figure 1 The overall process of the method provided in this embodiment;

[0044] Figure 2 A schematic diagram of a geological finite element diagram designed for the method provided in this embodiment;

[0045] Figure 3 A schematic diagram of the error between the current interface and the actual interface during the iteration of the method provided in this embodiment;

[0046] Figure 4 A calculation diagram in the method provided in this embodiment;

[0047] Figure 5 This is a schematic diagram of the general back propagation;

[0048] Figure 6 A comparison chart of the method results provided in this embodiment and the results of solving sparse linear equations using traditional methods;

[0049] Figure 7 The stability of the method provided in this embodiment under poor hyperparameters;

[0050] Figure 8 The drilling prior Figure 7 Correction of undesirable results. DETAILED DESCRIPTION

[0051] The present invention will be further described below in conjunction with the accompanying drawings and embodiments.

[0052] It should be noted that the following detailed descriptions are all illustrative and intended to provide further explanation of the present invention. Unless otherwise specified, all technical and scientific terms used herein have the same meanings as those commonly understood by those skilled in the art to which the present invention belongs.

[0053] It should be noted that the terms used herein are only for describing specific embodiments and are not intended to limit exemplary embodiments according to the present invention. As used herein, unless the context clearly indicates otherwise, the singular form is also intended to include the plural form. In addition, it should be understood that when the terms "comprising" and / or "including" are used in this specification, it indicates the presence of features, steps, operations, devices, components and / or combinations thereof.

[0054] In the absence of conflict, the embodiments in this application and the features in the embodiments may be combined with each other.

[0055] Embodiment 1:

[0056] Embodiment 1 of the present invention provides a tunnel seismic wave travel time tomography method based on a deep learning optimization paradigm. On the one hand, this method refers to the deep learning optimization paradigm to construct a forward calculation and reverse gradient update mode of seismic wave travel time imaging, realize the iterative update of the velocity model, and finally obtain good tunnel advanced geological prediction results.

[0057] First, a densely connected graph data structure is used to store geological model parameters to facilitate subsequent calculations.

[0058] Furthermore, the shortest path on the densely connected graph is calculated and combined with Fermat's theorem to determine the propagation path of the seismic wave.

[0059] The shortest path mentioned above refers to the sum of three paths from the transmitter to the current reflection surface, from the receiver to the current reflection surface, and from the current reflection surface to the real reflection surface.

[0060] Construct a seismic wave propagation calculation graph to obtain the travel time, and estimate the travel time deviation through the loss function. Use deep learning back propagation to obtain the gradients of all variables on the calculation graph, and then update the wave velocity model parameters and the possible geological interface positions.

[0061] The SGD optimization method is used to optimize the iterative update process, avoid falling into the local minimum and improve the quality of the results.

[0062] On the other hand, through the deep learning optimization paradigm, this method can add a variety of prior information to the optimization process, which improves the operator's control over the input information and improves the quality of imaging results.

[0063] By analyzing the geological sections, we can obtain the fault zones of different strata and determine the fault zones to which the nodes in each figure belong.

[0064] A linear smoothing loss function constraint is imposed on the nodes in the same fault zone, and the above-mentioned iterative optimization of seismic wave travel time is performed in parallel to smooth the wave velocity distribution within the fault while ensuring that the model is consistent with the seismic wave observation results.

[0065] For some embodiments, for possible borehole information, stricter borehole wave velocity constraints are added while optimizing the seismic wave travel time and smoothing the wave velocity within the fault, so as to ensure that the results are closer to the actual engineering situation while ensuring that the results meet the aforementioned conditions.

[0066] This embodiment introduces the specific details of each step of the above process in detail:

[0067] (1) Construction of a graph with tightly connected geological information and calculation of the shortest path

[0068] For unknown geological models, the finite element analysis method is used to divide the continuous bottom layer into a large number of discrete regular cells V. Each cell serves as a node in a tightly connected graph and records the seismic wave velocity v of its corresponding geological block.

[0069] For each cell node v t , an edge e connected to it will be constructed: Figure 2 As shown, if the node v i There is another node in its immediate neighborhood Edge ij will be constructed connecting these two nodes.

[0070] For each edge e ij , record the length of the side l ij (Euclidean distance between two nodes), edge weight s ij (The slowness of the seismic wave propagation at the edge, that is, the inverse of the seismic wave velocity ). The set of all edges is denoted by E.

[0071] The graph G = (V, E) constitutes a complete tightly connected graph describing geological information.

[0072] By executing the shortest path algorithm (such as Dijkstra algorithm or A* algorithm, etc.) on G, we can get the edge set of the shortest path between any two points. And the corresponding edge length set Edge weight set Edge Set represents the path before the seismic wave transmitter-receiver pair (m, n), and (a, b, ..., z) are the intermediate nodes on the path.

[0073] For the above path It can be further divided into two sections according to the reflection surface (fracture surface) r: from the transmitter m to the reflection point on the interface From the reflection point on the interface To receiver n. The total path is divided into:

[0074] In the iterative process, Figure 3 As shown, the current interface t There may be some error between the actual interface r and the actual interface r, so in the iterative calculation process, it is necessary to Add a Modification To compensate for this error:

[0075] To simplify the problem, Will share the slowness parameters of the node at the current position

[0076] The above corrections help to add the impact of interface movement in the shortest path calculation, preparing for subsequent deep learning optimization methods.

[0077] (2) Non-fitting calculation of seismic wave travel time

[0078] Seismic wave travel time tomography is an iterative optimization process, including multiple forward modeling and inversion steps. It is worth noting that in deep learning, the forward modeling step is consistent with the forward propagation process of the model, while the inversion step is consistent with the backward propagation process. Therefore, the deep learning method can efficiently and quickly iterate and optimize the model.

[0079] This project sets two learnable parameters in the forward modeling of seismic wave travel time tomography: velocity distribution v t and interface position r t .

[0080] (a) Calculation diagram of seismic wave travel time

[0081] By using the shortest path between each transmitter-receiver pair calculated in (1), we can easily perform the forward propagation of seismic wave travel time:

[0082]

[0083] In the formula That is the travel time of the seismic wave to be calculated. and Represents the shortest path edge set of the transceiver pair (m, n) The length and slowness of the i-th edge in . The summation process can be achieved by the shortest path segmentation method in (1): Split into three sub-sections.

[0084] By considering the transceiver-reflector propagation process as a whole, the above formula can be simplified to:

[0085]

[0086] In the formula That is, from the transmitter m to the reflection point The average slowness, From reflection point n to receiver The average slowness, It is the slowness between the current reflection point and the actual reflection point; From the transmitter m to the reflection point The total side length of From reflection point n to receiver The total side length of It is the total side length from the current reflection point to the actual reflection point.

[0087] calculate Figure 4 The forward calculation process and the back propagation process for the formula are demonstrated.

[0088] (b) Loss function

[0089] The loss function is a formula used to evaluate the deviation between the predicted value and the true value of the model, and can generate the basic gradient for back propagation. The accuracy of the model can be improved by iteratively minimizing the loss function. The present invention uses the mean square error (MSE) loss function as an evaluation index and training model. In the case of N transceiver pairs, the MSE loss function is expressed as follows:

[0090]

[0091] In the formula is the calculated seismic wave travel time of the nth transceiver pair, d n is the actual seismic wave travel time of the nth transceiver pair.

[0092] (3) Model parameter optimization

[0093] In order to optimize the model parameters velocity distribution v t and interface position r t , we need to calculate the updated gradients of both and Traditional methods obtain the gradient by solving a large sparse linear equation system, but the results are often unstable, computationally time-consuming, or stuck in local optimal solutions.

[0094] By using the automatic derivation and back-propagation mechanism of deep learning, we can avoid the calculation of sparse linear equations and directly update v t and r t Modern deep learning optimizers can also better alleviate the problem of local optimal solutions. The specific process is as follows:

[0095] (a) Obtaining parameter gradients through back propagation

[0096] The present invention uses modern deep learning tools (such as pytorch) to automatically record the calculation graph during forward calculation and calculate the gradient of each parameter through back propagation. Its basic concepts are as follows Figure 5 The back propagation principle is to apply the chain rule followed in derivation to the computational graph to find the derivative of the combination of two differentiable functions.

[0097] The back propagation process of the present invention has been Figure 3 The gradient generated from the loss function The velocity distribution and interface position of the model will be propagated along the computational graph: This completes the gradient calculation of the parameters.

[0098] In particular, in order to make the interface update more stable, the present invention resets the gradient of the interface update: Limiting the interface update to no more than one unit length each time makes the update process more continuous and stable.

[0099] The above method avoids solving large sparse linear equations.

[0100] (b) Update the model through the SGD optimizer

[0101] According to the parameter gradient obtained above, the model is updated as follows:

[0102] Node slowness:

[0103] Node wave speed:

[0104] Reflection point location:

[0105] Note that the reflection points only retain floating point values ​​when updating, and they need to be rounded to align with the grid map during the actual calculation process.

[0106] Compared with the traditional method that simultaneously calculates all propagation paths of all transmitter-receiver pairs and updates parameters at the same time, the present invention adopts the stochastic gradient descent (SGD) optimization method, which randomly selects a set of shortest paths for updating each time, thereby better alleviating the local minimum problem and obtaining better results.

[0107] Embodiment 2:

[0108] This embodiment is the result of smoothing embodiment 1, and adds more practical priors, specifically including:

[0109] (1) Smoothness loss constraint: constrain the predicted value through local linear estimation

[0110] Specifically, for each model node, its neighborhood U(p) is intercepted by calculating the distance to other nodes. Note that this area should not cross the fracture surface of the current model. The area size should cover a large range and should not be too small.

[0111] Define the domain U(p i ) for each node p j The relative slowness difference is: where d ij For p i to p j The Euclidean distance, s i 、s j is the corresponding current slowness.

[0112] The present invention proposes that the relative slowness difference can fully reflect the linear interpolation expectation of each point in the field to the center point (assuming that the infinite distance ). In actual constraints, the constraint effect of the entire neighborhood is approximately represented by the average of the maximum speed difference and the minimum speed difference:

[0113]

[0114] The program needs to calculate right The loss function (such as mean square error MSE) makes the calculation process differentiable, where Not included in the computation graph.

[0115] (2) (Optional) Drilling Priors: Constrain the predicted values ​​by presupposing the expected values ​​of nodes within a certain range

[0116] Specifically, if there is a borehole that samples part of the tunnel front address, the wave velocity at the corresponding location can be estimated by density measurement and other means. Then we get the estimated slowness

[0117] A priori constraints are placed on the estimated slowness and the model slowness:

[0118]

[0119] The invention adopts the hyperbolic tangent function, which can make the slowness model close to the estimated value still maintain a high update intensity, and at the same time can reduce the influence of the wrong estimation to a certain extent. The constant k in the formula is used to adjust the confidence of the estimated value.

[0120] Embodiment 3:

[0121] The method was tested by simulation model, and the results are as follows Figure 6 shown.

[0122] The input of this embodiment is an initial model for determining the number of fault zones, the initial seismic wave velocity is set to a uniform tunnel side detection result, and the model is iteratively optimized through the smoothness loss constraints of embodiment 1 and embodiment 2. The traditional linear equations model uses similar initial conditions.

[0123] The method provided by the present invention has a high tolerance for bad hyperparameters (learning rate). Figure 7 It shows that this method still has good stability under different hyper-parameters. Figure 8 The drilling priori constraint in application example 2 is Figure 7 The corrective effect of the adverse consequences produced.

[0124] Embodiment 4:

[0125] A tunnel seismic wave travel time tomography system, comprising:

[0126] a model building module configured to build a target geological model to form a tightly connected graph;

[0127] a propagation path determination module, configured to execute a shortest path algorithm on the densely connected graph to determine the shortest path between any two points, and determine the propagation path of the seismic wave in combination with Fermat's theorem;

[0128] A back propagation optimization module is configured to construct a computational graph of seismic wave propagation based on the propagation path, obtain travel time, estimate travel time deviation through a loss function, obtain the gradient of all variables on the computational graph through deep learning back propagation, and iteratively update the velocity distribution of the model and the possible geological interface position;

[0129] The prior information optimization module is configured to optimize the above-mentioned iterative update process using a stochastic gradient descent optimization method, and during the optimization process, a linear smoothing loss function constraint is imposed on the nodes of the same fault zone to smooth the wave velocity distribution within the fault while ensuring that the model conforms to the seismic wave observation results.

[0130] Those skilled in the art will appreciate that embodiments of the present invention may be provided as methods, systems, or computer program products. Therefore, the present invention may take the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware. Moreover, the present invention may take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0131] The present invention is described with reference to flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to embodiments of the present invention. It should be understood that each process and / or block in the flowchart and / or block diagram, as well as the combination of processes and / or blocks in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to produce a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 A process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in a block or multiple blocks.

[0132] These computer program instructions may also be stored in a computer-readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer-readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 A process or multiple processes and / or boxes Figure 1 A function specified in one or more boxes.

[0133] These computer program instructions can also be loaded onto a computer or other programmable data processing device so that a series of operating steps are executed on the computer or other programmable device to produce a computer-implemented process, thereby providing instructions for implementing the process in the computer or other programmable device. Figure 1 A process or multiple processes and / or boxes Figure 1 The steps for the functions specified in one or more boxes.

[0134] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. For those skilled in the art, the present invention may have various modifications and variations. Any modification, equivalent replacement, improvement, etc. made by those skilled in the art within the spirit and principle of the present invention without creative labor shall be included in the protection scope of the present invention.

Claims

1. A tunnel seismic wave travel time tomography method, characterized in that: The following steps are involved: Construct the target geological model to form a tightly connected map; Executing the shortest path algorithm on the densely connected graph to determine the shortest path between any two points, and combining Fermat's theorem to determine the propagation path of the seismic wave; Based on the propagation path, a calculation graph of seismic wave propagation is constructed to obtain travel time, and the travel time deviation is estimated through a loss function. The gradients of all variables on the calculation graph are obtained through deep learning back propagation, and the velocity distribution of the model and the possible geological interface positions are iteratively updated; The stochastic gradient descent optimization method is used to optimize the above iterative update process. During the optimization process, a linear smoothing loss function constraint is imposed on the nodes in the same fault zone to smooth the wave velocity distribution within the fault while ensuring that the model is consistent with the seismic wave observation results.

2. A tunnel seismic wave travel time tomography method as claimed in claim 1, characterized in that: The process of constructing the target geological model and forming a tightly connected graph includes: for an unknown geological model, using the finite element analysis method, dividing the continuous bottom layer into a plurality of regular cells, each cell being a node of the tightly connected graph, and recording the seismic wave velocity of the corresponding geological block; Construct the edges connected to each node, and for each edge, record the length and weight of the edge; A densely connected graph is constructed using nodes and edges.

3. A tunnel seismic wave travel time tomography method as claimed in claim 1, characterized in that The process of executing the shortest path algorithm on the densely connected graph to determine the shortest path between any two points includes: Execute the shortest path algorithm on the densely connected graph to obtain the edge set of the shortest path between any two points and the corresponding edge length set and edge weight set; The shortest path edge set is further divided into two segments according to the reflection surface: from the transmitter to the reflection point on the interface; from the reflection point on the interface to the receiver; In addition, a correction item is added to the edge set of the shortest path to compensate for the error between the current interface and the real interface, and the correction item will share the slowness parameter of the node at the current position.

4. A tunnel seismic wave travel time tomography method as claimed in claim 1, characterized in that: Combined with Fermat's theorem, the process of determining the propagation path of seismic waves includes determining the shortest path between each pair of transmitters and receivers based on the set of shortest path edges, and performing forward propagation of seismic wave travel time.

5. The tunnel seismic wave travel time tomography method according to claim 1, characterized in that: The loss function is a mean square error loss function.

6. A tunnel seismic wave travel time tomography method as claimed in claim 1, characterized in that: The process of obtaining the gradients of all variables on the calculation graph through deep learning back propagation and iteratively updating the wave velocity model parameters and the possible geological interface positions includes: recording the calculation graph during forward calculation and calculating the gradient of each parameter through back propagation. Back propagation applies the chain rule followed in derivation to the calculation graph to obtain the derivative of the combination of differentiable functions of velocity distribution and interface position; The gradients generated from the loss function propagate along the computational graph to the velocity distribution and interface positions of the model.

7. A tunnel seismic wave travel time tomography method as claimed in claim 6, characterized in that: Set a gradient for interface updates, which limits the interface to update only one unit length each time.

8. The tunnel seismic wave travel time tomography method according to claim 1, characterized in that: The process of applying linear smooth loss function constraints to nodes in the same fault zone includes: For each model node, its neighborhood is intercepted by calculating the distance to other nodes, and this area does not cross the fracture surface of the current model; Define the domain U(p i ) for each node p j The relative slowness difference is: where d ij For p i to p j The Euclidean distance, s i 、s j is the corresponding current slowness; The constraint of the entire neighborhood is represented by using the average of the maximum speed difference and the minimum speed difference.

9. A tunnel seismic wave travel time tomography method as claimed in claim 8, characterized in that: It also includes adding borehole prior constraints during the optimization process. If there is a borehole that samples part of the tunnel front address, the wave velocity at the corresponding position is estimated by density measurement and other means. Then we get the estimated slowness A priori constraints are placed on the estimated slowness and the model slowness.

10. A tunnel seismic wave travel time tomography system, characterized in that: include: a model building module configured to build a target geological model to form a tightly connected graph; a propagation path determination module, configured to execute a shortest path algorithm on the densely connected graph to determine the shortest path between any two points, and determine the propagation path of the seismic wave in combination with Fermat's theorem; A back propagation optimization module is configured to construct a computational graph of seismic wave propagation based on the propagation path, obtain travel time, estimate travel time deviation through a loss function, obtain the gradient of all variables on the computational graph through deep learning back propagation, and iteratively update the velocity distribution of the model and the possible geological interface position; The prior information optimization module is configured to optimize the above-mentioned iterative update process using a stochastic gradient descent optimization method, and during the optimization process, a linear smoothing loss function constraint is imposed on the nodes of the same fault zone to smooth the wave velocity distribution within the fault while ensuring that the model conforms to the seismic wave observation results.

Images