Tunnel earthquake advanced detection waveform inversion method and system

Abstract

The invention belongs to the technical field of tunnel advanced detection, and particularly relates to a tunnel earthquake advanced detection waveform inversion method and system, and the method comprises the steps: carrying out the forward simulation according to an initial speed model and an earthquake source, and obtaining a simulated seismic record of a simulated seismic wave; establishing an objective function by minimizing the difference between the observed seismic data and the simulated seismic record; and performing inversion optimization on the target function by adopting an enhanced whale optimization algorithm, and updating the initial velocity model through repeated iteration until the inversion reaches the minimum value of the target function, so that the problem that the inversion precision is reduced along with the increase of the problem dimension in the tunnel space waveform inversion is solved, and the tunnel space waveform inversion accuracy is improved. According to the method, the global search and local development capabilities are improved, and the inversion precision and the convergence speed are improved.

Description

Technical Field

[0001] This application belongs to the technical field of tunnel advanced detection, and specifically relates to a waveform inversion method and system for tunnel seismic advanced detection. Background Art

[0002] During the tunnel excavation process, complex geological conditions such as faults, karst caves, and soft surrounding rocks are often encountered, and these factors are likely to trigger geological disasters. Therefore, before tunnel construction, accurately grasping geological information is crucial for reducing risks. Full waveform inversion (FWI), as a high-resolution seismic imaging method, reconstructs the underground velocity structure by precisely matching the measured waveform and the simulated waveform. However, its application in tunnel detection has the defect of being highly sensitive to the initial model and easily falling into local minima, which limits its inversion accuracy in tunnel advanced detection.

[0003] The unscented hybrid simulated annealing algorithm (UHSA) has been successfully applied to the inversion research of tunnel velocity models. This algorithm combines the global optimization ability of the simulated annealing algorithm (SA) and the local optimization ability of the unscented Kalman filter (UKF), significantly improving the inversion accuracy. However, as the problem dimension increases, UHSA is prone to falling into local optimal solutions, especially when the number of velocity bodies increases, the inversion accuracy of the algorithm significantly decreases. Summary of the Invention

[0004] The technical problem to be solved by this application is to provide a waveform inversion method for tunnel seismic advanced detection, which solves the problem that the inversion accuracy decreases as the problem dimension increases in the waveform inversion of the tunnel space.

[0005] This application is implemented as follows: A waveform inversion method for tunnel seismic advanced detection includes: performing forward simulation based on an initial velocity model and a seismic source to obtain a simulated seismic record of the simulated seismic wave; establishing an objective function with the minimum difference between the observed seismic data and the simulated seismic record; using an enhanced whale optimization algorithm to perform inversion optimization on the objective function, and repeatedly iterating to update the initial velocity model until the inversion reaches the minimum value of the objective function. The using of the enhanced whale optimization algorithm to perform inversion optimization on the objective function includes: Initializing the parameters of the enhanced whale optimization algorithm, where the parameters include the initial value of the iteration number, the maximum iteration number, and the frequency factor; Generating an initial population, generating the initial population using a piecewise Chebyshev-logistic map, and adding the reverse solution of the current individual to the initial population; Adjusting the non-linear convergence factor; Updating the first coefficient vector and the second system vector and the first random number and the second random number ; Trigger an operation including a shrinking enclosure mechanism, a random search mechanism, a bubble net attack mechanism, and a frequency fluctuation mechanism according to the values of the first coefficient vector and the second random number ; Determine whether the maximum number of iterations or a preset termination condition is reached. If the termination condition is satisfied, output the optimal individual position and the fitness value of the optimal individual.

[0006] Further, for generating the initial population, a piecewise Chebyshev-logistic map is used to generate the initial population, and the reverse solution of the current individual is added to the initial population, including:[[]] Generating the reverse solution of the current individual to increase population diversity is expressed as:[[]] where represents the reverse solution of the current individual,[[]] represents the lower limit of the individual position variable,[[]] represents the upper limit of the individual position variable,[[]] represents the th th position variable of the th individual, and the calculation method is:[[]] represents the th chaotic variable,[[]] The update of is:[[]] with the perturbation factor representing the population size,[[]] is a parameter for controlling the behavior of the logistic map.[[]]

[0007] Further, adjusting the non-linear convergence factor includes:[[]] represents the non-linear convergence factor,[[]] represents the maximum number of iterations,[[]] represents the number of iterations.[[]]

[0008] Further, the shrinking enclosure mechanism includes: when is satisfied and the global optimal individual represents the prey space coordinates, and other individuals gradually approach the global optimal individual by adjusting their own positions. The population update formula is:[[]] where represents the weighted distance between the individual and the prey,[[]] represents the position of the individual at time,[[]] is the position of the global optimal individual of the current population,[[]] ​​​​Indicates the position of an individual at moment; The first coefficient vector and the second coefficient vector are calculated by the following formula: , , where, is the third random number between [0, 1], is the fourth random number between [0, 1], is the non-linear convergence factor.

[0009] Furthermore, the random search mechanism includes: when satisfying , in the iteration process, each individual randomly selects any feasible solution in the current population as the target position with a probability of and expands the search range by moving towards the target position. The population update formula is: , , where, is the target position randomly selected in the current population, represents the weighted distance between the individual and the prey, represents the position of the individual at moment, represents the position of the individual at moment.

[0010] Furthermore, the bubble net attack mechanism includes: when satisfying , the population update formula is: , , where, is set to 1, represents the Euclidean distance between the individual and the prey, represents the position of the individual at moment, is the position of the global optimal individual in the current population, represents the position of the individual at moment, is the modulo operation on the iteration number and the frequency factor .

[0011] Furthermore, the frequency fluctuation mechanism includes: when satisfying , At this time, the update of the population: , is the amplitude factor, defined as: , where is the angular frequency, calculated as ; is a random phase angle within the interval, is the oscillation factor, such that the amplitude factor varies periodically with time. By adjusting the angular frequency and the phase angle, the frequency and phase of the oscillation are controlled to guide the search agent to explore in the search space. The amplitude of the oscillation gradually decreases as the number of iterations increases. is the modulo operation on the number of iterations and the frequency factor , represents the maximum number of iterations.

[0012] Furthermore, it also includes using the low-altitude hovering mechanism of the red-tailed hawk to enhance its local development ability. The specific formula is expressed as: , , , , After normalization, replace and : , where, and represent the direction coordinates at time is expressed as the adaptive step size factor, is the average value of the positions, represents the initial value of the radius, represents the angle gain, is the random gain, is the control gain, represents the instantaneous radius value during the search, represents the instantaneous angle during the search, is the position of the global optimal individual of the current population, represents the position of the individual at time.

[0013] On the other hand, an embodiment of the present application provides a waveform inversion system for tunnel seismic forward detection, including: A forward modeling module for performing forward simulation based on an initial velocity model and a seismic source to obtain a simulated seismic record of the simulated seismic wave; An objective function establishment module for establishing an objective function with the minimum difference between the observed seismic data and the simulated seismic record; An inversion module for inversely optimizing the objective function by using an enhanced whale optimization algorithm, and repeatedly iterating to update the initial velocity model until the inversion reaches the minimum value of the objective function; The inversely optimizing the objective function by using the enhanced whale optimization algorithm includes: Initializing the parameters of the enhanced whale optimization algorithm, where the parameters include the initial value of the number of iterations, the maximum number of iterations, and the frequency factor; Generating an initial population, generating the initial population by using a piecewise Chebyshev-logistic map, and adding the reverse solution of the current individual to the initial population; Adjusting the non-linear convergence factor; Updating the first coefficient vector , the second system vector , the first random number and the second random number ; According to the values of the first coefficient vector and the second random number , triggering operations including a shrinking encircling mechanism, a random search mechanism, a bubble net attack mechanism, and a frequency fluctuation mechanism; Judging whether the maximum number of iterations or a preset termination condition is reached. If the termination condition is satisfied, output the optimal individual position and the fitness value of the optimal individual.

[0014] Compared with the prior art, the beneficial effects of the present application are as follows: The method of the present application improves the global search and local development capabilities. By initializing the population, the diversity and uniform distribution of the initial solutions are enhanced; a non-linear convergence factor is introduced, which can dynamically adjust the weights of exploration and development according to the search stage, better adapting to different optimization requirements; through the frequency fluctuation mechanism, the ability to jump out of the local optimum is improved; the local search accuracy is further optimized through the low-altitude hovering strategy. The inversion accuracy and convergence speed are improved. Description of the Drawings

[0015] Figure 1 is a flowchart of the method proposed in the embodiment of the present application; Figure 2 is a tunnel model parameter diagram provided in the embodiment of the present application; Figure 3The convergence curves of the whale optimization algorithm and its variants of the enhanced whale optimization algorithm in the Bent Cigar function for 500 iterations in the embodiments of the present application; Figure 4 The convergence curves of the whale optimization algorithm and its variants of the enhanced whale optimization algorithm in the Levy function for 500 iterations in the embodiments of the present application; Figure 5 The concrete residue tunnel model diagram in the embodiments of the present application; Figure 6 The concrete residue inversion result diagrams of the method (a) of the present application and the UHSA algorithm (b) in the embodiments of the present application; Figure 7 The convergence curves of the method of the present application and the UHSA algorithm for 100 iterations of the concrete residue in the embodiments of the present application; Figure 8 The transverse velocity profile diagrams at the longitudinal position of 20 m (a) and the transverse velocity profile diagram at 60 m (b) of the concrete residue of the method of the present application and the UHSA algorithm in the embodiments of the present application; Figure 9 The inclined layer change model diagram provided in the embodiments of the present application; Figure 10 The inclined layer change inversion result diagrams of the method (a) of the present application and the UHSA algorithm (b) provided in the embodiments of the present application; Figure 11 The convergence curves of the inclined layer change inversion for 100 iterations of the method of the present application and the UHSA algorithm provided in the embodiments of the present application; Figure 12 The transverse velocity profile diagrams at the longitudinal position of 32 m (a) and the transverse velocity profile diagram at 70 m (b) of the inclined layer change inversion of the method of the present application and the UHSA algorithm provided in the embodiments of the present application. Detailed implementation manners

[0016] In order to make the objectives, technical solutions and advantages of the present application clearer, the present application will be further described in detail below in conjunction with embodiments. It should be understood that the specific embodiments described herein are only used to explain the present application and are not used to limit the present application.

[0017] Refer to Figure 1 A waveform inversion method for tunnel seismic prediction ahead shown as follows, including: Performing forward simulation according to the initial velocity model and the seismic source to obtain the simulated seismic record of the simulated seismic wave; Establishing an objective function with the minimum difference between the observed seismic data and the simulated seismic record; Using the enhanced whale optimization algorithm to perform inversion optimization on the objective function, and repeatedly iterating to update the initial velocity model until the inversion reaches the minimum value of the objective function.

[0018] Among them, other inversion methods (such as travel time tomography, migration velocity analysis) can be used to establish the initial velocity model, and there is no limitation here. The process of forward simulation is to obtain the simulated seismic record of the simulated seismic wave according to the initial velocity model and the seismic source. Simulating the propagation process of the seismic wave through the given parameters of the initial velocity model and the seismic source is the forward simulation. The methods of forward simulation can include: the finite element method, the boundary element method, the pseudo-spectral method, the finite difference method, etc.

[0019] There are many types of objective functions defined for full waveform inversion, and there is no limitation here. For example, it can include: the least square form of the difference between the observed seismic data and the simulated seismic record, the cross-correlation weighted norm, the combination of cross-correlation and least squares, etc.

[0020] The optimization process of the objective function is the process of iterative updating of the initial velocity model.

[0021] In one embodiment, in order to solve the problem that the present application falls into a local optimal solution after optimization, especially when the number of velocity bodies increases, the inversion accuracy significantly decreases. The enhanced whale optimization algorithm is used to invert and optimize the objective function. It includes: Initializing the parameters of the enhanced whale optimization algorithm, the parameters include the initial value of the number of iterations, the maximum number of iterations, and the frequency factor; Generating an initial population, using a piecewise Chebyshev-Logistic mapping (Chebyshev-Logistic hybrid mapping) to generate the initial population, and adding the reverse solution of the current individual to the initial population to improve the uniformity of the population distribution by means of the piecewise Chebyshev-Logistic mapping. The chaotic variables generated by the piecewise Chebyshev-Logistic mapping. Introduce the Chebyshev mapping in the low value area to generate violent oscillations, and at the same time use the Logistic mapping to maintain the uniformity of traversal in the high value area, so as to combine the advantages of the two mappings and effectively avoid the problem of uneven coverage caused by random initialization. To further improve the coverage range of the population, a reverse learning strategy is adopted. This reverse learning strategy increases the population diversity by generating the reverse solution of the current individual.

[0022] Adjusting the non-linear convergence factor; the non-linear convergence factor optimizes the balance between exploration and development, and the behavior of dynamic adjustment enables it to widely explore the search space in the initial stage of the search, while being able to concentrate on developing potential optimal solutions in the later stage of the search. In the embodiment of the present application, in the early iteration, the convergence factor maintains a high steady state to maintain the global exploration persistence of the enhanced whale optimization algorithm; in the later iteration, the convergence factor value implements exponential decay to enhance the local development intensity to improve the convergence accuracy.

[0023] Updating the first coefficient vector The second system vector , the first random number and the second random number ; According to the first coefficient vector and the second random number The value of triggers the operations including shrinking and surrounding mechanism, random search mechanism, bubble net attack mechanism and frequency fluctuation mechanism; Among them, the first coefficient vector Used to control the moving distance of whales or individuals to the global optimal solution, the second system vector Used to control the whale's movement direction toward the global optimal solution; Determine whether the maximum number of iterations or the preset termination condition has been reached. If the termination condition is met, the optimal individual position and the fitness value of the optimal individual are output. The optimal individual position is the global optimal solution.

[0024] In one embodiment, a piecewise Chebyshev-logistic mapping is used to generate an initial population, and a reverse solution of the current individual is generated and added to the initial population, including: Generating the reverse solution of the current individual to increase population diversity is expressed as: ,in, represents the reverse solution of the current individual, represents the lower limit of the individual location variable, represents the upper limit of the individual position variable, Indicates The individual position variables, calculated as: , Indicates A chaotic variable, The update is: , disturbance factor , represents the population size, Parameters that control the behavior of the logistic mapping.

[0025] Then the fitness of the population and its anti-population is evaluated and sorted, and based on the sorting results, the top 50% of individuals with the best fitness values ​​are selected to form the initial population.

[0026] In one embodiment, adjusting the nonlinear convergence factor includes: , represents the nonlinear convergence factor, represents the maximum number of iterations, Indicates the number of iterations.

[0027] Higher nonlinear convergence factor The value makes The probability is significantly increased, making the enhanced whale optimization algorithm more inclined to global exploration. As the iteration progresses, in the middle and late stages, the value of the non-linear convergence factor rapidly decreases, enhancing the local development ability of the enhanced whale optimization algorithm. The enhanced whale optimization algorithm equipped with the non-linear convergence factor exhibits better global exploration performance in the early stage. This characteristic enables it to effectively break free from the limitation of local minima and avoid premature convergence.

[0028] In one embodiment, operations including a shrinking encircling mechanism, a random search mechanism, a bubble-net attack mechanism, and a frequency fluctuation mechanism are triggered according to the values of the first coefficient vector and the second random number .

[0029] Among them, the shrinking encircling mechanism. When , , the shrinking encircling mechanism is triggered. Under the framework of the enhanced whale optimization algorithm, the global optimal individual represents the prey space coordinates, and other individuals gradually approach the global optimal individual by adjusting their own positions. The population update formula is: , , where, represents the weighted distance between the individual and the prey, represents the position of the individual at time, is the position of the global optimal individual of the current population, represents the position of the individual at time; The first coefficient vector and the second coefficient vector are calculated by the following formula: , , where, is the third random number between [0,1], is the fourth random number between [0,1].

[0030] When , , the random search mechanism is triggered. During the iteration, each individual randomly selects any feasible solution in the current population as the target position with a probability of and expands the search range by moving towards this target position, thus enhancing the diversity and exploration ability of the algorithm. The position update equation is: , , Among them, is the target position randomly selected from the current population, represents the weighted distance between the individual and the prey, represents the individual at time, represents the individual at time.

[0031] Bubble net attack mechanism: When , is satisfied, the population update formula is: , , Among them, is set to 1, represents the Euclidean distance between the individual and the prey, represents the individual at time, is the position of the global optimal individual of the current population, represents the individual at time, is the modulo operation on the iteration number and the frequency factor . Humpback whales swim towards their prey along a spiral trajectory during predation and continuously release bubbles to form an encircling net, forcing the fish school to gather in the central area. This behavior is simulated through spiral motion.

[0032] When , is satisfied, the frequency fluctuation mechanism is triggered, and the population is updated: , is the amplitude factor, defined as: , Among them is the angular frequency, calculated as ; is a random phase angle within the interval, is the oscillation factor, making the amplitude factor vary periodically with time. By adjusting the angular frequency and the phase angle, the frequency and phase of the oscillation are controlled to guide the search agent to explore in the search space, making the amplitude of the oscillation gradually decrease as the number of iterations increases, is the modulo operation on the iteration number and the frequency factor , It represents the maximum number of iterations. To dynamically adjust the search rhythm of the population, a frequency factor is introduced , which controls the timing of vibration. A fluctuation is triggered every iterations

[0033] To enhance the local development mechanism, the low-altitude hovering mechanism of the red-tailed hawk is introduced into the whale optimization algorithm. After obtaining the range of the general global optimal individual position obtained by the whale optimization algorithm, the low-altitude hovering mechanism of the red-tailed hawk is used to approach the optimal solution to enhance the local development ability. The whale optimization algorithm combined with the low-altitude hovering mechanism of the red-tailed hawk is called the enhanced whale optimization algorithm. The specific formula is as follows: The specific formula is expressed as: , , , , After normalization, replace and : , Among them, and represent the direction coordinates at time is expressed as the adaptive step size factor, is the average value of the position, represents the initial value of the radius, represents the angle gain, is the random gain, is the control gain, represents the instantaneous radius value during search, represents the instantaneous angle during search, is the position of the global optimal individual of the current population, represents the position of the individual at time

[0034] This application embodiment also provides a tunnel seismic advanced detection waveform inversion system, including: A forward modeling module, which is used to perform forward modeling simulation according to the initial velocity model and the seismic source to obtain the simulated seismic record of the simulated seismic wave; A target function establishment module, which is used to establish a target function with the minimum difference between the observed seismic data and the simulated seismic record; An inversion module, which is used to perform inversion optimization on the target function by using the enhanced whale optimization algorithm, and repeatedly iterates to update the initial velocity model until the inversion reaches the minimum value of the target function; The inversion optimization of the objective function using the enhanced whale optimization algorithm includes: Initialize the parameters of the enhanced whale optimization algorithm, where the parameters include the initial value of the number of iterations, the maximum number of iterations, and the frequency factor; Generate an initial population. Use a piecewise Chebyshev - logistic map to generate the initial population, and add the reverse solution of the current individual to the initial population; Adjust the non - linear convergence factor; Update the first coefficient vector , the second system vector , the first random number and the second random number ; According to the values of the first coefficient vector and the second random number , trigger operations including a shrinking encircling mechanism, a random search mechanism, a bubble net attack mechanism, and a frequency fluctuation mechanism; Judge whether the maximum number of iterations or a preset termination condition is reached. If the termination condition is satisfied, output the optimal individual position and the fitness value of the optimal individual.

[0035] This application uses a simulation experimental environment to verify the effects of the embodiments of this application, including: Construction of a tunnel model. Refer to the tunnel model parameter diagram in Figure 2 . For acoustic wave propagation research using the tunnel model, the specific tunnel model parameters are as follows. The total length of the tunnel model is 150 meters, where the length of the tunnel section is 50 meters and the width is 15 meters. The upper and lower walls of the tunnel section are 40 meters and 45 meters away from the edges of the tunnel model respectively. To achieve high - quality acoustic wave signal acquisition, a total of 20 receivers are arranged in the tunnel model, including 8 on each of the upper and lower walls of the tunnel (both using an equally - spaced arrangement) and 4 at the tunnel face. The seismic source contains 3 excitation points, located at the uppermost, middle, and lowermost positions of the tunnel face respectively. In the setting of the acoustic wave propagation speed, the speed inside the tunnel is 340 m / s and the background speed is 3000 m / s.

[0036] Setting of forward - modeling parameters. Use a Ricker wavelet (Mexican hat wave) as the excitation source of the seismic wave, and its center frequency is set to 10 Hz for generating observed seismic data. The numerical simulation part is implemented based on the regular grid finite - difference method, and the grid size is 100×150. To ensure the simulation accuracy, the finite - difference operator uses a second - order accuracy in the time dimension and a fourth - order accuracy in the space dimension.

[0037] Setting of inversion parameters: The population size of the enhanced whale optimization algorithm is set to 300, and the fluctuation frequency is set to 5. The number of times to run the unscented Kalman filter in each iteration is set to 10. The objective function: , Among them, are model parameters, is the number of channels of the observed seismic data, and are respectively the th channel of the observed seismic data and the simulated seismic data.

[0038] Select the UHSA algorithm (Ultra-High-Speed Acoustic Full-Waveform Inversion method) as the comparative algorithm. The UHSA algorithm is a meta-heuristic algorithm specifically designed for the global inversion problem of tunnel space, with high efficiency and robustness, and can effectively handle complex geological inversion tasks. Through its application in tunnel scenarios, the UHSA algorithm has demonstrated its significant advantages in global optimization and parameter inversion, making it an ideal choice for comparison with the method of this application.

[0039] To verify the excellent performance of the method of this application embodiment in tunnel space waveform inversion, this application conducts a qualitative analysis of the synergistic effect of each enhancement mechanism based on the Bent Cigar function and the Levy function.

[0040] The Bent Cigar function presents a deep bowl-shaped structure, with relatively small function values near the global minimum point. As it moves away from the minimum point, the function values increase rapidly. The Levy function presents a complex multi-peak structure, with multiple local minimum points and one global minimum point, and the function surface has large fluctuations, making it difficult to optimize.

[0041] See Figure 3 The convergence curves of the Whale Optimization Algorithm and its variants of the enhanced Whale Optimization Algorithm of this application embodiment in the Bent Cigar function after 500 iterations and see Figure 4 The convergence curves of the Whale Optimization Algorithm and its variants of the enhanced Whale Optimization Algorithm of this application embodiment in the Levy function after 500 iterations. In the Bent Cigar function and the Levy function, conduct 500-iteration tests on each improved variant of the enhanced Whale Optimization Algorithm. Among them, CLM & RLI is to introduce the Whale Optimization Algorithm and the reverse learning initialization strategy; NCF is to add a non-linear convergence factor; FFM is to add a frequency fluctuation mechanism; LSS is to adopt a low-altitude hovering strategy; LSS & FFM is a variant of the Whale Optimization Algorithm that combines the low-altitude hovering strategy and the frequency fluctuation mechanism. The Whale Optimization Algorithm refers to WOA, and the variants of the Whale Optimization Algorithm refer to the method of this application. It can be seen that the fitness value of the method of this application changes significantly with the number of iterations.

[0042] For the inversion of concrete residues, see Figure 5The concrete residue tunnel model diagram in [reference], the real model is set as a model with a high-speed anomaly and a low-speed anomaly, which are 25 meters and 45 meters away from the tunnel face respectively. The speed of the high-speed anomaly is 5000 m / s, and the speed of the low-speed anomaly is 1000 m / s. In the figure, refers to the speed, and the speed here and in the following content all refer to the longitudinal wave speed. By comparing the inversion results of the method of this application and the UHSA algorithm, see Figure 6 (a) The concrete residue inversion result diagram of the method of this application in [reference] and Figure 6 (b) The concrete residue inversion result diagram of the UHSA algorithm in [reference], it can be seen that the method of this application locates the positions of the two speed anomalies very accurately, while the UHSA algorithm fails to correctly identify the positions of the anomalies. Further analyze the convergence curves of the method of this application and the UHSA algorithm in 100 iterations, see Figure 7 The convergence curves of the method of this application and the UHSA algorithm for concrete residue iteration 100 times provided in the embodiment of this application. It can be found that in the first 10 iterations, the fitness value of the method of this application rapidly drops from about 4.5 to nearly 3.5, showing the characteristic of rapid convergence; in contrast, the fitness value of the UHSA algorithm only slightly drops from about 5.5 to about 5.4 in the initial stage (the first 10 iterations), the change is very small, and then it tends to be stable, and the convergence speed is significantly slower.

[0043] See Figure 8 (a) The transverse velocity profile diagram at the longitudinal position of 20 m of the concrete residue of the method of this application and the UHSA algorithm in [reference], Figure 8 (b) The transverse velocity profile diagram at the longitudinal position of 60 m of the concrete residue of the method of this application and the UHSA algorithm in [reference]. By intercepting and comparing the transverse velocity profiles at the longitudinal positions of 20 m and 60 m, the speed and positioning of the method of this application are highly consistent with the actual speed of the real model, especially the low-speed anomaly is accurately located within the actual speed area. However, it can be seen from the UHSA speed that there is a large gap with the actual speed, and the existence of the low-speed anomaly cannot be identified. At the same time, the inversion result of the high-speed anomaly by the method of this application also shows high accuracy and better reflects the characteristics of the real model, while the UHSA algorithm still fails to accurately locate the high-speed anomaly. This result further verifies the superior performance of the method of this application in the inversion of complex geological models.

[0044] Inversion of inclined layer changes. See Figure 9 The inclined layer change model diagram in [reference], the real model is set as a model with a high-speed inclined layer and a low-speed inclined layer, which are 30 meters and 58 meters away from the tunnel face respectively. The speed of the high-speed anomaly is 5000 m / s, and the speed of the low-speed anomaly is 1000 m / s. From the inversion results of the method of this application and the UHSA algorithm, it can be known that Figure 10In (a), the inversion result of the inclined layer change of the method of the present application provided by the embodiment of the present application and Figure 10 in (b) is the inversion result diagram of the inclined layer change of the UHSA algorithm. It can be seen that the fitness value of the method of the present application decreases significantly with the number of iterations. Figure 11 The convergence curves of the method of the present application and the UHSA algorithm for the inversion of the inclined layer change iterated 100 times are shown. The speed of the method of the present application is highly consistent with the actual speed, and only a slight deviation appears at the boundary of the speed body. However, the UHS speed of the inversion result of the UHSA algorithm has a large gap with the actual speed, and the low-speed anomaly body cannot be identified. According to the comparison result of the convergence curves of the method of the present application and the UHSA algorithm in 100 iterations, the method of the present application drops rapidly in the first 10 iterations and then continuously converges to the optimal solution. In contrast, the UHSA algorithm falls into a local optimum from the beginning and fails to further converge.

[0045] See Figure 12 In (a) shown in, the transverse velocity profile at the longitudinal position of 32 m for the inversion of the inclined layer change of the method of the present application and the UHSA algorithm provided by the embodiment of the present application and Figure 12 in (b) shown in, the transverse velocity profile at 70 m. By intercepting the transverse velocity profiles at the longitudinal positions of 32 m and 70 m, the inversion result of the method of the present application fully reflects the velocity distribution and velocity values of the true model. In the intervals of 30 m to 100 m and 25 m to 100 m of the velocity profile, the velocity curves of the method of the present application basically fit the true velocity model. In contrast, there are large gaps between the velocity body inverted by the UHSA algorithm and the true model in terms of shape and velocity values.

[0046] In the embodiment of the present application, in population initialization, a piecewise Chebyshev-Logistic hybrid mapping is used to generate chaotic variables. By introducing the Chebyshev mapping in the low-value area to generate violent oscillations and using the Logistic mapping to maintain the uniformity of traversal in the high-value area, the advantages of the two mappings are combined, effectively avoiding the problem of uneven coverage caused by random initialization. To further improve the coverage range of the population, a reverse learning strategy is adopted, and the diversity of the population is increased by generating the reverse solution of the current individual. According to the qualitative analysis of each enhancement mechanism, it can be found that in the Bent Cigar function, the convergence value of the CLM & RLI curve is lower than that of the traditional WOA. It shows that the population initialization method adopted in the embodiment of the present application provides a proxy closer to the global optimal solution in the initial stage of iteration.

[0047] In the embodiment of the present application, non-linear convergence factor adjustment. A higher non-linear convergence factor value makes the probability increase significantly, so that the algorithm is more inclined to global exploration. As the algorithm iterates, in the middle and late stages The value drops rapidly, enhancing the local development ability of the algorithm. According to the qualitative analysis of each enhancement mechanism, it can be found that the algorithm equipped with a non-linear convergence factor (NCF) exhibits better global exploration performance in the early stage. This characteristic enables it to effectively break free from the limitation of local minima and avoid premature convergence. This performance advantage is fully verified in the tests of the BentCigar function and Levy function selected in this application.

[0048] In the embodiments of this application, the role of the frequency fluctuation mechanism is to enhance the ability of the whale optimization algorithm to jump out of local optimal solutions. In the traditional random search mechanism, an individual randomly selects any individual in the population and approaches it. This method lacks consideration of historical information and position distribution. Therefore, in multimodal functions, especially in cases like the Levy function where the number of local optimal solutions is large and the distribution is complex, it is easy to fall into local optima and difficult to find the global optimal solution. The frequency fluctuation mechanism proposed in the embodiments of this application improves this problem by adjusting the volatility during the search process. According to the qualitative analysis of each enhancement mechanism, the FFM curve shows a continuous downward trend throughout the iteration process. Compared with the standard whale optimization algorithm, this curve can find a solution better than the previous generation in each generation and always maintains a negative slope, thus effectively avoiding the problem of falling into local extreme points.

[0049] The method of the embodiments of this application combines the low-altitude hovering strategy of the red-tailed hawk with the spiral predation strategy of the whale optimization algorithm, thus achieving a stronger local development ability. According to the qualitative analysis of each enhancement mechanism, the convergence ability of the method of this application is particularly obvious in the convergence curve of the Bent Cigar function. The slope of the convergence curve is higher in the early stage, and the convergence speed is significantly accelerated. However, a stronger development ability is more likely to fall into local extreme points. In the Levy function, the LLS curve falls into the local optimal solution earlier than the original whale optimization algorithm. It should be noted that by combining the frequency fluctuation mechanism and the low-altitude hovering mechanism, the enhanced whale optimization algorithm can more effectively jump out of local extreme points while improving the convergence speed and accuracy. Compared with the original whale optimization algorithm, the whale optimization algorithm only containing the frequency fluctuation mechanism, and the whale optimization algorithm only containing the low-altitude hovering strategy, this combined strategy shows better convergence results (LSS & FFM curves) in the tests of the Bent Cigar function and Levy function.

[0050] Through in-depth analysis of the enhancement mechanism, the method of this application integrates the Chebyshev-Logistic hybrid mapping and the population initialization method of reverse learning, the non-linear convergence factor, the frequency fluctuation mechanism, and the low-altitude hovering strategy of the red-tailed hawk, thereby achieving optimization in terms of global search ability and local development ability. Specifically, the synergistic effect of the Chebyshev-Logistic hybrid mapping and reverse learning effectively improves the uniformity and diversity of the initial population distribution; the non-linear convergence factor can dynamically adjust the balance between exploration and development of the algorithm to adapt to the optimization requirements at different stages; the introduction of the frequency fluctuation mechanism enhances the ability of the algorithm to escape from local extrema; and the low-altitude hovering strategy of the red-tailed hawk further improves the accuracy and efficiency of local convergence. The experimental results show that the method of the embodiment of this application is superior to the UHSA algorithm in terms of inversion accuracy and convergence speed, fully verifying the effectiveness and superiority of the above enhancement mechanism.

[0051] The above are only the preferred embodiments of this application and are not intended to limit this application. Any modifications, equivalent replacements, and improvements made within the spirit and principle of this application shall be included within the protection scope of this application.

Claims

1. A waveform inversion method for tunnel seismic prediction ahead, characterized in that, Including: Performing forward modeling based on the initial velocity model and the seismic source to obtain the simulated seismic record of the simulated seismic wave; Establishing an objective function with the minimum difference between the observed seismic data and the simulated seismic record; Using the enhanced whale optimization algorithm to perform inversion optimization on the objective function, and repeatedly iterating to update the initial velocity model until the inversion reaches the minimum value of the objective function; The performing inversion optimization on the objective function using the enhanced whale optimization algorithm includes: Initializing the parameters of the enhanced whale optimization algorithm, where the parameters include the initial value of the iteration number, the maximum iteration number, and the frequency factor; Generating an initial population, generating the initial population using a piecewise Chebyshev-logistic map, and adding the reverse solution of the current individual to the initial population; Adjusting the non-linear convergence factor; Update the first coefficient vector , the second system vector , the first random number and the second random number ; According to the values of the first coefficient vector and the second random number , trigger operations including a contraction enclosure mechanism, a random search mechanism, a bubble net attack mechanism, and a frequency fluctuation mechanism; Judging whether the maximum iteration number or the preset termination condition is reached. If the termination condition is satisfied, output the position of the optimal individual and the fitness value of the optimal individual.

2. The tunnel seismic prediction waveform inversion method according to claim 1, characterized in that The generating an initial population, generating the initial population using a piecewise Chebyshev-logistic map, and adding the reverse solution of the current individual to the initial population includes: Generate the reverse solution of the current individual to increase population diversity, which is expressed as: , where represents the reverse solution of the current individual, represents the lower limit of the individual position variable, represents the upper limit of the individual position variable, represents the -th position variable of the -th individual, and the calculation method is: represents the -th chaotic variable, is updated as: , Perturbation factor , represents the population size, is a parameter for controlling the behavior of the logistic map.

3. The waveform inversion method for tunnel seismic prediction ahead according to claim 1, wherein Adjust the non - linear convergence factor, including: , represents the non - linear convergence factor, represents the maximum number of iterations, represents the number of iterations.

4. The waveform inversion method for tunnel seismic prediction ahead according to claim 1, characterized in that The contraction and encirclement mechanism includes: when it satisfies , and , the globally optimal individual represents the prey space coordinates, and other individuals gradually approach the globally optimal individual by adjusting their own positions. The population update formula is: ， ， Among them, represents the weighted distance between an individual and its prey, represents the position of the individual at time, is the position of the globally optimal individual in the current population, represents the position of the individual at time; The first coefficient vector and the second coefficient vector are calculated by the following formula: ， ， Among them, is the third random number between [0, 1], is the fourth random number between [0, 1], is the non-linear convergence factor.

5. The waveform inversion method for tunnel seismic advanced detection according to claim 1, characterized in that The random search mechanism includes: when satisfying , , in the iterative process, each individual randomly selects any feasible solution in the current population as the target position with a probability of , and expands the search range by moving towards the target position. The population update formula is: ， ， Among them, is the target position randomly selected from the current population, represents the weighted distance between the individual and the prey, represents the position of the individual at time, represents the position of the individual at time.

6. The waveform inversion method for tunnel seismic prediction ahead according to claim 1, characterized in that The bubble net attack mechanism includes: when satisfying , , the population update formula is: ， ， Among them, is set to 1, represents the Euclidean distance between the individual and the prey, represents the individual at time position, is the position of the global optimal individual of the current population, represents the individual at time position, is the modulo operation on the iteration number and the frequency factor performs the modulo operation.

7. The waveform inversion method for tunnel seismic advanced detection according to claim 1, characterized in that The frequency fluctuation mechanism includes: when satisfying , , the update of the population: ， is the amplitude factor, defined as: ， Among them is the angular frequency, calculated as ; is a random phase angle within the interval, is the oscillation factor, such that the amplitude factor varies periodically with time. By adjusting the angular frequency and the phase angle, the frequency and phase of the oscillation are controlled to guide the search agent to explore in the search space, such that the amplitude of the oscillation gradually decreases as the number of iterations increases, is the modulo operation on the number of iterations and the frequency factor , represents the maximum number of iterations.

8. The waveform inversion method for tunnel seismic prediction ahead according to claim 1, characterized in that, It also includes using the low-altitude hovering mechanism of the red-tailed hawk to enhance its local development ability, and the specific formula is expressed as: ， ， ， ， Replacement after normalization and : ， Among them, and represent the direction coordinates at a moment, is expressed as an adaptive step size factor, is the average value of the position, represents the initial value of the radius, represents the angle gain, is the random gain, is the control gain, represents the instantaneous radius value during the search, represents the instantaneous angle during the search, is the position of the globally optimal individual of the current population, represents the position of the individual at a moment.

9. A waveform inversion system for tunnel seismic advanced detection, characterized in that, Including: A forward modeling module for performing forward modeling based on the initial velocity model and the seismic source to obtain the simulated seismic record of the simulated seismic wave; An objective function establishment module for establishing an objective function with the minimum difference between the observed seismic data and the simulated seismic record; An inversion module for performing inversion optimization on the objective function using the enhanced whale optimization algorithm, and repeatedly iterating to update the initial velocity model until the inversion reaches the minimum value of the objective function; The performing inversion optimization on the objective function using the enhanced whale optimization algorithm includes: Initializing the parameters of the enhanced whale optimization algorithm, where the parameters include the initial value of the iteration number, the maximum iteration number, and the frequency factor; Generating an initial population, generating the initial population using a piecewise Chebyshev-logistic map, and adding the reverse solution of the current individual to the initial population; Adjusting the non-linear convergence factor; Update the first coefficient vector , the second system vector , the first random number and the second random number ; According to the values of the first coefficient vector and the second random number , trigger operations including a contraction enclosure mechanism, a random search mechanism, a bubble net attack mechanism, and a frequency fluctuation mechanism; Judging whether the maximum iteration number or the preset termination condition is reached. If the termination condition is satisfied, output the position of the optimal individual and the fitness value of the optimal individual.

Images